基于迭代誤差補償的發射點坐標精確求解方法

溫永祿，馬 駿，黃 蓓，李 強，許 諾

（中國運載火箭技術研究院，北京 100071）

0 引言

隨著現代化戰爭的動態化和復雜化發展，軍事信息系統將是提升現代化戰爭下指揮員指揮能力、發揮武器裝備作戰效能的重要抓手［1-2］，而軍事信息系統經常涉及到發射點精確計算問題，該問題屬于CGCS2000 坐標系下的大地主題解算。對于小射程的導彈武器，目前通常采用6°帶高斯平面坐標，并采用高斯平面坐標計算距離和方向等參數，但是由于高斯平面坐標存在分帶處理，對于大射程的導彈武器，在計算距離和方向等參數時，與真實情況存在較大偏差［3-4］。因此，為適應新形勢下的作戰和裝備發展需求，通常采用大地主題解算算法，解決不同射程下的發射點坐標計算問題。

目前針對橢球模型下的大地主題解算問題，還沒有可用于直接解算橢球面三角形的解析求解方法。國內外眾多學者一直關注該問題并提出了多種解算方法，分別適用于短距離、中距離和任何距離。對于短距離（小于120 km），比較典型是史賴伯（O. Schreiber）公式和高斯平均引數公式［5-6］，該公式是在勒讓德（Legendre）的基礎上，利用大地線終點在起始子午圈上的垂足，用以解算大地主題。此后，高斯（Gauss）提出把勒讓德技術改化成大地線起點和終點的平均緯度和平均方位角為引數的冪級數，這樣勒讓德級數式中所有的偶次冪都消失了，可使項數減少，從而大大加速了級數的收斂，形成高斯平均引數公式，最終達到大地主題問題解算的目的，但是其解算精度與距離有關，距離越長，收斂越慢，因此，只適用于較短的距離。對于中距離（小于400 km）大地主題解算，常用的是巴烏曼投影公式，巴烏曼采用橢球面對球面的等距離投影方法求解大地主題［7］。針對中長距離解算問題，貝塞爾（Bessel）提出將橢球面上的參數按一定條件投影到輔助球面上，然后由邊長和經差投影公式求出輔助球面參數表示的大地要素［8］，實現從球面到橢球面上的相應轉換。貝塞爾大地投影的方法，解算的距離相比高斯平均引數法要遠得多，依據白塞爾的這種解法，派生出許許多多的公式［9-10］，有的是逐漸趨近的解法，有的是直接解法，還有的是其簡化公式，例如導航使用的大地線長近似計算公式，白塞爾公式的直接解法——陳俊勇公式，保持緯度不變的大地投影——張志新公式［11］。韋森特（T.Vincenty）以貝塞爾公式為基礎，推導出嵌套系數公式［12-14］。1985年，張學廉以該公式為基礎，采用嵌層約化的方法對冪級數及系數進行改化，推導出3 個嵌套系數和兩個迭代改正數，進而求解距離和經差改正項，該方法統稱為韋森特公式，也稱為嵌套系數法，可適用于任何距離的大地主題解算，而且計算量小。韋森特公式能夠適用于線長從數厘米到近20 000 km的大地問題解算，大地線長的精度為毫米級，大地方位角的精度為0.000 1 s，并且我國于2008 年頒布了GJB6304-2008《2000 中國大地測量系統》，推廣了韋森特公式的使用。

針對軍事信息系統火力籌劃運用需求，為實現對目標的精確打擊，通常需要預先根據導彈武器的射程、發射點方位角（射向）、目標點坐標，精確計算發射點坐標，供載機航線規劃和作戰規劃使用，該問題是一種非典型的大地主題解算問題，目前尚未形成統一的計算方法解決該問題。因此，針對該問題，本文在大地主題正解的基礎上，提出一種基于迭代誤差補償的非典型大地主題解算方法，根據已知的發射點方位角、目標點坐標和射程，初步計算目標點方位角，采用射向迭代法，再調用大地主題正解模型求得新的發射點方位角，求取兩次發射點方位角之間的偏差，利用偏差對目標點方位角進行補償修正，進一步減小發射點方位角誤差，直至計算得到的發射點方位角與已知的發射點方位角誤差小于預設的極小值，從而解決發射點坐標精確計算問題，該方法計算適用射程范圍廣、精度高、計算量小，能夠滿足軍事信息系統工程化運用需求。

1 典型大地主題解算方法

假設發射點大地經緯度為B1，L1，發射點至目標點的大地方位角（射向）為A12和射程S，目標點大地經緯度為B2，L2，目標點至發射點的大地方位角為A21，各角度及距離示意圖如圖1 所示。

圖1 大地問題解算示意圖

CGCS2000 地球橢球參數表如表1 所示：

表1 CGCS2000 地球橢球參數表

1.1 大地主題正解

大地主題正解問題是根據發射點的大地經緯度B1、L1和發射點至目標點的大地方位角A12和射程S，推求目標點大地經緯度B2、L2和目標點至發射點的大地方位角A21。GJB6304-2008《2000 中國大地測量系統》［12］給出的標準算法計算步驟如下：

發射點（B1、L1）的地心緯度u1正切三角函數：

從赤道到發射點連線在球上的角距σ1的正切三角函數：

從發射點到目標點連線在球上的角距σ 近似值：

依次迭代計算如下3 個公式：

直至σ 無顯著變化之后（變化值小于MIN=10-10），計算：

目標點（B2、L2）的大地緯度B2正切三角函數：

發射點（B1、L1）和目標點（B2、L2）在輔助球上的經差正切三角函數：

發射點（B1、L1）和目標點（B2、L2）的經差ΔL

目標點（B2、L2）的大地經度L2：

目標點（B2、L2）的大地方位角A21的正切三角函數：

該模型同樣適用于已知目標點坐標、射程、目標點方位角，計算發射點坐標和發射點方位角的情況。

1.2 大地主題反解

大地主題反解問題是根據已知發射點、目標點的大地經緯度B1、L1、B2、L2，求兩點之間射程S 和正反大地方位角A12、A21。GJB6304-2008《2000 中國大地測量系統》［5］給出的標準算法如下。

發射點（B1、L1）的地心緯度正切三角函數：

目標點（B2、L2）的地心緯度正切三角函數：

依次迭代計算如下公式：

發射點（B1、L1）和目標點（B2、L2）之間的射程：

發射點（B1、L1）的大地方位角A12的正切三角函數：

目標點（B2、L2）的大地方位角A21的正切三角函數：

當兩點接近對趾時，反解公式可能無解。（對趾問題是指兩點的直線連線穿過參考橢球球心的情況，主要指射程在20 000 000 m 左右，具體解決方法見文獻［8］）。

2 基于迭代誤差補償的非典型大地主題解算方法

根據射程、射向（發射點方位角）和目標點坐標，計算發射點坐標、目標點方位角，該問題是非典型的大地正解問題［13］，已知目標點大地經緯度B2、L2和發射點至目標點的大地方位角A12和射程S，推求發射點大地經緯度B1、L1和目標點至發射點的大地方位角A21。具體求解流程如圖2 所示。

1）首先，根據A12，估算得出初步的A'21；

2）根據初步計算的A'21、目標點大地緯經度B2、L2、射程S，調用大地主題正算模型，計算發射點大地方位角A'12；

圖2 非典型的大地正解計算流程

4）采用迭代法進行求解，循環計算以下步驟，直至ΔA 小于給定的極小誤差值；

②根據更新后的A'21、輸入的目標點大地緯經度B2、L2和射程S，調用大地主題正算模型，計算發射點大地方位角A'12；

④判斷ΔA' 是否小于給定的極小誤差值（MIN=10-10）；如果滿足，則退出迭代；如果迭代次數超過最大迭代次數（TIMES=400），則逐步放大迭代誤差要求（乘以10 倍，MIN=10-9，10-8，…，10-3）直至迭代誤差不大于10-3；

5）根據最終計算得到的A'21、輸入的目標點大地緯經度B2、L2和射程S，調用大地主題正算［14］，得出發射點緯經度B1、L1。

3 算例分析

為驗證本文提出方法的正確性和適用性，在Win7 i7-6700 3.4 GHz Matlab 2012b 的環境下進行仿真測試，極小值取值10-10，計算時間約0.1 s。為驗證算法的正確性，選用國軍標GJB6304-2008《2000中國大地測量系統》中的算例進行計算對比；為驗證算法的適用性，在20 000 000 m 的射程范圍內，隨機在全球范圍內生成發射點和目標點坐標，形成計算樣例，對計算結果的適用性進行分析。

3.1 與GJB 中的算例對比

根據GJB6304-2008《CGCS2000 中國大地測量系統》中提供的解算算例，以目標點經度L2、目標點緯度B2、射程S、發射點方位角A12為輸入計算發射點經度L1、發射點緯度B1和目標點方位角A21，以標準中的算例結果為參考對計算結果進行對比分析。

通過與標準算例中的結果進行對比，計算偏差如圖3 所示，計算得到的發射點坐標緯度最大偏差為0.000 333''，經度最大偏差為0.000 121''，發射方位角偏差最大為0.000 196''，滿足發射點坐標的計算精度要求。

表2 仿真計算結果

圖3 與標準算例的計算偏差

3.2 隨機測試算例

在目標點坐標緯度范圍［-90，90］°，經度范圍［0，360）°，射程［0，20 000 000］m，射向［0，360）°范圍內，隨機生成測試樣例，測試次數設置為5 000 次，對結果進行統計與分析。首先，根據目標點坐標、射程、射向（發射點方位角），采用非典型的大地主題正解模型計算發射點坐標、目標點方位角［15］；然后再根據得到的發射點坐標，結合目標點坐標，采用國軍標中的大地主題反解公式計算射程、發射點和目標點方位角，并與已知射程、發射點方位角和非典型大地主題解算計算得到的發射點坐標進行對比，分析計算誤差如下頁圖4～圖6 所示。

根據仿真計算結果可看出，在設定的數據取值范圍內，隨機生成的5 000 個測試樣例中，射程計算偏差絕對值最大為6.08×10-3m，偏差量級為10-3m；發射點方位角偏差絕對值最大為7.23×10-7°m，偏差量級為10-7°（10-4''），目標點方位角偏差絕對值最大為7.78×10-7°，偏差量級為10-7°（10-4''），計算精度能夠滿足發射點精確計算要求，而且算法具有很好的魯棒性。

圖4 目標點坐標分布

圖5 射程偏差

圖6 目標點方位角偏差

圖7 發射點方位角偏差

4 結論

本文通過分析發射點精確計算特點，在大地主題正算（韋森特公式）、大地主題反算求解方法的基礎上，提出一種基于迭代誤差補償的非典型大地主題正解方法，采用國軍標GJB6304-2008 中的算例和隨機測試樣例，對算法的適用性和精確性進行驗證，通過仿真結果進行分析可知，本文提出的方法適用于20 000 km 內任意兩點的發射點精確計算，方法簡單實用，能夠滿足軍事指揮控制信息系統的發射點坐標精確計算要求。