CN110907993A - 一种最小二乘偏移成像方法及系统 - Google Patents
一种最小二乘偏移成像方法及系统 Download PDFInfo
- Publication number
- CN110907993A CN110907993A CN201811087883.9A CN201811087883A CN110907993A CN 110907993 A CN110907993 A CN 110907993A CN 201811087883 A CN201811087883 A CN 201811087883A CN 110907993 A CN110907993 A CN 110907993A
- Authority
- CN
- China
- Prior art keywords
- offset
- velocity field
- field
- imaging
- migration
- 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.)
- Pending
Links
- 238000003384 imaging method Methods 0.000 title claims abstract description 91
- 238000013508 migration Methods 0.000 claims abstract description 39
- 230000005012 migration Effects 0.000 claims abstract description 35
- 238000000034 method Methods 0.000 claims abstract description 24
- 230000008569 process Effects 0.000 claims abstract description 11
- 238000004088 simulation Methods 0.000 claims abstract description 9
- 230000009466 transformation Effects 0.000 claims description 6
- 230000001131 transforming effect Effects 0.000 claims description 6
- 230000006870 function Effects 0.000 claims description 4
- 238000004590 computer program Methods 0.000 claims description 3
- 238000004364 calculation method Methods 0.000 abstract description 21
- 238000010586 diagram Methods 0.000 description 9
- 230000000694 effects Effects 0.000 description 3
- 238000004422 calculation algorithm Methods 0.000 description 2
- 238000005516 engineering process Methods 0.000 description 2
- 238000001914 filtration Methods 0.000 description 2
- 238000005286 illumination Methods 0.000 description 2
- 238000004321 preservation Methods 0.000 description 2
- 230000001133 acceleration Effects 0.000 description 1
- 230000003044 adaptive effect Effects 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000002939 conjugate gradient method Methods 0.000 description 1
- 230000002950 deficient Effects 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 238000004519 manufacturing process Methods 0.000 description 1
- 238000013507 mapping Methods 0.000 description 1
- 238000007781 pre-processing Methods 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 230000000007 visual effect 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/34—Displaying seismic recordings or visualisation of seismic data or attributes
- G01V1/345—Visualisation of seismic data or attributes, e.g. in 3D cubes
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/24—Recording seismic data
- G01V1/247—Digital recording of seismic data, e.g. in acquisition units or nodes
-
- 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/36—Effecting static or dynamic corrections on records, e.g. correcting spread; Correlating seismic signals; Eliminating effects of unwanted energy
- G01V1/362—Effecting static or dynamic corrections; Stacking
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)
- Radar Systems Or Details Thereof (AREA)
Abstract
本发明涉及一种最小二乘偏移成像方法及系统,属于油气物探工程领域。本发明提到的最小二乘偏移成像方法为,将深度域偏移速度场变换为变换域偏移速度场计算得到偏移成像结果以及反偏移模拟数据;计算反偏移模拟数据与炮记录数据之间的残差;计算梯度,根据梯度迭代地更新偏移成像结果,直到残差值小于设定低值,输出最终成像结果。本发明通过将深度域偏移速度场变换为变换域偏移速度场,使得速度场的高速区域在变换域中被压缩,低速区域在变换域中被拉伸,可在保证计算精度的前提下,降低波场延拓的空间点数,提高计算效率,同时模型网格大大减小,有效减小了反演迭代过程中占用的计算内存。
Description
技术领域
本发明涉及一种最小二乘偏移成像方法及系统,属于油气物探工程领域。
背景技术
针对复杂探区中深部储层保幅成像问题,无倾角限制的逆时偏移技术被认为是较好的成像方法,然而单纯的逆时偏移技术在处理深部储层的成像问题时,还存在如下困难:(1)分辨率有限,由于逆时偏移仅仅是正演算子的共轭转置,而不是它的逆,另外由于采集孔径有限、地下模型复杂、地震波带宽有限,因此逆时偏移仅仅能够提供比较模糊的构造信息,无法对复杂油气藏进行精细成像;(2)保幅性差,受相同因素的制约,常规逆时偏移也无法提供准确的振幅信息,深部成像振幅弱且振幅不均衡,深部有效能量识别困难;(3)偏移噪声大,逆时偏移的偏移噪声主要包括采集脚印和低频噪声,常用的成像后滤波算法不保幅,在去除低频噪音的同时会引入高频噪声。
最小二乘偏移(Least-Squares Migration,LSM)是在最小二乘反演的理论框架下,通过构建反偏移模拟算子和偏移算子,将常规偏移成像结果作为最小二乘反演的第一步,通过不断更新成像结果,达到对深部储层高精度成像的目的。但计算成本太高是阻碍LSM偏移应用于实际生产的最主要因素。数据域最小二乘偏移每次迭代大概需要两次偏移的计算量,如果不采取一定的加速策略,想得到比较理想的反演结果往往需要进行几十次甚至上百次迭代。
一篇公布号为CN 107870363 A的中国发明专利申请文献与一篇公布号为CN107957591 A的中国发明专利申请文献都是利用最小二乘偏移成像方法进行储层保幅成像。这两种方法都是在常规深度域实现波场延拓的,其计算耗时长、效率低下,且计算占用内存大,对计算设备依赖过高,并且在大规模推广应用方面有缺陷,不利于推广应用。
发明内容
本发明的目的在于提供一种最小二乘偏移成像方法及系统,用于解决现有方法在常规深度域实现波场延拓,计算耗时长、效率低下,且计算占用内存大的问题。
为实现上述目的,本发明提出一种最小二乘偏移成像方法及系统,一种最小二乘偏移成像方法,包括以下步骤:
(1)采集野外的观测数据,计算得到炮记录数据,预设残差阈值为设定低值;
(2)利用子波估计方法计算震源子波函数;
(3)将深度域偏移速度场变换为变换域偏移速度场计算得到偏移成像结果以及反偏移模拟数据;
(4)计算所述反偏移模拟数据与所述炮记录数据之间的残差;;
(5)计算梯度,根据梯度迭代地更新偏移成像结果,直到所述残差值小于所述设定低值,输出最终成像结果。
本方法通过将深度域偏移速度场变换为变换域偏移速度场,使得速度场的高速区域在变换域中被压缩,低速区域在变换域中被拉伸,可在保证计算精度的前提下,降低波场延拓的空间点数,提高计算效率,同时模型网格大大减小,有效减小了反演迭代过程中占用的计算内存。
进一步的,最终成像结果为变换域偏移速度场的成像结果,或者为将变换域偏移速度场的成像结果变换回深度域偏移速度场的成像结果。
将变换域偏移速度场的成像结果变换回深度域偏移速度场的成像结果,结果更加直观,对比明确,与常规深度域偏移速度场的成像结果相比,大大提高成像精度。
进一步的,变换域偏移速度场为伪时间域偏移速度场,将深度域偏移速度场变换为伪时间域偏移速度场的变换公式为:
其中,ν为偏移速度,x为横向坐标,z为垂向深度,τ为单程垂直旅行时间。
由于时间与速度相互关系明确,计算简单,所以将变换域偏移速度场设定为
伪时间域。
进一步的,设定低值为1%。
最小二乘偏移成像方法是通过残差对成像结果进行修正,所以设定残差阈值为1%保证成像的精确。
进一步的,计算伪时间域偏移速度场的偏移成像结果过程中,伪时间域反偏移过程的波场延拓公式为:
通过此公式可以较为准确的获得伪时间域偏移速度场的偏移成像结果。
进一步的,在计算梯度过程中,伪时间域偏移过程的波场延拓公式为:
通过此公式可以较为准确的计算梯度以及更新成像结果。
进一步的,将伪时间域偏移速度场的成像结果变换回深度域偏移速度场的成像结果的变换公式为:
其中,ν为偏移速度,x为横向坐标,z为垂向深度,τ为单程垂直旅行时间。
通过此公式可以较为准确的获得成像结果。
一种最小二乘偏移成像系统,包括控制器,控制器包括处理器与存储器,其特征在于,处理器运行存储在所述存储器中的计算机程序,以实现上述方法。
附图说明
图1为本发明的流程示意图;
图2为本发明中深度域偏移速度场示意图;
图3为本发明中模拟野外放炮正演计算得到的炮记录数据示意图;
图4为本发明的震源子波示意图;
图5为本发明的伪时间域偏移速度场示意图;
图6为本发明第1次迭代得到的梯度更新量示意图;
图7为本发明伪时间域偏移速度场成像效果图;
图8为本发明反变换得到的最终的深度域偏移速度场成像效果图;
图9为现有技术中常规深度域偏移成像效果图;
图10为本发明深度域的迭代误差曲线。
具体实施方式
下面结合附图对本发明的实施方式作进一步说明。
最小二乘偏移成像方法实施例:
最小二乘偏移成像方法是基于反演理论的成像方法,算法核心是根据反偏移模拟数据与观测数据的匹配程度来判定成像的准确性,并根据残差对成像结果进行修正。
本发明的主要构思为将深度域偏移速度场变换为变换域偏移速度场进行计算成像结果,待反偏移模拟数据与观测数据的残差降到设定低值,再将变换域偏移速度场的成像结果变换回深度域偏移速度场的成像结果。
在此实施例中变换域偏移速度场为伪时间域偏移速度场,但是变换域偏移速度场不仅限于伪时间域偏移速度场,作为其他实施方式,变换域偏移速度场也可采用曲线坐标域偏移速度场、自适应网格域偏移速度场、小波域偏移速度场等,但是变换域偏移速度场中伪时间域偏移速度场成像最好,其他变换域偏移速度场改善不如伪时间域偏移速度场效果明显。
具体实施步骤如图1所示:
1)输入初始偏移成像结果、深度域偏移速度场,以及野外采集的观测数据,预设控制参数—残差阈值为1%(设定低值)。在实践中,根据需要,初始偏移成像结果可以是已有的深度域偏移剖面,在此假设没有,令其为0;复杂断块模型的深度域偏移速度场如图2所示,模型网格大小为1441*917,纵横向网格间距均为5m;图3为模拟野外放炮正演计算得到的炮记录数据。
2)从步骤1)中的观测数据中选取偏移距为200m以内的近偏移距数据,利用常规子波估计方法计算得到的震源子波函数如图4所示。
由于考虑了速度加权映射,原始深度域偏移速度场的纵向网格数目为917,新的伪时间域偏移速度场的纵向网格数目变为601,减少量达到30%,有效减小了反演迭代过程中占用的计算内存,且通过对比图2和图5,还可以发现,深度域偏移速度场中深部的高速区域在伪时间域偏移速度场被压缩,浅部的低速区域在伪时间域偏移速度场被拉伸,从而使整个空间采样更加均匀,因而可提高计算效率。
4)利用震源子波函数、伪时间域偏移速度场,通过伪时间域反偏移算子计算当前成像结果下的模拟记录,伪时间域反偏移过程的波场延拓公式为:
5)计算反偏移模拟数据与野外采集数据之间的残差。
6)实现数据残差的伪时间域偏移,计算对应的更新梯度,如图6所示,伪时间域偏移过程的波场延拓公式为;
7)通过共轭梯度法迭代更新成像结果,直到满足步骤1)中设定的残差阈值条件,即残差阈值降到1%以下,终止更新成像结果。最终伪时间域偏移速度场成像结果如图7所示。
其中,迭代更新过程为:
mk+1=mk-αkPkgk
gk=L*(Lmk-pobs)
其中,αk为第k次迭代的更新步长,gk为第k次迭代的梯度,Pk为预处理算子,在此为震源照明,mk为第k次迭代的成像结果,L为反偏移模拟算子,L*为L的伴随算子。
8)将最终伪时间域偏移速度场成像结果反变换回深度域偏移速度场成像结果,如图8所示,反变换公式为:
其中,ν为偏移速度,x为横向坐标,z为垂向深度,τ为单程垂直旅行时间。
步骤8)中又变回深度域偏移速度场,是为了与常规的深度域偏移速度场得出的成像结果进行比较,因此,也可以不变换回深度域偏移速度场,直接使用伪时间域偏移速度场的成像结果即可。
图9为现有技术中常规深度域偏移成像结果,通过对比图8和图9可以看出,本发明的反变换得到的最终的深度域偏移速度场成像结果(图8)相比现有技术中常规深度域偏移成像结果(图9)有较大改善,具体改善有以下几点:
(1)不仅低频噪音得到较好的压制,且不含常规偏移中由于滤波引入的浅层高频噪音;
(2)纵向分辨率有了较大提高,尤其是在复杂的逆掩断层附近;
(3)振幅更加均衡,在剖面两侧照明补偿较好。
图10为本发明深度域中数据残差以及模型残差随迭代次数的收敛曲线,其实深度域和伪时间域的收敛曲线结果是一样的,都是基于数据驱动,不受计算空间的影响,伪时间域的优点主要在于计算效率高、内存占用小,计算用时相比常规最小二乘偏移方法节省了近25%,更加高效率,所以伪时间域在保证结果相同的情况下提高计算效率,在一定程度上说明了方法的可靠性。
最小二乘偏移成像系统实施例:
最小二乘偏移成像系统,包括控制器,控制器包括处理器与存储器,处理器运行存储在存储器中的计算机程序,以实现上述最小二乘偏移成像方法,具体最小二乘偏移成像方法已经在上述方法实施例中介绍,这里不做赘述。
以上是对本发明的较佳实施方式进行了具体的说明,但本发明并不限于以上实例,熟悉本领域的技术人员在不违背本发明精神的前提下还可以做出种种的等同变形和替换,这些等同变形和替换均包含在本申请权利要求所限定的范围内。
Claims (8)
1.一种最小二乘偏移成像方法,其特征在于,包括以下步骤:
(1)采集野外的观测数据,计算得到炮记录数据,预设残差阈值为设定低值;
(2)利用子波估计方法计算震源子波函数;
(3)将深度域偏移速度场变换为变换域偏移速度场计算得到偏移成像结果以及反偏移模拟数据;
(4)计算所述反偏移模拟数据与所述炮记录数据之间的残差;
(5)计算梯度,根据梯度迭代地更新偏移成像结果,直到所述残差值小于所述设定低值,输出最终成像结果。
2.根据权利要求1所述的最小二乘偏移成像方法,其特征在于,所述最终成像结果为变换域偏移速度场的成像结果,或者为将变换域偏移速度场的成像结果变换回深度域偏移速度场的成像结果。
4.根据权利要求1所述的最小二乘偏移成像方法,其特征在于,所述设定低值为1%。
8.一种最小二乘偏移成像系统,包括控制器,所述控制器包括处理器与存储器,其特征在于,所述处理器运行存储在所述存储器中的计算机程序,以实现权利要求1-7中任意一项所述的方法。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201811087883.9A CN110907993A (zh) | 2018-09-18 | 2018-09-18 | 一种最小二乘偏移成像方法及系统 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201811087883.9A CN110907993A (zh) | 2018-09-18 | 2018-09-18 | 一种最小二乘偏移成像方法及系统 |
Publications (1)
Publication Number | Publication Date |
---|---|
CN110907993A true CN110907993A (zh) | 2020-03-24 |
Family
ID=69813534
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201811087883.9A Pending CN110907993A (zh) | 2018-09-18 | 2018-09-18 | 一种最小二乘偏移成像方法及系统 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN110907993A (zh) |
Cited By (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN112083493A (zh) * | 2020-08-19 | 2020-12-15 | 中国石油大学(华东) | 一种三维c-τ坐标系的圆锥波编码多震源最小二乘逆时偏移成像方法 |
Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN104360381A (zh) * | 2014-10-20 | 2015-02-18 | 李闯 | 一种地震资料的偏移成像处理方法 |
WO2016001750A2 (en) * | 2014-06-30 | 2016-01-07 | Cgg Services Sa | Seismic data processing using matching filter based cost function optimization |
US20160170059A1 (en) * | 2014-12-16 | 2016-06-16 | Pgs Geophysical As | Visco-acoustic reverse-time migration using pseudo-analytical method |
CN107356972A (zh) * | 2017-06-28 | 2017-11-17 | 中国石油大学(华东) | 一种各向异性介质的成像方法 |
CN107589443A (zh) * | 2017-08-16 | 2018-01-16 | 东北石油大学 | 基于弹性波最小二乘逆时偏移成像的方法及系统 |
-
2018
- 2018-09-18 CN CN201811087883.9A patent/CN110907993A/zh active Pending
Patent Citations (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2016001750A2 (en) * | 2014-06-30 | 2016-01-07 | Cgg Services Sa | Seismic data processing using matching filter based cost function optimization |
CN104360381A (zh) * | 2014-10-20 | 2015-02-18 | 李闯 | 一种地震资料的偏移成像处理方法 |
US20160170059A1 (en) * | 2014-12-16 | 2016-06-16 | Pgs Geophysical As | Visco-acoustic reverse-time migration using pseudo-analytical method |
CN107356972A (zh) * | 2017-06-28 | 2017-11-17 | 中国石油大学(华东) | 一种各向异性介质的成像方法 |
CN107589443A (zh) * | 2017-08-16 | 2018-01-16 | 东北石油大学 | 基于弹性波最小二乘逆时偏移成像的方法及系统 |
Non-Patent Citations (4)
Title |
---|
QINGYANG LI ET AL.: "Cross-correlation least-squares reverse time migration in the pseudo-time domain", 《JOURNAL OF GEOPHYSICS AND ENGINEERING》 * |
SUN XIAO-DONG ET AL.: "Least squares reverse-time migration in the pseudodepth domain and reservoir exploration", 《APPLIED GEOPHYSICS》 * |
孙郧松: "分频编码最小二乘偏移方法研究", 《中国优秀博硕士学位论文全文数据库(硕士) 基础科学辑》 * |
李振春等: "基于平面波加速的VTI介质最小二乘逆时偏移", 《地球物理学报》 * |
Cited By (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN112083493A (zh) * | 2020-08-19 | 2020-12-15 | 中国石油大学(华东) | 一种三维c-τ坐标系的圆锥波编码多震源最小二乘逆时偏移成像方法 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Yang et al. | Deep learning seismic random noise attenuation via improved residual convolutional neural network | |
CN112083482B (zh) | 基于模型驱动深度学习的地震超分辨反演方法 | |
CN106405651B (zh) | 一种基于测井匹配的全波形反演初始速度模型构建方法 | |
CN105277985A (zh) | 一种基于图像处理的ovt域地震数据规则化方法 | |
CN107368668A (zh) | 基于双重稀疏字典学习的地震数据去噪方法 | |
CN112364291B (zh) | 一种前置滤波极值点优化集合经验模态分解方法及装置 | |
CN109633752B (zh) | 基于三维快速Radon变换的海上拖缆资料自适应鬼波压制方法 | |
CN111290019B (zh) | 一种应用于最小二乘逆时偏移的l-bfgs初始矩阵求取方法 | |
CN104749631A (zh) | 一种基于稀疏反演的偏移速度分析方法及装置 | |
CN112435162B (zh) | 一种基于复数域神经网络的太赫兹图像超分辨重建方法 | |
Xue et al. | Airborne electromagnetic data denoising based on dictionary learning | |
CN110488354B (zh) | 一种q补偿的起伏地表棱柱波与一次波联合最小二乘逆时偏移成像方法 | |
CN111768349A (zh) | 一种基于深度学习的espi图像降噪方法及系统 | |
CN110907993A (zh) | 一种最小二乘偏移成像方法及系统 | |
Li et al. | A deep learning approach for acoustic FWI with elastic data | |
CN105929447B (zh) | 考虑地震子波拉伸效应的变顶点稀疏双曲线Radon变换方法 | |
US20230095632A1 (en) | Interpretive-guided velocity modeling seismic imaging method and system, medium and device | |
CN116068644A (zh) | 一种利用生成对抗网络提升地震数据分辨率和降噪的方法 | |
CN110261912A (zh) | 一种地震数据的插值和去噪方法及系统 | |
CN109188516A (zh) | Radon域能量扫描叠加的微地震事件定位方法 | |
CN109856673B (zh) | 一种基于优势频率迭代加权的高分辨Radon变换数据分离技术 | |
CN114509814A (zh) | 一种叠前地震资料随机噪音压制方法及系统 | |
CN112698389B (zh) | 一种地震资料反演成像方法及装置 | |
CN118033732B (zh) | 一种基于空域频域融合架构的地震数据重建方法 | |
CN114817854B (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 | ||
RJ01 | Rejection of invention patent application after publication | ||
RJ01 | Rejection of invention patent application after publication |
Application publication date: 20200324 |