楊 旭， 陳建國， 程潭武， 劉 銳

(1.中國地質(zhì)大學(xué)(武漢) 地質(zhì)過程與礦產(chǎn)資源國家重點實驗室， 湖北 武漢 430074;' 2.中國地質(zhì)大學(xué)(武漢) 資源學(xué)院， 湖北 武漢 430074)

1 研究背景

隨著空間信息技術(shù)的不斷發(fā)展，利用數(shù)字高程模型DEM提取流域水系已經(jīng)成為水系提取方法中熱議的話題，同時利用遙感影像提取多時相、多尺度的水系也已經(jīng)被廣泛應(yīng)用，如可以結(jié)合不同季節(jié)影像，進一步識別時令河與常年河等信息。由于二者在水系提取上存在缺陷，因此如何更為準確高效地提取水系仍是有待研究的問題[1]。

隨著GIS技術(shù)的發(fā)展，基于DEM提取水文要素在O′Callaghan等[2]和Tribe[3]提出的傳統(tǒng)理論基礎(chǔ)上得以快速發(fā)展。謝順平等[4]在傳統(tǒng)DEM水系提取的基礎(chǔ)上提出了面向分布式水文過程模擬和流域特征提取的數(shù)字水系模型，對復(fù)雜地形中所存在的洼地有效處理[5-7]使得對各洼地的填平一次完成。Martz等[8]提出了數(shù)字高程流域水系模型(DEDNM)用來表征地表水流特征，并能為計算柵格水流方向、流域分水線等提供依據(jù)。在此思想的基礎(chǔ)上，熊立華等[9]總結(jié)了根據(jù)流域的數(shù)字高程模型自動提取河網(wǎng)水系[10]的集水面積方法，并分析了臨界支撐面積[11-12]的取值對所提取的數(shù)字水系總長度以及平均坡降的影響?；贒EM不能準確提取面狀水系的問題，Turcotte等[13]在DEM的基礎(chǔ)上通過輸入數(shù)字河流和湖泊網(wǎng)絡(luò)(DRLN)來改善對湖泊和平原的處理能力。Hutchinson[14]提出并由Yamazaki等[15]改進的利用提取出的矢量河道與原始DEM進行疊加從而達到降低河道所在DEM高程，使得提取出的河網(wǎng)與實際水系相近。該方法雖然簡單高效，但在河道與河岸邊界之間容易產(chǎn)生不確定性高差。

目前基于遙感影像提取水體常用的方法主要有: 多波段譜間結(jié)構(gòu)關(guān)系法[16]、水體指數(shù)法[17]、決策樹法[18]、監(jiān)督分類法[19]以及歸一化差異性水體指數(shù)NDWI[20]等，其中NDWI的應(yīng)用最為廣泛。徐涵秋[21]對構(gòu)成指數(shù)的波長組合進行了修改,提出了改進的歸一化差異水體指數(shù)MNDWI，但由于MNDWI的構(gòu)成采用了中紅外波段,因此不適用于無中紅外波段的傳感器影像。為了提高遙感水系的提取精度[22]，Bhandari等[23]提出了一種利用歸一化植被指數(shù)(NDVI)和歸一化水體指數(shù)(NDWI)進行離散小波變換(DWT)和奇異值分解(SVD)的增強方法。

為了提高從空間數(shù)據(jù)中提取水系的精度，使得提取的水系分布與實際情況更加符合，使用NDWI和波段分析相結(jié)合的方法來提取面狀水域和干流，并使用數(shù)學(xué)形態(tài)學(xué)[24-25]對提出的二值圖進行除噪、細化等處理。通過DEM提取河網(wǎng)時，通過分析流域集水面積閾值與河網(wǎng)密度之間的關(guān)系確定河網(wǎng)集水面積閾值并對偽河道等進行合理剔除。最后以遙感影像提取出的干流為依據(jù)設(shè)置校正控制點，對DEM河網(wǎng)中水體的位置及面狀水域進行修正，從而完成最終流域水系的提取。

2 數(shù)字流域水系綜合提取思路

基于DEM與遙感信息的數(shù)字流域水系提取時技術(shù)路線圖如圖1，其主要內(nèi)容如下：

(1)對遙感影像進行校正、裁剪等預(yù)處理，并采用決策樹的方法提取出主干河道和面狀水域，并利用數(shù)學(xué)形態(tài)學(xué)中膨脹、腐蝕、擊中、邊緣提取及細化等運算方法將提取的水體結(jié)果進行斷線連接、去噪及細化等處理，從而得到連續(xù)的水系；

(2)基于DEM提取河網(wǎng)，在對原始DEM數(shù)據(jù)進行填洼等預(yù)處理后，結(jié)合流域集水面積閾值與河網(wǎng)密度之間的關(guān)系確定河網(wǎng)集水面積閾值，并通過對矢量水系運用疊置運算和鄰域分析等方法對提出的偽河道進行合理剔除；

(3)將遙感影像和DEM在提取水系的優(yōu)點進行互補，綜合遙感水系信息與DEM水系信息，并將兩者水系進行疊加。以遙感影像提取出的主干河流和面狀水域為依據(jù)，修正DEM提取中主干河道以及缺失的面狀水域，從而完成對研究區(qū)水系的提取。

圖1 數(shù)字流域水系綜合提取技術(shù)路線圖

3 數(shù)字流域水系綜合提取方法

3.1 遙感水系信息的提取

本研究區(qū)選取滄源流域，使用的是高分一號影像，拍攝于2015年2月6日。對高分遙感影像數(shù)據(jù)首先進行正射校正、圖像裁剪等圖像預(yù)處理，然后對遙感影像進行水系提取(見圖2)。

首先采用歸一化差異水體指數(shù)(NDWI)[26]對本研究區(qū)遙感影像進行水系提取，歸一化差異水體指數(shù)：

(1)

式中：ρGreen為綠色光波，在遙感影像中為第2波段；ρNir為近紅外光波，在遙感影像中為第4波段。

由于通過NDWI運算后，水體的值會很低，而其他地物的值會很高[27]，因此可以通過設(shè)置水體提取閾值完成對水系的初步提取，經(jīng)過反復(fù)實驗，選取水系提取閾值為0.1。但實際上通過NDWI提取后的水系圖上還夾雜著地物、陰影和細小水體等噪聲[28](圖3)。

圖2 遙感影像水體提取路線圖

圖3 滄源流域基于NDWI提取的水系圖

然后利用波段分析方法，通過設(shè)置限制性條件去除無用的噪聲。首先對4個波段的亮度直方圖進行統(tǒng)計(表1)。

表1 研究區(qū)4個波段統(tǒng)計特征表 m

由表1可知，Band 4的標準差最大，表明該波段內(nèi)地物的亮度取值距均值的離散程度最大，即地物間的差異表現(xiàn)最大，信息量最大。因此可以作為設(shè)置限制性條件從而剔除無用數(shù)據(jù)的最優(yōu)波段。經(jīng)過多次實驗，認為閾值設(shè)定為1100時剔除陰影等噪聲的效果最佳，見圖4。由圖4可看出，其雖然提取出了水系的基本輪廓，但仍舊存在少量陰影及水系斷線，從而造成提取的水系不連續(xù)的現(xiàn)象。

圖4 滄源流域基于NDWI和波段分析提取的水系圖

基于二值圖圖4中出現(xiàn)的水系不連續(xù)以及有噪聲的情況，采用數(shù)學(xué)形態(tài)學(xué)的方法對提取出的二值圖進行改進處理。數(shù)學(xué)形態(tài)學(xué)基本運算主要有膨脹和腐蝕：

膨脹運算：

(2)

腐蝕運算：

AΘB={x|(B)x?A}

(3)

式中：A為待處理圖像；B為結(jié)構(gòu)元素；z為B的映射位移；x為B的位移。

以上兩種形態(tài)學(xué)基本運算可以派生出擊中運算、邊界提取及細化處理運算：

擊中運算：

A?B=(AΘB1)∩(AcΘB2)

(4)

其中B=(B1,B2)

邊界提?。?/p>

C(A)=A-(AΘB)

(5)

細化處理：

C(A)=A-(A?B)

(6)

式中：C(A)為處理后圖像；Ac為A的補集。

基于以上數(shù)學(xué)形態(tài)學(xué)運算，對水系二值圖像進行膨脹處理，從而實現(xiàn)對斷線水系進行連接，從而形成連續(xù)水系；然后對膨脹后的圖像進行腐蝕處理，從而實現(xiàn)對圖中噪聲的剔除，同時在剔除噪聲后還可以很好地保持水系輪廓；再對腐蝕處理后的圖像進行細化處理，使得提取出的水系更為平滑，最終將處理后的二值圖進行邊緣檢測處理，得出最終的遙感水系圖,如圖5。

圖5 滄源流域基于數(shù)學(xué)形態(tài)學(xué)處理后的水系圖

3.2 基于DEM數(shù)字流域水系信息的提取

對于研究區(qū)ASTER DEM，首先需要對DEM中的洼地、平坦區(qū)進行預(yù)處理，處理流程見圖6。其主要步驟為：(1)計算柵格DEM中的水流流向柵格；(2)通過匯運算判定DEM中是否存在水流流向柵格的異常值，即判斷是否存在洼地、平坦地區(qū)等；(3)對柵格DEM中的洼地進行填洼平滑等預(yù)處理；(4)再次計算填洼后DEM中每個網(wǎng)格單元的水流方向。

圖6 DEM預(yù)處理流程圖

通過以上對DEM進行預(yù)處理的步驟后，需要進一步計算基于單元網(wǎng)格方向的流量，即計算每個網(wǎng)格單元的上游匯集水流面積。繼而通過設(shè)置集水面積閾值完成對流域河網(wǎng)水系的提取。本研究區(qū)主要位于溝壑區(qū)，可以認為河道起源于坡面由凸向凹的轉(zhuǎn)變[29]，隨著閾值的增加，河道的變化也會趨于平緩，因此可以通過河網(wǎng)密度與集水面積閾值之間的關(guān)系(見表2)進行曲線擬合確定合理的集水面積閾值。圖7為集水面積閾值與河網(wǎng)密度及其一階導(dǎo)數(shù)的關(guān)系曲線，由圖7(a)可知：隨著流域集水面積閾值的增加，河網(wǎng)密度持續(xù)遞減，由于圖7(a)的其臨界轉(zhuǎn)折點變化并不顯著，因此作出其一階導(dǎo)數(shù)曲線圖7(b)，在集水面積閾值設(shè)置大約為165 km2時表現(xiàn)為曲線遞減變緩的轉(zhuǎn)折點，即認為165 km2是坡面由凸向凹轉(zhuǎn)變的轉(zhuǎn)折點，也是該區(qū)域較為合理的流域集水閾值。

表2 研究區(qū)集水面積閾值與河網(wǎng)密度關(guān)系表

該流域地勢較陡，偽河道較少，但在下游出口處支流上可能存在偽河道，如圖8(a)。由于偽河道只存在于支流中，因此不考慮主流，此時可將支流河道DEM與原始DEM進行疊置運算，從而得到支流河道的高程，再對支流河道進行領(lǐng)域分析，即如果是位于平坦區(qū)的偽河道，其支流高程與鄰域單元高程差別較小，從而做出偽河道在支流中所在位置并剔除，從而得到圖8(b)。

4 基于DEM與遙感信息數(shù)字流域水系的綜合提取

通過DEM提取水系時，只能提取河網(wǎng)，不能對面狀水域進行精確地提取，并可能在面狀水域的存在的地方出現(xiàn)異常(圖9(a))；在通過遙感影像進行水系提取時，只能提取水系的主干河道以及面狀水域，對于細小水系并不能很好地提取。

利用DEM數(shù)據(jù)根據(jù)溝谷提取水系，對細小水系更為敏感；遙感影像主要根據(jù)水系在波段中的特征識別水系，對大面積水域更為敏感，綜合二者在水系提取上的優(yōu)點進行互補。由于遙感影像提取出的主干河流和面狀水系較為準確，因此以遙感影像提取出的主干河流和面狀水域為依據(jù)，對DEM提取出的水系進行修正。首先以遙感水系和DEM水系信息中分別以主干河道為基準，在河道轉(zhuǎn)彎處設(shè)置校正控制點，然后以遙感水系為基準對DEM水系中控制點圖層進行配準，完成對DEM主干河流的修正，即完成模擬地表徑流從沿岸流向主干河道，以確保DEM河網(wǎng)水系的主干河道與實際相符；在面狀水系處，綜合遙感水體信息，剔除DEM水系位于面狀水系處的異常，從而最終完成對研究區(qū)水系的提取(圖9(b))。

圖7 集水面積閾值與河網(wǎng)密度及其一階導(dǎo)數(shù)關(guān)系曲線圖

圖9 基于DEM與遙感影像提取的水系

5 結(jié) 論

(1)基于遙感影像提取水系時將波段運算作為NDWI的限制條件，同時利用數(shù)學(xué)形態(tài)學(xué)的方法對提取出的水系二值圖進行連接、去噪、邊緣檢測等處理，從而實現(xiàn)對遙感影像中水系的提取。

(2)對基于DEM數(shù)字高程提取水文要素的方法進行改進，通過流域集水面積閾值與河網(wǎng)密度之間的關(guān)系確定河網(wǎng)集水面積閾值，并通過疊置運算和鄰域分析對偽河道和過短支流進行合理剔除。

(3)基于DEM與遙感信息方法對數(shù)字流域水系進行提取，二者可以相互彌補各自在水系提取中的缺陷，從而使得提取出的成果可以較為準確地反映研究區(qū)水系的實際分布情況。通過與研究區(qū)地質(zhì)底圖中的水系圖進行對比，結(jié)果較為準確，從而表明采用該方法提取水文要素是合理、可行的。