李 森 王基福 林 彬
(大連海事大學(xué)信息科學(xué)技術(shù)學(xué)院 大連 116026)
輻射源的位置可以通過測(cè)量多個(gè)傳感器接收信號(hào)的到達(dá)時(shí)差(Time Difference Of Arrival, TDOA)來確定[1]。經(jīng)典的基于相關(guān)的TDOA估計(jì)器可以使用幾種廣義互相關(guān)技術(shù)中的一種來構(gòu)造[2]。在許多傳播環(huán)境中,每個(gè)傳感器的接收信號(hào)中除了包含直達(dá)輻射源信號(hào)外還包括一個(gè)或多個(gè)輻射源信號(hào)的多徑副本信號(hào)。由于受到多徑干擾的影響,傳感器間的互相關(guān)呈現(xiàn)多個(gè)峰值。互相關(guān)峰值之間的延遲分辨能力、相對(duì)幅度和相位偏移決定了多徑副本信號(hào)對(duì)直達(dá)輻射源信號(hào)TDOA估計(jì)精度的影響。多徑環(huán)境下的TDOA估計(jì)技術(shù)有很多種,在輻射源信號(hào)已知的情況下匹配濾波提供了估計(jì)多徑信號(hào)TDOA的近似最優(yōu)方法[3]。但是這種方法中多徑TDOA的分辨能力受到輻射源信號(hào)相關(guān)時(shí)間的限制,當(dāng)多徑TDOA之間的時(shí)間間隔小于輻射源信號(hào)的相關(guān)時(shí)間時(shí)就必須應(yīng)用高分辨率的多徑TDOA估計(jì)算法。
高分辨率的多徑TDOA估計(jì)算法主要可以分為以下幾類:第1類是在輻射源信號(hào)已知的情況下將多徑信號(hào)TDOA估計(jì)的多維優(yōu)化問題轉(zhuǎn)化為多個(gè)1維優(yōu)化問題,從而降低優(yōu)化問題求解的復(fù)雜度,這類算法主要有WRELAX (Weighted RELAXation)算法[4,5]和期望最大化(Expectation Maximization,EM)算法[6,7]。第2類是在輻射源信號(hào)未知的情況下,將時(shí)域的多徑TDOA估計(jì)問題轉(zhuǎn)化為對(duì)頻域多個(gè)正弦信號(hào)的頻率估計(jì)問題,然后通過子空間技術(shù)實(shí)現(xiàn)對(duì)TDOA的估計(jì)[8,9]。第3類是利用時(shí)延值與信道沖激響應(yīng)的稀疏支撐集之間的關(guān)系,將多徑信號(hào)的TDOA估計(jì)問題轉(zhuǎn)化為稀疏重構(gòu)問題,通過合理地構(gòu)造測(cè)量矩陣并結(jié)合不同的稀疏重構(gòu)算法來估計(jì)TDOA[10,11]。還有一些其他應(yīng)用于高分辨率多徑時(shí)間延遲估計(jì)的方法,如重要性采樣方法[12]和馬爾可夫鏈蒙特卡羅方法[13]等。
上述多徑信號(hào)TDOA估計(jì)算法都是假設(shè)接收傳感器處于高斯分布噪聲環(huán)境中,但在實(shí)際的無線通信和水聲通信環(huán)境中,自然現(xiàn)象或者人為因素的影響使得傳感器的輸入噪聲具有明顯的尖峰脈沖特性。這類尖峰脈沖特性噪聲的概率密度分布曲線與高斯分布的概率密度曲線類似,但比高斯分布概率密度曲線具有更厚的拖尾。對(duì)稱α-穩(wěn)定分布(Symmetric α-Stable distribution, SαS) 是描述這類具有尖峰脈沖性的非高斯噪聲的一種非常廣泛的數(shù)學(xué)模型[14]。SαS分布的概率密度可用3個(gè)參數(shù)來描述:特征指數(shù)α、分散系數(shù)γ和位置參數(shù)a。其中特征指數(shù)α( 0 <α ≤2)表示分布函數(shù)的拖尾厚度,α取值越小,表示SαS分布的拖尾越厚,噪聲脈沖性越強(qiáng);α取值越大,表示SαS分布的拖尾越薄,噪聲脈沖性越弱。當(dāng)α=2時(shí),SαS分布退化為高斯分布。當(dāng)0<α<2時(shí),SαS分布不具有α階及其以上階矩。由于SαS分布的這個(gè)特性,使得經(jīng)典的基于2階統(tǒng)計(jì)量提出的各種TDOA估計(jì)算法的性能發(fā)生退化,甚至不能得到正確的估計(jì)結(jié)果。于是Shao等人[15]提出了分?jǐn)?shù)低階統(tǒng)計(jì)量的概念并將傳統(tǒng)的2階和高階統(tǒng)計(jì)量信號(hào)處理方法拓展到分?jǐn)?shù)低階統(tǒng)計(jì)量信號(hào)處理方法。
為了解決SαS分布脈沖噪聲環(huán)境下的多徑信號(hào)時(shí)延估計(jì)問題,文獻(xiàn)[16]利用分?jǐn)?shù)低階統(tǒng)計(jì)量理論對(duì)高分辨率多徑時(shí)延估計(jì)EM算法進(jìn)行改進(jìn),提出一種基于分?jǐn)?shù)低階相關(guān)的P-EM算法,仿真結(jié)果表明P-EM算法對(duì)尖峰脈沖噪聲具有更好的適應(yīng)性。劉文紅等人[17]用分?jǐn)?shù)低階統(tǒng)計(jì)量理論中的最小分散系數(shù)準(zhǔn)則代替最小二乘準(zhǔn)則,并利用波形松弛思想將高維優(yōu)化問題分解為多個(gè)1維優(yōu)化問題,提出了一種韌性的高分辨率多徑時(shí)延估計(jì)方法P-WRELAX。文獻(xiàn)[18]將兩路接收信號(hào)的分?jǐn)?shù)低階互協(xié)方差譜序列當(dāng)作時(shí)間序列,利用子空間旋轉(zhuǎn)不變技術(shù)實(shí)現(xiàn)了脈沖噪聲環(huán)境下多徑時(shí)延的高分辨率估計(jì)?;诜?jǐn)?shù)低階統(tǒng)計(jì)量理論的信號(hào)處理算法要求選取的分?jǐn)?shù)低階統(tǒng)計(jì)量的階數(shù)p小于SαS分布的特征指數(shù)α,而在實(shí)際系統(tǒng)中SαS分布的特征指數(shù)是不易獲得的,從而限制了基于分?jǐn)?shù)低階統(tǒng)計(jì)量理論的信號(hào)處理算法的應(yīng)用。相關(guān)熵作為隨機(jī)變量間局部相似性的度量不僅具有良好的抗噪性能,而且不依賴于脈沖噪聲的先驗(yàn)知識(shí),近年來受到了國(guó)內(nèi)外學(xué)者的極大關(guān)注[19,20]。文獻(xiàn)[21]首先證明了在特定條件下最大相關(guān)熵準(zhǔn)則和分?jǐn)?shù)低階統(tǒng)計(jì)量理論中的最小分散系數(shù)準(zhǔn)則是等價(jià)的,然后提出了一種基于最大相關(guān)熵準(zhǔn)則的自適應(yīng)時(shí)延估計(jì)算法。文獻(xiàn)[22]將最大相關(guān)熵準(zhǔn)則分別與Farrow非整數(shù)時(shí)延濾波器和Lagrange非整數(shù)時(shí)延濾波器相結(jié)合解決了脈沖噪聲環(huán)境下的非整數(shù)時(shí)延估計(jì)問題。文獻(xiàn)[23]通過對(duì)兩路接收信號(hào)的互相關(guān)熵譜密度進(jìn)行白化加權(quán)處理提高了相關(guān)熵時(shí)延估計(jì)算法抗尖峰脈沖噪聲的能力。文獻(xiàn)[24]提出一種基于相關(guān)熵希爾伯特差值的時(shí)延估計(jì)算法,克服了窄帶射頻信號(hào)時(shí)延估計(jì)的精度受相對(duì)帶寬限制和脈沖噪聲影響的問題。文獻(xiàn)[25]基于最大相關(guān)熵準(zhǔn)則提出了一種針對(duì)時(shí)變線性系統(tǒng)狀態(tài)估計(jì)的擴(kuò)展遞歸最大相關(guān)熵自適應(yīng)算法,并利用所提出的算法解決了ADS-B報(bào)文的緩變時(shí)延估計(jì)問題。文獻(xiàn)[26]研究了脈沖噪聲環(huán)境下振幅衰減信號(hào)的時(shí)延估計(jì)問題,基于最大相關(guān)熵準(zhǔn)則和自適應(yīng)時(shí)延估計(jì)思想提出了兩步法時(shí)延估計(jì)算法;為了避免兩步估計(jì)的誤差累積和估計(jì)偏差,進(jìn)一步提出了一種一步法自適應(yīng)時(shí)延估計(jì)算法。受相關(guān)熵是隨機(jī)變量相似性度量思想的啟發(fā),本文對(duì)期望最大化高分辨率多徑TDOA估計(jì)算法進(jìn)行改進(jìn),提出了基于相關(guān)熵的期望最大化算法(Correntropy based Expectation Maximization, CEM)。仿真實(shí)驗(yàn)表明,該算法在不依賴于脈沖噪聲先驗(yàn)信息的情況下具有比PEM算法和P-WRELAX算法更好的脈沖噪聲適應(yīng)性和估計(jì)精度。
兩個(gè)隨機(jī)變量X和Y之間的相似性可以用相關(guān)熵來度量,其定義為[19]
相關(guān)熵具有對(duì)稱性和有界性,即Vσ(X,Y)=Vσ(Y,X)和Vσ(X,Y)≤Vσ(0)。由相關(guān)熵可以誘導(dǎo)出一個(gè)距離測(cè)度CIM(Correntropy Induced Metric)。與常規(guī)的距離測(cè)度不同,CIM距離測(cè)度表現(xiàn)出由核長(zhǎng)σ控制的混合范數(shù)特征。即當(dāng)核長(zhǎng)σ較大時(shí),CIM距離測(cè)度的L2范數(shù)距離區(qū)域也較大;當(dāng)核長(zhǎng)σ較小時(shí),其L2范數(shù)距離區(qū)域縮小而L1范數(shù)距離區(qū)域和L0范數(shù)距離區(qū)域擴(kuò)大?;贑IM距離測(cè)度,可以定義最大相關(guān)熵準(zhǔn)則為
在多徑信號(hào)傳播模型中,一般假設(shè)源信號(hào)s(k)是已知的或是可以通過某種方式獲得的,傳感器的接收信號(hào)可以表示為
期望最大化算法是一種迭代優(yōu)化策略,由于它的每一次迭代都分兩步,其中一個(gè)為期望步(Expectation step, E-step),另一個(gè)為極大步(Maximization step, M-step),所以算法簡(jiǎn)稱為EM算法。EM算法最初作為一種數(shù)據(jù)添加算法主要用于解決在數(shù)據(jù)量非常大的情況下存在數(shù)據(jù)缺失或者數(shù)據(jù)不可用的問題。隨著理論的發(fā)展,EM算法已經(jīng)不僅僅用于處理數(shù)據(jù)缺失問題,而成為一種用于解決復(fù)雜問題的簡(jiǎn)化策略。利用EM算法的基本思想,通過引入恰當(dāng)?shù)臐撛跀?shù)據(jù)可以將復(fù)雜的多維優(yōu)化問題分解為若干個(gè)獨(dú)立的1維優(yōu)化問題。
E-step:對(duì)i =1,2,···,D計(jì)算
式(8)是一個(gè)求解兩個(gè)參數(shù)的最優(yōu)化問題。它的求解過程可以首先根據(jù)期望步給定的m)求得m+1),然后再由m+1)求解新的m+1)。根據(jù)式(8)可得m+1)的估計(jì)值為
由于SαS分布脈沖噪聲不存在有限的2階矩,因此文獻(xiàn)[16]利用最小分散系數(shù)準(zhǔn)則代替了式(8)中的最小均方誤差準(zhǔn)則提出了一種適用于脈沖噪聲環(huán)境下的魯棒的EM算法,由于在該算法中迭代公式(9)和式(10)中的2階相關(guān)被p-階分?jǐn)?shù)低階相關(guān)所代替,因此稱該算法為P-EM算法,其中參數(shù)1
本文用最大相關(guān)熵準(zhǔn)則代替EM算法中的最小均方誤差準(zhǔn)則,即在第m次迭代的M-step中
對(duì)i=1,2,···,D計(jì)算
取源信號(hào) s(k)為B P S K 信號(hào),載波頻率為fc=0.25fs, fs為 采樣頻率且fs=105Hz, Ts=1/fs為采樣間隔,碼元持續(xù)時(shí)間為 20Ts,信號(hào)帶寬B =0.1/Ts。分別考慮多徑數(shù)目D =2 和 D =3兩種情況,在 D =2的情況下,兩路信號(hào)的TDOA和衰減系數(shù)分別是τ1=28Ts, a1=0.8和 τ2=32Ts, a2=0.6。D =3 是在D =2的 情況下增加了τ3=36Ts, a3=0.7的第3路多徑信號(hào)。各路多徑信號(hào)之間的TDOA之差為 ?τ =|τ3?τ2|=|τ2?τ1|=4Ts<1/B =10Ts,屬于基本相關(guān)法中的多徑不可分辨情況。
仿真中3路信號(hào)的TDOA和衰減系數(shù)初值分別設(shè)為0及0.5, β1, β2和β3分別設(shè)為0.25, 0.40和0.35,實(shí)驗(yàn)數(shù)據(jù)長(zhǎng)度為3000點(diǎn)。由于SαS分布不具有有限的2階矩,定義廣義信噪比 GSNR=10lg(σs2/γn),其中 γn為 SαS分布的分散系數(shù),σs2為信號(hào)的方差。采用均方根誤差(Root Mean Square Error, RMSE)和平均準(zhǔn)確率 PA來衡量EM算法、P-EM算法、P-WRELAX算法和CEM算法的估計(jì)性能。它們的定義分別為
實(shí)驗(yàn)1給定SαS分布的特征指數(shù)α=1.4,圖1給出了EM算法、P-EM算法、P-WRELAX算法和CEM算法分別在2徑和3徑情況下的均方根誤差和平均準(zhǔn)確率隨著廣義信噪比GSNR變化的曲線。從圖1可以明顯看出,4種算法在2徑情況下的估計(jì)性能都要優(yōu)于3徑情況下的估計(jì)性能,并且無論在2徑還是3徑情況下采用了脈沖噪聲抑制技術(shù)的P-EM算法、P-WRELAX算法和CEM算法的性能都要優(yōu)于沒有采用脈沖噪聲抑制技術(shù)的EM算法。特別是在2徑情況下P-EM算法、P-WRELAX算法和CEM算法在GSNR>–2 dB時(shí)平均準(zhǔn)確率為100%,均方根誤差為零,而EM算法要到GSNR>–8 dB時(shí)才會(huì)得到平均準(zhǔn)確率為100%,均方根誤差為零。與基于分?jǐn)?shù)低階矩的P-EM算法和P-WRELAX算法相比,基于相關(guān)熵的CEM算法在低信噪比下的性能更有優(yōu)勢(shì),例如在2徑情況下當(dāng)GSNR=–8 dB時(shí)P-EM算法的準(zhǔn)確率為58.5%, P-WRELAX算法的準(zhǔn)確率為64.5%,而CEM算法的準(zhǔn)確已經(jīng)達(dá)到90%以上。
實(shí)驗(yàn)2給定廣義信噪比GSNR=–3 dB,圖2給出了EM算法、P-EM算法、P-WRELAX算法和CEM算法的均方根誤差和平均準(zhǔn)確率隨著SαS分布噪聲的特征指數(shù)α變化的曲線。從圖中可以看出在2徑和3徑的情況下P-EM 算法、P-WRELAX算法和CEM算法的平均準(zhǔn)確率和均方根誤差在各種脈沖噪聲環(huán)境下的性能都要優(yōu)于EM算法,說明P-EM算法、P-WRELAX算法和CEM算法均具有對(duì)脈沖噪聲的抑制作用。但是基于相關(guān)熵的CEM算法對(duì)脈沖噪聲的適應(yīng)性要優(yōu)于基于分?jǐn)?shù)低階矩的P-EM算法和P-WRELAX算法,特別是在脈沖噪聲的特征指數(shù)α<1的強(qiáng)脈沖噪聲環(huán)境下,例如在3徑情況下當(dāng)α=0.6時(shí),P-EM算法和P-WRELAX算法的平均準(zhǔn)確率只有32.2%和60.9%,而CEM算法的平均準(zhǔn)確率為99.1%,這是因?yàn)檫@兩種算法中利用的p階分?jǐn)?shù)低階矩雖然在α<1的情況下可以通過求時(shí)間平均的方法來計(jì)算,但理論上是沒有定義的;而CEM算法中利用的相關(guān)熵在特征指數(shù) 0 <α ≤2的情況下都是有定義的。
圖1 不同廣義信噪比下各算法的估計(jì)性能
實(shí)驗(yàn)3給定廣義信噪比GSNR=–8 dB,圖3給出了CEM算法在SαS穩(wěn)定分布噪聲的特征指數(shù)α=0.5, 1.0, 1.5和2.0時(shí)的估計(jì)性能隨著核長(zhǎng)參數(shù)σ的變化曲線。從圖中可以看出,在某一給定特征指數(shù)的SαS穩(wěn)定分布噪聲環(huán)境下,為了獲得最佳的估計(jì)性能存在著最佳核長(zhǎng)選取問題,且最佳核長(zhǎng)值與SαS穩(wěn)定分布噪聲的特征指數(shù)密切相關(guān)。當(dāng)給定的核長(zhǎng)較小時(shí)(如σ=1),隨著噪聲脈沖性的增加,CEM算法的估計(jì)性能越來越好;當(dāng)給定的核長(zhǎng)較大時(shí)(如σ=14),隨著噪聲脈沖性的增加,CEM算法的估計(jì)性能卻越來越差。這是因?yàn)楫?dāng)核長(zhǎng)較小時(shí),相關(guān)熵誘導(dǎo)的CIM距離測(cè)度的L1和L0范數(shù)距離空間增大,L1和L0范數(shù)對(duì)脈沖噪聲具有很好的抑制作用,所以隨著噪聲脈沖性的增加,CEM算法的估計(jì)性能也越來越好。當(dāng)核長(zhǎng)較大時(shí),相關(guān)熵誘導(dǎo)的CIM距離測(cè)度的L1和L0范數(shù)距離空間被壓縮,與信號(hào)2階相關(guān)等價(jià)的L2范數(shù)距離空間增大,也就是說此時(shí)相關(guān)熵對(duì)脈沖噪聲的抑制能力減弱卻更多地聚集了信號(hào)中的2階相關(guān)信息,所以隨著噪聲脈沖性的增加,CEM算法的估計(jì)性能下降。CEM算法中核長(zhǎng)參數(shù)的選擇雖然與脈沖噪聲的特性有關(guān),但它與P-EM算法中分?jǐn)?shù)低階參數(shù)p(p<α)必須依賴于脈沖噪聲的先驗(yàn)信息是不同的。
圖2 不同特征指數(shù)α下各算法的估計(jì)性能
圖3 不同核長(zhǎng)σ和不同特征指數(shù)α下CEM算法的估計(jì)性能
圖4 不同核長(zhǎng)σ和不同廣義信噪比GSNR下CEM算法的估計(jì)性能
實(shí)驗(yàn)4給定SαS分布的特征指數(shù)α=0.5,圖4給出了CEM算法在廣義信噪比GSNR=–12 dB,–8 dB, –4 dB和0 dB時(shí)的估計(jì)性能隨著核長(zhǎng)參數(shù)σ的變化曲線。從圖中可以看出,在4種廣義信噪比的環(huán)境下,CEM算法能夠獲得最佳估計(jì)性能的核長(zhǎng)參數(shù)的取值范圍都在1.4~1.8范圍內(nèi),結(jié)合實(shí)驗(yàn)3可以說明相關(guān)熵函數(shù)中核長(zhǎng)參數(shù)的最優(yōu)選取與SαS分布的脈沖特性密切相關(guān),而與表征脈沖噪聲強(qiáng)弱的廣義信噪比無關(guān)。在相同脈沖特性不同脈沖強(qiáng)度的SαS分布噪聲環(huán)境下,當(dāng)核長(zhǎng)參數(shù)的取值在大數(shù)范圍內(nèi)時(shí),為了獲得相同的估計(jì)性能,對(duì)于強(qiáng)度較大的脈沖噪聲應(yīng)該選取較小的核長(zhǎng)。例如圖4中在σ>4的范圍內(nèi),為了獲得99.5%的平均準(zhǔn)確率,當(dāng)GSNR=0 dB時(shí),應(yīng)選取的核長(zhǎng)為σ=10,而當(dāng)脈沖噪聲的強(qiáng)度增加到GSNR=–4 dB時(shí),核長(zhǎng)參數(shù)降為σ=6。
為了解決基于EM算法的多徑信號(hào)TDOA估計(jì)算法在脈沖噪聲環(huán)境下性能退化的問題,本文基于最大相關(guān)熵準(zhǔn)則提出了基于相關(guān)熵的期望最大化高分辨率TDOA估計(jì)算法(CEM)。仿真實(shí)驗(yàn)結(jié)果表明CEM算法和基于分?jǐn)?shù)低階矩的P-EM算法和P-WRELAX算法一樣對(duì)脈沖噪聲具有很好的抑制作用,并且在強(qiáng)脈沖噪聲和低信噪比的環(huán)境下CEM算法的估計(jì)性能要優(yōu)于基于分?jǐn)?shù)低階矩的算法,而且CEM算法中參數(shù)的選取不依賴于脈沖噪聲的先驗(yàn)信息。