歐海軍, 馮騰飛, 沈月千, 張 堅(jiān)
(1.廣州市城市規(guī)劃勘測設(shè)計(jì)研究院, 廣州 510053; 2.同濟(jì)大學(xué)測繪與地理信息學(xué)院, 上海 200092;3.河海大學(xué)地球科學(xué)與工程學(xué)院, 南京 211100)
機(jī)載激光雷達(dá)系統(tǒng)(light detection and ranging, LiDAR)是一種以飛機(jī)為載體的集激光測距、全球定位系統(tǒng)(GPS)和慣性導(dǎo)航系統(tǒng)(INS)三種技術(shù)優(yōu)點(diǎn)于一體的系統(tǒng),通過三種技術(shù)的結(jié)合,能夠精準(zhǔn)地測定激光束落在激光腳點(diǎn)的三維空間坐標(biāo),從而獲得地物的三維坐標(biāo)點(diǎn)云數(shù)據(jù)[1]。但與其他數(shù)據(jù)類型不同,由于LiDAR點(diǎn)云數(shù)據(jù)分布反映著地物表面的空間變化,具有數(shù)據(jù)量大且分布較為離散等特點(diǎn),因此對點(diǎn)云進(jìn)行后續(xù)處理前,必須率先進(jìn)行三維點(diǎn)云數(shù)據(jù)濾波,以便實(shí)現(xiàn)更快速、高效的點(diǎn)云分類及地物提取[2]。
濾波算法性能的提升及其在不同地形的適用性一直備受中外學(xué)者的關(guān)注。Sithole等[3]提出了基于地形坡度的濾波算法,該算法通過對比目標(biāo)點(diǎn)與其臨近點(diǎn)間的高差值與給定閾值之間的大小關(guān)系,判斷目標(biāo)點(diǎn)是否為地面點(diǎn),但該方法地形坡度閾值設(shè)置為30%,并不適用于起伏較大的區(qū)域;Kraus等[4]提出了迭代最小二乘線性內(nèi)插的濾波算法,利用原始點(diǎn)云與最小二乘內(nèi)插的擬合高程求差,并根據(jù)差值是否服從正態(tài)分布來判斷待定點(diǎn)是否為非地面點(diǎn),但該算法只適用于地形起伏不大的機(jī)載點(diǎn)云區(qū)域,對于地形復(fù)雜的區(qū)域,擬合平面會隨著地面的起伏變化而變化,進(jìn)而造成較差的濾波效果,且數(shù)據(jù)插值固有的誤差會增加濾波過程的不確定性;Lindenmayer等[5]提出了基于圖像處理方式的數(shù)學(xué)形態(tài)學(xué)濾波算法,該方法在濾波區(qū)域選擇一個固定大小的窗口,并利用窗口中其他點(diǎn)與預(yù)先運(yùn)用數(shù)學(xué)形態(tài)運(yùn)算找到的最低點(diǎn)求差,通過與設(shè)定閾值的比較進(jìn)行判斷,若差值小于閾值則為地面點(diǎn),該方法也采用了固定閾值的方式,增加了濾波困難;張小紅[6]和劉經(jīng)南等[7]均基于激光腳點(diǎn)在一定程度上能夠反映地形表面的空間起伏變化這一原理提出了移動曲面擬合濾波算法,該方法利用曲面擬合的高程與實(shí)際點(diǎn)云數(shù)據(jù)的高差與固定的閾值比較,從而判別地面點(diǎn),但該方法不僅需要尋找種子區(qū)域,還采用了固定閾值,對濾波效果造成較大影響。此外,李鵬程等[8]基于高于地面3 m的點(diǎn)為非地面點(diǎn),采用區(qū)塊索引機(jī)制將整個點(diǎn)云區(qū)域分成很多區(qū)塊,提出基于距離圖像的分塊曲面濾波方法,該方法很大程度上依賴于區(qū)塊大小的確定,而區(qū)塊的確定需要人工干預(yù),大大增加了工作量。
對比現(xiàn)有的LiDAR點(diǎn)云數(shù)據(jù)濾波方法可知大部分方法存在一定的共同點(diǎn)及局限性。首先,一些算法需要在濾波之前對數(shù)據(jù)進(jìn)行插值處理,以此來減少存儲空間,提高運(yùn)算速度,但在數(shù)據(jù)內(nèi)插的過程中會破壞點(diǎn)云數(shù)據(jù)本身存在的內(nèi)在性質(zhì)以及點(diǎn)云之間的鄰近關(guān)系,最終使得濾波的效果并不佳[9-12];其次,部分算法需要選擇含地面點(diǎn)的種子區(qū)域,但選取方法缺乏科學(xué)的驗(yàn)證,可能會造成由于種子區(qū)域選取不當(dāng)帶來的濾波精度損失;此外,絕大部分的濾波都需要設(shè)定濾波閾值,但固定且單一的閾值使算法很大程度上不能自適應(yīng)地形的變化,造成較差的濾波效果?;谏鲜龇治?,提出一種改進(jìn)的自適應(yīng)閾值分塊濾波算法。一方面,對原始點(diǎn)云進(jìn)行不斷的分塊,降低計(jì)算時內(nèi)存的占用空間,加快算法的計(jì)算速度;另一方面,對分塊點(diǎn)云進(jìn)行既定坡度平面擬合,對超出自適應(yīng)濾波閾值的點(diǎn)進(jìn)行分類,并不斷進(jìn)行計(jì)算,最終達(dá)到地面點(diǎn)與非地面點(diǎn)分離的效果。
點(diǎn)云數(shù)據(jù)在獲取的過程中會產(chǎn)生噪聲,這些噪聲會嚴(yán)重干擾點(diǎn)云濾波的過程,甚至可能造成點(diǎn)云數(shù)據(jù)濾波的失敗。因此點(diǎn)云濾波之前,必須對點(diǎn)云數(shù)據(jù)進(jìn)行去噪處理。為了更好地得到無噪聲數(shù)據(jù),采用“兩步走”的方式進(jìn)行去噪,首先利用高斯濾波對點(diǎn)云進(jìn)行處理:根據(jù)點(diǎn)云中所有點(diǎn)的距離應(yīng)構(gòu)成高斯分布這一特性,計(jì)算每個點(diǎn)與其最鄰近k個點(diǎn)的平均距離,并取其平均值與方差作為參考值,進(jìn)而根據(jù)“3倍中誤差準(zhǔn)則”[13]剔除異常點(diǎn)。此外,由于試驗(yàn)數(shù)據(jù)為無序點(diǎn)云,高斯濾波去噪可能存在剔除空洞,因此利用K-D(K-dimensional)樹去噪的方式進(jìn)行二次去噪:對已經(jīng)進(jìn)行粗去噪的點(diǎn)云數(shù)據(jù)構(gòu)建K-D樹,并隨機(jī)取點(diǎn)求平均高差d,刪除所有大于2d的點(diǎn),實(shí)現(xiàn)點(diǎn)云去噪,為點(diǎn)云濾波提供更精確的基礎(chǔ)點(diǎn)云數(shù)據(jù)。
在傳統(tǒng)分塊曲面濾波算法中,通常將原始點(diǎn)云分塊后,對每一塊點(diǎn)云進(jìn)行濾波,從而分離地面點(diǎn)與非地面點(diǎn),其基本步驟及尚存問題分析[14-15]如下。
(1)建立區(qū)塊索引機(jī)制。分塊曲面濾波算法利用建立區(qū)塊索引的方式對離散的點(diǎn)云數(shù)據(jù)進(jìn)行組織,使離散點(diǎn)云可以有效的進(jìn)行索引定位和處理。
(2)利用距離圖像確定分塊數(shù)。根據(jù)一定的重采樣間距將離散點(diǎn)云數(shù)據(jù)組織為規(guī)則格網(wǎng)數(shù)據(jù),每個格網(wǎng)中點(diǎn)云的高程均值為該網(wǎng)格的高程信息數(shù)據(jù),并對所得數(shù)據(jù)進(jìn)行歸化處理,生成距離圖像,結(jié)合圖像中最大建筑物的跨度來確定區(qū)塊的高度與寬度,從而確定分塊數(shù),但是相對于海量的點(diǎn)云數(shù)據(jù),該操作效率非常低。
(3)點(diǎn)云濾波。對上面所建立的所有區(qū)塊進(jìn)行曲面擬合,選取塊內(nèi)一定數(shù)量的高程點(diǎn)作為曲面擬合初點(diǎn)進(jìn)行曲面擬合,并計(jì)算塊內(nèi)其他點(diǎn)云高程與擬合高程的差值,使其與固定閾值作比較,如果小于或等于給定的閾值,則判定該點(diǎn)為地面點(diǎn),否則為非地面點(diǎn)。這樣固定閾值的方式不僅不適用于所有的區(qū)塊,而且容易造成誤分現(xiàn)象。
傳統(tǒng)的分塊曲面濾波算法利用區(qū)塊索引的方式進(jìn)行組織,其分塊的方式取決于離散點(diǎn)云所生成的灰度圖像,與大部分濾波算法類似,一般根據(jù)經(jīng)驗(yàn)來選擇閾值,比較常用的是3 m或者5 m,但在實(shí)際中,固定的閾值往往使濾波結(jié)果出現(xiàn)不同程度的點(diǎn)云錯分現(xiàn)象(濾波算法誤差),從而造成點(diǎn)云濾波的兩類錯誤增加,不利于點(diǎn)云的后續(xù)處理。同時人工選取圖像最大建筑物跨度也加大了工作量。因此,通過考慮分塊面積以及塊域內(nèi)最大高程差兩個閾值影響因素建立自適應(yīng)閾值模型,提出了一種改進(jìn)的自適應(yīng)閾值分塊曲面濾波算法,其詳細(xì)過程如下。
Step1建立虛擬網(wǎng)格。由于曲面擬合至少需要6個點(diǎn)位數(shù)據(jù),而每次選取的點(diǎn)位均為塊內(nèi)最低高程點(diǎn)位,因此將離散點(diǎn)云區(qū)域分割成6塊,且以包含點(diǎn)云的最小矩形窗口為初始窗口。
Step2自動選取種子點(diǎn)云。在分塊完成之后,獲取每塊3個最低的點(diǎn)位并取平均值,降低每塊去噪過程中造成的影響,利用所獲得的點(diǎn)位作為擬合曲面的初始種子點(diǎn),避免人工選取種子區(qū)域。
Step3地面點(diǎn)的獲取。利用所選種子點(diǎn)位進(jìn)行曲面擬合,分別計(jì)算每個點(diǎn)位高程與擬合曲面高程的差值,差值小于自適應(yīng)閾值的點(diǎn)即為地面點(diǎn),否則為非地面點(diǎn)。其中自適應(yīng)閾值的確定如下。
在改進(jìn)算法實(shí)現(xiàn)過程中,隨著分塊次數(shù)的增加,分塊面積也隨之不斷的縮小,因此濾波閾值應(yīng)具有以下特點(diǎn):
(1)面積越大,閾值越大,當(dāng)面積增大的某種程度時,閾值不再隨面積增大而變化。
(2)閾值隨著面積的遞減而漸趨穩(wěn)定,當(dāng)面積減小到一定的值時,閾值不再變化。根據(jù)上述特點(diǎn),建立閾值與塊域面積的函數(shù)關(guān)系,利用控制變量法,對多塊點(diǎn)云進(jìn)行反復(fù)試驗(yàn),直到達(dá)到濾波效果的最優(yōu),并利用最小二乘方法,擬合建立所記錄閾值與對應(yīng)塊域面積的函數(shù)關(guān)系。
(3)盡管上述函數(shù)在整體上表現(xiàn)最優(yōu),但當(dāng)塊域起伏較大時,濾波的效果并不佳,因此除了面積這一因素之外,還需考慮與對應(yīng)面積中的地形起伏,而最大高程差能夠體現(xiàn)出地形起伏的特點(diǎn),且每塊點(diǎn)云數(shù)據(jù)中最大高程差也很大程度上反映了分塊點(diǎn)云的空間狀態(tài)。因此基于上述分析建立了顧及塊域面積以及塊域內(nèi)最大高差兩個因素的自適應(yīng)閾值模型,經(jīng)過點(diǎn)云數(shù)據(jù)初步修正和實(shí)驗(yàn)驗(yàn)證,得到自適應(yīng)閾值與兩因素的關(guān)系式為
y=0.076 1s2-0.355 1s+0.228 2Δhmax+
2.727 9
(1)
式(1)中:y為兩種因素影響下的閾值;s為分塊面積;Δhmax為區(qū)塊中高差最大值。
Step4重復(fù)Step 1—Step 3,迭代濾波直到?jīng)]有新的地面點(diǎn)產(chǎn)生。
由于點(diǎn)云塊域隨著算法的進(jìn)行,分塊區(qū)域內(nèi)點(diǎn)云數(shù)目不斷減小,且塊數(shù)也在不斷減小,當(dāng)塊數(shù)小于6時,不滿足曲面擬合條件,迭代終止。算法的主要流程如圖1所示。

圖1 改進(jìn)濾波算法流程圖
采用具有不同特性的三組機(jī)載激光雷達(dá)點(diǎn)云影像作為實(shí)驗(yàn)數(shù)據(jù),如圖2所示,試驗(yàn)區(qū)1包括植被、道路、建筑物(主要以廠房為主)以及建筑物附屬設(shè)施,整個區(qū)域包含195 990個點(diǎn)云數(shù)據(jù),區(qū)域內(nèi)最大高差約為20.27 m;試驗(yàn)區(qū)2地勢較為平緩,包含復(fù)雜低矮的建筑物,建筑物較多,且植被較少,區(qū)域內(nèi)最大高差約為30 m,包含162 426個離散點(diǎn)云;試驗(yàn)區(qū)3包含較多植被、新蓋的樓盤以及較高古塔,區(qū)域內(nèi)最大高差為40 m,包含174 417個離散點(diǎn)云。分別利用傳統(tǒng)分塊曲面濾波及自適應(yīng)閾值分塊曲面濾波對三個試驗(yàn)區(qū)點(diǎn)云影像進(jìn)行處理,通過國際攝影測量與遙感協(xié)會(ISPRS)的評價體系驗(yàn)證本文改進(jìn)算法的可行性及優(yōu)越性。該評價體系包括兩類誤差及總體誤差指標(biāo),其中第Ⅰ類誤差定義為將地面點(diǎn)誤分為非地面點(diǎn)的比例;第Ⅱ類誤差則是將非地面點(diǎn)誤分為地面點(diǎn)的比例;總體錯誤率即所有的類別錯分點(diǎn)占總點(diǎn)數(shù)比例,錯誤率越低代表整體濾波效果越好。具體如表1所示,第Ⅰ類誤差:b/e×100%;第Ⅱ類誤差:c/f×100%;總誤差:(b+c)/n×100%。

圖2 原始點(diǎn)云影像

表1 誤分率的定義
對上述實(shí)驗(yàn)數(shù)據(jù)進(jìn)行去噪,再分別利用傳統(tǒng)算法以及改進(jìn)算法進(jìn)行濾波,得到結(jié)果如圖3所示,其中圖3(a)和3(b)、圖3(e)和3(f)為傳統(tǒng)算法和本文改進(jìn)算法對原始點(diǎn)云進(jìn)行濾波得到的非地面點(diǎn)云對比圖,圖3(c)和圖3(d)為地面點(diǎn)云對比圖。就傳統(tǒng)方法濾波結(jié)果而言,分別存在地面點(diǎn)錯分為非地面點(diǎn)[圖3(c)紅圈所示]以及非地面點(diǎn)錯分為地面點(diǎn)[圖3(a)和圖3(e)紅圈所示]的現(xiàn)象。而對于改進(jìn)算法的濾波結(jié)果,可以明顯看出,圖3(b)和圖3(d)中的第II類錯誤已基本消除,圖3(f)中第I類錯誤也明顯減小,三塊點(diǎn)云的錯分現(xiàn)象均得到了極大的改善,表明改進(jìn)算法既能很好地適應(yīng)不同地形,也解決了閾值選擇問題對濾波帶來的精度損失,濾波效果明顯優(yōu)于傳統(tǒng)算法。

圖3 傳統(tǒng)算法與改進(jìn)算法濾波效果對比
進(jìn)一步利用誤差公式對傳統(tǒng)分塊濾波結(jié)果以及改進(jìn)濾波結(jié)果進(jìn)行定量評價,結(jié)果如表2和表3所示。由表2可知,I類誤差多于II類誤差,說明傳統(tǒng)算法更加側(cè)重于降低II類誤差。其中試驗(yàn)區(qū)1和試驗(yàn)區(qū)2的I類誤差都較II類誤差大,說明傳統(tǒng)算法在低矮建筑物較多的區(qū)域存在過濾現(xiàn)象,容易將地面點(diǎn)錯分為低矮建筑物。從表3可知,I類誤差明顯低于II類誤差,說明本文算法能很好地適應(yīng)不同的地形,具有較強(qiáng)的魯棒性。此外,對比表2和表3中三類誤差,可以看出改進(jìn)算法的濾波精度明顯提高。試驗(yàn)區(qū)2的I類誤差幅度最大,說明改進(jìn)算法更適用于建筑物較多、植被較少的地勢。與傳統(tǒng)分塊濾波算法相比,I類誤差值降低幅值較大,說明改進(jìn)算法有更強(qiáng)的適應(yīng)性。通過對II類誤差的對比分析,改進(jìn)算法相比于傳統(tǒng)分塊濾波算法II類誤差最大降低了3.34%,最小降低0.15%??傉`差最大降低了4.73%,最小降低了0.41%,說明了改進(jìn)算法比傳統(tǒng)分塊曲面濾波有更高的精度。

表2 傳統(tǒng)濾波算法誤差評價

表3 改進(jìn)濾波算法誤差評價
綜上所述,無論是在低矮建筑區(qū),還是在復(fù)雜多植被的城區(qū),改進(jìn)算法都能夠顯著減低I類誤差和II類誤差,驗(yàn)證了改進(jìn)算法的適用性以及有效性。
利用網(wǎng)格法對點(diǎn)云進(jìn)行逐級分塊,并以曲面擬合的方式自動獲取塊域種子點(diǎn),同時顧及分塊面積及塊域內(nèi)最大高程差兩個閾值影響因素建立了自適應(yīng)閾值模型,從而提出一種改進(jìn)的自適應(yīng)閾值分塊曲面濾波算法,很好地滿足了機(jī)載LiDAR點(diǎn)云濾波的要求。實(shí)驗(yàn)結(jié)果表明改進(jìn)算法在濾波時能夠?qū)崿F(xiàn)更有效的閾值選擇,且濾波精度較傳統(tǒng)方法有明顯提高。但對于閾值的影響因素還需進(jìn)一步研究,以期實(shí)現(xiàn)更高精度的自適應(yīng)閾值模型。