游 波,蔡志明
(1.武昌工學(xué)院電信與物聯(lián)網(wǎng)工程系,湖北武漢 430065;2.海軍工程大學(xué)電子工程學(xué)院水聲工程系,湖北武漢 430033)
水中聲納目標(biāo)進行戰(zhàn)術(shù)動作時會輻射瞬態(tài)信號.信號出現(xiàn)與消失時被檢測序列概率密度函數(shù)(Probability Density Function,PDF)的變化是其統(tǒng)計特性的基本變化,累積和檢測方法能夠利用該變化進行有效的信號檢測.在主動聲納檢測方面,累積和檢測后置處理利用匹配濾波的脈壓效應(yīng)將檢測問題轉(zhuǎn)化為瞬態(tài)信號檢測問題,依據(jù)脈沖出現(xiàn)前后概率密度特性的變化來實現(xiàn)檢測[1,2].在被動檢測方面,由于水聲瞬態(tài)信號未知結(jié)構(gòu)、未知長度、未知到達(dá)時間[3]和未知PDF模型等問題,檢測難度大于主動聲納.變門限累積和方法[4]解決了未知長度瞬態(tài)信號的上升沿檢測穩(wěn)健性問題,用非高斯模型t分布假設(shè)[5]作為描述瞬態(tài)信號的PDF形式也進一步提高了上升沿的檢測性能.實際上,下降沿出現(xiàn)時被檢測序列PDF的變化與上升沿相反,檢驗統(tǒng)計量、門限等檢測參數(shù)需要重新推導(dǎo)計算.原有累積和序貫檢測理論假設(shè)將0作為下降沿的判別門限,會帶來較大的檢測誤差,工程上采用的基于檢驗統(tǒng)計量均值和方差連續(xù)變化趨勢的下降沿判斷方法,需要反復(fù)調(diào)整經(jīng)驗系數(shù),精確度并不高.
低信噪比情況下瞬態(tài)信號檢測常采用基于變換域思想的方法,如基于高階譜的Power-law檢測器[6,7]、倒三譜[8]、短時傅里葉變換、Winger-Vill分布、小波變換、經(jīng)驗?zāi)B(tài)分解(Empirical Mode Decomposition,EMD)等.常規(guī)傅氏變換不能兼顧時間和頻率分辨率的要求,小波變換基函數(shù)的選擇對信號檢測分析結(jié)果有較大影響,基于EMD分解方法[9,10]能將本地干擾和信號分解到不同固有模態(tài)函數(shù)中,去除干擾后的重構(gòu)信號能提高檢測性能.有規(guī)瞬態(tài)信號的EMD分解是有效的,但對于比如魚雷出管這樣由寬頻噪聲填充的無規(guī)瞬態(tài)信號,EMD分解反而會增加模態(tài)混疊程度,效果并不理想.累積和檢測方法另辟蹊徑,依據(jù)瞬態(tài)信號出現(xiàn)與消失時PDF的變化進行檢測,因而建立PDF模型是關(guān)鍵環(huán)節(jié).區(qū)分上升沿和下降沿被檢測序列PDF變化,重新構(gòu)建下降沿的累積和檢測過程,提高下降沿檢測精度,對于獲取瞬態(tài)信號清晰的時頻結(jié)構(gòu)、進行有效的瞬態(tài)信號分類識別和作出及時的戰(zhàn)術(shù)判斷有著十分重要意義.
考慮到水聲瞬態(tài)信號的突發(fā)性和不穩(wěn)定性、聲信道傳播特性未知等因素,在一般意義上將聲納接收機輸入信號中的瞬態(tài)變化視作方差發(fā)生變化的高斯分布[11].給定獨立觀測值序列X={x1,x2,…,xn},n為被檢測序列長度.假設(shè)瞬態(tài)信號結(jié)束時,概率分布律從f0變化至f1,變化發(fā)生點θ未知,即:
式中f0(xi)=N(0,σ2)表示在H0假設(shè)(信號未結(jié)束)下被檢測序列服從均值為零、方差σ2的分布;通過對背景零均值化和功率歸一化處理,可將H1假設(shè)下(信號結(jié)束)背景噪聲的PDF分布設(shè)為f1(xi)=N(0,1).備擇假設(shè){H1:θ=ν(ν<n)}對原假設(shè){H0:θ=∞(無變點)}的對數(shù)似然比統(tǒng)計量[12]為:
式(2)中p1(X)表示H1假設(shè)下檢驗統(tǒng)計量X的多點聯(lián)合PDF,其中變點發(fā)生在i=ν的位置.p0(X)表示H0假設(shè)下X的PDF分布.令b=σ2·lnσ2/(σ2-1)為偏差,則對數(shù)似然比可寫成:
由于σ2為變化發(fā)生前瞬態(tài)信號加噪聲的總方差,因此σ2>1,(1-σ2)/2σ2<0.依據(jù)最大似然比準(zhǔn)則,對數(shù)似然比統(tǒng)計量等價于下述統(tǒng)計量:
由于
則
所以,Zn可遞推實現(xiàn),即得到下降沿檢測完整的檢驗統(tǒng)計量:
定義檢測器非線性量或更新量為yn=g(xn)=x2nb=x2n-σ2lnσ2/(σ2-1),是被檢測序列X的函數(shù).根據(jù)雅可比行列式,可得更新量的概率分布密度PDF:
圖1 原假設(shè)下E[g]隨方差σ2的變化關(guān)系
圖2 備擇假設(shè)下E[g]隨方差σ2的變化關(guān)系
常規(guī)累積和上升沿檢測理論模型可描述為序貫檢測[4,5]:
式(13)中0和h分別設(shè)為上、下門限,顯然將下降沿的檢測門限簡單設(shè)為0是不合適的,應(yīng)根據(jù)本節(jié)的結(jié)論重新推導(dǎo)和計算.完整的累積和瞬態(tài)檢測過程應(yīng)包括上升沿和下降沿的同時檢測,是兩者的結(jié)合.
下面分析在序貫檢測三種狀態(tài)下(如式12)第n時刻檢驗統(tǒng)計量Zn的概率,再進行N步狀態(tài)轉(zhuǎn)移和遞推,得到檢測概率和虛警概率的數(shù)學(xué)表達(dá).該分析模型為后續(xù)基于紐曼-皮爾遜準(zhǔn)則,在給定檢測和虛警概率條件下推算門限提供理論基礎(chǔ).該過程分析步驟與上升沿性能分析模型類似,下面重點給出下降沿分析模型的不同點.
關(guān)注具體的累積過程:Zn=Zn-1+g(xn),n時刻檢驗統(tǒng)計量Zn的概率分布律等于n-1時刻Zn在范圍(h,0)間的概率分布律fn-1(z)與更新量g(xn)概率分布律fg(z)的卷積[13,14]:
比較式(12)和式(13)可見,在推算三種狀態(tài)下檢驗統(tǒng)計量Zn的概率時,下降沿檢驗統(tǒng)計量Zn的范圍與上升沿顯著不同.定義π0(n)表示未檢測到下降沿的概率:
定義π1(n)表示在上述條件下在時刻n結(jié)束且檢測結(jié)果為檢測到下降沿的概率:
定義π2(n)表示檢測過程能持續(xù)到下一個時刻的概率,即檢驗統(tǒng)計量Zn在范圍(h,0)間的概率:
對結(jié)果進行歸一化處理,一次檢測過程持續(xù)到n時刻檢驗統(tǒng)計量Zn在范圍(h,0)間的概率分布密度fn(z)表為:
接下來的思路與上升沿檢測類似[13,14],由式(14)~(18)對累積過程進行狀態(tài)遞推,可得在n時刻之前未過門限或置零復(fù)位的假設(shè)下,檢測持續(xù)到時刻n+1的概率pon(n)=pon(n-1)π2(n)、持續(xù)到n時刻并判為未檢測到下降沿的概率p0(n)=pon(n-1)π0(n)、在n時刻結(jié)束且檢測到下降沿的概率p1(n)=pon(n-1)π1(n).
下降沿檢測包括正常檢測和延遲檢測兩部分之和[13,14],但與上升沿檢測不同的是,正常檢測是指初始檢測狀態(tài)為檢驗統(tǒng)計量進入穩(wěn)態(tài)(即瞬態(tài)信號已經(jīng)發(fā)生,H0假設(shè)),在H1假設(shè)下低于檢測門限,此時不必考慮0初始狀態(tài)的情形.延遲檢測是指H0假設(shè)下進入穩(wěn)態(tài),在H1假設(shè)下未低于門限,但在重新轉(zhuǎn)換成H0假設(shè)時低于門限的情形.檢測概率可表為:
式中p1(n)為H0假設(shè)下持續(xù)到n時刻,檢測到下降沿的概率.正常檢測概率Pdnorm包括初始狀態(tài)為穩(wěn)態(tài)時未經(jīng)置零低于門限和經(jīng)過k次置零低于門限的概率之和[14]:
式中ps1s(n)表示H0假設(shè)下穩(wěn)定初始狀態(tài)為fss(z)時檢驗過程一直持續(xù)到n時刻檢測判決為下降沿(記為1)的概率,N為被檢測序列的長度.p(0k)(m)為k次置零復(fù)位(未檢測到下降沿,記為0)的概率.p01(n)表示剩余(N-m)個序列點從0狀態(tài)開始,檢測到下降沿,記為1的概率.H0假設(shè)下進入穩(wěn)態(tài)時檢驗統(tǒng)計量Zn的概率密度函數(shù)[14]表為:
對累積過程進行狀態(tài)遞推的基礎(chǔ)是對更新量g(xn)進行量化,根據(jù)式(10)、式(11)更新量的PDF構(gòu)建馬爾科夫狀態(tài)轉(zhuǎn)移矩陣,這是進行虛警概率和檢測概率推算的必要條件,具體過程詳見之前的工作[14].虛警概率的計算考慮在H0假設(shè)下低于檢測門限的概率,詳細(xì)過程不再贅述.
前兩節(jié)性能分析模型給出了基于狀態(tài)遞推的檢測和虛警概率的計算過程,反過來,根據(jù)紐曼皮爾遜準(zhǔn)則,在方差變化高斯分布的PDF假設(shè),如式(1),在滿足虛警概率Pfa=10-4和檢測概率Pd=0.6條件下,可以反推累積和方法的門限h和方差σ2等檢測參數(shù).根據(jù)其隨瞬態(tài)變化長度N的變化規(guī)律,計算搜索過程可以進行優(yōu)化,結(jié)果如圖3和圖4所示.
圖3 門限隨著瞬態(tài)變化長度的結(jié)果圖
圖4 方差隨著瞬態(tài)變化長度的結(jié)果圖
在工程實現(xiàn)上,常規(guī)累積和檢驗方法的下降沿檢測是從當(dāng)前時刻開始,將檢驗統(tǒng)計量以L為長度分成M段,分別計算每段數(shù)據(jù)的方差V(xi)和均值J(xi),i=1~M.如果出現(xiàn)M段V(xi)的增大或者J(xi)的連續(xù)下降,即判為信號下降沿.M和L與采樣率、信噪比和瞬態(tài)信號特征有關(guān),均為經(jīng)驗值,檢測穩(wěn)健性較差.本文三個例子信號采樣率fs=1000 Hz,M=50,L=5.下面利用常規(guī)方法和改進方法對艦船轉(zhuǎn)舵仿真信號、魚雷出管仿真信號和水池重物落水實測信號進行下降沿檢測分析比較.
艦船轉(zhuǎn)舵輻射噪聲主要由機械撞擊和空化氣泡破滅噪聲組成.仿真信號模型采用多個起始時間不同的指數(shù)衰減正弦信號之和,同組合高斯包絡(luò)信號的疊加來表示,表達(dá)為[15]:
式中G(t)為組合高斯包絡(luò),設(shè)N=5,正弦信號幅度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,正 弦 信 號 頻 率f1=f4=150 Hz,f2=f3=120 Hz,f=50 Hz,f5=100 Hz,初相φk=0.仿真波形如圖5所示,時長1 s,采樣率為1000 Hz.
圖5 轉(zhuǎn)舵仿真信號
將信號以3 dB的信噪比疊加到高斯白噪聲背景1s處,2 s處信號結(jié)束,如圖6所示.轉(zhuǎn)舵信號由強變?nèi)酰谛旁氡茸內(nèi)鯐r常規(guī)累積和下降沿判定誤差較大,判為1.792 s,而新方法判定結(jié)果為2.034 s,絕對誤差為0.034 s,精度更高.表1為不同信噪比下進行1000次蒙特卡洛實驗得到的兩種方法平均絕對誤差,新方法下降沿檢測誤差有明顯下降.
圖6 轉(zhuǎn)舵仿真信號下降沿檢測效果比對
表1 不同信噪比下1000次蒙特卡洛實驗兩種方法平均絕對誤差
魚雷發(fā)射過程中的瞬態(tài)信號[16]主要包括三部分,蓋板聲、出管聲和啟動點火聲.蓋板聲是在發(fā)射管與海水保持連通后,前蓋板開啟過程中產(chǎn)生的機械沖擊類信號.采用包絡(luò)加權(quán)的多次諧波敲擊瞬態(tài)信號進行模擬,壓縮空氣膨脹將魚雷推出管的信號采用指數(shù)衰減包絡(luò)加權(quán)的帶限噪聲進行模擬,啟動點火信號為燃料點火燃燒和螺旋槳運動的噪聲,采用對數(shù)調(diào)頻和均勻分布的寬帶填充噪聲復(fù)合構(gòu)成,該噪聲相比前兩者信噪比更大,并且熱動力雷和電動力雷的啟動點火噪聲特性不同,因而本文針對前兩部分相對比較弱的噪聲進行仿真和檢測.
多次諧波敲擊瞬態(tài)信號的包絡(luò)設(shè)為上升沿陡、下降沿緩的形狀,可表為[16]:
式中A0為幅度,τ1和τ2為信號上升時刻和信號持續(xù)時間,仿真中τ1=0.2s,τ2=0.6s,a=1,b=0.5.被調(diào)制的是以平臺自身固有頻率為主、含有多次諧波分量的信號,可表為:
式中Ai為對應(yīng)諧波信號幅度,N為諧波次數(shù).取基頻為f0=100 Hz,諧波幅度分別為1、0.4、0.2、0.1.魚雷出管仿真結(jié)果如圖7所示,蓋板聲時長0.6 s,出管信號時長0.8 s,分別在1 s和3 s處以-3.7 dB和-3.5 dB的信噪比疊加到高斯白噪聲N(0,1)上,真實下降沿分別為1.6 s和3.8 s.檢測效果比對如圖8所示.可見當(dāng)瞬態(tài)信號信噪比較低時,常規(guī)下降沿檢測易受到隨機噪聲的干擾,出現(xiàn)誤判,出現(xiàn)多次下降沿.而新下降沿檢測算法的檢測結(jié)果更加準(zhǔn)確,絕對誤差為0.028 s和0.029 s.表2為不同信噪比下進行1000次蒙特卡洛實驗得到的兩種方法平均絕對誤差,可見新方法對蓋板聲和出管聲下降沿檢測精度更高.
圖7 魚雷出管信號仿真
圖8 魚雷出管信號下降沿檢測效果比對
表2 不同信噪比下1000次蒙特卡洛實驗兩種方法平均絕對誤差
將重物小鐵錘在2 m高度自由落入消聲水池,由于重物體積小,擊水聲弱,故采集的信號主要為脈動聲,如圖9所示.信號時長4 s,疏樣后的采樣率1000 Hz,主要頻率范圍為100~210 Hz,帶通濾波頻段[50 Hz,300 Hz],以-5 dB信噪比添加到帶限高斯白噪聲上,如圖10上圖所示.進行零均值化和背景噪聲功率歸一化處理,常規(guī)和改進下降沿的檢測結(jié)果如圖10中圖所示.依據(jù)均值和方差變化的常規(guī)方法在脈動信號信噪比逐漸減弱時出現(xiàn)較大的下降沿檢測誤差,精度明顯低于新方法.常規(guī)方法檢測到的信號下降沿為2.557 s,改進方法檢測到的信號下降沿為2.855 s.事實上信號的主要能量約分布在2.1 s至2.8 s的時間段,可見改進方法精度更高,穩(wěn)健性更好.圖11為不同信噪比下進行1000次蒙特卡洛實驗得到的兩種方法平均絕對誤差對比,可見隨著信噪比的降低,相對于常規(guī)方法,新方法對重物落水信號下降沿檢測誤差更小.
圖9 水池重物落水脈動信號
圖10 水池重物落水信號下降沿檢測效果比對
圖11 不同信噪比下1000次蒙特卡洛實驗兩種方法平均絕對誤差對比
在利用累積和算法對瞬態(tài)信號上升沿檢測的基礎(chǔ)上,本文依據(jù)信號結(jié)束時PDF的變化重新推導(dǎo)了檢驗統(tǒng)計量和更新量的概率密度函數(shù)表達(dá);在性能分析模型方面,仍然利用FFT方法,但嚴(yán)格區(qū)分了上升沿檢測和下降沿檢測的不同,推導(dǎo)了檢測概率的計算公式;依據(jù)紐曼皮爾遜準(zhǔn)則計算了在虛警概率Pfa=10-4和檢測概率Pd=0.6條件下,方差變化高斯過程假設(shè)下的檢測參數(shù).艦船轉(zhuǎn)舵仿真信號、魚雷出管仿真信號和水池重物落水實測脈動信號的驗證結(jié)果表明,低信噪比下,新方法對瞬態(tài)信號下降沿的判斷精確度更高,穩(wěn)健性更好.