邱 偉，包長(zhǎng)春

(國(guó)防科技大學(xué) 氣象海洋學(xué)院， 湖南 長(zhǎng)沙 410073)

信源到達(dá)角(Direction-Of-Amival, DOA)估計(jì)是陣列信號(hào)處理的重要研究方向之一，在雷達(dá)、聲吶、通信等領(lǐng)域都有廣泛的應(yīng)用[1-2]。傳統(tǒng)的DOA估計(jì)方法主要包括常規(guī)波束形成算法(Conventional BeamForming, CBF)、Capon算法、子空間類算法、最大似然估計(jì)算法等[2]，上述算法適用于均勻線性陣列，當(dāng)陣列出現(xiàn)壞道時(shí)，也就是部分陣元失效時(shí)，估計(jì)性能會(huì)嚴(yán)重下降。

利用信源在整個(gè)空域中的分布具有稀疏性的特點(diǎn)，近年來(lái)，基于壓縮感知的DOA估計(jì)方法層出不窮。該理論認(rèn)為，只要信號(hào)在某種變換域下具有稀疏性，那么就可以通過(guò)求解由該稀疏性約束的最優(yōu)化問(wèn)題，以遠(yuǎn)少于奈奎斯特采樣的采樣數(shù)據(jù)中高概率精確重構(gòu)原信號(hào)。將壓縮感知理論應(yīng)用到DOA估計(jì)問(wèn)題中，可以提高陣列存在壞道條件下的DOA估計(jì)性能，并且相對(duì)于傳統(tǒng)方法具有更好的角度分辨能力。文獻(xiàn)[3]提出了一種基于改進(jìn)正交匹配追蹤算法的DOA估計(jì)方法，利用信源功率譜主副瓣之間的差異設(shè)置合理的閾值來(lái)確定算法迭代的終止條件，在信源數(shù)未知的條件下仍能獲得高精度的估計(jì)結(jié)果。文獻(xiàn)[4]采用一種加權(quán)迭代最小范數(shù)(FOCal Underdetermined System Solver, FOCUSS)方法來(lái)求解稀疏約束下的DOA估計(jì)問(wèn)題，取得了較好的結(jié)果。針對(duì)多快拍數(shù)據(jù)，文獻(xiàn)[5]進(jìn)一步提出了多次快拍FOCUSS算法，通過(guò)求解多測(cè)量稀疏重構(gòu)問(wèn)題，實(shí)現(xiàn)DOA估計(jì)。文獻(xiàn)[6]提出了L1-SVD(L1-singular value decomposition)方法，利用奇異值分解降低數(shù)據(jù)量，并基于L1范數(shù)進(jìn)行稀疏重構(gòu)獲得信源的方位估計(jì)結(jié)果。文獻(xiàn)[7]針對(duì)在非均勻噪聲下DOA估計(jì)精度差以及分辨率低的問(wèn)題，利用矩陣補(bǔ)全理論引入彈性正則化因子重構(gòu)無(wú)噪條件下的協(xié)方差矩陣，然后利用稀疏重構(gòu)加權(quán)L1范數(shù)實(shí)現(xiàn)DOA參數(shù)估計(jì)，該方法能顯著抑制非均勻噪聲，并具有較好的DOA估計(jì)性能。針對(duì)基于L1范數(shù)優(yōu)化的DOA估計(jì)方法存在正則化參數(shù)選取困難、計(jì)算復(fù)雜度高等問(wèn)題，文獻(xiàn)[8-9]利用稀疏貝葉斯學(xué)習(xí)提出高效的DOA估計(jì)方法，通過(guò)優(yōu)化稀疏貝葉斯學(xué)習(xí)的基消除機(jī)制提高算法效率。文獻(xiàn)[10]針對(duì)有源欺騙干擾環(huán)境下的小樣本DOA估計(jì)問(wèn)題，提出自適應(yīng)極化濾波和聯(lián)合塊稀疏貝葉斯學(xué)習(xí)算法的DOA估計(jì)方法，在有效抑制干擾的同時(shí)能夠獲得高分辨率和高精度的DOA估計(jì)結(jié)果。

不同于上述直接基于陣列接收數(shù)據(jù)稀疏性進(jìn)行DOA估計(jì)的思想，文獻(xiàn)[11]提出了一種基于陣列協(xié)方差矩陣稀疏表示的DOA估計(jì)方法L1-SRACV(L1-sparse representation of array covariance vectors)，該方法利用接收陣的協(xié)方差矩陣估計(jì)誤差滿足漸近正態(tài)分布的特性來(lái)實(shí)現(xiàn)DOA估計(jì)，具有更強(qiáng)的抗噪聲性能，且無(wú)須確定正則化參數(shù)，但該方法在快拍數(shù)較小時(shí)由于協(xié)方差矩陣估計(jì)誤差較大，導(dǎo)致性能嚴(yán)重下降。針對(duì)該問(wèn)題，文獻(xiàn)[12]提出一種改進(jìn)方法，利用快速極大似然算法(Fast Maximum Likelihood, FML)穩(wěn)定協(xié)方差矩陣的小特征值，提升了算法在快拍數(shù)較少時(shí)的穩(wěn)定性。雖然基于協(xié)方差矩陣稀疏性的DOA估計(jì)方法通常具有更好的抗噪聲性能，但是傳統(tǒng)方法在求解過(guò)程中需要對(duì)矩陣進(jìn)行矢量化操作，并將二維稀疏重構(gòu)問(wèn)題轉(zhuǎn)化為常規(guī)稀疏性約束下的一維矢量稀疏優(yōu)化問(wèn)題，計(jì)算量大大增加，一定程度上限制了其實(shí)用性。

針對(duì)該問(wèn)題，本文提出了一種基于協(xié)方差矩陣稀疏性的高效DOA估計(jì)方法。該方法直接利用協(xié)方差矩陣的稀疏性，構(gòu)造二維稀疏表示模型，并通過(guò)直接求解二維稀疏性約束的優(yōu)化問(wèn)題實(shí)現(xiàn)DOA估計(jì)，避免了矩陣矢量化操作，從而大大提高了算法的效率，同時(shí)具有較高的DOA估計(jì)性能。

1 信號(hào)模型

考慮均勻線性陣列構(gòu)型，設(shè)其陣元個(gè)數(shù)為M，定義方位角θ為來(lái)波方向與x軸正方向的夾角，因此有0°<θ<180°。假設(shè)在遠(yuǎn)場(chǎng)條件下，有K個(gè)窄帶信源入射到該陣列，信源之間相互獨(dú)立。則第k個(gè)信源、來(lái)波方向?yàn)棣萲，入射到該均勻線性陣列(間距為d)的陣列響應(yīng)矢量為：

(1)

定義導(dǎo)向矢量矩陣為：

(2)

此時(shí)，陣列第m個(gè)陣元的接收信號(hào)可表示為：

(3)

(4)

(5)

(6)

此時(shí)，陣列的接收信號(hào)可以寫(xiě)成以下的矢量模型：

x(t)=AS(t)+N(t)

(7)

這里，導(dǎo)向矢量矩陣A與陣列的形狀以及信源的方位有關(guān)，而在一般的實(shí)際應(yīng)用中，陣列的形狀一般是固定的，所以矩陣A中任一列總是和某個(gè)信源的方位緊密聯(lián)系。

在實(shí)際處理中，陣列接收到的數(shù)據(jù)是在有限時(shí)間范圍內(nèi)的有限次快拍數(shù)據(jù)。在這段時(shí)間內(nèi)假定信源的方向不發(fā)生變化，信源信號(hào)的包絡(luò)雖然隨時(shí)間變化，但通常認(rèn)為它是一個(gè)平穩(wěn)隨機(jī)的過(guò)程，其統(tǒng)計(jì)特性不隨時(shí)間變化，這樣就可以定義陣列輸出信號(hào)x(t)的協(xié)方差矩陣為：

R=E[x(t)xH(t)]=ARsAH+σ2IM

(8)

2 本文所提DOA估計(jì)方法

(9)

(10)

首先定義高斯函數(shù)

(11)

(12)

圖1 不同ρ值下的高斯函數(shù)結(jié)果Fig.1 Values of Gaussian function under different values of ρ

Step1:初始化。

2)選擇遞減的ρ序列：ρ=[ρ0ρ1…ρr]，這里ρr=ηρr-1(r=1,2,…)，η為遞減因子。

Step2:迭代求解。

1)令ρ=ρr-1。

Step3:令r=r+1。重復(fù)Step2中的步驟2)和步驟3)，直至ρr≤ρt，算法終止，其中ρt為設(shè)定的某一閾值。

從上述算法描述可以看出，在求解重優(yōu)化問(wèn)題時(shí)，本文算法由于直接對(duì)二維協(xié)方差矩陣進(jìn)行處理，避免了傳統(tǒng)方法的矢量化操作，從而大大提升了算法的求解效率，增強(qiáng)了算法的實(shí)用性。

3 算法討論

3.1 算法計(jì)算復(fù)雜度分析

3.2 算法參數(shù)選取

4 仿真實(shí)驗(yàn)及結(jié)果分析

為了驗(yàn)證所提方法的有效性，采用仿真數(shù)據(jù)對(duì)算法性能進(jìn)行驗(yàn)證分析，并與L1-SRACV算法進(jìn)行比較。實(shí)驗(yàn)中，采用16陣元的均勻線性陣列，即M=16，陣元間距為半波長(zhǎng)，空域角度范圍設(shè)定為0°～180°，間隔為1°，因此，Q=181。

實(shí)驗(yàn)一：假設(shè)空域中存在兩個(gè)信源，來(lái)波方向?yàn)?0°和70°。陣列每個(gè)陣元的接收噪聲均為相互獨(dú)立的零均值復(fù)高斯白噪聲，信噪比設(shè)定為10 dB，并與信源互不相關(guān)，快拍數(shù)T=500。圖2(a)和圖2(b)分別給出了非相關(guān)信源和相關(guān)信源(相關(guān)系數(shù)為0.9)由不同方法得到的空域功率譜，對(duì)比可知，相對(duì)于常規(guī)Capon方法，稀疏重構(gòu)類方法包括本文方法、L1-SRACV方法和L1-SVD方法，其功率譜譜峰更尖銳、旁瓣更低，顯示了稀疏重構(gòu)類算法的優(yōu)越性，因此在后續(xù)的DOA估計(jì)結(jié)果分析中只比較三種稀疏重構(gòu)算法的DOA估計(jì)性能。另外，本文方法對(duì)于相關(guān)信源和非相關(guān)信源均能取得較好的估計(jì)效果，無(wú)須額外的去相關(guān)處理，且比L1-SRACV和L1-SVD方法結(jié)果的旁瓣更低，該結(jié)果表明本文方法對(duì)信源的相關(guān)性不敏感。

(a) 非相關(guān)信源DOA估計(jì)結(jié)果

此外，本實(shí)驗(yàn)還比較了本文方法、L1-SRACV方法和L1-SVD方法運(yùn)行所需的時(shí)間。算法在配置相同的計(jì)算機(jī)上運(yùn)行，計(jì)算機(jī)型號(hào)為聯(lián)想Thinkpad X230系列，CPU為Intel 酷睿i5-3210 2.50 GHz，內(nèi)存為4 GB，操作系統(tǒng)為Windows7 32位，三種方法的運(yùn)行時(shí)間分別為2.392 4 s、38.111 4 s和 4.641 8 s，可以看出，本文算法由于無(wú)須矩陣矢量化操作，算法運(yùn)行時(shí)間遠(yuǎn)低于L1-SRACV算法，同時(shí)也低于L1-SVD算法，顯示了本文算法的高效性。

實(shí)驗(yàn)二：考慮非相關(guān)信源，來(lái)波方向與實(shí)驗(yàn)一相同，快拍數(shù)設(shè)定為300，加入不同噪聲水平的復(fù)高斯白噪聲，每一信噪比下蒙特卡洛仿真次數(shù)設(shè)置為100，用以考察不同信噪比下不同方法的DOA估計(jì)成功率。圖3(a)和圖3(b)給出了信噪比分別為-10 dB和-5 dB時(shí)本文方法、L1-SRACV算法和L1-SVD算法的DOA估計(jì)結(jié)果，圖4給出了三種方法的DOA估計(jì)成功率。從圖3和圖4的結(jié)果可以看出，本文方法的DOA估計(jì)成功率相對(duì)于L1-SRACV和L1-SVD方法有較大提升，特別是在低信噪比條件下，成功率要明顯優(yōu)于L1-SRACV和L1-SVD方法，說(shuō)明了本文方法的穩(wěn)健性。

(a) 信噪比為-10 dB時(shí)的DOA估計(jì)結(jié)果

圖4 不同信噪比下DOA估計(jì)成功率Fig.4 DOA estimation successful rate as a function of SNR

實(shí)驗(yàn)三：本實(shí)驗(yàn)中，考察不同快拍數(shù)條件下的DOA估計(jì)成功率。仍然考慮兩個(gè)非相關(guān)信源，信源來(lái)波方向與實(shí)驗(yàn)一相同，信噪比固定為10 dB，圖5(a)和圖5(b)給出了快拍數(shù)分別為20和50時(shí)本文方法、L1-SRACV算法和L1-SVD算法得到的空域功率譜。從中可以看出，在低快拍數(shù)條件下，L1-SRACV算法和L1-SVD算法無(wú)法實(shí)現(xiàn)DOA估計(jì)，而本文方法仍能獲得高質(zhì)量的DOA估計(jì)結(jié)果。下面在不同快拍數(shù)條件下進(jìn)行100次蒙特卡洛仿真，計(jì)算DOA估計(jì)成功率。圖6給出了本文方法、L1-SRACV方法和L1-SVD方法在不同快拍數(shù)條件下的DOA估計(jì)成功率結(jié)果，從圖中可以看出，在高快拍數(shù)條件下，三種方法均能取得非常高的DOA估計(jì)成功率；在低快拍數(shù)條件下，L1-SRACV方法的DOA估計(jì)成功率明顯降低，L1-SVD方法性能優(yōu)于L1-SRACV，但是，在快拍數(shù)極低時(shí)估計(jì)成功率不高，而本文方法仍能保持較高的成功率，這顯示了本文方法在低快拍數(shù)條件下的優(yōu)越性。

(a) 快拍數(shù)為20時(shí)的DOA估計(jì)結(jié)果

圖6 不同快拍數(shù)下的DOA估計(jì)成功率Fig.6 DOA estimation successful rate under different snapshots

實(shí)驗(yàn)四：本實(shí)驗(yàn)中考察不同角度間隔下各種方法的DOA估計(jì)成功率。假設(shè)兩個(gè)非相關(guān)信源，一個(gè)信源的角度θ1固定為50°，另一個(gè)信源的角度θ2=θ1+Δθ，其中角度間隔Δθ的范圍為2°～20°，步長(zhǎng)為2°。信噪比固定為10 dB，快拍數(shù)設(shè)定為300。每一角度間隔步長(zhǎng)條件下進(jìn)行100次蒙特卡洛仿真。圖7給出了不同角度間隔條件下三種方法的DOA估計(jì)成功率。從結(jié)果可以看出，在角度間隔較大時(shí)，三種方法均能獲得較高的DOA估計(jì)成功率，而在角度間隔較小時(shí)，本文方法的DOA估計(jì)成功率低于L1-SRACV方法和L1-SVD方法。這表明本文方法在鄰近目標(biāo)分辨能力方面要弱于其他兩種方法，需要在下一步工作中進(jìn)行改進(jìn)。

圖7 不同角度間隔下的DOA估計(jì)成功率Fig.7 DOA estimation successful rate under different angular interval

實(shí)驗(yàn)五：本實(shí)驗(yàn)中，比較不同方法對(duì)DOA估計(jì)精度方面的性能，采用的性能指標(biāo)為均方根誤差(Root Mean Square Error, RMSE)，定義為：

(13)

(a) 不同信噪比下的DOA估計(jì)RMSE

實(shí)驗(yàn)六：本實(shí)驗(yàn)中，考察不同稀疏陣元數(shù)條件下的DOA估計(jì)成功率。仍然考慮兩個(gè)非相關(guān)信源，來(lái)波方向保持不變，快拍數(shù)設(shè)置為300，信噪比固定為10 dB。隨機(jī)選取若干個(gè)陣元的接收數(shù)據(jù)進(jìn)行DOA估計(jì)，各稀疏陣元條件下重復(fù)實(shí)驗(yàn)次數(shù)為100。圖9給出了不同稀疏陣元數(shù)條件下三種方法的DOA估計(jì)成功率。從結(jié)果可以看出，隨著稀疏陣元數(shù)的增加，兩種方法的DOA估計(jì)成功率均隨之提高，當(dāng)隨機(jī)選取的陣元數(shù)大于8時(shí)，三種方法的DOA估計(jì)成功率均在90%以上。該結(jié)果表明，基于信源稀疏特性的DOA估計(jì)方法能利用較少的有效陣元實(shí)現(xiàn)信源的測(cè)向，對(duì)于解決存在壞道情況下的信源DOA估計(jì)問(wèn)題意義重大。

圖9 不同稀疏陣元數(shù)下的DOA估計(jì)成功率Fig.9 DOA estimation successful rate under different numbers of sparse array sensors

實(shí)驗(yàn)七：本實(shí)驗(yàn)中，考慮陣列存在陣元位置偏差時(shí)的DOA估計(jì)性能。第一種情況為各陣元偏離程度滿足-0.05～0.05倍半波長(zhǎng)的均勻分布，第二種情況為各陣元位置偏離程度滿足-0.5～0.5倍半波長(zhǎng)的均勻分布，意味著第二種情況下的陣元位置誤差大于第一種情況。信噪比設(shè)定為10 dB，快拍數(shù)設(shè)定為300，兩種情況下的空域功率譜如圖10所示。從圖中可以看出，在第一種情況下，陣元位置偏差對(duì)本文方法和L1-SRACV方法的DOA估計(jì)性能影響較小，而第二種情況下，陣元位置偏移量增大，50°方向的信源功率譜峰值下降，且旁瓣水平升高，對(duì)于信源的檢測(cè)會(huì)有較大影響，特別是對(duì)于L1-SRACV算法，在陣元配置偏移量較大時(shí)，出現(xiàn)多個(gè)較高幅度水平的旁瓣，會(huì)導(dǎo)致信源的檢測(cè)概率下降。

(a) 第一種陣元位置偏離情況下DOA估計(jì)結(jié)果

保持仿真其他條件不變，加入不同功率水平的高斯白噪聲，每一信噪比下進(jìn)行100次蒙特卡洛仿真，三種方法在陣元位置存在第二種偏差情況下的成功率如圖11所示。從圖中可以看出，在低信噪比條件下，當(dāng)陣元位置偏差較大時(shí)，三種方法的DOA估計(jì)成功率相對(duì)于無(wú)偏差時(shí)下降較多，這表明陣列存在誤差時(shí)，DOA估計(jì)性能下降，需要進(jìn)一步對(duì)陣元位置進(jìn)行矯正。

圖11 存在第二種陣元位置誤差時(shí)不同信噪比下的DOA估計(jì)成功率Fig.11 DOA estimation successful rate in the presence of second type of sensor position error

實(shí)驗(yàn)八：最后，給出了不同方法對(duì)于多個(gè)信源的DOA估計(jì)結(jié)果。假設(shè)空間中有10個(gè)非相關(guān)窄帶遠(yuǎn)場(chǎng)信源，其DOA范圍為50°～140°，間隔為10°，信噪比設(shè)定為10 dB，快拍數(shù)設(shè)定為300。圖12給出了10個(gè)信源的DOA估計(jì)空域功率譜，從圖中可以看出，本文方法、L1-SRACV方法以及L1-SVD方法均能對(duì)多個(gè)信源實(shí)現(xiàn)準(zhǔn)確的DOA估計(jì)。

圖12 10個(gè)信源DOA估計(jì)結(jié)果Fig.12 DOA estimation result of 10 uncorrelated sources

5 結(jié)論

基于協(xié)方差矩陣稀疏特性的DOA估計(jì)方法具有良好的噪聲抑制性能，但是經(jīng)典的L1-SRACV算法計(jì)算復(fù)雜度高，并且在低快拍數(shù)條件下的估計(jì)性能嚴(yán)重下降。本文針對(duì)該問(wèn)題提出了一種直接基于協(xié)方差矩陣二維稀疏特性的DOA估計(jì)方法，避免了矩陣矢量化操作，從而大大提高了算法的運(yùn)行效率。同時(shí)，仿真實(shí)驗(yàn)結(jié)果表明該算法在低信噪比、低快拍數(shù)條件下仍然能取得良好的估計(jì)性能，顯示了本文方法相對(duì)于傳統(tǒng)算法的優(yōu)越性。下一步工作將對(duì)本文算法進(jìn)行改進(jìn)，提高其對(duì)于鄰近目標(biāo)的分辨能力以及陣列存在誤差下的DOA估計(jì)性能。