面向DEM地形复杂度分析的分形方法研究

2006年1月GeographyandGeo-InformationScienceJanuary2006

面向DEM地形复杂度分析的分形方法研究

周　侗1,2,龙　毅1,汤国安1,杨　昕1

(1.南京师范大学地理信息科学江苏省重点实验室,江苏南京210097;2.南通大学地理科学学院,江苏南通226007)

摘要:在基于栅格DEM数据的基础上,应用分形几何学方法,采用元分维模型理论,提出一种描述DEM地形复杂度的分形分析方法,并得出分析指标———地形分维指数(TFI)。实验证明,该指标可以有效描述栅格DEM数据反映的地形变化特征,分析窗口由小到大的规律性变化反映出地貌从微观到宏观的变化情况,因此其既可以作为对栅格DEM所描述地形进行评价的坡面因子,又可以作为基于DEM的地貌类型区自动划分的参照依据。关键词:DEM;元分维;地形;复杂度;地形分维指数

中图分类号:P208　　文献标识码:A　　文章编号:1672-0504(2006)01-0026-05

0　引言

20世纪50年代末提出数字高程模型(DEM)的

学研究地形描述,利用VC6610开发了相关的分析

软件,经过大量DEM数据验证,证明本实验采用元分维分析方法可以有效地通过对栅格DEM数据的分析来描述自然地形的复杂度状况,可以作为进一步自动划分地貌类型区的依据。概念后,DEM逐渐取代等高线地形图,成为地理数据库中的核心数据。它是进行地形分析的基础,广泛应用于测绘、遥感、地质、制图综合、地理形态、资源、环境等领域[1]。而基于DEM数据的地形特征分析问题一直是学术界的研究热点。描述地形特征的参数很多,不同学科和领域对于地形属性的分类也有差别[2],不同的地形要素反映了地形不同方面的特征,对应于不同的应用范围。目前对地形描述的坡面因子主要有坡度、坡向、坡长、坡度变率、坡向变率、平面曲率、剖面曲率、地形复杂度指数[3]等微观因子以及地形粗糙度、地形起伏度、高程变异系数、地表切割深度等宏观因子[4]。微观因子的作用范围大都集中在单一方向,或仅能评价地形细节的变化,无法评价地形整体状况;宏观地形因子能描述一个地区整体的地形状况,却无法揭示分析区域局部的地形差异。实际研究中多将若干因子进行组合[5],或者将地形变化分解到多个二维空间。一个能够综合反应地形信息状况的指标,对于栅格DEM数据的地形分析

1　实验方法分维估值方法种类繁多,但不同的分数维适用于不同的领域,因此在应用分形分析方法前应首先选取合适的分维估值方法。111　分维估值方法

经过大量研究证明,尺度—表面积法能够较好描述地表曲折度,能够综合地表局部的粗糙程度及起伏程度[10],因此本实验采用尺度—表面积法作为分维估值方法。其原理如下:对于所考虑的地形表面S={(x,y,z(x,y))|(x,y∈G)},若用正方形网格覆盖G,即将G分成n(r)×n(r)个边长为r的小正方形,则第i个小正方形Gi对应的S的表面积在r较小时,可由折平面面积Si来近似,从而S的面积A可写成:

A∝∑=n2(r)r2

j=12n(r)

(1)

及地貌状况的自动研究意义重大。

由量纲分析方法及线性回归分析原理[6]得到:

自曼德勃罗特(Mandelbrot)创立分形几何学以MMM

M∑[(logri)(logA(ri))-(∑logri)(∑logA(ri))]

i=1i=1来,分形方法在地学研究中得到广泛应用,如河流、D=2-i=1(2)MM

22[6]()()M∑logri-∑logri海岸线的分维计算、城镇结构分析、地图综合等,i=1i=1

分维扩展方法以及分形应用理论得到了较好的发式(2)即利用表面积—尺度关系求曲面分维值展,为分形方法应用于地形研究奠定了基础。国内一些研究成果也表明[7-9],分形理论可以有效地应用于地貌特征及DEM的研究。本文基于分形几何

的公式,D为曲面的分维。式中需要由不同的ri求出相应的面积A(ri)(i=1,2,…,M),由式(1)可知,只要对固定的码尺r求出相应的Si即可(j=1,2,

　　收稿日期:2005-09-29;　修订日期:2005-10-27　　基金项目:国家自然科学基金项目(40271089);南京师范大学重点科研基金(2004105XGQ2B55)　　作者简介:周侗(1978-),男,硕士研究生,从事GIS和DEM方面的研究。E-mail:zhoutong@ntu.edu.cn

第1期　　　　　　　　　　周　侗等:面向DEM地形复杂度分析的分形方法研究第27页

…,n2(r))。实验中Si采用三角平面法求得。如图1所示,对于给定的r,可以把Gi的4个顶点所对应的高程点作为空间三角形的顶点,形成两个空间三角形ΔAjBjCj和ΔDjBjCj,并以这两个三角形的面积之和作为Si,由于4个顶点的坐标可以由其行列号和格网单元尺寸计算得到,因而根据三角形面积计算公式(3)可求得空间三角形的面积,累加即可得到Gi及Si、Ai。

SΔ=

s(s-a)(s-b)(s-c)

(3)

指导方法进行分析验证。

在地理现象中,地理实体内的相邻部位总是在形态、结构和动力学成因等方面具有较大的相似性,同样地形对象间互相作用的过程也和距离有直接关系,元分维模型的思想基础正是基于这点,即对于分形对象的每一个空间单元,在相邻区域内和其他单元距离越近,则空间性质越相似。为了能够较好地反映这种空间特征,元分维理论采用“滑动窗

[11,12]口”的方法对DEM数据进行分析(图2):1)在滑动窗口覆盖的区域范围内,通过计算得到不同尺度下ri及Ai的值;2)通过线性回归方法得到分维值,将该值作为中心栅格点的分维值,命名为“元分维值”,这一参数包含中心点一定邻域范围内的地形情况,不像传统研究方法中仅对相邻栅格进行分析,因而一定程度上保证了分析结果的连续性;3)将滑动窗口按照一定方式(如图2c、图2d)遍历研究的

式中:a、b、c为三角形的三边,s=(a+b+c)/2。

Fig.1　Principleoftriangle-plane

图1　三角平面法原理

DEM数据,可得到与原数据一一对应的DTM数据,

112　实验模型

自然界的地表是连续的,栅格DEM是通过离散数据来实现对真实地表的模拟,如果仅用两相邻栅格的关系来描述真实地表地形的变化,如坡度、坡向的分析,一定程度上进一步夸大了这种离散性,而忽略了地形的连续性。为能够利用DEM数据描述地形在空间上的性质变化,既反映局部差异,又不割断目标内部的连续性,本实验采用元分维模型[10]作为

属性值为元分维值。由于实验得到的DTM数据与DEM源数据在栅格上一一对应,因而可以实现DEM源数据和DTM目标数据的叠加分析。

这种分析方法还可以通过对滑动窗口的大小进行调节,控制中心点的影响范围,得到的参数可以反映地形特征从宏观到微观的变化情况。这种基于尺度得到的元分维值,可以反映地形的复杂度状况,本文将此参数定义为地形分维指数(TerrainFractalDimensionIndex,TFI)。

Fig.2　Sizeandmodeoftheslidingwindow

图2　滑动窗口的大小及滑动方式

113　尺度的选取

这要求在实验时,必须对实验数据进行简单的分析,如果仅考虑到尺度—表面积法的估值过程,为

了能实现对一个单元的分维值估算,一般选取3×3栅格作为可能的最小滑动窗口。实验中,根据实际的分析数据可以考虑适当放大;最大的滑动窗口面积不应超过所研究对象中最小地貌单元的面积,否则将会出现地形分维指数反映不出地形变化的情况,失去了此研究方法的意义,因而应用此方法前必须选取合理的分析尺度。

地球信息科学理论中,空间尺度是描述地理空间的必要条件,是处理分析地理数据的前提,是空间分析与决策的基础[13]。尺度的选取也是实验的一个关键问题,滑动窗口选择过小,则割裂了目标内部的连续性,未能消除数据的离散特征;如果滑动窗口尺度过大,则类似于传统的对整个目标的分维估值研究,削弱了目标内部的局部差异,分析结果同样无效,所以必须选择合理的分析尺度。

第28页地理与地理信息科学　　　　　　　　　　　　　　　　第22卷

2　算法验证

211　模拟数据

本实验用大量的DEM数据验证地形分维指数

的有效性,前期利用理想地形的DEM数据进行算法有效性分析,其中选取有代表性的3组行数为40、列数为30的DEM数据来说明(图3)。图3a数据在坡度和坡向上均呈现出稳定的趋势;图3b数据是在坡度方向上有变化的DEM数据,其中有7处高程呈现转折;图3c数据为坡向上有变化的DEM数据,其中有5处高程转折。212　结果分析

对图3a～图3c进行基于3×3滑动窗口的分析,得到地形分维指数,将该值作为属性值建立图3d

～图3f的DTM数据。通过实验得到如下结论:1)地形分维指数的大小能够反映地形变化的程度,地表越复杂,指数越大,反之越小。将分析结果与原始数据对比分析可以看出,地形分维指数反映了原始地形的复杂度,当此区域DEM数据中地形变化剧烈时,得到的指数值较大,而地表趋于平缓时,指数变小,当地形未发生变化时,指数值为定值2。2)地形分维指数能够反映地形在各个方向上的变化,可以弥补平面曲率、剖面曲率等坡形因子受方向限制的不足。为了充分说明指数的有效性,将原始DEM数据与地形分维指数在剖面图中叠加,同时在便于观察的角度,将DEM数据的高程值和地形分维指数利用公式(4)进行极差标准化处理。

Sn(t)=(S(t)-Smin)/(Smax-Smin)

(4)

Fig.3　DEMsandtheDTMsofterrainFractal-DimensionIndexbasedontheDEMs

图3　一组DEM数据及3×3滑动窗口计算得到的地形分维指数DTM

　　结果如图4所示,其中图4a为图3b、图3e数据

沿坡度方向的剖面图,Y轴为经过极差标准化的高程数据和地形分维指数值。图4b为图3c、图3f数据沿坡向方向的剖面图,从图上反映出指数变化与

原始数据的高程变化吻合较好,证明该指数能够较好地反映原始DEM数据的地形变化情况。

Fig.4　TheelevationoftheoriginalDEMsdataandtheTFIintheTFI-DTMsbasedontheDEMs

图4　原始DEM数据的高程和对应TFI-DTM的地形分维指数剖面对照

3　应用分析

基于上述分析,实验选取渭北黄土高原陕西省

耀县1∶5万DEM数据。其地理位置为东经108°34′～109°06′,北纬34°50′～35°20′,东西宽约43km,南北

长约56km,地势西北高,东南低。实验结论如下:

第1期　　　　　　　　　　周　侗等:面向DEM地形复杂度分析的分形方法研究第29页

(1)地形分维指数DTM能够反映出原始DEM

地形变化情况,而不是简单反映原始地形高程分布情况(图5),即使原始数据的负地形,由于该地区地

形复杂,同样在相对应的地区地形分维指数较大,在DTM上呈现正地形,将TFI-DTM和原始数据等高线套合,在等高线密的地方,指数较大,等高线稀疏的地方,指数较小,吻合情况良好。

(2)滑动窗口的大小反映了该区域地形特征从宏观到微观的变化程度。滑动窗口的大小反映了该栅格邻域作用的范围,滑动窗口较大时,地形分维指数描述的是中心栅格较大邻域范围内的地形复杂度情况,表现为该区域的整体特征,趋向于宏观的态势;滑动窗口较小时,指数描述的则是小区域范围内的地貌破碎程度,表现为局部特征,趋向于微观的态势(图6)。

(a)耀县局部1∶5万DEM数据　　(b)TFI-DTM与等高线套合

图5　DEM数据对应TFI-DTM山体阴影图与等高线套合

Fig.5　TheDMEsandtheoverlapping

ofhill-shadesandcontoursoftheTFI-DTMs

(a)耀县局部1∶5万DEM数据　　　　　(b)200m尺度的TFI-DTM　　　　　(c)500m尺度的TFI-DTMFig.6　TheDEMsofYaoCountyof1∶50000scaleandthecorrespondingTFI-DTMsbasedondifferentscale

图6　耀县1∶5万DEM数据及对应的不同尺度的TFI-DTM

　　(3)按照传统的地貌类型划分方法,耀县从北向南大致分为土石山地峡谷、黄土高原丘陵沟壑和川

原平地3个地貌单元(图7a)。对3个区域分别进行地形分维指数计算,结果如图7b所示,可以看出在同一分析尺度下,北部山地区由于地形变化复杂,其地形分维指数明显大于中部残塬沟壑区,而南部台塬区由于地势平坦,指数最小。由于尺度的变化,分析单元地形分维指数反映的分析范围不同,造成指数的大小发生变化,使同一地区的地形分维指数平均值也呈规律性变化。

个滑动窗口指数的获得,即通过不同尺度ri得到相应的Ai值,然后通过线性回归的方法得到地形分维指数,因此每一个滑动窗口内都包含了一种分析尺度ri和Ai的对应关系图,可以作为进一步研究的分析工具;二是通过对滑动窗口尺寸的调节,可以控制中心栅格的影响区域,从而可以表述从微观到宏观不同的含义,并针对不同的研究用途加以利用。实验证明,该方法可以有效地描述地表的地形复杂度,即可以作为DEM数据的坡面研究因子来研究地表地形状态。同时,得到研究区域的地形分维指数后,可以将该指标结合单元的行列坐标计算单元间的空间距离,计算结果作为距离系数进行空间聚类分析,实现研究目标区域的初分类,从而为进一步自动划分地貌类型区提供依据。

　　　(a)地貌类型　　　(b)各地貌类型区地形分维指数比较

参考文献:

[1]　李志林,朱庆.数字高程模型[M].武汉:武汉大学出版社,

2001.16,260.

[2]　史文中.空间数据与空间分析不确定性原理[M].北京:科学出

图7　耀县地貌类型及对应的不同尺度的TFI-DTM

Fig.7　LandformtypeareasofYaoCounty

andtheTFI-DTMsbasedondifferentscale

4　结语

地形分维指数的研究方法借鉴了分形几何学的思想,并参照元分维模型理论,同时又是一种基于尺度的研究方法,这里的尺度有两层含义:一是指在每

版社,2005.171-173.

[3]　王雷,汤国安,刘学军,等.DEM地形复杂度指数及提取方法研

究[J].水土保持通报,2004,24(4):55-58.

[4]　汤国安,刘学军,闾国年.数字高程模型及地学分析的原理与方

法[M].北京:科学出版社,2005.208-233.

[5]　张婷,汤国安,王春,等.黄土高原地形因子间关联性的神经网

第30页

络分析[J].地球信息科学,2004,6(4):45-50.

地理与地理信息科学　　　　　　　　　　　　　　　　第22卷

[10]　龙毅.扩展分维模型在地图目标空间信息描述中的应用研究

[D].武汉大学,2002.118-120.

[11]　CHENGQ.Theglidingboxmethodsformulti-fractalmodeling

[J].Computers&Geosciences,1999,25:1073-1079.

[12]　PROKOPHA.Multi-fractalandslidingwindowcorrelationdi2

mensionanalysisofsedimentarytimeseries[J].Computers&Geosciences,1999,25:1009-1021.

[13]　李霖,应申.空间尺度基础性问题研究[J].武汉大学学报(信

[6]　王桥,毋河海.地图信息的分形描述与自动综合研究[M].武

汉:武汉测绘科技大学出版社,1998.

[7]　艾南山,陈嵘,李后强.走向分形地貌学[J].地理学与国土研

究,1999,15(1):92-96.

[8]　朱永清,李占斌,鲁克新,等.地貌形态特征分形信息维数与像

元尺度关系研究[J].水利学报,2005,36(3):333-338.

[9]　王光霞,崔凯.基于分形的DEM精度评估[J].测绘学院学报,

2005,22(2):107-109.

息科学版),2005,30(3):199-203.

AFractalMethodtoDescribetheTerrainComplexityReflectedbytheRasterDEM

ZHOUTong1,2,LONGYi1,TANGGuo-an1,YANGXin1

(1.JiangsuProvincialKeyLaboratoryofGeographicalInformationScience,NanjingNormalUniversity,Nanjing210097;

2.InstituteofGeographicalSciences,NantongUniversity,Nantong226007,China)

Abstract:Anewanalysismethodisputforwardinthispaper,whichisguidedbythemeta-fractaltheoryandcanbeappliedinde2scribingtheterraincomplexityoftherasterDEMs.Accordingtothetheory,ananalyticalguideline,whichiscalledterrainfractaldi2mensionindex(TFI),isdefined.Inthisexperimentation,agroupofDEMsofhighresolutionandhighprecisionwhosemapscaleis1∶10000and1∶50000inYaoCounty,ShanxiProvincearealsoselectedasthebasicdata.Theresultshowsthatthismethodcanef2fectivelyevaluatetheterraincomplexityoftheDEMsdata.Asanexperimentalresult,theTFIDTMsarehomologoustothesourceda2taofDEMspixelbypixel.Becausethescaleisadoptedintheexperimentation,itcanreflectthevarietyofterrainfromthemicroviewtothemacroviewthroughtheregularchangeoftheanalysisscale.Soitmaybeconsideredasatopographicalfactor,whichcanbeap2pliedforautomaticdivisionoflandformregion.

Keywords:DEM;meta-fractal;terrain;complexity;terrainfractaldimensionindex

(上接第25页)

MappingandQuantitativeAnalysisontheCorrelationbetweenLandslideandGeometric

AlignmentofStrikesofSlopeandStratumBasedonGIS

ZENGZhong-ping1,FUXiao-lin2,LIUXue-mei1,XUEChong-sheng3

(1.GraduateSchoolofChinaUniversityofGeosciences,Wuhan430074;2.ThreeGorgesGeologicalHazardMonitoringCenter,ChinaInstituteofGeologyandEnvironmentMonitoring,Yichang443002;

3.CollegeofEarthScience,ChinaUniversityofGeosciences,Wuhan430074,China)

Abstract:Thegeometricalignmentorconformitybetweenstrikesofslopeandgeologicbeddingplanessuchassoftorweakstratumsurfacesisoneofthecausativefactorsofhillslopelandslide.Topographic/bedding-planeintersectionangle(TOBIA)indexproducesspatiallydistributedfieldsofgeometricalignmentbetweentopographyandtheorientationofgeologicbeddingplanes.ComputationanddigitalmappingofTOBIAindexrequiresthederivationoffourspatiallydistributedvariables:topographicslope,slopeaspect,beddingdip,anddipazimuth.TOBIAindexisusefulforpredictionofpotentialslopefailureandastudyareaofThreeGorgesreservoirischo2sentoverifythecorrelationbetweenTOBIAindexanddipslopelandslides.Theapplicationdemonstratesthatthemethodherecanprovideanefficientmeansforestimatingslopetype.Basedonthepreparedlandslidesdatabase,twolandslidefactors,slopeandTOBIAindexhavebeencombinedintothebinarylogisticregressionmodel.Theregressionfunctionhasclassifiedsuccessfullyincreasedfrom55%to7115%oftheoverallslopepopulationafterindexisadded.TheresultsshowthatTOBIAindexcanbeaquantitativecausativefactordatasourceforlandslidehazardassessment.

Keywords:TOBIAindex;slopetype;ThreeGorgesreservoirarea;landslidehazardassessment