DEM尺度變換中直方圖相似度計算與應用*

牛亮,楊勤科

(西北大學城市與環境學院,西安710127)

研究表明,DEM及其基礎上提取的地形要素因受到水平分辨率的強烈影響,表現出明顯的尺度效應。認識理解尺度效應,是有效利用DEM數據資源,準確提取適應與應用尺度和目的地形數據的基礎[1]。基于多種比例尺地形圖生成的多種分辨率的DEM(原生DEM)為尺度效應研究提供了基礎,但由于不同比例尺地形圖之間并不具備完全相同的定位控制基礎,限制了研究。為此,研究者通過對高分辨率DEM進行變換方法生成多種分辨率(或者柵格)的DEM 數據。常見方法包括重采樣[2-4]、小波變換[5-6]、空間濾波[1]等。在這些方法中,如何將變換結果與原生DEM的分辨率建立聯系,也就是說如何確認變換結果的分辨率,是一個有待研究解決的問題,其中直方圖的相似性比較是一個有效方法。

基于1∶1萬數字地形圖數據插值建立多種柵格尺寸的 DEM,通過直方圖相似度比較找到與1∶5萬、1∶10萬、1∶25萬原生 DEM 等分辨率的DEM,為坡度的尺度效應分析研究和中低分辨率DEM坡度的變換提供依據。

1 數據基礎和研究方法

1.1 研究區域與數據基礎

本研究工作是在黃土高原中部延河流域的二級支流縣南溝流域內進行的。研究區為典型的黃土丘陵,地面起伏較大。基于10m分辨率DEM的坡度圖統計表明,坡度大于25°的地面約占總土地面積的50%。研究所用基礎數據為:(1)1∶1萬、1∶5萬、1∶10萬和1∶25萬四種不同比例尺的數字化地形圖,等高距分別為5 m、20 m、40m和100 m,與之對應的原生DEM,分辨率應分別為5 m、10 m、25 m 和50 m[7-9],標記為 DEM5a 、DEM 10a 、DEM 25a、DEM50a 。

1.2 研究方法

本研究基于1∶1萬數字化地形圖,生成多種柵格尺寸的DEM,然后通過與原生DEM 的對比,使分辨率得到匹配,方法如下。

1.2.1 多種柵格DEM生成 為了弄清地形形態及其與坡面流水的關系,本研究采用 ANUDEM[10-11]軟件對1∶1萬數字化地形圖進行空間插值,生成柵格尺寸為10m 、20 m、25 m、30m、40 m、50 m、60 m和70 m的 DEM,分別標記為DEM10b、DEM 20b、DEM25b、DEM30b、DEM40b、 DEM 50b、DEM 60b和DEM70b。所用參數設置為第二粗糙系數0.5,迭代次數40。

1.2.2 DEM相似性比較參數計算 通常比較兩個DEM之間的相似性和差異性是在比較圖像統計特征基礎上,通過觀察來實現的[1]。本文引入直方圖相似度指數和逼真度指數,以便客觀進行比較。

(1)統計特征:將DEM 看成簡單的矩陣,利用經典統計學統計其特征變量,包括:最小值、最大值、平均值、標準差、眾數等,并計算信息熵[12]。

(2)直方圖相似度:直方圖是數字圖像的重要特征之一,反映圖像(或其中一部分區域)的灰度的總體分布信息。在圖像變換中,變換前后2幅圖像的相似程度通常用直方圖相似性來度量。直方圖相似度指數包括:

①距離直方圖相似度指數:直方圖相似性度量通常是基于2幅直方圖相同窗口之間的距離或經過變換后的一些距離標準,如L1距離(曼哈頓距離)、L2(歐幾里得距離)距離、χ2距離等[13],這些方法側重于衡量兩幅圖像的差異,得到的距離越小表示兩幅圖像越相似。其中L1、L2和 χ2距離的計算公式為:

對于頻率直方圖,L1距離的取值范圍是[0,2],L2距離的取值范圍是,χ2距離的取值范圍是[0,1],值為0時最相似,值為最大值時最不相似。為了統一規范化相似性距離度量函數的取值范圍,將兩幅圖像的相似度定義在[0,1]之間,值越接近1則圖像越相似。對于累計頻率直方圖而言,統一規范化相似性距離度量函數的取值范圍比較困難,直接采用距離來判定相似性,距離越小越相似。

②相交直方圖相似度指數。另一種方法稱為直方圖的交[14],對于歸一化的直方圖,其具體計算公式如下:

直方圖的交運算就是取兩幅圖像的特征直方圖在各個數值級上的較小值,累加即表示圖像之間的相似程度。這種相似性實際上表示兩幅圖像在該直方圖特征意義下的公共部分。0≤D(X,Y)≤1,它的值越大,則兩幅圖像越相似,如果兩幅圖像完全一致 ,則 D(X,Y)=1。

(3)圖像逼真度指數。圖像逼真度,是描述待對比圖像是否存在較大誤差的一個重要技術指標[15],本文根據公式(5)計算,以逼真度接近數值1的程度作為處理前后圖像的近似性度量[15],越接近于1,逼真度越高,反之越小。

式中:I——原生DEM 的坡度和曲率;I′——多柵格DEM 的坡度和曲率;Iij——原生DEM 的坡度和曲率的第i行j列柵格單元的值;I′ij——多柵格DEM的坡度和曲率第i行j列柵格單元值。

1.3 試驗步驟

本研究的試驗只要包括DEM 生成、統計表、直方圖相似度和逼真度計算、分析比較等。工作流程如圖1示。

(1)DEM生成。將經過質量檢查編輯好的1∶1萬、1∶5萬、1∶10萬、1∶25萬數字地形圖分別輸入到ANUDEM軟件中,所采用的主要參數是第二粗糙系數0.5,迭代次數40。1∶1萬選擇的柵格尺寸分別為 10m 、20m 、25m 、30 m 、40m 、50m 、60m和 70 m,生成 標記為 DEM10b、DEM20b、DEM25b、DEM 30b、DEM 40b、DEM 50b、DEM 60b和 DEM 70b的多柵格DEM 。1∶5萬、1∶10萬和1∶25萬地形圖根據基于信息含量確定相應分辨率分別為10 m、25 m和50 m,生成對應分辨率 DEM,分別標記為:DEM 10a 、DEM 25a 、DEM 50a 。

(2)坡度和曲率的生成及統計參數比較。在ARCGIS軟件中分別對所有的DEM求取坡度和曲率,并記錄各坡度和曲率的空間統計特征參數,包括最小值、最大值、平均值、標準差、眾數和信息熵。通過比較各統計特征參數,找到統計特征參數最為相似的兩個坡度或者曲率。

(3)坡度和曲率直方圖相似度計算。生成坡度和曲率做頻率和累計頻率圖。頻率和累計頻率圖實質上也是直方圖的表現形式。通過直方圖相似度L1距離、L2距離、χ2距離和直方圖交四種方法分別計算基本比例尺1∶1萬地形圖生成DEM的坡度、曲率的頻率和累計頻率與 1∶5萬、1∶10萬和1∶25萬地形圖生成DEM的坡度、曲率的頻率和累計頻率的相似度,從而分別找到與1∶5萬、1∶10萬、1∶25萬比例尺DEM 的最為相似的1∶1萬比例尺插值DEM的柵格尺寸。

(4)坡度和曲率的直方圖逼真度計算。在A rcGIS軟件中,先使用重采樣工具對多柵格DEM重采樣生成目標分辨率DEM的柵格大小,使用柵格計算器進行柵格計算得出相應參數,帶入到公式(8)中得到逼真度值。

(5)分辨率匹配。根據統計參數、直方圖相似度、逼真度分別將多柵格DEM與多分辨率DEM進行匹配,根據結果得出相關結論。

圖1 研究技術路線圖

2 結果與分析

2.1 坡度和曲率的統計特征參數

坡度和曲率統計特征參數見表1和表2。從表1可以看出 DEM10b的坡度統計特征與目標數據DEM 10a上提取的坡度特征更為接近,DEM 25b的坡度統計特征與目標數據DEM 25a上提取的坡度特征更為接近。從表2中可以看到DEM 10b的曲率統計特征與目標數據DEM10a上提取的曲率特征更為接近,DEM25b的曲率統計特征與目標數據DEM25a上提取的曲率特征更為接近。綜合坡度和曲率的統計特征,可以認為DEM 10b、DEM25b的地形起伏反映能力分別相當于DEM10a、DEM25a。從表1中找不到與DEM 50a的坡度統計特征相近的由1∶1萬地形圖生成的DEM,從表2中可以看到DEM60b與目標數據DEM 50a的曲率統計特征較為相近。

表1 坡度的統計特征

表2 曲率的統計特征

2.2 坡度和曲率的相似度計算

坡度和曲率的頻率直方圖和累計頻率直方圖如下:以與DEM 10a的相似度比較為例(見圖2)。

對坡度和曲率的頻率和累計頻率分別使用直方圖L1距離、L2距離、χ2距離和直方圖交4種方法求取相關直方圖的相似度,與 DEM10a、DEM25a和DEM50a比較的相似度計算結果如表3-5所示(加黑數字為最大相似度值):

圖2 部分DEM坡度和曲率頻率分布

表3 與DEM 10a比較坡度、曲率直方圖相似度

表4 與DEM 25a比較坡度、曲率直方圖相似度

表5 與DEM 50a比較坡度、曲率直方圖相似度

從表3中可以看出DEM 10b與DEM 10a的坡度頻率、坡度累計頻率、曲率頻率和曲率累計頻率的所有相似度值都最高,因此可以認為DEM10b與DEM 10a最為相似。從表4中可以看出 DEM25b與DEM25a的坡度頻率和坡度累計頻率的各相似度值最高,DEM 20b和DEM25b與DEM 25a的曲率頻率和曲率累計頻率的相似度值接近并最為相似,綜合認為DEM25b與DEM 25a最為相似。從表5中可以看出隨著DEM柵格尺寸不斷變大,與DEM 50a的坡度頻率和累計頻率的相似度變大,是一個增函數,但是相似度是不是一直增大,或者有個峰值,還需生成更大柵格尺寸的DEM 進行驗證,再生成 DEM80b,DEM90b和 DEM100b進行驗證,見表6所示,結果顯示,DEM 80b的頻率各相似度值最高,而DEM 90b的累計頻率相似度值最高。同時也可以從表5看出,DEM60b與DEM50a的曲率頻率、曲率累計頻率的相似度最高,綜合來看,坡度和曲率的相似度不統一。

2.3 坡度和曲率的逼真度計算

采用公式(8),計算出的坡度和曲率的逼真度如表7-8所示:

表6 與DEM 50a比較坡度、曲率直方圖相似度

表7 坡度的逼真度值

表8 曲率的逼真度值

從表7中可以看出DEM 10b與DEM 10a,DEM 20b與DEM 25a,DEM 60b與DEM50a的逼真度值最高,都接近于1。從表8中得出的結果是DEM 10b與DEM 10a,DEM 20b與DEM25a,DEM 40b與 DEM 50a的逼真度值較高,但是曲率逼真度值總體較坡度的逼真度值小。

3 結論與討論

3.1 結論

在變換前后DEM統計特征比較的基礎上,通過引入直方圖相似度和逼真度指數兩個指標,可以對分辨率匹配提供依據。

(1)初步分析表明,基于1∶1萬地形圖插值生成的 10 m、20 m或 25 m柵格的 DEM,與基于1∶5萬、1∶10萬地形圖生成的原生DEM(分辨率為10 m、25 m)相似。從統計特征參數、直方圖相似度指數和逼真度指數這三種方法得不到統一的與1∶25萬原生分辨率DEM相似的基于1∶1萬地形圖插值柵格大小的DEM,并且各方法得到的相似度值并不高。也可以看出隨著比例尺變小,原生DEM與1∶1萬地形圖插值生成的柵格DEM的相似度變小。

(2)關于相似度指數,從本文的試驗結果看,四種直方圖相似度的計算方法對計算DEM統計特征的直方圖相似度同樣適用,取得基本相同的結果,但是L1距離、L2距離、χ2距離相比直方圖的交方法的優勢在于可以對累計直方圖進行計算。L1距離所反映的相似度值的幅度變化較大,所以更加直觀,因此L1距離是4種方法中最適合的相似度指數計算方法。

3.2 討論

(1)由于直方圖相似度和逼真度計算只是基于統計學原理,而DEM具有空間格局屬性,因此在直方圖相似度和逼真度計算基礎上,還必須充分考慮DEM的格局和地形信息含量等,才能形成全面準確的認識;

(2)隨著比例尺的變小,原生DEM 與1∶1萬地形圖插值生成的柵格DEM的相似度變小。可能的原因是,DEM地形分析的尺度推繹要考慮DEM格網單元異質性和尺度推繹范圍,由于尺度推繹非線性的根本原因,尺度推繹不能無限度的進行[16],隨著從高分辨率到低分辨率的轉換過程,尺度轉換的有效性在降低。

(3)利用1∶1萬地形圖插值,生成與1∶25萬地形圖原生DEM分辨率對應的DEM的有效性還需進一步探討。

[1] 楊勤科,Jupp D B L,郭偉玲,等.基于濾波方法的DEM尺度變換方法研究[J].水土保持通報,2008,28(6):58-62.

[2] 郝振純,池宸星,王玲,等,DEM空間分辨率的初步分析[J].地地球科學進展,2005,20(5):499-504.

[3] Wo lck D M,M cCabe G J.Differences in topographic characteristics computed from 100-m and 1000-m resolution digital elevation model data[J].hyd ro logical p rocesses,2000,14(2):987-1002.

[4] 吳險峰,劉昌明,王中根.柵格DEM的水平分辨率對流域特征的影響分析[J].自然資源學報,2003,18(2):168-173.

[5] 吳凡,祝國瑞.基于小波分析的地貌多尺度表達與自動綜合[J].武漢大學學報:信息科學版,2001,26(2):170-176.

[6] 于浩,楊勤科,張曉萍,等.基于小波多尺度分析的DEM 數據綜合研究[J].測繪科學,2008,33(3):93-95.

[7] 張彩霞,楊勤科,段建軍.高分辨率數字高程模型的構建方法[J].水利學報,2006,37(8):1009-1013.

[8] 師維娟,楊勤科,趙東波,等.中分辨率水文地貌關系正確DEM 建立方法研究:以黃土丘陵區為例[J].西北農林科技大學學報:自然科學版,2007,35(2):143-148.

[9] Yang QK,McVicar T R,Van Niel TG,etal.Imp roving terrain rep resentation of a digital elevationmodel by reducing source data errors and optim ising interpolation algorithm parameters:an example in the Loess Plateau[J].China International Journal of App lied Earth Observation and Geoinformation,2007,9(3):235-246.

[10] Hutchinson M F.ANUDEM Version 5.1.A.N.U.Centre for Resource and Environmental Studies[M].Canberra,2004.

[11] 楊勤科,M cvicar T R,李領濤,等.ANUDEM-專業化數字高程模型插值算法及其特點[J].干旱地區農業研究[J].2006,24(3):36-41.

[12] 朱偉,王東華,周曉光.基于信息熵的DEM最佳分辨率確定方法研究[J].遙感信息,2008(5):79-82.

[13] 馬煜,顧曉東,汪源源.基于平均窗口平移的直方圖相似性度量[J].應用科學學報,2008,26(1):28-33.

[14] 林麗惠.一種改進的顏色直方圖相似性度量算法[J].武夷學院學報,2009,28(2):58-61.

[15] 尹忠海,簡劍鋒,周利華,等.現有影像逼真度定義的缺陷及其改進[J].西安電子科技大學學報:自然科學版[J].2004,31(6):833-836,869.

[16] 劉學軍,盧華興,仁政,等.論DEM 地形分析中的尺度問題[J].地理研究,2007,26(3):433-442.