胡渴望 ， 吳根水 ， 屠 寧 ， 王艷奎

（1.中國空空導(dǎo)彈研究院 河南 洛陽 471009；2.航空制導(dǎo)武器航空科技重點(diǎn)實(shí)驗(yàn)室，河南 洛陽 471009）

飛行數(shù)據(jù)對于導(dǎo)彈性能的分析、改進(jìn)以及仿真模型的修正具有重要意義，由于外場試驗(yàn)環(huán)境的復(fù)雜，所測得的飛行數(shù)據(jù)往往具有噪聲。飛行數(shù)據(jù)本身具有非平穩(wěn)的特點(diǎn)，且導(dǎo)彈完整的飛行過程又分為中制導(dǎo)、末制導(dǎo)等階段，不同制導(dǎo)階段飛行數(shù)據(jù)具有明顯不同的變化特征，這使得傳統(tǒng)的降噪方法作用非常有限。

小波變換不同于傳統(tǒng)分析方法把時域信號表示為若干精確的頻率分量之和，而是將其表示為若干描述子頻帶的時域分量之和。小波分析具有自適應(yīng)的時—頻窗，可以聚焦到被分析信號的局部時域位置和局部頻帶，具有在時域、頻域表征信號局部特征的能力，小波降噪中對信號的多尺度分解，即通過不同中心頻率不同帶寬的濾波器對信號進(jìn)行濾波，然后根據(jù)噪聲和信號在各頻帶上的不同特性對小波系數(shù)進(jìn)行處理，再重構(gòu)實(shí)現(xiàn)信號的降噪[1]。

目前小波降噪方法主要有3種，即模極大值降噪、閥值降噪和空域相關(guān)濾波算法[2-3]。模極大值降噪基于信號和噪聲有不同的奇異性，提取信號在各個尺度上的模極大值，由信號模極大值重構(gòu)實(shí)現(xiàn)降噪。其優(yōu)點(diǎn)是無需估計噪聲水平，但缺點(diǎn)是模極大值的搜索困難、且投影算法過于復(fù)雜，計算量太大；而基于閥值的小波降噪由于其計算量小、易于實(shí)現(xiàn)得到了廣泛的應(yīng)用。但不論采用哪種閥值選取準(zhǔn)則，由閥值收縮后保留下來的小波系數(shù)，與原信號的小波系數(shù)存在一定偏差，使得降噪后信號過平滑；空域相關(guān)濾波算法利用相鄰尺度間小波系數(shù)的相關(guān)性來進(jìn)行濾波，需要設(shè)定處理相關(guān)系數(shù)的閥值，且可能會保留由噪聲引起的變換點(diǎn)，使得重構(gòu)的信號出現(xiàn)“毛刺”[4]，但空域相關(guān)濾波的最大優(yōu)點(diǎn)是能保留信號尖峰等細(xì)節(jié)信息。

空域相關(guān)濾波算法的應(yīng)用有兩個難點(diǎn)：即難以確定小波的分解層數(shù)和濾波閥值。這限制了空域相關(guān)濾波算法實(shí)現(xiàn)效果。本文將基于統(tǒng)計理論的白噪聲檢驗(yàn)方法[5]和相關(guān)空域?yàn)V波相結(jié)合，自適應(yīng)的確定最佳小波分解層數(shù)，并利用新的閥值準(zhǔn)則，使得信號的小波系數(shù)的保留更加精確，通過仿真和實(shí)測數(shù)據(jù)試驗(yàn)表明，本文算法實(shí)現(xiàn)了更好的濾波效果。

1 空域相關(guān)濾波

信號的突變點(diǎn)在不同分解尺度的同一位置都有較大的峰值出現(xiàn)，而噪聲并無此特性，基于此Xu提出了空域相關(guān)濾波算法[6]。即取相鄰尺度的小波系數(shù)的直接乘積進(jìn)行相關(guān)計算，根據(jù)相關(guān)性大小區(qū)分尺度上的噪聲和信號。這樣可以更好地保留信號邊緣信息。

空域相關(guān)濾波概括如下：對含噪信號進(jìn)行小波分解，設(shè)cd（l，n）為觀測數(shù)據(jù)在l尺度上n位置的小波系數(shù)，則相鄰尺度間相關(guān)系數(shù)為

為了使相關(guān)系數(shù)與小波系數(shù)具有可比性，將corr（l，n）的能量歸一化到 cd（l，n），即

上式中：

通過比較 newcorr（l，n）和 cd（l，n）的大小來判斷該位置小波系數(shù)是否信號能量占優(yōu)，即若 newcorr（l，n）≥cd（l，n），認(rèn)為 n 位置小波系數(shù)為信號系數(shù)，提取 cd（l，n）并保存為 cd′（l，n），將原位置的 cd（l，n）置零。再利用剩余的 cd（l，n）計算新的newcorr′（l，n）、Pw′、Pc′，再重復(fù)進(jìn)行上述計算，不斷更新cd′（l，n），直到該尺度上未被抽取的點(diǎn)的能量Pwl滿足該尺度上的噪聲能量閥值停止迭代計算。最后根據(jù)各個尺度上新的小波系數(shù)cd′（l，n）重構(gòu)濾波后的信號。

上述傳統(tǒng)空域?yàn)V波算法有不足之處：分解到哪個尺度上可得到最佳濾波效果全靠經(jīng)驗(yàn)；對于各尺度上認(rèn)定為噪聲能量占優(yōu)的小波系數(shù)完全舍棄，使得重構(gòu)信號在很多點(diǎn)失真較為嚴(yán)重；各尺度上噪聲能量閥值的設(shè)定會對濾波效果帶來很大的影響，如何精確設(shè)定閥值也稱為關(guān)鍵。

2 改進(jìn)的空域相關(guān)濾波

2.1 白噪聲檢驗(yàn)

為確定相關(guān)空域?yàn)V波的分解停止尺度，引入白噪聲檢驗(yàn)方法[4]，在此先做介紹，以便后文使用。由于有用信號和噪聲的非關(guān)聯(lián)性，二者各自的自相關(guān)函數(shù)也具有非關(guān)聯(lián)性，故噪聲的自相關(guān)函數(shù)反映了噪聲的分布情況。若小波系數(shù)序列為X=（x1，x2，…xi，…xN），定義噪聲自相關(guān)函數(shù)為：

定義序列

則判定小波系數(shù)序列X是否為噪聲序列就等價于判斷序列P的z個變量是否獨(dú)立并服從標(biāo)準(zhǔn)正態(tài)分布。若定義

則判斷序列P的z個變量是否獨(dú)立并服從標(biāo)準(zhǔn)正態(tài)分布就等價于判斷Q是否服從自由度為m的χ2分布，此處z取5。

2.2 改進(jìn)算法

信噪比低的信號，其小波分解最高尺度上往往噪聲能量絕對占優(yōu)，甚至完全為噪聲，這種情況下進(jìn)行相關(guān)空域?yàn)V波已失去意義，本文對分解所得的小尺度上的小波系數(shù)，在進(jìn)行相關(guān)空域?yàn)V波前，先進(jìn)行白噪聲檢驗(yàn)，若判斷是白噪聲，則全部舍去，再對下一層進(jìn)行檢驗(yàn)，直至某一層未通過白噪聲檢驗(yàn)，則認(rèn)為大于等于該層的尺度上，含有足夠的信號信息，應(yīng)進(jìn)行相關(guān)空域?yàn)V波。

對于各層噪聲能量閾值的設(shè)定，文獻(xiàn)[7]提出的作為停止相關(guān)性迭代計算的判定條件為：

其中nl為未被提取的小波系數(shù)的個數(shù)，cth為放大系數(shù)，于1～2之間依經(jīng)驗(yàn)取值。因?yàn)閘尺度上噪聲標(biāo)準(zhǔn)差σl由估計而得，閾值設(shè)定有誤差，使得相關(guān)空域?yàn)V波難以有效實(shí)施。本文為了更精確地提取信號系數(shù)，從相關(guān)性濾波的目的出發(fā)，經(jīng)過提取后剩余的小波系數(shù)，若噪聲能量絕對占優(yōu)，則較好的分離了信號占優(yōu)系數(shù)和噪聲占優(yōu)系數(shù)。而噪聲占優(yōu)系數(shù)大多處于同一能量級別，若設(shè)定閥值在臨界值變動時，則提取的小波系數(shù)的個數(shù)會有極其明顯的變化，而設(shè)定閾值在除臨界值外的相當(dāng)大的范圍內(nèi)變動時，提取的系數(shù)個數(shù)并無顯著變化。故本文采取由結(jié)果導(dǎo)向來確定閾值的新方法：在 σl和 2σl之間等間隔設(shè)置 m 個閾值（t1≤t2≤…≤tm），依此分別進(jìn)行相關(guān)性計算，得到m組cd（l，n），每組個數(shù)為Nm，由公式（8）可知，不等式右側(cè)閾值越大，提取的信號占優(yōu)小波系數(shù)個數(shù)越少，即Nm≤Nm-1。比較兩相鄰閾值提取的信號占優(yōu)小波系數(shù)的個數(shù)Nm，在某一組兩相鄰預(yù)設(shè)閾值情況下，若信號占優(yōu)系數(shù)的個數(shù)（Nm-1-Nm）變化最大時，則確定兩閾值中較大的閾值為最優(yōu)閾值，并確定提取出來的信號占優(yōu)小波系數(shù)為

本文中m取6。

針對分解到哪一尺度停止相關(guān)空域?yàn)V波的問題，對l尺度上提取 cd′（l，n）過之后剩余的小波系數(shù)序列（設(shè) cd″（l，n）為）非零項(xiàng)進(jìn)行相關(guān)性檢驗(yàn)，若判斷為白噪聲，則繼續(xù)對下一尺度上進(jìn)行相關(guān)性濾波，否則說明該尺度進(jìn)行相關(guān)濾波后剩余的小波系數(shù)已達(dá)到信號占優(yōu)，不宜通過相關(guān)性濾波進(jìn)行區(qū)分，上一尺度L即為最終相關(guān)濾波尺度。

如前所述，僅提取能量占優(yōu)位置小波系數(shù)，則舍棄了真實(shí)信號在其他位置上的系數(shù)，為了將未提取位置上的噪聲能量占優(yōu)系數(shù)調(diào)整到接近理想信號，考慮到不同尺度的噪聲方差沿分解尺度上近似為指數(shù)分布[5，8]，故取指數(shù)函數(shù)為萎縮函數(shù)，則最終的小波重構(gòu)系數(shù)算法為：

綜上所述，本文的濾波步驟如下：

1）對含噪信號進(jìn)行小波分解，得到近似系數(shù)和細(xì)節(jié)系數(shù)；

2）對前一步得到的細(xì)節(jié)系數(shù)用公式（7）進(jìn)行檢驗(yàn)，若cd（1）判斷為白噪聲，則cd （1）全賦零值，并對近似部分繼續(xù)進(jìn)行小波分解和判斷，直至某一層細(xì)節(jié)系數(shù)判斷為非白噪聲，進(jìn)入下一步；

3） 在預(yù)置的 6 個閥值 σl、1.2σl、1.4σl、1.6σl、1.8σl、2σl情況下，根據(jù)公式（1）、（2）進(jìn)行相關(guān)空域?yàn)V波，最后根據(jù)公式（9）確定最優(yōu)閥值，提取信號占優(yōu)小波系數(shù)；

4）對 l尺度剩余的小波系數(shù)序列 cd″（l）非零項(xiàng)用公式（7）進(jìn)行相關(guān)性檢驗(yàn)，若判斷為白噪聲，則在下一尺度上繼續(xù)進(jìn)行相關(guān)性濾波，直至分解到L+1尺度上，該尺度cd″（L+1）非零項(xiàng)判斷為非白噪聲，則停止分解，并確定L為最終分解尺度；

5）各尺度小波系數(shù)按式（10）進(jìn)行構(gòu)造，并進(jìn)行小波分解逆變換重構(gòu)信號。

3 試驗(yàn)結(jié)果與分析

為驗(yàn)證本文算法的有效性，分別采用仿真數(shù)據(jù)和實(shí)測數(shù)據(jù)進(jìn)行分析。仿真數(shù)據(jù)來自某彈道條件下全數(shù)字仿真。對產(chǎn)生的仿真數(shù)據(jù)施加均值為0標(biāo)準(zhǔn)差為5、10、15、20、25的高斯噪聲，仿真結(jié)果分別見表1和圖1所示。

表1 不同信噪比降噪結(jié)果SNR比較/dbTab.1 The denoised results of different SNR signals

由于篇幅的限制，在此僅給出噪聲標(biāo)準(zhǔn)差為10的仿真結(jié)果，如圖1所示。

圖1 噪聲標(biāo)準(zhǔn)差10時仿真結(jié)果Fig.1 The result of signal with 10 standard deviation white noise

由上圖可知，文中算法較好地去除了噪聲干擾，有效的提高了信噪比。

取某彈道外場飛行試驗(yàn)方位角數(shù)據(jù)，分別采用minimax準(zhǔn)則的可調(diào)整軟閥值小波降噪算法和本文算法進(jìn)行分析，結(jié)果如圖2所示，本文算法可有效改進(jìn)數(shù)據(jù)質(zhì)量，效果良好。

由圖2可知，文中提出的降噪方法，較之minimaxi閥值準(zhǔn)則更完善的保留了信號尖峰等細(xì)節(jié)信息，又克服了rigrsure閥值準(zhǔn)則過于保守的缺陷（去噪后依然含有過多噪聲），提高了信號的信噪比，改善了數(shù)據(jù)質(zhì)量。

4 結(jié) 論

本文針對飛行數(shù)據(jù)難以有效去噪這一實(shí)際問題，從相關(guān)空域?yàn)V波的目的出發(fā)，給出了一種空域相關(guān)濾波閥值確定方法，將小波系數(shù)分為能量占優(yōu)和噪聲占優(yōu)部分，對噪聲占優(yōu)系數(shù)進(jìn)行了萎縮重構(gòu)，并利用噪聲能量占優(yōu)系數(shù)來自適應(yīng)的確定相關(guān)空域?yàn)V波停止層數(shù)。給出了相關(guān)空域?yàn)V波的實(shí)用去噪流程，通過仿真信號和實(shí)測信號表明：本文提出的方法使得信噪比提高的同時完善的保留了信號的細(xì)節(jié)信息，提升了信號質(zhì)量。

圖2 飛行試驗(yàn)數(shù)據(jù)降噪結(jié)果Fig.2 The denoised result of real flight signal

[1]謬玲娟.小波分析在光纖陀螺信號濾波中的應(yīng)用研究[J].宇航學(xué)報，2000，21（1）:42-46.MIAO Ling-juan.Application of wavelet analysis in the signal processing of the fiber optic gyro[J].Journal of Astronautics.2000，21（1）:42-46.

[2]喬玉坤.小波域軟閥值降噪的地磁匹配導(dǎo)航仿真研究[J].兵工學(xué)報，2011，32（9）:1083-1088.QIAO Yu-kun.Simulation research on geomagnetic matching navigation based on soft-threshold wavelet denoising method[J].Acta Armamentarii，2011，32（9）:1083-1088.

[3]程海林.光譜測量數(shù)據(jù)降噪處理濾波器設(shè)計[J].天津理工大學(xué)學(xué)報，2012，28（1）:45-50.CHENG Hai-lin.Noise reducing filter design for the spectra data processing[J].Journal of Tianjin University of Technology，2012，28（1）:45-50.

[4]吳永超，馬秀榮，張媛.一種改進(jìn)的小波系數(shù)相關(guān)濾波算法[J].信息與電子工程，2012，10（5）；604-607.WU Yong-chao，MA Xiu-rong，ZHANG Yuan.An improved denoisingalgorithmusing correlation of wavelet coefficients[J].Information and Electronic Engineering，2012，10（5）:604-607.

[5]王斐，梁曉庚，鄭強(qiáng)，等.飛行數(shù)據(jù)改進(jìn)中值濾波和自適應(yīng)小波降噪[J].計算機(jī)仿真，2013，（7）:75-77.WANG Fei，LIANG Xiao-geng，ZHENG Qiang，et al.Flight data denoising algorithm using improved median filter and adaptive threshold shrinkage in wavelet transform[J].Computer Simulation，2013，（7）:75-77.

[6]王國林.小波分析在路面不平度信號降噪中的應(yīng)用研究[J].機(jī)械設(shè)計與制造，2011，4（4）；65-67.WANG Guo-lin.Research on wavelet application in road roughness denoising[J].Machinery Design and Manufacture，2011，4（4）:65-67.

[7]彭云輝，劉云峰.基于空域相關(guān)的激光陀螺信號濾波方法研究[J].紅外激光與工程，2007，36（4）:493-496.PENG Yun-hui，LIU Yun-feng.Spatial correlation denoising for signal of laser gyro[J].Infrared Laser Engineering，2007，36（4）:493-496.

[8]Yuan X H，Buckles B P.Subband noise estimation for adaptive wavelet shrinkage[C]//ICPR 2004 proceedings of the 17th International Conference Pattern Recognition，2004:858-888.