CN106932804A - 天文辅助的惯性/北斗紧组合导航系统及其导航方法 - Google Patents

Abstract

本发明公开一种天文辅助的惯性/北斗紧组合导航系统，惯性元件、北斗接收机和星敏感器将数据送入导航计算机模块，导航计算机模块进行处理后将导航解算结果送入上位机显示模块实时显示。本发明还公开一种天文辅助的惯性/北斗紧组合导航方法，首先计算速度、位置，用四元数描述载体的姿态运动；建立天文辅助的导航系统状态方程；根据各导航系统的不同工作特性，建立各子系统的量测方程；根据可见星及可见星的几何构型构造自适应卡尔曼滤波器，从而对惯性导航系统的导航误差进行修正。此种技术方案可以高效实时实现惯性/天文/北斗组合导航系统的子系统解算及信息融合处理，获取高精度、高可靠性的导航结果。

Description

天文辅助的惯性/北斗紧组合导航系统及其导航方法

技术领域

本发明属于飞行器组合导航技术领域，特别涉及种天文辅助的惯性/北斗紧组合导航系统及其导航方法。

背景技术

目前美国的GPS是使用最为广泛的卫星导航系统，卫星导航系统在国家的国防及经济上起着及其重要的作用。GPS系统的主导权在美国手里，这对我国的国防安全造成了一定的威胁，因此国内北斗导航系统的发展是必然趋势。由于我国北斗导航系统仍处于发展阶段，可见星少，且卫星信号容易受地形、建筑物等影响，也易受到电磁干扰而失去作用，因此惯性/BD紧组合导航算法的研究有着积极意义。我国近年空天技术发展迅速，空天飞行器具有跨空天、大机动、高超声速、长航时的飞行特性，惯性/北斗两组合导航系统无法满足其对高精度、高可靠性导航系统的需求。天文定姿、定位具有精度高且不受电磁干扰的优点，因此天文辅助的惯性/北斗紧组合导航系统是实现远程长航时、高机动、高超声速飞行器圆满完成任务的保证。

考虑到北斗卫星导航系统的可见星数目及卫星信息质量随着时间动态变化，因此,需要分析在天文辅助的惯性/北斗导航系统中北斗导航系统量测信息的可信赖性问题。在惯性、天文、北斗组合导航系统中的信息融合过程中主要是采取的卡尔曼滤波方法，系统参数的设置对整个系统有着极大的影响，参数设置不准确不仅影响系统的精确性，还易导致整个系统发散直至崩溃。鉴于北斗导航系统量测信息质量动态变化的情况，天文辅助的惯性/北斗紧组合导航系统的自适应滤波方法有待进一步研究。

发明内容

本发明的目的，在于提供一种天文辅助的惯性/北斗紧组合导航系统及其导航方法，其可解决高机动、高超声速、长航时飞行器对导航系统高精度、高可靠性的需求，可以高效实时实现惯性/天文/北斗组合导航系统的子系统解算及信息融合处理，获取高精度、高可靠性的导航结果，同时为天文辅助的惯性/北斗紧组合导航系统的工程化实现和应用提供了有效的支撑，具有突出的实用价值。

为了达成上述目的，本发明的解决方案是：

一种天文辅助的惯性/北斗紧组合导航系统，包括导航计算机模块、惯性元件、北斗接收机、星敏感器和上位机显示模块，惯性元件、北斗接收机和星敏感器将数据送入导航计算机模块，导航计算机模块进行处理后将导航解算结果送入上位机显示模块实时显示。

上述导航计算机模块包括惯导解算模块、量测转换模块和卡尔曼滤波模块，惯导解算模块和量测转换模块的输出端均连接卡尔曼滤波模块的输入端。

上述卡尔曼滤波模块包括SINS/CNS滤波器和SINS/BD自适应滤波器，SINS/CNS滤波器的输入端分别连接惯导解算模块和量测转换模块的输出端，SINS/BD自适应滤波器的输入端分别连接惯导解算模块和量测转换模块的输出端。

上述导航计算机模块采用PC机、工作站、基于WINDOWS系统或LINUX系统的嵌入式平台。

一种天文辅助的惯性/北斗紧组合导航方法，包括如下步骤：

步骤1，对于天文辅助的惯性/北斗紧组合导航系统，速度计算公式为：

式中，为载体相当于地理系的加速度在载体上的投影，为载体相当于地理系的速度在载体上的投影，为机体系到导航坐标系的姿态转移矩阵，为载体相当于地理坐标系的速度在导航坐标系上的投影，其中的地理坐标系即为导航坐标系；

位置计算公式为：

用四元数描述载体的姿态运动，计算公式为：

用矩阵表示为

式中，q为姿态四元数，为q对时间的导数，矢量为四元数形式，表示载体坐标系相对导航坐标系在载体坐标系上的投影；

步骤2，在上述基础上，建立天文辅助的SINS/BD紧组合导航系统状态方程为：

其中，状态变量为：

表示捷联惯性导航在地里坐标系下东、北、天方向平台误差角，δv_E,δv_N,δv_U表示惯性导航系统在地理坐标系下东、北、天方向的速度误差；δL,δλ,δh表示惯性导航系统在地里坐标系下纬度、经度、高度的误差；ε_bx,ε_by,ε_bz表示陀螺随机常数误差；ε_rx,ε_ry,ε_rz表示陀螺一阶马尔可夫过程随机误差；表示加速度计的一阶马尔可夫过程随机误差；δt_u,δt_ru表示北斗接收机和卫星之间的钟差和频差，分别以距离和速度体现；A(t)_20×20为系统的状态转移矩阵；G(t)_20×11为噪声系数矩阵；W(t)_11×1为系统的白噪声矢量；

步骤3，根据各导航系统的不同工作特性，建立各子系统的量测方程；

步骤4，根据可见星及可见星的几何构型构造自适应卡尔曼滤波器，从而对惯性导航系统的导航误差进行修正。

上述步骤3中，SINS/BD子系统的量测方程分为伪距量测部分和伪距率量测部分，其中伪距量测为利用惯导提供的位置信息计算出的伪距信息与北斗接收机提供的伪距信息之差；伪距率量测信息是通过惯导提供的速度信息、卫星多普勒频移信息以及卫星的速度信息计算出来的；

设接收机某时刻观测到n颗可用卫星，则伪距量测方程：

Z_ρ(t)＝H_ρ(t)X(t)+v_ρ(t)

其中H_ρ(t)＝[0_n×6 H_ρ1 0_n×6 H_ρ2]，e_i1、e_i2、e_i3为第i颗卫星与载体之间的方向余弦，v_ρ(t)为伪距量测噪声，R_N为地球参考椭球卯酉圈曲率半径；

伪距率量测方程为：

其中

为伪距率量测噪声。

上述步骤4的详细内容是：

(41)当可见星小于4时，设定量测噪声随着可见星的增多而逐渐减小，即可见星越多则越信赖北斗导航系统的量测信息，伪距和伪距率的量测噪声表示如下：

式中n代表可见星数目，v_ρ0，分别表示伪距、伪距率的初始量测噪声，α表示一个常系数；

(42)当北斗可见星大于4时，首先根据导航位置和卫星位置计算卫星的高度角和方位角，通过卫星的高度角和方位角计算GDOP值，最终对伪距和伪距率的量测噪声进行调节，过程如下：

a、卫星高度角、方位角计算

用户导航位置在地固坐标系下的位置为：

P_u＝[p_ux p_uy p_uz]

第i颗卫星在地固标系下的位置为：

P_s＝[p_sx p_sy p_sz]

设

D＝[d₁ d₂ d₃]＝(P_s-P_u)E

则高度角为：

方位角为：

b、几何精度因子GDOP的计算

设某时刻北斗卫星导航系统的可见星数目为n，设矩阵Q如下：

则GDOP的值由下式计算：

c、卫星到接收机伪距、伪距率量测噪声的设置

其中v_ρ0，分别表示伪距、伪距率的初始量测噪声，GDOP₀表示北斗导航系统的初始几何精度因子，GDOP_k表示北斗导航系统第k时刻的几何精度因子。

上述步骤3中，SINS/CNS子系统量测为惯导给出的姿态角与天文给出的姿态角信息之差，定义观测量为星敏感器测得的载体姿态角与捷联惯导的姿态角之差，量测方程如下：

Z_s(t)＝H_s(t)X(t)+v_s(t)

由于组合导航系统状态方程中的误差角状态量为惯性导航平台误差角，它描述的是数学平台与地里面坐标系之间的误差角关系，因此需要将惯性姿态误差角转化到平台误差角；在平台误差角为小量时，量测矩阵由以下公式表示：

因此，量测方程表示为：

式中，γ、θ、ψ分别为横滚角、俯仰角和航向角，v_s(t)为姿态误差角量测噪声。

上述步骤4中，滤波过程如下：

上式中，表示t_k-1时刻的状态对t_k时刻的状态的最优估计值，Φ_k/k-1表示t_k-1时刻至t_k时刻系统的状态转移矩阵，表示t_k-1时刻的系统状态估计值，K_k表示增益矩阵，H_k表示滤波模块t_k时刻的观测系数矩阵，P_k/k-1表示最优预测估计误差协方差阵，P_k/k表示滤波模块t_k时刻的系统误差协方差阵，Q_k-1表示t_k-1时刻的噪声方差矩阵，Γ_k/k-1表示t_k-1时刻的噪声矢量对t_k时刻状态矢量影响的噪声系数矩阵，R_k表示t_k时刻的量测方差矩阵，I为单位矩阵。

采用上述方案后，本发明的有益效果如下：

(1)在整个导航系统中，惯性与北斗导航系统采用的是紧组合方式，提高系统在可见星少时的适用性，天文信息的辅助明显提高了导航系统姿态角的精度。针对惯性/北斗紧组合模块，提出了基于卫星高度角和方位角的自适应滤波方法，大大增强了系统的自适应性。

(2)硬件层面上，整个导航系统移植性很强，可以在工作站、PC机、嵌入式等模块上工作。系统中的传感器逻辑上属于并列关系，传感器与导航计算机、导航计算机与上位机之间可以通过USB串口或RS232串口进行通信，提高了系统的灵活性，便于系统的扩展与裁剪。

(3)软件层面上，整个导航系统采用多线程模式，整个导航系统在不同任务间进行高速运行与轮转，任务并发进行保证了系统的实时性，利于系统在任务上的扩展与裁剪。

附图说明

图1是本发明的原理图；

图2是本发明导航系统硬件模块化构成图；

图3是有无天文辅助的姿态误差曲线对比图；

图4和图5分别是AKF与KF的位置误差曲线对比图、AKF与KF的速度误差曲线对比图；

在本次试验中，通过人为设定使卫星可见星每隔25秒实时变换一次，图4和图5的非自适应滤波的波形即能反映可见星变化的情况。

具体实施方式

以下将结合附图，对本发明的技术方案进行详细说明。

本发明提供一种天文辅助的惯性/北斗紧组合导航系统及其导航方法，采用可移植性很强的跨平台导航程序和导航硬件装置，导航程序由Qt Creator编写，考虑到导航系统的实时性及模块化配置，整个导航程序采取多线程方案，分别包括捷联惯导线程、SINS/BD紧组合滤波线程、SINS/天文滤波线程、串口数据读取解码及导航周期控制线程、上位机显示及导航模式配置线程；导航硬件装置由星敏感器仿真器、北斗接收机、MEMS惯性元件、导航计算机模块、上位机显示模块组成，其中导航计算机模块可以由PC机、工作站、基于WINDOWS系统或LINUX系统的嵌入式平台承担。

如图1所示，本发明一种天文辅助的惯性/北斗紧组合导航方法的原理是：对于惯性/北斗紧组合模块，根据可见星的高度角和方位角获得精度因子GDOP，再根据精度因子实时调整量测噪声，提高滤波模块的自适应性。除此之外，利用星光信息对整个导航系统的姿态角进行修正。具体实施方法如下：

1、天文辅助的惯性/北斗紧组合导航系统建模

1)组合导航系统惯导误差状态量方程建模

对于天文辅助的惯性/北斗紧组合导航系统，其捷联惯导解算过程如下，速度计算公式为：

式中，为载体相当于地理系的加速度在载体上的投影，为载体相当于地理系的速度在载体上的投影，为机体系到导航坐标系的姿态转移矩阵，为载体相当于地理坐标系的速度在导航坐标系上的投影，本模型的地理坐标系即为导航坐标系。

位置计算公式为：

用四元数描述载体的姿态运动，计算公式为：

用矩阵表示为

式中，q为姿态四元数，为q对时间的导数，矢量为四元数形式，表示载体坐标系相对导航坐标系在载体坐标系上的投影。

在上述基础上，建立天文辅助的SINS/BD紧组合导航系统状态方程为：

其中，状态变量为：

表示捷联惯性导航在地里坐标系下东、北、天方向平台误差角，δv_E,δv_N,δv_U表示惯性导航系统在地理坐标系下东、北、天方向的速度误差；δL,δλ,δh表示惯性导航系统在地里坐标系下纬度、经度、高度的误差；ε_bx,ε_by,ε_bz表示陀螺随机常数误差；ε_rx,ε_ry,ε_rz表示陀螺一阶马尔可夫过程随机误差；表示加速度计的一阶马尔可夫过程随机误差；δt_u,δt_ru表示北斗接收机和卫星之间的钟差和频差，分别以距离和速度体现。A(t)_20×20为系统的状态转移矩阵；G(t)_20×11为噪声系数矩阵；W(t)_11×1为系统的白噪声矢量；

2)组合导航系统量测方程建模

在此所叙述的天文辅助惯性/北斗紧组合导航系统，是利用天文的姿态信息对惯导的姿态信息进行滤波估计，利用北斗导航系统提供的伪距和多普勒信息对惯导的位置和速度信息进行滤波估计。根据各导航系统的不同工作特性，建立各子系统的量测方程。SINS/CNS子系统量测为惯导给出的姿态角与天文给出的姿态角信息之差。SINS/BD子系统的量测方程分为伪距量测部分和伪距率量测部分，其中伪距量测为利用惯导提供的位置信息计算出的伪距信息与北斗接收机提供的伪距信息之差；伪距率量测信息是通过惯导提供的速度信息、卫星多普勒频移信息以及卫星的速度信息计算出来的。

a、SINS/BD紧组合系统

设接收机某时刻观测到n颗可用卫星，则伪距量测方程：

Z_ρ(t)＝H_ρ(t)X(t)+v_ρ(t) (10)

其中H_ρ(t)＝[0_n×6 H_ρ1 0_n×6 H_ρ2]，e_i1、e_i2、e_i3为第i颗卫星与载体之间的方向余弦，v_ρ(t)为伪距量测噪声，R_N为地球参考椭球卯酉圈曲率半径。

伪距率量测方程为：

其中

为伪距率量测噪声。

b、SINS/CNS组合系统

在SINS/CNS子系统中，定义观测量为星敏感器测得的载体姿态角与捷联惯导的姿态角之差，量测方程如下：

Z_s(t)＝H_s(t)X(t)+v_s(t) (12)

由于组合导航系统状态方程中的误差角状态量为惯性导航平台误差角，它描述的是数学平台与地里面坐标系之间的误差角关系，因此需要将惯性姿态误差角转化到平台误差角。在平台误差角为小量时，量测矩阵可以由以下公式表示：

因此，量测方程可以表示为：

式中，γ、θ、ψ分别为横滚角、俯仰角和航向角，v_s(t)为姿态误差角量测噪声。

2、惯性/北斗紧组合导航系统自适应滤波方法

对于SINS/BD子系统，由于北斗导航系统的可见星数目随着时间变化，可见星数目及可见星的构型可导致伪距和伪距率信息质量的动态变化，因此根据可见星及可见星的几何构型构造自适应卡尔曼滤波器，从而对惯性导航系统的导航误差进行修正。

当北斗可见星大于4颗时，可以通过可见星的高度角和方位角计算出几何精度因子(GDOP)值，GDOP是衡量定位精度的很重要的一个系数，它代表北斗测距误差造成的接收机与空间卫星间的距离矢量放大因子。实际表征参与定位解的从接收机至空间卫星的单位矢量所勾勒的形体体积与GDOP成反比，GDOP越小则定位精度越高。在惯性/北斗紧组合导航系统中，可以根据导航位置和卫星位置计算各个卫星的高度角和方位角，进而根据高度角和方位角计算描述卫星定位精度的几何定位精度因子GDOP，最后通过GDOP的大小实时调整伪距、伪距率的量测噪声，从而使整个组合系统达到自适应调整的性能。其转换过程如下：

(1)当可见星小于4时，设定量测噪声随着可见星的增多而逐渐减小，即可见星越多则越信赖北斗导航系统的量测信息。伪距和伪距率的量测噪声表示如下：

式中n代表可见星数目(小于4)，v_ρ0，分别表示伪距、伪距率的初始量测噪声，α表示一个常系数，一般要根据系统的实测情况进行设定。

(2)当北斗可见星大于4时，首先根据导航位置和卫星位置计算卫星的高度角和方位角，通过卫星的高度角和方位角计算GDOP值，最终对伪距和伪距率的量测噪声进行调节，过程如下：

a、卫星高度角、方位角计算

用户导航位置在地固坐标系下的位置为：

P_u＝[p_ux p_uy p_uz] (17)

第i颗卫星在地固标系下的位置为：

P_s＝[p_sx p_sy p_sz] (18)

设

D＝[d₁ d₂ d₃]＝(P_s-P_u)E (20)

则高度角为：

方位角为：

b、几何精度因子GDOP的计算

设某时刻北斗卫星导航系统的可见星数目为n，设矩阵Q如下：

则GDOP的值可由下式计算：

c、卫星到接收机伪距、伪距率量测噪声的设置

其中v_ρ0，分别表示伪距、伪距率的初始量测噪声，GDOP₀表示北斗导航系统的初始几何精度因子，GDOP_k表示北斗导航系统第k时刻的几何精度因子。

(3)天文辅助的惯性/北斗紧组合导航系统中，对于CNS/SINS和SINS/BD子系统，各滤波模块过程如下：

上式中，表示t_k-1时刻的状态对t_k时刻的状态的最优估计值，Φ_k/k-1表示t_k-1时刻至t_k时刻系统的状态转移矩阵，表示t_k-1时刻的系统状态估计值，K_k表示增益矩阵，H_k表示滤波模块t_k时刻的观测系数矩阵，P_k/k-1表示最优预测估计误差协方差阵，P_k/k表示滤波模块t_k时刻的系统误差协方差阵，Q_k-1表示t_k-1时刻的噪声方差矩阵，Γ_k/k-1表示t_k-1时刻的噪声矢量对t_k时刻状态矢量影响的噪声系数矩阵，R_k表示t_k时刻的量测方差矩阵，I为单位矩阵。

综合上述，本发明公布了一种天文辅助的惯性/北斗紧组合导航系统及其导航方法，整个导航程序由Qt编写，采用多线程并发模式，以保证整个导航程序的实时性与可扩展性。导航程序由捷联惯导线程、SINS/BD紧组合滤波线程、SINS/天文滤波线程、串口数据读取解码及导航周期控制线程、上位机显示及导航模式配置线程组成。导航硬件装置由星敏感器仿真器、北斗接收机、MEMS惯性元件、导航计算机、上位机显示模块组成，其中导航计算机可以由PC机、工作站、基于WINDOWS系统或LINUX系统的嵌入式平台承担。导航方法包括以下步骤：首先在地理坐标系下分别建立惯性/北斗紧组合模型、惯性/天文组合模型，在惯性/北斗紧组合的基础上研究基于卫星高度角和方位角的自适应滤波方法。北斗卫星可见星数目的动态变化及卫星几何构型的变化，不仅影响整个系统的精确性，在极端情况下还会造成系统发散甚至崩溃。本发明可以有效解决这个问题，并获取高精度的导航结果，同时为天文辅助的惯性/北斗紧组合导航系统的工程实现和应用提供了有效的支撑。

以上实施例仅为说明本发明的技术思想，不能以此限定本发明的保护范围，凡是按照本发明提出的技术思想，在技术方案基础上所做的任何改动，均落入本发明保护范围之内。

Claims (9)

1.一种天文辅助的惯性/北斗紧组合导航系统，其特征在于：包括导航计算机模块、惯性元件、北斗接收机、星敏感器和上位机显示模块，惯性元件、北斗接收机和星敏感器将数据送入导航计算机模块，导航计算机模块进行处理后将导航解算结果送入上位机显示模块实时显示。

2.如权利要求1所述的天文辅助的惯性/北斗紧组合导航系统，其特征在于：所述导航计算机模块包括惯导解算模块、量测转换模块和卡尔曼滤波模块，惯导解算模块和量测转换模块的输出端均连接卡尔曼滤波模块的输入端。

3.如权利要求2所述的天文辅助的惯性/北斗紧组合导航系统，其特征在于：所述卡尔曼滤波模块包括SINS/CNS滤波器和SINS/BD自适应滤波器，SINS/CNS滤波器的输入端分别连接惯导解算模块和量测转换模块的输出端，SINS/BD自适应滤波器的输入端分别连接惯导解算模块和量测转换模块的输出端。

4.如权利要求1所述的天文辅助的惯性/北斗紧组合导航系统，其特征在于：所述导航计算机模块采用PC机、工作站、基于WINDOWS系统或LINUX系统的嵌入式平台。

5.一种天文辅助的惯性/北斗紧组合导航方法，其特征在于包括如下步骤：

步骤1，对于天文辅助的惯性/北斗紧组合导航系统，速度计算公式为：

$v_{en}^b=v_{en0}^b+{\∫}_0^t{\stackrel{\·}{v}}_{en}^{b}dt$

$v_{en}^n=C_b^n{v}_{en}^b$

式中，为载体相当于地理系的加速度在载体上的投影，为载体相当于地理系的速度在载体上的投影，为机体系到导航坐标系的姿态转移矩阵，为载体相当于地理坐标系的速度在导航坐标系上的投影，其中的地理坐标系即为导航坐标系；

位置计算公式为：

$L=L_0+{\∫}_0^t\frac{v_y^n}{R_m}dt$

$\mathrm{\λ}={\mathrm{\λ}}_0+{\∫}_0^t\frac{v_x^n}{R_n}\sec Ldt$

$h=h_0+{\∫}_0^t{v}_z^{n}dt$

用四元数描述载体的姿态运动，计算公式为：

用矩阵表示为

$\left[\begin{array}{c}{\stackrel{\·}{q}}_0\\ {\stackrel{\·}{q}}_1\\ {\stackrel{\·}{q}}_2\\ {\stackrel{\·}{q}}_3\end{array}\right]=0.5\×\left[\begin{array}{cccc}0& -{\mathrm{\ω}}_{nbx}^b& -{\mathrm{\ω}}_{nby}^b& -{\mathrm{\ω}}_{nbz}^b\\ {\mathrm{\ω}}_{nbx}^b& 0& {\mathrm{\ω}}_{nbz}^b& -{\mathrm{\ω}}_{nby}^b\\ {\mathrm{\ω}}_{nby}^b& -{\mathrm{\ω}}_{nbz}^b& 0& {\mathrm{\ω}}_{nbx}^b\\ {\mathrm{\ω}}_{nbz}^b& {\mathrm{\ω}}_{nby}^b& {\mathrm{\ω}}_{nbz}^b& 0\end{array}\right]$

式中，q为姿态四元数，为q对时间的导数，矢量为四元数形式，表示载体坐标系相对导航坐标系在载体坐标系上的投影；

步骤2，在上述基础上，建立天文辅助的导航系统状态方程为：

$\left[\begin{array}{c}{\stackrel{\·}{X}}_I\left(t\right)\\ {\stackrel{\·}{X}}_G\left(t\right)\end{array}\right]=\left[\begin{array}{cc}{F}_I\left(t\right)& 0\\ 0& F_G\left(t\right)\end{array}\right]\left[\begin{array}{c}{X}_I\left(t\right)\\ X_G\left(t\right)\end{array}\right]+\left[\begin{array}{cc}{G}_I\left(t\right)& 0\\ 0& G_G\left(t\right)\end{array}\right]\left[\begin{array}{c}{W}_I\left(t\right)\\ W_G\left(t\right)\end{array}\right]$

$\stackrel{\·}{X}{\left(t\right)}_{20\×1}=A{\left(t\right)}_{20\×20}X\left(t\right)+G{\left(t\right)}_{20\×11}W{\left(t\right)}_{11\×1}$

其中，状态变量为：

表示捷联惯性导航在地里坐标系下东、北、天方向平台误差角，δv_E,δv_N,δv_U表示惯性导航系统在地理坐标系下东、北、天方向的速度误差；δL,δλ,δh表示惯性导航系统在地里坐标系下纬度、经度、高度的误差；ε_bx,ε_by,ε_bz表示陀螺随机常数误差；ε_rx,ε_ry,ε_rz表示陀螺一阶马尔可夫过程随机误差；▽_x,▽_y,▽_z表示加速度计的一阶马尔可夫过程随机误差；δt_u,δt_ru表示北斗接收机和卫星之间的钟差和频差，分别以距离和速度体现；A(t)_20×20为系统的状态转移矩阵；G(t)_20×11为噪声系数矩阵；W(t)_11×1为系统的白噪声矢量；

步骤3，根据各导航系统的不同工作特性，建立各子系统的量测方程；

步骤4，根据可见星及可见星的几何构型构造自适应卡尔曼滤波器，从而对惯性导航系统的导航误差进行修正。

6.如权利要求5所述的天文辅助的惯性/北斗紧组合导航方法，其特征在于：所述步骤3中，SINS/BD子系统的量测方程分为伪距量测部分和伪距率量测部分，其中伪距量测为利用惯导提供的位置信息计算出的伪距信息与北斗接收机提供的伪距信息之差；伪距率量测信息是通过惯导提供的速度信息、卫星多普勒频移信息以及卫星的速度信息计算出来的；

设接收机某时刻观测到n颗可用卫星，则伪距量测方程：

Z_ρ(t)＝H_ρ(t)X(t)+v_ρ(t)

其中H_ρ(t)＝[0_n×6 H_ρ1 0_n×6 H_ρ2]，e_i1、e_i2、e_i3为第i颗卫星与载体之间的方向余弦，v_ρ(t)为伪距量测噪声，R_N为地球参考椭球卯酉圈曲率半径；

伪距率量测方程为：

$Z_{\stackrel{\·}{\mathrm{\ρ}}}\left(t\right)=H_{\stackrel{\·}{\mathrm{\ρ}}}\left(t\right)X\left(t\right)+v_{\stackrel{\·}{\mathrm{\ρ}}}\left(t\right)$

其中 为伪距率量测噪声。

7.如权利要求6所述的天文辅助的惯性/北斗紧组合导航方法，其特征在于：所述步骤4的详细内容是：

(41)当可见星小于4时，设定量测噪声随着可见星的增多而逐渐减小，即可见星越多则越信赖北斗导航系统的量测信息，伪距和伪距率的量测噪声表示如下：

$v_{\mathrm{\ρ}}\left(t\right)=\frac{\mathrm{\α}}{n}\×v_{\mathrm{\ρ}0}$

$v_{\stackrel{\·}{\mathrm{\ρ}}}\left(t\right)=\frac{\mathrm{\α}}{n}\×v_{\stackrel{\·}{\mathrm{\ρ}}0}$

式中n代表可见星数目，v_ρ0，分别表示伪距、伪距率的初始量测噪声，α表示一个常系数；

(42)当北斗可见星大于4时，首先根据导航位置和卫星位置计算卫星的高度角和方位角，通过卫星的高度角和方位角计算GDOP值，最终对伪距和伪距率的量测噪声进行调节，过程如下：

a、卫星高度角、方位角计算

用户导航位置在地固坐标系下的位置为：

P_u＝[p_ux p_uy p_uz]

第i颗卫星在地固标系下的位置为：

P_s＝[p_sx p_sy p_sz]

设

$E=\left[\begin{array}{ccc}-\frac{p_{ux}}{p_1}& \frac{p_{ux}}{p_1}& 0\\ -\frac{x_{ux}{x}_{uz}}{p_1{p}_2}& -\frac{x_{uy}{x}_{uz}}{p_1{p}_2}& \frac{p_1}{p_2}\\ \frac{x_{ux}}{p_2}& \frac{x_{uy}}{p_2}& \frac{x_{uz}}{p_2}\end{array}\right]$

D＝[d₁ d₂ d₃]＝(P_s-P_u)E

$s=\frac{d_2}{\sqrt{d_1^2+d_2^2+d_3^2}}$

则高度角为：

$El=\left\{\begin{array}{cc}\frac{\mathrm{\π}}{2},& s=1\\ a\tan \frac{s}{\sqrt{1-s^2}},& s\≠1\end{array}\right.$

方位角为：

b、几何精度因子GDOP的计算

设某时刻北斗卫星导航系统的可见星数目为n，设矩阵Q如下：

$Q=\left[\begin{array}{cccc}\sin {\mathrm{El}}_1& \cos {\mathrm{El}}_1\sin {\mathrm{Az}}_1& \cos {\mathrm{El}}_1\cos {\mathrm{Az}}_1& -1\\ \sin {\mathrm{El}}_2& \cos {\mathrm{El}}_1\sin {\mathrm{Az}}_1& \cos {\mathrm{El}}_1\cos {\mathrm{Az}}_1& -1\\ .& .& .& .\\ .& .& .& .\\ .& .& .& .\\ \sin {\mathrm{El}}_n& \cos {\mathrm{El}}_1\sin {\mathrm{Az}}_1& \cos {\mathrm{El}}_1\cos {\mathrm{Az}}_1& -1\end{array}\right]$

则GDOP的值由下式计算：

$GDOP=\sqrt{trace{\left(Q^{T}Q\right)}^{-1}}$

c、卫星到接收机伪距、伪距率量测噪声的设置

$v_{\mathrm{\ρ}}\left(t\right)=\frac{{\mathrm{GDOP}}_k}{{\mathrm{GDOP}}_0}\×v_{\mathrm{\ρ}0}$

$v_{\stackrel{\·}{\mathrm{\ρ}}}\left(t\right)=\frac{GDOP}{{\mathrm{GDOP}}_0}\×v_{\stackrel{\·}{\mathrm{\ρ}}0}$

其中v_ρ0，分别表示伪距、伪距率的初始量测噪声，GDOP₀表示北斗导航系统的初始几何精度因子，GDOP_k表示北斗导航系统第k时刻的几何精度因子。

8.如权利要求5所述的天文辅助的惯性/北斗紧组合导航方法，其特征在于：所述步骤3中，SINS/CNS子系统量测为惯导给出的姿态角与天文给出的姿态角信息之差，定义观测量为星敏感器测得的载体姿态角与捷联惯导的姿态角之差，量测方程如下：

Z_s(t)＝H_s(t)X(t)+v_s(t)

由于组合导航系统状态方程中的误差角状态量为惯性导航平台误差角，它描述的是数学平台与地里面坐标系之间的误差角关系，因此需要将惯性姿态误差角转化到平台误差角；在平台误差角为小量时，量测矩阵由以下公式表示：

$H_s=\frac{-1}{cos\mathrm{\θ}}\left[\begin{array}{ccc}sin\mathrm{\ψ}& cos\mathrm{\ψ}& 0\\ cos\mathrm{\ψ}cos\mathrm{\θ}& -sin\mathrm{\ψ}cos\mathrm{\θ}& 0\\ sin\mathrm{\ψ}sin\mathrm{\θ}& cos\mathrm{\ψ}sin\mathrm{\θ}& -\cos \mathrm{\θ}\end{array}\right]$

因此，量测方程表示为：

式中，γ、θ、ψ分别为横滚角、俯仰角和航向角，v_s(t)为姿态误差角量测噪声。

9.如权利要求5所述的天文辅助的惯性/北斗紧组合导航方法，其特征在于：所述步骤4中，滤波过程如下：

${\hat{X}}_{k/k-1}={\mathrm{\Φ}}_{k/k-1}{\hat{X}}_{k-1/k-1}$

$P_{k/k-1}={\mathrm{\Φ}}_{k/k-1}{P}_{k-1/k-1}{\mathrm{\Φ}}_{k/k-1}^T+{\mathrm{\Γ}}_{k/k-1}{Q}_{k-1}{\mathrm{\Γ}}_{k/k-1}^T$

$K_k=P_{k/k-1}{H}_k^T{(H_k{P}_{k/k-1}{H}_k^T+R_k)}^{-1}$

${\hat{X}}_{k/k}={\hat{X}}_{k/k-1}+K_k(Z_k-H_k{\hat{X}}_{k/k-1})$

$P_{k/k}=(I-K_k{H}_k)P_{k/k-1}{(I-K_k{H}_k)}^T+K_k{R}_k{K}_k^T$

上式中，表示t_k-1时刻的状态对t_k时刻的状态的最优估计值，Φ_k/k-1表示t_k-1时刻至t_k时刻系统的状态转移矩阵，表示t_k-1时刻的系统状态估计值，K_k表示增益矩阵，H_k表示滤波模块t_k时刻的观测系数矩阵，P_k/k-1表示最优预测估计误差协方差阵，P_k/k表示滤波模块t_k时刻的系统误差协方差阵，Q_k-1表示t_k-1时刻的噪声方差矩阵，Γ_k/k-1表示t_k-1时刻的噪声矢量对t_k时刻状态矢量影响的噪声系数矩阵，R_k表示t_k时刻的量测方差矩阵，I为单位矩阵。