一種顧及地形的點(diǎn)云自適應(yīng)坡度濾波方法

丁少鵬，劉如飛，蔡永寧2，王鵬

(1.山東科技大學(xué) 測(cè)繪科學(xué)與工程學(xué)院，山東 青島 266590；2.濟(jì)南市勘察測(cè)繪研究院，濟(jì)南 250101)

0 引言

三維激光掃描系統(tǒng)能快速獲取高精度的地理信息數(shù)據(jù)，近些年得到了廣泛的應(yīng)用。其中，車載移動(dòng)測(cè)量系統(tǒng)能夠精確地獲取公路兩側(cè)一些近地物信息，在智慧交通、道路管理和城市規(guī)劃等方面具有重要作用[1]。近些年許多專家學(xué)者圍繞著車載移動(dòng)數(shù)據(jù)中道路提取及目標(biāo)分類做了大量研究，但是如何能夠在地面濾波過(guò)程中顧及各種復(fù)雜地形特征，是一個(gè)亟待解決的問(wèn)題。

在點(diǎn)云地面濾波方面，許多學(xué)者做了相關(guān)的研究，其方法可以分為四類：首先，基于形態(tài)學(xué)的方法原理比較簡(jiǎn)單，如何選擇合適的濾波窗口是最關(guān)鍵的一步[1-3]，董保根等[4]在插值方面做了相應(yīng)的改進(jìn)，但是仍然無(wú)法在起伏的地形中取得滿意的效果。其次，基于不規(guī)則三角網(wǎng)的方法可以取得很好的濾波效果[5-6]，但是三角網(wǎng)的創(chuàng)建過(guò)程較復(fù)雜，占用內(nèi)存較大，濾波時(shí)間較長(zhǎng)，如何改進(jìn)算法提高濾波效率，減少內(nèi)存占用，是現(xiàn)在需要研究的主要問(wèn)題。而對(duì)于基于曲面擬合的方法[7-9]，需要選取足夠準(zhǔn)確的點(diǎn)進(jìn)行擬合運(yùn)算，參與擬合運(yùn)算的點(diǎn)越多，擬合出的平面就越準(zhǔn)確。所以很多算法都是先插值再擬合，但是插值運(yùn)算一般會(huì)造成一些精度損失。最后是基于坡度的濾波方法[10]，張皓等[11]提出通過(guò)設(shè)置坡度增量來(lái)判斷非平面地形，Susaki[12]提出利用數(shù)字地面模型(digital terrain model，DTM)計(jì)算坡度閾值，何培培等[13]又增加了不斷變化窗口，利用DTM自動(dòng)更新坡度閾值通過(guò)多次迭代進(jìn)行濾波。以上方法需要手動(dòng)輸入多個(gè)閾值參數(shù)，而不斷迭代計(jì)算坡度也增加了運(yùn)算成本。

總體來(lái)說(shuō)，相比前3類方法，基于坡度的方法原理簡(jiǎn)單，不需要插值運(yùn)算并能很好地保留原始地形特征信息，但如何選取合適的閾值進(jìn)行濾波是最主要問(wèn)題。本文提出一種顧及地形的點(diǎn)云自適應(yīng)坡度濾波方法，可動(dòng)態(tài)選取濾波窗口和坡度閾值，實(shí)現(xiàn)對(duì)多種地形的濾波。

1 坡度濾波算法改進(jìn)

基于坡度的濾波方法最早由Vosselman提出[10]，傳統(tǒng)的坡度法只能設(shè)置固定坡度閾值，對(duì)單一地形進(jìn)行濾波時(shí)可以取得不錯(cuò)的效果，無(wú)法顧及所有的地形。坡度閾值過(guò)大，無(wú)法過(guò)濾掉一些距離地面較近的非地面點(diǎn)，坡度閾值過(guò)小又會(huì)造成過(guò)度濾波，丟失一些地面點(diǎn)。除此之外，如果在處理不同的地形時(shí)使用單一尺度的窗口也不具有自適應(yīng)性。因此，傳統(tǒng)的坡度濾波方法主要存在無(wú)法自適應(yīng)選取窗口大小與坡度閾值這2個(gè)問(wèn)題。

算法主要分為數(shù)據(jù)組織與預(yù)處理、地形判斷與多尺度濾波、優(yōu)化處理與數(shù)據(jù)輸出3部分，對(duì)傳統(tǒng)坡度濾波方法存在的問(wèn)題進(jìn)行改進(jìn)，主要體現(xiàn)在地形判斷與多尺度濾波部分。具體方法如圖1所示。

圖1 算法流程圖

1.1 數(shù)據(jù)組織與預(yù)處理

1)創(chuàng)建格網(wǎng)索引。為了方便數(shù)據(jù)的組織，提高運(yùn)行效率，首先根據(jù)文獻(xiàn)[14]建立二維格網(wǎng)，將原始點(diǎn)云投影到xoy平面。

2)去噪與選點(diǎn)。傳統(tǒng)坡度濾波方法直接選取每個(gè)格網(wǎng)內(nèi)的最低點(diǎn)作為地面種子點(diǎn)。但是在采集點(diǎn)云的過(guò)程中，會(huì)因?yàn)橐恍┢≡诳諝庵须s物的遮擋，或掃描到低于地面的一些坑洞而產(chǎn)生噪點(diǎn)。這些噪點(diǎn)一般分布在整塊點(diǎn)云周邊，點(diǎn)數(shù)相對(duì)較少。若選取的最低點(diǎn)為噪點(diǎn)，則會(huì)對(duì)后期濾波造成干擾。為了準(zhǔn)確選取地面種子點(diǎn)，需要增加約束條件。

如圖2所示，本文在選取地面種子點(diǎn)前，先進(jìn)行高程分析統(tǒng)計(jì)去噪，即對(duì)每個(gè)格網(wǎng)中的點(diǎn)按高度分層，若某一層中點(diǎn)數(shù)小于平均點(diǎn)密度的1/10，則認(rèn)為該點(diǎn)為孤立的噪點(diǎn)，需要去除。同時(shí)為了避免所選取的最低點(diǎn)是位于地面以下的噪點(diǎn)，改為選取每個(gè)格網(wǎng)內(nèi)的次低點(diǎn)作為地面種子點(diǎn)[7]。

如圖3所示，若中心格網(wǎng)內(nèi)所包含的點(diǎn)全是高于地面的非地面點(diǎn)，所選取的地面種子點(diǎn)也會(huì)是非地面點(diǎn)。因此，可以利用高差進(jìn)行判斷，將高差閾值H設(shè)定為2倍格網(wǎng)邊長(zhǎng)，計(jì)算該中心點(diǎn)與8鄰域格網(wǎng)內(nèi)地面種子點(diǎn)的高差ΔHi(i=1，2…，8)，若ΔH大于H的鄰域格網(wǎng)點(diǎn)數(shù)N=8，則認(rèn)為該點(diǎn)是非地面點(diǎn)，需要去除。

圖2 格網(wǎng)分層孤立>點(diǎn)示意圖

圖3 中心點(diǎn)與8鄰域過(guò)高示意圖

1.2 地形判斷與多尺度濾波

如圖4所示，計(jì)算某窗口內(nèi)臨近點(diǎn)之間的坡度，并對(duì)坡度進(jìn)行分段，統(tǒng)計(jì)出坡度的頻率分布情況，由曲線可以看出窗口內(nèi)點(diǎn)云坡度值符合正態(tài)分布的統(tǒng)計(jì)規(guī)律。

圖4 坡度頻率分布圖

進(jìn)行濾波時(shí)，首先要計(jì)算較大窗口內(nèi)的坡度值，統(tǒng)計(jì)局部坡度分布情況，將地形分為3類；然后根據(jù)地形選用不同尺度的濾波窗口，計(jì)算出該窗口范圍內(nèi)的局部坡度閾值；最后利用計(jì)算出的局部坡度閾值，對(duì)相應(yīng)區(qū)域內(nèi)的點(diǎn)云進(jìn)行判斷，得到準(zhǔn)確的地面點(diǎn)。

1)地形判斷。進(jìn)行地形判斷時(shí)，所選擇的窗口大小應(yīng)該顧及周圍地形，不能過(guò)小，但是如果窗口過(guò)大，會(huì)同時(shí)包含多種地形信息，錯(cuò)判地形情況，導(dǎo)致濾波失敗?？紤]到以上問(wèn)題，根據(jù)所處理的數(shù)據(jù)，本文設(shè)定9×9格網(wǎng)作為地形判斷窗口Wb。

根據(jù)正態(tài)分布的統(tǒng)計(jì)規(guī)律，剔除一些可能有異常的坡度值，即計(jì)算出坡度均值μ與標(biāo)準(zhǔn)差σ，選取(μ-2σ，μ+2σ)內(nèi)的最大坡度作為該窗口內(nèi)判斷地形的坡度值S0?？紤]到實(shí)際道路路面會(huì)有起伏，所以選擇坡度值0.3作為路面與斜坡的分界值。若S0<0.3，則認(rèn)為是路面等平地；0.3≤S0<1則認(rèn)為是緩坡；S0≥1的部分可認(rèn)為是地形比較復(fù)雜的陡坡。

2)多尺度濾波。根據(jù)坡度值S0判斷出地形后，選擇合適的濾波窗口W0進(jìn)行濾波。路面等平坦地形坡度變化較小，采用較大窗口可以更準(zhǔn)確地判斷中心點(diǎn)坡度，所以直接采用地形判斷窗口Wb計(jì)算窗口內(nèi)中心點(diǎn)與鄰域格網(wǎng)內(nèi)點(diǎn)的坡度。山地等地形坡度變化較大，若窗口過(guò)大，計(jì)算出的坡度閾值可能無(wú)法準(zhǔn)確表示中心點(diǎn)坡度，所以要選擇較小的窗口?？紤]到數(shù)據(jù)要具有統(tǒng)計(jì)特性，數(shù)量不能過(guò)少，所以選擇5×5格網(wǎng)作為較小大的濾波窗口Ws。緩坡地形則選擇介于二者之間的7×7格網(wǎng)作為濾波窗口Wm。如式(1)選擇合適的濾波窗口。

(1)

選取合適的窗口后，需要計(jì)算出坡度閾值Sm來(lái)確定地面點(diǎn)。坡度閾值由中心點(diǎn)與窗口W0內(nèi)的所有鄰域格網(wǎng)點(diǎn)的坡度值決定。同樣選取(μ-2σ，μ+2σ)范圍內(nèi)所有坡度Si，剔除粗差或疑似異常值，按照式(2)計(jì)算坡度閾值Sm。然后以中心點(diǎn)與8鄰域格網(wǎng)點(diǎn)的最大坡度作為該點(diǎn)坡度，若小于閾值Sm則認(rèn)為是地面點(diǎn)，否則需去除。

(2)

式中：N為在(μ-2σ，μ+2σ)范圍內(nèi)的坡度數(shù)量。

1.3 優(yōu)化處理與數(shù)據(jù)輸出

1)去噪優(yōu)化。完成上述濾波后，可能還會(huì)存在少量噪點(diǎn)，如圖5紅點(diǎn)所示。這類噪點(diǎn)主要特征是周圍沒(méi)有臨近點(diǎn)，處于孤立的位置。本文基于鄰域判斷去噪，以待判點(diǎn)為中心點(diǎn)，以窗口W0為判斷范圍，統(tǒng)計(jì)該點(diǎn)的鄰域格網(wǎng)內(nèi)是否有點(diǎn)，如果有點(diǎn)的格網(wǎng)數(shù)小于窗口內(nèi)格網(wǎng)總數(shù)的一半，即周圍點(diǎn)數(shù)過(guò)少，則去除該中心點(diǎn)。

圖5 孤立點(diǎn)示意圖

2)提取所有地面點(diǎn)。以所得到的地面種子點(diǎn)提取所有地面點(diǎn)。對(duì)于中心格網(wǎng)內(nèi)各點(diǎn)，分別計(jì)算其與8鄰域格網(wǎng)內(nèi)地面種子點(diǎn)的坡度，求取均值S8m，以上文中的閾值Sm進(jìn)行判斷，若S8m

2 實(shí)驗(yàn)分析

算法由C++/QT語(yǔ)言實(shí)現(xiàn)。實(shí)驗(yàn)數(shù)據(jù)共2組，都以高程進(jìn)行渲染。如圖6(a)所示，數(shù)據(jù)1是泰安某一段公路，長(zhǎng)度約100 m，共有2 048 061個(gè)點(diǎn)。其中包括平坦路面、較陡的邊坡以及中間的緩坡等地形，并含有樹(shù)木、護(hù)欄、車輛等非地面物。公路較平坦，沒(méi)有明顯起伏。陡坡處地形較復(fù)雜，中間部分有明顯凸起，斜坡整體坡度在0.8左右，局部最大坡度在2以上。

圖6(c)為原始點(diǎn)云的側(cè)視圖，從圖中可更加直觀地看到所有地物以及地形變化情況。圖6(b)為濾波效果圖，與原數(shù)據(jù)相比，緩坡處的兩排樹(shù)已經(jīng)被過(guò)濾掉，整個(gè)地面保留較完整，無(wú)漂浮的噪點(diǎn)。同時(shí)，從圖6(d)中的濾波效果側(cè)視圖可以進(jìn)一步看到，平坦的公路和高低起伏的邊坡等地形都被完整地保留下來(lái)，公路上的汽車和公路右側(cè)的低矮灌木護(hù)欄已經(jīng)被過(guò)濾掉。

圖6 公路點(diǎn)云濾波效果對(duì)比圖

數(shù)據(jù)2為某礦山邊坡數(shù)據(jù)，如圖7(a)所示，相較于數(shù)據(jù)1，整個(gè)測(cè)區(qū)無(wú)植被，但是地形變化較大。藍(lán)色部分是坡度為1的邊坡，中部為平地，剩余部分是坡度近似垂直的陡坡，陡坡中部有凸起，其中地物主要包含1號(hào)區(qū)域及周圍平地上的卡車和2號(hào)區(qū)域內(nèi)的挖掘機(jī)。圖7(b)為濾波后的效果圖，從圖中可以看到卡車及挖掘機(jī)等地物都已被過(guò)濾，平地及緩坡部分的地面點(diǎn)保留比較完整。

圖7(c)為區(qū)域1濾波效果圖，白色為濾波后的點(diǎn)，紅色為原始點(diǎn)云。疊加后可發(fā)現(xiàn)，左側(cè)的A邊坡保留相對(duì)比較完整，B邊坡中的部分凸起坡度接近垂直，因?yàn)檫x取地面種子點(diǎn)時(shí)只保留了低點(diǎn)，格網(wǎng)內(nèi)的其他坡度較大點(diǎn)需要通過(guò)計(jì)算出的坡度閾值進(jìn)行比對(duì)，保留小于坡度閾值的點(diǎn)。而計(jì)算坡度閾值時(shí)，剔除了2σ以外的極值點(diǎn)，然后求取平均值，垂直部分的坡度的較大點(diǎn)就有可能大于坡度閾值，造成一部分不能保留。圖7(d)為區(qū)域2局部濾波效果圖，A邊坡是與地面垂直的陡坡部分，同樣也沒(méi)有被完整地保留下來(lái)，而B(niǎo)邊坡部分基本上保留完整。

分析2塊點(diǎn)云數(shù)據(jù)的實(shí)驗(yàn)結(jié)果，可以看出算法可以過(guò)濾掉2塊點(diǎn)云的絕大部分的地物點(diǎn)。在數(shù)據(jù)1中，平地、緩坡、陡坡都是連續(xù)的，且無(wú)近似垂直的陡坡，所以整體地形保留較完整。數(shù)據(jù)2中，除無(wú)法保留與地面近似垂直的邊坡外，其他邊坡都可以較完整地保留下來(lái)。

圖7 礦山邊坡點(diǎn)云濾波效果對(duì)比圖

為了進(jìn)一步驗(yàn)證濾波效果，將本文方法與傳統(tǒng)的坡度法和擬合法進(jìn)行對(duì)比。傳統(tǒng)的坡度法設(shè)置固定的8鄰域?yàn)V波窗口和固定坡度閾值0.2和2，擬合法則是在取最低點(diǎn)后去除明顯非地面點(diǎn)，局部擬合出二次曲面，計(jì)算所有點(diǎn)到曲面的距離，分別設(shè)定閾值0.5 m和1 m進(jìn)行濾波。手動(dòng)將2組數(shù)據(jù)的地面點(diǎn)和地物點(diǎn)進(jìn)行分類，GP表示地面點(diǎn)數(shù)，OP表示非地面點(diǎn)數(shù)，OFP表示被過(guò)濾的地面點(diǎn)的點(diǎn)數(shù)，IFP表示錯(cuò)判為地面點(diǎn)的點(diǎn)數(shù)。按公式(3)計(jì)算出一類誤差I(lǐng)e、二類誤差I(lǐng)Ie和總體誤差A(yù)e，對(duì)3種方法濾波結(jié)果進(jìn)行定量評(píng)價(jià)。一類誤差和二類誤差反映了適用性，總體誤差反映了算法的可行性，統(tǒng)計(jì)結(jié)果如表1所示。

(3)

表1 濾波誤差分析表

從表1中可以看出，對(duì)于公路數(shù)據(jù)，本文提出的顧及地形的點(diǎn)云自適應(yīng)坡度濾波方法相比傳統(tǒng)的坡度濾波方法在一類誤差和總體誤差方面有很大地提升。利用傳統(tǒng)的坡度濾波方法，設(shè)置坡度閾值0.2，二類誤差僅為0.02%，但因?yàn)橛蟹瞧矫娴牡匦?造成了大范圍的過(guò)度濾波,一類誤差高達(dá)54.15%。擬合法采用了二次曲面擬合，效果較好，但相對(duì)而言本文算法的3種誤差都更低，可以取得更好的濾波效果。

礦山邊坡數(shù)據(jù)中，傳統(tǒng)坡度法依然會(huì)大范圍過(guò)度濾波，設(shè)置坡度為2的情況下，雖保留下了大部分地面點(diǎn)，但非地面點(diǎn)也無(wú)法去除，造成二類誤差較高。擬合法雖然保留了整體的地形特征，但是二類誤差卻整體偏大，無(wú)法過(guò)濾掉卡車等地物。從表中的數(shù)據(jù)結(jié)合濾波效果圖來(lái)看，相對(duì)于傳統(tǒng)坡度法和擬合法而言，本文方法雖然一類誤差、總體誤差稍微偏高，但是可以保留下整體地形特征，算法具有一定的可行性。

3 結(jié)束語(yǔ)

本文提出的顧及地形的點(diǎn)云自適應(yīng)坡度濾波方法，可以根據(jù)不同地形自適應(yīng)地確定坡度閾值，完成地面濾波。算法無(wú)需迭代，提高了點(diǎn)云濾波的自動(dòng)化水平。最后的實(shí)驗(yàn)結(jié)果表明，算法在去除地物點(diǎn)的同時(shí)可以保留地形特征信息。但是在垂直的陡坡地形中，雖然能夠過(guò)濾地物，邊坡信息卻無(wú)法完整保留下來(lái)。所以，如何能夠在處理這種地形時(shí)保留更多的地面點(diǎn)和進(jìn)一步提高自動(dòng)化水平是后續(xù)研究重點(diǎn)。