楊 勇 張忠政 程 浩
(1.鞍鋼集團(tuán)礦業(yè)弓長(zhǎng)嶺有限公司露采分公司,遼寧 遼陽(yáng) 111007;2.東北大學(xué)資源與土木工程學(xué)院,遼寧 沈陽(yáng) 110819)
在金屬礦山開(kāi)采過(guò)程中,微震監(jiān)測(cè)作為礦山安全開(kāi)采的主要監(jiān)測(cè)手段之一,得到了廣泛的應(yīng)用[1]。通過(guò)監(jiān)測(cè)井下作業(yè)、巖體破裂等因素激發(fā)的彈性波信號(hào),利用拾取到的不同傳感器的信號(hào)進(jìn)行微震源定位,進(jìn)而為預(yù)警井下事故提供依據(jù)。微震數(shù)據(jù)有效信號(hào)的到時(shí)拾取是微震源定位處理過(guò)程的關(guān)鍵步驟之一,直接影響反演結(jié)果的準(zhǔn)確性。但金屬礦山開(kāi)采環(huán)境復(fù)雜,通常所采集到的微震數(shù)據(jù)中包括了大量的隨機(jī)噪聲,嚴(yán)重影響了微震數(shù)據(jù)有效信號(hào)到時(shí)拾取的精度。所以,對(duì)微震數(shù)據(jù)進(jìn)行隨機(jī)噪聲去除處理是提高有效信號(hào)到時(shí)拾取精度的關(guān)鍵步驟。
時(shí)頻峰值濾波(Time-Frequency Peak Filtering,TFPF)作為一種時(shí)頻分析方法,它通過(guò)求解信號(hào)的偽維格納—維利分布(Pseudo Wigner-Ville Distribution,PWVD)頻譜來(lái)估計(jì)瞬時(shí)頻率,可以克服變換域中基函數(shù)的影響,同時(shí)又無(wú)需任何的假定條件即可達(dá)到去除噪聲的目的[2]。憑借其獨(dú)特的優(yōu)越性和有效性,TFPF 已逐步地被廣泛應(yīng)用于多種類數(shù)據(jù)的隨機(jī)噪聲去除。林紅波等[3-8]對(duì)TFPF 去噪方法進(jìn)行了較為系統(tǒng)的研究,并將其應(yīng)用于實(shí)際地震數(shù)據(jù)隨機(jī)噪聲去除。但由于TFPF 較為適宜處理平穩(wěn)信號(hào),所以,固定窗長(zhǎng)制約了對(duì)于非平穩(wěn)信號(hào)隨機(jī)噪聲去除的效果。
上世紀(jì)90 年代,Huang 等[9-12]針對(duì)信號(hào)非平穩(wěn)問(wèn)題提出了經(jīng)驗(yàn)?zāi)B(tài)分解(Empirical Mode Decomposition,EMD)方法。該方法是一種自適應(yīng)分解方法,無(wú)需任何假定條件即可將信號(hào)分解成模態(tài)分量,突出信號(hào)的物理特性。但是這種方法會(huì)帶來(lái)模態(tài)混疊的問(wèn)題。針對(duì)上述問(wèn)題,變分模態(tài)分解方法(Variational Mode Decomposition,VMD)被提出,該方法通過(guò)估算各信號(hào)分量的信號(hào)主頻與帶寬,將原始信號(hào)分解成主頻不同帶寬不同的信號(hào)分量,解決了模態(tài)混疊的問(wèn)題,使信號(hào)分解成頻率不同的近似平穩(wěn)信號(hào)。基于VMD 分解方法的優(yōu)越性,劉沖等[13]和喬云等[14]將其結(jié)合于小波變換應(yīng)用于不同領(lǐng)域的去噪研究。
通過(guò)上述TFPF 與VMD 方法的特點(diǎn)分析,本研究提出了一種基于VMD 的TFPF 隨機(jī)噪聲壓制方法。首先,將含噪信號(hào)通過(guò)VMD 分解,獲取若干個(gè)不同主頻占優(yōu)的近似平穩(wěn)信號(hào)分量;然后,對(duì)不同頻率的信號(hào)分量靈活地選取不同的窗長(zhǎng)進(jìn)行TFPF 去噪,得到處理后的各信號(hào)分量;最后,進(jìn)行信號(hào)重構(gòu),獲取隨機(jī)噪聲去除以后的信號(hào)。并將通過(guò)理論與實(shí)際微震數(shù)據(jù)進(jìn)行該方法有效性的驗(yàn)證。
VMD 是一種基于模態(tài)分解的自適應(yīng)性較高的信號(hào)分解算法。VMD 通過(guò)將模態(tài)估計(jì)轉(zhuǎn)化為變分問(wèn)題,解決該變分問(wèn)題的構(gòu)造和求解來(lái)確定非平穩(wěn)信號(hào)中不同的信號(hào)主頻與帶寬,將原始信號(hào)分解成K個(gè)具有一定帶寬的信號(hào)模態(tài)分量,進(jìn)而將有效信號(hào)和噪音分為不同的模態(tài)分量區(qū)分開(kāi)來(lái)。VMD 分解的結(jié)果要求各模態(tài)分量的頻帶寬度之和最小,其分解出來(lái)的各階模態(tài)分量相加即可得到原始信號(hào)。
VMD 將信號(hào)看作是不同的具有一定帶寬的信號(hào)分量組成:

通過(guò)對(duì)信號(hào)分量進(jìn)行希爾伯特變換以得到其單邊頻譜,隨后將各信號(hào)分量頻譜調(diào)整到以估計(jì)中心頻率為中心的基頻帶上,并估計(jì)各信號(hào)分量的帶寬,由此引出約束性變分問(wèn)題。VMD 的約束方程式為

式中,{xk}為第k個(gè)分量;{ωk}為各信號(hào)分量的中心頻率;?t為Tikhonov 矩陣;δ(t)為脈沖函數(shù);*表示卷積;x(t)是原始信號(hào)。
為了解決式(2)約束性問(wèn)題,引入二次懲罰因子和Lagrange 乘法算子,將上述約束性變分問(wèn)題轉(zhuǎn)化為無(wú)約束優(yōu)化問(wèn)題:




式(6)和式(7)中,τ表示的是噪聲容限參數(shù),ε為給定的一個(gè)大于0 的判別精度。當(dāng)?shù)Y(jié)果滿足式(7)中的收斂條件,則迭代完成,原始信號(hào)被分解為K個(gè)信號(hào)分量。
TFPF 算法通過(guò)頻率調(diào)制因子μ對(duì)含噪信號(hào)進(jìn)行頻率調(diào)制,得到解析信號(hào),通過(guò)求解析信號(hào)的偽維格納—維利分布(PWVD)頻譜及其峰值,得到解析信號(hào)的瞬時(shí)頻率估計(jì),達(dá)到去噪的目的。其具體步驟如下:
(1)將含噪信號(hào)表示為有效信號(hào)和隨機(jī)噪聲的疊加:

式中,x(t)為有效信號(hào);r(t)為隨機(jī)噪聲。
(2)用一個(gè)頻率調(diào)制因子μ對(duì)包含隨機(jī)噪聲的信號(hào)進(jìn)行頻率調(diào)制:

(3)求解析信號(hào)的偽維格納—維利分布頻譜(PWVD):

(4)求解析信號(hào)的PWVD 分布頻譜的峰值,作為解析信號(hào)瞬時(shí)頻率估計(jì):

在TFPF 算法中,窗長(zhǎng)的選擇會(huì)對(duì)信號(hào)的去噪效果產(chǎn)生影響。通過(guò)理論數(shù)據(jù)對(duì)TFPF 固定窗長(zhǎng)的問(wèn)題進(jìn)行探討。對(duì)理論數(shù)據(jù)進(jìn)行加噪處理,得到含噪信號(hào)。再分別選取窗長(zhǎng)為7、11、15,對(duì)含噪信號(hào)進(jìn)行濾波,如圖1所示。從處理結(jié)果來(lái)看,各窗長(zhǎng)處理結(jié)果都能對(duì)噪聲進(jìn)行壓制,還原出原始信號(hào),說(shuō)明TFPF算法在噪聲去除的效果上能取得良好的效果。分別選取波峰處和低振幅信號(hào)處進(jìn)行局部顯示,可以觀察到在波峰處,短窗長(zhǎng)結(jié)果有更好的信號(hào)擬合度,長(zhǎng)窗長(zhǎng)則存在波峰擬合度較低的問(wèn)題;在低振幅信號(hào)處,長(zhǎng)窗長(zhǎng)對(duì)噪聲的去除效果更加明顯,而短窗長(zhǎng)去除噪聲的效果相比長(zhǎng)窗長(zhǎng)較差。

圖1 TFPF 固定窗長(zhǎng)影響Fig.1 Effect of fixed window length of TFPF
根據(jù)地震數(shù)據(jù)的經(jīng)驗(yàn)公式,TFPF算法中固定窗長(zhǎng)的選取:

式中,fs為地震波的采樣頻率;fd為主頻。
由式(12)可知,TFPF 方法中窗長(zhǎng)長(zhǎng)度與微震信號(hào)的采樣頻率成正比,與微震信號(hào)的主頻成反比,若fs增大或者fd減小,則窗長(zhǎng)相應(yīng)變大,且窗長(zhǎng)取值通常為奇數(shù)。因此在利用TFPF 處理高頻信號(hào)時(shí),應(yīng)采用短窗長(zhǎng)進(jìn)行處理,而處理低頻信號(hào)時(shí),應(yīng)采用長(zhǎng)窗長(zhǎng)進(jìn)行處理。這樣TFPF 算法在進(jìn)行信號(hào)去噪的過(guò)程中既能有效地壓制噪聲也能更好地保持有效信號(hào)的幅值。
本文方法的處理流程如下:
(1)對(duì)含噪信號(hào)進(jìn)行VMD 分解,將信號(hào)分解成具有不同主頻占優(yōu)的若干個(gè)近似平穩(wěn)信號(hào)分量。
(2)通過(guò)上述窗長(zhǎng)選取原則,對(duì)高頻信號(hào)分量使用短窗長(zhǎng)進(jìn)行TFPF 去噪,對(duì)低頻信號(hào)分量采用長(zhǎng)窗長(zhǎng)進(jìn)行TFPF 去噪。
(3)將TFPF 去噪后的各階信號(hào)分量進(jìn)行精確的信號(hào)重構(gòu),獲取隨機(jī)噪聲壓制以后的結(jié)果。
選用3 個(gè)主頻分別為20,25,30 Hz,帶寬分別為2,3,4 的Ricker 子波組合作為人工合成信號(hào)進(jìn)行理論數(shù)據(jù)實(shí)驗(yàn)。該信號(hào)采樣率為1 ms,共計(jì)1 000 個(gè)采樣點(diǎn)。向該信號(hào)中加入隨機(jī)噪聲,使信號(hào)的信噪比降低至5 dB 左右,此時(shí),合成信號(hào)被噪聲污染嚴(yán)重,低頻處有效信號(hào)基本被噪聲淹沒(méi),原始有效信號(hào)無(wú)法被識(shí)別出來(lái)。利用VMD 分解將原始合成信號(hào)和加噪信號(hào)進(jìn)行信號(hào)分解,取分解信號(hào)分量個(gè)數(shù)K=4,可以得到主頻不同帶寬不同的4 階信號(hào)分量,如圖2所示。

圖2 理論數(shù)據(jù)模型Fig.2 Theoretical data model
對(duì)比原始合成信號(hào)和含噪信號(hào)的VMD 分解結(jié)果,可以觀察到原始合成信號(hào)經(jīng)過(guò)加噪處理后,其VMD 分解得到的4 階分量中第一信號(hào)分量出現(xiàn)低振幅噪聲,剩余三階信號(hào)分量基本被噪聲覆蓋,無(wú)法識(shí)別出有效信號(hào)。根據(jù)TFPF 窗長(zhǎng)選擇規(guī)則,分別對(duì)含噪各階信號(hào)分量靈活地選取利用TFPF 去噪的自適應(yīng)窗長(zhǎng),進(jìn)行TFPF 降噪處理,再將經(jīng)過(guò)TFPF 降噪后各階信號(hào)分量進(jìn)行信號(hào)重構(gòu),得到壓噪后的信號(hào)。進(jìn)行TFPF去噪后,各階信號(hào)分量的噪聲都有明顯的壓制(圖3)。對(duì)上述壓噪處理后各階分量進(jìn)行信號(hào)重構(gòu),即可得到去噪信號(hào)(圖4)。

圖3 各階信號(hào)分量對(duì)比Fig.3 Comparison of signal components of different order

圖4 人工合成數(shù)據(jù)壓噪信號(hào)Fig.4 Denoising signal of theoretical data
本文方法是將VMD 分解與TFPF 相結(jié)合,通過(guò)分解原始信號(hào)并將各信號(hào)分量以合理的窗長(zhǎng)進(jìn)行TFPF 去噪,為了證明其有效性,下面分別應(yīng)用不同固定窗長(zhǎng)TFPF 去噪方法與本文VMD-TFPF 去噪方法對(duì)上述理論數(shù)據(jù)進(jìn)行隨機(jī)噪聲的去除結(jié)果對(duì)比。為了更充分直觀地展示實(shí)驗(yàn)成果,將原始數(shù)據(jù)不同固定窗長(zhǎng)(取窗長(zhǎng)長(zhǎng)度為7、11、15)TFPF 去噪方法和VMD-TFPF 去噪方法的去噪結(jié)果置于同一參考系中,并將其波峰、波谷處局部放大進(jìn)行保幅、保真、去噪效果的對(duì)比,如圖5所示。

圖5 人工合成數(shù)據(jù)壓噪對(duì)比Fig.5 Comparison of denoising result of theoretical data
觀察去噪結(jié)果對(duì)比圖,固定窗長(zhǎng)TFPF 去噪和VMD-TFPF 去噪結(jié)果都可以有效地壓制隨機(jī)噪聲,將原始信號(hào)還原出來(lái)。觀察波峰數(shù)據(jù)(圖5(b)),本文方法在去噪效果和保幅保真上達(dá)到平衡,與原始數(shù)據(jù)擬合度更高,而窗長(zhǎng)為7 的TFPF 去噪效果其有效信號(hào)遭到破壞,而窗長(zhǎng)為11、15 的TFPF 去噪效果其幅值有一定的衰減。觀察波谷數(shù)據(jù)(圖5(c)),亦可見(jiàn),本文的去噪結(jié)果與原始信號(hào)擬合度更高。
為了進(jìn)一步說(shuō)明本文方法的有效性,表1 給出了含噪信號(hào)與去噪所得信號(hào)均方誤差(MSE)、信噪比(SNR)和峰值信噪比(PSNR)的對(duì)比。

表1 合成數(shù)據(jù)實(shí)驗(yàn)的MSE、SNR和PSNR 對(duì)比Table 1 Comparison of MSE,SNR and PSNR for the synthetic-data experiment
理論數(shù)據(jù)實(shí)驗(yàn)結(jié)果表明,本文方法的MSE 值比固定窗長(zhǎng)TFPF 去噪方法更低,與原始信號(hào)擬合度更高;而在信噪比和峰值信噪比的指標(biāo)上,VMD-TFPF去噪方法和固定窗長(zhǎng)TFPF 去噪方法都有顯著的提高,說(shuō)明結(jié)果能夠壓制噪聲,達(dá)到還原出原始信號(hào)的目的,但本文方法的SNR和PSNR 比固定窗長(zhǎng)TFPF去噪方法更高,其去噪效果最佳。
選取某地實(shí)際礦山微震數(shù)據(jù)(圖6)進(jìn)行實(shí)際數(shù)據(jù)VMD-TFPF 方法的應(yīng)用測(cè)試,觀察可見(jiàn),實(shí)際信號(hào)中包含有大量的隨機(jī)噪聲,這對(duì)有效信號(hào)的識(shí)別產(chǎn)生了極大的干擾。對(duì)該微震數(shù)據(jù)進(jìn)行傅里葉變換,結(jié)果如圖7所示。經(jīng)過(guò)傅里葉變換后,分析其頻率特征發(fā)現(xiàn),該信號(hào)頻率主頻基本介于0 ~50 Hz 之間。由該地區(qū)微震數(shù)據(jù)的頻譜特征可知,數(shù)據(jù)的有效信號(hào)部分具有低頻特性,而在實(shí)際數(shù)據(jù)的傅里葉變換圖中可見(jiàn),該信號(hào)在高頻部分也有分布,說(shuō)明該部分是由隨機(jī)噪聲引起的,因而應(yīng)該壓制引起高頻分布的隨機(jī)噪聲。

圖6 某地實(shí)際礦山微震信號(hào)Fig.6 Actual micro-seismic data

圖7 實(shí)際礦山微震信號(hào)頻譜圖Fig.7 Frequency spectrum of the actual micro-seismic data
如圖8所示,對(duì)該實(shí)際數(shù)據(jù)進(jìn)行VMD 分解,取K值為8,由式(12)TFPF 算法窗長(zhǎng)選取原則,對(duì)8 個(gè)分量分別選取自適應(yīng)窗長(zhǎng)進(jìn)行TFPF 去噪處理。最后將各個(gè)處理后的信號(hào)分量組進(jìn)行信號(hào)重構(gòu)操作,以還原出有效信號(hào)。觀察利用VMD-TFPF 方法壓噪處理的結(jié)果,并與實(shí)際數(shù)據(jù)進(jìn)行對(duì)比(圖9),顯示本文方法能夠很好地去除隨機(jī)噪聲,而且還保護(hù)了信號(hào)的幅值,較為真實(shí)地還原出了有效信號(hào)。

圖8 實(shí)際礦山微震信號(hào)VMD 分解各階分量圖Fig.8 VMD decomposition of actual mine microseismic signals

圖9 實(shí)際去噪結(jié)果對(duì)比Fig.9 Comparison of actual denoising effects
微震數(shù)據(jù)中隨機(jī)噪聲的有效去除,對(duì)于準(zhǔn)確識(shí)別微震事件的到達(dá)時(shí)間,進(jìn)而提高微震源定位的精度具有重要的實(shí)際意義。針對(duì)金屬礦山微震數(shù)據(jù)中隨機(jī)噪聲的去除問(wèn)題進(jìn)行研究,提出了一種基于VMD 的TFPF 隨機(jī)噪聲去除方法,充分結(jié)合2 種方法的優(yōu)勢(shì),通過(guò)理論數(shù)據(jù)和實(shí)際數(shù)據(jù)的測(cè)試,可得以下結(jié)論:
(1)微震數(shù)據(jù)通過(guò)VMD 方法進(jìn)行分解,可獲得近似平穩(wěn)的多個(gè)信號(hào)分量,避免了直接利用TFPF 方法進(jìn)行非平穩(wěn)信號(hào)的處理。
(2)TFPF 方法可針對(duì)不同近似平穩(wěn)的多個(gè)信號(hào)分量選取對(duì)應(yīng)的窗長(zhǎng),進(jìn)而達(dá)到最大限度地去除隨機(jī)噪聲、保留有效信號(hào)。
(3)基于VMD 的TFPF 去噪方法能夠有效地去除隨機(jī)噪聲,提高微震數(shù)據(jù)的信噪比,同時(shí),更好地保留有效信號(hào)的幅值。