低信噪比微震P波震相初至自動拾取方法

賈瑞生,譚云亮,孫紅梅,洪永發

(1.山東科技大學礦山災害預防控制省部共建國家重點實驗室培育基地,山東青島 266590;2.山東科技大學信息科學與工程學院,山東青島266590;3.山東科技大學礦業與安全工程學院,山東青島 266590)

?

低信噪比微震P波震相初至自動拾取方法

賈瑞生1,2,譚云亮1,3,孫紅梅2,洪永發2

(1.山東科技大學礦山災害預防控制省部共建國家重點實驗室培育基地,山東青島 266590;2.山東科技大學信息科學與工程學院,山東青島
266590;3.山東科技大學礦業與安全工程學院,山東青島 266590)

摘 要:針對微震信號固有的低信噪比、非平穩性、隨機性等特征,發展了一種基于希爾伯特-黃變換(HHT)和AIC準則相結合的震相初至自動拾取方法,該方法首先對含噪信號進行經驗模態分解(EMD)及內蘊模態函數(IMF)重構,在保留微震信號固有特征基礎上實現降噪;其次,基于Hilbert變換計算出歸一化包絡信號,通過設置包絡閾值搜索震相初至的大致位置,并以該位置為基礎為AIC函數選擇合適的計算時窗以降低計算復雜度;最后,在選擇的時窗內應用AIC準則計算出P波震相初至。應用該方法對隨機選取的高、低信噪比兩類共80個微震信號進行震相初至拾取,若以人工拾取結果為基準、時差在10 ms以內視為準確拾取,結果顯示,高信噪比信號拾取準確率為100%,低信噪比信號拾取準確率為92%。

關鍵詞:微震P波;震相初至拾取;經驗模態分解;希爾伯特-黃變換;AIC準則

責任編輯:常 琛

賈瑞生,譚云亮,孫紅梅,等.低信噪比微震P波震相初至自動拾取方法[J].煤炭學報,2015,40(8):1845-1852.doi:10.13225/ j.cnki.jccs.2014.1122

微震監測是預警沖擊地壓、煤與瓦斯突出、礦井突水等煤礦災害事故的主要手段之一,目前已在微震信號頻譜特征提取、災害危險預警及防治等方面取得許多研究成果[1],這些成果的取得均建立在震源準確定位的基礎之上,而震源準確定位技術的關鍵環節之一則是微震P波初至到時的準確拾取[2-3],因此微震P波震相初至自動拾取是微震監測系統的關鍵技術之一,是震源定位及震源機制解釋的前提與基礎。

微震監測系統拾取的震動信號具有瞬態性、多樣性及不確定性等特點,信號中夾雜有機械震動、環境噪聲、電磁噪聲等多種干擾信號,因此,如何在低信噪比信號中準確辨識微震P波到時十分困難。對于地震事件,目前已發展多種震相自動識別及拾取方法, 如Allen等基于時間域信號能量變化構建特征函數,提出了長短時均值比方法(STA/ LTA)[4-6];Akaike等基于地震波可劃分為局部平穩段的假設,依據到達前后地震波形數據統計的差別,提出了判別震相初至的AIC準則[7-11];Saragiotis等基于地震波形的偏斜度和峰度提出了PAI-S/ K方法[12];常旭等應用分形理論對地震記錄進行了有效的初至自動拾取[13];馬強等綜合應用STA/ LTA方法及AIC準則,基于Delaunay三角剖分剔除干擾信號,提出了多步驟的P波自動拾取方法[14];王繼等應用單臺Akaike信息準則和多臺AIC最小二乘互相關方法,發展了震相自動精確檢測技術,實現了流動地震臺陣觀測震相初至的自動拾取[15];劉勁松等通過分析STA/ LTA,AIC,PAI-S/ K等幾種方法的原理及特點,提出了移動時窗峰度的快速算法和改進的峰度拾取初至算法[16]。以上方法均已在不同場合得到實際應用,能夠實現對較高信噪比地震波的到時拾取,但當信號信噪比較低時可能會出現拾取結果精度不夠甚至錯誤。對于微地震事件,其P波震相到時自動拾取的特殊性在于:相較于天然地震而言,微震信號頻率較高、信噪比低,自動拾取震相初至易受外部噪聲干擾。

為此,本文從信號分析理論入手,引入希爾伯特-黃變換對含噪微震信號進行降噪處理及包絡分析,通過設置包絡閾值為AIC算法選擇合適的計算時窗,再在時窗內應用AIC方法,實現了低信噪比微震P波震相初至的自動拾取。

1 微震P波震相自動拾取算法分析

1.1 長短時均值比(STA/ LTA)法

設x(t)為t時刻P波垂向幅值,f(t)為高通濾波后的垂向幅值數據,則有

其中,α為高通濾波因子,由于P波震相在垂直方向幅值較大,可以通過構造信號特征函數來放大垂直方向幅值,常用的P波拾取特征函數為

STA/ LTA方法的基本原理是根據微震P波特征函數值的長短時均值比判斷其初至[8-10],應用STA/ LTA的基本形式如下:

其中,t0為當前時刻;t1,t2為t0之前某時刻,且滿足t2

圖1　STA/ LTA法拾取P波到時Fig.1 Pick up P-arrival time based on STA/ LTA

STA/ LTA方法的突出優點是算法穩定可靠、拾取效率高,但其觸發點一般滯后于實際初至點(圖1),當微震信號中存在干擾信號時,通常會產生震相初至點的錯誤判斷[14],因而適用于高信噪比微震信號的震相初至拾取。

1.2 AIC方法

AIC方法的基本原理是求取地震信號AIC函數的局部最小值,Sleema提出了AR-AIC準則[10],它根據自回歸過程將地震波形數據分成2個局部統計時段(圖2),AR-AIC函數表示為

其中,k為兩個局部統計時段分界點;p為AR過程階數;l為地震波形數據長度;分別為2個局部統計時段的擬合誤差;C為一個常數。為了求出震相初至,必須求出該函數中AR模型的階數和系數,該方法計算復雜度較高,不利于震相初至的實時拾取。

圖2　AR-AIC法拾取P波初至Fig.2 Pick up P-arrival time based on AR-AIC

不同于AR-AIC模型,Maeda提出直接由地震波形數據計算AIC函數[11],求取AIC函數的局部最小值(圖3),該值對應的位置即為震相初至,AIC函數表示為

其中,x(i)(i=1,2,…,l)為地震波形離散數據;k的取值范圍是數據窗口內所有采樣點。

圖3　微震波形和時窗[500,2 500]內的AIC函數曲線Fig.3 A microseismic signal and AIC function curve whentime window between 500 ms to 2 500 ms

對比以上2種AIC方法,后者不需要計算AR模型的階數即可直接求取AIC值,在震相初至拾取實時性要求很高的情況下,是一種更加高效實用的算法,但是該算法需要在震相初至的附近尋找一個合適的時窗來計算AIC值,這是因為不同的時窗可能使AIC函數局部最小值出現的位置不同。圖3(a),4(a)為同一個微震波形數據在不同時窗內的波形,圖3(b),4(b)為對應的AIC函數曲線,其中圖4由于時窗設置不合理導致震相初至拾取錯誤,如何合理選擇時窗是AIC法準確拾取震相初至的關鍵問題之一。

圖4　微震波形和時窗[1 000,3 000]內的AIC函數曲線Fig.4 A micro-seismic signal and AIC function figure whentime window between 1 000 ms to 3 000 ms

此外,微震監測過程中拾取的信號易受外界噪聲干擾,使得AIC法震相初至拾取的精度大打折扣,有時還會造成微震事件的誤拾取[10-11],因此,如何降低環境噪聲對震相初至拾取的影響成為目前該領域研究的關鍵問題。由于微震信號及外部噪聲具有隨機性、非平穩性等特征,常規降噪方法(如Fourier方法)無法有效消除微震信號中夾雜的外部噪聲,而HHT則能對非平穩信號進行EMD分解并重構,在保留信號固有非平穩特征基礎上實現降噪。

2 對震相自動拾取算法的改進

2.1 EMD變換及信號降噪處理

Huang等提出了分析處理非平穩信號的EMD方法[17],即Huang變換,其核心是:任何非平穩信號都是由一些互不相同的、非正弦的、簡單的IMF函數組成的,因此可將信號分解成從高頻到低頻若干個基本時間序列的組合,每個基本時間序列即是一個IMF分量。基于HHT的核心思想,含噪微震信號的EMD分解及消噪算法如下:

(1)確定微震信號x(t)時間曲線上所有局部極大、極小值點,分別將其用3次樣條函數擬合為原數據序列的上、下包絡線,上、下包絡線對應的坐標均值連接成均值線m1(t),由式(8)計算出h1(t),即

(2)判斷h1(t)是否滿足IMF函數所需條件,令

當SD值介于[0.2,0.3]之間時,篩選過程終止,得到信號x(t)的第1個IMF分量;若SD值不在[0.2,0.3]之間,此時把h1(t)看作新的時間序列曲線,重復以上步驟繼續分解,則有

經過k次迭代后,當SD值介于[0.2,0.3]之間時,就得到了第1階IMF分量h1k(t),將其記為c1(t),它表示微震信號的最高頻部分。

(3)從原信號x(t)中減去該信號的高頻部分c1(t),得到頻率較低的信號殘差r1(t):

(4)將r1(t)看作一個新的信號序列,重復以上步驟,得到一系列的ci(t)和最后一個不可分解的信號殘差r(t),則原信號可表達如下:

(5)剔除高頻IMF分量,對剩余的IMF分量進行重構,得到降噪后的微震信號:

圖5　一個含噪微震信號的EMD分解Fig.5 Empirical mode decomposition of a micro seismic signal with noise

圖5顯示,含噪信號經EMD分解后得到7個IMF分量,這些IMF分量從高頻到低頻依次排列,由于微震頻率范圍一般介于0~200 Hz之間,觀察IMF分量的頻譜分布,可以判斷出imf1為高頻噪聲成分,故可將其剔除,重構剩余的IMF即得到降噪后的微震有效信號。圖6為EMD降噪前后微震信號的波形及頻譜,對比圖6(a),(c)可知,經過EMD降噪后較好地保留了信號的尖峰、突變及波形特征,說明EMD法能充分保留信號的瞬態非平穩特征;對比圖6(b),(d)可知,圖6(d)可清晰辨別出微震信號的時頻特征,無用高頻干擾噪聲被有效濾除,為進一步識別震相初至奠定基礎。

2.2 Hilbert變換及震相初至拾取算法

Hilbert變換是非平穩信號分析與處理的有效方法,實信號x(t)的Hilbert變換可定義[18]如下

Hilbert變換即是信號和時間倒數的卷積,實信號x(t)和它的Hilbert變換可組成解析信號:

其中, a(t)為解析信號的幅值;φ(t)為解析信號的相位,分別表示為

式中,a(t)為實信號x(t)的包絡。

將降噪后的微震信號進行Hilbert變換并進行包絡分析,即按時序將上下包絡的絕對值相加并將其歸一化處理,則圖6(c)的Hilbert包絡信號波形如圖7所示,對同一類微震信號使用同一閾值(本文設為0.3)沿時間軸順序搜索,當包絡值大于設定閾值時,即可得到P波震相初至的大致位置t0,本例中t0= 1 243 ms。在經過HHT降噪后的微震信號波形上以t0為基準向前及向后分別取500個采樣點作為計算時窗,在該時窗內應用AIC函數求解P波震相初至,求解結果如圖8所示。

圖6　去噪前后的微震信號波形及其頻譜Fig.6 Microseismic signal waveform and its spectrum before and after de-noising

圖7　Hilbert包絡信號波形Fig.7 Waveform of Hilbert envelope signal

圖8　時窗[744,1 743]內的微震波形和對應AIC函數曲線Fig.8 Microseismic signal when time window between 744 msand 1 743 ms and the corresponding AIC function curve

綜上,基于HHT及AIC求解震相初至的改進算法(簡稱HHT-AIC法)如下:①對含噪微震信號進行EMD分解,按照頻率從高到低獲得一系列IMF分量,根據外部噪聲特征對IMF分量進行剔除,并把剩余的IMF分量進行重構,實現微震信號的降噪;②對降噪后的微震信號進行Hilbert變換,求出其包絡信號,并對其進行歸一化處理;③設置包絡閾值,沿時間軸查找第1個大于包絡閾值的時刻,該時刻即為震相初至的粗略值,以該時刻為基準分別向前、向后各取500個采樣點作為計算時窗;④在選取的時窗內應用式(7)計算AIC函數值,AIC函數局部最小值對應的時刻即為震相初至。

3 算法檢驗及討論

3.1 實驗數據來源

為檢驗算法的有效性,我們分別選取一組高信噪比(標記為H組)及一組低信噪比(標記為L組)微震信號進行對比試驗,其中H組信號來源于我國西部某水電站微震監測數據,該水電站所在地主要為硬巖結構,巖石破裂產生的微震信號起跳明顯(圖9(a)),監測環境較為理想,信號中噪聲干擾較少,從中隨機抽取30組微震信號用于實驗;L組微震信號來源于我國西北某煤礦井下采場,屬軟巖地質構造,且機械振動、運輸、放炮等產生的外部噪聲較多,所采集的微震信號信噪比低,波形起跳不明顯(圖9(b)),從中隨機抽取50組微震信號用于實驗。

圖9　高信噪比和低信噪比微震波形Fig.9 Micro seismic waveform with high signal to noise ratio and Micro seismic waveform with low signal to noise ratio

分別應用AIC方法及本文方法對上述2組不同類型的微震信號進行處理,并將其與人工拾取的震相初至結果對比,從拾取結果的正確率、算法耗時等方面進行評價。

3.2 實驗結果及分析

應用Matlab編程實現AIC算法及HHT-AIC算法,對圖9(a)所示的高信噪比微震波形進行震相初至拾取,與人工拾取的結果對比如圖10所示。

圖10所示的AIC及HHT-AIC兩種拾取方法都獲得了同樣的震相初至時刻739 ms,在相同軟硬件配置的PC機上,AIC法運算耗時0.641 s,HHT-AIC法耗時1.892 s,AIC法耗時最少,但AIC法在選擇時窗大小及位置不合理時可能會出現圖4(b)出現的錯誤拾取結果,故AIC法是不穩定的;HHT-AIC法包括EMD降噪(耗時1.689 s)、Hilbert包絡分析及AIC函數計算(耗時0.203 s)等過程,計算復雜度較高,但當處理高信噪比信號時可以省去EMD降噪過程,僅保留Hilbert包絡分析及AIC函數計算過程(以下稱HT-AIC法),此時算法僅耗時0.203 s,大大低于AIC法的計算耗時,且能保持震相初至拾取精度不變,這是因為HHT-AIC法利用Hilbert變換進行包絡分析時為AIC函數選擇了合適的計算時窗(圖10(c)),時窗內僅有1 000個采樣點數據,大大降低了AIC函數的計算強度。

對H組中30個高信噪比微震信號進行處理,以人工拾取的震相初至作為參照,震相初至拾取結果見表1。其中AIC法的計算時窗為P波初至前750個至P波初至后2 250個采樣點數據,時窗內共3 000個采樣點數據;HHT-AIC時窗為P波初至前后各500個采樣點數據,時窗內共1 000個采樣點數據; HT-AIC法所采用的時窗位置及大小與HHT-AIC方法所采用的時窗相同。

從表1可以看出,對于H組中的高信噪比微震信號,3種方法均獲得了準確率很高的震相初至拾取結果,對比它們的算法平均耗時,可以看出HT-AIC算法耗時最少。

圖10　使用AIC法及HHT-AIC法拾取高信噪比微震P波初至結果對比Fig.10 Detection on micro seismic P-arrival time by using AIC and HHT-AIC method under high SNR

表1　高信噪比微震信號震相初至拾取結果Table 1 Result of detection on microseismic P-arrival time

分別應用AIC,HHT-AIC方法對圖9(b)所示的低信噪比微震信號進行震相初至拾取,與人工拾取的結果對比如圖11所示。圖11(a)為原始微震信號,起跳位置較為模糊,即使人工辨識信號起跳點也存在一定困難,通過信號局部放大,人工拾取的震相初至t=631 ms;圖11(b)為未經消噪處理直接應用AIC法求取的震相初至t=2 326 ms,拾取結果出現錯誤; 圖11(c)為經過EMD降噪后的微震信號,相較于圖11(a)而言信號波形較為清晰,且充分保留了微震信號的隨機非平穩特征;圖11(d)為在消噪基礎上應用HHT-AIC法計算得到的震相初至t=617 ms,計算結果與人工拾取結果相比誤差在20 ms之內,說明HHT-AIC法能有效處理低信噪比信號震相初至的自動拾取,算法具有較強的抗噪性能。

從算法耗時來看,AIC方法耗時0.635 s,但自動拾取的震相初至出現錯誤; HHT - AIC法耗時1.912 s,雖然算法耗時較長,但獲得了正確的拾取結果。若在圖11(a)所示的微震波形上直接應用HTAIC法,算法耗時0.204 s,但拾取的震相初至為927 ms,與人工拾取結果631 ms相比產生了較大誤差,這是因為該含噪信號波形起跳不明顯,HT-AIC算法通過包絡閾值搜尋震相初至的大致位置出現了偏差,導致為AIC函數的選取計算時窗不合理,從而使AIC函數計算得到的震相初至出現較大誤差。因此,HT-AIC算法在處理低信噪比微震信號震相初至拾取時也是不穩定。

對L組中50個低信噪比微震信號進行處理,以人工拾取的震相初至作為參照,震相初至拾取結果見表2。表2中HHT-AIC法震相初至的拾取準確率最高,而HT-AIC法耗時最少。綜合表1和2的拾取結果可以看出:對于高信噪比微震信號,HT-AIC法利用Hilbert變換對信號震相初至位置進行了粗略判定,為AIC函數選擇了位置及大小合適的計算時窗,有效降低了AIC函數的計算強度,因而能在保持較高識別準確率的前提下使算法耗時最少;但HT-AIC方法在處理低信噪比微震信號時,如果不對含噪信號進行降噪處理,仍會導致計算時窗位置選擇不合理,最終出現震相拾取結果錯誤;HHT-AIC法則首先應用EMD法對低信噪比微震信號進行降噪,提高了信號的信噪比,然后再應用HT-AIC算法實現震相初至的自動拾取,識別準確率比AIC法提高約10%,比HT-AIC法提高約6%。

圖11　使用AIC法及HHT-AIC法拾取低信噪比微震P波初至結果對比Fig.11 Detection on micro seismic P-arrival time by using AIC and HHT-AIC method under low SNR

表2　低信噪比微震信號震相初至拾取結果Table 2 Result of detection on microseismic P-arrival time

3.3 關于震相自動拾取方法的討論

AIC算法在時窗恰好包含有震相初至前后各一段波形的前提下,能夠得到較好的拾取結果,故AIC算法適用于已知震相初至大致位置的情況,工程上通常與STA/ LTA法組合使用,即首先使用STA/ LTA法獲取震相初至的大致位置,再以該位置為基礎向前及向后各取一段數據作為AIC函數的時窗計算出震相初至位置,但在處理低信噪比微震信號時,應用AIC 及STA/ LTA方法時震相初至辨識的準確率均顯著下降,算法抗噪性能較差。

HHT-AIC算法通過信號的EMD分解與IMF分量重構實現降噪,充分保留了微震信號固有的隨機非平穩特征,并在消噪濾波的上應用Hilbert變換求得歸一化包絡信號波形,通過設定包絡閾值搜索震相初至的大致位置后,再應用AIC法求解震相初至,因而HHT-AIC算法具有較強的抗噪性能,在相同信噪比情況下,震相初至辨識的準確率顯著提高。

4 結 論

(1)對于高信噪比微震信號,HT-AIC方法能依據Hilbert變換及包絡閾值為AIC函數確定計算時窗位置及大小,包絡閾值參數選取簡單,自動拾取的震相初至與人工拾取結果一致,誤差在10 ms內拾取占比為100%,且算法實時性強,克服了傳統AIC方法由于時窗選擇不合理造成初至錯誤拾取的不足。

(2)對于低信噪比微震信號,HHT-AIC方法基于EMD分解及IMF重構兩個過程對含噪信號進行降噪,充分保留了微震信號固有的瞬態非平穩特征,在此基礎上再應用HT-AIC方法實現震相初至的自動拾取,通過實驗與人工拾取結果相比,時差在10 ms內的微震信號占比為92%,自動拾取準確率高于AIC方法及HT-AIC方法,表現出較強的抗噪性能;不足之處是算法復雜度較高,但以犧牲時間為代價換取精度也是慣常的選擇。

參考文獻:

[1]齊慶新,竇林名.沖擊地壓理論與技術[M].徐州:中國礦業大學出版社,2008:32-38.

[2]鞏思園,竇林名,馬小平,等.煤礦礦震定位中異向波速模型的構建與求解[J].地球物理學報,2012,55(5):1757-1763.

Gong Siyuan,Dou Linming,Ma Xiaoping,et al.Study on the construction and solution technique of anisotropic velocity model in the location of coal mine tremor[J].Chinese Journal of Geophysics, 2012,55(5):1757-1763.

[3]朱權潔,姜福興,王存文,等.微震波自動拾取與多通道聯合定位優化[J].煤炭學報,2013,38(8):397-403.

Zhu Quanjie, Jiang Fuxing, Wang Cunwen, et al.Automated micro-seismic event arrival picking and multi-channel recognition and location[J].Journal of China Coal Society,2013,38(8):397-403.[4]Allen R V.Automatic earthquake recognition and timing from single traces[J].Bull.Seismol.Soc.Amer.,1978,68 (5):1521 -1532.

[5]Allen R V.Automatic phase pickers:Their present use and future prospects[J].Bull.Seismol.Soc.Amer.,1982,72(6):225-242.[6]Baer M,Kardolfer U.An automatic phase picker for local and tele-seismic events[J].Bull.Seismol.Soc.Amer.,1987,77 (4): 1437-1445.

[7]Akaike H.Information theory and an extension of the maximum likelihood principle[A].2ndInternational Symposium on Information Theory (Tsahkadsor)[C].1971:267-281.

[8]Takanami T,Kitagawa G.Multivariate time-series model to estimate the arrival times of S-waves[J].Computers and Geosciences, 1993,19(2):295-301.

[9]Leonard M,Kennett M B L N.Multi-component autoregressive techniques for the analysis of seismograms[J].Phys.Earth Planet.Interiors,1999,113(1-4):247-264.

[10]Sleeman R,Orild V E.Robust automatic P-phase picking:an online implementation in the analysis of broadband seismogram recordings [J].Phys.Earth Planet.Interiors,1999,113(1-4):265-275.

[11]Maeda N.A method for reading and checking phase times in autoprocessing system of seismic wave data[J].J.Seismol.Soc.Jpn., 1985,38(3):365-379.

[12]Saragiotis, Christos D.PAI-S/ K: A robust automatic seismic P phase arrival identification scheme[J].IEEE Transactions on Geosciences and Remote Sensing,2002,40(6):1395.

[13]常 旭,劉伊克.地震記錄的廣義分維及其應用[J].地球物理學報,2002,11(6):839-846.

Chang Xu,Liu Yike.The generalized fractal dimension of seismic records and its application[J].Chinese Journal of Geophysics, 2002,11(6):839-846.

[14]馬 強,金 星,李山有,等.用于地震預警的P波震相到時自動拾取[J].地球物理學報,2013,56(7):2313-2321.

Ma Qiang,Jin Xing,Li Shanyou,et al.Automatic P-arrival detection for earthquake early warning[J].Chinese Journal of Geophysics,2013,56(7):2313-2321.

[15]王 繼,陳九輝,劉啟元,等.流動地震臺陣觀測初至震相的自動檢測[J].地震學報,2006,28(1):42-51.

Wang Ji,Chen Jiuhui,Liu Qiyuan,et al.Automatic onset phase picking for portable seismic array observation[J].Acta Seismologica Sinica,2006,28(1):42-51.

[16]劉勁松,王 赟,姚振興.微地震信號到時自動拾取方法[J].地球物理學報,2013,56(5):1660-1666.

Liu Jinsong,Wang Yun,Yao Zhenxing.On micro-seismic first arrival identification:A case study[J].Chinese Journal of Geophysics, 2013,56(5):1660-1666.

[17]Huang N E,Zheng S,Steven R,et al.The empirical mode decomposition and the Hilbert spectrum for nonlinear non-stationary time series analysis[A].Proceedings:Mathematical,Physical and Engineering Sciences [C].London: The Royal Society Press,1998: 903-995.

[18]Rashid A.IIR discrete-time Hilbert transformers[J].IEEE Transactions on Acoustics,Speech,and Signal Processing,1987,35(8): 1116-1119.

Jia Ruisheng,Tan Yunliang,Sun Hongmei,et al.Method of automatic detection on micro-seismic P-arrival time under low signal to noise ratio[J].Journal of China Coal Society,2015,40(8):1845-1852.doi:10.13225/ j.cnki.jccs.2014.1122

Method of automatic detection on micro-seismic P-arrival time under low signal-to-noise ratio

JIA Rui-sheng1,2,TAN Yun-liang1,3,SUN Hong-mei2,HONG Yong-fa2

(1.State Key Laboratory of Mining Disaster Prevention and Control Co-founded by Shandong Province and Ministry of Science and Technology,Shandong University of Science and Technology,Qingdao 266590,China;2.College of Information Science and Engineering,Shandong University of Science and Technology,Qingdao 266590,China;3.College of Mining and Safety Engineering,Shandong University of Science and Technology,Qingdao 266590,China)

Abstract:According to the low signal-to-noise ratio,non-stationary,and inherent randomness of micro-seismic signals, a new method based on Hilbert-Huang Transform (HHT) and AIC criterion was proposed.Firstly,to implement denoising and filtering,the signals with noise were decomposed with empirical mode decomposition (EMD) and then reconstructed with intrinsic mode function (IMF);Secondly,based on the Hilbert transform,the normalized envelope signal was calculated,the coarse position of P-arrival was searched by setting the envelope threshold,and on the basis of the position,the time window of AIC function was determined;Finally,the P-arrival time was calculated by using AIC criterion in the selected time window.This method was applied to 80 random micro-seismic signals with both high SNR and low SNR.Considering the manual picking results as a benchmark and the time difference within 10 ms as an accurate result,the result shows that for the signals with high SNR,the picking accuracy is 100%.For the signals with low SNR,the accuracy is 92%.

Key words:micro-seismic P-arrival;seismic phase first arrival identification;empirical mode decomposition;Hilbert-Huang transform;Akaike information criterion

作者簡介:賈瑞生(1972—),男,安徽碭山人,教授,博士。E-mail:jrs716@163.com

基金項目:山東省自然科學基金資助項目(ZR2013EEM019);國家“十二五”科技支撐計劃資助項目(2012BAK04B06)

收稿日期:2014-08-30

中圖分類號:TD32

文獻標志碼:A

文章編號:0253-9993(2015)08-1845-08