基于GIS的黃土高塬溝壑區(qū)硯瓦川流域地形特征提取與分析*

方磊,劉文兆,2,李懷有

(1.西北農(nóng)林科技大學(xué)資源環(huán)境學(xué)院,陜西楊陵712100;2.中國科學(xué)院水利部水土保持研究所,陜西楊陵712100;3.黃河水利委員會西峰水土保持科學(xué)試驗站,甘肅 慶陽745000)

地形形態(tài)分布特征對地表徑流和土壤侵蝕、植被分布和土地利用、水土保持措施布設(shè)等具有重要影響[1],它直接影響著地表物質(zhì)的遷移與能量的轉(zhuǎn)換,標(biāo)志著地表基本環(huán)境格局的形成與發(fā)展。

基于數(shù)字高程模型(Digital Elevation M odel,簡稱DEM)的數(shù)字地形分析技術(shù)經(jīng)過40多年的發(fā)展,其諸多基礎(chǔ)理論問題得到了深入的研究[2],而利用GIS中DEM的自動提取技術(shù)已經(jīng)成為地形因子的主要提取手段。同時DEM的應(yīng)用也從簡單的地形要素提取擴展到更廣的范圍,如輔助建立水文模型、土壤侵蝕模型和洪水預(yù)報模型等。目前國內(nèi)學(xué)者對基于DEM的數(shù)字地形分析研究主要是集中在地形因子提取與定量分析以及DEM信息容量等方面。如湯國安、趙牡丹等[3-4]系統(tǒng)研究了坡度分級方法并對陜北不同地貌區(qū)進行坡度提取、分級,最終對不同地面坡譜進行了統(tǒng)計分析,劉學(xué)軍等[5-9]將DEM誤差分為DEM模型誤差以及空間要素提取算法誤差,并在考慮DEM誤差自相關(guān)的前提下對基于DEM 的空間因素(包括坡度、坡向、河網(wǎng)等)的不同提取算法進行了系統(tǒng)的誤差分析,陳楠等對不同尺度DEM做了信息容量比較,定量化研究了信息損失量的規(guī)律性[10-11]。但以上諸多研究多集中在陜北丘陵溝壑區(qū),對高塬溝壑區(qū)研究較少,尤其是對位于隴東,有天下第一塬號稱的董志塬相關(guān)研究更是鮮有報道。高塬溝壑區(qū)地貌復(fù)雜,塬高、溝深、坡陡,地勢起伏較大,是水土流失嚴重地區(qū)之一,其中水蝕和重力侵蝕強烈,所以如何正確快速掌握地形特征規(guī)律,為當(dāng)?shù)叵嚓P(guān)部門尤其是水保部門提供決策依據(jù)是亟待解決的問題。

基于以上情況,本研究選擇典型高塬溝壑區(qū)——硯瓦川流域為研究區(qū)域,在眾多空間因子中選取了3個最能突出反映該地區(qū)地形特征的指標(biāo)因子,指標(biāo)選取基本原則就是這些地形因子能從不同角度反映該區(qū)域地形特征,獨立性強,最終確定坡度、坡長和剖面曲率三個地形因子作為該區(qū)域地形特征的主要指標(biāo),因為坡度大小直接影響著地表物質(zhì)的流動及能量轉(zhuǎn)換的規(guī)模與強度,坡長能夠影響坡面徑流的流速和流量,進而影響土壤侵蝕強度,而剖面曲率能夠反映出坡度的變化率。最后結(jié)合GIS,以數(shù)字地形分析技術(shù)為手段,建立相應(yīng)地形因子的數(shù)字地面模型(Digital Terrain M odel簡稱DTM),并統(tǒng)計分析各因子分布特點與規(guī)律,望給當(dāng)?shù)亟?jīng)濟建設(shè)提供支持,為后續(xù)相關(guān)研究提供參考。

1 材料與方法

1.1 研究區(qū)域概況

硯瓦川流域系涇河水系,位于甘肅省慶陽市西峰區(qū)境內(nèi) ,107°37′-107°55′E,35°31′-35°44′N,流域中軸線呈NW-ES走向,總面積約368 km2,該區(qū)屬半濕潤季風(fēng)氣候區(qū),年均溫為8.1℃,年均降水量為579 mm,最多年降水量為852.8 mm(1964年),最少年降水量336.7 mm(1986年),7-9月約占全年降水量的60%左右,并多以暴雨形式出現(xiàn),大量降水以地表經(jīng)流的形式匯入江河,發(fā)生暴洪的機率大。流域以黑壚土為主,占流域面積的70%,主要分布在塬面,土層厚約250 cm,其上松下粘,耕性良好,適耕期長。流域地面海拔高度為 940～1 409m,平均高程 1 243 m,平均溝壑密度2.17 km/km2,溝谷是侵蝕最劇烈的部位,重力侵蝕頻繁發(fā)生,對區(qū)域經(jīng)濟發(fā)展、居住安全以及下游水利工程正常運行帶來不利影響,水土保持任務(wù)艱巨。

1.2 數(shù)據(jù)基礎(chǔ)及預(yù)處理

目前制作DEM的主要方法有兩種,即通過遙感立體像得到DEM,或是基于數(shù)字地形圖通過內(nèi)插的方法得到DEM,后者是迄今為止最為成熟和經(jīng)濟實用的DEM建立方法,本文綜合考慮該區(qū)域的數(shù)據(jù)可獲性,選擇后者建立本研究區(qū)域的DEM,所使用到的基礎(chǔ)數(shù)據(jù)為覆蓋研究區(qū)域的4幅1∶5萬國家標(biāo)準(zhǔn)地形圖,基本參數(shù)為:54北京坐標(biāo)系,56黃海高程系,等高距為20 m,使用投影為高斯-克呂格投影。

一般來說,地形圖掃描后都有一定變形,必須進行預(yù)處理,也就是掃描糾正,并對糾正后的柵格圖進行檢查,才能得到矢量化準(zhǔn)確的數(shù)字等高線。本研究采用A rcGIS 9.2提供的Georeferencing工具進行圖像糾正,先在每幅地形圖上選取12個均勻分布的控制點,為方便計算,本研究選取的都是格網(wǎng)交點,然后建立鏈接表,并逐個對控制點進行誤差檢查,合格后對圖像進行重采樣,最后得到經(jīng)過精確糾正的掃描圖。接下來使用Geoway 3.5對該掃描圖進行等高線、高程點、河流、湖泊以及流域邊界線的矢量化,通過要素拼接最終得到硯瓦川流域完整的數(shù)字地形圖。

1.3 建立水文增強型DEM

DEM是一組有序數(shù)據(jù)陣列對地表高度空間變化的定量表現(xiàn)。是一定范圍內(nèi)規(guī)則格網(wǎng)點的平面坐標(biāo)(X,Y)及其高程(Z)的數(shù)據(jù)集,它主要是描述區(qū)域地貌形態(tài)的空間分布,可通過等高線或相似立體模型進行數(shù)據(jù)采集(包括采樣和量測),然后進行數(shù)據(jù)內(nèi)插而形成的。根據(jù)DEM對水文特征和地貌關(guān)系的表現(xiàn)能力,可將DEM分為2類,即水文增強型DEM(hydrologically enhanced DEM s/drainage enforced DEM s)[12],也稱為水文地貌關(guān)系正確DEM(hydrologically correct DEM s,H c-DEM)[13],以及非水文增強型DEM。本研究采用澳大利亞國立大學(xué)Hutchinson教授開發(fā)的ANUDEM軟件[13-14]來建立水文增強型DEM(H c-DEM),考慮到本研究區(qū)域為塬區(qū),相對于丘陵溝壑區(qū)地形平坦,用ANUDEM軟件插值時設(shè)置的三個關(guān)鍵參數(shù)[15]分別為:分辨率10m,第二糙度系數(shù)0.5,迭代次數(shù)35,最終得到水文增強的DEM:H c-DEM,圖1是基于H c-DEM建立的本研究區(qū)域的三維示意圖。

1.4 地形指標(biāo)計算

坡度是指地球表面任一點通過該點的切平面與水平地面的夾角,表示地表面在該點的傾斜程度,是一個重要的地形因子。坡度大小直接影響著地表物質(zhì)的流動及能量轉(zhuǎn)換的規(guī)模與強度。坡度的計算如公式(1),其實現(xiàn)過程是在3×3的DEM 柵格窗口中進行的,窗口在DEM數(shù)據(jù)矩陣中連續(xù)移動后完成整幅圖的計算工作。

式中:fx,fy——南北和東西方向上的偏導(dǎo)數(shù)。目前對 fx,fy的計算有二階差分、三階不帶權(quán)差分(Sharp nack算法)、三階反距離平方權(quán)差分(Horn算法)、三階反距離權(quán)差分、Frame差分、簡單差分等6種常見算法,對于這6種算法的分析比較,國內(nèi)學(xué)者已經(jīng)做了大量的研究[5,16-17],最終,本研究確定使用三階反距離平方權(quán)差分(Horn算法),該算法已經(jīng)集成到A rcGIS中,利用該算法提取的坡度圖如附圖4。在實際應(yīng)用中,坡度有兩種表達方法,一種是坡度即水平面與地面表面之間的夾角;另一種是坡度百分比即高程增加量與水平增加量之比的百分數(shù),本研究采用第一種方法來表達坡度。剖面曲率是對地面坡度沿最大坡降方向地面高程變化率的度量,通過在A rcGIS下對坡度求“slope”來實現(xiàn)提取(附圖5)。

坡長通常是指在地面上沿水流方向到其流向起點間的最大地面距離在水平面上的投影長度。坡長因子是水土保持上的重要因子之一,坡長通過影響坡面徑流的流速和流量,影響水流挾沙力,進而影響土壤侵蝕強度,因而是定量計算土壤流失的重要指標(biāo),也是土壤侵蝕模型中的必要參數(shù)。在不考慮其他影響因素時,坡長越大,徑流速度就愈大,匯集的徑流也愈大,水力侵蝕也愈強[18]。研究坡長對產(chǎn)流產(chǎn)沙過程對建立土壤侵蝕預(yù)報模型、制訂有效的水土保持體系有重要的理論價值。本研究對坡長的提取是基于Hickey和Van Remortel等提出的基于DEM提取區(qū)域坡度坡長因子的算法,H ickey和Van Rem ortel等編制了相應(yīng)AM L語言程序[19-21],后又改進為利用C++語言編寫的程序,使運行效率大為提高[22]。本研究正是通過運行改進后的程序來實現(xiàn)坡長的提取(附圖6),其流程圖如圖2。

圖2 硯瓦川流域坡長提取流程圖

2 結(jié)果與分析

2.1 流域坡度特征分析

將坡度圖屬性表導(dǎo)出,利用Visual Foxp ro進行頻率統(tǒng)計得到圖6,可知流域的坡度分布主要集中在1°～5°,其在累積頻率曲線上表現(xiàn)為曲線增勢較陡,當(dāng)坡度＞5°時,累積頻率曲線的增勢趨于緩和,且5°以上的各坡度頻率分布均小于3%。流域16%的面積坡度在 1°以下,坡度在9°以下的約占流域面積的50%,＞30°的地區(qū)不到流域總面積的20%。流域坡度分布在 0～78.5°,平均坡度為14.9°,坡度標(biāo)準(zhǔn)差為 14.6,變異系數(shù)為 0.98,說明研究區(qū)的地面較為平坦,溝壑發(fā)育還不完全。

圖3 硯瓦川流域坡度頻率曲線圖

根據(jù)文獻資料記錄,在黃土高原丘陵溝壑區(qū)縣南溝流域[23],基于5m分辨率DEM提取的坡度最小值為0°,最大值為 72.9°,平均坡度為28.9°;榆林地區(qū)典型樣區(qū)[24]10 m分辨率DEM提取的坡度最小值為0°,最大值為 59.8°,平均坡度為20.3°,坡度標(biāo)準(zhǔn)差為10.48,變異系數(shù)為0.52;延川樣區(qū)的平均坡度為 32.27°,坡度標(biāo)準(zhǔn)差為 11.94,變異系數(shù)為0.37。比較結(jié)果知,相對于丘陵溝壑區(qū),本研究區(qū)的坡度分布范圍較廣,平均坡度小,說明區(qū)域內(nèi)平緩坡面積比較大,而坡度標(biāo)準(zhǔn)差和變異系數(shù)均大于丘陵溝壑區(qū),反映了研究區(qū)域的地形起伏波動程度要大于丘陵溝壑區(qū),具體統(tǒng)計數(shù)據(jù)見表1。

表1 研究區(qū)及丘陵溝壑區(qū)若干樣區(qū)地形特征值統(tǒng)計表[23-24]

2.2 流域剖面曲率的特征分析

導(dǎo)出流域的剖面曲率圖屬性表,并在V isual Foxp ro下進行頻率統(tǒng)計(圖4)。本研究區(qū)的剖面曲率分布主要集中在1～3,其中曲率為1的地區(qū)面積最廣,約占流域的29.2%,其在累積頻率曲線上的表現(xiàn)為曲線增勢很陡,當(dāng)曲率＞1時,累積頻率曲線的增勢趨于緩和,4以上各曲率的頻率分布均小于3%;由累積頻率曲線知,流域50%的區(qū)域坡面曲率小于7,＞24的地區(qū)不到流域總面積的20%。流域剖面曲率分布在0～75.4,平均剖面曲率為12,曲率標(biāo)準(zhǔn)差為13.6,變異系數(shù)為1.13。

圖4 硯瓦川流域剖面曲率頻率曲線圖

根據(jù)文獻資料記錄[24],在黃土高原丘陵溝壑區(qū),榆林地區(qū)典型樣區(qū)10 m分辨率DEM提取的剖面曲率最小值為0,最大值為67.8,平均剖面曲率為21.6,曲率標(biāo)準(zhǔn)差為13.2,變異系數(shù)為0.61。比較結(jié)果知,相對于丘陵溝壑區(qū),本研究區(qū)的剖面曲率分布稍廣,而平均剖面曲率不到丘陵溝壑區(qū)一半,說明陡坡區(qū)域相比丘陵溝壑區(qū)較少,曲率標(biāo)準(zhǔn)差非常接近,但變異系數(shù)是丘陵溝壑區(qū)的一倍,表明高塬溝壑區(qū)的坡度變化率的離散程度比較大,地貌起伏的波動程度要大于丘陵溝壑區(qū)。

2.3 流域坡長特征分析

通過導(dǎo)出坡長圖的屬性表,在Visual Foxp ro下進行頻率統(tǒng)計(圖5),本研究區(qū)的坡長分布主要集中在10～150 m,其在累積頻率曲線上表現(xiàn)為曲線增勢較陡,其中坡長為10 m的坡分布最廣,占流域的11.21%,150m以上各坡長的頻率分布均小于2%;由累積頻率曲線知,流域50%的區(qū)域坡長小于80m,當(dāng)坡長大于200m時,累積頻率曲線的趨勢變緩,坡長大于240 m的坡只占流域總面積的20%。坡長較大值均分布在海拔1 300 m以上的高塬區(qū),主要在流域的邊界處,溝壑發(fā)育不完全;坡長較小值主要分布在溝谷區(qū),這些地區(qū)地表破碎。流域坡長分布在0～9 583.1m,平均坡長為181.8m,標(biāo)準(zhǔn)差為350.8,變異系數(shù)為1.93,說明研究區(qū)的地面較為平坦,溝壑發(fā)育還不完全。

圖5 硯瓦川流域坡長頻率曲線圖

根據(jù)黃土高原丘陵溝壑區(qū)文獻資料記錄,在縣南溝流域[23],基于5 m分辨率DEM提取的坡長最小值為2.5 m,最大值為485.2m,平均坡長為33.7m;基于25 m分辨率DEM提取的坡長最小值為12.8 m,最大值為443.8 m,平均坡長為54.2 m;在延河流域10m分辨率的典型樣區(qū)[24],坡長最小值為5m,最大值為1 642.25m,平均坡長為64.16m,坡長標(biāo)準(zhǔn)差為65.30,變異系數(shù)為1.02;榆林樣區(qū)基于10m分辨率DEM提取的坡長,最小值為5 m,最大值為571.3 m,平均坡長為37.05 m,坡長標(biāo)準(zhǔn)差為38.2,變異系數(shù)為1.03。比較結(jié)果知,相對于丘陵溝壑區(qū),本研究區(qū)的坡長分布范圍非常廣,平均坡長是丘陵溝壑區(qū)的3～6倍,坡長標(biāo)準(zhǔn)差也是丘陵溝壑區(qū)的5～9倍,變異系數(shù)約是丘陵溝壑區(qū)的2倍,反映了本研究區(qū)地表在水平方向上的平均轉(zhuǎn)折變化水平較小,坡長平均值較大,說明研究區(qū)的連續(xù)平坦區(qū)域面積較大,這也與本區(qū)域塬面占主導(dǎo)的實際地形特征相一致。

3 結(jié)論

(1)流域的坡度分布主要集中在1°～5°,且5°以上的各坡度的頻率分布均小于3%。

(2)剖面曲率分布主要集中在1～3,其中曲率1約占流域的29.2%,4以上各曲率的頻率分布均小于3%。

(3)坡長分布主要集中在10～150 m,10 m的坡長坡占流域的11.21%,150 m以上各坡長的頻率分布均小于2%。

(4)相對丘陵溝壑區(qū),高塬溝壑區(qū)的坡度分布范圍較廣,平均坡度小,區(qū)域內(nèi)平緩坡面積比較大,但坡度變化率的離散程度比較大,地貌起伏的波動程度要大于丘陵溝壑區(qū)。

這些結(jié)果的獲取有助于區(qū)域水保措施部門因地制宜,制定出合理的水保措施,為流域規(guī)劃方案的制訂提供理論依據(jù)和技術(shù)指導(dǎo),本研究僅考慮了坡度、剖面曲率與坡長因子,并未考慮其他地形因子如坡向、高程極差、高程曲線積分等,需要今后進一步進行研究。

[1] 朱顯謨.中國黃土高原土地資源[M].西安:陜西科學(xué)技術(shù)出版社,1986.

[2] 李志林,朱慶.數(shù)字高程模型[M].武漢:武漢大學(xué)出版社,2001.

[3] 湯國安,宋佳.基于DEM坡度圖制圖中坡度分級方法的比較研究[J].水土保持學(xué)報,2006,20(2):157-160.

[4] 趙牡丹,湯國安,陳正江,等.黃土丘陵溝壑區(qū)不同坡度分級系統(tǒng)及地面坡譜對比[J].水土保持通報,2002,22(4):33-36.

[5] 劉學(xué)軍,任志峰,王彥芳,等.基于DEM的任意方向坡度計算方法[J].地域研究與開發(fā),2009,28(4):139-141.

[6] 劉學(xué)軍,晉蓓,王彥芳.DEM流徑算法的相似性分析[J].地理研究,2008,27(6):1347-1357.

[7] 劉學(xué)軍,卞璐,盧華興,等.顧及DEM誤差自相關(guān)的坡度計算模型精度分析[J].測繪學(xué)報,2008,37(2):200-206.

[8] 劉學(xué)軍,王永君,龔健雅,等.DEM流域網(wǎng)絡(luò)提取算法的誤差特性分析[J].測繪學(xué)報,2007,36(2):224-230.

[9] 劉學(xué)軍,龔健雅,周啟鳴,等.基于 DEM 坡度坡向算法精度的分析研究[J].測繪學(xué)報,2004,33(3):258-263.

[10] 陳楠,王欽敏.基于地形因子信息量的數(shù)字高程模型分辨率的選擇:以黃土高原的研究為例[J].武漢大學(xué)學(xué)報:信息科學(xué)版,2009,34(6):692-695.

[11] 陳楠,林宗堅,湯國安,等.數(shù)字高程模型的空間信息不確定性分析[J].測繪通報,2005(11):14-17.

[12] 楊勤科,M cVicar T R,VanN iel T G,等.ANUDEM和TIN兩種建立DEM方法的對比研究[J].水土保持通報,2006,26(6):84-88.

[13] Hutchinson M F.A new p rocedure for gridding elevation and stream line data with automatic removal of spurious pits[J].Journal of Hyd rology,1989,106(3/4):211-232.

[14] 周買春,黎子浩.數(shù)值地形圖的生成及其水文地貌特征評價[J].水利學(xué)報,2002(2):71-74.

[15] 楊勤科,M cVicar T R,VanNiel T G,等.用 ANUDEM建立水文地貌關(guān)系正確DEM的方法研究[J].測繪科學(xué),2006,31(6):155-157.

[16] 陳楠,王欽敏,湯國安,等.6種坡度提取算法的應(yīng)用范圍分析:以在黃土丘陵溝壑區(qū)的研究為例[J].測繪信息與工程,2006,31(4):20-22.

[17] 劉學(xué)軍,王葉飛,曹志東,等.基于 DEM 的坡度坡向誤差空間分布特征研究[J].測繪通報,2004(12):11-13.

[18] 唐克麗,陳永宗.黃土高原地區(qū)土壤侵蝕區(qū)域特征及綜合治理途徑[M].北京:中國科學(xué)技術(shù)出版社,1990.

[19] Van Remortel R,H am ilton M,Hickey R.Estimating the LS fac tor for RUSLE through iterative slope length p rocessing of digital elevation data[J].Cartography,2001,30(1):27-35.

[20] H ickey R.Slope angle and slope length solutions for GIS[J].Cartography,2000,29(1):1-8.

[21] H ickey R,Sm ith A,Jankow ski P.Slope Length Calcu lations from a DEM within ARC/INFO GRID[J].Computers,Environment and U rban Systems,1994,18(5):365-380.

[22] Van Remortel H M,H ickey R.Calcu lations LS Factor Using DEM dataaccording to RUSLE-based C riteria.www.on linegeographer.com/slope/rusleV4pc.am l2004.

[23] 謝紅霞.延河流域土壤侵蝕時空變化及水土保持環(huán)境效應(yīng)評價研究[D].西安:陜西師范大學(xué),2008.

[24] 張婷.基于DEM的陜北黃土高原多地形因子空間關(guān)聯(lián)特征研究[D].西安:西北大學(xué),2005.