CN111323830A - 一种基于大地电磁和直流电阻率数据的联合反演方法 - Google Patents

Abstract

本发明公开了一种基于大地电磁和直流电阻率数据的联合反演方法，包括：根据观测数据确定反演区域，并进行剖分得到初始网格，设定反演参数向量的初值；基于反演参数向量的当前值进行高斯‑牛顿反演迭代计算出模型改变量以及线性搜索步长，再计算出反演参数向量的迭代更新值；判断是否达到反演终止条件，若未达到，进行下一次迭代；若达到，判断是否达到渐进网格反演终止条件，若达到，反演结束；若未达到，对网格单元进行细化剖分并更新反演参数向量值，再返回进行高斯‑牛顿反演迭代，直至反演结束。本发明研发出反演网格自适应调整的联合反演技术，有效地降低大地电磁法和直流电阻率数据的联合反演解释的多解性难题，提高联合解释准确度。

Description

一种基于大地电磁和直流电阻率数据的联合反演方法

技术领域

本发明涉及到应用地球物理领域，尤其是涉及一种基于大地电磁和直流电阻率数据的联合反演方法。

背景技术

在应用地球物理勘探领域，大地电磁法法(MT)和直流电阻法(DC)被广泛应用到矿产资源勘探中，并且应用效果显著。但是，单一的大地电磁法法或直流电阻法方法难以避免反演多解性问题，影响资料解释的准确度，从而衍生出了联合反演解释技术。联合反演能够可有效降低反演的多解性，提高资料解释的精度和可靠性。除此之外，传统联合反演方式大多数都是在固定的正反演网格下进行，这技术给资料解释带来极其的不稳定性。这是因为高精度的正演算法是进行反演前提条件，要获取高精度的正演结果势必会造成网格密度的增加，从而导致反演的网格单元数急剧增多，从而加重了反演的多解性；其次，倘若降低反演网格的密度来降低反演的多解性，势必影响正演求解精度，从而影响反演结果的准确度。因此，传统联合反演方式中选用固定的正反演网格，其中，若设定的网格密度低会降低反演结果的准确度，若设定的网格密度过高，会加重反演的多解性，如何获取到合适的网格密度是亟待研究的。

发明内容

本发明的目的是提供了一种基于大地电磁和直流电阻率数据的联合反演方法，具体提供了一种反演网格自适应调整的联合反演技术，构建由粗网格到细网格的渐进网格的反演策略，极大改善了反演网格和正演网格之间平衡性，以及有效降低反演面临多解性难题。

其中，本发明提供的一种基于大地电磁和直流电阻率数据的联合反演方法，包括如下步骤：

步骤S1：采集观测点的观测数据，所述观测数据包括大地电磁法观测的视电阻率和相位以及直流电阻率法观测的视电阻率；

步骤S2：根据观测数据确定反演区域，并对所述反演区域进行剖分得到初始网格，并设定所述初始网格对应的反演参数向量的初值，所述反演参数向量的元素值与网格单元的电阻率值相关；

步骤S3：基于反演参数向量的当前值进行高斯-牛顿反演迭代计算出模型改变量以及线性搜索步长，再根据模型改变量以及线性搜索步长计算出反演参数向量的迭代更新值；

其中，所述高斯-牛顿反演迭代对应的方程、线性搜索步长、以及迭代更新值的公式如下：

m_k+1＝m_k+λδm_k

式中，

J_k为第k次对应的雅克比矩阵J，J_k(MT)为观测数据中大地电磁法数据对第k次迭代对应的反演参数向量的偏导数矩阵，J_k(DC)为观测数据中直流电阻率法数据对第k次迭代对应的反演参数向量的偏导数矩阵，m_k、m_k+1分别为第k次、k+1次迭代对应的反演参数向量，k＝1时，m_k为反演参数向量的初值；W_d为数据方差矩阵，W_m为模型误差矩阵，m_apr为先验模型，λ为线性搜索步长，δm_k为第k次迭代对应的模型改变量，μ为正则化因子，δd_k为第k次迭代对应的正演计算值与观测数据的差；

判断是否达到反演终止条件，若未达到，返回步骤S3进行下一次迭代；

若达到，判断是否达到渐进网格反演终止条件，若达到，反演结束；若未达到，对网格单元进行细化剖分并更新反演参数向量值，再返回步骤S3进行高斯-牛顿反演迭代，其中，基于细化后新增的网格单元个数更新反演参数向量中的元素个数，所述渐进网格反演终止条件为网格剖分终止条件。

在反演参数向量中细网格单元对应的元素值与被细化的网格单元的元素值相等。

本发明提供渐进网格反演技术，反演从粗网格到细网格逐步进行，有效且准确地实现自适应的网格划分，得到恰当的网格密度，从而降低反演的多解性，提高反演结果的准确度。其中，根据反演结束后得到反演参数向量可以计算出电阻率。

进一步优选，步骤S4中对网格单元再进行剖分时的剖分规则如下：

按照如下公式计算当前每个网格单元对应的模型灵敏度，并对模型灵敏度超过预设阈值的网格单元进行标记，再对标记的网格单元进行剖分；

式中，S_j为第j个网格单元对应的模型灵敏度，J_ij为当前雅克比矩阵J中第i个观测点第j个网格单元对应的元素，即第i行第j列元素或者说是第i个观测点对应第j个网格单元参数的偏导数，n为观测点的个数；

雅克比矩阵J如下所示：

式中，

分别表示为第1个、第2个、第n个观测点的观测数据，m₁、m₂、m_M分别表示为第1个、第2个、第M个网格单元对应在反演参数向量中的值。

正则化反演技术严重依赖于先验信息，当先验信息不足或者不正确时将会产生错误的反演结果，减少反演网格单元的数量与提高反演单元质量是实现稳定反演的重要途径。生成高质量网格的难度在于如何选择单元细化或组合的标准，最直接的方法是通过分析模型分辨率矩阵来实现。本发明的研究工作利用模型灵敏度先验信息，进行网格自适应生成。本发明网格细化以模型灵敏度为依据进行局部细化，通过S_j对网格进行优化实质是调整偏导数矩阵，从而改善偏导数性质，进而得到高质量的反演网格，同时避免了过度细化网格。

进一步优选，判断是否达到所述渐进网格反演终止条件的过程为：计算模型灵敏度的方差，若当前方差与前一次的模型灵敏度的方差的变化值在预设变化范围，当前剖分合理，达到渐进网格反演终止条件；否则，未达到渐进网格反演终止条件，继续选择网格单元进行细化剖分。

进一步优选，所述反演参数向量的元素个数与网格单元数相对应，所述初始网格对应的反演参数向量的初值获取过程如下：

首先，获取反演区域的电阻率的下限和上限值；

其次，分别在所述下限和上限值的区间范围取值作为网格单元的电阻率值，再按照如下公式计算反演参数向量初值中各个元素值；

式中，m_j为第j个网格单元对应在反演参数向量中的元素值，ρ_j为第j个网格单元对应的电阻率值，通常情况模型电阻率初值采用均匀半空间值，ρ_min、ρ_max分别为预设的反演区域的电阻率的下限和上限值。

进一步优选，所述数据方差矩阵如下所示：

式中，ε为确保分母不为零的一个任意正实数，χ_i为第i个观测点的观测数据方差，α_i为第i个观测点对应的数据权重系数。

进一步优选，数据权重系数α_i的取值如下：

式中，N_MT为大地电磁法参与反演的数据个数、N_DC为直流电阻率法参与反演的数据个数。

进一步优选，步骤S4中所述反演终止条件为：当前迭代对应的数据均方差误差小于预设阈值或者迭代次数达到预设迭代总数；

RMS为数据均方根误差，

表示第i个观测点的观测数据，

表示当前反演迭代后再进行正演计算得到的第i个观测点对应的数据，n为观测点的个数。

进一步优选，所述高斯-牛顿反演迭代对应的方程是依据目标函数最优化求解得到，所述目标函数如下：

式中，

表示目标函数，d_obs为反演数据向量；μ为正则化因子；

为数据误差函数；

为模型误差函数。

进一步优选，所述模型误差函数选用最小结构模型约束，如下所示：

式中，β_s，β_x，β_y为模型误差三部分之间的比例系数，Ω为反演区域或网格区域，x,y为x.y坐标方向。

进一步优选，随着反演迭代，所述正则化因子μ逐步减小。

有益效果

1、本发明提供的所述方法极大改善了反演网格和正演网格之间平衡性，以及有效降低反演面临多解性难题。其无需顾忌正演网格密度，通过迭代过程实现了网格自适应的划分，解决了网格划分固定带来的弊端。

2、本发明的优选方案中利用模型灵敏度为依据优化反演网格，对模型灵敏度高的网格单元进一步剖分，生成高质量的反演网格，构建了从粗网格到细网格的渐进网格的反演测量。

附图说明

图1为为根据本发明提供的一种基于大地电磁和直流电阻率数据的联合反演方法的流程图；

图2为均匀半空间下存在38个异常体，异常体均为截面为矩形的四棱柱，低阻体电阻率为10Ω·m，高阻体电阻率为1000Ω·m，背景电阻率为100Ω·m；

图3为理论模型(如图1所示)TM模式正演结果；

图4理论模型反演得到电阻率示意图；

图5渐进网格反演得到的每次格变化图。

具体实施方式

下面将结合实施例对本发明做进一步的说明。

本发明提供的一种基于大地电磁和直流电阻率数据的联合反演方法，实现了渐进网格反演技术，反演从粗网格到细网格逐步进行，有效且准确地实现自适应的网格划分，得到恰当的网格密度，从而降低反演的多解性，提高反演结果的准确度。其中，根据反演结束后得到反演参数向量可以计算出电阻率。

首先，根据地球物理反问题的病态性，构建正则化反演的目标函数

为：

式中：m为反演参数向量，d_obs为反演数据向量，μ为正则化因子；

为数据误差函数；

为模型误差函数。

由于DC与MT联合反演的变量为电阻率，采用变换函数以确保反演电阻率在实际岩、矿石的取值范围内。构建反演参数向量与电阻率的关系，由以下公式(2)可知，若电阻率已知，通过公式2可以计算出反演参数向量；反之，若已知反演参数向量，则可以计算出电阻率。

式中，m_j为反演空间参数，即第j个网格单元对应在反演参数向量中的元素值；ρ_j为反演单元的电阻率值；ρ_min、ρ_max分别为研究区电阻率的下限和上限。

反演数据由MT与DC数据组成，即观测数据

其中，ρ_TM,a、

为大地电磁TM(Transverse Magnetic)模式视电阻率与相位；ρ_TE,a、

为TE(Transverse Electric)模式视电阻率与相位；ρ_DC,s为直流电阻率法观测的视电阻率。数据误差函数可表示为如下统一形式：

或者表示为：

式中：G(m)为正演函数，分别由MT与DC两种方法观测数据组成，由于正演函数是本领域的公知常识，本发明对此不进行具体的阐述；W_d为n×n的数据方差矩阵，其通过下面的公式进行构造：

其中，ε为确保分母不为零的一个任意正实数，n为观测数据的个数；χ_i为观测数据方差；α_i为数据权重系数，其作用是平衡MT与DC数据在反演中的权重，以避免由于两种数据间误差范围的不同造成某一数据集的权重过大，而另一数据不起作用。在数据误差为高斯噪声的假设下，反演后数据误差函数期望值为N_MT+N_DC，N_MT为MT参与反演的数据个数、N_DC为DC参与反演的数据个数，此两个参数均是已知参数。如果希望MT与DC两个数据集在反演中取相同比重，数据权重系数可表示为：

先验信息通常根据模型空间分布特征进行构造，根据其反演后模型的特点可大致分为平滑与陡变两种约束。其中平滑模型一般利用模型参数空间梯度进行计算。支持陡变模型的约束方法有多种，如最大变化量支持模型、最小支持模型、最小梯度支持模型等；另外，通过模型误差的L1与L2范数形式可以实现模型分块与平滑形式的反演约束。本文实施例优选采用最小结构模型约束，其他可行的实施例中，可以选用其他的，本发明对此不进行具体的限定，本发明实施例中，模型误差函数如下：

式中：β_s，β_x，β_y为模型误差三部分之间的比例系数；m_apr为先验模型，Ω为反演区域或网格区域，x,y为x.y坐标方向，其中，模型误差可统一表示为：

式中，W_m为统一定义形式等价模型误差矩阵。

根据以上公式可开展反演目标函数最优化求解问题。将式(3)、(7)代入式(1)，可得第k次高斯-牛顿迭代公式为：

式中：

J_k为第k次对应的雅克比矩阵J，J_k(MT)为MT对第k次迭代对应的反演参数向量的偏导数矩阵，J_k(DC)为DC对第k次迭代对应的反演参数向量的偏导数矩阵，雅克比矩阵J的计算式如下：

J是偏导数矩阵，并且是满阵，其利用观测数据对模型的导数来进行计算，储存偏导数矩阵需要较大的内存，计算偏导数矩阵往往也是反问题求解中最耗时的部分。另外，公式(8) 中，数据变化量为：

其中，G_MT(m_k)表示第k模型的MT方法的正演解，G_DC(m_k)表示第k模型的DC的正演解，

表示MT的观测数据，d_{DC_obs}＝d_obs(ρ_DC,s)表示 DC的观测数据。将公式(10)记为第k次正演得到的计算数据与观测数据的差。为确保每一次迭代目标函数可以充分下降，对模型更新步长进行线性搜索，表达式为：

m_k+1＝m_k+λδm_k (11)

λ为线性搜索步长，m_k+1为k+1次的电阻率模型，m_k为k次的电阻率模型，δm_k为第k 次的电阻率模型改变量，满足：

则，

为了求解公式(8)本发明采用双共轭梯度迭代算法(BICGSTAB)，实现大地电磁和直流电阻率满足的高斯-牛顿方程高效、稳定正则化反演计算。

另外，正则化反演技术严重依赖于先验信息，当先验信息不足或者不正确时将会产生错误的反演结果，减少反演网格单元的数量与提高反演单元质量是实现稳定反演的重要途径。生成高质量网格的难度在于如何选择单元细化或组合的标准，最直接的方法是通过分析模型分辨率矩阵来实现。本发明的研究工作利用模型灵敏度先验信息，进行网格自适应生成，该模型灵敏度计算表达式为：

其中，i表示观测数据的个数，j表示离散单元的个数，模型灵敏度S_j是第j个模型m_j改变时对整个观测数据集的影响。不考虑数据方差与稳定因子的条件下，S_j实质上反映了高斯 -牛顿优化过程中偏导数矩阵。通过S_j对网格进行优化实质是调整偏导数矩阵的主对角元素，从而改善偏导数性质，进而得到高质量的反演网格。

基于研发的网格优化算法，进一步研究了一种渐进网格反演技术。本发明提出的渐进网格反演技术以模型灵敏度为依据进行网格细化，反演由粗网格到细网格逐步进行，从而降低反演的多解性，提高反演结果的准确度。

基于上述原理性内容的阐述，本发明提供的一种基于大地电磁和直流电阻率数据的联合反演方法，包括如下步骤：

步骤S1：采集观测点的观测数据，所述观测数据包括大地电磁法观测的视电阻率和相位以及直流电阻率法观测的视电阻率。

采集观测点的观测数据，即MT和DC的观测数据

步骤S2：根据观测数据确定反演区域，并对所述反演区域进行剖分得到初始网格，并设定所述初始网格对应的反演参数向量的初值。

其中，根据观测点的坐标来确定反演区域，反演区域包含存在数据的观测点。本实施例中，采用网格剖分器实现反演区域的非结构化三角单元初始剖分。本实施例中，还优选设定最小网格单元面积信息，即初始划分或者后续剖分过程中，要保证网格单元的面积不得小于最小网格单元面积，以免剖分过细。

步骤S3：基于反演参数向量的当前值进行高斯-牛顿反演迭代计算出模型改变量以及计算出线性搜索步长，再根据模型改变量以及线性搜索步长计算出反演参数向量的迭代更新值。

其中，求解区域的偏导数矩阵

数据方差矩阵W_d，再利用高斯-牛顿反演迭代的方程计算模型改变量以及计算出线性搜索步长，最终实现了反演参数向量的更新。

此外，正则化因子μ是一个平衡参数，其作用是平衡数据误函数和模型误函数之间的关系，通常情况下反演是先拟合模型误差函数，故本实施例中初始将μ设置一个较大值，在反演逐步进行中改变两种之间的关系，而减小μ直至反演结束。其他可行的实施例中，也可以选择固定值，本发明对此不进行具体的限定。

步骤S4：判断是否达到反演终止条件，若未达到，返回步骤S3进行下一次迭代；

若达到，判断是否达到渐进网格反演终止条件，若达到，反演结束；若未达到，对网格单元进行细化剖分并更新反演参数向量值，再返回步骤S3进行高斯-牛顿反演迭代，其中，基于细化后新增的网格单元个数更新反演参数向量中的元素个数，所述渐进网格反演终止条件为网格剖分终止条件。

其中，终止条件为：数据均方根误差

小于预设阈值或或达到设定的迭代次数。关于预设阈值或，通常会设置一个较小的值，

表示第i个观测点的数据，

表示每次反演迭代后得到的模型进行正演计算得到的观测点处的数据。应当理解，本发明中两终止条件中满足其中一个都停止反演。退出反演；若不满足，则还需要进行剖分，其中，规则如下：

不满足进入互换定理计算雅克比矩阵

并计算模型灵敏度S_j，分析模型灵敏度 S_j的平均值

方差

首先对模型灵敏度灵敏度大的单元进行标记，然后对标记好的单元进行细化处理，从而生成最优化反演网格单元，在这过程中需要评价每次迭代过程的方差，确保前后几次的方差的变化范围较小来整体评价反演网格单元剖分的好坏，如果结果合理然后循环跳步骤S3直至反演结束，得到最终的反演结果。

设计如图2所示的二维模型，背景电阻率为100Ω·m，共设置38个异常体，异常体均为截面为矩形的四棱柱，低阻体电阻率为10Ω·m，高阻体电阻率为1000Ω·m。其中，接近地表交替有14个低阻和14个高阻异常体，异常体截面宽10m、高30m、上顶埋深为20m，相邻两个异常体间隔为100m；稍深交替有4个低阻与3个高阻体，异常体截面宽100m、高100m、上顶埋深为100m，相邻两个异常体间隔为400m；再向下，位于断面中部有1个低阻体，异常体截面宽400m、高200m、上顶埋深为300m；最深部有1个低阻和1个高阻异常体，异常体截面宽500m、高500m，上顶埋深为500m，两个异常体间隔为600m。

近地表异常体用于模拟MT的静态效应，静态效应严重干扰视电阻率，通常采用校正方法进行压制；为避免校正方法引起的数据误差，可通过其它浅层方法进行补充，进而开展综合解释和联合反演以提高解释精度。

直流电阻率法采用二极装置，共布设电极320根，相邻电极距10m，最大极距200m，数据集点数6000个；AMT测点距40m，共80个测点，测量频段1～4096Hz，按2的指数共 13个频点，将正演数据加入3％随机噪声作为反演数据集。

图3为正演得到的TM模式卡尼亚电阻率与相位。二维条件下TM模式卡尼亚电阻率受静态效应影响严重，断面图出现条带状异常，无法识别深部异常体；静态效应对相位影响小，相较视电阻率断面图可以分辨更多的异常体。

采用本文介绍的渐进反演算法，分别对DC、AMT、DC与AMT数据进行渐进反演。反演过程对网格进行5次细化，反演结果如图4所示。

图4(a)为直接反演DC数据所得断面，反演结果没有得到深部异常体；近地表的28个异常体反演效果较好，但受稍深的7个异常体影响存在畸变。图4(b)为AMT反演结果，从反演断面可以判断出几乎所有异常体，但是很多异常相互联结无法分离；近地表异常反演效果不理想，异常体位置有高、低阻体出现，但反演所得电阻率值与真实值相差较大，说明模拟设定的点距与频率条件下AMT对浅层不均匀体有一定的反演能力。图4(c)为联合反演两个数据集的反演结果，反演有效地给出了浅层与深部所有异常的位置；相比于单独反演，近地表28个异常体的电阻率与空间位置均较单独反演更接近真实值；稍深的7个异常体形态与设定模型有一定的出入，但异常体的位置与反演电阻率值均较单独反演精确；只有联合反演可分辨断面中部独立的低阻异常体；深部低阻异常的反演效果优于高阻，分析原因为在本文所采用的模型设置下，低阻二度体相较高阻二度体更易形成沿走向方向的电流场，有利于反演。

为进一步分析渐进反演过程，现对MT反演进行分析，反演过程网格变化如图5所示。反演区单元数从7836到17932，逐渐增加，第1、2次网格反演，数据均方根误差(RMS) 值快速下降，后面优化网格RMS值下降量逐渐减少，趋于稳定；反演过程中，初始网格正则化因子出现增大现象，随后稳定下降，表明渐进网格方法对正则化因子的初始依赖不大，可以通过自动调整确保反演稳定进行；最后两重网格数据误差下降不多，更多的对模型进行优化；由于地表的灵敏度较大，网格总体由上到下逐渐加密，从图5看，由电阻率分布引起的网格变化过程清晰，体现了基于模型灵敏度的网格优化特点。

需要强调的是，本发明所述的实例是说明性的，而不是限定性的，因此本发明不限于具体实施方式中所述的实例，凡是由本领域技术人员根据本发明的技术方案得出的其他实施方式，不脱离本发明宗旨和范围的，不论是修改还是替换，同样属于本发明的保护范围。

Claims (10)

1.一种基于大地电磁和直流电阻率数据的联合反演方法，其特征在于：包括如下步骤：

步骤S1：采集观测点的观测数据，所述观测数据包括大地电磁法观测的视电阻率和相位以及直流电阻率法观测的视电阻率；

步骤S2：根据观测数据确定反演区域，并对所述反演区域进行剖分得到初始网格，并设定所述初始网格对应的反演参数向量的初值；

步骤S3：基于反演参数向量的当前值进行高斯-牛顿反演迭代计算出模型改变量以及线性搜索步长，再根据模型改变量以及线性搜索步长计算出反演参数向量的迭代更新值；

其中，所述高斯-牛顿反演迭代对应的方程、线性搜索步长、以及迭代更新值的公式如下：

m_k+1＝m_k+λδm_k

式中，

J_k为第k次对应的雅克比矩阵J，J_k(MT)为观测数据中大地电磁法数据对第k次迭代对应的反演参数向量的偏导数矩阵，J_k(DC)为观测数据中直流电阻率法数据对第k次迭代对应的反演参数向量的偏导数矩阵，m_k、m_k+1分别为第k次、k+1次迭代对应的反演参数向量，k＝1时，m_k为反演参数向量的初值；W_d为数据方差矩阵，W_m为模型误差矩阵，m_apr为先验模型，λ为线性搜索步长，δm_k为第k次迭代对应的模型改变量，μ为正则化因子，δd_k为第k次迭代对应的正演计算值与观测数据的差；

步骤S4：判断是否达到反演终止条件，若未达到，返回步骤S3进行下一次迭代；

若达到，判断是否达到渐进网格反演终止条件，若达到，反演结束；若未达到，对网格单元进行细化剖分并更新反演参数向量值，再返回步骤S3进行高斯-牛顿反演迭代，其中，基于细化后新增的网格单元个数更新反演参数向量中的元素个数，所述渐进网格反演终止条件为网格剖分终止条件。

2.根据权利要求1所述的方法，其特征在于：步骤S4中对网格单元再进行剖分时的剖分规则如下：

按照如下公式计算当前每个网格单元对应的模型灵敏度，并对模型灵敏度超过预设阈值的网格单元进行标记，再对标记的网格单元进行剖分；

式中，S_j为第j个网格单元对应的模型灵敏度；J_ij为当前雅克比矩阵J中第i个观测点第j个网格单元对应的元素，n为观测点的个数；

雅克比矩阵J如下所示：

式中，

分别表示为第1个、第2个、第n个观测点的观测数据，m₁、m₂、m_M分别表示为第1个、第2个、第M个网格单元对应在反演参数向量中的值。

3.根据权利要求2所述的方法，其特征在于：判断是否达到所述渐进网格反演终止条件的过程为：计算模型灵敏度的方差，若当前方差与前一次的模型灵敏度的方差的变化值在预设变化范围，当前剖分合理，达到渐进网格反演终止条件；否则，未达到渐进网格反演终止条件，继续选择网格单元进行细化剖分。

4.根据权利要求1所述的方法，其特征在于：所述反演参数向量的元素个数与网格单元数相对应，所述初始网格对应的反演参数向量的初值获取过程如下：

首先，获取反演区域的电阻率的下限和上限值；

其次，分别在所述下限和上限值的区间范围取值作为网格单元的电阻率值，再按照如下公式计算反演参数向量初值中各个元素值；

式中，m_j为第j个网格单元对应在反演参数向量中的元素值，ρ_j为第j个网格单元对应的电阻率值，ρ_min、ρ_max分别为预设的反演区域的电阻率的下限和上限值。

5.根据权利要求1所述的方法，其特征在于：所述数据方差矩阵如下所示：

式中，ε为确保分母不为零的一个任意正实数，χ_i为第i个观测点的观测数据方差，α_i为第i个观测点对应的数据权重系数。

6.根据权利要求5所述的方法，其特征在于：数据权重系数α_i的取值如下：

式中，N_MT为大地电磁法参与反演的数据个数、N_DC为直流电阻率法参与反演的数据个数。

7.根据权利要求1所述的方法，其特征在于：步骤S4中所述反演终止条件为：当前迭代对应的数据均方差误差小于预设阈值或者迭代次数达到预设迭代总数；

RMS为数据均方根误差，

表示第i个观测点的观测数据，

表示当前反演迭代后再进行正演计算得到的第i个观测点对应的数据，n为观测点的个数。

8.根据权利要求1所述的方法，其特征在于：所述高斯-牛顿反演迭代对应的方程是依据目标函数最优化求解得到，所述目标函数如下：

式中，

表示目标函数，d_obs为反演数据向量；μ为正则化因子；

为数据误差函数，

为模型误差函数。

9.根据权利要求8所述的方法，其特征在于：所述模型误差函数选用最小结构模型约束，如下所示：

式中，β_s，β_x，β_y为模型误差三部分之间的比例系数，Ω为反演区域或网格区域，x,y为x.y坐标方向。

10.根据权利要求1所述的方法，其特征在于：随着反演迭代，所述正则化因子μ逐步减小。

Images