宋傳峰,秘金鐘,黨亞民,薛樹強
(中國測繪科學研究院,北京 100830)
隨著GEO (Geostationary Earth orbit)與IGSO (Inclined Geo Synchronous orbit)衛(wèi)星在導航定位領(lǐng)域中的應(yīng)用,對其軌道的精度要求越來越高[1]。利用數(shù)值算法積分衛(wèi)星軌道是人造衛(wèi)星軌道預(yù)報與精密定軌中的基礎(chǔ)工作[2]。由于二體問題可以直接得出衛(wèi)星位置、速度和加速度的解析解,因此本文避開復(fù)雜的攝動力計算,在二體意義下分別對GEO和IGSO衛(wèi)星進行了軌道數(shù)值積分,并與其理論值進行了分析比較,得出了有益結(jié)論。若考慮攝動力,則不便于結(jié)果比較。
目前,在MEO定軌中使用的積分方法有RKF與Admas預(yù)測-校正系統(tǒng)結(jié)合的方法 (如Gamit軟件),Collocation方法 (如Bernese軟件)等。本文驗證了RKF與Admas預(yù)測-校正系統(tǒng)結(jié)合的方法在GEO與IGSO衛(wèi)星軌道積分時的適用性。
Runge-Kutta-fehlberg方法 (RKF)是一種嵌套的RK方法,即同時給出n階和n+1階兩組RK計算公式,用兩組公式計算出來的xi+1之差來估計截斷誤差,根據(jù)截斷誤差的大小來控制步長[3]。以下給出7階和8階的RKF公式
截斷誤差
式中,h為積分步長,ck,βkλ為積分算法的系數(shù)。
在衛(wèi)星運動方程的數(shù)值積分中,單步法只是用于多步法的起步算法,當采用單步法推出足夠的步點后,就可采用高精度的多步法往前推算[3]。本文使用變步長的RKF7(8)階積分法積分10步后,再使用定步長的10階Adams預(yù)測-校正系統(tǒng)進行積分。以下給出顯式的Adams-Bashforth公式和隱式的Adams-Moulton公式。
顯式
隱式
在二體意義下,GEO地球同步靜止軌道衛(wèi)星在地心地固坐標系中的坐標為一固定值,但軌道積分是在慣性坐標系下進行的,因此需要計算二體問題下地固坐標系與慣性坐標系間的相互轉(zhuǎn)換矩陣。
假設(shè)積分初始時刻t0(s),兩坐標系坐標軸完全重合,則任意t(s)時刻,地固系向慣性系轉(zhuǎn)換的矩陣為
在慣性系下,GEO衛(wèi)星繞Z軸旋轉(zhuǎn)的角速度與地球自轉(zhuǎn)角速度相等,因此易得下式
慣性系下GEO衛(wèi)星在X-Y平面內(nèi)運動,即Z≡0,由此可得:
式中,X、Y、Z表示衛(wèi)星在慣性系下的坐標。
在二體意義下,GEO衛(wèi)星的軌跡為圓,且其在地固坐標系下的速度為0,因此給定初始坐標可得任意時刻GEO衛(wèi)星在地固和慣性坐標系下的位置、速度和加速度。
IGSO衛(wèi)星軌道是特殊的地球同步圓軌道,軌道高度跟GEO一樣,均屬于地球同步軌道[4]。由于IGSO衛(wèi)星軌道在慣性系下也是圓軌道,因此可看做由GEO衛(wèi)星軌道繞X或Y軸或X-Y平面內(nèi)一條直線旋轉(zhuǎn)i角,i為IGSO衛(wèi)星軌道傾角,本文中按繞X軸旋轉(zhuǎn)i角計算。
若t0時刻,IGSO衛(wèi)星在其軌道直角坐標系(慣性系繞X軸旋轉(zhuǎn)i角)下坐標為(,,0)T,速度為(-ωe,ωe,0)T,則t時刻其在軌道直角坐標系中坐標為(X*,Y*,0)T=Rz(α),,0)T,速度為(-ωeY*,ωeX*,0)T。ωe,α與3.1節(jié)的意義相同。
t時刻IGSO衛(wèi)星在地心慣性坐標系下位置和速度:
式中,X、Y、Z、VX、VY、VZ為衛(wèi)星在慣性坐標系下位置和速度,i表示軌道傾角,X*、Y*為衛(wèi)星在軌道直角坐標系下平面坐標。
將慣性系下位置左乘矩陣RT,R同式 (6),即可得地固系下衛(wèi)星位置,再轉(zhuǎn)成大地坐標系 (B,L,H)即可得IGSO衛(wèi)星的星下點軌跡,圖1為經(jīng)度為45°,傾角為10°和30°的IGSO衛(wèi)星星下點軌跡 (二體意義下)。
對不同經(jīng)度的GEO衛(wèi)星、不同經(jīng)度和傾角的IGSO衛(wèi)星,先使用變步長的RKF7(8)階積分法積分10步,截斷誤差設(shè)定為10-12,再使用定步長75s的10階Adams預(yù)測-校正系統(tǒng)進行積分48h(2圈),每隔15min輸出一組結(jié)果。不同經(jīng)度的GEO衛(wèi)星在地心地固坐標系下的積分初值見表1,其速度均為0;不同經(jīng)度和傾角的IGSO衛(wèi)星在慣性坐標系下的積分初值見表2。
圖1 IGSO衛(wèi)星星下點軌跡 (二體意義)
表1 GEO衛(wèi)星積分初值
表2 IGSO衛(wèi)星積分初值
由于篇幅原因,僅給出了Y方向差異的統(tǒng)計信息和經(jīng)度為30°的GEO衛(wèi)星差異變化圖。圖2是經(jīng)度為30°的GEO衛(wèi)星軌道數(shù)值積分結(jié)果與二體意義下結(jié)果在地心地固坐標系下的差異,歷元間隔為15分鐘,表3是GEO衛(wèi)星軌道積分結(jié)果與真值 (二體意義)在地固系下Y方向差異的統(tǒng)計信息。
圖2 GEO衛(wèi)星軌道數(shù)值積分差異
表3 GEO衛(wèi)星軌道數(shù)值積分Y方向差異統(tǒng)計信息
由于篇幅原因,僅給出了X方向差異的統(tǒng)計信息和經(jīng)度為45°傾角為30°的IGSO衛(wèi)星差異變化圖。圖3是經(jīng)度為45°傾角為30°的IGSO衛(wèi)星軌道數(shù)值積分結(jié)果與二體意義下結(jié)果在地心地固坐標系下的差異,歷元間隔為15min,表4是IGSO衛(wèi)星軌道積分結(jié)果與真值 (二體意義)在地固系下X方向差異的統(tǒng)計信息。
表4 IGSO衛(wèi)星軌道數(shù)值積分X方向差異統(tǒng)計信息
圖3 IGSO衛(wèi)星軌道數(shù)值積分差異
由以上圖表可以得出,二體意義下GEO、IGSO衛(wèi)星軌道數(shù)值積分結(jié)果與真值相差甚小,均在亞毫米級或以內(nèi),與現(xiàn)有精密定軌精度相比可以忽略,適當減小積分步長可獲得更小的理論值差異,但這將明顯增加計算耗時。結(jié)果表明對MEO(如GPS)衛(wèi)星適用的軌道數(shù)值積分方法對GEO和IGSO衛(wèi)星也同樣適用。
本文僅在二體意義下驗證了數(shù)值積分方法的適用性,現(xiàn)實中導航衛(wèi)星攝動力復(fù)雜多變,因此建立合適的力學模型再結(jié)合數(shù)值積分方法進行軌道積分是統(tǒng)計定軌中的重要步驟。
[1]劉吉華,歐吉坤,鐘世民,等.GEO與IGSO衛(wèi)星星間差分的精密定軌[C]//中國衛(wèi)星導航學術(shù)年會組委會.第一屆中國衛(wèi)星導航學術(shù)年會文集.北京:出版社不詳,2010:1-7.
[2]李得海,袁運斌,歐吉坤,等.12階Runge-kutta2次算法的衛(wèi)星軌道積分研究[J].武漢大學學報:信息科學版,2010,35(11):1335-1338.
[3]趙齊樂.GPS導航星座及低軌衛(wèi)星的精密定軌理論和軟件研究[D].武漢:武漢大學,2004.
[4]楊 洋,范 麗,董緒榮.IGSO星座分析與優(yōu)化設(shè)計[J].兵工自動化,2008,27(11):53-55.
[5]ES-HAGH M.Step Variable Numerical Orbit Integration of a Low Earth Orbiting Satellite[J].Journal of the Earth &Space Physics,2005,31(1):1-12.
[6]XU Guo-chang.GPS Theory Algorithms and Applications[M].Berlin:Springer-Verlag,2003.
[7]IERS Convention Centre.IERS Conventions:2010[M].Frankfurt am Main:Verlag des Bundesamts für Kartographie und Geod?sie,2010.
[8]LEICK A.GPS Satellite Surveying[M].Hoboken:John Wiley &Sons Inc.,2004.