孔振興 鄒進(jìn)貴 邊鴻明

1 武漢大學(xué)測(cè)繪學(xué)院，湖北 武漢，430079

2 浙江省國土勘測(cè)規(guī)劃有限公司，浙江 杭州，310000

利用三維激光掃描技術(shù)能快速、高效地獲取物體表面的三維坐標(biāo)信息，被廣泛應(yīng)用于地鐵隧道形變監(jiān)測(cè)、滑坡監(jiān)測(cè)、文物保護(hù)、建筑物三維重建、自動(dòng)駕駛等領(lǐng)域［1?5］。在大量無規(guī)則點(diǎn)云數(shù)據(jù)中快速檢測(cè)出其中包含的基本形狀，有利于對(duì)物體進(jìn)行快速分類，可用于物體的非監(jiān)督分類［6］。現(xiàn)實(shí)中有很多物體包含圓柱體形，如何快速得到圓柱體的模型參數(shù)并將其用于三維重建或分類是一個(gè)熱門問題。利用幾何模型方法，用圓柱模型的數(shù)學(xué)表達(dá)式建立誤差方程，進(jìn)而求得圓柱模型參數(shù)［7?9］，不太適用于大數(shù)據(jù)量的點(diǎn)云模型檢測(cè)。利用圓柱體的高斯特性可將提取過程分為兩個(gè)階段進(jìn)行特征參數(shù)提取［10］；采用遺傳算法可以實(shí)現(xiàn)圓柱幾何特征信息全局最優(yōu)解的評(píng)定［11］；利用隨機(jī)抽樣一致（random sample consensus，RANSAC）算法擬合圓柱面，在點(diǎn)云中將兩點(diǎn)作為種子點(diǎn)進(jìn)行圓柱面的檢測(cè)，不能保證選取的兩個(gè)點(diǎn)的質(zhì)量，影響計(jì)算得到的模型參數(shù)［12］?；赗ANSAC的算法得到的模型參數(shù)有很好的魯棒性，能從包含大量局外點(diǎn)的數(shù)據(jù)集中估計(jì)出較高精度的參數(shù)?？紤]到點(diǎn)云數(shù)據(jù)量大且其中可能包含大量噪聲點(diǎn)，本文提出了一種基于改進(jìn)選點(diǎn)策略的RANSAC圓柱面檢測(cè)算法。

1 傳統(tǒng)基于RANSAC的點(diǎn)云圓柱面檢測(cè)

1.1 k維樹

k維樹（k?dimensional tree，KD?Tree）是k維的二叉樹索引，主要用于檢索多屬性的數(shù)據(jù)或多維數(shù)據(jù)。與二叉樹不同的是，KD?Tree的每個(gè)結(jié)點(diǎn)均表示k維空間的點(diǎn)。如當(dāng)k=3時(shí)，每個(gè)結(jié)點(diǎn)表示一個(gè)三維坐標(biāo)（X，Y，Z）。KD?Tree能像二叉樹一樣，按照一定的索引規(guī)則完成點(diǎn)的精確查找，由每個(gè)內(nèi)部結(jié)點(diǎn)決定搜索走向，最終搜索到所查詢的結(jié)點(diǎn)的塊。索引結(jié)構(gòu)中相似性查詢有兩種基本的方式：①R半徑查詢：通過提供查詢點(diǎn)、查詢距離的閾值，得到一個(gè)數(shù)據(jù)集合，其中每個(gè)數(shù)據(jù)都小于給定的閾值；②k最近鄰（k?nearest neighbor，KNN）個(gè)數(shù)查詢：通過KNN得到一個(gè)包含K個(gè)數(shù)據(jù)的，距離查詢點(diǎn)最近的點(diǎn)的集合［13］。

1.2 點(diǎn)云法向量估計(jì)

點(diǎn)云法向量估計(jì)最常用的方法是協(xié)方差分析法，對(duì)于點(diǎn)云中的點(diǎn)pi及其K鄰域點(diǎn)集Kn(pi)（鄰域點(diǎn)可通過KD?Tree空間查詢得到），構(gòu)造Kn(pi)的協(xié)方差矩陣Ci：

式中，Ci為對(duì)稱正定矩陣，其特征值為（λ0，λ1，λ2），對(duì)應(yīng)的3個(gè)特征向量（e0，e1，e2）相互正交，最小特征值λ0所對(duì)應(yīng)的特征向量e0即為最優(yōu)平面的法向量，可作為pi的法向量［14，15］。

1.3 圓柱面擬合

為了用任意兩點(diǎn)及其法線生成一個(gè)圓柱面，首先確定軸的方向a(l，m，n)，a=n1×n2，其中，n1、n2分別為圓柱體表面的任意兩點(diǎn)p1、p2的法向量。將兩條法線p1+tn1和p2+tn2沿軸線方向投影到a x=0平面上，兩條投影線相交得到點(diǎn)p12。由此可以得到中軸線的直線表達(dá)式。應(yīng)用閾值、點(diǎn)到圓柱體軸線的距離來判斷局內(nèi)點(diǎn)。統(tǒng)計(jì)滿足該模型參數(shù)的點(diǎn)云個(gè)數(shù)，在滿足局內(nèi)點(diǎn)個(gè)數(shù)閾值的圓柱體中尋找標(biāo)準(zhǔn)差最小的圓柱體，得到最佳擬合圓柱體。圓柱面示意圖見圖1。

圖1 圓柱面示意圖Fig.1 Diagram of Cylindrical Surface

2 改進(jìn)選點(diǎn)策略的RANSAC圓柱面檢測(cè)

2.1 改進(jìn)的種子點(diǎn)選取方式

傳統(tǒng)基于RANSAC算法的圓柱面參數(shù)通過空間中兩個(gè)點(diǎn)來確定，該算法的選點(diǎn)方式是從原始數(shù)據(jù)點(diǎn)中隨機(jī)選出兩個(gè)點(diǎn)作為種子點(diǎn)，獲取圓柱面模型參數(shù)的初始值，再根據(jù)初始值尋找數(shù)據(jù)中其他的模型內(nèi)點(diǎn)。該方法得到的參數(shù)模型穩(wěn)健性較差，受點(diǎn)云數(shù)據(jù)質(zhì)量的影響較大。在采樣次數(shù)一致情況下，滿足判斷準(zhǔn)則的圓柱面模型的數(shù)據(jù)集減少了，降低了得到最優(yōu)模型的概率。本文提出的改進(jìn)種子點(diǎn)選取方式是隨機(jī)選取3個(gè)點(diǎn)，通過KD?Tree建立的索引尋找每個(gè)點(diǎn)R半徑內(nèi)的鄰域點(diǎn)集Kn，采用協(xié)方差分析法來求得這3個(gè)點(diǎn)的法向量，并根據(jù)判斷準(zhǔn)則得到模型參數(shù)的初始值。

2.2 模型參數(shù)判斷準(zhǔn)則

依據(jù)兩點(diǎn)作為種子點(diǎn)方法得到的圓柱體參數(shù)模型絕大多數(shù)不滿足判斷準(zhǔn)則，考慮在模型參數(shù)初始值的選取上利用3個(gè)點(diǎn)的法向量來構(gòu)造初始值：

式中，n1、n2、n3表示過3點(diǎn)的法向量；a1、a2是由法向量n1、n2、n3構(gòu)造的軸線方向。對(duì)于規(guī)則光滑的圓柱體表面，理論上a1=λ×a2，考慮到獲取的圓柱體表面點(diǎn)云數(shù)據(jù)存在噪聲，引入兩條軸線的夾角這一判斷準(zhǔn)則：若θ<θ0（角度閾值），則軸線方向a=a1+a2，a為兩個(gè)向量的合向量。這樣可以抑制部分噪聲，優(yōu)化模型初始值的計(jì)算。

將3條法線p1+tn1、p2+tn2、p3+tn3沿軸線方向投影在a x=0平面上，理論上3條法線在該平面上的投影應(yīng)交于一點(diǎn)，實(shí)際得到兩兩相交的3個(gè)點(diǎn)。將得到的3點(diǎn)組成一個(gè)三角形，若三邊d12、d13、d23的距離均小于閾值d，則求三角形的重心c，否則重新選取3個(gè)點(diǎn)，得到新的模型參數(shù)。

設(shè)3個(gè)點(diǎn)在投影平面的投影點(diǎn)坐標(biāo)分別為（x1，y1，z1）、（x2，y2，z2）、（x3，y3，z3），則重心c的坐標(biāo)（Xc，Yc，Zc）計(jì)算公式如下：

2.3 算法流程

本文方法要對(duì)點(diǎn)云數(shù)據(jù)進(jìn)行多次抽樣，在抽樣結(jié)束后比較θ和dmax，將θ最小和dmax最大的那組數(shù)據(jù)得到的模型參數(shù)作為最佳模型參數(shù)。具體算法流程見圖2。

圖2 本文方法流程Fig.2 Flow Chart of the Proposed Method

3 仿真實(shí)驗(yàn)

以真實(shí)值為半徑（R=2 m），中軸線方向?yàn)長(0，0，1)，過點(diǎn)P（0，0，0）構(gòu)造的圓柱面點(diǎn)云個(gè)數(shù)CP為2 000、5 000、10 000、20 000、50 000、100 000、200 000的點(diǎn)云數(shù)據(jù)，并加入一定比例S的噪聲點(diǎn)（S=N/CP），N表示加入的噪聲點(diǎn)個(gè)數(shù)。圖3和圖4分別展示了未加入噪聲點(diǎn)和加入噪聲點(diǎn)后的完整圓柱面點(diǎn)云數(shù)據(jù)，考慮實(shí)際獲取的圓柱面點(diǎn)云數(shù)據(jù)可能只是整個(gè)圓柱面的一部分，本文構(gòu)造了圖5所示的點(diǎn)云數(shù)據(jù)。

圖3 未加入噪聲點(diǎn)的完整柱面Fig.3 Complete Cylindrical Surface Without Noise Points

圖4 加入噪聲點(diǎn)的完整柱面Fig.4 Complete Cylindrical Surface with Noise Points

圖5 加入噪聲點(diǎn)后的部分圓柱面點(diǎn)云數(shù)據(jù)Fig.5 Point Cloud Data of Partial Cylindrical Surface After Adding Noise Points

用基于RANSAC擬合圓柱面的方法和基于改進(jìn)選點(diǎn)策略的RANSAC圓柱面檢測(cè)方法分別處理仿真得到的點(diǎn)云數(shù)據(jù)，當(dāng)S<0.50時(shí)，得到了比較一致的實(shí)驗(yàn)結(jié)果。本文給出了當(dāng)S=0.21時(shí)，兩種方法處理得到的模型參數(shù)和處理所需時(shí)間，見表1。

表1 模型參數(shù)和耗時(shí)Tab.1 Model Parameters and Time Consumption

為了更直觀地比較本文方法和傳統(tǒng)基于RANSAC擬合圓柱面方法的效果和效率，本文利用Python的matplotlib包給出了不同點(diǎn)云數(shù)量下的處理時(shí)間，見圖6。不同點(diǎn)云數(shù)量下的擬合半徑圖見圖7。

圖6 耗時(shí)Fig.6 Time Consumption

圖7 擬合半徑Fig.7 Fitting Radius

在不同點(diǎn)云數(shù)量下，本文方法比傳統(tǒng)基于RANSAC的方法耗時(shí)要少，在點(diǎn)云數(shù)量較少時(shí)，兩者差距很小，但隨著點(diǎn)云數(shù)量的增加，兩者的檢測(cè)速度差距也逐漸增大，當(dāng)點(diǎn)云數(shù)據(jù)中包含多個(gè)需要檢測(cè)的圓柱面點(diǎn)云數(shù)據(jù)時(shí)，兩者的效率差距會(huì)更明顯。

存在噪聲時(shí)，利用兩種方法都可以擬合出圓柱體的半徑，且擬合出來的半徑大小具有較好的一致性，這說明兩種方法受噪聲的影響是基本一致的。但是傳統(tǒng)基于RANSAC擬合圓柱面半徑的方法不夠穩(wěn)健，在點(diǎn)云數(shù)量為12 100時(shí)，擬合出來的半徑和實(shí)際半徑相差34.4 mm，在點(diǎn)云數(shù)量為120 500時(shí)，擬合出來的半徑和實(shí)際半徑僅相差1.4 mm。同樣，本文方法在點(diǎn)云數(shù)量為12 100時(shí)擬合的半徑和實(shí)際半徑相差最大，這說明噪聲點(diǎn)對(duì)利用兩種方法得到的半徑都產(chǎn)生了影響，但本文方法穩(wěn)定性更好，在不同點(diǎn)云數(shù)量下都可以較好地接近真實(shí)值。

4 結(jié)束語

本文提出了一種基于改進(jìn)選點(diǎn)策略的RANSAC圓柱面檢測(cè)方法，在采樣次數(shù)一致的情況下，根據(jù)給定規(guī)則，提高了得到最優(yōu)模型參數(shù)的概率和效率。將隨機(jī)獲取的點(diǎn)云數(shù)據(jù)按照給定規(guī)則進(jìn)行篩選，進(jìn)而當(dāng)作種子點(diǎn)，用于計(jì)算模型參數(shù)，能有效地提高模型參數(shù)初始值的質(zhì)量。實(shí)驗(yàn)結(jié)果表明，該方法能提高模型檢測(cè)的速度和模型參數(shù)的質(zhì)量。