GIS的核心之一：數(shù)字地形分析

1、第9章 DEM與數(shù)字地形分析數(shù)字地面模型于1958年提出，特別是基于DEM的GIS空間分析方法的出現(xiàn)，使傳統(tǒng)的地形分析方法產(chǎn)生了革命性的變化，數(shù)字地形分析方法逐步形成和完善。目前，基于DEM的數(shù)字地形分析已經(jīng)成為GIS空間分析中最具特色的部分，在測繪、遙感及資源調(diào)查、環(huán)境保護(hù)、城市規(guī)劃、災(zāi)害防治及地學(xué)研究各方面發(fā)揮越來越重要的作用。本章首先介紹了數(shù)字高程模型的基本概念和建立步驟，然后從基本坡面因子、特征地形因子、水文因子和可視域等方面簡述數(shù)字地形分析的主要內(nèi)容和研究方法。9.1 基本概念9.1.1 數(shù)字高程模型數(shù)字高程模型（Digital Elevation Model，簡稱DEM）是通過有限

2、的地形高程數(shù)據(jù)實現(xiàn)對地形曲面的數(shù)字化模擬（即地形表面形態(tài)的數(shù)字化表示），它是對二維地理空間上具有連續(xù)變化特征地理現(xiàn)象的模型化表達(dá)和過程模擬。由于高程數(shù)據(jù)常常采用絕對高程（即從大地水準(zhǔn)面起算的高度），DEM也常常稱為DTM（Digital Terrain Model）?！癟errain”一詞的含義比較廣泛，不同專業(yè)背景對“Terrain”的理解也不一樣，因此DTM趨向于表達(dá)比DEM更為廣泛的內(nèi)容。從研究對象與應(yīng)用范疇角度出發(fā)，DEM可以歸納為狹義和廣義兩種定義。從狹義角度定義，DEM是區(qū)域表面海拔高程的數(shù)字化表達(dá)。這種定義將描述的范疇集中地限制在“地表”、“海拔高程”及“數(shù)字化表達(dá)”內(nèi)，觀念較為

3、明確。從廣義角度定義，DEM是地理空間中地理對象表面海拔高度的數(shù)字化表達(dá)。這是隨著DEM的應(yīng)用不斷向海底、地下巖層以及某些不可見的地理現(xiàn)象（如空中的等氣壓面等）延伸，而提出的更廣義的概念。該定義將描述對象不再限定在“地表面”，因而具有更大的包容性，有海底DEM、下伏巖層DEM、大氣等壓面DEM等。數(shù)學(xué)意義上的數(shù)字高程模型是定義在二維空間上的連續(xù)函數(shù)。由于連續(xù)函數(shù)的無限性，DEM通常是將有限的采樣點用某種規(guī)則連接成一系列的曲面或平面片來逼近原始曲面，因此DEM的數(shù)學(xué)定義為區(qū)域D的采樣點或內(nèi)插點Pj按某種規(guī)則連接成的面片M的集合： （9.1）DEM按照其結(jié)構(gòu)，可分為規(guī)則格網(wǎng)DEM、TIN、基于點的

4、DEM和基于等高線的DEM等。由于規(guī)則格網(wǎng)結(jié)構(gòu)簡單，算法設(shè)計明了，在實際運(yùn)用中被廣泛采用。本書中的DEM僅指規(guī)則格網(wǎng)DEM。9.1.2 數(shù)字地形分析數(shù)字地形分析（Digital Terrain Analysis, DTA），是指在數(shù)字高程模型上進(jìn)行地形屬性計算和特征提取的數(shù)字信息處理技術(shù)。DTA技術(shù)是各種與地形因素相關(guān)空間模擬技術(shù)的基礎(chǔ)。 地形屬性根據(jù)地形要素的關(guān)系特征和計算特征，可以歸納為地形曲面參數(shù)（parameters）、地形形態(tài)特征（features）、地形統(tǒng)計特征（statistics）和復(fù)合地形屬性（compound attributes）。地形曲面參數(shù)具有明確的數(shù)學(xué)表達(dá)式和物理定

5、義，并可在DEM上直接量算，如坡度、坡向、曲率等。地形形態(tài)特征是地表形態(tài)和特征的定性表達(dá)，可以在DEM上直接提取，其特點是定義明確，但邊界條件有一定的模糊性，難以用數(shù)學(xué)表達(dá)式表達(dá)，如在實際的流域單元的劃分中，往往難于確定流域的邊界。地形統(tǒng)計特征是指給定地表區(qū)域的統(tǒng)計學(xué)上的特征。復(fù)合地形屬性是在地形曲面參數(shù)和地形形態(tài)特征的基礎(chǔ)上，利用應(yīng)用學(xué)科（如水文學(xué)、地貌學(xué)和土壤學(xué)）的應(yīng)用模型而建立的環(huán)境變量，通常以指數(shù)形式表達(dá)。數(shù)字地形分析的主要內(nèi)容有兩方面，一是在復(fù)雜的現(xiàn)實世界地理過程中各影響因子和簡單、高效、精確、易于理解的抽象與計算機(jī)實現(xiàn)中找到平衡。簡單地說，就是提取描述地形屬性和特征的因子，并利用各

6、種相關(guān)技術(shù)分析解釋地貌形態(tài)、劃分地貌形態(tài)等。二是DTM的可視化分析。數(shù)字地形分析中可視化分析的重點在于地形特征的可視化表達(dá)和信息增強(qiáng)，以幫助傳達(dá)地形曲面參數(shù)、地表形態(tài)特征和復(fù)合地形屬性的信息。根據(jù)分析內(nèi)容，常用的數(shù)字地形分析的方法有以下幾種（圖9.1）：1.提取坡面地形因子地形定量因子是為有效地研究與表達(dá)地貌形態(tài)特征所設(shè)定的具有一定意義的參數(shù)或指標(biāo)。從地形地貌的角度考慮，地表是由不同的坡面組成的，而地貌的變化，完全源于坡面的變化。常用的坡面地形因子有坡度、坡向、平面曲率、坡面曲率、地形起伏度、粗糙度、切割深度等。2.提取特征地形要素（1）流域分析流域分析主要是根據(jù)地表物質(zhì)運(yùn)動的特性，特別是水流

7、運(yùn)動的特點，利用水流模擬的方法來提取水系、山脊線、谷底線等地形特征線，并通過線狀信息分析其面域特征。（2）可視域分析可視性分析包括兩方面內(nèi)容，一個是兩點之間的通視性(Intervisibility)，另一個是可視域(ViewShed)，即對于給定的觀察點所覆蓋的區(qū)域。3.地形統(tǒng)計特征分析地形統(tǒng)計分析是應(yīng)用統(tǒng)計方法對描述地形特征的各種可量化的因子或參數(shù)進(jìn)行相關(guān)、回歸、趨勢面、聚類等統(tǒng)計分析，找出各因子或參數(shù)的變化規(guī)律和內(nèi)在聯(lián)系，并選擇合適的因子或參數(shù)建立地學(xué)模型，從更深層次探討地形演化及其空間變異規(guī)律。數(shù)字地形分析方法坡面地形因子提取特征地形要素提取地形統(tǒng)計特征分析流域分析可視域分析提取坡度提取

8、坡面曲率圖9.1 數(shù)字地形分析常用方法9.2 DEM建立9.2.1 DEM建立的一般步驟數(shù)字高程模型的建立過程是一個模型建立過程。從模型論角度講，就是將源域（地形）表現(xiàn)在另一個域（目標(biāo)域或DEM）中的一種結(jié)構(gòu)，建模的目的是對復(fù)雜的客體進(jìn)行簡化和抽象，并把對客體（源域，DEM中為地形起伏）的研究轉(zhuǎn)移到對模型的研究上來。模型建立之初，首先要為模型構(gòu)造一個合適的空間結(jié)構(gòu)（spatial framework）?？臻g結(jié)構(gòu)是為把特定區(qū)域內(nèi)的空間目標(biāo)鑲嵌在一起而對區(qū)域進(jìn)行的劃分，劃分出的各個空間范圍稱為位置區(qū)域或空間域?？臻g結(jié)構(gòu)一般是規(guī)則的（如格網(wǎng)），或不規(guī)則的（如不規(guī)則三角網(wǎng)TIN）。建立在空間結(jié)構(gòu)基礎(chǔ)上

9、的模型是由n個空間域的有限集合組成。由于空間數(shù)據(jù)包含位置特征和屬性特征，而屬性特征是定義在位置特征上的，因此每一個空間域就是由空間結(jié)構(gòu)到屬性域的計算函數(shù)或域函數(shù)。模型的可計算性要求有兩點，一是空間域的數(shù)量、屬性域和空間結(jié)構(gòu)是有限的，二是域函數(shù)是可計算的。構(gòu)筑模型的一般內(nèi)容和過程為：采用合適的空間模型構(gòu)造空間結(jié)構(gòu)；采用合適的屬性域函數(shù)；在空間結(jié)構(gòu)中進(jìn)行采樣，構(gòu)造空間域函數(shù)；利用空間域函數(shù)進(jìn)行分析。當(dāng)空間結(jié)構(gòu)為歐幾里德平面，屬性域是實數(shù)集合時，模型為一自然表面。將歐幾里德平面充當(dāng)水平的XY平面，屬性域給出Z坐標(biāo)（或高程），模型即為數(shù)字高程模型。對于數(shù)字高程模型而言，空間結(jié)構(gòu)的構(gòu)造過程即為DEM的格

10、網(wǎng)化過程（形成格網(wǎng)），屬性值為高程，構(gòu)造空間域函數(shù)即為內(nèi)插函數(shù)的確定，利用空間域函數(shù)進(jìn)行分析就是求取格網(wǎng)點的函數(shù)值。9.2.2 規(guī)則格網(wǎng)DEM的建立DEM是在二維空間上對三維地形表面的描述。構(gòu)建DEM的整體思路是首先在二維平面上對研究區(qū)域進(jìn)行格網(wǎng)劃分（格網(wǎng)大小取決于DEM的應(yīng)用目的），形成覆蓋整個區(qū)域的格網(wǎng)空間結(jié)構(gòu)，然后利用分布在格網(wǎng)點周圍的地形采樣點內(nèi)插計算格網(wǎng)點的高程值，最后按一定的格式輸出，形成該地區(qū)的格網(wǎng)DEM（圖9.1）。x不規(guī)則分布點規(guī)則分布等高線分布 Y對每一格網(wǎng)點求取格網(wǎng)點高程圖9.2 格網(wǎng)DEM建立流程9.2.3 DEM內(nèi)插方法DEM建立過程中的關(guān)鍵環(huán)節(jié)是根據(jù)采樣點的值內(nèi)插計

11、算格網(wǎng)點上的高程值。內(nèi)插是指根據(jù)分布在內(nèi)插點周圍的已知參考點的高程值求出未知點的高程值，它是DEM的核心問題，貫穿于DEM的生產(chǎn)、質(zhì)量控制、精度評定、分析應(yīng)用的各個環(huán)節(jié)。隨著DEM的發(fā)展和完善，已經(jīng)提出了多種高程內(nèi)插方法。根據(jù)不同的分類標(biāo)準(zhǔn)，有不同的內(nèi)插方法分類，例如按數(shù)據(jù)分布規(guī)律分類，有基于規(guī)則分布數(shù)據(jù)的內(nèi)插方法、基于不規(guī)則分布的內(nèi)插方法和適合于等高線數(shù)據(jù)的內(nèi)插方法等；按內(nèi)插點的分布范圍，內(nèi)插方法分為整體內(nèi)插、局部內(nèi)插和逐點內(nèi)插法；從內(nèi)插函數(shù)與參考點的關(guān)系方面，又分為曲面通過所有采樣點的純二維插值方法和曲面不通過參考點的曲面擬合插值方法；從內(nèi)插曲面的數(shù)學(xué)性質(zhì)來講，有多項式內(nèi)插、樣條內(nèi)插、最小

12、二乘配置內(nèi)插等內(nèi)插函數(shù)；從對地形曲面理解的角度，內(nèi)插方法有克立金法、 多層曲面疊加法、加權(quán)平均法、分形內(nèi)插等。表9.1對各種DEM內(nèi)插分類方法進(jìn)行了簡要的總結(jié)和歸納。本小節(jié)僅從內(nèi)插點的分布范圍來看，簡要介紹整體內(nèi)插法、局部內(nèi)插法和逐點內(nèi)插法。詳細(xì)介紹參見第十章。表9.1 DEM內(nèi)插分類方法DEM內(nèi)插數(shù)據(jù)分布規(guī)則分布內(nèi)插方法不規(guī)則分布內(nèi)插方法等高線數(shù)據(jù)內(nèi)插方法內(nèi)插范圍整體內(nèi)插方法局部內(nèi)插方法逐點內(nèi)插方法內(nèi)插曲面與參考點關(guān)系純二維內(nèi)插曲面擬合內(nèi)插內(nèi)插函數(shù)性質(zhì)多項式內(nèi)插線性插值雙線性插值高次多項式插值樣條內(nèi)插有限元內(nèi)插最小二乘配置內(nèi)插地形特征理解克立金內(nèi)插多層曲面疊加內(nèi)插加權(quán)平均值內(nèi)插分形內(nèi)插傅立葉

13、級數(shù)內(nèi)插整體內(nèi)插是指在整個區(qū)域用一個數(shù)學(xué)函數(shù)來表達(dá)地形曲面。整體內(nèi)插函數(shù)通常是高次多項式，要求地形采樣點的個數(shù)大于或等于多項式的系數(shù)數(shù)目。整體內(nèi)插方法有整個區(qū)域上函數(shù)的唯一性、能得到全局光滑連續(xù)的DEM、充分反映宏觀地形特征等優(yōu)點。但由于整體內(nèi)插函數(shù)往往是高次多項式，它也有保凸性較差、不容易得到穩(wěn)定的數(shù)值解、多項式系數(shù)的物理意義不明顯、解算速度慢且對計算機(jī)容量要求較高、不能提供內(nèi)插區(qū)域的局部地形特征等缺點。在DEM內(nèi)插中，一般是與局部內(nèi)插方法配合使用，例如在使用局部內(nèi)插方法前，利用整體內(nèi)插去掉不符合總體趨勢的宏觀地物特征。另外也可用來進(jìn)行地形采樣數(shù)據(jù)中的粗差檢測。局部分塊內(nèi)插是將地形區(qū)域按一定

14、的方法進(jìn)行分塊，對每一分塊，根據(jù)其地形曲面特征單獨進(jìn)行曲面擬合和高程內(nèi)插。一般按地形結(jié)構(gòu)線或規(guī)則區(qū)域進(jìn)行分塊，分塊的大小取決于地形的復(fù)雜程度、地形采樣點的密度和分布。為保證相鄰分塊之間的曲面平滑連接，相鄰分塊之間要有一定寬度的重疊，或者對內(nèi)插曲面補(bǔ)充一定的連續(xù)性條件。這種方法簡化了地形的曲面形態(tài)，使得每一分塊可用不同的曲面表達(dá)，同時得到光滑連續(xù)的空間曲面。不同的分塊單元可以使用不同的內(nèi)插函數(shù)。常用的內(nèi)插函數(shù)有線性內(nèi)插、雙線性內(nèi)插、多項式內(nèi)插、樣條函數(shù)、多層曲面疊加法等。逐點內(nèi)插是以內(nèi)插點為中心，確定一個鄰域范圍，用落在鄰域范圍內(nèi)的采樣點計算內(nèi)插點的高程值。逐點內(nèi)插本質(zhì)上是局部內(nèi)插，但與局部分塊

15、內(nèi)插不同的是，局部內(nèi)插中的分塊范圍一經(jīng)確定，在整個內(nèi)插過程中其大小、形狀和位置是不變的，凡是落在該塊中的內(nèi)插點，都用該塊中的內(nèi)插函數(shù)進(jìn)行計算，而逐點內(nèi)插法的鄰域范圍大小、形狀、位置乃至采樣點個數(shù)隨內(nèi)插點的位置而變動，一套數(shù)據(jù)只用來進(jìn)行一個內(nèi)插點的計算。逐點內(nèi)插法要注意兩個問題，一是選擇合適的內(nèi)插函數(shù)，內(nèi)插函數(shù)決定著DEM精度、DEM連續(xù)性、內(nèi)插點鄰域的最小采樣點個數(shù)和內(nèi)插計算效率。二是確定內(nèi)插點鄰域，內(nèi)插點的鄰域大小和形狀、鄰域內(nèi)參加內(nèi)插計算的數(shù)據(jù)點的個數(shù)、采樣點的權(quán)重、采樣點的分布、附加信息等不僅會影響到DEM的內(nèi)插精度，也影響到內(nèi)插速度。逐點內(nèi)插方法計算簡單，內(nèi)插效率較高，應(yīng)用比較靈活，是

16、目前較為常用的一類DEM內(nèi)插方法。在建立DEM時，要根據(jù)情況選擇合適的、運(yùn)算效率高的方法。而眾多內(nèi)插方法并不是獨立的，而往往是相互結(jié)合使用，這在后續(xù)的章節(jié)里會講到。9.3 數(shù)字地形分析地形分析是地形環(huán)境認(rèn)知的一種重要手段，傳統(tǒng)的地形分析是基于二維平面地圖進(jìn)行的。從基于紙質(zhì)地圖的地形分析發(fā)展到到基于數(shù)字地圖的地形分析，計算機(jī)取代了大量的人工計算和繪制，地形分析的手段、功能發(fā)生了一次飛躍；可視化技術(shù)和虛擬現(xiàn)實技術(shù)的發(fā)展，使得建立三維實時、交互的仿真地形環(huán)境成為可能，同時也需要實現(xiàn)三維地形環(huán)境中的地形分析。特別是DEM的出現(xiàn)和大量應(yīng)用，使得從地形屬性中提取各類地形參數(shù)和特征因子更加的簡便和準(zhǔn)確。用來

17、描述地形特征和空間分布的地形參數(shù)很多，不同的應(yīng)用目的，不同的學(xué)科和領(lǐng)域?qū)Υ说睦斫夂头诸愐膊煌?。本章將綜合相關(guān)知識，著重介紹基本因子分析、地形特征提取、水文分析和可視域分析。9.3.1 基本因子分析本質(zhì)上講，DEM是地形的一個數(shù)學(xué)模型，可以看成是一個或多個函數(shù)的集合。實際上許多地形因子就是從這些函數(shù)進(jìn)行一階或二階推導(dǎo)出來的，也有的通過某種組合或復(fù)合運(yùn)算得到?；镜匦我蜃影ㄐ逼乱蜃樱ㄆ露取⑵孪?、坡度變化率、坡向變化率等）、面積因子（表面積、投影面積、剖面積）、體積因子（山體體積、挖填體積）和面元因子（相對高差、粗糙度、凹凸系數(shù)、高程變異等）。本節(jié)將闡述一些常用的基本地形因子，為了方便起見，并從實

18、際應(yīng)用角度考慮，本節(jié)這些地形因子的計算都是基于格網(wǎng)DEM。1. 坡度x (N)y (E)zPSlopeAspectn圖9.5 地表單元坡度示意圖o嚴(yán)格地講，地表面任一點的坡度是指過該點的切平面與水平地面的夾角。坡度表示了地表面在該點的傾斜程度，在數(shù)值上等于過該點的地表微分單元的法矢量與z軸的夾角（如圖9.5所示），即：Slope = （9.1） 當(dāng)具體進(jìn)行坡度提取時，常采用簡化的差分公式，完整的數(shù)學(xué)表示為： （9.2）式中，fx是X方向高程變化率，fy是Y方向高程變化率。地面坡度實質(zhì)是一個微分的概念，地面上每一點都有坡度，它是一個微分點上的概念，是地表曲面函數(shù)z = f(x,y)在東西、南北方

19、向上的高程變化率的函數(shù)。實際應(yīng)用中，坡度有兩種表示方式（如圖9.6）：l 坡度(degree of slope)：即水平面與地形面之間夾角。l 坡度百分比（percent slope）：即高程增量（rise）與水平增量（run）之比的百分?jǐn)?shù)。圖9.6 坡度的兩種表示方法擬合曲面法是解求坡度的最常用的方法。擬合曲面法，一般采用二次曲面，即在3×3的DEM柵格分析窗口中（如圖9.7）進(jìn)行，每個柵格中心為一個高程值，分析窗口在DEM數(shù)據(jù)矩陣中連續(xù)移動完成整個區(qū)域的計算工作。常用的計算fx、 fy的方法是三階反距離平方權(quán)，該算法也用于ArcView和ARC/INFO。其計算方法為： （9.3

20、）式（9.3），g為格網(wǎng)間距。2. 坡向坡向定義為：地表面上一點的切平面的法線矢量在水平面的投影與過該點的正北方向的夾角（如表7.1中的坡向示意圖所示，x軸為正北方向）。其數(shù)學(xué)表達(dá)公式為： （9.4）對于地面任何一點來說，坡向表征了該點高程值改變量的最大變化方向。在輸出的坡向數(shù)據(jù)中，坡向值有如下規(guī)定：正北方向為0°，順時針方向計算，取值范圍為0°360°。坡向可在DEM數(shù)據(jù)中用式9.4直接提取。但應(yīng)注意，由于式9.4求出坡向有與x軸正向和x軸負(fù)向夾角之分，此時就要根據(jù)fx和fy的符號來進(jìn)一步確定坡向值（如表9.2所示）圖9.7 3×3分析窗口表9.2 坡

21、向值的判斷fxfy = Aspect坡向示意=0>090y (E)360+180+180+x (N)=0-1<0270>0>0090=000<0-900360+<0>0-900180+=00180<0090180+注：上述情況假定所建立的DEM數(shù)據(jù)從南向北獲取的，且x軸與正北方向重合，否則上述公式求得的坡向值，還應(yīng)加上x軸偏離正北方向的夾角值。采用這種方法求取的坡向分級比較詳細(xì)，但實際應(yīng)用中往往需要給予歸并，在ArcView和ArcGIS軟件中，通常把坡向綜合成九種坡向：平緩坡（1）、北坡（0° - 22.5°, 337.5&

22、#176; - 360°）、東北坡（22.5° - 67.5°）、東坡（67.5° - 112.5°）、東南坡（112.5° - 157.5°）、南坡（157.5° - 202.5°）、西南坡（202.5° - 247.5°）、西坡（247.5° - 292.5°）、西北坡（292.5° - 337.5°）。圖9.8 原始DEM數(shù)據(jù)及實驗區(qū)等高線圖圖9.9 ARCVIEW軟件下提取的坡度圖圖9.10 由DEM提取的坡向圖3. 曲率曲率是對地形表面

23、一點扭曲變化程度的定量化度量因子，地面曲率在垂直和水平兩個方向上分量分別稱為平面曲率和剖面曲率。地形表面曲率反映了地形結(jié)構(gòu)和形態(tài)，同時也影響著土壤有機(jī)物含量的分布，在地表過程模擬、水文、土壤等領(lǐng)域有著重要的應(yīng)用價值和意義。剖面曲率是對地面坡度的沿最大坡降方向地面高程變化率的度量。數(shù)學(xué)表達(dá)式為： （9.5）圖9.11平面曲率示意圖x（N）y（E）zPV水平方向切地表面所得曲線平面曲率指在地形表面上，具體到任何一點P，指用過該點的水平面沿水平方向切地形表面所得的曲線在該點的曲率值（圖9.11所示）。平面曲率描述的是地表曲面沿水平方向的彎曲、變化情況，也就是該點所在的地面等高線的彎曲程度。從另一個角

24、度講，地形表面上一點的平面曲率也是對該點微小范圍內(nèi)坡向變化程度的度量。數(shù)學(xué)表達(dá)式為： （9.6）曲率數(shù)學(xué)表達(dá)式中，利用離散的DEM數(shù)據(jù)把地表曲面數(shù)學(xué)模擬為一個連續(xù)的曲面H(x,y)，x和y地面點的平面坐標(biāo)值，H(x,y)為地面點高程值，式中其它符號所表示的意義為：，是x方向高程變化率；，是y方向高程變化率；，對高程值在x方向上的變化率進(jìn)行同方向求算變化率，即x方向高程變化率的變化率；，對高程值在x方向上的變化率進(jìn)行y方向上求算變化率，即x方向高程變化率在y方向的變化率； ，對高程值在y方向上的變化率同方向上求算變化率，即y方向高程變化率的變化率。曲率因子的提取算法的基本原理為：在DEM數(shù)據(jù)的基

25、礎(chǔ)上，根據(jù)其離散的高程數(shù)值，把地表模擬成一個連續(xù)的曲面，從微分幾何的思想出發(fā)，模擬曲面上每一點所處的垂直于和平行于水平面的曲線，利用曲線曲率的求算方法的推導(dǎo)得出各個曲率因子的計算公式。利用公式求算出每一點的曲率值的關(guān)鍵在于確定得出式中各個參量的值，在DEM中求算高程的微分分量有一套獨特的算法，最常用是三階反距離平方權(quán)差分。對每一個柵格點都確定一個3×3的分析窗口，其過程如圖9.12所示。利用ArcView所提取的剖面曲率與平面曲率圖如圖9.13和9.14所示。圖9.13 ArcView提取的剖面曲率圖9.14 ArcView提取的平面曲率q值矩陣p值矩陣圖9.12 地面曲率提取步驟流

26、程圖abcdefghI高程值陣列a'xb'xc'xd'xe'xf'xg'xh'xI'xa'yb'yc'yd'ye'yf'yg'yh'yI'yKV剖面曲率Kh平面曲率4. 宏觀地形因子地形起伏度、地形表面粗糙度與地表切割深度等地形因子是描述和反映地形表面較大區(qū)域內(nèi)地形的宏觀特征，在較小的區(qū)域內(nèi)并不具備任何地理和應(yīng)用意義。這些參數(shù)對于在宏觀尺度上的水土保持、土壤侵蝕特征、地表發(fā)育、地貌分類等研究中具有重要的理論意義?；跂鸥馜EM計算宏觀地形因子時，關(guān)鍵在

27、于確定分析半徑的大小。不同地貌類型、不同分辨率的數(shù)據(jù)，計算宏觀地形因子所取的分析半徑大小是不一。因此，確定一個合適的分析窗口半徑或分析區(qū)域，使得求取的宏觀因子能夠準(zhǔn)確反映地面的起伏狀況與水土流失特征，是提取算法的核心步驟和決定信息提取效果與有效性的關(guān)鍵。 地線起伏度地形起伏度是指，在所指定的分析區(qū)域內(nèi)所有柵格中最大高程與最小高程的差?？杀硎緸槿缦鹿剑?(9.7)式中，RFi指分析區(qū)域內(nèi)的地面起伏度，Hmax指分析窗口內(nèi)的最大高程值，Hmin指分析窗口內(nèi)的最小高程值。地形的起伏是反映地形起伏的宏觀地形因子，在區(qū)域性研究中，利用DEM數(shù)據(jù)提取地形起伏度能夠直觀的反映地形起伏特征。在水土流失研究中

28、，地形起伏度指標(biāo)能夠反映水土流失類型區(qū)的土壤侵蝕特征，比較適合區(qū)域水土流失評價的地形指標(biāo)。 地形粗糙度地表粗糙度，一般定義為地表單元的曲面面積S曲面與其在水平面上的投影面積S水平之比。用數(shù)學(xué)公式表達(dá)為：R = S曲面 / S水平 (9.8)地表粗糙度能夠反映地形的起伏變化和侵蝕程度的宏觀地形因子。在區(qū)域性研究中，地表粗糙度是衡量地表侵蝕程度的重要量化指標(biāo)，在研究水土保持及環(huán)境監(jiān)測時研究地表粗糙度也有很重要的意義。實際應(yīng)用時，當(dāng)分析窗口為3×3時，可采用下面近似公式求解：R = 1/cos(S) (9.9)此時，基于DEM的地表粗糙度的提取主要分為以下兩個步驟： 根據(jù)DEM提取坡度因子

29、S； 根據(jù)公式R = 1/cos(S) 計算地表粗糙度。 地表切割深度地表切割深度是指地面某點的鄰域范圍的平均高程與該鄰域范圍內(nèi)的最小高程的差值。可用以下公式表示： (9.10)式中，Di指地面每一點的地表切割深度，Hmean指一個固定分析窗口內(nèi)的平均高程，Hmin指一個固定分析窗口內(nèi)的最低高程。地表切割深度直觀的反映了地表被侵蝕切割的情況，并對這一地學(xué)現(xiàn)象進(jìn)行了量化，是研究水土流失及地表侵蝕發(fā)育狀況時的重要參考指標(biāo)，其提取算法可參照地表起伏度的提取。9.3.2 地形特征分析雖然地表形態(tài)各式各樣，但地形點、地形線、地形面等地形結(jié)構(gòu)的基本特征構(gòu)成了地形的骨架，因此一般的地形特征提取主要是指地形特

30、征點、線、面的提取，并進(jìn)而通過基本要素的組合進(jìn)行地表形態(tài)分析。特征地形要素的提取更多地應(yīng)用較為復(fù)雜的技術(shù)方法，其中山谷線、山脊線的提取采用了全域分析法，成為數(shù)字高程模型地學(xué)分析中很具特色的數(shù)據(jù)處理內(nèi)容。1. 地形特征點提取地形特征點主要包括山頂點（peak）、凹陷點（pit）、脊點（ridge）、谷點（channel）、鞍點（pass），平地點（plane）等。利用DEM提取地形特征點，可通過一個3×3或更大的柵格窗口，通過中心格網(wǎng)點與8個鄰域格網(wǎng)點的高程關(guān)系來進(jìn)行判斷會獲取。即在一個局部區(qū)域內(nèi)，用x方向和y方向上關(guān)于高程z的二階導(dǎo)數(shù)的正負(fù)組合關(guān)系來判斷（見表9.3）。該方法假設(shè)DE

31、M表面為z = f(x,y),但由于真實地表與數(shù)學(xué)表面的差別，在利用該方法在DEM上提取特征點，結(jié)果常產(chǎn)生偽特征點。表9.3 地形特征點類型的判斷表名稱定 義領(lǐng)域高程關(guān)系山頂點（peak）是指在局部區(qū)域內(nèi)海拔高程的極大值點，表現(xiàn)為在各方向上都為凸起。，凹陷點（pit）是指在局部區(qū)域內(nèi)海拔高程的極小值點，表現(xiàn)為在各方向上都為凹陷。，脊點（ridge）是指在兩個相互正交的方向上，一個方向凸起，而另一個方向沒有凹凸性變化的點。，或，谷點（channel）是指在兩個相互正交的方向上，一個方向凹陷，而另一個方向沒有凹凸性變化的點。，或，鞍點（pass）是指在兩個相互正交的方向上，一個方向凸起，而另一個方

32、向凹陷的點。，或，平地點（plane）山頂點是在局部區(qū)域內(nèi)各方向上都沒有凹凸性變化的點。，(i,j)(i-1,j-1)(i-1,j)(i-1,j+1)(i,j+1)(i+1,j+1)(i+1,j)(i+1,j-1)(i,j-1)圖9.15 差分算法示意圖表9.3中的關(guān)于地形特征點的判斷是在局部區(qū)域內(nèi)利用x，y方向的凹凸性判斷的，該判斷法十分適合利用在DEM上判斷地形特征點。在DEM中可以利用差分的方法得到和的值。除上述算法外，在一個3×3的柵格窗口中，也可以直接利用中心格網(wǎng)點與8個鄰域格網(wǎng)點的高程關(guān)系來進(jìn)行判斷地形特征點。具體方法為：假設(shè)有一個如圖9.10所示的3×3窗口。

33、則：如果（Zi，j-1 - Zi，j）（Zi，j-1 - Zi，j）>0（1）當(dāng)Zi，j+1> Zi，j 則VR（i，j）= -1（2）當(dāng)Zi，j+1< Zi，j 則VR（i，j）= 1如果（Zi-1，j - Zi，j）（Zi+1，j - Zi，j）>0等高線山頂點鞍部圖例圖9.16 利用ArcView軟件及DEM數(shù)據(jù)提取的山頂、鞍部（3）當(dāng)Zi+1> Zi，j 則VR（i，j）= -1（4）當(dāng)Zi+1< Zi，j 則VR（i，j）= -1如果（1）和（4）或（2）和（3）同時成立，則VR（i，j）= 2如果以上條件都不成立，則VR（i，j）= 0其中，2.

34、 山脊線和山谷線提取山脊線和山谷線構(gòu)成了地形起伏變化的分界線（骨架線），因此它對于地形地貌研究具有重要的意義。另一方面，對于水文物理過程研究而言，由于山脊、山谷分別表示分水性與匯水性，山脊線和山谷線的提取實質(zhì)上也是分水線與匯水線的提取。這一特性又使得山脊線和山谷線在許多工程應(yīng)用方面有著特殊的意義。在對山脊線、山谷線的提取方法中，基于規(guī)則格網(wǎng)DEM的方法是主要。從原理上來分，主要分為以下四種：（1） 基于圖像處理技術(shù)的原理因為規(guī)則格網(wǎng)DEM數(shù)據(jù)事實上是一種柵格形式的數(shù)據(jù)，可以利用數(shù)字圖像處理中的技術(shù)來設(shè)計算法。利用數(shù)字圖像處理技術(shù)設(shè)計的算法大都采用各種濾波算子進(jìn)行邊緣提取。基于該原理有一種簡單移

35、動窗口的算法，其主要思路是： 計一個2×2窗口以對DEM格網(wǎng)陣列進(jìn)行掃描； 第一次掃描中，將窗口中的具有最低高程值的點進(jìn)行標(biāo)記，自始至終未被標(biāo)記的點即為山脊線上的點； 第二次掃描中，將窗口中的具有最高高程值的點進(jìn)行標(biāo)記，自始至終未被標(biāo)記的點即為山谷線上的點。以上方法存在兩個主要缺陷： 取特征點時必須排除DEM中噪聲的影響； 特征點連接成線時的算法設(shè)計較為困難。（2） 基于地形表面幾何形態(tài)分析原理基于地形表面幾何形態(tài)分析原理的典型算法就是斷面極值法。其基本思想就是地形斷面曲線上高程的極大值點就是分水點，而高程的極小值點就是匯水點。該方法的基本過程為：找出DEM的縱向與橫向的兩個斷面上的

36、極大、極小值點，作為地形特征線上的備選點；根據(jù)一定的條件或準(zhǔn)則將這些備選點劃歸各自所屬的地形特征線。這種算法存在兩個主要缺陷： 由于這種方法對地形特征線上的點的判定與其所屬的地形特征線的判定是分開進(jìn)行的，在確定地形特征線時，全區(qū)域采用一個相同的曲率閾值作為判定地形特征線上點的條件。因此它忽略了每條地形特征線必然存在的曲率變化現(xiàn)象。當(dāng)閾值選擇較大時，會丟失許多地形特征線上的點，導(dǎo)致后續(xù)跟蹤的地形特征線間斷且較短；如果選擇過小，會產(chǎn)生地形特征線上點的誤判，給后續(xù)地形特征線的跟蹤帶來困難。由于該方法只選擇縱、橫兩個斷面來去確定高程變化的極值點，因此它所確定的地形特征線具有一定的近似性，與實際的地形特

37、征線有一定的差異，有時候還會出現(xiàn)遺漏。（3） 基于地形表面流水物理模擬分析原理的算法這種算法的基本思想是：按照流水從高至低的自然規(guī)律，順序計算每一柵格點上的匯水量，然后按匯水量單調(diào)增加的順序，由高到低找出區(qū)域中的每一條匯水線。根據(jù)得到的匯水線，通過計算找出各自匯水區(qū)域的邊界線，就得到了分水線。算法采用了DEM的整體追蹤分析的思路與方法，分析結(jié)果具有系統(tǒng)性好，還便于進(jìn)行相應(yīng)的徑流成因分析，但是，該方法也存在以下兩個明顯的缺陷：由于該算法所計算的匯水量與高程有關(guān)，計算的結(jié)果必然是高程值大的地形特征線上的點的匯水量小，高程值小的地形特征線上的點的匯水量大。因此，可能導(dǎo)致低處非地形特征線上的點的匯水量

38、也較大而被誤認(rèn)為地形特征線上的點；而位于高處的地形特征線上的點會因為匯水量小而被排除；這就造成用該算法所確定的地形特征線（匯水線）的兩端效果很差。由于該算法降格匯水區(qū)域的公共邊界視為分水線，因此它所確定的分水線均為閉合曲線，這與實際的地形特征線（山脊線）不符。（4） 基于地形表面幾何形態(tài)分析和流水物理模擬分析相結(jié)合由于基于地形表面幾何形態(tài)分析原理和基于地形表面流水物理模擬的算法均存在一定的缺陷，因此可將兩者結(jié)合起來以實現(xiàn)地形特征線的提取。這種算法的基本思路是：首先采取較稀疏的DEM格網(wǎng)數(shù)據(jù)，按流水物理模擬算法去提取區(qū)域內(nèi)概略的地形特征線；然后用其引導(dǎo)，在其周圍鄰近區(qū)域?qū)Φ匦芜M(jìn)行幾何分析，來精確

39、的確定區(qū)域的地形特征線。這一算法的關(guān)鍵在于：求出已提取的概略的地形特征線與DEM格網(wǎng)線的交點，在該交點附近的一個小區(qū)域內(nèi)，對DEM數(shù)據(jù)進(jìn)行幾何分析，即找出該區(qū)域內(nèi)與概略的地形特征線正交方向地形斷面上高程變化的極值點，該點即為地形特征線的精確位置。這一算法的基本過程可歸納為：概略DEM的建立；地形流水物理模擬；概略地形特征線提?。坏匦螏缀畏治?；地形特征線精確確定。（5） 平面曲率與坡位組合法即首先利用DEM數(shù)據(jù)提取地面的平面曲率及地面的正負(fù)地形，取正地形上平面曲率的大值即為山脊，負(fù)地形上平面曲率的大值為山谷。該種方法提取的山脊、山谷的寬度可由選取平面曲率的大小來調(diào)節(jié)，方法簡便效果好。圖9.17

40、DEM提取的山脊線圖9.18 DEM提取的山谷線9.3.3 流域分析1. 流域定義降水匯集在地面低洼處，在重力作用下經(jīng)?；蛑芷谛缘匮亓魉旧硭斐傻牟坌喂鹊亓鲃樱纬伤^的河流。河流沿途接納很多支流，水量不斷增加。干流和支流共同組成水系。每一個河流或每一個水系都從一部分陸地面積上獲得補(bǔ)給，這部分陸地面積就是河流或水系的流域，也就是河流或水系在地面的集水區(qū)。把兩個相鄰集水區(qū)之間的最高點連接成的不規(guī)則曲線，就是兩條河流或水系的分水線，因此，流域也可以說是河流分水線以內(nèi)的地表范圍。高程格網(wǎng)和柵格數(shù)據(jù)運(yùn)算用于流域分析，以獲取流域和河網(wǎng)等在水文過程中非常重要的地形要素。2. 流域提取在格網(wǎng)DEM實現(xiàn)流域

41、地形分析，需要順序執(zhí)行如下的步驟：第一步：DEM洼地填充。由于數(shù)據(jù)噪音、內(nèi)插方法的影響，DEM數(shù)據(jù)中常常包含一些“洼地”，“洼地”將導(dǎo)致流域水流不暢，不能形成完整的流域網(wǎng)絡(luò)，因此在利用模擬法進(jìn)行流域地形分析時，要首先對DEM數(shù)據(jù)中的洼地進(jìn)行處理。填充洼地最常用的方法是之一是把 其單元值加高至周圍的最低單元值。第二步：水流方向確定（flow direction）。水流方向是指水流離開格網(wǎng)時的流向。流向確定目前有單流向和多流向兩種，但在流域分析中，常是在3´3局部窗口中找出八個周邊單元中一個最陡的坡度在流域分析中（如圖9.19），水流方向矩陣是一個基本量，這個中間結(jié)果要保存起來，后續(xù)的幾個環(huán)節(jié)都要用到水流方向矩陣。原始DEM與中心格網(wǎng)高差中心單元流向圖9.19 3