Voxelization 体素化
点云和三角形网格是非常灵活但不规则的几何类型。体素网格（voxel grid）是在 3D 网格上定义的 3D 中的另一种几何类型，体素可以被视为 2D 中像素（pixel）的 3D 对应物。Open3D 具有可用于处理体素网格的几何类型VoxelGrid。
From triangle mesh 从三角形网格
Open3D 提供了从三角形网格创建体素网格的方法create_from_triangle_mesh。它返回一个体素网格，其中所有与三角形相交的体素都设置为1 ，所有其他体素都设置为 0。参数voxel_size定义体素网格的分辨率。
print('input') bunny = o3d.data.BunnyMesh() mesh = o3d.io.read_triangle_mesh(bunny.path) # fit to unit cube mesh.scale(1 / np.max(mesh.get_max_bound() - mesh.get_min_bound()), center=mesh.get_center()) o3d.visualization.draw_geometries([mesh]) print('voxelization') voxel_grid = o3d.geometry.VoxelGrid.create_from_triangle_mesh(mesh, voxel_size=0.05) o3d.visualization.draw_geometries([voxel_grid])
From point cloud 从点云
体素网格也可以使用create_from_point_cloud方法从点云创建。如果点云的至少一个点在体素内，则体素被占用。体素的颜色是体素内所有点的平均值。参数voxel_size定义体素网格的分辨率。
print('input') armadillo = o3d.data.ArmadilloMesh() mesh = o3d.io.read_triangle_mesh(armadillo.path) N = 2000 pcd = mesh.sample_points_poisson_disk(N) # fit to unit cube pcd.scale(1 / np.max(pcd.get_max_bound() - pcd.get_min_bound()), center=pcd.get_center()) pcd.colors = o3d.utility.Vector3dVector(np.random.uniform(0, 1, size=(N, 3))) o3d.visualization.draw_geometries([pcd]) print('voxelization') voxel_grid = o3d.geometry.VoxelGrid.create_from_point_cloud(pcd, voxel_size=0.05) o3d.visualization.draw_geometries([voxel_grid])
Inclusion test 包含测试
体素网格还可用于测试点是否在占用的体素内。check_if_included方法将(n,3)数组作为输入并输出bool数组。
queries = np.asarray(pcd.points) output = voxel_grid.check_if_included(o3d.utility.Vector3dVector(queries)) print(output[:10])
Voxel carving 体素雕刻
方法create_from_point_cloud和create_from_triangle_mesh创建占用的体素仅在几何体的表面上。但是，可以从许多深度图（depth maps）或剪影(sihouettes)中雕刻体素网格。Open3D提供了体素雕刻的carve_depth_map方法和carve_silhouette方法。
下面的代码演示了首先从几何图形渲染深度图并使用这些深度图雕刻密集体素网格的用法。结果是给定形状的填充体素网格。
import open3d as o3d import numpy as np def xyz_spherical(xyz): x = xyz[0] y = xyz[1] z = xyz[2] r = np.sqrt(x * x + y * y + z * z) r_x = np.arccos(y / r) r_y = np.arctan2(z, x) return [r, r_x, r_y] def get_rotation_matrix(r_x, r_y): rot_x = np.asarray([[1, 0, 0], [0, np.cos(r_x), -np.sin(r_x)], [0, np.sin(r_x), np.cos(r_x)]]) rot_y = np.asarray([[np.cos(r_y), 0, np.sin(r_y)], [0, 1, 0], [-np.sin(r_y), 0, np.cos(r_y)]]) return rot_y.dot(rot_x) def get_extrinsic(xyz): rvec = xyz_spherical(xyz) r = get_rotation_matrix(rvec[1], rvec[2]) t = np.asarray([0, 0, 2]).transpose() trans = np.eye(4) trans[:3, :3] = r trans[:3, 3] = t return trans def preprocess(model): min_bound = model.get_min_bound() max_bound = model.get_max_bound() center = min_bound + (max_bound - min_bound) / 2.0 scale = np.linalg.norm(max_bound - min_bound) / 2.0 vertices = np.asarray(model.vertices) vertices -= center model.vertices = o3d.utility.Vector3dVector(vertices / scale) return model def voxel_carving(mesh, cubic_size, voxel_resolution, w=300, h=300, use_depth=True, surface_method='pointcloud'): mesh.compute_vertex_normals() camera_sphere = o3d.geometry.TriangleMesh.create_sphere() # setup dense voxel grid voxel_carving = o3d.geometry.VoxelGrid.create_dense( width=cubic_size, height=cubic_size, depth=cubic_size, voxel_size=cubic_size / voxel_resolution, origin=[-cubic_size / 2.0, -cubic_size / 2.0, -cubic_size / 2.0], color=[1.0, 0.7, 0.0]) # rescale geometry camera_sphere = preprocess(camera_sphere) mesh = preprocess(mesh) # setup visualizer to render depthmaps vis = o3d.visualization.Visualizer() vis.create_window(width=w, height=h, visible=False) vis.add_geometry(mesh) vis.get_render_option().mesh_show_back_face = True ctr = vis.get_view_control() param = ctr.convert_to_pinhole_camera_parameters() # carve voxel grid pcd_agg = o3d.geometry.PointCloud() centers_pts = np.zeros((len(camera_sphere.vertices), 3)) for cid, xyz in enumerate(camera_sphere.vertices): # get new camera pose trans = get_extrinsic(xyz) param.extrinsic = trans c = np.linalg.inv(trans).dot(np.asarray([0, 0, 0, 1]).transpose()) centers_pts[cid, :] = c[:3] ctr.convert_from_pinhole_camera_parameters(param) # capture depth image and make a point cloud vis.poll_events() vis.update_renderer() depth = vis.capture_depth_float_buffer(False) pcd_agg += o3d.geometry.PointCloud.create_from_depth_image( o3d.geometry.Image(depth), param.intrinsic, param.extrinsic, depth_scale=1) # depth map carving method if use_depth: voxel_carving.carve_depth_map(o3d.geometry.Image(depth), param) else: voxel_carving.carve_silhouette(o3d.geometry.Image(depth), param) print("Carve view %03d/%03d" % (cid + 1, len(camera_sphere.vertices))) vis.destroy_window() # add voxel grid survace print('Surface voxel grid from %s' % surface_method) if surface_method == 'pointcloud': voxel_surface = o3d.geometry.VoxelGrid.create_from_point_cloud_within_bounds( pcd_agg, voxel_size=cubic_size / voxel_resolution, min_bound=(-cubic_size / 2, -cubic_size / 2, -cubic_size / 2), max_bound=(cubic_size / 2, cubic_size / 2, cubic_size / 2)) elif surface_method == 'mesh': voxel_surface = o3d.geometry.VoxelGrid.create_from_triangle_mesh_within_bounds( mesh, voxel_size=cubic_size / voxel_resolution, min_bound=(-cubic_size / 2, -cubic_size / 2, -cubic_size / 2), max_bound=(cubic_size / 2, cubic_size / 2, cubic_size / 2)) else: raise Exception('invalid surface method') voxel_carving_surface = voxel_surface + voxel_carving return voxel_carving_surface, voxel_carving, voxel_surface armadillo_path = r'../data/ArmadilloMesh.ply' mesh = o3d.io.read_triangle_mesh(armadillo_path) visualization = True cubic_size = 2.0 voxel_resolution = 128.0 voxel_grid, voxel_carving, voxel_surface = voxel_carving( mesh, cubic_size, voxel_resolution) print("surface voxels") print(voxel_surface) o3d.visualization.draw_geometries([voxel_surface]) print("carved voxels") print(voxel_carving) o3d.visualization.draw_geometries([voxel_carving]) print("combined voxels (carved + surface)") print(voxel_grid) o3d.visualization.draw_geometries([voxel_grid])
