【小白深度教程 1.23】手把手教你使用 Open3D(6)对点云数据进行全局配准(RANSAC 算法)
1. 全局配准的概念
ICP 配准和彩色点云配准都属于局部配准方法,因为它们依赖于粗略对齐作为初始化。本教程展示了另一类配准方法,称为全局配准。该类算法不需要对齐作为初始化,通常生成较不紧密的对齐结果,并用于局部方法的初始化。
2. 可视化
以下辅助函数可视化变换后的源点云与目标点云:
def draw_registration_result(source, target, transformation):
source_temp = copy.deepcopy(source)
target_temp = copy.deepcopy(target)
source_temp.paint_uniform_color([1, 0.706, 0])
target_temp.paint_uniform_color([0, 0.651, 0.929])
source_temp.transform(transformation)
o3d.visualization.draw_geometries([source_temp, target_temp],
zoom=0.4559,
front=[0.6452, -0.3036, -0.7011],
lookat=[1.9892, 2.0208, 1.8945],
up=[-0.2779, -0.9482, 0.1556])
3. 提取几何特征
我们对点云进行降采样、估计法线,然后为每个点计算 FPFH 特征。FPFH 特征是一个 33 维的向量,用于描述点的局部几何属性。在 33 维空间中进行最近邻查询,可以返回具有相似局部几何结构的点。详情请参见 [Rasu2009]。
def preprocess_point_cloud(pcd, voxel_size):
print(":: Downsample with a voxel size %.3f." % voxel_size)
pcd_down = pcd.voxel_down_sample(voxel_size)
radius_normal = voxel_size * 2
print(":: Estimate normal with search radius %.3f." % radius_normal)
pcd_down.estimate_normals(
o3d.geometry.KDTreeSearchParamHybrid(radius=radius_normal, max_nn=30))
radius_feature = voxel_size * 5
print(":: Compute FPFH feature with search radius %.3f." % radius_feature)
pcd_fpfh = o3d.pipelines.registration.compute_fpfh_feature(
pcd_down,
o3d.geometry.KDTreeSearchParamHybrid(radius=radius_feature, max_nn=100))
return pcd_down, pcd_fpfh
4. 处理输入数据
下面的代码从两个文件中读取源点云和目标点云。它们以单位矩阵作为变换进行未对齐的初始化。
def prepare_dataset(voxel_size):
print(":: Load two point clouds and disturb initial pose.")
source = o3d.io.read_point_cloud("../../test_data/ICP/cloud_bin_0.pcd")
target = o3d.io.read_point_cloud("../../test_data/ICP/cloud_bin_1.pcd")
trans_init = np.asarray([[0.0, 0.0, 1.0, 0.0], [1.0, 0.0, 0.0, 0.0],
[0.0, 1.0, 0.0, 0.0], [0.0, 0.0, 0.0, 1.0]])
source.transform(trans_init)
draw_registration_result(source, target, np.identity(4))
source_down, source_fpfh = preprocess_point_cloud(source, voxel_size)
target_down, target_fpfh = preprocess_point_cloud(target, voxel_size)
return source, target, source_down, target_down, source_fpfh, target_fpfh
voxel_size = 0.05 # means 5cm for this dataset
source, target, source_down, target_down, source_fpfh, target_fpfh = prepare_dataset(
voxel_size)
对点云数据进行全局配准(RANSAC 算法)_files/60a457119e264d4cbf6058bc992ee8ec.png)
5. RANSAC
我们使用 RANSAC 进行全局配准。在每次 RANSAC 迭代中,从源点云中随机选取 ransac_n 个点。在目标点云中的对应点是通过在 33 维 FPFH 特征空间中查询最近邻来检测的。修剪步骤使用快速修剪算法快速拒绝早期的错误匹配。
Open3D 提供以下修剪算法:
CorrespondenceCheckerBasedOnDistance检查对齐的点云是否接近(小于指定的阈值)。CorrespondenceCheckerBasedOnEdgeLength检查源和目标对应关系中任意两条边(由两个顶点形成)的长度是否相似。本教程检查以下条件是否为真:(|edge_{source}| > 0.9 \times |edge_{target}|) 和 (|edge_{target}| > 0.9 \times |edge_{source}|)。CorrespondenceCheckerBasedOnNormal考虑任意对应关系的顶点法线相似性。它计算两个法向量的点积,并接受一个弧度值作为阈值。
只有通过修剪步骤的匹配才用于计算变换,并在整个点云上进行验证。核心函数是 registration_ransac_based_on_feature_matching。此函数的最重要超参数是 RANSACConvergenceCriteria,它定义了 RANSAC 迭代的最大次数和最大验证步骤。参数值越大,结果越准确,但算法耗时也越长。
我们根据 [Choi2015] 提供的经验值设置 RANSAC 参数。
def execute_global_registration(source_down, target_down, source_fpfh,
target_fpfh, voxel_size):
distance_threshold = voxel_size * 1.5
print(":: RANSAC registration on downsampled point clouds.")
print(" Since the downsampling voxel size is %.3f," % voxel_size)
print(" we use a liberal distance threshold %.3f." % distance_threshold)
result = o3d.pipelines.registration.registration_ransac_based_on_feature_matching(
source_down, target_down, source_fpfh, target_fpfh, distance_threshold,
o3d.pipelines.registration.TransformationEstimationPointToPoint(False),
4, [
o3d.pipelines.registration.CorrespondenceCheckerBasedOnEdgeLength(
0.9),
o3d.pipelines.registration.CorrespondenceCheckerBasedOnDistance(
distance_threshold)
], o3d.pipelines.registration.RANSACConvergenceCriteria(4000000, 500))
return result
result_ransac = execute_global_registration(source_down, target_down,
source_fpfh, target_fpfh,
voxel_size)
print(result_ransac)
draw_registration_result(source_down, target_down, result_ransac.transformation)
对点云数据进行全局配准(RANSAC 算法)_files/b737e00cce6249dea6e2e84deede568d.png)
Open3D 提供了更快的全局配准实现。请参考快速全局配准。
6. 局部精细化
出于性能考虑,全局配准仅在大幅降采样的点云上执行,结果也较不紧密。我们使用点到平面 ICP 进一步优化对齐效果。
def refine_registration(source, target, source_fpfh, target_fpfh, voxel_size):
distance_threshold = voxel_size * 0.4
print(":: Point-to-plane ICP registration is applied on original point")
print(" clouds to refine the alignment. This time we use a strict")
print(" distance threshold %.3f." % distance_threshold)
result = o3d.pipelines.registration.registration_icp(
source, target, distance_threshold, result_ransac.transformation,
o3d.pipelines.registration.TransformationEstimationPointToPlane())
return result
result_icp = refine_registration(source, target, source_fpfh, target_fpfh,
voxel_size)
print(result_icp)
draw_registration_result(source, target, result_icp.transformation)
对点云数据进行全局配准(RANSAC 算法)_files/b8a96b8d24ce4ee9842a446c3b20e153.png)
7. 快速全局配准
基于 RANSAC 的全局配准方案可能由于无数次模型提议和评估而耗费大量时间。[Zhou2016] 提出了一种更快的方法,快速优化少量对应关系的线过程权重。由于每次迭代中不涉及模型提议和评估,[Zhou2016] 提出的方法可以节省大量计算时间。
本教程比较了基于 RANSAC 的全局配准与 [Zhou2016] 实现的运行时间。
8. 处理输入数据
我们使用与上述全局配准示例中相同的输入。
voxel_size = 0.05 # means 5cm for the dataset
source, target, source_down, target_down, source_fpfh, target_fpfh = \
prepare_dataset(voxel_size)
对点云数据进行全局配准(RANSAC 算法)_files/28c2017fc82549019338fdd7b76562b3.png)
9. Baseline
在下面的代码中,我们对全局注册方法进行计时。
start = time.time()
result_ransac = execute_global_registration(source_down, target_down,
source_fpfh, target_fpfh,
voxel_size)
print("Global registration took %.3f sec.\n" % (time.time() - start))
print(result_ransac)
draw_registration_result(source_down, target_down, result_ransac.transformation)
对点云数据进行全局配准(RANSAC 算法)_files/6838c0569eda4018b3bde6bdb5441efa.png)
10. 快速全局配准
使用相同的基线输入,下面的代码调用 Zhou2016 的实现。
def execute_fast_global_registration(source_down, target_down, source_fpfh,
target_fpfh, voxel_size):
distance_threshold = voxel_size * 0.5
print(":: Apply fast global registration with distance threshold %.3f" \
% distance_threshold)
result = o3d.pipelines.registration.registration_fast_based_on_feature_matching(
source_down, target_down, source_fpfh, target_fpfh,
o3d.pipelines.registration.FastGlobalRegistrationOption(
maximum_correspondence_distance=distance_threshold))
return result
start = time.time()
result_fast = execute_fast_global_registration(source_down, target_down,
source_fpfh, target_fpfh,
voxel_size)
print("Fast global registration took %.3f sec.\n" % (time.time() - start))
print(result_fast)
draw_registration_result(source_down, target_down, result_fast.transformation)
对点云数据进行全局配准(RANSAC 算法)_files/c4180eb7bb6b4ae2b5acc2be4a2290b0.png)