汝小虎,王 聰,李 悅,姜文利
(1.國防科學(xué)技術(shù)大學(xué) 電子科學(xué)與工程學(xué)院, 長沙410073;2.海軍航空工程學(xué)院 電子信息工程系,山東 煙臺264000;3.解放軍61541 部隊, 北京100094)
線性調(diào)頻信號作為一種大時寬-帶寬積信號在脈沖壓縮雷達中有重要的應(yīng)用,準確估計出它的參數(shù)對雷達的干擾和個體識別等有重要意義。
已有很多文獻研究了LFM 信號的參數(shù)估計問題,包括T.J.Abatzoglou 提出的載頻和調(diào)制斜率聯(lián)合估計的ML 方法[1]、P.M.Djuric 等人提出的解相位方法[2]、S.Peleg 等人提出的離散多項相位變換(Discrete Polynomial-Phase Transform, DPT)方法[3],以及基于chirp 傅里葉變換方法[4]等。但是這些研究都只針對一個信號,當信號信噪比較低時,其性能會受到很大影響。如解相位方法的信噪比門限為8 dB左右,DPT 方法在信噪比為0 dB時估計方差會比克拉美-羅限(Cramér-Rao Lower Bound, CRLB)高出約60%。文獻[5]和[6]分別改進了解相位方法和DPT方法,但對調(diào)制斜率估計精度的提高程度有限。另外,這些LFM 信號參數(shù)估計方法常常要在離散的頻域或二維平面上進行極值搜索,其準確度和時間效率無法兼顧。
在實際環(huán)境中截獲到的具有相同調(diào)制波形的雷達脈沖不止一個,而且對LFM 信號感興趣的往往只是它的調(diào)制斜率。但是針對多個LFM 脈沖信號的調(diào)制斜率估計問題,尚未發(fā)現(xiàn)有文獻提出解決方法。本文提出的方法首先基于脈沖積累得到信噪比較高的調(diào)制波形估計,接著解出相位,之后用最小二乘法估計調(diào)制斜率。其關(guān)鍵步驟是估計調(diào)制波形,典型方法是Howard 提出的ML 估計方法[7],本文研究了其近似的快速計算方法,給出了其適用條件。本文方法的優(yōu)勢在于可以處理多個脈沖信號,不需要搜索操作,并且經(jīng)論證是一種最大似然估計方法,參數(shù)估計精度高。
由于數(shù)據(jù)采集的先后順序不同和多普勒效應(yīng)等原因,具有相同調(diào)制波形的雷達脈沖會有相對的時延和載頻的不一致,脈沖積累的前提是使得脈沖在時間和頻率上對齊。假設(shè)已實現(xiàn)了理想的時頻對齊,用下面的信號模型描述具有相同調(diào)制波形的LFM 脈沖:
其中,NS為采樣點數(shù),NP為脈沖數(shù);Ak=A ejk是幅度因子,它的模A 表示脈沖幅度(對于短時間內(nèi)截獲到的NP 個脈沖,可認為A 對所有脈沖都一樣),它的相位k是無用參數(shù),不進行估計;εk是復(fù)高斯噪聲,滿足E(εk)=0 ,E(εkεHk)=2σ2I,其中2σ2為噪聲功率,I 為單位矩陣;0、f 0 和K 分別為相位調(diào)制
的初相、載頻和調(diào)制斜率,其中TS為采樣間隔。
由于調(diào)制波形的ML 估計結(jié)果不可避免地會有初相存在,所以式(1)中的信號模型認為相位調(diào)制存在初相,它并不影響調(diào)制斜率的估計。該信號模型中,未知參數(shù)是{A,0,f0,K},本文關(guān)心的是調(diào)制斜率K 的估計。
其中A′k= NSAk。式(3)中未知參量是復(fù)幅度A′k和調(diào)制波形μ,NP個脈沖zk的聯(lián)合概率密度函數(shù)為
將式(4)中的求和項作如下變形:
由瑞利準則[8]可知,當μ是Z 矩陣最大特征值對應(yīng)的特征向量時,式(5)中μHZμ達到最大。所以調(diào)制波形的ML 估計 μ是Z 矩陣最大特征值對應(yīng)的特征向量。
直接計算矩陣的特征向量會使時間代價很大。文獻[7]提到,如果Z 矩陣的最大特征值超過0.5,可以通過迭代的方式快速求得 μ的近似結(jié)果,與直接求特征向量得到的結(jié)果非常接近。第4 節(jié)仿真實驗會證實,當脈沖信噪比高于0 dB時可以利用這種近似計算方法。
Tretter 指出,當信噪比較高(8 dB以上)時,觀測信號中的加性復(fù)高斯噪聲可等效為相位上的加性實高斯噪聲,并且在高斯噪聲條件下最小二乘估計等效于ML 估計[9]。所以,為了提高調(diào)制斜率的估計精度,得到調(diào)制波形的ML 估計后,用文獻[10]提到的基于相對無模糊相位重構(gòu)的方法解相位得到相位調(diào)制的估計 =[ (0), (1), …, (NS-1)]T,然后對它用最小二乘法進行二次曲線擬合估計調(diào)制斜率K。設(shè)相位調(diào)制的理論值為 (i)=a0 +a1 iTs +a2i2T2s,i=0,1, …,NS-1,令a =[ a0,a1,a2]T,對它的估計為
那么調(diào)制斜率K 的估計可由 a2/π得到。下面求 a的解析表達式。擬合誤差為
其中矩陣T 為
綜合以上分析,基于脈沖積累的LFM 脈沖調(diào)制斜率的估計方法總結(jié)如下。
(1)如果脈沖的平均信噪比低于0 dB,由zk,k =1,2, …, NP構(gòu)造脈沖積累矩陣Z,求該矩陣最大特征值對應(yīng)的特征向量,得到調(diào)制波形的估計 μ,然后執(zhí)行第3 步;否則,如果脈沖平均信噪比高于0 dB,執(zhí)行第2 步。
(2)設(shè)zmax為信噪比最大的脈沖,調(diào)制波形的初始估計為 μ(0)=norm[zmax] ,其中norm[·]表示歸一化。取m =0,重復(fù)計算
直到
其中β 是個很小的數(shù),比如0.001,調(diào)制波形的ML估計近似為 μ≈ μ(m+1)。
(3)用基于相對無模糊相位重構(gòu)的方法解出 μ的相位 。
(4)由式(8)、(9)計算二次曲線的擬合系數(shù) a,那么調(diào)制斜率的估計為 K = a2/π。
本節(jié)對3.2 節(jié)提出的調(diào)制斜率估計方法進行討論,證明該方法是最大似然的。
首先定義函數(shù)β=g(v),其中v 是一個幅度固定的信號矢量。函數(shù)g 的功能是先解出v 的相位,然后用最小二乘法擬合該相位,得到的二次項系數(shù)與π的比值為函數(shù)值β。根據(jù)3.2 節(jié)提出的調(diào)制斜率估計方法,有 K =g( μ)。另外,g 是一對一的函數(shù),這是因為如果v 的值發(fā)生改變,由于它的模是固定的,所以它的相位會發(fā)生變化,由式(8)可知所有的擬合系數(shù)都會改變,最后會得到不同的β;另一方面,如果β 發(fā)生變化,意味著二次曲線會發(fā)生變化,即v 的相位會改變,所以對應(yīng)不同的v。
如果v=μ,由g 函數(shù)的功能可知,調(diào)制斜率K=g(μ)。又因為 μ是調(diào)制波形μ的ML 估計,由ML 估計的不變性[11]可以推得 K =g( μ)是K 的ML估計,即本文提出的調(diào)制斜率估計方法是最大似然的,故將該方法稱為基于脈沖積累的調(diào)制斜率ML估計方法。
需要說明的是,上述結(jié)論成立的前提條件是能準確解出 μ的相位,這就要求 μ有較高的信噪比。本文基于脈沖積累估計調(diào)制波形,使得 μ可以滿足信噪比要求,第4 部分會對此進行仿真驗證。
本節(jié)先考查脈沖積累矩陣Z 的最大特征值隨信噪比的變化。仿真時設(shè)脈沖數(shù)NP=100,采樣率fS=250 MHz,脈沖寬度為2 μs,噪聲功率2σ2=2,式(1)中幅度因子Ak的相位在[0,2π] 間隨機分布,LFM 相位調(diào)制的初相0=π/3,載頻f0=40 MHz,調(diào)制斜率K=1 MHz/μs。
圖1 給出了脈沖積累矩陣Z 的最大特征值隨信噪比的變化,其中脈沖信噪比定義為脈沖內(nèi)信號與噪聲平均功率的比值,所以式(1)中脈沖的信噪比為A2/2σ2,圖中以dB 作為單位??梢钥闯?當信噪比高于0 dB時,Z 矩陣的最大特征值大于0.5,可以用ML 方法的近似迭代計算方法。此外,由圖1 可知,脈沖數(shù)降低會使Z 矩陣的最大特征值增大,而改變采樣率對最大特征值的影響很小。
圖1 脈沖積累矩陣的最大特征值與脈沖信噪比的關(guān)系Fig.1 The biggest eigenvalue of pulse accumulation matrix vs.the SNR of pulses
現(xiàn)在研究調(diào)制波形估計結(jié)果的信噪比與脈沖信噪比的關(guān)系。調(diào)制波形的估計為 μ=μ+ε,其中ε= μ-μ認為是估計結(jié)果中的噪聲,所以調(diào)制波形估計結(jié)果的信噪比為
由于 μ相對μ可能會有不同的初相,但是這對后續(xù)的調(diào)制參數(shù)估計沒有任何影響,所以在用式(10)研究信噪比之前應(yīng)該先消除這一初相上的不一致。
用ML 方法估計調(diào)制波形μ,脈沖數(shù)NP分別設(shè)為100、40 和10,其他設(shè)置不變,蒙特卡洛仿真200次,得到的調(diào)制波形估計結(jié)果的信噪比與脈沖信噪比的關(guān)系如圖2 所示。由圖可知,調(diào)制波形估計結(jié)果的信噪比相對脈沖的信噪比高出約10 lg NP(dB),說明脈沖積累可以大大改善信噪比,有利于準確解出相位調(diào)制。
圖2 調(diào)制波形估計結(jié)果的信噪比與脈沖信噪比的關(guān)系Fig.2 The SNR of the estimated modulation waveform vs.the SNR of pulses
本節(jié)研究不同脈沖數(shù)目條件下調(diào)制斜率K 的估計方差與信噪比的關(guān)系,分析提出的ML 估計方法的門限效應(yīng)和脈沖積累效果。先給出式(1)信號模型調(diào)制斜率估計的CRLB。當脈沖只有一個,即NP=1 時,文獻[12]給出了調(diào)制斜率估計的CRLB,約為90/(π2N5ST4S·SNR),其中SNR 為脈沖信噪比。對于式(1)中獨立同分布的NP個觀測數(shù)據(jù),參數(shù)估計的CRLB 是單次觀測的1/NP倍[11],所以脈沖積累條件下LFM 脈沖調(diào)制斜率估計的CRLB 為調(diào)制斜率的估計方差為
其中,M 為蒙特卡洛仿真次數(shù), Km為第m 次仿真對調(diào)制斜率K 的估計。脈沖數(shù)NP分別設(shè)為1、10 和40,蒙特卡洛仿真次數(shù)為200 次,其他設(shè)置與4.1 節(jié)相同。得到的調(diào)制斜率估計方差與脈沖信噪比的關(guān)系曲線如圖3 所示,其中縱坐標為10 lg(σ2K)??梢钥闯?本文提出的調(diào)制斜率估計方法存在門限效應(yīng),并且利用的脈沖數(shù)目越多信噪比門限越低。對于只有一個脈沖的情況,當信噪比在8 dB左右時調(diào)制斜率的估計方差就迅速偏離了CRLB;而對于40 個脈沖的情況,當信噪比在-5 dB左右時估計方差仍能達到CRLB。這說明本文提出的調(diào)制斜率估計方法可以在低信噪比的情況下大大降低估計方差,具有優(yōu)越的性能。
圖3 調(diào)制斜率估計方差與脈沖信噪比的關(guān)系Fig.3 The estimation variance of chirp rate vs.the SNR of pulses
現(xiàn)在進一步研究基于多個脈沖相對基于單一脈沖進行調(diào)制斜率估計的性能差異,即脈沖積累的效果。圖4 給出了單個脈沖和NP=10 個脈沖兩種情況下調(diào)制斜率估計方差的比值(設(shè)這個比值為γ)與脈沖信噪比的關(guān)系,其中縱坐標表示10 lg(γ)??梢钥闯霎斝旁氡容^高時, γ基本不隨信噪比變化,實際上此時γ應(yīng)該等于NP,即NP 個脈沖聯(lián)合估計相對單個脈沖的估計使得調(diào)制斜率的估計方差降低了NP倍;而當信噪比較低時,估計方差降低的倍數(shù)遠超過NP。比如對于4 dB的情況,相對于單個脈沖,由10 個脈沖對調(diào)制斜率進行估計使得估計方差降低了γ≈105倍,脈沖積累效果明顯。
圖4 單個和10 個脈沖條件下調(diào)制斜率估計方差的比值與脈沖信噪比的關(guān)系Fig.4 The ratio of the estimation variance of chirp rate under conditions of a single pulse and 10 pu lses vs.the SNR of pulses
本文利用多個脈沖對LFM 脈沖信號的調(diào)制斜率進行估計,提出了一種經(jīng)論證是最大似然的估計方法。該方法首先由脈沖積累估計出調(diào)制波形,大大改善了信噪比,進而準確地解出LFM 脈沖的相位,然后用最小二乘擬合的方法估計出調(diào)制斜率。提出的方法在低信噪比條件下可以使得調(diào)制斜率的估計方差達到CRLB,具有很低的信噪比門限。
需要指出的是,本文基于的信號模型假設(shè)脈沖已經(jīng)實現(xiàn)了理想的時頻對齊,但是在低信噪比條件下脈沖的時頻對齊會有較大的偏差,導(dǎo)致調(diào)制波形和調(diào)制參數(shù)估計精度的下降,對這一問題的解決還有待進一步的研究。
[ 1] Abatzoglou T J.Fast maximum likelihood joint estimation of frequency and frequency rate[J] .IEEE Transactions on Aerospace and Electronic Systems,1986,AES-22(6):708-715.
[ 2] Djuric P M, Kay S M.Parameter estimation of chirp signals[ J] .IEEE Transactions on Acoustics, Speech and Signal Processing, 1990, 38(12):2118-2126.
[ 3] Peleg S, Friedlander B.The discrete polynom ial-phase transform[ J] .IEEE Transactions on Signal Processing, 1995, 43(8):1901-1914.
[ 4] Xia X G.Discrete chirp-Fourier transform and its application to chirp rate estimation [ J] .IEEE Transactions on Signal Processing, 2000, 48(11):3122-3133.
[ 5] 金勝, 王峰, 鄧振淼,等.一種LFM 信號相位域快速高精度參數(shù)估計算法[J] .系統(tǒng)工程與電子技術(shù), 2011,33(2):264-267.
JIN Sheng, WANG Feng, DENG Zhen-miao, et al.Fast and accurate estimator on parameters of chirp signals in phase domain[ J] .Systems Engineering and Electronics, 2011, 33(2):264-267.(in Chinese)
[ 6] 周良臣,楊建宇,唐斌.一種高效的LFM 信號參數(shù)估計方法及性能分析[ J] .電子學(xué)報,2007,35(6):1128-1133.
ZHOU Liang-cheng, YANG Jian-yu, TANG B.An efficient parameter estimation and performance analysis for LFM signal[ J] .Acta Electronica Sinica, 2007, 35(6):1128-1133.(in Chinese)
[ 7] Howard S D.Estimation and correlation of radar pulse modulations for electronic support[ C]//Proceedings of the 2003 IEEE Aerospace Conference.Big Sky, MT:IEEE, 2003:2065-2072.
[ 8] Horn RA, Johnson, C R.Matrix analysis[M] .Cambridge:Cambridge University Press, 1985.
[ 9] Tretter S A.Estimating the frequency of a noisy sinusoid by linear regression[J] .IEEE Transactions on Information Theory, 1985, IT-31(6):832-835.
[10] 黃知濤, 周一宇, 姜文利.基于相對無模糊相位重構(gòu)的自動脈內(nèi)調(diào)制特性分析[ J] .通信學(xué)報, 2003, 24(4):153-160.
HUANG Zhi-tao,ZHOU Yi-yu, JIANG Wen-li.The automatic analysis of intrapu lse modulation characteristics based on the relatively non-ambiguity phase restoral[ J] .Journal on Communications,2003,24(4):153-160.(in Chinese)
[11] Kay S M.統(tǒng)計信號處理基礎(chǔ)——估計與檢測理論[M] .羅鵬飛, 張文明, 劉忠, 等, 譯.北京:電子工業(yè)出版社, 2006.
Kay S M.Fundamentals of Statistical Signal Processing:Estimation Theory and Detection Theory[M] .Translated by LUO Peng-fei,ZHANG Wen-ming,LIU Zhong,et al.Beijing:Publishing House of Electronic Industry,2006.(in Chinese)
[12] Peleg S, Porat B.The Cramér-Rao lower bound for signals with constant amplitude and polynom ial phase [ J] .IEEE Transactions on Signal Processing,1991,39(3):749-752.