3D点云开发框架Open3D实战教程（Python）

GitHub练习仓库：github.com/hyperplasma/open3d-practice

1 Open3D简介

Open3D（官网：www.open3d.org）是一个开源库，支持快速开发处理3D数据的软件。其后端用C++实现，经过高度优化并通过Python的前端接口公开。

API文档：www.open3d.org/docs/latest

Open3D提供了三种数据结构：点云（Point Cloud）、网格（Mesh）和RGB-D图像。对于每个表示，Open3D都实现了一整套基本处理算法，如I/O、采样、可视化和数据转换。此外，还包括一些常用的算法，如法线估计、ICP配准等。

内置模块及其功能如下：

模块	功能
Geometry 几何模块	数据结构和基本处理算法
Camera 相机模块	相机模型和相机轨迹
Odometry 里程计模块	RGB-D图像的跟踪与对齐
Registration 配准模块	全局和局部配准
Integration 积分模块	体积积分
I/O 输入输出模块	读写3维数据
Visualization 可视化模块	使用OpenGL呈现3D数据的可自定义GUI
Utility 辅助功能模块	辅助功能，如控制台输出、文件系统和特征包装器
Python模块	Open3D Python绑定和教程

本文将以Open3D为工具，进行对点云数据的相关操作。

---

2 点云基础操作

2.1 RGBD转PCD

import open3d as o3d
import matplotlib.pyplot as plt

color_raw = o3d.io.read_image("img/color_file_1.jpeg")
depth_raw = o3d.io.read_image("img/depth_file_1.png")
# 创建一个rgbd图像
rgbd_image = o3d.geometry.RGBDImage.create_from_color_and_depth(color_raw, depth_raw)
print(rgbd_image)

# 使用matplotlib显示图像
plt.subplot(1, 2, 1)
plt.title('Redwood grayscale image')
plt.imshow(rgbd_image.color)
plt.subplot(1, 2, 2)
plt.title('Redwood depth image')
plt.imshow(rgbd_image.depth)
plt.show()

# rgbd转pcd并显示
pcd = o3d.geometry.PointCloud.create_from_rgbd_image(
    rgbd_image,
    o3d.camera.PinholeCameraIntrinsic(
        o3d.camera.PinholeCameraIntrinsicParameters.PrimeSenseDefault
    )
)
# Flip it, otherwise the pointcloud will be upside down
pcd.transform([[1, 0, 0, 0], [0, -1, 0, 0], [0, 0, -1, 0], [0, 0, 0, 1]])
o3d.visualization.draw_geometries([pcd])
o3d.io.write_point_cloud("pointcloud/test.pcd", pcd, format='auto', write_ascii=False, compressed=False, print_progress=True)

color_file_1.jpeg（颜色数据）：

depth_file_1.png（深度数据）：

2.2 PLY转PCD

import open3d as o3d

def convert_ply_to_pcd(ply_file, pcd_file):
    point_cloud = o3d.io.read_point_cloud(ply_file) # 读取PLY文件
    o3d.io.write_point_cloud(pcd_file, point_cloud) # 保存为PCD文件

ply_file_path = "pointcloud/test_file_1.ply"
pcd_file_path = "pointcloud/test_file_1.pcd"
convert_ply_to_pcd(ply_file_path, pcd_file_path)

其中test_file_1.ply为通过大模型随机生成的PLY格式点云。

---

3 点云空间搜索

无论点云的来源如何，点云终归是离散了的点集，点与点的关联是孤立的。在诸多点云相关的应用中，恢复点云数据的拓扑关系，得到点与点之间的邻近关系，都是首要操作。

例如，通过激光雷达获取的点云数据，具有数据量大、分布不均匀等特点。作为三维领域中一个重要的数据来源，点云主要是表征目标表面的海量点的集合，并不具备传统网格数据的几何拓扑信息，所以点云数据处理中最为核心的问题就是建立离散点间的拓扑关系，实现基于邻域关系的快速查找。

3.1 KD-Tree

KD树（K-Dimensional Tree），是一种分割k维数据空间的数据结构。主要应用于多维空间关键数据的搜索（如：范围搜索和最近邻搜索）。KD树是二进制空间分割树的特殊的情况。 索引结构中相似性查询有两种基本的方式：

1. 范围查询（Range Searches）：给定查询点和查询距离的阈值，从数据集中找出所有与查询点距离小于阈值的数据。
2. K近邻查询（K-Neighbor Searches）：给定查询点及正整数K，从数据集中找到距离查询点最近的K个数据，当K = 1时就是最近邻查询（Nearest Neighbor Searches）。

相关函数：

• o3d.geometry.KDTreeFlann(pcd)：创建KD-Tree
• search_knn_vector_3d(search_pt, k)：K近邻搜索
• search_radius_vector_3d(search_pt，radius)：半径R近邻搜索
• search_hybrid_vector_3d(search_pt, radius, max_nn)：混合邻域搜索，返回半径radius内不超过max_nn个近邻点

3.1.1 原理

KD-Tree是一个二叉树，每个节点表示一个空间范围：

• node-data：数据矢量，数据集中具体的某个数据点。是k维矢量。
• range：空间矢量，该节点所代表的空间范围。
• split：垂直于分割超平面的方向轴序号。整数。
• left：由位于该节点分割超平面左子空间内所有数据点所构成的KD-Tree。
• right：由位于该节点分割超平面右子空间内所有数据点所构成的KD-Tree。
• parent：父节点。

以下通过一个简单直观的例子介绍KD-Tree算法——

【例】有6个二维数据点{(2, 3), (5, 4), (9, 6), (4, 7), (8, 1), (7, 2)}位于二维空间内。KD-Tree算法的目标即为确定如下所示的二维数据k-d树空间划分示意图中这些分割空间的分割线（多维空间即为分割平面，一般为超平面）。

以下将逐步展示KD-Tree是如何确定这些分割线的：

由于此例简单，数据维度只有2维，故可简单地给x、y两个方向轴编号为0、1，即split = {0, 1}。具体步骤如下——

1. 确定split域首先该取的值：分别计算x、y方向上数据的方差可知x方向上的方差最大，故split域值首先取0，即x轴方向。
2. 确定node-data的域值：根据x轴方向的值2, 5, 9, 4, 8, 7排序选出中值为7，故node-data = (7, 2)。由此可得，该节点的分割超平面即为通过(7, 2)并垂直于split = 0（x轴）的直线x = 7。
3. 确定左子空间和右子空间：分割超平面x = 7将整个空间分为两部分，如下图所示。x <= 7的部分为左子空间，包含3个节点{(2, 3), (5, 4), (4, 7)}；另一部分为右子空间，包含2个节点{(9, 6), (8, 1)}。

4. KD-Tree的构建是一个递归的过程，对左子空间和右子空间内的数据重复根节点的过程可得下一级子节点(5, 4)和(9, 6)（即左右子空间的根节点），同时将空间和数据集进一步细分。如此反复直到空间中只包含一个数据点。最后生成的KD-Tree如下图所示。

3.1.2 搜索算法

在KD-Tree中进行数据的查找是特征匹配的重要环节，其目的是检索在KD-Tree中与查询点距离最近的数据点。以下通过一个简单的实例来描述最近邻查找的基本思路——

【例】如下图所示，星号表示要查询的点(2.1, 3.1)。通过二叉搜索，顺着搜索路径很快就能找到最邻近的近似点，即叶子节点(2, 3)，但找到的叶子节点并不一定就是最邻近的（最邻近者肯定距离查询点更近）。为了找到真正的最近邻，还需进行回溯操作：算法沿搜索路径反向查找是否有距离查询点更近的数据点。

此例中先从(7, 2)点开始进行二叉查找，然后到达(5, 4)，最后到达(2, 3)，此时搜索路径中的节点为<(7, 2), (5, 4), (2, 3)>，首先以(2, 3)作为当前最近邻点，计算其到查询点(2.1, 3.1)的距离为0.1414，然后回溯至其父节点(5, 4)，并判断在该父节点的其他子节点空间中是否有距离查询点更近的数据点。以(2.1, 3.1)为圆心、以0.1414为半径作圆（如上图所示），发现该圆并不与超平面y = 4交割，因此无需进入(5, 4)节点右子空间中去搜索。

3.1.3 应用：利用KD-Tree评估点云分辨率

点云分辨率（点间距，密度）用于评价点云的疏密程度。

常用计算方法：针对点云中的每个点P，计算其各个最近点Q到它的欧式距离，再将所有点的计算结果取均值。

import open3d as o3d
import numpy as np

print("->正在加载点云... ")
pcd = o3d.io.read_point_cloud("pointcloud/res_test.pcd")
print(pcd)

# 将点云设置为灰色
pcd.paint_uniform_color([0.5, 0.5, 0.5])

# 建立KDTree
pcd_tree = o3d.geometry.KDTreeFlann(pcd)

# 将第1500个点设置为紫色
pcd.colors[1500] = [0.5, 0, 0.5]

# 使用K近邻，将第1500个点最近的5000个点设置为蓝色
print("使用K近邻，将第1500个点最近的5000个点设置为蓝色")
k = 2000   # 设置K的大小
[num_k, idx_k, _] = pcd_tree.search_knn_vector_3d(pcd.points[1500], k)    # 返回邻域点的个数和索引
np.asarray(pcd.colors)[idx_k[1:], :] = [0, 0, 1]  # 跳过最近邻点（查询点本身）进行赋色
print("k邻域内的点数为：", num_k)

# 使用半径R近邻，将第1500个点半径（0.02）范围内的点设置为红色
print("使用半径R近邻，将第1500个点半径（0.02）范围内的点设置为红色")
radius = 0.02   # 设置半径大小
[num_radius, idx_radius, _] = pcd_tree.search_radius_vector_3d(pcd.points[1500], radius)   # 返回邻域点的个数和索引
np.asarray(pcd.colors)[idx_radius[1:], :] = [1, 0, 0]  # 跳过最近邻点（查询点本身）进行赋色
print("半径r邻域内的点数为：", num_radius)

# 使用混合邻域，将半径R邻域内不超过max_num个点设置为绿色
print("使用混合邻域，将第1500个点半径R邻域内不超过max_num个点设置为绿色")
max_nn = 200   # 半径R邻域内最大点数
[num_hybrid, idx_hybrid, _] = pcd_tree.search_hybrid_vector_3d(pcd.points[1500], radius, max_nn)
np.asarray(pcd.colors)[idx_hybrid[1:], :] = [0, 1, 0]  # 跳过最近邻点（查询点本身）进行赋色
print("混合邻域内的点数为：", num_hybrid)

print("->正在可视化点云...")
o3d.visualization.draw_geometries([pcd])

示例结果：

3.2 Octree

八叉树（Octree）是一种树数据结构，其中每个内部节点有8个子节点，通常用于三维点云的空间划分。八叉树的非空叶节点包含属于同一空间细分的一个或多个点。八叉树是三维空间的有用描述，可用于快速査找附近的点。Open3D具有几何体类型的八叉树，可用于创建、搜索和遍历具有用户指定的最大树深度的八叉树max_depth。

如下图所示，左边是八叉树地图，右面是八叉树：

可以使用convert_from_point_cloud()从点云构造Octree。通过遵循从根节点到“最大深度”(max_depth)处的相应叶节点的路径，将每个点插入到树中。随着树深度的增加，内部（最终是叶子）节点表示三维空间的较小分区。

若点云具有颜色，则对应的叶节点将采用上次插入点的颜色。 可设置size_expand参数来增加根八叉树节点的大小，使其略大于原始点云边界以容纳所有点。

import open3d as o3d
import numpy as np

# 加载点云
print("->正在加载点云... ")
pcd = o3d.io.read_point_cloud("pointcloud/example.pcd")
print("原始点云：", pcd)

# 构建Octree
print("Octree分割")
octree = o3d.geometry.Octree(max_depth=5)
octree.convert_from_point_cloud(pcd, size_expand=0.1)
print("->正在可视化Octree...")
o3d.visualization.draw_geometries([octree])

max_depth的值越大，八叉树的节点越细腻。以下是max_depth = 5时的结果：

---

4 点云滤波

点云滤波的作用是去除数据中的离群点，并且平滑三维数据，使其更加容易处理。

体素（Voxel）是体积元素（Volume Pixel）的简称，包含体素的立体可以通过立体渲染或提取给定阈值轮廓的多边形等值面表现出来。体素是数字数据于三维空间分割上的最小单位，用于三维成像、科学数据与医学影像等领域。概念上类似二维空间的最小单位——像素，像素用在二维计算机图像的影像数据上。某些真正的三维显示器运用体素来描述它们的分辨率，例如可以显示512×512×512体素的显示器。

4.1 点云体素化

import open3d as o3d
import numpy as np

# 加载点云
print("->正在加载点云... ")
pcd = o3d.io.read_point_cloud("pointcloud/example.pcd")
print("原始点云：", pcd)

# 体素化点云
print("执行体素化点云")
voxel_grid = o3d.geometry.VoxelGrid.create_from_point_cloud(pcd, voxel_size=0.005)
print("->正在可视化体素...")
o3d.visualization.draw_geometries([voxel_grid])

4.2 体素下采样

import open3d as o3d
import numpy as np

# 点云体素化函数
def get_mesh(file_path):
    mesh = o3d.io.read_triangle_mesh(file_path)
    mesh.compute_vertex_normals()
    return mesh

# 点云体素化
print("->正在进行点云体素化...")
file_path = "pointcloud/example.ply"
N = 2000    # 将点划分为N个体素
pcd = get_mesh(file_path).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)))
print("体素下采样点云：", pcd)
print("->正在可视化体素下采样点云...")
o3d.visualization.draw_geometries([pcd])

print("执行体素化点云")
voxel_grid = o3d.geometry.VoxelGrid.create_from_point_cloud(pcd, voxel_size=0.05)
print("->正在可视化体素...")
o3d.visualization.draw_geometries([voxel_grid])

体素下采样使用常规体素网格从输入点云创建统一下采样的点云，通常用作许多点云处理任务的预处理步骤。该算法分两个步骤运行：

1. 将点云进行进行体素划分，点被存储到体素中。
2. 每个占用的体素通过平均内部的所有点来生成精确的一个点。一般而言取体素内点云的质心作为该体素的点的位置。

【例】另一个简单的体素下采样案例：

import open3d as o3d
import numpy as np

print("->正在加载点云... ")
pcd = o3d.io.read_point_cloud("pointcloud/example.ply")
print(pcd)

print("->正在可视化原始点云...")
# o3d.visualization.draw_geometries([pcd])

print("->正在体素下采样...")
voxel_size = 0.001
downpcd = pcd.voxel_down_sample(voxel_size)
print(downpcd)

原始点云图像可视化结果：

voxel_size = 0.01时的体素下采样结果：

voxel_size = 0.001时的体素下采样结果：

4.3 统计滤波

statistical outlier_removal()会移除距离相邻点更远的点，需要两个参数：

• num_neighbors：指定为了计算给定点的平均距离，需要考虑多少个邻居。
• std_ratio：K邻域的点数允许根据点云平均距离的标准偏差设置阈值水平。该数值越低，滤除的点数就越多。

import open3d as o3d
import numpy as np

print("->正在加载点云... ")
pcd = o3d.io.read_point_cloud("pointcloud/res_test.pcd")
print("原始点云：", pcd)

# 统计滤波
print("->正在进行统计滤波...")
num_neighbors = 200 # K邻域点的个数
std_ratio = 2.0 # 标准差乘数
# 执行统计滤波，返回滤波后的点云sor_pcd和对应的索引ind
sor_pcd, ind = pcd.remove_statistical_outlier(num_neighbors, std_ratio)
sor_pcd.paint_uniform_color([0, 0, 1])
print("统计滤波后的点云：", sor_pcd)
sor_pcd.paint_uniform_color([0, 0, 1])
# 提取噪声点云
sor_noise_pcd = pcd.select_by_index(ind, invert=True)
print("噪声点云：", sor_noise_pcd)
sor_noise_pcd.paint_uniform_color([1, 0, 0])

# 可视化统计滤波后的点云和噪声点云
o3d.visualization.draw_geometries([sor_pcd, sor_noise_pcd])

统计滤波后的可视化效果：

噪声点云可视化：

两者放到一起的效果：

4.4 半径滤波

radius_outlier_removal()移除给定球体中几乎没有邻居的点，需要两个参数:

• num_points：邻域球内的最少点数，低于该值的点为噪声点
• radius：为噪声点邻域半径的大小

import open3d as o3d
import numpy as np

print("->正在加载点云... ")
pcd = o3d.io.read_point_cloud("pointcloud/res_test_b.pcd")
print("原始点云：", pcd)

# 半径滤波
print("->正在进行半径滤波...")
num_points = 600  # 邻域球内的最少点数，低于该值的点为噪声点
radius = 0.5    # 邻域半径大小
# 执行半径滤波，返回滤波后的点云sor_pcd和对应的索引ind
sor_pcd, ind = pcd.remove_radius_outlier(num_points, radius)
sor_pcd.paint_uniform_color([0, 0, 1])
print("半径滤波后的点云：", sor_pcd)
sor_pcd.paint_uniform_color([0, 0, 1])
# 提取噪声点云
sor_noise_pcd = pcd.select_by_index(ind, invert=True)
print("噪声点云：", sor_noise_pcd)
sor_noise_pcd.paint_uniform_color([1, 0, 0])

# 可视化半径滤波后的点云和噪声点云
o3d.visualization.draw_geometries([sor_pcd, sor_noise_pcd])

---

5 点云特征提取：法线估计

平面的法线（Normal Line）是垂直于它的单位向量。点云表面的法线被定义为垂直于与点云表面相切的平面的向量。表面法线也可以计算点云中一点的法线，这被认为是一种十分重要的性质，在计算机视觉中应用甚广，例如用于推算出光源的位置。

estimate_normals()可用于计算每个点的法线。该函数查找相邻点并使用协方差分析计算相邻点的主轴。该函数的search_param参数需要KDTreeSearchParamHybrid类的实例，该类的两个关键构造参数radius和max_nn指定搜索半径和最大最近邻居。

以下代码中指定搜索半径为10厘米，最多考虑30个邻居，以节省计算时间。

import open3d as o3d
import numpy as np

print("->正在加载点云... ")
pcd = o3d.io.read_point_cloud("pointcloud/example.pcd")
print(pcd)

print("->正在估计法线并可视化...")
radius = 0.01   # 搜索半径
max_nn = 30     # 邻域内用于估算法线的最大点数
pcd.estimate_normals(search_param=o3d.geometry.KDTreeSearchParamHybrid(radius, max_nn)) # 执行法线估计
o3d.visualization.draw_geometries([pcd], point_show_normal=True)

print("->正在打印前10个点的法向量...")
print(np.asarray(pcd.normals)[:10, :])

---

6 点云分割

点云分割的目的是分块，从而便于单独处理。根据空间、几何和纹理等特征点进行划分，同一划分内的点云拥有相似的特征。

6.1 DBSCAN聚类分割

DBSCAN（Density-Based Spatial Clustering of Applications with Noise，基于密度的带噪声应用空间聚类）通过点的密度来识别聚类集群。聚类通常位于高密度区域，而异常值往往位于低密度区域。

根据该算法的先驱者的说法，DBSCAN具有3个主要优势：只需要最少的领域知识，可以发现任意形状的聚类，并且对于大型数据库来说非常高效。

Open3D通过cluster_dbscan()实现了DBSCAN，需要两个参数：

• eps：同一簇内的最大点间距。
• min_points：有效聚类的最小点数。

函数返回标签label，其中label = -1表示噪声。

import open3d as o3d
import numpy as np
import matplotlib.pyplot as plt

print("->正在加载点云... ")
pcd = o3d.io.read_point_cloud("pointcloud/test.pcd")
print(pcd)

print("->正在DBSCAN聚类...")
eps = 0.5           # 同一聚类中最大点间距
min_points = 50     # 有效聚类的最小点数
with o3d.utility.VerbosityContextManager(o3d.utility.VerbosityLevel.Debug) as cm:
    labels = np.array(pcd.cluster_dbscan(eps, min_points, print_progress=True))
max_label = labels.max()    # 获取聚类标签的最大值[-1,0,1,2,...,max_label]，label=-1为噪声，因此总聚类个数为max_label + 1
print(f"Point cloud has {max_label + 1} clusters.")
colors = plt.get_cmap("tab20")(labels / (max_label if max_label > 0 else 1))
colors[labels < 0] = 0  # labels = -1 的簇为噪声，以黑色显示
pcd.colors = o3d.utility.Vector3dVector(colors[:, :3])
o3d.visualization.draw_geometries([pcd])

对室内场景进行点云分割（运行耗时极长）：

6.2 RANSAC平面分割

RANSAC（Random Sample Consensus）是一种迭代的参数估计算法，主要用于从包含大量噪声数据的样本中估计模型参数。其核心思想是通过随机采样和模型验证来找到数据中最符合模型假设的点。

算法步骤：

1. 初始化：设置最大迭代次数max_iterations和内点阈值distance_threshold。
2. 随机采样：从数据集中随机选择最小数量的样本点来拟合模型（例如拟合平面需要3个点）。
3. 模型估计：使用选定的样本点计算模型参数（例如平面的方程）。
4. 模型验证：计算所有数据点到模型的距离，将距离小于distance_threshold的点标记为内点。
5. 评估模型：计算内点的数量，若内点数量超过预定的阈值并且模型质量优于之前的模型，则更新最佳模型。
6. 重复上述步骤，直至达到最大迭代次数或找到最优模型。

Open3D提供了segment_plane()方法，需要3个关键参数，方法返回模型系数plane_model和内点索引inliers：

• distance_threshold：用于定义点到平面的最大距离阈值，单位为米。如果点到平面的距离小于此阈值，则认为该点属于平面。默认值为0.01。
• ransac_n：用于RANSAC算法的迭代次数。RANSAC是一种随机采样一致性算法，用于估计平面模型参数。默认值为3。
• num_iterations：在RANSAC算法中每次迭代时使用的随机样本点的数量。默认值为10000。

RANSAC平面分割算法广泛应用于计算机视觉和点云处理领域，特别适用于以下场景：

• 点云平面分割：从三维点云数据中提取平面，如地面、墙壁等。
• 图像拼接：用于匹配图像中的特征点并估计变换矩阵。
• 3D物体识别：从点云数据中提取特定形状或结构。

import open3d as o3d

print("->正在加载点云... ")
pcd = o3d.io.read_point_cloud("pointcloud/res_test.pcd")
print(pcd)

print("->正在RANSAC平面分割...")
distance_threshold = 0.2    # 内点到平面模型的最大距离
ransac_n = 3        # 用于拟合平面的采样点数
num_iterations = 1000   # 最大迭代次数

# 返回模型系数plane_model和内点索引inliers
plane_model, inliers = pcd.segment_plane(distance_threshold, ransac_n, num_iterations)

# 输出平面方程
[a, b, c, d] = plane_model
print(f"Plane equation: {a:.2f}x + {b:.2f}y + {c:.2f}z + {d:.2f} = 0")

# 平面内点点云，标记为蓝色
inlier_cloud = pcd.select_by_index(inliers)
inlier_cloud.paint_uniform_color([0, 0, 1.0])
print(inlier_cloud)

# 平面外点点云，标记为红色
outlier_cloud = pcd.select_by_index(inliers, invert=True)
outlier_cloud.paint_uniform_color([1.0, 0, 0])
print(outlier_cloud)

# 可视化平面分割结果
o3d.visualization.draw_geometries([inlier_cloud, outlier_cloud])

原始点云可视化：

RANSAC平面分割效果：

然而事实上点云中有多个平面，上述方法只能分割出一个平面。要想分割多个平面，可以参照以下做法：

import open3d as o3d
import numpy as np

def plane_detection(pcd, tolerance=50):
    current_pcd = pcd
    planes = []
    while len(current_pcd.points) > tolerance:
        plane_model, inliers = current_pcd.segment_plane(distance_threshold=0.3, ransac_n=5, num_iterations=1000)
        if len(inliers) < tolerance:
            break
        inlier_indices = np.asarray(inliers)
        inlier_cloud = current_pcd.select_by_index(inlier_indices)
        current_pcd = current_pcd.select_by_index(inlier_indices, invert=True)

        normal_vector = plane_model[:3]
        point_in_plane = -normal_vector * plane_model[3] / np.linalg.norm(normal_vector) ** 2
        endpoint = point_in_plane + normal_vector * 2

        line = o3d.geometry.LineSet()
        line.points = o3d.utility.Vector3dVector([point_in_plane, endpoint])
        line.lines = o3d.utility.Vector2iVector([[0, 1]])

        planes.append(line)
        planes.append(inlier_cloud)

    return current_pcd, planes

def main(file_path):
    pcd = o3d.io.read_point_cloud(file_path)

    remain_pcd, planes = plane_detection(pcd)
    for plane in planes:
        plane.paint_uniform_color(np.random.rand(3))

    mesh_frame = o3d.geometry.TriangleMesh.create_coordinate_frame(size=5.0, origin=[0,0,0])
    # 可视化结果
    o3d.visualization.draw_geometries([remain_pcd, *planes, mesh_frame])

main("pointcloud/res_test.pcd")

6.3 隐藏点剔除

若希望从给定的视点渲染点云，但背景中的某些点会泄漏到前景中（它们不被其他点遮挡），则可以应用隐藏点移除算法。

在Open3D中，实现了Katz2007的方法，该方法从给定视图近似点云的可见性，无需曲面重建或法线估计。

import open3d as o3d
import numpy as np

print("->正在加载点云... ")
pcd = o3d.io.read_point_cloud("pointcloud/res_test.pcd")
print(pcd)

print("->正在剔除隐藏点...")
diameter = np.linalg.norm(np.asarray(pcd.get_max_bound()) - np.asarray(pcd.get_min_bound()))
print("定义隐藏点去除的参数")
camera = [0, 0, diameter]   # 视点位置
radius = diameter * 100 # 噪声点云半径（The radius of the sperical projection）
_, pt_map = pcd.hidden_point_removal(camera, radius)    # 获取视点位置能看到的所有点的索引pt_map

# 可视点点云
pcd_visible = pcd.select_by_index(pt_map)
pcd_visible.paint_uniform_color([0, 0, 1])  # 可视点为蓝色
print("可视点个数为：", pcd_visible)
# 隐藏点点云
pcd_hidden = pcd.select_by_index(pt_map, invert=True)
pcd_hidden.paint_uniform_color([1, 0, 0])   # 隐藏点为红色
print("隐藏点个数为：", pcd_hidden)
print("->正在可视化可视点和隐藏点点云...")
o3d.visualization.draw_geometries([pcd_visible, pcd_hidden])

---

7 点云曲面重建：Alpha Shapes

曲面重建技术在逆向工程、数据可视化、机器视觉、虚拟现实、医疗技术等领域中得到了广泛的应用。 例如在汽车、航空等工业领域中，复杂外形产品的设计仍需要根据手工模型，采用逆向工程的手段建立产品的数字化模型；根据测量数据建立人体以及骨骼和器官的计算机模型，在医学、定制生产等方面都有重要意义。

在许多情况下希望生成密集的三维几何体，即三角形网格。然而，从多视图立体方法或深度传感器中只能获得非结构化点云。要从非结构化输入中获得三角形网格，需要执行点云曲面重建。

Alpha Shapes（α-形状）是一种用于点云边界重建的算法。通过调节alpha参数，可以生成不同分辨率的边界形状。Alpha Shapes通过将点云中的点连成三角形，并移除长于alpha的边，从而生成符合数据特征的形状。

Open3D实现了create_from_point_cloud_alpha_shape()方法。

import open3d as o3d

# 加载点云
print("->正在加载点云... ")
pcd = o3d.io.read_point_cloud("pointcloud/res_test.pcd")
print("原始点云：", pcd)

# Alpha shapes
alpha = 0.03
print(f"alpha={alpha:.3f}")
mesh = o3d.geometry.TriangleMesh.create_from_point_cloud_alpha_shape(pcd, alpha)
mesh.compute_vertex_normals()
o3d.visualization.draw_geometries([mesh], mesh_show_back_face=True)

alpha = 0.03时效果如下：

alpha = 0.1时效果如下：