一階多次波聚焦變換成像

劉學(xué)建， 劉伊克， 胡昊， 謝宋雷

1 中國科學(xué)院地質(zhì)與地球物理研究所 中國科學(xué)院工程地質(zhì)力學(xué)重點(diǎn)實(shí)驗(yàn)室， 北京 10002 2 中國科學(xué)院大學(xué)， 北京 100049 3 中國科學(xué)院南海海洋研究所， 廣州 510301

?

一階多次波聚焦變換成像

劉學(xué)建1,2， 劉伊克1， 胡昊1,2， 謝宋雷3

1 中國科學(xué)院地質(zhì)與地球物理研究所 中國科學(xué)院工程地質(zhì)力學(xué)重點(diǎn)實(shí)驗(yàn)室， 北京 10002 2 中國科學(xué)院大學(xué)， 北京 100049 3 中國科學(xué)院南海海洋研究所， 廣州 510301

將多次波轉(zhuǎn)換成反射波并按傳統(tǒng)反射波偏移算法成像，是多次波成像的一種方法.聚焦變換能準(zhǔn)確的將多次波轉(zhuǎn)換為縱向分辨率更高的新波場記錄，其中一階多次波轉(zhuǎn)換為反射波.本文對聚焦變換提出了兩點(diǎn)改進(jìn)： 1)提出局部聚焦變換，以減小存儲(chǔ)量和計(jì)算量，增強(qiáng)該方法對檢波點(diǎn)隨炮點(diǎn)移動(dòng)的采集數(shù)據(jù)的適應(yīng)性； 2)引入加權(quán)矩陣，理論上證明原始記錄的炮點(diǎn)比檢波點(diǎn)稀疏時(shí)，共檢波點(diǎn)道集域的局部聚焦變換可以將多次波準(zhǔn)確轉(zhuǎn)換成炮點(diǎn)與檢波點(diǎn)有相同采樣頻率的新波場記錄.本文在第一個(gè)數(shù)值實(shí)驗(yàn)中對比了對包含反射波與多次波的原始記錄做局部聚焦變換和直接對預(yù)測的多次波做局部聚焦變換兩種方案，驗(yàn)證了第二種方案轉(zhuǎn)換得到的波場記錄信噪比更高且避免了第一個(gè)方案中切聚焦點(diǎn)這項(xiàng)比較繁雜的工作.第二個(gè)數(shù)值實(shí)驗(yàn)表明：在炮點(diǎn)采樣較為稀疏時(shí)，該方法能有效的將一階多次波轉(zhuǎn)換成反射波;轉(zhuǎn)換的反射波能提供更豐富的波場信息，成像結(jié)果更均衡、在局部有更高的信噪比,以及較高的縱向分辨率.

聚焦變換； 多次波成像； 多次波消除

1 引言

傳統(tǒng)的偏移算法往往只利用反射波對地下結(jié)構(gòu)成像，多次波要在偏移之前的預(yù)處理中盡可能的減掉.實(shí)際上多次波在地下比反射波傳播路徑更長、覆蓋范圍更廣(如圖1).近年來，很多的學(xué)者致力于多次波成像的研究，發(fā)展出多次波地震干涉成像、多次波直接偏移成像、多次波聚焦變換成像等方法(Berkhout, 1993; Berkhout and Verschuur, 1994; Schuster et al., 2004; Guitton, 2002; Berkhout and Verschuur, 2006; Vasconcelos et al.,2008).甚至有的學(xué)者嘗試?yán)枚啻尾ㄟM(jìn)行偏移速度分析(Manuel and Uren, 2001).多次波地震干涉成像、多次波聚焦變換成像都是先將多次波轉(zhuǎn)換成反射波，然后按傳統(tǒng)反射波偏移算法成像.

圖1 反射波與多次波的傳播路徑Fig.1 Propagating paths of the primary and the multiple

地震干涉法最早可以追溯到1968年，Claerbout(1968)在二維模型上進(jìn)行了被動(dòng)震源地震記錄的自相關(guān)實(shí)驗(yàn)，現(xiàn)已經(jīng)取得了廣泛的研究成果(Schuster, 2009).多次波地震干涉成像的簡單流程為：首先將共炮點(diǎn)道集內(nèi)(或共檢波點(diǎn)道集內(nèi))的道記錄兩兩之間做互相關(guān)生成虛擬炮記錄，疊加相同位置處的虛擬炮記錄，多次波轉(zhuǎn)換成反射波和噪聲(Wapenaar and Fokkema, 2006; Schuster, 2009)，然后再按傳統(tǒng)反射波偏移算法成像(Sheng, 2001; He et al., 2007; Jiang et al., 2007).地震干涉法的算法簡單、計(jì)算效率高.然而，地震干涉法不夠精確，生成的虛擬炮記錄有大量的串聲噪聲.另一方面，虛擬炮記錄同相軸的子波是原始地震記錄同相軸子波的互相關(guān).

Berkout等(1994)提出多次波直接偏移成像的思想.修改傳統(tǒng)的單程波偏移或逆時(shí)偏移方法，可以很容易的實(shí)現(xiàn)多次波的直接偏移成像，其原理為：將包含反射波與多次波的原始道記錄代替?zhèn)鹘y(tǒng)偏移算法中的點(diǎn)震源向地下正傳，將多次波(或原始道記錄)代替?zhèn)鹘y(tǒng)偏移算法中的反射波向地下反傳,同一位置處的正傳波場與反傳波場互相關(guān)成像(Guitton, 2002; Muijs, 2007; Liu et al., 2011a, 2011b; Lu et al., 2011; Wang et al., 2014).該方法不需要提前將多次波轉(zhuǎn)換為反射波，且有與傳統(tǒng)反射波偏移方法相同的偏移計(jì)算效率.然而其在成像域產(chǎn)生的串聲噪聲和結(jié)構(gòu)性假象較難消除.

聚焦變換從基于波動(dòng)方程的表面多次波消減方法SRME(surface-related multiples elimination) (Verschuur et al., 1992; Berkhout et al., 1997; Verschuur and Berkhout,1997; 李鵬等, 2007)發(fā)展而來，由Berkhout和Verschuur(2003)首次提出,已應(yīng)用到多次波聚焦變換消除、多次波成像和地震道插值等方面(Berkhout et al., 2004, 2006; Verschuur and Berkhout, 2005; Groenestijn and Verschuur, 2006).該方法將原始波場記錄中的反射波能量聚焦到時(shí)間為零的點(diǎn)周圍，準(zhǔn)確的將一階多次波轉(zhuǎn)換為反射波，二階多次波轉(zhuǎn)換為一階多次波，依次類推.在共炮點(diǎn)道集，反射波能量的聚焦點(diǎn)為炮點(diǎn)；在共檢波點(diǎn)道集，反射波能量的聚焦點(diǎn)為檢波點(diǎn).聚焦變換是最小二乘意義上的變換，多次波轉(zhuǎn)換的新波場記錄具有很高的信噪比.切去聚焦點(diǎn)周圍的能量，然后將新波場記錄中的多次波(由二階及更高階的多次波轉(zhuǎn)換而來)減去，得到的一階多次波轉(zhuǎn)換的反射波.將轉(zhuǎn)換的反射波成像，實(shí)質(zhì)上是對一階多次波的反射位置成像.

本文提出了局部聚焦變換，增強(qiáng)了該方法對檢波點(diǎn)隨炮點(diǎn)移動(dòng)的采集數(shù)據(jù)的適應(yīng)性.引入加權(quán)矩陣后，本文驗(yàn)證了：當(dāng)原始波場記錄的炮點(diǎn)采樣相對稀疏時(shí)，共檢波點(diǎn)道集域的局部聚焦變換能將多次波準(zhǔn)確轉(zhuǎn)換成炮點(diǎn)與檢波點(diǎn)有相同采樣頻率的新波場記錄.傳統(tǒng)反射波成像算法，我們采用的是對速度模型的橫向變化及陡傾角適應(yīng)性強(qiáng)的逆時(shí)偏移(Baysal et al., 1983; 胡昊等, 2013).

2 一階多次波聚焦變換成像原理

2.1 聚焦變換基本原理

在炮點(diǎn)與檢波點(diǎn)重合的觀測系統(tǒng)上，用單頻矩陣表示全波場的波場記錄，波場矩陣的列存儲(chǔ)共炮點(diǎn)或共檢波點(diǎn)道集記錄(如圖2).忽略算符細(xì)節(jié)，用ΔP表示反射波波場記錄，M表示多次波波場記錄，P= ΔP+M=ΔXS+ΔXR∩ΔXS

+(ΔXR∩)2ΔXS+…

圖2 頻率域數(shù)據(jù)矩陣的結(jié)構(gòu)(數(shù)據(jù)矩陣ΔP(xr, xs , ωi)的行列可以互換)Fig.2 Structure of data matrix in frequency domain, row and column of data matrix ΔP(xr, xs , ωi) are interchangeable

P表示包含反射波和多次波的原始波場記錄，根據(jù)反饋迭代模型(Verschuur et al., 1992; Berkhout and Verschuur, 2003, 2006) (圖3)：

圖3 反饋迭代模型Fig.3 The feedback model

ΔP=ΔXS,

(1)

(2)

(3)

其中：對角矩陣S為震源特性矩陣，對角元素為頻率域的震源子波；ΔX表示脈沖震源在地下的一次響應(yīng)；R∩為自由表面反射系數(shù)特性矩陣.

表面算符A記為：

(4)

引入表面算符，公式(2)、(3)改寫為：

M=ΔPAΔP+(ΔPA)2ΔP+…

=ΔPA(ΔP+ΔPAΔP+…)

=ΔPAP′,

(5)

P=ΔP+ΔPAΔP+(ΔPA)2ΔP+…

=ΔP+ΔPA(ΔP+ΔPAΔP+…)

=ΔP+ΔPAP′.

(6)

顯而易見，P′為缺少原始波場記錄P中最高階多次波的波場記錄.當(dāng)測線在自由表面之下時(shí)，自由表面的反射系數(shù)特性(R∩)包含在波場記錄中，此時(shí)表面算符A記為：

(7)

聚焦變換定義為：

(8)

聚焦變換最小二乘意義上的穩(wěn)定表達(dá)式為：

(9)

上標(biāo)H表示復(fù)數(shù)矩陣的共軛轉(zhuǎn)置，ε較為保證矩陣求逆計(jì)算穩(wěn)定的較小正常數(shù).

聯(lián)立公式(5)、(6)、(8)、(9)可以推出：

(10)

(11)

(12)

單頻矩陣I對應(yīng)波場記錄中聚焦點(diǎn)，反射波能量聚焦到該點(diǎn)周圍.在共炮點(diǎn)道集域，聚焦點(diǎn)對應(yīng)各炮的炮點(diǎn)；在共檢波點(diǎn)道集域，聚焦點(diǎn)為相應(yīng)的檢波點(diǎn).AP′為我們從多次波中還原得到的所有波場記錄；其中一階多次波轉(zhuǎn)換為反射波，二階多次波轉(zhuǎn)換為一階多次波，依次類推.A近似為震源特性矩陣的逆，因此聚焦變換具有去子波效應(yīng)，多次波轉(zhuǎn)換的新波場記錄比原波場記錄有更高的縱向分辨率.

2.2 聚焦變換的改進(jìn)

2.2.1 局部聚焦變換

在檢波點(diǎn)隨炮點(diǎn)移動(dòng)的觀測系統(tǒng)中，若用一個(gè)矩陣存儲(chǔ)全波場記錄，存儲(chǔ)量、計(jì)算量會(huì)相常的大.為了增強(qiáng)聚焦變換的適用性，我們提出局部聚焦變換：在逐個(gè)共炮點(diǎn)或共檢波點(diǎn)道集記錄上做聚焦變換.在炮點(diǎn)與檢波點(diǎn)重合的情況下，局部聚焦變換表示為：

(13)

相應(yīng)的：

(14)

2.2.2 引入加權(quán)矩陣

用全波場記錄矩陣的列存儲(chǔ)共檢波點(diǎn)道集記錄，行存儲(chǔ)共炮點(diǎn)道集記錄.引入對角矩陣Λ做加權(quán)矩陣，表示炮記錄的缺失：

(15)

Λk,l=1,k=l,炮存在；

Λk,l=0,k=l,炮缺失；

Λk,l=0,k≠l.

(16)

當(dāng)加權(quán)矩陣Λ的對角元素周期性的為1其余為0時(shí)，表示炮點(diǎn)比檢波點(diǎn)稀疏的觀測系統(tǒng).聯(lián)立公式(8)、(10)、(15)容易推導(dǎo)出：

Q=(ΛΔP)H[(ΛΔP)(ΛΔP)H+εI]-1(ΛP)

≈I+AP′

(17)

相應(yīng)的：

AP′≈(ΛΔP)H[(ΛΔP)(ΛΔP)H+εI]-1(ΛM),

(18)

在檢波點(diǎn)隨炮點(diǎn)移動(dòng)的觀測系統(tǒng)上：

Qj=(ΛΔPj)H[(ΛΔPj)(ΛΔPj)H+εI]-1(ΛPj),

(19)

(20)

2.3 一階多次波成像

(21)

在時(shí)間域?qū)⑥D(zhuǎn)換的反射波成像.將轉(zhuǎn)換的反射波成像，實(shí)質(zhì)上是對一階多次波的反射位置成像.本文采用對速度模型的橫向變化及陡傾角適應(yīng)性強(qiáng)的逆時(shí)偏移成像(互相關(guān)成像條件)，其基本原理為：

(22)

Ws(x,t)為時(shí)間正向延拓到地下的震源波場，Wg(x,t)為時(shí)間逆向延拓到地下的檢波點(diǎn)波場.

3 模型算例

3.1 Pluto1.5模型的完整數(shù)據(jù)

Pluto1.5聲波速度模型如圖4，共1387炮，炮間距和道間距為22.86 m，檢波點(diǎn)隨炮點(diǎn)移動(dòng).以炮點(diǎn)位于7985 m的第350炮為例，在共炮點(diǎn)道集域?qū)Ρ热鐖D5a、圖5b所示的兩種方案.如圖6a為包含反射波與多次波的原始炮記錄，圖6b、圖6c分別為通過SRME方法得到的炮記錄的反射波和多次波.

按照圖5a中的方案一，對如圖6a中包含反射波與多次波的原始炮記錄做局部聚焦變換(公式(13))得到新的波場記錄(圖7a).反射波聚焦到炮點(diǎn)周圍，如圖7a中水平箭頭所指位置；一階多次波轉(zhuǎn)換為反射波，二階多次波轉(zhuǎn)換為一階多次波，依次類推.減去炮點(diǎn)周圍的能量即可得到多次波轉(zhuǎn)換的新波場記錄(圖7b).按照圖5b中的方案二，直接對如圖6c中的多次波做局部聚焦變換(公式(14))得到多次波轉(zhuǎn)換的新波場記錄(圖8a).相比原始波場記錄，多次波轉(zhuǎn)換的新波場記錄(圖7b、圖8a)有更好的縱向分辨率，體現(xiàn)了聚焦變換的去子波效應(yīng).

由于SRME方法分離反射波與多次波過程中多次波能量的泄露(Berkhout and Verschuur, 2006)，方案一的優(yōu)勢在于能利用完全無損失的多次波信息.然而，當(dāng)采用方案一時(shí)，原始記錄的反射波能量并非最佳聚焦，部分能量泄露到多次波轉(zhuǎn)換的剖面中(如圖7b中傾斜箭頭所指位置).顯然，圖8a中的波場記錄比圖7b中波場記錄有更高的信噪比.這一現(xiàn)象或許與采樣數(shù)據(jù)的空間間隔不夠小以及引入的穩(wěn)定因子ε都使得聚焦變換不再完全準(zhǔn)確有關(guān).在方案一中經(jīng)常需要人工切除聚焦點(diǎn)及其周圍的能量，而這項(xiàng)工作非常繁雜.

圖4 Pluto1.5聲波速度模型Fig.4 Acoustic wave velocity of pluto1.5 model

圖5 局部聚焦變換(a)方案一和(b)方案二，(c)減去高階多次波轉(zhuǎn)換的多次波Fig.5 (a) First scheme and (b) second scheme for local focal transformation; (c) eliminate multiples transformed from high-order multiples

將多次波轉(zhuǎn)換的波場記錄(圖8a)中的高階多次波轉(zhuǎn)換的多次波減去(如圖5c)得到一階多次波轉(zhuǎn)換的反射波；對比圖8a與圖8b，圖8a中箭頭所指位置處的多次波較好的被消除.

3.2 Pluto1.5模型的抽稀數(shù)據(jù)

將原始炮集抽?。罕Ａ襞谔?hào)為奇數(shù)的原始記錄，炮號(hào)為偶數(shù)的原始記錄充零，此時(shí)炮間距是道間距的2倍.本例局部聚焦變換的工作流程采用如圖5b所示的方案二，直接對多次波做局部聚焦變換.與算例1不同之處在于：算例2首先在共炮點(diǎn)道集域分離反射波與多次波記錄；然后將反射波及多次波記錄分別抽到共檢波點(diǎn)道集域，炮記錄的缺失在共檢波點(diǎn)道集表現(xiàn)為數(shù)值為0的一個(gè)空道；多次波的局部聚焦變換在共檢波點(diǎn)道集域內(nèi)進(jìn)行；圖5c所示的消減高階多次波轉(zhuǎn)換的多次波在共檢波點(diǎn)道集內(nèi)進(jìn)行；最后將一階多次波轉(zhuǎn)換的反射波抽回到共炮域.共檢波點(diǎn)道集域的局部聚焦變換能夠從多次波中提取原始波場記錄中缺失的炮記錄信息.

利用SRME分離抽稀數(shù)據(jù)的反射波及多次波記錄，當(dāng)炮間距過大時(shí)應(yīng)當(dāng)采用其他方法分離反射波與多次波(Liu et al., 2009, 2010; 薛亞茹等, 2012).將反射波及多次波記錄抽到共檢波點(diǎn)道集域，按檢波器位置分別有1746個(gè)共檢波點(diǎn)道集.

以位于15324 m處的第850個(gè)共檢波點(diǎn)道集為例.為了方便與下文中的圖片對比，抽出第850個(gè)共檢波點(diǎn)道集的原始記錄(圖9a)在本文展示，圖9b為該共檢波點(diǎn)道集的反射波記錄.炮間距為道間距的2倍， 表現(xiàn)為每隔一道有一個(gè)數(shù)值為0的空道，反射波記錄的FK譜中出現(xiàn)空間假頻(圖9c).

圖6 原始炮集中第350炮的(a)原始記錄; (b)反射波估計(jì); (c)預(yù)測得到的多次波Fig.6 The 350th shot′s in original shot gather (a) Original wavefield record; b) Estimated primaries; (c) Predicted multiples

圖7 (a)原始記錄轉(zhuǎn)換的波場記錄; (b)切去圖(a)中聚焦點(diǎn)周圍的能量，得到多次波轉(zhuǎn)換的波場記錄Fig.7 (a)Wavefield record transformed from original record; (b) get wavefiled record transformed from multiples after muting energies around the focal point

圖8 (a)多次波直接轉(zhuǎn)換的波場記錄; (b)消減圖(a)中的多次波,得到一階多次波轉(zhuǎn)換的反射波Fig.8 (a) Wavefiled record transformed directly from multiples; (b) Primaries transformed from first-order multiples after subtracting multiples in figure(a)

圖9 抽稀原始炮集, 即將偶數(shù)炮充零并將奇數(shù)炮保留.然后在共炮域利用SRME分離反射波和多次波.將數(shù)據(jù)從共炮域整理到共檢波點(diǎn)道集域.圖中展示了抽稀數(shù)據(jù)第850個(gè)共檢波點(diǎn)道集的(a)原始記錄和(b)反射波記錄；(c)圖(b)中反射波的FK譜.Fig.9 Extract the original shot gather to be sparseness, i.e., fill the even number shots with zero and retain the odd number shots. And then, we separate primaries and multiples with SRME in the common shot domain. Rearrange data from common shot domain into common receiver domain. For the sparse data, the 850 th common receiver gather′s (a) Primaries and (b) Multiples are displayed in this figure; (c) Primaries′ FK spectrum in figure(b).

圖10 在共檢波點(diǎn)道集域，對多次波直接做局部聚焦變換(a) 多次波轉(zhuǎn)換的波場記錄，原始炮集中缺失的炮記錄信息得到還原; (b) 減去圖(a)中的多次波,得到一階多次波轉(zhuǎn)換的反射波記錄; (c) 圖(b)中反射波的FK譜.Fig.10 In the common receiver domain, do the local focal transformation directly with multiples(a) Wavefield record transformed from multiples, the information of missing shot records in original gather is retrieved; (b) Primaries transformed from first-order multiples got after subtracting multiples in figure(a); (c) Primaries′ FK spectrum in figure(b).

圖11 將一階多次波在共檢波點(diǎn)道集域轉(zhuǎn)換的反射波抽回到共炮域(a)抽稀炮集的反射波記錄；(b)一階多次波轉(zhuǎn)換的反射波記錄,原始炮集中缺失的炮記錄信息得到還原.圖中從左到右依次為第350、351、352、353炮.Fig.11 Rearrange primaries transformed from first order multiples in common receiver gather domain back into common shot domain(a)Primaries of shot gather extracted to be sparseness; (b)Primaries transformed from first order multiples, the information of missing shot record in original gather is retrieved. The 350th、351th、352th、353th shot are displayed from left to right.

圖12 灰白色與黑色曲線分別為圖11a、圖11b中第351炮反射波的平均振幅譜.兩條曲線被調(diào)整為有相同的最大值Fig.12 The gray and black curves are respectively average amplitude spectrums of 351th shots′ primaries in Fig.11a and Fig.11b. The maximum value of two curves are scaled to be identity

圖13 抽稀炮集中反射波的(如圖11a)逆時(shí)偏移剖面Fig.13 Reverse time migration section of primaries in shot gather extracted to be sparseness (as shown in Fig.11a)

圖14 一階多次波轉(zhuǎn)換的反射波(如圖11b)的逆時(shí)偏移剖面其有更均衡的能量分布、相對較高的縱向分辨率和局部更高的信噪比.Fig.14 Reverse time migration section of primaries transformed from first order multiples (as shown in Fig.11b), which has more balance energy distribution, relatively higher vertical resolution and locally has higher signal-to-noise ratio

直接對第850個(gè)共檢波點(diǎn)道集的多次波ΛMr850做局部聚焦變換(公式(20))，首先要將該道集覆蓋范圍內(nèi)的所有共檢波點(diǎn)道集的反射波估計(jì)按位置逐列存儲(chǔ)到矩陣ΛΔPr850中.從多次波中還原得到炮間距與道間距相同的新波場記錄(圖10a).將新波場記錄中高階多次波轉(zhuǎn)換的多次波減去，得到一階多次波轉(zhuǎn)換的反射波(圖10b).轉(zhuǎn)換反射波的FK譜如圖10c所示，其中沒有空間假頻.

原始抽稀炮集的反射波記錄如圖11a所示，將轉(zhuǎn)換的反射波記錄從共道集域(如圖10b)抽回到共炮域如圖11b所示.對比兩圖可知：從多次波還原得到了炮點(diǎn)采樣頻率與檢波點(diǎn)采樣頻率相同的新波場記錄，提供比原始記錄更豐富的疊前偏移信息.從圖12的振幅譜分析中可以看出，一階多次波轉(zhuǎn)換的反射波有比原始反射波更寬的頻帶范圍；由于聚焦變換的去子波效應(yīng)，一階多次波轉(zhuǎn)換的反射波比原始反射波有更高的縱向分辨率.

原始反射波、轉(zhuǎn)換反射波的逆時(shí)偏移結(jié)果分別如圖13、圖14所示，轉(zhuǎn)換反射波的逆時(shí)偏移成像即為一階多次波的聚焦變換逆時(shí)偏移成像.逆時(shí)偏移的成像條件是互相關(guān)成像條件；原始反射波逆時(shí)偏移的震源子波為15Hz雷克子波；為了與如圖12所示的振幅譜分析吻合，轉(zhuǎn)換反射波逆時(shí)偏移的震源子波采用20Hz雷克子波.對比圖13和圖14可知，一階多次波聚焦變換成像的優(yōu)勢體現(xiàn)在： 1)縱向分辨率更高，尤其在剖面的淺層，剖面能量整體分布更為均衡；2)在局部有更高的信噪比(如矩形虛線所圈位置)及更好的偏移噪聲壓制(如箭頭所指位置以及橢圓虛線所圈位置).

一階多次波成像剖面往往與傳統(tǒng)反射波成像剖面有近似的有效信息和不同的構(gòu)造假象，并在局部表現(xiàn)出更好的成像效果.因此： 1)一階多次波成像剖面可以成為反射波成像剖面解釋工作中的有益參考； 2)尋找使反射波成像剖面與一階多次波成像剖面最佳匹配的局部濾波器，壓制反射波成像剖面中的噪聲.此外，一階多次波成像提供了比傳統(tǒng)反射波成像更豐富的小角度信息，可以為速度分析作出貢獻(xiàn).

4 結(jié)論

本文敘述了聚焦變換的基本原理，并提出局部聚焦變換；理論上證明炮點(diǎn)比檢波點(diǎn)稀疏時(shí)，共檢波點(diǎn)道集域的局部聚焦變換可以將多次波準(zhǔn)確轉(zhuǎn)換成炮點(diǎn)與檢波點(diǎn)有相同采樣頻率的新波場記錄；數(shù)值算例很好地驗(yàn)證了上述理論.由于聚焦變換的去子波效應(yīng)，在兩個(gè)數(shù)值實(shí)驗(yàn)中，我們都可以看到多次波轉(zhuǎn)換的新波場記錄比原始波場記錄有更高的縱向分辨率.在第二個(gè)數(shù)值算例中，相比原始反射波成像，一階多次波聚焦變換成像的能量更為均衡、在局部表現(xiàn)出更好的偏移噪聲壓制和更高的信噪比、以及較高的縱向分辨率.通過理論分析和數(shù)值試驗(yàn)可以看出，多次波完全可以做為有效信號(hào)，提供更高的下地表覆蓋次數(shù)，對成像做出貢獻(xiàn).

致謝 感謝匿名審稿專家提出的寶貴意見.作者感謝國家自然科學(xué)基金項(xiàng)目、國家油氣重大專項(xiàng)以及中國科學(xué)院戰(zhàn)略性先導(dǎo)科技專項(xiàng)(B類)資助的聯(lián)合資助.

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

Berkhout A J. 1993. Migration of multiple reflections. ∥63th Ann. Mtg. and Intl. Expo., Soc. Expl. Geophys.. Expanded Abstracts, 1022-1025.Berkhout A J, Verschuur D J. 1994. Multiple technology: Part 2, migration of multiple reflections. ∥64th Ann. Internat. Mtg., Soc. Expl. Geophys.. Expanded Abstracts, 1497-1500. Berkhout A J, Verschuur D J. 1997. Estimation of multiple scattering by iterative inversion, Part I: Theoretical considerations.Geophysics, 62(5): 1586-1595, doi: 10.1190/1.1444261. Berkhout A J, Verschuur D J. 2003. Transformation of multiples into primary Reflections. ∥73rd Ann. Internat. Mtg., Soc. Expl. Geophys.. Expanded Abstracts, 1925-1928. Berkhout A J, Verschuur D J. 2004. Imaging multiple reflections, the concept. 74th Ann. Internat. Mtg., Soc. Expl. Geophys.. Expanded Abstracts, 1273-1276.

Berkhout A J, Verschuur D J. 2006. Imaging of multiple reflections.Geophysics, 71(4): SI209-SI220, doi: 10.1190/1.2215359.

Claerbout J F. 1968. Synthesis of a layered medium from its acoustic transmission response.Geophysics, 33(2): 264-269, doi: 10.1190/1.1439927.

Guitton A. 2002. Shot-profile migration of multiple reflections. ∥72nd Ann. Internat. Mtg., Soc. Expl. Geophys.. Expanded Abstracts, 1296-1299.

He R Q, Hornby B, Schuster G T. 2007. 3D wave-equation interferometric migration of VSP free-surface multiples.Geophysics, 72(5): S195-S203, doi: 10.1190/1.2743375. Hu H, Liu Y K, Chang X, et al. 2013. Analysis and application of boundary treatment for the computation of reverse-time migration.ChineseJ.Geophys. (in Chinese), 56(6): 2033-2042, doi: 10.6038/cjg20130624. Jiang Z Y, Sheng J M, Yu J H, et al. 2007. Migration methods for imaging different-order multiples.Geophys.Prospect., 55(1): 1-19, doi: 10.1190/1.2370257.

Li P, Liu Y K, Chang X, et al. 2007. Application of the equipoise pseudomultichannel matching filter in multiple elimination using wave equation method.ChineseJ.Geophys. (in Chinese), 50(6): 1844-1853.

Liu Y K, Jin D G, Chang X, et al. 2009. Multiple subtraction using statistically estimated inverse wavelets. ∥79th Ann. Internat. Mtg., Soc. Expl. Geophys.. Expanded Abstracts, 3098-3102.

Liu Y K, Jin D G, Chang X, et al. 2010. Multiple subtraction using statistically estimated inverse wavelets.Geophysics, 75(6): WB247-WB254, doi: 10.1190/1.3494082.

Liu Y K, Chang X. 2011a. Reverse time migration of multiples. ∥81th Ann. Internat. Mtg., Soc. Expl. Geophys.. Expanded Abstracts, 3326-3331.

Liu Y K, Chang X, Jin D G, et al. 2011b. Reverse time migration of multiples for subsalt imaging.Geophysics, 76(5): WB209-WB216, doi: 10.1190/1.3627888.

Lu S P, Whitmore N D, Valenciano A A et al. 2011. Imaging of primaries and multiples with 3D SEAM Synthetic. ∥81th Ann. Internat. Mtg., Soc. Expl. Geophys.. Expanded Abstracts, 3217-3221.

Manuel C D, Uren N F. 2001. Migration Velocity analysis using seismic multiples. ∥ASEG 15th Geophysical Conference and Exhibition, 1-4.

Muijs R, Robertsson J O A, Holliger K. 2007. Prestack depth migration of primary and surface-related multiple reflections: Part II—Identification and removal of residual multiples.Geophysics, 72(2): S71-S76, doi: 10.1190/1.2424544.

Schuster G T, Yu J, Sheng J, et al. 2004. Interferometric/daylight seismic imaging.Geophys.J.Int, 157(2): 838-852, doi: 10.1111/j.1365-246X.2004.02251.x.

Schuster G. 2008. Seismic Interferometry. Cambridge University Press.

Sheng J. 2001. Migrating multiples and primaries in CDP data by crosscorrelogram Migration. ∥71th Ann. Internat. Mtg., Soc. Expl. Geophys.. Expanded Abstracts, 1297-1300.

van Groenestijn G V, Verschuur D J. 2006. Reconstruction of missing data from multiples using the focal transform. ∥76th Ann. Internat. Mtg., Soc. Expl. Geophys.. Expanded Abstracts, 2737-2741.Vasconcelos I, Snieder R, Hornby B. 2008. Imaging internal multiples from subsalt VSP data—Examples of target-oriented interferometry.Geophysics, 73(4): S157-S168, doi: 10.1190/1.2944168.

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

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, doi: 10.1190/1.1444262.

Verschuur D J, Berkhout A J. 2005. Transforming multiples into Primaries: experience with field data. ∥75th Ann. Internat. Mtg., Soc. Expl. Geophys.. Expanded Abstracts, 2103-2106.

Wang Y B, Chang X, Hu H. 2014. Simultaneous reverse time migration of primaries and free-surface related multiples without multiple prediction.Geophysics, 79(1): S1-S9, doi: 10.1190/geo2012-0450.1.

Wapenaar K, Fokkema J. 2006. Green′s function representations for seismic interferometry.Geophysics, 71(4): SI33-SI46, doi: 10.1190/1.2213955.

Xue Y R, Chen X H, Ma J T. 2012. Multiples attenuation based on directional orthogonal polynomial transform.ChineseJ.Geophys. (in Chinese), 55(10): 3450-3458, doi: 10.6038/j.issn.0001-5733.2012.10.028.

附中文參考文獻(xiàn)

胡昊, 劉伊克, 常旭等. 2013. 逆時(shí)偏移計(jì)算中的邊界處理分析及應(yīng)用. 地球物理學(xué)報(bào), 56(6): 2033-2042, doi: 10.6038/cjg20130624.

李鵬, 劉伊克, 常旭等. 2007. 均衡逆多道匹配濾波法在波動(dòng)方程法壓制多次波中的應(yīng)用. 地球物理學(xué)報(bào), 50(6): 1844-1853.

薛亞茹, 陳小宏, 馬繼濤. 2012. 多方向正交多項(xiàng)式變換壓制多次波. 地球物理學(xué)報(bào), 55(10): 3450-3458, doi: 10.6038/j.issn.0001-5733.2012.10.028.

(本文編輯 汪海英)

Focal transformation imaging of first-order multiples

LIU Xue-Jian1,2, LIU Yi-Ke1, HU Hao1,2, XIE Song-Lei3

1KeyLaboratoryofEngineeringGeomechanics,InstituteofGeologyandGeophysics,ChineseAcademyofSciences,Beijing100029,China2UniversityofChineseAcademyofSciences,Beijing100049,China3SouthChinaSeaInstituteofOceanology,ChineseAcademyofSciences,Guangzhou510301,China

Surface-related multiples penetrate into the subsurface several times and contain abundant reflection information of small angles. Compared with primaries, they sometimes can provide higher fold and better illumination for subsurface. Instead of trashing multiples as noises during the seismic data processing, nowadays, a lot of methods have been proposed to image multiples. Conventional migration methods have been modified for directly imaging multiples, e.g., Kirchhoff migration of multiples, wave equation migration of multiples, or reverse-time migration (RTM) of multiples. Alternatively, the linear two-step procedures can be utilized for imaging multiples. Primaries can be extracted from multiples based on seismic interferometry or focal transformation, and then imaged by utilizing conventional migration methods. However, most methods for imaging multiples would generate so many artifacts due to the undesired interactions between seismic-events (i.e. primaries and different-order multiples), which are hard to be attenuated and seriously pollute the true-image. Moreover, the least-squares based focal transformation can transform first-order multiples into primaries, second-order multiples into first-order multiples, etc., with few noises generated. Multiples focal-transformed from higher-order multiples can be eliminated by surface-related multiples elimination (SRME). And then, primaries focal-transformed from first-order multiples can be imaged by any conventional migration methods, and we call the procedures focal transformation imaging of first-order multiples.The focal transformation is proposed by Berkhout and Verschuur, which is developed from SRME. Primaries are once subsurface response of sources, surface-related multiples include different-orders, and the raw data contain primaries and surface-related multiples. Focal transformation is manipulated with several matrices that contain the full data of mono-frequency, one column stores one shot gather or one common-receiver gather. Primaries need to be estimated prior to focal transformation, and the inverse of primaries matrix is used as focal operator. The full raw data will be transformed into the new wavefield records. In fact, the focal transformation is implemented with a stable form of least-squares sense. Primaries are focused around one point in the profile of zero-time and zero-offset, the so called focal-point; first-order multiples are transformed into primaries; second-order multiples are transformed into first-order multiples; etc. Focal transformation is a kind of least-squares transformation, which nearly doesn't generate noises. So, after muting the energies around focal point, SRME can be utilized to eliminate the multiples transformed from higher-order multiples, and transformed primaries are obtained. Alternatively, in order to utilize the information of multiples, we can directly implement the focal transformation of multiples that usually have been separated from primaries during the regular seismic data processing, and the muting work can be avoided. On the other hand, the focal transformation has the effect of deconvolution by utilizing the inverse of real source-signature matrix, so the transformed primaries have higher vertical resolution than the acquired primaries. However, when the receive array moves with the source position, the matrices will occupy large memory that are mostly wasted by off diagonal 0 elements. In this article, we put forward two improvements for focal transformation: (1) develop the local focal transformation for reducing storage and computation; (2) bring in the weighted matrix, and demonstrating that local focal transformation in common-receiver domain can transform multiples into new wavefield records retrieving the missing shot-gathers of acquired data. The local transformation is implemented specially for one shot gather or one common-receiver gather, so the local focal transformation has better adaptation to the acquired data whose traces move with corresponding source position. We introduce the diagonal matrix into the focal transformation. When source sampling rate is sparser than receiver sampling rate of acquired data and one column of the matrices stores one common-receiver gather, the diagonal elements of the weighted matrix are periodically 1 spaced by 0 diagonal elements. The focal transformation transforms multiples into new wavefield records where source sampling is the same with receiver sampling, when the source sampling is sparser than receiver sampling in the acquired data. The local focal transformation in common-receiver domain also can transform multiples into new wavefild records with denser source sampling, but the common-receiver gathers must be extracted from shot-gathers in advance.In the first numerical test, two workflows of local focal transformation of row data and direct local focal transformation of predicted multiples are compared. They are both implemented in common-shot domain. We verify the second work flow will generate wavefield records with higher signal to noise ratio, and the laborious task for muting the energies around focal point is avoided. With the first workflow, a few parts of primaries are not focused around the focal point and leaked into the profile transformed from multiples. With the second workflow, the new wavefiled records directly transformed from multiples nearly do not have noises, and multiples transformed from higher multiples can be successfully eliminated by SRME. Obviously, the wavefield records transformed from multiples have higher vertical resolution than the acquired data. The second numerical test demonstratesthat when source sampling is relatively sparse, the local focal transformation in common-receiver domain can effectively transform first-order multiples into primaries; the transformed primaries include the missing shot-gathers of the acquired data; the imaging result of transformed primaries is more balanced, locally has higher signal-to-noise ratio, and shows slightly higher vertical resolution. The separated primaries and multiples are both rearranged into common-receiver gathers in advance; the zero traces in common-receiver gathers show the missing shot records; that every other trace in common-receiver gathers are zero represents the sparse source sampling, and alias artifacts can be clearly seen in the FK domain. The direct local transformation of multiples on common-receiver gathers can retrieve the missing shot records from multiples, and the alias artifacts in the FK domain are avoided. In the common-receiver domain, multiples transformed from higher-order multiples are eliminated by SRME, and the transformed primaries are rearranged back into the common-shot domain. The transformed primaries have wider amplitude spectrum than acquired primaries, which demonstrate focal transformation has the effect of deconvolution. The transformed primaries and acquired primaries are both migrated by RTM.The focal transformation can transform first-order multiples into primaries with few noises, and it has the effect of deconvolution. The proposed local focal transformation is implemented specially for one shot or one common-receiver gather, so the local focal transformation can save computation and storage. The local focal transformation in common-receiver domain can retrieve the missing shot records of acquired data from multiples. When source sampling of acquired data is sparse, the RTM image of primaries transformed from first-order multiples is more balanced, locally has higher signal-to-noise ratio, and shows slightly higher vertical resolution.

Focal transformation; Imaging of multiples; Multiples elimination

10.6038/cjg20150614.

國家自然科學(xué)基金項(xiàng)目(41374138)，國家油氣重大專項(xiàng)(2011ZX05008-006)和中國科學(xué)院戰(zhàn)略性先導(dǎo)科技專項(xiàng)(B類)(XDB01020300)聯(lián)合資助.

劉學(xué)建，男，1987年生，在讀博士研究生，從事表面多次波消除、多次波成像方法以及最小二乘逆時(shí)偏移等方面的研究. E-mail:liuxuejian10@mails.ucas.ac.cn

10.6038/cjg20150614

P631

2014-05-08，2015-05-08收修定稿

劉學(xué)建，劉伊克，胡昊等. 2015. 一階多次波聚焦變換成像.地球物理學(xué)報(bào)，58(6):1985-1997,

Liu X J, Liu Y K, Hu H, et al. 2015. Focal transformation imaging of first-order multiples.ChineseJ.Geophys. (in Chinese),58(6):1985-1997,doi:10.6038/cjg20150614.