本發(fā)明屬于信號處理
技術(shù)領(lǐng)域：
，具體涉及一種基于VMD的礦山微震信號的降噪濾波方法。
背景技術(shù)：
：巖石發(fā)生破裂時誘發(fā)微震，形成微震數(shù)據(jù)，而煤礦井下噪聲污染嚴(yán)重，因此微震數(shù)據(jù)中包含了大量外部噪聲，需將微震有效信號從噪聲中分離出來。目前常用的巖石破裂微震信號的降噪濾波方法有經(jīng)驗?zāi)B(tài)分解(EMD)、集成經(jīng)驗?zāi)B(tài)分解(EEMD)，小波分析等，這些方法運算速度慢、抗噪性能差、誤判率高、拾取精度低、算法實時性不強(qiáng)。如EMD在分解過程中會產(chǎn)生模態(tài)混疊現(xiàn)象，即分解得到的一個或多個IMF中包含差異極大的特征時間尺度，信號和噪聲混疊在一個或多個IMF中，很難達(dá)到有效地降噪濾波效果。技術(shù)實現(xiàn)要素：針對現(xiàn)有技術(shù)中存在的上述技術(shù)問題，本發(fā)明提出了一種基于VMD的礦山微震信號的降噪濾波方法，變分模態(tài)分解(VMD)是一種新的信號分解方法，相比于其它模態(tài)分解技術(shù)，它有效地解決了模態(tài)混疊問題，具有很好的噪聲魯棒性，克服了現(xiàn)有技術(shù)的不足，具有良好的降噪效果。為了實現(xiàn)上述目的，本發(fā)明采用如下技術(shù)方案：一種基于VMD的礦山微震信號的降噪濾波方法，包括如下步驟：步驟1：讀取含噪微震信號x的時序序列x(t)，其中，t＝1，2，…，N，N為微震信號的采樣點個數(shù)；步驟2：對含噪微震信號x進(jìn)行VMD分解：對含噪微震信號x采用VMD分解為一系列變分模態(tài)分量，使各個模態(tài)的估計帶寬之和最小，約束條件為各模態(tài)之和等于輸入信號x，約束變分模型描述為式(1)和式(2)：s.t.∑kuk＝x(2)；式(1)中，{uk}:＝{u1,…,uK}為分解得到的K個有限帶寬的變分模態(tài)分量，{ωk}:＝{ω1,…,ωK}為各變分模態(tài)分量的頻率中心，δ(t)為狄拉克(Dirac)函數(shù)，*表示卷積，j2＝-1；式(2)中，x為含噪微震信號，為對所有的變分模態(tài)求和；為求解式(1)和式(2)的最優(yōu)解，引入擴(kuò)展的Lagrange將約束變分問題變?yōu)榉羌s束變分問題，其表達(dá)式為式(3)：其中，α為二次懲罰因子，λ(t)為拉格朗日乘法算子；求解該變分問題的具體步驟如下：步驟2.1：定義變分模態(tài)分量個數(shù)K值與懲罰因子α的值；步驟2.2：初始化步驟2.3：令n＝n+1，執(zhí)行整個循環(huán)；步驟2.4：執(zhí)行內(nèi)層第一個循環(huán)，根據(jù)式(4)更新uk；其中，為信號x(t)的傅立葉變換，步驟2.5：令k＝k+1，重復(fù)步驟2.4，直到k＝K，結(jié)束內(nèi)層第一個循環(huán)；步驟2.6：執(zhí)行內(nèi)層第二個循環(huán)，根據(jù)式(5)更新ωk；步驟2.7：令k＝k+1，重復(fù)步驟2.6，直到k＝K，結(jié)束內(nèi)層第二個循環(huán)；步驟2.8：執(zhí)行外層循環(huán)，根據(jù)式(6)更新λ；其中，τ為拉格朗日乘法算子λ(t)的更新步長參數(shù)；步驟2.9：重復(fù)步驟2.3至步驟2.8，直到滿足迭代停止條件如式(7)所示，結(jié)束整個循環(huán)，輸出結(jié)果，得到K個變分模態(tài)分量；其中，ε為求解精度；步驟3：計算原始微震信號x(t)與各變分模態(tài)分量uk的互相關(guān)系數(shù)；應(yīng)用式(8)計算含噪微震信號x的時序序列x(t)與各變分模態(tài)分量uk的互相關(guān)系數(shù)；其中，N為采樣點個數(shù)，且步驟4：將中心頻率大于200Hz且與原信號x的互相關(guān)系數(shù)小于0.3的變分模態(tài)分量作為噪聲濾除，對剩余的變分模態(tài)分量進(jìn)行重構(gòu)，得到降噪濾波后的微震信號。本發(fā)明原理如下：為實現(xiàn)微震信號的有效降噪濾波，本發(fā)明針對微震信號的隨機(jī)性、非平穩(wěn)性特點，讀取原含噪微震信號并進(jìn)行VMD分解，對分解后得到各個變分模態(tài)分量uk，計算原始微震信號x(t)與各變分模態(tài)分量uk的互相關(guān)系數(shù)，將中心頻率大于200Hz且與原信號互相關(guān)系數(shù)小于0.3的變分模態(tài)分量作為噪聲濾除，對剩余的變分模態(tài)分量進(jìn)行重構(gòu)，得到消噪濾波后的微震信號。本發(fā)明所帶來的有益技術(shù)效果：本發(fā)明依據(jù)VMD良好頻譜分解特征對低頻變分模態(tài)分量uk進(jìn)行重構(gòu)，能對實際信號進(jìn)行自適合分解，在充分保留微震信號隨機(jī)性和非平穩(wěn)的特征的基礎(chǔ)上，最大程度上對微震信號進(jìn)行濾波，該算法具有算法簡單、自適應(yīng)性和實時性強(qiáng)的特點，能對礦山微震監(jiān)測系統(tǒng)拾取的含噪微震信號進(jìn)行有效降噪濾波，具有很好的技術(shù)價值和應(yīng)用前景。附圖說明圖1為本發(fā)明一種基于VMD的礦山微震信號的降噪濾波方法的流程圖。圖2為含噪微震信號時序序列x(t)的示意圖。圖3為含噪微震信號x(t)經(jīng)VMD分解后得到的6個變分模態(tài)分量示意圖。圖4為6個變分模態(tài)分量對應(yīng)的頻譜圖。圖5為降噪后的微震信號時序序列示意圖。圖6為含噪微震信號與降噪后微震信號的對比圖。具體實施方式下面結(jié)合附圖以及具體實施方式對本發(fā)明作進(jìn)一步詳細(xì)說明：如圖1所示，一種基于VMD的礦山微震信號的降噪濾波方法，具體包括如下步驟：步驟1：讀取含噪微震監(jiān)測信號x的時序序列x(t)，其中，t＝1，2，…，N，N為微震信號的采樣點個數(shù)；步驟2：對含噪微震信號x進(jìn)行VMD分解：對含噪微震信號x采用VMD分解為一系列變分模態(tài)分量，使各個模態(tài)的估計帶寬之和最小，約束條件為各模態(tài)之和等于輸入信號x(t)，約束變分模型描述為式(1)和式(2)：s.t.∑kuk＝x(2)；式(1)中，{uk}:＝{u1,…,uK}為分解得到的K個有限帶寬的變分模態(tài)分量，{ωk}:＝{ω1,…,ωK}為各變分模態(tài)分量的頻率中心，δ(t)為狄拉克(Dirac)函數(shù)，*表示卷積，j2＝-1；式(2)中，x為原含噪微震信號，為對所有的變分模態(tài)求和。為求解式(1)和式(2)的最優(yōu)解，引入擴(kuò)展的Lagrange將約束變分問題變?yōu)榉羌s束變分問題，其表達(dá)式為式(3)：其中，α為二次懲罰因子，λ(t)為拉格朗日乘法算子。求解的具體步驟如下：步驟2.1：定義變分模態(tài)分量個數(shù)K值與懲罰因子α的值；步驟2.2：初始化步驟2.3：令n＝n+1，執(zhí)行整個循環(huán)；步驟2.4：執(zhí)行內(nèi)層第一個循環(huán)，根據(jù)式(4)更新uk；其中，為信號x(t)的傅立葉變換，步驟2.5：令k＝k+1，重復(fù)步驟2.4，直到k＝K，結(jié)束內(nèi)層第一個循環(huán)；步驟2.6：執(zhí)行內(nèi)層第二個循環(huán)，根據(jù)式(5)更新ωk；步驟2.7：令k＝k+1，重復(fù)步驟2.6，直到k＝K，結(jié)束內(nèi)層第二個循環(huán)；步驟2.8：執(zhí)行外層循環(huán)，根據(jù)式(6)更新λ；其中，τ為拉格朗日乘法算子λ(t)的更新步長參數(shù)；步驟2.9：重復(fù)步驟2.3至步驟2.8，直到滿足迭代停止條件如式(7)所示，結(jié)束整個循環(huán)，輸出結(jié)果，得到K個變分模態(tài)分量；其中，ε為求解精度；步驟3：計算原始微震信號x(t)與各變分模態(tài)分量uk的互相關(guān)系數(shù)；應(yīng)用式(8)計算原信號x(t)與各uk分量之間的互相關(guān)系數(shù)：其中，N為采樣點個數(shù)，且步驟4：將中心頻率大于200Hz且與原信號x互相關(guān)系數(shù)小于0.3的變分模態(tài)分量作為噪聲濾除，對剩余的變分模態(tài)分量進(jìn)行重構(gòu)，得到降噪濾波后的微震信號。如圖2所示，步驟1獲取以時間(ms)為橫軸，振幅為縱軸，微震信號時序序列表示為x(t),t＝1,2,…,3000，微震信號采樣點數(shù)據(jù)見表1。表1微震信號采樣點數(shù)據(jù)(可以存儲于Excel表中)序號時刻(時間/ms)振幅11‐1.83E‐0522‐0.12E‐0533‐0.68E‐05441.26E‐05550.34E‐05………29992999‐0.77E‐0530003000‐1.06E‐05按照步驟2的VMD算法對微震信號時序序列x(t)進(jìn)行變分模態(tài)分解，取K＝6，二次懲罰因子α＝5000，分解后的6個變分模態(tài)分量uk如圖3所示，6個變分模態(tài)分量對應(yīng)的頻譜圖如圖4所示。按照步驟3的方法，計算圖3中每個變分模態(tài)分量與原含噪微震信號x(t)的互相關(guān)系數(shù)，如表2所示。表2原信號與各uk的互相關(guān)系數(shù)R1R2R3R4R5R60.83960.41450.18320.19200.07350.0576圖3中每個變分模態(tài)分量的中心頻率值，如表3所示。表3各uk的中心頻率ω1ω2ω3ω4ω5ω625.583.6169.4238.4335.4418.6按照步驟4的方法，將中心頻率大于200Hz且與原信號互相關(guān)系數(shù)小于0.3的變分模態(tài)分量u4、u5、u6作為噪聲濾除，對剩余的變分模態(tài)分量u1、u2、u3進(jìn)行重構(gòu)，得到消噪濾波后的微震信號，如圖5所示。微震信號是隨機(jī)性、非平穩(wěn)的低頻信號，其頻率分布較為分散，主頻分布于0～200Hz頻帶，因此可將高頻部分的uk分量進(jìn)行剔除，將剩余的uk分量重構(gòu)，即可實現(xiàn)含噪微震信號的降噪濾波。本發(fā)明的優(yōu)點：通過VMD分解無需考慮模態(tài)混疊，在充分保留了微震信號本身的非平穩(wěn)特征的基礎(chǔ)上，減少和消除了信號中混雜的噪聲。如圖6所示，通過將含噪微震信號和降噪濾波后的微震信號對比可以看出，應(yīng)用VMD進(jìn)行微震信號的降噪濾波較完整地保留了原信號的波形特征，同時也較好地保留了原信號的尖峰與突變部分，得到了較好的降噪濾波效果。當(dāng)然，上述說明并非是對本發(fā)明的限制，本發(fā)明也并不僅限于上述舉例，本
技術(shù)領(lǐng)域：
的技術(shù)人員在本發(fā)明的實質(zhì)范圍內(nèi)所做出的變化、改型、添加或替換，也應(yīng)屬于本發(fā)明的保護(hù)范圍。當(dāng)前第1頁1 2 3