一种低信噪比地震初至信号自动增强方法及系统与流程

through sliding windows and fuzzy c-means",mathematics,vol.7,pp.221,2019.
11.[6]an,s.,hu,t.,cui,y.,duan,w.,and peng,g.,auto-pick first breaks with complex raypaths for undulate surface conditions.applied geophysics,2015,12(1):93
–
100.
[0012]
[7]宋龙龙,邹志辉,黄忠来,基于相邻虚拟道叠加的超虚折射干涉法及其在广角obs折射波增强中的应用.地球物理学报,2019,62(03):993
–
1006.


技术实现要素：

[0013]
本技术提供一种低信噪比地震初至信号的自动增强方法及系统，综合利用三维地震数据信息，实现针对低信噪比三维地震资料的初至自动增强，提高初至信噪比的提升效果。本发明的核心是：本发明对累积量干涉法进行了发展和扩充，通过高阶累积量和叠加得到两检波点间的虚折射记录，再通过多道褶积和叠加得到炮检点(炮点和检波点)间的折射波初至。本发明方法使用高阶累积量克服互相关对高斯噪声敏感的问题，提高初至信噪比的提升效果；相比svi方法所用的互相关使用两道数据进行道间时差估计，只能利用二维地震数据信息，本发明方法所用的高阶累积量引入多方向地震数据进行道间时差估计，综合利用了三维地震数据信息，并实现本发明方法应用于三维地震数据的处理。
[0014]
为解决上述技术问题，本技术实施例提供一种低信噪比地震初至信号自动增强方法及系统是这样实现的：
[0015]
采用一种基于累积量干涉法(cumulant seismic interferometry,ci)进行扩展的方法代替现有传统的svi方法，实现初至自动增强，并适用于低信噪比三维地震资料；所述低信噪比地震初至信号自动增强方法包括：通过高阶累积量得到检波点间的初至波，通过同相位叠加提高初至波信噪比，再通过多道褶积将检波点间的初至波恢复成炮检点间的初至波，再次通过同相位叠加提高信噪比。实现本发明所述方法的系统具体包括如下模块或步骤：
[0016]
a.去噪模块：对输入地震数据滤波去除异常噪声干扰；
[0017]
b.初至截取模块：利用时窗截取包含初至波的局部地震数据；
[0018]
c.初至增强模块：应用ci方法增强初至波；
[0019]
包括绕射波、折射波和直达波在内的初至波信号均可通过ci方法增强；反射波、背景噪声无法实现同相位叠加，被相对压制；应用ci方法增强初至波包括如下过程：
[0020]
c1.使用不同炮点，通过计算高阶累积量得到检波点间的初至波；
[0021]
c2.将使用不同炮点得到的检波点间的初至波叠加，可以实现同相位叠加；
[0022]
具体实施时，对于绕射波，通过高阶累积量和叠加，得到检波点间的虚绕射记录；再计算检波点间的虚绕射，得到叠加后的累积量函数。
[0023]
c3.通过多道褶积将不同检波点间的初至波恢复得到炮检点间的初至波；
[0024]
c4.将使用不同检波点得到的炮检点间的初至波叠加，可以实现同相位叠加；
[0025]
其中，应用ci方法增强初至波的初至增强模块中，具体以初至波中的绕射波为例，通过高阶累积量和叠加，得到检波点间的虚绕射记录；包括如下过程：
[0026]
c11.设{x0(t),x1(t),...,x
p
(t),,...,x
p-1
(t)}代表p个零均值的地震道，p阶累积量可表示如下：
[0027]
c(τ1,τ2,...,τ
p
,....,τ
p-1
)＝cum{x0(t),x1(t τ1),...,x
p
(t τ
p
),...,x
p-1
(t τ
p-1
)},
ꢀꢀꢀꢀꢀꢀ
(1)
[0028]
其中，cum{.}是累积量算子，c(.)是累积量函数，τ
p
(p＝1,
…
,p-1)表示地震道x
p
(t)与x0(t)之间的时延，p阶累积量包含p-1个时延；
[0029]
c12.计算检波点间的虚绕射，得到叠加后的累积量函数，记为c
1,2,...n
(τ2,τ3,...,τn)；
[0030]
将在炮点sk，k＝1,
…
,m，m表示炮点总数量，激发检波点ri(i＝1,
…
,n，n表示检波点总数量)接收的地震记录表示为u
k；i
(t)；针对地震记录包含绕射波和随机噪声的情况，即u
k；i
(t)＝g
k；i
(t) w(t)，对炮点sk在检波点r1,...,rn的地震记录计算得到高阶累积量[c
1,2,...,i,...n
(τ2,τ3,...,τi,...,τn)]k；
[0031]
时延τi表示为式(3)：
[0032][0033]
τi表示ri到r1之间的绕射波道间时差；和分别为炮点sk到d和d到ri的旅行时；对于不同检波点ri接收的绕射波信号，始终为共用的旅行时；使用不同炮点sk得到的[c
1,2,...,i,...n
(τ2,τ3,...,τi,...,τn)]k包含相同的时延；
[0034]
将[c
1,2,...,i,...n
(τ2,τ3,...,τi,...,τn)]k全部叠加，进行绕射波信号的同相位叠加，得到叠加后的累积量函数c
1,2,...n
(τ2,τ3,...,τn)，其中包含虚绕射记录；
[0035]
c2.通过多道褶积将虚绕射记录恢复成炮检点间的绕射波记录，包括如下过程：
[0036]
将累积量函数c
j,1,...,n
(τ1,...,τn)与原始地震记录u
s；j
(t)进行多道褶积，重构炮点s到检波点ri(i＝j,1,2,
…
,n)的绕射波记录，表示为式(4)：
[0037][0038]
其中，中的符号“^”表示该绕射波是重构得到的，区别于原始地震记录中的绕射波，下标j表示该重构信号是使用检波点rj的地震记录得到；c
j，1,2，...n
(t
1-τ1,...,t
n-τn)中包含n组检波点间的时延，使用n道原始地震记录u
s；j
(t)对每个时延进行褶积即多道褶积，将时延全部重构成炮检点间的旅行时；
[0039]
重构的地震记录中旅行时ti(i＝1,
…
,n)表示为式(5)：
[0040][0041]
即重构炮点s在检波点ri接收的绕射波记录；
[0042]
使用不同检波点得到的重构记录中包含相同的绕射波，通过进行同相位叠加，得到叠加后的重构绕射波记录；
[0043]
绕射波、折射波和直达波在内的初至波信号均可通过ci方法增强；反射波、背景噪声无法实现同相位叠加，被相对压制；
[0044]
d.结果检验模块：检查经步骤c4实现同相位叠加得到的初至波是否达到了满意的初至波信噪比提升效果，如果未达到，则以步骤c4实现同相位叠加得到的初至波作为步骤c新的输入数据，重复步骤c再次进行初至增强，得到加强初至波，直到实现满意的信噪比提升效果。
[0045]
具体通过基于能量比的长短窗平均(short term averaging/long term averaging,sta/lta)方法自动拾取初至；如果初至拾取精度达到后续的层析或折射静校正所要求的精度，则视为达到满意的初至波信噪比提升效果，否则迭代执行步骤c直到实现满意的信噪比提升效果。
[0046]
本发明的有益效果是：
[0047]
本发明提供一种基于高阶累积量的低信噪比地震资料初至自动增强技术，通过高阶累积量和叠加得到检波点间的初至波，再通过多道褶积和叠加得到增强后的炮检点间初至波。与现有技术相比，本发明具有以下技术优势：
[0048]
(一)本发明方法完全数据驱动，不依赖频率或视速度信息，以原始地震数据作为输入即可自动实现初至波增强。
[0049]
(二)本发明方法使用高阶累积量克服互相关对高斯噪声敏感的问题，提高初至信噪比的提升效果。
[0050]
(三)相比传统svi方法所用的互相关使用两道数据进行道间时差估计，只能利用二维地震数据信息，本发明方法所用的高阶累积量引入多方向地震数据进行道间时差估计，综合利用了三维地震数据信息，实现本发明方法应用于三维地震数据的处理。
附图说明
[0051]
为了更清楚地说明本技术实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本技术中记载的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动性的前提下，还可以根据这些附图获得其他的附图。
[0052]
图1是本发明系统的模块示意图。
[0053]
图2是通过高阶累积量和叠加得到检波点间的虚绕射记录；sk(k＝1,
…
,m)是炮点，ri(i＝1,
…
,n)是检波点，d为绕射点，符号cum表示计算高阶累积量；(a)为用于计算高阶累积量的地震记录，其中灰色粗线表示地震记录共用的旅行时，(b)为高阶累积量结果：包含了检波点间的虚绕射记录，与炮点有关的共用旅行时被抵消，实(虚)线表示正(负)的旅行时，(c)为叠加结果：将使用不同炮点sk得到的高阶累积量结果叠加，实现绕射波信号的同相位叠加，增强信噪比。
[0054]
图3是ci方法通过多道褶积和叠加恢复炮检点间的绕射波；图中，s是炮点，ri(i＝j,k,l,1,
…
,n)是检波点，d为绕射点，符号*是褶积符号；(a)是原始地震记录与累积量函数做多道褶积，灰色粗线代表共用的旅行时，实(虚)线表示正(负)的旅行时，(b)是多道褶积的结果，包含了相同的炮检点间虚绕波，(c)是叠加结果：将使用不同检波点得到的多道褶积结果叠加，实现绕射波的同相位叠加，进一步增强信噪比。
[0055]
图4是二维合成数据实例；(a)是截取了初至同相轴的原始数据，(b)是加入了非相关和相关的高斯噪声的数据，(c)是svi方法的结果，(d)是采用ci方法得到的结果。
[0056]
图5是三维实际数据实例；其中，白(黑)色表示波峰(谷)位置，虚线表示初至波截取时间的边界，(a)-(c)中的方框所圈的区域用于局部放大：(a)是原始数据，(b)是f-x滤波结果，(c)是ci方法的结果。
具体实施方式
[0057]
本技术实施例提供一种低信噪比地震初至信号增强方法及系统。
[0058]
为了使本技术领域的人员更好地理解本技术中的技术方案，下面结合附图，通过实施例进一步描述本发明，所描述的实施例仅仅是本技术一部分实施例，而不是全部的实施例，但不以任何方式限制本发明的范围。基于本技术中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都应当属于本技术保护的范围。
[0059]
本技术实施例的目的是提供一种低信噪比地震初至信号增强方法及系统，实现针对低信噪比三维地震资料的初至自动增强。
[0060]
本发明申请实施例提供了一种低信噪比地震初至信号自动增强方法及系统，基于高阶累积量干涉法(ci)实现初至增强，适用于低信噪比三维地震资料。所述方法包括，通过高阶累积量和叠加得到检波点间的初至波，再通过多道褶积和叠加得到增强后的炮检点间初至波。基于ci方法实现初至自动增强系统的具体实施方式如下：
[0061]
a.去噪模块：对输入地震数据进行高通滤波和中值滤波滤波，去除异常噪声干扰，包括超低频的背景噪声、高频异常噪声等。
[0062]
b.初至截取模块：利用时窗截取包含初至波的局部地震数据(需包含全部初至波)，减少全数据参与运算所需要的计算量，时窗长度一般为3-4个地震波波长。
[0063]
c.初至增强模块：应用ci方法增强初至波：
[0064]
c1.以绕射波为例说明ci方法的实施过程，通过高阶累积量和叠加得到两检波点间的虚绕射记录：
[0065]
累积量函数是可以用于描述确定性信号的协方差函数，函数值可以由累积量生成函数的泰勒展开系数确定[8](mendel,j.,tutorial on higher-order statistics(spectra)in signal processing and system theory:theoretical results and some applications.proceedings of the ieee,1991,79(3):278
–
305)，假设{x0(t),x1(t),...,x
p
(t),,...,x
p-1
(t)}代表p个零均值的地震道，它的p阶累积量可表示如下：
[0066]
c(τ1,τ2,...,τ
p
,....,τ
p-1
)＝cum{x0(t),x1(t τ1),...,x
p
(t τ
p
),...,x
p-1
(t τ
p-1
)},
ꢀꢀꢀꢀꢀꢀ
(1)
[0067]
其中，cum{.}是累积量算子，c(.)是累积量函数，τ
p
(p＝1,
…
,p-1)表示地震道x
p
(t)与x0(t)之间的时延，p阶累积量包含p-1个时延，理论上可以通过计算p-1次道集间的互相关函数得到[8]。
[0068]
计算两检波点间虚绕射的物理过程如下。将在炮点sk(k＝1,
…
,m，m表示炮点总数量)激发检波点ri(i＝1,
…
,n，n表示检波点总数量)接收的地震记录表示为u
k；i
(t)，只考虑地震记录包含绕射波和随机噪声的简单情况，即u
k；i
(t)＝g
k；i
(t) w(t)，对炮点sk在检波点r1,...,rn的地震记录计算高阶累积量的公式如下：(图2b)
[0069][0070]
其中，[c
1,2,...,i,...n
(τ2,τ3,...,τi,...,τn)]k的下标k表示该累积量函数使用炮点sk有关的地震记录得到，w(t)是高斯随机噪声，因此包含w(t)的高阶累积量分项为零。绕射波信号g
k；j
(t)是由地震子波s(t)与绕射波介质格林函数褶积得到的，即
和分别为sk到d和d到ri的旅行时。对于不同检波点ri接收的绕射波信号，为始终为共用的旅行时，而计算累积量函数的过程，相当于求取两点间的时延，共用旅行时被相互抵消。因此，累积量函数中的时延τi为
[0071][0072]
τi表示了ri到r1之间的绕射波道间时差。注意到τi与炮点无关，使用不同炮点sk得到的[c
1,2,...,i,...n
(τ2,τ3,...,τi,...,τn)]k包含了相同的时延，因此，将[c
1,2,...,i,...n
(τ2,τ3,...,τi,...,τn)]k(k＝1,
…
,m)全部叠加可以实现绕射波信号的同相位叠加，而非相关噪声却无法实现同相位叠加，因此被相对压制。最后得到叠加后的累积量函数c
1,2,...,i,...n
(τ2,τ3,...,τi,...,τn)，包含了虚绕射记录(图2c)，不再与炮点相关，因此下标k被删除。
[0073]
此步骤涉及累积量阶数的选取，使用待处理数据，选取累积量阶数为三阶、四阶和五阶的累积量函数计算虚绕射记录，综合评价结果的信噪比和时延估计的准确性，确定最优的累积量阶数。大多数情况下，四阶累积量可达到理想处理效果。
[0074]
c2.通过多道褶积将虚绕射记录恢复成炮检点间的绕射波记录：
[0075]
由于虚绕射记录仅包含了检波点间的绕射波信号，不具有真实的炮点位置，但在常规地震数据处理中经常需要使用炮检点间的地震记录，因此有必要进一步将虚绕射记录恢复为炮检点间的绕射波记录，实现的物理过程如下(图3)。将累积量函数c
j,1,2,...,i,...,n
(τ1,τ2,...,τi,,...,τn)与原始地震记录u
s；j
(t)进行多道褶积，可重构炮点s到检波点ri(i＝1,2,
…
,n)的绕射波记录，公式如下：
[0076][0077]
其中，中的符号“^”表示该绕射波是重构得到的，区别于原始地震记录中的绕射波，下标j表示该重构信号是使用检波点rj的地震记录得到。多道褶积是将待褶积地震数据进行傅里叶变换转换至频率域函数，在频率域乘积后，进行傅里叶反变换，得到时间域的多道褶积结果。由于c
j，1,2，...,i,...,n
(t
1-τ1,t
2-τ2,...,t
i-τi,...,t
n-τn)中包含了n组检波点间的时延，需要将这些时延全部重构成炮检点间的旅行时，因此需要使用n道原始地震记录u
s；j
(t)对每个时延进行褶积，即多道褶积中“多道”的含义。累积量函数中的绕射波与u
s；j
(t)共用d到rj的射线路径，但具有符号相反的旅行时和由于褶积对旅行时是相加的过程，因此，这对旅行时在多道褶积的过程中被抵消，由此重构的地震记录中旅行时ti(i＝1,2,
…
,n)为：
[0078][0079]
即重构了炮点s在检波点rj接收的绕射波记录(图3b)，该重构记录与检波点rj不再相关，因此使用不同检波点rj,rk,r
l
等得到的重构记录中包含了相同的绕射波，可以实现同相位叠加，而非相关噪声依然无法实现同相位叠加，再次被压制，最终得到叠加后的重构绕射波记录(图3c)。
[0080]
c3.除了绕射波，地震数据初至波还包含近偏移距的直达波和单点绕射波，均可使用ci方法增强：
[0081]
考虑炮检点均位于地表的地震观测系统，与绕射波的情况相同，同一炮点激发在不同检波点接收的远偏移距折射波共用了炮点到高速层顶的旅行时，共用旅行时在计算高阶累积量的过程中被抵消，得到了与炮点无关的虚折射记录，使用不同炮点得到相同的虚折射记录，可以实现同相位叠加增强信噪比。对于近偏移距的直达波是完全相同的情况。地震记录包含浅层反射波和背景噪声，会影响初至波的识别。对于反射波，相同炮点在不同检波点的反射波不共用射线路径，高阶累积量结果中的虚反射记录对应的时延会随着炮点位置的变化而变化，无法实现同相位叠加，因此相比初至波信号被相对压制。对于不相关的背景噪声也是相同的情况。因此，绕射波、折射波和直达波在内的初至波信号可以通过ci方法增强，而反射波、背景噪声因无法实现同相位叠加被相对压制。
[0082]
d.结果检验模块：检查输出道集是否达到了满意的初至波信噪比提升效果，如果未达到，则以步骤c2得到的重构绕射波记录为新的输入道集，即步骤c1中公式(2)所用的地震记录，重复步骤c再次加强初至波，直到实现满意的信噪比提升效果。
[0083]
判断输出道集是否达到满意的初至波信噪比提升效果，是通过基于能量比的长短窗平均(short term averaging/long term averaging,sta/lta)方法自动拾取初至；如果初至拾取精度达到后续的层析或折射静校正所要求的精度，则视为达到满意的初至波信噪比提升效果，否则迭代执行步骤c直到实现满意的信噪比提升效果。
[0084]
下面通过实例对本发明做进一步说明。
[0085]
实施例一：
[0086]
本发明通过有限差分方法求解三维声波方程得到三维合成数据，所用三维近地表模型包含上下两层，声波速度分别为1500m/s和2300m/s，地表为起伏地表，上下地层的分解面为存在一处断层的高速层顶结构，该近地表模型模拟中国西部地区的薄黄土或沙漠地区的结构。使用主频为40hz的雷克子波作为震源子波，在远偏移距，初至波路径可能包含沿速度分界面传播的折射波，经断层产生的绕射波或两者的组合。
[0087]
首先验证ci方法相比传统的svi方法具有更好的抗噪性。抽取远偏移距的原始二维数据，并用时窗进行初至波截取(图4a)。在原始数据中加入了非相干高斯噪声和视速度为
–
500和1500m/s的相干高斯噪声，含噪数据的初至波信噪比为
–
12db，初至同相轴几乎全部被噪声掩盖(图4b)。svi方法可以一定程度上压制噪声，初至波信噪比提高至
–
2db，但由于svi对相干噪声敏感，结果中仍有明显的相干噪声残留(图4c)。使用ci方法可以很好地压制相干噪声，初至同相轴变得清晰可见，信噪比提高到11db(图4d)。同时注意到，结果中的反射波被压制，一是由于反射波在ci方法中无法实现同相位叠加，相比初至波被相对压制，二是由于时窗截取数据覆盖了所有偏移距的初至波，但并不保证覆盖全部的反射波，因此，在ci方法的叠加步骤中，初至波的能力覆盖次数会高于反射波。
[0088]
实施例二：
[0089]
选取真实勘探工区的实际地震资料检验ci方法的应用效果，工区呈现起伏地表，部分原始数据道集存在强背景噪声，难以准确自动拾取初至，因此需要加强初至波的信噪比，保证初至拾取的精度，为后续的近地表速度建模提供可靠的输入数据。
[0090]
原始地震记录中的初至波同相轴被4至30hz的强相干噪声所掩盖(图5a)。商业软件中常用的去噪方法是f-x域滤波，在频率域使用梯形窗截取制定频带内的数据，经参数测试选取的最佳低频截断频率为7和12hz，高频截断频率为55和65hz。f-x滤波结果中(图5b)，
低频噪声被相对压制，但仍有较多噪声残留，这部分噪声与有效信号的频带范围接近，难以通过f-x域滤波分离。ci方法显著提高了噪声压制的效果(图5c)，初至波同相轴清晰可见，可以很好地保证初至自动拾取的准确性。取得如此好的初至增强效果，一是由于ci方法不依赖频率信息，可以更好地将频率接近的噪声同有效信号分离，二要归功于相邻区域内道集的初至波信噪比更优，ci方法实现了相邻区域内的道集初至波同相位叠加，高信噪比道集叠加在低信噪比道集中，将有助于初至信噪比的提升。
[0091]
需要注意的是，公布实施例的目的在于帮助进一步理解本发明，但是本领域的技术人员可以理解：在不脱离本发明及所附权利要求的精神和范围内，各种替换和修改都是可能的。因此，本发明不应局限于实施例所公开的内容，本发明要求保护的范围以权利要求书界定的范围为准。