空间直角坐标转换大地坐标的直接解法.docx

空间直角坐标转换大地坐标的直接解法.docx

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

空间直角坐标转换大地坐标的直接解法

空间直角坐标转换大地坐标的直接解法

第19卷第2期　　　　　　　　　　　　　　　　测　绘　工　程　　　　　　　　　　　　　　　　Vol.19№.22019年4月　　　　　　　　ENGINEERINGOFSURVEYINGANDMAPPING　　　　　　　　Apr.,2019

空间直角坐标转换大地坐标的直接解法

王仲锋1,杨凤宝2

（1.长春工程学院勘查与测绘学院,吉林长春130021;2.吉林水利水电勘测设计研究院,吉林长春130021）

摘　要:

根据地面点地心空间直角坐标与地心大地坐标的关系式,导出由空间直角坐标求解大地坐标的更为简洁的直接计算公式,并对计算误差进行分析讨论。

理论分析和实际验算结果表明,该直接解法引起的纬度误差不大于

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

关键词:

空间直角坐标;大地坐标;坐标转换;坐标系;大地测量

）022*******

Newstudyondirectcalculationmethodtoturnspacial

rectangularcoordinateintogeodeticcoordinate

WANGZhong2feng1,YANGFeng2bao2

（1.SchoolofProspectingandSurveying,ChangchunInstituteofTechnology,Changchun130021,China;2.SurveyDesignandResearchLnstituteofJilinWaterUtilities,Changchun130021,China）

Abstract:

Thispaperdeducesthesimplerdirectcalculatingformulafromtherelativeformulasofspacialrectangularcoordinateofpointwithitsgeodeticcoordinate;analyzesanddiscussesrelativecalculatinger2rors.Itisprovedbythetheoreticalanalysesandthepracticalcaculatingresultsthatthelatitudeerrorislessthansecond.Thisprecisioncansatisfytheneedofprecisiongeodeticsurvey.

Keywords:

spacialrectangularcoordinate;geodeticcoordinate;coordinatetransformation;coordinatesys2tem;geodeticsurvey

　　由地面点的地心空间直角坐标求解其地心大地坐标的应用越来越广泛。

例如,用GPS可直接测得地面点在WGS284下的地心空间直角坐标,若要计算地面点在WGS284下的高斯平面直角坐标,应首先将地心空间直角坐标转换为地心大地坐标。

由于地面点的地心空间直角坐标求解其地心大地坐标时“,大

[1]

地纬度B的计算比较复杂,通常采用迭代法”进行求解。

迭代法的公式虽然简明,但对一般用户而言,应用起来存在着诸多不便。

文献[2]虽然给出两种直接计算大地纬度B的公式,但所给公式亦比较复杂,而且不够直观。

从目前科技期刊中检索到的文献看,由地面点的地心空间直角坐标求解其地心大地坐标时,均使用迭代法。

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

其中

Δ=He≈Δ0+dΔ,

N+H

Δ=0He0,

N+H

00

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

2

2

2

2

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

2

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

（1）

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

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

2

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=

2

1-

Ne

N+H

tanB=

）+（e2-Ne

2

（1-e2N+H

）tanB=

（1-e2He

2

）+

N+H

tanB.

（8）

令

2

Δ=HeN+H

.

（9）

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

（2）中求解B

的严密公式。

但是,由于

N

=a/

1-e2sin2B,H=

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

这样避免迭代计算。

从理论上讲,上述

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

2　公式分析

2.1　

Δ的最大值由于Δ=He

2

2

H+N=

H+a/

1-e2

sin2

（B

10）

故

Δ2

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）]

2

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

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

2,

（13）

将式

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

2

2

（N+H）

=

1-e2

+

He

N+H

=（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

（B

.

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

ΔB″

ρ″

（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

sin2

B0

）],

即

N

=N

220

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

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

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

（19）

N

1-δ,

有

2

ΔB″tanB0-NeΔd=

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

2

2

其中

220

δ=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）

2

sinB-4（1-e2sin2B0）

2

2

220

Δ

（1-e2）

2

=

02-sin2B0-（.2

41-e2sin2B0）Y（1-e）

δ220==1/,20

ΔB″esin2B2

20

=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

cosBsinL

其中

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

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

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

dH≈-+

（1-e2）cosB0sinL

（23）

2

00

所以Δ0=

H0e2

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

020dΔ=-sin2B0-2

Y（1-e）

4（1-e2sin2B0）

220

=

-10

-5.4509966474×10=

-7

2200

N.2220

4（1-e）（1-esinB）

2.5　dH与dN对Δ的综合影响

对式（9）求微分得

Δ=d

Ne

dH-HdN.2

N（N+H）

2

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

）Y（1-e+Δ

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

B=arctan

=44°00′00.0000051504″,

将式（22）代入上式有

Ne

dΔ=

（N+H）2

2

sin2B

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

21-e2）

B+ΔB″=44°00′00.00000513″,

dB″-dN-HdN

（N+H）tanB=

ρ″N

2

N=a/

H=

1-e2sin2B=6388466.925631m,

Y

-N=159.999927m.

cosBsinL

（下转第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,

H0

=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×10-6s,完全可以满足精密大地测量的需要。

星历还是存在较大的误差,在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.[责任编辑:

李铭娜]