楊光亮
(1)中國地震局地震研究所,武漢 430071 2)地殼運動實驗室,武漢 430071)
一種去除地震背景噪音的新方法*
楊光亮1,2)
(1)中國地震局地震研究所,武漢 430071 2)地殼運動實驗室,武漢 430071)
結合希爾伯特-黃變換方法,根據信號背景噪音時間延續性假設,提出了一種自適應地震信號去噪算法。利用該算法實現了臺站地震監測信號的去噪分析。分析表明:1)該去噪算法能根據背景信號特性有效去除信號的低頻干擾;2)該算法可自適應完成信號分解去噪,計算效率高,無時間分辨率和頻率分辨率問題;3)該算法處理高頻干擾存在缺陷,有待進一步完善。
模態分解;固有模態函數;相關分析;特性時間延續;去噪
天然地震監測臺網接收到的有用信號淹沒在大量的干擾噪聲中,有些地震由于震級較小、信號信噪比低,將導致無法識別。微小地震反映殼內部的微破裂,如果能精確監測到微破裂的發展過程,對地震預警將起到重要作用。另外,由于大量的干擾噪音也影響對地震 P波初至的識別。因此,對地震監測弱信號的提取、去噪的研究具有重要實際意義。
地震監測信號是非平穩非各態歷經的隨機信號。目前廣泛用于地震信號去噪的方法,如傅立葉頻譜分析、功率譜分析、小波變換等方法,都是將地震信號近似為平穩隨機信號進行處理。非平穩信號的主要特征是其頻率是時間的函數。這種時變頻率雖然可以通過已有的時頻聯合分析方法,如短時傅立葉變換 (STFT)、W igner-Ville分布、小波變換(WT)等進行分析,但這些時頻分析法的最終理論依據都是傅立葉變換,在本質上還無法擺脫 Fourier變換的局限性,而且受 Heisenberg不確定原理的限制。
希爾伯特-黃變換 (HHT)是一種新的數據分析方法[1],它從根本上突破了傅立葉變換理論的限制,首次建立了一種基于瞬時頻率的信號分析方法,它可以比較精確地描繪出非平穩信號中頻率隨時間的變化規律,更容易避免假頻和多余信號分量等冗余問題的出現。它依據數據本身的時間尺度特征將信號分解為有限個固有模態函數 (I MF,然后對各模態分量進行變換從而得到信號能量在時間尺度上的分布規律,實現信號特征的量化提取。相對于傳統的數據分析方法,HHT有完全自適應性,能處理非線性非平穩數據,不受 Heisenberg測不準原理制約等優點,在客觀性和分辨率方面都具有明顯的優越性,能提取到更多、更接近實際的信號特性[2-5]。
地震信號具有短時、突變等特點,是一種典型的非平穩隨機信號。對于地震監測臺站接收到的信號,可以認為是背景噪音與有用信號的混合。假設背景噪音在一定時間范圍內是平穩的,并假設其信號特征在時間上有一定的延續性,據此可以將背景噪音作為信號基底 (圖 1),通過對基底和分析區信號的特征分析,可以尋找到兩者的內在關聯。

圖1 信號基底與分析區示意圖Fig.1 Schematic diagram of the backround signal and the analyzed area
HHT算法依據數據本身的時間尺度特征可以將信號分解為有限個固有模態函數 (I MF),經驗模態分解(EMD)的依據為:
1)整個信號中,零點數與極點數相等或至多相差 1。這一限制條件近似于傳統平穩高斯過程中關于窄帶的定義。

根據上述設想,我們設計了如下地震監測信號去噪算法,算法流程如圖 2所示。該算法需要解決3個關鍵問題:1)EMD分解時,包絡插值算法的選擇;2)EMD分解終止條件的選擇;3)相關閾值的設置。這里忽略了邊界效應的討論。Huang提出HHT時,采用三次樣條插值,但在實際使用過程中容易出現插值“過沖”和“欠沖”。這里我們采用線性插值算法計算包絡線,并采用 Rilling[6]提出的終止條件終止模態分解。

圖2 信號去噪流程Fig.2 Flowchart of denoising
相關閾值與信號長度密切相關。理論上,當分析信號和基底較長時,相互間的信號特征延續性較差,當信號較短時,則延續性較好。因此,相關閾值的選擇應根據信號長度及當地資料特性,經多次試驗來獲得。
圖3與圖 4中的信號是上海某地震監測臺一段連續垂直向地震信號,按上述算法分成基底區和分析區兩部分,作為分析數據,信號采樣頻率 100 Hz,信號總長度 600 s,分段后基底區與分析區均為 300 s。

圖3 基底原始信號Fig.3 Original background signal

圖4 待去噪的原始信號Fig.4 Orignal signal to be denoised
通過 EMD分解獲得了基底和分析信號的模態(圖 5)。基底分解后得到 13個基本模態和一個余量,即趨勢項。分析信號分解得到 11個模態和一個余量。原始信號經過模態分解后,得到從高頻到低頻的一系列模態分量 imf1,imf2,i mf3,imf4,…,這里并不是說 imf1的頻率一定比 imf2的高,而是 imf1中的某個局部的頻率比 imf2中相同局部的頻率要高,在相同時間的局部,imf1頻率一定大于 imf2。至于分解過程造成的誤差,主要是包絡方式的選取、邊界效應的處理和濾波停止條件的設計,會不斷累積到下一層分解中,并不一定是最后一個余量 (趨勢項)。
根據 EMD模態分解結果,各模態頻率成分由高到低排列,如果我們能從分析信號中識別噪音干擾,就可以據此對信號進行去噪。只需要將干擾信號模態直接剔出,然后將其他模態信號累加,就可得到去噪后的結果。在這種意義上說,模態分解相當于一個濾波器。

圖5 基底信號與分析信號模態分解結果Fig.5 Mode decomposition results of background signal and the analysis signal

圖6 基底與分析區信號相關系數Fig.6 Correlation coefficients between signals of background and analyzed area

圖7 噪音頻譜Fig.7 Noise spectrum
通過對分解出來的各固有模態作相關分析,計算得到相關系數關系如圖 6所示。相關系數有 4個大于 0.3的局部峰值,其中 2個集中在 0.3以上,另2個集中在 0.6以上,而其他峰值相對較小,因此,設置相關閾值分別為 0.3和 0.6時,可以得到 2組噪音模態,圖 7是兩組噪音模態的頻譜圖。相關閾值為 0.6時,剔除分析信號中最后 2個模態,其他模態通過模態累加就可以復合出去噪后信號,當相關閾值設置為 0.3時,剔除最后 4個模態,其他模態累加,可得到另一組復合后的分析信號(圖 8)。
該去噪算法是在假設背景噪音在一定時間范圍內是平穩的,并假設其信號特征在時間上有一定的延續性而作出的。一般認為時間序列越長,其信號特性延續性越差。為此,我們另外取了連續的兩段信號利用上述算法進行了計算,相關系數關系如圖8所示,x軸和 y軸分別為信號分解后的模態序數。從圖 8中可看到相關峰值明顯增多,相關平面不再平坦,其中相關系數大于 0.3的有 5個。該段分析信號去噪后時域和頻域結果如圖 9。通過實際資料的處理,表明該算法利用背景基底特征,能有效去除觀測信號中的部分噪音干擾。
從以上處理結果看,去除的干擾信息都集中在低頻端,這可能是低頻干擾信號更穩定,持續時間更長,更易被識別。由此可知,該去噪算法對高頻端的處理還有待進一步探索。
通過以上實例的分析,初步得出結論:1)該去噪算法能根據背景信號特性有效去除信號的低頻干擾;2)該算法去噪可自適應完成,計算效率高;3)該算法對高頻干擾存在缺陷,有待進一步完善。

圖8 去噪后基底與分析區信號相關系數Fig.8 Correlation coefficients between signals of background and analyzed area after denoising

圖9 10s信號去噪結果及頻譜圖Fig.9 Spectrum and denoising results of 10 s signal
1 HuangN E.The empiricalmode decomposition and the Hilbert spectrum for nonlinear and non-stationary time series analysis[J].Procedures of the Royal Cociety of London, 1998,A454(3):903-995.
2 Qiuhui Chen and Norden Huang.A B-spline approach for empirical mode decompositions[J].Advances in ComputationalMathematics,2006,24:171-195.
3 Norden E Huang Zhengshen and Steven R Long.A new view of nonlinearwaterwaves:the Hilbert spectrum[J].Annual Review of FluidMechanics,1999,31:417-457.
4 張立,等.基于 HHT提取重力固體潮的地震前兆信息[J].地震學報,2007,29(2):222-226.(ZhangLi,et al. HHT-based extraction of gravity tide of earthquake precursor information[J].Acta Seismologica Sinica,2007,29(2):222-226)
5 楊志華,等.一種基于 HHT的信號周期性分析方法及應用[J].中山大學學報 (自然科學版),2005,44(2):14-18.(Yang Zhihua,et al.Method and applications of periodic signal based on HHT[J].Journal of Sun Yat-sen University (science),2005,44(2):14-18)
6 Gabriel Rilling and Paulo Gon?alvés.Empirical mode decomposition as a filter bank[J].IEEE Signal ProcessingLetters,2004,11(2):112-114.
SEISM IC SIGNAL DENO ISING ALGORITHM BASED ON HHT AND ITS APPL ICATI ON
Yang Guangliang1,2)
(1)Institute of Seism ology,CEA,W uhan 430071 2)CrustalM ovem ent Laboratory,W uhan 430071)
On the basis of the Hilbert-Huang Transform(HHT)method and according to the assumption of the background noise signal continuity in time domain,an adaptive algorithm of seis mic signal denoising is presented. The algorithm,as used to analyze the seis mic wave,has certain positive results.Through the analysis following conclusions can be drawn.1)According to the background noise signal characteristics,the algorithm can effectively remove the low-frequency noise signal from the original signal.2)By the use of the adaptive algorithm we can decompose the signal and denoising,and avoid the ti me resolution and frequency resolution problem.3)The algorithm has the flaw for dealingwith high-frequency interferenc,and so needs to be improved.
EMD(EmpiricalMode Decomposition);I MF(Intrinsic Mode Function);correlation analysis;characteristics of the time extension;denoising
1671-5942(2011)04-0090-04
2011-03-17
中國地震局地震研究所重點基金(IS200916004)
楊光亮,男,博士,主要從事重力動態變化與重力殼幔結構反演研究工作.E-mail:vforyang@gmail.com
P207
A