CN110398773B - 一种针对部分缺失的地震数据的恢复重构方法 - Google Patents

一种针对部分缺失的地震数据的恢复重构方法 Download PDF

Info

Publication number
CN110398773B
CN110398773B CN201910474117.6A CN201910474117A CN110398773B CN 110398773 B CN110398773 B CN 110398773B CN 201910474117 A CN201910474117 A CN 201910474117A CN 110398773 B CN110398773 B CN 110398773B
Authority
CN
China
Prior art keywords
data
shot set
shot
seismic data
mapping
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
CN201910474117.6A
Other languages
English (en)
Other versions
CN110398773A (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.)
China University of Petroleum East China
Original Assignee
China University of Petroleum East China
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 East China filed Critical China University of Petroleum East China
Priority to CN201910474117.6A priority Critical patent/CN110398773B/zh
Publication of CN110398773A publication Critical patent/CN110398773A/zh
Application granted granted Critical
Publication of CN110398773B publication Critical patent/CN110398773B/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
    • G01V1/282—Application of seismic models, synthetic seismograms
    • 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
    • G01V1/30—Analysis
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations

Abstract

本发明公开了一种针对部分缺失的地震数据的恢复重构方法。本发明对地震数据采用坐标变换的方法，变换到映射域，减少直达的斜率，提高直达波的拉平度，提高了其在曲波域的凸集映射数据恢复精度，将恢复的数据再利用反映射公式变换到正常炮集。此外，利用先验约束信息约束进行进一步高精度恢复，从而进一步提高了数据恢复的信噪比和恢复精度。

Description

一种针对部分缺失的地震数据的恢复重构方法
技术领域
本发明涉及一种针对部分缺失的地震数据的恢复重构方法，属于油气物探工程领域。
背景技术
在当前的勘探领域中，地震数据采集的完整性关键，由于采集成本、坏道、噪声、地形等诸多因素的影响，实际采集到的数据往往会出现缺失道和野值道，这些不规则的或者稀疏的地震数据，给数据的处理与解释工作造成了一定的困扰。曲波域凸集映射方法可以较好的重建炮集，但是炮集中直达波的能量较强时，在迭代重建地震数据的过程中产生较强的干扰噪音，影响重建效果，可以在重建过程中将直达波进行切除处理，再进行数据恢复，然而有些直达波与反射波等波场耦合在一起，并且直达波及早至波在后期处理中，尤其在速度反演中具有重要作用。
因此，对于地震勘探中采集的地震数据缺失部分，进行高精度地震数据恢复是非常必要的。
发明内容
本发明的目的在于，提供一种针对部分缺失的地震数据的恢复重构方法，它可以解决当前技术中存在的问题，对含有较强直达波不完整地震数据，提高了地震数据恢复精度。
为解决上述技术问题，本发明采用如下的技术方案，一种针对部分缺失的地震数据的恢复重构方法，包括以下步骤：
S1，输入地震的观测数据D，并且拾取初至窗；
S2，利用初至窗，对观测数据利用坐标映射进行变换得到映射炮集，具体包括如下方式：
p2(x,z)＝p(x,z1+1)*(t2(x,z)-t(x,z1))/dt
+p(x,z1)*(t(x,z1+1)-t2(x,z))/dt
其中，(x,z)表示观测数据第x道第z个采样点，p(x,z)为观测数据D在(x,z)采样点处的波场值，p2(x,z)为映射炮集D2在(x,z)采样点处的波场值，t(x,z1)为观测数据在采样点(x,z1)处的采样时刻，t2(x,z)为映射炮集在(x,z)处的采样时刻，其中t2(x,z)、t(x,z1)、z1分别通过以下公式获得：
Figure BDA0002081663300000021
其中，gc(x)为第x道初至时间窗G(x)，dt为观测数据的采样时间间隔，dtxin为映射炮集的采样时间间隔，nt为炮集的每道的采样点数，int为取整函数；
S3，对映射炮集做曲波变换得到曲波域内的曲波系数矩阵；
S4，基于所述曲波系数矩阵，利用凸集映射方法对映射炮集进行迭代重构得到重构炮集；
S5，利用反坐标映射公式对重构炮集反变换得到恢复炮集，具体包括：
Figure BDA0002081663300000022
p1(x,z)为采样点处(x,z)反变换后的波场值，组成恢复炮集Dpre，mod为取余函数。
与现有技术相比，本发明对地震炮集采用坐标变换的方法，变换到映射域，减少直达的斜率，提高直达波的拉平度，提高了其在曲波域的凸集映射数据恢复精度，将恢复的数据再利用反映射公式变换到正常炮集。此外，利用先验约束信息约束进行进一步高精度恢复，从而进一步提高了数据恢复的信噪比和恢复精度。
附图说明
附图是用来提供对本发明的进一步理解，并且构成说明书的一部分，与下面的具体实施方式一起用于解释本发明，但并不构成对本发明的限制。在附图中：
图1为本发明的实施例的流程示意图；
图2为本常规技术中的完整地震数据；
图3为本说明书实施例所提供的缺失部分数据的观测数据；
图4为本发明所提供的映射炮集的示意图；
图5为本发明所提供的多次迭代得到的最终的恢复炮集的示意图；
图6为本发明所提供的对重构炮集进行反坐标变换的示意图；
图7为本发明所提供的基于约束所得到的目标地震数据的示意图；
图8为不采用重构约束的常规凸集映射地震数据的示意图；
图9a为常规方法所得到的地震数据与标准的数据做差所得到的示意图，图9b为本发明所得到的地震数据与标准的数据做差所得到的示意图。
下面结合附图和具体实施方式对本发明作进一步的说明。
具体实施方式
本发明的实施例1：一种提高缺失地震勘探数据恢复精度重构方法，如图1所示，包括以下步骤：
S1，输入地震的观测数据，并且拾取初至窗。
本发明可以针对不同观测数据进行不同选取，当炮点在中间时，选取炮点位置与两边的最远接收道的地震波场初至点组成的两条直线作为初至窗；当炮点在两端时，选取炮点位置与最远接收道的地震波场初至点一条直线作为初至窗。
理想的状态下，完成的地震数据常为以共炮方式存储地震数据表现为单炮地震数据，如图2所示。
在地震数据缺失时，则常常如图3所示。图3即为本说明书实施例所提供的缺失部分数据的观测数据。在测试时，可以基于图2中所示的数据进行一定比例的数据抽取(例如抽取70％)而得到
对图3选取炮点(0，1400)与(230，0)直线为变换G1，(230，0)与(461，1200)直线为变换G2，G1与G2构成映射变换的初至窗G，如图3中的直线G1，G2，五星为炮点，三角形为最远接收道的波场初至点。
S2，利用初至窗，对观测数据利用坐标映射进行变换得到映射炮集。
根据S1中拾取的初至窗G，对缺失的地震数据做映射变换，变换公式为：
Figure BDA0002081663300000041
其中，(x,z)表示观测数据第x道第z个采样点，p(x,z)为观测数据D在(x,z)采样点处的波场值，p2(x,z)为映射炮集Dy在(x,z)采样点处的波场值，t(x,z1)为观测数据在采样点(x,z1)处的采样时刻，t2(x,z)为变换后的映射炮集再(x,z)处的采样时刻。其中t2(x,z)、t(x,z1)、z1分别通过以下公式获得：
Figure BDA0002081663300000042
其中，gc(x)为第x道初至时间窗G(x)，在实际资料处理中可以近似的采用从炮点到检波点斜的直线，dt为观测数据的采样时间间隔，dtxin为映射炮集的采样时间间隔，nt为每道的采样点数，int为取整函数，映射炮集如图4所示。
S3，对映射炮集做曲波变换得到曲波域内的曲波系数矩阵。
曲线是稀疏表示地震波场的理想字典，依据反射系数在曲波域更为稀疏的特点，利用曲波变换实现凸集映射方法。考虑二维情况下的曲波变换，曲波系数矩阵可以表示为：
Figure BDA0002081663300000051
式中φm表示曲波基函数,m＝(j,kd,l)，其中j,kd,l分别为尺度、方向、位置参数，cm为曲波域内的系数矩阵，t1,t2为数据Dy的时间窗。
S4，基于所述曲波系数矩阵，利用凸集映射方法对映射炮集进行迭代重构得到重构炮集。
基于公式(3)，映射炮集Dy可以用曲波系数矩阵和基函数表示为：
Figure BDA0002081663300000052
具体而言，在曲波域地震的数据凸集映射插值方法的公式如下所示，第k次迭代的映射炮集在曲波域的硬阈值约束为：
Figure BDA0002081663300000053
本发明中阈值函数矩阵取值Sτ采用指数衰减的形式阈值函数:
τk＝τ0*e(1-k)/M,k＝1,2...num(6)
其中M为控制衰减快慢的系数。给定初始的最大阈值τ0，系数衰减项M和迭代次数k即可确定对应的阈值。大于阈值τk的系数保留，小于阈值的系数变为0，num为总迭代次数，e为自然常数。
S5，利用反坐标映射公式对重构炮集反变换得到恢复炮集。
对重构炮集在映射域重构炮集反坐标变换，具体公式如下：
Figure BDA0002081663300000054
其中，mod为取余函数，p1(x,z)为采样点处(x,z)反变换后的波场值，组成恢复炮集Dpre，如图6所示，此时的恢复炮集已经恢复了相当部分的缺失数据。
在一种实施例中，通过不断的迭代，阈值从大到小的选取，则缺失的地震数据可以被重构出来：
Figure BDA0002081663300000055
其中，Dy,obs为S2中映射地震数据Dy未缺失的部分，I为单位矩阵，
Figure BDA0002081663300000061
为第k次重构的重构炮集，P表示地震数据缺失矩阵，
Figure BDA0002081663300000062
为阈值函数矩阵取值，Φ,ΦT分别表示地震数据的曲波系数以及曲波域系数求取函数，重复S3到S5满足迭代次数，输出最终的恢复炮集，如图5所示。
进一步地，还可以对恢复炮集作为先验约束信息约束，迭代生成目标地震数据。
本发明中采用常规的凸集映射方法进行改进，迭代重构得到目标炮集，并且采用新的约束公式为：
Figure BDA0002081663300000063
其中
Figure BDA0002081663300000064
为第k+1次迭代的目标地震数据，其中Dk+1为采用常规凸集映射公式重构的地震数据，Lmk+1为加权矩阵因子，其值取决于对先验信息的准确性。本发明中对于窗函数G以上的部分取值为1，其余部分从上往下选取0.8到0.1递减的过程。
具体而言，Dk+1可以通过以下步骤获得，对第k次目标地震数据
Figure BDA0002081663300000065
做曲波变换得到曲波域内的曲波系数矩阵。当k＝1时，目标地震数据即为观测数据。
考虑二维情况下的曲波变换，对目标地震数据
Figure BDA0002081663300000066
可以用曲波系数矩阵和基函数表示为：
Figure BDA0002081663300000067
则地震数据Dk+1可以被重构出来：
Figure BDA0002081663300000068
其中，
Figure BDA0002081663300000069
为第k次迭代的得到目标地震数据，Dobs为观测炮集中未缺失的道集部分，I为单位矩阵，为第k次重构的重构炮集，P表示地震数据缺失矩阵，
Figure BDA00020816633000000610
为阈值函数矩阵取值，Φ,ΦT分别表示地震数据的曲波系数以及曲波域系数求取函数，如公式(6)所示。
当满足预设的迭代次数时，输出目标地震数据Dk+1，即随着迭代次数的增加，重构的误差逐渐减少，直到满足重构要求的地震数据，如图7所示。不采用重构约束的常规凸集映射地震数据如图8所示，通过两者比较可以发现图7对于直达波区域的重构更加精确。
将两者和标准地震数据作差比较，如图9所示，图9a为常规方法所得到的地震数据与标准的数据做差所得到的示意图，图9b为本发明所得到的地震数据与标准的数据做差所得到的示意图。
可以看出采用本说明书实施例所挺的方法所得到重构的地震数据误差更小。由于本发明在测试选用的缺失较大，缺失70％，因此常规方法数据存在较大误差，当缺失的地震数据的较少时，常规方法也可以取得较好的恢复效果，但是经过测试本发明方法误差也比常规方法恢复误差更小。
通过前述的方法，本发明对地震数据采用坐标变换的方法，变换到映射域，减少直达的斜率，提高直达波的拉平度，提高了其在曲波域的凸集映射数据恢复精度，将恢复的数据再利用反映射公式变换到正常炮集。此外，利用先验约束信息约束进行进一步高精度恢复，从而进一步提高了数据恢复的信噪比和恢复精度。

Claims (3)

1.一种针对部分缺失的地震数据的恢复重构方法，其特征在于，包括以下步骤：
S1，输入地震的观测数据D，并且拾取初至窗；
S2，利用初至窗，对观测数据利用坐标映射进行变换得到映射炮集，具体包括如下方式：
p2(x,z)＝p(x,z1+1)*(t2(x,z)-t(x,z1))/dt+p(x,z1)*(t(x,z1+1)-t2(x,z))/dt
其中，(x,z)表示观测数据第x道第z个采样点，p(x,z)为观测数据D在(x,z)采样点处的波场值，p2(x,z)为映射炮集在(x,z)采样点处的波场值，t(x,z1)为观测数据在采样点(x,z1)处的采样时刻，t2(x,z)为映射炮集在(x,z)处的采样时刻，其中t2(x,z)、t(x,z1)、z1分别通过以下公式获得：
Figure FDA0002713548430000011
其中，gc(x)为第x道的初至窗G(x)，dt为观测数据的采样时间间隔，dtxin为映射炮集的采样时间间隔，nt为炮集的每道的采样点数，int为取整函数；
S3，对映射炮集做曲波变换得到曲波域内的曲波系数矩阵；
S4，基于所述曲波系数矩阵，利用凸集映射方法对映射炮集进行迭代重构得到重构炮集；
S5，利用反坐标映射公式对重构炮集反变换得到恢复炮集，具体包括：
Figure FDA0002713548430000012
p1(x,z)为采样点处(x,z)反变换后的波场值，组成恢复炮集Dpre，mod为取余函数。
2.如权利要求1所述的方法，所述方法还包括：
以所述恢复炮集作为先验约束信息约束，迭代生成目标地震数据，具体包括如下方式：
Figure FDA0002713548430000021
其中
Figure FDA0002713548430000022
为第k+1次迭代的目标地震数据，其中Dk+1为采用常规凸集映射公式重构的地震数据，Lmk+1为加权矩阵因子，I为单位矩阵。
3.如权利要求1所述的方法，所述S4中，基于所述曲波系数矩阵，利用凸集映射方法对映射炮集进行迭代重构得到重构炮集，包括：
基于预设的阈值矩阵，对所述曲波系数矩阵进行滤波，生成滤波后的曲波系数矩阵，其中，所述阈值矩阵为τk＝τ0*e(1-k)/M,k＝1,2...num，其中，τ0为给定的初始最大阈值，M为系数衰减项，k为迭代次数，num为总迭代次数；
基于滤波后的曲波系数矩阵，利用凸集映射方法对映射炮集进行迭代重构得到重构炮集。
CN201910474117.6A 2019-06-02 2019-06-02 一种针对部分缺失的地震数据的恢复重构方法 Active CN110398773B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910474117.6A CN110398773B (zh) 2019-06-02 2019-06-02 一种针对部分缺失的地震数据的恢复重构方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910474117.6A CN110398773B (zh) 2019-06-02 2019-06-02 一种针对部分缺失的地震数据的恢复重构方法

Publications (2)

Publication Number Publication Date
CN110398773A CN110398773A (zh) 2019-11-01
CN110398773B true CN110398773B (zh) 2020-12-08

Family

ID=68322968

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910474117.6A Active CN110398773B (zh) 2019-06-02 2019-06-02 一种针对部分缺失的地震数据的恢复重构方法

Country Status (1)

Country Link
CN (1) CN110398773B (zh)

Family Cites Families (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103955904B (zh) * 2014-05-12 2017-05-24 东南大学 一种基于离散分数阶傅里叶变换相位信息的信号重建方法
GB2553257A (en) * 2016-03-28 2018-03-07 Seismic Apparition Gmbh De-aliased source separation method
CN108983287B (zh) * 2018-10-23 2020-08-25 东华理工大学 一种基于凸集投影算法的曲波变换反假频地震数据重建方法

Also Published As

Publication number Publication date
CN110398773A (zh) 2019-11-01

Similar Documents

Publication Publication Date Title
CN103492910B (zh) 时域中的同步小波提取和反卷积
CN111368247B (zh) 基于快速正交字典的稀疏表征正则化叠前avo反演方法
CN110895348B (zh) 一种地震弹性阻抗低频信息提取方法、系统及存储介质
CN112394414B (zh) 两步法地震绕射波场叠前分离的方法
WO2005026776A1 (en) Wide-offset-range pre-stack depth migration method for seismic exploration
GB2538804A (en) Improved method for inversion modelling
CN111290019B (zh) 一种应用于最小二乘逆时偏移的l-bfgs初始矩阵求取方法
CN111505713B (zh) 基于多点地质统计的叠前地震反演方法
CN113077386A (zh) 基于字典学习和稀疏表征的地震资料高分辨率处理方法
AU2014101578A4 (en) System and method for seismic adaptive optics
CN109541691B (zh) 一种地震速度反演方法
CN110398773B (zh) 一种针对部分缺失的地震数据的恢复重构方法
WO2021127382A1 (en) Full waveform inversion in the midpoint-offset domain
CN111175821A (zh) 一种vti介质的各向异性参数分步反演方法
CN114185095B (zh) 一种三维平面波域地震数据多次波压制的方法
CN109100789A (zh) 基于波场逆时传播的全波形反演子波提取方法
CN113031070A (zh) 一种深度域合成地震记录的制作方法
CN104570100B (zh) 多子波克希霍夫地震数据偏移方法
Zhao et al. Cost reduction of RTM surface offset gathers by local coherent events based interpolation
CN108387928B (zh) 基于地震特征变换空间的数据构建方法
CN111722287B (zh) 一种基于渐进数据同化方法的震相特征识别波形反演方法
CN113900140B (zh) 一种基于时空组合的地震数据优化的方法及装置
CN114624765B (zh) 一种相位域地震数据处理与重构方法、装置及可存储介质
CN113031058B (zh) 一种基于反射系数精确式的页岩vti储层叠前混合反演方法
CN113960660B (zh) 一种基于正演模拟的动校正畸变区域自动识别及切除方法

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