依據(jù)信噪比的雷達(dá)信號(hào)質(zhì)量估計(jì)

(空軍航空大學(xué) 航空作戰(zhàn)勤務(wù)學(xué)院,長(zhǎng)春 130022)

1 引 言

在戰(zhàn)場(chǎng)電磁環(huán)境中，微弱信號(hào)由于其功率低，甚至達(dá)不到接收機(jī)靈敏度，尤其是低截獲概率(Low Probability of Intercept,LPI)雷達(dá)采用功率管理技術(shù)，致使現(xiàn)有技術(shù)無(wú)法對(duì)此類(lèi)信號(hào)實(shí)現(xiàn)有效處理。信號(hào)增強(qiáng)技術(shù)是專(zhuān)門(mén)針對(duì)微弱信號(hào)偵察截獲難而進(jìn)行的研究，其中如何檢驗(yàn)增強(qiáng)技術(shù)是否有效是信號(hào)增強(qiáng)的重要部分，最直接的方法是對(duì)比增強(qiáng)前后信號(hào)質(zhì)量。在具體的驗(yàn)證中，本文通過(guò)對(duì)比信號(hào)的信噪比來(lái)衡量信號(hào)的質(zhì)量。

信噪比(Signal-to-Noise Ratio,SNR)一直是現(xiàn)代通信系統(tǒng)中一個(gè)非常重要的參數(shù)，在傳統(tǒng)電子戰(zhàn)領(lǐng)域研究得不多。許多學(xué)者在通信領(lǐng)域就信噪比估計(jì)這一問(wèn)題進(jìn)行了廣泛研究。文獻(xiàn)[1]研究了針對(duì)線(xiàn)性調(diào)頻信號(hào)在時(shí)間域上利用傅里葉級(jí)數(shù)擬合自相關(guān)主瓣達(dá)到信噪比估計(jì)的算法，主要分析了時(shí)寬和帶寬對(duì)結(jié)果的影響。文獻(xiàn)[2]為在載波頻率精確恢復(fù)前提高多進(jìn)制數(shù)字相位調(diào)制信號(hào)在低信噪比下的估計(jì)精度，提出了一種在符號(hào)定時(shí)恢復(fù)和幀同步后提取同步段符號(hào)，相關(guān)運(yùn)算后在頻域進(jìn)行信噪比估計(jì)的算法。文獻(xiàn)[3]研究了利用最大似然估計(jì)法解決在調(diào)頻差分混沌鍵控系統(tǒng)中的多徑瑞利衰落信道的信噪比估計(jì)問(wèn)題。文獻(xiàn)[4]將Gram-Schmidt正交化過(guò)程引入到緊縮投影近似子空間跟蹤(Projection Approximation Subspace Tracking deflation,PASTd)中，使計(jì)算得到的特征向量相互正交，從而保證算法具有更好的收斂性能。文獻(xiàn)[5]研究了在高動(dòng)態(tài)環(huán)境下基于數(shù)據(jù)輔助信噪比估計(jì)算法，其實(shí)質(zhì)是利用似然函數(shù)求微分的思想。

本文根據(jù)雷達(dá)信號(hào)復(fù)雜多變的特點(diǎn)，采用基于子空間分解的信噪比估計(jì)算法：首先針對(duì)傳統(tǒng)算法對(duì)自相關(guān)矩陣估計(jì)不準(zhǔn)的問(wèn)題進(jìn)行改進(jìn)，隨后利用粒子群優(yōu)化(Partical Swarm Optimization,PSO)算法[6]對(duì)自相關(guān)矩陣進(jìn)行特征值分解，然后利用最小描述長(zhǎng)度(Minimum Description Length,MDL)算法[7]確定出信號(hào)的維數(shù)，將信號(hào)分為信號(hào)子空間和噪聲子空間，進(jìn)而實(shí)現(xiàn)信噪比估計(jì)。本文在加性高斯白噪聲(Additive White Gaussian Noise，AWGN)信道條件下進(jìn)行理論推導(dǎo)和仿真討論，從標(biāo)準(zhǔn)偏差和方差的角度估計(jì)信噪比估計(jì)算法的性能。

2 信號(hào)模型

任何信號(hào)都可以表示為各種諧波之和，因此，任何信號(hào)估計(jì)問(wèn)題都可以等效為諧波估計(jì)問(wèn)題。在本文中假設(shè)信號(hào)模型如下：

(1)

3 信噪比估計(jì)計(jì)算

信號(hào)x(n)可以變形為

(2)

式中：αi=Aiejφi。

取時(shí)間長(zhǎng)度為M的時(shí)間窗，窗內(nèi)的信號(hào){x(n),n=n,n+1,…,n+M-1}可以表示為矢量的形式：

X(n)=[x(n),x(n+1),…,x(n+M-1)]T。

因此，

(3)

式中：

w(n)=[w(n),w(n+1),…,w(n+M-1)]T，
V(fi)=[1,ej2πfi,…,ej2π(M-1)fi]T。

X(n)的自相關(guān)矩陣為

(4)

式中：

V=[V(f1),V(f2),…,V(fP)]M×P，

自相關(guān)矩陣RX可以表示為

(5)

式中：λm是以降階排列的RX的特征值，λ1≥λ2≥…≥λM；qm是相應(yīng)的特征向量，而

(6)

式中：

Qs=[q1,q2,…,qP]M×P,

Qw=[qP+1,qP+2,…,qM]M×(M-P)。

在整個(gè)分析問(wèn)題過(guò)程中，關(guān)鍵是要能夠確定Λ的秩P，就能通過(guò)自相關(guān)矩陣RX的特征值分解實(shí)現(xiàn)區(qū)分開(kāi)信號(hào)子空間和噪聲子空間，由此信噪比估計(jì)值可以由以下公式計(jì)算：

(7)

3.1 自相關(guān)矩陣估計(jì)改進(jìn)

在進(jìn)行信噪比估計(jì)的過(guò)程中，如果采用式(4)計(jì)算出來(lái)的結(jié)果作為信號(hào)自相關(guān)矩陣，就會(huì)存在一定的偏差。

真實(shí)情況下X(n)的自相關(guān)矩陣為

RX=E[X(n)XH(n)]=

E[(S(n)+w(n))(S(n)+w(n))H]=

Rss+Rww+Rsw+Rws。

(8)

比較式(8)和式(4)，式(8)的前兩部分是式(4)中的兩項(xiàng)和，代表了信號(hào)和噪聲部分；后兩部分代表了信號(hào)和噪聲的互相關(guān)。在研究的模型當(dāng)中，假定的噪聲是零均值并且和信號(hào)是不相關(guān)的，因此式(4)中后兩部分應(yīng)該等于零。然而，在數(shù)據(jù)量有限的情況下，式(8)中后兩項(xiàng)就會(huì)變得重要，它直接導(dǎo)致了自相關(guān)矩陣估計(jì)偏離真實(shí)值。

本文就針對(duì)上述問(wèn)題在估計(jì)自相關(guān)矩陣時(shí)提出了改進(jìn)形式。改進(jìn)的信號(hào)矩陣形式為

(9)

式中：n=0,1,…,N-1，N是信號(hào)的長(zhǎng)度，又叫快拍數(shù)。

這樣改進(jìn)的目的是充分考慮到信號(hào)和噪聲之間的相關(guān)性，盡可能按照式(8)來(lái)構(gòu)造信號(hào)矩陣，相對(duì)于式(4)提高了估計(jì)精度。改進(jìn)后的信號(hào)自相關(guān)矩陣為

(10)

3.2 特征值計(jì)算

信號(hào)自相關(guān)矩陣的特征值計(jì)算是本文的一個(gè)重要組成部分，準(zhǔn)確快速地對(duì)自相關(guān)矩陣RX進(jìn)行特征值分解是基于子空間理論進(jìn)行信噪比估計(jì)的關(guān)鍵。普通的求解方法有迭代法和變換法，其中，迭代法是利用矩陣向量之間相乘進(jìn)而求得特征值和特征向量，常用的方法有Lanczos法[8]、Davidson法[9]等；變換法則是直接對(duì)矩陣進(jìn)行變換，使之變成容易求解特征值和特征向量的新矩陣。這些方法都是矩陣論中的經(jīng)典計(jì)算方法，但是普遍存在著存儲(chǔ)量大、精度低、收斂速度較慢等不足。另外可以采用一些分解矩陣的方法，例如Jacobi旋轉(zhuǎn)[10]和QR分解[11]，但其都存在計(jì)算量大、收斂速度慢的缺點(diǎn)。粒子群優(yōu)化算法是基于群體智能的進(jìn)化優(yōu)化算法，采用實(shí)數(shù)求解，需要調(diào)整的參數(shù)較少，且算法簡(jiǎn)單,易于實(shí)現(xiàn)。因此，本文選其作為RX矩陣特征值求解的方法。

3.2.1PSO算法

PSO初始化開(kāi)始時(shí)，會(huì)隨機(jī)產(chǎn)生一群粒子又叫做隨機(jī)解，在這個(gè)過(guò)程中粒子不斷跟蹤兩個(gè)極值，一個(gè)是個(gè)體極值，指的是粒子自身找到的最優(yōu)位置；另一個(gè)是全局極值，指的是迄今整個(gè)粒子群找到的最優(yōu)位置來(lái)更新速度和位置。

假設(shè)第i個(gè)粒子的位置是Xi=(xi1,xi2,…,xid)，速度是Vi=(vi1,vi2,…,vid)，其中d是粒子的維數(shù)；第i個(gè)粒子的歷史最優(yōu)位置為Pi=(pi1,pi2,…,pid)，其適應(yīng)度值記為fitness(Pi)；迄今為止粒子群搜索到的最好位置記為Pg=(pg1,pg2,…,pgd)，其適應(yīng)度記為fitness(Pg)。粒子按照下式調(diào)整自己的位置：

(11)

(12)

式中：1≤i≤N，N為粒子群規(guī)模；1≤d≤D，D為粒子維數(shù)；k為迭代次數(shù)(k≥0)；加速常數(shù)c1和c2是非負(fù)數(shù)；r1和r2是(0,1)區(qū)間的隨機(jī)數(shù)；慣性權(quán)重w一般取0.5～0.9之間的常數(shù)。在解空間內(nèi)粒子不斷跟蹤個(gè)體極值與全局極值進(jìn)行搜索，直到達(dá)到規(guī)定的迭代次數(shù)或者滿(mǎn)足規(guī)定的誤差標(biāo)準(zhǔn)為止。

3.2.2基于PSO的特征值求解

個(gè)體由粒子位置X和粒子速度V這兩部分組成，分別包含兩個(gè)分量，即(X,V)=((x1,x2),(v1,v2))，其中，(x1,x2)對(duì)應(yīng)矩陣特征值的實(shí)部和虛部，即特征值λ=x1+x2i；(v1,v2)表示特征值對(duì)應(yīng)的速度的實(shí)部和虛部。

把個(gè)體代入特征方程P(λ)=|A-λI|=0中，得P(λ)=A-(x1+x2i)I=0，設(shè)e=det(A-(x1+x2i)I)，對(duì)λ的求解轉(zhuǎn)化為e=det(A-(x1+x2i)I)的最小化問(wèn)題。即PSO的適應(yīng)度函數(shù)為

minf=|det(A-(x1+x2i)I)| 。

(13)

基于PSO的矩陣特征值求解步驟如圖1所示。

圖1 PSO算法求解矩陣特征值流程 Fig.1 PSO algorithm process for solving matrix eigenvalue

3.2.3PSO的改進(jìn)

為了防止PSO算法進(jìn)入局部最優(yōu)值“陷阱”，縮短算法收斂時(shí)間，本文對(duì)PSO算法進(jìn)行如下改進(jìn)：

(1)如果粒子的適應(yīng)度值不等于fitness(Pg)，利用式(11)、(12)更新其位置；

(2)如果粒子的適應(yīng)度值等于fitness(Pg)，利用如下公式更新其位置：

4 信號(hào)子空間維數(shù)的確定

含有噪聲的信號(hào)的維數(shù)可以由最終預(yù)測(cè)誤差(Final Prediction Error,FPE)[12]、赤池信息量(Akaike Information Criterion,AIC)準(zhǔn)則、自回歸傳遞函數(shù)或者最小描述長(zhǎng)度來(lái)計(jì)算，但是在接收雷達(dá)信號(hào)較少的情況下，最小描述長(zhǎng)度效果最好[13]。所以，本文選擇MDL原理來(lái)計(jì)算信號(hào)子空間的維數(shù):

(14)

則信號(hào)子空間的維數(shù)為

(15)

綜上理論分析，得出信噪比估計(jì)的過(guò)程如下：

(1)建立雷達(dá)信號(hào)的自相關(guān)矩陣RX；

(2)利用PSO算法對(duì)RX進(jìn)行特征值求解，求出其特征值序列；

(3)依據(jù)MDL原理表達(dá)式確定RX的維數(shù)P；

(4)按照式(7)求其信噪比估計(jì)值。

5 仿真與性能分析

為了測(cè)試算法的性能，我們選取雷達(dá)信號(hào)中最常見(jiàn)的常規(guī)雷達(dá)信號(hào)、線(xiàn)性調(diào)頻信號(hào)(Linear Frequency Modulation,LFM)和二相編碼信號(hào)(Binary Phase Shift Keying,BPSK)并進(jìn)行100次蒙特卡洛(Monte-Carlo)仿真實(shí)驗(yàn)，主要分析觀測(cè)長(zhǎng)度M和不同信號(hào)對(duì)實(shí)驗(yàn)結(jié)果的影響。在仿真實(shí)驗(yàn)中，采用信噪比的估計(jì)均值和估計(jì)均方根誤差(Root Mean Square Error,RMSE)來(lái)衡量算法的性能。RMSE定義如下：

5.1 仿真1

設(shè)置情景:接收機(jī)收到一個(gè)常規(guī)雷達(dá)信號(hào)，信號(hào)形式為單一載頻。信號(hào)長(zhǎng)度為256點(diǎn)，時(shí)間窗長(zhǎng)度為16，初始種群個(gè)數(shù)為50，慣性權(quán)重為0.8，空間維數(shù)為2，兩個(gè)加速常數(shù)都為0.5。在本實(shí)驗(yàn)中，設(shè)置信號(hào)的歸一化頻率為0.03。仿真結(jié)果如圖2所示。

圖2 一部常規(guī)雷達(dá)信號(hào)SNR估計(jì)Fig.2 SNR estimation of a conventional radar signal

5.2 仿真2

設(shè)置情景:接收機(jī)收到兩部常規(guī)雷達(dá)信號(hào)，信號(hào)形式均為單一載頻。在本實(shí)驗(yàn)中，設(shè)置兩個(gè)信號(hào)的歸一化頻率分別為0.03和0.07。仿真結(jié)果如圖3所示。

(a)SNR估計(jì)均值

(b)SNR估計(jì)均方根誤差圖3 兩部常規(guī)雷達(dá)信號(hào)SNR估計(jì)Fig.3 SNR estimation of two conventional radar signals

5.3 仿真3

設(shè)置情景:接雷收機(jī)收到兩部達(dá)信號(hào)，其中一部為常規(guī)雷達(dá)信號(hào)，信號(hào)形式為單一載頻，歸一化頻率值為0.03；另外一部為低截獲概率雷達(dá)，采用線(xiàn)性調(diào)頻的調(diào)制方式，歸一化頻率變化范圍0.01～0.02。仿真結(jié)果如圖4所示。

(a)SNR估計(jì)均值

(b)SNR估計(jì)均方根誤差圖4 常規(guī)和線(xiàn)性調(diào)頻樣式混合信號(hào)SNR估計(jì)Fig.4 The conventional and LFM mixed signal SNR estimation

5.4 仿真4

設(shè)置情景:接收機(jī)收到3部雷達(dá)信號(hào)，其中一部為常規(guī)雷達(dá)信號(hào)，信號(hào)形式為單一載頻，歸一化頻率值為0.03；另外一部為低截獲概率雷達(dá)，信號(hào)采用線(xiàn)性調(diào)頻的調(diào)制樣式，歸一化頻率變化范圍0.01～0.02；第三部也為低截獲概率雷達(dá)，信號(hào)采用二相編碼的調(diào)制樣式，歸一化頻率為0.025。仿真結(jié)果如圖5所示。

(b)SNR估計(jì)均方根誤差

5.5 仿真5

將本文的算法和EB方法、GT方法和M2M4方法[14-16]進(jìn)行性能對(duì)比，設(shè)置信號(hào)和仿真1的相同。仿真結(jié)果如圖6所示。

(a)不同方法的SNR估計(jì)均值

(b)SNR估計(jì)均方根誤差圖6 4種方法的SNR估計(jì)Fig.6 SNR estimation comparison among four methods

仿真1～4分別設(shè)置了不同的信號(hào)條件，從單一常規(guī)信號(hào)一直增加到多個(gè)復(fù)雜調(diào)制信號(hào)混疊輸入，從結(jié)果來(lái)看，經(jīng)過(guò)100次的蒙特卡洛實(shí)驗(yàn)，在-10～10 dB之間信噪比估計(jì)值與理論值基本相等，隨著SNR逐漸增大，估計(jì)的均方根誤差呈現(xiàn)指數(shù)衰減，屬于正常誤差范圍之內(nèi)。仿真5將本文算法同文獻(xiàn)[14-16]的信噪比估計(jì)算法進(jìn)行了性能對(duì)比，從仿真結(jié)果來(lái)看，當(dāng)SNR<-3 dB時(shí)，本文算法的估計(jì)值精度優(yōu)于其他算法，均方根誤差明顯小于其他算法。此外，在求解自相關(guān)矩陣的特征值時(shí)，本文運(yùn)算復(fù)雜度為O(n2)，而基于QR分解的算法復(fù)雜度為O(n3)。

6 結(jié) 論

本文研究了在復(fù)雜電磁環(huán)境中對(duì)接收雷達(dá)信號(hào)的信噪比估計(jì)，解決了在有限信號(hào)量下自相關(guān)矩陣估計(jì)偏差問(wèn)題，同時(shí)也解決了傳統(tǒng)算法求解特征值復(fù)雜度高、運(yùn)算速度慢的問(wèn)題。仿真結(jié)果表明，本文算法明顯優(yōu)于EB方法、GT方法和M2M4方法，有效提高了估計(jì)精度。此算法是一種進(jìn)行盲信噪比估計(jì)的有效方法，實(shí)現(xiàn)了將信噪比作為衡量信號(hào)質(zhì)量的指標(biāo)，為信號(hào)增強(qiáng)研究領(lǐng)域提供了準(zhǔn)確有效的檢驗(yàn)方法。在后續(xù)研究中，我們將考慮在更低信噪比的極端環(huán)境下的算法改進(jìn)。