楊志濤
(1.中國科學院國家天文臺,北京100012;2.中國科學院大學,北京100049;3.國家航天局空間碎片監(jiān)測與應用中心,北京100012)
初軌確定已有200多年的歷史,最初是針對小行星、彗星等自然天體,隨著天體力學的發(fā)展,特別是近代觀測技術的發(fā)展和計算條件的改善,初軌計算方法也有了進一步發(fā)展。人造地球天體上天后,跟蹤測量方式有所改變,初軌計算方法又根據(jù)衛(wèi)星運動規(guī)律進行諸多改進。但無論何種改進,基本上歸結為采用的動力學條件和表現(xiàn)形式的不同,而初軌確定的基本原理并無改變。雖然方法眾多,但從本質上講,初軌計算方法基本上可以分為兩大類[1-13]:一類以Laplace方法為代表,通過計算得到t0時刻的位置矢量r0和速度矢量,然后轉換成相應的軌道根數(shù);另一類以Gauss方法為代表,通常分為兩步,第一步通過計算得到兩個時刻的位置矢量r1和r2,第二步由r1和r2計算軌道根數(shù),實即蘭伯特 (Lambert)問題。
盡管初軌確定方法眾多,但其基本前提是一致的:建立在無任何先驗信息前提下的定軌方法[14]。短弧初軌確定的初衷就是在第一次獲得該“新”目標的短弧資料時,在沒有任何初始信息的前提下確定其軌道。在航天器發(fā)射初期和中途修正、對非合作目標和空間碎片的搜索發(fā)現(xiàn)或者編目中斷后的重新編目等工作中,初軌確定方法都有重要的應用價值。在具體的工程應用中,初軌確定的定軌誤差與資料誤差的關系是需要重點關注的一個問題。
本文首先介紹了改進的Laplace方法初軌確定的原理,然后利用仿真觀測數(shù)據(jù)對定軌誤差與觀測資料時長、資料誤差及資料密度等因素的關系進行了定量分析,最后總結給出一些有益的分析結論,期望能為相關工作提供技術參考。
在地心天球坐標系中,航天器相對中心天體(質心)的運動微分方程即:
式中,GM是中心天體引力常數(shù);加速度Fε表示中心天體質心引力以外的各種攝動力。在采用短弧情況下 (這也正對應初軌確定的實際背景)存在如下時間 (間隔)冪級數(shù)解:
式中,x,y分量的系數(shù)相同,在方程中即記作F,G,而z分量就是方程中的Fz,Gz。若記τ=Δt表示時間段長度,則上述系數(shù)F,G和Fz,Gz按τ的冪級數(shù)展開式如下:
護理床是行動不方便的病人在住院或居家護理時使用的病床.其主要目的是便于護理人員進行照顧,使病人康復.但是傳統(tǒng)護理床價格昂貴,而且功能模式單一,或是只能輔助翻身,或是只能輔助坐臥,而且病人一旦需要離開床體進行其他活動時,搬移病人便成了護理人員的一大難題[1-3].
雷達是空間碎片監(jiān)測中的常用設備,每次采樣可同時獲得測距和測角資料(t;ρ,A,h)。對于合作目標的監(jiān)測,可由衛(wèi)星導航系統(tǒng)給出每個采樣時刻的定位資料t,r(x,y,z)。這兩類資料在定軌方程中本質上相同,由測距和測角資料(t;ρ,A,h)結合測站坐標可轉換得到同一時刻定位資料t,r(x,y,z), 相應的定軌基本方程可直接表達為下列簡單形式:
對于一次測量采樣資料,上式對應的三個方程都是獨立的,原則上只需要兩次采樣即可定軌,給出定軌歷元t0時刻的位置r0(x0,y0,z0)和速度
關于上述線性方程組的求解,對于測角加測距或定位資料盡管只要兩次采樣就可定軌,但通常都是盡量利用多組資料定軌。k次采樣數(shù)據(jù)對應的基本方程個數(shù)是k×3,且往往有k×3≥n的狀態(tài),這里n=6是軌道量的維數(shù)。為了充分發(fā)揮多資料的作用,就需要選擇相應的最優(yōu)估計方法來求解定軌基本方程。但對于初軌確定而言,沒有必要采用過于復雜的估計方法,通常就采用簡單的等權最小二乘估計。
k組數(shù)據(jù)對應的基本方程組可整理為如下形式,其中
式中, 系數(shù)矩陣A3k×6和B3k×1由式 (5) 中的F,G,F(xiàn)z,Gz等量組成。采用等權最小二乘估計方法可得 到:Z=[(A3k×6)TA3k×6]-1(A3k×6)TB3k×1。由于系數(shù)矩陣的計算亦會涉及用到變量Z,因此實際計算是一個逐漸迭代收斂的過程。F,G等量可采用如下初值在定軌中進行迭代,并且這是一個較好的初始信息。
本節(jié)分析所用的初軌確定算法,即上文所述考慮J2項影響的改進的Laplace方法;而分析采用的觀測資料默認為等時間間隔的多組數(shù)值仿真的定位數(shù)據(jù)。由于初軌計算通常是半長徑的誤差Δa最大,因此本文以Δa的大小作為衡量定軌結果優(yōu)劣的標準,下文所述的定軌誤差均默認是指Δa。
除定軌模型外,初軌確定的誤差還和所用觀測資料的情況密切相關,例如資料時長、資料密度以及資料誤差的大小。有必要說明,此處的資料時長是指所用資料的最大時間間隔,以此表征資料弧段的大小;資料密度表征所用資料的稀疏程度,由該等時間間隔的長度表示。由上節(jié)所述原理可知,此初軌確定方法的核心是在動力學關系中采用了時間冪級數(shù)解 (這與通常的軌道分析解中采用的小參數(shù)冪級數(shù)解有所不同)的形式,表達式的收斂性及精度主要對時長τ(等效于弧長)比較敏感,這與軌道半長軸密切相關,而對其他軌道參數(shù)并不敏感,因此本文仿真分析的低地球軌道只涉及不同軌道高度的情況,并未對軌道偏心率和軌道傾角等參數(shù)增加算例。
本節(jié)主要針對定位資料,分析對低軌不同高度下近圓軌道的定軌誤差與資料時長、資料密度以及資料誤差大小的關系。分析所用的觀測資料由高精度數(shù)值法軌道預報軟件模擬計算給出:軌道外推的起始歷元為2007年7月1日12時0分0秒,終止歷元為2007年7月2日12時0分0秒,時間間隔為1s,模擬的軌道共有4組,軌道初值如表1所示;計算時考慮的攝動力模型包括:非球形引力 (20×20階次)、日月第三體質點引力、大氣阻尼、光壓和日月潮汐形變攝動等,其中大氣阻尼和光壓攝動的計算中所用的衛(wèi)星有效面質比均為0.02m2/kg,大氣模型為Jacchia-Roberts。
在下文分析中對觀測數(shù)據(jù)加入的隨機誤差,采用均值為0、標準差大小為σ的高斯分布的隨機數(shù),其中σ會給定不同的數(shù)值,并以此表征誤差大小。
定軌結果與所用觀測資料的弧段大小,即資料時長ΔT密切相關。此處取定軌歷元為模擬資料起始歷元開始的第500s,即2007年7月1日12時8分19秒;選取數(shù)據(jù)時固定資料密度為每5s取一點,同時保證定軌歷元為所選弧段的中間時刻,資料時長從10s(即3點)開始遞增,至定軌誤差為Δa=10km時停止,比較的標準是模擬資料中定軌歷元的軌道根數(shù)。對四條軌道均用上述方法計算,得到定軌誤差與資料時長的變化關系如圖1所示:
表1 模擬資料軌道初值列表Tab.1 Initial orbital values of simulated data
圖1 定軌誤差與資料時長關系圖Fig.1 The relationship between orbit determination error and data duration
由圖1可以看出,定軌誤差為Δa隨資料時長ΔT變化的總體趨勢相同:在ΔT較小 (100s內)時Δa會有小幅波動,隨著ΔT的遞增,Δa先會有一段平穩(wěn)期 (此時誤差在200m之內),之后Δa會隨之緩慢遞增至一峰值,然后開始快速遞減,當Δa減小至負值后,實際上定軌誤差 (Δa的絕對值)是在快速遞增。
顯然,誤差峰值及其對應的峰值時長是一對很有意義的量,且大小與軌道高度相關,現(xiàn)將二者分別記為和,其中h(單位默認為km)代表軌道高度。由上圖可見:誤差峰值及其對應的峰值時長均隨軌道高度遞增。
同時,考慮到定軌結果的實用性,當資料時長ΔT很大時,即使迭代計算可以收斂并給出定軌結果,仍有可能因為定軌誤差太大而導致定軌結果無意義。為此引進有效時長的概念:
假定當Δa大于某閾值Da時定軌結果即視為無效,由于Δa隨ΔT的增大而增大,當Δa達到閾值Da時對應的資料時長即稱為有效時長,并用TD表示。本文中取Da=10 km。
表2 誤差峰值及對應的資料時長Tab.2 Peak error and time length of corresponding data
以上結果均是由無資料誤差 (嚴格來講含有輸出數(shù)據(jù)時帶來的截斷誤差,其大小為1mm量級,在標準單位下是1.6×10-10)的數(shù)據(jù)計算得到,因此該定軌誤差主要反映的是模型誤差。對于峰值時長之后快速變化,顯然與時間冪級數(shù)展開式F,G的變化規(guī)律一致:計算F,G時舍棄的余項隨ΔT的增長以指數(shù)形式遞增,由此帶來的模型誤差亦會快速增大。詳細的誤差演化規(guī)律仍有待進一步研究。
表2中有效時長與軌道高度的關系如圖2所示,二者近似為線性關系,以此可作為選取定軌資料長度的簡易判據(jù)。
圖2 有效時長與軌道高度的關系Fig.2 The relationship between effective time length and orbital height
在3.1節(jié)的分析中,所用定軌數(shù)據(jù)均未加入誤差,但在實際的測量數(shù)據(jù)中隨機誤差是不可避免的。若在計算中對所選數(shù)據(jù)人為加入隨機誤差,對應的定軌結果必然會有不同。記資料誤差σ引起的定軌誤差為δaσ,a1、a2分別表示未加入、加入隨機誤差時半長徑的定軌結果,則δaσ=a2-a1。圖3是對300km高的軌道資料加入σ=20m的隨機誤差后,δaσ與ΔT的變化關系。
圖3 資料誤差引入的定軌誤差與資料時長的關系Fig.3 The relationship between orbit determination error caused by data errors and time length of data
圖4 300km高軌道δaσ隨ΔT變化關系Fig.4 Variation of δaσ with ΔT for 300km high orbit
由于δaσ是漸變趨于零,故嚴格定義理想時長Tσ為:使得δaσ小于閾值δa的最小資料時長。
圖5是取閾值δa=200m時對應的理想時長Tσ與資料誤差大小σ的關系圖,其中圖a對應資料間隔為5s,σ的步長為5m;圖b對應資料間隔為1s,σ的步長為2m的情況。定軌時刻均為所選弧段的中間時刻??梢钥闯鯰σ與σ近似成線性關系,對軌道高度的變化并不敏感,只在σ較大時有緩慢的遞增趨勢;同時通過兩圖中Tσ大小的比較可以看出,Tσ與資料間隔密切相關。
定軌誤差的大小與資料密度密切相關。此處以300km高的軌道為例分析定軌誤差δaσ與資料密度的關系,其中資料密度用相鄰數(shù)據(jù)的時間間隔表征。
圖5 理想時長Tσ與資料誤差大小σ關系(a圖:資料間隔為5s;b圖:資料間隔為1s)Fig.5 The relationship between ideal time Tσand data error size σ.(a: data interval is 5s;b:data interval is 1s)
首先計算在不加入隨機誤差的情況下,采用不同時間間隔的資料進行初軌確定,對應的定軌誤差與資料時長的關系如圖6所示。
圖6 無誤差300km高軌道不同資料密度定軌誤差圖Fig.6 Orbital determination error map of different data densities for error-free 300km high orbit
圖7 100s資料定軌誤差與時間間隔的關系Fig.7 The relationship between orbital determination error of 100s data and time intervals
圖8 200s資料定軌誤差與時間間隔的關系Fig.8 The relationship between orbital determination error of 200s data and time intervals
圖9 400s資料定軌誤差與時間間隔的關系Fig.9 The relationship between orbital determination error of 400s data and time intervals
由圖10、11可見,在總體趨勢上,δaσ隨δT的增大而增大,因為所用資料個數(shù)在隨之遞減。對于圖中跳變及遞減的部分,實際反映的是在資料組數(shù)固定的情況下,定軌誤差隨資料間隔的增大在減小 (詳見圖10、11)。如圖9中的最后一次跳變,對應的 δT為133s,之后直至 δT達到200s的區(qū)間之內,δaσ是隨之遞減的,分析可知,由于資料組數(shù)必為整數(shù),在最大資料時長固定為400s時,δT在區(qū)間[134s,200s]內時均只有3組數(shù)據(jù)。同時,資料誤差σ越大,定軌誤差δaσ對資料密度越敏感。
圖10、11依次是固定選用資料組數(shù)N為2組和5組,同時加入不同大小的隨機誤差時,定軌誤差隨時間間隔δT的變化關系圖。其中縱軸為統(tǒng)計計算的定軌誤差δaσ,計算方法同3.3中所述;橫軸為相鄰資料的時間間隔δT,由于選用資料組數(shù)已固定,因此δT越大對應的資料時長ΔT=NδT亦越大。計算時取δT的步長為2s,定軌歷元 (記為Tm)同3.1,選取資料使得Tm恰為其中間時刻。
圖10 2組資料定軌誤差與資料間隔關系Fig.10 The relationship between orbit determination errors of 2 groups of data and data intervals
圖11 5組資料定軌誤差與資料間隔關系Fig.11 The relationship between orbit determination errors of 5 groups of data and data intervals
由圖可見,在資料組數(shù)固定時,對于不同大小的隨機誤差,其定軌誤差隨時間間隔δT的變化趨勢是一致的,均會有一拐點。記拐點處對應的時間間隔為δTσ,定軌誤差為δaσ。在δT<δTσ時,增大時間間隔可降低隨機誤差造成的定軌誤差;δT>δTσ時,情況恰好相反。同時亦可發(fā)現(xiàn),在δT>δTσ的部分,幾條曲線很快便趨于一條,說明誤差的影響不再是主導因素 (實際已轉為模型誤差為主導,正如3.1中表現(xiàn)出來的變化趨勢)。以上兩圖中拐點對應的具體數(shù)據(jù)見表3:
表3 拐點處時長及誤差列表Tab.3 Time length and error list at the turning point
本節(jié)主要討論利用低軌近圓軌道的定位資料進行初軌確定的問題,分析了資料時長、資料密度以及資料誤差的大小對定軌精度的影響,并根據(jù)定性的分析得到以下幾點結論:
(1)初軌計算所用資料長度的選擇,并非只要程序收斂即可,本文建議引入有效時長 (與其閾值相關)的概念,以此作為資料長度的上限;有效時長隨軌道高度遞增,在初步估計時可考慮采用線性近似的公式;
(2)資料隨機誤差引起的定軌誤差會隨資料時長ΔT的增加而迅速減小,在具有足夠多觀測資料的情況下,應選取資料使其長度大于理想時長 (與其閾值及與資料密度相關)以消除資料隨機誤差的影響;
(3)在資料時長ΔT固定時,所用資料組數(shù)越大定軌誤差越??;在資料組數(shù)固定時,適當增大資料時間間隔δT可降低定軌誤差。
對于初軌確定誤差分析的工作,本文只對低軌近圓軌道的定位資料作了初步的分析,對于其余類型資料、其他各類軌道以及對應的理論分析仍有待進一步研究。