黃代軍

(浙江省測繪科學(xué)技術(shù)研究院 浙江杭州 310030)

機(jī)載LiDAR技術(shù)經(jīng)過20多年的發(fā)展，其硬件設(shè)備已經(jīng)取得了長足的進(jìn)步。作為主動式數(shù)據(jù)采集手段，機(jī)載LiDAR技術(shù)已經(jīng)廣泛應(yīng)用在場景數(shù)據(jù)采集中，如建筑物數(shù)據(jù)采集、電線巡檢數(shù)據(jù)采集和中小型區(qū)域地形植被數(shù)據(jù)采集等。相比于傳統(tǒng)的航空攝影測量，機(jī)載LiDAR不僅抗干擾能力強(qiáng)，而且采集的點(diǎn)云數(shù)據(jù)信息豐富，可以全天候作業(yè)。點(diǎn)云數(shù)據(jù)的優(yōu)點(diǎn)在于可以真實(shí)反映場景中地物的細(xì)節(jié)信息，如強(qiáng)度信息和位置信息等，但也存在一些缺點(diǎn)，如采樣粒度不均勻、冗余高和數(shù)據(jù)量大等。

點(diǎn)云數(shù)據(jù)的后處理問題一直是測繪行業(yè)研究的難題之一，相比于點(diǎn)云采集的硬件設(shè)備，點(diǎn)云數(shù)據(jù)的處理技術(shù)仍然相對滯后，其中點(diǎn)云數(shù)據(jù)的濾波更是點(diǎn)云數(shù)據(jù)處理中的關(guān)鍵點(diǎn)和難點(diǎn)之一[1]。經(jīng)過幾十年的發(fā)展，許多專家學(xué)者提出并驗(yàn)證了多種點(diǎn)云濾波算法。目前，常用的點(diǎn)云數(shù)據(jù)濾波算法主要包括基于不規(guī)則三角網(wǎng)濾波算法、移動窗口濾波算法、數(shù)學(xué)形態(tài)學(xué)濾波算法和移動曲面擬合濾波算法等。本文在傳統(tǒng)的移動曲面擬合濾波算法的基礎(chǔ)上，引入偏度概念，對原始算法進(jìn)行了改進(jìn)，并采用國際攝影測量與遙感學(xué)會提供的數(shù)據(jù)，通過試驗(yàn)驗(yàn)證改進(jìn)算法的可行性，同時與原始算法進(jìn)行對比，判斷兩種算法的效果。

1 常用的機(jī)載LiDAR點(diǎn)云數(shù)據(jù)濾波方法

1.1 基于不規(guī)則三角網(wǎng)濾波算法

Terrasolid軟件使用基于不規(guī)則三角網(wǎng)濾波算法?；诓灰?guī)則三角網(wǎng)濾波算法的機(jī)理：首先，對所有數(shù)據(jù)進(jìn)行計(jì)算得到初始參數(shù)，再通過選取地面種子點(diǎn)構(gòu)成不規(guī)則三角網(wǎng)；其次，設(shè)置條件，通過向不規(guī)則三角網(wǎng)中添加滿足條件的地面點(diǎn)來加密不規(guī)則三角網(wǎng)；最后，當(dāng)?shù)孛鏇]有滿足加入三角網(wǎng)的點(diǎn)集時，結(jié)束運(yùn)算。

基于不規(guī)則三角網(wǎng)濾波算法僅能計(jì)算一些連續(xù)表面環(huán)境的點(diǎn)云，具有一定的局限性，對于城市中的一些不是連續(xù)表面中的點(diǎn)會因?yàn)椴荒軡M足算法設(shè)置的條件而被遺漏。為了解決這一問題，李卉等[2]使用鏡像的方法將遺漏的點(diǎn)加入到不規(guī)則三角網(wǎng)中。此方法雖有一定的可行性，但是理論依據(jù)不足。為此，他們在算法中融入了一種融合區(qū)域增長思想，大大減少了漏選點(diǎn)。

1.2 移動窗口濾波算法

移動窗口濾波算法的機(jī)理：首先，在某個范圍內(nèi)取一個比較大的窗口，找到最低點(diǎn)；其次，設(shè)置一定的閾值，如果最低點(diǎn)與窗口中最高點(diǎn)的高差不大于此閾值，那么可將此點(diǎn)視為地面點(diǎn)。通過不斷減少窗口大小和閾值，尋找到更為精確的地面點(diǎn)。同時，對利用不同窗口大小得到的地面點(diǎn)進(jìn)行配權(quán)，以加權(quán)平均的方式得到DEM成果。

移動窗口濾液算法在使用時，如果設(shè)置的窗口太大，可能導(dǎo)致表面不連續(xù)的點(diǎn)被濾掉，窗口太小有可能會將建筑物等識別為地面點(diǎn)。設(shè)置的閾值過小，可能會將表面不連續(xù)的點(diǎn)識別為非地面點(diǎn)；設(shè)置的閾值過大，可能會將植被等點(diǎn)識別為地面點(diǎn)。所以使用移動窗口濾波算法的關(guān)鍵在于參數(shù)的選取是否合適。

1.3 數(shù)學(xué)形態(tài)學(xué)濾波算法

數(shù)學(xué)形態(tài)學(xué)濾波算法在數(shù)字圖像處理中的應(yīng)用相當(dāng)廣泛，其機(jī)理是通過具有一定形態(tài)的結(jié)構(gòu)實(shí)現(xiàn)圖像中對應(yīng)形狀的提取，達(dá)到圖像識別和圖像分析的目的。

數(shù)字形態(tài)學(xué)濾波算法包括兩個不同的運(yùn)算流程：開運(yùn)算與閉運(yùn)算，通過膨脹與腐蝕兩個基本運(yùn)算實(shí)現(xiàn)。先進(jìn)行腐蝕操作，再進(jìn)行膨脹操作的是開運(yùn)算；先進(jìn)行膨脹操作，再進(jìn)行腐蝕操作的是閉運(yùn)算。開運(yùn)算可以實(shí)現(xiàn)對機(jī)載LiDAR點(diǎn)云的數(shù)據(jù)濾波。

數(shù)學(xué)形態(tài)學(xué)濾波算法受窗口尺寸選取影響較大，在實(shí)際使用前，需要對點(diǎn)云數(shù)據(jù)所在環(huán)境具有一定的先驗(yàn)知識，才能設(shè)置合適的窗口尺寸，從而實(shí)現(xiàn)對機(jī)載LiDAR點(diǎn)云數(shù)據(jù)的濾波[3]。

1.4 移動曲面擬合濾波算法

移動曲面擬合濾波算法可以實(shí)現(xiàn)對數(shù)據(jù)的濾波，計(jì)算速度快且簡單明了。

移動曲面擬合濾波算法的機(jī)理：首先，在種子區(qū)域內(nèi)找到最低的3個角點(diǎn)作為初始地面點(diǎn)，對其進(jìn)行二維排序，通過求解平面方程的方式擬合初始平面，在平面方程中加入相鄰角點(diǎn)的平面坐標(biāo)，得到擬合高程值；其次，比較擬合高程與觀測高程的大小，如果二者之間的差值在給定的閾值范圍內(nèi)，則此相鄰角點(diǎn)為地面點(diǎn)，否則為非地面點(diǎn)；再次，將判斷為地面點(diǎn)的相鄰角點(diǎn)與組成初始平面的3個角點(diǎn)擬合為1個新的曲面；最后，按照上述方法繼續(xù)判斷相鄰角點(diǎn)是否為地面點(diǎn)，當(dāng)判斷得到的地面點(diǎn)數(shù)量達(dá)到6個時，在增加一個新地面點(diǎn)的同時去掉一個最遠(yuǎn)的角點(diǎn)。

重復(fù)上述操作，通過簡單的地形曲面擬合方式完成對整個測區(qū)的過濾[4]。

2 基于偏度平衡的移動曲面擬合濾波算法

本文在移動曲面擬合濾波算法的基礎(chǔ)上，引入統(tǒng)計(jì)學(xué)中的偏度概念。偏度(skewness)，也稱偏度系數(shù)，是對變量的概率分布密度函數(shù)曲線偏離對稱度的統(tǒng)計(jì)，用以表現(xiàn)概率密度曲線尾部的相對長度。

假設(shè)有一隨機(jī)變量X，同時擁有三階和四階矩陣，那么偏度Sk可表示為[5]

(1)

標(biāo)準(zhǔn)正態(tài)分布的偏度為0，當(dāng)Sk<0時稱為負(fù)偏態(tài)，Sk>0時稱為正偏態(tài)；Sk的絕對值越大，表明其偏離程度越大。本文引入偏度概念的目的在于剔除干擾點(diǎn)。

試驗(yàn)前，假設(shè)數(shù)據(jù)中的地面點(diǎn)概率密度分布函數(shù)服從標(biāo)準(zhǔn)正態(tài)分布，非地面點(diǎn)呈偏態(tài)分布，對數(shù)據(jù)進(jìn)行干擾，即通過不斷地從點(diǎn)云中剔除非地面點(diǎn)的方式使得偏態(tài)分布更接近于正態(tài)分布。同時，為了使樣本的分布逐漸趨于平衡，還需要不斷地對偏度值進(jìn)行校正。

移動曲面擬合濾波算法實(shí)現(xiàn)原理較為簡單，但對濾波數(shù)據(jù)的要求較高，只有在點(diǎn)云密度足夠高的前提下才能得到較好的效果。該算法通過判斷相鄰角點(diǎn)是否為地面點(diǎn)的方式不斷擬合新的曲面，同時不斷加入新點(diǎn)，擬合地形不斷變化，誤差不斷累積與傳遞就會導(dǎo)致濾波精度降低。

另外，該算法的閾值不可以隨著地形起伏的變化而自動調(diào)整，導(dǎo)致對坡度較大的地形特征點(diǎn)不能進(jìn)行有效分離；其次，該算法對場景中的所有點(diǎn)都是單獨(dú)進(jìn)行比較，會大大降低濾波效率。本文主要是對移動曲面擬合濾波算法中的誤差累積進(jìn)行改進(jìn)，改進(jìn)后的算法流程為[6]：

1)對機(jī)載LiDAR點(diǎn)云數(shù)據(jù)進(jìn)行預(yù)處理后，對點(diǎn)云區(qū)域進(jìn)行虛擬網(wǎng)格劃分，并為不同網(wǎng)格中的激光點(diǎn)建立高程序列索引。

2)確定初始窗口尺寸，在相鄰3×3窗口中選取6個高程最低點(diǎn)，構(gòu)成點(diǎn)集合；然后采用偏度校正的方式剔除虛擬網(wǎng)格中高程最低的非地面點(diǎn)，剩余的點(diǎn)組成種子點(diǎn)集。

3)使用最小二乘法通過種子點(diǎn)擬合點(diǎn)云數(shù)據(jù)區(qū)塊。

4)遍歷場景中所有點(diǎn)云，計(jì)算遍歷的激光角點(diǎn)的擬合高程與實(shí)際高程的差值Δhi，設(shè)置一定的閾值Δh。如果某個點(diǎn)的Δhi在閾值外，那么可認(rèn)定此點(diǎn)為非地面點(diǎn)，反之則為地面點(diǎn)，同時初步標(biāo)記地面點(diǎn)與非地面點(diǎn)。設(shè)置合適的閾值對點(diǎn)云濾波非常重要，直接影響濾波精度。本文閾值設(shè)定的方式為

Threshold=(Hmax-Hmin)×Scale

(2)

式中：Threshold為高差閾值，Hmax和Hmin分別為區(qū)域內(nèi)點(diǎn)云的最大高程與最小高程，Scale為窗口尺寸。

5)擴(kuò)大窗口尺寸，重復(fù)步驟3)和步驟4)，得到分類更加精細(xì)的激光角點(diǎn)。當(dāng)窗口尺寸達(dá)到閾值時，結(jié)束算法運(yùn)行，實(shí)現(xiàn)濾波。

基于偏度平衡的移動曲面擬合濾波算法的具體流程如圖1所示。

圖1 基于偏度平衡的移動曲面擬合濾波算法流程Fig.1 Flow Chart of Moving Surface Fitting Filtering Algorithm Based on Skewness Balance

3 試驗(yàn)驗(yàn)證與分析

3.1 數(shù)據(jù)來源

為了驗(yàn)證基于偏度平衡的移動曲面擬合濾波算法的可行性，本文采用國際攝影測量與遙感學(xué)會公布的CSite2首次回波數(shù)據(jù)及利用不同濾波方法得到的參考濾波數(shù)據(jù)作為試驗(yàn)數(shù)據(jù)對改進(jìn)后算法進(jìn)行驗(yàn)證與分析。

國際攝影測量與遙感學(xué)會提供的數(shù)據(jù)包含激光脈沖首次、末次的反射強(qiáng)度和反射位置，其中，Csite2的點(diǎn)密度可達(dá)到0.67點(diǎn)/平方米。

3.2 誤差指標(biāo)

濾波后的數(shù)據(jù)除了需要定性分析外，還需要使用誤差指標(biāo)對濾波結(jié)果進(jìn)行定量分析。目前常用的誤差指標(biāo)有兩種：

1)激光點(diǎn)的誤分率。誤分率可分為總誤差、一類誤差和二類誤差。其中，一類誤差主要是在進(jìn)行地面點(diǎn)分類時剔除部分地面點(diǎn)產(chǎn)生的誤差；二類誤差是在不能有效剔除地物點(diǎn)時產(chǎn)生的誤差；總誤差則是一類誤差與二類誤差通過加權(quán)求和得到的誤差。

改進(jìn)后算法的適用性可利用一類誤差和二類誤差進(jìn)行反映，通常是選擇有效的方法降低二類誤差，以得到更加全面的地形數(shù)據(jù)[7]。本文主要使用激光點(diǎn)的誤分率(表1)對改進(jìn)算法的效果進(jìn)行驗(yàn)證。

表1 濾波誤分率Tab.1 Filtering Error Rates參考結(jié)果濾波結(jié)果地面點(diǎn)個數(shù)地物點(diǎn)個數(shù)參考數(shù)據(jù)點(diǎn)個數(shù)誤分率/%地面點(diǎn)個數(shù)abe=a+bp= b / e地物點(diǎn)個數(shù)cdf=c+dq= c / f濾波后點(diǎn)數(shù)g=a+ch=b+dm=a+b+c+dr=(b+c)/m

表1中，a和d分別為在濾波過程中正確分類的地面點(diǎn)和地物點(diǎn)的個數(shù)，b為將地面點(diǎn)錯分為地物點(diǎn)的個數(shù)，c為將非地面點(diǎn)錯分為地面點(diǎn)的個數(shù)，e為正確地面點(diǎn)的個數(shù)，f為正確地物點(diǎn)個數(shù)，m為場景中所有的LiDAR點(diǎn)個數(shù)，p為一類誤差，q為二類誤差，r為總誤差。

2)在不同區(qū)域選取離散樣本點(diǎn)，計(jì)算濾波后地面模型的高程值與原始數(shù)據(jù)高程值的均方根誤差。將濾波前后高程值的均方根誤差和高程差值作為評價指標(biāo)對改進(jìn)算法的濾波效果進(jìn)行評價，均方根誤差越小，高程差值越小，表明濾波效果越好，精度越高[8]。

3.3 試驗(yàn)過程

為了更加直觀地反映基于偏度平衡的移動曲面擬合濾波算法的濾波效果，使用CSite2中的樣本數(shù)據(jù)Samp31和Samp23對改進(jìn)算法進(jìn)行驗(yàn)證。由于樣本數(shù)據(jù)Samp31和Samp23中既有原始數(shù)據(jù)又有濾波參考信息，因此可同時進(jìn)行定性和定量比對。

通過MATLAB軟件來運(yùn)行和實(shí)現(xiàn)改進(jìn)算法，并對傳統(tǒng)算法與改進(jìn)算法得到的試驗(yàn)結(jié)果進(jìn)行比對。在試驗(yàn)過程中，濾波效果受網(wǎng)格大小的影響較為嚴(yán)重，經(jīng)過試驗(yàn)比較，采用2 m×2 m的網(wǎng)格得到的效果較好。

在每個網(wǎng)格中選取3個高程最低點(diǎn)組成種子點(diǎn)集，并通過偏度平衡過濾掉網(wǎng)格中高程較高的點(diǎn)。由于網(wǎng)格面積較小，可能會覆蓋到建筑物表面，從而將建筑物點(diǎn)加入到地面點(diǎn)中，在進(jìn)行平面擬合時由于含有非地面點(diǎn)，導(dǎo)致誤差累積。剔除高程較大值后，得到所有最低點(diǎn)高程，對不同最低點(diǎn)密度進(jìn)行統(tǒng)計(jì)得到種子點(diǎn)直方圖，如圖2所示。

圖2 種子點(diǎn)高頻直方圖Fig.2 High Frequency Histogram of Seed Points

通過種子點(diǎn)進(jìn)行曲面擬合，以Samp31數(shù)據(jù)為例得到接近于地面點(diǎn)的擬合曲面如圖3所示。

圖3 擬合曲面Fig.3 Fitting Surface

使用改進(jìn)濾波算法與基于不規(guī)則三角網(wǎng)濾波算法對Samp31和Samp23試驗(yàn)數(shù)據(jù)進(jìn)行濾波，結(jié)果如圖4所示。

圖4 不同濾波方法試驗(yàn)結(jié)果Fig.4 Experimental Results of Different Filtering Methods

圖4(a)是采用TerraScan軟件的算法對Samp31數(shù)據(jù)進(jìn)行濾波處理得到的地面點(diǎn)集，可以看出其算法基本分離出了地面點(diǎn)，視覺上可以明顯看到地物點(diǎn)和地面點(diǎn)的區(qū)別，但實(shí)際是地面點(diǎn)中仍然含有少量地物點(diǎn)。

圖4(b)是采用基于偏度平衡的移動曲面擬合濾波算法對Samp31數(shù)據(jù)進(jìn)行濾波處理得到的地面點(diǎn)集，可以明顯看到空白區(qū)域?yàn)榻ㄖ锼趨^(qū)域，過濾掉存在的小型植被區(qū)域點(diǎn)云。圖4(c)是采用Terrasolid軟件的算法對Samp23數(shù)據(jù)進(jìn)行濾波處理得到的地面點(diǎn)集，圖4(d)是采用基于偏度平衡的移動曲面擬合濾波算法對Samp23數(shù)據(jù)進(jìn)行濾波處理得到的地面點(diǎn)集，可以看到提取效果明顯好于對Samp31數(shù)據(jù)進(jìn)行濾波處理的效果，但是仍存在一類誤差和二類誤差。其中，圖4(d)中的點(diǎn)云更能反映出地面區(qū)域特征，更明確地表達(dá)出建筑物及植被的信息。

3.4 結(jié)果分析

為了定量分析基于偏度平衡的移動曲面擬合濾波算法的濾波效果，本文使用誤分率這一評價指標(biāo)對濾波效果進(jìn)行評價。

由于點(diǎn)云數(shù)據(jù)量大，采用人工選取數(shù)據(jù)進(jìn)行濾波效果對比的方式不實(shí)際，因此，本文選擇與Csite樣本數(shù)據(jù)中標(biāo)記的參考濾波數(shù)據(jù)進(jìn)行對比的方式，分別使用移動曲面擬合濾波法和基于偏度平衡的移動曲面擬合濾波算法對Samp31數(shù)據(jù)和Samp23數(shù)據(jù)進(jìn)行濾波處理，并定量評價處理效果，結(jié)果如表2所示。

表2 濾波效果Tab.2 Filtering Effects 參考結(jié)果Samp31數(shù)據(jù)濾波結(jié)果地面點(diǎn)/個地物點(diǎn)/個參考數(shù)據(jù)點(diǎn)/個誤分率/%Samp23數(shù)據(jù)濾波結(jié)果地面點(diǎn)/個地物點(diǎn)/個參考數(shù)據(jù)點(diǎn)/個誤分率/%地面點(diǎn)/個15 23674513 1905.6516 34365713 1904.98地物點(diǎn)/個1 19912 18816 1057.441 05311 28416 1056.54濾波后點(diǎn)/個16 43512 93329 2956.6417 39611 94129 2955.84

根據(jù)表2可知，移動曲面擬合濾波法和基于偏度平衡的移動曲面擬合濾波算法的濾波效果都在10%以內(nèi)，一類誤差小于二類誤差，二類誤差和總誤差稍偏大，但都在可接受范圍內(nèi)，并且基于偏度平衡的移動曲面擬合濾波算法的一類誤差、二類誤差與總誤差均小于移動曲面擬合濾波算法。試驗(yàn)結(jié)果表明，基于偏度平衡的移動曲面擬合濾波算法適用于城市地區(qū)點(diǎn)云濾波，結(jié)果較為可靠，精度更高。

4 結(jié) 語

本文首先介紹了幾種常用的點(diǎn)云濾波方法，在分析移動曲面擬合濾波算法的基礎(chǔ)上，提出了一種基于偏度的移動曲面擬合濾波算法。該算法在確定種子點(diǎn)的過程中加入偏度概念，以降低種子點(diǎn)誤差，達(dá)到在上一次濾波完成之后繼續(xù)對非地面點(diǎn)和地面點(diǎn)進(jìn)行濾波，有效分離非地面點(diǎn)與地面點(diǎn)，從而提高濾波精度的目的。同時，使用國際攝影測量與遙感學(xué)會提供的數(shù)據(jù)對改進(jìn)算法進(jìn)行驗(yàn)證，結(jié)果表明，基于偏度的移動曲面擬合濾波算法的一類誤差和二類誤差都在10%以內(nèi)，且濾波效果優(yōu)于移動曲面擬合濾波算法，可以減少一類誤差、二類誤差和總誤差?；谄鹊囊苿忧鏀M合濾波算法適用于城市地區(qū)點(diǎn)云濾波，可取得較為可靠的濾波效果和更高的濾波精度。