CN110873897A - 基于方位弹性阻抗傅里叶级数展开的裂缝预测方法及系统 - Google Patents

基于方位弹性阻抗傅里叶级数展开的裂缝预测方法及系统 Download PDF

Info

Publication number
CN110873897A
CN110873897A CN201811024016.0A CN201811024016A CN110873897A CN 110873897 A CN110873897 A CN 110873897A CN 201811024016 A CN201811024016 A CN 201811024016A CN 110873897 A CN110873897 A CN 110873897A
Authority
CN
China
Prior art keywords
azimuth
angle
elastic impedance
wave velocity
series expansion
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
Application number
CN201811024016.0A
Other languages
English (en)
Other versions
CN110873897B (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 Petroleum and Chemical Corp
Sinopec Exploration Southern Co
Original Assignee
China Petroleum and Chemical Corp
Sinopec Exploration Southern Co
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 Petroleum and Chemical Corp, Sinopec Exploration Southern Co filed Critical China Petroleum and Chemical Corp
Priority to CN201811024016.0A priority Critical patent/CN110873897B/zh
Publication of CN110873897A publication Critical patent/CN110873897A/zh
Application granted granted Critical
Publication of CN110873897B publication Critical patent/CN110873897B/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
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/60—Analysis
    • G01V2210/61—Analysis by combining or comparing a seismic data set with other data
    • G01V2210/616—Data from specific type of measurement
    • G01V2210/6169—Data from specific type of measurement using well-logging
    • 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/6226—Impedance
    • 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/64—Geostructures, e.g. in 3D data cubes
    • G01V2210/646—Fractures

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

基于方位弹性阻抗傅里叶级数展开的裂缝预测方法及系统
技术领域
本发明涉及石油天然气勘探开发领域，更具体地，涉及一种基于方位 弹性阻抗傅里叶级数展开的裂缝预测方法及系统。
背景技术
储层发育天然微裂缝不仅可以连通孤立的孔隙，增加储层的储集空间， 而且可以为油气运移提供通道；同时在页岩气勘探中也能对水平井的工程 压裂起增产作用，水力压力作用使天然微裂缝相沟通形成复杂缝网，从而 实现商业开采，因此，天然微裂缝在油气勘探中尤为重要。
在地球物理领域，基于纵波方位各向异性信息识别天然微裂缝已经成 为主流技术(Liu,2012)。利用纵波方位各向异性技术检测裂缝的最佳资料 是宽方位地震数据。近几年来，OVT域数据处理技术的发展，对叠前天然 微裂缝预测起到推动作用，每一个OVT都有自己的方位角信息，进行偏移 后，能够保留方位角的信息，避免采用CMP道集资料进行叠后偏移，OVT 地震道集资料在叠前各向异性裂缝预测具有一定的优势。
裂缝介质的振幅随方位角发生变化，基于振幅随入射角与方位变化 (AVAAz)的特点可以进行裂缝预测。前人的研究表明，裂缝介质导致地 层在不同方向的振幅、衰减、弹性性质等属性会随着方位变化而变化，基 于AVAAz可以获得不同方位的属性(衰减、频散属性或弹性参数)，进而 可以进行裂缝发育强度与走向的预测。
近几十年来，大量基于波方位地震数据的微裂缝预测方法被研发。 Mallick和Craft(1998)将方位反射振幅近似为一条周期为π的余弦曲线，进 而预测裂缝强度与李峰走向。AL-Marzong(2004)将方位AVO梯度拟合成 椭圆，认为椭圆率可以指示裂缝密度，椭圆长半轴指示裂缝走向。Downton (2006)将HTI介质的Ruger近似反射系数公式线性化，通过反演得到各 向异性梯度指示裂缝密度。上述三种方法预测裂缝走向与实际情况均有90 度的差异。孙赞东等提出衰减各项异性的裂缝预测方法；张佳佳等将各项 异性参数与各项同性参数分离，从方位弹性阻抗中提取裂缝参数，但该方 法弹性阻抗反演是单独进行，没有考虑到弹阻抗之间的结构特征，且裂缝 参数较小通常掩盖与噪声中，提取稳定性不高。Downton等(2015)将Ruger (1998)近似方程进行傅里叶级数展开，并利用其二阶傅里叶系数进行各 向异性裂缝预测，该方法主要是利用地震反射特征，不能体现层内裂缝发 育特征。
目前的叠前天然微裂缝预测方法主要分为三类，一是基于方位振幅的 差异性，二是基于振幅的衰减与频散特征，三是各项异性参数反演法，综 合目前的裂缝预测方法，国内外的相关研究方法主要存在以下问题：
①基于方位振幅属性的差异性方法：但反映的是界面的各项异性裂缝 特征，并不能反映地层内部的各向异性裂缝特征；
②基于振幅的衰减与频散特征的方法：引起衰减和频散的因素很多， 且衰减和频散参数难以准确求取，反映的也是地层界面的信息，不是地层 内部信息；
③各项异性参数反演法：各项异性参数数量级相对于弹性参数很小， 且反演维度增加，反演的稳定较差，无法准确获得。
因此，有必要开发一种基于方位弹性阻抗傅里叶级数展开的裂缝预测 方法及系统。
公开于本发明背景技术部分的信息仅仅旨在加深对本发明的一般背景 技术的理解，而不应当被视为承认或以任何形式暗示该信息构成已为本领 域技术人员所公知的现有技术。
发明内容
本发明提出了一种基于方位弹性阻抗傅里叶级数展开的裂缝预测方 法，其能够通过方位弹性阻抗傅里叶级数展开，对叠前小尺度裂缝预测， 提升裂缝预测的稳定性与精度。
根据本发明的一方面，提出了一种基于方位弹性阻抗傅里叶级数展开 的裂缝预测方法。所述方法可以包括：根据原始道集，获得包含方位角信 息的角道集；根据方位角与入射角将所述角道集进行角度划分，获得不同 方位角、不同入射角的部分叠加地震数据；根据测井的纵波速度、横波速 度、密度、各向异性信息，计算与所述部分叠加地震数据对应的每一个方 位角、每一个入射角的方位弹性阻抗曲线；根据每一个方位角、每一个入 射角的方位弹性阻抗曲线开展正演，与对应的部分叠加地震数据进行匹配 标定，提取每一个方位角、每一个入射角对应的地震子波，进而建立每一 个方位角、每一个入射角对应的弹性阻抗模型；根据所述部分叠加地震资 料、所述弹性阻抗模型、所述地震子波，计算每一个方位角、每一个入射 角对应的方位弹性阻抗体，并对方位弹性阻抗体进行傅里叶级数展开，实 现裂缝预测。
优选地，通过公式(1)计算方位弹性阻抗曲线：
Figure BDA0001787780040000031
其中，EI(θ，φ)为方位弹性阻抗曲线，g为横波速度与纵波速度的比值 的平方，Vp为纵波速度，Vs为横波速度，ρ为密度，ε(v)、δ(v)、γ分别表 示为Thomsen弱各向异性参数，θ表示地震波入射角，φ表示为测线方向与 裂缝倾向的夹角，Vp0、Vs0、ρ0分别表示为纵波速度、横波速度、密度的 均值。
优选地，通过公式(2)计算方位弹性阻抗体：
Figure BDA0001787780040000041
其中，F(m)为目标函数，G为正演算子，m为模型参数向量，d为 观测数据向量，i＝1,2，…，N，N为采样点数量，Xn为地震数据的协方差 矩阵，C为积分矩阵，ξ为相对弹性阻抗，
Figure BDA0001787780040000042
表示第i组模型参数向量，
Figure BDA0001787780040000043
表示第i组模型参数向量的协方差矩阵，XCm-ξ表示模型误差的协 方差矩阵。
优选地，通过公式(3)进行傅里叶级数展开：
ln(EIA(θ，φ))＝A0+A2cos2φ+A4cos4φ (3)
其中，EIA(θ,φ)表示归一化后的方位弹性阻抗，A0、A2与A4分别为 零阶、二阶、四阶傅里叶系数，表达式分别为：
Figure BDA0001787780040000044
Figure BDA0001787780040000045
Figure BDA0001787780040000046
其中，g为横波速度与纵波速度的比值的平方，Vp为纵波速度，Vs为 横波速度，ρ为密度，ε(v)、δ(v)、γ分别表示为Thomsen弱各向异性参数， θ表示地震波入射角。
优选地，所述部分叠加地震数据的方位角与入射角分别为划分区域的 方位角与入射角的平均值。
根据本发明的另一方面，提出了一种基于方位弹性阻抗傅里叶级数展 开的裂缝预测系统，其特征在于，该系统包括：存储器，存储有计算机可 执行指令；处理器，所述处理器运行所述存储器中的计算机可执行指令， 执行以下步骤：根据原始道集，获得包含方位角信息的角道集；根据方位 角与入射角将所述角道集进行角度划分，获得不同方位角、不同入射角的 部分叠加地震数据；根据测井的纵波速度、横波速度、密度、各向异性信 息，计算与所述部分叠加地震数据对应的每一个方位角、每一个入射角的 方位弹性阻抗曲线；根据每一个方位角、每一个入射角的方位弹性阻抗曲 线开展正演，与对应的部分叠加地震数据进行匹配标定，提取每一个方位 角、每一个入射角对应的地震子波，进而建立每一个方位角、每一个入射 角对应的弹性阻抗模型；根据所述部分叠加地震资料、所述弹性阻抗模型、 所述地震子波，计算每一个方位角、每一个入射角对应的方位弹性阻抗体， 并对方位弹性阻抗体进行傅里叶级数展开，实现裂缝预测。
优选地，通过公式(1)计算方位弹性阻抗曲线：
Figure BDA0001787780040000051
其中，EI(θ，φ)为方位弹性阻抗曲线，g为横波速度与纵波速度的比值 的平方，Vp为纵波速度，Vs为横波速度，ρ为密度，ε(v)、δ(v)、γ分别表 示为Thomsen弱各向异性参数，θ表示地震波入射角，φ表示为测线方向与 裂缝倾向的夹角，Vp0、Vs0、ρ0分别表示为纵波速度、横波速度、密度的 均值。
优选地，通过公式(2)计算方位弹性阻抗体：
Figure BDA0001787780040000061
其中，F(m)为目标函数，G为正演算子，m为模型参数向量，d为 观测数据向量，i＝1,2，…，N，N为采样点数量，Xn为地震数据的协方差 矩阵，C为积分矩阵，ξ为相对弹性阻抗，
Figure BDA0001787780040000062
表示第i组模型参数向量，
Figure BDA0001787780040000063
表示第i组模型参数向量的协方差矩阵，XCm-ξ表示模型误差的协 方差矩阵。
优选地，通过公式(3)进行傅里叶级数展开：
ln(EIA(θ，φ))＝A0+A2cos2φ+A4cos4φ (3)
其中，EIA(θ,φ)表示归一化后的方位弹性阻抗，A0、A2与A4分别为 零阶、二阶、四阶傅里叶系数，表达式分别为：
Figure BDA0001787780040000064
Figure BDA0001787780040000071
Figure BDA0001787780040000072
其中，g为横波速度与纵波速度的比值的平方，Vp为纵波速度，Vs为 横波速度，ρ为密度，ε(v)、δ(v)、γ分别表示为Thomsen弱各向异性参数， θ表示地震波入射角。
优选地，所述部分叠加地震数据的方位角与入射角分别为划分区域的 方位角与入射角的平均值。
本发明具有其它的特性和优点，这些特性和优点从并入本文中的附图 和随后的具体实施方式中将是显而易见的，或者将在并入本文中的附图和 随后的具体实施方式中进行详细陈述，这些附图和具体实施方式共同用于 解释本发明的特定原理。
附图说明
通过结合附图对本发明示例性实施例进行更详细的描述，本发明的上 述以及其它目的、特征和优势将变得更加明显，其中，在本发明示例性实 施例中，相同的参考标号通常代表相同部件。
图1示出了根据本发明的基于方位弹性阻抗傅里叶级数展开的裂缝预 测方法的步骤的流程图。
图2a、2b、2c、2d、2e、2f分别示出了根据本发明的一个实施例的入 射角为25°时，方位角为25°、55°、85°、115°、145°、170°的部 分叠加地震剖面图。
图3a、3b、3c、3d、3e、3f分别示出了根据本发明的一个实施例的入 射角为25°时，方位角为25°、55°、85°、115°、145°、170°的方 位弹性阻抗体的示意图。
图4a、4b分别示出了根据本发明的一个实施例的DY2井成像测井资料 与裂缝预测结果的示意图。
具体实施方式
下面将参照附图更详细地描述本发明。虽然附图中显示了本发明的优 选实施例，然而应该理解，可以以各种形式实现本发明而不应被这里阐述 的实施例所限制。相反，提供这些实施例是为了使本发明更加透彻和完整， 并且能够将本发明的范围完整地传达给本领域的技术人员。
图1示出了根据本发明的基于方位弹性阻抗傅里叶级数展开的裂缝预 测方法的步骤的流程图。
在该实施例中，根据本发明的基于方位弹性阻抗傅里叶级数展开的裂 缝预测方法可以包括：步骤101，根据原始道集，获得包含方位角信息的角 道集；步骤102，根据方位角与入射角将角道集进行角度划分，获得不同方 位角、不同入射角的部分叠加地震数据；步骤103，根据测井的纵波速度、 横波速度、密度、各向异性信息，计算与部分叠加地震数据对应的每一个 方位角、每一个入射角的方位弹性阻抗曲线；步骤104，根据每一个方位角、每一个入射角的方位弹性阻抗曲线开展正演，与对应的部分叠加地震数据 进行匹配标定，提取每一个方位角、每一个入射角对应的地震子波，进而 建立每一个方位角、每一个入射角对应的弹性阻抗模型；步骤105，根据部 分叠加地震资料、弹性阻抗模型、地震子波，计算每一个方位角、每一个 入射角对应的方位弹性阻抗体，并对方位弹性阻抗体进行傅里叶级数展开， 实现裂缝预测。
在一个示例中，通过公式(1)计算方位弹性阻抗曲线：
Figure BDA0001787780040000091
其中，EI(θ，φ)为方位弹性阻抗曲线，g为横波速度与纵波速度的比值 的平方，Vp为纵波速度，Vs为横波速度，ρ为密度，ε(v)、δ(v)、γ分别表 示为Thomsen弱各向异性参数，θ表示地震波入射角，φ表示为测线方向与 裂缝倾向的夹角，Vp0、Vs0、ρ0分别表示为纵波速度、横波速度、密度的 均值。
在一个示例中，通过公式(2)计算方位弹性阻抗体：
Figure BDA0001787780040000092
其中，F(m)为目标函数，G为正演算子，m为模型参数向量，d为 观测数据向量，i＝1,2，…，N，N为采样点数量，C为积分矩阵，ξ为相对 弹性阻抗，
Figure BDA0001787780040000093
表示第i组模型参数向量，gr表示对m的分组，Xn为地震 数据的协方差矩阵，
Figure BDA0001787780040000094
表示第i组模型参数向量的协方差矩阵，XCm-ξ 表示模型误差的协方差矩阵，上标-1表示对矩阵的求逆。
在一个示例中，通过公式(3)进行傅里叶级数展开：
ln(EIA(θ，φ))＝A0+A2cos2φ+A4cos4φ (3)
其中，EIA(θ,φ)表示归一化后的方位弹性阻抗，A0、A2与A4分别为 零阶、二阶、四阶傅里叶系数，表达式分别为：
Figure BDA0001787780040000101
Figure BDA0001787780040000102
Figure BDA0001787780040000103
其中，g为横波速度与纵波速度的比值的平方，Vp为纵波速度，Vs为 横波速度，ρ为密度，ε(υ)、δ(v)、γ分别表示为Thomsen弱各向异性参数， θ表示地震波入射角。
在一个示例中，部分叠加地震数据的方位角与入射角分别为划分区域 的方位角与入射角的平均值。
具体地，根据本发明的基于方位弹性阻抗傅里叶级数展开的裂缝预测 方法可以包括：
根据原始道集，收集具有方位信息的叠前地震道集资料，并对原始道 集开展去噪、剩余时差校正等道集优化处理工作，提高资料的信噪比，并 利用收集的速度场优化后的地震道集转化为具有方位角信息的角道集。
分析地震道集资料的方位角和入射角分布范围以及信噪比情况，据此 将具有方位角信息角道集转化分为至少6个方位角、3个入射角，保持每个 角度范围的覆盖次数基本相当，当数据的覆盖次数较高时，可以划分为更 多的方位角与入射角，最终将划分的角度范围内的地震数据进行叠加即可 得到不同方位角、不同入射角的部分叠加地震数据，部分叠加地震数据的 方位角与入射角为划分区域的角度的平均值。
根据测井的纵波速度、横波速度、密度、各向异性信息，计算与部分 叠加地震数据对应的每一个方位角、每一个入射角的方位弹性阻抗曲线。
在HTI介质中，在非零入射角下，不同方位地震波速度表现出差异性， 因此地震振幅随方位变化。Rüger(1998)根据一阶扰动理论给出了HTI介 质的方位AVO的反射系数近似方程方程，具体形式如下式所示：
Figure BDA0001787780040000111
其中，g为横波速度与纵波速度的比值的平方，Zp、μ、Vp分别表示介 质的纵波阻抗、剪切模量与纵波速度，Δ·表示为上下介质参数的差异(下 层减去上层)，
Figure BDA0001787780040000112
表示介质界面两侧参数的平均，θ表示地震波入射角，φ表 示为测线方向与裂缝倾向的夹角。
由于Zp＝ρVp，
Figure BDA0001787780040000113
将其代入(4)式，将按照Vp，Vs，ρ，ε(v)， δ(v)，γ进行参数合并整理得到：
Figure BDA0001787780040000114
从上式可以看出，HTI各向异性介质的反射系数表达式由各向同性背 景反射系数Riso(θ)与各向异性扰动反射系数Rani(φ，θ)组成：
R(φ，θ)＝Riso(θ)+Rani(φ，θ) (6)
其中，Riso(θ)与Rani(φ，θ)的表达式分别为：
Figure BDA0001787780040000115
Figure BDA0001787780040000121
根据Connolly(1999)对叠前弹性阻抗的定义：
Figure BDA0001787780040000122
式(6)所表示的方位反射系数可以表示为：
Figure BDA0001787780040000123
对上式两侧进行积分可以将方位弹性阻抗表示为各向同性背景介质的 弹性阻抗EIiso(θ)与各向异性扰动介质的弹性阻抗EIani(θ，φ)的乘积：
EI(θ，φ)＝EIiso(θ)EIani(θ，φ) (11)
其中，EI(θ，φ)为方位弹性阻抗；EIiso(θ)为各向同性弹性阻抗，只与 地震波的入射角有关，表达式为：
Figure BDA0001787780040000124
EIani(θ，φ)表示为各向异性扰动的弹性阻抗，不仅与地震波入射角有 关，还与观测的方位有关，根据Connolly的弹性阻抗思想，可以推导得到：
Figure BDA0001787780040000125
将(13)式与(12)式代入(11)式可以得到方位弹性阻抗的完整表 达式，即方位弹性阻抗曲线为公式(1)。
根据每一个方位角、每一个入射角的方位弹性阻抗曲线开展正演，与 对应的部分叠加地震数据进行匹配标定，提取每一个方位角、每一个入射 角对应的地震子波，进而建立每一个方位角、每一个入射角对应的弹性阻 抗模型。
根据部分叠加地震资料、弹性阻抗模型、地震子波，计算每一个方位 角、每一个入射角对应的方位弹性阻抗体，并对方位弹性阻抗体进行傅里 叶级数展开，实现裂缝预测。
将模型参数看作服从一定分布的随机变量，将似然函数与先验信息相 结合求取最大后验概率解。贝叶斯定理描述了待反演模型参数m的后验概 率分布：
p(m|d)∝p(m)p(d|m) (14)
其中，p(m|d)为后验概率，p(d|m)为似然函数，描述的是利用反演结 果正演得到的理论数据与观测数据的匹配程度，p(m)为模型参数的先验概 率。
给定模型参数向量m后，假设地面地震的观测数据向量d的噪声n服 从均值为0、协方差矩阵为Xn的高斯分布，则地面观测叠前地震数据d与 模型参数m之间的似然函数可利用噪声的分布特征进行表示：
Figure BDA0001787780040000131
对于同一套地层而言，不同方位不同入射角弹性阻抗反射系数应在相 同的时间采样点位置。因此，假设待反演地震数据的采样点为N，方位角 个数为H，入射角个数为M，我们按照时间采样点将模型参数向量m分为N 组，第i(i＝1，…，N)组向量
Figure BDA0001787780040000132
包含第i(i＝1，…，N)个采样时刻的不 同方位不同入射角的反射系数，如公式(16)所示：
Figure BDA0001787780040000133
向量
Figure BDA0001787780040000134
内成员参数服从高斯分布，组与组之间具有稀疏性，服从修正 Cauchy分布，则先验概率密度函数为：
Figure BDA0001787780040000141
上式综合考虑待反演参数的稀疏特征与内部结构特征，利用上式作为 稀疏先验信息对方位弹性阻抗反演进行约束更为合理。在反演中，为了增 强反演的稳定性，补充低频信息，还需要考虑低频模型约束：
Figure BDA0001787780040000142
其中，算子C由H×M个积分算子P，
Figure BDA0001787780040000143
沿对角线排列组成， 为不同方位角不同入射角相对弹性阻抗组成的列向量：
Figure BDA0001787780040000144
Figure BDA0001787780040000145
其中，EI0(φj，θk)(j＝1，…，H；k＝1，…，M)为第j个方位角，第k个入射角 的弹性阻抗的初采样点的值。
待反演模型参数的先验分布由公式(17)与公式(18)的乘积得到：
p(m)＝pLFM(m)pmc(m) (21)
将公式(17)与公式 (21)代入公式(16)，在后验概率极大条件下求得目标函数，即利用组稀 疏正则化方法通过公式(2)反演计算每一个方位角、每一个入射角对应的 方位弹性阻抗体。
由于(1)式所表示的方位弹性阻抗方程是非线性的，求解较为复杂， 为方便求解，首先对(1)式中常数项进行合并，令：
Figure BDA0001787780040000151
令EIA(θ，φ)＝EI(θ，φ)/A(θ)，并对(1)式等号两端分别取对数，将其 进行线性化处理即可得到：
Figure BDA0001787780040000152
由于裂缝的方位各向异性信息与观测的方位有关，且上式中的方位信 息均表示为正弦与余弦的幂的形式，可以对其进行傅里叶技术展开，进而 可以得到方位弹性阻抗体的傅里叶级数展开式为(3)。
其中，A0为零阶傅里叶系数，背景项，与观测方位无关，也包含了各 向异性信息，但各向异性信息较弱掩盖在背景的纵波速度、横波速度、密 度中；A2与A4分别为二阶与四阶傅里叶系数，只与入射角和各向异性参数 有关，均反映了裂缝的各向异性特征，入射角越大反映的各向异性信息越 显著，因此在实际的裂缝预测中需要利用较大的入射角信息进行开展相关 研究。比较A2与A4，A2的数值明显远大于A4，表明利用A4进行裂缝预测 时会严重受到噪声的影响，因此，采用二阶傅里叶系数A2进行裂缝预测， 其值越大，表明裂缝发育强度越高。
本方法通过方位弹性阻抗傅里叶级数展开，对叠前小尺度裂缝预测， 提升裂缝预测的稳定性与精度。
应用示例
为便于理解本发明实施例的方案及其效果，以下给出一个具体应用示 例。本领域技术人员应理解，该示例仅为了便于理解本发明，其任何具体 细节并非意在以任何方式限制本发明。
根据本发明的基于方位弹性阻抗傅里叶级数展开的裂缝预测方法可以 包括：
根据原始道集，收集具有方位信息的叠前地震道集资料，并对原始道 集开展去噪、剩余时差校正等道集优化处理工作，提高资料的信噪比，并 利用收集的速度场优化后的地震道集转化为具有方位角信息的角道集。
分析地震道集资料的方位角和入射角分布范围以及信噪比情况，据此 将具有方位角信息角道集转化分为6个方位角、3个入射角，保持每个角度 范围的覆盖次数基本相当，得到不同方位角、不同入射角的部分叠加地震 数据，部分叠加地震数据的方位角与入射角为划分区域的角度的平均值。
根据测井的纵波速度、横波速度、密度、各向异性信息，利用公式(1) 计算与部分叠加地震数据对应的每一个方位角、每一个入射角的方位弹性 阻抗曲线。
根据每一个方位角、每一个入射角的方位弹性阻抗曲线开展正演，与 对应的部分叠加地震数据进行匹配标定，提取每一个方位角、每一个入射 角对应的地震子波，进而建立每一个方位角、每一个入射角对应的弹性阻 抗模型。
根据部分叠加地震资料、弹性阻抗模型、地震子波，利用组稀疏正则 化方法通过公式(2)反演计算每一个方位角、每一个入射角对应的方位弹 性阻抗体，通过公式(3)进行傅里叶级数展开，利用二阶傅里叶系数实现 裂缝预测。
图2a、2b、2c、2f、2d、2f分别示出了根据本发明的一个实施例的入 射角为25°时，方位角为25°、55°、85°、115°、145°、170°的部 分叠加地震剖面图，从图中可以观察，不同方位地震数据振幅存在明显的 差异。
图3a、3b、3c、3f、3d、3f分别示出了根据本发明的一个实施例的入 射角为25°时，方位角为25°、55°、85°、115°、145°、170°的弹 性阻抗的反演结果的示意图。
图4a、4b分别示出了根据本发明的一个实施例的DY2井成像测井资料 与裂缝预测结果的示意图，从图上可知，预测结果与成像测井资料相一致， 表明本方法对裂缝预测具有较高的精度。
综上所述，本发明通过方位弹性阻抗傅里叶级数展开，对叠前小尺度 裂缝预测，提升裂缝预测的稳定性与精度。
本领域技术人员应理解，上面对本发明的实施例的描述的目的仅为了 示例性地说明本发明的实施例的有益效果，并不意在将本发明的实施例限 制于所给出的任何示例。
根据本发明的基于方位弹性阻抗傅里叶级数展开的裂缝预测系统，其 特征在于，该系统包括：存储器，存储有计算机可执行指令；处理器，处 理器运行存储器中的计算机可执行指令，执行以下步骤：根据原始道集， 获得包含方位角信息的角道集；根据方位角与入射角将角道集进行角度划 分，获得不同方位角、不同入射角的部分叠加地震数据；根据测井的纵波 速度、横波速度、密度、各向异性信息，计算与部分叠加地震数据对应的 每一个方位角、每一个入射角的方位弹性阻抗曲线；根据每一个方位角、 每一个入射角的方位弹性阻抗曲线开展正演，与对应的部分叠加地震数据 进行匹配标定，提取每一个方位角、每一个入射角对应的地震子波，进而 建立每一个方位角、每一个入射角对应的弹性阻抗模型；根据部分叠加地 震资料、弹性阻抗模型、地震子波，计算每一个方位角、每一个入射角对 应的方位弹性阻抗体，并对方位弹性阻抗体进行傅里叶级数展开，实现裂 缝预测。
在一个示例中，通过公式(1)计算方位弹性阻抗曲线：
Figure BDA0001787780040000181
其中，EI(θ，φ)为方位弹性阻抗曲线，g为横波速度与纵波速度的比值 的平方，Vp为纵波速度，Vs为横波速度，ρ为密度，ε(v)、δ(v)、γ分别表 示为Thomsen弱各向异性参数，θ表示地震波入射角，φ表示为测线方向与 裂缝倾向的夹角，Vp0、Vs0、ρ0分别表示为纵波速度、横波速度、密度的 均值。
在一个示例中，通过公式(2)计算方位弹性阻抗体：
Figure BDA0001787780040000182
其中，F(m)为目标函数，G为正演算子，m为模型参数向量，d为 观测数据向量，i＝1,2，…，N，N为采样点数量，Xn为地震数据的协方差 矩阵，C为积分矩阵，ξ为相对弹性阻抗，
Figure BDA0001787780040000183
表示第i组模型参数向量，
Figure BDA0001787780040000184
表示第i组模型参数向量的协方差矩阵，XCm-ξ表示模型误差的协 方差矩阵。
在一个示例中，通过公式(3)进行傅里叶级数展开：
ln(EIA(θ，φ))＝A0+A2cos2φ+A4cos4φ (3)
其中，EIA(θ,φ)表示归一化后的方位弹性阻抗，A0、A2与A4分别为 零阶、二阶、四阶傅里叶系数，表达式分别为：
Figure BDA0001787780040000191
Figure BDA0001787780040000192
Figure BDA0001787780040000193
其中，g为横波速度与纵波速度的比值的平方，Vp为纵波速度，Vs为 横波速度，ρ为密度，ε(v)、δ(v)、γ分别表示为Thomsen弱各向异性参数， θ表示地震波入射角。
在一个示例中，部分叠加地震数据的方位角与入射角分别为划分区域 的方位角与入射角的平均值。
本系统通过方位弹性阻抗傅里叶级数展开，对叠前小尺度裂缝预测， 提升裂缝预测的稳定性与精度。
以上已经描述了本发明的各实施例，上述说明是示例性的，并非穷尽 性的，并且也不限于所披露的各实施例。在不偏离所说明的各实施例的范 围和精神的情况下，对于本技术领域的普通技术人员来说许多修改和变更 都是显而易见的。

Claims (10)

1.一种基于方位弹性阻抗傅里叶级数展开的裂缝预测方法，包括：
根据原始道集，获得包含方位角信息的角道集；
根据方位角与入射角将所述角道集进行角度划分，获得不同方位角、不同入射角的部分叠加地震数据；
根据测井的纵波速度、横波速度、密度、各向异性信息，计算与所述部分叠加地震数据对应的每一个方位角、每一个入射角的方位弹性阻抗曲线；
根据每一个方位角、每一个入射角的方位弹性阻抗曲线开展正演，与对应的部分叠加地震数据进行匹配标定，提取每一个方位角、每一个入射角对应的地震子波，进而建立每一个方位角、每一个入射角对应的弹性阻抗模型；
根据所述部分叠加地震资料、所述弹性阻抗模型、所述地震子波，计算每一个方位角、每一个入射角对应的方位弹性阻抗体，并对方位弹性阻抗体进行傅里叶级数展开，实现裂缝预测。
2.根据权利要求1所述的基于方位弹性阻抗傅里叶级数展开的裂缝预测方法，其中，通过公式(1)计算方位方位弹性阻抗曲线：
Figure FDA0001787780030000011
其中，EI(θ，φ)为方位弹性阻抗曲线，g为横波速度与纵波速度的比值的平方，Vp为纵波速度，Vs为横波速度，ρ为密度，ε(υ)、δ(υ)、γ分别表示为Thomsen弱各向异性参数，θ表示地震波入射角，φ表示为测线方向与裂缝倾向的夹角，Vp0、Vs0、ρ0分别表示为纵波速度、横波速度、密度的均值。
3.根据权利要求1所述的基于方位弹性阻抗傅里叶级数展开的裂缝预测方法，其中，通过公式(2)计算方位弹性阻抗体：
Figure FDA0001787780030000021
其中，F(m)为目标函数，G为正演算子，m为模型参数向量，d为观测数据向量，i＝1,2，…，N，N为采样点数量，Xn为地震数据的协方差矩阵，C为积分矩阵，ξ为相对弹性阻抗，
Figure FDA0001787780030000022
表示第i组模型参数向量，
Figure FDA0001787780030000023
表示第i组模型参数向量的协方差矩阵，XCm-ξ表示模型误差的协方差矩阵。
4.根据权利要求1所述的基于方位弹性阻抗傅里叶级数展开的裂缝预测方法，其中，通过公式(3)进行傅里叶级数展开：
ln(EIA(θ，φ))＝A0+A2cos2φ+A4cos4φ (3)
其中，A0、A2与A4分别为零阶、二阶、四阶傅里叶系数，表达式分别为：
Figure FDA0001787780030000031
Figure FDA0001787780030000032
Figure FDA0001787780030000033
其中，EIA(θ,φ)表示归一化后的方位弹性阻抗，g为横波速度与纵波速度的比值的平方，Vp为纵波速度，Vs为横波速度，ρ为密度，ε(v)、δ(v)、γ分别表示为Thomsen弱各向异性参数，θ表示地震波入射角。
5.根据权利要求1所述的基于方位弹性阻抗傅里叶级数展开的裂缝预测方法，其中，所述部分叠加地震数据的方位角与入射角分别为划分区域的方位角与入射角的平均值。
6.一种基于方位弹性阻抗傅里叶级数展开的裂缝预测系统，其特征在于，该系统包括：
存储器，存储有计算机可执行指令；
处理器，所述处理器运行所述存储器中的计算机可执行指令，执行以下步骤：
根据原始道集，获得包含方位角信息的角道集；
根据方位角与入射角将所述角道集进行角度划分，获得不同方位角、不同入射角的部分叠加地震数据；
根据测井的纵波速度、横波速度、密度、各向异性信息，计算与所述部分叠加地震数据对应的每一个方位角、每一个入射角的方位弹性阻抗曲线；
根据每一个方位角、每一个入射角的方位弹性阻抗曲线开展正演，与对应的部分叠加地震数据进行匹配标定，提取每一个方位角、每一个入射角对应的地震子波，进而建立每一个方位角、每一个入射角对应的弹性阻抗模型；
根据所述部分叠加地震资料、所述弹性阻抗模型、所述地震子波，计算每一个方位角、每一个入射角对应的方位弹性阻抗体，并对方位弹性阻抗体进行傅里叶级数展开，实现裂缝预测。
7.根据权利要求6所述的基于方位弹性阻抗傅里叶级数展开的裂缝预测系统，其中，通过公式(1)计算方位弹性阻抗曲线：
Figure FDA0001787780030000041
其中，EI(θ，φ)为方位弹性阻抗曲线，g为横波速度与纵波速度的比值的平方，Vp为纵波速度，Vs为横波速度，ρ为密度，ε(υ)、δ(v)、γ分别表示为Thomsen弱各向异性参数，θ表示地震波入射角，φ表示为测线方向与裂缝倾向的夹角，Vp0、Vs0、ρ0分别表示为纵波速度、横波速度、密度的均值。
8.根据权利要求6所述的基于方位弹性阻抗傅里叶级数展开的裂缝预测系统，其中，通过公式(2)计算方位弹性阻抗体：
Figure FDA0001787780030000051
其中，F(m)为目标函数，G为正演算子，m为模型参数向量，d为观测数据向量，i＝1,2，…，N，N为采样点数量，Xn为地震数据的协方差矩阵，C为积分矩阵，ξ为相对弹性阻抗，
Figure FDA0001787780030000055
表示第i组模型参数向量，
Figure FDA0001787780030000056
表示第i组模型参数向量的协方差矩阵，XCm-ξ表示模型误差的协方差矩阵。
9.根据权利要求6所述的基于方位弹性阻抗傅里叶级数展开的裂缝预测系统，其中，通过公式(3)进行傅里叶级数展开：
ln(EIA(θ，φ))＝A0+A2cos2φ+A4cos4φ (3)
其中，EIA(θ,φ)表示归一化后的方位弹性阻抗，A0、A2与A4分别为零阶、二阶、四阶傅里叶系数，表达式分别为：
Figure FDA0001787780030000052
Figure FDA0001787780030000053
Figure FDA0001787780030000054
其中，g为横波速度与纵波速度的比值的平方，Vp为纵波速度，Vs为横波速度，ρ为密度，ε(v)、δ(v)、γ分别表示为Thomsen弱各向异性参数，θ表示地震波入射角。
10.根据权利要求6所述的基于方位弹性阻抗傅里叶级数展开的裂缝预测系统，其中，所述部分叠加地震数据的方位角与入射角分别为划分区域的方位角与入射角的平均值。
CN201811024016.0A 2018-09-04 2018-09-04 基于方位弹性阻抗傅里叶级数展开的裂缝预测方法及系统 Active CN110873897B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201811024016.0A CN110873897B (zh) 2018-09-04 2018-09-04 基于方位弹性阻抗傅里叶级数展开的裂缝预测方法及系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201811024016.0A CN110873897B (zh) 2018-09-04 2018-09-04 基于方位弹性阻抗傅里叶级数展开的裂缝预测方法及系统

Publications (2)

Publication Number Publication Date
CN110873897A true CN110873897A (zh) 2020-03-10
CN110873897B CN110873897B (zh) 2021-07-27

Family

ID=69716249

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201811024016.0A Active CN110873897B (zh) 2018-09-04 2018-09-04 基于方位弹性阻抗傅里叶级数展开的裂缝预测方法及系统

Country Status (1)

Country Link
CN (1) CN110873897B (zh)

Cited By (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111897006A (zh) * 2020-07-15 2020-11-06 中国石油大学(华东) 一种基于方位弹性阻抗差异奇异值分解的裂缝密度及方向预测方法及系统与应用
CN113312583A (zh) * 2021-05-26 2021-08-27 中南大学 裂缝性页岩气储层性能参数预测方法及系统
CN113740896A (zh) * 2021-07-30 2021-12-03 中国石油大学(华东) 一种预测正交各向异性储层裂缝的方法及存储介质
CN113835119A (zh) * 2021-09-18 2021-12-24 中国石油大学(华东) 基于方位扩展弹性阻抗的页岩储层水平地应力预测方法
CN114002738A (zh) * 2021-11-02 2022-02-01 中海石油(中国)有限公司 高角度断裂约束的叠前裂缝预测方法、装置及存储介质
CN114019565A (zh) * 2021-11-02 2022-02-08 中海石油(中国)有限公司 快速叠前裂缝预测方法、装置、计算机存储介质及设备
US20230124730A1 (en) * 2021-10-14 2023-04-20 Halliburton Energy Services, Inc. Fracture Geometry And Orientation Identification With A Single Distributed Acoustic Sensor Fiber
CN117331123A (zh) * 2023-11-06 2024-01-02 成都理工大学 缝洞密度反演方法及缝洞储层预测系统

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104597491A (zh) * 2015-01-28 2015-05-06 中国石油大学(华东) 一种基于方位弹性阻抗差异的裂缝储层参数反演方法
CN106772578A (zh) * 2016-12-07 2017-05-31 中国矿业大学(北京) 一种合成地震记录的方法和装置
US20170160413A1 (en) * 2015-12-04 2017-06-08 Cgg Services Sas Method and apparatus for analyzing fractures using avoaz inversion
CN106842313A (zh) * 2015-12-04 2017-06-13 中国石油化工股份有限公司 基于方位叠前地震数据的各向异性参数反演方法
CN108363100A (zh) * 2018-01-16 2018-08-03 太原理工大学 基于封存单元与粗糙集的煤层气甜点区地震地质识别方法
US10816686B2 (en) * 2015-07-28 2020-10-27 Schlumberger Technology Corporation Seismic constrained discrete fracture network

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104597491A (zh) * 2015-01-28 2015-05-06 中国石油大学(华东) 一种基于方位弹性阻抗差异的裂缝储层参数反演方法
US10816686B2 (en) * 2015-07-28 2020-10-27 Schlumberger Technology Corporation Seismic constrained discrete fracture network
US20170160413A1 (en) * 2015-12-04 2017-06-08 Cgg Services Sas Method and apparatus for analyzing fractures using avoaz inversion
CN106842313A (zh) * 2015-12-04 2017-06-13 中国石油化工股份有限公司 基于方位叠前地震数据的各向异性参数反演方法
CN106772578A (zh) * 2016-12-07 2017-05-31 中国矿业大学(北京) 一种合成地震记录的方法和装置
CN108363100A (zh) * 2018-01-16 2018-08-03 太原理工大学 基于封存单元与粗糙集的煤层气甜点区地震地质识别方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
HUAIZHEN CHEN ET AL.: "Azimuthal elastic impedance inversion for fluid term and fracture weakness", 《CREWES RESEARCH REPORT》 *
HUAIZHEN CHEN ET AL.: "Estimation of modified fluid factor and dry fracture weaknesses using azimuthal elastic impedance", 《GEOPHYSICS》 *
王康宁等: "基于傅里叶级数展开的纵波方位各向异性裂缝预测", 《石油物探》 *

Cited By (12)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111897006A (zh) * 2020-07-15 2020-11-06 中国石油大学(华东) 一种基于方位弹性阻抗差异奇异值分解的裂缝密度及方向预测方法及系统与应用
CN113312583A (zh) * 2021-05-26 2021-08-27 中南大学 裂缝性页岩气储层性能参数预测方法及系统
CN113740896A (zh) * 2021-07-30 2021-12-03 中国石油大学(华东) 一种预测正交各向异性储层裂缝的方法及存储介质
CN113740896B (zh) * 2021-07-30 2022-06-03 中国石油大学(华东) 一种预测正交各向异性储层裂缝的方法及存储介质
CN113835119A (zh) * 2021-09-18 2021-12-24 中国石油大学(华东) 基于方位扩展弹性阻抗的页岩储层水平地应力预测方法
US20230124730A1 (en) * 2021-10-14 2023-04-20 Halliburton Energy Services, Inc. Fracture Geometry And Orientation Identification With A Single Distributed Acoustic Sensor Fiber
CN114002738A (zh) * 2021-11-02 2022-02-01 中海石油(中国)有限公司 高角度断裂约束的叠前裂缝预测方法、装置及存储介质
CN114019565A (zh) * 2021-11-02 2022-02-08 中海石油(中国)有限公司 快速叠前裂缝预测方法、装置、计算机存储介质及设备
CN114002738B (zh) * 2021-11-02 2023-08-15 中海石油(中国)有限公司 高角度断裂约束的叠前裂缝预测方法、装置及存储介质
CN114019565B (zh) * 2021-11-02 2023-10-13 中海石油(中国)有限公司 快速叠前裂缝预测方法、装置、计算机存储介质及设备
CN117331123A (zh) * 2023-11-06 2024-01-02 成都理工大学 缝洞密度反演方法及缝洞储层预测系统
CN117331123B (zh) * 2023-11-06 2024-04-02 成都理工大学 缝洞密度反演方法及缝洞储层预测系统

Also Published As

Publication number Publication date
CN110873897B (zh) 2021-07-27

Similar Documents

Publication Publication Date Title
CN110873897B (zh) 基于方位弹性阻抗傅里叶级数展开的裂缝预测方法及系统
CN104407378B (zh) 一种各向异性参数反演方法及装置
CN109425896B (zh) 白云岩油气储层分布预测方法及装置
CN106772604B (zh) 基于流体体积压缩系数的叠前地震反演方法
CN106556861B (zh) 一种基于全方位地震资料的方位avo反演方法
CN103487831B (zh) Avo地震正演计算方法
CN103713323B (zh) 一种全方位各向异性保幅成像与抽道集方法
CN103513277B (zh) 一种地震地层裂隙裂缝密度反演方法及系统
CN103424776A (zh) 一种碳酸盐岩油气储层裂缝地震检测方法
CN104459777A (zh) 基于流体体积模量avo反演的流体识别的方法及系统
CN111366975B (zh) 基于交叉梯度正则化约束的叠前地震ava反演方法
CN103675911B (zh) 一种基于纵波和转换横波联合反演截距和梯度的方法
CN104459778A (zh) 一种基于双相介质固液解耦的叠前地震反演方法及系统
CN105319589A (zh) 一种利用局部同相轴斜率的全自动立体层析反演方法
CN111506861B (zh) 一种目的层有利区域裂缝强度计算方法
CN109521469B (zh) 一种海底沉积物弹性参数的正则化反演方法
CN110618450B (zh) 基于岩石物理建模的致密储层智能化含气性预测方法
Thiel et al. Comparison of acoustic and elastic full‐waveform inversion of 2D towed‐streamer data in the presence of salt
CN109143351A (zh) 叠前各向异性特征参数反演方法及计算机可读存储介质
US9658354B2 (en) Seismic imaging systems and methods employing correlation-based stacking
CN109143352B (zh) 一种各向异性介质地震反射特征方程建立方法
CN112305602A (zh) 基于叠前多属性与古地貌融合技术的碳酸岩储层预测方法
Sun et al. Fractured reservoir modeling by discrete fracture network and seismic modeling in the Tarim Basin, China
CN102520445B (zh) 一种利用松弛因子叠前地震反演进行储层预测的方法
CN105487113B (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
GR01 Patent grant
GR01 Patent grant