S=arctg
2f2x+fy
A=270°+arctg(fy/fx)-90°fx/|fx|g为格网间距
fx=(z8-z2)/2g fy=(z6-z4)/2gfx=(z7-z1+z8-z2+z9-z3)/6g
fy=(z3-z1+z6-z4+z9-z7)/6gfx=(z7-z1+2(z8-z2)+z9-z3)/8gfy=(z3-z1+2(z6-z4)+z9-z7)/8g
fx=(z7-z1+(z8-z2)+z9-z3)/(4+
2)gfy=(z3-z1+(z6-z4)+z9-z7)/(4+2)gfx=(z7-z1+z9-z3)/4gfy=(z3-z1+z9-z7)/4gfx=(z5-z2)/g,fy=(z5-z4)/g
二阶差分、矢量算法、不完全四次曲面
三阶不带权差分、线性回归平面、非限制二次曲面、限制二次曲面三阶反距离平方权差分、带权限制二次曲面、带权非限制二次曲面三阶反距离权差分
Frame差分
简单差分
规则格网 坡度坡向算法的比较分析 DEM
干 旱 区 地 理 27卷40 0
可用类似方法导出。所用符号参见表1,并且(x,y)为3×3窗口中心格网点5的坐标。
由数值微分理论知,二阶差分是以二阶精度逼近其真值的〔25〕,其逼近总误差可表示为如下〔25,26〕:dfx=dfy=
2
mm
fx(ζ,y)+f(ν,y)m mA=aM+bm/tgS,
2222
仿上述过程,并参考表1,其余算法的坡度、坡向中误差如表2。
表2 DEM坡度、坡向中误差
Tab.2 ErrorsincalculatingslopesandaspectswithDEM
算法
三阶差分
三阶不带权差分三阶反距离平方权差分
三阶反距离权差分
Frame差分
6
2
6
(3),
2g2
m
()()+,
2g2
+
计算模型误差M
系数a
g2/6g2/6g2/6g2/6g2/6g/2
数据误差m系数b
0.707/g0.408/g0.433/g0.414/g0.500/g2.000/g
在(3)式中,第一项是由连续数据的离散化结果和公式的截断所引起的误差,可归结为由数学模型不精
m
ζνν确引起。其中ζx、变化的量,ζx∈y、x、y是依f
(x,x+8),ζy∈(y,y+8),νx∈(x-g,x),νy∈(y-g,y)。由于这些变量与x,y的相关关系一般并不清楚,具体数值难以确定,通常是给其一个上界。设Mx,My分别是fm关于x和y关于的上界,则(3)可改写成(4):
dfx=dfy
,
62g2
=My+6Mx+
2
简单差分
通过对(8)式和表2,坡度坡向误差,对地形曲面离;
(2)DEM误差,包括DEM数据采样误差和数据舍入误差;
(3)DEM格网分辨率;
(4)
),在本文即DEM。
在(4)式中,Mx,My分别是fm关于x和y在某一格网点上的上界,它是按最坏情况估计误差限,而且一般比实际误差大得多,这种保守的误差估计不反映实际误差的积累,考虑到误差的分布特性,此类误差也具有一定的随机性,可视为服从某种分布的
27〕
随机变量〔。这对在在3×3移动窗口中的操作是适宜的,由于每点处的Mx,My的大小、符号并不相同,而且不可知,因此Mx,My具有随机性。设Mx,My中误差相等且同为M,并顾及DEM中误差设为m,则通过误差传播定律有Mx、My的中误差:
mfx=mfy=
2
2
4 数据独立DEM分析方法与结果
(8),在格网分辨率一定的情况下,考察式(7)、
坡度坡向精度fx与fy和的计算模型误差M和DEM误差m有关,坡度坡向的精度取决于这两个
因素中哪个起支配作用,而这与分析方法有关。比
20〕21〕24〕
如Skidmore〔、Florinsky〔、Chang和Tsai〔等在
实际DEM上进行算法分析,这时DEM误差起主要作用,
而且DEM误差所引起的坡度坡向误差比算法产生的误差大得多,故主要考察的是数据误差对
22〕23〕坡度坡向的影响。Hodgson〔、Jones〔等在数学曲
2
6
2
M+
2
,2g2
2
(5)
对坡度、坡向公式微分,并考虑到S=arctg
f
2
x
+f
2
y
和tgS=f
2
x
+fy有:
2
面上对坡度坡向的分析,DEM数据本身无误差,坡度坡向误差主要来源于数学模型。由于没有区分误差来源与性质,他们的分析结论相互矛盾。本文将在数据独立DEM上对坡度坡向算法进行分析。4.1 数学曲面设计
dS=
,
dA=,(6)
(1+tg2S)tgStg2S
顾及(5)式,则得坡度中误差mS、坡向中误mA差为:
ms= mA=
22
6
4
M+2cosS,
2g22
M+2/tgS,2g
2
考虑到实际地形表面的复杂性,利用数学曲面进行DEM坡度坡向算法分析时,数学曲面应最大限度的接近于局部实际地形,简单曲面并不能很好地反映算法所具有的精度和相关参数的确定,因此在本文的研究中,选择了凹向半椭球和高斯合成曲面,并对此曲面在给定的区域边界内按一定分辨率格网化建立数学曲面DEM如图2、图3所示。其中
(7)
(7)为二阶差分的坡度、坡向中误差,令a=
2
g/6,b=1/g,则(7)可表示为一般公式:
mS=aM+bmcosS,
22222
(8)
规则格网 坡度坡向算法的比较分析 DEM
3期 李天文等:基于DEM坡度坡向算法精度的分析研究 401
凹向半椭球定义为:
x2/A2+y2/B2+z2/
C2=1, (z<0)。(9)高斯合成曲面定义为:
z=A1-(
22-(m)2
)e-(+1)mn35--B0.2()-()-()emmn
(
22
)-()mn
(10)。
上两式中A、B、C为地势起伏参数,(10)中m、n为范围控制参数。两种曲面上的坡度坡向真值计算公
-Ce
-〔(
22
)+1〕-()mn
(10)由(1)、(2)式导出。式可通过曲面表达式(9)、
4.2 实验方法与结果
对上述数学曲面按一定分辨率离散化后建立相应的DEM,在此DEM上通过所选算法计算值和理论值的比较则可定量描述坡度坡向算法精度。精度指标采用中误差RMSE(RootMeanSquareError,RMSE)。
直接对数学曲面的离散化而建立的DEM,DEM中的格网点数据并没有误差,其中仅包含了
图2
凹向半椭球示意图,其中A=400,B=300,C=300
Fig.2 Sketchfigureofconcaveellipticsphere
图3 高斯合成曲面示意图,其中A=3,B=10,C=1/3,范围:-500≤X,Y≤500
Fig.3 SketchfigureofGausscomposedcamber
表3 高斯合成曲面DEM(无误差)坡度坡向中