譚 麒，原瀚杰，陳 亮，張 雨，何 勇，董麗夢，黃達文，陳澤佳

（1.廣東電網(wǎng)有限責任公司肇慶供電局，廣東 肇慶 526000）

點云濾波技術是機載激光雷達（LiDAR）數(shù)據(jù)預處理的重要環(huán)節(jié)之一，也是基于LiDAR數(shù)據(jù)生成和制作數(shù)字高程模型（DEM）的重要方法。點云濾波的核心功能是區(qū)分離散點云中的地面點和非地面點。目前國內(nèi)外已有眾多學者針對該問題開展了大量研究工作，主要可分為基于形態(tài)學、曲面擬合、圖像分割、不規(guī)則三角網(wǎng)（TIN）和深度學習等技術的濾波方法?；谛螒B(tài)學的方法將點云轉(zhuǎn)換為柵格化高程數(shù)據(jù)，根據(jù)局部區(qū)域內(nèi)真實地形的幾何特征，利用相應的圖像處理方法進行濾波，具體包括漸進形態(tài)學[1]、多尺度表示[2]、圖像補全[3]、Hermite變換[4]、頂帽變換[5]等?；谇鏀M合的方法利用一個連續(xù)的參數(shù)表面對真實地形進行逼近，通過測量點與擬合曲面的接近程度確定該點的類別，根據(jù)函數(shù)模型的不同可分為線性預測[6]、薄板樣條[7]、二次曲面[8]等方法?；趫D像分割的方法同樣將點云進行柵格化處理，采用先聚類再判斷分割塊屬性的策略進行分類[9-10]，該過程中除直接利用高程值外，還可能采用回波數(shù)量[11]、最大回波高差[12]等信息作為額外判據(jù)?；赥IN的方法通常利用局部最低點構建初始TIN，然后通過逐步優(yōu)化TIN實現(xiàn)濾波，其中具有代表性的是Axelsson P[13]提出的漸進加密三角網(wǎng)（PTD）算法；在此基礎上，一系列改進的PTD算法得到發(fā)展，具體改進策略包括引入平滑約束[14]、增加山脊點以保留山體特征[15-16]、增加種子點[17]、引入支持向量機[18]或布料模擬技術[19]等?；谏疃葘W習的方法主要通過近年來計算機視覺領域發(fā)展迅速的卷積神經(jīng)網(wǎng)絡等技術實現(xiàn)濾波[20-22]，但該類方法通常對訓練樣本的標注質(zhì)量和數(shù)量有較高要求。

作為目前生產(chǎn)實踐中應用最廣泛的濾波方法之一，經(jīng)典PTD算法通常能在預先知曉數(shù)據(jù)地形類別、進而更改對應參數(shù)配置的前提下發(fā)揮良好性能；但對于平地/山區(qū)混合出現(xiàn)、地形類別難以明確定義的場景，利用同一套參數(shù)的適應性將可能下降，導致濾波結果存在山脊特征缺失或平坦地區(qū)凸起等問題。鑒于此，本文提出了一種基于地形預分類的改進PTD算法，通過二次濾波的思想，實現(xiàn)對混合地形點云數(shù)據(jù)的自適應高精度濾波。首先采用統(tǒng)一參數(shù)進行初次濾波，并利用濾波得到的地面點構建DEM，以DEM中的高程信息為引導進行地形預分類，得到平地和山區(qū)的地理范圍；然后針對不同地形區(qū)域的點云數(shù)據(jù)，分別采用不同的參數(shù)配置進行二次濾波；最后對不同地形下的二次濾波結果進行合并，得到最終濾波結果。

1 研究原理與方法

1.1 經(jīng)典PTD算法存在的問題

經(jīng)典PTD算法通過迭代加密三角網(wǎng)的方式，在地面點密度足夠的情況下，可保留地形的細節(jié)特征，在一定程度上能適應地形的復雜變化，但仍存在局限性。對于包含地勢起伏較大的山區(qū)和含有大型建筑物的平地區(qū)域等混合地形只使用一套閾值參數(shù)難以取得較好的結果，植被覆蓋的斜坡地形中一些小的起伏都可能導致地面點與非地面點的混淆，選取寬松的閾值可能導致山區(qū)無法檢測出完整的地形特征，而選取保守的閾值又可能使平地區(qū)域的建筑物被誤分為地面點。

1.2 基于地形預分類的改進PTD算法

針對經(jīng)典PTD算法的上述問題，本文提出了一種基于地形預分類的改進PTD算法。具體算法流程如圖1所示。

圖1 技術流程圖

1.2.1 地形預分類

本文將地形分為平地和山區(qū)，二者最顯著的區(qū)別是地勢起伏度、高程信息和坡度信息，因此以地勢起伏和高程信息為主，坡度信息為輔，對整體區(qū)域進行地形分類。地形預分類的具體步驟為：

1）初次濾波。首先利用適中的濾波閾值進行初次濾波，得到初始DEM；再根據(jù)DEM計算每個點的地勢起伏信息。地勢起伏可表示為在某點的一定范圍內(nèi)，高程最大值與最小值的差值。由于山區(qū)地勢起伏度大、高程大；平地區(qū)域地勢起伏度小、高程小，因此利用高程和地勢起伏信息點位相加后的數(shù)據(jù)可顯著區(qū)分因特征丟失而平整的山脊區(qū)域和平地區(qū)域。對于疊加后的數(shù)據(jù)，計算統(tǒng)計直方圖，遍歷整個直方圖，尋找兩個峰值之間的最低點獲取最優(yōu)閾值，并進行二值化。

2）數(shù)學形態(tài)學和區(qū)域生長。若僅根據(jù)上述數(shù)據(jù)進行分類，無法很好地區(qū)分平地和山區(qū)，且分類區(qū)域破碎，綜合利用數(shù)學形態(tài)學的腐蝕和膨脹可有效解決該問題。首先利用腐蝕形態(tài)學減緩邊界的擴張，有利于后期邊界處理；然后利用閉運算形態(tài)學填充物體內(nèi)的細小空洞，連接鄰近的物體，該過程中可調(diào)整窗口大小，以達到除去大塊區(qū)域中破碎分類區(qū)域的目的，有效解決大塊區(qū)域中的空洞問題。在山區(qū)的山谷部分，因其地勢平坦，會被誤分為平地，于是利用區(qū)域生長，對區(qū)域內(nèi)的小塊區(qū)域進行面積統(tǒng)計，小于閾值且其周圍存在大量不同類別數(shù)據(jù)的，則改變其類別，從而實現(xiàn)山區(qū)和平地的分類。

3）坡度數(shù)據(jù)輔助邊界處理。經(jīng)過上述流程，在平地與山區(qū)的交界處，仍存在嚴重的誤分現(xiàn)象。地形交界處坡度會急劇變化且山區(qū)坡度較大，因此可利用坡度數(shù)據(jù)減少這種誤分現(xiàn)象。坡度較大的點可認為是山區(qū)地形。

1.2.2 不同地形的濾波參數(shù)選擇

針對地形預分類后的點云，本文根據(jù)地形類型進行閾值設置：

1）平地區(qū)域。地面點選擇窗口大小不能太小，否則無法避免較大人工建筑物的影響；角度與距離閾值不能太大，對于平地而言，初始TIN就能較完整的表達，若角度與距離閾值過大會使原本平坦的地面出現(xiàn)凹陷與凸起現(xiàn)象，將誤分一些低矮地物，這些誤分點在平地上更為明顯。

2）山地區(qū)域。地面點選擇窗口大小不能太大，否則地面點密度不足，難以完整表達地形；角度與距離閾值不能太小，山區(qū)地形崎嶇，原本相鄰地面點的高程就相差較大，若閾值太小會使最終地面點密度過低，尤其是在大面積植被覆蓋的山區(qū)。

1.2.3 不同地形濾波結果的合成

分地形完成濾波后，需對濾波后的結果進行合并。由于分區(qū)域的點云是從同一個點云數(shù)據(jù)集中分離出來的，因此它們擁有相同維度，且不存在重疊區(qū)域，點云合并只是點云數(shù)量上的合并，且保留其地面點和非地面點屬性，也可只保留地面點的點云，便于后續(xù)DEM的生成。

2 實驗結果與分析

2.1 實驗數(shù)據(jù)與參數(shù)選擇

本文選取安徽省黃山市休寧縣的一個山地/居住地混合地區(qū)作為實驗區(qū)域，采集的LiDAR點云數(shù)據(jù)覆蓋面積約為1 km2，平均點密度約為2.0點/m2，高程區(qū)間為204～430 m。實驗區(qū)地勢起伏大、山脊陡而窄、植被覆蓋度高，植被下方地面點相較于裸露地表更稀疏；中間地帶地勢相對平坦，存在集中的居民地和河流，建筑密集但無單個超大面積建筑物。實驗數(shù)據(jù)的整體情況如圖2所示，圖2a中白色矩形內(nèi)為實驗區(qū)可見光遙感影像；圖2b為點云數(shù)據(jù)，激光點的高程由藍至紅逐漸升高，紅框內(nèi)為包含多個建筑物的居民區(qū)。本文對實驗數(shù)據(jù)進行人工編輯，生成了1 m分辨率的DEM，作為參考數(shù)據(jù)用于實驗的定量評價。

圖2 實驗數(shù)據(jù)

實驗數(shù)據(jù)的地形預分類結果如圖3所示，圖3b中白色區(qū)域表示山區(qū)、黑色區(qū)域表示平地。通過與圖2的比較發(fā)現(xiàn)，本文的預分類算法整體上取得了較好的地形區(qū)分效果。

圖3 去噪前后殘差值對比圖

圖3 原始數(shù)據(jù)與預分類結果

本文設置了3組對比實驗參數(shù)：經(jīng)典PTD平地參數(shù)、經(jīng)典PTD山地參數(shù)、本文方法。對平地和山區(qū)分別采用平地參數(shù)和山地參數(shù)，不同方案的參數(shù)設計如表1所示。

表1 不同方案的參數(shù)設置

2.2 濾波結果與分析

利用表1中的參數(shù)分別進行3次試驗，并利用濾波后的結果生成1 m分辨率的DEM，計算其與參考DEM之間的高程差異；再利用均方根誤差（RMSE）對結果進行評價。評價結果如表2所示，可以看出，采用本文方法的濾波精度比經(jīng)典PTD平地參數(shù)方法提升了0.551 m，而僅比經(jīng)典PTD山地參數(shù)方法提升了0.129 m，較小的提升很大程度上與試驗數(shù)據(jù)有關，該數(shù)據(jù)大部分區(qū)域為山區(qū)，因此利用本文方法精度提升不顯著。

表2 RMSE評價結果

3種濾波方法結果的TIN渲染圖如圖4所示，可以看出，在經(jīng)典PTD平地參數(shù)下，平地區(qū)域的濾波結果達到較好效果，但山地部分特征缺失嚴重，尤其是山脊部分缺失嚴重，導致山脊平整；在經(jīng)典PTD山地參數(shù)下，更好地保留了山地部分特征，但平地區(qū)域凸起現(xiàn)象明顯，尤其是建筑物區(qū)域；本文方法既保留了完整的山地特征，又兼顧了平地區(qū)域，達到了最好的效果。

圖4 濾波結果TIN渲染全局對比

3種濾波方法結果的局部TIN渲染結果對比如圖5所示，分別選取山區(qū)和平地區(qū)域進行對比，可以看出，經(jīng)典PTD山地參數(shù)方法和本文方法在山地區(qū)域的效果較好，山脊特征完整，而經(jīng)典PTD平地參數(shù)方法的山脊特征嚴重缺失，如圖5中紅框所示；經(jīng)典PTD平地參數(shù)方法和本文方法在平地區(qū)域相對平滑，沒有明顯的突起或凹陷，而經(jīng)典PTD山地參數(shù)方法存在明顯的突起現(xiàn)象，如圖5中藍框所示。

圖5 濾波結果TIN渲染局部對比

高程差值分布如圖6所示，可以看出，采用本文方法獲得的結果在平地區(qū)域和山地區(qū)域的濾波精度均有明顯提升，提高了PTD濾波的整體精度。

圖6 濾波結果高程差值圖

3 結語

本文在經(jīng)典PTD算法的基礎上，提出了一種基于地形預分類的改進PTD算法，通過利用和提取初次濾波結果中包含的高程信息實現(xiàn)地形預分類，進而有針對性地應用不同參數(shù)配置對不同地形進行濾波處理，以提高經(jīng)典PTD算法處理混合地形數(shù)據(jù)的能力。實驗結果表明，對于平地/山區(qū)混合地形數(shù)據(jù)，本文方法能較完整地保留不同地形的地表特征，所生成DEM的精度顯著優(yōu)于經(jīng)典PTD算法。