基于KD樹的規則格網DEM插值技術

黃 艷，盧學良

1.西安測繪研究所，陜西 西安，710054;

2.地理信息工程國家重點實驗室，陜西 西安，710054

1 引 言

數字高程模型(digital elevation model，DEM)作為一種表達地形的數學模型，自1958年由美國麻省理工學院Miller教授提出以來，由于其表達方式直觀多樣、便于更新等優勢，已經在測繪導航、農林規劃、地學分析等領域得到了廣泛的應用和研究，成為國家經濟建設以及國防建設的重要基礎數據之一[1]。DEM的生產途徑多種多樣，目前應用最為廣泛、生產效率最高的方法一般是通過航空航天攝影測量、合成孔徑雷達干涉測量或是機載激光掃描等方式獲取地形地物的點云數據，然后通過內插的方式獲得規則格網的DEM[2，3]。在這些方法中，插值是DEM生產的重要環節，常用的DEM插值算法有反距離加權插值算法、謝別德插值算法和徑向基函數插值算法等[4，5]，這些算法無不是利用待插值點鄰域范圍內的已知采樣點通過特定的加權方式完成插值。因此，鄰域采樣點的獲取是完成插值的前提，而為了高效獲取待插值點的鄰域采樣點，一般要對點云建立空間索引，常用的空間索引方法為柵格法。本文在深入分析柵格法索引的基礎上，引入了KD樹來實現鄰域點的快速獲取，并通過試驗對兩種空間索引的效果進行了對比分析。

2 DEM插值

DEM插值是根據已知采樣點的高程值去估計未知插值點高程值的過程[6]，其主要用于缺值估計、等值線內插和離散點云數據的格網化[7]。DEM插值可以通過已知采樣點來估計未知插值點，主要是因為地形具有空間異質性和空間相關性等基本特征，使得利用一些空間位置合理的采樣點(一般為鄰域采樣點)獲得對地形表面相對精確的描述成為可能。

2.1 鄰域采樣點的獲取

(1)柵格法索引及其存在的問題

在對大數據量的地面點云進行插值時，插值點鄰域采樣點的獲取一般不宜采用全局搜索的方法。為了提高搜索效率，首先要建立點云數據的空間索引，常用的空間索引方法為柵格法。

柵格法的基本思想如圖1所示，首先輸入點云數據，在確定點云數據的平面范圍后，將該平面區域按照設定的尺度柵格化，柵格單元一般為正方形，然后將各個小柵格與其內包含的點關聯起來。對于待插值點，在確定該點所在的柵格后，通常該柵格及其8鄰域柵格所包含的點即為待插值點的候選鄰域采樣點，最后經過距離約束或其他方式即可篩選出鄰域采樣點。

圖1 柵格索引建立及檢索流程

利用柵格法建立空間索引主要有兩個要點：一是柵格大小的確定;二是待插值點搜索范圍的確定，即待插值點所在柵格的鄰域柵格的范圍。其中柵格的大小與DEM的分辨率有關，一般是DEM分辨率的若干倍。而鄰域柵格的選擇范圍一般是其8鄰域，當然也可以向外擴充，原則上是柵格越大，鄰域的選擇范圍越小，柵格越小，鄰域的選擇范圍越大。

然而利用柵格法作為空間索引，可能存在某些點比部分鄰域采樣點距離待插值點更近的情況。如圖2所示，顯而易見的是doutside＜dinside，這樣很可能會導致計算出的待插值點的若干個最近鄰點不夠準確(并非最鄰近點)，進而會影響最終的插值精度。

圖2 鄰域柵格的內點與外點距離示意

目前有兩種辦法來解決上述情況。第一種方法是在DEM插值點設定某個距離閾值d，如果某些最鄰近采樣點與待插值點之間的距離大于d，那么它將會被排除。因為距離較遠的采樣點可能并不能提高待插值點的插值精度，甚至會有所損害，基于這點考慮，將柵格的邊長設置為d，便可以規避上述問題，但這樣會明顯降低空間索引的效率。另一種方法是如果理想中的鄰域柵格的選擇范圍是其8鄰域，在實際處理中把選擇范圍擴大到24鄰域，可以在很大程度上避免上述問題(并不能完全規避)，但同樣會降低空間索引的效率。

(2)KD樹索引

鑒于上述問題，本文在空間索引上使用了KD樹。KD樹(K-Dimension tree)最早由斯坦福大學的Jon Bentley教授在論文中提出[8]。它是一種應用于高維空間的數據索引結構，采用分而治之的思想將整個空間數據集在K維空間按照一定的規則順次劃分為若干個小空間。二維KD樹的構建如圖3所示。首先紅色的線條將點集分成兩個部分，然后每個子塊被藍色線條分別分成兩個部分，接著按照綠色線條、紫色線條、淺藍色線條的順序依次將子塊劃分，直至子塊中的點集規模小于局部規模閾值。KD樹實際上就是決策樹，每次找區分度最大的一條軸線進行二分，二分至葉子節點只包含一個點。在查找最鄰近點時，先找到葉子節點，然后回溯其祖先節點，以這個點的距離剪枝搜索，最終找到符合要求的若干個最鄰近點。

圖3 二維KD樹[9]

2.2 常用的插值算法

目前常用的插值算法有反距離加權插值算法、謝別德插值算法和徑向基函數插值算法等。

(1)反距離加權插值算法

反距離加權插值算法是根據相近相似的原理，對每個采樣點賦予一定的權重，權重隨著采樣點與插值點之間距離的增加而減小。任一插值點處的值是各采樣點加權之和，如式(1)所示。

其中，zp為插值點的高程值;λi為第i個點的權重;di為第i個采樣點到插值點的距離;d-u為距離衰減函數，冪指數u具有隨著距離的增加而減小其他位置影響的作用，通常取值為1或2。

(2)謝別德插值算法

謝別德插值算法本質上是一種標準的導數距離加權過程，權函數如式(2)所示。

其中，wi為權重;di為第i點距待插值點的距離;r為調整距離。

(3)徑向基函數插值算法

徑向基函數插值算法是一系列用于精確插值算子的統稱，其插值原理是任何一個表面都可以用多個曲面的線性組合去逼近。通常情況下，徑向基函數插值算法可以表述為兩部分之和，如式(3)所示。

其中，zp為插值點的高程值;λi為i第個點的權重;di為第i個采樣點到插值點的距離;φ(di)為徑向基函數;fi為“趨勢”函數，是次數小于m的基本多項式函數。

3 試驗與分析

3.1 試驗方案

試驗一：通過對前述兩種空間索引的對比分析，測試空間索引效果。試驗環境為Intel(R)Xeon(R)CPU E5-2630 v3(2.4GHz)，32G內存，操作系統是Window7 64位。試驗數據采用某區域5m分辨率的衛星影像通過密集匹配并濾波后得到的點云數據。點云數據有A、B兩塊，A區域的數據量為：356MB(文本文件)，含有10989051個點;B區域的數據量為：1.14GB(文本文件)，含有36551055個點。DEM的采樣分辨率設置為15m，最鄰近點的檢索數量設置為4個，距離閾值d設置為60m(即4倍的DEM采樣分辨率)。主要統計的指標有最鄰近點查詢的正確率、空間索引的構建時間以及DEM插值計算時間(鄰域點的搜索與插值)。對于柵格法由于不同參數的設定會產生不同的試驗結果，這里設定k為柵格大小與DEM分辨率的比率，n為待插值點所在柵格的鄰域柵格的范圍(n鄰域)。

試驗二：利用試驗一中表現效果較好的空間索引，通過反距離加權插值算法完成A區域的DEM的采樣，并通過目視評估采樣效果。

3.2 試驗分析與結論

試驗一的結果見表1。

表1 不同空間索引的統計指標

通過表中的指標統計，可以看出，點云數量和插值點數的比例并非理論上的9∶1(DEM分辨率和影像分辨率的比例為3∶1)，這是因為實際插值的范圍并非恰好是點云數據的平面分布范圍，而是矩形包圍盒的范圍，因此，插值點的數量會比理論上要多。對于柵格法，空間索引的構建時間與點云數量成正比，與柵格的大小沒有關系，而插值的計算時間與檢索的鄰域范圍成正相關關系，與插值點數成正比。當k=2、n=8時，正如前述的分析的那樣，并不能保證最鄰近點檢索的正確性，而對于k=2、n=24和k=4、n=8的情況，試驗結果與理論相符(與距離閾值d的設定有關)，能保證最鄰近點檢索的完全正確。對于KD樹，空間索引構建時間與點云的數量成正比，插值計算時間與插值點數成正比，并且KD樹能保證檢索的最鄰近點完全正確。對比兩種空間索引，在保證正確率的前提下，可以得出的結論是KD的效率要優于柵格法。

根據上述結論這里選取KD樹作為空間索引進行試驗二。試驗二的結果如圖4所示。

圖4 DEM采樣效果

圖4中左上角為原始點云，左下角為采樣后的DEM，右側分別為原始點云和采樣后DEM的局部放大圖，從圖中可以看出基于KD樹的DEM插值的效果較為理想。

4 結 語

本文對DEM插值中常用的空間索引進行了詳細的介紹和分析，指出了其中存在的問題，并在此基礎上引入了KD樹。試驗證明，KD樹在效率方面具有優勢，并且基于KD樹的DEM插值的效果很好。當然柵格法的特點是實現簡單、靈活，對于不同的數據特點、參數設置和成果要求，柵格法同樣可以加以運用。

[1]張亞南.DEM分辨率確定與尺度轉換方法研究[D].南京：南京師范大學，2014.

[2]張海霞，陳向陽，趙文普等.高分辨率遙感影像自動生成DEM的研究[J].測繪與空間地理信息，2015(7)：103-105.

[3]湯國安.我國數字高程模型與數字地形分析研究進展[J].地理學報，2014，69(9)：1305-1325.

[4]張錦明.DEM插值算法適應性研究[D].鄭州：信息工程大學，2012.

[5]王耀革.DEM建模與不確定性分析[D].鄭州：信息工程大學，2009.

[6]盧華興，劉學軍，湯國安.DEM誤差模型研究[C].鄭州：第三屆地理信息系統全國博士生學術論壇，2008.

[7]李新，程國棟，盧玲.空間內插方法比較[J].地球科學進展，2000，15(3)：260-265.

[8]Bentley JL.Multidimensional Binary Search Trees Used for Associative Searching[J].Communications of the ACM，1975，18(9)：509-517.

[9]Muja M，Lowe D G.Scalable Nearest Neighbor Algorithms for High Dimensional Data[J].IEEE Transactions on Pattern Analysis&Machine Intelligence，2014，36(11)：2227-2240.