最優(yōu)冗余線陣欠定信號雙重構(gòu)DOA估計方法

王 庚, 李 浩, 何 翼, 馮明月, 陳昌孝,*

(1. 空軍預(yù)警學(xué)院, 湖北 武漢 430019; 2. 海軍研究院, 北京 100161)

0 引 言

在陣列信號到達角(direction of arrival, DOA)估計中,當信號數(shù)超過陣元數(shù),稱為欠定DOA估計。稀疏線陣可以實現(xiàn)欠定DOA估計,因此稀疏線陣DOA估計問題成為了研究熱點。早在1968年,Moffet就提出了最小冗余線陣[1],相同陣元數(shù)量條件下,在所有差值虛擬陣列連續(xù)的稀疏線陣中,最小冗余線陣孔徑最大。盡管最小冗余線陣已經(jīng)被證明可有效增加陣列自由度,但由于限制條件過于苛刻,多數(shù)孔徑條件下不存在相應(yīng)最小冗余線陣結(jié)構(gòu),因此難以用于實際的陣列設(shè)計中。直到2010年,Ma提出了利用Khatri-Rao積對于欠定信號DOA估計問題的解決思路[2],基于此,有學(xué)者陸續(xù)提出多種具有嚴格解析式的稀疏線陣,包括嵌套類稀疏線陣[3-5]、互質(zhì)類稀疏線陣[6-8]以及它們的混合型[9]。但相同陣元數(shù)量下,上述稀疏線陣孔徑都小于最小冗余線陣,而陣列測向精度主要取決于陣列孔徑,因此如何設(shè)計一種性能接近最小冗余線陣且孔徑連續(xù)的稀疏線陣,是稀疏線陣設(shè)計的一個研究方向。

稀疏線陣可以利用陣列接收數(shù)據(jù)構(gòu)造Toeplitz協(xié)方差矩陣,進而實現(xiàn)更高精度的DOA估計?；诓罘痔摂M陣列,文獻[10]通過空間平滑重構(gòu)出Toeplitz協(xié)方差矩陣,該方法計算速度快,但估計精度低。為提高Toeplitz協(xié)方差矩陣估計精度,文獻[11-12]提出了基于低秩矩陣去噪模型的協(xié)方差矩陣重構(gòu)算法,文獻[13-14]提出了基于跡范數(shù)最小化的協(xié)方差矩陣插值方法,文獻[15-17]提出了基于協(xié)方差矩陣擬合的重構(gòu)方法。上述Toeplitz矩陣估計方法的估計精度較高,但計算量大,因而如何在保持估計精度的同時,降低計算量是協(xié)方差矩陣估計的關(guān)鍵。

稀疏重構(gòu)算法無需信號個數(shù)信息即可實現(xiàn)欠定DOA估計,因此成為新興研究熱點。在文獻[18]中,Malioutov提出了基于L1范數(shù)的奇異值分解(L1 based singular value decomposition, L1-SVD)算法。在文獻[19]中,Liu結(jié)合估計誤差限制提出了協(xié)方差矩陣稀疏表示算法。文獻[20]將采樣協(xié)方差矩陣進行向量化處理并進行稀疏表示,結(jié)合誤差限制提出了基于L1范數(shù)最小化的稀疏重構(gòu)類算法。以上3種稀疏重構(gòu)類DOA估計算法直接利用稀疏線陣采樣協(xié)方差矩陣,因此都受協(xié)方差估計精度不足的限制,難以實現(xiàn)高精度欠定DOA估計。

針對上述陣列設(shè)計和DOA估計方法中的問題,本文提出了最優(yōu)冗余陣結(jié)構(gòu),實現(xiàn)了在保持低冗余低耦合的同時,便于陣列設(shè)計。為提高DOA估計精度,本文提出了快速協(xié)方差向量稀疏表示(fast covariance vector sparse representation, FCVSR)算法。該方法利用凸優(yōu)化理論近似最優(yōu)解條件,推導(dǎo)出了Toeplitz協(xié)方差矩陣的解析計算式。進而構(gòu)造無冗余協(xié)方差向量稀疏表示模型,通過稀疏重構(gòu)稀疏信號功率向量,最終實現(xiàn)了低計算復(fù)雜度下的高精度欠定信號DOA估計。

1 最優(yōu)冗余陣列設(shè)計

稀疏線陣可以視為對等孔徑均勻線陣進行陣元稀疏選擇后所得的等效線陣結(jié)構(gòu),本文將與稀疏線陣孔徑相等的均勻線陣稱為初始均勻線陣。

1.1 最優(yōu)冗余陣

陣元數(shù)為M的初始均勻線陣的陣元位置為

PU={ud, 0≤u≤M-1}={pu1,pu2,…,puM}

(1)

式中:d為相鄰陣元的最小陣元間距,通常設(shè)定為信號載頻的半波長；pui表示為第i個陣元的位置坐標。

擁有i個物理陣元的稀疏線陣的陣元位置為

PS={ps1,ps2,…,psS}

(2)

稀疏線陣的差值虛擬陣列定義和文獻[10]中相同:

PD={psi-psj=md}, ?i,j=1,2,…，S

(3)

相應(yīng)的差值虛擬陣元數(shù)和權(quán)值函數(shù)分別為

(4)

w(m)=|M(m)|

(5)

如文獻[21]所述,陣列耦合系數(shù)矩陣Cc模型為

(6)

式中: 1=c0

最小冗余線陣的差值虛擬陣列連續(xù),且相同陣元數(shù)條件下陣列冗余度最小。陣列的冗余度定義為

(7)

式中:ρmax表示陣列孔徑與單位陣元間距的比值。

現(xiàn)有文獻通過窮舉法搜索出了陣元數(shù)不超過17的所有最小冗余線陣,具體如表1所示。

表1 不同陣元數(shù)下最小冗余陣的陣列孔徑

由表1可知,最小冗余線陣的陣列孔徑不連續(xù),大部分孔徑無對應(yīng)最小冗余陣結(jié)構(gòu),因此很難用于陣列設(shè)計中?；谏鲜鲎钚∪哂嚓嚨亩x,本文提出了一種固定孔徑條件下冗余度最小的稀疏線陣。

定義 1固定陣列孔徑下,擁有連續(xù)差值虛擬陣列且冗余度最小的稀疏線陣為廣義最小冗余陣。

相同孔徑下存在多個廣義最小冗余陣結(jié)構(gòu),盡管這些廣義最小冗余陣具有相同的冗余度,但由于物理陣元的分布不同,陣列耦合效應(yīng)也不同。考慮到陣列耦合效應(yīng)對DOA估計的不良影響,本文提出了耦合影響最低的最優(yōu)冗余線陣。

定義 2同一孔徑下陣元耦合影響最小的廣義最小冗余陣為最優(yōu)冗余線陣。

由式(6)可知,陣元間的耦合影響主要與陣元間距相關(guān),隨陣元間距擴大而迅速降低,因此可以按照相鄰陣元間距由小到大的順序,比較相應(yīng)陣元間距下的陣元對數(shù),選擇其中相同陣元間距下陣元對數(shù)最少的廣義最小冗余陣,即為最優(yōu)冗余線陣。

本文提出的最優(yōu)冗余線陣只限于陣列孔徑為101d以內(nèi)小型稀疏線陣,40d以內(nèi)的具體陣元位置如圖1所示。最優(yōu)冗余線陣包含了最小冗余線陣。

圖1 最優(yōu)冗余線陣的陣元位置示意圖Fig.1 Schematic diagram of array element positions of optimal redundancy linear arrays

1.2 接收信號模型

由稀疏線陣和相應(yīng)初始均勻線陣孔徑相等的內(nèi)在關(guān)系可得psS=puM。

基于稀疏線陣和初始均勻線陣的陣元位置關(guān)系,可得稀疏選擇矩陣?！蕒0,1}S×M,其各列的計算公式為

(8)

式中:ρj∈{0,1}S×1只有第j個元素為1,其他為0。

假定PS={0,d,4d,7d,9d},則對應(yīng)稀疏選擇矩陣Γ為

(9)

假定K個互不相關(guān)的遠場窄帶信號分別以-π/2<θ1,…,θK<π/2的入射角度照射到初始均勻線陣上,則相應(yīng)的陣列輸出為

(10)

式中:A=[a(θ1),…,a(θK)]為陣列流型矩陣;a(θk)=[1,…,e-jπ(M-1)sin θk]T為陣列方向向量;信號s(t)和噪聲n(t)為非相關(guān)零均值高斯白噪聲,滿足E[s(t)s(t)H]=Rs和E[n(t)n(t)H]=δn,其中δn為噪聲功率,信號功率Rs=diag(δ1,…,δK)。

若考慮陣列耦合,則初始均勻線陣陣列輸出為

x(t)=CcAs(t)+n(t)

(11)

相應(yīng)地,稀疏線陣的陣列輸出表示為

(12)

式中:aCS(θk)代表考慮陣列耦合影響后的稀疏線陣的陣列方向向量。

2 基于雙重構(gòu)的FCVSR方法

基于稀疏線陣和初始均勻線陣之間的對應(yīng)關(guān)系,可以將稀疏線陣的欠定DOA估計問題轉(zhuǎn)化為初始均勻線陣的DOA估計問題。初始均勻線陣的協(xié)方差矩陣具有Toeplitz性,本文將其稱為Toeplitz協(xié)方差矩陣。

2.1 Toeplitz協(xié)方差矩陣重構(gòu)

基于信號和噪聲相互獨立的假定,可以得到初始協(xié)方差矩陣Rx的理論值為

(13)

相應(yīng)的無噪初始協(xié)方差矩陣為

Tx=Rx-δnIM=ARsAH

(14)

式中:IM為M維單位矩陣,如文獻[15]中所述;Rx和Tx同為Toeplitz矩陣和Hermitian矩陣。

相應(yīng)的稀疏線陣的協(xié)方差矩陣為

Ry=ΓARsAHΓH+δnIS=ΓTxΓH+δnIS

(15)

基于L次快拍的協(xié)方差矩陣的估計值為

(16)

(17)

因此可以推得

(18)

式中: Asχ2(|S|2)代表服從自由度為|S|2的卡方分布。由卡方分布的性質(zhì)可得

(19)

綜上所得,可以將無噪初始協(xié)方差矩陣Tx的求解問題轉(zhuǎn)換為下面的范數(shù)最小優(yōu)化問題[22]:

(20)

由于L0范數(shù)是非凸的NP難問題,為了便于實際計算,將L0范數(shù)凸松弛為跡范數(shù)。

(21)

當快拍數(shù)超過天線陣元數(shù)的3倍時,估計得到的協(xié)方差矩陣在很大概率上為半正定矩陣,因此去掉該限制條件,對Tx的估計效果影響不大。

2.2 快速計算方法

將跡范數(shù)最小化模型進一步轉(zhuǎn)化為拉格朗日形式為

(22)

將E=Ry-ΓTxΓT代入式(22)并進一步推導(dǎo)可得

(23)

T[J]=T[JTxJ]

(24)

式中：T[·]表示對矩陣進行壓縮處理,按照Toeplitz矩陣結(jié)構(gòu)將矩陣CM×M壓縮為向量C(2M-1)×1。

式(24)實質(zhì)為含有2M-1個變量的線性方程組。為求得該線性方程組的解,首先將等式右邊轉(zhuǎn)化為

(25)

(26)

式中:t=[t1, …,tM-1,tM]T。

為將線性方程組中的矩陣轉(zhuǎn)化為分塊矩陣,將一個零向量插入到矩陣H的第M列,從而得到下列等式:

(27)

(28)

式中:R(·)表示復(fù)數(shù)的實數(shù)部分;I(·)表示虛數(shù)部分。進而可得

(29)

最終可求得向量t,基于Tx的Toeplitz性和Hermitian性,可以得到

(30)

2.3 協(xié)方差向量稀疏表示模型

根據(jù)式(13)和式(30),可得向量t的元素值為

(31)

進而可以推導(dǎo)出:

(32)

通過對協(xié)方差向量t進行稀疏表示可得

(33)

其中估計誤差向量來源于信號采樣數(shù)的有限性,空域離散化后的陣列流型矩陣AG由向量a(θGi)組成。

利用凸優(yōu)化方法,基于L1范數(shù)對初始協(xié)方差向量進行稀疏重構(gòu),求解后可得稀疏功率向量和相應(yīng)信號的DOA估計值。具體優(yōu)化模型如下:

(34)

式中:門限β是用于約束初始協(xié)方差向量的估計誤差,可以查表或利用Matlab數(shù)學(xué)軟件中的命令函數(shù)chi2inv(1-w,M)求得,本文中,置信參數(shù)w與文獻[23]中相同,采用參數(shù)值0.001。

通過求解凸優(yōu)化問題,可以得到基于信號功率的稀疏向量δG,通過對δG進行峰值搜索,可以得到目標信號功率和相應(yīng)DOA的估計值。

綜上所述,基于雙重構(gòu)的FCVSR算法的詳細步驟如下所示:

步驟 2基于式(29)求得協(xié)方差向量t;

步驟 3基于協(xié)方差向量稀疏表示模型式(33),對信號功率稀疏向量δG進行稀疏重構(gòu);

步驟 4對估計出的稀疏向量δG進行峰值搜索,得到的峰值功率和相應(yīng)網(wǎng)格對應(yīng)的DOA值,即為信號功率和DOA的估計結(jié)果。

3 仿真實驗分析

為了對最優(yōu)冗余線陣測向性能、FCVSR算法的估計性能進行仿真分析,本文設(shè)計了兩個仿真實驗。假定K個等功率非相關(guān)遠場窄帶信號照射到稀疏線陣天線上,其中DOA均勻分布于區(qū)間[-53°,53°]。實驗中,峰值搜索步長為0.1°,網(wǎng)格間距設(shè)定為1°。

本文對每種仿真條件進行了I=500次蒙特卡羅實驗,并計算出相應(yīng)的均方根誤差(root mean square error, RMSE),具體公式為

(35)

3.1 稀疏線陣性能對比分析

為驗證最優(yōu)冗余線陣相對于其他等陣元數(shù)的稀疏線陣的耦合效應(yīng)更小,耦合條件下的欠定估計能力更強,本節(jié)針對陣元數(shù)為10的4種稀疏線陣,利用經(jīng)典的基于空間平滑的多信號分類(spatial smoothing-multiple signal classification, SS-MUSIC)進行欠定信號DOA估計。與文獻[21]仿真實驗中的仿真條件相同,假定耦合系數(shù)c1=0.5ejπ/4,c2=0.5ej0.7π/2,c3=0.5ej0.7π/3,B=3。仿真實驗中4種對比稀疏線陣的陣元位置如圖2所示。

圖2 4種稀疏線陣的陣元位置示意圖Fig.2 Schematic diagram of array element positions of four sparse linear arrays

4種不同仿真條件下的空間譜如圖3所示,其中信號個數(shù)、信噪比(signal to noise ratio, SNR)和快拍數(shù)的具體取值標記在每幅圖的下面,圖中的黑色虛線對應(yīng)真實DOA值。

圖3 陣元耦合時不同條件下4種稀疏線陣的歸一化空間譜Fig.3 Normalizational spatial spectrum of four sparse linear arrays under different conditions with mutual coupling

觀察圖3可知:① 對比4種稀疏線陣的歸一化功率譜可發(fā)現(xiàn),綠色部分功率譜中產(chǎn)生了最多的偽峰,說明嵌套陣列的耦合效應(yīng)最強;② 其次,紅色部分的估計結(jié)果也產(chǎn)生了較多的偽峰,藍色部分的估計結(jié)果與紫色部分大致相近,但當信號個數(shù)K=15時,紫色部分的估計結(jié)果要遠好于藍色部分,說明最優(yōu)冗余線陣較超級嵌套陣列,受耦合影響的程度更小,耦合條件下的估計欠定信號DOA估計能力更強,遠超廣義互質(zhì)陣列和嵌套陣列;③ 4種不同條件下,紫色部分的峰值與黑色虛線的距離最近,說明4種稀疏線陣中,最優(yōu)冗余線陣的測向精度更高。綜上,仿真結(jié)果驗證了最優(yōu)冗余線陣的優(yōu)越性。

3.2 條件變化對估計效果的影響分析

前文的欠定DOA估計功率譜僅為一次仿真實驗結(jié)果,為進一步分析FCVSR算法的欠定估計能力,并與已有算法進行性能對比,本文設(shè)計了不同接收SNR、快拍數(shù)和角度間隔下的蒙特卡羅仿真實驗,分析提出的FCVSR算法與已有L1-SVD[18],基于L1范數(shù)的陣列協(xié)方差向量稀疏表示(L1 norm based sparse representation of array covariance vector, L1SRACV)[23]和基于L1范數(shù)的協(xié)方差矩陣稀疏重構(gòu)(L1 norm based covariance matrix sparse reconstruction, L1CMSR)[19]算法的DOA估計精度受條件變化的影響,以及算法在不同SNR條件下的單次DOA估計運算時間。仿真結(jié)果如圖4所示,其中,當SNR從-10 dB增加到10 dB時,信號個數(shù)K=11,快拍數(shù)L=100;當快拍數(shù)從20增加到220時,信號個數(shù)K=11,SNR=0 dB;當信號個數(shù)從2增加到22時,SNR=0 dB,快拍數(shù)L=100。

圖4 估計效果隨不同SNR、快拍數(shù)和信號個數(shù)的變化情況Fig.4 Estimation effect varies with different SNR, number of snapshots and angular interval

觀察圖4可知:① 提出的FCVSR算法相對于其他3種算法的RMSE更低,說明其估計精度更高;② 通過提高協(xié)方差估計精度,可以實現(xiàn)更高精度的DOA估計;③ FCVSR算法與已有的L1CMSR算法的計算時間相近,小于L1SVD和L1SRACV。

4 結(jié)束語

本文提出了陣列孔徑連續(xù)的最優(yōu)冗余陣結(jié)構(gòu)和一種高精度欠定DOA估計方法。該方法通過構(gòu)造矩陣重構(gòu)模型,實現(xiàn)了協(xié)方差矩陣的高精度估計,進而提高了DOA估計精度。仿真結(jié)果證明,所提陣列的陣列耦合更低,相同陣元數(shù)下的測向精度更好;本文提出的算法較現(xiàn)有方法具有更高的欠定信號DOA估計精度,同時保持了較低的計算復(fù)雜度。該方法無需信號個數(shù)的先驗信息,通過提高協(xié)方差的估計精度,實現(xiàn)了更高精度的DOA估計,后續(xù)可以將該解決思路應(yīng)用到稀疏貝葉斯學(xué)習(xí)等稀疏重構(gòu)類方法中。