CN104165640A - 基于星敏感器的近空间弹载捷联惯导系统传递对准方法 - Google Patents

Abstract

本发明公开了一种基于星敏感器的近空间弹载捷联惯导系统传递对准方法，包括以下步骤：1)以载体发射点处的惯性坐标系(简称发射点惯性坐标系)为导航坐标系，以待发射弹体上的捷联惯性导航系统(SINS)为子惯导，建立弹载捷联惯性导航系统传递对准状态方程；2)弹载捷联惯导系统导航信息和观测量的计算；3)量测方程的建立；4)根据建立的状态方程和量测方程，利用稀疏网格求积分卡尔曼滤波器估计弹体的数学平台失准角、速度误差、位置误差、安装误差和载体的挠曲变形，对子惯导系统进行修正，完成传递对准过程。

Description

基于星敏感器的近空间弹载捷联惯导系统传递对准方法

技术领域

本发明涉及组合导航传递对准技术领域，具体涉及一种基于星敏感器的近空间弹载捷联惯导系统传递对准方法。

背景技术

临近空间载体是21世纪航空航天技术领域新的制高点。捷联惯性导航系统是临近空间载体的主要导航设备，是一种完全自主的导航技术，具有短时精度高、输出连续、抗干扰能力强、导航信息全面等优势，但是缺点是导航误差随时间积累，临近空间载体的弹体子惯导系统在空中应急启动时需要利用主惯导或者其他辅助设备(如星敏感器)的输出信息进行传递对准，以达到自主、快速、高精度启动的目的。以星敏感器为观测手段的天文导航(CNS)，主要利用恒星来进行导航，具有隐蔽性好、自主性强、精度高和容易受天候影响的特点。SINS/CNS组合导航系统具有良好的定姿性能，在航空航天领域得到广泛应用，但目前主要存在以下问题：

(1)一般以当地地理坐标系为导航坐标系，利用星敏感器的测量平台坐标系相对地理坐标系的姿态角差值，作为组合系统滤波器的观测值，对惯导系统进行修正。此类方法的缺点是选取地理坐标系为导航坐标系，不便于考虑万有引力及地球非球形引力摄动影响，不适用于临近空间载体的运动学分析。

(2)选择当地地理系为作为星光矢量的参考坐标系，就需要将星敏感器测量得的载体坐标系的星光矢量转换到地理坐标系下，在此过程中不可避免的引入误差，使滤波器性能下降。

基于此，需要研究一种模型更简单直观、精度更高的姿态匹配传递对准方法。

发明内容

发明目的：本发明的目的在于解决现有技术中没有综合考虑载体运动学、星光测量参考系带来精度低问题，本发明基于星敏感器的近空间弹载捷联惯性导航系统传递对准方法，以发射点惯性坐标系为导航坐标系，直接利用星敏感器输出的高度角和方位角信息；建立任意失准角的状态方程和量测方程；综合考虑了系统的安装误差、杆臂和挠曲变形，建立了更全面的传递对准；利用稀疏网格求积分滤波器对子惯导系统的导航参数和惯性器件的误差进行修正和估计。

技术方案：本发明一种基于星敏感器的近空间弹载捷联惯导系统传递对准方法，包括以下步骤：

1)以载体发射点处的惯性坐标系(简称发射点惯性坐标系)为导航坐标系，以待发射的弹体上的捷联惯性导航系统(SINS)为子惯导，建立弹载捷联惯性导航系统传递对准状态方程；

2)弹载捷联惯导系统导航信息和观测量的计算。根据弹体SINS解算得到的姿态、位置以及星敏感器识别的导航星，查询导航星历中的天体格林时角和赤纬，得到弹体坐标系下子惯导解算的导航星高度角和方位角，与星敏感器输出的导航星高度角和方位角比较，得到导航星高度角误差和方位角误差；

3)量测方程的建立。利用子惯导解算的位置误差和姿态误差对星敏感器进行补偿，获得导航星高度角和方位角，建立星敏感器导航星高度角误差和方位角误差量测方程；

4)根据建立的状态方程和量测方程，利用稀疏网格求积分卡尔曼滤波器估计弹体的数学平台失准角、速度误差、位置误差、安装误差和载体的挠曲变形，对子惯导系统进行修正，完成传递对准过程。

进一步地，所述步骤1)基于星敏感器的近空间弹载捷联惯导系统传递对准方法，具体为：

状态变量X为

包括失准角φ_xφ_yφ_z、速度误差δv_xδv_yδv_z、位置误差δs_xδs_yδs_z、陀螺仪漂移误差ε_xε_yε_z、加速度计偏置误差安装误差μ_xμ_yμ_z、子惯导挠曲变形

24维的状态方程为 $\stackrel{\·}{X}=F\left(X\right)+\mathrm{Gw};$

系统状态方程的建立：

(1)数学平台失准角误差方程

${\stackrel{\·}{\mathrm{\φ}}}^n=-C_w^{-1}{C}_b^{\hat{i}}{\mathrm{\δ\ω}}_{\mathrm{ib}}^b$

其中：φ＝[φ_xφ_yφ_z]^T；

$C_w^{-1}=\frac{1}{\cos {\mathrm{\φ}}_y}\left[\begin{array}{ccc}\cos {\mathrm{\φ}}_y& \sin {\mathrm{\φ}}_x\sin {\mathrm{\φ}}_y& \cos {\mathrm{\φ}}_x\sin {\mathrm{\φ}}_y\\ 0& \cos {\mathrm{\φ}}_x{\cos \mathrm{\φ}}_y& -\sin {\mathrm{\φ}}_x\cos {\mathrm{\φ}}_y\\ 0& \sin {\mathrm{\φ}}_x& \cos {\mathrm{\φ}}_x\end{array}\right]$

i是发射点惯性坐标系，此处也是导航坐标系；

是惯导解算的发射点惯性坐标系，即数学平台坐标系；

b是弹体坐标系，即子惯导坐标系；

为子惯导解算的姿态矩阵，表示弹体坐标系b到数学平台坐标系的姿态转换矩阵；

为陀螺测量误差；

(2)速度误差方程

惯性坐标系下速度误差微分方程为，

$\mathrm{\δ}{\stackrel{\·}{V}}^i=(I-C_{\hat{i}}^i)C_b^{\hat{i}}{f}^b+C_{\hat{i}}^i{C}_b^{\hat{i}}\mathrm{\δ}{f}^b+\mathrm{\δ}{g}^i$

其中：δVⁱ＝[δv_x δv_y δv_z]^T；

f^b是载体子惯导IMU的比力测量值；

δf^b是子惯导IMU的比力测量误差；

是系至i系的变换阵；

δgⁱ是重力加速度误差；

(3)位置误差方程

惯性坐标系下位置误差δS微分方程为，

$\mathrm{\δ}{\stackrel{\·}{S}}^i={\mathrm{\δV}}^i$

其中：δSⁱ＝[δs_x δs_y δs_z]^T；

(4)姿态变换矩阵

式中：bm是载体坐标系即载机坐标系或星敏感器坐标系；

bs是子惯导坐标系；bh是子惯导水平坐标系；

是bh系到bm系的变换矩阵；是bs系到bh系的变换矩阵；

是主惯导的姿态矩阵；

μ是子惯导安装误差角，是机翼挠曲变形角

(5)安装误差与机翼挠曲变形误差

安装误差方程为

机翼挠曲变形引起的子惯导水平坐标系bh相对载机坐标系bm的挠曲变形角为：模型为：

其中：的方差为σ＝[σ_x σ_y σ_z]^T；η＝[η_x η_y η_z]^T为白噪声，其方差为Q_η＝[Q_ηx Q_ηy Q_ηz]^T，即η～N(0,Q_η)；β＝[β_x β_y β_z]^T为常数。

所述步骤2)具体为：

弹载捷联惯导系统导航信息和观测量的计算：

(1)发射点惯性坐标系下弹载捷联惯导系统导航信息计算；

惯性导航子系统利用惯性测量单元测量弹体的加速度信息和角速度信息，由SINS解算单元给出弹体的发射点惯性坐标系的位置信息和姿态信息；根据加速度信息、角速度信息、万有引力和地球非球形引力摄动，解算弹体的位置信息和姿态信息；结合导航星天体的格林时角(GHA)、赤纬(DEC)信息，转换成惯导解算的高度角(H_cb)和方位角(A_cb)。

(2)观测量的计算；

星敏感器观测设备固定在载体上，通过星敏感器对星光进行跟踪观测，输出导航星在星敏感器坐标系下(即载机坐标系)的高度角H_b和方位角A_b；

实际测量的对应某颗导航星的高度角H_ob和方位角A_ob；

H_ob＝H_b+v_h A_ob＝A_b+v_a

式中：v_h、v_a为星敏感器角度测量噪声；H_ob、A_ob为高度角方位角的测量值。

由于子惯导存在位置误差、速度误差、安装误差和杆臂挠曲效应，由子惯导解算的高度角H_cb和方位角A_cb中存在计算误差：

H_cb＝H_b+δh A_cb＝A_b+δa

观测量为子惯导解算的高度角方位角和星敏感器实际测量的高度角方位角之差：

δh＝H_cb-H_ob+v_h δa＝A_cb-A_ob+v_a

所述步骤3)具体为：

量测方程的建立；

量测信息来源于两部分：星敏感器输出的导航星高度角和方位角以及子惯导解算的导航星高度角和方位角；

考虑子惯导解算误差，星敏感器坐标系的高度角和方位角误差主要由子惯导解算的位置误差和姿态误差等造成；通过补偿获得星敏感器的高度角H_obl和方位角A_obl；

星敏感器高度角方位角传递对准量测方程为：

ΔH_b＝H_cb-H_obl

ΔA_b＝A_cb-A_obl；

所述步骤4)具体为：基于稀疏网格求积分卡尔曼滤波器的传递对准系统信息融合；

该步骤是利用子惯导的误差方程作为状态方程，利用高度角和方位角的系统测量方程，将天文导航星敏感器子系统和捷联惯导子系统高度角方位角输出的差值作为观测值，基于稀疏网格求积分卡尔曼滤波器对系统误差进行实时估计，并将估计误差发送到子惯导解算单元，对导航误差进行校正。

进一步地，所述步骤2)中星敏感器的近空间弹载捷联惯导系统传递对准方法子惯导解算高度角H_cb和方位角A_cb具体步骤为：

2.1)由导航星星历和子惯导系统解算的经度Longi_bs和纬度Latit_bs，获得当地地理坐标系下的解算高度角H_ct和方位角A_ct；

H_ct＝arcsin(sinDEC sinLatit_bs+cosDEC cosLatit_bs cost_Ec)

其中，DEC为赤纬，定义子午线角t_Ec，由下式计算得:

$A_{\mathrm{ct}}=\arccos \left(\frac{\sin \mathrm{DEC}-\sin H_{\mathrm{ct}}\sin {\mathrm{Latit}}_{\mathrm{bs}}}{\cos H_{\mathrm{ct}}\cos {\mathrm{Latit}}_{\mathrm{bs}}}\right)$

2.2)利用星光矢量与高度角和方位角之间关系及姿态转换矩阵获得弹体坐标系下的解算高度角H_cb和方位角A_cb；

其过程如下：

(1)由当地地理坐标系下的解算高度角H_ct和方位角A_ct得到当地地理坐标系下的解算星光矢量r_ct

r_ct＝[r_ctx r_ct_y r_ctz]^T＝[cosH_ctcosA_ct sinH_ct cosH_ctsinA_ct]^T

(2)利用姿态转换矩阵得到弹体坐标系的星光矢量r_cb＝[r_cbx r_cby r_cbz]^T

其中： $r_{\mathrm{cb}}=C_{\hat{t}}^b{r}_{\mathrm{ct}}$

姿态转换矩阵 $C_{\hat{t}}^b=C_{\hat{i}}^{\mathrm{bs}}{C}_e^i{C}_{\hat{t}}^e$

其中：是的转置矩阵，是子惯导解算的姿态矩阵；

是地球坐标系e到发射点惯性坐标系i的转换矩阵；

是子惯导解算北天东坐标系到地球坐标系e的转换矩阵，由当地经纬度确定。

(3)利用星光矢量与高度角方位角之间的几何关系及姿态转换矩阵获得弹体坐标系下解算高度角H_cb和方位角A_cb。

$H_{\mathrm{cb}}=\arcsin \frac{r_{\mathrm{cby}}}{\sqrt{r_{\mathrm{cbx}}^2+r_{\mathrm{cby}}^2+r_{\mathrm{cbz}}^2}}$

$A_{\mathrm{cb}}=\arctan \frac{r_{\mathrm{cbz}}}{r_{\mathrm{cbx}}}$ 其中：

进一步地，所述步骤3)中建立星敏感器导航星高度角误差和方位角误差的量测方程为：

3.1)星敏感器观测设备固定在载体上，星敏感器坐标系的高度角和方位角误差主要由惯导解算的位置误差和姿态误差造成；补偿状态变量中发射点惯性坐标系的位置误差，根据发射点惯性坐标系与地球坐标系转换关系求得当地经纬度；

发射点惯性坐标系的位置误差补偿表示为s_il＝s_i-δs_i其中s_il＝[s_xl s_yl s_zl]^Ts_i＝[s_x s_y s_z]^T

[s_xl s_yl s_zl]^T＝[s_x-δs_x s_y-δs_y s_z-δs_z]^T

补偿后地心惯性坐标系的位置其中 是发射点惯性坐标系i到地心惯性坐标系c的转换矩阵；

补偿后地球坐标系的位置其中 是地心惯性坐标系c到地球坐标系e的转换矩阵；

补偿后当地地理经纬度,

${\mathrm{Latit}}_{\mathrm{otl}}=\arctan \frac{s_{\mathrm{el}}^z}{\sqrt{{\left(s_{\mathrm{el}}^x\right)}^2+{\left(s_{\mathrm{el}}^y\right)}^2}}$ ${\mathrm{Longi}}_{\mathrm{otl}}=\arccos \frac{s_{\mathrm{el}}^x}{\sqrt{{\left(s_{\mathrm{el}}^x\right)}^2+{\left(s_{\mathrm{el}}^y\right)}^2}}$

由导航星星历和补偿后的经度Longi_ctl和纬度Latit_ctl，获得补偿后地理坐标系下的高度角H_otl方位角A_otl；

H_otl＝arcsin(sinDEC sinLatit_otl+cosDEC cosLatit_otl cost_Eotl)

其中定义子午线角t_Eotl，由下式计算得：

$A_{\mathrm{otl}}=\arccos \left(\frac{\sin \mathrm{DEC}-\sin H_{\mathrm{otl}}\sin {\mathrm{Latit}}_{\mathrm{otl}}}{\cos H_{\mathrm{otl}}\cos {\mathrm{Latit}}_{\mathrm{otl}}}\right)$

3.2)利用星光矢量r_otl与高度角H_otl方位角A_otl之间的关系及姿态转换矩阵获得载体系下的解算高度角H_obl和方位角A_obl。

其过程如下：

(1)考虑子惯导解算的位置误差，补偿后当地地理坐标系下的解算高度角H_otl和方位角A_otl得到当地地理坐标系下的星光矢量r_otl

r_otl＝[r_otlx r_otly r_otlz]^T＝[cosH_otl cosA_otl sinH_otl cosH_otl sinA_otl]^T

(2)考虑子惯导解算的姿态误差，利用补偿后姿态转换矩阵得到星敏感坐标系即载体坐标系的星光矢量r_obl＝[r_oblx r_obly r_oblz]^T

其中： $r_{\mathrm{obl}}=C_{\hat{\hat{t}}}^b{r}_{\mathrm{otl}}$

姿态转换矩阵 $C_{\hat{\hat{t}}}^b=C_{b\hat{h}}^{b\hat{m}}{C}_{b\hat{s}}^{b\hat{h}}{C}_{\hat{i}}^{\mathrm{bs}}{C}_i^{\hat{i}}{C}_e^i{C}_{\hat{\hat{t}}}^e$

是由滤波器估计的机翼挠曲变形角估计值给出；

是由滤波器估计的安装误差角估计值给出；

是子惯导解算的姿态矩阵；

是地理坐标系到地球坐标系e的转换矩阵，由子惯导解算补偿后的经度Longi_otl和纬度Latit_otl确定；

是惯导解算姿态误差补偿矩阵，表示导航坐标系i到数学平台坐标系之间转换矩阵，由惯导解算后的数学平台失准角确定。(3)利用星光矢量与高度角方位角转换关系，获得星敏感坐标系下高度角H_obl和方位角A_obl；

$H_{\mathrm{obl}}=\arcsin \frac{r_{\mathrm{obly}}}{\sqrt{r_{\mathrm{oblx}}^2+r_{\mathrm{obly}}^2+r_{\mathrm{oblz}}^2}}$

$A_{\mathrm{obl}}=\arctan \frac{r_{\mathrm{oblz}}}{r_{\mathrm{oblx}}}$ 其中：

本发明与现有技术相比，其有益效果是：(1)本发明采用基于星敏感器的高度角和方位角匹配传递对准方法进行姿态误差估计和修正，直接利用星敏感器测量给出载体坐标系下的高度角和方位角信息，模型简单直观，能提供高精度的姿态估计；

(2)本发明以发射点惯性坐标系为导航坐标系，综合考虑万有引力及地球非球形引力摄动影响解算出载体的位置信息和姿态信息，符合临近空间载体飞行高度，具有较好的适应性；

(3)本发明综合考虑了系统的安装误差杆臂挠曲变形等效应，建立任意失准角度的状态方程和量测方程，建立了比较全面的传递对准；

(4)本发明利用稀疏网格求积分卡尔曼滤波器对系统的导航误差进行校正，提高了惯性/星敏感器高度角方位角传递对准的精度。

附图说明

图1为本发明提出的一种基于星敏感器的近空间弹载捷联惯导系统传递对准系统结构框图；

图2为本发明星光矢量和高度角和方位角示意图；

图3为本发明一种基于星敏感器的近空间弹载捷联惯导系统传递对准方法的原理图；

图4为基于星敏感器的近空间弹载捷联惯导系统传递对准方法姿态误差仿真曲线图Ⅰ；

图5为基于星敏感器的近空间弹载捷联惯导系统传递对准方法姿态误差仿真曲线图Ⅱ。

具体实施方式

下面对本发明技术方案进行详细说明，但是本发明的保护范围不局限于所述实施例。

如图1所示，图1中：1—捷联惯导子系统  2—天文导航子系统  3—信息融合子系统

101—捷联惯导单元即弹载捷联惯性测量单元和SINS解算单元  102—弹载捷联子惯导解算的位置信息即经度和纬度

103—地理坐标系下解算高度角方位角

104—计算地理坐标系到弹体坐标系的姿态转换矩阵

105—星敏感器坐标系或载机坐标系与弹载坐标系之间的安装误差估计

106—星敏感器坐标系或载机坐标系与弹坐标系之间的挠曲变形角估计

107—弹体坐标系下解算高度角方位角

GHA—天体格林时角  DEC—天体赤纬  start_GHA—发射点时角

latit_bs—子惯导解算纬度  Longi_bs—子惯导解算经度

μ—子惯导安装误差—子惯导杆臂挠曲变形角

—解算地理坐标系到弹体坐标系转换矩阵

图2中：

o_b—坐标系原点，为载体质心

o_bx_b—为载体外壳纵轴对称轴，指向载体的头部

o_by_b—在载体纵向对称面内，垂直于纵轴向上

o_bz_b—按照右手规则确定指向载体的右侧

H_b—星敏感器载机载机坐标系的高度角

A_b—星敏感器载机坐标系的方位角

r_b—星敏感器载机坐标系的星光矢量

H_ob—星敏感器实际测量载机坐标系的高度角

A_ob—星敏感器实际测量载机坐标系的方位角

r_ob—星敏感器实际测量载机坐标系的星光矢量

H_cb—子惯导解算弹体坐标系的高度角

A_cb—子惯导解算弹体坐标系的方位角

r_cb—子惯导解算弹体坐标系的星光矢量

H_otl—补偿地理坐标系的高度角

A_otl—补偿地理坐标系的方位角

r_otl—补偿地理坐标系的星光矢量

H_obl—补偿载机坐标系的高度角

A_obl—补偿载机坐标系的方位角

r_obl—补偿星敏感器载体坐标系的星光矢量

各坐标系定义为：

下标i—发射点惯性坐标系

下标n—导航坐标系，此处是指发射点惯性坐标系

下标—子惯导解算发射点惯性坐标系，即数学平台坐标系

下标g—发射点重力坐标系

下标c—地心惯性坐标系

下标e—地球坐标系

下标b—子惯坐标系导即弹体坐标系

下标bs—子惯导坐标系即弹体坐标系

下标bm—载体坐标系即载机坐标系或星敏感器坐标系

下标bh—子惯导水平坐标系

下标t—北天东当地地理坐标系

下标—惯导解算北天东当地地理坐标系

下标—补偿北天东当地地理坐标系

说明书中提到的载体坐标系即为载机坐标系或星敏感器坐标系。

如图1所示，本发明提出一种基于星敏感器的近空间弹载捷联惯导系统传递对准系统，包括惯性导航子系统1、天文导航子系统2、信息融合子系统3。

导航子系统1包括捷联惯导单元即弹载捷联惯性测量单元和SINS解算单元101，SINS解算单元利用惯性测量单元的输出解算出载体的速度姿态和位置信息，惯导解算102位置经度longi_bs和纬度latit_bs与202星敏感器导航星天体的格林时角GHA、赤纬DEC和发射点时角start_GHA信息，获得地理坐标系解算高度角方位角103，子惯导利用解算地理坐标系到弹体坐标系的矩阵转换单元104、弹体坐标系与星敏感载体坐标系之间的安装误差角105和挠曲变形角106单元，给出解算弹体系下高度角方位角107。

天文导航星敏感器高度角方位角信息子系统2包括大视场星敏感器天文信息观测模块201、导航星星历计算模块202、星敏感器导航系统进行星敏感器天文定姿解算直接提供载体坐标系为参考的高度角方位角差定位解算部分203。该大视场星敏感器同时观测多颗恒星，输出载体的以载机坐标系为参考的星光高度角方位角信息，利用弹载捷联子惯导解算的弹体系下高度角方位角差，将其作为观测量，采用稀疏网格求容积卡尔曼滤波器对捷联子惯导系统的导航误差进行校正。

所述稀疏网格求积分卡尔曼滤波器信息融合子系统3，包括高度角方位角失准角计算单元301和稀疏网格求积分卡尔曼滤波器302，高度方位失准角计算单元301利用捷联惯性导航子系统1和天文导航子系统2提供给的高度方位角求得δH、δA，提供给稀疏网格求积分卡尔曼滤波器302；稀疏网格求容积卡尔曼滤波器以子惯导的误差方程作为状态方程，利用高度角和方位角的系统测量方程，将天文导航星敏感器子系统和捷联惯导子系统高度角方位角输出的差值作为观测值，基于稀疏网格求积分卡尔曼滤波器对系统误差进行实时估计，并将估计误差发送到子惯导解算单元，对导航误差进行校正。。

本发明提出一种基于星敏感器的近空间弹载捷联惯导系统传递对准方法，具体包括以下步骤：

步骤一：系统状态方程的建立

以载体的发射点处的惯性坐标系为导航坐标系，状态变量选取

状态变量X为

包括失准角φ_xφ_yφ_z、速度误差δv_xδv_yδv_z、位置误差δs_xδs_yδs_z、陀螺仪漂移误差ε_xε_yε_z、加速度计偏置误差安装误差μ_xμ_yμ_z、子惯导挠曲变形

24维的状态方程为 $\stackrel{\·}{X}=F\left(X\right)+\mathrm{Gw}.$

系统状态方程的建立：

(1)数学平台失准角误差方程

${\stackrel{\·}{\mathrm{\φ}}}^n=-C_w^{-1}{C}_b^{\hat{i}}{\mathrm{\δ\ω}}_{\mathrm{ib}}^b$

其中：φ＝[φ_x φ_y φ_z]^T；

$C_w^{-1}=\frac{1}{\cos {\mathrm{\φ}}_y}\left[\begin{array}{ccc}\cos {\mathrm{\φ}}_y& \sin {\mathrm{\φ}}_x\sin {\mathrm{\φ}}_y& \cos {\mathrm{\φ}}_x\sin {\mathrm{\φ}}_y\\ 0& \cos {\mathrm{\φ}}_x{\cos \mathrm{\φ}}_y& -\sin {\mathrm{\φ}}_x\cos {\mathrm{\φ}}_y\\ 0& \sin {\mathrm{\φ}}_x& \cos {\mathrm{\φ}}_x\end{array}\right]$

(2)速度误差方程

惯性坐标系下速度误差微分方程为，

$\mathrm{\δ}{\stackrel{\·}{V}}^i=(I-C_{\hat{i}}^i)C_b^{\hat{i}}{f}^b+C_{\hat{i}}^i{C}_b^{\hat{i}}\mathrm{\δ}{f}^b+\mathrm{\δ}{g}^i$

(3)位置误差方程

惯性坐标系下位置误差δS微分方程为，

$\mathrm{\δ}{\stackrel{\·}{S}}^i={\mathrm{\δV}}^i$

(4)姿态变换矩阵

(5)安装误差与机翼挠曲变形误差

安装误差方程为

机翼挠曲变形引起的子惯导水平坐标系bh相对载机坐标系bm的挠曲变形角为：模型为：

步骤二：弹载捷联惯导系统导航信息和观测量的计算：

(1)发射点惯性坐标系下弹载捷联惯导系统导航信息计算；

惯性导航子系统利用惯性测量单元测量弹体的加速度信息和角速度信息，由SINS解算单元给出弹体的发射点惯性坐标系的位置信息和姿态信息；根据加速度信息、角速度信息、万有引力和地球非球形引力摄动，解算弹体的位置信息和姿态信息；结合导航星天体的格林时角(GHA)、赤纬(DEC)信息，转换成惯导解算的高度角(H_cb)和方位角(A_cb)。

具体计算步骤如下：

A.由导航星星历和子惯导系统解算经度Longi_bs和纬度Latit_bs，获得当地地理坐标系下的解算高度角(H_ct)和方位角(A_ct)。

H_ct＝arcsin(sinDEC sinLatit_bs+cosDEC cosLatit_bs cost_Ec)

其中，DEC为赤纬，t_Ec为子午线角，

$A_{\mathrm{ct}}=\arccos \left(\frac{\sin \mathrm{DEC}-\sin H_{\mathrm{ct}}\sin {\mathrm{Latit}}_{\mathrm{bs}}}{\cos H_{\mathrm{ct}}\cos {\mathrm{Latit}}_{\mathrm{bs}}}\right)$

B.利用星光矢量与高度角和方位角之间关系及姿态转换矩阵获得弹体系下的解算高度角(H_cb)方位角(A_cb)。

其过程如下：

a.由当地地理坐标系下的解算高度角(H_ct)方位角(A_ct)得到当地地理坐标系下的解算星光矢量r_ct

r_ct＝[r_ctx r_cty r_ctz]^T＝[cosH_ct cosA_ct sinH_ct cosH_ctsinA_ct]^T

b.利用姿态转换矩阵得到弹体坐标系的星光矢量r_cb＝[r_cbx r_cby r_cbz]^T

其中： $r_{\mathrm{cb}}=C_{\hat{t}}^b{r}_{\mathrm{ct}}$

姿态转换矩阵 $C_{\hat{t}}^b=C_{\hat{t}}^{\mathrm{bs}}{C}_e^i{C}_{\hat{t}}^e$

c.利用星光矢量与高度角方位角之间几何关系及姿态转换矩阵获得弹体系下解算高度角(H_cb)方位角(A_cb)。

$H_{\mathrm{cb}}=\arcsin \frac{r_{\mathrm{cby}}}{\sqrt{r_{\mathrm{cbx}}^2+r_{\mathrm{cby}}^2+r_{\mathrm{cbz}}^2}}$

$A_{\mathrm{cb}}=\arctan \frac{r_{\mathrm{cbz}}}{r_{\mathrm{cbx}}}$ 其中：

惯导存在位置误差、速度误差、安装误差和杆臂挠曲效应，由子惯导解算的高度方位角(H_cb)和(A_cb)中存在计算误差：

H_cb＝H_b+δh A_cb＝A_b+δa

(2)观测量的计算。

星敏感器观测设备固定在载体上，通过星敏感器对星光进行跟踪观测，输出导航星在星敏感器坐标系下的高度角(H_b)和方位角(A_b)；

实际测量的对应某颗导航星的高度角H_ob和方位角A_ob；

H_ob＝H_b+v_h A_ob＝A_b+v_a

式中：v_h、v_a为星敏感器角度测量噪声；H_ob、A_ob为高度角方位角的测量值。

由于子惯导存在位置误差、速度误差、安装误差和杆臂挠曲效应，由子惯导解算的高度角(H_cb)和方位角(A_cb)中存在计算误差：

H_cb＝H_b+δh A_cb＝A_b+δa

观测量为子惯导解算的高度角方位角和星敏感器实际测量的高度角方位角之差：

δh＝H_cb-H_ob+v_h δa＝A_cb-A_ob+v_a

步骤三：量测方程的建立

量测信息来源于两部分：星敏感器输出的导航星高度角和方位角以及子惯导解算的导航星高度角和方位角；

(1)惯导解算的高度角方位角

由导航星星历和子惯导系统解算的经度(Longi_bs)纬度(Latit_bs)及姿态转换矩阵获得弹体坐标系下解算解算高度角(H_cb)方位角(A_cb)。

(2)星敏感器坐标系的高度角方位角量测方程

考虑子惯导解算误差，星敏感器坐标系的高度角和方位角误差主要由子惯导解算的位置误差和姿态误差造成；通过补偿获得星敏感器的高度角(H_obl)和方位角(A_obl)；

具体步骤如下：

A.补偿状态变量中发射点惯性坐标系的位置误差，根据发射点惯性坐标系与

地球坐标系转换关系求得当地经纬度。

发射点惯性坐标系的位置误差补偿表示为s_il＝s_i-δs_i其中s_il＝[s_xl s_yl s_zl]^Ts_i＝[s_x s_y s_z]^T

[s_xl s_yl s_zl]^T＝[s_x-δs_x s_y-δs_y s_z-δs_z]^T

补偿后地心惯性坐标系的位置其中 是发射点惯性坐标系i到地心惯性坐标系c的转换矩阵；

补偿后地球坐标系的位置其中 是地心惯性坐标系c到地球坐标系e的转换矩阵；

补偿后当地地理经纬度，

${\mathrm{Latit}}_{\mathrm{otl}}=\arctan \frac{s_{\mathrm{el}}^z}{\sqrt{{\left(s_{\mathrm{el}}^x\right)}^2+{\left(s_{\mathrm{el}}^y\right)}^2}}$ ${\mathrm{Longi}}_{\mathrm{otl}}=\arccos \frac{s_{\mathrm{el}}^x}{\sqrt{{\left(s_{\mathrm{el}}^x\right)}^2+{\left(s_{\mathrm{el}}^y\right)}^2}}$

由导航星星历和补偿后的经度(Longi_otl)和纬度(Latit_otl)，获得补偿后地理坐标系下的高度角(H_otl)和方位角(A_otl)，

H_otl＝arcsin(sinDEC sinLatit_otl+cosDEC cosLatit_otl cost_Eotl)

其中定义子午线角t_Eotl，由下式计算得:

$A_{\mathrm{otl}}=\arccos \left(\frac{\sin \mathrm{DEC}-\sin H_{\mathrm{otl}}\sin {\mathrm{Latit}}_{\mathrm{otl}}}{\cos H_{\mathrm{otl}}\cos {\mathrm{Latit}}_{\mathrm{otl}}}\right)$

B.利用星光矢量(r_otl)与高度角(H_otl)方位角(A_otl)之间关系及姿态转换矩阵获得载体系下的解算高度(H_obl)和方位角(A_obl)。

其过程如下：

a.考虑子惯导解算的位置误差，补偿后当地地理坐标系下的解算高度(H_otl)方位角(A_otl)得到当地地理坐标系下的星光矢量r_otl

r_otl＝[r_otlx r_otly r_otlz]^T＝[cosH_otl cosA_otl sinH_otl cosH_otl sinA_otl]^T

b.考虑子惯导解算的姿态误差，利用补偿后姿态转换矩阵得到星敏感坐标系即载体坐标系的星光矢量r_obl＝[r_oblx r_obly r_oblz]^T

其中： $r_{\mathrm{obl}}=C_{\hat{\hat{t}}}^b{r}_{\mathrm{otl}}$

姿态转换矩阵： $C_{\hat{\hat{t}}}^b=C_{b\hat{h}}^{b\hat{m}}{C}_{b\hat{s}}^{b\hat{h}}{C}_{\hat{i}}^{\mathrm{bs}}{C}_i^{\hat{i}}{C}_e^i{C}_{\hat{\hat{t}}}^e$

c.利用星光矢量与高度角方位角关系获得弹体坐标下高度角(H_obl)和方位角(A_obl)。

$H_{\mathrm{obl}}=\arcsin \frac{r_{\mathrm{obly}}}{\sqrt{r_{\mathrm{oblx}}^2+r_{\mathrm{obly}}^2+r_{\mathrm{oblz}}^2}}$

$A_{\mathrm{obl}}=\arctan \frac{r_{\mathrm{oblz}}}{r_{\mathrm{oblx}}}$ 其中：

星敏感器高度角方位角传递对准量测方程为

ΔH_b＝H_cb-H_obl

ΔA_b＝A_cb-A_obl

步骤四：基于的星敏感器的近空间弹载捷联惯导系统传递对准稀疏网格求积分卡尔曼滤波器信息融合子系统

信息融合子系统是基于稀疏网格求积分卡尔曼滤波器，利用弹载捷联子惯导的误差方程作为状态方程，利用高度角和方位角的系统测量方程，将星敏感器子系统和弹载捷联惯导子系统高度角方位角输出的差值作为观测值，基于稀疏网格求积分卡尔曼滤波算法对系统误差进行实时估计，并将估计误差发送到捷联惯导解算单元，对导航误差进行校正。

本发明的可行性通过如下仿真加以验证：

(1)陀螺仪随机常值漂移0.5°/h，随机白噪声0.5°/h，加速度计随机常值偏置0.1mg，随机白噪声0.1mg，星敏感器测量误差2″，杆臂长度3m/0.5m/1m，初始高度50KM，初始速度3Ma，初始速度误差0.5m/s，初始位置误差10m；

(2)初始姿态误差：俯仰角30°，航向角30°，横滚角30°；

(3)惯性传感器数据更新周期为5ms，滤波周期为0.1s，仿真时间30s；

通过计算机仿真，采用基于星敏感器的近空间弹载捷联惯导系统传递对准方法如附图4和附图5所示。附图4：安装误差为3′，采用星敏感器的高度角方位角匹配弹载捷联惯导传递对准30秒的姿态误差均值为0.4026′、1.0948′、0.7832′，方差分别为0.3050′、1.0014′、0.6882′；附图5：安装误差为5′，姿态误差均值分别为0.3587′、1.1071′、0.8703′，方差分别为0.3159′、1.0360′、0.7069。

附图4和附图5可见，本发明所提方法对姿态误差估计满足高超声速载体导航系统对姿态测量的高精度需求。

如上所述，尽管参照特定的优选实施例已经表示和表述了本发明，但其不得解释为对本发明自身的限制。在不脱离所附权利要求定义的本发明的精神和范围前提下，可对其在形式上和细节上作出各种变化。

Claims (4)

1.基于星敏感器的近空间弹载捷联惯导系统传递对准方法，其特征在于，包括以下步骤：

1)以载体发射点处的惯性坐标系为导航坐标系，以待发射弹体上的捷联惯性导航系统为子惯导，建立弹载捷联惯性导航系统传递对准状态方程；

2)弹载捷联惯导系统导航信息和观测量的计算：根据弹体SINS解算得到的姿态、位置以及星敏感器识别的导航星，查询导航星历中的天体格林时角和赤纬，得到弹体坐标系下子惯导解算的导航星高度角和方位角，与星敏感器输出的导航星高度角和方位角比较，得到导航星高度角误差和方位角误差；

3)量测方程的建立：利用子惯导解算的位置误差和姿态误差对星敏感器进行补偿，获得导航星高度角和方位角，建立星敏感器导航星高度角误差和方位角误差量测方程；

4)根据建立的状态方程和量测方程，利用稀疏网格求积分卡尔曼滤波器估计弹体的数学平台失准角、速度误差、位置误差、安装误差和载体的挠曲变形，对子惯导系统进行修正，完成传递对准过程。

2.根据权利要求1所述基于星敏感器的近空间弹载捷联惯导系统传递对准方法，其特征在于：

所述步骤1)具体为：

状态变量X为

包括失准角φ_xφ_yφ_z、速度误差δv_xδv_yδv_z、位置误差δs_xδs_yδs_z、陀螺仪漂移误差ε_xε_yε_z、加速度计偏置误差安装误差μ_xμ_yμ_z、子惯导挠曲变形

24维的状态方程为 $\stackrel{\·}{X}=F\left(X\right)+\mathrm{Gw};$

系统状态方程的建立：

(1)数学平台失准角误差方程

${\stackrel{\·}{\mathrm{\φ}}}^n=-C_w^{-1}{C}_b^{\hat{i}}{\mathrm{\δ\ω}}_{\mathrm{ib}}^b$

其中：φ＝[φ_x φ_y φ_z]^T；

$C_w^{-1}=\frac{1}{\cos {\mathrm{\φ}}_y}\left[\begin{array}{ccc}\cos {\mathrm{\φ}}_y& \sin {\mathrm{\φ}}_x\sin {\mathrm{\φ}}_y& \cos {\mathrm{\φ}}_x\sin {\mathrm{\φ}}_y\\ 0& \cos {\mathrm{\φ}}_x{\cos \mathrm{\φ}}_y& -\sin {\mathrm{\φ}}_x\cos {\mathrm{\φ}}_y\\ 0& \sin {\mathrm{\φ}}_x& \cos {\mathrm{\φ}}_x\end{array}\right]$

i是发射点惯性坐标系，此处也是导航坐标系；

是惯导解算的发射点惯性坐标系，即数学平台坐标系；

b是弹体坐标系，即子惯导坐标系；

为子惯导解算的姿态矩阵，表示弹体坐标系b到数学平台坐标系的姿态转换矩阵；

为陀螺测量误差；

(2)速度误差方程

惯性坐标系下速度误差微分方程为，

$\mathrm{\δ}{\stackrel{\·}{V}}^i=(I-C_{\hat{i}}^i)C_b^{\hat{i}}{f}^b+C_{\hat{i}}^i{C}_b^{\hat{i}}\mathrm{\δ}{f}^b+\mathrm{\δ}{g}^i$

其中：δVⁱ＝[δv_x δv_y δv_z]^T；

f^b是子惯导IMU的比力测量值；

δf^b是子惯导IMU的比力测量误差；

是系至i系的变换阵；

δgⁱ是重力加速度误差；

(3)位置误差方程

惯性坐标系下位置误差δS微分方程为，

$\mathrm{\δ}{\stackrel{\·}{S}}^i={\mathrm{\δV}}^i$

其中：δSⁱ＝[δs_x δs_y δs_z]^T；

(4)姿态变换矩阵

式中：bm是载体坐标系即载机坐标系或星敏感器载体坐标系；

bs是子惯导坐标系；bh是子惯导水平坐标系；

是bh系到bm系的变换矩阵；是bs系到bh系的变换矩阵；

是主惯导的姿态矩阵；

μ是子惯导安装误差角，是机翼挠曲变形角；

(5)安装误差与机翼挠曲变形误差

安装误差方程为

机翼挠曲变形引起的子惯导水平坐标系bh相对载机坐标系bm的挠曲变形角为：模型为：

其中：的方差为σ＝[σ_x σ_y σ_z]^T；η＝[η_x η_y η_z]^T为白噪声，其方差为Q_η＝[Q_ηx Q_ηy Q_ηz]^T，即η～N(0,Q_η)；β＝[β_x β_y β_z]^T为常数。

所述步骤2)具体为：

弹载捷联惯导系统导航信息和观测量的计算：

(1)发射点惯性坐标系下弹载捷联惯导系统导航信息计算；

惯性导航子系统利用惯性测量单元测量弹体的加速度信息和角速度信息，由SINS解算单元给出弹体的发射点惯性坐标系的位置信息和姿态信息；根据加速度信息、角速度信息、万有引力和地球非球形引力摄动，解算弹体的位置信息和姿态信息；结合导航星天体的格林时角GHA、赤纬DEC信息，转换成惯导解算的高度角H_cb和方位角A_cb；

(2)观测量的计算；

星敏感器观测设备固定在载体上，通过星敏感器对星光进行跟踪观测，输出导航星在星敏感器坐标系下(即载机坐标系)的高度角H_b和方位角A_b；

实际测量的对应某颗导航星的高度角H_ob和方位角A_ob；

H_ob＝H_b+v_h A_ob＝A_b+v_a

式中：v_h、v_a为星敏感器角度测量噪声；H_ob、A_ob为高度角方位角的测量值。

由于子惯导存在位置误差、速度误差、安装误差和杆臂挠曲效应，由子惯导解算的高度角H_cb和方位角A_cb中存在计算误差：

H_cb＝H_b+δh A_cb＝A_b+δa

观测量为子惯导解算的高度角方位角和星敏感器实际测量的高度角方位角之差：

δh＝H_cb-H_ob+v_h δa＝A_cb-A_ob+v_a

所述步骤3)具体为：

量测方程的建立；

量测信息来源于两部分：星敏感器输出的导航星高度角和方位角以及子惯导解算的导航星高度角和方位角；

考虑子惯导解算误差，星敏感器坐标系的高度角和方位角误差主要由子惯导解算的位置误差和姿态误差造成；通过补偿获得星敏感器的高度角H_obl和方位角A_obl；

星敏感器高度角方位角传递对准量测方程为：

ΔH_b＝H_cb-H_obl

ΔA_b＝A_cb-A_obl；

所述步骤4)具体为：基于稀疏网格求积分卡尔曼滤波器的传递对准系统信息融合；

该步骤是利用子惯导的误差方程作为状态方程，利用高度角和方位角的系统测量方程，将导航星敏感器子系统和捷联惯导子系统高度角方位角输出的差值作为观测值，基于稀疏网格求积分卡尔曼滤波器对系统误差进行实时估计，并将估计误差发送到子惯导解算单元，对导航误差进行校正。

3.根据权利要求1或2所述的基于星敏感器的近空间弹载捷联惯导系统传递对准方法，其特征在于，所述步骤2)中子惯导解算高度角H_cb和方位角A_cb具体步骤为：

2.1)由导航星星历和子惯导解算的经度Longi_bs和纬度Latit_bs，获得当地地理坐标系下的解算高度角H_ct和方位角A_ct；

H_ct＝arcsin(sinDEC sinLatit_bs+cosDEC cosLatit_bs cost_Ec)

其中，DEC为赤纬，t_Ec为子午线角，

$A_{\mathrm{ct}}=\arccos \left(\frac{\sin \mathrm{DEC}-\sin H_{\mathrm{ct}}\sin {\mathrm{Latit}}_{\mathrm{bs}}}{\cos H_{\mathrm{ct}}\cos {\mathrm{Latit}}_{\mathrm{bs}}}\right)$

2.2)利用星光矢量与高度角和方位角之间关系及姿态转换矩阵获得弹体坐标系下的解算高度角H_cb和方位角A_cb；

其过程如下：

(1)由当地地理坐标系下的解算高度角H_ct和方位角A_ct得到当地地理坐标系

下的解算星光矢量r_ct

r_ct＝[r_ctx r_cty r_ctz]^T＝[cosH_ctcosA_ct sinH_ct cosH_ctsinA_ct]^T

(2)利用姿态转换矩阵得到弹体坐标系的星光矢量r_cb＝[r_cbx r_cby r_cbz]^T

其中： $r_{\mathrm{cb}}=C_{\hat{t}}^b{r}_{\mathrm{ct}}$

姿态转换矩阵 $C_{\hat{t}}^b=C_{\hat{i}}^{\mathrm{bs}}{C}_e^i{C}_{\hat{t}}^e$

其中：是的转置矩阵，是子惯导解算的姿态矩阵；

是地球坐标系e到发射点惯性坐标系i的转换矩阵；

是子惯导解算北天东坐标系到地球坐标系e的转换矩阵，由当地经纬度确定；

(3)利用星光矢量与高度角方位角之间的几何关系及姿态转换矩阵获得弹体坐标系下解算高度角H_cb和方位角A_cb；

$H_{\mathrm{cb}}=\arcsin \frac{r_{\mathrm{cby}}}{\sqrt{r_{\mathrm{cbx}}^2+r_{\mathrm{cby}}^2+r_{\mathrm{cbz}}^2}}$

$A_{\mathrm{cb}}=\arctan \frac{r_{\mathrm{cbz}}}{r_{\mathrm{cbx}}}$ 其中：

4.根据权利要求1或2所述的基于星敏感器的近空间弹载捷联惯导系统传递对准方法，其特征在于，所述步骤3)中建立星敏感器导航星高度角误差和方位角误差量测方程为：

3.1)星敏感器观测设备固定在载体上，星敏感器坐标系的高度角和方位角误差主要由子惯导解算的位置误差和姿态误差造成；补偿状态变量中发射点惯性坐标系的位置误差，根据发射点惯性坐标系与地球坐标系转换关系求得当地经纬度；

发射点惯性坐标系的位置误差补偿表示为s_il＝s_i-δs_i其中s_il＝[s_xl s_yl s_zl]^Ts_i＝[s_x s_y s_z]^T

[s_xl s_yl s_zl]^T＝[s_x-δs_x s_y-δs_y s_z-δs_z]^T

补偿后地心惯性坐标系的位置其中 是发射点惯性坐标系i到地心惯性坐标系c的转换矩阵；

补偿后地球坐标系的位置其中 是地心惯性坐标系c到地球坐标系e的转换矩阵；

补偿后当地地理经纬度,

${\mathrm{Latit}}_{\mathrm{otl}}=\arctan \frac{s_{\mathrm{el}}^z}{\sqrt{{\left(s_{\mathrm{el}}^x\right)}^2+{\left(s_{\mathrm{el}}^y\right)}^2}}$ ${\mathrm{Longi}}_{\mathrm{otl}}=\arccos \frac{s_{\mathrm{el}}^x}{\sqrt{{\left(s_{\mathrm{el}}^x\right)}^2+{\left(s_{\mathrm{el}}^y\right)}^2}}$

由导航星星历和补偿后的经度Longi_otl和纬度Latit_otl，获得补偿后地理坐标系下的高度角H_otl和方位角A_otl；

H_otl＝arcsin(sinDEC sinLatit_otl+cosDEC cosLatit_otl cost_Eotl)

其中定义子午线角t_Eotl，由下式计算得:

$A_{\mathrm{otl}}=\arccos \left(\frac{\sin \mathrm{DEC}-\sin H_{\mathrm{otl}}\sin {\mathrm{Latit}}_{\mathrm{otl}}}{\cos H_{\mathrm{otl}}\cos {\mathrm{Latit}}_{\mathrm{otl}}}\right)$

3.2)利用星光矢量r_otl与高度角H_otl方位角A_otl之间的关系及姿态转换矩阵获得载体坐标系下的解算高度角H_obl和方位角A_obl；

其过程如下：

(1)考虑子惯导解算的位置误差，补偿后当地地理坐标系下的解算高度角H_otl和方位角A_otl得到当地地理坐标系下的星光矢量r_otl

r_otl＝[r_otlx r_otly r_otlz]^T＝[cosH_otl cosA_otl sinH_otl cosH_otl sinA_otl]^T

(2)考虑子惯导解算的姿态误差，利用补偿后姿态转换矩阵得到星敏感坐标系即载体坐标系的星光矢量r_obl＝[r_oblx r_obly r_oblz]^T

其中： $r_{\mathrm{obl}}=C_{\hat{\hat{t}}}^b{r}_{\mathrm{otl}}$

姿态转换矩阵 $C_{\hat{\hat{t}}}^b=C_{b\hat{h}}^{b\hat{m}}{C}_{b\hat{s}}^{b\hat{h}}{C}_{\hat{i}}^{\mathrm{bs}}{C}_i^{\hat{i}}{C}_e^i{C}_{\hat{\hat{t}}}^e$

是由滤波器对机翼挠曲变形角估计值给出；

是由滤波器对安装误差角估计值给出；

是子惯导解算的姿态矩阵；

是地理坐标系到地球坐标系e的转换矩阵，由子惯导解算补偿后的经度Longi_otl和纬度Latit_otl和确定；

是惯导解算姿态误差补偿矩阵，表示导航坐标系i到数学平台坐标系之间转换矩阵，由子惯导解算后的数学平台失准角确定；

(3)利用星光矢量与高度角方位角转换关系，获得星敏感坐标系下高度角H_obl和方位角A_obl；

$H_{\mathrm{obl}}=\arcsin \frac{r_{\mathrm{obly}}}{\sqrt{r_{\mathrm{oblx}}^2+r_{\mathrm{obly}}^2+r_{\mathrm{oblz}}^2}}$

$A_{\mathrm{obl}}=\arctan \frac{r_{\mathrm{oblz}}}{r_{\mathrm{oblx}}}$ 其中：