CN107621654A

CN107621654A - 一种基于最大相关熵的地震叠后波阻抗反演方法

Info

Publication number: CN107621654A
Application number: CN201710753538.3A
Authority: CN
Inventors: 王峣钧; 程三; 厍斌; 胡光岷; 张家树
Original assignee: University of Electronic Science and Technology of China; Southwest Jiaotong University
Current assignee: University of Electronic Science and Technology of China; Southwest Jiaotong University
Priority date: 2017-08-29
Filing date: 2017-08-29
Publication date: 2018-01-23

Abstract

本发明公开一种基于最大相关熵的地震叠后波阻抗反演方法，应用于地震勘探领域，采用地震反演模型为Robinson褶积模型，构建叠后波阻抗正演模型；并根据最大相关熵准则建立目标函数，最后通过迭代求解目标函数，实现地震叠后波阻抗反演，本申请的方法通过采用最大相关熵准则建立目标函数，使得目标函数能够更好的克服非高斯噪声，并且相对于传统方法鲁棒性更强，精度更高。

Description

一种基于最大相关熵的地震叠后波阻抗反演方法

技术领域

本发明属于地震勘探领域，特别涉及一种地震反演技术。

背景技术

地震反演是地震勘探中的关键技术。地震反演是根据地球物理观测数据，利用地震波在地下介质中的传播规律，来推测地球内部的结构、形态及物质成分，同时定量计算各种相关的地球物理参数的过程。地震反演为油气的勘探开发提供了重要的依据。

基于模型的反演是地震反演中一种重要且有效的方法。这种方法采用迭代的思想，从初始地质模型出发，利用合成记录与真实记录之间的残差计算出的模型修正量来更新模型，使模型正演地震资料与实际地震资料最佳吻合，最终的模型数据便是反演结果。

基于模型的反演通常采用L₂范数作为衡量误差的准则，这需要假设噪声服从高斯分布。但是当地震观测信号中包含脉冲噪声等非高斯噪声干扰时，由于L₂范数对大的脉冲噪声干扰比较敏感，导致算法失效。随着油气勘探的逐步深入，勘探目标逐渐由构造油气藏转为更为复杂的岩性油气藏，这就要求我们必须更加精细地描述油气藏，以便获取更可靠的储层信息，降低勘探开发的风险。对于某些复杂储层，噪声高斯假设很难成立。为了处理非高斯噪声，一些新的泛函准则开始逐渐应用在地震反演中。例如L₁范数，L₁/L₂混合范数(Bube and Nemeth 2007)、Huber范数(Guitton and Symes 2003)、自适应广义极值法(Zhang et al.2013)等。上述的诸多方法可以有效地降低脉冲噪声的影响，但是仍然不能完全消除脉冲噪声的干扰，因此当数据所含脉冲噪声的比例较大时，这些方法的鲁棒性往往较差。

发明内容

为解决上述技术问题，本申请提出了一种基于最大相关熵的地震叠后波阻抗反演方法，采用最大相关熵准则替换原有最小二乘准则，使得目标函数能够更好的克服非高斯噪声。

本申请采用的技术方案为：一种基于最大相关熵的地震叠后波阻抗反演方法，包括：

S1、构建叠后波阻抗正演模型；

S2、根据步骤S1建立的叠后波阻抗正演模型，根据最大熵准则建立目标函数；

S3、对目标函数进行迭代求解。

进一步地，步骤S1所述正演模型为：

s＝Gx

其中，G表示正演矩阵，x表示地层波阻抗系数的对数，s表示地震道。

进一步地，步骤S2所述的目标函数为：

其中，e_i表示残差向量e的第i个分量，s_i表示地震道s的第i个分量，G(i,：)为正演矩阵G的第i行；i＝1,2,3，…，n；N表示地震道的长度；σ为核宽，且σ＞0。

更进一步地，将一阶Tikhonov正则化算子作为约束项，则目标函数为：

其中，D为一阶差分算子，λ为正则化参数。

本发明的有益效果：本申请的一种基于最大相关熵的地震叠后波阻抗反演方法，采用地震反演模型为Robinson褶积模型，构建叠后波阻抗正演模型；并根据最大相关熵准则建立目标函数，使得目标函数能够更好的克服非高斯噪声；最后通过迭代求解目标函数，实现地震叠后波阻抗反演，本申请的方法相对于传统方法鲁棒性更强，精度更高，更有利于进行地质勘探中的资料处理。

附图说明

图1为本申请的方案流程图；

图2为本申请实施例提供的marmousi速度模型图和本申请方法生成叠后地震记录图；

其中，图2(a)为已知的marmousi速度模型图；图2(b)为通过本申请褶积方法生成叠后地震记录图；

图3为本申请实施例提供的含有5db高斯噪声情形下L₂范数反演结果和最大相关熵反演结果对比图；

其中，(a)为5db高斯噪声的合成地震记录图；(b)为L₂范数反演结果图；(c)为最大相关熵反演结果图；

图4为本申请实施例提供的在合成的叠后记录中添加5db拉普拉斯噪声效果图以及L₁范数、Huber范数和最大相关熵的反演结果对比图；

其中，(a)为含5db拉普拉斯噪声的合成地震记录图；(b)为L₂范数反演结果图；(c)为Huber范数反演结果图；(d)为最大相关熵反演结果图；

图5为本申请实施例提供的西南某页岩气工区叠后地震记录及L₁范数、Huber范数和最大相关熵的反演结果对比图；

其中，(a)为西南某页岩气工区叠后地震记录图；(b)为L₁范数反演结果图；(c)为L₂范数反演结果图；(d)为Huber范数反演结果图；(e)为最大相关熵反演结果图。

具体实施方式

为便于本领域技术人员理解本发明的技术内容，下面结合附图对本发明内容进一步阐释。

如图1所示为本申请的方案流程图，本申请的技术方案为：一种基于最大相关熵的地震叠后波阻抗反演方法，包括：

S1、构建叠后波阻抗正演模型；

本发明采用地震反演模型为Robinson褶积模型：

s(t)＝w(t)*r(t)+n(t) (1)

其中，s(t)为地震道；w(t)为地震子波；r(t)为地层反射系数；n(t)为随机噪声。令s＝[s₁,s₂,…,s_n-1]^T，r＝[r₁,r₂,…,r_n-1]^T，w＝[w₁,w₂,…,w_n-1]^T。(1)式可以写成矩阵形式：

根据Russell近似公式，当地层波阻抗系数之差远小于波阻抗的绝对值时，有：

因此反射系数可表示为：

将(4)式带入(2)式，可得到：

令：

则最终得到正演模型为：

s＝Gx (6)

通过最大相关熵准则建立目标函数：

其中，e_i＝s_i-G(i,:)·x，e_i表示残差向量e的第i个分量，s_i表示地震道s的第i个分量，G(i,：)为正演矩阵G的第i行；i＝1,2,3，…，n；N表示地震道的长度；σ为核宽，且σ＞0。将目标函数转换为取极小值形式：

为了提高反演的稳定性，我们将一阶Tikhonov正则化算子作为约束项添加到目标函数中，则新的目标函数可以表示为：

其中，D为一阶差分算子，λ为正则化参数，λ值的大小可为经验法或L曲线法确定。

S3、对目标函数进行迭代求解。从(9)式结构来看，该公式是高度非线性的，如果直接求导很难得到解析解，因此本申请采用迭代求解方式对该目标函数进行求解。

在处理病态优化问题时，共轭梯度算法和拟牛顿算法都是常用的算法。通常情况下，前者的收敛速度较后者低，但后者在计算海森矩阵逆时需要花费较大存储空间。本发明基于求解算法稳定性、低存储性、和迭代速率等综合考虑，采用Zhang et.al(2013)提出的拟共轭梯度反演算法(该反演算法结合了共轭梯度算法和拟牛顿算法)进行求解运算，此外，由于原始算法中需要根据不同的广义极值分布目标函数来计算步长不具有普适性，本申请采用Dianne(1991)提供的步长搜索的代码修改该步骤，由于该搜索方法只需要给出计算目标函数的值和梯度的两个函数便能够搜索得到步长，因此可以使得上述算法能适应任意目标函数情形。

以下通过两个实例对本申请的技术效果进行说明：

实施例1

在本实施例中采用经典的marmousi模型作为参考模型，生成叠后地震记录，并用提出的最大相关熵算法进行反演。

已知的marmousi速度模型如图2(a)所示，通过褶积方法生成叠后地震记录如图2(b)所示，地震子波是主频为40Hz的雷克子波，合成的地震记录的采样间隔为2ms。

为分析最大相关熵准则对含噪声记录的鲁棒性，首先分析在如图3(a)所示含有5db 高斯噪声情形下对比了如图3(b)所示的L₂范数反演结果和如图3(c)所示的最大相关熵反演结果。接下来，在合成的叠后记录中添加如图4(a)所示的5db拉普拉斯噪声(非高斯噪声)，由于L₁范数和Huber范数能压制住非高斯噪声，因此对比了L₁范数、Huber 范数和最大相关熵的反演结果如图3(b)、图4(c)、图4(d)所示。从含噪声结果来看，在高斯噪声条件下，最大相关熵的反演结果要比L₂范数的结果好；而对于非高斯噪声，最大相关熵的反演结果也明显好于比L₁范数和Huber范数的反演结果。由此可见，最大相关熵准则条件下的地震反演方法比传统的反演方法具有更强的鲁棒性。

图中的CMP表示共中心点(Common Middle Point)道集的道集号，Time表示时间，Amplitude表示振幅。

实施例2

本例中应用西部某页岩气工区实际资料进行了地震反演，叠后地震记录如图5(a)所示。对于页岩气，干酪根含量(TOC)是预测页岩中是否存在油气的一项关键指标。分析本工区交会图，密度和波阻抗与TOC都有着密切的联系。因此能否得到效果好的TOC依赖于地震波阻抗反演结果。

为了验证最大相关熵方法的有效性，以及相比于传统算法的优劣，如图5(b)至图5(e)所示本部分分别用L₁、L₂、Huber范数以及最大相关熵方法对实际数据进行了反演。四种方法得到的反演结果在井位置的波阻抗曲线与真实测井曲线；无论从反演剖面还是曲线对比来看，本申请提出的基于最大相关熵准则的地震反演方法比传统的反演方法精度更高，鲁棒性更强。

本领域的普通技术人员将会意识到，这里所述的实施例是为了帮助读者理解本发明的原理，应被理解为本发明的保护范围并不局限于这样的特别陈述和实施例。对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的权利要求范围之内。

Claims

1.一种基于最大相关熵的地震叠后波阻抗反演方法，其特征在于，包括：

S1、构建叠后波阻抗正演模型；

S3、对目标函数进行迭代求解。

2.根据权利要求1所述的一种基于最大相关熵的地震叠后波阻抗反演方法，其特征在于，步骤S1所述正演模型为：

s＝Gx

其中，G表示正演矩阵，x表示地层波阻抗的对数，s表示地震道。

3.根据权利要求1所述的一种基于最大相关熵的地震叠后波阻抗反演方法，其特征在于，步骤S2所述的目标函数为：

<mrow> <mi>&phi;</mi> <mrow> <mo>(</mo> <mi>x</mi> <mo>)</mo> </mrow> <mo>=</mo> <mi>m</mi> <mi>a</mi> <mi>x</mi> <mfrac> <mn>1</mn> <mi>N</mi> </mfrac> <munderover> <mo>&Sigma;</mo> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>n</mi> </munderover> <mfrac> <mn>1</mn> <mrow> <msqrt> <mrow> <mn>2</mn> <mi>&pi;</mi> </mrow> </msqrt> <mi>&sigma;</mi> </mrow> </mfrac> <mi>exp</mi> <mrow> <mo>(</mo> <mo>-</mo> <msup> <msub> <mi>e</mi> <mi>i</mi> </msub> <mn>2</mn> </msup> <mo>/</mo> <mn>2</mn> <msup> <mi>&sigma;</mi> <mn>2</mn> </msup> <mo>)</mo> </mrow> </mrow>

4.根据权利要求3所述的一种基于最大相关熵的地震叠后波阻抗反演方法，其特征在于，将一阶Tikhonov正则化算子作为约束项，则目标函数为：

<mrow> <mi>&phi;</mi> <mrow> <mo>(</mo> <mi>x</mi> <mo>)</mo> </mrow> <mo>=</mo> <mi>m</mi> <mi>i</mi> <mi>n</mi> <mo>&lsqb;</mo> <mn>1</mn> <mo>-</mo> <mfrac> <mn>1</mn> <mi>N</mi> </mfrac> <munderover> <mo>&Sigma;</mo> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>n</mi> </munderover> <mfrac> <mn>1</mn> <mrow> <msqrt> <mrow> <mn>2</mn> <mi>&pi;</mi> </mrow> </msqrt> <mi>&sigma;</mi> </mrow> </mfrac> <mi>exp</mi> <mrow> <mo>(</mo> <mo>-</mo> <msup> <msub> <mi>e</mi> <mi>i</mi> </msub> <mn>2</mn> </msup> <mo>/</mo> <mn>2</mn> <msup> <mi>&sigma;</mi> <mn>2</mn> </msup> <mo>)</mo> </mrow> <mo>+</mo> <mi>&lambda;</mi> <mo>|</mo> <mo>|</mo> <mi>D</mi> <mi>x</mi> <mo>|</mo> <msubsup> <mo>|</mo> <mn>2</mn> <mn>2</mn> </msubsup> <mo>&rsqb;</mo> </mrow>

其中，D为一阶差分算子，λ为正则化参数。