基于最大期望算法的多時變信號DOA估計方法

范保華, 左 樂,*, 唐 勇,2, 胡澤華

(1. 中國電子科技集團(tuán)公司第二十九研究所, 四川 成都 610036; 2. 電子科技大學(xué)電子科學(xué)與工程學(xué)院, 四川 成都 611731)

0 引 言

輻射源參數(shù)估計問題在雷達(dá)、聲吶、麥克風(fēng)陣列、射電天文以及移動通信中有廣泛應(yīng)用,因此其在陣列信號處理領(lǐng)域具有重要的研究價值[1-5]。微波測向通過在不同空間位置布置天線接收來波的電磁信號進(jìn)行。根據(jù)布置天線的方式可分為靜態(tài)陣列和時變合成陣列。靜態(tài)陣列每個單元均配有對應(yīng)的接收通道,通過不同陣元同時采集空間中來波信號的場分布,并通過參數(shù)估計等方法完成測向[5]。靜態(tài)陣的優(yōu)點是可以完成單脈沖測向,但缺點在于需通過較多的天線和相應(yīng)的接收通道來解模糊及高精度測向。

與同時采集來波信號場分布的靜態(tài)陣不同,時變陣列在測量過程中移動采樣位置,通過時間積累和空間掃描,利用簡化的采樣硬件獲取來波的空間場分布信息。將同時采集用分時采集替代,在處理時間與硬件資源之間進(jìn)行權(quán)衡。

無源合成陣列(passive synthetic array, PSA)既是典型的時變測向陣列。其通過較少的采樣單元數(shù)目,通過分時改變陣元位置來實現(xiàn)一維或二維場分布采集[6-7]。1996年,Friedlander首次將時變陣列用于測向,將靜態(tài)陣基于本征值的測向方法應(yīng)用于時變陣[6]。

無源合成陣列由于其簡單的系統(tǒng)架構(gòu)和較高的測向精度被廣泛應(yīng)用,但在實際工程應(yīng)用中易遇到以下難點:

(1) 信號分選問題。采用靜態(tài)陣測向時,由于采用多通道同時采樣,可同時測量每個信號的幅度、頻率、相位等信息,并以此計算得到該信號的入射角。然而采用時變陣測向時,每一次采樣僅能獲得入射場的幅度、相位、頻率等數(shù)據(jù),無法建立多個數(shù)據(jù)間的顯性關(guān)聯(lián)。當(dāng)多個輻射源在同一時間段內(nèi)交替輻射信號時,采樣數(shù)據(jù)相互交錯,但又無法事先獲得單個采樣數(shù)據(jù)與其輻射源的對應(yīng)關(guān)系。反觀入射角的估計,需首先完成信號分選,建立每個采樣信號與輻射源之間的關(guān)聯(lián)關(guān)系,才能采用來自同一個輻射源的數(shù)據(jù)進(jìn)行該輻射源信號參數(shù)的精確估計。

(2) 采樣過程中信號幅度時變問題。如對于來波為波束掃描的雷達(dá)信號[7],在不同時刻收到信號幅度不同,導(dǎo)致其接收的信噪比不穩(wěn)定,而不同信噪比(signal to noise ratio, SNR)的信號置信度不完全一樣,對于高精度測向需該考慮不同質(zhì)量信號對測向結(jié)果的貢獻(xiàn)程度。

(3) 輻射信號頻率時變問題。對于頻率間隔大于信號帶寬的信號,可以從頻域上將其分為多個輻射源。但對于頻率間隔小于信號帶寬的信號,如跳頻信號或頻率調(diào)制信號,此時信號的頻率時變特性就必須謹(jǐn)慎處理,因為若將同一輻射源產(chǎn)生的變頻信號分為多個輻射源進(jìn)行處理,則會造成信號分選錯誤。不僅如此,每組信號數(shù)據(jù)量的減少還會導(dǎo)致測向精度的下降。

(4) 信號的采樣位置取決于信號到達(dá)的時間,其分布并非均勻,空間采樣位置無周期性。

在眾多的波達(dá)方向(direction of arrival, DOA)估計方法中,最大似然(maximum likelihood, ML)方法被認(rèn)為性能穩(wěn)健且統(tǒng)計效率較高[5,8-9]。同時,ML方法容易結(jié)合最大期望(expectation maximization, EM)算法[10],進(jìn)行信號的分選以及多信號的參數(shù)估計。EM算法可將多信號進(jìn)行分離,其最大優(yōu)勢在于能將多目標(biāo)的優(yōu)化問題分解為幾個獨立的單目標(biāo)優(yōu)化問題,通過反復(fù)迭代將問題降維,從而實現(xiàn)運算量的縮減。Weinstein和Feder首先將EM算法應(yīng)用于測向領(lǐng)域,實現(xiàn)多徑時延和多輻射源定位[11]。隨后,Chung和Frenkel分別將EM算法應(yīng)用于多目標(biāo)點頻信號位置信息的獲取[12-13]。對于寬帶信號,Mada和Lu分別將EM算法應(yīng)用于信號位置參數(shù)和波形的聯(lián)合估計[14-15]。

在時變測向陣列的移動形式方面,文獻(xiàn)[16-19]將無源合成線陣用于輻射源的測向和定位。雖然無源合成線陣通過陣列的直線移動形成一個大的測向物理尺寸對空間來波進(jìn)行分時采集,但其存在移動路徑的左右無法分辨、移動姿態(tài)不易控制等問題。另一方面,由于陣列的旋轉(zhuǎn)容易控制,無源合成圓陣可通過兩傳感器的旋轉(zhuǎn)獲取各角度下的輻射源場分布信息,用于解算輻射源的二維入射角[20-32],并實現(xiàn)360°的全方位角覆蓋。文獻(xiàn)[21]提出采用兩天線旋轉(zhuǎn)形成一個無源合成圓陣,用于地球同步軌道衛(wèi)星測向。但該方法未考慮多輻射源以及相位模糊問題,因此僅適用于單輻射源以及窄帶、小角度的測向系統(tǒng)。文獻(xiàn)[22]基于信號子空間測向方法,但要求頻率和采樣間隔固定。文獻(xiàn)[23]采用粒子群優(yōu)化算法,進(jìn)行目標(biāo)參數(shù)計算,但是由于粒子群算法采用隨機(jī)搜索策略,其全局尋優(yōu)能力并不十分可靠[24]。文獻(xiàn)[25]所采用的陣元數(shù)較大,且只適用于單信號。文獻(xiàn)[26-27]采用多假設(shè)偽線性迭代最小二乘方法,但未考慮信號分選問題,故僅能處理單輻射源情形。對于多信號問題,文獻(xiàn)[28]僅利用相位差信息,提出采用正交譜配對追蹤方法解決了多目標(biāo)問題。但由于其譜值的計算過程要求頻率固定、采樣位置均勻,因此仍無法適應(yīng)數(shù)據(jù)采集過程中信號頻率、幅度、位置間隔改變等問題。除此之外,文獻(xiàn)[29-31]雖然實現(xiàn)了多目標(biāo)二維角度的無模糊估計,但均未考慮各采樣點信噪比各異的情況,會造成測向精度下降。

為同時解決上述4個無源合成陣列測向及分選所面臨的問題,本文提出采用EM算法,同時解決頻變多輻射源的測向以及信號分選問題。據(jù)筆者所知,尚未見將EM算法應(yīng)用于無源合成陣列多目標(biāo)測向的報道。對于每個輻射源的DOA估計,本文采用ML方法進(jìn)行二維測向,并考慮每個接收信號采樣位置間隔、頻率、SNR均不同的情形。但直接采用復(fù)數(shù)形式的接收信號對來波方位角及仰角的ML估計是一個非線性、多維的極值問題,需要全局多維搜索,其計算量較大。由于不同DOA產(chǎn)生的信號相位不同,相位干涉儀同樣適用于時變測向陣列。注意到接收信號的相位差與二維入射角的投影呈線性關(guān)系,本文采用相位數(shù)據(jù)推導(dǎo)了不同SNR條件下二維入射角的閉合形式解析表達(dá)式。其次,無源合成陣列與靜態(tài)陣列相同,當(dāng)涉及到采用干涉儀體制測量輻射源DOA時,會出現(xiàn)相位翻轉(zhuǎn)的模糊問題。為了解相位模糊的多值問題,通常需要增加傳感器數(shù)量來排除多值,這無疑會增加系統(tǒng)和算法的復(fù)雜度。通常的相位解模糊方法根據(jù)先計算出的一個粗角度來還原相位的真實值[32]。但已報道的相位解模糊方法只適合靜態(tài)線陣[32-34],不能直接用于無源合成圓陣。為了解決相位模糊問題,本文先采用無模糊的復(fù)數(shù)形式的ML方法,進(jìn)行粗測向,以減小譜峰搜索的運算量,再采用解模糊后的相位數(shù)據(jù)進(jìn)行入射角的解析運算。最后,數(shù)值實驗結(jié)果表明,本文提出的基于EM算法的能完成高精度測向,其測向誤差接近克拉米勞下限(Cramer Rao lower bound, CRLB),且測向與信號分選同時完成。

1 EM算法與模型建立

1.1 EM算法

EM算法最早于1977年由Dempster等人提出[10],其思路是將用于參數(shù)估計的數(shù)據(jù)進(jìn)行分組,然后在每組內(nèi)進(jìn)行參數(shù)的分別估計,再根據(jù)新的估計參數(shù)重新進(jìn)行數(shù)據(jù)的分組,以提高參數(shù)估計的質(zhì)量,并以此交替迭代,進(jìn)行數(shù)據(jù)的分解與參數(shù)的估計。EM算法的流程如圖1所示。隨著參數(shù)似然值的增加,分組結(jié)果最終也將收斂。EM算法的的一個顯著特點在于多個參數(shù)的估計不通過單次處理,而是每次參數(shù)估計后,用新參數(shù)重新劃分?jǐn)?shù)據(jù),再對劃分后的數(shù)據(jù)分別進(jìn)行新一輪的參數(shù)估計。EM算法采用迭代的方式將一個高維問題分解成多個低維問題,通過降維實現(xiàn)運算量的精簡。

結(jié)合多目標(biāo)測向與信號問題,由于輻射源的位置在旋轉(zhuǎn)周期內(nèi)相對固定,因此可依據(jù)輻射源的入射角參數(shù)進(jìn)行觀察數(shù)據(jù)的分解,即將來自于同一輻射源的接收數(shù)據(jù)分解入一組。而接收數(shù)據(jù)與輻射源位置的關(guān)系可通過ML方法求解的二維角度參數(shù)進(jìn)行表征?；诖?首先建立接收信號的模型。

1.2 接收信號模型

設(shè)間隔為d的兩天線在與XOY面內(nèi)旋轉(zhuǎn),旋轉(zhuǎn)中心為坐標(biāo)原點O,如圖2所示。在一個采樣周期內(nèi),該合成圓陣接收到來自P個位置固定的輻射源發(fā)出的N個信號。輻射源輻射的電磁信號頻率、幅度均不相同。

由于采樣時間極短,采樣周期內(nèi)的旋轉(zhuǎn)角度可忽略不計。第p個輻射源的第n個脈沖信號的接收信號時域表達(dá)式為

Vn(t)=Anexp(jωnt+jφp(n))+en(t),n=1,2，…,N

(1)

φp(n)=kp(n)dsinθpcos(φ(tn)-φp)

(2)

式中：旋轉(zhuǎn)角度φ(tn)∈[0,2π)為旋轉(zhuǎn)臂與x軸的夾角;仰角φp∈[0,2π)為來波方向與x軸夾角;方位角θp∈[0,π)為來波方向與z軸夾角;kp(n)=ωn/c為波數(shù),c為波的傳播速度。本文假設(shè)信號頻率測量無誤差。

此時,通過采樣點數(shù)為M的付氏變換,將時域信號轉(zhuǎn)化為頻域信號:

(3)

可得第無模糊相位測量值為

(4)

式中:ε(n)表示相位測量噪聲。由于不同采樣點間噪聲相互獨立,等效相位噪聲的方差[35]為

(5)

此即為頻域相位噪聲與時域高斯白噪聲的關(guān)系。其揭示了相位噪聲與信噪比的關(guān)系,即在采樣數(shù)一定的條件下,相位噪聲的方差與信噪比成反比。同時,等效相位噪聲服從高斯分布[36]。

當(dāng)天線間距d大于半個工作波長時,會出現(xiàn)相位模糊問題。由于相位的測量值在[-π,π)之內(nèi),為相位的真實值對2π取模后的余數(shù),故為

(6)

式中:h(n)為一整數(shù),也被稱為模糊數(shù)。由此可見,相位的測量值為信號的空間相差疊加相位噪聲后,對2π取模后的余數(shù)。

當(dāng)P個輻射源在一個采樣時間內(nèi)交替輻射信號,該合成圓陣的采樣信號相位為P個輻射源所產(chǎn)生的信號相位疊加,即

(7)

2 二維入射角估計

2.1 基于相位的二維測向方法

(8)

對條件概率密度函數(shù)的對數(shù)值取偏微分,并令其為零,即

(9)

和

(10)

有

(11)

及

(12)

根據(jù)三角函數(shù)的和差化積關(guān)系,以及每個相位差采樣相位誤差方差與SNR的反比關(guān)系,可得二維入射角的簡潔表達(dá)式為

(13)

(14)

式中:b1和b2為矩陣b=[b1,b2]T的兩個元素,可通過下式計算得到:

b=(det(A))-1[A][c1,c2]T,

(15)

式中:det(A)為矩陣A的行列式值。矩陣A=[(a11,a12)T·(a21,a22)T]的元素為

(16)

以及

(17)

式(16)和式(17)即通過無模糊相位值計算來波二維入射角的解析表達(dá)式。為得到無模糊相位值,還需進(jìn)行模糊解算。

2.2 基于復(fù)數(shù)響應(yīng)的相位模糊解算

(18)

式中:

w(n)=exp(jε(n))-1=x(n)+jz(n)

(19)

當(dāng)SNR較高時,w(n)的實部和虛部根據(jù)泰勒級數(shù)展開進(jìn)行一階近似,有

(20)

z(n)=sin(ε(n))?ε(n)

(21)

因此,與式(8)相對應(yīng)的yp=[yp(1),yp(2),…,yp(N)]的條件概率密度可寫為

(22)

(23)

(24)

式中:round為四舍五入運算。然后通過模糊數(shù)即可計算出相位真值,即

(25)

3 理論精度下限分析

采用信噪比分析時變陣對頻變目標(biāo)的測向誤差下限,作為評價時變陣測向性能的指標(biāo)。由式(1)可得關(guān)于方位角和仰角的Fisher矩陣元素為

(26)

式中：α1=θp;α2=φp;μn=cos(φp(n));νn=sin(φp(n))。再根據(jù)式(5)得Fisher矩陣的各元素為

(27)

求Fisher矩陣的逆,即可得到合成圓陣測量寬帶時變信號的方位角與仰角的克拉米勞精度下限,分別為

(28)

(29)

值得注意的是,若采樣點均勻,各采樣信號頻率、SNR相同,分別為k0c/2π和SNR0時,可得到仰角、方位角的克拉米勞下限表達(dá)式,分別為

var(Δφp)≥(NM×SNR0)-1(k0dsinθp)-2

(30)

var(Δθp)≥(NM×SNR0)-1(k0dcosθp)-2

(31)

可見,來波入射角遠(yuǎn)離陣面法向,仰角測向誤差增加,而方位角測向誤差減小。同時,測量相位的采樣點和場分布的采樣點數(shù)越多、SNR越高,測向誤差越小。

4 EM算法求解多信號參數(shù)

在多輻射源入射條件下,直接求解法既求下式的ML解：

(32)

式中：dis[x]=|mod(x,2π)|定義相位的距離,mod(x,2π)表示x對2π取余數(shù)運算。注意到,在信號未分選的情況下,尚無法確定相位值與輻射源的對應(yīng)關(guān)系,以及輻射源的入射角,因此對式(32)的直接求解是一個計算復(fù)雜度為非確定性多項式的問題[37]。

換個角度,若已知輻射源的入射角信息,則可依據(jù)角度信息計算得到該輻射源在合成圓陣處產(chǎn)生的理論相位差,并將與該相位差最接近的采樣數(shù)據(jù)來自該輻射源,以此進(jìn)行信號分選,從而建立輻射源與采樣數(shù)據(jù)之間的對應(yīng)關(guān)系。另一方面,若已知信號分選,即已建立各接收數(shù)據(jù)與其輻射源輻的對應(yīng)關(guān)系,則該P維問題變?yōu)镻個一維問題,通過接收數(shù)據(jù)進(jìn)行入射角解算。基于此將高維復(fù)雜問題降維的思想,本文先將觀察數(shù)據(jù)進(jìn)行隨機(jī)分組,根據(jù)每組頻率、幅度、相位差等數(shù)據(jù)分別進(jìn)行單輻射源的入射角估計,然后根據(jù)估計的入射角對觀測數(shù)據(jù)進(jìn)行重新歸類,并依次在估計和歸類間進(jìn)行迭代,直至算法收斂。合成陣列求解時變多目標(biāo)參數(shù)的EM算法流程如圖3所示,主要包含求期望值步驟和期望值最大化兩個步驟。

4.1 求期望值步驟

求期望值步驟的作用是進(jìn)行信號分選,即根據(jù)輻射源入射角計算該位置在合成圓陣處產(chǎn)生的相位差,并以采樣相位差距此理論相位差最小為依據(jù)進(jìn)行接收數(shù)據(jù)的重新分組,將的采樣數(shù)據(jù)歸為一組中。根據(jù)第p個輻射源第i次入射角估計值,其在合成圓陣處所產(chǎn)生的理論相位差為

(33)

根據(jù)每個觀測相位差與理論相位差的距離,將觀測相位差與理論相位差最接近的數(shù)據(jù)分為同一組,即

(34)

最終將觀測數(shù)據(jù)分為P組。

4.2 期望值最大化步驟

此步驟采用ML法,通過分組后的觀測數(shù)據(jù)分別計算每組數(shù)據(jù)對應(yīng)的二維入射角,即將求期望值步驟中每組觀測數(shù)據(jù)分別求解下式的ML解:

(35)

5 數(shù)值仿真實驗

本節(jié)通過數(shù)值實驗仿真分析無源合成圓陣的信號分選和測向性能,并驗證求解時變多目標(biāo)信號參數(shù)的EM算法的能力。

5.1 單輻射源情形

本節(jié)主要考察測向誤差的均方根值隨入射仰角的變化情況。兩天線間距為1.2 m,旋轉(zhuǎn)一周采樣點數(shù)為100。信號的中心頻率為1.5 GHz,頻率變化范圍為中心頻率的13%,即從1.4～1.6 GHz范圍內(nèi)隨機(jī)變化。SNR在-10～10 dB內(nèi)隨機(jī)分布。仿真采用的每個接收信號的SNR與來波頻率如圖4所示。

仰角從5°變化到75°,方位角固定為170°,分別進(jìn)行500次獨立計算,并將測向誤差的均方根與CRLB進(jìn)行對比,重點考察了考慮每個接收信號SNR不同的情形。計算結(jié)果對比如圖5所示?？紤]不同采樣點的SNR能提高測向精度,且考慮SNR的測向精度接近理論誤差下限。同時,入射仰角越大,方位角測向誤差越小,仰角測向誤差越大,與分析結(jié)論相符。

5.2 多輻射源情形

6 結(jié) 論

本文重點討論了采用無源合成陣列在多輻射源環(huán)境下對時變信號的參數(shù)估計問題,包括對多個輻射源的二維測向以及對其產(chǎn)生交錯信號的分選。建立了無源合成圓陣測向的數(shù)學(xué)模型,考慮了目標(biāo)頻變、測向位置非均勻、單次測量信噪比不同等因素。提出了采用無源合成陣列進(jìn)行時變多信號參數(shù)估計與分選的EM算法。該算法通過將高維多參數(shù)估計問題分解成多個低維參數(shù)估計問題,從而降低運算復(fù)雜度。根據(jù)基于相位數(shù)據(jù)的ML估計方法,推導(dǎo)出閉合形式且性能接近最優(yōu)估計的入射角計算解析解。此外,采用基于接收響應(yīng)復(fù)數(shù)的ML方法,利用旋轉(zhuǎn)帶來的基線多樣性實現(xiàn)了相位模糊的解算。通過兩種不同似然方法的有機(jī)結(jié)合,兼顧了入射角初值計算運算量與最終入射角的求解精度。通過理論下限推導(dǎo),分析了無源合成圓陣的測向精度,揭示了入射仰角、陣列口徑、采樣點數(shù)和信噪比對二維測向精度的影響。最后通過單輻射源和多輻射源的數(shù)值實驗驗證了分析結(jié)論和算法的有效性。