city_retrofit/city_model_structure/attributes/polyhedron.py

311 lines
9.8 KiB
Python
Raw Blame History

This file contains ambiguous Unicode characters

This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

"""
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=0.5, area_delta=0.001)
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, True):
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):
# print('wroooong', complementary_surface.type)
# print(complementary_surface.points)
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()