許亮華，郭永剛，杜修力

（1.中國水利水電科學研究院 工程抗震研究中心，北京 100048；2.北京工業(yè)大學 建筑工程學院，北京 100124）

1 研究背景

HHT（Hilbert-Huang transform）信號分析方法［1-2］是Norden E.Huang在1998年提出的一種全新的數(shù)據(jù)處理方法，方法核心由兩部分構成，即經(jīng)驗模式分解方法（EMD）和Hilbert變換。Hilbert變換得到的瞬時頻率、振幅以及時頻分布圖等可以反應出數(shù)據(jù)的非平穩(wěn)變化規(guī)律。雖然HHT方法在非平穩(wěn)數(shù)據(jù)分析中得到廣泛應用，但其方法本身仍然存在許多問題，如包絡線擬合以及邊界問題，沒有確切公式理論支持，模態(tài)混淆等。這些問題需要得到解決才能推動HHT方法的進一步發(fā)展和應用。大壩強震反應信號是研究大壩結構特性的珍貴數(shù)據(jù)［3～5］，可以從中獲得大壩的強震下非平穩(wěn)反應規(guī)律、特性變化規(guī)律和損傷信息等有價值的信息，從而檢驗和改進大壩抗震設計，驗證大壩抗震性能，保證大壩安全。大壩強震反應數(shù)據(jù)是一種非平穩(wěn)數(shù)據(jù)，作者在應用HHT方法分析處理中發(fā)現(xiàn)，常規(guī)的EMD方法對頻譜密集的數(shù)據(jù)分解不可避免產(chǎn)生模態(tài)混淆現(xiàn)象，導致HHT分析結果與實際偏差很大，分析結果無法作為結論進行應用。因此，應用HHT方法首先必須抑制模態(tài)混淆現(xiàn)象發(fā)生。

2 常規(guī)EMD方法

信號的EMD過程［1-2］是一個將非平穩(wěn)信號進行平穩(wěn)化處理的過程，分解目的是得到有近似窄帶特征且盡量滿足單分量函數(shù)［6］特征的固有模式函數(shù)（IMF），為下一步Hibert變換提供準備。

EMD過程是一個篩選獲得IMF的過程，通過時間序列上下包絡的平均值確定瞬時平衡位置，進而提取IMF。目前學者采用的EMD方法（為了與本文改進方法名稱上區(qū)分，本文稱之為常規(guī)EMD方法）步驟為：首先搜索查找出信號x（t）的局部極大值和極小值，利用三次樣條曲線插值連接局部極大值和極小值點，分別得到極大值包絡線和極小值包絡線；其次對極大值包絡線線和極小值包絡線取平均，得到瞬時平均值m（t）；最后用原始數(shù)列x（t）減去瞬時平均值m（t），得到一個去掉低頻的新數(shù)列 h（t）。

檢查h（t）是否滿足IMF的2個條件：（1）極值點數(shù)目和過零點數(shù)目相等或最多相差1個；（2）在任意點，上下包絡線的平均值為0。

若h（t）滿足上述兩個條件，則將h（t）作為1個IMF；若不滿足，將h（t）序列按照上述步驟進行重新“篩選”，直到滿足下式

滿足式（2）的兩個條件就得到第1個IMF，記為c1（t）。

將c1（t）從原始數(shù)列中分離出來。將r1（t）作為新的信號重復上述步驟處理。經(jīng)過處理后，原始數(shù)列x（t）依次被分解出n個IMF和一個殘余項rn（t）其中殘余項rn（t）代表了原始數(shù)列的趨勢。

通過EMD分解得到的IMF在特點上非常適合作Hilbert變換，從而構造解析信號得到瞬時頻率。以上是常規(guī)EMD分解過程，具體流程見圖1。

3 基于濾波的EMD分解

在非平穩(wěn)信號中，間歇信號、噪聲信號和奇異信號存在，系統(tǒng)模態(tài)密集等因素會使EMD處理的IMF產(chǎn)生模態(tài)混淆，原因在于篩選得到的IMF并不能有效滿足單分量函數(shù)特征。在提高信噪比，剔除奇異點，以及采取手段抑制端點效應的同時，如果信號保持窄帶特性，那么EMD篩選得到的IMF必然保持窄帶特性，這樣可以盡可能避免模態(tài)混淆或降低模態(tài)混淆程度。為了解決常規(guī)EMD分解會產(chǎn)生模態(tài)混淆問題，研究人員將濾波處理和EMD方法結合起來，根據(jù)信號的傅立葉譜頻帶特征進行濾波處理，使得EMD分解時信號具有窄帶特性。

3.1 頻帶濾波與EMD方法結合應用其他學者的思路 許多方法普遍都是采取預濾波處理的思路（圖3），對信號先行帶通濾波［7-8］。首先將信號根據(jù)不同頻帶參數(shù)帶通濾波獲得幾個帶通子信號，對每個子信號分別進行EMD分解，篩選出各自的IMF；最后從所有帶通信號的IMF中提取認為有用的IMF。

這種處理過程比較繁雜，生成了許多冗余IMF，失去了EMD方法的自適應性優(yōu)勢。這種處理思路對于批量數(shù)據(jù)的信號處理是不可行的。因為每個通道信號都要選擇許多帶通濾波的參數(shù)，再加上人為篩選判斷IMF的過程，使求解效率大大降低，無法滿足工程應用的要求。同時這樣處理結果可能無法保證分解信號的完備性和正交性。這種濾波方法特點都是EMD前進行信號濾波，未觸及EMD過程。

3.2 本文改進的EMD過程 基于前面濾波思路的不足，本文對于濾波思路進行改進，將濾波處理過程和EMD過程結合起來，在EMD分解的篩選過程中對信號進行濾波處理。通過改進，這種含濾波的EMD分解方法在保證了信號窄帶特性同時發(fā)揮了EMD方法自適應分解的優(yōu)點。

3.2.1 篩選流程的改進 首先對信號r（t）進行高通濾波獲得R（t），將R（t）為EMD處理的信號，獲取到一個IMF后，獲取殘余項時采用r（t）-c1（t）=r1（t）處理。將r1（t）項作為新的r（t）重復以上過程，直到滿足篩選停止條件。改進的EMD流程見圖2。

對于高頻噪聲較強且噪聲頻帶分布較寬的信號，可以先進行低通濾波將信號噪聲濾除提高信噪比后，再采用以上方法進行處理。

3.2.2 濾波參數(shù)設置過程信號r（t）高通濾波獲得R（t）過程中濾波參數(shù)的設定可以采用以下兩種方法設置。

（1）自適應設置。這種方法的濾波參數(shù)是在EMD分解過程中自適應設置并濾波，不需要人為確定濾波參數(shù)。根據(jù)每次篩選進程中剩余成分r（t）的傅立葉譜，程序自行調(diào)整濾波參數(shù)。為便于描述，以一個仿真信號的傅立葉譜圖（見圖4）為例進行說明。過程為：①r（t）傅立葉變換，求出傅立葉譜；②傅立葉譜平滑處理，得到平滑后的傅立葉譜；③查找出傅立葉譜中一系列的局部極大值和極小值頻率點，如例圖4中所示的A、B、C和D共4個局部極大值和a、b、c、d和e共5個局部極小值；④找出局部極大值中頻率最大的兩個點C和D點，通過C、D兩點確定出處兩點間的極小值點d點。將d點的頻率值作為本次高通濾波的最小頻率參數(shù)；⑤對r（t）采用Fourier濾波器高通濾波得到R（t）。

對傅立葉譜平滑處理是為了便于查找極大值極小值，平滑程度要根據(jù)具體信號處理調(diào)節(jié)。如果信號是密集頻率信號，譜平滑的次數(shù)要減少，否則，原本相鄰的兩個譜峰可能會合并成一個譜峰，濾波時便不能區(qū)分兩個頻帶。

（2）預先設定濾波參數(shù)。信號頻譜十分密集并且信噪比不高時，傅立葉譜平滑后劃分的頻帶可能也不滿足要求，EMD結果仍然有模態(tài)混淆現(xiàn)象，這時可以采用預先設定參數(shù)的方法。過程為：①對原始信號進行傅立葉分析，獲得信號的傅立葉譜圖。根據(jù)譜圖特點先劃分頻帶，獲得每個頻帶的頻率參數(shù)，按頻率大小順序依次設置濾波參數(shù)并保存；②在每次篩選IMF前，程序依次從設置的濾波參數(shù)中提取濾波參數(shù)對r（t）進行高通濾波。

4 實例分析

4.1 仿真信號分析 為討論改進前后的分析效果，本文用最簡單的簡諧波構建仿真信號（見圖5中No.1通道）。構成的各頻率成分是15Hz的間歇簡諧波及2、5和6Hz簡諧波（分別見圖5中No.2-5通道），構成的仿真信號含間歇信號具有密集模態(tài)的特點。

分別用常規(guī)EMD方法和本文改進的EMD方法對仿真信號進行分解。比較看出，常規(guī)EMD方法分解的imf1-4波形（見圖6）與仿真信號成分No.2-5波形（見圖5）明顯不同，說明已經(jīng)出現(xiàn)模態(tài)混淆。而本文改進的EMD方法分解得到的imf1-4（見圖7），只是在端點處有輕微端點效應，間歇信號兩端產(chǎn)生過渡信號。

對得到的EMD分解結果進一步進行Hilbert分析后，得到瞬時頻率曲線。常規(guī)EMD方法得到的瞬時頻率（見圖8）波動大，特別是構成信號中5和6Hz簡諧波的頻率特征因為模態(tài)混淆，出現(xiàn)了錯誤的結果。而改進的EMD方法得到的瞬時頻率（見圖9），反應了正確的頻率分布規(guī)律。這表明改進的EMD方法有效解決了模態(tài)混淆問題。

本文瞬時頻率分布圖中每個IMF的瞬時頻率曲線的顏色深淺與相應的IMF瞬時幅值相關，瞬時幅值最大時瞬時頻率用純黑色表示，瞬時幅值等于0時瞬時頻率用純白色表現(xiàn)。

4.2 應用分析 目前在我國在水工建筑物特別是大壩上獲得的強震動反應信號還不多，每個強震反應信號都是研究建筑物地震振動特征的珍貴資料。對于這些資料需要用各種方法，各個角度去分析應用，充分體現(xiàn)出這些數(shù)據(jù)的價值。

新豐江水電站［9］位于廣東省河源縣東江支流新豐江上，距廣州市東北約160km。1962年3月19日庫區(qū)發(fā)生了震級Ms為6.1級、震中烈度為8度的強烈地震。這次地震導致大壩右岸壩段頭部108m高程出現(xiàn)長達82m的水平裂縫，地震后對大壩裂縫進行了修復加固。為了研究大壩斷裂加固效果和裂縫上下兩部分在地震作用下反應差異，設立了大壩裂縫研究臺陣，以14、16、17號墩作為觀測對象，在108m高程裂縫上、下各20cm布設拾振器。

該大壩裂縫研究臺陣，在1989年11月26日記錄到震中距為2km的4.6級地震。本文提取這次地震中14壩段108m高程縫上和下、順河向的強震反應信號（圖10）進行HHT分析舉例。

4.2.1 數(shù)據(jù)預處理 通過信號的傅立葉譜判斷該大壩在此次地震中的卓越頻率在10Hz以下，為了減少高頻影響，對數(shù)據(jù)進行HHT處理之前先進行了10Hz以下的低通濾波，增加低頻的信噪比，濾波后波形見圖11。

4.2.2 HHT分析 對預處理后的強震反應數(shù)據(jù)采用改進的EMD分解得到每個通道信號的IMF結果（見圖12和圖13），對分解結果進行Hilbert分析，計算得到各個IMF的瞬時頻率和瞬時幅值，繪制得到瞬時頻率分布（為觀測強震階段，1～5s作了放大，分別見圖14和圖15）和邊際譜（見圖16）。

EMD分解方法分解間歇信號時，由于分解是利用信號的包絡，而包絡對突變處的處理會導致分解的各個IMF結果在間歇信號開始和結束兩個時間點處產(chǎn)生過渡信號。由IMF結果中的過渡信號獲得的瞬時頻率特征并非實際信息真實存在的瞬時頻率，得到的瞬時頻率在間歇信號兩端會有較大的波動變化。結構強震反應信號在其非平穩(wěn)性表現(xiàn)之一表現(xiàn)為各個頻率成分振動幅值都增大過程和衰減過程，在幅值開始增大和幅值衰減到脈動水平時的兩個時間點，相當于間歇信號的首尾兩個時間點，得到的IMF結果同樣會有過渡信號產(chǎn)生。研究瞬時頻率變化必須針對每個頻率成分的強震動階段，在強震動階段信噪比強，得到的瞬時頻率才能夠比較準確、接近實際值，從而反應出大壩結構強震動過程中的非平穩(wěn)變化規(guī)律。

4.2.3 HHT分析結論

（1）通常由于放大效應影響，同一次地震中高程高的監(jiān)測點監(jiān)測信號的幅值要比高程低的相同方向的監(jiān)測信號的幅值大。通過裂縫上下監(jiān)測點強震信號的波形振動幅度大小可以看出，裂縫上的振幅卻比裂縫下小（其它測點同樣是這種規(guī)律），這表明修復的大壩裂縫具有一定的緩沖吸能作用。放大效應曲線在裂縫處會有所突變，并不是沿高程單純放大。

（2）大壩強震動的卓越頻率。從傅立葉譜分析裂縫上下的卓越頻率分別是5.811和5.859Hz。從邊際譜分析本次強震的大壩卓越頻率分別是5.762和5.859Hz。邊際譜和傅立葉譜分別見圖16、17。兩種譜圖的統(tǒng)計結果相近。頻率分析結果與過去觀測分析結果變化不大，這說明修復后大壩整體特性尚屬穩(wěn)定。

分析裂縫上下監(jiān)測點強震信號的邊際譜和傅立葉譜，發(fā)現(xiàn)各中心頻率略有不同，表明裂縫對大壩震動的頻率變化有一定影響。且各中心頻率的譜幅值比值不成相似的比例關系，說明裂縫對不同頻率的影響效果并不相同。具體裂縫對大壩振動的不同頻率成分影響是否存在一定的規(guī)律有待繼續(xù)深入研究。

（3）分別觀察裂縫上下強震信號的瞬時頻率分布圖，可以看出裂縫上下強震信號各頻率成分在強震階段頻率的瞬時頻率曲線沒有異常突變，基本處于水平波動變化，裂縫上下強震信號的各成分瞬時頻率分布規(guī)律相似。這些表明大壩在這次強震中處于彈性振動狀態(tài)，大壩修復后的裂縫在這次地震中處于安全狀態(tài)。

5 結語

本文對常規(guī)的EMD分解過程進行改進，在EMD分解進程中采用濾波手段，這種改進方式可以有效抑制模態(tài)混淆。通過仿真數(shù)據(jù)和實際大壩強震數(shù)據(jù)的分析表明改進的EMD方法是有效可行的，可以用于大壩強震反應數(shù)據(jù)的分析處理中。

［1］Huang N E，Shen Z，Long S R，et al.The empirical mode decomposition and the Hilbert spectrum for nonlinear and non-stationary time series analysis［J］.Proceeding of the Royal Society of Londom，1998，454：903-995.

［2］Huang N E，Shen Z，Long S R.A new view of nonlinear water waves：the Hilbert spectrum［J］.Annu.Rev.Flu?id Mech.，1999，31：417-557.

［3］DL/T5416-2009，水工建筑物強震動安全監(jiān)測技術規(guī)范［S］.

［4］蘇克忠，張力飛，等.大壩強震安全監(jiān)測［M］.北京：中國水利水電出版社，1995.

［5］郭永剛.水利水電工程強震監(jiān)測和強震監(jiān)測儀器［J］.地球物理學進展，2005（2）：422-426.

［6］L.科恩.時—頻分析：理論與應用［M］.白居憲譯，西安：西安交通大學出版社，1998.

［7］Yangj N，Lei Ying，Pan Shu-wen，et al.System identification of linear structures based on Hilbert-Huang spec?tral analysis：part 1 normal modes［J］.Earthquake Engng Struct Dyn.，2003，32（9）：1443-1467.

［8］郭淑卿.EMD的頻帶濾波篩分方法［J］.中國民航大學學報，2008，26（4）：30-33，39.

［9］陳厚群，蘇克忠，等.中國水工結構重要強震數(shù)據(jù)及分析［M］.北京：地震出版社，2000.