結(jié)合EMD和功率譜熵的船舶軸頻電場(chǎng)線譜提取

程 銳，陳 聰，姜潤(rùn)翔

(1. 海軍工程大學(xué) 兵器工程系，湖北 武漢 430033；2. 海軍工程大學(xué) 理學(xué)院，湖北 武漢 430033)

結(jié)合EMD和功率譜熵的船舶軸頻電場(chǎng)線譜提取

程 銳1，陳 聰2，姜潤(rùn)翔1

(1. 海軍工程大學(xué) 兵器工程系，湖北 武漢 430033；2. 海軍工程大學(xué) 理學(xué)院，湖北 武漢 430033)

為實(shí)現(xiàn)強(qiáng)海洋背景噪聲中的微弱船舶軸頻電場(chǎng)信號(hào)檢測(cè)，提出了一種結(jié)合經(jīng)驗(yàn)?zāi)B(tài)分解（Empirical Mode Decomposition，EMD）和窄帶子區(qū)間功率譜熵的線譜提取新算法。首先，利用EMD方法從含噪信號(hào)中分解出一組有效固有模態(tài)函數(shù)（Intrinsic Mode Functions，IMFs），對(duì)各有效IMF的功率譜進(jìn)行子區(qū)間劃分；其次，定義并計(jì)算各子區(qū)間的能量峰值熵比（Energy Peak Entropy Ratio，EPER）特征；最后，通過對(duì)軸頻信號(hào)和環(huán)境噪聲物理特征差異的分析，結(jié)合K-均值聚類方法進(jìn)行特征量的篩選，實(shí)現(xiàn)線譜提取。海上實(shí)測(cè)數(shù)據(jù)的處理結(jié)果表明，相比于直接的功率譜分析，算法的線譜可提取下限降低了6.7 dB。

軸頻電場(chǎng)；經(jīng)驗(yàn)?zāi)B(tài)分解（EMD）；K-均值聚類；能量熵；線譜提取

0 引 言

不同金屬材料構(gòu)成的船舶在海水中，會(huì)在其周圍產(chǎn)生電場(chǎng)，其中由轉(zhuǎn)動(dòng)的螺旋槳對(duì)腐蝕和防腐電流的調(diào)制而激發(fā)出的軸頻電場(chǎng)[1]，具有隨距離衰減慢、線譜特征明顯的特點(diǎn)，通常基頻約為1～7 Hz，海洋環(huán)境電場(chǎng)噪聲在此頻段相對(duì)量級(jí)較低，因此，軸頻電場(chǎng)構(gòu)成了一種對(duì)船舶目標(biāo)的有效探測(cè)源和潛在攻擊源[2–3]。目前隨著軸頻電場(chǎng)抑制技術(shù)[4–5]的發(fā)展，低信噪比下的特征線譜提取方法，已成為軸頻電場(chǎng)應(yīng)用研究的一個(gè)熱點(diǎn)。

根據(jù)信號(hào)與噪聲在特征頻段內(nèi)的能量差異，文獻(xiàn)[6]提出了一種基于滑動(dòng)功率譜估計(jì)的檢測(cè)算法，本質(zhì)上是頻域的能量檢測(cè)，隨噪聲強(qiáng)度增加，算法性能下降迅速；文獻(xiàn)[7-9]基于小波良好的時(shí)頻局域化能力，在小波域中利用小波模極大值、小波尺度相關(guān)、小波包熵等特征進(jìn)行信號(hào)檢測(cè)，但小波類算法存在基函數(shù)選擇缺乏自適應(yīng)性、降噪閾值設(shè)定不易的缺陷；文獻(xiàn)[10]將經(jīng)驗(yàn)?zāi)B(tài)分解（Empirical Mode Decomposition，EMD）算法與4階累積量相結(jié)合，利用前者的信號(hào)自適應(yīng)分解能力和后者的抑制高斯過程的性質(zhì)，取得了不錯(cuò)的效果，但算法中部分參數(shù)需要多次試驗(yàn)才能確定。針對(duì)上述問題，本文利用軸頻信號(hào)和環(huán)境噪聲功率譜窄帶分布特點(diǎn)上的差異，基于EMD和定義的能量峰值熵比（Energy Peak Entropy Ratio，EPER）特征，結(jié)合K-均值聚類方法，實(shí)現(xiàn)了低信噪比下的軸頻電場(chǎng)信號(hào)特征頻率子區(qū)間的自動(dòng)確定。

1 基本方法簡(jiǎn)述

1.1 EMD

EMD方法[11]是一種高效的數(shù)據(jù)自適應(yīng)信號(hào)分解算法，不同于采用三角函數(shù)族的傅里葉分析和采用各種小波基函數(shù)的小波類算法，其基函數(shù)是采用逐級(jí)尺度篩選的方式，從數(shù)據(jù)本身后驗(yàn)地構(gòu)造出來，因此，本質(zhì)上來說EMD方法具有自適應(yīng)性。

EMD分解在頻域中具有二進(jìn)尺度分解特性[12]，不同階的IMF反映了不同的頻率分量。EMD方法的具體實(shí)現(xiàn)步驟參見文獻(xiàn)[11]。

1.2 功率譜熵

功率譜熵用于衡量功率譜中各頻率分量分布的集中程度。分布集中于少數(shù)頻點(diǎn)的功率譜，其功率譜熵越小，反之各分量的強(qiáng)度分布越均勻，則熵值會(huì)越大[13]。原始功率譜熵定義如下：

1.3 K-均值聚類

K-均值聚類是一種廣泛使用的迭代梯度分類算法，可無監(jiān)督地將樣本集進(jìn)行聚類，使得各聚類內(nèi)部的數(shù)據(jù)相關(guān)性比聚類之間的相關(guān)性要高，實(shí)施過程簡(jiǎn)單快速，通過尋找各簇的聚類中心，實(shí)現(xiàn)數(shù)據(jù)分類，具體流程可參考文獻(xiàn)[14]。

2 線譜提取算法

2.1 算法原理分析

軸頻電場(chǎng)信號(hào)在頻域上呈現(xiàn)分立多線譜結(jié)構(gòu)，且線譜頻率間存在諧波關(guān)系，線譜所在子區(qū)間功率譜熵值較小，而理論分析和大量實(shí)測(cè)數(shù)據(jù)都表明，海洋環(huán)境電場(chǎng)噪聲則沒有這樣的特征結(jié)構(gòu)，與軸頻信號(hào)的功率譜熵差異明顯。考慮到功率譜的能量幅值具有實(shí)際的物理意義，為克服采用原始功率譜熵進(jìn)行篩選帶來的偏差，算法定義了功率譜窄帶子區(qū)間的峰值熵比（EPER）特征：

其中Eimax為對(duì)頻率進(jìn)行等長(zhǎng)區(qū)間劃分后，子區(qū)間i的峰值能量，Hi為由（1）式計(jì)算出的相應(yīng)譜熵，由式（2）可知，EPER是線譜特性的有效度量。

采用EMD算法可依據(jù)信號(hào)的尺度特點(diǎn)自適應(yīng)地進(jìn)行尺度分解，由于海洋電場(chǎng)噪聲相對(duì)于軸頻信號(hào)來說通常頻率較低[1]，根據(jù)各IMF的頻率分布區(qū)間進(jìn)行有效IMF的選擇。

得到了各有效IMF分量的各子區(qū)間的EPER特征后，需要考慮如何選擇的問題，這也是算法的關(guān)鍵步驟，K-均值聚類算法可以根據(jù)樣本集自身的分布特點(diǎn)進(jìn)行聚類，避免了采用固定閾值的不自適應(yīng)性和可能導(dǎo)致的遺漏。

K-均值聚類初步篩選后，對(duì)于各候選子區(qū)間的進(jìn)一步處理，則需要根據(jù)EMD的分解特點(diǎn)和軸頻、噪聲間的物理特征差異，下面結(jié)合算例進(jìn)行介紹。

2.2 算法流程描述

1）采用EMD算法分解序列X，根據(jù)采樣率fs和軸頻特征頻段（1～7 Hz），選擇出L個(gè)有效IMF（文中取為前3個(gè)）；

2）對(duì)各有效IMF求取功率譜，進(jìn)行子區(qū)間等長(zhǎng)劃分，得到L×K個(gè)子區(qū)間；

3）按式（2）計(jì)算全部L×K個(gè)子區(qū)間的EPER值，采用K-均值聚類算法進(jìn)行聚類，選擇具有較大特征值的聚類內(nèi)的子區(qū)間進(jìn)行下一步處理；

4）對(duì)于聚類篩選后的候選子區(qū)間，舍棄出現(xiàn)于同一IMF中連續(xù)子區(qū)間，將剩下的分立子區(qū)間中EPER值最大的作為線譜頻率所在子區(qū)間，將該子區(qū)間內(nèi)的能量最大處頻率作為線譜頻率。

2.3 算例驗(yàn)證

利用自行研制的內(nèi)記式沉底電場(chǎng)測(cè)量系統(tǒng)，在東海某開放海域?qū)崪y(cè)大量船舶軸頻電場(chǎng)信號(hào)和海洋背景噪聲數(shù)據(jù)，為算法分析提供了有力的數(shù)據(jù)驗(yàn)證。下面采用不同情況下的實(shí)測(cè)數(shù)據(jù)進(jìn)行線譜提取，以驗(yàn)證算法的有效性。采樣頻率20 Hz。

2.3.1 “純”信號(hào)時(shí)

實(shí)測(cè)某船軸頻電場(chǎng)信號(hào)和背景電場(chǎng)噪聲的時(shí)域波形圖和歸一化功率譜分別如圖1所示。功率譜采用周期圖法估計(jì)。

由圖1可知，軸頻電場(chǎng)信號(hào)具有明顯的分立多倍頻線譜結(jié)構(gòu)，而背景電場(chǎng)噪聲則主要在某頻率附近呈連續(xù)性的能量遞減特點(diǎn)，峰值頻率和海水運(yùn)動(dòng)周期有關(guān)。利用所提算法對(duì)圖1（a）中“純”軸頻信號(hào)進(jìn)行處理，處理結(jié)果見表1。

可見，EMD將信號(hào)中倍頻分量分解到連續(xù)2個(gè)IMF中去，表明了EMD的二進(jìn)尺度分解能力，注意到子區(qū)間序號(hào)相隔較遠(yuǎn)，這是線譜分立性的體現(xiàn)。注意此時(shí)7.5 Hz和9.9 Hz處的線譜被K-均值聚類算法認(rèn)為是EPER值較小的一類，表1中未給出。

表 1 “純”信號(hào)時(shí)算法處理結(jié)果Tab. 1 Algorithm processing results with pure signal

2.3.2 “純”噪聲時(shí)

利用所提算法處理圖1（b）中環(huán)境噪聲數(shù)據(jù)，結(jié)果見表2。

表 2 “純”噪聲時(shí)算法處理結(jié)果Tab. 2 Algorithm processing results with pure noise

首先注意到此時(shí)這些子區(qū)間連續(xù)性地分布于同一個(gè)IMF中，這與圖1（b）中噪聲的功率譜分布特點(diǎn)也是相符的。實(shí)際中，環(huán)境電場(chǎng)主要來源于各種復(fù)雜形式的海水運(yùn)動(dòng)對(duì)地磁場(chǎng)的切割，這種運(yùn)動(dòng)會(huì)導(dǎo)致在相鄰的頻帶子區(qū)間上具有相近的能量和熵分布特性，即EPER值差別應(yīng)不大。反之，如果是軸頻信號(hào)，首先，相比于相鄰子區(qū)間，線譜所在子區(qū)間內(nèi)的能量峰值更大、熵值更小，導(dǎo)致EPER值很大，這種EPER值連續(xù)地取相近值很難出現(xiàn)；其次，考慮到EMD分解的二進(jìn)分解特性，多條特征線譜會(huì)被分解到不同的IMF中去，也不會(huì)出現(xiàn)子區(qū)間都位于一個(gè)IMF中的現(xiàn)象，將同一IMF中連續(xù)地取相近EPER值的子區(qū)間認(rèn)為是噪聲，是合理的。

2.3.3 混合信號(hào)時(shí)

由于測(cè)量期間未出現(xiàn)高海況，將圖1中信號(hào)和噪聲按指定信噪比進(jìn)行幅值線性疊加，以模擬生成相應(yīng)的低信噪比數(shù)據(jù)。設(shè)Es和En分別為信號(hào)和噪聲的能量，信噪比定義如下，

圖2分別為信噪比為–6.9 dB和–13.6 dB時(shí)的混合信號(hào)時(shí)域、頻域圖。從中可以看出，信號(hào)線譜都無法被分辨出。其中，–6.9 dB對(duì)應(yīng)信號(hào)的特征線譜強(qiáng)度開始小于噪聲強(qiáng)度的最大值，此時(shí)直接功率譜算法線譜提取開始失效。

對(duì)圖2中2種低信噪比下混合數(shù)據(jù)，采用文中算法進(jìn)行處理，結(jié)果分別如表3和表4所示。

表 3–6.9 dB時(shí)混合信號(hào)算法處理結(jié)果Tab. 3 Algorithm processing results with mixed signal under –6.9 dB SNR

表 4–13.6 dB時(shí)混合信號(hào)算法處理結(jié)果Tab. 4 Algorithm processing results with mixed signal under –13.6 dB SNR

2種信噪比下，EMD分解后的前4個(gè)IMF分量的功率譜圖分別如圖3所示，后面都是更低頻的分量，此處不予畫出。

表3和表4表明，文中算法都能夠正確地提取出特征線譜。由表4可知，即使噪聲子區(qū)間比信號(hào)子區(qū)間的EPER特征值要大，但在同一IMF中連續(xù)分布的特性并沒有改變，體現(xiàn)了EPER特征對(duì)噪聲強(qiáng)度變化的穩(wěn)健性，同時(shí)可注意到K-均值聚類方法成功選擇出了特征值較小的信號(hào)子區(qū)間，避免了采用閾值篩選可能的遺漏問題。

圖3（a）進(jìn)一步表明，雖然此時(shí)噪聲幅值較大，但EMD分離出了2.5Hz主導(dǎo)的信噪比較高的IMF2子帶，驗(yàn)證了EMD方法優(yōu)良的自適應(yīng)分解性能。

3 結(jié) 語

本文提出了一種結(jié)合EMD、功率譜子區(qū)間能量峰值熵比（EPER）特征、K-均值聚類方法的船舶軸頻電場(chǎng)線譜提取算法，通過對(duì)軸頻電場(chǎng)信號(hào)和環(huán)境電場(chǎng)噪聲物理特征差異的分析，確定了特征候選子區(qū)間的選擇方法。實(shí)測(cè)數(shù)據(jù)驗(yàn)證了算法在低信噪比下的良好性能。對(duì)算法應(yīng)用性能存在影響的參數(shù)分段數(shù)K、聚類數(shù)P，下一步將對(duì)其自動(dòng)確定方法進(jìn)行研究。

參考文獻(xiàn)：

[ 1 ]林春生, 龔沈光. 艦船物理場(chǎng)[M]. 北京: 兵器工業(yè)出版社,2007: 237–243.

[ 2 ]熊露, 姜潤(rùn)翔, 龔沈光. 淺海中船舶軸頻電場(chǎng)建模方法[J]. 國(guó)防科技大學(xué)學(xué)報(bào), 2014, 36(1): 98–103.XIONG Lu, JIANG Run-xiang, Gong Shen-guang. Ship modeling method of shaft-ELFE in shallow sea[J]. Journal of National University of Defense Technology, 2014, 36(1):98–103.

[ 3 ]DYMARKOWSKI K, UCZCIWEK J. The extremely low frequency electromagnetic signature of the electric field of the ship[C]//Proceedings of Undersea Defense Technology,Washington: Nexus Media Limited, 2001: 1–6.

[ 4 ]常明, 姜潤(rùn)翔, 張伽偉, 等. 基于主動(dòng)軸接地系統(tǒng)的船舶軸頻電場(chǎng)抵消方法[J]. 海軍工程大學(xué)學(xué)報(bào), 2015, 27(1): 64–67.CHANG Ming, JIANG Run-xiang, ZHANG Jia-wei, et al.Method of reducing ship's shaft-rate electric field based on ASG system[J]. Journal of Naval University of Engineering,2015, 27(1): 64–67.

[ 5 ]MILLER R C. Active shaft grounding and diagnostic system:U.S, 4, 873, 512[P]. 1989–10-10.

[ 6 ]李松. 船舶電場(chǎng)目標(biāo)特征提取與檢測(cè)方法研究[D]. 武漢: 海軍工程大學(xué), 2008.LI Song. Research on feature extraction and detection algorithms of ship electric field signal[D]. Wuhan: Naval University of Engineering, 2008.

[ 7 ]賈亦卓, 姜潤(rùn)翔, 龔沈光. 基于小波模極大值的船舶軸頻電場(chǎng)檢測(cè)算法研究[J]. 兵工學(xué)報(bào), 2013, 34(5): 579–584.JIA Yi-zhuo, JIANG Run-xiang, GONG Shen-guang. Research on wavelet modulus maximum-based detection algorithm of ship’s shaft-rate electric field[J]. Acta Ar-mamentarii, 2013,34(5): 579–584.

[ 8 ]胡鵬, 肖愛明, 龔沈光. 基于小波尺度相關(guān)的船舶軸頻電場(chǎng)信號(hào)檢測(cè)算法[J]. 探測(cè)與控制學(xué)報(bào), 2013, 35(2): 37–40.HU Peng, XIAO Ai-ming, GONG Shen-guang. Detection algorithm for ship shaft-rate electric field signal based on wavelet scale correlation [J]. Journal of Detection & Control,2013, 35(2): 37–40.

[ 9 ]胡鵬, 包中華, 龔沈光. 基于小波包熵的船舶軸頻電場(chǎng)信號(hào)消噪研究[J]. 應(yīng)用基礎(chǔ)與工程科學(xué)學(xué)報(bào), 2011, 19(4):672–678.HU Peng, BAO Zhong-hua, GONG Shen-guang. Denoising the signal of ship shaft-rate electric field based on wavelet packet entropy[J]. Journal of Basic Science and Engineering, 2011,19(4): 672–678.

[10]程銳, 姜潤(rùn)翔, 龔沈光. 基于EMD和4階累積量的船舶軸頻電場(chǎng)線譜提取[J]. 艦船科學(xué)技術(shù), 2016(1): 94–98.CHENG Rui, JIANG Run-xiang, GONG Shen-guang. Extraction of line spectrum of the ship shaft-rate electric field based on EMD and fourth-order cumulant[J]. Ship Science and Technology, 2016(1): 94–98.

[11]HUANG N E, SHEN Z, LONG S R, et al. The empirical mode decomposition and the Hilbert spectrum for nonlinear and nonstationary time series analysis[C]//Proceedings of the Royal Society of London A: Mathematical, Physical and Engineering Sciences, 1998, 454(1971): 903–995.

[12]FLANDRIN P, RILLING G, GONCALVES P. Empirical mode de-composition as a filter bank [J]. IEEE Signal Processing Letters, 2004, 11(2): 112–114.

[13]劉榮, 劉珩. 低信噪比下基于功率譜熵的語音端點(diǎn)檢測(cè)算法[J]. 計(jì)算機(jī)工程與應(yīng)用, 2009, 45(33): 122–124.LIU Rong, LIU Kan. Power spectrum entropy based voice activity detection algorithm in low signal-to-noise ratio conditions[J]. Computer Engineering and Applications, 2009,45(33): 122–124.

[14]HARTIGAN J A, WONG M A. Algorithm AS 136: A k-means clustering algorithm[J]. Journal of the Royal Statistical Society.Series C (Applied Statistics), 1979, 28(1): 100–108.

[15]張自力, 魏文博, 劉保華, 等. 海浪感應(yīng)電磁場(chǎng)的理論計(jì)算[J].海洋學(xué)報(bào), 2008, 30(1): 42–46.ZHANG Zi-li, WEI Wen-bo, LIU Bao-hua, et al. Theoretical calculation of electromagnetic field generated by ocean waves[J]. Acta Oceanologica Sinica, 2008, 30(1): 42–46.

Line spectrum extraction of ship shaft-rate electric field combining EMD and power spectra entropy

CHENG Rui1, CHEN Cong2, JIANG Run-xiang1

(1. Department of Weaponry Engineering, Naval University of Engineering, Wuhan 430033, China;2. College of Sciences, Naval University of Engineering, Wuhan 430033, China)

In order to implement detection of weak ship shaft-rate electric field signal buried in strong marine background noise, a new line spectrum extraction algorithm combining Empirical Mode Decomposition (EMD) and narrowband sub-interval power spectra entropy is proposed. Firstly, a set of effective Intrinsic Mode Functions (IMFs) were separated from noise-polluted signal by means of EMD method, of which the power spectra was then divided into subinterval sections.Furthermore, a new feature named Energy Peak Entropy Ratio (EPER) of each section was defined and computed. Finally,taking advantage of analyzing differences of physical properties between shaft-rate signal and ambient noise, together withK-means clustering method, line spectrum was extracted. Processing results of sea measured data indicate that, comparing to direct power spectra analysis, the line spectrum detection low bound is reduced by 6.7 dB.

shaft-rate electric field；empirical mode decomposition (EMD)；K-means clustering；energy entropy；line spectrum extraction

TP274.5

A

1672 – 7649(2017)09 – 0159 – 05

10.3404/j.issn.1672 – 7649.2017.09.032

2016 – 07 – 07；

2016 – 08 – 22

國(guó)家自然科學(xué)基金資助項(xiàng)目(51109215，51509252)

程銳(1988 – )，男，博士研究生，研究方向?yàn)槟繕?biāo)特性及其信息處理技術(shù)。