井震融合平面波分解地层倾角计算方法与应用

穆立华,杨国涛,高文中,张建坤,王方鲁

(中国石油天然气股份有限公司冀东油田分公司勘探开发研究院,河北唐山063004)

穆立华,杨国涛,高文中,等.井震融合平面波分解地层倾角计算方法与应用[J].石油物探,2017,56(6):-826

MU Lihua,YANG Guotao,GAO Wenzhong,et al.Formation dip calculation using plane-wave decomposition based on logging and seismic integration[J].Geophysical Prospecting for Petroleum,2017,56(6):-826

井震融合平面波分解地层倾角计算方法与应用

穆立华,杨国涛,高文中,张建坤,王方鲁

(中国石油天然气股份有限公司冀东油田分公司勘探开发研究院,河北唐山063004)

基于倾角测井得到的井口附近高精度倾角数据,提出了一种井震融合的平面波分解地层倾角计算方法。在有倾角测井的探区,根据Tikhonov正则化的思想在平面波分解滤波器中引入地层倾角测井数据约束,井震融合计算空间任意点的地层视倾角,通过三角变换将视倾角转化为真倾角和真倾向。理论模型测试及实际资料应用结果表明,该方法估计的地层倾角和地层厚度转换结果比常规平面波分解计算得到的结果更精确,具有一定的实际应用价值与推广意义。

倾角测井数据;平面波分解;井震融合;模型正则化;地层倾角;地层倾角测井信息;滤波器

地层倾角信息对于地震采集、处理、解释都具有重要意义。在采集阶段,地层倾角信息可以用于指导观测系统设计;在处理阶段,地层倾角对于速度分析、速度建模等有很大影响;在油藏描述阶段,真实的地层倾角信息有助于研究地质构造和沉积环境,进而追踪地下油气的分布情况。地层倾角的计算方法有局部倾斜叠加法[1]、复地震道分析法[2]、结构张量法[3-4]、平面波分解法[5]以及相干估计法[6-7]等。CLAERBOUT[5]提出了平面波分解滤波器,通过有限差分法求解平面波微分方程得到地震数据的局部斜率。一般情况下,局部平面波模型能够方便地表达地震数据,但是需要预测较多的滤波器系数,而这些系数需要人工交互处理,存在一定的局限性。FOMEL[8]对CLAERBOUT滤波器[5]进行了改进,结合适当的数据正则化方法,使滤波器系数预测减少到最少,局部平面波斜率预测结果具有清晰的物理意义。FOMEL[9]将改进的平面波分解滤波器应用于局部斜率估计、断层预测、数据插值和噪声衰减。杜婧等[10]使用平面波分解滤波器估计局部斜率属性,用于VSP波场分离研究。刘洋等[11]使用平面波分解滤波器得到地震数据的局部倾角属性,设计出预测构造数据体和构造走向的滤波器,该滤波器能在去噪的同时保护断层信息。刘海燕[12]利用平面波分解滤波器估计的地层倾角对地震数据进行调和,使其C3相干体算法的分辨率更高,有效消除了地层倾角导致的低相干干扰,突出了有效断层信息。陈鑫[13]在速度无关时差校正中使用平面波分解滤波器获得局部同相轴斜率信息,用于替代相关速度参数进行动校正处理。

分析发现,常规平面波分解方法只利用地震数据,未充分利用老油区的大量倾角测井数据。本文基于平面波分解滤波器计算局部斜率的方法,将地层倾角测井信息作为先验信息,结合Tikhonov正则化思想,提出了井震融合平面波分解地层倾角计算方法;利用理论数据和实际资料对方法进行了测试。

1 平面波分解方法原理

根据局部平面波的物理模型,可以定义局部平面波微分方程:

(1)

式中:P(t,x)代表局部平面波波场,σ代表局部斜率,两者都与时间t和炮检位置x有关。公式(1)表明,在x1处观测到的波场可以被x2处观测到的波场完美地预测,它们相差一个延迟时间σΔx,其中Δx=|x1-x2|。即平面波时延满足:

(2)

将方程(2)对时间和空间分别进行Z变换,得:

(3)

(4)

其中,

(5)

式中:N为空间滤波器的阶数;bkσ为局部斜率σ的函数。在低频处(N=1),通常采用三点中心滤波器,即:

(6)

方程(4)相当于一个二维滤波算子,将方程(4)代入方程(3),可以得到与地层倾角σ相关的方程:

(7)

方程(7)即利用平面波分解滤波器计算地层局部斜率的基本公式。

2 井震融合的地层倾角计算

(8)

方程(8)的本质是求解一个非线性问题。根据最优化理论,可将方程(8)转化为线性方程(9),然后使用高斯-牛顿法求解:

(9)

式中:Δσ表示斜率的增量,C′(σ)是滤波器C(σ)的导数。在迭代求解线性方程(9)的过程中,初始σ值被不断更新,逐步趋近于真实斜率。

在地震数据信噪比较低、分辨率不高的情况下,直接求解方程(9)往往不稳定,需要引入合适的正则化。FOMEL[7]给出正则化项:

(10)

式中:D为梯度算子,ε为正则化比例系数。通过增加方程(10)中正则化项,可以得到更加可靠的局部斜率。

在成熟探区,通常有各种测井数据,由井口附近的倾角测井数据能够得到较地震数据更为可靠的局部倾角信息。根据Tikhonov正则化的思想,在公式(8) 中加入先验测井信息约束,则得到井震融合的倾角估计方程:

(11)

式中:ξσ表示引入的测井约束,本文主要利用地层倾角测井信息。

由公式(11)得到的是沿主测线和联络测线方向的局部倾角,称之为视倾角,需要经过一定的三角换算转化为真倾角和真倾向。转换公式如下:

(13a)

(13b)

式中:θx和θy分别为沿主测线和联络测线的视倾角;θ为地层真倾角;φ为真倾向与y方向视倾向之间的夹角(正北方为y方向,正东方为x方向)。由公式(13)可以得到真实的地层倾角和倾向。

基于地层倾角信息约束的平面波分解地层倾角估计方法有两方面优势:①将经过预处理后的倾角测井信息作为约束加入到方程(11)中,有倾角测井值的网格不更新σ,在一定程度上减少了反演解的误差,降低了反演迭代次数,提高了反演效率。②利用倾角测井信息约束提高地层倾角估计的精度,有助于精细刻画薄互层、微幅构造等小尺度地质体。

3 测试与应用

3.1 理论模型测试

采用图1a所示合成3D偏移数据体[6]对本文方法进行了测试,图1b为合成数据的平面图以及4口合成地层倾角测井(垂直井)井口位置示意图。将图1b 中红五星周围3口井作为约束井,采用公式(11) 得到x方向和y方向的视倾角如图2所示。再利用公式(13)求出地层真倾角及倾向,如图3所示。为了进一步说明本文方法的正确性,我们对比了平面波分解方法未加入约束、加入井约束预测的倾角曲线以及合成地层倾角测井曲线,如图4所示。可以看出,加入井约束预测的倾角更加趋近于地层倾角测井的结果。

3.2 实际资料应用

3.2.1 地层倾角估计

南堡凹陷为经过多期构造运动形成的裂陷盆地,地质特征复杂。位于南堡凹陷南部的南堡3号(NP3)构造地层倾角变化大,严重影响地质认识程度。采用本文井震融合的地层倾角估计方法估算了南堡3号构造的地层倾角数据体。图5a为NP3号构造1599线偏移剖面,利用公式(11)迭代计算了地层视倾角/向,图5b和图5c为利用公式(13)由视倾角/向转换的地层真倾角/向。

根据南堡3号工区NP3-19井、NP306X1井和PG2井的轨迹,抽取未加入井约束地层倾角估计方法得到的地层倾角数据体中相应的倾角值,将其与地层倾角测井结果作对比,如图6所示。

图1 3D偏移数据体(引自参考文献[6])(a)及其井口位置(垂直井)(b)

图2 3D偏移数据体x方向(a)与y方向(b)的视倾角

图3 3D偏移数据体真倾角(a)与真倾向(b)

图4 未加入井约束、加入井约束预测的倾角与合成地层倾角测井曲线对比

图5 南堡3号构造1599线偏移剖面(a)及本文方法估计的真倾角(b)和真倾向(c)

由于地层倾角测井的采样率和井轨迹数据采样率不一致,故两者不可能一一对应,但是从图6可以清晰地看出,未加入井约束地层倾角估计方法得到的3口井地层倾角变化趋势与地层倾角测井结果基本吻合。

将NP3-19井、NP306X1井和PG2井测井数据作为约束加入公式(11)进行地层倾角估算,再抽取离NP3-19井较近的NP3-20井地层倾角信息进行对比,结果如图7所示。由图7可见,加入测井数据约束后估算出来的地层倾角更加趋近于地层倾角测井的结果,估算精度得到了提高。

图6 未加入井约束估算的地层倾角与测井结果对比 a NP3-19井; b NP306X1井; c PG2井

图7 NP3-20井地层倾角估算结果与测井结果对比

对比PG2,NP3-20以及P3002三口井的地层产状(表1)可以看出,加入测井数据约束后,地震数据体地层产状估算精度明显提升,与实际测量结果相对误差小于1%。P3002井与约束井距离较远(约800m左右),约束效果不明显,地层产状相对误差为1%左右,也在可接受范围内。

表1 毛庄组地层产状对比 (°)

3.2.2 地层厚度校正

为了说明井震融合平面波分解地层倾角估计在薄互层、微幅构造等小尺度地质体精细刻画中的作用,这里以地层厚度校正为例,对研究区内NP3-82井、PG2井、NP3-80井进行测井地层厚度校正,将地层视厚度换算成地层真厚度。

根据图5所示的地层倾角信息,利用测井曲线校正公式[14],将NP3-82井、PG2井、NP3-80井的垂直厚度转换成真厚度(图8)。图8a为NP3-82井—PG2井—NP3-80井地层垂直厚度对比图;图8b和图8c分别为未加入井约束和加入井约束方法估算地层倾角后,校正得到的NP3-82井—PG2井—NP3-80井地层真厚度对比图。从图8可以看出,PG2井和NP3-80井地层水平厚度变化较大,与海相地层沉积稳定的地质规律不符(图8a);采用未加入井约束方法估算地层倾角并校正地层厚度后,PG2井和NP3-80井地层水平厚度变化趋势有所减缓(图8b);而采用本文井震融合方法预测地层倾角并对地层厚度进行校正后,地层水平厚度变化趋势更加平稳,更加符合海相地层地质规律(图8c),小层划分精度也得到了提高。

图8 NP3-82井—PG2井—NP3-80井地层厚度对比 a 垂直厚度; b 未加入井约束方法得到的地层真厚度; c 本文井震融合方法得到的地层真厚度

4 结束语

本文提出一种井震融合的平面波分解地层倾角估算方法,根据Tikhonov正则化思想,融合地震和地层倾角测井信息实现地层倾角的估计。该方法引入测井信息约束后,不仅使平面波分解滤波器估计地层倾角的精度得到了提升,减小了与实测结果的误差,而且提高了迭代反演的收敛速率,因此具有一定的实际应用价值。

[1] OTTOLINI R.Signal/noise separation in dip space[R].San Francisco:Stanford Exploration Project,1983

[2] BARNES A E.A tutorial on complex seismic trace analysis[J].Geophysics,2007,72(6):W33-W43

[3] FEHMERS G C,HÖCKER C F W.Fast structural interpretation with structure-oriented filtering[J].Geophysics,2003,68(4):1286-1293

[4] HALE D.Local dip filtering with directional Laplacians[R].Golden:Center for Wave Phenomena,2007

[5] CLAERBOUT J F.Earth soundings analysis:processing versus inversion[M].Cambridge,Massachusetts,USA:Blackwell Scientific Publications,1992:94-98

[6] 陈双全,季敏.地震数据结构张量相干计算方法[J].石油物探,2012,51(3):233-238

CHEN S Q,JI M.Structure tensor coherence computation method of seismic data[J].Geophysical Prospecting for Petroleum,2012,51(3):233-238

[7] 乐靖,王晖,范廷恩,等.基于地震等时格架的倾角导向储层静态建模方法[J].石油物探,2017,56(3):449-458

LE J,W H,FAN T E,et al.A method of dip steering reservoir static modeling based on seismic isochronal stratigraphic framework[J].Geophysical Prospecting for Petroleum,2017,56(3):449-458

[8] FOMEL S.Plane wave prediction in 3-D[R].San Francisco:Stanford Exploration Project,2001

[9] FOMEL S.Applications of plane-wave destruction filters[J].Geophysics,2002,67(6):1946-1960

[10] 杜婧,王尚旭,刘国昌,等.基于局部斜率属性的 VSP 波场分离研究[J].地球物理学报,2009,52(7):1867-1872

DU J,WANG S X,LIU G C,et al.VSP wavefield separation using local slopes attribute[J].Chinese Journal of Geophysics,2009,52(7):1867-1872

[11] 刘洋,王典,刘财,等.局部相关加权中值滤波技术及其在叠后随机噪声衰减中的应用[J].地球物理学报,2011,54(2):358-367

LIU Y,WANG D,LIU C,et al.Weighted median filter based on local correlation and its application to poststack random noise attenuation[J].Chinese Journal of Geophysics,2011,54(2):358-367

[12] 刘海燕.C3相干算法的改进及其在断层识别中的应用[D].长春:吉林大学,2013

LIU H Y.The improvement of C3 coherence algorithm and the application in fault identification[D].Changchun:Jilin University,2013

[13] 陈鑫.基于平面波解构滤波的速度无关时差校正方法研究[D].长春:吉林大学,2015

CHEN X.Study on velocity-independent moveout correction using plane-wave destruction filtering[D].Changchun:Jilin University,2015

[14] 张卫江,孙广伯.容积法储量计算中定向井油层有效厚度的校正[J].石油勘探与开发,1998,25(5):51-52

ZHANG W J,SUN G B.Correction of effective thickness of directional wells in volumetric reserves calculation[J].Petroleum Exploration and Development,1998,25(5):51-52

(编辑：戴春秋)

Formationdipcalculationusingplane-wavedecompositionbasedonloggingandseismicintegration

MU Lihua,YANG Guotao,GAO Wenzhong,ZHANG Jiankun,WANG Fanglu

(ExplorationandDevelopmentResearchInstituteofJidongOilfieldBranchCompany,PetroChina,Tangshan063004,China)

Dip logging can obtain high precision dip angle data near the well trajectory.This paper proposes a method of calculating the formation dip angle by the plane-wave decomposition,based on logging and seismic integration.In the exploration area with dip logging,based on the idea of Tikhonov regularization,this method was used to calculate the apparent formation dip at any point in space by logging and seismic integration,through plane-wave decomposition filter with dip logging data as constraints.The calculated apparent dip was subsequently transferred to the true dip angle and tendency,using the triangle transform.Testing on both synthetic and field data showed that the proposed method estimated the formation dip angle and thickness with an accuracy higher than that of the common plane-wave decomposition dip calculation method.This indicates its practical application value and potential for widespread use.

dip logging data,plane-wave decomposition,logging and seismic integration,model regularization,formation dip,dip logging,filter

2016-09-14;改回日期2017-07-27。

穆立华(1962—),男,博士,教授级高级工程师,主要从事油气地质研究与管理工作。

P631

A

1000-1441(2017)06-0820-07

10.3969/j.issn.1000-1441.2017.06.007