CN101539634A - 通过防混淆、防泄漏傅立叶变换内插地震数据的方法 - Google Patents

通过防混淆、防泄漏傅立叶变换内插地震数据的方法 Download PDF

Info

Publication number
CN101539634A
CN101539634A CN200910128951A CN200910128951A CN101539634A CN 101539634 A CN101539634 A CN 101539634A CN 200910128951 A CN200910128951 A CN 200910128951A CN 200910128951 A CN200910128951 A CN 200910128951A CN 101539634 A CN101539634 A CN 101539634A
Authority
CN
China
Prior art keywords
frequency
fourier transform
obscuring
wavenumber
frequency component
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.)
Granted
Application number
CN200910128951A
Other languages
English (en)
Other versions
CN101539634B (zh
Inventor
M·A·肖纳维勒
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.)
PGS Geophysical AS
Original Assignee
PGS Geophysical AS
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 PGS Geophysical AS filed Critical PGS Geophysical AS
Publication of CN101539634A publication Critical patent/CN101539634A/zh
Application granted granted Critical
Publication of CN101539634B publication Critical patent/CN101539634B/zh
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/36Effecting static or dynamic corrections on records, e.g. correcting spread; Correlating seismic signals; Eliminating effects of unwanted energy
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00Details of seismic processing or analysis
    • G01V2210/50Corrections or adjustments related to wave propagation
    • G01V2210/57Trace interpolation or extrapolation, e.g. for virtual receiver; Anti-aliasing for missing receivers

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)
  • Complex Calculations (AREA)

Abstract

本发明涉及通过防混淆、防泄漏傅立叶变换内插地震数据的方法。通过将第一防泄漏傅立叶变换方法应用于时间变换后的地震数据中的混淆频率分量以及将第二防泄漏傅立叶变换方法应用于时间变换后的地震数据中的非混淆频率分量来产生估计频率-波数谱。第二防泄漏傅立叶变换方法应用从非混淆频率外推到混淆频率的绝对频率-波数谱来加权混淆频率的频率-波数分量。将逆时间和空间傅立叶变换应用于估计频率-波数谱,从而产生该地震数据的轨迹内插。

Description

通过防混淆、防泄漏傅立叶变换内插地震数据的方法
技术领域
本发明主要是关于地球物理勘探颌域。更具体地,本发明是关于在地震数据中内插轨迹或规则化轨迹的领域。
背景技术
在油气产业中,地球物理勘探通常用于帮助探索和评价地层。地球物理勘探技术产生地球的地下结构的知识,该知识有益于发现和提取有价值的矿物资源、尤其是诸如石油和天然气的碳氢化合物沉积物。地球物理勘探的公知技术是地震测量。在陆基地震测量中,地震信号产生在地球的表面上或在地球的表面附近,随后向下行进到地球的次表面内。在海洋地震测量中,地震信号也可向下行进通过覆盖地球的次表面的水体。地震能量源被用于产生地震信号,在传播进入地球后,该地震信号至少部分地由次表面地震反射器反射。这样的地震反射器通常是具有不同的弹性特性的地层之间的界面,该弹性特性特别是声波速度和岩石密度,不同的弹性特性导致在界面处的声阻抗的差异。被反射的地震能量在地球的表面处或地球的表面附近、在覆盖的水体中、或在钻孔中的已知深度处由地震传感器(也称为地震接收器)检测并被记录。
在陆地地震测量中用于产生地震信号的适当地震源可包括爆炸物或振动器。海洋地震测量通常采用由船牵引且被周期性地激活以产生声波场的水下地震源。产生该波场的地震源可以是几种类型,包括小型爆炸装填物、电火花或电弧、海洋振动器,以及通常包括枪。该地震源枪可以是水枪、蒸汽枪以及最通常地是气枪。通常,海洋地震源不是包括单一源元件,而是包括空间分布的源元件阵列。对于当前最通用的海洋地震源的形式-气枪,该布置尤为准确。
适当类型的地震传感器通常包括尤其用于陆地测量中的粒子速度传感器以及尤其用于海洋测量中的水压传感器。有时候代替粒子速度传感器或除了粒子速度传感器之外,使用粒子位移传感器、粒子加速传感器或压力梯度传感器。本领域中一般已知粒子速度传感器和水压力传感器分别作为地震检波器和水中地震检波器。地震传感器可单独配置,但更一般地被配置在传感器阵列中。此外,压力传感器和粒子速度传感器可以被一起配置在海洋测量中,成对并置或以阵列对并置。
在典型的海洋地震测量中,地震测量船在水面上通常以大约5节的速度行进,且包含地震采集装置、诸如导航控制、地震源控制、地震传感器控制以及记录装置。地震源控制装置使得在水体中被地震船牵引的地震源以所选择的次数开动。也称作地震电缆的地震拖缆为伸长的像电缆的结构,该结构在水体中由牵引地震源的地震测量船或由其他地震测量船牵引。通常,多个地震拖缆被牵引在地震船后。该地震拖缆包括传感器,以检测由地震源启动且从反射界面被反射的反射波场。常规地,地震拖缆包含诸如水中地震检波器的压力传感器,但是已经提出了这样的地震拖缆,除了水中地震检波器,其包含诸如地震检波器的水粒子速度传感器或诸如加速度计的粒子加速度传感器。压力传感器和粒子运动传感器可以很接近地配置,沿地震电缆成对并置或以阵列对并置。
处理在执行该测量中得到的结果地震数据以产生关于被测量的区域中的地层的特性和地质结构的信息。为了显示和分析这些地层的潜在的碳氢化合物含量,处理该被处理的地震数据。地震数据处理的目的是从该地震数据提取尽可能多的关于地层的信息以便充分地使地质次表面成像。为了识别在地球的次表面中有发现石油聚集的可能性的位置,在搜集、处理和解释地震数据方面花费大量的金钱。根据所记录的地震数据来构造反射器表面(该反射器表面限定感兴趣的地下地球层)的处理提供在深度或时间方面的地球的图像。
地球的次表面的结构的图像被生成以便使得解释者能选择最大可能具有石油聚集的位置。为了检验石油的存在,必须钻井。钻井以确定石油矿床存在与否是个极其昂贵且耗时的任务。由于该原因,存在改进地震数据的处理和显示的持续需求,以便生成地球的次表面的结构的图像(无论是由计算机还是人来进行解释,该图像将提高解释者的能力),以估定石油聚集存在于地球的次表面中特定位置的可能性。
无论在陆地上还是在水中,在地震数据的收集中可能出现两个问题。该数据可以被欠采样(under-sampled)(混淆)或该数据可以被非均匀地(不规则地)采样。在地震测量中的物理和经济限制通常导致地震数据作为欠采样或非均匀采样数据被采集。
欠采样数据通常被称为混淆数据。根据数据采样理论,要求嵌入数据中的波长不短于采样间隔的两倍。否则,相应于该嵌入波长的特征将是欠分解的(under-resolved)并且因此由于混淆而失真。
因此,混淆开始的时间混淆(temporal alias)频率、即Nyquist频率在频率-波数(f-k)域中为 f = 1 2 Δt . 混淆开始的空间混淆频率、即Nyquist波数在空间坐标中为 k = π Δx . 这里Δt是以毫秒为单位的采样时间间隔,而Δx是以优选单位、诸如米为单位的站间距。因此,在时间-空间(t-x)域中的大的采样间隔相应于相应f-k域中的小的Nyquist频率和Nyquist波数。
采用比理想要求的采样间隔大得多的采样间隔来记录地震数据会导致在随后的数据处理中的有害影响。然而,尤其是在3D测量的情况下,在地震测量期间以更精细的采样间隔收集数据大大增加地震数据采集的成本。因此,替代地,可以根据所采集的数据近似缺失的数据。因此,必须通过内插或外推空间上混淆的地震数据来寻求某些益处。
在非混淆的、均匀采样的地震数据中轨迹的内插是简明的。例如,可以通过在空间域中与sinc滤波器卷积或通过在傅立叶域中经由零填充来扩展带限信号的Nyquist波数来执行该内插。然而,这种较容易的轨迹内插假定该内插采用正交基函数来进行。当轨迹内插在不规则采样的格栅上执行时,数据中信号的能量向所有其他频率泄漏。该能量泄漏由采样的不规则性和边界效应引起。
因此,非均匀采样的数据需要被规则化为正交(规则)基本格栅。三维地震数据规则化要求在地震测量期间没有发生在源和接收器位置处的采集的地点产生地震轨迹。换句话说,来自不规则格栅上所采集的数据的地震轨迹被内插或外推到规则格栅。
因为包括电缆水平偏转、障碍回避、坏的轨迹的编辑以及经济的许多原因,特别是海洋地震数据通常被沿空间方向不规则地且稀疏地采样。然而,对于包括3D表面相关的多次波消除(multiple elimination)和3D波动方程偏移的若干地震处理应用来说,需要规则采样的数据。获得3D规则采样的数据的最佳方法是在交叉线方向以更多的冗余且以更宽的方位距离采集更多的数据,但是因为当前海洋采集技术的原因这一点要实现是困难且昂贵的。因此,数据规则化成为重要的处理工具。
因此,存在对用于在欠采样且非均匀采样的地震数据中内插轨迹的方法的需要。特别是,存在对轨迹内插方法的需要,该方法消弱由于混淆地震数据中的不规则采样造成的能量泄漏。
发明内容
本发明是用于在地震数据中内插轨迹的方法,该地震数据可以是欠采样且非均匀采样的。通过将第一防泄漏傅立叶变换方法应用于时间变换的地震数据中的混淆频率分量以及将第二防泄漏傅立叶变换方法应用于时间变换的地震数据中的非混淆频率分量来产生估计的频率-波数谱。第二防泄漏傅立叶变换方法应用从非混淆频率外推到混淆频率的绝对频率-波数谱来加权混淆频率的频率-波数分量。将逆时间和空间傅立叶变换应用于估计的频率-波数谱,从而产生该地震数据的轨迹内插。
附图说明
通过参考下列详细描述以及附图可以更容易地理解本发明及其优点,其中:
图1是说明用于在欠采样且非均匀采样的地震数据中内插轨迹的本发明的第一实施例的流程图;
图2是说明用于在欠采样且非均匀采样的地震数据中内插轨迹的本发明的第二实施例的初始部分的流程图;
图3示出说明在图2中开始的用于在地震数据中内插轨迹的本发明的第二实施例的中间部分的流程图。
图4示出说明在图2中开始的用于在地震数据中内插轨迹的本发明的第二实施例的最终部分的流程图。
图5是说明用于处理图2的非均匀采样的地震数据中的非混淆频率的本发明的实施例的流程图;以及
图6是说明用于处理图3的非均匀采样的地震数据中的混淆频率的本发明的实施例的流程图。
图7示出在较低非混淆频率处的频率-波数谱的图表;
图8示出用作用于较高混淆频率的加权函数的外推频率-波数谱的图表;
图9A-9D示出应用于合成数据例子的本发明的方法;
图10A-10D示出应用于图9A-9D中的合成数据例子的标准ALFT方法;
图11A-11D示出应用于场数据例子的本发明的方法;以及
图12A-12D示出应用于图11A-11D中的场数据例子的标准ALFT方法。
尽管将结合本发明的优选实施例描述本发明,但可以理解的是本发明不限于这些优选实施例。相反,本发明意图覆盖可被包括在由所附的权利要求限定的本发明的范围内的所有替换、修改和等价方案。
具体实施方式
最简明的轨迹内插方法为在每个时间样本处通过在波数(k)方向上的傅立叶变换、填充足够大数目的零样本以及逆傅立叶变换来在空间方向上进行1D内插。这种轨迹内插可以在频率-波数(f-k)域中针对数据的每个频率切片(slice)来等价地进行。该过程是确定性数据独立、sinc(或正弦基数)函数、轨迹内插器。只要地震数据不是欠采样的(混淆的)或非均匀采样的(在不规则格栅上),该过程就顺利执行。
作为地震数据规则化的轨迹内插是内插/外推问题,轨迹内插旨在利用所采集的非规则采样数据来估计在空间规则格栅上的地震轨迹。傅立叶变换起着估计频率-波数域中的频率分量的至关重要的作用,且其逆傅立叶变换在时间-空间域中在要求的规则格栅背面(back)再创建地震数据。基于傅立叶变换的数据重建的基本问题是诸如Sinc函数的全局基函数在非规则格栅上不正交。全局傅立叶基的非正交性导致能量从一个频率分量泄漏到其他分量上,被称为“谱泄漏”的现象。
本发明是一种用于在地震数据中内插轨迹的方法,该地震数据可以是欠采样的且非均匀采样的。本发明建立于用于内插缺失轨迹的“防泄漏(anti-leakage)傅立叶变换”(ALFT)方法。该ALFT被提出用于通过为在不规则格栅上采样的数据减少频率泄漏现象来解决地震数据规则化问题。ALFT在不规则采样的格栅上“再正交化(re-orthogonalize)”全局傅立叶基函数,这导致不规则格栅上的信号谱的良好估计。
对于标准傅立叶求和方法,求解每个傅立叶系数的次序对最终结果没有影响。然而,在ALFT方法中顺序是至关重要的,因为具有较大量值(能量)的傅立叶系数将比具有较小量值的系数贡献更多的能量泄漏。因此,为了减少泄漏,迭代地估计傅立叶系数,在每个点估计具有最大能量的系数。在每次估计后,通过更新输入数据,所计算的频率分量(傅立叶系数)将被重置为零。数学上,其等价于从原始输入地震数据中除去频率分量。
再次使用相同的最大能量标准,该新近被减的输入被用于求解下个傅立叶系数。重复该迭代过程直到所有傅立叶系数被解出,即,直到更新的输入中的所有值趋向于零(尤其是,低于阈值)。一般来说,全局傅立叶基函数仅仅在规则格栅上,即,对于规则采样的数据来说是正交的。该减法作为在不规则格栅上傅立叶基的正交化机制。换句话说,该傅立叶基函数被再正交化。这一方法导致用于使从一个频率到另一个频率的泄漏影响最小化的实际解决方案。
如果所使用的傅立叶分量的数目和范围是足够的,则在所有减法操作后不规则格栅上最终的更新的输入数据将趋向于零。在这种情况下,根据所获得的傅立叶系数重建的数据将与原始测量相匹配,对要求的内插方法的要求之一。
标准ALFT方法的问题是处置混淆数据(在有噪声的情况下)。混淆分量可以具有等于或大于非混淆分量的量值且可能被ALFT方法错误地从序列中挑出。本发明的方法通过使用来自较低非混淆时间频率的信息来详述ALFT方法,以帮助“解混淆(de-alias)”较高时间频率。尤其是,本发明使用非混淆较低频率来设计加权函数,该加权函数确定在ALFT过程中首先计算(并除去)哪些谱分量。通过在非混淆频率处将频率-波数谱外推到混淆频率来构造加权函数。
图1-6是说明用于轨迹内插的本发明的实施例的流程图。
图1和图2-4分别示出本发明的方法的两个实施例,而图5和6示出图2-4中所示的方法的进一步细节。图7-12说明在参考图1-6讨论的流程图中描述的技术中的一些。图6和7说明用于加权技术中的谱。图9A-9D和10A-10D就合成(synthetic)例子分别说明本发明的方法和标准ALFT的轨迹内插过程的结果。图11A-11D和12A-12D就场数据例子分别说明本发明的方法和标准ALFT的轨迹内插过程的结果。
图1示出用于在欠采样且非均匀采样的输入地震数据集中内插轨迹的本发明的第一实施例的流程图。
在框11,产生估计的频率-波数谱。将第一标准防泄漏傅立叶变换方法应用于时间变换后的(temporal-transformed)地震数据的非混淆频率分量,并且将第二非标准防泄漏傅立叶变换方法应用于时间变换后的地震数据的混淆频率分量。第二防泄漏傅立叶变换方法应用从(通常较低的)非混淆频率被外推到(通常较高的)混淆频率的绝对频率-波数谱以加权混淆频率。
在框12,将逆时间和空间傅立叶变换应用于来自框11的估计的频率-波数谱,从而产生估计的数据集。估计的频率-波数谱从频率-波数(f-k)域被逆变换回时间-空间(t-x)域,导致输入地震数据的要求的轨迹内插。
图2示出用于在欠采样且非均匀采样的地震数据中内插轨迹的本发明的第二实施例的初始部分的流程图。图2详述上面参考图1讨论的第一实施例的讨论。
在框21,获得输入地震数据。该地震数据被假定为在时间-空间(t-x)域中被表示。地震数据可以是混淆的和不规则采样的,因为本发明的方法适合于处置这两种情况。
在框22,时间傅立叶变换被应用于在框21中获得的输入地震数据。优选地,为了计算效率,所应用的傅立叶变换是快速傅立叶变换(FFT)。该输入数据从时间-空间(t-x)域被变换到频率-空间(f-x)域。
在框23,确定来自框22的变换后的输入地震数据的哪些频率分量为非混淆的以及哪些频率分量是混淆的。通常,较低的频率将是非混淆的而较高的频率将是混淆的。
在框24,为第一ALFT方法选择被指定为N1的第一迭代数目,该第一ALFT方法将被应用于在框23中确定的非混淆频率。根据经验为第一ALFT方法得出数目N1
在替换实施例中,代替第一迭代数目N1,为第一ALFT方法选择第一阈值ε1。随后,将第一ALFT方法迭代地应用于每个傅立叶系数,直到该系数减小到第一阈值ε1以下,而不是针对设定的迭代数目。仅仅为了说明的目的针对使用该迭代数目的实施例来说明本发明的方法,且该选择不应被认为是对本发明的限制。
在框25,从在框23中被确定为非混淆的频率分量中选择非混淆频率分量。
在框26,该过程以在框25中选择的非混淆频率分量进行到图5的框51。在那里,在图5中,第一(标准)ALFT方法将被应用于非混淆频率。
在框27,从图5的框59获得所选非混淆频率分量的估计频率-波数谱。
在框28,确定是否剩下更多的非混淆频率分量要被选择。如果剩下更多的非混淆频率分量,则该过程返回到框25。如果没有剩下更多的非混淆频率分量,则该过程继续到下个框29。
在框29,来自框27的所选非混淆频率分量的所有估计频率-波数谱被组合为非混淆估计频率-波数谱。
在框30,对于图2来说该过程结束且以来自框29的非混淆估计频率-波数谱进行到图3的框31。
图3示出说明在图2中开始的用于在地震数据中内插轨迹的本发明的第二实施例的中间部分的流程图。
在框31,该过程从图2的框29以非混淆估计频率-波数谱继续。
在框32,取来自框31的非混淆估计频率-波数谱的绝对值,从而在非混淆频率的频率-波数(f-k)域中产生绝对频率谱。
在框33,来自框32的绝对频率-波数谱被外推到混淆频率,从而在混淆频率的频率-波数(f-k)域中产生外推的绝对频率-波数谱。通常,该混淆频率包括较高的频率和波数。该外推的谱实际上(effectively)包含来自非混淆的较低频率的信息,且该信息将给出在混淆的较高频率处非混淆傅立叶分量的改进的选择。
图7-8说明将在框32和33中计算的示例谱。图7示出将在框32中计算的在非混淆的较低频率处的频率-波数谱的图表。图8示出将在框33中计算的外推的频率-波数谱的图表。该外推的谱被用作较高频率的加权函数,这些较高频率中的一些为混淆的。
原则上,较低频率被外推到较高频率且因此被外推到较大的带宽或波数。在实践中,在频率和波数值中均内插较低频率谱。该内插可包括平均(averaging)或平滑(smoothing)。图8示出2:1内插的结果。
在框34,为第二ALFT方法选择被指定为N2的第二迭代数目,该第二ALFT方法将被应用于在图2中的框23中确定的混淆频率。根据经验为第二ALFT方法得出数目N2。尽管通常N1>N2,但在一个实施例中,N1=N2。然而,为每个频率分量确定不同的迭代数目N,即迭代数目是频率相关的,这在本发明的范围内。
在替换实施例中,代替第二迭代数目N2,为第二ALFT方法选择第二阈值ε2。然后,第二ALFT方法被迭代地应用于每个傅立叶系数,直到该系数减小到第二阈值ε2以下,而不是针对设定的迭代数目。仅仅为了说明的目的针对使用该迭代数目的实施例来说明本发明的方法,且该选择不应被认为是对本发明的限制。
在框35,从在图2的框23中被确定为混淆的频率中选择混淆频率分量。
在框36,该过程以在框35中选择的混淆频率分量进行到图6的框61。在那里,在图6中,本发明的第二(非标准)ALFT方法将被应用于混淆频率。
在框37,从图6的框71获得所选混淆频率分量的估计频率-波数谱。
在框38,确定是否剩下更多的混淆频率分量要被选择。如果剩下更多的混淆频率,则该过程返回到框35。如果没有剩下更多的混淆频率,则该过程继续到下个框39。
在框39,图3中的过程结束,且以来自框31的非混淆估计频率-波数谱以及从框37获得的所选混淆频率分量的所有估计频率-波数谱进行到图4的框41。
图4示出说明在图2中开始且在图3中继续的用于在地震数据中内插轨迹的本发明的第二实施例的最终部分的流程图。
在框41,该过程以非混淆估计频率-波数谱和所选混淆频率分量的所有估计频率-波数谱从图3的框39继续。
在框42,来自框41的所选混淆频率分量的所有估计频率-波数谱被组合为混淆估计频率-波数谱。
在框43,分别来自框41和42的非混淆和混淆估计频率-波数谱被组合为总估计频率-波数谱。
在框44,逆空间傅立叶变换被应用于来自框43的总估计频率-波数谱,从而产生总估计频率-空间谱。逆空间傅立叶变换被设计为将轨迹变换到包括缺失轨迹的位置在内的要求的轨迹位置或规则(正交)格栅上的位置。优选地,为了计算效率,所应用的逆傅立叶变换为逆空间离散傅立叶变换(DFT)或非均匀快速傅立叶变换(NFFT)。总估计频率-波数谱从频率-波数(f-k)域被逆变换到频率-空间(f-x)域中的总估计频率-空间谱。
在框45,逆时间傅立叶变换被应用于来自框44的变换后的总估计频率-空间谱,从而产生总估计数据集。优选地,所应用的逆傅立叶变换为逆时间快速傅立叶变换(FFT)。总估计频率-空间谱进一步从频率-空间(f-x)域被逆变换到时间-空间(t-x)域中的总估计数据集。
替换地,在框44和45中应用的逆傅立叶变换为2D(时间和空间)快速傅立叶变换(FFT)。在任一情况下,最终结果为估计频率-波数谱从频率-波数(f-k)域被逆变换回时间-空间(t-x)域,从而产生输入地震数据的要求的轨迹内插。输入地震数据的该轨迹内插可实现许多目标,包括填充缺失的轨迹和使采样轨迹规则化。
针对在图2的框23中确定的非混淆频率应用标准ALFT。图5示出说明用于处理图2的非均匀采样的地震数据中的非混淆频率的本发明的实施例的流程图。
在框51,从图2的框26获得所选非混淆频率分量。
在框52,针对在框51中获得的所选非混淆频率分量建立频率-波数(f-k)域中的估计频率-波数谱。该估计频率-波数谱最初被设定为零。通过在下面的框55中增加所选波数分量将进一步建立该谱。
在框53,空间傅立叶变换被应用于在框51中选择的所选非混淆频率分量。优选地,为了计算效率,所应用的傅立叶变换为离散傅立叶变换(DFT)或非均匀快速傅立叶变换(NFFT)。频率分量从频率-空间(f-x)域被变换到频率-波数(f-k)域。
在框54,选择来自框53的变换后的频率分量中的最强波数分量。该最强波数分量是由在框53中计算的空间傅立叶变换所得到的傅立叶分量,其具有最大的量值(能量)。
在框55,在框54中选择的最强波数分量被增加到所选非混淆频率分量的在框52中建立并初始化的估计频率-波数谱。
在框56,逆空间傅立叶变换被应用于在框54中所选择的最强傅立叶分量。优选地,为了计算效率,所应用的逆傅立叶变换为逆离散傅立叶变换(DFT)或逆非均匀快速傅立叶变换(NFFT)。最强傅立叶分量从频率-波数(f-k)域被逆变换回频率-空间(f-x)域。
在框57,从在框51中获得的非混淆频率分量减去在框56中计算的逆变换后的最强分量。这些减法迭代地产生纠正的非混淆频率分量。
在框58,确定针对在框51中获得的非混淆频率分量是否发生了框53到57的N1次迭代。如果N1次迭代没有发生,则该过程返回到框53以进行另一迭代。如果发生了N1次迭代,则图5中的过程继续到框59。
在图2的框24中描述的替换实施例中,确定来自框57的纠正的非混淆频率分量是否已经在第一阈值ε1以下。如果没有在阈值ε1以下,则该过程返回到框53以进行另一迭代。如果在阈值ε1以下,则图5中的过程继续到框59。
在框59,对于图5来说该过程结束,且以在框55中迭代地构造的所选非混淆频率分量的估计频率-波数谱返回到图2中的框27。
针对在图2的框27中确定的混淆频率,应用本发明的非标准ALFT。图6示出用于处理来自图3的非均匀采样的地震数据中的混淆频率的本发明的实施例的流程图。
在框61,从图3的框36获得所选混淆频率分量。
在框62,针对在框51中获得的所选混淆频率分量建立频率-波数(f-k)域中的估计频率-波数谱。该估计频率-波数谱最初被设定为零。通过在下面的框67中增加所选波数分量将进一步建立该谱。
在框63,将空间傅立叶变换应用于在框61中选择的所选混淆频率分量。优选地,为了计算效率,所应用的傅立叶变换为离散傅立叶变换(DFT)或非均匀快速傅立叶变换(NFFT)。该频率分量从频率-空间(f-x)域被变换到频率-波数(f-k)域。
在框64,来自图3的框33的外推绝对频率-波数谱被应用于在框63中计算的变换后的频率分量,以在变换后的频率分量中加权波数分量。在图8中示出了该外推频率-波数谱的一个例子。
在框65,选择来自框64的变换和加权后的频率分量中的最强波数分量。该最强波数分量是由在框53中计算的空间傅立叶变换所得到的傅立叶分量,其具有最大的量值(能量)。
在框66,获得来自框63的变换后的频率分量中的未加权波数分量,其相应于在框65中确定的最强加权波数分量。
在框67,将在框66中获得的相应于最强波数分量的未加权波数分量增加到所选混淆频率分量的在框62中建立和初始化的估计频率-波数谱中。
在框68,将逆空间傅立叶变换应用于在框66中确定的最强未加权分量。优选地,为了计算效率,所应用的逆傅立叶变换为逆离散傅立叶变换(DFT)或逆非均匀快速傅立叶变换(NFFT)。该最强未加权分量从频率-波数(f-k)域被逆变换回频率-空间(f-x)域。
在框69,从在框61中获得的混淆频率分量减去在框68中计算的逆变换后的最强未加权分量。这些减法迭代地产生纠正的混淆频率分量。
在框70,确定针对在框51中获得的混淆频率分量是否发生了框63到69的N2次迭代。如果N2次迭代没有发生,则该过程返回到框63以进行另一迭代。如果发生了N2次迭代,则该过程继续到框71。
在图3的框34中描述的替换实施例中,确定来自框69的纠正的非混淆频率分量是否已经在第二阈值ε2以下。如果没有在阈值ε2以下,则该过程返回到框63以进行另一迭代。如果在阈值ε2以下,则图6中的过程继续到框71。
在框71,对于图6来说该过程结束,且以在框57中迭代地构造的所选混淆频率分量的估计频率-波数谱返回到图3中的框37。
图9A-12D示出本发明的方法和标准ALFT的轨迹内插过程的结果,以做比较。图9A-9D和10A-10D说明合成例子的轨迹内插过程的结果,而图11A-11D和12A-12D说明场数据例子的轨迹内插过程的结果。
图9A-9D示出应用于合成数据例子的本发明的方法。图9A示出原始合成地震数据。图9B示出除去轨迹以模拟混淆数据的输入数据。图9C示出采用本发明的方法的内插数据。图9D示出图9A中的原始数据和图9C中的内插数据之间的差异,其中小的差异指示紧密的一致性。
图10A-10D示出应用于图9A-9D中的合成数据例子的标准ALFT方法,以做比较。图10A示出如图9A中的原始合成地震数据。图10B示出如图9B中的除去轨迹以模拟混淆数据的输入数据。图10C示出采用标准ALFT方法的内插数据。图9D示出图10A中的原始数据和图10C中的内插数据之间的差异,其中较大的差异指示不像本发明的方法的上面图9D中那样紧密的一致性。
图11A-11D示出应用于场数据例子的本发明的方法。图11A示出原始场地震数据。图11B示出除去轨迹以模拟混淆数据的输入数据。图11C示出采用本发明的方法的内插数据。图11D示出图11A中的原始数据和图11C中的内插数据之间的差异。
图12A-12D示出应用于图11A-11D中的场数据例子的标准ALFT方法,以做比较。图12A示出如图11A中的原始场地震数据。图12B示出如图11B中的除去轨迹以模拟混淆数据的输入数据。图12C示出采用标准ALFT方法的内插数据。图12D示出图12A中的原始数据和图12C中的内插数据之间的差异。
图11A中的原始数据和图11C(本发明)中的内插数据之间的归一化均方根(NRMS)差为70%,而图12A中的原始数据和图12C(标准ALFT)中的内插数据之间的NRMS差为84%。因此,本发明的方法示出同样对于场数据例子来说比标准ALFT更好的结果。
本发明的防混淆(anti-alias)、防泄漏傅立叶变换方法提供比单独的标准ALFT更好的关于混淆地震数据的轨迹内插。估计计算的额外成本是非常有限的。如在此所公开的,本发明的方法可容易地扩展到多维实施例,包括具有两个空间维度加时间的3D,具有三个空间维度加时间的4D,以及具有四个空间维度加时间的5D。空间维度可包括源x,y和接收器x,y坐标的子集,或者等价地,纵测线(inline)中点、横测线中点、偏移和方位角坐标的子集。在替换实施例中,时间坐标可以是频率或深度坐标。
本发明的方法还可扩展到一和多维凸集投影(POCS)图像恢复算法。其他变型方案也是可能的,包括平滑频率-波数谱以产生更好的权重和应用其他加权方案。其他变型方案包括但不限于使用更高频率,或来自其他集合的信息,或使用不同方法来估计在较低频率处的频率-波数谱,诸如最小平方傅立叶变换。
应该理解的是前述仅仅是对本发明的具体实施例的详细描述且在不脱离本发明的范围的情况下可以根据本公开内容对所公开的实施例做出多种改变、修改和替换。因此,前述描述并不意为限制本发明的范围。更确切的说,本发明的范围仅由所附的权利要求及它们的等价确定。

Claims (22)

1.一种用于在地震数据中内插轨迹的方法,包括:
通过将第一防泄漏傅立叶变换方法应用于时间变换后的地震数据的非混淆频率分量以及将第二防泄漏傅立叶变换方法应用于时间变换后的地震数据的混淆频率分量来产生估计频率-波数谱,其中第二防泄漏傅立叶变换方法应用从非混淆频率外推到混淆频率的绝对频率-波数谱来加权混淆频率的频率-波数分量;以及
将逆时间和空间傅立叶变换应用于估计频率-波数谱,从而产生该地震数据的轨迹内插。
2.如权利要求1所述的方法,其中该第一防泄漏傅立叶变换方法为标准防泄漏傅立叶变换方法。
3.如权利要求1所述的方法,其中非混淆频率基本上为地震数据中较低的频率而混淆频率基本上为地震数据中较高的频率。
4.如权利要求2所述的方法,最初包括:
获得时间-空间域中的输入地震数据;
将时间傅立叶变换应用于输入地震数据,从而产生变换后的地震数据;
确定变换后的地震数据中的哪些频率分量相应于非混淆频率以及变换后的地震数据中的哪些频率分量相应于混淆频率;以及
将第一防泄漏傅立叶变换方法应用于来自变换后的地震数据的非混淆频率分量中的每一个。
5.如权利要求4所述的方法,其中时间傅立叶变换为快速傅立叶变换。
6.如权利要求4所述的方法,进一步包括:
为第一防泄漏傅立叶变换方法确定迭代数目N1
7.如权利要求6所述的方法,其中确定迭代数目N1包括:
为第一防泄漏傅立叶变换方法确定阈值ε1
8.如权利要求7所述的方法,进一步包括:
针对非混淆频率分量将估计频率-波数谱初始化为零;
将空间傅立叶变换应用于所选非混淆频率分量;以及
针对N1次迭代执行下列步骤:
在变换后的非混淆频率分量中选择最大波数分量;
将最大波数分量增加到非混淆频率分量的估计频率-波数谱中;
将逆空间傅立叶变换应用于所选最大波数分量;以及
从所选非混淆频率分量减去逆变换后的最大分量,从而产生纠正的频率分量。
9.如权利要求8所述的方法,其中空间傅立叶变换为离散傅立叶变换。
10.如权利要求8所述的方法,其中逆空间傅立叶变换为逆离散傅立叶变换。
11.如权利要求8所述的方法,其中空间傅立叶变换为非均匀快速傅立叶变换。
12.如权利要求8所述的方法,其中逆空间傅立叶变换为逆非均匀快速傅立叶变换。
13.如权利要求8所述的方法,进一步包括:
组合非混淆频率分量的估计频率-波数谱,从而产生非混淆估计频率-波数谱;
取非混淆估计频率-波数谱的绝对值,从而产生绝对频率-波数谱;
将绝对频率-波数谱外推到非混淆频率,从而产生外推绝对频率-波数谱;以及
将第二防泄漏傅立叶变换方法应用于来自变换后的地震数据的混淆频率分量中的每一个。
14.如权利要求13所述的方法,进一步包括:
在将绝对频率-波数谱外推到非混淆频率之前将绝对频率-波数谱平滑。
15.如权利要求13所述的方法,进一步包括:
为第二防泄漏傅立叶变换方法确定迭代数目N2
16.如权利要求15所述的方法,其中确定迭代数目N2包括:
为该第二防泄漏傅立叶变换方法确定阈值ε2
17.如权利要求8所述的方法,进一步包括:
针对混淆频率分量将估计频率-波数谱初始化为零;
在所选混淆频率分量处应用空间傅立叶变换;以及
针对N2次迭代执行下列步骤:
将外推绝对频率-波数谱应用于变换后的混淆频率分量,从而产生加权混淆频率-波数谱;
在加权混淆频率-波数谱中选择最大波数分量;
获得相应于最大加权波数分量的未加权波数分量;
将相应未加权波数分量增加到混淆频率分量的估计频率-波数谱中;
计算最大未加权波数分量的逆空间傅立叶变换;以及
从混淆频率分量减去所计算的逆变换后的最大波数分量,从而产生纠正的混淆频率分量。
18.如权利要求17所述的方法,进一步包括:
将逆空间傅立叶变换应用于所选混淆频率分量,从而返回频率-空间域产生内插的地震数据。
19.如权利要求18所述的方法,其中该逆空间傅立叶变换为逆离散傅立叶变换。
20.如权利要求19所述的方法,其中该逆空间傅立叶变换为逆非均匀快速傅立叶变换。
21.如权利要求18所述的方法,进一步包括:
将逆时间傅立叶变换应用于所选混淆频率分量,从而返回时间-空间域产生内插的地震数据。
22.如权利要求21所述的方法,其中该逆傅立叶变换为逆快速傅立叶变换。
CN200910128951.6A 2008-03-17 2009-03-17 通过防混淆、防泄漏傅立叶变换内插地震数据的方法 Expired - Fee Related CN101539634B (zh)

Applications Claiming Priority (3)

Application Number Priority Date Filing Date Title
US12/077,108 US7751277B2 (en) 2008-03-17 2008-03-17 Method for interpolating seismic data by anti-alias, anti-leakage Fourier transform
US12/077,108 2008-03-17
US12/077108 2008-03-17

Publications (2)

Publication Number Publication Date
CN101539634A true CN101539634A (zh) 2009-09-23
CN101539634B CN101539634B (zh) 2014-02-05

Family

ID=40848502

Family Applications (1)

Application Number Title Priority Date Filing Date
CN200910128951.6A Expired - Fee Related CN101539634B (zh) 2008-03-17 2009-03-17 通过防混淆、防泄漏傅立叶变换内插地震数据的方法

Country Status (11)

Country Link
US (1) US7751277B2 (zh)
EP (1) EP2103959B1 (zh)
CN (1) CN101539634B (zh)
AU (1) AU2009200673B2 (zh)
BR (1) BRPI0900890B1 (zh)
CA (1) CA2658300C (zh)
EA (1) EA014282B1 (zh)
EG (1) EG26391A (zh)
MX (1) MX2009002932A (zh)
MY (1) MY150168A (zh)
SG (1) SG155833A1 (zh)

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102636811A (zh) * 2012-04-10 2012-08-15 恒泰艾普石油天然气技术服务股份有限公司 一种海上二维地震资料中多次波的消除方法
CN105659111A (zh) * 2013-10-11 2016-06-08 雪佛龙美国公司 用于正则化地震数据的系统和方法
CN107787455A (zh) * 2015-06-23 2018-03-09 西门子股份公司 用于分析信号的方法以及用于执行该方法的装置
CN109001800A (zh) * 2018-07-20 2018-12-14 中国石油天然气股份有限公司 一种基于地震数据的时频分解与气藏检测方法及系统
CN109188535A (zh) * 2018-09-18 2019-01-11 中国科学院地质与地球物理研究所 地球物理数据处理的方法和装置

Families Citing this family (25)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8014950B2 (en) * 2008-06-19 2011-09-06 Chevron U.S.A. Inc. System and method for seismic trace analysis
US8321134B2 (en) * 2008-10-31 2012-11-27 Saudi Arabia Oil Company Seismic image filtering machine to generate a filtered seismic image, program products, and related methods
US8239135B2 (en) * 2009-05-07 2012-08-07 Pgs Geophysical As Method for calculation of seismic attributes from seismic signals
US8619498B2 (en) * 2010-09-24 2013-12-31 CGGVeritas Services (U.S.) Inc. Device and method for calculating 3D angle gathers from reverse time migration
US9043155B2 (en) 2010-10-07 2015-05-26 Westerngeco L.L.C. Matching pursuit-based apparatus and technique to construct a seismic signal using a predicted energy distribution
US8862408B2 (en) * 2011-09-28 2014-10-14 Westerngeco L.L.C. Determining one or more target positions in an acquisition domain for processing survey data
US9541659B2 (en) 2011-11-18 2017-01-10 Westerngeco L.L.C. Noise removal from 3D seismic representation
RU2488145C1 (ru) * 2012-01-10 2013-07-20 Министерство образования и науки РФ Федеральное государственное бюджетное образовательное учреждение высшего профессионального образования "Уральский государственный горный университет" Способ построения сейсмических изображений геологической среды
US9423518B2 (en) * 2012-02-09 2016-08-23 Pgs Geophysical As Method for processing dual-sensor streamer data with anti-alias protection
US20140121977A1 (en) * 2012-11-01 2014-05-01 Pgs Geophysical As Methods and systems for monitoring a petroleum reservoir
EP2784551A3 (en) * 2013-03-26 2015-10-28 CGG Services SA System and method for interpolating seismic data by matching pursuit in fourier transform
US20160084976A1 (en) * 2013-05-29 2016-03-24 Cgg Services Sa Processing of multi-sensor streamer data
CN104459770B (zh) * 2013-09-24 2017-06-16 中国石油化工股份有限公司 一种高维地震数据规则化方法
US20150276955A1 (en) * 2013-11-06 2015-10-01 Robert H. Brune Method and System for Extending Spatial Wavenumber Spectrum Of Seismic Wavefields On Land Or Water Bottom Using Rotational Motion
US10605941B2 (en) 2014-12-18 2020-03-31 Conocophillips Company Methods for simultaneous source separation
US10261207B2 (en) 2014-12-18 2019-04-16 Pgs Geophysical As Seismic noise mitigation system and method
US10267939B2 (en) 2015-09-28 2019-04-23 Conocophillips Company 3D seismic acquisition
US10809402B2 (en) 2017-05-16 2020-10-20 Conocophillips Company Non-uniform optimal survey design principles
EP3714294B1 (en) * 2017-11-20 2024-01-03 Shearwater Geoservices Software Inc. Offshore application of non-uniform optimal sampling survey design
WO2019100068A1 (en) * 2017-11-20 2019-05-23 Conocophillips Company Offshore application of non-uniform optimal sampling survey design
CN108802820B (zh) * 2018-05-28 2019-10-11 中国石油天然气股份有限公司 一种深度域反假频方法、装置及系统
US11481677B2 (en) 2018-09-30 2022-10-25 Shearwater Geoservices Software Inc. Machine learning based signal recovery
US11346971B2 (en) * 2019-06-26 2022-05-31 Saudi Arabian Oil Company Imaging subterranean features using Fourier transform interpolation of seismic data
US11215725B2 (en) * 2019-07-17 2022-01-04 Saudi Arabian Oil Company Seismic processing workflow for orthogonal wide azimuth 3D surveys
CN113341220B (zh) * 2021-08-05 2021-11-02 中国空气动力研究与发展中心设备设计与测试技术研究所 含噪多频衰减实信号频率估计方法

Family Cites Families (11)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US4594693A (en) * 1983-11-04 1986-06-10 Mobil Oil Corporation Seismic trace interpolation using f-k filtering
US4628492A (en) 1984-01-11 1986-12-09 Mobil Oil Corporation Method of avoiding aliasing in slant stacking of seismic data
US4922465A (en) 1989-05-30 1990-05-01 Geco A/S Interpolation of severely aliased events
US5235556A (en) 1992-01-10 1993-08-10 Halliburton Geophysical Services Inc. Interpolation of aliased seismic traces
GB9320540D0 (en) 1993-10-06 1993-11-24 Ensign Geophysics Ltd Seismic data acquisition
US5677892A (en) * 1996-08-14 1997-10-14 Western Atlas International, Inc. Unaliased spatial trace interpolation in the f-k domain
US5617372A (en) * 1996-08-14 1997-04-01 Western Atlas International, Inc. Unaliased spatial trace interpolation in the f-k domain
US6115726A (en) * 1997-10-03 2000-09-05 Kromos Technology, Inc. Signal processor with local signal behavior
US6943803B1 (en) 1998-09-21 2005-09-13 Evans & Sutherland Computer Corporation Anti-aliased, textured, geocentric and layered fog graphics display method and apparatus
US7027929B2 (en) * 2003-11-21 2006-04-11 Geo-X Systems Ltd. Seismic data interpolation system
US7239578B2 (en) 2005-03-03 2007-07-03 John M. Robinson Removal of noise from seismic data using radon transformations

Cited By (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102636811A (zh) * 2012-04-10 2012-08-15 恒泰艾普石油天然气技术服务股份有限公司 一种海上二维地震资料中多次波的消除方法
CN105659111A (zh) * 2013-10-11 2016-06-08 雪佛龙美国公司 用于正则化地震数据的系统和方法
CN105659111B (zh) * 2013-10-11 2018-02-06 雪佛龙美国公司 用于正则化地震数据的系统和方法
CN107787455A (zh) * 2015-06-23 2018-03-09 西门子股份公司 用于分析信号的方法以及用于执行该方法的装置
CN109001800A (zh) * 2018-07-20 2018-12-14 中国石油天然气股份有限公司 一种基于地震数据的时频分解与气藏检测方法及系统
CN109188535A (zh) * 2018-09-18 2019-01-11 中国科学院地质与地球物理研究所 地球物理数据处理的方法和装置

Also Published As

Publication number Publication date
EA200900300A1 (ru) 2009-10-30
EP2103959A2 (en) 2009-09-23
MY150168A (en) 2013-12-13
BRPI0900890A2 (pt) 2010-04-06
EG26391A (en) 2013-09-26
US7751277B2 (en) 2010-07-06
EP2103959A3 (en) 2011-01-19
AU2009200673B2 (en) 2013-08-29
US20090231956A1 (en) 2009-09-17
CN101539634B (zh) 2014-02-05
CA2658300C (en) 2014-10-28
EA014282B1 (ru) 2010-10-29
SG155833A1 (en) 2009-10-29
MX2009002932A (es) 2009-09-24
CA2658300A1 (en) 2009-09-17
AU2009200673A1 (en) 2009-10-01
BRPI0900890B1 (pt) 2020-10-13
EP2103959B1 (en) 2019-10-30

Similar Documents

Publication Publication Date Title
CN101539634B (zh) 通过防混淆、防泄漏傅立叶变换内插地震数据的方法
AU2012345565B2 (en) Separation of simultaneous source data
CN101556339B (zh) 对不规则接收机位置海洋地震拖缆数据进行消重影的方法
US8379481B2 (en) 3D deghosting of multicomponent or over/under streamer recordings using cross-line wavenumber spectra of hydrophone data
CN101487899B (zh) 三维双传感器拖缆数据中的波场分离方法
US10690793B2 (en) Simultaneous source acquisition and separation on general related sampling grids
EP3004942B1 (en) Systems and methods for de-noising seismic data
US9229123B2 (en) Method for handling rough sea and irregular recording conditions in multi-sensor towed streamer data
US11029432B2 (en) De-aliased source separation method
CN101881836B (zh) 用于根据地震信号来计算地震属性的方法
EP2360495A1 (en) DIP-based corrections for data reconstruction in three-dimensional surface-related multiple prediction
US20150066374A1 (en) Seismic data processing with frequency diverse de-aliasing filtering
US20150198729A1 (en) Regularization of spatially aliased seismic data
Huff et al. Near offset reconstruction for marine seismic data using a convolutional neural network
Gremillion Seismic data interpolation with shaping inversion to zero offset and least-squares flattening
Basu et al. Pre-conditioning of data before PZ summation in OBC survey–a case study
Cunnell et al. High Quality and Efficient Seismic Acquisition for Frontier Deepwater Areas

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
C14 Grant of patent or utility model
GR01 Patent grant
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20140205

Termination date: 20210317

CF01 Termination of patent right due to non-payment of annual fee