王 超,魏璽章
(國防科技大學電子科學與工程學院,湖南長沙410073)
空間信源定位是陣列信號處理中的一項重要研究內容,一直受到國內外學者的密切關注。早期的空間信源定位研究主要集中在遠場源的DOA估計上面,然而當信源靠近陣列而處于陣列的Fresnel區(qū)(近場區(qū))[1]時,波前的固有彎曲將不能忽略,即遠場平面波前的假設將不再成立,需要用球面波來精確描述,其形狀隨陣元位置具有非線性變化特性,可以通過信源的距離及其DOA來聯(lián)合確定,當信源的頻率未知時,近場源參數估計就成了頻率、距離和DOA的聯(lián)合估計問題。近場源定位在語音增強、聲源定位、雷達[2]、聲納、電子監(jiān)測和地震監(jiān)測等領域中具有廣泛的應用。
近年來,人們提出了許多近場源距離和DOA等多參數聯(lián)合估計問題[3-12],文獻[3]最早提出了最大似然估計方法,該方法具有最優(yōu)的估計性能,但需要多維搜索并且是高度非線性的。Huang等人[4]和Russell等人[5]把傳統(tǒng)的一維MUSIC算法推廣到二維MUSIC算法,實現(xiàn)距離和DOA的聯(lián)合估計,但是二維MUSIC算法仍然要作復雜度很大的二維搜索。Starer等人[6]提出的路徑跟蹤(path-following)進一步降低了二維MUSIC的運算量。為避免多維搜索,Challa和Shamsunder[7]提出了一種基于四階累積量矩陣的類ESPRIT方法,該方法理論上不需譜峰搜索且各參數自動配對,但由于二維參數分別由獨立的特征分解得到,所以實際上仍需要參數配對,而且陣元利用率較低。基于四階累積量的改進的MUSIC算法[8]也被提出,利用高階累積量帶來的虛擬陣元擴展和抑制色噪聲的作用,提高了傳統(tǒng)MUSIC算法的效率,但仍然需要非線性的多維搜索。Emmanuele和Karim[9]提出了一種加權線性預測方法,需要解多個約束最優(yōu)化方程和參數配對問題。陳建峰[10]、黃家才[11]和梁軍利[12]等人提出了基于高階累積量的頻率、距離和DOA聯(lián)合估計方法,通過選擇特定序號陣元輸出來構造多個特定結構的高階累積量,并利用不同累積量矩陣的結構關系來獲得待估計的參數,算法由于進行了多次特征分解導致需要進行參數匹配,或者是利用了線性陣列的對稱性而造成了孔徑損失。周祎[13]和胡增輝[14]利用聯(lián)合對角化在假設信源頻率已知時估計近場源距離和DOA時取得了一定的效果。
上述分析顯示在近場源定位場景中存在的問題主要集中在以下三個方面:(1)陣列孔徑的損失;(2)估計參數的配對;(3)需要計算量較大的搜索或者迭代運算;(4)信號源頻率、DOA、距離等多參數的聯(lián)合估計。為了克服上述缺點,本文提出了一種近場源信號頻率、到達角和距離的聯(lián)合估計算法,將聯(lián)合對角化方法用于近場源頻率、距離和DOA估計,首先基于陣列接收數據的協(xié)方差矩陣構造白化矩陣,再利用一組累積量矩陣的對角結構信息和聯(lián)合對角化矩陣得到陣列流形和信號源頻率的估計,進而由陣列導向矢量聯(lián)合對應的信號源頻率估計信號源的到達方向和距離。與基于四階累積量的類ESPRIT方法相比,本文所提算法能夠聯(lián)合估計頻率、到達角和距離,并且提高了陣元利用率,避免了參數匹配操作。
假定有M個近場窄帶信號源入射到間距為d的均勻線陣(如圖1)上,且均勻線陣由2N+1個陣元構成,從左至右陣元編號為-N,-N+1,…,0,…,N-1,N,中心陣元為編號為0的參考陣元,源信號經過解調到中頻并抽樣后,第m個源信號可表示為S m(t)=s m(t)ejωmt,s m(t)表示信號的復
包絡,ωm為信號頻率。故第l個陣元上的接收信號可表示為

式中,s m(t)=A mejφm,其中A m表示信源m的信號幅度,φm表示信源m的初始相位,n l(t)是陣元l接收到的加性復高斯白噪聲,用τlm表示第m個信源在陣元l和參考陣元0之間的相位差,可以表示為

式中,λm,θm,r m分別為信源m的波長、DOA和距離;d為均勻陣列陣元間距;f m=c/λm為信源m的中心頻率。如圖1所示,為避免模糊假設信源DOA滿足-π≤θm≤π/2。

圖1 均勻線陣近場源結構示意圖
利用近場源近似即Fresnel近似(本質為Taylor二階展開),可將式(2)近似為

假定M個信源可分別表示為s1(t)ejω1t,s2(t)ejω2t,…,s M(t)ejωMt。 式(1)可用矩陣形式 表示為

式中,

對信號模型作如下假設:
(1)源信號是相互統(tǒng)計獨立并且具有非零峰度的非高斯過程。
(2)每個陣元上的加性噪聲均是零均值的復高斯白噪聲,且噪聲與源信號以及各陣元噪聲之間是不相關的。
(3)陣元間距d≤λm/2;線陣陣元個數N S滿足N S>max[M,2]。
(4)基于具有2N+1個陣元的均勻對稱線陣的情況下,與文獻[7]中要求最大可容納信源數M<N不同,本文只要求M<2N+1即可。
(5)每個源信號距離和波達角不完全相同。即(θi,r i)≠(θj,r j),對任意1≤i≠j≤M成立。
除上述假設外,源信號數目M的確定屬于信號檢測問題,已有文獻詳細研究該問題[15],因此,文中假設M是已知的。
矩陣A與遠場情形一樣稱為近場信號的陣列流形矩陣,是DOA、距離以及頻率的函數。S是以各個源信號經過采樣得到的離散序列構成的行向量組成的矩陣,x是陣列接收數據矩陣。算法首先對白化后的陣列接收數據構造一組四階累積量,再由聯(lián)合對角化方法估計出陣列流形A和窄帶信號源頻率,最后根據陣列導向矢量a m和對應的源信號頻率f m的結構特點,利用最小二乘估計出DOA和距離。
(1)對混合觀測數據x進行白化處理,構造白化矩陣Q使得QA為酉矩陣,R x表示陣列接收信號的協(xié)方差矩陣,由假設中噪聲與信號不相關,協(xié)方差矩陣R x可表示為

對R x進行特征分解可得R x=VΣVH,其中Σ=diag(λ1,…,λM,σ2,…,σ2),且λ1≥…≥λM≥σ2=…=σ2,則根據文獻[16]白化矩陣Q可以取為

根據白化矩陣Q,進而得到白化數據(t)。

由上式,可以構造出M×(M+1)/2個四階累積量矩陣C kl。
(3)令W=QA,白化后的接收數據,由高階累積量的可加性以及高斯分布高階累積量為零的性質,根據四階累積量的定義,C kl可以寫成如下形式:

式中,W(i,k)表示矩陣W第i行第k列對應的元素;γ4sm為源信號s m對應的四階累積量,且

根據文獻[15]可知存在一個酉矩陣W使得C kl可以被同時對角化,用矩陣形式具體可表示如下:

式中,對角矩陣Γkl=diag{W(k,l)W?(l,1)γ4s1,…,W(k,M)W?(l,M)γ4sM}。
(4)得到酉矩陣W,利用近似對角化方法,把所有M×(M+1)/2個矩陣{C kl}聯(lián)合對角化,即可得到源信號的估計(t)WH以及混合矩陣的估計為,Q+為Q的偽逆。
由于聯(lián)合近似對角化算法理論可知,?AΛP,其中,Λ為一對角線元素非零的對角矩陣,P為一置換矩陣,它們分別對應于聯(lián)合正交化算法的幅度和順序不確定性,假設對角矩陣Λ=dia g{z1,…,z M},即令,則存在復常數z m(m=1,…,M),使得

但是這并不影響由陣列導向矢量即陣列流形的列矢量來確定相應信號源的二維參數,并且對角化酉矩陣W與對角矩陣Γkl是在一次特征值分解中得到,因此源信號與得到的陣列流形是行列相互對應的,即每個信號源的采樣離散序列在中的行位置與其對應的陣列導向矢量在中的列位置是順序對應的。對估計出的源信號每一行進行頻譜分析(文中采用離散傅里葉變換)得到M個信號源的頻率{f1,…,f m,…,f M}。
利用已經估計出來的信源頻率{f1,…,f m,…,f M},再用最小二乘方法得到陣列導向矢量中包含的距離和DOA參數,每個信源對應的頻率、距離和DOA都估計出來以后,可以將頻率、距離或DOA按升序或者降序排列,因此實際估計時可以忽略P帶來的影響。
由式(3)、(4)和(6)可以看出,導向矢量a m的相位是陣元位置的二次函數,已假定對角矩陣則

令

式中,~a m(i)表示~a m的第i個元素,所以c m是個2N×1維的矢量。
將式(9)代入式(10)可得

根據第m個源信號到陣元i-N和陣元i-N-1的距離r m,i-N和r m,i-N-1的距離差小于陣元間距(三角形兩邊之差小于第三邊)可以得到|τi-N,m-τi-N-1,m|=2π|r i-N,m-r|/λm<2πd/λm,由假設可知d≤min(λm/2),1<m<M,所以|τi-N,m-τi-N-1,m|<π。
若用φm(i)表示c m(i)的相位,將式(4)代入上式可得

用矩陣可表示為

對上式利用線性最小二乘即可得到ψm,?m的估計。
近場源的到達方向和距離的估計公式如下:

已知源信號頻率,再結合對應的陣列導向矢量即矩陣的列向量可以得到對應信號源的到達角和距離參數。
假定源信號是兩個等功率統(tǒng)計獨立的空間信號可以分別表示為ej2πf1t+jα1和ej2πf2t+jα2,其中α1和α2均為[0,2π]上的均勻分布。入射到九元均勻線陣上,且陣元間距d=min(λ1,λ2)/4。試驗中使用根均方誤差(RMSE)作為性能評價標準:

式中,為參數估計值,a為參數真值,N x為蒙特卡洛實驗次數。
實驗1 頻率、DOA和距離分別為(3 M Hz,30°,2λ1)和(4 M Hz,20°,λ2),控制信噪比(SNR)從10 dB到20 d B,數據長度為4096個快拍。蒙特卡洛仿真次數為300次。兩個信源頻率、DOA及距離的估計結果如表1所示。

表1 不同信噪比下本文所提方法與基于高階累積量的類ESPRIT方法的比較
表1中顯示的是信源各個參數的均方根誤差(RMSE),位于前面的數據是由基于高階累積量的類ESPRIT方法所得,后面的數據則是本文所提算法所得。從仿真結果中可以看出在快拍數固定不變,隨著信噪比的增加,各個參數估計的均方根誤差呈減小趨勢。而且從仿真結果中可以看出本文所提方法的估計效果明顯比基于高階累積量的類ESPRIT方法要好。
從圖2可以看出,在信噪比小于-10 dB時,頻率估計均方根誤差明顯減小;當信噪比大于-10 dB時,均方根誤差變化不再明顯,并逐漸趨于零。在SNR大于-10 dB時,頻率估計精度已經較高且比較穩(wěn)定。這恰好解釋了表1中本文所提算法的頻率估計RMSE基本不隨信噪比變化的原因。
實驗2 頻率、DOA和距離仍分別為(3 MHz,30°,2λ1)和(4 MHz,20°,λ2),固定信噪比為15 dB,快拍數從500個變化到2 500個后,300次蒙特卡洛仿真后兩個信源頻率、DOA及距離三維參數的估計結果如表2所示。

圖2 本文所提方法下頻率的估計結果

表2 不同快拍數下本文所提方法與基于高階累積量的類ESPRIT方法的比較
從圖3可以看出,頻率的估計精度和采樣點數成正比,采樣點數越多,頻率估計精度越好。這是因為頻率是根據估計得到的源信號離散序列經過FFT提取頻譜之后得到的,所以頻率的估計精度和采樣點數成正比。

圖3 頻率估計隨快拍數變化
表2中顯示的是信源各個參數的估計均方根誤差(RMSE),位于前面的數據是由基于高階累積量的類ESPRIT方法所得,后面的數據則是本文所提算法所得。可以看出,頻率、DOA和距離的估計精度隨采樣點數變化的趨勢一致,都是隨著采樣數目的增大而變好,并且在信噪比相同的情況下,同樣快拍數下本文所提算法表現(xiàn)更好。
實驗3 為了直觀表現(xiàn)本文方法的優(yōu)越性,將本文方法與文獻[17]中的基于高階累積量的類ESPRIT方法進行比較,試驗中信源為2個等功率的調幅信號,載波頻率分別為f1=3 M Hz,f2=4 M Hz,帶寬為25 k Hz,采樣頻率fs=12 MHz。兩個信源的到達角分別為30°和20°,距離分別為2λ1和λ2。波速為光速,陣元間距為最小波長的0.25倍。信噪比固定為15 dB,快拍數固定為2 048,做300次Monte-Carlo實驗,兩個信號的頻率、DOA和距離分布如圖4和圖5所示。
由上圖中可以直觀地看出,本文所提方法參數估計統(tǒng)計性能更好。
本文提出了一種基于高階累積量矩陣的聯(lián)合對角化的近場源頻率、DOA和距離三維參數的聯(lián)合估計方法。首先利用聯(lián)合對角化方法得到高階累積量矩陣的對角信息和對角化矩陣來估計陣列流形和各個源信號的離散序列構成的信號矩陣,對信號矩陣進行譜分析得到各獨立信號源的頻率,進而通過總體最小二乘方法估計出DOA和距離,所提算法在基于陣元數為2N+1的均勻對稱線陣時可估計的最大信源數目為2N,與基于高階累積量的類ESPRIT方法的仿真結果比較可以看出,本文所提方法在相同的信噪比和快拍數下表現(xiàn)更好,且從多次蒙特卡洛仿真的結果可以看出所提方法更穩(wěn)健,并且參數自動配對,不需要搜索或者迭代等過程。

圖4 基于四階累積量的類ESPRIT方法

圖5 本文所提方法
[1]Hoole P R P.Smart Antennas and Signal Processing for Communication,Biomedical and Radar Systems[M].Southampton:WIT Press,2001.
[2]郭藝奪,張永順,童寧寧,等.基于L型MIMO雷達二維DOA估計新算法[J].雷達科學與技術,2010,8(5):412-416.GUO Yi-duo,ZHANG Yong-shun,TONG Ningning,et al.New Algorithm for 2-D DOA Estimation Based on L Shaped MIMO Radar[J].Radar Science and Technology,2010,8(5):412-416.(in Chinese)
[3]Swindlehurst A L,Kailath T.Passive Direction of Arrival and Range Estimation for Near-Field Sources[C]∥4th Annual ASSP Workshop on Spectrum Estimation and Modeling,Minneapolis,Minnesota,USA:[s.n.],1988:123-128.
[4]Huang Y D,Barkat M.Near-Field Multiple Sources Localization by Passive Sensor Array[J].IEEE Trans on Antennas and Propagation,1991,39(7):968-975.
[5]Russell J,Bell K L,Van Trees H L.Broadband Passive Range Estimation Using MUSIC[C]∥2002 IEEE International Conference on Acoustics,Speech,and Signal Processing,Orlando,Florida,USA:[s.n.],2002:2920-2922.
[6]Starer D,Nehoraj A.Path-Following Algorithm for Passive Localization of Near-Field Sources[C]∥5th Annual ASSP Workshop on Spectrum Estimation and Modeling,Rochester,NY,USA:[s.n.],1990:322-326.
[7]Challa R N,Shamsunder S.Higher-Order Subspace Based Algorithm for Passive Localization of Near-Field Sources[C]∥29th Asilomar Conference on Signals,Systems and Computers,Pacific Grove,CA:[s.n.],1995:777-781.
[8]刁鳴,吳小強,李曉剛.基于四階累積量的測向方法研究[J].系統(tǒng)工程與電子技術,2008,30(2):226-228.
[9]Grosicki E,Abed-Meraim K,Hua Y.A Weighted Linear Prediction Method for Near-Field Source Localization[J].IEEE Trans on Signal Processing,2005,53(10):3651-3660.
[10]陳建峰,張賢達,吳云韜.近場源距離、頻率及到達角聯(lián)合估計算法[J].電子學報,2004,32(5):803-806.
[11]黃家才,石要武,陶建武.近場源DOA、距離、極化參數及頻率聯(lián)合估計算法[J].計算機工程與應用,2006(19):17-19.
[12]Liang Junli,Zeng Xianju,Ji Bangjie.A Computationally Efficient Algorithm for Joint Range-DOAFrequency Estimation of Near-Field Sources[J].Digital Signal Processing,2009,19(4):596-611.
[13]周祎,馮大政,劉建強.基于聯(lián)合對角化的近場源參數估計[J].電子與信息學報,2006,28(10):1766-1769.
[14]胡增輝,朱炬波,梁甸農.基于盲源分離的近場源參數估計[J].信號處理,2010,26(6):951-955.
[15]Shapo B,Bethel R.A Novel Passive Broadband Bayesian Detector/Tracker[C]∥Proceedings of the 2000 IEEE Sensor Array and Multichannel Signal Processing Workshop,Cambridge,MA,USA:[s.n.],2000:92-96.
[16]Belouchrani A,Abed-Meraim K,Cardoso J F,et al.A Blind Source Separation Technique Using Second-Order Statistics[J].IEEE Trans on Signal Processing,1997,45(2):434-444.
[17]吳建嵐,胡光波,季鋒.基于高階累積量的近場源三維參數估計算法[J].艦船電子工程,2010,30(11):50-53.