基于樣本協方差矩陣譜分離特性的波達方向估計方法

郭 拓， 王英民， 張立琛

(西北工業大學 航海學院，西安 710072)

波達方向估計(Direction of Arrival，DOA)是雷達、聲納、通信、無源定位等陣列信號處理中的一個重要研究內容之一。陣列波達方向估計經過近幾十年的發展，已經形成了兩類技術，即基于波束形成的技術和基于子空間分解的技術。基于波束形成技術的典型方法有常規波束形成器(Conventional Beam Forming，CBF)和最小方差無畸變響應(Minimum Variance Distortionless Response,MVDR)波束形成器；多重信號分類(Multiple Signal Classification，MUSIC)算法是經典的子空間類方法[1]。上述DOA估計算法都需要估計基陣接收到數據的統計協方差矩陣，然而在實際數據處理中都以樣本協方差矩陣代替統計協方差矩陣，在小快拍情況下估計出的統計協方差矩陣誤差很大，將會導致傳統DOA估計算法不能準確估計目標方位。

為了獲得與統計協方差矩陣更加相近的樣本協方差矩陣，通常情況下會對樣本協方差矩陣進行修正，主要的修正方法包括：對角加載方法[2]，收縮因子方法[3]等，通過對樣本協方差矩陣的特征值進行加權操作以改進其方位估計的分辨性能；前后向平均技術[4]，空間光滑技術[5]，通過移動窗口虛擬出更多的快拍以達到改進樣本協方差矩陣估計的目的；偽隨機噪聲再采樣技術[6-8]則是人為生成偽隨機噪聲信號加在接收到的信號上，從而對樣本協方差樣本協方差矩陣進行修正；而子空間泄漏分析方法[9]總共有兩個步驟：①采用根MUSIC算法估計目標方位；②根據估計所得的目標方位修正樣本協方差矩陣。

上述這些方法的思路都是設法改進樣本協方差矩陣，以達到其與大快拍數時的樣本協方差矩陣更加相近，而沒有考慮樣本協方差矩陣與快拍數的關系，以及不同快拍數時協方差矩陣本身所具有的不同特性。

本文將不沿著上文所述那樣通過改進樣本協方差矩陣來改進方位估計效果的思路，而是轉換研究思路，研究不同快拍數時樣本協方差矩陣本身具有的一些特性，比如譜的特性，具體的就是通過深入分析相干目標與獨立目標的樣本協方差矩陣在不同快拍下所具有的各異的譜分離(Spectral Separation)特性，提出一種基于小快拍條件下基于樣本協方差矩陣譜分離特性的DOA估計方法。

1 陣列接收信號模型

如圖1所示，水平線陣由N個陣元組成，假設有D個相互獨立的平面波入射，陣列第k次快拍X(k)的輸出如式(1)所示。

圖1 陣列信號模型Fig.1 Array signal model

(1)

式中：X(k)=[x1(k),x2(k), …,xN(k)]T，A(θ)=[a(θ1),a(θ2),…,a(θD)]；a(θi)=[1, ej2πd cos θi/λ, …, ej2π(N-1)dcos θi/λ]為第i個入射平面波的陣列流形；陣元間距為d，d為波長λ的一半；θi為第i個平面波入射方向與水平面的夾角；A(θ)為陣列流形矩陣；S(k)=[s1(k),s2(k), …,sD(k)]T,si(k)為第i個平面波信號的第k次快拍；其中T為轉置運算。N(k)=[n1(k),n2(k), …,nN(k)]為陣列接收到的加性高斯白噪聲。

當信號被假設具有零均值時，陣列接收信號的協方差矩R=E[X(k)·X(k)H]，即協方差矩陣需要通過求集總平均而得到。而實際應用中只能得到有限的采樣數據，如有M次觀測樣本X1,X2，…，XM用這M個樣本來估計樣本協方差矩陣如式(2)，H為復共軛轉置。

(2)

2 基于小快拍協方差矩陣的主特征空間波達方向估計方法

本節將首先深入分析當目標信號分別為相干信號與獨立信號時樣本協方差矩陣的譜分離特性，然后在此基礎上提出樣本協方差矩陣主特征空間的概念，最后提出基于小快拍的主特征空間波達方向估計方法。

2.1 相干信號樣本協方差矩陣特征的譜分離特性分析

從子空間的角度來解釋上述相干源情況下譜分離特性，就是在y較小時，樣本協方差矩陣得到較好的估計，相干信號的信息全部集中在一個最大特征值所對應的子空間上；而在y較大時，如y為1此時陣元數等于快拍數，樣本協方差矩陣沒得到很好的估計，此時信號子空間發生了泄漏，使得最大特征值之前出現一些連續的較大特征值，但是總體而言相干信號的主要信息還是集中于最大特征值處。本文將這個最大特征值所對應的子空間稱為主特征空間(Main Feature Space, MFS)，而不是MUSIC算法所假設的那樣，即信號子空間就是s個信源所對應的s個較大特征值所構成的空間，通過我們上述對相干源情況下不同快拍數時譜分離特性的仿真與分析，對于相干源來說顯然信號信息主要集中在主特征空間而不是MUSIC算法所使用的s個較大特征值所構成的空間，這就很好的解釋了MUSIC算法對相干目標DOA估計時為何失效。

(a) 相干源譜分離性

(b) 相干源譜分離性局部放大圖2 相干源譜分離性分析Fig.2 Analysis of spectra separation for coherent sources

2.2 獨立信號樣本協方差矩陣特征的譜分離特性分析

圖3 獨立源的譜分離性Fig.3 Spetral separation of independent source

由“2.1”節分析可知，對于相干源，主特征空間中總是有個信源的信息，只是在快拍數小時，主特征空間中有一定信號信息泄漏入非主特征空間中，做DOA估計時會導致旁瓣會比快拍數大時的高。

對于獨立源，在小快拍情況下，主特征空間中會有所有目標信源的信息；快拍數大時，主特征空間中會完全不存在一部分源的信息，導致DOA估計中沒有這些目標。

總之不管是相干源還是獨立源，小快拍時主特征空間中包含所有信號的信息(這里的小快拍是指快拍數與陣元數相當，具體是指快拍數是陣元數的2倍以下，即y>0.5的情況)，所以可以利用小快拍時樣本協方差矩陣的主特征空間進行方位估計。下文將根據上述分析的譜分離特性提出基于主特征空間的DOA估計方法。

2.3 基于主特征空間的DOA估計方法

(3)

(4)

根據子空間理論，信號子空間與噪聲子空間正交，且導向向量張成的空間與信號子空間張成的空間是同一空間，即Ena(θ)=0，Esa(θ)，其中En為噪聲子空間，Es為信號子空間。

主特征空間中包含的主要是信號信息，即主特征子空間Emain中主要是Es中的信息。根據分析，在小快拍情況下，不管是相干源還是獨立源，主特征空間中都包括所有信號的信息。這樣，根據導向向量與信號子空間形成的是同一空間，則Emaina(θ)→1。在極限情形下如果噪聲子空間En的信息沒有泄漏入Emain，則Emaina(θ)=1，那么arccos(Emaina(θ))=0。因此利用導向向量與信號子空間形成的是同一空間這一特性，定義arccos(Emaina(θ))的倒數為主特征空間方法的DOA估計幅度，且使用Emain與a(θ)的l1范數來歸一化，即

(5)

3 仿真實驗

3.1 相干信號DOA估計

由于在水池實驗中陣元數有限采用了被動合成孔徑技術進行了孔徑擴展，在仿真中采用被動合成孔徑算法ETAM[14](Extended Towed Array Measurements)。仿真采用一個半波長分布的6元均勻線陣，陣列運動速度為0.07 m/s，入射的信號為4個遠場窄帶信號，頻率為6.25 kHz，采樣頻率50 kHz，4個信號的入射角分別為-20°，-10°，0°，10°，信噪比SNR皆為0，參考聲速為1 500 m/s，噪聲為復高斯白噪聲，快拍數分別取30、12 000，ETAM算法經過8次合成，合成24個虛擬陣元，加6個物理陣元，則總線陣為30個陣元。

合成孔徑后30元陣取30個快拍和12 000個快拍時，所提出算法給出的目標DOA估計對比，如圖4所示。由圖4可知，無論快拍數大小，主特征空間方法都可以分辨出相干目標，只是在大快拍情況下旁瓣會低一些。圖5給出使用相同的數據，采用MUSIC算法獲得的目標DOA估計。由圖5可知，MUSIC算法在小快拍情況下，會出現偽目標，不能正確分辨目標。比較可知，本文提出的方法可以在小快拍情況下，對相干目標

圖4 相干源主特征空間DOA估計Fig.4 DOA of coherent source with main feature space

圖5 相干源MUSIC算法DOA估計Fig.5 DOA of coherent source with MUSIC

進行很好的分辨，且不需要像MUSIC算法那樣需預先擁有目標個數的先驗信息，這是相干源時，樣本協方差矩陣的譜分離特性帶來的好處，即信號信息全部集中于主特征空間，只有很少的信號泄漏于非主特征空間中。

3.2 獨立信號DOA估計

本節仿真條件與“3.1”節相同，快拍數分別取30， 12 000，使用ETAM算法經過8次合成，合成24個虛擬陣元，加6個物理陣元，則總線陣為30個陣元，信號源改為獨立源，頻率為6.24 kHz，6.245 kHz，6.25 kHz和6.254 kHz ，入射角還是-20°，-10°，0°，10°。主特征空間方法DOA估計，如圖6所示。

由圖6可知，當目標為獨立信源，小快拍(30快拍)情況下，主特征空間方法可以很好的估計出目標方位；而在大快拍(12 000快拍)情況下，會發現不能正確分辨出所有目標，有一個目標(0°方向的目標)沒能估計出。此仿真驗證了“2.2”節所分析的結論，即小快拍時，主特征空間中會有所有目標信源的信息，大快拍時，主特征空間中會完全不存在某些源的信息，導致DOA估計中沒有得到這些目標的方位。圖7為與圖6使用相同數據情況下，MUSIC算法在不同快拍情形下的DOA估計，發現與圖6的情況恰好相反，是在小快拍(30快拍)時不能準確估計出4個目標方位，而在多快拍時能準確估計出目標方位。

圖6 獨立源主特征空間DOA估計Fig.6 DOA of independent source with main feature space

圖7 獨立源MUSIC算法DOA估計Fig.7 DOA of independent source with MUSIC

4 水池驗證

為了驗證本文所提算法，在消聲水池(長20 m,寬8 m,深7 m，6面消聲)進行了相關實驗。由于可用陣元數有限，故使用合成孔徑方法進行孔徑的擴展，然后進行DOA估計。

實驗采用由6個水聽器組成的均勻線陣，設計頻率為6.25 kHz，水平布放于水下2.5 m深處，陣元間距為半波長0.12 m，陣列勻速運動速度為0.07 m/s，采用電機勻速拉動。相干源時，3個發射源發射頻率皆為6.25 kHz的連續CW(單頻)信號；獨立源時，3個源頻率分別為6.24 kHz，6.25 kHz和6.254 kHz；與接收陣在同一深度，3個源按遠場條件布放，與陣形成的角度分別為24°，33°，44° ；采樣頻率為50 kHz；參考聲速為1 500 m/s。ETAM算法合成6次，得到18個虛擬陣元，加6個物理陣元，共得到24個陣元，快拍數分別取24,12 000。實驗系統如圖8所示。

圖8 水池實驗系統Fig.8 Tank experiment system

4.1 相干信號DOA估計

圖9為水池實驗對相干信號使用主特征空間DOA估計的結果，正如“2.1”節相干信號樣本協方差矩陣譜分離特性的分析以及“3.1”節仿真驗證，無論快拍數大小，主特征空間方法都可以分辨出相干目標，只是在大快拍情況下旁瓣會低一些，圖10中快拍數為12 000時比24個快拍時DOA估計旁瓣低。圖10為水池實驗對相干信號使用MUSIC算法DOA估計的結果，由圖10可知，在小快拍情況下，MUSIC算法有偽目標出現；在大快拍情況下，才可以勉強分開相干目標。

圖9 水池實驗相干源主特征空間DOA估計Fig.9 Tank experiment DOA of coherent source with main feature space

圖10 水池實驗相干源MUSIC算法DOA估計Fig.10 Tank experiment DOA of coherent source with MUSIC

4.2 獨立信號DOA估計

水池實驗采用主特征方法對獨立信號DOA估計結果，從圖11可知，在小快拍情況下，3個目標都可以分辨開，這是由于在小快拍情況下主特征空間擁有所有目標信號的信息；而在大快拍情況下，主特征空間信號信息泄漏，從圖11虛線可知，3個目標有2個已很難分辨出。這與“2.2”節及“3.2”節仿真結論一致，進一步驗證了“2.2”節分析的結果。

圖11 水池實驗獨立源主特征空間DOA估計Fig.11 Tank experiment DOA of independent source with main feature space

5 海試驗證

海試數據采用1993-10-26/10-27，SACLANT研究中心在意大利西海岸北厄爾巴島附近進行2天海試的數據，海況等級為2～3級，該實驗使用的接收陣是48元的均勻線陣，采樣頻率1 kHz。參照文獻[15]，在335 Hz處有2個目標，本文對時域信號分段，每1 000數據點做為一個快拍，傅里葉變換后，取335 Hz處的數據構造樣本協方差矩陣。

圖12為使用本文提出方法得到的目標DOA估計，使用了48個快拍(快拍數等于陣元數)，與Qian等獲得的目標方位結果一致，即準確估計出目標的方位。圖13為與圖12使用相同的樣本協方差矩陣，通過MUSIC算法DOA估計所得到的結果，由圖13可以發現其分辨出3個目標，出現了偽目標，無法準確得到目標方位。

圖12 海試主特征空間DOA估計Fig.12 Sea trial DOA with main feature space

圖13 海試MUSIC算法DOA估計Fig.13 Sea trial DOA with MUSIC

6 結 論

本文針對MUSIC算法對相干目標分辨的失效以及在小快拍情況下傳統波達方向估計方法目標分辨能力的下降，通過分析相干目標與獨立目標的樣本協方差矩陣的譜分離特性，提出了基于主特征空間的波達方向估計方法。

(1)該方法在深入分析相干目標與獨立目標信號信息在特征空間分布的基礎上，發現在小快拍情況下，目標信息主要分布在樣本協方差矩陣的最大特征值及所對應的特征向量上，將其定義為主特征空間；利用導向向量與噪聲子空間正交，導向向量張成的空間與信號子空間為同一空間，定義了基于MFS的目標方位估計方法。

(2)通過仿真、水池實驗以及海試數據驗證表明，對于相干目標，不管快拍數大小，本文方法都可以實現目標方位的準確分辨，只是快拍數較大時旁瓣較低，彌補了MUSIC算法對相干目標方位估計性能下降的缺陷，且不需要先驗信息信源個數；對于獨立目標，在小快拍條件下也可以準確分辨目標，改善了在小快拍情況下，樣本協方差矩陣估計誤差導致的傳統DOA估計方法無法準確估計目標，甚至出現偽目標的情況?？傊?，本文所提方法在不需要信源個數等先驗信息情況下，使用小快拍樣本協方差矩陣，就能實現相干目標與獨立目標的估計，將會在UUV(Unmanned Underwater Vehicle)、魚雷等水下目標方位估計的小快拍應用場景中有廣闊的前景。

參 考 文 獻

[ 1 ] 王永良, 陳輝, 彭應寧, 等. 空間譜估計理論與算法[M]. 北京： 清華大學出版社, 2004.

[ 2 ] CARLSON B D. Covariance matrix estimation errors and diagonal loading in adaptive arrays[J]. IEEE Transactions on Aerospace Electronic Systems, 1988, 24(4): 397-401.

[ 3 ] CHEN Y, WIESEL A, ELDAR Y C, et al. Shrinkage algorithms for MMSE covariance estimation[J]. IEEE Transactions on Signal Processing, 2010, 58(10): 5016-

5029.

[ 4 ] PILLAI S U, KWON B H. Forward/backward spatial smoothing techniques for coherent signal identification[J]. IEEE Transactions on Acoustics Speech & Signal Processing, 1989, 37(1): 8-15.

[ 5 ] EVANS J E, SUN D F, JOHNSON J R. Application of advanced signal processing techniques to angle of arrival estimation in ATC navigation and surveillance systems[R]. Lexington MA: MIT Lincoln Laboratory, 1982.

[ 6 ] GERSHMAN A B, BOHME J F. Improved DOA estimation via pseudorandom resampling of spatial spectrum[J]. IEEE Signal Processing Letters, 1997, 4(2): 54 - 57.

[ 7 ] VASYLYSHYN V. Removing the outliers in root-MUSIC via pseudo-noise resampling and conventional beamformer[J]. Signal Processing, 2013, 93(12): 3423-3429.

[ 8 ] QIAN C, HUANG L, SO H C. Improved unitary root-MUSIC for DOA estimation based on pseudo-noise resampling[J]. IEEE Signal Processing Letters, 2014, 21(2): 140-144.

[ 9 ] SHAGHAGHI M, VOROBYOV S A. Subspace leakage analysis and improved DOA estimation with small sample size[J]. Signal Processing IEEE Transactions on, 2015, 63(12): 3251-3265.

[10] BAI Z D, SILVERSTEIN J W. No eigenvalues outside the support of the limiting spectral distribution of large-dimensional sample covariance matrices[J]. Annals of Probability, 1998, 26(1): 316-345.

[11] BAI Z D, SILVERSTEIN J W. Exact separation of eigenvalues of large dimensional sample covariance matrices[J]. Annals of Probability, 1999, 27(3): 1536-1555.

[12] 李華, 白志東, 肖玉山. 大維隨機矩陣的漸進特征[J]. 東北師大學報(自然科學版),2014,46 (4): 1-8.

LI Hua, BAI Zhidong, XIAO Yushan. The asymptotic properties of the large dimension random matrix[J]. Journal of Northeast Normal University (Natural Science), 2014 ,46(4): 1-8.

[14] STERGIOPOULOS S, SULLIVAN E J. Extended towed array processing by an overlap correlator[J]. The Journal of the Acoustical Society of America, 1989, 86(1): 158-171.

[15] CHAVALI V, WAGE K E, BUCK J R. Coprime processing for the Elba Island sonar data set[C]∥ 48th Asilomar Conference on Signals, Systems and Computers. California: IEEE, 2014.