CN111596346A - 弹性波速度反演方法和装置 - Google Patents
弹性波速度反演方法和装置 Download PDFInfo
- Publication number
- CN111596346A CN111596346A CN201910125374.9A CN201910125374A CN111596346A CN 111596346 A CN111596346 A CN 111596346A CN 201910125374 A CN201910125374 A CN 201910125374A CN 111596346 A CN111596346 A CN 111596346A
- Authority
- CN
- China
- Prior art keywords
- elastic wave
- depth
- inversion
- seismic
- parameter
- 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.)
- Granted
Links
- 238000000034 method Methods 0.000 title claims abstract description 55
- 238000003384 imaging method Methods 0.000 claims abstract description 83
- 238000013508 migration Methods 0.000 claims abstract description 79
- 230000005012 migration Effects 0.000 claims abstract description 79
- 238000012545 processing Methods 0.000 claims description 52
- 239000011159 matrix material Substances 0.000 claims description 48
- 238000004364 calculation method Methods 0.000 claims description 29
- 230000006870 function Effects 0.000 claims description 21
- 238000007781 pre-processing Methods 0.000 claims description 19
- 238000004590 computer program Methods 0.000 claims description 17
- 238000005259 measurement Methods 0.000 claims description 16
- 238000000926 separation method Methods 0.000 claims description 11
- 238000004422 calculation algorithm Methods 0.000 claims description 10
- 239000000284 extract Substances 0.000 claims description 4
- 238000000605 extraction Methods 0.000 claims description 4
- 230000008859 change Effects 0.000 abstract description 24
- 238000004519 manufacturing process Methods 0.000 description 19
- 238000010586 diagram Methods 0.000 description 14
- 238000005516 engineering process Methods 0.000 description 10
- 230000008569 process Effects 0.000 description 8
- 238000004891 communication Methods 0.000 description 5
- 238000012937 correction Methods 0.000 description 4
- 230000003287 optical effect Effects 0.000 description 4
- 239000000126 substance Substances 0.000 description 4
- 230000003068 static effect Effects 0.000 description 3
- 238000006243 chemical reaction Methods 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 239000002344 surface layer Substances 0.000 description 2
- 230000005540 biological transmission Effects 0.000 description 1
- 230000001413 cellular effect Effects 0.000 description 1
- 239000012530 fluid Substances 0.000 description 1
- 239000011521 glass Substances 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 239000004973 liquid crystal related substance Substances 0.000 description 1
- 239000000463 material Substances 0.000 description 1
- 238000012067 mathematical method Methods 0.000 description 1
- 230000005055 memory storage Effects 0.000 description 1
- 239000000203 mixture Substances 0.000 description 1
- 239000003208 petroleum Substances 0.000 description 1
- 230000000750 progressive effect Effects 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 239000000523 sample Substances 0.000 description 1
- 239000004065 semiconductor Substances 0.000 description 1
- 238000012546 transfer Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
- G01V1/30—Analysis
-
- Y—GENERAL 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
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
- Y02A90/30—Assessment of water resources
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
技术领域
本发明涉及石油地球物理勘探地震数据处理领域,尤其涉及一种弹性波速度反演方法和装置。
背景技术
油气地震勘探的根本任务是根据观测到的各种信息研究和提取有关地下介质的物性参数(如速度、密度等),并对储层的含油气性做出评价,完成这一任务主要有正演和反演两个途径。反演根据各种地球物理观测数据推测地球内部的结构、形态及物质成分,定量计算各种相关的地球物理参数,应用十分广泛。
油气勘探行业的许多观测结果在一定程度上依赖于反演结果的资料解释,究其原因,对于许多观测资料解释问题,没有一个可将多种观测结果(如有效信号及其能量衰减、噪声、其他各种相关误差等)关联起来的解析解,在这种情况下,必须求助于反演这一数学方法,先估算一个结果,然后对照观测数据检查反演结果,并进行适当修正,最终得到合理答案。
如今,叠前地震反演技术已成为油气勘探的常规技术,并在复杂储层精细预测、储层流体识别等领域展示了良好的应用前景。而叠前深度偏移(Pre-stack depthmigration,PSDM)是实现叠前地震反演的关键技术,PSDM是实现地质构造空间归位的一项处理技术,已知精确速度模型的情况下,叠前深度偏移被认为是精确地获得复杂构造内部映像最有效的手段,是一种真正的全三维叠前成像技术。但是,现有弹性波速度反演技术没有综合考虑光滑约束和非光滑约束对于反演速度的影响,导致反演出的弹性波速度不准确,对于复杂构造或速度横向变化较大的地震资料不能正确成像,不能有效修正陡倾地层和速度变化产生的地下图像畸变,进而降低地震数据解释精度,不能有效指导实际生产。
发明内容
有鉴于此,本发明提供一种弹性波速度反演方法和装置、电子设备、计算机可读存储介质,综合考虑光滑约束和非光滑约束对于反演速度的影响,精确反演出的弹性波速度,对于复杂构造或速度横向变化较大的地震资料能正确成像,能有效修正陡倾地层和速度变化产生的地下图像畸变,进而提高地震数据解释精度,有效指导实际生产,提高了生产效率。
为了实现上述目的,本发明采用如下技术方案:
第一方面,提供一种弹性波速度反演方法,包括:
获取一地质区域内地表的地震弹性波数据的共成像点道集和叠前深度偏移剖面;
从该叠前深度偏移剖面中提取反射界面深度;
根据该共成像点道集和该反射界面深度生成估计偏移深度;
根据该反射界面深度和该估计偏移深度建立弹性波正则化速度反演模型,该弹性波正则化速度反演模型包括:光滑约束项和非光滑约束项;
求解该弹性波正则化速度反演模型,得到该地质区域的弹性波速度。
进一步地,该根据该共成像点道集和该反射界面深度生成估计偏移深度,包括:
对该共成像点道集进行相似性扫描得到弹性波相干性衡量参数;
根据该弹性波相干性衡量参数和该反射界面深度计算该估计偏移深度。
进一步地,该共成像点道集的形状包括双曲部分和非双曲部分,该共成像点道集的横坐标为偏移距,
根据该弹性波相干性衡量参数计算该估计偏移深度时,采用如下公式:
进一步地,该弹性波正则化速度反演模型如下:
其中,J(Δλ)为目标函数,Δλ为模型更新量;A为(M1×M2)行N列的矩阵;矩阵A的元素为偏移深度对模型参数的导数xj为反演参数的横坐标位置,j代表反射界面上的一坐标点;λi为第i个反演参数,i=1表示纵波速度,i=2表示横波速度,i=3表示各向异性参数,该反演参数包括:弹性波速度和各向异性参数;M1为该共成像点道集的偏移距的数量;M2为用于更新反演参数的共成像点道集的数量;P代表纵波;S代表横波;向量b包含M1×M2个元素, 为不同偏移距估计偏移深度的平均值;α,β是正则化参数;|| ||代表求范数,系数γ1、γ2、γ3是相应项的加权因子,矩阵R表示纵波和横波偏移深度对纵横波速度导数的差值,y即纵波和横波偏移深度的差值,
DTV(Δλ)为非光滑约束,是全变分的离散化,可近似为Mζ(Δλ),
其中,dt指与估计速度和各向异性参数相关的参数的均匀变化,其中,ζ为大于零的常数;
Γ(Δλ)为光滑约束,且
其中,负拉普拉斯变换离散形式的矩阵D定义为:
进一步地,该求解该弹性波正则化速度反演模型,包括:
采用非线性化迭代算法,求解该弹性波正则化速度反演模型。
进一步地,采用非线性化迭代算法,求解该弹性波正则化速度反演模型,包括:
步骤1:设置初始模型更新量Δλ0、常数参数r、对称正定矩阵B0、迭代次数k=0;
步骤2:计算g(Δλk),其中,该g(Δλk)为J(Δλ)的梯度函数;
步骤3:判断||g(Δλk)||是否小于ξ×max{1,||g(Δλ0)||},其中ξ<0.1;
若是,执行步骤6,否则,执行步骤4;
J(Δλk+ωkdk)≤J(Δλk)+rωk[g(Δλk)]Tdk,
步骤5:令k=k+1,并利用下式更新Bk,返回步骤2;
其中,I为单位矩阵;
步骤6:停止迭代,得到最终的反演参数。
进一步地,该获取一地质区域内地表的地震弹性波数据的共成像点道集和叠前深度偏移剖面,包括:
获取一地质区域内地表的地震弹性波数据;
对该地震弹性波数据进行叠前深度偏移处理得到共成像点道集和叠前深度偏移剖面。
进一步地,该获取一地质区域内地表的地震弹性波数据,包括:
获取该地质区域内的地震炮集数据;
对该地震炮集数据进行预处理;
通过散度和旋度算子对预处理后的地震炮集数据进行波场分离,得到该地质区域内地表的地震弹性波数据。
进一步地,该对该地震弹性波数据进行叠前深度偏移处理得到共成像点道集和叠前深度偏移剖面,包括:
利用初始速度模型对该地震弹性波数据进行叠前深度偏移处理,得到共成像点道集;
沿着共成像点道集横轴方向进行求和,得到叠前深度偏移剖面。
第二方面,提供一种弹性波速度反演装置,包括:
数据获取模块,获取一地质区域内地表的地震弹性波数据的共成像点道集和叠前深度偏移剖面;
反射界面深度提取模块,从该叠前深度偏移剖面中提取反射界面深度;
估计偏移深度计算模块,根据该共成像点道集和该反射界面深度生成估计偏移深度;
建模模块,根据该反射界面深度和该估计偏移深度建立弹性波正则化速度反演模型,该弹性波正则化速度反演模型包括:光滑约束项和非光滑约束项;
模型求解模块,求解该弹性波正则化速度反演模型,得到该地质区域的弹性波速度。
进一步地,该估计偏移深度计算模块包括:
相似性扫描单元,对该共成像点道集进行相似性扫描得到弹性波相干性衡量参数;
深度计算单元,根据该弹性波相干性衡量参数和该反射界面深度计算该估计偏移深度。
进一步地,该弹性波正则化速度反演模型如下:
其中,J(Δλ)为目标函数,Δλ为模型更新量;A为(M1×M2)行N列的矩阵;矩阵A的元素为偏移深度对模型参数的导数xj为反演参数的横坐标位置,j代表反射界面上的一坐标点;λi为第i个反演参数,i=1表示纵波速度,i=2表示横波速度,i=3表示各向异性参数,该反演参数包括:弹性波速度和各向异性参数;M1为该共成像点道集的偏移距的数量;M2为用于更新反演参数的共成像点道集的数量;P代表纵波;S代表横波;向量b包含M1×M2个元素,
DTV(Δλ)为非光滑约束,是全变分的离散化,可近似为Mζ(Δλ),
其中,dt指与估计速度和各向异性参数相关的参数的均匀变化,其中,ζ为大于零的常数;
Γ(Δλ)为光滑约束,且
其中,负拉普拉斯变换离散形式的矩阵D定义为:
进一步地,该模型求解模块包括:
初始参数设置单元,设置初始模型更新量Δλ0、常数参数r、对称正定矩阵B0、迭代次数k=0;
梯度计算单元,计算g(Δλk),其中,该g(Δλk)为J(Δλ)的梯度函数;
判断单元,判断||g(Δλk)||是否小于ξ×max{1,||g(Δλ0)||},其中ξ<0.1;
J(Δλk+ωkdk)≤J(Δλk)+rωk[g(Δλk)]Tdk,
参数更新单元,令k=k+1,并利用下式更新Bk;
其中,I为单位矩阵;
反演参数输出单元,当判断单元的判断结果为是时,停止迭代,得到最终的反演参数。
进一步地,该数据获取模块包括:
地震弹性波数据获取单元,获取一地质区域内地表的地震弹性波数据;
数据处理单元,对该地震弹性波数据进行叠前深度偏移处理得到共成像点道集和叠前深度偏移剖面。
进一步地,该地震弹性波数据获取单元包括:
勘探数据获取子单元,获取该地质区域内的地震炮集数据;
预处理子单元,对该地震炮集数据进行预处理;
波场分离子单元,通过散度和旋度算子对预处理后的地震炮集数据进行波场分离,得到该地质区域内地表的地震弹性波数据。
进一步地,该数据处理单元包括:
叠前深度偏移处理子单元,利用初始速度模型对该地震弹性波数据进行叠前深度偏移处理,得到共成像点道集;
求和子单元,沿着共成像点道集横轴方向进行求和,得到叠前深度偏移剖面。
第三方面,提供一种电子设备,包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序,该处理器执行该程序时实现:
获取一地质区域内地表的地震弹性波数据的共成像点道集和叠前深度偏移剖面;
从该叠前深度偏移剖面中提取反射界面深度;
根据该共成像点道集和该反射界面深度生成估计偏移深度;
根据该反射界面深度和该估计偏移深度建立弹性波正则化速度反演模型,该弹性波正则化速度反演模型包括:光滑约束项和非光滑约束项;
求解该弹性波正则化速度反演模型,得到该地质区域的弹性波速度。
第四方面,提供一种计算机可读存储介质,其上存储有计算机程序,该计算机程序被处理器执行时实现:
获取一地质区域内地表的地震弹性波数据的共成像点道集和叠前深度偏移剖面;
从该叠前深度偏移剖面中提取反射界面深度;
根据该共成像点道集和该反射界面深度生成估计偏移深度;
根据该反射界面深度和该估计偏移深度建立弹性波正则化速度反演模型,该弹性波正则化速度反演模型包括:光滑约束项和非光滑约束项;
求解该弹性波正则化速度反演模型,得到该地质区域的弹性波速度。
本发明提供的弹性波速度反演方法和装置、电子设备、计算机可读存储介质,获取一地质区域内地表的地震弹性波数据的共成像点道集和叠前深度偏移剖面;从该叠前深度偏移剖面中提取反射界面深度;根据该共成像点道集和该反射界面深度生成估计偏移深度;根据该反射界面深度和该估计偏移深度建立弹性波正则化速度反演模型,该弹性波正则化速度反演模型包括:光滑约束项和非光滑约束项;求解该弹性波正则化速度反演模型,得到该地质区域的弹性波速度。其中,综合考虑了光滑约束和非光滑约束对于反演速度的影响,能够精确反演出的弹性波速度,提高地震成像的分辨率,对于复杂构造或速度横向变化较大的地震资料能正确成像,能有效修正陡倾地层和速度变化产生的地下图像畸变,进而提高地震数据解释精度,有效指导实际生产,提高了生产效率。
为让本发明的上述和其他目的、特征和优点能更明显易懂,下文特举较佳实施例,并配合所附图式,作详细说明如下。
附图说明
为了更清楚地说明本申请实施例或现有技术中的技术方案,下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图是本申请的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其他的附图。在附图中:
图1A为本发明实施例中的服务器S1与客户端设备B1之间的架构示意图;
图1B为本发明实施例中的服务器S1、客户端设备B1及数据库服务器S2之间的架构示意图;
图2是本发明实施例中的弹性波速度反演方法的流程示意图一;
图3示出了图2中步骤S300的具体步骤;
图4示出了图2中步骤S400的具体步骤;
图5示出了图2中步骤S100的具体步骤;
图6示出了图5中步骤S110的具体步骤;
图7示出了图5中步骤S120的具体步骤;
图8是本发明实施例中的弹性波速度反演装置的结构框图一;
图9示出了图8中估计偏移深度计算模块30的具体结构;
图10示出了图8中模型求解模块40的具体结构;
图11示出了图8中数据获取模块10的具体结构;
图12示出了图11中地震弹性波数据获取单元11的具体结构;
图13示出了图11中数据处理单元12的具体结构;
图14为本发明实施例电子设备的结构图。
具体实施方式
为了使本技术领域的人员更好地理解本申请方案,下面将结合本申请实施例中的附图,对本申请实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅仅是本申请一部分的实施例,而不是全部的实施例。基于本申请中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都应当属于本申请保护的范围。
本领域内的技术人员应明白,本发明的实施例可提供为方法、系统、或计算机程序产品。因此,本发明可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且,本发明可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等)上实施的计算机程序产品的形式。
需要说明的是,本申请的说明书和权利要求书及上述附图中的术语“包括”和“具有”以及他们的任何变形,意图在于覆盖不排他的包含,例如,包含了一系列步骤或单元的过程、方法、系统、产品或设备不必限于清楚地列出的那些步骤或单元,而是可包括没有清楚地列出的或对于这些过程、方法、产品或设备固有的其它步骤或单元。
需要说明的是,在不冲突的情况下,本申请中的实施例及实施例中的特征可以相互组合。下面将参考附图并结合实施例来详细说明本申请。
叠前深度偏移(Pre-stack depth migration,PSDM)是实现叠前地震反演的关键技术,PSDM是实现地质构造空间归位的一项处理技术,已知精确速度模型的情况下,叠前深度偏移被认为是精确地获得复杂构造内部映像最有效的手段,是一种真正的全三维叠前成像技术。但是,现有弹性波速度反演技术没有综合考虑光滑约束和非光滑约束对于反演速度的影响,导致反演出的弹性波速度不准确,对于复杂构造或速度横向变化较大的地震资料不能正确成像,不能有效修正陡倾地层和速度变化产生的地下图像畸变,进而降低地震数据解释精度,不能有效指导实际生产。
本发明实施例提供的弹性波速度反演方法和装置、电子设备、计算机可读存储介质,通过综合考虑光滑约束和非光滑约束对于反演速度的影响,能够精确反演出的弹性波速度,对于复杂构造或速度横向变化较大的地震资料能正确成像,能有效修正陡倾地层和速度变化产生的地下图像畸变,进而提高地震数据解释精度,有效指导实际生产,提高了生产效率。
有鉴于此,本申请提供了一种采用基于光滑约束项和非光滑约束项的弹性波正则化速度反演模型的弹性波速度反演装置,该装置可以为一种服务器S1,参见图1A,该服务器S1可以与至少一个客户端设备B1通信连接,所述客户端设备B1可以将地震弹性波数据发送至所述服务器S1,所述服务器S1可以在线接收所述地震弹性波数据。所述服务器S1可以在线或者离线对获取的地震弹性波数据进行预处理,获取一地质区域内地表的地震弹性波数据的共成像点道集和叠前深度偏移剖面;从该叠前深度偏移剖面中提取反射界面深度;根据该共成像点道集和该反射界面深度生成估计偏移深度;根据该反射界面深度和该估计偏移深度建立弹性波正则化速度反演模型,该弹性波正则化速度反演模型包括:光滑约束项和非光滑约束项;求解该弹性波正则化速度反演模型,得到该地质区域的弹性波速度。而后,所述服务器S1可以将该地质区域的弹性波速度在线发送至所述客户端设备B1。所述客户端设备B1可以在线接收该地质区域的弹性波速度。
另外,参见图1B,所述服务器S1还可以与至少一个数据库服务器S2通信连接,所述数据库服务器S2用于存储历史地震波数据与历史地震波速度。所述数据库服务器S2在线将所述历史地震波数据与历史地震波速度发送至所述服务器S1,所述服务器S1可以在线接收历史地震波数据与历史地震波速度,而后根据历史地震波数据与历史地震波速度得到初始速度模型,应用所述初始速度模型对所述地震波数据进行叠前深度偏移处理。
基于上述内容,所述客户端设备B1可以具有显示界面,使得用户能够根据界面查看所述服务器S1发送的该地质区域的弹性波速度。
可以理解的是,所述客户端设备B1可以包括智能手机、平板电子设备、网络机顶盒、便携式计算机、台式电脑、个人数字助理(PDA)、车载设备、智能穿戴设备等。其中,所述智能穿戴设备可以包括智能眼镜、智能手表、智能手环等。
在实际应用中,进行弹性波速度反演的部分可以在如上述内容所述的服务器S1侧执行,即,如图1A所示的架构,也可以所有的操作都在所述客户端设备B1中完成。具体可以根据所述客户端设备B1的处理能力,以及用户使用场景的限制等进行选择。本申请对此不作限定。若所有的操作都在所述客户端设备B1中完成,所述客户端设备B1还可以包括处理器,用于进行弹性波速度反演的具体处理。
所述服务器与所述客户端设备之间可以使用任何合适的网络协议进行通信,包括在本申请提交日尚未开发出的网络协议。所述网络协议例如可以包括TCP/IP协议、UDP/IP协议、HTTP协议、HTTPS协议等。当然,所述网络协议例如还可以包括在上述协议之上使用的RPC协议(Remote Procedure Call Protocol,远程过程调用协议)、REST协议(Representational State Transfer,表述性状态转移协议)等。
本申请采用基于光滑约束项和非光滑约束项的弹性波正则化速度反演模型进行弹性波速度反演,能够精确反演出的弹性波速度,对于复杂构造或速度横向变化较大的地震资料能正确成像,能有效修正陡倾地层和速度变化产生的地下图像畸变,进而提高地震数据解释精度,有效指导实际生产,提高了生产效率。具体通过下述实施例及应用场景进行具体说明。
为了能够精确反演出的弹性波速度,本申请实施例提供一种弹性波速度反演方法,参见图2,所述弹性波速度反演方法具体包括如下内容:
步骤S100:获取一地质区域内地表的地震弹性波数据的共成像点道集和叠前深度偏移剖面。
其中,所述地表的地震波数据包括横波数据和纵波数据。
具体地,可以分别针对横波数据和纵波数据,利用相应的弹性波初始速度模型,进行叠前深度偏移处理,生成横波数据和纵波数据的共成像点道集和叠前深度偏移剖面;所述共成像点道集的横坐标为偏移距;所述共成像点道集的纵坐标为深度;所述共成像点道集的数值为地震波振幅;其中,所述初始弹性波速度模型由常规地震弹性波速度扫描得到的均方根速度,经过时深转换后获得。
通常,地震波是由地震震源向四处传播的振动,指从震源产生向四周辐射的弹性波。地震波按传播方式可分为纵波(也称为P波)、横波(也称为S波)和面波(L波)三种类型;其中,纵波和横波均属于体波。地震发生时,震源区的介质发生急速的破裂和运动,这种扰动构成一个波源,由于地球介质的连续性,波动就向地球内部及表层各处传播开去,形成了连续介质中的弹性波。
步骤S200:从该叠前深度偏移剖面中提取反射界面深度。
其中,所述反射界面深度表征了反射界面位置,反射界面位置即其空间位置坐标。
步骤S300:根据该共成像点道集和该反射界面深度生成估计偏移深度。
其中,可以分别针对横波数据和纵波数据,得到对应的估计偏移深度。
步骤S400:根据该反射界面深度和该估计偏移深度建立弹性波正则化速度反演模型(也称地震波速度反演模型),该弹性波正则化速度反演模型包括:光滑约束项和非光滑约束项。
具体地,该弹性波正则化速度反演模型如下:
其中,J(Δλ)为目标函数,Δλ为模型更新量;A为(M1×M2)行N列的矩阵;矩阵A的元素为偏移深度对模型参数的导数由点到面射线追获得;xj为反演参数的横坐标位置,j代表反射界面上的一坐标点;λi为第i个反演参数,i=1表示纵波速度,i=2表示横波速度,i=3表示各向异性参数,该反演参数包括:弹性波速度和各向异性参数;M1为该共成像点道集的偏移距的数量;M2为用于更新反演参数的共成像点道集的数量;P代表纵波;S代表横波;向量b包含M1×M2个元素, 为不同偏移距估计偏移深度的平均值;α,β是正则化参数;|| ||代表求范数。系数γ1、γ2、γ3是相应项的加权因子,矩阵R表示纵波和横波偏移深度对纵横波速度导数的差值,y即纵波和横波偏移深度的差值,
DTV(Δλ)为非光滑约束,是全变分的离散化,可近似为Mζ(Δλ),
其中,dt指与估计速度和各向异性参数相关的参数的均匀变化,其中,ζ为大于零的常数;
Γ(Δλ)为光滑约束,且
其中,负拉普拉斯变换离散形式的矩阵D定义为:
为了简化符号,令
ψP(Δλ)的梯度表示为ψS(Δλ)的梯度表示为χy(Δλ)的梯度表示为DTV(Δλ)的梯度表示为Γ(Δλ)的梯度表示为J(Δλ)的梯度表示为▽J(Δλ);ψP(Δλ)、ψS(Δλ)、χy(Δλ)、DTV(Δλ)、Γ(Δλ)和J(Δλ)的海塞矩阵分别表示为HM,Δλ、HΓ,Δλ和HJ,Δλ,
采用非线性化迭代算法,求解所述弹性波正则化速度反演模型时,在每次迭代中,更新后的解表示为:
另外,Mζ(Δλ)的离散形式可以表示为:
Mζ(Δλ)的梯度表示为:
其中,
diag(φ′(Δλ))表示n×n的对角矩阵,它的第i个对角元素是φ′((LiΔλ)2),L是n×(n+1)的矩阵,第i行表示为Li,并且
M(Δλ)=htLTdiag(φ′(r))L.
所以,公式为:
因此,J(Δλ)的梯度为:
步骤S500:求解该弹性波正则化速度反演模型,得到该地质区域的弹性波速度。
具体地,通过对该弹性波正则化速度反演模型进行迭代循环求解,能够反演出该地质区域的最优弹性波速度,实现该地质区域内地下介质的地震波模型数据的反演。
其中,所述地下介质的地震波模型数据包括横波速度、纵波速度和各向异性参数。
通过上述技术方案可知,本发明实施例提供的弹性波速度反演方法,综合考虑了光滑约束和非光滑约束对于反演速度的影响,能够精确反演出的弹性波速度,提高地震成像的分辨率,对于复杂构造或速度横向变化较大的地震资料能正确成像,能有效修正陡倾地层和速度变化产生的地下图像畸变,进而提高地震数据解释精度,有效指导实际生产,提高了生产效率。
为了进一步提高反演出的弹性波速度精度,在一个可选的实施例中,参见图3,该步骤S300可以包括以下内容:
步骤S310:对该共成像点道集进行相似性扫描得到弹性波相干性衡量参数r1、r2;
其中,该弹性波相干性衡量参数是一种无量纲常数。
步骤S320:根据该弹性波相干性衡量参数和该反射界面深度计算该估计偏移深度。
具体地,该共成像点道集的形状包括双曲部分和非双曲部分,该共成像点道集的横坐标为偏移距,
根据该弹性波相干性衡量参数计算该估计偏移深度时,采用如下公式:
通过上述技术方案可知,本发明实施例提供的弹性波速度反演方法,能够精确计算估计偏移深度,进而提高地震成像的分辨率。
为了进一步提高反演出的弹性波速度精度,在一个可选的实施例中,参见图4,该步骤S400中,采用非线性化迭代算法,求解该弹性波正则化速度反演模型,具体包括以下内容:
步骤S410:设置初始模型更新量Δλ0、常数参数r、对称正定矩阵B0、迭代次数k=0;
步骤S420:计算g(Δλk),其中,该g(Δλk)为J(Δλ)的梯度函数;
步骤S430:判断||g(Δλk)||是否小于ξ×max{1,||g(Δλ0)||},其中ξ<0.1;
若是,执行步骤S460,否则,执行步骤S440;
J(Δλk+ωkdk)≤J(Δλk)+rωk[g(Δλk)]Tdk,
步骤S450:令k=k+1,并利用下式更新Bk,返回步骤S420;
其中,I为单位矩阵;
步骤S460:停止迭代,得到最终的反演参数。
通过上述技术方案可知,本发明实施例提供的弹性波速度反演方法,通过采用上述的非线性化迭代算法求解该弹性波正则化速度反演模型,能够以最快速度求得最优解,兼顾了计算速度和计算精度,进一步提高了反演出的弹性波速度精度。
为了进一步提高反演出的弹性波速度精度,在一个可选的实施例中,参见图5,该步骤S100可以包括以下内容:
步骤S110:获取一地质区域内地表的地震弹性波数据。
步骤S120:对该地震弹性波数据进行叠前深度偏移处理得到共成像点道集和叠前深度偏移剖面。
其中,在进行叠前深度偏移处理前,可以对地震弹性波数据进行预处理,生成可用于偏移的地震波数据,以有效提高数据精度。
为了进一步提高反演出的弹性波速度精度,在一个可选的实施例中,参见图6,该步骤S110可以包括以下内容:
步骤S111:获取该地质区域内的地震炮集数据。
具体地,利用地球物理方法,采集预设地质区域内的携带有地下地质信息的地震炮集数据;该地震炮集数据通过对三维地震炮集数据进行预处理后获得。其中,所述地下地质信息包括地质构造信息和地质岩性变化信息。
步骤S112:对该地震炮集数据进行预处理。
其中,考虑到通过传感装置采集的初始数据可处理性较差,该预处理包括对地震炮集数据进行噪声去除处理以及对地震炮集数据与预先保存的历史地震数据进行一一对应;进一步地,该预处理还可以包括观测系统加载。
步骤S113:通过散度和旋度算子对预处理后的地震炮集数据进行波场分离,得到该地质区域内地表的地震弹性波数据。
可以理解的是,上述方法通过将获取到的地震炮集数据进行初步处理,提高了数据后续的可处理性。
为了进一步提高反演出的弹性波速度精度,在一个可选的实施例中,参见图7,该步骤S120可以包括以下内容:
步骤S121:利用初始速度模型对该地震弹性波数据进行叠前深度偏移处理,得到共成像点道集。
其中,该初始速度模型根据历史地震波数据与历史地震波速度的对应关系建立。
另外,当不存在历史地震波数据与历史地震波速度时,可以采用常速度模型作为初始速度模型。
同时,对该地震弹性波数据进行叠前深度偏移处理之前,还可以对该地震弹性波数据进行预处理,该预处理包括对地震波数据进行噪声去除、静校正、以及对地震波数据与预先保存的历史地震波数据进行一一对应。
步骤S122:沿着共成像点道集横轴方向进行求和,得到叠前深度偏移剖面。
其中,沿着共成像点道集的横轴方向进行求和的方式也可以称为沿着共成像点道集的横轴方向进行叠加。
在一个可选的实施例中,反演地质区域内地下介质的地震波数据还可以包括如下步骤:
步骤1:对初始速度模型和反射界面位置进行网格离散化处理,生成速度网格点;
步骤2:通过两点射线跟踪的方式,分别计算炮点和检波点位置到反射点的射线;其中,该反射点位于反射界面位置上;该射线的信息至少包括射线走时和传播路径;
步骤3:将射线上的点分配至与点距离最近的速度网格点上;
步骤4:将速度网格点对应的参数设置为地质区域内地下介质的地震波数据;
步骤5:通过地震波速度反演模型,求解地质区域内地下介质的地震波数据。
其中,通过联合纵波和横波的共成像点道集的剩余深度校正量、射线路径网格上的速度值和各向异性参数值也可以建立地震波速度反演模型,通过非线性化迭代求解地震波速度反演模型,得到每个网格点的速度和各向异性参数值。
基于上述内容,本申请提供一种利用本发明实施例弹性波速度反演方法方法实现弹性波速度反演的场景,具体内容如下:
(1)获取地质区域内的地震炮集数据。
(2)对地震炮集数据进行预处理。
(3)通过散度和旋度算子对预处理后的地震炮集数据进行波场分离,得到地质区域内地表的地震弹性波数据。
(4)利用初始速度模型对所述地震弹性波数据进行叠前深度偏移处理,得到共成像点道集。
(5)沿着共成像点道集横轴方向进行求和,得到叠前深度偏移剖面。
(6)从所述叠前深度偏移剖面中提取反射界面深度。
(7)对所述共成像点道集进行相似性扫描得到弹性波相干性衡量参数;
(8)根据所述弹性波相干性衡量参数和所述反射界面深度计算所述估计偏移深度。
其中,共成像点道集的形状包括双曲部分和非双曲部分,所述共成像点道集的横坐标为偏移距,根据所述弹性波相干性衡量参数计算所述估计偏移深度时采用如下公式:
(9)根据所述反射界面深度和所述估计偏移深度建立弹性波正则化速度反演模型,所述弹性波正则化速度反演模型包括:光滑约束项和非光滑约束项;
其中,弹性波正则化速度反演模型如下:
其中,J(Δλ)为目标函数,Δλ为模型更新量;A为(M1×M2)行N列的矩阵;矩阵A的元素为偏移深度对模型参数的导数xj为反演参数的横坐标位置,j代表反射界面上的一坐标点;λi为第i个反演参数,i=1表示纵波速度,i=2表示横波速度,i=3表示各向异性参数,该反演参数包括:弹性波速度和各向异性参数;M1为所述共成像点道集的偏移距的数量;M2为用于更新速度模型的共成像点道集的数量;P代表纵波;S代表横波;向量b包含M1×M2个元素, 为不同偏移距估计偏移深度的平均值;α,β是正则化参数;||||代表求范数。系数γ1、γ2、γ3是相应项的加权因子,矩阵R表示纵波和横波偏移深度对纵横波速度导数的差值,y即纵波和横波偏移深度的差值,
DTV(Δλ)为非光滑约束,是全变分的离散化,可近似为Mζ(Δλ),
其中,dt指与估计速度和各向异性参数相关的参数的均匀变化,其中,ζ为大于零的常数;
Γ(Δλ)为光滑约束,且
其中,负拉普拉斯变换离散形式的矩阵D定义为:
(10)设置初始模型更新量Δλ0、常数参数r、对称正定矩阵B0、迭代次数k=0;
(11)计算g(Δλk),其中,所述g(Δλk)为J(Δλ)的梯度函数;
(12)判断||g(Δλk)||是否小于ξ×max{1,||g(Δλ0)||},其中ξ<0.1;
若是,执行步骤15,否则,执行步骤13;
J(Δλk+ωkdk)≤J(Δλk)+rωk[g(Δλk)]Tdk,
(14)令k=k+1,并利用下式更新Bk,返回步骤11;
其中,I为单位矩阵;
(15)停止迭代,得到最终的反演参数。
基于同一发明构思,本申请实施例还提供了一种弹性波速度反演装置,可以用于实现上述实施例所描述的方法,如下面的实施例所述。由于弹性波速度反演装置解决问题的原理与上述方法相似,因此弹性波速度反演装置的实施可以参见上述方法的实施,重复之处不再赘述。以下所使用的,术语“单元”或者“模块”可以实现预定功能的软件和/或硬件的组合。尽管以下实施例所描述的装置较佳地以软件来实现,但是硬件,或者软件和硬件的组合的实现也是可能并被构想的。
图8是本发明实施例中的弹性波速度反演装置的结构框图一。如图8所示,该弹性波速度反演装置包括:数据获取模块10、反射界面深度提取模块20、估计偏移深度计算模块30、建模模块40以及模型求解模块50。
数据获取模块10获取一地质区域内地表的地震弹性波数据的共成像点道集和叠前深度偏移剖面。
其中,所述地表的地震波数据包括横波数据和纵波数据。
具体地,可以分别针对横波数据和纵波数据,利用相应的弹性波初始速度模型,进行叠前深度偏移处理,生成横波数据和纵波数据的共成像点道集和叠前深度偏移剖面;所述共成像点道集的横坐标为偏移距;所述共成像点道集的纵坐标为深度;所述共成像点道集的数值为地震波振幅;其中,所述初始弹性波速度模型由常规地震弹性波速度扫描得到的均方根速度,经过时深转换后获得。
通常,地震波是由地震震源向四处传播的振动,指从震源产生向四周辐射的弹性波。地震波按传播方式可分为纵波(也称为P波)、横波(也称为S波)和面波(L波)三种类型;其中,纵波和横波均属于体波。地震发生时,震源区的介质发生急速的破裂和运动,这种扰动构成一个波源,由于地球介质的连续性,波动就向地球内部及表层各处传播开去,形成了连续介质中的弹性波。
反射界面深度提取模块20从所述叠前深度偏移剖面中提取反射界面深度。
其中,所述反射界面深度表征了反射界面位置,反射界面位置即其空间位置坐标。
估计偏移深度计算模块30根据所述共成像点道集和所述反射界面深度生成估计偏移深度。
其中,可以分别针对横波数据和纵波数据,得到对应的估计偏移深度。
建模模块40根据所述反射界面深度和所述估计偏移深度建立弹性波正则化速度反演模型,所述弹性波正则化速度反演模型包括:光滑约束项和非光滑约束项。
具体地,该弹性波正则化速度反演模型如下:
其中,J(Δλ)为目标函数,Δλ为模型更新量;A为(M1×M2)行N列的矩阵;矩阵A的元素为偏移深度对模型参数的导数由点到面射线追获得;xj为反演参数的横坐标位置,j代表反射界面上的一坐标点;λi为第i个反演参数,i=1表示纵波速度,i=2表示横波速度,i=3表示各向异性参数,该反演参数包括:弹性波速度和各向异性参数;M1为该共成像点道集的偏移距的数量;M2为用于更新反演参数的共成像点道集的数量;P代表纵波;S代表横波;向量b包含M1×M2个元素, 为不同偏移距估计偏移深度的平均值;α,β是正则化参数;|| ||代表求范数。系数γ1、γ2、γ3是相应项的加权因子,矩阵R表示纵波和横波偏移深度对纵横波速度导数的差值,y即纵波和横波偏移深度的差值,
DTV(Δλ)为非光滑约束,是全变分的离散化,可近似为Mζ(Δλ),
其中,dt指与估计速度和各向异性参数相关的参数的均匀变化,其中,ζ为大于零的常数;
Γ(Δλ)为光滑约束,且
其中,负拉普拉斯变换离散形式的矩阵D定义为:
为了简化符号,令
ψP(Δλ)的梯度表示为ψS(Δλ)的梯度表示为χy(Δλ)的梯度表示为DTV(Δλ)的梯度表示为Γ(Δλ)的梯度表示为J(Δλ)的梯度表示为ψP(Δλ)、ψS(Δλ)、χy(Δλ)、DTV(Δλ)、Γ(Δλ)和J(Δλ)的海塞矩阵分别表示为HM,Δλ、HΓ,Δλ和HJ,Δλ,
采用非线性化迭代算法,求解所述弹性波正则化速度反演模型时,在每次迭代中,更新后的解表示为:
另外,Mζ(Δλ)的离散形式可以表示为:
Mζ(Δλ)的梯度表示为:
其中,
diag(φ′(Δλ))表示n×n的对角矩阵,它的第i个对角元素是φ′((LiΔλ)2),L是n×(n+1)的矩阵,第i行表示为Li,并且
M(Δλ)=htLTdiag(φ′(r))L.
||APΔλ+bP||2=(APΔλ+bP,APΔλ+bP),||ASΔλ+bS||2=(ASΔλ+bS,ASΔλ+bS),
||RΔλ+y||2=(RΔλ+y,RΔλ+y),
所以,公式为:
因此,J(Δλ)的梯度为:
模型求解模块50求解所述弹性波正则化速度反演模型,得到所述地质区域的弹性波速度。
具体地,通过对该弹性波正则化速度反演模型进行迭代循环求解,能够反演出该地质区域的最优弹性波速度,实现该地质区域内地下介质的地震波模型数据的反演。
其中,所述地下介质的地震波模型数据包括横波速度、纵波速度和各向异性参数。
通过上述技术方案可知,本发明实施例提供的弹性波速度反演装置,综合考虑了光滑约束和非光滑约束对于反演速度的影响,能够精确反演出的弹性波速度,提高地震成像的分辨率,对于复杂构造或速度横向变化较大的地震资料能正确成像,能有效修正陡倾地层和速度变化产生的地下图像畸变,进而提高地震数据解释精度,有效指导实际生产,提高了生产效率。
图9示出了图8中估计偏移深度计算模块30的具体结构。如图9所示,该估计偏移深度计算模块30包括:相似性扫描单元31以及深度计算单元32。
相似性扫描单元31对所述共成像点道集进行相似性扫描得到弹性波相干性衡量参数r1、r2。
其中,该弹性波相干性衡量参数是一种无量纲常数。
深度计算单元32根据所述弹性波相干性衡量参数和所述反射界面深度计算所述估计偏移深度。
具体地,该共成像点道集的形状包括双曲部分和非双曲部分,该共成像点道集的横坐标为偏移距,
根据该弹性波相干性衡量参数计算该估计偏移深度时,采用如下公式:
通过上述技术方案可知,本发明实施例提供的弹性波速度反演装置,能够精确计算估计偏移深度,进而提高地震成像的分辨率。
图10示出了图8中模型求解模块40的具体结构。如图10所示,该模型求解模块40包括:初始参数设置单元41、梯度计算单元42、判断单元43、模型更新量计算单元44、参数更新单元45以及反演参数输出单元46。
其中,该模型求解模块40采用非线性化迭代算法进行模型求解。
初始参数设置单元41设置初始模型更新量Δλ0、常数参数r、对称正定矩阵B0、迭代次数k=0;
梯度计算单元42计算g(Δλk),其中,所述g(Δλk)为J(Δλ)的梯度函数;
判断单元43判断||g(Δλk)||是否小于ξ×max{1,||g(Δλ0)||},其中ξ<0.1;
J(Δλk+ωkdk)≤J(Δλk)+rωk[g(Δλk)]Tdk,
参数更新单元45令k=k+1,并利用下式更新Bk;
其中,I为单位矩阵;
反演参数输出单元46当判断单元的判断结果为是时,停止迭代,得到最终的反演参数。
通过上述技术方案可知,本发明实施例提供的弹性波速度反演装置,通过采用上述的非线性化迭代算法求解该弹性波正则化速度反演模型,能够以最快速度求得最优解,兼顾了计算速度和计算精度,进一步提高了反演出的弹性波速度精度。
图11示出了图8中数据获取模块10的具体结构。如图11所示,该数据获取模块10包括:地震弹性波数据获取单元11以及数据处理单元12。
地震弹性波数据获取单元11获取一地质区域内地表的地震弹性波数据。
数据处理单元12对所述地震弹性波数据进行叠前深度偏移处理得到共成像点道集和叠前深度偏移剖面。
其中,在进行叠前深度偏移处理前,可以对地震弹性波数据进行预处理,生成可用于偏移的地震波数据,以有效提高数据精度。
图12示出了图11中地震弹性波数据获取单元11的具体结构。如图12所示,该地震弹性波数据获取单元11包括:勘探数据获取子单元11a、预处理子单元11b以及波场分离子单元11c。
勘探数据获取子单元11a获取所述地质区域内的地震炮集数据。
具体地,利用地球物理方法,采集预设地质区域内的携带有地下地质信息的地震炮集数据;该地震炮集数据通过对三维地震炮集数据进行预处理后获得。其中,所述地下地质信息包括地质构造信息和地质岩性变化信息。
预处理子单元11b对所述地震炮集数据进行预处理。
其中,考虑到通过传感装置采集的初始数据可处理性较差,该预处理包括对地震炮集数据进行噪声去除处理以及对地震炮集数据与预先保存的历史地震数据进行一一对应;进一步地,该预处理还可以包括观测系统加载。
波场分离子单元11c通过散度和旋度算子对预处理后的地震炮集数据进行波场分离,得到所述地质区域内地表的地震弹性波数据。
可以理解的是,上述方法通过将获取到的地震炮集数据进行初步处理,提高了数据后续的可处理性。
图13示出了图11中数据处理单元12的具体结构。如图13所示,该数据处理单元12包括:叠前深度偏移处理子单元12a以及求和子单元12b。
叠前深度偏移处理子单元12a利用初始速度模型对所述地震弹性波数据进行叠前深度偏移处理,得到共成像点道集。
其中,该初始速度模型根据历史地震波数据与历史地震波速度的对应关系建立。
另外,当不存在历史地震波数据与历史地震波速度时,可以采用常速度模型作为初始速度模型。
同时,对该地震弹性波数据进行叠前深度偏移处理之前,还可以对该地震弹性波数据进行预处理,该预处理包括对地震波数据进行噪声去除、静校正、以及对地震波数据与预先保存的历史地震波数据进行一一对应。
求和子单元12b沿着共成像点道集横轴方向进行求和,得到叠前深度偏移剖面。
其中,沿着共成像点道集的横轴方向进行求和的方式也可以称为沿着共成像点道集的横轴方向进行叠加。
上述实施例阐明的装置、模块或单元,具体可以由计算机芯片或实体实现,或者由具有某种功能的产品来实现。一种典型的实现设备为电子设备,具体的,电子设备例如可以为个人计算机、膝上型计算机、蜂窝电话、相机电话、智能电话、个人数字助理、媒体播放器、导航设备、电子邮件设备、游戏控制台、平板计算机、可穿戴设备或者这些设备中的任何设备的组合。
在一个典型的实例中电子设备具体包括存储器、处理器以及存储在存储器上并可在处理器上运行的计算机程序,所述处理器执行所述程序时实现下述步骤:
获取一地质区域内地表的地震弹性波数据的共成像点道集和叠前深度偏移剖面;
从该叠前深度偏移剖面中提取反射界面深度;
根据该共成像点道集和该反射界面深度生成估计偏移深度;
根据该反射界面深度和该估计偏移深度建立弹性波正则化速度反演模型,该弹性波正则化速度反演模型包括:光滑约束项和非光滑约束项;
求解该弹性波正则化速度反演模型,得到该地质区域的弹性波速度。
从上述描述可知,本发明实施例提供的电子设备,可用于弹性波速度反演,综合考虑光滑约束和非光滑约束对于反演速度的影响,精确反演出的弹性波速度,对于复杂构造或速度横向变化较大的地震资料能正确成像,能有效修正陡倾地层和速度变化产生的地下图像畸变,进而提高地震数据解释精度,有效指导实际生产,提高了生产效率。
下面参考图14,其示出了适于用来实现本申请实施例的电子设备600的结构示意图。
如图14所示,电子设备600包括中央处理单元(CPU)601,其可以根据存储在只读存储器(ROM)602中的程序或者从存储部分608加载到随机访问存储器(RAM))603中的程序而执行各种适当的工作和处理。在RAM603中,还存储有系统600操作所需的各种程序和数据。CPU601、ROM602、以及RAM603通过总线604彼此相连。输入/输出(I/O)接口605也连接至总线604。
以下部件连接至I/O接口605:包括键盘、鼠标等的输入部分606;包括诸如阴极射线管(CRT)、液晶显示器(LCD)等以及扬声器等的输出部分607;包括硬盘等的存储部分608;以及包括诸如LAN卡,调制解调器等的网络接口卡的通信部分609。通信部分609经由诸如因特网的网络执行通信处理。驱动器610也根据需要连接至I/O接口606。可拆卸介质611,诸如磁盘、光盘、磁光盘、半导体存储器等等,根据需要安装在驱动器610上,以便于从其上读出的计算机程序根据需要被安装如存储部分608。
特别地,根据本发明的实施例,上文参考流程图描述的过程可以被实现为计算机软件程序。例如,本发明的实施例包括一种计算机可读存储介质,其上存储有计算机程序,该计算机程序被处理器执行时实现下述步骤:
获取一地质区域内地表的地震弹性波数据的共成像点道集和叠前深度偏移剖面;
从该叠前深度偏移剖面中提取反射界面深度;
根据该共成像点道集和该反射界面深度生成估计偏移深度;
根据该反射界面深度和该估计偏移深度建立弹性波正则化速度反演模型,该弹性波正则化速度反演模型包括:光滑约束项和非光滑约束项;
求解该弹性波正则化速度反演模型,得到该地质区域的弹性波速度。
从上述描述可知,本发明实施例提供的计算机可读存储介质,可用于弹性波速度反演,综合考虑光滑约束和非光滑约束对于反演速度的影响,精确反演出的弹性波速度,对于复杂构造或速度横向变化较大的地震资料能正确成像,能有效修正陡倾地层和速度变化产生的地下图像畸变,进而提高地震数据解释精度,有效指导实际生产,提高了生产效率。
在这样的实施例中,该计算机程序可以通过通信部分609从网络上被下载和安装,和/或从可拆卸介质611被安装。
计算机可读介质包括永久性和非永久性、可移动和非可移动媒体可以由任何方法或技术来实现信息存储。信息可以是计算机可读指令、数据结构、程序的模块或其他数据。计算机的存储介质的例子包括,但不限于相变内存(PRAM)、静态随机存取存储器(SRAM)、动态随机存取存储器(DRAM)、其他类型的随机存取存储器(RAM)、只读存储器(ROM)、电可擦除可编程只读存储器(EEPROM)、快闪记忆体或其他内存技术、只读光盘只读存储器(CD-ROM)、数字多功能光盘(DVD)或其他光学存储、磁盒式磁带,磁带磁磁盘存储或其他磁性存储设备或任何其他非传输介质,可用于存储可以被计算设备访问的信息。按照本文中的界定,计算机可读介质不包括暂存电脑可读媒体(transitory media),如调制的数据信号和载波。
为了描述的方便,描述以上装置时以功能分为各种单元分别描述。当然,在实施本申请时可以把各单元的功能在同一个或多个软件和/或硬件中实现。
本发明是参照根据本发明实施例的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器,使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。
这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中,使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品,该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。
这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上,使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理,从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。
还需要说明的是,术语“包括”、“包含”或者其任何其他变体意在涵盖非排他性的包含,从而使得包括一系列要素的过程、方法、商品或者设备不仅包括那些要素,而且还包括没有明确列出的其他要素,或者是还包括为这种过程、方法、商品或者设备所固有的要素。在没有更多限制的情况下,由语句“包括一个……”限定的要素,并不排除在包括所述要素的过程、方法、商品或者设备中还存在另外的相同要素。
本领域技术人员应明白,本申请的实施例可提供为方法、系统或计算机程序产品。因此,本申请可采用完全硬件实施例、完全软件实施例或结合软件和硬件方面的实施例的形式。而且,本申请可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等)上实施的计算机程序产品的形式。
本申请可以在由计算机执行的计算机可执行指令的一般上下文中描述,例如程序模块。一般地,程序模块包括执行特定任务或实现特定抽象数据类型的例程、程序、对象、组件、数据结构等等。也可以在分布式计算环境中实践本申请,在这些分布式计算环境中,由通过通信网络而被连接的远程处理设备来执行任务。在分布式计算环境中,程序模块可以位于包括存储设备在内的本地和远程计算机存储介质中。
本说明书中的各个实施例均采用递进的方式描述,各个实施例之间相同相似的部分互相参见即可,每个实施例重点说明的都是与其他实施例的不同之处。尤其,对于系统实施例而言,由于其基本相似于方法实施例,所以描述的比较简单,相关之处参见方法实施例的部分说明即可。
以上所述仅为本申请的实施例而已,并不用于限制本申请。对于本领域技术人员来说,本申请可以有各种更改和变化。凡在本申请的精神和原理之内所作的任何修改、等同替换、改进等,均应包含在本申请的权利要求范围之内。
Claims (18)
1.一种弹性波速度反演方法,其特征在于,包括:
获取一地质区域内地表的地震弹性波数据的共成像点道集和叠前深度偏移剖面;
从所述叠前深度偏移剖面中提取反射界面深度;
根据所述共成像点道集和所述反射界面深度生成估计偏移深度;
根据所述反射界面深度和所述估计偏移深度建立弹性波正则化速度反演模型,所述弹性波正则化速度反演模型包括:光滑约束项和非光滑约束项;
求解所述弹性波正则化速度反演模型,得到所述地质区域的弹性波速度。
2.根据权利要求1所述的弹性波速度反演方法,其特征在于,所述根据所述共成像点道集和所述反射界面深度生成估计偏移深度,包括:
对所述共成像点道集进行相似性扫描得到弹性波相干性衡量参数;
根据所述弹性波相干性衡量参数和所述反射界面深度计算所述估计偏移深度。
4.根据权利要求3所述的弹性波速度反演方法,其特征在于,所述弹性波正则化速度反演模型如下:
其中,J(Δλ)为目标函数,Δλ为模型更新量;A为(M1×M2)行N列的矩阵;矩阵A的元素为偏移深度对模型参数的导数xj为反演参数的横坐标位置,j代表反射界面上的一坐标点;λi为第i个反演参数,i=1表示纵波速度,i=2表示横波速度,i=3表示各向异性参数,该反演参数包括:弹性波速度和各向异性参数;M1为所述共成像点道集的偏移距的数量;M2为用于更新反演参数的共成像点道集的数量;P代表纵波;S代表横波;向量b包含M1×M2个元素, 为不同偏移距估计偏移深度的平均值;α,β是正则化参数;||||代表求范数,系数γ1、γ2、γ3是相应项的加权因子,矩阵R表示纵波和横波偏移深度对纵横波速度导数的差值,y即纵波和横波偏移深度的差值,
DTV(Δλ)为非光滑约束,是全变分的离散化,可近似为Mζ(Δλ),
其中,dt指与估计速度和各向异性参数相关的参数的均匀变化,其中,ζ为大于零的常数;
Γ(Δλ)为光滑约束,且
其中,负拉普拉斯变换离散形式的矩阵D定义为:
5.根据权利要求4所述的弹性波速度反演方法,其特征在于,所述求解所述弹性波正则化速度反演模型,包括:
采用非线性化迭代算法,求解所述弹性波正则化速度反演模型。
6.根据权利要求5所述的弹性波速度反演方法,其特征在于,采用非线性化迭代算法,求解所述弹性波正则化速度反演模型,包括:
步骤1:设置初始模型更新量Δλ0、常数参数r、对称正定矩阵B0、迭代次数k=0;
步骤2:计算g(Δλk),其中,所述g(Δλk)为J(Δλ)的梯度函数;
步骤3:判断||g(Δλk)||是否小于ξ×max{1,||g(Δλ0)||},其中ξ<0.1;
若是,执行步骤6,否则,执行步骤4;
J(Δλk+ωkdk)≤J(Δλk)+rωk[g(Δλk)]Tdk,
步骤5:令k=k+1,并利用下式更新Bk,返回步骤2;
其中,I为单位矩阵;
步骤6:停止迭代,得到最终的反演参数。
7.根据权利要求1所述的弹性波速度反演方法,其特征在于,所述获取一地质区域内地表的地震弹性波数据的共成像点道集和叠前深度偏移剖面,包括:
获取一地质区域内地表的地震弹性波数据;
对所述地震弹性波数据进行叠前深度偏移处理得到共成像点道集和叠前深度偏移剖面。
8.根据权利要求7所述的弹性波速度反演方法,其特征在于,所述获取一地质区域内地表的地震弹性波数据,包括:
获取所述地质区域内的地震炮集数据;
对所述地震炮集数据进行预处理;
通过散度和旋度算子对预处理后的地震炮集数据进行波场分离,得到所述地质区域内地表的地震弹性波数据。
9.根据权利要求7所述的弹性波速度反演方法,其特征在于,所述对所述地震弹性波数据进行叠前深度偏移处理得到共成像点道集和叠前深度偏移剖面,包括:
利用初始速度模型对所述地震弹性波数据进行叠前深度偏移处理,得到共成像点道集;
沿着共成像点道集横轴方向进行求和,得到叠前深度偏移剖面。
10.一种弹性波速度反演装置,其特征在于,包括:
数据获取模块,获取一地质区域内地表的地震弹性波数据的共成像点道集和叠前深度偏移剖面;
反射界面深度提取模块,从所述叠前深度偏移剖面中提取反射界面深度;
估计偏移深度计算模块,根据所述共成像点道集和所述反射界面深度生成估计偏移深度;
建模模块,根据所述反射界面深度和所述估计偏移深度建立弹性波正则化速度反演模型,所述弹性波正则化速度反演模型包括:光滑约束项和非光滑约束项;
模型求解模块,求解所述弹性波正则化速度反演模型,得到所述地质区域的弹性波速度。
11.根据权利要求10所述的弹性波速度反演装置,其特征在于,所述估计偏移深度计算模块包括:
相似性扫描单元,对所述共成像点道集进行相似性扫描得到弹性波相干性衡量参数;
深度计算单元,根据所述弹性波相干性衡量参数和所述反射界面深度计算所述估计偏移深度。
12.根据权利要求11所述的弹性波速度反演装置,其特征在于,所述弹性波正则化速度反演模型如下:
其中,J(Δλ)为目标函数,Δλ为模型更新量;A为(M1×M2)行N列的矩阵;矩阵A的元素为偏移深度对模型参数的导数xj为反演参数的横坐标位置,j代表反射界面上的一坐标点;λi为第i个反演参数,i=1表示纵波速度,i=2表示横波速度,i=3表示各向异性参数,该反演参数包括:弹性波速度和各向异性参数;M1为所述共成像点道集的偏移距的数量;M2为用于更新反演参数的共成像点道集的数量;P代表纵波;S代表横波;向量b包含M1×M2个元素, 为不同偏移距估计偏移深度的平均值;α,β是正则化参数;||||代表求范数,系数γ1、γ2、γ3是相应项的加权因子,矩阵R表示纵波和横波偏移深度对纵横波速度导数的差值,y即纵波和横波偏移深度的差值,
DTV(Δλ)为非光滑约束,是全变分的离散化,可近似为Mζ(Δλ),
其中,dt指与估计速度和各向异性参数相关的参数的均匀变化,其中,ζ为大于零的常数;
Γ(Δλ)为光滑约束,且
其中,负拉普拉斯变换离散形式的矩阵D定义为:
13.根据权利要求12所述的弹性波速度反演装置,其特征在于,所述模型求解模块包括:
初始参数设置单元,设置初始模型更新量Δλ0、常数参数r、对称正定矩阵B0、迭代次数k=0;
梯度计算单元,计算g(Δλk),其中,所述g(Δλk)为J(Δλ)的梯度函数;
判断单元,判断||g(Δλk)||是否小于ξ×max{1,||g(Δλ0)||},其中ξ<0.1;
J(Δλk+ωkdk)≤J(Δλk)+rωk[g(Δλk)]Tdk,
参数更新单元,令k=k+1,并利用下式更新Bk;
其中,I为单位矩阵;
反演参数输出单元,当判断单元的判断结果为是时,停止迭代,得到最终的反演参数。
14.根据权利要求13所述的弹性波速度反演装置,其特征在于,所述数据获取模块包括:
地震弹性波数据获取单元,获取一地质区域内地表的地震弹性波数据;
数据处理单元,对所述地震弹性波数据进行叠前深度偏移处理得到共成像点道集和叠前深度偏移剖面。
15.根据权利要求14所述的弹性波速度反演装置,其特征在于,所述地震弹性波数据获取单元包括:
勘探数据获取子单元,获取所述地质区域内的地震炮集数据;
预处理子单元,对所述地震炮集数据进行预处理;
波场分离子单元,通过散度和旋度算子对预处理后的地震炮集数据进行波场分离,得到所述地质区域内地表的地震弹性波数据。
16.根据权利要求14所述的弹性波速度反演装置,其特征在于,所述数据处理单元包括:
叠前深度偏移处理子单元,利用初始速度模型对所述地震弹性波数据进行叠前深度偏移处理,得到共成像点道集;
求和子单元,沿着共成像点道集横轴方向进行求和,得到叠前深度偏移剖面。
17.一种电子设备,包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序,其特征在于,所述处理器执行所述程序时实现权利要求1至9任一项所述的弹性波速度反演方法的步骤。
18.一种计算机可读存储介质,其上存储有计算机程序,其特征在于,该计算机程序被处理器执行时实现权利要求1至9任一项所述的弹性波速度反演方法的步骤。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910125374.9A CN111596346B (zh) | 2019-02-20 | 2019-02-20 | 弹性波速度反演方法和装置 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201910125374.9A CN111596346B (zh) | 2019-02-20 | 2019-02-20 | 弹性波速度反演方法和装置 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN111596346A true CN111596346A (zh) | 2020-08-28 |
CN111596346B CN111596346B (zh) | 2023-04-25 |
Family
ID=72181298
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201910125374.9A Active CN111596346B (zh) | 2019-02-20 | 2019-02-20 | 弹性波速度反演方法和装置 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN111596346B (zh) |
Citations (10)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20100135115A1 (en) * | 2008-12-03 | 2010-06-03 | Chevron U.S.A. Inc. | Multiple anisotropic parameter inversion for a tti earth model |
CN102176053A (zh) * | 2011-01-27 | 2011-09-07 | 中国科学院地质与地球物理研究所 | 提升波动方程叠前深度偏移成像效果的方法 |
CN102841375A (zh) * | 2012-09-06 | 2012-12-26 | 中国石油大学(华东) | 一种复杂条件下基于角度域共成像点道集的层析速度反演方法 |
CN102841376A (zh) * | 2012-09-06 | 2012-12-26 | 中国石油大学(华东) | 一种基于起伏地表的层析速度反演方法 |
CN104749631A (zh) * | 2015-03-11 | 2015-07-01 | 中国科学院地质与地球物理研究所 | 一种基于稀疏反演的偏移速度分析方法及装置 |
CN106959467A (zh) * | 2017-03-20 | 2017-07-18 | 中国科学院地质与地球物理研究所 | 地震波速度反演方法和装置 |
US9720117B1 (en) * | 2014-06-30 | 2017-08-01 | Pivotal Software, Inc. | Imaging subsurface properties using a parallel processing database system |
US20170285194A1 (en) * | 2016-03-30 | 2017-10-05 | Jaewoo Park | 2D Multiline Seismic Reflection Tomography With Seismic-Tie Constraint |
CN108037531A (zh) * | 2017-11-24 | 2018-05-15 | 电子科技大学 | 一种基于广义全变分正则化的地震反演方法及系统 |
CN108333628A (zh) * | 2018-01-17 | 2018-07-27 | 中国石油大学(华东) | 基于正则化约束的弹性波最小二乘逆时偏移方法 |
-
2019
- 2019-02-20 CN CN201910125374.9A patent/CN111596346B/zh active Active
Patent Citations (10)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20100135115A1 (en) * | 2008-12-03 | 2010-06-03 | Chevron U.S.A. Inc. | Multiple anisotropic parameter inversion for a tti earth model |
CN102176053A (zh) * | 2011-01-27 | 2011-09-07 | 中国科学院地质与地球物理研究所 | 提升波动方程叠前深度偏移成像效果的方法 |
CN102841375A (zh) * | 2012-09-06 | 2012-12-26 | 中国石油大学(华东) | 一种复杂条件下基于角度域共成像点道集的层析速度反演方法 |
CN102841376A (zh) * | 2012-09-06 | 2012-12-26 | 中国石油大学(华东) | 一种基于起伏地表的层析速度反演方法 |
US9720117B1 (en) * | 2014-06-30 | 2017-08-01 | Pivotal Software, Inc. | Imaging subsurface properties using a parallel processing database system |
CN104749631A (zh) * | 2015-03-11 | 2015-07-01 | 中国科学院地质与地球物理研究所 | 一种基于稀疏反演的偏移速度分析方法及装置 |
US20170285194A1 (en) * | 2016-03-30 | 2017-10-05 | Jaewoo Park | 2D Multiline Seismic Reflection Tomography With Seismic-Tie Constraint |
CN106959467A (zh) * | 2017-03-20 | 2017-07-18 | 中国科学院地质与地球物理研究所 | 地震波速度反演方法和装置 |
CN108037531A (zh) * | 2017-11-24 | 2018-05-15 | 电子科技大学 | 一种基于广义全变分正则化的地震反演方法及系统 |
CN108333628A (zh) * | 2018-01-17 | 2018-07-27 | 中国石油大学(华东) | 基于正则化约束的弹性波最小二乘逆时偏移方法 |
Non-Patent Citations (4)
Title |
---|
CAIXIA YU等: "Double parameterized regularization inversion method for migration velocity analysis in transversely isotropic media with a vertical symmetry axis", 《GEOPHYSICAL PROSPECTING》 * |
YANFEI WANG等: "Hybrid regularization methods for seismic reflectivity inversion", 《INT J GEOMATH》 * |
于彩霞等: "VTI介质快速偏移速度分析", 《石油地球物理勘探》 * |
卢昕婷等: "基于全变分原理的多震源混合数据直接偏移方法", 《地球物理学报》 * |
Also Published As
Publication number | Publication date |
---|---|
CN111596346B (zh) | 2023-04-25 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN105277978B (zh) | 一种确定近地表速度模型的方法及装置 | |
US9075159B2 (en) | System and method for seismic data inversion | |
CN108873066B (zh) | 弹性介质波动方程反射波旅行时反演方法 | |
AU2012268718B2 (en) | System and method for seismic data inversion by non-linear model update | |
US20120316844A1 (en) | System and method for data inversion with phase unwrapping | |
US20120316790A1 (en) | System and method for data inversion with phase extrapolation | |
CN105319589B (zh) | 一种利用局部同相轴斜率的全自动立体层析反演方法 | |
US8731838B2 (en) | Fresnel zone fat ray tomography | |
CN106959467B (zh) | 地震波速度反演方法和装置 | |
CN114167511B (zh) | 一种基于连分式展开向下延拓的位场数据快速反演方法 | |
EP3612863B1 (en) | Post-stack kirchhoff depth de-migration method for tilted transverse isotropic (tti) and heterogeneous media based on ray tracing on migrated data | |
CN104316961A (zh) | 获取风化层的地质参数的方法 | |
CN111665556A (zh) | 地层声波传播速度模型构建方法 | |
CN109581494B (zh) | 叠前偏移方法及装置 | |
CN111596346B (zh) | 弹性波速度反演方法和装置 | |
CN109655888B (zh) | 地震数据处理中光滑浮动基准面的定量选择方法及系统 | |
CN111665550A (zh) | 地下介质密度信息反演方法 | |
CN111665549A (zh) | 地层声波衰减因子反演方法 | |
US12032111B2 (en) | Method and system for faster seismic imaging using machine learning | |
US20220283329A1 (en) | Method and system for faster seismic imaging using machine learning | |
CN113495296B (zh) | 层析静校正量的确定方法、装置、设备及可读存储介质 | |
CN117950045A (zh) | 提高密度和速度反演精度的方法、装置、设备及介质 | |
CN108254787B (zh) | 波的走时获得方法及装置、成像方法及装置 | |
Weidle | Surface wave phase velocity maps from multiscale wave field interpolation | |
CN115032694A (zh) | 一种vsp初至旅行时层析成像方法及系统 |
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 |