基于DBIM的瞬变电磁电导率反演方法技术

基于DBIM的瞬变电磁电导率反演方法，涉及地球物理勘探。包括以下步骤：1)读取观测数据；2)建立初始模型；3)更新模型参数；4)计算模型场值；5)计算误差；6)计算Frechet导数；7)计算更新量；8)判断收敛条件。从瞬变电磁的频谱信息出发，只要提取接收信号的准确频谱信息，就可以进行反演。基于DBIM方法建立迭代反演过程，最后反演的结果能很好与实际数据相吻合，可以大大提高瞬变电磁系统的计算速度和反演精度。不仅适用于半航空瞬变电磁系统，而且适用于全航空瞬变电磁系统。

【技术实现步骤摘要】

本专利技术涉及地球物理勘探，尤其是涉及一种基于DBIM的瞬变电磁电导率反演方法。
技术介绍
瞬变电磁法(TransientElectromagneticMethods)又称时间域电磁方法(TimeDomainElectromagneticMethods)，简称TDEM或TEM。航空电磁(AEM,AirborneElectromagnetics)法是在地面电磁法基础上发展起来的一种空中测量的电磁法，被广泛用来进行资源勘探，并且其非常适合探测一些复杂地形和人员难以达到的地区。其中半航空瞬变电磁(Semi-TransientElectromagneticMethods)是为了克服全航空瞬变电磁的探测深度有限所发展出的新型航空电磁法，其是在地面铺设几公里长度的电性源，然后利用无人机或者直升机在空中接收磁场，通过分析接收到磁场的信息，来反演出地下介质的情况。航空电磁法是基于岩石电性及磁性差异，利用电磁感应原理，以固定翼飞机或直升机等飞行器作为运载工具，实时地球物理探测的勘探方法，其方法具有高效、经济、适应性强等特点，能够广泛应用于地面大面积的矿产普查、水工环境普查、详查和精细测量等工作。
技术实现思路
本专利技术的目的在于提供一种在瞬变电磁系统中能准确反演地层电导率，在频率域中直接进行反演，不需要进行频时转换，在数据处理过程比较简单，只需要准确提取相应频率的频谱即可，不仅适用于半航空瞬变电磁系统，同样适用于全航空瞬变电磁系统的基于DBIM的瞬变电磁电导率反演方法。本专利技术包括以下步骤：1)读取观测数据，具体步骤如下：将实际接收机接收数据作为输入数据。在步骤1)中，所述输入数据包括电压数据或磁场数据。2)建立初始模型，具体步骤如下：根据已有的信息建立初始模型，初步定义下列初始参数：ε＝(ε1,ε2...εn)；σ＝(σ1,σ2...σn)；h＝(h1,h2...hn)其中，ε：为地下分层模型的介电常数矩阵，ε1,ε2...εn：为地下第一层、第二层…第n层模型的介电常数；σ：为地下分层模型的电导率矩阵，σ1,σ2...σn：为地下第一层、第二层…第n层模型的电导率；h：为地下分层模型的电导率矩阵，h1,h2...hn：为地下第一层、第二层…第n层模型的厚度；在步骤2)中，所述初始模型包括分层信息、每层的厚度、每层介质电导率、介电常数等信息。3)更新模型参数，具体步骤如下：为迭代过程中更新模型的参数过程，每一次迭代的前一次都会得到一个模型的更新量，这里是δσ＝(δσ1,δσ2...δσn)，其中δσ为模型电导率更新矩阵，δσ1,δσ2...δσn为地下第一层、第二层…第n层的电导率更新量；得到电导率更新量之后，有：σk+1＝σk+δσk+v其中，σk+1为第k+1次迭代所需要的电导率参数矩阵；σk为第k次迭代过程中所需要的迭代矩阵；δσk+1为第k次迭代过程中所得到的模型电导率更新矩阵，其中第一次迭代，模型的更新量为0；4)计算模型场值，具体步骤如下：得到模型后，对模型的场值进行计算，采用频率域的计算方法：E＝<GEJ；J>H＝<GHJ；J>其中，E和H表示电场和磁场；GEJ和GHJ为电流源电场并矢Green函数和磁场Green函数，<；>为点乘积分和；5)计算误差，具体步骤如下：得到模型的场值后，对实际测量的场值与模型计算的场值做差，得到两者之间的误差；6)计算Frechet导数，具体步骤如下：Frechet导数矩阵由以下方程求出：其中为对应于地下介质所产生的等效源，为由源J'所产生的场值的变化量，得到去谱域的解后对其进行Fourier逆变换得到其频域的解。如下为在y方向源的作用下所产生的磁场的x分量：δHx＝δHx1+δHx2+δHx3δHx为磁场的x方向分量，δHx1、δHx2、δHx3分别为其分量的第一项、第二项与第三项，具体表达式如下：上式中Jy表示y方向的源的强度；y-l与yl表示沿y方向源的起点与终点；y为接收点的y坐标；x'表示源的x坐标，x、y表示接收点的x和y坐标；表示为沿着第n层的上表面积分到n层的下表面，即第n层厚度的积分；irip(r)Vrip(r)表示位于z'的1V串联脉冲电流源(或1A并联脉冲电压源)在深度z处的电流和电压，下标p＝e,h分别表示横磁波与横电波；J0、J1表示第0类Bessel函数和第1类Bessel函数；其中，kx和ky分别为x方向和y方向的波数；δσn为地下第n层的电导率参数变化量。对于上述的场值变化δHx＝F·δσ，其中F即为地层参数的Frechet导数。7)计算更新量，具体步骤如下：对于上述地层模型，建立如下成本函数，其中，C表示成本函数(CostFunction)；δf表示场值的误差；F为Frechet导数矩阵；δσk+1为第k次迭代所得到的k+1次的电导率更新量；fobs为测量到的场值；σk为第k次迭代所得到的电导率矩阵；γ2表示正则化系数；||||2表示矩阵的二范数。为了使成本函数取得最小值，等效于解如下方程：其中，为F矩阵的转置，若其为复数矩阵，则表示为共轭转置。由上式可以得到每一次迭代过程中电导率的更新量δσk+1。8)判断收敛条件，具体步骤如下：若未达到收敛条件，则返回步骤3)更新数据的参数；若达到收敛条件，则结束反演过程。本专利技术应用于瞬变电磁快速反演，是利用瞬变电磁系统接收到的信号的频谱信息在频率域基于变形波恩迭代方法(DBIMDistortedBorniterativemethods)的一种快速反演方法。在瞬变电磁系统中，发射信号为双极性方波或者是半正弦信号，同时利用接收线圈接收一次场或者是二次场的时域信息，本专利技术可以使用任意发射信号，只要准确记录发射信号和接收信号，就可以准确的求得对应发射信号下的理论接收波形；对于时域的接收信号对其进行简单的去直流，去噪等操作后做FFT运算，得到其对应频率的频谱信息，就可以用来进行反演运算。本专利技术不仅适用于半航空瞬变电磁系统，而且适用于全航空瞬变电磁系统。本专利技术从瞬变电磁的频谱信息出发，只要提取接收信号的准确频谱信息，就可以用本专利技术进行反演。本专利技术基于DBIM方法建立迭代反演过程，最后反演的结果能很好与实际数据相吻合，从而可以大大提高瞬变电磁系统的计算速度和反演精度。附图说明图1为本专利技术的系统反演流程图。图2为瞬变电磁地下模型分层图。图3为本专利技术理论模型5层的反演效果图。图4为本专利技术理论模型10层的反演效果图。图5为本专利技术实际数据效果图。具体实施方式以下将结合附图对本专利技术的技术方案及其突出效果作进一步说明。参见图1，本专利技术的具体实施步骤如下：(1)读取观测数据将实际接收机接收数据作为输入数据，这些数据可以是电压数据，也可以是接收到的磁场数据。(2)建立初始模型可根据已有的信息建立一个简单的初始模型，初始模型需要包括分层信息，每层的厚度，每层介质电导率、介电常数等信息，参见图2，在图2中，Z1、Z2、…Zn-1为层界面的高度，线源可采用电流源。初步定义下列初始参数：ε＝(ε1,ε2...εn)；σ＝(σ1,σ2...σn)；h＝(h1,h2...hn)ε：为地下分层模型的介电常数矩阵，ε1,ε2...εn：为地下第一层、第二层…第n层模型的介电常数；σ：为地下分层模型的电导率本文档来自技高网...

【技术保护点】
基于DBIM的瞬变电磁电导率反演方法，其特征在于包括以下步骤：1)读取观测数据，具体步骤如下：将实际接收机接收数据作为输入数据；2)建立初始模型，具体步骤如下：根据已有的信息建立初始模型，初步定义下列初始参数：ε＝(ε1,ε2...εn)；σ＝(σ1,σ2...σn)；h＝(h1,h2...hn)其中，ε：为地下分层模型的介电常数矩阵，ε1,ε2...εn：为地下第一层、第二层…第n层模型的介电常数；σ：为地下分层模型的电导率矩阵，σ1,σ2...σn：为地下第一层、第二层…第n层模型的电导率；h：为地下分层模型的电导率矩阵，h1,h2...hn：为地下第一层、第二层…第n层模型的厚度；3)更新模型参数，具体步骤如下：为迭代过程中更新模型的参数过程，每一次迭代的前一次都会得到一个模型的更新量，这里是δσ＝(δσ1,δσ2...δσn)，其中δσ为模型电导率更新矩阵，δσ1,δσ2...δσn为地下第一层、第二层…第n层的电导率更新量；得到电导率更新量之后，有：σk+1＝σk+δσk+1其中，σk+1为第k+1次迭代所需要的电导率参数矩阵；σk为第k次迭代过程中所需要的迭代矩阵；δσk+1为第k次迭代过程中所得到的模型电导率更新矩阵，其中第一次迭代，模型的更新量为0；4)计算模型场值，具体步骤如下：得到模型后，对模型的场值进行计算，采用频率域的计算方法：E＝<GEJ；J>H＝<GHJ；J>其中，E和H表示电场和磁场；GEJ和GHJ为电流源电场并矢Green函数和磁场Green函数，<；>为点乘积分和；5)计算误差，具体步骤如下：得到模型的场值后，对实际测量的场值与模型计算的场值做差，得到两者之间的误差；6)计算Frechet导数，具体步骤如下：Frechet导数矩阵由以下方程求出：δE~=<GEJ;J′>]]>δH~=<GHJ;J′>]]>其中为对应于地下介质所产生的等效源，为由源J'所产生的场值的变化量，得到去谱域的解后对其进行Fourier逆变换得到其频域的解；如下为在y方向源的作用下所产生的磁场的x分量：δHx＝δHx1+δHx2+δHx3δHx为磁场的x方向分量，δHx1、δHx2、δHx3分别为其分量的第一项、第二项与第三项，具体表达式如下：δHx1=-Jy2πy-ylρl∫0∞∫zn-1zn(irie(r)Vsie(s)J1(kρρl)dkρdzδσn+Jy2πy-y-lρ-l∫0∞∫zn-1zn(irie(r)Vsie(s)J1(kρρ-l)dkρdzδσn]]>δHx2=Jy2πy-ylρl∫0∞∫zn-1zn(irihVsih)J0(kρρl)dkρdzδσn-Jy2πy-y-lρ-l∫0∞∫zn-1zn(irihVsih)J0(kρρ-l)dkρdzδσn-Jy2π∫-ll∫0∞∫zn-1znkρ(irihVsih)J0(kρρ)dkρδσndy′dz]]>δHx3=-Jy2πy-ylρl∫0∞∫zn-1zn(irihVsih)J1(kρρl)dkρdzδσn+Jy2πy-y-lρ-l∫0∞∫zn-1zn(irihVsih)J1(kρρ-l)dkρdzδσn]]>上式中Jy表示y方向的源的强度；y‑l与yl表示沿y方向源的起点与终点；y为接收点的y坐标；x'表示源的x坐标，x、y表示接收点的x和y坐标；表示为沿着第n层的上表面积分到n层的下表面，即第n层厚度的积分；irip(r)Vrip(r)表示位于z'的1V串联脉冲电流源(或1A并联脉冲电压源)在深度z处的电流和电压，下标p＝e,h分别表示横磁波与横电波；J0、J1表示第0类Bessel函数和第1类Bessel函数；其中，kx和k...

【技术特征摘要】
1.基于DBIM的瞬变电磁电导率反演方法，其特征在于包括以下步骤：1)读取观测数据，具体步骤如下：将实际接收机接收数据作为输入数据；2)建立初始模型，具体步骤如下：根据已有的信息建立初始模型，初步定义下列初始参数：ε＝(ε1,ε2...εn)；σ＝(σ1,σ2...σn)；h＝(h1,h2...hn)其中，ε：为地下分层模型的介电常数矩阵，ε1,ε2...εn：为地下第一层、第二层…第n层模型的介电常数；σ：为地下分层模型的电导率矩阵，σ1,σ2...σn：为地下第一层、第二层…第n层模型的电导率；h：为地下分层模型的电导率矩阵，h1,h2...hn：为地下第一层、第二层…第n层模型的厚度；3)更新模型参数，具体步骤如下：为迭代过程中更新模型的参数过程，每一次迭代的前一次都会得到一个模型的更新量，这里是δσ＝(δσ1,δσ2...δσn)，其中δσ为模型电导率更新矩阵，δσ1,δσ2...δσn为地下第一层、第二层…第n层的电导率更新量；得到电导率更新量之后，有：σk+1＝σk+δσk+1其中，σk+1为第k+1次迭代所需要的电导率参数矩阵；σk为第k次迭代过程中所需要的迭代矩阵；δσk+1为第k次迭代过程中所得到的模型电导率更新矩阵，其中第一次迭代，模型的更新量为0；4)计算模型场值，具体步骤如下：得到模型后，对模型的场值进行计算，采用频率域的计算方法：E＝<GEJ；J>H＝<GHJ；J>其中，E和H表示电场和磁场；GEJ和GHJ为电流源电场并矢Green函数和磁场Green函数，<；>为点乘积分和；5)计算误差，具体步骤如下：得到模型的场值后，对实际测量的场值与模型计算的场值做差，得到两者之间的误差；6)计算Frechet导数，具体步骤如下：Frechet导数矩阵由以下方程求出：δE~=<GEJ;J′>]]>δH~=<GHJ;J′>]]>其中为对应于地下介质所产生的等效源，为由源J'所产生的场值的变化量，得到去谱域的解后对其进行Fourier逆变换得到其频域的解；如下为在y方向源的作用下所产生的磁场的x分量：δHx＝δHx1+δHx2+δHx3δHx为磁场的x方向分量，δHx1、δHx2、δHx3分别为其分量的第一项、第二项与第三项，具体表达式如下：δHx1=-Jy2πy-ylρl∫0∞∫zn-1zn(irie(r)Vsie(s)J1(kρρl)dkρdzδσn+Jy2πy-y-lρ-l∫0∞∫zn-1zn(irie(r)Vsie(s)J1(kρρ-l)dkρdzδσn]]>δHx2=Jy2πy-ylρl∫0∞∫zn-1zn(irihVsih)J...

【专利技术属性】
技术研发人员：柳清伙，梁冰洋，邱晨，朱春晖，韩峰，刘娜，刘海，
申请(专利权)人：厦门大学，
类型：发明
国别省市：福建;35