瞬变电磁信号叠前去噪方法与流程

本发明涉及地球物理信号处理与分析技术研究领域，尤其涉及一种瞬变电磁勘探的时间序列信号的叠前去噪技术。

背景技术：

长偏移距瞬变电磁阵列法(LOATEM)作为电磁法勘探中的一门新兴技术，愈来愈受到勘探地球物理界的重视。自上世纪90年代以来，国内学者开始从事该方法在石油勘探领域里的方法试验与解释方法研究，先后在南方碳酸盐岩地区和鄂尔多斯盆地成功地进行了方法试验，取得了一些令人鼓舞的结果，为该方法进一步发展奠定了基础。

尽管LOATEM法具有坚实的理论基础，尽管它在分辨率、信噪比等诸多因素方面优于频域电磁勘探方法，但在数据处理与解释方面仍存在一些十分突出的问题。首先，由于该方法记录的是时间序列，采样点密集(0.1至0.5毫秒采样间隔)，周期长(一般8至64秒)，阵列采集时数据量巨大，处理与反演任务相当繁重，现场一般不能及时给出电性剖面，从而影响了对采集数据质量以及勘探效果的评价。

随着工业化进程的加快，电磁法勘探即使在偏远的山区，其文化干扰依然存在，且到了无法忍受的地步，更不用说发达地区了。虽然LOATEM具有发射功率大、勘探深度深、抗干扰能力强、方便灵活等特点，但仍有许多测点由于受文化干扰的影响，信噪比很低，其数据往往无法解释。为了提高信解释质量，“去噪”是电磁法数据处理的前期基础工作。目前，去噪的方法很多。针对LOATEM 法时间采样密度高、单点时间序列记录数据量大、现场叠加处理困难等特点，本发明提出了一种瞬变电磁时间序列信号的叠前去噪技术。

瞬变电磁信号具有早期信号幅值高，衰减快，晚期信号弱，衰减慢，动态范围大等特点。常规的FFT方法因其主要是从全时段的频域特性分析出发，难以对瞬变电磁不同时窗的时间序列信号中的噪声较好地消除。也有学者研究与开发了三点指数逼近非线性平滑去噪和基于小波阈值的信号去噪技术，取得了一定的效果。目前瞬变电磁法提高信噪比的方法主要是通过叠加和加大发射源电流来。由于信噪比的提高与叠加次数的指数倍成正比，因而依靠叠加次数来提高信噪比有限。提高发射电流受仪器功率、器件及接地影响，也只能是一个有限的手段。开展信号识别与提取方法研究，在叠加与加大供电电流的前提下，开发新的去噪技术，是目前瞬变电磁法提高资料质量的重要途径。

技术实现要素：

本发明要解决的技术问题在于针对现有技术LOATEM法时间采样密度高、单点时间序列记录数据量大、现场叠加处理困难等的缺陷，提供一种可使得LOATEM方法资料处理快速有效，信噪比明显改善，资料畸变明显减小，资料品质明显提高，为视电阻率曲线计算，高分辨的精细反演奠定了坚实基础。

本发明解决其技术问题所采用的技术方案是：

提供一种瞬变电磁勘探的时间序列信号的叠前去噪方法，包括以下步骤：

S1、存贮瞬变电磁时间序列；

S2、对瞬变电磁测量的场分量时间序列的每个周期衰减曲线进行倒序进域递归滤波；

S3、对滤波后每个周期的场分量衰减曲线进行加窗处理，获取对数等间隔的场值衰减曲线；

S4、采取设定阈值的统计叠加处理方法，对对数等间隔的场值衰减曲线中的多个观测周期的同一时窗场值进行去飞点叠加。

本发明所述的方法中，步骤S1中，仅提取关断后的信号，且将每个周期内的两次判断后的二次场衰减曲线相加后再平均，形成一条衰减曲线，在一个测点进行K个周期观测时，则为K个周期衰减曲线的同一延迟时间的场值相加后平均。

本发明所述的方法中，步骤S3中进行加窗处理时，将关断后的时间分成N个对数等间隔时窗，对数域里的间隔宽度为Δ，每个时窗的中心时间点记为t₁,t₂......t_N，其所对应的场值分别计为EMF(t₁),EMF(t₂),.......EMF(t_N)，N为自然数；将每个时窗分成4份，依次记为t_i+Δ/4,t_i+Δ/2,t_i+Δ/2t_i+3Δ/4；为了计算中心时间点处的场分量值EMF(t_i)，将时窗宽度扩大到3Δ/2，以保证每个时窗有Δ/2的重叠；

将落入时窗中的K个采样点进行最小二乘拟合；设y是采样时间ts及参量的函数，(y_i,ts_i)是K对观测值，其中y_i是观测场值，ts_i是对应的采样时间，i＝1,......,K，K、M均为自然数；求参量使目标函数达到最小，其中σ_i为观测场值y_i误差的标准差；则时窗中心点所对应的场值为其中i＝1,......,N，N为自然数。

本发明所述的方法中，当时窗中观测场值点数少于3个时，采取中值或样条插值方法求取时窗中心点所对应的场值。

本发明所述的方法中，步骤S4中设瞬变电磁场分量衰减曲线的时窗为数为N，观测周期数为M，加窗后的第j个周期中的场值记为EMF^j(t_i)，其中i＝1,......,N；j＝1,......,M，M、N均为自然数；并设阈值为ε，则M个周期叠加处理如下：

④计算不同周期同一时窗中心点场值的平均值：

⑤设定阈值ε，踢除的场值，其中j＝1,......,M；

对余下的场值按①中公式重新计算时窗中心点场值。

本发明的有益效果：本发明使得LOATEM方法资料处理快速有效，信噪比明显改善，资料畸变明显减小，资料品质明显提高，为视电阻率曲线计算，高分辨的精细反演奠定了坚实基础。

附图说明

下面将结合附图及实施例对本发明作进一步说明，附图中：

图1是Z平面内极点与O点分布情况，其中+是极点，O是零点；

图2是陷波器的幅频特性；

图3是50Hz标准信号顺序与倒序滤波效果对比；

图4是实测衰减时间序列(工频干扰严重)；

图5是倒序和正序时域递归滤波效果对比；

图6是瞬变电磁资料的重叠加窗处理；

图7(a)为瞬变电磁磁资料叠前滤波视电阻率曲线图；

图7(b)为瞬变电磁磁资料叠后滤波视电阻率曲线图。

具体实施例

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅用以解释本发明，并不用于限定本发明。

叠前处理方法原理及流程

由于电偶极源瞬变电磁勘探方法具有勘探深度大，适应性强，其多分量观测时对高阻薄层有更高的分辨能力，在油气勘探中被广泛使用。目前，记录海量电偶极源瞬变电磁测深多分量时间序列成为可能，为增强时间序列资料予处理的信噪比提供了坚实的数据保证。

1.瞬变电磁时间序列存贮程

通常供电波型是双极性方波，占空比50％，周为T。在采集端接收到的信号中只提取关断后的信号。每个周期内的两次判断后的二次场衰减曲线加起来平均产，形成一条衰减曲线(注：每二次衰减曲线乘以负号与第一次衰减曲线相加除以2)。这样，在一个测点进行K个周期观测时，即是K个周期衰减曲线的同一延迟时间的场值加起来平均即可。对于叠前处理，不必保留时间序列数据，只是不断叠加，保留最后叠加结果即可。

2.叠前处理方法与技术

(1)叠前倒序时域递归滤波

①基本原理

设YUTEM采样时间段内记录的时间序列为{x_i}，i＝1,2,…,N，滤波因子为h_i,则滤波后的时间序列{y_i}为：

在实际运算中，滤波因子h_j只能取有限项，这就不可避免产生误差，为了减少误差，h_j要取成有限项但个数要相当多，YUTEM方式中最大项数N为4800，这样运算起来既占用内存大，又费机时，而递归滤波可以很好地解决这个问题。

递归滤波的思想是：认为输出值y_i之间彼此关联，因此在计算y_i时，要利用以前的计算结果y_i-1,y_i-2...。按照这种思想，常规的递归滤波公式为：

y_i＝a₀x_i+a₁x_i-1+…+a_nx_i-n-(b₁y_i-1+b₂y_i-2+…+b_my_i-m) (2)

其中m，n为自然数，a_j、b_j为递归滤波参数。

②Z平面法设计滤波器

LOATEM测量信号往往受到的是某些固定频率的干扰，因而滤波器的设计应是针对某个固定频率的陷波器。运用Z平面法可以方便地设计时域递归陷波器。所谓Z平面法设计滤波器，就是在Z平面上选择适当的0点和极点，按照对滤波器对振幅和相位的要求来设计滤波器。根据Z变换的基本理论，时域递归滤波器的幅频特性的Z变换可写为

Z平面法是利用单位园上的0点和其附近的极点对该点所对应幅频函数所具有陷波特性(如图1和图2所示)来快速确定W(z)，即在－f_r,f_r点处振幅值应为0，而在这两点之外幅值为1(归一化值)。要满足上述两个条件，则要求振幅函数必须同时含有0点和极点，且0点和极点应十分靠近。为此，(3)式可以写为

下面通过Z平面法来快速确定W(z)。

在Z平面内的零点与极点的复数形式可写为

z_z＝cosΩ_r±sinΩ_r＝R_z±iI_z (5)

z_p＝r_p cosΩ_p±r_p sinΩ_p＝R_p±iI_p (6)

其中r_p为滤波系数，Ω_r为要陷波的角频率，其数学表达式为：

Ω_r＝±180°·f_r/f_N

f_r为要陷波的频率，f_N为奈奎斯特频率，根据采样定理有：

f_N＝1/(2Δt)

将(5)、(6)式代入(4)式并通过简单的数学运算，则陷波器的频率响应可写为

G由奈奎斯特频率ω_N＝π时，令W(z)＝1而获取，即：

根据递归滤波器的性质，我们可以方便地写出时域递归滤波器的数学表达式：

⑥倒序滤波

在实际滤波过程中，首先将时间序列数据倒排，然后进行上述滤波处理，再重新倒排获取最后滤波结果，这样保证了时间序列前端数据的滤波质量，因为关断后的瞬变信号对资料处理与反演至关重要。以50Hz正弦波滤波为例，如图3所示。黑虚线为50Hz正弦波信号，红线为正序条件下用(8)式滤波效果，蓝线为倒序滤波效果。可以看出，正序滤波的前段有一段阻力衰减，信号明显畸变，效果较差，倒排滤波的前段数据滤波效果十分理想。基于倒序滤波技术，我们对实际数据进行了处理，处理效果如图4与图5所示，倒序滤波使前段时间序列信噪比大大提高，且保证了信号的真实性。

(2)重叠加窗处理

在每个周期的场分量衰减曲线进行滤波处理后，要求取对数域时窗曲线，形成对数域等间隔时窗数据。为了保证曲线连续、光滑，我们采取重叠加窗处理。其方法原理如图6所示。

将关断后的时间分成N个对数等间隔时窗，对数域里的间隔宽度为Δ，中心时间点记为t₁,t₂......t_N，其所对应的场值分别计为EMF(t₁),EMF(t₂),..EM.F(.t_N).。..将每个时窗分成4份，依次记为t_i+Δ/4,t_i+Δ/2,t_i+Δ/2t_i+3Δ/4。为了计算中心时间点处的场分量值EMF(t_i)，将时窗宽度扩大到3Δ/2，这样保证每个时窗有Δ/2的重叠，如图3所示。

将落入时窗中的K个采样点进行最小二乘拟合。设y是采样时间ts及参量的函数，(y_i,ts_i),(i＝1,......,K)是K对观测值，其中y_i是观测场值，ts_i是对应的采样时间。求参量使目标函数

达到最小，其中σ_i为观测场值y_i误差的标准差。则时窗中心点所对应的场值为

当时窗中观测场值点数少于3个时，可采取中值或样条插值方法求取时窗中心点所对应的场值。

(3)设定阈值的统计叠加处理

设瞬变电磁场分量衰减曲线的时窗为数为N，观测周期数为M，加窗后的第j个周期中的场值记为EMF^j(t_i),(i＝1,......,N；j＝1,......,M)，并设阈值为ε，则M个周期叠加处理按如下：

⑦计算不同周期同一时窗中心点场值的平均值：

⑧设定阈值ε，踢除的场值。

⑨对余下的场值按①式重新计算时窗中心点场值。

图7(a)和7(b)是单点叠前与叠后处理视电阻率曲线。可以看出，叠前滤波的视电阻率曲线品质明显提高。

本发明的大功率长偏移距瞬变电磁阵列勘探方法是一种用于深部勘探，特别是油气勘探中的电磁勘探方法。随着电子与计算机技术及仪器弱信号采集技术的迅猛发展，可同时海量存贮多分量多道观测场值的时间序列数据，为该方法的叠前信号处理提供了坚实的数据保证。本技术发明首先对瞬变电磁测量的场分量时间序列的每个周期衰减曲线进行倒序时域递归滤波，较好地消除了周期性的文化干扰，且具较高的保真性，其核心技术是倒序滤波。其次，进行对每个周期场分量衰减曲线加窗处理，获取对数等间隔的衰减曲线，其核心技术是时窗重叠与非线性最小二乘法求解时窗中心点场值。最后取设定阈值的多个周期统计叠加处理技术，对多个观测周期的同一时窗场分量值进行去飞点叠加处理，大大提高了场值估计的精度。

本发明使得LOATEM方法资料处理快速有效，信噪比明显改善，资料畸变明显减小，资料品质明显提高，为视电阻率曲线计算，高分辨的精细反演奠定了坚实基础。

应当理解的是，对本领域普通技术人员来说，可以根据上述说明加以改进或变换，而所有这些改进和变换都应属于本发明所附权利要求的保护范围。