王浙宇,韩孟纳,石永强,张哲浩,童明雷

(1.中国科学院国家授时中心,西安 710600;2.中国科学院大学天文与空间科学学院,北京100049;3.北京控制工程研究所,北京 100190)

0 引言

脉冲星是高速自转的致密星,具有超高温、超高压、超高密度、超强电磁场和超强引力场等极端的物理环境[1]。脉冲星在射电、红外、可见光、紫外、X 射线和γ射线等波段产生信号辐射,在X 射线频段辐射信号的脉冲星称为X 射线脉冲星。X射线是高能光子,集中了脉冲星的绝大部分辐射能量,易于小型化探测器进行观测,但难以穿过稠密的地球大气层,因此只能在空间进行观测。

X射线脉冲星能为航天器提供位置、速度、时间和姿态等高精度导航信息,在深空探测前沿研究领域具有重要的战略意义[2-3]。国外自20世纪70年代起开始将目光投向脉冲星导航,开展相应的理论和技术研究,并制定了一系列空间试验计划[4-6]。2011年,NASA的Goddard空间飞行中心(Goddard Space Flight Center,GSFC)联合美国大学空间研究联合会启动了空间站X射线计时与导航试验(the station explorer for X-ray timing and navigation technology,SEXTANT)项目[7]。2017年6月,中子星内部成分探测器(neutron star interior composition explorer,NICER)在国际空间站(international space station,ISS)部署,并开展SEXTANT 试验[8]。我国自2004年开始逐渐重视X射线脉冲星导航,近年来利用硬X射线调制望远镜(HXMT/慧眼卫星),天宫二号搭载的“天极”伽马暴偏振探测仪(POLAR)和脉冲星试验01星(XPNAV-1)的实测数据开展相应研究,取得了一系列成果[9-13]。

脉冲星试验01星是为建立脉冲星导航试验型数据库、检验脉冲星导航相关技术而研发,于2016年11月10日发射成功[14]。XPNAV-1搭载了两套不同的X射线探测器开展在轨试验,一套是大面阵微通道板(microchannel plate,MCP)探测器,另一套是Wolter-I掠入射聚焦型探测器。通过在轨技术试验,验证两种不同类型探测器的性能,为将来探测器的选型、改进等提供依据。此外,通过长期空间观测获得实测数据,为研究脉冲星物理特性[15-16]和脉冲星导航[17-18]提供数据支持。

蟹状星云脉冲星(Crab pulsar,简称Crab)是目前发现的X 射线波段辐射最强的脉冲星,它是公元1054年超新星爆发的产物,属于年轻的脉冲星,内部物理过程仍然活跃,导致其自转稳定性较差,难以用一组自转参数来准确地描述Crab的自转。因此,随着时间的推移,需要不断更新自转参数。本文采用了Jodrell Bank发布的Crab脉冲星星历[19],其中自转参数通常一个月更新一次,而当Crab出现周期跃变(glitch)时,星历的更新更加频繁。经过更新产生的各组星历参数内部自洽,相互之间存在自洽性问题,因此标准脉冲轮廓也需要更新。文献[20]提到了“使用不同星历的观测数据段分别构建各自的标准脉冲轮廓”的想法,这是针对星历更新带来的自洽性问题的一个考虑,不过该文献使用的Crab星历只经历了一次更新,因此没有开展进一步的研究。为了讨论星历更新引起的自洽性问题,本文将利用Wolter-I探测器观测Crab的两年数据进行计时分析,并据此改进Crab脉冲星长期(经历多次星历更新)观测数据的计时处理方法。

1 脉冲星试验01星两年观测数据计时分析

1.1 观测时长分析

脉冲星试验01 星两年观测数据涉及2016 年11月28日至2018年12月30日,总共观测了614次,其中观测时长约为50 min的数据文件近半。通过统计每次观测的时长和获得的光子数,可以得知观测时长与获得的光子数符合线性关系如图1 所示,平均光子计数率为14.8 counts/s。01星的观测数据存在少部分文件记录的光子数过少的情况,难以折叠出积分脉冲轮廓,体现为噪声信号。利用这些数据得到的计时残差及脉冲到达时间(time of arrival,TOA)不确定度非常大,可以通过3σ原则筛选,图1中的粗差数据即为筛选得到的结果。从粗差数据的分布反映出粗差的来源分为两种。一种是由于观测时间过短,探测的光子数较少导致折叠轮廓的信噪比(signal-to-noise ratio,SNR)过低;另一种是由于某些原因使得探测器的平均光子计数率降低,导致探测的光子数不如预期。在探测器正常工作且没有其他意外干扰下,为了保证观测质量,折叠信噪比较高的积分脉冲轮廓,需要限制观测时长的下限。图1 中红色虚线为观测时长等于35 min的界限,红色虚线右边部分的数据中,除了不在线性拟合曲线附近的粗差数据,其他数据中不存在粗差,说明观测时长大于35 min得到的光子数足以折叠出信噪比较高的积分脉冲轮廓。

图1 XPNAV-1观测时长和光子数的关系图Fig.1 The relationship between XPNAV-1 observation time and photon number

1.2 计时残差与定位精度

为了得到TOA,采用Taylor-FFT 法[21]在频域对标准脉冲轮廓与观测积分脉冲轮廓进行互相关计算。折叠脉冲轮廓时要选取恰当的相位子间隔数,根据文献[22]对XPNAV-1一年数据的分析,选取相位子间隔数Nbins=128 和Nbins=256 比较合适,并且两种选择的结果基本一致。为了更高的脉冲轮廓信噪比,避免引入更多噪声,本文选择相位子间隔数Nbins=128折叠脉冲轮廓。通过01星两年观测数据,根据3σ原则去除了74个数据点后,最终得到540个TOA。计时残差如图2所示,计时残差的均方根(root mean square,RMS)为70.90μs。

图2 XPNAV-1两年数据的TOA的残差,相位子间隔数N bins=128Fig.2 Residuals of TOAs of XPNAV-1 two-year data,the number of phase bins N bins=128

通常来说,航天器的定位精度σp与拟合前计时残差的RMS成正比,与TOA 个数NTOA的平方根成反比,即

因此,可以利用RMS和TOA 个数的平方根的比值来衡量定位精度。对于01星两年观测数据而言,该比值为3.05。

对计时残差进行线性拟合,其斜率可以估计星载原子钟的相对频率偏差kf,频率偏差的不确定度σkf衡量其估计的精度水平。01星的星载钟可以认为不存在频率偏差,因此估计其频率偏差kf没有意义,但可以评估估计水平。经计算,σkf=2.35×10-15,而对于目前的星载原子钟而言,其频率准确度通常在10-12量级[23],因此能够利用01星两年观测数据对星载原子钟进行频率修正。

2 脉冲星试验01星两年观测数据分组处理

2.1 观测数据分组方案

由于原始观测数据存在观测时长不均匀、观测日期不均匀的情况,为了讨论星历更新引起的自洽性问题对TOA 的影响,需要重新分配观测数据。将01星两年数据按照特定的观测时长分配为若干数据文件,每个数据文件的观测时长相同。由于Crab脉冲星星历在不断更新,考虑到不同脉冲星星历之间存在的自洽性问题,将使用相同星历的观测数据段看作集合进行分配,避免混合使用不同星历的观测数据。为了避免单一观测时长分配的选择性,对观测数据进行多种观测时长分配,每个数据文件的观测时长最低40 min,以5 min为步长,最高为140 min。

使用重新分配的数据文件,采用两种方法来获得TOA。方法一是将积分脉冲轮廓与所有数据折叠得到的标准脉冲轮廓互相关;方法二是将积分脉冲轮廓和与星历同步更新的标准脉冲轮廓互相关,标准脉冲轮廓的每次更新是利用更新星历对应时间段的所有数据折叠得到。方法一的标准脉冲轮廓由于使用了所有数据,因此无法排除星历更新带来的影响,而方法二的标准脉冲轮廓经过更新,其本身不包含其他时段星历信息,因此能够排除星历更新带来的影响。由于某些时间段观测数据不足,方法二中的标准脉冲轮廓可能信噪比过低,利用这些标准脉冲轮廓得到的计时结果存在问题。为了保证更新标准脉冲轮廓的精度,利用信噪比作为判断依据,筛选信噪比过低的标准脉冲轮廓,并以方法一使用的标准脉冲轮廓作为代替。本文定义信噪比为[18]

式中,Nsignal为脉冲信号成分的光子数,Ntotal为脉冲轮廓的总光子数。若是以信号成分和总体的平均光子计数率与观测时长T来表示信噪比,则有

式(3)中的平均光子计数率可以看作常数,因此能够对信噪比与观测时长进行拟合。图3是利用方法二中的标准脉冲轮廓的信噪比与观测时长得到的拟合曲线,系数为0.146 s-1/2。根据式(3)可以计算出利用20个观测时长为2 400 s的积分脉冲轮廓折叠出的标准脉冲轮廓信噪比为32.01。据此本文以SNR＜32作为判断依据来筛选信噪比过低的标准脉冲轮廓。

图3 SNR 与观测时长的关系图Fig.3 Relationship between SNR and observation time

2.2 两种方法得到的计时残差和定位精度

经过方法一和方法二两种方法处理,得到了以每组观测时长为自变量,计时残差的RMS和TOA个数为因变量的曲线如图4所示。与未经分组的处理结果相比,方法一和方法二在RMS为70μs时获得的TOA 数要远远大于540 个,而在TOA 数为540个时得到的RMS也要远远小于70μs,说明时长均匀分配的观测数据计时效果更好。对比图4(a)和图4(b)的RMS可知,方法二得到的RMS比方法一得到的RMS平均降低了14.2%,说明更新Crab脉冲星星历引起的自洽性问题影响了计时残差的精确度,而方法二能够有效降低计时残差的RMS,提高计时精度。因此,对于Crab脉冲星这种年轻的脉冲星,由于其自转参数变化相对较大,脉冲星星历需要不断地更新,如果要对其长期观测数据进行计时分析,方法二是比方法一更合理和有效的改进方法。

图4 两种方法处理得到的计时残差的RMS 和TOA的个数Fig.4 RMS of the timing residuals and number of the TOAs obtained through two methods processing

利用计算得到的RMS和TOA 个数,可以衡量不同观测时长下的方法一和方法二的定位精度。图5的结果反映了定位精度水平,与未分组的结果相比,分组后方法一和方法二的RMS和TOA 个数的平方根的比值皆低于3.05,说明均匀的观测时长能有效提高定位精度。虽然较短的观测时长折叠出的积分脉冲轮廓信噪比较低,计时残差的RMS偏高,但是获得了较多的TOA,整体的定位精度得到了提高。此外,与方法一相比,方法二优化了脉冲星星历的自洽性问题,定位精度平均提高了14.4%。

图5 RMS 和TOA个数的平方根的比值与每组观测时长的关系Fig.5 The relationship between the ratio of RMS to the square root of the number of TOAs and the observation time per group

图6反映了相对频率偏差的不确定度σkf与每组观测时长的关系。随着每组观测时长的增加,获得的TOA 数量减少,相对频率偏差的不确定度增大,说明获得更多的TOA 能够提高频率偏差的估计水平。与方法一相比,方法二将相对频率偏差的不确定度平均降低了26.2%,进一步体现了星历更新导致的自洽性问题和方法二处理长期观测数据的合理性和优越性。

图6 相对频率偏差的不确定度与每组观测时长的关系Fig.6 The relationship between the uncertainty of relative frequency deviation and the observation time per group

2.3 稳定度分析

当描述原子时的稳定性时,通常使用Allan方差或Hadamard方差。对于脉冲星计时而言,由于脉冲星观测数据间隔的随机性,不能直接使用原子时稳定度的分析方法。1997 年,Matsakis等[24]提出一种改进的σz方差,直接利用长度为τ的脉冲星计时残差序列来拟合一个三阶多项式,并利用三次项系数c3确定σz

式中,尖括号表示以c3不确定度的平方倒数为权做加权平均。利用01星的两年数据,分别给出了方法一和方法二对应的稳定度曲线如图7(a)和图7(b)所示,其中对照组为使用未进行数据分组的原始数据以方法一处理得到的结果。通过图7(a)和图7(b)可以看到,经过方法一和方法二处理得到的σz量级上并无改变,对稳定度的提高贡献不大。01星的2年观测数据的时间稳定度能达到10-12,但与一些稳定的X 射线毫秒脉冲星相比,其1年的观测就能达到这个水平[18]。因此,稳定的X 射线毫秒脉冲星是更具有潜力的研究对象。

图7 两种方法处理得到的时间稳定度Fig.7 Timing stability obtained through two methods processing

3 结论

本文通过研究,得到如下结论:

1)对于不同的观测时长分组情况,方法二得到的计时残差的RMS均明显低于方法一,平均降低了14.2%,说明方法二明显改进了计时残差的精确度。

2)与方法一相比,方法二不会降低获得的TOA数,因此方法二的定位精度均明显高于方法一,平均提高了14.4%。

3)方法二估计原子钟相对频率偏差的精度均高于方法一,平均提高了26.2%。

4)方法二对01星两年观测数据的时间稳定度没有量级上的改进,对标准轮廓更新的改进方法基本不影响时间稳定度。

以上几点反映了对于经历多次星历更新的Crab脉冲星X 射线观测数据,方法二的计时改进效果显著,有效提高了计时残差的精确度,显著提升了Crab脉冲星在航天器定位、星载原子钟频率校准等应用领域的性能,体现了改进方法在处理此类长期观测数据的合理性和优越性。尽管Crab脉冲星不是毫秒脉冲星,星历需要频繁更新,但是通过射电波段的密集观测可以准确测定其星历参数。而由于Crab脉冲星的X 射线波段非常明亮,易于被深空航天器携带的X 射线探测器探测到,那么通过注入射电波段观测的星历参数并结合X 射线波段的实测数据,可以校准深空航天器时间系统的频率偏差,提高其自主保持能力。