基于小波域模糊閾值算法的脈沖星微弱信號(hào)消噪

甘偉，袁冰妍

（中國電子科技集團(tuán)公司第二十研究所，西安 710068）

0 引言

脈沖星是一種具有超高壓、超強(qiáng)磁場(chǎng)、超高溫、高速自轉(zhuǎn)及超高穩(wěn)定周期的中子星被譽(yù)為自然界最穩(wěn)定的時(shí)鐘。在脈沖星發(fā)現(xiàn)后的半個(gè)世紀(jì)中，其研究主要集中在天文和物理等科研領(lǐng)域，并取得了極大成功。近三年來，隨著航天技術(shù)的發(fā)展，脈沖星在精確定時(shí)和深空導(dǎo)航中的應(yīng)用研究逐漸受到關(guān)注[1]。但是，脈沖星距離地球非常遙遠(yuǎn)，在傳播過程中受到色散延遲、散射、相對(duì)論效應(yīng)、地球運(yùn)動(dòng)及其他未知因素的影響[2]，在太陽系中接收到的信號(hào)非常微弱且湮沒在噪聲中。因此，脈沖星信號(hào)的降噪、檢測(cè)方法研究對(duì)于脈沖星后續(xù)的科研、應(yīng)用研究都具有極其重要的意義。

脈沖星信號(hào)降噪的目的有兩點(diǎn)：第一，提高脈沖輪廓的精度，保留微脈沖等細(xì)節(jié)信息。累積脈沖輪廓是脈沖星極冠區(qū)的輻射窗口[3]，在多個(gè)頻段獲取精確的脈沖輪廓對(duì)研究脈沖星特定高度上的物理輻射機(jī)制具有重要意義。第二，提高信噪比。脈沖到達(dá)時(shí)間（TOA）的測(cè)量精度正比于信噪比[4]。因此，提高信噪比對(duì)于脈沖星自身物理特征研究、引力波輻射的檢測(cè)、脈沖星定時(shí)、導(dǎo)航等依賴于TOA精確測(cè)量的科研和應(yīng)用領(lǐng)域具有重要意義。

傳統(tǒng)信號(hào)去噪方法主要包括線性濾波和非線性濾波方法，如中值濾波和Wiener濾波等。傳統(tǒng)去噪方法的不足在于信號(hào)變換后的熵增高，無法刻畫信號(hào)的非平穩(wěn)特性并且無法得到信號(hào)的相關(guān)性。為了克服上述缺點(diǎn)，人們開始使用小波變換來解決信號(hào)去噪問題。小波變換具有低熵性、多分辨率、去相關(guān)性，選基靈活性等特點(diǎn)。而通常的小波軟閾值和硬閾值去噪算法不能解決抑制噪聲和保留脈沖細(xì)節(jié)之間的矛盾。如通用閾值可以大幅度提高信噪比，但其重構(gòu)波形過于平滑，微脈沖等細(xì)節(jié)有較大損失；Stein無偏風(fēng)險(xiǎn)閾值可以保留部分細(xì)節(jié)，但卻削弱了對(duì)噪聲的抑制作用[5]。為了盡可能地在抑制噪聲的條件下保留信號(hào)有用細(xì)節(jié)信息，文中提出了模糊閾值小波降噪算法，實(shí)驗(yàn)仿真驗(yàn)證了該算法的有效性。

1 小波去噪及典型閾值算法

1.1 小波去噪原理

含噪一維信號(hào)可以表示成：

且其二進(jìn)制小波變換如下：

其中，φ為小波基。為真實(shí)信號(hào)，為噪聲。實(shí)際中，有用信號(hào)通常表現(xiàn)為低頻信號(hào)或是一些比較平穩(wěn)的信號(hào)，而噪聲則通常表現(xiàn)為高頻信號(hào)。

由于有用信號(hào)與噪聲有不同的Lipschitz指數(shù)，反映到小波域，表現(xiàn)為信號(hào)的小波系數(shù)隨分解層數(shù)的增大而變化緩慢，噪聲的小波系數(shù)隨分解層數(shù)的增大而迅速衰減。在第3尺度上真實(shí)脈沖星信號(hào)的小波系數(shù)占主要部分，而噪聲基本上為零。此時(shí)，只需要設(shè)定適當(dāng)?shù)拈撝担纯蓪⒃肼暈V除。所以去噪過程可以按如下方法處理：

首先對(duì)信號(hào)f進(jìn)行小波分解（如進(jìn)行 3層分解），則噪聲部分通常包含在高頻系數(shù)cD1，cD2，cD3中，因而可以采用門限閾值等方法對(duì)小波系數(shù)進(jìn)行處理，然后對(duì)信號(hào)進(jìn)行重構(gòu)即可以達(dá)到去噪的目的。對(duì)信號(hào)f進(jìn)行去噪的目的就是要抑制信號(hào)中的噪聲成份，恢復(fù)出有用信號(hào)。

圖1 信號(hào)的3層小波分解圖

圖1 中，cA1，cA2，cA3分別為分解后第1～3層的低頻系數(shù)，cD1，cD2，cD3分別為分解后第1～3層的高頻系數(shù)。

1.2 小波去噪步驟

一般來說，信號(hào)的去噪過程可分為三步進(jìn)行：

（1）信號(hào)的小波分解。選擇小波基并確定小波分解的層數(shù)N，然后對(duì)信號(hào)f進(jìn)行N層小波分解。

（2）小波分解后高頻系數(shù)的閾值處理。對(duì)第1到第N層的每一層高頻系數(shù)，選擇一種閾值法進(jìn)行閾值量化處理。

（3）小波重構(gòu)。根據(jù)小波分解的第N層的低頻系數(shù)和經(jīng)過閾值量化處理后的第1層到第N層的高頻系數(shù)，進(jìn)行信號(hào)的小波重構(gòu)。

在這三個(gè)步驟之中，最關(guān)鍵的就是如何選取閾值和如何進(jìn)行閾值的量化，它直接關(guān)系到信號(hào)去噪的質(zhì)量。

1.3 小波去噪典型閾值算法

（1）硬閾值函數(shù)法其閾值函數(shù)為：

式中，λ為預(yù)置閾值或門限值，通常?。ㄍㄓ瞄撝捣ǎ?，σn是噪聲的方差，N是信號(hào)的長度。硬閾值方法可以很好地保留信號(hào)邊緣等局部特征，軟閾值處理相對(duì)要平滑，但會(huì)造成邊緣模糊等失真現(xiàn)象。

（2）軟閾值函數(shù)法其閾值函數(shù)為：

軟閾值法由它估計(jì)出來的小波系數(shù)與被處理信號(hào)的小波系數(shù)存在恒定的偏差。

（3）Stein風(fēng)險(xiǎn)閾值法：是采用Stein無偏似然原理（SURE）來確定的閾值。其具體算法如下[6]：

（a）將每一層小波變換的分解系數(shù)wj,k先平方，然后按由小到大的順序進(jìn)行排列，得到一個(gè)排序 后 的 向 量W= [W1,W2,…,Wn]且W1≤W2≤,…,≤Wn，其中為小波系數(shù)的個(gè)數(shù)。

（b）計(jì)算風(fēng)險(xiǎn)向量R=[r1,r2,…,rn]，其相應(yīng)元素為：

通過此式的多次迭代比較，可以得到一個(gè)最小的ri。再以求出的該最小值ri作為風(fēng)險(xiǎn)值，最后由該風(fēng)險(xiǎn)值ri的下標(biāo)i，求出對(duì)應(yīng)的Wi。

2 模糊閾值小波降噪算法

可以通過選取不同的f(wj,k)，得到不同的偏大型模糊分布。如采用升嶺形隸屬度函數(shù)時(shí)，相應(yīng)的f(wj,k)為：

為簡(jiǎn)化算法，本文采用升半梯形函數(shù)即：

因此所形成的隸屬度函數(shù)為：

式中，k為隸屬度函數(shù)曲線調(diào)節(jié)因子。式中a和b可以在下限和上限之間通過為隸屬度函數(shù)賦合適的值來決定。 由于第j尺度中幅度小于等于的細(xì)節(jié)小波系數(shù)最有可能由噪聲產(chǎn)生，所以當(dāng)細(xì)節(jié)小波系數(shù)的幅度等于時(shí)，可將隸屬度函數(shù)值設(shè)為接近于0的值，文中用ε來表示該值。同理，由于第j尺度中幅度大于等于的小波系數(shù)最有可能是原始信號(hào)造成的，因此當(dāng)小波系數(shù)的幅度等于時(shí)，可將隸屬度函數(shù)值設(shè)為接近于1的值，文中用1-ε來表示該值。將這些參數(shù)帶入隸屬度函數(shù)中可得：

因此可以推出參數(shù)a和b的表達(dá)式為：

3 實(shí)驗(yàn)結(jié)果及分析

3.1 脈沖星標(biāo)準(zhǔn)累積輪廓信號(hào)的降噪

本實(shí)驗(yàn)采用脈沖星 B0540-69的標(biāo)準(zhǔn)累積輪廓作為原始信號(hào)，如圖2（a）所示。為驗(yàn)證本文提出的模糊閾值小波降噪算法對(duì)脈沖星信號(hào)的消噪效果，首先對(duì)原始標(biāo)準(zhǔn)輪廓信號(hào)進(jìn)行加噪處理，如圖2（b）所示，加噪后的信號(hào)其信噪比為6.0439dB。實(shí)驗(yàn)中信噪比作如下定義參考為：，f為原始信號(hào)，為含噪信號(hào)。

首先在小波域采用小波函數(shù)wavedec對(duì)含噪脈沖星信號(hào)進(jìn)行小波分解，然后采用本文提出的方法計(jì)算每個(gè)樣本點(diǎn)幅值的隸屬度并確定該樣本點(diǎn)否是信號(hào)或噪聲的小波分解系數(shù)，最后采用小波函數(shù)waverec對(duì)模糊閾值小波降噪算法處理后的小波系數(shù)進(jìn)行重構(gòu)。實(shí)驗(yàn)中只對(duì)含噪脈沖星信號(hào)進(jìn)行3個(gè)尺度的小波分解，其中隸屬度函數(shù)曲線調(diào)節(jié)因子k設(shè)為0.7。本實(shí)驗(yàn)在Matalb 7.6.0環(huán)境下進(jìn)行仿真，仿真結(jié)果如圖2所示。

圖2 （a） 脈沖星B0540-69的標(biāo)準(zhǔn)累積輪廓信號(hào)

圖2 （b） 脈沖星B0540-69的含噪累積輪廓信號(hào)

圖2 （c） 小波“db6”軟閾值法消噪后的信號(hào)

圖2 （d） 小波“sym4”軟閾值法消噪后的信號(hào)

圖2 （e） 模糊閾值小波降噪算法消噪后的信號(hào)

為了與模糊閾值小波降噪算法進(jìn)行對(duì)比，本文采用了目前最為經(jīng)典的小波降噪算法，將原始含噪信號(hào)在小波基“db6”和小波基“sym4”下展開并進(jìn)行軟閾值消噪處理，實(shí)驗(yàn)結(jié)果分別如圖2（c）和圖2（d）所示，模糊閾值小波降噪算法的實(shí)驗(yàn)結(jié)果如圖2（e）所示。

從圖2（a）可以看出：脈沖星B0540-69的標(biāo)準(zhǔn)累積輪廓有兩個(gè)脈沖尖峰，一個(gè)“大尖峰”，一個(gè)“小尖峰”。從圖2（c）和圖2（d）可以看出：小波“db6”軟閾值法和小波“sym4”軟閾值法消噪后的信號(hào)過于平滑，把重要信息“小尖峰”已經(jīng)平滑掉，損失了過多的細(xì)節(jié)信息。而如圖2（e）所示的模糊閾值小波降噪算法能夠很好的保留下“小尖峰”，說明該算法能保留一定的細(xì)節(jié)信息，不會(huì)出現(xiàn)過渡扼殺的現(xiàn)象。

表1 各算法降噪后比較

3.2 實(shí)測(cè)脈沖星累積輪廓信號(hào)的降噪

本實(shí)驗(yàn)是選取脈沖星 J0437-4715的輻射信號(hào)進(jìn)行實(shí)驗(yàn)，該脈沖星信號(hào)來自歐洲脈沖星網(wǎng)絡(luò)EPN數(shù)據(jù)庫（The European Pulsar Network Data Archive）。如圖3（a）是未經(jīng)消色散處理的該脈沖星累積輪廓信號(hào)。圖3（b）是經(jīng)消色散處理后的脈沖星累積輪廓含噪信號(hào)。圖3（c）是將消色散處理后的脈沖星含噪信號(hào)在小波基“db6”下展開，并進(jìn)行軟閾值消噪處理后的結(jié)果。圖3（d）是將消色散處理后的脈沖星含噪信號(hào)經(jīng)模糊閾值小波降噪算法處理后的結(jié)果。

從圖 3（c）可以看出小波“db6”軟閾值法消噪后的信號(hào)過于平滑，損失掉過多的細(xì)節(jié)信息，如部分脈沖頂過寬，勢(shì)必會(huì)引起周期測(cè)量不準(zhǔn)，且消噪后信號(hào)幅度減小。而模糊閾值小波降噪算法能避免這些情況。

圖3 （a）未經(jīng)消色散處理的脈沖星J0437-4715信號(hào)

圖3 （b）含噪脈沖星累計(jì)輪廓信號(hào)

圖3 （c）小波“Db6”軟閾值法消噪后的信號(hào)

圖3 （d）模糊閾值小波降噪算法消噪后的信號(hào)

4 結(jié)論

本文結(jié)合模糊理論，提出了一種模糊閾值小波降噪算法。實(shí)驗(yàn)結(jié)果表明：與小波“db6”和小波“sym4”軟閾值消澡算法相比，該算法能有效的提高脈沖星輻射信號(hào)的信噪比，能保留較多的有用細(xì)節(jié)信息，避免脈沖累積輪廓過于光滑而引起脈沖頂過寬。

[1] 孫景榮, 許錄平, 梁逸升 等. 中心差分Kalman 濾波方法在X射線脈沖星導(dǎo)航中的應(yīng)用[J]. 宇航學(xué)報(bào), 2008,29（6）:1829-1833.

[2] 仲崇霞,楊廷高. 小波域中的維納濾波在綜合脈沖星時(shí)算法中的應(yīng)用[J]. 物理學(xué)報(bào)，2007，56（10）:6157-6163.

[3] 徐軒彬,吳鑫基.脈沖星 PSR2111+46平均脈沖分析和譜特性[J]. 中國科學(xué)（A輯）,2002,32（12）: 134-1141.

[4] Sheikh S I. The use of variable celestial x-ray sources for spacecraft navigation[D]. Maryland University PhD, 2005.

[5] 閻迪, 許錄平, 謝振華. 脈沖星信號(hào)的模糊閾值小波降噪算法[J]. 西安交通大學(xué)學(xué)報(bào), 200741（10）:1193-1196.

[6] 周懷來,李錄明,羅省賢,李 枚.基于小波變換的改進(jìn)閾值函數(shù)自適應(yīng)去噪方法[J].物探化探技術(shù),2006,28（2）:173-177.