基于時頻面霍夫變換的諧波族檢測算法

唐玉志,劉 瑜,郭 艷,吳 娜

(1.西北工業(yè)大學,陜西西安 710072;2.西北核技術研究所,陜西西安 710024)

0 引言

直升機和無人機在現(xiàn)代戰(zhàn)爭中發(fā)揮著重要作用,可用于戰(zhàn)場的監(jiān)視、運輸、攻擊等,反潛直升機等配備的聲納系統(tǒng)的性能逐漸升級,對潛艇的威脅也越來越大。對直升機(包括螺旋槳驅動的反潛機)和無人機的探測和預警,是一項重要的戰(zhàn)術任務。低成本、可靠、探測距離遠、布設方便的探測系統(tǒng),部分國家已有設備正式列裝,其性能也在持續(xù)提高中。

直升機和大部分的無人機等由螺旋槳驅動,其發(fā)出的噪聲信號包含有一組或多組諧波信號,每一組諧波信號的譜線在頻域上呈等間隔分布,強度處于同一量級。將具有這種諧波關系的一組信號稱為諧波族信號(harmonic signal set)。在由聲學方法進行直升機等目標的探測時,以往多據(jù)其能量最高的譜線或者頻譜的整體形狀提取特征,結合模式識別算法進行檢測[1-5]。這類算法往往沒有有效利用其諧波性質,在強干擾或噪聲背景下檢測能力大幅降低。

本文使用諧波小波變換方法求得含噪聲諧波族信號(以模擬直升機和無人機等的諧波族信號)的時頻分布,在此基礎上進行Hough變換,基于頻率分集技術[6-8]中信號融合處理技術構建檢測量,并在給出檢測結果的同時,給出諧波族的基頻及其變化率。

1 信號的時頻分布及Hough變換

1.1 基于諧波小波的時頻分布

頻率范圍為m～n Hz(m,n為實數(shù),且0

式(1)表示的是一個強度在m ～n Hz范圍內為單位幅度的諧波小波。對于復信號 x(τ),加時間窗win(τ)以改善時間定位能力,其諧波小波變換系數(shù)為[10-11]:

hwt(m,n,t)代表了x(τ)與諧波小波w(t)的相關程度。hwt(m,n,t)2可作為信號在該帶寬范圍內能量的量度。在選定的時間和頻率范圍內變化(m,n,t)的組合,將 hwt(m,n,t)2排列在二維時頻平面上,即得到x(τ)的一個時頻分布[12-13],該二維平面即稱為時頻面。

1.2 二維圖像的Hough變換[12-13]

設有二維笛卡爾坐標系下的兩點 A(x,y)和B(x′,y′)。過兩點的直線為 AB,過原點的垂線長度為R,與橫坐標Ox的夾角為θ。過AB兩點的直線可以表示成y=kx+b,有:

代入直線AB的表達式,可以得到

圖1 Hough變換示意圖,(左為圖像空間,右為參數(shù)空間)Fig.1 Illustration of Hough transformation(lef t one is the image domain,right one is the parameter domian)

對于固定的點 A(x,y),當參數(shù)(R,θ)變化時,其軌跡為一條曲線,每一點代表了過點 A(x,y)的一條直線的R和θ。假設曲線上與直線AB相對應的點設為(R0,θ0),則直線 AB上所有點在參數(shù)空間的曲線都將通過點(R0,θ0),經(jīng)加權后必然在點(R0,θ0)處形成極大值。通過適當?shù)乃阉魉惴烧页鰳O值點(R0,θ0),它唯一的定義了直線 AB并可直接得到其截距和斜率。

2 諧波族信號時頻面的Hough變換

將使用諧波小波變換得到的信號時頻分布作為Hough變換的二維輸入圖形,以強度為變換域中各點的加權值。直升機、無人機等目標的諧波族噪聲包含多條譜線,在時頻面上,表現(xiàn)為多條直線,經(jīng)Hough變換以后,參數(shù)空間中出現(xiàn)多個極值點。極值點的個數(shù)對應于諧波信號的階數(shù)。

設諧波族信號基頻為 f 1(t),最高階數(shù)K,第n階諧波頻率f n(t),n=1,2,…,K。時頻分布像元數(shù)為M×N,時間和頻率分辨率分別為Δt和Δf。以時間中心和頻率0作為時頻分布的坐標原點,時頻分布都位于上半平面。

2.1 諧波信號頻率平穩(wěn)時的Hough變換

在直升機等目標飛行平穩(wěn)、且沒有多普勒效應時,噪聲譜線頻率不隨時間變化,時頻圖中各階諧波時頻線的斜率為0,在參數(shù)空間(R,θ)中,有:

(式中,θn為第n階諧波對應譜線的垂線與縱軸的夾應的(Rn,θn)各點等間距分布在縱坐標軸線上。

2.2 諧波族信號頻率變化時的Hough變換

諧波族信號頻率變化時,在某一時間長度內,可認為信號頻率是接近于線性變化的,變化率頻率的時變由多普勒示目標與測點的徑向相對速度。諧波譜線成為斜率與階數(shù)成正比的一組直線,其斜率為:

可在時頻分布計算中合理選擇參數(shù),使得k′n=tan(θn)?θn近似成立,則在Hough變換的參數(shù)空間(R,θ)中,θn落在一系列等間隔的點上。

[89]Maung Aung Myoe, In the Name of Pauk-Phaw: Myanmar’s China Policy Since 1948, Singapore: Institute of Southeast Asian Studies, 2011, p.190.

各階諧波譜線的垂線與縱軸的夾角分別為:θn,n=1,2,…,K 。在t=0的時刻,有 f 1=f 2/2=f 3/3=…,由幾何關系有:Rn=f n cosθn。理論上,參數(shù)空間(R,θ)中R的分布不再是等間距的。但當諧波族信號頻率變化率較低,k′n=tanθn?θn成立時,也有cosθn≈1成立。這時有:Rn=nf 1,分布成為等間距的。Hough變換中,Δt和Δf本身就有一定寬度,該近似帶來的誤差不致引起極值點大的偏移。

因此,在時頻面上譜線斜率較小的情況下,各階譜線在參數(shù)空間(R,θ)中對應的點就散布在一條直線附近。這條直線過參數(shù)空間的零點,可以表示為:

式中:

圖2 頻率時變諧波信號的時頻分布與Hough變換參數(shù)的示意圖Fig.2 Time-frequency distribution of frequency-shifting signal and parameter of the Hough transformation

式(7)中的k是Hough變換參數(shù)空間中多個極值點擬合為一條直線時的斜率,是一個平均值,用于進行諧波信號頻率變化率的計算時,可以得到比利用單一譜線更好的計算結果。

3 Hough變換參數(shù)空間中的諧波檢測算法

由2.2節(jié),要得到諧波基頻,需要先在(R,θ)空間的上半平面,即在射線與 R軸的角度 ψ∈取不同值的射線上進行搜。

3.1 ψ的搜索算法

對于某一確定的ψ,射線上各點的幅度表示了時頻面上各點經(jīng)Hough變換給予的能量的總貢獻。設射線上離散點總數(shù)為M 0個,第m個點的能量為p m,m=1,2,…,M 0。定義射線上的總能量Eψ和熵函數(shù)Sψ分別為:

熵函數(shù)Sψ具有以下性質:

3)當p m的分布位于以上兩種情況之間時,0

3.2 諧波基頻搜索算法

在確定ψ,從而確定了信號中各階諧波在參數(shù)空間中所在的射線以后,就可以在此射線上搜索諧波信號基頻的確切位置。在射線上定義函數(shù):

K為諧波族中的最高階數(shù)。使Pm取得最大值的點,就是諧波族信號基頻的時頻線在參數(shù)空間(R,θ)中的映射,據(jù)此可進一步由式(5)和式(6)得到目標的部分運動學參量。

3.3 Pm的意義

在無線通訊領域,由于電磁波在傳輸過程中的多途隨機衰落,導致信息在某些時刻不能準確檢測?，F(xiàn)常采取頻率分集的方法,將信息在多個載頻上同時發(fā)送,在檢測端采取某種融合方式將這些頻率攜帶的信息經(jīng)綜合后解碼,避免了單一頻率的隨機衰落引起的信息丟失。

直升機或無人機噪聲為諧波族信號。Pm反映了諧波族信號中各階分量的貢獻的和,將其用于對目標的檢測,可以得到比基于某單個頻率分量更好的檢測結果。Pm實際就是頻率分集技術中各頻率貢獻以等增益方式融合的結果。

3.4 算法的改進

Hough變換和參數(shù)空間中ψ的搜索算法計算量很大。為節(jié)省計算時間,可以先在時頻面上每一個時刻的譜分布上應用式(9)的諧波基頻搜索算法,將結果重新排列成一個關于t和基頻m的時頻面,其頻率范圍縮小為基頻可能的取值范圍。再將此時頻面進行Hough變換,找出極值點,就可以得到信號基頻的參數(shù)。經(jīng)這樣變換,可以使得Hough變換的計算量大幅降低為原來的1/K,并且在參數(shù)空間(R,θ)中只需要搜索極值點,而不需要進行一次二維搜索和一次一維搜索。

3.5 Hongh變換網(wǎng)格與檢測精度間的關系

在Hough變換的參數(shù)空間(R,θ)中,若網(wǎng)格劃分的較細密,則時頻面上某一點變換后對應的分布區(qū)域將減少,時頻面上的直線在參數(shù)空間中對應點的能量將更加集中,幅度更高,有利于進行檢測,也可以更加準確地定位基頻和頻率變化率。

網(wǎng)格加密帶來的代價是計算量的增加。

4 仿真計算實例及分析

4.1 平穩(wěn)諧波族信號的檢測算例

圖3中是兩組平穩(wěn)諧波族信號的時頻分布、Hough變換和據(jù)式(9)的基頻搜索結果。諧波族信號的基頻為10 Hz,最高階數(shù)6。時頻面中,Δt=1/15 s,Δf=0.2 Hz,Hough變換網(wǎng)格數(shù)500×500,θ范圍為-π/4～π/4。圖中上面一組是信噪比為14 d B時的結果,信號長度為5 s。下面一組是信噪比為-4 dB時的計算結果,信號長度為10 s。由于信號譜線斜率為零,諧波搜索只需在θ=0的軸線上進行,而不必再進行ψ的搜索。基頻搜索的結果為M=50,與時頻圖中譜線的位置一致,對應于基頻頻率M?Δf=10 Hz。

圖3 信噪比分別為14 dB和-4 dB時,諧波族信號的時頻分布、Hough變換和基頻搜索結果Fig.3 Time-frequency distribution,the Hough transformation,and base frequency result of the harmonic signal set while SNR=14 d Band SNR=-4 dB

4.2 頻率時變的諧波族信號的檢測算例

圖4 為一組頻率時變的諧波族信號的時頻分布及其Hough變換圖。Hough變換網(wǎng)格數(shù)為200×300(θ范圍為-π/4～π/4),信號頻率變化率引起的時頻面上基頻譜線的斜率為-0.045,截距為8.67。根據(jù)Hough變換參數(shù)空間的極值點,計算得到相應的結果為-0.049和8.11,誤差分別為8.9%和6.5%。

圖4 頻率時變諧波族信號的時頻圖及Hough變換Fig.4 Time-frequency distribution and Hough transformation of frequency-shifting harmonic signal set

4.3 仿真結果分析

在信噪比較高時,從圖3第一組和圖4的結果中,可以清楚看出信號時頻面強度分布與Hough變換參數(shù)空間中分布的對應關系,通過搜索算法得到了信號基頻的準確值(不可避免的,搜索過程中在基頻的整數(shù)倍處也有幅度較低的峰值出現(xiàn))。在圖3的第二組結果中,即使信噪比降到-4 d B,在時頻面上已很難分辨信號,但在諧波搜索曲線上,基頻和其整數(shù)倍位置處仍有較其他位置高的幅度可供檢測。

4.4 誤差分析

圖4中計算誤差主要有以下原因:1)由于時間-頻率分辨率互相制約,即使信號頻率為連續(xù)均勻變化,時頻面上的譜線也不是細的連續(xù)直線,而是具有一定寬度的階梯形狀,這是誤差的最主要來源,它與時頻不確定原理相對應,難以大幅度消除;2)Hough變換過程中,Δt和Δf具有一定寬度,網(wǎng)格劃分的疏密程度也有影響;3)多普勒效應引起的頻偏不是嚴格的線性關系,在時頻面的兩端帶來一定誤差;4)k′n=tanθn?θn的近似帶來少量誤差。圖4中基頻搜索結果與信號設定的參數(shù)基本無數(shù)值誤差,只存在由于第1)和第2)項帶來的小量不確定寬度。

5 結論

基于時頻面Hough變換的諧波檢測算法有以下優(yōu)點:1)符合無線通訊技術中可靠傳輸信息的頻率分集模型,有效利用信號的多階諧波分量;2)有效利用時頻面上的時間積累效應;3)在信號頻率隨時間變化時,可以利用多階信號的信息,求得頻率變化率,精度比利用單階信號時高;4)無論信號頻率變化與否,在Hough變換參數(shù)空間的特征檢測都較簡單。

仿真結果表明,此法在低信噪比下仍有較好的檢測性能,將會在螺旋槳驅動的直升機、無人機等軍事目標的檢測和識別領域中具有較好的應用前景。

[1]劉貫領.聲目標識別方法研究[D].南京:南京理工大學,2003.

[2]李在庭,高德勇.直升機聲信號特征提取和識別技術[J].兵工學報,1996,17(1):55-59.LI Zaiting,GAO Deyong.Research on feature extraction and recognition technology of helicopter acoustic signal[J].Acta Armamentarii,1996,17(1):55-59.

[3]呂國云,許學忠.戰(zhàn)場目標被動噪聲識別技術[J].探測與控制學報,2001,23(4):22-25.LV Guoyun,XU Xuezhong.Passive sound recognition technology for thetargetsin battlefield[J].Journal of Detection&Control,2001,23(4):22-25.

[4]祝龍石,莊志洪,張清泰.戰(zhàn)場聲目標噪聲特性分析[J].探測與控制學報(原現(xiàn)代引信),1996,18(2):57-61.ZHU Longshi,ZHUANG Zhihong,ZHANG Qingtai.Analyze on the character of the acoustic target in the war field[J].Journal of Detection&Control,1996,18(2):57-61.

[5]唐玉志,劉瑜,郭艷,等.直升機諧波信號的多頻循環(huán)平穩(wěn)檢測[J].探測與控制學報,2009,31(5):33-37.TANG Yuzhi,LIU Yu,GUO Yan,et al.Multi-frequency Cyclostationary Analysis of Harmonic Signal of Helicopter[J].Modern fuze,2009,31(5):33-37.

[6]韓崇昭,朱洪艷,段戰(zhàn)勝.多元信息融合[M].北京:清華大學出版社,2006.

[7]廖明,譚曉衡,張志華.頻率分集MC-CDMA在瑞利衰落信道下的性能分析[J].重慶郵電學院學報,2004,16(1):47-50.LIAO Ming,TAN Xiaoheng,ZHANG Zhihua.Performance of multi-carrier frequency diversity CDMA in Rayleigh fading channel[J].Journal of Chongqing University of Posts and Telecommunications,2004,16(1):47-50.

[8]唐秋玲,覃團發(fā),張淑儀.頻率選擇性衰落無線信道下的發(fā)送分集技術[J].南京大學學報(自然科學),2004,40(4):454-461.TANG Qiuling,QIN Tuanfa,ZHANG Shuyi.Transmit diversity scheme in the wireless channel of frequency-selective fading[J].Journal of Nanjing University(Natural Sciences),2004,40(4):454-461.

[9]David E Newland.Harmonic wavelet analysis[J].Mathematical and Physical Sciences,1993,443:203-225.

[10]李舜酩,許慶余.微弱振動信號的諧波小波頻域提取[J].西安交通大學學報,2004,38(1):51-55.LI Shunming,XU Qingyu.Harmonic wavelet extraction for weak vibration signal in frequency domain[J].Journal of Xi＇an Jiaotong University,2004,38(1):51-55.

[11]唐玉志,董永峰.噪聲和振動信號的諧波小波時頻分布[J].噪聲與振動控制,2009,29(4):42-45.TANG Yuzhi,DONG Yongfeng.Time-frequency distribution of acoustic and vibration signal based on the harmonic wavelet transformation[J].Noise and Vibration Control,2009,29(4):42-45.

[12]岡薩雷斯.數(shù)字圖像處理[M].北京:電子工業(yè)出版社,2003.

[13]徐毓,楊瑞娟,李鋒.Hough變換與數(shù)據(jù)融合[J].現(xiàn)代雷達,2001(6):16-18.XU Yu,YANG Ruijuan,LI Feng.Hough transform and data fusion[J].M odern Radar,2001(6):16-18.