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

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

Info

Publication number
CN112014875B
CN112014875B CN201910467123.9A CN201910467123A CN112014875B CN 112014875 B CN112014875 B CN 112014875B CN 201910467123 A CN201910467123 A CN 201910467123A CN 112014875 B CN112014875 B CN 112014875B
Authority
CN
China
Prior art keywords
model
inversion
energy error
iteration
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
CN201910467123.9A
Other languages
English (en)
Other versions
CN112014875A (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.)
Beijing Sunshine Geo Tech Co ltd
Original Assignee
Beijing Sunshine Geo Tech Co ltd
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 Beijing Sunshine Geo Tech Co ltd filed Critical Beijing Sunshine Geo Tech Co ltd
Priority to CN201910467123.9A priority Critical patent/CN112014875B/zh
Publication of CN112014875A publication Critical patent/CN112014875A/zh
Application granted granted Critical
Publication of CN112014875B publication Critical patent/CN112014875B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

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

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

叠前地震反演方法及装置
技术领域
本发明涉及地质勘探技术领域,尤其涉及一种叠前地震反演方法及装置。
背景技术
本部分旨在为权利要求书中陈述的本发明实施例提供背景或上下文。此处的描述不因为包括在本部分中就承认是现有技术。
储层物性特征参数反映了油气的空间分布和油气的储量估计,直接影响油气开采方案的设计和采收率,一直以来都是油气勘探开发研究的一个重点。随着油气田勘探开发难度的增大和技术的发展,近年来储层预测及流体检测地球物理方法与技术的研究取得了显著进展。地震反演是储层预测及流体检测技术系列之一,其是利用野外观测的地震数据,确定地质目标的物理参数。
传统的叠前地震反演是基于贝叶斯原理,其假设数据为高斯分布。这种假设只考虑大部分数据分布区域,不具有稀疏性,不能提高反演的分辨率;同时传统叠前地震反演收敛速度较慢。
因此,现有的叠前地震反演存在反演效率低、反演分辨率低的问题。
发明内容
本发明实施例提供一种叠前地震反演方法,用以提高反演效率和反演分辨率,该方法包括:
步骤101:采集目标储层的叠前地震数据;
步骤102:通过岩石物理建模确定叠前地震数据的先验模型和反演初始模型,将反演初始模型作为当前反演迭代模型;先验模型是反映目标储层的弹性参数的模型;反演初始模型是反映目标储层的初始弹性参数的模型;
步骤103:根据目标储层的地震角度道集及当前反演迭代模型确定能量误差;能量误差满足高斯概率分布;
若能量误差不小于预设能量误差,则执行步骤104;若能量误差小于预设能量误差,则执行步骤105;
步骤104:根据叠前地震数据的先验模型和当前反演迭代模型确定更新后的反演迭代模型,更新后的反演迭代模型满足柯西概率分布;将更新后的反演迭代模型作为当前反演迭代模型,返回执行步骤103进行迭代;
步骤105:将能量误差小于预设能量误差时的当前反演迭代模型作为叠前地震数据的反演结果。
本发明实施例还提供一种叠前地震反演装置,用以提高反演效率和反演分辨率,该装置包括:
数据采集模块,用于采集目标储层的叠前地震数据;
模型建立模块,用于通过岩石物理建模确定叠前地震数据的先验模型和反演初始模型,将反演初始模型作为当前反演迭代模型;先验模型是反映目标储层的弹性参数的模型;反演初始模型是反映目标储层的初始弹性参数的模型;
能量误差确定模块,用于根据目标储层的地震角度道集及当前反演迭代模型确定能量误差;
更新模块,用于若能量误差不小于预设能量误差,根据叠前地震数据的先验模型和当前反演迭代模型确定更新后的反演迭代模型,更新后的反演迭代模型满足柯西概率分布;将更新后的反演迭代模型作为当前反演迭代模型;
反演结果获取模块,用于将能量误差小于预设能量误差时的当前反演迭代模型作为叠前地震数据的反演结果。
本发明实施例还提供一种计算机设备,包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序,所述处理器执行所述计算机程序时实现上述叠前地震反演方法。
本发明实施例还提供一种计算机可读存储介质,所述计算机可读存储介质存储有执行上述叠前地震反演方法的计算机程序。
本发明实施例中,能量误差满足高斯概率分布,更新后的反演迭代模型满足柯西概率分布,同时在能量误差不小于预设能量误差时,将更新后的反演迭代模型作为当前反演迭代模型,迭代执行根据目标储层的地震角度道集及当前反演迭代模型确定能量误差,直至能量误差小于预设能量误差,因此,不仅能够提高叠前地震数据的反演效率,还能够提高叠前地震数据的反演分辨率。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案,下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本发明的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其他的附图。在附图中:
图1为本发明实施例提供的叠前地震反演方法的实现流程图;
图2为本发明实施例提供的叠前地震反演方法中步骤103的实现流程图;
图3为本发明实施例提供的叠前地震反演方法中步骤104的实现流程图;
图4为本发明实施例提供的叠前地震反演方法中步骤302的实现流程图;
图5为本发明实施例提供的叠前地震反演装置的功能模块图;
图6为本发明实施例提供的叠前地震反演装置中能量误差确定模块503的结构框图;
图7为本发明实施例提供的叠前地震反演装置中更新模块504的结构框图;
图8为本发明实施例提供的叠前地震反演装置中柯西概率分布确定单元702的结构框图。
具体实施方式
为使本发明实施例的目的、技术方案和优点更加清楚明白,下面结合附图对本发明实施例做进一步详细说明。在此,本发明的示意性实施例及其说明用于解释本发明,但并不作为对本发明的限定。
虽然本发明提供了如下述实施例或附图所示的方法操作步骤或装置结构,但基于常规或者无需创造性的劳动在所述方法或装置中可以包括更多或者更少的操作步骤或模块单元。在逻辑性上不存在必要因果关系的步骤或结构中,这些步骤的执行顺序或装置的模块结构不限于本发明实施例或附图所示的执行顺序或模块结构。所述的方法或模块结构的在实际中的装置或终端产品应用时,可以按照实施例或者附图所示的方法或模块结构进行顺序执行或者并行执行。
针对现有叠前地震反演存在的反演效率低、反演分辨率低的缺陷,本发明的申请人提出了一种叠前地震反演方法及装置,其在能量误差不小于预设能量误差时,根据叠前地震数据的先验模型和当前反演迭代模型确定更新后的反演迭代模型,将更新后的反演迭代模型作为当前反演迭代模型,迭代执行根据目标储层的地震角度道集及当前反演迭代模型确定能量误差,直至能量误差小于预设能量误差;同时能量误差满足高斯概率分布,更新后的反演迭代模型满足柯西概率分布,达到了提高叠前地震数据的反演效率和反演分辨率的目的。
图1示出了本发明实施例提供的叠前地震反演方法的实现流程,为便于描述,仅示出了与本发明实施例相关的部分,详述如下:
如图1所示,叠前地震反演方法,其包括:
步骤101,采集目标储层的叠前地震数据;
步骤102,通过岩石物理建模确定叠前地震数据的先验模型和反演初始模型,将反演初始模型作为当前反演迭代模型;先验模型是反映目标储层的弹性参数的模型;反演初始模型是反映目标储层的初始弹性参数的模型;
步骤103:根据目标储层的地震角度道集及当前反演迭代模型确定能量误差;
若能量误差不小于预设能量误差,则执行步骤104;若能量误差小于预设能量误差,则执行步骤105;
步骤104:根据叠前地震数据的先验模型和当前反演迭代模型确定更新后的反演迭代模型,更新后的反演迭代模型满足柯西概率分布;将更新后的反演迭代模型作为当前反演迭代模型,返回执行步骤103进行迭代;
步骤105:将能量误差小于预设能量误差时的当前反演迭代模型作为叠前地震数据的反演结果。
叠前地震数据,是指对地震数据处理中的偏移来说的,在地下介质的产状不是水平时,反射地震的同相轴会发生偏移,不能反映该处地下介质的真实产状,这时要求在地震资料处理时进行偏移归位,叠加之前偏移就称为叠前。叠前地震数据可以采用dr表示,在本发明的一实施例中,叠前地震数据dr至少包括共角度域叠前地震数据。
目标储层是所要研究的储层区域,例如可以是地层层位在1000米至2000米层位范围内的叠前地震数据,也可以是地层层位在2500米至4000米层位范围的叠前地震数据。本领域技术人员可以理解的是,还可以是上述层位范围之外的其他地层层位数据,本发明实施例不再详细赘述。
地球物理方法在页岩气勘探开发过程中起到了重要的技术支持作用,它在研究储层的复杂特征方面具有独特的优势。岩石物理建模是岩石物理研究的一个重要方面。岩石物理建模,是指为了模拟介质的弹性性质所提出的一种等效模型。
在采集到叠前地震数据dr后,可以基于岩石物理建模获得叠前地震数据dr的先验模型mp(x,t)=(vp,vsp)和反演初始模型mp0(x,t)=(vp0,vs0p0)。其中,x表示偏移距,t表示时间。先验模型mp(x,t)=(vp,vsp)和反演初始模型mp0(x,t)=(vp0,vs0p0)的大小可以用M表示。先验模型mp(x,t)反映目标储层的弹性参数的模型。其中,弹性参数至少包括纵波速度vp、横波速度vs及密度ρp,即叠前地震数据的先验模型mp(x,t)是反映目标储层的纵波速度vp、横波速度vs及密度ρp等弹性参数的模型。同理,反演初始模型mp0(x,t)是反映目标储层的初始弹性参数的模型,是先验模型的低频部分可以通过低通滤波获得,一般来讲反演初始模型mp0(x,t)是频率为10赫兹左右的低频弹性参数模型。其中,初始弹性参数至少包括初始纵波速度vp0、初始横波速度vs0及初始密度ρp0,即叠前地震数据的反演初始模型mp0(x,t)是反映目标储层的初始纵波速度vp0、初始横波速度vs0及初始密度ρp0等初始弹性参数的模型。
在基于岩石物理建模获得叠前地震数据dr的反演初始模型后,在初次迭代时将反演初始模型作为当前反演迭代模型mpi(x,t)=(vpi,vsipi),其中,i表示迭代次数。当前反演迭代模型mpi(x,t)是反映目标储层当前迭代反演时的纵波速度vpi、横波速度vsi及密度ρpi等弹性参数的模型。这里我们把包括反演初始模型在内的每次迭代的反演模型称为反演迭代模型,当前反演迭代模型为进行当前迭代的反演模型。
在将反演初始模型确定为当前反演迭代模型后,根据目标储层的地震角度道集角度道集及根据当前反演迭代模型正演的角度道集进行能量误差分析。能量误差是反映反演分辨率的指标,可以理解的是,能量误差越小反演分辨率越高,能量误差越大反演分辨率越低。
其中,预设能量误差为预先设定的能量误差,可以根据实际情况和具体需求预先设定。例如,可以预先设定该预设能量误差为10-2,或者0.8×10-2,又或者是1.2×10-2。本领域技术人员可以理解的是,还可以预先设定预设能量误差为除上述数值之外的其他数值,例如0.9×10-2,又或者是1.1×10-2;本发明实施例不再详细赘述。
在确定能量误差后判断能量误差与预设能量误差的关系,若能量误差不小于预设能量误差,则说明当前的反演分辨率不符合设定的分辨率的要求,此时则可以根据叠前地震数据dr的先验模型和当前反演迭代模型,对正在进行反演的反演模型(即当前反演迭代模型)进行更新,确定更新后的反演迭代模型,并将更新后的反演迭代模型作为当前反演迭代模型,返回执行步骤103进行迭代。更新后的反演迭代模型,与反演初始模型类似,都是基于纵波速度、横波速度及密度等弹性参数的模型。
在本发明实施例中,能量误差满足高斯分布,更新后的反演迭代模型满足柯西概率分布。高斯分布(又称为正态分布),高斯概率分布是指通过概率密度函数来定义高斯分布。在本发明实施例中,柯西概率分布是指与更新后的反演迭代模型相关的,反映目标储层的纵波速度、横波速度以及密度等弹性参数的概率分布。柯西概率分布,是一个数学期望不存在的连续性概率分布,当随机变量X满足它的概率密度函数时,称变量X满足柯西概率分布。在本发明实施例中,柯西概率分布是指与更新后的反演迭代模型相关的,反映目标储层的纵波速度、横波速度以及密度等弹性参数的概率分布。
在将更新后的反演迭代模型作为当前反演迭代模型迭代执行步骤103,即迭代执行根据目标储层的地震角度道集及当前反演迭代模型(在步骤104中更新的反演迭代模型)确定能量误差后,进一步继续判断能量误差与预设能量误差的大小关系,在满足能量误差不小于预设能量误差时继续执行步骤104更新反演迭代模型,并将更新后的反演迭代模型迭代执行步骤103得到能量误差,继续进一步判断能量误差与预设能量误差的大小关系,直至满足能量误差小于预设能量误差,即当前的反演满足了反演分辨率的要求,此时执行步骤105,即将能量误差小于预设能量误差时的当前反演迭代模型作为叠前地震数据dr的反演结果。
在本发明实施例中,能量误差满足高斯分布,更新后的反演迭代模型满足柯西概率分布,同时在能量误差不小于预设能量误差时,将更新后的反演迭代模型作为当前反演迭代模型,迭代执行根据目标储层的地震角度道集及当前反演迭代模型确定能量误差,直至能量误差小于预设能量误差,因此,不仅能够提高叠前地震数据的反演效率,还能够提高叠前地震数据的反演分辨率。
图2示出了本发明实施例提供的叠前地震反演方法中步骤103的实现流程,为便于描述,仅示出了与本发明实施例相关的部分,详述如下:
在本发明的一实施例中,为了进一步提高叠前地震数据的反演效率和反演分辨率,如图2所示,步骤103,根据目标储层的地震角度道集及当前反演迭代模型确定能量误差,包括:
步骤201,根据目标储层的地震角度道集及当前反演迭代模型确定叠前地震数据的反射系数;
步骤202,根据叠前地震数据的反射系数、子波数据及当前反演迭代模型确定合成地震记录;
步骤203,根据叠前地震数据及合成地震记录确定能量误差。
鉴于目标储层包含很多个地层介质段,在确定叠前地震数据dr的反射系数时,每一次反演过程都会得到很多个反射系数,再由得到的很多个反射系数依序组成一个反射序列。即本发明实施例中得到的反射系数为反射系数序列,可以以矩阵的形式表现。
其中,假设目标储层包含N个地层介质段,则根据下述公式可以得知,在一次迭代过程中可以获得N-1个反射系数,此处以第i次迭代过程为例进行说明。此时当前反演迭代模型表示为mpi(x,t)=(vpi,vsipi)。假设n表示在第i次迭代过程中的反射系数的序号,k表示目标储层中地层介质段的序号,即目标储层中第k个地层介质段,第i次迭代过程中第n个反射系数可以表示为Rpin。具体在第i次迭代过程中确定第n个反射系数Rpin时,可以通如下公式确定:
Figure BDA0002079770440000071
且满足:
Figure BDA0002079770440000072
Figure BDA0002079770440000073
Figure BDA0002079770440000074
Δvpin=vpi(k-1)-vpik; 公式(5)
Δvsin=vsi(k-1)-vsik; 公式(6)
Δρpin=ρpi(k-1)pik; 公式(7)
其中,θ表示地震角度道集,n=1,2,3…N-3,N-2,N-1,k=2,3…N-2,N-1,N,且满足k=n+1。N为大于2的正整数,N表示目标储层中地层介质段的总个数。
vpi(k-1)表示第i次迭代过程中目标储层内第k-1个地层介质段的纵波速度,vpik表示第i次迭代过程中目标储层内第k个地层介质段的纵波速度;vsi(k-1)表示第i次迭代过程中目标储层内第k-1个地层介质段的横波速度,vsik表示第i次迭代过程中目标储层内第k个地层介质段的横波速度;ρpi(k-1)表示第i次迭代过程中目标储层内第k-1个地层介质段的密度,ρpik表示第i次迭代过程中目标储层内第k个地层介质段的密度。
vapin表示第i次迭代过程中目标储层内第k-1个地层介质段的纵波速度vpi(k-1),与第i次迭代过程中目标储层内第k个地层介质段的纵波速度vpik的平均纵波速度;vasin表示第i次迭代过程中目标储层内第k-1个地层介质段的横波速度vsi(k-1),与第i次迭代过程中目标储层内第k个地层介质段的横波速度vsik的平均横波速度;ρapin表示第i次迭代过程中目标储层内第k-1个地层介质段的密度ρpi(k-1),与第i次迭代过程中目标储层内第k个地层介质段的密度ρpik的平均密度。
Δvpin表示第i次迭代过程中目标储层内第k-1个地层介质段的纵波速度vpi(k-1),与第i次迭代过程中目标储层内第k个地层介质段的纵波速度vpik的纵波速度差;Δvsin表示第i次迭代过程中目标储层内第k-1个地层介质段的横波速度vsi(k-1),与第i次迭代过程中目标储层内第k个地层介质段的横波速度vsik的横波速度差;Δρpin表示第i次迭代过程中目标储层内第k-1个地层介质段的密度ρpi(k-1),与第i次迭代过程中目标储层内第k个地层介质段的密度ρpik的密度差。
在利用上述公式(1)至公式(7)确定第i次迭代过程中的第n个反射系数Rpin后,可以依次确定第i次迭代过程中的N-1个反射系数Rpi,N-1个反射系数Rpi可以以序列形式表现。该N-1个反射系数Rpi即为第i次迭代过程中叠前地震数据dr的反射系数Rpi。进而利用叠前地震数据的反射系数Rpi,子波数据及当前反演迭代模型确定合成地震记录,具体的可以通过下述公式确定合成地震记录:
Figure BDA0002079770440000081
其中,ds表示合成地震记录,Rpi表示第i次迭代过程中得到的叠前地震数据dr的反射系数,w(t)表示叠前地震数据dr的子波数据,
Figure BDA0002079770440000082
表示褶积,L表示第i次迭代过程中得到的叠前地震数据dr的反射系数Rpi与叠前地震数据dr的子波数据w(t)的褶积,mpi(x,t)表示当前反演迭代模型。
在确定合成地震记录ds后,具体可以通过如下公式确定能量误差:
Figure BDA0002079770440000083
其中,E表示能量误差,dr表示叠前地震数据,ds表示合成地震记录。
在本发明实施例中,根据目标储层的地震角度道集及当前反演迭代模型确定叠前地震数据的反射系数,根据叠前地震数据的反射系数、子波数据及当前反演迭代模型确定合成地震记录,进而根据叠前地震数据及合成地震记录确定能量误差,可以进一步提高叠前地震数据的反演效率和反演分辨率。
图3示出了本发明实施例提供的叠前地震反演方法中步骤104的实现流程,为便于描述,仅示出了与本发明实施例相关的部分,详述如下:
在本发明的一实施例中,为了进一步提高叠前地震数据的反演效率和反演分辨率,如图3所示,步骤104中,根据叠前地震数据的先验模型和当前反演迭代模型确定更新后的反演迭代模型,包括:
步骤301,根据叠前地震数据的反射系数及子波数据确定叠前地震数据的高斯概率分布;
步骤302,根据当前反演迭代模型确定当前反演迭代模型的柯西概率分布;
步骤303,根据叠前地震数据的先验模型及当前反演迭代模型,与叠前地震数据的高斯概率分布及当前反演迭代模型的柯西概率分布,确定更新后的反演迭代模型。
在确定更新的反演迭代模型时,首先可以确定叠前地震数据的高斯概率分布,具体可以通过如下公式确定叠前地震数据的高斯概率分布:
Figure BDA0002079770440000091
其中,Gpi表示第i次迭代过程中得到的叠前地震数据dr的高斯概率分布,Rpi表示第i次迭代过程中得到的叠前地震数据dr的反射系数,w(t)表示叠前地震数据dr的子波数据,
Figure BDA0002079770440000092
表示褶积,L表示第i次迭代过程中得到的叠前地震数据dr的反射系数Rpi与叠前地震数据dr的子波数据w(t)的褶积。
在基于叠前地震数据dr的反射系数Rpi及子波数据w(t)确定叠前地震数据dr的高斯概率分布Gpi后,可以根据当前反演迭代模型mpi(x,t)确定当前反演迭代模型mpi(x,t)的柯西概率分布,假设当前反演迭代模型mpi(x,t)的柯西概率分布采用Cpi表示。
具体的,在确定更新后的反演迭代模型时,当前反演迭代模型为mpi(x,t),假设更新后的反演迭代模型采用mp(i+1)(x,t)表示,具体可以通过如下公式确定更新后的反演迭代模型采用mp(i+1)(x,t):
mp(i+1)(x,t)=(Gpi+μCpi+αCov-1)-1(Ldr+αCov-1mpi(x,t)); 公式(11)
其中,mp(i+1)(x,t)表示更新后的反演迭代模型,mpi(x,t)表示当前反演迭代模型,Gpi表示第i次迭代过程中得到的叠前地震数据dr的高斯概率分布,Cpi表示第i次迭代过程中得到的当前反演迭代模型mpi(x,t)的柯西概率分布;Cov-1表示叠前地震数据的先验模型mp(x,t)的协方差矩阵,μ表示当前反演迭代模型mpi(x,t)的稀疏性参数,用于调节当前反演迭代模型mpi(x,t)的稀疏性,α表示当前反演迭代模型mpi(x,t)的权重系数,用于调整当前反演迭代模型mpi(x,t)的权重。在本发明的一实施例中,满足0<μ<1;另外,0<α<1。
在本发明实施例中,根据叠前地震数据的反射系数及子波数据确定叠前地震数据的高斯概率分布,根据当前反演迭代模型确定当前反演迭代模型的柯西概率分布,根据叠前地震数据的先验模型及当前反演迭代模型,与叠前地震数据的高斯概率分布及当前反演迭代模型的柯西概率分布,确定更新后的反演迭代模型,更新后的反演迭代模型满足柯西概率分布,因此,能够进一步提高叠前地震数据的反演效率和反演分辨率。
在本发明的一实施例中,为了进一步提高叠前地震数据的反演效率和反演分辨率,在本发明实施例中,更新后的反演迭代模型处于区域岩石物理趋势上边界和下边界范围内。具体的,可以通过如下公式表示:
mlb≤mp(i+1)(x,t)≤mμb; 公式(12)
其中,mlb表示区域岩石物理趋势下边界,mμb表示区域岩石物理趋势上边界,可以通过实测测井资料统计获得,将将更新后的反演迭代模型mp(i+1)(x,t)(即反演结果)约束在区域岩石物理趋势下边界mlb与下边界mμb范围内,反演结果才符合岩石物理规律,可以进一步提高前地震数据的反演效率和反演分辨率。
图4示出了本发明实施例提供的叠前地震反演方法中步骤303的实现流程,为便于描述,仅示出了与本发明实施例相关的部分,详述如下:
在本发明的一实施例中,为了进一步提高叠前地震数据的反演效率和反演分辨率,如图4所示,步骤303,根据当前反演迭代模型确定当前反演迭代模型的柯西概率分布,包括:
步骤401,根据辅助算子及弹性参数协方差矩阵确定模型权重矩阵;
步骤402,根据模型权重矩阵及当前反演迭代模型确定当前反演迭代模型的柯西概率分布。
在确定第i次迭代过程中当前反演迭代模型mpi(x,t)的柯西概率分布时,可以引入一大小为3×3M(其中,当前反演迭代模型mpi(x,t)的大小,与先验模型mp(x,t)=(vp,vsp)和反演初始模型mp0(x,t)=(vp0,vs0p0)的大小一致,也可以用M表示)的辅助算子Di,该辅助算子Di是一个与当前反演迭代模型mpi(x,t)的纵波速度vpi、横波速度vsi及密度ρpi相关的矩阵。具体的,可以通过下述公式确定辅助算子Di
Figure BDA0002079770440000111
其中,Di表示辅助算子,用于提取某一时间点的弹性参数,i表示迭代次数,M表示当前反演迭代模型mpi(x,t)的大小,n1表示反演参数的个数,例如n1=3表示三参数反演,l表示辅助算子Di的列数,且1≤l≤3M。
进一步的,设置大小为3×3的弹性参数协方差矩阵Ψ,弹性参数协方差矩阵Ψ为当前反演迭代模型mpi(x,t)的弹性参数的协方差矩阵,反映了当前反演迭代模型mpi(x,t)的纵波速度vpi、横波速度vsi及密度ρpi等弹性参数的相互关系,可以提高反演结果的稳定性。在本发明的一实施例中,弹性参数协方差矩阵Ψ是大小为3×3的矩阵。
在确定辅助算子Di及弹性参数协方差矩阵Ψ后,可以基于辅助算子D及协方差矩阵弹性参数协方差矩阵Ψ确定模型权重矩阵Φi,模型权重矩阵Φi是指将协方差矩阵弹性参数协方差矩阵Ψ中当前反演迭代模型mpi(x,t)的纵波速度vpi、横波速度vsi及密度ρpi等弹性参数的权重分配到当前反演迭代模型mpi(x,t)后得到的矩阵。具体的,可以通过如下公式确定模型权重矩阵Φi
Φi=(Di)TΨ-1(Di); 公式(13)
其中,Φi表示模型权重矩阵,Di表示辅助算子,Ψ表示协方差矩阵弹性参数协方差矩阵。
在确定模型权重矩阵Φi后,具体可以基于模型权重矩阵Φi和当前反演迭代模型mpi(x,t),通过如下公式确定在第i次迭代过程中当前反演迭代模型mpi(x,t)的柯西概率分布Cpi
Figure BDA0002079770440000112
其中,Cpi表示在第i次迭代过程中当前反演迭代模型mpi(x,t)的柯西概率分布,其可以采用大小为3M×3M的矩阵形式表示。M表示当前反演迭代模型mpi(x,t)的大小,在本发明实施例中,根据辅助算子及弹性参数协方差矩阵确定模型权重矩阵,进而根据模型权重矩阵及当前反演迭代模型确定当前反演迭代模型的柯西概率分布,可以进一步提高叠前地震数据的反演效率和反演分辨率。
本发明实施例中还提供了一种叠前地震反演装置,如下面的实施例所述。由于这些装置解决问题的原理与叠前地震反演方法相似,因此这些装置的实施可以参见方法的实施,重复之处不再赘述。
图5示出了本发明实施例提供的叠前地震反演装置的功能模块,为便于说明,仅示出了与本发明实施例相关的部分,详述如下:
参考图5,所述叠前地震反演装置所包含的各个模块用于执行图1对应实施例中的各个步骤,具体请参阅图1以及图1对应实施例中的相关描述,此处不再赘述。本发明实施例中,所述叠前地震反演装置包括数据采集模块501、模型建立模块502、能量误差确定模块503、更新模块504及反演结果获取模块505。
数据采集模块501,用于采集目标储层的叠前地震数据。
模型建立模块502,用于通过岩石物理建模确定叠前地震数据的先验模型和反演初始模型,将反演初始模型作为当前反演迭代模型;先验模型是反映目标储层的弹性参数的模型;反演初始模型是反映目标储层的初始弹性参数的模型。
能量误差确定模块503,用于根据目标储层的地震角度道集及当前反演迭代模型确定能量误差,能量误差满足高斯概率分布。
更新模块504,用于若能量误差不小于预设能量误差,根据叠前地震数据的先验模型和当前反演迭代模型确定更新后的反演迭代模型,更新后的反演迭代模型满足柯西概率分布;将更新后的反演迭代模型作为当前反演迭代模型。
反演结果获取模块505,用于将能量误差小于预设能量误差时的当前反演迭代模型作为叠前地震数据的反演结果。
在本发明实施例中,能量误差确定模块503确定的能量误差满足高斯分布,更新模块504更新后的反演迭代模型满足柯西概率分布,同时更新模块504在能量误差不小于预设能量误差时,将更新后的反演迭代模型作为当前反演迭代模型,能量误差确定模块503迭代执行根据目标储层的地震角度道集及当前反演迭代模型确定能量误差,直至能量误差小于预设能量误差,因此,不仅能够提高叠前地震数据的反演效率,还能够提高叠前地震数据的反演分辨率。
图6示出了本发明实施例提供的叠前地震反演装置中能量误差确定模块503的结构示意,为便于说明,仅示出了与本发明实施例相关的部分,详述如下:
在本发明的一实施例中,为了进一步提高叠前地震数据的反演效率和反演分辨率,参考图6,所述能量误差确定模块503所包含的各个单元用于执行图2对应实施例中的各个步骤,具体请参阅图2以及图2对应实施例中的相关描述,此处不再赘述。本发明实施例中,所述能量误差确定模块503包括反射系数确定单元601、合成地震记录确定单元602及能量误差确定单元603。
反射系数确定单元601,用于根据目标储层的地震角度道集及当前反演迭代模型确定叠前地震数据的反射系数。
合成地震记录确定单元602,用于根据叠前地震数据的反射系数、子波数据及当前反演迭代模型确定合成地震记录。
能量误差确定单元603,用于根据叠前地震数据及合成地震记录确定能量误差。
在本发明实施例中,反射系数确定单元601根据目标储层的地震角度道集及当前反演迭代模型确定叠前地震数据的反射系数,合成地震记录确定单元602根据叠前地震数据的反射系数、子波数据及当前反演迭代模型确定合成地震记录,能量误差确定单元603根据叠前地震数据及合成地震记录确定能量误差,可以进一步提高叠前地震数据的反演效率和反演分辨率。
图7示出了本发明实施例提供的叠前地震反演装置中更新模块504的结构示意,为便于说明,仅示出了与本发明实施例相关的部分,详述如下:
在本发明的一实施例中,为了进一步提高叠前地震数据的反演效率和反演分辨率,参考图7,所述更新模块504所包含的各个单元用于执行图3对应实施例中的各个步骤,具体请参阅图3以及图3对应实施例中的相关描述,此处不再赘述。本发明实施例中,所述更新模块504包括高斯概率分布确定单元701、柯西概率分布确定单元702及更新单元703。
高斯概率分布确定单元701,用于根据叠前地震数据的反射系数及子波数据确定叠前地震数据的高斯概率分布。
柯西概率分布确定单元702,用于根据当前反演迭代模型确定当前反演迭代模型的柯西概率分布。
更新单元703,用于根据叠前地震数据的先验模型及当前反演迭代模型,与叠前地震数据的高斯概率分布及当前反演迭代模型的柯西概率分布,确定更新后的反演迭代模型。
在本发明实施例中,高斯概率分布确定单元701根据叠前地震数据的反射系数及子波数据确定叠前地震数据的高斯概率分布,柯西概率分布确定单元702根据当前反演迭代模型确定当前反演迭代模型的柯西概率分布,更新单元703根据叠前地震数据的先验模型及当前反演迭代模型,与叠前地震数据的高斯概率分布及当前反演迭代模型的柯西概率分布,确定更新后的反演迭代模型,更新后的反演迭代模型基于高斯概率分布及当前反演迭代模型的柯西概率分布,可以进一步提高叠前地震数据的反演效率和反演分辨率。
图8示出了本发明实施例提供的叠前地震反演装置中柯西概率分布确定单元702的结构示意,为便于说明,仅示出了与本发明实施例相关的部分,详述如下:
在本发明的一实施例中,为了进一步提高叠前地震数据的反演效率和反演分辨率,参考图8,所述柯西概率分布确定单元702所包含的各个子单元用于执行图4对应实施例中的各个步骤,具体请参阅图4以及图4对应实施例中的相关描述,此处不再赘述。本发明实施例中,所述柯西概率分布确定单元702的包括弹性参数协方差矩阵确定子单元801和柯西概率分布确定子单元802。
弹性参数协方差矩阵确定子单元801,用于根据辅助算子及弹性参数协方差矩阵确定模型权重矩阵。
柯西概率分布确定子单元802,用于根据模型权重矩阵及当前反演迭代模型确定当前反演迭代模型的柯西概率分布。
在本发明实施例中,弹性参数协方差矩阵确定子单元801根据辅助算子及弹性参数协方差矩阵确定模型权重矩阵,柯西概率分布确定子单元802,用于根据模型权重矩阵及当前反演迭代模型确定当前反演迭代模型的柯西概率分布,可以进一步提高叠前地震数据的反演效率和反演分辨率。
本发明实施例还提供一种计算机设备,包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序,所述处理器执行所述计算机程序时实现上述叠前地震反演方法。
本发明实施例还提供一种计算机可读存储介质,所述计算机可读存储介质存储有执行上述叠前地震反演方法的计算机程序。
综上所述,在本发明实施例中,能量误差满足高斯分布,更新后的反演迭代模型满足柯西概率分布,同时在能量误差不小于预设能量误差时,将更新后的反演迭代模型作为当前反演迭代模型,迭代执行根据目标储层的地震角度道集及当前反演迭代模型确定能量误差,直至能量误差小于预设能量误差,因此,不仅能够提高叠前地震数据的反演效率,还能够提高叠前地震数据的反演分辨率。
本领域内的技术人员应明白,本发明的实施例可提供为方法、系统、或计算机程序产品。因此,本发明可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且,本发明可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等)上实施的计算机程序产品的形式。
本发明是参照根据本发明实施例的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器,使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。
这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中,使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品,该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。
这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上,使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理,从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。
以上所述的具体实施例,对本发明的目的、技术方案和有益效果进行了进一步详细说明,所应理解的是,以上所述仅为本发明的具体实施例而已,并不用于限定本发明的保护范围,凡在本发明的精神和原则之内,所做的任何修改、等同替换、改进等,均应包含在本发明的保护范围之内。

Claims (7)

1.一种叠前地震反演方法,其特征在于,包括:
步骤101:采集目标储层的叠前地震数据;
步骤102:通过岩石物理建模确定叠前地震数据的先验模型和反演初始模型,将反演初始模型作为当前反演迭代模型;先验模型是反映目标储层的弹性参数的模型;反演初始模型是反映目标储层的初始弹性参数的模型;
步骤103:根据目标储层的地震角度道集及当前反演迭代模型确定叠前地震数据的反射系数;根据叠前地震数据的反射系数、子波数据及当前反演迭代模型确定合成地震记录;根据叠前地震数据及合成地震记录确定能量误差;能量误差满足高斯概率分布;
若能量误差不小于预设能量误差,则执行步骤104;若能量误差小于预设能量误差,则执行步骤105;
步骤104:根据叠前地震数据的反射系数及子波数据确定叠前地震数据的高斯概率分布;根据当前反演迭代模型确定当前反演迭代模型的的柯西概率分布;根据叠前地震数据的先验模型及当前反演迭代模型,与叠前地震数据的高斯概率分布及当前反演迭代模型的柯西概率分布,确定更新后的反演迭代模型,更新后的反演迭代模型满足柯西概率分布;将更新后的反演迭代模型作为当前反演迭代模型,返回执行步骤103进行迭代;
步骤105:将能量误差小于预设能量误差时的当前反演迭代模型作为叠前地震数据的反演结果。
2.如权利要求1所述的叠前地震反演方法,其特征在于,根据当前反演迭代模型确定当前反演迭代模型的柯西概率分布,包括:
根据辅助算子及弹性参数协方差矩阵确定模型权重矩阵;
根据模型权重矩阵及当前反演迭代模型确定当前反演迭代模型的柯西概率分布。
3.如权利要求1所述的叠前地震反演方法,其特征在于,更新后的反演迭代模型处于区域岩石物理趋势上边界和下边界范围内。
4.一种叠前地震反演装置,其特征在于,包括:
数据采集模块,用于采集目标储层的叠前地震数据;
模型建立模块,用于通过岩石物理建模确定叠前地震数据的先验模型和反演初始模型,将反演初始模型作为当前反演迭代模型;先验模型是反映目标储层的弹性参数的模型;反演初始模型是反映目标储层的初始弹性参数的模型;
能量误差确定模块,用于根据目标储层的地震角度道集及当前反演迭代模型确定叠前地震数据的反射系数;根据叠前地震数据的反射系数、子波数据及当前反演迭代模型确定合成地震记录;根据叠前地震数据及合成地震记录确定能量误差;能量误差满足高斯概率分布;
更新模块,用于若能量误差不小于预设能量误差,根据叠前地震数据的反射系数及子波数据确定叠前地震数据的高斯概率分布;根据当前反演迭代模型确定当前反演迭代模型的的柯西概率分布;根据叠前地震数据的先验模型及当前反演迭代模型,与叠前地震数据的高斯概率分布及当前反演迭代模型的柯西概率分布,确定更新后的反演迭代模型,更新后的反演迭代模型满足柯西概率分布;将更新后的反演迭代模型作为当前反演迭代模型;
反演结果获取模块,用于将能量误差小于预设能量误差时的当前反演迭代模型作为叠前地震数据的反演结果。
5.如权利要求4所述的叠前地震反演装置,其特征在于,所述更新模块包括弹性参数协方差矩阵确定子单元和柯西概率分布确定子单元;
所述弹性参数协方差矩阵确定子单元,用于根据辅助算子及弹性参数协方差矩阵确定模型权重矩阵;
所述柯西概率分布确定子单元,用于根据模型权重矩阵及当前反演迭代模型确定当前反演迭代模型的柯西概率分布。
6.一种计算机设备,包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序,其特征在于,所述处理器执行所述计算机程序时实现权利要求1至3任一所述方法。
7.一种计算机可读存储介质,其特征在于,所述计算机可读存储介质存储有执行权利要求1至3任一所述方法的计算机程序。
CN201910467123.9A 2019-05-31 2019-05-31 叠前地震反演方法及装置 Active CN112014875B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910467123.9A CN112014875B (zh) 2019-05-31 2019-05-31 叠前地震反演方法及装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910467123.9A CN112014875B (zh) 2019-05-31 2019-05-31 叠前地震反演方法及装置

Publications (2)

Publication Number Publication Date
CN112014875A CN112014875A (zh) 2020-12-01
CN112014875B true CN112014875B (zh) 2022-08-05

Family

ID=73501978

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910467123.9A Active CN112014875B (zh) 2019-05-31 2019-05-31 叠前地震反演方法及装置

Country Status (1)

Country Link
CN (1) CN112014875B (zh)

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102636812A (zh) * 2012-04-18 2012-08-15 中国石油天然气股份有限公司 一种获得碳酸盐岩储层储集空间体积的方法

Family Cites Families (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101556337B (zh) * 2008-04-10 2011-12-07 中国石油天然气集团公司 一种确定地下深层特殊岩性体的方法
CN102109616B (zh) * 2009-12-23 2012-08-15 中国石油天然气集团公司 一种沉积旋回约束的波阻抗反演方法
US9772415B2 (en) * 2011-08-05 2017-09-26 Saudi Arabian Oil Company Correcting time lapse seismic data for overburden and recording effects
CN107544091A (zh) * 2017-07-21 2018-01-05 北京阳光杰科科技股份有限公司 一种高精度储层孔隙度定量预测方法及其应用

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102636812A (zh) * 2012-04-18 2012-08-15 中国石油天然气股份有限公司 一种获得碳酸盐岩储层储集空间体积的方法

Also Published As

Publication number Publication date
CN112014875A (zh) 2020-12-01

Similar Documents

Publication Publication Date Title
CN101681394B (zh) 来自并发地球物理源的数据的迭代反演
CA2690373C (en) Method for velocity analysis using waveform inversion in laplace domain for geophysical imaging
CN103713315B (zh) 一种地震各向异性参数全波形反演方法及装置
US8406081B2 (en) Seismic imaging systems and methods employing tomographic migration-velocity analysis using common angle image gathers
EP2491516B1 (en) Method for optimization with gradient information
US9229123B2 (en) Method for handling rough sea and irregular recording conditions in multi-sensor towed streamer data
CN110031896B (zh) 基于多点地质统计学先验信息的地震随机反演方法及装置
Cherpeau et al. Stochastic structural modelling in sparse data situations
CA2964893C (en) Structure tensor constrained tomographic velocity analysis
WO2008157032A1 (en) Creating an absorption parameter model
US20130297274A1 (en) Methods and systems regarding models of underground formations
CN109581501B (zh) 用于沙漠区深度域速度建模的方法
CN111123359B (zh) 随钻测井与地层格架约束的井周地震成像探测方法及装置
CN106574980A (zh) 用于地下地质体的岩石性质估计的系统和方法
CN112014875B (zh) 叠前地震反演方法及装置
US10422898B2 (en) Seismic data processing
CN110646846B (zh) Vti介质各向异性参数确定方法、装置和设备
GB2502918A (en) Computer estimation method, and method for oil exploration and development using such a method
WO2021191722A1 (en) System and method for stochastic full waveform inversion
CN111751886A (zh) 一种基于微地震监测数据的页岩气藏裂缝建模方法
CN112649848A (zh) 利用波动方程求解地震波阻抗的方法和装置
CN113267810B (zh) 地震勘探全深度速度建模方法及装置
US20240168189A1 (en) Systems, devices, and methods for generating average velocity maps of subsurface formations
CN117950045A (zh) 提高密度和速度反演精度的方法、装置、设备及介质
CN115267901A (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