一種水聲瞬態(tài)信號(hào)的下降沿檢測(cè)方法

游 波，蔡志明

（1.武昌工學(xué)院電信與物聯(lián)網(wǎng)工程系，湖北武漢 430065；2.海軍工程大學(xué)電子工程學(xué)院水聲工程系，湖北武漢 430033）

1 引言

水中聲納目標(biāo)進(jìn)行戰(zhàn)術(shù)動(dòng)作時(shí)會(huì)輻射瞬態(tài)信號(hào).信號(hào)出現(xiàn)與消失時(shí)被檢測(cè)序列概率密度函數(shù)（Probability Density Function，PDF）的變化是其統(tǒng)計(jì)特性的基本變化，累積和檢測(cè)方法能夠利用該變化進(jìn)行有效的信號(hào)檢測(cè).在主動(dòng)聲納檢測(cè)方面，累積和檢測(cè)后置處理利用匹配濾波的脈壓效應(yīng)將檢測(cè)問(wèn)題轉(zhuǎn)化為瞬態(tài)信號(hào)檢測(cè)問(wèn)題，依據(jù)脈沖出現(xiàn)前后概率密度特性的變化來(lái)實(shí)現(xiàn)檢測(cè)［1，2］.在被動(dòng)檢測(cè)方面，由于水聲瞬態(tài)信號(hào)未知結(jié)構(gòu)、未知長(zhǎng)度、未知到達(dá)時(shí)間［3］和未知PDF模型等問(wèn)題，檢測(cè)難度大于主動(dòng)聲納.變門(mén)限累積和方法［4］解決了未知長(zhǎng)度瞬態(tài)信號(hào)的上升沿檢測(cè)穩(wěn)健性問(wèn)題，用非高斯模型t分布假設(shè)［5］作為描述瞬態(tài)信號(hào)的PDF形式也進(jìn)一步提高了上升沿的檢測(cè)性能.實(shí)際上，下降沿出現(xiàn)時(shí)被檢測(cè)序列PDF的變化與上升沿相反，檢驗(yàn)統(tǒng)計(jì)量、門(mén)限等檢測(cè)參數(shù)需要重新推導(dǎo)計(jì)算.原有累積和序貫檢測(cè)理論假設(shè)將0作為下降沿的判別門(mén)限，會(huì)帶來(lái)較大的檢測(cè)誤差，工程上采用的基于檢驗(yàn)統(tǒng)計(jì)量均值和方差連續(xù)變化趨勢(shì)的下降沿判斷方法，需要反復(fù)調(diào)整經(jīng)驗(yàn)系數(shù)，精確度并不高.

低信噪比情況下瞬態(tài)信號(hào)檢測(cè)常采用基于變換域思想的方法，如基于高階譜的Power-law檢測(cè)器［6，7］、倒三譜［8］、短時(shí)傅里葉變換、Winger-Vill分布、小波變換、經(jīng)驗(yàn)?zāi)B(tài)分解（Empirical Mode Decomposition，EMD）等.常規(guī)傅氏變換不能兼顧時(shí)間和頻率分辨率的要求，小波變換基函數(shù)的選擇對(duì)信號(hào)檢測(cè)分析結(jié)果有較大影響，基于EMD分解方法［9，10］能將本地干擾和信號(hào)分解到不同固有模態(tài)函數(shù)中，去除干擾后的重構(gòu)信號(hào)能提高檢測(cè)性能.有規(guī)瞬態(tài)信號(hào)的EMD分解是有效的，但對(duì)于比如魚(yú)雷出管這樣由寬頻噪聲填充的無(wú)規(guī)瞬態(tài)信號(hào)，EMD分解反而會(huì)增加模態(tài)混疊程度，效果并不理想.累積和檢測(cè)方法另辟蹊徑，依據(jù)瞬態(tài)信號(hào)出現(xiàn)與消失時(shí)PDF的變化進(jìn)行檢測(cè)，因而建立PDF模型是關(guān)鍵環(huán)節(jié).區(qū)分上升沿和下降沿被檢測(cè)序列PDF變化，重新構(gòu)建下降沿的累積和檢測(cè)過(guò)程，提高下降沿檢測(cè)精度，對(duì)于獲取瞬態(tài)信號(hào)清晰的時(shí)頻結(jié)構(gòu)、進(jìn)行有效的瞬態(tài)信號(hào)分類(lèi)識(shí)別和作出及時(shí)的戰(zhàn)術(shù)判斷有著十分重要意義.

2 檢驗(yàn)統(tǒng)計(jì)量推導(dǎo)

考慮到水聲瞬態(tài)信號(hào)的突發(fā)性和不穩(wěn)定性、聲信道傳播特性未知等因素，在一般意義上將聲納接收機(jī)輸入信號(hào)中的瞬態(tài)變化視作方差發(fā)生變化的高斯分布［11］.給定獨(dú)立觀測(cè)值序列X={x1，x2，…，xn}，n為被檢測(cè)序列長(zhǎng)度.假設(shè)瞬態(tài)信號(hào)結(jié)束時(shí)，概率分布律從f0變化至f1，變化發(fā)生點(diǎn)θ未知，即：

式中f0(xi)=N(0，σ2)表示在H0假設(shè)（信號(hào)未結(jié)束）下被檢測(cè)序列服從均值為零、方差σ2的分布；通過(guò)對(duì)背景零均值化和功率歸一化處理，可將H1假設(shè)下（信號(hào)結(jié)束）背景噪聲的PDF分布設(shè)為f1(xi)=N(0，1).備擇假設(shè)｛H1：θ=ν(ν＜n)｝對(duì)原假設(shè)｛H0：θ=∞（無(wú)變點(diǎn)）｝的對(duì)數(shù)似然比統(tǒng)計(jì)量［12］為：

式（2）中p1(X)表示H1假設(shè)下檢驗(yàn)統(tǒng)計(jì)量X的多點(diǎn)聯(lián)合PDF，其中變點(diǎn)發(fā)生在i=ν的位置.p0(X)表示H0假設(shè)下X的PDF分布.令b=σ2·lnσ2/(σ2-1)為偏差，則對(duì)數(shù)似然比可寫(xiě)成：

由于σ2為變化發(fā)生前瞬態(tài)信號(hào)加噪聲的總方差，因此σ2＞1，(1-σ2)/2σ2＜0.依據(jù)最大似然比準(zhǔn)則，對(duì)數(shù)似然比統(tǒng)計(jì)量等價(jià)于下述統(tǒng)計(jì)量：

由于

則

所以，Zn可遞推實(shí)現(xiàn)，即得到下降沿檢測(cè)完整的檢驗(yàn)統(tǒng)計(jì)量：

定義檢測(cè)器非線性量或更新量為yn=g(xn)=x2nb=x2n-σ2lnσ2/(σ2-1)，是被檢測(cè)序列X的函數(shù).根據(jù)雅可比行列式，可得更新量的概率分布密度PDF：

圖1 原假設(shè)下E[g]隨方差σ2的變化關(guān)系

圖2 備擇假設(shè)下E[g]隨方差σ2的變化關(guān)系

常規(guī)累積和上升沿檢測(cè)理論模型可描述為序貫檢測(cè)［4，5］：

式（13）中0和h分別設(shè)為上、下門(mén)限，顯然將下降沿的檢測(cè)門(mén)限簡(jiǎn)單設(shè)為0是不合適的，應(yīng)根據(jù)本節(jié)的結(jié)論重新推導(dǎo)和計(jì)算.完整的累積和瞬態(tài)檢測(cè)過(guò)程應(yīng)包括上升沿和下降沿的同時(shí)檢測(cè)，是兩者的結(jié)合.

3 檢測(cè)參數(shù)計(jì)算

3.1 下降沿檢測(cè)性能分析模型

下面分析在序貫檢測(cè)三種狀態(tài)下（如式12）第n時(shí)刻檢驗(yàn)統(tǒng)計(jì)量Zn的概率，再進(jìn)行N步狀態(tài)轉(zhuǎn)移和遞推，得到檢測(cè)概率和虛警概率的數(shù)學(xué)表達(dá).該分析模型為后續(xù)基于紐曼-皮爾遜準(zhǔn)則，在給定檢測(cè)和虛警概率條件下推算門(mén)限提供理論基礎(chǔ).該過(guò)程分析步驟與上升沿性能分析模型類(lèi)似，下面重點(diǎn)給出下降沿分析模型的不同點(diǎn).

關(guān)注具體的累積過(guò)程：Zn=Zn-1+g(xn)，n時(shí)刻檢驗(yàn)統(tǒng)計(jì)量Zn的概率分布律等于n-1時(shí)刻Zn在范圍(h，0)間的概率分布律fn-1(z)與更新量g(xn)概率分布律fg(z)的卷積［13，14］：

比較式（12）和式（13）可見(jiàn)，在推算三種狀態(tài)下檢驗(yàn)統(tǒng)計(jì)量Zn的概率時(shí)，下降沿檢驗(yàn)統(tǒng)計(jì)量Zn的范圍與上升沿顯著不同.定義π0(n)表示未檢測(cè)到下降沿的概率：

定義π1(n)表示在上述條件下在時(shí)刻n結(jié)束且檢測(cè)結(jié)果為檢測(cè)到下降沿的概率：

定義π2(n)表示檢測(cè)過(guò)程能持續(xù)到下一個(gè)時(shí)刻的概率，即檢驗(yàn)統(tǒng)計(jì)量Zn在范圍(h，0)間的概率：

對(duì)結(jié)果進(jìn)行歸一化處理，一次檢測(cè)過(guò)程持續(xù)到n時(shí)刻檢驗(yàn)統(tǒng)計(jì)量Zn在范圍(h，0)間的概率分布密度f(wàn)n(z)表為：

接下來(lái)的思路與上升沿檢測(cè)類(lèi)似［13，14］，由式（14）～（18）對(duì)累積過(guò)程進(jìn)行狀態(tài)遞推，可得在n時(shí)刻之前未過(guò)門(mén)限或置零復(fù)位的假設(shè)下，檢測(cè)持續(xù)到時(shí)刻n+1的概率pon(n)=pon(n-1)π2(n)、持續(xù)到n時(shí)刻并判為未檢測(cè)到下降沿的概率p0(n)=pon(n-1)π0(n)、在n時(shí)刻結(jié)束且檢測(cè)到下降沿的概率p1(n)=pon(n-1)π1(n).

下降沿檢測(cè)包括正常檢測(cè)和延遲檢測(cè)兩部分之和［13，14］，但與上升沿檢測(cè)不同的是，正常檢測(cè)是指初始檢測(cè)狀態(tài)為檢驗(yàn)統(tǒng)計(jì)量進(jìn)入穩(wěn)態(tài)（即瞬態(tài)信號(hào)已經(jīng)發(fā)生，H0假設(shè)），在H1假設(shè)下低于檢測(cè)門(mén)限，此時(shí)不必考慮0初始狀態(tài)的情形.延遲檢測(cè)是指H0假設(shè)下進(jìn)入穩(wěn)態(tài)，在H1假設(shè)下未低于門(mén)限，但在重新轉(zhuǎn)換成H0假設(shè)時(shí)低于門(mén)限的情形.檢測(cè)概率可表為：

式中p1(n)為H0假設(shè)下持續(xù)到n時(shí)刻，檢測(cè)到下降沿的概率.正常檢測(cè)概率Pdnorm包括初始狀態(tài)為穩(wěn)態(tài)時(shí)未經(jīng)置零低于門(mén)限和經(jīng)過(guò)k次置零低于門(mén)限的概率之和［14］：

式中ps1s(n)表示H0假設(shè)下穩(wěn)定初始狀態(tài)為fss(z)時(shí)檢驗(yàn)過(guò)程一直持續(xù)到n時(shí)刻檢測(cè)判決為下降沿（記為1）的概率，N為被檢測(cè)序列的長(zhǎng)度.p(0k)(m)為k次置零復(fù)位（未檢測(cè)到下降沿，記為0）的概率.p01(n)表示剩余(N-m)個(gè)序列點(diǎn)從0狀態(tài)開(kāi)始，檢測(cè)到下降沿，記為1的概率.H0假設(shè)下進(jìn)入穩(wěn)態(tài)時(shí)檢驗(yàn)統(tǒng)計(jì)量Zn的概率密度函數(shù)［14］表為：

對(duì)累積過(guò)程進(jìn)行狀態(tài)遞推的基礎(chǔ)是對(duì)更新量g(xn)進(jìn)行量化，根據(jù)式（10）、式（11）更新量的PDF構(gòu)建馬爾科夫狀態(tài)轉(zhuǎn)移矩陣，這是進(jìn)行虛警概率和檢測(cè)概率推算的必要條件，具體過(guò)程詳見(jiàn)之前的工作［14］.虛警概率的計(jì)算考慮在H0假設(shè)下低于檢測(cè)門(mén)限的概率，詳細(xì)過(guò)程不再贅述.

3.2 檢測(cè)參數(shù)的計(jì)算結(jié)果

前兩節(jié)性能分析模型給出了基于狀態(tài)遞推的檢測(cè)和虛警概率的計(jì)算過(guò)程，反過(guò)來(lái)，根據(jù)紐曼皮爾遜準(zhǔn)則，在方差變化高斯分布的PDF假設(shè)，如式（1），在滿(mǎn)足虛警概率Pfa=10-4和檢測(cè)概率Pd=0.6條件下，可以反推累積和方法的門(mén)限h和方差σ2等檢測(cè)參數(shù).根據(jù)其隨瞬態(tài)變化長(zhǎng)度N的變化規(guī)律，計(jì)算搜索過(guò)程可以進(jìn)行優(yōu)化，結(jié)果如圖3和圖4所示.

圖3 門(mén)限隨著瞬態(tài)變化長(zhǎng)度的結(jié)果圖

圖4 方差隨著瞬態(tài)變化長(zhǎng)度的結(jié)果圖

4 數(shù)據(jù)驗(yàn)證

在工程實(shí)現(xiàn)上，常規(guī)累積和檢驗(yàn)方法的下降沿檢測(cè)是從當(dāng)前時(shí)刻開(kāi)始，將檢驗(yàn)統(tǒng)計(jì)量以L為長(zhǎng)度分成M段，分別計(jì)算每段數(shù)據(jù)的方差V(xi)和均值J(xi)，i=1～M.如果出現(xiàn)M段V(xi)的增大或者J(xi)的連續(xù)下降，即判為信號(hào)下降沿.M和L與采樣率、信噪比和瞬態(tài)信號(hào)特征有關(guān)，均為經(jīng)驗(yàn)值，檢測(cè)穩(wěn)健性較差.本文三個(gè)例子信號(hào)采樣率fs=1000 Hz，M=50，L=5.下面利用常規(guī)方法和改進(jìn)方法對(duì)艦船轉(zhuǎn)舵仿真信號(hào)、魚(yú)雷出管仿真信號(hào)和水池重物落水實(shí)測(cè)信號(hào)進(jìn)行下降沿檢測(cè)分析比較.

4.1 艦船轉(zhuǎn)舵信號(hào)仿真和檢測(cè)結(jié)果

艦船轉(zhuǎn)舵輻射噪聲主要由機(jī)械撞擊和空化氣泡破滅噪聲組成.仿真信號(hào)模型采用多個(gè)起始時(shí)間不同的指數(shù)衰減正弦信號(hào)之和，同組合高斯包絡(luò)信號(hào)的疊加來(lái)表示，表達(dá)為［15］：

式中G(t)為組合高斯包絡(luò)，設(shè)N=5，正弦信號(hào)幅度A1=A2=A4=0.23，A3=0.25，A5=0.08，指數(shù)衰減系數(shù)λ1=0，λ2=λ4=3.3，λ3=1.96，λ5=5.2，正 弦 信 號(hào) 頻 率f1=f4=150 Hz，f2=f3=120 Hz，f=50 Hz，f5=100 Hz，初相φk=0.仿真波形如圖5所示，時(shí)長(zhǎng)1 s，采樣率為1000 Hz.

圖5 轉(zhuǎn)舵仿真信號(hào)

將信號(hào)以3 dB的信噪比疊加到高斯白噪聲背景1s處，2 s處信號(hào)結(jié)束，如圖6所示.轉(zhuǎn)舵信號(hào)由強(qiáng)變?nèi)?，在信噪比變?nèi)鯐r(shí)常規(guī)累積和下降沿判定誤差較大，判為1.792 s，而新方法判定結(jié)果為2.034 s，絕對(duì)誤差為0.034 s，精度更高.表1為不同信噪比下進(jìn)行1000次蒙特卡洛實(shí)驗(yàn)得到的兩種方法平均絕對(duì)誤差，新方法下降沿檢測(cè)誤差有明顯下降.

圖6 轉(zhuǎn)舵仿真信號(hào)下降沿檢測(cè)效果比對(duì)

表1 不同信噪比下1000次蒙特卡洛實(shí)驗(yàn)兩種方法平均絕對(duì)誤差

4.2 魚(yú)雷出管信號(hào)仿真和檢測(cè)結(jié)果

魚(yú)雷發(fā)射過(guò)程中的瞬態(tài)信號(hào)［16］主要包括三部分，蓋板聲、出管聲和啟動(dòng)點(diǎn)火聲.蓋板聲是在發(fā)射管與海水保持連通后，前蓋板開(kāi)啟過(guò)程中產(chǎn)生的機(jī)械沖擊類(lèi)信號(hào).采用包絡(luò)加權(quán)的多次諧波敲擊瞬態(tài)信號(hào)進(jìn)行模擬，壓縮空氣膨脹將魚(yú)雷推出管的信號(hào)采用指數(shù)衰減包絡(luò)加權(quán)的帶限噪聲進(jìn)行模擬，啟動(dòng)點(diǎn)火信號(hào)為燃料點(diǎn)火燃燒和螺旋槳運(yùn)動(dòng)的噪聲，采用對(duì)數(shù)調(diào)頻和均勻分布的寬帶填充噪聲復(fù)合構(gòu)成，該噪聲相比前兩者信噪比更大，并且熱動(dòng)力雷和電動(dòng)力雷的啟動(dòng)點(diǎn)火噪聲特性不同，因而本文針對(duì)前兩部分相對(duì)比較弱的噪聲進(jìn)行仿真和檢測(cè).

多次諧波敲擊瞬態(tài)信號(hào)的包絡(luò)設(shè)為上升沿陡、下降沿緩的形狀，可表為［16］：

式中A0為幅度，τ1和τ2為信號(hào)上升時(shí)刻和信號(hào)持續(xù)時(shí)間，仿真中τ1=0.2s，τ2=0.6s，a=1，b=0.5.被調(diào)制的是以平臺(tái)自身固有頻率為主、含有多次諧波分量的信號(hào)，可表為：

式中Ai為對(duì)應(yīng)諧波信號(hào)幅度，N為諧波次數(shù).取基頻為f0=100 Hz，諧波幅度分別為1、0.4、0.2、0.1.魚(yú)雷出管仿真結(jié)果如圖7所示，蓋板聲時(shí)長(zhǎng)0.6 s，出管信號(hào)時(shí)長(zhǎng)0.8 s，分別在1 s和3 s處以-3.7 dB和-3.5 dB的信噪比疊加到高斯白噪聲N(0，1)上，真實(shí)下降沿分別為1.6 s和3.8 s.檢測(cè)效果比對(duì)如圖8所示.可見(jiàn)當(dāng)瞬態(tài)信號(hào)信噪比較低時(shí)，常規(guī)下降沿檢測(cè)易受到隨機(jī)噪聲的干擾，出現(xiàn)誤判，出現(xiàn)多次下降沿.而新下降沿檢測(cè)算法的檢測(cè)結(jié)果更加準(zhǔn)確，絕對(duì)誤差為0.028 s和0.029 s.表2為不同信噪比下進(jìn)行1000次蒙特卡洛實(shí)驗(yàn)得到的兩種方法平均絕對(duì)誤差，可見(jiàn)新方法對(duì)蓋板聲和出管聲下降沿檢測(cè)精度更高.

圖7 魚(yú)雷出管信號(hào)仿真

圖8 魚(yú)雷出管信號(hào)下降沿檢測(cè)效果比對(duì)

表2 不同信噪比下1000次蒙特卡洛實(shí)驗(yàn)兩種方法平均絕對(duì)誤差

4.3 水池?cái)?shù)據(jù)驗(yàn)證

將重物小鐵錘在2 m高度自由落入消聲水池，由于重物體積小，擊水聲弱，故采集的信號(hào)主要為脈動(dòng)聲，如圖9所示.信號(hào)時(shí)長(zhǎng)4 s，疏樣后的采樣率1000 Hz，主要頻率范圍為100～210 Hz，帶通濾波頻段［50 Hz，300 Hz］，以-5 dB信噪比添加到帶限高斯白噪聲上，如圖10上圖所示.進(jìn)行零均值化和背景噪聲功率歸一化處理，常規(guī)和改進(jìn)下降沿的檢測(cè)結(jié)果如圖10中圖所示.依據(jù)均值和方差變化的常規(guī)方法在脈動(dòng)信號(hào)信噪比逐漸減弱時(shí)出現(xiàn)較大的下降沿檢測(cè)誤差，精度明顯低于新方法.常規(guī)方法檢測(cè)到的信號(hào)下降沿為2.557 s，改進(jìn)方法檢測(cè)到的信號(hào)下降沿為2.855 s.事實(shí)上信號(hào)的主要能量約分布在2.1 s至2.8 s的時(shí)間段，可見(jiàn)改進(jìn)方法精度更高，穩(wěn)健性更好.圖11為不同信噪比下進(jìn)行1000次蒙特卡洛實(shí)驗(yàn)得到的兩種方法平均絕對(duì)誤差對(duì)比，可見(jiàn)隨著信噪比的降低，相對(duì)于常規(guī)方法，新方法對(duì)重物落水信號(hào)下降沿檢測(cè)誤差更小.

圖9 水池重物落水脈動(dòng)信號(hào)

圖10 水池重物落水信號(hào)下降沿檢測(cè)效果比對(duì)

圖11 不同信噪比下1000次蒙特卡洛實(shí)驗(yàn)兩種方法平均絕對(duì)誤差對(duì)比

5 小結(jié)

在利用累積和算法對(duì)瞬態(tài)信號(hào)上升沿檢測(cè)的基礎(chǔ)上，本文依據(jù)信號(hào)結(jié)束時(shí)PDF的變化重新推導(dǎo)了檢驗(yàn)統(tǒng)計(jì)量和更新量的概率密度函數(shù)表達(dá)；在性能分析模型方面，仍然利用FFT方法，但嚴(yán)格區(qū)分了上升沿檢測(cè)和下降沿檢測(cè)的不同，推導(dǎo)了檢測(cè)概率的計(jì)算公式；依據(jù)紐曼皮爾遜準(zhǔn)則計(jì)算了在虛警概率Pfa=10-4和檢測(cè)概率Pd=0.6條件下，方差變化高斯過(guò)程假設(shè)下的檢測(cè)參數(shù).艦船轉(zhuǎn)舵仿真信號(hào)、魚(yú)雷出管仿真信號(hào)和水池重物落水實(shí)測(cè)脈動(dòng)信號(hào)的驗(yàn)證結(jié)果表明，低信噪比下，新方法對(duì)瞬態(tài)信號(hào)下降沿的判斷精確度更高，穩(wěn)健性更好.