2. 中国石油天然气集团公司油藏描述重点实验室, 兰州 730020
2. CNPC Reservoir Description Key Laboratory, Lanzhou 730020, China
高分辨率地震资料处理是地震精细解释的重要一环,现有的处理方法包括反褶积、谱白化和反Q滤波等,不同的方法具有各自的特点和优势,从不同的角度解决了地震资料处理中遇到的一些难题.随着岩性油气藏勘探及开发的不断深入,薄储层的问题日益凸显,亟需高精度的地震资料来进行描述和解释,而传统的方法暴露出在某些方面的局限性,迫切需要继续探索新的方法.
传统的反褶积方法假设地震信号是平稳信号、子波为最小相位及反射系数为白噪序列.但一般的地震资料难以满足以上条件.基于此,Rosa等首先提出了谱模拟反褶积方法,它对非白噪序列的反射系数具有很好的包容性,有效降低了反射系数非白噪成分对子波振幅谱求取的影响.国内学者在子波振幅谱的求取、反射序列特征及自适应方面进行了研究(赵波等,1996;唐博文等,2010;李振春等,2013).Rosa等人的方法对子波振幅谱的形态有明确限定;而反射系数的长度有限且非白噪声的特点导致二阶谱并不能完全分离子波成分,恢复的子波振幅谱在某些频段不准确,局限了方法的应用.
微分可以凸显信号的变化率,微分算子在频率域具有线性的单调递增特征,可以提高信号的高频成分,同时也压制了信号的低频成分.基于微分的这种特性,将地震记录作了少数几阶微分变换,在同相位条件下进行相加,获得了宽频高分辨率的地震记录,并且与谱白化算法和Gabor反褶积进行了对比,显示了此方法的良好效果(Sajid and Ghosh, 2014).从频谱上分析,Sajid的方法是对不同阶微分信号的振幅谱(信号微分谱)进行相加,从而拓宽了频带,但这种信号微分谱的直接相加缺乏依据,结果有待商榷..
本文提出了一种基于二阶谱及多阶微分融合的频谱拓展方法.这种方法不需要精确地求取子波的振幅谱,而通过不同阶微分信号的振幅谱来拟合期望记录的振幅谱,从而拓宽地震资料的频带.首先对地震记录作零相位化处理和有色补偿,然后求取不同阶微分信号的振幅谱(信号微分谱);在假设地震子波平滑的前提下,根据二阶谱的特性提取不同阶信号微分谱的包络,并与期望子波振幅谱建立线性方程组,采用L1正则化方法建立目标函数并求解,将所求系数用于信号微分谱的融合,得到宽频的地震资料.本方法对高频信号的恢复比传统的谱模拟反褶积方法更加准确.
1 方法原理 1.1 地震信号的二阶谱特征二阶谱算法原理简述如下:计算信号的振幅谱,并将结果作为输入信号再求取其振幅谱,即得到信号的二阶谱.
如果信号是白噪声,信号的功率谱密度在整个频域轴上为非零常数,即:
(1) |
则信号的二阶谱为单位脉冲响应πN0δ(ω),即信号能量在二阶谱上集中于零点附近.数值计算使用的是白噪声的一个样本,样本的数据越长越能体现白噪声的性质.
如果地震反射系数是白噪声序列,其二阶谱具有相同性质.但反射系数序列的频谱并非白谱,而是具有低频弱、高频强趋势的蓝谱,为此,反射系数系列首先要进行有色补偿.反射系数序列是白噪成分序列r w和非白噪成分序列r nw的褶积,公式为
(2) |
Walden等通过多次实验证实:一阶的ARMA模型(自回归滑动平均模型)可以较好地模拟反射系数的非白噪成分r nw,公式为
(3) |
其中L为滞后算子,ϕ1为AR参数,θ1为MA参数.在已知模型阶数的情况下,参数ϕ1和θ1可以通过矩估计法从声波测井资料中求取.反射系数白噪成分的振幅谱为
(4) |
其中,R为反射系数的振幅谱,R w为反射系数白噪成分的振幅谱,R nw为反射系数非白噪成分的振幅谱.
合成地震记录由反射系数与子波褶积所得,其频谱及二阶谱为
(5a) |
(5b) |
(5c) |
其中,S为地震记录s的振幅谱,W为子波w的振幅谱,S2、R w2、W 2分别表示信号s、r w、w的二阶谱.由二阶谱的性质可知:当反射系数序列足够长时,R w2近似为单位脉冲序列;子波w的振幅谱W为平滑序列(低频主导),即子波在二阶谱上位于低频端,由公式5c知:地震记录的低频成分基本反映了子波二阶谱的特征(如图 1a).当已知地震记录s时,可通过选取二阶谱S2的低值部分来恢复子波w的振幅谱(如图 1b).
具有平滑振幅谱的信号具有同样的特征,也可通过二阶谱的性质来恢复信号.由图 1a可知地震记录的二阶谱与子波的二阶谱并不完全重合,导致恢复的子波振幅谱在低幅值处求取不准确(箭头所指),给传统谱模拟反褶积求取子波振幅谱带来影响,但对下文中多阶微分融合的影响不大.
1.2 多阶微分融合微分表达信号的变化率,突出其细节信息;在频率域,微分虽不能拓宽频带,但可以提高主频,即具有分频属性.本文主要采用微分的这种性质来提取信号不同频率段信息,通过融合获得真实的宽频地震资料.为了准确利用二阶谱的性质,首先对信号进行有色补偿.
信号s(t)的n阶微分的振幅谱为
(6) |
其中n为微分阶数,S(f)为信号的振幅谱,S (n)(f)为n阶微分信号的振幅谱(简称信号的微分谱),n可以为任意实数.微分算子的频域表达式(2πif)n是关于f的线性函数,随着微分阶数n的增加,信号的高频逐渐被提升,而低频逐渐被压制.地震记录为离散信号,离散的微分算子为a t=(1, -1),则离散地震记录的n阶微分振幅谱为
(7) |
其中Δt为采样间隔,A(f)为微分算子的频谱.离散微分算子在频率域并非严格线性变化(如图 2),其低频近似线性变化,而高频趋于平缓.离散微分算子的这种性质避免了过高地提升高频噪声.
微分谱(公式7所求)为微分谱包络D n(低频主导)与反射系数谱R的乘积(公式8),且已知有色补偿后地震反射系数为白噪声序列,可根据2.1中二阶谱的性质来恢复微分谱包络D n.公式(8)为
(8) |
其中,S (n)为信号s(t)的n阶微分谱,S为地震记录振幅谱,W为地震子波振幅谱,R为反射系数振幅谱.
本文的目标是通过不同阶微分谱的融合来拟合出期望地震记录的振幅谱,相应地,如图 3所示,微分谱的包络(红色实线)须通过加权相加拟合为期望子波的振幅谱(红色虚线).对两者建立线性方程(微分阶数在频率域为实数,为了防止方程欠定,加入阶数间隔Δ,记:D n·Δ(f)=(2πf)n·Δ· W(f), n=1, 2…m),即
(9) |
可简写为
(10) |
其中D为多阶微分谱包络的矩阵,W q是期望子波振幅谱.为了简化运算的复杂性,采用L1正则化算法进行求解,该方法可以表达为Tikhonow正则化求泛函极小,建立目标函数为
(11) |
其中,λ为正则化因子,丨α丨为正则化约束.可求取拟合系数
设R=diag (R(f1), R(f2), …, R(fm)),其中序列R(f1), R(f2), …, R(fm)为反射系数振幅谱.由公式9得:
(12) |
上式可简写为
(13) |
其中,RD为多阶微分谱矩阵(公式8),F=RWq为期望记录的振幅谱.已知多阶微分谱矩阵及拟合系数
(14) |
颜色特征为反射系数的性质,所求得叠合谱需去掉颜色补偿.叠合谱与原始记录的相位谱联立得到新的频谱,反傅氏变换得到宽频地震记录.
2 模型验证建立地质模型(如图 4a),其中白色为泥岩,黄色为砂岩,假设地震波在两介质中的传播速度相同,均为4000 m/s,砂泥岩密度比为2:1;模型中不同地层厚度由深浅两端向中心逐渐变薄,分别为36 m,34 m,…,2 m,4 m,…,36 m,为便于与合成地震记录对比,地质模型的纵轴单位由深度(m)转化为双程走时间(ms).从地质模型中提取反射系数序列,与主频为50 Hz的雷克子波褶积得到合成地震记录(图 4b).为了测试噪声对处理结果的影响,在地震记录中加入S/N=4的随机噪声.
图 4c为采用传统谱模拟反褶积方法对合成地震记录处理后的剖面,此剖面在中心薄层处要比合成地震记录的分辨能力强,但振幅强度较弱;而且,某些道的噪声影响了地震波形(箭头所指的第4道),且反映厚层的同相轴响应信息较弱(横向箭头所指).图 4d为采用多阶微分融合方法处理后的剖面,此方法比传统谱模拟反褶积方法处理后资料分辨能力高,中心处可以分辨3 m左右的薄层;其次,反映厚层同相轴的响应信息增强;而且噪声对于整个剖面的影响不大.
图 5a为反射系数的振幅谱,图 5b为合成地震记录振幅谱,图 5c为传统谱模拟反褶积方法处理后的地震记录振幅谱,其频带得到拓展,但并不充分,而且某些高频段有异常抬升现象,图 5d是通过多阶微分融合得到的地震记录振幅谱,高频得到充分拓展,而且变化趋势更符合期望子波的振幅谱.传统谱模拟反褶积方法求取的子波振幅谱不准确(主要是低振幅段),存在剩余子波,导致地震资料的某些频段(特别是高频信息)补偿过度或不足.本文方法避免了这种缺陷,不仅保留了地震记录频率的局部变化特征,而且其各频率段的比例更加合理和充分,对降低噪声干扰也有一定优势.
实际资料采用柴达木盆地东坪三维区抽取的一条二维测线,工区内砂泥岩互层较多且比较薄,叠后地震数据经过了常规的叠前处理、叠前时间偏移及子波零相位化处理.子波零相位化消除了子波在叠后地震记录相位谱中的影响,在此基础上,多阶微分融合的方法对地震记录的振幅谱进行了有效拓展.实际地震资料中存在噪声,特别是高频段的随机噪声,所以期望子波要选择恰当(鉴于所选工区地震资料的信噪比相对较高,选择主频为35Hz的雷克子波作为其期望子波).
图 6a为原始的叠后地震记录剖面,分辨率明显不够,有些同相轴相互交错;图 6b为传统的谱模拟反褶积方法处理后的剖面,其分辨率得到了一定的提升;图 6c为采用多阶微分融合的方法处理后的剖面,与传统的方法对比,其分辨率明显提高(如图中箭头所指);而且同相轴的连续性更好,能量比较集中;噪声对图 6c剖面的影响要比图 6b剖面小;同时多阶微分融合的方法基本维持了原始叠后剖面的能量分布特征.处理的整体质量的提升得益于本文新方法对于地震记录频谱的合理变换,特别是对于高频信息的提取更加准确.
图 7是同一工区的另一剖面,原始叠加剖面(图 7a)分辨能力不够,同相轴连续性差;经过传统的谱模拟反褶积方法及多阶微分融合的方法处理后,剖面与合成地震记录都具有较好的标定,但在红线所圈区域,后者的分辨能力及同相轴连续性要优于前者.
图 8是图 7剖面经过多阶微分融合方法处理前后的振幅谱变化.处理后,信号频谱的高截频从45 Hz左右拓展到了70 Hz左右;而且因所选期望子波较理想,低频成分衰减较小,削弱了子波旁瓣的影响;同时,处理前信号振幅谱的细节变化得到了保留.
针对目前薄互层识别和分辨的难题,提出一种基于二阶谱及多阶微分融合的频谱拓展方法,经过研究得到以下结论:
(1)传统的谱模拟反褶积方法提取的子波振幅谱都存在一定的误差,这种误差严重影响了地震薄层响应的高频成分,造成对薄层的成像不准确,本文提出的方法不需要对子波振幅谱精确求取,通过多阶微分谱融合提取的高频信息更加准确.
(2)多阶微分谱融合过程采用L1正则化方法建立目标函数并求解,减小了运算的复杂程度,并获得了最优的拟合系数.
(3)地震记录振幅谱拓展后需要去掉有色补偿,还原原始反射系数序列.通过对模型数据及实际资料采用传统谱模拟反褶积方法和多阶微分融合的方法分别处理,证实了后者的优势及有效性.
致谢 感谢审稿专家提出的修改意见和编辑部的大力支持![] | Cai L F, Tian X M .2012. A blind seismic deconvolution method based on Particle Swarm Optimization[J]. Progress in Geophysics , 27 (3) : 1116–1122. DOI:10.6038/j.issn.1004-2903.2012.03.036 |
[] | Cai W T, Hu G Y, Fan T E, et al .2014. Quality factor extraction by spectral simulation method[J]. Progress in Geophysics , 29 (2) : 642–649. DOI:10.6038/pg20140223 |
[] | Chen S Q, Li X Y .2015. Seismic resolution enhancement based on the scale characteristics of Fourier transform[J]. Oil Geophysical Prospecting , 50 (2) : 213–218. |
[] | Ding J J, Dai Y S, Zhang Y N, et al .2013. Improving seismic data resolution based on high frequency compensation[J]. Progress in Geophysics , 28 (6) : 3214–3221. DOI:10.6038/pg20130648 |
[] | Dong E Q, Liu G Z, Zhang Z P .2001. Fast implementation technique of adaptive Kalman filtering deconvolution via dyadic wavelet transform[J]. Chinese Journal. of Geophysics. , 44 (2) : 255–262. DOI:10.3321/j.issn:0001-5733.2001.02.013 |
[] | Fan X D, Jing H, Liu Y Q .1995. Time-space variant spectrum whitening of seimic data[J]. Oil Geophysical Prospecting , 30 (4) : 550–555. |
[] | Gao J H, Zhang M, Zhang B .2011. An improved blind deconvolution based on mutual information rate[J]. Progress in Geophysics , 26 (2) : 540–548. DOI:10.3969/j.issn.1004-2903.2011.02.019 |
[] | Guo X, Yong X S, Gao J H, et al .2016. Seismic bandwidth broadening method based on multi-order differential fusion in frequency domain[J]. Geophysical Prospecting for Petroleum , 55 (2) : 271–279. |
[] | Hargreaves N, Treitel S, Smith M. 2013. Frequency extension, resolution, and sparse inversion[C].//83th Annual International Meeting. SEG, Expanded Abstracts, 3345-3349. |
[] | Jia L H, Wu C J, Luo Y X, et al .2002. High resolution seismic data processing techniques[J]. Geophysical Prospecting for Petroleum , 41 (4) : 484–488. |
[] | Li F L, Wang C S, Xu B X, et al .1990. LP norm deconvolution and application to seismic prospecting[J]. Acta Geophysica Sinica , 33 (5) : 585–592. |
[] | Li J S, Jiang X P, Wang X T .2011. The main energy spiking deconvolution[J]. Periodical of Ocean University of China , 41 (7-8) : 119–126. |
[] | Li Z C, Li D, Wang D Y, et al .2013. Adaptive spectral modeling deconvolution based on the constraint of signal to noise ratio spectrum[J]. Progress in Geophysics , 28 (1) : 0301–0309. DOI:10.6038/pg20130133 |
[] | Liu Z W, Wang Y C, Zhao H X, et al .2013. High-resolution processing methods of thin interbeds imaging[J]. Chinese Journal. of Geophysics. , 56 (4) : 1350–1359. DOI:10.6038/cjg20130429 |
[] | Rosa R, Ulrych T J .1991. Processing via spectral modeling[J]. Geophysics, 56 (8) : 1244–1251. DOI:10.1190/1.1443144 |
[] | Sajid M, Ghosh D .2014. A fast and simple method of spectral enhancement[J]. Geophysics, 79 (3) : V75–V80. DOI:10.1190/geo2013-0179.1 |
[] | Shang X M, Diao R, Feng Y P, et al .2014. The application of spectral modeling method to high resolution processing of seismic data[J]. Geophysical & Geochemical Exploration , 38 (1) : 75–80. |
[] | Sun C Y .2000. Spectrum modeling method and its application to seismic resolution improvement[J]. Oil Geophysical Prospecting , 35 (1) : 27–35. |
[] | Tang B W, Zhao B, and Wu Y H, et al .2010. A new way to realize spectral modeling deconvolution[J]. Oil Geophysical Prospecting , 45 (SZ1) : 66–70. |
[] | Ten Kroode F, Bergler S, Corsten C, et al .2013. Broadband seismic data-the importance of low frequencies[J]. Geophysics, 78 (2) : WA3–WA14. DOI:10.1190/geo2012-0294.1 |
[] | Walden A T, Hosken J W J .1985. An investigation of the spectral properties of primary reflection coefficients[J]. Geophysical Prospecting, 33 (3) : 400–435. DOI:10.1111/gpr.1985.33.issue-3 |
[] | Wan H, Fan X Y, Liu T, et al .2012. Methods and applications for improving pre-stack seismic data resolution[J]. Progress in Geophysics , 27 (1) : 304–311. DOI:10.6038/j.issn.1004-2903.2012.01.034 |
[] | Wang B L, Yin X Z, Ding L X, et al .2015. Study of fast stochastic inversion based on FFT-MA spectrum simulation[J]. Chinese Journal. of Geophysics. , 58 (2) : 664–673. DOI:10.6038/cjg20150227 |
[] | Wang D W, Liu Z, Chen X H, et al .2009. Discussion on characteristics of time-lapse seismic differential waveform of thin oil layer[J]. Progress in Geophysics , 24 (1) : 170–175. |
[] | Yun M H, Ji Z Y, Cui J T .1997. Can derivation improve the resolution of seismogram?[J]. Oil Geophysical Prospecting , 32 (5) : 740–748. |
[] | Yun M H, Tian J W, Zhang G F .1998. Theoretical study on multiple-order differential spectrum shaping filter[J]. Petroleum Instruments , 12 (6) : 22–24. |
[] | Zhang K, Li Y D, Liu G Z, et al .1999. Multiresolution seismic signal deconvolution[J]. Chinese Journal. of Geophysics. , 42 (4) : 529–535. DOI:10.3321/j.issn:0001-5733.1999.04.011 |
[] | Zhao B, Yu S P, Nie X B, et al .1996. Spectral-modeled deconvolution and its application[J]. Oil Geophysical Prospecting , 31 (1) : 101–116. |
[] | Zhao S L, Wang X M, Wang J H .1994. A simple frequency compensating method[J]. Oil Geophysical Prospecting , 29 (2) : 231–235. |
[] | 蔡连芳, 田学民.2012. 基于PSO的地震盲反褶积方法[J]. 地球物理学进展, 27 (3) : 1116–1122. DOI:10.6038/j.issn.1004-2903.2012.03.036 |
[] | 蔡文涛, 胡光义, 范廷恩, 等.2014. 谱模拟法提取地层品质因子[J]. 地球物理学进展, 29 (2) : 642–649. DOI:10.6038/pg20140223 |
[] | 陈双全, 李向阳.2015. 应用傅里叶尺度变换提高地震资料分辨率[J]. 石油地球物理勘探, 50 (2) : 213–218. |
[] | 丁进杰, 戴永寿, 张亚南, 等.2013. 基于高频补偿方法提高地震资料分辨率的初步研究[J]. 地球物理学进展, 28 (6) : 3214–3221. DOI:10.6038/pg20130648 |
[] | 董恩清, 刘贵忠, 张宗平.2001. 自适应Kalman滤波反褶积的快速实现方法[J]. 地球物理学报, 44 (2) : 255–262. DOI:10.3321/j.issn:0001-5733.2001.02.013 |
[] | 范小东, 景灏, 刘益强.1995. 地震资料时空变谱白化处理[J]. 石油地球物理勘探, 30 (4) : 550–555. |
[] | 高静怀, 张明, 张兵.2011. 一种改进的基于互信息率的盲反褶积方法[J]. 地球物理学进展, 26 (2) : 540–548. DOI:10.3969/j.issn.1004-2903.2011.02.019 |
[] | 郭欣, 雍学善, 高建虎, 等.2016. 基于频率域多阶微分融合的地震频带拓宽方法研究[J]. 石油物探, 55 (2) : 271–279. |
[] | 贾丽华, 吴长江, 罗焱鑫, 等.2002. 地震资料高分辨率处理技术[J]. 石油物探, 41 (4) : 484–488. |
[] | 李凤林, 王承曙, 徐伯勋, 等.1990. LP模反褶积及其在地震勘探中的应用[J]. 地球物理学报, 33 (5) : 585–592. |
[] | 李金山, 姜秀萍, 王修田.2011. 主能量脉冲反褶积[J]. 中国海洋大学学报, 41 (7-8) : 119–126. |
[] | 李振春, 李栋, 王德营, 等.2013. 信噪比谱约束的自适应谱模拟反褶积方法研究[J]. 地球物理学进展, 28 (1) : 301–309. DOI:10.6038/pg20130133 |
[] | 刘志伟, 王彦春, 赵会欣, 等.2013. 薄互层地震成像中高分辨率处理方法[J]. 地球物理学报, 56 (4) : 1350–1359. DOI:10.6038/cjg20130429 |
[] | 尚新民, 刁瑞, 冯玉苹, 等.2014. 谱模拟方法在高分辨率地震资料处理中的应用[J]. 物探与化探, 38 (1) : 75–80. |
[] | 孙成禹.2000. 谱模拟方法及其在提高地震资料分辨率中的应用[J]. 石油地球物理勘探, 35 (1) : 27–35. |
[] | 唐博文, 赵波, 吴艳辉, 等.2010. 一种实现谱模拟反褶积的新途径[J]. 石油地球物理勘探, 45 (S增刊1) : 66–70. |
[] | 万欢, 樊小意, 刘涛, 等.2012. 叠前地震资料提高分辨率处理方法及应用[J]. 地球物理学进展, 27 (1) : 304–311. DOI:10.6038/j.issn.1004-2903.2012.01.034 |
[] | 王保丽, 印兴耀, 丁龙翔, 等.2015. 基于FFT-MA谱模拟的快速随机反演方法研究[J]. 地球物理学报, 58 (2) : 664–673. DOI:10.6038/cjg20150227 |
[] | 王大伟, 刘震, 陈小宏, 等.2009. 薄油层时移地震差异波形特征探讨[J]. 地球物理学进展, 24 (1) : 170–175. |
[] | 云美厚, 纪阵宇, 崔景涛.1997. 求导能否提高地震记录的分辨率[J]. 石油地球物理勘探, 32 (5) : 740–748. |
[] | 云美厚, 田景文, 张国富.1998. 多级微分频谱整形滤波器的理论研究[J]. 石油仪器, 12 (6) : 22–24. |
[] | 章珂, 李衍达, 刘贵忠, 等.1999. 多分辨率地震信号反褶积[J]. 地球物理学报, 42 (4) : 529–535. DOI:10.3321/j.issn:0001-5733.1999.04.011 |
[] | 赵波, 俞寿朋, 聂勋碧, 等.1996. 谱模拟反褶积方法及其应用[J]. 石油地球物理勘探, 31 (1) : 101–116. |
[] | 赵圣亮, 王晓明, 王建华.1994. 一种简易的频率补偿方法[J]. 石油地球物理勘探, 29 (2) : 231–235. |