GPS廣播星歷參數擬合的雅可比矩陣數值導數計算方法①

邢志成，王解先

（同濟大學 土木工程學院 測量與國土信息工程系，上海200092）

0 引 言

用戶利用GPS廣播星歷參數可以計算出GPS衛星的位置，繼而解算用戶的位置；反過來，也可以利用GPS衛星的位置，擬合出廣播星歷的參數［1－2]。GPS廣播星歷就是由地面控制中心在確定衛星軌道的基礎上，以衛星的位置作為觀測量擬合出廣播星歷參數再上傳給衛星，然后，衛星再轉發這些星歷參數給用戶［3]。

利用最小二乘原理擬合廣播星歷參數的過程中，偏導數矩陣的求解有解析法和數值導數法兩種方式。前者需要推導位置矢量對待估參數的偏導數，形式復雜，若改變參數形式，偏導數的形式需要重新推導；與前者相比，數值導數計算方法形式簡單，計算方便，便于計算機操作［4]。采用數值導數計算方法，利用IGS精密星歷衛星坐標數據，根據最小二乘平差原理擬合出GPS廣播星歷的15個參數；并分析了計算參數的雅可比矩陣時，增量δ的選取對擬合精度和效率的影響，得出了一些結論。

1 GPS廣播星歷參數擬合算法

GPS廣播星歷參數以開普勒軌道根數為基礎，加上表示軌道攝動調和系數以及隨時間的變化率，一共有15個：6個開普勒軌道根數ω，Ω0，M）和9個攝動力參數（Δn，i，Ω，Cuc，Cus，Crc，Crs，Cic，Cis）.15個參數的具體意義參考文獻［2].

1.1 由廣播星歷計算衛星位置

若要計算時刻t衛星的空間坐標，根據廣播星歷參數，按如下步驟計算：

1）求長半軸A

2）計算平角速度n0

3）計算從需要時刻到參考時刻的時間差tk

4）改正平角速度n

5）計算平近點角Mk

6）按下式迭代計算偏近點角Ek

7）由下式計算真近點角

8）計算緯度參數φk

9）周期改正項

10）計算改正后的緯度參數uk

11）計算改正后的向徑rk

12）計算改正后的傾角ik

13）計算衛星在軌道平面內的坐標（xk，yk）

14）改正升交點的經度

15）最后計算衛星在地固坐標系中的坐標Rk＝（XkYkZk）

1.2 數值導數的概念及雅可比矩陣的確定

若y＝f（x1，x2，…，xi…），y到xi的導數的定義為［5]

當δ足夠小時，這一導數可以近似寫成數值導數的形式

結合1.1節，則一組衛星坐標Rk＝（XkYkZk）到Xi（i＝1，2，…15）的雅可比矩陣可表示為

若利用同一顆衛星的上下各5個歷元共11組坐標，則可以得到式（4）的參數的雅可比矩陣B.

1.3 平差模型

對于給定的參考時刻t，設參數為Xi（i＝1，2，…15），由一個歷元的衛星坐標Rk＝［XkYkZk]T，就可以根據式（1）建立如下的誤差方程

式（6）中的R0代表R的近似值，X0代表X的近似值，F（X0）＝R0.

則（6）式可以寫為

式（7）中L＝（R－R0）（根據最小二乘原理，利用間接平差求解，可得

通過式（5），對每個歷元的衛星坐標可列出3個誤差方程，若選取11個歷元，那么一共有33個方程，對于15個參數，進行間接平差迭代求解。迭代的收斂條件為：｜σi＋1－σi｜／σi＜ξ，其中，ξ是預先給定的閾值，取1.0e－05.是第i次迭代的單位權方差［2]。

2 算例分析

使用2010年1月1日從0時45分到3時15分第02顆GPS衛星的共11組IGS精密星歷數據（間隔15min）。15參數中六個軌道根數的初值如表2所示，其余的9個參數初值設為0.

表2 六個軌道根數初值

為了比較小量δi的取值對擬合精度和效率的影響，算例根據1.2節所述的方法，利用公式Δ＝f（x1，x2，…，xi＋δ…）－f（x1，x2…，xi…）得到的差值Δ所在的量級不同，采用12種方案確定增量δi的取值，分別為：Δ量級為千萬米級、百萬米級、十萬米級、萬米級、千米級；百米級、十米級、米級、分米級、厘米級、毫米級、微米級。當然，增量的取值對計算結果的影響和計算機的性能有很大關系，本算例是在32位操作系統下完成的。根據這個方法可以確定增量δi的取值范圍。表3列出了方案8米級的增量δi的取值。

表3 第8種方案“米級”增量δi的值

根據不同的方案，采用1.3節所述的算法進行計算。計算迭代的次數統計如表4所示。

表4 12種方案統計

表5給出了方案8擬合出的廣播星歷15參數與GPS廣播星歷15參數的差值。

表5 擬合參數與廣播星歷參數之差Δ

圖1給出了利用方案8擬合出來的廣播星歷參數計算11顆衛星的坐標與IGS精密星歷差值。

圖1 方案8米級小量計算坐標殘差

3 結 論

1）通過表4和圖1可知，在最小二乘原理的基礎上，選取適當的小量，利用數值導數計算方法可以有效地擬合出GPS廣播星歷，其擬合結果GPS位置與IGS精密星歷坐標僅在厘米范圍有差別。

2）當小量對衛星位置的影響為萬米到厘米時，可以擬合出廣播星歷15參數，擬合中迭代計算的次數較少，擬合出的參數與GPS廣播星歷給出的參數值基本一致。當小量超出上述范圍時，擬合精度和效率降低，甚至有可能導致擬合失敗。

3）采用方案4至9對擬合的效率和精度都沒有明顯的差異，說明方案4到方案9都是比較理想的預定小量的方法。

［1]郝金明.利用地面測軌資料擬合GPS廣播星歷［D].鄭州：解放軍測繪學院，1989.

［2]呂志偉，易維勇，曾志林.GPS廣播星歷參數擬合算法及其分析［J].測繪科學技術學報，2010，27（2）：83－86.

［3]陳劉成，韓春好，陳金平.廣播星歷參數擬合算法研究［J].北京：測繪科學，2007，32（3）：12－14.

［4]王解先，徐志京.三種坐標間轉換的雅克比矩陣數值導數計算方法［J].大地測量與地球動力學，2004，24（4）：19－23.