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

在閃電電場(chǎng)信號(hào)的采集、傳輸過(guò)程中不可避免地要受到噪聲的干擾[1]，如環(huán)境噪聲、儀器噪聲、信道噪聲等，這些噪聲的存在會(huì)直接影響到對(duì)一些微小閃電放電過(guò)程的認(rèn)識(shí)，也會(huì)給后期的數(shù)據(jù)處理和分析帶來(lái)不便，所以對(duì)閃電電場(chǎng)信號(hào)的去噪處理成為了亟待解決的問(wèn)題。

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

1 經(jīng)驗(yàn)?zāi)B(tài)分解和同步壓縮小波變換

1.1 經(jīng)驗(yàn)?zāi)B(tài)分解

經(jīng)驗(yàn)?zāi)B(tài)分解可將待分析信號(hào)自適應(yīng)地分解為若干個(gè)不同時(shí)間尺度的IMF分量，對(duì)信號(hào)x(t)進(jìn)行EMD分解的過(guò)程簡(jiǎn)介如下[6]。

首先，分別采用三次樣條插值函數(shù)來(lái)擬合原始信號(hào)x(t)的上下包絡(luò)線xmax(t)、xmin(t)，通過(guò)對(duì)得到的上下包絡(luò)線求均值得到原始信號(hào)的平均包絡(luò)線m1(t)

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

(1)

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

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

(2)

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

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

(3)

式中，m11(t)為h1(t)的上下包絡(luò)線的均值。繼續(xù)對(duì)h11(t)進(jìn)行循環(huán)判斷和篩選，重復(fù)以上的步驟。若上述過(guò)程循環(huán)了k次，記為h1k(t)

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

(4)

當(dāng)h1k(t)滿足結(jié)束條件時(shí)則終止篩選。

1.2 同步壓縮小波變換

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

首先對(duì)時(shí)域信號(hào)f(t)進(jìn)行連續(xù)小波變換得到小波系數(shù)wf(a,b)

(5)

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

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

(6)

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

(7)

通過(guò)對(duì)小波系數(shù)求偏導(dǎo)可估計(jì)瞬時(shí)頻率

(8)

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

(9)

式中，ak為離散尺度，k為尺度個(gè)數(shù)。

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

(10)

式中，c-1φ表示取實(shí)部。

2 仿真信號(hào)的去噪分析

2.1 組合去噪方法流程

所提出的組合算法利用了EMD算法能夠自適應(yīng)分解信號(hào)和SST算法可將噪聲壓縮為點(diǎn)狀噪聲或顆粒狀噪聲并集中分布的特點(diǎn)，從而利用中值濾波法可以實(shí)現(xiàn)對(duì)閃電電場(chǎng)信號(hào)的去噪處理。其具體步驟如下：

① 對(duì)含噪信號(hào)進(jìn)行EMD自適應(yīng)分解，得到各階IMF分量；

② 利用原始信號(hào)與各階IMF分量的相關(guān)系數(shù)大小，挑選出IMF的優(yōu)勢(shì)分量；

③ 利用IMF優(yōu)勢(shì)分量重構(gòu)出信號(hào)；

④ 將重構(gòu)后的信號(hào)進(jìn)行SST變換；

⑤ 對(duì)SST變換時(shí)頻圖進(jìn)行中值濾波；

⑥ 用SST逆變換重構(gòu)出去噪后的信號(hào)。

算法流程圖如圖1所示。

從圖1可知，本算法的關(guān)鍵流程是對(duì)信號(hào)進(jìn)行EMD分解和SST變換。

圖1 去噪流程圖

2.2 仿真信號(hào)實(shí)驗(yàn)結(jié)果

根據(jù)M.W.Wik等人的研究，可以把閃電電磁脈沖(LEMP)歸結(jié)為雙指數(shù)衰減型脈沖波形。因此，標(biāo)準(zhǔn)閃電波的時(shí)域波形可表示為[20]

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

(11)

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

圖2 原始信號(hào)和含噪信號(hào)

圖3 4種方法去噪效果圖

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

表1 4種去噪方法的濾波結(jié)果比較

3 實(shí)測(cè)閃電電場(chǎng)信號(hào)的去噪分析

3.1 數(shù)據(jù)來(lái)源及實(shí)驗(yàn)結(jié)果

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

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

圖4 地閃輻射場(chǎng)信號(hào)

在所提組合算法中，首先利用EMD算法將經(jīng)過(guò)去均值和歸一化后的地閃輻射場(chǎng)信號(hào)(如圖4)自適應(yīng)地分解為13個(gè)IMF分量，計(jì)算出原地閃輻射場(chǎng)信號(hào)與各階IMF分量相關(guān)系數(shù)的大小，挑選出相關(guān)系數(shù)大于0.2的各階IMF分量，利用這些優(yōu)勢(shì)分量重構(gòu)出地閃信號(hào)，然后對(duì)重構(gòu)后的信號(hào)進(jìn)行SST變換，將變換后的SST時(shí)頻面進(jìn)行中值濾波，之后再進(jìn)行SST逆變換，逆變換后的結(jié)果即為所提算法的去噪結(jié)果。4種方法去噪效果如圖5所示。

圖5 4種方法去噪效果圖

從圖5中可以看出，單獨(dú)用EMD算法和單獨(dú)用SST算法去噪效果并不明顯，信號(hào)還存在著大量噪聲，這會(huì)影響到對(duì)閃電信號(hào)的時(shí)頻分析和特征提取，進(jìn)而影響對(duì)閃電電場(chǎng)信號(hào)的自動(dòng)化識(shí)別。小波閾值去噪的結(jié)果要比另外3種方法去噪效果平滑，這是因?yàn)樾〔ㄩ撝捣ㄔ谌コ肼暤耐瑫r(shí)，也將信號(hào)中的細(xì)節(jié)部分平滑了，由于閃電信號(hào)的突變和尖峰處包含了大量信息，該方法在去噪的同時(shí)也將信號(hào)中的有用成分濾除掉了，因此這種方法對(duì)后期的數(shù)據(jù)處理是很不利的，而基于EMD和SST算法在有效去噪的同時(shí)也保留了信號(hào)的細(xì)節(jié)部分，便于對(duì)閃電微小放電過(guò)程如回?fù)舻鹊恼J(rèn)識(shí)，也便于后期數(shù)據(jù)處理。因此，組合算法在對(duì)閃電電場(chǎng)信號(hào)的去噪處理方面具有優(yōu)越性。

3.2 實(shí)驗(yàn)結(jié)果定量分析

分別選取5例云閃輻射場(chǎng)信號(hào)(IC)和5例地閃輻射場(chǎng)信號(hào)(CG)，計(jì)算在4種去噪方法下的相關(guān)系數(shù)、信噪比、均方誤差。其中，相關(guān)系數(shù)計(jì)算結(jié)果如圖6所示，信噪比計(jì)算結(jié)果如表2所示，均方誤差計(jì)算結(jié)果如圖7所示。相關(guān)系數(shù)和信噪比越大，均方誤差越小說(shuō)明去噪效果越好。

圖6 4種去噪方法下的相關(guān)系數(shù)

閃電信號(hào)小波閾值法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幅圖表可發(fā)現(xiàn)，最大的信噪比和相關(guān)系數(shù)、最小的均方誤差均是本文所提方法計(jì)算出的。因此，本文算法的去噪效果要優(yōu)于文中另外3種方法的去噪效果，為閃電電場(chǎng)信號(hào)的去噪提供新思路。

4 結(jié)束語(yǔ)

本文著眼于處理非平穩(wěn)信號(hào)的EMD算法和SST算法，將其結(jié)合起來(lái)發(fā)揮其各自的優(yōu)勢(shì)，并運(yùn)用于閃電電場(chǎng)信號(hào)的去噪研究，實(shí)驗(yàn)結(jié)果表明，相比于傳統(tǒng)小波閾值去噪和單獨(dú)使用EMD算法和單獨(dú)使用SST算法去噪，本文所提算法的去噪效果更好，這對(duì)閃電電場(chǎng)信號(hào)的后期處理，例如特征提取和分類識(shí)別具有重要的意義。