CN112578436A - 一种子波提取方法及系统 - Google Patents

Abstract

本发明提出了一种子波提取方法及系统。该方法包括：对地震道进行频谱分解，利用地震道的频谱特征，获取初始子波；将获取的初始子波看成多项式函数的系数，计算该函数的特征值；建立目标函数，利用遗传算法，求取最优解；根据求取的最优解组合，反向求解子波。本发明针对叠后反演技术过程中地震子波提取这一难点问题，将子波看成一个多项式的系数，通过遗传算法来优化该多项式，提供了一种优化子波提取的方法，为在地震勘探中解决地震反演等技术问题提供支持。

Description

一种子波提取方法及系统

技术领域

本发明属于油气勘探领域，涉及油气勘探中地震资料信号处理技术领域，更具体地，涉及一种子波提取方法及系统。

背景技术

地震子波为具有确定的起始时间、能量有限且具有一定连续长度的信号，它是地震记录中的基本单元，因此子波提取在地震勘探中占据着重要的位置。在地震反演问题中，子波提取的精度决定了最终反演的精度。

由于子波和反射系数序列“平滑度”的差别，地震子波和反射系数的复赛谱在时间轴上分布不同。因此在时间域进行低通滤波，可以在复赛谱中实现子波和反射系数的分离，达到子波提取的目的。

子波在频率域可以表示为它的振幅谱和相位谱。振幅谱可以通过地震记录振幅谱的包络计算出来，但是其相位谱却很难得到，这是因为子波的相位没有一定的规律，且随频率变化明显，不容易通过解析式进行描述，所以在实际应用中很难通过一个线性滤波器将两者完全分离。

发明内容

本发明针对叠后反演技术过程中地震子波提取这一难点问题，将子波看成一个多项式的系数，通过遗传算法来优化该多项式，提供了一种优化子波提取的方法。为在地震勘探中解决地震反演等技术问题提供支持。

根据本发明一个方面，提供一种子波提取方法，该方法包括：

对地震道进行频谱分解，利用地震道的频谱特征，获取初始子波；

将获取的初始子波看成多项式函数的系数，计算该函数的特征值；

建立目标函数，利用遗传算法，求取最优解；

根据求取的最优解组合，反向求解子波。

进一步地，通过复赛谱提取初始子波；

s(t)＝w(t)*r(t)+n(t)； (1)

其中，s(t)表示地震信号，w(t)表示子波，r(t)表示反射系数序列，n(t)表示附加干扰；

使用傅里叶变换，转换到频率域：s(ω)＝w(ω)r(ω)； (2)

取对数：ln s(ω)＝lnw(ω)+lnr(ω)； (3)

使用反傅里叶变换，转换回时间域：sp(t)＝wp(t)+rp(t)； (4)

则sp(t)、wp(t)、rp(t)是s(t)、w(t)、r(t)的复赛谱。

优选地，将初始子波看成多项式的一组系数，即：

其中，a_i，i∈[0，n]是子波的每个样点值；

设多项式的n个根为γ_i，i∈[0，n]，则多项式P(x)表示为：

令P(x)＝0，利用QR方法求得多项式方程的所有实根，得到不同的多项式系数。

进一步地，建立目标函数为：

Fitness＝MAX{R_(l)} (7)

其中相关系数R_(l)的公式为：

其中，x长度为n，为地震数据；y长度为m，为合成记录数据；l为是相关系数数组的下标；t为输入地震数据数组的下标。

进一步地，利用遗传算法，求取最优解包括：

步骤1：进行染色体编码并初始化；

步骤2：计算当前种群中的染色体的适应度；

步骤3：调整染色体的编码，产生新的染色体；

步骤4：根据染色体的编码调整根的数值；

步骤5：根据调整后的根生成调整子波，并计算新个体的适应度；

步骤6：将新的个体插入种群，并淘汰适应度差的个体；

步骤7：如果达到规定代数，返回最优染色体，否则返回步骤3。

进一步地，步骤1中染色体编码包括将染色体表示成如下向量：

C＝(c₀，c₁，c₂，c₃，c₄……c_i……c_n)，c_i∈[0，1]，取值0表示不进行调整，取值1表示进行调整。

进一步地，步骤3中调整染色体的编码包括：选择、交叉、变异操作。

(1)选择操作包括：

step1.计算当前种群关于适应度值的累积分布函数Pf_(fitness)；

step2.通过随机函数生成一个[0，1]区间的随机数ran；

step3.计算ran在累积分布函数Pf_(fitness)中的分布区间；

step4.选择该区间所对应的个体为新的个体；

step5.重复step2-step4，直到满足种群中个体的上限。

(2)交叉操作包括：

step1.给定交叉概率Pc_(cross)；

step2.利用随机函数生成一个[0，1]区间的随机数ran，判断ran是否超过交叉概率；

step3.如果通过，利用随机函数选择两个随机的进行交叉的个体A，B；

step4.利用随机函数选择交叉的染色体位置i，交换t＝A.c_i，A.c_i＝B.c_i，B.c_i＝t；

step5.重复step2-step4，直到满足种群中个体的上限。

(3)变异操作包括：

step1.给定变异概率Pm_(mutation)；

step2.利用随机函数生成一个[0，1]区间的随机数ran，判断ran是否超过变异概率；

step3.如果通过，利用随机函数选择一个随机的进行变异的个体A；

step4.利用随机函数选择变异的染色体位置i，进行取反操作，A.c_i＝～A.c_i；

step5.重复step2-step4，直到满足种群中个体的上限。

进一步地，根据染色体的编码调整根的数值，包括：

步骤4.1：遍历整个染色体；

步骤4.2：染色体编码为0，复制该位置处对应的根；

步骤4.3：染色体编码为1，调整该位置处对应的根；

设原有根为γ_k＝(a+bi)，调整为γ′_k＝(a+bi)/|(γ_k)|²，同时记录调整因子

步骤4.4：重复步骤4.2-步骤4.3，直到满足种群中个体的上限，调整常数项系数a₀＝a₀/|a′_k|。

进一步地，使用递推法，反向求解多项式系数，生成子波。

根据本发明的另一方面，提供一种子波提取系统，该系统包括：

存储器，存储有计算机可执行指令；

处理器，所述处理器运行所述存储器中的计算机可执行指令，执行以下步骤：

对地震道进行频谱分解，利用地震道的频谱特征，获取初始子波；

将获取的初始子波看成多项式函数的系数，计算该函数的特征值；

建立目标函数，利用遗传算法，求取最优解；

根据求取的最优解组合，反向求解子波。

本发明针对叠后反演技术过程中地震子波提取这一难点问题，将子波看成一个多项式的系数，通过遗传算法来优化该多项式，提供了一种优化子波提取的方法，为在地震勘探中解决地震反演等技术问题提供支持。

附图说明

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

图1是根据本发明实施例的一种子波提取方法的基本流程示意图。

图2是输入的地震数据和反射系数的示意图，包括该井在地震剖面上的投影。

图3是通过复赛谱提取的初始子波，该子波是一个零相位子波。

图4a-4d是多个通过本发明提取的混合相位子波。

具体实施方式

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

本发明公开一种综合利用井中反射系数和井旁地震数据进行子波提取的技术，属于石油地球物理勘探中的地震资料信号处理技术领域。本发明方法提出了一种通过多项式拟合进行子波提取的方法。如图1所示，本发明的方法包括以下步骤：

(1)对地震道进行频谱分解，利用地震道的频谱特征，获取初始子波；

(2)将获取的子波看成是一个多项式函数的系数，计算该函数的特征值；

(3)建立目标函数，利用遗传算法，求取最优解；

(4)根据特征值组合，反向求解子波。

所述方法包括以下过程：

步骤(1)对地震道进行频谱分解，利用地震道的频谱特征，获取初始子波。初始子波的提取有多种方式，可以通过确定性方法，或者通过统计性方法。因为本发明重点在后续的优化方法，因此，在此采取一种比较容易实现方法，通过复赛谱提取初始子波。

具体地，子波提取基于地震褶积模型这一物理基础，可以表示成：

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

其中，s(t)表示地震信号，w(t)表示子波，r(t)表示反射系数序列，n(t)表示附加干扰。

使用傅里叶变换，转换到频率域：

s(ω)＝w(ω)r(ω) (2)

取对数：

ln s(ω)＝lnw(ω)+lnr(ω) (3)

使用反傅里叶变换，转换回时间域：

sp(t)＝wp(t)+rp(t) (4)

则sp(t)、wp(t)、rp(t)是s(t)、w(t)、r(t)的复赛谱，由于子波和反射系数序列“平滑度”的差别，地震子波和反射系数的复赛谱在时间轴上分布不同。地震子波在时间轴原点附近，后者离原点远，即地震子波表现为低频，而反射系数表现为高频。因此，可以使用低通滤波，获取初始子波。

接下来，步骤(2)将获取的子波看成是一个多项式函数的系数，计算该多项式的根。

为了方便后面对子波形状进行调整，将子波看成多项式的一组系数，即：

其中，a_i，i∈[0，n]就是子波的每个样点值，对于这种表达形式的多项式，比较难以进行优化，因此，需要转换为连根式，方便后续调整。

设多项式的n个根为γ_i，i∈[0，n]，则多项式P(x)可以表示为：

可以令P(x)＝0，利用QR方法求得多项式方程的所有实根。

将多项式分解为连根式之后，则只需对多项式的根进行调整，就可以得到不同的多项式系数，方便后续求解。

接下来，步骤(3)中建立目标函数，利用遗传算法，求取最优解。

由于对上面的多项式进行优化，牵涉到大量的独立变量的值的调整，因此本发明使用遗传算法来进行优化。

本发明中，遗传算法是一种常用的寻优算法，基本操作包括选择，交叉，变异，基本过程如下：

step1.进行染色体编码并初始化；

step2.计算当前种群中的染色体的适应度；

step3.调整染色体的编码，产生新的染色体，调整的方法包括选择，交叉，变异等操作；

step4.根据染色体的编码调整根的数值；

step5.根据调整后的根生成调整子波，并计算新个体的适应度；

step6.将新的个体插入种群，并淘汰适应度差的个体；

step7.如果达到规定代数，返回最优染色体，否则返回step3。

对于本算法，需要进行以下说明：

首先，建立适应度函数(目标函数)，用于计算个体的适应度。根据公式(1)，已知初始子波或者调整后子波，结合给定的反射系数，即可以通过褶积计算出合成记录，有合成记录之后，就可以通过原始的井旁地震道数据对合成记录的准确度进行判断。本发明中通过计算合成记录和地震数据的最大相关系数作为适应度函数。

适应度函数(目标函数)：

Fitness＝MAX{R_(l)} (7)

相关系数R_(l)公式为：

其中，x长度为n，为地震数据；y长度为m，为合成记录数据；l为是相关系数数组的下标；t为输入地震数据数组的下标。

其次，本方法对染色体的编码主要用于从长度为n的根数组中，选择出要进行调整的根，因此，将染色体表示成如下向量：C＝(c₀，c₁，c₂，c₃，c₄……c_i……c_n)，c_i∈[0，1]，取值0表示不进行调整，取值1表示进行调整。

具体地，调整染色体的编码包括：

(1)选择操作：

step1.计算当前种群关于适应度值的累积分布函数Pf_(fitness)；

step2.通过随机函数生成一个[0，1]区间的随机数ran；

step3.计算ran在累积分布函数Pf_(fitness)中的分布区间；

step4.选择该区间所对应的个体为新的个体；

step5.重复step2-step4，直到满足种群中个体的上限。

Pf_(fitness)函数在计算过程中可根据给定的随机数自动生成。

(2)交叉操作：

step1.给定交叉概率Pc_(cross)；

step2.利用随机函数生成一个[0，1]区间的随机数ran，判断ran是否超过交叉概率；

step3.如果通过，利用随机函数选择两个随机的进行交叉的个体A，B；

step4.利用随机函数选择交叉的染色体位置i，交换t＝A.c_i，A.c_i＝B.c_i，B.c_i＝t；

step5.重复step2-step4，直到满足种群中个体的上限。

(3)变异操作：

step1.给定变异概率Pm_(mutation)；

step2.利用随机函数生成一个[0，1]区间的随机数ran，判断ran是否超过变异概率；

step3.如果通过，利用随机函数选择一个随机的进行变异的个体A；

step4.利用随机函数选择变异的染色体位置i，进行取反操作，A.c_i＝～A.c_i；

step5.重复step2-step4，直到满足种群中个体的上限。

通过遗传算子的调整，可以得到一组染色体，这组染色体的0/1编码决定了是否对该位置处的根的值进行调整。

接下来，根据染色体的编码调整根的数值包括：

step1.遍历整个染色体，；

step2.染色体编码为0，复制该位置处对应的根；

step3.染色体编码为1，调整该位置处对应的根，

设原有根为γ_k＝(a+bi)，调整为γ′_k＝(a+bi)/|(γ_k)|²，同时记录调整因子

step4.重复step2-step3，直到满足种群中个体的上限，调整常数项系数a₀＝a₀/|a′_k|。a、b为复数根γ_k的实部和虚部，a₀是原有系数的常数项。

最后，步骤(4)根据根组合，反向求解多项式系数，该系数就是子波。

调整子波的生成和最终优化后子波的生成是相同的，都是根据给定的根，反向求解出多项式系数，该系数就是调整后的子波。

在本发明中使用递推法：

设多项式的n个根为γ_i，i∈[0，n]，则多项式P(x)可以表示为：

设第k个根对应的多项式为：

则有：

4.1初始值：

对于k＝1，有P₁(x)＝-γ₁+x，有a₀＝-γ₁，a₁＝1。

4.2实数根的递推公式F1：

4.3复数根的递推公式F2：

因为是实系数方程，因此复数根以共轭复数根的方式成对出现，

设相邻两个复数根为γ_k+1＝(c+di)，γ_k+2＝(c-di)，则，

(x-γ_k+1)*(x-γ_k+2)＝(x-(c+di))*(x-(c-di))＝(x²-2cx+(c²+d²))

其中，c、d为复数根γ_k+1的实部和虚部。

带入递推公式F1，可得针对复数根的递推公式F2，

4.4求解多项式系数：

对于多项式来说，同时存在实数根和复数根，因此，可以先利用公式F1计算实数根的多项式系数，然后在利用公式F2继续计算复数根的多项式系数。在整个求解过程中，多项式的各个系数值是进行累加，因此计算先后顺序不会影响结果的正确性。

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

下面结合附图对本发明作进一步详细描述。图2是输入的地震数据和反射系数的示意图，包括该井在地震剖面上的投影。图3是通过复赛谱提取的初始子波，该子波是一个零相位子波。图4是多个通过本发明提取的混合相位子波，因为算法中使用到的随机函数，因此，生成的子波存在一定的差异，但是总体的趋势保持不变。

本实施例是本发明用于某地区的实际地震解释资料提取一些地震子波的实例，基于遗传算法的子波提取的步骤如下：

(1)如图2所示，输入地震数据和反射系数；

(2)利用复赛谱，提取初始子波，初始子波为零相位子波，如图3；

(3)分解该零相位子波为多项式的根，对根进行编码，然后利用遗传算法，对根进行迭代优化，通过目标函数，计算出迭代过程中适应度函数值最大值，利用该值对应的根，反向计算出多项式系数，该系数值作为新提取的子波。如图4a、4b、4c、4d所示，进行多次迭代计算，得到多个满足条件的解。

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

Claims (10)

1.一种子波提取方法，其特征在于，该方法包括：

对地震道进行频谱分解，利用地震道的频谱特征，获取初始子波；

将获取的初始子波看成多项式函数的系数，计算该函数的特征值；

建立目标函数，利用遗传算法，求取最优解；

根据求取的最优解组合，反向求解子波。

2.根据权利要求1所述的子波提取方法，其特征在于，通过复赛谱提取初始子波；

s(t)＝w(t)*r(t)+n(t)； (1)

其中，s(t)表示地震信号，w(t)表示子波，r(t)表示反射系数序列，n(t)表示附加干扰；

使用傅里叶变换，转换到频率域：s(ω)＝w(ω)r(ω)； (2)

取对数：ln s(ω)＝lnw(ω)+lnr(ω)； (3)

使用反傅里叶变换，转换回时间域：sp(t)＝wp(t)+rp(t)； (4)

则sp(t)、wp(t)、rp(t)是s(t)、w(t)、r(t)的复赛谱。

3.根据权利要求1所述的子波提取方法，其特征在于，将初始子波看成多项式的一组系数，即：

其中，a_i，i∈[0,n]是子波的每个样点值；

设多项式的n个根为γ_i，i∈[0,n]，则多项式P(x)表示为：

令P(x)＝0，利用QR方法求得多项式方程的所有实根，得到不同的多项式系数。

4.根据权利要求1所述的子波提取方法，其特征在于，建立目标函数：

Fitness＝MAX{R_(l)} (7)

其中相关系数R_(l)的公式为：

其中，x长度为n，为地震数据；y长度为m，为合成记录数据；l为是相关系数数组的下标；t为输入地震数据数组的下标。

5.根据权利要求1所述的子波提取方法，其特征在于，利用遗传算法，求取最优解包括：

步骤1：进行染色体编码并初始化；

步骤2：计算当前种群中的染色体的适应度；

步骤3：调整染色体的编码，产生新的染色体；

步骤4：根据染色体的编码调整根的数值；

步骤5：根据调整后的根生成调整子波，并计算新个体的适应度；

步骤6：将新的个体插入种群，并淘汰适应度差的个体；

步骤7：如果达到规定代数，返回最优染色体，否则返回步骤3。

6.根据权利要求5所述的子波提取方法，其特征在于，步骤1中染色体编码包括将染色体表示成如下向量：

C＝(c₀,c₁,c₂,c₃,c₄……c_i……c_n),c_i∈[0,1]，取值0表示不进行调整，取值1表示进行调整。

7.根据权利要求5所述的子波提取方法，其特征在于，步骤3中调整染色体的编码包括：选择、交叉、变异操作。

8.根据权利要求5所述的子波提取方法，其特征在于，步骤4中根据染色体的编码调整根的数值，包括：

步骤4.1：遍历整个染色体；

步骤4.2：染色体编码为0，复制该位置处对应的根；

步骤4.3：染色体编码为1，调整该位置处对应的根；

设原有根为γ_k＝(a+bi)，调整为γ′_k＝(a+bi)/|(γ_k)|²，同时记录调整因子

步骤4.4：重复步骤4.2-步骤4.3，直到满足种群中个体的上限，调整常数项系数a₀＝a₀/|a′_k|。

9.根据权利要求1所述的子波提取方法，其特征在于，使用递推法，反向求解多项式系数，生成子波。

10.一种子波提取系统，其特征在于，该系统包括：

存储器，存储有计算机可执行指令；

处理器，所述处理器运行所述存储器中的计算机可执行指令，执行以下步骤：

对地震道进行频谱分解，利用地震道的频谱特征，获取初始子波；

将获取的初始子波看成多项式函数的系数，计算该函数的特征值；

建立目标函数，利用遗传算法，求取最优解；

根据特征值组合，反向求解子波。

Images