基于稀疏陣列的電磁矢量傳感器多輸入多輸出雷達高分辨角度和極化參數聯合估計*

謝前朋 潘小義 陳吉源 肖順平

(國防科技大學, 電子信息系統復雜電磁環境效應國家重點實驗室, 長沙 410073)

針 對 雙 基 地 電 磁 矢 量 傳 感 器 多 輸 入 多 輸 出 (electromagnetic vector sensors multiple-input multipleoutput, EMVS-MIMO)雷達參數估計精度以及角度參數配對問題, 通過設計一種新的稀疏陣列和采用自動參數配對算法來實現高分辨的角度參數和極化參數聯合估計.首先, 通過設計稀疏的發射陣列和接收陣列來實現對EMVS-MIMO雷達陣列孔徑的擴展; 然后, 提出平行因子-三線性分解算法對接收數據的三階張量模型進行求解.所提出的平行因子-三線性分解算法能夠實現二維發射角、二維接收角、極化相位角和極化相位差的聯合參數自動配對; 且針對估計得到的發射導向矢量矩陣和接收導向矢量矩陣, 根據旋轉不變特性可以實現高精度的發射俯仰角和接收俯仰角測量.在得到精確的發射俯仰角和接收俯仰角之后, 相應的發射和接收方位角、極化角和極化相位差可以通過矢量叉積算法來進行估計.相比于現有算法, 所提出的算法能夠避免高維數據奇異值分解以及額外的參數配對過程; 且通過稀疏陣列設計, 角度參數估計精度能夠進一步地提升,仿真結果表明所提出的算法具有優良的角度參數估計性能.

1 引 言

多輸入多輸出(multiple-input multiple-output,MIMO)雷達參數估計是當前研究的一個熱點問題, 相比于使用相關波形的常規相控陣雷達, 其通過發射相互正交的信號波形能夠實現高的角度分辨率、靈活可控的波形設計以及陣列自由度的提升升[1?3].根據陣列的配置方式, MIMO 雷達可分為利用空間多樣性的統計MIMO雷達和利用波形多樣性的集中式MIMO雷達[4].并且, 集中式MIMO雷達可以進一步地劃分為單基地MIMO雷達和雙基地MIMO雷達.其中, 單基地MIMO雷達由于發射陣和接收陣距離較近, 相應的發射角(direction-of-departure, DOD)和接收角(directionof-arrival, DOA)是相同的; 而雙基地 MIMO 雷達由于發射陣和接收陣相距較遠, 其 DOD和DOA是不同的.本文主要針對雙基地MIMO雷達的角度參數估計展開研究.

近年來, 為了實現雙基地MIMO雷達的角度參數估計, 許多優良的算法被提出.在文獻[5,6]中, 譜峰搜索類算法二維Capon估計器和二維MUSIC估計器被提出.但是, 為了確保DOD和DOA的估計精度, 二維譜峰搜索類算法由于需要較小的搜索間隔, 具有較高的計算復雜度.為了降低計算代價, 相應的降維Capon和降維MUSIC算法在文獻[5,6]中也被提出.相比于二維譜峰搜索類算法, 降維類算法僅需要一維譜峰搜索過程.并且,降維類算法在降低計算復雜度的同時也能保持良好的估計精度.為了避免譜峰搜索過程, 文獻[7]利用旋轉不變技術 (estimated signal parameters via rotational invariance technique, ESPRIT) 來實現對雙基地MIMO雷達DOD和DOA的聯合估計.但是, 文獻[7]所提的算法需要額外的角度參數配對過程.在文獻[8]中, 一種修正的自動參數配對ESPRIT算法被提出.文獻[9]提出多項式求根MUSIC算法來實現對DOD和DOA的角度參數估計.文獻[10,11]提出聯合對角化算法來實現對雙基地MIMO雷達的角度參數的自動配對,并且所提出的算法相比于ESPRIT算法具有更高的估計精度.為了充分利用匹配濾波之后陣列接收數據的多維特性, 文獻[12]提出三階張量方法.文獻[13]利用高斯色噪聲的空時非相關特性, 通過構建延時相關的四階張量矩陣來實現對雙基地MIMO雷達中色噪聲的抑制.為了進一步地利用發射相關增益來提高角度參數估計精度, 文獻[14]利用波束空間變化技術來實現對發射波形進行約束, 從而把發射信號的能量和接收信號的能量集中到DOD和DOA所在的區域.針對沖擊噪聲背景下寬帶雙基地MIMO雷達的角度參數估計問題,文獻[15]利用Sigmoid變化來實現對聯合DOD和DOA角度參數估計性能的提升.針對雙基地MIMO雷達發射陣列以及接收陣列陣元損壞背景下的DOD和DOA估計問題, 文獻[16,17]利用圖像熵和低秩塊Hankel矩陣補全技術來實現對丟失數據的恢復.

盡管以上所提出的算法能夠實現良好的參數估計性能, 但它們主要針對標量發射陣列和標量接收陣列背景下的雙基地MIMO雷達角度參數估計問題.相比于標量陣列, 電磁矢量傳感器陣列(electromagnetic vector sensors, EMVS)不僅能夠提供角度信息, 同時也能夠提供極化信息.一個電磁矢量傳感器通常利用三個相互正交的電偶極子和三個相互正交的磁偶極子來實現對電場和磁場的測量[18].為了實現對雙基地EMVS-MIMO雷達的角度參數估計, 文獻[19]從發射EMVS陣列和接收EMVS陣列中提取出旋轉不變特性來實現對 2D-DOD和 2D-DOA的估計.但是, 文獻 [19]提出的ESPRIT-Like算法需要進行高維矩陣的奇異值分解.為了避免高維奇異值分解的計算復雜度, 文獻 [20]利用傳播算子 (propagator method,PM) 來實現對信號子空間的近似.文獻[21]進一步考慮利用EMVS-MIMO雷達陣列接收數據的多維特性, 提出基于協方差高階奇異值分解的2D-DOD和 2D-DOA 聯合估計算法.在文獻 [19?21]中, 面臨的共性問題在于為了實現2D-DOD和2D-DOA的參數配對, 需要進行構建額外的配對優化函數,并且, 在進行發射俯仰角和接收俯仰角參數估計的時候, 需要進行合適的參數選擇來實現旋轉不變關系的構建.如果不能選擇合適的旋轉不變參數, 那么估計得到的角度信息將會產生較大的誤差或者導致參數估計算法失效.因此, 為了進一步提升雙基地EMVS-MIMO的角度參數和極化參數的估計性能, 本文通過設計新型的稀疏發射EMVS陣列和稀疏接收EMVS陣列來實現對陣列孔徑的提升.同時為了避免2D-DOD和2D-DOA角度參數配對和額外旋轉不變參數的選擇, 本文提出利用平行因子-三線性分解 (PARAFAC-TALS) 算法來實現對角度參數和極化參數的估計.所提出的算法既能夠充分利用匹配濾波之后接收陣列數據的多維特性, 同時又能夠避免高維數據的奇異值分解.仿真實驗證明, 針對雙基地EMVS-MIMO雷達的角度參數和極化參數聯合估計, 所提出的算法具有較高的參數估計精度以及較低的計算復雜度.

2 信號模型

圖1 稀疏陣列 EMVS-MIMO 雷達系統Fig.1.EMVS-MIMO radar system with sparse linear array.

如圖1所示, 考慮一個包含M個稀疏EMVS發射陣列和N個稀疏EMVS接收陣列的雙基地EMVS-MIMO雷達系統, 其中發射稀疏陣列和接收稀疏陣列均包含兩個稀疏子陣.對于發射稀疏陣列, 第一個稀疏子陣的陣列個數為 M1且陣元間距為 M1λ , 第二個稀疏子陣的陣列個數為 M2且陣元間距為 M2λ , 其中l表示發射信號的波長.第一個發射子陣和第二個發射子陣的陣元間距為 M1λ.對于接收稀疏陣列, 第一個稀疏子陣的陣列個數為N1且陣元間距為 N1λ , 第二個稀疏子陣的陣列個數為 N2且陣元間距為 N2λ.第一個接收子陣和第二個接收子陣的陣元間距為 N1λ.因此, 發射稀疏陣列和接收稀疏陣列的陣元位置可以表示為

假設目標的個數為K, 則稀疏EMVS發射導向矢量和稀疏EMVS接收導向矢量為

對應接收俯仰角的導向矢量矩陣.ctk(θtk,?tk,γtk,ηtk)和 c rk(θrk,?rk,γrk,ηrk) 表示對應于發射陣列和接收陣列的電磁矢量傳感器的空間響應.

其 中 , Ftk(θtk,?tk) 和 Frk(θrk,?rk) 表 示 維 度 為6×2的空間角度位置矩陣, θtk,θrk∈ [0,π) 表示俯仰角, ?tk,?rk∈ [0,2π) 表示方位角; gtk(γtk,ηtk) 和grk(γrk,ηrk)表示維度為 2 ×1 的極化狀態矢量,γtk,γrk∈ [0,π/2)表 示 極 化 角 , ηtk,ηrk∈ [?π,π) 表示極化相位差.發射EMVS和接收EMVS的空間角度位置矩陣和極化狀態矢量可以表示為

于是, 根據矢量叉積算法, 歸一化波印廷矢量可以表示為

因此, 從公式 (11)和 (12)可以看出, 如果能夠得到歸一化波印廷矢量, 那么就可以對發射俯仰角、發射方位角、接收俯仰角和接收方位角的角度參數進行提取.并且, 從極化狀態矢量中可以提取相應的極化狀態角和極化相位差.

由于雙基地EMVS-MIMO雷達發射信號波形和接收信號波形的正交性, 匹配濾波之后的陣列接收數據可以表示為[19]

其中, At=[at1,at2,···,atK]和Ar=[ar1,ar2,···,arK]分別表示發射導向矢量矩陣和接收導向矢量矩陣, n (t) 表示加性高斯白噪聲矢量.對于L個采樣快拍, 總的陣列接收數據可以表示為

從公式(14)可以看出, 雙基地EMVS-MIMO雷達的陣列接收數據滿足多維張量結構, 如果直接對公式(14)進行協方差矩陣的求解, 則會破壞接收數據的空時特性.因此, 為了充分考慮發射陣列、接收陣列和采樣快拍之間的內在聯系, 這里采用張量結構來對陣列接收數據進行處理.

3 基于張量結構的PARAFAC算法

3.1 角度參數和極化參數聯合估計

為了實現對陣列接收數據的張量求解, 這里首先給出PARAFAC分解的定義[22]:

定義1[PARAFAC分解]一個維度為M × N ×L的三階張量X的PARAFAC分解定義為

其中, Dj(B) , Dk(C) 和 Di(A) 表示對角矩陣, 其對角線上的元素分別為加載矩陣A, B和C的第i, j和 k 行的元素.

因此, 根據 PARAFAC分解的定義, 公式(14)的陣列接收數據可以進一步地重新表示為

相應地, 關于 At和 Ar的聯立方程可以表示為

為了實現對發射導向矢量矩陣 At, 接收導向矢量矩陣 Ar和信號矩陣S的求解, 這里采用三線性迭代最小二乘算法.根據公式(17), (18)和(19),關于 At, Ar和S的最小二乘估計可以表示為

因此, 發射俯仰角 θtk,k=1,2,···,K 滿足如下的旋轉不變特性

表示第一個稀疏子陣的旋轉不變因子,

表示第二個稀疏子陣的旋轉不變因子.于是,Φt1(θt)和 Φt2(θt) 的估計可以表示為

進一步地, 估計得到的發射俯仰角的正弦值可以表示為

從公式(28)和(29)可以看出, 由于發射稀疏子陣的陣元間距大于半個波長, 所以以上估計得到的發射俯仰角是高精度周期模糊的.根據模糊周期與陣元間距的對應關系, 稀疏發射子陣1和稀疏發射子陣2的所有高精度稀疏模糊值可以表示為

因此, 為了得到高精度的無模糊發射俯仰角的正弦值, 對于

因此, 相應的高精度發射俯仰角可以被表示為

當利用解模糊算法得到所有目標的發射俯仰角之后, 可以利用矢量叉積算法來實現對歸一化波印廷矢量的求解.首先, 對于估計得到的發射導向矢量矩陣進行極化參數矩陣的提取.于是, 發射的電磁矢量傳感器陣列的空間響應可以表示為

對應于已估計得到的發射稀疏陣列的旋轉不變因子.根據公式(11), 發射EMVS陣列的歸一化波印廷矢量可以表示為

并且, 經過以上處理過程得到的發射角度估計參數和發射極化參數是自動配對的.因此, 額外的參數配對過程也能夠避免.

為了實現對接收角度參數和極化參數的聯合估計, 下面針對估計得到的接收導向矢量矩陣進行處理.首先, 為了得到高精度的接收俯仰角θrk,k=1,2,···,K, 定義如下的選擇矩陣:

從而可以得到滿足如下旋轉不變特性的接收俯仰角 θrk,k=1,2,···,K.

分別表示接收陣列的第一個稀疏子陣和接收陣列的第二個稀疏子陣的旋轉不變因子.因此, Φr1(θr)和 Φr2(θr) 的估計可以被表示為

進一步地, 估計得到的高精度接收俯仰角的正弦值可以表示為

同樣地, 由于接收稀疏子陣的陣元間距大于半個波長, 所以根據公式(46)和(47)得到的高精度發射俯仰角滿足周期模糊特性, 且針對稀疏子陣1和稀疏子陣2的所有的高精度稀疏模糊值可以表示為

因此, 為了解決以上稀疏陣列接收俯仰角的模糊問題, 對于

最終, 相應的高精度接收俯仰角為

當得到所有的接收俯仰角之后, 接收陣列的電磁矢量傳感器的空間響應可以表示為

最終, 經過以上的算法處理, 得到的對應于發射EMVS陣列和接收EMVS陣列的發射俯仰角、發射方位角、發射極化角、發射極化相位差和接收俯仰角、接收方位角、接收極化角、接收極化相位差能夠保證良好的估計精度.因此, 所提出的算法充分利用了稀疏陣列的孔徑擴展能力以及PARAFAC-TALS的自動參數配對特性.

3.2 克拉美羅界(Cramer-Rao Bound)以及計算復雜度分析

為了實現對所提算法性能的評價, 這里給出相應的CRB下界限以及所提算法的計算復雜度.CRB表示參數無偏估計的下界, 在衡量參數估計性能中具有重要的意義.對應于發射陣列和接收陣列 需 要 估 計 的 四 維 參 數 分 別 為 ( θt,?t,γt,ηt) 和(θr,?r,γr,ηr), 根據文獻 [23], 關于四維發射參數和四維接收參數的隨機CRB可以表示為

其 中, A =(At⊙Ar) 表 示 雙 基 地 EMVS-MIMO雷達的聯合發射接收矩陣,表示矩陣A的投影矩陣, ⊕ 表示Hadamard乘積,Rs表示多快拍信號協方差矩陣, 18×8表示維度為8×8的全1矩陣.D表示聯合發射接收矩陣A對四維發射參數 ( θt,?t,γt,ηt) 和四維接收參數(θr,?r,γr,ηr)的聯合導數矢量矩陣, 其具體形式如下:

因此, 通過以上詳細推導, 可以得到相應的CRB表達式.為了進一步驗證算法的性能, 所提算法的計算復雜度如表1所示.同時為了對比, 表1中也列出了文獻[19]中的ESPRIT-Like算法、文獻[20]中的PM-Like算法和文獻[21]中Tensor子空間算法的計算復雜度.四種算法的仿真軟件平臺為 MATLAB R2018a, 硬件平臺為 Lenovo ThinkPad X1 Carbon 筆記本, 其處理器為 Intel Core i7-6500U, CPU 為 2.5 GHz, 電 腦 內 存 為8 GB.在進行時間統計時, 信噪比設置為 1 0dB , 快拍數為100.運行時間能夠實現對各算法計算復雜度的直觀對比.從表1中給出的各算法的平均運行時間可以看出, 相比于其他算法, 本文算法具有較低的計算復雜度.如圖2所示, 進一步考慮10元發射陣列和10元接收陣列條件下, 目標數K=3時, 各種不同算法的計算復雜度隨快拍數L的變化關系.從圖2也可以看出, 所提算法相比于ESPRIT-Like算法、PM-Like算法和Tensor子空間算法具有較低的計算復雜度.Tensor子空間算法由于需要對四維矩陣進行奇異值分解, 因此在四種算法中具有最高的計算復雜度.

表1 不同算法的計算復雜度對比Table 1.Computational complexity comparison of different methods.

圖2 不同算法的計算復雜度隨快拍數的變化Fig.2.Comparison of computational complexity versus different snapshots number.

4 仿真實驗

為了驗證所提算法的有效性, 下面通過不同的仿真實驗來開展所提算法與文獻[19]中的ESPRITLike算法、文獻[20]中的PM-Like算法和文獻[21]中Tensor子空間算法的性能對比.如圖1所示,發射陣列和接收陣列均采用稀疏EMVS陣列.其中, 發射陣元個數 M =10(M1=3,M2=7) , 接收陣元個數 N =10(N1=3,N2=7).并且, 為了實現算法的有效對比, 利用ESPRIT-Like算法、PMLike算法和Tensor子空間算法進行參數估計時仍然采用文獻[19?21]中使用的均勻半波長EMVS陣列, 并設置相同的發射陣元數和接收陣元數.假設入射目標的個數 K =3 , 并且各個信號之間相互獨立, 相應的發射俯仰角、發射方位角、發射極化角、發射極化相位差和接收俯仰角、接收方位角、接收極化角、接收極化相位差如表2所示.其中表2中的第一個參數對應發射EMVS陣列, 第二個參數對應接收EMVS陣列.噪聲設置為相互獨立的零均值加性高斯白噪聲, 并且信號和噪聲之間相互獨立.

表2 目標回波參數表Table 2.Parameters of target signals.

4.1 算法的角度參數自動配對特性

首先, 通過星座圖來驗證所提算法的角度參數自動配對的有效性.在仿真中, 快拍數L設置為200, 信噪比設置為 1 0dB.利用100次蒙特卡洛仿真實驗結果來繪制如圖3所示的星座圖.從圖3(a)和圖3(b)可以看出, 所提出的PARAFAC-TALS算法能夠實現2D-DOD和2D-DOA的自動參數配對; 從圖3(c)到圖3(f)可以看出, 所提出的PARAFAC-TALS算法也能夠實現發射陣列四維參數和接收陣列四維參數的自動角度配對.因此, 所提算法的自動參數配對性能得到了有效的驗證.

4.2 信噪比對算法的影響

圖3 所提算法角度參數和極化參數估計星座圖 (a) 發射俯仰角和接收俯仰角; (b) 發射方位角和接收方位角; (c) 發射俯仰角和發射方位角; (d) 發射極化角和極化相位差; (e) 接收俯仰角和接收方位角; (f) 接收極化角和極化相位差Fig.3.Scatter plot of the angle parameters and polarization parameters by using the proposed method: (a) Scatter plot of the transmit elevation angle and receive elevation angle; (b) scatter plot of the transmit azimuth angle and receive azimuth angle;(c) scatter plot of the transmit elevation angle and azimuth angle; (d) scatter plot of the transmit polarization angle and polarization phase difference; (e) scatter plot of the receive elevation angle and azimuth angle; (f) scatter plot of the receive polarization angle and polarization phase difference.

下面進一步驗證信噪比對算法性能的影響, 其入射信號個數和快拍數與第一個實驗相同.信噪比的變化范圍為? 10dB 到 3 0dB , 變化的步長為 5 dB.在每個信噪比條件下蒙特卡洛仿真實驗次數為200.均方誤差的定義為其中表示估計得到的角度或極化參數, ? 表示真實的角度或極化參數, I表示蒙特卡洛仿真實驗次數.同時, 也給出了相應的檢測成功概率曲線.本文中檢測成功概率的準則和文獻[20,21]中給出的準則一致.從圖4(a)可以看出, 相比于ESPRITLike算法、PM-Like算法和Tensor子空間算法,所提算法具有最好的均方誤差估計性能, 其中的下標d表示角度參數, 下標p表示極化參數.從圖4(b)可以看出, 所提算法的角度參數檢測成功概率和極化參數檢測成功概率均好于其他三種算法, 且檢測概率基本上接近于1.

4.3 快拍數對算法的影響

為了評價采樣快拍數對所提算法性能的影響,入射信號個數和信噪比與第一個仿真實驗相同.快拍數的變化范圍為 1 00?1000 , 變化的步長為100.在每個信噪比條件下蒙特卡洛仿真實驗次數為200.同時檢測成功概率隨快拍數的變化如圖5(b)所示.從圖5可以看出, 本文所提算法的性能隨快拍數的增加表現出優良的估計精度和檢測概率.并且在快拍數比較小的情況下, 所提出的算法仍然能夠得到滿意的估計性能.

圖4 信噪比對算法的影響 (a) 均方誤差隨信噪比的變化; (b) 檢測概率隨信噪比的變化Fig.4.The effect of the SNR for different methods: (a) Curves of RMSE versus SNR; (b) curves of PSD versus SNR.

圖5 快拍數對算法的影響 (a) 均方誤差隨快拍數的變化; (b) 檢測概率隨快拍數的變化Fig.5.The effect of the snapshot for different methods: (a) Curves of RMSE versus snapshot; (b) curves of PSD versus snapshot.

4.4 角度分辨率分析

為了進一步驗證所提算法的角度分辨率性能,這里考慮兩個相鄰目標.第一個目標的發射四維參數和接收四維參數分別為 ( θt1,?t1,γt1,ηt1)= (40°,15°, 10°, 36°)和 ( θr1,?r1,γr1,ηr1)= (24°, 21°, 42°,17°), 第二個目標的發射四維參數和接收四維參數分別為(θt2,?t2,γt2,ηt2)=(40?+ ?,15?+ ?,10?+?,36?+?)和(θr2,?r2,γr2,ηr2)=(24?+?,21?+?,42?+?,17?+?).其中 D 的變化范圍為 1°—10° ,步長為 1°.快拍數為 200, 信噪比為 1 0dB , 每次角度間隔下的蒙特卡洛仿真實驗次數為200.從圖6可以看出, 本文所提算法具有較好的角度分辨性能.因此, 對于相鄰目標也能實現精確估計.

圖6 不同算法的目標分辨力比較 (a) 均方誤差隨角度間隔的變化; (b) 檢測概率隨角度間隔的變化Fig.6.Comparison of target resolution ability of different methods: (a) curves of RMSE versus angular separation; (b) curves of PSD versus angular separation.

5 結 論

本文的主要工作在于提出了一種新型稀疏雙基地EMVS-MIMO雷達系統, 有效地解決了雙基地EMVS-MIMO雷達的陣列孔徑擴展問題.同時利用平行因子-三線性分解算法解決了當前雙基地EMVS-MIMO雷達所面臨的2D-DOD和2DDOA的角度參數配對問題.本文利用稀疏發射EMVS陣列和稀疏接收EMVS陣列來構建針對于發射俯仰角和接收俯仰角的旋轉不變關系, 從而實現了高精度的角度參數求解.同時利用估計得到的發射導向矢量矩陣和接收導向矢量矩陣, 使得相應的方位角, 極化角和極化相位差可以通過歸一化波印廷矢量來求解.相比于ESPRIT-Like算法、PMLike算法和Tensor子空間算法, 所提出的算法能夠避免高維矩陣奇異值的分解以及角度參數配對過程中的高維譜峰搜索過程, 因此, 所提算法具有較低的計算復雜度.且通過仿真實驗可以發現, 本文所提算法在降低計算復雜度的同時能夠實現較高的參數估計精度和角度分辨率.在接下來的研究中, 為了進一步提升雙基地EMVS-MIMO雷達的角度參數和極化參數估計性能, 將主要圍繞設計靈活的稀疏發射EMVS陣列和稀疏接收EMVS陣列來展開.