技术领域
[0001] 本
发明主要涉及到导航方法领域，特指一种高轨卫星脉冲星/星光角距组合导航方法。
背景技术
[0002] 高轨卫星（轨道高度高于20000km）在
气候探测、灾难预警、数据中继等方面扮演了重要的角色。精确的
位置信息是保证高轨卫星有效工作的重要因素。目前，高轨卫星的轨道确定工作大都基于地面测控站实现。然而，随着高轨卫星数目的增多，地面测控站的负担逐渐增加。同时，空间环境的日益复杂也对高轨卫星的自主性提出了迫切需求。拥有一个高可靠、高稳定的自主
导航系统是实现高轨卫星自主性的必要前提。
[0003] 然而，由于轨道高度高于导航卫星，高轨卫星难以同时接收到四颗以上的导航卫星
信号。因此，无法直接利用现有的全球导航卫
星系统来完成高轨卫星的自主导航工作。此外，基于星光角距的天文导航系统(Celestial Navigation System，CNS)是一种成熟的自主导航系统。该方法通过测量地球和参考
恒星的视线夹角来实现导航
定位，且适用于整个近地空间。该方法还具有导航敏感器体积小、测量信息
采样周期短的优点。然而，该方法的性能随着轨道高度的增加而降低。对高轨卫星而言，CNS仅能获得优于1km的定位
精度。
[0004]
X射线脉冲星是一种高速旋转的
中子星，其远离地球并能产生X波段的电磁
辐射。X射线脉冲星的自转周期长期稳定，一些毫秒脉冲星的周期稳定度可媲美当前的
原子钟。欧美等国从上世纪70年代起开始从事利用X射线脉冲星实现地球卫星自主导航的方法。由于X射线脉冲星的信号全天区可见，X射线脉冲星导航系统的性能不会随着卫星轨道高度的变化而变化。因此，利用X射线脉冲星导航系统可自主提供可靠、稳定的高轨卫星的位置信息。然而，脉冲星信号微弱，且不连续。信号需要进行持续若干分钟的累积，才可获得可信的测量信息。在信号积累时间内的导航误差无法忽略，而这种误差也会大幅降低X射线
2
脉冲星导航系统的性能。当前，大多利用面积为1m 的X射线探测器来提高脉冲信号的
信噪比。然而，大面积探测器带来的体积和功耗负担，是一般卫星难以承受的。
[0005] 鉴于X射线脉冲星导航系统和CNS可互补，
申请号为200710191527.9的中国
专利文件公布了一种基于多种信息融合的组合导航方法（以下简记为方法1）。方法1融合了四种测量信息，其中就包含了星光角距和脉冲星测量信息。但方法1未使用小面积的探测器，且未使用非线性动静态
滤波器。申请号为200910063267.4的中国专利文件公布了一种基于脉冲星/CNS的组合导航方法（以下简记为方法2）。该方法基于相邻历元的观测数据做差方法，采用无迹卡尔曼滤波（Unscented Kalman filter，UKF）/H∞滤波器来削弱脉冲星星表位置误差和钟差的影响。方法2利用星光仰角进行组合导航，并利用联邦滤波融合X射线脉冲星导航系统和CNS的导航结果。
发明内容
[0006] 本发明要解决的技术问题为：针对当前的高轨卫星难以实现自主定位的问题，采用基于动静态非
线性滤波器的脉冲星/星光角距组合导航方法来完成高轨卫星的自主导航定位工作。
[0007] 为解决上述技术问题，本发明采用以下技术方案：
[0008] 一种高轨卫星脉冲星/星光角距组合导航方法，步骤为：
[0010] B、建立脉冲星导航观测方程；
[0011] C、建立星光角距观测方程；
[0012] D、利用动态滤波器处理动力学模型和星光角距信息。具体步骤为：
[0013] ①初始化状态量 和方差阵P0
[0014]
[0015] ②构建状态采样点和权重
[0016] 在 附近建立一系列的采样点，这些采样点的均值和协方差分别为 和Pk-1。状态变量为6×1维，所以13个采样点χ0，k-1，χ1，k-1，…，χ13，k-1及均值权重和方差权重 如下：
[0017]
[0018]
[0019] 在式(9)(10)中，n为状态量的维数,ξ=α2(n+κ)-n，其中α用来控制采样点分布，其取值在0,1之间，κ=3-n。β为与状态先验分布有关的参数，对高斯分布，β取为2。T T
当Pk-1＝AA时， 取A的第i行；当Pk-1＝AA 时， 取A的第i列。
[0020] ③动态滤波器时间更新
[0021] 采样点的一步预测值χk|k-1为
[0022] χk|k-1＝f(χk-1) (11)
[0023] 所有采样点的一步预测值加权后的结果 为
[0024]
[0025] 状态量估计的一步估计方差阵 为
[0026]
[0027] ④动态滤波器测量更新
[0028] 采样点一步预测的星光角距测量值Zk|k-1st为
[0029] Zk|k-1st＝hst(χk|k-1) (14)
[0030] 所有采样点的一步预测星光角距测量值加权后的结果 为
[0031]
[0032] 星光角距测量方差阵 为
[0033]
[0034] 星光角距测量值与状态量的协方差阵 为
[0035]
[0036] 滤波增益Kst,k为
[0037]
[0038] 动态滤波器的状态估计值xst,k+和方差阵Pst,k+分别为
[0039]
[0040]
[0041] E、利用静态滤波器融合动态滤波器结果和脉冲星观测信息。具体步骤为[0042] ①静态滤波器初始化
[0043] 静态滤波器的初始状态预估值xp,k-和方差阵Pp,k-为
[0044] xp,k-＝xst,k+Pp,k-＝Pst,k+ (21)
[0045] ②构建静态滤波器采样点
[0046] 静态滤波器的采样点ε0,k-1，ε1,k-1，…,ε13,k-1为
[0047]
[0048] 其中，当Pp,k-＝ATA时， 取A的第i行；当Pp,k-＝AAT时， 取A的第i列。ξ的定义与式(10)相同。
[0049] ③静态滤波器测量更新
[0050] 采样点一步预测的脉冲星测量值Zk|k-1p为
[0051] Zk|k-1p＝hp(εk|k-1) (23)
[0052] 所有采样点的一步预测脉冲星测量值加权后的结果 为
[0053]
[0054] 脉冲星测量方差阵 为
[0055]
[0056] 星光角距测量值与状态量的协方差阵 为
[0057]
[0058] 滤波增益Kp,k为
[0059]
[0060] 静态滤波器的状态估计量 和方差阵 分别为
[0061]
[0062]
[0063] 将获得的状态量 和方差阵 返回动态滤波器，用于k+1时刻，k=1,2，…。最终得到各个导航时刻的状态估计值 和方差阵 每个导航时刻的状态估计值即为高轨卫星自主获得的其在该时刻的位置信息，从而实现了高轨卫星的自主导航功能。
[0064] 优选地，本发明步骤A包括：
[0065] 在地心惯性系J2000.0中，建立近地卫星的轨道动力学方程为
[0066]
[0067] 其中，x=[rT,vT]T是
航天器的状态矢量，w=[wrT wvT]T为动力学模型噪声，可建模为零均值高斯白噪声。该噪声的方差阵为Q。a＝aTB+aNS+aT+aH.O.T是航天器受到的
加速度，包含以下几项。
[0068] 1）aTB＝-μEr/|r|3是航天器受到的地球二体引力加速度，其中μE是地球的引力常数。
[0069] 2） 是地球的非球形引力摄动。UNSE可以表示为
[0070]
[0071] 其中，Re是地球的半径，φ和λ分别是经度和纬度，Pn和Pnm是勒让德多项式，Jn是带谐项系数，λn,m是田谐项，Jn,m为田谐项系数。
[0072] 3） 是三体摄动加速度。μi是第i个天体的引力常数，ri是第i个天体相对于地球的位置矢量。
[0073] 4）aH.O.T是影响航天器加速度的高阶项。相对于已建模的摄动加速度，这些高阶项的影响可以忽略。
[0074] 式(1)可简记为
[0075] x＝f(x)+w (3)
[0076] 优选地，本发明步骤B包括：
[0077] 脉冲星导航的观测方程为
[0078] Zp＝hP(x)+Vp (4)1 j N T j
[0079] 其中，Vp为观测噪声，hP(x)＝[hP(x) L hP(x) L hP(x)],hP(x)为观测第j个脉冲星的观测方程，其表达式为
[0080]j
[0081] 其中， 为第j个脉冲星方向矢量的测量值， 为由星历预报的地球位置，D0 为第j个脉冲星相距太阳质心的距离，b为太阳质心相对于
太阳系质心的位置，μS为太阳的引力常数，c为光速。
[0082] 在式(4)中，Vp可看做零均值高斯白噪声，其标准方差由下式决定。
[0083]
[0084] 其中，W为脉冲宽度，BX为
宇宙背景辐射流量，FX为脉冲星流量，pf为脉冲星的脉冲比，d为脉冲宽度与脉冲周期的比，A为X射线探测器的面积，tm为脉冲星测量周期。
[0085] 优选地，本发明步骤C包括：
[0086] 星光角距观测方程为
[0087]
[0088] 其中，s为参考恒星的方向矢量。vα可看做零均值高斯白噪声，其标准方差由星敏感器和光学相机的精度决定。
[0090] （1）采用脉冲星/星光角距组合导航方法，能够在发挥脉冲星导航系统和CNS优势的同时，弥补两种导航系统单一运行的不足。由式(6)可知，脉冲星导航的测量精度与航天器轨道无关，仅取决于脉冲星的物理特性、星载探测器面积和脉冲星观测周期。因此，脉冲星导航对高轨卫星仍可获得理想的导航效果，但相应的，脉冲星的观测周期需延长，且需要大面积的X射线探测器。由式(7)可知，CNS的精度仅取决于导航敏感器的精度，与观测时间无关。因此，CNS的观测周期远小于脉冲星导航系统。此外，根据当前的工艺
水平，星敏感器和光学相机的体积功耗均远小于X射线探测器(王鹏.基于星敏感器的卫星自主导航及
姿态确定方法研究[D].哈尔滨工业大学,2008.)。但由于CNS测量的角度信息与航天器的轨道高度有关，CNS对高轨卫星的导航精度偏低。本专利提出的组合导航系统融合利用了脉冲星导航系统、CNS的测量信息。相比于脉冲星导航系统，该导航系统可在一个脉冲星观测周期内采用CNS测量信息进行导航，从而降低了整个导航系统的测量周期。再者，由于融合利用了CNS的测量信息，在导航精度允许的前提下，X射线探测器的面积也可适当降低，从而可以降低整个导航系统的体积和功耗。另一方面，该组合导航系统利用相对高精度的脉冲星测量信息来融合CNS的低精度测量信息，可获得优于CNS的导航精度。
[0091] （2）采用动静态非线性滤波器可实现脉冲星原始观测信息和星光角距信息的最优融合。现有的联邦滤波器等传统信息融合
算法首先利用动力学模型信息分别同脉冲星原始观测信息，星光角距信息进行导航计算。然后利用两组导航结果进行信息融合。可以看出，现有该种做法会导致动力学模型的信息被重复使用，从而导致了最终融合的两组导航结果之间具有相关性，信息融合结果次优。然而，从式（9）-（21）可以看出，在动静态非线性滤波器中，动力学模型的信息仅在动态滤波器中使用，静态滤波器是基于动态滤波器的结果进行工作。所以，动静态非线性滤波器中不存在信息重复利用的情况，可实现脉冲星原始观测信息和星光角距信息的最优融合。
[0092] （3）在动静态非线性滤波器中嵌入无迹变换（unscented transformation）降低线性化误差的影响。传统的滤波器是针对线性系统推导的（杨元喜.多源
传感器动、静态滤波融合导航[J].武汉大学学报(信息科学版),2003,28(4):386-388.），虽然作者利用
扩展卡尔曼滤波将动静态滤波器应用于非线性情况，但扩展卡尔曼滤波存在的线性化误差会影响动静态滤波器的性能（Chang,L.B.,Hu,B.Q.,Chang,G.B,et al.,Multiple Outliers Suppression Derivative-Free Filter Based on Unscented Transformation[J],Journal of Guidance,Control,and Dynamics,2012,35(6),1902-1906.）。本专利提出的动静态滤波器在式(9)和式(22)中采用了无迹变换的思想，不存在线性化过程，从而降低了线性化误差的影响。
[0093] 由于高轨卫星在日常生活中扮演着十分重要的角色，本发明实现的高轨卫星的自主导航定位功能，至少具有以下方面的优势：1、可大幅降低地面测控系统的人力投入和资金投入；2、若地面测控系统失效，高轨卫星仍能通过自主定位获得自身的位置信息，从而保证了高轨卫星功能实现。
附图说明
[0094] 图1是本发明具体应用实例的流程示意图。
[0095] 图2是本发明具体应用实例的效果对比图。
具体实施方式
[0096] 如图1所示，一种基于脉冲星和星光角距的中高轨自主导航算法，具体的实施步骤如下：
[0097] （1）建立近地卫星的轨道动力学方程；
[0098] 在地心惯性系J2000.0中，建立近地卫星的轨道动力学方程为
[0099]
[0100] 其中，x=[rT,vT]T是航天器的状态矢量，w=[wrT wvT]T为动力学模型噪声，可建模为零均值高斯白噪声。该噪声的方差阵为Q。a＝aTB+aNS+aT+aH.O.T是航天器受到的加速度，包含以下几项。
[0101] 1）aTB＝-μEr/|r|3是航天器受到的地球二体引力加速度，其中μE是地球的引力常数。
[0102] 2） 是地球的非球形引力摄动。UNSE可以表示为
[0103]
[0104] 其中，Re是地球的半径，φ和λ分别是经度和纬度，Pn和Pnm是勒让德多项式，Jn是带谐项系数，λn,m是田谐项，Jn,m为田谐项系数。
[0105] 3） 是三体摄动加速度。μi是第i个天体的引力常数，ri是第i个天体相对于地球的位置矢量。
[0106] 4）aH.O.T是影响航天器加速度的高阶项。相对于已建模的摄动加速度，这些高阶项的影响可以忽略。
[0107] 式(1)可简记为
[0108] x＝f(x)+w (3)
[0109] （2）建立脉冲星导航观测方程；
[0110] 脉冲星导航的观测方程为
[0111] Zp＝hP(x)+Vp (4)
[0112] 其中，Vp为观测噪声，hP(x)＝[hP1(x) L hPj(x) L hPN(x)]T,hPj(x)为观测第j个脉冲星的观测方程，其表达式为
[0113]
[0114] 其中， 为第j个脉冲星方向矢量的测量值， 为由星历预报的地球位置，D0j为第j个脉冲星相距太阳质心的距离，b为太阳质心相对于太阳系质心的位置，μS为太阳的引力常数，c为光速。
[0115] 在式(4)中，Vp可看做零均值高斯白噪声，其标准方差由下式决定。
[0116]
[0117] 其中，W为脉冲宽度，BX为宇宙背景辐射流量，FX为脉冲星流量，pf为脉冲星的脉冲比，d为脉冲宽度与脉冲周期的比，A为X射线探测器的面积，tm为脉冲星测量周期。
[0118] （3）建立星光角距观测方程；
[0119] 星光角距观测方程为
[0120]
[0121] 其中，s为参考恒星的方向矢量。vα可看做零均值高斯白噪声，其标准方差由星敏感器和光学相机的精度决定。
[0122] （4）利用动态滤波器处理动力学模型和星光角距信息。
[0123] 具体步骤为
[0124] ①初始化状态量 和方差阵P0
[0125]
[0126] ②构建状态采样点和权重
[0127] 在 附近建立一系列的采样点，这些采样点的均值和协方差分别为 和Pk-1。状态变量为6×1维，所以13个采样点χ0,k-1，χ1,k-1，…，χ13,k-1及均值权重和方差权重 如下：
[0128]
[0129]
[0130] 在式(9)(10)中，n为状态量的维数,ξ=α2(n+κ)-n，其中α用来控制采样点分布，其取值在0,1之间，κ=3-n。β为与状态先验分布有关的参数，对高斯分布，β取为2。T T
当Pk-1＝AA时， 取A的第i行；当Pk-1＝AA 时，
[0131] ③动态滤波器时间更新
[0132] 采样点的一步预测值χk|k-1为
[0133] χk|k-1＝f(χk-1) (11)
[0134] 所有采样点的一步预测值加权后的结果 为
[0135]
[0136] 状态量估计的一步估计方差阵 为
[0137]
[0138] ④动态滤波器测量更新
[0139] 采样点一步预测的星光角距测量值Zk|k-1st为
[0140] Zk|k-1st＝hst(χk|k-1) (14)所有采样点的一步预测星光角距测量值加权后的结果为
[0141]
[0142] 星光角距测量方差阵 为
[0143]
[0144] 星光角距测量值与状态量的协方差阵 为
[0145]
[0146] 滤波增益Kst,k为
[0147]
[0148] 动态滤波器的状态估计值xst,k+和方差阵Pst,k+分别为
[0149]
[0150]
[0151] （5）利用静态滤波器融合动态滤波器结果和脉冲星观测信息。具体步骤为[0152] ①静态滤波器初始化- -
[0153] 静态滤波器的初始状态预估值xp,k 和方差阵Pp,k 为- + - +
[0154] xp,k ＝xst,kPp,k ＝Pst,k (21)
[0155] ②构建静态滤波器采样点
[0156] 静态滤波器的采样点ε0,k-1，ε1,k-1，…,ε13,k-1为
[0157]
[0158] 其中，当Pp,k-＝ATA时， 取A的第i行；当Pp,k-＝AAT时， 取A的第i列。ξ的定义与式(10)相同。
[0159] ③静态滤波器测量更新
[0160] 采样点一步预测的脉冲星测量值Zk|k-1p为
[0161] Zk|k-1p＝hp(εk|k-1) (23)
[0162] 所有采样点的一步预测脉冲星测量值加权后的结果 为
[0163]
[0164] 脉冲星测量方差阵 为
[0165]
[0166] 星光角距测量值与状态量的协方差阵 为
[0167]
[0168] 滤波增益Kp,k为
[0169]
[0170] 静态滤波器的状态估计量 和方差阵 分别为
[0171]
[0172]
[0173] 将获得的状态量 和方差阵 返回动态滤波器，用于k+1时刻，k=1,2，…。最终得到各个导航时刻的状态估计值 和方差阵
[0174] 以某一地球高轨卫星为例，初始轨道根数为：半长轴26570.08km，偏心率0.0047，轨道倾角63.15°，升交点赤经136.31°，近交点角距211.64°，平近点角304.82°。导航时间为7天。图2给出了本专利所提出的组合导航系统、星光角距导航系统、脉冲星导航系统的性能对比情况。从图2可以看出，随着时间的增加，组合导航系统、脉冲星导航系统、星光角距导航系统（CNS）的位置估计误差曲线均逐渐收敛到0附近。但可以看出，组合导航系统的位置估计误差曲线收敛速度最快，且可获得的最终位置估计误差最小，即导航精度最高。相对的，脉冲星导航的误差曲线收敛速度优于星光角距导航系统，且脉冲星导航系统的最终导航精度也优于星光角距导航系统。因此，由图2可知，组合导航系统的性能优于单纯的星光角距导航系统和脉冲星导航系统。
[0175] 以上所述仅是本发明的优选实施方式，本发明的保护范围并不仅局限于上述
实施例，凡属于本发明思路下的技术方案均属于本发明的保护范围。应该提出，对于本技术领域的普通技术人员来说，在不脱离本发明原理前提下的改进和润饰，这些改进和润饰也应视为本发明的保护范围。
