侯志強(qiáng) 尹文筍 胡 偉 王曉培 張 巖 劉慶文

(中海石油(中國(guó))有限公司上海分公司 上海 200335)

海底節(jié)點(diǎn)(Ocean Bottom Node，OBN)地震勘探是近年發(fā)展起來的一種新的海底地震勘探技術(shù)，相比于常規(guī)海洋地震勘探技術(shù)，OBN勘探具有如下優(yōu)點(diǎn)：①接收點(diǎn)定位準(zhǔn)確,避免了接收點(diǎn)位置漂移對(duì)勘探結(jié)果帶來的誤差；②接收點(diǎn)位置穩(wěn)定，適宜進(jìn)行油藏地震監(jiān)測(cè)；③不受海底地形限制，可在海上有密集生產(chǎn)平臺(tái)和其他障礙物的區(qū)域開展地震采集工作；④接收點(diǎn)沉放于海底,且OBN節(jié)點(diǎn)具有很強(qiáng)的抓地能力，可很好地與海底耦合,采集的地震信號(hào)噪音小、頻帶寬，具有較高的矢量保真度；⑤不受接收電纜約束，可以靈活設(shè)計(jì)觀測(cè)系統(tǒng)，有利于實(shí)現(xiàn)長(zhǎng)偏移距、寬方位、全方位地震數(shù)據(jù)采集；⑥四分量接收，有利于利用橫波提高特殊目標(biāo)(如裂縫區(qū)、氣云帶等)的勘探精度，同時(shí)縱橫波聯(lián)合反演還能進(jìn)一步提高勘探精度，降低勘探風(fēng)險(xiǎn)。

基于上述認(rèn)識(shí)，近年來殼牌、雪佛龍、BP、CGG以及道達(dá)爾等公司紛紛開展了OBN技術(shù)的研究和應(yīng)用試驗(yàn)，并分別在Dalia油田[1]、尼日尼亞邦達(dá)地區(qū)[2]和墨西哥灣等探區(qū)取得了效果[3-4]。目前，業(yè)界對(duì)OBN技術(shù)的研究主要集中在采集裝備研發(fā)[5]、采集方案設(shè)計(jì)[6-8]和資料成像處理[4,9-12]等3個(gè)方面。在OBN資料成像處理方面，目前的解決方案可歸納為兩類：①只處理縱波分量，不處理橫波分量[9-13]；②縱橫波分量都處理[13-14]。第1種方案浪費(fèi)了橫波信息，且假定OBN的Z分量只接收到縱波，這與實(shí)際情況不符。第2種方案雖然考慮了2個(gè)水平分量中的橫波信息，但它將多分量地震波場(chǎng)看作幾個(gè)標(biāo)量波場(chǎng)的簡(jiǎn)單疊加，忽略了彈性波的矢量特征，這種思路同樣隱含了Z分量只接收到縱波和水平分量而只接收橫波的假設(shè)。

由于OBN資料的特殊性，現(xiàn)有技術(shù)并不能完全滿足OBN多分量資料的逆時(shí)偏移成像需求，主要有2個(gè)原因：①假設(shè)條件與實(shí)際情況偏離過大；②數(shù)值頻散會(huì)嚴(yán)重影響偏移質(zhì)量[15-16]。本文利用通量校正傳輸算法壓制多分量地震波場(chǎng)延拓過程中的數(shù)值頻散，提升OBN多分量地震資料縱橫波聯(lián)合逆時(shí)偏移的精度，具體思路為：在炮點(diǎn)和接收點(diǎn)波場(chǎng)延拓的每個(gè)時(shí)間步，依據(jù)相鄰2個(gè)時(shí)刻的波場(chǎng)以及漫射因子和反漫射因子構(gòu)建漫射通量和反漫射通量，并對(duì)每一時(shí)刻的波場(chǎng)進(jìn)行校正，最后利用校正后的波場(chǎng)進(jìn)行縱橫波成像。

1 基于一階速度-脹縮-旋轉(zhuǎn)方程的彈性波逆時(shí)偏移實(shí)現(xiàn)方法

基于一階速度-脹縮-旋轉(zhuǎn)方程[17]的多分量彈性波逆時(shí)偏移一般有3個(gè)步驟：①炮點(diǎn)波場(chǎng)的正向延拓；②接收點(diǎn)波場(chǎng)的逆時(shí)延拓；③炮、檢波場(chǎng)的互相關(guān)成像。由于輸入的三分量數(shù)據(jù)中每個(gè)分量都同時(shí)包含縱波與橫波，為消除偏移結(jié)果中的波場(chǎng)干擾，一般還要在炮、檢波場(chǎng)延拓過程中進(jìn)行縱、橫波解耦，以獲取純縱波和純橫波的偏移結(jié)果；另外，為壓制逆時(shí)偏移中的低頻噪聲干擾，還需要求取所有成像點(diǎn)在各個(gè)時(shí)刻的縱、橫波傳播方向，以便在后續(xù)成像過程中利用炮、檢波場(chǎng)的傳播方向選取部分波場(chǎng)進(jìn)行成像。

1.1 波場(chǎng)延拓方法

各向同性介質(zhì)中的一階速度-脹縮-旋轉(zhuǎn)方程為[17]

(1)

在交錯(cuò)網(wǎng)格空間中采用高階有限差分技術(shù)[18]對(duì)式(1)進(jìn)行差分離散，可得該方程炮點(diǎn)波場(chǎng)正向延拓和接收點(diǎn)波場(chǎng)逆時(shí)延拓的差分格式。二維情況下，以接收點(diǎn)波場(chǎng)逆時(shí)延拓的θ分量為例，有

θt(i,j)=θt+1(i,j)-

(2)

采用與文獻(xiàn)相同的推導(dǎo)方法,得到上式的穩(wěn)定性條件為

(3)

利用式(2)進(jìn)行炮點(diǎn)和接收點(diǎn)波場(chǎng)的延拓時(shí),還需要在截?cái)噙吔缣幉捎锰厥馑惴?以消除計(jì)算空間的突然截?cái)鄬?duì)波場(chǎng)延拓和成像的影響。本文在接收點(diǎn)波場(chǎng)延拓過程中采用PML吸收邊界條件[19]消除邊界反射，炮點(diǎn)波場(chǎng)延拓過程中則采用隨機(jī)邊界技術(shù)[20]解決截?cái)噙吔鐔栴}，以實(shí)現(xiàn)炮、檢波場(chǎng)的同步計(jì)算，提升效率。

相比于常規(guī)彈性波方程，方程(1)的優(yōu)勢(shì)在于該方程自動(dòng)實(shí)現(xiàn)了縱橫波的解耦，因此利用該方程進(jìn)行逆時(shí)偏移無(wú)需進(jìn)行波場(chǎng)分離工作；同時(shí)，一階速度-脹縮-旋轉(zhuǎn)方程能夠求取純縱波和純橫波的坡印廷矢量，獲得單一類型波的傳播方向信息。依據(jù)李凱瑞 等[21]的思路，在方程(1)延拓過程中可通過求取各成像點(diǎn)在不同時(shí)刻的坡印廷矢量來得到縱橫波傳播方向?；谑?1)的縱橫波坡印廷矢量計(jì)算公式為

(4)

式(4)中：Ep為縱波的坡印廷矢量；Es為橫波的坡印廷矢量?？v橫波的坡印廷矢量即代表了縱橫波的傳播方向。

1.2 波場(chǎng)成像方法

采用基于行波分離的歸一化互相關(guān)成像條件[22]進(jìn)行縱橫波成像。本文所用的二維成像公式為

(5)

(6)

式(5)、(6)中炮點(diǎn)波場(chǎng)和接收點(diǎn)波場(chǎng)的各個(gè)方向行波的求取方法為：在波場(chǎng)延拓過程中首先利用式(4)求取炮檢波場(chǎng)中的縱橫波傳播方向，然后據(jù)此將縱波波場(chǎng)和橫波波場(chǎng)分解為上、下、左、右4個(gè)方向的行波；橫波極性校正因子sign_S的確定方法參考文獻(xiàn)[21]。

2 OBN多分量地震數(shù)據(jù)逆時(shí)偏移中的數(shù)值頻散問題和校正方法

2.1 多分量地震數(shù)據(jù)逆時(shí)偏移中的數(shù)值頻散問題

OBN資料的特點(diǎn)是地震節(jié)點(diǎn)與海底的耦合性好,且避開了海水層對(duì)上行波的影響，因此地震數(shù)據(jù)高頻成分非常豐富，地震資料頻帶寬。這一特點(diǎn)為基于有限差分的彈性波逆時(shí)偏移技術(shù)提出了更高的要求，主要表現(xiàn)在：多分量地震波中的高頻成分要求在逆時(shí)偏移中必須采用更小的剖分網(wǎng)格，否則會(huì)在波場(chǎng)延拓過程中產(chǎn)生嚴(yán)重的數(shù)值頻散現(xiàn)象。這種頻散不是介質(zhì)物理性質(zhì)的反映，而是由波動(dòng)方程的數(shù)值求解方法引起的，本質(zhì)是一種誤差。這種誤差會(huì)使逆時(shí)偏移剖面中出現(xiàn)同相軸錯(cuò)斷或虛假同相軸等現(xiàn)象，降低偏移精度，因此必須采用適當(dāng)?shù)姆椒ㄟM(jìn)行消除。

利用圖1a所示的二維各向同性介質(zhì)模型及其合成兩分量記錄(圖1b、c)，說明高頻地震數(shù)據(jù)逆時(shí)偏移時(shí)所產(chǎn)生的數(shù)值頻散及其對(duì)偏移結(jié)果的影響。圖1b、c是采用虛譜法得到的炮點(diǎn)位于頂界面中間位置處的兩分量記錄，正演模擬該記錄時(shí)所用的震源為主頻80 Hz的Ricker子波。

利用式(1)～(6)對(duì)該兩分量記錄進(jìn)行逆時(shí)偏移，偏移所用的空間網(wǎng)格步長(zhǎng)為5 m(為分析問題方便，本文均采用正方體網(wǎng)格)，炮、檢波場(chǎng)延拓時(shí)均采用時(shí)間2階、空間12階差分格式(為方便說明問題，本文算例均基于時(shí)間2階、空間12階差分格式)，圖2a、b分別為炮點(diǎn)波場(chǎng)正向延拓過程中得到的300 ms時(shí)的縱橫波波場(chǎng)快照，圖2c、d為接收點(diǎn)波場(chǎng)逆時(shí)延拓過程中得到的縱橫波波場(chǎng)快照，可見炮、檢波場(chǎng)均可以觀察到明顯的數(shù)值頻散現(xiàn)象，如果利用這種具有嚴(yán)重?cái)?shù)值頻散的波場(chǎng)信息進(jìn)行互相關(guān)成像，必然會(huì)給偏移結(jié)果帶來誤差和假象(圖3a、b)。因此，當(dāng)多分量地震記錄的主頻較高或高頻成分豐富時(shí)，基于常規(guī)有限差分法的彈性波逆時(shí)偏移技術(shù)無(wú)法取得滿意的成像效果，必須采用適當(dāng)?shù)募夹g(shù)壓制數(shù)值頻散，降低偏移噪聲。

圖1 速度模型和x、z分量地震剖面Fig .1 Velocity model and seismic sections of x & z component

圖2 基于常規(guī)有限差分法的波場(chǎng)延拓快照Fig .2 Snapshots of wave field continuation based on conventional finite difference method

2.2 基于FCT的多分量地震數(shù)據(jù)逆時(shí)偏移中的數(shù)值頻散壓制方法

FCT技術(shù)最初由Boris等[23]提出，用于求解流體力學(xué)方程。之后，地球物理學(xué)者對(duì)其進(jìn)行了改進(jìn)，并應(yīng)用于地震波場(chǎng)模擬中的數(shù)值頻散壓制，取得了明顯效果[24-27]。本文將FCT技術(shù)用于一階速度-漲縮-旋轉(zhuǎn)方程的逆時(shí)偏移領(lǐng)域，以解決高頻數(shù)據(jù)延拓時(shí)炮、檢波場(chǎng)中的數(shù)值頻散壓制問題。

FCT技術(shù)的基本原理本文不贅述，以下僅給出利用FCT壓制炮點(diǎn)波場(chǎng)分量數(shù)值頻散的基本流程與實(shí)現(xiàn)方法，炮點(diǎn)波場(chǎng)其他分量和接收點(diǎn)波場(chǎng)各分量的校正可采用類似方法實(shí)現(xiàn)。

2.2.1漫射運(yùn)算

1) 求解t時(shí)刻的漫射通量。

式(7)中：Pt、Qt分別為時(shí)刻網(wǎng)格點(diǎn)(i,j)處x方向和z方向的漫射通量；η1(0<η1<1)為漫射因子，根據(jù)數(shù)值頻散程度，η1可以設(shè)置為常量或者變量。

2) 利用式(2)求取t+1時(shí)刻的分量θi+1(i,j)值。

3) 使用漫射通量修正分量。

Pt(i-1/2,j)]+[Qt(i,j+1/2)-Qt(i,j-1/2)]

(8)

2.2.2反漫射運(yùn)算

1) 計(jì)算時(shí)刻的漫射通量。

(9)

式(8)中：Pt+1、Qt+1為t+1時(shí)刻對(duì)應(yīng)網(wǎng)格點(diǎn)(i,j)處x方向和z方向的漫射通量；η2為反漫射因子。反漫射因子的設(shè)置準(zhǔn)則同漫射因子，實(shí)際運(yùn)算中由于反漫射修正要同時(shí)補(bǔ)償漫射運(yùn)算造成的損失和有限差分運(yùn)算過程中的振幅損失，因此其值一般略大于漫射因子。

2) 求解t+1時(shí)刻的反漫射通量。

(10)

3) 修正反漫射通量并校正分量，實(shí)現(xiàn)數(shù)值頻散壓制。

(11)

其中

(12)

(13)

(14)

(15)

以上過程即為炮點(diǎn)波場(chǎng)延拓過程中的FCT數(shù)值頻散校正方法。采用類似的方法可進(jìn)行炮點(diǎn)波場(chǎng)的其他分量以及接收點(diǎn)波場(chǎng)的所有分量的校正處理。圖4為本文方法得到的炮、檢波場(chǎng)快照，除在延拓過程中利用式(7)～(15)進(jìn)行FCT校正外，所用的波場(chǎng)延拓參數(shù)與圖2完全相同，炮、檢波場(chǎng)延拓中的數(shù)值頻散得到了明顯壓制。這樣，利用校正后的延拓?cái)?shù)據(jù)進(jìn)行波場(chǎng)成像得到了更準(zhǔn)確的縱橫波成像結(jié)果，如圖5所示。

圖4 基于FCT校正的波場(chǎng)延拓快照Fig .4 Snapshots of wave field continuation based on FCT correction

圖5 基于FCT校正的逆時(shí)偏移成像剖面Fig .5 Sections of RTM based on FCT correction

3 模型測(cè)試與實(shí)測(cè)資料處理試驗(yàn)

3.1 模型算例

分別采用常規(guī)有限差分彈性波逆時(shí)偏移算法和本文基于FCT校正的有限差分逆時(shí)偏移算法對(duì)該模型的正演數(shù)據(jù)進(jìn)行逆時(shí)偏移處理，所用的主要參數(shù)為：空間網(wǎng)格大小5 m×5 m，時(shí)間網(wǎng)格大小0.5 ms，炮點(diǎn)波場(chǎng)采用主頻80 Hz的Ricker子波進(jìn)行正向延拓，以炮孔徑作為偏移孔徑。圖7、8分別為上述兩種方法的縱橫波偏移結(jié)果，可以看出兩種算法均能準(zhǔn)確恢復(fù)該模型的速度層位，但常規(guī)有限差分算法偏移剖面中出現(xiàn)虛假同向軸，本文算法偏移剖面虛假同相軸不存在；本文算法的疊加剖面中同相軸比常規(guī)有限差分算法更加清晰，剖面整體精度比常規(guī)有限差分算法高。

圖6 縱橫波速度模型Fig .6 P and S velocity models

圖7 常規(guī)方法的逆時(shí)偏移成像剖面Fig .7 Sections of RTM based on conventional methods

圖8 本文算法的逆時(shí)偏移成像剖面Fig .8 Sections of RTM in this paper

3.2 實(shí)測(cè)資料處理試驗(yàn)

利用本文方法對(duì)某工區(qū)的OBN數(shù)據(jù)進(jìn)行逆時(shí)偏移處理。原始數(shù)據(jù)為四分量二維海底地震數(shù)據(jù)，本文以各向同性為假設(shè)，只取其中的x分量和z分量進(jìn)行逆時(shí)偏移，偏移所用的深度域縱波層速度模型由反射縱波的層析成像得到，橫波層速度模型則由轉(zhuǎn)換波速度分析結(jié)果與縱波模型換算得到，采用本文方法得到的縱橫波逆時(shí)偏移結(jié)果見圖9。

從圖9可以看出縱、橫波剖面中的構(gòu)造層位基本一致，但在相位和構(gòu)造細(xì)節(jié)方面存在一些差異。分析認(rèn)為，造成這些差異的可能原因有：①逆時(shí)偏移前建立的深度域縱橫波層速度模型存在較大誤差；②對(duì)于相同地層，反射縱波和轉(zhuǎn)換橫波具有不同響應(yīng)；③橫波速度小于縱波速度，在深度域中表現(xiàn)為橫波波長(zhǎng)小于縱波波長(zhǎng)，即深度域中橫波分辨率大于縱波分辨率。

圖9 實(shí)測(cè)OBN資料逆時(shí)偏移疊加剖面Fig .9 Measured OBN RTM stack sections

4 結(jié)論

1) 數(shù)值頻散是差分離散求解波動(dòng)方程產(chǎn)生的虛假波場(chǎng)，是一種計(jì)算誤差而非介質(zhì)本身的彈性響應(yīng)，相同條件下地震波的高頻信息越豐富，這種誤差越大。OBN勘探由于其特殊的接收方式而含有豐富的高頻信息，對(duì)OBN數(shù)據(jù)進(jìn)行基于有限差分的逆時(shí)偏移處理必須采用適當(dāng)?shù)姆椒ㄏ@種誤差。

2) 本文提出的基于FCT校正的OBN資料彈性波逆時(shí)偏移算法能夠消除由數(shù)值頻散產(chǎn)生的縱橫波逆時(shí)偏移成像誤差，且算法在頻散產(chǎn)生后再進(jìn)行波場(chǎng)校正，沒有改變現(xiàn)有的彈性波逆時(shí)偏移技術(shù)框架，只須在已有的軟件模塊中插入對(duì)應(yīng)的頻散校正代碼即可。因此，本文算法不僅適用于一階速度-漲縮-旋轉(zhuǎn)彈性波方程，而且適用于其他類型的地震波方程。

3) 本文提出的基于FCT校正的OBN資料彈性波逆時(shí)偏移算法屬于矢量波場(chǎng)處理方法范疇，在對(duì)多分量地震資料處理前無(wú)需進(jìn)行縱橫波解耦，避開了波場(chǎng)分離的難題，理論上具有更高的振幅保真性和成像精度。