CN112363222A - 叠后自适应宽带约束波阻抗反演方法及装置 - Google Patents

叠后自适应宽带约束波阻抗反演方法及装置 Download PDF

Info

Publication number
CN112363222A
CN112363222A CN202011175057.7A CN202011175057A CN112363222A CN 112363222 A CN112363222 A CN 112363222A CN 202011175057 A CN202011175057 A CN 202011175057A CN 112363222 A CN112363222 A CN 112363222A
Authority
CN
China
Prior art keywords
wave impedance
seismic
data
model
inversion
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
CN202011175057.7A
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 National Petroleum Corp
BGP Inc
Original Assignee
China National Petroleum Corp
BGP Inc
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 National Petroleum Corp, BGP Inc filed Critical China National Petroleum Corp
Priority to CN202011175057.7A priority Critical patent/CN112363222A/zh
Publication of CN112363222A publication Critical patent/CN112363222A/zh
Pending legal-status Critical Current

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/40Seismology; Seismic or acoustic prospecting or detecting specially adapted for well-logging
    • G01V1/44Seismology; Seismic or acoustic prospecting or detecting specially adapted for well-logging using generators and receivers in the same well
    • G01V1/48Processing data
    • G01V1/50Analysing data
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00Details of seismic processing or analysis
    • G01V2210/60Analysis
    • G01V2210/61Analysis by combining or comparing a seismic data set with other data
    • G01V2210/616Data from specific type of measurement
    • G01V2210/6169Data from specific type of measurement using well-logging

Landscapes

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

Abstract

本发明公开了一种叠后自适应宽带约束波阻抗反演方法及装置,该方法包括:获取目标储层的叠后地震数据、测井数据和地震层位解释数据;分析叠后地震数据的频带范围,利用测井数据和地震层位解释数据,建立该频带范围内初始的波阻抗模型;遍历每个地震道,根据每个地震道的地震数据和波阻抗模型,确定自适应阻尼因子;根据确定的自适应阻尼因子,确定对波阻抗模型进行迭代更新的波阻抗摄动量;根据确定的波阻抗摄动量,对每个地震道的波阻抗模型进行迭代更新,直到合成地震记录与实际地震数据的相关系数达到预设阈值,得到每个地震道的波阻抗反演结果:获取每个地震道的地震数据和波阻抗模型。本发明能够获得具有高分辨率和高稳定性的反演结果。

Description

叠后自适应宽带约束波阻抗反演方法及装置
技术领域
本发明涉及地质勘探技术领域,尤其涉及一种叠后自适应宽带约束波阻抗反演方法及装置。
背景技术
本部分旨在为权利要求书中陈述的本发明实施例提供背景或上下文。此处的描述不因为包括在本部分中就承认是现有技术。
波阻抗作为重要的岩石和地层属性,能够将地震数据与地质数据、测井数据关联起来。波阻抗反演是利用地震资料反演地层波阻抗的地震处理技术,相比于叠前地震数据,叠后地震数据具有高信噪比、低计算成本等优点,因而,基于叠后地震资料的波阻抗反演方法,被广泛应用于地下储层的定性或定量表征。
宽带约束波阻抗反演方法,是一种基于叠后地震资料进行波阻抗反演的方法,以叠后地震资料为基础,以地质数据、测井数据为约束条件,提出一个初始的波阻抗模型,并使用随机迭代算法对波阻抗模型进行修正,以获得一个具有较高分辨率的宽带波阻抗模型。目前,常用的宽带约束波阻抗反演方法,根据是否需要计算反射系数中间变量,可分为如下两类:
第一类是稀疏脉冲波阻抗反演方法,采用“两步法”来实现:第一步,在反射系数稀疏准则的约束下,对地震资料进行反褶积处理,得到宽带反射系数;第二步,利用递推算法、90度相位转换法或全变分正则法,获得绝对阻抗。常用的稀疏约束方法有:最大似然约束、稀疏范数(包括但不限于L1范数、柯西范数、LP范数等)约束、自回归模型谱拓展及稀疏层约束等。受反射系数“稀疏性”假设的限制,稀疏脉冲波阻抗反演方法只能识别具有较大反射系数的地层,从而丢失一些小反射系数信息,降低反演结果的分辨率。且在后续递推求解波阻抗的过程中,只利用了初始模型中起始点的信息,会导致反演结果的分辨率进一步降低,从而造成反演结果误差的传递和积累。
第二类是基于模型的波阻抗反演方法,采用“一步法”来实现:直接对初始模型进行迭代扰动,当正演模拟数据与观测数据达到一定的匹配标准时停止迭代,输出最终的绝对阻抗。在求解策略上,通常采用局部寻优的梯度法(例如,最速下降法、共轭梯度法、牛顿法)求解广义线性反问题,或采用全局寻优算法(例如,模拟退火、遗传算法、粒子群算法)求解非线性反问题,也可在上下介质模型参数差异较小的情况下,直接采用最小二乘法求解对数域模型参数的线性反问题。
分析可知,无论是采用“两步法”实现的稀疏脉冲波阻抗反演方法,还是采用“一步法”实现的模型波阻抗反演方法,都面临阻尼因子(也称正则化因子、惩罚项权重)主观选取的问题。由于不同的阻尼因子会出现区别较大的反演结果,因而,稀疏脉冲波阻抗反演方法和模型波阻抗反演方法,都会因主观选取阻尼因子而导致反演结果的不确定性。如果由井旁道的反演试验确定阻尼因子,单一的常数阻尼因子不利于整个叠后地震数据体的外推计算。
针对上述问题,目前尚未提出有效的解决方案。
发明内容
本发明实施例中提供了一种叠后自适应宽带约束地震波阻抗反演(AdaptiveBroadband Constraint Inversion,简称ABCI)方法,用以解决现有解决传统稀疏脉冲反演方法分辨率低,而基于模型的反演方法主观选取阻尼因子导致反演结果不确定性的技术问题,该方法包括:获取目标储层的叠后地震数据、测井数据和地震层位解释数据,叠后地震数据包括:多个地震道的地震数据;分析叠后地震数据的频带范围,利用测井数据和地震层位解释数据,建立频带范围内的每个地震道的波阻抗模型;遍历每个地震道,循环执行如下步骤,直到合成地震记录与实际地震数据的相关系数达到预设阈值,得到每个地震道的波阻抗反演结果:获取每个地震道的地震数据和波阻抗模型;根据每个地震道的地震数据和波阻抗模型,确定自适应阻尼因子;根据确定的自适应阻尼因子,确定对波阻抗模型进行迭代更新的波阻抗摄动量;根据确定的波阻抗摄动量,对每个地震道的波阻抗模型进行迭代更新。
本发明实施例中还提供了一种叠后自适应宽带约束地震波阻抗反演装置,用以解决现有解决传统稀疏脉冲反演方法分辨率低,而基于模型的反演方法主观选取阻尼因子导致反演结果不确定性的技术问题,该装置包括:数据获取模块,用于获取目标储层的叠后地震数据、测井数据和地震层位解释数据,叠后地震数据包括:多个地震道的地震数据;初始波阻抗模型建立模块,用于分析叠后地震数据的频带范围,利用测井数据和地震层位解释数据,建立频带范围内的每个地震道的波阻抗模型;波阻抗模型迭代更新模块,包括:数据选取模块、自适应阻尼因子计算模块、波阻抗摄动量确定模块和数据更新模块,用于循环执行数据选取模块、自适应阻尼因子计算模块、波阻抗摄动量确定模块和数据更新模块的功能,直到合成地震记录与实际地震数据的相关系数达到预设阈值,得到每个地震道的波阻抗反演结果;其中,数据选取模块,用于获取每个地震道的地震数据和波阻抗模型;自适应阻尼因子计算模块,用于根据每个地震道的地震数据和波阻抗模型,确定自适应阻尼因子;波阻抗摄动量确定模块,用于根据确定的自适应阻尼因子,确定对波阻抗模型进行迭代更新的波阻抗摄动量;数据更新模块,用于根据确定的波阻抗摄动量,对每个地震道的波阻抗模型进行迭代更新。
本发明实施例中还提供了一种计算机设备,用以解决现有解决传统稀疏脉冲反演方法分辨率低,而基于模型的反演方法主观选取阻尼因子导致反演结果不确定性的技术问题,该计算机设备包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序,处理器执行计算机程序时实现上述叠后自适应宽带约束波阻抗反演方法。
本发明实施例中还提供了一种计算机可读存储介质,用以解决现有解决传统稀疏脉冲反演方法分辨率低,而基于模型的反演方法主观选取阻尼因子导致反演结果不确定性的技术问题,该计算机可读存储介质存储有执行上述叠后自适应宽带约束波阻抗反演方法的计算机程序。
本发明实施例中,在获取到目标储层的叠后地震数据、测井数据和地震层位解释数据后,分析叠后地震数据的频带范围,并利用测井数据和地震层位解释数据,建立频带范围内的每个地震道的波阻抗模型(即初始波阻抗模型),进而遍历每个地震道,根据每个地震道的地震数据和波阻抗模型,确定自适应阻尼因子;根据确定的自适应阻尼因子,确定对波阻抗模型进行迭代更新的波阻抗摄动量;根据确定的波阻抗摄动量,对每个地震道的波阻抗模型进行迭代更新,直到合成地震记录与实际地震数据的相关系数达到预设阈值,得到每个地震道的波阻抗反演结果:获取每个地震道的地震数据和波阻抗模型,与现有技术中利用稀疏脉冲反演方法和基于模型的反演方法进行波阻抗反演的技术方案相比,本发明实施例中通过基于三维地震资料中噪声水平大小自动调节的自适应阻尼因子,应用于地震反演计算过程中,能够获得具有高分辨率和高稳定性的波阻抗反演结果。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案,下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本发明的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其他的附图。在附图中:
图1为本发明实施例中提供的一种叠后自适应宽带约束波阻抗反演方法流程图;
图2为本发明实施例中提供的一种可选的叠后自适应宽带约束波阻抗反演方法流程图;
图3为本发明实施例中提供的一维理论模型示意图;
图4为本发明实施例中提供的二维理论模型示意图;
图5为本发明实施例中提供的实际资料反演结果的剖面效果对比示意图;
图6为本发明实施例中提供的与图5对应的指控参数示意图;
图7为本发明实施例中提供的实际资料反演结果的平面效果对比示意图;
图8为本发明实施例中提供的一种叠后自适应宽带约束波阻抗反演方法的具体实现流程图;
图9为本发明实施例中提供的一种叠后自适应宽带约束波阻抗反演装置示意图。
具体实施方式
为使本发明实施例的目的、技术方案和优点更加清楚明白,下面结合附图对本发明实施例做进一步详细说明。在此,本发明的示意性实施例及其说明用于解释本发明,但并不作为对本发明的限定。
本发明实施例中提供了一种叠后自适应宽带约束地震波阻抗反演方法,通过基于三维地震资料中噪声水平大小进行自动调节的自适应阻尼因子,应用于地震反演计算过程中,能够获得具有高分辨率和高稳定性的波阻抗反演结果。
图1为本发明实施例中提供的一种叠后自适应宽带约束波阻抗反演方法流程图,如图1所示,该方法包括如下步骤:
S101,获取目标储层的叠后地震数据、测井数据和地震层位解释数据,叠后地震数据包括:多个地震道的地震数据。
需要说明的是,本发明实施例中目标储层可以是任意一个待反演研究区的储层。上述S101获取的叠后地震数据中,每个地震道对应一个二维的地震剖面记录。
S102,分析叠后地震数据的频带范围,利用测井数据和地震层位解释数据,建立频带范围内的每个地震道的波阻抗模型。
在通过分析确定目标储层叠后地震数据的频带范围后,可以利用目标储层的测井数据和地震层位解释数据,建立一个宽频带的波阻抗模型,作为迭代更新的初始波阻抗模型。
S103,遍历每个地震道,循环执行S1031~S1034,直到合成地震记录与实际地震数据的相关系数达到预设阈值,得到每个地震道的波阻抗反演结果。
由于不同地震道的地震噪声不同,因而,本发明实施例中,针对不同的地震道,提供不同的自适应阻尼因子,对不同的地震道进行波阻抗反演,能够大大提高反演结果的准确性。
S1031,获取每个地震道的地震数据和波阻抗模型。
需要说明的是,在通过上述S102建立各个地震道的初始波阻抗模型后,针对待反演的各个地震道,可以任意选取一个地震道的地震数据和波阻抗模型,以确定该地震道的自适应阻尼因子,对该地震道进行波阻抗迭代反演。
S1032,根据每个地震道的地震数据和波阻抗模型,确定自适应阻尼因子。
需要说明的是,上述S1032具体可以通过如下步骤来实现:根据每个地震道的地震数据和波阻抗模型,确定地震数据对波阻抗的梯度矩阵;根据梯度矩阵,确定波阻抗摄动量协方差矩阵和地震噪声协方差矩阵;根据波阻抗摄动量协方差矩阵和地震噪声协方差矩阵,确定自适应阻尼因子。
可选地,本发明实施例中提供的叠后自适应宽带约束波阻抗反演方法,确定梯度矩阵的表达式如下:
Figure BDA0002748475970000051
其中,G表示梯度矩阵;s表示合成地震记录;z表示波阻抗模型;r表示反射系数;W表示子波褶积矩阵;
Figure BDA0002748475970000052
为反射系数对波阻抗的偏导数。
在确定梯度矩阵后,可以通过如下公式(2)和(3)分别确定波阻抗摄动量协方差矩阵和地震噪声协方差矩阵:
Figure BDA0002748475970000061
Figure BDA0002748475970000062
其中,CΔz表示波阻抗摄动量协方差矩阵;
Figure BDA0002748475970000063
表示波阻抗摄动量的方差;Ce表示地震噪声协方差矩阵;
Figure BDA0002748475970000064
表示地震噪声的方差;I表示单位矩阵。
可选地,在确定自适应阻尼因子的时候,可以通过如下公式(4)来计算自适应阻尼因子:
Figure BDA0002748475970000065
其中,λ表示自适应阻尼因子。
S1033,根据确定的自适应阻尼因子,确定对波阻抗模型进行迭代更新的波阻抗摄动量。
本发明实施例中提供了一种计算波阻抗摄动量的表达式,如公式(5)所示,可以计算出每次对波阻抗模型进行迭代更新的波阻抗摄动量:
Figure BDA0002748475970000066
其中,Δzk表示第k次对波阻抗模型进行迭代更新时的波阻抗摄动量;Gk-1表示第k-1次对波阻抗模型进行迭代更新时的梯度矩阵;
Figure BDA0002748475970000067
表示Gk-1的转置;dobs表示实际地震数据;sk-1表示第k-1次对波阻抗模型进行迭代更新后得到的合成地震记录;λk-1表示第k-1次对波阻抗模型进行迭代更新时确定的自适应阻尼因子。
S1034,根据确定的波阻抗摄动量,对每个地震道的波阻抗模型进行迭代更新。
在具体实施时,可以利用公式(6)所示的波阻抗反演的迭代表达式,对每个地震道的波阻抗模型进行迭代更新:
zk=zk-1+Δzk (6)
其中,zk表示第k次迭代更新后得到的波阻抗模型;zk-1表示第k-1次迭代更新后得到的波阻抗模型。
在一个实施例中,如图2所示,本发明实施例中提供的叠后自适应宽带约束波阻抗反演方法还可以包括如下步骤:
S104,根据每个地震道的波阻抗反演结果,生成目标储层的波阻抗反演体;
S105,输出目标储层的波阻抗反演体。
通过遍历每个地震道,对每个地震道进行波阻抗反演后,可以根据各个地震道的波阻抗反演结果合成为目标储层的波阻抗反演体,进而输出该目标储层的波阻抗繁衍体。
本发明实施例中,利用泰勒级数展开方法构建波阻抗摄动量的正演方程,构建出了波阻抗摄动量的拟线性正演方程;基于贝叶斯理论推导出了自适应阻尼因子和波阻抗反演的解析迭代表达式,以直接反演绝对阻抗,获得高精度反演结果。
下面具体论述其理论推导过程:
假设平面界面两侧为各向同性、弹性介质,对于垂直入射的反射波,根据褶积模型可将观测地震数据dobs表示为:
dobs=s+e=Wr+e (7)
其中,
r=[r1,r2,...,rn-1]T (8)
ri=(zi+1-zi)/(zi+1+zi) (9)
其中,s表示合成地震记录;r表示反射系数序列;zi表示第i层介质波阻抗,W表示子波褶积矩阵,e为误差项,下标n为波阻抗模型参数个数。
由于合成地震记录s与波阻抗z之间为非线性关系,本发明实施例中,采用泰勒级数展开方法求解此类问题。将合成地震记录s在初始波阻抗z0处进行一阶泰勒展开,并且省略二次及以上的高阶项,可得:
Figure BDA0002748475970000073
其中,
Figure BDA0002748475970000071
表示Jacobi矩阵,即上述公式(1)所示的梯度矩阵。将公式(8)代入公式(7),并且令Δz=z-z0,Δs=dobs-s(z0),可得波阻抗摄动量的正演方程为:
Δs=GΔz+e (11)
由于地球物理反问题通常具有不适定性,直接通过最小化目标函数
Figure BDA0002748475970000072
来求取Δz会造成反演结果的强烈不稳定性。通常采用Tikhonov正则化理论获得稳定的反演解,比如最小长度解、最光滑解和最稀疏解等。然而正则化参数通常为固定的常数,选择时受人为主观因素影响大,在实际三维地震资料反演时不能保证各道反演结果的分辨率和反演结果都达到最佳折中。
贝叶斯反演方法能够综合多种先验信息,从而改善反问题求解的不适定性问题。假设波阻抗摄动量和地震噪声均服从零均值高斯分布,则先验信息和似然函数的概率密度函数可以分别表示为:
Figure BDA0002748475970000081
Figure BDA0002748475970000082
其中,CΔz和Ce分别为波阻抗摄动量和地震噪声的协方差矩阵,
Figure BDA0002748475970000083
Figure BDA0002748475970000084
分别为波阻抗摄动量和地震噪声的方差,c1和c2为归一化常系数。
根据贝叶斯理论,由公式(12)和公式(13)可得波阻抗摄动量的后验概率密度函数为:
Figure BDA0002748475970000085
求解公式(14)的最大后验概率解,相当于最小化目标函数:
Figure BDA0002748475970000086
对公式(15)求取关于波阻抗摄动量的导数并令其为零,可得最终波阻抗反演的迭代表达式为:
Figure BDA0002748475970000087
由公式(4)可得出
Figure BDA0002748475970000088
分析可知,当地震资料噪声水平较高时,
Figure BDA0002748475970000089
增大,因此阻尼因子也会增大,从而压制噪声获得相对稳定的反演结果;相反地,当地震资料几乎不含噪声时,
Figure BDA00027484759700000810
趋近零,此时阻尼因子也接近零,反演结果分辨率高。可见,本发明实施例中提供的自适应阻尼因子,在地震反演计算的过程中,能够基于三维地震资料中的噪声水平大小进行自动调节,以获得具有最佳分辨率和稳定性的反演结果。
为了验证本发明实施例的有效性,分别开展一维、二维理论模型,以及三维实际资料的方法测试。图3示出的为一维理论模型,理论波阻抗曲线采用无条件序贯高斯模拟方法获得。在模型参数设置中,设定波阻抗参数的期望为7200m/s·g/cm3,标准差为100m/s·g/cm3,并且模型空间结构设置为20ms的指数变差函数。图3中(a)所示为使用图3中(b)、(c)、(d)中所示的理论波阻抗曲线与30Hz雷克地震子波计算得到的合成地震记录,并且分别添加了不同程度的随机噪声。图3中(b)、(c)、(d)为以不同噪声含量的合成地震数据作为输入时,使用本发明实施例中提供的叠后自适应宽带约束波阻抗反演方法得到的反演结果对比,其中实线为理论模型,点划线为均值初始模型,点线为反演结果。可以看出随着地震数据中噪声含量的增加,反演结果的分辨率逐渐降低,阻尼因子也逐渐增大,表明了阻尼因子的自适应变化特征。同时,当地震噪声含量为30%时,可以看出反演结果与理论模型的趋势仍然较为一致,也表明了本发明实施例中提供的叠后自适应宽带约束波阻抗反演方法,具有较好的稳定性。
图4示出的为二维理论模型,图4中(a)所示为波阻抗的真实理论模型,包含了楔形体、透镜体、薄层和断层等复杂地质体信息,两条曲线为从中抽取的第30道和第90道波阻抗数据,代表测井数据。图4中(b)所示为合成地震记录,该合成记录由图4中(a)所示的波阻抗模型与30Hz雷克地震子波计算得到,并且添加了10%的随机地震噪声。图4中(c)所示为两口井沿层插值得到的波阻抗初始模型,其中测井数据进行了0-120Hz的带通滤波。可以看出,在此初始模型中,无法识别楔形体、透镜体和薄层的边界位置信息,并且断层的位置也发生了错乱。图4中(d)所示为采用本发明实施例中提供的叠后自适应宽带约束波阻抗反演方法获得的波阻抗反演结果,可以看出反演结果接近真实参考模型,能够准确反映出初始模型中所没有的楔形体、透镜体和薄层的边界位置信息,并且在一定程度上更加准确的反映了断层的位置,表明了反演结果的高分辨率优势。由于初始模型仅依赖井和层位建立,从反演结果对初始模型的改变可以看出,本发明实施例中提供的叠后自适应宽带约束波阻抗反演方法,获得的反演结果主要依赖地震资料,反演结果较为可靠。
图5为实际资料反演结果的剖面效果对比,数据来自伊拉克某研究区,利用研究区中的测井数据进行沿层插值,建立0-120Hz初始模型,然后分别利用本发明实施例中提供的叠后自适应宽带约束波阻抗反演方法和传统稀疏脉冲反演方法进行波阻抗反演。图5中(a)所示为反演时输入的地震数据,地震资料的主频约为35Hz,剖面中两口井经过了精细的井震标定,合成记录与井旁地震道匹配较好。图5中(b)和(c)所示为两种方法的反演结果,可以看出两种方法反演结果在整体趋势上一致,从浅到深波阻抗逐渐增加。但是本发明实施例中提供的叠后自适应宽带约束波阻抗反演方法得到的反演结果,其分辨率明显高于传统稀疏脉冲反演方法,如图5中箭头所示,本发明实施例中提供的叠后自适应宽带约束波阻抗反演方法获得的反演结果,能够清晰指示薄层位置,而在传统稀疏脉冲反演方法获得的反演结果中,这些薄层仅有微弱显示,分辨率与地震数据相当。从反演结果与地震波形特征的关系来看,如图5中椭圆圈所示,本发明实施例中提供的叠后自适应宽带约束波阻抗反演方法获得的反演结果,能更加清晰地展示复波内的薄层特征。
图6为图5中地震数据反演时的质控参数分析,这里反演过程中迭代阈值设置为85%,即当反演结果的合成地震记录与实际地震道相关系数达到85%则停止迭代,继续反演下一道。图6中(a)所示为各道迭代终止时的阻尼因子,可以看出,在反演的过程中各道的阻尼因子各不相同,其主要受地震噪声的影响,用于自动调节反演结果的分辨率和稳定性。图6中(b)和(c)所示分别为各道的最终迭代次数和相关系数,可以看出经过3次迭代修正,反演结果合成记录与实际地震道的相关系数均能达到95%以上。当然,在研究过程中也发现,相关系数的阈值并不是设置得越高,反演结果越好。当地震资料噪声含量较高时,如果相关系数的阈值设置过高,会造成部分地震道达不到此相关系数,迭代次数增加,反演结果与相邻地震道结果差异较大,从而产生地震反演中较为常见的“挂面条”现象。
图7为实际资料反演结果的平面效果对比,数据来自中国山西某煤层气研究区,目的层为一套具有明显低阻抗特征的煤系地层。相似地,针对本发明实施例中提供的叠后自适应宽带约束波阻抗反演方法和传统稀疏脉冲反演方法获得的反演结果进行对比。图7中(a)、(b)、(c)所示分别为初始模型数据体、利用本发明实施例中叠后自适应宽带约束波阻抗反演方法得到的反演数据体和利用稀疏脉冲反演方法得到的稀疏脉冲反演数据体,沿目的层顶界面上下各10ms时窗计算的波阻抗均方根切片,切片中的数值代表实际钻遇的煤层厚度。对比图7中(a)和(b)中箭头所指区域,可以看出本发明实施例获得的反演结果切片上,钻井显示的较厚煤层与低阻抗值对应关系较好,相反,在初始模型切片上,这些区域均对应高阻抗值,可以说明本发明实施例中叠后自适应宽带约束波阻抗反演方法对初始模型的依赖程度小,主要来源于地震。对比图7中(b)和(c)中所示的椭圆区域,同样可以看出本发明实施例获得的反演结果切片与实钻井上煤层厚度更加符合,并且对于工区东南侧煤层的横向变化刻画更加精细。
通过对一维理论模型、二维理论模型和实际资料的测试结果,可以验证本发明实施例中叠后自适应宽带约束波阻抗反演方法,具有阻尼因子自适应特征和反演结果高分辨率等优势,同时说明了本发明实施例中叠后自适应宽带约束波阻抗反演方法,对初始模型的依赖性较小,并且反演结果与实钻井资料更加吻合,这表明本发明实施例中叠后自适应宽带约束波阻抗反演方法,在薄层勘探中具有较好的应用前景。需要注意的是,由于本发明实施例中叠后自适应宽带约束波阻抗反演方法中宽带初始模型的建立主要依赖测井和地震解释层位,因此在构造复杂或地层厚度变化巨大的地质条件下,应用本发明实施例中叠后自适应宽带约束波阻抗反演方法会存一定的限制。
图8为本发明实施例中提供的一种叠后自适应宽带约束波阻抗反演方法的具体实现流程图,如图8所示,具体包括如下步骤:
(1)读取输入数据,包括地震数据、测井数据和层位解释结果;
(2)分析地震数据的频带范围,利用测井数据和层位解释结果,建立宽带波阻抗初始模型,并且提取地震子波;
(3)选择当前道的地震数据dobs和波阻抗初始模型z0
(4)根据泰勒级数展开,计算地震数据对波阻抗的梯度矩阵(Jacobi矩阵)G;
(5)计算波阻抗摄动量协方差矩阵
Figure BDA0002748475970000111
和地震噪声协方差矩阵
Figure BDA0002748475970000112
从而计算自适应阻尼因子λ;
(6)反演波阻抗摄动量Δz,并修正波阻抗模型z;
(7)重复步骤(3)~(6),迭代更新波阻抗模型,直至合成记录与实际地震数据的相关系数达到设定的阈值,获得当前道的波阻抗反演结果;
(8)重复步骤(3)~(7),逐道循环,直至所有地震道反演结束,输出最终的波阻抗反演体。
基于同一发明构思,本发明实施例中还提供了一种叠后自适应宽带约束波阻抗反演装置,如下面的实施例所述。由于该装置解决问题的原理与叠后自适应宽带约束波阻抗反演方法相似,因此该装置的实施可以参见叠后自适应宽带约束波阻抗反演方法的实施,重复之处不再赘述。
图9为本发明实施例中提供的一种叠后自适应宽带约束波阻抗反演装置示意图,如图9所示,该装置包括:数据获取模块91、初始波阻抗模型建立模块92和波阻抗模型迭代更新模块93;波阻抗模型迭代更新模块93包括:数据选取模块931、自适应阻尼因子计算模块932、波阻抗摄动量确定模块933和数据更新模块934。
其中,数据获取模块91,用于获取目标储层的叠后地震数据、测井数据和地震层位解释数据,叠后地震数据包括:多个地震道的地震数据;初始波阻抗模型建立模块92,用于分析叠后地震数据的频带范围,利用测井数据和地震层位解释数据,建立频带范围内的每个地震道的波阻抗模型;波阻抗模型迭代更新模块93,用于循环执行数据选取模块931、自适应阻尼因子计算模块932、波阻抗摄动量确定模块933和数据更新模块934的功能,直到合成地震记录与实际地震数据的相关系数达到预设阈值,得到每个地震道的波阻抗反演结果;数据选取模块931,用于获取每个地震道的地震数据和波阻抗模型;自适应阻尼因子计算模块932,用于根据每个地震道的地震数据和波阻抗模型,确定自适应阻尼因子;波阻抗摄动量确定模块933,用于根据确定的自适应阻尼因子,确定对波阻抗模型进行迭代更新的波阻抗摄动量;数据更新模块934,用于根据确定的波阻抗摄动量,对每个地震道的波阻抗模型进行迭代更新。
在一个实施例中,上述自适应阻尼因子计算模块932还用于:根据每个地震道的地震数据和波阻抗模型,确定地震数据对波阻抗的梯度矩阵;根据梯度矩阵,确定波阻抗摄动量协方差矩阵和地震噪声协方差矩阵;根据波阻抗摄动量协方差矩阵和地震噪声协方差矩阵,确定自适应阻尼因子。
可选地,上述自适应阻尼因子计算模块932还用于公式(1)所示的梯度矩阵的表达式,确定地震数据对波阻抗的梯度矩阵。
进一步地,上述自适应阻尼因子计算模块932还用于通过公式(2)所示的的波阻抗摄动量协方差矩阵确定波阻抗摄动量协方差矩阵,以及根据公式(3)所示的地震噪声协方差矩阵的表达式确定地震噪声协方差矩阵。
更进一步地,上述自适应阻尼因子计算模块932还用于根据公式(4)所示的自适应阻尼因子的表达式为,确定自适应阻尼因子。
在一个实施例中,上述波阻抗摄动量确定模块933还用于根据公式(5)所示的波阻抗摄动量的表达式,确定对波阻抗模型进行迭代更新的波阻抗摄动量。
在一个实施例中,上述数据更新模块934还用于利用公式(6)所示的波阻抗反演的迭代表达式,对每个地震道的波阻抗模型进行迭代更新。
在一个实施例中,本发明实施例中提供的叠后自适应宽带约束波阻抗反演装置还可以包括:波阻抗反演结果合成模块94,用于根据每个地震道的波阻抗反演结果,生成目标储层的波阻抗反演体;波阻抗反演体输出模块95,用于输出目标储层的波阻抗反演体。
基于同一发明构思,本发明实施例中还提供了一种计算机设备,用以解决现有解决传统稀疏脉冲反演方法分辨率低,而基于模型的反演方法主观选取阻尼因子导致反演结果不确定性的技术问题,该计算机设备包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序,处理器执行计算机程序时实现上述叠后自适应宽带约束波阻抗反演方法。
基于同一发明构思,本发明实施例中还提供了一种计算机可读存储介质,用以解决现有解决传统稀疏脉冲反演方法分辨率低,而基于模型的反演方法主观选取阻尼因子导致反演结果不确定性的技术问题,该计算机可读存储介质存储有执行上述叠后自适应宽带约束波阻抗反演方法的计算机程序。
综上所述,本发明实施例中了一种叠后自适应宽带约束波阻抗反演方法、装置、计算机设备及计算机可读存储介质,在获取到目标储层的叠后地震数据、测井数据和地震层位解释数据后,分析叠后地震数据的频带范围,并利用测井数据和地震层位解释数据,建立频带范围内的每个地震道的波阻抗模型(即初始波阻抗模型),进而遍历每个地震道,根据每个地震道的地震数据和波阻抗模型,确定自适应阻尼因子;根据确定的自适应阻尼因子,确定对波阻抗模型进行迭代更新的波阻抗摄动量;根据确定的波阻抗摄动量,对每个地震道的波阻抗模型进行迭代更新,直到合成地震记录与实际地震数据的相关系数达到预设阈值,得到每个地震道的波阻抗反演结果:获取每个地震道的地震数据和波阻抗模型,与现有技术中利用稀疏脉冲反演方法和基于模型的反演方法进行波阻抗反演的技术方案相比,本发明实施例中通过基于三维地震资料中噪声水平大小自动调节的自适应阻尼因子,应用于地震反演计算过程中,能够获得具有高分辨率和高稳定性的波阻抗反演结果。
本领域内的技术人员应明白,本发明的实施例可提供为方法、系统、或计算机程序产品。因此,本发明可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且,本发明可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等)上实施的计算机程序产品的形式。
本发明是参照根据本发明实施例的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器,使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。
这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中,使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品,该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。
这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上,使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理,从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。
以上所述的具体实施例,对本发明的目的、技术方案和有益效果进行了进一步详细说明,所应理解的是,以上所述仅为本发明的具体实施例而已,并不用于限定本发明的保护范围,凡在本发明的精神和原则之内,所做的任何修改、等同替换、改进等,均应包含在本发明的保护范围之内。

Claims (12)

1.一种叠后自适应宽带约束波阻抗反演方法,其特征在于,包括:
获取目标储层的叠后地震数据、测井数据和地震层位解释数据,所述叠后地震数据包括:多个地震道的地震数据;
分析所述叠后地震数据的频带范围,利用所述测井数据和所述地震层位解释数据,建立所述频带范围内的每个地震道的波阻抗模型;
遍历每个地震道,循环执行如下步骤,直到合成地震记录与实际地震数据的相关系数达到预设阈值,得到每个地震道的波阻抗反演结果:
获取每个地震道的地震数据和波阻抗模型;
根据每个地震道的地震数据和波阻抗模型,确定自适应阻尼因子;
根据确定的自适应阻尼因子,确定对波阻抗模型进行迭代更新的波阻抗摄动量;
根据确定的波阻抗摄动量,对每个地震道的波阻抗模型进行迭代更新。
2.如权利要求1所述的方法,其特征在于,根据每个地震道的地震数据和波阻抗模型,确定自适应阻尼因子,包括:
根据每个地震道的地震数据和波阻抗模型,确定地震数据对波阻抗的梯度矩阵;
根据所述梯度矩阵,确定波阻抗摄动量协方差矩阵和地震噪声协方差矩阵;
根据所述波阻抗摄动量协方差矩阵和所述地震噪声协方差矩阵,确定自适应阻尼因子。
3.如权利要求2所述的方法,其特征在于,所述梯度矩阵的表达式为:
Figure FDA0002748475960000011
其中,G表示梯度矩阵;s表示合成地震记录;z表示波阻抗模型;r表示反射系数;W表示子波褶积矩阵;
Figure FDA0002748475960000012
为反射系数对波阻抗的偏导数。
4.如权利要求3所述的方法,其特征在于,所述波阻抗摄动量协方差矩阵和所述地震噪声协方差矩阵的表达式分别为:
Figure FDA0002748475960000013
Figure FDA0002748475960000014
其中,CΔz表示波阻抗摄动量协方差矩阵;
Figure FDA0002748475960000015
表示波阻抗摄动量的方差;Ce表示地震噪声协方差矩阵;
Figure FDA0002748475960000016
表示地震噪声的方差;I表示单位矩阵。
5.如权利要求4所述的方法,其特征在于,所述自适应阻尼因子的表达式为:
Figure FDA0002748475960000021
其中,λ表示自适应阻尼因子。
6.如权利要求5所述的方法,其特征在于,所述波阻抗摄动量的表达式为:
Figure FDA0002748475960000022
其中,Δzk表示第k次对波阻抗模型进行迭代更新时的波阻抗摄动量;Gk-1表示第k-1次对波阻抗模型进行迭代更新时的梯度矩阵;
Figure FDA0002748475960000023
表示Gk-1的转置;dobs表示实际地震数据;sk-1表示第k-1次对波阻抗模型进行迭代更新后得到的合成地震记录;λk-1表示第k-1次对波阻抗模型进行迭代更新时确定的自适应阻尼因子。
7.如权利要求5所述的方法,其特征在于,利用如下公式,对每个地震道的波阻抗模型进行迭代更新:
zk=zk-1+Δzk
其中,zk表示第k次迭代更新后得到的波阻抗模型;zk-1表示第k-1次迭代更新后得到的波阻抗模型。
8.如权利要求1至7任一项所述的方法,其特征在于,所述方法还包括:
根据每个地震道的波阻抗反演结果,生成所述目标储层的波阻抗反演体;
输出目标储层的波阻抗反演体。
9.一种叠后自适应宽带约束波阻抗反演方法装置,其特征在于,包括:
数据获取模块,用于获取目标储层的叠后地震数据、测井数据和地震层位解释数据,所述叠后地震数据包括:多个地震道的地震数据;
初始波阻抗模型建立模块,用于分析所述叠后地震数据的频带范围,利用所述测井数据和所述地震层位解释数据,建立所述频带范围内的每个地震道的波阻抗模型;
波阻抗模型迭代更新模块,包括:数据选取模块、自适应阻尼因子计算模块、波阻抗摄动量确定模块和数据更新模块,用于循环执行所述数据选取模块、所述自适应阻尼因子计算模块、所述波阻抗摄动量确定模块和所述数据更新模块的功能,直到合成地震记录与实际地震数据的相关系数达到预设阈值,得到每个地震道的波阻抗反演结果;
其中,数据选取模块,用于获取每个地震道的地震数据和波阻抗模型;自适应阻尼因子计算模块,用于根据每个地震道的地震数据和波阻抗模型,确定自适应阻尼因子;波阻抗摄动量确定模块,用于根据确定的自适应阻尼因子,确定对波阻抗模型进行迭代更新的波阻抗摄动量;数据更新模块,用于根据确定的波阻抗摄动量,对每个地震道的波阻抗模型进行迭代更新。
10.如权利要求9所述的装置,其特征在于,所述装置还包括:
波阻抗反演体输出模块,用于根据每个地震道的波阻抗反演结果,生成所述目标储层的波阻抗反演体;以及输出目标储层的波阻抗反演体。
11.一种计算机设备,包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序,其特征在于,所述处理器执行所述计算机程序时实现权利要求1至8任一项所述叠后自适应宽带约束波阻抗反演方法。
12.一种计算机可读存储介质,其特征在于,所述计算机可读存储介质存储有执行权利要求1至8任一项所述叠后自适应宽带约束波阻抗反演方法的计算机程序。
CN202011175057.7A 2020-10-28 2020-10-28 叠后自适应宽带约束波阻抗反演方法及装置 Pending CN112363222A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011175057.7A CN112363222A (zh) 2020-10-28 2020-10-28 叠后自适应宽带约束波阻抗反演方法及装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011175057.7A CN112363222A (zh) 2020-10-28 2020-10-28 叠后自适应宽带约束波阻抗反演方法及装置

Publications (1)

Publication Number Publication Date
CN112363222A true CN112363222A (zh) 2021-02-12

Family

ID=74511288

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011175057.7A Pending CN112363222A (zh) 2020-10-28 2020-10-28 叠后自适应宽带约束波阻抗反演方法及装置

Country Status (1)

Country Link
CN (1) CN112363222A (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117434604A (zh) * 2023-02-20 2024-01-23 中国石油化工股份有限公司 地震数据处理方法、装置、存储介质及电子设备

Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO1996007935A1 (en) * 1994-09-02 1996-03-14 Exxon Production Research Company Method for deriving reservoir lithology and fluid content from pre-stack inversion of seismic data
CN104122581A (zh) * 2013-04-28 2014-10-29 中国石油化工股份有限公司 一种叠后声波阻抗反演方法
WO2016008105A1 (zh) * 2014-07-15 2016-01-21 杨顺伟 一种基于柯西分布的叠后波阻抗反演方法
US20160131781A1 (en) * 2014-11-12 2016-05-12 Chevron U.S.A. Inc. Creating a high resolution velocity model using seismic tomography and impedance inversion
CN106291682A (zh) * 2015-06-01 2017-01-04 中国石油化工股份有限公司 一种基于基追踪方法的叠后声波阻抗反演方法
CN106291677A (zh) * 2015-05-22 2017-01-04 中国石油化工股份有限公司 一种基于匹配追踪方法的叠后声波阻抗反演方法
WO2017024702A1 (zh) * 2015-08-11 2017-02-16 深圳朝伟达科技有限公司 一种射线弹性参数的反演系统
CN109061764A (zh) * 2018-09-07 2018-12-21 中国石油化工股份有限公司 一种分频融合波阻抗反演方法

Patent Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO1996007935A1 (en) * 1994-09-02 1996-03-14 Exxon Production Research Company Method for deriving reservoir lithology and fluid content from pre-stack inversion of seismic data
CN104122581A (zh) * 2013-04-28 2014-10-29 中国石油化工股份有限公司 一种叠后声波阻抗反演方法
WO2016008105A1 (zh) * 2014-07-15 2016-01-21 杨顺伟 一种基于柯西分布的叠后波阻抗反演方法
US20160131781A1 (en) * 2014-11-12 2016-05-12 Chevron U.S.A. Inc. Creating a high resolution velocity model using seismic tomography and impedance inversion
CN106291677A (zh) * 2015-05-22 2017-01-04 中国石油化工股份有限公司 一种基于匹配追踪方法的叠后声波阻抗反演方法
CN106291682A (zh) * 2015-06-01 2017-01-04 中国石油化工股份有限公司 一种基于基追踪方法的叠后声波阻抗反演方法
WO2017024702A1 (zh) * 2015-08-11 2017-02-16 深圳朝伟达科技有限公司 一种射线弹性参数的反演系统
CN109061764A (zh) * 2018-09-07 2018-12-21 中国石油化工股份有限公司 一种分频融合波阻抗反演方法

Non-Patent Citations (5)

* Cited by examiner, † Cited by third party
Title
张玉芬等: "多参数约束反演及其应用", 《地球科学-中国地质大学学报》, vol. 22, no. 02, 28 February 1997 (1997-02-28), pages 219 - 222 *
林小竹等: "有井多道反演", 《石油物探》, vol. 38, no. 04, 31 December 2009 (2009-12-31), pages 44 - 50 *
梁光河: "多井约束随机逆波阻抗反演研究", 《石油地球物理勘探》, vol. 33, no. 02, 30 April 1998 (1998-04-30), pages 151 - 160 *
王西文等: "地震波阻抗反演方法研究", 《岩性油气藏》, vol. 19, no. 03, 30 September 2007 (2007-09-30), pages 80 - 88 *
黄捍东等: "地震相控混沌反演在岩性油气藏储层预测中的应用——以 RY 坳陷 M 井区为例", 《CT理论与应用研究》, vol. 28, no. 5, 31 October 2019 (2019-10-31), pages 549 - 557 *

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117434604A (zh) * 2023-02-20 2024-01-23 中国石油化工股份有限公司 地震数据处理方法、装置、存储介质及电子设备

Similar Documents

Publication Publication Date Title
US8923093B2 (en) Determining the quality of a seismic inversion
CN108535775B (zh) 非平稳地震资料声波阻抗反演方法及装置
CN109254324B (zh) 全频保幅地震数据处理方法和装置
CN111123354B (zh) 基于频变反射振幅衰减预测致密气层的方法及设备
GB2414802A (en) Seismic event correlation and Vp-Vs estimation
US10310117B2 (en) Efficient seismic attribute gather generation with data synthesis and expectation method
CN110187384B (zh) 贝叶斯时移地震差异反演方法及装置
CN103616723A (zh) 基于avo特征的crp道集真振幅恢复方法
CN113093272A (zh) 基于卷积编码的时间域全波形反演方法
CN108957532B (zh) 储层预测方法及装置
WO2005026776A1 (en) Wide-offset-range pre-stack depth migration method for seismic exploration
CN110687597B (zh) 一种基于联合字典的波阻抗反演方法
CN111856566A (zh) 湖相滩坝砂体中薄储层预测方法及装置
CN111060961B (zh) 基于多信息约束反演的品质因子确定方法、装置及系统
US11340366B2 (en) Accurate velocity model estimation and imaging in the presence of localized attenuation (Q) anomalies
CN112711065A (zh) 叠前地震反演方法及装置
CN112363222A (zh) 叠后自适应宽带约束波阻抗反演方法及装置
CN114428343A (zh) 基于归一化互相关的Marchenko成像方法及系统
CN107678065B (zh) 提高地震分辨率的保构造井控空间反褶积方法和装置
CN111239805B (zh) 基于反射率法的块约束时移地震差异反演方法及系统
CN110988991B (zh) 一种弹性参数反演方法、装置及系统
CN114114421A (zh) 基于深度学习的导向自学习地震数据去噪方法及装置
CN114137606A (zh) 一种稳健的谱模拟反褶积方法
CN113806674A (zh) 古河道纵向尺度的量化方法、装置、电子设备及存储介质
CN111077573A (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