用于Nested陣列的多測量矢量DOA估計算法

王 璐，陶海紅，李 靖，裴 悅，智開宇

(1.西安電子科技大學 雷達信號處理國家重點實驗室，陜西 西安 710071；2.中國電子科技集團公司第五十四研究所，河北 石家莊 050081)

0 引言

波達方向 (Direction of Arrival，DOA) 估計是陣列信號處理領域的熱門方向之一[1]，其基本問題是確定同時處在空間某一區(qū)域內(nèi)多個感興趣的信號到達陣列參考陣元的方向角[2]。對于傳統(tǒng)均勻排布的陣元，采用子空間類[3]算法測向其自由度(Degree of Freedom，DOF)是有限的。對于陣元個數(shù)為N的均勻線性陣列(Uniform Linear Arrays，ULAs)，用傳統(tǒng)的子空間類算法測向其DOF為N-1，即DOF受到陣列孔徑的限制[4]。針對此情況，文獻[5]提出嵌套陣列，即Nested陣列，大大增加了自由度。

由于壓縮感知[6-7]技術可以突破Nyquist采樣定理的限制，僅通過少量快拍便可恢復出原始信號[8]，將壓縮感知技術引入DOA估計是大勢所趨[9-13]。文獻[14]提出正交匹配追蹤算法，文獻[15]提出基于奇異值分解的1范數(shù)重構(1Reconstruction after Singular Value Decomposition，L1-SVD)算法等。但是，匹配追蹤類算法的超分辨能力稍弱，1范數(shù)的稀疏約束能力大于MP算法，且基于1范數(shù)類方法在低信噪比、小快拍的現(xiàn)代復雜電磁信號環(huán)境中也能適應。

Nested陣列利用Khatri-Rao積操作[4]將協(xié)方差矩陣向量化后,可認為是單測量矢量(Single Measurement Vector，SMV)模型，傳統(tǒng)追蹤算法和正交匹配追蹤算法在該SMV模型下估計多個信源精度不高，由此本文在空間平滑算法[16]的基礎上提出了一種用于Nested陣列的多測量矢量(Multiple Measurement Vector，MMV)DOA估計算法，將SMV模型轉為MMV模型并采用L1-SVD算法完成DOA估計，最后實驗仿真分析比較3種算法在信源數(shù)多于物理陣元個數(shù)的情況下的測向性能。首先介紹了Nested陣列的信號模型，然后分別給出基于SMV模型和MMV模型的二級嵌套陣的稀疏表示，仿真比較了基于SMV模型稀疏表示的基追蹤算法和正交匹配追蹤算法以及所提算法的測向性能和均方根誤差，證明了在Nested陣列下所提算法相比于其他算法的優(yōu)勢。

1 Nested陣列的信號模型

1.1 Nested陣列

(1)

K級嵌套陣下第i級示意如圖1所示。

圖1 K級嵌套陣下第i級示意Fig.1 Schematic diagram of level i under K-level nested array

(2)

當嵌套級數(shù)大于2時，得到的差協(xié)同陣列將不再連續(xù)，即陣列間距不再連續(xù)，出現(xiàn)孔。故考慮二級嵌套陣的模型。

假設二級嵌套陣列總陣元數(shù)為N，第一級間距為d，陣元個數(shù)設置為N1，第二層間距為(N1+1)d，陣元個數(shù)為N2，則有N=N1+N2。二級嵌套陣物理陣列的陣元位置為：

S2={n1d,n2(N1+1)d},n1=1,2,…,N1,n2=1,2,…,N2。

(3)

通過Khatri-Rao積，二級嵌套陣列可以生成N2個陣元，其中不冗余的虛擬陣元個數(shù)為(N2-2)/2+N，間隔為d。N1=N2=3的二級嵌套陣列的陣元位置排布如圖2所示，通過Khatri-Rao積后形成的虛擬陣列及每個陣元出現(xiàn)次數(shù)如圖3所示。

圖2 物理陣元排布Fig.2 Physical array element arrangement

圖3 二級嵌套陣的加權函數(shù)Fig.3 Weighting function of two-level nested array

1.2 陣列模型

假設有D個遠場窄帶信號源以角度θi(i≤D)入射到N個陣元組成的線陣上，接收信號可表示為：

x[k]=As[k]+n[k]，

(4)

式中，A=[a(θ1)a(θ2)…a(θD)]為陣列流形矩陣，a(θi)為N×1維，為θi角度對應的導向矢量，第i個陣元的值為exp(j2πdisinθ/λ)；s[k]=[s1[k],s2[k],…,sD[k]]T表示源信號矢量，k=1,2,…,Ksnap,稱為快拍數(shù)，Ksnap表示陣元接收數(shù)據(jù)的總快拍數(shù)。

假設信號獨立且噪聲為加性高斯白噪聲。則有：

(5)

將Rxx向量化，可得：

(6)

2 基于MMV模型的L1-SVD算法

2.1 SMV向MMV模型轉化方法研究

在構造多快拍數(shù)據(jù)時，有2個基本要素需要注意：

(1) 轉換后的多快拍數(shù)據(jù)應該含有轉換前單快拍數(shù)據(jù)的所有信息，可以冗余但是盡量不要缺失。

(2) 因為多快拍數(shù)據(jù)是由同一組陣元上多次采樣得到的，導致多快拍數(shù)據(jù)的陣列流形是一致的。相對應的，在構造多快拍矩陣時也要注意多個快拍數(shù)據(jù)之間物理結構相同或者近似，方便之后選取測量矩陣。

針對這2點要素，可以開始陣元的選取。陣元的選取可以采用間隔選取或者平移選取等方式，只要各個快拍的陣元個數(shù)和位置都符合(2)的要求。借鑒空間平滑[16]的思想，將虛擬陣列接收數(shù)據(jù)劃分為重疊子陣，每個子陣數(shù)據(jù)看作一個快拍，構成新的接收數(shù)據(jù)矩陣的一列，但是這樣做是建立在規(guī)則陣列的基礎上的。與空間平滑算法目的是為解決協(xié)方差矩陣不滿秩不同的是，構建重疊子陣的目的更多是為了盡可能準確地構造多快拍的數(shù)據(jù)。首先，對式(6)的接收數(shù)據(jù)進行去冗余重排的處理，得到的陣列模型[4]為：

(7)

得到虛擬陣元數(shù)據(jù)后，可以將數(shù)據(jù)進行不重疊劃分或者重疊劃分子陣2種操作。2種操作理論上都是可行的，但是為了保證每個快拍的陣元之間的相位中心盡量近，避免模糊的同時又不丟失原虛擬接收陣列的信息，考慮重疊子陣的劃分方式。子陣之間盡可能重疊，即每次只向右滑一位直至最后一個陣元。陣元平滑示意如圖4所示。

圖4 陣元平滑示意Fig.4 Schematic diagram of array element smoothing

已經(jīng)確定好了平滑方式，下一步確定每一個子陣的陣元個數(shù)Q和構造的快拍數(shù)F。每個子陣陣元個數(shù)過少會影響陣列自由度，而快拍數(shù)在很多時候又會影響算法精度。在實際應用中，設計模型的快拍數(shù)和每個子陣的陣元個數(shù)可根據(jù)需要做調(diào)整。

(8)

(9)

(10)

(11)

式(9)中，Φ=diag(e-j(2π/λ)dsin θ1,e-j(2π/λ)dsin θ2,…,e-j(2π/λ)dsin θD)。Φf-1即是對Φ求f-1次冪。

(12)

(13)

至此，得到構造的多快拍數(shù)據(jù)模型。接下來，將應用MMV模型的經(jīng)典壓縮感知算法——L1-SVD算法來完成DOA估計。

2.2 Nested陣列下基于MMV模型的L1-SVD算法

Y=AobsS+Nnoise。

(14)

為了解決S在行項不存在稀疏性的問題，對S的每一行進行二范數(shù)處理，得到的結果為：

(15)

(16)

利用文獻[15]方法求解式(16)。二級嵌套陣下基于矩陣奇異值分解的1范數(shù)重構算法步驟如下：

步驟1：得到Ksnap次有限快拍下的物理陣列接收數(shù)據(jù)，x(k),k=1,2,…,Ksnap；

步驟2：求x(k)的采樣協(xié)方差矩陣Rxx；

步驟3：將Rxx向量化，得到z；

步驟6：按照式(13)構建稀疏字典Aobs；

步驟7：將Y按照文獻[15]進行矩陣奇異值分解；

步驟8：按照文獻[15]所示方法，利用CVX工具箱求解式(16)，所求結果行最大值對應角度即為所求角。

3 實驗仿真及性能分析

將式(7)模型稀疏化，令L為空間柵格個數(shù)，對p填充零元素，得到SMV模型為：

(17)

通過幾組仿真實驗驗證3種算法性能。比較在式(17)下采用文獻[13]的OMP算法、采用CVX工具箱求解的基追蹤算法、所提MMV模型下的L1-SVD算法3種算法的性能。定義均方根誤差(Root Mean Square Error，RMSE)為：

(18)

式中，D為信源個數(shù)；Mc表示蒙特卡羅試驗次數(shù)。

仿真實驗一：設置9個窄帶遠場非相干點頻信號源，分別以-60°，-45°，-30°，-15°，0°，15°，30°，45°，60°入射到陣面。信噪比均為10 dB，頻率分別為2.2，2.3，1.8，1.5，2.6，2.9，1.9，3.6，3.1 MHz，網(wǎng)格劃分間隔為0.5°，快拍數(shù)取1 024。噪聲假設為功率為1 W的高斯白噪聲。驗證信源個數(shù)多于物理陣元數(shù)目時算法估計性能，N=6時3種算法DOA估計結果如圖5所示。圖5中估計角度數(shù)值如表1所示。

圖5 N=6時3種算法DOA估計結果Fig.5 DOA estimation results of the three algorithms when N=6

表1 信源個數(shù)多于物理陣元時3種算法角度估計Tab.1 Angle estimation of three algorithms when the number of signal sources is more than that of physical array elements 單位：(°)

對比圖5和表1，當估計信源為9時，在上述仿真條件下，3種算法能夠對9個角度進行分辨，但是BP算法和OMP算法均存在誤差，而L1-SVD算法可以實現(xiàn)精確測向。統(tǒng)計3種算法在9個信號源下的均方誤差，OMP算法均方誤差約為1.013 8°,而BP算法均方誤差約為1.201 8°，L1-SVD算法均方誤差為0°。相較而言，在上述仿真條件下，3種算法測量精度排序為L1-SVD算法>OMP算法>BP算法。從信源估計個數(shù)分析，在同一仿真條件下所提算法對多個信源和單個信源均能實現(xiàn)準確估計，性能好于OMP算法和BP算法。

實驗仿真二：設置2個遠場窄帶點頻非相干信號源，入射角度分別為-5°，5°，頻率分別為1.3，2.4 MHz，快拍數(shù)為1 024，網(wǎng)格劃分間隔為0.5°。噪聲為功率為1 W的高斯白噪聲，信噪比變化范圍為-14∶1∶20。每個信噪比下進行500次蒙特卡羅實驗，3種算法RMSE隨信噪比變化如圖6所示。

圖6 3種算法RMSE隨信噪比變化Fig.6 Variation of RMSE of the three algorithms with the signal-to-noise ratio

由圖6可以看出，3種算法的曲線均呈現(xiàn)隨信噪比增加測向估計的RMSE降低的趨勢。其中，OMP算法在6 dB趨于穩(wěn)定，角度測向誤差趨于0.353 6°；BP算法在10 dB趨于穩(wěn)定，角度測向誤差趨于0；L1-SVD算法在0 dB趨于穩(wěn)定，且角度測向誤差趨于0。在上述條件下基于MMV模型的L1-SVD算法優(yōu)于BP算法和OMP算法，其中OMP算法在性能穩(wěn)定時還有0.353 6°誤差，性能為三者中最差。

考慮到上述仿真條件的SMV模型是Nested陣列對x(k)求取協(xié)方差矩陣向量化后去除冗余和重排得到，該模型噪聲將不再是高斯白噪聲且陣列數(shù)據(jù)可近似認為單快拍的相干數(shù)據(jù)，考慮到算法本身的因素，故SMV模型下2種算法測向精度不高。所提算法首先利用空間平滑算法，構造多快拍數(shù)據(jù)的同時對相干數(shù)據(jù)源做了處理，同時利用矩陣奇異值分解降低了單快拍轉為多快拍數(shù)據(jù)后的運算復雜度并削弱了部分噪聲，且考慮L1-SVD算法本身性能的因素，故性能相較于SMV模型下的BP算法和OMP算法得到提升。

實驗仿真三：設置2個窄帶遠場非相干點頻信號源，入射角度分別為-5°，5°，頻率分別為1.3，2.4 MHz，網(wǎng)格劃分間隔為0.5°，2個信源信噪比均為10 dB，噪聲為功率為1 W的高斯白噪聲。快拍數(shù)5∶50∶1 000變化。每個快拍數(shù)下進行500次蒙特卡羅實驗。3種算法RMSE隨快拍數(shù)變化如圖7所示。

圖7 3種算法RMSE隨快拍數(shù)變化Fig.7 Variation of RMSE of the three algorithms with the number of snapshots

由圖7可以看出，3種算法的曲線均呈現(xiàn)隨快怕數(shù)增加測向估計的RMSE降低的趨勢。其中，OMP算法在快拍數(shù)為255時趨于穩(wěn)定，角度測向RMSE趨于0.353 6°；BP算法快拍數(shù)從605開始角度測向RMSE在0值附近震蕩；L1-SVD算法快拍數(shù)從155開始角度測向RMSE在0值附近震動。在上述條件下所提算法優(yōu)于基于SMV模型的BP算法和OMP算法，其中OMP算法在性能穩(wěn)定時還有0.353 6°誤差，性能為三者中最差。

實驗仿真四：設置7個窄帶遠場非相干點頻信號源，入射角度分別為-43°，-28°，-14°，0°，18.5°，32°，44.5°，頻率分別為2.2，2.3，1.8，1.5，2.6，2.9，1.9 MHz，網(wǎng)格劃分間隔為0.5°，7個信源信噪比均為10 dB，噪聲為功率為1 W的高斯白噪聲?？炫臄?shù)50∶100∶2 000變化。每個快拍數(shù)下進行500次蒙特卡羅實驗。2種算法RMSE隨快拍數(shù)變化如圖8所示。

圖8 2種算法RMSE隨快拍數(shù)變化Fig.8 Variation of RMSE of the two algorithms with the number of snapshots

實驗仿真五：設置7個窄帶遠場非相干點頻信號源，入射角度分別為-43°，-28°，-14°，0°，18.5°，32°，44.5° ，頻率分別為2.2，2.3，1.8，1.5，2.6，2.9，1.9 MHz，網(wǎng)格劃分間隔為0.5°，快拍數(shù)取512，噪聲為功率為1 W的高斯白噪聲。信噪比-10∶2∶10變化。每個快拍數(shù)下進行500次蒙特卡羅實驗。2種算法RMSE隨信噪比變化如圖9所示。

圖9 2種算法RMSE隨信噪比變化Fig.9 Variation of RMSE of the two algorithms with the signal-to-noise ratio

MUSIC算法的復雜度體現(xiàn)在協(xié)方差矩陣的子空間分解，復雜度約為O(Q3) ，L1-SVD算法的復雜度約為O(D×L)3，復雜度高于MUSIC算法。但對比實驗仿真四和實驗仿真五可發(fā)現(xiàn),在快拍數(shù)較低如快拍數(shù)為50的情況下，L1-SVD算法的誤差約為3.412°，而MUSIC算法的誤差約為9.481°，證明低快拍數(shù)下L1-SVD算法性能好于MUSIC算法，但低信噪比的情況下，從圖9可知,L1-SVD算法性能弱于MUSIC算法。

4 結束語

本文從壓縮感知求解欠定問題的角度出發(fā)，將嵌套陣模型稀疏化，并用壓縮感知的經(jīng)典算法求解。本文比較了嵌套陣下基于SMV模型的BP算法、OMP算法和基于MMV模型的L1-SVD算法估計多個信號源角度情況，仿真表明，基于MMV模型的L1-SVD算法可準確估計信源個數(shù)多于物理陣元個數(shù)的目標信號的角度，且分辨率優(yōu)于其他2種算法。同時，從3種算法的信噪比和快拍數(shù)的RMSE曲線分析可知，基于MMV模型的L1-SVD算法性能優(yōu)于其他2種算法，且在小快拍數(shù)下估計多個信號源方向時，角度估計性能優(yōu)于MUSIC算法。本文方法完全可以擴展到二維Nested陣列，能夠實現(xiàn)對二維角度的估計。