1 引言

近些年来,三维(3D)点云数据被广泛应用于三维模型重建、文化遗产管理、机器人导航定位等领域^[1]。快速、准确的三维点云配准是关键技术和研究重点。三维点云配准的目的是找到能使两个输入点云对准一个共同坐标系的最佳刚体变换。在实际应用中,数据可能会受到严重的遮挡,两个点云之间重叠区域很小,这使得寻找最佳刚体变换的过程充满挑战。因此,针对重叠范围更小的点云,寻找快速、准确、鲁棒的配准算法是当今一个活跃的研究主题。

目前,研究最热的三维点云配准算法可以分为两类:基于特征的配准方法和基于无特征的配准方法。前者是对关键点进行提取^[2-3]并利用能够在刚体变换下保持不变的特征描述子^[4-6],找到两片点云的对应关系,进而进行配准。但是当重叠率较低时,此类方法提取特征的质量和速度会受到限制。后者是在原始点云的基础上直接进行配准。此类方法依靠底层采样策略识别无异常子集,从而建立两个子集的对应关系。随机采样一致性(RANSAC)算法在存在噪声和离群值的情况下,很难选出迭代所需要的三对优质对应点作为样本子集,配准效果会受到很大影响。此方法可以使用局部优化技术来改进^[7-9]。四点同余集(4PCS)算法采取的是枚举4个共面点所有一致子集的策略,计算代价比较高。Super4PCS算法对其进行改进,引入智能索引,使得其复杂度从二次降为线性。Le等^[10]通过采用图匹配问题的特殊实例,提出了新的采样策略。通过解决松弛凸问题的方式快速搜索对应关系,用于估计和验证假设,该方法在运行时间和处理含大量噪声和离群值的情况上有着很大的优势。

依据文献[ 11-13],当源点云与目标点云的重叠率低于60%时,两片点云的重叠程度低。对于低重叠率的点云数据,由于它们的重叠部分点云数据量很少,可以提取的特征有限,而直接利用原始点云配准会产生很多的离群子集,处理起来耗时而且容易匹配错误。因此,低重叠率的点云配准问题是配准中的难点。近年来,一些研究者将上述两大类方法各自的优势结合在一起,提出了综合性的方法,这为低重叠率点云的配准研究提供了新的方向。

文献[ 14]将多尺度稀疏特征(MSSF)嵌入到4PCS中;文献[ 15]先建立点邻域关系,再利用四点法拼接求取转换参数,它们均在一定程度上实现了低重叠率点云的高效全局配准。文献[ 11,16]均提出了基于关键点的Surper4PCS算法,文献[ 11]利用内在形状签名(ISS)提取关键点,而文献[ 16]利用三维尺度不变特征变换(3D-SIFT)算法提取关键点,两种方法均缩小了Super4PCS算法的检索范围,解决了其对于低重叠率点云配准耗时而且精度差的问题。

为了进一步提高低重叠率三维点云的配准精度和时间效率,本文提出一种将聚类区域分块和凸优化问题相结合的点云配准方法。对源点云和目标点云建立多尺度特征描述符,利用描述符的差异进行聚类分块,获取重叠区域以及它们的对应关系,最后利用半定规划(SDP)凸求解问题和迭代最近点(ICP)算法进行点云配准。以期能够获取到初始点云数据中非冗余的并且完整的点云数据,在不建立二次描述符、节省计算量的情况下获取两片点云的重叠区域,将基于对应关系和聚类分块的配准方法应用于最新的基于无特征的图匹配问题,很好地去除离群值和优化对应关系,提高点云配准尤其是低重叠率点云配准的精度和效率。

2 获取重叠区域

由于低重叠率的点云数据重叠区域很小,难以提取特征,容易出现错配的问题。因此本文方法将源点云与目标点云的重叠区域提取之后,再进行配准,由此来放大重叠区域在全局配准中的作用。

与文献[ 17-19]获取重叠区域的方式不同,本文方法在利用曲率特征建立多尺度描述符的基础上,直接构造该描述符的角度差异,衡量区域之间的相似程度,减少了重新建立描述符的计算量;而且本文方法在进行聚类分块的同时保留了点之间的对应关系,为下一步配准奠定了基础。

2.1 构建多尺度描述符

对于点云数据来说,曲率特征和法向量是最基本的特征,可以很好地反映点云的凹凸程度,利用它们可以在简化点云数据的基础上保证其几何特性的完整性。因此,本文借鉴文献[ 20]来构建多尺度描述符。构建方式如下。

1) 通过多次改变邻域半径r的大小来实现多尺度描述符。假设要考虑L个不同的邻域半径,它们是r₁<r₂<…<r_L。对于每个查询点I,在每个邻域半径r_l内构建协方差矩阵表示为

Cl=1Sl∑Ii∈Sl(Ii-I)(Ii-I)Τ,(1)

式中:l=1,…,L;S_l= Ii|Ii-I≤rl。

2) 利用奇异值分解(SVD)算法分解(1)式,可以获得三个特征值λ_l1≥λ_l2≥λ_l3及其对应的特征向量n_l1、n_l2、n_l3。其中,最小的特征值所对应的特征向量n_l3即为平面的法向量,记为n_l。为避免后续法向量的方向影响区域分块,将区域法向量的方向设置为一致,将不指向视点方向的法向量进行取反处理。

3) 将特征值进行归一化,得到向量d_l,表示为

dl=(λl1λl2λl3)Τλl1+λl2+λl3。(2)

为了充分利用特征值的差异来表现点云块的变化,定义Δd_l=d_l+1-d_l,将要考虑的L个不同邻域半径的Δd_l连接起来构建基于特征值的描述符D,表示为

D=Δd1⋮ΔdL-1。(3)

合并L个不同邻域半径的法向量n_l得到基于法线的描述符N,表示为

N=(n1,…,nL)。(4)

将(3)式和(4)式综合起来得到基于特征值和法线的多尺度描述符(N,D)。

2.2 聚类区域分块

建立多尺度特征描述符之后,就可以依据其对源点云与目标点云建立对应关系。在源点云与目标点云重叠区域小的情况下,利用多尺度描述符的差异对有待配准的点云在建立对应关系的同时进行聚类分块,获取源点云与目标点云的重叠区域,达到充分利用重叠部分作用的目的。文献[ 17-19]对于低重叠率点云的配准问题也采用了区域分块的策略。文献[ 17-18]均是利用曲率差值和法线差值建立区域分块之后,再建立快速点特征直方图(FPFH)描述符,并利用FPFH的方差来衡量点云块之间的相似程度。文献[ 19]将特征点的区域分割问题转换为谱空间的聚类问题并利用SVD进行区域配准。为减少计算量,本文直接利用上一节得到的多尺度描述符建立角度差异来衡量区域之间的相似程度,进行聚类分块^[21]。聚类分块的流程如下。

1) 依据上述得到的多尺度描述符对下采样之后的源点云P与目标点云Q执行最近邻搜索,形成种子匹配,并对它们进行排序,记为(p_i,q₁),…,(p_j,q_n),其中p_i代表源点云中的种子点,q_n代表目标点云中的种子点。

2) 对于第一个种子匹配(p_i,q₁),通过以下方式形成其聚类点云分块。

①定义近邻搜索范围

对于近邻搜索,将搜索范围定义为ε₁。即当源点云中的点p在搜索范围α= p-pi时,目标点云中的点q的搜索范围为(β= q-q1)∩( β-α<ε₁)。

②定义角度差异ε₂

对于源点云P中近邻搜索范围的任意点p,用向量θ表示其描述符与点p_i的描述符之间的角度,θ=(θ₁,θ₂,…,θ_L)^Τ,L表示要考虑到的多尺度邻域半径个数,对于其中每个邻域半径r_l,

θl=arccos<nlpi,nlp>。(5)

同理,对于目标点云Q中近邻搜索范围的任意点q,用向量ψ表示其描述符与点q₁的描述符之间的角度,ψ=(ψ₁,ψ₂,…,ψ_L)^Τ,对于其中每个邻域半径r_l,

ψl=arccos<nlq1,nlq>。(6)

将点p和点q之间的角度差异定义为

Δθ=1L∑l=1Lθl-ψl,allθl-ψl<ε2inf,otherwise。(7)

求出最小的Δθ对应的点p和点q。

③判断点p和点q是否满足条件 Dp-Dq<ε₃,ε₃为两点之间基于描述符D的距离差异的阈值。

将同时满足以上条件的对应关系(p,q)全部聚合为第一个对应关系点云分块M₁。

3)对于后续的每一个种子匹配(p_i,q₂),…,(p_j,q_n),都执行步骤2,最终得到n个对应关系点云分块,即M₁,M₂,…,M_n。

至此,通过对应关系聚类区域分块,可以得到源点云与目标点云中相似的点云块,以及它们之间的对应关系。

3 改进的三维点云配准算法

3.1 算法设计

基于凸优化问题的配准算法是Le等^[7]于2017年提出的新方法,该方法来源于图匹配,作者将其命名为SDRSAC(Semidefinite-Based Randomized Approach for Robust Point Cloud Registration without Correspondences)。它与4PCS算法一样,是一种基于无特征的随机点云配准算法。但其又与4PCS算法以及其同系列的改进方法不同,它采用的是SDP凸求解问题搜索对应关系而不是枚举的方式,大大提高了搜索的时间效率,增强了处理离群值的能力。本文对SDRSAC算法进行了改进,将其用于基于特征和对应关系的点云配准,同时,达到了减小配准搜索范围、节省计算量、提高配准精度和时间效率的效果。

对于SDRSAC算法,作者在每次迭代时的采样策略是从源点云P中随机选取N_sample个采样点组成点集P_N,再从目标点云中随机选取N_sample个采样点组成点集Q_N,然后对于点集P_N和点集Q_N进行匹配潜力的判断,去除离群值,进而求解对应的最优关系集,通过不断迭代选择最终评分中最高的、最适用的最优对应关系集,求解最佳变换并完成最后的配准。

本文对其进行改进,将从上一节得到的对应关系点云分块M₁,M₂,…,M_n在源点云P中的点云集记为P_c,在目标点云Q中的点云集记为Q_c,把它们以及它们的对应关系代入SDRSAC算法中,利用其求得的最优变换进行粗配准,并利用ICP算法进行细配准。与原始点云P和Q相比,P_c和Q_c数据量虽然小,但是重叠率较高,而且包含高质量的对应关系,为配准工作减少了很大的计算量,并提高了精确度。图1给出了本文算法的主要流程图。

图 1. 本文算法主要流程图

Fig. 1. Main flow chart of proposed algorithm

下载图片 查看所有图片

3.2 算法实现

配准算法的实现步骤如下。

1) 选择每次迭代时的采样点。

每次迭代时,从P_c中随机选取N_sample个采样点组成点集P_c1,利用上一节得到的对应关系聚合,得到目标点云Q_c中与P_c1相对应的点集Q_c1。

2) 计算点集的匹配潜力,去除离群值。

为了表述方便,A=P_c1,B=Q_c1,利用下列函数求解点集的匹配潜力,匹配潜力由H表示,H是对称矩阵。

Hab,cd=f(Aa,Ac,Bb,Bd),ifδ(Aa,Ac)-δ(Bb,Bd)≤γ0,otherwise,(8)

式中:a、b、c、d表示点的索引;ab和cd表示H的行和列的索引;A_a与A_c属于源点集A,B_b与B_d属于目标点集B,且(A_a,B_b)和(A_c,B_d)是由上一节得到的对应关系;γ>0是预定义的阈值;δ(A_a,A_c)表示两个三维点云之间的欧氏距离;f=exp -δ(Aa,Ac)-δ(Bb,Bd)用于评判两个线段的长度差。通过求解,如果A_aA_c和B_bB_d的长度相差小于γ,继续将(A_a,B_b)和(A_c,B_d)视为匹配候选,否则去掉它们的对应关系。

3) 求解优化的对应关系集。

为了求解优化的对应关系集 CAB*,假设置换矩阵X∈{0,1}^N×N,如果p_i∈A和q_i∈B属于对应集,X中第i行j列的元素的值为1,否则为0。

假设最优解包含m对的对应关系,而且m<N,将X堆叠起来得到X₁,X₁={ XΤi,1Xi,2Τ… XΤi,N}^Τ。

搜索最佳对应关系的问题可以转化为

maxXX1ΤHX1。(9)

同时满足以下条件: Xi,j∈{0,1} ∀i,j∈{1,…,N};∑j=1NXi,j≤1 ∀i∈{1,…,N};∑i=1NXi,j≤1 ∀j∈{1,…,N};∑i=1N∑j=1NXi,j=m。

令Y=X₁X1Τ∈ RN2×N2,则(9)式可以转化为

maxX,Ytr(HY)。(10)

(10)式满足以下条件,tr(Y)=m。

将上述非凸问题凸优化,得到

maxX,Ytr(HY),(11)

式中:tr(Y)=m;Y-X₁X1Τ≥0。

同时,在(11)式的基础上加上限制条件,确保源点云A中的一个点p_i只能与目标点云中的0个或者1个点q_i相匹配,而不能出现一个以上相匹配的情况;与此同时,要满足步骤2中的匹配潜力的要求,如果A_aA_c和B_bB_d的长度相差大于定义的阈值γ,则不允许其进行匹配。

通过SDP凸求解器以及从步骤2中求得的矩阵H来求解上述加上限制的(11)式,可以求得X₁。接下来利用线性分配问题(LP)将X₁投影到置换矩阵的空间X,得到优化的对应关系集 CAB*。

4) 迭代选取最优的对应关系集完成配准。

上述步骤3只是一次迭代中求得的最优的对应关系集,利用ICP将其进行细化。重复进行步骤1、2和3的迭代过程,求解并进行比较,得出最终评分中最高的、最适用的最优对应关系集,求得最佳变换(R^*,t^*),并进行配准。

4 实验结果与分析

4.1 实验环境及数据

实验数据来源于斯坦福大学公布的三维点云模型和geometryhub网站的实体三维扫描数据。实验环境为Windows10的64位操作系统 ,MATLAB 2018b开发平台。

依据文献[ 11-13],当源点云与目标点云的重叠率低于60%时,两片点云的重叠程度低,当它们的重叠率高于60%时,重叠程度较高。为了证明本文方法的有效性,实验分为重叠率高和重叠率低的两种情况进行。

对于点云模型,使用不同视角下获取的点云数据来体现它们的重叠率。实验选取有代表性的几组数据进行结果展示。重叠率较高的点云数据:Bunny000、Bunny045和Dragon000、Dragon024,它们存在明显的重叠,重叠率约为75%。重叠率较低的点云数据:Bunny045、Bunny090和Dragon000、Dragon048,它们的重叠率较低,重叠率约为50%;Happy48 和 Happy96的重叠率约为40%。除此之外,为了更好地体现本文实验方法的普适性,特意选取了从geometryhub网站获得的两组实体三维扫描数据,它们分别是Chair1、Chair2和Motobike1、Motobike2,重叠率约为40%。

4.2 点云配准评价标准

本文采用最大公共点集(LCP)作为点云配准精度的评价指标。LCP与刚体变换T相关,表示对于刚体变换之后,源点云中的点能够在目标点云中找到与之距离在一定容错范围内的点集。LCP的计算式为m/n,其中n表示目标点云中点的数量,m表示在目标点云中能够找到的与源点云中点相对应的点的数量。LCP的值越高,表示源点云与目标点云最佳对齐的程度越高,配准的精度越高。

4.3 配准结果分析

对选取的几组点云数据,按照本文所提的将聚类区域分块和凸优化问题相结合的点云配准方法进行配准。即对下采样之后的源点云和目标点云先构建多尺度协方差矩阵求解特征值和法向量,再构建多尺度特征描述符,依据描述符的差异进行对应关系聚类分块,而后代入SDP凸求解问题中求解最佳对应集,并进行最佳变换,最后利用ICP细化配准。

为了展现本文方法的优越性,本文设计了对比实验,与经典的Go-ICP(Globally optimal ICP)、Super4PCS算法以及最新的配准算法^[11、21]进行了对比。其中,文献[ 11]是一种基于关键点的融合Super4PCS和ICP的点云配准算法。文献[ 21]是一种基于快速描述符和对应传播的全局点云配准算法,本文提出的算法在文献[ 21]算法的基础上进行了优化改进。

4.3.1 重叠率高的点云配准结果分析

对于重叠率较高的两组点云数据Bunny000、Bunny045和Dragon000、Dragon024,它们的配准结果分别如图2、3所示。从图2、3中可以直观地看出,与最新的文献[ 11]和文献[ 21]算法相比,本文方法得到的配准结果中两种颜色的贴合程度更优,效果更好。

图 2. Bunny000与Bunny045配准结果。(a)初始位置;(b)利用文献[ 11]方法配准;(c)利用文献[ 21]方法配准;(d)利用本文方法配准

Fig. 2. Registration results of Bunny000 and Bunny045. (a) Initial position; (b) registering using method of Ref. [11]; (c) registering using method of Ref. [21]; (d) registering using proposed method

下载图片 查看所有图片

图 3. Dragon000与Dragon024配准结果。(a)初始位置;(b)利用文献[ 11]方法配准;(c)利用文献[ 21]方法配准;(d)利用本文方法配准

Fig. 3. Registration results of Dragon000 and Dragon024. (a) Initial position; (b) registering using method of Ref. [11]; (c) registering using method of Ref. [21]; (d) registering using proposed method

下载图片 查看所有图片

为了更准确地评价本文配准方法的准确性,使用4.2节的评价标准计算LCP,并记录配准算法的运行时间,与最新的文献[ 11]、[21]算法以及经典的配准算法Go-ICP和Super4PCS进行对比,对比结果如表1所示。

依据4.2节的点云配准评价标准,LCP的值越高,点云配准的精度越高。从表1可以看出,对于重叠程度较高的两组点云数据,本文方法的精度值均高于其他四种方法,配准时间较Go-ICP、Super4PCS和文献[ 21]三种方法更快。对于Bunny000与Bunny045这组数据,文献[ 11]的配准方法运行时间更具优势,原因是文献[ 11]方法在提取特征点时速度比本文的多尺度描述符快很多,而且特征点的提取缩短了后续Super4PCS算法的检索时间。

4.3.2 重叠率低的点云配准结果分析

对重叠率低的三组斯坦福公共点云数据进行配准的结果如图4~6所示。从图中可以直观地看出,本文方法得到的配准结果中两片点云的错位程度最低,配准效果最好。

图 4. Bunny045与Bunny090配准结果。(a)初始位置;(b)利用文献[ 11]方法配准;(c)利用文献[ 21]方法配准;(d)利用本文方法配准

Fig. 4. Registration results of Bunny045 and Bunny090. (a) Initial position; (b) registering using method of Ref. [11]; (c) registering using method of Ref. [21]; (d) registering using proposed method

下载图片 查看所有图片

图 5. Dragon000与Dragon048配准结果。(a)初始位置;(b)利用文献[ 11]方法配准;(c)利用文献[ 21]方法配准;(d)利用本文方法配准

Fig. 5. Registration results of Dragon000 and Dragon048. (a) Initial position; (b) registering using method of Ref. [11]; (c) registering using method of Ref. [21]; (d) registering using proposed method

下载图片 查看所有图片

图 6. Happy048与Happy096配准结果。(a)初始位置;(b)利用文献[ 11]方法配准;(c)利用文献[ 21]方法配准;(d)利用本文方法配准

Fig. 6. Registration results of Happy048 and Happy096. (a) Initial position; (b) registering using method of Ref. [11]; (c) registering using method of Ref. [21]; (d) registering using proposed method

下载图片 查看所有图片

使用4.2节的评价标准计算LCP,并记录配准算法的运行时间,与最新的文献[ 11]、[21]以及经典的配准算法Go-ICP和Super4PCS进行对比,对比结果如表2所示。由表2可以看出,对于重叠程度较低的三组点云数据,本文方法的配准精度均高于其他四种方法,而且与表1相比,本文方法在重叠率低的点云数据配准方面,较其他四种方法的配准精度优势更大。对于低重叠率的点云数据,本文方法

较Go-ICP、Super4PCS和文献[ 11]三种方法运行速度较快,较文献[ 21]方法运行速度稍慢。

由表1和表2可以得出,本文方法能够有效提高三维点云配准的精度和时间效率,尤其是在低重叠率三维点云配准方面,本文方法较其他四种方法表现出更多的精度优势。Go-ICP、Super4PCS两种算法是基于无特征的全局点云配准算法,对低重叠率点云配准的计算代价比较高,时间较慢。本文方法有效减少了图匹配问题的配准搜索范围,节省了计算量,比最新的适用于低重叠率的两种方法文献[ 11]、[21]具有更高的配准效果。

除此之外,为了更好地展现本文算法的普适性,对从geometryhub网站获得的两组实际三维扫描数据进行配准验证。这两组数据的初始重叠率大约为40%。实验结果如图7、8所示。其中,图(a)、(b)分别是从两个不同的视角扫描得到的初始点云数据,图(c)为利用本文方法得到的点云配准结果。从图中可以看出,利用本文方法都可以将其很好地配准。上述点云数据的配准性能如表3所示。

图 7. Chair1与Chair2配准结果。(a)初始点云1;(b)初始点云2;(c)利用本文方法配准

Fig. 7. Registration results of Chair1 and Chair2. (a) Initial point cloud 1; (b) initial point cloud 2; (c) registering using proposed method

下载图片 查看所有图片

图 8. Motobike1与Motobike2配准结果。(a)初始点云1;(b)初始点云2;(c)利用本文方法配准

Fig. 8. Registration results of Motobike1 and Motobike2. (a) Initial point cloud 1; (b) initial point cloud 2; (c) registering using proposed method

下载图片 查看所有图片

综上所述,从配准精度和配准时间上综合考虑,本文方法的可行性更好,对初始重叠率较低的点云配准具有普适性。

5 结论

本文提出了一种新的针对低重叠率的三维点云配准算法。将点云的特征值和法线信息综合起来构建多尺度描述符,依据描述符的差异直接进行对应关系点云聚类分块,达到放大点云重叠区域作用的同时,保持着点之间的优质对应关系。与此同时,将得到的对应关系代入原先没有对应关系的凸优化点云配准算法中,缩小了凸优化配准的搜索范围,并为其提供了良好的初值,减少了计算量。实验结果表明,本文算法能够有效提高三维点云配准的精度和时间效率,尤其是在初始重叠率较低的点云数据配准中其精度优势更为明显,相比于最新的适用于低重叠率的两种方法具有更高的配准效果;在配准效率方面,本文方法对于低重叠率的点云配准比基于关键点的融合Super4PCS和ICP的点云配准算法效率更高,是基于快速描述符和对应传播的全局点云配准算法的有效补充。后续将研究重叠率更低和表面更为复杂的实际物体点云模型的配准。