小波包節(jié)點分段閾值降噪在水聲監(jiān)聽中的應(yīng)用?

趙 杰 楊 英123 惠 力123 王 志123 初士博123 劉茂科123

(1齊魯工業(yè)大學(xué)(山東省科學(xué)院) 濟南 250353)

(2山東省科學(xué)院海洋儀器儀表研究所 山東省海洋環(huán)境監(jiān)測技術(shù)重點實驗室 青島 266061)

(3國家海洋監(jiān)測設(shè)備工程技術(shù)研究中心 青島 266061)

0 引言

目前，水聲監(jiān)聽在海洋生物調(diào)查研究、艦艇監(jiān)聽監(jiān)測、深海勘探等方面需求量具大且應(yīng)用前景十分廣闊，其監(jiān)聽信號主要為寬帶中低頻信號和窄帶高頻信號，如魚類的摩擦聲頻率一般在100 Hz～8 kHz，魚鰾發(fā)聲多數(shù)在75～100 Hz；艦船噪聲譜主要集中在1 kHz以下；水聲場測量、海洋地震測量、淺地層剖面測量等主要采用寬帶中低頻進行測量，頻率范圍在幾百Hz到十幾kHz；側(cè)掃聲吶及地貌掃描聲吶等采用窄帶高頻信號，頻率范圍在十幾kHz到幾百kHz。在實際的水聲探測、接收和識別過程中，海洋環(huán)境噪聲信息會弱化相關(guān)聲目標信號監(jiān)聽的準確度[1-2]，嚴重影響目標的后續(xù)識別與處理，特別是寬帶中低頻目標信號更容易摻雜大量復(fù)雜的噪聲信息。因此，如何將目標信號與噪聲有效分離與提取成為水聲監(jiān)聽研究的重中之重。

傳統(tǒng)的數(shù)字信號處理及去噪方法多數(shù)是以傅里葉變換為基礎(chǔ)的，如司新新等[3]利用短時傅里葉變換進行時變和非平穩(wěn)信號的處理和分析，姚家駿等[4]在地震波特征量分析中采用短時傅里葉變換，曲麗榮[5]選用合適的短時傅里葉窗口函數(shù)對信號頻率特性分析。以上研究均以傅里葉變換為基礎(chǔ)的傳統(tǒng)數(shù)字信號去噪法屬于全局變換，無法同時表征時頻雙域特性，只能反映信號整體特性，時域無法局部化，難以檢測局部突變信號，去噪局限性較大。小波變換去噪方法具有良好的時頻局部性、多分辨率、去相關(guān)性、快速算法和選基靈活等特點，與傅里葉變換去噪方法相比具有明顯的優(yōu)越性，成為國內(nèi)外許多學(xué)者專家應(yīng)用研究的去噪方法。Kalpana等[6]對水聲通信中提高信噪比的去噪技術(shù)進行研究，利用Gabor小波變換來實現(xiàn)去噪提高信噪比，并與香農(nóng)小波進行了對比研究。Hill等[7]使用雙數(shù)復(fù)小波系數(shù)法實現(xiàn)圖像噪聲處理研究。普利達大學(xué)研究人員Ghribi等[8]提出用小波變換分解的雙線性結(jié)構(gòu)，通過雙通道對稱方法實現(xiàn)小波自適應(yīng)去噪聲。周建等[9]利用一種新的小波閾值去噪方法對振動信號進行去噪分析。張振鳳等[10]利用改進閾值隨實際分解層數(shù)而變化的方法進行去噪，通過信噪比和均方差的比較，該方法優(yōu)于傳統(tǒng)固定閾值、stein無偏似然估計閾值、極大極小值的小波去噪方法。曹建華等[11]利用小波分析法對振動信號進行去噪處理，對去噪后信號進行功率譜密度分析。余本富等[12]基于自適應(yīng)分解層數(shù)和閾值的小波去噪算法降低電能質(zhì)量擾動信號中的噪聲。李戰(zhàn)明等[13]對小波分析中4種去噪方法分析比較，通過對比分析采用改進閾值函數(shù)提升小波變化閾值法用于心電信號去噪分析。以上小波變換的研究方法在閾值和閾值函數(shù)選取等方面具有較高的參考價值，但其主要針對頻率變化緩慢和頻帶較窄的非平穩(wěn)振動、心電等生物醫(yī)學(xué)、圖像處理等方面信號，不能同時兼顧低頻和高頻目標，無法實現(xiàn)寬頻帶信號分頻段聚焦處理，因此，對頻段變換范圍較大的水聲監(jiān)聽信號消噪及提取處理適應(yīng)性較差。

水聲信號屬于非平穩(wěn)、非線性時變信號，受噪聲及混響等干擾的影響，探測獲得的水聲信號往往會存在很多突變的奇異點。小波變換具有多分辨和時頻局部分析的能力，與短時傅里葉去噪法相比，能更好地去噪，小波和小波包降噪處理的核心是閾值和閾值函數(shù)的選取，常用的小波硬閾值和軟閾值函數(shù)達到了一定的去噪效果，但硬閾值處理后的小波系數(shù)在正負閾值處不連續(xù)，軟閾值可有效避免間斷連續(xù)性較好，但小波系數(shù)較大時，系數(shù)絕對值與系數(shù)之間存在偏差，易造成高頻信息損失，重構(gòu)信號易產(chǎn)生誤差。小波包分解系數(shù)按照不同的閾值規(guī)則能更有效去除噪聲，然而不同頻段的小波包分解系數(shù)對目標信號和噪聲信號反映不同，用同一閾值和閾值規(guī)則對整個頻段進行處理無法達到較為理想的效果。因此，本文采用小波包相對能量法進行最優(yōu)分解層判別，考慮各個節(jié)點頻段有用信息的比例，按照比例大小確定閾值準則及函數(shù)，進一步提高分辨率，優(yōu)化高頻信號分析，去噪提取結(jié)果更加精細優(yōu)化。利用小波包節(jié)點系數(shù)分段閾值處理的方法對監(jiān)聽獲得的100 Hz～50 kHz范圍內(nèi)水下聲音數(shù)據(jù)降噪處理分析，對不同頻段范圍內(nèi)的水聲信號進行提取，噪聲信號分離，結(jié)果表明：該方法既可以對1 kHz上下范圍內(nèi)的較強信號進行降噪提取，也可對較高頻率相對較弱的信號實現(xiàn)降噪分離。

1 小波包算法基本原理

小波包算法為最佳子帶樹結(jié)構(gòu)(Optimal subband tree structuring)，包括分解算法和重構(gòu)算法。用子帶樹來表示小波包分解，利用多次迭代的小波轉(zhuǎn)換對信號的低頻和高頻部分進行多層次劃分，根據(jù)信號特征，自適應(yīng)選擇頻帶范圍，分析信號細節(jié)部分。小波包克服了小波分析的高頻分辨率低的缺點，對信號在全頻帶范圍內(nèi)進行正交分解，提高時-頻分辨率[14]。子帶樹結(jié)構(gòu)圖如圖1所示。

圖1 子帶樹結(jié)構(gòu)Fig.1 Subband tree structure

將尺度空間Vj和小波子空間Wj統(tǒng)一用U(n,j)新的子空間來表示：

定義子空間U(n,j)為函數(shù)Un(t)的閉包空間，Un(t)是函數(shù)U2n(t)的閉包空間，令Un(t)滿足雙尺度方程：

式(2)系數(shù)為多分辨率分析中的濾波器系數(shù)，滿足正交關(guān)系，g(k)=(-1)kh(1-k)。當(dāng)n=0，

在多分辨率分析中，U0(t)和U1(t)分別退化為尺度函數(shù)φ(t)和小波基函數(shù)ψ(t)：

小波包變換的分解與合成中，存在由多種小波基組成的基底庫，因此需要根據(jù)不同場合、不同信號并且在綜合考慮小波基的正交性、緊支性、對稱性等特性的基礎(chǔ)上選擇合適的小波基。從小波庫中尋找使代價函數(shù)最小的小波包基，選擇關(guān)鍵是代價函數(shù)[15]。

式(5)中，s代表信號，si表示信號在小波包基上的分解系數(shù)。求出使代價函數(shù)最小的小波包序列，從而求出最優(yōu)基。代價函數(shù)稱為熵或平均信息量，代價函數(shù)自變量是信號在某個基向量組上分解系數(shù)，函數(shù)實部因變量即為代價。采用二叉樹自下向上搜索最優(yōu)小波包基，小波包基的選取是后續(xù)數(shù)據(jù)處理的關(guān)鍵。

1.1 小波包相對能量確定最優(yōu)分解層數(shù)

小波包分解需要獲得一個最優(yōu)的分解層數(shù)。分解層數(shù)直接影響信號的去噪效果，分解層數(shù)過多對各層小波系數(shù)閾值處理過程中造成信息量丟失嚴重、信噪比下降、運算量增加。分解層數(shù)過少，信噪比提高不多，去噪效果不明顯[16]。原始信號和噪聲信號的小波系數(shù)在不同層次有不同特性，原始信號的小波系數(shù)會隨著分解層次增加而增大，噪聲信號小波系數(shù)減小。

通過最優(yōu)小波包中節(jié)點能量判斷的方法確定小波包分解層數(shù)。根據(jù)最優(yōu)小波基節(jié)點能量的大小，判斷該節(jié)點有用信號能量大小。能量越高，表示該頻段的有效成分越多，反之，噪聲成分越多。據(jù)此，根據(jù)能量值判斷有效成分，確定最優(yōu)分解層數(shù)。

如圖1所示，U(i,j)表示第j層的第i個節(jié)點，每個節(jié)點系數(shù)表示某一頻段內(nèi)的信號特征，其中j=0,1,···,n，i=0,1,···,2j-1。能量計算公式為[17]

式(6)中，s(i,j)為U(i,j)節(jié)點的小波包系數(shù)。由此可得出每層分解后子帶的相對能量公式為

其中，k(i,j)為節(jié)點i在第j層的相對能量，代表有效成分的含量，對每層的節(jié)點相對能量k(i,j)由大到小進行排序，k(i,j)的最大值記為C(1,j)，最小值記為C(n,j)，即C(1,j)≥C(2,j)≥···≥C(n,j)，其中j代表層數(shù)。首先對信號第一層分解相對能量進行初判：C(1,1)≥0.9、C(2,1)＜0.1信號能量頻段相對集中的情況，這種情況下對每層的k(2,j)進行對比，當(dāng)?shù)趈層的C(2,j)大于第j-1層C(2,j-1)時，將j-1作為最優(yōu)分解層數(shù)。

其他情況下，對信號逐層分解，將C(1,j)首次小于0.2的層數(shù)j進行記錄，如果j＞1，最優(yōu)分層數(shù)為m=(j-1)；如果j=1，最優(yōu)分解層數(shù)為1。

1.2 小波包分段閾值處理法

水聲信號的頻率范圍比較廣，不能簡單地將低頻段作為有用信息、高頻段作為噪聲信號處理，按小波相對能量法實現(xiàn)分段閾值降噪處理，該方法與單一閾值處理相比充分考慮信號頻率段與噪聲分布情況。針對不同的分解系數(shù)選取最合適的閾值函數(shù)和閾值th，最大程度去除噪聲，保留有效成分[18]。目前，常用的閾值估計方法主要有以下幾種：

(1)固定閾值估計法(sqtwolog)，在正態(tài)高斯噪聲模型下，針對多維獨立正態(tài)變量分布得出的結(jié)論。

(2)無偏似然閾值估計法(rigrsure)，小波分解系數(shù)的平方記為λ，從小到大的順序排列，即λ1≤λ2≤···≤λn。設(shè)風(fēng)險向量R，其元素為ri，

風(fēng)險向量R的最小值rmin作為風(fēng)險值，對應(yīng)求出對應(yīng)的λmin，則閾值為

(3)啟發(fā)式閾值法(heursure)，是前兩種閾值的綜合，設(shè)S為N個小波系數(shù)的平方和。令

(4)極大極小閾值法(minimaxi)，也是固定閾值的一種，產(chǎn)生最小均方誤差的極值。

以上各式中σ=MAD/0.6745為噪聲標準方差，MAD為各高頻子帶系數(shù)的中值；N為信號長度。

如1.1節(jié)所述，若C(1,1)≥ 0.9、C(2,1)＜0.1信號能量頻段相對集中，只對C(1,1)≥0.9進行小波去噪。其他情況下，信號最優(yōu)分解層m上，節(jié)點相對能量k(i,j)≥0.2的節(jié)點系數(shù)采用無偏似然估計閾值估計方法；k(i,j)≥0.1、k(i,j)＜0.2的節(jié)點系數(shù)采用極大極小閾值法；k(i,j)＜0.1的節(jié)點系數(shù)采用統(tǒng)一閾值法選擇。多閾值準則選取表如1所示。

表1 多閾值準則選取表Table1 Multiple threshold criteria selection table

目前常用的閾值函數(shù)為硬閾值函數(shù)和軟閾值函數(shù)。硬閾值將絕對值不大于閾值的元素置零，存在不連續(xù)現(xiàn)象，數(shù)學(xué)表達式為

軟閾值函數(shù)在式(13)的基礎(chǔ)上實現(xiàn)連續(xù)點收縮為零，有效避免間斷，連續(xù)性較好，但小波系數(shù)較大時，系數(shù)絕對值與系數(shù)之間存在偏差，易造成高頻信息損失，重構(gòu)信號易產(chǎn)生誤差，其表達式為

考慮到硬閾值和軟閾值函數(shù)處理信號細節(jié)成分的局限性，利用文獻[19]中的新型閾值函數(shù)，公式如下：

式(15)中，x表示小波包分解后的節(jié)點系數(shù)，y為閾值函數(shù)處理后的小波節(jié)點系數(shù)，th為選定的閾值。式(14)和式(15)中sign(x)為符號函數(shù)，當(dāng)x＞0，sign(x)=1；當(dāng)x=0，sign(x)=0；當(dāng)x＜0，sign(x)=-1。閾值函數(shù)處理結(jié)果對比如圖2、圖3所示。

圖2 軟、硬閾值函數(shù)結(jié)果對比Fig.2 Comparison of soft and hard threshold function results

圖3 新、軟閾值函數(shù)結(jié)果對比Fig.3 Comparison of new and soft threshold function results

2 監(jiān)聽水聲信號去噪分析

監(jiān)聽水聲信號來源于阿拉斯加費班克大學(xué)Erin Pettit和華盛頓大學(xué)JeffNystuen等在阿拉斯加州冰灣利用自主水聽器(被動水聽器(PAL))獲得的被動水聲記錄[20]，該海域的水深為160 m，水聽器的布放水深為90 m。所有數(shù)據(jù)采樣頻率為100 kHz，從100 Hz～50 kHz共分為64個頻段，其中100 Hz～3 kHz的范圍內(nèi)的頻率分辨率為200 Hz，3 kHz～50 kHz范圍的為分辨率為1 kHz。讀取2009年某一時間段內(nèi)的原始音頻文件，其時頻特性標度圖如圖4所示。

由圖4的時-頻域分布標度圖可以看出原始水聲監(jiān)聽信號頻率隨時間變化的強弱分布，本文的主要工作是對不同頻段范圍、不同強弱的信號進行消噪并且提取。將最原始水聲監(jiān)聽信號按db6進行小波包分解計算，計算每層節(jié)點相對能量并由大到小排序(每層節(jié)點僅排前四位)，如表2所示。根據(jù)每層節(jié)點相對能量的大小按照1.1節(jié)判別方法獲得最優(yōu)分解層數(shù)為5層。按小波包5層將水聲信號分解，進行最優(yōu)樹結(jié)構(gòu)和節(jié)點系數(shù)頻率大小排序，分別如圖5和圖6所示。的方式實現(xiàn)消噪分離，處理后時-頻域標度分布如圖7所示。依此頻段來進行本文方法與全局單一閾值對比，將小波包硬閾值、小波包軟閾值在5層分解的基礎(chǔ)上對信號進行消噪處理，獲得的時-頻域信號分布圖，分別如圖8和圖9所示。由圖7、圖8和圖9三幅圖片對比可以看出，本文方法具有良好的優(yōu)越性。

表2 節(jié)點相對能量排序表Table2 Nodes relative energy ranking table

圖4 監(jiān)聽水聲信號時-頻分布標度圖Fig.4 Time-frequency distribution scale diagram for monitoring underwater acoustic signals

圖5 最優(yōu)樹結(jié)構(gòu)Fig.5 Optimal tree structure

圖6 節(jié)點系數(shù)頻率排序圖Fig.6 Node coefficient frequency sorting diagram

圖7 分段閾值降噪信號時-頻分布標度圖Fig.7 Time-frequency distribution scale diagram of segmental threshold denoising signal

圖6中橫軸為時間，左邊縱軸為頻率大小順序，右邊縱軸為小波包分解第五層節(jié)點系數(shù)31-62按照頻率順序由下往上進行排序。由圖6可看出監(jiān)聽水聲信號頻段主要分布在31、32、34三個節(jié)點上，頻率范圍在0～5 kHz。通過本文方法，首先對頻段在節(jié)點31處監(jiān)聽水聲信號，頻率范圍大約在0.1 kHz～1.6 kHz，按照節(jié)點系數(shù)分頻段閾值處理

圖8 硬閾值降噪信號時-頻分布標度圖Fig.8 Time-frequency distribution scale diagram of hard threshold denoising signal

圖9 軟閾值降噪信號時-頻分布標度圖Fig.9 Time-frequency distribution scale diagram of soft threshold denoising signal

圖7為利用本文方法的消噪提取結(jié)果，信號提取結(jié)果較為理想，噪聲基本消除。圖8中為小波包硬閾值5層分解處理后的監(jiān)聽水聲信號時-頻域標度分布圖，與圖4原始信號相比，信號降噪較為明顯，但依然存在噪聲。圖9為小波包軟閾值降噪處理，與圖8相比，噪聲頻率信號進一步減少，但降噪不夠徹底，仍存在小部分噪聲。小波包分段閾值降噪分離的關(guān)鍵在于在未知真實有效的水聲信號前提下，節(jié)點系數(shù)按頻率排序，節(jié)點系數(shù)頻率由低到高、信號強度由強到弱的排序為31、32、34、33、37，其頻率范圍大約在0.1 kHz～1.6 kHz、1.7 kHz～3.3 kHz、 3.4 kHz～5 kHz、 5.1 kHz～6.7 kHz、6.8 kHz～8.4 kHz，按照本文方法對剩余4個節(jié)點頻率信號進行提取，其消噪提取結(jié)果較為理想，結(jié)果如圖10～13所示。由圖12、圖13可以看出頻率范圍相對較高、信號強度較弱的節(jié)點33、節(jié)點37的信號消噪提取亦可通過本文方法實現(xiàn)。

由于監(jiān)聽的真實水聲信號未知，參照文獻[21]平滑度指標、信噪比增益和峰值信噪比三種方法對監(jiān)聽水聲信號進行去噪評價。三種方法的定義及公式表示如下：

(1)平滑度指標是指去噪后的差分數(shù)的方差根和原始信號的差分數(shù)的方差根之比，記為r，

其中，f′(n)為去噪后函數(shù)，f(n)為原始信號。

(2)信噪比增益是指去噪后的信噪比和去噪前的信噪比的比值，記為GSNR，

其中，SNRdn為去噪后信噪比，SNRn為原始信號信噪比。

(3)峰值信噪比信號最大可能功率和影響它的表示精度的破壞性噪聲功率的比值，記為PSNR，

其中，f′(n)為去噪后函數(shù)，f(n)為原始信號。

圖10 節(jié)點32頻段降噪信號時-頻分布標度圖Fig.10 Node 32 time-frequency distribution scale diagram of segmental threshold denoising signal

圖11 節(jié)點34頻段降噪信號時-頻分布標度圖Fig.11 Node 34 time-frequency distribution scale diagram of segmental threshold denoising signal

圖12 節(jié)點33頻段降噪信號時-頻分布標度圖Fig.12 Node 33 time-frequency distribution scale diagram of segmental threshold denoising signal

圖13 節(jié)點37降噪信號時-頻分布標度圖Fig.13 Node 37 time-frequency distribution scale diagram of segmental threshold denoising signal

平滑度指標數(shù)值越小，信噪比增益和峰值信噪比越大，去噪效果越好。表3是采用上述三種方法分別對小波包硬閾值，小波包軟閾值和小波包分段閾值的降噪處理結(jié)果進行去噪評價，由表3可以看出本文方法是可靠的。

表3 評價結(jié)果Table3 The evaluation results

3 結(jié)論

通過對水下監(jiān)聽信號的去噪分析可以得出本文方法的關(guān)鍵在于(1)按照每層小波包節(jié)點系數(shù)相對能量的變化確定小波包最優(yōu)分解層數(shù)；(2)按相對能量法實現(xiàn)分頻段閾值降噪處理。根據(jù)節(jié)點相對能量大小，判斷該節(jié)點有用信號能量大小。能量越高，表示該頻段的有效成分越多，反之，噪聲成分越多。針對節(jié)點系數(shù)相對能量大小選取最合適的閾值函數(shù)和閾值th，最大程度去除噪聲，保留有效成分。本文方法去噪效果優(yōu)于傳統(tǒng)小波、小波包全局單一閾值處理等方法，可在信號噪聲的處理中進行廣泛的推廣和應(yīng)用，特別在100 Hz～50 kHz頻率范圍內(nèi)的水聲監(jiān)聽信號降噪處理中具有重要參考價值。