范 林 峰，胡 瑞 林，張 小 艷，王 珊 珊，袁 延 西

基于GIS和DEM的水系三維分形計(jì)盒維數(shù)的計(jì)算

范 林 峰1，胡 瑞 林1，張 小 艷1，王 珊 珊1，袁 延 西2

（1.中國(guó)科學(xué)院工程地質(zhì)力學(xué)重點(diǎn)實(shí)驗(yàn)室，中國(guó)科學(xué)院地質(zhì)與地球物理研究所，北京 100029;2.西安科技大學(xué)地質(zhì)與環(huán)境學(xué)院，陜西 西安 710054）

實(shí)現(xiàn)具有三維立體特征的河網(wǎng)水系分形維數(shù)的準(zhǔn)確計(jì)算，是目前分形地貌學(xué)面臨的主要問(wèn)題之一。該文依據(jù)三維分形計(jì)盒維數(shù)的基本原理，提出基于DEM和ArcGIS的水系三維分形計(jì)盒維數(shù)的計(jì)算方法，并以湖北省恩施市境內(nèi)水系為例，驗(yàn)證了方法的可靠性，說(shuō)明三維計(jì)盒維數(shù)比二維計(jì)盒維數(shù)更能反映水系的空間分布特征和內(nèi)部結(jié)構(gòu)。

三維計(jì)盒維數(shù);分形;DEM;GIS

分形理論逐漸成為處理復(fù)雜性系統(tǒng)的有力工具，在地貌學(xué)領(lǐng)域尤其是在水系研究中應(yīng)用廣泛［1－10］。在眾多分形維數(shù)中，計(jì)盒維數(shù)由于簡(jiǎn)單實(shí)用，得到了廣泛應(yīng)用。何隆華［1］提出根據(jù)水系的計(jì)盒維數(shù)劃分流域地貌侵蝕發(fā)育階段的方法;鄒謹(jǐn)敞等［2］利用計(jì)盒維數(shù)法分析甘肅中部及鄰區(qū)水系的分形特征及其與滑坡分布的關(guān)系;盤(pán)小東等［3］計(jì)算了渾江斷裂帶水系的計(jì)盒維數(shù)，并分析了水系分形特征與構(gòu)造活動(dòng)性的關(guān)系;馬宗偉等［4］分析了長(zhǎng)江中下游河流的分形特性，并在此基礎(chǔ)上分析了不同水系特征下洪水的特點(diǎn);孫祝友等［5］應(yīng)用GIS技術(shù)提取了萊州灣東岸河流信息，并運(yùn)用計(jì)盒維數(shù)法對(duì)萊州灣東岸的黃水河和王河水系進(jìn)行分形分析;韓杰等［6］在滑坡泥石流災(zāi)害區(qū)劃中應(yīng)用水系分形維數(shù)的空間變化表征水系空間形態(tài)特征的差異，從而反映云南省滑坡泥石流災(zāi)害發(fā)生的程度。目前對(duì)水系分形維數(shù)的計(jì)算大多是在二維平面上進(jìn)行的，然而，天然溝谷中的水系具有三維立體特征，要得到關(guān)于水系空間形態(tài)和內(nèi)部結(jié)構(gòu)真正客觀(guān)、科學(xué)的研究結(jié)論，必須對(duì)三維空間中的水系分形特征進(jìn)行深入探討。本文根據(jù)分形理論，利用ArcGIS軟件的空間分析功能，計(jì)算水系的三維分形計(jì)盒維數(shù)。

1 利用DEM生成水系

河網(wǎng)水系分布圖是計(jì)算水系三維分形維數(shù)的必備資料之一，可以由數(shù)字高程模型（DEM）自動(dòng)提取。目前最常用的提取方法是坡面流模擬法，依據(jù)水總是沿斜坡坡度最大方向流動(dòng)的原理，首先確定DEM中每個(gè)柵格單元的水流方向，據(jù)此計(jì)算每個(gè)柵格單元的上游給水區(qū)，再選擇合適的集流閾值確定河網(wǎng)。由于該方法依據(jù)水文學(xué)匯流概念判別水流路徑，被認(rèn)為是一種較好的方法;但仍存在一定問(wèn)題，即DEM中不能存在局部凹陷和平坦的柵格單元，否則會(huì)使水流無(wú)法確定。因此，在提取河網(wǎng)水系之前，需要對(duì)DEM進(jìn)行預(yù)處理，包括洼地的確定、填充和平地的抬升，再利用經(jīng)過(guò)洼地填充處理的DEM源數(shù)據(jù)生成水系。上述操作過(guò)程可在ArcGIS軟件中實(shí)現(xiàn)，由于篇幅所限，在此不再贅述。

2 三維計(jì)盒維數(shù)的計(jì)算模型

目前研究中使用的計(jì)盒維數(shù)一般是二維計(jì)盒維數(shù)。計(jì)算時(shí)，首先用邊長(zhǎng)為r的正方形網(wǎng)格覆蓋二維分形圖形，得到與圖形相交的非空網(wǎng)格數(shù)Nr，逐步改變網(wǎng)格的邊長(zhǎng)r，得到一系列Nr值。二維計(jì)盒維數(shù)可用下式表示：

式中：Df為研究對(duì)象的計(jì)盒維數(shù)。

本文的研究對(duì)象是三維空間中的水系，不能用二維正方形網(wǎng)格覆蓋，而應(yīng)采用三維小立方體對(duì)空間水系進(jìn)行覆蓋，計(jì)算被水系占據(jù)的小立方體的數(shù)目（圖1）;不斷變換小立方體的邊長(zhǎng)，得到一系列被水系占據(jù)的小立方體數(shù)目，進(jìn)而求得水系的三維計(jì)盒維數(shù)。這里，核心的問(wèn)題是如何計(jì)算被水系占據(jù)的小立方體的數(shù)目。

要解決這個(gè)問(wèn)題，可首先在XOY二維平面上用不同邊長(zhǎng)的正方形網(wǎng)格覆蓋二維水系分布圖，確定非空網(wǎng)格;再由水系高程，確定XOY平面內(nèi)非空網(wǎng)格對(duì)應(yīng)的水系在Z軸上占據(jù)的小立方體數(shù)目（圖2）。

圖1 立方體覆蓋示意Fig.1 Schematic diagram of cube covering

圖2 立方體個(gè)數(shù)計(jì)算Fig.2 Counting of covering cubes

XOY平面內(nèi)每個(gè)非空網(wǎng)格對(duì)應(yīng)的水系在Z軸上占據(jù)的小立方體數(shù)目可由下式計(jì)算：

式中：Ni，j為用邊長(zhǎng)ri的正方形網(wǎng)格覆蓋時(shí)，第j個(gè)非空網(wǎng)格對(duì)應(yīng)的水系在Z軸上占據(jù)的立方體數(shù)目;ri為小正方形網(wǎng)格（即小立方體）的邊長(zhǎng);max（hj）、min（hj）分別為第j個(gè)非空網(wǎng)格對(duì)應(yīng)的水系的最大、最小高程;INT為向下取整運(yùn)算。

以邊長(zhǎng)為ri的小立方體對(duì)水系進(jìn)行覆蓋時(shí)，被水系占據(jù)的小立方體的總數(shù)Ni為：

3 三維計(jì)盒維數(shù)的計(jì)算過(guò)程

基于上述計(jì)算模型，可以根據(jù)研究區(qū)的DEM和水系分布圖，利用ArcGIS的空間分析功能，計(jì)算水系的三維計(jì)盒維數(shù)。計(jì)算步驟如下：1）在XOY平面內(nèi)生成邊長(zhǎng)為ri的Fishnet網(wǎng)格，并將其由線(xiàn)狀要素轉(zhuǎn)化為面狀要素，生成面域網(wǎng)格;2）將面域網(wǎng)格與二維水系分布圖求交，確定二維平面上被水系占據(jù)的非空網(wǎng)格;3）利用ArcGIS的Zonal Statistics功能，求出每個(gè)非空網(wǎng)格內(nèi)水系的最大高程max（hj）和最小高程 min（hj）;4）由式（2）和式（3）計(jì)算用邊長(zhǎng)為ri的小立方體對(duì)三維水系進(jìn)行覆蓋時(shí)，被水系占據(jù)的小立方體數(shù)目Ni;5）改變Fishnet網(wǎng)格邊長(zhǎng)ri，重復(fù)上述步驟，得到一系列Ni值;6）將Ni和ri值繪入雙對(duì)數(shù)坐標(biāo)系，由式（4）求得研究區(qū)水系的三維計(jì)盒維數(shù)。

4 計(jì)算實(shí)例

本文以湖北省恩施市境內(nèi)的清江水系（圖3、圖4）為例，按照上述模型計(jì)算出水系的二維和三維分形計(jì)盒維數(shù)（圖5、表1）分別為1.18和1.26。

圖3 恩施二維水系分布Fig.3 2D river system in Enshi

圖4 恩施三維水系分布 Fig.4 3D river system in Enshi

圖5 水系二維和三維計(jì)盒維數(shù)計(jì)算結(jié)果 Fig.5 Ln-ln plot of Ni and ri for 3D and 2D fractal dimension

表1 水系二維和三維計(jì)盒維數(shù)計(jì)算Table 1 Calculation of 3D and 2D fractal dimension of river system

5 結(jié)果分析與討論

（1）歐式幾何中，直線(xiàn)的維數(shù)是1，三維實(shí)體的維數(shù)是3。三維空間中的水系可以看做是空間曲線(xiàn)段的集合，其計(jì)盒維數(shù)應(yīng)該在1～3之間。本文的計(jì)算結(jié)果與上述分析一致，且與文獻(xiàn)［11］的結(jié)果相近，因此可以認(rèn)為本文提出的方法是可行的。

（2）從表1和圖5可以看出，三維空間內(nèi)lnNi值大于二維平面的lnNi值，且二者的差值隨著ri的增加有逐漸減小的趨勢(shì)。當(dāng)ri≥3 200 m時(shí)，二者的差值為0，即覆蓋水系所需的小立方體數(shù)目與正方形網(wǎng)格數(shù)目相等。這是由于當(dāng)非空正方形網(wǎng)格內(nèi)三維水系的最大高程和最小高程接近到一定程度時(shí)，在Z軸上僅需要1個(gè)立方體即可將此網(wǎng)格對(duì)應(yīng)的三維水系覆蓋（式（2））。由此可以看出，對(duì)于地形起伏較小的地區(qū)，水系的三維計(jì)盒維數(shù)與二維計(jì)盒維數(shù)將十分接近，三維水系分形維數(shù)更適用于地形起伏較大的地區(qū)。

（3）利用本文方法計(jì)算水系的三維分形維數(shù)受原始數(shù)據(jù)的精度影響較大。由于本文采用的DEM分辨率是90 m，因此，本文將小立方體的最小尺寸設(shè)置為100 m。當(dāng)DEM的數(shù)據(jù)精度提高時(shí)，應(yīng)該采用更小尺寸的立方體對(duì)研究區(qū)的水系進(jìn)行覆蓋，從而獲得更為精確可靠的計(jì)算結(jié)果。

（4）對(duì)水系的分形研究大多采用計(jì)盒維數(shù)法，主要是因?yàn)樵摲椒ê?jiǎn)單方便。但是，該方法僅考慮被水系占據(jù)的正方形網(wǎng)格和小立方體數(shù)目，而沒(méi)有考慮每個(gè)正方形網(wǎng)格和小立方體中水系的條數(shù)或長(zhǎng)度，由此計(jì)算出的分形維數(shù)不能完全反映空間水系的分形結(jié)構(gòu)特征。信息維由于考慮了每個(gè)正方形網(wǎng)格和小立方體中水系的條數(shù)，因而能彌補(bǔ)這方面的不足。此外，還有一些學(xué)者認(rèn)為空間水系具有多重分形的特征。多重分形分析比單分形分析（如計(jì)盒維數(shù)、信息維）更能真實(shí)地反映空間水系的復(fù)雜性和本質(zhì)特征，所獲得的標(biāo)度指數(shù)和多重分維可以從整體上反映出水系的多重分形結(jié)構(gòu)特征，而奇異指數(shù)和分維譜函數(shù)則能反映出水系多重分形的局部特征。因此，需進(jìn)一步研究空間三維水系的多重分形特征。

［1］何隆華.水系的分形維數(shù)及其含義［J］.地理科學(xué)，1996，12（2）：124－128.

［2］鄒謹(jǐn)敞，邵順妹.甘肅中部及鄰區(qū)水系分形研究結(jié)果與滑坡分布的關(guān)系［J］.地震研究，1994，17（4）：383－388.

［3］盤(pán)曉東，康力.渾江斷裂帶及水系的分形特征和構(gòu)造活動(dòng)性研究［J］.東北地震研究，2002，18（1）：7－15.

［4］馬宗偉，許有鵬，李嘉峻.河流形態(tài)的分維及與洪水關(guān)系的探

討——以長(zhǎng)江中下游為例［J］.水科學(xué)進(jìn)展，2005，16（4）：530－534.［5］孫祝友，杜國(guó)云，李德一，等.基于GIS技術(shù)的萊州灣東岸河流

分形研究［J］.測(cè)繪科學(xué)，2007，32（3）：120－121.

［6］韓杰，陸桂華，李海濤.水系分維在滑坡泥石流災(zāi)害區(qū)劃中的應(yīng)用［J］.自然災(zāi)害學(xué)報(bào)，2009，18（4）：63－71.

［7］朱俊江，詹文歡，唐誠(chéng)，等.水系的分形特征與紅河斷裂帶活動(dòng)性關(guān)系初探［J］.華南地震，2002，22（1）：1－7.

［8］連建發(fā)，劉傳正，高尚嶸.三峽庫(kù)區(qū)溝谷水系的分形幾何特征研究［J］.勘察科學(xué)技術(shù)，2006（6）：25－28.

［9］馬宗偉，許有鵬，鐘善錦.水系分形特征對(duì)流域徑流特性的影響——以贛江中上游流域?yàn)槔跩］.長(zhǎng)江流域資源與環(huán)境，2009，18（2）：163－169.

［10］杜尚海，蘇小四，朱琳.松花江流域地表水系分形維及其影響因素分析［J］.水文，2009，29（5）：30－35.

［11］朱曉華.地理空間信息的分形與分維［M］.北京：測(cè)繪出版社，

2007.77－84.

Calculation of 3D Box Dimension of River System Based on GIS and DEM

FAN Lin－feng1，HU Rui－lin1，ZHANG Xiao－yan1，WANG Shan－shan1，YUAN Yan－xi2
（1.KeyLaboratoryofEngineeringGeomechanics，InstituteofGeologyandGeophysics，CAS，Beijing100029;2.CollegeofGeologyandEnvironment，Xi′anUniversityofScienceandTechnology，Xi′an710054，China）

Calculation of the fractal dimension of 3D river system is one of the major problems in fractal geomorphology.Based on the principle of box counting method，a method to calculate the 3D fractal dimension of river system based on ArcGIS and DEM is proposed.Enshi was used as an example to verify the reliability of this method.The results indicate that 3D box dimension can reflect more information about the spatial distribution and inner structure of river system than 2D box dimension.

3D box dimension;fractal theory;DEM;GIS

P 64;P208

A

1672－0504（2012）06－0028－03

2011－12- 11;

2012－03－29

國(guó)家科技支撐計(jì)劃課題“重大滑坡、泥石流災(zāi)害綜合風(fēng)險(xiǎn)評(píng)價(jià)技術(shù)”（2008BAK50B04－3）;中國(guó)科學(xué)院知識(shí)創(chuàng)新工程項(xiàng)目“地震地質(zhì)災(zāi)害的內(nèi)外動(dòng)力耦合作用機(jī)制”（KZCX2－YW－Q03－2）;湖北省交通運(yùn)輸廳科技項(xiàng)目“宜巴高速公路巖溶隱患辨識(shí)與風(fēng)險(xiǎn)防范研究”

范林峰（1988－），男，碩士研究生，主要從事地質(zhì)災(zāi)害風(fēng)險(xiǎn)評(píng)價(jià)方面的研究。E－mail：fanlinfeng1988＠163.com