CN111257938A - 基于小波互相关时移地震虚拟震源波场重构方法和系统 - Google Patents

基于小波互相关时移地震虚拟震源波场重构方法和系统 Download PDF

Info

Publication number
CN111257938A
CN111257938A CN202010216844.5A CN202010216844A CN111257938A CN 111257938 A CN111257938 A CN 111257938A CN 202010216844 A CN202010216844 A CN 202010216844A CN 111257938 A CN111257938 A CN 111257938A
Authority
CN
China
Prior art keywords
wavelet
time
cross
correlation
seismic
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
CN202010216844.5A
Other languages
English (en)
Inventor
周阳田丽
赵杨
左佳卉
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
China University of Petroleum Beijing
Original Assignee
China University of Petroleum Beijing
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by China University of Petroleum Beijing filed Critical China University of Petroleum Beijing
Priority to CN202010216844.5A priority Critical patent/CN111257938A/zh
Publication of CN111257938A publication Critical patent/CN111257938A/zh
Pending legal-status Critical Current

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/36—Effecting static or dynamic corrections on records, e.g. correcting spread; Correlating seismic signals; Eliminating effects of unwanted energy
    • G01V1/362—Effecting static or dynamic corrections; Stacking
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/36—Effecting static or dynamic corrections on records, e.g. correcting spread; Correlating seismic signals; Eliminating effects of unwanted energy
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/36—Effecting static or dynamic corrections on records, e.g. correcting spread; Correlating seismic signals; Eliminating effects of unwanted energy
    • G01V1/364—Seismic filtering
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/30—Noise handling
    • G01V2210/32—Noise reduction
    • G01V2210/324—Filtering
    • G01V2210/3246—Coherent noise, e.g. spatially coherent or predictable
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/50—Corrections or adjustments related to wave propagation
    • G01V2210/51—Migration
    • G01V2210/512—Pre-stack

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Acoustics & Sound (AREA)
  • Environmental & Geological Engineering (AREA)
  • Geology (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • General Physics & Mathematics (AREA)
  • Geophysics (AREA)
  • Geophysics And Detection Of Objects (AREA)

Abstract

本发明及一种基于小波互相关时移地震虚拟震源波场重构方法和系统，包括以下步骤：S1.将垂直地震剖面VSP地震数据进行上行波和下行波分离；S2.将上行波和下行波都进行小波变换，选择一检波器作为虚拟震源，把检波器的下行波小波系数，与其它检波器的上行波小波系数互相关，得到小波互相关系数；S3.对小波互相关系数进行时间尺度域、时间尺度波数域或时间空间域滤波；S4.将经过滤波的小波互相关系数逆变换回时空域，得到虚拟震源道集；S5.将虚拟震源道集进行共中心点叠加，得到虚拟震源激发波场。其能够有效压制噪声，并对非平稳信号实现高分辨率分离，无需近地表模型即可有效地抑制伪像，优化叠加成像的结果。

Description

基于小波互相关时移地震虚拟震源波场重构方法和系统
技术领域
本发明是关于一种基于小波互相关时移地震虚拟震源波场重构方法和系统，属于地震数据波处理技术领域。
背景技术
时移地震技术是利用不同时间观测的地震有效信息的差异，进行储层监测、完善油藏管理方案、提高油气采收率的技术。时移地震的可重复性是长期油藏监测最重要的影响因素。时移地震的可重复性与采集中的仪器因素，环境、背景噪声等采集因素密切相关。其中仪器因素较易控制，而环境、背景噪声往往很难控制，是造成监测结果不可重复的主要原因。环境、背景噪声主要源于由震源、检波器和近地表条件的变化导致的时间延迟和相位、振幅的畸变。现有技术中，通常通过将震源和检波器埋入地下的方式来提升时移地震的可重复性。但对于沙漠环境，近地表发育溶岩结构的地区的成像结果尤为复杂，会引入多次波、面波、散射噪声，以及不同波型互相干扰，使上下行波场发生畸变。
为了解决上述问题，引入虚拟震源法(VS)对时移地震进行检测，虚拟震源法不依赖于近地表速度模型，其将上下行波场进行互相关重构，将地面震源重构到地下观测系统的位置，压制来自近地表仪器耦合和观测系统变化形成的噪声，在储层监测和时移地震处理中，提高了可重复性。然而，由于有限的震源孔径和野外预处理的结果不能完全满足虚拟震源的假设条件，因此仍然存在重复性差和图像质量差的问题。尤其，在强烈非均质的上覆层中，下行波场往往被近地表反射和自由表面多次波的上行波干扰，使二者难以分离。一旦下行的P波受到污染，就无法有效抑制包含多次波、散射和面波的噪声和串扰，导致传统虚拟震源叠加的信噪比低于实际震源叠加，使成像质量下降。
发明内容
针对上述现有技术的不足，本发明的目的是提供了一种基于小波互相关时移地震虚拟震源波场重构方法和系统，在传统虚拟震源方法处理流程中，加入小波域虚拟震源波场重构，对小波系数进行时间尺度域TS或时间尺度波数域TSK滤波，从而有效压制噪声，并对非平稳信号实现高分辨率分离。该方法无需近地表模型即可有效地抑制伪像，优化叠加成像的结果。
为实现上述目的，本发明提供了一种基于小波互相关时移地震虚拟震源波场重构方法，包括以下步骤：S1.通过若干检波器获得垂直地震剖面VSP地震数据，将VSP地震数据进行上行波和下行波分离；S2.将分离后的上行波和下行波都进行小波变换，选择其中一个检波器作为虚拟震源，把检波器的下行波小波系数与其它检波器的上行波小波系数互相关，得到小波互相关系数；S3.对小波互相关系数进行时间尺度域TS或时间尺度波数域TSK滤波；S4.将经过滤波的小波互相关系数逆变换回时空域，得到虚拟震源道集；S5.将虚拟震源道集进行共中心点叠加，得到虚拟震源激发波场。
进一步，步骤S1中垂直地震剖面VSP地震数据采用下式表示：
Figure BDA0002424697860000021
其中，
Figure BDA0002424697860000022
表示时间互相关，rA、rB和rS分别是检波器A、检波器B和震源的空间坐标，V(rB|rA；τ)是将rA视为虚拟震源时在检波器rB处的干涉结果，τ表示时间延迟，D(rB|rS；t)代表检波器rA在近偏移距位置记录的下行波场，U(rB|rS；t)是检波器rB接收到的上行波场。
进一步，步骤S1中垂直地震剖面VSP地震数据分别通过去噪和加窗获得上行波和下行波，上行波和下行波的公式为：
U(rB|rS；t)＝Up(rB|rS；t)+UM(rB|rS；t)+USh(rB|rS；t)...
D(rB|rS；t)＝Dp(rB|rS；t)+DM(rB|rS；t)+DSh(rB|rS；t)...
其中，D(rB|rS；t)是检波器rB接收到的下行波场，Up(rB|rS；t)、Dp(rB|rS；t)分别是直达P波的上行波场和下行波场，UM(rB|rS；t)、DM(rB|rS；t)分别为多次波近地表反射的上行波场和下行波场，USh(rB|rS；t)、DSh(rB|rS；t)分别为S波的上行波场和下行波场。
进一步，步骤S2中进行小波变换为：对于一个给定的地震信号D(t)采用以下公式进行小波变换：
Figure BDA0002424697860000023
Figure BDA0002424697860000024
是小波系数，φ*(t)表示满足适应条件且均值为零的母子波，α为小波尺度，τ是局部时间。
进一步，步骤S2中小波互相关系数采用计算公式为：
Figure BDA0002424697860000025
其中，WV(α，τ)为小波互相关系数，
Figure BDA0002424697860000026
和
Figure BDA0002424697860000027
分别是上行波和下行波的小波系数，
Figure BDA0002424697860000028
是下行波的小波系数的共轭函数，T是信号周期。
进一步，小波互相关系数满足如下条件：
WV(α，ω)＝αV(ω)|φ(αω)|2
|WV(α，ω)|＝α|V(ω)||φ(αω)|2
∠WV(α，ω)＝∠V(ω)
WV(α，ω)是小波尺度和频域中的小波互相关系数，V(ω)是信号互相关谱，φ(αω)是母小波的频谱，∠表示相位谱，ω是角频率，α为小波尺度。
进一步，步骤S3中时间尺度空间域的三维滤波方程为：
WVo(α，τ，x)＝H(α，τ，x)*α，τWVi(α，τ，x)
其中，WV0(α，τ，x)是时间尺度空间域中小波互相关经过滤波的小波互相关系数；H(α，τ，x)是时间尺度空间域中的滤波函数；α，τWVi(α，τ，x)是时间尺度空间域中的小波互相关原始信号。
进一步，对步骤S3中时间尺度空间域的三维滤波方程进行二维傅立叶变换，得到时间频率波数域，所述时间频率波数域的三维滤波方程为：
Figure BDA0002424697860000031
其中，WV0(α，ω，k)是滤波结果在时间频率波数域中小波互相关经过滤波的小波互相关系数，VV(ω，k)是频率波数域中虚拟震源道集的互相关谱，p(α，ω，k)是依赖于尺度因子α的滤波函数；
Figure BDA0002424697860000035
是滤波函数中只与ω，k相关的部分；
Figure BDA0002424697860000032
是滤波函数中只与尺度因子α相关的部分；Hkα(α，k)是滤波函数中只与α，k相关的部分；φ(αω)是母小波的频谱。
进一步，小波互相关系数的逆变换公式为：
Figure BDA0002424697860000033
其中，D(t)为经过处理后的地震信号，
Figure BDA0002424697860000034
是小波系数，φ(t)表示满足适应条件且均值为零的母子波，α为小波尺度，τ是局部时间；Cφ是比例系数。
基于相同的发明构思，本发明还公开了一种基于小波互相关时移地震虚拟震源波场重构系统，包括：上下行波分离模块，用于将垂直地震剖面VSP地震数据进行上行波和下行波分离；小波互相关系数获取模块，用于将上行波和下行波都进行小波变换，选择一检波器作为虚拟震源，把检波器的下行波小波系数，与其它检波器的上行波小波系数互相关，得到小波互相关系数；滤波模块，用于对小波互相关系数进行时间尺度域、时间尺度波数域或时间空间域滤波；逆变换模块，用于将经过滤波的小波互相关系数逆变换回时空域，得到虚拟震源道集；共中心点叠加模块，用于将虚拟震源道集进行共中心点叠加，得到虚拟震源激发波场。
本发明由于采取以上技术方案，其具有以下优点：
1、在传统虚拟震源方法处理流程中，加入小波域虚拟震源波场重构，对小波系数进行时间尺度域TS或时间尺度波数域TSK滤波，从而有效压制噪声，并对非平稳信号实现高分辨率分离，该方法无需近地表模型即可有效地抑制伪像，优化叠加成像的结果。
2、本发明将地震干涉方法从一维时间域拓展到了多维TS域和TSK域，提高信噪比，使依赖于频率和波数的噪声和信号能够更好地分离。
3、相比于常规虚拟震源方法，本发明能够有效抑制与近地表结构相关的噪声，大幅提升时移地震的可重复性。
附图说明
图1是本发明一实施例中基于小波互相关时移地震虚拟震源波场重构方法的流程图；
图2是经过共道集噪声处理的单炮信号图，其中，上方的虚线上方为P波的下行波场，下方的虚线下方是P波的上行波场；
图3中图3(a)是某一虚拟单炮信号的下行波场；图3(b)是该下行波场对应的频率-波数谱；
图4中图4(a)是经过TS滤波后的FK域的下行波场；图4(b)是过TS滤波后的下行波场对应的频率-波数谱，图4中箭头标出的部分与图3中相应位置的图像相比，图4(a)的图像清晰度提高，图4(b)的波形噪音较低；
图5中图5(a)是仅经过TS滤波的虚拟单炮信号图；图5(b)是经过TS滤波和TSX三维滤波的虚拟单炮信号图；
图6是对TSX三维滤波进行二维傅立叶变换后得到的TSK三维数据体；
图7是进行SK滤波后得到的三维数据体；
图8中图8(a)是采用传统虚拟震源法经过共中心点叠加处理的波场剖面图，图8(b)是采用本发明一实施例中方法经过共中心点叠加处理的波场剖面，将图8(a)和图8(b)中箭头标出的部分进行比较，可见图8(b)中图像更加清晰，噪音较低；
图9是13次时移地震的二维剖面形成的三维数据体，图9(a)非虚拟震源法得到的三围数据体；图9(b)现有技术中虚拟震源法得到的三围数据体；图9(c)本实施例中虚拟震源法得到的三围数据体。
具体实施方式
为了使本领域技术人员更好的理解本发明的技术方向，通过具体实施例对本发明进行详细的描绘。然而应当理解，具体实施方式的提供仅为了更好地理解本发明，它们不应该理解成对本发明的限制。在本发明的描述中，需要理解的是，所用到的术语仅仅是用于描述的目的，而不能理解为指示或暗示相对重要性。
实施例一
本实施例公开了一种基于小波互相关时移地震虚拟震源波场重构方法，如图1所示，包括以下步骤：
S1.通过若干检波器获得垂直地震剖面VSP地震数据，并将垂直地震剖面VSP地震数据进行上行波和下行波分离。
其中，垂直地震剖面VSP地震数据采用下式表示：
Figure BDA0002424697860000051
其中，
Figure BDA0002424697860000052
表示时间互相关，rA、rB和rS分别是检波器A、检波器B和震源的空间坐标，V(rB|rA；τ)是将rA视为虚拟震源时在检波器rB处的干涉结果，τ表示时间延迟，D(rB|rS；t)代表检波器rA在近偏移距位置记录的下行波场，U(rB|rS；t)是检波器rB接收到的上行波场。
垂直地震剖面VSP地震数据分别通过去噪和加窗获得上行波和下行波，上行波和下行波的公式为：
U(rB|rS；t)＝Up(rB|rS；t)+UM(rB|rS；t)+USh(rB|rS；t)...
D(rB|rS；t)＝Dp(rB|rS；t)+DM(rB|rS；t)+DSh(rB|rS；t)...
其中，D(rB|rS；t)是检波器rB接收到的下行波场，Up(rB|rS；t)、Dp(rB|rS；t)分别是直达P波的上行波场和下行波场，UM(rB|rS；t)、DM(rB|rS；t)分别为多次波近地表反射的上行波场和下行波场，USh(rB|rS；t)、DSh(rB|rS；t)分别为S波的上行波场和下行波场。
S2.将上行波和下行波都进行小波变换，选择一检波器作为虚拟震源，把检波器的下行波小波系数，与其它检波器的上行波小波系数互相关，得到小波互相关系数。
其中，进行小波变换步骤具体为：对于一个给定的地震信号D(t)采用以下公式进行小波变换：
Figure BDA0002424697860000053
其中，
Figure BDA0002424697860000054
是小波系数，φ*(t)表示满足适应条件且均值为零的母子波，α为小波尺度，τ是局部时间。
通过上行波和下行波的小波系数互相关计算小波互相关系数，小波互相关系数采用如下公式计算：
Figure BDA0002424697860000055
其中，WV(α，τ)为小波互相关系数，
Figure BDA0002424697860000061
和
Figure BDA0002424697860000062
分别是上行波和下行波的小波系数，
Figure BDA0002424697860000063
是下行波的小波系数的共轭函数，T是信号周期。
为了保证小波尺度和频域中的小波的同相，小波互相关系数满足如下条件：
WV(α，ω)＝αV(ω)|φ(αω)|2
|WV(α，ω)|＝α|V(ω)||φ(αω)|2
∠WV(α，ω)＝∠V(ω)
WV(α，ω)是小波尺度和频域中的小波互相关系数，V(ω)是信号互相关谱，φ(αω)是母小波的频谱，∠表示相位谱，ω是角频率，α为小波尺度。该方法仅对振幅谱滤波，保留原始相位，并从虚拟震源中提取出格林函数的运动学特征，是一种保相技术。
S3.对小波互相关系数进行时间尺度域TS或时间尺度波数域TSK滤波。
其中，时间尺度空间域的三维滤波方程为：
WVo(α，τ，x)＝H(α，τ，x)*α，τWVi(α，τ，x)
其中，WV0(α，τ，x)是时间尺度空间域中小波互相关经过滤波的小波互相关系数；H(α，τ，x)是时间尺度空间域中的滤波函数；α，τWVi(α，τ，x)是时间尺度空间域中的小波互相关原始信号。
其中，滤波函数H(α，τ，x)旨在改善同相轴识别、处理非平稳特性并压制虚拟震源法处理过程中产生的相干噪声。三维滤波函数H(α，τ，x)表示为：
Figure BDA0002424697860000064
Figure BDA0002424697860000065
滤波函数中只与α，k相关的部分，可以是简单的α，k带通滤波器，允许特定波速的信号通过。
对时间尺度空间域的三维滤波方程进行二维傅立叶变换，得到时间频率波数域，所述时间频率波数域的三维滤波方程为：
Figure BDA0002424697860000066
其中，WV0(α，ω，k)是滤波结果在时间频率波数域中小波互相关经过滤波的小波互相关系数，V(ω，k)是频率波数域中虚拟震源道集的互相关谱，p(α，ω，k)是依赖于尺度因子α的滤波函数；
Figure BDA0002424697860000067
是滤波函数中只与ω，k相关的部分；
Figure BDA0002424697860000068
是滤波函数中只与尺度因子α相关的部分；Hkα(α，k)是滤波函数中只与α，k相关的部分；φ(αω)是母小波的频谱。
S4.将经过滤波的小波互相关系数逆变换回时空域，得到虚拟震源道集。
其中，小波互相关系数的逆变换公式为：
Figure BDA0002424697860000071
其中，D(t)为经过处理后的地震信号，
Figure BDA0002424697860000072
是小波系数，φ(t)表示满足适应条件且均值为零的母子波，α为小波尺度，τ是局部时间，Cφ是比例系数。
S5.虚拟震源道集进行共中心点叠加，得到虚拟震源激发波场。
实施例二
本实施例以沙漠环境中的时移地震实际数据对实施例一中的方法进行了测试。具体过程如下：
在19个月内进行了13次重复二维勘探，前6次(S1-S6)在3个月完成，中断了17个月后又在一周内完成了后7次勘探(S7-S13)。13次勘探均使用了Mertz 26可控震源，震源重复激发的精度达1m。平均间隔30m的80个检波器和水听器固定在最深达50m的地下，组成一条二维测线。为了有效去除线性和散射噪声，实现最佳照明，进行了密集的三维面激发，将炮距设置为7.5×7.5m。偏移距范围为0到2400m。目标层位于地下约2000米处，上覆有多套速度差异较大的地层。近地表为厚度从几米到几十米不等的砂层，在近地表以下是简单的层状结构，反射界面的倾角均小于5°。
地震数据预处理时，先自动拾取初至的波形，再通过初至的波形的位置截取直达波的时窗用于下行波场近似。图2表示经过共道集噪声处理的单炮信号，上方虚线上面的数据无需处理，可直接用于p波的下行波场近似，对应30m或100ms以内的小偏移距。下方虚线下面的数据反射到达较晚，大于300ms，对应大偏移距，用于p波的上行波场近似。由于该区域地面情况复杂，下行波场往往包含来自近地表的干扰。
如图3所示，图3(a)是某一虚拟单炮信号的频率-波数域即FK域的下行波场；图3(b)是该下行波场对应的频率-波数谱。图3(a)中浅色区域为散射噪声，被浅色区域包裹的深色区域表示来自目标层的强反射。在FK域中的下行波场中的能量明显失真，因此，需进行适当滤波。与常规虚拟震源法相比，小波互相关和多通道组合滤波能将多维空间中的反射波信号从背景散射噪声中分离，提高信噪比，增加了反射波信号。本实施例中对每个炮集的上行波和下行波进行小波互相关。进行小波互相关时，由于在有效带宽内，目标储层的反射信号出现在约1.2-1.4s处。首先通过匹配目标储层的反射信号设置参数，本实施例中有效带宽优选为20-50Hz，参数优选α＝20Hz，τ＝0.2s和η＝5％的滑动二维窗口内的最大值。然后，将滤波函数与图3a中所示的小波系数振幅谱进行卷积。
经过小波互相关滤波的结果如图4所示，图4(a)是经过TS滤波后的FK域的下行波场；图4(b)是过TS滤波后的下行波场对应的频率-波数谱。在互相关计算之前，去除散射噪声能够进一步提升整体的信噪比。从图4(a)和图4(b)中可以看出经过滤波后，背景噪声明显减弱，但是浅层(＜0.6s)的高频成分(＞60Hz)的信号被衰减，从而降低分辨率。
随后，采用多分量滤波器处理被面波污染的虚拟单炮信号。将经过TS滤波生成的虚拟单炮信号变换到时间尺度空间域上，进行TSX三维滤波，TSX三维滤波的结果如图5所示。图5(a)是仅经过TS滤波的虚拟单炮信号图；图5(b)是经过TS滤波和TSX三维滤波的虚拟单炮信号图。从图5(a)和图5(b)的对比中可以看出，经过TSX三维滤波的虚拟单炮信号强度明显提高，细节特征也更加明显，信噪比得到有效提高。
图6为对TSX三维滤波进行二维傅立叶变换后得到的TSK三维数据体。通过在TSK三维数据体中调整小波互相关系数可以获得特定的时间和频率分辨率。TSK滤波器包含时间信息，因而能够处理小波的非平稳性。当滤波器时间不变时，此方法简化为SK滤波，如图7所示。由此可见，即使使用简单的TSK滤波器，小波互相关也能够消除残留的面波，这是传统虚拟震源法所不能实现的。
图8中，图8(a)是采用传统虚拟震源法经过共中心点叠加处理的波场剖面图，图8(b)是采用本实施例中方法经过共中心点叠加处理的波场剖面。从图8(a)和图8(b)中可以看出，图8(b)的波场信号的信噪比和同相轴连续性的得到了显著提升。
图9是13次时移地震的二维剖面形成的三维数据体，图9(a)非虚拟震源法得到的三围数据体；图9(b)现有技术中虚拟震源法得到的三围数据体；图9(c)本实施例中虚拟震源法得到的三围数据体。从图9(a)、图9(b)和图9(c)中可以看出，经过本实施例中方法处理的时移地震波场具有更好的重复性。
本实施例中方法的优点在于能够有效消除强振幅噪声，特别是大偏移距数据中常出现的面波和陆上地震数据中的近地表溶岩产生的散射。不同的信号分解方法相互组合可以有效去除噪声。例如，小倾角的简单层状结构横向变化较小，TS滤波加TSX阈值是理想的选择。而对于包含断层或盐丘的地质结构，各向异性较强，存在潜在散射体，则需将信号转换到自适应TSK域或TK进行滤波。通用三维滤波比二维滤波更贵，但能更有效的解决高维度地震数据中的非平稳性。
实施例三
基于相同的发明构思，本实施例公开了一种基于小波互相关时移地震虚拟震源波场重构系统，包括：
上下行波分离模块，用于将垂直地震剖面VSP地震数据进行上行波和下行波分离；
小波互相关系数获取模块，用于将上行波和下行波都进行小波变换，选择一检波器作为虚拟震源，把检波器的下行波小波系数，与其它检波器的上行波小波系数互相关，得到小波互相关系数；
滤波模块，用于对小波互相关系数进行时间尺度域、时间尺度波数域或时间空间域滤波；
逆变换模块，用于将经过滤波的小波互相关系数逆变换回时空域，得到虚拟震源道集；共中心点叠加模块，用于将虚拟震源道集进行共中心点叠加，得到虚拟震源激发波场。
上述内容仅为本申请的具体实施方式，但本申请的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本申请揭露的技术范围内，可轻易想到变化或替换，都应涵盖在本申请的保护范围之内。因此，本申请的保护范围应以权利要求的保护范围为准。

Claims (10)

1.一种基于小波互相关时移地震虚拟震源波场重构方法，其特征在于，包括以下步骤：
S1.通过若干检波器获得垂直地震剖面VSP地震数据，将所述VSP地震数据进行上行波和下行波分离；
S2.将分离后的上行波和下行波都进行小波变换，选择其中一个检波器作为虚拟震源，把所述检波器的下行波小波系数与其它检波器的上行波小波系数互相关，得到小波互相关系数；
S3.对小波互相关系数进行时间尺度域TS或时间尺度波数域TSK滤波；
S4.将经过滤波的小波互相关系数逆变换回时空域，得到虚拟震源道集；
S5.将所述虚拟震源道集进行共中心点叠加，得到虚拟震源激发波场。
2.如权利要求1所述的基于小波互相关时移地震虚拟震源波场重构方法，其特征在于，所述步骤S1中垂直地震剖面VSP地震数据采用下式表示：
Figure FDA0002424697850000012
其中，
Figure FDA0002424697850000013
表示时间互相关，rA、rB和rS分别是检波器A、检波器B和震源的空间坐标，V(rB|rA；τ)是将rA视为虚拟震源时在检波器rB处的干涉结果，τ表示时间延迟，D(rB|rS；t)代表检波器rA在近偏移距位置记录的下行波场，U(rB|rS；t)是检波器rB接收到的上行波场。
3.如权利要求2所述的基于小波互相关时移地震虚拟震源波场重构方法，其特征在于，所述步骤S1中垂直地震剖面VSP地震数据分别通过去噪和加窗获得上行波和下行波，上行波和下行波的公式为：
U(rB|rS；t)＝Up(rB|rS；t)+UM(rB|rS；t)+USh(rB|rS；t)...
D(rB|rS；t)＝Dp(rB|rS；t)+DM(rB|rS；t)+DSh(rB|rS；t)...
其中，D(rB|rS；t)是检波器rB接收到的下行波场，Up(rB|rS；t)、Dp(rB|rS；t)分别是直达P波的上行波场和下行波场，UM(rB|rS；t)、DM(rB|rS；t)分别为多次波近地表反射的上行波场和下行波场，USh(rB|rS；t)、DSh(rB|rS；t)分别为S波的上行波场和下行波场。
4.如权利要求1-3任一项所述的基于小波互相关时移地震虚拟震源波场重构方法，其特征在于，所述步骤S2中进行小波变换为：对于一个给定的地震信号D(t)采用以下公式进行小波变换：
Figure FDA0002424697850000011
其中，
Figure FDA0002424697850000021
是小波系数，φ*(t)表示满足适应条件且均值为零的母子波，α为小波尺度，τ是局部时间。
5.如权利要求4所述的基于小波互相关时移地震虚拟震源波场重构方法，其特征在于，步骤S2中小波互相关系数计算公式为：
Figure FDA0002424697850000022
其中，WV(α，τ)为小波互相关系数，
Figure FDA0002424697850000023
和
Figure FDA0002424697850000024
分别是上行波和下行波的小波系数，
Figure FDA0002424697850000025
是下行波的小波系数的共轭函数，T是信号周期。
6.如权利要求5所述的基于小波互相关时移地震虚拟震源波场重构方法，其特征在于，所述小波互相关系数满足如下条件：
WV(α，ω)＝αV(ω)|φ(αω)|2
|WV(α，ω)|＝α|V(ω)||φ(αω)|2
∠WV(α，ω)＝∠V(ω)
WV(α，ω)是小波尺度和频域中的小波互相关系数，V(ω)是信号互相关谱，φ(αω)是母小波的频谱，∠表示相位谱，ω是角频率，α为小波尺度。
7.如权利要求5或6所述的基于小波互相关时移地震虚拟震源波场重构方法，其特征在于，所述步骤S3中时间尺度空间域滤波采用三维滤波方程：
WVo(α，τ，x)＝H(α，τ，x)*α，τWVi(α，τ，x)
其中，WV0(α，τ，x)是时间尺度空间域中小波互相关经过滤波的小波互相关系数；H(α，τ，x)是时间尺度空间域中的滤波函数；α，τWVi(α，τ，x)是时间尺度空间域中的小波互相关原始信号。
8.如权利要求7所述的基于小波互相关时移地震虚拟震源波场重构方法，其特征在于，对步骤S3中时间尺度空间域的三维滤波方程进行二维傅立叶变换，得到时间频率波数域，所述时间频率波数域的三维滤波方程为：
Figure FDA0002424697850000026
其中，WV0(α，ω，k)是滤波结果在时间频率波数域中小波互相关经过滤波的小波互相关系数，VV(ω，k)是频率波数域中虚拟震源道集的互相关谱，p(α，ω，k)是依赖于尺度因子α的滤波函数；
Figure FDA0002424697850000027
是滤波函数中只与ω，k相关的部分；
Figure FDA0002424697850000028
是滤波函数中只与尺度因子α相关的部分；Hkα(α，k)是滤波函数中只与α，k相关的部分；φ(αω)是母小波的频谱。
9.如权利要求7所述的基于小波互相关时移地震虚拟震源波场重构方法，其特征在于，所述步骤S4中小波互相关系数的逆变换公式为：
Figure FDA0002424697850000031
其中，D(t)为经过处理后的地震信号，
Figure FDA0002424697850000032
是小波系数，φ(t)表示满足适应条件且均值为零的母子波，α为小波尺度，τ是局部时间，Cφ是比例系数。
10.一种基于小波互相关时移地震虚拟震源波场重构系统，其特征在于，包括：
上下行波分离模块，用于将垂直地震剖面VSP地震数据进行上行波和下行波分离；
小波互相关系数获取模块，用于将所述上行波和下行波都进行小波变换，选择一检波器作为虚拟震源，把所述检波器的下行波小波系数，与其它检波器的上行波小波系数互相关，得到小波互相关系数；
滤波模块，用于对所述小波互相关系数进行时间尺度域、时间尺度波数域或时间空间域滤波；
逆变换模块，用于将经过滤波的所述小波互相关系数逆变换回时空域，得到虚拟震源道集；
共中心点叠加模块，用于将所述虚拟震源道集进行共中心点叠加，得到虚拟震源激发波场。
CN202010216844.5A 2020-03-25 2020-03-25 基于小波互相关时移地震虚拟震源波场重构方法和系统 Pending CN111257938A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010216844.5A CN111257938A (zh) 2020-03-25 2020-03-25 基于小波互相关时移地震虚拟震源波场重构方法和系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010216844.5A CN111257938A (zh) 2020-03-25 2020-03-25 基于小波互相关时移地震虚拟震源波场重构方法和系统

Publications (1)

Publication Number Publication Date
CN111257938A true CN111257938A (zh) 2020-06-09

Family

ID=70946179

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010216844.5A Pending CN111257938A (zh) 2020-03-25 2020-03-25 基于小波互相关时移地震虚拟震源波场重构方法和系统

Country Status (1)

Country Link
CN (1) CN111257938A (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113740906A (zh) * 2021-10-19 2021-12-03 中国地质大学（北京） 一种水下垂直缆地震波干涉成像方法及装置
CN113933894A (zh) * 2020-07-13 2022-01-14 中国石油天然气股份有限公司 虚拟源成像的方法及装置

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107402405A (zh) * 2016-05-18 2017-11-28 中国石油化工股份有限公司 静相位虚源道集构建方法
US20180217284A1 (en) * 2017-01-27 2018-08-02 Saudi Arabian Oil Company Virtual source redatuming using radiation pattern correction
CN110178056A (zh) * 2016-11-17 2019-08-27 沙特阿拉伯石油公司 利用小波互相关进行虚拟源去噪
CN110907989A (zh) * 2018-09-17 2020-03-24 中国石油化工股份有限公司 重建拟地面地震反射波成像方法及系统

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107402405A (zh) * 2016-05-18 2017-11-28 中国石油化工股份有限公司 静相位虚源道集构建方法
CN110178056A (zh) * 2016-11-17 2019-08-27 沙特阿拉伯石油公司 利用小波互相关进行虚拟源去噪
US20180217284A1 (en) * 2017-01-27 2018-08-02 Saudi Arabian Oil Company Virtual source redatuming using radiation pattern correction
CN110907989A (zh) * 2018-09-17 2020-03-24 中国石油化工股份有限公司 重建拟地面地震反射波成像方法及系统

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
YANG ZHAO 等: "WAVELET-CROSSCORRELATION-BASED INTERFEROMETRIC REDATUMING IN 4D SEISMIC", 《GEOPHYSICS》 *
曾靖雯 等: "基于小波变换的虚震源信号去噪研究", 《地球物理学进展》 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113933894A (zh) * 2020-07-13 2022-01-14 中国石油天然气股份有限公司 虚拟源成像的方法及装置
CN113740906A (zh) * 2021-10-19 2021-12-03 中国地质大学（北京） 一种水下垂直缆地震波干涉成像方法及装置

Similar Documents

Publication Publication Date Title
Mousavi et al. Automatic noise-removal/signal-removal based on general cross-validation thresholding in synchrosqueezed domain and its application on earthquake data
US6987706B2 (en) Removal of noise from seismic data using high resolution radon transformations
US7953556B2 (en) Geophone noise attenuation and wavefield separation using a multi-dimensional decomposition technique
US9250341B2 (en) Device and method for extrapolating specular energy of reverse time migration three dimensional angle gathers
US20210311218A1 (en) Surface wave prediction and removal from seismic data
CN113156515A (zh) 一种声波远探测成像降噪处理方法和装置
CN111257939B (zh) 一种时移地震虚拟震源双向波场重构方法和系统
CN111257938A (zh) 基于小波互相关时移地震虚拟震源波场重构方法和系统
CN109581481B (zh) 一种便携式高频可控震源地震信号谐波干扰消除方法
Verschuur et al. Multiple removal results from Delft University
CA2497296C (en) Removal of noise from seismic data using improved radon transformations
Ventosa et al. Window length selection for optimum slowness resolution of the local-slant-stack transform
Zhao et al. Wavelet-crosscorrelation-based interferometric redatuming in 4D seismic
Wang et al. Multicomponent seismic noise attenuation with multivariate order statistic filters
CN112505782B (zh) 四维地震中辐射模式校正的干涉基准面重建方法及系统
Staring et al. R-EPSI and Marchenko equation-based workflow for multiple suppression in the case of a shallow water layer and a complex overburden: A 2D case study in the Arabian Gulf
Yu et al. VSP and DAS data noise attenuation using complex wavelet transform
US9014985B2 (en) System and method for compensating time and offset varying near-surface effects in seismic data background
Backus et al. Multiple reflections as an additive noise limitation in seismic reflection work
Telling et al. Evaluation of a broadband marine source
Zhou et al. Applying wavelet transform to suppress ghost in ocean-bottom node dual-sensor technology
Lu et al. Receiver deghosting in the tx domain based on super-Gaussianity
CN114002740A (zh) 一种低信噪比地震初至信号自动增强方法及系统
Akuhara et al. Application of Inverse Water-Layer Filter Method
Hu et al. Coherent noise attenuation for passive seismic data based on iterative two-dimensional model shrinkage

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
RJ01 Rejection of invention patent application after publication
RJ01 Rejection of invention patent application after publication

Application publication date: 20200609