林新平 張小飛 沈金清

（南京航空航天大學電子信息工程學院，南京，211106）

引 言

近年來，信號波達方向(Direction of arrival,DOA)估計被廣泛應用于雷達、無線通信、聲吶和生物醫(yī)學等領域，是陣列信號處理中的一個研究熱點[1-3]。經典的超分辨率DOA估計算法，例如多重信號分類(Multiple signal classification,MUSIC)算法[4]、借助旋轉不變性進行信號參數估計(Estimation of signal parameters via rotational invariance technique,ESPRIT)算法[5]等，都是基于特征結構的子空間類方法。子空間類方法最初是針對傳統(tǒng)滿陣提出的，要求陣元間距嚴格小于等于電磁波半波長，以避免出現角度估計模糊。然而，傳統(tǒng)滿陣由于陣元間距過小，存在陣元間互耦嚴重、陣列孔徑不足以及測向精度不高等缺點，已經無法滿足許多航空航天以及軍事民用領域的應用需求。

文獻[6]提出了互質陣列的概念。它是一種陣元間距大于半波長的非均勻陣列，由兩個陣元數與陣元間距存在互質關系的均勻子陣穿插拓撲構成。相比于傳統(tǒng)滿陣，互質陣列具有陣元互耦更低，陣列孔徑更大，定位測向精度更高等優(yōu)點。文獻[6]中證明了一個具有M+N-1個陣元的互質線陣，能夠獲得O{MN}的空間自由度。因此，基于互質陣列的空間譜估計研究成為當下信號處理領域的熱點問題之一。

文獻[7]中提出了一種互質線陣下基于矢量化協(xié)方差矩陣的DOA估計方法，稱之為虛擬化方法。該方法通過矢量化接收信號協(xié)方差矩陣進行數據重構，從而得到一個通過虛擬陣列接收到的單快拍信號。特別地，虛擬陣列能夠提供比物理陣列更大的陣列孔徑和更高的空間自由度。文獻[8]中則提出了基于互質特性的聯合MUSIC方法，該方法利用互質線陣的兩個子陣單獨估計DOA，通過比對子陣間的估計結果消除測向模糊，稱之為解模糊方法。相比于虛擬化方法，解模糊方法實現簡單，且DOA估計性能更優(yōu)[9]。然而文獻[8]中的MUSIC算法需要全局譜搜索，算法復雜度較高。文獻[10]和文獻[11]中分別提出了MUSIC部分空間譜搜索的方法和求根MUSIC算法，復雜度相對于文獻[8]中的方法具有較大改善。相對于MUSIC算法，互質線陣下的ESPRIT算法[12]顯然在算法復雜度上要低得多。文獻[13-15]中通過對接收數據矩陣重構和實值化的方法達到降低互質陣列中ESPRTI算法角度估計復雜度的目的。然而，上述DOA估計方法只利用了接收信號的空域信息，并未對時域信息加以利用[16]。

為了能夠充分利用互質線陣接收信號的時域與空域信息，本文提出了一種基于共軛增廣(Conjugate augmented)的DOA估計算法。算法首先通過利用不同時長間隔下陣元接收信號間互相關函數的共軛對稱特性，構造一個共軛增廣虛擬陣列接收信號，虛擬接收信號具有和實際信號相同的信號形式，但相比于原有的互質陣列，虛擬陣列中增加了近一倍的鏡像虛擬陣元，陣列孔徑更大，空間自由度更高。然后采用基于互質特性的聯合酉ESPRIT(Unitary ESPRIT)算法[17]對虛擬信號實現DOA估計，稱該算法為基于共軛增廣的酉ESPRIT(CA-UESPRIT)算法。相比于傳統(tǒng)互質線陣下的聯合UESPRIT算法，CA-UESPRIT算法自由度更高，DOA估計性能更優(yōu)。CA-UESPRIT算法通過引入酉變換矩陣降低了ESPRIT算法中的復乘運算次數，算法復雜度較低。

1 陣列結構與數據模型

考慮如圖1所示互質線陣，它由兩個均勻子陣重合第1個陣元拓撲而成，子陣結構如圖2所示，其中子陣1和子陣2陣元數分別為M和N,陣元間距為d1=Nλ/2和d2=Mλ/2，λ為電磁波波長且M和N互為質數?；ベ|線陣總陣元數為M+N-1。

假設互質線陣置于x軸坐標系上，以圖1虛線位置為坐標原點，那么互質線陣陣元位置與坐標原點距離可表示為

圖1 互質線陣Fig.1 Coprime linear array

圖2 子陣拓撲圖Fig.2 Topology of subarray

假設K個遠場窄帶獨立平面波信號以角度θ=[θ1,θ2,…,θK]入射到上述互質陣列中，信號快拍數為J。在t時刻，假設第k個入射信號表示為sk(t)=Akejωkt，其中Ak為信號幅度，ωk表示信號頻率，k=1,2,…,K。那么互質線陣中第p(p=1,2,…,M+N-1)個陣元接收信號可以表示為[18]

式中：np(t)表示加性噪聲，lp表示集合L的第p個元素，根據式(1)可知l1=0。

陣列接收信號用矩陣形式可以表示為[19]

式中：s(t)=[s1(t),s2(t),…,sK(t)]T為信源矢量，n(t)為均值為0的加性高斯白噪聲，協(xié)方差均為，且噪聲與信號間互不相關。為互質陣列的方向矩陣，第k列表示為

2 基于共軛增廣的DOA估計算法

CA-UESPRIT算法實現DOA估計主要包含兩個步驟，先是利用不同時長間隔下陣元接收信號間互相關函數的共軛對稱特性，構造一個共軛增廣虛擬陣列接收信號，隨后通過基于互質特性的聯合UESPRIT算法[17]求解虛擬陣列的DOA估計結果。

2.1 共軛增廣過程

本節(jié)通過構造共軛增廣虛擬陣列以擴展互質線陣的陣列孔徑及提高空間自由度。根據式(2)，互質線陣中第1個陣元接收信號與第p個陣元接收信號關于時長τ(τ≠0)的互相關函數可以表示為[16]

式中，Rsksk(τ)表示成[16]

(τ)可以表示為[16]

根據式(7)，可以將式(5)重寫為

構造矢量Rs(τ)=[Rs1s1(τ),Rs2s2(τ),…,Rsksk(τ)]T，R(τ)=[Rx1x1(τ),Rx2x1(τ),…,RxM+N-1x1(τ)]T，可以得到

根據式(6)可知 Rs(τ)=Rs*(-τ)，因此可以得到 (R(-τ))*=A*Rs(τ)。令矩陣 A(2)和矢量 R(2)(τ)分別表示矩陣A和矢量R(τ)的第2～M+N-1行，那么R(2)(τ)=A(2)Rs(τ)，構造為

圖3 虛擬互質線陣Fig.3 Virtual coprime linear array

設圖3中虛線位置為坐標原點，可以看出，虛擬陣列關于原點鏡像對稱，稱之為共軛增廣虛擬陣列。其中，虛擬陣列陣元總數為2(M+N-1)-1。根據式(10)，構造偽快拍數據矩陣

式中：Ts為偽采樣間隔，偽快拍數為偽信號矩陣。那么y可以看作是虛擬互質陣列的虛擬接收信號值。

2.2 DOA估計

本節(jié)利用基于互質特性的聯合UESPRIT算法求解虛擬信號y的DOA估計結果。將圖3中虛擬陣列拆分為兩個均勻子陣，如圖4所示。

圖4中子陣1和子陣2陣元數分別為2M-1和2N-1,陣元間距分別為d1=Nλ/2和d2=Mλ/2。設y1和y2分別為子陣1和子陣2的接收信號值，那么y1和y2可以表示成

式中：Ae,1=[ae,1(θ1),ae,1(θ2),…,ae,1(θK)]∈ C(2M-1)×K和 Ae,2=[ae,2(θ1),ae,2(θ2),…,ae,2(θK)]∈C(2N-1)×K分別表示子陣1和子陣2的方向矩陣，Ae,1和Ae,2第k列分別表示為

圖4 虛擬陣列子陣拓撲圖Fig.4 Subarray topology of virtual array

式中：k=1,2,…,K?？梢缘玫教摂M子陣1和子陣2的接收信號協(xié)方差矩陣分別為R1=E{y1yH1}和

虛擬子陣中陣元間距均大于半波長，因此先利用UESPRIT算法分別得到兩個子陣的模糊DOA估計結果，再利用互質特性聯合子陣DOA估計值消除測向模糊。

2.2.1 UESPRIT算法

兩個虛擬子陣具有相似的陣列結構，因此，以子陣1為例詳細闡述UESPRIT算法過程，然后推廣至子陣2同理得到子陣2的DOA估計結果。

對于子陣1，引入酉矩陣為

式中：ΠM-1為反對角矩陣。酉矩陣滿足對協(xié)方差矩陣R1進行酉變換

式中：Re(·)表示取實值操作。酉變換后虛擬子陣1協(xié)方差矩陣變?yōu)閷崝稻仃?，方向矩陣由Ae,1變?yōu)?/p>

可知方向矩陣Ae,1滿足旋轉不變關系式

式中：φx=diag(e-j2πd1sin(θ1)/λ,e-j2πd1sin(θ2)/λ,…,I(2M-2)×(2M-2)]。J1Ae,1和J2Ae,1分別表示選取方向矩陣Ae,1的前(2M-2)行和后(2M-2)行。根據文獻[20]，酉變換后，方向矩陣由 Ae,1變?yōu)?Ω1，因此式(18)變?yōu)?/p>

對式(17)中實數協(xié)方差矩陣進行特征分解，取最大K個特征值對應的特征矢量構成信號子空間Es1，可知Es1與實數方向矩陣Ω1間滿足關系式Es1=Ω1T，其中T是一個維度為K×K的非奇異矩陣，將Es1=Ω1T代入式(19)可得

對式(20)中T-1φuT=(K1Es1)+K2Es1進行特征分解，可得特征值為φu的對角元素，記第k個特征值為αk，k=1,2,…,K，那么虛擬子陣1的角度估計值可以表示為

同理，對于虛擬子陣2，對協(xié)方差矩陣R2進行酉變換后得到R2,u，參考式(18—21)，可以得到虛擬子陣2的角度估計結果，k=1,2,…,K。注意到虛擬子陣1和虛擬子陣2陣元間距均大于半波長，因此得到的角度估計值均為模糊值，需要解模糊處理。

2.2.2 角度解模糊

對于第k個入射信源θk，假設虛擬子陣1和子陣2中分別存在估計模糊角度值θ'k和θ″k，并有以下關系式[8]

式中：θk為第 k個入射信源的真實值，k=1,2,…,K。根據 d1=Nλ/2和d2=Mλ/2，將式(14)和(15)分別代入式(22)和(23)，可以得到式中：k1，k2為整數，因為|sin(θk)-sin(θ'k)|＜2，|sin(θk)-sin(θ″k)|＜2，所以k1和k2的取值范圍分別為k1=-(N-1),…,0,…,(N-1)和k2=-(M-1),…,0,…,(M-1)。考慮式(24)和(25)中θ'k，θ″k和θ″k可互換，因此虛擬子陣1和子陣2關于θk的模糊值個數分別為N和M(包含真實值)。其中，聯立式(24)和(25)可以得到

因為M和N互為質數，當且僅當k1=k2=0時上式成立，也就是說，當且僅當θk=θ'k=θ″k時，式(26)成立?？梢缘贸鼋Y論，虛擬子陣1和子陣2中關于θk的所有模糊值中只有一個是相等的，也就是真實估計值。

2.3 所提DOA估計算法步驟

本文提出的CA-UESPRIT算法步驟如下：

(1)得到互質線陣中第p(p=1,2,…,M+N-1)個陣元接收信號xp(t)，利用不同時長間隔τ求陣元接收信號間的互相關函數Rxpx1(τ)。

(2)根據式(8—11)得到共軛增廣虛擬陣列接收信號y。

(3)利用UESPRIT算法分別對虛擬子陣1和子陣2接收信號y1和y2進行DOA估計，步驟如式(16—21)，得到兩個子陣的模糊角估計值

3 性能分析

3.1 最大可辨識信源數

根據文獻[1,2],一個具有N個傳感器的陣列最多能辨識N-1信源。因此，互質陣列解模糊方法下最大可辨識信源數取決于傳感器數較小的子陣。令T=Min(M,N)，其中Min(M,N)表示取M和N中的最小值。根據文獻[17]，傳統(tǒng)互質線陣下基于互質特性的聯合UESPRIT算法，最大可辨識信源數為T-1。然而，根據2.2節(jié)描述，本文提出的CA-UESPRIT算法通過共軛增廣特性使得虛擬子陣陣元數分別達到2M-1和2N-1，因此CA-UESPRIT算法的最大可辨識信源數為Min(2M-1,2N-1)-1，也就是2T-2。

3.2 復雜度分析

以復乘次數來評估DOA估計算法的復雜度，CA-UESPRIT算法主要運算復雜度包括：構造虛擬陣列接收信號y所需復雜度為O{Jw(2(M+N-1)-1)}，其中Jw為偽快拍數；得到虛擬子陣1和子陣2協(xié)方差矩陣R1和R2的復雜度分別為O{Jw(2M-1)2}和O{Jw(2N-1)2}；得到實數協(xié)方差矩陣R1,u和R2,u的復雜度分別為O{2(2M-1)3}和O{2(2N-1)3}。因此，總運算復雜度為O{Jw(2M+2N-3+(2M-1)2+(2N-1)2)+2(2M-1)3+2(2N-1)3}。本文利用UESPRIT算法實現DOA估計，不需要計算協(xié)方差矩陣特征分解以及后續(xù)總體最小二乘法的復雜度，因此復雜度相比于ESPRIT算法要低。

3.3 優(yōu)點總結

本文提出的CA-UESPRIT算法具有如下優(yōu)點：

(1)相比于傳統(tǒng)互質線陣下基于互質特性的聯合UESPRIT算法，CA-UESPRIT算法空間自由度更高。

(2)CA-UESPRIT算法具有比聯合UESPRIT算法更優(yōu)的DOA估計性能。

(3)CA-UESPRIT算法通過引入酉變換矩陣降低了ESPRIT算法中的復乘運算次數，算法復雜度較低。

4 仿真結果

本文采用1000次蒙特卡洛仿真的求根均方誤差(Root mean squares error,RMSE)評估DOA估計算法性能，定義RMSE為

圖5為本文DOA估計算法在信噪比SNR=0dB條件下角度估計圖。圖5中橫坐標為仿真次數，縱坐標為估計出的角度值。仿真參數設置為陣元數M=3,N=5，快拍數J=200。其中圖5(a)入射信源數 K=3, θ =[25°,30°,35°],圖 5(b)入射信源數 K=4, θ=[25°,30°,35°,40°]。從圖 5中可以看出，本文提出的CA-UESPRIT算法可以準確估計出信源角度值，且仿真結果驗證了3.1節(jié)中提出的本文算法可辨識信源數能達到2T-2=4，其中T=Min(M,N)。而根據文獻[17]，傳統(tǒng)互質線陣下基于互質特性的聯合UESPRIT算法，可辨識信源數僅能達到T-1=2。此外，對比圖5(a)和圖5(b),可以看出隨著信源數的增加，算法估計角度精確性變差，圖5(a)中RMSE值對應于不同仿真次數分別為0.0842,0.0831,0.0838,0.0866,0.0827，如圖 5(b)所示,RMSE 值分別為 0.1330,0.1295,0.1372，0.1276,0.1372?？梢缘贸鲂旁磾翟黾邮沟肦MSE性能降低。

圖6為互質線陣下不同DOA估計算法隨信噪比SNR變化下的RMSE性能對比圖，仿真參數均設置為陣元數M=3,N=5和K=2，θ=[20°,40°]。算法包括基于共軛增廣的CA-UESPRIT算法、CA-ESPRIT算法和CA-PM算法以及傳統(tǒng)互質線陣下基于互質特性的聯合UESPRIT算法[17]、ESPRIT算法[9]和PM算法。從圖6可以看出，隨著信噪比的增加，所提CA-UESPRIT算法的RMSE性能逐步變優(yōu)，也就是說，信噪比的增加帶來角度估計精度的提升。對比不同的DOA估計算法性能，可以看出，所提的CA-UESPRIT算法DOA估計性能與CA-ESPRIT算法十分相近，但CA-UESPRIT算法將協(xié)方差矩陣從復數域轉化到實數域，所以算法復雜度更低，且CA-UESPRIT算法DOA估計性能要優(yōu)于CA-PM算法。此外，相比于互質線陣下基于互質特性的聯合DOA估計算法[9，17],本文提出基于共軛增廣特性DOA估計算法由于構造的虛擬陣列信號具有更大的空間自由度，因此在RMSE性能全面占優(yōu)。

圖5 DOA估計結果散點圖Fig.5 Scatter of DOA estimation results

圖7 為互質線陣下不同DOA估計算法隨快拍數J變化下的RMSE性能對比圖，仿真參數均設置為陣元數M=3,N=5和K=2，θ=[20°,40°]，算法類型與圖6中相同。從圖7可以看出，隨著快拍數的增加，所提的CA-UESPRIT算法DOA估計性能提升。此外，對比不同的DOA估計算法，可以看出基于共軛增廣類的DOA估計算法具有更優(yōu)的RMSE性能。此外，所提的CA-UESPRIT算法DOA估計性能優(yōu)于CA-PM算法且接近于CA-ESPRIT算法，但CA-UESPRIT算法復雜度較低(具體分析參考圖6分析結果)。

圖6 互質陣列下不同算法RMSE性能對比圖(J=200)Fig.6 Comparison of RMSE performance with different algorithms for coprime array(J=200)

圖7 互質陣列下不同算法RMSE性能對比圖(SNR=0dB)Fig.7 Comparison of RMSE performance with different algorithms for coprime array(SNR=0dB)

圖8 為本文提出的CA-UESPRIT算法下RMSE隨信噪比SNR和快拍數J變化的性能對比圖，仿真參數均設置為陣元數M=3,N=5和K=2，θ=[20°,40°]。從圖8中可以看出，隨著信噪比SNR和快拍數J的增大，所提的DOA估計算法角度估計性能變優(yōu)。

5 結束語

本文提出了互質線陣中基于共軛增廣的DOA估計算法。該算法首先利用不同時長間隔下物理陣元接收信號間互相關函數的共軛對稱特性，構造虛擬陣列接收信號以擴展陣列孔徑和提高空間自由度。然后采用基于互質特性的聯合UESPRIT算法實現DOA估計。相比于傳統(tǒng)互質線陣下聯合的UESPRIT算法，本文提出的基于共軛增廣的UESPRIT(CA-UESPRIT)算法角度估計性能更優(yōu)，空間自由度更高。此外，相比于ESPRIT算法，CA-UESPRIT算法通過引入酉變換矩陣降低復乘運算次數，算法復雜度較低。探索更高效的DOA估計算法將會是下一步的研究方向。

圖8 不同快拍數條件下本文DOA估計算法RMSE性能對比圖Fig.8 Comparison of RMSE performance of DOA estimation algorithm under different snapshots