具体实施方式
现将详细参照本发明的示例性实施例,所述实施例的示例在附图中示出,其中,相同的标号始终指的是相同的部件。以下将通过参照附图来说明所述实施例,以便解释本发明。
图1示出根据本发明示例性实施例的提取叠前校正道集地震子波的方法的流程图。这里,作为示例,所述方法可由用于提取叠前校正道集地震子波的设备来实现,也可以完全通过计算机程序来实现。
如图1所示,在步骤S100,对采集的地震记录和测井数据分别进行预处理,以获取地震叠前道集数据和测井道集数据。
这里,所述预处理是采用本领域的常规方法,例如将原始采集的地震记录和测井数据经过静校正、去噪、振幅补偿、动校正、叠加以及偏移等处理,最终形成地震叠前道集数据和测井道集数据,在此不再赘述。
在步骤S200,对地震叠前道集数据和测井道集数据分别进行非线性滤波来获取待校正地震叠前道集数据和待处理测井道集数据。
具体说来,由于提取高精度的地震子波对数据的信噪比和保真性要求较高,因此,为了获取精度高、抗噪能力强的地震子波,需要分别对地震叠前道集数据和测井道集数据进行非线性滤波,以提高地震叠前道集数据或测井道集数据的信噪比,为后续获取精度高、抗噪能力强的地震子波提供高信噪比的地震资料。
作为示例,可通过如下方式对地震叠前道集数据进行非线性滤波:
首先,可将地震叠前道集数据按道集顺序重采样到一维空间域,得到一维空间域的地震叠前道集数据,并利用预先设定的高通滤波器对所述一维空间域的地震叠前道集数据进行滤波处理,获得高频脉冲序列,例如,可按照下面的等式(2-1)对地震叠前道集数据进行滤波处理:
B(i)=A(i)*HF,等式(2-1)
其中,i表示地震叠前道集数据的道号,A(i)表示一维空间域的地震叠前道集数据,HF表示预先设定的高通滤波器,B(i)表示高频脉冲序列,*表示褶积。
然后,利用通过B(i)来获取脉冲干扰E(i),其中,脉冲干扰E(i)为|B(i)|的局部极大值所对应的B(i)的样值,接下来,从一维空间域的地震叠前道集数据A(i)去除脉冲干扰E(i),可进一步获取去除噪声后的一维空间域的地震叠前道集数据,并通过重采样将一维空间域的地震道集数据采样到二维空间域中,从而获取去除噪声后的待校正地震叠前道集数据。这里,应注意,为了更好地实现噪声压制,可对上述步骤进行多次迭代来提高地震叠前道集数据的信噪比,从而获取待校正地震叠前道集数据。此外,还可基于类似的方法来获取待处理测井道集数据。
图2示出根据本发明示例性实施例的非线性滤波前后的地震叠前道集数据的示意图。
图2中的(a)示出非线性滤波前的地震叠前道集数据的示图,图2中的(b)示出非线性滤波后的地震叠前道集数据的示图。这里,由图2中的(a)和(b)比较可以看出,经非线性滤波处理后的地震叠前道集数据(即,待校正地震叠前道集数据)信噪比较高,同相轴的连续性更好,因此,可通过对地震叠前道集数据或测井道集数据进行非线性滤波处理来为后续获取精度高、抗噪能力强的地震子波提供高信噪比的地震资料。
再次参照图1,在步骤S300,利用混合算法确定每道待校正地震叠前道集数据的相位校正角,其中,所述混合算法为最大方差模法和解析法。
具体说来,由于山地勘探中近地表降速带风化层的横向厚度和速度的变化,地下介质的非均一性和地下构造的复杂性等造成了相位特性的差异较大,使得动校正后道集质量受到了严重影响,因此,为了提高地震子波的质量,必须进行相位校正处理。而相位校正的关键是依据一定的判别准则求取相位校正角。下面,将结合图3来具体描述根据本发明示例性实施例的确定待校正地震叠前道集数据的相位校正角的步骤的流程图。
如图3所示,在步骤S310,可利用待校正地震叠前道集数据来确定待校正地震叠前道集数据的初始模型道。
具体说来,在步骤310,可利用下面的等式(3-1)来确定待校正地震叠前道集数据的初始模型道:
等式(3-1)
其中,m(t)表示第t时刻待校正地震叠前道集数据的初始模型道,Xi(t)表示第i道第t时刻的待校正地震叠前道集数据,N表示待校正地震叠前道集数据的道集个数且为正整数。
接下来,在步骤S320,可基于所述初始模型道,利用最大方差模法来获取待校正地震叠前道集数据的参考模型道,其中,所述参考模型道是指相位为零时的初始模型道。
作为示例,可通过如下方式获取待校正地震叠前道集数据的参考模型道:
这里,可利用降次幂的最大方差模准则获取关于所述参考模型道的目标函数:
其中,P表示关于参考模型道的目标函数,i、t和N的定义可参照上述示例,z(α,t)表示待校正地震叠前道集数据的参考模型道,其中,z(α,t)=m(t)cosα-H[m(t)]sinα,其中,α为未知常数,m(t)表示第t时刻待校正地震叠前道集数据的初始模型道,H[m(t)]表示第t时刻参考模型道的希尔伯特数据。
这里,使得P取得极小值时的常数α为初始模型道的相位校正角,在确定了初始模型道的相位校正角后,可进一步获取待校正地震叠前道集数据的参考模型道z(α,t)。
接下来,在步骤S330,针对待校正地震叠前道集数据中的每道数据,基于参考模型道,可利用解析法来确定所述每道数据的相位校正角。
作为示例,可通过如下方式获取每道地震叠前道集数据的相位校正角:
这里,在待校正地震叠前道集数据的参考模型道已知的情况下,可利用解析法获取关于所述相位校正角的目标函数:
其中,Q表示关于第i道待校正地震叠前道集数据的相位校正角的目标函数,M为第i道待校正地震叠前道集数据中地震数据的传播时间,i、t、Xi(t)、H[Xi(t)]和z(α,t)的定义可参照上述示例。
这里,使得Q取得极大值时的βi为第i道待校正地震叠前道集数据的相位校正角。
再次参照图1,在步骤S400,利用所述相位校正角分别对相应道的待校正地震叠前道集数据进行相位校正,以获取待处理地震叠前道集数据。通过所述处理,能够提高动校正后地震叠前道集的质量。
在步骤S500,基于待处理地震叠前道集数据和待处理测井道集数据,利用各道待处理地震叠前道集数据的预定匹配滤波算子来确定待处理地震叠前道集数据的最终匹配滤波算子。
图4示出根据本发明示例性实施例的确定待处理地震叠前道集数据的最终匹配滤波算子的步骤的流程图。
具体说来,如图4所示,在步骤S510,将每道待处理地震叠前道集数据分别与相应的预定匹配滤波算子进行匹配滤波,使得匹配滤波后的每道待处理地震叠前道集数据的实际输出趋近于期望输出的相应道的待处理测井道集数据,并确定匹配滤波后每道数据的实际输出与期望输出的误差。
具体说来,在步骤S510,可通过下面的等式(5-1)确定匹配滤波后每道数据的实际输出与期望输出的误差:
ei=Yi(t)*pi(t)-Li(t),等式(5-1)
其中,ei表示匹配滤波后第i道数据的实际输出与期望输出的误差,Yi(t)表示第i道第t时刻的待处理地震叠前道集数据(即,实际输出),Li(t)表示第i道第t时刻的待处理测井道集数据(即,期望输出),pi(t)表示第i道第t时刻的待校正地震叠前道集数据的预定匹配滤波算子。
在步骤S520,可利用每道数据的误差来获取每道数据的误差能量,基于所述误差能量应用最小二乘法来确定每道数据的匹配滤波算子。
这里,可通过下面的等式(5-2)确定每道数据的误差能量:
等式(5-2)
其中,Ei表示匹配滤波后第i道数据的实际输出与期望输出的误差能量,Yi(t)表示第i道第t时刻的待处理地震叠前道集数据,Li(t)表示第i道第t时刻的待处理测井道集数据,pi(t)表示第i道第t时刻的待处理地震叠前道集数据的预定匹配滤波算子。
其中,应用最小二乘法原理,可令误差能量对所述预定滤波算子的偏导数等于零,从而获得该道的匹配滤波算子。
在步骤S530,将需要进行处理的各道数据的匹配滤波算子进行平均处理,以得到待处理地震叠前道集数据的最终匹配滤波算子。这里,作为示例,可根据各道待处理叠前道集数据的相关性和信噪比来确定需要进行处理的匹配滤波算子。例如,可选择相关性好或信噪比高的待处理叠前道集数据的匹配滤波算子求平均,以获得待处理地震叠前道集数据的最终匹配滤波算子。
具体说来,在步骤S530,可通过下面的等式(5-3)确定待处理地震叠前道集数据的最终匹配滤波算子:
等式(5-3)
其中,p(t)表示待处理地震叠前道集数据的最终匹配滤波算子,i、t、N和pi(t)的定义可参照上述示例。
再次返回图1,在步骤S600,将最终匹配滤波算子与待处理地震叠前道集数据进行匹配滤波,以获取待处理地震叠前道集数据和待处理测井道集数据的匹配道集数据。
具体说来,可通过下面的等式(6-1)获取待处理地震叠前道集数据和待处理测井道集数据的匹配道集数据:
yi(t)=Yi(t)*p(t),等式(6-1)
其中,yi(t)表示待处理地震叠前道集数据和待处理测井道集数据的匹配道集数据的匹配道集,Yi(t)表示第i道第t时刻的待处理地震叠前道集数据样点值,p(t)表示最终匹配滤波算子。
图5示出根据本发明示例性实施例的待处理叠前道集数据和待处理测井道集数据进行匹配滤波前后的示意图。
图5中的(a)示出未进行匹配滤波前的待处理地震叠前道集数据和待处理测井道集数据的示图,图5中的(b)示出进行匹配滤波后的待处理地震叠前道集数据和待处理测井道集数据的示图。这里,由图5中的(a)可以看出,待处理地震叠前道集数据和待处理测井道集数据在匹配滤波前,存在明显的拼接痕迹(如图中虚线框中所示),频率不同且同相轴间存在时差;由图5中的(b)可以看出,在待处理地震叠前道集数据和待处理测井道集数据进过匹配滤波后,不仅消除了同相轴间的时差,而且衔接处的振幅、频率和相位也得到了较好的匹配,此外,深浅层的反射波数据也能较好地拼接。
通过上述处理,可实现待处理地震叠前道集数据和待处理测井道集数据的无缝匹配。
接下来,再次参照图1,在步骤S700,基于匹配道集数据和待处理测井道集数据来确定地震叠前道集数据的地震子波频谱,并利用所述地震子波频谱来获取最优地震子波。
图6示出根据本发明示例性实施例的获取最优地震子波的步骤的流程图。
如图6所示,在步骤S710,利用所述匹配道集数据和待处理测井道集数据来确定地震叠前道集数据的地震子波频谱。
下面,将结合图7来详细说明如何确定地震叠前道集数据的地震子波频谱。
如图7所示,在步骤S711,利用所述匹配道集数据和待处理测井道集数据确定地震叠前道集数据的初始子波频谱和谱一致函数。
具体说来,为了获取所述地震子波频谱,可利用所述匹配道集数据或待处理测井道集数据的相邻道集数据的相似性,来获取与地震子波频谱相关的量,例如,可基于所述匹配道集数据获取其信噪比和功率谱信息,同样,基于类似的方法,可基于待处理测井道集数据获取其信噪比和功率谱信息。
作为示例,可通过如下的方式确定地震叠前道集数据的初始子波频谱和谱一致函数:
首先,匹配道集数据基于褶积原理可被表示为:
y(t)=r(t)*w(t),等式(7-1)
在等式(7-1)中,y(t)表示匹配道集输出记录,r(t)表示地层反射系数序列(输入记录),w(t)表示地震子波。这里,作为示例,所述地震子波可作为一个滤波器。
这里,由输入输出的线性褶积模型原理,观测到的输入输出可被表示为:
input(t)=r(t)+ε(t),等式(7-2)
output(t)=y(t)+η(t)
其中,input(t)表示观测的输入,output(t)表示观测的输出,r(t)和y(t)的定义可参照上述示例,ε(t)表示输入噪声,η(t)表示输出噪声。
对等式(7-2)两边求取自相关,再做傅里叶变换:
Sinput(f)=Sr(f)+Sε(f),等式(7-3)
Soutput(f)=Sy(f)+Sη(f)
其中,Sinput(f)表示观测输入input(t)的功率谱,Soutput(f)为观测输出的output(t)的功率谱,Sr(f)表示有效信号的功率谱,Sy(f)表示匹配道集的功率谱,Sε(f)表示输入噪声的功率谱,Sη(f)表示输出噪声的功率谱。
这里,由于输入噪声和输出噪声是相互独立的高斯白噪,因此,可有:
Sinput,output(f)=Sry(f),等式(7-4)
其中,Sinput,output(f)为观测输入和观测输出的互功率谱,Sry(f)表示输入记录和输出记录的互功率谱。
这里,由于待处理叠前地震道集数据在每个频点处都有噪声,因此,地震子波频谱可被表示为:
等式(7-5)
其中,W(f)表示子波频谱,SNRy表示匹配道集的信噪比,表示地震子波频谱的无偏项,Sinput,output(f)为观测输入和观测输出的互功率谱,Sinput(f)表示观测输入input(t)的功率谱。
这里,由于观测系统的单位冲激响Winput,output(f)可应用来获取因此,通过等式(7-5)可知,只需通过使每个频率与无偏项相乘来修正观测系统,而不需要常规方法的最小相位假设。因此,通过这种方式,可以提高地震子波的可靠性。
另一方面,定义谱一致函数为:
等式(7-6)
其中,表示所述匹配道集的谱一致函数,SNRl表示待处理测井数据的信噪比,SNRy表示匹配道集数据的信噪比。
在步骤S712,利用所述谱一致函数对所述初始子波频谱进行无偏项替代转换,来获取地震子波频谱。
具体说来,可将无偏项用谱一致函数来表示,并代入等式(7-5),即得到地震子波频谱的表达式:
其中,W(f)表示地震子波频谱,表示匹配道集的谱一致函数,SNRl表示测井叠前道集的信噪比,SInput,Output(f)表示观测输入input(t)和观测输出output(t)的互功率谱,Sinput(f)表示观测输入input(t)的功率谱。
再次参照图6,在步骤S720,将所述地震子波频谱进行傅里叶反变换,得到时间域的地震叠前道集地震子波。
在步骤S730,将时间域的地震叠前道集地震子波进行加窗平滑处理来获取最优地震子波。
具体说来,可将时间域的地震叠前道集地震子波进行常规的中值或均值滤波处理来获取最优地震子波。
图8示出根据本发明示例性实施例的最优地震子波的示意图。
由图8可以看出,通过本方法获取的地震子波波形变化自然、稳定、主瓣能量集中,旁瓣震荡较少,且随着时间轴快速衰减。
综上所述,在根据本发明示例性实施例的提取叠前校正道集地震子波的方法中,能够提高地震子波的精度和抗噪能力,并且不需要最小相位的假设,提高了地震子波的可靠性,为后续反褶积和反演处理打下了良好的基础。
尽管已经参照其示例性实施例具体显示和描述了本发明,但是本领域的技术人员应该理解,在不脱离权利要求所限定的本发明的精神和范围的情况下,可以对其进行形式和细节上的各种改变。