CN102313901B - 一种初至波迭代拾取的方法 - Google Patents

Abstract

本发明是石油物探资料处理的初至波迭代拾取的方法，计算每一道每个采样点处的能量比值，得到能量比值曲线，确定初至预测值，以此为起点找到其中最大振幅值，计算相邻道的初至时间差的绝对值和的平均值，取一道和前后两道初至时间，排列计算可信度系数，对可信度系数为0的初至利用前后相邻4个可靠初至进行拉格朗日多项式拟合插值，得到初至预测值。本发明对于初至波形变化较大、各种波相互干扰，同时外界干扰等的影响的复杂地区，可以获得准确的初至时间，解决了初至拾取方法的难题。

Description

一种初至波迭代拾取的方法

技术领域

本发明涉及石油地球物理勘探资料处理技术，是一种初至波迭代拾取的方法。

背景技术

石油地球物理勘探资料处理中的初至拾取，早期通过手工拾取各道的初至时间，这种方法工作量往往非常庞大，道数甚至超过一万道，这样全部用手工拾取就成为一项费时、费工、十分繁杂的工作。随着地震勘探的不断深入和数据处理精度的不断提高，准确快速地拾取地震波初至是非常重要的。在折射波静校正、VSP解释、浅层折射波勘探和层析成像勘探中，常用初至波的旅行时来反演表层速度和厚度等介质参数，特别是在一些表层较为复杂的地区，如黄土塬、沙漠等，初至拾取是各种初至波静校正的基础。早期的初至拾取方法已不能满足需求。

近年来，公开了不同的初至波自动拾取方法。早期Peraldi和Clement公开了拐点校正法，将地震道与参考道进行互相关，归一化的互相关函数极大值对应的时移即初至时间；1980年Hatherrly公开了线性最小平方预测技术和拐点校正相结合的方法，提出首先识别第一个峰值和拐点，然后估算二者的统计差值；1983年Gelchinsky和Shtivelman公开了一种相关法，假设各道的脉冲形状不发生变化，然而实际情况与假设不符；1987年黄成之等采用统计方法将地震初至波记录分成信号和噪音两个部分，并使这两部分统计特征之间的差别为很大；1990年陈遵德等采用模糊检测方法确定初至波的隶属函数进行初至波检测；1991年Lin Zhang采用相似度准则，建立地震波同相轴的目标函数，用最短路径算法求取目标函数极小，达到自动拾取地震剖面同相轴的目的，但此方法建立目标函数的工作量较大；1994年庄东海等采用将地震记录初至拾取看做一个模式识别过程，充分利用地震记录的瞬时(局部)特征和整体特征，用神经网络方法进行地震记录初至拾取；1995年江玉乐等提出同极性能量比值法，即改进的能量比值法；1996年Fabio Boschetti等提出了一种基于分形维的初至检测算法，该方法是基于在地震道随着信号的出现其分形维就发生变化的特征来确定地震道初至；新近发展起来的BP神经网络法利用多参数特征进行模式识别，充分利用地震记录的瞬时(局部)特征和整体特征来确定地震道初至。

2006年潘树林等利用单道边界检测公式和样条插值进行拾取初至，即：

$R\left(r\right)=\left|\frac{N}{M}(N-M)\right|---\left(1\right)$

$M=\underset{p=r}{\overset{r-n}{\mathrm{\Σ}}}{x}_p$ $N=\underset{p=r}{\overset{r+n}{\mathrm{\Σ}}}{x}_p$

式中：R(r)为第r个样点的边界特征，R(r)越大，代表r点的边界特征越强，最大边界特征值所对应的点就是能量边界。该能量比值公式对初至波到达前的干扰波能量变化过于敏感，对于较差的资料初至拾取误差较大。

综上所述，初至波自动拾取的方法有很多，对于地震记录背景较平静且初至波形简单地区的资料，用一般的初至拾取方法都能满足要求，但对于复杂地区，由于初至波形变化较大，各种波相互干扰，同时外界干扰等的影响，常常难以获得准确的初至时间，一直是初至拾取方法的难题。

发明内容

本发明目的是提供一种消除初至波到达前随机干扰波的影响，有效地避免误拾取续至波的初至波迭代拾取的方法。

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

1)激发并采集记录单炮地震数据。

2)利用以下公式计算每一道每个采样点处的能量比值：

$R\left(r\right)=\frac{{\left(\frac{N}{M+A^{2}W}\right)}^4\×{(M-N)}^2}{r^2}---\left(2\right)$

$M=\underset{p=1}{\overset{r}{\mathrm{\Σ}}}{x}_p^2$ $N=\underset{p=r}{\overset{r+W}{\mathrm{\Σ}}}{x}_p^2$ $A=\frac{1}{n}\underset{r=1}{\overset{n}{\mathrm{\Σ}}}\left|x_r\right|$

式中：每一道的第一个采样点到当前采样点的长度为长时窗，r为长时窗的结束样点，n为一道的采样点数，N为当前点之后W个采样点的能量之和，W为短时窗的计算点数，A为所用道的所有点振幅绝对值和的平均值，M为所用道的长时窗内所有采样点的能量之和；

步骤2)所述的计算中短时窗的计算点数W为20个采样点，采样点为短时窗的长度。

3)取任一道，利用步骤2)中公式(2)计算该道每个采样点处的能量比值，得到能量比值曲线，在能量比值曲线上确定门槛值对应的时间点，在对应的原始地震道上该时间点为初至预测值；

上述门槛值为大于0小于1，缺省为0.5，初至越弱，门槛值取的越小，初至越强，门槛值取的越大。

4)以初至预测值为起点，长度缺省30个采样点为第一时窗，找到其中最大振幅值。

修正波峰点的能量参数与最大振幅的乘积，作为振幅门槛值。

以初至预测值为起点，长度是第一时窗的1/2为第二时窗，第二时窗内第一个振幅大于振幅门槛值的位置，就是所要找的初至波峰点，如果在第二个时窗内找不到满足条件的振幅大于振幅门槛值的位置，则该道初至不存在；

步骤4)所述的修正波峰点的能量参数大于0小于1，缺省能量参数取0.5，初至越弱，该能量参数取的越小，反之，初至越强，该能量参数取的越大；

步骤4)所述的第一时窗长度，如果实际初至波峰点在确定的初至波峰点的下方，增加第一时窗的长度，反之，减少第一时窗的长度；

步骤4)所述的初至波峰点的确定，对于三维地震资料，每个排列都按照相同的方法处理，得到每一炮的初至波峰点。

5)计算相邻道的初至时间差的绝对值和的平均值为：

式中：t₁，t₂，...，t_m是第一道到最后一道的初至时间，m为一炮的总道数；

以五道为一组，任一道初至时间t_i，i为道号，i大于等于1小于等于m，满足下面条件中一个，该道初至的可信度系数为1，如果均不满足，该道初至的可信度系数为0：

第一个条件：j＝i-4，i-3，i-2，i-1；

第二个条件：

j＝i-3，i-2，i-1，i+1；

第三个条件：

j＝i-2，i-1，i+1，i+2；

第四个条件：

j＝i-1，i+1，i+2，i+3；

第五个条件：

j＝i+1，i+2，i+3，i+4；

t_j为t_i相邻道的初至时间，j为道号，j大于等于1小于等于m，门槛值k为大于0小于等于2，缺省情况下为1。

步骤5)中所述的门槛值k，对信噪比高的资料，该门槛值参数取的大一些，反之，该门槛值参数取的小一点；

步骤5)所述的参数对于三维地震资料，每个排列利用相邻道的初至时间差的绝对值和的平均得到一个平均值，参数

的值为将所有排列得到的平均值相加除以排列个数，利用该参数计算每道初至的可信度系数。

6)每炮除了前两道和后两道，任取一道初至时间t_i，对该道和相邻的前后四道初至时间t_i-2，t_i-1，t_i，t_i+1，t_i+2排序，如果该道排序后的下标与排序前下标的差的绝对值等于2，该道初至的可信度系数为0，否则该道初至的可信度系数为1。

步骤6)所述的排序对于三维地震资料，每个排列都按照上述方法计算可信度系数。

7)通过步骤5)和步骤6)，每一炮除了前两道和最后两道外每道初至都得到两个可信度系数，对于任一道，如果得到的两个可信度系数相同，该道初至可信度系数不变，如果两次得到的可信度系数不同，则该道初至的可信度系数为0，对于每炮的前两道和后面两道取步骤5)得到可信度系数，每一道可信度系数等于1的初至为可靠初至。

8)对可信度系数为0的初至利用前后相邻4个可靠初至进行拉格朗日多项式拟合插值，得到初至预测值；

对于边缘的可信度系数为0的初至，则利用二次曲线拟合插值，得到边缘的初至预测值；

初至预测值向上找最近的波峰点，向下找第一个波峰点，如果向上最近波峰点的振幅值小于向下第一个波峰点的振幅值和步骤4)中修正波峰点的能量参数的乘积，初至预测值下方的第一个波峰点为初至波峰点，否则初至预测值向上最近的波峰点为初至波峰点。

本发明对于初至波形变化较大、各种波相互干扰，同时外界干扰等的影响的复杂地区，可以获得准确的初至时间，解决了初至拾取方法的难题。

附图说明

本发明附图说明如下：

图1是陆上某地区原始单炮记录；

图2(a)为图1炮第5道原始地震道；

图2(b)为第5原始地震道用能量比值公式得到的能量比值曲线；

图3为图1所示的单炮记录初步拾取波峰点的效果图；

图4是对图3所示初至波峰点通过确定可靠初至效果图；

图5是图1第7道原始地震道重新寻找初至波峰点示意图；

图6是对图4二次拾取效果图。

具体实施方式

本发明是根据能量比值公式和初至波频率随炮检距、检波点和炮点的变化，利用能量参数和门槛值确定初至波峰点的位置，然后对找到的初至波峰点利用相邻道之间关系和排序计算可信度系数，从而确定可靠初至，最后对可信度系数为0的初至利用可靠初至作为约束条件进行拟合插值，重新拾取可信度系数为0的初至的方法。

本发明提供的技术主要是利用能量比值公式计算每一道能量比值曲线、利用能量参数和门槛值确定初至波峰点的位置、计算可信度系数、确定可靠初至和重新拾取。

1)激发并采集记录单炮地震数据，如图1，该炮220道。

2)利用以下公式计算每一道每个采样点处的能量比值：

$R\left(r\right)=\frac{{\left(\frac{N}{M+A^{2}W}\right)}^4\×{(M-N)}^2}{r^2}---\left(2\right)$

$M=\underset{p=1}{\overset{r}{\mathrm{\Σ}}}{x}_p^2$ $N=\underset{p=r}{\overset{r+W}{\mathrm{\Σ}}}{x}_p^2$ $A=\frac{1}{n}\underset{r=1}{\overset{n}{\mathrm{\Σ}}}\left|x_r\right|$

式中：每一道的第一个采样点到当前采样点的长度为长时窗，r为长时窗的结束样点，M为所用道的长时窗内所有采样点的能量之和；n为一道的采样点数；N为当前点之后W个点的能量之和，W为短时窗的计算点数；A为所用道的所有点振幅绝对值和的平均值；(N-M)²是压制初至波到达前的随机干扰；

加强初至边界能量；A²W是在实际计算中，初至点以上少数样点的幅值有可能为一些接近零的点，从而导致M趋于0，提高拾取的稳定性；

降低了续至波的能量比值，从而避免误拾取续至波；

计算中短时窗的计算点数W取20个采样点。

3)利用步骤2)中公式(2)，该炮的每一道都得到一个能量比值曲线，门槛值取0.5。对第5道，图2(a)为第5道原始地震道，图2(b)为第5道通过能量比值公式(2)得到的能量比曲线，在该道能量比曲线上门槛值0.5对应的时间点为1142ms，在对应的原始地震道上该时间点为初至预测值。

4)以初至预测值1142ms为起点，向后搜索30个采样点，找出最大振幅值0.25，取修正波峰点的能量参数为0.5，初至预测值1142ms向后15个采样点中第一个振幅大于0.25和0.5乘积的点为该道上第290个采样点，即0.163＞0.125，该采样点就是第5道所要找的初至波峰点的位置，该炮的每一道都用此方法进行寻找初至波峰点，得到如图3所示的拾取初至波峰点示意图。

5)计算相邻道的初至时间差的绝对值和的平均值

为：

该炮第一道到最后一道的初至时间为：

t₁＝1172ms，t₂＝1168ms，t₃＝1164ms，...，t₂₂₀＝1140ms则

从第一道开始，以每五道为一组，判断每一道与前后相邻道之间的关系，取第5道，t₅前后8道的初至时间为：

t₁＝1172ms，t₂＝1168ms，t₃＝1164ms，t₄＝1160ms，t₅＝1156ms，

t₆＝1148ms，t₇＝1184ms，t₈＝1084ms，t₉＝1132ms

门槛值k为取0.4，满足下面条件中一个，第5道初至的可信度系数为1，如果均不满足，第5道初至的可信度系数为0：

第一个条件：|t₅-t_j|＜11.288，j＝1，2，3，4；

第二个条件：|t₅-t_j|＜11.288，j＝2，3，4，6；

第三个条件：|t₅-t_j|＜11.288，j＝3，4，6，7；

第四个条件：|t₅-t_j|＜11.288，j＝4，6，7，8；

第五个条件：|t₅-t_j|＜11.288，j＝6，7，8，9；

t_j为t₅相邻的初至时间，j为道号。

通过计算知，第一个条件和第二个条件均满足，这样该道初至的可信度系数为1。

利用门槛值参数k＝0.4对每一道初至确定可信度系数。

6)该炮除了前两道和后两道，每道与该道相邻的前面两道和后面两道共五道为一组，判断该道与它相邻这四道之间的顺序关系确定可信度系数。对于第8道初至时间t₈，t₈和前后相邻四道为：

t₆＝1148ms，t₇＝1184ms，t₈＝1084ms，t₉＝1132ms，t₁₀＝1124ms，

初至时间t₈排序前与前后道下标顺序为：6，7，8，9，10，对初至时间t₆，t₇，t₈，t₉，t₁₀进行由小到大排序，排序后五道初至时间顺序为t₈，t₁₀，t₉，t₆，t₇，得到的下标顺序为：8，10，9，6，7。

这样第8道排序前下标为8，排序后下标在排序前下标6的位置，排序前后下标差的绝对值等于2，该道初至的可信度系数为0，其它道用类似的方法进行判断，除了前两道和后两道每道都又得到一个可信度系数。

7)通过步骤5)和步骤6)，该炮除了前两道和后两道外每道初至都得到两个可信度系数，对于任一道，如果得到的两个可信度系数相同，该道的可信度系数不变，如果得到的可信度系数不同，则该道的可信度系数为0。前两道和最后两道通过步骤5)知可信度系数均为1，这样该炮每道初至都得到一个可信度系数为：

1，1，1，1，1，1，0，0，1，1，…，1，1，1

每一道可信度系数等于1的初至为可靠初至。

8)从第一道开始寻找可信度系数为0的初至，出现可信度系数为0的道为第7道和第8道，后面的可靠初至为第9道和第10道，利用相邻的可靠初至第5道、第6道、第9道和第10道对应的初至时间t₅，t₆和t₉，t₁₀构造拉格朗日多项式，对第7道和第8道进行拟合插值，得到初至预测值。

第7道初至预测值是该道的第282个采样点处，该位置向上最近的波峰点为第259个采样点，振幅为0.064，向下第一个波峰点为第285个采样点，振幅为0.159，步骤4)中的修正波峰点的能量参数为0.5，由于0.064小于0.5×0.159＝0.0785，该道初至波峰点的位置在第285个采样点处，如图5所示。第8道按照同样方法确定初至波峰点的位置。

按照上面的方法对整炮进行拟合插值，重新确定可信度系数为0的初至的位置，得到如图6所示。

Claims (9)

1.一种初至波迭代拾取的方法，特点是通过以下技术方案实现：

1）激发并采集记录单炮地震数据；

2）利用以下公式计算每一道每个采样点处的能量比值：

$R\left(r\right)=\frac{{\left(\frac{N}{M+A^{2}W}\right)}^4\×{(M-N)}^2}{r^2},$

$M=\underset{p=1}{\overset{r}{\mathrm{\Σ}}}{x}_p^2,$ $N=\underset{p=r}{\overset{r+W}{\mathrm{\Σ}}}{x}_p^2,$ $A=\frac{1}{n}\underset{r=1}{\overset{n}{\mathrm{\Σ}}}\left|x_r\right|,$

式中：每一道的第一个采样点到当前采样点的长度为长时窗，r为长时窗的结束样点，n为一道的采样点数，N为当前点之后短时窗所有采样点的能量之和，W为短时窗的计算点数，A为所用道的所有点振幅绝对值和的平均值，M为所用道的长时窗内所有采样点的能量之和；

3）取任一道，利用步骤2）中公式计算该道每个采样点处的能量比值，得到能量比值曲线，在能量比值曲线上确定门槛值对应的时间点，在对应的原始地震道上该时间点为初至预测值；

4）以初至预测值为起点，长度缺省30个采样点为第一时窗，找到其中最大振幅值，修正波峰点的能量参数与最大振幅的乘积，作为振幅门槛值，以初至预测值为起点，长度是第一时窗的1/2为第二时窗，第二时窗内第一个振幅大于振幅门槛值的位置，就是所要找的初至波峰点，如果在第二时窗内找不到满足条件的振幅大于振幅门槛值的位置，则该道初至不存在；

5）计算相邻道的初至时间差的绝对值和的平均值

式中：t₁,t₂,...,t_m是第一道到最后一道的初至时间，m为一炮的总道数；

以五道为一组，任一道初至时间t_i，i为道号，i大于等于1小于等于m，满足下面条件中一个，该道初至的可信度系数为1,如果均不满足，该道初至的可信度系数为0：

第一个条件：

j＝i-4,i-3,i-2,i-1；

第二个条件：

j＝i-3,i-2,i-1,i+1；

第三个条件：

j＝i-2,i-1,i+1,i+2；

第四个条件：

j＝i-1,i+1,i+2,i+3；

第五个条件：

j＝i+1,i+2,i+3,i+4；

t_j为t_i相邻道的初至时间，j为道号，j大于等于1小于等于m，门槛值k为大于0小于等于2，缺省情况下为1；

6）每炮除了前两道和后两道，任取一道初至时间t_i，对该道和相邻的前后四道初至时间t_i-2,t_i-1,t_i,t_i+1,t_i+2排序，如果该道排序后的下标与排序前下标的差的绝对值等于2，该道初至的可信度系数为0，否则该道初至的可信度系数为1；

7）通过步骤5）和步骤6），每一炮除了前两道和最后两道外每道初至都得到两个可信度系数，对于任一道，如果得到的两个可信度系数相同，该道初至可信度系数不变，如果两次得到的可信度系数不同，则该道初至的可信度系数为0，对于每炮的前两道和后面两道取步骤5）得到可信度系数，每一道可信度系数等于1的初至为可靠初至；

8）对可信度系数为0的初至利用前后相邻4个可靠初至进行拉格朗日多项式拟合插值，得到初至预测值；

对于边缘的可信度系数为0的初至，则利用二次曲线拟合插值，得到边缘的初至预测值；

初至预测值向上找最近的波峰点，向下找第一个波峰点，如果向上最近波峰点的振幅值小于向下第一个波峰点的振幅值和步骤4）中修正波峰点的能量参数的乘积，初至预测值下方的第一个波峰点为初至波峰点，否则初至预测值向上最近的波峰点为初至波峰点。

2.根据权利要求1所述的方法，特点是步骤2)所述的计算中短时窗的计算点数W为20个。

3.根据权利要求1所述的方法，特点是步骤3)所述的门槛值为大于0小于1，缺省为0.5，初至越弱，门槛值取的越小，初至越强，门槛值取的越大。

4.根据权利要求1所述的方法，特点是步骤4）所述的修正波峰点的能量参数大于0小于1，缺省能量参数取0.5，初至越弱，该能量参数取的越小，反之，初至越强，该能量参数取的越大。

5.根据权利要求1所述的方法，特点是步骤4）所述的第一时窗长度，如果实际初至波峰点在确定的初至波峰点的下方，增加第一时窗的长度，反之，减少第一时窗的长度。

6.根据权利要求1所述的方法，特点是步骤4）所述的初至波峰点的确定，对于三维地震资料，每个排列都按照相同的方法处理，得到每一炮的初至波峰点。

7.根据权利要求1所述的方法，特点是步骤5）中所述的门槛值k，对信噪比高的资料，该门槛值参数取的大一些，反之，该门槛值参数取的小一点。

8.根据权利要求1所述的方法，特点是步骤5）所述的参数对于三维地震资料，每个排列利用相邻道的初至时间差的绝对值和的平均得到一个平均值，参数

的值为将所有排列得到的平均值相加除以排列个数，利用该参数

计算每道初至的可信度系数。

9.根据权利要求1所述的方法，特点是步骤6）所述的排序对于三维地震资料，每个排列都按相同的方法计算可信度系数。

Images