畢麗飛，曾志毅，張建中，黃忠來，芮擁軍，刁 瑞

（1.中國石油化工股份有限公司勝利油田分公司科技處，山東東營257001；2.海洋科學(xué)與探測技術(shù)教育部重點實驗室，中國海洋大學(xué)海洋地球科學(xué)學(xué)院，山東青島266100；3.中國石油化工股份有限公司勝利油田分公司物探研究院，山東東營257001）

水力壓裂微震監(jiān)測技術(shù)是一項通過觀測注水壓裂過程中所誘發(fā)的微震事件來監(jiān)測和評估裂縫發(fā)育情況的技術(shù)[1-4]。微震事件識別和到時拾取是微地震數(shù)據(jù)處理的關(guān)鍵環(huán)節(jié)，直接影響對裂隙成像的可靠性[5]。通常情況下，水力壓裂產(chǎn)生的微震事件能量十分微弱，微震觀測資料的信噪比較低。低信噪比資料的微震事件識別及其到時拾取是微震監(jiān)測技術(shù)的重要研究課題。

微震監(jiān)測記錄的數(shù)據(jù)量巨大，人工識別并拾取微震事件初至到時的工作量極大，且容易引入人為誤差。微震事件到時拾取是根據(jù)有效信號與背景噪聲在能量、偏振、頻譜和統(tǒng)計等方面的特性差異來實現(xiàn)的。常見的初至到時拾取方法基于單道記錄的算法，其中以長短時窗能量均值比（STA/LTA）算法[6-9]及其改進方法為代表[10-12]。該類方法計算速度快，可以滿足實時處理要求，但是對低信噪比微震資料的應(yīng)用效果不佳。近年來，一些學(xué)者綜合利用微震信號與環(huán)境噪聲的多種特征差異，提高了識別和拾取效果。如：宋維琪等[13]利用小波分解與Akaike信息準(zhǔn)則相結(jié)合，有效壓制噪聲引起的局部效應(yīng)，提高微震事件初至拾取的準(zhǔn)確性；吳治濤等[14]聯(lián)合小波變換與偏振分析實現(xiàn)微震事件的全自動識別和拾取P 波初至；譚玉陽等[15]通過綜合地震信號與環(huán)境噪聲在能量、偏振和統(tǒng)計方面存在的差異，提出一種微震事件初至拾取算法，有效提高了抗噪能力和拾取精度。

GIBBONS等[16]的研究表明，具有相似震源位置及破裂機理的地震事件通常在記錄上表現(xiàn)出相似的波形特征。一些學(xué)者發(fā)展了基于事件多道記錄相似特征的初至拾取算法[17-19]，如：譚玉陽等[20]通過多道相似系數(shù)實現(xiàn)微震事件的識別，避免了由于某些道數(shù)據(jù)的信噪比過低而遺漏微震事件的情況；VANDECAR 等[21]利用多道互相關(guān)與最小二乘法相結(jié)合，確定遠震震相的相對到時，提升了拾取效率和抗噪性；DE MEERSMAN 等[22]和AKRAM 等[23]分別基于互相關(guān)的迭代算法優(yōu)化了微震事件的初至拾取結(jié)果；魏夢祎等[24]提出一種基于波形互相關(guān)的微震事件自動拾取方法，通過疊加獲取高信噪比參考道，提高了初至拾取精度。在拾取微震事件的過程中，以上方法需要對各道微震事件記錄進行疊加得到信噪比更高的疊加參考道，但微震事件信號的初至極性與震源機制相關(guān)，各道微震事件的初至極性會不一致，即初至極性具有正負(fù)之分[25-28]，致使通過波形直接疊加得到的疊加參考道的信噪比反而降低，影響了對微震事件的識別和到時拾取精度。

本文提出了一種適于初至極性反轉(zhuǎn)的微震到時拾取方法。通過計算任意兩道之間的延時互相關(guān)，獲得各道微地震的相對到時；利用相鄰道相乘后再疊加得到參考道，避免了道間波形相似系數(shù)（Semblance）和直接波形疊加得到參考道時無法有效克服極性反轉(zhuǎn)影響的困難。本文首先介紹了微震事件識別和初至拾取的方法原理和技術(shù)流程，然后分別利用合成數(shù)據(jù)和實際資料對方法的性能進行了測試和分析。

1 方法原理

1.1 相對走時計算

取一時窗內(nèi)的微震記錄，利用互相關(guān)函數(shù)計算任意兩道微震記錄之間的互相關(guān)系數(shù)序列。時窗內(nèi)應(yīng)包含同一微震事件在各道記錄的P 波或S波初至，時窗長度應(yīng)大于微地震波的時長以及最短和最長偏移距檢波器之間的到時差，可根據(jù)射孔信號時長及最短與最長偏移距檢波器的到時差長度進行估計?；ハ嚓P(guān)系數(shù)計算公式為：

式中：xi和xj分別為地震道i和j的微震記錄；W為時窗長度；k為xi和xj的時間間隔；ci，j（k）為互相關(guān)系數(shù)?？紤]到初至極性存在反轉(zhuǎn)的情況，對互相關(guān)系數(shù)ci，j（k）取絕對值，其最大值cmaxij對應(yīng)的k為xi和xj的到時差，記為Δti，j。兩個地震道i和j包含的同一微震事件的到時差Δti，j可表示為：

式中：ti和tj分別為地震道i和j的微震初至到時。對所有地震道記錄進行兩兩互相關(guān)處理，可得到兩兩記錄中微震事件的到時差。以4個地震道（i＝1，2，3，4）為例，可得到6個到時差Δt1，2，Δt1，3，…，Δt3，4。為了避免方程組求解過程的不適定性，增加一個約束方程∑ti＝0，使相對到時的均值為0。這樣有下列方程：

若參與計算的微震記錄有N道，則需利用（1）式進行N（N－1）/2次互相關(guān)函數(shù)計算，得到N（N－1）/2個兩道間的到時差數(shù)據(jù)。（3）式寫成矩陣形式為：

式中：A為0，1和－1組成的稀疏系數(shù)矩陣；Δt為到時差向量；t為到時向量，也是待求向量。（4）式的最小二乘解為：

利用（5）式求出t后，再用t中的各個分量對各道記錄進行時差校正。時差校正后，時窗內(nèi)各道記錄的微震初至波在時間上將趨于一致。

1.2 參考道獲取

利用（5）式求出的是相對于時差校正后各道記錄平均到時的相對到時。求微震事件在各道記錄上的絕對到時，還需要求出各道的平均到時，即微震事件在疊加道或參考道上的初至到時。求參考道及其微震事件的絕對到時的常用做法是：先將時窗內(nèi)時差校正后的各道記錄波形直接進行疊加，形成高信噪比的參考道，再在參考道上拾取絕對到時。但是，當(dāng)時窗內(nèi)微震事件的初至極性存在反轉(zhuǎn)時，對校正后的各道微震記錄進行直接疊加，會使微震記錄振幅相互抵消，反而降低了參考道記錄的信噪比，從而影響微震事件絕對到時的拾取精度。有的方法通過對各道的初至極性進行識別或通過震源機制反演，將初至極性校正為相同極性后，再疊加得到參考道[23]，這樣做將會增加計算量。為了解決這個問題，我們利用相鄰道零延遲相乘再疊加的方法獲取參考道（圖1），計算式為：i＝1

式中：N是相鄰地震道數(shù)的配對數(shù)；W是時窗長度；w是時窗內(nèi)采樣點序號；x′i（w）表示時差校正后的地震道i在采樣點w處的幅值；X（w）為獲得的參考道在采樣點w處的幅值。圖1a中，時差校正后微震初至被“拉平”，其中左邊記錄（黑色地震道）與右邊記錄（紅色地震道）的初至極性相反。對時差校正后的微震記錄進行相鄰道的零延遲相乘，初至極性相同（同為正或負(fù)極性）的相鄰道記錄的零延遲相乘后的幅值均為正值；而初至極性相反（一正一負(fù)）的情況僅出現(xiàn)在不同初至極性區(qū)域的分界處，如圖1b中第8道（紅色）所示。幅值為負(fù)的道的數(shù)量很少，對整個相乘道的疊加值影響很小，而隨機噪聲的極性無一定規(guī)律，疊加后振幅減弱，從而得到信噪比更高的參考道，如圖1c所示。

圖1 相鄰道零延遲相乘再疊加得到參考道示意a存在初至極性相反的地震道；b 相鄰道零延遲相乘道；c參考道

1.3 絕對到時拾取

由相鄰道零延遲相乘再疊加得到的參考道，雖然沒有完整的微震信號波形，但是在微震信號初至?xí)r刻附近有較大的振幅或能量變化，因此，本文利用STA/LTA 算法判斷參考道是否含有有效微震信號以及確定其初至到時。STA/LTA 算法的計算式為：

式中：X（w）為參考道記錄；WL和WS分別為長、短時窗的長度；m為采樣時刻；STA 和LTA 分別為短時窗和長時窗內(nèi)的能量均值。在初至到時附近，STA 比LTA 變化更快，得到時窗內(nèi)的STA/LTA會出現(xiàn)明顯的極大值，一般近似地將STA/LTA 曲線的極大值點作為初至到時點。若時窗內(nèi)的STA/LTA 的極大值超過所設(shè)閾值[9]，則認(rèn)為該時窗內(nèi)含微震事件，并且同時將STA/LTA 曲線的極大值點所對應(yīng)的時刻取為該微震事件在參考道的初至到時，記為T0；否則認(rèn)為該時窗內(nèi)不含微震事件。本文以各個時窗的長短窗能量比值的均值為依據(jù)，將該均值的某個倍數(shù)（以下稱為長短窗能量比值均值倍數(shù)）作為微震事件的識別閾值，這樣，就可在每個時窗內(nèi)自適應(yīng)地設(shè)置閾值，從而避免固定閾值所帶來的識別錯誤和遺漏。

獲得參考道的初至到時T0，再結(jié)合（5）式中各道的相對到時ti，則各道微震事件的初至絕對到時Ti可由下式得到：

1.4 算法流程

本文的微震事件初至拾取算法流程如圖2所示，關(guān)鍵步驟概括如下。

1）確定滑動時窗大小及滑動步長，設(shè)置用于判別微震時間的長短窗能量比值均值的倍數(shù)b。

2）提取時窗內(nèi)的微震記錄，計算該時窗內(nèi)的道間互相關(guān)系數(shù)，確定微震事件的到時差Δti，j，并利用（5）式計算各道微震事件的相對到時ti。

3）利用得到的相對到時，對各道微震信號進行時差校正，對時差校正后的微震記錄進行相鄰道零延遲相乘再疊加處理，得到參考道。

4）利用STA/LTA 算法得到參考道的長短時窗能量比曲線，并計算能量比均值與設(shè)置的倍數(shù)b的乘積，即微震事件識別閾值R。判斷該時窗能量比極大值F是否大于微震事件識別閾值R，若是，進行步驟5）；否則返回步驟2），進行下個時窗的處理。

5）拾取STA/LTA 能量比曲線的極大值對應(yīng)時刻T0，再利用（8）式求取時窗內(nèi)各道微震事件初至的絕對到時Ti。

圖2 微震事件初至到時拾取算法流程

2 合成數(shù)據(jù)測試

首先利用合成數(shù)據(jù)來檢驗本文方法的應(yīng)用效果。圖3 為合成微震記錄，時間長度為3s，采樣率為1 ms。圖3a為無噪聲微震記錄，含有3個有效微震事件，從第1個到第3個微震事件的振幅不斷減小，其中第1個微震事件初至極性在第8～16道上發(fā)生反轉(zhuǎn)。圖3b為添加均值為0，標(biāo)準(zhǔn)差為0.4的高斯白噪的微震記錄，其中第3個微震事件信號已經(jīng)完全淹沒在噪聲中。分別采用常規(guī)方法和本文方法進行測試，選用的滑動時窗長度和滑動步長分別為120個采樣點和15個采樣點，取各個時窗參考道的長短窗能量比均值的3.5倍，作為微震事件的識別閾值。

圖4a和圖4b分別為采用各道記錄波形線性疊加參考道和本文的相鄰道零延遲相乘疊加參考道得到的STA/LTA 能量比曲線（藍色曲線）及自適應(yīng)閾值曲線（紅色曲線）?？梢钥闯?，圖4a所示的微地震記錄直接疊加的STA/LTA 能量比曲線只有1個明顯的峰值（紅色箭頭所指），指示了第2個微震事件，而對能量最強但存在初至極性反轉(zhuǎn)的第1個微震事件和信噪比較低的第3個微震事件卻沒反映出來；圖4b所示的曲線有3 個明顯的峰值（紅色箭頭所指），3個微震事件都被識別出來。可見，本文方法不但適于初至極性存在反轉(zhuǎn)的情況，對低信噪比資料也有很強的適應(yīng)性。

圖3 合成微震記錄

圖4 采用不同方法識別參考道微震事件的STA/LTA 能量比曲線（藍色）及自適應(yīng)閾值曲線（紅色）

下面用本文方法對微震事件初至到時進行拾取，關(guān)鍵步驟的處理結(jié)果如圖5 所示。圖5a 為根據(jù)圖3b中3個微震信號所在的時窗分別提取出來的經(jīng)時差校正后的含微震時間記錄，與理論微震事件的到時一致，說明微震事件相對到時計算的準(zhǔn)確性。圖5b是相鄰道零延遲相乘結(jié)果，在一定程度上壓制了不相干的背景噪聲。圖5c是對圖5b各道進行疊加的結(jié)果，即參考道，具有明顯的振幅極值，消除了初至極性反轉(zhuǎn)對參考道的不利影響。圖5d 為利用STA/LTA 算法得到的參考道的STA/LTA 能量比曲線及其極大值對應(yīng)的參考道的初至到時T0（紅色虛線所示）。

圖5 用本文方法對微震事件初至到時拾取過程的相關(guān)結(jié)果（左、中、右分別對應(yīng)圖3中的第1、第2和第3個微震事件）

利用得到的各道相對到時ti和參考道初至到時T0并用（8）式計算出微震事件在各道的初至到時Ti。圖6為采用本文方法拾取微震事件初至到時的誤差（拾取到時與理論到時之差）。從圖6可以看出，整體拾取誤差較小，說明采用本文方法可以得到較為準(zhǔn)確的微震事件初至的絕對到時。同時我們采用基于波形線性疊加參考道的識別方法拾取第2個微震事件初至到時，其相關(guān)結(jié)果和拾取誤差如圖7所示，與采用本文方法得到的第2個微震事件初至到時拾取結(jié)果僅相差1 ms（1個采樣點長度），說明本文方法的拾取精度與常用拾取方法的拾取精度相當(dāng)。但是，采用常用方法只能識別出不存在極性反轉(zhuǎn)的第2個微震事件，未能識別出含有極性反轉(zhuǎn)的第1個微震事件和信噪比很低的第3個微震事件，而采用本文方法可以識別出這3個微震事件?？梢?，采用本文方法不僅可以較為準(zhǔn)確地識別和拾取低信噪比微震事件，且能較好地克服初至極性反轉(zhuǎn)所帶來的不利影響，具有更強的適用性。

圖6 采用本文方法拾取微震事件初至到時的誤差（拾取到時減去理論到時）

圖7 采用基于波形線性疊加參考道的識別方法拾取第2個微震事件到時的相關(guān)結(jié)果及拾取誤差

3 實際資料處理

圖8 實際野外檢波器布設(shè)位置

采用本文方法對某區(qū)實際地面微震監(jiān)測資料進行處理和分析。野外布設(shè)42個檢波器，位置如圖8所示。圖9為一段經(jīng)過帶通濾波后的地面微震監(jiān)測記錄，時間長度為4 s，時間采樣間隔為1 ms。采用本文方法和常用方法處理該段實際資料時選取的參數(shù)為：滑動時窗為200個采樣點，滑動步長為15個采樣點，判別微震事件的閾值設(shè)置為該時窗內(nèi)長短窗能量比均值的3.5倍。圖10a和圖10b分別為采用微地震記錄波形線性疊加方法和本文方法得到的參考道STA/LTA 能量比曲線（藍色）及自適應(yīng)閾值曲線（紅色），可以看出，微地震記錄波形線性疊加參考道的STA/LTA 能量比曲線沒有出現(xiàn)明顯的極大值，而本文方法的STA/LTA 能量比曲線出現(xiàn)了2個超過識別閾值的明顯的峰值（圖10b中紅色箭頭所示），檢測出該段微震記錄中的2個微震事件。

圖9 野外實際微震監(jiān)測記錄

利用本文方法對該段微震記錄的2個微震事件的初至到時進行了拾取。圖11為實際微震信號的拾取過程的相關(guān)結(jié)果。圖11a為時差校正后的微震信號記錄，可以看出，微震事件波形基本上被校正拉平，其中第1個微震事件中的第5～12道和第35～40道的微震初至極性發(fā)生了反轉(zhuǎn)，波形線性疊加得到的參考道振幅會相互抵消，這是利用常規(guī)波形線性疊加參考道方法無法識別出微震事件的主要原因。第2個微震事件記錄的信噪比較低，個別道初至同相軸未被完全拉平，也導(dǎo)致波形線性疊加得到的參考道振幅較小，漏識了該微震事件。本文方法利用相鄰道零延遲相乘再疊加得到參考道，突出了被拉平的各道初至同相軸的作用，受部分未拉平道的影響相對較小。圖11b為利用相鄰道零延遲相乘的結(jié)果，有效壓制了背景噪聲，提高了信噪比。圖11c為由相鄰道零延遲相乘再疊加形成的參考道。圖11d為利用STA/LTA 算法對參考道進行處理得到的STA/LTA 能量比曲線，以及由該曲線的極大值點得到參考道的初至到時T0（圖中紅色虛線所示）。最后把各道的相對到時ti和參考道的絕對到時T0相加，得到了微震事件在各道的初至絕對到時Ti，如圖12所示（紅點為各道微震事件的初至到時時刻）。

圖10 采用不同方法得到的參考道STA/LTA 能量比曲線（藍色）及自適應(yīng)閾值曲線（紅色）

圖11 實際微震事件初至到時拾取過程的相關(guān)結(jié)果（左列對應(yīng)915～1 115 ms時窗記錄；右列對應(yīng)2 565～2 765 ms時窗記錄）

為進一步說明本文方法對實際資料中微震事件的識別率及初至拾取精度，以震源掃描定位算法[28]得到微震事件數(shù)目作為標(biāo)準(zhǔn)，分別采用本文方法和常用的基于波形線性疊加的微震事件識別方法對同一段壓裂的微地震監(jiān)測數(shù)據(jù)進行識別和拾取，結(jié)果如表1所示。從表1可以看出，采用本文方法的誤拾率和漏拾率都小于采用波形線性疊加參考道方法的結(jié)果，尤其是漏拾率遠小于后者，說明采用本文方法能識別出更多的有效事件，具有更高的微震事件識別能力。

圖12 采用本文方法對某段監(jiān)測記錄的微震事件到時的拾取結(jié)果（紅點為初至到時時刻）

由于采用掃描定位方法確定的初至到時受到速度模型精度的影響，因此本文在采用掃描定位法確定到時的基礎(chǔ)上進行人工拾取，再以人工拾取結(jié)果為準(zhǔn)確值，計算本文方法和波形線性疊加參考道方法拾取微震初至到時的誤差。波形線性疊加參考道方法共識別出26個微震事件，這里分別計算了采用本文方法和波形線性疊加參考道法拾取這26個微震事件的到時誤差，統(tǒng)計的拾取誤差的頻率如圖13所示。從圖13可以看出，采用波形線性疊加參考道方法的拾取誤差主要分布在2～21 ms，而采用本文方法的拾取誤差主要集中在2～12 ms?？傮w上，采用本文方法的拾取誤差要低于采用波形線性疊加參考道方法的誤差。

表1 采用不同方法識別微震事件的結(jié)果

圖13 采用不同方法拾取初至到時的誤差率

4 結(jié)論

常用的多道相似系數(shù)法由于在識別有效微震事件時受初至極性反轉(zhuǎn)的影響，因而會導(dǎo)致微震事件的遺漏。本文根據(jù)微震事件的波形相似特征，并用相鄰道零延遲相乘再疊加的方法取代常用的波形線性疊加方法，獲取更高信噪比的參考道，提出了一種從多道微震監(jiān)測記錄中識別和拾取微震到時的方法。該方法不但能夠壓制不相干噪聲，從低信噪比監(jiān)測資料中獲取有效微震事件信息，而且能夠克服初至極性反轉(zhuǎn)的不利影響。合成數(shù)據(jù)和實際資料應(yīng)用結(jié)果表明，本文方法具有較強的微震事件識別能力和較高的微震到時的拾取精度，對低信噪比資料和存在極性反轉(zhuǎn)資料有較強的適應(yīng)性。

在判別參考道各時窗是否存在微震事件時，以各個時窗的長短窗能量比均值為依據(jù)，將該均值的某個倍數(shù)作為判別的閾值，這樣自適應(yīng)地設(shè)置閾值，避免了固定閾值帶來的識別誤差。

采用本文方法能較好地區(qū)別具有波形相似性的微震記錄與不規(guī)則噪聲干擾。但在實際資料中，除背景噪聲外，還可能包括一些相干噪聲，如反射波和折射波及續(xù)至波等，這些噪聲會對互相關(guān)時差估計和參考道產(chǎn)生影響，導(dǎo)致對微震事件識別困難及對到時拾取產(chǎn)生誤差，需要通過分析壓裂區(qū)域微震事件的初至曲線特征來進一步研判和去除錯誤的拾取結(jié)果。