基于水庫氣槍震源數(shù)據(jù)的頻域水準(zhǔn)比例因子反褶積方法分析

游秀珍 李軍 林彬華 黃艷丹 巫立華 郭陽

福建省地震局，福州市華鴻路7號 350003

0 引言

自20世紀(jì)60年代以來，氣槍震源在海洋石油勘探中得到廣泛應(yīng)用，積累了大量寶貴經(jīng)驗。隨著氣槍震源技術(shù)的日益成熟，為使人工氣槍震源能在地球物理探測領(lǐng)域也發(fā)揮作用，研究人員發(fā)展了更低頻率的震源——大容量氣槍震源。大容量氣槍震源以激發(fā)能量大、低頻成分豐富、可操控、可重復(fù)、成本低、破壞性小、綠色無污染等特點受到越來越多的關(guān)注。近些年來，國內(nèi)地震科學(xué)研究人員利用大容量氣槍震源開展了一系列深部結(jié)構(gòu)探測科學(xué)實驗，如福建陸域及中國臺灣海峽西部地區(qū)陸海聯(lián)合探測實驗和地學(xué)長江計劃實驗等，并在云南賓川、新疆呼圖壁和甘肅張掖建設(shè)了地震信號發(fā)射臺用以監(jiān)測地下介質(zhì)性質(zhì)的變化（陳颙等，2016；王寶善等，2016）。

水庫氣槍實驗中，由于水庫容積小，高壓氣泡在水庫中產(chǎn)生劇烈振動，同時出現(xiàn)長時間的振蕩子波，因而導(dǎo)致激發(fā)的氣槍震源信號不完全重復(fù)，進(jìn)而影響地下介質(zhì)變化的監(jiān)測精度。因格林函數(shù)可反映研究區(qū)域地下介質(zhì)的信息，地震學(xué)家們提出利用反褶積方法獲取氣槍震源至接收臺站之間的近似格林函數(shù)以消除震源的影響（王寶善等，2012；欒奕等，2016）。劉自鳳等（2015）、魏蕓蕓等（2016）在利用氣槍震源識別地殼介質(zhì)波速變化之前都先進(jìn)行反褶積處理以排除震源變化的干擾。王寶善等（2012）通過比較交叉干涉與互相關(guān)方法認(rèn)為，反褶積的結(jié)果雖然信噪比較低，但是能較好地去除激發(fā)環(huán)境的影響。翟秋實等（2016）對比了頻率域和時間域反褶積的效果，指出頻率域反褶積計算效率有優(yōu)勢，但信噪比不如時間域反褶積。

前人的相關(guān)工作側(cè)重利用反褶積直接處理數(shù)據(jù)，或是與其它方法的效果進(jìn)行比較，而對頻率域水準(zhǔn)反褶積方法計算存在的問題鮮有研究。本文采用頻率域水準(zhǔn)比例因子反褶積方法和互相關(guān)時延檢測技術(shù)，研究反褶積方法引入的水準(zhǔn)比例因子取不同值時對反褶積計算結(jié)果的影響，以及臺站的背景噪聲水平對計算結(jié)果的影響，并分析反褶積法消除震源影響的效果。

1 水庫氣槍實驗概況

2016年12月，福建省地震局在漳州南靖南一水庫開展主動源探測實驗。實驗使用福建省地震局自主設(shè)計組裝的大容量氣槍震源系統(tǒng)，選用4支Bolt-1500LL型氣槍組合陣列，其單支氣槍容量為2000 in3，槍陣總?cè)萘窟_(dá)8000 in3。槍陣尺寸為7m×7m，沉放深度為12m，每一工況實驗定點激發(fā)9次。實驗觀測系統(tǒng)由陸海臺陣構(gòu)成：陸域臺陣包括88個福建測震臺、福建周邊（廣東、江西、浙江、湖南和中國臺灣）部分測震臺站、40個實時傳輸?shù)牧鲃优_和100個PDS（非實時傳輸?shù)囊巴饬鲃优_）；海域臺陣由57套OBS形成的網(wǎng)格狀臺陣組成。本文選取資料時考慮到氣槍記錄信號的清晰度，選擇震中距為200m的岸邊測震臺L3583為近場參考臺，震中距為18.53km的測震臺NJJS和震中距為19.05km的測震臺NJNK為遠(yuǎn)場臺站（圖1）。

圖1 南一水庫氣槍激發(fā)位置及實驗觀測臺站分布

2 方法原理

2.1 頻率域水準(zhǔn)比例因子反褶積方法

地震波數(shù)據(jù)包含震源函數(shù)、震源與地震儀器之間的格林函數(shù)、儀器響應(yīng)和地面運動噪聲，可將其表示為

式中，u（t）為地震儀器記錄的數(shù)據(jù)；s（t）為震源時間函數(shù)；g（t）為格林函數(shù)；i（t）為儀器響應(yīng)；n（t）為噪聲。

我們把靠近氣槍震源的近場臺站信號當(dāng)作近似的震源時間函數(shù)，由于實驗期間觀測設(shè)備沒有發(fā)生變化，因此，不考慮儀器響應(yīng)的影響，同時忽略噪聲影響。將遠(yuǎn)場臺站和近場臺站記錄分別進(jìn)行傅里葉變換，再相除即為頻率域反褶積

式中，G（ω）為近似格林函數(shù)的頻譜；U（ω）為遠(yuǎn)場臺站氣槍記錄的頻譜；S（ω）為近場臺站氣槍記錄的頻譜。由于式（2）中分母可能出現(xiàn)為零或者趨于零的情況，以致頻率域的除法不穩(wěn)定。采用Helmberger等（1971）提出的頻率域因子反褶積方法，將式（2）改進(jìn)為

式中，S*（ω）為S（ω）的復(fù)共軛；max為求最大值運算符；c為水準(zhǔn)比例因子。

2.2 互相關(guān)時延檢測方法

互相關(guān)時延檢測方法是一項被廣泛應(yīng)用的成熟技術(shù)，其原理是以不同時延計算兩列相似波形的互相關(guān)系數(shù)，互相關(guān)系數(shù)最大時所對應(yīng)的時間延遲，即是兩列波形的延時?；ハ嚓P(guān)時延測量是以一個采樣間隔作為步進(jìn)的，然而相關(guān)函數(shù)的最大峰值通常并不落在采樣點上而是偏離采樣點，這可以通過擬合的方式解決以獲得更高的時間測量精度（Céspedes et al，1995）。為此，本文采用拋物線擬合方法重建相關(guān)函數(shù)，得到新相關(guān)函數(shù)的峰值和其對應(yīng)的位置（圖2）。擬合前后相關(guān)函數(shù)峰值對應(yīng)位置的時間延遲計算公式為

圖2 拋物線擬合曲線

式中，y1為擬合前的最大相關(guān)系數(shù)；y0和y2分別為擬合前相關(guān)函數(shù)峰值相鄰采樣點對應(yīng)的相關(guān)系數(shù)；ΔT為采樣時間間隔。

3 計算結(jié)果

以2016年12月南一水庫A09工況為例，先根據(jù)氣槍的激發(fā)時間將原始波形進(jìn)行分段截取，再將垂向分量記錄進(jìn)行去均值、2～8Hz濾波、歸一化等預(yù)處理，圖3是參考臺L3583和遠(yuǎn)場臺NJNK預(yù)處理后的波形。由于工況實驗持續(xù)時間不足半小時，其間未發(fā)生有影響的地震，且水庫水位基本穩(wěn)定，故可認(rèn)為地下介質(zhì)是不變的。在實驗過程中裝載氣槍的浮臺會隨氣槍激發(fā)發(fā)生偏移，因此，福建省地震局的工作人員采用懸掛重物、繩索固定等方式對浮臺進(jìn)行固定，根據(jù)浮臺上安裝的GPS監(jiān)測結(jié)果可知，激發(fā)位置的變化較小，基本控制在1m以內(nèi)。我們首先分析反褶積與疊加的先后順序?qū)Y(jié)果的影響。圖4為NJNK臺的計算結(jié)果，將不同順序計算的近似格林函數(shù)做互相關(guān)，得到互相關(guān)系數(shù)峰值為0.999，可見反褶積與疊加的先后順序?qū)ψ罱K結(jié)果的影響較小，因此，進(jìn)行數(shù)據(jù)處理時先疊加再反褶積計算。

圖3 L3583、NJNK臺預(yù)處理后的波形

圖4 不同順序計算的近似格林函數(shù)對比

3.1 水準(zhǔn)比例因子的影響

頻率域水準(zhǔn)反褶積方法引入水準(zhǔn)比例因子以獲得合理的水準(zhǔn)閾值，計算時使分母中頻譜振幅小于閾值的提升至閾值，從而提高反褶積計算結(jié)果的穩(wěn)定性。為探究水準(zhǔn)因子對反褶積結(jié)果的影響，在相同情況下，設(shè)定水準(zhǔn)比例因子為 0.1、0.01、0.001、0.0001、0.00001和0.000001共6個不同值，對反褶積后格林函數(shù)的振幅信噪比、格林函數(shù)相對于反褶積前波形的走時差及兩者的互相關(guān)系數(shù)峰值進(jìn)行比較。振幅信噪比采用有效信號幅值絕對值的最大值與噪聲振幅的均方根之比。

表1為A09工況NJNK臺不同水準(zhǔn)比例因子得到的格林函數(shù)振幅信噪比、走時差和互相關(guān)系數(shù)峰值，表2為 A13工況 NJJS、NJNK臺的計算結(jié)果。由表1、2可見，隨著水準(zhǔn)因子值減小，格林函數(shù)的振幅信噪比呈減小趨勢；格林函數(shù)相對于反褶積處理前波形的走時差有差別，不同的水準(zhǔn)比例因子可能引起數(shù)毫秒的走時變化；不同工況、不同臺站計算的互相關(guān)系數(shù)峰值變化情況也不同。圖5為A09工況NJNK臺不同水準(zhǔn)比例因子反褶積后的結(jié)果。由圖5可見，水準(zhǔn)比例因子越小，背景噪聲越大，小到一定值時有效信號則湮沒在背景噪聲中。因此，在使用頻率域水準(zhǔn)反褶積方法時，需要根據(jù)實際情況選取合適的水準(zhǔn)比例因子。由圖5還可見，反褶積后的波形相對反褶積前有一個較為明顯的走時變化，約為0.3s，其原因可能是，氣槍激發(fā)信號主要由氣槍氣體釋放所產(chǎn)生的壓力脈沖和由于氣泡在水中的振蕩所產(chǎn)生的氣泡脈沖構(gòu)成，壓力脈沖屬于高頻信號，主要頻率約為數(shù)十到數(shù)百赫茲，衰減較快，遠(yuǎn)場測震臺站幾乎無法記錄到壓力脈沖；而氣泡脈沖則屬于低頻信號，主要頻率為3～8Hz，遠(yuǎn)場測震臺站記錄到的主要是氣泡脈沖傳播引起的地面振動，根據(jù)水庫氣槍浮臺下方庫底OBS的記錄以及氣槍上方2m處水聽器的記錄，氣泡脈沖約在壓力脈沖后0.2s左右產(chǎn)生，同時，由于進(jìn)行本次氣槍實驗時，浮臺距岸邊的測震儀相對較遠(yuǎn)，約200m，因此，測震儀的記錄并不能完全等效于震源時間函數(shù)。根據(jù)浮臺與岸邊測震儀的記錄，及波在水中的傳播速度（1.5km／s，由于氣槍與岸邊測震臺站間距離僅200m，因此，我們認(rèn)為波在水中或地表較淺的上層傳播）可知，氣槍激發(fā)產(chǎn)生的壓力脈沖傳播至岸邊測震臺約需0.1s，氣泡脈沖則需0.3s左右，即本文中利用震中距為200m的岸邊測震記錄作為近似震源時間函數(shù)時，本身就引入了0.3s左右的時間誤差，這可能是本研究反褶積計算后波形走時出現(xiàn)0.3s變化的主要原因。這也說明，將距氣槍震源200m的臺站記錄作為震源時間函數(shù)可能并不完全適用。但由于本研究均使用同一臺站的觀測數(shù)據(jù)，因此，并不影響本文后續(xù)分析結(jié)論。

表1 A 09工況NJNK臺不同水準(zhǔn)比例因子格林函數(shù)振幅信噪比、走時差和互相關(guān)系數(shù)峰值比較

表2 A 13工況NJJS與NJNK臺不同水準(zhǔn)比例因子格林函數(shù)振幅信噪比、走時差和互相關(guān)系數(shù)峰值比較

圖5 A09工況NJNK臺不同水準(zhǔn)比例因子反褶積后的結(jié)果

3.2 背景噪聲水平對結(jié)果的影響

在使用頻率域水準(zhǔn)反褶積方法監(jiān)測地殼介質(zhì)變化前，有必要了解反褶積計算結(jié)果的影響因素。本文采用仿真數(shù)值模擬方法，通過改變遠(yuǎn)場臺或參考臺的背景噪聲水平以及不同水準(zhǔn)比例因子，探討背景噪聲對反褶積計算的影響。

以往研究中，多采用隨機白噪聲數(shù)值來模擬噪聲的變化，而實際上，地震儀記錄的天然噪聲并不是絕對的隨機白噪聲，而是隨時間、空間具有一定變化規(guī)律的噪聲，因此，本文選取氣槍未激發(fā)時段連續(xù)半個月測震臺站記錄的噪聲數(shù)據(jù)，以仿真噪聲變化對分析結(jié)果的影響。對所選噪聲數(shù)據(jù)每隔10min抽取一段，并作基線校正、2～8Hz濾波和兩頭尖滅處理，再計算每段噪聲的標(biāo)準(zhǔn)偏差，結(jié)果如圖6所示，由圖6可見，噪聲呈規(guī)律性日變化。

圖6 噪聲日變化

計算步驟如下：①把遠(yuǎn)場臺記錄數(shù)據(jù)預(yù)處理后的噪聲部分置換為零，再將抽取的噪聲加到該臺的氣槍記錄上，將生成的新信號當(dāng)作遠(yuǎn)場臺記錄；②把新的遠(yuǎn)場臺記錄與參考臺記錄做反褶積得到近似格林函數(shù)（根據(jù)前文所述，水準(zhǔn)比例因子取0.001）；③由于反褶積會引入高頻噪聲，因此需要再次濾波，頻帶范圍為2.5～6Hz；④采用互相關(guān)時延方法，計算根據(jù)不同水準(zhǔn)比例因子得到的近似格林函數(shù)P波段（0.5s窗長）走時誤差隨噪聲的日變化；⑤改變遠(yuǎn)場臺記錄的信噪比，即改變抽取噪聲的幅值，計算格林函數(shù)走時誤差隨噪聲的日變化；⑥對參考臺按照同樣的流程進(jìn)行處理。

圖7、8分別為A09工況遠(yuǎn)場臺NJNK和參考臺L3583加噪聲反褶積后格林函數(shù)P波段走時識別誤差隨噪聲的日變化。圖7（a）為遠(yuǎn)場臺氣槍記錄的信噪比約為49時不同水準(zhǔn)比例因子的計算結(jié)果，圖8（a）為參考臺記錄的信噪比約為50時不同水準(zhǔn)比例因子的結(jié)果，由圖7（a）、8（a）可見，當(dāng)因子相同，且遠(yuǎn)場臺、參考臺記錄的信噪比相差不多時，前者計算的走時誤差只有數(shù)毫秒，而后者識別的走時誤差明顯大于前者，可見參考臺記錄的信噪比對結(jié)果的影響遠(yuǎn)大于遠(yuǎn)場臺；另外，水準(zhǔn)比例因子的不同對反褶積計算后波形的走時有影響，當(dāng)水準(zhǔn)比例因子取值過小時，走時誤差較大。本研究只是改變臺站的背景噪聲水平，并未更改氣槍信號，因此，不同信噪比的氣槍記錄反褶積計算后波形的走時變化理論上應(yīng)該相同，但結(jié)果顯示，由不同的信噪比計算得到的走時變化并不相同（圖7（b）、8（b）），這說明臺站的背景噪聲水平對結(jié)果有較大影響，氣槍記錄的信噪比越大，所得結(jié)果精度越高。遠(yuǎn)場臺不同信噪比的氣槍記錄所得走時誤差變化趨勢一致，說明識別誤差均由噪聲引起，這就要求在正式的應(yīng)用中需要花費較大精力去排除噪聲影響。當(dāng)遠(yuǎn)場臺記錄的信噪比大于10時，走時誤差一般小于6ms。對參考臺具有不同信噪比的記錄，由于信噪比不同會影響反褶積中的近似震源函數(shù)，因此得到的走時誤差變化趨勢一致性較差。此外，參考臺記錄的信噪比對走時識別結(jié)果影響更大，信噪比為30左右的參考臺記錄，走時誤差可能達(dá)到20ms左右，即使對于信噪比100左右的參考臺記錄，走時誤差也可能達(dá)到幾毫秒。此外，不同水準(zhǔn)比例因子及臺站記錄的信噪比計算的走時誤差日變化與背景噪聲日變化間沒有明顯關(guān)聯(lián)。

圖7 遠(yuǎn)場臺加噪聲反褶積后P波段走時誤差隨噪聲的日變化

圖8 參考臺加噪聲反褶積后P波段走時誤差隨噪聲的日變化

3.3 消除震源影響的效果

為了檢驗利用頻率域水準(zhǔn)比例因子反褶積方法消除氣槍震源影響的效果，我們對工作壓力為1000、1200、1500、1800、2000psi時NJNK測震臺記錄進(jìn)行反褶積計算（水準(zhǔn)比例因子取0.001），以槍壓為1000psi的A09工況作為參考，比較不同槍壓激發(fā)的觀測記錄反褶積前后波形的時間延遲。

圖9（a）、9（b）分別為不同工作壓力的氣槍激發(fā)信號反褶積前后的波形。由于氣槍震源信號重復(fù)性較高，因此以豎線為參照，比較各波形初始震相波段內(nèi)峰值的到時位置。由圖9可見，隨著槍壓的增加，反褶積前波形峰值位置的到時存在一定差異，且槍壓為2000psi時較明顯，而經(jīng)反褶積計算后，波形峰值位置的到時差異變小。設(shè)時間窗長為 0.3s，步長為0.2s，

圖9 不同工作壓力的氣槍信號反褶積前后波形對比

圖10 NJNK臺不同槍壓氣槍信號相對A09反褶積前后波形的時間延遲

計算不同槍壓的激發(fā)信號在反褶積前后分別與A09工況不同波段間的互相關(guān)時間延遲（圖10）。A9、A10工況的槍壓差為200psi，反褶積前兩列波形不同波段的時間延遲較小，反褶積后時間延遲稍大，但都小于0.01s（即1個采樣點）；反褶積前后初始震相波段計算的應(yīng)為第2個時間窗，該波段算得的時間延遲為時間延遲圖中第2個圓點，反褶積后的時間延遲基本比反褶積前小，特別是較A09工況而言，A13工況反褶積前時間延遲達(dá)2個采樣點，經(jīng)反褶積后時間延遲顯著縮短。由此可見，利用頻率域水準(zhǔn)反褶積方法，可以有效消除震源差別較大的影響。

4 結(jié)論與討論

本文以2016年福建南一水庫氣槍實驗資料為研究對象，通過分析水準(zhǔn)比例因子及臺站背景噪聲對反褶積結(jié)果的影響，研究利用反褶積法消除震源影響的效果，得出以下結(jié)論與認(rèn)識。

（1）水準(zhǔn)比例因子影響反褶積后波形的信噪比及走時，水準(zhǔn)比例因子越小，反褶積后波形的信噪比越小，當(dāng)水準(zhǔn)比例因子取值過小時，走時識別誤差較大。實際應(yīng)用時，需進(jìn)行分析后選取合適的水準(zhǔn)比例因子。

（2）臺站的背景噪聲對反褶積計算結(jié)果有影響。氣槍記錄的信噪比越大，走時識別精度越高，參考臺記錄的信噪比對結(jié)果的影響遠(yuǎn)大于遠(yuǎn)場臺；當(dāng)遠(yuǎn)場臺記錄的信噪比大于10時，走時誤差一般小于6ms；當(dāng)參考臺記錄的信噪比為30左右時，走時誤差可能達(dá)到20ms左右，即使對于信噪比100左右的參考臺數(shù)據(jù)，走時誤差也可能達(dá)幾毫秒。

（3）利用頻域水準(zhǔn)比例因子反褶積方法，可以消除不同工作壓力引起的震源影響，驗證了該方法去除震源效應(yīng)的有效性。該方法在震源差異較大時去除震源效應(yīng)的效果較好，在震源差異較小時，反而可能引入方法誤差。

本文對頻域水準(zhǔn)反褶積計算結(jié)果的影響因素及反褶積法消除震源影響的效果進(jìn)行了初步分析，研究結(jié)果可為反褶積方法應(yīng)用于氣槍監(jiān)測地殼介質(zhì)變化中的應(yīng)用提供參考。