310 lines
9.7 KiB
Python
310 lines
9.7 KiB
Python
"""
|
||
Polyhedron module
|
||
SPDX - License - Identifier: LGPL - 3.0 - or -later
|
||
Copyright © 2020 Project Author Guille Gutierrez guillermo.gutierrezmorote@concordia.ca
|
||
"""
|
||
import numpy as np
|
||
from trimesh import Trimesh
|
||
from helpers.geometry_helper import GeometryHelper
|
||
from helpers.configuration_helper import ConfigurationHelper
|
||
from city_model_structure.attributes.surface import Surface
|
||
|
||
|
||
class Polyhedron:
|
||
"""
|
||
Polyhedron class
|
||
"""
|
||
|
||
def __init__(self, surfaces):
|
||
self._surfaces = list(surfaces)
|
||
self._polygons = [s.polygon for s in surfaces]
|
||
self._polyhedron = None
|
||
self._volume = None
|
||
self._faces = None
|
||
self._vertices = None
|
||
self._mesh = None
|
||
self._centroid = None
|
||
self._max_z = None
|
||
self._max_y = None
|
||
self._max_x = None
|
||
self._min_z = None
|
||
self._min_y = None
|
||
self._min_x = None
|
||
self._geometry = GeometryHelper(delta=5.0, area_delta=0.01)
|
||
print(f'[{self._geometry._area_delta} {self._geometry._delta}]')
|
||
|
||
def _position_of(self, point, face):
|
||
vertices = self.vertices
|
||
for i in range(len(vertices)):
|
||
# ensure not duplicated vertex
|
||
if i not in face and self._geometry.almost_equal(vertices[i], point):
|
||
return i
|
||
return -1
|
||
|
||
@property
|
||
def vertices(self) -> np.ndarray:
|
||
"""
|
||
Polyhedron vertices
|
||
:return: np.ndarray(int)
|
||
"""
|
||
if self._vertices is None:
|
||
vertices, self._vertices = [], []
|
||
_ = [vertices.extend(s.points) for s in self._surfaces]
|
||
for vertex_1 in vertices:
|
||
found = False
|
||
for vertex_2 in self._vertices:
|
||
found = False
|
||
if self._geometry.almost_equal(vertex_1, vertex_2):
|
||
found = True
|
||
break
|
||
if not found:
|
||
self._vertices.append(vertex_1)
|
||
self._vertices = np.asarray(self._vertices)
|
||
return self._vertices
|
||
|
||
@staticmethod
|
||
def _point(coordinates):
|
||
return coordinates[0], coordinates[1], coordinates[2]
|
||
|
||
@staticmethod
|
||
def _get_regions(point_index, points_list):
|
||
if point_index == 0:
|
||
# first point in the polygon so the triangle is the points n-1, n, 0
|
||
triangle_left = ' '.join(str(e) for e in [*points_list[len(points_list) - 6:], *points_list[0:3]])
|
||
# remove point n
|
||
rest_points_left = ' '.join(str(e) for e in [*points_list[:len(points_list) - 3]])
|
||
elif point_index == 3:
|
||
# second point in the polygon so the triangle is the points n, 0, 1
|
||
triangle_left = ' '.join(str(e) for e in [*points_list[len(points_list) - 3:], *points_list[0:6]])
|
||
# remove point 0
|
||
rest_points_left = ' '.join(str(e) for e in [*points_list[3:]])
|
||
else:
|
||
# normal point index-2¸index-1, index
|
||
triangle_left = ' '.join(str(e) for e in [*points_list[point_index - 6:point_index + 3]])
|
||
# remove middle point (index - 1)
|
||
rest_points_left = ' '.join(str(e) for e in [*points_list[0:point_index - 3], *points_list[point_index:]])
|
||
if point_index < len(points_list) - 6:
|
||
# normal point index, index+1, index+2
|
||
triangle_right = ' '.join(str(e) for e in [*points_list[point_index:point_index + 9]])
|
||
rest_points_right = ' '.join(str(e) for e in [*points_list[0:point_index + 3], *points_list[point_index + 6:]])
|
||
elif point_index == (len(points_list) - 6):
|
||
# last two points in the polygon so the triangle is the points n-1, n, 0
|
||
triangle_right = ' '.join(str(e) for e in [*points_list[point_index:], *points_list[0:3]])
|
||
rest_points_right = ' '.join(str(e) for e in [*points_list[:len(points_list - 3)]])
|
||
else:
|
||
# last point in the polygon so the triangle is n, 0, 1
|
||
triangle_right = ' '.join(str(e) for e in [*points_list[len(points_list) - 3:], *points_list[0:6]])
|
||
rest_points_right = ' '.join(str(e) for e in [*points_list[3:]])
|
||
return (Surface(triangle_left, remove_last=False), Surface(rest_points_left, remove_last=False)), \
|
||
(Surface(triangle_right, remove_last=False), Surface(rest_points_right, remove_last=False))
|
||
|
||
def _triangulate(self, surface):
|
||
triangles = []
|
||
complementary_surface = surface
|
||
triangles_count = len(surface.points) - 2
|
||
points_list = surface.points_list
|
||
point_index = 0
|
||
area = surface.area
|
||
while len(triangles) < triangles_count:
|
||
# get triangles and regions in both direction to find ears
|
||
left_direction, right_direction = Polyhedron._get_regions(point_index, points_list)
|
||
# todo: use enum to describe triangle or rest instead 0 1
|
||
right_area = right_direction[0].area + right_direction[1].area
|
||
left_area = left_direction[0].area + left_direction[1].area
|
||
if self._geometry.almost_same_area(area, left_area) and self._geometry.almost_same_area(area, right_area):
|
||
# Both seems to be an ear, choose the more precise
|
||
if np.abs(left_area-area) < np.abs(right_area-area):
|
||
area = left_direction[1].area
|
||
point_index = 0
|
||
triangles.append(left_direction[0])
|
||
points_list = left_direction[1].points_list
|
||
complementary_surface = left_direction[1]
|
||
else:
|
||
area = right_direction[1].area
|
||
point_index = 0
|
||
triangles.append(right_direction[0])
|
||
points_list = right_direction[1].points_list
|
||
complementary_surface = right_direction[1]
|
||
elif self._geometry.almost_same_area(area, left_area):
|
||
area = left_direction[1].area
|
||
point_index = 0
|
||
triangles.append(left_direction[0])
|
||
points_list = left_direction[1].points_list
|
||
complementary_surface = left_direction[1]
|
||
elif self._geometry.almost_same_area(area, right_area):
|
||
area = right_direction[1].area
|
||
point_index = 0
|
||
triangles.append(right_direction[0])
|
||
points_list = right_direction[1].points_list
|
||
complementary_surface = right_direction[1]
|
||
else:
|
||
point_index = point_index + 3
|
||
if point_index >= len(points_list):
|
||
return triangles
|
||
if len(points_list) == 9:
|
||
# the rest point's are already a triangle
|
||
triangles.append(complementary_surface)
|
||
|
||
return triangles
|
||
|
||
@property
|
||
def faces(self) -> [[int]]:
|
||
|
||
"""
|
||
Polyhedron faces
|
||
:return: [[int]]
|
||
"""
|
||
|
||
if self._faces is None:
|
||
self._faces = []
|
||
|
||
for surface in self._surfaces:
|
||
|
||
face = []
|
||
points = surface.points
|
||
if len(points) != 3:
|
||
sub_surfaces = self._triangulate(surface)
|
||
for sub_surface in sub_surfaces:
|
||
face = []
|
||
points = sub_surface.points
|
||
for point in points:
|
||
face.append(self._position_of(point, face))
|
||
self._faces.append(face)
|
||
else:
|
||
for point in points:
|
||
face.append(self._position_of(point, face))
|
||
self._faces.append(face)
|
||
return self._faces
|
||
|
||
@property
|
||
def _polyhedron_mesh(self):
|
||
if self._mesh is None:
|
||
self._mesh = Trimesh(vertices=self.vertices, faces=self.faces)
|
||
return self._mesh
|
||
|
||
@property
|
||
def volume(self):
|
||
"""
|
||
Polyhedron volume in cubic meters
|
||
:return: float
|
||
"""
|
||
if self._volume is None:
|
||
if not self._polyhedron_mesh.is_volume:
|
||
self._volume = np.inf
|
||
else:
|
||
self._volume = self._polyhedron_mesh.volume
|
||
return self._volume
|
||
|
||
@property
|
||
def max_z(self):
|
||
"""
|
||
Polyhedron maximal z value
|
||
:return: float
|
||
"""
|
||
if self._max_z is None:
|
||
self._max_z = ConfigurationHelper().min_coordinate
|
||
for surface in self._surfaces:
|
||
for point in surface.points:
|
||
if self._max_z < point[2]:
|
||
self._max_z = point[2]
|
||
return self._max_z
|
||
|
||
@property
|
||
def max_y(self):
|
||
"""
|
||
Polyhedron maximal y value
|
||
:return: float
|
||
"""
|
||
if self._max_y is None:
|
||
self._max_y = ConfigurationHelper().min_coordinate
|
||
for surface in self._surfaces:
|
||
for point in surface.points:
|
||
if self._max_y < point[1]:
|
||
self._max_y = point[1]
|
||
return self._max_y
|
||
|
||
@property
|
||
def max_x(self):
|
||
"""
|
||
Polyhedron maximal x value
|
||
:return: float
|
||
"""
|
||
if self._max_x is None:
|
||
self._max_x = ConfigurationHelper().min_coordinate
|
||
for surface in self._surfaces:
|
||
for point in surface.points:
|
||
self._max_x = max(self._max_x, point[0])
|
||
return self._max_x
|
||
|
||
@property
|
||
def min_z(self):
|
||
"""
|
||
Polyhedron minimal z value
|
||
:return: float
|
||
"""
|
||
if self._min_z is None:
|
||
self._min_z = self.max_z
|
||
for surface in self._surfaces:
|
||
for point in surface.points:
|
||
if self._min_z > point[2]:
|
||
self._min_z = point[2]
|
||
return self._min_z
|
||
|
||
@property
|
||
def min_y(self):
|
||
"""
|
||
Polyhedron minimal y value
|
||
:return: float
|
||
"""
|
||
if self._min_y is None:
|
||
self._min_y = self.max_y
|
||
for surface in self._surfaces:
|
||
for point in surface.points:
|
||
if self._min_y > point[1]:
|
||
self._min_y = point[1]
|
||
return self._min_y
|
||
|
||
@property
|
||
def min_x(self):
|
||
"""
|
||
Polyhedron minimal x value
|
||
:return: float
|
||
"""
|
||
if self._min_x is None:
|
||
self._min_x = self.max_x
|
||
for surface in self._surfaces:
|
||
for point in surface.points:
|
||
if self._min_x > point[0]:
|
||
self._min_x = point[0]
|
||
return self._min_x
|
||
|
||
@property
|
||
def center(self):
|
||
"""
|
||
Polyhedron centroid
|
||
:return: [x,y,z]
|
||
"""
|
||
x = (self.max_x + self.min_x) / 2
|
||
y = (self.max_y + self.min_y) / 2
|
||
z = (self.max_z + self.min_z) / 2
|
||
return [x, y, z]
|
||
|
||
def stl_export(self, full_path):
|
||
"""
|
||
Export the polyhedron to stl given file
|
||
:param full_path: str
|
||
:return: None
|
||
"""
|
||
self._polyhedron_mesh.export(full_path, 'stl_ascii')
|
||
|
||
def obj_export(self, full_path):
|
||
"""
|
||
Export the polyhedron to obj given file
|
||
:param full_path: str
|
||
:return: None
|
||
"""
|
||
self._polyhedron_mesh.export(full_path, 'obj')
|
||
|
||
def show(self):
|
||
self._polyhedron_mesh.show()
|