CN103713317B - 基于对时变子波进行分频段处理的反褶积方法 - Google Patents

Abstract

本发明涉及一种提高地震记录分辨率的同时，压制噪声，以恢复反射系数为目标的反褶积方法，尤其是一种基于时变子波的分频段处理的反褶积方法，以解决常用反褶积方法提高分辨率的同时却降低信噪比，和求解不稳定的问题。该方法包括：采用井控多道统计自相关法提取地震记录中的时变子波；依据井资料提取高、低频算子，用于构造对应原始地震记录的高、低频信息，作为约束条件；对分频信息扩散滤波处理；原始地震记录及其高、低频信息联立反褶积。本发明使地震资料与测井资料有机的结合起来，避免了常规反褶积方法提取的不可控的反褶积因子，压制噪音的同时能有效提高分辨率，考虑了地震信号的衰减，实现了时域点谱白化的功能。

Description

基于对时变子波进行分频段处理的反褶积方法

技术领域

本发明涉及一种提高地震记录分辨率的同时，压制噪声，以恢复反射系数为目标的反褶积方法，尤其是一种基于对时变子波进行分频段处理的反褶积方法。

背景技术

提高地震资料分辨率一直是地震资料数字处理工作者的追求目标。高分辨率地震资料数据处理的关键环节就是压缩地震子波，或者去除地震波在地下传播过程中干涉、调谐等效应对地下地层的影响，拓宽有效地震信号的频带范围，特别是较为准确地拓宽高频成分。薄层的识别和厚度成为当前地震石油勘探领域的主要研究方向。薄层反射系数反演方法也就应运而生。

谱反演的原理就是根据时间域褶积模型，从地震记录中去除地震子波的影响，进而得到反射系数序列。在时间域内一个脉冲对的表达为：

g(t)=r₁δ(t-t₁)+r₂δ(t-t₁-T) (1)

其中r₁为层顶部反射系数，r₂为层底部反射系数，t为时间位置，t₁为顶部反射的时间位置，T为层厚度。

将分析点放在层的中心点位置，进行傅立叶变换后，用三角法则进行化简约去t，得到：

g(f)=2r_ecos(πfT)+2r_osin(πfT) (2)

基于褶积原理，s(t,f)是地震数据，w(t,f)是已知的子波。

$s(t,f)=w(t,f){\∫}_{-t_w}^{t_w}\left\{r_e\left(t\right)\cos \right[\mathrm{\πfT}\left(t\right)]+i{r}_o\left(t\right)\sin [\mathrm{\πfT}\left(t\right)\left]\right\}\mathrm{dt}---\left(3\right)$

最后整理得到单层目标函数：

$O(r_e,r_o,T,t)=$

$\underset{f_l}{\overset{f_H}{\∫}}\left[\begin{array}{c}{a}_e\left\{\mathrm{Re}\right[S(t,f)/W(t,f)]-\underset{-t_w}{\overset{t_w}{\∫}}{r}_e\left(t\right)\cos [\mathrm{\πfT}\left(t\right)\left]\mathrm{dt}\right\}+\\ a_o\left\{\mathrm{Im}\right[S(t,f)/W(t,f)]-\underset{t_w}{\overset{t_w}{\∫}}{r}_o\left(t\right)\sin [\mathrm{\πfT}\left(t\right)\left]\mathrm{dt}\right\}\end{array}\right]\mathrm{df}---\left(4\right)$

目标函数中，S(t,f)/W(t,f)可以认为是反褶积过程，在这里我们假设子波是已知。子波的正确与否将很大程度上影响反演的结果，因此子波的提取是一个非常重要的过程。但是我们对于子波往往是不知道的，如果子波知道了，反射系数也就知道了。因此我们做的是尽可能的使子波接近于真实子波。厚度T,我们在反射系数中可以认为是第一层和最后一层的距离，第二层和倒数第二层的距离，依次类推。在信号的奇偶分解中我们可以认为是进行奇偶分解的信号对之间的距离。

以上的推导都是基于窗口中心位置的坐标为零，但实际上数据都是从零开始的，而窗口的中心位置坐标不可能为零，而是半个窗口的位置。因此我们实际应用中应该对上面的推导加一时移。

实部和虚部时移后的谱可以表示为：

Im[e^2iπfΔtg(f)]=2r_osin(πfT)cos(2πfΔt)+2r_ecos(πfT)sin(2πfΔt) (5)

Re[e^2iπfΔtg(f)]=2r_ecos(πfT)cos(2πfΔt)-2r_osin(πfT)sin(2πfΔt) (6)

多层的目标函数变为：

$O(r_e,r_o,T,t)=$

$\underset{f_l}{\overset{f_H}{\∫}}\left[\begin{array}{c}{a}_e\left\{\mathrm{Re}\right[S(t,f)/W(t,f)]-\underset{0}{\overset{t_w}{\∫}}\{2r_e\left(t\right)\cos \left[\mathrm{\πfT}\left(t\right)\right]\cos 2\mathrm{\πf\Δt}+2r_o\left(t\right)\sin \left[\mathrm{\πfT}\left(t\right)\right]\sin 2\mathrm{\πf\Δt}\left\}\mathrm{dt}\right\}+\\ a_o\left\{\mathrm{Im}\right[S(t,f)/W(t,f)]-\underset{0}{\overset{t_w}{\∫}}\{2r_o\left(t\right)\sin \left[\mathrm{\πfT}\left(t\right)\right]\cos 2\mathrm{\πf\Δt}-2r_e\left(t\right)\cos \left[\mathrm{\πfT}\left(t\right)\right]\sin 2\mathrm{\πf\Δt}\left\}\mathrm{dt}\right\}\end{array}\right]\mathrm{df}---\left(7\right)$

其中，t_w为半个窗口的宽度；a_e，a_o为平衡参数。

求解该目标函数等价的矩阵形式为：

任何计算反射系数的反演方法实质都是一个反褶积方法，都是一个系统工程，许多小环节的处理手段都会影响到反演的结果，例如，子波提取、噪声及初始模型选取都对反演结果有影响。由于多层模型目标函数推导的假设条件，所有反射层位都是关于同一个时窗中心点对称的；但是，在实际情况下很难满足这样的条件，通常认为时窗内的每个点都有反射系数存在（实际上是稀疏的），才能应用目标函数公式（7）进行求解，反演出反射系数的大小。从理论上，在没有反射系数存在的地方，反演出的大小几乎趋于零。但是，实际问题的求解，时窗中心位置的选择很难判断，再者目标函数的解非常依赖于反演算法的全局优化性。直接求解目标函数公式（7）现阶段还不能取得很好的效果。通常是利用其他有效的方法去检测反射系数的位置，确定每层的厚度，再反演得到反射系数的大小。

总的来讲，此类反演方法的三个技术难点：

1、确定反射系数的初始位置；

2、噪声对反演结果的影响；

3、地震衰减对反演结果的影响。

发明内容

针对上述反演技术的不足之处，借鉴上述反演目标函数奇偶分解的思想，即对原始信号分解的思想，本发明提供一种可以提高地震记录分辨率的同时，压制噪声，以恢复反射系数为目标的基于对时变子波进行分频段处理的反褶积方法。

为实现上述目的，本发明提供一种基于对时变子波进行分频段处理的反褶积方法，包括以下步骤：

S1、采用井控下多道统计自相关法提取时变地震子波；

S2、根据井资料构造高、低频滤波因子，提取地震记录的高、低频信息；

S3、对步骤S2中的各频段信息做扩散滤波，去噪的同时增强同相轴的连续性，形成高信噪比和产状保持良好的高、中、低频三个数据体；

S4、对步骤S3中的高、中、低频数据体进行反褶积，以得到提频剖面。

在步骤S1中，采用井控下多道统计自相关法提取地震记录中的时变地震子波的具体步骤为：

（101）、提取井旁邻近道，多道均衡处理；

（102）、等间隔取时窗中心；

（103）、各时窗信号逐次处理；

（104）、各时窗子波提取结束，对子波做能量归一化处理，第一个时窗直接归一化，其余时窗，以第一个时窗为窗函数，加窗归一化；

（105）、各时窗子波光滑处理。

在步骤（101）中，将井旁邻近道资料与通过检测后得出的井资料进行对比，以确定井旁邻近道资料中时变地震子波是否准确。

在步骤（103）中，各时窗信号逐次处理的具体步骤为：

（a）、加窗处理；

（b）、求自相关；

（c）、该时窗信号处理，保证其中心点位最大振幅值；

（d）、带通滤波；

（e）、滤波后的信号变换到频率域，去噪处理；

（f）、反变换到时间域，再做一次带通滤波，得到该时窗的子波。

在步骤S2中，井控提取高、低频滤波因子P_L、P_H，构造高、低频信息S_L、S_H，

${\stackrel{\‾}{W}}^{-1}{S}_L=R_L,S_L=P_{L}S,R_L=P_{L}R,---\left(11\right)$

${\stackrel{\‾}{W}}^{-1}{S}_H=R_H,S_H=P_{H}S,R_H=P_{H}R.---\left(12\right)$

在步骤S3中，去除高、低频信息中的随机噪音，以形成高、中、低频三个数据体。

在步骤S4中，在去噪后的高、中、低频三个数据体进行反褶积，利用下式得出提频剖面，

$\min J={\left|\right|{\stackrel{\‾}{S}}_L-W_L{R}_L\left|\right|}_2+{\left|\right|\stackrel{\‾}{S}-\mathrm{WR}\left|\right|}_2+{\left|\right|{\stackrel{\‾}{S}}_H-W_H{R}_H\left|\right|}_2.---\left(13\right)$

与现有技术相比，本发明具有以下优点：

1、考虑了地震信号的能量衰减问题，并在井控下精细提取时变子波；

2、井控下提取分频滤波因子，并做分频去噪，有效压制了噪声的影响；同时有效地保持了高、低频信息能量特征，具有更高的时间分辨率和横向连续性；

3、在高、中、低频数据体的约束下的分频联合反褶积，增加了求解问题的稳定性，实现了时域点谱白化的功能。

附图说明

图1为本发明的方法流程图；

图2为地震道时频谱；

图3为由图2中提取的时变地震子波。

具体实施方式

下面结合附图和实施例对本发明作进一步详细说明。

如图1至图3所示，本发明提供一种基于对时变子波进行分频段处理的反褶积方法，包括以下步骤：

步骤S1、采用井控下多道统计自相关法提取地震记录中的时变地震子波的具体步骤为。

具体步骤为：

（101）、提取井旁邻近道，多道均衡处理；

（102）、等间隔取时窗中心；

（103）、各时窗信号逐次处理；

（104）、各时窗子波提取结束，对子波做能量归一化处理，第一个时窗直接归一化，其余时窗，以第一个时窗为窗函数，加窗归一化；

（105）、各时窗子波光滑处理。

在步骤（101）中，将井旁邻近道资料与通过检测后得出的井资料进行对比，以确定井旁邻近道资料中时变地震子波是否准确。

在步骤（103）中，各时窗信号逐次处理的具体步骤为：

（a）、采用以下公式进行加窗处理，

$w_i=\frac{1}{2}[1+\cos (-\mathrm{\π}+i*\frac{2\mathrm{\π}}{n-1})]$ i=0,1,…,n-1， (9)

其中，n为时窗信号长度；

（b）、求自相关；

（c）、该时窗信号处理，保证其中心点位最大振幅值；

（d）、采有以下公式进行带通滤波，

${\mathrm{flt}}_i=\frac{1}{2}[1+\cos ((i-n)*\frac{\mathrm{\π}}{N})]$ i=1,2,…,2n-1， (10)

其中，n为时窗长度，N为半子波长度；

（e）、滤波后的信号变换到频率域，去噪处理；

（f）、反变换到时间域，再做一次带通滤波，得到该时窗的子波。

地震记录是地震子波与反射系数褶积的结果，反褶积过程即可以看作是去子波的过程；而且地震数据由浅到深主频趋势也是由高到低变化的，因此为了保证反演结果的正确性，采取分时窗提取子波的方法，提取一系列时变子波。

步骤S2，根据井资料构造高、低频滤波因子，提取地震记录的高、低频信息。

具体的，井控提取高、低频滤波因子P_L、P_H，构造高、低频信息S_L、S_H，

${\stackrel{\‾}{W}}^{-1}{S}_L=R_L,S_L=P_{L}S,R_L=P_{L}R,---\left(11\right)$

${\stackrel{\‾}{W}}^{-1}{S}_H=R_H,S_H=P_{H}S,R_H=P_{H}R;---\left(12\right)$

由于时变地震子波带有不同频带区域的两种信号，而地震记录的频谱相对于反射系数的频谱，只是缺少一部分低频和高频成分的事实，因此，通过上式可增加反射系数的高频信号和低频信号比重。

步骤S3，去除步骤S2中各频段信息的噪音。

具体的，对各频段信息做扩散滤波，去噪的同时增强同相轴的连续性，形成高信噪比和产状保持良好的高、中、低频三个数据体。

S4、对步骤S3中的高、中、低频数据体进行反褶积，以得到提频剖面。

具体的，在对去噪后的高、中、低频数据体约束下做反褶积，利用下式求解提频剖面（近似反射系数剖面） ，

$\min J={\left|\right|{\stackrel{\‾}{S}}_L-W_L{R}_L\left|\right|}_2+{\left|\right|\stackrel{\‾}{S}-\mathrm{WR}\left|\right|}_2+{\left|\right|{\stackrel{\‾}{S}}_H-W_H{R}_H\left|\right|}_2.---\left(13\right)$

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

Claims (5)

1.一种基于对时变子波进行分频段处理的反褶积方法，其特征在于，包括以下步骤：

S1、采用井控下多道统计自相关法提取时变地震子波；

S2、根据井资料构造高、低频滤波因子，提取地震记录的高、低频信息；

S3、对步骤S2各频段信息进行扩散滤波，以得到高信噪比和增强同相轴连续性的高、中、低频数据体；

S4、对步骤S3中的高、中、低频数据体进行反褶积，以得到提频剖面。

2.根据权利要求1所述基于对时变子波进行分频段处理的反褶积方法，其特征在于，在步骤S1中，采用井控下多道统计自相关法提取地震记录中的时变地震子波的具体步骤为：

(101)、提取井旁邻近道，多道均衡处理；

(102)、等间隔取时窗中心；

(103)、各时窗信号逐次处理；

(104)、各时窗子波提取结束，对子波做能量归一化处理，第一个时窗直接归一化，其余时窗，以第一个时窗为窗函数，加窗归一化；

(105)、各时窗子波光滑处理。

3.根据权利要求2所述基于对时变子波进行分频段处理的反褶积方法，其特征在于，在步骤(101)中，将井旁邻近道资料与通过检测后得出的井资料进行对比，以确定井旁邻近道资料中时变地震子波是否准确。

4.根据权利要求3所述基于对时变子波进行分频段处理的反褶积方法，其特征在于，在步骤(103)中，各时窗信号逐次处理的具体步骤为：

(a)、加窗处理；

(b)、求自相关；

(c)、该时窗信号处理，保证其中心点位最大振幅值；

(d)、带通滤波；

(e)、滤波后的信号变换到频率域，去噪处理；

(f)、反变换到时间域，再做一次带通滤波，得到该时窗的子波。

5.根据权利要求1所述基于对时变子波进行分频段处理的反褶积方法，其特征在于，在步骤S3中，对各频段信息做扩散滤波，去噪的同时增强同相轴的连续性，形成高信噪比和产状保持良好的高、中、低频三个数据体。