基于兩次傅里葉變換的時(shí)域MUSIC波達(dá)方向估計(jì)

(江西科技學(xué)院 信息工程學(xué)院，南昌 330098)

1 引 言

寬帶信號(hào)具有較強(qiáng)的抗干擾能力，在聲吶、雷達(dá)和通信等領(lǐng)域被廣泛應(yīng)用和分析[1-2]。到目前為止，已有許多學(xué)者對(duì)基于寬帶信號(hào)的被動(dòng)檢測(cè)和方位估計(jì)問題進(jìn)行了深入研究，并取得了一定的成果[3-5]。近年來，具有高分辨能力的多重信號(hào)分類(Multiple Signal Classification,MUSIC)方法得到了快速發(fā)展，并在聲吶、雷達(dá)和通信等領(lǐng)域得到了廣泛應(yīng)用[6-9]。由于頻域和波束域 MUSIC方法對(duì)數(shù)據(jù)平穩(wěn)性的要求限制了其在被動(dòng)檢測(cè)和方位估計(jì)中的有效應(yīng)用。如何在付出代價(jià)較小的條件下放寬頻域和波束域MUSIC方法對(duì)數(shù)據(jù)平穩(wěn)性的要求，是MUSIC方法應(yīng)用于寬帶目標(biāo)被動(dòng)檢測(cè)和方位估計(jì)需要解決的關(guān)鍵問題[10-11]。

針對(duì)頻域和波束域MUSIC波達(dá)方向估計(jì)方法在快拍數(shù)較少時(shí)難以穩(wěn)定實(shí)現(xiàn)目標(biāo)波達(dá)方向估計(jì)問題，司偉建等人[12]利用延時(shí)相關(guān)函數(shù)所包含的角度信息，并結(jié)合混沌優(yōu)化思想實(shí)現(xiàn)空間譜和波達(dá)方向(Direction of Arrival，DOA)估計(jì)；張志剛等人[13]利用觀測(cè)矩陣的結(jié)構(gòu)信息，通過自適應(yīng)迭代加權(quán)，在秩缺失情況下，提高了DOA估計(jì)精度；文獻(xiàn)[14-16]采用了壓縮感知與MUSIC相結(jié)合，提高了MUSIC算法在快拍數(shù)較少情況下的DOA估計(jì)性能；張濤濤等人[17]通過構(gòu)造基于最佳準(zhǔn)則的聚焦矩陣，然后利用滑動(dòng)平均實(shí)現(xiàn)對(duì)寬帶MUSIC方法改進(jìn)。以上方法雖然能在快拍數(shù)不能滿足MUSIC方法要求時(shí)提高波達(dá)方向估計(jì)性能，但均沒有給出如何只利用時(shí)頻處理技術(shù)在快拍數(shù)不能滿足MUSIC方法要求時(shí)提高其對(duì)DOA的估計(jì)性能。為此，本文結(jié)合時(shí)域常規(guī)最佳處理和MUSIC基本思想，提出了一種基于兩次傅里葉變換的時(shí)域MUSIC波達(dá)方向估計(jì)方法(本文稱之為FTMUSIC方法)。該方法是基于時(shí)域與頻域相結(jié)合的陣列信號(hào)處理方法，在一次有效快拍數(shù)情況下，通過多個(gè)時(shí)間點(diǎn)的累積獲得良好的協(xié)方差矩陣估計(jì)，獲得與導(dǎo)向權(quán)向量穩(wěn)定正交的噪聲子空間，降低快拍數(shù)對(duì)MUSIC方法的影響，且具有較好的噪聲抑制效果，進(jìn)一步拓寬了MUSIC方法在工程領(lǐng)域的應(yīng)用。

2 MUSIC方法

MUSIC波達(dá)方向估計(jì)方法是在噪聲(包括干擾)與信號(hào)不相關(guān)的假設(shè)下，對(duì)協(xié)方差矩陣Rx作特征分解得到信號(hào)特征向量Us和噪聲特征向量Uv。依據(jù)噪聲向量與導(dǎo)向權(quán)向量的正交性[10-11]，可獲得來自掃描角度θ上的空間譜為

(1)

式中：W(θ)=[ej2πfτ1,ej2πfτ2,…,ej2πfτN]為導(dǎo)向權(quán)向量;τn=(n-1)dcos(θ)/c，c為聲速，d為陣間距。

從MUSIC方法輸出空間譜過程可知，MUSIC方法實(shí)現(xiàn)時(shí)應(yīng)該滿足：信號(hào)與噪聲不相關(guān)，即Rx中不應(yīng)含有信號(hào)與噪聲相關(guān)成分的貢獻(xiàn)；數(shù)學(xué)上要求Rx滿秩；權(quán)向量應(yīng)為復(fù)權(quán)，以便在目標(biāo)波達(dá)方向獲得與噪聲子空間的正交性。

(2)

3 FTMUSIC方法

3.1 理論分析

對(duì)于時(shí)域數(shù)據(jù)，通常一次快拍即可獲得具有信號(hào)與噪聲不相關(guān)和滿秩條件的協(xié)方差矩陣[18]。若能引入復(fù)導(dǎo)向權(quán)向量，則可在單次快拍條件下實(shí)現(xiàn)DOA估計(jì)，進(jìn)而降低有效快拍數(shù)對(duì)MUSIC方法估計(jì)性能的影響。FTMUSIC方法的基本思路為：通過兩次快速傅里葉變換將各陣元時(shí)域?qū)崝?shù)據(jù)轉(zhuǎn)換為相移后的復(fù)數(shù)解析數(shù)據(jù)，用復(fù)數(shù)數(shù)據(jù)引入復(fù)導(dǎo)向權(quán)向量，再結(jié)合時(shí)域常規(guī)最佳處理求取經(jīng)相移后協(xié)方差矩陣，最后利用特征分解思想求取具有正交特性的噪聲子空間，獲得來波方向估計(jì)值。

假設(shè)單個(gè)目標(biāo)信號(hào)從θ1方向輻射到N元的陣列上，則第n號(hào)陣元接收的實(shí)數(shù)信號(hào)可表示為

xn(t)=s(t-τn(θ1))+vn(t) 。

(3)

式中：τn(θ1)，1≤n≤N為目標(biāo)信號(hào)到達(dá)第n號(hào)陣元相對(duì)于參考陣元的時(shí)間延遲，其只與目標(biāo)信號(hào)相對(duì)陣列所處的方位θ1有關(guān)。

由FFT變換的性質(zhì)可將式(3)寫為頻域離散化形式，即

Xn(fm)=S(fm)·exp(-j2πfmτn(θ1)+Vn(fm)。

(4)

由式(4)所示的頻域形式也可發(fā)現(xiàn)，對(duì)陣列接收時(shí)域復(fù)解析數(shù)據(jù)進(jìn)行的時(shí)延τn可在頻域通過相移的方式來實(shí)現(xiàn)，即插入相移因子exp(j2πfmτn(θ))改變各陣元接收數(shù)據(jù)在掃描角度θ上的相位不一致性。此時(shí)，可將式(4)用矩陣形式表達(dá)：

XF(fm)=S(fm)·exp(-jφ1)+V(fm) 。

(5)

式中：φ1=[2πfmτ1(θ1)，2πfmτ2(θ1)，…，2πfmτN(θ1)]。

為了求取掃描角度θ方向的空間譜，此時(shí)可對(duì)處理頻帶內(nèi)的頻域數(shù)據(jù)按頻率單元乘以相移因子exp(jφ)，改變各陣元接收數(shù)據(jù)在掃描角度θ上的相位差異，φ=[2πfmτ1(θ),2πfmτ2(θ),…,2πfmτN(θ)]T。此時(shí)式(5)可轉(zhuǎn)換為

YF(fm)=XF(fm)·exp(jφ)=XF(fm)·

[exp(j2πfmτ1(θ)),…,exp(j2πfmτN(θ))]T=

XF(fm)Wτ,fl≤fm≤fh。

(6)

式中：fs為系統(tǒng)采樣頻率，fl為處理頻帶下限，fh為處理頻帶上限，Wτ=exp(jφ)。

為了構(gòu)造時(shí)域復(fù)解析數(shù)據(jù)，接下來對(duì)式(6)如下處理：

(7)

(8)

(9)

式中：I1×N=[1,1,…,1]1×N為加法器。

(10)

與MUSIC方法相比，F(xiàn)TMUSIC方法所需快拍次數(shù)較少。在數(shù)據(jù)長度滿足式(10)情況下，只需一次快拍(包含L采樣點(diǎn)數(shù))即可穩(wěn)定得到具有與目標(biāo)波達(dá)方向?qū)?yīng)導(dǎo)向權(quán)向量正交特性的噪聲子空間。

由于引入了時(shí)域解析數(shù)據(jù)，F(xiàn)TMUSIC方法放寬了對(duì)快拍數(shù)的條件限制，且不需做子帶分解，計(jì)算量大大減小，使得FTMUSIC方法有更寬的適用范圍。不過，在一個(gè)掃描角度上，F(xiàn)TMUSIC方法只有一組子空間，而MUSIC方法相應(yīng)于每個(gè)子帶均有一組子空間，子空間個(gè)數(shù)數(shù)倍于FTMUSIC方法。因此，需要數(shù)值仿真和實(shí)測(cè)數(shù)據(jù)進(jìn)一步對(duì)比FTMUSIC方法和MUSIC方法的噪聲抑制及抗干擾能力。

3.2 實(shí)現(xiàn)流程

根據(jù)上面分析可知，基于FTMUSIC方法的DOA估計(jì)方法可分為如下步驟實(shí)現(xiàn)：

Step1 對(duì)N元陣列采集的時(shí)域離散數(shù)據(jù)xn(m)進(jìn)行濾波處理，并對(duì)濾波后的N路離散數(shù)據(jù)做一次M點(diǎn)長度的FFT處理，得到相應(yīng)的頻域數(shù)據(jù)。

Step2 在掃描角度θ上，按式(7)對(duì)處理頻帶內(nèi)的頻域數(shù)據(jù)乘以相移因子exp(j2πfmτn(θ))。

Step3 對(duì)相移補(bǔ)償后的N路頻域數(shù)據(jù)的共軛形式做FFT，得到時(shí)域復(fù)解析數(shù)據(jù)，此時(shí)的復(fù)解析數(shù)據(jù)為經(jīng)過相移后的時(shí)域復(fù)解析數(shù)據(jù)。

Step5 按式(9)對(duì)噪聲特征向量進(jìn)行處理，可獲得FTMUSIC方法在掃描角度θ上的空間譜估計(jì)值PFTMUSIC(θ)。

Step6 改變掃描角度θ，重復(fù)Step 2～ 5，可獲得FTMUSIC方法在不同掃描角度上的空間譜PFTMUSIC(θ),θ∈[θl,θh]，θl為掃描角度下限，θh為掃描角度上限。

Step7 通過對(duì)PFTMUSIC(θ),θ∈[θl,θh]進(jìn)行峰值篩選，可對(duì)目標(biāo)波達(dá)方向?qū)崿F(xiàn)估計(jì)。

3.3 數(shù)值仿真分析

為了進(jìn)一步驗(yàn)證FTMUSIC方法對(duì)有效快拍數(shù)的寬容性，對(duì)FTMUSIC方法進(jìn)行數(shù)值仿真分析。仿真條件：接收陣為32陣元的等間隔水平直線陣，相鄰陣元間距為1 m，模擬信號(hào)為700～800 Hz寬帶高斯白噪聲，信號(hào)長度0.1 s，SLR(Spectrum Level Ratio)為信號(hào)與背景噪聲譜級(jí)比，背景噪聲為加性高斯白噪聲，模擬信號(hào)輸入方向?yàn)?0°。系統(tǒng)采樣率為20 kHz，一次采樣長度為1 s，F(xiàn)TMUSIC方法由一次快拍實(shí)現(xiàn)，MUSIC方法所分子帶數(shù)為M=100，每一子帶由76次快拍(每一次快拍包含512個(gè)采樣樣本)實(shí)現(xiàn)。

圖1為SLR=-25～0 dB情況下，由MUSIC方法與FTMUSIC方法通過200次獨(dú)立統(tǒng)計(jì)所得正確檢測(cè)概率。圖2為SLR=-25～0 dB情況下，由MUSIC方法與FTMUSIC方法通過200次獨(dú)立統(tǒng)計(jì)所得DOA估計(jì)均方誤差(Root Mean Squared Error,RMSE)。圖3為SLR=-15 dB情況下，由MUSIC方法與FTMUSIC方法所得單一時(shí)刻空間譜。

圖1 MUSIC與FTMUSIC檢測(cè)目標(biāo)成功率Fig.1 Detect target success rate of MUSIC and FTMUSIC

圖2 MUSIC與FTMUSIC波達(dá)方向估計(jì)均方根誤差Fig.2 RMSE of DOA estimation of MUSIC and FTMUSIC

圖3 2種方法輸出空間譜

由圖2和圖3仿真結(jié)果可知，在有效快拍數(shù)較少情況下，對(duì)于相同的正確檢測(cè)概率和DOA估計(jì)精度，相比MUSIC方法，F(xiàn)TMUSIC方法輸出空間譜對(duì)最低信噪比要求上降低了5 dB，背景噪聲級(jí)和旁瓣級(jí)得到了3 dB以上的改善，數(shù)值仿真結(jié)果進(jìn)一步證明了FTMUSIC方法對(duì)有效快拍數(shù)的寬容性得到了有效改善。

在該數(shù)值仿真數(shù)據(jù)處理中，F(xiàn)TMUSIC方法與頻域MUSIC方法計(jì)算量比較是在Intel(R) Core(TM) i7-7500U CPU@2.70 GHz 2.90 GHz的計(jì)算機(jī)上利用Matlab2014a的CPU TIME測(cè)出的，如表1所示。該實(shí)驗(yàn)采樣率較高，由于FTMUSIC方法只對(duì)一次快拍數(shù)據(jù)進(jìn)行處理，且并未進(jìn)行子帶分解，其計(jì)算量小于頻域MUSIC方法。

表1 兩種方法計(jì)算時(shí)間比較Tab.1 Comparison of processing time of two methods

4 實(shí)測(cè)數(shù)據(jù)處理

為進(jìn)一步考核FTMUSIC方法性能，下面分別利用FTMUSIC方法和MUSIC方法對(duì)某次試驗(yàn)數(shù)據(jù)進(jìn)行處理并作對(duì)比。試驗(yàn)中，接收陣為均勻分布32陣元的水平直線陣，陣元間距為8 m；處理數(shù)據(jù)帶寬為60～120 Hz，處理數(shù)據(jù)時(shí)間長度為100 s，該時(shí)間段內(nèi)目標(biāo)方位在33°、50°、55°、69°、77°(運(yùn)動(dòng))、82°、109°和153°(運(yùn)動(dòng))附近。MUSIC方法具體處理過程如下：首先將數(shù)據(jù)分塊，每塊數(shù)據(jù)為1 024個(gè)，記為一次快拍樣本長度，數(shù)據(jù)塊之間重疊512個(gè)，快拍數(shù)為76，做FFT后選取60～120 Hz頻段，對(duì)每個(gè)頻點(diǎn)分別作協(xié)方差估計(jì)；而FTMUSIC方法通過512階帶通濾波器選取60～120 Hz頻段，采用一次快拍數(shù)據(jù)進(jìn)行分析處理。圖4和圖5分別為MUSIC方法和FTMUSIC方法處理數(shù)據(jù)所得時(shí)間歷程圖，圖6為單一時(shí)刻兩種方法輸出空間譜。

圖4 MUSIC方法所得方位歷程圖Fig.4 Bearing time/record of MUSIC method

圖5 FTMUSIC方法所得方位歷程圖Fig.5 Bearing time/record of FTMUSIC method

圖6 兩種方法輸出空間譜(t=60 s)Fig.6 Output spatial spectrum of two methods(t=60 s)

由圖4～5可知，F(xiàn)TMUSIC方法所得方位歷程圖顯示目標(biāo)航跡清晰，目標(biāo)方位明晰可辨，進(jìn)一步證明了FTMUSIC方法在理論上和實(shí)測(cè)數(shù)據(jù)處理上的正確性；而MUSIC方法在0～100 s時(shí)間段內(nèi)無法對(duì)55°附近目標(biāo)實(shí)現(xiàn)有效檢測(cè)，對(duì)77°附近目標(biāo)檢測(cè)效果較差，且無目標(biāo)處譜級(jí)較大。

同樣，由圖6也可以看出，F(xiàn)TMUSIC方法較MUSIC方法有更低的背景噪聲級(jí)和旁瓣級(jí)。

5 結(jié)束語

本文首先介紹了頻域MUSIC波達(dá)方向估計(jì)方法基本原理，并分析了頻域MUSIC方法對(duì)快拍數(shù)的需求，然后提出了一種基于兩次傅里變換的時(shí)域MUSIC波達(dá)方向估計(jì)方法——FTMUSIC方法。該方法是基于時(shí)域與頻域相結(jié)合的陣列信號(hào)處理方法，利用兩次傅里葉變換實(shí)現(xiàn)了時(shí)域復(fù)解析數(shù)據(jù)構(gòu)造和相移補(bǔ)償，在一次有效快拍數(shù)情況下通過多個(gè)時(shí)間點(diǎn)的累積可獲得良好的協(xié)方差矩陣估計(jì)，可獲得與導(dǎo)向權(quán)向量穩(wěn)定正交的噪聲子空間，降低了快拍數(shù)對(duì)MUSIC方法的影響，且具有較好的噪聲抑制效果，進(jìn)一步拓寬了MUSIC方法在工程領(lǐng)域的應(yīng)用。數(shù)值仿真及實(shí)測(cè)數(shù)據(jù)處理結(jié)果均表明，在有效快拍數(shù)較少情況下，相比頻域MUSIC方法，對(duì)于同樣的正確檢測(cè)概率和DOA估計(jì)精度，F(xiàn)TMUSIC方法輸出空間譜對(duì)最低信噪比要求上降低了5 dB，背景噪聲級(jí)和旁瓣級(jí)得到了3 dB以上的改善，目標(biāo)檢測(cè)和方位估計(jì)性能得到了明顯提高。

如何在FTMUSIC方法基礎(chǔ)上進(jìn)一步降低一次快拍包含的采樣點(diǎn)數(shù)還有待進(jìn)一步研究。