快速回轉波近地表速度建模方法

郭振波 孫鵬遠 錢忠平 李培明 唐博文 熊定鈺

(東方地球物理公司物探技術研究中心，河北涿州 072751)

0 引言

地震勘探由于激發和接收均在地表，因此需要消除地表起伏及近地表速度橫向變化對反射信號的影響[1]，其主要處理方式有兩種：①通過靜校正技術將炮點、檢波點校正到統一基準面上[2-3]，在消除近地表對反射波影響的同時也可使一些基于水平地表假設的方法(如水平疊加等)可行；②直接進行基于起伏地表的處理，在算法內部考慮近地表的影響(如起伏地表條件下的速度建模、偏移成像等)[4-5]。不管采用哪種方式，都需要進行近地表建模，因此近地表建模方法的研究對油氣地震勘探具有重要意義。

近地表建模主要利用地震記錄中的初至波，可分為基于波動理論和射線理論的層析，也可分為波形層析和旅行時層析[6]。由于觀測誤差、噪聲等因素的影響，很難直接進行波形匹配，因此目前以旅行時層析為主[7]。由于計算能力的限制，旅行時層析以基于射線理論的方法為主，而基于波動理論的旅行時層析還未得到大規模應用[8]。

常規的基于射線理論的旅行時層析方法(常規射線層析)[9]是基于射線追蹤方程或程函方程、通過迭代求解正演旅行時與觀測旅行時的最優擬合，可獲取較高精度的近地表模型。但由于需要大量的正演及反演計算，計算量較大，特別是對于目前的高密度、寬方位地震數據，限制了方法的大規模應用。除了常規射線層析反演，還有一類方法是通過對速度分布做進一步假設構建快速近地表建模方法，該類方法在天然地震研究中廣泛應用，通常用來構建一維參考模型[10]。Diebold等[11]采用截斷時間求和方法利用折射波與反射波旅行時構建一維速度模型； Rühl等[12]利用最大深度方法構建近地表模型并將其應用于靜校正，實際數據測試驗證了由該方法計算的長波長靜校正量的準確性； Nowack[13]、Novotn等[14]、胡自多等[15]、徐濤等[16]利用 Herglotz-Wiechert方法構建近地表模型； Shi等[17]在線性速度分布假設下采用層剝離的方法進行快速近地表建模。

本文通過假設速度隨深度線性變化發展了一種快速的回轉波近地表建模方法。通過在算法內部采用多基準面校正降低地表起伏對反演結果的影響，提高反演精度；應用基于局部加權的穩定射線參數估計方法降低初至拾取誤差對反演結果的影響，增強反演的穩定性。理論模型數據與實際數據測試驗證了方法的有效性及高效性。

1 方法原理

1.1 基本原理

假設在最大炮檢距范圍內近地表速度橫向不變但隨深度線性變化，此時近地表的速度變化趨勢可表述為

v(z)=v0+gz

(1)

式中：v0為地表速度；z為深度；g為速度梯度。在該介質中傳播的地震波稱為回轉波。給定射線入射角度θ、起始坐標(x0，z0)可導出其射線路徑、炮檢距、初至時間的理論表達式(圖1)。具體來說，其射線路徑為一圓弧，圓心位于點[18]

(2)

半徑為

(3)

式中p=sinθ/v為射線參數。相應地可求得對應的炮檢距H及初至時間t

(4)

對于速度反演，已知炮檢距、初至時間，求取地下不同位置處的速度參數。假設速度水平橫向不變，地震波的傳播射線參數保持不變，可據此求得回轉點處的速度為

圖1 地震波在速度隨深度線性變換介質中的傳播路徑

(5)

由式(2)和式(3)可求得回轉點的深度為

(6)

射線參數p可由初至時間估算； 地表速度v0由近炮檢距數據計算或由先驗信息給定。速度隨深度變化的梯度通過求解目標函數

(7)

獲得。式中W為加權因子，調整炮檢距項(第一項)與時間項(第二項)的權重。當炮檢距單位為m、時間單位為s時，加權因子W設為1000，以均衡由于不同測量單位固有的加權效應。

利用式(5)～式(7)可獲取一個炮檢距—初至時間對所對應回轉點的速度及其深度。通過在由近及遠的炮檢距范圍內進行初至反演，可獲取由淺到深的速度參數；通過對不同空間位置點進行相同的處理，可獲取整個三維空間的速度參數。給定離散的空間網格，將相應的速度參數投影到對應的網格點上，然后進行對應的速度內插、平滑等處理，可得到規則網格點上最終的速度參數。

1.2 多基準面初至時間校正

由于上述方法基于水平地表假設，起伏地表情況下需要將其校正到水平基準面上。如圖2所示，本文通過兩方面的策略保持校正精度：①將數據抽至CMP道集后再進行后續的處理，保證炮檢點的空間局部性；②根據地表起伏的情況，針對不同的炮檢距范圍利用不同的基準面進行高程校正，盡量縮小校正時間以減少高程校正帶來的誤差。一般情況下，近炮檢距范圍采用CMP點處的高程作為反演的基準面，遠炮檢距范圍內采用炮檢點的平均高程作為反演基準面。

圖2 多基準面校正示意圖

1.3 基于局部加權的穩定射線參數估計

射線參數又稱為視慢度，可通過估算炮檢距—初至時間曲線不同點的斜率獲取。射線參數估計的精度及穩定性直接影響最終的反演結果。本文提出了一種基于局部加權的穩定射線參數估計方法，可有效降低初至拾取誤差、異常值的影響，實現相對穩定的射線參數估計。

對于一個CMP道集，利用一個固定寬度的窗口沿炮檢距—初至時間進行滑動，對平滑窗內的數據進行直線擬合[19]，所擬合直線的斜率即為該平滑窗中心位置對應的射線參數，如圖3a所示。利用這種直線擬合的方式可消除小的拾取誤差對射線參數估計的影響，但是對于大的拾取誤差容錯能力較弱。借鑒Cleveland[20]局部加權平滑的相關思想，本文在直線擬合的過程中引入局部加權減少大的拾取誤差對射線參數估計的影響。

基于局部加權的直線擬合，可以表述為求解目標函數[20]

(8)

式中：N為平滑窗內的點數；tk、Hk為第k個點處的初至時間與炮檢距；t0、p分別為所擬合直線的截距時間與斜率(即射線參數)；wk為第k個點處的加權系數。

穩定的參數估計方法主要分為兩步：①將所有點處的加權系數wk設為1，首先進行常規線性擬合；②計算所有點與擬合直線的距離，根據給定的可容許的最大時間差Δtmax，重新設置其加權系數為

(9)

利用新的加權系數再進行一次線性擬合，消除異常初至對射線參數估計的影響，如圖3b所示。

圖3 射線參數估計示意圖

2 理論數據測試

為了說明本文方法的有效性，選用Amoco靜校正基準測試模型[21](圖4)進行理論數據的測試。該模型包含了大部分常見的近地表地質構造，如高速層出露(區域A)、局部高速、低速異常體(區域B)、淺層低速層(區域C)、近地表復雜構造(區域D)及極淺層低速體(區域E)等，可在一定程度上說明該方法對不同近地表地質構造的適應性。

圖4 Amoco靜校正基準測試模型

原有模型橫向范圍為50km，縱向深6km，由于本文只進行近地表模型的建立，因此只截取了淺部3.2km的部分模型進行測試。模型縱、橫向采樣間隔均為5m，橫向共10001個采樣點，縱向共647個采樣點，最大速度為5760 m/s，最小速度約為800m/s。觀測數據共1998炮，炮點以25m的間隔均勻分布于地表以下10m處，第一個炮點位于橫向位置10m處。按照陸上采集方式布設觀測系統，最大炮檢距為7.5km，檢波點間距為10m。為了避免初至拾取巨大的工作量以及拾取誤差對分析的影響，利用有限差分求解程函方程計算得到各個檢波點處的初至時間作為反演的輸入。

本文方法反演的速度模型如圖5所示。對比圖4與圖5可知，本文方法除了在復雜構造區域(圖4中D標識區域)反演效果較差外，其余部分反演結果與真實模型形態均相似。對于高速層出露區域(圖4中A標識區域)反演結果具有與真實相似的構造形態。高速、低速體標識區域(圖4中B標識區域)，近地表低速異常均得到有效的反演，部分高速層的構造形態較為清晰。由于該方法基于速度線性遞增假設，在反演之后對沿深度出現速度反轉的部分進行人為處理。B標識區域中近地表高速異常構造形態上的差異均是由該處理引起的。淺層低速層區域(圖4中C標識區域)是最常見的地質構造，反演結果與真實模型結果較為接近，反演精度較高。近地表復雜構造區域(圖4中D標識區域)，利用簡單的線性速度假設已不能對其進行很好的近似，因此反演結果與真實模型差距較大，但極淺層的速度變化還是得到了一定的反映。極淺層低速體區域(圖4中E標識區域)類似常規山地地質情況，地表存在一極淺低速層，下面緊鄰致密高速巖石。通過對比圖4與圖5中E標識區域可知，在這種地質條件下該反演方法的反演結果同樣是可接受的。

圖5 原始數據回轉波層析結果

將圖5所示反演結果作為初始模型，然后利用射線層析做進一步的反演優化，得到最終的射線層析反演結果(圖6)。對比圖6與圖4可知，反演結果與真實模型非常接近，說明本文方法得到的模型可作為射線類初至旅行時層析的初始模型。相對于常規的梯度類初始模型，本文方法具有如下優勢：①結果與真實模型更為接近，將其作為初至旅行時層析的初始模型可減少迭代次數，對于此次試驗，相對于梯度類初始模型，可節省約5次迭代；②計算效率極高，該理論數據處理僅耗時25s，但反演結果卻反應了近地表速度變化的大體結構，可對進一步初至旅行時層析反演的參數設定起到參考作用。

圖6 將圖5所示結果作為初始模型的常規射線層析結果

對比圖5與圖6可知，除了A、D標識區域因為構造過于復雜不滿足假設條件而引起較大差別外，其余部分反演結果均可接受。對于近地表復雜構造區域，由于反演分辨率的限制，即使采用常規射線層析方法也很難獲取較高精度的反演模型。分別抽取真實模型、本文方法反演模型、初至旅行時層析反演模型地表以下50和100m處的速度曲線，如圖7所示。對比分析圖5～圖7可知：①本文方法以極小的計算成本快速獲取了一個精度較高近地表模型；②該近地表模型在未出現速度反轉區域反演精度較高，在存在速度反轉區域可靠性變差；③該速度模型可用作初至旅行時層析反演的初始模型以提高其收斂速率，在構造相對簡單區域，該模型可直接用于層析靜校正等處理。

圖7 地表以下不同深度處速度曲線對比

為了驗證拾取誤差對反演結果的影響，在初至中加入噪聲進行相應的數值測試。加入的噪聲主要分為兩類，一類加在所有數據上，為幅值在[-12ms，12ms]范圍內的隨機噪聲，模擬小的拾取誤差； 第二類為在隨機選取50%炮上加入幅值在[-0.3s，0.3s]范圍內的隨機噪聲，模擬大的異常初至。圖8為其中三炮加噪前、后的初至對比。

圖9為利用加入噪聲之后的初至得到的反演結果。對比圖9與圖5可知，加入隨機噪聲對反演結果基本沒有影響，驗證了基于局部加權的穩定射線參數估計方法的抗噪性和有效性。

圖8 原始初至與加入隨機噪聲后的初至對比

圖9 利用加入噪聲初至進行回轉波層析的結果

3 實際資料測試

為了驗證本文所述方法在實際資料上的應用效果，選取中國西部M工區三維實際資料中的兩束線進行測試。該工區地表最大高差約為400m，地表相對較為平緩。由數據分析可知，近地表存在較為明顯的低速異常體以及低速層厚度變化。試驗所用數據共867炮、8352個檢波點，道間距為30m。定義觀測系統后，共68條CMP線，CMP線間距為30m；每條CMP線共800個CMP，CMP間距為15m。采用自動拾取與人工修正相結合的方法進行初至拾取，利用拾取后的初至分別進行本文方法以及常規射線類層析的近地表速度反演(圖10)。

對比圖10a與圖10b可知，不管是地表附近的低速異常還是高速層的基本形態，本文方法與常規射線層析方法反演得到的速度模型基本一致。分炮檢距疊加是驗證是否解決長波長靜校正的一種常用

圖10 實際數據本文方法(a)和常規射線層析(b)反演結果

質控手段，若不同炮檢距段疊加剖面橫向連續性好，沒有因為疊加所用炮檢距范圍的變化而出現疊加剖面上同相軸的錯動，則說明長波長靜校正量得到了有效解決。圖11為利用本文方法進行層析靜校正之后的分炮檢距疊加結果，可以看出，同相軸連續性好，從另一個側面說明了本文所述方法反演得到的模型能夠滿足層析靜校正的需求。

由表1的效率對比可見，本文方法比常規的射線類層析方法快了兩個數量級，大大降低了近地表建模對大規模計算資源的要求，可用于大數據量處理，實現快速近地表速度建模。

圖11 本文方法層析靜校正后分炮檢距疊加結果

本文方法常規方法(10次迭代)計算資源1個節點(32線程/節點)15個節點(32線程/節點)計算時間/min147等價單節點耗時/min1705

4 討論

本文方法之所以能夠進行快速的近地表建模，主要依賴于對近地表速度分布局部橫向不變及線性遞增假設或者對射線路徑的圓弧狀假設?？蓮囊韵聨讉€方面討論這種假設的合理性。

(1)關于速度局部橫向均勻假設。地質構造的尺度通常要比地震觀測最大炮檢距的尺度大得多，在最大炮檢距范圍內將其假設為橫向均勻是可接受的，也是地震資料處理過程中常用的一種假設條件。本文方法與常規射線類層析方法的關系可類比于疊前時間偏移與疊前深度偏移的關系。

(2)關于速度縱向線性遞增假設。由于壓實的作用，地層的速度通常是隨深度逐漸增加的，特別是對于距地表幾百米以內的速度分布(近地表建模主要關注的區域)這種現象更為明顯[22]。近地表速度分布的線性遞增假設符合大部分地區地質規律。

(3)本文方法基于射線路徑圓弧狀假設，換一個角度來說，等價于利用圓弧狀射線路徑去擬合實際觀測旅行時。由于射線類層析分辨率的限制(不小于第一菲涅耳帶)[23]，這種假設與常規射線類層析反演方法同樣是非常接近的。圖12為利用圖6所示的反演結果對不同位置的三炮進行射線追蹤的結果，由射線路徑分布可見，除了速度結構特別復雜的區域外，大部分射線路徑與圓弧較為接近。

對于大部分地區，以上假設是成立的，也就是說本文所述方法可取得較為理想的效果。對于復雜地區，本文所述方法反演的精度受到一定限制，但由于其極高的計算效率，仍然具有較強實際應用意義：①本文方法的反演結果可以作為常規射線層析方法的初始模型，相對于簡單的梯度模型，使用精度更高的初始模型可減少迭代次數；②以極少的計算成本獲取了近地表速度分布的基本規律、不同地區的復雜程度，可對隨后的常規射線層析類方法等處理起到參考作用。

圖12 利用圖6所示反演模型計算的三炮的射線路徑

鑒于本文所述方法主要利用回轉波，因此在滿足回轉波假設的地區將獲得更佳的效果。可通過兩條途徑來判定是否滿足假設：①基于對該探區的地質認識，比如黃土塬等探區由于壓實作用，近地表速度橫向變化不大、隨深度逐漸增加[24]；②基于炮集上初至隨炮檢距的變化形態，若在消除地表高程的影響之后，初至波視速度隨炮檢距線性增加[25]，則可判定該地區適用于本文方法。

5 結束語

本文發展了一種利用回轉波進行快速近地表建模的方法，通過采用多基準面高程校正方法增強了對起伏地表的適應性，通過采用基于局部加權的穩定射線參數估計方法提高了抗噪性。相對于常規的射線類層析反演方法，在大部分情況下，該方法能夠以極小的計算成本提供了一個與之相當的近地表模型，具有極大的實際應用價值。

與常規的射線類層析反演方法類似，本文所述方法存在無法對速度反轉進行反演、復雜構造區域反演效果較差的問題。一種解決策略是將本文所述方法作為初始模型采用精度更高的反演方法(如波動方程旅行時反演、波形反演等)進行近地表建模；另一種可能的解決策略是采用新的反演思路，如利用深度神經網絡等機器學習算法，這也是下一步研究方向。