稀疏圓陣的解相干求根MUSIC算法

， ， ，

(空軍預警學院， 湖北武漢 430019)

0 引言

均勻圓陣(Uniform Circular Array, UCA)在雷達、聲納和無線通信等領域應用日益廣泛，這是因為相比于常規的均勻線陣，UCA具有可以同時提供俯仰角和方位角信息、方向圖在任何方位上具有基本相同的波束形狀和圓陣天線更易實現共形等優勢[1-2]。UCA的超分辨DOA(Direction of Arrival)估計算法主要包括陣元域直接處理和波束域處理兩類。在陣元域直接處理時，由于UCA導向矢量不具備范德蒙德結構，因此無法使用root-MUSIC, ESPRIT等算法減少運算量，以及平滑處理、矩陣重構等方法解相干。文獻[3]研究了兩種不同放置方式的圓陣的分辨特性，文獻[4-5]提出了直接利用圓陣的特殊陣列設置實現對相干源的DOA估計，不過這類算法的解相干性能和陣列設置的結構密切相關，而且沒有考慮陣列稀疏的情況。

波束域處理通常采用波束變換的方法得到具有范德蒙德結構的波束域導向矢量。陣列內插變換[6]和相位模式激勵方法[7]是UCA波束變換的常用方法，前者適用于較小方位區域的處理，增加處理方位的范圍會導致誤差增大；后者利用空間離散傅里葉變換，適用于全方位處理，在陣元數目足夠多的情況下，波束域導向矢量可以看作具有共軛對稱和范德蒙德結構。以相位模式激勵方法為基礎，文獻[8]提出了UCA-RB-MUSIC，UCA-ESPRIT算法，文獻[9-10]提出了Unitary root-MUSIC算法，這些算法在波束域進行處理，避免了譜搜索過程，降低了計算量；文獻[11-14]提出了針對相干源的圓陣波束域DOA估計算法，這些算法利用波束域導向矢量的共軛對稱和范德蒙德結構，通過矩陣重構實現了對相干源的解相干。但是，以上算法的適用范圍是陣元數大于2K(K是最大模式值，與圓陣半徑和信號頻率的乘積成正比)，當陣元數減少，即陣列稀疏時會產生很大誤差，導致上述算法失效。為了解決稀疏圓陣的DOA估計問題，文獻[15]定量分析了相位模式激勵類算法的誤差，并且提出了一種改進的波束變換方法，通過循環迭代的方式可以消除誤差的主要部分，從而使波束域導向矢量更加接近理想的共軛對稱和范德蒙德結構。文獻[16]提出的sparse UCA root-MUSIC(SR)算法則是不改變相位模式激勵的變換矩陣，而對波束域的導向矢量進行誤差補償，從而消除誤差的主要部分。文獻[15]中的方法復雜度較高，需要信號方位的先驗信息，而且在陣元數小于K時只能消除部分誤差，導致算法性能急劇下降，因此只適用于陣元數大于K的情況。SR算法理論上適用于任意陣元數的UCA，但是得到的波束域導向矢量不再具有共軛對稱或者范德蒙德結構，無法使用平滑、前后向平均等解相干方法。

本文在SR算法的基礎上，首先改進了波束變換矩陣，對變換矩陣進行相位校正，從而使波束域導向矢量具備了共軛對稱結構，進而提出了稀疏圓陣的解相干求根MUSIC算法(Sparse UCA Decorrelation Root-MUSIC , SDR)，這種算法適用于任意陣元數的UCA，而且進行了前后向平均處理，可以實現對相干信號源的DOA估計，在低信噪比和低快拍數情況下有更好的估計性能。

1 陣列模型

這里只考慮入射信號與陣列在同一個平面的情況，假設均勻圓陣由M個各向同性的陣元組成，圓陣水平放置，半徑為r。以圓心為參考點，則陣列的導向矢量可以表示為

a(θ)=[e-j2πr/λ cos(θ-γ1), e-j2πr/λ cos(θ-γ2),…,e-j2πr/λ cos(θ-γM)]T

(1)

陣元的接收數據可以表示為

x(t) =As(t)+n(t)

(2)

式中，A=[a(θ1),a(θ2),…,a(θN)]為N個信號的陣列流形，s(t)為N×1維的信號矢量，n(t)為M×1維的噪聲矢量，x(t)為M×1維的陣元接收數據矢量。

bkejkθ+Δ

(3)

b(θ)≈FaB(θ)=

[b-Ke-jKθ,b-(K-1)e-j(K-1)θ,…,bKejKθ]T

(4)

式中，F=diag([b-K,b-(K-1),…,bK])，aB(θ)=[e-jKθ,e-j(K-1)θ,…,ejKθ]T。

常規相位模式激勵法以F-1W=F-1[w-K,w-(K-1),…,wK]H作為波束形成器，那么aB(θ)即為波束域導向矢量，具備共軛對稱和范德蒙德結構，可以使用平滑處理及矩陣重構等方法解相干[12-15]。

SR算法以W作為波束形成器，考慮陣元數為奇數且M<2K的陣列稀疏情況，記P=(M-1)/2，取|k|≤P，此時Δ中存在P<|k±qM|

|k±qM|

co(θ)=Hob(θ)

(5)

式中，Ho=[I1?I?Ir]，I為(2P+1)×(2P+1)維的單位矩陣，Il為I的后K-P列，Ir為I的前K-P列。顯而易見，這種方法得到的波束域維數為2P+1，將P<|k±qM|

對于補償后的導向矢量仍然可以進行root-MUSIC運算，避免了計算量較大的譜搜索過程。上面提到的補償矩陣Ho的形式只適用于K<3P+1的情形，當陣元數繼續減少時，文獻[17]給出了補償矩陣相應的形式，在3P+1

2 改進的波束變換和平滑方法

2.1 改進的波束變換方法

SR算法改進了傳統的相位模式激勵方法，適用于陣元稀疏的情況，而且可以采用root-MUSIC方法減少計算量。但是它的缺陷在于不適用于相干信號源的DOA估計，而實際應用中，由于多徑等情況的存在，相干信號源是不可避免的。SR算法不能解相干的原因是它的波束域導向矢量co(θ)不再具有共軛對稱或者范德蒙德結構，因此，首先考慮改進SR算法，得到具有共軛對稱結構的波束域導向矢量。

這里考慮陣元數目為偶數的情況，即M=2P，同樣地，P

He=[0l?I?0r]+[Il?0?Ir]

(6)

式中，0l和0r都為(2P+1)×(K-P)維的零矩陣，I為(2P+1)×(2P+1)維的單位矩陣，Il為I的后K-P+1列，Ir為I的前K-P+1列，0為(2P+1)×(2P-1)維的零矩陣。

可以看出，當陣元數目為偶數時，補償矩陣有別于陣元數目為奇數的情況。下面討論這兩種情況下波束域的特征。為了直觀地說明問題，這里取K=5，P=3。首先得出陣元數目M=6和M=7兩種情況下的補償矩陣。

(7)

相應的波束域導向矢量為

(8)

觀察以上兩種情況下的波束域導向矢量可以看出：

1)ce(θ)的首末兩項是重復的，其中一個是冗余的，co(θ)則無冗余，也就是ce(θ)的有效維數為M。這是由于補償后的波束域是2P+1維的，當陣元數為2P時，少于波束域維數。

3)ce(θ)的元素中相加項的下標同為奇數或者同為偶數，而co(θ)的元素中相加項的下標為一奇數一偶數。

在上述觀察結果的基礎上，為了使波束域導向矢量具有共軛對稱結構，再考慮bk=jkJk

T=diag([(-1)P,(-1)P-1,…,(-1)0,1,…,1])

(9)

綜上所述，新的波束形成器是TW，對應的波束域導向矢量d(θ)滿足共軛對稱結構，表示為

d(θ)=Tce(θ)=THeFaB(θ)

(10)

2.2 波束域前后向平滑處理

通過前一小節所述的改進波束變換方法，可以得到具有共軛對稱結構的波束域導向矢量，從而可以利用前后向平均處理的方法，實現對相干源的DOA估計。需要說明的是，這種方法針對的是時域相位不一致的相干信號，即信號的相位差非零且非π的情況，實際中的相干源基本上都滿足這一要求。

首先，介紹波束域前后向平均處理的方法，采用前一節提出的波束變換方法后，可以得到波束域數據：

y(t)=TWx(t)=TWAs(t)+TWn(t)=

Ds(t)+TWn(t)

(11)

式中，D=[d(θ1)d(θ2)…d(θN)]。相應的波束域數據協方差矩陣可以表示為

R= E [y(t)yH(t)]

(12)

然后，可以進行前后向平均處理：

(13)

式中，J為(2P+1)×(2P+1)維的反對角線單位矩陣。

(14)

令z=ejθ，可得出關于z的2P階方程，通過對方程求根即可得出相干信號源的DOA。

R=DE [s(t)s*(t)]DH+TWE [n(t)n*(t)]WHTH=

(15)

(16)

ρ*ρT)，下面對它進一步求解：

(17)

式中，P=[ρρ*]，P又可以表示為

(18)

很明顯可以得到，在|angle(ρ(1))-angle(ρ(2))|

3 仿真分析

實驗中，圓陣半徑r=λ，

最大模式值為K=

圖1 r=λ情況下|bk|與k的關系圖

實驗1：不同陣元數時的解相干性能

實驗仿真比較本文算法與文獻[13]中MODETOEP算法及MODEFBSS算法在不同陣元數情況下的解相干性能。實驗中陣元數目分別為16，10和6，其中，陣元數16是陣列不稀疏的情況，陣元數10和6分別是陣元數目大于和小于最大模式K的陣列稀疏的情況， 3種陣元數目分別對應陣元間距為0.39λ，0.63λ和1.05λ。兩個相干信號入射方位分別為-20°和20°，衰落系數分別為1，0.9ej1.3。其中，MODETOEP算法和MODEFBSS算法需要采用傳統的相位模式激勵法，波束域維數為15，波束域導向矢量采用aB(θ)，MODEFBSS算法子陣陣元數為10；本文算法波束域維數分別為17，11和7，波束域導向矢量采用d(θ)。實驗結果由200次蒙特卡洛仿真得出，快拍次數為100。圖2是16元陣3種算法的解相干性能，圖3是10元陣和6元陣使用本文算法的解相干性能。

當陣元數為16時，3種算法均有效，但是MODETOEP算法估計偏差始終為11°左右，估計精度較差，這是因為MODETOEP等矩陣重構類算法適用范圍是ρ為實數矢量，當ρ為復數矢量時，會引起估計角度的偏移。隨著信噪比增大，本文算法和MODEFBSS算法估計偏差逐漸減小。在信噪比小于0 dB時，很明顯本文算法估計精度好于MODEFBSS算法，這是由于本文算法考慮了|k+qM|<11的所有誤差，而MODEFBSS算法需要忽略9<|k±qM|<11的誤差。

當陣元數為6或者10，即陣列稀疏時，MODETOEP算法和MODEFBSS算法失效，無法正確估計出信號源方位，因為這兩種算法要求波束域導向矢量滿足范德蒙德結構，然而陣列稀疏時相位模式激勵有較大誤差，波束域導向矢量并不滿足范德蒙德結構。本文算法由于考慮了模式激勵的誤差，在陣列稀疏時仍有很好估計性能，當信噪比大于-5 dB時，估計偏差小于1°，估計精度好于16元陣的MODEFBSS算法。

實驗2：稀疏陣列對非相關信號源估計性能

實驗仿真對于陣元數為6的稀疏UCA，比較不同快拍數和不同信噪比時本文算法與SR算法對于信號的估計性能。實驗中，采用兩個不相關信號，入射角度分別為-10°和10°，比較信噪比的影響時快拍數為50，比較快拍數影響時，信噪比為-3 dB。當估計方位與目標真實方位相差小于0.1倍的目標方位間距時，認為估計成功。實驗結果由200次蒙特卡洛仿真得出。

從圖4中可以看出，在快拍數大于150或者信噪比大于0 dB時，本文算法與SR算法估計性能基本一致，但是當快拍數小于150或者信噪比小于0 dB時，本文算法成功概率和RMSE都要優于SR算法。這是因為本文算法通過改進波束變換矩陣，增加了平滑處理，具有平均的意義。

(a)成功概率與信噪比的關系

(b)估計偏差與信噪比的關系圖2 16元陣3種算法解相干性能

(a)成功概率與信噪比的關系

(b)估計偏差與信噪比的關系圖3 10陣和6元陣使用本文算法的解相干性能

(a)成功概率與信噪比的關系

(b)RMSE與信噪比的關系

(c)成功概率與快拍數的關系

(d)RMSE與快拍數的關系圖4 稀疏陣列的估計性能

實驗3：稀疏陣列的解相干能力

實驗仿真陣元數為6和8的兩種稀疏圓陣，其中針對6元陣仿真了1組和2組相干源的情況，針對8元陣仿真了2組和3組相干源的情況。其中各組相干源入射角度分別是(10°,-10°)，(50°,-50°)和(100°,-100°)。實驗快拍數為100，信噪比為0 dB。圖5和圖6為利用本文算法求解所得根的分布，其中橫軸表示根的實部，縱軸表示根的虛部，與橫軸正方向的夾角即為估計的入射信號方位角。

(a)1組相干源

(b)2組相干源圖5 6元陣根分布圖

(a)1組相干源

(b)2組相干源

(c)3組相干源圖6 8元陣根分布圖

從圖5和圖6可以看出，6元陣和8元陣分別可以估計出4個相干源和6個相干源的入射方位。理論上，對于M元陣，在波束域有M個有效的自由度，可以估計M-1個相干源，然而實際中由于噪聲等因素的影響，在對M-1個相干源估計時會出現多個虛假方位，導致估計失敗，最多只能估計出M-2個相干源。

4 結束語

本文首先在SR算法的基礎上，對波束變換過程中的誤差進行補償，得到了波束域導向矢量的正確表達形式，相比傳統的相位模式激勵方法要求陣元數大于2K，該方法適用于陣元數小于2K的稀疏圓陣。然后，通過分析新的波束域導向矢量的結構特征，提出對波束域進行相位校正的方法，從而實現了波束域導向矢量的共軛對稱結構。最后，利用這種共軛對稱結構進行前后向平均，實現了對相干信號源的解相干。

相比其他的圓陣波束域解相干算法，該方法考慮了波束變換的誤差，從而不再需要陣元數目大于2K，在陣列稀疏的情況下仍然具有較高的估計精度。與陣元域處理算法相比，本文算法可以采用root-MUSIC算法，避免了譜搜索過程，減少了計算量。

[1] JACKSON B R, RAJAN S, LIAO B J, et al. Direction of Arrival Estimation Using Directive Antennas in Uniform Circular Arrays[J]. IEEE Trans on Antennas and Propagation, 2015, 63(2):736-747.

[2] WANG P, LI Y H, VUCETIC B. Millimeter Wave Communications with Symmetric Uniform Circular Antennas Arrays[J]. IEEE Communications Letters, 2014, 18(8):1307-1310.

[3] 吳垚,陳輝. 兩種均勻圓陣特性對比分析研究[J]. 雷達科學與技術, 2011, 9(1):62-66.

WU Yao, CHEN Hui. Comparison and Analysis of Characteristics of Two Different Uniform Circular Arrays[J]. Radar Science and Technology, 2011, 9(1):62-66.(in Chinese)

[4] WU Y, CHEN H, CHEN Y C. A Method of 2-D DOA Estimation of Coherent Signals Based on Uniform Circular Array via Spatial Smoothing[C]∥CIE International Radar Conference, Chengdu:IEEE, 2011:312-314.

[5] CHEN H, ZHAO M. A Modified Method of the Wideband Coherent Signals DOA Estimation on the Uniform Circular Arrays[C]∥2nd International Conference on Electronics, Communications and Control, Zhoushan:IEEE, 2012:3166-3169.

[6] 蔣柏峰,呂曉德. 一種基于導向矢量變換的DOA估計預處理方法[J]. 電子與信息學報, 2012, 34(7):1552-1557.

[7] DAVIES D. A Transformation Between the Phasing Techniques Required for Linear and Circular Aerial Arrays[J]. IET Proceedings of the Institution Electrical Engineers, 1965, 112(11):2041-2045.

[8] MATHEWS C P, ZOLTOWSKI M D. Eigenstructure Techniques for 2-D Angle Estimation with Uniform Circular Arrays[J]. IEEE Trans on Signal Processing, 1994, 42(9):2395-2407.

[9] BELLONI F, KOIVUNEN V. Unitary Root-MUSIC Technique for Uniform Circular Array[C]∥IEEE International Symposium on Signal Processing and Information Technology, Darmstadt, Germany:IEEE, 2003:451-454.

[10] 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.

[11] HAN Yong, GAO Bingxia, DONG Wenyi, et al. An Improved Mode-MUSIC Algorithm for DOA Estimation of Coherent Sources Based on Hybrid Array[C]∥International Conference on Signal Processing, Beijing:IEEE, 2014:358-362.

[12] GAO Shuyan, CHEN Hui, WANG Yongliang, et al. A Novel Algorithm for Estimating DOA of Coherent Signals on Uniform Circular Array[C]∥CIE International Conference on Radar, Shanghai:IEEE, 2006:1-4.

[13] WANG Ju, GU Xiaopeng, LIU Yan. The Coherent DOA Estimation of Uniform Circular Array Based on SVD Algorithm[C]∥International Conference on Electrical and Control Engineering, Yichang:IEEE, 2010:6147-6151.

[14] ZHAO Hui, WANG Hongjin. MODE-SUMWE Algorithm for Estimating DOA of Coherent Signals on Uniform Circular Array[C]∥International Conference on Networks Security, Wireless Communications and Trusted Computing, Wuhan:IEEE, 2010:318-320.

[15] BELLOIN F, KOIVUNNEN V. Beamspace Transform for UCA:Error Analysis and Bias Reduction[J]. IEEE Trans on Signal Processing, 2006, 54(8):3078-3089.

[16] GOOSSENS R, ROGIER H, WERBROUCK S. UCA Root-MUSIC with Sparse Uniform Circular Arrays[J]. IEEE Trans on Signal Processing, 2008, 56(8):4095-4099.

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