趙筱彤, 周建江
(南京航空航天大學雷達成像與微波光子學教育部重點實驗室, 江蘇 南京 211106)
隨著信息技術的發展,單部雷達已經不能滿足復雜戰場環境的需求,多部雷達協同工作成為了必然的發展趨勢[1]。學者們結合了無線通信系統中多輸入多輸出(multiple input multiple output,MIMO)的概念以及綜合脈沖孔徑雷達 (synthetic impulse and aperture radar, SIAR)技術,在2004年IEEE雷達會議上正式提出了MIMO雷達[2]。隨著對MIMO雷達的深入研究,將MIMO雷達分為共址MIMO雷達和統計MIMO雷達[3]。本文選取單基地共址MIMO雷達進行研究,其在接收端采用匹配濾波方式,產生虛擬陣列孔徑,相當于擴大了陣元數目,使得空間分辨增益以及角度分辨力都得到了提高[4]。近年來利用MIMO雷達進行目標定位[5]受到廣泛的關注,現有的MIMO雷達多信號分類(multiple signal classification,MUSIC)算法大多是針對平穩信號進行研究[6]。但是在軍事領域中,機載雷達的探測對象為機動目標,其運動會產生多普勒效應,故接收信號為非平穩信號。
MUSIC算法自1979年的Schmidt創立開始,一直在創新和改進[7-8],時至今日,仍有很多學者致力于MUSIC算法的研究[9]。針對非平穩信號的MUSIC算法,一般采用時頻分析[10]。現有的MIMO雷達時頻MUSIC算法僅僅是利用空間時頻分布(spatial time frequency distribution, STFD)矩陣代替接收信號的協方差矩陣,完成對非平穩信號的處理[11-12]。本文在此基礎上對非均勻線陣單基地MIMO雷達展開降維預處理,降低算法運算量,同時參考盲源分離技術[13],通過白化處理、時頻分析、時頻點篩選、正交聯合對角化以及Givens旋轉等信號處理,對MIMO雷達時頻MUSIC算法進行改進,提高算法的估計精度,降低算法的均方根誤差(root mean squared error, RMSE),使算法具有低信噪比(signal to noise ratio,SNR)、低快拍數的優勢。
在電子對抗戰中,隨著電子偵察技術的快速發展,必須要提高雷達的生存能力。雷達在保證精確探測到遠方弱小目標的同時,還要防止被偵察設備、反輻射導彈截獲[14]。因此雷達工作的SNR環境通常較低,而現有的MIMO雷達時頻MUSIC算法所適用的最低SNR多為0 dB左右[15],而仿真實驗證明本文提出的改進MUSIC算法將SNR降低至-2 dB,驗證了MIMO雷達使用本文改進MUSIC算法進行角度估計時,具有低截獲性能。
MIMO雷達自身具有寬發窄收的特點,具有較好的反偵察、低截獲能力[16]。針對測向MIMO雷達,發射天線陣元間距增大有助于支持目標方面的空間分集,而分集會導致目標分量在接收信號中更有利的分布[17]。同時,綜合考慮MIMO雷達的發射和接收陣,使其可以在小規模陣列條件下實現更高的角度分辨率。本文參考文獻[18]中分布式發射和密布式接收的MIMO陣列,此陣列可以有效抑制目標閃爍,同時實現波達方向估計。
考慮到單基地MIMO雷達系統的特點[19]以及本文研究的目的,采用圖1所示單基地MIMO雷達系統模型。假設發射天線為M個陣元組成的非均勻線陣,接收天線為N個陣元組成的均勻線陣,線陣布陣方式如圖1所示。考慮到資源利用率問題,發射陣列中前(k+1)個陣元排布與接收陣列完全一致,為陣元間距d的密布式陣列,故接收陣列陣元數N=k+1;而后(k-1)個陣元為了達到抑制目標閃爍,提高角度估計精度的目的,本文將發射天線間距擴大為(k+1)d,故整個發射陣列陣元數M=2k。圖1中陣元間距d=λ/2,λ為信號波長。
假設圖1中的探測目標共有Q個,第q個目標的探測角度為θq,且目標相對于雷達滿足遠場條件。因此接收信號的矩陣形式[20]如下所示:
(1)
式中:取快拍數為L,則t=1,2,…,L。S(t)=[S1(t),S2(t),…,SM(t)]為天線發射信號矢量,根據MIMO雷達的特點,SmSmH=I,w(t)為均值為0的高斯白噪聲。
機動目標運動速度的變化會產生多普勒頻移,所以目標的回波信號為非平穩信號。故bq(t)=βqexp{j2π(fqt+(μq/2)t2)} 為第q個目標的散射信號,其中βq為第q個目標的雷達散射截面積(radar cross section, RCS)的復振幅,fq和μq為第q個目標的多普勒頻移和調頻斜率。將式(1)中的發射天線導向矢量at(θq)和接收天線的導向矢量ar(θq)分別記為
at(θq)=[at1(θq);at2(θq)]=
[1,e-j2π(d/λ)sin θq,…,e-j2π(k(d/λ)sin θq,
e-j2πk(k+k+1)(d/λ)sinθq,…,e-j2π(M-1)(d/λ)sin θq]T
(2)
ar(θq)=[1,e-j2π(d/λ)sin θq,…,e-j2π(N-1)(d/λ)sin θq]T
(3)
MIMO雷達在信號處理的過程中,在接收陣列處設置M個匹配濾波器,對接收信號進行一維時域匹配濾波,形成MN個虛擬陣列[21]。濾波后的信號y(t)如下所示:
y(t)=A(θ)b(t)+W(t)
(4)
L個快拍下的y(t)形成矩陣Y,如下所示:
Y=AB+W
(5)
式中:Y=[y(1),y(2),…,y(L)]為MN×L維信號矩陣;A(θ)=[a1(θ1),a2(θ2),…,aq(θq)]為MN×Q維的虛擬陣列聯合導向矩陣,其中aq(θq)=at(θq)?ar(θq)是由發射天線導向矢量at(θq)與接收天線導向矢量ar(θq)進行Kronecker積運算得到的;矩陣W為MN×L維、均值為0、方差為δ2的高斯白噪聲矩陣;B=[b(1),b(2),…,b(L)]為Q×L維散射系數,其中b(t)的表達式如下所示:
(6)
通過文獻[22]的研究可知,當發射線陣與接收線陣均為陣元間距是半波長的均勻陣列時,有效的虛擬陣元數目為發射陣元數與接收陣元數的和減1。由圖1可知,本文采用的單基地MIMO雷達非均勻發射線陣的前(k+1)個陣元是間距為d的均勻線陣,而后(k-1)個陣元可視為每間隔(k+1)個間距為d的陣元取一個作為發射陣元。結合空間卷積原理[23],該模型的有效虛擬陣元數目中發射陣元數可以理解為是發射陣列中工作和未工作的,即范圍內所有陣元數目為(2k+1)+(k+1)(k-2)+1。由于接收線陣為均勻陣列,因此該模型有效虛擬陣元數是
Ne=[(2k+1)+(k+1)(k-2)+1]+N-1=
k2+k+N-1
(7)
這種方法填補了非均勻發射線陣中陣元間距大于半波長的部分。根據有效虛擬陣元數目Ne,可將MN×1維的導向矢量a(θ)通過線性變換為g(θ)=[1,ejπ sin θ,…,ejπ(Ne-1)sin θ]T即Ne×1維的導向矢量,線性變換過程如下所示:
a(θ)=at(θ)?ar(θ)=Fg(θ)
(8)
由于非均勻發射陣列具有前半部分和后半部分均為均勻線陣的特點,因此針對這一特殊性,將變換矩陣F調整為
F=[F1;F2]
(9)
式中:F1是針對發射線陣中前(k+1)段均勻陣列的變換矩陣;F2是針對后(k-1)段均勻陣列的變換矩陣,表達式為
(10)
(11)
因此,F1、F2組成了MN×Ne維的變換矩陣F。根據式(8)將導向矩陣A變為
A=FG
(12)
其中,虛擬線陣的導向矢量矩陣G=[g(θ1),g(θ2),…,g(θq)],將式(12)代代入式(1)中,得到
X=FGBS+W
(13)
由式(13)可以看出,信號位于由G構成的低維空間之中,因此,可以將信號轉換為低維信號。由于目標信息保留在轉換后的低維信號中,還可以在低維空間中進行角度估計,降低了算法的運算量。假設一個Ne×MN維的降維變換矩陣D,則式(13)轉變為
XD=DFGBS+DW
(14)
為了保證降維處理后的噪聲仍是均值為0、方差為δ2的高斯白噪聲,變換矩陣就必須滿足DDH=INe,因此,變換矩陣D為
(15)
對匹配濾波處理之后的接收信號Y,也進行相應的降維處理,得到低維空間下的匹配信號Z,其表達式如下:
Z=DY
(16)
通過以上所述,將MN×L維矩陣Y變換為Ne×L維矩陣Z。降低了信號的維度空間,同時保留了目標探測信息。這種針對非均勻陣列的降維操作不僅僅大大降低了運算量,還不影響后續的MUSIC算法。
正常情況下,雷達獲取數據都具有一定的相關性[24],本文采用白化處理方式降低信號之間的相關性,提高信號的收斂性[25],為后續時頻分析中時頻點的篩選以及低SNR下的MUSIC算法提供了更高的精度。
降維處理后接收信號Z(t)的自相關函數為
(17)
如果沒有噪聲,處理后信號Z(t)的自相關函數為
RD=(DAB)(DAB)H
(18)
因此,由式(17)和式(18)可得
(DAB)(DAB)H=RZ-δ2I
(19)
假設對應的白化矩陣為W,定義為WRDWH=W(DAB)·(DAB)HWH=I,結合酉矩陣定義UUH=I,即可得
DAB=W#U
(20)
上式中W#為W的摩爾-彭若斯廣義逆矩陣,將式(20)代入式(19)得

(21)
令K=diag[(λ1-δ2)1/2,(λ2-δ2)1/2,…,(λNe-δ2)1/2],那么W#=UK,進而可以得到白化矩陣W的表達式為
W=K-1UH
(22)
因此,可以得到白化處理后的信號H為
H=WZ
(23)
WVD是一種二次型時頻分布的方式[26]。信號矢量H(t)基于WVD分布的STFD矩陣[27]表達式如下:

(24)
式中的STFD矩陣是由信號自身的時頻分布矩陣和信號之間的時頻分布矩陣組成的,分別稱為自項和交叉項。由于信號矩陣H(t)中包含Ne個信號信息,因此會產生交叉項。假設Ne=2,即H=[H1;H2],那么其相應的STFD矩陣為
(25)
其中,主對角線元素對應的是非平穩信號自身時頻點的WVD分布結果,而其他部分元素對應的則是非平穩信號間形成時頻點的WVD分布結果。

(26)
由于自項與矩陣范數之間的比率在0~1之間,因此設置閾值如下:
(27)
通過現有知識[28]大量仿真得到,當ζ>0.1時,可判定該時頻點為自項點。
時頻分析中篩選出來的n個時頻點所對應的STFD矩陣按列串聯為一組復矩陣Rjd=[DH1;DH2;…;DHn],本文采用基于Givens旋轉的正交聯合對角化[29]的方法,來計算復矩陣的近似實值公共正交基。由聯合對角化原理可知,聯合矩陣具有如下性質:
Rjd=VCVH
(28)
式中:C=[C1;C2;…;Cn]是擬對角矩陣級聯;V是正交矩陣。因此式(28)可以拆分為Rjd1=VC1VH,…,Rjdn=VCnVH。該算法的目的是通過迭代運算找到一個正交矩陣V使得矩陣C盡可能是對角陣,且與矩陣Rjd共有一種“平均特征結構”。式(28)等價于求解約束最大化問題:
(29)
本文對每一個STFD矩陣進行Givens旋轉,令c2+s2=1,其中c=cos(θ),s=sin(θ),θ為旋轉角度,旋轉矩陣如下:
(30)
由式(30)可知,通過對n個STFD矩陣的一系列Givens旋轉迭代運算,即可實現對該矩陣組的正交聯合對角化,所有Givens旋轉矩陣的乘積即是聯合對角化矩陣V,迭代收斂時最后一步的C,即為最優對角陣組。
由于V中含有信號子空間和噪聲子空間全部的信息,因此,直接利用V和C進行空間譜估計。通過對角陣列C中的每一組對角矩陣得到矩陣的聯合特征值,而后找到對應的特征向量UN,最終改進MUSIC算法的空間譜函數為
(31)
綜合本節內容,本文所提出的MIMO雷達改進MUSIC算法信號處理流程,如圖2所示:
機載雷達在電子戰中通常采用射頻隱身技術來對雷達進行保護。因此,提高MIMO雷達的射頻隱身能力,對抗截獲接收機是本文研究MIMO雷達改進MUSIC算法的主要目的。
由雷達距離方程可以得到,MIMO雷達最大截獲距離為
(32)
式中:GT和GR為發射和接收天線增益;M是陣列數目;PT為雷達發射功率;η為信號占空比;tB為信號駐留時間;δ為目標RCS;N0為雷達接收機處的噪聲功率譜密度;SNRmin為雷達間最小可探測SNR。由于第3節對MIMO雷達接收信號進行了相關的信號處理,故式(32)中陣列數目M為第2.1節中有效虛擬陣元數目Ne。
射頻隱身技術中用施里海爾截獲因子衡量雷達的射頻隱身性能,施里海爾截獲因子定義為截獲接收機最大截獲距離與雷達最大探測距離之比。結合式(32)通過推導,可以得到下式所示的MIMO雷達截獲因子表達式:
(33)
式中:Gi是截獲接收機天線增益;γ∈(0,1)為非相參積累損失,本文中令Ni與N0相等;Bi為截獲接收機帶寬;SNRr,min為雷達最小可檢測信噪比;SNRi,min為截獲接收機最小可檢測SNR,∑(θ)為雷達在截獲接收機方向歸一化方向圖因子。
由文獻[30]可知,雷達的功率域被截獲概率為
p=1-e-0.5α
(34)
由式 (33)可知,MIMO雷達的射頻隱身能力受到信號發射功率PT和駐留時間tB的影響。首先在環境參數確定的情況下,SNR與發射功率PT成正比,即SNR越小,PT就越小,相應的α就越小,p就越小,雷達射頻隱身能力就越好;其次由文獻[31]可知,信號駐留時間tB與快拍數L成正比,即L越小,tB越小,相應的α就越小,p就越小,雷達射頻隱身能力就越好。
仿真實驗模擬在某一時刻單基地MIMO雷達探測運動目標方位角度的場景。雷達模型和圖1所示,模型是發射陣列為分布式線陣、接收陣列為密布式線陣的單基地MIMO雷達,令k=4,則M=8,N=5,式(7)中的有效虛擬陣列Ne=24,陣元間距d為半波長。假設有3個運動目標,在某一時刻散射信號起止頻率分別為fs1=0.12、fe1=0.32、fs2=0.15、fe2=0.35、fs3=0.18、fe3=0.38。
本文共設計3組實驗,通過第1節中介紹的MIMO雷達時頻MUSIC算法與本文改進MUSIC算法對比實驗的方式,分析兩種算法的探測性能、臨近角度分辨能力以及射頻隱身能力。
假設3個運動目標回波信號到達角分別為θ1=-10°、θ2=3°、θ3=10°,快拍數L為100,SNR為10 dB。MIMO雷達時頻MUSIC算法與本文改進MUSIC算法對比仿真實驗結果如圖3所示。圖3中,MIMO雷達使用兩種算法均能探測出3個目標角度。圖3(a)中時頻MUSIC算法波峰的高度為30 dB;圖3(b)中改進MUSIC算法波峰的高度為40 dB。顯然,本文改進MUSIC算法的波束更為尖銳,指向精度有所提高。
假設3個運動目標回波信號到達角分別為θ1=-10°,θ2=9°,θ3=10°,快拍數L為100,SNR為5 dB。MIMO雷達時頻MUSIC算法與本文改進MUSIC算法對比仿真實驗結果如圖4所示。
仿真目標2和目標3的角度相差僅為1°,圖4(a)中目標2、目標3探測角度重合,探測角度為9.5°,無法分辨;圖4(b)中探測角度為9°、10°,可以分辨出臨近目標,證明了本文提出的改進MUSIC算法具有較好的分辨能力。
通過第4節的分析可知,SNR和快拍數越小,射頻隱身性能就越好。因此實驗3仿真兩種MUSIC算法的估計精度隨著SNR和快拍數的變化過程,并且證明了本文算法的低截獲性能。
本文采用RMSE來衡量兩種MUSIC算法的估計精度。利用蒙特卡羅實驗使仿真過程趨近于真實過程。基于蒙特卡羅實驗的RMSE表達式如下所示:
(35)
式中,Q為目標數目;P為蒙特卡羅次數,本次實驗選取蒙特卡羅數為P=200。
假設3個運動目標,其信號到達角分別為θ1=-10°、θ2=3°、θ3=10°、SNR=5 dB。兩種MUSIC算法RMSE隨快拍數L變化的對比仿真實驗結果如圖5所示。
由圖5可見,當快拍數低于80時,本文提出的MIMO雷達改進MUSIC算法的RMSE低于時頻MUSIC算法,證明改進MUSIC算法在低快拍數的條件下具有較好的估計精度,同時也證明了在相同的角度估計精度下,改進MUSIC算法所需的快拍數更小,具有更好的射頻隱身能力。
假設3個運動目標的信號到達角不變,快拍數L為100。兩種MUSIC算法的RMSE隨著SNR變化的對比仿真實驗結果如圖6所示,其中圖6(a)為MIMO雷達時頻MUSIC算法的RMSE,圖6(b)為MIMO雷達本文改進MUSIC算法的RMSE。如圖6所示,在SNR=-15 dB的狀態下,圖6(b)的估計誤差比圖6(a)小20°,證明在具有相同射頻隱身能力時,本文改進MUSIC算法具有更好的估計精度。除此之外,當2種MUSIC算法估計精度達到最高即RMSE為0時,圖6(a)中的最小SNR為0 dB,圖6(b)中的最小SNR為-2 dB,證明了本文改進MUSIC算法所需的最小SNR降低了2 dB,具有更好的射頻隱身能力。
功率域低截獲概率技術主要是降低雷達發射機的峰值功率。設雷達與截獲接收機具有相同的噪聲功率譜密度,噪聲溫度為290 K,噪聲系數為2,信號占空比為0.01,目標RCS為1 m2。參考文獻[30]中雷達和截獲接收機的參數。根據本文陣列排布特點,陣元數M=Ne=24,令駐留時間為26 ms,根據式(33)和式(34)得到如圖7所示的雷達被截獲概率隨SNR變化曲線圖。
圖6仿真實驗中兩種算法估計精度最高時,所需最低SNR分別為0 dB和-2 dB,代入圖7可知時頻MUSIC算法的被截獲概率為1.74%,而本文改進算法的被截獲概率為1.55%。圖6(b)中改進算法在-5 dB時,才開始有較大的估計誤差,此時的被截獲概率為1.31%。因此本文改進算法保證了雷達低截獲性能。
本文提出一種低截獲的MIMO雷達改進MUSIC算法。本文在解決單基地MIMO雷達非均勻陣列探測運動目標問題時,采用降維預處理、白化處理、時頻分析、時頻點篩選、正交聯合對角化幾個信號處理過程,提高了算法的估計精度以及射頻隱身能力。與MIMO雷達時頻MUSIC算法的對比仿真實驗結果表明:① 本文提出的MIMO雷達改進MUSIC算法指向精度有所提高;② 本文改進MUSIC算法可以分辨角度差為1°的臨近目標,分辨能力有所提高;③ 本文改進MUSIC算法在低SNR、低快拍數條件下估計精度更高,所適用的SNR環境降低了2 dB。本文算法的研究僅考慮了MIMO雷達對運動目標的定位問題,在后續研究中可以將算法應用在雷達對運動目標角度跟蹤之中。