脉冲星导航试验卫星观测数据处理与分析

易韦韦，偶晓娟,，许静文，李晶，李冰

（1. 宇航动力学国家重点实验室，西安 710043；2. 63769部队，西安 710043）

0 引 言

脉冲星是高速自转的中子星，具有极其稳定的周期性，被誉为自然界最精准的天文时钟。因此，脉冲星能够成为人类在宇宙中航行的“灯塔”，为近地轨道、深空[1]和星际空间飞行[2]的航天器提供自主导航信息服务。导航卫星信号容易受到人为干扰，星座还可能被摧毁，作用具有局限性。脉冲星信号有更强的抗干扰能力，几乎无法被摧毁。因此，各国都十分重视对此技术的探索和研究。

美国针对脉冲星定位导航进行了一系列研究计划和试验。1995年，美国航空航天局（National Aeronautics and Space Administration，NASA）发射了罗希X射线计时探测器，实现了对Crab脉冲星的观测，验证了脉冲星导航和轮廓折叠算法。1999年，斯坦福大学线性加速器中心和美国海军研究实验室共同设计研制的非常规天体试验（Unconventional Stellar Aspect，USA）系统搭载在先进研究和对地观测卫星（Advanced Research and Global Observation Satellite，ARGOS）上升空，试验结果表明了X射线脉冲星导航是可行的。

2004年，美国国防高级研究计划局的战术技术办公室启动了XNAV（X-ray Source-based Navigation for Autonomous Position Determination）计划，其目标是验证利用X射线源进行航天器导航的可行性，研制开展空间试验所需要的有效载荷。XNAV计划仅实施了概念可行性阶段，即脉冲星特征描述、导航算法开发、探测器原型设计、应用效果评估。

2011年，NASA的戈达德空间飞行中心（Goddard Space Flight Center）联合美国大学空间研究联合会启动了“空间站X射线计时与导航技术试验”（Station Explorer for X-ray Timing and Navigation Technology，SEXTANT）项目，目标是通过观测2～3颗脉冲星实现低轨航天器1 km的实时定规精度，项目正在实施过程中。

1999年，欧洲空间局（European Space Agency，ESA）发射了多镜面X射线观测卫星。2004年，启动了“ESA深空探测器脉冲星导航”计划，开展了脉冲星导航理论方法和仿真验证研究工作。2009年，俄罗斯科学院宇宙空间研究所宣称已在国际天体物理学天文台和国际γ射线天体物理实验卫星上进行了脉冲星观测试验，并取得了良好效果。2015年，印度发射了多波段天文观测卫星Astrosat实现了对Crab的观测。

2016年11月10日，我国在酒泉卫星发射中心发射了脉冲星导航试验卫星（X-ray Pulsar-based Navigation-1，XPNAV-1），目的是开展X射线脉冲星在轨观测，验证脉冲星导航技术。卫星质量约270 kg，轨道高度500 km，倾角97.4°。XPNAV-1搭载了两个X射线探测器：一个为掠入射Wolter-I聚焦型X射线探测器（简称Wolter-I探测器），另一个为准直型微通道板（Microchannel Plates，MCP）探测器（简称MCP探测器）[3-4]。

2017年5月，XPNAV-1卫星的首批在轨试验数据发布，这批数据是Wolter-I探测器的探测数据，数据下载地址为http://www.beidou.gov.cn/xpnavdata.rar。

文献[3]～[4]中钱学森空间技术实验室针对这批数据介绍了脉冲到达时间转换涉及的相关技术，给出了处理数据的软件、处理流程及方法，指出XPNAV-1卫星能够实现“看得见”，对推动脉冲星导航技术向前发展有着重要意义。

本文在此基础上，引入核回归法对轮廓进行了降噪处理，并分析这一方法对折叠周期所分格数对计时残差影响的改变，得出了一些有意义的结论。

1 光子到达时间转换

XPNAV-1卫星发布的数据采用协调世界时（Coordinated Universal Time，UTC）记录光子到达时间[3]，需要将之转换为太阳系质心力学时（Barycentric Dynamical Time，TDB）或太阳系质心坐标时（Barycentric Coordinate Time，TCB）以进行修正和后续处理。

本文选择将其转化为TDB，转换次序为UTC-TAITT-TDB国际原子时（International Atomic Time，TAI）、地球时（Terrestrial Time，TT）。通过调用国际天文协会（International Astronomical Union，IAU）的标准基本天文程序库SOFA函数，能够更加准确快捷地完成这一任务。需要用到的SOFA函数包括：iauCal2jd、iauUtctai、iauTaitt、iauDtdb、iauTttdb。转换过程中要注意的是，XPNAV-1星的数据文件中的UTC秒是相对UTC2008年1月1日0时0分0秒的秒数。

时间转化完成后，还需要对其进行修正，消除各类影响，得到准确的光子到达太阳系质心时间。这些影响包括：Roemer延迟、Shapiro延迟、时钟修正、Einstrin延迟等[4-6]。依据数据情况和影响大小，对Roemer延迟和Shapiro延迟进行修正[4-6]。

Roemer延时计算公式为[6]

表1 Crab自转参数[8]Table 1 Crab rotation parameters

Shapiro延时是由第三体引力导致的，第i颗行星带来的影响计算公式为[4]

JPL会及时更新各天体的星历数据，利用拉格朗日插值和JPL星历，可得到需要时刻的行星星历。结合式（2），可计算各行星带来的Shapiro延迟。

2 脉冲积分轮廓生成

光子到达时间转换完成后，使用Jodrell Bank天文台利用Crab射电观测数据得到的脉冲频率（）及其一、二阶导数（、）作为脉冲星自转模型初始值，对所有观测时刻的脉冲频率进行外推，继而利用历元折叠可获得各时段观测数据的折叠轮廓。

1）基于历元折叠的脉冲轮廓生成

历元折叠是从光子到达时刻（Time of Arrival，TOA）中恢复X射线脉冲星密度函数的方法[1-2，4，7]。其过程是：将一个脉冲周期P分为nbin个等间隔的格，简称bin，然后依据到达时间将每个光子折叠到对应的bin中，设总的光子数为，折叠后每个bin内的光子数为那么最后，将计算得到的结果进行标准化处理，即可得出脉冲轮廓。设总的观测时长为，那么第个bin的标准化光子数计算公式为

折叠采用的Crab自转参数见表1。折叠的部分结果见图1所示，bin数为128。图1中能够清晰地看到Crab的脉冲轮廓，表明了光子到达时间的转换和折叠方法是正确的。也证明了脉冲星试验01星对Crab的“可见”。

图1 脉冲轮廓结果Fig. 1 Results of pulse profiles

2）基于核回归的轮廓降噪

核回归是一种非参数回归方法，对观测量的分布具有较少限制，比参数回归有更大的适应性。采用核回归对得到的原始轮廓进行降噪处理[4]，其原理如下。

设非线性模型为

高斯核函数具有相当高的灵活性，当特征数远小于样本数时，一般使用高斯核函数。

因此核函数选用高斯核函数，即

图1是原始折叠及核回归处理后的部分脉冲轮廓。bin数为128。图1（a）图为折叠轮廓，图1（b）图为核回归处理后的轮廓。从上到下分别为公开数据的第5、15、25组的处理结果，有效观测时长分别为485.8 s、2 554.9 s和5 997.0 s。显然，核回归法能够平滑轮廓，且观测时间短，对信噪比提高更为明显。

3 计时残差估计

在参数完全精确的前提下，可准确获得任意时刻的脉冲星信号的相位。然而，局限于当前的观测手段和理论水平，难以准确标定出脉冲星自转频率的各阶导数。由此，导致了计时噪声，通常用计时残差表示。生成脉冲轮廓的目的实际上是为后续计时残差处理，继而更新脉冲星星历和标准轮廓作铺垫。首先需要进行相位延迟估计。

1）相位延迟估计

目标函数为

式（10）指出相位延迟估计需要标准轮廓这一先验知识，根据表1参数和所有光子到达时间，推算出标准轮廓，如图2所示（bin数为128）。

图2 标准轮廓结果Fig. 2 Result of normal pulse profile

2）计时残差估计

通过相位延迟，可计算脉冲到达时间残差[6]为

表2是按照表1自转参数得到的计时残差结果。结果显示，bin数的改变导致计时残差出现较大变化，核回归算法有效降低了噪声，改善了bin数对计时残差的影响，XPNAV-1卫星发布的首批近一个月的数据处理得到的计时残差均方根约为40 μs。

文献[1]处理了XPNV-1星3.5个月数据，得到计时残差均方根为55.94 μs。其值较大的主要原因是Crab周期变化快，观测时间越长，计时残差将越大。因公开数据信息已被处理简化，且实际处理结果较好（Jodrell Bank天文台公布自转参数准确），本文不再拟合计时残差和更新Crab星历参数。

表2 计时残差估计结果Table 2 Time residual results

4 结 论

本文针对我国脉冲星导航试验卫星公开发布的第一批观测数据进行初步分析，主要成果和结论总结如下：

1）文中方法实现了脉冲星试验01星数据处理，得到了各组脉冲轮廓，得出了多组观测的计时残差，能够为后续星历更新、脉冲星导航提供支持。Crab脉冲轮廓的获得，表明了脉冲星试验01星具备“看得见”的能力。

2）核回归算法能够提高轮廓信噪比，有效改善bin数对计时残差的影响。35组数据拟合前计时残差均方根约40 μs，折算为距离约12 km，因此要想达到千米级甚至更高的深空导航精度，需进一步提高观测精度。