基于傅里葉梅林變換的寬帶非線(xiàn)性調(diào)頻信號(hào)時(shí)差/尺度差估計(jì)方法

張 亮， 于洪波， 孫殿星， 程樹(shù)強(qiáng)， 張翔宇， 王國(guó)宏

(1.海軍航空大學(xué)，山東 煙臺(tái) 264000； 2.中國(guó)人民解放軍94326部隊(duì)，濟(jì)南 250000；3.中國(guó)人民解放軍94201部隊(duì)，濟(jì)南 250000)

0 引言

電子戰(zhàn)是信息化戰(zhàn)爭(zhēng)的重要組成部分，對(duì)目標(biāo)的精確定位是電子戰(zhàn)的一項(xiàng)重要內(nèi)容。目標(biāo)定位系統(tǒng)通常可分為有源定位系統(tǒng)和無(wú)源定位系統(tǒng)，典型的有源定位系統(tǒng)為雷達(dá)，由于需要對(duì)外輻射電磁信號(hào)以獲取目標(biāo)回波，戰(zhàn)場(chǎng)上容易被敵方偵察而遭受反輻射武器打擊。相較而言，無(wú)源定位系統(tǒng)自身不輻射電磁信號(hào)，隱蔽性好，典型體制為雙星無(wú)源定位[1]。寬帶條件下，雙星無(wú)源定位系統(tǒng)根據(jù)兩路接收信號(hào)估計(jì)時(shí)差(TDOA)和尺度差(SDOA)，再利用衛(wèi)星、目標(biāo)和地球的幾何關(guān)系對(duì)輻射源進(jìn)行精確定位。為估計(jì)時(shí)差/尺度差，傳統(tǒng)方法[2-3]是計(jì)算兩路接收信號(hào)寬帶互模糊函數(shù)(Wideband Cross Ambiguity Function,WBCAF)，該方法能夠同時(shí)估計(jì)2個(gè)參數(shù)，且適用于多種源信號(hào)形式，而存在的問(wèn)題是計(jì)算量大。WBCAF計(jì)算過(guò)程包含兩步，首先計(jì)算一路信號(hào)的尺度版本，然后將其與另一路信號(hào)進(jìn)行相關(guān)處理，由于第2步可利用FFT，IFFT在頻域快速實(shí)現(xiàn)，因此WBCAF快速實(shí)現(xiàn)的難點(diǎn)為第1步，針對(duì)該問(wèn)題，文獻(xiàn)[2]提出多速率采樣重構(gòu)濾波方法，但是當(dāng)尺度因子為多位小數(shù)時(shí)，該方法無(wú)法獲得精確的尺度版本。部分研究考慮僅在變換域計(jì)算WBCAF，以規(guī)避接收信號(hào)尺度版本獲取難題，其中，文獻(xiàn)[3]利用小波變換內(nèi)積定理，提出基于交互小波變換的WBCAF實(shí)現(xiàn)方法，但母小波的選取對(duì)計(jì)算結(jié)果影響較大，文獻(xiàn)[4]進(jìn)一步提出基于梅林變換的WBCAF實(shí)現(xiàn)方法，但梅林變換的數(shù)值計(jì)算方法有多種[5-9]，所用方法計(jì)算精度不高。另外，如果源信號(hào)形式已知，時(shí)差/尺度差估計(jì)問(wèn)題會(huì)變得相對(duì)容易，例如針對(duì)線(xiàn)性調(diào)頻(LFM)信號(hào)，文獻(xiàn)[10-11]利用分?jǐn)?shù)階傅里葉變換(Fractional Fourier Transform,FRFT)測(cè)量?jī)陕方邮招盘?hào)調(diào)頻斜率，根據(jù)兩路信號(hào)調(diào)頻斜率比例關(guān)系估計(jì)尺度差，再對(duì)第1路信號(hào)“拉伸”或者“壓縮”，與第2路信號(hào)進(jìn)行相關(guān)處理估計(jì)時(shí)差。該改進(jìn)方法通過(guò)兩步級(jí)聯(lián)處理降低了運(yùn)算整體復(fù)雜度，低信噪比(SNR)條件下2個(gè)參數(shù)估計(jì)均方誤差逼近克拉美羅下界(Cramer-Rao Lower Bound,CRLB)，遺憾的是當(dāng)源信號(hào)為非線(xiàn)性調(diào)頻(NLFM)信號(hào)時(shí)，該方法不適用。

針對(duì)上述問(wèn)題，本文引入時(shí)間尺度(TS)、尺度估計(jì)(SE)基本概念，利用兩路信號(hào)頻譜包絡(luò)是彼此尺度版本的特點(diǎn)，提出基于傅里葉梅林變換的寬帶非線(xiàn)性調(diào)頻信號(hào)時(shí)差/尺度差估計(jì)方法。

1 傳統(tǒng)與改進(jìn)方法

寬帶條件下，利用雙星無(wú)源定位系統(tǒng)對(duì)地面一個(gè)固定輻射源進(jìn)行定位，以第1路信號(hào)為參考基準(zhǔn)，則2個(gè)獨(dú)立接收機(jī)接收到的信號(hào)可表示為[10-11]

(1)

式中：x1(t)，x2(t)分別為第1路和第2路接收信號(hào)；A1，A2為源信號(hào)s(t)到達(dá)2個(gè)接收機(jī)的路徑增益；Δt為兩路信號(hào)時(shí)差(TDOA)；Δα為尺度差(SDOA)，與衛(wèi)星相對(duì)輻射源的徑向速度及光速有關(guān)；n1(t)，n2(t)分別為2個(gè)接收機(jī)基底復(fù)高斯白噪聲。輻射源相對(duì)衛(wèi)星的徑向速度體現(xiàn)為接收信號(hào)尺度的變化，傳統(tǒng)估計(jì)方法[2-3]是利用寬帶互模糊函數(shù)(WBCAF)實(shí)現(xiàn)。連續(xù)信號(hào)z(t)與x(t)的WBCAF定義為

(2)

式中：γ為WBCAF尺度因子；τ為時(shí)移。令z(t)=x2(t)，x(t)=x1(t)，易知，當(dāng)γ=1+Δα，τ=Δt時(shí)，|Xzx(γ,τ)|

取得最大值，峰值搜索可得到時(shí)差和尺度差(需減1)。WBCAF可表示成卷積形式頻域快速實(shí)現(xiàn)，即

Xzx(γ,τ)=z(t)?x*(-γt)=F-1{F[z(t)]F*[x(γt)]}

(3)

式中：?為卷積符號(hào)；F[·]，F(xiàn)-1[·]分別為傅里葉變換與逆變換表示符號(hào)。由式(3)可知，為計(jì)算兩路接收信號(hào)WBCAF，需要設(shè)定尺度因子搜索區(qū)間和間隔，然后根據(jù)不同的尺度因子對(duì)第1路信號(hào)進(jìn)行“拉伸”或“壓縮”，搜索區(qū)間的擴(kuò)大會(huì)導(dǎo)致運(yùn)算量成倍增加。由式(1)可知，如果源信號(hào)為單頻信號(hào)，則尺度上的差異反映為接收信號(hào)頻率上的不同，如果源信號(hào)為L(zhǎng)FM信號(hào)，則反映為調(diào)頻斜率上的不同，通過(guò)估計(jì)兩路信號(hào)頻率或者調(diào)頻斜率可得SDOA，然后根據(jù)估計(jì)的SDOA對(duì)x1(t)進(jìn)行“拉伸”或者“壓縮”，再與x2(t)進(jìn)行時(shí)域相關(guān)處理可得TDOA，該改進(jìn)方法雖然能夠顯著減少運(yùn)算量，但當(dāng)源信號(hào)為NLFM信號(hào)時(shí)則不適用。雷達(dá)、通信系統(tǒng)常用NLFM信號(hào)包括二次調(diào)頻(QFM)、指數(shù)調(diào)頻(EFM)、雙曲調(diào)頻(HFM)和正弦調(diào)頻(SFM)等[12]，復(fù)數(shù)形式依次為

(4)

(5)

(6)

(7)

式中：T為時(shí)寬，f1,f2分別為4種信號(hào)的起始頻率和截止頻率;T0為HFM延時(shí)時(shí)間。對(duì)式(4)～(7)相位求導(dǎo)可得QFM,EFM,HFM和SFM瞬時(shí)頻率分別為

(8)

(9)

(10)

(11)

2 時(shí)間尺度與尺度估計(jì)

為介紹所提方法，本章引入2個(gè)基本概念，即時(shí)間尺度(TS)和尺度估計(jì)(SE)，然后利用梅林變換(MT)的特例尺度變換(ST)給出TS和SE快速實(shí)現(xiàn)方案。

2.1 基本概念

所謂時(shí)間尺度，即連續(xù)信號(hào)x(t)到y(tǒng)(t)的映射過(guò)程[13]，可表示為

(12)

(13)

由式(12)可知，當(dāng)0<α<1時(shí)，信號(hào)時(shí)域擴(kuò)張，當(dāng)α>1時(shí)，信號(hào)時(shí)域收縮。根據(jù)傅里葉變換尺度性質(zhì)可知

(14)

式中，Y(f)，X(f)分別為y(t)，x(t)的傅里葉變換，f為頻率自變量。可以看出，對(duì)信號(hào)進(jìn)行TS還會(huì)導(dǎo)致帶寬成比例減小(0<α<1)或增大(α>1)，但信號(hào)時(shí)寬-帶寬積不變[14]。設(shè)連續(xù)信號(hào)z(t)=TSα1[x(t)]，α1為尺度因子。所謂尺度估計(jì)，即z(t)和x(t)已知時(shí)對(duì)α1的估計(jì)，可表示為

(15)

(16)

式中：◇為尺度互相關(guān)符號(hào)；*表示共軛轉(zhuǎn)置。由式(2)可知，當(dāng)τ=0時(shí)，WBCAF與尺度互相關(guān)函數(shù)近似，區(qū)別在于不同的積分區(qū)間。

2.2 實(shí)現(xiàn)方法

尺度變換(ST)是梅林變換(MT)復(fù)自變量實(shí)部取0.5的特例，表達(dá)式[13]為

(17)

式中：S{·}為ST表示符號(hào)；Df(c)為信號(hào)f(t)的ST，c為尺度自變量。逆尺度變換(Inverse Scale Transform,IST)為

(18)

式中,S-1{·}為IST表示符號(hào)。令t=eu，代入式(17)，得到

(19)

(20)

IST同樣可利用IFFT快速實(shí)現(xiàn)，得到快速逆尺度變換(Inverse Fast Scale Transform,IFST)，計(jì)算過(guò)程與FST相反，先計(jì)算Df(c)的逆傅里葉變換，再進(jìn)行對(duì)數(shù)采樣。為解決TS快速實(shí)現(xiàn)問(wèn)題，引用ST尺度不變性，即

Dz(c)=eicln α1Dx(c)

(21)

式中,Dx(c),Dz(c)分別為x(t),z(t)的ST，兩者包絡(luò)相同、相位不同。計(jì)算式(21)的IST，易得

z(t)=S-1{eicln α1Dx(c);t} 。

(22)

利用式(22)可得x(t)在某一尺度因子下的TS后信號(hào)。SE同樣利用ST快速實(shí)現(xiàn)[15]，具體為

(23)

3 基于傅里葉梅林變換的時(shí)差/尺度差估計(jì)

3.1 算法原理

為介紹本文所提方法，暫不考慮接收機(jī)噪聲，重寫(xiě)式(1)為

(24)

利用傅里葉變換時(shí)移不變性消除時(shí)差影響，容易得到

(25)

式中，X1(f)，X2(f)分別為x1(t)，x2(t)的傅里葉變換。易知

(26)

可以看出，|X2(f)|實(shí)際為|X1(f)|的尺度版本，尺度因子為1/(1+Δα)，當(dāng)|X1(f)|和|X2(f)|均已知時(shí)，對(duì)1/(1+Δα)的估計(jì)屬SE概念，利用式(23)可估計(jì)1/(1+Δα)，進(jìn)而得到尺度差Δα，然后對(duì)x1(t)進(jìn)行時(shí)間尺度(尺度因子取1+Δα)，與x2(t)進(jìn)行時(shí)域相關(guān)處理可進(jìn)一步估計(jì)TDOA。為直觀(guān)顯示，圖1給出了本文所提方法實(shí)現(xiàn)流程。

圖1 時(shí)差/尺度差估計(jì)流程

3.2 計(jì)算復(fù)雜度

首先，分析時(shí)間尺度的計(jì)算復(fù)雜度。TS實(shí)現(xiàn)過(guò)程包含“FST+相位修正+IFST”3個(gè)環(huán)節(jié)，計(jì)算復(fù)雜度為O[(NlnN)lb(NlnN)]，總共需要進(jìn)行3NlnN+NlnN·lb(NlnN)次復(fù)數(shù)乘法。其次，分析SE計(jì)算復(fù)雜度。SE實(shí)現(xiàn)的關(guān)鍵在于計(jì)算尺度互相關(guān)函數(shù)，具體需要進(jìn)行2次FST，1次復(fù)數(shù)乘積和1次IFST，計(jì)算復(fù)雜度同樣為O[(NlnN)lb(NlnN)]，總共需要進(jìn)行4NlnN+1.5·NlnNlb(NlnN)次復(fù)數(shù)乘法。再次，分析傳統(tǒng)基于寬帶互模糊函數(shù)的時(shí)差/尺度差估計(jì)方法計(jì)算復(fù)雜度。對(duì)于1個(gè)尺度因子，為計(jì)算式(3)需要進(jìn)行1次TS，2次FFT，1次 IFFT和1次點(diǎn)乘，總共需要進(jìn)行MN+1.5MN·lbN+3MNlnN+MNlnNlb(NlnN)次復(fù)乘運(yùn)算，計(jì)算復(fù)雜度為O[(MNlnN)lb(NlnN)]，M為尺度因子取值個(gè)數(shù)。最后，分析所提方法計(jì)算復(fù)雜度。為估計(jì)尺度差需要進(jìn)行2次FFT和1次SE，而估計(jì)時(shí)差需要進(jìn)行1次TS和1次時(shí)域相關(guān)處理(2次FFT，1次IFFT和1次點(diǎn)乘)，因此，總共需要進(jìn)行N+2.5NlbN+7NlnN+2.5NlnN·lb(NlnN)次復(fù)數(shù)乘法運(yùn)算，計(jì)算復(fù)雜度為O[(NlnN)·lb(NlnN)]。為直觀(guān)顯示，設(shè)M分別為10，100和1000，圖2給出了復(fù)乘次數(shù)隨信號(hào)采樣點(diǎn)數(shù)變化曲線(xiàn)。所提方法計(jì)算量明顯低于傳統(tǒng)方法。上述分析未考慮峰值搜索環(huán)節(jié)，由于傳統(tǒng)方法需要進(jìn)行二維搜索，實(shí)際運(yùn)算量會(huì)更大。

圖2 復(fù)乘次數(shù)隨信號(hào)采樣點(diǎn)數(shù)變化曲線(xiàn)

4 仿真結(jié)果與分析

4.1 參數(shù)設(shè)置

輻射源發(fā)射N(xiāo)LFM信號(hào)，時(shí)寬2.5 μs，帶寬250 MHz，采樣頻率1 GHz，雙星無(wú)源定位系統(tǒng)接收的兩路NLFM信號(hào)時(shí)差/尺度差分別為3.75 μs和0.1，QFM，EFM，HFM和SFM起始頻率和截止頻率分別為10 MHz和260 MHz，HFM延遲時(shí)間為1個(gè)時(shí)域采樣間隔。4種信號(hào)時(shí)頻分布如圖3所示。

圖3 NLFM信號(hào)時(shí)頻分布

4.2 可行性分析

以SFM信號(hào)時(shí)差/尺度差估計(jì)為例，利用2個(gè)仿真試驗(yàn)依次對(duì)第1章所提WBCAF實(shí)現(xiàn)方法和第3章所提基于傅里葉梅林變換的時(shí)差/尺度差估計(jì)方法的可行性進(jìn)行驗(yàn)證。為展示細(xì)節(jié)，接收信號(hào)中均不添加噪聲(4.3節(jié)將進(jìn)行抗噪效能分析)。

4.2.1 仿真試驗(yàn)1

設(shè)接收信號(hào)時(shí)長(zhǎng)為10 μs，兩路接收信號(hào)時(shí)域波形和頻譜如圖4所示。第2路信號(hào)中的SFM信號(hào)相比第1路信號(hào)中的SFM信號(hào)存在時(shí)域收縮和頻譜擴(kuò)張現(xiàn)象，與2.1節(jié)中的TS概念一致。

圖4 接收信號(hào)時(shí)域波形和頻譜

設(shè)WBCAF尺度因子取值0.5～1.5，間隔0.01，利用第1章方法，計(jì)算兩路接收信號(hào)寬帶互模糊函數(shù)，結(jié)果如圖5所示。

圖5 接收信號(hào)寬帶互模糊函數(shù)

由圖5可以看出，兩路信號(hào)WBCAF呈明顯沖擊特性，峰值搜索可得時(shí)差和尺度差分別為3.745 μs，0.096，與仿真所設(shè)3.75 μs，0.1基本一致。

需要注意的是，寬帶互模糊函數(shù)X軸實(shí)際為尺度因子，與尺度差存在減1的關(guān)系，圖5尺度差坐標(biāo)已進(jìn)行減1處理。

4.2.2 仿真試驗(yàn)2

按照?qǐng)D1流程，對(duì)所提基于傅里葉梅林變換的時(shí)差/尺度差估計(jì)方法可行性進(jìn)行驗(yàn)證。首先，利用2.2節(jié)方法計(jì)算圖4(b)兩路信號(hào)頻譜尺度互相關(guān)函數(shù)(SCF)結(jié)果,如圖6所示。可以看出，對(duì)頻譜取包絡(luò)時(shí)，估計(jì)得到的頻域尺度因子為0.909 2，對(duì)應(yīng)時(shí)域尺度因子為1.099 9，尺度差為0.099 9，與仿真所設(shè)尺度差一致，而不對(duì)頻譜取包絡(luò)，估計(jì)結(jié)果完全錯(cuò)誤。

圖6 接收信號(hào)頻譜尺度互相關(guān)函數(shù)

根據(jù)估計(jì)的時(shí)域尺度因子，對(duì)第1路信號(hào)進(jìn)行時(shí)間尺度，消除尺度差對(duì)時(shí)差估計(jì)的影響，接收信號(hào)時(shí)域波形和頻譜如圖7所示。

圖7 時(shí)間尺度后的接收信號(hào)時(shí)域波形和頻譜

兩路接收信號(hào)中的SFM時(shí)寬、帶寬基本一致，不再存在收縮、擴(kuò)展現(xiàn)象。以第1路信號(hào)為匹配信號(hào)，對(duì)第2路信號(hào)進(jìn)行相關(guān)處理，結(jié)果如圖8所示，峰值搜索估計(jì)得到時(shí)差為3.75 μs，與仿真所設(shè)參數(shù)一致。

圖8 相關(guān)處理結(jié)果

4.3 抗噪效能分析

對(duì)傳統(tǒng)基于寬帶互模糊函數(shù)時(shí)差/尺度差估計(jì)方法[2-3]和所提基于傅里葉梅林變換的時(shí)差/尺度差估計(jì)方法的抗噪效能進(jìn)行驗(yàn)證，使用的評(píng)估指標(biāo)為平均相對(duì)誤差(MRE)[16-17]，同時(shí),為驗(yàn)證所提算法的適用范圍，選取LFM信號(hào)(脈寬、初始頻率、截止頻率同4.1節(jié)),4種NLFM信號(hào)(QFM,EFM,HFM,SFM)和零均值高斯白噪聲信號(hào)(脈寬同4.1節(jié))作為源信號(hào)。SNR取值-10～10 dB，間隔1 dB，運(yùn)行蒙特卡羅仿真500次，得到尺度差估計(jì)平均相對(duì)誤差隨信噪比變化曲線(xiàn)如圖9所示，圖10為時(shí)差估計(jì)平均相對(duì)誤差隨信噪比變化曲線(xiàn)。

圖9 尺度差估計(jì)平均相對(duì)誤差隨信噪比變化曲線(xiàn)

圖10 時(shí)差估計(jì)平均相對(duì)誤差隨信噪比變化曲線(xiàn)

傳統(tǒng)方法與所提方法對(duì)不同信號(hào)形式時(shí)差/尺度差估計(jì)效能不盡相同，總體上看，傳統(tǒng)方法更適合低SNR條件，當(dāng)SNR大于-4 dB時(shí)，對(duì)LFM，QFM，EFM和SFM這4種信號(hào)形式尺度差估計(jì)MRE小于5%，時(shí)差估計(jì)MRE小于1%，而所提方法則要求SNR大于5 dB，需要注意的是，高SNR條件下所提方法估計(jì)精度要高于傳統(tǒng)方法，原因是計(jì)算WBCAF需人為設(shè)定尺度因子取值范圍和間隔(本節(jié)仿真中，尺度因子取0.5～1.5，間隔0.01)，當(dāng)取值范圍一定時(shí)，降低間隔能夠提高估計(jì)精度，但同時(shí)也意味著更大的計(jì)算量，所提方法不存在此問(wèn)題。

同時(shí)可以看出，傳統(tǒng)方法對(duì)HFM尺度差估計(jì)效果極不理想，不同SNR條件下尺度差估計(jì)MRE大于13%，估計(jì)誤差明顯高于所提方法，分析原因是HFM具備“多普勒不變性”，其寬帶模糊函數(shù)脊線(xiàn)近似為線(xiàn)性，脊線(xiàn)幅度較為接近[18-19]，峰值搜索無(wú)法準(zhǔn)確估計(jì)時(shí)差/尺度差，而所提方法是對(duì)信號(hào)頻譜包絡(luò)進(jìn)行尺度估計(jì)，實(shí)現(xiàn)原理不同。傳統(tǒng)方法與所提方法對(duì)噪聲信號(hào)尺度差估計(jì)MRE小于5%的臨界SNR分別為2 dB和3 dB，時(shí)差估計(jì)MRE小于1%的臨界SNR分別為1 dB和2 dB，兩種方法估計(jì)效能相差不大。

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

針對(duì)雙星無(wú)源定位系統(tǒng)中的時(shí)差/尺度差估計(jì)問(wèn)題，提出基于傅里葉梅林變換的寬帶非線(xiàn)性調(diào)頻信號(hào)時(shí)差/尺度差估計(jì)方法，核心思想是利用傅里葉變換“時(shí)移不變性”、梅林尺度變換“尺度不變性”對(duì)時(shí)差/尺度差進(jìn)行級(jí)聯(lián)估計(jì)。仿真表明，當(dāng)SNR大于5 dB時(shí)，所提方法能夠準(zhǔn)確估計(jì)4種寬帶非線(xiàn)性調(diào)頻信號(hào)時(shí)差和尺度差，時(shí)差估計(jì)平均相對(duì)誤差小于1%，尺度差估計(jì)平均相對(duì)誤差小于5%；相比傳統(tǒng)基于寬帶互模糊函數(shù)的時(shí)差/尺度差估計(jì)方法，所提方法計(jì)算復(fù)雜度更低；相比基于FRFT的時(shí)差/尺度差估計(jì)方法，所提方法可適用LFM，QFM，EFM，HFM和SFM，甚至噪聲脈沖，理論上，所提方法可適用任意信號(hào)形式，限于篇幅不做進(jìn)一步仿真。所提方法要求SNR大于5 dB，低SNR條件下的時(shí)差/尺度差估計(jì)問(wèn)題是下一步研究工作的重點(diǎn)。