空间直角坐标转换大地坐标的直接解法Word文档下载推荐.docx

空间直角坐标转换大地坐标的直接解法Word文档下载推荐.docx

《空间直角坐标转换大地坐标的直接解法Word文档下载推荐.docx》由会员分享，可在线阅读，更多相关《空间直角坐标转换大地坐标的直接解法Word文档下载推荐.docx（18页珍藏版）》请在冰点文库上搜索。

鉴于以上原因,本文导出直接计算大地纬度B的公式更为简单、直观。

其中

Δ=He≈Δ0+dΔ,

N+H

Δ=0He0,

00

02dΔ=2sin2B0-4（1-e2sin2B0）Y（1-e）

2

X=（N+H）cosBcosL,Y=（N+H）cosBsinL,

Z=[N（1-e）+H]sinB.

（1）

当已知地面点的空间直角坐标X、Y、Z求其大地坐标L、B及大地高H时,有以下直接计算公式:

L=arctan（Y/X）,B=arctanH=

Y（1-e+Δ）

（2）

XY

-N=-N.

cosBcosLcosBsinL

（3）

1　公式推导

空间直角坐标与大地坐标的关系为

收稿日期:

2019204217

作者简介:

王仲锋（1962-）,男,教授,博士.

（4）

N

=a/

1-e2sin2B0,

H0

=

cosB0

sinL

-N0

（5）B

=arctan{ZsinL/[Y（1-e2

）]},

（6）e2=（a2-b2）/a2.

（7）

在式

（2）中,可由式

（1）直观导出L和H,对于求B可做如下推导:

由式

（1）知

Z

Zsin2

Y/sinL=LY=（）（N+H）

tanB=

1-

Ne

）+（e2-Ne

（1-e2N+H

）tanB=

（1-e2He

）+

tanB.

（8）

令

Δ=HeN+H

.

（9）

将式（9）代入式（8）整理,即可得到式

（2）中求解B

的严密公式。

但是,由于

1-e2sin2B,H=

（Y/cosBsinL）-N,N和H均是B的函数,用式（3）求解Δ时做近似处理,即将其中的N换成N0,H换成H0

这样避免迭代计算。

从理论上讲,上述

近似处理必将对B的解造成影响,但这种影响是微不足道的,基本上属于计算取舍误差的范畴。

2　公式分析

2.1　

Δ的最大值由于Δ=He

H+N=

H+a/

1-e2

sin2

（B

10）

故

Δ2

max

=HeHe

H+N=+a

.（11）

minH取西安80坐标系的参数,即a=6378140m,e2=

[1**********]049取大地高H=10000m（在地球表面上,这一高度实际上是少见的,它已超过了珠峰的高度）时,算得Δmax=1.05×

10-5。

而在珠峰峰顶

B=27°

59′16.94241″雪面大地高H=8821.4016m

时,Δ=9.239193×

10-6。

可见Δ是一个小量,在多

数地区Δ远远小于10-5。

2.2　B0与B的差值

由于Δ是一个小量,故将

B=arctanY（1-e2

+Δ）

在Δ=0处用泰勒级数展开,并取一次项得

B=arctanZsinL

Y（1-e2

）

-

（ZsinL）2+[Y（1-e2）]

Δ=B0+ΔB,（12）其中

ΔB=-（ZsinL）2+[Y（1-e2）]

2,

（13）

将式

（1）代入式（13）,可得

（N+H）

+

He

=（1-e2

+Δ）,

（14）

整理得

ΔB=-（2（1-e2+Δ）2sin2B+（1-e2）2cos2B

（15）

考虑到Δ是一个小量,故在上式右端的分子和分母中取（1-e2+Δ）≈（1-e2）,于是有（秒值）

ΔB″≈-ρ″sin2B2（1-e2）Δ≈-ρ″sin2B22（1-e2）

Δ

（16）

当取大地高H=10000m和B=45°

（此时ΔB″=

max）时,代入

ΔB″≈-ρ″sin2B2（1-e2）

Δ=

-ρ″sin2B2

2（1-e2

）H+a/

1-e2sin2

17）

算得ΔB″max≈-1.08″（此时,B=45°

在中国的这一地区,H远远小于10000m）。

在珠峰峰顶,ΔB″=0.79506″。

可见,就地面点而言,B0与B的差值在极端情况下约差1″,而多数地区均远远小于1″。

2.3　N0与N的差值

由于

sinB=sin（B0+ΔB）≈sinB0+cosB0

ΔB″

ρ″

sin2

B=sin2

B0

+sin2B2

ΔB″ρ″+cosB0ΔB

″22ρ″

≈

sin2B0+sin2B0ΔρB″20

″

=sinB+ε.

ε=sin2B0

（18）

N=a/

1-e2sin2B=a/

1-e2（sin2B0+ε

）=a/

1-e2sin2B0-e2ε=

a/

（1-e2sin2B0）（1-[e2ε/（1-e2sin2B0

）]=

N0

/

（1-[e2

ε/（1-e2

）],

即

=N

220

（1-[e2ε/（1-esinB）]=

再将式（21）、式（16）和式

（2）及dB″=ΔB″代入上式

（19）

1-δ,

有

ΔB″tanB0-NeΔd=

（N+H）ρ″2（1-e0sin2B2）

δ=e2ε/（1-esinB）=2

eΔB″sin2B0/[（1-e2sin2B0）ρ″].

（20）

20

Ne2sin2BtanB-Δ=-（N+H）2（1-e2sin2B0）2（1-e2）

由于δ是一个很小的量,故将式（19）在δ=0处用泰勒级数展开,并取一次项可得

N=N（1-δ/2）=N-（Nδ）/2,即

N-NN

Ne-（N+H）

sinB-4（1-e2sin2B0）

（1-e2）

02-sin2B0-（.2

41-e2sin2B0）Y（1-e）

δ220==1/,20

ΔB″esin2B2

=N≈2（1-e2sin2B0）ρ″

（24）

或

N-N

取Δ=Δ0,B=B0代入上式得

002022020dΔ≈-sinB-.2

Y（1-e）4（1-e2sin2B0）

（25）

（21）

02200

N.2220

4（1-e）（1-esinB）

取B0=45°

及ΔB″=1″时,（N-N的相对误差很小。

可知,式（24）便是式（3）中的第2项。

N）/N=

3　算　例

例1设L=124°

B=44°

H=160m取西安80坐标系的参数,利用式

（1）算得X=-2569823.337900m,Y=3809919.776743m,

Z=4408204.814268m现用文中给出的公式直接

1/60000000。

可见,以B0代替B计算出的N0与2.4　H0与H的差值

由式

（2）知

H=

Y

-N,

cosBsinL

反算B,并与已知值进行比较。

解:

因为

L=arctan（Y/X）=120°

00′00.0000000167″,

故有

02

YYdB″B=arctan[ZsinL/Y（1-e）]=44°

00′00.0174018213″,dH=d-N=tanB-dN=cosBsinLcosBsinLρ″0

N=a/1-e2sin2B0=6388466.92743915m,

（N+H）tanBdB″（22）-dN,00

ρ″H=-N=160.51860893m.0

dH=H-H0,dN=N-N.

　　由式（22）看出,由于N的值较大,dH受dB″影响亦较大。

考虑到式（16）、式（21）和式（22）可近似为

dH≈-+

（1-e2）cosB0sinL

（23）

所以Δ0=

H0e2

10-7,00=1.[1**********]3×

H+N

020dΔ=-sin2B0-2

Y（1-e）

4（1-e2sin2B0）

-10

-5.4509966474×

10=

-7

2200

2.5　dH与dN对Δ的综合影响

对式（9）求微分得

Δ=d

dH-HdN.2

N（N+H）

-0.[1**********]74×

102

）Y（1-e+Δ

Δ=Δ0+dΔ=1.[1**********]6×

10-7,

B=arctan

=44°

00′00.0000051504″,

将式（22）代入上式有

dΔ=

（N+H）2

sin2B

ΔB″Δ=-0.01739669″≈-ρ″（,

21-e2）

B+ΔB″=44°

00′00.00000513″,

dB″-dN-HdN

（N+H）tanB=

ρ″N

1-e2sin2B=6388466.925631m,

-N=159.999927m.

（下转第12页）

NedN

tanBdB″-.

（N+H）ρ″N

表3　插值坐标结果比较

插值时间/s

卫星坐标ΔX/m

卫星坐标ΔY/m

卫星坐标ΔZ/m

356410-6.30186-2.56248-3.[1**********]-6.27646-2.55479-3.[1**********]-6.25107-2.54709-3.[1**********]-6.22570-2.53938-3.[1**********]-6.20194-2.53167-3.3566356460-6.17500-2.52395-3.[1**********]

-6.14966

-2.51623

-3.345

36

图1　PRN2精密星历广播星历比较图

4　结　论

在实验过程中,比较成功获取广播星历的所有卫星,时间持续一周。

由于卫星分布、地域限制等因素,部分位于精密星历中的卫星接收不到广播星历,故没

有参与比较。

1）通过比较可以看出,内插后的精密星历和广播（上接第9页）

　　由以上计算过程可以看出,本例直接反算出的B与已知值仅差5.15×

10-6s,B0+ΔB″与已知值仅差5113×

10-6s,完全可以满足精密大地测量的需要。

例2已知珠穆朗玛峰峰顶B=27°

59′16194241″,L=86°

55′31.72137″,以及雪面大地高H=8821.4016m,算得X=302726.854413m,Y=5636102.390135m,Z=2979527.619433m,其坐标系为北京54坐标系,现用本文给出的公式直接反算B,并与已知值进行比较。

利用直接法算得L=86°

55′[1**********]54″,B0

=27°

59′17.7373370011″,N0

=6382951.34371611m,

=8834.42452880m,Δ2

=9.[1**********]7×

10-6,dΔ=-0.[1**********]8×

10-6,Δ=9.2376938839×

10

-6

B=27°

59′16.94241162″,N=

6382951.275382m,H=8821.401627m。

可以看出,本例直接反算出B与已知值差1.62×

星历还是存在较大的误差,在X方向、Y方向、Z方向、

都产生了m级的误差。

由此可以看出,广播星历精度

并没有很大的改进,对于高精度的测量工作还是选取精密星历为益。

2）对于拉格朗日多项式,在插值广播星历时,阶数跟所选取得样本点个数相同,在一个时段内其插值精度受阶数影响,同时还受样本点间的时间间隔影响。

插值阶数不宜过大,当阶数过大时会在插值区间两端产生较大跳动,中间部分则比较吻合,即所谓的龙格现象,一般在5阶左右较为合适。

3）对于切比雪夫多项式,在插值精密星历时,要达到cm级的定位精度,至少需要10阶。

10阶以上对于坐标精度没有显著提高。

参考文献

[1]崔先强,焦文海,秦显平.GPS广播星历参数拟合算法的探讨[J].测绘科学,2019,31

（1）:

25226.

[2]李征航,黄劲松.GPS测量与数据处理[M].武汉:

武汉大学

出版社,2019.

[3]魏子卿,葛茂荣.GPS相对定位的数学模型[M].北京:

测绘

出版社,1998.

[4]蔡艳辉,程鹏飞,李夕银.卫星坐标的内插与拟合[J].全球定

位系统,2019（3）:

10213.

[5]李庆扬,关　治,白峰杉.数值计算原理[M].北京:

清华大学

出版社,20001.

[责任编辑:

李铭娜]

4　结束语

根据地面点地心空间直角坐标与地心大地坐标的关系式,导出由空间直角坐标求解大地坐标的简明直

接计算公式,并对计算误差进行分析讨论。

理论分析和实际验算结果表明,该直接解法引起的纬度误差不大于10-5s,可以满足精密大地测量需要。

[1]孔祥元,郭际明,刘宗泉.大地测量学基础[M].武汉:

武汉大

学出版社,2019.

[2]周忠谟,易杰军,周　琪.GPS卫星测量原理与应用[M].北

京:

测绘出版社,1999.

[3]束婵方,李　斐,沈　飞.空间直角坐标向大地坐标转换的新

算法[J].武汉大学学报:

信息科学版,2019,34（5）:

14216.

[4]李　宏.空间直角坐标至大地坐标的直接严密互换[J].测绘

与空间地理信息,2019（6）:

48250.

[5]张　锐.坐标转换中大地高对平面坐标和高程的影响[J].测

绘工程,2019,14（4）:

17220.[责任编辑: