基于提升靜態小波變換的風廓線雷達地雜波抑制方法

張雯雯，劉黎平，阮征，葛潤生

(中國氣象科學研究院 災害天氣國家重點實驗室，北京，100081)

風廓線雷達探測技術發展于20世紀60年代，它以凈空湍流作為探測目標，但由于空氣的運動速度比較小，產生的多普勒頻移也較小，而且大氣散射回波信號一般都很弱，所以，很容易受到地物返回信號等低頻雜波的干擾，形成地雜波[1]。當湍流回波與地雜波同處于一個脈沖寬度內時，強地雜波有可能淹沒湍流回波，這給信號檢測帶來不便。尤其是在監測大氣邊界層時．風廓線雷達回波受到地雜波的影響更為明顯，可能導致嚴重的多普勒頻率偏移[2]。因此，抑制地雜波一直是風廓線雷達需要解決的問題。國內的風廓線雷達地雜波抑制方法有非相參 MTI對消器、IIR濾波器和頻域濾波法等[3]。前 2種濾波器設備量要求大，會產生多普勒頻率偏移，因此，限制了它的推廣和應用。而頻域濾波法則是考慮到地物雜波在頻域中的特點提出的，它包括零頻附近回波譜對稱消法[4]、NOAA風廓線雷達業務運行網采用的均值替代中心對稱譜法[5]等。這些方法對于風速較大時(地物雜波譜與大氣譜分離的情況下)雜波的抑制效果較好；而在風速較小(雜波譜與大氣回波譜交疊)情況下，往往容易將零頻附近的部分大氣譜消去，導致多普勒頻移計算產生偏差，從而影響了雷達的測量精度。Jordan等[6]提出了基于小波變換去地物雜波干擾的方法。把經過相干積分后的時域信號變換到小波域中，在小波域中可以區分雜波和晴空大氣回波，將雜波成分去除，再將信號重構即得到已抑制雜波的大氣回波。但在實際應用中發現，這個方法有以下幾個缺點：難以給出確定的信號分解層數；小波系數截取閾值的設定不理想；由于采用的是離散小波變換，不具有平移不變這一特性，在重構時會帶來較大的重構誤差，所以，在大部分情況下去地雜波效果不明顯。針對這些現象，本文作者提出了一種基于提升靜態小波變換(SWT)的去風廓線雷達地雜波技術，并通過理論推導和仿真實驗驗證此方法的合理性及有效性。

1 風廓線雷達回波信號特征分析

風廓線雷達信號包括了大氣湍流的散射回波、均勻背景噪聲和或多或少的地雜波和間歇性雜波。由于本文主要研究地雜波信號的去除，所以，暫時不考慮間歇性雜波。

地雜波是雷達附近的靜止目標如樹、建筑物或電力線等散射的信號進入雷達天線旁瓣造成的，地雜波產生的多普勒功率譜在零頻移附近[3]。當觀察雷達接收到回波信號時，可以看到回波的多普勒功率譜不僅僅包含具有均勻噪聲背景的凈空信號的譜峰，還附加有雜波譜峰，而且這種雜波譜很有可能與凈空信號交疊在一起。地雜波能量常常比湍流回波信號能量高幾個數量級，當雜波譜與凈空信號譜交疊時，就會淹沒凈空信號，在這種情況下，難以區分多普勒功率譜中凈空信號和雜波。這里假設地雜波和大氣湍流的散射回波的多普勒功率譜都服從高斯分布，即滿足：

其中：σf為多普勒頻率譜寬；f為多普勒平均頻率，通常對于地雜波信號，f=0。

2 SWT的提升分解實現

2.1 DWT的提升分解實現

DWT的多相[7]表示如圖1所示。

圖1 一次離散小波分解和重構的多相表示Fig.1 Polyphase representation of decomposition and reconstruction of DWT

利用Laurent多項式的Euclidean算法將P(z)分解為有限步提升來實現，如式(1)所示。

式中：P(z)為母小波的多相表達式；R為常數；si(z)和ti(z)為分解表達式。

傳統小波一次分解的提升框架實現結構見圖2。

圖2 一次分解的提升框架實現結構Fig.2 Lifting scheme frame representation of decomposition

2.2 SWT的提升分解算法

DWT和SWT的區別僅僅在于SWT在分解端沒有對信號進行下采樣[8]，重構端也沒有上采樣，所以，將圖1中分解端的下采樣算子與、重構端的上采樣算子與P(z)等效易位，再分別去掉分解端的下采樣算子與重構端的上采樣算子，則可得 SWT的多相表示，如圖3所示。

由式(2)和式(3)得：

圖3 SWT的多相表示Fig.3 Polyphase representation of SWT

3 基于提升SWT的地雜波抑制方法

3.1 小波基的選取和最佳分解層數的確定

濾波的第1步是選擇合適的小波基，使得信號能量能分布在少數幾個基底上。通常為了兼顧濾波的實時性和效果，希望所選的小波能同時具有下列性質[7]：① 對稱性或反對稱性，以避免信號失真；② 較短的支撐，以減小運算時間；③ 正交性；④ 較高的消失矩，以更好地匹配待分析的信號。

然而，實際上同一個小波基不可能同時具有上述所有性質，就性質①來說，在所有小波基中，僅僅只有Haar小波是滿足條件的，但它由于消失矩階數過低而阻止了其在實際中的應用；而且較短的支撐和較高的消失矩本身就是一對矛盾，所以，在選取小波基時應折衷地考慮上述性質。考慮到db小波的高階原點矩等于零，并綜合計算量以及上述因素，這里選擇 db3小波。其一次SWT分解的提升實現如式(7)所示：

式中：a=-2.425 497 25；b=-0.079 339 46；c=0.352 387 66；d=2.895 347 45；e=-0.561 414 91；f=0.019 750 53；g=2.315 458 04；h=0.431 879 99。

在提升框架中，根據實際情況分別對信號進行延拓處理，采取對稱處理的延拓方式。

為了能夠更好地去除地雜波信號，最好能夠使得分解后信號的最低層低頻系數中不僅包含地雜波信號的所有成分，而且還盡可能少地包含回波信號成分，也即能夠保證將地雜波信號和回波信號很好地分離的同時，使得地雜波信號集中在小波分解域的最低層，以便于將其集中去除。一般來說，地雜波信號強度較回波信號強很多，所以，其經分解后在各層小波系數上的值也要大很多，為了確定最佳分解層數，可以計算信號在各層小波系數上的2階原點矩：

式中：D為2階原點矩的值；Xji為第j層分解的小波系數；N為SWT點數。

依次增加信號的分解層數，并計算各層系數的 2階原點矩，找出能夠使信號能量剛好在最低層發生明顯躍變的最大分解層數，即為最佳分解層數。表1給出了對信號進行不同層的 SWT處理后，其各自對應的每一層的原點矩。

表1 對信號進行不同層的SWT處理后各層對應的D Table 1 D after SWT to different layers

這里，a1～a5和 d1～d5分別表示信號經過 SWT分解后的低頻和高頻系數。從表1可以看出：當信號分解層數小于4層時，雖然能保證地雜波信號能量均集中在其對應的最低層低頻系數上，但是，由于分解層數過少，所以，不能很好地將回波信號與地雜波信號分離，也就是說，其最低層也含有大量的有用回波。為了達到信雜分離的目的，就只能增加分解層數，但是，如果對其分解過多，如進行5層小波分解，則又會由于分解層數過多而導致雜波信號被同時分解到信號最低層的高頻和低頻系數上，這樣，若繼續利用3.2節方法來確定閾值的話，則很可能會導致閾值過大而保留過多的地雜波，甚至出現錯誤的結果，同時也會造成計算量增大。所以，當分解層數為4層時效果最佳，此時，不僅能最大程度地將地雜波信號和回波信號分解在不同的小波層上，而且可以保證地雜波信號均在第4層低頻系數上，有利于后續的處理。

3.2 小波閾值的確定與小波系數的截取

去除地雜波就是將小波域中已經確定為地雜波成分的大值系數予以截取，但在以往的研究中，都是選取后一半小波系數的最大值作為截取閾值，且這種處理針對所有小波系數是不科學的。因為若信號回波頻率也比較低，則可能出現后一半小波系數都是噪聲的情況；若繼續取后一半系數的最大值作為截取閾值，則可能會出現由于閾值過小，而使得有用信號也被濾除，導致信號失真的現象。而即便后一半小波系數中也包含回波信號，由于截取閾值取的是后一半系數的最大值，若回波信號比較大，即這個最大值比較大，也可能會由于閾值過大而保留了過多的地雜波，導致信號失真。

這里選取了最鄰近地雜波頻段的2個頻段的平均值作為閾值，這是因為最鄰近的那幾層所含信號頻段與它最類似，可以最真實的還原最低層信息。且由于經過3.1節的處理，雜波信號僅存在于SWT的最低層低頻信號中；所以，這里僅針對 SWT的最低層系數進行處理即可，這樣可以更好地濾除地雜波信號。

4 仿真分析

根據上述方法，對風廓線雷達信號進行仿真，仿真均是在信雜比為-50 dB，信噪比為20 dB，f=200 Hz，σf=62.5，c=0.01，發射機重復頻率為1 280 Hz的條件下進行。這里的雜波僅考慮地雜波情況，并假設其σf′=12.5。圖4所示為對信號進行SWT處理前后的對比。從圖4(a)和(b)可以看出：由于地雜波信號很大，風廓線雷達的回波信號已經完全被其掩蓋，此時，信號的特征主要表現為地雜波的特征，若對其進行SWT處理，使得地雜波信號均集中在SWT的最低層的低頻信號上，并選擇合適的閾值對其濾波，此時可以將大部分地雜波信號剔除，并還原風廓線雷達回波信號(圖 4(c)和 4(d))。

圖5所示為信號在相同情況下進行DWT處理后小波系數和去雜波后的信號的頻譜。從圖5可以看出：經DWT重構后的信號會發生頻譜的折疊，而且由于其處理的點數較少，并有少量失真，其去雜波效果明顯不如 SWT的去雜效果，這也說明了本文方法的優越性。

圖4 信號在進行SWT處理前后的對比Fig.4 Comparison chart before and after SWT

圖5 信號進行DWT后小波系數和去雜波后的頻譜Fig.5 Wavelet coefficients after DWT and de-clutter spectrum

5 結論

根據DWT的提升實現方法，給出了SWT的提升實現方法，并將其應用于風廓線雷達信號的地雜波去除，提出了一種新的確定最優分解層數和閾值的方法；同時，通過大量的實驗結果和圖例驗證了這種方法的有效性。相對于傳統的基于離散小波變換去地雜波的方法，本文提出的方法具有以下幾個特點：

(1) 采用SWT代替了DWT，克服了DWT不具有平移不變性的缺點，提高了信號的重構精度,也使得其能應用于更低信雜比的情況。

(2) 與傳統的基于DWT去地雜波方法不同的是，SWT的提升實現方法是一種明確有效的確定最優分解層數和閾值的方法，避免了實踐中由于分解層數或閾值選取不恰當而造成算法失效的問題。

(3) 由于提升小波變換本身運算簡單,運算量比較小，所以，這種方法在大大提高了算法性能的同時，其計算量相對于傳統的基于離散小波變換去地雜波的方法并沒有很大的提高，易于實現。

[1] 何平. 相控陣風廓線雷達[M]. 北京: 氣象出版社, 2006:23-70.HE Ping. The phased array profiler radar[M]. Beijing:Meteorological Press, 2006: 23-70.

[2] 丁敏, 黃登才, 宋金雷. 小波變換去風廓線雷達地雜波技術[J]. 中南大學學報: 自然科學版, 2007, 38(1): 949-953.DING Min, HUANG Deng-cai, SONG Jin-lei. Removing ground clutter from wind profiler data using wavelet[J]. Journal of Central South University: Science and Technology, 2007,38(1): 949-953.

[3] 劉艷. 多普勒天氣雷達地物雜波時域和頻域抑制研究[D]. 成都: 成都信息工程學院電子工程系, 2007: 33-40.LIU Yan. The research on ground clutter suppression of Doppler weather radar in time domain and frequency domain[D].Chengdu: Chengdu University of Information Technology.Electron Engineering Department, 2007: 33-40.

[4] Passarelli R E, Romanik’P, Geotis S G, et al. Ground clutter rejection in the frequency domain[C]//Proceedings of 20th Conf on Radar Meteorology. Boston, MA: Amer Meteor Soc, 1981:295-300.

[5] Barth M F, Chadwick R B, van de Kamp D W. Data processing algorithms used by NOAA’s wind profiler demonstration network[J]. Ann Geophys, 1994, 12(5): 518-528.

[6] Jordan J R, Lataitis R J, Carter D A. Removing ground and intermittent clutter contamination from wind profiler signals using wavelet transforms[J]. J Atmos Oceanic Techonl, 1997,14(3): 1280-1297.

[7] 成禮智, 王紅霞, 羅永. 小波的理論與應用[M]. 北京: 科學出版社, 2004: 121-133.CHENG Li-zhi, WANG Hong-xia, LUO Yong. The theory and application of Wavelet[M]. Beijing: Science Press, 2004:121-133.

[8] 孟晉麗. 基于鄰域相關性的小波域濾波算法研究[D]. 西安:西北工業大學自動化系, 2006: 87-95.MENG Jin-li. Wavelet denoising based on adjacent dependencies[D]. Xi’an: Northwestern Polytechnical University.Department of Automation Control, 2006: 87-95.