X射线源/地心矢量观测的卫星自主导航应用

李璟璟, 张玉兔, 王文从, 胡慧君, 邵思霈

(山东航天电子技术研究所, 山东 烟台 264000)

0 引 言

自20世纪60年代X射线源天体被发现以来,基于X射线源的天文导航技术成为实现航天器自主导航的又一研究热点。2016年11月,我国首次发射了以X射线导航为目的的X射线脉冲星探测卫星,获取了宝贵的在轨数据。2017年6月,长四乙发射的“慧眼”卫星,也搭载了X射线探测载荷,通过在轨获取X射线数据,同样能够为X射线导航提供有力支持。X射线脉冲星导航方法主要采用脉冲到达时间(time of arrival, TOA)观测的手段,该方法具有完全自主、全信息、高精度等优点,但是也存在脉冲折叠时间长、探测器面积大、TOA估计精度受限等缺点。此外,要实现完全的自主导航,还需要通过相当长时间的脉冲星观测,建立脉冲星数据库[1-4]。相比于TOA观测技术,X射线成像技术不需要光子累积时间,探测器面积小。此外,与传统光学波段敏感器相比,X射线成像探测器成本低、重量轻,星图识别速度快,既可以作为现有导航手段的补充,也可独立提供高精度的导航信息[5-7]。

结合X射线新型导航手段的优点,将X射线探测器应用于自主天文导航[8]。提出基于X射线源/地心矢量观测的航天器自主天文导航方法。建立了自主导航观测模型,并针对系统观测残差进行了残差分析和建模,利用UKF滤波器,实现了航天器自主导航信息的解算。

1 基于X射线源/地心矢量观测的自主导航方法

基于X射线源/地心矢量观测的自主导航原理如图1所示。依据卫星、X射线导航星以及地球三者的几何关系,可以得到的角度观测量有星光角距α、地球视角β以及星光仰角γ。星光角距是从卫星上观测到的导航X射线源的星光矢量方向与地心矢量方向之间的夹角。地球视角是从卫星上观测地球边缘构成的圆锥角。星光仰角是从卫星上观测的导航X射线源与地球边缘的切线方向之间的夹角。结合卫星轨道动力学模型,并运用最优状态估计方法,可以得到航天器自主导航信息。

图1 X射线源/地心矢量观测自主导航原理Fig.1 Autonomous navigation based on X-ray/geocentricvector observation

基于X射线源/地心矢量观测的航天器自主导航系统主要由探测器单元和导航计算单元构成。其中探测器单元包括X射线探测器和地球敏感器,用于获取二者矢量信息。星载导航计算单元包括观测量合成模块、星历数据库和最优滤波模块。系统基本方案如图2所示。

图2 自主导航系统方案Fig.2 Autonomous navigation system scheme

两种探测器方案可以实现X射线源矢量方向的测量。第1种是编码孔径盘像素探测器。通过对天空X射线源进行成像,获取视野中的X射线源特征、轮廓等信息,进而获取X射线源在惯性坐标系下的矢量方向。研究表明,利用有效面积为650 cm2的X射线探测器观测SC0-X1射线源,只需1 ms便可以获取可识别的图像[5]。第2种方案类似于传统的星扫描仪。通过准直器和机械扫描的方法扫描已知X射线源,通过一定的搜索算法确定X射线源单位矢量方向[9]。

2 观测方程建模

2.1 观测方程的构建

以X射线源为观测对象,通过星图识别,获得X射线源的单位方向矢量。其方位信息来自像平面坐标系中的X射线源矢量Sp。设(xc,yc)为像平面坐标系中的像点坐标,fc为成像焦距,则

(1)

(2)

通过地球敏感器获取地平信息,并结合地球物理参数信息,计算得到地心在地球敏感器基准坐标系中的矢量。利用红外地球敏感器输出的地平扫入脉冲、地平扫出脉冲以及磁基准脉冲。通过几何关系,可以得到地球敏感器坐标系中地心方向的矢量为

(3)

卫星通过X射线源探测器观测得到X射线源方向矢量,通过地球敏感器敏感地平并计算获得地心矢量。从而进一步合成得到关于两者夹角的实际观测量。为了进行状态估计,还需要对该实际观测量进行数学上的表达。由图1所示的几何关系可知,X射线源方向矢量和地心矢量的夹角α为位置的函数,可建立观测方程为

(4)

式中,rc为J2000地心惯性系下卫星的位置矢量;s为卫星到X射线源的单位矢量;va为测量噪声,通常为高斯白噪声。

如果以星光仰角为观测矢量,则可直接由地球敏感器获得地平方向矢量,与星光方向矢量共同构建观测方程

(5)

式中,Re为地球半径。

2.2 测量残差分析

将系统观测残差建模为加性误差,则式(4)～式(5)表示的观测方程可以统一写为以下形式,即

Yob=Y+ΔY+vob

(6)

式中,Yob为实际观测量;Y为真实量;ΔY为X射线源方向矢量观测误差和地球敏感器方向矢量观测误差引起的观测量残差的和;vob为测量噪声,服从零均值高斯分布。

根据式(4),令

M=-(rc·s)/rc=-r·s

(7)

式中,r为探测器在地心惯性坐标系中的单位位置矢量。由X射线源方向矢量观测误差和地球敏感器方向矢量观测误差引起的点积偏差可表示为ΔM,有

ΔM=ΔMs+ΔMe

(8)

式中,ΔMs为X射线源方向矢量观测误差引起的点积偏差,ΔMe为地球敏感器方向矢量观测误差引起的点积偏差。由ΔM引起的观测量偏差为ΔY,令观测量Y′为

Y′=Y+ΔY=arccos(M+ΔM)

(9)

真实量为

Y=arccos(M)

(10)

可得

ΔM=cos(Y+ΔY)-cos(Y)

(11)

式(11)进行麦克劳林展开,有

(12)

略去高阶小项,整理得

(13)

根据式(7),对s求偏导数,可得到X射线源观测误差引起的观测量偏差ΔMs为

ΔMs=-r·Δs

(14)

即由X射线源方向矢量角位置偏移Δs引起的误差传播函数。当Δs为角秒量级时,ΔMs极小。上式对时间求导,并取绝对值,得

(15)

近地圆轨道卫星速率不超过8 km/s,轨道高度一般为200～800 km,X射线源天体角位置误差一般为角秒量级。由此估算,|ΔMs|一小时变化率约为10-5量级。可见,X射线源角位置误差引起的ΔMs极小且变化极为缓慢。

同理,式(7)两端对r求偏导数,得

ΔMe=-Δr·s

(16)

即由地心矢量角位置偏移Δr引起的误差传播函数。上式对时间求导,并在等式两边取绝对值,可得地心方向矢量观测误差引起的单位方向矢量点积偏差ΔMe为

(17)

单圆锥扫描地球敏感器在轨测量误差约为0.25°～1°,双圆锥扫描地球敏感器在轨测量误差约为0.05°～0.07°,面阵地球敏感器在轨测量误差约为0.02°～0.05°,当测量误差为常值误差或缓变量时[10],|ΔMe|的变化率极小。对于近地圆轨道卫星,当采用敏感器光轴在轨道平面内与地平切线垂直的最佳UKF滤波安装方式时[11,12],真实观测量Y约为90°+arcsin(Re/rc),由于ΔM为缓变量,根据式(13),可将ΔY建立为缓变模型,即

(18)

式中,ws为高斯白噪声。

3 状态方程建模

卫星轨道运动学采用牛顿二体引力模型[13],同时考虑地球非球形引力及其他摄动力的作用,卫星运动学模型建立为

(19)

式中,v为J2000地心惯性系下的卫星速度矢量;μ为地球引力常数;P(rc,v)为受摄运动项。将上式展开,并利用勒让德多项式化简,整理后即可得到卫星运动状态方程

(20)

式中,状态变量取为x=[x,y,z,vx,vy,vz]T,分别为卫星在J2000地心惯性系中3个坐标轴上的位置和速度,J2为地球引力二阶带谐项系数。在这里,考虑到J2摄动项相对其他摄动项为大项,因而在方程中单独提出,ΔFx、ΔFy、ΔFz为地球非球形摄动的高阶摄动项和日、月摄动以及太阳光压摄动和大气摄动等其他摄动力影响。wo为过程噪声,假定为高斯白噪声。

为了补偿观测量偏差引起的系统估计误差,提高导航参数估计精度,根据系统测量残差的分析,将X射线源矢量观测误差和地心矢量观测误差引起的系统残差扩展补充到系统状态方程中,选取状态量X=[xΔY]T,建立扩展状态方程为

(21)

式中,w=[wows]T为扩展状态噪声,设为零均值高斯白噪声。

4 导航算法

式(4)、式(5)和式(21)构成的自主导航系统是复杂的非线性系统。针对这一类系统的最优滤波问题,Julier等提出了UKF滤波算法。该算法通过设计合理的采样Sigma点逼近非线性系统状态概率密度分布函数,在提高滤波精度的同时规避了雅克比矩阵的计算[14]。经过Unscented变换,UKF滤波可以获得三阶准确度的非线性系统函数统计量,在非线性动态滤波中得到越来越多的应用[15-16]。因此,本文应用UKF算法对自主导航系统状态量进行滤波计算,非线性系统的状态方程式(21)和观测方程式(6)经离散化后,可整理为以下形式，即

Xk+1=F(Xk)+Wk

(22)

Yob,k=G(Xk)+Vk

(23)

式中,k代表离散化采样点序列,状态向量为Xk∈RL,测量向量为Yob,k∈RM,L和M分别为状态向量维数和观测向量维数,Wk和Vk分别为状态噪声和测量噪声,满足

(24)

式中,Qk≥0为状态噪声方差阵;Rk>0为测量噪声方差阵;δkj为Kronecker符号。

，i=0

(25)

(26)

(27)

Sigma点在非线性系统中传播,并对其结果进行加权计算,得到状态变量和估计均方误差阵的一步预测,即

χi,k/k-1=F(χi,k-1)

(28)

(29)

(30)

重新生成Sigma点并预测观测向量为

(31)

(32)

(33)

计算增益矩阵为

(34)

(35)

(36)

对状态变量和估计均方误差阵进行更新,可得到k时刻导航信息为

(37)

(38)

(39)

式中,β是状态分布参数,β=2 对于高斯分布是最优的。

由于观测量获取时间间隔与导航计算周期往往是不同步的,因而需要基于导航周期最近估计点进行导航时刻的位速估计,从而获得满足导航与控制需求的实时位速信息。因此在导航周期中,进行两步独立的UKF计算,一方面根据观测量到达周期进行UKF导航估计,从而为下一周期的导航滤波提供状态估值;另一方面根据当前时刻t进行UKF导航估计,用于飞行器制导与控制决策。导航周期计算流程如图3所示。其中tm

图3 导航周期计算流程Fig.3 Calculation flow of navigation periods

5 仿真分析

为了验证本文提出方法在航天器自主导航应用中的有效性,以实现近地圆轨道卫星自主导航为目标进行仿真。

卫星的初始轨道六根数为

半长轴a=7 208.72 km,偏心率e=0.011,轨道倾角i=65°,升交点赤经Ω=30°,近地点角距ω=32°。

J2000地心惯性坐标系下初始三轴位置估计误差为1 km,初始三轴速度估计误差为10 m/s,初始系统残差估计误差为0.5°,滤波器初始位置、速度估计方差阵为

仿真中所采用的系统过程噪声矩阵为

Qk=diag[1 m 1 m 1 m 2×10-3m/s

2×10-3m/s 2×10-3m/s]2

测量残差初始估计方差为(7°)2,测量残差状态噪声为(0.001°)2,观测量噪声取为(4×10-6)°。

仿真中导航周期取为5 s。

参照实际敏感器能力,取X射线源矢量方向观测精度3″,地心矢量方向观测精度为0.02°,采用标准UKF方法进行导航计算。卫星在J2000地心惯性坐标系下的位置、速度估计误差随时间变化如图4和图5所示。

图4 位置估计误差Fig.4 Position estimation error

图5 速度估计误差Fig.5 Velocity estimation error

在经过一定时间的收敛过程后,滤波结果最终达到稳定,图6为测量残差的估计误差。

图6 测量残差估计误差Fig.6 Observation residual estimation error

从局部放大图可见,测量残差估计误差可达到10-4量级,这表明滤波器能够实现对测量残差的有效估计。考虑到滤波收敛时间,需要在导航方案设计中,合理设置自主导航滤波起始点,从而保证任务时间点的滤波精度。

按照以上条件进行蒙特卡罗仿真,并对稳定后的滤波结果求取均方根误差均值,三轴位置和速度估计误差如表1所示。可见,三轴位置估计误差约为100～200 m,三轴速度估计误差约为0.1～0.15 m/s。该精度水平能够满足航天器在轨自主导航的需求。

表1 均方根误差均值

Table 1 Mean of RMSE

6 结 论

本文以近地航天器自主天文导航为应用背景,提出一种基于X射线源/地心矢量观测的航天器自主导航方法,将空间X射线探测器拓展应用于航天器自主天文导航。进行了系统建模与分析,将X射线源矢量观测误差和地心矢量观测误差引起的系统观测残差建模为缓变过程,并扩展补充到状态方程中。利用UKF算法进行了导航信息解算。仿真表明,该方法能够满足航天器自主天文导航需求,具有一定的工程应用参考价值。