摘 要:流形模糊發生的原因是陣列流形上的多個矢量線性相關,子空間類算法無法靠自身解模糊。為解決子空間類算法中的模糊問題,提出了一種新的基于協方差陣擬合的解模糊算法。該方法用預估的導向矢量與入射波功率來擬合協方差陣,用最陡下降法尋找最優的功率參數,認定功率較大的入射波為真實來波。與子空間算法結合,可解決測向模糊問題。數值仿真表明,此方法在陣元數量較多時仍能有效工作。
關鍵詞:陣列; 測向; 模糊; 協方差陣; 擬合; 子空間
中圖分類號:TP301.6文獻標志碼:A
文章編號:1001-3695(2009)09-3293-03
doi:10.3969/j.issn.1001-3695.2009.09.026
Approach of ambiguity resolution based on covariance matrix fitting
LIU Hong-shenga,b, XIAO Xian-cib
(a.School of Communication Information Engineering, b.School of Electronic Engineering, University of Electronic Science Technology of China, Chengdu 610054, China)
Abstract:The reason of the manifold ambiguity is the linear dependence of steering vectors on manifold, this kind of ambiguity can not be resolved by subspace algorithms themselves. To resolve the manifold ambiguity, this paper proposed an approach based on covariance matrix fitting. The implementation of the approach was to fit the covariance matrix by pre-estimated manifold vectors and corresponding powers, then found optimal solution of powers by method of steepest descent. At last, regarded the DOAs with bigger powers as true DOAs. Combining with subspace algorithms, the approach could identify the true DOAs. Simulation verify its effectiveness and robustness when there are more sensors in array
Key words:array; direction finding; ambiguity; covariance; fitting; subspace
如果將子空間類測向算法應用于某些陣列(如稀疏線陣),可能出現測向模糊問題[1~7]。在子空間類算法中,測向模糊表現為空間譜上出現了沒有真實來波對應的偽峰值,因而無法確定真實的波達方向。Schmidt[8]在提出子空間算法的同時也研究了測向模糊問題,發現模糊的產生是由于多個陣列流形矢量線性相關造成的,他根據流形矢量間的相關性,定義了模糊的“階”數。
Manikas等人[1~3]對這類模糊問題進行了進一步的研究。他們從陣列流形的幾何特性入手,發現了陣列流形很多特性,如線陣的陣列流形是一條嵌于N維復空間的超螺旋曲線(其各階曲率為常數)。在這些特性的基礎上,提出了一套確定陣列測向模糊發生角度的方法,根據這些方法,當陣列幾何結構確定后,可找到一些可能發生模糊的角度集合。
由于流形模糊是子空間類測向算法特有的,可以通過附加方法進行解模糊。Abramovich等人利用Pillai等人[9]提出的陣列加強技術,對稀疏線陣的空間協方差陣進行加強,加強后的陣列維數為稀疏線陣的相伴陣(co-array),在加強陣的基礎上應用子空間類算法即可得到無模糊的DOA[10]。這種解模糊方法簡單實用且速度快,缺點是與陣列結構高度相關,且目前得到結論只適用于一類特殊線陣(如全加強陣)。另外一種解模糊算法則與陣列結構關系不大,該方法先假定子空間方法測得的所有DOA為真實的入射波方向,然后對各入射波的功率進行估計,如果在某個DOA入射功率為零,則該DOA為模糊方向。基于第二種方法的原理,Abramovich等人[10]采用協方差陣對角擬合準則,利用線性規劃理論得出一種解模糊方法(LP法)。該方法由于利用了現有的線性規劃算法,速度快,但是當陣列規模稍大(大于4)時,線性規劃算法的敏感性失效。本文提出一種基于協方差陣直接擬合準則的解模糊方法,直接利用相關陣擬合方法解決模糊問題,該方法直接對入射波功率進行估計。由于采用了最陡下降法尋優,相對于線性規劃方法更為穩健,可適用于較大規模的陣列。
1 子空間類算法的測向模型與模糊生成集
設有M個遠場輻射源S(t)入射由N個陣元構成的陣列,陣列的輸出矢量X(t)可由下式描述:
X(t)=A(q)S(t)+N(t)(1)
其中:S(t)∈CM×1是入射波信號復包絡;N(t)∈CN×1是零均值加性高斯白噪聲。導向矩陣A(q)∈CN×M具有以下形式:
A(q)=[a(q1)…a(qM)](2)
其中:q=[q1…qM]是波達方向(DOA);a(qi)為導向矢量,導向矢量的軌跡為陣列流形[a(q):q∈Ω],Ω是參數空間。
觀測數據協方差陣為
R=E{X(t)XH(t)}=A(q)PAH(q)+σ2I(3)
其中:P=diag(p1,…,pM),是由各入射波功率構成的對角陣。
對于線陣,DOA q只是方位角θ∈[0,π),陣列流形a(θ)具有如下形式:
a(θ)=exp(-jπ(r cos θ))
其中:r是線陣陣元坐標,為一N維矢量。可見a(θ)是嵌于N維復空間CN上的曲線。如果這條曲線上有多個線性相關的矢量,那么就會發生流形模糊。Manikas等人[1,3]發現這條曲線是一條超螺旋(常曲率)曲線,對這條曲線按一定規則分割,分割點上的矢量構成線性相關的矢量組(模糊生成集)。
對于面陣,DOA q是二維參數:方位角θ∈[0,2π),仰角∈[0,π/2)。陣列流形為
a(θ,)=exp(-jπ(rx cos θ+ry sin θ)cos )
其中:矢量(rx,ry)是陣元的位置坐標,面陣的流形隨著(θ,φ)變化而在N維復空間中展開成一個復曲面。該曲面可看做由兩簇(θ,φ)曲線覆蓋而成。當θ確定時,曲線也是一條超螺旋線,利用線陣結論,可對其分割從而找到模糊生成集[2]。
2 流形模糊與可辨識性
對于子空間分解類算法,模糊現象是由于陣列流形上矢量的線性相關而發生的,這種模糊被稱為流形模糊。
為克服流形模糊,一種直觀的方法是改造陣列流形,使之不存在任何線性相關的矢量組合。由于均勻線陣(ULA)的導向矩陣是范得蒙陣,不同DOA對應的導向矢量必定不相關,是目前惟一可確定的無模糊線陣,若要避免線陣流形模糊,可將模糊陣列流形改造成ULA的流形,顯然只有一些特殊幾何布局的線陣可以實現這樣的改造。全加強陣就是這樣一種陣,它的伴隨陣與ULA相同,可以通過加強技術構造虛擬陣元,對協方差陣完型(complete),完型后的協方差陣等效于ULA的協方差陣,對其應用標準的子空間類算法,算法中搜索的流形就是ULA流形,因此可實現無模糊測向。由于平面陣、立體陣的復雜性,尚不能找到類似的方法來實現解模糊。
另一種解流形模糊的方法是通過估計輻射源的其他特征參數(如功率)來確定模糊DOA。這種方法的可行性取決于該測向系統是否可辨識。
對于如式(1)的統計參數模型,如果入射信號S(t)是零均值的隨機窄帶信號,且各信號分量彼此不相關。可知陣列輸出服從N維復高斯分布:X(t)~CN(0,A(q)PAH(q)+σ2I)。那么不可辨識作如下定義:
如果存在兩組不完全相同的參數q≠q,P≠P使
CN(0,A(q)PAH(q)+σ2I)=CN(0,A(q)PAH(q)+σ2I)(4)
則稱參數q與q不可辨識。
因為在實踐中是采用測向模型式(1)的樣本來估計模型里邊的參數(q,P),一旦出現不可辨識情形,參數(q,P)和(q,P)都對應著相同的樣本,無法由這些樣本來分辨兩者。更為具體地說,如果(q,q)是兩組不完全相同的角度,總有相應的非零功率(P,P)使式(4)成立,這將無法通過統計模型式(1)產生的樣本來區分兩者。所以在研究流形模糊的解法前,先要確定流形模糊是否可辨識。
式(4)可以簡化為A(q)PAH(q)=A(q)PAH(q)。
當各入射波彼此互不相關時,P=diag(p1,…,pn),P=diag(p1,…,pn),都是對角陣,可進一步簡化為
∑Mi=1pia(qi)aH(qi)=∑Di=1pia(qi)aH(qi)(5)
這就是測向模型式(1)不可辨識的充分條件。式(5)中的D是模糊集里的角度個數,M是真實的入射波個數,D>M,模糊方向個數k=D-M。
從統計模型角度定義的不可辨識與流形模糊不是一個概念,雖然流形模糊發生了,但是如果可辨識就可設計相應的方法辨識哪些方向是真實的來波方向。
要確定一個流形模糊集里的模糊角度是否可辨識,需要考察它們是否滿足式(5)。根據Abramovich等人[10]的研究,在輻射源不相干的情況下,用P-M(Proukakis和Manikas)技術找到的線陣模糊集,如果其秩小于陣元數,那么其中的流形模糊角度都是可辨識的。本文進一步研究了面陣P-M流形模糊集,證實也是可辨識的。
上邊的討論說明,當陣列的流形模糊可辨識時,沒有任何兩組不完全相同的角度參數及相應的非零功率參數滿足式(5)。比較式(5)與式(3)可知,兩者就差一個對角陣,這個對角陣各元素為噪聲功率,式(5)兩邊加上該對角陣就構成了協方差陣。實踐中協方差陣已經估計出來,而模型式(1)又是可辨識的,若想由含有流形模糊角度的D個波達方向來擬合協方差陣,那么k個模糊角度所對應的功率必然為零。利用這個特性,就可以確定子空間算法所測的方向中,哪些是真實的波達方向,哪些是模糊方向。
3 基于協方差陣擬合的最小二乘擬合的解模糊算法
當由觀測數據估計出協方差陣R^后,子空間類算法需要對其作特征值分解,得出噪聲子空間及信號子空間,最后在陣列流形上搜尋距信號子空間最近的點,這些點被認為是需要估計的DOA。但是,當模糊出現時,在其中的一些點所代表的方向上并沒有真實的入射波,或者可以說其入射波的功率為零,所以可以根據這一特點來確定真實的DOA。
設由子空間類算法已經估計出一組DOA q^,采用P-M方法對其分析,確定q^中有分量屬于一模糊集,且滿足秩條件,也就是說可能出現模糊。
設對應于q^i的入射波功率為p^i,那么可以由這些入射波功率及流形矢量A(q^)擬合協方差陣:C=A(q^)P^AH(q^)。其中:P^=diag(p^1,p^2,…,p^M)。擬合誤差為ε=‖R^-C‖2。這里R^=(1/L)∑Lt=1X(t)XH(t),L是快拍數。
求ε對各入射波功率的偏導數:
ε/p^i=tr((R^-C)(R^-C)H)/p^i=2tr((C-R^)C/p^i)
令上式為零,取i=1,…,M可得一組方程:
aH(q^i)Ca(q^i)=aH(q^i)R^a(q^i)
各入射波功率可由下式求出:
diag(P^)=((AHA)⊙(AHA))-1 diag(AHR^A)
其中:diag(P^)表示一個由括號中矩陣對角線元素構成的矢量;⊙表示Hadamard乘,也就是參與運算的兩個等維數矩陣具有相同下標的元素相乘,結果是由乘積構成的一個等維數的新矩陣。
鑒于求逆的困難,實際中可采用最陡下降法進行迭代運算。Assumpcao[11]給出了算法的詳細步驟。
在已知有M個輻射源的情況下,將Manikas的模糊定位方法與本章的入射波功率估計方法結合,可以得到一個解流形模糊的方法:采用子空間類算法估計DOA,得q^=[q^1…q^D],也就是說空間譜上有D個峰值。如果D=M,則無模糊發生,子空間算法估計出的DOA正確,測向完成;如果D>M,則發生測向模糊,保存子空間類算法中得到的協方差陣R^并計算βi=aH(q^i) R^a(q^i),進入下面的解模糊步驟。
對應于各方向q^=[q^1…q^D],設各次迭代的入射波功率初始值為p^=[p^…p^D],計算C=A(q^)P^A(q^),ε=‖C-R^‖2。
a)計算各方向梯度:
di=2(βi-aH(q^i)R^a(q^i))
b)計算迭代步長:
μ=(-∑d2i)/(2‖A(q^)diag(d1,d2,…,dM)AH(q^)‖2)
c)迭代后各入射波功率:
p^′i=max(p^i+μd^i,0)
計算C′=A(q^)diag(p^′1,…,p^′M)AH(q^),ε′=‖C′-R^‖2。
比較本次迭代前后的擬合誤差:|ε′-ε|。如果大于容差,用迭代后的各功率p^′i取代c)中的各假設功率p^i,繼續迭代;如果小于容差,迭代結束,把p^i當做各入射波功率的估計值。考察估計的各入射波功率的大小,按大小次序排列,取前M個入射波的DOA作為真實來波方向,其他方向舍棄。
4 數值仿真
在數值仿真中首先確定陣列的模糊集,得出模糊出現的角度,然后按這些角度設定入射波,人為產生子空間算法的測向模糊。然后按協方差陣擬合算法算出真正的入射波方向。
已經在多個陣列上完成上述數值仿真,結果表明算法可行。因為本算法的重點在于解模糊,所以在仿真中都采用了較高的信噪比(20 dB)、各仿真入射波等功率和較多的快拍數(500)。當擬合誤差ε<0.01時停止迭代求解。
1)線陣
采用文獻[2]中的例3,陣元坐標矢量為
r=[-3,-2.2,1.1,1.9,2.2]
任取該陣的一個秩為4的模糊集為
[29.69°,61.04°,84.29°,106.57°,132.05°]
實驗1 任取[29.69°,61.04°,84.29°,106.57°]為真實的入射方向,四個入射波的功率都是1 W。利用解模糊算法求各方向上的功率,獨立隨機實驗50次,求得功率的均值及均方根誤差(RMSE)列于表1。該實驗的目的是考察本算法能否達到解模糊的目的。
實驗2 在實驗1的條件下,利用Abramovich的線性規劃(LP)方法估計各波達方向的功率,重復50次,取均值及RMSE列于表1。
由表1中數據可知,當模糊方向沒有真實入射波時,估得的功率也很小,比真實入射波方向的估計功率小兩個數量級以上,因此可認為模糊方向上的功率是由噪聲及擬合算法中的迭代允許誤差產生的。采用協方差陣擬合的功率估計方法,可以將產生模糊的DOA辨識出來。相對應的是,LP方法在本例中失效,在更多更大規模的陣列解模糊實驗中,LP方法同樣失效,經分析是由于線性規劃在大規模優化問題應用中的敏感性造成的。
2)面陣
采用文獻[2]中的例子,陣元坐標矢量為
rx=[-2.5,-2,-2.5,2.5,2,2.5]
ry=[1.5,0,-1.5,1.5,0,-1.5]
任取該陣一個秩為4模糊集為
[(30°,28.66°),(30°,65.44°),(210°,87.35°),
(210°,59.46°),(210°,14.06°)]
也就是說,在該模糊集中任取四個方向作為真實來波方向,那么子空間算法一定會在余下的那個角度上產生一個虛假的峰值。利用解模糊算法可判定這個虛假DOA。
實驗3 為考察本算法能否解面陣模糊,任取[(30°,28.66°),(30°,65.44°),(210°,87.35°),(210°,59.46°)]作為真實入射波方向,設定入射波功率皆為1 W。將由本算法估出功率的均值及均方誤差值列于表2。
實驗4 在實驗3的條件下利用LP方法解模糊,發現此方法仍無效。由此方法估計出的功率均值及RMSE列于表2中的最后兩行。
考察表2中的實驗數據可以得出與線陣實驗同樣的結論:在本例中,基于協方差陣擬合方法有效而LP方法失效。
5 結束語
本文中提出用協方差陣擬合的方法來解決子空間測向算法的流形模糊問題,仿真結果表明:該方法相對于Abramovich提出的線性規劃方法適用于更多陣列,但是本方法所解的模糊僅限于Manikas模糊集中的模糊角度,當模糊角度多于陣元個數時(如文獻[9]中式(70)所示例子),由于不可辨識,無法利用本方法來解模糊。
參考文獻:
[1]MANIKAS A, PROUKAKIS C. Modeling and estimation of ambiguities in linear arrays[J]. IEEE Trans on Signal Processing,1998, 46(8):2166-2179.
[2]MANIKAS A, PROUKAKIS C, LEFKADITIS V. Investigative study of planar array ambiguities based on “hyperhelical” parameterization [J]. IEEE Trans on Signal Processing, 1999, 47(6): 1532-1541.
[3]PROUKAKIS C, MANIKAS A. Study of ambiguities of linear arrays[C]//Proc of IEEE ICASSP’94.[S.l.]:IEEE Press, 1994:549-552.
[4]劉洪盛,肖先賜.一種基于最小流形長度的高精度線陣設計[J]. 航空學報,2008,29(2):462-466.
[5]TAN C M, FOO S E, BEACH M A, et al. Ambiguity in MUSIC and ESPRIT for direction of arrival estimation[J].Electronic Letter, 2002,38(24): 1598-1600.
[6]司偉建.MUSIC算法多值模糊問題研究[J].系統工程與電子技術,2004,26(7):960-962.
[7]司偉建,崔冬槐,司錫才.一種新的解模糊方法研究[J].制導與引信,2007,28(1):44-47.
[8]SCHMIDT R O. Multiple emitter location and signal parameter estimation[J].IEEE Trans on Antennas and Propagation,1986,34(3):276-280.
[9]PILLAI S,BAR-NESS Y, HABER F. A new approach to array geo-metry for improved spatial spectrum estimation[C]//Proc of IEEE International Conference on Acoustics, Speech, and Signal Processing. 1985:1816-1819.
[10]ABRAMOVICH Y I, SPENCER N K, GOROKHOV A Y. Resolving manifold ambiguities in direction-of-arrival estimation for nonuniform linear antenna arrays[J].IEEE Trans on Singal Processing, 1999,47(10):2629-2643.
[11]D’ASSUMPCAO H A. Some new signal processors for arrays of sensors[J].IEEE Trans on Information Theory, 1980,26(4):441-453.