蔡杰雄
(中國石油化工股份有限公司石油物探技術(shù)研究院,江蘇南京211103)
地震勘探所利用的數(shù)據(jù)主要是反射波。反射波場由背景速度擾動和高波數(shù)速度擾動共同構(gòu)成[1-2],因此基于反射波的成像分為疊前深度偏移成像(確定高波數(shù)擾動)和反射波層析成像(確定背景速度擾動)兩個步驟。偏移成像用于定位地下反射位置,其實現(xiàn)了從數(shù)據(jù)空間到成像空間的映射;層析成像用于反演傳播路徑上的速度場,其實現(xiàn)了從成像空間到模型空間的映射。當然,這兩個步驟也可以耦合在一起進行數(shù)據(jù)域的迭代反演成像(如全波形反演)。相比于常規(guī)射線走時層析,結(jié)合偏移成像在成像域進行波動方程線性化走時層析速度反演是當前比較實用有效且精度較高的技術(shù)組合。
基于傍軸近似的常規(guī)射線類成像方法是目前應(yīng)用于地震數(shù)據(jù)層析及偏移的重要方法。傳統(tǒng)的射線追蹤方法一般局限于對射線路徑及走時的描述,其實現(xiàn)靈活高效、沒有傾角限制且容易拓展到起伏地表。但是,高頻近似下的常規(guī)射線追蹤認為中心射線代表著地震波的主能量,在實現(xiàn)上僅僅利用中心射線來描述地震波傳播,這樣的近似處理只能反映地震波的高頻運動學特征。且在數(shù)值計算上對于復(fù)雜介質(zhì)可能存在焦散及多到達時問題,因此常規(guī)射線類成像方法(包括偏移和層析)在復(fù)雜構(gòu)造情況下的應(yīng)用效果并不十分理想。
高斯束傳播算子是射線傳播算子的發(fā)展,HILL[3-4],HALE[5]以及POPOV等[6-7]將高斯束方法應(yīng)用到偏移處理中并取得了較好的成像效果。Hill于1990年首先提出高斯束疊后深度偏移方法,并于2001年將該方法推廣到疊前深度偏移,利用共偏移距道集中的對稱性解決了高斯束疊前深度偏移中的執(zhí)行效率問題,非常適合海上拖纜地震數(shù)據(jù)的成像處理。對于陸上地震數(shù)據(jù)的成像處理,GRAY[8]于2005年提出了快速精確的、適合于陸地起伏地表數(shù)據(jù)的炮域高斯束疊前深度偏移方法,解決了高斯束偏移方法中的局部平面波分解因地表高程變化及地表速度橫向變化帶來的問題。自HILL奠定了高斯束偏移方法的理論基礎(chǔ)以來,衍生出了一系列實用化的束偏移技術(shù)(如控制束偏移和快速束偏移及自適應(yīng)束偏移等)。高斯束疊前深度偏移技術(shù)的實現(xiàn)主要包括單個獨立的高斯束傳播的求解及所有高斯束疊加成像兩個步驟。單個獨立的高斯束傳播分兩步求得,即通過運動學射線追蹤求取中心射線的路徑及走時,通過動力學射線追蹤獲取中心射線附近的走時及高頻能量分布。利用相互獨立的高斯束描述波傳播,既保持了射線方法的高效性和靈活性,又考慮了波場的動力學特征。高斯束偏移利用相互獨立的高斯束疊加并成像,解決了射線類方法中的多路徑問題,兼具了初至波到達時Kirchhoff積分偏移的靈活性和波動方程偏移的精確性,是一種精致、精確且實現(xiàn)上靈活高效的深度偏移方法。
與成像方法類似,在層析速度反演中也有介于常規(guī)射線理論和波動理論之間的方法,如胖射線層析[9]、菲涅爾體層析[10-12]、高斯束層析[13]、波路徑層析[14]等。胖射線層析從射線層析向波動層析走近了一步,但其理論基礎(chǔ)不夠嚴格,僅僅是對波動理論感性認識的結(jié)果?;诜颇鶢栿w的波形層析方法是近年來發(fā)展較快的重要方法。它利用射線周圍的波動性質(zhì)來同時擬合走時和波場的振幅,其本質(zhì)是有限頻層析的一種。SEMTCHENOK等[13]最早提出了高斯束層析方法,利用高斯束展布范圍內(nèi)的波場分布建立層析方程,對于一個炮檢對通過建立多個層析方程來減小層析反演的病態(tài)性,在一定程度上提高了反演穩(wěn)定性。但是SEMTCHENOK在建立層析敏感度核函數(shù)(以下簡稱為核函數(shù))時未給出詳細的推導(dǎo)過程,僅僅是根據(jù)直觀認識利用高斯束的局部波場振幅作為加權(quán)系數(shù)進行求解計算,其理論基礎(chǔ)不嚴格。上面所述的層析速度反演方法大多在數(shù)據(jù)域?qū)崿F(xiàn),在面對復(fù)雜構(gòu)造低信噪比數(shù)據(jù)時,受疊前數(shù)據(jù)復(fù)雜波場及低信噪比的影響,實際應(yīng)用較為困難。邵榮峰等[15]將高斯束層析方法應(yīng)用到成像域,并通過自動拾取等技術(shù)實現(xiàn)了偏移速度分析,萬弘等[16]將構(gòu)造約束等實用化方法加入到高斯束層析中以提高反演的穩(wěn)健性,但兩者的理論基礎(chǔ)仍然與文獻[13]相似,同樣是簡單采用了高斯束振幅作為加權(quán)系數(shù)構(gòu)建層析核函數(shù)。李輝等[17]從高斯束作為正演工具的角度重新推導(dǎo)了成像域高斯束層析核函數(shù),推導(dǎo)得到的層析核函數(shù)仍然是以高斯束振幅作為加權(quán)系數(shù)進行計算。為了更好地將偏移成像與層析反演結(jié)合起來,本文將直接從偏移成像條件出發(fā),推導(dǎo)建立成像走時擾動與速度擾動的線性關(guān)系,更直觀地體現(xiàn)層析與偏移的關(guān)系,給出新的成像域走時層析核函數(shù)。鑒于高斯束傳播算子在疊前深度偏移領(lǐng)域中具有較強的實用性[18-22],且能夠高效提取三維共方位-反射角度成像道集[23],本文從高斯束偏移角度道集出發(fā),在波動方程的一階Born近似和Rytov近似下,推導(dǎo)成像域反射波走時層析方程及其核函數(shù),從而發(fā)展一種新的成像域的波動方程線性化走時層析反演方法,并利用高斯束傳播算子計算該走時層析核函數(shù)(簡稱基于高斯束傳播算子的成像域波動方程線性化走時層析方法為高斯束層析),從而形成結(jié)合高斯束偏移及高斯束層析的迭代反演建模及成像方法和流程。
結(jié)合疊前深度偏移成像提取的角度域共成像點道集,成像域的反射波層析成像可以將反投影路徑分解成兩個透射波分支,類似于高斯束或單程波偏移成像的處理方式。其基本思想是認為成像點的走時擾動由炮點到成像點的走時擾動以及成像點到檢波點的走時擾動之和決定。從而,成像域反射波走時層析成像可以建立起包含兩個核函數(shù)分支的走時擾動對速度擾動的線性關(guān)系。這兩個核函數(shù)分別代表從炮點到成像點以及從成像點到檢波點的透射波路徑??紤]到陸上地震數(shù)據(jù)采樣時的地表非水平、空間采樣不規(guī)則和數(shù)據(jù)信噪比低,為了更方便高效地提取方位-反射角成像道集,將偏移成像與層析成像更好地配套起來,本文選擇高斯束算子作為波傳播算子用于計算格林函數(shù),發(fā)展起高斯束偏移和高斯束層析配套方法。相對于常規(guī)射線層析,由于高斯束層析核函數(shù)具有更嚴格的理論推導(dǎo)、更逼近波傳播實際而反演精度更高,同時由于構(gòu)建的層析矩陣稀疏性較低而使得反演穩(wěn)定性更強?;诟咚故鴤鞑ニ阕拥某上裼蜃邥r層析方法與高斯束偏移相結(jié)合,可形成具有典型特征波成像特點、適應(yīng)復(fù)雜構(gòu)造低信噪比數(shù)據(jù)的成像與建模工具。
高斯束偏移從實現(xiàn)上講,主要包括運動學和動力學射線追蹤及高斯束波場疊加兩個部分。前者給出了單個獨立高斯束的計算,后者則表達了空間中任意一點的格林函數(shù)可由其鄰域內(nèi)所有高斯束的疊加得到,如(1)式所示。該格林函數(shù)可用于實現(xiàn)高斯束疊前深度偏移,亦可用于求解高斯束層析核函數(shù)。
(1)
(1)式表示從震源點x=(x,y,z)出發(fā),在目標點x′=(x′,y′,z′)處的格林函數(shù)由每個相互獨立的高斯束在該點的積分得到。其中,ω為圓頻率;p為慢度矢量,代表傳播方向;uGB表示從震源點x=(x,y,z)出發(fā),在目標點x′=(x′,y′,z′)接收的單個高斯束波場,其表達為(2)式,表示高斯束的構(gòu)建通過運動學射線追蹤確定其中心射線的軌跡和走時τ(s),通過動力學射線追蹤獲取中心射線附近的動力學參數(shù)Q和M,從而得到單個高斯束波場。該表達式在文獻[3]中已經(jīng)給出并討論了其具體計算方法,本文僅給出表達式,不再贅述。
(2)
給出了高斯束計算格林函數(shù)的表達式之后,通過頻率域波場相關(guān)的成像條件即可實現(xiàn)高斯束疊前深度偏移。進一步,根據(jù)高斯束疊前深度偏移成像條件,可以推導(dǎo)成像域高斯束走時層析核函數(shù),并通過公式(1) 利用高斯束積分表達格林函數(shù)實現(xiàn)該核函數(shù)的計算,從而將高斯束偏移與高斯束層析有機結(jié)合起來,實現(xiàn)復(fù)雜介質(zhì)的偏移與層析迭代。
對于炮域?qū)崿F(xiàn)的高斯束疊前深度偏移,需要從炮點和接收點分別進行高斯束正演傳播,以求得從炮點出發(fā)的下行高斯束波場和從反射點出發(fā)到達檢波器的上行高斯束波場。然后類似于單程波波動方程偏移,使用上行波場和下行波場的互相關(guān)成像條件:
(3)
(3)式給出的是單炮成像結(jié)果。其中,G(x,xs,ω)與G(x,xr,ω)分別表示從炮點xs和從檢波點xr到成像點x的格林函數(shù),“*”代表共軛;Ds(xs,xr,pr,ω)表示基于炮道集的加高斯窗局部傾斜疊加。將(1)式代入(3)式,進一步表示成:
(4)
(4)式說明了任一單炮的成像結(jié)果是由從炮點出發(fā)的所有下行方向高斯束波場與從檢波點出發(fā)的所有上行方向高斯束波場的相關(guān)得到。當速度存在誤差導(dǎo)致成像存在誤差時,可以從(4)式出發(fā)推導(dǎo)建立速度誤差與成像誤差的線性關(guān)系,給出基于高斯束偏移成像道集的成像域高斯束走時層析核函數(shù),將在下一節(jié)詳細介紹。
在疊前深度偏移過程中提取方位-反射角度道集是進行后續(xù)成像域走時層析速度反演的必備過程。計算地下成像點的方位-反射角的一種有效且高效的方法是分別估算從炮點和檢波點出發(fā)到達成像點處的波場傳播方向。一旦獲得兩者的傳播方向,即可通過簡單的向量代數(shù)運算得到方位角和反射角(張角)。對于射線類偏移方法而言,這個過程相對容易,只需通過旅行時場的空間導(dǎo)數(shù)分別計算震源波場的射線慢度ps和檢波點波場射線慢度pr即可。對高斯束函數(shù)式((2)式)解析求解,即可方便高效地計算旅行時場的空間導(dǎo)數(shù)。
假設(shè)高斯束鄰域內(nèi)任意一點Q,其對應(yīng)的中心射線上的點R,根據(jù)(2)式所示的高斯束函數(shù),則點Q的旅行時可以由點R的旅行時及動力學射線追蹤參量表示為:
(5)
得到高斯束的旅行時可進一步通過旅行時場的空間導(dǎo)數(shù)得到點Q的高斯束波場慢度矢量P(Q)=?T(Q)/?x,因此可以得到慢度矢量的解析表達式:
(6)
k=1,2,3
分別得到炮點波場和檢波點波場的慢度矢量后,根據(jù)反射角和方位角的定義,即可給出具體求解方法,詳見參考文獻[21],此處不再贅述。
層析反演方法可以統(tǒng)一表示成正問題的一階近似:
(7)
仍然基于層析方程式(7),本文從高斯束偏移成像條件出發(fā),推導(dǎo)出基于高斯束積分計算格林函數(shù)的成像域波動方程線性化走時層析(簡稱為高斯束層析)核函數(shù)。從而實現(xiàn)基于高斯束偏移成像道集的成像域高斯束層析方法。
從公式(4)出發(fā),去掉所有的積分項,僅考慮任意一個炮檢對、任一傳播方向的成像結(jié)果(所有炮檢對、所有方向高斯束的成像結(jié)果疊加即為單炮成像結(jié)果),將角度域高斯束偏移(GBM)成像條件重寫為(8)式。為了后續(xù)表達易于區(qū)分,將從震源出發(fā)的下行高斯束波場表示為S(x,ps,ω;xs,xr),將從檢波點出發(fā)的上行高斯束波場表示為R(x,pr,ω;xs,xr),有:
(8)
(8)式是對應(yīng)頻率域任一炮檢對某一傳播方向p的高斯束成像結(jié)果(共方位-反射角成像道集)。其中,x=(x,y,z)表示成像點坐標;θ,φ分別表示成像點的反射張角及方位角。
在波動方程的一階Born近似下,波場U可以分解為背景波場U0和擾動波場ΔU:
(9)
因此從成像條件(8)式可以近似得到擾動像:
(10)
式中:ΔS,ΔR分別為一階Born近似散射場;S0,R0分別為從炮點和檢波點出發(fā)到成像點的背景波場。該成像條件表明,成像點x處像的擾動來自炮點端和檢波點端兩個分支的影響。
進一步地,根據(jù)文獻[24]和文獻[25],一階Born近似散射場ΔS與ΔR可以表示為:
式中:VS和VR分別表示從炮點和從檢波點到成像點的Born波路徑;k0=ω/v0表示背景模型波數(shù);Δv為待反演的速度擾動;G(x,p,ω;y)表示由高斯束積分計算的從點y到點x的格林函數(shù);S0(y,ps,ω;xs,xr)和R0(y,pr,ω;xs,xr)分別表示在背景速度模型中從炮點和檢波點傳播到空間任意一點y的波場。
如果將(11)式代入(10)式,形式上可以給出像的擾動(左端項)與速度擾動(右端項)的關(guān)系式。但實際操作時,顯然不能直接利用該式進行層析反演。對比數(shù)據(jù)域?qū)游龇囱菘梢苑治?數(shù)據(jù)域反演利用正演數(shù)據(jù)與實測數(shù)據(jù)的差在某種范數(shù)下(一般是二范數(shù))最小作為誤差泛函,其實測數(shù)據(jù)是客觀的,可直接用來做逼近標準;如果直接利用(10)式進行成像域反演,則由于客觀上無法得到真實像IGBM(x,θ,φ,ω)從而無法得到擾動像ΔIGBM(x,θ,φ,ω),因此其本質(zhì)是利用一個未知的中間量來估計另一未知量Δv。直接利用像的擾動這個概念來進行反演缺乏嚴格的判斷標準。從另一個角度分析,像的擾動是一個綜合概念,其實質(zhì)包括了走時(位置)擾動和振幅擾動。實際計算時需要將像的擾動退化為走時擾動(深度擾動)或振幅擾動,從而分別建立與速度擾動的關(guān)系式。由于像域的振幅擾動影響因素太復(fù)雜,實際層析反演一般退化為僅利用走時擾動。
考慮退化到僅利用走時擾動進行層析反演,進一步地,將(10)式重寫為:
(12)
整理得:
(13)
在波動方程的Rytov近似下,波場可以表示為u=(A0+ΔA)ei(φ0+Δφ)。成像值是兩個波場相關(guān)得到,因此同樣可以表示I=(A0+ΔA)ei(φ0+Δφ),I0(x,θ,φ,w)=A0eiφ0,由于ΔA?A0,則(13)式可以進一步表示成:
(14)
由于Δφ趨于零,則(14)式左端項近似等于iΔφ,兩邊取虛部得:
(15)
進一步地,根據(jù)文獻[26],單頻相位擾動與單頻走時擾動有近似關(guān)系:
(16)
同時將(11)式及(16)式代入(15)式整理,最終得到成像域單頻走時擾動與速度擾動的關(guān)系式:
(17)
其中,KF為頻率域走時層析核函數(shù),其兩個分支分別表示為:
(18a)
(18b)
從(18)式可以看出,成像域走時層析核函數(shù)的表現(xiàn)形式與LIU等[12]給出的數(shù)據(jù)域菲涅爾體走時層析核函數(shù)的表現(xiàn)形式類似,其本質(zhì)都是有限頻核函數(shù),關(guān)鍵是背景速度下格林函數(shù)的計算。
需要說明的是,在背景模型中的波場可以表示成子波與格林函數(shù)的乘積,因此在推導(dǎo)得到(18)式的過程中約去了子波項(假設(shè)子波不變),僅剩下格林函數(shù)項。
由于實際操作時,走時擾動Δt在時空域(角度域成像道集)測量得到,與頻率無關(guān)。因此,定義帶限地震信號的走時擾動可以用單頻走時擾動加權(quán)疊加得到:
(19)
(20)
最終得到成像域帶限走時擾動與速度擾動關(guān)系式:
(21)
其中,KT為帶限走時層析核函數(shù),其兩個分支分別表示為:
(22a)
(22b)
至此,基于高斯束角度道集、在波動方程的一階Born近似和Rytov近似下導(dǎo)出了成像域帶限走時層析方程(21)式及其對應(yīng)的核函數(shù)表達式(22)式。從(22)式可以看出,該核函數(shù)的本質(zhì)是有限頻核函數(shù),其求解主要是背景波場中格林函數(shù)的計算。利用(1)式的高斯束積分計算格林函數(shù)是一種精度較高且計算量較小的實用化計算方式。
需要說明的是,盡管最終得到的有限頻核函數(shù)在推導(dǎo)過程引入了多次線性化處理,但有限頻核函數(shù)相對于高頻近似的射線層析核函數(shù)仍然更接近于波動傳播的實際情況,是對波動層析的一種折中近似。用此方法替換常規(guī)射線層析能提高精度,且能借鑒常規(guī)射線層析流程實現(xiàn)技術(shù)的實用化。為進一步提高對復(fù)雜地質(zhì)體(復(fù)雜構(gòu)造、強橫向變速等)的成像效果,結(jié)合高斯束成像的優(yōu)勢,有必要在后續(xù)研究中加強非線性項的合理近似和算法優(yōu)化的方法研究。
下面通過數(shù)值試驗分析單頻和帶限走時層析核函數(shù)的特征。設(shè)計背景模型參數(shù):網(wǎng)格個數(shù)Nx=Ny=Nz=201,網(wǎng)格間距dx=dy=dz=10m,速度為3000m/s,震源和檢波點坐標分別為(500,1000,1000)和(1500,1000,1000),帶限走時層析核函數(shù)的頻率范圍取0~40Hz,單頻核函數(shù)計算取中心頻率20Hz。圖1給出了對應(yīng)的單頻及帶限走時層析核函數(shù)。從圖1g和圖1h上可以看出,核函數(shù)在炮檢連線的中心線上的敏感度為零,這與常規(guī)射線層析的核函數(shù)僅分布在射線上的假設(shè)完全相反。劉玉柱[10]詳細分析了這種差異,并指出了常規(guī)射線層析之所以仍然能取得較好效果的原因。
圖1 單頻和帶限走時層析核函數(shù)a 對應(yīng)單頻10Hz的走時層析核函數(shù)(XZ切片); b 對應(yīng)單頻10Hz的走時層析核函數(shù)(YZ切片); c 對應(yīng)單頻20Hz的走時層析核函數(shù)(XZ切片); d 對應(yīng)單頻20Hz的走時層析核函數(shù)(YZ切片); e 對應(yīng)單頻30Hz的走時層析核函數(shù)(XZ切片); f 對應(yīng)單頻30Hz的走時層析核函數(shù)(YZ切片); g 帶限走時層析核函數(shù)(XZ切片); h 帶限走時層析核函數(shù)(YZ切片)
圖1反映的是單頻和帶限走時層析核函數(shù)。根據(jù)(22)式,該核函數(shù)僅僅是成像域走時層析核函數(shù)的一個分支,兩個分支(炮點到成像點及成像點到檢波點)則構(gòu)成了成像域走時層析反演核函數(shù)的完整形態(tài),如圖2所示。用此帶有一定寬度的高斯束層析核函數(shù)替換常規(guī)射線層析核函數(shù)能更準確逼近波實際傳播方式,提高層析反演精度和穩(wěn)定性。
圖2 成像域高斯束走時層析反演核函數(shù)
成像域高斯束層析反演的實際操作流程及具體計算方法可以完全借鑒常規(guī)射線層析反演的框架和算法,包括成像剖面層位及成像道集RMO的自動拾取,矩陣求解,層位約束正則化方法等等。兩者的區(qū)別僅僅是在核函數(shù)的表達與計算上,因此該方法的實用性較強。
需要進一步說明的是,成像域走時層析的具體實現(xiàn)與其走時誤差的具體計算形式有關(guān)。XIE等[27-29]給出了基于炮偏移后按照炮索引排列的成像道集的剩余時差(RMO)所對應(yīng)的核函數(shù)的表現(xiàn)形式,該核函數(shù)的兩個分支并不對稱,這與成像道集上每一個RMO代表一炮的成像走時誤差有關(guān)。本文所推導(dǎo)的核函數(shù)與XIE等給出的核函數(shù)表達形式一致,但
由于是針對高斯束偏移所提取的方位角度域成像道集,所提取的RMO僅與方位角及反射角有關(guān)[30],核函數(shù)則表現(xiàn)為關(guān)于剖面上反射軸法方向?qū)ΨQ的兩個分支,這一點并不會影響理論上的精度,但會使得層析的實現(xiàn)更加自然方便。
根據(jù)某地區(qū)實際地質(zhì)構(gòu)造設(shè)計的速度模型如圖3a 所示,斷層發(fā)育,地層高陡,速度橫向變化較大。橫縱向采樣點為731×550,橫向采樣間隔為10m,縱向采樣間隔為5m。利用聲波正演得到疊前炮記錄(圖3b),炮間距和道間距都是10m,中心放炮,每炮361道。利用等梯度速度模型作為層析初始模型(圖4a)。對初始模型進行高斯束疊前深度偏移,輸出偏移剖面(圖4b)及角度道集(圖5),可以看出初始模型對應(yīng)的偏移剖面上的各層位成像深度并不準確,繞射波也沒有完全收斂。圖5對比了CDP400位置處提取的初始模型角度道集和層析后角度道集,由于速度偏低,初始角度道集同相軸上翹;經(jīng)過層析反演更新速度模型,偏移后角度道集得到了拉平。圖6是高斯束層析經(jīng)過5次迭代后得到的速度模型及其相應(yīng)的高斯束偏移剖面。圖7是常規(guī)射線層析經(jīng)過8次迭代得到的速度模型及其高斯束偏移剖面。與常規(guī)射線層析對比可知,本文發(fā)展的高斯束層析方法能夠提供更加豐富的速度信息,迭代次數(shù)也小于常規(guī)射線層析。從抽取的單道速度曲線對比(圖8)可以看出,高斯束層析反演的速度值分辨率更高,更逼近真實速度。對比圖6b和圖7b的偏移深度可以看出,高斯束層析偏移結(jié)果更加逼近真實深度,反射界面歸位到正確的深度位置,繞射波收斂,斷層位置聚焦更好,更新后偏移提取的角度道集更加接近真實速度模型下偏移提取的角度道集(圖5),說明高斯束層析對復(fù)雜構(gòu)造模型的速度反演精度優(yōu)于常規(guī)射線層析方法。
圖3 某地區(qū)復(fù)雜速度模型(a)及其炮集記錄(b)
圖4 初始速度模型(a)及其高斯束疊前深度偏移剖面(b)
圖5 CDP400角度道集a 初始角度道集; b 正確速度模型角度道集; c 高斯束層析角度道集; d 常規(guī)射線層析角度道集
圖6 高斯束層析更新速度場(a)及其層析后高斯束偏移剖面(b)
實際資料來自中國東北某復(fù)雜構(gòu)造工區(qū)。該區(qū)基底之下是上古生界變質(zhì)巖,變質(zhì)巖之上發(fā)育火山巖。目的層波場復(fù)雜,存在多組斷裂系統(tǒng),速度變化大,深度域建模困難,精確成像難度較大。處理要求縱向上能刻畫火山機構(gòu)特征和噴發(fā)期次,橫向上能分辨接觸關(guān)系。合理突出與火山巖接觸的地層反射。前期進行了時間域速度建模,以此作為深度域初始模型進行后續(xù)高斯束層析反演建模。
本文的研究工作是在前期處理人員進行了預(yù)處理、疊前時間偏移及深度域初始建模及偏移成像的基礎(chǔ)上開展的,主要通過高斯束層析反演后的速度模型進行高斯束疊前深度偏移來體現(xiàn)本文方法的實際應(yīng)用效果。從圖9所示的三維體模型上可以看出,高斯束層析反演的速度模型分辨率明顯高于常規(guī)射線層析反演結(jié)果。進一步對比過井速度剖面及偏移剖面可以更清楚地反映該效果。
從圖10所示的層析反演的速度模型上看,高斯束層析結(jié)果體現(xiàn)了速度模型的更多細節(jié)信息,進一步,通過與井曲線的對比可以看出,高斯束層析反演速度模型的精度更高(圖11)。
圖7 常規(guī)射線層析更新速度場(a)及其層析后高斯束偏移剖面(b)
圖8 真實速度模型、初始速度模型、射線層析反演及高斯束層析反演速度模型中抽取CDP400處的速度對比
圖9 某實際資料常規(guī)射線層析(a)與高斯束層析(b)速度模型體
圖10 某過井線射線層析速度模型(a)及高斯束層析速度模型(b)(綠線為井位置)
圖12為利用射線層析模型及高斯束層析模型對實際數(shù)據(jù)進行高斯束偏移的結(jié)果,圖13為圖12中藍框的放大顯示。從圖12及圖13可以看出,同樣利用高斯束偏移方法,高斯束層析反演的速度模型對應(yīng)的偏移結(jié)果斷裂成像質(zhì)量更高,斷點更干脆,同相軸更連續(xù),整體成像質(zhì)量更高。從該實際資料處理可以認識到,結(jié)合高斯束層析與高斯束偏移技術(shù)的處理方式,可以更好地適應(yīng)復(fù)雜構(gòu)造資料的深度域速度建模,提供更高質(zhì)量的成像結(jié)果。
圖11 射線層析模型、高斯束層析模型及實鉆井速度曲線對比
圖12 利用射線層析模型(a)及高斯束層析模型(b)對實際數(shù)據(jù)進行了高斯束偏移的結(jié)果
圖13 利用射線層析模型(a)及高斯束層析模型(b)對實際數(shù)據(jù)進行高斯束偏移的結(jié)果(局部放大顯示)
本文從高斯束積分表達格林函數(shù)出發(fā),引出了高斯束偏移及提取方位反射角成像道集方法?;诮嵌扔蚋咚故瞥上駰l件,在波動方程的一階Born近似和Rytov近似下,推導(dǎo)給出了成像域波動方程線性化走時層析成像方程及其核函數(shù)表達式,并利用高斯束傳播算子計算該核函數(shù)。利用高斯束層析核函數(shù)替代常規(guī)射線層析核函數(shù)(該核函數(shù)為常數(shù)1),可以改進層析反演精度,加快反演收斂,從而形成了基于高斯束算子的偏移成像與層析成像的聯(lián)合迭代速度建模與偏移處理方法,相對于常規(guī)方法,該方法具有更高的精度和更強的實用性。
本文發(fā)展的高斯束層析反演方法利用高斯束傳播算子計算成像域走時層析核函數(shù),提供了一種新的成像域波動方程線性化近似走時層析反演方法。該方法主要利用反射數(shù)據(jù)的走時信息,對初始模型的依賴性較低。但該方法僅能反演速度模型的低波數(shù)成分,因而主要用于背景速度建模,為偏移成像服務(wù)及為更高精度的反演方法(如FWI)提供初始模型。高斯束層析方法與高斯束偏移技術(shù)相結(jié)合,可形成具有典型特征波成像特點的、適應(yīng)低信噪比數(shù)據(jù)的成像與建模工具,真正體現(xiàn)偏移成像與速度建模一體化的處理思想。
[1]王華忠,馮波,王雄文,等.特征波反演成像理論框架[J].石油物探,2017,56(1):38-49
WANG H Z,FENG B,WANG X W,et al.The theoretical framework of characteristic wave inversion imaging[J].Geophysical Prospecting for Petroleum,2017,56(1):38-49
[2]王華忠,馮波,王雄文,等.地震波反演成像方法與技術(shù)核心問題分析[J].石油物探,2015,54(2):115-125
WANG H Z,FENG B,WANG X W,et al.Analysis of seismic inversion imaging and its technical core issues[J].Geophysical Prospecting for Petroleum,2015,54(2):115-125
[3]HILL N R.Gaussian beam migration[J].Geophysics,1990,55(11):1416-1428
[4]HILL N R.Prestack Gaussian-beam depth migration[J].Geophysics,2001,66(4):1240-1250
[5]HALE D.Migration by the Kirchhoff,slant stack,and Gaussian beam methods[J].Center for Wave Phenomena,1992:CWP-126
[6]POPOV M M,SEMTCHENOK N M,VERDEL A R,et al.Reverse time migration with Gaussian beams and velocity analysis applications[J].Expanded Abstracts of 70thEAGE Annual Conference,2008:F048
[7]POPOV M M,SEMTCHENOK N M,VERDEL A R,et al.Depth migration by the Gaussian beam summation method[J].Geophysics,2010,75(2):S81-S93
[8]GRAY S H.Gaussian beam migration of common-shot records[J].Geophysics,2005,70(1):953-959
[9]VASCO D W,PETERSON J E,MAJER E L.Beond ray tomography:wavepaths and Fresnel volumes[J].Geophysics,1995,60(6):1790-1804
[10]劉玉柱.菲涅爾體地震層析成像理論與應(yīng)用研究[D].上海:同濟大學,2011
LIU Y Z.Theory and applications of Fresnel volume seismic tomography[D].Shanghai:Tongji University,2011
[11]劉玉柱,謝春,楊積忠.基于Born波路徑的高斯束初至波波形反演[J].地球物理學報,2014,57(9):2900-2909
LIU Y Z,XIE C,YANG J Z.Gaussian beam first-arrival waveform inversion based on Born wavepath[J].Chinese Journal of Geophysics,2014,57(9):2900-2909
[12]LIU Y Z,DONG L G,WANG Y M,et al.Sensitivity kernels for seismic Fresnel volume tomography[J].Geophysics,2009,74(5):U35-U46
[13]SEMTCHENOK N M,POPOV M M,VERDEL A R.Gaussian beam tomography[J].Expanded Abstracts of 71stEAGE Annual Conference,2009:U32
[14]BAKKER P,GERRITSEN S,CAO Q.3D RTM-based wave path tomography tested at a realistic scale[J].Expanded Abstracts of 77thEAGE Annual Conference,2015:WS05-C02
[15]邵榮峰,方伍寶,蔡杰雄,等.高斯束層析偏移速度建模方法及應(yīng)用[J].石油物探,2016,55(1):91-99
SHAO R F,FANG W B,CAI J X,et al.A method of migration velocity analysis based on Gaussian beam tomography and its application[J].Geophysical Prospecting for Petroleum,2016,55(1):91-99
[16]萬弘,楊勤勇,蔡杰雄,等.地質(zhì)構(gòu)造約束高斯束層析反演方法與應(yīng)用[J].石油物探,2017,56(5):707-717
WAN H,YANG Q Y,CAI J X,et al.A method of geological structure constrained tomographic inversion based on Gaussian beam and its application[J].Geophysical Prospecting for Petroleum,2017,56(5):707-717
[17]李輝,王華忠,劉守偉.基于高斯束的速度層析方法研究[J].石油物探,2017,56(1):116-125
LI H,WANG H Z,LIU S W.A velocity tomography algorithm based on Gaussian beam[J].Geophysical Prospecting for Petroleum,2017,56(1):116-125
[18]蔡杰雄,方伍寶,王華忠.高斯束深度偏移的實現(xiàn)與應(yīng)用研究[J].石油物探,2012,51(5):469-475
CAI J X,FANG W B,WANG H Z.Realization and application of Gaussian beam depth migration[J].Geophysical Prospecting for Petroleum,2012,51(5):469-475
[19]黃建平,楊繼東,李振春,等.基于有效鄰域波場近似的起伏地表保幅高斯束偏移[J].地球物理學報,2016,59(6):2245-2256
HUANG J P,YANG J D,LI Z C,et al.An amplitude-preserved Gaussian beam migration based on wave field approximation in effective vicinity under irregular topographical conditions[J].Chinese Journal of Geophysics,2016,59(6):2245-2256
[20]岳玉波,李振春,錢忠平,等.復(fù)雜地表條件下保幅高斯束偏移[J].地球物理學報,2012,55(4):1376-1383
YUE Y B,LI Z C,QIAN Z P,et al.Amplitude-preserved Gaussian beam migration under complex topographic conditions[J].Chinese Journal of Geophysics,2012,55(4):1376-1383
[21]劉強,張敏,李振春,等.各向異性介質(zhì)共炮域高斯束偏移[J].石油地球物理勘探,2016,51(5):930-937
LIU Q,ZHANG M,LI Z C,et al.Common-shot domain Gaussian beam migration in anisotropic media[J].Oil Geophysical Prospecting,2016,51(5):930-937
[22]代福材,黃建平,李振春,等.角度域黏聲介質(zhì)高斯束疊前深度偏移方法[J].石油地球物理勘探,2017,52(2):283-293
DAI F C,HUANG J P,LI Z C,et al.Angle domain prestack Gaussian beam migration for visco-acoustic media[J].Oil Geophysical Prospecting,2017,52(2):283-293
[23]蔡杰雄,王華忠,王立歆.基于三維高斯束算子解析的方位反射角道集提取技術(shù)研究[J].石油物探,2016,55(1):76-83
CAI J X,WANG H Z,WANG L X.Azimuth-opening angle domain common-image gathers from 3D Gaussian beam migration[J].Geophysical Prospecting for Petroleum,2016,55(1):76-83
[24]WOODWARD M J.Wave-equation tomography[J].Geophysics,1992,57(1):15-26
[25]JEROEN J,SPETZLER J,SMEULDERS D,et al.Validation of first-order diffraction theory for the traveltimes and amplitudes of propagating waves[J].Geophysics,2006,71(6):167-177
[26]SPETZLER G,SNIEDER R.The fresnel volume and transmitted waves[J].Geophysics,2004,69(3):653-663
[27]XIE X B,YANG H.A migration velocity updating method based on the shot index common image gather and finite-frequency sensitivity kernel[J].Expanded Abstracts of 77thAnnual Internat SEG Mtg,2007:2767-2771
[28]XIE X B,YANG H.The finite-frequency sensitivity kernel for migration residual moveout and its applications in migration velocity analysis[J].Geophysics,2008,73(6):S241-S249
[29]XIE X B,YANG H.A wave-equation migration velocity analysis approach based on the finite-frequency sensitivity kernel[J].Expanded Abstracts of 78thAnnual Internat SEG Mtg,2008:3093-3097
[30]蔡杰雄,王華忠,陳進,等.基于高斯束傳播算子的成像域走時層析成像方法[J].地球物理學報,2017,60(9):3539-3554
CAI J X,WANG H Z,CHEN J,et al.Traveltime tomography in the image domain based on the Gaussian-beam-propagator[J].Chinese Journal of Geophysics,2017,60(9):3539-3554