激光雷達技術和攝影測量方法在無人機數字測圖中的精度對比分析

李 濤，常 江,胡東升，廉旭剛，呂俊沛

(1.華陽新材料科技集團有限公司，山西 陽泉 045000；2.太原理工大學 礦業工程學院，太原 030024)

目前，無人機搭載激光雷達模塊采集點云數據，獲取地表高精度三維坐標，相比數碼航空攝影測量有獨特優勢[1]。首先，激光雷達直接獲取地物表面高精度三維坐標，不受山體、建筑物陰影的影響，單純數碼航空攝影測量在陰影區域判讀困難。其次，激光雷達的多次回波對植被具有較好的穿透性，而植被覆蓋一直是數碼航空攝影不能解決的問題。最后，攝影測量提高精度需經過后差分技術[2-3]和SFM技術[4-5]，即通過設計控制點的增多與布設位置的調整，通過拍攝多視角照片結合特征匹配算法，DEM精度可以達到±9 cm，而原始激光點云即可達到±10 cm的高程精度。同時,有多篇文章提出Lidar獲取高精度DEM的關鍵技術，大致可以分為兩類，即通過設計改進濾波算法和通過調整采集方案。例如,提出過綠減過紅指數提取植被[6-7]、距離限制濾波提取地面點[8]等改進算法；提出從設備選擇、點云密度設計、植被覆蓋密集山區數據獲取方法、點云數據分類組合算法、空白區處理等方面[9-10]進行探討,并提出改進方法。

目前，對一種點云數據(Lidar點云或者攝影測量點云)精度研究比較多，但是缺乏Lidar點云和攝影測量點云對同一片地區DEM的對比誤差研究，定量化描述Lidar點云和攝影測量點云的誤差缺乏方法也缺乏實際應用案例。本文提出一種點云對比方法，并結合某地區的實際數據進行案例演示，為后續此類型研究提供參考方案。

1 研究區概況和數據準備

1.1 研究區概況

研究區位于呂梁市孝義市驛馬鄉下荊封村的某個廠區。地形表面復雜，溝谷縱橫，屬大陸性半干旱氣候。研究區大部分地區植被茂密，除了一條山路和廠區外，只有少量地表裸露，可以滿足不同地貌條件之間高程數據的對比。

地理位置東經111.604 2°—111.608 1°，北緯36.957 2°—36.960 0°之間，高程最低947 m，最高992 m，研究區面積大約0.055 km2。如圖1所示。

圖1 研究區概況

1.2 數據準備

數據采集使用飛馬D2000無人機航測平臺，一次飛行搭載D-LIDAR2000激光雷達模塊，重疊度設置為25%。另一次飛行搭載索尼 A6000相機，有效像素2 430萬，航向和旁向重疊分別為80%和60%，采用網絡RTK/PPK高精度融合POS+免像控方案。數據采集時間是2021年9月20日，搭載索尼 A6000相機采集的數據，可滿足高精度1∶500地形圖精度需求。搭載D-LIDAR2000激光雷達模塊采集的數據可以直接用于本次的研究。雷達模塊參數如表1所示。航測模塊參數如表2所示。

表1 D-LIDAR2000激光雷達模塊的參數

表2 SONY a6000航測模塊的參數

2 點云數據處理

2.1 數據修正

在本次研究中，首先需要清除數據在采集過程中或者使用相關軟件處理過程中產生的異常點、孤立點[11]。使用相關濾波算法時，要注意清除低點(不是正常的低矮地物)。使用點云分類工具，設置相關參數分離低點和孤立點，為了保證清除完整，可以多次執行該命令。特殊噪點需要手動清除，使用分類工具classify above line，從側視圖中觀察發現噪點。

2.2 數據分層

對已經清除噪點的兩份點云數據進行點云分類，兩份數據分別劃分為植被層、建筑物層、地面層(地面層生成DEM)，主要使用點云分類工具，設置不同的參數，保存對應層的數據。兩份數據提取相同層的時候盡量保證參數一致，使得兩份數據相同層的點云在數量和質量上保持一致。每層數據的保存需要規定格式，方便后續的分析處理。

1)攝影測量點云提取植被層。使用點云分類工具的by vegetation index功能，因為攝影測量點云是帶有顏色(RGB值)的，其原理[6]主要是通過過綠減過紅指數(EG-ER)把RGB值當參數，通過計算結果和閾值來判斷是否是綠色，再將密度值和高程差值作為輸入特征變量[7]，使用支持向量機(SVM)算法做點云分類，來識別是否是植被。設置的min value和max value含義是規定綠色的范圍，目標找到這范圍內的綠色點。如圖2所示。

(a)植被提取前

EG-ER作為一種改進的顏色指數通過將EG指數圖像與ER指數圖像相減，發現基于閾值的EG-ER指數可以較好地將植被與背景分離，計算公式如下：

EG-ER=3g-2.4r-b.

(1)

其中ER為過紅指數，計算方法如下：

ER=1.4r-g.

(2)

上式的R,G,B值是歸一化后的值，即0～255歸一到0～1之間。

2)點云提取地面點層制作DEM。使用點云分類工具的ground功能，該功能的原理是發現點云之間的距離和角度的關系[12]。因為不同地物構建的三角網在角度、邊長、高度上有差異，通過設置相關參數，排除地表的植被和建筑，只留下地面點。分離地面點的原理有兩個，距離限制濾波[8]分類地面點和非地面點，角度濾波去除地物的側面信息。然后依據局部地形設置動態閾值，以表面擬合區域生長算法擴充地面種子點，循環迭代逐漸逼近真實地面[13]。

濾除起伏突然比較大的點，比如建筑物點，植被點等，使用角度濾除，一般計算θ值的大小，計算公式如下：

獲得地面種子點之后，便可以進行表面擬合區域生長算法擴充地面種子點，其流程可以描述為：

a.利用地面種子點構建TIN；

b.選擇一個未分類點，尋找該點所在的三角形；

c.計算該點到三角形平面的距離S和該點與三角形3個角點的夾角，并選出最大角度，其幾何意義如圖3所示；

圖3 S和θ的幾何意義

d.根據人工設置的最大距離閾值Smax和最大角度閾值θmax，再規定區間，將該點加入地面點；

e.重復步驟b～d直到所有未分類的點全部判斷完畢。

上述原理在提取過程中體現為設置不同的參數，Max building size設置需要根據實地情況填寫，必須大于當地的最大建筑面積。Iteration angle和Iteration distance兩個參數是三角網的參數，代表角度和距離。城區設置一般偏小，比如4°和0.6 m；山區設置一般偏大，比如6°和1.4 m。如果得到過濾點云密度太低，需要選擇加密點云ground層，因為ground層的點云密度越大，使用Export lattice Model工具插值得到規則格網DEM越精確。

2.3 DEM點云數據處理

按照2.2描述進行數據分層，一共得到6個文件，分別是Lidar點云的building.txt，tree.txt，DEM.txt和攝影測量點云的building.txt，tree.txt，DEM.txt。然后主要使用python語言，寫代碼實現數據的分析處理。

主要思想是，由于兩份數據的坐標系相同，同一地區范圍得到的同一種地物類型的點云，必定存在相同坐標點，兩幅點云的相同坐標點的高程差應該是差異較小。但是由于數據獲取方式不同、數據處理方式不同和數據類型不同，兩幅點云的高程存在不容忽視的差異，所以需要通過以下方法量化分析二者差異。如圖4所示。

圖4 算法實現方法圖

使用Lidar點云和影像點云的DEM.txt來舉例：通過上述處理，得到Lidar點云大約20萬個點，攝影測量點云17萬個點，使用python的pandas庫可以實現找相同坐標點的操作，并且可以最大程度保留原始數據[14-15]。經過上述操作，保留大約16萬對相同坐標點，是原始數據Lidar點云數量的80%，是原始數據攝影測量點云的95%。然后對新表的列名重新命名，進而計算Lidar點云高程減去攝影測量點云高程的差值。如表3所示。

表3 兩份DEM的交集數據以及高程差

從Lidar-photo這一列中，可以看到，兩份DEM的高程差是普遍存在的。由于軟件處理植被點云可能不完整，或者由于無人機飛行初始和結尾狀態的不穩定性[9]，導致出現的高程差比較大，稱之為粗差點，這些點是需要排除的。

選擇剔除的方法一般是認為兩倍標準差是數據的合理范圍，大于兩倍標準差的可以認為是粗差，選擇剔除。如圖5所示。

圖5 箱型圖觀察交集DEM排除大于兩倍標準差的點云變化

箱型圖中黑色圓圈點表示大于1.5倍標準差的數據，矩形方框表示數據的75%集中的區間，綠色的線表示平均值。下表是處理前后的關鍵統計量的變化，平均值變化較小，進一步說明排除大于兩倍標準差的方法是合理的。

經過上述處理，剔除5%的粗差點云，得到的點云數量占原來的95%。如表4所示。

表4 統計量觀察交集DEM排除大于兩倍標準差的點云變化

從數據中可以看出，已經把處理前最大高差19 m的數據排除，但存在最大高差5 m或者4 m，而選擇保留是因為在實際中可能存在5 m或4 m的高差。現在兩份DEM高差數據大部分的高差有幾十厘米甚至更小，從第一四分位(25%)為-0.583 75 m和第三四分位(75%)為0.147 m可以看出。因此，待研究的數據整體是比較合理的，研究結果是具有參考意義的。

2.4 植被和建筑點云數據處理

數據處理辦法和2.3一樣，但是需要注意因為Lidar點云的穿透性，是存在同一個x、y坐標有多個高程值，所以第一步需要對保存的植被層刪除重疊點，保留重疊的多個點的最高值，把最高值認為是植被高度。對于建筑物點云的重疊情況，選擇重疊點云的平均值作為建筑物的高程。

可以看到,在處理重疊點后，兩次輸出的點云數量是不一樣的。比如，植被層，第一次輸出9 429 398個，第二次輸出9 352 880個，重疊率0.81%。

3 數據對比分析

3.1 平均絕對誤差(EMA)的比較

平均絕對誤差，全稱是Mean Absolute Error，EMA的值越小，說明影像高程數據越靠近Lidar高程數據。EMAP[17]是EMA的變形，它是一個百分比值，表示偏離真實值的程度，其值越小則兩份數據之間差異越小。把之前數據處理得到的植被的2 040個點、建筑物的205個點、DEM的153 238個點輸入下列公式中，計算結果如表5所示。

表5 植被、建筑、DEM的EMA和EMAP的對比表

數據結果說明，不論影像點云還是Lidar點云，建筑物的精度是最高的。植被和DEM差異比較小，但是由于DEM是地面點合成的，從物理的角度來說確實應該比植被的精度要高，所以可以說明這份數據處理過程是合理的，其對比分析具有一定參考價值。

3.2 均方根誤差(ERMS)的比較

均方根誤差(ERMS)反映測量數據偏離真實值的程度，ERMS越小，表示測量精度越高。把Lidar點云數據當作真值，把影像點云數據當作觀測值。同樣把之前數據處理得到的植被的2 040個點、建筑物的205個點、DEM的153 238個點輸入到上述代碼中，得到以下計算結果，如表6所示。

表6 植被、建筑、DEM的ERMS的對比表

從ERMS結果中仍然可以看出，建筑物點云的數據精度最高，和其他兩種數據的差異是比較大的，植被和DEM比較接近，DEM稍優于植被。所以，再次證明不論Lidar點云還是影像點云，在做精度對比的時候一定要先將地物進行分類。

3.3 Lidar點云和影像點云的DEM驗證

可以明顯看到在平坦或者規則地區，Lidar的DEM和影像的DEM高度重合，但是在山地或者植被覆蓋的地方二者差異還是比較明顯的，如圖6所示。

(a)硬化地表DEM點云

4 結論

1)當不做點云分類時得到的均方根誤差(ERMS)達到724 mm，而做過點云分類之后，植被層的ERMS可以達到1 138 mm，相比不做點云分類的724 mm漲幅57%；但是建筑物層的ERMS計算得到108 mm，降幅達到85%；DEM通過計算為627 mm，降幅13%。可見精度對比分析第一步，必須先分類，再對比分析，使用文中點云數據處理方法，得到相關結論。

2)經過傳統攝影測量得到的影像點云，按照文中的分類方法，分類得到的硬化地表地區的點云，其精度是可以認為比較接近Lidar點云，所以在建筑不太密集的城區，可以繼續發揮傳統攝影測量在成本和效率方面的優勢。

3)在實際應用案例中，經過該方法處理后，DEM點云數量損失為5%，而植被和建筑物點云數量損失稍大，可以說明該對比方法的合理性；而且數據處理過程中，在同一坐標系下多次找兩份點云數據的相同坐標點，確保對比方法的嚴密性，得到的結論具有實際參考意義。