CN105319594A - 一种基于最小二乘参数反演的傅里叶域地震数据重构方法 - Google Patents

Abstract

本发明提供了一种基于最小二乘参数反演的傅里叶域地震数据重构方法，属于油气及煤层气地震勘探与开发领域。所述方法基于最小二乘参数反演，对不含假频的非均匀采样的地震数据进行插值，从而重构地震数据。本发明方法不仅可以适用于不规则地震数据而且可以避免谱泄漏，并且计算效率很高。

Description

一种基于最小二乘参数反演的傅里叶域地震数据重构方法

技术领域

本发明属于油气及煤层气地震勘探与开发领域，具体涉及一种基于最小二乘参数反演的傅里叶域地震数据重构方法。

背景技术

地震勘探是利用地下介质弹性和密度的差异,通过观测和分析大地对人工激发地震波的响应,推断地下岩层的性质和形态的地球物理勘探方法。在这个过程中，地震波场是一个时间和空间都连续的多维信号。而实际记录的是这些连续信号的离散采样结果。根据信号采样理论，当采样频率大于Nyquist频率时，原始的连续信号可以由采样后的离散记录重构。但是在实际地震资料采集过程中，由于技术、经济以及地表条件的限制，对连续地震信号的空间采样并不是完美的。例如，在空间上，密集布置激发点(炮点)和接收点(检波点)不仅在技术和逻辑上无法实现，而且经济上也不允许。又例如，为了避开河流、断崖、村镇等因素时，激发点(炮点)和接收点(检波点)的布置就会变得不是很规则。而且在实际施工过程中还会有哑炮等情况。这些情况造成采集到的地震资料不仅在空间上采样不均匀，而且还可能产生假频。这会对地震资料后续的多道处理带来问题，例如基于多道的噪音压制、偏移、反演等。以用于岩性和流体识别的AVA分析为例，在的AVA分析中，为了使AVA技术可以用于复杂构造，需要利用叠前偏移来产生用于AVA分析的道集。波动方程偏移算法要求波场的空间采样是均匀的，Kirchhoff偏移虽然没有此要求，但是如果输入波场在空间采样是不均匀的，那么Kirchhoff偏移后的结果就不能保证相对振幅信息不被破坏，那么它产生的道集也就无法用于AVA分析。对于这个问题，通常的做法是面元规则化，但是面元规则化也只能使振幅相对保持，而不能解决空间方位上地震数据稀疏的问题。更好地解决这一问题的方法是，对地震数据进行插值或者外推，也就是重构。重构地震数据可以使用延拓算子进行，但是该方法需要准确的速度模型，如果速度模型不准确，则会得到错误的结果。另外一种方法是用空间预测滤波算子，但是这种方法和其他空间差值算法一样运算量很大。

发明内容

本发明的目的在于解决上述现有技术中存在的难题，提供一种基于最小二乘参数反演的傅里叶域地震数据重构方法，使得空间上采用不均匀的地震数据规则化，从而使得地震数据可以更好地满足后续地震资料处理和解释的需求。

本发明是通过以下技术方案实现的：

一种基于最小二乘参数反演的傅里叶域地震数据重构方法，基于最小二乘参数反演，对不含假频的非均匀采样的地震数据进行插值，从而重构地震数据。

所述方法包括：

(1)对原始地震数据即共炮道集进行预处理，得到NMO校正后的数据D；

(2)对于NMO校正后的数据D，计算：b＝A^H·W·D，

其中 $A=A_{\mathrm{nm}}=\frac{{\mathrm{\ΔS}}_F}{{4\mathrm{\π}}^2}{e}^{-i(m_{\mathrm{kx}}\·\mathrm{\Δ}{k}_x\·x_n+m_{\mathrm{ky}}\·\mathrm{\Δ}{k}_y\·y_n)},$ A^H是A的共轭转置矩阵，W是加权值，它是一个与局部采样密度成反比的值，m_kx，m_ky为空间波数采样序号，Δk_x，Δk_y为空间波数间隔，ΔS_F＝(X_i+1-X_i)(Y_i+1-Y_i)，(X_i，Y_i)、(X_i+1，Y_i+1)为空间上两点；

(3)计算：H＝A^H·W·A+λ²·I，然后求其逆H^-1，其中A、W与第(2)步中的相同，σ是先验模型变差，c是常数；I是单位对角矩阵；

(4)由第(2)、第(3)步得到的b和H^-1计算

(5)计算插值重构结果：其中A是与第(2)步中的相同，是由第(4)步计算所得。

所述步骤(1)中的预处理包括：置道头、去噪、振幅补偿、静校正、反褶积、分选得到CMP道集并对CMP道集进行NMO校正。

与现有技术相比，本发明的有益效果是：现有的方法是空间域进行的，例如邻近点插值或者多项式拟合。邻近点插值等空间域的插值方法，要求地震数据是规则采样的，缺少的只是个别点，而空间域拟合方法虽然不要求地震数据是规则的，并可以计算很多点，但是可能会出现假频并且计算效率不高，毕竟空间点的个数远远大于空间波数，而本发明基于最小二乘参数反演的傅里叶域地震数据重构方法不仅可以适用于不规则地震数据而且可以避免谱泄漏，并且计算效率很高。

附图说明

图1本发明方法的步骤框图。

具体实施方式

下面结合附图对本发明作进一步详细描述：

对于不含假频的均匀采样的数据D(p·Δx,q·Δy,ω)，对其进行空间2D离散Fourier变换：

$\tilde{D}(m_{\mathrm{kx}}\·\mathrm{\Δ}{k}_x,m_{\mathrm{ky}}\·\mathrm{\Δ}{k}_y,\mathrm{\ω})=\underset{p}{\mathrm{\Σ}}\underset{q}{\mathrm{\Σ}}D(p\·\mathrm{\Δx},q\·\mathrm{\Δy},\mathrm{\ω}){\·e}^{i(m_{\mathrm{kx}}\·\mathrm{\Δ}{k}_x\·p\·\mathrm{\Δx}+m_{\mathrm{ky}}\·\mathrm{\Δ}{k}_y\·q\·\mathrm{\Δy})}\·\mathrm{\Δx}\·\mathrm{\Δy}---\left(1\right)$

这样就可以得到然后利用进行离散Fourier反变换进行插值，得到空间任意点(x,y)上的值D(x,y,ω)：

$D(x,y,\mathrm{\ω})=\frac{{\mathrm{\Δk}}_x\·\mathrm{\Δ}{k}_y}{4{\mathrm{\π}}^2}\underset{m_{\mathrm{kx}}}{\mathrm{\Σ}}\underset{m_{\mathrm{ky}}}{\mathrm{\Σ}}\tilde{D}(m_{\mathrm{kx}}\·\mathrm{\Δ}{k}_x,m_{\mathrm{ky}}\·\mathrm{\Δ}{k}_y,\mathrm{\ω})e^{-i(m_{\mathrm{kx}}\·\mathrm{\Δ}{k}_x\·x+m_{\mathrm{ky}}\·\mathrm{\Δ}{k}_y\·y)}---\left(2\right)$

本发明借鉴上述思路，对于不含假频的非均匀采样的地震数据D(x_n,y_n,ω)进行如下处理，以得到

$\tilde{D}(m_{\mathrm{kx}}\·\mathrm{\Δ}{k}_x,m_{\mathrm{ky}}\·\mathrm{\Δ}{k}_y,\mathrm{\ω})=\underset{n}{\mathrm{\Σ}}D(x_n,y_n,\mathrm{\ω}){\·e}^{i(m_{\mathrm{kx}}\·\mathrm{\Δ}{k}_x\·x_n+m_{\mathrm{ky}}\·\mathrm{\Δ}{k}_y\·y_n)}\·\mathrm{\Δ}{S}_n---\left(3\right)$

在得到以后，再利用(2)插值得到任意点(x,y)或者规则网格点(pΔx,qΔy)上的数据。

其中，p，q为空间样点序号，Δx，Δy为空间采样间隔，ω为时间频率，Δk_x，Δk_y为空间波数间隔，m_kx，m_ky为空间波数采样序号。

但是，实际上对于不规则采样数据，利用(3)式进行分别对坐标x和y进行求和是不可行的。另外，假如非均匀采样的地震数据D(x_n,y_n,ω)是来自均匀采样的数据D(p·Δx,q·Δy,ω)，那么由(1)式和(3)式得到的Fourier谱应该是一致的，但是实际上并非如此，也就是说(3)式产生了谱泄漏，尤其是非均匀采样的地震数据D(x_n,y_n,ω)采样间隔比较大或者含有噪音的时候谱泄漏严重。为此，本发明采用最小二乘参数反演的方法来得到而不是直接利用(3)式进行的计算。具体方式如下：

对于(2)式，写成矩阵形式有：

对于含有噪音的情况，(4)式的更一般形式是：

那么利用最小二乘参数估计可以得到的一个估计值：

$\hat{\tilde{D}}={(A^H\·W\·A+{\mathrm{\λ}}^2\·I)}^{-1}{A}^H\·W\·D---\left(6\right)$

其中A^H是A的共轭转置矩阵。W是加权值，它是一个与局部采样密度成反比的值。σ是先验模型变差，c是常数。I是单位对角矩阵。

如图1所示，实现本发明的主要步骤如下：

①原始地震数据即共炮道集进行预处理：置道头、去噪、振幅补偿、静校正、反褶积、分选得到CMP道集并对CMP道集进行NMO校正，得到NMO校正后的数据D，即D(x_n,y_n,ω)，三个参数分别表示X坐标、Y坐标和频率)，这样做的目的是为了提高本发明的应用效果，也就是说本发明应用于NMO校正后的数据会效果更好；

②对于NMO校正后的数据D，计算：b＝A^H·W·D，

其中 $A=A_{\mathrm{nm}}=\frac{{\mathrm{\ΔS}}_F}{{4\mathrm{\π}}^2}{e}^{-i(m_{\mathrm{kx}}\·\mathrm{\Δ}{k}_x\·x_n+m_{\mathrm{ky}}\·\mathrm{\Δ}{k}_y\·y_n)},$ A^H是A的共轭转置矩阵，W是加权值，它是一个与局部采样密度成反比的值；

③计算：H＝A^H·W·A+λ²·I，然后求其逆H^-1，其中A、W与第②步中的相同，σ是先验模型变差，c是常数。I是单位对角矩阵；

④由②、③步分别得到b、和H^-1以后，

⑤计算插值重构结果：其中A是与第②步中的相同，是由第④步计算所得。

本发明借鉴利用不含假频的均匀采样数据和Fourier变换进行插值的思路，基于最小二乘参数反演，提出了一种利用不含假频的非均匀采样的地震数据进行插值的方法，从而实现地震数据规则化或者重构。

上述技术方案只是本发明的一种实施方式，对于本领域内的技术人员而言，在本发明公开了应用方法和原理的基础上，很容易做出各种类型的改进或变形，而不仅限于本发明上述具体实施方式所描述的方法，因此前面描述的方式只是优选的，而并不具有限制性的意义。

Claims (3)

1.一种基于最小二乘参数反演的傅里叶域地震数据重构方法，其特征在于：所述方法基于最小二乘参数反演，对不含假频的非均匀采样的地震数据进行插值，从而重构地震数据。

2.根据权利要求1所述的基于最小二乘参数反演的傅里叶域地震数据重构方法，其特征在于：所述方法包括：

(1)对原始地震数据即共炮道集进行预处理，得到NMO校正后的数据D；

(2)对于NMO校正后的数据D，计算：b＝A^H·W·D，

其中 $A=A_{\mathrm{nm}}=\frac{{\mathrm{\ΔS}}_F}{{4\mathrm{\π}}^2}{e}^{-i(m_{\mathrm{kx}}\·\mathrm{\Δ}{k}_x\·x_n+m_{\mathrm{ky}}\·\mathrm{\Δ}{k}_y\·y_n)},$ A^H是A的共轭转置矩阵，W是加权值，它是一个与局部采样密度成反比的值，m_kx，m_ky为空间波数采样序号，Δk_x，Δk_y为空间波数间隔，ΔS_F＝(X_i+1-X_i)(Y_i+1-Y_i)，(X_i，Y_i)、(X_i+1，Y_i+1)为空间上两点；

(3)计算：H＝A^H·W·A+λ²·I，然后求其逆H^-1，其中σ是先验模型变差，c是常数；I是单位对角矩阵；

(4)由第(2)、第(3)步得到的b和H^-1计算

(5)计算插值重构结果：

3.根据权利要求2所述的基于最小二乘参数反演的傅里叶域地震数据重构方法，其特征在于：所述步骤(1)中的预处理包括：置道头、去噪、振幅补偿、静校正、反褶积、分选得到CMP道集并对CMP道集进行NMO校正。