CN104102814A - 一种基于大地电磁数据反演电阻率和磁化率的方法及系统 - Google Patents

Abstract

本发明涉及一种基于大地电磁数据反演电阻率和磁化率的方法，所述方法包括：步骤1，在TM极化模式下测量地面的视电阻率实测值，并建立初始地质模型，获得初始地质模型各层的电阻率和磁化率；步骤2，在水平层状介质条件下，基于初始地质模型的电阻率和磁化率，采用一维正演方法计算地面的视电阻率理论值；步骤3，逐层修正初始地质模型的电阻率和磁化率；步骤4，基于修正后的电阻率和磁化率，重新计算地面的视电阻率理论值，若视电阻率实测值和重新计算的视电阻率理论值的拟合差小于设定的误差期望值，则停止对电阻率和磁化率的修正。本发明的方法不用进行繁复偏导数矩阵计算，迭代过程稳定收敛，可逐层判断地层的地电阻率参数和磁化率参数。

Description

一种基于大地电磁数据反演电阻率和磁化率的方法及系统

技术领域

本发明涉及电法勘探领域，特别是涉及一种基于大地电磁数据反演电阻率和磁化率的方法及系统，其适用于探测隐伏的地下磁性矿体或隐伏的含磁性岩体的解释，

背景技术

大地电磁法作为一种重要的勘探地球物理方法，应用领域广泛，分支众多。自首次提出之后，20世纪60年代相继提出了音频大地电磁法(AMT)，1971年，加拿大多伦多大学的D.W.Strangbway教授和他的研究生MyronGoldstein又提出了可控源音频大地电磁法(CSAMT)。

目前，大地电磁在数据处理方面，其主流方法依然是电阻率反演，2002年苏东刘给出了大地电磁电阻率一维“正演修正”法反演。而对于磁导率参数，Ward(1961)和Fraser(1973)通过观测高磁化率薄板地质体附近的电磁数据，表明磁化率的影响是可以被观测到的。Zhiyi Zhang和Douglas W.Oldenburg(1997)在假定电阻率分布已知的情况下，从一维模型的基础上恢复出了磁化率的分布。

Les P.Beard和Jonathan E.Nyquist(1998)在航空电磁法的研究中指出，在地下地质体磁导率远大于真空中磁导率的地方，磁导率的影响在航空电磁勘探中产生了一个与频率无关的相位偏移，而这种偏移随磁导率的增加而增大，随磁性物质的规模增加而增大，随感应器接近地球表面距离的减小而增大。在实际应用中表明，当忽略磁化率的影响时，反演得到的电阻率分布不符合实际目标体情况。

从大地电磁数据中反演电阻率和磁化率参数，由于二者具有非线性相关性，所需要的理论和算法都比较复杂，要准确和精细的反演恢复出来，具有相当的难度。在众多反演方法中比较经典的反演方法如：Occam反演、非线性共轭梯度反演(NLCG)等，已经相当完善，虽然前者不依赖初始模型，便可获得较好的反演结果，但由于反演需要繁复计算偏导数矩阵和反复调用正演模型计算，反演速度极慢；后者则除了需要繁复计算偏导数矩阵，并严重依赖初始模型，初始模型的选择正确与否，直接影响最终的反演结果。因此，目前取得的研究成果还很少，主要处于研究阶段。

博斯蒂克(bostick)反演法是一种具有代表性的近似反演技术，尽管反演不够精确，但运算简便，能够直观的给出地下电阻率随深度的变化形式，所以得到广泛的应用。在大地电磁测深数据实时处理和现场处理的系统中大多配备了该反演程序。

本发明从博斯蒂克(bostick)反演的基本理论出发，利用大地电磁法TE极化模式和TM极化模式及其观测数据与地下电阻率和磁化率的相关关系和特点，提出了一种基于电阻率和磁化率曲线拟合分解改正，从大地电磁数据中反演电阻率和磁化率参数的方法及系统。

发明内容

本发明所要解决的技术问题是提供一种基于大地电磁数据反演电阻率和磁化率参数的方法及系统，进行电阻率和磁化率曲线拟合分解改正，用于解决现有技术中从大地电磁数据中反演电阻率和磁化率参数的方法存在的计算繁琐、精度不高的问题。

本发明解决上述技术问题的技术方案如下：一种基于大地电磁数据反演电阻率和磁化率的方法，包括：

步骤1，在TM极化模式下测量地面视电阻率实测值ρ_ai，并基于该视电阻率实测值ρ_ai建立初始地质模型，获得该初始地质模型各层的电阻率ρ_i和磁化率x_i；

步骤2，在水平层状介质条件下，基于初始地质模型的电阻率ρ_i和磁化率x_i，采用一维正演方法计算地面的视电阻率理论值ρ_ati；

步骤3，逐层修正初始地质模型的电阻率ρ_i，修正公式为：

${\mathrm{\ρ}}_i^k={\mathrm{\ρ}}_i^{k-1}+\mathrm{\α}\frac{\mathrm{\ξ}{|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k|}^2}{1+\mathrm{\ξ}|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k|};$

逐层修正初始地质模型的磁化率x_i，修正公式为：

$x_i^k=x_i^{k-1}+\mathrm{\β}\frac{|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k|}{{\mathrm{\ρ}}_{\mathrm{ati}}^k(1+\mathrm{\ξ}|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k\left|\right)};$

其中，i为地质模型的层数，k为模型修正次数；α和β为设定的步长；ξ为控制电阻率与磁化率所占比例的权参数，其取值范围为0至∞；

步骤4，基于修正后的各层电阻率和磁化率重新计算地面视电阻率理论值若地面视电阻率实测值ρ_ai和重新计算的地面视电阻率理论值的拟合差小于预先设定的误差期望值，则停止步骤3中对电阻率和磁化率的修正。

在上述技术方案的基础上，本发明还可以做如下改进。

进一步，基于视电阻率实测值ρ_ai进行Bostick反演方法，计算出电性层层数及各层的厚度参数和电阻率，同时取各层初始的磁化率x_i均为0，得到初始地质模型。

进一步，计算出的电性层中包括有一个虚拟电性层，且该虚拟电性层为所述初始地质模型的最底层电性层。

进一步，所述步骤3中，当ξ取值为0时，不需要进行电阻率的修正，且在测量地面的电磁场分量时仅需考虑磁化率的影响。

进一步，所述步骤3中，当ξ取值为无穷大时，不需要进行磁化率的修正，且在测量地面的电磁场分量时仅需考虑电阻率的影响。

本发明的技术方案还包括一种基于大地电磁数据反演电阻率和磁化率的系统，包括以下模块：

模型建立模块，用于在TM极化模式下测量地面视电阻率实测值ρ_ai，并基于该视电阻率实测值ρ_ai建立初始地质模型，获得该初始地质模型各层的电阻率ρ_i和磁化率x_i；

电阻率理论值计算模块，用于在水平层状介质条件下，基于初始地质模型的电阻率ρ_i和磁化率x_i，采用一维正演方法计算地面的视电阻率理论值ρ_ati；

修正模块，用于逐层修正初始地质模型的电阻率ρ_i，修正公式为：

${\mathrm{\ρ}}_i^k={\mathrm{\ρ}}_i^{k-1}+\mathrm{\α}\frac{\mathrm{\ξ}{|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k|}^2}{1+\mathrm{\ξ}|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k|};$

用于逐层修正初始地质模型的磁化率x_i，修正公式为：

$x_i^k=x_i^{k-1}+\mathrm{\β}\frac{|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k|}{{\mathrm{\ρ}}_{\mathrm{ati}}^k(1+\mathrm{\ξ}|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k\left|\right)};$

其中，i为地质模型的层数，k为模型修正次数；α和β为设定的步长；ξ为控制电阻率与磁化率所占比例的权参数，其取值范围为0至∞；

判断模块，用于基于修正后的各层电阻率和磁化率调用电阻率理论值计算模块重新计算地面的视电阻率理论值若地面视电阻率实测值ρ_ai和重新计算的地面视电阻率理论值的拟合差小于预先设定的误差期望值，则停止所述修正模块对电阻率和磁化率的修正。

进一步，所述模型建立模块包括：

反演模块，其用于根据视电阻率实测值ρ_ai进行Bostick反演方法，计算出电性层层数及各层的厚度参数和电阻率；

磁化率设置模块，其用于将各层初始的磁化率x_i取为0；

模型生成模块，其用于根据电性层层数、各层的厚度参数、电阻率及磁化率，生成初始地质模型。

进一步，所述反演模块计算出的电性层中包括有一个虚拟电性层，且该虚拟电性层为所述初始地质模型的最底层电性层。

进一步，所述修正模块中，当ξ取值为0时，不需要进行电阻率的修正，且在测量地面的电磁场分量时仅需考虑磁化率的影响。

进一步，所述修正模块中，当ξ取值为无穷大时，不需要进行磁化率的修正，且在测量地面的电磁场分量时仅需考虑电阻率的影响。

本发明的有益效果是：本发明的方法不用进行繁复偏导数矩阵计算，迭代过程稳定收敛，且可逐层判断地层的地电阻率参数和磁化率参数。本发明可用于地层、岩浆岩体和磁性矿体的勘探工作，通过区分地下岩层在构造变动和岩石蚀变产生的大地电磁响应，划分出多种有用的物理参数，进行隐伏矿体的深部空间定位。

附图说明

图1为本发明所述反演电阻率和磁化率的方法的流程示意图；

图2为进行大地电磁测量的装置及模式示意图；

图3为本发明实施例中的地质模型示意图；

图4为基于图3的地质模型的视电阻率-磁化率曲线图，图4(a)为A点TE模式下的曲线图，图4(b)为B点TE模式下的曲线图，图4(c)为A点TM模式下的曲线图，图4(d)为B点TM模式下的曲线图；

图5为应用例中实测获得的大地电磁数据曲线图；

图6为应用例中建立的初始地质模型的数据示意图；

图7为某沉积变质型铁矿分布区进行音频大地电磁测量给出的TM模式视电阻率拟断面图；

图8为沉积变质型铁矿分布区应用音频大地电磁TM模式观测数据，通过曲线拟合分解改正法反演获得的地层电阻率反演断面图；

图9为沉积变质型铁矿分布区应用音频大地电磁TM模式观测数据，通过曲线拟合分解改正法反演获得的地层磁化率反演断面图。

具体实施方式

以下结合附图对本发明的原理和特征进行描述，所举实例只用于解释本发明，并非用于限定本发明的范围。

如图1所示，本实施例给出了一种从大地电磁数据中反演电阻率和磁化率的方法，包括：

步骤1，在TM极化模式下测量由N个周期点T_i(i＝1,2,3.…N)在地面的视电阻率实测值ρ_ai，并基于该视电阻率实测值ρ_ai进行Bostick反演,获得N个电性层和一个虚拟电性层N+1(该虚拟电性层为反演需要而设立的最底层)，取第i层厚度和电阻率值分别为

$\begin{array}{c}{H}_i=\sqrt{\frac{{\mathrm{\ρ}}_{\mathrm{ai}}{T}_i}{2{\mathrm{\π\μ}}_0}}-\sqrt{\frac{{\mathrm{\ρ}}_{a(i-1)}{T}_{i-1}}{2{\mathrm{\π\μ}}_0}}\\ {\mathrm{\ρ}}_i={\mathrm{\ρ}}_{\mathrm{ai}}\frac{1+(\lg {\mathrm{\ρ}}_{\mathrm{ai}}-\lg {\mathrm{\ρ}}_{a(i-1)})/(\lg T_i-\lg T_{i-1})}{1-(\lg {\mathrm{\ρ}}_{\mathrm{ai}}-\lg {\mathrm{\ρ}}_{a(i-1)})/(\lg T_i-\lg T_{i-1})}\end{array}$

其中第1层和第N+1层的厚度和电阻率分别为

$H_1\sqrt{\frac{{\mathrm{\ρ}}_{a1}{T}_1}{{2\mathrm{\π\μ}}_0}}$

ρ＝ρ_a1

H_N+1＝∞

ρ_N+1＝ρ_N

需知，每个电性层(除最底层外)与频率一一对应。得到地下N+1层的各层厚度参数H_i及各层的电阻率参数ρ_i，同时令各层磁化率x_i均为零，从而获得初始的地质模型。

步骤2，在水平层状介质条件下，基于初始地质模型的电阻率ρ_i和磁化率x_i，采用一维正演方法计算地面的视电阻率理论值ρ_ati；

步骤3，逐层修正初始地质模型的电阻率ρ_i，修正公式为：

${\mathrm{\ρ}}_i^k={\mathrm{\ρ}}_i^{k-1}+\mathrm{\α}\frac{\mathrm{\ξ}{|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k|}^2}{1+\mathrm{\ξ}|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k|};---\left(1\right)$

逐层修正初始地质模型的磁化率x_i，修正公式为：

$x_i^k=x_i^{k-1}+\mathrm{\β}\frac{|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k|}{{\mathrm{\ρ}}_{\mathrm{ati}}^k(1+\mathrm{\ξ}|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k\left|\right)};---\left(2\right)$

其中，i为地质模型的层数(i＝1,2,3.…N)，k为模型修正次数(k＝1,2,…k)；α和β为设定的步长；ξ为控制电阻率与磁化率所占比例的权参数，其取值范围为0至∞；

步骤4，基于修正后的各层电阻率和磁化率重新计算地面的视电阻率理论值若地面视电阻率实测值ρ_ai和重新计算的地面视电阻率理论值的拟合差小于预先设定的误差期望值，则停止步骤3中对电阻率和磁化率的修正。

对应地，本实施例也给出了相应的反演电阻率和磁化率的系统，包括：

模型建立模块，用于在TM极化模式下测量地面的视电阻率实测值ρ_ai，并基于该视电阻率实测值ρ_ai建立初始地质模型，获得该初始地质模型的各层厚度参数H_i及各层电阻率ρ_i参数，同时令各层磁化率x_i参数均为零。

电阻率理论值计算模块，用于在水平层状介质条件下，基于初始地质模型的电阻率ρ_i和磁化率x_i，采用一维正演方法计算地面的视电阻率理论值ρ_ati。

修正模块，用于逐层修正初始地质模型的电阻率ρ_i，修正公式为上述的公式(1)；用于逐层修正初始地质模型的磁化率x_i，修正公式为上述的公式(2)。

判断模块，用于基于修正后的电阻率和磁化率调用电阻率理论值计算模块重新计算地面的视电阻率理论值若地面视电阻率实测值ρ_ai和重新计算的地面视电阻率理论值的拟合差小于设定的误差期望值，则停止修正模块对电阻率和磁化率的修正。

此外，所述模型建立模块包括：

反演模块，其用于根据视电阻率实测值ρ_ai进行Bostick反演方法，计算出电性层层数及各层的厚度参数和电阻率；

磁化率设置模块，其用于将各层初始的磁化率x_i取为0；

模型生成模块，其用于根据电性层层数、各层的厚度参数、电阻率及磁化率，生成初始地质模型。

上述各个模块的功能及具体实施过程与相应步骤一致，下面进一步阐述相应步骤具体的实施过程。

一、初始地质模型的建立

(一)在TM极化模式对大地电磁场进行测量。

本实施例以天然交变电磁场为场源，当交变电磁场以波的形式在地下介质中传播时，通过麦克斯韦方程组，可以把相应电磁场分为两组，一组包括分量E_y、H_x、H_z为H偏振波(TM极化)，另一组包括分量E_x、H_y、E_z为E偏振波(TE极化)，两组彼此独立。如图2所示，其为大地电磁法野外矢量测量模式，其中构造线方向沿y轴方向，Z_yx＝|E_y|/|H_x|为TM测量模式阻抗张量，Z_xy＝|E_x|/|H_y|为TE测量模式阻抗张量。

图3示意了一个地质模型，其为一个矩形低阻棱柱体，埋深为100m，规模为200m×200m，A测点在棱柱体上方，B点在棱柱体外，异常体电阻率为10Ωm，磁化率分别取0、0.01、0.05、0.5、1，围岩电阻率为100Ω.m，磁化率为0。图4(a)至图4(b)显示当该矩形低阻棱柱体具有不同磁化率时，引起的正演响应在测点A和B上的不同变化曲线。从图4可以看出，磁化率的变化对于TE模式的正演响应影响较小，即便磁化率值增大到1单位，对于正演响应也是没有明显变化。与TE模式相比，磁化率对TM模式的影响就相对比较明显。当磁化率等于1个单位时，个别频率上，相对于磁化率为零时，引起的视电阻率的偏差可达5％。另外与层状介质相似，只有当κ>0.1单位时，才会引起视电阻率有比较明显的变化，一般认为κ＝0.01单位时是磁化率产生变化的门限值。

因此，通过对大地电磁的理论分析和数值模拟结果分析，可知对于相同地质模型，磁化率的变化对大地电磁TM模式响应的影响明显大于TE模式，即TE观测模式受磁化率参数的影响很小，甚至可以忽略不计，而TM观测模式受磁化率参数的影响较大，其中磁化率门限值一般认为是0.01。从而，从TM模式观测数据中可以利用有效方法分解出磁性地层的磁化率信息，本实施例则主要是对大地电磁场在TM模式下进行测量，以分析电阻率和磁化率的变化特性。

(二)获得初始地质模型及对应参数

本实施例在地质表面测量的周期性输入函数E_y、H_x,在TM极化模式下测得视电阻率实测值并进一步获取-f频率曲线。基于-f频率曲线，本实施例采用Bostick反演方法，以视电阻率实测值为参数计算电性层厚度、电性层层数和电阻率，基于电性层厚度、电性层层数和电阻率建立初始地质模型，并令该初始地质模型的各层的初始磁化率为0。

另外，为便于后续描述，本文中采用ρ_ai表示视电阻率实测值采用分别表示第k次修正后的地质模型各层的电阻率、磁化率，当测量有N个周期点T_i(i＝1,2….N)地面视电阻率时，则初始地质模型共有N+1层，参考上述对于步骤1的描述，其中第1层和第N+1层的厚度和电阻率分别为

$H_1\sqrt{\frac{{\mathrm{\ρ}}_{a1}{T}_1}{{2\mathrm{\π\μ}}_0}}$

ρ＝ρ_a1

H_N+1＝∞

ρ_N+1＝ρ_N

以此类推，获得初始的地质模型。

二、视电阻率理论值的计算

在水平层状介质条件下，基于初始地质模型的电阻率ρ_i和磁化率x_i，采用一维正演方法计算地面的视电阻率理论值ρ_ati，具体的计算过程如下所述。

一维正演计算水平层状介质条件下地面理论观测数据或参数,由下列公式求解：

$Z\left(T\right)=\frac{E\left(T\right)}{H\left(T\right)}$

式中E(T)和H(T)分别是水平正交的电场强度和磁场强度的谱,其中T＝2π/ω，T为周期,ω为圆频率。则求取地面的视电阻率为:

${\mathrm{\ρ}}_a\left(T\right)=\frac{1}{\mathrm{\ω}{\mathrm{\μ}}_1}{\left|Z_{1,N}\left(T\right)\right|}^2$

Z_1,N表示N层介质情况下第一层顶面处的波阻抗，μ₁为第一层介质的磁导率，则递推公式为：

$Z_m\left(T\right)=Z_{0m}\frac{1-L_{m+1}\exp (-{2k}_m{H}_m)}{1+L_{m+1}\exp (-{2k}_m{H}_m)}$

式中m＝1,2,3,...,N-1,且：

$L_{m+1}=\frac{Z_{0m}-Z_{m+1}\left(T\right)}{Z_{0m}+Z_{m+1}\left(T\right)}$

其中Z_0m＝-iωμ_m/k_m为第m层的特征阻抗，μ₀为真空中的磁导率，μ_m＝μ₀(1+x_m)为第m层的磁导率，x_m为第m层的磁化率，为第m层的复波数，σ_m为第m层的电导率(电阻率ρ_m倒数)。其中，第N层顶面处的阻抗为

Z_N(T)＝Z_0N＝-iωμ_N/K_N

三、初始地质模型相关参数的修正

利用曲线拟合分解改正法进行反演计算，通过下式来修正初始地质模型各地层的电阻率值：

${\mathrm{\ρ}}_i^k={\mathrm{\ρ}}_i^{k-1}+\mathrm{\α\Δ}{\mathrm{\ρ}}_i={\mathrm{\ρ}}_i^{k-1}+\mathrm{\α}\frac{\mathrm{\ξ}{|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k|}^2}{1+\mathrm{\ξ}|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k|}$

其中为反演的电阻率，k为修正次数，ρ_ai为TM极化模式下的视电阻率实测值，为一维正演计算的视电阻率的理论值，α为步长,通常取一绝对值较小的数,若则α＞0，若则α＜0；ξ为控制电阻率与磁化率所占比例的权参数，其取值范围为0至∞，当ξ取值为0时，不需要进行电阻率的修正，且在测量地面的电磁场分量时仅需考虑磁化率的影响；当ξ取值为无穷大时，不需要进行磁化率的修正，且在测量地面的电磁场分量时仅需考虑电阻率的影响。

通过下式来修正反演的初始地质模型各层的磁化率值：

$x_i^k=x_i^{k-1}+\mathrm{\β\Δ}{x}_i=x_i^{k-1}+\mathrm{\β}\frac{|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k|}{{\mathrm{\ρ}}_{\mathrm{ati}}^k(1+\mathrm{\ξ}|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k\left|\right)}$

β为步长,通常取一绝对值较小的数。

四、判断是否满足误差期望值

基于修正后的电阻率磁化率采用一维正演方法重新计算地面的视电阻率理论值具体过程如第二部分所述。

事先给定一个拟合误差期望值ε₀，当视电阻率实测值与一维正演计算的视电阻率理论值的之差在所有周期点上总体相对误差ε(即拟合差)小于预先指定的一个很小的误差期望值ε₀时,即：

$\mathrm{\ϵ}=\sqrt{\frac{1}{n}\underset{1}{\overset{n}{\mathrm{\Σ}}}\left(\frac{{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k}{{\mathrm{\ρ}}_{\mathrm{ai}}}\right)}\≤{\mathrm{\ϵ}}_0,$

进行数据存盘，结束反演过程，否则继续步骤3的修正。

以构建的初始地质模型为基础，应用前述曲线拟合分解改正法反演。如图5所示，在一个含有磁性地质体二维地质单元中，矢量测量获得TE和TM两条视电阻率曲线，以TM视电阻率数据为基础，利用Bostick直接反演方法建立初始模型，并令模型各层初始磁化率为零。从图5可知，测量频率分别为：8192、4096、2048、1024、512、256、128、64、32、16、8、4、2、1,共14频点，初始模型为15层，各层厚度(1,2,3….N,N+1)分别为：23.27，7.197，8.866，12.848，23.272，42.287，71.06，99.59，114.57，138.41，203.26，381.17，689.5，1690.9，∞(米)；各层电阻率分别为：31，28，26，21，20，35，80,，120，93，30，50，68，85，83，83(欧姆米)；各层磁化率均为0。

反演结果如图6所示，拟合步长参数α＝0.5、β＝0.5，电阻率与磁化率的权参数ξ＝0.1，迭代20次，拟合误差小于0.5％。从图6可知，反演获得的地层电阻率和磁化率结果，除最底层(N+1)共14层，各层厚度不变，各层电阻率反演结果为：28，25，22，16，18，25，63，77，77，93，42，50，59，74(欧姆米)；各层磁化率为：0.06，0.08，0.11，0.24，0.15，0.07，0.06，0.05，0.14，0.35，0.15，0.11，0.09，0.06(SI单位)；最终的拟合误差小于0.5％。

实践中，应用对电阻率和磁化率的曲线拟合分解改正技术，在含磁性地层的地质与地球物理条件下开展各种实验。

举例说明，如图7所示，选择某沉积变质型铁矿分布区进行音频大地电磁测量，其中围岩主要为太古代片麻岩，目标地质体为沉积变质型铁矿，通过测量给出的TM模式视电阻率拟断面图，测线点距40m，共75个测点，采集频率从10HZ到100000HZ，共40个频点。

通过曲线拟合分解改正法对75个测深点单独进行反演，将各测点获得的地层电阻率反演结果绘制成断面图形式，分析断面图中的电阻率特征，表现为地下存在一个向形的电性构造特征，这与该地区沉积变质型铁矿的控矿构造特点相吻合，其效果如图8所示。

通过曲线拟合分解改正法对75个测深点单独进行反演，将各测点获得的地层磁化率反演结果绘制成断面图形式，分析断面图中的磁化率特征，表现为在地下200-400米空间内和600-1500米空间内存在两个中等磁化率异常带，推断为含一定品位的沉积变质型铁矿存在所引起，其效果如图9所示。

根据模型计算和现场试验，当地下存在含磁化率异常地质体时，应用本实施例提出的方法，所获得的电场与磁场观测数据包含有磁化率信息因素的影响，当忽略磁化率的影响时，反演得到的电阻率分布不符合实际目标体情况，通过曲线拟合分解改正法对观测数据进行一维反演，不用繁复计算偏导数矩阵，迭代过程稳定收敛，在获得地下地层电阻率参数和地层厚度参数的同时，获得地下地层磁化率参数。

以上所述仅为本发明的较佳实施例，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (10)

1.一种基于大地电磁数据反演电阻率和磁化率的方法，其特征在于，包括：

步骤1，在TM极化模式下测量地面视电阻率实测值ρ_ai，并基于该视电阻率实测值ρ_ai建立初始地质模型，获得该初始地质模型各层的电阻率ρ_i和磁化率x_i；

步骤2，在水平层状介质条件下，基于初始地质模型的电阻率ρ_i和磁化率x_i，采用一维正演方法计算地面的视电阻率理论值ρ_ati；

步骤3，逐层修正初始地质模型的电阻率ρ_i，修正公式为：

${\mathrm{\ρ}}_i^k={\mathrm{\ρ}}_i^{k-1}+\mathrm{\α}\frac{\mathrm{\ξ}{|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k|}^2}{1+\mathrm{\ξ}|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k|};$

逐层修正初始地质模型的磁化率x_i，修正公式为：

$x_i^k=x_i^{k-1}+\mathrm{\β}\frac{|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k|}{{\mathrm{\ρ}}_{\mathrm{ati}}^k(1+\mathrm{\ξ}|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k\left|\right)};$

其中，i为地质模型的层数，k为模型修正次数；α和β为设定的步长；ξ为控制电阻率与磁化率所占比例的权参数，其取值范围为0至∞；

步骤4，基于修正后的各层电阻率和磁化率重新计算地面视电阻率理论值若地面视电阻率实测值ρ_ai和重新计算的地面视电阻率理论值的拟合差小于预先设定的误差期望值，则停止步骤3中对电阻率和磁化率的修正。

2.根据权利要求1所述的方法，其特征在于，所述步骤1中建立初始地质模型具体包括：基于视电阻率实测值ρ_ai进行Bostick反演方法，计算出电性层层数及各层的厚度参数和电阻率，同时取各层初始的磁化率x_i均为0，得到初始地质模型。

3.根据权利要求2所述的方法，其特征在于，计算出的电性层中包括有一个虚拟电性层，且该虚拟电性层为所述初始地质模型的最底层电性层。

4.根据权利要求1所述的方法，其特征在于，所述步骤3中，当ξ取值为0时，不需要进行电阻率的修正，且在测量地面的电磁场分量时仅需考虑磁化率的影响。

5.根据权利要求1所述的方法，其特征在于，所述步骤3中，当ξ取值为无穷大时，不需要进行磁化率的修正，且在测量地面的电磁场分量时仅需考虑电阻率的影响。

6.一种基于大地电磁数据反演电阻率和磁化率的方法，其特征在于，包括：

模型建立模块，用于在TM极化模式下测量地面视电阻率实测值ρ_ai，并基于该视电阻率实测值ρ_ai建立初始地质模型，获得该初始地质模型各层的电阻率ρ_i和磁化率x_i；

电阻率理论值计算模块，用于在水平层状介质条件下，基于初始地质模型的电阻率ρ_i和磁化率x_i，采用一维正演方法计算地面的视电阻率理论值ρ_ati；

修正模块，用于逐层修正初始地质模型的电阻率ρ_i，修正公式为：

${\mathrm{\ρ}}_i^k={\mathrm{\ρ}}_i^{k-1}+\mathrm{\α}\frac{\mathrm{\ξ}{|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k|}^2}{1+\mathrm{\ξ}|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k|};$

用于逐层修正初始地质模型的磁化率x_i，修正公式为：

$x_i^k=x_i^{k-1}+\mathrm{\β}\frac{|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k|}{{\mathrm{\ρ}}_{\mathrm{ati}}^k(1+\mathrm{\ξ}|{\mathrm{\ρ}}_{\mathrm{ai}}-{\mathrm{\ρ}}_{\mathrm{ati}}^k\left|\right)};$

其中，i为地质模型的层数，k为模型修正次数；α和β为设定的步长；ξ为控制电阻率与磁化率所占比例的权参数，其取值范围为0至∞；

判断模块，用于基于修正后的各层电阻率和磁化率调用电阻率理论值计算模块重新计算地面的视电阻率理论值若地面视电阻率实测值ρ_ai和重新计算的地面视电阻率理论值的拟合差小于预先设定的误差期望值，则停止所述修正模块对电阻率和磁化率的修正。

7.根据权利要求6所述的系统，其特征在于，所述模型建立模块包括：

反演模块，其用于根据视电阻率实测值ρ_ai进行Bostick反演方法，计算出电性层层数及各层的厚度参数和电阻率；

磁化率设置模块，其用于将各层初始的磁化率x_i取为0；

模型生成模块，其用于根据电性层层数、各层的厚度参数、电阻率及磁化率，生成初始地质模型。

8.根据权利要求7所述的系统，其特征在于，所述反演模块计算出的电性层中包括有一个虚拟电性层，且该虚拟电性层为所述初始地质模型的最底层电性层。

9.根据权利要求6所述的系统，其特征在于，所述修正模块中，当ξ取值为0时，不需要进行电阻率的修正，且在测量地面的电磁场分量时仅需考虑磁化率的影响。

10.根据权利要求6所述的系统，其特征在于，所述修正模块中，当ξ取值为无穷大时，不需要进行磁化率的修正，且在测量地面的电磁场分量时仅需考虑电阻率的影响。