吳明東，付建紅，孟慶祥

(1.河南理工大學(xué) 測繪與國土信息工程學(xué)院，河南 焦作 454000；2.武漢大學(xué) 遙感信息工程學(xué)院，武漢 430079)

0 引言

機載激光掃描(light detection and ranging，LiDAR)集全球定位系統(tǒng)、慣性導(dǎo)航系統(tǒng)與激光測量3種技術(shù)于一體，能夠快速高效地獲取地物表面每個采樣點的精確三維坐標(biāo)，廣泛應(yīng)用于數(shù)字城市、災(zāi)害監(jiān)測、海岸工程、林業(yè)等方面[1]。在對一個區(qū)域進行機載激光掃描時，往往需要多航帶或多次掃描才能完成，而多次掃描的激光點云需要通過配準(zhǔn)轉(zhuǎn)換到統(tǒng)一的坐標(biāo)系中[2-3]，這也是利用點云進行三維重建中最重要且最困難的步驟，直接影響著三維重建的精度與效率，是后續(xù)點云重建的技術(shù)基礎(chǔ)[4]。

目前，點云配準(zhǔn)的算法有很多，主要分為2類。一類是基于點云數(shù)據(jù)特性建立的特征匹配法，該方法利用物體表面上的幾何特征信息，對特征進行提取和描述，通過確定特征之間的對應(yīng)關(guān)系來實現(xiàn)點云配準(zhǔn)。劉劍等[5]提出一種快速點特征直方圖(FPFH)特征描述與Delaunay三角剖分相結(jié)合的特征點配準(zhǔn)方法，有效地提升了點云配準(zhǔn)的精度和效率。李健等[6]結(jié)合當(dāng)下最熱門的深度學(xué)習(xí)技術(shù)，將深度學(xué)習(xí)用于三維點云局部特征的提取，融合點云配準(zhǔn)算法解決了多維度大角度變化的點云配準(zhǔn)問題，使深度學(xué)習(xí)在點云配準(zhǔn)領(lǐng)域?qū)崿F(xiàn)了很好的應(yīng)用。盛敏等[7]針對RGBD點云數(shù)據(jù)的配準(zhǔn)問題，利用待配準(zhǔn)的RGBD點云模型的曲率和顏色特征度定義了一種基于二者特征相似性的點云配準(zhǔn)方法。但是這類方法需要點云包含豐富的顏色信息，同時對點云片的曲率有要求。機載LiDAR數(shù)據(jù)覆蓋范圍廣、點云密度小，覆蓋區(qū)內(nèi)地物數(shù)據(jù)組成單一(主要為地表和建筑物屋頂)，同時缺少顏色信息，對曲率的計算和特征的提取相對困難。第二類是基于點云集之間相互迭代計算的方法。其中比較主流的是迭代最近點 (iterative closest point，ICP)算法，該算法是Besl等[8]于1992年提出的，其基本思想是通過旋轉(zhuǎn)、平移使2個點集之間的歐氏距離最小。文獻[9-10]分析了ICP算法的基本原理并用于對地面掃描數(shù)據(jù)的配準(zhǔn)。喬世權(quán)等[11]通過定義源數(shù)據(jù)與目標(biāo)數(shù)據(jù)的曲面距離，利用PCA算法計算旋轉(zhuǎn)矩陣的初始值改進了傳統(tǒng)的ICP算法，并用于點云配準(zhǔn)。針對點云配準(zhǔn)中的尺度和收斂速度問題，趙夫群等[12]通過在ICP算法中加入帶邊界的尺度矩陣，提出了一種改進的尺度迭代最近點算法，解決點云配準(zhǔn)中尺度變換的問題。為了解決強噪聲及密度不均勻點云高效配準(zhǔn)的難題，彭真等[13]提出了一種基于關(guān)鍵點提取與優(yōu)化的改進的 ICP點云配準(zhǔn)算法?；谔卣鞯钠ヅ浜突贗CP算法的匹配在地面激光掃描數(shù)據(jù)，特別是工業(yè)測量等方面得到了廣泛的應(yīng)用。對于機載LiDAR條帶的配準(zhǔn)問題，這些算法的研究和應(yīng)用相對較少。

結(jié)合機載LiDAR點云數(shù)據(jù)特點，本文提出一種基于機載LiDAR點云高程相關(guān)的配準(zhǔn)新方法。該方法首先在基準(zhǔn)點云中選取某一目標(biāo)點，通過目標(biāo)點與周圍若干點進行曲面擬合形成擬合曲面；然后在另一片點云建立的搜索區(qū)中遍歷采樣點，比較每個采樣點和其鄰近點組成的點集在擬合曲面上的擬合高程與點集的真實高程數(shù)據(jù)之間相關(guān)系數(shù)，選擇相關(guān)系數(shù)最大的采樣點作為同名點；最后通過2片點云中若干同名點對完成點云配準(zhǔn)。

1 基于高程相關(guān)的機載點云配準(zhǔn)

1.1 點云預(yù)處理

對于LiDAR點云數(shù)據(jù)而言，點云數(shù)據(jù)量越大將增加點云配準(zhǔn)的時間、影響點云配準(zhǔn)的效率；點云噪聲會對點云同名點的尋找產(chǎn)生干擾，導(dǎo)致錯誤匹配點對[14]。因此，為保證點云配準(zhǔn)的可靠性，在對點云進行配準(zhǔn)前先對點云進行預(yù)處理。點云的預(yù)處理包括點云的去噪和精簡。K-D樹在點云去噪中有著廣泛的運用，在進行點云去噪中，對全部點云進行遍歷，并逐一查找各點的鄰域，計算鄰域內(nèi)各點到中心點的平均距離，當(dāng)平均距離過大時，可將該點視為噪聲點。點云精簡的過程中，在整體降低數(shù)據(jù)量的同時盡可能多地保留一些能夠反映環(huán)境的表面特征的點位數(shù)據(jù)。

1.2 曲面擬合模型建立

對于任意機載LiDAR三維點集M中的一點p，定義p的一個鄰域R，將鄰域R內(nèi)的所有點進行曲面擬合，擬合的曲面作為這些點所代表的局部真實場景的逼近，稱為最鄰近曲面。常見的曲面擬合算法有很多[15-19]，本文采用二次曲面擬合模型算法來擬合曲面，如式(1)所示。

(1)

式中：(Xi，Yi，Zi)為鄰域R中第i個點的三維空間坐標(biāo)值；a、b、c、d、e、f為相應(yīng)的二次曲面擬合系數(shù)。

本文主要思想是在基準(zhǔn)點云上選擇局部點擬合曲面，計算出待配準(zhǔn)點云上的局部點云在擬合曲面上的高程，求解該局部點云在擬合曲面上的高程和其真實高程之間的相關(guān)系數(shù)，并在搜索區(qū)域內(nèi)移動擬合曲面，取相關(guān)系數(shù)最大時的局部點云作為基準(zhǔn)點云的對應(yīng)點云，進而確定同名點，完成點云配準(zhǔn)。因此，如何快速、準(zhǔn)確地選取合適的局部擬合曲面對于后續(xù)的配準(zhǔn)工作有決定意義。選取的擬合曲面應(yīng)遵循以下原則。

1)參與擬合的點應(yīng)選擇地形起伏較明顯或與周圍地面有明顯區(qū)分的地物點。

2)曲面的范圍選取不宜過大，必須保證所擬合的曲面為局部小范圍曲面。

3)原始點云擬合后的曲面偏移距離要小。參與擬合的點云為局部真實的地面或地物點，若某一點或多點與擬合后的曲面偏移距離較大，則說明該局部點云不滿足曲面模型，將重新選擇點云。

圖1 局部點云查找示意圖

局部點云選擇如圖1所示。用于擬合曲面的局部點云選擇過程為:假設(shè)點p為點云數(shù)據(jù)中滿足曲面選取原則中的一個點，以該點為中心尋找在周圍點云中距離點p一定閾值范圍內(nèi)的k個點，利用選取的k個點擬合二次曲面，并計算每個點的擬合殘差，若擬合殘差超限，則移動中心點位置，重新擬合曲面，直到生成合適的最鄰近曲面。

1.3 高程相關(guān)系數(shù)求解

在本文中，求2組變量Z和Z′的相關(guān)系數(shù)ρ(Z，Z′)。相關(guān)系數(shù)的取值范圍為[-1，1]，當(dāng)ρ(Z，Z′)取正時，表示正相關(guān)；當(dāng)ρ(Z，Z′)取負時表示負相關(guān)；ρ(Z，Z′)取0時，表示不相關(guān)。相關(guān)系數(shù)越大，表示2組變量相關(guān)性越大，當(dāng)ρ(Z，Z′)=1時，說明2組變量完全相關(guān)。

1.4 點云配準(zhǔn)

設(shè)M、N分別為存在交集的2個相鄰機載LiDAR點云片S1和S2上的點集，S1、S2重合區(qū)域C=S1∩S2。記M落入C區(qū)域的采樣點所構(gòu)成的集合為CM={Mi|Mi∈R3}，N中落入C區(qū)域的采樣點所構(gòu)成的集合為CN={Ni|Ni∈R3}。具體關(guān)系如圖2所示。

圖2 點云數(shù)據(jù)集關(guān)系圖

將點集M固定，點集N作為浮動點集參與配準(zhǔn)。然后在CM中按照1.2節(jié)的方法選擇一個初始目標(biāo)點作為待定點p，以該點為中心搜索范圍l以內(nèi)的鄰近點，建立擬合曲面ρ。為了在CN上搜索點p的同名點，須估計出該同名點可能存在的范圍，建立一個搜索區(qū)域，如圖3所示。在搜索區(qū)域中選擇一點q作為局部采樣點集的中心點，以l為半徑查找q的鄰接點組成采樣點集。同時，為了方便求取搜索區(qū)中采樣點集的平面坐標(biāo)落在擬合曲面的坐標(biāo)范圍內(nèi)，須對目標(biāo)點集和采樣點集進行坐標(biāo)中心化處理，即將點集中心點平面坐標(biāo)設(shè)置為(0，0)，點集中其余各點依次用相對于該中心點的坐標(biāo)表示。通過遍歷規(guī)則在搜索區(qū)域中選擇不同的采樣點，對比不同采樣點集中各點在擬合曲面上的高程和各點真實的高程之間的相關(guān)系數(shù)，選擇相關(guān)系數(shù)最大的一組采樣點集合作為目標(biāo)點集的對應(yīng)點集，從而尋找出同名點對。

圖3 搜索區(qū)域內(nèi)查找點示意圖

點云配準(zhǔn)的具體過程如下。

1)2片相鄰點云S1和S2，將S1固定作為基準(zhǔn)點云，S2作為浮動點云與S1配準(zhǔn)。

2)在CM中選擇采樣點p作為基準(zhǔn)點云的目標(biāo)點，并以p點為中心、l為半徑搜索p點的鄰近點，組成目標(biāo)點集。以p點的平面坐標(biāo)為(0，0)更改目標(biāo)點集的平面坐標(biāo)。

3)將p點及鄰近點按照式(1)生成擬合曲面ρ，計算點集中每一點到擬合曲面的偏移距離，若偏移距離大于設(shè)定閾值則重新選取采樣點，重復(fù)以上步驟，直到找到合適的曲面。

4)在CN中預(yù)估p的同名點可能出現(xiàn)的范圍，建立搜索區(qū)。

5)在搜索區(qū)選擇采樣點q，并以q為中心、l為半徑查找q點的鄰近點，組成搜索點集。以q點的平面坐標(biāo)為(0，0)更改搜索點集的平面坐標(biāo)。

6)設(shè)搜索點集中每一個點的高程為Zi，并計算每一個點在擬合曲面上對應(yīng)的高程Z′i。

7)計算Zi和Z′i的集合Z、Z′的相關(guān)系數(shù)ρ(Z，Z′)。

8)在搜索區(qū)中選擇其他采樣點，循環(huán)步驟5)、步驟6)、步驟7)，直到搜索區(qū)中的采樣點遍歷完。

9)計算相關(guān)系數(shù)ρmax=max{ρi}，對應(yīng)搜索點集的采樣點qi即為p點的同名點。

10)循環(huán)步驟2)至步驟9)，在CM中選擇若干目標(biāo)點，并在CN中尋找對應(yīng)的同名點，生成同名點對(pk，qk)。

11)遍歷所有的同名點對(pk，qk)，根據(jù)同名點之間的坐標(biāo)關(guān)系計算出旋轉(zhuǎn)和平移矩陣，反復(fù)迭代計算，直到達到最佳配準(zhǔn)效果，完成配準(zhǔn)。

點云配準(zhǔn)的具體實施流程如圖4所示。

圖4 點云配準(zhǔn)總流程圖

2 實驗及結(jié)果

機載LiDAR數(shù)據(jù)覆蓋范圍較廣，且主要反映地表及地表各種建筑物信息。為驗證本文方法的正確性和有效性，選取了4組機載LiDAR點云數(shù)據(jù)進行配準(zhǔn)實驗。實驗系統(tǒng)環(huán)境為內(nèi)存8 GB、Windows10操作系統(tǒng)，數(shù)據(jù)處理的開發(fā)平臺為Visual Studio 2012，編程語言為C#語言，點云的可視化軟件為Cloud Compare 3D點云處理軟件。

第1組實驗為某地表A的相鄰2條航帶的機載LiDAR數(shù)據(jù)，地面分辨率10 cm。圖5(a)為該地表配準(zhǔn)前的原始點云數(shù)據(jù)(紅色部分表示基準(zhǔn)點云，白色部分表示待配準(zhǔn)點云)，可以看出2片點云在配準(zhǔn)前存在明顯的位置偏差。圖5(b)為傳統(tǒng)ICP算法配準(zhǔn)結(jié)果，配準(zhǔn)精度RMS(root mean square)為0.204 69 m。圖5(c)為本文方法的配準(zhǔn)結(jié)果，配準(zhǔn)精度RMS為0.218 24 m。

圖5 地表A配準(zhǔn)圖

第2組實驗為某地表B的相鄰2條航帶的機載LiDAR數(shù)據(jù)，地面分辨率10 cm。圖6(a)為該地表配準(zhǔn)前的原始點云數(shù)據(jù)(紅色部分表示基準(zhǔn)點云，白色部分表示待配準(zhǔn)點云)。第2組數(shù)據(jù)與第1組數(shù)據(jù)均為略有起伏的地表，與第1組數(shù)據(jù)的區(qū)別在于第2組數(shù)據(jù)地表B中存在一棟規(guī)則建筑物。圖6(b)為傳統(tǒng)ICP算法配準(zhǔn)結(jié)果，配準(zhǔn)精度RMS為0.412 20 m。圖6(c)為本文方法的配準(zhǔn)結(jié)果，配準(zhǔn)精度RMS為0.398 42 m。

第3組實驗為某建筑物的相鄰2條航帶的機載LiDAR數(shù)據(jù)，地面分辨率10 cm。圖7(a)為該建筑物配準(zhǔn)前的原始點云數(shù)據(jù)(紅色部分表示基準(zhǔn)點云，白色部分表示待配準(zhǔn)點云)。圖7(b)為傳統(tǒng)ICP算法配準(zhǔn)結(jié)果，配準(zhǔn)精度RMS為0.426 53 m。圖7(c)為本文方法的配準(zhǔn)結(jié)果，配準(zhǔn)精度RMS為0.385 73 m。

第4組實驗為某建筑物屋頂?shù)南噜?條航帶的機載LiDAR數(shù)據(jù)。圖8(a)為該建筑物屋頂配準(zhǔn)前的原始點云數(shù)據(jù)(紅色部分表示基準(zhǔn)點云，白色部分表示待配準(zhǔn)點云)。圖8(b)為傳統(tǒng)ICP算法配準(zhǔn)結(jié)果，配準(zhǔn)精度RMS為0.493 54 m。圖8(c)為本文方法的配準(zhǔn)結(jié)果，配準(zhǔn)精度RMS為0.437 11 m。

圖6 地表B配準(zhǔn)圖

圖7 建筑物配準(zhǔn)圖

圖8 建筑物屋頂配準(zhǔn)圖

在上述幾組實驗中，對點云預(yù)處理后，使用了傳統(tǒng)ICP算法與本文所述算法進行了對比實驗，實驗相關(guān)參數(shù)及結(jié)果見表1所示。對比以上4組實驗的結(jié)果以及表1數(shù)據(jù)中的配準(zhǔn)精度可以看出：就機載LiDAR數(shù)據(jù)而言，對于地表或含有少量建筑物的地表點云數(shù)據(jù)的配準(zhǔn)方面，傳統(tǒng)ICP算法和本文提出的算法區(qū)別不大，均能達到較好的配準(zhǔn)效果；對于單個地面建筑物的配準(zhǔn)而言，本文提出的算法通過多次迭代后的配準(zhǔn)精度與傳統(tǒng)ICP算法的配準(zhǔn)精度相當(dāng)，甚至略優(yōu)于ICP算法的配準(zhǔn)精度。

表1 不同算法的配準(zhǔn)結(jié)果

機載LiDAR數(shù)據(jù)為大范圍內(nèi)的條帶點云數(shù)據(jù)，點云密度較小。以上幾組實驗數(shù)據(jù)均為局部小范圍的數(shù)據(jù)，其實驗結(jié)果存在一定局限。為驗證本文算法對大范圍地表點云的配準(zhǔn)效果，進行了2組實驗(實驗5、實驗6)。實驗5為某城區(qū)相鄰2條航帶的機載LiDAR點云數(shù)據(jù)，覆蓋范圍約為600 m×200 m。區(qū)域內(nèi)多為規(guī)則建筑物屋頂、馬路及景觀樹。圖9(a)為該城區(qū)2條相鄰航帶的原始LiDAR點云數(shù)據(jù)(紅色部分、白色部分分別表示不同的條帶數(shù)據(jù))。圖9(c)為該城區(qū)點云數(shù)據(jù)利用傳統(tǒng)ICP算法配準(zhǔn)結(jié)果。圖9(e)為該城區(qū)使用本文方法的配準(zhǔn)結(jié)果。實驗6為某地表相鄰2條航帶的機載LiDAR點云數(shù)據(jù)，覆蓋范圍約為300 m×400 m，區(qū)域內(nèi)為坡度較緩的地表(地面有少量樹林、一條馬路及路燈若干)。圖9(b)為該地表2條相鄰航帶的原始LiDAR點云數(shù)據(jù)(紅色部分、白色部分分別表示不同的條帶數(shù)據(jù))。圖9(d)為該地表數(shù)據(jù)利用傳統(tǒng)ICP算法配準(zhǔn)結(jié)果。圖9(f)為該地表使用本文方法的配準(zhǔn)結(jié)果。

圖9 實驗5、實驗6機載LiDAR點云配準(zhǔn)圖

實驗5、實驗6的配準(zhǔn)精度結(jié)果如表2所示，表中還列出了2組實驗的點云配準(zhǔn)變換參數(shù)。從圖9和表2中的實驗結(jié)果可以看出，對于大范圍的機載LiDAR條帶點云數(shù)據(jù)的配準(zhǔn)，本文所提出的基于機載LiDAR局部點云擬合高程相關(guān)的配準(zhǔn)算法與經(jīng)典的ICP算法精度相差無幾，2組實驗都達到了較高的精度水平。而本文提出的方法，通過在基準(zhǔn)點云中選取特征點，在待配準(zhǔn)點云中使用本文提出的算法計算特征點對應(yīng)的同名點，通過所計算的若干同名點對反復(fù)迭代完成2片點云的配準(zhǔn)，在2片待配準(zhǔn)點云初始位置不是很好的情況下也能完成點云的高精度配準(zhǔn)。

表2 實驗5、實驗6不同方法配準(zhǔn)參數(shù)及精度表

3 結(jié)束語

針對機載LiDAR點云的特點，本文提出了一種基于機載LiDAR點云高程相關(guān)的配準(zhǔn)方法，通過將基準(zhǔn)點云片內(nèi)局部范圍內(nèi)的點擬合曲面作為特征曲面，在待配準(zhǔn)點云中對比求取在擬合曲面的擬合高程，計算擬合高程與實際高程的相關(guān)系數(shù)，選擇待配準(zhǔn)點云中相關(guān)系數(shù)最大位置作為配準(zhǔn)的同名曲面。通過建立機載LiDAR三維點云的局部小范圍擬合曲面，將離散機載三維點云配準(zhǔn)問題簡化為實際高程與擬合高程的相關(guān)性對比問題。實驗表明，本文算法與傳統(tǒng)ICP算法相比,精度相當(dāng)。但是，在實驗過程中也發(fā)現(xiàn)，對于擬合曲面點的選取、擬合半徑的設(shè)定還有值得深入研究的地方，對于不同地形、地物、不同密度的點云設(shè)定不同的擬合半徑和搜索半徑往往得出不一樣的結(jié)論。針對不同地形、地物，不同密度的點云該如何選擇擬合半徑這一問題將是后期研究工作的重點。