目录
一、概述
1.1原理
1.2实现步骤
1.3应用场景
二、代码实现
2.1关键函数
2.1.1步骤
2.1.2函数代码
2.2完整代码
三、实现效果
3.1原始点云
3.2重建后点云
Open3D点云算法汇总及实战案例汇总的目录地址:
Open3D点云算法与点云深度学习案例汇总(长期更新)-CSDN博客
一、概述
行进立方体算法(Marching Cubes)是一种用于从三维标量场数据(如体素网格)生成等值面(isosurface)的算法。其基本原理是通过遍历三维标量场中的每个体素单元,根据每个体素的顶点值,构建出相应的三角网格。
Open3D目前并没有直接的行进立方体算法实现,需要结合其他库来完成,例如scikit-image中的marching_cubes方法。结合 Open3D 和 scikit-image编写行进立方体算法,可以实现从点云数据生成高质量的等值面
1.1原理
- 遍历三维标量场中的每个体素单元。
- 根据体素单元的八个顶点值,确定该单元内部是否包含等值面。
- 使用预定义的查找表,根据顶点值的组合生成对应的三角形面片。
- 将所有生成的三角形面片合并,构成完整的等值面。
1.2实现步骤
- 加载点云数据:读取点云数据文件。
- 估计法向量:计算点云中每个点的法向量。
- 生成体素网格:从点云数据生成体素网格。
- 使用行进立方体算法:从体素网格生成等值面。
- 可视化和保存结果:可视化生成的等值面,并将其保存到文件中。
1.3应用场景
- 医学成像:如CT、MRI数据的三维重建。
- 计算机图形学:生成多边形网格模型。
- 科学可视化:如气象数据、地质数据的三维可视化
二、代码实现
2.1关键函数
- 这段代码使用行进立方体算法(Marching Cubes)从体素网格生成三维等值面。具体来说,它将输入的体素网格转换为三维体积数据,并利用
- scikit-image 库的 marching_cubes 方法生成等值面。生成的等值面被转换为 Open3D 的三角网格对象 TriangleMesh,可以用于可视化和进一步处理。
2.1.1步骤
- 获取体素索引:从体素网格中提取所有体素的索引,并将其转换为 NumPy 数组。
- 获取边界和体素大小:获取体素网格的最小和最大边界,以及体素的大小。
- 计算体素网格形状:根据边界和体素大小计算体素网格的形状。
- 创建体素数组:创建一个空的体素数组,用于存储体积数据。
- 填充体素数组:将体素网格中的每个体素标记为1,生成体积数据。
- 行进立方体算法:使用 scikit-image 的 marching_cubes 方法生成等值面。
- 转换为三角网格:将生成的顶点、面片和法向量转换为 Open3D 的 TriangleMesh 对象。
- 返回生成的三角网格:返回生成的 TriangleMesh 对象,可以用于可视化和进一步处理。
2.1.2函数代码
def marching_cubes_from_voxel_grid(voxel_grid, level):
"""
使用行进立方体算法从体素网格生成等值面。
参数:
voxel_grid (open3d.geometry.VoxelGrid): 输入的体素网格。
level (float): 等值面的阈值。
返回:
open3d.geometry.TriangleMesh: 生成的三角网格。
"""
# 获取体素网格中的所有体素索引
voxel_indices = np.asarray([v.grid_index for v in voxel_grid.get_voxels()])
# 获取体素网格的最小和最大边界
min_bound = voxel_grid.get_min_bound()
max_bound = voxel_grid.get_max_bound()
# 获取体素大小
voxel_size = voxel_grid.voxel_size
# 计算体素网格的形状
volume_shape = ((max_bound - min_bound) / voxel_size).astype(int) + 1
# 创建一个空的体素数组
volume = np.zeros(volume_shape, dtype=np.float32)
# 将体素网格中的每个体素标记为1
for voxel in voxel_indices:
# 计算体素在体素数组中的索引
grid_index = ((voxel * voxel_size - min_bound) / voxel_size).astype(int)
# 确保索引在有效范围内
if np.all(grid_index >= 0) and np.all(grid_index < volume_shape):
volume[tuple(grid_index)] = 1
# 使用 scikit-image 的 marching_cubes 方法生成等值面
verts, faces, normals, _ = measure.marching_cubes(volume, level=level)
# 将结果转换为 Open3D 的 TriangleMesh
mesh = o3d.geometry.TriangleMesh()
mesh.vertices = o3d.utility.Vector3dVector(verts * voxel_size + min_bound)
mesh.triangles = o3d.utility.Vector3iVector(faces)
mesh.vertex_normals = o3d.utility.Vector3dVector(normals)
return mesh
2.2完整代码
import open3d as o3d
import numpy as np
from skimage import measure
def create_voxel_grid_from_point_cloud(pcd, voxel_size):
"""
将点云转换为体素网格。
参数:
pcd (open3d.geometry.PointCloud): 输入的点云数据。
voxel_size (float): 体素的大小。
返回:
open3d.geometry.VoxelGrid: 生成的体素网格。
"""
voxel_grid = o3d.geometry.VoxelGrid.create_from_point_cloud(pcd, voxel_size=voxel_size)
return voxel_grid
def marching_cubes_from_voxel_grid(voxel_grid, level):
"""
使用行进立方体算法从体素网格生成等值面。
参数:
voxel_grid (open3d.geometry.VoxelGrid): 输入的体素网格。
level (float): 等值面的阈值。
返回:
open3d.geometry.TriangleMesh: 生成的三角网格。
"""
# 将体素网格转换为3D numpy数组
voxel_indices = np.asarray([v.grid_index for v in voxel_grid.get_voxels()])
min_bound = voxel_grid.get_min_bound()
max_bound = voxel_grid.get_max_bound()
voxel_size = voxel_grid.voxel_size
volume_shape = ((max_bound - min_bound) / voxel_size).astype(int) + 1
volume = np.zeros(volume_shape, dtype=np.float32)
for voxel in voxel_indices:
grid_index = ((voxel * voxel_size - min_bound) / voxel_size).astype(int)
if np.all(grid_index >= 0) and np.all(grid_index < volume_shape):
volume[tuple(grid_index)] = 1
# 使用scikit-image的marching_cubes方法生成等值面
verts, faces, normals, _ = measure.marching_cubes(volume, level=level)
# 将结果转换为Open3D的TriangleMesh
mesh = o3d.geometry.TriangleMesh()
mesh.vertices = o3d.utility.Vector3dVector(verts * voxel_size + min_bound)
mesh.triangles = o3d.utility.Vector3iVector(faces)
mesh.vertex_normals = o3d.utility.Vector3dVector(normals)
return mesh
# 加载点云数据
pcd = o3d.io.read_point_cloud("hand.pcd")
o3d.visualization.draw_geometries([pcd], window_name="Cloud", width=800, height=600)
# 估计法向量
pcd.estimate_normals(search_param=o3d.geometry.KDTreeSearchParamHybrid(radius=0.1, max_nn=30))
# 生成体素网格
voxel_size = 0.05
voxel_grid = create_voxel_grid_from_point_cloud(pcd, voxel_size)
# 定义等值面阈值
isosurface_level = 0.5
# 使用行进立方体算法从体素网格生成等值面
mesh = marching_cubes_from_voxel_grid(voxel_grid, isosurface_level)
# 可视化生成的等值面
o3d.visualization.draw_geometries([mesh], window_name="Marching Cubes", width=800, height=600)
# 保存生成的等值面
# o3d.io.write_triangle_mesh("marching_cubes_mesh.ply", mesh)