白塞尔大地主题解算课程设计

大地测量学课程设计

设计题目:白塞尔大地主题解算 学 院: 矿业学院 专 业: 测绘工程 班 级: 测绘122 学 号: 1208010255 学生姓名: 郑英智 指导教师: 张 俊 2014年12月26日 目录

1. 白塞尔法大地主题正算步骤 1 2. 白塞尔法大地主题反算步骤 3 3. 同一平行圈弧长、子午线弧长与大地线比较大小 4 4. 程序代码 5 5. 演算示例 9 6. 参考文献 12 7. 心得体会 13

教师评语 14

白塞尔大地主题解算 摘要:白塞尔法解算大地主题的基本思想是将椭球面上的大地元素按照白塞尔投影条件(,

椭球面大地线投影到球面上为大圆弧;,大地线和大圆弧上相应点的方位角相等;,球面上

任意一点的纬度等于椭球面上相应点的归化纬度。)投影到辅助球面上，继而在球面上进行

大地主题解算，最后再将球面上的计算结果换算到椭球面上。 1.白塞尔法大地主题正算步骤

2已知 、、、、()、S，计算、、()。 BLABLAaAAe1112221221 (1)将椭球面元素投影到球面上 2utanu,1,etanB求:,由B1111 A,计算辅助量和 ,01

sinA,cosusinAtan,,tanusecA , 011111 , ,计算球面长度，将S化为

,,,S，,sin,cos(2,，,)，,sin2,cos(4,，2,) 11式中系数分别为: ''''246k7k15k,,,,(1,，,，？) , bAb4256

''246,Bkk37k'',,,,(,，，？) A48512 ''46,Ckk'',,,,(,，？) A128128 22'2k,ecosA 0

24635kkk ,1，,,，？A4256246kk15kB,,，，？416512 46kk3C,,？128512

,,,,S上式右端含有，因此需要迭代计算。第一次迭代取近似值，第二次计算取 0

,,,S，,sin,cos(2,，,)，,sin2,cos(4,，2,) 01010 - 1 -

,以后计算用代换,代入上式迭代计算，直到所要求的精度为止。一般取0 ''|,-,|,0.001。 0

(2)解算球面三角形 ,计算A2

，，cosusinA11 tanA,2,,cosucos,cosAsinusin,,111，， ,计算 u2

sinu,sinucos,，cosucosAcos,或 tanu,,cosAtan(,，,)2111221 ,计算 ,

，，,sinAsin1 ,tan,,,cosucos,sinusin,cosA,111，， (3)将球面元素换算到椭球面上 ,由u求B 22

sinu22'tanB, 或tanB,1，etanu 222221,e1,sinu2 ,将球面经差化为椭球面经差l，求L ,2

''',,,,-sinA,,，,sin,cos(2,，,)，,sin2,cos(4,，2,) l011 式中

2462eeee234'2'2',(，，),(1，e)k，ek, 281616128 2224ee'2'',(1，e)k,k, 1632 24e'',,k 256 2'22k,ecosA0

''',0.0002式中的最大值为，故在计算时通常可以略去不计。 ,象限的判定 - 2 -

符号 sinA+ + - - 1 tan,符号 + - + - |,|,-|,||,|,,- |,| ,

l L,L，21

符号 sinA- - + + 1 符号 tanA+ - + - 2

|A|,,|A|,，|A|2,,|A|A2222 2 其中、为锐角。 |A||,|2 2.白塞尔法大地主题反算步骤

2已知 、、、、，计算()、()、S。 BLBLAAAAae1122121221、 (1)将椭球面元素投影到球面上 ,由B求u

22tanu,1,etanBtanu,1,etanB ，，l= L,L112221 , a,sinusinua,cosucosu112212 , b,cosusinub,sinucosu112212 ,采用逐次趋近法，由l计算 ,

在反算中，已知椭球面上经差l，球面经差上的对应经差未知，为了由l求，由下,,

,A式可知还需计算、、,，计算又还需量，故需要进行迭代计算。 ,,01 ,,第一次趋近，取l; ， p,cosusin,q,b,bcos,212

，，,cosusinp2tanA,tanA或 ,11,,qcosusinusinucosusin,,1212，， - 3 -

判断的象限 A1 p符号 + + - - q符号 + - + -

A,|A|,,|A|,，|A|2,,|A|11111

sin,,cosusin,sinA，(cosusinu,sinucosucos,)cosA2112121 cos,,sinusinu，cosucosucos,1212 ,sin，，,,arctan, ,,cos，， 判断象限 , cos,+ - ,, |,|,-|,| sinA,cosusinA011

,,tantanusecA,111, sinu,1tan,,tanA11,sinA0,

'',,,, l+sinA,,'，,sin,cos(2,，,)，,sin2,cos(4,，2,) 011仿照上述计算步骤迭代计算，直到为止。 |,-,|,,i，1i

(2)将球面元素换算到椭球面上

1 ,,S,,-,sin,cos(2,，,),,sin2,cos(4,，2,)11, ,，，cosusin1 tanA ,2,,cosusinucos,sinucosu,1212，， 或

，，cosusinA11tanA ,2,,cosucos,cosAsinusin,,111，， A象限的判断与前面一致 2 - 4 -

3.同一平行圈弧长、子午线弧长与大地线比较大小

子午线弧长计算公式: M,a,acos2B，acos4B,acos6B，acos8B02468 3535m,2a,，，，，，mmmm？00468,2816128,157mm,24a,，，，mm268,223216,37m,4a,，m，m式中: ,46881632,

mm,68a,，6,3216,m,8a,8,128,

''Nl''bcosB, 平行圈弧长公式:() S,NcosB,bl11'''',,

不同纬度对应的一些弧长的数值

子午线弧长 平行圈弧长 B 1′ 1″ 1′ 1″

0? 110 576m 1 842.94m 30.716m 111 321m 1 855.36m 30.923m 15? 110 656 1 844.26 30.738 107 552 1 792.54 29.876 30? 110 863 1 847.71 30.795 96 488 1 608.13 26.802 45? 111 143 1 852.39 30.873 78 848 1 314.14 21.902 60? 111 423 1 857.04 30.951 55 801 930.02 15.5 75? 111 625 1 860.42 31.007 28 902 481.71 8.028 90? 111 696 1 861.6 31.027 0 0 0 利用白塞尔大地主题反算求解大地线长S

纬度为30?，经差为1?的平行圈弧长S=96 488m，两点间大地线长为96 487.595m

- 5 -

经度为30?，纬度差为1?的子午线弧长X=110 863m，两点间大地线长为110 862.869m

通过比较可知，同一平行圈或同一子午线两点间大地线长度与对应的平行圈弧长或子午线弧长相等。

4.程序代码

仅体现计算过程，不包含定义、输入输出部分 正算:

eps = e2 / (1 - e2) b = a / Sqr(1 + eps)

u1 = Atn(Sqr(1 - e2) * Tan(B1)) sinA0 = Cos(u1) * Sin(A1) cosA0 = Sqr(1 - sinA0 * sinA0) sgm1 = Atn(Tan(u1) / Cos(A1)) xk2 = eps * cosA0 * cosA0 xk4 = xk2 * xk2 xk6 = xk4 * xk2

alpha = (1 - xk2 / 4 + 7 * xk4 / - 15 * xk6 / 256) / b beta = xk2 / 4 - xk4 / 8 + 37 * xk6 / 512 gamma = xk4 / 128 - xk6 / 128

sgm = alpha * S Do sgm0 = sgm

sgm = alpha * S + beta * Sin(sgm0) * Cos(2 * sgm1 + sgm0) sgm = sgm + gamma * Sin(2 * sgm0) * Cos(4 * sgm1 + 2 * sgm0) Dsgm = Abs(sgm - sgm0) * p Loop While Dsgm > 0.0001 '计算反方位角A2 - 6 -

sinA2 = Cos(u1) * Sin(A1)

cosA2 = Cos(u1) * Cos(sgm) * Cos(A1) - Sin(u1) * Sin(sgm) tanA2 = sinA2 / cosA2 A2 = Abs(Atn(sinA2 / cosA2)) sinA1 = Sin(A1)

If sinA1 < 0 And tanA2 > 0 Then A2 = A2 If sinA1 < 0 And tanA2 < 0 Then A2 = pi - A2 If sinA1 > 0 And tanA2 > 0 Then A2 = pi + A2 If sinA1 > 0 And tanA2 < 0 Then A2 = 2 * pi - A2 '计算大地纬度B2

sinu2 = Sin(u1) * Cos(sgm) + Cos(u1) * Cos(A1) * Sin(sgm) B2 = Atn(sinu2 / Sqr(1 - e2) / Sqr(1 - sinu2 * sinu2)) '计算大地经度L2

sinl = Sin(A1) * Sin(sgm)

cosl = Cos(u1) * Cos(sgm) - Sin(u1) * Sin(sgm) * Cos(A1) tanlambda = sinl / cosl lambda = Abs(Atn(sinl / cosl))

If tanlambda > 0 And sinA1 > 0 Then lambda = lambda If tanlambda < 0 And sinA1 > 0 Then lambda = pi - lambda If tanlambda < 0 And sinA1 < 0 Then lambda = -lambda If tanlambda > 0 And sinA1 < 0 Then lambda = lambda - pi e4 = e2 * e2 e6 = e4 * e2

xk2 = e2 * cosA0 * cosA0 xk4 = xk2 * xk2 xk6 = xk4 * xk2

alpha1 = (e2 / 2 + e4 / 8 + e6 / 16) - e2 * (1 + e2) * xk2 / 16 + 3 * xk4 * e2 / 128

beta1 = e2 * (1 + e2) * xk2 / 16 - e2 * xk4 / 32 gamma1 = e2 * xk4 / 256

xx = alpha1 * sgm + beta1 * Sin(sgm) * Cos(2 * sgm1 + sgm) xx = xx + gamma1 * Sin(2 * sgm) * Cos(4 * sgm1 + 2 * sgm) l = lambda - sinA0 * xx L2 = L1 + l 反算:

eps = e2 / (1 - e2) b = a / Sqr(1 + eps)

u1 = Atn(Sqr(1 - e2) * Tan(B1))

u2 = Atn(Sqr(1 - e2) * Tan(B2)) DL = L2 - L1

aa1 = Sin(u1) * Sin(u2) - 7 -

aa2 = Cos(u1) * Cos(u2) bb1 = Cos(u1) * Sin(u2) bb2 = Sin(u1) * Cos(u2)

'获取数据逐次趋近法计算经差lambda,大地方位角A1,球面长度 lambda = DL Do

lambda0 = lambda

p1 = Cos(u2) * Sin(lambda0) q = bb1 - bb2 * Cos(lambda0) A1 = Atn(p1 / q)

If p1 > 0 And q > 0 Then A1 = Abs(A1) If p1 > 0 And q < 0 Then A1 = pi - Abs(A1) If p1 < 0 And q < 0 Then A1 = pi + Abs(A1) If p < 0 And q > 0 Then A1 = 2 * pi - Abs(A1) sinsgm = p1 * Sin(A1) + q * Cos(A1) cossgm = aa1 + aa2 * Cos(lambda0) sgm = Atn(sinsgm / cossgm) If cossgm > 0 Then sgm = Abs(sgm) If cossgm < 0 Then sgm = pi - Abs(sgm)

sinA0 = Cos(u1) * Sin(A1) sgm1 = Atn(Tan(u1) / Cos(A1)) cosA0 = Sqr(1 - sinA0 * sinA0) e4 = e2 * e2 e6 = e4 * e2

xk2 = e2 * cosA0 * cosA0 xk4 = xk2 * xk2 xk6 = xk4 * xk2

alpha1 = (e2 / 2 + e4 / 8 + e6 / 16) - e2 * (1 + e2) * xk2 / 16 + 3 * xk4 * e2 / 128

beta1 = e2 * (1 + e2) * xk2 / 16 - e2 * xk4 / 32 gamma1 = e2 * xk4 / 256

dll = alpha1 * sgm + beta1 * Sin(sgm) * Cos(2 * sgm1 + sgm) dll = dll + gamma1 * Sin(2 * sgm) * Cos(4 * sgm1 + 2 * sgm) lambda = dll * sinA0 + DL

Dlambda = Abs(lambda - lambda0) * p Loop While Dlambda > 0.0001 cosA0 = Sqr(1 - sinA0 * sinA0) xk2 = eps * cosA0 * cosA0 xk4 = xk2 * xk2 xk6 = xk4 * xk2

alpha = (1 - xk2 / 4 + 7 * xk4 / - 15 * xk6 / 256) / b - 8 -

beta = xk2 / 4 - xk4 / 8 + 37 * xk6 / 512 gamma = xk4 / 128 - xk6 / 128

Y = gamma * Sin(2 * sgm) * Cos(4 * sgm1 + 2 * sgm) S = (sgm - beta * Sin(sgm) * Cos(2 * sgm1 + sgm) - Y) / alpha

sinA2 = Cos(u1) * Sin(A1)

cosA2 = Cos(u1) * Cos(sgm) * Cos(A1) - Sin(u1) * Sin(sgm) tanA2 = sinA2 / cosA2

A2 = Abs(Atn(sinA2 / cosA2)) sinA1 = Sin(A1)

If sinA1 < 0 And tanA2 > 0 Then A2 = A2 If sinA1 < 0 And tanA2 < 0 Then A2 = pi - A2 If sinA1 > 0 And tanA2 > 0 Then A2 = pi + A2 If sinA1 > 0 And tanA2 < 0 Then A2 = 2 * pi - A2 5.演算示例(教材151页) 用克式椭球进行计算

- 9 -

用1975国际椭球进行计算

正算计算结果 30?30′00.00″ -0.505 246 786

114?20′00.00″ -0.610 0 770

225?00′00.00″ -0.707 105 781

10000000.00 0.356 537 888 S

0.530858149 0.003 350 562

-0.6097472 1.766 987 874E-06 0.792563224 4.622 109 698E-10 -0.692911606 -1.103 777 280

1.571 441 354 5.270 424 445E-03 1.055 977 930 -37?43′44.1351″

1.393 840 997 51?16′32.4976″ 1.572 478 988 50?21′22.46″ -0.609 7 747

- 10 -

反算计算结果

0.530 858 149 -0.656 869 486

-0.309 151 284 0.682 919 788

-0.526 600 070 0.400 921 954 趋近次数 计算 值 1 3 4 2

-1.100 563 429 -1.103 768 055 -1.103 777 253 -1.103 777 280 -0.705 956 271 -0.707 102 492 -0.707 105 771 -0.707 105 781 p -0.708 255 369 -0.707 109 082 -0.707 105 7 -0.707 105 780 q

224?54′24.6737″ 224?59′59.0390″ 224?59′59.9973″ 225?00′00.0002″

0.999 999 996 0.999 998 594 0.999 998 584 0.999 998 584

2.755 182 939E-04 -1.677 035 942E-04 -1.682 4 780E-04 -1.682 660 878E-04

1.570 520 809 1.572 473 3 1.572 478 973 1.572 478 9 sinA0.608 797 601 0.609 786 905 0.609 7 739 0.609 7 747 0

-0.692 114 035 -0.692 909 315 -0.692 911 599 -0.692 911 606

3.350 558 481E-03 3.350 561 868E-03 3.350 561 878E-03 3.350 561 878E-03

1.770 381 674E-6 1.766 997 603E-6 1.766 987 901E-6 1.766 987 874E-6

4.639 919 221E-10 4.622 160 703E-10 4.622 109 841E-10 4.622 109 695E-10 ,1

5.263 861 5 5.270 405 594 5.270 424 392 5.270 424 446 xx

-1.103 768 055 -1.103 777 253 -1.103 777 280 -1.103 777 280

1.571 441 354 -0.6097747

1.055 977 929 -0.505246786

1.393 840 996 50?21′22.47″ -4.367 3 449E-10 S 10000000.003 Y 6.参考文献 - 11 -

[1]郭际明,丁士俊,苏新洲,刘宗泉.大地测量学基础实践教程.武汉:武汉大学出版社,2009

[2]裴连磊. 用C语言实现大地主题解算[J]. 价值工程,2013,20:235-236. [3]丁士俊,杨艳梅,史俊波,程新明. 大地主题解算几种不同算法在计算中应注意的问题[J]. 黑

龙江工程学院学报(自然科学版),2013,03:1-5.

[4]周振宇,郭广礼,贾新果. 大地主题解算方法综述[J]. 测绘科

学,2007,04:190-191+174+200. [5]王建强,胡明庆. 贝赛尔大地主题解算分析[J]. 测绘科学,2012,01:30-31. [6]徐晓晗,谢云开,李亚军. 大地主题解算实用算法[J]. 科学技术与工程,2012,09:2062-2068. [7]许厚泽. 关于正反大地主题解算方法的

综合研究[J]. 测量制图学报,1958,04:274-288. [8]史国友,周晓明,贾传荧. 贝塞尔大地主题正解的改进算法[J]. 大连海事大学学

报,2008,01:15-19. 7.心得体会 - 12 -

本次课程设计在一开始我就有感觉到并不容易，因为平时上课翻开书就是满篇公式，完全不知所云，尽管知道只是一门基础课，但我还是深深感受到这门课的深奥与晦涩，学起来确实有些困难。可是不管怎么难，既然是基础课，为了学习以后更为专业的知识，我们还是应该去学习、了解、掌握它。

我们的设计题目是白塞尔大地主题解算，设计的第一步便是认真看书，理解白塞尔投影的三个条件，然后跟着书上的步骤一步一步往下进行。知道了原理，有了设计思路，工作才正式开始。

这次设计最让我头疼也是让我收获最多的是程序设计。虽然之前学过C语言，但是只是懂得一些基础，而程序设计中有很多变量，编写时比较麻烦，写起来十分不方便。因此，我决定用其他语言来设计程序。相比C语言来说，用VB来编写较为简便。VB是面向对象的，可以一边写一边定义，可视化的用户界面可以给人简单直观的感受。为此，我决定先学习VB再进行课程的设计。

我在图书馆借了相关的教材，配合网上的教学视频，很快就入了门。我一边学习VB，一边进行课程设计，不断地尝试与修改，一遍又一遍的调试，最后总算完成了本次设计。从基本的语法，到用户界面设计，再到程序的调试，最后完成程序设计，我从一无所知渐渐掌握知识。每一次做出成果后我都会想办法进一步改善，思考它有哪些漏洞和不足，不断完善，使程序更简练、更全面、更准确。

在这次课程设计中，我发现教材上有一些印刷错误，虽然问题不大，但对于还不熟悉这门课的我们来说，在学习过程中难免会造成疑惑，多走弯路。此外，教材上的算法虽然简单易懂，但是它在计算系数时是用椭球参数直接代入后化简的公式，一个程序只能依据一个参考椭球，不够灵活。我在此基础上稍微做了修改，可以以不同参考椭球为依据，使程序更具普遍性，让课程设计更完善。

此次设计让我明白了思路即出路，有什么不懂不明白的地方要及时请教或上网查询，只要认真钻研，动脑思考，动手实践，就没有弄不懂的知识，收获颇丰。

回顾起此次课程设计，至今我仍感慨颇多，从理论到实践，在这段日子里，可以说得是苦多于甜，但是可以学到很多很多的东西，同时不仅可以巩固了以前所学过的知识，而且学到了很多在书本上所没有学到过的知识。通过这次课程设计使我懂得了理论与实际相结合是很重要的，只有理论知识是远远不够的，只有把所学的理论知识与实践相结合起来，从理论中得出结论，才能真正学以致用，从而提高自己的实际动手能力和思考的能力。在设计的过程中遇到问题，可以说得是困难重重，但可喜的是最终都得到了解决。

本次课程不仅让我们利用大地测量学基本知识完成了椭球大地主题正、反算以及椭球上子午线弧长和平行圈弧长的计算工作，培养了我们学生思考、动手操作的能力。还让我们进一步了解和掌握了大地线是椭球面上最短程曲线的特性以及大地主题解算在测绘学科中的地位和在实际应用中的意义，让我们对理论知识有了更深层次的理解。也使我们利用计算机解决问题的基本技能得到训练和提高。更重要的是，在学习实践中，我们学会了很多学习的方法。而这是日后最实用的，真的是受益匪浅。要面对社会的挑战，只有不断的学习、实践，再学习、再实践，才能为走向工作岗位打下一定实践基础。这对于我们的将来会有很大的帮助。

8.教师评语 - 13

- - 14 -