CN112711065A - 叠前地震反演方法及装置 - Google Patents

叠前地震反演方法及装置 Download PDF

Info

Publication number
CN112711065A
CN112711065A CN201911023165.XA CN201911023165A CN112711065A CN 112711065 A CN112711065 A CN 112711065A CN 201911023165 A CN201911023165 A CN 201911023165A CN 112711065 A CN112711065 A CN 112711065A
Authority
CN
China
Prior art keywords
seismic
inversion
stack
dimensionless
determining
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
CN201911023165.XA
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 CN201911023165.XA priority Critical patent/CN112711065A/zh
Publication of CN112711065A publication Critical patent/CN112711065A/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/28Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/282Application of seismic models, synthetic seismograms
    • 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/30Analysis
    • 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
    • G01V1/362Effecting static or dynamic corrections; Stacking

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

本发明提供了一种叠前地震反演方法及装置,其中,该方法包括以下步骤:根据地震数据,确定地震中噪声类别和叠前地震正演模型;根据叠前地震正演模型和最小梯度支持正则化项,确定叠前反演残差模型;根据权重参数和叠前反演残差模型,构建自适应混合范数反演目标函数;根据自适应混合范数反演目标函数、叠前正演模型和地震中噪声类别,确定叠前地震反演结果。本发明通过判断地震数据中噪声类别来选择范数从而可以自适应压制噪声,解决了叠前反演中噪声的自适应性技术问题,又通过加入最小梯度正则化项使反演得到稳定且唯一的近似解,得到更聚焦的、边界更为清晰地质构造特征,提高了反演的稳定性和有效性。

Description

叠前地震反演方法及装置
技术领域
本发明涉及地球物理勘探解释技术领域,特别涉及一种叠前地震反演方法及装置。
背景技术
地震反演受频带宽度和噪声等的影响很难获得稳定且相对唯一的结果,通常我们都假设地震噪声符合高斯分布,但实际地震记录中,噪声分布非常复杂,远非简单的高斯分布能够描述。尤其对于地震反演来说,噪声类型应该是多种分布的叠加。如Omid Karimi在Arild Buland研究基础上,讨论了地震噪声的分布情况,认为地震噪声服从偏高斯分布,发展了偏高斯分布的贝叶斯AVO反演(Karimi and Omre,2010)。Liu和Zhang等讨论了基于L1、L2混合范数的三参数非高斯反演方法(Liu et al,2015)。以上方法虽然较好地解决了地震数据的非高斯问题,但都不能根据地震噪声分布特征进行自适应反演。
由于地震数据的带限性和不完备性,地震反演问题涉及到的***方程组呈现出病态性,带来的后果就是观测数据中很小的噪声干扰,都会造成方程组的解出现很大的波动。因此在实际地震记录存在噪声的情况下,直接求解目标函数,得到的解往往是不稳定的、偏离真实情况的。
针对上述问题,目前尚未提出有效的解决方案。
发明内容
本发明实施例提供了一种叠前地震反演方法及装置,解决了现有技术中地震数据不能根据地震噪声分布特征进行自适应反演以及地震反演结果易被噪声干扰导致数据不准确的技术问题。
本发明实施例提供了一种叠前地震反演方法,该方法包括:
根据地震数据,确定地震中噪声类别和叠前地震正演模型;
根据叠前地震正演模型和最小梯度支持正则化项,确定叠前反演残差模型;
根据权重参数和叠前反演残差模型,构建自适应混合范数反演目标函数;
根据自适应混合范数反演目标函数、叠前正演模型和地震中噪声类别,确定叠前地震反演结果。
本发明实施例还提供了一种叠前地震反演装置,该装置包括:
噪声分类模块,用于根据地震数据,确定地震中噪声类别;
叠前地震正演模型构建模块,用于根据地震数据,确定叠前正演模型;
叠前反演残差模型构建模块,用于根据叠前地震正演模型和最小梯度支持正则化项,确定叠前反演残差模型;
反演目标函数构建模块,用于根据权重参数和叠前反演残差模型,构建自适应混合范数反演目标函数;
反演数据获取模块,用于根据自适应混合范数反演目标函数、叠前正演模型和地震中噪声类别,确定叠前地震反演结果。
本发明实施例还提供了一种计算机设备,包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序,所述处理器执行所述计算机程序时实现上述所述方法。
本发明实施例还提供了一种计算机可读存储介质,所述计算机可读存储介质存储有执行上述所述方法的计算机程序。
本发明实施例提供技术方案的有益技术效果是:
本发明实施例的方案通过判断地震数据中噪声类别来选择范数从而可以自适应压制噪声,解决了叠前反演中噪声的自适应性技术问题,又通过加入最小梯度正则化项使反演得到稳定且唯一的近似解,得到更聚焦的、边界更为清晰地质构造特征,提高了反演的稳定性和有效性。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案,下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本发明的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其他的附图。
图1是本发明实施例提供的叠前地震反演方法的流程示意图;
图2是本发明实施例提供的一种常规反演结果的示意图;
图3是本发明实施例提供的叠前地震反演方法的反演结果示意图;
图4是本发明实施例提供的叠前地震反演装置的结构示意图。
具体实施方式
下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅是本发明一部分实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都属于本发明保护的范围。
发明人发现:地震反演受频带宽度和噪声等的影响很难获得稳定且相对唯一的结果,需要在反演过程中加入一定的先验信息。这种先验信息有的是为了减少反演的不适定性,有的是体现了地质学家的先验知识,有的是为了提高反演算法对地震数据噪声的压制效果等。通常我们都假设地震噪声符合高斯分布,但实际地震记录中,噪声分布非常复杂,远非简单的高斯分布能够描述。尤其对于地震反演来说,噪声类型应该是多种分布的叠加。如Omid Karimi在Arild Buland研究基础上,讨论了地震噪声的分布情况,认为地震噪声服从偏高斯分布,发展了偏高斯分布的贝叶斯AVO反演(Karimi and Omre,2010)。Liu和Zhang等讨论了基于L1、L2混合范数的三参数非高斯反演方法(Liu et al,2015)。以上方法虽然较好地解决了地震数据的非高斯问题,但都不能根据地震噪声分布特征进行自适应反演。
由于地震数据的带限性和不完备性,地震反演问题涉及到的***方程组呈现出病态性,带来的后果就是观测数据中很小的噪声干扰,都会造成方程组的解出现很大的波动。因此在实际地震记录存在噪声的情况下,直接求解目标函数,得到的解往往是不稳定的、偏离真实情况的。解决地震反演问题病态性的常用方法就是在目标函数中加入正则化项,这样不仅能使反演得到稳定且唯一的近似解,还可以整合一些诸如测井、地震或地质构造解释等先验信息,从而提高反演的准确性。此外,正则化参数的选取也是正则化约束反演的一个重要环节。如果正则化参数选取得过小,则反演的不适定性得不到较好地改善,正则化方法无法发挥作用;如果正则化参数选取得过大,则会出现反演结果过拟合而偏离真实情况的现象。常用的正则化项有Tikhonov正则化,全变差正则化(Rudin and Osher 1992),广义交叉验证法(Wahba 1990)等。
由于发明人发现了上述技术问题,提出提出一种自适应混合范数的最小梯度正则化的叠前反演方法,该方法能自适应压制噪声,又能得到更聚焦的、边界更为清晰地质构造特征。
基于以上所述可知,本发明实施例提供方案的主要目的在于实现一种基于混合范数的自适应地震三参数叠前反演方法,提高反演的抗噪性,降低多解性,提高反演结果精度。该发明所描述方法能够有效克服不同类型地震噪声的影响,提高成像边界处的精度,可应用于实际油气勘探地震资料反演中,提高油气储层预测精度。
图1是本发明实施例中在本发明实施例提供的叠前地震反演方法的流程示意图,如图1所示,该方法包括:
步骤101:根据地震数据,确定地震中噪声类别和叠前地震正演模型;
步骤102:根据叠前地震正演模型和最小梯度支持正则化项,确定叠前反演残差模型;
步骤103:根据权重参数和叠前反演残差模型,构建自适应混合范数反演目标函数;
步骤104:根据自适应混合范数反演目标函数、叠前正演模型和地震中噪声类别,确定叠前地震反演结果。
本发明实施例提供技术方案的有益技术效果是:
本发明实施例的方案通过判断地震数据中噪声类别来选择范数从而可以自适应压制噪声,解决了叠前反演中噪声的自适应性技术问题,又通过加入最小梯度正则化项使反演得到稳定且唯一的近似解,得到更聚焦的、边界更为清晰地质构造特征,提高了反演的稳定性和有效性。
下面结合附图2至图3,对本发明实例涉及的各个步骤进行详细介绍如下。
一、首先介绍步骤101。
在一个实施例中,根据地震数据,确定地震中噪声类别,包括:
按如下公式计算峭度参数:
Figure BDA0002247859810000041
其中:χ为峭度参数,无量纲;Xi为地震信号值,无量纲;
Figure BDA0002247859810000042
为地震信号均值,无量纲;σ为地震信号标准差,无量纲;N为采样长度,无量纲;
在峭度参数为0时,确定地震中噪声类别为高斯分布;
在峭度参数为负值时,确定地震中噪声类别为亚髙斯分布;
在峭度参数为正值时,确定地震中噪声类别为超高斯分布。
具体实施时,已知地震数据中噪声主要存在高斯分布、超高斯分布和亚高斯分布噪声。我们可以采用峭度对噪声的高斯、亚高斯和超高斯类型进行区分。满足高斯分布的随机变量二阶及二阶以上的高阶统计量均为零。因此:
①当峭度应等于零时,为高斯分布;
②当峭度与零相差较大时,我们可以判断该随机变量偏离高斯分布,具体为:
峭度为负的信号称作亚髙斯信号;
峭度为正的信号称作超高斯信号。
在一个实施例中,根据地震数据,确定叠前地震正演模型,包括:
根据地震数据中叠后地震道子波与反射系数的卷积,确定地震道矩阵;
根据地震道矩阵,确定实际地震道反射率;
根据地震道矩阵和实际地震道反射率,确定合成地震记录;
根据合成地震记录,确定叠前正演模型。
在一个实施例中,根据地震数据中叠后地震道子波与反射系数的卷积,确定地震道矩阵,包括:
按如下公式确定地震数据中叠后地震道子波与反射系数的卷积表达式:
Figure BDA0002247859810000051
按如下公式确定地震道矩阵的表达式:
Figure BDA0002247859810000052
其中:tn代表地震道的第n个样本,无量纲;wk代表提取的地震子波的第k个项,无量纲;S是地震道;W是子波函数矩阵;Rp是P波反射率;D是正演差分矩阵函数矩阵;
Figure BDA0002247859810000053
Figure BDA0002247859810000054
为矩阵算子,无量纲;IP为P波阻抗的自然对数,无量纲;ln代表第n个样本反射系数,无量纲。
在一个实施例中,按如下公式根据地震道矩阵,确定实际地震道反射率:
Rpp(θ)=c1Rp+c2Rs+c3RD,
其中:c1=1+tan2θ,c2=-8α2tan2θ,c3=-0.5tan2θ+2α2sin2θ,α=Vs/Vp;α为纵横波速比,无量纲;Vs为横波速度,单位:m/s;Vp为纵波速度,单位:m/s;Rpp代表实际地震道反射率;Rs是S波反射率;RD是密度反射率;c1,c2,c3为方程系数,无量纲。在一个实施例中,按如下公式根据地震道矩阵和实际地震道反射率,确定合成地震记录:
Figure BDA0002247859810000061
其中:s(θ)为任意角度的地震道,无量纲;IS为S波阻抗的自然对数,无量纲;ID表示密度的自然对数,无量纲;W(θ)是任意角度子波函数,无量纲;θ是地震反射波入射角,单位是度。
在一个实施例中,如权利要求6中所述的方法,按如下公式根据合成地震记录,确定叠前正演模型:
s=Gx
Figure BDA0002247859810000062
其中:s为实际地震道,无量纲;x是输入模型参数矩阵,包括纵波阻抗自然对数、横波阻抗自然对数、密度自然对数,无量纲;G为正演算子矩阵,无量纲;s(θ)是s的一种表达形式,表示任意角度的地震道,以cdp形式表示,即地震数据共深度点道集矩阵函数,无量纲。
具体实施时,根据Robinson褶积模型,叠后地震道可以用子波与反射系数的卷积来表示,具体如下:
Figure BDA0002247859810000063
其中,tn代表地震道的第n个样本,无量纲;wk代表提取的地震子波的第k个项,无量纲;
Figure BDA0002247859810000071
Figure BDA0002247859810000072
为矩阵算子,无量纲;ln代表第n个样本反射系数,无量纲。
将(1)表示简写成矩阵地震道矩阵表达式:
S=WRp=W(DIp) (2)
其中,S是地震道;W是子波函数矩阵;Rp是P波反射率;D是正演差分矩阵函数矩阵;IP为P波阻抗的自然对数,无量纲。
将方程(2)扩展到叠前地震数据表达式:
Rpp(θ)=c1Rp+c2Rs+c3RD, (3)
其中,c1=1+tan2θ,c2=-8α2tan2θ,c3=-0.5tan2θ+2α2sin2θ,α=Vs/Vp,α为纵横波速比,无量纲;Vs为横波速度,单位:m/s;Vp为纵波速度,单位:m/s;Rpp代表实际地震道反射率;Rs是S波反射率;RD是密度反射率;c1,c2,c3为方程系数,无量纲。
为简单起见,我们假设α=0.5的恒定背景比。结合式(2)、式(3),于是对于任意角度的地震道,基于褶积模型的合成地震记录有如下表达式:
Figure BDA0002247859810000073
其中,s(θ)为任意角度的地震道,无量纲;IS为S波阻抗的自然对数,无量纲;ID表示密度的自然对数,无量纲;IP为P波阻抗的自然对数,无量纲。W(θ)是任意角度子波函数,无量纲;θ是地震反射波入射角,单位是度;W是子波函数矩阵;D是正演差分矩阵函数矩阵;c1,c2,c3为方程系数,无量纲。
将式(4)写成矩阵形式得到叠前正演模型,可表示为:
s=Gx
Figure BDA0002247859810000074
式(5)中,s为实际地震道,有不同的表达形式,无量纲;s(θ)是s的一种表达形式,表示任意角度的地震道,在实际生产中常以cdp形式表示,即地震数据共深度点道集(common depth point)矩阵函数,无量纲;G是正演算子矩阵,无量纲,D是差分矩阵函数,无量纲;W(θ)是任意角度子波函数,无量纲;c1,c2,c3为方程系数,无量纲;IS表示S波阻抗的自然对数,无量);ID表示密度的自然对数,无量纲;x是输入模型参数矩阵,包括纵波阻抗自然对数、横波阻抗自然对数、密度自然对数,无量纲;
二、其次介绍步骤102.
在一个实施例中,根据叠前地震正演模型和最小梯度支持正则化项,确定叠前反演残差模型,包括:
根据叠前地震正演模型,确定叠前地震反演目标函数;
根据最小梯度支持正则化项和叠前地震反演目标函数,确定叠前反演残差模型。
在一个实施例中,按如下公式根据叠前正演模型,确定叠前地震反演目标函数:
Figure BDA0002247859810000081
其中:J(X)为叠前地震反演目标函数,无量纲;s为实际地震道,无量纲;G为正演算子矩阵,无量纲;x为输入模型参数矩阵,无量纲。
在一个实施例中,按如下公式根据最小梯度支持正则化项和叠前地震反演目标函数,确定叠前反演残差模型:
Figure BDA0002247859810000082
其中:Φ(x)为反演目标函数,无量纲;矩阵D为一阶差分矩阵,无量纲;T为矩阵转置操作,无量纲;σ为阈值因子,无量纲;λ为正则化参数,无量纲。
具体实施时,全文所提到的“叠前地震反演目标函数”均指的是“叠前地震反演目标”作为该函数的定语,该函数的全称是“叠前地震反演目标函数”,而不是“反演”作为动词去“反演”目标函数;
全文所提到的“叠前反演残差模型”均指的是“叠前反演残差”作为该函数的定语,该函数的全称是“叠前反演残差模型”,而不是“反演”作为动词去“反演”残差模型;
全文所提到的“自适应混合范数反演目标函数”均指的是“自适应混合范数反演目标”作为该函数的定语,该函数的全称是“自适应混合范数反演目标函数演”,而不是“反演”作为动词去“反演”目标函数。
具体实施时,为了更准确地求解反演模型,将直接求解矩阵方程s=Gx,转化为最小化如下的目标函数:
Figure BDA0002247859810000083
解决地震反演问题病态性的常用方法就是在式(6)所示的目标函数中加入正则化项。J(X)为目标函数,无量纲;s是实际地震道,有不同的表达形式,无量纲;G是正演算子矩阵,无量纲;x是输入模型参数矩阵,无量纲;
①Tikhonov正则化建立
将Tikhonov泛函以拉格朗日乘数法的形式引入到模型的目标函数中:
Figure BDA0002247859810000091
上式的第一部分表示数据拟合项,第二部分表示正则化项。其中,λ为正则化参数,用来调节着数据拟合项与正则项的权重相对大小;矩阵D为平滑算子,通常取为一阶差分矩阵。Tikhonov正则化选取L2范数作为约束项,它可以有效减少原问题解的震荡性,从而给出了平滑的近似解。J(X)为目标函数,无量纲;s是实际地震道,有不同的表达形式,无量纲;G是正演算子矩阵,无量纲;x是输入模型参数矩阵,无量纲。
②全变差正则化建立
但是当地层波阻抗分层界面本身是不光滑的时候,Tikhonov正则化会导致界面被过分模糊,从而使反演结果的分辨率降低。提出了全变差正则化(Total variationregularization,简称TV正则化),全变差正则化泛函可表示为:
Figure BDA0002247859810000092
其中,矩阵D同样为一阶差分矩阵。与Tikhonov正则化不同之处在于,全变差正则化选取L1范数作为约束项。J(X)为目标函数,无量纲;s是实际地震道,有不同的表达形式,无量纲;G是正演算子矩阵,无量纲;x是输入模型参数矩阵,无量纲;λ为正则化参数,用来调节着数据拟合项与正则项的权重相对大小。
③最小梯度正则化建立
TV正则化只能根据整体梯度变化控制块状特征,即层内部和边界施加同样约束,显然不能合理突出层边界的效应。为了使反演结果可调,改进Tikhonov正则化、全变差正则化,定义最小梯度支持正则化函数,
最小梯度支持正则化加入到公式(4)中叠前反演模型后的形式如下:
Figure BDA0002247859810000093
其中:Φ(x)为自适应混合范数反演目标函数,无量纲;矩阵D为一阶差分矩阵,无量纲;矩阵G为正演算子矩阵,无量纲;x为输入模型参数矩阵,无量纲;T为矩阵转置操作,无量纲;σ为阈值因子,无量纲,用来控制解的锐利程度。在反演过程中那些大于阈值σ的模型扰动予以保留,而将小于阈值σ的模型扰动则进行了平滑处理。λ为正则化参数,无量纲,用来调节着数据拟合项与正则项的权重相对大小。
三、接着介绍步骤103.
在一个实施例中,按如下公式根据权重参数和叠前反演残差模型,构建自适应混合范数反演目标函数:
Figure BDA0002247859810000101
其中:Φ(x)为自适应混合范数反演目标函数,无量纲;γ(nb)为权重参数,无量纲;nb表示背景噪声,无量纲;s为实际地震道,无量纲;矩阵G为正演算子矩阵,无量纲;x为输入模型参数矩阵,无量纲
具体实施时,构建自适应混合范数反演目标函数。
自适应混合范数反演目标函数可定义为:
Figure BDA0002247859810000102
其中权重参数γ(nb),无量纲,控制L2范数和L4范数的权值大小;nb表示背景噪声,无量纲;
四、最后介绍步骤104。
在一个实施例中,根据自适应混合范数反演目标函数、叠前正演模型以及地震中噪声类别,确定叠前地震反演结果,包括:
在地震中噪声类别为高斯分布或超高斯分布时,权重参数γ(nb)≈0,自适应混合范数反演目标函数采用L2范数,得出此时的自适应混合范数反演目标函数的结果;
当地震中噪声类别为亚高斯噪声时,权重参数γ(nb)≈1,自适应混合范数反演目标函数采用L4范数,得出此时的自适应混合范数反演目标函数的结果;
通过计算地震数据与叠前地震正演模型的数据的均方根误差,当此误差与自适应混合范数反演目标函数的结果相差最小时,此时的自适应混合范数反演目标函数的结果即为反演数据结果。
具体实施时,可以看出:
γ(nb)=0时,自适应混合范数反演目标函数退化为L2范数;当γ(nb)≈1时,自适应混合范数反演目标函数退化为L4范数;自适应准则应满足以下条件:
自适应权重参数应是反演误差的函数,从而对权重进行自适应调节;
当噪声为高斯或是超高斯噪声时,权重参数γ(nb)≈0,使用L2范数;反之,当噪声为亚高斯噪声时,权重参数γ(nb)≈1,使用L4范数。
具体实施时,还包括对正则化参数优化。如果正则化参数选取得过小,则反演的不适定性得不到较好地改善,正则化方法无法发挥作用;如果正则化参数选取得过大,则会出现反演结果过拟合而偏离真实情况的现象,需要根据具体情况进行调整测试。
具体实施时,根据公式(5)~(10)迭代更新模型。从初始模型出发,利用正演合成记录与叠前地震记录的残差(实际叠前地震记录与正演合成地震记录的均方根误差)计算出的模型修正量来不断更新模型,使模型正演合成记录与实际地震资料达到最佳吻合,最后输出反演数据结果。
具体实施时,图2至图3为常规反演和自适应混合范数最小梯度支持正则化的反演方法对比结果,可以发现,常规反演不能有效克服噪声造成的抖动,并且地层轮廓不能清晰地体现出来,而自适应混合范数反演结果(图3)信噪比和分辨率要明显高于常规反演方法(图2)。
本发明实施例的方案通过判断地震数据中噪声类别来选择范数从而可以自适应压制噪声,解决了叠前反演中噪声的自适应性技术问题,又通过加入最小梯度正则化项使反演得到稳定且唯一的近似解,得到更聚焦的、边界更为清晰地质构造特征,提高了反演的稳定性和有效性。
本发明是一种基于自适应混合范数的最小梯度支持正则化叠前反演技术,通过峭度构建高斯噪声和非高斯噪声的自适应权重,引入最小梯度支持正则化使反演结果边缘更清晰,轮廓更明显,为油气藏的精确识别奠定基础,指导后续油气检测及钻井安全。
本发明适用于在油气勘探前期的探井井位选择及后期的油气识别。这一方法为油气藏的精确识别提供了依据。
基于同一发明构思,本发明实施例中还提供了一种叠前地震反演装置,如下面的实施例所述。由于叠前地震反演装置解决问题的原理与叠前地震反演方法相似,因此叠前地震反演装置的实施可以参见叠前地震反演方法的实施,重复之处不再赘述。以下所使用的,术语“单元”或者“模块”可以实现预定功能的软件和/或硬件的组合。尽管以下实施例所描述的装置较佳地以软件来实现,但是硬件,或者软件和硬件的组合的实现也是可能并被构想的。
图4是本发明实施例的叠前地震反演装置的一种结构框图,如图4所示,包括:
噪声分类模块01,用于根据地震数据,确定地震中噪声类别;
叠前地震正演模型构建模块02,用于根据地震数据,确定叠前正演模型;
叠前反演残差模型构建模块03,用于根据叠前地震正演模型和最小梯度支持正则化项,确定叠前反演残差模型;
反演目标函数构建模块04,用于根据权重参数和叠前反演残差模型,构建自适应混合范数反演目标函数;
反演数据获取模块05,用于根据自适应混合范数反演目标函数、叠前正演模型和地震中噪声类别,确定叠前地震反演结果。
本发明实施例还提供了一种计算机设备,包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序,所述处理器执行所述计算机程序时实现上述所述方法。
本发明实施例还提供了一种计算机可读存储介质,所述计算机可读存储介质存储有执行上述所述方法的计算机程序。
本领域内的技术人员应明白,本发明的实施例可提供为方法、***、或计算机程序产品。因此,本发明可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且,本发明可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等)上实施的计算机程序产品的形式。
本发明是参照根据本发明实施例的方法、设备(***)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器,使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。
这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中,使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品,该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。
这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上,使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理,从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。
以上所述仅为本发明的优选实施例而已,并不用于限制本发明,对于本领域的技术人员来说,本发明实施例可以有各种更改和变化。凡在本发明的精神和原则之内,所作的任何修改、等同替换、改进等,均应包含在本发明的保护范围之内。

Claims (20)

1.一种叠前地震反演方法,其特征在于,包括:
根据地震数据,确定地震中噪声类别和叠前地震正演模型;
根据叠前地震正演模型和最小梯度支持正则化项,确定叠前反演残差模型;
根据权重参数和叠前反演残差模型,构建自适应混合范数反演目标函数;
根据自适应混合范数反演目标函数、叠前正演模型和地震中噪声类别,确定叠前地震反演结果。
2.如权利要求1中所述的方法,其特征在于,根据地震数据,确定地震中噪声类别,包括:
按如下公式计算峭度参数:
Figure FDA0002247859800000011
其中:χ为峭度参数,无量纲;Xi为地震信号值,无量纲;
Figure FDA0002247859800000012
为地震信号均值,无量纲;σ为地震信号标准差,无量纲;N为采样长度,无量纲;
在峭度参数为0时,确定地震中噪声类别为高斯分布;
在峭度参数为负值时,确定地震中噪声类别为亚髙斯分布;
在峭度参数为正值时,确定地震中噪声类别为超高斯分布。
3.如权利要求1中所述的方法,其特征在于,根据地震数据,确定叠前地震正演模型,包括:
根据地震数据中叠后地震道子波与反射系数的卷积,确定地震道矩阵;
根据地震道矩阵,确定实际地震道反射率;
根据地震道矩阵和实际地震道反射率,确定合成地震记录;
根据合成地震记录,确定叠前正演模型。
4.如权利要求3中所述的方法,其特征在于,根据地震数据中叠后地震道子波与反射系数的卷积,确定地震道矩阵,包括:
按如下公式确定地震数据中叠后地震道子波与反射系数的卷积表达式:
Figure FDA0002247859800000021
按如下公式确定地震道矩阵的表达式:
S=WRp=W(DIp)
其中:tn代表地震道的第n个样本,无量纲;wk代表提取的地震子波的第k个项,无量纲;S是地震道;W是子波函数矩阵;Rp是P波反射率;D是正演差分矩阵函数矩阵;
Figure FDA0002247859800000022
Figure FDA0002247859800000023
为矩阵算子,无量纲;IP为P波阻抗的自然对数,无量纲;ln代表第n个样本反射系数,无量纲。
5.如权利要求4中所述的方法,其特征在于,按如下公式根据地震道矩阵,确定实际地震道反射率:
Rpp(θ)=c1Rp+c2Rs+c3RD,
其中:c1=1+tan2θ,c2=-8α2tan2θ,c3=-0.5tan2θ+2α2sin2θ,α=Vs/Vp;α为纵横波速比,无量纲;Vs为横波速度,单位:m/s;Vp为纵波速度,单位:m/s;Rpp代表实际地震道反射率;Rs是S波反射率;RD是密度反射率;c1,c2,c3为方程系数,无量纲。
6.如权利要求5中所述的方法,其特征在于,按如下公式根据地震道矩阵和实际地震道反射率,确定合成地震记录:
Figure FDA0002247859800000024
其中:s(θ)为任意角度的地震道,无量纲;IS为S波阻抗的自然对数,无量纲;ID表示密度的自然对数,无量纲;W(θ)是任意角度子波函数,无量纲;θ是地震反射波入射角,单位是度。
7.如权利要求6中所述的方法,其特征在于,按如下公式根据合成地震记录,确定叠前正演模型:
s=Gx
Figure FDA0002247859800000031
其中:s为实际地震道,无量纲;x是输入模型参数矩阵,包括纵波阻抗自然对数、横波阻抗自然对数、密度自然对数,无量纲;G为正演算子矩阵,无量纲;s(θ)是s的一种表达形式,表示任意角度的地震道,以cdp形式表示,即地震数据共深度点道集矩阵函数,无量纲。
8.如权利要求1中所述的方法,其特征在于,根据叠前地震正演模型和最小梯度支持正则化项,确定叠前反演残差模型,包括:
根据叠前地震正演模型,确定叠前地震反演目标函数;
根据最小梯度支持正则化项和叠前地震反演目标函数,确定叠前反演残差模型。
9.如权利要求8中所述的方法,其特征在于,按如下公式根据叠前正演模型,确定叠前地震反演目标函数:
Figure FDA0002247859800000032
其中:J(X)为叠前地震反演目标函数,无量纲;s为实际地震道,无量纲;G为正演算子矩阵,无量纲;x为输入模型参数矩阵,无量纲。
10.如权利要求9中所述的方法,其特征在于,按如下公式根据最小梯度支持正则化项和叠前地震反演目标函数,确定叠前反演残差模型:
Figure FDA0002247859800000033
其中:Φ(x)为自适应混合范数反演目标函数,无量纲;矩阵D为一阶差分矩阵,无量纲;T为矩阵转置操作,无量纲;σ为阈值因子,无量纲;λ为正则化参数,无量纲。
11.如权利要求1中所述的方法,其特征在于,按如下公式根据权重参数和叠前反演残差模型,构建自适应混合范数反演目标函数:
Figure FDA0002247859800000034
其中:Φ(x)为自适应混合范数反演目标函数,无量纲;γ(nb)为权重参数,无量纲;nb表示背景噪声,无量纲;s为实际地震道,无量纲;矩阵G为正演算子矩阵,无量纲;x为输入模型参数矩阵,无量纲。
12.如权利要求1中所述的方法,其特征在于,根据自适应混合范数反演目标函数、叠前正演模型以及地震中噪声类别,确定叠前地震反演结果,包括:
在地震中噪声类别为高斯分布或超高斯分布时,权重参数γ(nb)≈0,自适应混合范数反演目标函数采用L2范数,得出此时的自适应混合范数反演目标函数的结果;
当地震中噪声类别为亚高斯噪声时,权重参数γ(nb)≈1,自适应混合范数反演目标函数采用L4范数,得出此时的自适应混合范数反演目标函数的结果;
通过计算地震数据与叠前地震正演模型的数据的均方根误差,当此误差与自适应混合范数反演目标函数的结果相差最小时,此时的自适应混合范数反演目标函数的结果即为反演数据结果。
13.一种叠前地震反演装置,其特征在于,包括:
噪声分类模块,用于根据地震数据,确定地震中噪声类别;
叠前地震正演模型构建模块,用于根据地震数据,确定叠前正演模型;
叠前反演残差模型构建模块,用于根据叠前地震正演模型和最小梯度支持正则化项,确定叠前反演残差模型;
反演目标函数构建模块,用于根据权重参数和叠前反演残差模型,构建自适应混合范数反演目标函数;
反演数据获取模块,用于根据自适应混合范数反演目标函数、叠前正演模型和地震中噪声类别,确定叠前地震反演结果。
14.如权利要求13中所述的装置,其特征在于,噪声分类模块具体用于:
按如下公式计算峭度参数:
Figure FDA0002247859800000041
其中:χ为峭度参数,无量纲;Xi为信号值,无量纲;
Figure FDA0002247859800000042
为信号均值,无量纲;σ为信号标准差,无量纲;N为采样长度,无量纲;
在峭度参数为0时,确定地震中噪声类别为高斯分布;
在峭度参数为负值时,确定地震中噪声类别为亚髙斯分布;
在峭度参数为正值时,确定地震中噪声类别为超高斯分布。
15.如权利要求13中所述的装置,其特征在于,正演模型构建模块具体用于:
根据地震数据中叠后地震道子波与反射系数的卷积表达式,确定地震道矩阵;
根据地震道矩阵,确定实际地震道反射率;
根据地震道矩阵和实际地震道反射率,确定合成地震记录;
根据合成地震记录,确定叠前正演模型。
16.如权利要求13中所述的装置,其特征在于,叠前反演残差模型构建模块具体用于:
根据叠前地震正演模型,确定叠前地震反演目标函数;
根据最小梯度支持正则化项和叠前地震反演目标函数,确定叠前反演残差模型。
17.如权利要求13中所述的装置,其特征在于,反演目标函数构建模块具体用于按以下公式根据权重参数和叠前反演残差模型,构建自适应混合范数反演目标函数:
Figure FDA0002247859800000051
其中:Φ(x)为自适应混合范数反演目标函数,无量纲;γ(nb)为权重参数,无量纲;nb表示背景噪声,无量纲;s为实际地震道,无量纲;矩阵G为正演算子矩阵,无量纲;x为输入模型参数矩阵,无量纲。
18.如权利要求13中所述的装置,其特征在于,反演数据获取模块具体用于:
当地震中噪声类别为高斯分布或超高斯分布时,权重参数γ(nb)≈0,自适应混合范数反演目标函数采用L2范数,得出此时的自适应混合范数反演目标函数的结果;
当地震中噪声类别为亚高斯噪声时,权重参数γ(nb)≈1,自适应混合范数反演目标函数采用L4范数,得出此时的自适应混合范数反演目标函数的结果;
通过计算地震数据与叠前地震正演模型的数据的均方根误差,当此误差与自适应混合范数反演目标函数的结果相差最小时,此时的自适应混合范数反演目标函数的结果即为反演数据结果。
19.一种计算机设备,包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序,其特征在于,所述处理器执行所述计算机程序时实现权利要求1至12任一所述方法。
20.一种计算机可读存储介质,其特征在于,所述计算机可读存储介质存储有执行权利要求1至12任一所述方法的计算机程序。
CN201911023165.XA 2019-10-25 2019-10-25 叠前地震反演方法及装置 Pending CN112711065A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201911023165.XA CN112711065A (zh) 2019-10-25 2019-10-25 叠前地震反演方法及装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201911023165.XA CN112711065A (zh) 2019-10-25 2019-10-25 叠前地震反演方法及装置

Publications (1)

Publication Number Publication Date
CN112711065A true CN112711065A (zh) 2021-04-27

Family

ID=75540902

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201911023165.XA Pending CN112711065A (zh) 2019-10-25 2019-10-25 叠前地震反演方法及装置

Country Status (1)

Country Link
CN (1) CN112711065A (zh)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114460654A (zh) * 2022-02-22 2022-05-10 成都理工大学 基于l1l2混合范数的半航空瞬变电磁数据反演方法及装置
CN116381793A (zh) * 2023-04-13 2023-07-04 中国矿业大学(北京) 结构tv正则化联合道间差异约束的叠前反演方法及装置
CN116381787A (zh) * 2023-04-12 2023-07-04 中国矿业大学(北京) 叠前反演方法、装置、电子设备及介质

Citations (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103293552A (zh) * 2013-05-24 2013-09-11 中国石油天然气集团公司 一种叠前地震资料的反演方法及***
US20140350903A1 (en) * 2011-09-01 2014-11-27 Geofoft Inc. Methods And Systems For The Inversion Of Magnetic Data From Remnant And Induced Sources In Geophysical Exploration
CN104237937A (zh) * 2014-07-28 2014-12-24 中国石油化工股份有限公司 叠前地震反演方法及其***
CN104977618A (zh) * 2014-04-09 2015-10-14 中国石油集团东方地球物理勘探有限责任公司 一种评价页岩气储层及寻找甜点区的方法
CN106501872A (zh) * 2016-09-26 2017-03-15 中国石油天然气股份有限公司 一种裂缝储层地应力特征的计算方法及装置
CN106646625A (zh) * 2016-09-27 2017-05-10 中国科学院电子学研究所 一种锐边界模型的瞬变电磁反演方法
CN106842313A (zh) * 2015-12-04 2017-06-13 中国石油化工股份有限公司 基于方位叠前地震数据的各向异性参数反演方法
CN107831543A (zh) * 2017-09-29 2018-03-23 中国石油化工股份有限公司 叠前地震反演方法及***
CN109143356A (zh) * 2018-08-29 2019-01-04 电子科技大学 一种自适应混合范数字典学习地震波阻抗反演方法

Patent Citations (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20140350903A1 (en) * 2011-09-01 2014-11-27 Geofoft Inc. Methods And Systems For The Inversion Of Magnetic Data From Remnant And Induced Sources In Geophysical Exploration
CN103293552A (zh) * 2013-05-24 2013-09-11 中国石油天然气集团公司 一种叠前地震资料的反演方法及***
CN104977618A (zh) * 2014-04-09 2015-10-14 中国石油集团东方地球物理勘探有限责任公司 一种评价页岩气储层及寻找甜点区的方法
CN104237937A (zh) * 2014-07-28 2014-12-24 中国石油化工股份有限公司 叠前地震反演方法及其***
CN106842313A (zh) * 2015-12-04 2017-06-13 中国石油化工股份有限公司 基于方位叠前地震数据的各向异性参数反演方法
CN106501872A (zh) * 2016-09-26 2017-03-15 中国石油天然气股份有限公司 一种裂缝储层地应力特征的计算方法及装置
CN106646625A (zh) * 2016-09-27 2017-05-10 中国科学院电子学研究所 一种锐边界模型的瞬变电磁反演方法
CN107831543A (zh) * 2017-09-29 2018-03-23 中国石油化工股份有限公司 叠前地震反演方法及***
CN109143356A (zh) * 2018-08-29 2019-01-04 电子科技大学 一种自适应混合范数字典学习地震波阻抗反演方法

Non-Patent Citations (9)

* Cited by examiner, † Cited by third party
Title
YAOJUN WANG ET AL: "《Sharp and laterally constrained multitrace impedance inversion based on blocky coordinate descent》", ACTA GEOPHYSICA, 6 June 2018 (2018-06-06), pages 2 - 7 *
张志勇等: "《基于最小梯度支撑的2.5D井地电位法正则化聚焦反演》", 中国有色金属学报, vol. 25, no. 11, 30 November 2015 (2015-11-30), pages 3182 - 3189 *
杨亚华等: "《叠前同时反演在灰质发育区识别储层及流体的应用研究》", 地球物理学进展, vol. 32, no. 1, pages 332 - 338 *
梁立锋等: "《利用叠前密度反演预测岩性的应用研究》", 海洋石油, vol. 31, no. 2, pages 53 - 58 *
王兴建等: "《纵横波联合反演在苏里格含气性检测中的应用》", 天然气工业, vol. 28, no. 10, pages 44 - 45 *
王峣钧等: "《基于自适应混合范数的快速多道波阻抗反演》", GPS/SEG北京2018国际地球物理会议暨展览电子论文集, pages 1068 - 1069 *
程三: "《基于构造约束的地震信号多道反演方法研究》", 中国优秀硕士学位论文全文数据库 基础科学辑, no. 2018, 15 September 2018 (2018-09-15), pages 12 - 17 *
程三: "《基于构造约束的地震信号多道反演方法研究》", 中国优秀硕士学位论文全文数据库 基础科学辑, pages 12 - 13 *
霍凤斌等: "《叠前同时反演在储层识别和流体预测中的应用――以A油气田B构造C组为例》", 油气地球物理, vol. 12, no. 1, pages 9 - 13 *

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114460654A (zh) * 2022-02-22 2022-05-10 成都理工大学 基于l1l2混合范数的半航空瞬变电磁数据反演方法及装置
CN116381787A (zh) * 2023-04-12 2023-07-04 中国矿业大学(北京) 叠前反演方法、装置、电子设备及介质
CN116381787B (zh) * 2023-04-12 2023-09-01 中国矿业大学(北京) 叠前反演方法、装置、电子设备及介质
CN116381793A (zh) * 2023-04-13 2023-07-04 中国矿业大学(北京) 结构tv正则化联合道间差异约束的叠前反演方法及装置
CN116381793B (zh) * 2023-04-13 2023-09-01 中国矿业大学(北京) 结构tv正则化联合道间差异约束的叠前反演方法及装置

Similar Documents

Publication Publication Date Title
US8923093B2 (en) Determining the quality of a seismic inversion
US10310113B2 (en) Q-compensated full wavefield inversion
US9470811B2 (en) Creating a high resolution velocity model using seismic tomography and impedance inversion
US11294087B2 (en) Directional Q compensation with sparsity constraints and preconditioning
CN108535775B (zh) 非平稳地震资料声波阻抗反演方法及装置
CN112711065A (zh) 叠前地震反演方法及装置
US20150293245A1 (en) Method and device for the generation and application of anisotropic elastic parameters in horizontal transverse isotropic (hti) media
US11175421B2 (en) Device and method for mitigating cycle-skipping in full waveform inversion
US9322943B2 (en) Method and apparatus for pre-stack deghosting of seismic data
CN108693555B (zh) 智能化时变盲反褶积宽频处理方法及装置
US10310117B2 (en) Efficient seismic attribute gather generation with data synthesis and expectation method
AU2014280832B2 (en) Seismic data spectrum restoring and broadening
CN110895348B (zh) 一种地震弹性阻抗低频信息提取方法、***及存储介质
WO2022232572A1 (en) Method and system for high resolution least-squares reverse time migration
US20180017692A1 (en) Device and method for estimating pre-stack wavelet model from seismic gathers
US20160320507A1 (en) Time lapse seismic data processing
CN110687597B (zh) 一种基于联合字典的波阻抗反演方法
Li et al. Robust Q-compensated multidimensional impedance inversion using seislet-domain shaping regularization
CN114114421B (zh) 基于深度学习的导向自学习地震数据去噪方法及装置
CN112363222A (zh) 叠后自适应宽带约束波阻抗反演方法及装置
US9784869B2 (en) Noise models by selection of transform coefficients
CN113495296B (zh) 层析静校正量的确定方法、装置、设备及可读存储介质
CN117388921B (zh) 弹性参数的叠前反演方法、装置和电子设备
CN113985480A (zh) 一种基于角度校正的avo反演方法及装置
Xin et al. A blind seismic inversion based balanced combination of L2 and total variation regularizations

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