顧怡鳴,宮在曉,李整林
(1. 中國科學(xué)院聲學(xué)研究所聲場聲信息國家重點(diǎn)實(shí)驗(yàn)室,北京 100190;2. 中國科學(xué)院大學(xué),北京 100190)
線性調(diào)頻(Linear Frequency Modulation, LFM)是主動聲吶常用的探測信號,通過分析 LFM 信號目標(biāo)回波可估計(jì)目標(biāo)與探測平臺的相對速度。傳統(tǒng)的LFM測速方法主要有多普勒頻移補(bǔ)償法和組合脈沖法。多普勒頻移補(bǔ)償法利用長脈寬 LFM 信號的多普勒敏感性,對探測信號進(jìn)行多普勒頻移搜索來估計(jì)相對速度[1]。組合脈沖法利用LFM信號模糊函數(shù)的時(shí)延-多普勒耦合性,通過正負(fù)調(diào)頻信號匹配濾波峰值的時(shí)延差估計(jì)速度,也被擴(kuò)展應(yīng)用于雙曲調(diào)頻信號[2]。但主動聲吶的發(fā)射占空比通常較小,組合脈沖增加了信號的發(fā)射周期。
分?jǐn)?shù)階傅里葉變換(Fractional Fourier Transform, FrFT)是一種時(shí)頻旋轉(zhuǎn)變換,近年來已經(jīng)被廣泛用于 LFM 信號的檢測和參數(shù)估計(jì)[3-13]。文獻(xiàn)[14]首先將 FrFT應(yīng)用于主動聲吶回波檢測,并利用海上實(shí)驗(yàn)數(shù)據(jù)證明了 FrFT在混響背景下可以獲得比匹配濾波器更高的信噪比和速度估計(jì)精度。但FrFT較高的計(jì)算復(fù)雜度限制了其實(shí)際應(yīng)用,針對該問題,學(xué)者們提出了一些快速 FrFT算法[15-21]。采樣型算法對時(shí)間和頻率進(jìn)行尺度歸一化,可以借助快速傅里葉變換(FFT)算法極大地降低計(jì)算量[15-16];特征向量分解型算法[12-13]通過設(shè)計(jì)離散傅里葉變換矩陣求解對應(yīng)的特征向量;文獻(xiàn)[21]中利用莫比烏斯反演公式證明了分?jǐn)?shù)階傅里葉級數(shù)的存在,并以此為基礎(chǔ)提出了低復(fù)雜度的離散FrFT算法。FrFT方法檢測 LFM 信號主要利用了其解調(diào)頻的性質(zhì),舍棄一些冗余的性質(zhì),對變換形式進(jìn)行簡化以降低計(jì)算量,例如二次相位變換[22-24](Quadratic Phase Transform, QPT)或簡明分?jǐn)?shù)階傅里葉變換[7-8],以及具有高次變換核的多項(xiàng)式時(shí)頻變換[25-26]。
現(xiàn)有的 FrFT估計(jì)速度方法需要對變換參數(shù)進(jìn)行搜索,在自主水下機(jī)器人(Autonomous Underwater Vehicle, AUV)等計(jì)算力有限的平臺難以實(shí)時(shí)處理。本文針對該問題提出了一種基于二次相位變換的線性調(diào)頻信號速度估計(jì)方法。QPT參數(shù)失配時(shí)能量無法聚焦,幅度譜存在一定展寬。本文推導(dǎo)了 QPT參數(shù)與信號參數(shù)失配時(shí)變換域幅度譜的解析解,分析了二次相位變換幅度譜展寬與多普勒頻移之間的關(guān)系,提出了利用單次QPT后幅度譜寬度的速度估計(jì)方法?;趩蝹€(gè)線性調(diào)頻信號進(jìn)行QPT后幅度譜的展寬信息實(shí)現(xiàn)了目標(biāo)相對速度估計(jì),避免了對變換參數(shù)的掃描,減少了一個(gè)維度的計(jì)算量。針對單次QPT結(jié)果無法確定速度符號的問題,提出了對稱參數(shù)QPT符號判斷方法。通過數(shù)值仿真驗(yàn)證了方法的有效性,并與 FrFT參數(shù)掃描方法進(jìn)行了對比。
FrFT的變換參數(shù)cotα與信號的調(diào)頻斜率k匹配時(shí),變換核中的 ( t2/2)cotα項(xiàng)對信號實(shí)現(xiàn)解調(diào)頻,變換結(jié)果的幅度取得最大值。文獻(xiàn)[14]依據(jù)該性質(zhì)給出了通過掃描 FrFT變換參數(shù)搜索最大峰值的測速方法。以長度為20 s、頻帶400~500 Hz的LFM信號為例,根據(jù)式(5)得到FrFT對LFM信號u-多普勒模糊(見圖1)。由圖1可以看出,當(dāng)變換參數(shù)與信號速度匹配時(shí)變換結(jié)果的幅度取得最大值。FrFT測速方法的原理是在變換域-多普勒平面內(nèi)搜索最大值,找到最大值對應(yīng)的變換參數(shù),進(jìn)而通過變換參數(shù)與多普勒頻移之間的關(guān)系得到目標(biāo)相對速度的估計(jì)值。
圖1 FrFT對LFM信號的u-多普勒模糊圖Fig.1 Schematic diagram of the u-Doppler ambiguity of FrFT for LFM signal
二次相位變換可以看作是 FrFT的簡化,信號x(t)的QPT定義為[25]
其中,u和κ是變換參數(shù)。當(dāng)變換參數(shù)和線性調(diào)頻信號的調(diào)頻斜率k匹配,即k =κ 時(shí),信號能量會在變換域聚集為峰值,此類情況已經(jīng)在現(xiàn)有文獻(xiàn)中廣泛討論;當(dāng)變換參數(shù)與信號參數(shù)失配時(shí)k≠κ ,信號的能量在變換域表現(xiàn)為一定展寬,已有文獻(xiàn)尚未發(fā)現(xiàn)有關(guān)變換參數(shù)與信號參數(shù)失配時(shí)特征的研究。以下給出變換參數(shù)與 LFM 信號調(diào)頻斜率不匹配情況下變換結(jié)果的解析解。
多普勒頻移后的回波信號二次相位變換可表示為
變換后信號的幅度譜由兩種形式的 Fresnel積分(如圖2所示)決定,式(10)對應(yīng)的幅度譜大部分區(qū)域?yàn)榱悖鐖D3所示,能量集中在一定的寬度,定義幅度譜寬度Bα為幅度譜峰值-3 dB對應(yīng)的u域?qū)挾?,能量集中的范圍的寬度?/p>
圖2 宗量為x的Fresnel積分Fig.2 Fresnel integral with argument x
圖3 多普勒頻移后的回波信號QPT的幅度譜示意圖Fig.3 QPT amplitude spectrum of the echo with Doppler frequency shift
以上推導(dǎo)表明,參數(shù)失配時(shí)LFM信號QPT的幅度譜寬度與信號多普勒頻移存在確定關(guān)系,利用該特征可估計(jì)目標(biāo)相對速度。若變換參數(shù)為原始信號參數(shù)時(shí)κ= k,LFM信號QPT的幅度譜的寬度為
式中:c為聲速;k和T為LFM信號的調(diào)頻斜率和脈沖寬度。只需要通過估計(jì)幅度譜寬度Bk即可得到相對速度。QPT幅度譜的寬度信息可通過閾值檢測獲取,即:
其中:d為檢測門限,根據(jù)閾值檢測幅度譜的前后沿u0和u1可估計(jì)幅度譜寬度的估計(jì)值目標(biāo)靜止時(shí)QPT幅度譜為sinc函數(shù),幅度譜寬度為。需要注意的是,式(10)并不是完美的矩形窗函數(shù),檢測門限較高時(shí)估計(jì)的寬度會偏小,需要通過修正系數(shù)進(jìn)行補(bǔ)償。
根據(jù) QPT幅度譜的寬度可以估計(jì)相對速度 v的大小,相對速度的符號則可以根據(jù)幅度譜前沿與的關(guān)系來確定:
假設(shè) LFM 信號長度為 20 s,頻帶 400~500 Hz,對相對速度為 0和 30 m·s-1的信號進(jìn)行QPT后得到的幅度譜分別如圖4(a)和圖4(b)所示,可以看出,當(dāng)相對速度為0時(shí),信號能量在變換域聚焦,信號QPT幅度譜幾乎沒有展寬。相對速度為30 m·s-1時(shí),信號QPT幅度譜展寬明顯,并且與式(11)理論分析的寬度一致。圖4(c)為QPT對LFM信號u-多普勒模糊,可以看出,隨著相對速度的增加,QPT幅度譜逐漸展寬。同時(shí),當(dāng)相對速度為正時(shí),幅度譜分布位于f0- kτ0右側(cè),相對速度為負(fù)時(shí)位于f0- kτ0左側(cè)。由于信號的能量一定,QPT幅度譜的幅度隨展寬而降低。經(jīng)過閾值檢測后的寬度圖像如圖4(d)所示,根據(jù)寬度可以進(jìn)一步估計(jì)目標(biāo)相對速度。需要注意,單次QPT方法依賴式(14)判定速度的正負(fù),實(shí)際應(yīng)用中不能準(zhǔn)確獲取時(shí)延τ0,因此單次 QPT提供的信息在判斷速度的正負(fù)方面存在缺陷。
圖4 單次QPT的速度估計(jì)方法Fig.4 Velocity estimation method with single QPT
單次 QPT譜寬度僅提供了相對速度的大小信息,不能確定相對速度的符號,這里給出對稱參數(shù)QPT速度符號判斷方法。圖 5為符號判斷方法的示意圖,給定一個(gè)初始寬度B0,若信號不存在多普勒頻移時(shí),B0對應(yīng)的兩組變換參數(shù)分別為分別對回波信號進(jìn)行變換參數(shù)為κ+和κ-的 QPT,分別可以得到兩組變換的幅度譜寬度Bκ+和Bκ-,兩次對稱參數(shù) QPT幅度譜寬度為
圖5 對稱參數(shù)QPT速度符號判斷示意圖Fig.5 Schematic diagram of speed symbol judgment with symmetrical parameter QPT
由式(15)可以看出,當(dāng)信號相對速度為正時(shí),伸縮因子η> 0,Bκ+< Bκ-;當(dāng)信號相對速度為負(fù)時(shí),伸縮因子η<0,Bκ+< Bκ-。通過對比兩次對稱參數(shù)QPT幅度譜寬度可判斷相對速度的符號。
本節(jié)通過數(shù)值仿真實(shí)驗(yàn)驗(yàn)證 QPT速度估計(jì)方法的有效性,并將QPT方法與FrFT參數(shù)掃描方法的性能進(jìn)行對比。仿真環(huán)境為無限空間,海水聲速為 1 500 m·s-1。仿真實(shí)驗(yàn)的部署如圖 6(a)所示,藍(lán)圈靜止位于坐標(biāo)原點(diǎn)(0, 0)的收發(fā)合置的主動聲吶平臺,黑圈代表目標(biāo)起始位置(-10 km, 5 km),虛線代表目標(biāo)自東向西以速度 200 m·s-1勻速運(yùn)動的航線。根據(jù)目標(biāo)運(yùn)動速度和相對平臺的方位可以得到目標(biāo)與平臺的相對速度如圖6(b)所示,正速度表示目標(biāo)接近平臺。探測信號為長度 20 s、頻帶 400~500 Hz的LFM信號,信號發(fā)射周期為60 s。在整個(gè)20 km的運(yùn)動距離內(nèi),總計(jì)發(fā)射16個(gè)脈沖。
圖6 仿真實(shí)驗(yàn)部署Fig.6 Deployment of simulation experiment
仿真接收信號的采樣頻率為5 000 Hz,假設(shè)背景噪聲為高斯白噪聲,回波信號頻帶內(nèi)的信噪比為20 dB。仿真回波信號直接經(jīng)過QPT后的幅度譜如圖7(a)所示,回波QPT后幅度譜能量集中的位置與回波的時(shí)延有關(guān)。為了更加直觀地表現(xiàn)相對速度對 QPT幅度譜的影響,將 u坐標(biāo)變換為u′坐標(biāo),結(jié)果如圖 7(b)所示。對比真實(shí)相對速度(圖 6(b))可以觀察到 QPT幅度譜寬度與相對速度大小的對應(yīng)關(guān)系,設(shè)置檢測門限為,進(jìn)行寬度檢測并依據(jù)式(12)的函數(shù)關(guān)系即可估計(jì)目標(biāo)的相對速度。圖8為對稱參數(shù)QPT得到的幅度譜寬度檢測結(jié)果,其中藍(lán)線和紅線分別對應(yīng)正參數(shù)和負(fù)參數(shù)變換對應(yīng)的幅度譜寬度,通過對比兩者之間的大小可以判斷速度的正負(fù)。
圖7 仿真數(shù)據(jù)QPT譜寬度估計(jì)Fig.7 Estimation of the of QPT spectrum width of simulation data
圖8 對稱參數(shù)QPT幅度譜寬度估計(jì)結(jié)果Fig.8 Estimation of the amplitude spectrum width by using QPT of symmetric parameters
在實(shí)際應(yīng)用中,F(xiàn)rFT測速方法不需要對0~π所有的變換參數(shù)進(jìn)行搜索,通過速度和變換參數(shù)的關(guān)系限定參數(shù)掃描的范圍可以降低計(jì)算量。速度掃描范圍確定后,計(jì)算量與掃描的分辨率相關(guān),掃描分辨率也直接決定了 FrFT方法的速度估計(jì)精度。仿真中,F(xiàn)rFT方法的速度掃描范圍限定在-30~30 m·s-1范圍內(nèi),掃描間隔分別采用3、1和0.3 m·s-1,對應(yīng)的FrFT點(diǎn)數(shù)Nα分別為10、60和100點(diǎn)。分別利用 FrFT和 QPT兩種方法對相對速度進(jìn)行估計(jì),得到速度誤差如圖9所示。由圖9可以看出,F(xiàn)rFT方法的速度估計(jì)精度主要與參數(shù)掃描的分辨率有關(guān),仿真數(shù)據(jù)單次QPT方法的速度估計(jì)精度基本可以保持在1 m·s-1以內(nèi),與60點(diǎn)的FrFT方法精度相當(dāng)。
圖9 FrFT參數(shù)掃描方法與單次QPT方法速度估計(jì)誤差對比Fig.9 Comparison of velocity estimation error between FRFT parameter scanning method and QPT method
FrFT類算法的計(jì)算復(fù)雜度分析在許多文獻(xiàn)中已有討論[7-8,15,26],盡管QPT簡化了變換核,但計(jì)算復(fù)雜度仍與 FrFT的量級相當(dāng)。忽略計(jì)算復(fù)雜度的低階項(xiàng),兩種變換基于 FFT的快速算法計(jì)算復(fù)雜度量級均為O(NlgN),其中N表示單次處理的時(shí)域數(shù)據(jù)點(diǎn)數(shù)。遠(yuǎn)程主動探測回波處理的時(shí)間窗較長,處理點(diǎn)數(shù)的增加使算法的時(shí)間復(fù)雜度急劇上升。以第2節(jié)仿真實(shí)驗(yàn)中采用的參數(shù)(窗長度40 s,采樣點(diǎn)數(shù)20 000點(diǎn))為例,評估不同參數(shù)掃描量的計(jì)算量,結(jié)果如圖10所示,Nα=1對應(yīng)的虛線表示單次QPT的計(jì)算量。如果只估計(jì)相對速度大小,QPT計(jì)算量僅為精度相當(dāng)?shù)腇rFT方法的1/60;增加判斷相對速度符號,則計(jì)算量為FrFT方法的1/30。
圖10 不同參數(shù)掃描量的計(jì)算量估計(jì)Fig.10 The computational costs of different scanning parameters
參數(shù)失配時(shí)能量在變換域并非聚焦于一點(diǎn),因此 QPT方法的抗噪聲性能差于 FrFT參數(shù)掃描方法。假設(shè)回波信號的帶內(nèi)信噪比為 20 dB,不同相對速度回波信號的QPT幅度譜波形如圖11所示,相對速度為0時(shí),信號能量完全聚焦,背景起伏很低,隨著相對速度提高,背景起伏逐漸升高。定義峰值信噪比為輸出信噪比,通過 20次蒙特卡洛仿真估計(jì)QPT和FrFT兩種方法的輸出信噪比隨相對速度的變化,結(jié)果如圖12所示,藍(lán)線和紅線分別表示QPT方法和FrFT方法的輸出信噪比??梢钥闯?,QPT方法的輸出信噪比在相對速度為0時(shí)輸出信噪比與FrFT方法相同,隨著目標(biāo)相對速度提高,信號能量在變換域逐漸分散,輸出信噪比逐漸降低。上述數(shù)值分析表明,QPT方法的抗噪聲性能比FrFT方法差,該方法計(jì)算效率提高的代價(jià)是輸出信噪比的降低。雖然經(jīng)典的 FFT譜分析方法通過分析信號頻偏也可獲得目標(biāo)速度,但主動聲吶接收的目標(biāo)回波信號通常信噪比很低,F(xiàn)FT完全無法聚焦信號的能量,輸出信噪比遠(yuǎn)低于以上兩種方法。
圖11 不同相對速度QPT幅度譜的波形Fig.11 Waveforms of QPT spectrums corresponding to different relative velocities
圖12 QPT方法和FrFT方法的輸出信噪比對比Fig.12 Comparison of output signal to noise ratio between QPT method and FrFT method
本文提出了一種基于二次相位變換的主動聲吶線性調(diào)頻信號的速度估計(jì)方法,通過減少變換次數(shù)降低了 FrFT類算法速度估計(jì)的計(jì)算量。單次QPT速度估計(jì)方法利用LFM變換后幅度譜寬度與信號多普勒頻移的關(guān)系對速度進(jìn)行估計(jì),速度符號的判斷依賴準(zhǔn)確的到達(dá)時(shí)刻。針對實(shí)際情況下難以準(zhǔn)確估計(jì)時(shí)延的問題,提出了對稱參數(shù)QPT速度符號判斷方法,通過對比對稱參數(shù)變換的幅度譜寬度判斷速度符號。由于不需要對變換參數(shù)進(jìn)行搜索,相較于 FrFT參數(shù)掃描的方法,所提速度估計(jì)方法可以降低一個(gè)維度的計(jì)算量。仿真實(shí)驗(yàn)結(jié)果表明,該方法的速度估計(jì)精度可以達(dá)到1 m·s-1,如果只估計(jì)相對速度大小,計(jì)算量僅為精度相當(dāng)?shù)腇rFT參數(shù)掃描方法的1/60。但該方法計(jì)算效率提高的代價(jià)是輸出信噪比的降低,在實(shí)際應(yīng)用中需要進(jìn)行權(quán)衡??紤]到復(fù)雜海洋環(huán)境條件下主動聲吶接收回波信號存在畸變,本文的方法還有待復(fù)雜信道下真實(shí)主動聲吶實(shí)驗(yàn)數(shù)據(jù)的驗(yàn)證。