吉同元 章定文 毛志遠 徐孝彬 李鵬飛

(1東南大學交通學院，南京 211189)(2華設(shè)設(shè)計集團股份有限公司，南京 210014)(3河海大學機電工程學院，常州 213022)

隨著科技的發(fā)展，聲吶已經(jīng)廣泛應用于水下目標探測、缺陷檢測和地形重構(gòu)[1-5]等，成為認識、開發(fā)和利用水下資源的重要手段.當聲吶用于水下樁基檢測時，一般將聲吶安裝在支架上并放入水中進行掃測，直接獲得水下樁基及周邊地形的空間三維點云數(shù)據(jù)，有效顯示樁基外觀狀態(tài)和損傷情況.但是由于聲吶支架穩(wěn)定性不高，聲吶頭俯仰角變化時支架會發(fā)生晃動，使得不同俯仰角的數(shù)據(jù)拼接區(qū)域誤差較大；并且聲吶設(shè)備獲取的點云存在不同程度的缺陷，如噪聲、空洞和尖銳特征的丟失等.因此，需要提出處理聲吶點云的方法.

在聲吶點云處理方面的關(guān)鍵技術(shù)包括去噪、配準和特征提取等，由于聲吶采集圖像是通過二自由度云臺獲取，因此在同一測量點下有多幅點云數(shù)據(jù)圖需要進行配準.Biber等[6]在2003年第一次提出了正態(tài)分布變換(NDT)算法，主要用于二維平面點云數(shù)據(jù)之間的匹配，該算法根據(jù)參考數(shù)據(jù)來構(gòu)建多維變量的正態(tài)分布，然后用優(yōu)化的方法求出使得概率密度之和最大的變換參數(shù)，從而確定2個點云間的最優(yōu)匹配.Magnusson[7]2009年首次在三維空間中使用NDT算法，并提出了3D-NDT算法，使得NDT算法突破了維度的限制，能夠適用于利用機器人進行采集的點云數(shù)據(jù)配準.經(jīng)過十幾年的發(fā)展，NDT算法已較為成熟，被廣泛應用于點云的配準中.對點云配準后，需進行特征提取，由于樁基是圓柱體，處理點云數(shù)據(jù)時要進行圓柱聚類分割.Schnabel等[8]在2007年提出了基于隨機抽樣一致(RANSAC)方法的形狀三維點云分割方法，從一組包含“局外點”的觀測數(shù)據(jù)集中，通過迭代方式估計數(shù)學模型的參數(shù).但是在實際處理中，由于缺乏對各類點云進行判斷，容易出現(xiàn)過擬合現(xiàn)象.基于區(qū)域分割的方法分為自下而上[9]和自上而下的非種子區(qū)域方法[10]，通過選取種子點，以法向量夾角或曲率作為閾值進行區(qū)域生長.基于聚類的方法廣泛應用于無監(jiān)督分割任務[11]，劉慧等[12]把獲取的植物深度圖像和彩色圖像轉(zhuǎn)換為顯著性點云，通過均值漂移(mean-shift)算法生成超體素,最后利用局部凸連接打包一波帶走(LCCP)算法融合聚類，對樹木實現(xiàn)了精確的枝葉分離.歐式聚類是一種基于歐式聚類度量的聚類算法[13]，其基于KD-Tree 的近鄰查詢算法，計算鄰域點到該點的歐氏距離，將在閾值范圍內(nèi)的點聚為一類，然后反復迭代.模糊c-均值聚類算法[14](FCM)是基于劃分的聚類算法，主要根據(jù)各個對象之間的相似性來完成目標對象的聚類.基于密度的噪點空間聚類(DBSCAN)算法是一種基于密度的空間數(shù)據(jù)聚類，可以發(fā)現(xiàn)某一密度條件下任意形狀的聚類[15-16].李海倫等[17]提出把點云的坐標、法向量、平均曲率和高斯曲率組成8維特征向量,然后使用遺傳算法聯(lián)合模糊聚類算法對汽車零件進行分割,獲得了合理的分類.Qi等[18]提出Point Net++點云分割網(wǎng)絡(luò)框架，并對點云數(shù)據(jù)集添加標簽后用此網(wǎng)絡(luò)進行訓練分割，最終取得了精確的分割結(jié)果.RANSAC、區(qū)域生長和聚類分割算法有各自的優(yōu)缺點，因此研究者通常將3種算法融合起來以得到最好的效果.隨著深度學習的發(fā)展及計算機計算能力的提高，深度學習逐漸成為點云分割領(lǐng)域的研究熱點.除了分割之外，對水下樁基的檢測極為重要，而關(guān)于樁基的缺陷分析和尺寸測量的文獻較少.

本文搭建了聲學可視化檢測平臺，并通過該平臺獲取水下樁基的三維點云數(shù)據(jù).采用NDT匹配方法實現(xiàn)同一站點不同俯仰角下的配準，利用RANSAC方法實現(xiàn)水面擬合，然后采用多種濾波算法濾除雜散點云；提出了基于圓柱特征的分割方法，實現(xiàn)樁基分割與擬合，完成樁基尺寸測量和缺陷檢測.

1 水下聲學可視化檢測系統(tǒng)

本文搭建了水下樁基聲學可視化檢測系統(tǒng)，如圖1所示.將Teledyne BlueView生產(chǎn)的BV5000 S3聲吶系統(tǒng)通過定位銷安裝在云臺上，采用云臺掃描實現(xiàn)水下三維成像，其他組件包含S3接線盒、聲吶云臺纜線、USB B到A纜線以及網(wǎng)線.

圖1 水下樁基聲學可視化檢測平臺

2 水下聲學目標尺寸檢測方法

本文提出了水下樁基聲學可視化檢測方法，該方法由單站點云匹配、水面校準、點云濾波和樁基擬合與分析4部分組成.算法流程如下：①利用NDT配準算法實現(xiàn)同一站點下的不同俯仰角間掃描點云的配準拼接；②對水面進行平面分割及擬合，利用擬合后的水面參數(shù)進行矯正，從而實現(xiàn)對整個點云的矯正；③采用直通濾波、球形區(qū)域濾波、高斯統(tǒng)計濾波及半徑濾波去除其余的雜散點；④通過獲取的點云數(shù)據(jù)進行樁基圓柱分割和擬合，分析尺寸變化.

2.1 基于NDT算法的單站點云配準

在同一站點下，聲吶采集數(shù)據(jù)時采用球面掃描，根據(jù)選定的傾斜角度(45°、15°、-15°、-45°)自動掃描，不同角度的掃描聲吶圖像如圖2所示.由于在同一位置進行測量時，位置偏移較小，旋轉(zhuǎn)角度較大，本文結(jié)合已知的旋轉(zhuǎn)角度，采用NDT算法進行精確配準.采用兩兩點云配準，將第1幅聲吶點云作為基準，不斷迭代兩幀點云獲得變換矩陣，再通過變換矩陣將后續(xù)的聲吶點云疊加到第1幅聲吶點云的基準上.

(a) 45°

NDT配準算法將參考點云網(wǎng)格化，根據(jù)網(wǎng)格內(nèi)的點計算概率密度函數(shù)，并解算多維正態(tài)分布參數(shù)：均值和協(xié)方差矩陣.對于要配準的點云，通過變換將其轉(zhuǎn)換到參考點云的網(wǎng)格中.計算點的概率密度函數(shù)如下：

(1)

式中，μ表示包括點云x在內(nèi)體素網(wǎng)格內(nèi)的均值向量；C表示包括點云x在內(nèi)體素網(wǎng)格內(nèi)的協(xié)方差矩陣；b為常量.

均值和協(xié)方差矩陣，即各體素網(wǎng)格內(nèi)的μ和C可以定義為

(2)

(3)

式中，xi(i=1,2,…,h)為網(wǎng)格中的全部點云.

NDT配準目標函數(shù)由每個網(wǎng)格計算出的概率密度相加得到，在對各映射點概率分布求和過程中，計算s(p)的公式如下：

(4)

式中，x′i為x通過轉(zhuǎn)換矩陣得到的在目標點云坐標系中的位置信息；μ′i表示x′i對應的均值向量；C′表示x′i對應網(wǎng)格的協(xié)方差矩陣.三維變換矩陣T(p,x)可表示為

(5)

式中，初始變換參數(shù)p=[t|r|φ]，t為平移偏移量，r和φ分別為旋轉(zhuǎn)軸和旋轉(zhuǎn)角，t={tx,ty,tz}，r={rx,ry,rz}；s=sinφ；c=cosφ；o=1-c.

采用牛頓優(yōu)化算法對目標函數(shù)進行優(yōu)化，即尋找變換參數(shù)p使得目標函數(shù)達到最小值.令f=s(p)，為了使函數(shù)f最小，對于每次操作都要進行以下方程處理：

HΔp=-g

(6)

式中，Δp為變換參數(shù)的變化量；g為f的轉(zhuǎn)置梯度，其元素可表示為

(7)

其中，pi為變換參數(shù)元素.H為f的Hessian矩陣，其元素可表示為

(8)

2.2 點云平面分割及水平面矯正

由于聲吶系統(tǒng)利用三腳架座底安裝，因此聲吶獲取的三維點云數(shù)據(jù)與z軸不垂直，與水面呈一定的角度.此外，用聲吶獲取的初始點云數(shù)據(jù)如圖3所示，其存在鏡像噪聲，不利于后續(xù)處理，需完成水面點云分割及水平面矯正.

圖3 聲吶獲得的初始點云數(shù)據(jù)

采用RANSAC方法對水平面進行分割，并擬合獲取水平面參數(shù)方程.對于水平面的矯正，先通過獲取兩向量之間的旋轉(zhuǎn)矩陣，再根據(jù)旋轉(zhuǎn)矩陣將點云進行變換，從而矯正點云.RANSAC是根據(jù)一組包含異常數(shù)據(jù)的樣本數(shù)據(jù)集，計算出數(shù)據(jù)的數(shù)學模型參數(shù)，得到有效樣本數(shù)據(jù)的算法.

平面的模型為

Ax+By+Cz+D=0

(9)

采用RANSAC擬合平面后，得到了該平面的法向量，即

M={A,B,C}T

(10)

取水平面的法向量為

Q={q1,q2,q3}T={0,0,1}T

(11)

繞Q旋轉(zhuǎn)角度為β，向量表示為Q=RM，M通過旋轉(zhuǎn)矩陣R旋轉(zhuǎn)到Q.根據(jù)羅德里格旋轉(zhuǎn)公式求出旋轉(zhuǎn)矩陣：

(12)

2.3 點云濾波算法

平面分割及矯正后，三維點云數(shù)據(jù)中還混雜著大量的噪聲點、粗差點和冗余點，會對樁基的真實形態(tài)造成干擾，因此需要采用濾波算法將其濾除.本文首先采用直通濾波器，過濾掉在指定維度方向上取值不在給定值域內(nèi)的點；然后采用球形區(qū)域濾波，在給定球半徑r，得到一個以原點(0,0,0)為圓心、r為半徑的球，選擇剔除球內(nèi)的點和剔除球外的點進行球形區(qū)域濾波；其次，采用基于連通分析的點云半徑濾波保留原始點云中符合假設(shè)條件的點云；最后，采用體素濾波對點云壓縮，減少點云數(shù)量.

2.4 樁基分割與缺陷檢測

為了實現(xiàn)水下樁基的檢測與尺寸測量，采用RANSAC方法提取圓柱形狀.但是在實際使用中，不是圓柱特征的目標也會被提取出來，因此為了準確提取樁基，還需要判斷分割后點云是否為圓柱特征.在分割出多類點云數(shù)據(jù)后，計算每一類點云的曲率和法向量，提出基于圓柱形狀特征的二次判斷方法，獲取多個真正樁基點云.計算點集F={Fi}(i=1,2,…,k)，點Fi的法向量di可以表示為

(13)

曲率δ可以表示為

(14)

式中，λ1、λ2、λ3(λ1≥λ2≥λ3)為點云的特征值.

圓柱特征為：假設(shè)圓柱表面無噪聲，法線協(xié)方差的最小特征值總是零，最大特征值和最小特征值之比無窮大.但是實際情況下點云總是含有噪聲，因此提出采用最大特征值與最小特征值的比值σ作為圓柱判斷條件之一，其公式為

(15)

另外根據(jù)樁基主要特征為垂直地面，其法向量也需要滿足一定的要求，即其法向量的方向與底面呈現(xiàn)一定角度.判斷分割后每一類點云是否為圓柱點云的具體公式如下：

(16)

式中，σmax=150；ath=0.6.a1、a2、a3為點云的特征向量.其中，a3為法向量.

因此，提出基于圓柱特征的RANSAC圓柱提取方法，具體步驟如下：

①先采用RANSAC方法進行圓柱分割得到k類.

②計算點云的法向量a3和最大特征值與最小特征值的比值σ.

③判斷分割后每一類點云是否為圓柱點云，若是，則保存為圓柱點云；否則，返回步驟②.

④分割的k類計算判斷完成后，提取圓柱的原始點云.

⑤自動選取圓柱點云分析.

⑥基于RANSAC擬合圓柱，通過圓柱半徑分析樁基的缺陷.

空間圓柱示意圖如圖4所示，圓柱的幾何特征為：圓柱面上的點到其軸線的距離恒等于半徑r.其中，中軸線方程如圖4所示，P(x,y,z)為圓柱面上任意一點，P0(x0,y0,z0)為圓柱軸線上一點，L={l,m,n}為圓柱軸線向量，r為圓柱底圓半徑，θ為PP0與軸線夾角.因此，根據(jù)圓柱軸線上的某起始點P0、圓柱的軸線法向量、圓柱底面半徑，可以確定圓柱方程如下：

(x-x0)2+(y-y0)2+(z-z0)2-r2=

(17)

圖4 空間圓柱示意圖

3 實驗與分析

混凝土灌注樁由于施工質(zhì)量、樁身材料以及地質(zhì)條件等因素影響，容易出現(xiàn)縮徑和擴徑現(xiàn)象.縮徑和擴徑是指澆筑混凝土后的樁身局部直徑小于或大于設(shè)計尺寸，縮徑現(xiàn)象會導致單樁承載力不足，對上部結(jié)構(gòu)穩(wěn)定性造成影響；擴徑會造成混凝土的浪費和間接產(chǎn)生局部縮徑的現(xiàn)象.因此，針對圖5所示某橋墩樁基縮徑現(xiàn)象和圖6所示某靠船墩樁基擴徑現(xiàn)象進行實驗分析.

圖5 某橋墩樁基縮徑

圖6 某靠船墩樁基擴徑

3.1 縮徑實驗分析

3.1.1 點云配準

聲吶在45°、15°、-15°、-45°這4個角度采得的點云數(shù)據(jù)如圖7所示.從圖中可知，4個點云數(shù)據(jù)特征明顯不同，在已知旋轉(zhuǎn)角度情況下結(jié)合NDT匹配算法效果最佳，設(shè)置NDT網(wǎng)格為2 m，線搜索步長為0.1，最小轉(zhuǎn)換差異為0.001，配準后的結(jié)果如圖8所示.由配準結(jié)果可知，當使用NDT算法配準時，點云模型多處樁基部位的配準結(jié)果都較好，歐幾里得適合度得分為0.013 676 3，精度較高.由于在配準過程中未采用特征點匹配計算，因此時間也較快，配準時間僅為15 s.

(a) 45°

(a) 俯視圖

3.1.2 點云平面分割及水平面矯正

水面校準前點云如圖9(a)所示，由于聲吶坐底傾斜等原因?qū)е曼c云數(shù)據(jù)與z軸不垂直，水面的點云相較于水底點云更平坦且變形較少，可用作參考平面.通過RANSAC擬合水面點云，對NDT配準后的點云數(shù)據(jù)進行校準，效果如圖9(b)所示，水下目標物到水面的距離得到修正，也能真實反映目標物情況，通過水面擬合出的平面方程為

-0.054x-0.102y+0.993z-2.126=0

(18)

該平面法向量為n1={-0.054,-0.102,0.993}T，將n1與水平面法向量n2={0,0,1}T代入下式：

(19)

解算出的校準度數(shù)α=6.630°，然后計算其旋轉(zhuǎn)矩陣R1為

(20)

運算時間僅為2 s，為后續(xù)樁基處理和分析節(jié)省大量時間.

(a) 校準前點云數(shù)據(jù)

3.1.3 點云濾波

由于校準后的點云數(shù)據(jù)含有各類雜散點，因此需對其進行濾波.首先對校準后的點云數(shù)據(jù)進行直通濾波，截取z值在-1.5～1.9之間的點云，濾波時間僅需1 s，效果如圖10(a)所示.然后進行球外濾波，選擇球的球心為(0,0,0)，半徑為2.5；再進行球內(nèi)濾波，選擇球的球心為(0,-1,1)，半徑為1.5，此過程耗時48 s，效果如圖10(b)所示.其次進行半徑濾波，選擇半徑r為0.1，近鄰個數(shù)為100，耗時40 s，效果如圖10(c)所示.最后采用體素網(wǎng)格濾波方法，設(shè)置體素網(wǎng)格的濾波網(wǎng)格大小為0.02 m×0.02 m×0.02 m，該過程較快，耗時1 s，得到濾波后的結(jié)果如圖10(d)所示.與濾波前作對比，濾波前點云數(shù)量為2 011 143，濾波后點云數(shù)量為43 144，圖中的離散點被移除，樁基的輪廓也更加清晰.

(a) 直通濾波效果

3.1.4 樁基分割

通過RANSAC擬合分割后的目標點云如圖11所示，不同顏色顯示為不同類別，分為6類.從圖中可明顯看到有非圓柱體目標被認為是圓柱體，這是過擬合現(xiàn)象，因此需要結(jié)合提出的基于圓柱特征判斷準則進行再次判斷.

圖11 RANSAC方法縮徑實驗點云分割結(jié)果

為了提取出有效樁基點云，采用提出的圓柱特征判斷準則，對點云進行二次判斷，其結(jié)果如表1所示.如果只采用曲率特征的話，滿足條件的點云類別有2、3、4、5，再增加法向量判斷準則后剩下點云類別3、4、5，結(jié)合圖11可知，這3類是對應圓柱點云.因此，證明了提出的基于圓柱特征判斷準則的RANSAC圓柱點云提取方法有效性.

表1 縮徑實驗點云曲率與法向量特征值

為了驗證算法先進性，采用歐式聚類和區(qū)域生長算法進行分析，其結(jié)果如圖12所示.從圖中可知，不同顏色表示不同類別，歐式聚類分割算法和區(qū)域生長算法無法有效提取圓柱特征，同時分割出的類別過多，無法應用在水下聲吶點云處理，因此本文提出的算法效果最佳.

(a) 本文算法

完成圓柱擬合后再進行圓柱填充，得到填充結(jié)果如圖13所示.樁基分割及圓柱擬合速度快，整個過程僅需1 s，由分割結(jié)果可知該樁基上半部分參數(shù)x=1.756，y=0.823，z=0.534，l=0.037，m=0.012，n=0.999，r=0.768，下半部分參數(shù)x=1.709，y=0.874，z=-0.929，l=0.022，m=-0.011，n=1，r=0.594.擬合得到的樁基上半部分直徑為1.536 m，該樁基的實測值為1.5 m，誤差為0.036 m，精度為97.6%.擬合得到樁基的下半部分直徑為1.188 m，存在嚴重的縮徑病害，應及時采取有效措施.

(a) 俯視圖

3.2 擴徑實驗分析

3.2.1 點云配準

聲吶在45°、15°、-15°、-45°這4個角度采得的點云數(shù)據(jù)如圖14所示.從圖中可知，4個點云數(shù)據(jù)特征明顯不同，在已知旋轉(zhuǎn)角度情況下結(jié)合NDT匹配算法效果最佳，設(shè)置NDT網(wǎng)格為2 m，線搜索步長為0.01，最小轉(zhuǎn)換差異為0.001，配準后的結(jié)果如圖15所示.由配準結(jié)果可知，當使用NDT算法配準時，點云模型多處部位如樁基、船體的配準結(jié)果都較好，歐幾里得適合度得分為0.278 975，精度較高.由于在配準過程中未采用特征點匹配計算，因此時間也較快，配準時間僅為17 s.

(a) 45°

3.2.2 點云平面分割及水平面矯正

水面校準前點云如圖16(a)所示，由于聲吶坐底傾斜等原因?qū)е曼c云數(shù)據(jù)與z軸不垂直，水面的點云相較于水底點云更平坦且變形較少，可用作參考平面.通過RANSAC擬合水面點云，對NDT配準后的點云數(shù)據(jù)進行校準，效果如圖16(b)所示，水下目標物到水面的距離得到修正，也能真實反映目標物情況，通過水面擬合出的平面方程為

0.049x+0.064y+0.997z-5.398=0

該平面法向量為n3={0.049,0.064,0.997}T，解算出的校準度數(shù)為4.622°，然后計算其旋轉(zhuǎn)矩陣R2為

(21)

運算時間僅為2 s.

(a) 校準前點云數(shù)據(jù)

3.2.3 點云濾波

由于校準后的點云數(shù)據(jù)含有各類雜散點，因此需對其進行濾波.首先對校準后的點云數(shù)據(jù)進行直通濾波，截取z值在-3～1.7之間的點云，濾波時間僅需1 s，效果如圖17(a)所示.然后進行球外濾波，選擇球的球心為(0,0,0)，半徑為6；再進行球內(nèi)濾波，選擇球的球心為(0,0,0)，半徑為3.5，此過程耗時58 s，效果如圖17(b)所示.其次進行半徑濾波，選擇半徑r為0.1，近鄰個數(shù)為50，耗時26 s，效果如圖17(c)所示.最后采用體素網(wǎng)格濾波方法，設(shè)置體素網(wǎng)格的濾波網(wǎng)格大小為0.05 m×0.05 m×0.05 m，該過程較快，耗時1 s，得到濾波后的結(jié)果如圖17(d)所示.與濾波前作對比，濾波前點云數(shù)量為2 675 132，濾波后點云數(shù)量為34 591，圖中的離散點被移除，樁基的輪廓也更加清晰.

(a) 直通濾波效果

3.2.4 樁基分割

通過RANSAC擬合分割后的目標點云如圖18所示，不同顏色顯示為不同類別，分為5類.結(jié)合提出的基于圓柱特征判斷準則進行再次判斷.

為了提取出有效樁基點云，本文采用提出的圓柱特征判斷準則，對點云進行二次判斷，其結(jié)果如表2所示，采用曲率特征和法向量進行判斷，并結(jié)合圖18可知，點云類別4、5這2類是對應圓柱點云，證明了本文提出的基于圓柱特征判斷準則的RANSAC圓柱點云提取方法的有效性.

表2 擴徑實驗點云曲率與法向量特征值

為了驗證本文算法先進性，采用歐式聚類和區(qū)域生長算法進行分析，其結(jié)果如圖19所示.歐式聚類分割和區(qū)域生長算法無法有效提取圓柱特征，同時分割出類別過多，無法應用在水下聲吶點云處理，因此本文算法效果最佳.

(a) 本文算法

進行圓柱擬合后再進行圓柱填充，得到填充結(jié)果如圖20所示.樁基分割及圓柱擬合速度快，整個過程僅需1 s，由分割結(jié)果可知該樁基上半部分參數(shù)x=-3.574，y=3.704，z=-1.323，l=0.11，m=0.162，n=0.981，r=0.526；下半部分參數(shù)x=-3.571，y=3.412，z=-0.827，l=0.167，m=0.019，n=0.986，r=1.127.擬合得到的樁基上半部分直徑為1.052 m，該樁基的直徑實測值為1 m，測量誤差為0.052 m，精度為 94.8%.擬合得到樁基的下半部分直徑為2.254 m，存在明顯的擴徑現(xiàn)象.

(a) 俯視圖

4 結(jié)論

1)設(shè)計了水下樁基聲學可視化檢測系統(tǒng)，實現(xiàn)水下成像，解決了水下樁基三維量化數(shù)據(jù)難以獲取的問題.

2)提出了基于水面校準和圓柱特征分割的水下樁基尺寸檢測方法.將提出的算法與歐式聚類算法和區(qū)域生長算法進行對比，本文算法能夠精確地分割出樁基點云數(shù)據(jù)，證明了提出的分割算法的有效性及準確性.

3)縮徑實驗中樁徑擬合值與實測值誤差為 0.036 m，精度達97.6%；擴徑實驗中樁徑擬合值與實測值誤差為0.052 m，精度達94.8%.