CN112327359B - 一种基于成像能流矢量的弹性逆时偏移方法 - Google Patents

一种基于成像能流矢量的弹性逆时偏移方法 Download PDF

Info

Publication number
CN112327359B
CN112327359B CN202011096231.9A CN202011096231A CN112327359B CN 112327359 B CN112327359 B CN 112327359B CN 202011096231 A CN202011096231 A CN 202011096231A CN 112327359 B CN112327359 B CN 112327359B
Authority
CN
China
Prior art keywords
imaging
elastic wave
time
field
elastic
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
CN202011096231.9A
Other languages
English (en)
Other versions
CN112327359A (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.)
Institute of Oceanographic Instrumentation Shandong Academy of Sciences
Original Assignee
Institute of Oceanographic Instrumentation Shandong Academy of Sciences
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 Institute of Oceanographic Instrumentation Shandong Academy of Sciences filed Critical Institute of Oceanographic Instrumentation Shandong Academy of Sciences
Priority to CN202011096231.9A priority Critical patent/CN112327359B/zh
Publication of CN112327359A publication Critical patent/CN112327359A/zh
Application granted granted Critical
Publication of CN112327359B publication Critical patent/CN112327359B/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. analysis, for interpretation, for correction
    • G01V1/36Effecting static or dynamic corrections on records, e.g. correcting spread; Correlating seismic signals; Eliminating effects of unwanted energy
    • 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. analysis, for interpretation, for correction
    • 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. analysis, for interpretation, for correction
    • G01V1/34Displaying seismic recordings or visualisation of seismic data or attributes
    • G01V1/345Visualisation of seismic data or attributes, e.g. in 3D cubes
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00Details of seismic processing or analysis
    • G01V2210/50Corrections or adjustments related to wave propagation
    • G01V2210/51Migration
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01VGEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00Details of seismic processing or analysis
    • G01V2210/70Other details related to processing
    • G01V2210/74Visualisation of seismic data

Abstract

本发明属于勘探地球物理学领域,涉及一种基于成像能流矢量的弹性逆时偏移方法。该方法包括:基于弹性波正向延拓算子,获得震源端弹性波质点振动速度矢量场;基于弹性波反向延拓算子,获得检波端弹性波应力张量场;对所述的震源端质点振动速度矢量场和检波端应力张量场进行点积运算,构建得到各网格节点在检查点时刻的成像能流矢量;根据能流矢量散度成像条件,对所构建的成像能流矢量进行成像,获得弹性波反射能量成像结果。本发明从矢量化弹性波场中实现了标量化弹性波反射能量成像结果提取,得到的成像结果物理意义明确;而且仅利用震源端弹性波质点振动速度矢量场和检波端弹性波应力张量场进行成像,存储消耗较低,存储效率高。

Description

一种基于成像能流矢量的弹性逆时偏移方法
技术领域
本发明属于勘探地球物理学领域,涉及一种基于成像能流矢量的弹性逆时偏移方法。
背景技术
弹性逆时偏移方法,可以利用多分量地震数据对复杂地球介质进行构造成像。合理的成像条件,作为弹性逆时偏移方法中重要组成部分,可以从震源端和检波端弹性波场中实现标量成像结果提取,从而实现复杂地球介质的构造成像。其中,基于能量密度的弹性能量互相关成像条件,可以从弹性质点振动速度-应力矢量场中提取出弹性波反射能量成像结果,物理意义明确。但是,以质点振动速度- 应力为自变量表征的弹性能量互相关成像条件存储效率较低且海底数据处理时计算稳定性差,这就限制了基于弹性能量互相关成像条件的弹性逆时偏移方法的实际应用。
中国发明专利申请CN201810967971.1,公开一种基于lowrank有限差分的弹性逆时偏移方法及系统,该方法通过利用平面波理论及lowrank分解方法构建了时空兼顾的有限差分算子,是基于波数域实现的,计算量大,运算复杂。
发明内容
本发明的目的是提供一种基于成像能流矢量的弹性逆时偏移方法,可以解决基于弹性能量互相关成像条件的弹性逆时偏移方法存储效率低且海底数据处理稳定性差的问题,提高弹性能量逆时偏移成像方法在实际应用中效率。
为实现上述目的,本发明采用的技术方案是:一种基于成像能流矢量的弹性逆时偏移方法,包括:
基于弹性波正向延拓算子,获得震源端弹性波质点振动速度矢量场;
基于弹性波反向延拓算子,获得检波端弹性波应力张量场;
对所述的震源端质点振动速度矢量场和检波端应力张量场进行点积运算,构建得到各网格节点在检查点时刻的成像能流矢量;
根据能流矢量散度成像条件,对所构建成像能流矢量进行成像,获得弹性波反射能量成像结果。
作为本发明的一种优选方式,所述成像能流矢量的构建方法为:
根据设定的检查点时刻,提取震源端弹性波质点振动速度矢量场
Figure RE-GDA0002844072110000021
根据设定的检查点时刻,提取检波端弹性应力张量场
Figure RE-GDA0002844072110000022
以所述的相同检查点时刻的震源端质点振动速度矢量场
Figure RE-GDA0002844072110000023
和检波端应力张量场
Figure RE-GDA0002844072110000024
为自变量,根据下式进行点积运算构建成像能流矢量:
Figure RE-GDA0002844072110000025
式中,
Figure RE-GDA0002844072110000026
表示各个网格点处第i个检查点时刻Ti的成像能流失量,
Figure RE-GDA0002844072110000027
为成像能流矢量沿x方向的分量,
Figure RE-GDA0002844072110000028
为成像能流矢量沿y方向的分量,
Figure RE-GDA0002844072110000029
为成像能流矢量沿z方向的分量,
Figure RE-GDA00028440721100000210
表示各个检查点处第i个检查点时刻Ti的检波端弹性应力张量场,
Figure RE-GDA00028440721100000211
表示各个检查点处第i个检查点时刻Ti的震源端弹性波质点振动速度矢量场;
成像能流矢量采用分量形式,表示为:
Figure RE-GDA00028440721100000212
保存每一检查点时刻的成像能流密度矢量。
进一步优选地,对所构建的成像能流矢量进行成像的方法包括:
根据设定检查点时刻,确定检查点时间间隔:
ΔTi=Ti-Ti-1,i=1,2,...,M (2)
其中,ΔTi为第i个检查点时刻处检查点时间间隔,Ti为任意一个检查点时刻时间,Ti+1表示相邻上一个检查点时刻时间;
针对每一检查点时刻的成像能流密度矢量及检查点时间间隔,利用如下所示的基于能流矢量的弹性能量成像条件进行成像,获得最终弹性能量成像结果:
Figure RE-GDA0002844072110000031
其中,Imenergy表示基于能流矢量散度成像条件求取的弹性波反射能量成像结果,i表示检查点时刻序号,共计N个检查点时刻,ishot表示炮点序号,共计Ishot 炮,符号
Figure RE-GDA0002844072110000032
表示散度运算,
Figure RE-GDA0002844072110000033
表示第i个检查点时刻成像能流矢量散度,采用分量形式表示为
Figure RE-GDA0002844072110000034
ΔTi表示检查点时间间隔。
进一步优选地,利用交错网格有限差分格式对弹性波动方程进行离散,获得弹性波延拓算子。
进一步优选地,对所述的弹性波延拓算子加载给定震源子波及介质模型,沿时间正向传播方向进行弹性波场延拓,获得弹性波正向延拓算子。
进一步优选地,所述的给定震源子波加载到弹性波场延拓算子中弹性波正应力或剪切应力上。
进一步优选地,对所述的弹性波延拓算子利用预先给定的介质模型,将多分量地震记录作为边界条件,沿时间反向传播方向进行弹性波场延拓获得弹性波反向延拓算子。
相对于现有技术,本发明的有益效果如下:提供了一种基于成像能流矢量的弹性逆时偏移方法,该方法构建了成像能流矢量,实现震源端弹性波场和检波端弹性波场相关矢量化关系表征,并利用所构建的能流矢量散度成像条件实现弹性波反射能量成像结果提取,物理意义明确,避免了所有检查点时刻震源端应力张量场存储消耗,提高了存储效率,而且解决了基于弹性能量互相关成像条件的弹性逆时偏移方法在海底多分量地震数据处理中计算不稳定问题。
附图说明
图1是基于成像能流矢量的弹性逆时偏移方法实施的流程示意图;
图2是预先给定的震源子波波形;
图3a是构建的凹陷模型的纵波速度场;
图3b是构建的凹陷模型的横波速度场;
图3c是构建的凹陷模型的密度场;
图4a是正向外推至1s时刻构建的震源端弹性波质点振动速度波场沿x方向的分量;
图4b是正向外推至1s时刻构建的震源端弹性波质点振动速度波场沿z方向的分量;
图5a是施加震源后地表采集到的时间采样时长2.7s的多分量地震记录沿x 方向的分量;
图5b是施加震源后地表采集到的时间采样时长2.7s的多分量地震记录沿z 方向的分量;
图6a是反向外推至1s时刻构建的震源端弹性波应力场沿xx方向的分量;
图6b是反向外推至1s时刻构建的震源端弹性波应力场沿xz方向的分量;
图6c是反向外推至1s时刻构建的震源端弹性波应力场沿zz方向的分量;
图7a是检查点1s时刻构建的成像能流矢量沿x方向的分量;
图7b是检查点1s时刻构建的成像能流矢量沿z方向的分量;
图8a是采用能流矢量散度成像条件获得的单炮弹性波能量成像剖面;
图8b是最终采用能流矢量散度成像条件获得的多炮叠加的弹性波能量成像剖面。
具体实施方式
为了使本技术领域的人员更好地理解本说明书中的技术方案,下面将结合本说明书实施例中的附图,对本说明书实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅仅是本说明书中的一部分实施例,而不是全部的实施例。基于本说明书中的一个或多个实施例,本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例,都应当属于本说明书实施例保护的范围。
本发明提供的其中一个实施例是:一种基于成像能流矢量的弹性逆时偏移方法,通过在一阶弹性波动方程离散基础上构建正向延拓算子及反向延拓算子,可以在实现弹性波场正向延拓及反向延拓的同时,实现检查点时刻震源端质点振动速度矢量场和检波端应力张量场构建并保存,从而可以有效减少内存消耗。通过利用所构建的成像能流矢量,实现震源端质点振动速度矢量场和检波端弹性应力场相关矢量化特征表征,并结合所构建的能流矢量散度成像条件,可以从弹性矢量波场中提取得到标量化的弹性波反射能量成像结果,并解决海底多分量地震数据处理时计算不稳定问题。
下面以一个具体的应用场景为例对本实施例进行说明。一种基于成像能流矢量的弹性逆时偏移方法,流程如图1所示,具体包括以下步骤:
S1:基于弹性波正向延拓算子,获得震源端弹性波质点振动速度矢量场。
在本发明的一个实施例中,利用交错网格有限差分格式对弹性波动方程进行离散,获得弹性波延拓算子;对所述的弹性波延拓算子加载给定震源子波及介质模型,沿时间正向传播方向延拓得到弹性波正向延拓算子,实现震源端弹性波场延拓;根据设定的检查点时刻,保存所有检查点时刻震源端弹性波质点振动速度矢量场。
例如一些实施场景中,利用交错网格有限差分格式对弹性波动方程进行离散,获得弹性波延拓算子。
其中,一阶弹性波动方程表示为:
Figure RE-GDA0002844072110000051
其中
Figure RE-GDA0002844072110000052
是弹性波应力张量场,考虑到σij=σji,i,j=x,y,z,可以表示为σ=(σxxyyzzxyxzyz)T,上标“T”表示转置,
Figure RE-GDA0002844072110000053
表示应力张量的时间导数,C是弹性介质刚度矩阵,
Figure RE-GDA0002844072110000054
是微分矩阵, lx,ly和lz分别表示沿x,y和z方向上的空间导数,v=(vx,vy,vz)T表示弹性波质点振动速度矢量场,vx,vy和vz分别表示质点振动速度矢量场沿x,y和z方向的分量,
Figure RE-GDA0002844072110000055
表示质点振动速度矢量场的时间导数,ρ是密度。
其中,对于各向同性介质而言,刚度矩阵可以表示为:
Figure RE-GDA0002844072110000061
其中,λ和μ是拉梅系数。
利用交错网格有限差分格式及PML边界条件对一阶弹性波动方程(1)进行离散,得到弹性波延拓算子:
Figure RE-GDA0002844072110000062
其中,上标“S”表示震源端弹性波场,(n+1)Δt和nΔt表示整时间点, (n+1/2)Δt和(n-1/2)Δt为半时间节点,n=1,2,...,N,T0=NΔt表示总的地震记录接收时长,(σS)(n+1)Δt和(σS)nΔt分别表示(n+1)Δt和nΔt时刻震源端弹性波应力张量场,(vS)(n+1/2)Δt和(vS)(n-1/2)Δt分别表示(n+1/2)Δt和(n-1/2)Δt时刻震源端弹性波质点振动速度矢量场,η是边界吸收系数,在目标区域内吸收系数η=0,在边界吸收区域内吸收系数η=200(0.5-0.5cos(πr/R)),r=1,2,...,R,R是吸收层的厚度,π表示圆周率,Δt为时间采样间隔,
Figure RE-GDA0002844072110000063
Figure RE-GDA0002844072110000064
分别表示交错网格有限差分算子。
其中l表示交错网格差分格式,下标“x”,“y”和“z”分别表示沿着x,y 和z方向的空间差分格式,上标“+”表示向前差分格式,上标“-”表示向后差分格式,具体可以表示为:
Figure RE-GDA0002844072110000071
其中,Δx、Δy和Δz分别表示沿x、y和z方向上的采样间隔,
Figure RE-GDA0002844072110000072
表示
Figure RE-GDA0002844072110000073
阶交错网格有限差分系数,u(xi,yj,zk)表示任意空间网格点(xi,yj,zk)的弹性波场函数,既可以是弹性质点振动速度矢量场,也可以是弹性应力张量场。
有限差分系数
Figure RE-GDA0002844072110000074
通过以下方程确定:
Figure RE-GDA0002844072110000075
一些实施场景中,对所述的弹性波延拓算子加载给定震源子波及介质模型,沿时间正向传播方向延拓构建弹性波正向延拓算子,实现震源端弹性波场延拓。其中,震源子波可以表示为:
Fp=h(t)Φ(x,y,z) (5)
其中,Fp表示震源子波,h(t)表示地震子波函数,Φ(x,y,z)表示震源空间衰减函数。一些实施场景中,h(t)可以为雷克子波 h(t)=[1-2(πfpt)2]exp[-(πfpt)2],也可以为其他子波,本说明书对此不作限定,其中,此处雷克子波表达式中π表示圆周率,fp表示主频,t表示时间,exp 表示指数函数。
一些实施场景中,Φ(x,y,z)=exp{-α2[(x-xs)2+(y-ys)2+(z-zs)2]},α表示衰减参数,(xs,ys,zs)表示震源位置。附图3为预先给定的总时长为2700ms,主频为30Hz的雷克子波波形的震源子波。
一些实施场景中,将震源子波加载到弹性波场延拓算子中弹性波正应力上,也可以将震源子波加载到弹性波场延拓算子中弹性波剪切应力上,本发明对此不作限定。
一些实施场景中,给定介质纵波速度模型、横波速度模型及密度模型,在弹性波延拓算子(式(2))中施加震源子波,就可以得到弹性波正向延拓算子:
Figure RE-GDA0002844072110000081
式中,上标“S”表示震源端弹性波场函数。
一些实施场景中,利用上述的弹性波正向延拓算子,从0时刻NΔt到时刻沿着时间正向方向进行波场延拓,就可以得到所有时刻震源端弹性波质点振动速度矢量-应力张量场。
一些实施场景中,预先设定M个检查点,第i,i=1,2,...,M个检查点对应时间节点为Ti,根据设定的检查点时刻Ti,保存所有检查点时刻震源端弹性质点振动速度矢量场
Figure RE-GDA0002844072110000082
其中,交错网格有限差分格式中,质点振动速度张量定义在半时间节点上,所以检查点时刻对应为半时间节点。
附图2为给定的震源子波波形,附图3a为给定的纵波速度模型、图3b为给定的横波速度模型、图3c为给定的密度模型。附图4a是根据附图2给定的震源子波、附图3a给定的纵波速度模型、图3b给定的横波速度模型和图3c给定的密度模型采用弹性波正向延拓算子沿正向时间外推至1s时刻构建的震源端弹性波质点振动速度矢量场沿x方向的分量;图4b是根据给定的根据附图2给定的震源子波、附图3a给定的纵波速度模型、图3b给定的横波速度模型和图3c给定的密度模型采用弹性波正向延拓算子沿正向时间外推至1s时刻构建的震源端弹性波质点振动速度矢量场沿z方向的分量。从附图4a至附图4b可以明显看出弹性波质点振动速度矢量场各分量之间具有不同传播特征,说明弹性质点振动速度矢量场具有明显矢量特征。
S2:基于弹性波反向延拓算子,获得检波端弹性波应力张量场。
本发明的一个实施例中,基于弹性波场延拓算子,可以利用给定的介质模型、多分量地震记录,沿时间反向传播方向延拓构建弹性波反向延拓算子,实现检波端弹性波场延拓,获得检波端弹性波应力张量场。
在一些实施场景中,所述弹性波反向延拓算子,可以包括:利用预先给定的介质模型,将多分量地震记录作为边界条件,沿时间反向传播方向进行弹性波场延拓,即弹性波反向延拓算子。其中,弹性波反向延拓算子为:
Figure RE-GDA0002844072110000091
其中,上标“R”表示检波端波场函数。利用上述的弹性波反向延拓算子,从NΔt时刻到0时刻沿着时间反向方向进行波场延拓,就可以得到所有整时间点时刻检波端弹性波质点振动速度矢量-应力张量场。
图5为给定的多分量地震记录,其中,图5a是地表采集到的2.7s的多分量地震记录沿x方向的分量示意图,图5b是地表采集到的2.7s的多分量地震记录沿z 方向的分量记录图。附图6a是根据附图5给定的多分量地震记录、附图3a给定的纵波速度模型、图3b给定的横波速度模型和图3c给定的密度模型采用弹性波反向延拓算子沿反向时间外推并进行时间插值得到的1s时刻构建的检波端弹性波应力张量沿xx方向的分量;图6b是根据附图5给定的多分量地震记录、附图3a给定的纵波速度模型、图3b给定的横波速度模型和图3c给定的密度模型采用弹性波反向延拓算子沿反向时间外推并进行时间插值得到的1s时刻构建的检波端弹性波应力张量沿xz方向的分量;图6c是根据附图5给定的多分量地震记录、附图3a给定的纵波速度模型、图3b给定的横波速度模型和图3c给定的密度模型采用弹性波反向延拓算子沿反向时间外推并进行时间插值得到的1s时刻构建的检波端弹性波应力张量沿zz方向的分量。从附图6a至附图6c可以明显看出弹性波应力张量场各分量之间具有不同传播特征,说明弹性波应力张量场具有明显矢量特征。
S3:根据设定的检查点时刻,对所述的震源端质点振动速度矢量场和检波端应力张量场进行点积运算,构建成像能流矢量,并保存检查点时刻成像能流矢量。
本发明的一个实施例中,构建成像能流矢量,包括:根据设定的检查点时刻,提取震源端弹性波质点振动速度矢量场
Figure RE-GDA0002844072110000101
根据设定的检查点时刻Ti,对检波端弹性应力张量场(σR)nΔt进行时间插值,提取检查点时刻检波端应力张量场
Figure RE-GDA0002844072110000102
基于检查点时刻时间一致性原则,对震源端质点振动速度矢量场
Figure RE-GDA0002844072110000103
和检波端应力张量场
Figure RE-GDA0002844072110000104
进行点积运算,可以构建得到该检查点时刻处根据下式构建成像能流矢量
Figure RE-GDA0002844072110000105
保存检查点时刻成像能流矢量。
在一些实施场景中,根据预先给定的检查点时刻Ti,i=1,2,...,M,当弹性波场反向延拓到达第i个检查点对应的检查点时刻(n+1/2)Δt时,从所述的保存的所有检查点时刻的震源端弹性质点振动速度矢量场中,提取出第i个检查点时刻 Ti的震源端弹性质点振动速度矢量场
Figure RE-GDA0002844072110000106
在一些实施场景中,根据预先给定的检查点时刻Ti,当弹性波场反向延拓到达第i个检查点对应的检查点时刻(n+1/2)Δt时,对所述的当前整时间节点nΔt 时刻对应的弹性波应力张量场(σR)nΔt和下一个整时间节点(n+1)Δt时刻对应的弹性波应力场(σR)(n+1)Δt进行时间插值1/2[(σR(nΔt+(σR)(n+1)Δt],可以提取当前检查点时刻Ti对应的检波端应力张量场
Figure RE-GDA0002844072110000107
在一些实施场景中,基于检查点时刻时间一致性原则,对相同的给定的检查点时刻Ti下的震源端弹性波质点振动速度矢量场和检波端弹性波应力张量场进行点积运算,可以构建得到该检查点时刻下成像能流矢量
Figure RE-GDA0002844072110000108
Figure RE-GDA0002844072110000109
其中,
Figure RE-GDA00028440721100001010
为成像能流矢量沿x方向的分量,
Figure RE-GDA0002844072110000111
为成像能流矢量沿y方向的分量,
Figure RE-GDA0002844072110000112
为成像能流矢量沿z方向的分量,符号“.”为点积。成像能流矢量采用分量形式,可以表示为:
Figure RE-GDA0002844072110000113
在一些实施场景中,在沿着时间反向方向进行弹性波场延拓时,遍历所有的检查点时刻Ti,i=1,2,...,M,求取所有检查点时刻的成像能流密度:
Figure RE-GDA0002844072110000114
附图7a为根据附图4a中1s时刻震源端弹性波质点振动速度矢量场沿x方向分量、附图4b中震源端弹性波质点振动速度矢量场沿z方向分量、附图6a中1s时刻检波端弹性波应力张量沿xx方向的分量、附图6b中1s时刻检波端弹性波应力张量沿xz方向的分量,进行点积构建得到1s时刻成像能流矢量
Figure RE-GDA0002844072110000115
沿x方向分量
Figure RE-GDA0002844072110000116
附图7b为根据附图4a中1s时刻震源端弹性波质点振动速度矢量场沿x 方向分量、附图4b中震源端弹性波质点振动速度矢量场沿z方向分量、附图6b中 1s时刻检波端弹性波应力张量沿xz方向的分量、附图6c中1s时刻检波端弹性波应力张量沿zz方向的分量,进行点积构建得到1s时刻成像能流矢量
Figure RE-GDA0002844072110000117
沿z 方向分量
Figure RE-GDA0002844072110000118
S4:根据构建的能流矢量散度成像条件,对所构建成像能流矢量进行成像,获得弹性波反射能量成像结果。
本发明的一个实施例中,可以利用给定的每个检查点时刻和相邻上一个检查时刻,求取两者之差作为检查点时间间隔ΔTi;利用得到的每一个检查点时刻的成像能流矢量
Figure RE-GDA0002844072110000119
及检查点时间间隔ΔTi,根据构建的能流矢量散度成像条件进行成像,获得最终弹性能量成像结果。
例如一些实施场景中,利用给定的任意一个检查点时刻Ti和相邻上一个检查时刻Ti-1,利用下述公式可以确定检查点时间间隔:
ΔTi=Ti-Ti-1,i=1,2,...,M (9)
其中,Ti为任意一个检查点时刻时间,Ti-1表示相邻上一个检查点时刻时间,两者之差即为第i个检查点时刻处检查点时间间隔ΔTi。在一些实施场景中,对于第一个检查点时刻而言,没有相邻上一个检查点时刻,此时取上一个检查点时刻 Ti-1=0,i=1。在一些实施场景中,如果相邻检查点时间间隔是相等的,那么任意一个检查点时间间隔ΔTi=ΔT,其中ΔT为第一个检查点时刻T1和第二个检查点时刻T2之间的时间间隔,即ΔT=T2-T1
一些实施场景中,针对每一个检查点时刻的成像能流矢量
Figure RE-GDA0002844072110000121
及检查点时间间隔ΔTi,利用构建的能流矢量散度成像条件进行成像。其中,所述构建的能流矢量散度成像条件,可以表示为:
Figure RE-GDA0002844072110000122
其中,Imenergy表示基于能流矢量散度成像条件求取的弹性波反射能量成像结果,i表示检查点时刻序号,共计N个检查点时刻,ishot表示炮点序号,共计Ishot 炮,符号
Figure RE-GDA0002844072110000123
表示散度运算,
Figure RE-GDA0002844072110000124
表示第i个检查点时刻成像能流矢量散度,采用分量形式表示为
Figure RE-GDA0002844072110000125
ΔTi表示检查点时间间隔。
附图8a为利用所构建的所有检查点时刻的成像能流矢量采用能流矢量散度成像条件得到的单炮的弹性能量成像结果,附图8b为利用所构建的所有检查点时刻的成像能流矢量采用能流矢量散度成像条件得到的多炮的弹性能量成像结果。
本发明提供的基于成像能流矢量的弹性逆时偏移方法,通过在弹性波场延拓基础上发展的弹性波正向延拓算子及弹性波反向延拓算子,可以在实现弹性波场正向延拓及反向延拓,获得弹性波质点振动速度矢量场和应力张量场,从而可以有效保持弹性波矢量特征。通过构建的成像能流矢量,只需要消耗用以存储所有检查点时刻的震源端弹性波质点振动速度矢量场和当前检查点时刻的检波端弹性波应力张量场的内存,解决了弹性能量互相关成像条件中内存消耗较高问题。通过利用构建的能流矢量散度成像条件从矢量弹性波场中获得标量化的弹性能量成像结果,提高成像效率。

Claims (5)

1.一种基于成像能流矢量的弹性逆时偏移方法,其特征在于,包括:
基于弹性波正向延拓算子,获得震源端弹性波质点振动速度矢量场;
基于弹性波反向延拓算子,获得检波端弹性波应力张量场;
对所述的震源端质点振动速度矢量场和检波端应力张量场进行点积运算,构建得到各网格节点在检查点时刻的成像能流矢量;
根据能流矢量散度成像条件,对所构建的成像能流矢量进行成像,获得弹性波反射能量成像结果;
所述成像能流矢量的构建方法为:
根据设定的检查点时刻,提取震源端弹性波质点振动速度矢量场
Figure FDA0003612005290000011
根据设定的检查点时刻,提取检波端弹性应力张量场
Figure FDA0003612005290000012
以所述的相同检查点时刻的震源端质点振动速度矢量场
Figure FDA0003612005290000013
和检波端应力张量场
Figure FDA0003612005290000014
为自变量,根据下式进行点积运算构建成像能流矢量:
Figure FDA0003612005290000015
式中,
Figure FDA0003612005290000016
表示各个网格点处第i个检查点时刻Ti的成像能流失量,
Figure FDA0003612005290000017
为成像能流矢量沿x方向的分量,
Figure FDA0003612005290000018
为成像能流矢量沿y方向的分量,
Figure FDA0003612005290000019
为成像能流矢量沿z方向的分量,
Figure FDA00036120052900000110
表示各个检查点处第i个检查点时刻Ti的检波端弹性应力张量场,
Figure FDA00036120052900000111
表示各个检查点处第i个检查点时刻Ti的震源端弹性波质点振动速度矢量场;
成像能流矢量采用分量形式,表示为:
Figure FDA0003612005290000021
保存每一检查点时刻的成像能流密度矢量;
对所构建的成像能流矢量进行成像的方法包括:
根据设定检查点时刻,确定检查点时间间隔:
ΔTi=Ti-Ti-1,i=1,2,...,M (2)
其中,ΔTi为第i个检查点时刻处检查点时间间隔,Ti为任意一个检查点时刻时间,Ti+1表示相邻上一个检查点时刻时间;
针对每一检查点时刻的成像能流密度矢量及检查点时间间隔,利用如下所示的基于能流矢量的弹性能量成像条件进行成像,获得最终弹性能量成像结果:
Figure FDA0003612005290000022
其中,Imenergy表示基于能流矢量散度成像条件求取的弹性波反射能量成像结果,i表示检查点时刻序号,共计N个检查点时刻,ishot表示炮点序号,共计Ishot炮,符号
Figure FDA0003612005290000023
表示散度运算,
Figure FDA0003612005290000024
表示第i个检查点时刻成像能流矢量散度,采用分量形式表示为
Figure FDA0003612005290000025
ΔTi表示检查点时间间隔。
2.根据权利要求1所述的基于成像能流矢量的弹性逆时偏移方法,其特征在于,利用交错网格有限差分格式对弹性波动方程进行离散,获得弹性波延拓算子。
3.根据权利要求2所述的基于成像能流矢量的弹性逆时偏移方法,其特征在于,对所述的弹性波延拓算子加载给定震源子波及介质模型,沿时间正向传播方向进行弹性波场延拓,获得弹性波正向延拓算子。
4.根据权利要求3所述的基于成像能流矢量的弹性逆时偏移方法,其特征在于,所述的给定震源子波加载到弹性波场延拓算子中弹性波的正应力或剪切应力上。
5.根据权利要求4所述的基于成像能流矢量的弹性逆时偏移方法,其特征在于,对所述的弹性波延拓算子利用预先给定的介质模型,将多分量地震记录作为边界条件,沿时间反向传播方向进行弹性波场延拓获得弹性波反向延拓算子。
CN202011096231.9A 2020-10-14 2020-10-14 一种基于成像能流矢量的弹性逆时偏移方法 Active CN112327359B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011096231.9A CN112327359B (zh) 2020-10-14 2020-10-14 一种基于成像能流矢量的弹性逆时偏移方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011096231.9A CN112327359B (zh) 2020-10-14 2020-10-14 一种基于成像能流矢量的弹性逆时偏移方法

Publications (2)

Publication Number Publication Date
CN112327359A CN112327359A (zh) 2021-02-05
CN112327359B true CN112327359B (zh) 2022-06-14

Family

ID=74313175

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011096231.9A Active CN112327359B (zh) 2020-10-14 2020-10-14 一种基于成像能流矢量的弹性逆时偏移方法

Country Status (1)

Country Link
CN (1) CN112327359B (zh)

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112904426B (zh) * 2021-03-27 2022-09-30 中国石油大学(华东) 一种解耦弹性波逆时偏移方法、系统及应用
CN115469362B (zh) * 2022-09-15 2023-10-10 中山大学 一种地震勘探中的能流密度矢量计算方法

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102156296A (zh) * 2011-04-19 2011-08-17 中国石油大学(华东) 地震多分量联合弹性逆时偏移成像方法
CN109557582A (zh) * 2018-12-17 2019-04-02 中国石油大学(华东) 一种二维多分量地震资料偏移成像方法及系统

Family Cites Families (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5999488A (en) * 1998-04-27 1999-12-07 Phillips Petroleum Company Method and apparatus for migration by finite differences
RU2570825C2 (ru) * 2010-06-02 2015-12-10 Эксонмобил Апстрим Рисерч Компани Эффективное вычисление миграции по волновому уравнению угловых сейсмограмм
KR101219746B1 (ko) * 2010-08-24 2013-01-10 서울대학교산학협력단 탄성 매질에서의 주파수 영역 역시간 구조보정을 이용한 지하구조의 영상화 장치 및 방법
US10520618B2 (en) * 2015-02-04 2019-12-31 ExxohnMobil Upstream Research Company Poynting vector minimal reflection boundary conditions
KR101564094B1 (ko) * 2015-07-02 2015-10-29 한국지질자원연구원 지하구조 영상의 품질 향상을 위해 절대값 함수가 적용된 탄성 매질에서의 역시간 구조보정 장치 및 방법
CN105807315B (zh) * 2016-03-14 2018-03-20 中国石油大学(华东) 弹性矢量逆时偏移成像方法
CN109143339B (zh) * 2018-08-14 2020-06-09 中国石油天然气集团有限公司 基于横波应力不变量的弹性逆时偏移成像方法及装置
CN110857998B (zh) * 2018-08-23 2021-11-05 中国石油化工股份有限公司 一种基于lowrank有限差分的弹性逆时偏移方法及系统
CN111221037B (zh) * 2020-01-21 2021-07-23 中国石油大学(华东) 解耦弹性逆时偏移成像方法和装置
CN111239804B (zh) * 2020-02-12 2021-07-02 中国石油大学(华东) 一种弹性能量逆时偏移成像方法、装置、设备及系统

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102156296A (zh) * 2011-04-19 2011-08-17 中国石油大学(华东) 地震多分量联合弹性逆时偏移成像方法
CN109557582A (zh) * 2018-12-17 2019-04-02 中国石油大学(华东) 一种二维多分量地震资料偏移成像方法及系统

Also Published As

Publication number Publication date
CN112327359A (zh) 2021-02-05

Similar Documents

Publication Publication Date Title
CN108459350B (zh) 一种深度域地震子波提取与地震记录合成的一体化方法
CN112327359B (zh) 一种基于成像能流矢量的弹性逆时偏移方法
JP5379163B2 (ja) 地震探査データのスペクトルシェーピングインバージョン法及びマイグレーション法
EP2497043B1 (en) Seismic imaging systems and methods employing a 3d reverse time migration with tilted transverse isotropy
US9063245B2 (en) Efficient wavefield compression in seismic imaging
AU2010238561B2 (en) Method for full-bandwidth deghosting of marine seismic streamer data
US11360224B2 (en) Inversion, migration, and imaging related to isotropic wave-mode- independent attenuation
EP2330443B1 (en) Method for full-bandwidth source deghosting of marine seismic streamer data
CA2947410A1 (en) Fast viscoacoustic and viscoelastic full-wavefield inversion
Chen et al. A normalized wavefield separation cross-correlation imaging condition for reverse time migration based on Poynting vector
CN112882099B (zh) 一种地震频带拓宽方法、装置、介质及电子设备
CN109143339A (zh) 弹性逆时偏移成像方法及装置
WO2018102813A2 (en) Seismic acquisition geometry full-waveform inversion
Zhao et al. Frequency‐domain double‐plane‐wave least‐squares reverse time migration
CN108445539A (zh) 一种消除地震子波旁瓣干扰的方法、设备及系统
CN111505714B (zh) 基于岩石物理约束的弹性波直接包络反演方法
Ma et al. Stable absorption compensation with lateral constraint
Gao et al. Deep learning vertical resolution enhancement considering features of seismic data
Slob et al. Unified elimination of 1D acoustic multiple reflection
Zwartjes et al. Optimizing reconstruction for sparse spatial sampling
US10338250B2 (en) Method of removing incoherent noise
US20080232195A1 (en) Apparatus and Method for Processing Geophysical Information
Sollberger et al. A reconstruction algorithm for temporally aliased seismic signals recorded by the insight Mars lander
Dondurur An approximation to sparse-spike reflectivity using the gold deconvolution method
AU2015377943A1 (en) Method, system and non-transitory computer-readable medium for forming a seismic image of a geological structure

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