(空軍預警學院, 湖北武漢 430019)
武裝直升機由于其本身獨特的優(yōu)點,已經成為了現(xiàn)代空襲作戰(zhàn)中的“殺手锏”,而直升機旋翼作為典型的微動目標,已經受到國內外學者的廣泛研究[1-6]。文獻[7-8]通過電磁仿真和實測數(shù)據指出旋翼目標回波中存在閃爍現(xiàn)象,但是當前窄帶雷達對旋翼特征進行提取主要是在旋翼回波不發(fā)生閃爍情況下進行的[9-13],如時頻分析方法[9]、正交匹配追蹤算法[9]、Hough變換法[9]、高階矩函數(shù)分析法[10-11]以及函數(shù)構建法[12]和自相關法[13]。其中正交匹配追蹤算法、高階矩函數(shù)分析法需要用到回波的相位信息,而Hough變換存在著計算量大,無法快速進行轉速估計,函數(shù)構建法雖然計算量較小,但是在旋翼散射點均勻分布時會失效,而自相關法必須在葉片個數(shù)已知的情況下才能準確估計旋翼轉速。因此,大部分現(xiàn)有的旋翼目標特征提取方法并不符合旋翼回波中存在閃爍的實際,閃爍情況下的旋翼目標特征提取方法仍需進一步研究。
針對上述問題,本文提出了基于閃爍現(xiàn)象的旋翼微動特征提取方法,根據旋翼回波中閃爍的時域特點和時頻域特點,通過自相關處理估計出回波的閃爍間周期,利用時頻分析和骨架提取方法判斷葉片的奇偶數(shù)及估計旋翼目標的最大微多普勒頻率,然后通過閃爍周期與葉片數(shù)和轉速的關系,估計出目標轉速的多種可能值,并同時計算不同轉速情況下的葉片尺寸估計值,然后通過先驗區(qū)間對轉速和尺寸的估計值進行初步選擇,利用選擇結果對信號在相同噪聲大小的情況下進行重構,將重構信號能量與原始信號能量最接近的估計值作為最終的特征提取結果,進而完成對閃爍條件下旋翼微動特征的提取。
基于窄帶雷達,在觀測時間較短的情況下,旋翼目標不會出現(xiàn)距離單元走動,若觀測時間較長,旋翼目標出現(xiàn)距離單元走動現(xiàn)象,此時對旋翼目標進行運動補償后,旋翼目標的運動狀態(tài)可等效為懸停狀態(tài)[14]。為簡化分析,假設目標與雷達處于同一個平面,如圖1所示,雷達到旋翼目標旋轉中心的距離為RC,假設旋翼微動目標為散射點模型,旋翼某一個葉片上的散射點P到旋翼中心的距離為r(0≤r≤l,l為旋翼葉片長度),P點到雷達的距離為RP,同時P點以角速度ω繞旋轉中心C旋轉,初始時刻P點的初始旋轉角為θ。
圖1 旋翼目標模型
旋翼目標經過脈沖壓縮后的視頻回波信號為[15]
(1)
Rij(tm)=RC+rijcos(2πfrottm+θij)
(2)
式中,rij為葉片上各散射點到旋轉中心的長度,frot為旋翼目標的轉動頻率,θij為葉片上各散射點的初始相位。
由式(1)、式(2)可知,散射點的微多普勒呈現(xiàn)余弦函數(shù)形式,而余弦函數(shù)的周期與旋翼旋轉周期是一致的。且不同葉片上散射點引起的微多普勒曲線有著不同的初相,而同一葉片散射點初相是一致的,差異在于多普勒頻率峰值大小的不同,旋翼目標瞬時微多普勒頻率為
(3)
當旋翼目標處于懸停狀態(tài)時,從推力平衡性考慮,多個旋翼葉片的轉速是相同的,轉速誤差不會超過1 r/min,表明此時旋翼回波只有一個調制周期特征[17]。理論和實際都證明,當雷達視線方向與旋翼葉片垂直時,回波的幅度最強,當偏離垂直方向時,回波的幅度銳減,從而形成閃爍現(xiàn)象[7-8]。
圖2 時域閃爍現(xiàn)象示意圖
當旋翼葉片采用式(1)所示的散射點模型時,還需要滿足單個葉片上的散射點間隔小于波長的一半,否則不會出現(xiàn)閃爍現(xiàn)象[6]。當出現(xiàn)閃爍時,其閃爍脈沖信號的形狀是由慢時間軸上的sinc函數(shù)給定的,其脈沖閃爍周期Tint[18]滿足
(4)
由于旋翼目標產生的這種具有周期性的閃爍脈沖信號,可以將慢時間維的旋翼目標回波信號看作周期平穩(wěn)或循環(huán)平穩(wěn)信號,這種特征為旋翼目標的檢測與特征提取提供了有利條件。
本文提出的基于閃爍現(xiàn)象的旋翼微動特征提取方法共分四步:第一步是尋找含有旋翼目標信號的單個距離單元,對其進行傅里葉變換估計出旋翼目標的閃爍周期;第二步是利用時頻分析和圖像骨架[3]提取的方法估計出目標的最大瞬時多普勒頻率以及旋翼目標葉片的奇偶數(shù);第三步是利用前兩步的結果,找出葉片數(shù)、轉動頻率和尺寸三者可能的組合;第四步是估計出回波信號中噪聲的功率,利用第三步的結果對回波進行重構,并在重構信號中添加與回波中相同功率大小的噪聲,然后對比原始回波信號與重構信號的能量,能量最接近的重構信號所包含的特征即為旋翼目標的特征。
在窄帶雷達回波中,旋翼目標不會出現(xiàn)距離單元走動,那么含有旋翼目標信號的通常只有一個距離單元,在噪聲中尋找到這個距離單元是本文的基本前提。因此,本文利用能量門限選擇該距離單元,門限η定義如下:
(5)
式中,Ed表示第d個距離單元的能量,Em表示所有距離單元的平均能量。當某個距離單元的能量低于門限η,則表示該距離單元只包含噪聲,當能量高于門限η,則表示該距離單元是旋翼目標所在的距離單元。對該距離單元進行自相關處理,表達式可寫為
(6)
式中,s(tm)為旋翼目標的回波方位向信號,n(tm)為噪聲,x(tm)=s(tm)+n(tm),從式中可以看出進行自相關處理后噪聲的幅度變?yōu)樵瓉淼?/M,能夠有效地提升回波的信噪比,使用該方法估計閃爍周期具有良好的抗噪性能。此時閃爍周期為
(7)
式中,τ1為R(τ)取最大值時的時延值,τ2為R(τ)取次最大值時的時延值。
采用局部分析方法即時頻分析來對旋翼目標所在的距離單元進行處理,可以得到其慢時間與微多普勒頻率的關系。STFT是最常用的基于匹配濾波的時頻分析方法,當STFT選取的窗函數(shù)為高斯窗時,稱之為Gabor變換,根據測不準原理,Gabor變換具有最佳的時頻分辨率。
對式(1)進行Gabor變換可以得到微動目標信號的時頻圖,其表達式[9]為
(8)
對具有時域閃爍現(xiàn)象的旋翼微動回波信號進行時頻分析,同樣會在時頻域出現(xiàn)閃爍現(xiàn)象,且閃爍的峰值對應旋翼目標的最大微多普勒頻率。時頻閃爍示意圖如圖3所示。
(a) 偶葉片
(b) 奇葉片圖3 時頻閃爍現(xiàn)象示意圖
由于時頻分析受到其分辨率的影響,直接提取參數(shù)的誤差較大,為了減少后續(xù)處理的運算量以及提高參數(shù)提取的準確度,從圖像域的角度出發(fā),對時頻圖圖像進行閾值處理、平滑處理和二值化處理,而后提取出閃爍位置的時頻曲線骨架,此時的圖像矩陣中只包含“0”和“1”兩個值,其中“1”出現(xiàn)的位置即為時頻曲線的骨架,此時可以通過時頻圖骨架準確提取出微動目標的最大瞬時微動頻率,并且通過判斷閃爍在正負頻率圖像上是否對稱來估計旋翼目標葉片的奇偶數(shù),若正頻率圖像與負頻率圖像對稱,即可判斷為旋翼目標的葉片數(shù)為偶數(shù),若不對稱,則為奇數(shù)。
熱點分析法屬于局部自相關分析方法,根據在一定分析規(guī)模內的所有要素,計算每個要素統(tǒng)計值,得到每個要素的z值和p值[35],通過熱點分析,可以識別出老年人口高、低值在空間上聚類的區(qū)域,公式如下[31]:
在估計出旋翼葉片的奇偶數(shù)后,可以進一步估計得到目標的尺寸:
(9)
估計得到尺寸之后,可以完整地得到關于葉片尺寸、轉動頻率以及葉片數(shù)的多種組合,這些組合具有相同的時域閃爍和時頻域閃爍特征,根據現(xiàn)有公開資料[19-20],直升機旋翼的葉片數(shù)基本是 2~8個,葉片的尺寸處于5~12 m,旋翼轉動頻率在3~7 Hz之間。利用這個先驗條件,對多種可能組合進行初步的篩選,可以得到符合實際情況的幾種組合,如圖4所示,當葉片數(shù)為偶數(shù)時,相同閃爍頻率情況下,最多只有3種可能的轉動頻率,并分別對應不同葉片數(shù),當葉片數(shù)為奇數(shù)時,最多只有2種可能的轉動頻率。因此,經篩選后可能組合最多只有3個。
(a) 偶葉片
(b) 奇葉片圖4 多值區(qū)間
利用式(5)中得到的門限η,選擇出只包含噪聲的距離單元,計算這些距離單元內噪聲的平均功率:
(10)
式中,M為脈沖數(shù),Dn為只包含噪聲的距離單元,Sn為這Dn個距離單元的噪聲。
在估計出回波信號中噪聲的平均功率后,將2.2節(jié)初步估計得到的關于葉片尺寸、轉動頻率以及葉片數(shù)的多種組合作為先驗知識,利用式(1)對不同組合的回波進行時域重構,且添加噪聲功率的大小與真實回波信號中包含的噪聲功率大小保持一致,此時重構信號表達式可寫為
(11)
對信號在相同噪聲功率大小下進行重構后,由于這些重構后的信號具有相同的時域閃爍和時頻域閃爍特征,所以依然無法找出目標組合,但是由于每個組合中的葉片數(shù)是不相同的,若每個散射點反射系數(shù)相同,散射點間隔d=λ/3(一定會出現(xiàn)閃爍現(xiàn)象),其時域回波的能量是不相同的,因此,可以根據重構信號能量與回波信號能量的相似程度來準確估計葉片數(shù)目,從而找出目標組合。
基于閃爍現(xiàn)象的旋翼微動特征提取方法具體流程如圖5所示。
圖5 算法流程圖
為了更加清晰得到本文方法對微動目標轉動頻率(frot=ω/2π)、尺寸和葉片數(shù)的提取,仿真選擇了幾種典型直升機的真實數(shù)據作為目標參數(shù),仿真信噪比為脈壓后10 dB。為了滿足方位向采樣條件PRF=4ωl/λ,以及在觀測時間內至少獲得2次時域閃爍,仿真參數(shù)設置如表1所示,目標參數(shù)如表2所示。
表1 仿真參數(shù)
表2 目標參數(shù)
(a) 時域閃爍結果
(b) 自相關處理結果
(c) 時頻域閃爍結果
(d) 骨架提取結果圖6 AH-1W特征提取結果
利用AH-1W直升機真實參數(shù)進行仿真,其仿真結果如圖6所示。從圖6(a)可以看出,此時回波中出現(xiàn)了明顯的時域閃爍現(xiàn)象,對其進行自相關處理結果如圖6(b)所示,根據圖6(b)中峰值位置與次峰值位置的時延可以估計出此時閃爍間隔為0.102 s,對其進行時頻分析可以得到圖6(c)所示結果,其閃爍出現(xiàn)的位置與時域閃爍位置相同。與圖3不同的是,由于在回波中添加了噪聲,能夠體現(xiàn)出葉片數(shù)量的正弦調制曲線完全被淹沒在了噪聲背景中,因此此時無法直接從曲線的條數(shù)對葉片數(shù)進行直觀的判斷,對圖6(c)所示結果進行骨架提取得到圖6(d)所示結果,此時可以判斷此時旋翼葉片數(shù)量為偶數(shù)并估計出目標最大瞬時微多普勒頻率為1 488.4 Hz。根據閃爍間隔、葉片數(shù)(偶葉片數(shù)只考慮2,4,6,8個)與閃爍周期的關系,可以得到閃爍周期可能值為0.204 3 s,0.408 5 s,0.612 7 s,0.817 0 s,那么對應的轉速可能值為4.896 0 Hz,2.448 0 Hz,1.632 0 Hz,1.224 0 Hz,根據式(9)可以計算出目標尺寸可能值為7.238 0 m,14.476 0 m,21.713 9 m,28.951 9 m,此時通過判斷葉片尺寸和旋翼轉速是否滿足設置區(qū)間來對可能值進行篩選,符合區(qū)間僅有一組值,葉片長度為7.238 0 m,旋翼轉速為4.896 0 Hz,葉片數(shù)為2個。由于符合設置區(qū)間的僅有一組數(shù)據,因此不需要利用估計值對信號進行重構來進一步進行判斷,此時尺寸估計誤差為0.849%,轉速估計誤差為0.082%。
利用SA341直升機真實參數(shù)進行仿真,其仿真結果如圖7所示。根據圖7(b)中峰值位置與次峰值位置的時延可以估計出此時閃爍間隔為0.028 s,從圖7(d)可看出,此時旋翼葉片數(shù)量為奇數(shù)并估計出目標最大瞬時微多普勒頻率為1 306.6 Hz,根據閃爍間隔、葉片數(shù)(奇葉片數(shù)只考慮3,5,7個)與閃爍周期的關系,可以得到閃爍周期可能值為0.168 s,0.280 s,0.392 s,那么對應的轉速可能值為5.952 4 Hz,3.571 4 Hz,2.551 0 Hz,根據式(3)可以計算出目標尺寸可能值為5.240 5 m,8.734 2 m,12.227 9 m,此時通過判斷葉片尺寸和旋翼轉速是否滿足設置區(qū)間來對可能值進行篩選,符合區(qū)間有兩組估計組,葉片長度為5.240 5 m,旋翼轉速為5.952 4 Hz,葉片數(shù)為3個以及葉片長度為8.734 2 m,旋翼轉速為3.571 4 Hz,葉片數(shù)為5個。由于符合設置區(qū)間的有兩組數(shù)據,還需要進行進一步的篩選,因此分別利用兩組估計值對信號進行重構,重構信號的能量與原始回波的能量相近的一組估計值即為目標特征的估計值。第一組估計值重構信號與原信號能量偏差為4.96%,第二組估計值重構信號與原信號能量偏差為176.66%,因此選擇第一組估計值作為輸出,此時與真實參數(shù)相比,尺寸估計誤差為0.181%,轉速估計誤差為0.793%。
(a) 時域閃爍結果
(b) 自相關處理結果
(c) 時頻域閃爍結果
(d) 骨架提取結果圖7 SA341特征提取結果
由于其余3種直升機特征提取的仿真方法相同,僅僅只有仿真參數(shù)設置的不同,為簡化分析,5種直升機的參數(shù)估計結果如表3所示,估計誤差如表4所示。
表3 參數(shù)估計結果
表4 估計誤差
從表3和表4可以看出,本文方法能夠準確估計出目標的葉片數(shù)目,并且葉片尺寸誤差和旋翼轉速誤差都較小,具有較強的穩(wěn)定性,并且利用本文方法仿真時原始信號與重構信號的初相均為隨機相位,因此本文方法可以在不知道目標初相情況下對旋翼目標特征進行提取。
本文提出了一種基于閃爍現(xiàn)象的旋翼微動特征提取方法,該方法通過自相關處理估計出旋翼目標回波的閃爍間隔,并利用閃爍間隔與葉片數(shù)的關系估計得到多種可能的閃爍周期,通過對其進行時頻分析并提取相應的時頻圖骨架,判斷葉片的奇偶數(shù)及估計旋翼目標的最大微多普勒頻率,并利用旋翼目標葉片數(shù)量、葉片尺寸和旋翼轉速所處的大致范圍作為先驗知識對估計值進行初步篩選,然后在保持與原始信號中噪聲功率相同情況下,利用估計值對回波信號進行重構,從能量差異的角度出發(fā)對估計值進行第二次篩選,將原始信號能量與重構信號能量最接近的估計值輸出作為最終的估計值。通過仿真結果可以看出,本文方法可以準確地估計出旋翼目標的葉片數(shù),旋翼轉速和尺寸誤差均在可接受范圍內。