CN110579800B - 一种基于高精度同步挤压变换的地震数据数字处理方法 - Google Patents

一种基于高精度同步挤压变换的地震数据数字处理方法 Download PDF

Info

Publication number
CN110579800B
CN110579800B CN201910996601.5A CN201910996601A CN110579800B CN 110579800 B CN110579800 B CN 110579800B CN 201910996601 A CN201910996601 A CN 201910996601A CN 110579800 B CN110579800 B CN 110579800B
Authority
CN
China
Prior art keywords
frequency
sampling interval
time
time sampling
seismic data
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.)
Active
Application number
CN201910996601.5A
Other languages
English (en)
Other versions
CN110579800A (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.)
Southwest Petroleum University
Original Assignee
Southwest Petroleum University
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 Southwest Petroleum University filed Critical Southwest Petroleum University
Priority to CN201910996601.5A priority Critical patent/CN110579800B/zh
Publication of CN110579800A publication Critical patent/CN110579800A/zh
Application granted granted Critical
Publication of CN110579800B publication Critical patent/CN110579800B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

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. analysis, for interpretation, for correction
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/40—Transforming data representation
    • G01V2210/48—Other transforms

Abstract

本发明公开了一种基于高精度同步挤压变换的地震数据数字处理方法，包括以下步骤：S1.输入待处理的所有地震数据；S2.选择一道地震数据，进行傅里叶变换，根据其最大动态范围，找到其对应的有效频带的截止频率；S3.根据有效频带的截止频率反演时间采样间隔；S4.利用反演的时间采样间隔对该道数据进行重采样；S5.对重采样的数据进行高精度同步挤压变换；S6.对于每一道地震数据，重复步骤S2～S5直至所有地震数据处理完成。本发明根据有效频带的截止频率计算出适合的时间采样间隔以保证同步挤压变换中瞬时频率计算的精度，然后对输入数据进行重采样，再进行同步挤压变换即可得到高精度的同步挤压变换的时频谱。

Description

一种基于高精度同步挤压变换的地震数据数字处理方法
技术领域
本发明属于地震资料数字处理领域，特别是涉及一种基于高精度同步挤压变换的地震数据数字处理方法。
背景技术
同步挤压变换具有很好的时频分辨率，已经被广泛应用于数字信号处理的各个方面。同步挤压变换能够将非线性和非平稳信号分解为一系列的固有模态函数，其中，固有模态函数可被认为是一系列有着准确数学定义的近似简谐波成分的组合。同步挤压变换与短时傅里叶变换，S变换，广义S变换，小波变换，改进的短时傅里叶变换等时频分析方法相比，同步挤压变换能够沿着频率方向对时频谱进行挤压实现时频谱重排，使得时频能量聚焦到信号的真实瞬时频率上，提高了时频谱的时频分辨率。因此，同步挤压变换更适合非线性和非平稳信号的分析和处理。
但是，同步挤压变换中瞬时频率的计算精度随着时间采样间隔和实际频率的增大而减小，导致同步挤压变换的时频谱的精度降低，这将影响后续的地震资料处理与解释。
发明内容
本发明的目的在于克服现有技术的不足，提供一种基于高精度同步挤压变换的地震数据数字处理方法，根据有效频带的截止频率计算出适合的时间采样间隔以保证同步挤压变换中瞬时频率计算的精度，然后对输入数据进行重采样，再进行同步挤压变换即可得到高精度的同步挤压变换的时频谱。
本发明的目的是通过以下技术方案来实现的：一种基于高精度同步挤压变换的地震数据数字处理方法，包括以下步骤：
S1.输入待处理的所有地震数据；
S2.选择一道地震数据，进行傅里叶变换，根据其最大动态范围，找到其对应的有效频带的截止频率；
S3.根据有效频带的截止频率反演时间采样间隔；
S4.利用反演的时间采样间隔对该道数据进行重采样；
S5.对重采样的数据进行高精度同步挤压变换；
S6.对于每一道地震数据，重复步骤S2～S5直至所有地震数据处理完成。
进一步地，所述步骤S2包括：
S201.对输入的单道地震数据s(t)以时间采样间隔Δt进行傅里叶变换，并提取傅里叶变换得到的复数信号X(f)+j·Y(f)的振幅谱G(f)和峰值振幅G(fp)；
Figure BDA0002239918530000021
其中，j是虚数单位，Max[ ]表示取最大值，fp为峰值振幅G(fp)对应的峰值频率；
S202.计算地震资料有效频带的截止频率fd；
Figure BDA0002239918530000022
其中，d是G(f)的最大动态范围，单位分贝；L(f)表示所有满足
Figure BDA0002239918530000023
的频率成分。
进一步地，所述步骤S3中包括：
S301.根据有效频带的截止频率反演得到的时间采样间隔Δτ1：
Figure BDA0002239918530000024
其中,K为瞬时频率的计算精度，根据用户需求设定；
根据反演得到的时间采样间隔Δτ1和原始的时间采样间隔Δt，确定最终反演的时间采样间隔Δτ：
Figure BDA0002239918530000025
S302.若根据有效频带的截止频率反演得到的时间采样间隔Δτ1大于原始的时间采样间隔Δt，说明原始的时间采样间隔Δt已经能够满足处理过程的精度要求，因此，最终反演的时间采样间隔Δτ取原始的时间采样间隔Δt，反之，则取Δτ1。
进一步地，所述步骤S4包括：
对输入信号s(t)以时间采样间隔Δt进行傅里叶变换得到的复数信号，并在频率域对输入信号s(t)以最终反演的时间采样间隔Δτ进行重采样，得到S(f)：
Figure BDA0002239918530000026
-[2Δτ]-1≤f＜-[2Δt]-1为反演的时间采样间隔Δτ对应的最大负频率与原始时间采样间隔Δt对应的最大负频率之间的负频率成分，-[2Δt]-1≤f＜[2Δt]-1为原始时间采样间隔Δt对应的频率部分，[2Δt]-1≤f＜[2Δτ]-1为原始时间采样间隔Δt对应的最大频率与反演的时间采样间隔Δτ对应的最大频率间的其他频率成分。
现在将时间采样间隔变为Δτ，必然会引入其他频率成分，为了只保留原始数据包含的-[2Δt]-1≤f＜[2Δt]-1对应的频率成分的信息，让除-[2Δt]-1≤f＜[2Δt]-1对应的频率成分之外的其他频率成分的复信号等于0。
对重采样以后的频域信号S(f)进行反傅里叶变换即可得到重采样以后的时域信号z(τ)，τ是时间变量，z(τ)的时间采样间隔为根据有效频带的截止频率反演的时间采样间隔Δτ，且z(τ)包含的采样点数为输入信号s(t)包含的采样点数的
Figure BDA0002239918530000031
倍。
进一步地，所述步骤S5包括：
S501.把重采样以后的数据z(τ)作为输入数据，进行时频分析，得到z(τ)的时频谱H(τ,f)；
S502.计算时频谱H(τ,f)的瞬时频率p(t,f)；
Figure BDA0002239918530000032
其中，hr(τ,f)为时频谱H(τ,f)的实部，hi(τ,f)为时频谱H(τ,f)的虚部；
S503.根据时频谱H(τ,f)和瞬时频率p(t,f)计算输入信号s(t)的高精度同步挤压变换N(t,ω)；
N(t,ω)＝∫H(t,f)df,if p(t,f)∈[ω,ω+Δω]
其中，ω表示频率，Δω代表频率增量。
本发明的有益效果是：本发明根据有效频带的截止频率和要求的瞬时频率的计算精度，可以反演出更为适合的时间采样间隔；用反演的时间采样间隔对输入数据进行重采样，再进行高精度同步挤压变换。常规的同步挤压变换的瞬时频率的计算精度容易受到时间采样间隔和瞬时频率的大小的影响。本发明提高了同步挤压变换过程中的瞬时频率的计算精度，获得的同步挤压变换的时频谱具有很好的精度和很高的时频分辨率，这对于后续地震资料的处理及解释具有重要意义。
附图说明
图1为本发明的方法流程图；
图2为实施例中，改进的短时傅里叶变换的时频谱、不同时间采样间隔的同步挤压变换结果以及高精度同步挤压变换的结果；
图3为实施例中合成信号的示意图；
图4为实施例中合成信号的时频谱示意图；
图5为实施例中改进的短时傅里叶变换效果示意图；
图6为实施例中不同属性的时间切片示意图；
图7为实施例中短时傅里叶变换的时频谱中不同谐波成分的时间切片示意图；
图8为实施例中高精度同步挤压变换的时频谱中不同谐波成分的时间切片。
具体实施方式
下面结合附图进一步详细描述本发明的技术方案，但本发明的保护范围不局限于以下所述。
考虑到输入数据的频率越高，要求时间采样间隔越小，若时间采样间隔太小，将导致计算量的增加；若时间采样间隔太大，将导致计算误差增加。本发明先通过傅里叶变换获得有效频带的截止频率，然后根据有效频带的截止频率和要求的瞬时频率的计算精度反演出时间采样间隔，然后根据反演的时间采样间隔对输入数据进行重采样，再进行同步挤压变换，即可获得高精度的同步挤压变换的时频谱，具体地：
如图1所示，一种基于高精度同步挤压变换的地震数据数字处理方法，包括以下步骤：
S1.输入待处理的所有地震数据；
S2.选择一道地震数据，进行傅里叶变换，根据其最大动态范围，找到其对应的有效频带的截止频率；
S3.根据有效频带的截止频率反演时间采样间隔；
S4.利用反演的时间采样间隔对该道数据进行重采样；
S5.对重采样的数据进行高精度同步挤压变换；
S6.对于每一道地震数据，重复步骤S2～S5直至所有地震数据处理完成。
其中，本文案例中得到时频谱H(t,f)使用的时频分析方法是改进的短时傅里叶变换，可以表示为：
Figure BDA0002239918530000041
其中，τ为时移因子，sa(t)为输入信号s(t)的解析信号，w(t)为窗函数；
在本申请的实施例中，图2为改进的短时傅里叶变换的时频谱、不同时间采样间隔的同步挤压变换结果以及高精度同步挤压变换的结果，在该实施例中，同步挤压变换是在包括改进的短时傅里叶变换在内的时频分析方法产生的时频谱的基础上进行的，把时频谱进行挤压重排，以此达到提高分辨率的目的，这里给出改进的短时傅里叶变换的时频谱的目的是由于该方法产生的时频谱是输入信号的一种表征方式，基本上对应着输入信号的频率，同步挤压变换的结果与输入信号的频率相比存在误差，因此用改进的短时傅里叶变换的时频谱为同步挤压变换的结果存在多大的误差提供一种参考，与高精度同步挤压变换相比，来体现高精度同步挤压变换的效果；具体地，其中图2(a)为改进的短时傅里叶变换的时频谱，图2(b)、图2(c)、分别为1ms、2ms时间采样间隔下的同步挤压变换示意图，可以看出，同步挤压变换的结果与输入信号的频率存在误差，且时间采样间隔越大，误差越大；图2(d)为高精度同步挤压变换示意图，高精度同步挤压变换的结果基本达到了输入信号的频率，提高了同步挤压变换的精度；
在本申请的实施例中，考虑到自然界中存在的大多数信号都可以认为是非平稳信号，故在测试、比较高精度同步挤压变换和同步挤压变换结果时合成一种具有普遍性和代表性的非平稳信号，合成信号如图3所示，合成信号频率随着时间的变化而变化；图4为图3所示的合成信号的时频谱，图4可以看出，同步挤压变换的结果与输入信号的频率存在误差，特别是高频成分；高精度同步挤压变换的结果基本与输入信号的频率相一致，提高了同步挤压变换时频谱的精度；图5改进的短时傅里叶变换效果示意图，可以看出在改进的短时傅里叶变换产生的时频谱中，弱能量成分被强能量成分所掩盖；在高精度同步挤压变换的时频谱上，这种现象得到有效改善，有效的凸显了弱能量成分，且时频能量更聚焦，提高了同步挤压变换的时频谱的分辨率和精度；图6为不同属性的时间切片示意图；图7为在短时傅里叶变换的时频谱中不同谐波成分的时间切片示意图；图8为在高精度同步挤压变换的时频谱中不同谐波成分的时间切片；从图6～8可以看出，图6、7、8主要是表现箭头所示的河道，在高精度同步挤压变换的时频谱上，河道清晰可见，分辨率有显著的提升。
经过以上处理过程以后，得到了高精度的同步挤压变换的时频谱，和常规同步挤压变换的时频谱相比，该方法得到的时频谱具有更高的精度，能够准确的展现输入数据的时频关系和精准的刻画了该区域的构造特征，对于后续地震资料的处理及解释具有重要意义。
以上所述是本发明的优选实施方式，应当理解本发明并非局限于本文所披露的形式，不应该看作是对其他实施例的排除，而可用于其他组合、修改和环境，并能够在本文所述构想范围内，通过上述教导或相关领域的技术或知识进行改动。而本领域人员所进行的改动和变化不脱离本发明的精神和范围，则都应在本发明所附权利要求的保护范围内。

Claims (5)

1.一种基于高精度同步挤压变换的地震数据数字处理方法，其特征在于：包括以下步骤：
S1.输入待处理的所有地震数据；
S2.选择一道地震数据，进行傅里叶变换，根据其最大动态范围，找到其对应的有效频带的截止频率；
S3.根据有效频带的截止频率反演时间采样间隔；
S4.利用反演的时间采样间隔对该道数据进行重采样；
S5.对重采样的数据进行高精度同步挤压变换；
S6.对于每一道地震数据，重复步骤S2～S5直至所有地震数据处理完成。
2.根据权利要求1所述的一种基于高精度同步挤压变换的地震数据数字处理方法，其特征在于：所述步骤S2包括：
S201.对输入的单道地震数据s(t)以时间采样间隔Δt进行傅里叶变换，并提取傅里叶变换得到的复数信号X(f)+j·Y(f)的振幅谱G(f)和峰值振幅G(fp)；
Figure FDA0002804261370000011
其中，j是虚数单位，Max[ ]表示取最大值，fp为峰值振幅G(fp)对应的峰值频率；
S202.计算地震资料有效频带的截止频率fd；
Figure FDA0002804261370000012
其中，d是G(f)的最大动态范围，单位分贝；L(f)表示所有满足
Figure FDA0002804261370000013
的频率成分。
3.根据权利要求2所述的一种基于高精度同步挤压变换的地震数据数字处理方法，其特征在于：所述步骤S3中包括：
S301.根据有效频带的截止频率反演得到的时间采样间隔Δτ1：
Figure FDA0002804261370000014
其中,K为瞬时频率的计算精度，根据用户需求设定；
根据反演得到的时间采样间隔Δτ1和原始的时间采样间隔Δt，确定最终反演的时间采样间隔Δτ：
Figure FDA0002804261370000021
S302.若根据有效频带的截止频率反演得到的时间采样间隔Δτ1大于原始的时间采样间隔Δt，说明原始的时间采样间隔Δt已经能够满足处理过程的精度要求，因此，最终反演的时间采样间隔Δτ取原始的时间采样间隔Δt，反之，则取Δτ1。
4.根据权利要求2所述的一种基于高精度同步挤压变换的地震数据数字处理方法，其特征在于：所述步骤S4包括：
对输入信号s(t)以时间采样间隔Δt进行傅里叶变换得到的复数信号，并在频率域对输入信号s(t)以最终反演的时间采样间隔Δτ进行重采样，得到S(f)：
Figure FDA0002804261370000022
对重采样以后的频域信号S(f)进行反傅里叶变换即得到重采样以后的时域信号z(τ)，τ是时间变量，z(τ)的时间采样间隔为根据有效频带的截止频率反演的时间采样间隔Δτ，且z(τ)包含的采样点数为输入信号s(t)包含的采样点数的
Figure FDA0002804261370000023
倍。
5.根据权利要求3所述的一种基于高精度同步挤压变换的地震数据数字处理方法，其特征在于：所述步骤S5包括：
S501.把重采样以后的数据z(τ)作为输入数据，进行时频分析，得到z(τ)的时频谱H(τ,f)；
S502.计算时频谱H(τ,f)的瞬时频率p(t,f)；
Figure FDA0002804261370000024
其中，hr(τ,f)为时频谱H(τ,f)的实部，hi(τ,f)为时频谱H(τ,f)的虚部；
S503.根据时频谱H(τ,f)和瞬时频率p(t,f)计算输入信号s(t)的高精度同步挤压变换N(t,ω)；
N(t,ω)＝∫H(t,f)df,if p(t,f)∈[ω,ω+Δω]
其中，ω表示频率，Δω代表频率增量。
CN201910996601.5A 2019-10-19 2019-10-19 一种基于高精度同步挤压变换的地震数据数字处理方法 Active CN110579800B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910996601.5A CN110579800B (zh) 2019-10-19 2019-10-19 一种基于高精度同步挤压变换的地震数据数字处理方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910996601.5A CN110579800B (zh) 2019-10-19 2019-10-19 一种基于高精度同步挤压变换的地震数据数字处理方法

Publications (2)

Publication Number Publication Date
CN110579800A CN110579800A (zh) 2019-12-17
CN110579800B true CN110579800B (zh) 2021-06-01

Family

ID=68815053

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910996601.5A Active CN110579800B (zh) 2019-10-19 2019-10-19 一种基于高精度同步挤压变换的地震数据数字处理方法

Country Status (1)

Country Link
CN (1) CN110579800B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111623986A (zh) * 2020-05-19 2020-09-04 安徽智寰科技有限公司 基于同步压缩变换与时频匹配的信号特征提取方法及系统

Family Cites Families (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7826968B2 (en) * 2007-04-30 2010-11-02 Mediatek Inc. Method, device and vehicle utilizing the same
CN103048684B (zh) * 2011-10-11 2015-10-28 中国石油化工股份有限公司 一种多分量地震资料面波压制方法
WO2014191011A1 (en) * 2013-05-27 2014-12-04 Statoil Petroleum As High resolution estimation of attenuation from vertical seismic profiles
CN105467451B (zh) * 2016-01-13 2018-05-15 中国石油集团东方地球物理勘探有限责任公司 基于全变差最小化约束的地震反射系数反演方法
CN109239780A (zh) * 2017-07-10 2019-01-18 中国石油化工股份有限公司 基于同步挤压小波变换去除面波的方法
CN107390267B (zh) * 2017-07-27 2019-03-01 西安交通大学 一种同步挤压变换域的地震资料衰减补偿方法
CN107576981B (zh) * 2017-08-31 2019-02-12 大连理工大学 一种基于监测位移和截止频率的层间位移修正方法

Also Published As

Publication number Publication date
CN110579800A (zh) 2019-12-17

Similar Documents

Publication Publication Date Title
Jenkins General considerations in the analysis of spectra
CN110687595B (zh) 一种基于时间重采样和同步挤压变换的地震数据处理方法
CN110907827B (zh) 一种马达瞬态失真测量方法及系统
CN110579800B (zh) 一种基于高精度同步挤压变换的地震数据数字处理方法
JPH0619390B2 (ja) デイジタル・フ−リエ変換の後処理方法
Zhang et al. Compound fault extraction method via self-adaptively determining the number of decomposition layers of the variational mode decomposition
CN106324342A (zh) 一种基于查表的谐波检测方法
CN108801296B (zh) 基于误差模型迭代补偿的传感器频响函数计算方法
CN116046968A (zh) 一种液相色谱工作站数据处理方法、系统及可存储介质
CN111273345B (zh) 一种基于高精度时频瞬时相位的地震数据时频谱处理方法
CN112505413B (zh) 一种时频分析方法和系统
CN112067927B (zh) 中高频振荡检测方法及装置
CN112702687A (zh) 一种快速确认喇叭或者整机失真的方法
CN106980043A (zh) 一种基于汉宁窗的改进相位差校正法
Wolf et al. Amplitude and frequency estimator for aperiodic multi-frequency noisy vibration signals of a tram gearbox
Junhong et al. Bearing fault diagnosis based on improved LMD
CN117169590B (zh) 一种基于软件变采样率的电力谐波分析的方法和装置
CN115691537B (zh) 一种耳机音频信号的分析与处理系统
CN109709397B (zh) 一种加连续Hanning窗的电网谐波非同步压缩感知检测方法
CN108120498B (zh) 用于变压器声学在线监测的双通道输入信号合成方法
CN114088077B (zh) 一种改进的半球谐振陀螺信号去噪方法
CN114142853B (zh) 一种基于插值dft信号同步的数字锁相放大处理方法
WO2022000656A1 (zh) 马达系统失真的测量方法及设备、计算机可读存储介质
CN115146221A (zh) 基于多重同步压缩的mstct时频分析方法及系统
Yang et al. Frequency-domain correction of sensor dynamic error for step response

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