針對橢圓形隧道底部的橫斷面分割及濾波

牛曉潔 季 民 徐 飛

(1. 山東科技大學 測繪與空間信息學院, 山東 青島 266590;2. 福州市勘測院, 福建 福州 350000)

0 引言

隧道是國家現代化建設的重要組成部分之一,由于隧道在施工及運營過程中,受附近林木、山體、人造建筑物等的影響會產生一定的變形,進而對隧道安全產生不利影響,因此對隧道進行定期檢查、測量、維修是隧道工程施工運營中的必不缺少的工作。隧道斷面變形的傳統方法主要使用全站儀、斷面儀、經緯儀等工具,雖有較高的測量精度,但測量工作時間長、斷面數據少、測量效率低、且隧道內部環境惡劣時測量人員的危險性較高。近二十年來地面激光掃描技術得到了迅猛發展,地面激光掃描技術被廣泛應用于橋梁、古建筑、城市街景等的三維測量工作中。相比于傳統測量手段,地面激光掃描技術具有高效率、高斷面點密度、斷面數據量多、測量范圍大、成果用途廣泛的特點。

目前,針對三維激光點云的隧道橫斷面提取、變形分析的方法已經有了較為充足的發展,國內外諸多學者在該方面進行了深入研究。

托雷等人[1]提出一種雙向投影提取中軸線的方法,即將圓形地鐵隧道點云旋轉使隧道法向與Y軸正交,然后將點云分別投影至XOY、YOZ平面結合隨機抽樣一致算法(random sample consensus,RANSAC)提取出兩平面內的中線,然后通過兩條平面曲線來表述空間曲線從而對隧道橫斷面進行提取。朱寧寧等人[2]使用雙向投影的方式提取出隧道中軸線,然后依據中軸線提取隧道橫斷面,最后使用橢圓模型對隧道進行橫斷面進行擬合及濾波。程效軍等人[3]通過對隧道點云雙向投影獲取隧道在XOY、YOZ平面的姿態變化,使用多項式方程對兩條平面曲線進行擬合并插值,然后計算各中軸點的法平面以對隧道點云進行分割,最后計算各切片內的點云到中軸線的距離,使用給定的距離閾值以對隧道內部的點云噪聲進行濾除。盧小平等人[4]使用雙向投影的方法提取出隧道中軸線,然后依據中軸線提取隧道橫斷面切片,接著對橫斷面點云進行旋轉使法線平行于Y軸,然后使用橢圓柱模型對橫斷面進行擬合及濾波。藍秋萍等人[5]在隧道點云上建立高斯球,選取高斯球上的最大圓作為隧道橫斷面,使用最小二乘法對大圓進行擬合將最大圓所在的平面的法向量作為隧道中線;然后根據隧道中線提取橫斷面并提取特征線。陳林等人[6]針對盾構隧道首先使用基于距離差特征的方法識別并提取出環片連接處的環縫點云,接著對環縫點云進行最小二乘擬合;然后將環片兩端的環縫曲線中心連接作為隧道中軸線。杜榮武等人[7]使用移動三維激光掃描系統對地鐵進行掃描工作,根據以隧道軌跡作為基準對隧道櫻花斷面進行截取,然后對橫斷面離散點云進行擬合并與隧道設計半徑相比較從而對隧道管片錯臺、凈空收斂等進行計算。徐飛等人[8]根據隧道邊界特征提取水平邊界并提取出水平中線并以此提取出隧道橫斷面,然后對橫斷面進行擬合及濾波,最后結合隧道設計參數對隧道變形量進行可視化分析。趙亞波等人[9]使用移動三維激光掃描系統對盾構隧道進行測量作業,使用三維激光點云生成正射影像,根據影像中的環片點云信息計算每一環片的錯臺量。林景峰等人[10]使用Delaunay三角網提取隧道水平邊界進而提取水平中線并以此提取橫斷面,然后擬合出空間中軸線。王鑫等人[11]使用三維激光掃描儀對對盾構隧道進行掃描作業,提取出隧道中線及橫斷面并與全站儀數據進行對比分析驗證了其精度。

以上文獻對隧道的橫斷面提取進行了較為深入的研究,然而目前的學者專家主要是使用規則的圓形、橢圓形、圓柱形模型對隧道進行研究和分析;然而除一些盾構隧道外,隧道并非絕對圓形,且隧道襯砌施工后往往隧道底部已不再屬于圓形,如圖1所示,如不對隧道底部進行精細化直接進行橢圓、圓柱擬合的話誤差往往較大,因此對隧道底部點云進行精細處理是有必要的。本文提出一種針對隧道底部為矩形的隧道橫斷面提取及濾波方法,該方法首先將隧道點云進行XOY平面投影從而提取出水平中線,接著使用水平中線提取出隧道點云切片,然后使用RANSAC算法結合平面約束模型對隧道底部點云進行濾波分離,從而實現半圓形半矩形隧道橫斷面的圓形部分與矩形部分的分割,最后對圓形隧道橫斷面的擬合及濾波。

圖1 隧道點云

1 理論方法

1.1 隧道腰部中線提取

本文提出一種先依據隧道底部特征提取隧道中軸線的方法,如圖2所示,隧道底部平面與隧道壁面往往形成略大于90°的夾角α、約等于90°的β、γ。而在隧道頂部及中部位置,則不存在這一特征,因此可以利用這一特征提取隧道底部邊界線。

圖2 隧道底部示意圖

首先將隧道點云進行旋轉使軸向方向與Y軸平行,如圖3所示。

圖3 點云旋轉俯視圖

旋轉公式如式(1)所示,式中(X,Y,Z)1為擬合坐標系中的點的三維坐標,(X,Y,Z)0為原始坐標系中的點的三維坐標,Rω為旋轉矩陣。

(1)

隨機選取種子點,以種子點為原點,Y軸方向為基準,建立局部空間直角坐標系,如圖4所示,此時以種子點為圓心建立kd樹,計算種子點R鄰域內的點云在XOZ平面內所形成的夾角,若該角度約等于閾值α,則該點作為合格點提取出來,若該角度與閾值α相差較多,則認為該點為非合格點。遍歷所有點云,提取出所有合格點即可提取邊界線。

圖4 隧道腰部邊界線提取示意圖

獲取邊界AB與CD的點云之后,即可獲取AB與CD的中線EF,此時的EF與隧道走向相同,因此可將EF作為隧道中軸線并使用多項式方程(2)進行最小二乘擬合,作為提取隧道橫斷面的基準線，

(2)

式中，a為二次項系數；b為一次項系數；c為常數；x為多項式自變量。

1.2 隧道橫斷面提取

橫斷面是隧道形變分析、超欠挖分析、開挖控制、混凝土厚度控制、侵界檢測、錯臺檢測、斷面間體積計算等作業的重要參考,因此對隧道橫斷面進行提取是隧道工程中必不可少的工作。

獲得隧道中軸線后可依據中軸線為基準對隧道橫斷面進行分割。如圖5所示,以中軸線上一點G為基準建立切向量τ,再建立切向量法平面ψ作為切割隧道的基準線,將ψ前后厚度為d/2的點云進行切割從而獲得隧道點云切片,d為切片厚度。

圖5 隧道切割示意圖

獲取點云切片后,由于隧道底部點云呈平面特征,上半部點云呈圓柱特征,因此使用RANSAC算法對平面進行檢核,從而對隧道底部點云進行分割。傳統RANSAC檢測點云平面時首先從點云切片中隨機選取三個點并進行平面模型擬合;然后計算其余點到平面模型的距離li與設定閾值λ進行比較,若li小于閾值,則認為該距離對應的點為平面點,反之則認為是非平面點;設置迭代次數n,對所有點云進行隨機抽取并迭代計算,選取出最佳模型作為最終平面模型。

若直接使用傳統RANSAC方法進行檢測,則平面模型的選取僅依靠隨機選取的點云與平面厚度即閾值λ來決定;這就導致在平面檢核過程中檢核出的平面并非底部平面點云,而是符合最佳模型的“類平面點云”。

為避免這一現象,本文加入平面法向量約束對RANSAC算法進行改進,具體步驟如下：

(1)從點云切片中隨機選取三個點構建函數模型為

(3)

(2)計算各點li與設定閾值λ進行比較,判斷模型優劣;

(3)使用最小二乘原理對閾值λ內的模型點云進行擬合,求取法向量τ(υ,ω,ξ)；

(4)計算各點云至平面函數模型Ax+By+Cz+D=0的向量(ai,bi,ci);

(5)對比計算各向量(ai,bi,ci)與法向量(A′,B′,C′)的夾角μi,若μi的數值大于設定閾值μ0,則該平面模型內的點云為符合最佳模型的“類平面點云”;反之則該平面模型內的點云為底部平面點云。以此實現對隧道底部點云的分割。

獲得隧道切片后(圖6)，將切片內的點云Gi(X,Y,Z)向平面ψ進行投影,得到投影點Qi(X′,Y′,Z′),從而獲得隧道橫斷面輪廓線。

圖6 隧道切片投影示意圖

1.3 橫斷面擬合及濾波

獲得剔除底部點云的隧道橫斷面輪廓線后,即可使用橢圓模型對隧道橫斷面進行擬合。本文使用間接平差對對橢圓標準方程進行擬合。首先使用對橫斷面進行坐標旋轉使橫斷面輪廓線平行于XOZ平面,將三維點云轉化為二維點云從而使用橢圓模型進行擬合。具體步驟如下：

(1)建立間接平差函數模型為

(4)

式中，(x0,y0)，a，b分別代表橢圓幾何中心，長半軸，短半軸；xi,yi分別表示點云在平面ψ的橫、縱坐標。

(2)轉化為誤差方程式為

(5)

根據泰勒公式將方程線性化為

(6)

式中

式中，(x0,y0)，a0，b0分別為參數橢圓幾何中心坐標、長半軸、短半軸的初值；Δx，Δy，Δa，Δb為各參數增量。將式(5)轉化為誤差方程的矩陣形式為

(9)

式中

(10)

根據最小二乘原理,將二倍中誤差2σ作為橫斷面濾波的參數,即若橫斷面點至擬合橢圓曲線的誤差vi>2σ時將該點作為噪點進行濾除;反之,則將該點作為合格點保留并重新進行橢圓曲線的最小二乘擬合,直至橢圓參數的改正值Δx,Δy,Δa,Δb趨于收斂。然后將濾波后的點云逆向旋轉至原始坐標系中,旋轉公式為

(11)

3 實驗分析

使用FARO X130三維掃描儀對某隧道進行掃描作業,隧道內壁半徑RS為5.5 m。掃描儀垂直分辨率為0.009°,水平分辨率為0.009°,掃描最遠距離為130 m。選取130 m隧道作為實驗區,得到隧道點云1 900 553個。在VC++環境下使用文中方法對隧道進行橫斷面提取及擬合。

3.1 隧道邊界及中線提取

首先對隧道進行預處理得到符合實驗要求的隧道點云,接著根據隧道底部邊界特征對隧道底部邊界進行提取,如圖7所示,共得到邊界點云8 153個。

圖7 隧道腰部邊界及中線

接著遍歷兩條邊界點,依據中點公式(12)求取中線點,式中(xleft,yleft,zleft)、(xright,yright,zright)分別表示兩邊的邊界點坐標,(xmid,ymid,zmid)表示中線點坐標。

(12)

求取中線如圖7所示,得到中線共有4 030個點,然后使用多項式方程對中線進行最小二乘擬合。由于該中線為三維曲線難以檢核,因此本文使用設計樁點的二維數據進行精度檢核,獲得中線對比數據如表1所示。

表1 中線坐標對比

從表1中可觀察到ΔX、ΔY最大值為11 mm、最小值分別為5 mm和2 mm,ΔD最大值為14.9 mm、最小值為8.6 mm,中線提取精度較高。擬合后的中線即可作為提取隧道橫斷面的基準線。

3.2 隧道橫斷面分割及擬合

3.2.1 橫斷面分割

提取并擬合出隧道中線后,以此為基準構建橫斷面方程,切割橫斷面切片,本文設置橫斷面切片厚度為0.15 m,共有650個橫斷面切片,切片厚度與中線點間隔一致為0.15 m。由于橫斷面切片互相連接不宜展示,為便于展示,本文實驗部分采用間隔4 m,切片厚度為0.11 m的橫斷面數據進行展示,如圖8所示為隧道切片。

圖8 橫斷面切片

獲得點云切片后,傳統方法往往直接使用RANSAC算法結合平面模型Ax+By+Cz+D=0對隧道底部點云進行分割,分割情況如下：

(1)設置迭代次數n為500,λ為0.01 m,如圖9(a)所示為提取到的最優平面模型內的點云,可以觀察到由于底面點云并非在同一平面內,而是有一定厚度的點云面片,因此當λ較小時RANSAC算法提取到的最優平面點云只是底面。

(2)設置迭代次數n為500,λ為0.2 m;如圖9(b)所示為提取到的最優平面模型內的點云,可以觀察到當λ較大時RANSAC算法提取到的最優平面點云為部分切片,而不是底部平面點云。

(3)設置迭代次數n為500,λ為0.06 m,如圖9(c)所示為提取到的最優平面模型內的點云,可以觀察到RANSAC算法提取到的最優平面點云為底面點云。

(a) λ為0.01 m (b) λ為0.2 m

(c) λ為0.6 m圖9 橫斷面分割

從圖9中可以看到直接使用平面模型進行選取最優參數時,計算機經常會選擇出符合平面模型參數的非底面點云,這是因為雖然該部分并非平面或者并非全部底面,但在數學模型內卻包含最多點云。此外測量作業中還會出現點云缺失、漏測、點云密度不均、點云解算軟件的內插操作等各種情況,往往造成平面模型的不確定性,進而使得選取參數時產生不確定性。因此就需要工作人員不斷調節誤差閾值λ,以使RANSAC模型提取出的是平面點云;然而這就使得工作效率較低,且可能出現同樣閾值λ在部分切片分割出了平面點云,部分切片分割出了非平面點云,進而需要對不同切片使用不同閾值λ,這就進一步增加了操作的復雜性。

針對該種情況,本文加入法向量τ(υ,ω,ξ)約束,利用法向量自動選取出符合底面特征的最優平面模型。

由圖7中的中線可求得隧道坡度在15°上下波動,因此設置法向量τ(υ,ω,ξ)為(0,-0.26,1),將該法向量作為平面模型內點云的約束,約束角θ為15°,設置迭代次數n為500;由于底面矩形高h為0.355 m,部分高度略有偏差,因此λ應略微大于h;且由于有法向量作約束,因此λ可取較大值,本文設置為0.36 m。將切片內的點云向橫斷面方程進行投影從而獲取橫斷面輪廓線,由圖10中可觀察到隧道切片得到了完整分割。

圖10 橫斷面分割

3.2.2 橫斷面擬合及濾波降噪

混凝土厚度檢測、隧道表面起伏差檢測、隧道表面的變形監測都需要橫斷面作為重要參數。因此對橫斷面進行提取及擬合是三維激光技術中的重要工作。

(1)將隧道切片進行旋轉,從原始坐標系中旋轉至擬合坐標系,如圖11所示,圖中旋轉后的切片點云平行于XOY平面,從而實現三維數據向二維數據的轉化,為橫斷面擬合做準備。

圖11 橫斷面切片旋轉

(2)直接使用最小二乘結合橢圓模型擬合,如圖12所示。從圖中可以觀察到擬合橫斷面形狀與隧道原始形狀差距較大。

圖12 直接擬合

(3)使用剔除掉底部點云后的橫斷面進行擬合,如圖13所示。從圖中可觀察到擬合橫斷面與原始橫斷面形狀接近,更能反映隧道原始橫斷面的真實情況。

圖13 橫斷面切割底面后擬合

(4)橫斷面濾波。對橫斷面進行最小二乘擬合后即可根據橫斷面參數對橫斷面輪廓點云進行降噪,本文取局部隧道共334 000個點進行濾波分析。如圖14所示為降噪前后的誤差分布圖,由圖可觀察到濾波前誤差分布在-0.3～0.15 m;濾波后誤差分布在-0.1～0.1 m,濾波前后誤差分布改變明顯,非點得到明顯剔除。

(a)濾波前

(b)濾波后圖14 濾波前后誤差分布圖

對濾波前后的噪點進行統計分析,從表2中可觀察到剔除噪點37 598個;平均值從-0.009 m降低至0 m;均方差從0.041 m降低至0.019 m,濾波前后誤差明顯降低。綜上所述本文降噪效果較好,降噪方法較優。

表2 去噪前后誤差分析

4 結束語

本文針對隧道底部為矩形的隧道橫斷面分割、擬合及濾波提出了一種利用隧道腰部邊界特征提取隧道腰部中線,接著依據隧道腰部中線提取隧道橫斷面,最后剔除出隧道底部點云并擬合橫斷面。得出以下結論：

(1)針對半圓形隧道上下不對稱、中線難以提取的特點,本文依據邊界特征提取隧道腰部邊界,進而提取腰部中線,并使用設計中樁進行檢驗;實驗表明本文提取中線與設計中線樁點ΔX、ΔY均在11 mm以內,ΔZ在15 mm以內,提取精度較高,效果良好。

(2)針對半圓形隧道橫斷面提取問題,依據底面坡度數據計算法向量,在平面檢核中作為約束條件進行參數求取;接著根據隧道下部矩形設計高度進行分割。相比于直接使用RANSAC進行平面檢核,該方法針對性更強,分割效果更好。

(3)使用剔除掉底部點云的隧道橫斷面進行擬合,相比于直接使用最小二乘法橫斷面擬合該方法形狀較為符合橫斷面形狀。根據最小二乘原理對橫斷面點云進行濾波并以局部隧道進行實驗驗證,濾波前誤差分布在-0.3～0.15 m,濾波后誤差分布在-0.1～0.1 m;此外濾波前后均方根誤差分別為0.041 m和0.019 m,綜上所述本文方法擬合精度較高。