CN110070567A - 一种地面激光点云配准方法 - Google Patents

Abstract

本发明公开了一种地面激光点云配准方法，通过原始点云数据自动估计点云的平均间距和长度，获取建立点云金字塔的体素格网大小。在金字塔点云顶层，建立所有点的FPFH特征，然后通过双向一致性检测获得初始匹配，可以快速有效实现点云的初始配准。在此基础上，利用初始匹配结果邻域内的的三维点进行由粗到精的ICP优化，在大幅减少参与ICP点云数量的同时，同名点邻域内的点云信息更丰富，有利于ICP算法迭代快速收敛。

Description

一种地面激光点云配准方法

技术领域

本发明属于三维重建技术领域，涉及两站带有适当重叠区域的点云数据之间的配准。

背景技术

近年来，激光扫描技术发展迅速。激光扫描技术具有高精度，高效率的特点，能够快速并精确的获取目标相对于扫描仪的距离，从而能够获取被扫描目标表明的点相对于扫描仪的精确的三维点坐标。随着技术的发展，点云获取越来越容易，成本也越来越低。点云处理技术在逆向工程、人体扫描、文物保护、虚拟现实、增强现实、三维城市重建、自动驾驶以及机器人视觉等方面应用广泛。由于对象或者场景的复杂性、激光传感器扫描范围的局限性等原因，往往需要对对象多次扫描或将一个大的场景分为多个区域分别扫描。多次扫描的数据具有不同的坐标系统，点云处理技术中的基础和至关重要的技术就是点云配准，将不同的扫描结果统一到同一个坐标系统下面。主流的点云配准算法采用由粗到精的配准流程，首先完成初始配准，通过人机交互旋转三个以上的同名点从而计算转换参数，或者是计算点的特征(如persist analysis FPFH，ISS FPFH)或者是直接的随机采样(如Keypoint-based4 Points Congruent Sets,Super4points Congruent Sets等)完成初始配准，接着在初始配准结果的基础上采用ICP(Iterative Closest Point)等方法完成精确配准，最终实现点云的配准。无论是persist analysis FPFH还是ISS FPFH都需要对大量点和点的邻域计算和分析，计算量大，计算效率不高，Keypoint-based4 Points CongruentSets依赖于特征点提取，对特征点提取要求较高，需要能够提取适当数目的特征点,并且它还和Super4points Congruent Sets一样需要设置先验的点云重叠参数，否则会导致配准错误或者是配准花费大量的时间。这些方法在抗噪性能、配准的鲁棒性和配准效率方面都存在不足。

发明内容

本发明所要解决的技术问题是，针对现有技术不足，提供一种地面激光点云配准方法，提高配准速度和配准精度。

为解决上述技术问题，本发明所采用的技术方案是：步骤1：读入由点云获取设备采集的具有不同视角的源点云S₀和目标点云T₀两视角点云；

步骤2：自动估计点云的平均点间距d_mean和平均长度L_mean；

步骤3：点云下采样，构建源点云与目标点云金字塔S₀、S₁、S₂和T₀、T₁、T₂；

步骤4：计算点云S₂与T₂中所有点的法向量和所有点的FPFH(Fast Point FeatureHistogram)描述子；

步骤5：双向匹配S₂与T₂中所有点的FPFH描述子，得到匹配对应关系C；

步骤6：利用SAC-IA(随机采样一致性)算法从C中将内点C₂筛选出来，同时得到初始配准转换参数Tr₂；

步骤7：在点云S₁与T₁上利用步骤6中得到的内点C₂和初始转换参数Tr₂，提取C₂邻域内的点云，用于第一次ICP优化，得到第一层优化的转换参数Tr₁；

步骤8：在点云S₀与T₀上利用步骤6中的内点C₂和步骤9中初步ICP优化后的转换参数Tr₁，采用最近邻点代替原始对应点，仍然采用点邻域来完成第二次ICP精配准优化，得到最终配准结果。

与现有技术相比，本发明所具有的有益效果为：本发明配准速度快，配准精度高，本发明利用平均点间距与平均点云长度确定相关参数、建立金字塔结构点云、所有点的FPFH特征描述子直接双向匹配和特征点引导的ICP精配准方法，确保了快速有效地配准类似地面激光扫描仪获取的海量的点云。

附图说明

图1是配准流程图；

图2是输入数据示意图；

图3是金字塔点云示意图；

图4初始匹配对示意图；(a)双向匹配结果；(b)经过SAC-IC保留下的对应关；

图5由粗到精的配准过程示意图；(a)Tr₂配准结果；(a)Tr₀配准结果。

具体实施方式

步骤1：读入由点云获取设备采集的具有不同视角的源点云S⁰和目标点云T⁰两视角点云；

步骤2：自动估计点云的平均点间距d_mean和平均长度L_mean；

步骤2.1：在两个点云中分别随机采样K₀(K₀＝1000)个采样点，对每一个采样点P_i的k₁邻域Q进行分析，k₁取经验值4，计算邻域中每一个点P_j到点P_i的距离d_ij，计算点P_i的邻域点到P_i的平均距离作为P_i点附近的点云平均间距的估计值d_i；

步骤2.2：用同样的方法计算出每一个采样点附近的点云平均点间距d_i,然后将所有采样点附近的点云平均点间距的平均值作为整个点云的平均点间距d，

步骤2.3：用相同的方法分别计算出源点云和目标点云的平均点间距d_src d_tgt，以它们的平均值作为最终的平均点间距d_mean，

步骤2.4：在两个点云中分别随机选取k₂(k₂＝2000)个采样点，，判断采样点中x,y,z坐标的上下界，x_max,x_min,y_min,y_max,z_max,z_min，

步骤2.5：计算由采样点所代表子点云的所在三维空间立方体的斜对角线长，作为点云的长度L,

步骤2.6：分别计算源点云和目标点云的长度L_src,L_tgt，以他们的平均值作为最终的平均点云长度L_mean。

步骤3：点云下采样，构建源点云与目标点云金字塔S₀、S₁、S₂和T₀、T₁、T₂；

步骤3.1：获取输入的点云的点数n，当点云点数在100万以上时，采用等间隔采样的方法进行下采样，计算间隔因子Δ，每隔Δ个点保留一个点，当点云的点数在100万以下时略过这一步，Δ采用以下式子计算：

步骤3.2：由步骤2.6中计算得到平均点云长度L_mean可以确定接下来要进行体素化下采样的体素格网大小δ₁,δ₂，将原始点云S⁰,T⁰以所在的空间细分为长度为δ₁的体素格网，体素内的所有点用这些点的重心替代，实现点云的压缩，得到压缩后的点云S₁和T₁。

上式中，Grid₁，为第1层金字塔点云格网的数目，取值在(200,600)之间，本实例取400，Grid₂为第2层金字塔点云格网的数目，取值在(40,120)之间，本实例取80。

步骤3.3：在压缩后的点云S₁和T₁上采用体素格网大小为δ₂的体素下采样得到进一步压缩的点云S₂和T₂，从而生成金字塔点云S₀、S₁、S₂和T₀、T₁、T₂结果如图4所示。

步骤4：计算点云S₂和T₂中所有点的法向量和所有点的FPFH(Fast Point FeatureHistogram)描述子；

步骤4.1：给定点云中任意一点P_i,找到距离点P_i的距离小于r_norm(r_norm＝2*δ₂)的所有邻域点集Q，假设Q中含有的点个数为k₁,质心为计算正定的方差协方差矩阵M：

步骤4.2：接着对矩阵M进行SVD分解，则最小特征值对应的特征向量即为点P_i点所在表面的法向量；

步骤4.3：按以上方法对点云S₂和T₂中每一个点进行法向量计算，针对由扫描仪所获取的点云，连接点云设备中心与每一个点，判断点云法向量与设备中心到该点的矢量方向的夹角，当该夹角小于90度时，将该点的法向量反向，从而实现法向量的定向。

步骤4.4：给定点云中任意一点P_i,利用Kd树搜索距离点P_i的距离小于r_FPFH(r_FPFH＝3*δ₂)的所有邻域点集合Q1(Q1不包括点P_i本身),设Q1中含有的点个数为k₁,连接Q1中所有的点P_ij和参考点P_i，每一对点定义一个局部坐标系UVW(选择两个点的法向方向和两点连接线之间的夹角较小的点对应得法向量作为U轴,该点用P_s表示而另一个点表示为P_t，其对应的法向量分别表示为U轴上的单位矢量为U轴和两点之间的连线的叉乘方向作为V轴，V轴上的单位矢量为U轴和V轴的矢量的叉乘作为W轴，W轴上的单位矢量为)。步骤4.5：由步骤4.4建立局部坐标系后，计算出用来表示两点之间法向量的差异特征的三个值α、θ(α和的取值范围为(-1,1),θ的取值范围为计算Q1内所有点与P_i两两法向差异特征，形成33维的特征描述子SPFH_i，其中前11维统计α的分布，中间11维统计的分布，剩下的11维统计θ的分布，

步骤4.6：计算完所有点的特征描述子SPFH_i，给定任意点P_i，找到距离点P_i的距离小于r_FPFH(r_FPFH＝3*δ₂)的所有邻域点集合Q1(Q1不包括点P_i本身),设Q1中含有的点个数为k₁，通过反距离加权的方式求所有邻域点的SPFH_i之和，得到P_i的未归一化的特征描述子FPFH′，将特征描述子FPFH′的每一维变量除以所有维度的值之和使得结果中所有维度的值之和等于1，得到该点的归一化特征描述子FPFH，按以上方法对点云S₂与T₂中所有点进行FPFH描述子计算。

上式中l为点P_ij到点P_i的距离，为权函数。

步骤5：双向匹配S₂和T₂中所有点的FPFH描述子，得到匹配对应关系C：

步骤5.1：给定点云S₂中任意一点P_i，其特征描述子为FPFH_i，利用KD-树查找在点云T₂中是否存在点P_j，它的特征描述子FPFH_j与FPFH_i的欧式距离D最小并且小于阈值ε₁，则以P_j为参考点，查找在点云S₂中是否存在点P_i′,它的特征描述子FPFH_i′与FPFH_j的欧式距离D最小并且小于阈值ε₁，当P_i′与P_i为同一个点时，认为点云S₂中的点P_i与点云T₂中点P_j匹配。

步骤5.2：遍历点云S₂中所有点，得到匹配对应关系C，利用直线连接C中的所有队员关系，如图4的(a)示意。

步骤6：运用SAC-IA(随机采用一致性)算法从C中将内点筛选出来，同时得到初始配准转换参数；

由步骤5已经计算出到匹配对应关系C，随机选取3对匹配点，一对匹配点由源点云S₂中一点和目标点云T₂中一点组成，3对匹配点含有源点云S₂中三点和和目标点云T₂中三点和给定阈值r_sac＝3*δ₂，源点云和目标点云中三点任意两点之间的距离必须大于阈值，否则进行下一次采样，接着由这三对对应点计算出转换参数Tr_rand，C中所有的匹配点应用Tr_rand统一到相同坐标系，计算匹配点对在相同坐标系下的点位误差RMSE，当点位误差小于阈值时，认为该点为内点，遍历所有匹配点对，得到内点个数，用内点个数来评价转换参数Tr_rand，随机采样足够多次之后，留下内点数最多的转换参数Tr_best作为初始转换参数Tr₂，以其对应的内点作为内点C₂。

T_目标＝R*S_待转+t

由SAC-IA算法剔除掉错误匹配的匹配对应关系见图4的(b)所示。

步骤7：在点云S₁与T₁上利用步骤8中得到的内点和初始转换参数，将内点邻域内的点云用来完成第一次ICP精配准优化，得到第一层优化的转换参数Tr₁；

步骤7.1：由步骤8得到初始转换参数Tr₂，和初始匹配对应关系C₂，选取点云S₁与T₁中属于C₂中对应点云的点的半径为r_sac的邻域S′₁、T′₁，

步骤7.2：对点云S′₁应用转换参数Tr₂得到点云S′_1Tr，点云S′_1Tr与T′₁应用ICP(Iterative Closest Point)进行第一次配准优化得到转换参数Tr₁和点位优化后的匹配对应点关系C₁。

步骤8：在点云S₀与T₀上利用步骤6中的内点C₂和步骤7中初步ICP优化后的转换参数，采用最近邻点代替原始对应点，仍然采用点邻域来完成第二次ICP精配准优化，得到最终配准结果。

步骤8.1：由步骤7得到转换参数Tr₁，和匹配对应关系C₂，选取点云S₀与T₀中属于中对应点云C₂的点的半径为δ₁的邻域S′₀、T′₀

步骤8.2：对点云S′₀应用转换参数Tr₁得到点云S′_0Tr，点云S′_0Tr与T′₀应用ICP(Iterative Closest Point)进行第二次次配准优化得到转换参数Tr₀，

步骤8.3：利用Tr₀对源点云进行计算，得到配准结果，最后配准结果示意图为图5中的(b)。

对本发明进行测试，可以看出本发明具有参数自动设置，配准速度快，配准精度高的特点，本发明中的平均点间距与平均点云长度确定相关参数、建立金字塔结构点云、所有点的FPFH特征描述子直接双向匹配和特征点引导的ICP精配准方法，确保了快速有效地配准类似地面激光扫描仪获取的海量的点云，具有一定实用价值。

Claims (9)

1.一种地面激光点云配准方法，其特征在于，包括以下步骤：

1)读入由点云获取设备采集的具有不同视角的源点云S₀和目标点云T₀两视角点云；

2)估计点云的平均点间距d_mean和平均长度L_mean；

3)构建源点云金字塔S₀、S₁、S₂与目标点云金字塔T₀、T₁、T₂；

4)计算点云S₂与T₂中所有点的法向量和所有点的FPFH描述子；

5)双向匹配S₂与T₂中所有点的FPFH描述子，得到匹配对应关系C；

6)利用SAC-IA方法从C中将内点C₂筛选出来，同时得到初始配准转换参数Tr₂；

7)在点云S₁与T₁上利用步骤6)中得到的内点C₂和初始转换参数Tr₂，提取C₂邻域内的点云，用于第一次ICP优化，得到第一层优化的转换参数Tr₁；

8)在点云S₀与T₀上利用内点C₂和转换参数Tr₁，采用最近邻点代替原始对应点，采用点邻域完成第二次ICP优化，得到最终配准结果。

2.根据权利要求1所述的地面激光点云配准方法，其特征在于，点云的平均点间距d计算公式为：

其中，k₁为常数；K₀为采样点个数；d_ij为邻域中每一个点P_j到点P_i的距离。

3.根据权利要求2所述的地面激光点云配准方法，其特征在于，K₀＝1000；k₁＝4。

4.根据权利要求1所述的地面激光点云配准方法，其特征在于，利用平均点云长度L_mean确定体素格网大小δ₁,δ₂，将原始点云S₀，T₀所在的空间细分为长度为δ₁的体素格网，体素内的所有点用这些点的重心替代，实现点云的压缩，得到压缩后的点云S₁和T₁；

其中，Grid₁为第1层金字塔点云格网的数目，Grid₂为第2层金字塔点云格网的数目；L_mean为平均点云长度。

5.根据权利要求4所述的地面激光点云配准方法，其特征在于，Grid₁取值为(200,600)；Grid₂取值为(40,120)。

6.根据权利要求4所述的地面激光点云配准方法，其特征在于，步骤4)的具体实现过程包括：

1)给定点云中任意一点P_i,找到距离点P_i的距离小于r_norm的所有邻域点集Q，设Q中含有的点个数为n₁,质心为计算正定的方差协方差矩阵M：其中P_ij为Q的第j点；

2)对矩阵M进行SVD分解，则最小特征值对应的特征向量即为点P_i所在表面的法向量；

3)按以上步骤1)和步骤2)对点云S₂和T₂中每一个点进行法向量计算，连接点云设备中心与每一个点，判断点云法向量与设备中心到该点的矢量方向的夹角，当该夹角小于90度时，将该点的法向量反向，从而实现法向量的定向；

4)利用Kd树搜索距离点P_i的距离小于r_FPFH的所有邻域点集合Q1,设Q1中含有的点个数为k₂，连接Q1中所有的点P’_ij和参考点P_i，每一对点定义一个局部坐标系UVW；

5)在局部坐标系UVW下计算出用来表示两点之间法向量的差异特征的三个值α、θ，计算Q1内所有点与P_i两两法向差异特征，形成33维的特征描述子SPFH_i，其中前11维统计α的分布，中间11维统计的分布，剩下的11维统计θ的分布；

6)通过反距离加权的方式求所有邻域点的SPFH_i之和，得到P_i的未归一化的特征描述子FPFH′，将特征描述子FPFH′的每一维特征除以所有特征之和，得到点P_i的归一化特征描述子FPFH，重复步骤4)～6)，对点云S₂与T₂中所有点进行FPFH描述子计算。

7.根据权利要求4所述的地面激光点云配准方法，其特征在于，步骤5)的具体实现过程包括：

1)给定点云S₂中任意一点P_i，其特征描述子为FPFH_i，利用KD-树查找在点云T₂中是否存在点P_j，它的特征描述子FPFH_j与FPFH_i的欧式距离D最小并且小于阈值ε₁，则以P_j为参考点，查找在点云S₂中是否存在点P_i′,它的特征描述子FPFH_i′与FPFH_j的欧式距离D最小并且小于阈值ε₁，当P_i′与P_i为同一个点时，认为点云S₂中的点P_i与点云T₂中点P_j匹配；

2)遍历点云S₂中所有点，得到匹配对应关系C。

8.根据权利要求4所述的地面激光点云配准方法，其特征在于，步骤7)的具体实现过程包括：

1)由初始转换参数Tr₂和初始匹配对应关系C₂，选取点云S₁与T₁中属于C₂中对应点云的点的半径为r_sac的邻域S′₁、T′₁；

2)对点云S′₁应用转换参数Tr₂得到点云S′_1Tr，点云S′_1Tr与T′₁应用ICP进行第一次配准优化得到转换参数Tr₁和点位优化后的匹配对应点关系C₁。

9.根据权利要求4所述的地面激光点云配准方法，其特征在于，步骤8)的具体实现过程包括：

1)由转换参数Tr₁和匹配对应关系C₂，选取点云S₀与T₀中属于中对应点云C₂的点的半径为δ₁的邻域S′₀、T′₀；

2)对点云S′₀应用转换参数Tr₁得到点云S′_0Tr，点云S′_0Tr与T′₀应用ICP进行第二次次配准优化得到转换参数Tr₀；

3)利用Tr₀对源点云进行计算，得到配准结果。