何兵壽,武雪嶠,高琨鵬

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

各向異性理論的研究始于19世紀(jì)50年代[1],但早期的研究工作主要由數(shù)學(xué)物理學(xué)家主導(dǎo),其應(yīng)用范圍也局限在晶體光學(xué)和晶體彈性方面。1911年,RUDZKI[2]強(qiáng)調(diào)了地震各向異性的重要性并第一次實(shí)現(xiàn)了各向異性彈性波的定量計(jì)算。隨后該領(lǐng)域的研究進(jìn)入停滯期,直到20世紀(jì)50年代,人們才開始研究地震各向異性對(duì)勘探地震的意義[3-9],在這個(gè)時(shí)期取得的主要成果包括:各向異性介質(zhì)中地震波垂向和橫向傳播速度不一致的原因[4]、薄互層組合各向異性的產(chǎn)生機(jī)理[5]、方位各向異性的提出[6]、橫波分裂現(xiàn)象的發(fā)現(xiàn)[7]與觀測(cè)[8]及各向異性本構(gòu)方程的建立[9]等。

從20世紀(jì)80年代起,地震各向異性的研究進(jìn)一步細(xì)化,開始針對(duì)不同類型的各向異性介質(zhì)展開研究,并提出了許多現(xiàn)在仍然廣泛應(yīng)用的理論。其中,最具代表性的是THOMSEN[10]提出的基于弱各向異性假設(shè)利用5個(gè)新參數(shù)描述橫向各向同性(TI)介質(zhì)中縱、橫波的各向異性特性的理論,這5個(gè)參數(shù)分別為:vP0,vS0,ε,γ和δ,前兩個(gè)參數(shù)代表P波和S波在對(duì)稱軸方向的速度,后3個(gè)參數(shù)統(tǒng)稱為THOMSEN參數(shù)(無(wú)量綱),分別代表縱、橫波的各向異性強(qiáng)度和縱波的變異系數(shù)。由于這5個(gè)參數(shù)具有更直觀的物理意義且便于應(yīng)用,因此,該表征方式一經(jīng)提出就受到了高度關(guān)注。

THOMSEN的工作是各向異性理論與技術(shù)的有力推動(dòng),幾十年來(lái),業(yè)界在各向異性領(lǐng)域開展了大量的研究工作并取得了諸多成果。主要包括各向異性地震波場(chǎng)模擬方法、橫波分裂和縱橫波慢度的研究、快慢橫波的提取方法、裂縫或各向異性參數(shù)反演方法、TI介質(zhì)中準(zhǔn)P波方程的推導(dǎo)以及準(zhǔn)P波或彈性波各向異性疊前偏移方法等。其中:疊前偏移技術(shù)是各向異性地震勘探領(lǐng)域的關(guān)鍵技術(shù)之一,理論上說,地震波的高精度疊前偏移不僅能實(shí)現(xiàn)地下構(gòu)造的準(zhǔn)確成像,而且能為疊前反演提供準(zhǔn)確的共成像點(diǎn)道集。因此,業(yè)界對(duì)此項(xiàng)技術(shù)的研究尤為重視。

TI介質(zhì)模型是當(dāng)前最具代表性和應(yīng)用潛力的各向異性理論,本文中的分析與論述僅限于該類介質(zhì)。TI介質(zhì)中縱波的偏振方向一般不與其傳播方向平行,橫波的偏振方向一般不與其傳播方向垂直,因此它們都不是嚴(yán)格理論意義上的縱波與橫波,但又具有與理論縱、橫波相似的物理性質(zhì),因此稱為準(zhǔn)縱(quasi-P,qP)波和準(zhǔn)橫(quasi-S,qS)波。由于TI介質(zhì)中qP波和qS波是耦合在一起的,不存在理論嚴(yán)謹(jǐn)?shù)膓P波方程和qS波方程,因此,只有采用彈性波方程對(duì)三分量數(shù)據(jù)進(jìn)行地震各向異性逆時(shí)偏移才具有理論的嚴(yán)謹(jǐn)性。目前地震勘探工作仍以單分量縱波勘探為主,許多探區(qū)缺乏三分量地震數(shù)據(jù),工業(yè)界迫切需要一種能夠解決單分量qP波成像問題的逆時(shí)偏移方法,該需求推動(dòng)了TI介質(zhì)中qP波方程逆時(shí)偏移技術(shù)的產(chǎn)生和發(fā)展,具體思路為:在一定的假設(shè)條件下,對(duì)TI介質(zhì)中的彈性波方程進(jìn)行簡(jiǎn)化與推導(dǎo),得到TI介質(zhì)中的近似qP波方程,并由此結(jié)合逆時(shí)偏移的基本思路實(shí)現(xiàn)TI介質(zhì)中的qP波逆時(shí)偏移。

本文主要分析當(dāng)前基于TI介質(zhì)理論的qP波逆時(shí)偏移技術(shù),首先給出TI介質(zhì)中qP波方程逆時(shí)偏移的基本流程,然后分析該流程中各個(gè)環(huán)節(jié)的研究現(xiàn)狀和存在的主要問題,并給出下一步研究建議。

1 TI介質(zhì)中qP波方程逆時(shí)偏移的實(shí)現(xiàn)思路與基本流程

TI介質(zhì)中qP波方程逆時(shí)偏移一般通過數(shù)值求解各成像點(diǎn)的炮點(diǎn)波場(chǎng)和接收點(diǎn)波場(chǎng),并利用時(shí)間一致性原理對(duì)炮、檢波場(chǎng)進(jìn)行互相關(guān)來(lái)實(shí)現(xiàn),典型的計(jì)算流程如圖1所示。

圖1中各向異性參數(shù)建模不屬于本文討論范圍,炮記錄為采集得到的單分量單炮記錄,其它黑框部分是中間結(jié)果,紅框部分組成了整個(gè)TI介質(zhì)中qP波方程逆時(shí)偏移技術(shù)體系,本文只論述紅框部分的各個(gè)環(huán)節(jié)。qP波方程是逆時(shí)偏移的理論基礎(chǔ),激發(fā)源設(shè)置方法、qP波方程的數(shù)值解法、吸收邊界條件和隨機(jī)邊界組成qP波方程正向和逆時(shí)延拓技術(shù);波場(chǎng)傳播方向優(yōu)選、偏移噪聲壓制、成像條件應(yīng)用和波場(chǎng)拉伸校正組成波場(chǎng)成像技術(shù)。本文對(duì)TI介質(zhì)中qP波逆時(shí)偏移研究現(xiàn)狀的分析也從這3個(gè)方面展開。

圖1 TI介質(zhì)中qP波方程逆時(shí)偏移的計(jì)算流程

2 TI介質(zhì)中qP波方程的研究現(xiàn)狀

TI介質(zhì)中的qP波方程一般由具有垂直對(duì)稱軸的橫向各向同性(transversely isotropic media with vertical symmetry axis,VTI)介質(zhì)中的qP波方程進(jìn)行坐標(biāo)旋轉(zhuǎn)得到。VTI介質(zhì)中qP波方程的推導(dǎo)主要有弱各向異性近似[10-12]和聲學(xué)近似[13-21]兩種思路。弱各向異性近似的研究源于THOMSEN[10],他假定3個(gè)THOMSEN參數(shù)均遠(yuǎn)小于1,然后對(duì)VTI介質(zhì)中的qP波相速度表達(dá)式進(jìn)行泰勒展開并只保留線性項(xiàng),得到弱各向異性近似下的qP波相速度,再利用相速度與頻率波數(shù)之間的關(guān)系公式得到弱各向異性近似頻散關(guān)系方程,將其變換到時(shí)間域即可得到VTI介質(zhì)的qP波方程。TSVANKIN等[11-12]得到了與之相同的結(jié)果。

聲學(xué)近似的研究始于ALKHALIFAH[13-15],他假定沿對(duì)稱軸方向的橫波速度為零,然后利用VTI介質(zhì)中的精確頻散關(guān)系得到聲學(xué)假設(shè)條件下的頻散關(guān)系式,最后經(jīng)過傅里葉反變換得到VTI介質(zhì)中的qP波方程,再對(duì)該方程進(jìn)行坐標(biāo)旋轉(zhuǎn)得到TI介質(zhì)中的qP波方程。ALKHALIFAH方程是一個(gè)同時(shí)包含時(shí)間和空間四階偏導(dǎo)數(shù)的偏微分方程,方程中還存在混合偏導(dǎo)數(shù)項(xiàng),不便于求解與工程應(yīng)用。為解決此問題,學(xué)者們通過引入不同的輔助函數(shù)對(duì)該方程進(jìn)行了降階處理[16-21],將其轉(zhuǎn)換為更容易求解的二階方程組,比較典型的方程主要有:ZHOU提出的方程[16]、DUVENECK等提出的方程[18]、FLETCHER等提出的方程[20]和KANG等[21]提出的方程,這些方程均可以較好地描述TI介質(zhì)中qP波的運(yùn)動(dòng)學(xué)特征。

沈銘成的研究[23]表明,對(duì)于各向異性介質(zhì)qP波,ALKHALIFAH聲學(xué)近似假設(shè)的精度明顯高于弱各向異性近似。因此目前業(yè)界傾向于采用聲學(xué)近似方程解決TI介質(zhì)的qP波偏移問題。但是聲學(xué)近似的qP波方程存在以下缺陷:①存在菱形偽橫波干擾(圖2中黑色箭頭所指處);②當(dāng)THOMSEN參數(shù)δ>ε時(shí)難以穩(wěn)定求解;③當(dāng)對(duì)稱軸傾角急劇變化時(shí),方程不能穩(wěn)定求解[22]。

圖2 由ZHOU提出的方程得到的TI介質(zhì)中的波場(chǎng)快照(ε=0.3,δ=0.1,θ=45°)

上述偽橫波干擾不是TI介質(zhì)的地震響應(yīng),而是方程本身的誤差解,該誤差會(huì)給逆時(shí)偏移帶來(lái)假象并產(chǎn)生嚴(yán)重的數(shù)值頻散問題(圖2中紅色箭頭所指處),在偏移過程中必須壓制。目前的壓制思路主要有:①依據(jù)偽橫波的產(chǎn)生機(jī)理,在數(shù)值求解過程中采用適當(dāng)?shù)姆椒ㄟM(jìn)行壓制[14,22-25];②推導(dǎo)不存在偽橫波干擾的純qP波方程[26-33]。第一種思路不改造方程本身,只在數(shù)值求解過程中采用適當(dāng)方法壓制偽橫波,主要有各向異性介質(zhì)參數(shù)匹配法[14,24]、波場(chǎng)濾波法[25]和輔助波場(chǎng)壓制方法[22]等。第二種思路側(cè)重從方程的角度解決問題,通過推導(dǎo)純qP波方程來(lái)消除偽橫波干擾,由于方程只存在qP波解,偽橫波干擾自動(dòng)消除。

VTI介質(zhì)中純qP波方程的推導(dǎo)思路總結(jié)為:首先用一階泰勒展開式近似表示VTI介質(zhì)精確相速度公式中的平方根項(xiàng),得到相速度近似公式,然后由此得到VTI介質(zhì)中qP波的頻散關(guān)系式,最后利用該頻散關(guān)系得到不同形式的純qP波方程。如對(duì)該頻散關(guān)系等式的兩邊同時(shí)乘以等式右邊的分母項(xiàng)并進(jìn)行時(shí)間和空間方向的傅里葉反變換,得到一個(gè)四階的偏微分方程[33],或者通過引入輔助變量將其轉(zhuǎn)化為一個(gè)二階偏微分方程組[26,34]。也可以對(duì)頻散關(guān)系式只進(jìn)行時(shí)間方向的傅里葉反變換得到時(shí)間波數(shù)域方程,然后用偽譜法求解。但這樣得到的方程分母項(xiàng)中包含了各向異性參數(shù),不利于偽譜法求解,HARLAN[35]將相速度公式的分母項(xiàng)中所包含的各向異性參數(shù)舍去,得到相速度的HARLAN近似式,并進(jìn)一步得到純qP波方程。由于HARLAN近似式存在較大誤差,郭成鋒等[36]對(duì)該近似式進(jìn)行了修正,在非橢圓項(xiàng)前添加了一個(gè)與ε有關(guān)的修正系數(shù),得到了3種改進(jìn)型HARLAN公式,并以近似精度最高的改進(jìn)式為基礎(chǔ)推導(dǎo)了TI介質(zhì)純qP波方程。此外、黃翼堅(jiān)等[32],杜啟振等[26]和ZHAN等[27]均利用不同的近似方法推導(dǎo)了TI介質(zhì)純qP波方程,這些方程都能在運(yùn)動(dòng)學(xué)上較好地?cái)M合TI介質(zhì)中qP波的傳播過程,且不存在偽橫波干擾,但許多方程中存在分?jǐn)?shù)形式的偽微分算子,增加了求解的難度與計(jì)算成本。

為克服上述局限,XU等[30]通過在各向同性介質(zhì)聲波方程中加入一個(gè)標(biāo)量算子實(shí)現(xiàn)對(duì)相速度的控制,并由此導(dǎo)出VTI介質(zhì)的純qP波方程,該方程的優(yōu)點(diǎn)是能在各向異性參數(shù)δ>ε情況下穩(wěn)定求解,缺點(diǎn)是描述的相速度只在波前面梯度的方向與VTI介質(zhì)相速度相符,在其它方向則存在明顯誤差,當(dāng)ε較大時(shí),這種相速度的誤差會(huì)導(dǎo)致波前能量分布的不均衡并產(chǎn)生噪聲。XU等[30]進(jìn)一步對(duì)該方程進(jìn)行了改正,他通過在聲波方程的空間微分項(xiàng)中加入各向異性參數(shù),使方程的相速度曲線更接近VTI介質(zhì)實(shí)際相速度,改進(jìn)后方程的相速度不僅在波前面的傳播方向上與VTI介質(zhì)的相速度相等,在其它方向上也更準(zhǔn)確,且不存在波前能量分布不均問題。ZHANG等[37]和張建敏[38]以該方程為基礎(chǔ),通過投影變換將該方程推廣至TTI介質(zhì)。ZHANG等[37]提出的TI介質(zhì)純qP波方程波前能量分布均勻(圖3),不存在偽微分算子,結(jié)構(gòu)簡(jiǎn)單,便于求解。張建敏[38]還通過引入輔助變量得到TI介質(zhì)中的一階純qP波方程,該方程能夠采用旋轉(zhuǎn)交錯(cuò)網(wǎng)格有限差分技術(shù)求解,且相同條件下,一階方程的計(jì)算精度高于二階方程。

圖3 由ZHANG等[37]提出的方程得到的TI介質(zhì)中的波場(chǎng)快照(ε=0.3,δ=0.1,θ=45°)

3 TI介質(zhì)中qP波方程的延拓方法

本節(jié)主要介紹qP波方程逆時(shí)偏移中的激發(fā)源設(shè)置方法、qP波方程的數(shù)值解法、吸收邊界條件和炮點(diǎn)波場(chǎng)重構(gòu)策略等方面的研究現(xiàn)狀。

3.1 激發(fā)源設(shè)置方法

逆時(shí)偏移成像一般通過炮點(diǎn)正向延拓波場(chǎng)和接收點(diǎn)逆時(shí)延拓波場(chǎng)的互相關(guān)來(lái)實(shí)現(xiàn),其中炮點(diǎn)波場(chǎng)的計(jì)算需要給定各向異性參數(shù)模型和震源子波,震源子波的類型、形態(tài)、主頻以及頻寬等決定炮點(diǎn)波場(chǎng)的諸多屬性,炮點(diǎn)波場(chǎng)的屬性又會(huì)影響逆時(shí)偏移的成像結(jié)果。不合理的震源子波會(huì)弱化炮點(diǎn)波場(chǎng)與接收點(diǎn)波場(chǎng)之間的相關(guān)性,降低成像準(zhǔn)確性和分辨率[39]。因此,設(shè)置合理的震源子波是獲取準(zhǔn)確逆時(shí)偏移成像結(jié)果的前提。但相對(duì)于qP波逆時(shí)偏移領(lǐng)域的其它技術(shù),激發(fā)源設(shè)置方面的研究相對(duì)滯后,特別是震源子波對(duì)逆時(shí)偏移的影響機(jī)理尚未明確,許多工作只得到了一些感性認(rèn)識(shí),未上升到理論高度。

目前,業(yè)界在該領(lǐng)域的研究主要集中在理論子波的特性分析方面[39],而實(shí)際地面記錄的炮點(diǎn)子波遠(yuǎn)非理論子波,將震源子波設(shè)置成理論子波不合適,而設(shè)定為野外近源道記錄的子波或激發(fā)震源的理論子波則會(huì)給逆時(shí)偏移帶來(lái)嚴(yán)重的頻散問題,并降低偏移剖面質(zhì)量。因此,劉學(xué)義[39]直接從實(shí)測(cè)資料中提取子波并將其設(shè)置為震源子波進(jìn)行逆時(shí)偏移,這種思路可以降低由于震源子波與記錄子波不一致時(shí)對(duì)偏移結(jié)果的影響,它比直接將震源子波設(shè)置為理論子波更具合理性。問題的關(guān)鍵是如何從炮記錄中提取出子波,目前大都采用反演領(lǐng)域或其它處理技術(shù)中的子波求取方法提取逆時(shí)偏移子波[40-41],還未見到針對(duì)qP波逆時(shí)偏移這一特定問題研發(fā)的子波提取技術(shù)。

3.2 qP波方程的數(shù)值解法

qP波方程的數(shù)值求解方法很多,每種qP波方程都有對(duì)應(yīng)的求解方法。歸納起來(lái),主要包括有限差分法、偽譜法以及有限差分與偽譜法的混合法3種,而有限差分法又包括交錯(cuò)網(wǎng)格法、旋轉(zhuǎn)網(wǎng)格法、旋轉(zhuǎn)交錯(cuò)網(wǎng)格法和緊致交錯(cuò)網(wǎng)格法等,幾乎每種各向同性領(lǐng)域的波動(dòng)方程數(shù)值求解技術(shù)都可以直接或稍加改動(dòng)后用于解決qP波方程的數(shù)值求解問題。關(guān)于qP波方程的數(shù)值求解問題已有多人加以綜述[42-43],本文不贅述。

3.3 吸收邊界條件

接收點(diǎn)波場(chǎng)逆時(shí)延拓過程中,需要采用邊界條件來(lái)壓制計(jì)算邊界的偽反射。目前常用的邊界條件有:①基于單程波方程的吸收邊界[44-45]。這類方法用擬微分算子分解波動(dòng)算子并只保留向外傳播部分,然后用Pade近似得到不同精度的邊界微分方程,其效果取決于邊界微分方程的精度,已有的低階單程波方程均只能在小入射角前提下取得良好近似,當(dāng)入射角較大時(shí)方程存在明顯誤差。提高單程波方程的階數(shù)可以減小這種誤差[46],但會(huì)增大算法復(fù)雜度和計(jì)算量,同時(shí)帶來(lái)穩(wěn)定性隱患。②基于波場(chǎng)衰減的吸收邊界[47-48]。這類方法通過對(duì)邊界附近的波場(chǎng)進(jìn)行衰減消除邊界反射,其效果取決于衰減函數(shù)的選擇,模型復(fù)雜時(shí)很難找到一個(gè)合理的衰減函數(shù)。③基于完全匹配層(Perfectly Matched Layer,PML)技術(shù)的吸收邊界條件[32-37,49-53]。此類方法的基本思路為,在計(jì)算區(qū)域周圍設(shè)置吸收衰減帶并沿坐標(biāo)軸方向分解波場(chǎng),在不同邊界區(qū)域按照不同方式對(duì)分解后的波場(chǎng)進(jìn)行衰減,達(dá)到消除邊界反射的目的。

PML邊界條件最早出現(xiàn)于電磁波模擬領(lǐng)域[49-50],隨后被引入地震波數(shù)值模擬與逆時(shí)偏移領(lǐng)域[32-37,51-53],之后許多研究工作拓展了PML的使用范圍[52-56],目前該邊界條件已成為TI介質(zhì)qP波逆時(shí)偏移領(lǐng)域應(yīng)用最為廣泛的邊界處理方法。張巖[52]將一個(gè)二階qP波方程變換到頻率域并進(jìn)行坐標(biāo)變換,在頻率域?qū)Σ▓?chǎng)進(jìn)行分裂處理,然后逆變換至?xí)r間域得到該方程的二階PML邊界條件并在模型數(shù)據(jù)的qP波逆時(shí)偏移中取得了良好效果。黃建平等[53]給出了TI介質(zhì)一階速度-應(yīng)力qP波方程的PML邊界條件,在交錯(cuò)網(wǎng)格中推導(dǎo)了該邊界條件的高階有限差分格式并用模型數(shù)據(jù)加以驗(yàn)證;ZHANG等[37]推導(dǎo)出了二階純qP波方程的PML吸收邊界及對(duì)應(yīng)的有限差分格式,此算法能夠很好地吸收TI介質(zhì)中入射到截?cái)噙吔绲耐庑胁?ZHANG等[37]還給出了TI介質(zhì)中一階純qP波方程的PML吸收邊界條件,該邊界條件同樣具有良好的外行波吸收性能。

3.4 炮點(diǎn)波場(chǎng)重構(gòu)策略

地震波逆時(shí)偏移常采用以計(jì)算換存儲(chǔ)的思路降低臨時(shí)文件存儲(chǔ)量和硬盤讀寫量,即進(jìn)行兩次炮點(diǎn)波場(chǎng)延拓,第一次延拓時(shí)存儲(chǔ)部分波場(chǎng)信息,然后利用這些信息對(duì)炮點(diǎn)波場(chǎng)進(jìn)行逆時(shí)重構(gòu)。與此同時(shí),接收點(diǎn)波場(chǎng)的逆時(shí)延拓也并行處理,這樣就能實(shí)現(xiàn)炮點(diǎn)和接收點(diǎn)波場(chǎng)的同步延拓,避免逆時(shí)偏移過程中的臨時(shí)文件存儲(chǔ)和讀寫工作,發(fā)揮并行計(jì)算機(jī)的計(jì)算優(yōu)勢(shì),提升偏移效率。

目前的炮點(diǎn)波場(chǎng)重構(gòu)方法主要有4種:①基于檢查點(diǎn)技術(shù)的重構(gòu)方法[54],在正傳過程中每隔一定時(shí)間間隔設(shè)置一個(gè)檢查點(diǎn)并儲(chǔ)存該點(diǎn)波場(chǎng),反傳過程中從檢查點(diǎn)出發(fā)通過反復(fù)遞推重構(gòu)任意時(shí)刻的波場(chǎng);②基于抽樣插值法的重構(gòu)方法[55],在滿足采樣定理的前提下對(duì)正傳波場(chǎng)進(jìn)行規(guī)則抽樣,并通過數(shù)據(jù)壓縮算法對(duì)抽樣波場(chǎng)進(jìn)行存儲(chǔ),采用適當(dāng)?shù)牟逯捣椒▽?duì)信號(hào)重采樣,插值重構(gòu)任意時(shí)刻的有效邊界層波場(chǎng),通過最后兩個(gè)時(shí)刻波場(chǎng)的反傳并加入插值得到的有效邊界層重構(gòu)炮點(diǎn)波場(chǎng);③基于有效邊界存儲(chǔ)的波場(chǎng)重構(gòu)方法[56-58],該方法存儲(chǔ)每次正傳波場(chǎng)的四周有效邊界層,并利用最后兩次時(shí)間采樣的波場(chǎng)做逆時(shí)傳播,將每個(gè)時(shí)間采樣點(diǎn)對(duì)應(yīng)的有效邊界加入反傳波場(chǎng)中從而重構(gòu)出該時(shí)間點(diǎn)的炮點(diǎn)波場(chǎng)值;④基于隨機(jī)邊界條件波場(chǎng)重構(gòu)方法[56-60],在計(jì)算區(qū)域外設(shè)置隨機(jī)彈性參數(shù)層,炮點(diǎn)波場(chǎng)正向延拓時(shí)記錄下最后若干時(shí)刻所有網(wǎng)格點(diǎn)的波場(chǎng)值,并以該波場(chǎng)值為初始條件,進(jìn)行炮點(diǎn)波場(chǎng)的逆時(shí)重構(gòu)。

前3種方法一般結(jié)合吸收邊界實(shí)現(xiàn)炮點(diǎn)波場(chǎng)重構(gòu),具有較高的計(jì)算精度,能夠很好地解決二維偏移問題,但在三維情況下,受計(jì)算機(jī)內(nèi)存限制,這些信息都必須存儲(chǔ)到硬盤中,頻繁訪問硬盤不利于發(fā)揮并行機(jī)的計(jì)算能力優(yōu)勢(shì),因此,前3種策略很難用于解決三維逆時(shí)偏移問題。隨機(jī)邊界策略能夠避免頻繁訪問硬盤,但不合理的隨機(jī)彈性參數(shù)無(wú)法完全破壞人工邊界反射波場(chǎng)的相干性,致使隨機(jī)邊界產(chǎn)生的噪聲對(duì)接收點(diǎn)附近地層的成像產(chǎn)生嚴(yán)重影響[57]。因此,邊界處隨機(jī)彈性參數(shù)的構(gòu)造方式是隨機(jī)邊界策略的核心技術(shù)。

解決這一問題的代表性思路是利用隨機(jī)介質(zhì)模型構(gòu)建隨機(jī)邊界[61-65],該隨機(jī)介質(zhì)模型由大、小兩種尺度的非均勻介質(zhì)組成,大尺度非均質(zhì)介質(zhì)是模型的均值速度,小尺度非均質(zhì)介質(zhì)是加在均值速度上的隨機(jī)擾動(dòng)?；谶@種思路,奚先等[64]根據(jù)隨機(jī)過程的譜分解定理,給出了一種利用局部自相關(guān)函數(shù)構(gòu)造非平穩(wěn)隨機(jī)介質(zhì)模型的方法,陳可洋等[65]則利用改進(jìn)的混合型自相關(guān)函數(shù)構(gòu)建各種尺度因子條件下的多層隨機(jī)介質(zhì)模型,張麗美等[66]基于自相關(guān)函數(shù)類型、相關(guān)長(zhǎng)度和速度擾動(dòng)標(biāo)準(zhǔn)差3個(gè)參數(shù)描述的隨機(jī)介質(zhì)表達(dá)方式構(gòu)建參數(shù)化的隨機(jī)邊界,并通過模型試驗(yàn)給出了3個(gè)參數(shù)的取值原則。此外,SHEN等[67]采用形態(tài)隨機(jī)的大尺度隨機(jī)速度顆粒構(gòu)建隨機(jī)邊界,實(shí)現(xiàn)了對(duì)隨機(jī)邊界散射波場(chǎng)低頻分量的壓制。

上述思路在各向同性聲波方程的隨機(jī)邊界構(gòu)建方面取得了效果,但由于TI介質(zhì)的qP波方程包含多個(gè)彈性參數(shù),這些彈性參數(shù)共同決定方程數(shù)值求解的穩(wěn)定性,各向異性參數(shù)的微小波動(dòng)極有可能引起各向異性視速度的很大波動(dòng),使偏移過程發(fā)散。因此,TI介質(zhì)qP波方程隨機(jī)邊界構(gòu)建中各參數(shù)的隨機(jī)值設(shè)置是相互影響的,不能簡(jiǎn)單套用各向同性聲波方程的隨機(jī)邊界構(gòu)建方法,而應(yīng)當(dāng)針對(duì)qP波方程逆時(shí)偏移這一特殊問題,研發(fā)具有完善理論的隨機(jī)邊界彈性參數(shù)賦值方法。劉國(guó)峰等[68]在這方面進(jìn)行了嘗試,對(duì)比分析了3種隨機(jī)參數(shù)設(shè)置方法,即:邊界外擴(kuò)區(qū)域的各彈性參數(shù)互不影響,獨(dú)立處理;外擴(kuò)區(qū)域設(shè)置為各向同性區(qū)域和外擴(kuò)區(qū)域設(shè)置為均勻各向異性區(qū)域。認(rèn)為將外擴(kuò)區(qū)域設(shè)置為均勻各向異性區(qū)域能取得最優(yōu)的模型數(shù)據(jù)成像效果。劉國(guó)峰等[68]將各彈性參數(shù)的隨機(jī)設(shè)置作為一個(gè)整體對(duì)待,其研究思路具有很強(qiáng)的針對(duì)性和理論嚴(yán)謹(jǐn)性,但未給出各個(gè)參數(shù)的具體設(shè)置方法,故將各彈性參數(shù)組合為數(shù)組,在一定的約束條件下產(chǎn)生隨機(jī)數(shù)組,并以該數(shù)組實(shí)現(xiàn)qP波逆時(shí)偏移隨機(jī)邊界的設(shè)置或許是一種可行的解決方案,但目前業(yè)界還缺乏這方面的研究工作。

4 TI介質(zhì)中qP波方程逆時(shí)偏移的波場(chǎng)成像方法

本節(jié)主要介紹qP波方程逆時(shí)偏移中的成像條件、偏移噪聲壓制、波場(chǎng)傳播方向求取和波場(chǎng)拉伸校正等方面的研究現(xiàn)狀。

4.1 成像條件

地震波逆時(shí)偏移領(lǐng)域常用激發(fā)時(shí)間成像條件、互相關(guān)成像條件和能量比成像條件等對(duì)地震波場(chǎng)進(jìn)行成像。其中激發(fā)時(shí)間成像條件無(wú)法解決多路徑情況,速度場(chǎng)復(fù)雜時(shí)容易丟失波場(chǎng)信息;能量比成像條件常用于估算反射系數(shù),在非彈性界面處不穩(wěn)定;互相關(guān)成像條件在各向同性逆時(shí)偏移中最常用,在實(shí)測(cè)資料的處理中取得了良好效果。為提高深部地層的成像效果或壓制逆時(shí)偏移中的低波數(shù)干擾,又在此基礎(chǔ)上提出了炮點(diǎn)波場(chǎng)歸一化的互相關(guān)成像條件[69]和基于行波分離的互相關(guān)成像條件[70],前者可以提升深部地層的照明效果,后者需要依據(jù)傳播方向進(jìn)行波場(chǎng)分解,二者結(jié)合能夠獲得更好的成像效果。因此,常采用公式(1)所示的成像條件實(shí)現(xiàn)波場(chǎng)成像。

(1)

式中:x,y,z分別為三維直角坐標(biāo)系的3個(gè)坐標(biāo);t為時(shí)間;tmax為記錄長(zhǎng)度;I為成像結(jié)果;s為炮點(diǎn)波場(chǎng);r為接收點(diǎn)波場(chǎng);Φ為波場(chǎng)傳播方向;N為用于成像的波場(chǎng)傳播方向個(gè)數(shù)。

公式(1)所示的成像條件本質(zhì)上是一種成像思路或準(zhǔn)則,沒有限定具體的介質(zhì)類型或方程形式,同樣可用于解決qP波的逆時(shí)偏移成像問題。但前提是首先要獲得各成像點(diǎn)在不同時(shí)刻的波場(chǎng)傳播方向并依據(jù)該方向進(jìn)行波場(chǎng)分解,因此如何準(zhǔn)確地分離出沿各個(gè)方向傳播的qP波是該成像條件的關(guān)鍵。郭旭等[71]通過求取VTI介質(zhì)中qP波的坡印廷矢量指示波場(chǎng)傳播方向并據(jù)此實(shí)現(xiàn)不同傳播方向的波場(chǎng)分解,其思路同樣可用于TTI介質(zhì),但坡印廷矢量的求取又要求對(duì)應(yīng)的qP波方程中必須包含計(jì)算坡印廷矢量所需的應(yīng)力分量和不同波型的質(zhì)點(diǎn)震動(dòng)速度分量?，F(xiàn)有的許多qP波方程沒有顯式包含應(yīng)力分量和質(zhì)點(diǎn)震動(dòng)速度分量,必須研究并采用其它技術(shù)求取qP波的傳播方向,目前此方面的研究成果很少。

4.2 偏移噪聲壓制

qP波方程的逆時(shí)偏移噪聲主要指由層間反射導(dǎo)致的低頻噪聲,這種低頻噪聲往往具有很強(qiáng)的能量,污染地層的成像結(jié)果(圖4)。偏移噪聲的本質(zhì)是偏移算法本身的誤差。

圖4 TI介質(zhì)逆時(shí)偏移中的偏移噪聲

目前的偏移噪聲壓制思路有:①依據(jù)逆時(shí)偏移噪聲的特點(diǎn),研究并采用適當(dāng)方法對(duì)已經(jīng)產(chǎn)生的噪聲進(jìn)行壓制處理,如拉普拉斯濾波法和帶通濾波法[72]等;②依據(jù)偏移噪聲的產(chǎn)生機(jī)理,通過改進(jìn)波場(chǎng)延拓方法或成像條件避免產(chǎn)生噪聲,如基于行波分離的互相關(guān)成像方法[70-72]等。前者屬于先產(chǎn)生后壓制的思路,本質(zhì)上是一種信號(hào)處理技術(shù),很難保證不傷害有效信號(hào);后者注重從根源上消除其產(chǎn)生的“土壤”,一般依據(jù)炮、檢波場(chǎng)的傳播方向信息進(jìn)行角度加權(quán)互相關(guān)成像,理論上具有更好的壓制效果。

上述兩種思路同樣沒有介質(zhì)或方程的限制,均可用于TI介質(zhì)qP波逆時(shí)偏移領(lǐng)域,楊富森等[73]和QIN[74]通過模型數(shù)據(jù)驗(yàn)證了這些思路對(duì)各向異性逆時(shí)偏移的適應(yīng)性。

4.3 波場(chǎng)傳播方向求取

求取qP波傳播方向的必要性在于,波場(chǎng)分解或角度加權(quán)互相關(guān)成像條件的應(yīng)用前提是要根據(jù)波的傳播方向?qū)ε邳c(diǎn)和接收點(diǎn)波場(chǎng)進(jìn)行分類,然后只選擇對(duì)最終成像有意義的波場(chǎng)參與互相關(guān)運(yùn)算達(dá)到消除偏移噪聲或提高照明均衡度的效果。目前一般以坡印廷矢量來(lái)指示波的傳播方向,YOON等[75]和郭鵬等[76]給出的地震波坡印廷矢量E的計(jì)算公式為:

E=pv

(2)

式中:p為聲壓;v為質(zhì)點(diǎn)振動(dòng)速度矢量。

由于各向同性介質(zhì)的一階聲波方程本身就包含了p和v這兩個(gè)量,因此其坡印廷矢量的求取是非常方便的,但對(duì)于TI介質(zhì),二階或高階qP波方程很難直接得到p和v,此時(shí)無(wú)法在波場(chǎng)延拓過程中求得坡印廷矢量,只有將TI介質(zhì)中的qP波方程寫為一階應(yīng)力-速度方程組的形式才方便得到qP波的坡印廷矢量,進(jìn)而獲得波的傳播方向信息并據(jù)此對(duì)炮點(diǎn)和接收點(diǎn)波場(chǎng)進(jìn)行分類,因此,采用一階應(yīng)力-速度qP波方程進(jìn)行逆時(shí)偏移有利于提高成像精度。

需要指出的是,地震波的坡印廷矢量是指單位時(shí)間內(nèi)穿過與能量流動(dòng)方向垂直單位截面的能量。當(dāng)同一時(shí)刻在某成像點(diǎn)處存在多個(gè)具有不同傳播方向的波時(shí)(在彈性分界面處這種現(xiàn)象一定存在),公式(2)求得的坡印廷矢量不代表某一特定的波的傳播方向,而是該時(shí)刻所有波的傳播方向的加權(quán)平均。由于逆時(shí)偏移中只需要部分波場(chǎng)參與成像,因此需要求出某一個(gè)或幾個(gè)特定波的傳播方向,坡印廷矢量無(wú)法實(shí)現(xiàn),說明坡印廷矢量法不適用于多組波的情況。

如圖5所示,當(dāng)下行波I在t時(shí)刻傳播至彈性分界面上的A點(diǎn)處時(shí)(假定入射角小于臨界角),在該時(shí)刻會(huì)瞬間產(chǎn)生一個(gè)反射波R和一個(gè)透射波T,也就是說,A點(diǎn)在t時(shí)刻同時(shí)存在3個(gè)不同傳播方向的波,而此時(shí)該點(diǎn)的坡印廷矢量P和上述3個(gè)波的傳播方向都不相同,它其實(shí)是這3個(gè)波傳播方向的加權(quán)平均,說明在多組波情況下,坡印廷矢量法求得的波傳播方向存在較大誤差。這一誤差會(huì)傳遞給后續(xù)的波場(chǎng)分解或角度域權(quán)系數(shù)求取過程,導(dǎo)致最終成像精度的下降,且該誤差在傳遞過程中的收斂性也難以確定,它對(duì)成像結(jié)果的影響目前也難以定量估算。

理論上,采用射線追蹤類的方法可以消除這種誤差,但用射線追蹤法解決該問題的缺陷在于:①?gòu)?fù)雜構(gòu)造條件下,TI介質(zhì)的射線追蹤技術(shù)難以解決多路徑問題;②射線追蹤法是一種“模型驅(qū)動(dòng)”類算法,它的計(jì)算精度嚴(yán)重依賴彈性參數(shù)模型的精度,實(shí)際數(shù)據(jù)處理中給定的彈性參數(shù)模型存在較大誤差,模型的誤差會(huì)傳遞給射線追蹤過程,最終降低波場(chǎng)成像精度,并且該誤差在傳遞過程中的收斂性同樣難以確定;③在波場(chǎng)延拓過程中必須同時(shí)進(jìn)行射線追蹤運(yùn)算,增加了計(jì)算量與算法復(fù)雜度。這說明在TI介質(zhì)逆時(shí)偏移領(lǐng)域,用射線追蹤法解決波場(chǎng)傳播方向求取問題具有不適應(yīng)性,很難取代坡印廷矢量法。

圖5 多路徑情況下的坡印廷矢量示意

盡管坡印廷矢量法存在上述缺陷,但目前仍然是地震波逆時(shí)偏移領(lǐng)域的主流算法。進(jìn)一步提升波場(chǎng)分解精度的關(guān)鍵在于如何降低多組波情況下計(jì)算誤差,目前業(yè)界在這方面沒有可行的解決思路并缺乏深入的研究工作。

4.4 波場(chǎng)拉伸校正

逆時(shí)偏移的子波拉伸效應(yīng)是指偏移后反射波在垂向上的子波延續(xù)長(zhǎng)度和頻率隨地震波速度、地層傾角和反射角的變化而變化,造成共成像點(diǎn)道集(Common Imaging Gather,CIG)或角度域共成像點(diǎn)道集(Angle Domain Common Imaging Gather,ADCIG)中不同偏移距的地震道具有不同的頻譜和分辨率的現(xiàn)象。圖6為均勻TI介質(zhì)中qP波方程逆時(shí)偏移脈沖響應(yīng)分別用不同方式顯示的結(jié)果(該脈沖響應(yīng)的炮點(diǎn)位于地表2000m位置處),采用灰度顯示時(shí)(圖6a)看不出拉伸現(xiàn)象,但對(duì)該數(shù)據(jù)抽稀并采用波形顯示時(shí)(圖6b)則可以看到明顯的拉伸現(xiàn)象,偏移距越大,波場(chǎng)在深度方向上的持續(xù)長(zhǎng)度越大。說明逆時(shí)偏移的子波拉伸效應(yīng)不是偏移算法本身的誤差,而是一種視角效應(yīng),圖6a的視角是波前的傳播方向,是偏移過程的視角,圖6b的視角是垂直向下方向,是偏移剖面的視角,說明偏移過程和偏移結(jié)果視角的不同是造成逆時(shí)偏移波場(chǎng)拉伸效應(yīng)的原因。

目前解決逆時(shí)偏移子波拉伸問題的思路有3種:①先偏移后壓縮[77-78],先用常規(guī)偏移方法得到拉伸后的偏移結(jié)果,再研發(fā)相關(guān)技術(shù)對(duì)已經(jīng)拉伸的數(shù)據(jù)進(jìn)行壓縮;②先壓縮后偏移[79-80],在偏移前先對(duì)地震炮記錄進(jìn)行時(shí)間上的壓縮,然后再偏移成像;③在偏移過程中壓縮[81-82],即在偏移過程中依據(jù)拉伸效應(yīng)的產(chǎn)生機(jī)理采用適當(dāng)方法避免產(chǎn)生拉伸。第一種思路的本質(zhì)是試圖恢復(fù)已經(jīng)被破壞或損失的信息,很難達(dá)到理想效果,對(duì)于第二種思路,由于在偏移前壓縮地震記錄會(huì)破壞地震數(shù)據(jù)原來(lái)的動(dòng)力學(xué)特征,而逆時(shí)偏移本身并不具備恢復(fù)這些動(dòng)力學(xué)信息的能力,必然以降低偏移結(jié)果的保真性為代價(jià),同時(shí)這種思路在操作上也存在巨大困難;第三種思路具有理論基礎(chǔ)的先天優(yōu)勢(shì),在理論上能夠減弱子波拉伸對(duì)偏移結(jié)果的影響。

圖6 qP波方程逆時(shí)偏移的脈沖響應(yīng)

ZHU等[81]依據(jù)逆時(shí)偏移子波拉伸的產(chǎn)生機(jī)理給出了各向同性聲波方程逆時(shí)偏移的子波拉伸校正方法,算法的本質(zhì)是在偏移過程中進(jìn)行視角轉(zhuǎn)換,只要地震波入射角和反射角求取準(zhǔn)確,這種算法就能消除逆時(shí)偏移中的波場(chǎng)拉伸現(xiàn)象。楊佳佳等[82]改進(jìn)了反射角和入射角的計(jì)算方法,將該算法擴(kuò)展到彈性波領(lǐng)域,并在模型數(shù)據(jù)和實(shí)測(cè)資料的處理中均取得了效果。

ZHU等[81]和楊佳佳等[82]的工作對(duì)于qP波逆時(shí)偏移的波場(chǎng)拉伸校正雖然具有很強(qiáng)的借鑒意義,但不能完全照搬,qP波逆時(shí)偏移中波場(chǎng)拉伸校正的特殊性表現(xiàn)在:①算法只適用于激發(fā)時(shí)間成像條件,而激發(fā)時(shí)間成像條件對(duì)深部地層的成像能力不足;②地震波入射角和反射角的求取需要已知波的傳播方向信息,在各向同性聲波領(lǐng)域,可以方便地利用坡印廷矢量來(lái)指示波的傳播方向,但目前的許多qP波方程不能提供求取坡印廷矢量必需的振動(dòng)速度和應(yīng)力信息,導(dǎo)致坡印廷矢量求取失敗。因此,針對(duì)具體qP波方程本身的特點(diǎn),研發(fā)適用于各向異性qP波方程的逆時(shí)偏移波場(chǎng)拉伸校正技術(shù)極為必要,但截止目前,還未見到該方程的研究成果與文獻(xiàn)報(bào)道。

5 存在的主要問題及未來(lái)發(fā)展展望

綜上所述,本文認(rèn)為在TI介質(zhì)中的qP波逆時(shí)偏移領(lǐng)域尚存在以下問題需要解決:

1) 對(duì)qP波方程的動(dòng)力學(xué)精度缺乏深入研究。由于qP波方程是彈性波方程的近似,其近似效果的評(píng)價(jià)應(yīng)該同時(shí)包括運(yùn)動(dòng)學(xué)和動(dòng)力學(xué)精度兩個(gè)方面,而目前的研究工作均停留在方程運(yùn)動(dòng)學(xué)精度的分析層面,對(duì)方程動(dòng)力學(xué)精度的研究尚未展開。運(yùn)動(dòng)學(xué)特征的近似程度決定了構(gòu)造成像的準(zhǔn)確度,動(dòng)力學(xué)特征的近似程度決定了偏移結(jié)果中的能量信息能否用于反演或解釋巖性,如果對(duì)方程動(dòng)力學(xué)精度的認(rèn)識(shí)不足,則會(huì)導(dǎo)致對(duì)偏移結(jié)果中地震波能量信息的利用率不足或應(yīng)用失敗。

2) 對(duì)qP波方程逆時(shí)偏移中的隨機(jī)邊界技術(shù)缺乏深入研究。因?yàn)殡S機(jī)邊界技術(shù)的炮點(diǎn)波場(chǎng)逆時(shí)重構(gòu)是提升逆時(shí)偏移效率的關(guān)鍵,在各向異性情況下,由于qP波方程包含多個(gè)各向異性參數(shù),這些參數(shù)共同決定求解算法的穩(wěn)定性,所以隨機(jī)邊界中各向異性參數(shù)的設(shè)置是相互關(guān)聯(lián)的。如果用各向同性領(lǐng)域的隨機(jī)值設(shè)置方法分別獨(dú)立設(shè)置邊界處的各向異性參數(shù),會(huì)導(dǎo)致延拓或重構(gòu)過程出現(xiàn)不穩(wěn)定現(xiàn)象。因此,各向同性逆時(shí)偏移中的隨機(jī)邊界技術(shù)不能直接用于解決各向異性逆時(shí)偏移問題,而應(yīng)當(dāng)研發(fā)具有完善理論的隨機(jī)邊界各向異性參數(shù)賦值方法。

3) 在qP波方程逆時(shí)偏移中波場(chǎng)傳播方向的求取方面需要展開更深入研究。不同時(shí)刻波場(chǎng)傳播方向信息的準(zhǔn)確求取是波場(chǎng)成像、偏移噪聲壓制以及波場(chǎng)拉伸校正的關(guān)鍵,實(shí)際求解中可以用坡印廷矢量指示波的傳播方向,但坡印廷矢量的求取要求首先已知各成像點(diǎn)在不同時(shí)刻的應(yīng)力矢量和質(zhì)點(diǎn)振動(dòng)速度矢量,對(duì)qP波方程提出了以下要求:方程中最好能包含應(yīng)力矢量與振動(dòng)速度矢量;或便于求解不同時(shí)刻各質(zhì)點(diǎn)的應(yīng)力矢量與振動(dòng)速度矢量;最好是一階方程且便于快速穩(wěn)定求解。同時(shí),現(xiàn)有的坡印廷矢量技術(shù)還存在當(dāng)同一點(diǎn)在同一時(shí)刻有不同傳播方程的多個(gè)波時(shí),坡印廷矢量只能得到一個(gè)混合波場(chǎng)傳播方向的缺點(diǎn)。因此,必須深入研究更準(zhǔn)確的波場(chǎng)傳播方向求取技術(shù)。

6 結(jié)束語(yǔ)

本文初步總結(jié)了TI介質(zhì)中qP波逆時(shí)偏移技術(shù)的研究現(xiàn)狀,分析了現(xiàn)有技術(shù)的實(shí)現(xiàn)思路和存在的主要問題,在此基礎(chǔ)上給出了進(jìn)一步研究建議。

各向異性理論的qP波逆時(shí)偏移是一個(gè)龐大的技術(shù)體系,受文章篇幅和作者水平的限制,本文無(wú)法就這一技術(shù)體系中每個(gè)細(xì)節(jié)的研究現(xiàn)狀與存在的問題展開討論,只能討論該體系中部分重點(diǎn)環(huán)節(jié),即使在這些重點(diǎn)討論的環(huán)節(jié)中也可能因?yàn)樽髡咚絾栴}遺漏其中的重要方法和文獻(xiàn)。

本文在總結(jié)過程中沒有考慮各向異性建模問題,而這一問題又是逆時(shí)偏移技術(shù)應(yīng)用前必須解決的,各向異性建模中的相關(guān)技術(shù)和難點(diǎn)問題需要另文討論。

地震各向異性資料的成像處理是一個(gè)系統(tǒng)工程,不能只靠逆時(shí)偏移一個(gè)環(huán)節(jié)解決問題,還應(yīng)當(dāng)包括去噪、靜校正、反褶積和彈性參數(shù)建模等諸多環(huán)節(jié),只有解決了各個(gè)環(huán)節(jié)的問題才有可能最終解決各向異性地層的地震波成像問題。