拖纜與OBN資料聯(lián)合成像域最小二乘逆時(shí)偏移成像

楊華臣， 張建中,2*

1 海底科學(xué)與探測(cè)技術(shù)教育部重點(diǎn)實(shí)驗(yàn)室， 中國(guó)海洋大學(xué)海洋地球科學(xué)學(xué)院， 青島 266100 2 青島海洋科學(xué)與技術(shù)國(guó)家實(shí)驗(yàn)室， 海洋礦產(chǎn)資源評(píng)價(jià)與探測(cè)技術(shù)功能實(shí)驗(yàn)室， 青島 266061

0 引言

拖纜地震是海洋油氣勘探的主要技術(shù)(Rickett，2003；劉學(xué)建和劉伊克，2016；葉月明等，2019；張瑞等，2020)，具有施工效率高、采集成本低和覆蓋次數(shù)高等優(yōu)點(diǎn).在地質(zhì)結(jié)構(gòu)較為簡(jiǎn)單的地區(qū)，僅使用拖纜地震資料就可以獲得高質(zhì)量的偏移成像剖面(韓復(fù)興等，2015；張力起等，2019；段心標(biāo)等，2020).然而，在存在高陡構(gòu)造、古潛山等深部地質(zhì)構(gòu)造復(fù)雜的區(qū)域，地震波場(chǎng)復(fù)雜，因拖纜的長(zhǎng)度有限，難以記錄深部地層的反射波信號(hào)，使得深部地層的偏移成像效果不佳，無法正確刻畫深部的復(fù)雜地質(zhì)構(gòu)造(Wong et al.，2012；Zhang and Schuster，2013；Li et al.，2017；張振波和羅偉，2021).

海底地震是新發(fā)展起來的海洋地震勘探技術(shù)(Sollid and Ursin，2003；Zhang，2014；趙維娜等，2019).震源在近海面被拖拽著激發(fā)，海底地震節(jié)點(diǎn)(OBN)放在海底接收，其偏移距可以達(dá)到數(shù)十千米，能夠記錄到深部高陡構(gòu)造的反射波信號(hào)，對(duì)深部的復(fù)雜地質(zhì)構(gòu)造進(jìn)行有效成像(Zelt et al.，1998；Charvis and Operto，1999；Wong et al.，2011；Wang et al.，2012；Cheng et al.，2018).然而，OBN資料采集的施工周期較長(zhǎng)、成本高.為了減少OBN資料采集的成本，通常布設(shè)OBN的間距較大(Wong et al.，2012；鐘廣見等，2014；趙維娜等，2019，2020)，導(dǎo)致OBN資料對(duì)海底介質(zhì)的覆蓋次數(shù)較低，且在空間上的照明分布非常不均勻(Zhang et al.，2019)，在某些淺部地層甚至沒有覆蓋，從而使得OBN資料的偏移成像剖面，特別是對(duì)淺層的成像質(zhì)量較差(Zelt et al.，1998；Wong et al.，2011，2012；鐘廣見等，2014).

一般地，拖纜地震資料對(duì)淺部地層的覆蓋次數(shù)較高，OBN資料對(duì)深部地層的覆蓋次數(shù)相對(duì)較高，這兩種地震資料具有優(yōu)勢(shì)互補(bǔ)性.我們聯(lián)合使用拖纜與OBN資料進(jìn)行波形反演，獲得了比使用其中一種資料更好的效果(Yang and Zhang，2019).Wong等(2012)利用一個(gè)理論模型合成資料，說明了存在鉆井平臺(tái)等障礙物時(shí)拖纜與OBN資料聯(lián)合偏移成像的優(yōu)勢(shì).可見，聯(lián)合使用拖纜與OBN資料進(jìn)行偏移成像，可以有機(jī)地發(fā)揮兩種資料的優(yōu)點(diǎn)，克服缺點(diǎn)，獲得更高質(zhì)量的偏移成像效果.因此，本文研究一種拖纜與OBN資料聯(lián)合的最小二乘逆時(shí)偏移成像方法(LSRTM).

此外，常規(guī)數(shù)據(jù)域LSRTM需要通過多次偏移與反偏移來估計(jì)Hessian矩陣(Ji，2009；Di et al.，2016；王曉毅等，2021)，計(jì)算量大，聯(lián)合使用拖纜與OBN資料時(shí)該問題更嚴(yán)重.本文的拖纜與OBN資料聯(lián)合的LSRTM，通過在成像域計(jì)算點(diǎn)擴(kuò)散函數(shù)來估計(jì)Hessian矩陣，無需進(jìn)行多次偏移和反偏移運(yùn)算，極大地提高了LSRTM的效率.將本文方法應(yīng)用于復(fù)雜構(gòu)造模型的合成資料，驗(yàn)證了聯(lián)合LSRTM對(duì)于諸如高陡構(gòu)造等復(fù)雜地質(zhì)構(gòu)造的成像優(yōu)勢(shì).

圖1 傾斜界面反射波射線路徑示意圖五角星表示震源，三角形表示檢波點(diǎn)，圓點(diǎn)表示反射點(diǎn). α、β、h和x分別表示反射界面傾角、入射角、反射點(diǎn)深度和偏移距.Fig.1 The ray path of the reflection wave with a sloping reflectorThe star, triangle, and dot indicate the shot, receiver, and the reflection point, individually. The α, β, h, x indicate the dip angle of the reflector, incident angle, depth of the reflection point, and offset between the shot and receiver, respectively.

1 拖纜資料偏移成像的問題

因拖纜的長(zhǎng)度有限，拖纜地震采集系統(tǒng)的偏移距不夠大，難以記錄高陡構(gòu)造的反射波信號(hào).圖1是一個(gè)傾斜界面及反射波射線路徑示意圖.S為炮點(diǎn)，R為接收點(diǎn)，反射點(diǎn)位于O點(diǎn)，界面傾角為α，偏移距為x，反射點(diǎn)深度為h，入射角或反射角為β.這些參數(shù)之間存在如下關(guān)系：

x=h[tan(β+α)+tan(β-α)]，

(1)

當(dāng)α=β=0時(shí)，x=0，此時(shí)界面水平，反射波垂直入射到反射界面并返回，為水平反射面的自激自收情形；當(dāng)α≠0，β=0時(shí)，x=0，為傾斜反射面的自激自收情形；當(dāng)α+β=90°時(shí)，x趨于無窮大，反射波從反射點(diǎn)沿水平方向向左傳播，在觀測(cè)面不能接收到反射波信號(hào).因此，在觀測(cè)面接收到的反射波的反射角β必介于0°與90°-α之間，且α越大，x將越大.

圖2是用式(1)計(jì)算的偏移距x隨反射界面傾角α、反射點(diǎn)深度h以及反射角β的變化曲線.圖2a表明，當(dāng)h=4 km，β=30°，α<40°時(shí)，偏移距x隨α增加而逐漸從約4.6 km增加到約10.3 km；當(dāng)α>40°，x隨α增加而急劇增加.α=25°對(duì)應(yīng)的偏移距x=6.06 km.這說明，當(dāng)?shù)貙觾A角大于25°時(shí)，偏移距x小于6 km的拖纜檢波器接收不到反射波.從圖2b看出，當(dāng)α和β一定時(shí)，x隨h增加而線性地增加.當(dāng)α=50°，β=30°時(shí)，要接收到來自深度h=1.2 km的反射波，偏移距x不能小于6.4 km.圖2c表明，當(dāng)α=50°，h=4 km時(shí)，x隨β增加而增加，且增加的速率逐漸變大.β越大，炮點(diǎn)離反射點(diǎn)水平距離越大，接收到反射波的偏移距也越大.當(dāng)α=50°，h=4 km時(shí)，β=16°對(duì)應(yīng)的偏移距x=6.29 km.這表明，反射角大于16°的反射波不能被長(zhǎng)度為6 km的拖纜中的檢波器記錄到.綜上，反射界面傾角α、反射點(diǎn)深度h或反射角β增加都會(huì)使偏移距x增加.因此，只有在偏移距達(dá)到一定程度時(shí)，才能記錄到深部高陡構(gòu)造的反射波信號(hào).然而，常規(guī)的拖纜長(zhǎng)度通常在6 km左右(馬光克等，2019；張振波和羅偉，2021)，難以記錄深部高陡構(gòu)造的反射波信號(hào)，導(dǎo)致其偏移成像孔徑小、成像效果差.

2 聯(lián)合成像方法

2.1 基本原理

為了聯(lián)合使用拖纜與OBN資料進(jìn)行LSRTM，獲得地下介質(zhì)的反射率成像剖面，定義聯(lián)合LSRTM的目標(biāo)函數(shù)為：

(2)

為了獲得目標(biāo)函數(shù)的極小值，令其對(duì)反射率模型的一階導(dǎo)數(shù)為零，可得：

(3)

其中，LT表示偏移算子，LTd表示觀測(cè)地震記錄的逆時(shí)偏移成像剖面，LTL為Hessian矩陣.

采用聲波方程進(jìn)行偏移與反偏移(劉夢(mèng)麗等，2018；柯璇等，2019).為了壓制逆時(shí)偏移成像剖面中的低頻背景噪聲，采用上下行波分離的互相關(guān)成像條件(Liu et al.，2011)：

+pu(x,t|xs)qd(x,t|xs)]dt,

(4)

其中，I表示逆時(shí)偏移成像剖面，s表示震源序號(hào)，T表示地震記錄最大采樣時(shí)間，pu和pd分別表示震源正向傳播波場(chǎng)的上、下行波，qu和qd分別表示檢波點(diǎn)反向傳播波場(chǎng)的上、下行波，x表示空間坐標(biāo)，xs表示震源的空間坐標(biāo)，t表示時(shí)間.在波數(shù)域進(jìn)行上下行波分離，具體方法請(qǐng)參見Liu等(2011).

將式(3)改寫為下列形式：

Hm=mmig，

(5)

式中：

(6)

(7)

(8)

(9)

通過式(6)和式(7)分別求得H(H的計(jì)算方法在下一節(jié)介紹)和mmig后，采用阻尼最小二乘算法求解式(5)，得到反射率模型m.

2.2 計(jì)算Hessian矩陣

式(5)中Hessian矩陣很大，直接計(jì)算和存儲(chǔ)Hessian矩陣是非常困難的(Tang，2009；李振春等，2014).傳統(tǒng)的數(shù)據(jù)域LSRTM通過多次偏移和反偏移來估計(jì)式(5)中的Hessian矩陣(劉夢(mèng)麗等，2018；柯璇等，2019).例如，劉夢(mèng)麗等(2018)和柯璇等(2019)分別在模型試驗(yàn)中進(jìn)行了30次和50次偏移和反偏移運(yùn)算，從而使傳統(tǒng)數(shù)據(jù)域LSRTM極其耗費(fèi)時(shí)間.本文通過計(jì)算點(diǎn)擴(kuò)散函數(shù)來求取Hessian矩陣，僅需要進(jìn)行1次偏移和反偏移運(yùn)算.

式(5)表明，如果已知反射率模型和對(duì)應(yīng)的逆時(shí)偏移成像剖面，則Hessian矩陣可以作為未知量求出.因此，如果給定一個(gè)反射率模型，并通過逆時(shí)偏移成像方法獲得其對(duì)應(yīng)的像，則可以通過阻尼最小二乘算法等求解式(5)，獲得Hessian矩陣.

特別地，當(dāng)反射率模型m中僅第i點(diǎn)存在擾動(dòng)，即反射率模型中僅第i個(gè)元素為非零值，逆時(shí)偏移成像后該點(diǎn)對(duì)應(yīng)的像為Hessian矩陣中第i列的元素.此時(shí)，式(5)可以表示為：

(10)

其中，n表示反射率模型中元素的總個(gè)數(shù)，mi表示僅第i個(gè)點(diǎn)存在擾動(dòng)時(shí)的反射率模型，即：

(11)

式中，η是一個(gè)給定的小量，本文設(shè)置為0.01.

將式(11)代入式(10)，進(jìn)行矩陣相乘，可得：

(12)

式中，等號(hào)右端是mi對(duì)應(yīng)的像，可以使用逆時(shí)偏移成像方法得到.因此，通過式(12)就可以得到Hessian矩陣H中第i列的元素.該列元素體現(xiàn)了第i個(gè)點(diǎn)的反射率與其像之間的映射關(guān)系，即該點(diǎn)的點(diǎn)擴(kuò)散函數(shù).為了簡(jiǎn)潔，將式(12)記為：

(13)

由于地震波的頻帶有限，第i個(gè)點(diǎn)的像僅存在于該點(diǎn)鄰域內(nèi)(以第i個(gè)點(diǎn)的位置為中心兩倍波長(zhǎng)為半徑的區(qū)域內(nèi))，在遠(yuǎn)離該點(diǎn)的區(qū)域其成像值為零.因此，可以在反射率模型m中同時(shí)設(shè)置若干個(gè)空間上具有一定間距的散射點(diǎn)，并給定其擾動(dòng)，再通過逆時(shí)偏移成像方法得到它們的像，即它們的點(diǎn)擴(kuò)散函數(shù).最后，通過插值方法得到反射率模型向量m中任意一點(diǎn)的點(diǎn)擴(kuò)散函數(shù).求得每一點(diǎn)的點(diǎn)擴(kuò)散函數(shù)后，就獲得了Hessain矩陣中每一列的元素.此時(shí)，式(5)中的Hessian矩陣可以表示為：

(14)

當(dāng)速度模型存在速度突變界面時(shí)，逆時(shí)偏移成像剖面中通常會(huì)存在一定的低頻背景噪聲(楊仁虎，2021；諸峰等，2022).這些低頻背景噪聲與點(diǎn)擴(kuò)散函數(shù)存在明顯的差異.點(diǎn)擴(kuò)散函數(shù)在空間上分布相對(duì)較為稀疏，表現(xiàn)出高波數(shù)的特點(diǎn).因此，通過對(duì)逆時(shí)偏移成像剖面進(jìn)行高通濾波，就可以實(shí)現(xiàn)對(duì)低頻背景噪聲的壓制，得到準(zhǔn)確的點(diǎn)擴(kuò)散函數(shù).

3 模型資料試驗(yàn)分析

3.1 資料合成

基于Marmousi2模型，建立了2個(gè)速度結(jié)構(gòu)相同但模型尺寸不同的速度模型，如圖3所示(模型左右兩側(cè)12.5 km區(qū)域內(nèi)拖纜與OBS資料的覆蓋次數(shù)較少，故僅顯示了模型中間部分).兩個(gè)模型在Z方向的大小均為5.375 km，在X方向的大小分別為33.75 km和112.5 km.模型的頂部為海水層，其厚度為0.25 km、速度為1.5 km·s-1.海水層之下0.25～1.5 km深度區(qū)間分布5個(gè)水平勻速層，速度分別為1.7 km·s-1、1.9 km·s-1、2.4 km·s-1、2.8 km·s-1和3 km·s-1.在Z方向1.5 km以下區(qū)域存在從Marmousi2模型截取的一系列高陡地層和斷層.整體上，圖3a模型中深部(Z>1.5 km)地層傾角主要介于12°～46°，圖3b模型對(duì)應(yīng)的地層傾角大致介于1.2°～6°.模型的最小和最大速度分別為1.5 km·s-1和5.56 km·s-1.

圖3 速度模型(a) 小尺寸的理論速度模型；其中三角形和圓點(diǎn)分別代表圖4的OBN和炮點(diǎn)的位置； (b) 大尺寸的理論速度模型； (c)和(d) 分別是圖(a)和(b)所示模型的平滑速度模型.Fig.3 Velocity models(a) Theoretical velocity model with smaller size； The triangle and dot indicate the OBN and shot location in Fig.4, respectively； (b) Theoretical velocity model with bigger size； (c) and (d) Indicate the smooth version of velocity models shown in (a) and (b), respectively.

在合成拖纜地震資料時(shí)，將拖纜長(zhǎng)度設(shè)置為12 km，道間距為12.5 m，共960道.震源以100 m的間距從模型左側(cè)X=12.5 km處開始激發(fā)直到模型的最右側(cè).拖纜位于震源的左側(cè)，并且隨著震源的移動(dòng)而同步移動(dòng).拖纜中的水聽器與震源之間的最小偏移距為50 m，最大偏移距為12.05 km.

在合成OBN資料時(shí)，將OBN分別以0.4 km、0.8 km和1.6 km的間距均勻布設(shè)在海底(Z=0.25 km)，震源以100 m的間距從模型最左側(cè)X=0 km處開始激發(fā)直到模型的最右側(cè).

在合成上述兩種數(shù)據(jù)時(shí)，均采用邊長(zhǎng)12.5 m的正方形網(wǎng)格離散速度模型、1 ms的時(shí)間采樣率和14 s的時(shí)間采樣長(zhǎng)度.使用時(shí)間二階、空間十二階精度的有限差分方法模擬地震波場(chǎng)(Alford et al.，1974；Liu and Sen，2009).采用PML邊界條件壓制模型邊界處的反射波場(chǎng)(王維紅等，2013).考慮到實(shí)際資料的頻帶范圍不同，分別采用20 Hz和15 Hz的雷克子波作為震源合成拖纜與OBN資料.

圖4a是震源位于海面X=17.4 km處(圖3a中圓點(diǎn)所示位置)的拖纜資料炮集記錄.圖4b是位于海底X=16.8 km處(圖3a三角形所示位置)的OBN資料的共接收點(diǎn)道集記錄.可以看出，拖纜與OBN資料存在兩個(gè)重要的差別.第一，拖纜只能記錄震源一側(cè)的地震波場(chǎng)，而OBN可以記錄震源兩側(cè)的地震波場(chǎng).第二，拖纜只能記錄偏移距在12.05 km以內(nèi)的地震波場(chǎng)，而OBN能夠記錄到更大偏移距的地震波場(chǎng)信息.這使得OBN資料具有更大的觀測(cè)孔徑，更容易接收到深部諸如高陡構(gòu)造等復(fù)雜地質(zhì)構(gòu)造的反射波信號(hào).

圖4 合成地震記錄(a) 炮點(diǎn)位于圖3a中圓點(diǎn)所示位置處的拖纜地震記錄； (b) OBN位于圖3a中三角形所示位置處的地震記錄.Fig.4 Synthetic seismic data(a) Towed-streamer seismic data with the shot located at the dot in Fig.3a； (b) OBN data recorded at the triangle in Fig.3a.

3.2 拖纜資料成像

圖3a、b所示速度模型具有相同的地質(zhì)結(jié)構(gòu)，但中深部的地層傾角不同.分別用125 m×125 m和1250 m×125 m窗口對(duì)這兩個(gè)速度模型進(jìn)行10次平滑處理，得到平滑速度模型，如圖3c、d所示.基于該平滑速度模型，分別對(duì)合成的偏移距6 km以內(nèi)的拖纜地震資料進(jìn)行逆時(shí)偏移成像，結(jié)果如圖5a、b所示.可以看出，圖5a、b中淺部水平同相軸都非常連續(xù)、清晰.然而，圖5a對(duì)中深部的高陡地層和斷層反映不清晰、不連續(xù)，甚至在局部區(qū)域不成像.相比之下，圖5b對(duì)中深部的地層和斷層反映得更清晰、同相軸連續(xù)性更好.該測(cè)試表明，采用常規(guī)的拖纜地震勘探觀測(cè)系統(tǒng)，對(duì)傾斜程度較低的速度界面可以進(jìn)行較好地成像，而對(duì)中深部高陡地層和斷層的成像質(zhì)量較差.

圖5 拖纜資料的逆時(shí)偏移成像剖面(a)和(b) 分別使用圖3a、b所示模型合成的偏移距6 km以內(nèi)的拖纜地震記錄； (c) 使用圖3a所示模型合成的偏移距12.05 km以內(nèi)的拖纜地震記錄.Fig.5 RTM images of the towed-streamer seismic data(a) and (b) Using the towed-streamer seismic data with the offset between 0～6 km generated with velocity models shown in Fig.3a,b, respectively； (c) Using the towed-streamer seismic data with the offset between 0～12.05 km generated with the velocity model shown in Fig.3a.

為了分析增加拖纜長(zhǎng)度對(duì)中深部高陡地層和斷層偏移成像的效果，對(duì)圖3a所示模型合成的偏移距12.05 km以內(nèi)全部拖纜地震資料進(jìn)行逆時(shí)偏移成像，結(jié)果如圖5c所示.對(duì)比圖5a、c可以看出，偏移距長(zhǎng)度從6 km增加到12.05 km，中深部高陡地層和斷層的成像質(zhì)量有一定提高，特別是在黑色圓圈所示范圍內(nèi).然而，相比于圖5b，模型中深部高陡地層和斷層反映得仍然不夠清晰、同相軸不夠連續(xù)，如黑色箭頭所指.這表明增加有限的拖纜長(zhǎng)度，對(duì)中深部高陡地層和斷層成像質(zhì)量的改善很有限.

3.3 OBN資料成像

基于圖3a對(duì)應(yīng)的平滑速度模型，對(duì)合成的OBN間距為0.4 km、0.8 km和1.6 km的資料分別進(jìn)行逆時(shí)偏移成像，結(jié)果如圖6所示.可以看出，中深部的高陡地層和斷層均得到了良好地成像.這表明OBN資料能夠記錄到模型中深部大角度傾斜界面的反射波信號(hào)，從而能對(duì)復(fù)雜的高陡地層和斷層進(jìn)行成像.而且，適度增加OBN間距對(duì)中深部高陡構(gòu)造成像質(zhì)量的影響相對(duì)較小，說明利用稀疏的OBN資料對(duì)深部復(fù)雜構(gòu)造進(jìn)行成像的潛力.但是，OBN間距的增加會(huì)降低OBN資料對(duì)海底介質(zhì)的覆蓋次數(shù)，其逆時(shí)偏移成像剖面中有效同相軸的振幅會(huì)隨著OBN間距的增加而減小，導(dǎo)致偏移成像剖面的信噪比降低.整體上，OBN資料對(duì)淺部地層成像效果不佳，且OBN間距越大，偏移噪聲越嚴(yán)重.這種偏移噪聲主要是淺層介質(zhì)大角度的反射波和折射波在逆時(shí)延拓過程中與震源正傳波場(chǎng)互相關(guān)產(chǎn)生的.

圖6 基于圖3a所示速度模型合成的OBN資料逆時(shí)偏移成像剖面(a) OBN間距為0.4 km； (b) OBN間距為0.8 km； (c) OBN間距為1.6 km.Fig.6 RTM images of the OBN data generated using the velocity model shown in Fig.3a(a) The interval of OBN is 0.4 km； (b) The interval of OBN is 0.8 km； (c) The interval of OBN is 1.6 km.

3.4 聯(lián)合成像

采用本文方法基于圖3a所示速度模型合成的偏移距6 km以內(nèi)的拖纜地震記錄和接收間距分別為0.4 km、0.8 km和1.6 km的OBN資料進(jìn)行聯(lián)合偏移成像.

圖7a、b分別是計(jì)算的拖纜和接收間距為0.4 km的OBN資料對(duì)應(yīng)的點(diǎn)擴(kuò)散函數(shù).圖中的點(diǎn)狀成像區(qū)表示一個(gè)點(diǎn)擴(kuò)散函數(shù)，反映了該位置處的散射點(diǎn)與逆時(shí)偏移成像后的像之間的映射關(guān)系，即Hessian矩陣中的一列元素.基于這些離散點(diǎn)對(duì)應(yīng)的點(diǎn)擴(kuò)散函數(shù)，通過插值方法求得空間任意點(diǎn)對(duì)應(yīng)的點(diǎn)擴(kuò)散函數(shù)值，從而得到Hessian矩陣.需要注意的是，在圖7a的淺部和圖7b的底部存在明顯的低頻背景噪聲.這些噪聲會(huì)嚴(yán)重影響點(diǎn)擴(kuò)散函數(shù)的幅值與形態(tài).因此，對(duì)點(diǎn)擴(kuò)散函數(shù)進(jìn)行了高通濾波，去除低頻的背景噪聲.圖7c、d分別表示噪聲壓制后拖纜和OBN資料對(duì)應(yīng)的點(diǎn)擴(kuò)散函數(shù).

圖7 點(diǎn)擴(kuò)散函數(shù)(a) 拖纜資料的點(diǎn)擴(kuò)散函數(shù)； (b) 接收間距為0.4 km的OBN資料的點(diǎn)擴(kuò)散函數(shù)； (c) 去噪后拖纜資料的點(diǎn)擴(kuò)散函數(shù)； (d) 去噪后OBN資料的點(diǎn)擴(kuò)散函數(shù).Fig.7 Point spread functions(a) Point spread functions of the towed-streamer seismic data； (b) Point spread functions of the OBN data with an interval of 0.4 km； (c) Denoised point spread functions of the towed-streamer seismic data； (d) Denoised point spread functions of the OBN data.

圖8是兩種資料聯(lián)合偏移與單一資料偏移結(jié)果的對(duì)比.其中，圖8a是理論垂直反射率模型，采用式(15)進(jìn)行計(jì)算(朱光，2016；方修政等，2018)：

(15)

圖8 不同資料LSRTM成像剖面對(duì)比(a) 理論垂直反射率模型；其中黑線代表圖9中成像曲線的位置； (b) 僅使用拖纜地震資料； (c) 僅使用接收間距為0.4 km的OBN資料； (d) 聯(lián)合使用拖纜和接收間距為0.4 km的OBN資料； (e) 聯(lián)合使用拖纜和接收間距為0.8 km的OBN資料； (f) 聯(lián)合使用拖纜和接收間距為1.6 km的OBN資料.Fig.8 Comparison of LSRTM images using different seismic data(a) Theoretical vertical reflectivity model；The black line indicates the location of image curves shown in Fig.9； (b) Using only the towed-streamer seismic data； (c) Using only the OBN data with and interval of 0.4 km； (d) Using the towed-streamer and OBN data with an interval of 0.4 km； (e) Using the towed-streamer and OBN data with an interval of 0.8 km； (f) Using the towed-streamer and OBN data with an interval of 1.6 km.

其中，m表示反射率，x表示空間坐標(biāo)，δv表示偏移速度與理論速度之差，v0表示偏移速度.圖8b、c分別是僅使用拖纜地震記錄和僅使用接收間距為0.4 km的OBN記錄的LSRTM成像剖面，圖8d、e、f分別是拖纜資料與接收間距為0.4 km、0.8 km和1.6 km的OBN資料聯(lián)合LSRTM成像剖面.對(duì)比圖8b和圖5a可以看出，拖纜資料LSRTM成像剖面對(duì)中深部大角度傾斜界面成像仍然不清楚.對(duì)比圖8c和圖6a可以看出，使用OBN資料LSRTM成像剖面的中深部同相軸的振幅得到了相對(duì)的增強(qiáng).然而，淺部的偏移噪聲整體上仍然非常明顯.

對(duì)比單一地震資料(圖8b、c)和兩種資料(圖8d)聯(lián)合后的LSRTM成像剖面可以看出，聯(lián)合偏移成像剖面中，淺部地層的成像質(zhì)量與拖纜資料偏移成像剖面的質(zhì)量相當(dāng)，深部地層的成像質(zhì)量與OBN資料偏移成像的成像質(zhì)量相當(dāng).對(duì)于模型中部高陡地層和斷層，聯(lián)合成像剖面同相軸的連續(xù)性和清晰度高于僅使用拖纜資料的成像剖面，其噪聲弱于僅使用OBN資料的成像剖面.

對(duì)比圖8e和圖6b可以看出，僅使用OBN資料逆時(shí)偏移成像剖面中存在明顯的偏移噪聲，聯(lián)合使用兩種資料對(duì)偏移噪聲有著很好地壓制效果，其LSRTM成像剖面中幾乎看不到明顯的偏移噪聲.對(duì)比圖8e、f可以發(fā)現(xiàn)，當(dāng)OBN間距由0.8 km增加到1.6 km后，雖然LSRTM成像剖面中可以看到較弱的偏移噪聲，但是高陡地層和斷層的成像質(zhì)量幾乎相當(dāng).這表明聯(lián)合偏移成像對(duì)稀疏的OBN資料逆時(shí)偏移成像噪聲有著很好地壓制效果.

綜上所述，拖纜與OBN資料聯(lián)合成像整體上優(yōu)于單一資料成像，可以得到更好的從淺到深的偏移成像效果.

圖9是X=16.25 km處(圖8a中黑線)OBN資料的RTM、LSRTM成像曲線和理論反射率曲線及其頻譜.RTM成像值(圖9c)在淺部(Z<2 km)相對(duì)較大，而在中深部相對(duì)較小.相比于RTM成像值，LSRTM在深部的成像值(圖9e)相對(duì)更大.在Z=3.35 km和5.14 km左右，理論反射率曲線(圖9a)具有一個(gè)明顯的波谷，RTM成像曲線在對(duì)應(yīng)位置的波動(dòng)非常微弱，而LSRTM成像曲線在對(duì)應(yīng)位置具有一個(gè)明顯的波谷.然而，在Z=4 km左右，理論反射率曲線存在一個(gè)明顯的波谷，而RTM成像曲線和LSRTM成像曲線對(duì)其均沒有明顯的反映.這主要是由于該深度處地層傾角極大，RTM成像值非常小，如圖9c所示.圖9b、d、f分別是理論反射率曲線、RTM成像曲線和LSRTM成像曲線的波數(shù)譜.RTM成像曲線和理論反射率曲線的波數(shù)譜的相關(guān)系數(shù)為0.72，而LSRTM成像曲線和理論反射率曲線的波數(shù)譜的相關(guān)系數(shù)為0.81.這表明LSRTM成像曲線與理論反射率曲線更接近.

圖9 X=16.25 km(圖8a中黑線)處的成像曲線及其波數(shù)譜(a) 理論反射率曲線； (b)(a) 所示反射率曲線的波數(shù)譜； (c) OBN間距0.4 km時(shí)逆時(shí)偏移成像值； (d)(c) 所示成像值曲線的波數(shù)譜； (e) OBN間距0.4 km時(shí)最小二乘逆時(shí)偏移成像值； (f)(e) 所示成像值曲線的波數(shù)譜.Fig.9 Image curves and their corresponding wave-number spectrums at X=16.25 km indicated by the black line in Fig.8a(a) Theoretical vertical reflectivity curve； (b) The wave-number spectrum of the curve in (a)； (c) The RTM image curve of the OBN data with an interval of 0.4 km； (d) The wave-number spectrum of the curve in (c)； (e) The LSRTM image curve of the OBN data with an interval of 0.4 km； (f) The wave-number spectrum of the curve in (e).

4 結(jié)論

受限于拖纜的長(zhǎng)度，海洋拖纜地震資料對(duì)中深部高陡構(gòu)造的偏移成像質(zhì)量較差，甚至難以成像.稀疏OBN資料對(duì)深部和高陡構(gòu)造成像效果較好，但在中淺部引入偏移噪聲.OBN間距越大，偏移噪聲越嚴(yán)重.本文提出了一種拖纜與OBN資料聯(lián)合成像域LSRTM方法，其中，在成像域通過計(jì)算點(diǎn)擴(kuò)散函數(shù)來估計(jì)Hessian矩陣，避免了常規(guī)數(shù)據(jù)域LSRTM通過多次偏移和反偏移來估計(jì)Hessian矩陣極其耗時(shí)的缺點(diǎn).理論模型測(cè)試表明，該方法把拖纜資料與OBN資料有機(jī)地融合在一起，發(fā)揮了這兩種資料優(yōu)勢(shì)互補(bǔ)的作用，能夠獲得比僅使用單一資料更好的偏移成像效果，是深部復(fù)雜構(gòu)造成像的有效方法.

由于缺乏同一測(cè)線的拖纜和OBN實(shí)際資料，本文僅使用模型的合成資料進(jìn)行實(shí)驗(yàn)和分析.對(duì)于實(shí)測(cè)資料還存在各種噪聲的影響及壓制等問題，需要進(jìn)一步研究.