基于DSP的快速MUSIC測角算法

金芳曉,崔 偉

(1. 中國電子科技集團公司第三十八研究所, 安徽 合肥 230088)(2. 合肥移瑞通信技術有限公司, 安徽 合肥 230000)

0 引 言

基于陣列信號處理的MUSIC算法[1],從空間分辨率上突破了瑞利限制,實現了高分辨波達方向(DOA) 估計向超分辨DOA估計的過渡。其核心原理是以信號子空間和噪聲子空間的正交性為基礎,劃分空間來進行DOA估計[2],廣泛用于雷達目標識別跟蹤、無線電監測、電子對抗、參數估計等領域。但MUSIC算法在超分辨DOA估計時,存在運算復雜度高,不利于硬件實現等問題[3]。

自從MUSIC 算法提出以來,大量學者對算法的實時硬件實現、降低其運算量進行了深入研究。MUSIC算法的運算量主要集中在矩陣特征值分解和空間譜搜索兩部分。針對特征值分解部分,文獻[4]利用接收數據近似快速估計信號子空間和噪聲子空間,簡化了特征值分解問題。文獻[5]中提出了一種采用多級維納濾波求解信號特征矢量和信號空間維數,有效地降低了算法運算量。文獻[6]為進一步提高DOA估計精度,提出了基于降維子空間的多重信號分類快速算法,利用最大噪聲子空間對估計結果修正,提高測角精度。對于MUSIC算法中空間譜搜索部分的簡化,文獻[7]和文獻[8]分別采用陣元降維和有限域搜索等方法研究了空間譜的快速搜索問題。文獻[9]空間相移因子構建方向向量,進而利用數字傅里葉變換(DFT)實現MUSIC譜函數的快速計算。這些算法都從MUSIC算法中關鍵子步驟入手降低算法運算量,但是距離工程實現還有一定的差距,且算法中用于空間譜搜索的所有搜索角度的導向矢量需占用較大的內存資源進行存儲,也在一定程度上影響了工程應用。

為此,本文結合TI1843 DSP處理平臺,從降低空間譜搜素運算量和減少陣列導向矢量占用內存空間兩方面出發,對基于數字信號處理(DSP)的快速MUSIC算法的硬件實現進行研究。在空間譜搜索方面,考慮MUSIC算法硬件工程實現中,在測角范圍內存在有多目標,為保證測角角精度和角分辨率,分為粗搜和精搜兩步減少譜峰搜索的次數。以角分辨率為步進進行測角范圍內的粗搜;在多個目標的方向,以粗搜索的步進作為搜索范圍區間,采用黃金分割算法[10]進行精搜,在保證測角精度的同時,大大降低了搜索次數。在空間譜搜索方面,本文利用陣列陣元間的關系,采用迭代的方式構建陣列導向矢量,大大縮短了導向矢量的構建時長,結合上述粗搜+精搜方式進行譜峰搜索,搜索次數大大減少,不再需要提前存儲導向矢量,大大減少了內存空間。綜上所述,本文算法針對MUSIC算法處理時長和內存資源消耗最大的兩大部分,從工程實現上入手,大幅減少了空間譜譜峰的搜索次數,迭代生成導向矢量不再占用內存資源,具有較大的實際應用價值。

1 MUSIC算法介紹

1.1 算法模型

考慮K個遠場的窄帶信號s(t)=[s1(t),s2(t), …,sK(t)]T,從不同方位角θ1,θ2, …,θK入射到M個陣元組成的均勻線陣,相鄰陣元間距為d,則t時刻接收到的數據x(t)=[x1(t),x2(t), …,xM(t)]T為

X(t)=AS(t)+N(t)

(1)

式中:A=[a(θ1),a(θ2), …,a(θK)]為M×K維的陣列流形,也稱方向矩陣;a(θk)為M×1維的方向矢量,a(θk)=[1,ej2πdsinθk/λ, …, ej2π(M-1)dsinθk/λ]T,θk表示第k個信號的方向角,λ為波速;n(t)為M×1維的高斯噪聲矢量n(t)=[n1(t),n2(t), …,nM(t)]T。

陣列數據的協方差矩陣為

R=E[XXH]=AE[SSH]AH+σ2I=ARsAH+σ2I

(2)

由于信號與噪聲相互獨立,數據協方差矩陣可分解為與信號、噪聲相關的兩部分,其中Rs是信號的協方差矩陣,ARsAH是信號部分。

對R進行特征值分解有

(3)

式中:US是由大特征值對應的特征矢量張成的子空間即信號子空間;UN是由小特征值對應的特征矢量張成的子空間即噪聲子空間;ΣS為大特征值組成的對角陣;ΣN為小特征值組成的對角陣。理想條件下,數據空間中的信號子空間與噪聲子空間是相互正交的,即信號子空間中的導向矢量與噪聲子空間正交,如下所示。

aH(θ)UN=0

(4)

經典的MUSIC算法正是基于上述這個性質提出的,但考慮到實際接收數據是有限長的,數據協方差的最大似然估計為

(5)

(6)

MUSIC算法的譜估計公式為

(7)

1.2 算法步驟及運算量

根據上面介紹,整理MUSIC算法的步驟歸納如下。

(2) 對得到協方差矩陣進行特征值分解 (EVD):R=UΣUH。

(3) 變化搜索角度θ,按照式(7)空間譜計算譜函數并進行譜峰搜素,譜峰所對應的角度便是DOA估計值。

進一步,將算法移植到TI1843S DSP的硬件系統上,進行了100次Monte Carlo隨機試驗。試驗條件設定如下:8陣元均勻線陣,單快拍,角度搜索次數100。表1給出了MUSIC算法各個子步驟DOA估計所需要運行時間。從表1可以看出,在三個子步驟中,耗時最多的便是空間譜構建及譜峰搜索(表1中的離散譜估計部分),是影響MUSIC算法實時性的關鍵因素。為此,本文提出了簡化空間譜構建和搜索的方法。

表1 MUSIC算法子步驟運算量評估Tab. 1 Evaluation of sub-step computation of MUSIC algorithm

2 導向矢量搜索和構建的簡化

為減少空間譜譜峰搜索次數,同時保證測角精度,本文采用黃金分割算法對空間譜搜索進行改進。該方法可以在滿足算法精度要求的前提下,成倍的降低運算量,且工程實現方便。

另一方面,傳統MUSIC算法在硬件實現時,一般根據系統測角范圍和角精度要求,需要提前構建好空間譜矩陣,存儲在系統中,用于測角時候進行各個角度的空間譜譜峰搜索。但是,這種方式需要較大內存,且測角范圍越大,角精度越高,所需的內存越大。而且,采用本文提出的黃金分割算法對空間譜進行搜索,并不需要調用測角范圍內的所有角度對應的空間譜,所以采用傳統遍歷的方式進行DOA估計很大程度上浪費了存儲空間。為此,本文利用陣元間相位關系迭代構建空間譜,大大較少了空間譜的構建時間。進一步,結合黃金分割算法僅構建黃金分割點所對應角度的空間譜,從而實現在空間譜搜索過程中構建相應角度的空間譜,這樣在保證系統實時性的前提和測角精度的前提下,大大地減少了內存資源的浪費。

2.1 黃金分割空間譜搜索算法

為減少搜索次數,本文采用一維黃金分割精搜算法取代傳統空間譜搜索。首先簡單介紹黃金分割算法原理[11],具體如下。

函數f(x):設該函數在區間(a,c)上有且僅有一個極大值,其中b點位于(a,c)區間的黃金分割點,x點位于(b,c)區間的黃金分割點,如圖1所示。最大值的搜索過程如下:區間(a,c)上計算的b和x兩點對應的函數值,當f(b)>f(x)時,即滿足f(b)>f(a)且f(b)>f(x),則可判定極大值點位于三元點af(c)且f(x)>f(c),則可判定極大值點位于三元點組b

圖1 黃金分割函數示意圖Fig.1 Schematic diagram of golden section function

本文將上述介紹的黃金分割應用于MUSIC算法譜峰搜索中,在保證搜索精度的基礎上,大大地減少了搜索次數?？紤]到實際情況,在測角范圍內存在有多目標,而黃金分割算法只適用于的搜索范圍內有且僅有一個極大值的設定下,那么如何設定搜索的步進并采用黃金分割對多個目標進行逐一搜索,是算法工程應用的關鍵。為此,在已知信源個數的情況下,將空間譜搜索分為粗搜和細搜兩個部分。其中,粗搜是為了檢測到多個目標,精搜則是為了保證測角的角精度。具體的采用黃金分割精搜對多目標譜峰搜索的方式如下。

根據工程需要,設定測角搜索范圍為[-θ,θ],所需角度分辨率為Δθ,角精度為ε,目標個數為K。陣列的波束寬度的近似計算公式為

采用上述粗搜+精搜的方式,在保證角精度和角分辨率的基礎上,減少了搜索范圍。傳統測角算法,為保證測角的精度,應以角精度ε為搜索步進,對測角范圍內全角度搜索。而采用本文算法全角度粗搜是以角分辨率Δθ為步進,由ε≈Δθ/5可知,全角度粗搜次數約為傳統全角度搜索的1/5。而精搜部分只是在目標角度局部,以角精度為步進進行搜索,且結合黃金分割方法,進一步減少了局部搜索的次數,與局部遍歷搜索相比搜索次數約為其1/3,粗搜和精搜總搜索次數較傳統遍歷的方式相比大大減少。下面以具體參數為例,對比本文搜索方式和傳統遍歷的方式搜索次數上減少的程度。設定測角搜索范圍為[-30°,30°],所需角度分辨率為1°,角精度為0.2°,目標個數為3。采用傳統遍歷搜索的方式,為保證角分辨和角精度,搜索步進設為0.2°,搜索次數約為300次。而采用本文搜索方式,第一步粗搜以角分辨率1°為步進,搜索次數為60次;第二步精搜,對于第k個目標,在[θk-1°,θk+1°]的2°范圍內以角精度0.2°為步進,采用黃金分割算法進行精搜,搜索次數為4次,4個目標共需12次。由此可得,本文搜索精搜+粗搜一共72次,對比傳統遍歷的方式,搜索次數大大減少,且根據分析可得,搜索范圍越大,角精度要求越高,傳統遍歷方式搜索次數越多,而本文所提出的搜索方式優勢越明顯。

2.2 導向矢量構建

下面基于均勻線陣導向矢量的模型,推導導向矢量構建的簡化方法。在DSP中復數采用實部和虛部分開表示,故將式(1)中的導向矢量a(θk)=[1,ej2πsinθkd/λ,…,ej2π(M-1)sinθkd/λ],a(θk)拆分為實部和虛部兩部分,即

a(θk)=areal(θk)+j·aimag(θk)

(8)

其中,

areal(θk)=[1,cos(2πsinθkd/λ,…,cos(2π(M-1)

sinθkd/λ)]

aimag(θk)=[0,sin(2πsinθkd/λ),…,sin(2π(M-1)

sinθkd/λ)]

便可利用areal(θk)和aimag(θk)來構建導向矢量a(θk)。此方式雖然簡單直接,但是其運算量卻較為復雜。以M個陣元為例,構建θk方向的導向矢量a(θk)需進行4(M-1)次三角函數運算。為此,常規方法一般為提前存儲好所需的導向矢量,方便角度搜索直接調用,且角精度要求越高,所需的存儲空間越大。

為簡化構建,減少存儲空間,基于導向矢量a(θk)的數學表達式,本文利用三角函數和差化積和推導陣元間導向元素的遞進關系,大大減少運算量。具體推導如下:

取areal(θk)和aimag(θk)中的任意元素,令

qreal(m,θk)=cos(2π(m-1)sinθkd/λ),

m=2,3,…,M-1

(9)

qimag(m,θk)=sin(2π(m-1)sinθkd/λ),

m=2,3,…,M-1

(10)

利用和差化積,可得

qreal(m,θk)=cos(2π(m-2)sinθkd/λ+2πsinθkd/λ)=

cos(2π(m-2)sinθkd/λ)cos(2πsinθkd/λ)-

sin(2π(m-2)sinθkd/λ)sin(2πsinθkd/λ)=

qreal(m-1,θk)qreal(1,θk)-qimag(m-1,θk)·

qimag(1,θk)

(11)

進一步推導,可得

qimag(m,θk)=sin(2π(m-2)sinθkd/λ+2πsinθkd/λ=

sin(2π(m-2)sinθkd/λ)cos(2πsinθkd/λ)+

cos(2π(m-2)sinθkd/λ)sin(2πsinθkd/λ)=

qimag(m-1,θk)qreal(1,θk)+qreal(m-1,θk)·

qimag(1,θk)

基于此可以看出,只需求解qreal(2,θk)和qimag(2,θk),后面m-2個導向矢量元素qreal(m,θk)和qimag(m,θk),便可利用式(11)和式(12)迭代得到。采用此方法,僅需要4次三角函數運算來構建qreal(2,θk)和qimag(2,θk),大大地降低了導向矢量構建的運算復雜度。

采用此方式迭代生成對應角度的導向矢量,計算量小。同時,配合前文提出的粗搜+精搜的搜索方式,空間譜無需按角精度為步進遍歷搜索,搜索次數減少。為此,采用本文方法空間譜搜索的次數少,且每次搜索導向矢量構建的運算量少,故每次搜索可根據所需的角度生成相應的導向矢量即可,無需提前存儲導向矢量,大大減少了存儲空間,具有較強的實際應用價值。

3 試驗結果

為了驗證本文提出的簡化的導向矢量構建和空間譜搜索算法,在實際DSP中運用的有效性,對傳統方式和本文方式從空間譜搜索、導向矢量構建時長和所占內存大小三個方面進行工程化對比分析。試驗的硬件平臺采用TI1843 DSP環境,其主頻為600 MHz進行了100次Monte Carlo隨機試驗。試驗測試條件采用8陣元均勻線陣,目標個數為2,為滿足角分辨率要求,設定兩個目標的角度分別為5°和10°。針對導向矢量構建時長和內存占用情況以及空間譜搜索時長的驗證分析,設定如下試驗:一、導向矢量構建時長和內存占用情況試驗,設定100個搜索角度,分別采用傳統和本文兩種方式構建100個角度的導向矢量,并統計總的構建時長和內存占用情況進行對比;二、分析對比多目標下空間譜搜索時長,對角精度(也就是搜索步進)分別為1°和0.1°兩種情況,采用傳統遍歷搜索和本文粗搜+精搜兩種方式進行空間譜搜索的時長對比。

試驗對比結果見表2,由表2可以看出:一、采用本文所提出的遞進的導向矢量構建方式,在構造100個方向角度導向矢量的相同條件下,與傳統構建方式相比構建時間大大縮短。進一步分析內存資源的占用情況,由于本文采用陣元間迭代構建導向矢量的方式,構建時間大大縮短,且黃金分割譜峰搜索也不需遍歷所有角度,故可采用搜索時生成對應角度的空間譜的方式進行角度估計,無需再設置內存用來存儲100個搜索角度8陣元的導向矢量,大大地減少了內存空間的占用。二、空間譜搜索時長方面,傳統的方式是以角精度為步進遍歷所有角度進行搜索,由表2可以看出,采用傳統遍歷方式,以0.1°為角精度與1°為角精度相比時長約為其10倍,而采用本文粗搜+精搜的方式,時長上大幅縮短,且隨著角精度的提高,黃金分割譜峰搜索的優勢越明顯。

表2 空間譜搜索、構建時長和內存對比分析Tab.2 Comparative analysis of spatial specturm search time, construction time and memory

4 結束語

為提高MUSIC算法的DSP硬件實現技術,本文從MUSIC算法運算量和所占資源最大的子步驟入手,提出了簡化的空間譜構建和譜搜索方法。試驗中采用TI1843 DSP硬件平臺,驗證了方法的可行性和有效性,對MUSIC算法工程實現具有較高的實際應用意義。