數據驅動的控制束偏移方法

呂考考 徐基祥 張 才 李凌高 孫夕平

(中國石油勘探開發研究院，北京 100083)

0 引言

射線束偏移是一種靈活、高效且精度較高的深度域成像方法。它克服了傳統基于射線的Kirchhoff偏移方法的某些缺陷，如高斯束偏移使用復值的初始束參數解決了常規射線追蹤振幅在焦散區奇異性問題[1-3]；高斯束偏移在局部平面波假設前提下使每條高斯束獨立傳播，實現了多波至成像[4-5]，并能對陡傾構造成像。

在高斯束理論研究的基礎上[6-8]，Hill[9-10]先后提出疊后和疊前高斯束偏移算法。后來，Hill的成果被拓展到共炮點道集[11-12]、各向異性介質[13-14]以及彈性介質[15-17]。國內學者在高斯束偏移方面也做了大量研究，如李瑞忠等[18]利用局部傾斜疊加法實現了高斯束疊后偏移；李振春等[19]、岳玉波等[20-21]實現了角度域、炮域和復雜地表條件的保幅高斯束深度偏移。此外，黃建平等[22-25]、韓冰凱等[26]、吳娟等[27]和代福材等[28]也針對高斯束偏移技術進行了深入研究并發表了相應著述。

經過幾十年的快速發展，高斯束偏移在理論上已經得到很好的證實，并作為Kirchhoff偏移的一種替代方法逐漸應用于實際地震數據處理。射線束偏移中一個重要步驟就是對疊前數據局部傾斜疊加后形成τ-p域道集上的同相軸進行選擇，該同相軸最終貢獻于最后的成像。盡管高斯束偏移克服了Kirchhoff偏移的一些缺點，但傳統的高斯束偏移在成像原理上與Kirchhoff類似，也是將所有反射能量投影到旅行時橢圓等時線上，這就意味著將τ-p域道集上的所有點(同相軸)進行偏移，不僅增大了計算量，而且還會產生一些畫弧噪聲和偏移假象，尤其是在偏移疊加效果不好的復雜構造區域。

針對這個問題，很多學者做了有益的研究。Gao等[29]發展了快速束偏移方法，顯著提高了偏移效率，但該方法在處理過程中會造成有效信號能量的缺失，很難解決復雜構造的精確成像問題。Vetle等[30]詳細介紹了控制束偏移方法的優勢，該方法是由CGG Veritas開發的高斯束偏移增強版，但未給出其具體原理和實現方法。Sherwood等[31]和Ting等[32]對控制束偏移方法進行了大量的研究和實際應用分析。Hu等[33]將慢度信息與高斯束偏移方法結合實現了一種慢度驅動的高斯束疊前深度偏移方法，該方法雖在某種程度能壓制一些相干噪聲，但它僅在共炮點域計算了檢波點的水平慢度信息，忽略了對共檢波點域炮點慢度信息的使用。為充分利用疊前地震數據炮點和檢波點慢度信息，Yang等[34-35]提出一種基于優化策略的數據驅動的高斯束偏移方法，可有效地壓制傳統高斯束偏移所產生的噪聲，且具實用性。

本文基于控制束偏移的思想，結合疊前數據的炮點和檢波點水平慢度信息，并設計一個偏移質量控制因子進行偏移同相軸的選取，對偏移中引起的噪聲(畫弧和相干噪聲等)和假象進行壓制從而提高成像質量。具體實現上，首先對局部共炮點和共檢波點道集進行局部傾斜疊加，并利用相干分析在τ-p域拾取炮點和束中心位置水平慢度信息；然后設計一個偏移質量控制因子對同相軸進行篩選(去除非一次反射同相軸)；最后采用經典高斯束偏移的算法進行成像。三個理論模型數據試驗和實際資料應用證明了本文方法的可行性、有效性以及適用性，并且對比了本文方法和傳統高斯束方法在Mar-mousi模型上的計算效率。

1 方法原理

1.1 射線參數拾取

最常用的提取水平射線參數信息的方法有三種：平面波分解[36]，多道互相關和相干分析[37]。Chopra等[38]指出，即使在噪聲嚴重的數據中，相干性分析依然是一種穩健并能提供高分辨率同相軸連續性剖面的方法。所以，本文采用相干性分析在τ-p域道集來估計水平射線參數信息，除了使用實值道集以外，還利用數據的解析道集計算相干系數，其公式可寫為

(1)

式中：τ表示截距時間；p表示水平慢度或射線參數；x0為中心道的坐標，xi表示第i個局部道集的坐標；n表示一個時窗寬度半徑；u(t,xi)和uH(t,xi)分別表示地震數據及其希爾伯特變換后數據；Δt為采樣間隔。隨著時窗寬度的增大，雖然相干性計算越穩定，但是局部同相軸監測的分辨率會越來越低，所以本文使用的時窗寬度為1。

本文使用炮點水平射線參數(ps)、檢波點水平射線參數(pr)及反射時間(t)三個參數表征地震數據中的一個局部同相軸。根據炮檢點互換原理，需拾取相同的同相軸，即在共檢波點域道集的τ-p域和共炮點域道集的τ-p域上分別拾取具有相同反射時間的ps和pr。圖1是使用相干分析拾取上述三個參數(x或z方向分量)的簡單示例，紅叉指示局部最大值位置，即相干能量最強的參數位置。

1.2 偏移質量控制因子

實際中有許多因素影響射線束偏移成像質量，本文假設速度模型和拾取的射線參數等都是準確的，唯一考慮的因素是成像的聚焦程度，所以使用兩條射線的距離作為偏移質量控制因子，如圖2所示。在射線束偏移中，射線末端交點所估計的偏移位置，在該深度上炮點旅行時和檢波點旅行時之和ts+tr等于反射時間t。如果兩個射線末端點之間的距離太大，那么這個結果不能用于疊加成像。由于射線束偏移不是點對點成像，而是使用一個偏移“波形”代替這個成像點。所以，為了不損失有效信號能量，本文設置一個距離范圍對偏移的同相軸進行選擇，即成像距離在這個范圍內的同相軸。控制因子的表達式可寫為

t-αΔt

|ds-dr|<βΔx

(2)

式中：ts、tr和t分別表示炮點、檢波點射線旅行時和反射時間；Δt表示時間間隔；Δx表示速度模型水平間隔；ds和dr分別表示兩條射線末端位置坐標；α和β是兩個可控參數，通過設置這兩個參數控制時間和距離大小范圍。

圖1 射線參數拾取

圖2 偏移質量控制因子

本文將α的值設為4、β的值設為6，已知反射時間t，利用拾取的炮點射線參數ps和檢波點射線參數pr發射射線得到兩條射線交點位置附近的炮點射線旅行時ts和檢波點射線旅行時tr，以及兩條射線末端位置坐標ds和dr，最后判斷兩條射線的旅行時之和以及末端坐標位置之差是否滿足式(2)，據此篩選待成像的炮點和檢波點射線參數。

1.3 高斯束偏移公式

根據Clearbout[39-40]提出的互相關成像原理，反射界面存在于地下這樣一些點上，在這些點上，下行波的波前到達或產生與上行波的波前到達或產生在時間上是一致的。互相關成像公式可寫為

(3)

式中：ω為角頻率；xs表示炮點坐標；D(x,xs,ω)表示下行波場；U(x,xs,ω)表示上行波場；上標*表示復數共軛。

上行波場可由地表接收的地震波場向下延拓得到，據Hill[10]的研究可知

(4)

式中：xr表示檢波點坐標；G(x,xs,ω)表示從xr到x的格林函數；U(xr,xs,ω)表示地表接收到的地震記錄。

下行波場可近似表示為格林函數，即D(x,xs,ω)≈G(x,xs,ω)，則共炮域的疊前偏移成像公式表示為

(5)

在高斯束偏移中，格林函數由一系列從不同角度出射的高斯束的疊加積分求得

(6)

式中：uGB(x,xs,ω)表示頻率域高斯束函數；A(x,xs)和T(x,xs)分別表示高斯束的復值振幅和旅行時間；psx和psz分別表示水平慢度和垂直慢度。

為了使接收到的地震波場與高斯束表示的波場相匹配，需要進行加窗處理，根據Hill[10]給出的高斯窗函數

(7)

式中：w0為初始束寬；ωr為參考角頻率；L為束中心位置；ΔL表示束中心間隔。

將式(6)、式(7)代入式(5)，并引入相位校正因子，消除檢波點和束中心位置不一致的影響，得到

(8)

式中：A(xr,xs)=A(x,xr)A(x,xs)、T(xr,xs)=T(x,xr)+T(x,xs)，均為復值；Ds(L,prx,ω)表示地震記錄的局部傾斜疊加，其表達式為

(9)

根據Hale的成像方法[41-42]，對式(8)進行化簡，最終的成像公式為

[Ar(xr,xs)D(L,prx,t)-

Ai(xr,xs)DH(L,prx,t)]

(10)

式中：D(L,prx,t)為時間域的局部傾斜疊加道集；DH(L,prx,t)為其對應的希爾伯特變換道集。

綜上所述，本文數據驅動的控制束偏移原理和算法實現流程分別如圖3和圖4所示。

圖3 數據驅動的控制束偏移原理

圖4 數據驅動的控制束偏移實現流程

2 方法驗證

為了說明本文方法的正確性、適用性以及計算效率，分別使用三個理論模型和某工區的實際資料對本文算法效果進行驗證，并且對比了本文方法與傳統高斯束偏移方法在Marmousi模型上的計算效率。

第一個理論模型是簡單的三層層狀模型，其單炮記錄如圖5a所示。對比使用傳統高斯束偏移 (圖5b)和使用本文控制束偏移(圖5c)成像結果，在后者中未見畫弧(圖5b紅色箭頭所示)等噪聲。

第二個理論模型是洼陷模型(圖6a)。模型尺寸為640m×375m，水平網格間距為15m，垂向網格間距為8m。正演記錄由二階有限差分計算得到。觀測系統為中間激發、兩邊接收，共121炮，炮間距為30m，每炮接收道數為121，道間距為30m。每道采樣點數為750，采樣間隔為4ms，其單炮記錄如圖6b所示。圖6c和圖6d是對第61炮記錄分別使用傳統高斯束偏移和本文控制束偏移的結果。圖6e和圖6f是對所有炮記錄分別使用傳統高斯束和本文控制束偏移結果。在圖6c和圖6e中可發現主要包括廣角反射疊加噪聲和畫弧噪聲等(紅色箭頭)大量偏移噪聲。使用本文控制束偏移方法進行成像，這些噪聲可被很好地壓制，得到較清晰的偏移剖面(圖6d和圖6f)。

圖5 層狀模型的偏移成像

圖6 洼陷模型的偏移成像

使用Marmousi模型分析本文方法對復雜構造的成像效果。圖7a是Marmousi速度模型，圖7b是該模型的第一炮地震記錄。圖7c和圖7d分別是該單炮的傳統高斯束偏移結果和本文的控制束偏移結果，從兩圖對比可以看出在本文方法成像結果中偏移噪聲和一些假象得到了很好的壓制(紅色箭頭)。圖7e和圖7f是所有炮偏移的結果，由圖可見，本文方法可以有效地壓制傳統高斯束偏移中所產生的噪聲和偏移假象。然而，由于成像中舍棄了一些弱能量信號，在大傾角斷面成像方面，控制束成像結果能量較弱，表明控制束偏移方法目前還有一定的局限性。

為了檢驗本文方法對實際資料的適用性及有效性，選用采集于M山地起伏地表的地震數據進行偏移成像。圖8a是通過速度建模得到的該工區的速度模型。圖8b是某一炮的地震記錄，可以看出信噪比較低。圖8c和8d分別是該單炮的傳統高斯束偏移結果和本文的控制束偏移結果，圖8d中沒有圖8c中出現的偏移噪聲和假象，并且圖8d中的有效反射能量(紅色箭頭)更加明顯。圖8e和圖8f是所有炮的偏移結果，圖8f相比于圖8e具有更高的信噪比和更強的同相軸連續性。

圖8 實際資料的偏移成像

3 結束語

本文充分利用炮點和檢波點的水平慢度信息并設計一個偏移質量控制因子，在傳統高斯束偏移的基礎上，發展了一種數據驅動的控制束偏移算法流程。從共檢波點和共炮點的τ-p域道集拾取水平慢度信息，利用相干分析提高了拾取精度；通過偏移質量控制因子對偏移同相軸進行篩選。這些做法不僅可壓制傳統高斯束偏移所產生的噪聲和偏移假象，且可提高計算效率。

模型試驗和實際資料應用結果表明，本文數據驅動的控制束偏移算法能有效壓制低信噪比地震數據成像噪聲，提高復雜山地低信噪比數據的深度域成像質量，特別是能顯著改善干擾嚴重的淺層的成像效果。

炮點和檢波點道集優勢同相軸提取是本文算法的關鍵之一，后續將深入研究基于壓縮感知技術的稀疏平面波數據分解方法，以進一步提高數據分解和提取精度。

感謝科羅拉多礦業學院CWP提供的SU軟件平臺支持。