井约束非稳态相位校正方法与流程

本发明属于油气勘探地震资料处理领域，是基于井约束的非稳态相位校正方法，是地震数据波形一致性校正的重要内容。

现有技术

关于子波的零相位化，国外研究起步较早，Levy和Oldenburg、Longbottom等、White在Wiggins提出的盲反褶积基础上，提出使得输出地震数据非高斯性最强的常相位校正方法。相对Wiggins方法，他们提出的方法自由参数变少，算法更加稳定。为了考虑地震子波相位随时间变化的情况，Baan采用滑动时窗方法对子波相位进行估计，但是该方法有效的前提是相位变化分段平稳，为了避免这一限制条件，Baan和Fomel进一步将相位估计看做一个最优化反演问题，利用局部方差模最大准则进行非稳态相位估计，从而提高相位估计的精度，改善相位校正的稳定性。Fomel指出在缺少测井资料的情况下，可以采用地震道与其包络局部相似度最大的准则，估算非稳态数据的子波相位，该方法相对于局部方差模方法具有更大的动态范围，因此拾取的相位信息更加准确。

国内学者在相位校正方面也做了很多创新性的研究。周兴元证明，常相位校正只是在不改变地震道包络的前提下改变子波的波形。白治信基于常相位假设在频率域实现了最佳相位校正。陈必远等针对子波相位受复杂因素的影响，提出了一种随时间、空间及频率变化的相位校正方法。郭向宇等在分析已有各种子波处理方法优缺点的基础上，提出了一种新的估计混合相位子波相位的算法。李合群等同时对相位和时差进行了局部调整。宋宗平等在叠前道集中应用常相位校正技术，提高叠加剖面质量。国九英和周兴元基于最大叠加能量法，提出一种地表一致性相位校正方法，采用迭代算法求取炮点、检波点的相位校正算子，对地震记录进行校正。高少武等采用纯相位滤波器来逼近地表变化引起的相位变化，对地表一致性振幅处理和静校正处理后的地震数据进行了地表一致性相位校正处理。单联瑜等提出一种新的相位校正判别准则，先用最大方差模方法对模型的剩余相位进行校正，然后用解析法计算各道的相位校正角，从而实现精确的相位校正。

在地震勘探资料处理中，子波的零相位化有助于提高地震资料的分辨率、改善叠加剖面的质量。常规的相位校正方法是利用测井合成记录对井旁地震记录进行相位估计，然后对整条剖面进行常相位校正，该方法没有考虑地震子波 相位的非稳态性(相位随时间和空间变化)。虽然通过局部相似度方法利用最大方差模准则或包络最大相似度准则可以估计出随时间和空间变化的相位属性，但是由于零相位判别准则本身具有一定的局限性，因此精度有限。

由于近地表低、降速带风化层的横向厚度和速度的变化，以及地震排列上炮点耦合、检波点耦合条件的差异等近地表因素以及地下介质的非均一性和地下构造的复杂性等造成了地震子波相位随时间和空间变化(李振春等，2008)。将地震子波校正到零相位，不仅可以改善叠加剖面的质量，而且能够有效提高资料的分辨率。因此，在常规地震数据处理中，有必要进行相位校正方法方面的研究。

技术实现要素：

本发明的目的是针对上述问题，在局部地震属性和局部平面波模型下，提出了一种井约束的非稳态相位校正方法，该方法不仅考虑了子波相位的非稳态性，在相位校正时充分考虑其随时间和空间的变化,而且充分利用了测井合成地震记录进行相位校正精度高的优点。

本发明的技术方案是：

一种井约束非稳态相位校正方法，结合测井数据和井旁地震记录包络形成模型道，利用局部相似度准则估计出井旁地震数据的局部相位属相，再利用局部倾角信息将井旁地震道的相位信息插值到整个地震剖面，从而得到随时间和空间变化的相位信息，利用该信息即可对整条剖面进行相位校正。

上述方案进一步包括包括：

(1)利用测井数据和井旁地震道包络构建模型道，利用测井数据记录时间短于地震数据记录时间，将测井数据记录时间之外的数据用井旁地震道的包络来代替，合成与原地震数据记录时间相同的模型道；

(2)根据模型道，以相似度准则对井旁地震道进行非稳态相位估计；

(3)利用平面波预测算子将井旁地震道的局部相位延拓到其他位置；

(4)对所有测井数据位置处进行非稳态相位估计，然后对所有局部相位属性剖面取平均，得到随时间和空间变化的局部相位属性；

(5)根据估计出来的局部相位属性对地震剖面进行非稳态相位校正，并输 出相位校正之后的地震剖面。

上述方案还包括：

(2)估计井旁地震道的局部相位，根据模型道，以相似度准则对井旁地震道进行非稳态相位估计，两个时间序列a_i和b_i的相似系数定义为：

公式(1)看作两个最小二乘反演问题解的乘积，

c²＝c₁c₂ (2)

c₁＝(a^Ta)^-1(a^Tb) (3)

c₂＝(b^Tb)^-1(b^Ta) (4)

其中a、b分别是信号a(t)、b(t)的向量表示，x^Ty表示向量x和y点乘，由向量a、b的元素分别构造对角矩阵A、B，对方程局部化后相当于对反演加入正则化条件，标量c₁和c₂变为向量c₁和c₂，利用整型正则化，方程(3)、(4)变为：

c₁＝[λ²I+S(A^TA-λ²I)]^-1SA^Tb (5)

c₂＝[λ²I+S(B^TB-λ²I)]^-1SB^Tb (6)

其中S为平滑算子，向量c₁和c₂各元素对应相乘后得到随时间变化的局部相似度c，其结果采用整形正则化反演进行迭代近似得到上式的解。

(3)利用平面波预测算子将井旁地震道的局部相位延拓到其他位置，采用平面波预测滤波器算子将井数据位置处的相位信息插值到其他位置处，平面波预测滤波器用下面的微分方程来描述：

其中P(t,x)为地震波场，σ为随时间和空间变化的局部倾角；

假定地震剖面s由一系列地震道组成，s＝[s₁ s₂ ... s_N]^T，平面波解构滤波器(PWD)算子根据相邻道预测该道数据，并将原始道与预测道相减，数学表达式为：

r＝Ds (8)

其中r为剩余量，D为PWD算子，定义为：

其中I为单位算子，P为预测算子，通过沿同相轴倾角方向移动原始道来预测下一道地震数据，P_i,j表示由第i道预测第j道的算子；

式(8)通过正则化共轭梯度法最小化预测误差r估计出主要的局部倾角信息，PWC算子是PWD算子的逆，可以通过递归算法快速实现PWC算子：

c₁＝s₁,c_k＝s_k+P_k-1,kc_k-1,k＝2,3,...,N (10)。

(5)根据估计出来的局部相位属性对地震剖面进行非稳态相位校正，并输出相位校正之后的地震剖面，地震道x(t)经过相位旋转后变为：

x_rot(t)＝x(t)cosφ(t)+H[x(t)]sinφ(t) (11)

式中，φ(t)为随时间变化相位旋转角度，H[x(t)]为地震道x(t)的Hilbert变换，采用不同的零相位判别准则，得到不同的相位校正方法。

发明效果

井约束非稳态相位校正方法考虑了子波相位的非稳态性，而且充分利用了测井合成地震记录进行相位校正精度高的优点。其具体优势和特点表现在以下几个方面：

一、方法效果的可靠性。该方法综合不仅考虑了随时间和空间变化的相位信息，而且考虑了地震剖面的局部构造信息，因此更为合理。

二、操作简单易实现。该方法流程及参数设置简单，运算速度快，不依据任何假设，不受常规相位校正方法应用的限制。

三、相位校正精度高。该方法不受零相位判别准则局限性的限制，并且充分利用了测井合成地震记录进行相位校正精度高的优点

附图说明

图1为本发明的井约束非稳态相位校正技术流程图；

图2为对一维合成数据采用不同相位校正方法的效果对比图；其中，图2(a)为波阻抗测井曲线，图2(b)为反射系数曲线，图2(c)为与零相位Ricker 子波褶积之后的地震记录，图2(d)为对合成地震记录进行-90°到0度线性变化相位旋转之后的地震记录，图2(e)、图2(f)、图2(g)分别校正结果；图2(h)为这三种方法估算的相位角对比图。

图3是二维实际资料单道相位校正；其中，图3(a)原始叠加剖面，图3(b)为应用PWD估计的局部倾角场。图3(c)为井旁地震数据，图3(d)为井旁地震数据包络和测井合成数据合成的模型道，图3(e)为估计的相位角，图3(f)为校正后的地震道。

图4是二维实际资料相位校正；其中图4(a)整条地震剖面的相位信息，图4(b)校正结果，图4(c)、图4(d)分别为图4(a)、图4(b)的局部放大图。

具体实施方式

本实施例以某探区二维地震资料为目标靶区，应用本方法对该资料进行处理，以验证本方法的效果，流程图见图1。

1)利用测井数据和井旁地震道包络构建模型道；

2)估计井旁地震道的局部相位；

3)利用平面波预测算子将井旁地震道的局部相位延拓到其他位置；

4)对所有测井数据位置处进行非稳态相位估计，然后对所有局部相位属性剖面取平均，得到随时间和空间变化的局部相位属性；

5)根据估计出来的局部相位属性对地震剖面进行非稳态相位校正，并输出相位校正之后的地震剖面。

发明的详细操作简图如图1所示。主要技术关键点为以下四个：①、构建模型道；②、估计井旁地震道的局部相位；③、局部相位延拓；④、对剖面进行非稳态相位校正。

井约束的非稳态相位校正方法实现过程如下：

(1)利用测井数据和井旁地震道包络构建模型道

当有测井资料时，首先采用测井资料合成地震记录，由于测井数据记录时间 往往短于地震数据记录时间，为此，测井数据记录时间之外的数据用井旁地震道的包络来代替，这样就合成了一道与原地震数据记录时间相同的模型道。

(2)估计井旁地震道的局部相位

根据模型道，可以对井旁地震道进行非稳态相位估计。相位估计依赖于相似度准则。两个时间序列a_i和b_i的相似系数定义为：

公式(1)可以看作两个最小二乘反演问题解的乘积，

c²＝c₁c₂ (2)

c₁＝(a^Ta)^-1(a^Tb) (3)

c₂＝(b^Tb)^-1(b^Ta) (4)

其中a、b分别是信号a(t)、b(t)的向量表示，x^Ty表示向量x和y点乘。由向量a、b的元素分别构造对角矩阵A、B，对方程局部化后相当于对反演加入正则化条件，标量c₁和c₂变为向量c₁和c₂，利用整型正则化，方程(3)、(4)变为：

c₁＝[λ²I+S(A^TA-λ²I)]^-1SA^Tb (5)

c₂＝[λ²I+S(B^TB-λ²I)]^-1SB^Tb (6)

其中S为平滑算子，向量c₁和c₂各元素对应相乘后得到随时间变化的局部相似度c，其结果可以采用整形正则化反演进行迭代近似得到上式的解(Claerbout，1992)。

(3)利用平面波预测算子将井旁地震道的局部相位延拓到其他位置。

为了将井数据位置处的相位信息插值到其他位置处，本专利采用平面波预测滤波器算子来解决这一问题。

平面波预测滤波器来源于描述地震数据的局部平面波模型，可以用下面的微分方程来描述：

其中P(t,x)为地震波场，σ为随时间和空间变化的局部倾角。

假定地震剖面s由一系列地震道组成，s＝[s₁ s₂ ... s_N]^T。平面波解构滤波 器(PWD)算子(Fomel，2002)根据相邻道预测该道数据，并将原始道与预测道相减，数学表达式为：

r＝Ds (8)

其中r为剩余量，D为PWD算子，定义为：

其中I为单位算子，P为预测算子，通过沿同相轴倾角方向移动原始道来预测下一道地震数据，P_i,j表示由第i道预测第j道的算子。式(8)可以通过正则化共轭梯度法最小化预测误差r估计出主要的局部倾角信息。该反演问题的约束条件通常是使得数据空间倾角平滑变化(Fomel，2007)。

PWC算子是PWD算子的逆，可以通过递归算法快速实现PWC算子：

c₁＝s₁,c_k＝s_k+P_k-1,kc_k-1,k＝2,3,...,N (10)

该公式的具体推导过程请参考文献(Fomel等人，2006)。

(4)对所有测井数据位置处进行非稳态相位估计，然后对所有局部相位属性剖面取平均，得到随时间和空间变化的局部相位属性。

(5)根据估计出来的局部相位属性对地震剖面进行非稳态相位校正，并输出相位校正之后的地震剖面。

地震道x(t)经过相位旋转后变为：

x_rot(t)＝x(t)cosφ(t)+H[x(t)]sinφ(t) (11)

式中，φ(t)为随时间变化相位旋转角度，H[x(t)]为地震道x(t)的Hilbert变换，采用不同的零相位判别准则，可以得到不同的相位校正方法。

相对于局部方差模和包络局部相似度准则，基于井数据的非稳态相位校正方法具有最高的精度，为了证明这一点，我们采用某探区的测井数据进行测试。图2(a)为波阻抗测井曲线，图2(b)为反射系数曲线，图2(c)为与零相位Ricker子波褶积之后的地震记录，图2(d)为对合成地震记录进行-90°到0度线性变化相位旋转之后的地震记录，对该记录分别在最大方差模准则、包络 最大相似度准则、测井约束最大相似度准则下进行非稳态相位校正，校正结果分别如图2(e)、图2(f)、图2(g)所示，对比可以看出，由于最大方差模有效的前提是反射系数序列符合随机超高斯分布，因此该方法将各个波峰对应的位置作为一个反射层，将其校正为零相位，但是存在薄层时，无法将其校正到正确的相位。包络最大相似度准则受包络计算的影响，同样无法校正薄层地震响应的相位，而井约束相位校正方法的精度最高。图2(h)为这三种方法估算的相位角对比图，从图中可以看出，0.02s位置左右三种方法估计的相位角比较接近，而其他位置处利用最大方差模和包络最大相似度准则无法估计出正确的相位角。

为了进一步验证本专利方法的有效性，对某探区实际资料进行相位校正处理，图3(a)原始叠加剖面，图中五角星表示井所在的位置，实线部分为测井的时间范围。图3(b)为应用PWD估计的局部倾角场。图3(c)为井旁地震数据，图3(d)为井旁地震数据包络和测井合成数据合成的模型道。下面采用局部相似度对井旁地震道进行相位估计，估计的相位角如图3(e)所示，利用该相位角对井旁地震道进行相位校正，校正后的地震道如图3(f)所示，从图中可以看出，该地震道子波基本校正到零相位，尤其是箭头所示子波相位校正效果更为明显。

利用测井合成数据对井旁数据进行相位估计以后，采用平面波预测算子将该相位信息插值到其他位置，从而得到整条地震剖面的相位信息，如图4(a)所示，利用插值后的相位信息对原始叠加剖面进行相位校正，校正结果如图4(b)所示，图4(c)、图4(d)分别为图4(a)、图4(b)的局部放大图。从图中可以看出，采用本专利所提出的相位校正方法可以有效将地震剖面的子波校正到零相位，从而改善叠加剖面的质量，提高叠加剖面的分辨率。

相对于线性插值等数学插值方法，平面波预测方法考虑了地震剖面的局部构造信息，因此更为合理。