大地線極點(diǎn)歸化緯度的迭代求解法

李 鑫，姚德新，金立新

(1.蘭州交通大學(xué) 測(cè)繪與地理信息學(xué)院，蘭州 730070；2.地理國情監(jiān)測(cè)技術(shù)應(yīng)用國家地方聯(lián)合工程研究中心，蘭州 730070；3.甘肅省地理國情監(jiān)測(cè)工程實(shí)驗(yàn)室，蘭州 730070；4.中鐵第一勘察設(shè)計(jì)院集團(tuán)有限公司，西安 710043；5.甘肅鐵道綜合工程勘察院有限公司,蘭州 730000；6.海軍工程大學(xué)，武漢 430043)

對(duì)大地線的研究在橢球大地測(cè)量學(xué)中是非常重要的內(nèi)容，它是經(jīng)典大地主題解算基礎(chǔ)。目前國外學(xué)者對(duì)大地線的研究主要方向?yàn)榭巳R勞定理，利用其他方法對(duì)其重新論證、對(duì)廣義克萊勞方程研究等方面，國內(nèi)學(xué)者目前研究方向主要是解決疆域確定等諸多需要繪制大地線的問題[1]，文中主要研究大地線極點(diǎn)歸化緯度的求解問題。早期白塞爾提出一種需要復(fù)雜迭代計(jì)算的解算長(zhǎng)距離大地問題的公式，但當(dāng)時(shí)沒有電子計(jì)算機(jī)，所以解算起來十分復(fù)雜繁瑣，因此國內(nèi)外許多學(xué)者都致力于非迭代解法的研究，赫爾默特在正解中消除了迭代；索達(dá)諾等研究了反解的非迭代計(jì)算；紀(jì)兵借助Mathemaica代數(shù)系統(tǒng)進(jìn)行了重新推導(dǎo)，得到了形式簡(jiǎn)單、便于實(shí)用的正反解直接解形式[2]。文中主要研究基于白塞爾方程的大地線極點(diǎn)歸化緯度求解，為大地主題正反算提供一種新的思路。白塞爾微分方程將大地線映射變形為平面曲線、將大地線映射變形為橢圓弧。性質(zhì)是保持歸化緯度、大地線方位角、類歸化緯度不變。深入研究發(fā)現(xiàn)，兩個(gè)白塞爾微分方程有規(guī)律性、關(guān)聯(lián)性，有內(nèi)在的本質(zhì)聯(lián)系，都可以表達(dá)為歸化緯度的函數(shù)，也可以表達(dá)為大地緯度的函數(shù)，同樣可以表達(dá)為類歸化緯度即球面弧長(zhǎng)的函數(shù)。在當(dāng)前橢球大地測(cè)量學(xué)研究中，以白塞爾微分方程為基礎(chǔ)所表達(dá)的系列函數(shù)來研究大地線新的映射拓?fù)湫再|(zhì)，對(duì)當(dāng)前大地線的理論研究有很關(guān)鍵的承接作用，并且求解大地線極點(diǎn)歸化緯度后，可以直接計(jì)算大地主題反解中大地正反方位角，因此以白塞爾微分方程表達(dá)的系列函數(shù)為基礎(chǔ)，對(duì)大地線極點(diǎn)緯度進(jìn)行求解十分必要。

1 克萊勞定理正弦表達(dá)

由文獻(xiàn)[3]可知大地線長(zhǎng)度S與大地經(jīng)度L、大地緯度B、大地方位角A的一階微分關(guān)系式為：

(1)

式中：N為卯酉圈曲率半徑;M為子午圈曲率半徑[4-6],根據(jù)式(1)可得克萊勞方程[3]：

r·sinA=C.

(2)

式中：r為平行圈半徑，根據(jù)文獻(xiàn)[7],可將式(2)改變形式得到新形式：

(3)

式中：u為歸化緯度，un為此條大地線上最高點(diǎn)的歸化緯度。其計(jì)算方法使用克氏方程[7]。

(4)

2 重要定義推導(dǎo)

為得出大地線極點(diǎn)歸化緯度迭代式，首先必須推演得出白塞爾球面弧長(zhǎng)和球面經(jīng)差的具體表達(dá)，此表達(dá)式要以歸化緯度u為自變量且包含大地線極點(diǎn)歸化緯度un。結(jié)合二者的表達(dá)式可以推演以歸化緯度為自變量的橢球面和白塞爾球面經(jīng)度之差的具體表達(dá)式，進(jìn)一步可以得到迭代式。

2.1 白塞爾球面弧長(zhǎng)

在單位圓球面上易知大圓弧微分方程[8]：

du=cosαdσ.

(5)

式中：u為歸化緯度;α為球面方位角;σ為球面弧長(zhǎng)。根據(jù)白塞爾投影條件，大地線和大圓弧上相應(yīng)點(diǎn)的方位角相等，則式(5)可寫為：

du=cosAdσ.

(6)

由式(3)克萊勞定理易得：

(7)

對(duì)式(6)求積分可得到球面弧長(zhǎng)的歸化緯度表達(dá)式為：

(8)

2.2 白塞爾球面經(jīng)差(歸化經(jīng)度)

在單位圓球面上易知大圓弧微分方程[8]：

(9)

式中：ω為球面經(jīng)差。同理，根據(jù)白塞爾投影條件，大地線和大圓弧上相應(yīng)點(diǎn)的方位角相等，則式(9)可寫為：

(10)

式(10)結(jié)合式(3)則可得出：

(11)

將式(8)推演得到的球面弧長(zhǎng)的歸化緯度表達(dá)形式代入式(11)，對(duì)其求積分可得球面經(jīng)差的歸化緯度表達(dá)式為：

(12)

式(12)推演過程如下：

(13)

證畢。

3 白塞爾球面直角三角形定理

(14)

根據(jù)上述數(shù)值，得到白塞爾球面直角三角形示意圖如圖1所示，圖中N為球面極點(diǎn)；u0為球面赤道；大地線上三點(diǎn)的投影點(diǎn)位分別為：P0(u0,ω0)，Pi(ui,ωi)和Pn(un,ωn)，各點(diǎn)位之間和極點(diǎn)之間球面弧長(zhǎng)均在圖1表示。

圖1 白塞爾球面大地線元素示意圖

結(jié)合圖1，根據(jù)球面直角三角形Napier通用規(guī)則[9]可得：

(15)

結(jié)合三角函數(shù)對(duì)式(15)進(jìn)行形式變換，則可得到式(3)、式(8)和式(12)，對(duì)于其克萊勞定理正弦形式、球面弧長(zhǎng)和球面經(jīng)差的歸化緯度表達(dá)式的正確性得到驗(yàn)證，并且由此可知，大地線投影的球面弧長(zhǎng)σ和球面經(jīng)差ω從大地線升交點(diǎn)起算。

4 經(jīng)度縮量

經(jīng)度縮量是指大地線在橢球面和白塞爾球面上的經(jīng)度差異。由文獻(xiàn)[10]，白塞爾微分方程中橢球面上經(jīng)差L與白塞爾球面上經(jīng)差ω微分關(guān)系式為：

(16)

結(jié)合式(11)，則白塞爾經(jīng)差微分方程變換為：

(17)

式(17)即為大地線經(jīng)差與歸化緯度、球面弧長(zhǎng)的微分關(guān)系式。對(duì)式(17)求積分，將其等式右邊級(jí)數(shù)展開：

(18)

式(18)中e為橢圓第一偏心率[11]，又根據(jù)式(8)歸化緯度和球面弧長(zhǎng)的關(guān)系可得：

(19)

結(jié)合(18)、(19)兩式并且整理形式后，式(17)可得到如下表達(dá)：

(20)

對(duì)式(20)第一項(xiàng)求積分可得：

(21)

積分過程如式(22)。

(22)

式(20)后三項(xiàng)中正弦函數(shù)部分通過倍角函數(shù)公式轉(zhuǎn)化一次冪余弦函數(shù)：

(23)

結(jié)合上文推演，對(duì)式(20)求積分得：

(24)

將式(24)整理可得：

(25)

將式(25)簡(jiǎn)寫：

(26)

式(26)中：

k2=

(27)

其中c=cosun。

式(26)進(jìn)一步可寫為；

(28)

式(28)中l(wèi)為經(jīng)度之差，且式中?具體表達(dá)式為：

(29)

式(28)中等式右邊第二項(xiàng)為經(jīng)度縮量，定義為ε，具體表達(dá)式為：

(30)

(31)

一個(gè)基本弧長(zhǎng)的歸化經(jīng)度，亦即極點(diǎn)un的橢球面經(jīng)差(歸化經(jīng)度)定義為ln，具體表達(dá)式為：

(32)

基本弧長(zhǎng)的橢球面經(jīng)差(歸化經(jīng)度)小于90°，亦即經(jīng)度有縮量，稱為經(jīng)度縮量。經(jīng)度縮量大小取決于兩個(gè)元素，即偏心率和大地線極點(diǎn)緯度。大地線經(jīng)差(歸化經(jīng)度)從大地線升交點(diǎn)起算。大地線升交點(diǎn)至極點(diǎn)的經(jīng)差為一個(gè)基本弧長(zhǎng)的歸化經(jīng)度。大地線升交點(diǎn)經(jīng)度加上兩個(gè)基本弧長(zhǎng)的歸化經(jīng)度等于大地線降交點(diǎn)經(jīng)度。

5 大地線極點(diǎn)歸化緯度求解

利用(8)、(12)兩式，或者白塞爾球面直角三角形的Napier法則，可得球面經(jīng)差余弦表達(dá)：

(33)

根據(jù)上文，結(jié)合式(33)，則可得到：

cosΔσ12=

sinu2sinu1+cosu2cosu1cos(Δl12+Δε12).

(34)

結(jié)合上文，證明過程如下：

cos(σ2-σ1)=cosσ2cosσ1+sinσ2sinσ1=

cosu2cosu1cosω2cosω1+cosu2cosu1sinω2sinω1+

sinu2sinu1=cosu2cosu1cos(ω2-ω1)+sinu2sinu1=

sinu2sinu1+cosu2cosu1cosΔω12，

(35)

證畢。

同理結(jié)合式(33)亦可得：

(36)

式(36)中：

Δω12=Δl12+Δε12.

(37)

其中Δε12是橢球面和白塞爾球面經(jīng)度差之差，定義為經(jīng)度縮量之差。結(jié)合經(jīng)度縮量部分的內(nèi)容，經(jīng)過形式變換，得到經(jīng)度縮量之差的表達(dá)式為：

(38)

結(jié)合上部分內(nèi)容，若大地線上的兩點(diǎn)不跨越大地線極點(diǎn)，則此大地線上兩點(diǎn)各元素如圖2所示：l1,l2和ln分別是大地線上點(diǎn)1、2和此大地線極點(diǎn)的大地經(jīng)度與此大地線升交點(diǎn)大地經(jīng)度L0之間的經(jīng)度之差。

圖2 兩點(diǎn)不跨越大地線極點(diǎn)示意圖

根據(jù)式(34)、式(36)、式(37)，則迭代式確定為：

(39)

又因?yàn)棣う?2是微小量，故初值可選定為：

(40)

(39)、(40)兩式即為大地線極點(diǎn)的歸化緯度計(jì)算式。

6 算例分析

大地線極點(diǎn)歸化緯度最直接的應(yīng)用是大地坐標(biāo)方位角的求解，為驗(yàn)證大地線極點(diǎn)歸化緯度正確性，根據(jù)大地線長(zhǎng)短，選取短距離(100 000 m以內(nèi))和超長(zhǎng)距離(5 000 000 m～10 000 000 m以內(nèi))的大地線進(jìn)行正反大地方位角的求解，并且將求解結(jié)果與傳統(tǒng)大地主題反算得到的正反大地方位角對(duì)比，從而驗(yàn)證大地線極點(diǎn)歸化緯度的正確性。

首先利用短距離大地線進(jìn)行計(jì)算，已知大地線上起點(diǎn)和終點(diǎn)的大地坐標(biāo)分別為：

(41)

選擇克拉索夫斯基橢球參數(shù)，按照白塞爾大地主題解算方法反算可得大地線長(zhǎng)度和正反大地方位角為：

(42)

然后利用大地線極點(diǎn)的歸化緯度求解正反大地方位角，首先根據(jù)式(40)求出大地線極點(diǎn)歸化緯度的迭代初值，計(jì)算初值所需數(shù)據(jù)和初值計(jì)算結(jié)果如表1所示。

表1中W為輔助函數(shù)，表達(dá)式為[12]：

(43)

式中：橢球第一偏心率e和白塞爾大地主題反算一致，均采用克拉索夫斯基橢球參數(shù)。

結(jié)合輔助函數(shù)，由文獻(xiàn)[3]，可得到大地緯度B和歸化緯度u的關(guān)系式為：

(44)

由式(44)可得表1中兩點(diǎn)歸化緯度的正弦和余弦值，結(jié)合表中所有數(shù)據(jù)，得到迭代初值。

得到迭代初值后，將初值代入式(39)，利用Wolfram Mathematica編程，進(jìn)行迭代計(jì)算，迭代結(jié)果如表2所示。

表1 大地線極點(diǎn)歸化緯度的迭代初值

表2 大地線極點(diǎn)歸化緯度的迭代結(jié)果

由表中迭代計(jì)算結(jié)果可以看出迭代計(jì)算4次，精度已經(jīng)足夠高，用第4次迭代結(jié)果結(jié)合表1數(shù)據(jù)，代入式(3)得：

(45)

根據(jù)式(45)，可得此大地線正反方位角為：

(46)

式(46)結(jié)果與式(42)白塞爾大地主題解算結(jié)果一致，大地線極點(diǎn)歸化緯度求解正確，說明短距離適用。接下來對(duì)長(zhǎng)距離進(jìn)行驗(yàn)證，與短距離同理，已知大地線上起點(diǎn)和終點(diǎn)的大地坐標(biāo)分別為：

(47)

同理，選擇克拉索夫斯基橢球參數(shù)，按照白塞爾大地主題解算方法反算可得大地線長(zhǎng)度和正反大地方位角為：

(48)

與短距離計(jì)算極點(diǎn)歸化緯度過程相同，過程數(shù)據(jù)和迭代初值計(jì)算結(jié)果如表3所示。

表3各數(shù)據(jù)計(jì)算方法與表1同理，得到迭代初值后，結(jié)合式(39)，利用Wolfram Mathematica編程，進(jìn)行迭代計(jì)算，迭代結(jié)果如表4所示。

表3 大地線極點(diǎn)歸化緯度的迭代初值

表4 大地線極點(diǎn)歸化緯度的迭代結(jié)果

與短距離迭代結(jié)果類似，收斂速度快，再次印證公式合理性。同理，利用第四次迭代結(jié)果計(jì)算可得：

(49)

根據(jù)式(49)，利用反正弦函數(shù)計(jì)算結(jié)果后化為度分秒，結(jié)合象限判斷，可以得到大地線正反方位角為：

(50)

式中結(jié)果與式(48)白塞爾大地主題解算結(jié)果相同，證明大地線極點(diǎn)歸化緯度求解正確，說明長(zhǎng)距離適用。

7 結(jié)束語

利用大地線克萊勞定理和白塞爾微分方程推演得到了大地線在橢球面和白塞爾球面經(jīng)差的函數(shù)關(guān)系式，得到兩球面之間經(jīng)度縮量之差的具體表達(dá)式，最后結(jié)合三角函數(shù)巧妙地分離了大地線流動(dòng)點(diǎn)和最高點(diǎn)的歸化緯度，得到大地線極點(diǎn)歸化緯度的迭代求解式，代入實(shí)際數(shù)據(jù)計(jì)算迭代結(jié)果，將其運(yùn)用在大地主題反算中大地線短距離(不跨越赤道)和長(zhǎng)距離(跨越赤道)正反方位角的求解，驗(yàn)證了大地線極點(diǎn)歸化緯度迭代求解法的可靠性。在橢球測(cè)量學(xué)中，根據(jù)已知大地線極點(diǎn)歸化緯度，可為大地主題直接解算提供新思路，并且在大地線新的理論研究中，大地線映射拓?fù)渌鶚?gòu)建的新型橢球的定位等問題也需要求出大地線極點(diǎn)的歸化緯度，因此在目前大地線新的理論研究中，大地線極點(diǎn)的歸化緯度求解有著承上啟下的重要意義。