靳宇航，王海涌，刘 涛，贾平会，王永海

(1.北京航空航天大学宇航学院，北京100191；2.北京航天万源科技有限公司，北京100176；3.北京航天长征飞行器研究所，北京100076)

0 引言

在导弹的导航系统中,惯性/天文/GPS三者或两两之间的组合导航是目前组合导航的主要研究范畴。其中,惯性/天文组合导航系统、惯性/GPS组合导航系统已经得到应用,但仍有一定的局限性。星敏感器成本高、体积大、价格昂贵,应用时需要主动控制以保障姿态稳定,避免成像拖尾。GPS及北斗导航信号易受干扰,在实战中难以使用。地磁导航具有自主性,长期工作无积累误差,随着器件和算法的发展,地磁导航越来越引起了人们的关注。因此,某些型号应用自主导航、转向惯性和地磁的组合寻求技术解决方案。地磁场的磁强范围在50000nT～60000nT之间,目前已有的三轴磁强计分辨率可达0.1nT,MEMS磁强计的分辨率可以达到10nT左右,可以用于地磁辅助导航。惯性/地磁组合导航主要有基于批处理算法的地磁匹配辅助导航和地磁序贯滤波辅助导航两种方式。前者是将地磁匹配定位算法得到的载体位置估计作为量测输入,再利用滤波技术估计惯导系统的误差,该算法需要存储大量的地磁数据,并且累积一段路径才能完成定位,实时性较差；后者直接利用磁强计的测量值作为滤波器的量测输入,需要的存储空间较小,同时结合递推滤波技术,直接把地磁测量值作为观测量,每获得一个测量值便对系统进行实时更新,具有连续修正能力,导弹和卫星通常采用后者。

雷达高度表(Radar Altimeter,RA)作为一种主动式距离传感器,与气压高度表和激光高度表相比,雷达高度表具有测高范围大、精度高、功耗低、体积小等优点。雷达测距技术已成功应用于轨道测量、弹道修正、卫星定位、交会对接、巡航制导、航天着陆、大地测量、微波遥感等方面。本文利用雷达高度表的测量信息,建立了一种估计导弹位置的方法,结合磁强计的输出进行姿态矫正,用一个观测式就可以包含姿态和位置信息,实现对导弹的组合导航。

1 地磁场模型

地磁场是一个矢量场,它是时间和位置的函数,但是地磁场随着时间变化缓慢,目前普遍采用的全球地磁模型有国际地磁参考场模型(IGRF)和世界地磁模型(WMM)。本文取用IGRF作为研究的背景场,地磁场势函数V表示为:

式中,Re表示地球半径,r表示载体的地心距,λ表示地理经度,θ代表了地理纬度,gmn(t)和hmn(t)是Gauss系数,IGRF模型提供相应的Gauss系数表,Pmn(cosθ)为Schmidt准归一化n次m阶Legendre函数。

地磁场强度B在地磁球面坐标系f与地球坐标系e中的转换关系如图1所示。

从而得到地磁场矢量B在地球坐标系e下的矢量Be,表示为:

2 捷联惯导/地磁/高度表组合系统模型

本文涉及的坐标系有:导弹发射点惯性坐标系l系、地球坐标系e系、地心惯性系i系、本体坐标系b系。 文中以发射点惯性系l系为导航坐标系,建立捷联惯导/地磁/雷达高度表组合导航系统的状态方程和量测方程。

2.1 组合导航系统状态方程

系统状态模型为惯性导航系统的误差模型,惯性器件加速度计、陀螺建模零偏和白噪声的叠加,建立15维状态方程:

其中,状态变量为:

其中,ϕx、ϕy、ϕz为捷联惯导数学平台失准角,δvx、δvy、δvz为发射点惯性系3个轴上的速度误差, δx、 δy、 δz为位置误差,εx、εy、εz为陀螺仪常值漂移,为加速度计常值偏置。系统噪声为:

式中,各项分别为三轴陀螺仪误差和三轴加速度计误差,均建模为Guass白噪声。

2.2 雷达高度表量测模型

导弹在飞行过程中,通过向地球反射表面(海平面或地表平面)发射一系列雷达脉冲,由弹上雷达信号接收机测得回波信号后,可根据时间间隔Δt得到雷达高度表测量高度H′为:

其中,c为电磁波传播速度。为了获得雷达高度表测得的海拔观测高度Ho,需在实测值H′的基础上进行高度补偿,即:

补偿值Δh由当地地表海拔等因素决定。

接下来推导高度表的输出和捷联惯导解算高度的关系,导弹的高度示意图如图2所示。其中,Re为当前所在位置的地球半径,Re0为发射点至地心距离。

在l系Ol－XlYlZl下,位置矢量OlP与发射点水平面夹角θ0满足:

式中,x、y、z分别为惯导位置解算结果在l系下的坐标,则该时刻惯导测得的海拔计算高度Hc为:

Hc与Ho之差ΔH满足如下关系:

其中,

wh为误差项,受导航误差及量测噪声影响。δx、δy、δz分别为捷联解算与标称弹道的位置误差。

2.3 地磁量测模型

(1)姿态误差角与平台失准角的关系

在捷联惯性导航系统中,载体姿态角是通过姿态矩阵计算出来的。理想情况下,导航计算机计算的导航坐标系(î系)应和理想的导航坐标系(l系)一致,然而,实际计算的导航系与理想的导航系之间将产生偏差,对应的误差角为平台失准角,用角矢量表示。 在小角度假设下,l系与î系之间满足如下关系:

式中,[φ×]为平台失准角φ的反对称矩阵。定义姿态误差角矢量η为:

式中,Cbî为捷联解算的姿态矩阵。忽略2阶小量,可得l系下姿态误差角与平台失准角关系:

式中,M为误差角转换矩阵:

(2)地磁观测模型

三轴磁强计测量值Bm可表示为:

磁强计在载体预估状态下,根据地磁场模型计算的地磁场强度值Bc可表示为:

定义上述磁强计测量值Bm与根据地磁场模型计算的地磁场强度值Bc之差为量测Z1:

地磁场矢量Bb与其估计值之间满足如下数学关系:

将式(22)和式(23)代入式(21)中,可得:

由式(25)可知,磁强计测量误差主要由两部分组成:姿态误差引起的测量误差和位置误差引起的测量误差,本文重点研究姿态部分引起的误差。

式(25)右边第一项可以简写成:

根据微分运算法则:

式中,Bl为地磁场矢量在l系下的表示。定义为:

得到:

由于组合系统状态变量中代表姿态信息的是平台失准角φ,而式(31)中的量测是用姿态误差角η表示的,将前面推导的它们的关系代入式(31)中,得到:

因此,地磁量测模型为:

其中,量测矩阵H1满足:

2.4 系统观测方程

在推导了雷达高度表和磁强计的观测模型之后,得到如下的系统观测方程:

式中,[η×]为姿态误差角η的元素组成的反对称矩阵。将式(28)和式(29)代入式(26),可得:

Z(t)包含地磁误差Z1和高度误差Z2两部分,。 其中,H1的 表达如式(34)所示,

V为观测噪声,仿真设定参考实际器件误差。

2.5 滤波器设计

本文推导得到的测量方程是非线性方程,采用工程上常用的扩展Kalman滤波方法,并采用反馈矫正滤波器。图3为组合导航系统原理图。H2表达式如下所示:

扩展Kalman滤波方程为以下形式:

3 仿真分析

3.1 仿真条件

1)导弹发射点位于北纬 39.98°,东经116.34°,采用垂直向东发射,主动段飞行时间为160s,全程飞行时间为1000s。

2)初始俯仰角为90°,偏航角为0°,横滚角为0°；俯仰角误差为3′,偏航角误差为6′,横滚角误差为 3′。

3)采样周期:陀螺及加速度计为0.01s,磁强计和高度表为0.1s,Kalman滤波周期为0.1s。

4)元器件参数:陀螺常值漂移为0.003(°)/s,陀螺测量白噪声标准差为0.0003(°)/s；加速度计常值偏置为1mg,加速度计白噪声标准差为0.1mg；磁强计测量误差为100nT；雷达高度表测量误差为50m。

3.2 仿真结果及分析

捷联惯性导航及与磁强计高度表组合导航相应的位置误差、速度误差、姿态角误差如图4～图8所示。经过1000s的飞行,纯捷联解算的位置误差为44.58km,速度误差为52.40m/s,俯仰角误差为3.21°,偏航角误差为3.07°,横滚角误差为2.94°。利用磁强计和高度表进行组合导航之后,上述误差的发散均得到有效抑制,位置误差为2.68km,速度误差为8.4m/s,俯仰角误差为0.31°,偏航角误差为0.25°,横滚角误差为0.19°。

由仿真结果可知,加上磁强计和高度表之后进行组合,导航精度得到明显的提高。因为测量方程中既包含了姿态误差,又包含了位置误差,引入了磁强计和高度表的测量值,采用滤波估计的方法能够较好地估计系统的误差。

4 结论

本文基于雷达高度表和磁强计的测量信息,提出一种导弹捷联惯性/地磁/高度表组合导航方法。以捷联惯导误差方程为基础建立系统的状态模型,以磁强计测量值、地磁场模型计算的地磁场强度之差和高度表测量与捷联惯导解算高度之差共同作为量测,只用一个观测表达式即同时包含载体的位置和姿态信息。引入状态反馈,通过混合校正的Kalman滤波进行信息融合,得到系统导航信息的最优估计。仿真结果表明,该算法能有效抑制捷联解算误差的发散。值得一提的是,惯性器件和磁强计可以选用低成本、中等精度的MEMS器件,经过组合导航之后可以到达较高的精度,该方法完全自主,具有一定工程应用价值。