楊高朝

YANG Gaozhao

國(guó)家測(cè)繪地理信息局 海南基礎(chǔ)地理信息中心，?？?570203

Hainan Geomatics Center,National Administration of Surveying,Mapping and Geoinformation,Haikou 570203,China

1 引言

地面三維激光掃描技術(shù)在測(cè)繪領(lǐng)域的應(yīng)用越來越廣泛[1-4]，打破了傳統(tǒng)測(cè)量數(shù)據(jù)獲取與處理模式，推動(dòng)了點(diǎn)云數(shù)據(jù)處理理論的發(fā)展。由于受到測(cè)量設(shè)備和環(huán)境的限制，物體表面完整數(shù)據(jù)的獲得往往需要通過多站掃描才能完成。因此，為了得到物體完整的點(diǎn)云數(shù)據(jù)，需要對(duì)這些局部點(diǎn)云數(shù)據(jù)進(jìn)行坐標(biāo)轉(zhuǎn)換，最終拼接到同一個(gè)坐標(biāo)系里。商業(yè)中應(yīng)用最廣泛最成熟的方法就是通過在重建物體上貼標(biāo)志點(diǎn)[5]來實(shí)現(xiàn)點(diǎn)集的配準(zhǔn)，但有些文物古跡的恢復(fù)和逆向工程存在無法貼標(biāo)簽甚至無法接觸的問題，因此研究點(diǎn)云自動(dòng)配準(zhǔn)算法不但具備重要的實(shí)用價(jià)值，還具有重要的商業(yè)價(jià)值。

目前，許多國(guó)內(nèi)外學(xué)者提出了多種點(diǎn)云配準(zhǔn)算法。到20世紀(jì)90年代中后期比較流行的基于幾何特征的匹配算法，如幾何Hush[6]算法、Spine image等，其中最常見的就是基于迭代最近點(diǎn)（Iterative Closest Point，ICP）的方法[7]。但I(xiàn)CP算法可能存在局部最優(yōu)缺陷，如果兩站點(diǎn)云沒有較高的重疊度或沒有很好的初姿，很可能陷入局部最優(yōu)，最后得到錯(cuò)誤的匹配結(jié)果。Blais等對(duì)點(diǎn)云進(jìn)行隨機(jī)采樣以提高速度，算法精度受到影響。Okatani提出了通過考慮平移和旋轉(zhuǎn)錯(cuò)位的特性來獲得最佳矩陣的方法，該方法易陷入局部收斂。朱延娟等[8]提出了基于點(diǎn)云曲率和鄰域曲率相似度的粗配準(zhǔn)。張廣鵬等[9]提出了主元素貼合法的全局配準(zhǔn)。劉艷豐[10]提出了基于四點(diǎn)算法的全局配準(zhǔn)。

本文在總結(jié)分析國(guó)內(nèi)外學(xué)者專家理論研究的基礎(chǔ)上，針對(duì)上述點(diǎn)云數(shù)據(jù)配準(zhǔn)中所存在的問題，從研究初始匹配入手，結(jié)合改進(jìn)的ICP算法，設(shè)計(jì)出一種利用點(diǎn)云法向量夾角變化度及特征曲率相結(jié)合的自動(dòng)配準(zhǔn)方法。

2 自動(dòng)配準(zhǔn)的主要流程

點(diǎn)云自動(dòng)配準(zhǔn)一般分兩個(gè)過程：粗配準(zhǔn)、精配準(zhǔn)[11-14]。粗配準(zhǔn)是通過特征提取[15]、信號(hào)濾波等手段，來獲得不同視角數(shù)據(jù)集間相對(duì)應(yīng)的控制點(diǎn)，從而估計(jì)出配準(zhǔn)變換參數(shù)。精確配準(zhǔn)采用迭代方式逐漸逼近最佳結(jié)果，由于迭代方式不能保證得到最佳收斂結(jié)果，精確配準(zhǔn)需要粗配準(zhǔn)為其提供一個(gè)迭代的起始位置，來增大收斂到最佳結(jié)果的概率。具體步驟如圖1。

圖1 點(diǎn)云配準(zhǔn)流程圖

首先P、Q兩站點(diǎn)云有一定的重疊度[16]（一般特征比較明顯的點(diǎn)云超過60%即可得到正確的配準(zhǔn)結(jié)果），根據(jù)每個(gè)點(diǎn)云的法矢量大小重采樣點(diǎn)云，為每個(gè)點(diǎn)云定義一個(gè)N維向量，向量值為點(diǎn)云與鄰近點(diǎn)云的法向量夾角。然后根據(jù)這些法向量夾角進(jìn)行對(duì)應(yīng)點(diǎn)搜索，進(jìn)行粗配準(zhǔn)，再進(jìn)行精匹配，求取轉(zhuǎn)換參數(shù)。具體詳細(xì)步驟見下文。

2.1 粗配準(zhǔn)

為了縮小點(diǎn)云之間的旋轉(zhuǎn)和平移錯(cuò)位，使得精確配準(zhǔn)不致得到錯(cuò)誤的結(jié)果，需要進(jìn)行點(diǎn)云粗配準(zhǔn)。在研究上述算法的基礎(chǔ)上，提出了法向量鄰域夾角搜索的粗配準(zhǔn)算法。首先對(duì)每個(gè)點(diǎn)構(gòu)造一個(gè)不變特征向量，即每個(gè)點(diǎn)的法向量與其K-近鄰點(diǎn)法向量[17]的夾角組成的向量，以此向量作為輸入量構(gòu)建k[θ1,θ2,…,θk]T，其中θk為弧度值，在兩點(diǎn)集中搜索最近點(diǎn)。法向量夾角如圖2、圖3。

圖2 模型點(diǎn)云法向量夾角

圖3 目標(biāo)點(diǎn)云法向量夾角

在計(jì)算法向量夾角的過程中，法向量的估計(jì)是比較重要的一步。由于采集到的原始數(shù)據(jù)往往很龐大，為了加快速度，先對(duì)原始數(shù)據(jù)進(jìn)行點(diǎn)云抽希[18]，然后采用基于局部表面擬合的方法進(jìn)行法向量估計(jì)[19]（為了使配準(zhǔn)速度更快，本文選擇平面擬合方案，但曲面擬合可能效果更優(yōu)，但費(fèi)時(shí)較多）。為此，對(duì)于點(diǎn)云中的每個(gè)掃描點(diǎn)pi，搜索與其最相近的K個(gè)相鄰點(diǎn)，然后計(jì)算這些點(diǎn)最小二乘意義上的局部平面P，此平面可以表述如下：

式中，n為平面P的法向量；d為P到坐標(biāo)原點(diǎn)的距離。

計(jì)算該點(diǎn)與最近點(diǎn)的法向量夾角余弦值，如圖2和圖3所示。當(dāng)法向量夾角值越大時(shí)，余弦值就越小，說明該點(diǎn)所在區(qū)域起伏變化較大，特征較為明顯，反之則較為平坦。因?yàn)槠教共糠值狞c(diǎn)的幾何特征較為接近，不適合作為特征點(diǎn)集進(jìn)行配準(zhǔn)，所以選取適當(dāng)?shù)拈撝郸?0.99（經(jīng)過實(shí)驗(yàn)求得），保留法向量夾角余弦值小于ε的點(diǎn)。

設(shè)待配準(zhǔn)的兩個(gè)點(diǎn)云數(shù)據(jù)分別為P和Q，利用上述特征點(diǎn)的提取方法，分別對(duì)兩個(gè)點(diǎn)云進(jìn)行特征點(diǎn)提取，得到P的特征點(diǎn)集為Q的特征點(diǎn)集為，其中m和n分別為P和Q的特征點(diǎn)的個(gè)數(shù)。

2.2 獲取匹配點(diǎn)對(duì)

2.2.1 對(duì)應(yīng)點(diǎn)對(duì)的初步獲取

幾何特征信息種類越多，計(jì)算就越復(fù)雜，但是特征信息量少，特征識(shí)別度低，就很難收斂。經(jīng)過多種特征點(diǎn)云的實(shí)驗(yàn)，最終決定選取兩種基本幾何特征搜索匹配點(diǎn)[20]。下面詳細(xì)介紹點(diǎn)云幾何特征的建立過程。

（1）對(duì)于點(diǎn)集Pt中的每個(gè)點(diǎn)pti，以該點(diǎn)鄰域法向量夾角加權(quán)值k[θ1,θ2,…,θk]T的特征距離作為該點(diǎn)的特征量：

（2）對(duì)于點(diǎn)集Pt中的每個(gè)點(diǎn)pti，利用鄰域點(diǎn)集中k個(gè)最近點(diǎn)，求出鄰域的重心O(pti)，以該點(diǎn)到其鄰域重心的距離值作為該點(diǎn)的特征量：

對(duì)于點(diǎn)集Pt中的每個(gè)點(diǎn)pti，為其在點(diǎn)集Qt中搜索對(duì)應(yīng)點(diǎn)qtj，若qtj是pti的對(duì)應(yīng)點(diǎn)或近似對(duì)應(yīng)點(diǎn)，則對(duì)兩點(diǎn)來說，以上提出的兩種幾何特征應(yīng)該相同或近似相同。因此對(duì)于點(diǎn)集Pt中的一點(diǎn)pti和點(diǎn)集Qt中的一點(diǎn)qtj，只要滿足以下4個(gè)條件則認(rèn)為它們具有對(duì)應(yīng)關(guān)系：

閾值ε1、ε2取值過大，會(huì)導(dǎo)致產(chǎn)生很多錯(cuò)誤的匹配點(diǎn)對(duì)，取值過小，無法獲取足夠的匹配點(diǎn)對(duì)，根據(jù)其他學(xué)者的經(jīng)驗(yàn)以及本人反復(fù)試驗(yàn)，式中取ε1=ε2=0.01，使用以上關(guān)系進(jìn)行過濾并初步建立點(diǎn)對(duì)應(yīng)關(guān)系。由于點(diǎn)云中本身就存在很多特征相似區(qū)域，使用多個(gè)特征關(guān)系進(jìn)行對(duì)應(yīng)點(diǎn)的建立可以很好地避免出現(xiàn)大量一對(duì)多的情形，再通過合理閾值的選擇可以為點(diǎn)集Pt中的每個(gè)點(diǎn)pti初步找到其在Qt上的對(duì)應(yīng)點(diǎn)。最后建立初始匹配點(diǎn)對(duì)數(shù)組并記為：

式中，N為匹配點(diǎn)對(duì)的數(shù)量。

2.2.2 錯(cuò)誤對(duì)應(yīng)點(diǎn)對(duì)的剔除

采用上述方法確立的對(duì)應(yīng)點(diǎn)集中會(huì)存在錯(cuò)誤對(duì)應(yīng)點(diǎn)對(duì)（或稱為不可靠、不兼容點(diǎn)對(duì)），或者會(huì)產(chǎn)生多對(duì)互相矛盾的點(diǎn)對(duì)。為此，必須引入評(píng)價(jià)標(biāo)準(zhǔn)或約束準(zhǔn)則去除錯(cuò)誤對(duì)應(yīng)點(diǎn)對(duì)，這是許多ICP改進(jìn)算法的關(guān)鍵，也是國(guó)內(nèi)外諸多學(xué)者研究的重點(diǎn)。區(qū)別于深度圖像可利用圖像自身的灰度、顏色、紋理等輔助信息[21-23]，對(duì)于點(diǎn)云，能夠利用的信息主要是點(diǎn)云自身的曲面幾何特征。

如果以上獲取的初始匹配點(diǎn)對(duì)都是正確匹配點(diǎn)對(duì)，依據(jù)幾何剛性變換原理可知，對(duì)于任意兩對(duì)匹配點(diǎn)對(duì)和應(yīng)滿足點(diǎn)對(duì)間距離的偏差在極小的允許范圍之內(nèi)，即可表示選取閾值ε=0.02，對(duì)每個(gè)點(diǎn)對(duì)計(jì)算P中除其以外與其符合剛性距離約束條件的點(diǎn)的數(shù)目Qi，若P中一個(gè)點(diǎn)對(duì)滿足：

2.2.3 匹配點(diǎn)對(duì)搜索精度檢驗(yàn)

為了檢驗(yàn)匹配點(diǎn)搜索精度，引入度量ε5：

2.2.4 計(jì)算初始匹配參數(shù)

在求得目標(biāo)點(diǎn)集P和Q匹配點(diǎn)對(duì)后，需要計(jì)算兩點(diǎn)云之間的平移變換矩陣，比較典型的方法有最小二乘法、四元數(shù)法[24]、矩陣的奇異值分解法（SVD）[25]等。經(jīng)過對(duì)比發(fā)現(xiàn)，最小二乘法求解收斂速度更快。根據(jù)旋轉(zhuǎn)矩陣和平移向量對(duì)目標(biāo)點(diǎn)云進(jìn)行初步轉(zhuǎn)換。

2.3 精配準(zhǔn)

2.3.1 改進(jìn)后的ICP算法

進(jìn)行初始配準(zhǔn)后，兩站點(diǎn)云的初始姿態(tài)基本上接近，可以作為ICP算法的初值進(jìn)行精匹配。

根據(jù)經(jīng)典的ICP算法，本文在以下方面對(duì)ICP算法進(jìn)行改進(jìn)，以提高算法的效率：

（1）使用基于法向夾角的點(diǎn)云重采樣[26]；

（2）應(yīng)用法向投影選取對(duì)應(yīng)點(diǎn)，該搜索過程中使用kd-tree或八叉樹可以提高速度；

（3）權(quán)重設(shè)置選取距離的平方加權(quán)平均值；

（4）誤匹配點(diǎn)剔除選取距離和法向夾角相兼容的原則；

（5）目標(biāo)函數(shù)選取點(diǎn)到面的誤差最小。

法向投影-鄰域搜索將源數(shù)據(jù)投影到目標(biāo)點(diǎn)云上，然后在目標(biāo)點(diǎn)云上進(jìn)行最近點(diǎn)搜索（點(diǎn)周圍的幾何特征信息也作為判斷的標(biāo)準(zhǔn)如曲率等）。

誤差準(zhǔn)則函數(shù)選擇：

（1）如圖4，對(duì)應(yīng)點(diǎn)對(duì)直接距離的平方和最小[27]，這一形式是最常用的誤差準(zhǔn)則。

其中，pi是目標(biāo)點(diǎn)云的坐標(biāo)向量；qki是模型點(diǎn)云的坐標(biāo)向量；R是變換矩陣；t是平移矩陣；ε是對(duì)應(yīng)點(diǎn)對(duì)直接距離的平方和。當(dāng)ε達(dá)到極值最小時(shí)，此時(shí)的R和t就是要求的值。

圖4 最近點(diǎn)搜索算法

（2）如圖5，每個(gè)源數(shù)據(jù)點(diǎn)到以其對(duì)應(yīng)點(diǎn)法向量建立的局部平面的距離平方和最小。

其中，ni是pi的法向量。其他變量的含義同上。表示R?qki+t投影到pi點(diǎn)的切平面的距離的平方。采用此目標(biāo)函數(shù)，只要兩個(gè)點(diǎn)在與切平面平行的方向上滑動(dòng)，距離都不變，因此在優(yōu)化過程中允許滑動(dòng)，這是此目標(biāo)函數(shù)獲得較快收斂速度的直觀原因。

圖5 法向量垂直搜索算法

根據(jù)圖6實(shí)驗(yàn)結(jié)果，可知第二種方案無論是收斂速度上還是精度上都優(yōu)于第一種方案，本文選取第二種方案作為誤差準(zhǔn)則函數(shù)。

圖6 精配準(zhǔn)誤差圖（比較點(diǎn)到點(diǎn)與點(diǎn)到面的差別）

2.3.2 精匹配實(shí)現(xiàn)步驟

運(yùn)用改進(jìn)的ICP算法進(jìn)行精細(xì)配準(zhǔn)步驟如下：

步驟1構(gòu)造參與算法的兩個(gè)初始點(diǎn)集，分別記作P和X；設(shè)定閾值Δ＞0，把點(diǎn)云的特征加入到配準(zhǔn)策略上來，比如點(diǎn)云的曲率、法向量夾角等。

步驟2初始化p0=p,k=0。

步驟4計(jì)算Rk和Tk,使

最小，ni是目標(biāo)點(diǎn)云的法向量。

步驟5對(duì)p0進(jìn)行三維變換，得到pk+1=Rkp0+Tk。

步驟6k=k+1,如果k=1,轉(zhuǎn)到步驟3；如果k≥2,且則終止迭代并轉(zhuǎn)到步驟7，否則轉(zhuǎn)到步驟3。

步驟7最后所得到的Rk、Tk即為所要求的精細(xì)配準(zhǔn)參數(shù)，利用所得的參數(shù)，對(duì)第一片點(diǎn)云進(jìn)行轉(zhuǎn)換，配準(zhǔn)結(jié)束。

2.3.3 誤差評(píng)定

為了評(píng)價(jià)點(diǎn)云自動(dòng)匹配的精度，由文獻(xiàn)[28]可得：

其中，mp為配準(zhǔn)中誤差；pi、p′i為公共點(diǎn)對(duì)轉(zhuǎn)換前后向量值。

3 實(shí)驗(yàn)及結(jié)果分析

在Windows平臺(tái)下利用MFC和OpenGL軟件實(shí)現(xiàn)了算法和配準(zhǔn)界面，如圖7。采用MFC單文檔結(jié)構(gòu)，利用分割視圖法，管理3個(gè)視圖。第二個(gè)視圖顯示model（模型）點(diǎn)集，第三個(gè)是data（目標(biāo)）點(diǎn)集，然后點(diǎn)擊菜單打開配準(zhǔn)窗口，可在第一個(gè)視圖同時(shí)顯示模型點(diǎn)集、目標(biāo)點(diǎn)集和配準(zhǔn)后的點(diǎn)集。提供手動(dòng)配準(zhǔn)，自動(dòng)配準(zhǔn)，Mesh（或點(diǎn)云）旋轉(zhuǎn)平移、點(diǎn)云刪除、重選、縮放以及點(diǎn)云數(shù)據(jù)簡(jiǎn)化等功能。

圖7 基于MFC和OpenGL實(shí)現(xiàn)的配準(zhǔn)軟件

3.1 實(shí)驗(yàn)結(jié)果

為了驗(yàn)證算法的有效性，本實(shí)驗(yàn)利用網(wǎng)上下載的不同視角的兩片點(diǎn)云，點(diǎn)云個(gè)數(shù)29876個(gè)。首先對(duì)點(diǎn)云進(jìn)行壓縮，保留一些特征比較明顯的點(diǎn)，壓縮后的點(diǎn)云為7652個(gè)，如圖8所示。接著對(duì)兩站點(diǎn)云進(jìn)行粗配準(zhǔn)，結(jié)果如圖9（a），最后再對(duì)點(diǎn)云進(jìn)行精匹配，結(jié)果如圖9（b），配準(zhǔn)結(jié)束后還可以查看本次實(shí)現(xiàn)的配準(zhǔn)精度以及配準(zhǔn)所花費(fèi)的時(shí)間。為了更能顯示本文方案的可用性，選取一個(gè)現(xiàn)今比較成熟的方案，應(yīng)用同樣的數(shù)據(jù)同時(shí)實(shí)驗(yàn)，以便和本文方案進(jìn)行比較。經(jīng)過對(duì)比，選取了文獻(xiàn)[20]配準(zhǔn)算法對(duì)掃描點(diǎn)云進(jìn)行迭代配準(zhǔn)，其配準(zhǔn)效果如圖10所示。為更形象地查看迭代過程及結(jié)果，根據(jù)每次迭代的殘差，結(jié)合MATLAB生成線狀圖，如圖11。

圖8 兩組點(diǎn)云初始位置

圖9 本文配準(zhǔn)效果

圖10 文獻(xiàn)[20]配準(zhǔn)效果

圖11 文獻(xiàn)[20]和本文算法迭代配準(zhǔn)誤差圖

為了更直觀地分析問題，對(duì)本次實(shí)驗(yàn)進(jìn)行統(tǒng)計(jì)，如表1。

表1 本文算法和文獻(xiàn)[20]算法的比較

3.2 結(jié)果分析

為了定量比較上述兩組實(shí)驗(yàn)結(jié)果，通過計(jì)算多視點(diǎn)云重疊區(qū)域的最鄰近點(diǎn)對(duì)距離平方和的平均值的平方根，作為整體配準(zhǔn)后的距離中誤差，用于評(píng)估多視點(diǎn)云間的配準(zhǔn)精度，同樣的數(shù)據(jù)運(yùn)用式（12）求得本文算法配準(zhǔn)后的中誤差為30.4 μm,文獻(xiàn)[20]配準(zhǔn)后的中誤差為38.7 μm。已知原始點(diǎn)云的密度為0.15 mm，實(shí)驗(yàn)表明，本文算法要優(yōu)于文獻(xiàn)[20]的算法。另外從實(shí)驗(yàn)結(jié)果還可以看出，本文算法配準(zhǔn)中誤差近似為0.3倍點(diǎn)云密度，因此本文配準(zhǔn)算法的精度比較理想。

在算法迭代收斂方面，使用本文算法進(jìn)行配準(zhǔn)共迭代9次滿足平差終止條件。其中3個(gè)平移參數(shù)的改正數(shù)收斂較快，3次迭代時(shí)已收斂；3個(gè)旋轉(zhuǎn)角參數(shù)的改正數(shù)迭代后期收斂相對(duì)較慢，第5次迭代時(shí)3個(gè)旋轉(zhuǎn)角參數(shù)的改正數(shù)絕對(duì)值都小于6×10-5，而本文中設(shè)定的限差為4×10-5，直到第9次滿足限差要求，迭代結(jié)束。使用文獻(xiàn)[20]算法同樣對(duì)兩站點(diǎn)云進(jìn)行配準(zhǔn)，一直迭代13次才滿足迭代終止條件，從收斂效率上可以看出，本文算法要優(yōu)于文獻(xiàn)[20]的算法。

在算法運(yùn)行速度方面，本文方案粗配準(zhǔn)和精配準(zhǔn)一共用了6.04 s，文獻(xiàn)[20]一共用了8.09 s。由于本文方案盡管使用的是基于點(diǎn)云特征的配準(zhǔn)方法，但并沒有直接計(jì)算擬合曲面曲率，直接通過法向量夾角作為三維點(diǎn)云剛體不變的特征來實(shí)現(xiàn)自動(dòng)配準(zhǔn)。因此和傳統(tǒng)利用點(diǎn)云擬合曲率實(shí)現(xiàn)點(diǎn)云自動(dòng)配準(zhǔn)的方案相比，本文方案提高了運(yùn)行速度。

綜上所述，本文算法無論是從配準(zhǔn)精度上還是從迭代效率上都優(yōu)于文獻(xiàn)[20]。實(shí)驗(yàn)配準(zhǔn)效果圖和實(shí)驗(yàn)獲得的對(duì)比充分說明了本文算法的有效性。

4 結(jié)論與展望

本文將利用鄰域法向夾角搜索進(jìn)行的初始配準(zhǔn)作為初始值，并利用改進(jìn)的ICP算法對(duì)提供的初始值進(jìn)行精配準(zhǔn)。實(shí)驗(yàn)結(jié)果表明，改進(jìn)后的ICP算法在對(duì)應(yīng)點(diǎn)的搜索時(shí)間和配準(zhǔn)精度等方面都優(yōu)于經(jīng)典的ICP算法，既解決了傳統(tǒng)的ICP算法進(jìn)行配準(zhǔn)，算法陷入局部最小，從而得到錯(cuò)誤的配準(zhǔn)結(jié)果的問題，又提高了點(diǎn)云配準(zhǔn)的效率。點(diǎn)云中本身就存在很多特征相似區(qū)域，單一的特征匹配識(shí)別度低，在實(shí)驗(yàn)中發(fā)現(xiàn)，很多特征較少的點(diǎn)云也會(huì)陷入局部收斂，但采用多個(gè)幾何特征來搜索對(duì)應(yīng)點(diǎn)所花費(fèi)的時(shí)間也會(huì)相應(yīng)增加。因此，如何在兩者間找到一個(gè)平衡點(diǎn)并在實(shí)際生產(chǎn)中廣泛應(yīng)用是接下來需要解決的問題。