李 勰,滿海鈞,劉舒蒔,陳光明
(1.中國科學院國家空間科學中心,北京100190;2.北京航天飛行控制中心,北京100094;3.航天飛行動力學技術重點實驗室,北京100094;4.中國科學院大學,北京100049)
熱層大氣密度的不確定性是影響低軌航天器定軌預報精度的最關鍵因素,而描述熱層大氣密度變化的經驗模式受建模數據以及對于密度變化自身物理機制認識的局限,雖歷經半個世紀的發展,依然無法突破15%的精度瓶頸。經驗密度模式依據建立時采用的數據源,主要分為兩大類:①基于衛星阻力數據,典型代表是美國史密松天文臺發布的Jacchia系列和法國宇航局發布的DTM系列;②基于地面非相干散射雷達和衛星原位探測數據,典型代表是美國的MSIS系列。上述3個系列的經驗密度模式在國內外航天工程任務中得到了廣泛應用。
進入21世紀以來,隨著各類工程、科學和軍事應用對航天器軌道精度的要求越來越高,傳統的經驗密度模式成為制約精度提升的一個關鍵因素。為此在原有模式基礎上,通過引入新的觀測數據和輻射指數,對部分方程進行修改,美國和歐洲相繼對上述不同模式進行了升級。Picone等在MSIS90基礎上融合新的觀測數據,并考慮500 km以上不規則原子氧和離子氧的影響,開發了NRLMSISE00模式,成為當前科學和工程應用領域使用最廣泛的模式;Bowman等在Jacchia71基礎上,引入新的太陽輻射指數和地磁指數,建立了Jacchia-Bowman模式,作為北美防空司令部軌道計算的標準模式;法國宇航局在歐盟第七研究計劃支持下,通過融合GRACE、GOCE衛星數據,推出了DTM2013模式,并宣稱該模式是目前誤差最小的模式。
由于密度模式精度對軌道計算,尤其是軌道預報影響顯著,因此模式精度的評估往往是工程應用部門關注的焦點。Marcos分別對Jacchia70/71/77、MSIS77/79/83/86 和 Jacchia71、DTM2000、NRLMSISE00和JB2006進行了系統分析比較,Vallado從應用層面比較了Jacchia60/70/71、MSIS86/90和NRLMSISE00等模式。由于驅動JB2008模式輸入參數的發布存在約3個月延遲,且最新的DTM2013模式不對外公開發布,因此國內工程領域實用的密度模式仍是數十年前的產品,除汪宏波等利用CAHMP衛星加速度計數據分析了不同輻射指數對Jacchia71、DTM94、NRLMSISE00和JB2006模式精度的影響外,國內對于新近發布模式的精度及使用建議關注較少。
中國在熱層大氣建模領域尚沒有自主的經驗模式,也沒有類似于馬歇爾工程熱層模型(Marshall Engineering Thermosphere Model,MET)的工程專用模式,因此,在實際工程應用前如何評估模式精度及計算效率,并據此選擇合適的經驗模式具有重要的理論和現實意義。本文利用CHAMP衛星的歷史數據和中國自主衛星探測數據,結合不同類型航天器的軌道預報效果,綜合評估Jacchia-Robert(JR)、NRLMSISE00(MSISE00)、JB2008和DTM2013 4種不同模式的精度,分析它們在不同應用場景下的優勢和特點,為后續的工程應用以及關聯的空間環境探測提供參考。
綜合考慮數據的精度和可獲得性,采用2類數據作為基準開展模式的精度評估:①使用廣泛且精度得到業內認可的CHAMP衛星數據;②中國自主熱層大氣密度探測與精密定軌試驗衛星(Atmospheric density detecting and Precise Orbit Determination,APOD)的觀測數據。上述數據覆蓋了400~500 km高度的大部分區域,也包含了不同太陽輻射水平和地磁條件,具有較好的代表性。
CHAMP衛星是德國地學研究中心研制的重力測地衛星。由于CHAMP衛星的軌道特性,5個月左右就能覆蓋全部緯度范圍內的所有地方,有利于考察模式誤差的全球分布特征。利用其搭載的高精度加速度計,Bruinsma等和Sutton等均開展了密度反演工作。本文采用Sutton反演的密度數據開展工作(密度下載地址為http://sisko.colorado.edu/suttion/data.htm l)。數據時間區間為2002~2008年,分辨率為45 s,涵蓋第23太陽活動周的大部分時段。
APOD衛星是中國首顆以大氣密度探測為主要目標的微納衛星,在國內首次獲取了500 km高度以下的長期連續觀測數據。由于晨昏軌道特性,APOD衛星不同緯度的地方時基本不變,數據主要集中在地方時早晚6點30分左右,這也是APOD衛星有別于其他大氣密度探測衛星的一個顯著特點。本文采用自2015年12月1日至2016年12月31日共計397天的觀測數據,數據采樣頻率為1 s。
將上述實測結果作為密度真值,對4種不同模式計算,結果與真值進行比較。定義變量R為模式計算值相對于觀測真值的誤差,如式(1)所示:



其中,N為數據個數。而模式相對真實密度數據偏差的彌散程度(隨機誤差)σ如式(3)所示:

按照CHAMP衛星的軌跡和實測空間環境參數,采用JR、MSISE00、JB2008和DTM2013 4種不同模式計算,結果與實測密度進行比較。由于數據量較大,為提高計算效率且不失統計特性,對數據進行日平均降采樣處理,其相對誤差的均值和標準偏差如圖1和表1所示。

圖1 不同模式計算結果與CHAMP衛星數據的比較Fig.1 Com parisons between differentmodels and CHAMP observations

表1 基于CHAMP數據的不同模式相對誤差統計結果Table 1 Statistical results of relative error for different models based on CHAMP
可以看出:JR與MSISE00模式存在比較明顯的非線性系統偏差,隨著時間推移,系統偏差明顯增加,2008年相對于2002年誤差擴大超過1個數量級。圖2給出了數據時間對應的不同波段太陽輻射流量,太陽輻射的整體水平呈逐年下降趨勢,2007和2008年處于一個比較低的水平。綜上,JR和MSIS00模式計算結果隨著太陽輻射的減弱而逐漸偏高;JR和MSIS00模式呈現出了非常明顯的周年/半年變化特征,且變化的幅度較大,成為顯著的誤差來源;而JB2008和DTM2013模式與觀測結果符合相對較好,相對誤差均值分別為0.4%和5.2%,除2007、2008以外,基本沒有明顯的系統偏差,其隨機偏差也均在10%左右,較JR和MSISE00模式降低約50%;JB2008模式由于引入了新的半年變化,其模式誤差基本沒有呈現周年變化的現象。

圖2 F10.7 cm太陽輻射流量指數Fig.2 F10.7 Solar radiation flux
中國自主原位探測載荷主要是質譜儀和真空計,相較于加速度計反演結果,原位探測數據精度略低。為了確保探測數據的準確性和可靠性,李勰等提出了一種基于動力學的在軌標校方法,標校后數據的系統誤差優于15%,因此作為評估模式精度的基準數據是合理可行的。需要說明的是由于DTM2013模式沒有正式對外公開發布,驅動模式運行的數據更新截止到2014年,本文只給出了JR、JB2008和MSISE00的比較結果。圖3為3種模式計算日均值與觀測日均值的比較,表2為數據統計結果。可以看出:JR和JB2008模式的系誤差均值分別為2.2%和-1.1%,系統偏差較小。而MSISE00模式誤差均值為12.4%,系統偏差略大。3個模式相對誤差的標準偏差基本相當,在15%~22%之間。相比日均值結果可以發現,通過平均也可以有效降低數據的隨機誤差。

圖3 3種模式相對APOD觀測數據日均值誤差Fig.3 Daily average error of 3 different models vs APOD observations

表2 基于APOD數據的模式相對誤差統計結果Tab le 2 Statistical results of relative error for different models based on APOD observations
由于APOD數據覆蓋時間有限,選取與APOD數據時間(2016年)范圍內太陽輻射流量(70~125 sfu)水平相當的2005年CHAMP衛星數據的比較結果,可以看出不同模式基于2類數據的評估結論基本一致:JB2008的表現均為最優,JR和MSISE00模式精度相當,較JB2008而言,有比較明顯的半年變化的系統偏差。此外,JB2008還呈現除了與其他2個模式不同的誤差特性,即JB2008計算結果較觀測值整體偏低,而其他2個模式結果偏高。
針對APOD衛星晨昏軌道特性,圖4給出了3種模式相對誤差的緯度、地方時分布,其中緯度按照5°,地方時按照1.5 h間隔進行數據平均。從圖中可以看出:3種模式的主要誤差均集中在南北半球的中高緯地區,且北半球多表現為正向誤差,南半球表現為負向誤差。此外,晨側和昏側的誤差特性也截然相反,晨側模式結果低于觀測結果,而昏側模式結果普遍高于觀測結果,MSISE00模式在昏側的誤差尤為顯著。從統計數據來看,JB2008模式的誤差均值和標準偏差均最小,分別為-0.4%和18.7%(圖4)。
將不同模式應用于軌道計算,通過軌道確定和軌道預報的效果進行分析比較是精度評估的另一種重要和有效手段,可以為工程應用提供技術依據,具有較強的現實意義。
大氣阻尼加速度計算公式如式(4)所示:



圖4 3種模式相對誤差的緯度地方時分布Fig.4 Relative error of latitude vs local time for 3 differentmodels
軌道預報精度主要取決于動力學模型誤差。目前常用的地球重力場、大行星等保守力攝動模型本身精度較高,且不會產生長期變化;而作為非保守力的大氣阻力攝動,其誤差傳播隨時間平方發散,是預報誤差的主要來源。因此,可以通過比較軌道預報的誤差來評估所使用的密度模式的精度。
以航天器對外發布的精密軌道為基準,采用不同模式進行軌道計算,將計算軌道與基準軌道進行比較。對于數據弧段內,通過比較解算阻尼系數評估不同模式的定軌精度差異;對于預報弧段,統計預報軌道與基準軌道誤差的均值和最大值來評估各個模式的預報精度。評估采用的軌道數據分別為CHAMP衛星2003年1月~6月的快速科學軌道(Rapid Scientific Orbit,RSO)數據和天宮一號2012年1月~4月的工程軌道數據,其標稱精度分別為10 cm和10 m。定軌預報過程中,為了降低其他誤差因素的影響,空間環境參數均采用實測值,面質比為固定值,阻尼系數采用定軌解算結果。
天宮一號目標飛行器是中國首個自主研制的載人空間試驗平臺,于2011年9月29日發射入軌。本文選取2012年1月~4月的測量數據,利用1天的弧段進行軌道計算,分別采用JR、MSISE00、JB2008和DTM2013 4種模式求解大氣阻尼系數,結果如圖5所示。由于模式存在誤差,在軌道解算過程中,阻尼系數會吸收模式的誤差,確保解算軌道與觀測符合最好。從圖中可以看出,JR、MSISE00、JB2008和DTM2013 4種不同模式解算阻尼系數的均值分別為2.576、2.308、2.074和2.169,反映不同模式之間的系統誤差水平。由于天宮一號目標飛行器的構型比較復雜,其阻尼系數的真實值較難確定,因此無法直接給出哪個模式最優,只能進行相對比較。如果模式計算結果高于實際密度,為了彌補模式誤差對軌道的影響,解算阻尼系數將低于真實值,反之亦然。因此,從不同模式解算結果來看,DTM2013模式計算密度整體偏高,而JR模式計算密度整體偏低,2個模式平均偏差可達24.2%,而相對偏差最小的2個模式是DTM2013和JB2008,平均偏差為4.58%,這體現了新近開發模式與早期發布模式的差異。
以不同模式解算阻尼系數的平均值為基準,其相對均值的相對偏差可以反映不同模式的隨機偏差,JR、MSISE00、JB2008和DTM2013 4種模式結果分別為23.7%、22.7%、21.1%和26.7%。可以看出,4種模式的隨機偏差大致相當,在21%~26%之間,隨機偏差最小的是DTM2013模式。
進一步利用解算軌道和阻尼系數進行24 h軌道預報,將預報軌道與基準軌道進行比較,用以統計最大位置和速度誤差。

圖5 不同模式定軌解算大氣阻尼系數結果Fig.5 Estimation of drag coefficient using different models
為了更直觀反映預報效果,將位置誤差按照500 m間隔進行頻度統計,總共121 d數據樣本,結果如圖6所示。統計結果表明:JR、MSISE00、JB2008和DTM2013 4種模式預報24 h位置偏差小于500 m的天數分別為22 d、38 d、36 d和22 d,MSISE00和DTM2013模式明顯優于Jacchia系列;預報誤差小于1 km的天數分別為51 d、64 d、63 d和59 d,仍然是MSISE00和DTM2013模式較優;而預報誤差大于3 km的天數分別為15 d、20 d、16 d和13 d,MSISE00模式效果略差。從天宮一號軌道預報效果來看,基于探測數據評估中表現不佳的MSISE00模式反而表現最好,說明實測密度數據反映的模式誤差不能完全體現在軌道預報誤差中。另一方面,由于模式誤差隨預報時間的平方發散,如果預報超過24 h,模式誤差的影響可能會更顯著增加。
利用CHAMP衛星2003年1月至6月(共計181 d)的快速科學軌道進行定軌,將解算軌道預報24 h與RSO軌道進行比較,仍然按照500m的間隔統計最大位置誤差的頻度分布,結果見圖7。結果表明:絕大多數時間的軌道預報精度均優于1.5 km,主要原因是CHAMP衛星軌道高度較高,受大氣阻力的影響比天宮一號顯著減小,因此受密度模式誤差的影響也相應減小。JR、MSISE00、DTM2013和JB2008模式預報24 h位置誤差優于1 km的天數分別為139 d、149 d、148 d和138 d,仍然是MSISE00和DTM2013表現優于其他2個模式,這與天宮一號軌道預報評估結果完全一致。

圖6 不同模式24 h軌道預報誤差頻率直方圖Fig.6 Frequency distribution of 24 h orbit p rediction error using differentmodels

圖7 不同模式24 h軌道預報誤差頻率直方圖(CHAMP)Fig.7 Frequency distribution of 24 h orbit prediction error using differentmodels(CHAMP)
2003年為太陽活動高年,為了考察不同模式在太陽地磁活動較強情況下的響應,選取全年地磁指數Ap大于50 nT的時段(共計14 d)進行軌道解算結果比較,其中超過100 nT的有6 d,2003年10月29至31日太陽輻射指數和地磁指數均達到了較強活動水平。
將所選14 d CHAMP衛星RSO軌道作為觀測量進行軌道解算,不同模式解算軌道位置的均方根誤差如圖8所示。從圖中可以看出:在地磁擾動最強烈的時段(圖8中橫軸第9~11個算例),所有模式解算軌道的位置均方根誤差均呈現增長趨勢,但區別不大,而11月20日的地磁擾動沒有10月底劇烈,但4個模式的表現卻出現了比較顯著的差別,JB2008和DTM2013顯著優于另外2個模式。

圖8 地磁擾動情況下不同模式解算位置均方根誤差Fig.8 RM S of estimated position error using differentmodels during geomagnetic disturbance
對圖8中位置均方根誤差再進行均值和標準偏差統計,結果見表3。由于數據樣本有限,誤差的統計特性不足,但依然反映出JB2008和DTM2013模式的誤差相對較低。

表3 不同模式定軌位置均方根誤差的均值和標準偏差Table 3 M ean and standard deviation of position RMS error for differentmodels
本文從應用支持決策的角度,對當前主流的JR、MSISE00、JB2008和DTM2013經驗密度模式進行了全面的精度評估。
1)利用國外權威發布和國內自主探測的密度數據,計算不同模式相對觀測誤差。結果表明:DTM2013和JB2008的系統偏差和隨機誤差較早期發布的JR和MSISE00模式有明顯的改善,整體表現更為優異。
2)從應用角度出發,以天宮一號目標飛行器和CHAMP衛星精密軌道數據為基準,計算4種模式軌道預報相對精密軌道的偏差。結果表明:MSISE00和DTM2013的效果略優于JR和JB2008模式,主要原因之一是24 h的預報時間對模式誤差傳遞的敏感度有限;在地磁擾動情況下,JB2008和DTM2013的軌道解算精度更高。
3)新發布模式由于同化了新近獲取的觀測數據,在模式精度改善方面取得了明顯進步。盡管在文中軌道預報的應用效果不明顯,但通過實測數據的比對,仍然證明了DTM2013和JB2008等新模式的精度優勢。
綜合后續的應用需求,給出以下建議:
1)高精度觀測數據對于改善現有模式精度具有顯著作用。在當前微納衛星平臺日益發展的條件下,利用微納平臺開展分布式探測,不斷提升大氣密度觀測能力和水平,將有助于開展模式修正工作。
2)為了改善對太陽輻射和地磁活動響應的不足,新近發布的模式采用了不同的指數驅動模式運行,包括很多衛星觀測數據。而目前中國在這方面的能力還十分有限,提高對太陽的多譜段觀測能力和地磁活動的觀測能力對于改善現有密度模式和提高預報能力具有重要意義。
3)與新發布的模式比較,目前國內工程中廣泛應用的MSISE00模式已經沒有明顯的優勢。針對后續應用對精度需求的不斷提高,應著手開展對現有密度模式的改進,并建立專門服務航天應用的工程模型。