基于諧波集檢測的飛行目標水下聲探測算法研究

穆森， 李京華， 張恒， 陳果

(1.西北工業大學 電子信息學院， 陜西 西安 710072; 2.陸軍研究院 近地面探測技術實驗室， 江蘇 無錫 214035;3.宜昌測試技術研究所， 湖北 宜昌 443003)

0 引言

反潛直升機對于水下潛艇具有很大的威脅，為了提高水下潛艇對抗反潛機的能力，有必要盡快解決水下對空探測的技術問題。反潛直升機飛行時，旋翼輻射噪聲的聲波穿越空- 海界面入水后形成的頻譜具有諧波特性，為水下探測空中反潛飛機提供了可行途徑。

目前關于被動聲探測的研究多集中在同一種介質內，如空氣介質中的聲探測或海水介質中的水下聲探測。魏麗萍等[1]給出了諧波集(HS)檢測算法，并成功地對實測的空氣中直升機聲信號進行了HS檢測；馬寧等[2]提出了一種基于小波包和神經網絡的檢測算法，可以較好地檢測識別不同型號的空氣中直升機信號；Salloum等[3]開發并建立了一個聲學空氣檢測(AAD)系統，實現了空中目標的檢測和跟蹤，并對小型飛機、直升機和超輕飛機進行了分類；周關林等[4]提出隨機共振水聲信號的檢測，為檢測水下微弱信號提出了新思路。關于空中飛行目標的水下聲探測，即空氣- 水跨界聲探測問題的研究大多集中在空中聲源激發的水下聲場理論建模方面，而空中目標的水下跨界檢測算法則研究很少。例如：以色列研制了潛艇反直升機和低空慢速飛行目標的聲探測浮標[5]，用于檢測海面上空目標的聲信號；韓建輝等提出在水下檢測空中聲源的平滑偽維格納·維爾分布- 高度門限(SPWVD-HT)算法[6]、雙譜冪次法[7]，但所檢測聲信號均采用空中運動聲源激發的水下聲場建模仿真信號。

本文提出基于HS檢測與最強譜線檢測相結合的檢測算法，根據背景噪聲進行自適應門限計算，并對功率譜中HS進行檢測，若HS不存在，則將最大線譜譜峰與門限進行比較，判斷目標是否出現，以實現空中移動聲源的水下跨界檢測。應用該算法對水下單水聽器接收的螺旋槳飛機飛行噪聲實測信號進行了目標檢測仿真，并對該算法的檢測性能進行了分析。

1 飛行目標噪聲及入水傳播特性分析

帶有螺旋槳的飛機在飛行時會有較大的噪聲，尤其是直升機飛行時，其輻射噪聲主要由主槳、尾槳和發動機產生。噪聲可分為氣動噪聲和機械噪聲兩大類[8-9]，機械噪聲屬于高頻噪聲，傳播到空氣中很快被吸收，其噪聲源只對飛機近場噪聲起作用；氣動噪聲主要集中在頻譜的低頻部分，對遠場噪聲起決定性作用，功率較強且持續不斷。由于直升機旋翼在結構上有周期轉動機制，其遠場噪聲功率譜是由具有諧波關系的離散線譜疊加于寬帶連續譜上構成的，呈諧波關系的頻率分量構成一個HS。例如，直升機由于有主槳和尾槳，通常由兩個不同基頻的HS組成，某直升機的時域波形和功率譜如圖1所示。圖1(b)中的慢變曲線是分離出的趨勢項波形，其中*和°表示頻譜中離散線譜構成的兩個HS，分別對應主槳槳葉通過頻率(BPF)及其各次諧波頻率，以及尾槳BPF及其各次諧波頻率。

螺旋槳BPF用fBPF表示，它對于某一種螺旋槳飛機而言是一個固定常數，可由(1)式計算[10]：

fBPF=iwr/60，

(1)

式中：i為諧波次數；w為槳葉個數；r為轉速。當飛機改變速度和狀態時，螺旋槳旋轉速度不變，fBPF的值也不變。

一般檢測問題中假設噪聲為平穩高斯隨機過程，然而對于水下探測，由于目標和傳感器周圍的環境條件(尤其是風速和海浪)影響，有可能改變噪聲的統計特性，成為非平穩非高斯過程。為便于處理，基于微分法思想，假設噪聲在一個非常短的時間內保持平穩，可將其視為平穩信號，故信號處理時將數據分為若干個短數據段進行。

空氣中點聲源發出的聲波傳播入水過程會有較大衰減。按照波動理論[11]，在滿足空氣中點聲源高度要求和聲傳播路徑的相位變化條件基礎上，空氣中高度為h的點聲源可以等效為水中深度為d的點聲源，同時其聲壓級衰減20lg(kαd)，其中kα為聲波在空氣中的波數。圖2所示為空氣中點聲源的水下等效源[11]示意圖，直觀地表明了空氣中點聲源可等效為水下聲源加上空氣聲源到水面的傳播衰減。

因此可將聲源空氣- 水跨界傳播問題轉換為水下聲傳播問題，進而分析水下聲場。根據對波動方程解法的不同，本文選取基于拋物方程的RAMGEO模型進行仿真計算，該模型改進了拋物方程模型，并與穩定的自身初始場結合[12]，適用于淺海低頻噪聲環境。本文具體仿真參數如下：聲源深度10 m，聲源頻率200 Hz，接收水聽器在水下呈垂直線陣排布，海水深度250 m. 利用RAMGEO模型仿真計算得到的聲源水下傳播衰減趨勢如圖3所示。

由圖3(a)可以看出，聲源在水下環境水平距離200 m以內的范圍大致呈線性衰減，水平距離每增加約10 m，聲壓級降低2 dB，在水平距離超出200 m范圍后，聲壓級衰減明顯減小。這表明只要聲源聲壓級足夠大，聲波就可以實現在海水介質中的遠距離傳播。聲源在深度范圍上的衰減如圖3(b)所示。由圖3(b)可以看出，在淺海環境下水深100 m之后，聲波傳播衰減趨勢發生了較大變化，這可能是因為海底聲波反射以及底層泥沙沉積物引起的聲能損耗造成的。

由于海上直升機目標往往帶有較大的噪聲聲壓級，按照等效源理論完全可以滿足聲波的遠距離傳播，從而為海上戰場聲目標的探測提供了有力手段。

2 基于自適應門限的HS檢測算法

2.1 算法原理

目標檢測基于以下兩種假設：1)H0-沒有目標；2)H1-目標出現。檢測算法原理框圖如圖4所示。

由圖4可知，本文提出的檢測算法是一種在頻域中采用HS和最強譜線檢測相結合的檢測方法，首先需要計算采集信號功率譜和背景噪聲功率譜擬合曲線，并根據擬合曲線計算檢測門限(DT)，在功率譜中搜尋線譜，將線譜按譜峰高低進行降序存儲，若在搜尋的若干個線譜中有2個以上存在諧波關系，則認為有HS并判斷為目標出現，輸出狀態為H1；若沒有HS，則將譜線中最強的譜線峰值與DT進行比較，該門限由背景噪聲功率譜擬合曲線確定：若大于DT，則判斷為目標出現、狀態為H1，否則判斷為沒有目標、狀態為H0.

2.2 HS檢測

HS檢測時，必須首先對水聽器接收信號時間序列采用Welch法[13-14]計算功率譜：

(2)

式中：PSD(f)為接收信號的功率譜；M為數據窗長度；L為數據段總數；xb(n)為第b個數據段中的第n個數據，b=1,2,3，…；d(n)為第n個數據的窗值。

在HS檢測過程中，線譜的判斷是算法中至關重要的一步，如第1節所述，螺旋槳飛機聲信號的特征頻率(基頻)及其諧波集中在中低頻段，因此線譜檢測過程中首先需要選擇相應的頻率范圍，然后進行線譜搜索和判定。在選定頻段內以滑動窗搜索所有譜峰，按譜峰值大小降序排列，選用前幾個線譜譜峰做HS檢測。線譜的搜索及判定方法如下：

1) 線譜值必須是功率譜中的局部極大值Lmax，該極大值對應的頻點為fk(k=1,2,…,K，K為線譜數量)，局部極大值判斷過程是：將當前檢測點k的功率譜值分別與其左、右相鄰點k-1和k+1的功率譜值求1階差分Δl和Δr，若滿足Δl<0且Δr>0，則Lmax暫定為此鄰域內的局部極大值；

2) 若步驟1中Lmax的大小在k點的左、右峰寬Lsize鄰域范圍內均為最大值，則Lmax暫定為線譜峰值；

3) 類似于步驟2，若在Lmax點的左、右峰寬Lsize鄰域范圍內，分別存在極小值Lminl和Lminr，且比值2×Lmax/(Lminl+Lminr)大于某個閾值，則可認為譜峰Lmax是顯著的，可將fk頻點處的Lmax作為線譜；

4)Lmax譜線的峰高定義為

hs=|Lmax-(Lminl+Lminr)/2|，

(3)

從而線譜的選取兼顧了譜峰左右邊界、峰寬和峰高等特征判據，對線譜特征的辨識度較高。

HS檢測按照以下方程[15]判斷：

|fk/fm-imk|≤ε，

(4)

式中：fm為假定基頻，m為主槳和尾槳基頻個數，m∈(1,2)；imk為諧波次數，imk=1,2,3；ε為選定的閾值。

若檢測到若干線譜對應的頻率與該假定基頻有明顯的諧波關系，則認為檢測到一個HS，該假定基頻即為HS的基頻。HS構建策略如下：

1) 根據線譜搜索及判定方法進行線譜檢測；

2) 根據線譜峰高hs的大小，選出若干峰高較大線譜的頻點fk；

3) 將K個線譜的頻點依次作為假定基頻fm，例如將第m個頻點作為基頻，其他K-1個頻點fk用(4)式進行諧頻檢測，即只要(4)式中imk取值為2或3時fk滿足(4)式，即可認為fk與fm具有諧頻關系，它們對應的線譜具有諧波關系。

在檢測中對于某個基頻，有時會出現不止一個HS，這時需要用權函數決定取舍，對每個HS內所有頻率的能量進行計算并作為權值，選擇具有最大能量的HS，與其對應的基頻作為提取的螺旋槳飛機特征頻率。

在HS檢測中，由于諧波信號特性不同，例如頻率檢測范圍、線譜個數、線譜譜峰寬度、HS內諧波數目以及閾值等參數，需要對實際信號進行分析后確定這些參數的最佳值，從而使諧波檢測效果最佳。

2.3 自適應DT的計算

由2.1節的算法原理可知，在進行目標檢測時首先進行HS檢測，若沒有HS則尋找線譜中最強的譜峰，并與DT比較，若目標出現，則由目標發出的噪聲級再加上背景噪聲的噪聲級就會大于該DT. DT計算公式為

TDT=α+P(f)，

(5)

式中：TDT為實際計算檢測門限DT；α為門限調節因子；P(f)為無目標時背景噪聲信號功率譜的多項式擬合曲線。在實際環境中，背景噪聲的大小是變化的，為了使檢測門限能自適應地隨之調整，本文提出采用α來動態調節門限，它可以根據無目標時不同環境下測得的P(f)來計算，使DT成為隨環境變化的自適應門限。

一般情況下，實際背景噪聲發生變化的主要原因是環境變化使得某一部分頻率分量的噪聲能量增加或減少，從而造成背景噪聲功率譜曲線[16]上某些頻率分量的功率譜幅度增大或減小，因此這種變化是一種非線性變化，在缺乏先驗知識的條件下，對這種變化進行定量描述是很困難的。因此本文提出采用最小二乘法對背景噪聲功率譜進行擬合，從而對背景噪聲級進行近似定量描述，α和DT的計算過程如下：

1) 計算第l幀背景噪聲nl(t)的功率譜，并取分貝數后為

Pn(f,l)=10lg(PSDn(f,l))，

(6)

式中：Pn(f,l)為第l幀背景噪聲的功率譜分貝數；PSDn(f,l)為第l幀背景噪聲的功率譜幅值。

2) 對Pn(f,l)采用最小二乘法進行曲線擬合，設擬合多項式的階數為k，第l幀的擬合公式為

P(f,l)=p1fk+p2fk-1+…+pkf+pk+1，

(7)

式中：P(f,l)為第l幀背景噪聲的功率譜擬合曲線；pk為擬合系數，擬合曲線如圖5中背景噪聲中的連續實線。

3) 門限調節因子α的計算示意如圖5所示，圖5中背景噪聲功率譜最大值對應的頻點為fα，則α取值為背景噪聲功率譜的最大值(圖5中黑色實心圓點)與其擬合曲線P(f,l)在fα頻點取值的差值，第l幀門限調節因子α的計算公式為

α(l)=max[Pn(f,l)]-P(fα,l) ，

(8)

式中：max[Pn(f,l)]為第l幀背景噪聲功率譜的最大值；fα為第l幀背景噪聲功率譜最大值對應的頻率值；P(fα,l)為第l幀背景噪聲功率譜頻率值為fα時對應的擬合曲線值。

4) 根據第l幀背景噪聲計算的判決門限DT為

TDT(f,l)=α(l)+P(f,l).

(9)

從(9)式可知，這樣計算出的門限值是一個隨頻率變化的動態曲線，且隨著時間的變化而變化。

3 實測信號檢測實驗及結果分析

3.1 實測信號檢測實驗

采用實測數據對本文算法進行驗證。實測數據的時域波形和時頻譜圖如圖6所示。圖6中3組數據文件Data1～Data3的記錄時間長度分別為622 s、974 s、651 s，3組數據為水下5 m單水聽器接收的某中型螺旋槳飛機水面上空飛行時的飛行噪聲信號，飛機勻速直線往返飛行，飛行高度約為150 m，數據記錄采樣頻率為10 kHz.

從時頻譜圖中可以同時觀察信號在時域和頻域上的信息，包括信號的頻譜特征以及頻率隨時間的變化。以實測信號Data1為例，分析圖6所示時頻譜圖中在200 s左右目標出現的時刻，縱軸所示信號的能量急劇增強并在200 Hz上下伴有明顯的線譜特征，線譜之間存在諧頻特征，因此對該實測信號進行諧波檢測是完全可行的。

對實測信號做預處理后進行諧波檢測，每次取1 024個點作為1個數據段，將其分為互相重疊50%的子段，加窗作快速傅里葉變換(FFT)，所得功率譜的頻率分辨率約為1 Hz. 選取子段中前10個最大的譜峰進一步作HS檢測，若出現3個以上諧波即1個HS時，則認為成功檢測到目標信號HS. 這里選取實測信號Data1的子段信號，檢測到的目標信號諧波關系如圖7所示，具有諧波關系的線譜在圖中用圓點示出。對3次實測信號Data1～Data3結合自適應門限調節因子作諧波檢測，發現目標時檢測到的HS如表1所示。

Tab.1 Harmonic set of measured signals Hz

3.2 實測信號檢測結果分析

以第1次實測信號Data1為例進行分析，如表1所示，HS檢測的線譜頻點有：192.9、249.0、273.4、397.9、500.5、556.6、1 018.1，其中存在諧波關系的HS有{192.9，397.9}，{249.0，500.5，1 018.1}，{273.4，556.6}。在實測過程中存在船舶噪聲干擾，如圖8所示為采集的一段無飛機的純船舶噪聲段，由其頻譜分析可知此船噪聲主要集中在245 Hz左右，對其單獨進行HS檢測，檢測到的HS{244,488,967}與Data1中的HS{249.0,500.5,1 018.1}基本一致，從而目標檢測時可用此HS排除船舶噪聲的干擾，并以此檢測識別飛機目標。

由于實測信號受風浪、船只等外界因素的影響，對3次實測信號檢測到的諧波線譜不完全相同。但由于是同一架飛機，HS檢測所得諧波基頻相差不大，表明諧波檢測算法在不斷變化的環境條件下仍能保持穩定的檢測性能。

根據圖4的檢測算法原理，若某幀信號檢測的線譜間不存在諧波關系，則尋找線譜中最強的譜峰，與DT進行比較來判斷有無目標，若大于DT則輸出H1，否則判斷為H0，繼續觀測。為了保證檢測效果，若連續檢測10次都判斷為H1，則認為發現目標，其中DT為發現目標前最后一幀信號計算所得的門限值，檢測概率的計算公式為

(10)

使用該方法進行檢測的結果如圖9所示，根據圖9的檢測結果可以獲得實測信號的相關預警信息，在檢測到目標后即發出預警信息，同時若檢測到HS，則可利用HS的基頻對目標進行輔助識別，判斷目標是否為飛機。按照實驗中該型號飛機的正常巡航速度v=36 m/s計算，預警距離為

(11)

式中：t1為發現目標的時刻；t2為目標過頂的時刻。

由(11)式可得預警距離如表2所示。從圖9和表2可知，從發現目標開始直至目標到達最靠近點(CPA)[17]，可以保持較遠的預警距離，并在檢測到目標后保持較高的檢測概率。

表2 實測信號預警信息

4 結論

本文針對螺旋槳飛機的水下聲探測，根據其輻射噪聲信號的諧波特性，提出了HS檢測與最強譜線檢測相結合的檢測算法；采用門限調節因子α自適應地調整DT，用實測信號進行了目標檢測試驗驗證。主要得出以下結論：

1) 螺旋槳飛機噪聲信號具有諧波特性，通過基于拋物方程的RAMGEO模型進行仿真計算表明，可將聲源空氣- 水跨界傳播問題轉換為水下聲傳播問題，從而按照等效源理論進行空中目標的水下探測。

2) 在HS檢測過程中，線譜檢測是至關重要的一個環節，可以通過選取兼顧譜峰左右邊界、峰寬和峰高等特征來提取線譜。

3) 若某幀信號檢測的線譜間不存在諧波關系，則可尋找線譜中最強的譜峰，從而利用基于自適應門限調節因子α的最強譜線檢測方法增強HS檢測方法的穩定性，有效提高目標的檢測距離。