黎昕婷，鐘舜聰，鐘劍鋒

（福州大學(xué)機(jī)械工程及自動(dòng)化學(xué)院，福州 350108）

0 概述

基于傳感器陣列的波達(dá)方向（Direction of Arrival，DOA）估計(jì)是聲納、雷達(dá)、通信、語(yǔ)音處理等領(lǐng)域的研究熱點(diǎn)，然而現(xiàn)代化技術(shù)發(fā)展的需求使DOA 估計(jì)不再局限于處理傳統(tǒng)窄帶信號(hào)。在工程實(shí)踐中，大部分信號(hào)是非平穩(wěn)或譜時(shí)變的。在窄帶信號(hào)源假設(shè)條件下，由于源信號(hào)的導(dǎo)向矢量具有頻率不一致性，因此傳統(tǒng)子空間算法，如多重信號(hào)分類算法（Multiple Signal Classification，MUSIC）、旋轉(zhuǎn)不變子空間算法（ESPRIT）等，并不適用于工程實(shí)踐［1-3］。

近年來(lái)，研究人員對(duì)寬帶信號(hào)源條件下的DOA估計(jì)進(jìn)行了大量研究［4-6］。寬帶信號(hào)DOA 估計(jì)算法原理上［7］包含非相干信號(hào)子空間方法（Incoherent Signals-subspace Method，ISM）［8］與相干信號(hào)子空間方 法（Coherent Signals-subspace Method，CSM）［9］。非相干信號(hào)子空間方法基于頻率分解原理，將寬帶源信號(hào)分解為一系列窄帶信號(hào)后進(jìn)行DOA 估計(jì)。相干信號(hào)子空間方法基于頻率聚焦原理，其中最為典型的是雙邊相關(guān)變換［10］（Two-sided Correlation Transformation，TCT）方法，利用聚焦矩陣將分解得到的子信號(hào)變換到參考頻率點(diǎn)上，進(jìn)而使用窄帶的處理方法實(shí)現(xiàn)DOA 估計(jì)。但是該方法需要CSM 預(yù)估信號(hào)源角度，而聚焦矩陣對(duì)預(yù)估角度的依賴導(dǎo)致最終DOA 估計(jì)結(jié)果產(chǎn)生偏差。文獻(xiàn)［11］提出一種聚焦的FTOPS 算法，利用參考頻點(diǎn)的信號(hào)子空間與陣列方向矢量投影矩陣間的正交性對(duì)DOA 進(jìn)行估計(jì)。文獻(xiàn)［12］對(duì)平滑自相關(guān)矩陣進(jìn)行特征分解，再根據(jù)特征向量空間之間的過(guò)渡性構(gòu)建聚焦矩陣，從而實(shí)現(xiàn)完美聚焦的目的。但是該算法僅針對(duì)環(huán)境噪聲為高斯噪聲，并且需要預(yù)設(shè)參考頻率的情況。文獻(xiàn)［13］基于壓縮感知理論，利用陣列協(xié)方差矩陣稀疏迭代估計(jì)的方法實(shí)現(xiàn)寬帶信號(hào)DOA 估計(jì)，但是該方法在信源數(shù)目預(yù)估出現(xiàn)錯(cuò)誤時(shí)，其空間譜結(jié)果易產(chǎn)生偽峰。文獻(xiàn)［14］對(duì)信號(hào)子空間聚焦法進(jìn)行改進(jìn)，采用奇異值分解重構(gòu)協(xié)方差矩陣，通過(guò)Root-正交傳播算子實(shí)現(xiàn)DOA 估計(jì)，改進(jìn)方法雖然降低了運(yùn)算量，但是仍需要預(yù)估參考頻點(diǎn)子頻帶。文獻(xiàn)［15］提出一種頻域時(shí)延補(bǔ)償方法，該方法無(wú)需對(duì)角度進(jìn)行預(yù)估，但是運(yùn)算復(fù)雜度高，難以滿足信號(hào)實(shí)時(shí)處理的要求。

針對(duì)寬帶非平穩(wěn)信號(hào)，研究人員嘗試從時(shí)頻域角度進(jìn)行研究，根據(jù)寬帶信號(hào)的時(shí)頻信息獲得更準(zhǔn)確的DOA 估計(jì)結(jié)果［16-17］。基于可調(diào)窗函數(shù)的S 變換和小波變換時(shí)頻分析方法能進(jìn)行多分辨率分析。文獻(xiàn)［18］構(gòu)建一種時(shí)頻域的陣列數(shù)據(jù)模型，根據(jù)信號(hào)的時(shí)頻信息來(lái)提高非平穩(wěn)信號(hào)的DOA 估計(jì)性能。文獻(xiàn)［19］利用小波包變換對(duì)信號(hào)進(jìn)行分解，再使用MUSIC 算法對(duì)每個(gè)子帶進(jìn)行空間譜估計(jì)。文獻(xiàn)［20］將S 變換應(yīng)用于MUSIC 算法，并對(duì)跳頻及交叉chirp 陣列信號(hào)進(jìn)行分析，然而該方法仍存在需要預(yù)知信號(hào)源個(gè)數(shù)的問(wèn)題。文獻(xiàn)［21］提出基于小波變換的多重信號(hào)分類改進(jìn)算法，根據(jù)信號(hào)的時(shí)頻域特征來(lái)提高算法的分辨率。

本文提出一種基于改進(jìn)MUSIC 算法的寬帶信號(hào)DOA 估計(jì)。通過(guò)對(duì)接收信號(hào)進(jìn)行S 變換，獲得多分辨的時(shí)頻譜矩陣，同時(shí)構(gòu)建時(shí)頻域陣列信號(hào)模型，根據(jù)頻率段不同時(shí)刻的功率譜矩陣呈聯(lián)合對(duì)角化結(jié)構(gòu)的特點(diǎn)，設(shè)計(jì)一種新的空間譜，從而實(shí)現(xiàn)寬帶信號(hào)的DOA 估計(jì)。

1 基于S 變換的時(shí)頻陣列數(shù)據(jù)模型

假設(shè)M個(gè)陣元等間隔線性排布，陣元間距為d，已知P個(gè)寬帶信號(hào)源sp(t)從不同方向入射，入射角度分別為()θ1，θ2，…，θP，不考慮傳播過(guò)程中的信號(hào)幅值衰減，第m個(gè)陣元的接收信號(hào)如式（1）所示：

其 中：τmp=(m-1)dsinθp/v為第p個(gè)信號(hào)到達(dá)第m個(gè)陣元產(chǎn)生的時(shí)延；v為信號(hào)的傳播速度；nm(t)為第m個(gè)陣元的加性噪聲。

使用S 變換對(duì)信號(hào)xm(t)進(jìn)行處理，變換結(jié)果如式（2）所示：

其中：τ為時(shí)間因子；fi(i=1，2，…，I)為頻率分量。將上式寫(xiě)成傅里葉頻譜形式進(jìn)行推導(dǎo)，如式（3）所示：

其 中：xfi，p(τ-τmp)為信號(hào)分量；nfi，m(τ)為噪聲分量。假設(shè)S 變換的頻寬和中心頻率滿足窄帶信號(hào)條件，即變換得到的信號(hào)分量為窄帶信號(hào)，如式（4）所示：

在向量形式下，陣列接收信號(hào)的時(shí)頻域模型表示如式（5）所示：

2 ST_MUSIC 算法

2.1 算法原理

根據(jù)經(jīng)典MUSIC 算法原理［22］，ST(τ，fi)的協(xié)方差矩陣如式（7）所示：

其中：RX(τ，fi)為(τ)的協(xié)方差矩陣，由于加性噪聲與源信號(hào)不相關(guān)，且自身互不相關(guān)，則變換后的噪聲方差可用表示；σi為矩陣奇異值。

由于實(shí)際接收數(shù)據(jù)長(zhǎng)度有限，因此數(shù)據(jù)協(xié)方差矩陣取其最大似然估計(jì)，如式（8）所示：

對(duì)于第n(n=1，2，…，P)個(gè)信號(hào)源，矢量bn(f)的定義需要滿足式（9）：

不考慮噪聲項(xiàng)，根據(jù)功率譜矩陣的聯(lián)合對(duì)角化結(jié)構(gòu)性質(zhì)［23］，建立如下等式：

當(dāng)RY=RY(τk，f)，k=1，2，…，K時(shí)，代價(jià)函 數(shù)［23］的簡(jiǎn)化結(jié)果如式（12）～式（14）所示：

2.2 算法步驟

本文提出的ST_MUSIC 算法主要分為以下5 個(gè)步驟：1）根據(jù)信號(hào)的有效頻段設(shè)置S 變換的頻率分量fi(i=1，2，…，I)；2）利用式（5）對(duì)接收信號(hào)進(jìn)行S 變換得ST(τ，fi)，根據(jù)式（8）計(jì)算得到協(xié)方差矩 陣3）根據(jù)式（6）構(gòu)建時(shí)頻域?qū)蛳蛄縣(θ，fi)；4）設(shè)置搜索范圍和搜索步進(jìn)Δθ，根據(jù)式（15）計(jì)算空間譜估計(jì)P(θ)；5）通過(guò)譜搜索獲得P(θ)所有譜峰所在的位置，從而得到DOA 的最大估計(jì)值。

2.3 算法分析

根據(jù)上述算法原理，ST_MUSIC 算法滿足如下條件：

1）本文算法要求符合所有寬帶信號(hào)的窄帶分量互不相關(guān)條件；

2）算法可以直接擴(kuò)展到二維源定位的情況，此時(shí)，算法的偏轉(zhuǎn)角表現(xiàn)為方位角與俯仰角的結(jié)合，平面范圍的譜搜索轉(zhuǎn)變?yōu)榭臻g譜搜索，峰值點(diǎn)坐標(biāo)為DOA 估計(jì)結(jié)果。

2.4 算法復(fù)雜度分析

本文對(duì)ST_MUSIC 算法和后續(xù)仿真中使用的雙邊相關(guān)變換方法（TCT）、基于壓縮感知理論的算法（CS_TCT）［13］及基于小波變換的MUSIC 算 法［21］（CWT_MUSIC）的計(jì)算復(fù)雜度進(jìn)行分析。本文設(shè)空間譜估計(jì)的觀測(cè)范圍的搜索點(diǎn)數(shù)為N，信號(hào)數(shù)據(jù)長(zhǎng)為L(zhǎng)，一個(gè)M×M維的數(shù)據(jù)協(xié)方差矩陣做特征分解，其計(jì)算量為O(M3)。文獻(xiàn)［13］已計(jì)算TCT 算法的計(jì)算復(fù)雜度為O(M3+2N3M3+N2M3)。在k個(gè)多頻點(diǎn)采用 CS_TCT 算法，其運(yùn)算復(fù)雜度為O(kN3M3+M2)。相 比TCT、CS_TCT 算 法，由 于ST_MUSIC 與CWT_MUSIC 算法中包含的I個(gè)尺度參數(shù)帶來(lái)較大的計(jì)算量，使算法復(fù)雜度顯著提升，CWT_MUSIC 算法的總計(jì)算量約為O(2IMLlbL+IM3+N)［21］。相 比CWT_MUSIC 算 法，ST_MUSIC算法在譜估計(jì)運(yùn)算中增加了O(IM)的運(yùn)算量，由于其不需要做特征值分解，總計(jì)算量仍小于CWT_MUSIC 算法，為O(2IMLlbL+N+IM)。

3 仿真結(jié)果分析

3.1 算法有效性分析

在多聲源場(chǎng)景下，本文對(duì)ST_MUSIC 算法進(jìn)行仿真實(shí)驗(yàn)，并與TCT、CS_TCT 以及CWT_MUSIC 這3 種算法進(jìn)行對(duì)比。本文采用16 元均勻線性陣列，構(gòu)建信噪比為0 dB、頻率范圍為165～300 Hz 的兩不相干線性調(diào)頻信號(hào)，設(shè)置信號(hào)入射角度分別為-20°和20°，采樣頻率為4 kHz，信號(hào)總長(zhǎng)度為1 024 個(gè)數(shù)據(jù)點(diǎn)。當(dāng)陣元數(shù)為16 時(shí)，CS_TCT、TCT、CWT_MUSIC、ST_MUSIC 算法的空間譜圖如圖1 所示。從圖1 可以看出，4 種算法均可以較準(zhǔn)確地估計(jì)信號(hào)角度，但是估計(jì)效果存在差異，當(dāng)信源數(shù)估計(jì)不準(zhǔn)確時(shí)，CS_TCT 算法的偽峰抑制效果受限，而偽峰抑制效果最優(yōu)的TCT 算法在精度上較其他3 種算法略差，CWT_MUSIC 算法與本文算法具有較優(yōu)的分辨率和精度。

圖1 不同算法的空間譜圖對(duì)比（陣元數(shù)為16）Fig.1 Spatial spectrogram comparison among different algorithms（the number of array elements is 16）

為進(jìn)一步驗(yàn)證算法的有效性，本文將陣元數(shù)減少至12，其他仿真條件不變，進(jìn)行二次仿真實(shí)驗(yàn)。不同算法的空間譜圖對(duì)比如圖2 所示。從圖2 可以看出，TCT 算法在陣元數(shù)減少后出現(xiàn)多個(gè)明顯偽峰，結(jié)果出現(xiàn)錯(cuò)亂，其他3 種算法對(duì)陣元數(shù)敏感度低，結(jié)果更穩(wěn)定。CS_TCT 算法估計(jì)結(jié)果準(zhǔn)確度為-20.5°和20.4°。雖然陣元數(shù)的減少不影響CS_TCT 算法最終估計(jì)結(jié)果，但是在算法仿真的零度位置出現(xiàn)較高偽峰。CWT_MUSIC 算法估計(jì)結(jié)果為-20.8°和20.9°。本文算法估計(jì)結(jié)果為-20.6°和20.6°，驗(yàn)證了本文算法的有效性。

圖2 不同算法的空間譜圖對(duì)比（陣元數(shù)為12）Fig.2 Spatial spectrogram comparison among different algorithms（the number of array elements is 12）

為進(jìn)一步探究陣元數(shù)對(duì)算法有效性產(chǎn)生的影響，在同等仿真條件下，ST_MUSIC 算法在不同陣元數(shù)下的空間譜圖如圖3 所示。從圖3 可以看出，陣元數(shù)的增加使主峰更尖銳，在提高估計(jì)精度的同時(shí)也會(huì)帶來(lái)更多的旁瓣，但其對(duì)估計(jì)結(jié)果沒(méi)有影響，而陣元數(shù)過(guò)少降低估計(jì)結(jié)果的精度，當(dāng)陣元數(shù)為4 時(shí)，DOA 估計(jì)結(jié)果偏差3°和4°，相比陣元數(shù)16，陣元估計(jì)結(jié)果的誤差保持在±0.2°。

圖3 在不同陣元數(shù)下ST_MUSIC 算法的空間譜圖Fig.3 Spatial spectrogram of ST_MUSIC algorithm under different numbers of array elements

3.2 分辨性能對(duì)比

本文采用16元均勻線性陣列，信噪比設(shè)為-10 dB，分別在信號(hào)入射角度相距120°（-45°和75°）、90°（-45°和45°）、60°（-45°和15°）、30°（-20°和10°）、10°（-5°和5°）條件下，進(jìn)行50 次隨機(jī)重復(fù)實(shí)驗(yàn)并計(jì)算4 種算法的平均分辨率。分辨率參數(shù)采用ρ=[P(θ1)+P(θ2)]/2-min{P(θ1)，P(θ1+Δθ)，…，P(θ2)}［15］，其中θ1、θ2為信號(hào)入射角，Δθ為搜索步進(jìn)，分辨率單位為dB。

在仿真實(shí)驗(yàn)中，當(dāng)信噪比降低為-10 dB 時(shí)，TCT算法的仿真結(jié)果不穩(wěn)定，并出現(xiàn)較為嚴(yán)重的偽峰，因此，本節(jié)僅對(duì)其他3 種算法進(jìn)行分辨率對(duì)比，如表1所示。從表1 可以看出，分辨率隨兩信號(hào)的入射角距增大而增高，基于壓縮感知的CS_TCT 算法相較于TCT 算法增大了角度分辨率，且分辨率結(jié)果較穩(wěn)定。本文算法在60°和120°處出現(xiàn)較高的分辨值，但是這3 種算法分辨率差距不大。

表1 不同算法的分辨率對(duì)比Table 1 Resolution comparison among different algorithms

3.3 誤差性能對(duì)比

本文仿真條件與3.1 節(jié)設(shè)置一致。本文仿真信號(hào)的信噪比范圍設(shè)置為-15～10 dB，分別通過(guò)50 次重復(fù)隨機(jī)實(shí)驗(yàn)對(duì)比4 種算法DOA 估計(jì)結(jié)果的均方根誤差，如圖4 所示。在信噪比為-15 dB 與-10 dB 時(shí)，TCT 算法估計(jì)的均方根誤差大于1°，分別為3.33°、1.60°，因此圖4 中未顯示其結(jié)果，CS_TCT 算法的DOA 估計(jì)結(jié)果在高信噪比條件下表現(xiàn)更佳。此外，在實(shí)驗(yàn)過(guò)程中，當(dāng)信噪比為-15 dB 與-10 dB 時(shí)，TCT算法的估計(jì)成功率低于50%，而CWT_MUSIC 算法與ST_MUSIC 算法在整個(gè)信噪比范圍內(nèi)估計(jì)成功率始終保持90%以上。

圖4 不同算法的均方根誤差對(duì)比Fig.4 Root mean square error comparison among different algorithms

3.4 計(jì)算量對(duì)比

本節(jié)仿真條件的設(shè)置與3.1 節(jié)保持一致。本文分別設(shè)置2 種仿真情形，分別為2 個(gè)信號(hào)（-20°和20°）和3 個(gè)信號(hào)（-60°、-20°和20°），通過(guò)50 次隨機(jī)重復(fù)實(shí)驗(yàn)對(duì)比算法的運(yùn)算時(shí)間。不同算法的平均運(yùn)算時(shí)間如圖5 所示，CS_TCT 算法的運(yùn)算時(shí)間最短，與復(fù)雜度分析結(jié)果對(duì)應(yīng)，基于頻率聚焦的TCT 算法明顯比基于頻率分解的CWT_MUSIC 算法與ST_MUSIC 算法運(yùn)算時(shí)間短，但是這4 種算法的運(yùn)算時(shí)間均在正?？山邮芊秶鷥?nèi)，ST_MUSIC 算法略優(yōu)于CWT_MUSIC 算法。

圖5 不同算法的平均運(yùn)算時(shí)間Fig.5 Average computing time of different algorithms

3.5 二維定位仿真實(shí)驗(yàn)

本文對(duì)所提算法進(jìn)行二維聲源定位仿真實(shí)驗(yàn)，其他仿真條件不變，在信噪比為0 dB 條件下，采用真實(shí)語(yǔ)音信號(hào)進(jìn)行實(shí)驗(yàn)，在不考慮環(huán)境混響的情況下，設(shè)置2 個(gè)語(yǔ)音聲源的方位角和俯仰角，分別為40°和40°、40°和-20°，定位結(jié)果如圖6 所示，從中可以看出，該算法在二維定位中能夠得到準(zhǔn)確的DOA 估計(jì)結(jié)果。

圖6 ST_MUSIC 算法的二維DOA 估計(jì)譜圖Fig.6 Two-dimensional DOA estimation spectrogram of ST_MUSIC algorithm

4 結(jié)束語(yǔ)

針對(duì)被動(dòng)探測(cè)系統(tǒng)中的寬帶信號(hào)DOA 估計(jì)問(wèn)題，本文提出一種基于S 變換且無(wú)需預(yù)估信源數(shù)的多重信號(hào)分類改進(jìn)算法。根據(jù)S 變換的多分辨率特性，通過(guò)構(gòu)建時(shí)頻陣列信號(hào)模型，提高多源DOA 估計(jì)的空間分辨率，利用譜搜索實(shí)現(xiàn)DOA 估計(jì)，實(shí)現(xiàn)多源寬帶信號(hào)的聲源定位。二維語(yǔ)音定位仿真實(shí)驗(yàn)驗(yàn)證了該算法的有效性。仿真結(jié)果表明，該算法具有較優(yōu)的分辨率性能和估計(jì)性能。后續(xù)將從S 變換參數(shù)的角度優(yōu)化本文算法，同時(shí)通過(guò)增強(qiáng)信號(hào)分量、弱化噪聲分量，降低算法運(yùn)算復(fù)雜度并提升其在復(fù)雜噪聲情況下的估計(jì)性能。