地震叠前三参数同步反演系统

资 源 简 介

本项目旨在基于MATLAB平台开发一套完整的地震叠前AVO（振幅随偏移距变化）反演工具，核心目标是利用叠前地震道集数据同时反演地层的纵波速度（Vp）、横波速度（Vs）和密度（Density）三个关键弹性参数。系统主要依据Zoeppritz方程的线性近似（如Aki-Richards近似或Fatti近似）建立反射系数与弹性参数之间的数学关系，结合地震子波和低频模型构建正演褶积矩阵。在反演策略上，项目采用基于贝叶斯理论的统计反演框架或阻尼最小二乘法（Damped Least Squares），引入柯西分布、高斯分布或全变分（TV）正则化作为先验约束，以解决反演过程中的不适定问题并提高抗噪能力。功能模块包括：1. 数据加载与预处理模块，支持读取SEGY格式的叠前角道集或偏移距道集，并进行角度叠加与平滑；2. 子波提取与标定模块，利用测井数据提取统计子波或确定性子波；3. 低频模型构建模块，基于测井插值和层位约束生成背景低频趋势；4. 反演核心算法模块，实现目标函数的构建与迭代优化求解；5. 结果可视化模块，能够并排显示反演的Vp、Vs、密度剖面以及合成地震记录与原始数据的残差对比，用于直接评估反演质量。该系统适用于油气勘探开发中的岩性识别和流体检测，特别是在致密砂岩储层预测和含气性检测方面提供重要技术支撑。

详 情 说 明

地震叠前三参数同步反演系统 (MATLAB)

项目简介

本项目是一个基于MATLAB平台开发的地震叠前AVO（振幅随偏移距变化）反演演示系统。该系统的核心功能是利用叠前地震角道集数据，通过物理模型驱动的方式，同时反演地层的纵波阻抗（Zp）、横波阻抗（Zs）和密度（Density）三个关键弹性参数。

代码目前实现了一个完整的闭环流程：从合成地质模型的构建、基于Fatti近似的叠前正演模拟，到构建大型稀疏矩阵方程进行同步反演，最后对反演结果进行可视化评估。该工具主要用于演示基于线性近似的地震反演原理，特别是在处理砂岩、含气砂岩及致密灰岩等复杂岩性时的反演效果。

功能特性

• 合成地质模型生成：能够自动生成包含泥岩背景、高速砂岩、含气砂岩（低Vp/低密度）和致密灰岩的多层状地质模型。
• 低频模型构建：实现了基于滑动平均滤波的低频趋势提取，模拟实际反演中缺失低频信息的情况，并将其作为贝叶斯反演的先验约束。
• 叠前AVO正演模拟：

* 生成雷克（Ricker）子波。 * 基于Fatti方程计算不同角度（近、中、远）的AVO反射系数。 * 实现褶积模型并添加指定信噪比的随机噪声。

• 三参数同步反演算法：

* 构建基于阻尼最小二乘法（Damped Least Squares）的目标函数。 * 构建大型稀疏核矩阵 G，将反演问题转化为线性方程组求解。 * 支持纵波阻抗、横波阻抗和密度的加权正则化约束，以解决密度反演的不稳定性。

• 结果验证与分析：并排展示真实模型、先验低频模型与最终反演结果，并计算合成地震记录与观测数据的残差。

系统要求

• MATLAB R2016a 或更高版本。
• Signal Processing Toolbox（用于 filtfilt 等信号处理函数）。

使用方法

1. 确保MATLAB当前工作目录包含 main.m 文件。
2. 直接运行 main 函数。
3. 程序将依次输出初始化状态、正演模拟进度、反演求解状态。
4. 运行结束后，系统将自动生成包含地层参数剖面对比和地震道集残差分析的图表。

代码逻辑与实现细节

main.m 是整个系统的入口和核心控制器，其内部逻辑严格遵循以下步骤：

1. 参数设置与地质建模

程序首先定义了采样间隔（2ms）、样点数（300）、三个叠前角度（10°, 20°, 30°）以及信噪比。 随后调用 generate_earth_model 函数构建真实的弹性参数模型（Vp, Vs, Density）。该模型硬编码了四种岩性：

• 背景泥岩
• 高速砂岩
• 含气砂岩（表现为由于流体导致的Vp和密度下降）
• 致密灰岩（表现为极高的速度和密度）

2. 先验模型与去噪

为了模拟真实反演场景，代码通过 build_low_freq_model 函数利用 filtfilt 对真实模型进行平滑处理，生成缺乏高频细节的低频模型（Prior Model），这将在反演方程中作为初始猜测和背景约束。

3. 基于Fatti近似的正演模拟

正演过程模拟了地震波在不同入射角下的反射特征：

• 子波生成：生成指定主频的零相位雷克子波。
• AVO系数计算：计算 Fatti 近似下的三个系数（c1, c2, c3）。计算过程依赖于背景Vp/Vs比率，反映了纵波阻抗、横波阻抗和密度项对反射系数的贡献权重。
• 褶积生成：分别计算纵波、横波和密度的反射率分量，利用加权叠加生成不同角度的反射系数序列，并与子波进行“相同（same）”模式的褶积，最后通过RMS标准添加随机噪声。

4. 反演矩阵方程构建

反演模块的核心在于构建线性方程组 $d = G cdot m$。

• 模型参数 ($m$)：待求解参数为纵波阻抗、横波阻抗和密度的自然对数（ln(Zp), ln(Zs), ln(Rho)）。
• 核矩阵 ($G$)：这是一个大型稀疏矩阵，物理意义涵盖了从阻抗域到反射系数域的差分操作，以及从反射系数到地震数据域的褶积操作。代码通过 build_inversion_matrix 函数构建了差分矩阵 $D$ 和褶积矩阵 $W$，并将它们与Fatti系数结合，组装成针对近、中、远三个角度子矩阵的复合矩阵。
• 正则化 ($W_m$)：为了克服反演的不适定性，特别是密度项的不稳定性，代码引入了加权矩阵。通过设置 $lambda_{zp}, lambda_{zs}, lambda_{rho}$，对不同参数施加不同程度的约束（例如代码中对密度项给予了更高的权重 0.5）。

5. 求解与后处理

系统采用阻尼最小二乘法求解模型扰动量，即求解 $Delta m$ 使得目标函数 $J(m) = ||d - Gm||^2 + ||W_m(m - m_{prior})||^2$ 最小化。 求解得到的扰动量叠加回先验模型（m_inverted = m_prior + perturbation），最终恢复出绝对的弹性参数。最后一步是利用反演结果重新进行正演，计算预测数据与观测数据的残差，用于质量控制（QC）。

关键算法说明

• Fatti AVO 方程：代码使用了 Fatti (1994) 提出的近似公式，将反射系数线性化为纵波阻抗反射率、横波阻抗反射率和密度反射率的加权和。
• 差分算子：在构建反演矩阵时，利用 $R_i approx 0.5 times (ln(Z_i) - ln(Z_{i-1}))$ 的关系，将层状阻抗参数转化为界面反射系数。
• 稀疏矩阵优化：为了提高计算效率，系统利用MATLAB的 sparse 和 spdiags 函数构建核矩阵，避免了存储大量的零元素。