王 凌 李國林 謝 鑫 齊 率

(海軍航空工程學(xué)院 煙臺 264001)

1 引言

近年來利用目標(biāo)空間中信源的信息冗余特性來提高空間譜估計(jì)的性能得到了廣泛關(guān)注。信號的非圓性作為時域上的信息擴(kuò)充，Galy首先將其引入空間譜估計(jì)研究中，提出了非圓信號測向的MUSIC算法[1](MUSIC algorithm for NonCircular signals,NC-MUSIC),Abeida和Delmas進(jìn)一步對非圓信號MUSIC類測向算法的Cramer-Rao Bound (CRB)進(jìn)行了分析[2]。近年來國內(nèi)外不少文獻(xiàn)對相關(guān)問題進(jìn)行了研究[3－12]。目前，研究較多的一類非圓信號為直線信號，例如幅度鍵控信號和二進(jìn)制相移鍵控信號等，其相位星座圖沿直線分布，非圓率ρ=1。文獻(xiàn)[3]中利用非圓信號為實(shí)值信號的特點(diǎn)，通過每次陣列快拍構(gòu)造了一組數(shù)據(jù)矩陣，實(shí)現(xiàn)了非圓相干信號的完全解相干。文獻(xiàn)[4,5]在構(gòu)造的 Toeplitz矩陣中均利用非圓信號的共軛信息，降低了陣列孔徑損失，提高了陣列自由度。文獻(xiàn)[6,7]利用非圓信號橢圓協(xié)方差E[ss]≠0特性，構(gòu)造虛擬陣元，提高了角度估計(jì)的精度。但上述算法模型背景中均假設(shè)非圓信號的初始相位為零，這在實(shí)際環(huán)境中幾乎不可能滿足。容易證明，非圓信號s(t)和其任意旋轉(zhuǎn)s(t)ejφ具有不同的概率分布，因此對于非圓信號來說，初始相位不能忽略。如果考慮初始相位，上述基于非圓信號特性的算法均會存在估計(jì)性能下降的情況，特別是文獻(xiàn)[3－5]中直接利用非圓信號的實(shí)值特性構(gòu)造數(shù)據(jù)矩陣，當(dāng)考慮到初始相位時，信號模型將不再滿足實(shí)值特性，對算法的影響是嚴(yán)重的。近幾年，部分文獻(xiàn)對上述算法進(jìn)行了修正[8,9]。

以上文獻(xiàn)都是針對1維DOA估計(jì)的非圓算法，針對非圓信號背景下的2維DOA估計(jì)問題，公開文獻(xiàn)中研究很少。僅有Haardt[13,14]等提出將2維酉ESPRIT用于非圓信號測向，以及劉劍[15,16]提出的基于雙平行線結(jié)構(gòu)的擴(kuò)展 MUSIC算法。但是上述算法運(yùn)算量都較大，需要構(gòu)造維數(shù)為33MM×的數(shù)據(jù)矩陣，并進(jìn)行特征分解。2維DOA和初相的聯(lián)合估計(jì)可以看作是在3維參數(shù)空間進(jìn)行點(diǎn)估計(jì)，基本思路是通過3維譜峰搜索實(shí)現(xiàn)，但運(yùn)算量巨大?；蛘吒鶕?jù)二次型的性質(zhì)，將多維搜索轉(zhuǎn)化為1維搜索或求根。

針對上述方法計(jì)算量大，且算法復(fù)雜的問題，本文提出了一種基于垂直陣列結(jié)構(gòu)的任意初始相位非圓信號2維DOA和初相聯(lián)合估計(jì)方法，利用陣列結(jié)構(gòu)特點(diǎn)將3維參數(shù)估計(jì)轉(zhuǎn)化為3個子陣上的2維參數(shù)估計(jì)問題，然后在每一個子陣列上，同時利用信號子空間的旋轉(zhuǎn)不變性和噪聲子空間的正交特性將2維估計(jì)轉(zhuǎn)化為1維估計(jì)，且得到的2維DOA和初始相位能實(shí)現(xiàn)自動配對。該算法只需對 3個22MM×維的數(shù)據(jù)矩陣進(jìn)行一次特征分解，且可在3個維度并行處理，相較于文獻(xiàn)[15]中算法，運(yùn)算量得到降低。由于利用了非圓信號特性，本文算法相較于具有相同陣列結(jié)構(gòu)的文獻(xiàn)[17]中算法，估計(jì)精度和可估計(jì)信源數(shù)都得到提升，且可估計(jì)多于陣元數(shù)的信源，最多可估計(jì)2(M?1)個信源。

2 陣列信號模型

考慮采用如圖1所示的傳感器陣列系統(tǒng)模型，該模型由3個相互垂直的均布線性子陣組成，記為，各子陣列陣元均為全向的，且陣列位于遠(yuǎn)場。第1個子陣由陣元0,1,…,M?1組成，第2個子陣由陣元0,1,…,N?1組成，第3個子陣由陣元0,1,…,P?1組成，各個子陣的陣元間距分別為dx,dy,dz?？臻g有D個相同載頻的不相關(guān)窄帶信源入射該陣列系統(tǒng)，信源數(shù)D假設(shè)為已知或已估計(jì)得到，波達(dá)方向矢量角分別為 (θ1,θ2,…,θD)，其中θk=(αk,βk,γk),αk,βk,γk分別為第k個源信號入射方向與x軸，y軸和z軸的夾角。顯然有如下關(guān)系成立

可知3個角度變量中只有兩個獨(dú)立，因此，波達(dá)方向矢量角可以用1×2的矢量表示，即θk=(αk,βk)。稱αk,βk分別為方位角和俯仰角。以陣元0為參考陣元，則子陣列Xa中第i個陣元的輸出信號可表示為

其中sk(t)為第k個直線信號，,λ為其中心波長，nxi(t)為子陣列Xa中第i個陣元中的復(fù)圓高斯白噪聲，噪聲功率為σ2，且噪聲與信源不相關(guān)。直線信號，其中φk為第k個直線信號入射到參考陣元的初始相位。文獻(xiàn)[3－7]中均假定φk為零，然后利用實(shí)值性進(jìn)行處理。本文考慮較為合理的情形，即各直線信號的非圓相位不全為零。

則子陣列Xa的輸出信號矢量可表示為

式中

為M×D維陣列流型矩陣，a(αk)為對應(yīng)的方向向量，且有

同理可以得到子陣列Ya和Za的輸出信號矢量

式中

3 本文算法描述

本文算法主線是利用空間信號入射到垂直陣列后的陣列流型特點(diǎn)，將2維DOA和初相聯(lián)合估計(jì)問題轉(zhuǎn)化為3個可并行處理的2維參數(shù)估計(jì)，然后在每一子陣列上，同時利用信號子空間的旋轉(zhuǎn)不變性和噪聲子空間的正交性，通過一次特征分解得到的信號和噪聲子空間將1維DOA和初始相位的2維搜索轉(zhuǎn)換為1維估計(jì)問題，這樣就把3維參數(shù)估計(jì)問題全部轉(zhuǎn)化為1維估計(jì)，使算法運(yùn)算量得到大幅的降低。

3.1 算法描述

受文獻(xiàn)[13,18]啟發(fā)，并利用非圓信號橢圓協(xié)方差非零特性，首先在x軸向構(gòu)造

式中，ΠΜ為M×M維置換矩陣，其反對角線上元素為1，其余元素為0。,從bxk(αk,φk)的表達(dá)式中可以看出，此時已將(α,β,φ)的3維估計(jì)轉(zhuǎn)化為(α,φ)的2維估計(jì)問題。改寫bxk如下式

式中⊙為哈達(dá)瑪積。從bxk的表達(dá)式不難發(fā)現(xiàn)，通過利用信號非圓特性構(gòu)造式(14)等效于在x軸向的負(fù)向端擴(kuò)展了一倍的陣列孔徑，且擴(kuò)展虛擬陣元的幅相增益為e?2jφk。按照上文思想，在y軸和z軸向同樣構(gòu)造擴(kuò)展的觀測矢量Gy(t)和Gz(t)，則經(jīng)過擴(kuò)展后的虛擬陣列結(jié)構(gòu)如圖2所示。

圖2 擴(kuò)展后的虛擬陣列系統(tǒng)模型

先計(jì)算x軸向的擴(kuò)展協(xié)方差矩陣

式中UxS為D個大特征值對應(yīng)的特征矢量張成的信號子空間，UxN為2M?D個小特征值對應(yīng)特征矢量張成的噪聲子空間。如果直接利用噪聲子空間的正交性進(jìn)行估計(jì)，即x軸向的方向向量bxk與噪聲子空間UxN正交，得到

由于噪聲影響，式(18)并不完全接近于0，此時初始相位和1維DOA的估計(jì)就轉(zhuǎn)化為最小化下式來實(shí)現(xiàn)

此時需要對α,φ進(jìn)行2維搜索，為了避免復(fù)雜搜索，本文考慮在利用噪聲子空間特性的同時使用信號子空間的旋轉(zhuǎn)不變性來降維。以x軸向?yàn)槔?，擴(kuò)展后的虛擬陣列如圖3所示。經(jīng)典的子陣劃分方法為圖中虛線所示。但由于非圓信號初始相位的影響，導(dǎo)致按圖中虛線劃分的子陣間將不再滿足旋轉(zhuǎn)不變性。由于虛擬陣元幅相增益的存在，可以考慮按實(shí)線方式進(jìn)行子陣劃分，此時

圖3x 軸向子陣劃分

其中V1和V2為選擇矩陣，同樣令

因?yàn)锽x2=Bx1Ψx，其中Ψx=diag[u(α1),u(α2),…,u(αD)]，根據(jù)信號子空間旋轉(zhuǎn)不變原理，易得：，其中Tx為一個非奇異轉(zhuǎn)換矩陣，可得

利用最小二乘或總體最小二乘方法求解式(24)得到兩子陣信號子空間之間的旋轉(zhuǎn)不變關(guān)系矩陣Θx，通過對Θx特征分解得到的Ψx即可求得空間信源的方位角αk，然后將得到的αk代入式(19)，通過1維搜索即可求得對應(yīng)于第k個信源的初始相位φk，這樣就把(α,φ)的2維估計(jì)問題轉(zhuǎn)化為兩個1維估計(jì)。

同理，對另外兩個軸向的旋轉(zhuǎn)不變關(guān)系Θy和Θz進(jìn)行特征分解可得到Ψy和Ψz，進(jìn)而可計(jì)算出各個信號與y軸和z軸之間的夾角βi和γj。2維角度參數(shù)的配對可以通過計(jì)算d(k,i,j)的最小值來實(shí)現(xiàn)，d(k,i,j)的定義式為

3.2 算法步驟

由上文分析，可得到本文算法步驟如下：

(1)通過陣列天線獲取垂直陣列的觀測矢量X(t),Y(t)和Z(t)；

(2)根據(jù)式(14)形式構(gòu)造擴(kuò)展觀測矢量Gx(t),Gy(t)和Gz(t)；

(3)計(jì)算擴(kuò)展協(xié)方差矩陣Rxx,Ryy和Rzz；

(4)對擴(kuò)展協(xié)方差矩陣進(jìn)行特征值分解，得到各維度上相應(yīng)的信號和噪聲子空間；

(5)利用式(24)中信號子空間的旋轉(zhuǎn)不變性計(jì)算信源與各軸之間夾角α,β和γ，并用式(25)進(jìn)行配對；

(6)利用上步得到結(jié)果代入式(19)中，根據(jù)噪聲子空間的正交性得到對應(yīng)的初始相位φ。

4 數(shù)值仿真

為驗(yàn)證本文算法的有效性，從3個方面進(jìn)行數(shù)值仿真驗(yàn)證，仿真中對比了具有相同陣列結(jié)構(gòu)的文獻(xiàn)[17]中算法。第1步驗(yàn)證本文算法對空間信源2維DOA和初始相位的聯(lián)合估計(jì)性能；第2步驗(yàn)證本文算法對快拍數(shù)和信噪比的敏感度。最后驗(yàn)證本文算法的抗過載能力，即多于陣元數(shù)的信源估計(jì)能力。

4.1 2維DOA和初始相位聯(lián)合估計(jì)驗(yàn)證

仿真1利用圖1中的陣列系統(tǒng)模型，取陣元數(shù)M,N和P均為8，為避免估值模糊，陣元間距dx,dy和dz均取為λ/2。3個空間信號的2維到達(dá)方向分別為(30°,90°,60°),(60°,60°,45°),(85°,75°,1 5 .9°)。其附加初始相位分別為ejπ/6,ej2π/3和ejπ/3。信噪比為15 dB，快拍數(shù)為500，進(jìn)行100次Monte-Carlo試驗(yàn)。圖4(a)和4(b)就是采用本文算法得到2維DOA估計(jì)的星座圖和初始相位的 MUSIC譜圖。從仿真結(jié)果可知，本文算法能夠準(zhǔn)確地估計(jì)得到2維DOA和初始相位，并能實(shí)現(xiàn)自動配對。

圖4 2維DOA和初相聯(lián)合估計(jì)

4.2 信噪比和快拍數(shù)對本文算法的影響

仿真1只是定性說明了本文算法的估計(jì)性能，下面通過兩個仿真試驗(yàn)驗(yàn)證本文算法的估計(jì)均方根誤差大小隨信噪比和快拍數(shù)的變化趨勢，并加入文獻(xiàn)[17]中算法作為比較。

仿真2陣列和信號模型與仿真1相同，文獻(xiàn)算法采用同樣模型參數(shù)?？炫臄?shù)固定為1000，改變信噪比，在0 dB到20 dB分別進(jìn)行50次的Monte-Carlo統(tǒng)計(jì)試驗(yàn)，得到如圖5所示的2維DOA估計(jì)性能與信噪比變化曲線。從圖中可以看出，本文算法對每個信號的估計(jì)均方根誤差都優(yōu)于文獻(xiàn)[17]中算法，并且隨著信噪比的提升，算法估計(jì)誤差逐漸減小。定義第k個信號的均方根誤差為

仿真3仿真參數(shù)同上。固定信噪比為5 dB，改變快拍數(shù)，在400次到1400次分別進(jìn)行50次的Monte-Carlo統(tǒng)計(jì)試驗(yàn)，得到如圖6所示的估計(jì)誤差與快拍數(shù)變化曲線。從仿真結(jié)果可得，隨著快拍數(shù)的增加，均方根誤差逐漸減小，在快拍數(shù)高于900次后，誤差趨于收斂?？傮w來看，本文算法在低信噪比和短快拍環(huán)境下，對2維DOA估計(jì)的均方根誤差能控制在1°以內(nèi)，且優(yōu)于文獻(xiàn)中算法。

4.3 陣列抗過載性能檢驗(yàn)

由于使用了利用非圓信號特性構(gòu)造的擴(kuò)展觀測矢量，因此本文能夠估計(jì)的空間信源數(shù)為2(M? 1 )，能夠?qū)Χ嘤陉囋獢?shù)的空間信號進(jìn)行估計(jì)，此時文獻(xiàn)[17]中算法已經(jīng)失效。

圖5 估計(jì)誤差隨信噪比變化曲線

仿真4取陣元數(shù)M,N和P均為4，4個空間信號的2維到達(dá)方向分別為(30°,90°,60°)，(60°,60°,4 5 °),(85°,7 5 °,1 5 .9°),(70°,2 5°,7 5 .6°)。其附加初始相位分別為ejπ/6,ej2π/3和ejπ/3,和ej5π/6。信噪比為20 dB，快拍數(shù)為500，進(jìn)行100次Monte-Carlo試驗(yàn)。圖7(a)和7(b)就是采用本文算法得到2維DOA估計(jì)的星座圖和初始相位的MUSIC譜圖。此時陣列處于過載狀態(tài)，傳統(tǒng)的子空間估計(jì)方法失效，本文算法仍能準(zhǔn)確得到結(jié)果。圖7(a)和7(b)中不同來波方向信源的2維DOA和初始相位估計(jì)誤差有區(qū)別的原因是仿真中采用了不同的信號形式。

5 結(jié)論

以往的2維DOA和初始相位聯(lián)合估計(jì)問題，一方面需要復(fù)雜的3維搜索，另一方面需要對高維度數(shù)據(jù)矩陣進(jìn)行特征分解運(yùn)算，針對傳統(tǒng)算法的復(fù)雜性，本文提出了一種基于垂直陣列結(jié)構(gòu)的任意初始相位非圓信號2維DOA和初相聯(lián)合估計(jì)方法，利用垂直陣列特點(diǎn)，將3維參數(shù)估計(jì)問題轉(zhuǎn)化為可并行處理的3個2維參數(shù)估計(jì)，在每一個子陣上，同時使用噪聲子空間正交性和信號子空間旋轉(zhuǎn)不變性，將2維參數(shù)估計(jì)進(jìn)一步轉(zhuǎn)化為1維估計(jì)問題，最終只需要對協(xié)方差矩陣進(jìn)行一次特征分解即可實(shí)現(xiàn) 2維DOA和初相的聯(lián)合估計(jì)及自動配對。由于使用擴(kuò)展觀測矢量構(gòu)造的協(xié)方差矩陣，本文算法適用于空間信源處于過載的情形和低信噪比、短快拍環(huán)境，因此適合多目標(biāo)和干擾下對2維DOA和初相估計(jì)要求較高的機(jī)載、彈載陣列天線系統(tǒng)和小型地面陣列天線系統(tǒng)。數(shù)值仿真驗(yàn)證了本文算法的上述優(yōu)良性能。

圖6 估計(jì)誤差隨快拍數(shù)變化曲線

圖7 過載環(huán)境下聯(lián)合估計(jì)