脉冲星导航变封装段的轮廓恢复方法

1.本发明涉及导航方法领域，具体涉及脉冲星导航变封装段的轮廓恢复方法。


背景技术：

2.x射线脉冲星导航是一种新兴的航天器自主导航方法，具有脉冲星导航适用范围广、导航稳定性高、导航精度不会因所处位置的不同而改变等传统天文导航不具备的优势，且其导航精度与探测器相关，在未来还有很大的发展空间。同时，脉冲星信号还具有周期稳定度高，在提供导航信息的同时还可以提供时间基准信息的特点。
3.脉冲轮廓不仅是脉冲星的辨识标准，其作为脉冲星信号处理过程的重要步骤，可以直接应用于脉冲到达时间估计及求解脉冲星导航的基本观测量，也可通过长时间观测构建高精度的脉冲星轮廓模板。历元折合作为目前基本的脉冲轮廓恢复方法已被广泛应用，利用航天器轨道信息对光子到达时间序列进行修正后，计算每个光子到达时间相对于参考历元(轮廓折叠基点)的相位，并将其折合到一个脉冲周期内所划分的各个封装段(bin)中，并统计每个封装段中的光子数，即可恢复累积的脉冲轮廓。
4.脉冲星轮廓恢复是脉冲星导航及信号处理过程中的一项关键技术。对于现有方法，由于历元折合设置的封装段是不发生变化的，想要恢复高精度的脉冲轮廓，提高轮廓信噪比，需要延长脉冲星的观测时间以获取更多的光子到达时间，相应折合到一个脉冲周期内每个封装段里的光子数也会增多，但太长观测时间将影响脉冲星导航的实时性和导航精度，还会增加系统的计算负担。在观测时间较短时，观测数据量少，恢复的脉冲轮廓信噪比低，则会影响脉冲到达时间的估计精度，进而影响导航精度。


技术实现要素：

5.因此，针对上述问题，提出一种能以较短的脉冲星观测时间恢复得到高精度脉冲轮廓的方法具有重要意义。为了实现上述目的，本发明提供了一种脉冲星导航变封装段的轮廓恢复方法，具体技术方案如下：
6.具体的，所述方法步骤如下：
7.步骤s1：光子数据质心修正，具体是：基于时间转换模型，将脉冲星观测间隔内的全部光子到达时间序列转换为光子传播到达太阳系质心处的时间；
8.步骤s2：脉冲星信号周期搜索，具体是：光子数据质心修正之后，变化不同的试验周期进行历元折合，通过检验函数，得到脉冲星信号的最佳折合周期；
9.步骤s3：光子数据相位信息计算，具体是：根据脉冲星信号的相位传播模型，计算观测间隔内全部光子数据的相位信息，并取余数，保留各光子相位信息的小数部分；
10.步骤s4：变bin轮廓折叠，得到恢复后的脉冲轮廓，具体是：光子相位信息排序；bin信息计算，设置轮廓恢复bin数，计算每个bin中所包含的光子数，利用各bin中的首尾光子的相位信息计算bin宽度；恢复轮廓坐标取各bin的中值点，并计算相应的轮廓高度系数；根据计算得到的bin信息恢复变bin脉冲轮廓。
11.进一步地，所述封装段为在一个脉冲周期内划分而成在历元折合过程中的光子存储空间，通过统计封装段中的光子数，恢复累积的脉冲轮廓。
12.优选的，所述步骤s1中的时间转换模型如下：
[0013][0014]
其中，t
ssb
是光子传播到达太阳系质心的时间，t
sc
是令航天器处探测到光子的时间，n是脉冲星位置矢量，r
sc
是航天器相对于太阳系质心的位置矢量，c是光速，d0是脉冲星在基准传播时间t0时的位置，b是太阳系质心相对于太阳质心的位置矢量，pk是航天器相对于第k颗行星的位置矢量，μk为第k颗行星的引力常数，ns为考虑影响的行星总数。
[0015]
优选的，所述步骤s2中的检验函数根据卡方检验构建，所述检验函数如下：
[0016][0017]
其中，为各bin的平均光子数，nb为轮廓恢复的bin数，χ2为检验函数；
[0018]
χ2的最大值对应的折合周期为脉冲星信号的最佳折合周期，表达式如下：
[0019][0020]
其中为脉冲星信号的最佳折合周期。
[0021]
优选的，所述步骤s3中的相位传播模型如下：
[0022][0023]
其中，t表示任意光子到达时间，t0表示恢复基点时刻，f0、f1、f2分别表示t0时刻脉冲星信号自转频率、自转频率一阶导数和自转频率二阶导数，在仅考虑周期的情况下，
[0024]
进一步地，所述步骤s4包括以下步骤：
[0025]
步骤s401，光子相位信息排序，具体是：将所述步骤s3获得的光子相位信息取余数操作，并按照从小到大的顺序对取余数后的光子相位信息进行排序；
[0026]
步骤s402，bin信息计算，具体是：通过信息准则估计脉冲星信号轮廓折叠的理论最佳bin数，对每个bin中所包含的光子数进行计算，并从排序后的光子相位信息序列里将相应光子填充进每个bin中，利用填充入各bin中的首尾光子的相位信息计算封装段宽度；
[0027]
步骤s403，计算恢复脉冲轮廓坐标，具体是：根据步骤s402计算得到的bin信息计算恢复轮廓坐标，所述恢复轮廓坐标取各bin内光子相位信息的中值点，并计算相应的轮廓高度系数统计各bin内累计光子数；
[0028]
步骤s404，绘制恢复后的脉冲轮廓，具体是：根据步骤s403计算出的恢复轮廓坐标绘制恢复后的脉冲轮廓图。
[0029]
优选地，所述步骤s402中的信息准则为贝叶斯信息准则，其基本原理表达式如下：
[0030][0031]
其中，表示恢复轮廓，nb表示bin的个数，tb表示每个bin的长度，ci是第i个封装段内的光子总数，ti为第i个封装段的中间时刻，b表示理论最佳bin数函数；
[0032]
当b为最小值时的bin数作为轮廓折叠的理论最佳bin数，形如下式：
[0033][0034]
其中，为轮廓折叠的理论最佳bin数。
[0035]
优选的，所述步骤s402中计算bin宽度的公式如下：
[0036][0037]
其中，wi为首尾光子的相位信息计算封装段宽度；
[0038]
为第i个封装段末尾填充光子的到达时间；为第i个封装段起始填充光子的到达时间；其中，i为封装段的编号，i＝1，2，...，nb；φ()表示利用s3所述的相位传播模型对相应光子的相位信息进行计算。
[0039]
优选的，所述步骤s403中计算恢复轮廓坐标的表达式如下：
[0040][0041]
其中，xi为恢复后的脉冲轮廓横坐标；yi为恢复后的脉冲轮廓纵坐标。
[0042]
应用本发明的技术方案，具有以下有益效果：
[0043]
本发明提出的轮廓恢复方法不固定历元折合过程中封装段的宽度，按照光子序列的相位信息将光子分配到不同的封装段内，每个封装段内的光子数量相同。因此，恢复的轮
廓光子数量较少的部分，其封装段的宽度较大，反之，光子数越多的部分，封装段宽度越小，即能使恢复轮廓的波峰特征更加明显。
[0044]
本发明提出的轮廓恢复方法利用计算得到的轮廓高度系数代替传统方法中对各封装段内累计光子数的统计，来恢复脉冲轮廓。
[0045]
本发明提出的轮廓恢复方法恢复的脉冲轮廓可以使后续脉冲到达时间估计的结果脱离轮廓封装段的限制，对比相位差的宽度可以直接利用光子相位信息进行计算，其精度与光子相位信息的计算精度相关。
[0046]
除了上面所描述的目的、特征和优点之外，本发明还有其它的目的、特征和优点。下面将参照图，对本发明作进一步详细的说明。
附图说明
[0047]
构成本技术的一部分的附图用来提供对本发明的进一步理解，本发明的示意性实施例及其说明用于解释本发明，并不构成对本发明的不当限定。在附图中：
[0048]
图1是历元折合原理图；
[0049]
图2是本发明实施例中轮廓恢复方法的步骤流程图；
[0050]
图3是本发明实施例中
[0051]
图4是本发明实施例中的恢复轮廓对比图；
[0052]
图5是本发明实施例中的基于历元折合的恢复轮廓模型图(bin数为1000)；
[0053]
图6是本发明实施例中的基于历元折合的恢复轮廓模型图(bin数为10000)；
[0054]
图7是本发明实施例中的相位差理论计算精度结果对比图。
具体实施方式
[0055]
以下结合附图对本发明的实施例进行详细说明，但是本发明可以根据权利要求限定和覆盖的多种不同方式实施。
[0056]
作为常用的脉冲轮廓恢复方法，传统历元折合的基本原理如下：
[0057]
设观测时长为t
obs
，含有n
p
个完整的脉冲星周期，即n
p
＝floor(t
obs
/p)，其中，p为脉冲星周期，所述p可划分为nb个封装段，每个封装段的长度为tb。具体为：将记录在后续周期的光子到达时间序列折合到第一个周期；计算每个封装段中的光子数ci(i＝1，2，...，nb)；通过归一化光子数，得到恢复轮廓
[0058]
归一化公式如下：
[0059][0060]
[0061]
其中，表示恢复轮廓；n
p
表示完整的脉冲周期个数，其表达式为n
p
＝floor(t
obs
/p)，其中，t
obs
为观测时长，p为脉冲星周期，所述p可划分为nb个封装段；tb表示每个封装段的长度；cj(ti)为第j个周期折合入第i个封装段的光子数；ci是第i个封装段内的光子总数；ti为第i个封装段的中间时刻。
[0062]
如图1所示的历元折合原理图，传统历元折合方法，其封装段的宽度是固定的，在轮廓折叠bin数确定后，其宽度就不再发生变化。脉冲轮廓在各封装段的高度是由累积在封装段中的光子数量来体现的。
[0063]
针对现有脉冲轮廓恢复技术方案所存在的问题，本实施例提出了脉冲星导航变封装段的轮廓恢复方法。其中，所述封装段为在一个脉冲周期内划分而成在历元折合过程中的光子存储空间，通过统计封装段中的光子数，恢复累积的脉冲轮廓。如图2所示，所述轮廓恢复方法具体包括以下步骤：
[0064]
步骤s1：光子数据质心修正。基于时间转换模型，将脉冲星观测间隔内的全部光子到达时间序列转换为到达太阳系质心处的时间；
[0065]
步骤s2：脉冲星信号周期搜索。变化不同的试验周期进行历元折合，依据统计学原理，通过相应的检验函数，得到脉冲星信号的最佳折合周期；
[0066]
步骤s3：光子数据相位信息计算。根据脉冲星信号的相位传播模型，计算观测间隔内全部光子数据的相位信息，并取余数，保留各光子相位信息的小数部分；
[0067]
步骤s4：变bin轮廓折叠，得到恢复后的脉冲轮廓，具体是：光子相位信息排序；bin信息计算，设置轮廓恢复bin数，计算每个bin中所包含的光子数，利用各bin中的首尾光子的相位信息计算bin宽度；恢复轮廓坐标取各bin的中值点，并计算相应的轮廓高度系数代替传统方法中对各bin内累计光子数的统计；根据计算得到的bin信息恢复变bin脉冲轮廓。
[0068]
本实施例提出的方法的具体实施过程如下：
[0069]
具体的，所述步骤s1为光子数据质心修正。实际在轨运行时，航天器的精确轨道信息未知，且已知探测器接收到的脉冲星信号是它在轨道不同位置处接收到的，因此需要利用预估轨道信息将航天器探测的光子到达时间全部转换到光子到达太阳系质心处的时间。
[0070]
进一步地，考虑到光子到达时间在传播过程中所受到的roemer延迟和shapiro延迟，假设脉冲星的视线方向位置在太阳系中不变，并忽略脉冲星自行及公式中存在高阶项的影响，相应的时间转换模型可如下式表示：
[0071][0072]
其中，t
ssb
是光子传播到达真空太阳系质心的时间；t
sc
是航天器处探测到光子的时间；n是脉冲星位置矢量；r
sc
是航天器相对于太阳系质心的位置矢量；c是光速；d0是脉冲星在基准传播时间t0时的位置；b是太阳系质心相对于太阳质心的位置矢量；pk是航天器相对于第k颗行星的位置矢量；μk为第k颗行星的引力常数；ns为考虑影响的行星总数。
[0073]
通过计算式(3)所示的时间转换模型，可以得到t
ssb
，即完成光子数据的质心修正。
[0074]
具体的，步骤s2为脉冲星信号周期搜索。在光子数据质心修正之后，为恢复得到准确的脉冲观测轮廓，需要确定脉冲星信号的最佳折合周期，即通过周期搜索获取航天器接收脉冲星信号的准确周期信息。
[0075]
进一步地，周期搜索的基本原理如下：
[0076]
基于脉冲星的先验周期信息设置周期搜索范围，利用范围内不同的试验周期p进行历元折合，得到每个封装段内的光子数ci，其中，i为bin的编号，i＝1，2，...，nb，根据卡方检验建立如下式所示的检验函数：
[0077][0078]
其中，ci为每个bin内的光子数，其中，i为bin的编号，i＝1，2，...，nb；为各封装段的平均光子数；nb为轮廓恢复的bin数；χ2为检验指标函数。
[0079]
当检验指标函数达到最大时，即χ2最大时，其对应的即为估计得到的脉冲星信号的最佳折合周期，表达式如下：
[0080][0081]
具体的，步骤s3为光子数据相位信息计算。恢复基点时刻的脉冲轮廓需要计算所有基于该点的光子到达时间所在的相位。步骤s1中已经将航天器处的光子到达时间转换到太阳系质心处，相应的相位信息可以利用太阳系质心处的时间相位模型外推获得，所述时间相位模型具体可表示为如下形式。
[0082][0083]
其中，t表示任意光子到达时间，t0表示恢复基点时刻，f0、f1、f2分别表示t0时刻脉冲星信号自转频率、自转频率一阶导数和自转频率二阶导数，在仅考虑周期的情况下，
[0084]
进一步地，通过式(6)计算观测时间内全部光子数据基于轮廓恢复基点的相位信息。
[0085]
具体的，步骤s4变bin轮廓折叠，恢复得到脉冲轮廓；如图3所示，具体步骤如下：
[0086]
步骤s401，光子相位信息排序。
[0087]
具体的，将步骤s3获得的光子相位信息取余数操作，并按照从小到大的顺序对取余数后的光子相位信息进行排序；
[0088]
步骤s402，bin信息计算。
[0089]
具体的，可以通过经验直接设置轮廓恢复的bin数；也可以通过信息准则来估计脉冲星信号轮廓折叠的理论最佳bin数，如一种基于贝叶斯信息准则的方法，其基本原理表达
式如下：
[0090][0091]
其中，b为理论最佳bin数函数。
[0092]
进一步地，将使式(7)达到最小值的bin数作为轮廓折叠的理论最佳bin数，即
[0093][0094]
其中，为轮廓折叠的理论最佳bin数。
[0095]
进一步地，对每个bin中所包含的光子数进行计算，并从排序后的光子相位信息序列里将相应光子填充进每个bin中。
[0096]
cb＝floor(n/nb)
ꢀꢀꢀ
(9)
[0097]
其中，n为光子总数，cb为填充进每个bin的光子数。
[0098]
进一步地，利用填充入各bin中的首尾光子的相位信息计算封装段宽度wi，其中，i为bin的编号，i＝1，2，...，nb。
[0099][0100]
其中，wi为首尾光子的相位信息计算bin宽度；为第i个封装段末尾填充光子的到达时间；为第i个封装段起始填充光子的到达时间；其中，i为封装段的编号，i＝1，2，...，nb；φ()表示利用s3所述的相位传播模型对相应光子的相位信息进行计算；i为bin的编号，i＝1，2，...，nb。
[0101]
步骤s403，计算恢复脉冲轮廓坐标。根据计算得到的bin信息，即可恢复变bin脉冲轮廓。
[0102]
具体的，恢复轮廓坐标取各bin内光子相位信息的中值点，并计算相应的轮廓高度系数代替传统方法中对各bin内累计光子数的统计，具体可表示为如下形式：
[0103][0104]
其中，xi为恢复后的脉冲轮廓横坐标；yi为恢复后的脉冲轮廓纵坐标。
[0105]
步骤s404，根据计算出的横纵坐标即可绘制出恢复后的脉冲轮廓图。
[0106]
本实施例提出的轮廓恢复方法不固定历元折合过程中封装段的宽度，能以较少的数据量恢复得到信噪比较高的脉冲轮廓，且无需额外的降噪处理，原理及步骤简单，可以直
接应用于在轨信号处理过程。在相同光子到达时间序列的情况下，本方法不固定封装段的数量，选择最佳封装段能避免在轮廓恢复过程中因设置封装段数过大而加大历元折合噪声的方差，分散光子的分布；封装段数过小则会出现恢复轮廓过拟合的问题。本方法能充分利用观测数据，恢复轮廓的波峰特征更加明显，最终提高脉冲到达时间的估计精度。另外，利用本发明恢复的脉冲轮廓可以使后续计算得到的相位差脱离封装段的限制，可以有效提高相位差估计值的精度，在数据量越少的情况下，提升的效果越明显。
[0107]
为了说明本实施例的目的和优势，本实施例选择nicer探测器2018年11月的在轨实际观测数据，具体观测号为(obs_id：1013010143)，观测目标为脉冲星crab。
[0108]
选取其中一段持续时间仅为800s左右的观测数据，分别应用本实施例提出的轮廓恢复方法以及传统历元折合方法折叠脉冲轮廓，2种方法恢复的轮廓对比如图4所示。
[0109]
从图4中可以看出，本实施例方法恢复的轮廓特征更加明显，峰值更加尖锐；在主波峰处以及次波峰处，轮廓相比与传统历元折合方法恢复的轮廓相比均偏窄。由于脉冲轮廓是多周期折叠得到的，这一效果在图中表现并不十分明显。在光子数据较少的情况下，恢复轮廓的波峰宽度能减少约5％-10％左右，进而可以提高脉冲到达时间的理论估计精度。
[0110]
对于传统历元折合方法，在轮廓折叠过程中，设置不同的bin数，会对恢复的轮廓造成较大的影响。
[0111]
如图5所示的bin数为1000的轮廓模型以及如图6所示的bin数为10000的轮廓模型，对于传统历元折合方法，后续轮廓与模板之间相位差的计算结果与相位信息的分辨率相关，即与封装段的宽度相关。因此，增大轮廓折叠bin数，可以提高相位差计算的理论精度。在光子数量相同的情况下，bin数的增加会加大历元折合噪声的方差，分散光子的分布，使恢复的轮廓信噪比较低，同样会造成相位差计算精度的下降，导致相位差的计算精度无法达到更高的水平。
[0112]
本实施例提出的轮廓恢复方法以及传统历元折合方法折叠脉冲轮廓分别与脉冲轮廓模板进行对比，相应的相位差理论计算精度结果如图7所示。
[0113]
以取光子相位信息计算结果小数点后的4为有效数字为例，从图7中可以看出，应用轮廓恢复方法的后续相位差结果的计算精度仅与利用相位传播方程计算得到单光子相位精度相关，不受折叠bin数的设置。
[0114]
传统历元折合的方法轮廓恢复的效果与人为设置的bin数直接相关，受其制约；在光子数量充足的情况下，在bin数设置为10000时，其后续相位差计算的理论精度才可以与本实施例方法相当。
[0115]
因此，在光子数量较少的情况下，本实施例提出的轮廓恢复方法的效果要远优于传统历元折合方法。
[0116]
以上所述仅为本发明的优选实施例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。