利用GECAM衛星Crab脈沖星觀測數據的定軌分析

韓大煒，鄭世界，庹攸隸，葛明玉，宋黎明，李新喬，文向陽，熊少林

中國科學院 高能物理研究所 粒子天體物理重點實驗室，北京 100049

脈沖星被稱為天然的GPS衛星，其脈沖信號的長期時間穩定度很高［1］。空間飛行器利用其脈沖信息可實現自主導航，從而為太陽系內距離地球較遠的行星乃至星際旅行的航天器進行導航定軌、守時等提供幫助。自1967年人類發現第一顆脈沖星［2］、1974年Downs提出射電脈沖星導航概念［3］、1981年Chester和Butman提出X射線脈沖星星際導航［4］以來，X射線脈沖星導航技術一直被世界各國爭相研究。

國際上脈沖星導航的研究經歷了從理論提出到算法研究、空間試驗等階段。目前開展在軌試驗仍是脈沖星導航研究的主要手段。1999年美國在ARGOS衛星上開展了非常規恒星定位試驗（Unconventional Stellar Aspect，USA）進行脈沖星導航研究［5-7］。2017年NASA實施了國際空間站上的X射線計時導航技術探索（The Sta?tion Explorer for X-ray Timing and Navigation Technology，SEXTANT）項目，使用的載荷稱為中子星內部組成探測器（Neutron Star Interior Composition Explorer，NICER），從而開展了脈沖星自主導航試驗；NICER由56個相同的X射線望遠鏡組成，總有效面積約為1800 cm2，每個望遠鏡都由一個集中的X射線光學元件和一個像素硅漂移探測器（SDD）組成；SEXTANT的試驗結果顯示通過觀測4顆毫秒脈沖星獲得的定位精度可達5 km內（1σ，其中σ為標準差）［8］。

21世紀以來中國的X射線脈沖星導航研究蓬勃發展，涉及時空坐標體系、導航脈沖星源［9］，導航理論、算法與誤差分析及空間X射線探測器研制等各個方面，研究者們利用空間高能探測設備積極開展導航理論研究、脈沖星導航體制探索［10］及空間在軌驗證。迄今為止開展的空間試驗包括天宮二號空間實驗室上的“天極望遠鏡”——伽馬射線暴偏振探測器（Gamma-Ray Bursts Polarimeter，POLAR），主要科學目標是用于伽馬暴的偏振測量，探測能段為15～500 keV。其具有很大的視場（超過2π立體角）且在該能段具有較大的有效面積（約200 cm2）、較高的時間分辨率（80 ns），因此也能探測到Crab脈沖星的X射線輻射信號，從而成功進行了脈沖星導航定軌試驗，定軌精度為10 km量級［11］。

脈沖星導航試驗衛星01星（X-Ray Pulsar Navigation Test Satellite No.1，XPNAV-1）使用Wolter-I型望遠鏡［12-13］和SDD探測器觀測Crab給出了初步的導航精度，利用25次觀測的數據計算出脈沖周期的均方根誤差（RMS）為0.00264，對應的脈沖到達時間（Time of Arrival，TOA）誤 差 為89.2 μs，距 離 誤 差 為26.8 km［14］。

中國首顆X射線天文衛星——“慧眼”衛星（Insight， Hard X-ray Modulation Telescope）使用3臺望遠鏡、通過對蟹狀脈沖星的觀測成功確定了軌道。試驗分別使用“慧眼”衛星3臺望遠鏡（高能X射線望遠鏡（High Energy X-Ray Tele?scope，HE）、中能X射線望遠鏡（Medium Energy X-Ray Telescope，ME）和低能X射線望遠鏡（Low Energy X-Ray Telescope，LE））各5天的觀測數據，均成功確定了軌道參數，綜合所有數據可知定位精度在10 km（3σ）［15］。

本文通過分析GECAM衛星對Crab脈沖星的觀測數據說明基于軌道動力學的脈沖輪廓波形典型顯著性分析方法進行脈沖星定軌的原理和步驟，并給出計算結果和與其他衛星的對比。

1 GECAM衛星

懷柔一號引力波暴高能電磁對應體全天監測器（Gravitational Wave High-Energy Elec?tromagnetic Counterpart All-Sky Monitor，GECAM）衛星是中國科學院“空間科學”（二期）戰略性先導科技專項的天文衛星，其科學目標是對引力波伽馬暴、快速射電暴高能輻射、特殊伽馬暴和磁星爆發等全天隨機爆發的高能暫現源的監測研究。GECAM衛星采用600 km高度低軌和29°低傾角的天文觀測軌道［16］，包含2顆相同的運行于同一軌道面的微小衛星，分別為GECAM-A和GECAM-B。伽馬射線探測器（Gamma-Ray Detector，GRD）是GECAM的主要探測器載荷，其主要任務和功能是實現對空間伽馬射線暴發源的寬能段能譜和時變觀測，并通過多個探頭實現對全天的監測及對伽馬暴的定位。

根據GECAM衛星的構型，GRD（圖1中的圓形探測器）近于均勻地分布在橢球形的穹頂結構上，每顆衛星共有25個GRD模塊，均安裝在衛星的穹頂艙，每個GRD模塊均能探測從前部入射的伽馬射線光子，視場約為2π。25個GRD分別指向不同方向，在立體角上基本均勻分布［17］。平均情況下，約12個探頭可觀測到Crab脈沖星。某個方向上所有探頭理論上均可接收X射線光子，有效面積合計約250 cm2［18］。GRD的主要技術指標如表1［17］所示。

圖1 GECAM有效載荷的構型設計Fig.1 Configuration design of GECAM payload

表1 GECAM-GRD的主要技術指標［17］Table 1 Main technical parameters of GECAM-GRD［17］

總之GECAM衛星的探測能區寬、時間分辨率和時間精度都很高［19］，觀測天區覆蓋廣。因此GECAM能長時間接收到Crab脈沖星發射出的X射線光子，可用于開展脈沖星導航試驗。

2 SEPO算法

基于軌道動力學的脈沖輪廓波形典型顯著性分析方法（Significance Enhancement of Pulseprofile with Orbit-dynamics，SEPO）是結合脈沖輪廓和軌道動力學顯著性分析的方法。經典導航領域常用的卡爾曼濾波方法已在中國廣泛應用于X射線脈沖星的仿真分析和在軌導航試驗中［20-21］。筆者采用與之不同的SEPO算法，此方法已在天極望遠鏡（POLAR）和“慧眼”衛星（Insight-HXMT）中進行了試驗［11，15］。

SEPO算法的基本原理是軌道動力學模型能在短時間內得到高精度的軌道預報，但具有長期性漂移，特別是低軌衛星由于大氣模型的變化導致軌道的長期預報偏差很大。另外當初始軌道存在偏差時預報軌道與真實軌道自然會出現偏離，且該偏離隨著時間的延長越來越大。根據真實軌道能計算出脈沖輪廓，而軌道偏離會導致脈沖輪廓發生變形，導致輪廓信號的顯著性降低。與“慧眼”衛星類似，GECAM衛星工作在約600 km高度的低地球圓軌道，其軌道預報的變化效應和對輪廓的影響比較類似，同樣可采用SEPO算法進行定軌。

SEPO算法的整個計算過程包括利用軌道動力學外推、計算脈沖波形的顯著性和基于網格搜索的軌道擬合3大步驟，詳見文獻［11，15］。

SEPO算法的優勢是僅使用一顆脈沖星即可實現衛星定軌。目前的試驗結果表明這種方法更適用于軌道根數變化劇烈且不存在完全“平行”的軌道，這樣這些變化可很好地顯示在脈沖輪廓中。另外由于脈沖星在不同能段的輪廓形狀會有所差異，所以不同類型的空間探測器由于工作原理和能段不同無法套用“脈沖標準輪廓”，SEPO算法使用自身探測的輪廓突破了該限制。

相較于最初在POLAR試驗中的應用，使用的SEPO算法做了進一步的完善，體現在利用軌道動力學加上偏差計算軌道時由POLAR時的向后外推更改為往前推算。因此計算結果可得到“現在”的軌道根數，而非POLAR計算時得到的“過去”的某個起始點的軌道根數，因而本文的計算具有了更大的工程應用價值。

3 數據篩選和預處理

通常情況下使用X射線脈沖星的觀測數據進行導航定軌的處理步驟分為4步：事例挑選、太陽系質心修正、自轉頻率搜索和脈沖輪廓疊加、脈沖TOA的計算。

對GECAM衛星的觀測數據處理利用了前3個步驟，使用SEPO算法并未用到TOA的計算。具體步驟如下。

3.1 事例挑選

1） 剔除地球遮擋數據

根據Jodrell Bank射電天文臺提供的Crab坐標值可知其在J2000坐標系下為赤經RA=05d34m31.972s、赤緯DEC=22h00m52.07s，可將該RA、DEC轉換到直角坐標系（xcrab，ycrab，zcrab）中。衛星的姿態數據中提供了衛星在J2000的直角坐標系中每秒的位置（x，y，z），這樣可先計算衛星的位置向量與Crab方向向量之間的夾角，從而可計算GECAM衛星與Crab的連線與地心之間的夾角。若Crab完全被地球遮擋（夾角小于67°即被地球遮擋），則剔除該部分的光子事例。

2） 探測器篩選

根據衛星的姿態及各GRD探頭的位置可計算Crab脈沖信號分別入射到25個GRD探測器的夾角。在衛星的姿態文件中給出了逐秒的衛星姿態四元數。根據姿態四元數及衛星載荷信息中各個GRD探頭與衛星主軸的位置關系可逐秒計算得到Crab脈沖星的位置矢量與各個GRD探頭的夾角。對于各個GRD探測器探測到的光子事例，選擇Crab與各探頭夾角小于70°時間段內的光子。

3） 能段篩選

選擇GECAM的能道不大于100的高增益事例，對應的是約10～100 keV的低能段事例。這樣能使來自于非目標源Crab的光子事例占比較小，增加Crab脈沖信號的顯著性。

選取GECAM-B星的數據，時間段為2021年7月23日—8月31日，總光子數為20.98億個。對所有光子進行折疊并扣除本底后統計脈沖輪廓的光子數，脈沖輪廓的光子數約為800萬個，約占總光子數的0.36%。

3.2 太陽系質心修正

太陽系質心坐標系是一個較好的慣性參考系，為得到準確的脈沖星脈沖到達時間，需將探測到的脈沖光子到達時間轉換到太陽系質心坐標系中。這主要包括真空延遲修正、相對論修正及不同時間系統之間的轉換等。使用太陽系星歷DE405進行太陽系質心修正。

3.3 脈沖輪廓折疊

由于直接利用了Fermi衛星的Crab脈沖星星歷，因此未作自轉頻率搜索。根據星歷計算出所有觀測的X射線脈沖光子的相位，統計所有光子的相位信息得到Crab脈沖星的X射線脈沖輪廓。

4 計算結果與分析

利用SEPO算法計算時首先給出估計的初始軌道參數和參數空間，利用軌道動力學模型進行30天的軌道預報，計算得到基于該預報軌道的觀測脈沖輪廓；然后通過對此觀測脈沖輪廓的顯著性分析搜索最優軌道參數。具體流程如圖2所示。

圖2 脈沖星定軌算法的處理流程Fig.2 Processing flow of pulsar orbit determination

進行軌道外推時使用的積分器為Runge-Kutta-Fehlberg7（8），積分步長為60 s，在此情況下可滿足軌道外推的需要，考慮的主要動力學模型參數如表2所示。

表2 軌道動力學模型中的參數與設置Table 2 Parameters and settings in orbital dynamic models

根據軌道動力學模型通過起始點向前推算了30天軌道位置與GECAM衛星真實測量的位置存在偏差，計算結果如圖3所示。可見偏差的均值為?5.10、?3.17、?0.77 km，1σ標準差分別為7.56、9.28、5.22 km。在計算的30天內偏差總的趨勢是逐漸放大，呈現比較復雜的變化趨勢，平均在km級。計算Crab脈沖星的輪廓分為500個相位，每個相位的周期約為66 μs，折算到距離為20 km。可見軌道計算的偏差對脈沖星輪廓的影響在一個脈沖輪廓相位的范圍內，未影響計算結果。

圖3 GECAM軌道起點向前外推30天的位置與真實測量值的偏差Fig.3 Deviations between positions extrapolated for?ward for 30 days from orbit starting point of GECAM and true measured values

GECAM衛星的軌道狀態用J2000坐標系下的6個軌道根數表示，分別是半長軸（a）、偏心率（e）、軌道傾角（i）、升交點赤經（ω）、近地點幅角（Ω）和平近地點角（w）。

首先選取GECAM衛星在觀測結束時刻的位置作為初始值，兼顧軌道參數覆蓋范圍和計算效率，選用軌道根數初始值的偏移量如表3所示。

表3 軌道參數初始值的偏移量Table 3 Bias of initial value of orbital elements

在此參數空間內的每個格點處分別計算得觀測脈沖輪廓。圖4為觀測Crab的原始輪廓和半長軸分別改變?520 m和+480 m的模擬輪廓，可見其發生了畸變。

圖4 GECAM觀測Crab脈沖星輪廓和半長軸變化的仿真輪廓Fig.4 Profile of Crab pulsar observed with GECAM and profile simulated with semimajor axis

直接對觀測脈沖輪廓進行處理。定義脈沖輪廓的顯著性χ2為

式中：p（i）為經過歸一化的脈沖輪廓；pˉ為p（i）的平均值；k為脈沖輪廓的相位數目。

計算每個軌道參數格點對應脈沖輪廓的顯著性，結果如圖5所示，圖中計算結果為不同參數的帶誤差卡方計算結果，擬合結果為最優值擬合結果（高斯擬合），0值點為實際軌道根數值。可見脈沖輪廓的顯著性χ2隨軌道參數的偏離有顯著變化，在真實值（0值）附近χ2具有極大值。通過高斯函數曲線擬合得到的最大值可分別給出6個軌道參數的最優解。以半長軸為例，當起始時刻的半長軸偏離±50 m時計算得到的脈沖輪廓會發生畸變，其顯著性χ2也從3.10×104下降到2.85×104（偏離±50 m）。

圖5 χ2 隨軌道參數偏移量的變化結果Fig.5 χ2 variation results with orbit elements bias

采用Sheikh等計算脈沖TOA誤差的方法［6，22］計算軌道參數誤差，每個參數的誤差σ為

式中：W為參數最優值的高速擬合過程中的半高全寬（FWHM）；S/N為信噪比。采用相同的計算方法［6］：

式中：S/Nfiltered為加了濾波限幅的信噪比。為限制S/N可能出現極大的異常值，采用限幅公式對式（3）進行限制，最大為1000。從而得到每個軌道根數的最優值和3σ誤差（99.7%置信度），如表4［11，15］所示。

表4 軌道參數最優值和誤差［11，15］Table 4 Best estimated values of elements and errors［11，15］

通過比較表4［11，15］可見3次試驗中半長軸的誤差比較接近；分析原因是半長軸的輕微改變就會導致衛星運動的軌道周期發生變化，使衛星位置的變化隨時間不斷累積，從而使脈沖輪廓發生較大的改變，因此SEPO方法對半長軸有較好的約束。其余5個軌道參數的誤差中“慧眼”衛星是最優的，GECAM衛星與POLAR望遠鏡相比偏心率、傾角和升交點赤經的誤差較小，而近地點幅角和平近地點角的誤差較大。

5 總結與展望

利用“懷柔一號”GECAM衛星對Crab脈沖星的在軌觀測數據進行了脈沖星導航定軌的原理驗證。通過建立基于軌道動力學模型和輪廓分析的導航算法初步實現了利用脈沖星探測進行空間飛行器的定軌。針對GECAM-B星40天的在軌觀測數據得到了初始歷元時刻6個軌道參數的最優值與誤差，其最優值均與真實值在3倍誤差范圍內吻合。3σ的定軌誤差與“慧眼”衛星（Insight-HXMT）脈沖星導航的試驗結果［15］相比除半長軸和傾角相近外，其他4個軌道根數的誤差是“慧眼”衛星的2～4倍。證明了GECAM衛星同樣可用于確定軌道根數，但由于本底誤差遠高于“慧眼”衛星，故導航誤差也大于“慧眼”衛星的10 km。

另外值得注意的是由于GECAM作為具有極寬視場的全天監測的天文衛星，盡管選取了70°的視場本底仍然遠高于“慧眼”衛星，因此選取40天的數據后軌道根數變化顯著性的擬合才有比較明顯的概率分布，且最優值偏離在3倍標準差之內。可見GECAM用于定軌所需的時間遠超過“慧眼”衛星的5天，而與POLAR試驗中的31天接近。這個結果也說明本底會對計算造成較大的影響。

需要說明的是，截至目前的3次試驗均未考慮航天器上的時鐘時間偏移及脈沖星內在的紅噪聲和躍變屬性等的影響。

從望遠鏡本身的研制而言，GECAM衛星的特點在于載荷對衛星平臺的資源需求較低。相較于“慧眼”衛星觀測脈沖星時指向特定天區的方式，GECAM采用大視場的設計，指向精度和姿軌控制的要求低。定軌數據處理過程中通過衛星姿態和觀測角度的方法排除一部分非Crab的本底事例。另外單顆GECAM衛星的質量只有180 kg，單星發射包絡尺寸不大于580 mm（X）×1050 mm（Y）×1364 mm（Z）［16］，與“慧眼”衛星的約2500 kg相比GECAM衛星的質量大幅度減少。

而安裝在中國天宮二號上的POLAR載荷盡管其本身質量較小，包含電控箱只有30 kg，但依賴于空間實驗室提供能源、姿軌控、測控等支持。GECAM衛星實現了載荷衛星一體化研制。

總之，本文在“慧眼”等傳統的天文衛星和POLAR等空間實驗室的研究基礎上，采用衛星小型化和載荷輕型化之后的GECAM衛星成功實現了定軌，為中國未來脈沖星定軌和導航空間試驗的開展提供了新思路。