CN114647003B - 转换波角度域共成像点道集剩余时差计算方法和设备 - Google Patents

转换波角度域共成像点道集剩余时差计算方法和设备 Download PDF

Info

Publication number
CN114647003B
CN114647003B CN202210272563.0A CN202210272563A CN114647003B CN 114647003 B CN114647003 B CN 114647003B CN 202210272563 A CN202210272563 A CN 202210272563A CN 114647003 B CN114647003 B CN 114647003B
Authority
CN
China
Prior art keywords
wave
point
gather
imaging point
horizontal
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
CN202210272563.0A
Other languages
English (en)
Other versions
CN114647003A (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.)
Sun Yat Sen University
Original Assignee
Sun Yat Sen University
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 Sun Yat Sen University filed Critical Sun Yat Sen University
Priority to CN202210272563.0A priority Critical patent/CN114647003B/zh
Publication of CN114647003A publication Critical patent/CN114647003A/zh
Application granted granted Critical
Publication of CN114647003B publication Critical patent/CN114647003B/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/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. analysis, for interpretation, for correction
    • G01V1/30Analysis
    • G01V1/303Analysis for determining velocity profiles or travel times
    • YGENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02ATECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/30Assessment of water resources

Abstract

本申请属于地震波油气勘探领域,具体涉及一种转换波角度域共成像点道集剩余时差计算方法和设备,其方法包括:S1、获取待勘测区域的地震数据;S2、通过预先建立的剩余时差解析式计算转换波角度域共成像点道集的剩余时差;其中,剩余时差解析式的建立步骤包括:S10、根据斯涅尔定律和地震波走时关系,将得到转换点的炮检距和地震波走时作为已知参数的位置表达式;其中,转换点的位置包括反射界面为水平或倾斜时转换点到检波点的水平距离和转换点深度;S20、基于位置表达式建立反射界面为水平或倾斜时的转换波角度域共成像点道集的剩余时差解析式。本申请的方法可准确地计算转换波角道集剩余时差,从而实现基于角道集的转换波速度精确建模。

Description

转换波角度域共成像点道集剩余时差计算方法和设备
技术领域
本申请属于地震波油气勘探领域,具体涉及一种转换波角度域共成像点道集剩余时差计算方法和设备。
背景技术
随着地震波油气勘探目标复杂程度的不断提高,多分量地震数据得到越来越多的应用。但受限于转换波速度场的精度,多分量地震数据中的转换波叠前深度偏移成像仍没有获得广泛的应用,因此如何获取高精度的转换波速度场是转换波叠前深度偏移成像得以实际应用的关键问题之一。
随着地震勘探技术的发展,作为地震波的速度分析方法,偏移速度分析方法从基于射线理论发展到基于波动方程层析成像理论,后者可以为复杂地质环境提供更为准确的理论基础。在众多的波动方程层析理论中,基于角度域共成像点道集(角道集)的层析方法得到广泛的应用。其中,基于弹性波角道集的层析成像主要包括计算高精度的弹性波角道集和使用准确的剩余时差公式更新速度。
对于计算弹性波角道集而言,目前已有多种理论方法可以计算获得高精度的弹性波角道集;对于使用剩余时差公式而言,目前已获得了较为准确的纵波剩余时差公式。而对于转换波剩余时差公式的推导而言,虽然也有地球物理学者对其进行了相关的研究,例如赵杨等在2020年发表的论文《来自弹性逆时偏移的表面偏移距道集及速度分析》,但这些研究均从偏移距域进行推导,因此获得的转换波角道集剩余时差公式仍然不够准确。
如何准确地计算转换波角道集剩余时差,成为亟待解决的技术问题。
发明内容
(一)要解决的技术问题
鉴于现有技术的上述缺点、不足,本申请提供一种转换波角度域共成像点道集剩余时差计算方法、设备和可读存储介质。
(二)技术方案
为达到上述目的,本申请采用如下技术方案:
第一方面,本申请实施例提供一种转换波角度域共成像点道集剩余时差计算方法,该方法包括:
S1、获取待勘测区域的地震数据,所述地震数据包括目标检波点上检测得到的炮检距和地震波走时;
S2、通过预先建立的剩余时差解析式计算转换波角度域共成像点道集的剩余时差,所述转换波角度域共成像点道集基于所述地震数据生成;其中,所述剩余时差解析式的建立步骤包括:
S10、根据斯涅尔定律和地震波走时关系,将得到转换点的炮检距和地震波走时作为已知参数的位置表达式;其中,转换点的位置包括反射界面为水平或倾斜时转换点到检波点的水平距离和转换点深度;
S20、基于所述位置表达式建立反射界面为水平或倾斜时的转换波角度域共成像点道集的剩余时差解析式。
可选地,当反射界面为水平时,所述位置表达式为:
Figure GDA0004149826310000021
其中,h为炮检距,t为地震波走时,xS为转换点到检波点的水平距离,VP和VS为纵波和横波的速度值,下标P和S表示纵波和横波,z为转换点的深度,a、b、A、p为简化所述位置表达式定义的符号,无实际含义,具体表示如下:
Figure GDA0004149826310000031
可选地,当反射界面为倾斜时,所述位置表达式为:
Figure GDA0004149826310000032
其中,xS为转换点到检波点的水平距离,z为转换点深度,h为炮检距,tT为坐标变换后的地震波走时,xT为坐标变换后转换点到检波点的水平位置,θ为地层倾角,zT为坐标变换后的转换点深度,下标T表示以水平地层为x轴的笛卡尔坐标变换至以倾斜地层为x轴的笛卡尔坐标,α为地震波入射角,β为反射角,d为笛卡尔坐标原点至炮点的水平距离。
可选地,S20包括:
基于反射界面为水平时的位置表达式建立反射界面为水平时的转换波角度域共成像点道集的剩余时差解析式;
基于反射界面为倾斜时的位置表达式建立反射界面为倾斜时的转换波角度域共成像点道集的剩余时差解析式。
可选地,基于反射界面为水平时的位置表达式建立反射界面为水平时的转换波角度域共成像点道集的剩余时差解析式,包括:
基于炮检距、成像点到检波点的水平位置、地震波入射角确定角道集成像点深度计算公式,所述角道集成像点深度计算公式为:
Figure GDA0004149826310000033
其中,z′为成像点深度,x′S为成像点到检波点的水平距离,利用所述转换点位置表达式求得,α为入射角。
基于成像点深度和地下反射界面的深度确定反射界面为水平时的转换波角度域共成像点道集的剩余时差解析式,所述剩余时差解析式为:
Figure GDA0004149826310000041
其中,z0为地下反射界面的深度,Δz为转换波角度域共成像点道集的剩余时差。
可选地,当反射界面为倾斜时,所述剩余时差解析式为:
Figure GDA0004149826310000042
其中,h为炮检距,θ为地层倾角,z0为地下反射界面的深度,Δz为转换波角度域共成像点道集的剩余时差,x′T为坐标变换后成像点到检波点的水平位置,其利用所述转换点位置表达式求得,α为入射角,V′P和V′S分别为角度域共成像点道集偏移速度。
可选地,S1之后,S2之前还包括:基于所述地震数据,通过弹性波叠前深度偏移技术生成反射界面为水平或倾斜时的转换波角度域共成像点道集。
可选地,通过弹性波叠前深度偏移技术生成反射界面为水平时的转换波角度域共成像点道集,包括:
利用弹性波解耦延拓方程构建正向时间延拓的纵波质点振动速度矢量和正向时间延拓的纵波应力场;
基于所述正向时间延拓的纵波质点振动速度矢量和所述纵波应力场构建正向时间延拓的纵波Poynting矢量;
利用所述弹性波解耦延拓方程构建逆向时间延拓的纵波质点振动速度矢量和逆向时间延拓的横波应力场;
基于所述逆向时间延拓的纵波质点振动速度矢量和所述横波应力场构建正向时间延拓的横波Poynting矢量;
基于所述纵波Poynting矢量和所述横波Poynting矢量确定转换波角度域共成像点道集中的入射角;
基于所述纵波质点振动速度矢量、所述纵波质点振动速度矢量和所述入射角,利用以下公式获得水平反射界面模型的转换波角度域共成像点道集:
Figure GDA0004149826310000051
其中,
Figure GDA0004149826310000052
为正向时间延拓的纵波质点振动速度场沿x方向的分量,
Figure GDA0004149826310000053
为逆向时间延拓的横波质点振动速度场沿x方向的分量,
Figure GDA0004149826310000054
为正向时间延拓的纵波质点振动速度场沿z方向的分量,
Figure GDA0004149826310000055
为逆向时间延拓的横波总质点振动速度场沿z方向的分量,α为入射角,αk为离散角,ishot为炮序列号,shot为最大炮数,Tmax为最大偏移时间。
第二方面,本申请实施例提供一种电子设备,包括:存储器、处理器及存储在所述存储器上并可在所述处理器上运行的计算机程序,所述计算机程序被所述处理器执行时实现如上第一方面任一项所述的转换波角度域共成像点道集剩余时差计算方法的步骤。
第三方面,本申请实施例提供一种计算机可读存储介质,所述计算机可读存储介质上存储有计算机程序,所述计算机程序被处理器执行时实现如上第一方面任一项所述的转换波角度域共成像点道集剩余时差计算方法的步骤。
(三)有益效果
本申请的有益效果是:本申请提出了一种转换波角度域共成像点道集剩余时差计算方法、设备和可读存储介质,其中的方法包括:S1、获取待勘测区域的地震数据;S2、通过预先建立的剩余时差解析式计算转换波角度域共成像点道集的剩余时差;其中,剩余时差解析式的建立步骤包括:S10、根据斯涅尔定律和地震波走时关系,将得到炮检距和地震波走时作为已知参数的转换点位置表达式;其中,转换点的位置包括反射界面为水平或倾斜时转换点到检波点的水平距离和转换点深度;S20、基于位置表达式建立反射界面为水平或倾斜时的转换波角度域共成像点道集的剩余时差解析式。本申请的方法可准确地计算转换波角道集剩余时差,从而实现基于角道集的转换波速度精确建模。
附图说明
本申请借助于以下附图进行描述:
图1为本申请一个实施例中的转换波角度域共成像点道集剩余时差计算方法流程示意图;
图2为本申请一个实施例中的水平反射界面射线路径图;
图3为本申请一个实施例中的倾斜反射界面射线路径图;
图4为本申请一个实施例中的弹性介质水平层速度模型图;
图5为本申请一个实施例中的弹性介质倾斜层速度模型图;
图6为本申请一个实施例中的偏移速度等于真实速度的转换波角道集的数值解图;
图7为本申请一个实施例中的在横波速度高于真实横波速度5%的水平反射界面的转换波角道集的数值解图;
图8为本申请一个实施例中的在纵波速度低于真实纵波速度5%的水平反射界面的转换波角道集的数值解图;
图9为本申请一个实施例中的在纵波速度高于真实纵波速度3%且横波速度低于真实横波速度5%的水平反射界面的转换波角道集的数值解图;
图10为本申请一个实施例中的在横波速度高于真实横波速度5%的倾斜反射界面的转换波角道集的数值解图;
图11为本申请一个实施例中的在纵波速度低于真实纵波速度5%的倾斜反射界面的转换波角道集的数值解图;
图12为本申请一个实施例中的在纵波速度高于真实纵波速度5%且横波速度高于真实横波速度5%的倾斜反射界面的转换波角道集的数值解图;
图13为本申请一个实施例中的转换波角道集成像深度解析值与转换波角道集对比图;
图14为本申请再一实施例中的电子设备的架构示意图。
具体实施方式
为了更好的解释本发明,以便于理解,下面结合附图,通过具体实施方式,对本发明作详细描述。可以理解的是,以下所描述的具体的实施例仅仅用于解释相关发明,而非对该发明的限定。另外还需要说明的是,在不冲突的情况下,本申请中的实施例及实施例中的特征可以相互组合;为了便于描述,附图中仅示出了与发明相关的部分。
实施例一
图1为本申请一个实施例中的转换波角度域共成像点道集剩余时差计算方法流程示意图,如图1所示,本实施例的转换波角度域共成像点道集剩余时差计算方法包括:
S1、获取待勘测区域的地震数据,地震数据包括目标检波点上检测得到的炮检距和地震波走时;
S2、通过预先建立的剩余时差解析式计算转换波角度域共成像点道集的剩余时差,转换波角度域共成像点道集基于所述地震数据生成;其中,剩余时差解析式的建立步骤包括:
S10、根据斯涅尔定律和地震波走时关系,将得到转换点的炮检距和地震波走时作为已知参数的位置表达式;其中,转换点的位置包括反射界面为水平或倾斜时转换点到检波点的水平距离和转换点深度;
S20、基于位置表达式建立反射界面为水平或倾斜时的转换波角度域共成像点道集的剩余时差解析式。
本实施例的转换波角度域共成像点道集剩余时差计算方法,准确地计算转换波角道集剩余时差,从而实现基于角道集的转换波速度精确建模,为基于弹性波角道集的层析成像、揭示勘测区域的内部地质构造特征提供保证。
为了更好地理解本发明,以下对本实施例中的各步骤进行展开说明。
本实施例S1中,获取待勘测区域的地震数据记录,提取地震数据记录的地震波走时。采集的地震数据记录中可以包括多个地震道数据的集合,每一地震道数据记录可以通过对地震信号进行采样得到,其中,对地震信号进行采样时,以时间间隔Δt对采样地区的地震信号进行采样,获得多道地震数据,构成地震记录。
本实施例S2中,通过预先建立的剩余时差解析式计算转换波角度域共成像点道集的剩余时差的步骤包括:通过输入采集的地震数据中的炮检距信息和走时信息到转换波角道集剩余时差公式可以解析计算出成像点的深度与角度(解析解),之后可以利用深度与角度的对应关系映射到利用偏移技术求得的角道集中(数值解),从而得到剩余时差。
本实施例中,根据斯涅尔定律和地震波走时关系,将推导得到炮检距和地震波走时作为已知参数的转换点位置表达式,以下对推导过程进行详细介绍。
需要说明的是,下文中提及的角道集是指角度域共成像点道集,为便于描述而将其简称为角道集。
图2为本申请一个实施例中的水平反射界面射线路径图,以图2所示的水平反射界面为例,在地下反射界面为水平时,斯涅尔(Snell)定律及地震波走时关系分别表示为:
Figure GDA0004149826310000081
Figure GDA0004149826310000082
其中,h为炮检距,t为地震波走时,xS为转换点到检波点的水平距离,VP和VS为纵波和横波的速度值,下标P和S表示纵波和横波,z为地表与转换点的距离,即转换点的深度。
联立Snell定律及地震波走时表达式,可以得到关于xS的一元三次方程:
Figure GDA0004149826310000083
其中,
Figure GDA0004149826310000091
a,b,c,d为简化所述转换点位置表达式定义的符号,无实际含义。
由判别式Δ=B2-4AC<0可知xS的一元三次方程有三个不相等的实数根,同时xS需要满足0<xS<h/2。因此,xS的解析表达式为:
Figure GDA0004149826310000092
其中,
Figure GDA0004149826310000093
其中,A,B和p为简化所述转换点位置表达式定义的符号,无实际含义。
其次,利用Snell定律表达式,转换点的深度表示为:
Figure GDA0004149826310000094
综合上述推导,水平反射界面时的转换点位置表达式为:
Figure GDA0004149826310000095
其中,
Figure GDA0004149826310000101
在本申请一实施例中,图3为本申请一个实施例中的倾斜反射界面射线路径图;以图3所示的倾斜反射界面为例,在地下反射界面为倾斜时,Snell定律及地震波走时关系分别表示为:
Figure GDA0004149826310000102
Figure GDA0004149826310000103
其中,h为炮检距,tT为坐标变换后的地震波走时,xT为坐标变换后转换点到检波点的水平位置,zT为坐标变换后的转换点深度,下标T表示坐标变换,VP和VS为纵波和横波的速度值,下标P和S表示纵波和横波,α为地震波入射角,β为反射角。
联立炮检距表达式和地震波走时表达式,可以得到关于xT的表达式:
Figure GDA0004149826310000104
其中,
Figure GDA0004149826310000105
其次,利用Snell定律表达式,变换后转换点的深度表示为:
Figure GDA0004149826310000111
综合上述推导,图3所示的关系,倾斜反射界面时的转换点位置表达式为:
Figure GDA0004149826310000112
其中,
Figure GDA0004149826310000113
h为炮检距,t为地震波走时,θ为地层倾角,hT为坐标变换后的炮检距,tT为坐标变换后的地震波走时,xS为转换点到检波点的水平位置,xT为坐标变换后转换点到检波点的水平位置,d为笛卡尔坐标原点至炮点的水平距离,下标T表示以水平地层为x轴的笛卡尔坐标变换至以倾斜地层为x轴的笛卡尔坐标,VP和VS为纵波和横波的速度值,下标P和S表示纵波和横波。m,n,u,v,D,E和q为简化所述转换点位置表达式定义的符号,无实际含义。
本实施例中,S20包括:
基于反射界面为水平时的位置表达式建立反射界面为水平时的转换波角道集的剩余时差解析式;
基于反射界面为倾斜时的位置表达式建立反射界面为倾斜时的转换波角道集的剩余时差解析式。
以下分别对基于位置表达式建立反射界面为水平和倾斜时的转换波角道集的剩余时差解析式的步骤进行具体说明。
在本申请一实施例中,以图2所示的水平反射界面为例,在地下反射界面为水平时,地震波入射角可以表示为:
Figure GDA0004149826310000121
因此,在从地震数据中获取炮检距h以及地震波走时t的前提下,使用偏移速度V′P和V′S获得的角道集成像点深度公式为:
Figure GDA0004149826310000122
则水平反射界面的转换波角道集剩余时差公式为:
Figure GDA0004149826310000123
其中,z′为成像点深度,z0为地下反射界面的深度,Δz为转换波角道集的剩余时差。x′S为成像点到检波点的水平位置,利用所述转换点位置表达式求得,α为入射角。
在本申请一实施例中,以图3所示的倾斜反射界面为例,在地下反射界面为倾斜时,地震波入射角可以表示为:
Figure GDA0004149826310000124
因此,在从地震数据中获取炮检距h以及地震波走时t的前提下,使用偏移速度V′P和V′S获得的角道集成像点深度公式为:
Figure GDA0004149826310000125
则水平反射界面的转换波角道集剩余时差公式为:
Figure GDA0004149826310000126
其中,z′为成像点深度,z0为地下反射界面的深度,Δz为转换波角道集的剩余时差。x′S为成像点到检波点的水平位置,利用所述转换点位置表达式求得,α为入射角。
Figure GDA0004149826310000127
因此,在从地震数据中获取炮检距h以及地震波走时t的前提下,使用偏移速度VP′和VS′获得的角道集成像点深度公式为:
Figure GDA0004149826310000131
则倾斜反射界面的转换波角道集剩余时差公式为:
Figure GDA0004149826310000132
其中,z′为成像点深度,z0为地下反射界面的深度,Δz为转换波角道集的剩余时差。x′T为坐标变换后成像点到检波点的水平位置,其利用所述转换点位置表达式求得,α为入射角。
在其他一些可选地实施方式中,S1之后,S2之前还包括:基于所述地震数据,通过弹性波叠前深度偏移技术生成反射界面为水平或倾斜时的转换波角度域共成像点道集。
以下以弹性逆时偏移成像为例,对通过弹性波叠前深度偏移技术生成反射界面为水平时的转换波角度域共成像点道集的过程进行具体说明。
首先,以图4所示的弹性介质水平层速度模型为例,在二维情况下,可利用如下的弹性波解耦延拓方程构建正向时间延拓的纵波质点振动速度矢量:
Figure GDA0004149826310000133
Figure GDA0004149826310000134
Figure GDA0004149826310000141
其中,
Figure GDA0004149826310000142
为正向时间延拓的弹性波总质点振动速度场沿x方向的分量,
Figure GDA0004149826310000143
为正向时间延拓的弹性波总质点振动速度场沿z方向的分量,
Figure GDA0004149826310000144
为正向时间延拓的纵波质点振动速度场沿x方向的分量,
Figure GDA0004149826310000145
为正向时间延拓的纵波质点振动速度场沿z方向的分量,
Figure GDA0004149826310000146
为正向时间延拓的纵波应力场。
Figure GDA0004149826310000147
为正向时间延拓的横波质点振动速度场沿x方向的分量,
Figure GDA0004149826310000148
为正向时间延拓的横波总质点振动速度场沿z方向的分量,
Figure GDA0004149826310000149
为正向时间延拓的横波应力场的xx分量,
Figure GDA00041498263100001410
为正向时间延拓的横波应力场的zz分量,
Figure GDA00041498263100001411
为正向时间延拓的横波应力场的xz分量。ρ表示弹性介质密度,λ和μ表示弹性介质的杨氏模量,其与纵横波速度的关系为
Figure GDA00041498263100001412
Figure GDA00041498263100001413
上标F表示正向时间延拓。
之后,保存正向时间延拓的纵波质点振动速度矢量
Figure GDA00041498263100001414
以及正向时间延拓的纵波应力场
Figure GDA00041498263100001415
因此,正向时间延拓的纵波Poynting矢量为:
Figure GDA00041498263100001416
其中,
Figure GDA00041498263100001417
为正向时间延拓的纵波Poynting矢量沿x方向的分量,
Figure GDA00041498263100001418
为正向时间延拓的纵波Poynting矢量沿z方向的分量。
其次,以图4所示的弹性介质水平层速度模型为例,在二维情况下,可利用如下的弹性波解耦延拓方程构建逆向时间延拓的纵波质点振动速度矢量:
Figure GDA0004149826310000151
Figure GDA0004149826310000152
Figure GDA0004149826310000153
其中,
Figure GDA0004149826310000154
为逆向时间延拓的弹性波总质点振动速度场沿x方向的分量,
Figure GDA0004149826310000155
为逆向时间延拓的弹性波总质点振动速度场沿z方向的分量。
Figure GDA0004149826310000156
为逆向时间延拓的纵波质点振动速度场沿x方向的分量,
Figure GDA0004149826310000157
为逆向时间延拓的纵波质点振动速度场沿z方向的分量,
Figure GDA0004149826310000158
为逆向时间延拓的纵波应力场。
Figure GDA0004149826310000159
为逆向时间延拓的横波质点振动速度场沿x方向的分量,
Figure GDA00041498263100001510
为逆向时间延拓的横波总质点振动速度场沿z方向的分量,
Figure GDA00041498263100001511
为逆向时间延拓的横波应力场的xx分量,
Figure GDA00041498263100001512
为逆向时间延拓的横波应力场的zz分量,
Figure GDA00041498263100001513
为逆向时间延拓的横波应力场的xz分量。ρ表示弹性介质密度,λ和μ表示弹性介质的杨氏模量,其与纵横波速度的关系为
Figure GDA00041498263100001514
Figure GDA00041498263100001515
此时,保存逆向时间延拓的横波质点振动速度矢量
Figure GDA00041498263100001516
以及逆向时间延拓的横波应力场
Figure GDA00041498263100001517
因此,正向时间延拓的横波Poynting矢量为:
Figure GDA0004149826310000161
其中,
Figure GDA0004149826310000162
为逆向时间延拓的横波Poynting矢量沿x方向的分量,
Figure GDA0004149826310000163
为逆向时间延拓的横波Poynting矢量沿z方向的分量。因此,转换波角道集中的入射角α可由下式计算获得:
Figure GDA0004149826310000164
最后,利用以下的转换波角道集数值解的求取公式:
Figure GDA0004149826310000165
获得水平反射界面模型的转换波角道集。
其中,I为转换波角道集,Tmax为偏移时间,ishot为炮序列号,shot为最大炮数,α为入射角,αk为离散角。
需要说明的是,上述的弹性逆时偏移成像仅仅是示例性的说明,并不构成对偏移成像方法的具体限定。
为了验证本申请提出的剩余时差解析式的准确性,基于上述方法,得到如图4所示的水平反射界面模型位于模型4km处的转换波角道集,分别是:
如图6所示的偏移速度等于真实速度的转换波角道集;
如图7所示的在横波速度高于真实横波速度5%的水平反射界面的转换波角道集;
如图8所示的在纵波速度低于真实纵波速度5%的水平反射界面的转换波角道集;
如图9所示的在纵波速度高于真实纵波速度3%且横波速度低于真实横波速度5%的水平反射界面的转换波角道集;
以及如图5所示的倾斜反射界面模型图位于模型3km处的转换波角道集,分别是:
如图10所示的在横波速度高于真实横波速度5%的倾斜反射界面的转换波角道集;
如图11所示的在纵波速度低于真实纵波速度5%的倾斜反射界面的转换波角道集;
如图12所示在纵波速度高于真实纵波速度5%且横波速度高于真实横波速度5%的倾斜反射界面的转换波角道集。
利用剩余时差公式计算解析值,并与所述求取的角道集进行深度上的对比。具体地,将不同偏移速度值带入水平地下反射界面转换波角道集成像深度解析值,并与水平反射界面转换波角道集进行深度对比。
图13为本申请一个实施例中的转换波角道集成像深度解析值与转换波角道集对比图,图13中(a)为与图6所示的转换波角道集的对比图,(b)为与图7所示的转换波角道集的对比图,(c)为与图8所示的转换波角道集的对比图,(d)为与图9所示的转换波角道集的对比图,(e)为与图10所示的转换波角道集的对比图,(f)为与图11所示的转换波角道集的对比图,(g)为与图12所示的转换波角道集的对比图,图中的灰色粗实线为解析解,表示角道集中成像点的深度。
在本申请一示例性实施例中,基于上述方法,计算偏移速度均与真实速度相等时的转换波角道集成像深度解析值,并与图6所示的转换波角道集进行对比,对比结果如图13中(a)所示。
计算纵波速度等于真实纵波速度而横波速度存在误差时的转换波角道集成像深度解析值,并于例如图7所示的转换波角道集进行对比,对比结果如图13中(b)所示;
计算横波速度等于真实横波速度而纵波速度存在误差时的转换波角道集成像深度解析值,并于例如图8所示的转换波角道集进行对比,对比结果如图13中(c)所示;
计算纵波速度和横波速度均不等于真实速度时的转换波角道集成像深度解析值,并于例如图9所示的转换波角道集进行对比,对比结果如图13中(d)所示。
在本申请一示例性实施例中,将不同偏移速度值带入倾斜地下反射界面转换波角道集成像深度公式计算解析值,并与倾斜反射界面转换波角道集进行深度对比。
在本申请一示例性实施例中,基于上述方法,计算纵波速度等于真实纵波速度而横波速度存在误差时的转换波角道集成像深度解析值,并于例如图10所示的转换波角道集进行对比,对比结果如图13中(e)所示;
计算横波速度等于真实横波速度而纵波速度存在误差时的转换波角道集成像深度解析值,并于例如图11所示的转换波角道集进行对比,对比结果如图13中(f)所示;
计算纵波速度和横波速度均不等于真实速度时的转换波角道集成像深度解析值,并于例如图12所示的转换波角道集进行对比,对比结果如图13中(g)所示。
从对比结果可以看出,用解析式求得的解析解可以完全匹配上角道集中成像点的深度。
实施例二
本申请第三方面通过实施例四提供了一种电子设备,包括:存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序,计算机程序被处理器执行时实现如上实施例中任意一项所述的转换波角度域共成像点道集剩余时差计算方法的步骤。
图14为本申请再一实施例中的电子设备的架构示意图。
图14所示的电子设备可包括:至少一个处理器101、至少一个存储器102、至少一个网络接口104和其他的用户接口103。电子设备中的各个组件通过总线系统105耦合在一起。可理解,总线系统105用于实现这些组件之间的连接通信。总线系统105除包括数据总线之外,还包括电源总线、控制总线和状态信号总线。但是为了清楚说明起见,在图14中将各种总线都标为总线系统105。
其中,用户接口103可以包括显示器、键盘或者点击设备(例如,鼠标,轨迹球(trackball)或者触感板等。
可以理解,本实施例中的存储器102可以是易失性存储器或非易失性存储器,或可包括易失性和非易失性存储器两者。其中,非易失性存储器可以是只读存储器(Read-OnlyMemory,ROM)、可编程只读存储器(Programmable ROM,PROM)、可擦除可编程只读存储器(Erasable PROM,EPROM)、电可擦除可编程只读存储器(Electrically EPROM,EEPROM)或闪存。易失性存储器可以是随机存取存储器(Random Access Memory,RAM),其用作外部高速缓存。通过示例性但不是限制性说明,许多形式的RAM可用,例如静态随机存取存储器(Static RAM,SRAM)、动态随机存取存储器(Dynamic RAM,DRAM)、同步动态随机存取存储器(Synchronous DRAM,SDRAM)、双倍数据速率同步动态随机存取存储器(Double Data RateSDRAM,DDRSDRAM)、增强型同步动态随机存取存储器(Enhanced SDRAM,ESDRAM)、同步连接动态随机存取存储器(Synchlink DRAM,SLDRAM)和直接内存总线随机存取存储器(DirectRambus RAM,DRRAM)。本文描述的存储器102旨在包括但不限于这些和任意其它适合类型的存储器。
在一些实施方式中,存储器102存储了如下的元素,可执行单元或者数据结构,或者他们的子集,或者他们的扩展集:操作系统1021和应用程序1022。
其中,操作系统1021,包含各种系统程序,例如框架层、核心库层、驱动层等,用于实现各种基础业务以及处理基于硬件的任务。应用程序1022,包含各种应用程序,用于实现各种应用业务。实现本发明实施例方法的程序可以包含在应用程序1022中。
在本发明实施例中,处理器101通过调用存储器102存储的程序或指令,具体的,可以是应用程序1022中存储的程序或指令,处理器101用于执行第一方面所提供的方法步骤。
上述本发明实施例揭示的方法可以应用于处理器101中,或者由处理器101实现。处理器101可能是一种集成电路芯片,具有信号的处理能力。在实现过程中,上述方法的各步骤可以通过处理器101中的硬件的集成逻辑电路或者软件形式的指令完成。上述的处理器101可以是通用处理器、数字信号处理器、专用集成电路、现成可编程门阵列或者其他可编程逻辑器件、分立门或者晶体管逻辑器件、分立硬件组件。可以实现或者执行本发明实施例中的公开的各方法、步骤及逻辑框图。通用处理器可以是微处理器或者该处理器也可以是任何常规的处理器等。结合本发明实施例所公开的方法的步骤可以直接体现为硬件译码处理器执行完成,或者用译码处理器中的硬件及软件单元组合执行完成。软件单元可以位于随机存储器,闪存、只读存储器,可编程只读存储器或者电可擦写可编程存储器、寄存器等本领域成熟的存储介质中。该存储介质位于存储器102,处理器101读取存储器102中的信息,结合其硬件完成上述方法的步骤。
另外,结合上述实施例中的转换波角度域共成像点道集剩余时差计算方法,本发明实施例可提供一种计算机可读存储介质,计算机可读存储介质上存储有计算机程序,计算机程序被处理器执行时实现如上方法实施例中的任意一种转换波角度域共成像点道集剩余时差计算方法。
应当注意的是,在权利要求中,不应将位于括号之间的任何附图标记理解成对权利要求的限制。词语“包含”不排除存在未列在权利要求中的部件或步骤。位于部件之前的词语“一”或“一个”不排除存在多个这样的部件。本发明可以借助于包括有若干不同部件的硬件以及借助于适当编程的计算机来实现。词语第一、第二、第三等的使用,仅是为了表述方便,而不表示任何顺序。可将这些词语理解为部件名称的一部分。
此外,需要说明的是,在本说明书的描述中,术语“一个实施例”、“一些实施例”、“实施例”、“示例”、“具体示例”或“一些示例”等的描述,是指结合该实施例或示例描述的具体特征、结构、材料或者特点包含于本发明的至少一个实施例或示例中。在本说明书中,对上述术语的示意性表述不必须针对的是相同的实施例或示例。而且,描述的具体特征、结构、材料或者特点可以在任一个或多个实施例或示例中以合适的方式结合。此外,在不相互矛盾的情况下,本领域的技术人员可以将本说明书中描述的不同实施例或示例以及不同实施例或示例的特征进行结合和组合。
尽管已描述了本发明的优选实施例,但本领域的技术人员在得知了基本创造性概念后,则可对这些实施例作出另外的变更和修改。所以,权利要求应该解释为包括优选实施例以及落入本发明范围的所有变更和修改。
显然,本领域的技术人员可以对本发明进行各种修改和变型而不脱离本发明的精神和范围。这样,倘若本发明的这些修改和变型属于本发明权利要求及其等同技术的范围之内,则本发明也应该包含这些修改和变型在内。

Claims (10)

1.一种转换波角度域共成像点道集剩余时差计算方法,其特征在于,该方法包括:
S1、获取待勘测区域的地震数据,所述地震数据包括目标检波点上检测得到的炮检距和地震波走时;
S2、通过预先建立的剩余时差解析式计算转换波角度域共成像点道集的剩余时差,所述转换波角度域共成像点道集基于所述地震数据生成;其中,所述剩余时差解析式的建立步骤包括:
S10、根据斯涅尔定律和地震波走时关系,将得到转换点的炮检距和地震波走时作为已知参数的位置表达式;其中,转换点的位置包括反射界面为水平或倾斜时转换点到检波点的水平距离和转换点深度;
S20、基于所述位置表达式建立反射界面为水平或倾斜时的转换波角度域共成像点道集的剩余时差解析式。
2.根据权利要求1所述的转换波角度域共成像点道集剩余时差计算方法,其特征在于,当反射界面为水平时,所述位置表达式为:
其中,h为炮检距,t为地震波走时,xS为转换点到检波点的水平距离,VP和VS为纵波和横波的速度值,下标P和S表示纵波和横波,z为转换点的深度,a、b、A、p为简化所述位置表达式定义的符号,无实际含义,具体表示如下:
3.根据权利要求1所述的转换波角度域共成像点道集剩余时差计算方法,其特征在于,当反射界面为倾斜时,所述位置表达式为:
其中,xS为转换点到检波点的水平距离,z为转换点深度,h为炮检距,tT为坐标变换后的地震波走时,xT为坐标变换后转换点到检波点的水平位置,θ为地层倾角,zT为坐标变换后的转换点深度,下标T表示以水平地层为x轴的笛卡尔坐标变换至以倾斜地层为x轴的笛卡尔坐标,α为地震波入射角,β为反射角,d为笛卡尔坐标原点至炮点的水平距离。
4.根据权利要求1所述的转换波角度域共成像点道集剩余时差计算方法,其特征在于,S20包括:
基于反射界面为水平时的位置表达式建立反射界面为水平时的转换波角度域共成像点道集的剩余时差解析式;
基于反射界面为倾斜时的位置表达式建立反射界面为倾斜时的转换波角度域共成像点道集的剩余时差解析式。
5.根据权利要求4所述的转换波角度域共成像点道集剩余时差计算方法,其特征在于,基于反射界面为水平时的位置表达式建立反射界面为水平时的转换波角度域共成像点道集的剩余时差解析式,包括:
基于炮检距、成像点到检波点的水平位置、地震波入射角确定角道集成像点深度计算公式,所述角道集成像点深度计算公式为:
其中,z′为成像点深度,x′S为成像点到检波点的水平距离,利用所述转换点位置表达式求得,α为入射角;
基于成像点深度和地下反射界面的深度确定反射界面为水平时的转换波角度域共成像点道集的剩余时差解析式,所述剩余时差解析式为:
其中,z0为地下反射界面的深度,Δz为转换波角度域共成像点道集的剩余时差。
6.根据权利要求4所述的转换波角度域共成像点道集剩余时差计算方法,其特征在于,当反射界面为倾斜时,所述剩余时差解析式为:
其中,h为炮检距,θ为地层倾角,z0为地下反射界面的深度,Δz为转换波角度域共成像点道集的剩余时差,x′T为坐标变换后成像点到检波点的水平位置,其利用所述转换点位置表达式求得,α为入射角,VP′和VS′分别为角度域共成像点道集偏移速度。
7.根据权利要求1所述的转换波角度域共成像点道集剩余时差计算方法,其特征在于,S1之后,S2之前还包括:基于所述地震数据,通过弹性波叠前深度偏移技术生成反射界面为水平或倾斜时的转换波角度域共成像点道集。
8.根据权利要求7所述的转换波角度域共成像点道集剩余时差计算方法,其特征在于,通过弹性波叠前深度偏移技术生成反射界面为水平时的转换波角度域共成像点道集,包括:
利用弹性波解耦延拓方程构建正向时间延拓的纵波质点振动速度矢量和正向时间延拓的纵波应力场;
基于所述正向时间延拓的纵波质点振动速度矢量和所述纵波应力场构建正向时间延拓的纵波Poynting矢量;
利用所述弹性波解耦延拓方程构建逆向时间延拓的纵波质点振动速度矢量和逆向时间延拓的横波应力场;
基于所述逆向时间延拓的纵波质点振动速度矢量和所述横波应力场构建正向时间延拓的横波Poynting矢量;
基于所述纵波Poynting矢量和所述横波Poynting矢量确定转换波角度域共成像点道集中的入射角;
基于所述纵波质点振动速度矢量、所述纵波质点振动速度矢量和所述入射角,利用以下公式获得水平反射界面模型的转换波角度域共成像点道集:
其中,为正向时间延拓的纵波质点振动速度场沿x方向的分量,为逆向时间延拓的横波质点振动速度场沿x方向的分量,为正向时间延拓的纵波质点振动速度场沿z方向的分量,为逆向时间延拓的横波总质点振动速度场沿z方向的分量,α为入射角,αk为离散角,ishot为炮序列号,shot为最大炮数,Tmax为最大偏移时间。
9.一种电子设备,其特征在于,包括:存储器、处理器及存储在所述存储器上并可在所述处理器上运行的计算机程序,所述计算机程序被所述处理器执行时实现如上权利要求1至8任一项所述的转换波角度域共成像点道集剩余时差计算方法的步骤。
10.一种计算机可读存储介质,其特征在于,所述计算机可读存储介质上存储有计算机程序,所述计算机程序被处理器执行时实现如上权利要求1至8任一项所述的转换波角度域共成像点道集剩余时差计算方法的步骤。
CN202210272563.0A 2022-03-18 2022-03-18 转换波角度域共成像点道集剩余时差计算方法和设备 Active CN114647003B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210272563.0A CN114647003B (zh) 2022-03-18 2022-03-18 转换波角度域共成像点道集剩余时差计算方法和设备

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210272563.0A CN114647003B (zh) 2022-03-18 2022-03-18 转换波角度域共成像点道集剩余时差计算方法和设备

Publications (2)

Publication Number Publication Date
CN114647003A CN114647003A (zh) 2022-06-21
CN114647003B true CN114647003B (zh) 2023-05-09

Family

ID=81995642

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210272563.0A Active CN114647003B (zh) 2022-03-18 2022-03-18 转换波角度域共成像点道集剩余时差计算方法和设备

Country Status (1)

Country Link
CN (1) CN114647003B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115951407B (zh) * 2022-09-15 2023-09-29 中山大学 多次波成像角度域共成像点道集计算方法及计算设备

Family Cites Families (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US6546339B2 (en) * 2000-08-07 2003-04-08 3D Geo Development, Inc. Velocity analysis using angle-domain common image gathers
US20100118652A1 (en) * 2008-11-13 2010-05-13 Jorg Friedrich Schneider Determination of depth moveout and of residual radii of curvature in the common angle domain
US9482770B2 (en) * 2011-03-22 2016-11-01 Exxonmobil Upstream Research Company Residual moveout estimation through least squares inversion
US9442207B1 (en) * 2015-06-12 2016-09-13 Chevron U.S.A. Inc. System and method for computing residual moveout from seismic images
CN112462427B (zh) * 2020-11-13 2022-02-18 中国石油大学(华东) 多分量地震资料保幅角度域共成像点道集提取方法及系统
CN112305615B (zh) * 2020-11-13 2022-02-18 中国石油大学(华东) 一种地震资料角度域共成像点道集提取方法及系统

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
角度域弹性波Kirchhoff叠前深度偏移速度分析方法;杜启振;李芳;侯波;秦童;毕丽飞;;地球物理学报;第54卷(第05期);第1327-1339页 *
面向目标的叠前角度道集提取策略;袁江华;刘洪;首皓;秦月霜;;石油物探;第46卷(第04期);第334-338页 *

Also Published As

Publication number Publication date
CN114647003A (zh) 2022-06-21

Similar Documents

Publication Publication Date Title
Taylor et al. Deconvolution with the ℓ 1 norm
AU783065B2 (en) Extraction of P-wave and S-wave velocities from multi- component seismic data by joint velocity inversion
Hale A method for estimating apparent displacement vectors from time-lapse seismic images
CN108983285B (zh) 一种基于矩张量的多种震源波场模拟方法及装置
CN115373024B (zh) 基于地层记录沉降反演被动陆缘地壳结构的方法及装置
US6226595B1 (en) Method and apparatus using multi-target tracking to analyze borehole images and produce sets of tracks and dip data
Delph et al. Constraining crustal properties using receiver functions and the autocorrelation of earthquake‐generated body waves
US11733413B2 (en) Method and system for super resolution least-squares reverse time migration
CN114647003B (zh) 转换波角度域共成像点道集剩余时差计算方法和设备
Wang et al. Seismic velocity inversion transformer
EP3978961B1 (en) System and method for quantitative seismic integration modeling workflow
Wang et al. Full‐waveform inversion of high‐frequency teleseismic body waves based on multiple plane‐wave incidence: Methods and practical applications
Ren et al. Seismic data inversion with acquisition adaptive convolutional neural network for geologic forward prospecting in tunnels
Konuk et al. Tensorial elastodynamics for anisotropic media
CN109188522B (zh) 速度场构建方法及装置
CN115469362B (zh) 一种地震勘探中的能流密度矢量计算方法
US20230184972A1 (en) Method and system for seismic processing using virtual trace bins based on offset attributes and azimuthal attributes
US11940585B2 (en) System and method for estimating one-way propagation operators
Liu et al. Automatic estimation of traveltime parameters in VTI media using similarity-weighted clustering
CN110764146B (zh) 一种基于声波算子的空间互相关弹性波反射波形反演方法
CN109581521B (zh) Tti各向异性的局部层析方法及系统
Yao et al. An empirical attenuation model of the peak ground acceleration (PGA) in the near field of a strong earthquake
Raeisdana et al. Oriented NMO correction of VTI data in the absence of near-offset traces
Cabrera et al. 3-D prestack depth migration: Implementation and case history
CN113607920B (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