一种三维点云数据自动配准方法与流程

本发明属于三维测量领域，更具体地，涉及一种三维点云数据自动配准方法。

背景技术：

三维模型在工业检测、文物保护、生物医学等领域都具有重要的意义。随着三维测量技术的发展，对点云模型的处理技术已成为近年来研究的热点。点云处理中一个重要步骤就是将同一场景中不同视点获得的点云数据统一到同一坐标系，即点云配准。

点云配准方法主要分为以下几类：(1)手动配准：手动选取两片点云之间的对应点实现拼合，该方法需要人工干预，拼合的效果往往依赖操作者的经验；(2)基于外部辅助的配准方法：利用转台或者在被测物上粘贴标志点的方式实现辅助配准，此类方法对环境要求较高，配准的精度易受外部辅助的影响；(3)基于自身形貌特征的配准方法：利用被测物体自身的形貌特征实现点云之间的匹配进而实现点云的自动配准，无需使用额外的设备或者标志点，能够满足不同测量场合的需求，但算法的稳定性和精度与有辅助的点云配准技术相比还有待加强。

针对基于自身形貌特征的配准方法，已有研究单位进行相关研究，并取得了一定的成果；Rusu利用特征点与周围邻域点的法向量夹角作为特征构建了一种快速点特征直方图(Fast Point Feature Histograms，FPFH)用于特征点的查找和匹配，完成初始配准得到较理想的初始位置之后，采用精确配准算法提高配准精度；Chen和Besl提出了迭代最近点算法(Iterative closest point,ICP)来进行两幅点云数据的精确配准，该方法通过迭代寻找对应点集中的最近点作为对应点，不断优化刚性变换矩阵，最终获得精确配准。但该方法对两片点云的初始位置要求较高，需要点云之间初始相对位置大致相近。

现有的点云配准方法通常采用粗配准和精确配准相结合的方式，其中，粗配准一般通过计算点的特征寻找对应点，估计得到两幅点云间的初始位置关系；精配准则是采用ICP及其改进算法对粗配准结果进行进一步优化，实现点云数据的精确配准；但现有的利用点的特征进行点云配准的方法对于点云噪声、体外孤立点以及点云密度变化等干扰因素较为敏感，易出现大量的误匹配点，降低了配准的精度和稳定性。

技术实现要素：

本发明的目的在于针对现有的技术问题，提供一种三维点云数据自动配准方法，配准精度高并且具有很强的稳定性，完成配准结果正确性的自动判断。

为实现上述目的，本发明的技术解决方案是基于旋转不变特征描述子的点云自动配准方法，首先对待配准的两片点云进行采样得到特征点，分别计算这些特征点的旋转不变特征子，对两片点云中特征点的旋转不变特征进行匹配搜索，得到特征点之间的初始对应关系；然后采用随机采样一致算法(Random Sample Consensus，RANSAC)对初匹配点集中存在的误匹配点进行判断和去除，得到优化的特征点对应关系，计算得到两片点云之间的大致刚性变换关系，实现粗配准；并提供一种刚性变换一致性的检测算法，利用匹配特征点间的局部坐标系变换关系对粗配准结果进行约束性检测，进行粗配准结果正确性的验证；最后采用迭代最近点算法优化点云数据间的刚性变换关系，最终实现点云的自动精确配准。

具体包括如下步骤：

步骤1：读入采集到的待配准的源点云P和目标点云Q；

步骤2：分别计算源点云P和目标点云Q的密度，并从所述源点云P中随机选取若干个点构成源特征点集S₁，从目标点云Q中随机选取若干个点构成目标特征点集S₂；

步骤3：根据源点云与目标点云的密度，分别计算源特征点集S₁与目标特征点集S₂中每个特征点的局部旋转平移不变坐标系；

步骤4：根据特征点的局部旋转平移不变坐标系计算每个特征点的高维特征描述；对目标点云Q与源点云P进行特征点匹配，获得初匹配点集C；

步骤5：利用随机采样一致算法去除所述步骤4中获得的所述初匹配点集中的误匹配，利用基于奇异值分解法的刚性变化估计算法计算获得旋转矩阵R和平移矩阵T，获得所述目标点云和所述源点云之间的粗配准结果；

步骤6：采用基于刚性变换一致性的配准错误检测算法判定步骤5获得的粗配准结果是否正确，结果正确则转入步骤7，结果不正确则返回配准失败的结果；

步骤7：采用迭代最近点算法迭代优化上述估计得到的旋转矩阵R和平移矩阵T，实现所述目标点云和源点云的自动精确配准。

总体而言，通过本发明所构思的以上技术方案与现有技术相比，能够取得下列有益效果：

(1)本发明提供的三维点云数据自动配准方法，通过从点云数据中选取特征点集，为每个特征点构建局部旋转平移不变坐标系，并计算得到在该特征点上的高维特征描述，实现了对传统人工标志点的替代，可以用于诸如文物数字化扫描等应用场景当中，可以减少人工粘贴、去除标志点的工作量；

(2)本发明提供的三维点云数据自动配准方法，利用所得到的特征点集和高维特征描述，可以有效实现不同点云数据之间的对应点查找与匹配，实现任意初始姿态下点云数据的自动配准；

(3)本发明提供的三维点云数据自动配准方法，通过将利用局部旋转不变坐标系计算得到的局部旋转矩阵结果与利用随机抽样一致性算法求解得到的全局旋转平移矩阵结果进行比较，判断其之间是否具有旋转平移一致性，可以进而有效判断粗配准结果是否成功，实现了精确配准。

附图说明

图1是按照本发明提供的三维点云数据自动配准方法的整体流程示意图。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅用以解释本发明，并不用于限定本发明。此外，下面所描述的本发明各个实施方式中所涉及到的技术特征只要彼此之间未构成冲突就可以相互组合。

本发明实施例所提供的三维点云数据自动配准方法，其流程如图1所示，包括如下步骤：

步骤1：读入由三维测量设备采集得到的目标点云和源点云；

步骤2：分别随机选取源点云和目标点云中的若干个点构成两组特征点集，具体为：

根据预设的特征点采样比例h₁和h₂，分别从源点云P与目标点云Q中随机采样得到对应的特征点集

其中，n₁是指源点云P的特征点集的特征点数目，n₂是指目标点云Q的特征点集的特征点数目，n₁＝h₁·N₁，n₂＝h₂·N₂。

步骤3：对采样得到的特征点集，计算其中每个特征点的局部旋转平移不变坐标系；具体包括以下子步骤：

步骤3.1：对于源点云中的每个点p_i与目标点云Q中的每个点q_i，利用最近邻搜索方法计算点p_i的最近邻点p_i′，以及点q_i的最近邻点q_i′

获取点云p_i的点云密度d_p，以及点云q_i的点云密度d_q；

其中，||p_i-p_i′||是指点p_i与p_i′之间的欧式距离，||q_i-q_i′||是指之间点q_i与点q_i′之间的欧式距离；i＝1,2,3,…,N，N为点云数据中点的个数；源点云P中包含N₁个点，对源点云P而言，N＝N₁；目标点云中Q包含N₂个点，对目标点云Q而言，N＝N₂；

步骤3.2：任意选取特征点集中一点p_i，构造点p_i的协方差矩阵：

其中，w_i＝1/||p_j-p_i||,半径r_loc为p_i的邻域计算半径；实施例中，取r_loc＝15d_p；

步骤3.3：采用下式求解协方差矩阵的特征值和特征向量；

COV(p_i)V＝EV

其中，E是由构成的3×3对角阵

是特征值对应的特征向量，l＝1,2,3；

步骤3.4：将分别设为x⁺,y⁺,z⁺，的反向向量设为x^-,y^-,z^-；

局部坐标系中z轴的方向由下式确定：

步骤3.5：采用步骤3.4的方法确定局部坐标系中x轴方向；局部坐标系中y轴方向由z×x确定，由此建立起以p_i为原点的局部旋转平移不变坐标系F_i＝{x,y,x}，；采用与上述步骤3.2～3.4同样的方法，建立以点q_i为原点的局部旋转平移不变坐标系。

步骤4：根据计算得到的特征点的局部旋转平移不变坐标系，计算特征点的高维特征描述，对两幅点云进行特征点匹配，得到初匹配点集；具体包括如下子步骤：

步骤4.1：对源点云特征点集中任一特征点p_i，其高维特征描述为高维向量形式在目标点云特征点集所对应的特征几何所构成的高维空间中搜索与距离最近的特征及次近的特征

在本发明中，最近邻特征的搜索是通过快速近似最近邻搜索方法(Fast Approximate Nearest Neighbor Search)进行加速计算；

步骤4.2：计算特征分别到最近的特征次近的特征的欧式距离和计算与比值，并判断该比值与预设的距离比阈值τ的大小；

具体地，根据下式判定得到与之间是否具有正确对应关系

若则表明与匹配成功，否则匹配失败；若与匹配成功，则其分别对应的源点云与目标点云中的特征点也正确对应，该点对构成一个匹配点对；

对源点云的其他特征点也进行上述处理，由所有的匹配点对构成初匹配点集C。

步骤5：采用随机采样一致算法(Random Sample Consensus，RANSAC)去除初匹配点集中的误匹配，利用奇异值分解法(SVD)计算旋转矩阵R和平移矩阵T，得到源点云和目标点云之间的粗配准关系，具体包括如下子步骤：

步骤5.1：设定随机采样的次数S_num和初始残余误差E_rr，从初匹配点集C中随机采样，选取3对匹配点作为初始点；

若选取的某个特征点的匹配点有多个，则从相应点集中随机选取其中一点作为该查找点的对应点；

步骤5.2：利用SVD算法求解旋转矩阵R和平移矩阵T，具体如下：

步骤5.21：设由步骤5.1得到匹配点集为P＝{p₁,p₂,p₃,},P′＝{p′₁,p′₂，p′₃}，利用下式计算点集的质心；

其中，p_k和p′_k为任一对匹配点的三维坐标；

步骤5.22：利用下式将点集P和P′相对于各自质心做平移，得到新点集Q＝{q₁,q₂,q₃,}和Q′＝{q′₁,q′₂,q′₃}；

q_k＝p_k-g,q′_k＝p′_k-g′,(k＝1,2,3)；

步骤5.23：利用下式计算3×3的矩阵M：

步骤5.24：对M矩阵进行奇异值分解M＝UΛV^T；

其中，上标T为矩阵的转置，U、V为3×3的酉矩阵，Λ为3×3的对角阵，定义3×3的对角阵A为：

步骤5.25：利用下式计算3×3的旋转矩阵R和3×1的平移矩阵T：

R＝UAV^T，T＝g′-Rg；

步骤5.3：对初匹配点集中所有的匹配点对，根据下式计算旋转平移后的距离误差d_err；

d_err＝||p_τ-(R·q_τ+T)||₂；

若距离误差d_err小于给定的初始残余误差E_rr，则判定该匹配点是内点；按照上述方法求出所有内点并统计内点的个数m′；

实施例中，残余误差即为距离阈值；(p_τ,q_τ)为初匹配点集中的一个匹配点对，τ＝1,2,…r，r为初匹配点集中的匹配点对的数目；

步骤5.4：如果内点个数m′大于给定的内点个数m，则利用该内点数据集中的点采用SVD算法重新求解旋转矩阵R、平移矩阵T；

根据下式计算残余误差E′_rr；

若残余误差小于初始残余误差，则将此时计算得到的旋转平移矩阵作为最佳估计目标模型参数，更新旋转平移矩阵R、T和初始残余误差E′_rr；否则转到步骤5.5；若内点的个数小于给定内点的个数，则直接转到步骤5.5；

步骤5.5：重复随机采样S_num次，重复步骤5.2～5.4，得到该组抽样所对应的内点数，对所有抽样的内点数进行排序，选取内点数最多的抽样结果作为最佳抽样；

利用该最佳抽样下得到的内点数据集C，根据步骤5.2求解旋转矩阵R、平移矩阵T，作为最佳旋转平移矩阵R_ran和T_ran；其中C₁中包含S对匹配点对。

步骤6：采用刚性变换一致性的检测算法校验粗配准关系是否正确；

步骤6.1：利用下式计算初匹配点集C中的匹配点对(p_loc,q_loc)其匹配点之间的旋转平移矩阵(R_loc,T_loc)：

其中，R_loc为两个匹配点之间的局部旋转矩阵，T_loc为两个匹配点之间的局部平移矩阵；

F_ploc、F_qloc分别是在点p_loc,q_loc建立的局部坐标系，loc＝1,2,…S，S为匹配点集C中匹配点对的数目；

步骤6.2：将求解得到的旋转矩阵R_loc转换为欧拉角表示，即

R_loc→(α_loc,β_loc,γ_loc)；将步骤5.6中得到的旋转矩阵R_ran转化为欧拉角表示，即R_ran→(α_ran,β_ran,γ_ran)；

步骤6.3：根据下式计算欧拉角之间的角度差值d_a：

其中，Δ(η₁,η₂)²＝(η₁-η₂)²,η₁,η₂为欧拉角；

步骤6.4：根据下式计算与T_ran＝(t_x,t_y,t_z)^T之间的距离差值d_t：

步骤6.5：根据d_a和d_t是否均小于给定的阈值σ_a和σ_t，判定两种方法求解出的旋转平移矩阵是否一致；

实施例中，σ_a取0.5236，σ_t取15D_den，D_den为点云密度，即旋转矩阵相差小于30°，平移矩阵距离小于15倍的点云密度时，则表明本步骤所计算的旋转平移矩阵与步骤5所获得的旋转平移矩阵是一致的；

步骤6.6：对最佳抽样下得到的内点数据集C内的其他匹配点对也执行步骤6.1到步骤6.5，得到所有匹配点对之间的局部旋转平移矩阵与上述最佳旋转平移矩阵(R_ran,T_ran)之间的一致性关系；

最终统计得到最佳抽样下得到的内点数据集C中所有满足一致性关系的匹配点对的数目s；

步骤6.7：计算满足一致性关系的匹配点对的数目s与匹配点集C中匹配点对的数目S的一致性比值λ＝s/S；

若λ≥阈值τ_λ，则表明利用RANSAC求解得到的旋转矩阵和利用局部旋转不变坐标系计算得到的旋转矩阵是一致的，则判定为配准成功；否则，判定配准失败；实施例中，阈值τ_λ＝0.7。

步骤7：采用改进的ICP算法优化点云之间的刚性变换关系，实现点云的自动精确配准；刚性变换关系包括旋转平移矩阵；具体步骤如下：

步骤7.1：设定距离阈值ω作为迭代终止的条件；其中，ω>0；距离阈值ω根据源点云P的点云密度d_p；实施例中，取ω＝5d_p；

步骤7.2：在源点云中随机选取若干点作为待匹配点；

步骤7.3：用逆向投影法在目标点云中查找源点云中待匹配点的对应点；

步骤7.4：采用基于点到面距离度量作为ICP算法所需求解的目标函数，不断迭代计算源点云到目标点云的刚性变换关系；

步骤7.5：当目标函数值小于距离阈值ω时，停止迭代；并将此时求解得到的刚性变换关系作为最终结果，完成点云的自动精确配准。

本领域的技术人员容易理解，以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。