張辰銳 王劍書 雷 濤
1(呂梁學院物理系 山西 呂梁 033000) 2(西北工業大學電子信息工程學院 陜西 西安 710129) 3(陜西科技大學電氣與信息工程學院 陜西 西安710021) 4(蘭州交通大學電子與信息工程學院 甘肅 蘭州 730070)
波達方向估計(DOA)是陣列信號處理的一項重要課題,并且廣泛應用于雷達、聲吶和通信等領域[1-3]。從多重信號分類(MUSIC)算法[4]提出開始,國內外學者對平穩信號的子空間高分辨DOA算法進行了深入的研究。本文的興趣在于準平穩信號QSS(Quasi-stationary signals)的DOA估計算法研究。準平穩信號是統計量在一定短時間內保持不變但會隨著時間片段變化的一類信號[5],例如:語音、音樂和一些機械信號通常被認為是準平穩信號。機場鳥聲監測,麥克風陣列和地面聲目標監測等環境都有QSS的DOA估計的實際應用,從而受到廣大學者的重視。
2010年,Ma等[5]提出一種針對QSS的Khatri-Rao(KR)子空間的方法,該方法利用QSS信號的短時平穩的特點,構造了新的數據模型,獲得了比普通MUSIC算法更高的角度分辨率。方法中提出的協陣列(co-array)的概念、對于非均勻線陣nested array[6-7]、co-preme array[8-9]等以及基于高階累積量(HOC)的DOA估計[7]都具有重要的意義。然而,當信號源相干時,源信號的協方差矩陣不再為對角矩陣,KR子空間模型的構造將失效。對于相干源的DOA估計問題,前人也進行了深入的研究。20世紀80年代,Shan等[10]提出了空間平滑(SS)方法,通過縮小均勻線陣有效陣列孔徑,解決了相干源導致的協方差矩陣秩缺失的問題。前后向空間平滑(FBSS)的發展擴大了SS方法的陣列孔徑,提升了自由度[11-12]。2005年,Han等[13]提出一種通過構造Toeplitz矩陣并利用類旋轉不變技術(ESPIRIT)的方法求解的DOA估計方法。2007年,Ye等[14]對此進行改進,重新構造了矩陣模型,提升了自由度。2010年,Choi[15]通過特征向量方法(EVM)和相關向量方法(CVM)分別構造前向后向向量,并使用類ESPIRIT方法求解DOA,其中EVM方法在全相干源環境下具有很好的效果。上述相干源DOA估計方法一般假設源信號為平穩信號,噪聲為不相關的高斯白噪聲,實際環境中,噪聲可能包含部分干擾源信號,可能具有一定的相關性,這將導致上述相干源DOA估計算法的性能下降。
本文針對準平穩相干源和相關噪聲的環境,提出一種基于KR子空間和前后向空間平滑(FBSS)的DOA估計方法。該方法首先利用接收信號各幀的協方差矩陣構造Toeplitz矩陣,然后將這些矩陣利用KR子空間方法構造一個新的模型,最后利用FBSS方法進行DOA估計。同時,該方法將非相干子空間處理(ISSM)方法[15]進行了擴展,可進行寬帶陣列信號模型的準平穩相干信號的DOA估計。仿真實驗表明,本文方法具有較高的檢測概率和較低的均方根誤差(RMSE),同時也驗證了本文方法對于寬帶陣列信號模型的有效性。
考慮由2N+1個各向同性的陣元組成的均勻線陣,其陣元位置記為-N,-N+1,…,0,…,N。假設存在K個遠場窄帶信號源入射到上述陣列,源信號記為s(t)=[s1(t),s2(t),…,sK(t)]T,陣元接收到的信號記為x(t)=[x-N(t),x-N+1(t),…,xN(t)]T,陣元接收的噪聲為v(t)=[v-N(t),v-N+1(t),…,vN(t)]T,則接收信號模型為:
x(t)=As(t)+v(t)
(1)
式中:A=[a(θ1),a(θ2),…,a(θK)]為陣列流型矩陣,a(θk)=[e-j2π(-N)d sinθk/λ,e-j2π(-N+1)d sinθk/λ,…,e-j2πNd sinθk/λ]T,θk為第k個信號的來向,λ為源信號的波長,d為均勻線陣的陣元間距,且等于源信號的半波長。
這里對源信號與噪聲作如下假設:
1) 假設源信號的數目已知,且每個信號的來向不一樣。
2) 假設源信號為準平穩信號,且各個源信號之間可能相干。
3) 假設噪聲信號為平穩信號,與各個源信號不相關,但各陣元接收的噪聲信號可能具有相關性。
由于源信號為準平穩信號且相干,噪聲相關,傳統的相干源DOA估計算法的性能可能會下降甚至失效。
對接收信號x(t)按采樣后的數據長度L進行分幀,認為幀內信號為平穩信號,則第m幀的信號的協方差矩陣為:
Rx,m=E{xm(t)xm(t)H}=ARs,mAH+C
(2)
式中:Rs,m為第m幀內源信號的協方差矩陣,C為噪聲信號的協方差矩陣。使用Rx,m的第n行元素構造Toeplitz矩陣:
(3)
式中:n=-N,-N+1,…,N,則上式可以寫為:
Wn,m=ArDn,mArH+Cn
(4)
式中:Ar=[ar(θ1),ar(θ2),…,ar(θK)]為陣元0,1,…,N組成的子陣列流型,ar(θk) = [1,e-j2πd sinθk/λ,e-j2πNd sinθk/λ]T,Dn,m為矩陣ARs,m的第n行元素構成的對角矩陣,Cn如下:
(5)
式中:C(n1,n2)為C在(n1,n2)位置的元素,n1,n2=-N,-N+1,…,N。
將協方差矩陣各行構造的Toeplitz矩陣累加起來,構造成新的Toeplitz矩陣:
(6)

通過KR子空間和FBSS方法進行DOA估計。首先對Wm向量化:
ym
(7)
式中:vec(·)為按列向量化運算,⊙為KR乘積[5]運算,dm=diag(Dm)。令Y[y1,y2,…,yM],則:
(8)
式中:Ψ=[d1,d2,…,dM],1M=[1,…,1]T∈M。
對于準平穩信號,一個普遍合理的假設:矩陣[Ψ,1M]∈M×(K+1)是列滿秩的[5]。

(9)
令Φ則:
Z=BΦ
(10)
式(10)與信號模型式(1)有類似的結構,將Φ的每一行視作一路源信號,將B視作陣列流型,則Z為該陣列流型下的接收信號。式(10)通過計算消除了噪聲項,所以該模型不僅能處理不相關的高斯白噪聲,而且對未知的相關噪聲也具有很好的抑制效果。由于B其行有冗余,則Z的行也有冗余。消除Z中重復的行,并將各行按陣列流型A排列可以得到:

(11)


式中:E=diag[e-j2πd sinθ1/λ,e-j2πd sinθ2/λ,…,e-j2πd sinθK/λ]。則有:
(12)


準平穩信號在現實中多為寬帶信號,這里提出寬帶陣列信號模型的DOA方法。第n個陣元接收的信號用各個源信號的時延混合表示:

t=0,1,2,…
(13)
(14)
式中:NSTFT為短時傅里葉變換的窗口長度,f∈[-1/2,1/2]為歸一化頻率。對xn(t)進行STFT,并令


(15)


(16)
式中:f1和f2分別為處理的最小歸一化頻率和最大歸一化頻率。搜索空間譜函數式(15)譜峰位置則可得寬帶源信號的來向。
使用檢測概率和均方根誤差(RMSE)來對算法進行評價。定義DOA的“正確檢測”:一次實驗中檢測到的所有的DOA值分別與其對應的真實DOA值誤差不超過2°。檢測概率定義為:蒙特卡洛實驗中正確檢測的實驗次數與總實驗次數的比值。RMSE定義為蒙特卡洛實驗中正確檢測實驗的DOA估計的均方根誤差:
(17)

實驗一:信噪比(SNR)與噪聲相關系數對算法的影響。使用9元均勻線陣,N=4,源信號使用準平穩的拉普拉斯分布的信號,來向分別為-40°、-20°、0°、15°、30°、50°。各陣元添加的噪聲為具有相關性的高斯白噪聲。源信號每一幀長度為150,幀數為150,前5個源信號互為相干信號,第6個信號為獨立信號。對于各信噪比與噪聲相關系數分別進行500次蒙特卡洛實驗,本文方法(記為KR-FBSS)與FBSS方法、改進的Toeplitz結構方法[13](記為TBD)、EVM和CVM方法的實驗結果如圖1所示。從圖1可以看出,在各組相關噪聲條件下,本文方法具有最好的檢測概率,隨著噪聲相關系數增加,對比方法的檢測概率下降較明顯,而本文方法檢測概率下降并不大。

(a) 噪聲相關系數0.01

(b) 噪聲相關系數0.05

(c) 噪聲相關系數0.1

(d) 噪聲相關系數0.2 圖1 不同信噪比和噪聲相關系數下檢測概率
實驗二:相近的DOA值的估計性能。與實驗一使用一樣的陣列結構,使用兩個相干的準平穩拉普拉斯分布的信號源,幀數和幀長度均為150,來向分別為3.2°和8.6°。各陣元添加噪聲為相關的高斯白噪聲,相關系數為0.1。在各信噪比條件下分別進行500次蒙特卡洛實驗,實驗結果如圖2所示??梢钥闯?,本文算法具有更高的檢測概率和更好的RMSE表現。

(a) 不同SNR下的檢測概率

(b) 不同SNR下的RMSE 圖2 相近的DOA值估計性能
實驗三:分幀長度與幀數對DOA估計的影響。與實驗二使用一樣的陣列結構和源信號,在各陣元添加相關系數0.1的高斯白噪聲噪聲,SNR為6 dB,在不同相關系數和不同幀長度與不同幀數條件下分別進行500次蒙特卡洛實驗,實驗結果如圖3所示??梢钥闯?,隨著幀長度和幀數增加,RMSE更小,為保證RMSE盡量小,幀長度和幀數目均不宜太小,同時也可以看出噪聲相關系數對本文算法影響較小。

(a) 幀長度變化,幀數為150

(b) 幀數變化,幀長度為150 圖3 幀長度和幀數對RMSE的影響
陣列與4.1節一致,使用4路相干的語音信號作為源信號,來向分別為-40°、-20°、0°、30°。源信號采樣頻率為8 000 Hz,各陣元添加信噪比為6 dB的高斯白噪聲,相關系數為0.1。對比的方法為相干子空間處理方法[16](CSSM)和文獻[17]的方法(CSSM-PM)。選取的頻率范圍為500~2 000 Hz,并離散為25個頻帶,兩種方法使用的預估計DOA值都采用真實DOA值,聚焦矩陣的構造均使用RSS方法[16]??臻g譜實驗結果見圖4,可以看出本文的方法能分辨所有源信號的DOA值,CSSM-PM方法只能勉強分辨源信號的DOA,并且存在偽峰,CSSM方法則不能分辨。其中CSSM-PM方法要求相關噪聲的協方差矩陣具有Toeplitz結構,然而這一條件多數情況下也很難滿足。

圖4 寬帶DOA估計空間譜
在準平穩信號源和相關噪聲環境下,本文研究了采用Khatri-Rao(KR)子空間和前后向空間平滑理論的波達方向估計方法。與常見前后向空間平滑方法及Khatri-Rao子空間等方法相比,具有更高的檢測概率和更低的均方根誤差。同時,由于實際環境中的噪聲可能包含部分干擾源信號,可能具有一定的相關性。本文算法通過利用KR子空間方法將各幀Toeplitz矩陣組合成新的信號模型,因此對相關噪聲有較強的適應性,比常見相關算法表現更好,且在噪聲相關系數變化下均有接近或較好的性能表現。此外,算法的RMSE分別隨著SNR、幀長度和幀數的增加而變小,同時算法對寬帶陣列信號模型也具有一定的有效性。