CN106556865A - 一种串联型地震信号优化时频变换方法 - Google Patents

一种串联型地震信号优化时频变换方法 Download PDF

Info

Publication number
CN106556865A
CN106556865A CN201611055163.5A CN201611055163A CN106556865A CN 106556865 A CN106556865 A CN 106556865A CN 201611055163 A CN201611055163 A CN 201611055163A CN 106556865 A CN106556865 A CN 106556865A
Authority
CN
China
Prior art keywords
frequency
time
seismic
signal
input signal
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
CN201611055163.5A
Other languages
English (en)
Other versions
CN106556865B (zh
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.)
Chengdu Univeristy of Technology
Original Assignee
Chengdu Univeristy of Technology
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 Chengdu Univeristy of Technology filed Critical Chengdu Univeristy of Technology
Priority to CN201611055163.5A priority Critical patent/CN106556865B/zh
Publication of CN106556865A publication Critical patent/CN106556865A/zh
Application granted granted Critical
Publication of CN106556865B publication Critical patent/CN106556865B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

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
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00Details of seismic processing or analysis
    • G01V2210/40Transforming data representation

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

一种串联型地震信号优化时频变换方法是一种石油勘探数据处理与解释技术,它利用完备集经验模态分解先将原始地震道x(t)分解成若干个不同频段的固有模态函数(IMF)信号分量,并将其作为新的输入信号imfi(t),再利用能量最优准则目标函数求取最优频域窗参数,并利用优化时频分布计算公式次求取每个新输入信号的时频分布TFRi(t,f),最后将得到的每个时频分布进行叠加,即可得到时频分辨率明显提高的地震信号时频谱TFR(t,f)。这一技术既充分利用了地震勘探数据的有效信息,具有很高的时频分辨率,还可以成倍地提高计算速度,不仅提高了地震勘探资料中信息的利用率,有利于提高地震资料解释的可靠性。

Description

一种串联型地震信号优化时频变换方法
技术领域
本发明涉及石油地震勘探数据处理和解释领域,是通过一种串联型地震信号优化时频变换方法来解决其他常规时频变换方法存在的分辨率和能量聚集性较低的问题,改善时频分辨率和时频能量在时频平面中的聚集性,有效提高从时频域地震信号中提取地层结构的能力,可作为地震信号处理中可靠的时频变换技术。
背景技术
时频变换方法是直接应用地震资料进行烃类监测的重要手段之一,也是地震信号处理的重要步骤。时频变换方法旨在通过构造一种时间和频率的密度函数,将一个一维的时间信号以二维的函数形式表示出来,以揭示地震信号中所包含的频率分量,使我们不但能够同时了解到地震信号在时间-频率域的信息,而且可以精准地刻画地震信号随时间变化的特征。因此,时频变换方法对提高地震资料处理和解释的精度以及可靠性有着非常重要的意义。
现有技术中地震信号时频变换方法很多,例如:1)分数阶变换方法,是在短时傅里叶变换等基础上提出的一种对非平稳信号分析的方法,将原来的变换方法从时间-频率域推广到时间-分数阶域,其具有短时傅里叶变换的优点,增强了处理非平稳信号的灵活性。2)同步挤压变换,是在其它变换(如连续小波变换)基础上根据信号各组成分量的时频特征、压缩时频曲线,达到提高时频分辨率的方法。它属于一种时频重排算法,与原始的连续小波变换不同的是,它在提高时频分辨率的同时,也支持信号的重构,主要用于检测地震信号中的低频微弱信号。具体实现的方法有同步挤压小波变换。3)自适应稀疏时频分析方法,它首先通过交替分裂Bregman迭代算法自适应计算给定地震信号的最优高斯窗,使其在时频平面的能量紧凑性达到最大值;其次,利用生成的最优窗在凸约束下进行逆分解,使其混合范数的转置系数最小化,通过以上两个步骤达到提高时频分辨率的目的。
随着地震勘探的不断发展和完善,人们对时频分析的精度要求也大大提高。常规的时频分析方法往往只利用了通频带以内的地震信息,对截频和限频以外的信息很少利用,尤其是地震信号的振幅信息变化较大时,对于那些高、低频弱振幅信息的识别和利用更少。然而,由于传播路径的复杂性及不同时刻不同岩性的地层对地震波吸收、衰减程度的差异使得地震信号呈现明显的时变、非平稳等特征的这类地震信号,需要在时间-频率域同时精确的揭示信号的时频分布特征,这对地震信号的时频分析方法的精度和可靠性提出了更高的要求。
发明内容
本发明的目的是提供一种通过地震信号的自适应固有模态函数分解,并用优化时频变换求解每个固有模态函数的时频谱,通过叠加时频谱来提高时频分辨率的串联型地震信号优化时频变换方法。这种新方法不仅可以提高地震信号的时频分辨率,而且充分利用了原始地震道中的有效信息,同时具有很高的计算效率。
本发明所用的时频分析方法,摒弃了短时傅里叶变换、小波变换等时频变换方法的窗函数在整个时间信号上滑动,汲取了上述时频变换方法以及经验模态分解的优点并进行了改进和发展,可以精确地刻画地震信号中隐藏的高、低频弱振幅信息,其主要优点在于:构建了由信号分量主频自动控制最优参数的频域窗函数,频率覆盖范围广且时频分辨率高,同时具有很高的计算效率,适用于实际地震数据的处理。
本发明充分利用了地震勘探信号中的有效信息,可将其中的有效信息精确地展示在时频平面中,细致刻画了地震信号中不同频率信号分量随时间变化的特征,为后期地震资料的解释提供了有力的保障。
一种串联型地震信号优化时频变换方法,其具有以下优越性:
(1)本发明由于对输入的原始地震道进行完备集经验模态分解(CEEMD),使输入的原始地震道中各信号分量的频率特征逐渐趋于窄带化,利用寻优算法快速找到对应的最优频域窗函数,以便充分利用原始地震道中的高、低频弱振幅信息。
(2)本发明使用的最优窗函数,消除了其他时频变换的窗函数参数存在的人为主观影响因素,并利用构建的能量最优准则和优化时频分布计算公式,使时频平面中的能量聚集性达到最优,使地震信号的时频分辨率得到了明显提高。
(3)常规的时频变换方法是窗函数沿着整个时间信号滑动,而本发明使最优窗函数与输入信号在频率域快速匹配,以便于最优窗可以在很窄的有限带宽内快速滑动,因此大大减少了时频变换所需的时间,使计算效率得到数倍的提升。
本发明的具体实现原理如下:
首先输入三维地震勘探数据,从其中抽出某一道地震记录进行完备集经验模态分解。假设输入的原始地震道为x(t),根据经验模态分解原理找到其对应的第一个分量信号h1(t),并根据经验模态函数(IMF)的判定条件检验h1(t)是否为IMF。
重复上述操作N次,得到N个IMF。对N个IMF求均值,作为完备集经验模态分解的第一个IMF。重复上述操作依次得到第i个IMF,即:
式中,ωk是均值为零且方差为1的高斯白噪声,ε是固定系数,i表示进行完备集经验模态分解的次数。
以每个IMF分量为新的输入信号,此处,将其表示为imti(t),则输入信号imti(t)的主频fd如下计算:
式中,表示快速傅里叶变换,表示求取最大值的运算。
然后,依据主频fd来确定频域窗函数G(f,fd)的位置,因此新的输入信号imti(t)的主频与频域窗函数的主频同步移动,此时的频域窗函数G(f,fd)如下定义:
其中,η和γ为控制G(f,fd)的参数。
利用主频为fd的频域窗函数G(f,fd)和新输入信号imti(t)构建能量最优关系准则的目标函数,得
式中,IMFi(f)表示新输入信号imti(t)的傅里叶变换,即:
IMFi(f)=F[imti(t)]
利用上述能量最优关系准则目标函数,求取最优窗参数γopt如下:
然后,构建最优频域窗函数如下:
利用上式的最优频域窗函数,建立时频分辨率明显提高的优化时频分布计算公式,则输入信号imti(t)的优化时频分布TFRi(t,f)的计算公式如下:
式中,b=2fd,a=fd/5,b-a表示新输入信号imti(t)的带宽。由于该公式在很窄的有限带宽内计算,因此大大减少了时频变换的计算时间,使计算效率得到数倍的提升。
最后,可计算原始地震道x(t)的优化时频分布TFR(t,f),其表达式如下:
附图说明
图1是合成的Cross-chirp信号的时频谱对比图。其中,横坐标为时间,单位为秒,图1(a)纵坐标为振幅,图1(b)-1(d)纵坐标为频率,单位为赫兹。图1(b)是采用本发明的技术得到的时频谱,图1(c)是采用常规S变换得到的时频谱,图1(d)是采用小波变换得到的时频谱。
图2是某地区采集的一道地震数据的时频谱对比图。其中,横坐标为时间,单位为秒,图2(a)纵坐标为振幅,图2(b)-2(d)纵坐标为频率,单位为赫兹。图2(b)是采用本发明的技术得到的时频谱,图2(c)是采用常规S变换得到的时频谱,图2(d)是采用小波变换得到的时频谱。
图3是某地区采集的二维地震数据的瞬时谱剖面对比图,横坐标为道号,纵坐标为时间,单位为秒。其中,图3(b)是采用本发明的技术得到的25Hz的瞬时谱剖面图,图3(c)是采用常规S变换得到的25Hz瞬时谱剖面图。
具体实施方式
为了实现上述目标,本发明采取以下技术方案:一种串联型地震信号优化时频变换方法,其包括以下步骤:1)输入原始地震道x(t);2)利用完备集经验模态分解将原始地震道分解成多个不同频段的固有模态函数(IMF)窄带信号分量,并将每个固有模态函数窄带信号分量作为新的输入信号imti(t);3)求取每个新的输入信号imti(t)的主频fd,并依据能量最优准则目标函数求取每个频域窗函数的位置及其最优参数γopt;4)利用优化时频分布计算公式依次计算每个新输入信号imti(t)的时频分布谱TFRi(t,f);5)将4)得到的每个时频分布谱进行叠加,得到原始地震道x(t)的时频分布谱,从而得到了时频分辨率明显提高的地震道时频分布谱TFR(t,f)。
本发明的实施实例说明:
图1是合成的Cross-chirp信号的时频谱对比图。在图1(b)中可以明显看出,信号的频率分辨率在低频端和高频端都很高,能量聚集性均很好,表明时频分辨率高。图1(c)是常规S变换的时频谱,虽然较好地反映了信号的时频分布特征,但存在一些缺陷:在低频端,信号的能量聚集性和频率分辨率较好;在高频端,能量聚集性较差,且随着频率的增高,频率分辨率和收敛性变得更差。图1(d)是小波变换的时频谱,但是存在明显的缺陷:整个时频平面的能量聚集性比图1(a)差,且低频端的时频谱出现了畸变。
图2是某地区采集的一道地震数据的时频谱对比图。在图2(c)-(d)的时频谱中,70Hz以上的高频信息消失,另外,图2(d)的其时间-频率关系不能很好地对应,增加了后期解释的不确定性,但在图2(b)中,可以精确描述70-150Hz隐藏的高频信息,较好地反映了原始地震道中的高频信息。因此,图2(b)弥补了图2(c)-(d)的缺点,既挖掘了其他时频变换方法不能刻画的隐藏的高频信息,又显示了更高的时频分辨率。
图3是某地区采集的二维地震数据的瞬时谱剖面对比图。图3(b)是采用本发明得到的25Hz瞬时谱剖面,可以清晰地刻画目的地层的能量分布情况,整个瞬时谱剖面的分辨率较强,在0.60s、0.70s及0.95s处能量连续性较好且比较集中,能更好地对应能量变化与地层的关系。在图3(c)的常规S变换的25Hz瞬时谱剖面中,其能量连续性较差,只能看到0.60s及0.95s处的能量较连续,整个瞬时谱剖面的分辨率和能量聚集性较差,不能清晰地刻画25Hz的信号能量随地层的分布规律。通过比较可以看出:采用本发明得到的瞬时谱剖面,其能量的横向连续性、聚焦性以及分辨率等都得到了较大的提升。
上述各实施例仅用于说明本发明,其中方法的各实施步骤等都是可变化的,凡在本发明技术方法的基础上进行的类似变换和改进,均不排除在本发明的保护范围之外。

Claims (3)

1.一种串联型地震信号优化时频变换方法,其包括以下步骤:
1)输入原始地震道x(t);
2)利用完备集经验模态分解将原始地震道x(t)分解成多个不同频段的固有模态函数窄带信号分量,并将每个固有模态函数窄带信号分量作为新的输入信号imfi(t);
3)求取每个新的输入信号imti(t)的主频fd,并依据能量最优准则目标函数求取每个频域窗函数的位置及其最优参数γopt
4)利用优化时频分布计算公式依次计算每个新的输入信号imfi(t)的时频分布谱TFRi(t,f);
5)将4)得到的每个新的输入信号的时频分布谱叠加,得到原始地震道x(t)的时频分布谱,从而得到了时频分辨率明显提高的地震道时频分布谱TFR(t,f)。
2.根据权利要求1所述的一种串联型地震信号优化时频变换方法,其特征在于:所述步骤3)中利用主频为fd的频域窗函数G(f,fd)和新输入信号imti(t)构建了能量最优关系准则目标函数:
式中,IMFi(f)表示新的输入信号imti(t)的傅里叶变换,即:
IMFi(f)=F[imti(t)]
频域窗函数G(f,fd)如下定义:
其中,η和γ为控制G(f,fd)的参数。
利用上述能量最优关系准则目标函数,求取最优频域窗参数γopt如下:
3.根据权利要求1所述的一种串联型地震信号优化时频变换方法,其特征在于:所述步骤4)利用了能量最优准则目标函数得到的最优频域窗参数γopt,建立了优化时频分布计算公式:
式中,b=2fd,a=fd/5,b-a表示新的输入信号imti(t)的带宽;最优频域窗函数如下:
CN201611055163.5A 2016-11-25 2016-11-25 一种串联型地震信号优化时频变换方法 Active CN106556865B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201611055163.5A CN106556865B (zh) 2016-11-25 2016-11-25 一种串联型地震信号优化时频变换方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201611055163.5A CN106556865B (zh) 2016-11-25 2016-11-25 一种串联型地震信号优化时频变换方法

Publications (2)

Publication Number Publication Date
CN106556865A true CN106556865A (zh) 2017-04-05
CN106556865B CN106556865B (zh) 2018-09-07

Family

ID=58444940

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201611055163.5A Active CN106556865B (zh) 2016-11-25 2016-11-25 一种串联型地震信号优化时频变换方法

Country Status (1)

Country Link
CN (1) CN106556865B (zh)

Cited By (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107132576A (zh) * 2017-07-05 2017-09-05 西安交通大学 基于二阶挤压小波变换的地震资料时频分析方法
CN107505652A (zh) * 2017-07-26 2017-12-22 山东科技大学 一种基于能量分布特征的矿山微震信号辨识方法
CN107632320A (zh) * 2017-08-21 2018-01-26 西安交通大学 基于同步提取s变换的地震资料时频分析方法
CN108760294A (zh) * 2018-06-13 2018-11-06 苏州大学 旋转机械设备自动诊断系统及方法
CN110333071A (zh) * 2019-06-28 2019-10-15 华北电力大学 一种利用窄带倒谱变换的机械振动信号处理方法
CN110941013A (zh) * 2018-09-21 2020-03-31 中国石油化工股份有限公司 时间频率域能量聚焦方法和储层预测方法

Cited By (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107132576A (zh) * 2017-07-05 2017-09-05 西安交通大学 基于二阶挤压小波变换的地震资料时频分析方法
CN107132576B (zh) * 2017-07-05 2018-10-30 西安交通大学 基于二阶同步挤压小波变换的地震资料时频分析方法
CN107505652A (zh) * 2017-07-26 2017-12-22 山东科技大学 一种基于能量分布特征的矿山微震信号辨识方法
CN107505652B (zh) * 2017-07-26 2018-11-13 山东科技大学 一种基于能量分布特征的矿山微震信号辨识方法
WO2019019565A1 (zh) * 2017-07-26 2019-01-31 山东科技大学 一种基于能量分布特征的矿山微震信号辨识方法
CN107632320A (zh) * 2017-08-21 2018-01-26 西安交通大学 基于同步提取s变换的地震资料时频分析方法
CN108760294A (zh) * 2018-06-13 2018-11-06 苏州大学 旋转机械设备自动诊断系统及方法
CN110941013A (zh) * 2018-09-21 2020-03-31 中国石油化工股份有限公司 时间频率域能量聚焦方法和储层预测方法
CN110333071A (zh) * 2019-06-28 2019-10-15 华北电力大学 一种利用窄带倒谱变换的机械振动信号处理方法
CN110333071B (zh) * 2019-06-28 2021-09-10 华北电力大学 一种利用窄带倒谱变换的机械振动信号处理方法

Also Published As

Publication number Publication date
CN106556865B (zh) 2018-09-07

Similar Documents

Publication Publication Date Title
CN106556865A (zh) 一种串联型地震信号优化时频变换方法
CN106597532B (zh) 一种结合井资料与层位资料的叠前地震数据频带拓展方法
CN107272062B (zh) 一种数据驱动的地下介质q场估计方法
CN106842298A (zh) 一种基于匹配追踪的不整合强反射自适应分离方法
CN110187388B (zh) 一种基于变分模态分解的稳定地震品质因子q估计方法
Sun et al. Cross-correlation analysis and time delay estimation of a homologous micro-seismic signal based on the Hilbert–Huang transform
CN103728660A (zh) 基于地震数据的多道匹配追踪方法
CN106707334B (zh) 一种提高地震资料分辨率的方法
CN108020863A (zh) 一种基于地震奇偶函数的碳酸盐岩薄储层孔隙度预测方法
CN107132579A (zh) 一种保地层结构的地震波衰减补偿方法
CN107703546A (zh) 一种基于小波变换的新阈值函数地震资料去噪方法
CN104730576A (zh) 基于Curvelet变换的地震信号去噪方法
CN107918146B (zh) 一种基于非线性挤压s时频变换的弱信号检测方法
CN106483563A (zh) 基于互补集合经验模态分解的地震能量补偿方法
CN104391324A (zh) 依赖频率的avo反演前的地震道集动校拉伸校正预处理技术
CN105334532A (zh) 一种地震子波估计方法
CN102854530B (zh) 基于对数时频域双曲平滑的动态反褶积方法
CN113281809B (zh) 一种地震信号的谱分析方法
Wu‐Yang et al. Research and application of improved high precision matching pursuit method
CN107728213A (zh) 一种小波新阈值函数地震资料去噪方法
CN108957529B (zh) 一种基于属性的无井区子波估计方法
Li et al. Suppression of strong random noise in seismic data by using time-frequency peak filtering
Xin et al. Seismic high-resolution processing method based on spectral simulation and total variation regularization constraints
CN110673211B (zh) 一种基于测井与地震数据的品质因子建模方法
Zhao et al. Automatic events extraction in pre-stack seismic data based on edge detection in slant-stacked peak amplitude profiles

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
GR01 Patent grant
GR01 Patent grant