基于M-FIPM的無網格DOA估計算法

陳 濤, 史 林, 申夢雨

(哈爾濱工程大學信息與通信工程學院, 黑龍江 哈爾濱 150001)

0 引 言

波達方向(direction of arrival, DOA)估計是陣列信號處理的一個重要分支,在雷達以及通信等領域有著廣泛的應用[1-2]。傳統的子空間類DOA估計算法,無法直接處理相關信源,并且對信噪比以及快拍數的魯棒性較差[3-4]。而由于壓縮感知理論的完善與快速發展,以l1-svd算法[5]為代表的稀疏類DOA估計算法,近年來逐漸受到國內外眾多學者的廣泛關注與研究[6-8]。

稀疏類DOA估計算法依據空間角度域的離散化處理,將入射信源DOA參數的估計問題轉變為一個稀疏信號重構問題,并利用凸優化等數學手段獲得最終的參數估計結果,在低信噪比以及少快拍的情況下依然可以保持良好的估計性能,并且能夠直接處理相干信源,相比于子空間類算法具有顯著的優勢。但是,稀疏類DOA估計算法內在的網格失配問題會引起不可忽視的算法模型誤差,嚴重影響稀疏DOA估計算法的性能[9-10]。

離網稀疏貝葉斯學習(off-grid sparse Bayestan learning, OGSBL)算法[11-14]以及網格細化[8]技術,分別利用網格點處一階泰勒級數展開與局部細密網格,來解決網格失配問題。雖然上述兩類算法能夠在一定程度上提高網格失配情況下稀疏類DOA估計算法的性能,但依然沒有擺脫空間角度域離散化處理,不能從根本上解決網格失配問題。

直到原子范數概念[15-17]的提出,才為徹底解決網格失配問題提供了一條行之有效的道路,進而形成了一類基于原子范數最小化[18-21](atomic norm minimization, ANM)的無網格DOA估計算法,簡記為ANM-DOA算法。

ANM-DOA算法首先在連續的角度空間上建立以導向矢量為原子的無限維原子集合,然后以天線陣列接收單快拍數據向量的原子范數來表征該向量的稀疏性,最后建立一個以原子范數與擬合誤差的加權和為目標函數的凸優化問題[22-24]。依據原子范數與半正定規劃(semi-definite programming, SDP)問題的等價性,ANM-DOA算法成功地將入射信源DOA參數的估計問題轉化為一個SDP問題,并最終利用Toeplitz矩陣的Vandermonde分解,獲得DOA參數的估計結果。

由此看見,ANM-DOA算法不再依賴于空間角度域的離散化處理,進而從根本上解決了網格失配問題。同時,ANM-DOA算法的核心便在于SDP問題的求解,這同時也是該算法中復雜度最高的運算?，F有的ANM-DOA算法利用內點法[25](interior point method, IPM)或交替方向乘子算法[26-27](alternating direction multiplier method, ADMM)來求解SDP問題。內點法每次迭代過程的復雜度為O(M6),而ADMM算法每次迭代過程的復雜度僅為O(M3)。不過,ADMM算法的收斂速度慢,達到收斂時所需要的迭代次數要遠高于內點法。

所以,基于上述兩種算法的SDP問題求解過程均不理想,如何能夠快速求解SDP問題已然成為當前ANM-DOA算法領域的一種重要研究方向。而快速內點法[28](fast IPM, FIPM)的提出,將SDP問題的維度由O(M2)將低為O(M),從而每次迭代過程的復雜度變為O(M3)。這樣,FIPM不僅保留了內點法收斂速度快的特點,同時又綜合了ADMM算法每次迭代過程算法復雜度低的優勢,所以在相同的SDP問題維度,即陣元數目相同時,FIPM算法的求解速度最快。但是,目前的FIPM僅僅適用于單快拍情況下ANM-DOA算法中的SDP問題。

針對FIPM無法處理多快拍數據的問題,本文提出了一種基于多快拍FIPM(multiple snapshots FIPM, M-FIPM)的ANM-DOA算法。該算法首先對天線陣列接收多快拍數據的協方差矩陣進行特征值分解,然后依據以特征值為系數的特征向量加權和建立滿足FIPM算法模型的SDP問題,最終再利用FIPM算法與Toeplitz矩陣的Vandermonde分解[29-30]獲得DOA參數的估計結果。

在估計精度方面,所提出的基于M-FIPM的ANM-DOA算法,不受網格失配問題的制約,所以與稀疏類DOA估計算法相比,估計精度有顯著的提升。另外,在運行時間方面,基于M-FIPM的ANM-DOA算法每次迭代過程的復雜度為O(M3),達到收斂時所需的迭代次數也要遠小于ADMM算法,所以具有較高的運算效率。

1 信號模型

考慮空間中K個遠場窄帶信號,入射到具有M陣元的均勻線陣上,陣元間距為λ/2。信源角度θk∈[-π/2,π/2],k=1,2,…,K。

陣列天線接收信號的多快拍數據為

Y=AS+N

(1)

式中：S=[s1,s2,…,sK]T∈CK×T,為信號矩陣,T為快拍數,sk=[sk(1),sk(2),…,sk(T)]表示第k個信號。N=[n1,n2,…nm,nm+1,…,nM]T∈CM×T為加性零均值高斯白噪聲矩陣,nm∈CT為第m個陣元的噪聲數據,m=1,2,…,M。

A=[a(θ1),a(θ2),…,a(θk),a(θk+1),…,a(θK)]∈CM×K為信號的導向矢量矩陣,具體地,a(θk)∈CM可以被表示為

a(θk)=[1,e-jπsin θk,…,e-jπ(M-1)sin θk]T

(2)

式中： j為虛數單位。

因此,當入射信源之間獨立時,天線陣列接收多快拍數據的協方差矩陣RY∈CM×M為

RY=E[YYH]=AE[SSH]AH+σ2IM=ARSAH+σ2IM

(3)

式中：E表示求數學期望;RS∈CK×K表示信號協方差矩陣;σ2為噪聲方差;IM是M維的單位矩陣。

在實際應用當中,由于快拍數T有限,所以實際協方差的計算常采用下面的公式進行計算:

(4)

2 多快拍數據的降維處理

首先,對式(3)中的陣列協方差矩陣進行特征值分解,可以得到:

(5)

在式(5)左右兩端分別右乘US,得到:

(6)

對式(3)兩端也同時右乘US,則有:

RYUS=ARSAHUS+σ2IMUS

(7)

綜合式(6)和式(7)可以得到:

(8)

將式(8)的左邊記為

(9)

式中：amk表示導向矢量a(θk)中的第k個元素;vki是矩陣V=RSAHUS中第k行第i列的數據。

類似地,式(8)右邊也可以被表示為

(10)

式中：eMi為矩陣US中第M行第i列的元素,即ei=[e1i,e2i,…,eMi]T為協方差矩陣RY的特征值λi所對應的特征向量,i=1,2,…,K。

由式(9)和式(10)可以得到:

(11)

然后,令:

(12)

因此,有如下的等式關系成立:

(13)

所以有以下的等式成立:

(14)

式中：q∈LM即為新構建的觀測向量。

從式(14)可以得出,新的觀測向量與單快拍下基于ANM的DOA估計算法模型具有相同的形式,所以可以利用FIPM算法進行最優解的求解。

3 M-FIPM算法

將式(14)中的pk記為pk=|pk|ejφk,其中φk表示pk的相位,|pk|表示pk的模,則新觀測向量q的具體形式可以被表示為

(15)

(16)

形如式(16)所示的原子范數可以等價為如下SDP問題的最優解:

(17)

傳統的求解SDP問題的算法那包括IPM與ADMM。其中,IPM中每次迭代過程的復雜度為O(M6),而ADMM每次迭代過程的復雜度僅為O(M3)。不過,與IPM相比,ADMM的收斂速度慢,達到收斂時所需要的迭代次數要遠大于IPM。因此,IPM和ADMM的運算量均不理想。而本文利用FIPM求解式(17)中的SDP問題,提出基于多快拍快速內點法的ANM-DOA估計算法,記為M-FIPM。該算法每次迭代過程的復雜度也為O(M3),而收斂速度也與IPM算法相當,是目前運算效率最高的SDP問題求解算法。

形如式(17)所示的半定規劃過程又可以表示為如下的形式:

(18)

(19)

(20)

根據上述對偶錐的定義,可以將式(18)的Lagrangian函數寫為如下形式:

(21)

則可以得到式(18)的增廣KKT條件為

(22)

lg(t-qHToep-1(u)q)

(23)

式中：Toep-1(u)為Toep(u)的逆。

M-FIPM算法的具體流程如下所示。

初始值閾值ε≥0,迭代次數i=1,初始的障礙因子s1>0,初始參數u0∈intZ。

輸出入射信號的角度值θk;

步驟 1利用式(5)對接收信號協方差矩陣進行特征值分解,得到特征值和對應的特征向量;

步驟 2利用式(14)構建新的觀測向量q;

步驟 3利用M-FIPM算法求解式(17)對應的SDP問題的最優解u*;

步驟 4根據u*構建Toeplitz矩陣,并對該矩陣進行Vandermonde分解,得到入射信號角度參數θk的估計結果。

4 算法仿真實驗

本節對算法的有效性進行仿真驗。采用的陣列是8陣元的均勻線陣,對空間中4個獨立的遠場窄帶信源進行DOA估計。

在M-FIPM算法以及ANM算法中,正則化參數選取τ=0.5,對偶偏差設置為10-7,縮放因子β=10。而對于l1-svd算法,正則化參數取值不變,空間角度離散網格的間距為0.5°。同樣,MUSIC算法的空間譜峰搜索步長也為0.5°。

首先,對多重信號分類(multiple signal classification, MUSIC),l1-svd,ANM以及M-FIPM 4種算法在不同信噪比(signal to noise ratio, SNR)以及不同快拍數下的估計精度進行仿真對比。

通過均方根誤差(root mean square error, RMSE)來衡量DOA估計算法的估計精度,其計算方法如下:

(24)

式中：L∈R+是Monte-Carlo的實驗次數;θl是第l次實驗中信源DOA參數的估計結果;θ為真實的入射角度。

實驗 1令SNR在0～20 dB之間均勻變化,步長為5 dB,快拍數固定為200。信源的角度在[-60°,60°]之間隨機選取,最小角度間隔為10°。在每個SNR下做200次Monte-Carlo實驗,上述4種算法的RMSE統計結果如圖1所示。

實驗 2令快拍數分別取10,20,50,100與200,保持SNR固定為10 dB,其他仿真條件保持不變。在每個快拍數下做200次Monte-Carlo實驗,上述4種算法的RMSE統計結果如圖2所示。

從圖1與圖2中所示的仿真結果中可以看出,上述4種DOA估計算法隨SNR和快拍數的增大,估計精度都有所提升。不過,MUSIC算法的估計精度受SNR與快拍數的影響最大,在SNR較低以及快拍數較小時,MUSIC算法的估計精度最低。而隨著快拍數和SNR的增大,由于網格失配問題的存在,導致MUSIC算法的估計精度要優于l1-svd算法。

ANM算法以及本文所提出的M-FIPM算法,對SNR以及快拍數的魯棒性更好,同時能夠解決網格失配問題,所以在相同的仿真條件下,估計精度要高于MUSIC算法以及l1-svd算法。而M-FIPM算法在降維過程中由于舍棄了協方差矩陣中的小特征值分量,對噪聲進行了抑制,所以在估計精度方面還要優于ANM算法。

然后,通過仿真實驗來驗證IPM、ADMM以及M-FIPM 3種SDP問題求解算法的求解效率。而算法的求解效率主要由收斂速度以及每次迭代的算法復雜度所共同決定。

在算法復雜度方面,IPM算法每次迭代的復雜度為O(M6),而ADMM算法以及M-FIPM算法每次迭代的復雜度均為O(M3)。因此,本文所提出的M-FIPM算法在復雜度方面具有顯著優勢。

在算法收斂速度方面,SDP問題求解算法的收斂速度主要體現在算法達到收斂時所需要的迭代次數。因此,本文通過對比在不同的SDP問題維度下,上述3種求解算法達到收斂時的迭代次數來驗證所提出M-FIPM算法在收斂速度方面的性能。

實驗 3固定SNR為10 dB,快拍數為200,陣元數由10均勻變化至50,步長為10。其他仿真條件與前文的仿真實驗相同,在每個陣元數下做200次Monte-Carlo實驗,收斂時上述3種算法平均迭代次數的統計結果如圖3所示。

從圖3中所示的仿真結果中可以看出,ADMM算法達到收斂時所需要的迭代次數遠高于IPM算法以及M-FIPM算法。而M-FIPM算法達到收斂時所需要的迭代次數雖然要略高于IPM算法,但根據前文所分析的結果,M-FIPM算法每一次迭代過程的算法復雜度要遠低于IPM算法,所以M-FIPM算法的求解效率依然要優于IPM算法。

在圖4中,給出了在不同陣元數下,上述3種SDP問題求解算法的平均運行時間。而該實驗結果也驗證了本文之前的結論,即M-FIPM算法具有最高的求解效率,在SDP問題的維度相同時,所需的運算時間最短。

最后,對比MUSIC、l1-svd、ANM以及M-FIPM共4種DOA估計算法的平均運行時間。在不同陣元數下,做200次Monte-Carlo實驗,統計上述4種算法的平均運行時間,統計結果如圖5所示。從圖5的仿真結果中可以看出,在相同的仿真條件下,MUSIC算法的運算速度最快,l1-svd算法的運算速度最慢,而本文所提出的M-FIPM算法,在平均運行時間方面,要優于ANM算法。

5 結 論

針對現有快速內點法僅適用于單快拍模型的缺陷,本文提出了一種多快拍快速內點法。該算法利用天線陣列接收多快拍數據協方差矩陣的特征向量加權和來構造適用于原子范數最小化理論的觀測模型,然后將DOA估計問題轉化為SDP問題,并最終通過Toeplitz矩陣的Vandermonde分解獲得最終的DOA估計結果。仿真實驗驗證了該算法在估計精度與運行時間方面的優勢。