基于協(xié)方差擬合的MIMO雷達高精度DOA估計算法

竇慧晶, 謝金鑫, 孫 璐，楊 帆

(北京工業(yè)大學信息學部，北京 100124)

多輸入多輸出(multiple input and multiple out，MIMO)雷達系統(tǒng)是近年來提出的一種新體制雷達系統(tǒng)，MIMO技術給雷達系統(tǒng)的應用性能帶來了全新的突破，與傳統(tǒng)的雷達相比，MIMO雷達在參數目標估計、目標參數檢測、估計性能分析、空時自適應處理抑制、雷達干擾波形設計等方面有潛在的優(yōu)勢[1]. MIMO雷達在發(fā)射端能夠利用發(fā)射陣列的各個陣元發(fā)射彼此相互正交的信號，并且在接收端通過利用匹配濾波器進行濾波，實現各發(fā)射信號波形分離. 因而，MIMO雷達相較于傳統(tǒng)雷達具有能產生大量虛擬陣元的優(yōu)勢，并且有效拓展了MIMO雷達系統(tǒng)陣列孔徑，極大地提高了DOA參數估計的精度.

DOA估計的研究一直以來都被作為提高MIMO雷達性能研究中的一項重要任務. MIMO雷達DOA估計算法多由陣列信號處理算法進行改進而來，基本上可以分為參數化估計方法、非參數化估計方法和稀疏類估計方法. 參數化估計方法如多重信號分類算法(multiple signal classification, MUSIC)、旋轉不變子空間信號參數估計(estimation of signal parameters by rotational invariance technique,ESPRIT)算法；非參數化估計方法中經典方法有單頻率最小二乘算法、自適應迭代法等；稀疏類方法包括壓縮感知算法等[2-4].

MUSIC算法在應用于DOA估計時雖然有精度高的特點，但仍存在一些不可被忽視的問題，如計算量大，對相干及小信噪比信號無法進行有效估計等. 文獻[5]根據MUSIC算法得到功率譜函數，利用其具有在譜峰處出現斜率突變的特性，對MUSIC功率譜進行3次求導，以得到更高分辨率的譜峰，但在低快拍、低信噪比條件下性能較差，實際應用性不高. 文獻[6]利用接收數據陣列同接收數據的共軛陣列同步進行虛擬拓展，以形成陣列的二次虛擬拓展，然后利用拓展后形成的陣列導向矢量和協(xié)方差矩陣共同進行DOA估計，通過構造更多陣元的方法，極大地拓展了陣列孔徑，從而能夠獲得更高分辨率，但是利用不斷拓展矩陣的方法存在運算量大的缺點. 文獻[7]提出利用降維變換思想結合共軛ESPRIT算法的思想，構造數據矩陣，其原理是在最小二乘ESPRIT方法的基礎上進行角度估計，從而實現信號之間相干源解相關，獲得較好的DOA估計精度，該方法對干擾信號下的DOA估計有較好抗干擾效果，但并沒有顯著地提高DOA估計精度. 在實際情況中噪聲大部分是非高斯的，且一般具有顯著信號沖擊性，再加上相干信號源的存在，也能造成信號子空間的向量會散發(fā)到信號噪聲子空間去，影響估計精度.

本文結合擬合估計算法和矩陣預處理技術的特點，在此基礎上提出基于MIMO雷達擬合高精度DOA估計算法. 首先，構建十字型接收信號陣列得到回波信號，對協(xié)方差矩陣進行擬合，利用凸優(yōu)化進行求解，得到更接近真實值的接收信號協(xié)方差矩陣，從而提高MIMO雷達DOA估計精度；然后，針對相干信號和噪聲干擾問題，利用矩陣預處理技術實現縮小矩陣的特征值分布，使得有用信號能量更加集中，在相干信號和噪聲環(huán)境下也能夠獲得較好的估計值；最后，采用酉ESPRIT算法將復數域內的復雜計算轉換為實數域內的計算，完成DOA估計，達到降低算法復雜度的目的.

1 信號模型

對MIMO雷達系統(tǒng)進行仿真，收發(fā)陣列利用十字型陣列模型[8]，其中xOy平面上的參考點O為由交叉點2個相互垂直子陣共同構成；x軸和y軸各個半軸等間距地分布著Mx和My個陣元，并滿足dx=dy，發(fā)射陣列同時由2Mx+2My-3個陣元發(fā)射載頻和帶寬相同的一組正交信號S，即S=[s1,s2,…,s2Mx+2My-3]，從中選取的第n個發(fā)射信號為sn=[sn(1),sn(2),…,sn(L)]，式中n的取值為n=1,2,…,2Mx+2My-3，L代表一個脈沖周期內包含子脈沖數. 假設位于遠場空域存在k個目標，信號功率為p，其中第k個目標的角度位置信息(θk,φk)為與坐標軸投影后與x軸和y軸的夾角，(αk,βk)為其對應的空間俯仰角和方位角. 經過目標反射后可以在接收端得到2Mx+2My-3個正交發(fā)射信號波形，則接收端的信號可表示為

xq(n)=A(θk，φk)diag(χq)AT(θk，φk)S(n)+wq(n)

(1)

(2)

(3)

其中

(4)

yq=G(θ,φ)χq+nq

(5)

對yq進行列堆棧，即可由此得到Q次脈沖累加下接收到的回波數據Y，可表示為

Y=[y1,y2,…,yQ]=G(θ,φ)η+N

(6)

2 擬合高精度DOA估計算法

2.1 擬合估計

通過已知回波信號Y，可求得協(xié)方差為

R=E[YYH]=
G(θ,φ)diag(p)GH(θ,φ)+diag(σ)=
D+diag(σ)

(7)

其中：σ=[σ1,σ2,…,σM]T∈CM×1為噪聲功率矢量；D=G(θ,φ)diag(p)GH(θ,φ)為Hermitian Topelitz矩陣，在實際中考慮到采樣點數有限這一條件，通常對信號協(xié)方差矩陣采用近似估計

(8)

vec(ΔR)～CN(0,ΣR)

(9)

式中：vec(ΔR)為ΔR的矢量化運算；ΣR(RT?R)/L. 用替代R會存在較大的估計誤差，這也是子空間類算法直接在的基礎上求取信號子空間進行DOA估計效果不如稀疏類算法的主要原因. 同時對相干信源，會出現秩的虧損，需要額外進行空間平滑或去相干處理才能在的基礎上進行DOA估計[10]. 為了能夠進行更準確的協(xié)方差估計，利用文獻[9]中提出的協(xié)方差擬合準則

(10)

由式(7)可知，由于協(xié)方差矩陣R的理論值是未知由參數θ、φ、p、σ構成的非線性函數，利用協(xié)方差擬合準則的性質，可得

(θ*，φ*，p*,σ*)=arg minf(θ,φ,p,σ)

(11)

式中(θ*，φ*，p*,σ*)為未知參數的估計值.

以采樣點數大于陣元數為例，將式(11)代入式(10)，對未知參數的估計可以等效為

(12)

因為式(12)為不能直接求解非線性函數，可通過將其轉換為半正定規(guī)劃問題，即

(13)

2.2 矩陣預處理

矩陣預處理技術是在解決大型線性方程組和特征值很分散的矩陣線性方程組時，采用的一種較為方便的解法，其原理是通過對系數矩陣中特征值分布進行集中，實現降低矩陣的條件數的目的，達到調整矩陣的病態(tài)特性、加快收斂的效果. 把這一技術引入到DOA估計中，可以使矩陣內有用信號更加集中，即使在面對相干信號和噪聲情況下依然能夠準確估計出波達方向[11-13]. 其方法為：

對于方程Ay=b，對矩陣A進行特征值集中. 設λm(m=1,2,…,i)是A的i個特征值，λmin和λmax分別代表A的最小和最大特征值，A的條件數即為cond(A)=λmax/λmin，條件數的大小能夠反映矩陣病態(tài)程度，為得到更準確結果，通常使條件數的值越小越優(yōu). 假設要處理矩陣方程表達式為f(A)=a×A+b×I，式中a、b代表常數，原方程組表達式可變?yōu)?g(A)y=f(A)b. 此時可得處理后A的系數矩陣：g(A)=a×A2+b×A.

利用二次函數的對稱性，把x軸上較大的區(qū)間通過映射變?yōu)閥軸上的一個較小區(qū)間. 為將[λmin,λmax]經過二次函數g(y)=a×y2+b×y的作用能夠映射到較小區(qū)間上，這里將(λmax+λmin)/2作為二次函數的中心對稱軸，令

(14)

則

其中有

(15)

(16)

對式(16)進行計算求解得到

(17)

對條件數的改善程度為

(18)

(19)

2.3 基于酉ESPRIT算法的二維空間聯(lián)合估計

(20)

(21)

(22)

(23)

此時令

(24)

定義一組新的矩陣

(25)

將式(25)2個矩陣合并為

P=?u+j?v

(26)

將新生成的P矩陣進行特征值分解，取其特征值組成新的對角矩陣Σ，進而

(27)

根據(θk，φk)與(αk,βk)存在對應關系[16]

(28)

可以求出入射波的方向角與俯仰角的估計，即

(29)

2.4 算法性能分析

比較本文算法與引入的以下3種算法的運算復雜度. 本文所采用算法復雜度在進行擬合估計中較高，算法的復雜度近似為O[(1/4)(QM2+KM4)+2K(M5+4K2)]，文獻[8]RD-ESPRIT算法復雜度可以表示成O[(QM3+M6+2K2(M2-M+2K2)]，文獻[17] Unitary-ESPRIT算法的復雜度可以表示為O[(1/4)(QM4+M6)+2K[M2+6K3]，文獻[18] 中所采用RD-MUSIC算法的復雜度可以表示為O[QM6+(1/4)K2(2QM2+M4+K2)]. 其中Q為快拍數，K為目標個數. 本文算法通過采用酉變換，將復數域運算變?yōu)閷崝涤蜻\算，降低了算法的復雜度，與通過降維運算降低復雜度的RD-ESPRIT算法相比具有相近的復雜度，低于另外2種算法，可得本文算法在改進了估計精度的前提下，并沒有因此產生較高的算法運算復雜度.

3 仿真結果分析

實驗一為驗證本文所提出算法DOA估計的效果，實驗條件為在采樣點數為100，信噪比為5 dB，陣元數為Mx=4、My=4的條件下，設3個目標，與接收陣列俯仰角分別為70°、30°、60°，方位角對應為30°、50°、45°. 進行實驗仿真. 結果如圖2所示，易知估計值已經基本上接近真實值，可證明本文所提出的算法能夠實現對目標進行有效估計.

實驗二為驗證采樣點、陣元數、信噪比等因素變化下對本文算法的DOA估計精度的影響，這里引入文獻[8] RD-ESPRIT算法、文獻[17] Unitary-ESPRIT算法及文獻[18] RD-MUSIC算法3種算法作為對比，同等條件下對比4種算法的RMSE仿真結果. 仿真條件設為：1) 陣元數為Mx=4、My=4，信噪比為5 dB，當采樣點數在20～200變化時；2) 采樣點數為100，信噪比為5 dB，陣元數為Mx=4、My=4，信噪比變化范圍為-10～10時；3) 采樣點數為100，信噪比為5 dB，陣元數變化范圍為5～20時. 設置3個目標，與接收陣列的俯仰角和方位角分別為(70°,30°)、(30°,50°)、(60°,45°)，在每個變化的條件下做100次MontE-Carlo實驗. 取每個算法得到的俯仰角和水平角RMSE均值繪制得到圖3～5的仿真結果.

由圖3可得，當采樣點數增加時，4種算法RMSE曲線都會降低，代表 DOA估計精度提高；相比較于其他3種算法，本文采用的算法總體RMSE曲線較低，即DOA估計精度較高. 由圖4可得，算法的估計精度與信噪比成正相關，本文算法在不同的信噪比條件下都優(yōu)于其他3種算法，且在低信噪比的情況下估計效果改善明顯. 由圖5知，隨著陣元數變化，當陣元數大于4時4種算法的RMSE曲線相近，但本文算法估計精度略高于其他3種. 綜合分析可知，在非相干信號源的DOA估計中，當信噪比、采樣點數、陣元數都達到比較理想情況時，4種DOA估計算法精度相差不大，本文最優(yōu)，而且當仿真條件不在理想情況時本文所提出的算法DOA估計精度始終表現最好.

實驗三為了進一步驗證在相干信源情況下本文算法的DOA精度，引入3個發(fā)射信號源，設定第1個信號源和第3個信號源為相干信源，在其他仿真條件和參數與實驗二中相同的情況下，進行100次MontE-Carlo實驗的仿真得到圖6～8的仿真結果.

同實驗二中的實驗結果相比，相干信號源的環(huán)境下，4種算法的DOA估計精度都有所降低，但本文算法仍可以取得相對較好的DOA估計結果. 從圖6中可得出當采樣點數高于60時，本文算法即開始優(yōu)于其他3種比較算法；當采樣點數較低時，根據式(12)，協(xié)方差擬合準則應用于采樣點數大于2倍陣元數，所以估計結果誤差較大. 從圖7中得到，在信噪比較低時本文算法表現依舊優(yōu)于其他3種算法；圖8中陣元數在Mx=My=2～12變化時，可看到隨著陣元數的不斷增加，4種實驗算法的 DOA估計精度都在不斷提升，當陣元數大于等于4時本文提出的算法即可實現最優(yōu)；但陣元數小于4時，十字陣列性能未能發(fā)揮出來，導致估計性能較差.

由實驗結果可知，本文提出的算法在不同仿真環(huán)境下，都可以實現較高精度的DOA估計，并且在相干信源下相比其他算法也有較好的表現，從理論和仿真中證明了所提算法的可行性. 實際應用中考慮成本和環(huán)境問題均無法達到信噪比、陣元數和采樣點數最高的理論值，本文算法將在實際應用中DOA估計精度占優(yōu). 未來將會繼續(xù)研究在不影響估計精度的前提下如何降低算法復雜度，尋求在低信噪比和低快拍數的情況下如何獲取更高估計精度的方法，并對低角度情況下產生的多徑效應影響估計精度的問題進行研究.

4 結論

1) 相對于傳統(tǒng)MIMO雷達DOA估計算法，通過引入擬合估計算法，獲得較高的DOA估計精度.

2) 利用矩陣預處理技術，使得改進的算法在相干信號和噪聲環(huán)境下均能實現DOA的高精度估計.

3) 與稀疏類算法和MUSIC類算法相比，本文算法有很大優(yōu)勢，在提高估計精度的前提下，并沒有增加運算量.