余 海，張廷龍，劉文杰，潘家利，陳 潔

(1.海南省南海氣象防災(zāi)減災(zāi)重點實驗室，海南 ?？?570203；2.海南省氣象災(zāi)害防御技術(shù)中心，海南 ?？?570203；3.海南大學(xué) 生態(tài)與環(huán)境學(xué)院，海南 ?？?70228)

大氣電場儀常被用來監(jiān)測晴天大氣電場和雷暴活動引起的電場，晴天電場強度通常為100 V·m-1左右，而雷暴電場要遠(yuǎn)大于晴天電場.因此，雷暴電場通常被認(rèn)為是疊加于晴天電場上的擾動[1-3].但雷暴擾動電場顯著高于背景的晴天電場，在具體討論雷暴電場時，通常不考慮晴天電場的貢獻(xiàn).組網(wǎng)大氣電場儀由于探測范圍更大、多站信號可綜合分析等優(yōu)勢，常常用來對雷暴過程或雷電的監(jiān)測和預(yù)警.隨著雷暴過程的靠近和遠(yuǎn)離，地面電場不僅會出現(xiàn)強度增強和減弱，還會發(fā)生電場極性的反轉(zhuǎn)，另外，閃電發(fā)生時，地面電場表現(xiàn)為快速的跳變，這些信息在一定程度上都可作為預(yù)警信息參量.受大氣電場儀自身的靈敏度和安裝環(huán)境等因素的影響，實測的大氣電場信號中常夾雜著一些其他非雷暴電場信號，也就是噪聲信號.如外場架設(shè)的大氣電場儀易受周邊植被的影響，當(dāng)雷暴過境時，隨風(fēng)擺動的植被給電場信號造成不同程度的信號“污染”.此外，還有來自工頻噪聲、地物反射噪聲以及其他輻射源信號等的干擾，這給后期資料的分析及應(yīng)用造成很大的困難.在具體應(yīng)用大氣平均電場時，需要對原始信號進(jìn)行去噪處理，保留雷暴大氣電場的變化趨勢和閃電引起的電場變化.目前，小波去噪應(yīng)用較為廣泛，其原理是首先將原始信號進(jìn)行小波分解，然后對分解得到的系數(shù)選擇一個合適的閾值，最后對降噪處理后的系數(shù)進(jìn)行小波逆變換，得到降噪后的信號.茍學(xué)強[4]等對地面電場資料的分析發(fā)現(xiàn)，閃電引起的較大的電場變化具有強烈的奇異性和明顯的多重分形性.余蓉[5]等以雙指數(shù)衰減脈沖信號上疊加高斯白噪聲來仿真大氣電場信號，并計算去噪前后仿真信號與實測大氣電場信號的均方誤差，對比分析不同小波函數(shù)、閾值分析方法和小波分解層數(shù)的去噪效果.李鵬[6]等對比了加權(quán)平均濾波、FIR數(shù)字低筒濾波器和小波閾值去噪方法對閃電瞬態(tài)電場信號的去噪效果，發(fā)現(xiàn)小波閾值去噪方法優(yōu)于其他2種方法，并且確定小波分解層數(shù)應(yīng)為5～7層.李艷[7]等、陳紅兵[8]等、李振亞[9]等都利用了小波技術(shù)對地面大氣電場數(shù)據(jù)進(jìn)行濾波或去噪處理.

雷暴云產(chǎn)生的閃電可以分為云閃和地閃，而云閃所占比例約為總閃電的2/3[3].由于地閃活動給社會生產(chǎn)生活造成的危害更大，目前多省份建立的閃電定位系統(tǒng)都為地閃監(jiān)測系統(tǒng)，而對雷暴進(jìn)行全閃電(云閃和地閃)監(jiān)測的較少.對于雷暴而言，云閃通常要早于地閃發(fā)生且云閃的頻率要高于地閃[3].因此，在缺少云閃觀測的前提下，利用地面電場信號提取總閃電頻數(shù)特征有利于雷電的臨近預(yù)警.筆者將采用小波分析方法開展地面電場觀測資料去噪處理，并在此基礎(chǔ)上，研究地面大氣電場波形中閃電信號包括閃電頻次及極性等提取方法，實現(xiàn)對雷暴總閃電頻數(shù)的評估.目前，國內(nèi)利用大氣電場儀進(jìn)行雷電預(yù)警業(yè)務(wù)應(yīng)用時[10]，一般采用設(shè)定多級電場強度閾值來劃分雷電預(yù)警等級，在實際業(yè)務(wù)化運行中容易出現(xiàn)漏報及誤報.利用小波技術(shù)處理地面大氣電場數(shù)據(jù)并提取相關(guān)閃電信息，結(jié)合地閃定位網(wǎng)，可以得到雷暴的云閃的時間演變特征，以此表征雷暴在不同發(fā)展階段的閃電活躍度，較為準(zhǔn)確地對地閃活動進(jìn)行預(yù)警.海南島對流活動頻繁[11]，閃電活動頻發(fā)[12-13]，雷電預(yù)警預(yù)報技術(shù)的提升十分必要.海南島目前已布設(shè)50多部地面大氣電場儀，以期通過本研究能夠為組網(wǎng)電場儀觀測業(yè)務(wù)應(yīng)用及雷電預(yù)警預(yù)報提供的關(guān)鍵技術(shù)和指標(biāo).

1 資料與方法

1.1 資料本文分析的資料為2013年和2014年6～9月自行布設(shè)的大氣電場儀獲取的雷暴電場部分個例.該大氣電場儀采樣頻率為125 Hz，電場強度探測范圍為±100 kV·m-1，探測范圍為10～15 km.其中，電場的正、負(fù)極性定義為：當(dāng)?shù)孛骐妶鍪芾妆┰苾?nèi)的正電荷控制時，地面電場為正，反之為負(fù).

1.2 小波分析方法小波分析方法用于電場信號的去噪原理，對原始時間序列的電場信號作小波變換，通過伸縮平移運算對信號逐步進(jìn)行多尺度細(xì)化，達(dá)到高頻處時間細(xì)分，低頻處頻率細(xì)分，而噪聲信號一般處于高頻處，而電場信號處于低頻處.通過多尺度的頻率細(xì)化，基本將噪聲細(xì)分至高頻部分，電場信號細(xì)分至低頻部分，然后選擇合適的閾值對各部分進(jìn)行閾值處理，完成電場信號與噪聲信號的區(qū)分.其具體過程可分為3個部分:1)分解過程.選定一種小波，對原始信號進(jìn)行指定層數(shù)的小波分解，分解得到低頻近似信號和高頻細(xì)節(jié)信號.2)閾值過程.對分解得到的各層系數(shù)選擇一個合適閾值，并對細(xì)節(jié)系數(shù)作閾值處理.3)重建過程.將降噪處理后的系數(shù)進(jìn)行小波逆變換，重建恢復(fù)原始信號[14].

1.2.1 小波分解層數(shù)對大氣電場信號進(jìn)行小波分析時，首先應(yīng)確定小波分解的層數(shù)，以便合理地對大氣電場數(shù)據(jù)波形進(jìn)行分析.計算電場信號的功率譜，能夠更好地找出電場信號的能量主要集中在哪一頻段范圍內(nèi).具體方法是對電場信號進(jìn)行快速傅里葉變換，然后取其幅值的平方除以信號樣本長度，其值即為該信號樣本的功率譜.選取2014年8月的一次雷暴地面大氣電場信號為樣本，經(jīng)計算得出如圖1所示的功率譜圖.

圖1 大氣地面電場數(shù)據(jù)序列功率譜

圖1所示，該信號功率譜能量主要集中在0.06 Hz以下，處于非常低的頻段范圍內(nèi).小波分解是將原始信號分解成一低頻近似信號和一高頻細(xì)節(jié)信號，多層分解的過程就是不斷的將低頻近似信號再分解成一低頻近似信號和一高頻細(xì)節(jié)信號.本文中的電場信號采集頻率為125 Hz，根據(jù)采樣定理可知小波分解的頻段區(qū)間，通過7層小波分解后低頻近似信號可以逼近至0.5 Hz以下并結(jié)合多層分解中不同層低頻近似與原始信號的對比情況發(fā)現(xiàn)，第七層低頻近似信號與原始信號變化特征基本相似.因此，擬對原始電場信號作7層小波分解，以期提取其電場波形的主要變化特征.

1.2.2 小波函數(shù)選取考慮到計算的功率譜以及選定的分解層數(shù)情況，在選取小波函數(shù)方面采取了較為簡單的對比方法，即選取不同小波函數(shù)，對比分析7層分解后的第七層低頻近似信號與原始電場信號的異同性，以此來驗證不同小波函數(shù)在分析電場信號方面的去噪效果.

隨機選取了7種小波函數(shù)(db5，db10，sym5，coif4，bior3.5，bior4.4和haar)進(jìn)行對比實驗.為了分析函數(shù)的處理效果，引入均方誤差(MSE)以及信噪比(SNR)參數(shù)來衡量電場信號處理前后的誤差及效果，MSE及SNR計算公式

(1)

(2)

其中，N為信號樣本數(shù)，s1和s2為小波處理前后的信號，k表示信號樣本中第k個樣本.MSE值越小，表示分解后的低頻近似信號與原始信號相似性越高；SNR值越大，效果越好.

表1給出了7種小波的信噪比、均方誤差計算值.由表1可以看出，均方誤差最小的小波函數(shù)和信噪比最大的小波函數(shù)均為db10小波函數(shù)，表明原始電場信號經(jīng)db10小波處理后，信號均方誤差較小，數(shù)據(jù)的離散程度較小，數(shù)據(jù)的可靠性較大，另外信噪比值最大，系統(tǒng)的隨機噪聲處理較好.因此，地面大氣電場信號分析可以采用db10函數(shù).

表1 選取的7種小波的信噪比、均方誤差計算值

1.2.3 閾值函數(shù)選取為進(jìn)一步探討選用函數(shù)的適用性，選用5種信噪比較大和均方誤差較小的小波函數(shù)(db10，db5，sym5，coif4和bior 4.4)和4種閾值分析法(Rigrsure，Sqtwolog，Heursure和Minimaxi)對各自的去噪效果進(jìn)行對比分析.

Rigrsure閾值法是一種基于stein的無偏似然估計原理的自適應(yīng)閾值選擇.對一個給定的閾值t，得到似然估計，再將非似然t最小化，得到了所選的閾值.Sqtwolog閾值法采用的是固定的閾值形式，產(chǎn)生的閾值大小是信號長度的對數(shù)的開方.Heursure閾值法是Rigrsure閾值法和Sqtwolog閾值法的綜合，如果信噪比很小似然估計有很大的噪聲，出現(xiàn)此情況，就采用固定的閾值.Minimaxi閾值法采用的也是一種固定的閾值，產(chǎn)生一個最小均方誤差的極值，而不是無誤差[14].

圖2 不同小波函數(shù)與閾值分析方法去噪后電場信號均方誤差

由圖2可知，在4 種閾值分析法中，Rigrsure閾值分析法的均方誤差是最小的，即降噪效果是最好的，其次是Minimaxi閾值法.Rigrsure閾值分析法下小波變換的信號均方誤差最小是sym5小波函數(shù)，其次是db10小波函數(shù).綜上所述，采用db10和sym5小波函數(shù)都能取得較好的去噪效果，但結(jié)合閾值分析方法，sym5小波配合Rigrsure閾值方法，去噪效果更優(yōu).

2 實例結(jié)果分析

選用sym5函數(shù)作為用于小波變換的小波函數(shù)，并采用Rigrsure閾值分析方法實現(xiàn)地面電場信號的去噪.此次選取2014年8月20日的地面電場數(shù)據(jù)，該過程閃電活動較為頻繁，且信號上疊加有一定的噪聲(圖3a)，采用上述去噪方法進(jìn)行測試，得到去噪后的波形(圖3b).

硬件實力支撐的同時，新科空調(diào)也沒有忽視軟實力的提升。舒樂華透露，在2019年，新科空調(diào)將繼續(xù)堅持“人才”、“創(chuàng)新”和“品牌”三大戰(zhàn)略，堅持“渠道為王，客戶為尊”的經(jīng)營理念，不斷開拓進(jìn)取。

圖3 地面電場波形圖

由圖3可見，電場信號在進(jìn)行去噪處理后，波形重疊部分明顯減少，波形更為平滑且突出顯示了閃電活動的主要變化趨勢.表明采用該方法可以較為顯著的降低信號中的噪聲，為數(shù)據(jù)的后期分析和處理起到質(zhì)量控制的作用.該雷暴過程地面電場波形從樣本從1.6×105序列點開始振蕩，電場強度急劇變化的時間大約持續(xù)了33 min.地面電場急劇變化過程中，電場強度主要集中在-20～25 kV·m-1，這一時段內(nèi)閃電活動頻發(fā)，電場波形主體為正極性，閃電引起的電場變化主要為負(fù)極性.

為了進(jìn)一步驗證sym5 小波函數(shù)在Rigrsure閾值方法下對地面大氣電場數(shù)據(jù)的處理效果，還選取了3 次伴有閃電活動的地面電場數(shù)據(jù)進(jìn)行小波去噪處理,如圖4所見.

圖4 地面電場原始波形及去噪后的電場波形

由圖4可見，3個樣本個例原始波形的噪聲相對較為明顯，經(jīng)去噪處理后，波形的重疊度明顯降低，即噪聲信號顯著得到抑制，電場變化更為清晰.對于有閃電活動的電場信號而言，經(jīng)sym5小波函數(shù)和Rigrsure閾值方法處理后，不僅能較好地保持原信號中電場波形強度陡增或者電場極性突變的部分，且能夠較為清晰地辨識出不同極性閃電放電次數(shù)等信息，對于后續(xù)開展大氣電場相關(guān)研究提供了有利的條件.

上述觀測實例的分析結(jié)果表明，采用Rigrsure閾值分析方法和sym5小波函數(shù)對地面大氣電場波形進(jìn)行去噪處理，在保留閃電波動信息的同時對信號也進(jìn)行了較好的平滑處理.

3 閃電信息提取方法

大氣電場信號去噪是為了更為精確地分析該信號的電場特征，并在分析該信號特征的基礎(chǔ)上，分析其代表的物理意義及閃電活動特征.地面電場波形的快速跳變能夠指示當(dāng)前時段閃電活動的規(guī)律，利用小波分析方法，對地面電場波形進(jìn)行去噪處理并分析其大氣電學(xué)層面的特征，包括電場強度、閃電頻數(shù)以及極性等閃電特征.

一次完整的閃電活動(地閃或者云閃)持續(xù)時間從幾十毫秒至1秒左右[3]，尤其是電場快速跳變的地閃回?fù)舴烹姾驮苾?nèi)電荷擊穿過程時間更為短暫(微秒級).另外文中電場數(shù)據(jù)的采樣率和數(shù)據(jù)記錄存儲速率分別為125 Hz和10 Hz(即1 s內(nèi)存儲采樣數(shù)據(jù)中的10個電場值)，意味著記錄的電場波形中快速跳變部分可以表示為一次閃電.快速跳變的部分在波形上表現(xiàn)為脈沖形式，識別這些脈沖需要計算電場波形中峰值和谷值點，然后求解該峰值或谷值點附近波形的上升沿或下降沿的斜率k.因此，本文采用差分計算方法來識別電場波形中的跳變部分，具體公式

ΔE=Ei+1-Ei，

(3)

其中，ΔE為相鄰電場數(shù)據(jù)的差值.

根據(jù)式(3)對整列電場樣本作差分計算，得出電場樣本的差分值.根據(jù)已識別的電場波形中峰值和谷值點，分析該點前后各5個差分值，計算該段上升沿或下降沿的斜率k(k=ΔE/Δt)，根據(jù)k值的變化情況確定該點是否是閃電放電引起的電場變化，如果是，則記為1次閃電，并同時記錄極性.以此類推，識別出整列電場樣本中閃電頻次及相應(yīng)極性.

在自動識別閃電的算法中，判斷k值的變化是其中的重點，對應(yīng)到式(3)中即差分值ΔE的取值.ΔE的絕對值的取值大于某個閾值時，是判斷該脈沖是否為閃電的關(guān)鍵.一次閃電持續(xù)時間大約在1 s以內(nèi)，而考慮到本文中電場數(shù)據(jù)的記錄頻率為10 Hz(即1 s記錄10個樣本點)，電場波形序列相鄰差分值的顯著差異可以表征為閃電脈沖.以2014年8月20日的地面大氣電場資料為例，詳細(xì)分解了電場閾值的取值方法.樣本序列在1.6×105開始振蕩，電場強度急劇變化的時間大約持續(xù)了33 min，在此期間電場強度在-30～+30 kV·m-1.

圖5 2014年8月20日電場波形電場快速反轉(zhuǎn)變化部分(紅色圓點為閃電脈沖尖峰或脈沖起始點)

圖5是當(dāng)天雷暴天氣下的大氣電場波形(晴天電場部分未在圖中給出).從圖5中可以看出，地面電場主要受正電荷區(qū)控制，當(dāng)電荷積聚到一定程度時，閃電通道被擊穿，通道內(nèi)正負(fù)電荷瞬間中和，引起地面電場波形出現(xiàn)極性反轉(zhuǎn)變化.電場波形中出現(xiàn)了脈沖尖峰，共計41個尖峰，其中電場負(fù)極性變化的有38個，正極性的3個.通過分析41個尖峰持續(xù)時間、上升沿和下降沿斜率推斷電場儀在其探測范圍內(nèi)采集了41次閃電的信息，閃電頻數(shù)為1.2 fl·min-1.此處正、負(fù)極性僅僅是針對波形中脈沖尖峰的變化方向定義，與正負(fù)極性閃電的定義無關(guān).本文中為與自動識別閃電方法區(qū)分開，將這種通過觀察電場波形、逐一分析尖峰脈沖持續(xù)時間與上升沿斜率的方法稱為人工識別方法.

設(shè)計自動識別電場波形中閃電頻數(shù)算法時，對比計算了不同電場判斷閾值下，算法自動識別的結(jié)果與人工識別的結(jié)果的誤差.

圖6 不同閾值下自動識別的總閃電頻數(shù)與人工識別的閃電頻數(shù)間誤差

圖6給出了不同閾值下算法自動識別出與人工識別的總閃電頻數(shù)的誤差計算結(jié)果.由圖6可知，不同閾值下，兩者的誤差計算結(jié)果先減小，至ΔE=1.5 kV·m-1時，誤差降至最低.隨著閾值的增大誤差隨之增大，表明電場判斷閾值越小，算法識別的脈沖就越多，會導(dǎo)致部分非閃電脈沖的電場波動也被識別為閃電脈沖，從而導(dǎo)致誤差增大；電場判斷閾值越大，算法識別的脈沖就越少，導(dǎo)致部分電場幅值較小的閃電脈沖被過濾，從而導(dǎo)致誤差增大.僅考慮總閃電頻數(shù)的誤差計算結(jié)果，電場閾值ΔE=1.5 kV·m-1時，算法自動識別閃電頻數(shù)的誤差最小.

為確保閾值選擇的合理性，還計算了人工識別與不同閾值下算法自動識別下每5 min時段內(nèi)閃電頻數(shù)的統(tǒng)計結(jié)果，如圖7所示.

圖7 每5 min內(nèi)人工識別閃電頻數(shù)與不同閾值下閃電頻數(shù)對比情況

圖7為有閃電過程的電場波形區(qū)段(晴天電場部分波形未顯示)人工識別與不同閾值(ΔE分別為1.0,1.3,1.5，1.7,2.0,2.5,3.0,4.0,5.0,6.0)下算法自動識別的結(jié)果.每5 min的閃電頻數(shù)的曲線總體變化特征相似，但在具體時段閃電頻數(shù)有所不同.為進(jìn)一步查看具體的差異，對不同時段的閃電頻數(shù)與人工識別的閃電頻數(shù)做了相關(guān)性計算.根據(jù)圖7的統(tǒng)計結(jié)果和選取的10組電場閾值ΔE，分別計算不同閾值下每5 min時間序列下自動識別結(jié)果與人工識別結(jié)果的相關(guān)系數(shù)，計算結(jié)果如表2所示.

表2 每5 min內(nèi)不同閾值下算法自動識別的閃電頻數(shù)與人工識別閃電頻數(shù)相關(guān)性計算結(jié)果

由表2可知，每5 min內(nèi)不同閾值下算法自動識別的閃電頻數(shù)與人工識別閃電頻數(shù)相關(guān)性計算結(jié)果整體上呈現(xiàn)閾值越大、相關(guān)性越差的規(guī)律特征，且閾值越小，算法識別的閃電頻數(shù)相對偏高(ΔE分別為1.0和1.3時)，閾值越大，部分時段漏掉的閃電過程可能越多(ΔE分別為5.0和6.0時).因此，在選定閾值時，應(yīng)當(dāng)綜合考慮總閃電頻數(shù)誤差以及分時段內(nèi)兩者間的相關(guān)性結(jié)果.在本文中綜合兩者計算結(jié)果，選定電場閾值為1.5 kV·m-1.

此自動識別方法在實際應(yīng)用中，應(yīng)當(dāng)對地面電場數(shù)據(jù)進(jìn)行多樣本計算，以此來確定此類地面電場儀獲取的電場波形的電場判斷閾值.由于云閃先于地閃發(fā)生[3]，通過上述方法自動識別電場波形中的總閃電信息，再結(jié)合地閃定位監(jiān)測網(wǎng)的數(shù)據(jù)，能夠獲取地閃發(fā)生前的閃電活動情況(主要是閃電發(fā)生時間及頻次)，這將大大提前雷電預(yù)警的時效.

4 結(jié)束語

利用小波分析技術(shù)研究了地面大氣電場資料的去噪技術(shù)以及提取去噪電場信號中閃電頻次、極性等信息的方法，分析方法與結(jié)果能夠在分析不同天氣系統(tǒng)下地面大氣電場特征以及獲取雷暴天氣系統(tǒng)下閃電活動信息等方面起到關(guān)鍵性作用，為后期開展海南島電場儀組網(wǎng)觀測及基于此網(wǎng)絡(luò)的雷電預(yù)警技術(shù)研發(fā)奠定了良好的技術(shù)基礎(chǔ).

通過算法自動識別并提取電場波形中的閃電頻數(shù)等信息時，電場儀數(shù)據(jù)記錄存儲頻率的高低以及電場強度幅值大小都會直接影響閃電頻數(shù)識別誤差.此外，電場儀的探測范圍有限，部分發(fā)生距離較遠(yuǎn)的閃電或雷暴云靠近或遠(yuǎn)離電場儀站點時引起的幅值較小的脈沖波形也會導(dǎo)致算法識別的誤差增大.下一步工作將著重研究多部大氣電場儀組網(wǎng)觀測結(jié)果.綜合考慮不同電場儀站點環(huán)境及地形等因素，利用本文方法對電場儀觀測網(wǎng)內(nèi)不同電場資料進(jìn)行小波分析，分析不同測站的差分閾值，提取同一天氣系統(tǒng)下不同電場儀識別的閃電頻數(shù)等信息，綜合多站結(jié)果確定該天氣系統(tǒng)下總閃電頻數(shù)等信息，為雷電預(yù)警提供關(guān)鍵參數(shù)及提前量.