基于局部重力场模型的惯性/重力组合导航方法

迟凤阳，宋海岩，王秋滢

(1.黑龙江工程学院 电气与信息工程学院，哈尔滨 150050；2.哈尔滨工程大学 水声工程学院，哈尔滨 150001)

对于水下航行器而言，导航定位精度是其重要的技术指标之一，目前水下航行器的导航方式主要以惯性导航为主。然而，惯性导航系统的定位误差随时间累积而发散，单独使用惯性导航系统很难完成长距离、高精度的导航任务。因此，如何修正惯性导航的误差成为关键问题[1-2]，无论是天文导航、卫星导航或者是声学导航等方式，都在一定程度上减弱了水下航行器的隐蔽性。随着水下导航技术的不断进步，许多国家都在研制以无源导航作为惯性导航辅助手段的水下组合导航系统，其中，重力异常辅助导航是一种常用的无源辅助导航方式[3-4]。重力辅助导航不需要接近水面工作，不受外界的干扰，能够满足水下航行器长时间水下作业的需求。文献[5-6]的研究结果表明,在有重力辅助导航修正的情况下，水下潜器的导航定位误差有所降低。

文中以高精度的地球重力场模型为基础，利用局部实测重力异常数据进行融合，构建局部精确的地球重力场模型，进而得到高精度的重力异常分布图。同时，针对惯性/重力组合导航系统的量测方程无法精确表达的问题，提出利用粒子滤波算法对惯导的定位误差进行修正。为了保证滤波算法的精度，引入正则化方法，增加粒子的多样性。最后，通过仿真验证了正则粒子滤波算法(RPF)在惯性/重力组合导航系统中的有效性。

1 基本原理

1.1 移去-恢复技术

移去-恢复技术[7-8]有3个阶段：第一阶段，利用实测的重力异常值与重力场模型计算得到的重力异常值做差，得到残差重力异常；第二阶段，将残差重力异常通过修正得到新的重力异常数据；第三阶段，利用更新后的重力异常数据得到新的重力场模型。

移去过程：Δg=ΔgF-ΔgG M，

变换过程：Δg→Δg′，

(1)

1.2 重力位系数数值积分算法

(2)

(3)

(4)

(5)

(6)

(7)

(8)

2 局部重力场模型的构建

根据移去-恢复技术，残差重力异常表示为

(9)

(10)

(11)

(12)

(13)

根据式(13)得到新的重力位系数，再利用新的位系数计算得到重力异常，重复上述过程，直到残差重力异常不再发生明显的变化，就可以得到适合于局部地区的重力场模型。图1为重力位系数迭代更新过程。

图1 重力位系数迭代更新过程

3 惯性/重力组合导航系统

3.1 工作原理

在惯性/重力组合导航系统中，需要预先绘制出高精度的重力异常基准图，利用惯导输出的经纬度信息和经纬度误差预测值，在高精度的重力异常基准图中得到重力异常的预测值。与此同时，通过重力测量仪实时获得水下航行器所在位置处的重力异常量测信息，将两种重力异常信息的差值作为量测信息，并通过组合导航滤波算法完成对经纬度误差的估计[11-12]，从而校正惯导输出的位置信息。惯性/重力组合导航系统的原理如图2所示。

图2 惯性/重力组合导航系统原理

3.2 组合导航系统模型

1)系统状态方程。在组合导航系统中，通常以惯导系统的速度误差方程、位置误差方程和姿态误差方程作为系统的状态方程，这里假设姿态误差角都为小角度，且在推导过程中忽略二阶小量，则可以得到惯导系统的线性状态方程。假设惯导系统采用东北天地理坐标系作为导航坐标系，不考虑高度通道，φx，φy，φz为失准角;δVx，δVy为速度误差;δλ，δL为经纬度误差;εx，εy，εz为陀螺漂移;▽x,▽y为加速度计零偏。

惯导系统的速度误差方程可以表示为

(14)

姿态误差方程可以表示为

(15)

位置误差方程可以表示为

(16)

(17)

式中:F为系统状态转移矩阵;G为系统噪声矩阵，具体形式见式(14)～(16)。

2)系统量测方程。将重力仪输出的重力异常值作为该组合导航系统的量测量，重力异常值与载体位置某种函数关系为

Gg(φk,λk)=Gt(φk,λk)+vg.

(18)

式中：Gg(φk,λk)为重力仪输出的重力异常量测值；(φk,λk)为载体真实位置处的纬度和经度；Gt(φk,λk)为真实位置处的重力异常值；vg为重力仪量测噪声。

然而，式(18)之间的函数关系无法用解析式进行具体表达，可以考虑采用粒子滤波算法来解决组合导航系统的位置估计问题。由于重力异常基准图存在误差，则由重力异常基准图估算(φk,λk)处的重力异常Gm(φk,λk)可以表示为

Gm(φk,λk)=Gt(φk,λk)+vm.

(19)

(20)

由式(18)～(20)可以得到量测方程

(21)

这里假设一维量测量可以表示为z，量测噪声v=vg-vm，则式(21)可以表示为

(22)

4 改进的粒子滤波算法

标准的粒子滤波算法在重采样过程中容易出现样本枯竭的现象，最终影响到算法的估计精度。因此，可以考虑采用正则化的方法来解决粒子滤波所存在的样本枯竭问题[13-14]。

正则粒子滤波(RPF)算法流程归纳如下:

3)量测更新。在k时刻更新粒子权值。

5)时刻k=k+1，转到1)。

将RPF算法作为组合导航系统的滤波算法，在时间更新阶段可以得到位置误差预测值，进而对惯导解算的位置信息进行修正，利用修正后的预测位置在重力异常图中确定相应位置的重力异常预测值。在量测更新阶段，取重力异常测量值和重力异常预测值之间的差值作为量测值，更新各粒子的权重，对系统的状态变量进行估计。在重采样阶段，利用一个连续型概率分布去平滑地拟合一个离散型分布序列，然后从这个连续型分布中抽取Ns个粒子，这样避免重采样过程中高权值的粒子被多次复制而产生样本枯竭问题。

5 仿真分析

5.1 基于局部数据改进重力场模型

美国国家地理空间情报局(NGA)经过多年的研究和准备工作，于2008年4月公布了新的全球重力场模型EGM2008，该模型利用新的全球重力测量数据和先进的数据构建方法，对以往的重力场模型进行较大的修正，模型的空间分辨率达到5′×5′，是目前公认精度较好的地球重力场模型。利用EGM2008作为初始模型[15]，按照局部重力异常数据融合方法计算出局部重力场模型。其中，实测重力异常数据由NOAA提供，范围在115°E～135°E、20°N～40°N之间，重力异常变化值在-265.8～266.1 mGal之间。表1给出了局部改进重力场模型及EGM2008模型与实测重力异常数据比较的统计结果。其中，局部改进重力场模型的计算值与实测值的误差大大减小，标准差由原来的16.47 mGal下降到5.26 mGal，精度提高了近3倍，误差起伏由原来的223.9 mGal下降到86.7 mGal。

表1 局部改进重力场模型及EGM2008模型

图3表示范围在115°E～135°E、20°N～40°N之间地区的实测重力异常等值线，图4表示该地区局部重力场模型计算的重力异常等值线，图5表示由EGM2008地球重力场模型计算的重力异常等值线。其中，图3与图4的等值线形状十分接近，局部区域得到一定程度的细化，而图5的等值线图形与图3的等值线图形差距较大，不能完全反应该地区的重力异常分布情况。

图3 实测重力异常等值线

图4 局部重力场模型重力异常等值线

图5 EGM2008模型重力异常等值线

从表1和图3～5可以看出，利用局部精确数据改进得到的重力场模型能够精确地反映出该区域的重力异常分布情况，进而可以为惯性/重力组合导航系统提供精确的重力异常基准图。

5.2 RPF算法在组合导航系统中的应用

在组合导航的仿真分析中，惯性导航的输出数据由惯导模拟器产生，设置水下航行器的初始位置纬度φ=31.2°、经度λ=122°，三轴姿态角分别为0°、0°、225°，以艏艉向3 m·s-1的速度匀速航行。陀螺常值漂移为0.05°·h-1，加速度计零偏为10-4g。惯导数据更新率为0.1 s，粒子滤波周期为10 s，粒子数为500个。

由组合导航系统的工作原理可知，重力异常预测值可以在重力异常基准图中进行查找，而仿真中所用到的重力异常基准图通过前文仿真中得到的局部重力场模型计算获得，通过双线性插值将基准重力异常图分辨率提高为1′×1′，不在网格上的重力异常值采用双线性插值获得。而重力异常量测值通过水下航行器运行轨迹所对应的重力异常值进行表示。

惯性/重力组合导航系统定位误差通过式(23)均方根误差(RMSE)进行计算。

(23)

式中：t为滤波时间；δλk为经度估计误差；δφk为纬度估计误差。

为了验证算法对惯导系统初始位置误差的鲁棒性，分别使惯导单独运行0 min、20 min、30 min和40 min后，再进入组合导航系统工作模式，得到的仿真结果如图6～13所示。

仿真条件1：惯导单独运行0 min,即仿真开始阶段就进入组合导航系统工作模式，无初始位置误差。

从图6和图7可以看出，惯导系统单独解算的位置误差随时间累积，无法满足高精度导航的需求，通过RPF粒子滤波算法可以准确地估计出位置误差，经过校正后组合导航系统输出的轨迹更接近真实轨迹。根据式(23)可以计算出组合导航系统的定位误差为252.63 m。

图6 位置误差估计

图7 校正后输出轨迹

仿真条件2：惯导单独工作20 min后进入组合导航系统工作模式，初始位置误差约为0.7 n mile。

从图8～9可以看出，RPF粒子滤波算法同样可以准确估计出位置误差，组合导航系统的定位误差为365.34 m，相比于无初始位置误差的情况，系统的定位误差要大一些。

图8 位置误差估计

图9 校正后输出轨迹

仿真条件3：惯导单独工作30 min后进入组合导航系统工作模式，初始位置误差约为1.8 n mile。

图11 校正后输出轨迹

从图10～11可以看出，RPF粒子滤波算法依然可以准确地估计出位置误差，组合导航系统的定位误差为376.28 m，相比于仿真条件2的情况，系统的定位误差没有太大的区别。

图10 位置误差估计

仿真条件4：惯导单独工作40 min后进入组合导航系统工作模式，初始位置误差约为3.2 n mile。

图13 校正后输出轨迹

从图12～13可以看出，组合导航系统的定位误差为382.21 m，相比于仿真条件2和仿真条件3的情况，系统的定位误差同样没有太大的变化，说明RPF算法不受初始位置误差影响，实现了对组合导航系统位置误差的精确估计。

图12 位置误差估计

6 结束语

文中研究了一种基于局部重力场模型的惯性/重力组合导航方法。首先，利用移去-恢复技术将已知重力场模型的重力位系数进行修正，可以得到适应局部地区的重力场模型，为组合导航系统提供精确的重力异常基准库。其次，将正则粒子滤波算法作为组合导航的滤波估计算法，解决了系统量测模型非线性和无具体解析式的问题。最后，仿真结果说明正则化粒子算法不受初始位置误差的影响，具有很强的鲁棒性，能够有效地提高惯性/重力组合导航系统的定位精度。