CN104714257A - 一种基于逐步插值校正的多重分形克里金插值的局部地磁图构建方法 - Google Patents

Abstract

本发明属于地磁导航技术领域，具体涉及采用基于逐步插值校正的多重分形克里金插值法对稀疏实测的地磁数据进行加密重构的一种基于逐步插值校正的多重分形克里金插值的局部地磁图构建方法。本发明包括：提取实测数据中的地磁异常场数据；计算变异函数表达式及权重系数，利用克里金插值法将实测地磁数据插值为过渡分辨率；选定多重分形估值区域，确定盒子数；确定奇异系数；推算多重分形克里金插值方程；采用交叉验证法验证构图精度并构建地磁基准图；采用逐步校正的方式改善插值结果。本发明利用克里金插值描述地磁异常场信息的趋势，多重分形理论还原原始信息中的局部奇异特征，采用逐步插值校正的方式，构建高精度、高分辨率的局部地磁基准图。

Description

一种基于逐步插值校正的多重分形克里金插值的局部地磁图构建方法

技术领域

本发明属于地磁导航技术领域，具体涉及采用基于逐步插值校正的多重分形克里金插值法对稀疏实测的地磁数据进行加密重构的一种基于逐步插值校正的多重分形克里金插值的局部地磁图构建方法。

背景技术

地磁异常场是由地壳表面分布和局部地理特征如矿产、岩石、人造磁场等产生的，具有长期的稳定性，几乎不随时间的变化而变化，并且地磁异常场在空间上分布频率较丰富，位置变化对应的场强变化明显，反映在地磁基准图上就是细节信息丰富，局部特征明显。可以认为在较长时间内，地磁异常值的变化只与空间地理位置相关，非常适合作为地磁导航的参考场。

李庆谋在《多维分形克里格方法》中将克里金插值法与多重分形理论有机结合，研究了多重分形克里金法在非规则分布的时空信号插值中的应用。克里金法是一种无偏、线性且最优的估计，它的求解过程是建立在空间相关系数矩阵极小意义下的，这就决定了它是一个空间滑动平均或低通滤波过程。由于高频、局部及弱信号在相关系数矩阵的方差中所占比例很小，因此不可避免地造成这些信号被抑制或消除。而多重分形理论可以弥补这一不足，它能够有效地还原原始信号中的局部奇异特征，保留更多的高频、局部及弱信号。常帅在《一种基于多重分形克里金法的局部地磁图构建方法》中证明了多重分形克里金插值法在构建局部地磁基准图方面具有一定的优越性，但插值区域内必须含有充足的测点才能对待插值点的地磁异常场值进行准确估计。在现有数据较稀疏的情况下，用局部估计整体的过程中会导致部分分形特征丢失，影响插值效果。

发明内容

本发明的目的在于提出一种构建高精度、高分辨率的局部地磁基准图的基于逐步插值校正的多重分形克里金插值的局部地磁图构建方法。

本发明的目的是这样实现的：

本发明包括如下步骤：

(1)提取实测数据中的地磁异常场数据：利用地磁场模型计算研究区域内各点的主磁场强度，在实测数据中减去主磁场强度，得到地磁异常场数据；

(2)计算变异函数表达式及权重系数，利用克里金插值法将实测地磁数据插值为过渡分辨率：利用实测地磁异常场数据计算不同间隔距离对应的变异函数值，进而得到间隔距离与变异函数值序列对，利用选定理论模型对其进行拟合运算，得到各已知数据点之间以及已知数据点与待插值点之间的变异函数表达式；在保证插值的估计是无偏的且估计误差的方差最小的前提下，通过求解空间相关系数方程得到各实测数据点的权重系数，完成对克里金插值法中各未知量的求解；选定过渡分辨率，利用克里金插值法将稀疏实测数据插值为过渡分辨率的数据，并将其作为多重分形克里金插值的基准数据；

(3)选定多重分形估值区域，确定盒子数：在将基准数据插值为目标数据的过程中，保证测度不变原则，将克里金插值圆域等效为多重分形插值的正方形邻域，并将这一正方形边长作为待插值点邻域的最大尺度；利用盒计数法，将最大尺度不断缩小奇数倍，得到待插值点邻域的最小尺度，最大尺度与最小尺度的比值即为盒子数；

(4)确定奇异系数：在待插值点邻域内，计算测度与尺度序列对，在测度与尺度的双对数坐标系中，利用最小二乘拟合法得到待插值点的奇异系数；

(5)推算多重分形克里金插值方程：利用待插值点较大尺度与小尺度内的测度表达式及步骤(3)、(4)中计算的盒子数和奇异系数，推导待插值点在其小邻域内的多重分形克里金插值方程；

(6)采用交叉验证法验证构图精度并构建地磁基准图：将多重分形克里金插值结果与相同分辨率的已知数据进行比较，并采用平均估计误差百分比、相对均方差、均方根预测误差、平均绝对离差四个有效性评价指标对构图精度进行评估；

(7)采用逐步校正的方式改善插值结果：若插值结果无法满足要求，调整过渡分辨率的值，返回步骤(2)，重新进行插值。

步骤(1)为：

地磁异常场强度占地磁场总强度的4％-5％，值在2000nT-4000nT之间，主磁场强度值在50000nT-70000nT之间，通过实测地磁场强度值的大小，判断出其中是否含有地磁主磁场，若含有，则利用国际地磁参考场模型IGRF计算地磁主磁场的变化趋势，并在实测总磁场值中将其除去，得到地磁异常场数据。

步骤(2)包括以下几个步骤：

(2.1)计算克里金插值的变异函数表达式及权重系数：

变异函数是利用实测地磁异常场数据和理论模型进行拟合得到的，利用实测数据计算出研究区域的变异函数值：

$\mathrm{\γ}\left(h\right)=\frac{1}{2n}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{[Z_i^{\′}-Z_i]}^2;$

其中h为间隔距离，Z′_i和Z_i是间隔为h的一组采样数据，n为被h分隔的数据对数目，通过调整h值，得到γ(h)～h序列对，拟合为高斯模型：

$\mathrm{\γ}\left(h\right)=C_0+C_1(1-e^{{(-\frac{h}{a})}^2});$

拟合后得到未知参数C₀、C₁、a的值，得到变异函数的计算公式；

权重系数λ_i：

$\left\{\begin{array}{c}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\λ}}_i\·\mathrm{\γ}(x_i,x_j)+\mathrm{\μ}=\mathrm{\γ}(x_0,x_j),1\≤j\≤n\\ \underset{i=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\λ}}_i=1\end{array}\right.,$

式中γ(x_i,x_j)是样本x_i与x_j之间的变异函数值；γ(x₀,x_j)是待插值点x₀与样本x_j之间的变异函数值，μ为拉格朗日乘数，求解空间相关系数方程可得到各待插值点的权重系数λ₁,λ₂,......,λ_n；

(2.2)利用克里金插值法将稀疏实测数据插值为过渡分辨率数据：

利用克里金估值将稀疏实测数据插值为过渡分辨率的数据：

$Z^{*}=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\λ}}_i\·Z_i$

Z^*待插值点x₀处的属性值；Z_i为估计区域内1～n个采样点的实测地磁异常场强度值。

步骤(3)包括：

(3.1)选定多重分形估值区域：

若克里金插值法的估值半径为R，则将圆域等效为多重分形插值的正方形区域时，为保证测度不变，取正方形边长为将这一边长作为多重分形插值中待插值点的最大尺度；

(3.2)确定盒子数：

插值区域内的数据点是按照等步长均匀分布的，按盒计数法原理，以待插值点为中心，将最大尺度N·δ不断缩小奇数倍，直到取到最小尺度δ，N为最大尺度与最小尺度的比值，也称盒子数。

步骤(4)为计算奇异系数：

根据分维定义，拟合最佳逼近时，尺度δ与测度m(s,δ)的关系为：

m(s,δ)＝bδ^a；

式中：a，b为常数。a为奇异系数，表征待插值点小邻域内的凹凸特性；

两边取对数：

lnm(s,δ)＝lnb+alnδ；

δ尺度下的测度m(s,δ)是以s点为中心，δ为边长的正方形邻域的质量，在待插值点的最小尺度到最大尺度邻域内计算测度值，得到测度与尺度序列对，利用最小二乘法拟合得到待插值点的奇异系数:

$m(s,\mathrm{\δ})=\underset{i=s-\frac{\mathrm{\δ}}{2}}{\overset{s+\frac{\mathrm{\δ}}{2}}{\mathrm{\Σ}}}\mathrm{\ρ}(s,\mathrm{\δ})\·(i\·\mathrm{\Δh})=\mathrm{\ρ}(s,\mathrm{\δ})\·{\mathrm{\δ}}^2$

式中ρ(s,δ)正方形邻域的面密度，为δ邻域所包含的已知点的地磁异常场强度均值。

步骤(5)包括：

(1)确定待插值点的最大尺度与最小尺度内测度

尺度为N·δ邻域及尺度为δ邻域内的测度为：

Z^*(N·δ)²＝b(N·δ)^a

Z_δδ²＝bδ^a；

(2)推算多重分形克里金插值方程

重分形克里金插值方程：

$Z_{\mathrm{\δ}}=N^{2-a}\·Z^{*}=N^{2-a}\·\underset{i=1}{\overset{n}{\mathrm{\Σ}}}({\mathrm{\λ}}_i\·Z_i)$

Z^*为N·δ邻域内地磁异常场强度的均值，将基准数据插值为目标分辨率数据的估值结果；Z_δ为δ邻域内地磁异常场强度的均值，也是待插值点的地磁异常场估计值。

本发明的有益效果在于：在对数据进行加密重构的基础上，利用克里金插值描述地磁异常场信息的大体趋势，多重分形理论还原原始信息中的局部奇异特征，并采用逐步插值校正的方式，构建高精度、高分辨率的局部地磁基准图。

附图说明

图1本发明的方法流程图；

图2本发明中区域4053三维比较图。

具体实施方式

下面结合附图对本发明作进一步地描述。

该方法是在稀疏实测的磁异常数据的基础上，利用克里金插值法将数据插值为过渡分辨率，再利用多重分形克里金插值法实现数据的加密重构。采用逐步校正的方式，不断细化过渡分辨率的值，可使插值结果得到明显改善。

本发明提出的基于逐步插值校正的多重分形克里金插值法，包括以下步骤：

步骤一、提取实测数据中的地磁异常场数据

利用地磁场模型计算研究区域内各点的主磁场强度，并在实测数据中将其减去，最终得到地磁异常场数据。

步骤二、计算变异函数表达式及权重系数，利用克里金插值法将实测地磁数据插值为过渡分辨率

利用实测地磁异常场数据计算不同间隔距离对应的变异函数值，进而得到间隔距离与变异函数值序列对，利用选定理论模型对其进行拟合运算，最终得到各已知数据点之间以及已知数据点与待插值点之间的变异函数表达式。在保证插值的估计是无偏的且估计误差的方差最小的前提下，通过求解空间相关系数方程得到各实测数据点的权重系数，完成对克里金插值法中各未知量的求解。

选定过渡分辨率，利用克里金插值法将稀疏实测数据插值为过渡分辨率的数据，并将其作为多重分形克里金插值的基准数据。

步骤三、选定多重分形估值区域，确定盒子数

在将基准数据插值为目标数据的过程中，保证测度不变原则，将克里金插值圆域等效为多重分形插值的正方形邻域，并将这一正方形边长作为待插值点邻域的最大尺度。

利用盒计数法，将最大尺度不断缩小奇数倍，得到待插值点邻域的最小尺度，最大尺度与最小尺度的比值即为盒子数。

步骤四、确定奇异系数

在待插值点邻域内，计算测度与尺度序列对，在测度与尺度的双对数坐标系中，利用最小二乘拟合法得到待插值点的奇异系数。

步骤五、推算多重分形克里金插值方程

利用待插值点较大尺度与小尺度内的测度表达式及步骤三、四中计算的盒子数和奇异系数，推导待插值点在其小邻域内的多重分形克里金插值方程。

步骤六、采用交叉验证法验证构图精度并构建地磁基准图

将多重分形克里金插值结果与相同分辨率的已知数据进行比较，并采用平均估计误差百分比、相对均方差、均方根预测误差、平均绝对离差四个有效性评价指标对构图精度进行评估。

步骤七、采用逐步校正的方式改善插值结果

若插值结果无法满足要求，调整过渡分辨率的值，返回步骤二，重新进行插值。

本发明提出一种采用逐步插值校正的方式将多重分形理论与克里金插值法相结合的地磁异常场数据插值方法，具体包括以下几个步骤：

步骤一、提取实测数据中的地磁异常场数据

地磁异常场强度占地磁场总强度的4％-5％，值在2000nT-4000nT之间，主磁场强度值在50000nT-70000nT之间，通过实测地磁场强度值的大小，可以判断出其中是否含有地磁主磁场，若含有，则利用国际地磁参考场模型(IGRF)计算地磁主磁场的变化趋势，并在实测总磁场值中将其除去，得到地磁异常场数据。

本发明选取编号为4053的研究区域，其数据来源于NOAA发布的航空磁测数据。由于研究区域的边界缺少参考点，以至于边界处的插值精度较低。本发明考虑了边界效应的影响，并通过实验验证了第一、二周参考点对边界处的构图精度有着重要的影响，不可缺少，本发明将研究区域的最外围五周作为边界区域。表1列出了研究区域及中心区域的经纬度范围，中心区域共65×65个待插值点，经纬度分辨率均为0.004°。选取过渡分辨率为0.002°，利用克里金插值法将中心区域插值为过渡分辨率的数据。

表1区域4053的研究区域及中心区域的经纬度范围

步骤二、计算变异函数表达式及权重系数，利用克里金插值法将实测地磁数据插值为过渡分辨率

变异函数是利用实测地磁异常场数据和理论模型进行拟合得到的。根据式(1)可利用实测数据计算出研究区域的变异函数值。

$\mathrm{\γ}\left(h\right)=\frac{1}{2n}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{[Z_i^{\′}-Z_i]}^2---\left(1\right)$

式中：h为间隔距离，Z_i'和Z_i是间隔为h的一组采样数据，n为被h相分隔的数据对数目。通过调整h值，可得到γ(h)～h序列对，利用式(2)将其拟合为高斯模型，

$\mathrm{\γ}\left(h\right)=C_0+C_1(1-e^{{(-\frac{h}{a})}^2})---\left(2\right)$

拟合后可得到未知参数C₀、C₁、a的值，进而得到变异函数的计算公式。

为了保证克里金插值的估计是无偏估计，且估计的方差小于采样值的其他线性组合产生的方差，权重系数λ_i应满足：

$\left\{\begin{array}{c}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\λ}}_i\·\mathrm{\γ}(x_i,x_j)+\mathrm{\μ}=\mathrm{\γ}(x_0,x_j),1\≤j\≤n\\ \underset{i=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\λ}}_i=1\end{array}\right.---\left(3\right)$

式中：γ(x_i,x_j)是样本x_i与x_j之间的变异函数值；γ(x₀,x_j)是待插值点x₀与样本x_j之间的变异函数值，可利用上述变异函数公式计算得到；μ为拉格朗日乘数。求解式(3)，即空间相关系数方程，可得到各待插值点的权重系数λ₁,λ₂,......,λ_n。

利用克里金估值公式(4)，将稀疏实测数据插值为过渡分辨率的数据(基准数据)，以便于增加数据密度，缓解在用局部去估计整体的过程中，部分分形特征丢失的问题。

$Z^{*}=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\λ}}_i\·Z_i---\left(4\right)$

式中：Z^*待插值点x₀处的属性值；Z_i为估计区域内1～n采样点的实测地磁异常场强度值。

步骤三、选定多重分形估值区域，确定盒子数

若克里金插值法的估值半径为R，则将圆域等效为多重分形插值的正方形区域时，为保证测度不变，取正方形边长为将这一边长作为多重分形插值中待插值点的最大尺度。

插值区域内的数据点是按照等步长均匀分布的，依据盒计数法原理，以待插值点为中心，将最大尺度N·δ不断缩小奇数倍，直到取到最小尺度δ。N为最大尺度与最小尺度的比值，也称盒子数，其最佳取值由实验测得，本发明将N值取为2。

步骤四、确定奇异系数

根据分维定义，拟合最佳逼近时，尺度δ与测度m(s,δ)的关系可表示为：

m(s,δ)＝bδ^a    (5)

式中：a，b为常数。a为奇异系数，表征了待插值点小邻域内的凹凸特性。将式(5)

两边取对数：

lnm(s,δ)＝lnb+alnδ    (6)

在双对数坐标系中，式(6)为一条直线，a为该直线的斜率。

对于曲面而言，δ尺度下的测度m(s,δ)是以s点为中心，δ为边长的正方形邻域的质量。在待插值点的最小尺度到最大尺度邻域内利用式(7)计算测度值，得到测度与尺度序列对，利用最小二乘法将其拟合为式(6)形式，最终得到待插值点的奇异系数。

$m(s,\mathrm{\δ})=\underset{i=s-\frac{\mathrm{\δ}}{2}}{\overset{s+\frac{\mathrm{\δ}}{2}}{\mathrm{\Σ}}}\mathrm{\ρ}(s,\mathrm{\δ})\·(i\·\mathrm{\Δh})=\mathrm{\ρ}(s,\mathrm{\δ})\·{\mathrm{\δ}}^2---\left(7\right)$

式中ρ(s,δ)正方形邻域的面密度，可理解为δ邻域所包含的已知点的地磁异常场强度均值。

步骤五、推算多重分形克里金插值方程

根据式(5)、(7)，将尺度为N·δ邻域及尺度为δ邻域内的测度表示为：

Z^*(N·δ)²＝b(N·δ)^a    (8)

Z_δδ²＝bδ^a    (9)

(8)式与(9)式联立得多重分形克里金插值方程：

$Z_{\mathrm{\δ}}=N^{2-a}\·Z^{*}=N^{2-a}\·\underset{i=1}{\overset{n}{\mathrm{\Σ}}}({\mathrm{\λ}}_i\·Z_i)---\left(10\right)$

式中：Z^*为N·δ邻域内地磁异常场强度的均值，也是按照步骤二，将基准数据插值为目标分辨率数据的估值结果，本发明中目标分辨率取为0.001°；Z_δ为δ邻域内地磁异常场强度的均值，也是待插值点的地磁异常场估计值。

步骤六、采用交叉验证法验证构图精度并构建地磁基准图

将多重分形克里金插值结果与相同分辨率的已知数据进行比较，并采用平均估计误差百分比、相对均方差、均方根预测误差、平均绝对离差四个有效性评价指标对构图精度进行评估。

(1)平均估计误差百分比(PAEE)

$\mathrm{PAEE}=\frac{1}{n\·\stackrel{\‾}{Z}}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{[\hat{Z}\left(X_i\right)-Z\left(X_i\right)]}^2\×100\%---\left(11\right)$

式中：n为测点的个数；为所有测点地磁要素平均值；为位置X_i处地磁要素随机变量的估计值；Z(X_i)为位置X_i处地磁要素随机变量的真实值。如果PAEE趋近于0，则说明估计是无偏的，说明所选的插值方法具有一定的可靠程度。

(2)相对均方差(RMSE)

$\mathrm{RMSE}=\frac{1}{n\·s^2}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{[\hat{Z}\left(X_i\right)-Z\left(X_i\right)]}^2---\left(12\right)$

式中：s²为所有测点地磁要素的样本方差。相对均方差反映了随机变量之间的离散程度，可以用来给指标的权重赋值。

(3)均方根预测误差(RMSPE)

$\mathrm{RMSPE}=\sqrt{\frac{1}{n}\underset{i=1}{\overset{1}{\mathrm{\Σ}}}{[\hat{Z}\left(X_i\right)-Z\left(X_i\right)]}^2}---\left(13\right)$

均方根预测误差越小，预测值就越接近它们的真实值，因此可以利用均方根预测误差来选择最优的插值方法。

(4)平均绝对离差(ADD)

它是实测数据与插值估算数据之间的差值的绝对值加权后再取平均值所得，有效地避免了残差计算中正负值抵消的问题。

基于逐步插值校正的多重分形克里金插值法的仿真结果见表2，图2。图2(a)为原始数据中心区域三维图；图2(b)克里金插值构建的地磁基准图；图2(c)基于逐步插值校正的多重分形克里金插值法所构地磁基准图。

表2基于逐步插值校正的多重分形克里金法的构图精度

步骤七、采用逐步校正的方式改善插值结果

若构图精度无法满足要求，调整过渡分辨率的值，返回步骤二，重新进行插值。

Claims (6)

1.一种基于逐步插值校正的多重分形克里金插值的局部地磁图构建方法，其特征在于，包括如下步骤：

(1)提取实测数据中的地磁异常场数据：利用地磁场模型计算研究区域内各点的主磁场强度，在实测数据中减去主磁场强度，得到地磁异常场数据；

(2)计算变异函数表达式及权重系数，利用克里金插值法将实测地磁数据插值为过渡分辨率：利用实测地磁异常场数据计算不同间隔距离对应的变异函数值，进而得到间隔距离与变异函数值序列对，利用选定理论模型对其进行拟合运算，得到各已知数据点之间以及已知数据点与待插值点之间的变异函数表达式；在保证插值的估计是无偏的且估计误差的方差最小的前提下，通过求解空间相关系数方程得到各实测数据点的权重系数，完成对克里金插值法中各未知量的求解；选定过渡分辨率，利用克里金插值法将稀疏实测数据插值为过渡分辨率的数据，并将其作为多重分形克里金插值的基准数据；

(3)选定多重分形估值区域，确定盒子数：在将基准数据插值为目标数据的过程中，保证测度不变原则，将克里金插值圆域等效为多重分形插值的正方形邻域，并将这一正方形边长作为待插值点邻域的最大尺度；利用盒计数法，将最大尺度不断缩小奇数倍，得到待插值点邻域的最小尺度，最大尺度与最小尺度的比值即为盒子数；

(4)确定奇异系数：在待插值点邻域内，计算测度与尺度序列对，在测度与尺度的双对数坐标系中，利用最小二乘拟合法得到待插值点的奇异系数；

(5)推算多重分形克里金插值方程：利用待插值点较大尺度与小尺度内的测度表达式及步骤(3)、(4)中计算的盒子数和奇异系数，推导待插值点在其小邻域内的多重分形克里金插值方程；

(6)采用交叉验证法验证构图精度并构建地磁基准图：将多重分形克里金插值结果与相同分辨率的已知数据进行比较，并采用平均估计误差百分比、相对均方差、均方根预测误差、平均绝对离差四个有效性评价指标对构图精度进行评估；

(7)采用逐步校正的方式改善插值结果：若插值结果无法满足要求，调整过渡分辨率的值，返回步骤(2)，重新进行插值。

2.根据权利要求1所述的一种基于逐步插值校正的多重分形克里金插值的局部地磁图构建方法，其特征在于，所述的步骤(1)为：

地磁异常场强度占地磁场总强度的4％-5％，值在2000nT-4000nT之间，主磁场强度值在50000nT-70000nT之间，通过实测地磁场强度值的大小，判断出其中是否含有地磁主磁场，若含有，则利用国际地磁参考场模型IGRF计算地磁主磁场的变化趋势，并在实测总磁场值中将其除去，得到地磁异常场数据。

3.根据权利要求1所述的一种基于逐步插值校正的多重分形克里金插值的局部地磁图构建方法，其特征在于，所述的步骤(2)包括以下几个步骤：

(2.1)计算克里金插值的变异函数表达式及权重系数：

变异函数是利用实测地磁异常场数据和理论模型进行拟合得到的，利用实测数据计算出研究区域的变异函数值：

$\mathrm{\γ}\left(h\right)=\frac{1}{2n}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{[Z_i^{\′}-Z_i]}^2;$

其中h为间隔距离，Z′_i和Z_i是间隔为h的一组采样数据，n为被h分隔的数据对数目，通过调整h值，得到γ(h)～h序列对，拟合为高斯模型：

$\mathrm{\γ}\left(h\right)=C_0+C_1(1-e^{{(-\frac{h}{a})}^2});$

拟合后得到未知参数C₀、C₁、a的值，得到变异函数的计算公式；

权重系数λ_i：

$\left\{\begin{array}{c}\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\λ}}_i\·\mathrm{\γ}(x_i,x_j)+\mathrm{\μ}=\mathrm{\γ}(x_0,x_j),1\≤j\≤n\\ \underset{i=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\λ}}_i=1\end{array}\right.,$

式中γ(x_i,x_j)是样本x_i与x_j之间的变异函数值；γ(x₀,x_j)是待插值点x₀与样本x_j之间的变异函数值，μ为拉格朗日乘数，求解空间相关系数方程可得到各待插值点的权重系数λ₁,λ₂,......,λ_n；

(2.2)利用克里金插值法将稀疏实测数据插值为过渡分辨率数据：

利用克里金估值将稀疏实测数据插值为过渡分辨率的数据：

$Z^{*}=\underset{i=1}{\overset{n}{\mathrm{\Σ}}}{\mathrm{\λ}}_i\·Z_i$

Z^*待插值点x₀处的属性值；Z_i为估计区域内1～n个采样点的实测地磁异常场强度值。

4.根据权利要求1所述的一种基于逐步插值校正的多重分形克里金插值的局部地磁图构建方法，其特征在于，所述的步骤(3)包括：

(3.1)选定多重分形估值区域：

若克里金插值法的估值半径为R，则将圆域等效为多重分形插值的正方形区域时，为保证测度不变，取正方形边长为将这一边长作为多重分形插值中待插值点的最大尺度；

(3.2)确定盒子数：

插值区域内的数据点是按照等步长均匀分布的，按盒计数法原理，以待插值点为中心，将最大尺度N·δ不断缩小奇数倍，直到取到最小尺度δ，N为最大尺度与最小尺度的比值，也称盒子数。

5.根据权利要求1所述的一种基于逐步插值校正的多重分形克里金插值的局部地磁图构建方法，其特征在于，所述的步骤(4)为计算奇异系数：

根据分维定义，拟合最佳逼近时，尺度δ与测度m(s,δ)的关系为：

m(s,δ)＝bδ^a；

式中：a，b为常数，a为奇异系数，表征待插值点小邻域内的凹凸特性；

两边取对数：

lnm(s,δ)＝lnb+alnδ；

δ尺度下的测度m(s,δ)是以s点为中心，δ为边长的正方形邻域的质量，在待插值点的最小尺度到最大尺度邻域内计算测度值，得到测度与尺度序列对，利用最小二乘法拟合得到待插值点的奇异系数:

$m(s,\mathrm{\δ})=\underset{i=s-\frac{\mathrm{\δ}}{2}}{\overset{s+\frac{\mathrm{\δ}}{2}}{\mathrm{\Σ}}}\mathrm{\ρ}(s,\mathrm{\δ})\·(i\·\mathrm{\Δh})=\mathrm{\ρ}(s,\mathrm{\δ})\·{\mathrm{\δ}}^2$

式中ρ(s,δ)正方形邻域的面密度，为δ邻域所包含的已知点的地磁异常场强度均值。

6.根据权利要求1所述的一种基于逐步插值校正的多重分形克里金插值的局部地磁图构建方法，其特征在于，所述的步骤(5)包括：

(1)确定待插值点的最大尺度与最小尺度内测度

尺度为N·δ邻域及尺度为δ邻域内的测度为：

Z^*(N·δ)²＝b(N·δ)^a

Z_δδ²＝bδ^a；

(2)推算多重分形克里金插值方程

重分形克里金插值方程：

$Z_{\mathrm{\δ}}=N^{2-a}\·Z^{*}=N^{2-a}\·\underset{i=1}{\overset{n}{\mathrm{\Σ}}}({\mathrm{\λ}}_i\·Z_i)$

Z^*为N·δ邻域内地磁异常场强度的均值，将基准数据插值为目标分辨率数据的估值结果；Z_δ为δ邻域内地磁异常场强度的均值，也是待插值点的地磁异常场估计值。