一種基于dwt_sta/lta的含噪信號p波初至峰度拾取方法
【技術領域】
[0001] 本發明屬于信號處理技術領域，尤其是涉及一種礦山微震信號P波初至時刻聯合 拾取方法。
【背景技術】
[0002] P波初至拾取是礦山微震識別、震源定位及震源機制解釋中最為基礎和重要的一 步。目前主要采用人工拾取P波到時，但其易受個人經驗、情緒等因素影響，且拾取數據量 大時非常耗時。鑒于上述原因，國內外提出了一系列P波初至自動拾取方法，但到目前為止 還沒有任何一種算法能夠在不同的震源環境、震中距范圍、噪聲背景下進行微震初至的一 致性檢測，更不用說對后續震相初至的有效檢測和拾取。
[0003] 常見的P波拾取方法包括長短時窗平均值比法、分形維數法、自回歸模型、人工神 經網絡法、高階統計量法等。長短時窗平均值比法利用短時窗平均值與長時窗平均值之比 來反映信號幅值的變化，并取比值大于某一設定閾值時為P波到時。長短時窗平均值比法 具有算法簡單、計算速度快、平均值特征函數多樣等特點，但其拾取低信噪比信號困難，且 閥值較大時，則可能拾取不到P波初至或增大拾取誤差；閥值較小時，則可能過早的被噪音 觸發。分形維數法認為噪音與信號疊加時分形維數發生變化，并以此作為P波拾取的依據。 該算法抗噪性能好，拾取精度較高，但其對插值的準確性、時窗和步長的依賴性很高，稍有 不慎就會嚴重影響拾取結果，且其計算速度較慢。自回歸模型將微震信號分為兩個局部統 計時段，并取自相關最小值點作為P波到時。該算法拾取精度較高，但其對低信噪比和尾部 震蕩信號拾取不穩定，且其采用最小值點作為P波到時，因此對純噪音信號仍會拾取到時 (HaijiangZhang等，2003)。人工神經網絡法采用峰值振幅、時窗內均方根振幅比、峰值與 其前后峰值的包絡斜率及噪聲與信號的比值等作為神經網絡輸入因素，人工神經網絡法是 一種較為綜合的P波初至拾取方法，對微震波形適應性強，但其參數計算工作量大，學習速 度慢，實現復雜。高階統計量法（PAI-S/K法）由Saragiotis等（2002,2004)提出，取峰度 (或偏度）最大值點作為P波初至時刻，該算法借助于信號只包含噪音時，峰度值和偏度值 趨于0,而當P波到達時峰度和偏度值開始增加，并當滑動窗口內包含部分微震信號時達到 最大值。該算法拾取精度高，但未考慮尖刺和尾部震蕩信號對拾取的影響，由此可能產生較 大的拾取誤差，且其采用最大值拾取P波到時，對純噪音信號仍會拾取到時。
[0004] 可見現有的礦山微震信號P波初至拾取方法存在很大的局限，需要研究一種拾取 精度高、穩定性好的自動拾取方法。
【發明內容】
[0005] 本發明所要解決的技術問題是提供一種改進的含噪信號P波初至峰度拾取方法， 解決峰度法（PAI-K法）拾取低信噪比、刺突、尾部震蕩及純噪音信號不穩定的技術問題，該 礦山微震信號P波初至拾取方法適用性強、準確性高。
[0006] 針對礦山微震信號P波初至PAI-K法拾取精度較高，但對低信噪比、刺突、尾部震 蕩及純噪音信號拾取不穩定的技術問題，本發明提出了W-S/L-K拾取法，該方法借助于小 波分解高頻細節信號能夠較好的保留P波初至信息以及STA/LTA法能很好地克服刺突和尾 部震蕩的影響并確定P波大致到時，進而采用權重分析排除上述干擾，極大地降低了錯誤 拾取率，增強了P波拾取的準確性。
[0007] -種基于DWT_STA/LTA的含噪信號P波初至峰度拾取方法，包括以下幾個步驟：
[0008] 步驟1 :提取微震信號X(n)
[0009] 從礦山微震信號中提取待拾取信號X(n)，n= 1，2,…，N，其中，N為所述微震信號 的采樣點個數，取N= 400~600,微震信號采樣頻率600Hz;
[0010] 步驟2 :對微震信號采用離散小波變換DWT得到細節信號D1,i= 1，2, 3, 4 ;
[0011] 步驟3 :依據STA/LTA算法按照以下公式確定細節信號P波大致到時k1;
[0013] 其中，STA(k，D1)和LTA(k，D1)分別為短、長時窗的振幅平均值；x(n，D1)為細節信 號D1的第n個采樣信號；k為采樣的第k個點，k=WWA，Wm+1，…，N;WSTJPWm分別為短、 長時窗的長度，分別取10和60個米樣點；
[0014] 若細節信號DJt應的短時窗和長時窗的振幅平均值的比值A(k，D1)的最大值大 于設定閥值a，則取首個觸發點對應的采樣點序號為該細節信號P波大致到時ki，否 貝1J，記該細節信號P波大致到時k 0。其中，i= 1，2, 3, 4 ;
[0015] 步驟4 :計算細節信號的P波大致到時Ic1= 0的個數L，i= 1，2, 3, 4 ;
[0016] 步驟5 :利用權重分析計算P波初至準確到時；
[0017] 若步驟4中L>2,則認為該信號為純噪音信號，反之則采用以下權重公式計算P波 準確到時k' ：
[0019] 其中，1/8、3/8、3/8和1/8分別為04、03、0 2和01的權重系數洫1和1^1'分別為1?1的 最大振幅和局部峰度拾取到時，R1SD1在區間[ki-b^i+b]的信號，b為確定局部峰度拾取 區間的常量，取b= 10~15 ;
[0020]
計算X(n,Di)的滑動峰度值K(k，Di)，k為采樣的第k個點， 動時窗的長度，取M= 20~30,并取K(I^D1)最大值點作為局部峰度拾取到時k/ ;
[0021] 步驟6 :依據步驟5輸出拾取結果，該信號為純噪音信號或者P波初至準確到時為 k'。
[0022]所述設定閾值a= 2~3。
[0023] 利用MalIat算法快速實現離散小波變換。
[0024] 并采用低通濾波器h和高通濾波器g獲得原始信號c。（等同于本發明中的X(n)) 的小波系數{山，d2，…，七}，從而將微震信號X(n)分解得到細節信號D1 (取i= 1，2, 3, 4)。
[0027] 式中：c]+lik為第k個采樣點的近似信號值c]+1，d]+1，k為第k個采樣點的細節信號 值dj+1，k為采樣的第k個點，1為濾波器確定指標，j為尺度參數（j= 0, 1，…，J-l)，J為 最大分解層數。為簡化說明，記Cj+1= {(^+1』}和1+1= {d_j+lik}。
[0028] 設定Ici= 0時，A;= 0,k/ = 0,i= 1，2, 3, 4,簡化計算過程，提高計算速度。
[0029] 有益效果
[0030] 本發明提供了一種基于DWT_STA/LTA的含噪信號P波初至峰度拾取方法，包括如 下步驟：提取礦山微震信號；離散小波變換（DWT)得到微震細節信號；計算各細節信號的 STA/LTA值；判定STA/LTA最大值是否大于閥值a，若是則取首個觸發點對應的采樣點序號 為該細節信號的大致到時，否則記該細節信號到時為0 ;再對到時非0的細節信號進行局部 峰態拾取；最后采用權重分析確定P波最終到時。本發明拾取純噪音信號時，細節信號STA/ LTA觸發數通常小于2,設定STA/LTA觸發數大于等于2即可排除噪音信號；拾取低信噪比 信號時，STA/LTA對含噪音較大的細節信號無拾取，且細節信號無拾取的數目通常不大于 2,因此可用其余細節信號拾取初至，從而排除低信噪比的影響；通常刺突信號對STA/LTA 值影響較小，設定合適的STA/LTA觸發值即可排除刺突的影響；STA/LTA通常在P波初至時 已觸發，不會受尾部震蕩的影響。由此，通過一種基于DWT_STA/LTA的含噪信號P波初至峰 度拾取方法有效排除了低信噪比、刺突、尾部震蕩及純噪音信號的干擾，極大地降低了錯誤 拾取率，增強了P波拾取的準確性。
【附圖說明】
[0031] 圖1是本發明所述方法流程圖；
[0032] 圖2是典型波形W-S/L-K法實現過程圖，其中，（a)為礦山微震原始信號x(n)(n =1，2,…，400)，虛線為人工拾取P波初至點；(b)為原始信號對應的峰度值K(k)(k= 21，22,…，400)，括號中數值分別為PAI-K法拾取P波初至點及其對應的峰度值；（c)為小 波分解所得細節信號D1Q= 1，2, 3, 4) ;(d)為細節信號D1Q= 1，2, 3, 4)的STA(I^D1)/ LTA(I^D1)值A(I^D1) ;(e)為細節信號局部峰態拾取k/(i= 1，2,3,4);
[0033] 圖3是PAI-K和W-S/L-K法工程拾取效果對比圖，其中，（a)為PAI-K拾取誤差圖； (b)為W-S-L-K拾取誤差圖；
[0034] 圖4是PAI-K和W-S/L-K法拾取典型案例圖，其中，虛線、點線和實線分別對應人 工、PAI-K法和W-S/L-K法拾取到時，（a)為含刺突純噪音信號，（b)為低信噪比信號，（c) 為起震明顯信號，（d)為低信噪比且含刺突信號，（e)為極低信噪比、含刺突且尾部存在震 蕩信號，（f)為起震不明顯且尾部存在震蕩的低信噪比信號。
【具體實施方式