CN110187382B - 一种回折波和反射波波动方程旅行时反演方法 - Google Patents

Abstract

本发明公开了一种回折波和反射波波动方程旅行时反演方法，包括地震资料回折波和反射波分离；回折波和反射波旅行时残差求取；旅行时残差目标函数构建；背景速度与反射系数梯度计算；反射系数模型更新迭代；背景速度模型更新迭代等技术环节。本发明通过采用回折波和反射波进行波动方程旅行时反演来获取可靠的地下速度初始模型，进而提高后续全波形反演和逆时偏移的反演和成像精度。

Description

一种回折波和反射波波动方程旅行时反演方法

技术领域

本发明属于地震波反演技术领域，具体涉及一种回折波和反射波波动方程 旅行时反演方法。

背景技术

随着勘探程度的深入，地震勘探的目标已从大型构造油气藏向中小型复杂 隐蔽型油气藏转变。面对日趋复杂的勘探目标，地震勘探的难度也随之增加， 对当前的偏移方法提出新的挑战。逆时偏移是一种叠前、深度域、双程波波动 方程偏移方法，假设条件少、精度高、无倾角限制，已成为复杂构造条件下最 为常用的成像方法。地震反演通过地表观测到的数据估算速度、密度等地球物 理参数，进而反推地下结构和构造，可以有效地识别地质构造、预测自然灾害 和勘探油气藏。与其他反演方法相比，全波形反演利用地震波的动力学和运动 学信息，是地球深部构造成像、浅层环境调查及油气勘探速度建模的有效工具。

逆时偏移和全波形反演方法对初始速度模型有较强的依懒性。当初始模型 与真实模型相差较远时，很难得到满意的成像和反演结果。偏移速度分析、走 时层析、Laplace域波形反演、包络波形反演等可用于初始速度建模，但都有一 定局限性。回折波本质上是一种直达波，波形反演中正向和反向传播波场的夹 角接近180度(对应低波数分量)，可以获得较好的低频速度模型。但回折波的 穿透深度有限，需要较大的偏移距，实际中往往很难满足。深层速度建模还需 要依靠反射波场。反射波波形反演采用反射波场来获取地下模型低频分量的新 方法。反射波波形反演中，模型参数被拆分为背景参数和扰动参数/反射系数， 反射波的振幅和相位/走时分别由反射系数和背景速度场决定，通过交替进行逆 时偏移/最小二乘逆时偏移和波形反演来更新反射系数和背景速度模型。与常规 全波形反演相比，反射波波形反演的敏感核/梯度更加光滑，能获得可靠的低频 速度模型。但现有方法大多单独采用回折波和反射波进行浅层和深部速度建模。

发明内容

本发明的目的在于是获取纵波速度低频模型，为后续地震资料弹性波全 波形反演和逆时偏移提供可靠的初始模型。

一种回折波和反射波波动方程旅行时反演方法，包括以下步骤：

步骤(1)：地震资料回折波和反射波分离；

步骤(2)：采用动态图像规整法获取回折波和反射波模拟数据与观测数据 的旅行时残差；

步骤(3)：构建旅行时残差目标函数；

步骤(4)：计算背景速度与反射系数梯度；

步骤(5)：反射系数模型更新迭代；

步骤(6)：背景速度模型更新迭代；

步骤(7)：重复上面的步骤，直到满足收敛条件；

其中，步骤(3)中旅行时残差目标函数构建方法如下：

采用如下的目标函数：

其中，τ₁和τ₂分别表示回折波和反射波模拟和观测数据的旅行时残差；T为最 大记录时间；

采用动态图像规整法(dynamic image warping)来估算局部旅行时残差；极 小化函数为：

其中，p_obs和Δp_obs为观测回折波和反射波，p和Δp为模拟回折波和反射波；

方程(3)关于自变量求偏导数等于0得：

其中，变量上方的“点”表示对时间的偏导数；

基于方程(4)和隐函数的求导法则可以得到：

同理可得：

伴随方程的虚源为目标函数对模拟波场的偏导数，则

进一步地，如上所述的回折波和反射波波动方程旅行时反演方法，步骤(4) 中背景速度与反射系数梯度计算公式如下：

目标函数对背景参数的梯度为：

目标函数对反射系数的梯度为：

通过方程(17)和方程(18)可得到背景速度和反射系数的梯度；

其中，Δp为反射波场，p为背景波场；v为背景速度；m为反射系数；q和Δq为 拉格朗日乘子函数。

进一步地，如上所述的回折波和反射波波动方程旅行时反演方法，步骤(5) 反射系数反演迭代步骤如下：

只采用反射波场更新反射系数，目标函数变为：

反射系数反演迭代步骤如下：

a.求解方程)得到正向背景波场p和扰动波场Δp，初始条件为：

b.求解伴随方程得到反向背景波场q，终值条件为：

c.通过方程(18)计算目标函数关于反射参数的梯度；

d.对反射系数梯度进行预条件处理；

e.采用抛物线拟合法计算迭代步长，对反射系数模型进行更新迭代。

进一步地，如上所述的回折波和反射波波动方程旅行时反演方法，步骤(6) 背景速度模型更新迭代方法如下：

a.求解方程得到正向背景波场p和扰动波场Δp，初始条件为：

b.求解伴随方程得到反向背景波场q，终值条件为：

c.通过方程(17)计算由背景场决定的背景速度梯度

d.求解伴随方程得到反向背景波场q和扰动波场Δq，终值条件为：

e.通过方程(17)计算由反射场决定的背景速度梯度

f.将步骤c和e中的梯度加权求和得到总的背景速度梯度；

g.对背景速度梯度进行预条件处理；

h.采用抛物线拟合法计算迭代步长，对背景速度模型进行更新迭代。

有益效果：

本发明通过回折波和反射波波动方程旅行时反演方法，来获取速度模型 的低频分量，进而提高后续地震资料全波形反演和逆时偏移的反演和成像精 度。

附图说明

图1为为本发明回折波和反射波波动方程旅行时反演的流程图；

图2为实施例1中双层介质模型；

图3(a)为实施例1中双层介质模型的反射波旅行时反演方法的敏感核；

图3(b)为实施例1中双层介质模型的回折波旅行时反演方法的敏感核；

图3(c)为实施例1中双层介质模型的回折波和反射波波动方程旅行时反 演方法的敏感核；

图4为实施例2中Sigsbee2A模型

图5(a)为实施例2中Sigsbee2A初始模型为随深度线性变化的一维模型；

图5(b)为实施例2中Sigsbee2A模型回折波和反射波波动方程旅行时反 演方法的反演结果；

图5(c)为实施例2中Sigsbee2A模型常规反射波波形反演结果。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚，下面本发明中的技术方 案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例， 而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有 作出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

反射波波形反演是一种有效的初始模型建模手段，而回折波在建立浅层速 度模型方面具有一定优势。因此，同时采用反射波和回折波进行波形反演可以 改善低频模型的反演精度。本发明将给出一种新的回折波和反射波波动方程旅 行时反演方法，来获取速度模型的低频分量，进而提高后续地震资料全波形反 演和逆时偏移的反演和成像精度。

如图1所示，为本发明实施回折波和反射波波动方程旅行时反演的流程图， 具体包括：

(1)地震资料回折波和反射波分离；

(2)回折波和反射波旅行时残差求取；

(3)旅行时残差目标函数构建；

(4)背景速度与反射系数梯度计算；

(5)反射系数模型更新迭代；

(6)背景速度模型更新迭代；

(7)重复上面的步骤，直到满足收敛条件。

步骤(1)中地震资料回折波和反射波分离方法如下：

为避免全波形反演中高频偏移等时线，必须将反射波和回折波分离。实际 地震资料中回折波震相清晰、且多以初至波形式出现，可采用现有的初至拾取 方法进行识别，基于初至时间将回折波和反射波分离。在Born近似理论框架下， 模拟数据的回折波和反射波分别位于背景波场和扰动波场中(自然分离)。模拟 反射波中只包含一次反射波，实际资料需要进行多次反射波去除。

步骤(2)回折波和反射波旅行时残差求取方法如下：

窗口互相关或动态图像规整法(dynamic image warping)(Ma和Hale，2013) 可以拾取模拟与观测数据的旅行时时差。与互相关方法相比，动态图像规整法 在时差变化剧烈时精度更高。本发明中采用动态图像规整法获取回折波和反射 波模拟数据与观测数据的旅行时残差。

步骤(3)中旅行时残差目标函数构建方法如下：

采用如下的目标函数：

其中，τ₁和τ₂分别表示回折波和反射波模拟和观测数据的旅行时残差。T为最 大记录时间。

采用动态图像规整法(dynamic image warping)来估算局部旅行时残差。极 小化函数为：

其中，p_obs和Δp_obs为观测回折波和反射波，p和Δp为模拟回折波和反射波。

方程(3)关于自变量求偏导数等于0得：

其中，变量上方的“点”表示对时间的偏导数。

基于方程(4)和隐函数的求导法则可以得到：

同理可得：

伴随方程的虚源为目标函数对模拟波场的偏导数，则

步骤(4)中背景速度与反射系数梯度计算方法如下：

二维声波方程为：

其中，v为传播速度，s为震源项。

对于背景速度v，背景波场p(包括回折波)可以通过求解方程(9)得到。当 存在模型扰动Δv时，扰动/反射波场为Δp，且满足：

方程(10)与方程(9)相减，化简并忽略高阶微小量可得：

其中，

可近似表示反射系数。

对给定的背景速度v和反射系数m，求解方程(9)和方程(11)可以得到反射波 Δp，即为反偏移过程。

方程(1)给出基于旅行时的目标函数，背景波场p满足方程(9)，反射波场Δp 满足方程(11)。采用拉格朗日乘子法求解该约束优化问题。目标函数变为：

其中，q和Δq为拉格朗日乘子函数，

和

分部积分方程(12)并令

和

得：

和

其中，背景场残差反传时s′＝-s₁；反射场残差反传时s′＝-s₂。

目标函数对背景参数的梯度为：

目标函数对反射系数的梯度为：

通过方程(17)和方程(18)可得到背景速度和反射系数的梯度。为得到较好 的反演效果，本发明中反射系数和背景速度交替进行反演迭代。

步骤(5)反射系数模型更新迭代方法如下：

只采用反射波场更新反射系数，目标函数变为：

反射系数反演迭代步骤如下：

a.求解方程(9)和方程(11)得到正向背景波场p和扰动波场Δp，初始条件为：

b.求解伴随方程(15)(s′＝-s₂)得到反向背景波场q，终值条件为：

c.通过方程(18)计算目标函数关于反射参数的梯度；

d.对反射系数梯度进行预条件(逆海森矩阵)处理；

e.采用抛物线拟合法计算迭代步长，对反射系数模型进行更新迭代。

步骤(6)背景速度模型更新迭代方法如下：

为了避免高频偏移等时线的影响，分别进行背景场(包含回折波)和反射场残 差反向延拓，背景速度反演迭代步骤如下：

a.求解方程(9)和方程(11)得到正向背景波场p和扰动波场Δp，初始条件为：

b.求解伴随方程(15)(s′＝-s₁)得到反向背景波场q(背景场残差反向延拓)， 终值条件为：

c.通过方程(17)计算由背景场决定的背景速度梯度

d.求解伴随方程(15)和(16)(s′＝-s₂)得到反向背景波场q和扰动波场Δq(反 射场残差反向延拓)，终值条件为：

e.通过方程(17)计算由反射场决定的背景速度梯度

f.将步骤c和e中的梯度加权求和(17)得到总的背景速度梯度；

g.对背景速度梯度进行预条件(逆海森矩阵)处理；

h.采用抛物线拟合法计算迭代步长，对背景速度模型进行更新迭代。

步骤(7)中重复上面的步骤，直到满足收敛条件。本发明中采用目标函数小 于某极小值(如初始残差的1/10)和迭代次数大于30次两种迭代终止条件。

下面通过几个例子来分析本发明中提出的回折波和反射波波动方程旅行时 反演方法的有效性。

实施例1：

首先以一个炮点-检波点对为例来说明本发明的优势。双层介质模型如图2 所示，炮点和检波点位于三角形指示的位置。时间步长1ms，空间间隔10m， 震源为15Hz的雷克子波。图3为不同反演方法的敏感核，由图3(a)-3(b)可知， 回折波和反射波波动方程旅行时反演可以结合常规反射波反演和回折波反演的 优点，得到更好的低频敏感核，可同时进行近地表和深层初始速度模型建模。

实施例2：

下面采用复杂的Sigsbee2A模型(如图4所示)对新提出的回折波和反射波 波动方程旅行时反演方法进行测试。时间步长1.5ms，空间间隔16m，炮点(36 炮)和检波点均匀分布于地表。震源为10Hz的雷克子波。初始模型为随深度线 性变化的一维模型(与真实模型相差较远)，如图5(a)所示。图5(b)给出Sigsbee2A 模型回折波和反射波波动方程旅行时反演方法的反演结果。图5(c)给出常规反射 波波形反演结果。由图可知，通过本发明中提出的反射波旅行时反演方法可以 获取较好的速度初始模型，进而采用常规全波形反演来恢复模型的高频分量。

最后应说明的是：以上实施例仅用以说明本发明的技术方案，而非对其限 制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员 应当理解：其依然可以对前述各实施例所记载的技术方案进行修改，或者对其 中部分技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的 本质脱离本发明各实施例技术方案的精神和范围。

Claims (2)

1.一种回折波和反射波波动方程旅行时反演方法，其特征在于，包括以下步骤：

步骤(1)：地震资料回折波和反射波分离；

步骤(2)：采用动态图像规整法获取回折波和反射波模拟数据与观测数据的旅行时残差；

步骤(3)：构建旅行时残差目标函数；

步骤(4)：计算背景速度与反射系数梯度；

步骤(5)：反射系数模型更新迭代；

步骤(6)：背景速度模型更新迭代；

步骤(7)：重复上面的步骤，直到满足收敛条件；

其中，步骤(3)中旅行时残差目标函数构建方法如下：

采用如下的目标函数：

其中，τ₁表示模拟回折波和观测回折波的旅行时残差；τ₂表示模拟反射波和观测反射波的旅行时残差；T为最大记录时间；

采用动态图像规整法(dynamic image warping)来估算局部旅行时残差；极小化函数为：

其中，d_obs和Δd_obs分别为观测回折波和观测反射波，d和Δd分别为模拟回折波和模拟反射波；

方程(3)关于自变量求偏导数等于0得：

其中，变量上方的“点”表示对时间的偏导数；

基于方程(4)和隐函数的求导法则可以得到：

同理可得：

伴随方程的虚源为目标函数对模拟波场的偏导数，则

步骤(4)中背景速度与反射系数梯度计算公式如下：

目标函数对背景速度的梯度为：

目标函数对反射系数的梯度为：

通过方程(17)和方程(18)可得到背景速度和反射系数的梯度；

其中，Δp为正向扰动波场，p为正向背景波场；v为背景速度；m为反射系数；q为反向背景波场；和Δq为伴随扰动波场；

步骤(6)背景速度模型更新迭代方法如下：

a.求解方程得到正向背景波场p和正向扰动波场Δp，初始条件为：

x表示空间坐标的水平分量，z表示空间坐标的垂直分量；

b.求解伴随方程得到反向背景波场q，终值条件为：

c.通过方程(17)计算由背景波场决定的背景速度梯度

d.求解伴随方程得到反向背景波场q和伴随扰动波场Δq，终值条件为：

e.通过方程(17)计算由反射波场决定的背景速度梯度

f.将步骤c和e中的梯度加权求和得到总的背景速度梯度；

g.对总的背景速度梯度进行预条件处理；

h.采用抛物线拟合法计算迭代步长，对背景速度模型进行更新迭代。

2.根据权利要求1所述的回折波和反射波波动方程旅行时反演方法，其特征在于，步骤(5)反射系数反演迭代步骤如下：

只采用反射波场更新反射系数，目标函数变为：

反射系数反演迭代步骤如下：

a.求解方程得到正向背景波场p和正向扰动波场Δp，初始条件为：

b.求解伴随方程得到反向背景波场q，终值条件为：

c.通过方程(18)计算目标函数关于反射系数的梯度；

d.对反射系数梯度进行预条件处理；

e.采用抛物线拟合法计算迭代步长，对反射系数模型进行更新迭代。

Images