陳珂，黃小羽，羅明良

(西華師范大學(xué) 國(guó)土資源學(xué)院，四川 南充 637009)

1 引 言

起伏度是一定面積內(nèi)最高和最低點(diǎn)之高差[1]，其研究源于前蘇聯(lián)學(xué)者提出的切割深度，現(xiàn)已成為劃分地貌類型的重要指標(biāo)，也是地貌圖編制的基本劃分依據(jù)[2]。利用數(shù)字高程模型作為基本信息源對(duì)區(qū)域地勢(shì)起伏度的研究逐漸興盛起來(lái)[3]，業(yè)已完成不同比例尺起伏度最佳統(tǒng)計(jì)單元求解[4-7]。起伏度研究在諸多方面應(yīng)用廣泛，如區(qū)域人口分布和經(jīng)濟(jì)發(fā)展[8-9]、區(qū)域土地時(shí)空變化和土地利用/土地覆蓋空間結(jié)構(gòu)[10]，以及用于水土流失定量評(píng)價(jià)及區(qū)域滑坡災(zāi)害評(píng)價(jià)等[5，11]。

基于DEM提取地形起伏度的關(guān)鍵在于確定一個(gè)合適的分析窗口，使得在宏觀范圍內(nèi)求取的起伏度能夠準(zhǔn)確反映地面的起伏狀況[1]。用地形起伏度分析使用窗口遞增方法，得到平均起伏度趨于平穩(wěn)的曲線，通過(guò)定性判別或均值變點(diǎn)等方法得到最佳分析單元[3]?，F(xiàn)有研究已經(jīng)表明窗口遞增法適用于大區(qū)域、較低分辨率的DEM[12]；面臨縣域尺度、較高分辨率的DEM，平均起伏度趨向不斷增加，如何確定窗口遞增法的最大分析窗口？ASTER GDEM作為全球公開共享數(shù)據(jù)，具有約30m水平分辨率[13]，在地貌研究中具有越來(lái)越廣泛的應(yīng)用空間[13-14]，探索基于30m×30m柵格分辨率DEM提取縣域尺度地形起伏度的方法，并結(jié)合實(shí)際地形對(duì)其結(jié)果進(jìn)行地學(xué)解釋和分析顯得尤為必要。

本文以ASTER GDEM數(shù)據(jù)，以川東平行嶺谷區(qū)的華鎣市為研究區(qū)域，使用窗口遞增方法，考察最大起伏度隨分析窗口變化趨勢(shì)，確定最大分析窗口，得到最大起伏度隨分析窗口趨于平穩(wěn)的曲線；基于均值變點(diǎn)法選取適宜分析窗口，得到華鎣市地形起伏度。研究對(duì)于深入認(rèn)識(shí)華鎣市地貌形態(tài)有借鑒意義。

2 研究樣區(qū)與數(shù)據(jù)

2.1 研究樣區(qū)概況

研究樣區(qū)為華鎣市，部分位于國(guó)家級(jí)地質(zhì)公園內(nèi)，地處四川盆地川東平行嶺谷區(qū)的華鎣山中段西緣、渠江東岸，海拔介于206m～1704m；高登山為最高點(diǎn)。襄渝鐵路沿華鎣山腳縱貫?zāi)媳保瑢⑷凶匀环譃榇笾孪嗟鹊臇|西兩大部分，東西兩部分的地形地貌有著明顯差異。西部低丘廣布，地勢(shì)偏低，以丘陵和低丘為主，深丘較少；平壩、臺(tái)地散布其間；多紫色土，土壤肥沃，灌溉便利，是主要的農(nóng)耕區(qū)。東部山地為主，地勢(shì)較高。東南部是華鎣山中段的西翼部分，低山為主，海拔超過(guò)1000m的華鎣山山脊一帶出現(xiàn)中山地貌。東部山區(qū)石灰?guī)r廣布，喀斯特地貌突出，有成片的石林和眾多溶洞，如瓦店的石林、溪口的仙鶴洞等。土壤為黃壤，土層薄，肥力低，農(nóng)耕條件差，自然資源豐富。

2.2 實(shí)驗(yàn)數(shù)據(jù)

本文所采用的數(shù)據(jù)是ASTER GDEM(先進(jìn)星載熱發(fā)射和反射輻射儀全球數(shù)字高程模型)數(shù)據(jù)，該數(shù)據(jù)是由美國(guó)國(guó)家航空航天局(NASA)與日本經(jīng)濟(jì)產(chǎn)業(yè)省(METI)于2009年6月共同推出的，是目前較精確的地球陸地海拔地形數(shù)據(jù)。空間分辨率為1×1弧秒(約30m×30m)，每個(gè)分片包含3601行×3601列。全球范圍內(nèi)，在置信度為95%時(shí)，ASTER GDEM的垂直精度約為20m，水平精度達(dá)30 m[13]。

3 研究方法與步驟

3.1 起伏度求解

基于ArcMap空間分析(Spatial Analyst)提取地形起伏度，使用鄰域分析(Neighborhood Statistics)實(shí)現(xiàn)。分析中使用n×n(n=3，5，7，…，83)分析窗口，計(jì)算窗口內(nèi)最大高程與最小高程之差；并求取正方形窗口下的平均地形起伏度及最大起伏度，并依此數(shù)據(jù)對(duì)其進(jìn)行擬合，得到網(wǎng)格單元面積與平均地形起伏度的關(guān)系。

3.2 最佳統(tǒng)計(jì)單元分析

利用窗口遞增法可以得到地形起伏度隨分析窗口變化曲線，曲線形態(tài)變化有助于推斷起伏度最佳統(tǒng)計(jì)單元。為了避免目視解譯帶來(lái)的主觀性，沿用了均值變點(diǎn)分析法，以合理計(jì)算出曲線上由陡變緩的點(diǎn)，即最佳統(tǒng)計(jì)單元。均值變點(diǎn)分析法是一種數(shù)理統(tǒng)計(jì)學(xué)上的對(duì)非線性數(shù)據(jù)進(jìn)行處理的方法，其計(jì)算步驟如下[3]：

(1)令i=2，…，N，對(duì)每個(gè)i將樣本分為兩段：X1，X2，Xi-1和Xi，Xi+1，…，XN。計(jì)算每段樣本的算術(shù)平均值Xi1和Xi2及統(tǒng)計(jì)量Si：

(1)

(2)計(jì)算統(tǒng)計(jì)量：

(2)

(3)計(jì)算期望值：E(S-Si)，i=2，3，…，N；其中S為原始樣本的統(tǒng)計(jì)量，Si為樣本分段后的統(tǒng)計(jì)量。

(3)

基于上述方法，作出不同窗口下S與Si差值的變化曲線，找到曲線上S與Si差值的最大值點(diǎn)，通過(guò)對(duì)比分析，得到適合本區(qū)域的最佳統(tǒng)計(jì)單元面積。

4 實(shí)驗(yàn)結(jié)果

4.1 不同窗口與平均、最大地形起伏度的對(duì)應(yīng)關(guān)系

最大分析窗口的確定是窗口遞增法的重要步驟之一。華鎣市面積僅470km2，區(qū)域較?。谎芯渴褂玫腄EM分辨率較高；屬于典型的小區(qū)域、較高分辨率條件下的地形起伏度研究。研究表1及圖1(a)可以發(fā)現(xiàn)，在分析窗口遞增條件下，地形起伏度均值由陡變緩趨勢(shì)不明顯。

按照地貌發(fā)育的基本理論，該區(qū)域存在使最大高差達(dá)到相對(duì)穩(wěn)定的最佳分析面積[4]?；诖?，研究擯棄了根據(jù)平均起伏度判定終止分析窗口的一貫做法，轉(zhuǎn)而考察區(qū)域最大起伏度何時(shí)趨于穩(wěn)定。從表1和圖1(b)可以看出，分析窗口擴(kuò)大至79～83時(shí)，最大地形起伏度保持1182m不變。據(jù)此，分析使用的最大分析窗口為83×83，面積約為6.2km2。圖1(a)為分析窗口與平均地形起伏度對(duì)應(yīng)關(guān)系曲線；圖1(b)為分析窗口與最大地形起伏度對(duì)應(yīng)關(guān)系曲線。

圖1 分析窗口與平均、最大起伏度對(duì)應(yīng)關(guān)系圖

表1分析窗口與平均起伏度(m)、最大起伏度(m)及窗口面積(104m2)對(duì)應(yīng)關(guān)系

窗口大小平均起伏度最大起伏度窗口面積窗口大小平均起伏度最大起伏度窗口面積3×316.08200.000.8145×45199.66936.00182.255×531.10289.002.2547×47205.41958.00198.817×744.88353.004.4149×49211.05974.00216.099×957.52396.007.2951×51216.58988.00234.0911×1169.18431.0010.8953×53222.01995.00252.8113×1380.03466.0015.2155×55227.351005.00272.2515×1590.18518.0020.2557×57232.631046.00292.4117×1799.75552.0026.0159×59237.811061.00313.2919×19108.81589.0032.4961×61242.911075.00334.8921×21117.43628.0039.6963×63247.941087.00357.2123×23125.67703.0047.6165×65252.881096.00380.2525×25133.57756.0056.2567×67257.741110.00404.0127×27141.16768.0065.6169×69262.521137.00428.4929×29148.47789.0075.6971×71267.221151.00453.6931×31155.53814.0086.4973×73271.851163.00479.6133×33162.35833.0098.0175×75276.411175.00506.2535×35168.97847.00110.2577×77280.901181.00533.6137×37175.40860.00123.2179×79285.331182.00561.6939×39181.68872.00136.8981×81289.711182.00590.4941×41187.80894.00151.2983×83294.021182.00620.0143×43193.79916.00166.41

4.2 變點(diǎn)法尋求最佳統(tǒng)計(jì)單元

從圖1(b)可見，擬合曲線存在由陡變緩的拐點(diǎn)，該點(diǎn)所對(duì)應(yīng)的區(qū)域面積即所要求取的最佳分析區(qū)域面積。為進(jìn)一步精確判斷拐點(diǎn)的大小，研究采用均值變點(diǎn)法對(duì)其進(jìn)行分析(圖2)。從圖2可以看出，系列分析窗口3×3、5×5，…，83×83對(duì)應(yīng)的S與Si差值最大點(diǎn)為第13個(gè)點(diǎn)，其窗口大小為27×27，正方形窗口邊長(zhǎng)約為810m，對(duì)應(yīng)的區(qū)域面積約為0.65km2。這一結(jié)果略小于朱紅春等基于1∶1萬(wàn)DEM得到的黃土高原綏德樣區(qū)的1km2[12]，大于程維明等基于1∶10萬(wàn)DEM得到的0.45km2[7]。

4.3 華鎣山地形起伏度分析

根據(jù)數(shù)字地貌圖制圖規(guī)范，地形起伏度可劃分為7級(jí)：<30m、30m～70m、70m～200m、200m～500m、500m～1000m、1000m～2500m、>2500m。本文按此標(biāo)準(zhǔn)對(duì)華鎣山地形起伏度進(jìn)行分級(jí)，用上述分析結(jié)果，以27×27分析窗口作為最佳統(tǒng)計(jì)單元，得到華鎣山地形起伏度分級(jí)圖(圖3)。

從華鎣山整體范圍來(lái)看東部為山區(qū)，西部為阡-陵平壩。中低山地集中在東部地區(qū)，中東部地區(qū)存在巖溶中山、低山及巖溶槽谷組成；西部多丘陵、平壩等。其中東北部低起伏區(qū)為華鎣山天池湖所在地；中部及南部起伏度超過(guò)500m的地方為華鎣山山脊部分。華鎣山處于川東帚狀構(gòu)造帶，由一系列近于平行的狹長(zhǎng)不對(duì)稱的箱狀高背斜組成，軸部由于斷裂破碎帶存在，強(qiáng)烈差異侵蝕形成陡峻山體和深窄谷底，在圖3中體現(xiàn)出較大的起伏度，最大起伏度768m，其中中低山部分約為全區(qū)面積的40.85%；這與華鎣市全市470km2，中低山面積198km2較為相符[15]。從圖3中還可以看出，華鎣市起伏度在70m以下的部分約為總面積的17%，這與實(shí)際情況有一定差異。考察這部分區(qū)域，大多為植被覆蓋的山麓地帶或城市建成區(qū)。由于ASTER GDEM利用遙感影像立體像對(duì)建立DEM，在一定程度上不能完全消除建筑物、水系、云層及植被覆蓋等的高度影響[16-17]，導(dǎo)致部分地區(qū)高程垂直精度異常。

圖2 S與Si差值曲線

圖3 華鎣山地形起伏度分級(jí)圖

5 結(jié)束語(yǔ)

地形起伏度的量測(cè)是重要的地形因子之一，利用ASTER GDEM計(jì)算研究區(qū)的起伏度，并將結(jié)果與常規(guī)測(cè)量結(jié)果進(jìn)行比較也是衡量DEM可用性和可靠性的重要方法。本文以川東平行嶺谷華鎣市為例，基于ARSTER GDEM數(shù)據(jù)，通過(guò)窗口遞增分析法、結(jié)合均值變點(diǎn)，對(duì)縣域尺度華鎣市地形起伏度進(jìn)行分析，結(jié)果表明，縣域尺度、較高分辨率條件下，平均地形起伏度曲線趨向增大；但最大起伏度可尋找到趨于穩(wěn)定的區(qū)間；均值變點(diǎn)探測(cè)的華鎣市的最佳柵格分析窗口邊長(zhǎng)約為27×27，對(duì)應(yīng)最佳統(tǒng)計(jì)單元為0.65km2，以這一分析單元得到的華鎣市中低山面積約為全區(qū)的40.85%；其結(jié)果與建立在常規(guī)測(cè)量基礎(chǔ)上的地理統(tǒng)計(jì)結(jié)果基本吻合，證明了分辨率為30m×30m的國(guó)際共享數(shù)據(jù)ASTER GDEM可以完成類似川東華鎣山脈及谷地的地貌起伏度計(jì)算。但研究區(qū)內(nèi)起伏度小于70m的區(qū)域面積偏大，除受ASTER GDEM自身精度的影響，分析使用的窗口分析法也存在較大的邊界效應(yīng)，在邊界區(qū)域得到的地形起伏度存在一定誤差。后續(xù)研究應(yīng)當(dāng)側(cè)重從地貌基礎(chǔ)認(rèn)知出發(fā)，考慮基于DEM的山地地貌實(shí)體表達(dá)方法[12，18]，以求更加客觀地反映區(qū)域地形起伏度。

參考文獻(xiàn)：

[1] 湯國(guó)安，李發(fā)源，劉學(xué)軍.數(shù)字高程模型教程[M].北京：科學(xué)出版社，2010.

[2] 周成虎，程維明，錢金凱.數(shù)字地貌遙感解析與制圖[M].北京：科學(xué)出版社，2009.

[3] 韓海輝，高婷，易歡，等.基于變點(diǎn)分析法提取地勢(shì)起伏度——以青藏高原為例[J].地理科學(xué)，2012，32(1)：101-104.

[4] 涂漢明，劉振東.中國(guó)地勢(shì)起伏度研究[J].測(cè)繪學(xué)報(bào)，1991，20(4)：311-319.

[5] 劉新華，楊勤科，湯國(guó)安.中國(guó)地形起伏度的提取及在水土流失定量評(píng)價(jià)中的應(yīng)用[J].水土保持通報(bào)，2001，21(1)：57-59+62.

[6] 劉愛(ài)利，湯國(guó)安.中國(guó)地貌基本形態(tài)DEM的自動(dòng)劃分研究[J].地球信息科學(xué)，2006，8(4)：8-14.

[7] 郎玲玲，程維明，朱啟疆，等.多尺度DEM提取地勢(shì)起伏度的對(duì)比分析——以福建低山丘陵區(qū)為例[J].地球信息科學(xué)，2007，9(6)：1-6+135-136.

[8] 封志明，唐焰，楊艷昭，等.中國(guó)地形起伏度及其與人口分布的相關(guān)性[J].地理學(xué)報(bào)，2007，62(10)：1073-1082.

[9] 周自翔，李晶，任志遠(yuǎn).基于GIS的關(guān)中-天水經(jīng)濟(jì)區(qū)地形起伏度與人口分布研究[J].地理科學(xué)，2012，32(8)：951-957.

[10] 谷建立，張海濤，陳家贏，等.基于DEM的縣域土地利用空間自相關(guān)格局分析[J].農(nóng)業(yè)工程學(xué)報(bào)，2012，28(23)：216-224.

[11] 郭芳芳，楊農(nóng)，孟暉，等.地形起伏度和坡度分析在區(qū)域滑坡災(zāi)害評(píng)價(jià)中的應(yīng)用[J].中國(guó)地質(zhì)，2008，35(1)：131-143.

[12] 朱紅春，陳楠，劉海英，等.自1∶10000比例尺DEM提取地形起伏度——以陜北黃土高原的實(shí)驗(yàn)為例[J].測(cè)繪科學(xué)，2005，30(4)：86-88.

[13] 康曉偉，馮鐘葵.ASTER GDEM數(shù)據(jù)介紹與程序讀取[J].遙感信息，2011，36(6)：69-72.

[14] 張學(xué)儒，官冬杰，牟鳳云，等.基于ASTER GDEM數(shù)據(jù)的青藏高原東部山區(qū)地形起伏度分析[J].地理與地理信息科學(xué)，2012，28(3)：11-14.

[15] 周申立.華鎣山滑坡研究[M].成都：西南交通大學(xué)出版社，2008.

[16] 郭笑怡，張洪巖，張正祥，等.ASTER-GDEM與SRTM3數(shù)據(jù)質(zhì)量精度對(duì)比分析[J].遙感技術(shù)與應(yīng)用，2011，26(3)：334-339.

[17] 朱圣軍.ASTER GDEM數(shù)據(jù)精度分析及其在石油勘探中的應(yīng)用[J].物探裝備，2011，21(1)：54-59.

[18] 羅明良，湯國(guó)安.地貌認(rèn)知及空間剖分的山頂點(diǎn)提取[J].測(cè)繪科學(xué)，2010，35(5)：126-127+253.