縱橫波彈性阻抗聯(lián)合反演方法

張廣智，杜炳毅，陳懷震，高建虎，李 超，李 遠(yuǎn)

1.中國石油大學(xué)（華東）地球科學(xué)與技術(shù)學(xué)院，山東 青島 266580

2.中國石油勘探開發(fā)研究院西北分院，蘭州 730020

3.東方地球物理公司新興物探開發(fā)處，河北 涿州 072750

0 引言

為了充分利用實(shí)際地震資料非零炮檢距的特點(diǎn)，近年來，許多學(xué)者對疊前地震反演做了大量的研究工作［1－5］。Connolly［6］根據(jù)聲阻抗的思想提出了彈性阻抗（elastic impedance）的概念，并且推導(dǎo)了彈性阻抗公式。Larsen［7］利用PP波和PS波疊前地震數(shù)據(jù)同步反演方法獲取了較好的縱波阻抗和橫波阻抗的估計(jì)值。Milos［8］將彈性阻抗反演應(yīng)用到實(shí)際資料中，充分利用小角度道集和大角度道集資料，提高了反演的可靠性。Duffaut等［9］在Connolly研究的基礎(chǔ)上提出橫波彈性阻抗的概念，將轉(zhuǎn)換波地震數(shù)據(jù)與測井?dāng)?shù)據(jù)聯(lián)系起來。由于彈性阻抗值隨著角度變化沒有統(tǒng)一的量綱，不利于不同角度的彈性阻抗進(jìn)行對比，Whitecombe［10］對彈性阻抗做了標(biāo)準(zhǔn)化處理。Gonzalez［11］將轉(zhuǎn)換波反射系數(shù)表示成縱波入射角的函數(shù)，進(jìn)而得到轉(zhuǎn)換波彈性阻抗公式，用不同角度轉(zhuǎn)換波彈性阻抗的交會來識別巖性和流體。馬勁風(fēng)［12］、王保麗等［13－15］對彈性阻抗反演都做了大量的研究工作。崔杰等［16］研究了各向異性介質(zhì)中的轉(zhuǎn)換波彈性阻抗，并且分析了各向異性參數(shù)對轉(zhuǎn)換波彈性阻抗的影響。陳天勝等［17］對縱橫波聯(lián)合反演方法做了相應(yīng)的研究。

筆者選取合適的轉(zhuǎn)換波彈性阻抗公式，對選取轉(zhuǎn)換波彈性阻抗進(jìn)行標(biāo)準(zhǔn)化處理；在此基礎(chǔ)上探索研究縱橫波彈性阻抗聯(lián)合反演方法，并且選用二維Marmousi模型和實(shí)際井模型，在合成角度道集添加信噪比（S/N）為3的隨機(jī)噪聲，對該反演算法進(jìn)行試算；最后從反演的彈性阻抗和轉(zhuǎn)換波彈性阻抗體中提取了精確的彈性參數(shù)。

1 縱橫波彈性阻抗聯(lián)合反演

1.1 基本原理

1.1.1 縱橫波彈性阻抗定義

Connolly［6］在1999年率先提出了彈性阻抗的概念，并在Aki－Richards近似式［2］的基礎(chǔ)上推導(dǎo)了彈性阻抗的表達(dá)式，將彈性阻抗EI表示為入射角θ、縱波速度vP、橫波速度vS及密度ρ的函數(shù)：

其中：

式中：vP1、vS1為分界面上部介質(zhì)的縱波速度和橫波速度；vP2、vS2為分界面下部介質(zhì)的縱波速度和橫波速度（圖1）。

圖1 彈性波在彈性界面上的反射與透射Fig.1 Reflection and transmission on elastic interface

Duffaut［9］在轉(zhuǎn)換波反射系數(shù)近似的基礎(chǔ)上，推導(dǎo)了轉(zhuǎn)換波彈性阻抗公式：

式中：

2006 年，Gonzalez［11］利用 Aki－Richards［2］近似將轉(zhuǎn)換波彈性阻抗公式表示成入射角θ的函數(shù)形式：

式中：

利用Goodway（1997）三層含氣砂巖與頁巖模型［1］（參數(shù)見表1），將由Zoeppritz方程得到的精確轉(zhuǎn)換波反射系數(shù)、由Aki－Richards近似得到的轉(zhuǎn)換波反射系數(shù)及由轉(zhuǎn)換波彈性阻抗得到的反射系數(shù)進(jìn)行對比，圖2是轉(zhuǎn)換波反射系數(shù)曲線對比圖。可以看出，在小角度情況下，由Duffaut轉(zhuǎn)換波彈性阻抗得到的反射系數(shù)比由Gonzalez轉(zhuǎn)換波彈性阻抗得到的反射系數(shù)更加接近Zoeppritz方程得到的反射系數(shù)，且方程形式更加簡單。因此，采用Duffaut推導(dǎo)的轉(zhuǎn)換波彈性阻抗公式。

表1 Goodway模型參數(shù)Table 1 Parameters of Goodway model

1.1.2 縱橫波彈性阻抗標(biāo)準(zhǔn)化

Whitecombe［10］于2002年對彈性阻抗進(jìn)行標(biāo)準(zhǔn)化處理，消除了入射角對彈性阻抗尺度的影響。彈性阻抗標(biāo)準(zhǔn)化公式如下：

其中：

式中：vPi，vSi和ρi（i＝1，…，N）分別是第i個采樣點(diǎn)的縱波速度、橫波速度及密度。

同理，筆者對轉(zhuǎn)換波彈性阻抗也進(jìn)行標(biāo)準(zhǔn)化處理：

圖2 Goodway模型不同近似的轉(zhuǎn)換波反射系數(shù)Fig.2 Different PS wave approximate reflection coefficients of Goodway model

1.2 反演流程

王保麗等［15］提出了彈性阻抗反演的基本流程。筆者針對縱橫波彈性阻抗聯(lián)合反演，提出了以下反演流程（圖3），包括縱波角度道集和轉(zhuǎn)換波角度道集的提取、彈性阻抗和轉(zhuǎn)換彈性阻抗的計(jì)算、地震子波的提取和縱橫波彈性阻抗反演等步驟。

1）縱波角度道集和轉(zhuǎn)換波角度道集的提取。在進(jìn)行縱橫波彈性阻抗聯(lián)合反演之前，需提取不同角度的縱波地震數(shù)據(jù)體和轉(zhuǎn)換波地震數(shù)據(jù)體。對于縱波地震數(shù)據(jù)體，把共中心點(diǎn)道集轉(zhuǎn)化為共角度道集［18］；對轉(zhuǎn)換波地震數(shù)據(jù)體，把共轉(zhuǎn)換點(diǎn)道集轉(zhuǎn)換為共角度道集［19］。

2）彈性阻抗和轉(zhuǎn)換波彈性阻抗的計(jì)算。根據(jù)測井資料的縱波速度、橫波速度以及密度利用公式（4）和（5）計(jì)算井位置處的彈性阻抗曲線和轉(zhuǎn)換波彈性阻抗曲線。彈性阻抗及轉(zhuǎn)換波彈性阻抗的計(jì)算一方面會影響地震子波的提取，另一方面可為縱橫波彈性阻抗聯(lián)合反演提供低頻信息，從而可以提高反演的分辨率。

3）地震子波的提取。地震子波的提取條件是使觀測地震數(shù)據(jù)與合成地震記錄的誤差能量差達(dá)到最小。分別提取不同角度的縱波和轉(zhuǎn)換波角度道集的地震子波。地震子波的提取是反演的關(guān)鍵問題之一，會直接影響地震反演的精度［20］。

4）縱波彈性阻抗體和轉(zhuǎn)換波彈性阻抗體反演。獲取地震子波和測井處的彈性阻抗曲線和轉(zhuǎn)換波彈性阻抗曲線之后，對不同角度的縱波角度道集和轉(zhuǎn)換波角度道集進(jìn)行反演，得到彈性阻抗體和轉(zhuǎn)換波彈性阻抗體［13］。

圖3 縱橫波彈性阻抗聯(lián)合反演流程Fig.3 Workflow of joint elastic impedance inversion of PP and PS wave

1.3 彈性參數(shù)提取

為了將彈性阻抗體和轉(zhuǎn)換波彈性阻抗體與實(shí)際的儲層預(yù)測聯(lián)系起來，需從彈性阻抗體和轉(zhuǎn)換波彈性阻抗體中提取縱波速度、橫波速度及密度等彈性參數(shù)。由于每一個采樣點(diǎn)處的K值不同，在提取彈性參數(shù)時需計(jì)算每個采樣點(diǎn)處的a（θ），b（θ，K），c（θ，K），m（θ，K），n（θ，K）。對某一采樣點(diǎn)處的彈性阻抗和轉(zhuǎn)換波彈性阻抗作線性化處理，即對式（4）和（5）兩邊取對數(shù)，有

式中：EI0＝vP0ρ0，EIPS0＝vS0ρ0。

在不同角度入射的情況下，所有采樣點(diǎn)線性化處理的彈性阻抗和轉(zhuǎn)換波彈性阻抗可以寫成矩陣形式，即：

其中：

式中：θn為第n個入射角；diag表示對角矩陣。

通過求解上述矩陣，可以從彈性阻抗體和轉(zhuǎn)換波彈性阻抗體中提取縱波速度、橫波速度和密度等彈性參數(shù)，進(jìn)一步識別儲層特征。

2 二維模型及實(shí)際井模型測試

2.1 2D模型縱橫波彈性阻抗聯(lián)合反演

根據(jù)縱橫波彈性阻抗聯(lián)合反演的基本原理，利用Marmousi2模型得到合成角度道集，并在PP和PS角度道集分別加入了SNR為3的隨機(jī)噪音，對該反演算法進(jìn)行模型測試。圖4是CDP 42（common depth point，CDP）正演得到的縱波角度道集和轉(zhuǎn)換波角度道集（用PP波時間表示）；圖中圓圈內(nèi)是氣層在不同角度下的地震特征，較其他層段反射振幅更加明顯。圖5是抽取CDP42處不同角度的彈性阻抗和轉(zhuǎn)換波彈性阻抗的反演結(jié)果；從圖中可以看出，不同角度的彈性阻抗和轉(zhuǎn)換波彈性阻抗反演結(jié)果與原始值吻合程度較高。圖6是由縱橫波彈性阻抗聯(lián)合反演和由單純縱波彈性阻抗反演提取的彈性參數(shù)與真實(shí)值的對比；結(jié)果顯示，縱橫波彈性阻抗聯(lián)合反演得到的彈性參數(shù)比單純縱波彈性阻抗反演精度高，特別是密度曲線更加穩(wěn)定，提取的彈性參數(shù)與真實(shí)值吻合更好。

圖7是不同角度的縱波和轉(zhuǎn)換波角度疊加剖面（用PP波時間表示）；在縱波角度疊加剖面中，藍(lán)色圈內(nèi)地震反射振幅明顯增加，為主要的含氣儲層，但是該儲層在轉(zhuǎn)換波角度疊加剖面中沒有響應(yīng)；這是因?yàn)檗D(zhuǎn)換波對流體的變化不敏感。圖8是縱橫波彈性阻抗聯(lián)合反演得到的不同角度的彈性阻抗剖面和轉(zhuǎn)換波彈性阻抗剖面；在彈性阻抗剖面和轉(zhuǎn)換波彈性阻抗剖面中均能清楚地顯示氣層所在的位置，這為下一步彈性參數(shù)的提取提供了數(shù)據(jù)基礎(chǔ)。圖9是經(jīng)縱橫波彈性阻抗聯(lián)合反演提取的彈性參數(shù)（縱波速度、橫波速度和密度）剖面與模型剖面對比；結(jié)果顯示，縱橫波彈性阻抗聯(lián)合反演的彈性參數(shù)估計(jì)值與真實(shí)值基本吻合，能夠準(zhǔn)確反映巖層信息和流體信息。

2.2 實(shí)際井模型的縱橫波彈性阻抗聯(lián)合反演

對某實(shí)際工區(qū)A井建立模型，在PP和PS角度道集加入S/N為3的隨機(jī)噪音進(jìn)行反演，驗(yàn)證該反演方法的有效性。圖10是正演得到的縱波和轉(zhuǎn)換波角度道集（用PP波時間表示）。從該圖中可以看出，時間為2.0～2.1s時，振幅隨入射角的變化有明顯的變化，為主要的儲層段；表明利用該特征進(jìn)行疊前縱橫波彈性阻抗聯(lián)合反演可以較好地識別儲層。圖11是不同角度彈性阻抗和轉(zhuǎn)換波彈性阻抗的反演結(jié)果，圖中顯示估算值與真實(shí)值匹配程度較高。圖12是從彈性阻抗和轉(zhuǎn)換波彈性阻抗中提取的彈性參數(shù)與真實(shí)值的對比。可見，復(fù)雜的實(shí)際井資料模型，也能得到較好的反演結(jié)果。

圖4 正演得到的角度道集（用PP波時間表示）Fig.4 Pre－stack synthetic seismogram（represented by PP wave time）

圖5 CDP 42處彈性阻抗（a，b，c）及轉(zhuǎn)換波彈性阻抗（d，e，f）反演結(jié)果Fig.5 Inversion result of EI（a，b，c）and EIPS（d，e，f）in CDP 42

圖6 CDP 42縱橫波彈性阻抗聯(lián)合反演（紅色曲線）及單純縱波彈性阻抗反演（綠色曲線）提取的彈性參數(shù)與真實(shí)值（藍(lán)色曲線）對比Fig.6 Comparison of joint PP and PS EI inversion（red curve）and single PP EI inversion（green curve）with original data（blue curve），elastic parameters extracted fromEI and EIPSin CDP 42

3 結(jié)語

在選取精確的轉(zhuǎn)換波彈性阻抗公式的基礎(chǔ)上進(jìn)行縱橫波彈性阻抗聯(lián)合反演方法研究，建立縱橫波彈性阻抗聯(lián)合反演流程，該方法充分利用了縱波和轉(zhuǎn)換波的地震數(shù)據(jù)，既能較好地利用疊前信息，又具有較高的反演分辨率。

二維模型以及實(shí)際井資料模型反演結(jié)果表明，縱波彈性阻抗和轉(zhuǎn)換波彈性阻抗反演精度較高，從中提取的縱波速度、橫波速度及密度等彈性參數(shù)的估測值與真實(shí)值之間吻合程度很好，尤其是密度反演更加精確，克服了常規(guī)地震反演中密度反演精度不高的缺點(diǎn)，為儲層預(yù)測和流體識別提供了可靠的地球物理依據(jù)。

在反演中依然存在以下問題需要進(jìn)一步研究：轉(zhuǎn)換波由偏移距道集轉(zhuǎn)化為角度道集的算法；就實(shí)際資料而言，反演之前橫縱波速度比K的確定。

圖7 縱波角度疊加剖面（a）和轉(zhuǎn)換波角度疊加剖面（b）Fig7 PP wave angle stacked sections（a）and PS wave angle stacked sections（b）

圖8 彈性阻抗（a，b，c）及轉(zhuǎn)換波彈性阻抗（d，e，f）反演結(jié)果Fig.8 Inversion result of EI（a，b，c）and EIPS（d，e，f）

圖9 彈性參數(shù)真實(shí)值（a，b，c）與反演值（d，e，f）Fig.9 Real elastic parameters（a，b，c）and inverted elastic parameters（d，e，f）

圖10 正演得到的角度道集（用PP波時間表示）Fig.10 Pre－stack synthetic seismogram（representecl by PP wave time）

圖11 彈性阻抗（a，b，c）及轉(zhuǎn)換波彈性阻抗（d，e，f）反演Fig.11 Inversion result of EI（a，b，c）and EIPS（d，e，f）

圖12 縱橫波彈性阻抗聯(lián)合反演提取的彈性參數(shù)（紅色曲線）與真實(shí)值（藍(lán)色曲線）對比Fig.12 Comparison of joint PP and PS EIinversion（red curve）with original data（blue curve），elastic parameters extracted fromEI and EIPS

（References）：

［1］印興耀，張繁昌，孫成禹.疊前地震反演［M］.青島：中國石油大學(xué)出版社，2010：34－56.Yin Xingyao， Zhang Fanchang， Sun Chengyu.Prestack Seismic Inversion［M］.Qingdao：China Universtiy of Petroleum Press，2010：34－56.

［2］Aki K，Richards P.Quantitative Seismology：Theory and Methods［M］.Mill Valley：University Science Books，1980.

［3］張廣智，王丹陽，印興耀，等.基于 MCMC疊前地震反演方法研究［J］.地球物理學(xué)報(bào)，2011，54（11）：2926－2932.Zhang Guangzhi，Wang Danyang，Yin Xingyao，et al.Study on Prestack Seismic Inversion Using Markov Chain Monte Carlo［J］.Chinese Journal of Geophysics，2011，54（11）：2926－2932.

［4］陳建江，印興耀，張廣智.層狀介質(zhì)AVO疊前反演［J］.石油地球物理勘探，2006，41（6）：656－662.Chen Jianjiang， Yin Xingyao， Zhang Guangzhi.Prestack AVO Inversion of Layered Medium［J］.Oil Geophysical Prospecting，2006，41（6）：656－662.

［5］張廣智，陳懷震，印興耀，等.基于各向異性AVO的裂縫彈性參數(shù)疊前反演方法［J］.吉林大學(xué)學(xué)報(bào)：地球科學(xué)版，2012，42（3）：845－851.Zhang Guangzhi，Chen Huaizhen，Yin Xingyao，et al.Method of Fracture Elastic Parameter Inversion Based on Anisotropic AVO［J］.Journal of Jilin University：Earth Science Edition，2012，42（3）：845－851.

［6］Connolly P.Elastic Impedance［J］.The Leading Edge，1999，18（4）：438－452.

［7］Larsen J.AVO Inversion by Simultaneous PP and PS Inversion［D］.Calgary：University of Calgary，1999：9－37.

［8］Milos S.Elastic Impedance Inversion in Practice［C］//SEG Technical Program Expanded Abstracts.Calgary：［s.n.］，2000：689－692.

［9］Duffaut K，Martin L.Shear－Wave Elastic Impedance［J］.The Leading Eadge，2000，19（11）：1222－1229.

［10］Whitecombe D.Elastic Impedance Normalization［J］.Geophysics，2002，67（1）：60－62.

［11］Gonzalez E.Physical and Quantitative Interpretation of Seismic Attribute for Rocks and Fluids Identification［D］.Alto：Stanford University，2006：6－35.

［12］馬勁風(fēng).地震勘探中廣義彈性阻抗的正反演［J］.地球物理學(xué)報(bào)，2003，46（1）：118－124.Ma Jinfeng.Forward Modeling and Inversion Method of Generalized Elastic Impedance in Seismic Exploration［J］.Chinese Journal of Geophysics，2003，46（1）：118－124.

［13］王保麗，印興耀，張繁昌.彈性阻抗反演及應(yīng)用研究［J］.地球物理學(xué)進(jìn)展，2005，20（1）：89－92.Wang Baoli，Yin Xingyao，Zhang Fanchang.Elastic Impedance Inversion and Its Application［J］.Progress in Geophysics，2005，20（1）：89－92.

［14］王保麗，印興耀，張繁昌，等.基于Fatti近似的彈性阻抗方程反演［J］.地球物理學(xué)進(jìn)展，2008，23（1）：192－197.Wang Baoli，Yin Xingyao，Zhang Fanchang，et al.Elastic Impedance Equation Based on Fatti Approximation and Inversion［J］.Progress in Geophysics，2008，23（1）：192－197.

［15］王保麗，印興耀，張繁昌，等.基于Gray近似的彈性波阻抗方程及反演［J］.石油地球物理勘探，2007，42（4）：435－439.Wang Baoli，Yin Xingyao，Zhang Fanchang，et al.Gray Approximation Based on Elastic Wave Impedance Equation and Inversion［J］.Oil Geophysical Prospecting，2007，42（4）：435－439.

［16］Cui Jie，Han Liguo，Liu Qiankun，et al.P－SV Wave Elastic Impedance and Fluid Identification Factor in Weakly Anisotropic Media［J］.Applied Geophysics，2010，7（2）：135－142.

［17］陳天勝，劉洋，魏修成.縱波和轉(zhuǎn)換波聯(lián)合AVO反演方法研究［J］.中國石油大學(xué)學(xué)報(bào)：自然科學(xué)版，2006，26（1）：33－37.Chen Tiansheng，Liu Yang，Wei Xiucheng.Joint Amplitude Versus Offset Inversion of P－P and P－SV Seismic Data［J］.Journal of China University of Petroleum：Edition of Natural Science，2006，26（1）：33－37.

［18］張廣智，王文龍，印興耀，等，基于射線理論的角度道集提取方法綜述［J］.地球物理學(xué)進(jìn)展，2011，26（5）：1696－1670.Zhang Guangzhi，Wang Wenlong，Yin Xingyao，et al，Overview of Extraction Methods of Angle Gathers Based on Ray Theory［J］.Progress in Geophysics，2011，26（5）：1696－1707.

［19］Tessmer G.Processing of PS－Reflection Data Applying a Common Conversion Point Stacking Technique［J］.Geophysical Prospecting，1990，38：267－286.

［20］張廣智，劉洪，印興耀.井旁道地震子波精細(xì)提取方法［J］.石油地球物理勘探，2005，40（2）：158－162.Zhang Guangzhi，Liu Hong，Yin Xingyao，Method for Fine Picking up Seismic Wavelet at Uphole Trace［J］.Oil Geophysical Prospecting，2005，40（2）：158－162.