子午線長度正反算公式及相應系數(shù)值

兩種子午線長度正反算公式及相應系數(shù)值楊國清李曉(鄭州測繪學校河南鄭州450015)TwoFormulasAndRelevantCoefficietsForComputingTheLenthOfMeridian

YangGuoqingLiXiao摘要：本文給出了兩種高精度的子午線長度直接正反算公式，以及公式中相應于四種橢球的系數(shù)值。關鍵詞：子午線長度；正反算；公式系數(shù)在大地測量的高斯投影計算中，需要計算中央子午線上緯度為B的點到赤道的子午線弧長S，或者說縱坐標X的值（常記作Xb）；而在高斯投影反算中，需要計算中央子午線0上坐標為X的點所對應的緯度（常記作Bf）。所以，子午線弧長計算及其反算是經(jīng)常要進行的計算工作。這個問題其實是個老問題，已有解決方法；但老問題的解決方法也常有需要完善和補充的地方。例如，一些書上提供了B一S的正算公式，而S—B的反算則是告訴要用迭代的方法計算；而迭代計算給人一種不直觀的感覺，不便手工計算，即使是用計算機編程計算，相對于直接計算公式，迭代計算編程時也更加麻煩。文獻［2］,［3］,［4］各給出了一種直接的反算公式。文獻［3］未給出公式推導過程；文獻［2］和［4］都是采用某種迭代法來獲得式中所需的系數(shù)值。三文獻都只給出了公式中所需的與54坐標系橢球和80坐標系橢球相應的系數(shù)值，而現(xiàn)在常用的WGS84的橢球和我國2000坐標系橢球的相應系數(shù)值則不見公開發(fā)表。本文給出了兩種形式的直接計算公式，提出了有別于前述文獻的公式系數(shù)計算方法。筆者計算了兩種公式中常用橢球的子午線弧長正反算有關系數(shù)值，借此文奉獻給讀者。一、公式推導倍角型公式的推導（1）正算公式（B—S）的推導緯度為B的點至赤道間的子午線弧長S的積分計算公式為TOC\o"1-5"\h\zS=a（1-e2）IB（1-e2sin2B）-32dB （1）0此式無法直接積分，將被積函數(shù)按泰勒級數(shù)展開(1一e2sin2B)-32=1+3e2sin2B+15e4sin4B+35e6sin6B+315e8sin8(1一e2sin2B)-32 8 16 128693 e10sin10B+— （2）256為了積分方便，把sinB的各偶次乘方化為倍角的余弦函數(shù)，再進行合并整理，可把上式化成下面的形式(1-e2sin2B)-32=A-Bcos2B+Ccos4B-Dcos6B+Ecos8B-Fcos10B+—⑶式中的A,B,C,D,E,F為與橢球元素有關的常系數(shù)，稍后再給出其表達式。將（3）式代入（1）式

S=a(1-e2)jB(A-Bcos2B+Ccos4B-Dcos6B+Ecos8B-Fcos10B+…)dB0進行積分，并將結果寫成下面的形式TOC\o"1-5"\h\zS=AB一Asin2B+Asin4B一Asin6B+Asin8B一Asin10BH— (4)0 1 2 3 4 5式中，/ 3 45 175 11025 43659 、人、式中，A=(1+—e2+——e4+ e6+ e8+ e10+…)a(1-e2)0 4 64 256 16384 65536A= (3e2+15e4+竺e6+2205e8+72765e10+--)a(1-e2)/24 16 512 2048 65536A= (He4+些e6+至e8+冬e10+ 2)/464 256 4096 16384A=(史e6512+空e8+工eA=(史e6512(315(163842048 131072e8+獨5e10+...)a(1-e(315(16384A=5(A=5(693(131072e10+…)a(1-e2)/10上面各式右邊第一個括號內的值便分別是前面（3）式中的A,B,C,D,E,F。式⑷中"喝一項已總是小于0式⑷中"喝一項已總是小于0.03mm.，因而一般子午線長度計算取下式精度已足夠（5）S=AB-Asin2B+Asin4B-A（5）01 2 3（2）反算公式（S—B）的推導將⑸式兩邊同除以%，經(jīng)移項得B=—+—^sin2B一冬sin4B+少sin6BA0 A0 A0 A0再令E=S冗，很明顯，因氣牝a（橢球長半軸），所以可認為E是B的一個近似值。事實上，E和B兩者之差在任何時候都小于0.0026弧度。如果將上式右邊的B改寫成E,并將各系數(shù)換成%氣,R3，從而寫成(6)B=E+%sin2E+R2sin4E+R3sin(6)則合適調整一下R1,R2,R3各系數(shù)的大小，就能使公式兩邊相等。下面的問題就是如何確定R『R2,R3。取B=45°，用正算公式精確算出子午線長S45再算出E45=S45/A0，代入（6）式可得方程兀?45/180=E45+R]sin2E45+R2sin4E45+R3sin6E45如果我們分布均勻的取B=22.5°,45°,67.5°，利用上述方法，則可用（6）式列出3個方程，將3個方程聯(lián)立解算，可唯一的確定R1,R2,R3三個待定系數(shù)。經(jīng)實算檢驗，這樣確定的系數(shù)，用（6）式計算，在0°至90。的范圍內，誤差均小于0.000005〃。順便提及，筆者是將（4）式取至A4sin8B項，并計算了相應的反算系數(shù)，達到了正算誤差小于0.0000001mm，反算誤差小于0.000000015〃的精度，并將這樣的計算值作為標準值檢驗（5），（6）兩式的實算誤差。乘方型公式的推導將（2）式代入（1）式，得S=jB(1+-e2sin2B+15e4sin4B+35e6sin6B+315e8sin8B+竺ei0sini0)dB0 2 8 16 128 256直接積分，經(jīng)合并整理后，便可得到后面的式（10）。再利用類似于從（5）式推出（6）式的方法，可從（10）式推出（12）式。并同樣用解方程的方法可確定式中所需系數(shù)。二、正弦函數(shù)倍角型的公式（正反算聯(lián)合誤差小于0.2mm）正算公式（B—S）：TOC\o"1-5"\h\zS=AB-Asin2B+Asin4B-Asin6B （7）0 1 2 3用本文提供的系數(shù)計算，此式誤差小于0.035mm。反算公式（S一B）：E=S/A0 （8）B=E+%sin2E+R2sin4E+R3sin6E （9）上面的公式中，B 點的緯度，S——赤道至該點的子午線弧長，A,A,A,A，R,R,R——系數(shù)，其大小與橢球參數(shù)有關。（注意此行的改動）0 12 3 1 2 3用本文提供的系數(shù)計算，此式誤差小于0.000005〃，換算為長度，小于0.14mm。計算中應注意：式（7）中，與A0相乘的B應以弧度為單位；式（8）和式（9）的計算結果也是以弧度為單位的；度化弧度和弧度化度用到的n為3.1415926535898n的取位應與結果要求的精度相應，要取10位或更多位有效數(shù)字。常用橢球的正弦函數(shù)倍角型子午線弧長正反算公式中的系數(shù)值見表1和表2。正如前面提及的，以上公式取至了6倍角項，如果取至8倍角項則能達到更高的精度。筆者也算出了相應的全套系數(shù)，只是實踐中大約用不上，因而未在此提供。三、正弦函數(shù)乘方型的公式（正反算誤差小于0.0045mm）有人喜歡用正弦函數(shù)乘方型的公式，這里再給出一個這種類型的公式和相應參數(shù)。正算公式(8-S)：TOC\o"1-5"\h\zS=CB-cosB（CsinB+Csin3B+Csin5B+Csin7B） （10）用本文提供的系數(shù)計算，此式誤差小于0.0045mm。反算公式（S—B）：E=S/C0 (11)B=E+cosE(ksinE-ksin3E+ksin5E-ksin7E) (12)12 3 4上面的公式中，B點的緯度，S——赤道至該點的子午線弧長，C,C,C,C,C，k,k,k,k——系數(shù)，其大小與橢球參數(shù)有關。0 1 2 3 4 1 2 3 4用本文提供的系數(shù)計算，此式誤差小于0.00000015〃，換算為長度，小于0.0045mm。計算中的注意事項同前面的正弦函數(shù)倍角型公式。常用橢球的正弦函數(shù)乘方型子午線弧長正反算公式中的系數(shù)值見表1和表2。表11954坐標系和1980坐標系橢球基本數(shù)據(jù)和子午線弧長計算參數(shù)參考橢球名1954北京坐標系橢球1980西安坐標系橢球橢球長半軸a63782456378140基本扁率a1/298.31/298.257數(shù)據(jù)短半軸b6356863.01877316356755.2881575第一偏心率平方e26.693421622966E-036.694384999588E-03第二偏心率平方e'26.738525414683E-036.739501819473E-03極曲率半徑C6399698.9017836399596.651988子午A 0 6367558.496876367452.13279線弧A 1 16036.4802716038.52823長倍A 2 16.8280716.83265角型A 3 0.021980.02198公式R 1 2.518464777E-032.518828476E-03參數(shù)R 2 3.6998859E-063.7009546E-06R 3 7.4449E-097.4655E-09子午C 0 6367558.4968756367452.132788線弧C32005.77987432009.857529長乘C 2 133.923808133.960153方型C 3 0.6972630.697553公式C 4 0.0039380.003940參數(shù)K 1 5.0517739022E-035.0525055926E-03K 2 2.98386764E-052.98473300E-05K 3 2.414996E-072.416048E-07K2.1866E-092.1879E-09兀=3.1415926535898 p〃=206264.8062471〃

表2WGS84坐標系和2000國家坐標系橢球基本數(shù)據(jù)和子午線弧長計算參數(shù)參考橢球名WGS84坐標系橢球2000國家坐標系橢球橢球長半軸a63781376378137基本扁率a1/298.2572235631/298.257222101數(shù)據(jù)短半軸b6356752.31424526356752.3141404第一偏心率平方e26.694379990141E-036.694380022901E-03第二偏心率平方e'26.739496742276E-036.739496775479E-03極曲率半徑C6399593.6257586399593.625864子午A 0 6367449.145826367449.14577線弧A 1 16038.5086616038.50874長倍A 2 16.8326116.83261角型A 3 0.021980.02198公式R2.518826585E-032.518826597E-03參數(shù)R 2 3.7009504E-063.7009491E-06R 3 7.4655E-097.4481E-09子午C 0 6367449.1458236367449.145771線弧C 1 32009.81853032009.818687長乘C 2 133.959889133.959891方型C 3 0.6975510.697551公式C0.0039400.003940參數(shù)K 1 5.0525017879E-035.0525018127E-03K 2 2.98472850E-052.98472853E-05K 3 2.416042E-072.416042E-07K2.1878E-092.1878E-09兀=3.1415926535898 p〃=206264.8062471〃筆