基于多源DEM數(shù)據(jù)的都江堰灌區(qū)地形因子提取及對比分析

肖文全,曹依帆,秦 濤,羅 尚,翟 星,馮丹禹,趙 輝

(四川水發(fā)勘測設(shè)計(jì)研究有限公司,成都 610072)

0 引言

20世紀(jì)50年代,數(shù)字高程模型(DEM)的概念被提出,它被定義為一種以數(shù)字?jǐn)?shù)組形式顯示地面標(biāo)高的實(shí)體模型[1]。近年來,隨著計(jì)算機(jī)技術(shù)、地理信息技術(shù)(GIS)和遙感技術(shù)(RS)的迅速發(fā)展,DEM數(shù)據(jù)的種類逐漸增多,分辨率也逐步提高?；贒EM數(shù)據(jù)提取的地形因子可以有效地描述地貌的形態(tài)特征,從而推動(dòng)地貌形態(tài)分類的研究進(jìn)展[2]。因此,DEM數(shù)據(jù)在眾多基礎(chǔ)地理信息數(shù)據(jù)中占有極為重要的地位。針對DEM的研究也成為當(dāng)下的熱點(diǎn),高藝偉等[3]基于ASTER GDEM2數(shù)據(jù)提取陜南地區(qū)地形起伏度,并采用均值變點(diǎn)法確定最佳鄰域分析單元;韓萍等[4]基于多源DEM數(shù)據(jù),提取青島市部分地區(qū)小流域水文特征并進(jìn)行對比分析;胡最等[5]依據(jù)1∶25萬DEM數(shù)據(jù),開展湖南省地貌形態(tài)特征分類研究。雖然各專家學(xué)者已對各式DEM展開不同探討,但利用多源DEM數(shù)據(jù)進(jìn)行地形因子提取和對比分析的研究還略顯不足,且已有研究表明,不同分辨率的DEM數(shù)據(jù)所表達(dá)的地形信息的容量和精度等均存在顯著差異[6]。鑒于此,針對同一研究區(qū),利用不同分辨率DEM數(shù)據(jù)提取地形因子并進(jìn)行定量分析具有一定研究價(jià)值。本文以都江堰灌區(qū)作為研究區(qū),首先,基于三種不同分辨率的DEM數(shù)據(jù),提取高程、坡度、坡向、地形起伏度、地表粗糙度、地表切割深度等六個(gè)地形因子;接著,對提取結(jié)果進(jìn)行定性和定量分析;最后,利用皮爾遜相關(guān)性分析篩選出研究區(qū)最佳地形因子。研究結(jié)果可為都江堰灌區(qū)的地貌類型精細(xì)化劃分、災(zāi)害風(fēng)險(xiǎn)評估及土壤侵蝕等提供依據(jù)和參考。

1 研究區(qū)概況與數(shù)據(jù)來源

1.1 研究區(qū)概況

都江堰灌區(qū)位于四川盆地西部,介于東經(jīng)103°29′～105°24′、北緯29°24′～31°29′之間。從地形地貌條件來看,灌區(qū)地貌類型以丘陵為主,丘陵面積占總面積的比例高達(dá)71.1%[7];海拔高度介于210～1 004 m,以龍泉山脈為分界,具有典型的西北高、東南低特征。從氣候水文方面來看,灌區(qū)屬亞熱帶濕潤氣候類型,全年氣候適宜,年均降雨量900～1 240 mm,降雨量受季節(jié)影響較大,呈現(xiàn)春旱、夏洪、秋澇、冬干的特點(diǎn)[8]。目前,都江堰灌區(qū)實(shí)際農(nóng)田灌溉面積達(dá)75.37萬hm2,涉及成都、綿陽等7市40縣(市、區(qū))。研究區(qū)概況如圖1所示。

圖1 研究區(qū)概況

1.2 數(shù)據(jù)來源

本文利用的三種DEM數(shù)據(jù)分別為ASTR GDEM、SRTM DEM及ALOS DEM。其中,ALOS DEM空間分辨率為12.5 m,來源于NASA地球數(shù)據(jù)中心(https：//search.asf.alaska.edu/#/);ASTR GDEM空間分辨率為30 m,SRTM DEM空間分辨率為90 m,均來源于地理空間數(shù)據(jù)云(http：//www.gscloud.cn/)。基礎(chǔ)數(shù)據(jù)的準(zhǔn)確性是提取地形因子的關(guān)鍵,基于GIS軟件平臺,對三種原始DEM數(shù)據(jù)進(jìn)行鑲嵌、裁剪及投影變換等預(yù)處理,以獲取研究區(qū)同一坐標(biāo)系下準(zhǔn)確的基礎(chǔ)數(shù)據(jù)。

2 研究方法

2.1 歸一化

不同的地形因子具有不同的原始值,并且對量綱和量綱單位也具有不同的要求,為了使各地形因子具有可比性,需要采取適當(dāng)?shù)氖侄螌Ω鞯匦我蜃舆M(jìn)行處理,從而將各地形因子的值限定在一定的范圍內(nèi)[9]。因此,為了便于開展地形因子間的相關(guān)性分析,本文對各地形因子的原始值進(jìn)行歸一化處理,使處理結(jié)果全部介于[0,1]。

(1)

式中,X為歸一化值;x為各地形因子原始值;xmin為各地形因子的最小值;xmax為各地形因子的最大值。

2.2 皮爾遜相關(guān)性分析

皮爾遜相關(guān)系數(shù)又叫皮爾遜積矩相關(guān)系數(shù),它是一種可以準(zhǔn)確量算兩個(gè)變量之間的相關(guān)程度的經(jīng)典統(tǒng)計(jì)學(xué)方法,被廣泛地應(yīng)用于國民消費(fèi)結(jié)構(gòu)[10]及災(zāi)害敏感性評估[11]等方面。以變量x和變量y為例,那么x與y之間的相關(guān)系數(shù)計(jì)算公式為：

(2)

式中,R表示變量x和變量y之間的相關(guān)系數(shù);n表示變量x和變量y的個(gè)數(shù)。相關(guān)系數(shù)與具體的相關(guān)程度如表1所示。

表1 皮爾遜相關(guān)系數(shù)與相關(guān)程度

3 地形因子提取與分級

3.1 高程

高程是地貌形態(tài)的最基本要素,DEM所提供的原始數(shù)據(jù)就是高程數(shù)據(jù)[12]。對三種高程數(shù)據(jù)按分辨率從高到低排序,獲得研究區(qū)的高程值分別介于210～1 004 m、192～1 014 m及254～1 042 m,考慮前人對高程的分級標(biāo)準(zhǔn)并結(jié)合研究區(qū)實(shí)際高程分布情況,將高程值分為<400 m、400～600 m、600～800 m、800～1 000 m及>1 000 m等五個(gè)等級,如圖2所示。

(a)12.5 m分辨率 (b)30 m分辨率 (c)90 m分辨率

3.2 坡度

坡度定量反映了地表的陡緩程度,是常見的微觀地形因子之一[13]。利用Arc GIS 10.2軟件的表面分析工具,從三種不同DEM數(shù)據(jù)中提取坡度因子,按分辨率從高到低排序,最終提取結(jié)果為0～69.36°、0～71.28°、0～47.97°。依據(jù)曾佩楓等[1]對坡度的等級劃分,本文將坡度分為<3°、3°～5°、5°～10°、10°～15°、15°～20°、20°～25°、25°～30°及>30°等八個(gè)等級,如圖3所示。

(a)12.5 m分辨率 (b)30 m分辨率 (c)90 m分辨率

3.3 坡向

坡向是地表任意一點(diǎn)切平面的法線矢量在水平面的投影與過該點(diǎn)正北方向的夾角,可通俗地理解為自高而低的方向[1]。與坡度的提取及排序方式相似,將研究區(qū)坡向劃分為平坦-1、北(0～22.5°及337.5°～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°)等共九個(gè)類別,如圖4所示。

(a)12.5 m分辨率 (b)30 m分辨率 (c)90 m分辨率

3.4 地形起伏度

地形起伏度(Rf)是指某一特定區(qū)域內(nèi)最高點(diǎn)高程與最低點(diǎn)高程的差值,能夠較為直觀地反映地表形態(tài)[14],是地貌類型劃分研究中最常用的因子之一,其具體計(jì)算公式為：

Rf=Hmax-Hmin

(3)

最佳分析窗口的確定是影響地形起伏度準(zhǔn)確性的最關(guān)鍵因素,基于前人[2,15-16]研究成果,對于ALOSDEM、ASTR GDEM及SRTM DEM,分別采取13×13窗口,15×15窗口,17×17窗口為最佳分析窗口,其面積分別為0.026 4 km2、0.202 5 km2及2.34 km2。利用Arc GIS 10.2軟件平臺的焦點(diǎn)統(tǒng)計(jì)工具,從三種DEM數(shù)據(jù)中提取地形起伏度因子,按分辨率從大到小排序,最終結(jié)果為0～216 m、0～363 m、0～509 m。依據(jù)以往研究成果的數(shù)字格局地貌研究分類方法,將三種提取結(jié)果劃分為五類,分別為<30 m(平原)、30～50 m(臺地)、50～200 m(丘陵)、200～500 m(小起伏山地)及>500 m(中起伏度山地),如圖5所示。

(a)12.5 m分辨率 (b)30 m分辨率 (c)90 m分辨率

3.5 地表粗糙度

地表粗糙度(R)是指某一固定區(qū)域內(nèi)地球表面積與其投影面積的比值,是反映地表形態(tài)的一個(gè)重要宏觀因子[2]。在提取的坡度數(shù)據(jù)基礎(chǔ)上,利用Arc GIS 10.2軟件的柵格計(jì)算器工具,獲取研究區(qū)地表粗糙度數(shù)據(jù)。結(jié)果分別為1～2.87、1～3.12及1～1.49。參照賈麗娜[17]的分級標(biāo)準(zhǔn),將三種分辨率的地表粗糙度數(shù)據(jù)分為1～1.03、1.03～1.08、1.08～1.14、1.14～1.21、1.21～1.29、1.29～1.39、1.39～1.57、1.57～1.93及>1.93等共九個(gè)等級,如圖6所示。

(a)12.5 m分辨率 (b)30 m分辨率 (c)90 m分辨率

3.6 地表切割深度

地表切割深度(D)定義為某一特定區(qū)域內(nèi)平均高程與最小高程的差值[18]。在提取地表切割深度的過程中,最佳分析窗口的選擇與提取地形起伏度的窗口相同,提取結(jié)果分別為0～112.23、0～204.94及0～285.27。參照梁宏艷等[2]的劃分標(biāo)準(zhǔn),將三種提取結(jié)果分為三個(gè)等級,分別為<30 m(微切割山地)、30～100 m(淺切割山地)及100～500 m(中等切割山地),如圖7所示。

(a)12.5 m分辨率 (b)30 m分辨率 (c)90 m分辨率

4 結(jié)果分析

4.1 地形因子對比分析

為了更直觀地分析三種不同分辨率地形因子的等級劃分結(jié)果,對各地形因子不同等級的柵格數(shù)量、面積、柵格數(shù)量占比及面積占比進(jìn)行統(tǒng)計(jì)。不同分辨率下各等級高程數(shù)據(jù)的柵格數(shù)量占比如表2所示。

表2 研究區(qū)不同分辨率下各等級柵格數(shù)量占比(高程)

由表2可知,對于30 m和90 m分辨率的高程因子來說,各等級所占比例接近,且400～600 m等級的占比最高,達(dá)70%以上。而對于12.5 m分辨率的高程因子,400～600 m所占比例卻不足50%,更多的區(qū)域被劃分為<400 m等級。

不同分辨率下各等級坡度數(shù)據(jù)的柵格數(shù)量占比如表3所示。

表3 研究區(qū)不同分辨率下各等級柵格數(shù)量占比(坡度)

由表3可知,三種分辨率的坡度數(shù)據(jù)在柵格數(shù)量占比上未表現(xiàn)出明顯規(guī)律。除<3°分級的差異較大外,其余各等級所占比例相差不大。

不同分辨率下各等級坡向數(shù)據(jù)的柵格數(shù)量占比如表4所示。

表4 研究區(qū)不同分辨率下各等級柵格數(shù)量占比(坡向)

由表4可知,三種分辨率的坡向數(shù)據(jù),除在平面所占比例相差較大外,其余八個(gè)分級所占比例接近,且均在10%上下波動(dòng)。

不同分辨率下各等級地形起伏度數(shù)據(jù)的柵格數(shù)量占比如表5所示。

表5 研究區(qū)不同分辨率下各等級柵格數(shù)量占比(地形起伏度)

由表5可知,30 m和90 m分辨率的地形起伏度數(shù)據(jù)在對平原和丘陵的劃分結(jié)果相似,而12.5 m和90 m分辨率的地形起伏度數(shù)據(jù)對臺地的劃分結(jié)果相似。此外,90 m分辨率的地形起伏度數(shù)據(jù)中,有極少數(shù)區(qū)域被劃分為中起伏度山地。

不同分辨率下各等級地表粗糙度數(shù)據(jù)的柵格數(shù)量占比如表6所示。

表6 研究區(qū)不同分辨率下各等級柵格數(shù)量占比(地表粗糙度)

由表6可知,三種分辨率的地表粗糙度數(shù)據(jù)中,各等級柵格數(shù)量所占比例差別不大,90%以上的柵格屬于1～1.08之間。

不同分辨率下各等級地表切割深度數(shù)據(jù)的柵格數(shù)量占比如表7所示。

表7 研究區(qū)不同分辨率下各等級柵格數(shù)量占比(地表切割深度)

由表7可知,對于30 m分辨率和90 m分辨率的地表切割深度劃分結(jié)果,各級占比差異不大。其中,70%以上的區(qū)域被劃分為微切割山地,20%以上的區(qū)域被劃分為淺切割山地,僅有極少區(qū)域被劃分為中等切割山地。而對于12.5 m分辨率的地形切割深度分級結(jié)果來說,僅有微切割山地和淺切割山地兩種類型,且微切割山地的占比高達(dá)94.78%。

4.2 最佳地形因子確定

為了確定研究區(qū)的最佳地形因子,利用皮爾遜相關(guān)系數(shù)對各因子進(jìn)行相關(guān)性分析,具體步驟為：①利用Arc GIS 10.2軟件的隨機(jī)選點(diǎn)功能,生成2 000個(gè)獨(dú)立的隨機(jī)點(diǎn);②提取各地形因子的原始值至隨機(jī)點(diǎn);③利用公式(1),將提取的原始值進(jìn)行歸一化處理;④基于SPSS 25軟件,對歸一化過后的數(shù)據(jù)進(jìn)行皮爾遜相關(guān)分析。各因子的相關(guān)性分析結(jié)果如表8所示。

表8 地形因子相關(guān)性分析結(jié)果

由表8的分析結(jié)果并結(jié)合表1可知,坡向與其他地形因子間的相關(guān)系數(shù)均小于0.2,即極低相關(guān)性,表明坡向與其他地形因子相比,對研究區(qū)地形特征差異表達(dá)不明顯,應(yīng)予以剔除。除此之外,切割深度與起伏度、切割深度與坡度的相關(guān)系數(shù)均大于0.8,即極強(qiáng)相關(guān),表明三者表達(dá)的地形信息有所重復(fù)?？紤]到起伏度和坡度已被大量研究確定為地形劃分的主要因子,因此,選擇保留起伏度和坡度而刪除切割深度。最終,研究區(qū)最佳地貌類型劃分的地形因子為高程、坡度、地形起伏度及地表粗糙度。

5 結(jié)論與討論

以都江堰灌區(qū)為研究區(qū),首先,基于三種不同分辨率的DEM數(shù)據(jù),提取高程、坡度、坡向、地形起伏度、地表粗糙度及地表切割深度等六個(gè)地形因子;接著,對提取結(jié)果進(jìn)行等級劃分和統(tǒng)計(jì);然后,利用皮爾遜相關(guān)性分析確定研究區(qū)地貌分類的最佳地形因子。研究結(jié)果如下：

(1)由不同分辨率提取的研究區(qū)地形因子結(jié)果存在明顯差異。對高程而言,90%以上的柵格分布于<600 m等級;對坡度而言,除<3°等級相差較大外,其余等級差距較小;對坡向而言,平面等級的占比相差在10%以內(nèi),其余等級差距不大;對地形起伏度而言,12.5 m分辨率DEM的提取結(jié)果與其他兩種提取結(jié)果差異較大,所有的區(qū)域都被劃分為平原、臺地及丘陵。呈現(xiàn)出分辨率越低,提取結(jié)果種類越多的趨勢;對地表粗糙度而言,各等級柵格數(shù)量所占比例接近;對地表切割深度而言,12.5 m分辨率DEM的提取結(jié)果僅劃分出微地形切割和淺地形切割兩種,且微地形切割的比例高達(dá)95%。

(2)由皮爾遜相關(guān)性分析結(jié)果可知,都江堰灌區(qū)的最佳地形因子為高程、坡度、地形起伏度及地表粗糙度。