地震子波提取方法与流程

本发明涉及油气及煤层气地震勘探与开发领域，更具体地，涉及一种地震子波提取方法。

背景技术：

随着油气勘探开发的需要，储层预测和精细描述越来越引起重视，围绕这个目的而开展的探索和研究也越来越多，地震反演就是最重要的一点。当反演的目标是波阻抗时，我们称之为波阻抗反演。目前的地震波阻抗反演按照所使用的资料分为叠后反演和叠前反演两大类，无论哪种反演，其实质都是想去除子波的影响，从而将地震剖面转化为能与钻井、地质等资料直接对比的形式，因此，反演在许多情况下提高了常规地震的分辨率，并提高了油藏参数研究的水平，这对于研究复杂油气储层的空间分布、开展复杂油气藏精细描述都具有十分重要的意义。

而想要去除子波的影响，首先要提取地震子波。目前，地震子波提取的方法主要分为两类：确定性子波提取和统计性子波提取。这两种子波提取方法都是基于同一个物理基础，即地震褶积模型，无噪声褶积模型表示如下：

s(t)＝w(t)*r(t)(1)

其中，s(t)表示地震信号，w(t)表示子波，r(t)表示反射系数序列，*表示褶积。

统计性子波提取是通过地震数据本身来估计子波。假设地震子波是不变的，反射系数序列是具有白噪谱的随机序列，则地震道的自相关就等同于地震子波的自相关估计，从而获得地震子波的振幅谱(子波二阶统计量)，对于相位谱，则一般假设为最小相位或者零相位。由于地震子波实际上是混合相位，因此出现了包含子波相位谱的高阶统计量方法。这种方法的优点是不需要井数据，但是它是基于统计的方法，因此需要大量统计，另外，基于高阶统计量的方法往往需要求解方程，使用最优化技术，这都增加了计算量，并且求解过程并不确保稳定。

确定性子波提取方法是利用测井资料计算出反射系数序列，然后结合井旁地震道，由褶积理论求出地震子波，通常使用的计算方法有最小平方维纳滤波、谱除法等。它的优点是不需要对反射系数序列的分布作任何假设，能得到较为准确的子波。但是在计算方法上没有太多的选择余地，同时现有的方法都是每次利用一口井，获得一个子波，然后对多个子波进行平均，而不能一次同时利用多口井。因此，有必要开发一种精确简便的地震子波提取方法。

公开于本发明背景技术部分的信息仅仅旨在加深对本发明的一般背景技术的理解，而不应当被视为承认或以任何形式暗示该信息构成已为本领域技术人员所公知的现有技术。

技术实现要素：

本发明提出了一种地震子波提取方法，其能够采用字典学习方法，可以同时对多个井进行精确高效的地震子波提取。

本发明提出了一种地震子波提取方法。所述方法可以包括：基于井资料与地震资料，获得多个井的反射系数序列与对应的井旁地震道，其中所述多个井的反射系数序列以反射系数矩阵表示，所述对应的井旁地震道以井旁地震道矩阵表示；基于所述反射系数矩阵与所述井旁地震道矩阵，建立关于地震子波矩阵的字典学习方程；基于所述反射系数矩阵、所述井旁地震道矩阵与所述地震子波矩阵，建立目标函数，计算使所述目标函数的值最小的地震子波矩阵，作为最终的地震子波矩阵；对所述最终的地震子波矩阵进行提取与处理，获得平均子波。

优选地，基于井资料与地震资料，获得多个井的反射系数序列与对应的井旁地震道，其中所述多个井的反射系数序列以反射系数矩阵表示，所述对应的井旁地震道以井旁地震道矩阵表示包括：基于所述井资料与所述地震资料，计算每口井的反射系数序列，并将所述每口井的反射系数序列转换到时间域；对每口井进行井-震标定，获得每口井对应的井旁地震道；将时间域的每个反射系数序列作为列向量，获得所述反射系数矩阵，将每口井对应的井旁地震道作为列向量，获得所述井旁地震道矩阵。

优选地，所述字典学习方程为：

y＝d·x(2)

其中，y表示井旁地震道矩阵，x表示反射系数矩阵，d表示地震子波矩阵。

优选地，所述目标函数为：

其中，f(x)表示目标函数，x表示反射系数矩阵、y表示井旁地震道矩阵，d表示地震子波矩阵，γ表示权重因子，||·||p表示p范数。

优选地，对所述最终的地震子波矩阵进行提取与处理，获得平均子波包括：分别提取所述最终的地震子波矩阵的每一列作为一个地震子波序列，获得多个地震子波序列对应的起始位置与长度；针对每个地震子波序列，基于其对应的起始位置与长度，对其进行傅立叶变换，获得所述长度内的所述地震子波序列的振幅谱与相位谱；针对多个地震子波序列的振幅谱与相位谱进行平均计算，获得平均振幅谱与平均相位谱，进而获得所述平均子波。

优选地，所述平均子波为：

其中，a(ω)表示平均振幅谱，表示平均相位谱，i表示虚数单位，w(ω)表示频率域的平均子波。

优选地，还包括：将所述频率域的平均子波通过傅里叶反变换转换至时间域，获得时间域的平均子波。

本发明的有益效果在于：采用字典学习方法，可以同时利用多口井及相应的井旁地震道数据，同时对多个井进行地震子波提取，并且设定了规则化项，能够更加精确地提取地震子波。

本发明的方法具有其它的特性和优点，这些特性和优点从并入本文中的附图和随后的具体实施方式中将是显而易见的，或者将在并入本文中的附图和随后的具体实施方式中进行详细陈述，这些附图和具体实施方式共同用于解释本发明的特定原理。

附图说明

通过结合附图对本发明示例性实施方式进行更详细的描述，本发明的上述以及其它目的、特征和优势将变得更加明显，其中，在本发明示例性实施方式中，相同的参考标号通常代表相同部件。

图1示出了根据本发明的地震子波提取方法的步骤的流程图。

具体实施方式

下面将参照附图更详细地描述本发明。虽然附图中显示了本发明的优选实施方式，然而应该理解，可以以各种形式实现本发明而不应被这里阐述的实施方式所限制。相反，提供这些实施方式是为了使本发明更加透彻和完整，并且能够将本发明的范围完整地传达给本领域的技术人员。

图1示出了根据本发明的地震子波提取方法的步骤的流程图。

根据本发明的地震子波提取方法可以包括：

步骤101，基于井资料与地震资料，获得多个井的反射系数序列与对应的井旁地震道，其中多个井的反射系数序列以反射系数矩阵表示，对应的井旁地震道以井旁地震道矩阵表示；

步骤102，基于反射系数矩阵与井旁地震道矩阵，建立关于地震子波矩阵的字典学习方程；字典学习方法主要应用在稀疏信号表示，即假定信号是可以稀疏表示的，稀疏就是在许多采样点上的信号都是零。这一假设在地震子波提取中是可以得到满足的，因为通常在实际工作中，将测井资料得到的深度域波阻抗曲线方波化，然后转换到时间域时并计算出时间域反射系数序列时，许多样点值是零，因此，可以采用字典学习方法提取地震子波；

步骤103，基于反射系数矩阵、井旁地震道矩阵与地震子波矩阵，建立目标函数，计算使目标函数的值最小的地震子波矩阵，作为最终的地震子波矩阵；

步骤104，对最终的地震子波矩阵进行提取与处理，获得平均子波。

在一个示例中，基于井资料与地震资料，获得多个井的反射系数序列与对应的井旁地震道，其中多个井的反射系数序列以反射系数矩阵表示，对应的井旁地震道以井旁地震道矩阵表示包括：基于井资料与地震资料，计算每口井的反射系数序列，并将每口井的反射系数序列转换到时间域；对每口井进行井-震标定，获得每口井对应的井旁地震道；将时间域的每个反射系数序列作为列向量，获得反射系数矩阵，将每口井对应的井旁地震道作为列向量，获得井旁地震道矩阵。

在一个示例中，字典学习方程为：

y＝d·x(2)

其中，y表示井旁地震道矩阵，x表示反射系数矩阵，d表示地震子波矩阵。

在一个示例中，目标函数为：

其中，f(x)表示目标函数，x表示反射系数矩阵、y表示井旁地震道矩阵，d表示地震子波矩阵，γ表示权重因子，一般取5-10，||·||p表示p范数，p一般取1。目标函数的优化方法可以选用最优方向法、k-svd方法或者梯度下降法，当选择合适的步长时，梯度下降法可能优于最优方向法、k-svd方法。

在一个示例中，对最终的地震子波矩阵进行提取与处理，获得平均子波包括：分别提取最终的地震子波矩阵的每一列作为一个地震子波序列，获得多个地震子波序列对应的起始位置与长度；针对每个地震子波序列，基于其对应的起始位置与长度，对其进行傅立叶变换，获得长度内的地震子波序列的振幅谱与相位谱；针对多个地震子波序列的振幅谱与相位谱进行平均计算，获得平均振幅谱与平均相位谱，进而获得平均子波。

在一个示例中，平均子波为：

其中，a(ω)表示平均振幅谱，表示平均相位谱，i表示虚数单位，w(ω)表示频率域的平均子波。

在一个示例中，还包括：将频率域的平均子波通过傅里叶反变换转换至时间域，获得时间域的平均子波。

具体地，最终的地震子波矩阵包含了地震子波的信息，是由地震子波构成的。这些地震子波的长度都很长，但是起始点不同。接下来需要将这些地震子波进行融合或重构，形成一个地震子波。将最终的地震子波矩阵中的每一列提取出来作为地震子波序列，并用图形形式显示出来，以确定每个地震子波序列的起始点及长度，然后利用傅里叶变换求得该长度内的地震子波序列的振幅谱和相位谱，最后对每一列的地震子波序列对应的振幅谱和相位谱按频率分别进行平均计算，求得平均振幅谱和平均相位谱，利用公式(4)，求得频率域的平均子波，通过傅里叶反变换至时间域，获得时间域的平均子波。

本发明采用字典学习方法，可以同时利用多口井以及相应的井旁地震道数据，同时对多个井进行地震子波提取，并且设定了规则化项，能够更加精确地提取地震子波。

应用示例

为便于理解本发明实施方式的方案及其效果，以下给出一个具体应用示例。本领域技术人员应理解，该示例仅为了便于理解本发明，其任何具体细节并非意在以任何方式限制本发明。

输入工区内的所有井资料和地震资料，计算每口井上的反射系数序列，并将每口井上的反射系数序列转换到时间域，时间采样间隔定位1毫秒，然后进行井-震标定，获得每口井所对应的井旁地震道，即获得与井上的反射系数相应的起始点与终止点，将时间域的每个反射系数序列作为列向量，获得反射系数矩阵，将每口井对应的井旁地震道作为列向量，获得井旁地震道矩阵；基于反射系数矩阵与井旁地震道矩阵，建立关于地震子波矩阵的字典学习方程为公式(2)。基于反射系数矩阵、井旁地震道矩阵与地震子波矩阵，建立目标函数为公式(3)，其中，γ＝7，p＝1，采用最优方向法计算使目标函数的值最小的地震子波矩阵，作为最终的地震子波矩阵。

将最终的地震子波矩阵中的每一列提取出来作为地震子波序列，并用图形形式显示出来，以确定每个地震子波序列的起始点及长度，然后利用傅里叶变换求得该长度内的地震子波序列的振幅谱和相位谱，最后对每一列的地震子波序列对应的振幅谱和相位谱按频率分别进行平均计算，求得平均振幅谱和平均相位谱，利用公式(4)，求得频率域的平均子波，通过傅里叶反变换转换至时间域，获得时间域的平均子波。

综上所述，本发明采用字典学习方法，可以同时利用多口井以及相应的井旁地震道数据，同时对多个井进行地震子波提取，并且设定了规则化项，能够更加精确地提取地震子波。

本领域技术人员应理解，上面对本发明的实施方式的描述的目的仅为了示例性地说明本发明的实施方式的有益效果，并不意在将本发明的实施方式限制于所给出的任何示例。

以上已经描述了本发明的各实施方式，上述说明是示例性的，并非穷尽性的，并且也不限于所披露的各实施方式。在不偏离所说明的各实施方式的范围和精神的情况下，对于本技术领域的普通技术人员来说许多修改和变更都是显而易见的。本文中所用术语的选择，旨在最好地解释各实施方式的原理、实际应用或对市场中的技术的改进，或者使本技术领域的其它普通技术人员能理解本文披露的各实施方式。