CN105954803A - 叠后地震反演方法及装置 - Google Patents

Abstract

本发明提供一种叠后地震反演方法及装置，其中方法包括：根据波阻抗的初始模型与实际地震数据得到合成地震数据；根据实际地震数据与合成地震数据优化波阻抗的初始模型，得到第一最优波阻抗；对第一最优波阻抗进行高通滤波，得到波阻抗的高频反演结果；根据实际地震数据得到实际地震数据包络，并根据合成地震数据得到合成地震数据包络；根据实际地震数据包络与合成地震数据包络优化波阻抗的初始模型，得到第二最优波阻抗；对第二最优波阻抗进行低通滤波，得到波阻抗的低频反演结果；根据波阻抗的高频反演结果和波阻抗的低频反演结果得到波阻抗的最终反演结果。通过实施本发明，提高了反演精度，在储层预测中实现定量化描述。

Description

叠后地震反演方法及装置

技术领域

本发明涉及地球物理勘探技术领域，具体地，涉及一种叠后地震反演方法及装置。

背景技术

地震勘探是在地表激发人工震源，由震源引起的震动以地震波的形式向地下传播，并在一定的条件下向上反射传回地表，然后由检波器数据反射回地震波，得到地震记录(也叫地震数据)。地震数据是经由地下介质而得到的，必然受到地下介质物理性质(如岩性、密度、孔隙度等)的影响，其中不仅包含有反射波，还包含有噪声、直达波、面波等。要利用地震数据进行地下储层的油气预测，首先需要进行数据处理，去除各种噪声、直达波和面波。地震数据中的面波能量非常强，频率仅在十几赫兹以内，通常在频率域对其进行去除。处理后的地震数据，其频率一般从十几赫兹到几十赫兹，十几赫兹以下的频率量值很小。

地震反演是储层地球物理勘探领域的重要课题，叠后地震反演是现今储层预测最常用的一项技术，通常是把地震数据转换成波阻抗剖面以进行储层预测。由于地震数据十几赫兹以下的数据量值很小，在叠后地震反演中通常仅能得到和地震频带匹配的波阻抗数据，无法得到地下波阻抗参数的背景值(也即为低频值)，从而难以进行储层的定量化预测。在实际反演中，通常的做法是利用地震工区内的测井数据进行井震标定，在地震层位的约束下，利用经过标定的测井数据进行数学插值以建立地下介质物性参数的低频模型(也叫初始模型)，再利用初始模型进行反演以得到相对精确的地下介质物性参数。但无论插值方式多先进，都难以完全反应地下构造在空间上的真实变化，尤其在测井数量稀少的勘探阶段更是如此。因此利用测井数据插值得到初始模型进行反演，难以在储层预测中实现定量化描述。

发明内容

本发明提供一种叠后地震反演方法及装置，以解决现有技术难以在储层预测中实现定量化描述的问题。

为了实现上述目的，本发明实施例提供一种叠后地震反演方法，包括：根据波阻抗的初始模型与实际地震数据得到合成地震数据；根据实际地震数据与合成地震数据优化波阻抗的初始模型，得到第一最优波阻抗；对第一最优波阻抗进行高通滤波，得到波阻抗的高频反演结果；根据实际地震数据得到实际地震数据包络，并根据合成地震数据得到合成地震数据包络；根据实际地震数据包络与合成地震数据包络优化波阻抗的初始模型，得到第二最优波阻抗；对第二最优波阻抗进行低通滤波，得到波阻抗的低频反演结果；根据波阻抗的高频反演结果和波阻抗的低频反演结果得到波阻抗的最终反演结果。

在其中一种实施例中，根据波阻抗的初始模型与实际地震数据得到合成地震数据，包括：在实际地震数据目的层段提取对应的子波；根据初始模型得到反射系数；根据反射系数与子波得到合成地震数据。

在其中一种实施例中，按如下公式根据反射系数与子波得到合成地震数据：d＝R*w；其中，矢量d为合成地震数据；矢量R为反射系数，其中，反射系数R的第i个采样点数据R_i(i＝1,2,3,...,n-1)如下：最后一个采样点数据R_n＝R_n-1；Z_i为波阻抗的初始模型第i个采样点数据；Z_i+1为波阻抗的初始模型第i+1个采样点数据；矢量w为实际地震数据目的层段对应的子波。

在其中一种实施例中，根据实际地震数据与合成地震数据优化波阻抗的初始模型，得到第一最优波阻抗，包括：计算合成地震数据与实际地震数据之间的误差；当误差处于预设范围内时，波阻抗的初始模型为第一最优波阻抗；当误差不处于预设范围内时，设定波阻抗的初始模型的修改次数为j，j＝1，并按照如下方式进行迭代处理：采用非线性全局寻优算法第j次修改波阻抗的初始模型，得到波阻抗的第j次修改模型；根据波阻抗的第j次修改模型与实际地震数据得到第j次修改合成地震数据；计算第j次修改合成地震数据与实际地震数据之间的误差；当误差处于预设范围内时，波阻抗的第j次修改模型为第一最优波阻抗；当误差不处于预设范围内时，将迭代处理中的j替代为j+1。

在其中一种实施例中，按照如下方式计算合成地震数据与实际地震数据之间的误差：其中，J_d为合成地震数据与实际地震数据之间的误差；

矢量d₀为实际地震数据，如下；d₀＝[d₀₁ d₀₂ d₀₃ … d_0n]^T；其中，d_0i(i＝1,2,3...n)为实际地震数据d₀第i个采样点数据；矢量d为合成地震数据，如下：d＝[d₁ d₂ d₃ …d_n]^T；其中，d_i(i＝1,2,3...n)为合成地震数据d第i个采样点数据。

在其中一种实施例中，按照如下方式根据实际地震数据得到实际地震数据包络：e₀＝[e₀₁ e₀₂ e₀₃ … e_0n]^T，其中，矢量e₀为实际地震数据包络，e_0i(i＝1,2,3...n)为实际地震数据包络e₀第i个采样点数据；d_0i的希尔伯特变换如下：按照如下方式根据合成地震数据得到合成地震数据包络：e＝[e₁ e₂ e₃ … e_n]^T，其中，矢量e为合成地震数据包络，e_i(i＝1,2,3...n)为合成地震数据包络e第i个采样点数据；d_i的希尔伯特变换如下：

在其中一种实施例中，根据实际地震数据包络与合成地震数据包络优化波阻抗的初始模型，得到第二最优波阻抗，包括：计算合成地震数据包络与实际地震数据包络之间的误差；当误差处于预设范围内时，波阻抗的初始模型为第二最优波阻抗；当误差不处于预设范围内时，设定波阻抗的初始模型的修改次数为k，k＝1，并按照如下方式进行迭代处理：采用非线性全局寻优算法第k次修改波阻抗的初始模型，得到波阻抗的第k次修改模型；根据波阻抗的第k次修改模型与实际地震数据得到第k次修改合成地震数据包络；计算第k次修改合成地震数据包络与实际地震数据包络之间的误差；当误差处于预设范围内时，波阻抗的第k次修改模型为第二最优波阻抗；当误差不处于预设范围内时，将迭代处理中的k替代为k+1。

在其中一种实施例中，按照如下方式计算合成地震数据包络与实际地震数据包络之间的误差：其中，J_e为合成地震数据包络与实际地震数据包络之间的误差。

在其中一种实施例中，按照如下方式根据波阻抗的高频反演结果和波阻抗的低频反演结果得到波阻抗的最终反演结果：Z＝Z_1h+Z_2l；其中，矢量Z为波阻抗的最终反演结果，矢量Z_1h为波阻抗的高频反演结果，矢量Z_2l为波阻抗的低频反演结果。

本发明实施例还提供一种叠后地震反演装置，包括：地震数据合成模块，用于根据波阻抗的初始模型与实际地震数据得到合成地震数据；高频反演模块，用于根据实际地震数据与合成地震数据优化波阻抗的初始模型，得到第一最优波阻抗；对第一最优波阻抗进行高通滤波，得到波阻抗的高频反演结果；低频反演模块，用于根据实际地震数据得到实际地震数据包络，并根据合成地震数据得到合成地震数据包络；根据实际地震数据包络与合成地震数据包络优化波阻抗的初始模型，得到第二最优波阻抗；对第二最优波阻抗进行低通滤波，得到波阻抗的低频反演结果；叠后地震反演模块，用于根据波阻抗的高频反演结果和波阻抗的低频反演结果得到波阻抗的最终反演结果。

在其中一种实施例中，地震数据合成模块具体用于：在实际地震数据目的层段提取对应的子波；根据初始模型得到反射系数；根据反射系数与子波得到合成地震数据。

在其中一种实施例中，高频反演模块具体用于按如下方式得到第一最优波阻抗：计算合成地震数据与实际地震数据之间的误差；当误差处于预设范围内时，波阻抗的初始模型为第一最优波阻抗；当误差不处于预设范围内时，设定波阻抗的初始模型的修改次数为j，j＝1，并按照如下方式进行迭代处理：采用非线性全局寻优算法第j次修改波阻抗的初始模型，得到波阻抗的第j次修改模型；根据波阻抗的第j次修改模型与实际地震数据得到第j次修改合成地震数据；计算第j次修改合成地震数据与实际地震数据之间的误差；当误差处于预设范围内时，波阻抗的第j次修改模型为第一最优波阻抗；当误差不处于预设范围内时，将迭代处理中的j替代为j+1。

在其中一种实施例中，低频反演模块具体用于按如下方式得到第二最优波阻抗：根据实际地震数据得到实际地震数据包络，并根据合成地震数据得到合成地震数据包络；计算合成地震数据包络与实际地震数据包络之间的误差；当误差处于预设范围内时，波阻抗的初始模型为第二最优波阻抗；当误差不处于预设范围内时，设定波阻抗的初始模型的修改次数为k，k＝1，并按照如下方式进行迭代处理：采用非线性全局寻优算法第k次修改波阻抗的初始模型，得到波阻抗的第k次修改模型；根据波阻抗的第k次修改模型与实际地震数据得到第k次修改合成地震数据包络；计算第k次修改合成地震数据包络与实际地震数据包络之间的误差；当误差处于预设范围内时，波阻抗的第k次修改模型为第二最优波阻抗；当误差不处于预设范围内时，将迭代处理中的k替代为k+1。

在其中一种实施例中，叠后地震反演模块具体用于按照如下方式根据波阻抗的高频反演结果和波阻抗的低频反演结果得到波阻抗的最终反演结果：Z＝Z_1h+Z_2l；其中，矢量Z为波阻抗的最终反演结果，矢量Z_1h为波阻抗的高频反演结果，矢量Z_2l为波阻抗的低频反演结果。

在本发明实施例提供的技术方案中，首先基于初始模型进行反演，对反演结果高通滤波得到和地震频带匹配的高频反演结果，然后计算地震数据的包络，再次利用初始模型进行反演，对反演结果低通滤波得到低频反演结果，合并高频反演结果和低频反演结果为最终反演结果，突出了地震数据十几赫兹以下的低频量级，提高了反演精度，在储层预测中实现定量化描述。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1是本发明实施例中叠后地震反演方法的流程图；

图2是本发明实施例中步骤101的具体流程图；

图3是本发明实施例中步骤102的具体流程图；

图4是本发明实施例中步骤105的具体流程图；

图5是本发明实施例中所用检验井的实际波阻抗数据示意图；

图6是本发明实施例中所用检验井的实际波阻抗数据的低频部分和波阻抗初始模型的对比示意图；

图7是本发明实施例中所用检验井的实际波阻抗数据的高频部分示意图；

图8是本发明实施例中所用检验井坐标对应处的一道实际地震数据示意图；

图9是本发明实施例中所用检验井坐标对应处的一道实际地震数据的频率谱；

图10是本发明实施例中实际地震数据对应的子波示意图；

图11是本发明实施例中第一最优波阻抗数据和所用检验井的实际波阻抗数据的对比示意图；

图12是本发明实施例中的第一最优波阻抗数据的高频部分和实际波阻抗数据的高频部分的对比示意图；

图13是本发明实施例中计算得出的实际地震数据的包络和所用检验井坐标对应处的一道实际地震数据的对比示意图；

图14是本发明实施例中计算得出的实际地震数据的包络的频率谱和所用检验井坐标对应处的一道实际地震数据的频率谱的对比示意图；

图15是本发明实施例中第二最优波阻抗数据和所用检验井的实际波阻抗数据的对比示意图；

图16是本发明实施例中第二最优波阻抗数据的低频部分和所用检验井的实际波阻抗数据的低频部分的对比示意图；

图17是本发明实施例中波阻抗的最终反演结果和所用检验井的实际波阻抗数据的对比示意图；

图18为本发明实施例中叠后地震反演装置的结构示意图。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

本发明提供了一种叠后地震反演方法，以解决现有技术中利用测井数据插值得到初始模型进行反演，难以在储层预测中实现定量化描述的问题。

图1是本发明实施例中叠后地震反演方法的流程图。如图1所示，叠后地震反演方法包括：

步骤101：根据波阻抗的初始模型与实际地震数据得到合成地震数据。

步骤102：根据实际地震数据与合成地震数据优化波阻抗的初始模型，得到第一最优波阻抗。

步骤103：对第一最优波阻抗进行高通滤波，得到波阻抗的高频反演结果。

步骤104：根据实际地震数据得到实际地震数据包络，并根据合成地震数据得到合成地震数据包络。

步骤105：根据实际地震数据包络与合成地震数据包络优化波阻抗的初始模型，得到第二最优波阻抗。

步骤106：对第二最优波阻抗进行低通滤波，得到波阻抗的低频反演结果。

步骤107：根据波阻抗的高频反演结果和波阻抗的低频反演结果得到波阻抗的最终反演结果。

具体实施时，先建立波阻抗的初始模型。建立波阻抗的初始模型有多种方式，比如可以：采集经过保幅处理的实际地震数据和地震工区内相应的测井数据，根据测井数据进行井震标定，对经过标定的测井数据进行插值。

下面结合具体的实施例详细介绍本发明实施例提供的一种叠后地震反演方法，以四川潼南须家河组一口检验井及其对应的一道地震数据为例进行说明。

通过如下方式建立波阻抗的初始模型：

采集经过保幅处理的地震数据、地震工区内相应的测井数据，根据已知先验地质层位与钻井层位，进行井震标定，把深度域测井的声波时差曲线标定为时间域，利用经过标定的测井数据进行插值以建立波阻抗的初始模型。

波阻抗的初始模型可以用矢量Z₀来表示：

Z₀＝[Z₁ Z₂ Z₃ … Z_n]^T，

其中，Z_i(i＝1,2,3,...,n)为第i个采样点数据，矢量Z₀共有n个元素。

图5是本发明实施例中所用检验井的实际波阻抗数据示意图。如图5所示，图5中实线为本发明中所用检验井的实际波阻抗数据，数据共有192个采样点。图6是本发明实施例中所用检验井的实际波阻抗数据的低频部分和波阻抗初始模型的对比示意图。如图6所示，图6中实线为图5中实际波阻抗数据通过0Hz-10Hz的低通滤波后得到的低频部分，虚线为利用地震工区中其它测井数据通过插值建立的0Hz-10Hz的波阻抗初始模型Z₀。图7是本发明实施例中所用检验井的实际波阻抗数据的高频部分示意图。如图7所示，图7中实线为图5中实际波阻抗数据减去图6中低频部分(实线)后得到的高频部分。对比图6可以看出，通过插值建立的波阻抗初始模型和实际波阻抗数据的低频部分存在很大的差距。

得到合成地震数据有多种方式，比如可以按如下方式得到合成地震数据：

图2是本发明实施例中步骤101的具体流程图；如图2所示，步骤101具体包括：

步骤201：在实际地震数据目的层段提取对应的子波。

步骤202：根据初始模型得到反射系数。

步骤203：根据反射系数与子波得到合成地震数据。

具体实施时，先执行步骤201；实际地震数据可以用矢量d₀来表示，如下：

d₀＝[d₀₁ d₀₂ d₀₃ … d_0n]^T，

其中，d_0i(i＝1,2,3...n)为实际地震数据d₀第i个采样点数据，矢量d₀共有n个元素。

对应的子波可以用矢量w来表示。

图8是本发明实施例中所用检验井坐标对应处的一道实际地震数据示意图。如图8所示，图8中的实线为本发明中所用检验井坐标对应处的一道实际地震数据，数据共有192个采样点；图9是本发明实施例中所用检验井坐标对应处的一道实际地震数据的频率谱。如图9所示，频率谱缺失了10Hz以下的频率。图10是本发明实施例中实际地震数据对应的子波示意图。

然后，执行步骤202，得到反射系数可以有多种方式，比如可以按如下方式得到反射系数：

由初始模型计算得到反射系数，反射系数可以用矢量R来表示，其中，反射系数R的第i个采样点数据R_i(i＝1,2,3,...,n-1)表述如下：

$R_i=\frac{Z_{i+1}-Z_i}{Z_{i+1}+Z_i},$

最后一个采样点数据R_n＝R_n-1；

其中，Z_i为波阻抗的初始模型第i个采样点数据；Z_i+1为波阻抗的初始模型第i+1个采样点数据；R_n为反射系数第n个采样点数据。

最后，执行步骤203得到合成数据；具体实施时，可按如下公式得到合成地震数据：

d＝R*w；

上述公式中的符号“*”为卷积符号，合成地震数据为反射系数与子波的卷积。

其中，矢量d为合成地震数据，如下：

d＝[d₁ d₂ d₃ … d_n]^T；

d_i(i＝1,2,3...n)为合成地震数据d第i个采样点数据，矢量d共有n个元素。

在得到合成地震数据之后，执行步骤102；图3是本发明实施例中步骤102的具体流程图；如图3所示，步骤102具体包括：

步骤301：计算合成地震数据与实际地震数据之间的误差。

步骤302：判断误差是否处于预设范围内。

步骤303：当误差处于预设范围内时，波阻抗的初始模型为第一最优波阻抗。

步骤304：当误差不处于预设范围内时，设定波阻抗的初始模型的修改次数为j。

步骤305：设定j＝1。

步骤306：采用非线性全局寻优算法第j次修改波阻抗的初始模型，得到波阻抗的第j次修改模型。

步骤307：根据波阻抗的第j次修改模型与实际地震数据得到第j次修改合成地震数据。

步骤308：计算第j次修改合成地震数据与实际地震数据之间的误差。

步骤309：判断误差是否处于预设范围内。

步骤310：当误差处于预设范围内时，波阻抗的第j次修改模型为第一最优波阻抗。

步骤311：当误差不处于预设范围内时，j＝j+1，返回步骤306。

具体实施时，先执行步骤301，可以按照如下方式得到合成地震数据与实际地震数据之间的误差：

$J_d=\mathrm{\Σ}{(d_0-d)}^2=\underset{i=1}{\overset{n}{\Σ}}{(d_{0i}-d_i)}^2,$

其中，J_d为合成地震数据与实际地震数据之间的误差。

实施例中，第一最优波阻抗可以用矢量Z₁来表示，它具有与波阻抗的初始模型Z₀相同的元素个数。图11是本发明实施例中第一最优波阻抗数据和所用检验井的实际波阻抗数据的对比示意图。如图11所示，图11中虚线是本实施例中经过反演得到的第一最优波阻抗Z₁的数据，实线即为图5中本发明实施例中所用检验井的实际波阻抗数据；对比可以看出，第一最优波阻抗和实际波阻抗数据在某些部分还是存在较大误差。

在得到第一最优波阻抗之后，执行步骤103；实施例中，波阻抗的高频反演结果(高频部分)可以用矢量Z_1h来表示，它具有与第一最优波阻抗Z₁相同的元素个数。图12是本发明实施例中的第一最优波阻抗数据的高频部分和实际波阻抗数据的高频部分的对比示意图。如图12所示，图12中虚线为第一最优波阻抗Z₁进行高通滤波后得到的高频反演结果Z_1h，实线即为图7中本发明中所用检验井的实际波阻抗数据的高频部分；对比可以看出，第一最优波阻抗的高频部分和实际波阻抗数据的高频部分误差很小。

在得到波阻抗的高频反演结果之后，执行步骤104；具体实施时，按照如下方式根据实际地震数据得到实际地震数据包络：

e₀＝[e₀₁ e₀₂ e₀₃ … e_0n]^T，

$e_{0i}=\sqrt{d_{0i}^2+{\hat{d}}_{0i}^2};$

其中，矢量e₀为实际地震数据包络，e_0i(i＝1,2,3...n)为实际地震数据包络e₀第i个采样点数据；

d_0i的希尔伯特变换如下：

${\hat{d}}_{0i}=d_{0i}*h_i,$

$h_i=\frac{1-{(-1)}^i}{i\mathrm{\π}}.$

具体实施时，按照如下方式根据合成地震数据得到合成地震数据包络：

e＝[e₁ e₂ e₃ … e_n]^T，

$e_i=\sqrt{d_i^2+{\hat{d}}_i^2};$

其中，矢量e为合成地震数据包络，e_i(i＝1,2,3...n)为合成地震数据包络e第i个采样点数据；

d_i的希尔伯特变换如下：

${\hat{d}}_i=d_i*h_i.$

图13是本发明实施例中计算得出的实际地震数据的包络和所用检验井坐标对应处的一道实际地震数据的对比示意图。如图13所示，图13中虚线是本发明中计算得到的实际地震数据的包络e₀，实线即为图8中检验井坐标对应处的一道实际地震数据；对比可以看出，包络可以很好地反应实际地震数据反射特征及趋势。图14是本发明实施例中计算得出的实际地震数据的包络的频率谱和所用检验井坐标对应处的一道实际地震数据的频率谱的对比示意图。如图14所示，图14中虚线为是本发明中计算得到的实际地震数据的包络e₀的频率谱，实线即为图9中实际地震数据的频率谱；对比可以看出，包络的频率谱很好地弥补了实际地震数据频率谱中的低频部分。

在得到实际地震数据包络与合成地震数据包络之后，执行步骤105；图4是本发明实施例中所用检验井的实际波阻抗数据示意图；如图4所示，步骤105具体包括：

步骤401：计算合成地震数据包络与实际地震数据包络之间的误差。

步骤402：判断误差是否处于预设范围内。

步骤403：当误差处于预设范围内时，波阻抗的初始模型为第二最优波阻抗。

步骤404：当误差不处于预设范围内时，设定波阻抗的初始模型的修改次数为k。

步骤405：设定k＝1。

步骤406：采用非线性全局寻优算法第k次修改波阻抗的初始模型，得到波阻抗的第k次修改模型。

步骤407：根据波阻抗的第k次修改模型与实际地震数据得到第k次修改合成地震数据包络。

步骤408：计算第k次修改合成地震数据包络与实际地震数据包络之间的误差。

步骤409：判断误差是否处于预设范围内。

步骤410：当误差处于预设范围内时，波阻抗的第k次修改模型为第二最优波阻抗。

步骤411：当误差不处于预设范围内时，k＝k+1，返回步骤406。

具体实施时，先执行步骤401，可以按照如下方式得到合成地震数据包络与实际地震数据包络之间的误差：

$J_e=\mathrm{\Σ}{(e_0-e)}^2=\underset{i=1}{\overset{n}{\Σ}}{(e_{0i}-e_i)}^2,$

其中，J_e为合成地震数据包络与实际地震数据包络之间的误差。

实施例中，第二最优波阻抗可以用矢量Z₂来表示，它具有与波阻抗的初始模型Z₀相同的元素个数。图15是本发明实施例中第二最优波阻抗数据和所用检验井的实际波阻抗数据的对比示意图。如图15所示，图15中虚线是本发明中经过包络反演得到的第二最优波阻抗Z₂的数据，实线即为图5中所用检验井的实际波阻抗数据；对比可以看出，第二最优波阻抗和实际波阻抗数据在某些部分还是存在较大误差。

在得到第二最优波阻抗之后，执行步骤106；实施例中，波阻抗的低频反演结果(低频部分)可以用矢量Z_2l来表示，它具有与第二最优波阻抗Z₂相同的元素个数。图16是本发明实施例中第二最优波阻抗数据的低频部分和所用检验井的实际波阻抗数据的低频部分的对比示意图。如图16所示，图16中虚线为第二最优波阻抗Z₂进行低通滤波后得到的低频反演结果Z_2l，实线即为图6中所用检验井的实际波阻抗数据的低频部分(实线)；对比可以看出，第二最优波阻抗的低频部分和实际波阻抗数据的低频部分误差很小。

在得到波阻抗的高频反演结果之后，执行步骤107；其中，波阻抗的最终反演结果可以用矢量Z来表示：

Z＝Z_1h+Z_2l。

图17是本发明实施例中波阻抗的最终反演结果和所用检验井的实际波阻抗数据的对比示意图。如图17所示，图17中虚线为波阻抗的最终反演结果Z，实线即为图5中所用检验井的实际波阻抗数据；对比可以看出，综合常规反演高频部分和包络反演低频部分的最终反演结果很大程度上还原了真实的波阻抗。

基于同一发明构思，本发明实施例中还提供了一种叠后地震反演装置，由于该装置解决问题的原理与叠后地震反演方法相似，因此该装置的实施可以参见方法的实施，重复之处不再赘述。

图18是本发明实施例中叠后地震反演装置的结构示意图，如图18所示，装置中可以包括：

地震数据合成模块181，用于根据波阻抗的初始模型与实际地震数据得到合成地震数据；

高频反演模块182，用于根据实际地震数据与合成地震数据优化波阻抗的初始模型，得到第一最优波阻抗；对第一最优波阻抗进行高通滤波，得到波阻抗的高频反演结果；

低频反演模块183，用于根据实际地震数据得到实际地震数据包络，并根据合成地震数据得到合成地震数据包络；根据实际地震数据包络与合成地震数据包络优化波阻抗的初始模型，得到第二最优波阻抗；对第二最优波阻抗进行低通滤波，得到波阻抗的低频反演结果；

叠后地震反演模块184，用于根据波阻抗的高频反演结果和波阻抗的低频反演结果得到波阻抗的最终反演结果。

在本发明的一个具体实施例中，地震数据合成模块具体用于：在实际地震数据目的层段提取对应的子波；根据初始模型得到反射系数；根据所述反射系数与所述子波得到合成地震数据。

在本发明的一个具体实施例中，高频反演模块具体用于按如下方式得到第一最优波阻抗：计算合成地震数据与实际地震数据之间的误差；当误差处于预设范围内时，波阻抗的初始模型为第一最优波阻抗；当误差不处于预设范围内时，设定波阻抗的初始模型的修改次数为j，j＝1，并按照如下方式进行迭代处理：采用非线性全局寻优算法第j次修改波阻抗的初始模型，得到波阻抗的第j次修改模型；根据波阻抗的第j次修改模型与实际地震数据得到第j次修改合成地震数据；计算第j次修改合成地震数据与实际地震数据之间的误差；当误差处于预设范围内时，波阻抗的第j次修改模型为第一最优波阻抗；当误差不处于预设范围内时，将迭代处理中的j替代为j+1。

在本发明的一个具体实施例中，低频反演模块具体用于按如下方式得到第二最优波阻抗：

根据实际地震数据得到实际地震数据包络，并根据合成地震数据得到合成地震数据包络；

计算合成地震数据包络与实际地震数据包络之间的误差；

当误差处于预设范围内时，波阻抗的初始模型为第二最优波阻抗；

当误差不处于预设范围内时，设定波阻抗的初始模型的修改次数为k，k＝1，并按照如下方式进行迭代处理：

采用非线性全局寻优算法第k次修改波阻抗的初始模型，得到波阻抗的第k次修改模型；

根据波阻抗的第k次修改模型与实际地震数据得到第k次修改合成地震数据包络；

计算第k次修改合成地震数据包络与实际地震数据包络之间的误差；

当误差处于预设范围内时，波阻抗的第k次修改模型为第二最优波阻抗；

当误差不处于预设范围内时，将迭代处理中的k替代为k+1。

在本发明的一个具体实施例中，叠后地震反演模块具体用于按照如下方式根据波阻抗的高频反演结果和波阻抗的低频反演结果得到波阻抗的最终反演结果：Z＝Z_1h+Z_2l；其中，矢量Z为波阻抗的最终反演结果，矢量Z_1h为波阻抗的高频反演结果，矢量Z_2l为波阻抗的低频反演结果。

由上述实施例可见，在本发明实施例提供的技术方案中，首先基于初始模型进行反演，对反演结果高通滤波得到和地震频带匹配的高频反演结果，然后计算地震数据的包络，再次利用初始模型进行反演，对反演结果低通滤波得到低频反演结果，合并高频反演结果和低频反演结果为最终反演结果，突出了地震数据十几赫兹以下的低频量级，解决了因地震数据低频量级小导致反演结果难以定量化的问题，很大程度上提高了反演精度，在地下复杂储层的预测中有明显的效果。

本领域内的技术人员应明白，本发明的实施例可提供为方法、系统、或计算机程序产品。因此，本发明可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且，本发明可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等)上实施的计算机程序产品的形式。

本发明是参照根据本发明实施例的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器，使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。

这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中，使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品，该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。

这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上，使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理，从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。

以上所述的具体实施例，对本发明的目的、技术方案和有益效果进行了进一步详细说明，所应理解的是，以上所述仅为本发明的具体实施例而已，并不用于限定本发明的保护范围，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (14)

1.一种叠后地震反演方法，其特征在于，包括：

根据波阻抗的初始模型与实际地震数据得到合成地震数据；

根据实际地震数据与合成地震数据优化波阻抗的初始模型，得到第一最优波阻抗；

对第一最优波阻抗进行高通滤波，得到波阻抗的高频反演结果；

根据实际地震数据得到实际地震数据包络，并根据合成地震数据得到合成地震数据包络；

根据实际地震数据包络与合成地震数据包络优化波阻抗的初始模型，得到第二最优波阻抗；

对第二最优波阻抗进行低通滤波，得到波阻抗的低频反演结果；

根据波阻抗的高频反演结果和波阻抗的低频反演结果得到波阻抗的最终反演结果。

2.根据权利要求1所述的叠后地震反演方法，其特征在于，所述根据波阻抗的初始模型与实际地震数据得到合成地震数据，包括：

在实际地震数据目的层段提取对应的子波；

根据初始模型得到反射系数；

根据所述反射系数与所述子波得到合成地震数据。

3.根据权利要求2所述的叠后地震反演方法，其特征在于，按如下公式根据所述反射系数与所述子波得到合成地震数据：

d＝R*w；

其中，矢量d为合成地震数据；

矢量R为反射系数，其中，反射系数R的第i个采样点数据R_i(i＝1,2,3,...,n-1)如下：

$R_i=\frac{Z_{i+1}-Z_i}{Z_{i+1}+Z_i};$

最后一个采样点数据R_n＝R_n-1；

Z_i为波阻抗的初始模型第i个采样点数据；Z_i+1为波阻抗的初始模型第i+1个采样点数据；

矢量w为实际地震数据目的层段对应的子波。

4.根据权利要求1所述的叠后地震反演方法，其特征在于，所述根据实际地震数据与合成地震数据优化波阻抗的初始模型，得到第一最优波阻抗，包括：

计算合成地震数据与实际地震数据之间的误差；

当误差处于预设范围内时，波阻抗的初始模型为第一最优波阻抗；

当误差不处于预设范围内时，设定波阻抗的初始模型的修改次数为j，j＝1，并按照如下方式进行迭代处理：

采用非线性全局寻优算法第j次修改波阻抗的初始模型，得到波阻抗的第j次修改模型；

根据波阻抗的第j次修改模型与实际地震数据得到第j次修改合成地震数据；

计算第j次修改合成地震数据与实际地震数据之间的误差；

当误差处于预设范围内时，波阻抗的第j次修改模型为第一最优波阻抗；

当误差不处于预设范围内时，将迭代处理中的j替代为j+1。

5.根据权利要求4所述的叠后地震反演方法，其特征在于，按照如下方式计算合成地震数据与实际地震数据之间的误差：

$J_d=\mathrm{\Σ}{(d_0-d)}^2=\underset{i=1}{\overset{n}{\Σ}}{(d_{0i}-d_i)}^2;$

其中，J_d为合成地震数据与实际地震数据之间的误差；

矢量d₀为实际地震数据，如下；

d₀＝[d₀₁ d₀₂ d₀₃ … d_0n]^T；

其中，d_0i(i＝1,2,3...n)为实际地震数据d₀第i个采样点数据；

矢量d为合成地震数据，如下：

d＝[d₁ d₂ d₃ … d_n]^T；

其中，d_i(i＝1,2,3...n)为合成地震数据d第i个采样点数据。

6.根据权利要求5所述的叠后地震反演方法，其特征在于，按照如下方式根据实际地震数据得到实际地震数据包络：

e₀＝[e₀₁ e₀₂ e₀₃ … e_0n]^T，

$e_{0i}=\sqrt{d_{0i}^2+{\hat{d}}_{0i}^2};$

其中，矢量e₀为实际地震数据包络，e_0i(i＝1,2,3...n)为实际地震数据包络e₀第i个采样点数据；

d_0i的希尔伯特变换如下：

${\hat{d}}_{0i}=d_{0i}*h_i,$

$h_i=\frac{1-{(-1)}^i}{i\mathrm{\π}};$

按照如下方式根据合成地震数据得到合成地震数据包络：

e＝[e₁ e₂ e₃ … e_n]^T，

$e_i=\sqrt{d_i^2+{\hat{d}}_i^2};$

其中，矢量e为合成地震数据包络，e_i(i＝1,2,3...n)为合成地震数据包络e第i个采样点数据；

d_i的希尔伯特变换如下：

${\hat{d}}_i=d_i*h_i.$

7.根据权利要求6所述的叠后地震反演方法，其特征在于，所述根据实际地震数据包络与合成地震数据包络优化波阻抗的初始模型，得到第二最优波阻抗，包括：

计算合成地震数据包络与实际地震数据包络之间的误差；

当误差处于预设范围内时，波阻抗的初始模型为第二最优波阻抗；

当误差不处于预设范围内时，设定波阻抗的初始模型的修改次数为k，k＝1，并按照如下方式进行迭代处理：

采用非线性全局寻优算法第k次修改波阻抗的初始模型，得到波阻抗的第k次修改模型；

根据波阻抗的第k次修改模型与实际地震数据得到第k次修改合成地震数据包络；

计算第k次修改合成地震数据包络与实际地震数据包络之间的误差；

当误差处于预设范围内时，波阻抗的第k次修改模型为第二最优波阻抗；

当误差不处于预设范围内时，将迭代处理中的k替代为k+1。

8.根据权利要求7所述的叠后地震反演方法，其特征在于，按照如下方式计算合成地震数据包络与实际地震数据包络之间的误差：

$J_e=\mathrm{\Σ}{(e_0-e)}^2=\underset{i=1}{\overset{n}{\Σ}}{(e_{0i}-e_i)}^2;$

其中，J_e为合成地震数据包络与实际地震数据包络之间的误差。

9.根据权利要求1所述的叠后地震反演方法，其特征在于，按照如下方式根据波阻抗的高频反演结果和波阻抗的低频反演结果得到波阻抗的最终反演结果：

Z＝Z_1h+Z_2l；

其中，矢量Z为波阻抗的最终反演结果，矢量Z_1h为波阻抗的高频反演结果，矢量Z_2l为波阻抗的低频反演结果。

10.一种叠后地震反演装置，其特征在于，包括：

地震数据合成模块，用于根据波阻抗的初始模型与实际地震数据得到合成地震数据；

高频反演模块，用于根据实际地震数据与合成地震数据优化波阻抗的初始模型，得到第一最优波阻抗；对第一最优波阻抗进行高通滤波，得到波阻抗的高频反演结果；

低频反演模块，用于根据实际地震数据得到实际地震数据包络，并根据合成地震数据得到合成地震数据包络；根据实际地震数据包络与合成地震数据包络优化波阻抗的初始模型，得到第二最优波阻抗；对第二最优波阻抗进行低通滤波，得到波阻抗的低频反演结果；

叠后地震反演模块，用于根据波阻抗的高频反演结果和波阻抗的低频反演结果得到波阻抗的最终反演结果。

11.根据权利要求10所述的叠后地震反演装置，其特征在于，所述地震数据合成模块具体用于：

在实际地震数据目的层段提取对应的子波；

根据初始模型得到反射系数；

根据所述反射系数与所述子波得到合成地震数据。

12.根据权利要求10所述的叠后地震反演装置，其特征在于，所述高频反演模块具体用于按如下方式得到第一最优波阻抗：

计算合成地震数据与实际地震数据之间的误差；

当误差处于预设范围内时，波阻抗的初始模型为第一最优波阻抗；

当误差不处于预设范围内时，设定波阻抗的初始模型的修改次数为j，j＝1，并按照如下方式进行迭代处理：

采用非线性全局寻优算法第j次修改波阻抗的初始模型，得到波阻抗的第j次修改模型；

根据波阻抗的第j次修改模型与实际地震数据得到第j次修改合成地震数据；

计算第j次修改合成地震数据与实际地震数据之间的误差；

当误差处于预设范围内时，波阻抗的第j次修改模型为第一最优波阻抗；

当误差不处于预设范围内时，将迭代处理中的j替代为j+1。

13.根据权利要求10所述的叠后地震反演装置，其特征在于，所述低频反演模块具体用于按如下方式得到第二最优波阻抗：

根据实际地震数据得到实际地震数据包络，并根据合成地震数据得到合成地震数据包络；

计算合成地震数据包络与实际地震数据包络之间的误差；

当误差处于预设范围内时，波阻抗的初始模型为第二最优波阻抗；

当误差不处于预设范围内时，设定波阻抗的初始模型的修改次数为k，k＝1，并按照如下方式进行迭代处理：

采用非线性全局寻优算法第k次修改波阻抗的初始模型，得到波阻抗的第k次修改模型；

根据波阻抗的第k次修改模型与实际地震数据得到第k次修改合成地震数据包络；

计算第k次修改合成地震数据包络与实际地震数据包络之间的误差；

当误差处于预设范围内时，波阻抗的第k次修改模型为第二最优波阻抗；

当误差不处于预设范围内时，将迭代处理中的k替代为k+1。

14.根据权利要求10所述的叠后地震反演装置，其特征在于，所述叠后地震反演模块具体用于按照如下方式根据波阻抗的高频反演结果和波阻抗的低频反演结果得到波阻抗的最终反演结果：

Z＝Z_1h+Z_2l；

其中，矢量Z为波阻抗的最终反演结果，矢量Z_1h为波阻抗的高频反演结果，矢量Z_2l为波阻抗的低频反演结果。