基于多級維納濾波的ESPRIT算法

張 超，黎仁剛，顧 軍

(中國船舶重工集團公司第七二三研究所，江蘇 揚州 225001)

基于多級維納濾波的ESPRIT算法

張 超，黎仁剛，顧 軍

(中國船舶重工集團公司第七二三研究所，江蘇 揚州 225001)

因奇異值分解引入了龐大的計算量，子空間分解類的測向算法在工程中難以實現實時處理。多級維納濾波(MSWF)算法，避免了奇異值分解運算，有效地減少了運算量。通過對多級維納濾波算法進行改進，提出了一種信源數目的估計方法。與傳統多級維納濾波算法相比，該方法提高了測向精度。通過仿真證明了該算法的有效性。

陣列信號處理；多級維納濾波；信源數目估計；測向精度

0 引 言

波達方向估計(DOA)在通信、雷達、電子對抗、導航等領域中有著舉足輕重的地位[1]。隨著信息技術的不斷發展，電磁環境變得越來越復雜，傳統的測向技術無法滿足當下的需求，以MUSIC[2]和ESPRIT[3]算法為代表的超分辨測向算法應運而生，并受到廣大研究人員的青睞。但其算法需要進行子空間分解，計算量龐大，不易實現實時性。為了減少算法的運算量，提高測向速度，許多學者也提出了一些方案，文獻[4]通過協方差矩陣實值化處理，將復矩陣分解轉換成矩陣分解，文獻[5]通過將陣元空間轉換到波束空間，降低矩陣維度。

以上所提出的算法都需要進行子空間分解，文獻[6]提到了傳播算子(PA)，此方法利用線性算子從陣列數據中提取信號子空間與噪聲子空間，以矩陣相乘代替特征值分解，減小了運算量，但此方法無法對信源數目進行估計，需要預先知道信源數目才能構造分塊矩陣。Goldstein等人在1998年提出了基于多級維納濾波的方法[7]，此方法在自適應波束形成、自適應干擾抑制等領域廣泛應用。文獻[8]將多級維納濾波技術應用在子空間分解類算法中，避免了特征值分解，大大降低了運算量。文獻[9]通過協方差矩陣來初始化數據，不需要訓練信號就能進行測向，但在低信噪比、小快拍數的情況下，誤差變大甚至可能失效。文獻[10]通過對初始化數據進行優化，使其在低信噪比、小快拍數的條件下，估計性能得到改善。文獻[11]提出了一種信號源數目與波達方向聯合估計的方法，但此方法需要進行M(陣元數)次迭代，在陣元數較多的情況下，運算量也較為龐大。本文通過改進多級維納濾波算法，優化初始矩陣，提出一種高效的信源數目估計方案，在陣元數較多的條件下，與文獻[11]中方法相比，有效降低了計算量。

1 改進算法

1.1 信號模型

對于均勻圓陣，N個全向天線均勻分布在圓周上，如圖1所示，D個不相關的遠場窄帶信號入射，方位角與俯仰角分別為θi和βi(i=1,2,…,D)。陣元的輸出數據模型為：

X=As+n

(1)

A=[α(θ1,β1),…,α(θD,βD)]

(2)

α(θ,β)=[ej2πcosθsinβ,ej2πcos(θ-φ1)sinβ,…, ej2πcos(θ-φN-1)sinβ]T

(3)

φi=2π(i-1)/N，i=1,2,…,N

(4)

式中:A為均勻圓陣的方向矢量；s為入射信號；n為陣列上的噪聲，均值為0，方差為σ2。

噪聲相互獨立，且與信號互不相關。

對于均勻圓陣而言，其方向矢量的表達式相對較復雜，且不具備均勻線陣的范得蒙矩陣形式。因此，首先要將其進行波束空間轉換變成虛擬線陣。

設空間離散傅里葉變換(DFT)矩陣:

(5)

(6)

(7)

(8)

定義變換矩陣C:

C=diag{j-M,…,j-1,j0,j-1,…,j-M}

(9)

FH=VHCWH

(10)

根據X=As+n，求出X的協方差矩陣RX，再將矩陣RX進行變換，得到:

RF=FHRXF

(11)

式中：RF為在波束空間下對應的協方差矩陣，取RF的實部代替RF。

1.2 多級維納濾波優化

維納濾波技術是在最小均方誤差(MMSE)準則下的一種線性濾波，在已知信號與噪聲相關函數的條件下，通過Winner-Holf方程RxWx=r=As+nxd對信號進行最優濾波。Goldstein等人提出的多級維納濾波技術的思想是通過對觀測數據進行多次正交投影，每次投影得到2個子空間，其中一個子空間平行于參考信號與上次觀測數據的互相關矢量，另一個子空間正交于這個子空間，然后再對這個正交于互相關矢量的子空間用同樣的方法繼續進行正交投影。其結構如圖1所示。

多級維納濾波算法需要參考信號di(k)與觀測數據Xi(k)，而實際工程中，并不一定能夠獲得，因此采用任意單元的數據來代替參考信號。信號之間相互獨立，假設為高斯白噪聲。觀測信號每次投影，分別得到與觀測數據平行的分量和正交的分量，這些代表的就是信號子空間與噪聲子空間。具體步驟如下所示：

第1步，初始化觀測數據和陣列接收數據：

(12)

(13)

式中：R為協方差矩陣。

第2步，求互相關矢量hi(i=1,2,3，…,D)：

(14)

(15)

(16)

(17)

第3步，估計信號源數目，確定信號子空間：

(18)

(19)

式中：S代表信號子空間；G代表噪聲子空間。

數據初始化使用協方差矩陣來代替，在文獻[10]中提出的初始數據優化方法，通過將協方差矩陣行取平均，得到：

(20)

此方法應用于均勻線陣，在低信噪比、小快拍數條件下，適當提高了算法的估計性能。

本文算法是基于均勻圓陣的，經過波束空間映射到虛擬均勻線陣，模式數為2M+1。本文提出的優化方法是選用協方差矩陣第M行初始化數據：

(21)

針對均勻圓陣，本文所提出的優化方法比起文獻[10]中提出的方法，在低信噪比、小快拍數下的精度更高，下文仿真實驗可見。

1.3 信號源數目估計

在常規的子空間分解測向算法中，信號源數目的估計一般采用Akaike信息論(AIC)準則和最小描述長度(MDL)準則，這2個準則具有相同的似然函數，只是懲罰函數不一樣。但這2個準則的使用前提需要對陣列協方差矩陣進行特征值分解。而在多級維納濾波算法中，通過正交投影取代了特征值分解。陣列協方差矩陣分解所得到的特征值分別代表信號或噪聲功率的大小。與特征值分解類似，通過多級維納濾波算法得到的正交投影信號d1(k),d2(k),…,dM(k)，可以用σi來代表其功率：

(22)

表1 不同信噪比下σi與特征值λi的比較

從表1可以看出，多級維納濾波算法在低信噪比時，信號能量會泄露到噪聲子空間中，導致無法判斷信號子空間。文獻[7]中提出了一種信號源數目的估計方法，判斷信號子空間。但該方法在陣元數比較多的時候，計算量較大。本文提出一種估計方法，適合多陣元數條件下的信號源數目估計。

第1步：初始化后進行4次投影得到正交投影信號，求出代表其功率的σi。

第2步：構造判決函數，判斷信號子空間需要濾波的級數：

fori=1:M-4

δi=(σi+1+2σi+2)/3

Δi=(σi-δi)/δi/2

ifΔi<1

ThenD=i

elseifσi<1

ThenD=i

endif

endfor

根據判決，得到信號源數目為i。

第3步：構造信號子空間S=[h1,h2,…,hD]。

1.4ESPRIT算法

經過多級維納濾波正交投影，得到的信號子空間S=[h1,h2,…,hD]相似于B：

(23)

式中：Jm(ζ)為第一類貝塞爾函數，m為階數，ζ為變量。

根據貝塞爾函數的性質可得：

ΓZ0=Z-1Φ+Z1Φ*

(24)

Γ=(λ/πr)diag{-(M-1),…,0,…,M-1}

(25)

Φ=diag(sinβ1ejθ1,sinβ2ejθ2,…,sinβDejθD)

(26)

式中：Z-1，Z0和Z1分別為信號子空間Z前2M-1行，中2M-1行和后2M-1行。

最后求解得到Φ，即可估計出信號源方位角與俯仰角。

2 仿真與分析

本文采用陣元數為16的均勻圓陣，模式數2M+1為9，2個不相干的信號源入射角分別為(75°,20°)與(60°,55°)，對比常規波束空間ESPRIT算法(算法1)、文獻[10]中所提出的算法(算法2)以及本文算法(算法3)3種算法的性能。定義聯合均方根誤差(RMSE)為：

(27)

式中：P為蒙特卡洛試驗次數；θi與βi分別為第i次試驗的方位角與俯仰角的估計值；θ與β分別為方位角與俯仰角的真實值。

試驗1：選擇信噪比為4dB與6dB，在不同的快拍數下，進行200次蒙特卡洛試驗，仿真結果如圖3所示。

從圖3可以看出，小快拍數的情況下，本文算法的估計性能優于常規的多級維納濾波算法，接近協方差矩陣分解算法的估計性能。隨著快拍數的增加，各算法的估計性能也隨之提高。

試驗2：選擇快拍數為40次與80次，在不同的信噪比下，進行200次蒙特卡洛試驗，仿真結果如圖4所示。

從圖4可以看出，在低信噪比的情況下，本文算法的估計性能優于常規的多級維納濾波算法。隨著信噪比的提高，各算法的估計性能也隨之提高，并且均方根誤差趨于一致。

試驗3：對于本文提出的信號源數目估計方法，分別在低信噪比、小快拍數下進行200次蒙特卡洛試驗，觀察其成功概率。

從圖5、圖6中可以看出，在低信噪比時隨著快拍數增加和信噪比的提高，信號源數目估計的成功概率都將提高。本文提出的信號源數目估計方法中引入了σi<1這一判決條件，可避免在大信噪比時噪聲功率偏小，導致噪聲功率之間的相對偏差較大，引起估計錯誤。

傳統的多級維納濾波算法的運算量為O(M3)，而本文所提出的信號源數目估計方法的運算量為O[M2(K+3)]，其中M為陣元數目，K為信號源數目。相比文獻[11]中的信號源數目估計方法，本文方法在陣元數較多而信號源數目較少的情況下，能有效減少運算量。

3 結束語

根據多級維納濾波的特點，本文提出了一種基于多級維納濾波的改進算法，用協方差矩陣代替參考信號，并通過優化初始數據來提高估計精度。與常規的子空間分解算法相比，其避免了特征值分解，大大減小了運算量。此外，在多級維納濾波投影過程中，低信噪比條件下存在子空間泄露問題。本文提出了一種信號源數目的估計方法，能夠在低信噪比時有效估計出信號源的數目。尤其是在大陣元數、信號源數目不多的情況下，本文方法比常規多級維納濾波算法運算量更低，更適用于工程實現。

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

[2]SCHMIDTRO.Multipleemitterlocationandsignalparameterestimation[J].IEEETransactionsonAntennasandPropagation,1986,34(3):276-280.

[3]ROYR,PAULRAJA,KAILATHT.ESPRIT——asubspacerotationapprochtoeatimationofParameterofcisoidsinnoise[J].IEEETransactionsonAcoustics,SpeechandSignalProcessing,1986,34(5):1340- 1342.

[4] 孫繼紅,楊科,張瑾.一種實值MUSIC算法[J].空間電子技術,2006(3):36-38.

[5] 李紹濱,趙宜楠,胡航.均勻圓陣波束空間二維DOA算法[J].哈爾濱工業大學學報,2004,36(6):796- 798.

[6]MARCOSS,BEBIDIRM.Onahighresolutionarrayprocessingmethodnon-basedoneigenanalysisApproach[C]//ETTERDM.InternationalConferenceonAcoustics,speechandSignalProcessing.NewMexico,USA:IEEEPress,1990:2955-2958.

[7]GOLDSTAINJS,REEDIS,ScharfLL.Amultistagerepresentionofthewinnerfilterbasedonorthogonalprojection[J].IEEETransactionsonInformationTheory,1998,44(7):2943-2959.

[8] 黃志忠,李萬玉,強勇.多級維納濾波在子空間類測向算法中的應用[J].火控雷達技術,2009(1):11-14.

[9] 于紅旗,黃知濤,周一宇,徐欣.一種不需要特征值分解的ESPRIT方法[J].信號處理,2008(3):357-360.

[10]李堰,宋愛民,劉劍.多級維納濾波測向算法的參考信號優化方法[J].電訊技術,2010(12):75-80.

[11]包志強,吳順君,張林讓.一種信源個數與波達方向聯合估計的新算法[J].電子學報,2006(12):2170- 2174.

ESPRIT Algorithm Baesd on Multi-stage Winner Filtering

ZHANG Chao,LI Ren-gang,GU Jun

(The 723 Institute of CSIC,Yangzhou 225001,China)

Because of the enormous calculating quantity introduced by singular value decomposition,the direction finding algorithm based on subspace decomposition is difficult to be realized in engineering in real time.The multi-stage Wiener filtering (MSWF) algorithm avoids the singular value decomposition calculation,reduces the calculating quantity effectively.This paper puts forward an estimation method of source number by improving MSWF algorithm.The algorithm improves the direction-finding accuracy compared with traditional MSWF algorithms.The validity of this algorithm is validated through simulation.

array signal processing;multi-stage Wiener filtering;estimation of source number;direction-finding accuracy

2016-09-01

TN971.1

A

CN32-1413(2017)01-0074-05

10.16426/j.cnki.jcdzdk.2017.01.016