第３１卷第３期 计算机仿真 ２０１４年３月 文章编号：１００６—９３４８（２０１４）０３—００５５—０５
基 于高 斯样条函数的水下重力辅助惯性导航
王志刚ｌＩ ．顾雪峰 （１．中国人民解放军９１５５０部队，辽宁大连１１６０２３；２．海军工程大学兵器工程系，湖北武汉４３００３３； ３．海军工程大学兵器科研部，湖北武汉４３００３３）
摘要：在水下重力辅助惯性导航定位问题的研究中，利用卡尔曼滤波实现水下重力辅助惯性导航，需要对重力量测值及其误 差建模．为此提出一种基于二维高斯样条函数逼近的水下重力辅助导航新方法，将测量重力与惯导指示重力之差表示为连 续的解析形式，解决了利用卡尔曼滤波实现水下重力辅助惯性导航，需要对重力量测值及其误差精确建模的问题。在分辨 率为２ ×２ 的某区域重力异常数据基础上进行仿真分析，局部重力场解析重构后的平均误差小于０．１９ｍＧａｌ，潜器的平均经、 纬定位误差分别小于０．８２和０．２０ｎｍｉｌｅ，为提高惯性导航精度提供了依据。 关键词：重力辅助惯性导航；高斯样条函数；局部重力场；扩展卡尔曼滤波；惯性导航系统 中图分类号：ＴＮ９１ｌ 文献标识码：Ｂ
Ｕｎｄｅｒｗａｔｅｒ Ｇｒａｖｉｔｙ Ａｉｄｅｄ Ｉｎｅｒｔｉａｌ Ｎａｖｉｇａｔｉｏｎ Ｂａｓｅｄ ｏｎ Ｇａｕｓｓ Ｓｐｌｉｎｅ Ｆｕｎｃｔｉｏｎ
ＷＡＮＧ Ｚｈｉ－ｇａｎｇ ，ＧＵ Ｘｕｅ－ｆｅｎｇ （１．Ｔｈｅ ９１５５０ｔｈ Ｔｒｏｏｐ ｏｆ ＰＬＡ，Ｄａｌｉａｎ ｌ１６０２３，Ｃｈｉｎａ； ２．Ｄｅｐａｒｔｍｅｎｔ ｏｆ Ｗｅａｐｏｎ，Ｎａｖａｌ Ｕｎｉｖｅｒｓｉｔｙ ｏｆ Ｅｎｇｉｎｅｅｒｉｎｇ，Ｗｕｈａｎ Ｈｕｂｅｉ ４３００３３，Ｃｈｉｎａ； ３．Ｄｅｐａｒｔｍｅｎｔ ｏｆ Ｓｃｉｅｎｃｅ，Ｎａｖａｌ Ｕｎｉｖｅｒｓｉｔｙ ｏｆ Ｅｎｇｉｎｅｅｒｉｎｇ，Ｗｕｈａｎ Ｈｕｂｅｉ ４３００３３，Ｃｈｉｎａ）
ＡＢＳＴＲＡＣＴ：Ｉｍｐｌｅｍｅｎｔｉｎｇ ｕｎｄｅｒｗａｔｅｒ ｇｒａｖｉｔｙ ａｉｄｅｄ ｉｎｅ￣ｉａｌ ｎａｖｉｇａｔｉｏｎ ｏｎ Ｋａｌｍａｎ ｆｉｌｔｅｒ ｒｅｑｕｉｒｅｓ ｔｈｅ ｍｏｄｅｌｉｎｇ ｏｆ ｇｒａｖｉｔｙ ｍｅａｓｕｒｅｍｅｎｔｓ ａｎｄ ｔｈｅｉｒ ｅｒｒｏｒｓ．Ｄｕｅ ｔｏ ｔｈａｔ ｂａｃｋｇｒｏｕｎｄ，ａ ｎｅｗ ｐａｔｔｅｒｎ ｏｆ ｇｒａｖｉｔｙ ａｉｄｅｄ ｎａｖｉｇａｔｉｏｎ ｂａｓｅｄ ｏｎ ２
－Ｄ Ｇａｕｓｓ ｓｐｌｉｎｅ ｆｕｎｃｔｉｏｎ ｗａｓ ｇｉｖｅｎ ｉｎ ｔｈｉｓ ｐａｐｅｒ．Ｉｎ ｏｕｒ ｎｅｗ ｐａｔｔｅｒｎ，ｔｈｅ ｄｉｆｆｅｒｅｎｃｅ ｂｅｔｗｅｅｎ ｍｅａｓｕｒｅｄ ｇｒａｖｉｔｙ ａｎｄ ｉｎｄｉｃａｔｅｄ ｇｒａｖｉｔｙ ｗａｓ ｔｈｅｎ ｅｘｐｒｅｓｓｅｄ ａｓ ａ ｃｏｎｔｉｎｕａｌ ａｎａｌｙｔｉｃａｌ ｅｑｕａｔｉｏｎ ｔｏ ｍｏｄｅ ｔｈｅ ｇｒａｖｉｔｙ ｍｅａｓｕｒｅｍｅｎｔｓ ａｎｄ ｔｈｅｉｒ ｅｒ． ｒ０ｒｓ．Ａｆｔｅｒ ｔｈａｔ，ｔｈｅ ｅｘｔｅｎｄｅｄ Ｋａｌｍａｎ ｆｉｌｔｅｒ ｃａｎ ｂｅ ｉｍｐｌｅｍｅｎｔｅｄ，ｗｉｔｈ ｔｈｅ ｇｒａｖｉｔｙ ｄｉｆｅｒｅｎｃｅ ｕｓｅｄ ａｓ ｍｅａｓｕｒｅｍｅｎｔ，ｔｏ ｅｓｔｉｍａｔｅ ｔｈｅ ｐｏｓｉｔｉｏｎ ｅｒｒｏｒ ｏｆ ＩＮＳ．Ｆｉｎａｌｌｙ，ａ ｓｉｍｕｌａｔｉｏｎ ｗａｓ ｄｏｎｅ ｏｎ ２ ｘ２ ｇｒａｖｉｔｙ ａｎｏｍａｌｙ ｄａｔａｂａｓｅ，ａｎｄ ｆｒｏｍ ｔｈｅ ｒｅ— ｓｕｈｓ ｗｅ ｃａｎ ｓｅｅ ｔｈａｔ ｔｈｅ ｍｅａｎ ｅｒｒｏｒ ｏｆ ｔｈｅ ｒｅｃｏｎｓｔｒｕｃｔｅｄ ｌｏｃａｌ ｇｒａｖｉｔｙ ｆｉｅｌｄ ｉｓ ｌｅｓｓ ｔｈａｎ ０．１９ｍＧａｌ，ａｎｄ ｔｈｅ ｍｅａｎ ｌｏｃａ・ ｔｉｏｎ ｅｒｒｏｒ ｉｎ ｌｏｎｇｉｔｕｄｅ ａｎｄ ｌａｔｉｔｕｄｅ ｉｓ ０．８２ ａｎｄ ０．２０ ｎａｕｔｉｃａｌ ｍｉｌｅ ｒｅｓｐｅｃｔｉｖｅｌｙ． ＫＥＹＷＯＲＤＳ：Ｇｒａｖｉｔｙ ａｉｄｅｄ ｉｎｅ￣ｉａｌ ｎａｖｉｇａｔｉｏｎ；Ｇａｕｓｓ ｓｐｌｉｎｅ ｆｕｎｃｔｉｏｎ；Ｌａｏｃａｌ ｇｒａｖｉｔｙ ｆｉｅｌｄ；Ｅｘｔｅｎｄｅｄ Ｋａｌｍａｎ ｆｉｌｔｅｒ； Ｉｎｅｒｔｉａｌ ｎａｖｉｇａｔｉｏｎ ｓｙｓｔｅｍ
ｌ 引言 惯性导航系统（ＩＮＳ）的误差随时间积累发散，无法长时 间保持高精度．因此，必须要通过其它导航方式实时或定期 修正ＩＮＳ。出于军事上的隐蔽性要求，水下军用潜器很难利 用卫星导航或无线电导航信息，需要借助另外一种无源的辅 助导航手段。重力匹配辅助惯性导航＿１。］由于具有完全自主 性和高隐蔽性。是很好的选择，因而获得了各国的广泛关注。 目前，实现重力辅助惯性导航主要有两类方法。一类是 目前国内普遍研究较多的相关匹配方法ｒ４ ｊ。这类方法以 收稿Ｅｔ期：２０１３—０９—０３修回日期：２０１３—１０—１２ ＩＣＣＰ算法（迭代最近点等值线算法）为主要代表。另一类是 已经在美国走向实用的卡尔曼滤波方法［３］。文章即是以卡 尔曼滤波方法为主要对象，重点解决应用卡尔曼滤波算法需 要对重力量测值及其误差精确建模的问题。目前常用的各 类全球重力场模型［７］无论在阶次和精度方面都已经达到了 很高的水平，它本质上属于球坐标系下的调和分析，因而应 用过程中需要计算大量的Ｌｅｇｅｎｄｒｅ系数，使得模型计算量巨 大，而这对于定位的实时性是一大考验，另外现有的大部分 全球重力场模型非线性比较严重，也不太利于滤波估计时的 线性化处理。在此情况之下，希望能有一种更简单的地球重 力场模型去描述局部区域的重力变化。同时这种模型的误差
一５５— 统计特性是可知的。 近年来，高斯样条函数在医学影像、数字地形构建、计算 机图形学及地球物理等领域有较为广泛的应用，文献［８］讨 论了一维高斯样条在协方差函数代数确定中的应用，文献 ［９］采用高斯函数作为样条基函数对计算区域重力异常进行 二维整体逼近，文献［１０］讨论了高斯样条函数逼近局部重力 异常的相关影响因素。采用高斯样条函数逼近局部重力异 常运算简单且最终解析式是统一的，文章在该样条逼近函数 的基础上建立了重力差量测方程并结合具体的动基座惯导 编排方程，成功运用扩展Ｋａｌｍａｎ滤波实现了对潜器惯导位 置误差的估计。 ２重力辅助惯性导航原理 重力辅助惯性导航的原理如图１所示，辅助导航系统是 利用高分辨率的重力背景场数据组成的基本特征信息数据 库（数字重力图），结合惯性导航系统给出的导航参数（主要 是位置信息）在数字重力图中找到相应的指示重力值；随后， 海洋重力仪器（或重力梯度仪）实时测量重力数据，测量重力 数据再根据测深测潜仪以及ＩＮＳ提供的导航参数进行相应 的改正处理” ，处理后的实测重力数据与指示重力数据以 及ＩＮＳ导航参数一并送人到中心匹配计算机中进行组合导 航计算。 利用卡尔曼滤波实现重力辅助导航，关键是建立重力观 测值模型，由引言分析知道，目前常用的全球重力场模型，由 于其复杂性，不适合辅助导航的局部应用要求，为此如何建 立一种新的、精确、实用的重力观测值模型是文章将要解决 的问题。 为解决以上问题，文章引入了重力高斯样条函数，这一 函数的引人使得可以将重力差量测值表示为连续的解析形 式，因而可以建立滤波估计所需的量测方程，有了量测方程 再结合惯导误差编排方程就可通过扩展卡尔曼滤波算法对 导航误差进行最优估计。最后，计算的匹配位置信息与ＩＮＳ 导航误差一道送人ＩＮＳ的中心解算计算机中对其进行修正。 图１重力辅助惯性导航原理示意图 ３二维高斯样条函数逼近 ３．１ 利用二维高斯样条函数建立局部重力场模型 一５６一 设高斯样条函数表达式为Ｇ（ ）：ｅｘｐ（一 ），Ｐ，￣－［ＩＸＹ 平面上等分离散格网点集合｛（ Ｙｉ）ｌ ｉ＝１，…，ｍ；ｊ＝１， …， ｝及其对应函数值 （ｉ＝１，…，ｍ；ｊ＝１，…，ｎ），Ｘ方向 一维高斯样条函数解析式为Ｌ ＝ｓｐａｎ｛Ｇｘ（（ 一 ）ｌａｘ）｝，Ｙ 方向一维高斯样条函数解析式为Ｌｙ＝ｓｐａｎ｛Ｇｙ（（Ｙ— Ｙ１）／ａｙ）｝，则ＸＹ平面的二维高斯样条函数可以写成 与Ｌ 的张量积形式： Ｌ（ ，Ｙ）＝ ｏ Ｌｙ ＝ｓｐａｎ｛Ｇ （（ｘ— ｆ）／Ａｘ）Ｇ （（Ｙ一，， ）／／ｔｙ）｝ （Ｉ） 即有
（ ，ｙ）＝∑∑ｃｉｄ ｃｘ（（ — ）／ａｘ） （（ｙ一 ）／△ｙ） （２） 不妨设
Ｐ＝ －Ｘ１ ）， 、（ｘ －ｘｚ ），…， （ ］
＝ Ｙ －Ｙｌ －Ｙ２ ，Ｇｙ（（ｙ －Ｙｍ ）］
（３） 则上式可以简化为矩阵形式：
，ｙ）＝ （４） ＬＣ＝（ｃ“）ｍ 则根据插值条件｛Ｌ（ｘ ，乃）＝ｚＩ ｌ ｉ＝１…．，ｍ ＝ １…．，ｎ｝，有如下线性方程：
＝ Ｙ＝
ＸＣ ｙ７＝Ｚ Ｇ （０） Ｇ （１） … （ｍ一１） （一１）Ｇｘ（０） … （ｍ一２） ！ ｊ ｊ ｊ （１一ｍ）Ｇｘ（２一ｍ） …Ｇｘ（０）
Ｇｙ（０） Ｇｙ（１） …Ｇｙ（ｎ一１）。 Ｇｙ（一１）Ｇｙ（０） …Ｇｙ ｎ一２） ； ； ； Ｇｙ（１一ｎ） （２一ｎ） …Ｇｙ（０）．
Ｚ＝（２ ）
（５）
易知Ｘ，ｌ，均为非奇异矩阵，即（５）式有唯一解，解方程即 得到系数矩阵Ｃ，将系数矩阵代入式（２）便得到该局部重力 异常基准图二维高斯样条函数逼近解析式。 由上式可以看出系数矩阵Ｃ由Ｘ、Ｙ与ｚ唯一确定，而 、 Ｙ由矩阵的阶数和ａ ，ａ 决定。文献［９］曾对求逆计算误差 与ａ和计算阶数的关系进行过详尽的探讨。经讨论可知：系 数矩阵Ｃ的求解过程中如 ，ｙ矩阵阶数过高且ａ ，。 过小或 阶数过低且ａ ，ａ 过大时， ，ｙ矩阵的求逆存在较大误差。 而当Ｘ，Ｙ矩阵阶数小于５或ａ ，ａ 值小于５时，矩阵求逆运算