CN108549102B - 联合梯度结构张量和多窗分析的地层结构曲率估计方法 - Google Patents

联合梯度结构张量和多窗分析的地层结构曲率估计方法 Download PDF

Info

Publication number
CN108549102B
CN108549102B CN201810273768.4A CN201810273768A CN108549102B CN 108549102 B CN108549102 B CN 108549102B CN 201810273768 A CN201810273768 A CN 201810273768A CN 108549102 B CN108549102 B CN 108549102B
Authority
CN
China
Prior art keywords
analysis
dimensional
window
calculating
dip angle
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
CN201810273768.4A
Other languages
English (en)
Other versions
CN108549102A (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.)
Xian Jiaotong University
Original Assignee
Xian Jiaotong University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Xian Jiaotong University filed Critical Xian Jiaotong University
Priority to CN201810273768.4A priority Critical patent/CN108549102B/zh
Publication of CN108549102A publication Critical patent/CN108549102A/zh
Application granted granted Critical
Publication of CN108549102B publication Critical patent/CN108549102B/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/306—Analysis for determining physical properties of the subsurface, e.g. impedance, porosity or attenuation profiles

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

本发明公开一种联合梯度结构张量和多窗分析的地层结构曲率估计方法，包括：步骤01：采集三维叠后地震数据并变换获得每道地震信号的复地震道；步骤02：计算得到三维瞬时相位数据体和瞬时相位的梯度向量；步骤03：每个样点附近选取相邻分析窗；步骤04：在每个相邻分析窗内，计算视倾角和相似度；步骤05：在所有相邻分析窗内选择相似度最大的分析窗，将其对应的视倾角指定为分析样点的视倾角；步骤06：判断所有样点是否都计算得到了视倾角；步骤07：计算三维视倾角数据体的偏导数；步骤08：计算地层的结构曲率。本发明解决了在三维低吼地震数据地层结构曲率估计的稳定性问题及抗断层干扰问题，且实现简单易于操作。

Description

联合梯度结构张量和多窗分析的地层结构曲率估计方法
技术领域
本发明属于勘探地球物理领域，特别涉及一种地层几何曲率估计方法。
背景技术
石油及天然气是一种重要的战略资源，在我国国民经济及国家安全中占有举足轻重的地位。地震属性是对地震勘探数据的一种度量，能够在视觉上增强具有解释意义的特征，对地层特征以及储层特性有直接反映。因此，地震属性可在油气藏勘探开发中发挥重要的作用。几何属性是地震资料解释中最为重要的地震属性之一，主要用于增强和显示地震层位的几何形态，主要包括倾角\方位角、连续性和地层曲率等。倾角\方位角属性体包含有重要的地震地层学和地理学信息，不但可直接用于勘探工区的解释，还可为后续处理及解释提供基础数据(例如结构曲率计算、结构导向滤波。而地层结构曲率可通过对视倾角求取一阶导数而获取,可用于描述地质体的几何形态变化，不但对地层的弯曲、褶皱及断层等结构反映敏感，而且还可用于监控地震资料处理的质量。因此，准确稳定地从三维叠后数据中估计地层曲率具有重要的理论价值及应用价值。
目前已有多项技术可用于地层曲率估计，主要技术如下所示：
现有技术1：基于层位拾取的地层曲率估计方法
该类方法首先利用软件人工或者自动拾取目的层位，然后计算拾取到层位的二阶空间导数。基于拾取层位的二阶空间导数可以估计地层的曲率。
现有技术1的特点：
优点：简单易行。
缺点：1、拾取层位需要软件的协助。
2、只能针对目的层位计算地层曲率，不能针对这个三维地震数据体计算。
3、受人为干扰因素影响较大。
现有技术2：基于地震道互相关的地层曲率估计方法
该类方法首先计算三维地震数据各地震道与相邻道之间的互相关，进而基于互相关值估计相邻道之间的时差以获得视倾角，最后对视倾角计算空间偏导数以估计地层曲率。
现有技术2的特点：
优点：实现简单、计算量小且不受人为因素干扰。
缺点：1、由于仅采用相邻道计算互相关，受噪声影响较大。
2、振幅横向变化及跨越断层时误差较大。
现有技术3：基于倾角扫描的地层曲率估计方法
该类方法首先按照一定的采样间隔以及范围设定一系列视倾角，然后沿每个视倾角计算相似度或者相干值，选取具有最大的相似度(或相干值)的视倾角做为分析点的倾角，最后对倾角计算空间偏导数以估计地层曲率。
现有技术3的特点：
优点：较为稳定、抗噪性能好且不受人为因素干扰。
缺点：1、由于每个设定的视倾角下都要计算相似度或，导致计算量较大。
2、倾角估计的精度受限于视倾角的采样间隔。
3、分析窗跨越断层时误差较大。
现有技术4：基于梯度结构张量的地层曲率估计方法
该类方法首先在分析窗内计算所有地震数据的梯度向量，然后由梯度向量构建梯度协方差矩阵，接着对协方差矩阵进行特征分解以得到主特征向量，利用主特征向量估计地层的倾角，最后对倾角计算空间偏导数以估计地层倾角。
现有技术4的特点：
优点：估计准确、抗噪性能好且不受人为因素干扰。
缺点：振幅横向变化及跨越断层时曲率估计误差较大。
发明内容
本发明的目的在于提供一种联合梯度结构张量和多窗分析的地层结构曲率估计方法，以解决上述技术问题。
为了实现上述目的，本发明采用如下技术方案：
联合梯度结构张量和多窗分析的地层结构曲率估计方法，包括以下步骤：
步骤01：采集三维叠后地震数据并利用Hilbert变换获得每道地震信号的复地震道；
步骤02：在复地震道的基础之上计算得到三维瞬时相位数据体，同时计算瞬时相位的梯度向量；
步骤03：针对三维数据体的每个样点附近选取众多相邻分析窗；
步骤04：在每个相邻分析窗内构建相位梯度的协方差矩阵并对其进行特征分解以得到三个特征值及特征向量，利用主特征向量计算视倾角，同时利用三个特征值计算相似度；
步骤05：在所有相邻分析窗内选择相似度最大的分析窗，将其对应的视倾角指定为分析样点的视倾角；
步骤06：判断三维数据中的所有样点是否都计算得到了视倾角：若没有返回步骤03；若完成进入步骤07；
步骤07：计算三维视倾角数据体的偏导数；
步骤08：利用视倾角的偏导数计算地层的结构曲率。
进一步的，步骤01具体包括：
若用f(x,y,t)表示三维地震数据，其中x和y表示三维地震数据在两个空间方向的指标，t表示三维地震数据在时间方向的指标；固定x和y，分别对每道信号f(x,y,t)进行Hilbert变换得到h(x,y,t)为：
Figure BDA0001613170190000041
其中τ为积分临时变量；则对应的复地震道c(x,y,t)为
c(x,y,t)＝f(x,y,t)+jh(x,y,t)＝A(x,y,t)ejφ(x,y,t)，
其中j为虚数单位，A(x,y,t)为瞬时振幅，φ(x,y,t)为瞬时相位。
进一步的，步骤02具体包括：
利用复地震道c(x,y,t)的实部f(x,y,t)和虚部h(x,y,t)计算得到三维瞬时相位体φ(x,y,t)如下：
Figure BDA0001613170190000042
同时可由三维瞬时相位体φ(x,y,t)计算得到瞬时相位的梯度向量如下：
Figure BDA0001613170190000043
进一步的，步骤03具体包括：
假设当前的分析点为(x,y,t)，在其附近选取众多相邻分析窗；假设分析窗大小为Mx×My×Mt；其中Mx，My及Mt均为奇数，且Mx＝2Wx+1，My＝2Wy+1，Mt＝2Wt+1；得到MxMyMt个相邻分析窗，且这些相邻分析窗均必须包含当前分析点(x,y,t)。
进一步的，步骤04具体包括：
在步骤03得到的每个相邻分析窗内利用步骤02得到的相位梯度构建协方差矩阵，构建当前分析点为(x,y,t)的相位梯度协方差矩阵C(x,y,t)：
Figure BDA0001613170190000051
其中mx，my及mt均为求和临时指标；
对相位梯度协方差矩阵C(x,y,t)进行特征值分解得到三个特征值μ1，μ2，μ3及对应的三个特征向量v1，v2，v3，并有μ1≥μ2≥μ3：
Figure BDA0001613170190000052
其中μ1及v1为协方差矩阵C(x,y,t)的主特征值及主特征向量；利用主特征向量v1的三个元素v1x，v1y及v1t估计在此相邻分析窗内地层沿x方向的视倾角px(x,y,t)及沿y方向的视倾角py(x,y,t)：
Figure BDA0001613170190000053
Figure BDA0001613170190000054
利用三个特征值μ1，μ2及μ3计算相似度similarity(x,y,t)如下：
Figure BDA0001613170190000055
在每个相邻分析窗内均计算视倾角和相似度；若分析点(x,y,t)共有M个相邻分析窗，记第m个相邻分析窗内估计的沿x及y方向的视倾角分别为px(m,x,y,t)及py(m,x,y,t)，记第m个相邻分析窗内估计的相似度为similarity(m,x,y,t)。
进一步的，步骤05具体包括：
分析点(x,y,t)周围共有M个相邻分析窗，选出相似度最大的相邻分析窗序号为win_opt(x,y,t)：
Figure BDA0001613170190000061
则分析点(x,y,t)处的沿x方向及y方向的视倾角分别为Px(x,y,t)及Py(x,y,t)：
Px(x,y,t)＝px(win_opt(x,y,t),x,y,t)，
Py(x,y,t)＝py(win_opt(x,y,t),x,y,t)。
进一步的，步骤07具体包括：
针对沿x方向的视倾角数据体Px(x,y,t)分别计算其沿x方向偏导数Pxx(x,y,t)及y方向的偏导数Pxy(x,y,t)：
Figure BDA0001613170190000062
Figure BDA0001613170190000063
针对沿y方向的视倾角数据体Py(x,y,t)分别计算其沿x方向偏导数Pyx(x,y,t)及y方向的偏导数Pyy(x,y,t)：
Figure BDA0001613170190000064
Figure BDA0001613170190000065
Δx和Δy为三维地震数据体在x和y方向上的采样间隔。
进一步的，步骤08具体包括：
利用步骤07中计算得到的三维视倾角数据体的偏导数Pxx(x,y,t)、Pxy(x,y,t)、Pyx(x,y,t)及Pyy(x,y,t)计算各点的最大正曲率MP(x,y,t)和最大负曲率MN(x,y,t)：
Figure BDA0001613170190000071
Figure BDA0001613170190000072
其中
a(x,y,t)＝0.5×Pxx(x,y,t)，
b(x,y,t)＝0.5×Pyy(x,y,t)，
c(x,y,t)＝0.5×[Pxy(x,y,t)+Pyx(x,y,t)]。
相对于现有技术，本发明具有以下有益效果：
本发明可更加稳定地估计三维地震数据的结构曲率属性，且减少断层等不连续结构对曲率属性估计的干扰。
附图表说明
图1为本发明流程图；
图2为分析点附近的相邻分析窗在空间方向及空间-时间方向示意图；其中，图2(a)为沿二维空间方向观察相邻分析窗示意图；图2(b)为沿空间-时间方向观察相邻分析窗示意图；
图3为分析点附件的所有相邻分析窗示意图；
图4为利用商业软件计算的最大正曲率沿层切片及利用本发明方法计算的的最大正曲率沿层切片；其中，图4(a)为利用商业软件计算的最大正曲率沿层切片；图4(b)为利用本发明方法计算的的最大正曲率沿层切片；
图5为利用商业软件计算的最大负曲率沿层切片及利用本发明方法计算的的最大负曲率沿层切片；其中，图5(a)为利用商业软件计算的最大负曲率沿层切片；图5(b)为利用本发明方法计算的的最大负曲率沿层切片。
具体实施方式
下面结合附图及表格对本发明做进一步详细的说明。
本发明为针对三维叠后地震数据利用梯度结构张量和多窗分析估计地层结构曲率的方法。本发明首先对地震数据进行复地震道分析以获得瞬时相位，然后在瞬时相位数据体上构建梯度结构张量并结合多窗分析技术进行地层倾角估计，进而对倾角计算空间导数，最终利用倾角的空间导数估计三维叠后地震数据的结构曲率，为后续地震资料解释人员提供更为可靠的结构曲率属性。
请参阅图1所示，本发明一种联合梯度结构张量和多窗分析的地层结构曲率估计方法，包括以下步骤：
步骤01：采集三维叠后地震数据并利用Hilbert变换获得每道地震信号的复地震道；
步骤02：在复地震道的基础之上计算得到三维瞬时相位数据体，同时计算瞬时相位的梯度向量；
步骤03：针对三维数据体的每个样点附近选取众多相邻分析窗；
步骤04：在每个相邻分析窗内构建相位梯度的协方差矩阵并对其进行特征分解以得到三个特征值及特征向量，利用主特征向量计算视倾角，同时利用三个特征值计算相似度；
步骤05：在所有相邻分析窗内选择相似度最大的分析窗，将其对应的视倾角指定为分析样点的视倾角；
步骤06：判断三维数据中的所有样点是否都计算得到了视倾角：若没有返回步骤03；若完成进入步骤07；
步骤07：计算三维视倾角数据体的偏导数；
步骤08：利用视倾角的偏导数计算地层的结构曲率。
步骤01中采集三维叠后地震数据并利用Hilbert变换获得每道地震信号的复地震道，具体包括：
若用f(x,y,t)表示三维地震数据，其中x和y表示三维地震数据在两个空间方向的指标，t表示三维地震数据在时间方向的指标。固定x和y，分别对每道信号f(x,y,t)进行Hilbert变换得到h(x,y,t)为：
Figure BDA0001613170190000091
其中τ为积分临时变量；则对应的复地震道c(x,y,t)为
c(x,y,t)＝f(x,y,t)+jh(x,y,t)＝A(x,y,t)ejφ(x,y,t)，
其中j为虚数单位，A(x,y,t)为瞬时振幅，φ(x,y,t)为瞬时相位。
步骤02中在复地震道的基础之上计算得到三维瞬时相位数据体并计算瞬时相位的梯度向量，具体包括：
由步骤01得到的复地震道c(x,y,t)的实部f(x,y,t)和虚部h(x,y,t)计算得到三维瞬时相位体φ(x,y,t)如下：
Figure BDA0001613170190000092
其中arctan表示反正切函数。由三维瞬时相位体φ(x,y,t)计算得到瞬时相位的梯度向量如下：
Figure BDA0001613170190000093
而φ(x,y,t)的偏导数由下列方法计算：
Figure BDA0001613170190000094
Figure BDA0001613170190000095
Figure BDA0001613170190000096
其中max表示取最大值，ε表示阻尼因子，取值较小(一般取0.001-0.01之间)。另外，计算f(x,y,t)以及h(x,y,t)的偏导数通过中心差分实现：
Figure BDA0001613170190000101
Figure BDA0001613170190000102
Figure BDA0001613170190000103
Figure BDA0001613170190000104
Figure BDA0001613170190000105
Figure BDA0001613170190000106
Δx，Δy及Δt为三维地震数据体在x，y及t方向上的采样间隔。
步骤03中针对三维数据体的每个样点附件选取众多相邻分析窗，具体包括：
假设当前的分析点为(x,y,t)，在其附近选取众多相邻分析窗。假设分析窗大小为Mx×My×Mt(其中Mx，My及Mt均为奇数，且Mx＝2Wx+1，My＝2Wy+1，Mt＝2Wt+1)，可以得到MxMyMt个相邻分析窗，且这些相邻分析窗均必须包含当前分析点(x,y,t)。假设分析窗大小为3*3*3，若用黑点表示当前分析点，则沿二维空间方向观察相邻分析窗如图2(a)所示，共有1个中心点分析窗和8个非中心点分析窗；沿空间-时间方向观察相邻分析窗如图2(b)所示，共有1个中心点分析窗和2个非中心点分析窗。综上，分析窗大小为3*3*3时，共有27个相邻分析窗(如图3所示1个中心点分析窗和26个非中心点分析窗，其中第1个为中心点分析窗)。
步骤04中在每个相邻分析窗内构建相位梯度的协方差矩阵并对其进行特征分解以得到三个特征值及主特征向量，利用主特征向量计算视倾角，同时利用三个特征值计算相似度，具体包括：
在步骤03得到的每个相邻分析窗内利用步骤02得到的相位梯度构建协方差矩阵，以中心点分析窗为例(即图3中的第一个窗)，构建当前分析点为(x,y,t)的相位梯度协方差矩阵C(x,y,t)：
Figure BDA0001613170190000111
其中mx，my及mt均为求和临时指标；
对相位梯度协方差矩阵C(x,y,t)进行特征值分解得到三个特征值(μ1，μ2及μ3，并有μ1≥μ2≥μ3)及对应的三个特征向量(v1，v2及v3)：
Figure BDA0001613170190000112
其中μ1及v1为协方差矩阵C(x,y,t)的主特征值及主特征向量。利用主特征向量v1的三个元素v1x，v1y及v1t估计在此相邻分析窗内地层沿x方向的视倾角px(x,y,t)及沿y方向的视倾角py(x,y,t)：
Figure BDA0001613170190000113
Figure BDA0001613170190000114
利用三个特征值μ1，μ2及μ3计算相似度similarity(x,y,t)如下：
Figure BDA0001613170190000115
在每个相邻分析窗内均计算视倾角和相似度。若分析点(x,y,t)共有M个相邻分析窗，记第m个相邻分析窗内估计的沿x及y方向的视倾角分别为px(m,x,y,t)及py(m,x,y,t)，记第m个相邻分析窗内估计的相似度为similarity(m,x,y,t)。
步骤05：在所有相邻分析窗内选择相似度最大的分析窗，将其对应的视倾角指定为分析样点的视倾角，具体包括：
分析点(x,y,t)周围共有M个相邻分析窗，选出相似度最大的相邻分析窗序号为win_opt(x,y,t)：
Figure BDA0001613170190000121
则分析点(x,y,t)处的沿x方向及y方向的视倾角分别为Px(x,y,t)及Py(x,y,t)：
Px(x,y,t)＝px(win_opt(x,y,t),x,y,t)，
Py(x,y,t)＝py(win_opt(x,y,t),x,y,t)。
步骤06：判断三维数据中的所有样点是否都计算得到了视倾角：若没有返回步骤03；若完成进入步骤07，具体包括：
对三维地震数据体中每一点均进行步骤3至步骤6的处理，可得到沿x方向的视倾角数据体Px(x,y,t)及沿y方向的视倾角数据体Py(x,y,t)。
步骤07中计算三维视倾角数据体的偏导数，具体包括：
针对沿x方向的视倾角数据体Px(x,y,t)分别计算其沿x方向偏导数Pxx(x,y,t)及y方向的偏导数Pxy(x,y,t)：
Figure BDA0001613170190000122
Figure BDA0001613170190000123
同样地，针对沿y方向的视倾角数据体Py(x,y,t)分别计算其沿x方向偏导数Pyx(x,y,t)及y方向的偏导数Pyy(x,y,t)：
Figure BDA0001613170190000131
Figure BDA0001613170190000132
上述求偏导数过程中，Δx和Δy为三维地震数据体在x和y方向上的采样间隔。
步骤08中利用视倾角的偏导数计算地层的结构曲率，具体包括：
为简化计算，利用步骤07中计算得到的三维视倾角数据体的偏导数Pxx(x,y,t)、Pxy(x,y,t)、Pyx(x,y,t)及Pyy(x,y,t)计算各点的各种曲率(最大正曲率MP(x,y,t)，最大负曲率MN(x,y,t))：
Figure BDA0001613170190000133
Figure BDA0001613170190000134
其中：
a(x,y,t)＝0.5×Pxx(x,y,t)，
b(x,y,t)＝0.5×Pyy(x,y,t)，
c(x,y,t)＝0.5×[Pxy(x,y,t)+Pyx(x,y,t)]。
以某油田的三维叠后地震数据为例，分别采用商业软件及本发明方法计算最大正曲率和最大负曲率。商业软件计算得到的最大正曲率和最大负曲率的时间切片如图4(a)及图5(a)所示，可以看出受噪声影响较大。利用本发明方法计算得到的最大正曲率及最大负曲率的时间切片如图4(b)及图5(b)所示，可以看出相对于商业软件的结果，本发明的结果反映地质结构更清晰且抗噪性能较好。
最后需要说明的是，以上算例对本发明的目的，技术方案以及有益效果提供了进一步的验证，这仅属于本发明的具体实施算例，并不用于限定本发明的保护范围，在本发明的精神和原则之内，所做的任何修改，改进或等同替换等，均应在本发明的保护范围内。

Claims (5)

1.联合梯度结构张量和多窗分析的地层结构曲率估计方法，其特征在于，包括以下步骤：
步骤01：采集三维叠后地震数据并利用Hilbert变换获得每道地震信号的复地震道；
步骤02：在复地震道的基础之上计算得到三维瞬时相位数据体，同时计算瞬时相位的梯度向量；
步骤03：针对三维数据体的每个样点附近选取众多相邻分析窗；
步骤04：在每个相邻分析窗内构建相位梯度的协方差矩阵并对其进行特征分解以得到三个特征值及特征向量，利用主特征向量计算视倾角，同时利用三个特征值计算相似度；
步骤05：在所有相邻分析窗内选择相似度最大的分析窗，将其对应的视倾角指定为分析样点的视倾角；
步骤06：判断三维数据中的所有样点是否都计算得到了视倾角：若没有返回步骤03；若完成进入步骤07；
步骤07：计算三维视倾角数据体的偏导数；
步骤08：利用视倾角的偏导数计算地层的结构曲率；
步骤01具体包括：
若用f(x,y,t)表示三维地震数据，其中x和y表示三维地震数据在两个空间方向的指标，t表示三维地震数据在时间方向的指标；固定x和y，分别对每道信号f(x,y,t)进行Hilbert变换得到h(x,y,t)为：
Figure FDA0002179567550000011
其中τ为积分临时变量；则对应的复地震道c(x,y,t)为
c(x,y,t)＝f(x,y,t)+jh(x,y,t)＝A(x,y,t)ejφ(x,y,t)，
其中j为虚数单位，A(x,y,t)为瞬时振幅，φ(x,y,t)为瞬时相位；
步骤02具体包括：
利用复地震道c(x,y,t)的实部f(x,y,t)和虚部h(x,y,t)计算得到三维瞬时相位体φ(x,y,t)如下：
Figure FDA0002179567550000021
同时可由三维瞬时相位体φ(x,y,t)计算得到瞬时相位的梯度向量如下：
Figure FDA0002179567550000022
φ(x,y,t)的偏导数由下列方法计算：
Figure FDA0002179567550000023
Figure FDA0002179567550000024
Figure FDA0002179567550000025
其中max表示取最大值，ε表示阻尼因子，取值为0.001-0.01；
步骤05具体包括：
分析点(x,y,t)周围共有M个相邻分析窗，选出相似度最大的相邻分析窗序号为win_opt(x,y,t)：
Figure FDA0002179567550000026
则分析点(x,y,t)处的沿x方向及y方向的视倾角分别为Px(x,y,t)及Py(x,y,t)：
Px(x,y,t)＝px(win_opt(x,y,t),x,y,t)，
Py(x,y,t)＝py(win_opt(x,y,t),x,y,t)。
2.如权利要求1所述的联合梯度结构张量和多窗分析的地层结构曲率估计方法，其特征在于，步骤03具体包括：
假设当前的分析点为(x,y,t)，在其附近选取众多相邻分析窗；假设分析窗大小为Mx×My×Mt；其中Mx，My及Mt均为奇数，且Mx＝2Wx+1，My＝2Wy+1，Mt＝2Wt+1；得到MxMyMt个相邻分析窗，且这些相邻分析窗均必须包含当前分析点(x,y,t)。
3.如权利要求2所述的联合梯度结构张量和多窗分析的地层结构曲率估计方法，其特征在于，步骤04具体包括：
在步骤03得到的每个相邻分析窗内利用步骤02得到的相位梯度构建协方差矩阵，构建当前分析点为(x,y,t)的相位梯度协方差矩阵C(x,y,t)：
Figure FDA0002179567550000031
其中mx，my及mt均为求和临时指标；
对相位梯度协方差矩阵C(x,y,t)进行特征值分解得到三个特征值μ1，μ2，μ3及对应的三个特征向量v1，v2，v3，并有μ1≥μ2≥μ3：
Figure FDA0002179567550000032
其中μ1及v1为协方差矩阵C(x,y,t)的主特征值及主特征向量；利用主特征向量v1的三个元素v1x，v1y及v1t估计在此相邻分析窗内地层沿x方向的视倾角px(x,y,t)及沿y方向的视倾角py(x,y,t)：
Figure FDA0002179567550000041
Figure FDA0002179567550000042
利用三个特征值μ1，μ2及μ3计算相似度similarity(x,y,t)如下：
Figure FDA0002179567550000043
在每个相邻分析窗内均计算视倾角和相似度；若分析点(x,y,t)共有M个相邻分析窗，记第m个相邻分析窗内估计的沿x及y方向的视倾角分别为px(m,x,y,t)及py(m,x,y,t)，记第m个相邻分析窗内估计的相似度为similarity(m,x,y,t)。
4.如权利要求1所述的联合梯度结构张量和多窗分析的地层结构曲率估计方法，其特征在于，步骤07具体包括：
针对沿x方向的视倾角数据体Px(x,y,t)分别计算其沿x方向偏导数Pxx(x,y,t)及y方向的偏导数Pxy(x,y,t)：
Figure FDA0002179567550000044
Figure FDA0002179567550000045
针对沿y方向的视倾角数据体Py(x,y,t)分别计算其沿x方向偏导数Pyx(x,y,t)及y方向的偏导数Pyy(x,y,t)：
Figure FDA0002179567550000046
Figure FDA0002179567550000047
Δx和Δy为三维地震数据体在x和y方向上的采样间隔。
5.如权利要求1所述的联合梯度结构张量和多窗分析的地层结构曲率估计方法，其特征在于，步骤08具体包括：
利用步骤07中计算得到的三维视倾角数据体的偏导数Pxx(x,y,t)、Pxy(x,y,t)、Pyx(x,y,t)及Pyy(x,y,t)计算各点的最大正曲率MP(x,y,t)和最大负曲率MN(x,y,t)：
Figure FDA0002179567550000051
Figure FDA0002179567550000052
其中
a(x,y,t)＝0.5×Pxx(x,y,t)，
b(x,y,t)＝0.5×Pyy(x,y,t)，
c(x,y,t)＝0.5×[Pxy(x,y,t)+Pyx(x,y,t)]。
CN201810273768.4A 2018-03-29 2018-03-29 联合梯度结构张量和多窗分析的地层结构曲率估计方法 Active CN108549102B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201810273768.4A CN108549102B (zh) 2018-03-29 2018-03-29 联合梯度结构张量和多窗分析的地层结构曲率估计方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201810273768.4A CN108549102B (zh) 2018-03-29 2018-03-29 联合梯度结构张量和多窗分析的地层结构曲率估计方法

Publications (2)

Publication Number Publication Date
CN108549102A CN108549102A (zh) 2018-09-18
CN108549102B true CN108549102B (zh) 2020-03-17

Family

ID=63517466

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201810273768.4A Active CN108549102B (zh) 2018-03-29 2018-03-29 联合梯度结构张量和多窗分析的地层结构曲率估计方法

Country Status (1)

Country Link
CN (1) CN108549102B (zh)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109782339A (zh) * 2019-01-14 2019-05-21 西安交通大学 一种基于3D-DnCNN网络的叠后三维地震资料随机噪声压制方法
CN110095813B (zh) * 2019-05-17 2020-04-24 中国石油大学(北京) 用于确定复地震道的瞬时相位差的方法及装置
CN111239818B (zh) * 2020-02-12 2022-02-25 成都理工大学 一种基于三维倾角属性体校正的古地貌分析方法

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104199109A (zh) * 2014-09-18 2014-12-10 中国石油集团川庆钻探工程有限公司地球物理勘探公司 确定钻井目标层位的视倾角的方法和设备
CN105093300A (zh) * 2015-07-27 2015-11-25 中国石油天然气股份有限公司 一种地质体边界识别方法及装置
CN106199710A (zh) * 2016-06-29 2016-12-07 中国石油化工股份有限公司 基于混合倾角扫描振幅变化率的潜山储层地震识别方法

Family Cites Families (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103364828A (zh) * 2013-07-01 2013-10-23 西安交通大学 一种基于递推和乘幂法的第三代相干体算法快速实现方法
CN104020492B (zh) * 2013-07-01 2015-10-28 西安交通大学 一种三维地震资料的保边滤波方法
CN104777513B (zh) * 2015-05-11 2017-07-07 西南石油大学 地震数据梯度信息不连续性边界检测方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104199109A (zh) * 2014-09-18 2014-12-10 中国石油集团川庆钻探工程有限公司地球物理勘探公司 确定钻井目标层位的视倾角的方法和设备
CN105093300A (zh) * 2015-07-27 2015-11-25 中国石油天然气股份有限公司 一种地质体边界识别方法及装置
CN106199710A (zh) * 2016-06-29 2016-12-07 中国石油化工股份有限公司 基于混合倾角扫描振幅变化率的潜山储层地震识别方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
Robust estimates of 3D reflector dip and azimuth;Kurt J.Marfurt;《geophysics》;20060831;第71卷(第4期);第P29-P40页 *
基于GST算法的储层裂缝识别技术在某页岩气工区的应用;朱颜;《中国优秀硕士学位论文全文数据库 基础科学辑》;20170315;正文第22、26-31页 *

Also Published As

Publication number Publication date
CN108549102A (zh) 2018-09-18

Similar Documents

Publication Publication Date Title
Wu et al. Least-squares horizons with local slopes and multigrid correlations
CA2683618C (en) Inverse-vector method for smoothing dips and azimuths
US7952960B2 (en) Seismic imaging with natural Green's functions derived from VSP data
CN108139498B (zh) 具有振幅保持的fwi模型域角度叠加
US20100131205A1 (en) Method for identifying and analyzing faults/fractures using reflected and diffracted waves
US8154951B2 (en) Model-based relative bearing estimation of three-component receivers
NO330675B1 (no) Fremgangsmate for prosessering av et seismisk datasett
CN108549102B (zh) 联合梯度结构张量和多窗分析的地层结构曲率估计方法
US20150168574A1 (en) Seismic trace attribute
Wang et al. Robust seismic volumetric dip estimation combining structure tensor and multiwindow technology
Wu et al. 3D seismic image processing for unconformities
Poliannikov et al. Interferometric hydrofracture microseism localization using neighboring fracture
CN113779812A (zh) 一种大斜度井环境中利用随钻测井数据的薄夹层识别方法
Ding et al. Reliability analysis of seismic attribute in the detection of fault-karst
Yan et al. Multidirectional eigenvalue-based coherence attribute for discontinuity detection
Zeng et al. Two methods for determining geophone orientations from VSP data
US10705241B2 (en) Determining sea water resistivity
US10816688B2 (en) Method and apparatus for measuring seismic data
Domenzain et al. Joint full-waveform ground-penetrating radar and electrical resistivity inversion applied to field data acquired on the surface
CN109143398B (zh) 一种自动网格层析深度域速度的建模方法
Yan et al. Improved eigenvalue-based coherence algorithm with dip scanning
Sui et al. Seismic coherence attribute based on eigenvectors and its application
CN116068644A (zh) 一种利用生成对抗网络提升地震数据分辨率和降噪的方法
Hirabayashi Real-time event location using model-based estimation of arrival times and back azimuths of seismic phases
Wu et al. Toward accurate seismic flattening: Methods and applications

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