CN115793064A - 一种改进的半航空瞬变电磁数据中激电信息的提取方法 - Google Patents

Abstract

本发明公开了一种改进的半航空瞬变电磁数据中激电信息的提取方法，包括：对反演的参数施加合理的对数域边界范围约束；从数据中截取早期段响应进行常规电阻率反演，为整体数据进行四参反演提供零频电阻率初始模型；反演过程中采取在前几次迭代中只反演零频电阻率和充电率，后续反演中再进行四个参数同时更新的策略；在进行后续的四参同时反演时，对时间常数和频率相关系数施加一个小幅度变化范围的约束；在四参反演时，令满足拟合差要求的测点不再进行反演更新来提高反演效率和稳定性。本发明具有逻辑简单、准确可靠等优点，在地球物理勘探技术领域具有很高的实用价值和推广价值。

Description

一种改进的半航空瞬变电磁数据中激电信息的提取方法

技术领域

本发明涉及地球物理勘探技术领域，尤其是一种改进的半航空瞬变电磁数据中激电信息的提取方法。

背景技术

半航空瞬变电磁法结合了地面发射源功率大、探测深以及空中飞行器快速测量、克服地形影响的优势，在矿产勘查、工程地质方面得到广泛的应用。而当地下存在极化体时，受激电效应影响，测得出现快速衰减和符号反转(负值)的异常响应。尤其近些年来，瞬变电磁响应中的激电效应已经被视为一种有效的寻找极化矿产的探测工具，如何从含激电效应的数据中准确的提取电阻率信息和激电信息是近些年该领域研究的热点和难点。

自上世纪80年以后，瞬变电磁响应中的激电效应开始广泛受到关注，早期通常对激电效应所引起的符号反转的晚期响应数据进行剔除，对剩余响应采用常规方法进行反演处理，这样的话必然会丢失大量的深部响应，导致反演结果不准确。现如今，普遍采用Cole-Cole电阻率模型对该现象进行模拟，但引入该模型也就意味着多个参加加入到反演当中，随之也带来了反演多解性更加严重的问题。为了改善多解性严重带来的反演不稳定问题，有的学者采取将极化响应去总响应中去除从而只反演电阻率的策略，该方法虽然在一定程度上可以提高电阻率反演的准确性，但不可避免的失去了地下极化信息，同时如何得到准确的极化响应也是一个难点。也有学者采取从含激电效应的数据中只反演零频电阻率和充电率的策略，这些方法虽然可以一定程度上能拟合异常响应，但时间常数和频率相关系数带来的影响就没有考虑，从而影响结果的准确性的同时，也失去了二者的地下分布信息。

另外，在“专利公开号为CN110673218A、名称为一种接地导线源瞬变电磁响应中IP信息的提取方法”的中国专利申请中，其包括：利用受IP效应影响较小的垂直磁场反演获取地下电阻率信息；基于获得的地下电性结构正演得到未受IP效应影响的电场响应；在观测的响应中去除IP效应的影响，得到纯IP响应；对获得的IP响应进行反演，获得极化率、频率相关系数和时间常数的IP信息。该技术的缺点在于，该技术的缺点在于垂直磁场依然受激电效应影响，导致反演的电阻率信息不准确，进而影响了后续的IP响应反演的准确性。

因此，急需要提出一种逻辑简单、提取可靠的改进的半航空瞬变电磁数据中激电信息的提取方法。

发明内容

针对上述问题，本发明的目的在于提供一种改进的半航空瞬变电磁数据中激电信息的提取方法，本发明采用的技术方案如下：

一种改进的半航空瞬变电磁数据中激电信息的提取方法，其包括以下步骤：

步骤S01，采集前T时刻的半航空瞬变电磁响应数据，并进行常规电阻率反演，得到初始的零频电阻率；所述T为大于0的自然数；

步骤S02，预设初始状态的充电率、时间常数、频率相关系数并和前T时刻的数据反演得到零频电阻率共同构成初始激电反演模型；

步骤S03，分别对零频电阻率、充电率、时间常数、频率相关系数进行对数域边界范围约束；

步骤S04，预设反演阈值K；所述K为大于0的自然数；预设初始状态的测点向量ψ为零向量；

步骤S05，若当前反演次数n小于反演阈值K时，则对充电率、零频电阻率进行反演；当反演次数n大于或等于阈值K时，对充电率、时间常数、频率相关系数和零频电阻率进行同步反演，同时对时间常数和频率相关系数施加一个小幅度变化范围的约束；

步骤S06，根据步骤S05中反演的充电率、时间常数、频率相关系数和零频电阻率，计算相应的雅可比矩阵；

同时，若当前的测点向量ψ为非零向量，则获取向量ψ中记录的测点位置所对应雅可比矩阵值置零；

步骤S07，根据反演计算结果更新横向约束反演模型；

步骤S08，求得更新后的横向约束反演模型的正演响应数据与实际数据的拟合差RMS和相邻次迭代拟合差相对变化ΔRMS；

步骤S09，若拟合差RMS小于5％，或相邻次迭代拟合差相对变化ΔRMS小于1％，则输出提取的激电信息；否则，则返回步骤S05，同时，求得横向约束反演模型中任一测点的拟合差RMS_ψ和相邻次迭代拟合差相对变化ΔRMS_ψ；若测点的拟合差RMS_ψ大于或等于5％，或相邻次迭代拟合差相对变化ΔRMS_ψ大于或等于1％，则不记录该测点位置信息，否则将该测点信息位置信息记录于测点向量ψ中。

进一步地，所述步骤S03中，分别对零频电阻率、充电率、时间常数、频率相关系数进行对数域边界范围约束，其对数域边界约束转换函数的表达式为：

其中，M为常规的参数表达，M_max和M_min为参数的上限和下限，

表示参数对数域的表达式。

更进一步地，当包括模拟激电效应时，引入Cole-Cole模型代替原本的电阻率模型模拟激电效应，其表达式为：

其中，

表示含激电效应的频散电阻率，ρ₀表示初始的零频电阻率；m₀表示初始的充电率；τ表示时间常数；c表示频率相关系数。

进一步地，所述步骤S07中，对横向约束反演模型施加约束条件，并与数据拟合项进行合并，在第(k+1)次迭代时，反演的表达式为：

其中，w_d表示数据加权因子；w_m表示模型加权因子；

表示参数对数域的表达式；

表示第k次迭代的雅可比矩阵；M_k+1表示第k+1次迭代的模型参数向量；M_k表示第k次迭代的模型参数向量；d_k表示第k次迭代的数据向量；d_obs表示实际数据向量；L表示模型横向光滑因子；e_obs表示数据项残差；e_r表示模型项残差。

进一步地，所述模型横向光滑因子的表达式为：

进一步地，所述拟合差RMS的表达式为：

其中，d_i′为第i^th个测点计算的模型响应；d_{obs_i}表示第i^th个测点的实际响应；N_S表示测线的测点个数；TN表示数据时间道数。

进一步地，对时间常数和频率相关系数施加一个相对于初始激电反演模型±50％的幅度变化范围。

(1)本发明巧妙地采用四个参数同时反演，使得四个参数同时对数据进行拟合，能够得到丰富的地下电阻率信息和极化信息。另外，本发明根据初始反演次数，采用充电率和零频电阻率反演，以解决多解性严重的问题，以避免四参反演的难度急剧增大的问题。

(2)本发明通过零频电阻率、充电率、时间常数和频率相关系数进行对数域边界范围约束，使得参数更新在一个合理的范围当中，避免无意义的解，保证信息的可靠性。

(3)本发明巧妙地对横向约束反演模型施加约束条件，并与数据拟合项进行合并，更进一步地减少反演的多解性。由于零频电阻率和充电率对于响应的灵敏度相对较高，因此，在反演前几次迭代中只反演这两个参数，使这两个参数能得到较好的拟合，在之后的反演中再引入灵敏度相对较低的时间常数和频率相关系数，来获得更多的地下极化信息。

(4)本发明巧妙地对时间常数和频率相关系数施加一个小变化范围约束，使得反演主要作用在灵敏度相对较高的零频电阻率和充电率上，让二者的结果更加理想，虽然时间常数和频率相关系数被约束在一个小范围内，但二者的异常值同样也能体现地下激电异常体的分布。

(5)本发明有效改善多参反演所带来的多解性严重的问题，不仅能够得到反映地下电阻率分布的电阻率信息，同时也能得到反映地下极化异常的极化信息，为后续地质解释增加了更多的地电信息。

(6)本发明采用拟合差和相邻次迭代拟合差相对变化进行判断；由于反演方式采用采用拟二维反演，即基于一维模型，对整个剖面数据进行同时反演，这样在整体反演过程中就可以“监视”每一个测点的拟合差变化，令拟合差值低于5％的或者相邻两次迭代拟合差相对变化低于1％的测点停止迭代，令其对应的雅可比矩阵置零。

综上所述，本发明具有逻辑简单、准确可靠等优点，并且能准确提取四个激电参数在地下的分布情况，在地球物理勘探技术领域具有很高的实用价值和推广价值。

附图说明

为了更清楚地说明本发明实施例的技术方案，下面将对实施例中所需使用的附图作简单介绍，应当理解，以下附图仅示出了本发明的某些实施例，因此不应被看作是对保护范围的限定，对于本领域技术人员来说，在不付出创造性劳动的前提下，还可以根据这些附图获得其他相关的附图。

图1为本发明的流程示意图。

图2为发明中模型一的参数及数据示意图。

图3为模型一的反演结果。

具体实施方式

为使本申请的目的、技术方案和优点更为清楚，下面结合附图和实施例对本发明作进一步说明，本发明的实施方式包括但不限于下列实施例。基于本申请中的实施例，本领域普通技术人员在没有做出创造性劳动的前提下所获得的所有其他实施例，都属于本申请保护的范围。

实施例

如图1至图3所示，本实施例提供了一种改进的半航空瞬变电磁数据中激电信息的提取方法，其包括以下步骤

第一步，采集1ms之前的半航空瞬变电磁相应数据，并进行常规电阻率反演，得到初始的零频电阻率。

第二步，预设充电率、时间常数、频率相关系数并上一步得到的零频电阻率结果共同构成激电反演的初始模型；

第三步，分别零频电阻率、充电率、时间常数、频率相关系数进行对数域边界范围约束。具体来说：

对数域边界约束转换函数的表达式为：

其中，M为常规的参数表达，M_max和M_min为参数的上限和下限，

表示参数对数域的表达式。在本实施例中，四个激电参数的取值范围：ρ₀＝10^-4～10⁴Ωm，m₀＝0～0.98，τ＝10^-3～5×10³s，c＝0～0.6，根据实际地质资料可进一步缩小四个激电参数的取值范围。

模拟激电效应时，引入Cole-Cole模型代替原本的电阻率模型模拟激电效应，其表达式为：

其中，

表示含激电效应的频散电阻率，ρ₀表示初始的零频电阻率；m₀表示初始的充电率；τ时间常数；c频率相关系数。

在本实施例中，为了进一步减少反演的多解性，对模型空间也施加约束条件，并与数据拟合项进行合并，在第k^th次迭代时，反演的表达式为：

Δl＝LM

其中，w_d表示数据加权因子；w_m表示模型加权因子；

表示参数对数域的表达式；

表示第k次迭代的雅可比矩阵；M_k+1表示第k+1次迭代的模型参数向量；M_k表示第k次迭代的模型参数向量；d_k表示第k次迭代的数据向量；d_obs表示实际数据向量；L表示模型横向光滑因子；e_obs表示数据项残差；e_r表示模型项残差。

模型横向光滑因子的表达式为：

第四步，预设反演阈值为5。

第五步，若当前反演次数n小于反演阈值5时，则对充电率、零频电阻率进行反演；否则，对时间常数和频率相关系数施加一个幅度变化范围的约束，并对充电率、时间常数、频率相关系数和零频电阻率进行同步反演。

第六步，根据反演参数计算相应参数的雅可比矩阵。同时读取向量ψ中记录的测点位置所对应雅可比矩阵值置零。

第七步，更新横向约束反演模型；

第八步，求得横向约束反演模型中拟合差RMS和相邻次迭代拟合差相对变化ΔRMS；

在本实施例中，拟合差RMS的表达式：

其中，d_i′为第i^th个测点计算的模型响应；d_{obs_i}表示第i^th个测点的实际响应；N_S表示测线的测点个数；TN表示数据时间道数。

第九步，若拟合差RMS大于5％，或相邻次迭代拟合差相对变化ΔRMS大于1％，则返回第五步，继续进行反演计算，同时求得横向约束反演模型中每个测点的拟合差RMS_ψ和相邻次迭代拟合差相对变化ΔRMS_ψ；若相应测点的拟合差RMS_ψ大于或等于5％，或相邻次迭代拟合差相对变化ΔRMS_ψ大于或等于1％，则不记录该测点位置信息，否则将该测点信息位置信息记录于ψ中。；否则，输出提取的激电信息。

实施例2

本实施例提供了一种改进的半航空瞬变电磁数据中激电信息的提取方法，其半航空瞬变电磁***采用1400米长的导线源，发射电流为20A。接收线圈有效面积为50平方米，飞行器飞行高度为30米，测线平行于长导线远，距离线源400米。

设计一个如图2(a)所示的低阻极化矿脉地质模型；图2(d)为该模型的具体模型参数；图2(b)(c)分别为该模型的无噪声的剖面数据和施加模拟噪声后的剖面数据。图2(e)测点400米测点处的响应，其中展示了不含噪声的模拟响应，加噪声后的模拟响应和施加的背景噪声分布。

对四个激电参数施加一个合理的约束范围：零频电阻率约束在[0,5000],充电率范围为[0,0.9],时间常数范围为[0,0.1]s,频率相关系数范围为[0,0.6]。

在本实施例中，对前1ms的相应数据进行常规电阻率反演，所得结果作为整体数据四个激电参数反演的初始零频电阻率模型，另外三个激电参数的初始模型均为均匀半空间模型，其中充电率的选值在0.1-0.3范围之间，时间常数选值为0.001s，频率相关系数选值为0.3。

反演过程中，采取在前几次迭代中固定时间常数和频率相关系数，只反演零频电阻率和充电率，在后续反演中再引入时间常数和频率相关系数，进行四个参数同时更新的策略，在此实施例中选择在前5次迭代中只反演零频电阻率和充电率，在5次以后的迭代中同时进行四个参数的反演。

在进行四参同时反演时，对时间常数和频率相关系数施加一个小幅度变化范围的约束。令时间常数约束在[0.0005，0.0015]s的变化范围内，频率相关系数约束在[0.25,0.35]的范围内。

在四参反演的同时令满足拟合差要求的测点不再进行反演更新来提高反演效率和稳定性。

如图3所示，采用传统一维阻尼最小二乘反演和本申请中的改进的反演提取方法对数据进行反演后，拟合差均达到了一个较低的水平，与传统一维阻尼最小二乘反演结果相比，采用改进后的反演提取方法后，对于零频电阻率，其反演效果得到明显的改善，中间低阻矿脉的形态得到有效恢复。对于充电率，时间常数和频率相关系数结果，使用改进方法后极大的提高了反演效果，对于地下激电异常的形态刻画与实际地电分布匹配，且很好的与实际低阻极化矿脉相对应，不仅提供了较为准确的地下电阻率分布信息，也得到了反映地下激电异常的极化信息。

上述实施例仅为本发明的优选实施例，并非对本发明保护范围的限制，但凡采用本发明的设计原理，以及在此基础上进行非创造性劳动而作出的变化，均应属于本发明的保护范围之内。

Claims (7)

1.一种改进的半航空瞬变电磁数据中激电信息的提取方法，其特征在于，包括以下步骤：

步骤S01，采集前T时刻的半航空瞬变电磁响应数据，并进行常规电阻率反演，得到初始的零频电阻率；所述T为大于0的自然数；

步骤S02，预设初始状态的充电率、时间常数、频率相关系数并和前T时刻的数据反演得到零频电阻率共同构成初始激电反演模型；

步骤S03，分别对零频电阻率、充电率、时间常数、频率相关系数进行对数域边界范围约束；

步骤S04，预设反演阈值K；所述K为大于0的自然数；预设初始状态的测点向量ψ为零向量；

步骤S05，若当前反演次数n小于反演阈值K时，则对充电率、零频电阻率进行反演；当反演次数n大于或等于阈值K时，对充电率、时间常数、频率相关系数和零频电阻率进行同步反演，同时对时间常数和频率相关系数施加一个小幅度变化范围的约束；

步骤S06，根据步骤S05中反演的充电率、时间常数、频率相关系数和零频电阻率，计算相应的雅可比矩阵；

同时，若当前的测点向量ψ为非零向量，则获取向量ψ中记录的测点位置所对应雅可比矩阵值置零；

步骤S07，根据反演计算结果更新横向约束反演模型；

步骤S08，求得更新后的横向约束反演模型的正演响应数据与实际数据的拟合差RMS和相邻次迭代拟合差相对变化ΔRMS；

步骤S09，若拟合差RMS小于5％，或相邻次迭代拟合差相对变化ΔRMS小于1％，则输出提取的激电信息；否则，则返回步骤S05，同时，求得横向约束反演模型中任一测点的拟合差RMS_ψ和相邻次迭代拟合差相对变化ΔRMS_ψ；若测点的拟合差RMS_ψ大于或等于5％，或相邻次迭代拟合差相对变化ΔRMS_ψ大于或等于1％，则不记录该测点位置信息，否则将该测点信息位置信息记录于测点向量ψ中。

2.根据权利要求1所述的一种改进的半航空瞬变电磁数据中激电信息的提取方法，其特征在于，所述步骤S03中，分别对零频电阻率、充电率、时间常数、频率相关系数进行对数域边界范围约束，其对数域边界约束转换函数的表达式为：

其中，M为常规的参数表达，M_max和M_min为参数的上限和下限，

表示参数对数域的表达式。

3.根据权利要求2所述的一种改进的半航空瞬变电磁数据中激电信息的提取方法，其特征在于，当包括模拟激电效应时，引入Cole-Cole模型代替原本的电阻率模型模拟激电效应，其表达式为：

其中，

表示含激电效应的频散电阻率，ρ₀表示初始的零频电阻率；m₀表示初始的充电率；τ表示时间常数；c表示频率相关系数。

4.根据权利要求1所述的一种改进的半航空瞬变电磁数据中激电信息的提取方法，其特征在于，所述步骤S07中，对横向约束反演模型施加约束条件，并与数据拟合项进行合并，在第(k+1)次迭代时，反演的表达式为：

其中，w_d表示数据加权因子；w_m表示模型加权因子；

表示参数对数域的表达式；

表示第k次迭代的雅可比矩阵；M_k+1表示第k+1次迭代的模型参数向量；M_k表示第k次迭代的模型参数向量；d_k表示第k次迭代的数据向量；d_obs表示实际数据向量；L表示模型横向光滑因子；e_obs表示数据项残差；e_r表示模型项残差。

5.根据权利要求4所述的一种改进的半航空瞬变电磁数据中激电信息的提取方法，其特征在于，所述模型横向光滑因子的表达式为：

6.据权利要求1所述的一种改进的半航空瞬变电磁数据中激电信息的提取方法，其特征在于，所述拟合差RMS的表达式为：

其中，d′_i为第i^th个测点计算的模型响应；d_{obs_i}表示第i^th个测点的实际响应；N_S表示测线的测点个数；TN表示数据时间道数。

7.据权利要求1所述的一种改进的半航空瞬变电磁数据中激电信息的提取方法，其特征在于，对时间常数和频率相关系数施加一个相对于初始激电反演模型±50％的幅度变化范围。

Images