基于EMD和SST算法的閃電電場信號去噪研究

(1.西北師范大學 物理與電子工程學院，甘肅 蘭州 730070； 2.甘肅省智能信息技術與應用工程研究中心，甘肅 蘭州 730070； 3.國網蘭州供電公司，甘肅 蘭州 730070)

在閃電電場信號的采集、傳輸過程中不可避免地要受到噪聲的干擾[1]，如環境噪聲、儀器噪聲、信道噪聲等，這些噪聲的存在會直接影響到對一些微小閃電放電過程的認識，也會給后期的數據處理和分析帶來不便，所以對閃電電場信號的去噪處理成為了亟待解決的問題。

目前，閃電信號的去噪主要依賴于硬件濾波電路，但這種技術存在局限性，去噪效果并不理想，因此對閃電信號去噪算法的研究越來越受到人們的重視。李鵬等人應用傳統數字濾波和小波閾值法進行了閃電信號去噪的比較研究[2]。高太長等人針對回擊過程中電場變化情況，提出了一種基于多小波變換的閃電信號去噪方法[3]。劉志田等人提出了基于提升小波變換的閃電電場變化信號去噪方法[4]。火元蓮等人提出了一種基于雙密度雙樹小波變換的閃電信號去噪方法[5]。近年來，由Huang[6]等人提出的經驗模態分解算(Empirical Mode Decomposion,EMD)和由Daubechies[7]等人提出的同步壓縮小波變換算法(Synchrosqueezing Wavelet Transform,SST)在非平穩信號的處理中得到了蓬勃發展。EMD主要依據原始時間序列自身的時間尺度特征分析序列，而不必預先設定任何基函數[8]。其主要思想是它能將復雜信號分解為若干個本征模態函數(Intrinsic Mode Function,IMF)，由于分解是基于信號序列時間尺度的局部特性[9]，因此具有自適應性，相比于短時傅里葉變換、小波分解等方法其在處理非平穩信號上有明顯優勢，所以一經提出就在不同工程領域得到了廣泛的應用[10-12]。SST是基于小波變換提出的一種新的時頻分析方法，可將一維信號變換到二維空間，即可將時域信號轉化成高分辨率的時頻譜，再結合時頻譜重排的思想，通過壓縮任一中心頻率附近區間值，從而得到同步壓縮小波變換量值[13]。本文將EMD算法和SST算法結合起來，利用EMD算法能夠自適應分解信號和SST算法可將噪聲壓縮為點狀噪聲或顆粒狀噪聲并集中分布的特點[14]，從而選用中值濾波可達到對噪聲的抑制[15-16]，實現對閃電電場信號的去噪。

1 經驗模態分解和同步壓縮小波變換

1.1 經驗模態分解

經驗模態分解可將待分析信號自適應地分解為若干個不同時間尺度的IMF分量，對信號x(t)進行EMD分解的過程簡介如下[6]。

首先，分別采用三次樣條插值函數來擬合原始信號x(t)的上下包絡線xmax(t)、xmin(t)，通過對得到的上下包絡線求均值得到原始信號的平均包絡線m1(t)

m1(t)=(xmax(t)+xmin(t))/2

(1)

用x(t)減去m1(t)得到h1(t)

h1(t)=x(t)-m1(t)

(2)

通過判斷得到的h1(t) 分量是否滿足構成IMF分量的2個條件，來決定其是否為第一階IMF分量。一般情況下，第一次分解得到的h1(t)信號分量很難滿足判定條件，因此需要將h1(t)信號分量作為新的原始信號，重復上述步驟，即可以得到

h11(t)=h1(t)-m11(t)

(3)

式中，m11(t)為h1(t)的上下包絡線的均值。繼續對h11(t)進行循環判斷和篩選，重復以上的步驟。若上述過程循環了k次，記為h1k(t)

h1k(t)=h1(k-1)(t)-m1k(t)

(4)

當h1k(t)滿足結束條件時則終止篩選。

1.2 同步壓縮小波變換

同步壓縮小波變換[17-18]是在連續小波變換(CWT)的基礎上提出的一種新的時頻分析方法[19]。通過對CWT結果進行壓縮，降低了噪聲的能量，一定程度上壓制了噪聲，且SST也是一種可逆變換，即通過反變換可完全恢復原始信號，其過程如下。

首先對時域信號f(t)進行連續小波變換得到小波系數wf(a,b)

(5)

式中，a為尺度因子;b為時間因子;φ*為母小波的共軛。

令諧波信號f(t)=Acos(wt)，根據Plancherel定理，式(5)的頻率域變換為

(6)

將f(t)的傅里葉變換代入式(6)得

(7)

通過對小波系數求偏導可估計瞬時頻率

(8)

通過式(8)，將小波系數wf(a,b)從時間-尺度平面映射到時間-頻率平面，即wf?wf(a,b),b」。

(9)

式中，ak為離散尺度，k為尺度個數。

同步壓縮小波變換的反變換為

(10)

式中，c-1φ表示取實部。

2 仿真信號的去噪分析

2.1 組合去噪方法流程

所提出的組合算法利用了EMD算法能夠自適應分解信號和SST算法可將噪聲壓縮為點狀噪聲或顆粒狀噪聲并集中分布的特點，從而利用中值濾波法可以實現對閃電電場信號的去噪處理。其具體步驟如下：

① 對含噪信號進行EMD自適應分解，得到各階IMF分量；

② 利用原始信號與各階IMF分量的相關系數大小，挑選出IMF的優勢分量；

③ 利用IMF優勢分量重構出信號；

④ 將重構后的信號進行SST變換；

⑤ 對SST變換時頻圖進行中值濾波；

⑥ 用SST逆變換重構出去噪后的信號。

算法流程圖如圖1所示。

從圖1可知，本算法的關鍵流程是對信號進行EMD分解和SST變換。

圖1 去噪流程圖

2.2 仿真信號實驗結果

根據M.W.Wik等人的研究，可以把閃電電磁脈沖(LEMP)歸結為雙指數衰減型脈沖波形。因此，標準閃電波的時域波形可表示為[20]

E(t)=E0(e-αt-e-βt)

(11)

式中，E0為脈沖波形的幅值系數；α、β為波前衰減系數和波尾衰減系數。在仿真過程中E0=30 V/m，α=2.0×107s-1，β=2.0×106s-1，采樣頻率fs=60 MHz。在雙指數衰減脈沖(如圖2(a))上疊加一信噪比(SNR)為30 dB 的高斯白噪聲(如圖2(b))。分別利用小波閾值法(小波函數取db5小波，分解層數取為6)、EMD算法、SST算法和所提出的組合方法對含噪信號進行消噪處理。在組合算法中，先用EMD算法自適應地將含噪信號分解為7個IMF分量，計算出各階IMF分量與原始信號的相關系數大小，選擇相關系數大于0.3的各階IMF分量重構出信號(EMD算法去噪結果)，然后將重構信號進行SST變換，將變換后的時頻圖進行中值濾波，之后進行SST逆變換，逆變換后的結果即為本文去噪方法的去噪結果。去噪結果如圖3所示，圖中依次為小波閾值去噪、EMD去噪、SST去噪、組合方法去噪的結果。

圖2 原始信號和含噪信號

圖3 4種方法去噪效果圖

為了定量地比較組合方法和另外3種方法的去噪效果，分別計算了信噪比(SNR)、相關系數(CC)和均方誤差(MSE)。這3種參數均能從不同角度來度量去噪效果，其中信噪比和相關系數越大，均方誤差越小，說明去噪效果越好。計算結果如表1所示，結果表明基于EMD和SST算法的信噪比和相關系數比其余3種方法要大，均方誤差較另外3種方法要小，說明組合去噪方法去噪效果優于其他3種算法的去噪效果。

表1 4種去噪方法的濾波結果比較

3 實測閃電電場信號的去噪分析

3.1 數據來源及實驗結果

本文選用的實驗數據是2009年在青海大通地區用快電場變化測量儀記錄的觀測資料。大通地區海拔較高，屬于雷電多發區，每個測站安裝有用于閃電輻射脈沖三維定位的閃電VHF輻射源到達時間差(TOA)定位系統和GPS同步的高精度時鐘(±25 ns) ，快電場變化探測儀帶寬為100 Hz～5 MHz，時間常數為1 ms、采樣速率為2.5 MHz。

在對閃電電場信號進行去噪之前，先進行去均值和歸一化處理。去均值即對閃電電場信號進行零均值化處理。由于采集到的閃電電場信號的距離大小不同，閃電放電強度不同，為了減少數據的分散性，對閃電電場信號進行歸一化，以得到幅值范圍統一的信號。考慮到EMD算法需要在整個信號長度范圍內作樣條插值，當信號采樣點數較多時，特別是極值點數目多的情況下，采用三次樣條擬合法計算量特別大，而采集到的一維閃電電場信號記錄時間超過800 ms，包含了2096000個采樣點，因此，直接對閃電信號進行EMD分解是很難實現的。為了提高EMD算法的計算速率，首先對閃電電場信號以1∶1000進行了重采樣，對重采樣后的信號分別用小波閾值法、EMD算法、SST算法和組合算法對地閃個例信號(如圖4所示)進行了去噪處理。

圖4 地閃輻射場信號

在所提組合算法中，首先利用EMD算法將經過去均值和歸一化后的地閃輻射場信號(如圖4)自適應地分解為13個IMF分量，計算出原地閃輻射場信號與各階IMF分量相關系數的大小，挑選出相關系數大于0.2的各階IMF分量，利用這些優勢分量重構出地閃信號，然后對重構后的信號進行SST變換，將變換后的SST時頻面進行中值濾波，之后再進行SST逆變換，逆變換后的結果即為所提算法的去噪結果。4種方法去噪效果如圖5所示。

圖5 4種方法去噪效果圖

從圖5中可以看出，單獨用EMD算法和單獨用SST算法去噪效果并不明顯，信號還存在著大量噪聲，這會影響到對閃電信號的時頻分析和特征提取，進而影響對閃電電場信號的自動化識別。小波閾值去噪的結果要比另外3種方法去噪效果平滑，這是因為小波閾值法在去除噪聲的同時，也將信號中的細節部分平滑了，由于閃電信號的突變和尖峰處包含了大量信息，該方法在去噪的同時也將信號中的有用成分濾除掉了，因此這種方法對后期的數據處理是很不利的，而基于EMD和SST算法在有效去噪的同時也保留了信號的細節部分，便于對閃電微小放電過程如回擊等的認識，也便于后期數據處理。因此，組合算法在對閃電電場信號的去噪處理方面具有優越性。

3.2 實驗結果定量分析

分別選取5例云閃輻射場信號(IC)和5例地閃輻射場信號(CG)，計算在4種去噪方法下的相關系數、信噪比、均方誤差。其中，相關系數計算結果如圖6所示，信噪比計算結果如表2所示，均方誤差計算結果如圖7所示。相關系數和信噪比越大，均方誤差越小說明去噪效果越好。

圖6 4種去噪方法下的相關系數

閃電信號小波閾值法EMD算法SST算法本文算法GC115.123.420.226.8GC217.320.419.223.5GC314.716.915.421.3GC419.521.620.126.9GC520.123.521.228.2IC119.123.421.626.3IC220.423.522.928.4IC313.615.914.820.3IC49.413.511.218.9IC52.93.13.08.4

圖7 4種去噪方法下的均方誤差

觀察這3幅圖表可發現，最大的信噪比和相關系數、最小的均方誤差均是本文所提方法計算出的。因此，本文算法的去噪效果要優于文中另外3種方法的去噪效果，為閃電電場信號的去噪提供新思路。

4 結束語

本文著眼于處理非平穩信號的EMD算法和SST算法，將其結合起來發揮其各自的優勢，并運用于閃電電場信號的去噪研究，實驗結果表明，相比于傳統小波閾值去噪和單獨使用EMD算法和單獨使用SST算法去噪，本文所提算法的去噪效果更好，這對閃電電場信號的后期處理，例如特征提取和分類識別具有重要的意義。