趙曉靜 高大治 孫 凱 李玉征 郭豆豆

(中國(guó)海洋大學(xué)信息科學(xué)與工程學(xué)部 青島 266100)

0 引言

線譜具有能量高、傳播距離遠(yuǎn)等特點(diǎn)[1]，被廣泛應(yīng)用于目標(biāo)運(yùn)動(dòng)參數(shù)估計(jì)。線譜多普勒頻移包含著目標(biāo)的運(yùn)動(dòng)速度、到達(dá)最近點(diǎn)時(shí)間、最近距離等運(yùn)動(dòng)參數(shù)信息。1994年，F(xiàn)erguson 等[2]提出了利用線譜多普勒頻移估計(jì)目標(biāo)運(yùn)動(dòng)參數(shù)的方法，該方法利用時(shí)頻分析獲取瞬時(shí)頻率，與多普勒頻移模型進(jìn)行擬合，獲取目標(biāo)的運(yùn)動(dòng)參數(shù)。而后，Quinn[3]推導(dǎo)了接收器和運(yùn)動(dòng)聲源運(yùn)動(dòng)方向不共線時(shí)多普勒頻移表達(dá)式。鄒紅星等[4]提出了Dopplerlet 變換，該變換可將多普勒頻移去除。吳國(guó)清等[5]利用鄒紅星等提出的Dopplerlet 變換和匹配投影算法，提取線譜多普勒中包含的目標(biāo)速度和最近距離等運(yùn)動(dòng)信息，而后又提出了一種利用多線譜的多普勒頻移估計(jì)運(yùn)動(dòng)目標(biāo)最近距離的方法[6]，先用Wigner-Ville分布提取瞬時(shí)頻率，再定義多普勒頻移的基函數(shù)，實(shí)現(xiàn)目標(biāo)的運(yùn)動(dòng)參數(shù)估計(jì)。該方法的缺點(diǎn)是需要大量的先驗(yàn)知識(shí)。高偉等[7]提出了一種基于目標(biāo)噪聲輻射強(qiáng)度和線譜多普勒頻移估計(jì)目標(biāo)最近距離的方法，該方法不需要大量的先驗(yàn)知識(shí)，且搜索范圍易確定。徐靈基[8]針對(duì)水下目標(biāo)線譜多普勒進(jìn)行深入研究，通過(guò)引入了一種時(shí)頻變換–匹配Wigner變換，獲取多普勒頻移，從而實(shí)現(xiàn)目標(biāo)運(yùn)動(dòng)參數(shù)估計(jì)。Liang[9]等提出了Doppler-Chirplet 變換，該變換直接利用觀測(cè)到的瞬時(shí)頻率與源速度的關(guān)系作為Chirplet 變換的核心，可直接估測(cè)運(yùn)動(dòng)目標(biāo)的速度。同時(shí)，還有多尺度Chirplet 路徑追蹤方法[10]等。劉凱悅等[11]提出了一種基于線譜多普勒的水下對(duì)空中聲源的運(yùn)動(dòng)參數(shù)獲取方法。這些算法不僅提高了計(jì)算精度，也縮短了計(jì)算時(shí)間?，F(xiàn)有的基于多普勒頻移獲取運(yùn)動(dòng)參數(shù)大都基于瞬時(shí)頻率。高德洋等[12]提出了一種不需要提取瞬時(shí)頻率的線譜目標(biāo)運(yùn)動(dòng)參數(shù)獲取方法，他首先提出Doppler-warping變換，將warping 變換的思想運(yùn)用到多普勒頻移的去除上，又利用頻譜求熵的方法獲取目標(biāo)的運(yùn)動(dòng)參數(shù)，相較于Ferguson 等的算法，不需要提取瞬時(shí)頻率，因此適用性更高，計(jì)算精度也更高。

分形維數(shù)與熵一樣，都是表征系統(tǒng)復(fù)雜程度的物理量。因此，分形維數(shù)在信號(hào)處理方面有廣泛的應(yīng)用。現(xiàn)階段，分形維數(shù)被廣泛應(yīng)用于圖像分類[13-14]、故障診斷[15]、地形分析[16]、頻譜感知[17-18]等方面。截至目前，鮮有學(xué)者將分形維數(shù)應(yīng)用于線譜信號(hào)目標(biāo)運(yùn)動(dòng)參數(shù)估計(jì)中。本文將頻譜分形維數(shù)應(yīng)用到目標(biāo)運(yùn)動(dòng)速度獲取中。首先利用Doppler-warping變換對(duì)接收信號(hào)重采樣，然后通過(guò)計(jì)算頻譜的分形維數(shù)作為衡量多普勒頻移消除程度的代價(jià)函數(shù)，最后通過(guò)獲取代價(jià)函數(shù)最小值對(duì)應(yīng)的速度值獲取目標(biāo)的運(yùn)動(dòng)速度。

1 理論

1.1 多普勒頻移及Doppler-warping變換

當(dāng)聲源沿不經(jīng)過(guò)接收器的路徑做勻速直線運(yùn)動(dòng)時(shí)，接收器接收信號(hào)存在多普勒頻移，其運(yùn)動(dòng)示意圖如圖1 所示。其中v為目標(biāo)聲源的運(yùn)動(dòng)速度，t0為到達(dá)最近點(diǎn)時(shí)間，r0為最近距離，接收信號(hào)瞬時(shí)頻率由Quinn[3]推導(dǎo)的表達(dá)式得出，如式(1) 所示：

圖1 聲源與接收器運(yùn)動(dòng)示意圖Fig.1 Geometrical model of the wayside acoustic testing model

其中，f0為目標(biāo)聲源發(fā)射信號(hào)的頻率。本文中，設(shè)置t0、r0已知，搜索目標(biāo)勻速直線運(yùn)動(dòng)速度v。由式(1)可以看出，由于多普勒效應(yīng)的存在，接收信號(hào)的頻率不再為單頻，而是隨時(shí)間變化的。

Doppler-warping 變換由中國(guó)海洋大學(xué)水聲研究團(tuán)隊(duì)[12]提出，warping 算子的表達(dá)式如式(2)所示：

使用式(2)對(duì)接收信號(hào)進(jìn)行重采樣，即將式(2)代入式(1)的t中，經(jīng)整理后f(t)=f0，具體推導(dǎo)過(guò)程見(jiàn)文獻(xiàn)[12]。此時(shí)多普勒頻移被消除，原始信號(hào)的單頻特性恢復(fù)。接收信號(hào)經(jīng)過(guò)Doppler-warping 變換前后，LOFAR圖如圖2所示。

圖2 接收信號(hào)Doppler-warping 變換前后Fig.2 LOFAR figure of single-frequency signal and LOFAR figure after Doppler-warping transform

1.2 Sevcik維數(shù)

分形維數(shù)是衡量事物分形特性的一種度量參數(shù)，也是衡量事物復(fù)雜程度的物理量。分形維數(shù)主要可分為Katz 分形維數(shù)、Hausdorff維數(shù)、計(jì)盒維數(shù)、Sevcik 分形維數(shù)和自相似分形維數(shù)等幾種。其中，Sevcik 維數(shù)多用于衡量信號(hào)的波形的分形特性，且計(jì)算過(guò)程較簡(jiǎn)單，計(jì)算時(shí)間短。因此，本文采用Sevcik 維數(shù)來(lái)衡量信號(hào)頻譜的分形特性。Sevcik 維數(shù)的表達(dá)式如式(3)所示：

其中，N為信號(hào)點(diǎn)個(gè)數(shù)，L為信號(hào)波形長(zhǎng)度，其表達(dá)式如式(4)所示：

對(duì)于一個(gè)時(shí)域信號(hào)s(t)，若其時(shí)長(zhǎng)為T(mén)，采樣率為fs，信號(hào)點(diǎn)個(gè)數(shù)為N，則其經(jīng)過(guò)快速傅里葉變換(Fast Fourier transform,FFT)后，頻譜S(f)的xk對(duì)應(yīng)信號(hào)頻率點(diǎn)fk。此時(shí)其信號(hào)波形長(zhǎng)度L的表達(dá)式經(jīng)整理后如式(6)所示：

1.3 算法步驟

由于多普勒的存在，線譜信號(hào)的頻譜能量不再集中在原始頻率上，而是在中心頻率周?chē)幸欢ǖ恼箤?，此時(shí)，信號(hào)頻譜的分形維數(shù)會(huì)比不存在多普勒時(shí)明顯升高，通過(guò)Doppler-warping 變換后，可將多普勒頻移去除，此時(shí)信號(hào)恢復(fù)單頻特性，其頻譜的分形維數(shù)會(huì)明顯降低。Doppler-warping 算子中包含最近點(diǎn)時(shí)間、最近距離、速度等運(yùn)動(dòng)參數(shù)，因此，Doppler-warping算子中運(yùn)動(dòng)參數(shù)的準(zhǔn)確性影響多普勒的去除程度?？赏ㄟ^(guò)線譜和連續(xù)譜結(jié)合的方式，估計(jì)目標(biāo)到達(dá)最近點(diǎn)的時(shí)間和最近距離，本文假設(shè)以上兩參數(shù)已知，利用Doppler-warping 變換結(jié)合頻譜Sevcik 維數(shù)單獨(dú)進(jìn)行速度估計(jì)，其步驟可概括為以下幾步：

(1)設(shè)置速度搜索網(wǎng)格[v1,v2,···,vk]。應(yīng)用時(shí)，可根據(jù)需要設(shè)置速度搜索步長(zhǎng)，在常見(jiàn)目標(biāo)速度估計(jì)應(yīng)用中可將搜索步長(zhǎng)設(shè)為0.1 m/s。

(2) 最近點(diǎn)時(shí)間和最近距離兩參數(shù)由其他方法估計(jì)，設(shè)為已知，計(jì)算每個(gè)速度值下的Dopplerwarping算子。

(3) 將不同算子代入接收時(shí)域信號(hào)s(t)，對(duì)接收信號(hào)進(jìn)行重采樣，獲得信號(hào)s′(t)。

(4) 對(duì)s′(t)進(jìn)行傅里葉變換后獲得頻譜S′(f)，取模并計(jì)算每一個(gè)S′(f)的Sevcik分形維數(shù)D(vk)。

(5) 找到D(vk)取值最小時(shí)對(duì)應(yīng)的速度v，即為速度估計(jì)值。

算法流程圖如圖3所示。

圖3 算法流程Fig.3 Flow chart of the algorithm

2 仿真及算法性能

2.1 仿真及結(jié)果

基于以上理論，進(jìn)行仿真驗(yàn)證，本文算法對(duì)空氣中目標(biāo)和水中目標(biāo)同樣適用，這里以空氣中目標(biāo)為例。聲源信號(hào)為600 Hz 單頻信號(hào)，采樣率為10000 Hz，信號(hào)時(shí)長(zhǎng)10 s。聲源運(yùn)動(dòng)速度為10 m/s，最近距離為10 m，最近點(diǎn)時(shí)間為5 s。聲速為340 m/s。接收信號(hào)及處理結(jié)果如圖4所示。

圖4 仿真結(jié)果Fig.4 Simulation results

設(shè)置速度搜索網(wǎng)格為[1 : 0.2 : 20]，搜索速度v，結(jié)果如圖4(d)所示。

利用本文方法進(jìn)行搜索速度，從圖4(d)可以看出，在真實(shí)速度值10 m/s 處，Sevcik 維數(shù)明顯低于錯(cuò)誤值處，因此能準(zhǔn)確搜索出速度結(jié)果，證明本文提出的方法有效。

2.2 算法性能

為了體現(xiàn)本文提出算法的性能，與用熵作為代價(jià)函數(shù)的搜索結(jié)果進(jìn)行比較。在相同的仿真條件下，分別用熵作為代價(jià)函數(shù)和Sevcik 維數(shù)作為代價(jià)函數(shù)，進(jìn)行速度估計(jì)。估計(jì)結(jié)果如圖5所示。

圖5 熵和Sevcik 搜索結(jié)果對(duì)比Fig.5 Comparision of estimation result of entropy and Sevcik

由于本文算法是尋找代價(jià)函數(shù)最小值，因此定義最小代價(jià)函數(shù)值為0 dB，定義3 dB 處對(duì)應(yīng)的兩速度值的差為主瓣寬度。根據(jù)以上定義，用Sevcik維數(shù)作為代價(jià)函數(shù)時(shí)，主瓣寬度為0.76 m/s，而用熵作為代價(jià)函數(shù)，主瓣寬度為3.81 m/s?？梢?jiàn)用Sevcik維數(shù)和熵作為代價(jià)函數(shù)時(shí)，雖然都可以搜索出正確速度值，但是本文提出的算法，搜索結(jié)果在正確值位置處起伏更明顯，瓣更窄，性能更好。

在實(shí)際應(yīng)用中，往往接收信號(hào)存在噪聲，有時(shí)甚至信噪比(Signal to noise ratio,SNR) 很低，因此，本文利用仿真分析不同SNR下兩種算法的性能對(duì)比。仿真條件為：信號(hào)時(shí)長(zhǎng)20 s，聲源運(yùn)動(dòng)速度為10 m/s，最近距離為10 m，最近點(diǎn)時(shí)間為10 s。分別在仿真的接收信號(hào)中加入高斯白噪聲，SNR計(jì)算方法采用譜SNR。兩種算法搜索結(jié)果對(duì)比如圖6 和表1所示。

表1 典型SNR 下兩方法速度估計(jì)誤差對(duì)比Table 1 Comparison of speed estimation errors between the two methods under typical SNR

圖6 SNR 對(duì)兩種算法搜索結(jié)果的影響Fig.6 Influence of SNR on the two algorithms

從圖6 和表1 中可以看出，在仿真條件相同的情況下，SNR 過(guò)低時(shí)，Sevcik 維數(shù)和熵作為代價(jià)函數(shù)效果均不佳，這是因?yàn)楫?dāng)SNR 較低時(shí)，噪聲對(duì)信號(hào)頻譜有較大影響，因此兩種方法均會(huì)失效。但是，從圖6 可以看出，Sevcik 維數(shù)作為代價(jià)函數(shù)抗噪性提高5 dB 左右。主要原因是兩種算法都基于頻譜的混亂程度，即使在SNR 低于-5 dB 的情況下，多普勒的存在與否也會(huì)對(duì)頻譜的Sevcik 維數(shù)產(chǎn)生較大影響，即存在多普勒時(shí)頻譜的Sevcik 維數(shù)明顯大于不存在時(shí)，因此本文算法可行。而只有在SNR大于-5 dB 以后，用熵值作為代價(jià)函數(shù)時(shí)，多普勒的存在與否對(duì)頻譜的影響較大，熵值作為代價(jià)函數(shù)的算法才會(huì)有效。SNR 對(duì)多普勒是否存在時(shí)兩種不同代價(jià)函數(shù)的影響如圖7(a)、圖7(b)所示。在相同SNR下，多普勒的存在對(duì)頻譜Sevcik維數(shù)的影響明顯大于頻譜熵(如圖7(c)所示)。

圖7 SNR 對(duì)多普勒是否存在時(shí)兩種代價(jià)函數(shù)的影響Fig.7 The effect of the SNR on the presence or absence of Doppler on the two cost functions

3 實(shí)驗(yàn)及結(jié)果

3.1 實(shí)驗(yàn)

針對(duì)上述理論和仿真結(jié)果，進(jìn)行實(shí)驗(yàn)。陸面目標(biāo)采用某卡車(chē)運(yùn)動(dòng)過(guò)程中輻射的線譜噪聲，數(shù)據(jù)來(lái)源為acousticstoday.org。該信號(hào)時(shí)長(zhǎng)30 s，采樣頻率為12000 Hz。聲速為347 m/s。最近距離為35 m，速度為20 km/h。其到達(dá)最近點(diǎn)時(shí)間為15.2 s，本文在處理過(guò)程中將其設(shè)為已知量。接收信號(hào)LOFAR圖以及搜索結(jié)果如圖8所示。

圖8 卡車(chē)實(shí)驗(yàn)結(jié)果Fig.8 Truck experiment results

圖8(a)為接收器接收信號(hào)LOFAR 圖，該圖中線譜多普勒十分明顯。設(shè)置速度網(wǎng)格為[1:0.1:15]，單獨(dú)搜索速度的結(jié)果如圖8(b)所示，結(jié)果為5.8 m/s，相對(duì)誤差僅為4.4%，該結(jié)果證明了本文算法的有效性。本文方法與熵值方法搜索結(jié)果對(duì)比如圖8(c)所示。計(jì)算兩種算法速度搜索結(jié)果的主瓣寬度，熵值算法的主瓣寬度為3.15 m/s，本文算法的主瓣寬度僅為0.74 m/s，結(jié)果表明本文算法主瓣寬度明顯更窄。

為了進(jìn)一步證明本文算法的有效性，選取2022年6 月在青島近海進(jìn)行的海試實(shí)驗(yàn)數(shù)據(jù)，實(shí)驗(yàn)聲源為航道上經(jīng)過(guò)的貨輪，該目標(biāo)為非合作目標(biāo)，貨輪的運(yùn)動(dòng)軌跡經(jīng)接收船船載AIS(船舶自動(dòng)識(shí)別系統(tǒng))記錄如圖9(a)所示，距離變換如圖9(b)所示。接收信號(hào)為貨輪在運(yùn)動(dòng)過(guò)程中輻射線譜噪聲。實(shí)驗(yàn)地聲速為1510 m/s，選取的信號(hào)時(shí)長(zhǎng)為1200 s，信號(hào)采樣頻率為80000 Hz。由于信號(hào)采樣率過(guò)高，而處理信號(hào)所選線譜頻段僅為980～1000 Hz 左右，因此對(duì)接收信號(hào)進(jìn)行降采樣以減少計(jì)算量。接收信號(hào)頻譜圖和LOFAR 圖如圖9(c)、圖9(d) 所示。貨輪運(yùn)動(dòng)平均速度為5.3 m/s，到達(dá)水聽(tīng)器最近距離為2010 m，選取的信號(hào)中，目標(biāo)到達(dá)最近點(diǎn)的時(shí)間為49 s。利用本文提出方法進(jìn)行運(yùn)動(dòng)速度搜索，搜索結(jié)果如圖9(e)所示。設(shè)置速度網(wǎng)格為[1:0.1:10]，速度搜索結(jié)果為5.1 m/s，相對(duì)誤差僅為3.8%。本文算法和熵值算法速度搜索對(duì)比結(jié)果如圖9(f)所示。計(jì)算兩種算法的主瓣寬度，熵值算法的主瓣寬度為1.31 m/s，而本文算法的主瓣寬度僅為0.28 m/s，因此本文算法的性能較熵值算法有明顯提高。

圖9 海試實(shí)驗(yàn)結(jié)果Fig.9 Marine experiment results

以上兩組實(shí)驗(yàn)結(jié)果誤差分別為4.4%和3.8%，均能較為準(zhǔn)確估計(jì)目標(biāo)運(yùn)動(dòng)速度值，速度搜索相對(duì)誤差均小于5%。兩組實(shí)驗(yàn)結(jié)果均證明了本文提出算法的有效性。

4 結(jié)論

本文提出了Doppler-warping 變換結(jié)合頻譜分形維數(shù)的方法來(lái)獲取目標(biāo)的運(yùn)動(dòng)速度。傳統(tǒng)方法通過(guò)迭代算法提取接收信號(hào)的瞬時(shí)頻率，計(jì)算時(shí)間較長(zhǎng)，本文提出的算法過(guò)程簡(jiǎn)單，計(jì)算時(shí)間較短。相對(duì)于用熵作為衡量頻譜復(fù)雜度的代價(jià)函數(shù)，本文利用Sevcik維數(shù)作為代價(jià)函數(shù)，經(jīng)仿真檢驗(yàn)，抗噪能力較熵值方法提高5 dB。經(jīng)仿真和實(shí)驗(yàn)檢驗(yàn)，本文算法代價(jià)函數(shù)主瓣寬度更窄，具有更好的性能。空氣和水中聲目標(biāo)實(shí)驗(yàn)均可準(zhǔn)確搜索出目標(biāo)運(yùn)動(dòng)速度，兩實(shí)驗(yàn)誤差均小于5%，證明了本文所提算法的有效性。下一步擬將本文方法用于寬帶噪聲目標(biāo)和雙目標(biāo)的運(yùn)動(dòng)參數(shù)獲取中，實(shí)現(xiàn)算法應(yīng)用的擴(kuò)展。