白塞尔大地主题解算

白塞尔大地主题解算。

方向：

学号：

姓名：

\

一．基本思路：

；

基本思想：将椭球面上的大地元素按照白塞尔投影条件投影到辅助球面上，继而在球面上进行大地主题解算，最后在将球面上的计算结果换算到椭球面上。其关键问题是找出椭球面上的大地元素与球面上相应元素之间的关系式，同时解决在球面上进行大地主题解算的方法。

正算流程：

1.计算起点的归化纬度

2.计算辅助函数值，解球面三角形

3.按公式计算相关系数A,B,C 以及α,β

4.计算球面长度

#

5.计算纬度差改正数

6.计算终点大地坐标及大地方位角

、

011122S B C A

{sin (cos )}

σσσ=-+101010

22222sin ()sin sin cos cos σσσσσσ+=+101010

22222cos ()cos cos sin sin σσσσσσ+=-001101

522B C A

[cos ()]sin ()

σσσσσσ=++++010122L A sin [(sin ()sin )]

λδασβσσσ-==++-2111u u u A sin sin cos cos cos sin σσ

=

+222

222

2222

2222221111u B u B W u B u B B arctan e u sin cos cos tan sin tan cos ⎧==⎪⎪⎨⎪=⎪⎩⎡⎤==--1

111

A arctan u u A sin sin []

cos cos sin sin cos σλσσ=-21L L λδ

=+

-

反算流程：

1.辅助计算

2.用逐次趋近法同时计算起点大地方位角、球面长度及经差，第一次趋近时，取δ＝０。

!

计算下式,重复上述计算过程2.

3.计算大地线长度S

;

4.

计算反方位角

二．已知数据

L λδ=+211212u p

A u u u u q

sin cos tan cos sin sin cos cos λλ=

=

-2121p p u q b b A arctan

q

sin cos cos λλ==-=11

p A q A sin cos tan cos σσ

+=

11

p A q A sin sin cos σ=+12a a cos cos σλ

=+arctan sin cos σσσ⎛⎫

=

⎪⎝⎭

011

A u A sin cos sin =111u A tan tan sec σ=21+σσσ

=02122L A sin [(sin sin )]

λδασβσσ-==+-L λδ

=+11222222S A B C B C sin (cos )sin (cos )

σσσσσ=++-+

…

三．源代码：

#include <>

#include <>

#define e //克拉索夫斯基椭球体第一偏心率

void main()

{

int k,B10,B11,L10,L11,A10,A11,B20,B21,L20,L21,A20,A21;

double B12,L12,A12,B22,L22,A22;

double B1,L1,A1,S,B2,L2,A2,L,pi;

【

double A,B,C,afa,beta;

double a1,a2,b1,b2,p,q,x,y;

double

W1,W2,sinu1,sinu2,cosu1,cosu2,sinA0,cotsigma1,sin2sigma1,cos2sigma1,sigma0,sin2,cos2,sigma ,sins,coss,delta0,delta,lamda;

pi=4*atan(1);

printf("白塞尔大地主题正算请输入1\n白塞尔大地主题反算请输入2\n");

scanf("%d",&k);

if(k==1)

{

printf("请输入大地线起点纬度B经度L，大地方位角A及大地线长度S:\n");

scanf("%d%d%lf%d%d%lf%d%d%lf%lf",&B10,&B11,&B12,&L10,&L11,&L12,&A10,&A11,&A1 2,&S);

`

B1=(B10+(float)B11/60+B12/3600)*pi/180;

L1=(L10+(float)L11/60+L12/3600)*pi/180;

A1=(A10+(float)A11/60+A12/3600)*pi/180;

W1=sqrt(1-e*e*sin(B1)*sin(B1)); //计算起点规划纬度

sinu1=sin(B1)*sqrt(1-e*e)/W1; //计算起点规划纬度

cosu1=cos(B1)/W1; //计算起点规划纬度

sinA0=cosu1*sin(A1); //计算辅助函数值

cotsigma1=cosu1*cos(A1)/sinu1; //计算辅助函数值

sin2sigma1=2*cotsigma1/(cotsigma1*cotsigma1+1); //计算辅助函数值

cos2sigma1=(cotsigma1*cotsigma1-1)/(cotsigma1*cotsigma1+1); //计算辅助函数