趙文豪,閆 利,張云生,陳斯煬

(1.武漢大學(xué) 測繪學(xué)院, 湖北 武漢 430079; 2.國家基礎(chǔ)地理信息中心, 北京 100830； 3.中南大學(xué) 地球科學(xué)與信息物理學(xué)院, 湖南 長沙 410083 )

隨著激光掃描技術(shù)和多視影像密集匹配方法的發(fā)展，獲取高密度的點云越來越容易，為了減少冗余和提高數(shù)字高程模型(DEM)處理、網(wǎng)絡(luò)傳輸?shù)葢?yīng)用的效率，有必要對DEM進行綜合以獲取多尺度DEM數(shù)據(jù)庫[1]。目前的DEM綜合方法按照DEM存儲方式可以大致分為2類：面向三角網(wǎng)(TIN)的綜合方法和面向格網(wǎng)DEM的綜合方法[2]。其中3維Douglas-Pucker (DP) 算法是由費立凡等人根據(jù)經(jīng)典2維DP算法發(fā)展而來的DEM綜合算法，最早用于處理離散點，也就是TIN表達形式的DEM[3-4]。由于原始算法處理海量數(shù)據(jù)時效率較低，何津在利用3維DP算法綜合規(guī)則格網(wǎng)DEM時，利用行列掃描的方式來減少計算量[5]。該方法大幅提高了效率，但可能會出現(xiàn)漏掉重要地形點的情況。根據(jù)朱雪堅等人的研究發(fā)現(xiàn)，基準(zhǔn)面的選擇和掃描方向的選擇對結(jié)果存在較大影響[6]。Ma等人認(rèn)為高程剖面線可以在一定程度上反映DEM的全局特征，因此提出一種4方向剖面線DP簡化方法來綜合DEM[7]。該方法生成剖面線時依賴于VIP特征[8]，若VIP特征保留較少，4個方向的掃描難以確保掃描到全局的DEM特征。

針對以上問題，本文提出一種基于多方向剖面線簡化的DEM綜合方法。首先基于多方向剖面線DP簡化算法提取DEM關(guān)鍵特征點，然后構(gòu)建三角網(wǎng)再內(nèi)插規(guī)則格網(wǎng)DEM，以達到DEM綜合的目的。

1 實驗數(shù)據(jù)

為了驗證本文算法，本文采用如圖1所示的兩組由航空激光點云數(shù)據(jù)內(nèi)插的格網(wǎng)DEM進行實驗，內(nèi)插過程采用surfer 8.0實現(xiàn)，選擇克里金插值算法生成間距為2 m格網(wǎng)DEM。如圖1(a)所示的數(shù)據(jù)1，格網(wǎng)大小為1 416像素×943像素，高程范圍從1 125～1 623 m；如圖1(b)所示的數(shù)據(jù)2，格網(wǎng)大小為913像素×555像素，高程范圍是從1 500～2 115 m，兩組數(shù)據(jù)地形起伏較大。

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

2 基于DP算法的DEM簡化原理及其改進

2.1 DP算法簡介

DP算法，是Douglas和Pucker在1973年提出的一種將二維曲線化簡的算法[3]。DP算法主要有2種形式，其一為“錨點前進法”，另外一種為“分而治之法”，后一種方法效率更高[4]。因此本文采用后種方法，圖2是其原理示意圖。如圖2所示的一條曲線，將這個序列中的首點A和末點J連成一條直線，然后計算序列中間各個點到這條直線的距離di。如果最大距離大于預(yù)先設(shè)置的閾值T(該閾值與簡化后DEM要保留的精度相關(guān)，閾值越小，保留點越多，本文實驗取經(jīng)驗值50)，設(shè)最大距離對應(yīng)的點為S，然后對端點AS和SJ之間的點序列進行相同迭代操作，如果最大距離小于T則停止。最后將插入的端點與A、J連成折線替代原始曲線，從而起到簡化的作用。

圖2 二維DP算法原理圖

2.2 基于多方向剖面線的特征提取

文獻[7]以VIP特征為起點，向東南西北4個方向發(fā)射射線與DEM有效邊界相交形成剖面線，然后進行DP簡化。如果為了更好地覆蓋整個DEM區(qū)域，最直接的方法是將DEM有效區(qū)域的每行的端點看成一條剖面線，然后采用2維DP算法簡化，這類似于傳統(tǒng)3維DP算法按照行列掃描的方式。為了獲取更多方向的剖面線，本文不針對原始DEM直接運算，而是圍繞原始DEM中心位置，按照一定角度間隔對原始DEM進行旋轉(zhuǎn)。對于給定的角度θi，采用式(1)計算每個點旋轉(zhuǎn)以后的坐標(biāo)。

(1)

(2)

其中(X,Y)和(X′,Y′)分別是旋轉(zhuǎn)前后的DEM水平坐標(biāo)，(XC,YC)為旋轉(zhuǎn)前DEM中心位置水平坐標(biāo)。

提取特征后，采用式(2)對特征提取結(jié)果進行變換，以便于將特征提取結(jié)果重新納入到原始DEM所在的坐標(biāo)系，也便于融合所有旋轉(zhuǎn)DEM提取的特征結(jié)果。圖3(a)DEM為圖1(a)所示DEM逆時針旋轉(zhuǎn)15°的結(jié)果。進行DEM旋轉(zhuǎn)后，只需要沿著新DEM水平方向生成剖面線即可實現(xiàn)多方向的剖面線生成。剖面線生成過程中，依次掃描新DEM的每一行，其起始點和終點為DEM有效區(qū)域。圖3(b)為圖3(a) 所示剖面線的高程剖面圖。對于生成的每一條剖面線，采用如圖2所示的2維DP算法進行簡化以獲得特征點，并保留每個特征點對應(yīng)的距離，作為每個特征點的興趣值R，以用于后續(xù)處理。

圖3 DEM旋轉(zhuǎn)結(jié)果與對應(yīng)剖面線圖

當(dāng)旋轉(zhuǎn)間隔較小時，相鄰剖面線獲取的特征點可能重疊，因此本文在完成多方向DP簡化獲得特征序列后，根據(jù)每個特征點記錄的興趣值，在局部進行極值抑制。具體過程為循環(huán)選取每個特征點，然后尋找其水平距離半徑r(本文所有r取經(jīng)驗值3 m)范圍的所有特征點集S，如果當(dāng)前特征點的興趣值R為S中的最小值或者最大值，則保留當(dāng)前點，否則舍棄。經(jīng)過抑制局部極值，局部范圍內(nèi)只保留一個特征明顯的點。由于提取的特征點為離散3維點，不能用原始DEM搜索其鄰域點。為了快速搜索每個點周圍鄰域內(nèi)的點，本文采用K-D樹進行搜索[9]。

綜上所述，本文基于多方向剖面線簡化的特征提取過程可以歸納為圖4所示的流程。

如圖5所示為按照30°角度間隔對圖1(a)所示DEM進行多方向特征提取的結(jié)果。圖5(a) 為多方向剖面線DP簡化后直接疊加的結(jié)果，一共有32 739個點，進行局部抑制極值后，保留6 458個點，結(jié)果如圖5(b)所示。

圖4 多方向剖面線簡化特征提取流程圖

圖5 特征點提取結(jié)果示意

獲得DEM的特征點后，將這些點構(gòu)建三角網(wǎng)，然后再采用線性內(nèi)插的方式生成10 m大小格網(wǎng)的DEM，即為本文方法DEM綜合的結(jié)果。

3 實驗結(jié)果及分析

3.1 評價指標(biāo)

獲得綜合后的規(guī)則格網(wǎng)DEM，分別采用高程中誤差、平均坡度對結(jié)果進行定量評價；采用等高線套合對綜合以后的DEM進行可視化定性分析[10]，指標(biāo)定義如下：

1) 高程中誤差：通過比較簡化以后DEM上每個格網(wǎng)點高程與原始DEM對應(yīng)位置的高程計算而來。

2) 坡度均值差：指簡化以后的DEM坡度平均值，與原始DEM的坡度平均值進行比較的差值。差值越大說明與原始的DEM吻合越低，則DEM簡化過程中，特征保留效果越差。

3) 等高線套合：分別對比簡化后DEM生成等高線與原始等高線套合在原始DEM上的吻合程度，從整體上定性分析綜合前和綜合后的等高線。吻合程度越高，則表明綜合后特征保留越好。

3.2 結(jié)果與分析

利用DP算法進行DEM特征點提取，傳統(tǒng)做法多采用0°和90°方向聯(lián)合提取，本文通過一定間隔角度旋轉(zhuǎn)，可以進行多方向剖面線DP簡化。本文對不同間隔角度進行了實驗，最后選取15°、30°、45°間隔的結(jié)果進行定量和定性評價。

1) 定量評價結(jié)果。 高程中誤差和坡度均值差如表1所示。從表1結(jié)果可以看出高程中誤差隨著旋轉(zhuǎn)角度間隔的增大而增大。結(jié)果表明多方向的DP算法的精度比傳統(tǒng)的算法要高。表1中坡度均值差也表明多方向簡化結(jié)果比傳統(tǒng)兩個方向的結(jié)果精度高。

表1 定量評價結(jié)果

2) 定性評價結(jié)果。 圖6是傳統(tǒng)方法和本文方法重構(gòu)DEM等高線套合結(jié)果；圖7為等高線套合局部放大結(jié)果(藍(lán)色為傳統(tǒng)方法重構(gòu)DEM的等高線，紅色為原始的DEM等高線，綠色、紫色、淡紫色分別為45°、30°、15°多方向的DP算法重構(gòu)DEM等高線)。從所得的等高線套合情況可以看出藍(lán)色等高線偏離紅色等高線最多，其余顏色的等高線與紅色等高線走勢更加吻合。從圖7的局部放大結(jié)果中也可以明顯看出藍(lán)色與紅色的吻合沒有其余顏色與紅色吻合得好。

圖6 等高線套合效果

圖7 等高線套合局部放大結(jié)果

實驗結(jié)果表明本文所提出的多方向DP算法內(nèi)插的DEM比傳統(tǒng)的DP算法內(nèi)插的DEM的精度更好。傳統(tǒng)的DP算法內(nèi)插時只考慮了2個方向，而本文方法通過格網(wǎng)DEM的旋轉(zhuǎn)可以獲取多方向剖面線的整合結(jié)果。故而用多方向算法得到的重構(gòu)DEM比起傳統(tǒng)算法得到的重構(gòu)DEM更接近原始DEM。

4 結(jié) 論

針對DEM綜合的問題，本文提出一種利用多方向剖面線DP算法簡化提取特征點的方法，從實驗結(jié)果可看出，不管是從高程誤差統(tǒng)計還是坡度均值上，都可以看出本文方法得到的DEM與原始DEM地形形態(tài)更加吻合，能更好地保持地形形態(tài)特征。

本文在DEM綜合時，只考慮了特征點，進一步將加入線特征，使得綜合后的DEM能夠更好地保持原始DEM特征。

[1] ZHOU Q, CHEN Y. Generalization of DEM for terrain analysis using a compound method [J]. ISPRS Journal of Photogrammetry and Remote Sensing, 2011,66(1),38-45.

[2] 董有福, 湯國安.利用地形信息強度進行DEM地形簡化研究[J].武漢大學(xué)學(xué)報(信息科學(xué)版), 2013,38(3):353-357.

[3] DOUGLAS D H, PEUCKER T K., Algorithms for the reduction of the number of points required to represent a digitized line or its caricature [J]. The Canadian Cartographer, 1973, 10(2):112-122 .

[4] 費立凡,何津,馬晨燕,等.3維Douglas-Peucker算法及其在DEM自動綜合中的應(yīng)用研究[J].測繪學(xué)報, 2006(3): 278-284.

[5] 何津,費立凡.再論三維Douglas-Peucker算法及其在DEM綜合中的應(yīng)用[J]. 武漢大學(xué)學(xué)報(信息科學(xué)版), 2008 (2):160-163.

[6] 朱雪堅,葉遠(yuǎn)智,湯國安.運用三維Douglas-Peucker算法提取DEM地形特征[J].測繪通報,2014(3):118-121

[7] MA T, CHEN Y, HUA Y, et al. DEM generalization with profile simplification in four directions[J]. Earth Science Informatics, 2017, 10(1): 29-39.

[8] Chen Z T, GUEVARA J A. Systematic selection of very important points (VIP) from digital terrain model for constructing triangular irregular networks[C]. In Auto-carto, 1987, 8, pp. 50-56.

[9] ARYA S, MOUNT D M. Approximate Nearest Neighbor Searching [C], Proc. 4th Annu. ACM-SIAM Sympos. on Discrete Algorithms(SODA’93), 1993, 271-280.

[10] 吳艷蘭,胡海,胡鵬,等.數(shù)字高程模型誤差及其評價的問題綜述[J]. 武漢大學(xué)學(xué)報(信息科學(xué)版), 2011, 36(5): 568-574.