CN107229071A

CN107229071A - 一种地下构造反演成像方法

Info

Publication number: CN107229071A
Application number: CN201710376504.7A
Authority: CN
Inventors: 国运东; 黄建平; 崔超; 李振春; 李庆洋; 李闯
Original assignee: China University of Petroleum East China
Current assignee: China University of Petroleum East China
Priority date: 2017-05-25
Filing date: 2017-05-25
Publication date: 2017-10-03
Anticipated expiration: 2037-05-25
Also published as: CN107229071B

Abstract

本发明公开了一种地下构造反演成像方法，包括以下步骤：S1，输入背景速度场、观测数据、震源，获取初始反射系数模型；S2，对所述观测数据进行编码，获得超炮集；对震源进行编码，获得编码后的震源；S3，根据编码后的震源和初始反射系数模型，获得正演模拟数据；S4，获得正演模拟数据与超炮集之间的数据残差；S5，根据数据残差确定梯度，基于梯度进行反射系数模型迭代，获得最终的反射系数模型。本发明通过对观测数据和震源采用分段阶梯编码，充分加大多震源炮数据之间的差异性，减少临炮之间的干扰，在多震源炮数据成像时较为好的压制串扰噪音，提高了多震源炮数据成像反演时的成像精度。

Description

一种地下构造反演成像方法

技术领域

本发明涉及一种地下构造反演成像方法，属于油气物探工程领域。

背景技术

在当前的勘探领域中，地下数据采集技术尤为关键，这其中多震源数据采集方式开始得到关注，所谓多震源数据是多震源同时激发，波场在地下介质中传播，在地面接受得到的数据，其是地下介质对多震源的响应。多震源数据采集方式改善了传统采集方式照明不足的问题，也大大降低了三维数据采集成本。但是当前在对多震源数据进行处理时，由于不相关炮的波场进行互相关而产生的串扰噪声，由于不能有效压制这些串扰噪声，最终影响了成像质量。因此，在对这样的多震源数据做处理的时候，压制串扰噪音以提高其成像精度是非常必要的。

发明内容

本发明的目的在于，提供一种地下构造反演成像方法，它可以解决当前技术中存在的问题，有效压制多震源数据中之间的串扰噪声，提高反演成像精度。

为解决上述技术问题，本发明采用如下的技术方案：一种地下构造反演成像方法，包括以下步骤：

S1，输入背景速度场、观测数据和震源，获取初始反射系数模型；

S2，对所述观测数据进行编码，获得超炮集；对震源进行编码，获得编码后的震源；

S3，根据编码后的震源和初始反射系数模型，获得正演模拟数据；

S4，获得正演模拟数据与超炮集之间的数据残差；

S5，根据数据残差确定梯度，基于梯度进行反射系数模型迭代，获得最终的反射系数模型。

与现有技术相比，本发明对观测数据和震源采用分段阶梯编码，充分加大多震源炮数据之间的差异性，减少临炮之间的干扰，在多震源数据成像时较为好的压制串扰噪音，提高了多震源数据成像反演时的成像精度。此外，利用分段阶梯编码可以直接对多震源数据进行偏移，而不需要对每炮数据进行单独偏移成像再叠加，从而提高了反演的效率。

附图说明

附图是用来提供对本发明的进一步理解，并且构成说明书的一部分，与下面的具体实施方式一起用于解释本发明，但并不构成对本发明的限制。在附图中：

图1为本发明的编码示意图；

图2为本发明实施例的流程示意图；

图3为本发明实施例的背景速度场示意图；

图4为本发明实施例观测数据中的单炮观测记录示意图；

图5为本发明实施例中超炮集的记录示意图；

图6为本发明实施例第5次迭代得到的反射系数模型示意图；

图7为本发明实施例最终的反射系数模型示意图；

图8为常规方法得到最终的反射系数模型示意图；

图9为图7和图8的局部放大后的对比，其中图(a)和(c)为常规方法所得，图(b)(d)为本发明所用的方法所得；

图10为常规方法与本发明方法的迭代误差曲线对比；

下面结合附图和具体实施方式对本发明作进一步的说明。

具体实施方式

本发明的实施例1：一种地下构造反演成像方法，如图2所示，包括以下步骤：

S1，输入背景速度场、观测数据和震源。在实践中，根据需要，观测数据常为多震源所采集得到的多震源数据，例如，将本发明的方法应用到国际标准盐丘模型上，图3为所述盐丘模型的背景速度场；使用波恩正演得到单炮数据如图4，根据背景速度场、观测数据、震源获得初始反射系数模型。

S2，对所述观测数据进行编码，获得超炮集；对震源进行编码，获得编码后的震源。也就是对炮集利用分段阶梯平面波进行编码形成编码的超炮集D_obs，在本实施例中观测数据有 160炮，其中每段平面波含有12炮数据，多角度平面波编码经验公式来计算第m炮第n次迭代的编码激发时间延迟量，采用如下计算方法：

其中，tdely_m,n为第m炮第n次编码的激发时间延迟量，mod为取余函数，N_s是编码的总炮数，N_T是每个分段阶梯编码的炮数，int为取整函数，Δx为相邻震源在地面上的间隔距离，v为表层速度场，P(n)为第n次编码的射线参数，即图1的角P所在的虚线，其垂直方向(箭头c)为主轴波场传播方向，A(n)为第n次编码的分段平面波的次轴传播方向(即箭头d)，G代表每段编码的炮数，这里当G越大则主轴传播方向(图1箭头c)能量越弱，反之越强。

因为分段阶梯编码方式，加大了临炮之间的相干性，另外与常规方法相比较，其包含多个波场传播方向信息，因此其可以更好的压制串扰噪音，更精确的对复杂区域成像。所得编码结果如图1所示，图1中P为分段阶梯的主轴方向，其垂直方向(箭头c)为主轴波场传播方向，横坐标为距离，纵坐标为延迟时间，星号代表震源位置，其纵坐标值代表激发延迟时间，箭头d表示分段平面波的次轴传播方向，虚线表示本发明方法可以包含不同的波场传播方向。

将多震源数据最小二乘偏移的正算子用L表示，偏移算子用L^*表示，观测数据总共有 Ns炮，而我们将观测数据编码成M_ga个超炮集，当N_s能够被M_ga整除时，满足N_s＝M_gan_s， n_s就是每个超炮集中编码的炮数。不能整除的时候也就是最后最后一次编码的炮数就少了，但是对整体并无影响。例如是5炮作为分成两个超道集，第一个超道集就是3，第二个就是2，此时n_s＝3。通常可以涉及让N_s能够被M_ga整除，例如若只用一个超炮集进行计算，Ns＝n_s。

对单炮数据利用分段阶梯编码的公式进行编码，其中得到超炮集，如图5所示。

S3，根据编码后的震源和初始反射系数模型，获得正演模拟数据。也就是根据编码后的震源进行多炮编码的正演模拟得到正演模拟炮数据D_cal；

S4，获得正演模拟数据与超炮集之间的数据残差；反演的过程即为通过不断迭代最小化一个目标函数，这里的目标函数即数据残差，表示为

其中，k表示迭代次数，γ表示进行偏移的编码后的超炮道集的序数。δ^(k,γ)表示第γ个超炮道集在第k次迭代时的残差即为：

δ^(k,γ)＝L^(k,γ)m-d^(k,γ)， (3)

S5，根据数据残差确定梯度，基于梯度进行反射系数模型迭代，获得最终的反射系数模型。具体的说，对数据残差进行反传得到第k次的更新梯度值g^(k)，用线性搜索方法或者抛物拟合方法求取更新步长α，然后更新反射系数模型:

m^(k+1)＝m^(k)+α^(k)g^(k)， (4)

其中，α^(k)表示步长，g^(k)表示梯度，他们的表达式分别为：

当误差小于预设值时，所得的反射系数模型即为最终的反射系数模型。在对前述盐丘模型的反演中，经过反复迭代直至满足误差条件，其中第5次迭代的结果如图6所示，经过50 次迭代后输出最终的反射系数模型如图7所示。传统全波形反演第50次迭代的反演结果如图 8，可以看出图7中对盐丘下成像及复杂构造成像更好，可以较为精确的刻画复杂构造，减少串扰噪音的影响。为了更好的比较对图7和图8的，将它们的典型区域进行放大比较，放大都的局部对比如图9所示。其中图9(a)和图9(c)为传统方法所得，图9(b)和图9(d)为本申请方法所得。

如图9(b)相对于图9(a)中箭头所指的位置，其在高陡构造处噪音更少，反射轴更加清楚并且对细小的构造也刻画出来，如图9(d)相对于图9(c)中箭头所指的位置，对于盐下的断层反射轴可以更清楚的刻画出来。

相比于传统的方法，采用本发明方法的反演误差曲线(图10)整体的收敛更快，例如误差降低到0.1的时候，本文方法迭代36次，传统方法需要迭代46次，计算效率更高。

通过前述的方法，本发明对观测数据和震源采用分段阶梯编码，充分加大多震源炮数据之间的差异性，减少临炮之间的干扰，在多震源炮数据成像时较为好的压制串扰噪音，提高了多震源数据成像反演时的成像精度。此外，利用分段阶梯编码可以直接对多震源数据进行偏移，而不需要对每炮数据进行单独偏移成像再叠加，从而提高了反演的效率。

Claims

1.一种地下构造反演成像方法，其特征在于，包括以下步骤：

S4，获得正演模拟数据与超炮集之间的数据残差；

2.根据权利要求1所述的一种地下构造反演成像方法，其特征在于，所述步骤S2中的编码，具体包括：采用分段平面波编码。

3.根据权利要求2所述的一种地下构造反演成像方法，其特征在于，所述采用分段阶梯平面波编码，具体包括：采用分段平面波编码计算第m炮第n次编码的激发时间延迟量，所述公式如下:

<mrow> <msub> <mi>tdely</mi> <mrow> <mi>m</mi> <mo>,</mo> <mi>n</mi> </mrow> </msub> <mo>=</mo> <mfrac> <mrow> <mi>&Delta;</mi> <mi>x</mi> </mrow> <mi>v</mi> </mfrac> <mo>{</mo> <mi>tan</mi> <mi> </mi> <mi>P</mi> <mrow> <mo>(</mo> <mi>n</mi> <mo>)</mo> </mrow> <mo>*</mo> <mo>&lsqb;</mo> <mi>m</mi> <mo>-</mo> <mi>mod</mi> <mrow> <mo>(</mo> <mi>m</mi> <mo>,</mo> <mi>G</mi> <mo>)</mo> </mrow> <mo>+</mo> <mfrac> <mi>G</mi> <mn>2</mn> </mfrac> <mo>&rsqb;</mo> <mo>-</mo> <mi>tan</mi> <mi> </mi> <mi>A</mi> <mrow> <mo>(</mo> <mi>n</mi> <mo>)</mo> </mrow> <mo>*</mo> <mo>&lsqb;</mo> <mi>mod</mi> <mrow> <mo>(</mo> <mi>m</mi> <mo>,</mo> <mi>G</mi> <mo>)</mo> </mrow> <mo>-</mo> <mfrac> <mi>G</mi> <mn>2</mn> </mfrac> <mo>&rsqb;</mo> <mo>,</mo> </mrow>

其中，tdely_m,n为第m炮第n次编码的激发时间延迟量，mod为取余函数，N_s是编码的总炮数，N_T是每个分段阶梯编码的炮数，int为取整函数，Δx为相邻震源在地面上的间隔距离，v为表层速度场，P(n)为第n次编码的射线参数，A(n)为第n次编码的分段平面的次轴方向，这里当G代表每段编码的炮数。

4.根据权利要求1所述的一种地下构造反演成像方法，其特征在于，所述步骤S5中的根据数据残差确定梯度，具体包括：对数据残差进行反传计算得到梯度。

5.根据权利要求1所述的一种地下构造反演成像方法，其特征在于，所述骤S5中的基于梯度进行反射系数模型迭代，具体包括：用线性搜索方法或者抛物拟合方法求取更新步长α，从而更新反射系数模型；m^(k+1)＝m^(k)+α^(k)g^(k)，其中α^(k)表示步长，g^(k)表示梯度，m^(k)表示反射系数模型。