CN115267901A - 动态坐标系弹性波逆时偏移方法、电子设备及介质 - Google Patents
动态坐标系弹性波逆时偏移方法、电子设备及介质 Download PDFInfo
- Publication number
- CN115267901A CN115267901A CN202110475407.XA CN202110475407A CN115267901A CN 115267901 A CN115267901 A CN 115267901A CN 202110475407 A CN202110475407 A CN 202110475407A CN 115267901 A CN115267901 A CN 115267901A
- Authority
- CN
- China
- Prior art keywords
- coordinate system
- wave
- frequency
- calculating
- longitudinal
- 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 37
- 238000013508 migration Methods 0.000 title claims abstract description 31
- 230000005012 migration Effects 0.000 title claims abstract description 31
- 238000003384 imaging method Methods 0.000 claims abstract description 57
- 238000004088 simulation Methods 0.000 claims description 28
- 238000006243 chemical reaction Methods 0.000 claims description 20
- 238000001514 detection method Methods 0.000 claims description 16
- 238000004590 computer program Methods 0.000 claims description 11
- 238000004364 calculation method Methods 0.000 abstract description 7
- 238000010586 diagram Methods 0.000 description 10
- 238000013507 mapping Methods 0.000 description 6
- 230000008901 benefit Effects 0.000 description 5
- 238000012545 processing Methods 0.000 description 5
- 230000000694 effects Effects 0.000 description 2
- 238000002474 experimental method Methods 0.000 description 2
- 239000002245 particle Substances 0.000 description 2
- 239000003208 petroleum Substances 0.000 description 2
- 239000004215 Carbon black (E152) Substances 0.000 description 1
- 230000001427 coherent effect Effects 0.000 description 1
- 238000004891 communication Methods 0.000 description 1
- 238000011156 evaluation Methods 0.000 description 1
- 239000012530 fluid Substances 0.000 description 1
- 230000006870 function Effects 0.000 description 1
- 229930195733 hydrocarbon Natural products 0.000 description 1
- 150000002430 hydrocarbons Chemical class 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- 230000008569 process Effects 0.000 description 1
- 230000000644 propagated effect Effects 0.000 description 1
- 230000009466 transformation Effects 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
- G01V1/30—Analysis
- G01V1/306—Analysis for determining physical properties of the subsurface, e.g. impedance, porosity or attenuation profiles
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/60—Analysis
- G01V2210/62—Physical property of subsurface
- G01V2210/624—Reservoir parameters
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
本申请公开了一种动态坐标系弹性波逆时偏移方法、电子设备及介质。该方法可以包括：步骤1：将观测的地震记录分为低、中、高频尺度数据；步骤2：更新低频纵横波速度的成像结果，判断残差是否小于第一阈值，若是，则进入步骤3；步骤3：更新中频纵横波速度的成像结果，判断残差是否小于第二阈值，若是，则进入步骤4；步骤4：更新高频纵横波速度的成像结果，判断残差是否小于第三阈值，若是，则进入步骤5；步骤5：输出低、中、高频尺度数据的纵横波速度的成像结果。本发明将地震记录分为低、中、高频尺度数据，自适应选取不同尺度的网格及坐标系进行反演，提升成像精度，提高计算效率。
Description
技术领域
本发明涉及石油地球物理勘探技术领域，更具体地，涉及一种动态坐标系弹性波逆时偏移方法、电子设备及介质。
背景技术
石油勘探在处理实际问题时会遇到许多复杂问题，实际采集的数据质量不高，计算效率和精度无法保持平衡，考虑到这些，进行地震数据处理时基于声波介质假设考虑会简化问题。若把地下介质看作声波介质，波场的弹性特性被视为相干噪声，而不是一个信息来源。但简化问题的同时又出现了新问题，由于地下介质同时存在两种形式的弹性性质，地震波在中间传播时纵横波同时存在，建立在声波介质假设上的问题只能主要围绕纵波速度和密度进行，从而忽略了横波存在，为此把地下介质看作声波介质是相对不准确的，不符合精细勘探的要求。横波和纵波在精细勘探中有着无法比拟的作用，在介质中同时传播，但他们又有不同的特性，为了勘探进一步深入，要充分利用地震波的各种信息，横波对流体检测方面的优势不可忽略。纵波及横波除了以上用途外，在岩性勘探中还可以用来进行反演多参数，从而得到更多的地球物理参数，同时得到的结果和岩石物理实验比照，两者结合起来可通过经验公式得到介质孔隙度等对储层评价极其重要的参数。综合上述考虑，针对精细岩性油气藏勘探，发展弹性介质假设下的成像方法迫在眉睫。弹性波最小二乘逆时偏移逐渐得到了发展与应用，但存在的最大问题就是计算量巨大的问题，现有的计算机难以对高密度采集的海量地震数据进行成像。剧烈起伏地表给地震成像带来巨大挑战，目前最好的解决方法是采用贴体网格，但贴体网格坐标系下的波动方程增加了更多分量的计算项，进一步增加计算负担。
因此，有必要开发一种基于数据残差驱动的动态坐标系多重网格弹性波最小二乘逆时偏移方法、电子设备及介质。
公开于本发明背景技术部分的信息仅仅旨在加深对本发明的一般背景技术的理解，而不应当被视为承认或以任何形式暗示该信息构成已为本领域技术人员所公知的现有技术。
发明内容
本发明提出了一种动态坐标系弹性波逆时偏移方法、电子设备及介质，其将地震记录分为低、中、高频尺度数据，自适应选取不同尺度的网格及坐标系进行反演，提升成像精度，提高计算效率。
第一方面，本公开实施例提供了一种动态坐标系弹性波逆时偏移方法，包括：
步骤1：将观测的地震记录分为低、中、高频尺度数据；
步骤2：在笛卡尔坐标系下更新纵横波速度的成像结果，计算低频尺度数据的模拟记录与低频地震记录的残差，判断所述残差是否小于第一阈值，若是，则进入步骤3；
步骤3：进行垂向曲网格坐标系转换，在所述垂向曲网格坐标系下更新纵横波速度的成像结果，计算中频尺度数据的模拟记录与中频地震记录的残差，判断所述残差是否小于第二阈值，若是，则进入步骤4；
步骤4：进行贴体网格坐标系转换，在所述贴体网格坐标系下更新纵横波速度的成像结果，计算高频尺度数据的模拟记录与高频地震记录的残差，判断所述残差是否小于第三阈值，若是，则进入步骤5；
步骤5：输出低、中、高频尺度数据的纵横波速度的成像结果。
优选地，所述步骤2包括：
采用大尺度网格，在笛卡尔坐标系下计算正向延拓的低频地震弹性波震源波场与逆时延拓的低频地震弹性波检波点波场；
计算低频尺度拉梅常数的梯度；
利用链式法则计算纵横波速度的梯度与更新步长，更新纵横波速度的成像结果；
基于弹性波伯恩近似计算低频扰动波场及模拟记录；
计算低频尺度数据的模拟记录与低频地震记录的残差，判断所述残差是否小于第一阈值，若是，则进入步骤3，若否，则继续进行步骤2迭代，直至所述残差小于第一阈值。
优选地，所述步骤3包括：
进行垂向曲网格坐标系转换；
采用中尺度网格，在所述垂向曲网格坐标系下计算正向延拓的中频地震弹性波震源波场与逆时延拓的中频地震弹性波检波点波场；
计算中频尺度拉梅常数的梯度；
利用链式法则计算纵横波速度的梯度与更新步长，更新纵横波速度的成像结果；
基于弹性波伯恩近似计算中频扰动波场及模拟记录；
计算中频尺度数据的模拟记录与中频地震记录的残差，判断所述残差是否小于第二阈值，若是，则进入步骤4，若否，则继续进行步骤3迭代，直至所述残差小于第二阈值。
优选地，通过公式(1)进行垂向曲网格坐标系转换：
其中，zi-1(ξ)表示第i层的顶界面的真实高程，zi(ξ)表示第i层的底界面的真实高程，ηi-1(ξ)和ηi(ξ)分别表示相应的垂向曲网格坐标系下第i层顶底界面的网格点数，ηmax是垂向曲网格坐标系的纵轴最大值，(x(ξ,η),z(ξ,η))是基于垂向曲网格坐标系映射的笛卡尔坐标系内网格节点的坐标值，(ξ,η)是垂向曲网格坐标系内网格节点的坐标值。
优选地，所述步骤4包括：
进行所述贴体网格坐标系转换；
采用小尺度网格，在所述贴体网格坐标系下计算正向延拓的高频地震弹性波震源波场与逆时延拓的高频地震弹性波检波点波场；
计算高频尺度拉梅常数的梯度；
利用链式法则计算纵横波速度的梯度与更新步长，更新纵横波速度的成像结果；
基于弹性波伯恩近似计算高频扰动波场及模拟记录；
计算高频尺度数据的模拟记录与高频地震记录的残差，判断所述残差是否小于第三阈值，若是，则进入步骤5，若否，则继续进行步骤4迭代，直至所述残差小于第三阈值。
优选地，通过公式(2)进行贴体网格坐标系转换：
其中，(x,z)是笛卡尔坐标系内网格节点的坐标值，(E,H)是贴体网格坐标系内网格节点的坐标值。
优选地，通过公式(3)计算纵横波速度的梯度：
其中，
为纵波速度的梯度，ρ为密度，vp为纵波速度，
为拉梅常数λ的梯度，
为横波速度的梯度，vs为横波速度，
为拉梅常数μ的梯度，j为频率参量，j＝L,M,H，L为低频尺度数据标识，M为中频尺度数据标识，H为高频尺度数据标识。
优选地，通过公式(4)更新纵横波速度的成像结果：
其中，
为第k次迭代纵波速度分量的成像结果，αvp为纵波速度分量的更新步长，
为纵波速度的梯度，
为第k次迭代横波速度分量的成像结果，αvs为横波速度分量的更新步长，
为横波速度的梯度，j为频率参量，j＝L,M,H，L为低频尺度数据标识，M为中频尺度数据标识，H为高频尺度数据标识。
作为本公开实施例的一种具体实现方式，
第二方面，本公开实施例还提供了一种电子设备，该电子设备包括：
存储器，存储有可执行指令；
处理器，所述处理器运行所述存储器中的所述可执行指令，以实现所述的动态坐标系弹性波逆时偏移方法。
第三方面，本公开实施例还提供了一种计算机可读存储介质，该计算机可读存储介质存储有计算机程序，该计算机程序被处理器执行时实现所述的动态坐标系弹性波逆时偏移方法。
本发明的方法和装置具有其它的特性和优点，这些特性和优点从并入本文中的附图和随后的具体实施方式中将是显而易见的，或者将在并入本文中的附图和随后的具体实施方式中进行详细陈述，这些附图和具体实施方式共同用于解释本发明的特定原理。
附图说明
通过结合附图对本发明示例性实施例进行更详细的描述，本发明的上述以及其它目的、特征和优势将变得更加明显，其中，在本发明示例性实施例中，相同的参考标号通常代表相同部件。
图1示出了根据本发明的一个实施例的笛卡尔坐标系到垂向曲网格坐标系的映射示意图。
图2示出了根据本发明的一个实施例的笛卡尔坐标系到贴体网格坐标系的映射示意图。
图3示出了根据本发明的一个实施例的动态坐标系弹性波逆时偏移方法的步骤的流程图。
图4a与图4b分别示出了根据本发明的一个实施例的加拿大逆掩断层速度模型中纵波速度与横波速度的示意图。
图5a与图5b分别示出了根据本发明的一个实施例的加拿大逆掩断层纵波速度反射系数模型与横波速度反射系数模型的示意图。
图6a与图6b分别示出了根据本发明的一个实施例的低频大尺度网格笛卡尔坐标系下纵波速度与横波速度的成像结果的示意图。
图7a与图7b分别示出了根据本发明的一个实施例的中频中尺度网格笛卡尔坐标系下纵波速度与横波速度的成像结果的示意图。
图8a与图8b分别示出了根据本发明的一个实施例的高频小尺度网格笛卡尔坐标系下纵波速度与横波速度的成像结果的示意图。
具体实施方式
下面将更详细地描述本发明的优选实施方式。虽然以下描述了本发明的优选实施方式，然而应该理解，可以以各种形式实现本发明而不应被这里阐述的实施方式所限制。
本发明提供一种动态坐标系弹性波逆时偏移方法，包括：
步骤1：将观测的地震记录分为低、中、高频尺度数据；
步骤2：在笛卡尔坐标系下更新纵横波速度的成像结果，计算低频尺度数据的模拟记录与低频地震记录的残差，判断残差是否小于第一阈值，若是，则进入步骤3；
步骤3：进行垂向曲网格坐标系转换，在垂向曲网格坐标系下更新纵横波速度的成像结果，计算中频尺度数据的模拟记录与中频地震记录的残差，判断残差是否小于第二阈值，若是，则进入步骤4；
步骤4：进行贴体网格坐标系转换，在贴体网格坐标系下更新纵横波速度的成像结果，计算高频尺度数据的模拟记录与高频地震记录的残差，判断残差是否小于第三阈值，若是，则进入步骤5；
步骤5：输出低、中、高频尺度数据的纵横波速度的成像结果。
在一个示例中，步骤2包括：
采用大尺度网格，在笛卡尔坐标系下计算正向延拓的低频地震弹性波震源波场与逆时延拓的低频地震弹性波检波点波场；
计算低频尺度拉梅常数的梯度；
利用链式法则计算纵横波速度的梯度与更新步长，更新纵横波速度的成像结果；
基于弹性波伯恩近似计算低频扰动波场及模拟记录；
计算低频尺度数据的模拟记录与低频地震记录的残差，判断残差是否小于第一阈值，若是，则进入步骤3，若否，则继续进行步骤2迭代，直至所述残差小于第一阈值。
在一个示例中，步骤3包括：
进行垂向曲网格坐标系转换；
采用中尺度网格，在垂向曲网格坐标系下计算正向延拓的中频地震弹性波震源波场与逆时延拓的中频地震弹性波检波点波场；
计算中频尺度拉梅常数的梯度；
利用链式法则计算纵横波速度的梯度与更新步长，更新纵横波速度的成像结果；
基于弹性波伯恩近似计算中频扰动波场及模拟记录；
计算中频尺度数据的模拟记录与中频地震记录的残差，判断残差是否小于第二阈值，若是，则进入步骤4，若否，则继续进行步骤3迭代，直至所述残差小于第二阈值。
在一个示例中，通过公式(1)进行垂向曲网格坐标系转换：
其中，zi-1(ξ)表示第i层的顶界面的真实高程，zi(ξ)表示第i层的底界面的真实高程，ηi-1(ξ)和ηi(ξ)分别表示相应的垂向曲网格坐标系下第i层顶底界面的网格点数，ηmax是垂向曲网格坐标系的纵轴最大值，(x(ξ,η),z(ξ,η))是基于垂向曲网格坐标系映射的笛卡尔坐标系内网格节点的坐标值，(ξ,η)是垂向曲网格坐标系内网格节点的坐标值。
在一个示例中，步骤4包括：
进行贴体网格坐标系转换；
采用小尺度网格，在贴体网格坐标系下计算正向延拓的高频地震弹性波震源波场与逆时延拓的高频地震弹性波检波点波场；
计算高频尺度拉梅常数的梯度；
利用链式法则计算纵横波速度的梯度与更新步长，更新纵横波速度的成像结果；
基于弹性波伯恩近似计算高频扰动波场及模拟记录；
计算高频尺度数据的模拟记录与高频地震记录的残差，判断残差是否小于第三阈值，若是，则进入步骤5，若否，则继续进行步骤4迭代，直至所述残差小于第三阈值。
在一个示例中，通过公式(2)进行贴体网格坐标系转换：
其中，(x,z)是笛卡尔坐标系内网格节点的坐标值，(E,H)是贴体网格坐标系内网格节点的坐标值。
在一个示例中，通过公式(3)计算纵横波速度的梯度：
其中，
为纵波速度的梯度，ρ为密度，vp为纵波速度，
为拉梅常数λ的梯度，
为横波速度的梯度，vs为横波速度，
为拉梅常数μ的梯度，j为频率参量，j＝L,M,H，L为低频尺度数据标识，M为中频尺度数据标识，H为高频尺度数据标识。
在一个示例中，通过公式(4)更新纵横波速度的成像结果：
其中，
为第k次迭代纵波速度分量的成像结果，αvp为纵波速度分量的更新步长，
为纵波速度的梯度，
为第k次迭代横波速度分量的成像结果，αvs为横波速度分量的更新步长，
为横波速度的梯度，j为频率参量，j＝L,M,H，L为低频尺度数据标识，M为中频尺度数据标识，H为高频尺度数据标识。
具体地，步骤1：输入纵横波偏移速度场、密度场、观测记录、地表高程文件、观测系统文件，将观测的地震记录分为低、中、高频尺度数据。
步骤2：采用大尺度网格，在笛卡尔坐标系下计算正向延拓的低频地震弹性波震源波场与逆时延拓的低频地震弹性波检波点波场；采用一阶速度-应力方程计算笛卡尔坐标系下正向延拓的低频地震弹性波震源波场：
其中，vx和vz是水平分量和垂直分量的质点速度，τxx和τzz为正应力，τxz为切应力，x(x,z)为空间网格坐标，其中x为水平坐标，z为垂直坐标，t为时间，λ和μ是拉梅常数，ρ是密度。
采用笛卡尔坐标系下的伴随状态方程计算逆时延拓的低频地震弹性波检波点波场：
其中，
分别为vx,vz,τxx,τzz,τxz的伴随变量，xr为检波点空间坐标，
和
分别为水平分量和垂直分量的观测低频记录与模拟低频记录的残差。
通过公式(7)计算低频尺度拉梅常数的梯度：
其中，vx和vz是水平分量和垂直分量的质点速度，τxx和τzz为正应力，τxz为切应力。
利用链式法则通过公式(3)计算纵横波速度的梯度，计算纵横波速度的更新步长，通过公式(4)更新纵横波速度的成像结果；
基于弹性波伯恩近似，通过笛卡尔坐标系反偏移公式计算低频扰动波场及模拟记录：
其中，δvx,δvz,δτxx,δτzz,δτxz为vx,vz,τxx,τzz,τxz的扰动波场，m(λ)＝δλ和m(μ)＝δμ分别为λ和μ的参数扰动。
计算低频尺度数据的模拟记录与低频地震记录的残差，判断残差是否小于第一阈值，若是，则进入步骤3，若否，则继续进行步骤2迭代，直至残差小于第一阈值。优选地，第一阈值为0.4。
图1示出了根据本发明的一个实施例的笛卡尔坐标系到垂向曲网格坐标系的映射示意图。
步骤3：通过公式(1)进行垂向曲网格坐标系转换，笛卡尔坐标系到垂向曲网格坐标系的映射如图1所示；
采用中尺度网格，在垂向曲网格坐标系下计算正向延拓的中频地震弹性波震源波场与逆时延拓的中频地震弹性波检波点波场；采用垂向曲网格坐标系一阶速度-应力方程计算正向延拓的中频地震弹性波震源波场：
采用垂向曲网格坐标系下的伴随状态方程计算逆时延拓的中频地震弹性波检波点波场：
通过下式计算中频尺度拉梅常数的梯度：
利用链式法则通过公式(3)计算纵横波速度的梯度，计算纵横波速度的更新步长，通过公式(4)更新纵横波速度的成像结果；
基于弹性波伯恩近似，通过垂向曲网格坐标系反偏移公式计算中频扰动波场及模拟记录：
计算中频尺度数据的模拟记录与中频地震记录的残差，判断残差是否小于第二阈值，若是，则进入步骤4，若否，则继续进行步骤3迭代，直至所述残差小于第二阈值。优选地，第二阈值为0.25。
图2示出了根据本发明的一个实施例的笛卡尔坐标系到贴体网格坐标系的映射示意图。
步骤4：通过公式(2)进行贴体网格坐标系转换，笛卡尔坐标系到贴体网格坐标系的映射如图2所示；
采用小尺度网格，在贴体网格坐标系下计算正向延拓的高频地震弹性波震源波场与逆时延拓的高频地震弹性波检波点波场；采用贴体网格坐标系一阶速度-应力方程计算正向延拓的高频地震弹性波震源波场：
采用贴体坐标系下的伴随状态方程计算逆时延拓的高频地震弹性波检波点波场：
通过下式计算高频尺度拉梅常数的梯度：
利用链式法则通过公式(3)计算纵横波速度的梯度，计算纵横波速度的更新步长，通过公式(4)更新纵横波速度的成像结果；
基于弹性波伯恩近似，通过贴体网格坐标系反偏移公式计算高频扰动波场及模拟记录：
计算高频尺度数据的模拟记录与高频地震记录的残差，判断残差是否小于第三阈值，若是，则进入步骤5，若否，则继续进行步骤4迭代，直至所述残差小于第三阈值。优选地，第三阈值为0.15。
步骤5：输出低、中、高频尺度数据的纵横波速度的成像结果。
本发明还提供一种电子设备，电子设备包括：存储器，存储有可执行指令；处理器，处理器运行存储器中的可执行指令，以实现上述的动态坐标系弹性波逆时偏移方法。
本发明还提供一种计算机可读存储介质，该计算机可读存储介质存储有计算机程序，该计算机程序被处理器执行时实现上述的动态坐标系弹性波逆时偏移方法。
为便于理解本发明实施例的方案及其效果，以下给出三个具体应用示例。本领域技术人员应理解，该示例仅为了便于理解本发明，其任何具体细节并非意在以任何方式限制本发明。
实施例1
图3示出了根据本发明的一个实施例的动态坐标系弹性波逆时偏移方法的步骤的流程图。
如图3所示，该动态坐标系弹性波逆时偏移方法包括：步骤1：将观测的地震记录分为低、中、高频尺度数据；步骤2：在笛卡尔坐标系下更新纵横波速度的成像结果，计算低频尺度数据的模拟记录与低频地震记录的残差，判断残差是否小于第一阈值，若是，则进入步骤3；步骤3：进行垂向曲网格坐标系转换，在垂向曲网格坐标系下更新纵横波速度的成像结果，计算中频尺度数据的模拟记录与中频地震记录的残差，判断残差是否小于第二阈值，若是，则进入步骤4；步骤4：进行贴体网格坐标系转换，在贴体网格坐标系下更新纵横波速度的成像结果，计算高频尺度数据的模拟记录与高频地震记录的残差，判断残差是否小于第三阈值，若是，则进入步骤5；步骤5：输出低、中、高频尺度数据的纵横波速度的成像结果。
图4a与图4b分别示出了根据本发明的一个实施例的加拿大逆掩断层速度模型中纵波速度与横波速度的示意图。
图5a与图5b分别示出了根据本发明的一个实施例的加拿大逆掩断层纵波速度反射系数模型与横波速度反射系数模型的示意图。
采用国际标准的起伏地表加拿大逆掩断层模型验证本方法，图4a、图4b、图5a、图5b在这里用于作为成像对比。
图6a与图6b分别示出了根据本发明的一个实施例的低频大尺度网格笛卡尔坐标系下纵波速度与横波速度的成像结果的示意图。
图7a与图7b分别示出了根据本发明的一个实施例的中频中尺度网格笛卡尔坐标系下纵波速度与横波速度的成像结果的示意图。
图8a与图8b分别示出了根据本发明的一个实施例的高频小尺度网格笛卡尔坐标系下纵波速度与横波速度的成像结果的示意图。
从图6a、图6b、图7a、图7b、图8a、图8b中可以看出，最终的成像结果与真实成像结果同相轴位置一致，低频噪音得到了很好地压制、成像能量均衡，分辨率高、信噪比高，得到了很高质量的成像结果。除此之外，相比于现有方法，本方法计算效率提高了两倍以上。
实施例2
本公开提供一种电子设备包括，该电子设备包括：存储器，存储有可执行指令；处理器，处理器运行存储器中的可执行指令，以实现上述动态坐标系弹性波逆时偏移方法。
根据本公开实施例的电子设备包括存储器和处理器。
该存储器用于存储非暂时性计算机可读指令。具体地，存储器可以包括一个或多个计算机程序产品，该计算机程序产品可以包括各种形式的计算机可读存储介质，例如易失性存储器和/或非易失性存储器。该易失性存储器例如可以包括随机存取存储器(RAM)和/或高速缓冲存储器(cache)等。该非易失性存储器例如可以包括只读存储器(ROM)、硬盘、闪存等。
该处理器可以是中央处理单元(CPU)或者具有数据处理能力和/或指令执行能力的其它形式的处理单元，并且可以控制电子设备中的其它组件以执行期望的功能。在本公开的一个实施例中，该处理器用于运行该存储器中存储的该计算机可读指令。
本领域技术人员应能理解，为了解决如何获得良好用户体验效果的技术问题，本实施例中也可以包括诸如通信总线、接口等公知的结构，这些公知的结构也应包含在本公开的保护范围之内。
有关本实施例的详细说明可以参考前述各实施例中的相应说明，在此不再赘述。
实施例3
本公开实施例提供一种计算机可读存储介质，该计算机可读存储介质存储有计算机程序，该计算机程序被处理器执行时实现所述的动态坐标系弹性波逆时偏移方法。
根据本公开实施例的计算机可读存储介质，其上存储有非暂时性计算机可读指令。当该非暂时性计算机可读指令由处理器运行时，执行前述的本公开各实施例方法的全部或部分步骤。
上述计算机可读存储介质包括但不限于：光存储介质(例如：CD－ROM和DVD)、磁光存储介质(例如：MO)、磁存储介质(例如：磁带或移动硬盘)、具有内置的可重写非易失性存储器的媒体(例如：存储卡)和具有内置ROM的媒体(例如：ROM盒)。
本领域技术人员应理解，上面对本发明的实施例的描述的目的仅为了示例性地说明本发明的实施例的有益效果，并不意在将本发明的实施例限制于所给出的任何示例。
以上已经描述了本发明的各实施例，上述说明是示例性的，并非穷尽性的，并且也不限于所披露的各实施例。在不偏离所说明的各实施例的范围和精神的情况下，对于本技术领域的普通技术人员来说许多修改和变更都是显而易见的。
Claims (10)
1.一种动态坐标系弹性波逆时偏移方法，其特征在于，包括：
步骤1：将观测的地震记录分为低、中、高频尺度数据；
步骤2：在笛卡尔坐标系下更新纵横波速度的成像结果，计算低频尺度数据的模拟记录与低频地震记录的残差，判断所述残差是否小于第一阈值，若是，则进入步骤3；
步骤3：进行垂向曲网格坐标系转换，在所述垂向曲网格坐标系下更新纵横波速度的成像结果，计算中频尺度数据的模拟记录与中频地震记录的残差，判断所述残差是否小于第二阈值，若是，则进入步骤4；
步骤4：进行贴体网格坐标系转换，在所述贴体网格坐标系下更新纵横波速度的成像结果，计算高频尺度数据的模拟记录与高频地震记录的残差，判断所述残差是否小于第三阈值，若是，则进入步骤5；
步骤5：输出低、中、高频尺度数据的纵横波速度的成像结果。
2.根据权利要求1所述的动态坐标系弹性波逆时偏移方法，其中，所述步骤2包括：
采用大尺度网格，在笛卡尔坐标系下计算正向延拓的低频地震弹性波震源波场与逆时延拓的低频地震弹性波检波点波场；
计算低频尺度拉梅常数的梯度；
利用链式法则计算纵横波速度的梯度与更新步长，更新纵横波速度的成像结果；
基于弹性波伯恩近似计算低频扰动波场及模拟记录；
计算低频尺度数据的模拟记录与低频地震记录的残差，判断所述残差是否小于第一阈值，若是，则进入步骤3，若否，则继续进行步骤2迭代，直至所述残差小于第一阈值。
3.根据权利要求1所述的动态坐标系弹性波逆时偏移方法，其中，所述步骤3包括：
进行垂向曲网格坐标系转换；
采用中尺度网格，在所述垂向曲网格坐标系下计算正向延拓的中频地震弹性波震源波场与逆时延拓的中频地震弹性波检波点波场；
计算中频尺度拉梅常数的梯度；
利用链式法则计算纵横波速度的梯度与更新步长，更新纵横波速度的成像结果；
基于弹性波伯恩近似计算中频扰动波场及模拟记录；
计算中频尺度数据的模拟记录与中频地震记录的残差，判断所述残差是否小于第二阈值，若是，则进入步骤4，若否，则继续进行步骤3迭代，直至所述残差小于第二阈值。
4.根据权利要求3所述的动态坐标系弹性波逆时偏移方法，其中，通过公式(1)进行垂向曲网格坐标系转换：
其中，zi-1(ξ)表示第i层的顶界面的真实高程，zi(ξ)表示第i层的底界面的真实高程，ηi-1(ξ)和ηi(ξ)分别表示相应的垂向曲网格坐标系下第i层顶底界面的网格点数，ηmax是垂向曲网格坐标系的纵轴最大值，(x(ξ,η),z(ξ,η))是基于垂向曲网格坐标系映射的笛卡尔坐标系内网格节点的坐标值，(ξ,η)是垂向曲网格坐标系内网格节点的坐标值。
5.根据权利要求1所述的动态坐标系弹性波逆时偏移方法，其中，所述步骤4包括：
进行所述贴体网格坐标系转换；
采用小尺度网格，在所述贴体网格坐标系下计算正向延拓的高频地震弹性波震源波场与逆时延拓的高频地震弹性波检波点波场；
计算高频尺度拉梅常数的梯度；
利用链式法则计算纵横波速度的梯度与更新步长，更新纵横波速度的成像结果；
基于弹性波伯恩近似计算高频扰动波场及模拟记录；
计算高频尺度数据的模拟记录与高频地震记录的残差，判断所述残差是否小于第三阈值，若是，则进入步骤5，若否，则继续进行步骤4迭代，直至所述残差小于第三阈值。
6.根据权利要求5所述的动态坐标系弹性波逆时偏移方法，其中，通过公式(2)进行贴体网格坐标系转换：
其中，(x,z)是笛卡尔坐标系内网格节点的坐标值，(E,H)是贴体网格坐标系内网格节点的坐标值。
7.根据权利要求2、3、5中任意一项所述的动态坐标系弹性波逆时偏移方法，其中，通过公式(3)计算纵横波速度的梯度：
其中，
为纵波速度的梯度，ρ为密度，vp为纵波速度，
为拉梅常数λ的梯度，
为横波速度的梯度，vs为横波速度，
为拉梅常数μ的梯度，j为频率参量，j＝L,M,H，L为低频尺度数据标识，M为中频尺度数据标识，H为高频尺度数据标识。
8.根据权利要求2、3、5中任意一项所述的动态坐标系弹性波逆时偏移方法，其中，通过公式(4)更新纵横波速度的成像结果：
其中，
为第k次迭代纵波速度分量的成像结果，αvp为纵波速度分量的更新步长，
为纵波速度的梯度，
为第k次迭代横波速度分量的成像结果，αvs为横波速度分量的更新步长，
为横波速度的梯度，j为频率参量，j＝L,M,H，L为低频尺度数据标识，M为中频尺度数据标识，H为高频尺度数据标识。
9.一种电子设备，其特征在于，所述电子设备包括：
存储器，存储有可执行指令；
处理器，所述处理器运行所述存储器中的所述可执行指令，以实现权利要求1-8中任一项所述的动态坐标系弹性波逆时偏移方法。
10.一种计算机可读存储介质，其特征在于，该计算机可读存储介质存储有计算机程序，该计算机程序被处理器执行时实现权利要求1-8中任一项所述的动态坐标系弹性波逆时偏移方法。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110475407.XA CN115267901B (zh) | 2021-04-29 | 2021-04-29 | 动态坐标系弹性波逆时偏移方法、电子设备及介质 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110475407.XA CN115267901B (zh) | 2021-04-29 | 2021-04-29 | 动态坐标系弹性波逆时偏移方法、电子设备及介质 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN115267901A true CN115267901A (zh) | 2022-11-01 |
| CN115267901B CN115267901B (zh) | 2024-08-30 |
Family
ID=83745268
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202110475407.XA Active CN115267901B (zh) | 2021-04-29 | 2021-04-29 | 动态坐标系弹性波逆时偏移方法、电子设备及介质 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN115267901B (zh) |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN104572575A (zh) * | 2015-01-27 | 2015-04-29 | 南京航空航天大学 | 一种特大变形动态网格生成方法 |
| CN109188519A (zh) * | 2018-09-17 | 2019-01-11 | 中国石油大学(华东) | 一种极坐标下的弹性波纵横波速度反演系统及方法 |
| CN110488354A (zh) * | 2019-07-19 | 2019-11-22 | 中国石油大学(华东) | 一种q补偿的起伏地表棱柱波与一次波联合最小二乘逆时偏移成像方法 |
| US20200233112A1 (en) * | 2019-01-22 | 2020-07-23 | Saudi Arabian Oil Company | Avo imaging condition in elastic reverse time migration |
-
2021
- 2021-04-29 CN CN202110475407.XA patent/CN115267901B/zh active Active
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN104572575A (zh) * | 2015-01-27 | 2015-04-29 | 南京航空航天大学 | 一种特大变形动态网格生成方法 |
| CN109188519A (zh) * | 2018-09-17 | 2019-01-11 | 中国石油大学(华东) | 一种极坐标下的弹性波纵横波速度反演系统及方法 |
| US20200233112A1 (en) * | 2019-01-22 | 2020-07-23 | Saudi Arabian Oil Company | Avo imaging condition in elastic reverse time migration |
| CN110488354A (zh) * | 2019-07-19 | 2019-11-22 | 中国石油大学(华东) | 一种q补偿的起伏地表棱柱波与一次波联合最小二乘逆时偏移成像方法 |
Non-Patent Citations (2)
| Title |
|---|
| 李文杰;魏修成;宁俊瑞;张建伟;王红旗;: "叠前弹性波逆时深度偏移及波场分离技术探讨", 物探化探计算技术, no. 06, 15 November 2008 (2008-11-15) * |
| 陈可洋;陈树民;李来林;吴清岭;范兴才;刘振宽;: "弹性波联合叠前逆时偏移数值试验", 石油物探, no. 01, 25 January 2014 (2014-01-25) * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN115267901B (zh) | 2024-08-30 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN103713315B (zh) | 一种地震各向异性参数全波形反演方法及装置 | |
| CN103238158B (zh) | 利用互相关目标函数进行的海洋拖缆数据同时源反演 | |
| CN101688926B (zh) | 基于拉普拉斯域波形反演的地下结构成像装置和分析方法 | |
| KR101948509B1 (ko) | 지구 물리학적 데이터의 반복 반전의 아티팩트 감소 | |
| WO2020123084A1 (en) | Machine learning-augmented geophysical inversion | |
| CN106526674A (zh) | 一种三维全波形反演能量加权梯度预处理方法 | |
| BRPI0918020B1 (pt) | métodos de exploração sísmica | |
| WO2018102813A2 (en) | Seismic acquisition geometry full-waveform inversion | |
| CA2987521A1 (en) | Method for improved geophysical investigation | |
| CN110954950B (zh) | 地下横波速度反演方法、装置、计算设备及存储介质 | |
| Hu | An improved immersed boundary finite-difference method for seismic wave propagation modeling with arbitrary surface topography | |
| CN110879412A (zh) | 地下横波速度反演方法、装置、计算设备及存储介质 | |
| Shragge et al. | Tensorial elastodynamics for isotropic media | |
| CN109459787A (zh) | 基于地震槽波全波形反演的煤矿井下构造成像方法及系统 | |
| WO2022153984A1 (ja) | 学習データ生成方法、モデル生成方法および学習データ生成装置 | |
| Chen | A nested grid, nonhydrostatic, elastic model using a terrain-following coordinate transformation: The radiative-nesting boundary conditions | |
| Ha et al. | 3D Laplace-domain waveform inversion using a low-frequency time-domain modeling algorithm | |
| CN108508481B (zh) | 一种纵波转换波地震数据时间匹配的方法、装置及系统 | |
| US11119233B2 (en) | Method for estimating elastic parameters of subsoil | |
| CN111914609B (zh) | 井震联合叠前地质统计学弹性参数反演方法及装置 | |
| CN105403919B (zh) | 一种逆时偏移成像方法及装置 | |
| CN115267901B (zh) | 动态坐标系弹性波逆时偏移方法、电子设备及介质 | |
| CN111175822B (zh) | 改进直接包络反演与扰动分解的强散射介质反演方法 | |
| US10338253B2 (en) | Method of suppressing spectral artefacts of wavefield decomposition caused by imperfect extrapolation | |
| Oosterlee et al. | Shifted-Laplacian preconditioners for heterogeneous Helmholtz problems |
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 |