利用测地坐标计算椭球面上凸多边形面积的算法研究

(s
x1
, s y1 , s x2 , s y2 , s x3 , s y3 。过 P1 , P2 , P3 三点分别作测地格网线（图 2）, 若把两内角均为直角的
) (
) (
)
三个曲边“梯形” 各分成曲边“矩形” 和“直角三角形” 两部分，而在“直角三角形” 面积的计 算中，将其看成所相应的曲边“矩形” 的一半左右。因这两个曲边“直角三角形” 仅有一条边不等 长，故其面积之比应为长度归化因子之比。则由式（5）可得出该椭球面三角形的面积为
目（03-04-02） 作者简介：施一民(1942), 男 ,浙江宁波人, 教授，博士生导师。.E-mail：yimshi@ 1
中国科技论文在线 ___________________________________________________________________________
(
)(
)
2 2 2 3 2 − s3 y1 + s y1 s y 2 − 2 s y1 s y 2 + s y1 s y 2 + s y 2 − 2 s y1 s y 2 s x2 − s x1 x3
)( )(s
− s x2
) )]
+{
=
3 3 s3 y1 s x3 − s x2 + s y 2 s x1 − s x3 + s y3 s x2 − s x1
(
)(
)
(
)(
)
[ (
)
(
)
(
)]
+
1 s y s y s y + s y3 s x1 − s x3 + s y2 s y1 s y2 + s y1 s x2 − s x1 + s y3 s y2 s y3 + s y2 s x3 − s x2 8R02 1 3 1
[
(
)(
)
(
)(
)
(
[7] 由式（7）可知,利用三顶点的测地坐标所算得的椭球面三角形面积可分为主项 F0 以及含 小项 δF 两部分。主项 F0 的表示式完全等同于用平面坐标计算面积的常用公式。
Hale Waihona Puke 摘 要：本文推导出用三顶点的测地坐标计算地球椭球面上三角形面积的公式，公式表明，其 主项的表示式与按平面坐标求面积的计算式完全一致，而附加项的表示式亦有规律可循，因此 该公式的适用范围可由椭球面三角形推广至椭球面上任意凸多边形。与高斯平面上计算的面积 相比，由于它不再蒙受投影变形的影响，更接近于实际面积,这就为至今难以实施的椭球面上面 积计算开拓了一个新的途径, 有利于作出更客观的 GIS 空间量度与分析。实际数据的验算充分 证实了该算法的正确性和有效性。 关键词：测地坐标 椭球面三角形 凸多边形 面积公式
由（3） 、 （4）两式，得出
(4)
A = ∫ EG ds y ds x = ∫ nds y ds x = ∫
D D
s x2
s x1
∫
s y2
s y1
⎛ s2 ⎜1 − y ⎜ 2R 2 0 ⎝
⎞ ⎟ds y ds x ⎟ ⎠
=∫
s x2
s x1
3 ⎡ s3 y 2 − s y1 ⎢ s y2 − s y1 − 6 R02 ⎢ ⎣
F=
1 m 1 s yi + s yi +1 s xi − s xi +1 + ∑ 2 i =1 24 R02
(
)(
)
∑ s 3yi s xi −2 − s xi −1 +
i =1
m
(
)
1 8R02
∑s
i =1
m
yi
s yi −1 s yi + s yi −1 s xi − s xi −1
（8）
(
)(
(
)
(
)
(
)
( )( ) )( )]
1 的 R02
24 R
2 0
1 1 1 s y1 + s y2 s x1 − s x2 + s y2 + s y3 s x2 − s x3 + s y3 + s y1 s x3 − s x1 2 2 2 1 3 3 s3 + y s x − s x2 + s y 2 s x1 − s x3 + s y3 s x2 − s x1 24 R02 1 3
)
P1 s x 1 , s y1
(
)
(B0 , L 0 )
P0
图 2：测地坐标系中的椭球面三角形 Fig.2 Ellipsoidal triangle in geodesic coordinate system
设椭球面上由三条大地线所围成的三角形 P 1 P2 P 3 的 三 顶 点 的 测 地 坐 标 分 别 为
[1]
1 利用测地坐标计算椭球面上测地格网曲边“矩形”的面积
按一定间隔的一族大地线及一族测地平行线构成了地球椭球面上的正交参数曲线格网（测地格 网）
[2 ]
,借此可度量椭球面上由两对正交的坐标曲线所围成的四内角均为直角的四边形面积。在图 1
所示的测地格网的一个微分网格 PQTL 中，当 Q，L 点无限接近于 P 点，则有
− s y2 s x2 − s x1 +
(
)
s3 y2 6R
2 0
(s
x2
⎡ − s x1 − ⎢ s y1 − s y2 s x2 − s x1 ⎢ ⎣
) (
)(
)
⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦
− s y2 s x3 − s x2 +
(
)
s
3 y2
6 R02
(s
x3
⎡ − s x2 − ⎢ s y3 − s y2 s x3 − s x2 ⎢ ⎣
→ ⎛ → ⎜ ∂r ∂r × A = ∫∫ dσ = ∫∫ ⎜ ⎜ ∂s y ∂s x D ⎝
⎞ ⎟ ds y ds x ⎟ ⎟ ⎠
（3）
---------------------------------------------------------------------------------------------------------------基金项目：国家自然科学基金资助项目(40471114) ,地球空间环境与大地测量教育部重点实验室开放基金资助项
)
若将上式分成主项及次要项两部分，则有
F0 =
1 m ∑ s y + s yi +1 s xi − s xi +1 2 i =1 i
∂r ∂r × ∂s y ∂s x
→ →
→
→
2
∂r ∂r ∂r ∂r =( × )⋅( × ) ∂s y ∂s x ∂s y ∂s x
→ → → → → →
→
→
→
→
∂r ∂r ∂r ∂r ∂r ∂r ∂r ∂r )⋅( )−( )⋅( ) = EG − F 2 = n 2 ⋅ ⋅ ⋅ ⋅ =( ∂s x ∂s y ∂s y ∂s x ∂s x ∂s x ∂s y ∂s y
(
)(
)
3
中国科技论文在线 ___________________________________________________________________________
2 ⎡ sy 1 ⎢ 1− 3 2 ⎤ − s3 s 2 R y y 0 − 1 2 2 ⎥ ÷ ⎢1 + 2 6 R0 ⎦ sy ⎥ ⎢ ⎢ 1 − 22 2 R0 ⎢ ⎣ 2 ⎡ sy 3 − 1 ⎢ 3 ⎤ − s3 s 2 R02 y y − 3 2 2 ⎥ ÷ ⎢1 + 2 ⎢ 6 R0 ⎥ ⎦ ⎢ 1 − s y2 2 R02 ⎢ ⎣
2 利用测地坐标计算椭球面三角形的面积
2
中国科技论文在线 ___________________________________________________________________________
P3
(s x
3
,sy3
)
F
P2
(s x
2
,sy2
(
)(
)
)(
)
( (
) ( ) (
)(
)
)
(
)
( − (s
1 2 2 2 3 2 [ s3 y1 + s y1 s y3 − 2 s y1 s y3 + s y1 s y3 + s y3 − 2 s y1 s y3 s x3 − s x1 8 R02
3 y2 2 2 2 2 s y3 + s y s y3 + s 3 + s y2 s y − 2s y y3 − 2 s y 2 s y 3 3 2 2
L
r ′ d sx
T
dσ
r ( s x , s y ) d sy P
Q
(B0 , L 0 )
P0
图 1：椭球面上基于测地格网的微分面积
Fig.1 Differential area based on geodesic Grid on ellipsoidal surface
由向量运算中的 Lagrange 恒等式得到：
∂r ∂r ds x ds y , PL = PQ = ∂s x ∂s y
→ ⎛ → ⎞ ∂r ⎜ ∂r ⎟ dσ = ⎜ ds y × ds x ⎟ ∂s x ⎜ ∂s y ⎟ ⎝ ⎠
→
→
(1)
微分四边形 PQTL 的面积可用来代表该测地格网在椭球面上的微分面积 dσ ，即有