基于多延遲相關(guān)參數(shù)的信號處理算法

， ,2， ，

(1.解放軍理工大學(xué)氣象海洋學(xué)院， 江蘇南京 211101；2.南京氣象雷達(dá)開放實(shí)驗(yàn)室， 江蘇南京 210008)

0 引言

天氣雷達(dá)是大氣探測的重要工具之一，其回波數(shù)據(jù)是進(jìn)行災(zāi)害性天氣預(yù)警和短時(shí)天氣預(yù)報(bào)的重要參考資料[1]。但是雷達(dá)回波信號中總是摻雜有部分噪聲干擾信號，有時(shí)在回波信噪比較低的情況下，噪聲信號甚至?xí)把蜎]”有用的降水回波信號。這嚴(yán)重制約回波信號的譜矩估計(jì)質(zhì)量，影響雷達(dá)輸出回波數(shù)據(jù)的應(yīng)用潛力[2]。因此，對天氣雷達(dá)信號處理算法進(jìn)行改進(jìn)，以期消除噪聲信號對估計(jì)得到的回波數(shù)據(jù)質(zhì)量的影響，是天氣雷達(dá)研究領(lǐng)域非常有意義的一件工作。

天氣雷達(dá)信號處理發(fā)展過程中，曾出現(xiàn)過多種信號處理算法[3]，包括快速傅里葉變換(FFT)算法、矢量相位變換(VPC)算法、標(biāo)量相位變換(SPC)算法、脈沖對處理(PPP)算法等。但在天氣雷達(dá)實(shí)際信號處理中，不僅要考慮處理得到的回波數(shù)據(jù)質(zhì)量，另外為了保證天氣雷達(dá)能夠?qū)崟r(shí)處理，還要兼顧所選信號處理算法的運(yùn)算速度。因此，目前天氣雷達(dá)信號處理器中脈沖對處理算法應(yīng)用最為廣泛。在實(shí)際應(yīng)用過程中，在弱回波區(qū)域，脈沖對處理算法得到的速度譜寬數(shù)據(jù)經(jīng)常會(huì)存在一些異常數(shù)據(jù)，這嚴(yán)重影響回波數(shù)據(jù)的正常使用[4]。因此，為提高天氣雷達(dá)數(shù)據(jù)質(zhì)量，本文對傳統(tǒng)信號處理算法進(jìn)行改進(jìn)，引入多階延遲自相關(guān)參數(shù)。

1 傳統(tǒng)信號處理算法介紹

文中信號處理主要是指時(shí)域I/Q(即同相/正交)數(shù)據(jù)處理得到譜矩參數(shù)的過程，其中時(shí)域I/Q數(shù)據(jù)是目標(biāo)回波信號中一組相關(guān)的時(shí)域采樣脈沖序列，如式(1)所示：

sn=In+jQn，n= 1,2,3,…,N

(1)

式中，In和Qn分別表示I路和Q路信號，N表示采樣脈沖的數(shù)量。

目前天氣雷達(dá)信號處理器主要使用脈沖對處理(Pulse Pair Processing, PPP)算法，其主要是根據(jù)回波信號的自相關(guān)參數(shù)進(jìn)行計(jì)算[5-6]。主要處理過程為：首先，對時(shí)域I/Q數(shù)據(jù)進(jìn)行濾波，消除地物雜波信號；然后，根據(jù)濾波后的時(shí)域I/Q數(shù)據(jù)計(jì)算各自相關(guān)參數(shù)，如式(2)所示；最后，根據(jù)自相關(guān)參數(shù)計(jì)算回波譜矩參數(shù)，如式(3)所示：

(2)

(3)

式中，R(mTs)表示m階自相關(guān)參數(shù)，Ts表示脈沖重復(fù)周期，N表示平均噪聲功率， *表示復(fù)數(shù)的共軛。

2013年，Ivic和Torres等對美國新一代天氣雷達(dá)網(wǎng)中WSR-88D多普勒天氣雷達(dá)平均噪聲功率作了系統(tǒng)研究，結(jié)果表明，同一部雷達(dá)在正常工作情況下，不同時(shí)間、不同方位、不同仰角以及不同徑向距離，其平均噪聲功率通常存在較大差異，有時(shí)甚至超過2 dB[7-8]。然而，目前應(yīng)用在業(yè)務(wù)觀測中的天氣雷達(dá)都是以較高仰角時(shí)晴空回波功率作為雷達(dá)平均噪聲功率，應(yīng)用該平均噪聲功率輸出的回波數(shù)據(jù)質(zhì)量必將受到影響，尤其是回波信噪比較低的區(qū)域。當(dāng)設(shè)置的平均噪聲功率低于實(shí)際值時(shí)，速度譜寬估值偏大；反之，速度譜寬估值偏小，甚至?xí)沟没夭ㄗV寬數(shù)據(jù)區(qū)域縮小(當(dāng)R(0)

2 信號處理算法改進(jìn)

Janssen和Spek采用5.56 GHz相控陣?yán)走_(dá)對氣象目標(biāo)進(jìn)行了長期觀測，并對大量回波信號進(jìn)行分析，研究結(jié)果表明：超過75%的回波信號功率譜都具有高斯形狀[9]。因此，實(shí)際應(yīng)用中通常將雷達(dá)接收的天氣回波信號頻譜近似為高斯分布，并用高斯函數(shù)表示天氣回波信號功率譜S(v)[10]：

(4)

由于功率譜S(f)與自相關(guān)函數(shù)R(l)互為傅里葉變換，因此，根據(jù)離散傅里葉變換(DFT)可以計(jì)算回波信號延遲mTs時(shí)相應(yīng)的自相關(guān)參數(shù)R(mTs)：

(5)

(6)

式中，Ts表示脈沖重復(fù)周期，mTs表示自相關(guān)函數(shù)延遲時(shí)間，δm表示沖擊函數(shù)，當(dāng)且僅當(dāng)m為零時(shí)，自相關(guān)參數(shù)R(mTs)包含噪聲分量。

根據(jù)式(6)可知，在脈沖重復(fù)周期Ts和雷達(dá)波長λ已知的前提下，只需要兩個(gè)不同自相關(guān)參數(shù)就可以推導(dǎo)出速度譜寬σv的估計(jì)值。另外，為了消除傳統(tǒng)信號處理算法中平均噪聲功率對回波數(shù)據(jù)的影響，避免使用R(0)計(jì)算速度譜寬。又考慮到天氣雷達(dá)回波脈沖樣本較少，且脈沖相關(guān)性容易受延遲時(shí)間影響，因此，通常嚴(yán)格控制自相關(guān)函數(shù)的延遲時(shí)間，通常m≤3。本文提出Lag1&2算法(即使用一階延遲自相關(guān)參數(shù)R(Ts)和二階延遲自相關(guān)參數(shù)R(2Ts))和Lag1&3算法(即使用一階延遲自相關(guān)參數(shù)R(Ts)和三階延遲自相關(guān)參數(shù)R(3Ts))兩種速度譜寬估計(jì)算法：

(7)

(8)

式中，σv2和σv3分別表示采用Lag1&2算法和Lag1&3算法得到的速度譜寬估計(jì)值，R(Ts)，R(2Ts)，R(3Ts)分別表示一階、二階、三階延遲自相關(guān)參數(shù)。

3 實(shí)測回波數(shù)據(jù)分析

本文將選取一次強(qiáng)降水過程的回波數(shù)據(jù)進(jìn)行分析，比較傳統(tǒng)信號處理算法與Lag1&2算法、Lag1&3算法得到的速度譜寬數(shù)據(jù)差異。由于是進(jìn)行信號處理算法研究，實(shí)驗(yàn)中所有處理過程都是對時(shí)域I/Q數(shù)據(jù)進(jìn)行處理。

3.1 回波數(shù)據(jù)選取

圖1(a)、(b)分別顯示的是一次強(qiáng)降水過程回波反射率因子圖像和回波平均徑向速度圖像。其中橫縱坐標(biāo)分別表示距離雷達(dá)站的徑向距離，(0,0)位置表示雷達(dá)站所在位置。為濾除噪聲回波，本文對參數(shù)估計(jì)得到的回波數(shù)據(jù)設(shè)置信噪比門限，即當(dāng)回波信噪比(SNR)大于-10 dB時(shí)，輸出該距離庫回波數(shù)據(jù)。

3.2 譜寬數(shù)據(jù)圖像比較

圖2(a)、(b)和(c)分別顯示的是采用傳統(tǒng)信號處理算法(如式(3)所示)、Lag1&2算法(如式(7)所示)、Lag1&3算法(如式(8)所示)處理得到的速度譜寬數(shù)據(jù)圖像。比較上述圖像，可以得出以下結(jié)論：

1) 在弱回波區(qū)域，采用傳統(tǒng)信號處理算法得

到的速度譜寬數(shù)據(jù)更容易受噪聲干擾，如圖2(a)所示，在近雷達(dá)站區(qū)域以及降水回波邊緣區(qū)域(由圖1(a)可知，主要是回波強(qiáng)度較弱區(qū)域)，其速度譜寬估計(jì)值明顯大于圖2(b)、(c)，而且速度譜寬數(shù)據(jù)連續(xù)性較差；

2) 采用Lag1&2算法和Lag1&3算法得到的速度譜寬數(shù)據(jù)較為接近，但是在回波強(qiáng)度較弱的區(qū)域，Lag1&3算法得到的速度譜寬數(shù)據(jù)更加平滑，如圖2(b)、(c)中近雷達(dá)站區(qū)域。

(a)反射率因子數(shù)據(jù)PPI圖像

(b)徑向速度數(shù)據(jù)PPI圖像圖1回波強(qiáng)度和速度數(shù)據(jù)圖像

(a)傳統(tǒng)算法 (b)Lag1&2算法(c)Lag1&3算法圖2速度譜寬數(shù)據(jù)圖像

3.3 譜寬散點(diǎn)數(shù)據(jù)分析

為定量分析傳統(tǒng)信號處理算法、Lag1&2算法、Lag1&3算法速度譜寬數(shù)據(jù)的差異，選取圖2所示PPI圖像中部分區(qū)域距離庫數(shù)據(jù)進(jìn)行分析。數(shù)據(jù)區(qū)間如下(共包含9 576個(gè)距離庫的回波數(shù)據(jù))：

區(qū)域一，方位角240°～270°，徑向距離60～90 km區(qū)間內(nèi)速度譜寬數(shù)據(jù)；

區(qū)域二，方位角330°～360°，徑向距離30～60 km區(qū)間內(nèi)速度譜寬數(shù)據(jù)。

圖3(a)、(b)和(c)分別顯示的是Lag1&2算法與傳統(tǒng)信號處理算法、Lag1&3算法與傳統(tǒng)信號處理算法、Lag1&3算法與Lag1&2算法速度譜寬偏差隨信噪比的分布圖像。分析上述圖像，可以得出以下結(jié)論：

1) 當(dāng)回波信噪比小于5 dB時(shí)，與Lag1&2算法、Lag1&3算法相比，傳統(tǒng)信號處理算法得到的速度譜寬數(shù)據(jù)普遍偏大，部分?jǐn)?shù)據(jù)偏差甚至超過-10 dB；

2) 當(dāng)回波信噪比大于10 dB時(shí)，Lag1&2算法、Lag1&3算法與傳統(tǒng)算法速度譜寬偏差均勻分布在-1.5～+1.5 m/s區(qū)間內(nèi)；

3) 當(dāng)回波信噪比大于10 dB時(shí)，Lag1&2算法與Lag1&3算法得到速度譜寬數(shù)據(jù)基本一致，其偏差大都不超過±0.5 m/s；

4) 當(dāng)回波信噪比小于5 dB時(shí)，Lag1&2算法與Lag1&3算法速度譜寬偏差分布范圍變大。

3.4 討論

(a)Lag1&2算法與傳統(tǒng)算法

(b)Lag1&3算法與傳統(tǒng)算法

當(dāng)回波信噪比大于10 dB時(shí)，回波信號功率要遠(yuǎn)大于噪聲功率，此時(shí)應(yīng)用傳統(tǒng)信號處理算法計(jì)算速度譜寬(即式(3))，噪聲功率對計(jì)算結(jié)果的影響可以忽略不計(jì)。當(dāng)回波信噪比小于5 dB時(shí)，回波信號功率接近甚至小于噪聲功率，此時(shí)噪聲功率對傳統(tǒng)信號算法速度譜寬數(shù)據(jù)影響較大，由于設(shè)置的平均噪聲功率偏低，導(dǎo)致傳統(tǒng)信號處理算法速度譜寬數(shù)據(jù)偏大，而且隨信噪比降低，噪聲功率的影響愈加嚴(yán)重。另外，天氣回波信號功率譜形狀可能受回波信噪比影響，當(dāng)回波信噪比小于5 dB時(shí)，回波信號功率譜與高斯函數(shù)偏差較大，導(dǎo)致Lag1&2算法與Lag1&3算法計(jì)算得到的速度譜寬數(shù)據(jù)存在一定誤差，因此出現(xiàn)圖3(c)中回波信噪比小于5 dB時(shí)譜寬偏差波動(dòng)較大的現(xiàn)象。

(c)Lag1&3與Lag1&2算法圖3速度譜寬偏差隨信噪比分布圖像(圖中速度譜寬偏差分別指σv2 -σv、σv3-σv、σv3-σv2)

4 結(jié)束語

對傳統(tǒng)信號處理算法計(jì)算速度譜寬式(3)進(jìn)行分析，發(fā)現(xiàn)平均噪聲功率是影響輸出譜寬數(shù)據(jù)質(zhì)量的一個(gè)不穩(wěn)定因素，為了消除噪聲功率的影響，就要避免使用零階延遲自相關(guān)參數(shù)。因此，本文提出采用多階延遲自相關(guān)參數(shù)進(jìn)行計(jì)算，如采用一階、二階、三階延遲自相關(guān)參數(shù)。由于該信號處理算法不涉及平均噪聲功率，因此其可以提高低信噪比區(qū)域的回波數(shù)據(jù)質(zhì)量。通過對傳統(tǒng)信號處理算法、Lag1&2算法、Lag1&3算法輸出的速度譜寬數(shù)據(jù)進(jìn)行比較，發(fā)現(xiàn)：當(dāng)信噪比小于5 dB時(shí)，傳統(tǒng)信號處理算法輸出的譜寬數(shù)據(jù)中摻雜有大量變化異常大的數(shù)值，而Lag1&2算法和Lag1&3算法輸出的譜寬數(shù)據(jù)與真實(shí)情況更接近；另外，信噪比小于5 dB時(shí)，Lag1&2算法和Lag1&3算法輸出的譜寬數(shù)據(jù)偏差波動(dòng)較大，作者認(rèn)為可能是由于低信噪比回波信號功率譜與高斯分布差異較大的原因。

天氣雷達(dá)實(shí)際信號處理過程中需要考慮各種因素的影響，如回波信噪比、譜寬等因素，不同信號處理算法在不同條件下處理性能也是有差異的，接下來需要系統(tǒng)論證 Lag1&2算法和Lag1&3算法在不同條件下的處理性能。另外，后面將開展研究應(yīng)用多階延遲自相關(guān)參數(shù)算法改善偏振回波數(shù)據(jù)質(zhì)量。

[1] 趙暢. 多普勒雷達(dá)及多源資料在局地短臨預(yù)報(bào)中的應(yīng)用[D]. 南京：南京信息工程大學(xué), 2014.

[2] FANG M, DOVIAK R J, MELNIKOV V. Spectrum Width Measured by WSR-88D:Error Sources and Statistics of Various Weather Phenomena[J]. Journal of Atmospheric and Oceanic Technology, 2004, 21(6):888-904.

[3] SIRMANS D, BUMGARNER B. Numerical Comparison of Five Mean Frequency Estimators[J]. Journal of Applied Meteorology, 1975, 14(6):991-1003.

[4] MONAKOV A, MONAKOV Y. Error Analysis of Pulse-Pair Estimates[J]. IEEE Trans on Aerospace and Electronic Systems, 2011, 47(3):2222-2230.

[5] SKOLNIK M I. 雷達(dá)手冊[M].3版. 南京電子技術(shù)研究所,譯. 北京：電子工業(yè)出版社, 2010:781-784.

[6] 盧小佳. 多普勒天氣雷達(dá)信號處理的研究[D]. 合肥：安徽大學(xué), 2014.

[7] IVIC I R, CURTIS C, TORRES S M. Radial-Based Noise Power Estimation for Weather Radars[J]. Journal of Atmospheric and Oceanic Technology, 2013, 30(12):2737-2757.

[8] IVIC I R. On the Use of a Radial-Based Noise Power Estimation Technique to Improve Estimates of the Correlation Coefficient on Dual-Polarization Weather Radars[J]. Journal of Atmospheric and Oceanic Technology, 2014, 31(9):1867-1880.

[9] JASSEN L H, VAN DER SPEK G A. The Shape of Doppler Spectra from Precipitation[J]. IEEE Trans on Aerospace and Electronic Systems, 1985, 21(2):208-219.

[10] TORRES S M. Estimation of Doppler and Polarimetric Variables for Weather Radars[D]. Norman, OK:University of Oklahoma, 2001.

[11] 郗蘊(yùn)天,陳春. 一種補(bǔ)償?shù)椭蓄l接收機(jī)I/Q不平衡的方法[J]. 雷達(dá)科學(xué)與技術(shù), 2015, 13(6):604-608.

XI Yuntian, CHEN Chun. A Compensation Method for I/Q Imbalance of Low-IF Receiver[J]. Radar Science and Technology, 2015, 13(6):604-608.(in Chinese)