基于重磁梯度数据张量不变量的目标深度估计方法及系统与流程

本发明属于重力勘探领域，尤其涉及基于重磁梯度数据张量不变量的目标深度估计方法及系统。

背景技术：

本部分的陈述仅仅是提供了与本公开相关的背景技术信息，不必然构成在先技术。

重力勘探起初是用来研究地球形状的重力学测量方法，是人类最早于用来研究地球的手段，如今重力勘探仍然是地球物理勘探中一个重要的分支，其目的是测量地下由密度分布不均和差异引起的重力或重力梯度异常，进而用来研究和推断地下地质体空间和密度分布属性等，更好地认识地下空间并服务于人类社会。重力勘探本身具有快速有效、轻便和成本低的优点，被广泛地应用于研究大地构造、勘探矿产和油气源分布，以及工程与环境调查领域。

发明人在研究中发现，地质目标体的深度位置计算是重磁数据解释的重点问题，而位场数据具有较低的深度分辨率特点且易受磁倾角和磁偏角的影响，深度估计误差较大，是一个难点问题；而若采用二维剖面数据进行处理解释时，其具有局限性且易受噪声干扰的影响；而通常采用重力异常或垂向重力异常数据处理解释，其信息量少，数据利用率较低。

技术实现要素：

为克服上述现有技术的不足，本发明提供了基于重磁梯度数据张量不变量的目标深度估计方法，对于单独地质体或距离较远的多个地质体情况下的深度估计问题，可以获得较高精度的计算结果。

为实现上述目的，本发明的一个或多个实施例提供了如下技术方案：

基于重磁梯度数据张量不变量的目标深度估计方法，包括：

构建引力位在球体模型三个方向上的二阶梯度数据平面图，根据重力梯度张量数据，计算其张量矩阵的特征值及特征值组合的张量不变量；

针对球体模型，根据特征值组合后的张量不变量与总水平导数，计算张量不变量与总水平导数比值的反正切值，建立水平距离与目标深度之间的关系；

根据反正切值平面图上的两条45°等值线之间的距离，对目标进行深度估计。

进一步的技术方案，重力梯度张量数据具体为：重力场是引力位在笛卡尔坐标系x,y,z方向上的矢量分量，重力梯度张量为地质体的引力位的二阶导数，即重力三分量gx,gy,gz在x,y,z三个方向上的导数，其中gx在x，y，z方向的一阶梯度分别为gxx,gxy,gxz，gy在x，y，z方向的一阶梯度分别为gyx,gyy,gyz，gz在x，y，z方向的一阶梯度分别为gzx,gzy,gzz，可构成梯度张量矩阵：

其中矩阵t为对称矩阵，且重力位满足拉普拉斯方程，即有gxy＝gyx，gxz＝gzx，gyz＝gzy，gxx+gyy+gzz＝0，即矩阵t有5个独立的分量。

进一步的技术方案，由于重力梯度张量矩阵t为一个实对称阵，通过计算特征值将其对称化，即张量矩阵t可以通过坐标转换成对称阵λ：

e^tte＝λ

其特征值λ和特征向量e分别为：

e＝[e1e2e3]

其中，特征向量e表示新笛卡尔坐标系的正交基底，在任何坐标转换下，矩阵t将包含三个不变量要素：

i0＝trace(t)＝gxx+gyy+gzz＝0

可以通过三个不变量要素，得到求解特征值的参数方程：

λ³-i0λ²+i1λ-i2＝0

求其方程的根，得到三个特征值λ1,λ2,λ3为

λ1＝a+b(1)

其中，

即三个特征值可利用重力梯度张量矩阵中的各个分量表示。

进一步的技术方案，计算三个特征值组合后的不变量μ，即

其中不变量μ可利用重力梯度张量的各个分量表示，且为不随坐标转换而变化的标量，标量值与坐标系的旋转无关，因此也与参考基无关。

进一步的技术方案，针对球体模型，假设均匀密度的球体在地面的投影坐标为(x0,y0,0)，球体中心埋深为h，半径为r，万有引力常数为g，剩余密度(与周围岩体密度差)为ρ，则剩余质量m＝4πr³ρ/3，在地面观测点产生的各个重力梯度张量分量异常为：

进一步的技术方案，将球体重力梯度张量分量公式(5)～(10)带入三个特征值求解公式(1)～(3)中，并计算其特征值组合的张量不变量μ，有

则

进一步的技术方案，令总水平导数为tdx，并带入公式(5)、(6)，即

构造公式et并根据公式(11)、(12)简化，

简化后：

令即表示测量点与目标中心的水平距离，则

通过公式et建立水平距离与埋深之间的关系，且与目标剩余质量无关。

进一步的技术方案，当et为45°时，即r²+h²＝3rh。

其方程解r1,r2为

即

利用et的两条45°等值线之间距离的0.4472倍进行深度估计。

以上一个或多个技术方案存在以下有益效果：

本发明针对地质体中的球体模型，提供了基于重磁梯度数据张量不变量的目标深度估计方法，对于单独地质体或距离较远的多个地质体情况下的深度估计问题，可以获得较高精度的计算结果。

利用重力梯度张量数据处理解释时，其数据量更多，蕴含信息更丰富；采用特征值组合的张量不变量，其与坐标系的旋转无关，也与参考基无关，能够减少磁倾角和磁偏角的影响，且特征值在数据挖掘中也具有重要作用。

相对于二维剖面曲线数据，利用平面数据进行处理解释时，其数据更具有规律性，准确率更高且抗噪能力更强。

附图说明

构成本发明的一部分的说明书附图用来提供对本发明的进一步理解，本发明的示意性实施例及其说明用于解释本发明，并不构成对本发明的不当限定。

图1(a)-图1(f)为本公开实施例子的球体模型重力梯度张量分量图；

图2(a)-图2(b)为本公开实施例子的模型体特征值组合的不变量μ和反正切比值et图；

图3为本公开实施例子的反正切比值et的45°等值线图；

图4为本公开实施例子的方法的流程示意图。

具体实施方式

应该指出，以下详细说明都是示例性的，旨在对本发明提供进一步的说明。除非另有指明，本文使用的所有技术和科学术语具有与本发明所属技术领域的普通技术人员通常理解的相同含义。

需要注意的是，这里所使用的术语仅是为了描述具体实施方式，而非意图限制根据本发明的示例性实施方式。如在这里所使用的，除非上下文另外明确指出，否则单数形式也意图包括复数形式，此外，还应当理解的是，当在本说明书中使用术语“包含”和/或“包括”时，其指明存在特征、步骤、操作、器件、组件和/或它们的组合。

在不冲突的情况下，本发明中的实施例及实施例中的特征可以相互组合。

实施例一

本实施例公开了基于重磁梯度数据张量不变量的目标深度估计方法。

参见附图4为本公开方法的流程示意图，包括以下步骤：

步骤(1)：根据球体模型对象的重力梯度张量数据，计算其梯度矩阵的特征值及特征值组合的张量不变量；

步骤(2)：利用张量不变量与总水平导数比值的反正切值，建立水平距离与埋深之间的关系。

步骤(3)：通过反正切比值图上的两条45°等值线之间的距离，对目标深度进行估计。

所述步骤(1)中，具体内容包括：

(1-1)重力场是引力位在笛卡尔坐标系x,y,z方向上的矢量分量，重力梯度张量为地质体的引力位的二阶导数，即重力三分量gx,gy,gz在x,y,z三个方向上的导数，其中gx在x，y，z方向的一阶梯度分别为gxx,gxy,gxz，gy在x，y，z方向的一阶梯度分别为gyx,gyy,gyz，gz在x，y，z方向的一阶梯度分别为gzx,gzy,gzz，可构成梯度张量矩阵：

其中矩阵t为对称矩阵，且重力位满足拉普拉斯方程，即有gxy＝gyx，gxz＝gzx，gyz＝gzy，gxx+gyy+gzz＝0，即矩阵t有5个独立的分量。

(1-2)由于重力梯度张量矩阵t为一个实对称阵，可以通过计算特征值将其对称化，即张量矩阵t可以通过坐标转换成对称阵λ：

e^tte＝λ

其特征值λ和特征向量e分别为：

e＝[e1e2e3]

其中，特征向量e表示新笛卡尔坐标系的正交基底，在任何坐标转换下，矩阵t将包含三个不变量要素：

i0＝trace(t)＝gxx+gyy+gzz＝0

可以通过三个不变量要素，得到求解特征值的参数方程：

λ³-i0λ²+i1λ-i2＝0

求其方程的根，得到三个特征值λ1,λ2,λ3为

λ1＝a+b(1)

其中，

即三个特征值可利用重力梯度张量矩阵中的各个分量表示。

(1-3)计算三个特征值组合的不变量μ，即

其中不变量μ可利用重力梯度张量的各个分量表示，且为不随坐标转换而变化的标量，标量值与坐标系的旋转无关，因此也与参考基无关。

所述步骤(2)中，具体包括：

(2-1)针对球体模型，假设均匀密度的球体在地面的投影坐标为(x0,y0,0)，球体中心埋深为h，半径为r，万有引力常数为g，剩余密度(与周围岩体密度差)为ρ，则剩余质量m＝4πr³ρ/3，在地面观测点产生的各个重力梯度张量分量异常为：

(2-2)将球体重力梯度张量分量公式(5)～(10)带入三个特征值求解公式(1)～(3)中，并计算其特征值组合的张量不变量μ，有

则

(2-3)令总水平导数为tdx，并带入公式(5)、(6)，即

构造公式et并根据公式(11)、(12)简化，

简化后：

令即表示测量点与目标中心的水平距离，则

通过公式et建立水平距离与埋深之间的关系，且与目标剩余质量无关。

所述步骤(3)中，具体包括：

(3-1)当et为45°时，即r²+h²＝3rh

其方程解r1,r2为

即

利用et的两条45°等值线之间距离的0.4472倍进行深度估计。

上述实施例子中，对含有球体模型的目标体空间，得到其引力位在x,y,z三个方向上的二阶导数，构成重力梯度张量数据平面图。

通过各个分量的重力梯度张量数据，计算其梯度矩阵的特征值及特征值组合的张量不变量，其值与坐标系的旋转无关，也与参考基无关。

利用张量不变量与总水平导数数据，通过反正切运算，建立水平距离与埋深之间的关系，得到其数据平面图。

通过反正切运算的数据平面图，利用图上两条45°等值线之间的距离，对目标体的深度进行估计。

本发明针对地质体中的球体模型，可近似代表自然界中一些轴状地质体，包括矿巢、岩株、穹隆、空洞等地质结构，提供了基于重磁梯度数据张量不变量的目标深度估计方法，对于单独地质体或距离较远的多个地质体情况下的深度估计问题，可以获得较高精度的计算结果。

为了验证申请的构思的方案的效果，下面进行仿真：

对一球体模型进行仿真，设置其各个参数：球体深度h为30m，其水平坐标(x0,y0)为坐标原点(0,0)，万有引力常数g＝6.67×10^-11n·m²/kg²，半径r为1m,剩余密度ρ为500kg/m³，单位转换有1m/s²＝10⁵mgal，测量范围为-100m到+100m，构成各重力梯度张量数据gxx,gxy,gxz,gyy,gyz,gzz平面数据图，如图1(a)-图1(f)。

计算其重力梯度张量数据的特征值及特征值组合的张量不变量μ，其值与坐标系的旋转无关，也与参考基无关，如图2(a)所示，并计算与总水平导数比值的反正切值et，建立水平距离与埋深之间的关系，如图2(b)所示。

通过反正切比值图et上的两条45°等值线之间的距离，计算其目标深度：如图3所示。

其等值线之间的距离为：

δr＝78.5417-11.4763＝67.0654m

目标的深度估计为：

h≈0.4472δr＝29.9916m

这种方法对于单独地质体或距离较远的多个地质体情况下的深度估计问题，可以获得较高精度的计算结果。

实施例二

本实施例的目的是提供一种计算装置,包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述程序时实现以下步骤，包括：

步骤(1)：根据球体模型对象的重力梯度张量数据，计算其梯度矩阵的特征值及特征值组合的张量不变量；

步骤(2)：利用张量不变量与总水平导数比值的反正切值，建立水平距离与埋深之间的关系。

步骤(3)：通过反正切比值图上的两条45°等值线之间的距离，对目标深度进行估计。

实施例三

本实施例的目的是提供一种计算机可读存储介质。

一种计算机可读存储介质，其上存储有计算机程序，该程序被处理器执行时执行以下步骤：

步骤(1)：根据球体模型对象的重力梯度张量数据，计算其梯度矩阵的特征值及特征值组合的张量不变量；

步骤(2)：利用张量不变量与总水平导数比值的反正切值，建立水平距离与埋深之间的关系。

步骤(3)：通过反正切比值图上的两条45°等值线之间的距离，对目标深度进行估计。

实施例四

本实施例的目的是提供基于重磁梯度数据张量不变量的目标深度估计系统，包括服务器，所述服务器被配置为：

步骤(1)：根据球体模型对象的重力梯度张量数据，计算其梯度矩阵的特征值及特征值组合的张量不变量；

步骤(2)：利用张量不变量与总水平导数比值的反正切值，建立水平距离与埋深之间的关系。

步骤(3)：通过反正切比值图上的两条45°等值线之间的距离，对目标深度进行估计。

以上实施例二、三和四的装置中涉及的各步骤与方法实施例一相对应，具体实施方式可参见实施例一的相关说明部分。术语“计算机可读存储介质”应该理解为包括一个或多个指令集的单个介质或多个介质；还应当被理解为包括任何介质，所述任何介质能够存储、编码或承载用于由处理器执行的指令集并使处理器执行本发明中的任一方法。

本领域技术人员应该明白，上述本发明的各模块或各步骤可以用通用的计算机装置来实现，可选地，它们可以用计算装置可执行的程序代码来实现，从而，可以将它们存储在存储装置中由计算装置来执行，或者将它们分别制作成各个集成电路模块，或者将它们中的多个模块或步骤制作成单个集成电路模块来实现。本发明不限制于任何特定的硬件和软件的结合。

以上所述仅为本发明的优选实施例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

上述虽然结合附图对本发明的具体实施方式进行了描述，但并非对本发明保护范围的限制，所属领域技术人员应该明白，在本发明的技术方案的基础上，本领域技术人员不需要付出创造性劳动即可做出的各种修改或变形仍在本发明的保护范围以内。