1、航空学报Acta Aeronautica et Astronautica SinicaFeb.15 2023 Vol.44 No.3ISSN 1000-6893 CN 11-1929/V526597-1X射线脉冲星导航的最优观测周期确定苏剑宇1，方海燕1，*，高敬敬1，赵良21西安电子科技大学 空间科学与技术学院，西安 7101262北京临近空间飞行器系统工程研究所 空间物理重点实验室，北京 100101摘 要：最优观测周期的确定对脉冲星导航计算具有重要意义。首先推导了航天器处观测脉冲相位估计方差的理论下界，然后给出了观测脉冲相位估计误差与脉冲星观测周期之间的关系式，并以观测脉冲相位估计的均方
2、误差最小为准则，给出了最优观测周期的近似计算公式，最后利用脉冲星 PSR B0531+21的实测数据验证了该关系式的正确性。仿真结果表明所给最优观测周期计算公式对脉冲星 PSR B0531+21的预测误差为 44 s，证明了所给公式的正确性，为导航中脉冲星观测周期的确定提供了理论基础。关键词：脉冲星导航；观测脉冲相位；观测周期；理论下界；均方误差中图分类号：V249.32+3 文献标识码：A 文章编号：1000-6893（2023）03-526597-12X 射 线 脉 冲 星 导 航（X-ray Pulsar-based Navigation，XPNAV）是一种新型的自主导航技术，可为航天器
3、提供位置、时间等导航信息1-4，实现航天器高精度自主导航。目前，国内外相继开展了 X 射线脉冲星导航试验，如中国空间实验室天宫二号（Tiangong-2，TG-2）的 暴偏振探测科学实验5，中国首颗射线脉冲星导航试验卫星（XPNAV-1）在轨开展的射线脉冲星的探测与脉冲星导航体制的验证6-7，以及国内首颗空间X 射 线 天 文 硬 X 射 线 调 制 望 远 镜 卫 星（Hard X-ray Modulation Telescope，HXMT）Insight-HXMT 的脉冲星定轨精度验证实验8。如美国NICER（Neutron star Interior Composition Explor
4、er）项 目 的 SEXTANT（Station Explorer X-ray Timing And Navigation Technology）搭载国际空间站（International Space Station，ISS）开展的定轨精度验证工作9。可见，X 射线脉冲星观测与导航验证工作正进入蓬勃发展的时期。XPNAV 的基本原理可描述为10：在脉冲星观测周期内，航天器上安装的 X 射线探测器会记录到一串光子到达时间（Time of Arrivals，TOA），利用观测的光子 TOA，通过一定的算法提取出在这段观测时间内的某一时刻处航天器所接收的脉冲星信号的相位以及多普勒频率，由于航天器在任
5、意时刻观测的脉冲信号相位和频率都可以用该时刻的位置、速度以及太阳系质心（Solar System Barycenter，SSB）处的脉冲星信号模型准确表示，因此以估计得到的航天器在某一时刻的观测信号相位和多普勒频率作为导航观测量，并将其表示为有关航天器在该时刻位置和速度的关系式作为导航观测方程，获得任意时刻的观测相位即可得到相应时刻的位置信息。高精度http:/ 引用格式：苏剑宇，方海燕，高敬敬，等.X 射线脉冲星导航的最优观测周期确定 J.航空学报，2023，44（3）：526597.SU J Y，FANG H Y，GAO J J，et al.Determination of optimal
6、 observation period for X-ray pulsar-based navigation J.Acta Aeronautica et Astronautica Sinica，2023，44（3）：526597（in Chinese）.doi：10.7527/S1000-6893.2022.26597收稿日期：2021-11-02；退修日期：2021-11-30；录用日期：2022-03-14；网络出版时间：2022-03-2311：40网络出版地址：https：/ 内 外 学 者 也 相 继 提 出 许 多 脉 冲 相 位 估 计算法11-18。在上述导航算法中，为获得高精度
7、的位置与速度估计，需要实现观测相位和频率的高精度估计。文献 19 提出在星际转移轨道段可将航天器的运动轨迹简化为匀加速模型，基于此建模了航天器轨道状态误差与其所引起的脉冲TOA 偏差之间的关系，并对导航观测方程进行了修正以消除 TOA 偏差的影响。但线性化过程或产生高阶截断误差，且该截断误差随脉冲星观测周期的增大而增大。所以通常为了将轨道位置误差进行线性化，选择压缩脉冲星观测周期或将观测周期分段。如文献 20 提出将观测时段划分为短时段，在每个较短时段内利用针对恒定多普勒频移模型的相位估计算法估计相位，再借助相位跟踪的方法提高单段相位估计精度并获得动态的多普勒频移估计。文献21 提出了由轮廓形
8、变来求解速度误差引起的多普勒频移的新思路，定义了轮廓熵以衡量累积轮廓的不同形变程度，推导了存在速度误差时，频率偏差与轮廓熵之间的关系，利用最优化方法通过使熵最小（轮廓最为清晰）从而求解出多普勒频移。该方法本质上仍需要假定在观测时段内的观测频率恒定，即航天器的速度恒定，所以观测周期不能过大。Anderson 和 Pines22研究了基于锁相环的脉冲相位跟踪方法，并给出了针对小流量脉冲星的相位跟踪方法。文献23 将较短观测时段内的预报位置误差建模为随时间近似线性变化，且将这种线性变化的位置误差对光子到达时间修正的影响建模为恒定的相位偏差加上恒定的多普勒频移，在此基础上，提出了一种基于相位和多普勒频
9、移联合观测的 X 射线脉冲星导航方法，并进行了实验验证与分析。以上算法在进行观测相位与多普勒频率估计时，未对脉冲星观测周期的选取进行深入研究。而最优观测周期的确定对脉冲星导航计算及观测计划制定等具有重要的意义，本文对该问题进行研究与讨论，以期为观测方案的制定提供参考。1观测相位估计误差最优脉冲星观测周期是为得到高精度的观测相位估计，所以将观测相位估计误差最小时对应的脉冲星观测周期称为最优脉冲星观测周期。航天器在(t0，tf)的脉冲星观测周期内，任意时刻的观测相位10为sc(t)=?sc(t)+?()t+n r?(t)c n r(t)c（1）式 中：?sc(t)为 观 测 相 位 的 预 测 值
10、；定 义fs=?()t0+n r?(t0)c为脉冲星源频率；n为所观测脉冲星的单位方向矢量；c为光速；r为位置误差。sc(t)的定义中包含了航天器位置信息，所以对sc(t)的估计误差最小时对应的脉冲星观测周期即为最优脉冲星观测周期。根据式（1），对sc(t)的估计误差取决于对 r的估计误差。任意时刻的位置误差r为r=r0+m=1(r0)(m)m！Tm（2）式中：T为脉冲星观测周期长度，T=tf-t0；r0为航天器初始时刻t0的位置误差。则r引起的相位误差为=0+fsnm=1(r0)(m)m！Tmc（3）式 中：0为 对 应t0时 刻 的 初 始 预 测 相 位。将式（3）展开到 2阶，表示为=
11、0+f0T+f?02T2+f?()6T3 t0 t0+T（4）式中：0=fsnr0c为由r0引起的初始相位偏差；初 始 频 率 偏 差 为f?0=fsna0c；f?()6T3为高 阶 截 断 误 差，可 近 似 表 示 为f?()fSna?(t0)c；v0、a0分别为初始时刻t0时的速度误差与加速度误差。记=(0，f0，f?0)，则对的估计误差由两部分组成，分别为的估计误差与截断误差。截断误差中的a?可利用加速度计航空学报526597-3测 量 相 邻 两 时 刻 的 加 速 度，再 进 行 差 分 运 算求解。根据上述分析，对的估计将同时存在方差与偏差，采用均方误差准则，并用mse()表示的
12、估计均方误差，则均方误差mse()为mse()=var()+b2()（5）式中：var()为的估计方差；b2为的估计偏差，b2()=f?2()36T6。方差是由的估计方差 引 起 的，所 以 求 出 估 计 方 差，即 可 得 到 使mse()最小的T的表达式。根据文献 12，利用探测器探测到的光子到达时间序列，使用最大似然估计可实现对脉冲相位延迟的估计，利用最大似然估计法对的估计可表示为?=arg max()i=1Mkln(b+sh(i)（6）式中：s、b分别为脉冲星源的流量与背景流量；h()为相位 0，1)的脉冲星的标准轮廓，满足01h()d=1，h(n+)=h()（n为整数）12。Ema
13、dzadeh 证明了利用式（6）的最大似然估计对脉 冲 相 位 延 迟 的 估 计 方 差 即 为 克 拉 美 罗 界（Cramer-Rao Lower Bound，CRLB）。所 以 下面推导对估计的 CRLB。2最优脉冲星观测周期确定2.1的 CRLB推导任意时刻的相位误差为(t)=0+f0(t-t0)+f?02(t-t0)2（7）将的 CRLB表示13为CRLB()=I-1()（8）式中：I()为 Fisher 信息矩阵，根据 Emadzadeh和 Speyer13的结论，I()的第 i行第 j列元素Iij的表达式为Iij=t0tf()t；j()t；i1()t；dt（9）其中：(t)=b
14、+sh()（10）由Iij的表达式可得Iij=Iji，所以 I 为对称矩阵，且I13=12t0tf s(t-t0)h?()(t)2b+sh()(t)dt=12I22（11）为确定I()矩阵，需要求I11、I12、I13、I23、I33 5个元素。Emadzadeh 给出了I11、I12、I13的表达式，附录 A 给出了I23、I33的推导过程，结果如式（A19）所示。得到I11、I12、I13、I23、I33的表达式，即可求得的CRLB()。附 录 B 给 出 了CRLB()的 推 导 过程。直接求解化解CRLB()难以化解，这里首先化解CRLB(0)I11=(2I223-4I13I33)I1
15、1a（12）式中：a为I()的行列式，将I11、I12、I13、I23、I33代入式（12），化解得到(2I223-4I13I33)I11a=(2I223-4I13I33)I112I33I212-4I12I13I23+4I313-4I11I33I13+2I11I223=9N4-450N3+1 077N2-852N+204N4-106N3+105N2+100N-112 9 （13）所以CRLB(d0)=9I11=92sI0T。同理得到CRLB(f0)I22=(2I213-2I11I33)I22a=64N4-520N3+1 032N2-782N+194N4-106N3+105N2+100N-112
16、64（14）所以CRLB(f0)=64I22=64N(N-1)(2N-1)6P3112sI01922sI0T3（15）CRLB(f?0)I33=(2I212-4I11I13)I33a=36N3-162N2+222N-102N3-107N2+212N-112 36（16）所以CRLB(f?0)=36I337202sI0T5。综上可得航空学报526597-4CRLB()=92sI0T 1922sI0T3 7202sI0T5 T（17）式（17）即为对参数估计的 CRLB。2.2最优观测周期确定得到=(0，f0，f?0)的 CRLB 后，可确定的估计误差。根据式（17）、式（4）中的估计方差var()可表示为var()=var(0)+var(f0T)+var()f?0T22=92sI0T+1922sI0T3T2+1802sI0T5T4=3812sI0T（18）式中：2sI0为脉冲星参数对估计方差的影响；T为脉冲星观测周期对估计方差的影响，如果只考虑估计方差，则观测周期越大，估计方差越小，此时不存在最优观测周期。根据式（5），对的估计中会产生偏差，所以在计算对的估计均方误差时需要考虑估计偏差的
