張　華，龐燕青，賀振華，楊　波，何光明，巫芙蓉，陳愛萍 ，謝小翠

(1.成都理工大學(xué) 油氣藏地質(zhì)及開發(fā)工程國家重點實驗室，四川 成都 610059；　2.中國石油 川慶鉆探工程有限公司地球物理勘探公司，四川 成都 610213；　3.中國石油 物探技術(shù)試驗基地 山地地震技術(shù)試驗基地，四川 成都 610213；4.中國石油 新疆油田分公司，新疆 克拉瑪依 834000；　5. 中國石油 新疆油田分公司 陸梁油田作業(yè)區(qū)，新疆 克拉瑪依，834000)

?

高階統(tǒng)計約束稀疏盲反褶積技術(shù)研究與應(yīng)用

張華1,2，3，龐燕青4，賀振華1，楊波5，何光明2,3，巫芙蓉2，3，陳愛萍2,3，謝小翠2,3

(1.成都理工大學(xué) 油氣藏地質(zhì)及開發(fā)工程國家重點實驗室，四川 成都 610059；2.中國石油 川慶鉆探工程有限公司地球物理勘探公司，四川 成都 610213；3.中國石油 物探技術(shù)試驗基地 山地地震技術(shù)試驗基地，四川 成都 610213；4.中國石油 新疆油田分公司，新疆 克拉瑪依 834000；5. 中國石油 新疆油田分公司 陸梁油田作業(yè)區(qū)，新疆 克拉瑪依，834000)

以貝葉斯基本理論為基礎(chǔ)，利用高階雙譜域估算的高精度地震子波為約束，在修正柯西判別準(zhǔn)則下，實現(xiàn)地震數(shù)據(jù)的稀疏盲反褶積拓頻處理，獲取比常規(guī)反褶積提頻方法分辨率更高的可靠數(shù)據(jù)。模擬數(shù)值實驗和多波實際數(shù)據(jù)處理結(jié)果表明：雙譜域估算的子波精度較高，以該子波為先驗約束的稀疏盲反褶積能極大提高縱波和轉(zhuǎn)換波數(shù)據(jù)的分辨率，滿足精細(xì)解釋能分辨更小薄層的要求，可作為多波高分辨率處理的有效手段。

子波估計；提頻；反褶積；反射系數(shù)；地震勘探

反褶積是實際地震資料進(jìn)行提高分辨率處理的主要手段。隨著精細(xì)石油勘探的不斷深入，對高分辨率地震資料處理的要求也越來越高，而以高斯白噪和子波最小相位為假設(shè)的傳統(tǒng)反褶積方法無法滿足地震資料分辨更小薄層的需要。

自Wiggins (1978)[1]提出了稀疏度準(zhǔn)則——最小熵反褶積以來，隨后出現(xiàn)了大量關(guān)于反演子波和反射系數(shù)的一大類盲反褶積技術(shù)。Canadas (2002)[2]在盲反褶積框架下，實現(xiàn)了地震反射系數(shù)與地震子波同時反演算法；劉財(2008)[3]等人利用自然梯度盲反褶積法估算地震子波,得到更接近于實際情況的混合相位子波,用于反褶積后的地震資料分辨率得到了提高；孟大江(2009)[4]對稀疏反褶積和高精度地震子波提取進(jìn)行了深入研究；李力民等(2011)[5]采用預(yù)條件共軛梯度迭代反演地震子波和反射系數(shù)的方法。盡管這些方法反演準(zhǔn)則、求解算法存在差異，但理論架構(gòu)基本相似，都不需要白噪反射系數(shù)和最小相位假設(shè)，以通過對地震記錄的子波反演后，再反演出稀疏反射系數(shù)，達(dá)到盲反褶積提頻的目的，但是反射系數(shù)的求取很大程度依賴子波的精度和效率，如何高效、準(zhǔn)確地估算出地震子波從而提高盲反褶積過程中反射系數(shù)估算精度和可信度，是在實際資料應(yīng)用一個難題。

以高階累積量和高階譜為基礎(chǔ)的理論從20世紀(jì)80年代發(fā)展以來，無論是從頻率域還是從時間域，基于高階累計量理論的子波估計為我們提供了新的思路，利用高階累計量能同時估算子波的振幅信息和相位信息，通過近幾年大量的學(xué)者研究，能從多解的相位中估算出精度較高的地震子波[6-7]。

本文在前人研究的基礎(chǔ)上，提出一種高階統(tǒng)計約束稀疏盲反褶積技術(shù)，把在高階雙譜域中進(jìn)行多相位疊加估算出的高精度地震子波，作為貝葉斯稀疏盲反褶積的高階約束，進(jìn)行地震反射系數(shù)的高精度估算，從而在多波高分辨率處理資料中達(dá)到能分辨薄層的目的，為后續(xù)多波聯(lián)合解釋與反演[8-10]提供高質(zhì)量的資料。

1　雙譜域多相位疊加估計子波

利用地震資料的雙譜進(jìn)行地震子波估計，是地震資料處理中行之有效的方法。用褶積模型表示地震記錄可為：

(1)

式中：s(t)，w(t)，r(t)和n(t)分別為地震記錄數(shù)據(jù)、地震子波數(shù)據(jù)、反射系數(shù)數(shù)據(jù)和噪聲數(shù)據(jù)；t為時間，s；*為褶積算子。

根據(jù)高階累積量的定義和傅里葉變換公式得到雙譜域地震子波估計表達(dá)式：

(2)

式中，Bs(ω1,ω2)為地震記錄s(t)的雙譜域數(shù)據(jù)；W(ω)為地震子波w(t)的傅里葉譜域數(shù)據(jù)；H為共軛算子；ω表示角頻率，rad/s；γ3r是反射系數(shù)的三階累計量,常數(shù)。

根據(jù)公式(2)可得到子波相位譜估計表達(dá)式：

(3)

式中：ψ(ω1,ω2)為子波雙譜域相位譜數(shù)據(jù)；φ(ω)為子波的傅里葉域相位譜數(shù)據(jù)。

子波相位估算方法很多，如Brillinger[11]提出的相位重構(gòu)法以及Rosenblatt[12]等提出的BMU遞歸算法等，但該遞歸算法忽略了雙譜相位的平均，對初值依賴較大，不利于噪聲的壓制。

為了更好的解決雙譜域中多值多解和抗噪性問題，提高相位估算精度，必須充分考慮雙譜的平均特性，這里采用雙譜域多相位疊加估算子波相位的辦法。

以雙譜域中兩個角頻率ω1和ω2的相位值為縱橫坐標(biāo)建立坐標(biāo)系。為了便于統(tǒng)一計算，可把ψ(ω1,ω2)寫成φ(ω1,ω2)形式，在第一象限取值，且奇偶分開，則公式(3)可表示為：

(4)

式中：M[]表示相位主值算子，將[]中的值限定在模為π的(-π,π]之內(nèi)。

以此類推，橫坐標(biāo)相位上每個相位的相位值可通過公式(6)得到：

(6)

由上式估算出的相位沒有受到模2π的影響，考慮雙譜平均性，將他們疊加求平均不會引起相位多解性問題，得到最終的地震子波相位。

2　高階統(tǒng)計約束的稀疏盲反褶積

基于貝葉斯理論的稀疏盲反褶積提頻技術(shù)，由于其采用反演的手段重構(gòu)高分辨率地震剖面，不受傳統(tǒng)反褶積提頻頻帶的限制的要求，可以最大程度的提高資料分辨率，滿足多波等資料分辨更小薄層的要求。但該方法在求解目標(biāo)函數(shù)時，需要通過多次循環(huán)迭代同時求取未知的地震子波和反射系數(shù)，在實際生產(chǎn)中效率過低，精度也不高，如果在雙譜域中先求出高精度的地震子波，再反演目標(biāo)函數(shù)求取反射系數(shù)，勢必提高其處理效率的同時，得到的精度更高的反射系數(shù)。

根據(jù)褶積模型和貝葉斯理論公式可得近似式[13]：

(7)

(8)

(9)

對上式(9)取對數(shù)得式(10)：

(10)

(11)

為了能在提高分辨率的同時，有效保護(hù)數(shù)據(jù)的弱反射信號，需要對目標(biāo)函數(shù)進(jìn)行修正的柯西準(zhǔn)則判斷，即：

(12)

把雙譜域估計的地震子波矩陣W=[w0,…,wp]T帶入目標(biāo)函數(shù)(11)中[14]，那么對其求取反射系數(shù)r的偏導(dǎo)數(shù)：

(13)

圖1　理論模型數(shù)值試驗Fig.1　Numerical test of theoretical modela.理論實驗子波;b.理論實驗反射系數(shù);c.褶積加噪合成記錄;d.估算子波(紅色)與理論子波(藍(lán)色)對比; e.脈沖反褶積反射系數(shù),紅色箭頭為反射系數(shù)位置;f.本文方法反射系數(shù)，紅色箭頭為反射系數(shù)位置

圖3　 轉(zhuǎn)換波數(shù)據(jù)處理前、后效果圖對比Fig.3　Comparison of PS-wave data before and after processinga. 處理前轉(zhuǎn)換波剖面(左)和振幅譜(右);b. 常規(guī)譜白化提頻后(左)和本文方法提頻后(右)轉(zhuǎn)換波剖面; c. 常規(guī)譜白化提頻后振幅譜(左)和本文方法提頻后(右)振幅譜;d. 處理前(左)、常規(guī)方法(中)和本文方法(右)轉(zhuǎn)換波儲層段最大波峰振幅

(14)

(15)

3　數(shù)值實驗及實際資料應(yīng)用

3.1理論模型數(shù)值試驗

圖1c為圖1a理論子波與圖1b理論反射系數(shù)褶積并加噪聲的合成記錄，圖1d為從圖1c估算的地震子波與圖1a理論子波的對比圖。從圖中可知：雙譜域估算的子波與理論子波相似度高，尤其是子波的主峰和旁瓣與理論子波吻合度較高。圖1f為利用本文方法得到圖1c的反射系數(shù)，該反射系數(shù)與圖1b的理論反射系數(shù)相似度較高，且與圖1e脈沖反褶積得到的反射系數(shù)相比更加可靠、準(zhǔn)確，證明了該方法的正確性。

3.2實際資料應(yīng)用

為驗證上述理論和方法的實際應(yīng)用效果，選取四川盆地某地區(qū)多波資料目的層進(jìn)行應(yīng)用測試。

圖2為縱波數(shù)據(jù)處理前后效果圖對比。圖2a與圖2b，2c比較可知：常規(guī)譜白化和本文方法都大幅提高了剖面的分辨率，構(gòu)造細(xì)節(jié)凸顯，而且頻譜也比處理前所有展寬。進(jìn)一步比較圖2b和圖2c可知：本文方法對提高剖面深層同相軸的分辨率效果更好，與測井合成記錄的吻合度更高，而且有效拓展高頻振幅譜的同時保留了低頻段信息。

圖3為轉(zhuǎn)換波數(shù)據(jù)處理前后效果圖對比。2b與圖2b，2c比較知：處理前轉(zhuǎn)換波分辨率較低，采用常規(guī)譜白化處理后同相軸連續(xù)性變差，斷層發(fā)育被掩蓋，而利用本文方法不但大幅提高分辨率，而且同相軸連續(xù)性較好，斷層清晰可見，與測井合成記錄吻合度高，同時本文方法處理后的振幅譜比常規(guī)方法保留了更加豐富的低頻信息和高頻拓展信息。

多波多分量勘探與單一縱波勘探相比，它所能提供的地震屬性信息將會成倍的增加，并能衍生出各種組合參數(shù)。利用這些參數(shù)估算地層巖性與流體等，將比僅依靠縱波的可靠性更高，可以最大限度消除利用單純縱波進(jìn)行儲層預(yù)測的不唯一性。

在對縱波和轉(zhuǎn)換波分量資料分別進(jìn)行多種屬性實驗中，相對而言，振幅類屬性較為穩(wěn)健，因而計算了縱波和轉(zhuǎn)換波最大波峰振幅類屬性。從圖2中縱波剖面與圖3中轉(zhuǎn)換波剖面中儲層段T3x2—T3x1中進(jìn)行最大波峰振幅分布圖分析知：本文方法處理后最大波峰振幅能量更聚焦、收斂性更好，而且可見氣井位于縱波中強(qiáng)振幅能量區(qū)(黑圈處)，并同時位于轉(zhuǎn)換波弱振幅區(qū)(黑圈處)，這與儲層相面含氣預(yù)測結(jié)果吻合：縱波波為中強(qiáng)反射，轉(zhuǎn)換波為弱反射，縱波強(qiáng)于轉(zhuǎn)換波反射。

4　結(jié)論

1) 本文提出一種高階統(tǒng)計約束稀疏盲反褶積技術(shù)，在高階累計量域中利用多相位疊加，提取的混合相位子波精度高，通過子波約束稀疏盲反褶積求取的反射系數(shù)精度高、強(qiáng)弱能量分明。

2) 利用本文方法在實際多波資料處理中，提高分辨率明顯，尤其是低分辨的轉(zhuǎn)換波層序內(nèi)幕更加明顯，斷層清晰，與測井有較好的吻合度，同相軸強(qiáng)弱能量更加分明，為后續(xù)研究地層巖石解釋提供更好的地震資料。

[1]Wiggings R A.Minimum entropy deconvolution[J].Geoexploration,1978,16(1):21-35.

[2]Canadas G.A mathematical framework for blind deconvolution inversion problems[J].The 72nd SEG Annual Meeting,2002,2202-2205.

[3]劉財,謝金娥,郭全仕,等.一種基于自然梯度的地震盲反褶積方法[J].石油物探，2008，47(5):439-443.

Liu Cai,Xie Jin'e ,Guo Quanshi,et al.Seismic blind deconvolution based on natural gradient[J].Geophysical Prospecting Petroleum,2008，47(5):439-443.

[4]孟大江.稀疏變換反褶積和高精度子波提取方法研究[D].吉林大學(xué)，2009.

Meng Dajiang.Deconvolution based on sparse transform and high precision wavelet extraction methods[D].Jilin University,2009.

[5]李力民,姚泳汐.基于Cauchy 稀疏約束的同時迭代反演反射系數(shù)與地震子波的研究[J].重慶科技學(xué)院學(xué)報(自然科學(xué)版),2011，1673-1980.

LI Limin,YAO Yongxi.Iteration Inversion of Reflectivity Coefficient and Seimic Wavelet Based on the Cauchy Sparse Constraint[J].Journal of Chongqing University of Science and Technology(Natural Sciences Edition),2011，1673-1980.

[6]李大衛(wèi),田繼東,成竹,等.基于高階統(tǒng)計的混合相位子波估計及其處理[J].勘探地球物理進(jìn)展，2005，1671-8585.

Li Dawei,Tian Jidong,Cheng Zhu,et al.Estimation of mixed-phase wavetet base on higher-order statistics[J].Progress in Exploration Geophysics ,2005，1671-8585.

[7]陳建,戴永壽,張亞南,等.基于高階統(tǒng)計量的地震子波提取方法評價[J].石油地球物理勘探，2013,48(3):497-503.

Chen Jian,Dai Yongshou,Zhang Yanan.Evaluation approaches for wavelet pickup based on high-order statistics[J].Oil Geophysical Prospecting,2013,48(3):497-503.

[8]楊克明,張虹.地震三維三分量技術(shù)在致密砂巖裂縫預(yù)測中的應(yīng)用-以川西新場氣田為例[J].石油與天然氣地質(zhì),2008,29(5):683-689.

Yang Keming,Zhang Hong.Application of 3D-3C seismic technique to the prediction of fractures in tight sandstone-an example from Xinchang gas field in the Western Sichuan Depression[J].Oil & Gas Geology,2008,29(5):683-689.

[9]葉泰然,付順,呂其彪,等.多波地震聯(lián)合反演預(yù)測相對優(yōu)質(zhì)儲層——以川西深層致密碎屑砂巖為例[J].石油與天然氣地質(zhì),2009,30(3):357-362,369.

Ye Tairan,Fu Shun,Lü Qibiao,et al.Application of multiwave joint inversion to the prediction of relatively high-quality reservoirs—an example from the prediction of deep tight clastic sandstone reservoirs in the western Sichuan Basin[J].Oil & Gas Geology,2009,30(3):357-362,369.

[10]袁書坤,陳開遠(yuǎn),Bob A Hardage,等.縱波-轉(zhuǎn)換橫波聯(lián)合反演技術(shù)在頁巖氣藏勘探中的應(yīng)用[J].石油與天然氣地質(zhì),2015,36(1):168-174.

Yuan Shukun,Chen Kaiyuan,Bob A Hardage,et al.Application of joint PP and PS inversion to the exploration of shale-gas reservoirs[J].Oil & Gas Geology,2015,36(1):168-174.

[11]Brillinger D.R.An introduction topolyspectra[J].Ann Math Statist,1965,36(1):1351-1374.

[12]Rosenblatt M,Van Ness J W.Estimation of the bispectrum[J].Ann.Math.Statist,1965,36 (1):1120-1136.

[13]李國福.盲源反褶積方法與應(yīng)用研究[J].石油化工應(yīng)用，2010,1673-5285.

LI Guofu.The method and applied research of blind deconvolution[J].Petrochemical Industry Application,2010,1673-5285.

[14]孫建業(yè),高偉,劉懷山,等.基于高階統(tǒng)計約束實現(xiàn)同態(tài)盲反褶積[J].石油地球物理勘探,2013,48(1):31-36.

Sun Jianye,Gao Wei,Liu Huaishan,et al.Seismic homomorphic blind deconcolution based on high order statistics constraint[J].Oil Geophysical Prospecting,2013,48(1):31-36.

(編輯張玉銀)

Study and application of high-order statistics constrained sparse blind deconvolution technique

Zhang Hua1，2，3,Pang Yanqing4，He Zhenhua1,Yang Bo5,He Guangming2，3,Wu Furong2,3,Chen Aiping2，3,Xie Xiaocui2,3

(1.StateKeyLaboratoryofOilandGasReservoirGeologyandExploitation,ChengduUniversityofTechnology,Chengdu,Sichuan610059，China；2.GeophysicalExplorationCompany,ChuanqingDrillingEngineeringCo.Ltd.,CNPC,Chengdu，Sichuan610213，China；3.MountainGeophysicalTechnologyTestCenter,CNPC,Chengdu,Sichuan610213，China；4.XinjiangOilfieldCompany,CNPC,Karamay,Xinjiang834000，China;5.LuliangFieldOperationDistrict,XinjiangOilfieledCompany,Karamay,Xinjiang834000，China)

Based on the basic theory of Bayesian and using the high precision wavelet estimated in high-order double spectrum domain as constraints,this paper realized sparse blind deconvolution frequency-expanding processing under the modified Cauchy criteria and obtained data that are more reliable and have higher frequency than the conventional deconvolution frequency-enhancing method.The results of numerical experiments and real data processing show that the wavelet estimated in the dual spectrum domain has higher accuracy,under the priori constraint of which the sparse blind deconvolution can greatly enhance the resolution of both P-wave and PS-wave.It can meet the requirement of resolute very thin layers,thus can be used as an effective method for high resolution processing.

wavelet estimation,frequency improving,deconvolution,reflectance,seismic exploration

2015-03-02;

2016-07-10。

張華(1981—)，男，博士研究生、工程師，地球物理資料處理和反演解釋。E-mail:zhanghua_sc@cnpc.com.cn。

中國石油天然氣集團(tuán)公司科學(xué)研究與技術(shù)開發(fā)項目(2013E-38-08)。

0253-9985(2016)04-0591-07

10.11743/ogg20160417

P631

A