CN112050799B - 一种基于磁梯度张量缩并量之比的测距定位方法 - Google Patents

Abstract

本发明公开了一种基于磁梯度张量缩并量之比的测距定位方法，属于磁探测与定位技术领域。所述定位方法包括以下步骤：将十个同类型三轴磁强计构成三轴磁强计平面阵列并将它们的各个敏感轴相互对齐；采集平面阵列的十个三轴磁强计输出；计算阵列四个面心处的磁性体磁梯度张量的独立分量值；由每两面心与磁性体之间的距离值的平方比为固定数值构成关于磁性体位置坐标的线性方程组，计算方程组得出磁性体位置坐标值。本发明相比于三轴磁强计阵列的空间立体结构，平面结构易于多个三轴磁强计组阵，特别是对于芯片化的微型三轴磁强计更具优势，也减少了阵列的安放空间。

Description

一种基于磁梯度张量缩并量之比的测距定位方法

技术领域

本发明涉及一种基于磁梯度张量缩并量之比的测距定位方法，属于磁探测与定位技术领域。

背景技术

磁梯度张量是磁场分量的空间变化率。磁梯度张量数据在磁性体的不同方向上呈现出的信息量极为丰富，有利于磁性体的识别和定位，是磁性体探测与定位的一种无源方式，具有很好的隐蔽性与广阔的应用前景，例如隐伏地质体勘探、UXO探测定位、反潜追踪等。近年来，磁梯度张量场的数据解释理论与应用得到人们的重视。

采用单一标量磁探测仪进行磁定位通常需要移动探测仪进行多次测量，定位结果还存在多解问题。目前，空间立体构型的三轴磁强计阵列可以实现单次测量定位，但这种阵列占用的空间较大。本发明提出的三轴磁强计阵列的平面结构易于多个三轴磁强计组阵，特别是对于芯片化的微型三轴磁强计更具优势，也减少了阵列的安放空间；提出的基于磁梯度张量缩并量之比的测距定位算法能通过单次测量获得磁性体位置坐标的唯一解，无需事先知道磁性体的磁矩，也消除了磁导率带来的定位误差。磁性体磁矩的最小二乘反演法抑制了磁场梯度张量测量噪声的干扰。

发明内容

本发明的目的是提出一种基于磁梯度张量缩并量之比的测距定位方法，已解决现有的磁探测方法需事先知道磁性体的磁矩、磁导率会带来定位误差等问题。

一种基于磁梯度张量缩并量之比的测距定位方法，所述方法包括以下步骤：

步骤一、以三轴磁强计平面阵列所在的平面为z＝0建立空间直角坐标系，在面z＝0上设置四个磁梯度张量测量点，它们的空间坐标分别为P₁(x₁,y₁,0)、P₂(x₂,y₂,0)、P₃(x₃,y₃,0)和P₄(x₄,y₄,0)；

步骤二、以所述四个磁梯度张量测量点为局部十字形平面阵列的交叉点，由多个局部十字形阵列构成三轴磁强计平面阵列；

步骤三、将磁梯度张量缩并量C_T的定义式为

式中，T_αβ为磁梯度张量G的第α行、第β列元素，

磁梯度张量测量点相对于磁性体的位置矢量与磁性体磁矩矢量之间夹角φ的余弦值为：

其中，λ₁、λ₂、λ₃分别为G的特征值λ的第一个值、第二个值和第三个值，

对于磁偶极子磁性体，其磁梯度张量的缩并量C_T的计算式为：

式中，系数

为磁梯度张量测量点相对于磁性体的位置矢量与磁性体磁矩矢量之间的夹角，μ为介质磁导率，m为磁性体磁矩的大小，r为磁梯度张量测量点与磁性体之间的距离；

步骤四、利用三轴磁强计平面阵列测量阵列所在平面四个点处的磁梯度张量值，由式(1)和式(2)分别计算得到平面四个点处的磁梯度张量缩并量及夹角余弦值，由于磁性体磁矩的大小是不变的，根据式(3)，由不同点处的磁梯度张量缩并量之比消除未知的磁性体磁矩的大小，消除介质磁导率的影响，因此，得到不同测量点距磁性体距离之比的平方ξ_ij为：

式中，j,k＝1,2,3,4且j≠k；

步骤五、式(4)包含有三个独立方程，所述三个独立方程构成关于磁性体位置坐标的方程组，如式(5)所示：

式中，x_Q、y_Q和z_Q分别为磁性体位置坐标，式(5)中的四个磁梯度张量测量点空间坐标是已知的，通过求解式(5)得到磁性体的位置坐标值，

磁性体的磁梯度张量只有五个独立分量，选为T₁₁、T₂₂、T₁₂、T₂₃和T₁₃，则有

式中，x_R、y_R和z_R分别为测量点相对于磁性体的位置坐标三分量值，m_x、m_y和m_z分别为磁性体磁矩的三分量值，

步骤六、将磁性体的位置坐标三分量计算值和不同测量点处的磁场梯度张量五个独立分量代入式(6)，得到关于m_x、m_y和m_z的超定线性方程组，利用最小二乘法计算磁性体磁矩的三分量值。

本发明的主要优点是：

(1)现有的十字形三轴磁强计阵列只能测量阵列中心处的磁梯度张量，本发明提出一种三轴磁强计平面阵列是同时测量同一平面内多点处的磁梯度张量的阵列结构，易于多个三轴磁强计组阵，特别是对于芯片化的微型三轴磁强计更具优势，也减少了阵列的安放空间。

(2)根据磁梯度张量缩并量与距离之间的关系，再求出不同测量点间的距离之比，从而求解位置定位，本发明所提的磁性体定位算法无需事先知道磁性体的磁矩，也消除了磁导率带来的定位误差。

(3)由不同测量点处的磁场梯度张量五个独立分量构成超定线性方程组，利用最小二乘法计算磁性体磁矩的三分量值，减少了磁场梯度张量测量噪声对解算磁性体磁矩的影响。

附图说明

图1是三轴磁强计平面阵列的结构示意图；

图2是基于三轴磁强计平面阵列的磁梯度张量定位算法流程图；

图3是磁性体位置坐标三分量的平均绝对误差值随三轴磁强计噪声标准差的变化曲线，其中，图3(a)为位置坐标x_Q的平均绝对误差值；图3(b)为位置坐标y_Q的平均绝对误差值；图3(c)为位置坐标z_Q的平均绝对误差值；

图4是磁性体磁矩三分量的平均相对误差值随三轴磁强计噪声标准差的变化曲线，其中，图4(a)为磁矩分量m_x的平均相对误差值；图4(b)为磁矩分量m_y的平均相对误差值；图4(c)为磁矩分量m_z的平均相对误差值。

具体实施方式

下面将结合本发明实施例中的附图对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

参照图2所示，一种基于磁梯度张量缩并量之比的测距定位方法，所述方法包括以下步骤：

步骤一、以三轴磁强计平面阵列所在的平面为z＝0建立空间直角坐标系，在面z＝0上设置四个磁梯度张量测量点，它们的空间坐标分别为P₁(x₁,y₁,0)、P₂(x₂,y₂,0)、P₃(x₃,y₃,0)和P₄(x₄,y₄,0)；

步骤二、以所述四个磁梯度张量测量点为局部十字形平面阵列的交叉点，由多个局部十字形阵列构成三轴磁强计平面阵列；

步骤三、将磁梯度张量缩并量C_T的定义式为

式中，T_αβ为磁梯度张量G的第α行、第β列元素，

磁梯度张量测量点相对于磁性体的位置矢量与磁性体磁矩矢量之间夹角φ的余弦值为：

其中，λ₁、λ₂、λ₃分别为G的特征值λ的第一个值、第二个值和第三个值，

对于磁偶极子磁性体，其磁梯度张量的缩并量C_T的计算式为：

式中，系数

为磁梯度张量测量点相对于磁性体的位置矢量与磁性体磁矩矢量之间的夹角，μ为介质磁导率，m为磁性体磁矩的大小，r为磁梯度张量测量点与磁性体之间的距离；

步骤四、利用三轴磁强计平面阵列测量阵列所在平面四个点处的磁梯度张量值，由式(1)和式(2)分别计算得到平面四个点处的磁梯度张量缩并量及夹角余弦值，由于磁性体磁矩的大小是不变的，根据式(3)，由不同点处的磁梯度张量缩并量之比消除未知的磁性体磁矩的大小，消除介质磁导率的影响，因此，得到不同测量点距磁性体距离之比的平方ξ_ij为：

式中，j,k＝1,2,3,4且j≠k；

步骤五、式(4)包含有三个独立方程，所述三个独立方程构成关于磁性体位置坐标的方程组，如式(5)所示：

式中，x_Q、y_Q和z_Q分别为磁性体位置坐标，式(5)中的四个磁梯度张量测量点空间坐标是已知的，通过求解式(5)得到磁性体的位置坐标值，

磁性体的磁梯度张量只有五个独立分量，选为T₁₁、T₂₂、T₁₂、T₂₃和T₁₃，则有

式中，x_R、y_R和z_R分别为测量点相对于磁性体的位置坐标三分量值，m_x、m_y和m_z分别为磁性体磁矩的三分量值，

步骤六、将磁性体的位置坐标三分量计算值和不同测量点处的磁场梯度张量五个独立分量代入式(6)，得到关于m_x、m_y和m_z的超定线性方程组，利用最小二乘法计算磁性体磁矩的三分量值。

下面提出一具体实施例：

步骤1、如图1所示，一种三轴磁强计平面阵列由十个三轴磁强计组成，这十个三轴磁强计分别位于点M₁、M₂、…、M₁₀，三轴磁强计的编号随之从1到10。将十个三轴磁强计对应的敏感轴相互对齐，并按图1所示的配置方式放置在安装平面上，并在这个安装面上固定测量坐标系Oxyz，安装面的方程为z＝0。点M_i(i＝1,2,…,10)的平面位置坐标为(x_Mi,y_Mi)。点M₁、M₂、M₅、M₆所构成的矩形的中心点为P₁，其平面位置坐标为(x₁,y₁)；点M₂、M₃、M₆、M₇所构成的矩形的中心点为P₂，其平面位置坐标为(x₂,y₂)；点M₃、M₄、M₇、M₈所构成的矩形的中心点为P₃，其平面位置坐标为(x₃,y₃)；磁力计M₇、M₈、M₉，M₁₀所构成的矩形的中心点为P₄，其平面位置坐标为(x₄,y₄)。

步骤2、由多通道数据采集器同步采集十个三轴磁强计输出的磁场分量测量值，再按式(7)、(8)、(9)和(10)分别计算点P₁、P₂、P₃和P₄处的磁梯度张量矩阵

和

式中，

和

为第i(i＝1,2,…,10)个三轴磁强计输出的磁场分量测量值。

步骤3、由点P_j(j＝1,2,3,4)处的磁梯度张量矩阵

按式(11)计算几何不变量

式中，

为磁梯度张量矩阵

的第α行、第β列元素。

步骤4、计算

的特征值

根据式(12)计算P_j点相对于磁性体的位置矢量与磁性体磁矩矢量之间夹角的余弦值，

式中，

分别为

的第一个值、第二个值和第三个值。

再根据式(13)计算系数γ_j。

步骤5、根据式(14)计算不同测量点与磁性体之间距离之比的平方值ξ_jk，

式中，r_j为点P_j与磁性体之间的距离，j,k＝1,2,3,4且j≠k，ξ_jk存在三个独立的值。

步骤6、由ξ_jk的三个独立值和式(14)，得到如式(15)所示的关于磁性体位置坐标x_Q、y_Q和z_Q的三个独立方程，

当ξ₁₂、ξ₂₃和ξ₃₄存在等于1的情况时，式(15)的三个方程中有线性方程，此时将此方程代入其他方程，就可较容易地求出磁性体位置坐标，因此不再赘述此时的磁性体位置坐标计算公式。

而通常情况下，ξ₁₂≠1、ξ₂₃≠1且ξ₃₄≠1，由于z₁＝z₂＝z₃＝z₄＝0，则有

由式(16)中的第一个方程和第二方程分别减去第三个方程，以消除坐标变量z_Q，得到关于坐标变量x_Q和y_Q的方程组为

其中

步骤7、由式(18)得到磁性体位置坐标的计算值

和

为

再将式(18)代入式(19)计算得到磁性体位置坐标的计算值

式中，±是根据三轴磁强计阵列是位于所探测定位的磁性体的上方或下方以及测量坐标轴的取向来确定。

步骤8、由磁梯度张量矩阵的测量值

和磁性体位置坐标的计算值

和

按式(21)计算磁矩矢量，得到磁矩矢量的计算值

式中，μ为磁导率，

和

和

分别为矩阵

的五个独立分量。

步骤9、由P_j处的

计算出的磁矩矢量

按式(21)计算磁性体的磁矩矢量，得到磁性体的磁矩矢量计算值

这样可通过磁性体的磁矩矢量值判断磁性体的状态和类型。

为了表征三轴磁强计平面阵列反演磁性体位置及磁矩的精度，定义磁性体定位的绝对误差分别为

和

其中x_Q、y_Q和z_Q为磁性体位置坐标的真值；定义磁性体磁矩三分量反演的相对误差分别为

和

m_x、m_y和m_z为磁性体磁矩分量的真值。

在Matlab数值仿真平台下，对该方法进行仿真实验。位于点Q处磁性体的磁矩三分量分别为m_x＝9×10⁵A·m²、m_y＝8×10⁶A·m²和m_z＝7×10⁶A·m²，点Q的坐标为(60m，50m，10m)。采用图1所示的三轴磁强计平面阵列测量点P1、P2、P3和P4处的磁梯度张量值，阵列平面方程为z＝0。第1个三轴磁强计在阵列平面中的位置坐标为x_M1＝-0.1m和y_M1＝0.2m，第2个三轴磁强计在阵列平面中的位置坐标为x_M2＝0.1m和y_M2＝0.4m，第3个三轴磁强计在阵列平面中的位置坐标为x_M3＝0.3m和y_M3＝0.6m，第4个三轴磁强计在阵列平面中的位置坐标为x_M4＝0.5m和y_M4＝0.8m，第5个三轴磁强计在阵列平面中的位置坐标为x_M5＝0.1m和y_M5＝0m，第6个三轴磁强计在阵列平面中的位置坐标为x_M6＝0.3m和y_M6＝0.2m，第7个三轴磁强计在阵列平面中的位置坐标为x_M7＝0.5m和y_M7＝0.4m，第8个三轴磁强计在阵列平面中的位置坐标为x_M8＝0.7m和y_M8＝0.6m，第9个三轴磁强计在阵列平面中的位置坐标为x_M9＝0.7m和y_M9＝0.2m，第10个三轴磁强计在阵列平面中的位置坐标为x_M10＝0.9m和y_M10＝0.4m。这样，点P1的位置坐标为x₁＝0.1m和y₁＝0.2m，点P2的位置坐标为x₂＝0.3m和y₂＝0.4m，点P3的位置坐标为x₃＝0.5m和y₃＝0.6m，点P4的位置坐标为x₄＝0.7m和y₄＝0.4m。

三轴磁强计在每个轴的噪声是相互独立的高斯过程，其均值为0，标准差为σ_magnetometer，在不同的三轴磁强计噪声标准差的情况下，基于三轴磁强计平面阵列的磁张量缩并量比的磁性体定位算法进行50次蒙特卡洛仿真实验。当σ_magnetometer从0nT增加到0.45nT时，反演出的磁性体位置坐标三分量的平均绝对误差值随三轴磁强计噪声标准差的变化曲线如图3所示，由图3可知，三个位置坐标分量反演的平均绝对误差值随σ_magnetometer的增加而增大。反演出的磁性体磁矩三分量的平均绝对误差值随三轴磁强计噪声标准差的变化曲线如图4所示，由图4可知，磁性体磁矩三分量的平均相对误差值随σ_magnetometer的增加而增大。

相比于现有的三轴磁强计阵列及其磁梯度张量定位方法，本发明提出的一种三轴磁强计平面阵列可使用十只三轴磁强计同时测量平面内四个点的磁梯度张量，相比于三轴磁强计阵列的空间立体结构，平面结构易于多个三轴磁强计组阵，特别是对于芯片化的微型三轴磁强计更具优势，也减少了阵列的安放空间。提出的基于磁梯度张量缩并量之比的测距定位算法是根据磁梯度张量缩并量与距离之间的关系，求出不同测量点间的距离之比，由三个独立的距离比值解算磁性体位置，这种磁梯度张量定位算法无需事先知道磁性体的磁矩，也消除了磁导率带来的定位误差。在磁性体磁矩的反演计算中，利用最小二乘法求解超定线性方程组获得磁性体磁矩的三分量值，这样就减少了磁场梯度张量测量噪声对解算磁性体磁矩的影响。

Claims (1)

1.一种基于磁梯度张量缩并量之比的测距定位方法，其特征在于，所述方法包括以下步骤：

步骤一、以三轴磁强计平面阵列所在的平面为z＝0建立空间直角坐标系，三轴磁强计平面阵列由十个三轴磁强计组成，这十个三轴磁强计分别位于点M₁、M₂、…、M₁₀，三轴磁强计的编号随之从1到10，将十个三轴磁强计对应的敏感轴相互对齐，并按L型阵列构型的配置方式放置在安装平面上，并在所述安装平面上固定测量坐标系Oxyz，安装面的方程为z＝0，点M_i(i＝1,2,…,10)的平面位置坐标为(x_Mi,y_Mi)，点M₁、M₂、M₅、M₆所构成的矩形的中心点为P₁，其平面位置坐标为(x₁,y₁)；点M₂、M₃、M₆、M₇所构成的矩形的中心点为P₂，其平面位置坐标为(x₂,y₂)；点M₃、M₄、M₇、M₈所构成的矩形的中心点为P₃，其平面位置坐标为(x₃,y₃)；磁力计M₇、M₈、M₉，M₁₀所构成的矩形的中心点为P₄，其平面位置坐标为(x₄,y₄)，在面z＝0上设置四个磁梯度张量测量点，所述四个磁梯度张量测量点的空间坐标分别为P₁(x₁,y₁,0)、P₂(x₂,y₂,0)、P₃(x₃,y₃,0)和P₄(x₄,y₄,0)；

步骤二、以所述四个磁梯度张量测量点为局部十字形平面阵列的交叉点，由多个局部十字形阵列构成三轴磁强计平面阵列；

步骤三、将磁梯度张量缩并量C_T的定义式为：

式中，T_αβ为磁梯度张量G的第α行、第β列元素，

磁梯度张量测量点相对于磁性体的位置矢量与磁性体磁矩矢量之间夹角φ的余弦值为：

其中，λ₁、λ₂、λ₃分别为G的特征值λ的第一个值、第二个值和第三个值，

对于磁偶极子磁性体，其磁梯度张量的缩并量C_T的计算式为：

式中，系数

为磁梯度张量测量点相对于磁性体的位置矢量与磁性体磁矩矢量之间的夹角，μ为介质磁导率，m为磁性体磁矩的大小，r为磁梯度张量测量点与磁性体之间的距离；

步骤四、利用三轴磁强计平面阵列测量阵列所在平面四个点处的磁梯度张量值，由式(1)和式(2)分别计算得到平面四个点处的磁梯度张量缩并量及夹角余弦值，由于磁性体磁矩的大小是不变的，根据式(3)，由不同点处的磁梯度张量缩并量之比消除未知的磁性体磁矩的大小，消除介质磁导率的影响，因此，得到不同测量点距磁性体距离之比的平方ξ_ij为：

式中，j,k＝1,2,3,4且j≠k；

步骤五、式(4)包含有三个独立方程，所述三个独立方程构成关于磁性体位置坐标的方程组，如式(5)所示：

式中，x_Q、y_Q和z_Q分别为磁性体位置坐标，式(5)中的四个磁梯度张量测量点空间坐标是已知的，通过求解式(5)得到磁性体的位置坐标值，

磁性体的磁梯度张量只有五个独立分量，选为T₁₁、T₂₂、T₁₂、T₂₃和T₁₃，则有：

式中，x_R、y_R和z_R分别为测量点相对于磁性体的位置坐标三分量值，m_x、m_y和m_z分别为磁性体磁矩的三分量值，

步骤六、将磁性体的位置坐标三分量计算值和不同测量点处的磁场梯度张量五个独立分量代入式(6)，得到关于m_x、m_y和m_z的超定线性方程组，利用最小二乘法计算磁性体磁矩的三分量值。

Images