CN105182407A - 一种三维地震数据道内插方法 - Google Patents

Abstract

本发明提供一种三维地震数据道内插方法，该方法包括：沿着地表二维平面的第一方向和第二方向，对预处理后的三维地震数据排序；将排序后的三维地震数据布散于三维空间中，并对排序后的三维地震数据划分若干窗口，使得相邻的窗口具有重叠部分；选定第一方向作为褶积方向，沿时间方向顺序遍历各个窗口，并对遍历到的每一窗口进行傅里叶变换、褶积计算、反傅立叶变换等处理，从而得到所有窗口的第一窗内时间-空间插值数据，其中需要对重叠部分的第一窗内时间-空间插值数据进行均值处理；最终输出插值结果。本发明是将二维f-x域预测扩展到三维f-x,y域，不受倾角限制和空间假频的影响，充分利用三维数据信息，内插后的信噪比和精度都获得提高。

Description

一种三维地震数据道内插方法

技术领域

本发明涉及地震勘探数据处理领域，具体地，涉及一种三维地震数据道内插方法。

背景技术

随着油气勘探的地质目标越来越复杂，对提高成像精度和纵横向分辨率的要求越来越高。近年来高密度地震勘探技术发展较快，其显著点是面元尺度小、面元属性均匀、空间采样密和接收道数多。但高密度地震勘探增加了炮密度和道密度，明显增加了地震采集成本，经济技术可行性是制约其广泛应用的主要原因。另外，常规地震勘探的老资料炮线距通常倍于检波线距，空间采样密度不够。老资料挖潜增效，也需要新的技术支撑，以改善成像的质量和品质。

地震道内插处理技术可以用来加密空间采样，降低采集成本。虽然不能增加信息量，但可以克服空间假频，提高偏移成像质量。地震道内插是以往叠后室内增加空间采样，提高成像精度和分辨率的主要方法。目前，地震道内插方法主要有：f-x域预测道内插法、f-k域谱扩展法、多项式拟合法、最小平方倾角分解法、以及稀疏反演傅里叶重建法等。其中f-x域预测道内插法基于线性同相轴在域可预测的理论，利用相移因子中f与道间时差t的对偶关系实现道内插，不需要相干同相轴的倾角信息，可以正确内插具有空间假频的地震道，处理效果较好。

发明内容

但是，现有的f-x域预测道内插法主要为二维地震数据道内插，对于三维地震数据应用二维道内插算法，应用的是局部二维数据而不是局部三维数据，会影响三维地震数据道内插的质量，通常是内插精度低、效果不稳定。

为此，非常需要一种三维地震数据频率-空间域道内插精度高、效果稳定的方法，能适应三维地震数据插值，增加空间采样的密度，并用于改善地震成像的质量和品质。

本发明实施例的主要目的在于提供一种三维地震数据道内插方法，以解决现有的三维地震数据道内插精度低、效果不稳定的问题。

为了实现上述目的，本发明实施例提供一种三维地震数据道内插方法，包括：

步骤1，采集三维地震数据；

步骤2，对所述三维地震数据进行预处理；

步骤3，沿着地表二维平面的相互正交的第一方向和第二方向，对预处理后的所述三维地震数据排序；其中，所述第一方向为所述三维地震数据的共中心点号的分布方向，所述第二方向为所述三维地震数据的线号的分布方向；或，所述第一方向为所述三维地震数据的线号的分布方向，所述第二方向为所述三维地震数据的共中心点号的分布方向；

步骤4，以时间方向以及所述第一方向和第二方向为轴建立三维空间，将排序后的所述三维地震数据布散于所述三维空间中，在所述三维空间中对排序后的所述三维地震数据划分若干三维的窗口，使得每个窗口对应于排序后的所述三维地震数据的一部分，各个窗口对应的数据量相等，且相邻的窗口具有重叠部分；

步骤5，选定所述第一方向作为褶积方向，沿时间方向顺序遍历各个窗口，并对遍历到的每一窗口执行如下处理：

步骤51，按照如下公式，将当前窗口的三维地震数据沿所述时间方向做傅里叶变换，得到第一频率-空间域数据：

S_i,j(f)＝FFT(S_i,j(t))

其中，i、j分别表示所述三维地震数据在所述第一方向、所述第二方向上的编号；

t表示所述三维地震数据在所述时间方向上对应的时间；

f表示频率；

S_i,j(f)表示对应编号i、j的第一频率-空间域数据；

S_i,j(t)表示对应编号i、j的三维地震数据；

步骤52，按照如下公式，在所述第一方向的正向进行褶积，计算所述第一方向的预测插值算子：

$\underset{x,y}{\Σ}{S}_{i+1-2x,j-1+y}\left(f\right)p_{x,y}\left(f\right)=S_{i,j}\left(f\right)$

其中，x＝1,2,…,lx；

y＝1,2,…,ly；

p_x,y(f)表示对应频率f的所述第一方向的预测插值算子；

lx表示预测插值算子p_x,y(f)在所述第一方向上的长度；

ly表示预测插值算子p_x,y(f)在所述第二方向上的长度；

S_i+1-2x,j-1+y(f)表示对应编号i+1-2x_、j-1+y的第一频率-空间域数据；

S_i,j(f)表示对应编号i、j的第一频率-空间域数据；

步骤53，按照如下公式，将所述第一方向的预测插值算子以及所述第一频率-空间域数据在所述第一方向的正向进行褶积，得到第一窗内频率-空间域正向预测的插值数据：

${G^{\′}}_{i,j}\left(f\right)=\underset{x,y}{\Σ}{S}_{i+1-x,j-1+y}\left(f\right)p_{x,y}(f/2)$

其中，p_x,y(f/2)表示对应频率f/2的所述第一方向的预测插值算子；

S_i+1-x,j-1+y(f)表示对应编号i+1-x、j-1+y的第一频率-空间域数据；

G'_i,j(f)表示对应编号i、j的第一窗内频率-空间域正向预测的插值数据；

步骤54，计算所述第一方向的预测插值算子的共轭，并按照如下公式将该共轭与所述第一频率-空间域数据在所述第一方向的反向进行褶积，得到第一窗内频率-空间域反向预测的插值数据：

${G^{\′\′}}_{i,j}\left(f\right)=\underset{x,y}{\Σ}{S}_{i-1+x,j+1-y}\left(f\right)p{*}_{x,y}(f/2)$

其中，p*_x,y(f/2)为对应频率f/2的第一方向的预测插值算子p_x,y(f/2)的共轭；

S_i-1+x,j+1-y(f)表示对应编号i-1+x、j+1-y的第一频率-空间域数据；

G”_i,j(f)表示对应编号i、j的第一窗内频率-空间域反向预测的插值数据；

步骤55，按照如下公式，对所述第一窗内频率-空间域正向预测的插值数据和所述第一窗内频率-空间域反向预测的插值数据求平均，得到第一窗内频率-空间域插值数据：

G_i,j(f)＝(G”_i,j(f)+G'_i,j(f))/2

其中，G_i,j(f)表示对应编号i、j的第一窗内频率-空间域插值数据；

步骤56，按照如下公式，对所述第一窗内频率-空间域插值数据进行反傅立叶变换，得到第一窗内时间-空间插值数据：

G_i,j(t)＝FFT^-1(G_i,j(f))；

其中，G_i,j(t)表示对应编号i、j的第一窗内时间-空间插值数据；

步骤57，针对当前窗口与前一窗口的重叠部分，将该重叠部分在当前窗口内对应的第一窗内时间-空间插值数据与该重叠部分在前一窗口内对应的第一窗内时间-空间插值数据进行均值处理，并将该重叠部分在当前窗口内对应的第一窗内时间-空间插值数据更新为该均值处理的结果；

步骤6，将所有窗口的三维地震数据和第一窗内时间-空间插值数据作为所述第一方向的插值结果进行输出。

借助于上述技术方案，本发明在计算预测插值算子的过程中使用的是局部三维地震数据而不是二维地震数据，更能适应局部三维的构造变化，并且由于信息量的增加，抗干扰能力也得到增强。相比于现有技术，本发明是将二维f-x域预测扩展到三维f-x,y域，不受倾角限制和空间假频的影响，充分利用三维数据信息，内插后的信噪比和精度都获得提高。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1是本发明示例性方法仅在第一方向进行插值的流程示意图；

图2是本发明示例性方法在第一方向和第二方向均进行插值的流程示意图；

图3是本发明实施例内插前的三维地震数据剖面；

图4是本发明实施例内插后的三维地震数据剖面。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

本领域技术技术人员知道，本发明的实施方式可以实现为一种系统、装置、设备、方法或计算机程序产品。因此，本公开可以具体实现为以下形式，即：完全的硬件、完全的软件(包括固件、驻留软件、微代码等)，或者硬件和软件结合的形式。

根据本发明的实施方式，提出了一种三维地震数据道内插方法。

在本文中，附图中的任何元素数量均用于示例而非限制，以及任何命名都仅用于区分，而不具有任何限制含义。

示例性方法

本发明提供一种示例性的三维地震数据道内插方法，如图1所示，该方法包括：

步骤S1，采集三维地震数据。

步骤S2，对三维地震数据进行预处理。

具体实施时，可选地，对三维地震数据进行包括加载观测系统、静校正、叠前去噪、振幅恢复、动校正和叠加等预处理。

步骤S3，沿着地表二维平面的相互正交的第一方向和第二方向，对预处理后的三维地震数据排序，其中，第一方向为三维地震数据的共中心点号的分布方向，第二方向为三维地震数据的线号的分布方向，或，第一方向为三维地震数据的线号的分布方向，第二方向为三维地震数据的共中心点号的分布方向。

根据排序结果，分别在第一方向和第二方向上对三维地震数据编号，例如，排序后的三维地震数据在第一方向上具有编号i，在第二方向上具有编号j。

该步骤的目的是在地表二维平面上对三维地震数据进行排序，使得排序后的三维地震数据满足三维频率-空间域道的规则化要求。

具体地，该步骤可以包括如下所示的过程：

S31，确定三维地震数据对应的三维地震记录道的地下反射点在地表二维平面的垂直投影点。

S32，根据垂直投影点沿第一方向和第二方向的分布情况，对三维地震记录道排序。

S33，根据三维地震记录道的排序情况，对预处理后的三维地震数据排序。

步骤S4，以时间方向、第一方向和第二方向为轴建立三维空间，将排序后的三维地震数据布散于三维空间中，在该三维空间中对排序后的三维地震数据划分若干三维的窗口，使得每个窗口对应于排序后的三维地震数据的一部分，各个窗口对应的数据量相等，且相邻的窗口具有重叠部分(即具有重叠的数据)。

该步骤是在时间(时间方向)和空间上(地表二维平面的第一方向和第二方向)将三维地震数据划分成多组三维地震数据，每组三维地震数据对应一个窗口。

具体地，该步骤可以按照如下的标准，对排序后的所述三维地震数据划分窗口：每个窗口沿时间方向的跨度为第一设定值(例如为2秒、3秒等)、沿第一方向的跨度为第二设定值(例如为15道、20道等)、沿第二方向的跨度为第三设定值(例如为15道、20道等)、并且相邻的窗口具有预设比例的重叠(重叠部分的数据量占到窗口所含全部数据量的一预设比例，该预设比例大于0且小于100％)。

本发明中，之所以令相邻的窗口具有重叠，是因为利用重叠部分的地震道可获得更好的内插效果，而如果相邻的窗口之间不重叠，则内插结果会出现不连续的情况。

步骤S5，选定第一方向作为褶积方向，然后沿时间方向顺序遍历所有窗口，并对遍历到的每一窗口执行如下步骤S501～步骤S507的处理：

步骤S501，将当前窗口的三维地震数据沿时间方向做傅里叶变换，得到第一频率-空间域数据：

S_i,j(f)＝FFT(S_i,j(t))

其中，i、j分别表示三维地震数据在第一方向、第二方向上的编号；

t表示三维地震数据在时间方向上对应的时间；

f表示频率；

S_i,j(f)表示对应编号i、j的第一频率-空间域数据；

S_i,j(t)表示对应编号i、j的三维地震数据。

步骤S502，按照如下公式，在第一方向的正向进行褶积，计算第一方向的预测插值算子：

$\underset{x,y}{\Σ}{S}_{i+1-2x,j-1+y}\left(f\right)p_{x,y}\left(f\right)=S_{i,j}\left(f\right)$

其中，x＝1,2,…,lx；

y＝1,2,…,ly；

p_x,y(f)表示对应频率f的第一方向的预测插值算子；

lx表示预测插值算子p_x,y(f)在第一方向上的长度；

ly表示预测插值算子p_x,y(f)在第二方向上的长度；

S_i+1-2x,j-1+y(f)表示对应编号i+1-2x、j-1+y的第一频率-空间域数据；

S_i,j(f)表示对应编号i、j的第一频率-空间域数据。

步骤S503，按照如下公式，将对应频率f/2的第一方向的预测插值算子p_x,y(f/2)与对应编号i+1-x_、j-1+y的第一频率-空间域数据S_i+1-x,j-1+y(f)在第一方向的正向进行褶积，得到第一窗内频率-空间域正向预测的插值数据：

${G^{\′}}_{i,j}\left(f\right)=\underset{x,y}{\Σ}{S}_{i+1-x,j-1+y}\left(f\right)p_{x,y}(f/2)$

其中，p_x,y(f/2)表示对应频率f/2的第一方向的预测插值算子；

S_i+1-x,j-1+y(f)表示对应编号i+1-x、j-1+y的第一频率-空间域数据；

G'_i,j(f)表示对应编号i、j的第一窗内频率-空间域正向预测的插值数据。

步骤S504，计算第一方向的对应频率f/2的预测插值算子p_x,y(f/2)的共轭，并将该共轭与对应编号i-1+x、j+1-y的第一频率-空间域数据S_i-1+x,j+1-y(f)在第一方向的反向进行褶积，得到第一窗内频率-空间域反向预测的插值数据：

${G^{\′\′}}_{i,j}\left(f\right)=\underset{x,y}{\Σ}{S}_{i-1+x,j+1-y}\left(f\right)p{*}_{x,y}(f/2)$

其中，p*_x,y(f/2)为对应频率f/2的第一方向的预测插值算子p_x,y(f/2)的共轭；

S_i-1+x,j+1-y(f)表示对应编号i-1+x、j+1-y的第一频率-空间域数据；

G”_i,j(f)表示对应编号i、j的第一窗内频率-空间域反向预测的插值数据。

步骤S505，按照如下公式，对对应编号i、j的第一窗内频率-空间域正向预测的插值数据和第一窗内频率-空间域反向预测的插值数据求平均，即得到对应编号i、j的第一窗内频率-空间域插值数据：

G_i,j(f)＝(G”_i,j(f)+G'_i,j(f))/2

其中，G_i,j(f)表示对应编号i、j的第一窗内频率-空间域插值数据；

G'_i,j(f)表示对应编号i、j的第一窗内频率-空间域正向预测的插值数据；

G”_i,j(f)表示对应编号i、j的第一窗内频率-空间域反向预测的插值数据。

步骤S506，按照如下公式，对对应编号i、j的第一窗内频率-空间域插值数据G_i,j(f)进行反傅立叶变换，得到对应编号i、j的第一窗内时间-空间插值数据：

G_i,j(t)＝FFT^-1(G_i,j(f))

其中，G_i,j(t)表示对应编号i、j的第一窗内时间-空间插值数据；

G_i,j(f)表示对应编号i、j的第一窗内频率-空间域插值数据。

步骤S507，针对当前窗口与前一窗口的重叠部分，将该重叠部分在当前窗口内对应的第一窗内时间-空间插值数据与该重叠部分在前一窗口内对应的第一窗内时间-空间插值数据进行均值处理，并将该重叠部分在当前窗口内对应的第一窗内时间-空间插值数据更新为该均值处理的结果。

可选地，该步骤所采用的均值处理是线性插值处理或者三次函数插值处理。

具体实施时，当期望均值处理的速度较快时，可以采用线性插值法，当期望均值处理的结果更精确时，可以采用三次函数插值法。

循环执行步骤S501～步骤S507，直至遍历完所有窗口。

步骤S6，将所有窗口的三维地震数据和第一窗内时间-空间插值数据作为第一方向的插值结果进行输出。

本发明提供的示例性方法，在计算预测插值算子的过程中使用的是部分三维地震数据而不是二维地震数据，更能适应局部三维的构造变化，并且由于信息量的增加，抗干扰能力也得到增强。本发明提供的三维地震道内插方法是将二维f-x域预测扩展到三维f-x,y域，不受倾角限制和空间假频的影响，充分利用三维数据信息，内插后的信噪比和精度都获得提高。

本发明提供的示例性方法中，第一方向可以是三维地震数据的共中心点号的分布方向，或三维地震数据的线号的分布方向，以上步骤S5是选定第一方向为褶积方向进行插值，即，步骤S5可以是沿共中心点号的分布方向进行插值，也可以是沿线号的分布方向进行插值。

进一步的，由于第一方向的预测插值算子和第二方向的预测插值算子是相互独立的，可以分别求取，因此为了获取更密集的插值效果，该示例性方法还可以在沿第一方向插值之后，基于第一方向的插值结果，继续对第二方向进行插值。

可选地，如图2所示，本示例性方法在步骤S6之后，还包括：

步骤S7，继续选定第二方向作为褶积方向，然后沿时间方向顺序遍历所有窗口，并对遍历到的每一窗口执行如下步骤S701～步骤S707的处理：

步骤S701，按照如下公式，将当前窗口的三维地震数据和第一窗内时间-空间插值数据沿时间方向做傅里叶变换，得到第二频率-空间域数据：

I_m,n(f)＝FFT(I_m,n(t))

其中，t表示三维地震数据和第一窗内时间-空间插值数据在时间方向上对应的时间；

m、n分别表示在第二方向、第一方向上按三维地震数据和第一窗内时间-空间插值数据在时间方向上对应的时间，对三维地震数据和第一窗内时间-空间插值数据重新排序后的编号；

f表示频率；

I_m,n(f)表示对应编号m、n的第二频率-空间域数据；

I_m,n(t)表示对应编号m、n的三维地震数据或第一窗内时间-空间插值数据。

步骤S702，按照如下公式，在第二方向的正向进行褶积，计算第二方向的预测插值算子：

$\underset{y,x}{\Σ}{I}_{m+1-2y,n-1+x}\left(f\right)p_{y,x}\left(f\right)=I_{m,n}\left(f\right)$

其中，x＝1,2,…,lx'；

y＝1,2,…,ly'；

p_y,x(f)表示对应频率f的第二方向的预测插值算子；

lx'表示预测插值算子p_y,x(f)在第一方向上的长度；

ly'表示预测插值算子p_y,x(f)在第二方向上的长度；

I_m+1-2y,n-1+x(f)表示对应编号m+1-2y、n-1+x的第二频率-空间域数据；

I_m,n(f)表示对应编号m、n的第二频率-空间域数据。

步骤S703，按照如下公式，将第二方向的对应频率f/2的预测插值算子p_y,x(f/2)与对应编号m+1-y、n-1+x的第二频率-空间域数据I_m+1-y,n-1+x(f)在第二方向的正向进行褶积，得到第二窗内频率-空间域正向预测的插值数据：

${G^{\′}}_{m,n}\left(f\right)=\underset{y,x}{\Σ}{I}_{m+1-y,n-1+x}\left(f\right)p_{y,x}(f/2)$

其中，p_y,x(f/2)表示对应频率f/2的第二方向的预测插值算子；

I_m+1-y,n-1+x(f)表示对应编号m+1-y、n-1+x的第二频率-空间域数据；

G'_m,n(f)表示对应编号m、n的第二窗内频率-空间域正向预测的插值数据。

步骤S704，按照如下公式，计算对应频率f/2的第二方向的预测插值算子p_y,x(f/2)的共轭，并将该共轭与对应编号m-1+y_、n+1-x的第二频率-空间域数据I_m-1+y,n+1-x(f)在第二方向的反向进行褶积，得到第二窗内频率-空间域反向预测的插值数据：

${G^{\′\′}}_{m,n}\left(f\right)=\underset{y,x}{\Σ}{I}_{m-1+y,n+1-x}\left(f\right)p{*}_{y,x}(f/2)$

其中，p*_y,x(f/2)表示对应频率f/2的第二方向的预测插值算子p_y,x(f/2)的共轭；

I_m-1+y,n+1-x(f)表示对应编号m-1+y、n+1-x的第二频率-空间域数据；

G”_m,n(f)表示对应编号m、n的第二窗内频率-空间域反向预测的插值数据。

步骤S705，按照如下公式，对对应编号m、n的第二窗内频率-空间域正向预测的插值数据和第二窗内频率-空间域反向预测的插值数据求平均，即得到对应编号m、n的第二窗内频率-空间域插值数据：

G_m,n(f)＝(G”_m,n(f)+G'_m,n(f))/2

其中，G_m,n(f)表示对应编号m、n的第二窗内频率-空间域插值数据；

G'_m,n(f)表示对应编号m、n的第二窗内频率-空间域正向预测的插值数据；

G”_m,n(f)表示对应编号m、n的第二窗内频率-空间域反向预测的插值数据。

步骤S706，按照如下公式，对对应编号m、n的第二窗内频率-空间域插值数据G_m,n(f)进行反傅立叶变换，得到对应编号m、n的第二窗内时间-空间插值数据：

G_m,n(t)＝FFT^-1(G_m,n(f))

其中，G_m,n(t)表示对应编号m、n的第二窗内时间-空间插值数据；

G_m,n(f)表示对应编号m、n的第二窗内频率-空间域插值数据。

步骤S707，针对当前窗口与前一窗口的重叠部分，将该重叠部分在当前窗口内对应的第二窗内时间-空间插值数据与该重叠部分在前一窗口内对应的第二窗内时间-空间插值数据进行均值处理，并将该重叠部分在当前窗口内对应的第二窗内时间-空间插值数据更新为该均值处理的结果。

可选地，该步骤所采用的均值处理是线性插值处理或者三次函数插值处理。

具体实施时，当期望均值处理的速度较快时，可以采用线性插值法，当期望均值处理的结果更精确时，可以采用三次函数插值法。

循环执行步骤S701～步骤S707，直至遍历完所有窗口。

步骤S8，将所有窗口的三维地震数据、第一窗内时间-空间插值数据和第二窗内时间-空间插值数据作为所述第一方向和所述第二方向最终的插值结果进行输出。

本发明提供的示例性方法中，第二方向可以是三维地震数据的共中心点号的分布方向，或三维地震数据的线号的分布方向，以上步骤S7是选定第二方向为褶积方向进行插值，即，步骤S7可以是沿共中心点号的分布方向进行插值，也可以是沿线号的分布方向进行插值。

步骤S7完成之后，继续执行步骤S8，输出先后经过第一方向和第二方向插值的结果。

实施例一

本实施例为应用本发明示例性方法的一具体实施例，该实施例仅对第一方向进行插值。

如下所示，该实施例包括如下过程：

1)采集三维地震数据。

2)对采集的三维地震数据进行加载观测系统、静校正、叠前去噪、振幅恢复、动校正和叠加等处理，获得叠后三维地震数据。

3)沿着三维地震数据的共中心点号的分布方向和线号的分布方向，对叠后三维地震数据排序。设定三维地震数据的共中心点号的分布方向为第一方向，线号的分布方向为第二方向，在第一方向上，排序后的三维地震数据具有编号i，在第二方向上，排序后的三维地震数据具有编号j。

4)将排序后的三维地震数据散布于以时间方向、第一方向和第二方向为轴的三维空间中，并以每个窗口在时间方向的跨度为2秒，第一方向的跨度为15道，第二方向的跨度也为15道，相邻的两个窗口重叠部分的数据量占到窗口全部数据量的50％为划分标准，在该三维空间中对排序后的三维地震数据划分窗口。

5)选定第一方向为褶积方向，然后顺序遍历所有窗口，并对遍历到的每一窗口执行如下处理：

51)将当前窗口内的三维地震数据S_i,j(t)沿时间方向做傅立叶变换，得到第一频率-空间域数据S_i,j(f)，其中，S_i,j(f)＝FFT(S_i,j(t))。

52)按照如下公式在第一方向的正向进行褶积，计算第一方向的预测插值算子p_x,y(f)：

$\underset{x,y}{\Σ}{S}_{i+1-2x,j-1+y}\left(f\right)p_{x,y}\left(f\right)=S_{i,j}\left(f\right)$

其中，预测插值算子p_x,y(f)在第一方向的长度为5，x＝1,2,…,5，预测插值算子p_x,y(f)在第二方向的长度为3，y＝1,2,…,3；

53)按照如下公式将第一方向的对应频率f/2的预测插值算子p_x,y(f/2)与对应编号i+1-x、j-1+y的第一频率-空间域数据S_i+1-x,j-1+y(f)在第一方向的正向进行褶积，得到第一窗内频率-空间域正向预测的插值数据G'_i,j(f)；其中， ${G^{\′}}_{i,j}\left(f\right)=\underset{x,y}{\Σ}{S}_{i+1-x,j-1+y}\left(f\right)p_{x,y}(f/2).$

54)计算第一方向的对应频率f/2的预测插值算子p_x,y(f/2)的共轭p*_x,y(f/2)，并将该共轭与对应编号i-1+x、j+1-y的第一频率-空间域数据S_i-1+x,j+1-y(f)在第一方向的反向进行褶积，得到第一窗内频率-空间域反向预测的插值数据G”_i,j(f)；其中， ${G^{\′\′}}_{i,j}\left(f\right)=\underset{x,y}{\Σ}{S}_{i-1+x,j+1-y}\left(f\right)p{*}_{x,y}(f/2).$

55)对对应编号i、j的第一窗内频率-空间域正向预测的插值数据和第一窗内频率-空间域反向预测的插值数据求平均，即得到对应编号i、j的第一窗内频率-空间域插值数据G_i,j(f)；其中，G_i,j(f)＝(G”_i,j(f)+G'_i,j(f))/2。

56)对对应编号i、j的第一窗内频率-空间域插值数据G_i,j(f)进行反傅立叶变换，得到对应编号i、j的第一窗内时间-空间插值数据G_i,j(t)；其中，G_i,j(t)＝FFT^-1(G_i,j(f))。

57)针对当前窗口与前一窗口的重叠部分，将该重叠部分在当前窗口内对应的第一窗内时间-空间插值数据与该重叠部分在前一窗口内对应的窗内时间-空间插值数据进行均值处理，并将该重叠部分在当前窗口内对应的第一窗内时间-空间插值数据更新为该均值处理的结果。

循环执行51)～57)，直至遍历完所有窗口。

6)将所有窗口的三维地震数据、第一窗内时间-空间插值数据作为所述第一方向的插值结果进行输出。

如图3所示，为本实施例内插前的三维地震数据剖面；如图4所示为本实施例内插后的三维地震数据剖面；对比图3、图4可知内插后的信噪比和精度都获得提高。

本发明预测插值算子使用的是局部三维数据体而不是二维数据，更能适应局部三维的构造变化，并且由于信息量的增加，抗干扰能力也得到增强。三维地震道内插方法是将二维f-x域预测扩展到三维f-x,y域，不受倾角限制和空间假频的影响，充分利用三维数据信息，内插后的信噪比和精度都获得提高。

本发明能应用于叠前和叠后三维地震数据，关键是通过排序满足三维f-x,y域道内插的规则化要求。在三维地震数据应用中，第一方向和第二方向的预测插值算子是相互独立的，两个方向的预测插值算子可分别求取，因此可先沿其中一个方向完成插值后，再进行另一个方向的插值，以获得更加密集的内插效果。本发明实现简便、能够适应大数据量的处理。

以上所述的具体实施例，对本发明的目的、技术方案和有益效果进行了进一步详细说明，所应理解的是，以上所述仅为本发明的具体实施例而已，并不用于限定本发明的保护范围，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (9)

1.一种三维地震数据道内插方法，其特征在于，包括：

步骤1，采集三维地震数据；

步骤2，对所述三维地震数据进行预处理；

步骤3，沿着地表二维平面的相互正交的第一方向和第二方向，对预处理后的所述三维地震数据排序；其中，所述第一方向为所述三维地震数据的共中心点号的分布方向，所述第二方向为所述三维地震数据的线号的分布方向；或，所述第一方向为所述三维地震数据的线号的分布方向，所述第二方向为所述三维地震数据的共中心点号的分布方向；

步骤4，以时间方向以及所述第一方向和第二方向为轴建立三维空间，将排序后的所述三维地震数据布散于所述三维空间中，在所述三维空间中对排序后的所述三维地震数据划分若干三维的窗口，使得每个窗口对应于排序后的所述三维地震数据的一部分，各个窗口对应的数据量相等，且相邻的窗口具有重叠部分；

步骤5，选定所述第一方向作为褶积方向，沿时间方向顺序遍历各个窗口，并对遍历到的每一窗口执行如下处理：

步骤51，按照如下公式，将当前窗口的三维地震数据沿所述时间方向做傅里叶变换，得到第一频率-空间域数据：

S_i,j(f)＝FFT(S_i,j(t))

其中，i、j分别表示所述三维地震数据在所述第一方向、所述第二方向上的编号；

t表示所述三维地震数据在所述时间方向上对应的时间；

f表示频率；

S_i,j(f)表示对应编号i、j的第一频率-空间域数据；

S_i,j(t)表示对应编号i、j的三维地震数据；

步骤52，按照如下公式，在所述第一方向的正向进行褶积，计算所述第一方向的预测插值算子：

$\underset{x,y}{\Σ}{S}_{i+1-2x,j-1+y}\left(f\right)p_{x,y}\left(f\right)=S_{i,j}\left(f\right)$

其中，x＝1,2,…,lx；

y＝1,2,…,ly；

p_x,y(f)表示对应频率f的所述第一方向的预测插值算子；

lx表示预测插值算子p_x,y(f)在所述第一方向上的长度；

ly表示预测插值算子p_x,y(f)在所述第二方向上的长度；

S_i+1-2x,j-1+y(f)表示对应编号i+1-2x、j-1+y的第一频率-空间域数据；

S_i,j(f)表示对应编号i、j的第一频率-空间域数据；

步骤53，按照如下公式，将所述第一方向的预测插值算子以及所述第一频率-空间域数据在所述第一方向的正向进行褶积，得到第一窗内频率-空间域正向预测的插值数据：

${G^{\′}}_{i,j}\left(f\right)=\underset{x,y}{\Σ}{S}_{i+1-x,j-1+y}\left(f\right)p_{x,y}(f/2)$

其中，p_x,y(f/2)表示对应频率f/2的所述第一方向的预测插值算子；

S_i+1-x,j-1+y(f)表示对应编号i+1-x、j-1+y的第一频率-空间域数据；

G'_i,j(f)表示对应编号i、j的第一窗内频率-空间域正向预测的插值数据；

步骤54，计算所述第一方向的预测插值算子的共轭，并按照如下公式将该共轭与所述第一频率-空间域数据在所述第一方向的反向进行褶积，得到第一窗内频率-空间域反向预测的插值数据：

${G^{\′\′}}_{i,j}\left(f\right)=\underset{x,y}{\Σ}{S}_{i-1+x,j+1-y}\left(f\right)p{*}_{x,y}(f/2)$

其中，p*_x,y(f/2)为对应频率f/2的第一方向的预测插值算子p_x,y(f/2)的共轭；

S_i-1+x,j+1-y(f)表示对应编号i-1+x、j+1-y的第一频率-空间域数据；

G”_i,j(f)表示对应编号i、j的第一窗内频率-空间域反向预测的插值数据；

步骤55，按照如下公式，对所述第一窗内频率-空间域正向预测的插值数据和所述第一窗内频率-空间域反向预测的插值数据求平均，得到第一窗内频率-空间域插值数据：

G_i,j(f)＝(G”_i,j(f)+G'_i,j(f))/2

其中，G_i,j(f)表示对应编号i、j的第一窗内频率-空间域插值数据；

步骤56，按照如下公式，对所述第一窗内频率-空间域插值数据进行反傅立叶变换，得到第一窗内时间-空间插值数据：

G_i,j(t)＝FFT^-1(G_i,j(f))；

其中，G_i,j(t)表示对应编号i、j的第一窗内时间-空间插值数据；

步骤57，针对当前窗口与前一窗口的重叠部分，将该重叠部分在当前窗口内对应的第一窗内时间-空间插值数据与该重叠部分在前一窗口内对应的第一窗内时间-空间插值数据进行均值处理，并将该重叠部分在当前窗口内对应的第一窗内时间-空间插值数据更新为该均值处理的结果；

步骤6，将所有窗口的三维地震数据和第一窗内时间-空间插值数据作为所述第一方向的插值结果进行输出。

2.根据权利要求1所述的三维地震数据道内插方法，其特征在于，在所述步骤7之后，还包括：

步骤7，继续选定所述第二方向作为褶积方向，沿时间方向顺序遍历所有窗口，并对遍历到的每一窗口执行如下处理：

步骤71，按照如下公式，将当前窗口的三维地震数据和第一窗内时间-空间插值数据沿时间方向做傅里叶变换，得到第二频率-空间域数据：

I_m,n(f)＝FFT(I_m,n(t))

其中，t表示三维地震数据和第一窗内时间-空间插值数据在时间方向上对应的时间；

m、n分别表示在所述第二方向、所述第一方向上按三维地震数据和第一窗内时间-空间插值数据在时间方向上对应的时间，对三维地震数据和第一窗内时间-空间插值数据排序后的编号；

I_m,n(f)表示对应编号m、n的第二频率-空间域数据；

I_m,n(t)表示对应编号m、n的三维地震数据或第一窗内时间-空间插值数据；

步骤72，按照如下公式，在所述第二方向的正向进行褶积，计算所述第二方向的预测插值算子：

$\underset{y,x}{\Σ}{I}_{m+1-2y,n-1+x}\left(f\right)p_{y,x}\left(f\right)=I_{m,n}\left(f\right)$

其中，x＝1,2,…,lx'；

y＝1,2,…-ly'；

p_y,x(f)表示对应频率f的所述第二方向的预测插值算子；

lx'表示预测插值算子p_y,x(f)在所述第一方向上的长度；

ly'表示预测插值算子p_y,x(f)在所述第二方向上的长度；

I_m+1-2y,n-1+x(f)表示对应编号m+1-2y、n-1+x的第二频率-空间域数据；

I_m,n(f)表示对应编号m、n的第二频率-空间域数据；

步骤73，按照如下公式，将所述第二方向的预测插值算子以及所述第二频率-空间域数据在所述第二方向的正向进行褶积，得到第二窗内频率-空间域正向预测的插值数据：

${G^{\′}}_{m,n}\left(f\right)=\underset{y,x}{\Σ}{I}_{m+1-y,n-1+x}\left(f\right)p_{y,x}(f/2)$

其中，p_y,x(f/2)表示对应频率f/2的所述第二方向的预测插值算子；

I_m+1-y,n-1+x(f)表示对应编号m+1-y、n-1+x的第二频率-空间域数据；

G'_m,n(f)表示对应编号m、n的第二窗内频率-空间域正向预测的插值数据；

步骤74，按照如下公式，计算所述第二方向的预测插值算子的共轭，并按照如下公式将该共轭与所述第二频率-空间域数据在所述第二方向的反向进行褶积，得到第二窗内频率-空间域反向预测的插值数据：

${G^{\′\′}}_{m,n}\left(f\right)=\underset{y,x}{\Σ}{I}_{m-1+y,n+1-x}\left(f\right)p{*}_{y,x}(f/2)$

其中，p*_y,x(f/2)表示对应频率f/2的第二方向的预测插值算子p_y,x(f/2)的共轭；

I_m-1+y,n+1-x(f)表示对应编号m-1+y_、n+1-x的第二频率-空间域数据；

G”_m,n(f)表示对应编号m、n的第二窗内频率-空间域反向预测的插值数据；

步骤75，按照如下公式，对所述第二窗内频率-空间域正向预测的插值数据和所述第二窗内频率-空间域反向预测的插值数据求平均，得到第二窗内频率-空间域插值数据：

G_m,n(f)＝(G”_m,n(f)+G'_m,n(f))/2

其中，G_m,n(f)表示对应编号m、n的第二窗内频率-空间域插值数据；

步骤76，按照如下公式，对所述第二窗内频率-空间域插值数据进行反傅立叶变换，得到第二窗内时间-空间插值数据：

G_m,n(t)＝FFT^-1(G_m,n(f))

其中，G_m,n(t)表示对应编号m、n的第二窗内时间-空间插值数据；

步骤77，针对当前窗口与前一窗口的重叠部分，将该重叠部分在当前窗口内对应的第二窗内时间-空间插值数据与该重叠部分在前一窗口内对应的第二窗内时间-空间插值数据进行均值处理，并将该重叠部分在当前窗口内对应的第二窗内时间-空间插值数据更新为该均值处理的结果；

步骤8，将所有窗口的三维地震数据、第一窗内时间-空间插值数据和第二窗内时间-空间插值数据作为所述第一方向和所述第二方向最终的插值结果进行输出。

3.根据权利要求1所述的三维地震数据道内插方法，其特征在于，所述步骤2包括：

对所述三维地震数据进行加载观测系统、静校正处理、叠前去噪处理、振幅恢复处理、动校正处理和叠加处理。

4.根据权利要求1所述的三维地震数据道内插方法，其特征在于，所述步骤3包括：

确定所述三维地震数据对应的三维地震记录道的地下反射点在所述地表二维平面的垂直投影点；

根据所述垂直投影点沿所述第一方向、所述第二方向的分布情况，对所述三维地震记录道排序；

根据所述三维地震记录道的排序情况，对预处理后的所述三维地震数据排序。

5.根据权利要求1所述的三维地震数据道内插方法，其特征在于，所述步骤4包括：

按照每个窗口沿所述时间方向的跨度为第一设定值、沿所述第一方向的跨度为第二设定值、沿所述第二方向的跨度为第三设定值，且相邻的窗口的重叠部分的数据量占所述窗口内全部数据量的一预设比例为划分标准，对排序后的所述三维地震数据划分窗口。

6.根据权利要求1所述的三维地震数据道内插方法，其特征在于，所述步骤57中的均值处理为线性插值处理。

7.根据权利要求1所述的三维地震数据道内插方法，其特征在于，所述步骤57中的均值处理为三次函数插值处理。

8.根据权利要求2所述的三维地震数据道内插方法，其特征在于，所述步骤77中的均值处理为线性插值处理。

9.根据权利要求2所述的三维地震数据道内插方法，其特征在于，所述步骤77中的均值处理为三次函数插值处理。