王維紅， 井洪亮

東北石油大學 地球科學學院， 大慶 163318

?

基于稀疏反演三維表面多次波壓制方法

王維紅， 井洪亮

東北石油大學 地球科學學院， 大慶 163318

三維表面多次波壓制是海洋地震資料預處理中的重要研究課題，基于波動理論的三維表面多次波壓制方法(3D SRME)是數(shù)據(jù)驅(qū)動的方法，理論上來說，可有效壓制復雜構(gòu)造地震數(shù)據(jù)表面多次波.但該方法因?qū)υ嫉卣饠?shù)據(jù)采集要求高而很難在實際資料處理中廣泛應用.本文基于貢獻道集的概念，將稀疏反演方法引入到表面多次波壓制中，應用稀疏反演代替橫測線積分求和，無需對橫測線進行大規(guī)模重建，進而完成三維表面多次波預測，這樣可有效解決實際三維地震數(shù)據(jù)橫測線方向稀疏的問題.基于縱測線多次波積分道集為拋物線的假設，為保證預測后三維表面多次波和全三維數(shù)據(jù)預測的多次波在運動學和動力學特征上基本一致，文中對預測數(shù)據(jù)實施基于穩(wěn)相原理的相位校正.理論模型和實際數(shù)據(jù)的測試結(jié)果表明，本文基于稀疏反演三維表面多次波壓制方法可在橫測線稀疏的情況下，有效壓制三維復雜介質(zhì)地震資料中的表面多次波，從而更好地提高海洋地震資料的信噪比，為高分辨率地震成像提供可靠的預處理數(shù)據(jù)保障.

稀疏反演; 表面多次波; 貢獻道集; 多次波預測

1 引言

多次反射波常出現(xiàn)于海洋地震勘探資料中，而且通常能量很強，從而對一次波的處理和成像造成強烈干擾，所以對當前廣泛應用的成像方法而言，多次波壓制預處理尤為重要(Jakubowicz，1998；李翔和胡天躍，2009；王保麗等，2014).在眾多的多次波壓制方法中，基于數(shù)據(jù)驅(qū)動的3D SRME(Surface-related multiple elimination)方法被認為是最有可能成功預測和壓制表面多次波的方法(Dragoset et al.，2010).基于Berkhout(1980)提出的數(shù)據(jù)矩陣的思想，該方法最初由Verschuur等提出(1992).依據(jù)惠更斯原理(Berkhout and Verschuur,1997; Verschuur and Berkhout, 1997; Sun, 1999)，SRME方法分為兩步來實現(xiàn)：一是多次波預測；二是多次波自適應相減.而多次波預測過程又可分為兩步：一是計算多次波的貢獻道集(Kabir and Verschuur,1995；van Dedem and Verschuur, 2001；Hokstad and Sollie, 2003)；二是對貢獻道集進行疊加，形成預測的多次波道集(van Dedem and Verschuur, 2005).SRME方法雖基于數(shù)據(jù)驅(qū)動，可處理復雜構(gòu)造地震資料，但有效預測多次波的一個重要假設是輸入全地震數(shù)據(jù)(Verschuur，1991；Dragoset et al.，2006)，也就是說每個檢波點位置都有炮點(炮點間隔等于檢波點間隔).常規(guī)觀測系統(tǒng)采集的海洋二維地震資料，在進行簡單的數(shù)據(jù)規(guī)則化后，上述假設條件就基本得到滿足，所以在二維地震資料多次波壓制處理中，SRME方法得到成功應用.

從原理上來講，2D SRME方法可直接擴展到3D，而且三維 SRME方法避免了二維算法的所有假設(Matson and Abma，2005；Dragoset et al.，2010).然而，當前海洋三維地震資料采集的縱測線稀少，因此橫測線方向的采樣間隔較大，這就存在計算成本高、存儲以及空間采樣不足等問題，致使全三維SRME技術(shù)一直未能得到有效應用.

Verschuur(1991)和Ross 等(1999)等將三維地震資料按照二維方法進行多次波壓制預處理，也就是說對不同的檢波線依次進行處理，在完成簡單的數(shù)據(jù)規(guī)則化后，可完成依據(jù)波動理論的多次波壓制處理.很顯然，該方法不能充分考慮三維介質(zhì)的傾角和方位角等信息，對于地層構(gòu)造簡單的地區(qū)可以應用，但是對于稍微復雜的地層，該方法在多次波預測的走時、相位和振幅等諸多方面都存在很大誤差，難以取得理想的多次波預處理結(jié)果(van Dedem and Verschuur，2005).因此，進行全三維SRME實用方法研究是地震勘探的迫切需求.van Dedem和 Verschuur(2001)、Hokstad和Sollie(2003)提出貢獻道集概念，為3D SRME的有效實現(xiàn)探索新的思路.本文借鑒貢獻道集的思想，利用基于稀疏反演方法實現(xiàn)3D SRME表面多次波壓制.

依據(jù)實際海洋地震資料采集的特點計算三維地震資料表面多次波的貢獻道集.先對縱測線進行道集記錄的簡單外推和規(guī)則化后，形成可按照二維方法求取貢獻道集的地震數(shù)據(jù)，橫測線方向數(shù)據(jù)很稀疏，但是橫測線方向地震同相軸可視為雙曲線或拋物線，因此可將該同相軸沿給定的曲線路徑進行積分，進而可進行反演計算，以代替橫測線方向的求和過程，降低3D SRME對運算數(shù)據(jù)的要求，避免數(shù)據(jù)重建的處理，減少貢獻道集菲涅爾帶求和過程中的空間假頻.由于拋物Radon變換求解矩陣小，頻率域解耦和計算效率高等優(yōu)點，本文針對三維地震數(shù)據(jù)橫測線稀疏的情況，并依據(jù)其拋物線型同相軸的假設，進行稀疏反演三維表面多次波壓制方法研究，并對理論和實際模型進行了試算，顯示了多次波壓制方法的效果和實用性.

2 三維SRME表面多次波壓制基本原理

Verschuur(1992)等提出了壓制表面多次波的2D SRME級數(shù)展開法，即：

P(ω)=D(ω)-A(ω)D2(ω)+A2(ω)D3(ω)-…，

(1)

其中，P(ω)為只含有效波的地震數(shù)據(jù)，D(ω)為含多次波的地震數(shù)據(jù)，A(ω)為與表面反射系數(shù)r0和震源子波S(ω)相關(guān)的表面算子，公式為

(2)

其中，r0為表面反射系數(shù)，S-1(ω)為逆源子波.

從式(1)和(2)可知，2DSRME方法依賴于原始地震數(shù)據(jù)、表面反射系數(shù)以及震源子波信息.但是實際上震源子波未知，Verschuur等(1992)假設多次波壓制后地震數(shù)據(jù)能量最小，利用最小二乘方法實現(xiàn)了多次波的自適應衰減.

理論上講，級數(shù)展開法要求預測所有階表面多次波，而且要求數(shù)據(jù)矩陣的多次相乘計算.基于多次波形成的物理原理，Berkhout和Verschuur(1997)等給出了多次波壓制的SRME迭代計算法，迭代形式如下：P(i+1)(ω)=D(ω)-A(ω)P(i)(ω)D(ω)，

(3)

其中，i代表迭代次數(shù).一般來說，首次迭代時，取P(1)(ω)=D(ω)，類似于級數(shù)展開法，通過預測多次波與總波場的最小二乘匹配可估計出震源子波，進而實現(xiàn)多次波壓制處理.

在公式(3)中，定義IR(ω)代表第i次迭代的多次波壓制結(jié)果，即IR(ω)=P(i)(ω)，定義M(ω)為預測多次波.則多次波預測公式可表述為

(4)

對任一頻率，依據(jù)式(4)預測多次波需要計算矩陣乘法，也就是對矩陣IR(ω)第i行和矩陣D(ω)的第j列進行內(nèi)積計算.

由公式(4)，單道二維地震數(shù)據(jù)的多次波預測可進一步表示為

(5)

將二維多次波預測方法直接擴展到三維，得到3DSRME多次波預測公式為

×D(xk,yk,ω;xs,ys)，

(6)

其中，xr、yr、xs、ys分別表示檢波點與炮點的空間坐標，xk、yk分別表示求和空間范圍內(nèi)震源和檢波點位置的坐標變量.

依據(jù)貢獻道集的概念，式(6)中的求和可進一步用貢獻道集作以說明，若令Mxy為多次波貢獻道集，則有

Mxy(xr,yr,xs,ys,ω;xk,yk)=IR(xr,yr,ω;xk,yk)

×D(xk,yk,ω;xs,ys)，

(7)

上式也可寫成褶積表達式為

mxy(xr,yr,xs,ys,t;xk,yk)=ir(xr,yr,t;xk,yk)

*d(xk,yk,t;xs,ys)，

(8)

因此，應用貢獻道集的思想，地震數(shù)據(jù)多次波預測的數(shù)學描述為

(9)

其中，mxy為多次波貢獻道集的時空域表示.不難發(fā)現(xiàn)，SRME方法壓制多次波可分為三個步驟：即貢獻道集計算、貢獻道集疊加和多次波自適應相減.

設計如圖1a所示的三維地質(zhì)模型，測試分析二維和三維SRME方法的多次波預測和壓制效果.模型最上層為海水，下面為含斷層、陡傾角和尖滅等復雜構(gòu)造的地層.深度域速度模型三個空間方向的網(wǎng)格數(shù)均為400，網(wǎng)格間距為5m.設計Inline方向和Crossline方向均100個檢波器接收的觀測系統(tǒng)，也就是道間距和炮間距均為20m，震源子波為主頻25Hz的零相位Ricker子波，時間方向采樣間隔為4ms，記錄道長為2s.按傳統(tǒng)三維SRME多次波預測對地震數(shù)據(jù)的要求，在每個檢波點都布設了炮點，共模擬10000炮地震記錄.圖1b為正演的三維某單炮地震記錄，可見能量較強的表面多次波.

文中應用公式(7)對圖1所示的復雜模型數(shù)據(jù)進行三維多次波的預測與壓制試算，同時也應用公式(5)對抽取三維數(shù)據(jù)中一條二維測線(Inline50)進行二維SRME的計算.2D和3DSRME算法預測和壓制多次波效果如圖2所示.圖2a為正演模擬的含多次波Inline50線共零偏移距剖面，從圖中可以看出含有能量較強的多次波.圖2b和2c分別為2DSRME和3DSRME算法的多次波預測結(jié)果，對比分析可知，3DSRME算法預測結(jié)果的同相軸更加連續(xù)，特別需要說明的是圖中箭頭所示位置處，適應三維復雜介質(zhì)的全三維SRME算法能夠清晰預測出多次波.如前所述，當?shù)貙觾A角較大或者較為復雜時，2DSRME對多次波預測的各種假設已經(jīng)不符合實際地質(zhì)條件，也就是說對于上述地震數(shù)據(jù)體，2DSRME算法預測的多次波出現(xiàn)較為明顯的走時誤差.3DSRME應用全三維的預測算法，其預測的多次波位置符合實際數(shù)據(jù)情況，因此預測復雜構(gòu)造數(shù)據(jù)多次波時，3DSRME算法精度要明顯高于2DSRME算法.圖2d和2e分別為采用相同自適應相減方法(石穎等，2013)得到的2DSRME和3DSRME多次波壓制結(jié)果，由上兩圖中箭頭所示位置可知，利用2DSRME方法多次波壓制不徹底，多次波能量有較大殘余，相比而言，3DSRME方法壓制多次波的效果較為理想.

圖1 三維復雜模型與正演模擬炮記錄 (a) 三維速度模型； (b) 正演模擬三維單炮記錄.Fig.1 Three dimensional complex geological model and modeling shot record (a) 3D velocity model; (b) The 3D shot record modeled by wave-equation method.

3 基于稀疏反演三維表面多次波壓制

對圖1所示的理論模型，設計觀測系統(tǒng)的炮間隔與檢波點間隔相等，因此模擬數(shù)據(jù)為全波場數(shù)據(jù)，符合SRME方法預測表面多次波的數(shù)據(jù)假設條件.對實際地震數(shù)據(jù)來說，其Inline方向采樣較密，稍做預處理后，可滿足多次波預測計算的要求，但是由于觀測系統(tǒng)的限制，其Crossline方向通常較為稀疏，如果直接按照公式(9)所述方法預測多次波，將會產(chǎn)生嚴重的假頻，通常采用數(shù)據(jù)插值和規(guī)則化的方法，獲得全波場數(shù)據(jù)，再利用SRME方法計算，但是數(shù)據(jù)插值與規(guī)則化后，將產(chǎn)生更大的地震數(shù)據(jù)體，因此計算效率和存儲成本是制約SRME方法有效應用的兩個重要問題.不論是二維還是三維SRME方法，多次波預測的計算效率可通過應用CPU/GPU異構(gòu)平臺加速計算來解決(石穎等，2013)，對海量地震數(shù)據(jù)的規(guī)則化預處理而言，波場存儲則具有很大的挑戰(zhàn)性.

由公式(9)可知，多次波預測由貢獻道集計算與疊加兩步組成.對于實際的三維地震資料而言，預測單道多次波的貢獻道集為三維數(shù)據(jù)體，可將貢獻道集計算分為兩個步驟.為此，定義如下數(shù)據(jù)向量：

(10)

公式(10)表示沿主線方向進行疊加求和后得到的數(shù)據(jù)向量.理論上，再按照Crossline方向進行疊加求和就得到最終預測的多次波.實際上，Crossline方向數(shù)據(jù)通常較為稀疏，無法進行求和運算.vanDedem和Verschuur(2001)證明公式(10)表示的數(shù)據(jù)的時距曲線可近似表示為如下雙曲線形式：

(11)

其中，τ為截距時間，ya為雙曲線頂點位置坐標，v為速度，滿足v=vnmo/2，vnmo為正常時差校正速度.

將式(11)進行Taylor展開，可得到拋物線方程為

(12)

其中，q=1/2τv2為曲率參數(shù)，公式(12)的假設條件是勘探目標層的埋藏深度遠大于Crossline方向的寬度，對多數(shù)實際三維地震數(shù)據(jù)，上述假設是可以滿足的.

圖2 2D和3D SRME預測及壓制多次波對比 (a) 含多次波的原始共偏移距剖面； (b) 2D SRME方法預測的多次波； (c) 3D SRME方法預測的多次波； (d) 基于2D SRME方法預測的多次波壓制結(jié)果； (e) 基于3D SRME方法預測的多次波壓制結(jié)果.Fig.2 Comparison of multiple prediction and suppression using 2D and 3D SRME (a) The modeled common offset section; (b) Predicted surface-related multiple by 2D SRME; (c) Predicted surface-related multiple by 3D SRME; (d) Multiple suppression based on 2D SRME multiple prediction; (e) Multiple suppression based on 3D SRME multiple prediction.

對于頻率域的數(shù)據(jù)，公式(12)的拋物Radon變換形式可表示為

(13)

其中，qi(i=1,2,…,Nq)是曲率參數(shù)，yaj(j=1,2,…,Nya)為拋物線頂點位置，上述方程用矩陣形式表示為

M=LHD，

(14)

其中，L為算子矩陣，上角標H代表復共軛轉(zhuǎn)置，D為數(shù)據(jù)矩陣.其反變換可表示為

(15)

公式(14)模型空間的估計可通過最小二乘反演問題的計算來實現(xiàn).假設模型空間是稀疏的，在模型求解時應用Cauchy概率密度函數(shù)和Gauss-Cauchy準則，可得到稀疏的高分辨率的解(SacchiandUlrych,1995;Tradetal.,2003;HokstadandSollie, 2003).模型空間的最大后驗估計，也就是其解的形式為

(16)

公式(16)為用于求解超定方程的加權(quán)最小二乘解(Wang, 2003).但本文求解的是數(shù)據(jù)空間點數(shù)小于模型空間的情況，也就是所謂的欠定問題，根據(jù)Sacchi和Ulrych(1995)以及Hokstad和Sollie(2003)等給出的矩陣形式，可將公式(16)寫為

(17)

(18)

式中，mj為模型空間矩陣的第j個分量.由式(18)可知，式(17)所表示的為模型空間的非線性稀疏解，該方程可應用迭代法進行求解(SacchiandUlrych,1995;HokstadandSollie,2003).

應用非線性拋物稀疏反演替代多次波貢獻道集稀疏的Crossline方向的疊加，可大幅減少三維數(shù)據(jù)的疊前規(guī)則化預處理.同時這里需要指出的是在得到預測的多次波數(shù)據(jù)時，不需進行模型空間數(shù)據(jù)的反變換和求和運算，而直接對稀疏的模型空間數(shù)據(jù)疊加即得到預測的多次波，但是這種三維多次波預測結(jié)果存在相位誤差.依據(jù)穩(wěn)相近似理論(Bleistein,1984;Wapenaar, 1992)，可應用如下公式進行相位校正：

(19)

模型空間數(shù)據(jù)求和后，應用公式(19)進行校正，就得到最終三維表面多次波預測數(shù)據(jù)為

(20)

上述基于稀疏反演的三維表面多次波預測主要包括四個步驟，即：沿Inline方向貢獻道集求和，Crossline方向貢獻道集數(shù)據(jù)的稀疏反演，反演模型空間數(shù)據(jù)相位校正以及校正后疊加，得到如式(20)所示的三維表面多次波預測結(jié)果.為了對比分析預測效果，文中也進行了多次波的自適應相減計算(BerkhoutandVerschuur, 1997;Abmaetal., 2005; 石穎等, 2013).不同于常規(guī)的SRME方法，基于稀疏反演的表面多次波預測的方法，可避免Crossline方向的數(shù)據(jù)規(guī)則化處理，有效緩解數(shù)據(jù)存儲的壓力，降低三維SRME方法對地震采集數(shù)據(jù)的要求.

4 數(shù)據(jù)測試分析

4.1 理論數(shù)據(jù)試算

對所形成的基于稀疏反演三維表面多次波預測方法，文中也應用CPU/GPU并行加速技術(shù)提高三維數(shù)據(jù)體多次波預測和壓制的計算效率.為驗證稀疏反演3D SRME方法，設計了含四層水平層狀介質(zhì)的理論模型，正演模擬時采用全3D地震數(shù)據(jù)采集形式，觀測系統(tǒng)包括56條炮線，在x方向和y方向分別放置56個檢波器，每條炮線等間隔放置56個震源點，每炮3136道接收，道間距和炮間距均為15 m，時間采樣間隔為4 ms.抽取地震數(shù)據(jù)中的5炮進行顯示和分析(圖3)，圖3a為原始含多次波模型數(shù)據(jù)，圖3b為采用傳統(tǒng)2D SRME方法預測的多次波，圖3c為利用傳統(tǒng)3D SRME方法對全波場數(shù)據(jù)預測的多次波，圖3d為抽取全波場數(shù)據(jù)中四條測線(線間距為75 m)，利用稀疏反演方法預測的多次波.分析可知，相比于2D SRME方法，傳統(tǒng)3D SRME和基于稀疏反演的3D SRME方法改善了預測效果.此外，拋物稀疏反演的3D SRME方法避免了Crossline方向的大規(guī)模數(shù)據(jù)規(guī)則化，極大地降低了計算成本，節(jié)約了大部分存儲空間.

為更好理解本文給出的稀疏反演表面多次波預測原理，圖4給出了三維地震數(shù)據(jù)中單道多次波預測的對比分析.該方法可有效避免數(shù)據(jù)規(guī)則化，特別是稀疏Crossline方向數(shù)據(jù)的大規(guī)模重建處理.基于圖3抽出的四條橫測線數(shù)據(jù)(線間距為75 m)，在完成沿主線方向的貢獻道集計算與疊加后，圖4a給出了某一接收點處 Crossline貢獻道集剖面，在稀疏反演三維表面多次波預測計算中，該道集就是稀疏反演的輸入數(shù)據(jù)，根據(jù)公式(17)計算，稀疏反演的模型空間數(shù)據(jù)如圖4b所示.圖4c為單道多次波預測數(shù)據(jù)的對比，重點對比密集Crossline貢獻道集疊加(傳統(tǒng)3D SRME方法)和模型道集疊加(稀疏反演方法)多次波預測結(jié)果.其中第1道為正演模擬的含多次波的原始道數(shù)據(jù)；基于全波場數(shù)據(jù)，利用傳統(tǒng)3D SRME方法預測的多次波如圖4c中第2道所示；第3道為輸入的數(shù)據(jù)Inline網(wǎng)格不變，但只有4條主線(橫測線方向稀疏)，采用本文稀疏反演三維SRME方法預測的多次波.與第2道對比可知，稀疏反演方法預測多次波在能量級別與走時上和傳統(tǒng)3D SRME方法基本一致，但是相位出現(xiàn)了反轉(zhuǎn)現(xiàn)象，需要應用公式(20)進行穩(wěn)態(tài)相位校正，圖4c第4道為相位校正后的結(jié)果，和第2道對比可知，稀疏反演多次波預測結(jié)果經(jīng)相位校正后和傳統(tǒng)全數(shù)據(jù)SRME方法預測得到的多次波無明顯差別，相位校正后的數(shù)據(jù)體就是多次波自適應減的輸入數(shù)據(jù).

圖3 表面多次波預測與對比 (a) 原始炮數(shù)據(jù)； (b) 2D SRME方法預測； (c) 傳統(tǒng)3D SRME方法預測； (d) 拋物稀疏反演3D SRME方法預測.Fig.3 Surface-related multiple prediction comparison for different prediction approaches (a) Shot records with multiple； (b) Predicted multiple shot gathers using 2D SRME； (c) Predicted multiple shot gathers using traditional 3D SRME； (d) Predicted multiple using the proposed sparse inversion 3D SRME prediction algorithm.

4.2 實際數(shù)據(jù)計算

在理論模型驗證的基礎上，本文也將研究形成的稀疏反演三維表面多次波預測與壓制算法應用于中國某海洋探區(qū)的實際地震資料.該數(shù)據(jù)有6條Inline線(拖纜)，每條拖纜有144個檢波器，線間距為100 m，道間距為12.5 m，炮間距為50 m，炮線和檢波線平行，時間方向采樣間隔為4 ms.根據(jù)稀疏反演多次波預測對地震數(shù)據(jù)的要求，對Inline線方向數(shù)據(jù)進行了規(guī)則化處理，規(guī)則化后炮間距也為12.5 m，計算所用炮數(shù)為204炮.規(guī)則化后數(shù)據(jù)Inline測線方向道間距為12.5 m，而Crossline方向的道間距為100 m，由于橫測線方向間距過大，傳統(tǒng)的3D SRME方法是無法預測準確的多次波數(shù)據(jù)，而且會帶來嚴重的假頻干擾.

應用稀疏反演三維表面多次波壓制方法對主線規(guī)則化后的海洋地震數(shù)據(jù)進行了多次波預測和壓制

圖4 稀疏反演三維SRME方法預測多次波單道數(shù)據(jù)分析 (a) Crossline多次波貢獻道集(作為稀疏反演的輸入數(shù)據(jù))； (b) 稀疏反演的模型空間域剖面； (c) 原始單道數(shù)據(jù)(1)、密集道集疊加(2)、模型道疊加(3)與相位校正道(4)對比.Fig.4 Single trace analysis of multiple prediction by three-dimensional SRME using sparse inversion (a) The Crossline multiple contributions used as input for the sparse inverse processing; (b) Sparse inversion data in model space domain; (c) The single trace of original data (trace1), stack of dense gathers result (trace2), stack of model gathers result(trace3), and the phase correction trace (trace4).

圖5 三維實際海洋數(shù)據(jù)表面多次波壓制 (a) 原始地震數(shù)據(jù)共偏移距剖面； (b) 2D SRME方法預測多次波； (c) 稀疏反演3D SRME方法預測多次波； (d) 2D SRME方法多次波壓制； (e) 稀疏反演3D SRME方法多次波壓制.Fig.5 3D marine dataset for multiple suppression (a) The common offset section of original seismic data; (b) The multiple predicted by 2D SRME; (c) The multiple predicted by sparse inversion 3D SRME; (d) The demultiple section by 2D SRME; (e) The demultiple section by sparse inversion 3D SRME.

處理，為對比分析計算效果，文中也應用二維SRME方法對該實際數(shù)據(jù)進行了處理.圖5a為三維海洋數(shù)據(jù)體中的某條Inline線的含多次波的原始共偏移距剖面，可以看出該數(shù)據(jù)含有較強的海底一階、二階甚至高階的表面多次波.圖5b為2D SRME方法預測的表面多次波，圖5c為基于稀疏反演三維方法預測的表面多次波.對比可知，三維稀疏反演方法預測的表面多次波的同相軸連續(xù)性更好，與原始多次波數(shù)據(jù)在能量對比關(guān)系上具有更好的一致性，也就是說可更好的與原始剖面的多次波同相軸進行擬合，為多次波自適應相減的有效處理提供了數(shù)據(jù)基礎.

針對兩種方法的預測結(jié)果，文中采用前述理論模型測試的多次波壓制算法進行處理，2D SRME及稀疏反演3D SRME多次波壓制結(jié)果分別如圖5d和5e所示.對圖中所示的雙程走時約為0.55 s位置處箭頭所示同相軸分析表明，2D SRME方法壓制多次波后，其海底一階表面多次波同相軸的能量仍很強，而3D稀疏反演的SRME方法中多次波能量得到很好的壓制，同樣，在雙程走時約為1.5 s和1.8 s附近的黑色箭頭也顯示3D稀疏反演SRME方法多次波壓制效果明顯優(yōu)于2D算法.雙程走時約為0.8 s處箭頭所示同相軸為一有效波同相軸，這一點從原始共偏移距剖面和預測的多次波剖面對比可以判斷出，由于2D SRME預測的多次波存在振幅和相位的誤差，使得自適應多次波匹配相減無法有效實現(xiàn)，有效波的振幅受到一定程度的損傷，與2D SRME方法相比，基于全三維假設的3D稀疏反演的SRME方法可較為精確的預測出多次波走時，進而可與原始數(shù)據(jù)中表面多次波有效擬合，可實現(xiàn)保幅的表面多次波壓制處理.

5 結(jié)論

相比于二維方法壓制三維地震數(shù)據(jù)的多次波，三維SRME方法考慮了地震波在地下介質(zhì)中傳播機理和實際情況，更為適合復雜構(gòu)造海洋地震資料的多次波壓制處理，但是該算法對數(shù)據(jù)的采集要求很高.實際采集地震數(shù)據(jù)的橫測線方向通常較為稀疏，為提高計算效率、避免大數(shù)據(jù)體的規(guī)則化處理以及規(guī)則化數(shù)據(jù)的存儲，本文基于貢獻道集的概念，將三維表面多次波貢獻道集求取分兩步實現(xiàn)，也就是在完成縱測線方向的多次波貢獻積分后，應用高分辨率拋物Radon變換稀疏反演代替橫測線的多次波道集求和，進而完成三維表面多次波的預測.應用CPU/GPU并行加速技術(shù)提高三維地震資料多次波預測和壓制的計算效率，可大幅縮短海洋地震資料預處理周期.對于稀疏反演的三維多次波預測數(shù)據(jù)和全三維全數(shù)據(jù)預測的多次波存在的相位差異，文中應用穩(wěn)相原理近似完成預測多次波的相位校正.理論模擬和實際海洋勘探三維數(shù)據(jù)的試算表明，稀疏反演三維表面多次波壓制方法具有計算精度高、適用范圍廣和易于實現(xiàn)的特點.

Abma R, Kabir N, Matson K H, et al. 2005. Comparisons of adaptive subtraction methods for multiple attenuation.TheLeadingEdge, 24(3): 277-280.

Berkhout A J. 1980. Seismic Migration: Imaging of Acoustic Energy by Wave Field Extrapolation. Part A: Theoretical Aspects. Elsevier.

Berkhout A J, Verschuur D J. 1997. Estimation of multiple scattering by iterative inversion, part I: Theoretical considerations.Geophysics, 62(5): 1586-1595.Bisley R, MacNeill M. 2008. Application of true-azimuth 3D SRME in the Northwest Shelf, Australia. ∥78th Annual International Meeting, SEG, Expanded Abstracts, 2436-2440.Bleistein N. 1984. Mathematical Methods for Wave Phenomena. Academic Press.

Dragoset B, Moore I, Kostov C. 2006. The impact of field-survey characteristics on surface-related multiple attenuation.GeophysicalProspecting, 54(6): 781-791. Dragoset B, Verschuur D J, Moore I, et al. 2010. A perspective on 3D surface-related multiple elimination.Geophysics, 75(5): 75A245-75A261.

Hokstad K, Sollie R. 2003. 3D surface-related multiple elimination using parabolic sparse inversion. ∥73rd Annual International Meeting, SEG, Expanded Abstracts, 1961-1964.

Jakubowicz H. 1998. Wave equation prediction and removal of interbed multiples. ∥68th Annual International Meeting, SEG, Expanded Abstracts, 1527-1530.

Kabir M N, Verschuur D J. 1995. Restoration of missing offsets by parabolic Radon transform.GeophysicalProspecting, 43(3): 347-368.

Kelamis P G, Verschuur D J. 2000. Surface-related multiple elimination on land seismic data—Strategies via case studies.Geophysics, 65(3): 719-734.

Li X, Hu T Y. 2009. Surface-related multiple removal with inverse scattering series method.ChineseJ.Geophy. (in Chinese), 52(6): 1633-1640, doi: 10.3969/j.issn.0001-5733.2009.06.026.Matson K H, Abma R. 2005. Fast 3D surface-related multiple elimination using azimuth moveout for multiples. ∥75th Annual International Meeting, SEG, Expanded Abstracts, 2064-2067.

Ross W S, Yu Y, Gasparotto F A. 1999. Travel time prediction and suppression of 3-D multiples.Geophysics, 64(1): 261-277.

Sacchi M D, Ulrych T J. 1995. High-resolution velocity gathers and offset space reconstruction.Geopyhsics, 60(4): 1169-1177.

Shi Y, Wang W H, Li Y, et al. 2013. 3D surface-related multiple prediction approach investigation based on wave equation.ChineseJ.Geophys. (in Chinese), 56(6): 2023-2032, doi: 10.6038/cjg20130623.

Sun Y. 1999. Anti-aliasing multiple prediction beyond 2-D. ∥69th Annual International Meeting SEG, Expanded Abstracts, 1338-1341

Trad D, Ulrych T, Sacchi M. 2003. Latest views of the sparse Radon transform.Geophysics, 68(1): 386-399.

van Dedem E J, Verschuur D J. 2001. 3D surface multiple prediction using sparse inversion. ∥71st Annual International Meeting, SEG, Expanded Abstracts, 1285-1288.

van Dedem E J, Verschuur D J. 2005. 3D surface-related multiple prediction: A sparse inversion approach.Geophysics, 70(3): V31-V43.

Verschuur D J. 1991. Surface-related multiple elimination, an inversion approach[Ph. D. thesis]. Netherlands: Delft University of Technology.Verschuur D J, Berkhout A J. 1997. Estimation of multiple scattering by iterative inversion, Part II: Practical aspects and examples.Geophysics, 62(5): 1596-1611.

Verschuur D J, Berkhout A J, Wapenaar C P A. 1992. Adaptive surface-related multiple elimination.Geophysics, 57(9): 1166-1177.

Wang B L, Sacchi M D, Yin X Y, et al. 2014. Multiple attenuation

based on amplitude preserving Radon transform.ChineseJ.Geophys. (in Chinese), 57(6): 1924-1933, doi: 10.6038/cjg20140623.

Wang Y H. 2003. Sparseness-constrained least-squares inversion: Application to seismic wave reconstruction.Geophysics, 68(5): 1633-1638.

Wapenaar C P A. 1992. The stationary phase method and its application in seismics and acoustics. Lecture notes: Akoestische Beeldvorming, Delft University of Technology.

附中文參考文獻

李翔, 胡天躍. 2009. 逆散射級數(shù)法去除自由表面多次波. 地球物理學報, 52(6): 1633-1640, doi: 10.3969/j.issn.0001-5733.2009.06.026.

石穎, 王維紅, 李瑩等. 2013. 基于波動方程三維表面多次波預測方法研究. 地球物理學報, 56(6): 2023-2032, doi: 10.6038/cjg20130623.

王保麗, Sacchi M D, 印興耀等. 2014. 基于保幅拉東變換的多次波衰減. 地球物理學報, 57(6): 1924-1933, doi: 10.6038/cjg20140623.

(本文編輯 張正峰)

3D surface-related multiple elimination based on sparse inversion

WANG Wei-Hong, JING Hong-Liang

CollegeofEarthSciences,NortheastPetroleumUniversity,Daqing163318,China

Three dimensional surface-related multiple elimination (SRME) is one of the important topics in the processing of seismic data from marine exploration,theoretically, the data-driven SRME based on wave-equation, can suppress all surface-related multiples from complex structure, both in 2D and 3D sense. But actually, because of the high requirement for seismic data acquisition, it is usually difficult to apply 3D SRME for field data demultiple processing.3D multiple suppression approach based on sparse inversion is analyzed.We classify multiple suppression method into two categories, filter and SRME method respectively. For the seismic data from complex geological structure, filter approach doesn′t work well. However, the data-driven SRME approach based on wave equation, can suppress multiple better, which has no requirement for velocity information.For SRME, the full wavefield data requirement is an important disadvantage, which cannot be meet for almost all marine field data. Therefore the data reconstruction is necessary for traditional multiple suppression using SRME.For current marine acquisition geometries, the data is densely sampled in the inline direction, but very sparsely in the crossline direction, we introduce contribution gather concept, and calculate sparsely sampled crossline multiple contribution by means of sparse inversion algorithm.Therefore, the data reconstruction is unnecessary before demultiple,comparing with traditional 3D SRME algorithm, which decrease storage cost greatly. Based on the assumption that the crossline time-distance curves are hyperbolic or parabolic, i.e., after integrated along inline direction after finished the first step processing of contribution gathers, we applied the phase correction algorithm which based on the principle of stationary phase approximation, predicted multiples by the proposed method and result by the full 3D SRME method are basically the same on kinematics and dynamics characteristics.Multiple prediction and adaptive subtraction are two crucial steps for 3D multiple suppression using SRME. We simulate complex model data to test the proposed 3D multiple prediction algorithm, and 10 thousands shot records are modelled, the shot interval and trace interval are both 25 m. The result comparison show the 3D multiple prediction approach can predict the multiple′s amplitude and phase correctly, and also, the subtraction result is superior than 2D algorithm. The horizontal four-layered media is also designed to test the sparse inversion 3D multiple suppression algorithm, there are 3136 shot records in total, the trace and shot intervals are both 25 m, the line interval is 75 m. The single trace and common-offset result show that the proposed approach can predict the multiple′s amplitude and traveltime correctly, and they are very close to the full data circumstance. The test on field data from some area in China show that the proposed sparse inversion method is applicable and effective, where the trace interval is 12.5 m, shot interval is 50 m, and the line interval is 100 m. After summation along inline direction, the partial integration data is transformed to Radon domain using apex-shifted Radon approach based on the assumption of hyperbolic or parabolic events, stacked the Radon imaging, and also applied the phase correction to the sparse inversion solution, the predicted multiple is acquired.After theoretical investigation and data tests, we have the following conclusions, (1) the proposed method is suitable for simple and complex 3D model data, (2) for real seismic data, the inline direction reconstruction is needed, (3) because only forward Radon transform, not inverse Radon transform is used,phase correction is demanded, and (4) the algorithm do not rely on the assumption of full data using sparse inversion approach, and also, widen the application extent of 3D multiple suppression.The result show that the proposed 3D multiple suppression algorithm can improve the S/N ratio during the course of preprocessing, and provide the high quality data for the subsequent high resolution imaging.

Sparse inversion; Surface-related multiple; Contribution gather; Multiple prediction

10.6038/cjg20150725.

國家自然科學基金項目(41474118),國家高技術(shù)研究發(fā)展計劃(863)(2012AA061202)和大連理工大學海岸和近海工程國家重點實驗室開放基金項目(LP1509)聯(lián)合資助.

王維紅，男，1975年生，博士，副教授，主要從事地震資料數(shù)字處理方面的研究. E-mail: wwhsy@sina.com

10.6038/cjg20150725

P631

2015-04-02，2015-07-06收修定稿

王維紅， 井洪亮. 2015. 基于稀疏反演三維表面多次波壓制方法.地球物理學報，58(7):2496-2507,

Wang W H, Jing H L. 2015. 3D surface-related multiple elimination based on sparse inversion.ChineseJ.Geophys. (in Chinese),58(7):2496-2507,doi:10.6038/cjg20150725.