CN111538084B - Ovt域数据转换成方位角度域成像道集的方法及系统 - Google Patents

Ovt域数据转换成方位角度域成像道集的方法及系统 Download PDF

Info

Publication number
CN111538084B
CN111538084B CN202010529547.6A CN202010529547A CN111538084B CN 111538084 B CN111538084 B CN 111538084B CN 202010529547 A CN202010529547 A CN 202010529547A CN 111538084 B CN111538084 B CN 111538084B
Authority
CN
China
Prior art keywords
domain
imaging
ovt
space
velocity
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
CN202010529547.6A
Other languages
English (en)
Other versions
CN111538084A (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.)
China University of Petroleum East China
Original Assignee
China University of Petroleum East China
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 China University of Petroleum East China filed Critical China University of Petroleum East China
Priority to CN202010529547.6A priority Critical patent/CN111538084B/zh
Publication of CN111538084A publication Critical patent/CN111538084A/zh
Application granted granted Critical
Publication of CN111538084B publication Critical patent/CN111538084B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

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/301—Analysis for determining seismic cross-sections or geostructures
    • 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/301—Analysis for determining seismic cross-sections or geostructures
    • G01V1/302—Analysis for determining seismic cross-sections or geostructures in 3D data 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/28—Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/30—Analysis
    • G01V1/303—Analysis for determining velocity profiles or travel times
    • 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/307—Analysis for determining seismic attributes, e.g. amplitude, instantaneous phase or frequency, reflection strength or polarity
    • 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
    • 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
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/30—Noise handling
    • G01V2210/32—Noise reduction
    • G01V2210/322—Trace stacking
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/50—Corrections or adjustments related to wave propagation
    • G01V2210/51—Migration
    • G01V2210/512—Pre-stack
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/50—Corrections or adjustments related to wave propagation
    • G01V2210/58—Media-related
    • G01V2210/586—Anisotropic media
    • 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/622—Velocity, density or impedance
    • G01V2210/6222—Velocity; travel time
    • 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/626—Physical property of subsurface with anisotropy
    • 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/63—Seismic attributes, e.g. amplitude, polarity, instant phase
    • G01V2210/632—Amplitude variation versus offset or angle of incidence [AVA, AVO, AVI]

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)
  • Image Processing (AREA)

Abstract

本发明公开了一种OVT域数据转换成方位角度域成像道集的方法及系统，所述方法，首先，计算Dix层速度，并利用Dix层速度设置初始的深度域空间的层速度，并采用优化目标函数的方式，确定优化后的深度域空间的层速度；根据优化后的深度域空间的层速度，计算最大反射角，确定反射角范围；对方位角范围和反射角范围的区域进行网格划分，获得方位角度域网格；根据优化后的深度域空间的层速度和OVT域成像道集数据中的叠前时间偏移均方根速度，计算每个方位角对应的反射角，得到方位角域成像道集。本发明采用求解目标函数的方式，实现了横向变速介质中的OVT域数据转换成方位角域成像道集以提高AVA/AVAZ反演精度。

Description

OVT域数据转换成方位角度域成像道集的方法及系统
技术领域
本发明涉及地震数据处理技术领域，特别涉及一种OVT域数据转换成方位角度域成像道集的方法及系统。
背景技术
石油地球物勘探中，人工地震是最广泛应用的方法，通过人工激发震源产生的地震波在地球介质中传播，并根据不同几何形态布置的检波器接收地震波反射信号，获取不同维度的数据体，进一步通过现代计算机软件系统对数据体进行成像处理，得到地下介质的构造形态和特征，即地震像。根据地震像的维度不同，地震勘探可以分为二维和三维勘探；由于现代计算机处理能力的提高和实际勘探目标的需要，三维勘探已经成为工业界最常用的勘探方法。三维勘探观测系统，包括在地表二维平面上配置的震源点和检波点四维坐标，加上旅行时方向，使得最终接收到的地震数据具有五个维度；这个五维坐标有不同的空间表示形式，上述震源点、检波点和旅行时坐标空间数据通过变换，可以得到中心点、偏移距和旅行时坐标下的数据，其中中心点表示炮点和检波点连线的中心位置，偏移距表示从震源点指向检波点的矢量，在极坐标系中可以使用方位角和绝对偏移距表示。振幅随偏移距变化(AVO，Amplitude versus Offset)或振幅随反射角变化(AVA，AmplitudeversusAngle)关系(由于水平层状介质假设情况下偏移距和反射角度有等价关系，因此以下描述只用AVA术语)，是重要的岩性和油气特征分析技术手段，过去AVA分析和反演并不考虑方位特性(即偏移距的矢量特征)，所以有一定的局限性。如今宽方位采集数据已经被广泛应用，考虑方位AVA分析技术是比较现实的需求，但方位AVA反演的前提需要随方位和偏移距(绝对偏移距)变化的数据体，因此五维OVT(OffsetVectorTile)数据体技术出现对方位AVO/AVA反演有着重要的意义。
五维OVT域处理技术已经成为国内外广泛应用的宽方位地震勘探技术(印兴耀等，2018；詹仕凡等，2015)，对国内“两宽一高”技术的发展有着重要的意义。上世纪末在宽方位地震数据采集观测系统设计研究中开始出现了OVT概念(Vermeer，1998；Cary,1999)，随后数位学者(Jenner et al.,2001；Williams&Jenner,2001,2002)实际应用时发现方位道集在研究方位各向异性速度分析以及地震振幅随入射角、方位角(AVA/AVAZ)变化特性方面具有重要优势。五维OVT域叠前道集数据用于AVA/AVAZ分析能够充分利用地震数据方位特性，有效提高对地下复杂介质的分析能力，比如用OVT域道集的属性分析或者AVA/AVAZ反演估计裂缝密度和方位信息等，在储层预测和油藏描述中发挥着重要的作用。基于OVT域道集的分析和反演技术使得地震解释真正进入五维空间去分析和识别储层和油藏特性，显然，精确的方位角度域成像道集生成方法是OVT域叠前AVA/AVAZ反演的基础。
从OVT域数据生成方位角域成像道集有两种途径：一是直接法，即OVT域地震数据通过深度域或者时间域偏移成像方法直接输出方位角度域成像道集，由于输入空间和输出空间差异造成并行设计较难，该方法精度高但效率较低，尤其是难以在宽方位高密度地震资料中实现广泛的应用；二是间接法，即OVT域地震数据通过时间域偏移成像方法生成OVT域成像道集，然后基于均匀或者层状介质假设把成像道集从方位偏移距空间投影到方位角度(反射角)空间，生成最终的方位角度域成像道集，该方法精度相对直接法比较低，但效率较高且有传统处理流程的支撑，如果使用保幅的偏移成像方法，成像后倾角、偏移距时差精确校正并且像点正确归位，每个成像点基于当前位置速度做层状介质假设进行偏移距到反射角转换也比较合理，获得的成像道集振幅精度能够满足生产需要。间接法的偏移距-反射角转换基于Ostrander(1984)、Todd&Backus(1985)、Walden(2006)和Resnick(1993)的理论，主要包含几个方面的信息：介质是均匀或者横向均匀(层状)介质；旅行时公式满足双曲时差公式；某个反射角的地震数据不是来自单一偏移距反射数据，而是来自一定范围偏移距具有一定带宽的反射数据叠加结果。Bale(2001)、Mukhopadhyay&Mallick(2011)为了更精确表示具有垂直对称轴的横向各向同性介质的偏移距-反射角转换关系，基于非双曲时差或者射线法得到相关的映射方法。在实际应用，中等和小偏移距道集AVA/AVAZ分析仍然是主流，所以考虑大偏移非双曲时差或者具有垂直对称抽的横向各向同性介质并不是重点问题，而横向变速介质问题的考虑会对提高AVA/AVAZ反演精度提高更有帮助。如何实现横向变速介质中的OVT域数据转换成方位角域成像道集以提高AVA/AVAZ反演精度成为一个亟待解决的技术问题。
发明内容
本发明的目的是提供一种OVT域数据转换成方位角度域成像道集的方法及系统，以实现横向变速介质中的OVT域数据转换成方位角域成像道集以提高AVA/AVAZ反演精度。
为实现上述目的，本发明提供了如下方案：
一种OVT域数据转换成方位角度域成像道集的方法，所述方法包括如下步骤：
采用时间域偏移成像方法将OVT域数据转换成OVT域成像道集；
根据OVT域成像道集中的叠前时间偏移均方根速度和OVT域成像道集中的成像射线单程旅行时，计算Dix层速度；
利用Dix层速度设置初始的深度域空间的层速度，并采用优化目标函数
Figure BDA0002534733240000031
的方式，确定优化后的深度域空间的层速度；其中，vDIX-INT为Dix层速度，vINT为深度域空间的层速度，r0'为深度域空间的成像出射点；
根据优化后的深度域空间的层速度和OVT域成像道集数据中的叠前时间偏移均方根速度，计算最大反射角，确定反射角范围(0，θmax)，其中，θmax为最大反射角；
对方位角
Figure BDA0002534733240000032
范围为(0，360°)和反射角θ范围为(0，θmax)的区域进行网格划分，获得方位角度域网格；
根据优化后的深度域空间的层速度和OVT域成像道集数据中的叠前时间偏移均方根速度，计算方位角度域网格中每个方位角对应的反射角，得到面元化结果，作为方位角域成像道集。
可选的，所述根据OVT域成像道集中的叠前时间偏移均方根速度和OVT域成像道集中的成像射线单程旅行时，计算Dix层速度，具体包括：
根据OVT域成像道集中的叠前时间偏移均方根速度和OVT域成像道集中的成像射线单程旅行时，利用公式
Figure BDA0002534733240000041
计算Dix层速度vDIX-INT(t)；其中，t表示时间，vRMS(t,r0)为OVT域成像道集中的叠前时间偏移均方根速度，r0为OVT域成像道集中的的成像射线出射点。
可选的，所述利用Dix层速度设置初始的深度域空间的层速度，并采用优化目标函数
Figure BDA0002534733240000042
的方式，确定优化后的深度域空间的层速度，具体包括：
利用公式vINT＝(z,r)设置初始的深度域空间的层速度vINT，其中，
Figure BDA0002534733240000043
r＝r0，r0为OVT域成像道集中的的成像射线出射点，t表示时间，z表示(x，y，z)空间的z轴的取值，r表示(x，y，z)空间的(x，y)平面的位置值；
利用深度域空间的层速度vINT，采用求解程函方程
Figure BDA0002534733240000044
的方式确定深度域空间的成像射线旅行时τ2；
根据深度域空间的成像射线旅行时τ2，利用公式
Figure BDA0002534733240000045
计算深度域空间的成像出射点r0'；
根据深度域空间的成像出射点计算目标函数
Figure BDA0002534733240000046
的目标函数值；
判断所述目标函数值是否小于目标函数阈值，得到判断结果；
若所述判断结果表示否，则利用共轭梯度法更新深度域空间的层速度，返回步骤“利用深度域空间的层速度vINT，采用求解程函方程
Figure BDA0002534733240000051
的方式确定深度域空间的成像射线旅行时τ2”；
若所述判断结果表示是，则输出更新后的深度域空间的层速度，作为优化后的深度域空间的层速度。
可选的，所述根据优化后的深度域空间的层速度和OVT域成像道集数据中的叠前时间偏移均方根速度，计算方位角度域网格中每个方位角对应的反射角，得到面元化结果，作为方位角域成像道集，具体包括：
根据优化后的深度域空间的层速度和OVT域成像道集数据中的叠前时间偏移均方根速度，利用公式
Figure BDA0002534733240000052
计算每个方位角对应的反射角，得到方位角域成像道集。
可选的，所述根据优化后的深度域空间的层速度和OVT域成像道集数据中的叠前时间偏移均方根速度，计算方位角度域网格中每个方位角对应的反射角，得到面元化结果，作为方位角域成像道集，之后还包括：
根据地震道与面元化结果的面元中心的距离，对所述面元化结果进行加权插值和均化处理。
一种OVT域数据转换成方位角度域成像道集的系统，所述系统：
时间域偏移成像转换模块，用于采用时间域偏移成像方法将OVT域数据转换成OVT域成像道集；
Dix层速度计算模块，用于根据OVT域成像道集中的叠前时间偏移均方根速度和OVT域成像道集中的成像射线单程旅行时，计算Dix层速度；
深度域空间的层速度优化模块，用于利用Dix层速度设置初始的深度域空间的层速度，并采用优化目标函数
Figure BDA0002534733240000053
的方式，确定优化后的深度域空间的层速度；其中，vDIX-INT为Dix层速度，vINT为深度域空间的层速度，r0'为深度域空间的成像出射点；
最大反射角计算模块，用于根据优化后的深度域空间的层速度和OVT域成像道集数据中的叠前时间偏移均方根速度，计算最大反射角，确定反射角范围(0，θmax)，其中，θmax为最大反射角；
网格划分模块，用于对方位角
Figure BDA0002534733240000061
范围为(0，360°)和反射角θ范围为(0，θmax)的区域进行网格划分，获得方位角度域网格；
面元化模块，用于根据优化后的深度域空间的层速度和OVT域成像道集数据中的叠前时间偏移均方根速度，计算方位角度域网格中每个方位角对应的反射角，得到面元化结果，作为方位角域成像道集。
可选的，所述Dix层速度计算模块，具体包括：
Dix层速度计算子模块，用于根据OVT域成像道集中的叠前时间偏移均方根速度和OVT域成像道集中的成像射线单程旅行时，利用公式
Figure BDA0002534733240000062
计算Dix层速度vDIX-INT(t)；其中，t表示时间，vRMS(t,r0)为OVT域成像道集中的叠前时间偏移均方根速度，r0为OVT域成像道集中的的成像射线出射点。
可选的，所述深度域空间的层速度优化模块，具体包括：
初始化子模块，用于利用公式vINT＝(z,r)设置初始的深度域空间的层速度vINT，其中，
Figure BDA0002534733240000063
r＝r0，r0为OVT域成像道集中的的成像射线出射点，t表示时间，z表示(x，y，z)空间的z轴的取值，r表示(x，y，z)空间的(x，y)平面的位置值；
深度域空间的成像射线旅行时求解子模块，用于利用深度域空间的层速度vINT，采用求解程函方程
Figure BDA0002534733240000064
的方式确定深度域空间的成像射线旅行时τ2；
深度域空间的成像出射点计算子模块，用于根据深度域空间的成像射线旅行时τ2，利用公式
Figure BDA0002534733240000065
计算深度域空间的成像出射点r0'；
目标函数值计算子模块，用于根据深度域空间的成像出射点计算目标函数
Figure BDA0002534733240000071
的目标函数值；
判断子模块，用于判断所述目标函数值是否小于目标函数阈值，得到判断结果；
深度域空间的层速度更新子模块，用于若所述判断结果表示否，则利用共轭梯度法更新深度域空间的层速度，返回步骤“利用深度域空间的层速度vINT，采用求解程函方程
Figure BDA0002534733240000072
的方式确定深度域空间的成像射线旅行时τ2”；
优化后的深度域空间的层速度输出子模块，用于若所述判断结果表示是，则输出更新后的深度域空间的层速度，作为优化后的深度域空间的层速度。
可选的，所述面元化模块，具体包括：
面元化子模块，用于根据优化后的深度域空间的层速度和OVT域成像道集数据中的叠前时间偏移均方根速度，利用公式
Figure BDA0002534733240000073
计算每个方位角对应的反射角，得到方位角域成像道集。
可选的，所述系统还包括：
加权插值和均化处理模块，用于根据地震道与面元化结果的面元中心的距离，对所述面元化结果进行加权插值和均化处理。
根据本发明提供的具体实施例，本发明公开了以下技术效果：
本发明公开了一种OVT域数据转换成方位角度域成像道集的方法及系统，所述方法包括如下步骤：采用时间域偏移成像方法将OVT域数据转换成OVT域成像道集；根据OVT域成像道集中的叠前时间偏移均方根速度和OVT域成像道集中的成像射线单程旅行时，计算Dix层速度；利用Dix层速度设置初始的深度域空间的层速度，并采用优化目标函数的方式，确定优化后的深度域空间的层速度；根据优化后的深度域空间的层速度和OVT域成像道集数据中的叠前时间偏移均方根速度，计算最大反射角；对方位角范围为
Figure BDA0002534733240000081
和反射角范围为(0，θmax)的区域进行网格划分，获得方位角度域网格；根据优化后的深度域空间的层速度和OVT域成像道集数据中的叠前时间偏移均方根速度，计算方位角度域网格中每个方位角对应的反射角，得到面元化结果，作为方位角域成像道集。本发明采用求解目标函数的方式，实现了横向变速介质中的OVT域数据转换成方位角域成像道集以提高AVA/AVAZ反演精度。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动性的前提下，还可以根据这些附图获得其他的附图。
图1为本发明提供的层状介质射线示意图；
图2为本发明提供的横向变速介质的成像射线示意图；
图3为本发明提供的一种OVT域数据转换成方位角度域成像道集的方法的流程图；
图4为本发明提供的用于AVA/AVAZ反演的方位和反射角面元示意图；
图5为本发明提供的用于AVA/AVAZ反演的方位和反射角面元中道集分布；
图6为本发明提供的成像道集在方位偏移距空间的分布情况示意图；
图7为本发明具体的实施例提供的成像道集转换示意图，图7(a)为OVT域叠前时间偏移成像道集，图7(b)为层状介质中偏移距到反射角的转换关系转换得到的方位角度域成像道集示意图，图7(c)为采用本发明的方法得到的方位角度域成像道集示意图；
图8为本发明具体的实施例提供的单一方位角度的方位角度域成像道集的剖面图，其中，图8(a)为方位角为30°反射角为4.5°的剖面图，图8(b)为方位角为90°反射角为4.5°的剖面图，图8(c)为方位角为30°反射角为13.5°的剖面，图8(d)为方位角为90°反射角为13.5°的剖面图。
具体实施方式
本发明的目的是提供一种OVT域数据转换成方位角度域成像道集的方法及系统，实现横向变速介质中的OVT域数据转换成方位角域成像道集以提高AVA/AVAZ反演精度。
为使本发明的上述目的、特征和优点能够更加明显易懂，下面结合附图和具体实施方式对发明作进一步详细的说明。
本发明的目的是实现横向变速介质中从OVT域叠前时间偏移成像道集到方位角度域成像道集的转换。首先，研究横向变速介质中精确的层速度和叠前时间偏移速度之间的复杂关系，引入了反射角和偏移距非线性转换公式；然后，实现了从OVT域叠前时间偏移成像道集到方位角度域成像道集的转换；最后，本发明所得到的精确的方位角度域成像道集生成方法能为叠前AVA/AVAZ反演提供有效的输入数据。
本发明的研究内容包括：1.从传统的层状介质中偏移距和反射角的转换关系，研究横向变速介质中的偏移距和反射角之间的非线性转换关系；2.使用该非线性转换关系，确定已知均方根速度和OVT域叠前时间偏移成像道集的反演问题的目标函数；3.根据输入的均方根速度使用Dix公式求解初始层速度，并求解非线性方程组得到最终的层速度，把OVT域叠前时间偏移道集映射到方位角度域成像道集；4.本发明在方位角度域成像道集映射过程中采用了面元化算法。
1、层状介质中偏移距到反射角的转换关系：
偏移距到反射角的转换关系是传统AVA/AVAZ反演算法的基本步骤，简单的推导方式可以从均匀介质或者层状介质(如图1)入手，基于斯奈尔定律的几何关系逐步建立反射角的求解公式。用射线参数表示反射角的公式如下：
Figure BDA0002534733240000101
其中，θ是入射角或者反射角，t是双程旅行时，r是三维观测系统中沿某炮点到检波点方位的水平距离，
Figure BDA0002534733240000102
为射线参数,vINT表示层速度。
本发明在已知均方根速度vRMS作为正常时差校正速度，τ是零偏移距自激自收的双程旅行时时间，根据下列双曲旅行时方程：
Figure BDA0002534733240000103
然后得到
Figure BDA0002534733240000104
式(1)是从偏移距到反射角的转换最常用的基本公式。本发明考虑大偏移距或者各向异性对应的反射角，使用非双曲时差的射线方程组通过求解得到双曲旅行时方程(2)的比较复杂的公式(3)；鉴于实际生产中大多数AVA/AVAZ使用仍集中在分析中小偏移对应的振幅随偏移变化特征，下面重点考虑横向速度变化介质中偏移距-角度转换问题。
2、横向变速介质中偏移距到反射角的转换关系：
间接法偏移距-反射角的转换涉及时深转换的问题，成像射线是建立时深转换的关键；成像射线是垂直于地表的地震射线。如图2所示，成像射线表达了叠前时间偏移成像道集空间(X，Y，T)和原模型空间(X，Y，Z)之间的对应关系。在层状横向均匀介质中，成像射线是垂直于地表的直线，因此建立在成像射线基础上的时深转换及偏移距-角度转换关系相对比较简单的；而在横向变速介质中，成像射线是一个弯曲射线，如图2所示，其层速度vINT与叠前时间偏移所采用的均方根速度vRMS之间是一个非常复杂的非线性关系。Dix公式是常用的均方根速度和层速度转换的公式，广义Dix公式所描述的层速度和叠前时间偏移的均方根速度之间的转换关系如下：
Figure BDA0002534733240000111
其中，vDIX-INT是Dix公式计算的层速度，rC表示成像射线地表位置(x0，y0),
Figure BDA0002534733240000112
表示几何扩散量。方程(4)说明，Dix层速度在横向变速介质中并不能完全表示真实的层速度，与真实的层速度之间的关系需要考虑几何扩散的影响。该广义Dix公式与一般的Dix区别是对自激自收双程旅行时的导数变成了偏导数，Li&Fome根据方程(4)、几何扩散表达式和程函方程表达了成像射线的非线性偏微分方程的边值问题。该非线程方程组是一个椭圆方程的柯西问题，是病态的，本发明直接从方程(4)出发，把上述问题变成一个优化问题，代价函数如下：
Figure BDA0002534733240000113
观测到的成像射线体现在叠前时间偏移成像道集的规则坐标空间中，而Dix层速度可以通过已知的叠前时间偏移均方根速度得到，然后该目标函数通过未知的真实层速度vINT预测成像射线与观测到的成像射线进行对比，越接近，估计的vINT越正确。也就是说，如果预测的rC不能满足方程(5)，前一次估计的vINT需要通过优化算法进行合理更新，然后进一步降低目标函数方程(5)的误差，最终通过多次迭代求解得到合理的层速度vINT值。
如前叙述，保幅的叠前时间偏移距已经校正了正常时差和倾角时差，与深度偏移相比，只是没有完全考虑介质的横向变速问题，但在仅考虑中小偏移距对应的反射角转换问题，双曲旅行时公式(2)近似也是足够的。观察公式(3)，在求解反射角的过程中，均方根速度vRMS、偏移距r和旅行时t是已知的，那么vINT和sinθ是线性映射关系。变换方程(3)得到
Figure BDA0002534733240000114
公式(6)代入公式(5)，得到：
Figure BDA0002534733240000121
上式表示固定偏移距下反射角求取的优化问题代价函数，从公式中可知，波前几何扩散可以近似表示为
Figure BDA0002534733240000122
公式(8)几何扩散表达从形式上仍然有一定的物理意义：tsinθvRMS是水平传播距离，r是偏移距，两者的对比关系加上纵横向速度比(层速度和均方根速度)的校正来近似表示几何扩散量。该近似进一步弱化了原公式(5)的非线性特征，有利于公式(7)代价函数表达的优化问题求解。由于反射角与偏移距相关，与通过公式(5)计算vINT不同的是：如果直接使用公式(7)计算反射角时，首先需要固定偏移距大小，在vINT对应的空间(z，r)计算每个网格点的反射角大小，根据成像射线再映射到(τ，r0)空间，且需要根据公式(6)生成层速度场，因此方程(7)表达的反射角求解方程直接求解比较复杂，仍然可以通过公式(5)首先求取vINT(z，r)，然后根据成像射线得到每个样点(τ，r0)所需的层速度vINT(τ，r0)，再使用公式(3)根据不同的偏移距得到对应的反射角。
如图3所示，本发明提供的一种OVT域数据转换成方位角度域成像道集的方法包括如下步骤：
步骤301，采用时间域偏移成像方法将OVT域数据转换成OVT域成像道集；
步骤302，根据OVT域成像道集中的叠前时间偏移均方根速度和OVT域成像道集中的成像射线单程旅行时，计算Dix层速度。
步骤302具体包括：根据OVT域成像道集中的叠前时间偏移均方根速度和OVT域成像道集中的成像射线单程旅行时，利用公式
Figure BDA0002534733240000123
计算Dix层速度vDIX-INT(t)；其中，t表示时间，vRMS(t,r0)为OVT域成像道集中的叠前时间偏移均方根速度，r0为OVT域成像道集中的的成像射线出射点。
步骤303，利用Dix层速度设置初始的深度域空间的层速度，并采用优化目标函数
Figure BDA0002534733240000131
的方式，确定优化后的深度域空间的层速度；其中，vDIX-INT为Dix层速度，vINT为深度域空间的层速度，r0'为深度域空间的成像出射点。
横向变速介质中的层速度求取过程，使用了快速行进方法(FMM,Fast MarchingMethod)法求解程函方程并用共轭梯度法求解优化问题。
步骤303具体包括：利用公式vINT＝(z,r)设置初始的深度域空间的层速度vINT，其中，
Figure BDA0002534733240000132
r＝r0，r0为OVT域成像道集中的的成像射线出射点，t表示时间，z表示(x，y，z)空间的z轴的取值，r表示(x，y，z)空间的(x，y)平面的位置值；利用深度域空间的层速度vINT，采用求解程函方程
Figure BDA0002534733240000133
的方式确定深度域空间的成像射线旅行时τ2；根据深度域空间的成像射线旅行时τ2，利用公式
Figure BDA0002534733240000134
计算深度域空间的成像出射点r0'；根据深度域空间的成像出射点计算目标函数
Figure BDA0002534733240000135
的目标函数值；判断所述目标函数值是否小于目标函数阈值，得到判断结果；若所述判断结果表示否，则利用共轭梯度法更新深度域空间的层速度，返回步骤“利用深度域空间的层速度vINT，采用求解程函方程
Figure BDA0002534733240000136
的方式确定深度域空间的成像射线旅行时τ2”；若所述判断结果表示是，则输出更新后的深度域空间的层速度，作为优化后的深度域空间的层速度。
步骤304，根据优化后的深度域空间的层速度和OVT域成像道集数据中的叠前时间偏移均方根速度，计算最大反射角，确定反射角范围(0，θmax)，其中，θmax为最大反射角。
具体的，根据公式
Figure BDA0002534733240000141
和需要分析的目标层τ，已知叠前时间偏移道集最大偏移距|r|、均方根速度和层速度，求取最大反射角θmax。
步骤305，对方位角
Figure BDA0002534733240000142
范围为(0，360°)和反射角θ范围为(0，θmax)的区域进行网格划分，获得方位角度域网格。
根据范围(0,θmax)和道集方位范围
Figure BDA0002534733240000143
分布情况面元化AVA/AVA所需的方位角度域网格
Figure BDA0002534733240000144
一般方位角间隔为60°，反射角间隔为10～20°；如果最大反射角为60°，间隔为20°，那么就是通常所用的6×3面元化形式。
步骤306，根据优化后的深度域空间的层速度和OVT域成像道集数据中的叠前时间偏移均方根速度，计算方位角度域网格中每个方位角对应的反射角，得到面元化结果，作为方位角域成像道集。
方位角度道集面元化算法：
OVT域叠前时间偏移成像道集进行一些去噪等处理校正后，具有一些较好的特性，它们携带方位和偏移距信息，单个成像点道集内偏移距分布均匀、整体的振幅能量分布比较均衡，成像点道集之间的方位和偏移分布信息基本相同；但是OVT域成像点道集数据量较大，信噪比不高。面向AVA/AVAZ反演分析的OVT域成像道集需要进一步处理：OVT域成像道集转换到方位角度域成像道集，此时需要保持原有的振幅分布特性，也需要考虑方位合理分割保留振幅方位上差异性减少反演的信息冗余提高收敛性，同时因同样偏移距浅层到深层反射角范围减小，所以根据目标层位保持深度方向反射角分割的一致性；在方位和反射角合理分割采样后，需要对每个面元的内相关原始道集进行叠加或者加权插值。
以方位角度域坐标系中6×3面元分割方式(如图4)为例对某一成像点道集某一目标层方位和反射角进行道集分布检测(如图5)。
步骤106具体包括：根据优化后的深度域空间的层速度和OVT域成像道集数据中的叠前时间偏移均方根速度，利用公式
Figure BDA0002534733240000151
计算每个方位角对应的反射角，得到方位角域成像道集。
步骤106之后还包括：根据地震道与面元化结果的面元中心的距离，对所述面元化结果进行加权插值和均化处理。
在图5中，虽然是宽方位角数据，但纵横向测线面元内道集分布并不均匀，需要对反射角的范围进一步缩小或者对道集进行外插，用以保证AVA/AVAZ输入数据分布的均衡性。在面元内需要叠加或者插值得到面元中心点的振幅值，该问题涉及到散乱点插值问题，这里选择使用反距离插值算法。
反距离插值算法是利用已知数据散点与插值点距离不同进行加权插值的方法，权值稀疏与距离成反比，重点是权值求取问题。反距离插值函数和权值表示如下：
Figure BDA0002534733240000152
其中，
Figure BDA0002534733240000153
是方位角度域空间
Figure BDA0002534733240000154
的振幅值，
Figure BDA0002534733240000155
和An分别是第n个散点的权值和振幅值，
Figure BDA0002534733240000156
是第n个散点与
Figure BDA0002534733240000157
点距离dn的倒数；根据近似AVA公式，振幅与sin2θ成线性关系，令ρ＝sin2θ，在极坐标
Figure BDA0002534733240000158
空间内求取距离dn，公式如下：
Figure BDA0002534733240000159
公式(10)的度量距离可以写成平方形式，即
Figure BDA00025347332400001510
这里使用ρ∈[0，1]，因为公式(9)反距离插值仍然是线性插值，该极坐标空间建立的振幅与坐标的关系模型满足线性形式，如果进一步考虑方位角与振幅的AVAZ近似关系，公式(11)中包含
Figure BDA00025347332400001511
和
Figure BDA00025347332400001512
同样能够体现振幅与方位之间的这种关系。
本发明还提供一种OVT域数据转换成方位角度域成像道集的系统，所述系统：
时间域偏移成像转换模块，用于采用时间域偏移成像方法将OVT域数据转换成OVT域成像道集。
Dix层速度计算模块，用于根据OVT域成像道集中的叠前时间偏移均方根速度和OVT域成像道集中的成像射线单程旅行时，计算Dix层速度。
所述Dix层速度计算模块，具体包括：Dix层速度计算子模块，用于根据OVT域成像道集中的叠前时间偏移均方根速度和OVT域成像道集中的成像射线单程旅行时，利用公式
Figure BDA0002534733240000161
计算Dix层速度vDIX-INT(t)；其中，t表示时间，vRMS(t,r0)为OVT域成像道集中的叠前时间偏移均方根速度，r0为OVT域成像道集中的的成像射线出射点。
深度域空间的层速度优化模块，用于利用Dix层速度设置初始的深度域空间的层速度，并采用优化目标函数
Figure BDA0002534733240000162
的方式，确定优化后的深度域空间的层速度；其中，vDIX-INT为Dix层速度，vINT为深度域空间的层速度，r0'为深度域空间的成像出射点。
所述深度域空间的层速度优化模块，具体包括：初始化子模块，用于利用公式vINT＝(z,r)设置初始的深度域空间的层速度vINT，其中，
Figure BDA0002534733240000163
r＝r0，z表示(x，y，z)空间的z轴的取值，r表示(x，y，z)空间的(x，y)平面的位置值；深度域空间的成像射线旅行时求解子模块，用于利用深度域空间的层速度vINT，采用求解程函方程
Figure BDA0002534733240000164
的方式确定深度域空间的成像射线旅行时τ2；深度域空间的成像出射点计算子模块，用于根据深度域空间的成像射线旅行时τ2，利用公式
Figure BDA0002534733240000165
计算深度域空间的成像出射点r0'；目标函数值计算子模块，用于根据深度域空间的成像出射点计算目标函数
Figure BDA0002534733240000171
的目标函数值；判断子模块，用于判断所述目标函数值是否小于目标函数阈值，得到判断结果；深度域空间的层速度更新子模块，用于若所述判断结果表示否，则利用共轭梯度法更新深度域空间的层速度，返回步骤“利用深度域空间的层速度vINT，采用求解程函方程
Figure BDA0002534733240000172
的方式确定深度域空间的成像射线旅行时τ2”；优化后的深度域空间的层速度输出子模块，用于若所述判断结果表示是，则输出更新后的深度域空间的层速度，作为优化后的深度域空间的层速度。
最大反射角计算模块，用于根据优化后的深度域空间的层速度和OVT域成像道集数据中的叠前时间偏移均方根速度，计算最大反射角，确定反射角范围(0，θmax)，其中，θmax为最大反射角。
网格划分模块，用于对方位角
Figure BDA0002534733240000173
范围为(0，360°)和反射角θ范围为(0，θmax)的区域进行网格划分，获得方位角度域网格。
面元化模块，用于根据优化后的深度域空间的层速度和OVT域成像道集数据中的叠前时间偏移均方根速度，计算方位角度域网格中每个方位角对应的反射角，得到面元化结果，作为方位角域成像道集。
所述面元化模块，具体包括：面元化子模块，用于根据优化后的深度域空间的层速度和OVT域成像道集数据中的叠前时间偏移均方根速度，利用公式
Figure BDA0002534733240000174
计算每个方位角对应的反射角，得到方位角域成像道集。
所述系统还包括：加权插值和均化处理模块，用于根据地震道与面元化结果的面元中心的距离，对所述面元化结果进行加权插值和均化处理。
为了说明本发明的技术效果，本发明提供了如下具体的实施例。
某探区内50平方公里目标成像面元：横测线Inline号范围1751-2447，采样间隔为12.5m，联络测线Crossline号范围2701-3623，采样间隔为6.25m，纵向时间采样样点数2501，采样间隔为2ms。目标数据量巨大，如果采用直接法生成方位角度域成像道集用于AVA/AVAZ，该较小的勘探面积也需要较长的地震资料处理时间而采用间接法生成方位角度域成像道集则比较实用。叠前时间偏移宽方位数据处理流程提供了每个中心点有320道的OVT域宽方位成像道集数据，本发明采用上述方法来生成AVA/AVAZ所需要的不同方位和不同反射角的地震数据。
如图6所示，OVT成像道集在方位-偏移距域的分布情况，首先考虑图3中方位的面元化方案，即360°按60°间隔分割形成的面元中，地震道分布极其不均匀，需要在面元化时利用反距离插值方法和面元均化方法得到合理的角道集数据；如图6中的偏移距方向所示，不同方位之间在偏移距方向道集分布不一致，导致不同方位之间的反射角求取后也是不均匀分布的，即使偏移距均匀分布，由于纵向速度的变化，反射角在同一方位不同时间方向的反射角分布也是不均匀。鉴于此，首先运用本文提出的偏移距转换反射角的方法，对每个成像道集的每一方位每一深度切片做反射角的转换，然后利用每个面元中分布的振幅散点数据进行插值(等同于加权叠加)得到最终的方位角度域成像道集。
如图7所示，其中，纵坐标为time(时间)，纵坐标为方位角和反射角对应的trace(数据分布)。图7a显示原始OVT域叠前时间偏移成像道集，通过传统的偏移距-角度转换方式，如公式(3)，可以得到方位角度域成像道集(如图7b)，选取了6个方位即30°、90°、150°、210°、270°和330°，每个方位的反射角范围为0-72°，采样间隔为6°；输入的OVT域成像道集仅320道，分布极不均匀，某些方位角度点(如图6)缺少成像道集，需要通过面元均化借道方法来弥补方位角度域内缺道的现象；图7c是用本发明的方法，即基于真实层速度与叠前偏移所用均方根速度之间的非线性关系，把偏移距映射到角度域空间得到的成像道集，该方法结果整体表现比传统的方法能量均衡和清晰，尤其浅层也得到了较好的处理，角道集分布情况更加合理。
实际AVA/AVAZ反演由于目标层有一定深度，足够大的偏移距下只有浅层反射角才能达到0～72°的范围，且为了反演的稳定性，一般会采用较大的反射角采样间隔，从而有效体现反射振幅随角度变化的差异性，减少反演的冗余信息。在图7b和7c中，单一方位的角道集数据由于原始输入道集稀疏仍然有不连续现象，在更大反射角采样间隔情况下，能够更加符合AVA/AVAZ输入数据的要求，由于反射角大采样间隔仅有稀疏几道(如6×3方案，每个方位就3道)。图8显示方位角度域成像道集沿方位角度向量点的切片，即单一方位角度向量点的剖面，用以说明本发明的方法对该实际资料三维工区五维OVT地震成像道集数据的处理情况，这里显示了四个方位角度向量的情况，方位角度向量坐标分别是(30°,4.5°)、(90°,4.5°)、(30°,13.5°)和(90°,13.5°)；每个剖面是一个三维数据体，上面是横向切片，下面是两个纵向方向的切片。从图8a与图8b或图8c与图8d分别比较，反射角相同，方位不同时，振幅的差异；从图8a与图8c或图8b与图8d分别比较，方位相同，反射角不同时，振幅的差异更为明显，总之，该数据能够用于AVA/AVAZ反演分析，相比于传统的分方位AVA反演会更稳定。
本说明书中各个实施例采用递进的方式描述，每个实施例重点说明的都是与其他实施例的不同之处，各个实施例之间相同相似部分互相参见即可。
本文中应用了具体个例对发明的原理及实施方式进行了阐述，以上实施例的说明只是用于帮助理解本发明的方法及其核心思想，所描述的实施例仅仅是本发明的一部分实施例，而不是全部的实施例，基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

Claims (10)

1.一种OVT域数据转换成方位角度域成像道集的方法，其特征在于，所述方法包括如下步骤：
采用时间域偏移成像方法将OVT域数据转换成OVT域成像道集；
根据OVT域成像道集中的叠前时间偏移均方根速度和OVT域成像道集中的成像射线单程旅行时，计算Dix层速度；
利用Dix层速度设置初始的深度域空间的层速度，并采用优化目标函数
Figure FDA0003382981370000011
的方式，确定优化后的深度域空间的层速度；其中，vDIX-INT为Dix层速度，vINT为深度域空间的层速度，r0'为深度域空间的成像出射点；
根据优化后的深度域空间的层速度和OVT域成像道集中的叠前时间偏移均方根速度，计算最大反射角，确定反射角范围(0，θmax)；其中，θmax为最大反射角；
对方位角
Figure FDA0003382981370000012
范围为(0，360°)和反射角θ范围为(0，θmax)的区域进行网格划分，获得方位角度域网格；
根据优化后的深度域空间的层速度和OVT域成像道集中的叠前时间偏移均方根速度，计算方位角度域网格中每个方位角对应的反射角，得到面元化结果，作为方位角度域成像道集。
2.根据权利要求1所述的OVT域数据转换成方位角度域成像道集的方法，其特征在于，所述根据OVT域成像道集中的叠前时间偏移均方根速度和OVT域成像道集中的成像射线单程旅行时，计算Dix层速度，具体包括：
根据OVT域成像道集中的叠前时间偏移均方根速度和OVT域成像道集中的成像射线单程旅行时，利用公式
Figure FDA0003382981370000013
计算Dix层速度vDIX-INT(t)；其中，t表示时间，vRMS(t,r0)为OVT域成像道集中的叠前时间偏移均方根速度，r0为OVT域成像道集中的成像射线出射点。
3.根据权利要求2所述的OVT域数据转换成方位角度域成像道集的方法，其特征在于，所述利用Dix层速度设置初始的深度域空间的层速度，并采用优化目标函数
Figure FDA0003382981370000021
的方式，确定优化后的深度域空间的层速度，具体包括：
利用公式vINT＝(z,r)设置初始的深度域空间的层速度vINT，其中，
Figure FDA0003382981370000022
r＝r0，z表示(x，y，z)空间的z轴的取值，r表示(x，y，z)空间的(x，y)平面的位置值；
利用深度域空间的层速度vINT，采用求解程函方程
Figure FDA0003382981370000023
的方式确定深度域空间的成像射线旅行时τ2；
根据深度域空间的成像射线旅行时τ2，利用公式
Figure FDA0003382981370000024
计算深度域空间的成像出射点r′0；
根据深度域空间的成像出射点计算目标函数
Figure FDA0003382981370000025
的目标函数值；
判断所述目标函数值是否小于目标函数阈值，得到判断结果；
若所述判断结果表示否，则利用共轭梯度法更新深度域空间的层速度，返回步骤“利用深度域空间的层速度vINT，采用求解程函方程
Figure FDA0003382981370000026
的方式确定深度域空间的成像射线旅行时τ2”；
若所述判断结果表示是，则输出更新后的深度域空间的层速度，作为优化后的深度域空间的层速度。
4.根据权利要求3所述的OVT域数据转换成方位角度域成像道集的方法，其特征在于，所述根据优化后的深度域空间的层速度和OVT域成像道集中的叠前时间偏移均方根速度，计算方位角度域网格中每个方位角对应的反射角，得到面元化结果，作为方位角度域成像道集，具体包括：
根据优化后的深度域空间的层速度和OVT域成像道集中的叠前时间偏移均方根速度，利用公式
Figure FDA0003382981370000031
计算每个方位角对应的反射角，得到方位角度域成像道集；其中，θ为反射角。
5.根据权利要求1所述的OVT域数据转换成方位角度域成像道集的方法，其特征在于，所述根据优化后的深度域空间的层速度和OVT域成像道集中的叠前时间偏移均方根速度，计算方位角度域网格中每个方位角对应的反射角，得到面元化结果，作为方位角度域成像道集，之后还包括：
根据地震道与面元化结果的面元中心的距离，对所述面元化结果进行加权插值和均化处理。
6.一种OVT域数据转换成方位角度域成像道集的系统，其特征在于，所述系统：
时间域偏移成像转换模块，用于采用时间域偏移成像方法将OVT域数据转换成OVT域成像道集；
Dix层速度计算模块，用于根据OVT域成像道集中的叠前时间偏移均方根速度和OVT域成像道集中的成像射线单程旅行时，计算Dix层速度；
深度域空间的层速度优化模块，用于利用Dix层速度设置初始的深度域空间的层速度，并采用优化目标函数
Figure FDA0003382981370000032
的方式，确定优化后的深度域空间的层速度；其中，vDIX-INT为Dix层速度，vINT为深度域空间的层速度，r′0为深度域空间的成像出射点；
最大反射角计算模块，用于根据优化后的深度域空间的层速度和OVT域成像道集中的叠前时间偏移均方根速度，计算最大反射角，确定反射角范围(0，θmax)，其中，θmax为最大反射角；
网格划分模块，用于对方位角
Figure FDA0003382981370000041
范围为(0，360°)和反射角θ范围为(0，θmax)的区域进行网格划分，获得方位角度域网格；
面元化模块，用于根据优化后的深度域空间的层速度和OVT域成像道集中的叠前时间偏移均方根速度，计算方位角度域网格中每个方位角对应的反射角，得到面元化结果，作为方位角度域成像道集。
7.根据权利要求6所述的OVT域数据转换成方位角度域成像道集的系统，其特征在于，所述Dix层速度计算模块，具体包括：
Dix层速度计算子模块，用于根据OVT域成像道集中的叠前时间偏移均方根速度和OVT域成像道集中的成像射线单程旅行时，利用公式
Figure FDA0003382981370000042
计算Dix层速度vDIX-INT(t)；其中，t表示时间，vRMS(t,r0)为OVT域成像道集中的叠前时间偏移均方根速度，r0为OVT域成像道集中的成像射线出射点。
8.根据权利要求7所述的OVT域数据转换成方位角度域成像道集的系统，其特征在于，所述深度域空间的层速度优化模块，具体包括：
初始化子模块，用于利用公式vINT＝(z,r)设置初始的深度域空间的层速度vINT，其中，
Figure FDA0003382981370000043
r＝r0，r0为OVT域成像道集中的成像射线出射点，t表示时间，z表示(x，y，z)空间的z轴的取值，r表示(x，y，z)空间的(x，y)平面的位置值；
深度域空间的成像射线旅行时求解子模块，用于利用深度域空间的层速度vINT，采用求解程函方程
Figure FDA0003382981370000051
的方式确定深度域空间的成像射线旅行时τ2；
深度域空间的成像出射点计算子模块，用于根据深度域空间的成像射线旅行时τ2，利用公式
Figure FDA0003382981370000052
计算深度域空间的成像出射点r′0；
目标函数值计算子模块，用于根据深度域空间的成像出射点计算目标函数
Figure FDA0003382981370000053
的目标函数值；
判断子模块，用于判断所述目标函数值是否小于目标函数阈值，得到判断结果；
深度域空间的层速度更新子模块，用于若所述判断结果表示否，则利用共轭梯度法更新深度域空间的层速度，返回步骤“利用深度域空间的层速度vINT，采用求解程函方程
Figure FDA0003382981370000054
的方式确定深度域空间的成像射线旅行时τ2”；
优化后的深度域空间的层速度输出子模块，用于若所述判断结果表示是，则输出更新后的深度域空间的层速度，作为优化后的深度域空间的层速度。
9.根据权利要求8所述的OVT域数据转换成方位角度域成像道集的系统，其特征在于，所述面元化模块，具体包括：
面元化子模块，用于根据优化后的深度域空间的层速度和OVT域成像道集中的叠前时间偏移均方根速度，利用公式
Figure FDA0003382981370000055
计算每个方位角对应的反射角，得到方位角度域成像道集；其中，θ为反射角。
10.根据权利要求6所述的OVT域数据转换成方位角度域成像道集的系统，其特征在于，所述系统还包括：
加权插值和均化处理模块，用于根据地震道与面元化结果的面元中心的距离，对所述面元化结果进行加权插值和均化处理。
CN202010529547.6A 2020-06-11 2020-06-11 Ovt域数据转换成方位角度域成像道集的方法及系统 Active CN111538084B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010529547.6A CN111538084B (zh) 2020-06-11 2020-06-11 Ovt域数据转换成方位角度域成像道集的方法及系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010529547.6A CN111538084B (zh) 2020-06-11 2020-06-11 Ovt域数据转换成方位角度域成像道集的方法及系统

Publications (2)

Publication Number Publication Date
CN111538084A CN111538084A (zh) 2020-08-14
CN111538084B true CN111538084B (zh) 2022-02-11

Family

ID=71980944

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010529547.6A Active CN111538084B (zh) 2020-06-11 2020-06-11 Ovt域数据转换成方位角度域成像道集的方法及系统

Country Status (1)

Country Link
CN (1) CN111538084B (zh)

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103543466A (zh) * 2012-07-17 2014-01-29 中国石油化工股份有限公司 一种时间域地震层速度反演方法
WO2015053811A1 (en) * 2013-10-11 2015-04-16 Chevron U.S.A. Inc. System and method for regularizing seismic data
CN105093299A (zh) * 2015-07-24 2015-11-25 中国石油天然气集团公司 一种基于炮检距向量片技术优化观测系统的方法及装置
CN106094029A (zh) * 2016-08-24 2016-11-09 中国石油集团川庆钻探工程有限公司地球物理勘探公司 利用偏移距矢量片地震数据预测储层裂缝的方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US10191165B2 (en) * 2015-01-13 2019-01-29 Cgg Services Sas Using an offset vector tile gather to image a subsurface

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103543466A (zh) * 2012-07-17 2014-01-29 中国石油化工股份有限公司 一种时间域地震层速度反演方法
WO2015053811A1 (en) * 2013-10-11 2015-04-16 Chevron U.S.A. Inc. System and method for regularizing seismic data
CN105093299A (zh) * 2015-07-24 2015-11-25 中国石油天然气集团公司 一种基于炮检距向量片技术优化观测系统的方法及装置
CN106094029A (zh) * 2016-08-24 2016-11-09 中国石油集团川庆钻探工程有限公司地球物理勘探公司 利用偏移距矢量片地震数据预测储层裂缝的方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
Azimuthal velocity analysis of 3D seismic for fractures:Altamont-Bluebell field;Khaled Al Dulaijan et al.;《CREWES Research Report》;20161231;第28卷;第1-15页 *
均值迭代法求取层速度的低频趋势;戴丹青等;《石油地球物理勘探》;20140831;第49卷(第4期);第672-677页 *

Also Published As

Publication number Publication date
CN111538084A (zh) 2020-08-14

Similar Documents

Publication Publication Date Title
US8619498B2 (en) Device and method for calculating 3D angle gathers from reverse time migration
US6859734B2 (en) Method and system for limited frequency seismic imaging
Zhang et al. Angle gathers from reverse time migration
US20100118652A1 (en) Determination of depth moveout and of residual radii of curvature in the common angle domain
CN110456417B (zh) 一种地震数据多次波压制方法
CN109765616B (zh) 一种保幅波场延拓校正方法及系统
US20120010820A1 (en) Fresnel Zone Fat Ray Tomography
CN102901984B (zh) 真地表地震数据倾角道集构建方法
CN102866426A (zh) 一种利用avo大角度道集分析岩体油气信息的方法
Ivanov et al. Traveltime parameters in tilted orthorhombic medium
Zhang et al. Target-oriented Gaussian beam migration using a modified ray tracing scheme
CN111103620A (zh) 一种岩巷超前探测三维偏移成像方法
Gong et al. Combined migration velocity model-building and its application in tunnel seismic prediction
CN111538084B (zh) Ovt域数据转换成方位角度域成像道集的方法及系统
Cheng et al. 3D Kirchhoff prestack time migration in average illumination-azimuth and incident-angle domain for isotropic and vertical transversely isotropic media
Yan et al. Analysis of converted-wave extended images for migration velocity analysis
CN113917533B (zh) Ti介质双联动全方位成像的系统性实现方法
CN109490964B (zh) 一种改进的高精度avo弹性参数快速反演方法
Liu et al. Reverse-time migration and amplitude correction in the angle-domain based on Poynting vector
Zhang et al. Research on joint inversion method of seismic compressional and shear waves from ocean-bottom nodes
XU et al. A fast and efficient common conversion point stacking technique for converted waves
CN111158053B (zh) 裂缝预测方法及装置
Khaksar et al. Angle-and azimuth-domain common-image gathers by reverse time migration for 3D elastic isotropic media
US12000971B2 (en) Method and system for seismic processing using virtual trace bins based on offset attributes and azimuthal attributes
Zhang et al. Polarization-based wave-equation migration velocity analysis in acoustic media

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