From 59cb3f9e98ca3dae2c1f494d31f13a61b60c517a Mon Sep 17 00:00:00 2001 From: pilar Date: Thu, 4 Jun 2020 12:50:10 -0400 Subject: [PATCH] Geometrical methods for dividing buildings by plane --- city_model_structure/plying_with_geometry.py | 81 ++++++++++++++++++++ 1 file changed, 81 insertions(+) create mode 100644 city_model_structure/plying_with_geometry.py diff --git a/city_model_structure/plying_with_geometry.py b/city_model_structure/plying_with_geometry.py new file mode 100644 index 00000000..bd8e93ca --- /dev/null +++ b/city_model_structure/plying_with_geometry.py @@ -0,0 +1,81 @@ +import trimesh +from trimesh import intersections as inter +import open3d as o3d +import numpy as np + + +def extract_points(segment_list): + point_list = np.asarray(segment_list[0]) + for segment in segment_list: + found = False + for new_point in segment: + for point in point_list: + comparison = new_point == point + if comparison.all(): + found = True + break + if not found: + point_list = np.concatenate((point_list, [new_point])) + return point_list + + +def point_cloud_to_mesh(point_list, normal_list): + pcd = o3d.geometry.PointCloud() + pcd.points = o3d.utility.Vector3dVector(point_list) + pcd.normals = o3d.utility.Vector3dVector(normal_list) + distances = pcd.compute_nearest_neighbor_distance() + avg_dist = np.mean(distances) + radius = 3 * avg_dist + bpa_mesh = o3d.geometry.TriangleMesh.create_from_point_cloud_ball_pivoting( + pcd, + o3d.utility.DoubleVector([radius, radius * 2])) + mesh_result = trimesh.Trimesh(vertices=np.asarray(bpa_mesh.vertices), faces=np.asarray(bpa_mesh.triangles)) + return mesh_result + + +def merge_meshes(mesh1, mesh2): + v_1 = mesh1.vertices + f_1 = mesh1.faces + v_2 = mesh2.vertices + f_2 = mesh2.faces + length = len(v_1) + v_merge = np.concatenate((v_1, v_2)) + f_merge = np.asarray(f_1) + + for item in f_2: + point1 = item.item(0) + length + point2 = item.item(1) + length + point3 = item.item(2) + length + surface = np.asarray([point1, point2, point3]) + f_merge = np.concatenate((f_merge, [surface])) + + mesh_merge = trimesh.Trimesh(vertices=v_merge, faces=f_merge) + + return mesh_merge + + +# example +vertices = [[0, 0, 0], [0, 0, 1], [0, 1, 0], [0, 1, 1], [1, 0, 0], [1, 0, 1], [1, 1, 0], [1, 1, 1]] +faces = [[0, 1, 2], [1, 3, 2], [1, 5, 3], [5, 3, 7], [0, 4, 1], [4, 5, 1], [4, 6, 5], [5, 6, 7], [0, 2, 4], [2, 6, 4], + [3, 7, 6], [2, 3, 6]] +mesh = trimesh.Trimesh(vertices=vertices, faces=faces) +normal_plane = [0, 0, 1] +point_plane = [0, 0, 0.5] +normal_plane_opp = [0, 0, -1] + +# Step 1: division and extraction of segments +mesh_1 = inter.slice_mesh_plane(mesh, normal_plane, point_plane) +mesh_1_segments = inter.mesh_plane(mesh, normal_plane, point_plane) + +# Step 2: create mesh from point cloud +points = extract_points(mesh_1_segments) +normals = np.empty((len(points), 3)) +for i in range(0, len(normals)): + normals[i] = normal_plane_opp +mesh_2 = point_cloud_to_mesh(points, normals) + +# Step 3: merge both meshes +mesh_final = merge_meshes(mesh_1, mesh_2) + + +print(mesh_final.vertices, mesh_final.faces)