基于泛克里金法的重力场插值方法
【技术领域】
[0001] 本发明属于惯性/重力组合导航系统技术领域，具体涉及一种基于泛克里金法的 重力场插值方法。
【背景技术】
[0002] 导航用背景图是重力匹配导航的基础，背景图的精度及分辨率直接影响匹配精 度。目前大部分海域的重力异常信息，都采用卫星反演的手段获得。但该方法的分辨率较 低，最高只有Γ*Γ。在匹配过程中，为了保证重力异常的测量精度，载体必须以低速航行， 从实用性的角度考虑，匹配时间不宜过长，在长度有限的航迹上，卫星反演背景图包含的信 息量太少，无法满足匹配定位需求，因此必须对其进行插值处理。
[0003] 目前常用的背景图插值算法包括距离倒数插值法及三次样条插值法。距离倒数插 值法首先搜索最靠近待插值点的四个背景图采样点（东北、西北、东南、西南方向），之后计 算待插值点到四个采样点的距离，以距离的倒数作为权重，取四个采样点重力异常值的加 权和，作为待插值点的重力异常。三次样条插值法对采样点之间的空间区域，采用连续、二 阶可导，且导数连续的三次函数拟合重力异常分布，其插值结果更加光滑。这两种算法计算 量较小，但均没有针对重力场的特点进行针对性的建模描述。
[0004] 由于地质结构千变万化，重力异常随之剧烈起伏，这会导致在不同区域内，重力异 常的均值和方差并非常数，也即重力场是一个非平稳场，上述两种算法人为设定非平稳的 重力场在空间上按照特定的模型分布，因此其插值结果不可避免的存在偏差。此外，地质结 构在空间上是存在连续性的，这种连续性表现在重力异常值上，就是重力异常在空间上存 在相关性，而这种相关性在不同空间范围内，又是各异的，上述两种算法均没有针对重力场 的相关性进行建模描述，因此无法准确表现出地质结构的连续性。最后，上述两种插值算 法，其插值的权重系数仅由间隔距离计算得到，亦不是最优的。综上所述，上述两种插值算 法，既无法达到无偏最优的插值目标，又没有描述重力场的相关性，因此插值精度较低。
【发明内容】
[0005] 本发明的目的在于克服现有技术的不足，提供一种设计合理且插值精度高的基于 泛克里金法的重力场插值方法。
[0006] 本发明解决其技术问题是采取以下技术方案实现的：
[0007] -种基于泛克里金法的重力场插值方法，于包括以下步骤：
[0008] 步骤1、针对空间中的重力异常模型建立回归模型；
[0009] 步骤2、针对空间中的重力异常模型建立相关模型；
[0010] 步骤3、计算回归模型和相关模型的权重系数，代入重力异常模型，便可获得待插 值位置上的重力异常空间分布表达式，再将待插值位置代入重力异常模型，便可计算出待 插值位置的重力值。
[0011] 而且，所述重力异常模型为：
[0012] yx = f Ix…fpx* β i，1+zx
[0013] 式中，x为空间位置，是一个1*2维位置矩阵，第一列保存经度，第二列保存炜度； yx为模型解算出的重力值。
[0014] 而且，所述步骤1建立回归模型的方法包括以下步骤：
[0015] (1)建立由全部因子xl, x2,…，xn组成的自变量集合，设自变量集合包含η个因 子，并设定因子被选入回归模型的显著性水平为α 1，被剔除出回归模型的显著性水平为 α 2,0〈α 1 彡 α 2〈1 ;
[0016] (2)在背景图上抽取m组经炜度与重力值数据，作为回归模型的样本；
[0017] (3)从全部因子中，抽取出任意两个因子xi，xj，依据m对样本值，利用最小二乘法 计算其参数向量τ i, jl，得到用这两个因子描述的空间重力场ykl = xi, xj* τ i, jl，并计 算这个模型的剩余平方和Si, jl = k = lmyk-ykl2 ;
[0018] (4)重复执行步骤（3)，遍历全部因子，计算任意两个因子所组成的模型的剩余平 方和Si, jl，组成离差矩阵Sl ;
[0019] (5)计算各个因子在离差矩阵中的贡献Vjl= (Sj，yl)2Sj，jl，设贡献最大的因子 编号为kl，对其做显著性水平为a 1的显著性F检验：
[0020] F = VjlSEl(M-I-I)
[0021] 其中，SEl = k = lmyk-y2_Vjl，y = k = Imykm ;
[0022] 当F>Fa ll，m-1-l时，该因子选入回归模型；
[0023] (6)在kl因子选入回归模型后，调整离差矩阵S2 :
[0024] 按照步骤（5)的方法重新计算剩余因子在新离差矩阵中的贡献，并对贡献最大的 因子k2做显著性水平为a 1，自由度为m-2-l的显著性F检验，作为其是否被引入回归模型 的判据；
[0025] (7)在选入两个因子后，首先按照步骤（5)的方法，将离差矩阵S2调整为S3， 在选入新因子之前，先对已经选入的两个因子，按照被剔除出回归模型的显著性水平 a 2做显著性F检验，以判断是否需要将其从回归模型中剔除，F>F a 21，m-3-l时保留， F〈Fa 21，m-3-l时剔除，之后再引入新的因子；
[0026] (8)按照步骤（7)筛选的方法，对自变量集合中的因子进行遍历，直到新的因子不 再被引入回归模型，同时回归模型内的因子也不再被剔除时，即可获得重力场的回归模型。
[0027] 而且，所述步骤2建立相关模型的方法包括以下步骤：
[0028] (1)对重力矩阵Y进行平差处理，得到Y' ；
[0029] (2)计算Y'的相关矩阵，对于空间中任意两个采样点位置Xi与xj，Y'的相关矩 阵为：
[0030] EY，xi*Y，xj = Ri，j = σ 2R Θ，xi，xj
[0031] 上式中，σ 2为背景图重力值的方差；
[0032] (3)使用高斯相关模型对Z进行描述，高斯相关模型的相关函数形式如下：
[0033] CO = σ 2 = Cl
[0034] Ch = Cle_h Θ 2
[0035] Chh -① 0
[0036] (4)计算相关尺度Θ的最优解：
[0037] min θ = c Θ E Rlm 〇 2
[0038] 其中，R为矩阵Ri，j所对应的行列式值，σ 2为背景图重力值的方差；通过调整相 关尺度Θ，结合高斯模型，来逼近Z的相关矩阵Ri，j，构造和描述相关模型。
[0039] 本发明的优点和积极效果是：
[0040] 1、本发明采用逐步回归法对重力场的非平稳性及相关性分别建模描述，使用多元 逐步回归法，结合拉格朗日乘数法，完整的描述了重力场的非平稳性，保证了插值结果的无 偏性，具有较高的插值精度。
[0041] 2、本发明拉格朗日乘数法，保证了插值结果的最优性。
[0042] 3、本发明采用回归模型的引入方式，使得构建的重力异常模型完整表述重力场的 相关性。
【附图说明】
[0043] 图1为逐步回归法的计算步骤；
[0044] 图2为某重力异常序列的相关部分；
[0045] 图3为重力异常相关部分的相关函数；
[0046] 图4为高斯模型协方差函数；
[0047] 图5为4#测线插值结果；
[0048] 图6为4#测线采用三次样条插值误差变异函数；
[0049] 图7为4#测线采用本专利算法插值误差变异函数；
【具体实施方式】
[0050] 以下结合附图对本发明实施例做进一步详述：
[0051] -种基于泛克里金法的重力场插值方法，是基于泛克里金法来设计的。泛克里金 法是一种针对非平稳场的插值算法，该算法是建立在变异函数理论分析基础上，对有限区 域内变量取值进行无偏最优估计的一种方法。本发明的算法假定空间中的重力异常符合如 下模型：
[0052] yx = fix…fpx* β i，1+zx (I)
[0053] 式（I)中，x为空间位置，是一个1*2维位置矩阵，第一列保存经度，第二列保存炜 度；yx为模型解算出的重力值。
[0054] 将m个背景图采样点位置xl…xm及与之相对应的m个重力异常值yxl"^yxm代入 式（1)，可以得到：
[0062] 于是，式⑵可简写为：Y = F* β +Z (3)
[0063] 一种基于泛克里金法的重力场插值方法，包括以下步骤：
[0064] 步骤1、针对空间中的重力异常模型建立回归模型
[0065]由于重力场是非平稳场，在不同的区域内，重力异常的均值并非常数。针对这个特 点，本算法采用以经度、炜度为自变量，重力异常值为因变量的多项式，来拟合回归模型，并 利用逐步回归法，确定组成回归模型的多项式因子（下面简称为因子）及其阶数。
[0066] 本步骤是对公式（1)中的fix…fpx建模，