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

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

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

地形形態分布特征對地表徑流和土壤侵蝕、植被分布和土地利用、水土保持措施布設等具有重要影響[1],它直接影響著地表物質的遷移與能量的轉換,標志著地表基本環境格局的形成與發展。

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

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

1 材料與方法

1.1 研究區域概況

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

1.2 數據基礎及預處理

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

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

1.3 建立水文增強型DEM

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

1.4 地形指標計算

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

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

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

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

2 結果與分析

2.1 流域坡度特征分析

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

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

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

表1 研究區及丘陵溝壑區若干樣區地形特征值統計表[23-24]

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

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

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

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

2.3 流域坡長特征分析

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

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

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

3 結論

(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)相對丘陵溝壑區,高塬溝壑區的坡度分布范圍較廣,平均坡度小,區域內平緩坡面積比較大,但坡度變化率的離散程度比較大,地貌起伏的波動程度要大于丘陵溝壑區。

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

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

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

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

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

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

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

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

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

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

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

[11] 陳楠,林宗堅,湯國安,等.數字高程模型的空間信息不確定性分析[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] 周買春,黎子浩.數值地形圖的生成及其水文地貌特征評價[J].水利學報,2002(2):71-74.

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

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

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

[18] 唐克麗,陳永宗.黃土高原地區土壤侵蝕區域特征及綜合治理途徑[M].北京:中國科學技術出版社,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] 謝紅霞.延河流域土壤侵蝕時空變化及水土保持環境效應評價研究[D].西安:陜西師范大學,2008.

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