盧 祥 袁運彬 蔣振偉

1)大地測量與地球動力學國家重點實驗室，武漢 430077

2)中國科學院大學，北京 100049

基于龍格庫塔法的GLONASS軌道仿真研究*

盧 祥1，2)袁運彬1)蔣振偉1，2)

1)大地測量與地球動力學國家重點實驗室，武漢 430077

2)中國科學院大學，北京 100049

將龍格庫塔法引入GLONASS軌道仿真，研究步長參數(shù)設(shè)置對結(jié)果的影響，定量分析仿真結(jié)果與廣播星歷的差異和算法效率，并進行長時間仿真試驗分析。結(jié)果表明，步長參數(shù)設(shè)置為1～300 s時能兼顧算法效率與結(jié)果精度;仿真結(jié)果與廣播星歷的差異隨時間延長而增大，仿真300 min時X、Y、Z方向差異分別為266.43、246.13、336.06 m，每仿真一個歷元24顆星數(shù)據(jù)耗時約2 ms。龍格庫塔法適用于GLONASS軌道仿真。

龍格庫塔法;GLONASS;軌道仿真;積分步長;積分時長

GNSS(GPS、BD、GALLEO、GLONASS)仿真系統(tǒng)在接收機測試、算法驗證等工作中有著重要作用［1］。與 GPS、BD、GALILEO 系統(tǒng)仿真情況不同，國內(nèi)外對GLONASS仿真研究較少。Spirent公司于2007年實現(xiàn)了GLONASS的系統(tǒng)仿真，但其方法及代碼均未公開。國內(nèi)開展GNSS系統(tǒng)仿真較早的單位實現(xiàn)了GPS、GALILEO、BD 3大系統(tǒng)的仿真，但鮮見公開報道GLONASS系統(tǒng)仿真研究。趙春梅等運用動力學定軌法進行了GALILEO軌道仿真，但由于GLONASS軌道采用衛(wèi)星坐標、速度、加速度方式發(fā)布，不同于其他3大系統(tǒng)，傳統(tǒng) GPS、BD、GALILEO中的動力學定軌方法不能直接應(yīng)用于GLONASS軌道仿真。傳統(tǒng)動力學定軌方法公式繁雜、計算量大，通常難以滿足 GLONASS高性能仿真的要求［2－6］。本文首次將龍格庫塔法引入到GLONASS軌道仿真，分析步長、仿真時間對仿真結(jié)果的影響，定量分析算法效率，并測試算法的穩(wěn)定性。

1 基于龍格庫塔法的GLONASS軌道仿真方法

采用簡化GLONASS軌道力學模型進行龍格庫塔積分，簡化力學模型如下:

根據(jù)以上微分方程，采用四階龍格庫塔積分對GLONASS軌道進行積分，積分模型為［7］:

根據(jù)式(7)～(11)積分模型，可按以下步驟進行積分:

由初始位置X0=X(t0=tb)計算時刻t的衛(wèi)星位置，必須以積分步長h重復(fù)積分，t1=t0+h，t2=t1+h，…，tm=tm－1+h，其中 t－ h≤tm＜ t0如果 tm=th，最后一次積分步長保持h不變;如果tm＞t－h(huán)，則最后一次積分步長h=t－tm。本文在計算過程中，將x··LS、y··LS、z··LS固定為常數(shù)。

2 仿真試驗及分析

2.1 步長參數(shù)選擇

利用龍格庫塔積分法仿真GLONASS軌道時，步長參數(shù)的選擇直接影響到仿真效率與仿真精度。較短的步長仿真精度較高，但積分效率較低;較長的步長仿真精度較差，但積分效率較高。選擇合適的積分步長可以在保證積分效率的同時達到較好的精度。為分析龍格庫塔積分方法中步長選擇對仿真結(jié)果產(chǎn)生的影響，選取2012-02-20由IGS提供的GLONASS廣播星歷進行試驗。

以廣播星歷中1號星為例進行試驗，以0:15:00 數(shù)據(jù)為起始數(shù)據(jù)(圖 1)，分別以 1、10、30、60、90、120、300、900 s為積分步長，積分 30 min，得到衛(wèi)星位置 Xn、速度值(即 0:45:00 的數(shù)據(jù))，將 Xn、與廣播星歷中0:45:00的衛(wèi)星位置、速度作差，計算得到不同積分步長下的仿真結(jié)果與廣播星歷的差異如圖2所示。

圖1 試驗示意圖Fig.1 Experimental sketch

圖2顯示，積分步長在1～300 s變化時，仿真30 min軌道數(shù)據(jù)，與真實軌道差異均保持在1 m左右;但當積分步長升至900 s時，誤差顯著增長，X、Y、Z 方向分別為 29.31、26.48、18.92 m。因此，在GLONASS軌道仿真中，宜將積分步長選擇在1～300 s之間，具體選擇還要結(jié)合軌道仿真要求的數(shù)據(jù)頻率確定。

圖2 不同步長30 min積分結(jié)果與廣播星歷差異Fig.2 Difference between integration results and broadcast ephemeris with different step length

2.2 仿真結(jié)果與廣播星歷的差異

為分析利用龍格庫塔法仿真GLONASS軌道的效果，設(shè)計試驗定量研究仿真結(jié)果與廣播星歷的差異。

圖3 試驗示意圖Fig.3 Experimental sketch

圖4 不同積分時長積分結(jié)果與廣播星歷差異Fig.4 Difference between integration results and broadcast ephemeris with different time duration

在2012-02-20 IGS發(fā)布的GLONASS廣播星歷中以1號星0:15:00為起點(圖3)，以60 s為積分步長分別積分 30 min、60 min、90 min、…、300 min，以廣播星歷中對應(yīng)時刻的軌道位置作為參考值，計算得到不同積分時長下的軌道差異如圖4所示。

隨著時間的延長，仿真結(jié)果與廣播星歷差異逐漸增大。30 min時，X、Y、Z三個方向的差異分別為1.25、0.75、1.21 m，優(yōu)于文獻［3］采用 4、5、6、7 階多項式仿真GPS系統(tǒng)軌道30 min的精度(10 m級)，與文獻［3］采用20階多項式擬合法仿真30 min GPS系統(tǒng)的精度相當(m級)。90 min仿真結(jié)果與廣播星歷的差異增至10 m級，到300 min時X、Y、Z 三個分量差異分別為 266.4、246.1、336.0 m，相比GLONASS衛(wèi)星軌道高度19 100 km，這些差異非常微小，可認為與廣播星歷基本吻合，表明龍格庫塔法適用于GLONASS軌道仿真。

2.3 算法效率

接收機進行野外測試前，需進行大量的仿真測試，特別是檢測50 Hz甚至更高頻的數(shù)據(jù)時，對軌道仿真算法效率有很高的要求。定量研究龍格庫塔法的效率，選取不同時長，積分步長設(shè)置為1 s，仿真GLONASS星座(24顆星)軌道數(shù)據(jù)(表1)。圖5為算法效率圖，龍格庫塔法每仿真一個歷元24顆星數(shù)據(jù)耗時約2 ms，可滿足50 Hz甚至更高頻率GLONASS的仿真要求。

2.4 長時間仿真穩(wěn)定性

使用龍格庫塔法進行長時間軌道仿真時，仿真結(jié)果發(fā)散性、算法穩(wěn)定性都是關(guān)鍵問題。為分析龍格庫塔法進行GLONASS軌道長時間仿真的穩(wěn)定性，以2012-02-20 00:15:00數(shù)據(jù)為起點進行GLONASS星座(24顆星)軌道仿真，得到 1 d、2 d、20 d的仿真軌道。選取1號星為例，繪制廣播星歷、仿真軌道1 d、2 d、20 d軌道軌跡圖(圖6)，將20天1號星的仿真軌道與廣播星歷作差得到20 d仿真軌道與廣播星歷X、Y、Z方向差異(圖7)。

由圖6與7可知，利用龍格庫塔法仿真所得的GLONASS軌道軌跡與廣播星歷軌跡高度相似;使用龍格庫塔法進行長時間GLONASS軌道仿真時，盡管仿真軌道形狀與廣播星歷類似，但仿真軌道與廣播星歷的差值隨時間延長逐漸增大。仿真時間到達480 h(20 d)時，X、Y、Z三個方向與廣播星歷的差異分別為 44 906.55、43 388.11、110 309.45 m。

表1 算法效率Tab.1 Efficiency of the algorithm

圖5 算法效率Fig.5 Efficiency of the algorithm

3 結(jié)論

系統(tǒng)研究了利用龍格庫塔法實現(xiàn)GLONASS軌道仿真的性能與涉及的關(guān)鍵技術(shù)，給出龍格庫塔法的數(shù)學模型。研究結(jié)果顯示:

1)步長參數(shù)設(shè)置為1～300 s時能兼顧算法效率與結(jié)果精度;

2)以60 s為積分步長，30 min仿真軌道與廣播星歷 X、Y、Z 三個方向差異分別為 1.25、0.75、1.21 m，300 min仿真軌道與廣播星歷X、Y、Z分量差異分別為 266.4、246.1、336.0 m;

3)利用龍格庫塔法每仿真一個歷元24顆星數(shù)據(jù)耗時約為2 ms，可滿足接收機仿真測試工作中常用的50 Hz高性能數(shù)據(jù)仿真要求;

4)使用龍格庫塔法仿真20 d GLONASS軌道，仿真軌道形狀與廣播星歷類似，但仿真軌道與廣播星歷的差值隨時間的延長逐漸增大，20 d時X、Y、Z三個方向與廣播星歷的差異分別為44 906.55、43 388.11、110 309.45 m。

1 趙春梅，歐吉坤，文援蘭.基于GALILEO及GPS-GALILEO組合系統(tǒng)的仿真分析［J］.系統(tǒng)仿真學報，2005，17(4):1 008 －1 011.(Zhao Chunmei，Ou Jikun，Wen Yuanlan.Simulation analysis based on the combination of GALILEO and integrated GPS-GALILEO［J］.Journal of System Simulation，2005，17(4):1 008 －1 011)

2 李得海，等.基于GALELIO/GPS仿真系統(tǒng)的單點定位研究［J］.系統(tǒng)仿真學報，2009，21(10):2 828 －2 831.(Li Dehai，et al.Research on point positioning based on simulated GALILEO and GPS system［J］.Journal of System Simulation，2009，21(10):2 828 －2 831)

3 白文娟.基于OpenGL的GPS衛(wèi)星軌道仿真與可視化實現(xiàn)［D］.哈爾濱:哈爾濱工程大學，2009.(Bai Wenjuan.Simulation and visualization realization of GPS orbit based on PpenGL［D］.Harbin:Harbin Engineering University，2009)

4 戴沖.衛(wèi)星導(dǎo)航系統(tǒng)空間段仿真關(guān)鍵技術(shù)研究［D］.長沙:國防科技大學，2011.(Dai Chong.Study on key technology of space segment of satellite navigation system［D］.Changsha:National University of Defense Technology，2011)

5 謝杰，等.GPS信號仿真器方案設(shè)計與實現(xiàn)［J］.計算機仿真，2012，29(2):36 －39.(Xie Jie，et al.Design and implementation of GPS signal simulator based on ARM［J］.Computer Simulation，2012，29(2):36 －39)

6 于清德，等.GNSS幾何性能仿真分析［J］.計算機仿真，2010，27(6):87 －91.(Yu Qingde，et al.simulation analysis of GNSS geometric performance［J］.Computer Simulation，2010，27(6):87 －91)

7 李慶楊，王能超，易大義.數(shù)值分析［M］.北京:清華大學出版社，2008.(Li Qingyang，Wang Nengchao，Yi Dayi.Numerical analysis［M］.Beijing:Tsinghua University Press，2008)

RESEARCH ON GLONASS ORBIT SIMULATION BASED ON RUNGE-KUTTA METHOD

Lu Xiang1，2)，Yuan Yunbin1)and Jiang Zhenwei1，2)
1)State Key Laboratory of Geodesy and Earth’s Dynamics，Wuhan 430077
2)University of Chinese Academy of Sciences，Beijing100049

Runge-Kutta method was introduced to GLONASS orbit simulation to analyze how the step length parameter effects on simulated results，the efficiency of Runge-Kutta method and the difference between simulated results and broadcast orbit.The result shows that step the balance between algorithm efficiency and accuracy can be keeped though setting step length as 1 s－300 s.The difference between simulated results and broadcast ephemeris increases with time，and reaches about 250m in X，Y，Z directions for 300 minutes simulation time duration.It takes 2 ms to simulate 1 epoch 24 satellites orbit data.According the research，Runge-Kutta method can be used for GLONASS orbit simulation.

Runge-Kutta;GLONASS;orbit simulation;integration step length;integration time duration

P228

A

1671-5942(2014)03-0137-05

2013-10-13

國家自然科學基金重點項目(41231064);國家863計劃項目(2012AA121803)。

盧祥，男，1989年生，碩士生，主要研究方向為GNSS精密定位。E－mail:luxiang.whigg@gmail.com。