基于解析時間波場外推與波場分解的逆時偏移方法研究

胡江濤, 王華忠

同濟大學(xué)海洋與地球科學(xué)學(xué)院波現(xiàn)象與反演成像研究組， 上海 200092

?

基于解析時間波場外推與波場分解的逆時偏移方法研究

胡江濤, 王華忠

同濟大學(xué)海洋與地球科學(xué)學(xué)院波現(xiàn)象與反演成像研究組， 上海 200092

雙程波方程逆時深度偏移是復(fù)雜介質(zhì)高精度成像的有效技術(shù)，但其結(jié)果中通常包含成像方法引起的噪音和假象，一般的濾波方法會破壞成像剖面上的振幅，其中的假象也會給后續(xù)地質(zhì)解釋帶來困擾.將波場進行方向分解然后實現(xiàn)入射波與反射波的相關(guān)成像能夠有效地消除這類成像噪音，并提高逆時偏移成像質(zhì)量.波傳播方向的分解通常在頻率波數(shù)域?qū)崿F(xiàn)，它會占用大量的存儲和計算資源，不便于在沿時間外推的逆時深度偏移中應(yīng)用.本文提出解析時間波場外推方法，可以在時間外推的每個時間片上實現(xiàn)波傳播方向的顯式分解，逆時深度偏移中利用分解后的炮檢波場進行對應(yīng)的相關(guān)運算，實現(xiàn)成像噪音和成像信號的分離.在模型和實際數(shù)據(jù)上的測試表明，相比于常規(guī)互相關(guān)逆時偏移成像結(jié)果，本文方法能夠有效地消除低頻成像噪音和特殊地質(zhì)構(gòu)造導(dǎo)致的成像假象.

解析時間波場外推； 波場分解； 噪音消除； 逆時深度偏移

1 引言

基于雙向波波動方程的逆時偏移方法(Baysal et al.,1983; McMechan, 1983; Whitmore, 1983)克服了Kirchhoff偏移的焦散和單走時問題，單向波波動方程的陡傾角成像限制問題，成為復(fù)雜介質(zhì)高精度成像的有力工具.借助GPU和多核CPU等高性能計算設(shè)備的發(fā)展(劉紅偉等，2010； 李博等， 2010； 劉守偉等， 2013)，逆時偏移方法已經(jīng)被廣泛應(yīng)用于工業(yè)生產(chǎn)中.但其成像結(jié)果通常存在一定的、由成像方法引起的成像噪音和假象，這是由于不對應(yīng)的波場進行相關(guān)成像而引起的.

低頻成像噪音是逆時偏移成像中最常見的問題.在速度模型中存在高波數(shù)的速度擾動或反射界面時，正、反傳波場中存在較強的背向散射(或反射)，低頻成像噪音非常嚴重.Loewenthal 等(1987)采用平滑速度場的方式來消除背向散射波場進而消除低頻噪音.Fletcher 等(2006)采用壓制反射界面背向散射的方式來消除低頻噪音.Yoon et和Marfurt (2006)、康瑋和程玖兵(2012)采用波印廷矢量方法進行波場傳播方向預(yù)測，實現(xiàn)成像結(jié)果的角度切除，進而消除低頻噪音.Guitton 等(2007)采用最小平方濾波的方式實現(xiàn)低頻噪音消除.Zhang和Sun(2009)采用Laplace濾波方法來進行低頻噪音消除，并通過頻率波數(shù)域關(guān)系實現(xiàn)了Laplace濾波的振幅校正.Liu 等 (2007， 2011)提出采用波場分解方式來構(gòu)造逆時偏移成像條件，實現(xiàn)入射波和反射波的相關(guān)成像，進而達到消除低頻噪音的目的.杜啟振等(2013)對上述低頻噪音壓制方法做了詳細分類并對比了各種方法的優(yōu)缺點，指出基于波場分解的成像方式在復(fù)雜介質(zhì)中能更好地消除低頻噪音和實現(xiàn)振幅保真.

除了低頻成像噪音，F(xiàn)ei 等(2010，2014)分析了逆時偏移對偏移速度的敏感性，并指出在速度模型精度較差時，當向斜構(gòu)造存在強速度梯度時，逆時偏移會產(chǎn)生較強的假象.產(chǎn)生這類假象的構(gòu)造通常會出現(xiàn)在巖丘邊界，潛山頂面及火成巖與圍巖的接觸面等地方.他指出在這類構(gòu)造區(qū)域，采用波場分解成像方式能夠有效地消除成像假象.雖然波場分解成像能夠有效地消除成像噪音，但它通常需要在頻率波數(shù)域進行實現(xiàn).在逆時深度偏移過程中，波場分解通常需要將炮點和檢波點波場存盤并進行關(guān)于時間的傅里葉變換，而在沿時間的外推過程中，時間維處于數(shù)據(jù)的最慢維，因此一般的波場分離成像條件的I/O量和計算量非常大.

針對上述問題，本文提出了解析時間波場外推及波場分解方法.在逆時深度偏移過程中，僅利用外推過程中的、每個時間步的解析時間波場進行不同傳播方向的波場的顯式分離，避免了時間域向頻率域轉(zhuǎn)換的I/O量和計算量.方向分解后的波場進行對應(yīng)的波場相關(guān)可以有效地壓制成像噪音和假象.模型和實際數(shù)據(jù)的測試表明了所提出方法的有效性.

2 逆時深度偏移及波場分解成像條件

(2)其中，dR為檢波點端接收的炮波場，gx=(gx,gy)為檢波點空間位置，uR為延拓后的檢波點波場.

雙向波沿時間的正向和反向外推可以模擬聲波方程(1)和(2)所能模擬的所有波現(xiàn)象，在偏移速度場比較精確的情況下，理論上可以對一次和多次散射(和反射)波進行偏移成像.這是逆時深度偏移成為復(fù)雜介質(zhì)成像首選方法技術(shù)的原因.存在的問題是計算代價大和不對應(yīng)的波場相關(guān)產(chǎn)生的噪音和假象.一般地，后者需要更改成像條件來解決.

本文采用有限差分方法求解公式(1)和公式(2)來實現(xiàn)波場的延拓.將延拓的波場應(yīng)用如下的零延遲相關(guān)成像條件即可提取地下界面的成像結(jié)果，

(3)

其中，tmax為延拓的最大時間值.

在速度模型比較準時，炮點正向外推和檢波點反向外推波場都會產(chǎn)生前向和背向散射(和反射)，炮點和檢波點波場中同方向傳播的波場相關(guān)產(chǎn)生低頻噪音，某些反方向傳播的波場會形成成像假象(Liuetal. 2011;Feietal.，2010).因此可以利用波場的方向分解來實現(xiàn)相對應(yīng)波場的相關(guān)，進而壓制低頻噪音和/或成像假象.

Liu等(2007)提出將炮點波場和檢波點波場進行上下行波分解：

(4a)

(4b)

其中，上標up代表上行波，down代表下行波.需要說明的是，Liu等(2007)所指炮檢端的上下行波是從炮點端來看的，其上下行波的定義與單向波方程的上下行波定義一致.

將公式(4)帶入公式(3)，可得出相關(guān)成像的結(jié)果為

=I1(x)+I2(x)+I3(x)+I4(x),

(5)

其中，I1(x)為炮點下行波場與檢波點上行波場相關(guān)成像結(jié)果，I2(x)為炮點上行波場與檢波點下行波場相關(guān)成像結(jié)果，I3(x)為炮點下行波場與檢波點下行波場相關(guān)成像結(jié)果，I4(x)為炮點上行波場與檢波點上行波場相關(guān)成像結(jié)果.在速度模型精度較高時，I1(x)和I2(x)分別是炮檢波場中反向傳播的波場的相關(guān)，對應(yīng)為反射界面(或散射點)的像.Fei等(2014)指出在速度模型精度不是太高時，I2(x)在特殊地質(zhì)構(gòu)造區(qū)域可能導(dǎo)致成像假象.I3(x)和I4(x)分別是炮檢波場中同向傳播波場的相關(guān)，表現(xiàn)為低頻噪音.由公式(5)也可看出，采用波場分解成像條件可以有效地實現(xiàn)逆時偏移成像噪音和成像信號的分離.

3 解析時間波場外推及波場分解

雖然波場分離成像條件能夠有效地分離逆時偏移中的信號和噪音，但是波場分解方法是其實現(xiàn)的核心.下面以炮點波場為例來說明波場分解過程，檢波點波場的方向分解與其類似，不再贅述.波場方向分解通常借助傅里葉變換來實現(xiàn)，波場的時間和空間傅里葉變換可以表示為

×e-i(ω t+xkx+yky+zkz)dtdxdydz.

(6a)

波場的上下行波的定義為

(7a)

(7b)

同理可以定義左右行波和前后行波.由公式(7)可以看出波場的方向由頻率和波數(shù)的符號共同決定，也就是說常規(guī)波場分離需要在頻率域進行實現(xiàn).時間域逆時偏移實現(xiàn)常規(guī)波場分離，需要將地下波場存盤，并進行關(guān)于最慢維的時間傅里葉變換，所以其I/O量和計算量很大.但是時間域外推中比較便于在波場時間片上獲取空間波數(shù)信息，所以為了提高計算效率，Liu等(2011)根據(jù)頻率和波數(shù)的符號關(guān)系推導(dǎo)出

(8)

其中，

(9a)

(9b)

(9c)

(9d)

由公式(7)和公式(8)可以看出，波場的顯式分離需要由頻率和波數(shù)的符號共同決定，而為了實現(xiàn)的高效，需要盡量利用空間波數(shù)的符號來進行波場分離.那么如果能夠獲取一個時間波場，其頻譜僅包含正或負頻率，那么波場的方向分解就可以依靠空間波數(shù)的符號來實現(xiàn).這樣的時間信號在信號處理中廣泛應(yīng)用，這就是復(fù)信號，也即解析時間信號.解析時間信號的實部為信號本身，而虛部為信號的Hilbert變換.將解析時間信號概念擴展到波場上，我們將其稱為解析時間波場.通過常規(guī)波場傳播方程可以獲得波場本身(即解析時間波場的實部)，但是無法獲得波場的Hilbert變換(即解析時間波場的虛部).可以對地下波場進行Hilbert變換來實現(xiàn)解析時間波場虛部的構(gòu)造，但Hilbert變換為一個時間方向的卷積，這就需要將波場存盤并進行卷積運算，I/O量同樣會增加.Hilbert變換為一個90°相移運算，在頻率域是一個相位濾波器.而波場傳播算子為線性算子，也即傳播方程(1)和(2)的源項和波場成線性關(guān)系(Liuetal., 2006).那么參考波動方程相位編碼的實現(xiàn)(Liuetal., 2006)，可以對波動方程的源項進行Hilbert變換.由于波場傳播算子的線性特征，波場也進行了90°相位移動，這樣就獲得了解析時間波場的虛部.炮點的解析時間波場外推可以表示為

(10)

(11)

(12a)

(12c)

(12d)其中，Re代表取實運算.同樣，可以寫出檢波點的上下行波：

(13a)

4 數(shù)值試驗

為了測試方法的有效性，下面采用向斜模型進行測試，如圖1所示.該模型第一層，第二層和第三層速度值分別為2000m·s-1, 2200m·s-1和4500m·s-1.第二層和第三層之間在實際情況可能為巖丘邊界，火成巖邊界及潛山頂面等地質(zhì)構(gòu)造.模擬數(shù)據(jù)一共100炮，炮間隔為50m，最大偏移距為4000m.采用常規(guī)相關(guān)逆時偏移成像方法得到的結(jié)果如圖2所示，由圖可以看出，成像結(jié)果受到嚴重的低頻噪音影響.

圖1 準確速度模型Fig.1 Exact velocity model

圖2 傳統(tǒng)逆時偏移成像結(jié)果Fig.2 Result of conventional RTM

下面利用例子來驗證本文分解方法的有效性.圖3為采用本文解析時間波場外推方法得到的CDP號為380位置處的解析VSP記錄，其中炮點布設(shè)在CDP號為300位置處.采用常規(guī)波場外推方法僅能得到該解析VSP記錄的實部.圖4為常規(guī)實波場和解析波場在0.5s時間片上的垂直波數(shù)譜對比.由圖3可以看出該記錄在0.5s時間片上僅存在下行波記錄.而常規(guī)實波場的波數(shù)譜并不能指示這一信息，所以常規(guī)波場分離方法通常在頻率域?qū)崿F(xiàn).但在本文解析時間波場的波數(shù)譜上則可以清楚地識別該下行波記錄.利用本文波場分解方法在每個時間片上進行上下行波分解可以得到分離的上下行波記錄，如圖5所示.

在成像空間應(yīng)用上述波場分解方法將波場分解為各個方向矢量，然后進行成像，得到的炮點下行和檢波點上行成像結(jié)果如圖6a所示，炮點上行和檢波

圖3 解析VSP記錄(a)解析VSP記錄實部； (b)解析VSP記錄虛部.Fig.3 Analytic VSP record(a) Real part of analytic VSP record; (b) Imaginary part of analytic VSP record.

圖4 波數(shù)譜對比圖(a)常規(guī)實波場波數(shù)譜；(b)解析波場波數(shù)譜.Fig.4 Comparison of wavenumber spectrums(a) Wavenumber spectrum of conventional real wavefield; (b) Wavenumber spectrum of analytic wavefield.

圖5 分離的上下行波(a)下行波記錄；(b)上行波記錄.Fig.5 Separated up and down going wavefield(a) Down-going wavefield; (b) Up-going wavefield.

圖6 波場分離成像結(jié)果(a)炮點下行與檢波點上行成像結(jié)果；(b)炮點上行與檢波點下行成像結(jié)果；(c)炮點下行與檢波點下行成像結(jié)果；(d)炮點上行與檢波點上行成像結(jié)果.Fig.6 Results of wavefield separation imaging(a) Result from down going source and up going receiver wavefield; (b) Result from up going source and down going receiver wavefield; (c) Result from down going source and down going receiver wavefield; (d) Result from up going source and up going receiverwavefield.

點下行成像結(jié)果如圖6b所示，炮點下行和檢波點下行成像結(jié)果如圖6c所示，炮點上行和檢波點上行成像結(jié)果如圖6d所示.由圖可以看出，圖6a和圖6b為有效信號，而圖6c和圖6d為低頻噪音.由圖也可以看出，本文方法可以實現(xiàn)各個成像分量的分離.在這種速度模型比較精確的情況，可以疊加圖6a和圖6b結(jié)果作為最終成像結(jié)果，也即等效波場分離成像公式結(jié)果.在這種情況下，等效波場分離成像條件的結(jié)果和本文方法結(jié)果等價.

在實際情況下，由于速度建模方法對速度結(jié)構(gòu)的分辨能力限制，通常只能得到一個相對平滑的速度模型.將準確速度模型進行平滑來模擬實際情況，平滑速度模型如圖7所示.采用常規(guī)相關(guān)逆時偏移成像，結(jié)果如圖8a所示.由圖可以看出，其結(jié)果存在一定的低頻噪音，還存在一定的成像假象.采用本文方法提取炮點下行和檢波點上行成像，炮點上行和檢波點下行成像的疊加結(jié)果如圖8b所示，該結(jié)果等價于等效波場分離成像公式得到的結(jié)果.由圖8b可以看出，相比于圖8a，成像假象明顯減少，但是在向斜的中心存在一個較強的成像假象，且該假象呈現(xiàn)與向斜構(gòu)造相似的形態(tài).這個假象產(chǎn)生的主要原因為：由于第二層和第三層之間存在強速度分界面，因此平滑速度模型中存在強梯度速度過渡帶，當炮點和檢波點波場經(jīng)過這一區(qū)域時，都會產(chǎn)生向上回傳的波場，而由于向斜獨特的地質(zhì)形狀，導(dǎo)致炮點上行波和檢波點下行波場都匯聚到向斜中心，進而在應(yīng)用成像條件時產(chǎn)生假象.由此，可以看出這類假象是由炮點上行和檢波點下行波場進行相關(guān)所導(dǎo)致的.采用本文方法分別提取炮點下行和檢波點上行成像如圖8c所示，炮點上行和檢波點下行成像結(jié)果如圖8d所示.由圖可以看出，圖8c為有效信號，而圖8d則是圖8b所包含的假象.在實際情況下，這種假象可能導(dǎo)致巖丘邊界、潛山面及火成巖邊界的錯誤解釋.由圖8c和圖8d也可以看出在這種情況下，本文方法相比于等效波場分離成像條件能夠更好地分離各個成像分量，進而消除假象.

圖7 平滑速度模型Fig.7 Smooth velocity model

圖8 成像結(jié)果對比(a)傳統(tǒng)成像結(jié)果；(b)等效波場分離成像結(jié)果；(c)炮點下行與檢波點上行成像結(jié)果；(d)炮點上行與檢波點下行成像結(jié)果.Fig.8 Imaging results comparison(a) Result of conventional RTM; (b) Result from equivalent wavefield imaging condition; (c) Result from down going source and up going receiver wavefield; (d) Result from up going source and down going receiver wavefield.

圖9 不同偏移方法成像結(jié)果的對比(a)相關(guān)成像結(jié)果；(b)Laplacian濾波方法；(c)波場方向分解方法.Fig.9 Comparison of imaging results from different methods(a) Result of correlation imaging condition; (b) Result from Laplacian filter; (c) Result from wavefield separation method.

下面對比本文方法與常規(guī)濾波方法在噪音消除方面的表現(xiàn).Laplacian濾波方法由于其高效性而被廣泛應(yīng)用于低頻噪音消除中，因此，這里將本文方法與該方法進行對比.測試的模型為一個兩層速度模型，速度分界面在1500m處，第一層速度為2000m·s-1，第二層速度為2500m·s-1.由該模型速度分界面的速度差異可以得出：在該模型上，低頻成像噪音強度較弱，同時成像噪音僅分布在成像子波之上.這樣，利用相關(guān)成像方法得到的成像子波就可近似地作為其他成像方法在成像振幅保真度方面的參考.圖9展示了相關(guān)成像方法、Laplacian濾波方法以及本文波場分解方法的結(jié)果.由圖可以看出，Laplacian濾波方法和本文波場分解方法都較好地對低頻背景噪音進行了去除.但是對比二者能量關(guān)系可以發(fā)現(xiàn)，Laplacian濾波方法顯著地改變了成像振幅.為了更清楚地進行對比，下面抽取各種方法在CDP號為300位置處的單道成像結(jié)果進行對比，如圖10所示.由圖10可以看出，常規(guī)相關(guān)方法由于在反射層之上受到低頻背景噪音的影響，其成像振幅小于零.Laplacian濾波方法和本文方法的結(jié)果在反射層之上的成像振幅趨近于零，由此可以看出它們都較好地實現(xiàn)了噪音的消除.但是進一步對比成像子波可以發(fā)現(xiàn)，本文波場分解方法的成像子波與相關(guān)方法更好地貼合，而Laplacian濾波方法的成像子波能量則存在顯著差異.由此可以看出本文方法在成像子波保真性方面具有較好的表現(xiàn).

圖10 不同成像結(jié)果的單道對比(紅色代表相關(guān)成像結(jié)果；綠色代表Laplacian濾波結(jié)果；藍色代表波場方向分解結(jié)果)Fig.10 Trace comparison of imaging result from different methods

圖11 實際數(shù)據(jù)偏移速度模型Fig.11 Migration velocity model of real data

圖12 實際數(shù)據(jù)成像結(jié)果對比(a) 傳統(tǒng)逆時偏移結(jié)果； (b) 波場分離成像結(jié)果.Fig.12 Imaging results comparison of real data(a) Result of conventional RTM; (b) Result of wavefield separation.

下面采用某地區(qū)實際數(shù)據(jù)進行方法測試，偏移采用的速度模型為層析速度反演結(jié)果，如圖11所示.在速度模型中可以看出較強的速度分界面，因此在常規(guī)逆時偏移結(jié)果(如圖12a)上存在較強的低頻噪音.采用本文方法提取炮點下行和檢波點上行成像，炮點上行和檢波點下行成像的結(jié)果，如圖12b所示.由圖可以看出，低頻噪音得到了很好的消除，同時產(chǎn)生低頻噪音區(qū)域的構(gòu)造更加清楚.

5 結(jié)論

有效的炮點及檢波點波場分解及相應(yīng)波場的相關(guān)成像是高質(zhì)量逆時深度偏移的核心.本文提出了解析時間波場外推及波場分解方法，能夠在逆時深度偏移外推過程中的每個時間片上僅應(yīng)用空間傅里葉變換實現(xiàn)波傳播方向的有效分解，進而實現(xiàn)相應(yīng)波場的相關(guān)成像，得到高質(zhì)量的逆時深度偏移成像結(jié)果，低頻噪音和成像假象被有效壓制.相比于頻率波數(shù)域波場分離方法，該方法具有較高的實現(xiàn)效率，與常規(guī)的相關(guān)+低頻濾波的逆時深度偏移方法相比，增加的計算量較少.當然，該方法也可以應(yīng)用于方位角度成像道集的生成，為后續(xù)的基于逆時深度偏移的速度層析反演和AVA分析與反演提供高質(zhì)量的工具.

Baysal E, Kosloff D D, Sherwood J W C. 1983. Reverse time migration.Geophysics, 48(11): 1514-1524.

Du Q Z, Zhu Y T, Zhang M Q, et al. 2013. A study on the strategy of low wavenumber noise suppression for prestack reverse-time depth migration.ChineseJ.Geophys. (in Chinese), 56(7): 2391-2401, doi: 10.6038/cjg20130725.Fei T W, Luo Y, Aramco S, et al. 2010. De-blending reverse-time migration. ∥SEG Expanded Abstracts, 3130-3134.

Fei T W, Luo Y, Qin F H. 2014. An endemic problem in reverse-time migration. ∥SEG Expanded Abstracts, 3811-3815.

Fletcher R P, Fowler P J, Kitchenside P, et al. 2006. Suppressing unwanted internal reflections in prestack reverse-time migration.Geophysics, 71(6): E79-E82.Guitton A, Kaelin B, Biondi B. 2007. Least-squares attenuation of reverse-time migration artifacts.Geophysics, 72(1): S19-S23.

Kang W, Cheng J B. 2012. Methods of suppressing artifacts in prestack reverse time migration.ProgressinGeophys. (in Chinese), 27(3): 1163-1172.Li B, Liu H W, Liu G F, et al. 2010. Computational strategy of seismic pre-stack reverse time migration on CPU/GPU.ChineseJ.Geophys. (in Chinese), 53(12): 2938-2943, doi: 10.3969/j.issn.0001-5733.2010.12.017.

Liu F Q, Hanson D W, Whitmore N D, et al. 2006. Toward a unified analysis for source plane-wave migration.Geophysics, 71(4): S129-S139.

Liu F Q, Zhang G Q, Morton S A, et al. 2007. Reverse-time migration using one-way wavefield imaging condition. ∥SEG Expanded Abstracts, 26: 2170-2174.

Liu F Q, Zhang G Q, Morton S A, et al. 2011. An effective imaging condition for reverse-time migration using wavefield decomposition.Geophysics, 76(1): S29-S39.

Liu H W, Li B, Liu H, et al. 2010. The algorithm of high order finite difference pre-stack reverse time migration and GPU implementation.ChineseJ.Geophys. (in Chinese), 53(7): 1725-1733, doi: 10.3969/j.issn.0001-5733.2010.07.024.

Liu S W, Wang H Z, Chen S C, et al. 2013. Implementation strategy of 3D reverse time migration on GPU/CPU clusters.ChineseJ.Geophys. (in Chinese), 56(10): 3487-3496, doi: 10.6038/cjg20131023.

Loewenthal D, Stoffa P L, Faria E L. 1987. Suppressing the unwanted reflections of the full wave equation.Geophysics, 52

(7): 1007-1012.

McMechan G A. 1983. Migration by extrapolation of time-dependent boundary values.GeophysicalProspecting, 31(3): 413-420.

Whitmore N D. 1983. Iterative depth migration by backward time propagation. ∥SEG Expanded Abstracts, 382-385.

Yoon K, Marfurt K J. 2006. Reverse-time migration using the Poynting vector.ExplorationGeophysics, 37(1): 102-107.

Zhang Y, Sun J. 2009. Practical issues in reverse time migration: true amplitude gathers, noise removal and harmonic source encoding.FirstBreak, 26(1): 29-35.

附中文參考文獻

杜啟振, 朱釔同, 張明強等. 2013. 疊前逆時深度偏移低頻噪聲壓制策略研究. 地球物理學(xué)報, 56(7): 2391-2401, doi: 10.6038/cjg20130725.

康瑋, 程玖兵. 2012. 疊前逆時偏移假象去除方法. 地球物理學(xué)進展, 27(3): 1163-1172, 10.6038/j.issn.1004-2903.2012.03.041.

李博, 劉紅偉, 劉國峰等. 2010. 地震疊前逆時偏移算法的CPU/GPU實施對策. 地球物理學(xué)報, 53(12): 2938-2943, doi: 10.3969/j.issn.0001-5733.2010.12.017.

劉紅偉, 李博, 劉洪等. 2010. 地震疊前逆時偏移高階有限差分算法及GPU實現(xiàn). 地球物理學(xué)報, 53(7): 1725-1733, doi: 10.3969/j.issn.0001-5733.2010.07.024.

劉守偉, 王華忠, 陳生昌等. 2013. 三維逆時偏移GPU/CPU機群實現(xiàn)方案研究. 地球物理學(xué)報, 56(10): 3487-3496, doi: 10.6038/cjg20131023.

(本文編輯 胡素芳)

Reverse time migration using analytical time wavefield extrapolation and separation

HU Jiang-Tao, WANG Hua-Zhong

WPI,SchoolofOceanandEarthScience,TongjiUniversity,Shanghai200092，China

Two-way wave equation based reverse time migration is powerful in imaging complex area. However, its result usually contains imaging noise which poses challenge to seismic interpretation. The imaging amplitude would be distorted by simple noise filter. Wavefield separation imaging condition is an effective way to suppress the imaging noise and enhance the imaging quality. Conventionally, the wavefield separation is achieved in frequency wavenumber domain which cannot be conveniently implemented in time domain reverse time migration. We propose the analytical wavefield extrapolation method which can efficiently separate wavefield into its directional components during the time domain extrapolation. We use this method in reverse time migration and separate both the source and receiver wavefield. And the imaging noise and signal can be separated after applying the imaging condition to the separated wavefield.The wavefield propagation direction is usually defined by temporal and spatial Fourier transform. In the Fourier domain, the wavefield propagation direction is defined by the sign of frequency and spatial wavenumber. In order to efficiently apply the wavefield separation imaging condition in reverse time migration, we extend the analytical time signal into the wavefield and call it analytical time wavefield. Since the analytical time wavefield contains only the positive frequency component, the wavefield propagation direction can be defined by the sign of spatial wavenumber. To avoid the I/O cost in generating the analytical wavefield at every imaging point, we propose an analytical wavefield propagation equation based on the linear relation between the source term and wavefield in wave equation. We solve the proposed equation by finite difference method. Then we separate the source and receiver wavefield into their up and down going components and apply the imaging condition to the separated wavefield. Four imaging components (i.e., imaging component from up going source and down going receiver wavefield, imaging component from down going source and down going receiver wavefield, imaging component from down going source and up going receiver wavefield and imaging component from up going source and up going receiver wavefield) are effectively separated. Then the imaging noise and signal are separated.We use the syncline model and real data to test the proposed method. Numerical example on syncline model shows that the proposed method can effectively separate imaging component from down going source and up going receiver wavefield, imaging component from up going source and down going receiver wavefield, imaging component from down going source and down going receiver wavefield and imaging component from up going source and up going receiver wavefield. When the migration velocity is not accurate, the correlation imaging condition and equivalent wavefield separation imaging condition would generate imaging artifact in central area of the syncline. And the proposed method can eliminate the artifact. Numerical example on real data shows that this method can generate imaging result free from low frequency noise.Wavefield separation imaging condition can effectively separate the imaging noise and signal in reverse time migration and enhance the imaging quality. The proposed analytical wavefield extrapolation method can separate the source and receiver wavefield into their up and down going components. The imaging noise can be suppressed by applying the imaging condition to the separated wavefield. We want to further use this method to generate common angle imaging gather in reverse time migration and provide input data for AVA and migration velocity inversion.

Analytical time wavefield extrapolation; Wavefield separation; Noise elimination; Reverse time migration

973項目(2011CB201002)，國家自然科學(xué)基金(41374117)以及國家重大專項(2011ZX05005-005-008HZ,2011ZX05006-002,2011ZX05023)資助.

胡江濤，男，1987年生，博士研究生，主要研究方向為地震波傳播與成像方法及應(yīng)用.E-mail:huidaojiaxiang@126.com

10.6038/cjg20150822.

10.6038/cjg20150822

P631

2015-01-08，2015-07-22收修定稿

胡江濤, 王華忠.2015.基于解析時間波場外推與波場分解的逆時偏移方法研究.地球物理學(xué)報，58(8):2886-2895,

Hu J T, Wang H Z. 2015. Reverse time migration using analytical time wavefield extrapolation and separation.ChineseJ.Geophys. (in Chinese),58(8):2886-2895,doi:10.6038/cjg20150822.