1# ----------------------------------------------------------------------------
2# - Open3D: www.open3d.org -
3# ----------------------------------------------------------------------------
4# Copyright (c) 2018-2023 www.open3d.org
5# SPDX-License-Identifier: MIT
6# ----------------------------------------------------------------------------
7
8import open3d as o3d
9import numpy as np
10
11
12def xyz_spherical(xyz):
13 x = xyz[0]
14 y = xyz[1]
15 z = xyz[2]
16 r = np.sqrt(x * x + y * y + z * z)
17 r_x = np.arccos(y / r)
18 r_y = np.arctan2(z, x)
19 return [r, r_x, r_y]
20
21
22def get_rotation_matrix(r_x, r_y):
23 rot_x = np.asarray([[1, 0, 0], [0, np.cos(r_x), -np.sin(r_x)],
24 [0, np.sin(r_x), np.cos(r_x)]])
25 rot_y = np.asarray([[np.cos(r_y), 0, np.sin(r_y)], [0, 1, 0],
26 [-np.sin(r_y), 0, np.cos(r_y)]])
27 return rot_y.dot(rot_x)
28
29
30def get_extrinsic(xyz):
31 rvec = xyz_spherical(xyz)
32 r = get_rotation_matrix(rvec[1], rvec[2])
33 t = np.asarray([0, 0, 2]).transpose()
34 trans = np.eye(4)
35 trans[:3, :3] = r
36 trans[:3, 3] = t
37 return trans
38
39
40def preprocess(model):
41 min_bound = model.get_min_bound()
42 max_bound = model.get_max_bound()
43 center = min_bound + (max_bound - min_bound) / 2.0
44 scale = np.linalg.norm(max_bound - min_bound) / 2.0
45 vertices = np.asarray(model.vertices)
46 vertices -= center
47 model.vertices = o3d.utility.Vector3dVector(vertices / scale)
48 return model
49
50
51def voxel_carving(mesh, cubic_size, voxel_resolution, w=300, h=300):
52 mesh.compute_vertex_normals()
53 camera_sphere = o3d.geometry.TriangleMesh.create_sphere(radius=1.0,
54 resolution=10)
55
56 # Setup dense voxel grid.
57 voxel_carving = o3d.geometry.VoxelGrid.create_dense(
58 width=cubic_size,
59 height=cubic_size,
60 depth=cubic_size,
61 voxel_size=cubic_size / voxel_resolution,
62 origin=[-cubic_size / 2.0, -cubic_size / 2.0, -cubic_size / 2.0],
63 color=[1.0, 0.7, 0.0])
64
65 # Rescale geometry.
66 camera_sphere = preprocess(camera_sphere)
67 mesh = preprocess(mesh)
68
69 # Setup visualizer to render depthmaps.
70 vis = o3d.visualization.Visualizer()
71 vis.create_window(width=w, height=h, visible=False)
72 vis.add_geometry(mesh)
73 vis.get_render_option().mesh_show_back_face = True
74 ctr = vis.get_view_control()
75 param = ctr.convert_to_pinhole_camera_parameters()
76
77 # Carve voxel grid.
78 centers_pts = np.zeros((len(camera_sphere.vertices), 3))
79 for cid, xyz in enumerate(camera_sphere.vertices):
80 # Get new camera pose.
81 trans = get_extrinsic(xyz)
82 param.extrinsic = trans
83 c = np.linalg.inv(trans).dot(np.asarray([0, 0, 0, 1]).transpose())
84 centers_pts[cid, :] = c[:3]
85 ctr.convert_from_pinhole_camera_parameters(param)
86
87 # Capture depth image and make a point cloud.
88 vis.poll_events()
89 vis.update_renderer()
90 depth = vis.capture_depth_float_buffer(False)
91
92 # Depth map carving method.
93 voxel_carving.carve_depth_map(o3d.geometry.Image(depth), param)
94 print("Carve view %03d/%03d" % (cid + 1, len(camera_sphere.vertices)))
95 vis.destroy_window()
96
97 return voxel_carving
98
99
100if __name__ == "__main__":
101
102 armadillo_data = o3d.data.ArmadilloMesh()
103 mesh = o3d.io.read_triangle_mesh(armadillo_data.path)
104 cubic_size = 2.0
105 voxel_resolution = 128.0
106
107 carved_voxels = voxel_carving(mesh, cubic_size, voxel_resolution)
108 print("Carved voxels ...")
109 print(carved_voxels)
110 o3d.visualization.draw([carved_voxels])