基于自適應震源子波提取與校正的瑞利波波形反演

邵廣周 獨 婷 吳 華

(①長安大學地質(zhì)工程與測繪學院，陜西西安 710054；②長安大學理學院，陜西西安 710064)

0 引言

在淺層介質(zhì)結(jié)構(gòu)探測中，瑞雷波的能量約占波場總能量的三分之二[1-2]，因而地震記錄中最突出的波形是直達瑞利波及其頻散波組。瑞利波是由P波和SV波干涉形成，其頻散波組中包含豐富的S波速度信息，利用瑞利波反演可以得到地下結(jié)構(gòu)的S波速度。全波形反演技術(shù)因具有精度高、可進行多參數(shù)建模的特點，越來越受到勘探地球物理學家的重視。然而，全波形反演存在以下不足：其大多數(shù)應用集中于體波，面波常被視為噪聲而去除；常規(guī)全波形反演通常缺乏低頻信息，多數(shù)研究以反演P波速度為主。而瑞利波的低頻特性恰好可以不受低頻信息缺失的限制，且瑞利波相速度與S波速度密切相關(guān)，因此，利用瑞利波波形反演淺層S波速度結(jié)構(gòu)較常規(guī)全波形反演具有一定優(yōu)勢。

早期瑞利波波形反演主要用于研究全球或區(qū)域大尺度的地球結(jié)構(gòu)[3-6]，由于研究中使用的面波波長大于淺層結(jié)構(gòu)的尺寸，使該方法在淺層分辨率方面存在不足。隨著研究的進展，出現(xiàn)了適用于淺層構(gòu)造成像的面波全波形反演方法。Zeng等[7]證明了瑞利波波形反演用于探測淺層速度結(jié)構(gòu)的可行性；Tran等[8]基于高斯—牛頓法在時域進行彈性全波形反演，成功重構(gòu)了20m范圍內(nèi)的淺層S波速度結(jié)構(gòu)；Groos等[9-10]研究了介質(zhì)黏彈性品質(zhì)因子對反演結(jié)果的影響，隨后實現(xiàn)了淺層瑞利波全波形反演，并通過模型試驗分析其穩(wěn)定性；Borisov等[11]利用基于包絡(luò)的目標函數(shù)反演淺層大尺度、不均勻結(jié)構(gòu)，證明了基于包絡(luò)的目標函數(shù)能有效減少周期跳躍問題。雖然瑞利波波形反演研究已經(jīng)取得了很大進展，但仍存在震源子波估計問題，需要進一步研究。

由于野外資料的震源子波未知，在反演中使用不準確的地震子波將不可避免地影響波形反演的結(jié)果[12]，因此，精確的子波估計是波形反演成功的關(guān)鍵。處理地震子波的方法主要有兩類：①在反演的迭代過程中進行震源子波校正濾波[9,13-16]，即在正演模擬時首先采用預先選取的地震子波(如Ricker子波)獲得合成數(shù)據(jù)(通常為共炮點道集)，然后利用最優(yōu)校正濾波器對合成數(shù)據(jù)進行震源子波校正濾波?；诰€性阻尼最小二乘理論，通過使觀測數(shù)據(jù)和波形校正濾波后的合成數(shù)據(jù)的傅里葉系數(shù)的殘差最小，在頻域構(gòu)造校正濾波器。雖然這種子波反演(或校正)方法可以獲得較好的反演結(jié)果，但反演過程需要多次迭代，增加了運行時間。Park等[17]采用阻尼函數(shù)對共中心點(CMP)道集進行子波校正，提高了計算效率。②將零炮檢距地震道作為經(jīng)驗子波或采用虛擬震源法，通過特征道褶積消除震源子波的影響[18-21]。

Schuster[22]提出一種反射波加窗提取地震子波的方法，可快速獲得近似震源子波，但窗口大小是人為確定的，窗口選取不合適會影響子波提取精度。為了更準確估計震源子波，本文采用自適應子波提取與子波校正聯(lián)合的方法。首先，對觀測數(shù)據(jù)進行自適應窗口時移疊加提取震源子波，子波提取時可自動確定窗口大??；其次，將提取的子波作為初始模型進行子波反演校正；再次，建立一組斷層模型并采用交錯網(wǎng)格有限差分法獲得其對應的合成記錄，進行波形反演試驗，采用不同的子波處理方法得到各自的反演結(jié)果，對這些反演結(jié)果進行對比分析以驗證各方法的應用效果；最后，通過一組實際數(shù)據(jù)進一步驗證聯(lián)合窗口自適應子波提取和子波校正波形反演方法的實際應用效果。

1 方法原理

1.1 瑞利波波形反演基本理論

全波形反演的基本原理是尋找一個能夠?qū)τ^測數(shù)據(jù)進行合理解釋的“最優(yōu)模型”，該模型能夠擬合特定地震波場(直達、反射、折射、面波等)的相位和振幅。

設(shè)炮點xs處激發(fā)、檢波點xr處接收到的t時刻地震波的第i個位移分量為ui(xs,xr,t)。如果已知地下模型參數(shù)分布，則可用有限差分方法求解波動方程，得到的模型正演數(shù)據(jù)記為umod。將理論模型數(shù)據(jù)與實際觀測數(shù)據(jù)uobs進行比較，得到數(shù)據(jù)殘差為δu=umod-uobs。如果δu很小，則表明該模型能夠很好地解釋觀測數(shù)據(jù)。波場的能量殘差E可用殘差向量的L2范數(shù)表示

(1)

瑞利波波形反演的基本步驟為：①給定一個初始模型，然后通過波場數(shù)值模擬計算初始模型或迭代修正模型的瑞利波場；②計算觀測波場與模擬波場之間的波場殘差即目標函數(shù)；③進行反演迭代，直到獲得滿足迭代精度要求的地下速度結(jié)構(gòu)成像。

在進行初始模型或迭代修正模型的波場數(shù)值模擬計算時，震源子波的估計或校正是波形反演中的一個關(guān)鍵步驟，需給定與觀測波場一致的震源子波，否則將會對波形反演結(jié)果帶來負面影響。

1.2 窗口自適應子波提取

通過自適應設(shè)定一個窗口，在窗口內(nèi)對瑞利波記錄進行震源子波提取(圖1)。具體步驟如下。

圖1 窗口自適應子波提取流程

(1)統(tǒng)計每一道峰值大小以及對應的時刻，根據(jù)最小時刻和最大時刻之間的時間寬度各向兩端擴大一定比例(如5%)，確定自適應窗口的時間長度T；再根據(jù)最大峰值所在的道號以及第一個峰值小于某個閾值(如最大峰值的20%)所在的道號確定自適應窗口的位置長度R(R為道數(shù))，可以自適應地獲得一個大小為T×R的局部窗口。

(2)在局部窗口內(nèi)通過第一道與其他道的互相關(guān)計算時間偏移，通過時間偏移對局部窗口中的所有道進行平移。

(3)對窗口內(nèi)所有道進行疊加以提高子波的信噪比。

(4)根據(jù)需要調(diào)整子波長度。

1.3 子波校正濾波

波形反演中觀測數(shù)據(jù)與合成數(shù)據(jù)的殘差主要由真實模型與當前計算模型的差異造成，也有一些是由于觀測數(shù)據(jù)對應的真實子波與合成記錄使用的假定子波之間存在差異造成的。為了消除子波帶來的影響，Groos等[9]提出一種校正濾波方法處理波形反演中的子波問題，采用如下基本思路。

首先對當前模型采用同一個假定子波得到各炮合成記錄，然后將合成記錄與觀測數(shù)據(jù)進行反褶積，得到子波校正濾波器，利用這些濾波器對前面正演步驟中使用的假定子波進行濾波，得到每一炮對應的最優(yōu)子波。在多尺度波形反演過程中，每個頻段只進行一次震源子波校正濾波。子波校正濾波器是在頻率域構(gòu)建的，利用線性阻尼最小二乘法使校正濾波后的觀測數(shù)據(jù)與合成數(shù)據(jù)的傅里葉系數(shù)的殘差最小，進而構(gòu)造一個最優(yōu)校正濾波器。

采用目標函數(shù)

(2)

(3)

式(2)中

(4)

表示用fk加權(quán)后的數(shù)據(jù)殘差。

(5)

表示對數(shù)據(jù)進行歸一化處理。其中

(6)

(7)

可得濾波器的傅里葉系數(shù)

(8)

式中“*”號表示復共軛。

2 模型試算

為了驗證子波提取與校正方法的有效性，建立如圖2所示斷層模型。模型尺寸為40m×12m，斷層上、下盤埋深分別為8、4m。上層介質(zhì)P、S波速度分別為800、200m/s，下層介質(zhì)P、S波速度分別為1200、400m/s。上、下層介質(zhì)的密度均為2000kg/m3。數(shù)值模擬采用交錯網(wǎng)格有限差分法。網(wǎng)格單元尺寸為0.2m×0.2m，時間步長為0.05ms。炮點位置分別為0、10、20、30、40m，炮間距為10m。每炮記錄36道接收，道間距為1m，記錄長度為0.6s，采樣率為0.5ms。每炮的炮點、接收點位置固定。

圖2 斷層模型S波速度分布圖

采用圖3所示的S波速度遞增模型作為初始模型，分別以真實震源子波、錯誤子波、錯誤子波聯(lián)和子波校正、子波提取得到的子波、提取子波聯(lián)和子波校正等5種子波測試瑞利波波形反演效果。為驗證本文方法不受真實子波的影響，對于前三種情況，模型對應的合成記錄均采用30Hz 正弦衰減子波。對于模型試算中的其他情況，模型對應的合成記錄則采用 30Hz 的Ricker子波。

圖3 反演初始模型

第一種情況，當實際資料的震源子波已知時，正、反演過程均采用同一個真實子波(30Hz正弦衰減子波)，反演結(jié)果如圖4所示。圖4a為經(jīng)過79次迭代后的斷層模型速度結(jié)構(gòu)反演結(jié)果?？梢钥闯?，斷層界面在地下4、8m處有清晰的速度分層，基本可以準確重建斷層模型。圖4b為真實模型、反演使用的初始模型和反演結(jié)果對應的第3炮記錄的波形對比圖(每隔2道顯示1道)。為便于分析，圖4c中只顯示了第10道。從圖4c可以看出，真實模型對應的波形與初始模型對應的波形差異較大，但與反演結(jié)果基本吻合。由此可見，當震源已知時，瑞利波波形反演結(jié)果是可靠的。

圖4 采用真實震源子波的斷層模型瑞利波波形反演結(jié)果

第二種情況，采用錯誤的地震子波進行波形反演，即真實子波為30Hz的正弦衰減子波，反演采用30Hz的Ricker子波。圖5為經(jīng)過92次迭代得到的反演結(jié)果。可以看出，斷層位置明顯錯誤，分層邊界不清晰，上、下層速度誤差較大。由此可見，使用錯誤震源子波得到的反演結(jié)果明顯不可靠。因此，為提高反演精度，有必要在反演之前先進行震源子波提取。

第三種情況，反演采用錯誤的震源子波，但同時進行了子波校正，亦是真實子波為30Hz的正弦衰減子波，反演采用30Hz的Ricker子波。圖6為經(jīng)90次迭代得到的反演結(jié)果。從反演結(jié)果可以識別出斷層構(gòu)造的基本輪廓。與未進行子波校正濾波的反演結(jié)果(圖5)相比，反演精度得到了一定程度的提高，但仍存在斷層位置錯誤、分層邊界不明顯、速度反演結(jié)果誤差大等不足。由此可見，在真實震源子波未知的情況下，基于震源子波校正的波形反演方法可以在一定程度上提高反演精度。然而，對于更為復雜的橫向非均勻模型，精度改進有限，無法滿足高精度反演的要求。

圖5 采用錯誤子波的斷層模型瑞利波波形反演結(jié)果

圖6 采用錯誤子波同時進行子波校正后的斷層模型瑞利波波形反演結(jié)果

第四種情況，真實子波為30Hz的Ricker子波，反演采用從不同炮記錄中提取的子波，提取結(jié)果如圖7所示，其中紅色曲線為真實的Ricker子波，其他顏色的曲線分別對應不同炮記錄的子波提取結(jié)果。可以看出，雖然每一炮提取的子波有所不同，但都與真實子波比較接近，可以直接用于波形反演。

圖7 從斷層模型模擬數(shù)據(jù)不同炮提取的子波與真實子波對比

圖8為僅使用提取的子波進行反演對應的斷層模型反演結(jié)果。圖8a為經(jīng)97次迭代后的S波速度反演結(jié)果，反演使用的初始模型與圖3相同。可以看出，斷層成像準確，速度分層較為清晰，上、下各層速度接近真實值。但在水平方向25～40m范圍內(nèi)，下層速度約為300m/s，與真實值400m/s有較大差異。這可能是提取子波與真實子波之間的誤差引起的。圖8b為真實模型、初始模型和反演結(jié)果對應的第3炮記錄的波形對比(每隔2道顯示1道)，圖8c為圖8b中的第28道。從圖中可以看出，真實模型對應的波形與初始模型相差較大，但與反演結(jié)果基本吻合。說明僅利用提取的子波進行波形反演，反演結(jié)果基本準確。

第五種情況，為了進一步提高反演精度，在子波提取的基礎(chǔ)上加入子波校正，進一步對地震子波進行約束。聯(lián)合使用提取的子波和子波校正進行反演，圖9為經(jīng)過89次迭代后的反演結(jié)果。從反演結(jié)果可以看出，斷層模型的成像效果較好，速度分層清晰(圖9a)。值得注意的是，8m以下的斷層上盤速度比僅用提取的子波反演的速度(圖8a)更準確(接近真實值400m/s)。雖然淺部速度誤差較大，且存在局部速度異常，但經(jīng)過子波校正后，反演精度明顯提高。圖9b為真實模型、初始模型和反演結(jié)果對應的第3炮記錄的波形對比(每隔2道顯示1道)，圖9c為圖9b中的第34道?？梢钥闯?，真實模型對應的波形與初始模型相差較大，但與反演結(jié)果基本吻合。因此，聯(lián)合使用提取的子波和子波校正可以進一步提高反演精度，特別是提高深層的速度反演精度。另外，反演迭代次數(shù)比僅使用提取的子波進行反演相比減少了8次，說明聯(lián)合使用子波提取與校正的方法可以節(jié)省一定的計算時間。

圖8 只采用提取子波的斷層模型瑞利波波形反演結(jié)果

圖9 采用提取子波聯(lián)合子波校正濾波后的斷層模型瑞利波波形反演結(jié)果

本文采用多尺度反演策略，在給定的頻帶內(nèi)對相

應的觀測數(shù)據(jù)和提取的子波從低頻到高頻逐帶進行反演。在反演過程中，每次使用一個新的帶寬更大的數(shù)據(jù)子集時，子波校正濾波器都會估計出一個最優(yōu)的震源子波用于波形反演，直到數(shù)據(jù)帶寬再次增加。圖10為第3炮提取的子波在相應頻帶寬度內(nèi)經(jīng)過子波校正后得到的地震子波?？梢钥闯觯Ｕ蟮淖硬S著頻帶寬度的變化而變化。當頻帶寬度為0～70Hz時，采用提取的子波聯(lián)合子波校正得到的子波(藍色虛線)與真實子波(紅色實線)及提取的子波(藍色實線)基本重合。圖11是每一炮提取的子波經(jīng)子波校正濾波后得到的最終子波。可以看出，最終的子波與真實子波相似。

圖10 不同頻帶范圍內(nèi)提取子波經(jīng)校正濾波后對應的子波與真實子波、提取子波對比

圖11 不同炮提取的子波經(jīng)過子波校正濾波后得到的最終子波與真實子波對比

由圖4～圖9不同子波處理的斷層模型反演測試結(jié)果表明，基于子波提取與校正的波形反演可有效提高反演精度，特別是可以提高深層速度結(jié)構(gòu)的反演精度。

3 實際資料測試

選取了一個人工隧道進行實際資料瑞利波探測試驗。Wang 等[23]和Smith 等[24]開展的類似研究證明，二維或三維彈性全波形反演方法在淺層隧道探測和成像方面具有廣闊的應用前景。本次試驗利用基于二維彈性全波形反演方法驗證震源子波對隧道探測效果的影響?，F(xiàn)場測試點如圖12所示，隧道寬約6m，頂板埋深約3m，底板埋深約6m。測線垂直于隧道走向，位于隧道上方的平臺上。從西向東0～47m之間放置48個間隔1m的低頻(4.5Hz)檢波器。炮點位于12、18、24、30、36m處，炮間距為6m。采用15磅錘擊震源，記錄長度為512ms，采樣率為0.5ms。隧道位于測線水平位置27～33m處。

圖12 試驗點概貌

反演計算區(qū)域為59m×25m，網(wǎng)格單元尺寸為0.25m×0.25m，采用0.05ms時間步長。波場計算采用了30Hz的正弦衰減子波，合成記錄長度為512ms，采樣率為0.5ms。每炮記錄共48道，道間距為1m，炮間距為6m。采用自由表面邊界和完全匹配層吸收邊界(厚度為5m)進行計算。除去吸收邊界厚度，與實際數(shù)據(jù)匹配的模型尺寸為49m×20m。

圖13為現(xiàn)場實際數(shù)據(jù)的5個單炮記錄，從圖13a和圖13b可以看出，由于隧道空腔的存在，第27～第33道之間(紅色圓圈為空洞位置)同相軸的能量及相位發(fā)生了明顯的變化。由于野外采集的面波數(shù)據(jù)具有較高的信噪比，對數(shù)據(jù)僅進行了5～70Hz的帶通濾波預處理。圖14為初始模型(圖14a)僅使用提取的子波(圖14b)和聯(lián)合使用子波提取與校正經(jīng)過20次迭代得到的反演結(jié)果(圖14c)。從圖14b、圖14c可以看出，在27～33m處(紅色圓圈為空洞實際位置)存在明顯的低速異

圖13 隧道探測試驗野外實際數(shù)據(jù)單炮記錄圖(Z分量)

圖14 實際數(shù)據(jù)波形反演采用的初始模型及反演結(jié)果

常，其位置與隧道實際位置幾乎一致。隧道頂板埋深可以很好地成像，但底部界面反演深度略深于實際埋深，總體成像效果令人滿意。圖15為第4炮的觀測波形與預測波形對比。圖16為相應的反演擬合誤差收斂曲線對比。對比可以看出，僅使用子波提取方法得到的隧道寬度反演誤差比聯(lián)合使用子波提取與校正方法的誤差大得多。實測數(shù)據(jù)反演結(jié)果表明，聯(lián)合子波提取與校正可以進一步提高瑞利波波形反演精度。

圖15 第4炮觀測數(shù)據(jù)與反演預測數(shù)據(jù)對比

圖16 歸一化誤差收斂曲線對比

4 結(jié)論

震源子波的選取是影響地震波場正演精度的重要因素之一，進而影響瑞利波全波形反演結(jié)果的精度。理論模型的測試結(jié)果表明如下認識。

(1)采用正確子波對應的反演結(jié)果最好，采用錯誤子波對應的反演結(jié)果最差。如果震源子波未知，直接進行瑞利波波形反演是不可靠的。

(2)窗口自適應子波提取方法可以從觀測記錄中提取子波，并將提取的子波直接用于波形反演，能夠獲得較好的反演效果。

(3)在子波提取方法基礎(chǔ)上，對提取的子波進行子波校正濾波，進一步約束震源子波。試驗結(jié)果表明，聯(lián)合子波提取與子波校正能夠有效地提高深部速度的反演精度，使反演速度更接近實際值。