周毅恒, 楊 軍, 夏賽強(qiáng), 呂明久

(空軍預(yù)警學(xué)院, 湖北 武漢 430019)

0 引 言

直升機(jī)飛行過程中,其旋翼高速轉(zhuǎn)動,會使雷達(dá)回波信號產(chǎn)生額外的頻率調(diào)制,稱為微多普勒現(xiàn)象,此時回波信號中會包含直升機(jī)的旋翼長度、旋轉(zhuǎn)頻率等信息,如果能夠較為精確地提取出這些物理特征,則可實(shí)現(xiàn)對機(jī)型進(jìn)行穩(wěn)健識別,對國土安全具有重要意義。

對于旋翼類目標(biāo)的微動參數(shù)估計問題,目前主要有以下幾種方法:一是基于變換域的微動參數(shù)估計方法,在時頻域通過逆Radon變換、Hough變換等方法將邊緣檢測問題轉(zhuǎn)換為峰值提取問題,此類方法估計精度較高,魯棒性強(qiáng),但是運(yùn)算量較大;二是基于圖像域的估計方法,通過對時頻圖進(jìn)行骨架提取的操作獲得目標(biāo)微動曲線,進(jìn)而實(shí)現(xiàn)對微動目標(biāo)的參數(shù)估計,該方法易受噪聲影響,隨著信噪比(signal to noise ratio,SNR)下降,估計精度將會越來越差。實(shí)測數(shù)據(jù)表明,閃爍現(xiàn)象在旋翼回波中普遍存在,上述方法在該條件下均難以有效估計出目標(biāo)的微動參數(shù),主要問題在于當(dāng)旋翼目標(biāo)微動回波出現(xiàn)閃爍時,回波信號幅度項(xiàng)增加了辛格(sinc(·))函數(shù)調(diào)制,時頻域中會出現(xiàn)余弦包絡(luò)、零頻帶和閃爍帶,雖然余弦包絡(luò)包含目標(biāo)的微動信息,但是易受其余微多普勒特征的影響,檢測難度很大。

近年來,深度學(xué)習(xí)在圖像識別、信號處理領(lǐng)域得到了廣泛運(yùn)用,為解決上述問題提供了一種新的思路。本文在研究閃爍現(xiàn)象下旋翼回波時頻圖特點(diǎn)的基礎(chǔ)上,提出一種有效的旋翼微動參數(shù)估計方法,針對零頻帶、閃爍帶以及大量無規(guī)律分布的噪點(diǎn)影響余弦包絡(luò)檢測的問題,將其都等效為時頻圖像中的噪聲,利用去噪卷積神經(jīng)網(wǎng)絡(luò)(denosing convolutional neural network,DnCNN)架構(gòu)分別訓(xùn)練去噪、去閃爍網(wǎng)絡(luò),消除噪聲、零頻帶和閃爍帶的影響。針對傳統(tǒng)逆Radon變換采用遍歷法進(jìn)行微動參數(shù)搜索,計算量較大的問題,采用黃金分割法對其搜索過程進(jìn)行改進(jìn),提高運(yùn)算速度,進(jìn)而完成對旋翼目標(biāo)微動參數(shù)的估計。

1 旋翼目標(biāo)信號建模及分析

不失一般性,假設(shè)雷達(dá)旋翼回波已經(jīng)經(jīng)過平動補(bǔ)償,建立以雷達(dá)所在位置為原點(diǎn)的OXYZ空間直角坐標(biāo)系,如圖1所示。該直升機(jī)擁有一個旋翼,多個葉片,以轉(zhuǎn)動角頻率正向旋轉(zhuǎn)。各個旋翼葉片長度均為,雷達(dá)到直升機(jī)旋翼中心(葉轂)的距離為,旋翼中心、雷達(dá)以及地面所形成的夾角為,方位角為。

圖1 雷達(dá)和旋翼葉片的幾何關(guān)系圖Fig.1 Geometric relationship between radar and rotor blades

設(shè)雷達(dá)發(fā)射信號為線性調(diào)頻信號:

(1)

散射點(diǎn)模型下直升機(jī)旋翼葉片可等效為多個散射點(diǎn),如圖2所示。

圖2 旋翼葉片散射點(diǎn)模型圖Fig.2 Scattering point model of rotor blades

旋翼葉片上散射點(diǎn)的回波信號進(jìn)行脈沖壓縮后,可表示為

(2)

式中:為葉片數(shù)目;為旋翼上的第個葉片;為單個葉片上的散射點(diǎn)數(shù)目;為葉片上的第個散射點(diǎn); 為對應(yīng)散射點(diǎn)的散射系數(shù);為觀測目標(biāo)的時間;sinc(·)表示辛格函數(shù);為線性調(diào)頻信號的帶寬; ()表示旋翼上對應(yīng)散射點(diǎn)到雷達(dá)的距離;c為光速; 是旋翼葉片對應(yīng)散射點(diǎn)的初相角;為發(fā)射信號的波長;′()表示脈壓后的噪聲。對應(yīng)圖1,令方位角=0,有：

()=+ cos( +)cos

(3)

當(dāng)旋翼葉片上各個散射點(diǎn)之間的距離 <2時,旋翼回波將出現(xiàn)閃爍現(xiàn)象,如圖3所示，單旋翼雙葉片散射點(diǎn)模型的閃爍現(xiàn)象仿真結(jié)果圖。圖3(a)表示=0、=0時葉片與雷達(dá)相對位置圖,圖3(b)為=03 m、=01 m、=120(每個葉片上的散射點(diǎn)數(shù)目)、=6 m、=12π rad/s、=0512 s、兩個葉片初始角度為0°和180°時通過短時傅里葉變換后的仿真結(jié)果。從圖3(b)可以看出,不加入噪聲時,時頻域中出現(xiàn)余弦包絡(luò)、零頻帶和閃爍帶,余弦包絡(luò)能量相對較弱。圖3(c)表示脈壓后SNR=10 dB時旋翼回波通過短時傅里葉變換后的仿真結(jié)果。從圖3(c)可以看出,當(dāng)回波中存在噪聲時,信號中的余弦包絡(luò)被噪點(diǎn)覆蓋,難以辨識。

圖3 雙葉片散射點(diǎn)模型閃爍現(xiàn)象示意圖Fig.3 Schematic diagram of flashing phenomenon of double blades’ scattering point model

余弦包絡(luò)是各個葉片葉尖散射點(diǎn)回波的時頻結(jié)果,其微多普勒頻率可表示為

(4)

(5)

2 旋翼目標(biāo)微動參數(shù)估計

圖4 微動參數(shù)估計流程圖Fig.4 Flow chart of micro-motion parameters estimation

2.1 基于DnCNN網(wǎng)絡(luò)的時頻圖預(yù)處理

211 網(wǎng)絡(luò)結(jié)構(gòu)

Zhang等提出一種用于圖像去噪的卷積神經(jīng)網(wǎng)絡(luò),稱為DnCNN神經(jīng)網(wǎng)絡(luò),該網(wǎng)絡(luò)包含17層,如圖5所示。輸入為帶噪灰度圖片=+,其中表示原始圖像,表示噪聲,第1層包含卷積層(Conv)和ReLu層,該卷積層包含64個卷積核且卷積核大小為3×3×1,ReLu代表一種非線性激活函數(shù)ReLu()=max(0,),用于避免因網(wǎng)絡(luò)層數(shù)增加造成的梯度消失等問題。第2層到第16層具有相同結(jié)構(gòu),包含卷積層(Conv)、批標(biāo)準(zhǔn)化(batch normalization,BN)層和ReLu層,卷積層包含64個卷積核且卷積核大小為3×3×64,在機(jī)器學(xué)習(xí)領(lǐng)域,有一個重要的獨(dú)立同分布(IID)假設(shè),即訓(xùn)練數(shù)據(jù)與測試數(shù)據(jù)滿足相同分布,BN層處于卷積層和ReLu層之間,可以使神經(jīng)網(wǎng)絡(luò)每一層的輸入都保持相同的分布特性,較好地保證訓(xùn)練出的模型在測試集上有一個很好的效果。第17層僅為一個卷積層,包含64個卷積核且卷積核大小為3×3×64。輸出為殘差圖,即中噪聲的分布圖片。整個網(wǎng)絡(luò)的損失函數(shù)定義為

圖5 DnCNN結(jié)構(gòu)Fig.5 DnCNN structure

(6)

DnCNN最大的特點(diǎn)就是引入了BN層和殘差學(xué)習(xí),BN層已經(jīng)在前文介紹,這里不再贅述,所謂的殘差學(xué)習(xí)如圖6所示。其思想是不再通過神經(jīng)網(wǎng)絡(luò)擬合潛在的映射(;)=(+;)=,而是通過神經(jīng)網(wǎng)絡(luò)擬合出(;)=(+;)=。依據(jù)文獻(xiàn)[22],直接擬合函數(shù)(;)的難度更大,而且隨著網(wǎng)絡(luò)層數(shù)的增加,不僅會出現(xiàn)梯度消失和梯度彌散的問題,網(wǎng)絡(luò)模型也會退化,而擬合函數(shù)(;)的難度相對而言更小,網(wǎng)絡(luò)模型更加魯棒。

圖6 神經(jīng)網(wǎng)絡(luò)映射方式Fig.6 Neural network mapping

2.1.2 網(wǎng)絡(luò)訓(xùn)練方式

去噪和去閃爍網(wǎng)絡(luò)通過串行方式進(jìn)行訓(xùn)練,如圖7所示。本文不直接利用DnCNN訓(xùn)練一個既可去噪又可去閃爍的網(wǎng)絡(luò),因?yàn)楦鶕?jù)文獻(xiàn)[22],輸入的訓(xùn)練對如果表現(xiàn)特征差距較大,那么網(wǎng)絡(luò)將不易訓(xùn)練。顯然,集合與集合的特征差距相比集合與集合的特征差距要小,所以分別訓(xùn)練兩種模型。

圖7 網(wǎng)絡(luò)訓(xùn)練方式示意圖Fig.7 Schematic diagram of network training

2.2 基于改進(jìn)逆Radon變換的微動參數(shù)估計方法

逆Radon變換是由平面中點(diǎn)的Radon變換推導(dǎo)出的,常用于微動目標(biāo)參數(shù)的估計,其數(shù)學(xué)原理如下。

如圖8所示平面內(nèi)任意一點(diǎn)(,)=(-)(-),記為,∠=,為該平面內(nèi)的一條直線,表達(dá)式為

圖8 xoy平面內(nèi)任意點(diǎn)示意圖Fig.8 Schematic diagram of any point in the xoy plane

=cos+sin

(7)

點(diǎn)沿該直線的Radon變換可表示為

(8)

(9)

式中：=cos；=sin。可見,經(jīng)過逆Radon變換,余弦曲線聚焦為一點(diǎn),該點(diǎn)坐標(biāo)為(,)。

(10)

(11)

(12)

(13)

將以(,)為圓心,半徑為10個單位內(nèi)的|(,)|置為0,計算|(,)|并確定其坐標(biāo)(,),由式(10)和式(11)估計出第二個葉片的微動參數(shù)。

將以(,)為圓心,半徑為10個單位內(nèi)的|(,)|置為0,計算|(,)+1|以及||(,)+1|-|(,)||。

(1) ||(,)|-|(,)+1||≤delta·|(,)|,那么求出|(,)+1|的坐標(biāo)(+1,+1),并依據(jù)式(10)和式(11)得出第+1個葉片的微動參數(shù),更新值為=+1,重復(fù)步驟3。

2.3 性能分析

利用黃金分割法的頻率估計方法具有估計精度高、迭代次數(shù)少、計算量小的優(yōu)點(diǎn)。下面對估計精度與迭代次數(shù)的關(guān)系進(jìn)行分析。

(14)

通過式(1)即可得到黃金分割法所需的迭代次數(shù)為

(15)

若利用遍歷搜索,對于轉(zhuǎn)動頻率的搜索區(qū)間[,]和搜索精度,需要的搜索次數(shù)為

(16)

若設(shè)置頻率搜索區(qū)間為[2π,14π],估計精度為0001時,利用式(15)可以得到黃金分割法搜索需要的迭代次數(shù)僅為22次。當(dāng)估計精度要求相同,通過遍歷搜索時,利用式(16)可以得到需要的搜索次數(shù)為37 699次?？梢钥闯?基于黃金分割法的頻率搜索方法在保證估計精度的同時有效減少了搜索次數(shù),運(yùn)算速度大大提高。

3 仿真結(jié)果與分析

3.1 雷達(dá)信號仿真參數(shù)設(shè)置及數(shù)據(jù)集準(zhǔn)備

訓(xùn)練集和測試集數(shù)據(jù)在64位Window10系統(tǒng),英特爾Core i7-8750H 2.2 GHz CPU,內(nèi)存8 GB設(shè)備上,通過Matlab 2019生成;神經(jīng)網(wǎng)絡(luò)的訓(xùn)練在同樣設(shè)備配置下,采用Pytorch框架,GPU版本實(shí)現(xiàn);旋翼目標(biāo)微動參數(shù)的估計在同樣設(shè)備配置下,通過Matlab 2019得到。

本文所采用的雷達(dá)信號是模擬信號,帶寬=1 MHz、采樣頻率=2 MHz、脈沖重復(fù)頻率PRF=8 000 Hz、脈沖寬度=100 μs、觀測時間=0.512 s、載頻=1 GHz、雷達(dá)反射截面積RCS=1、光速c=3×10m/s,考慮實(shí)際情況,設(shè)置脈壓后SNR為5～15 dB之間任意數(shù)。在散射點(diǎn)模型下生成旋翼目標(biāo)的雷達(dá)回波并作短時傅里葉變換獲得訓(xùn)練集、、各400張,測試集40張,每張圖片的像素值為180×180,圖片類型為灰度圖,圖片塊的像素值大小為40×40,每次投入網(wǎng)絡(luò)訓(xùn)練數(shù)為128個。為確保無重復(fù)參數(shù)作用,導(dǎo)致數(shù)據(jù)集失效,任意數(shù)的生成采用洗牌算法,具體仿真參數(shù)如表1和表2所示。

表1 微動參數(shù)

表2 短時傅里葉變換參數(shù)

3.2 神經(jīng)網(wǎng)絡(luò)訓(xùn)練步驟

采用正交初始化的方式設(shè)置初始權(quán)值,優(yōu)化器采用Adam優(yōu)化器。在訓(xùn)練過程中,前15輪學(xué)習(xí)率為0.000 1,后15輪學(xué)習(xí)率為0.000 01。為了使數(shù)據(jù)得到增強(qiáng),會隨機(jī)對圖片塊作翻轉(zhuǎn)、平移等操作。

3.3 仿真結(jié)果

基于深度學(xué)習(xí)的圖片預(yù)處理結(jié)果

圖9(a)、圖10(a)和圖11(a)為隨機(jī)選取的3張測試集圖片,其雷達(dá)直升機(jī)旋翼回波模型仿真參數(shù)如表3所示。圖9(b)、圖10(b)和圖11(b)分別為圖9(a)、圖10(a)和圖11(a)經(jīng)過去噪網(wǎng)絡(luò)以后的結(jié)果?？梢钥闯?去噪網(wǎng)絡(luò)有效去除了噪點(diǎn)的影響,恢復(fù)出了圖片中潛在的時頻特征,最外側(cè)余弦包絡(luò)也被很好地保留。圖9(c)、圖10(c)和圖11(c)分別為圖9(b)、圖10(b)和圖11(b)經(jīng)過去閃爍網(wǎng)絡(luò)以后的結(jié)果?？梢钥闯?閃爍帶和零頻帶被消除,而且余弦包絡(luò)能量被增強(qiáng)。

圖9 預(yù)處理結(jié)果(1)Fig.9 Preprocessing results (1)

圖10 預(yù)處理結(jié)果(2)Fig.10 Preprocessing results (2)

圖11 預(yù)處理結(jié)果(3)Fig.11 Preprocessing results (3)

表3 雷達(dá)直升機(jī)旋翼回波模型仿真參數(shù)

常規(guī)逆Radon變換的結(jié)果

圖12的3個分圖分別表示圖9(a)、圖10(a)和圖11(a)不做任何處理直接通過逆Radon變換得出的結(jié)果?？梢钥闯?無論旋翼目標(biāo)的微動參數(shù)如何變化,整個圖片僅會中心位置處出現(xiàn)一個聚焦點(diǎn),因此不對閃爍條件下的時頻結(jié)果進(jìn)行預(yù)處理,直接使用逆Radon變換是無法估計出旋翼的微動參數(shù)的。

圖12 預(yù)處理結(jié)果(4)Fig.12 Preprocessing results (4)

基于改進(jìn)逆Radon變換的參數(shù)估計結(jié)果

表4 旋翼長度估計結(jié)果

圖13 預(yù)處理結(jié)果(5)Fig.13 Preprocessing results (5)

表5 旋翼目標(biāo)轉(zhuǎn)動角頻率估計結(jié)果及誤差

表6 旋翼葉片初相估計結(jié)果

表7 旋翼葉片數(shù)估計結(jié)果

表8 不同算法運(yùn)算速度比較

以旋翼葉片長度的均值、轉(zhuǎn)動角頻率的誤差、各個葉片初相角的估計值、葉片數(shù)量的估計值作為度量標(biāo)準(zhǔn),從表4～表7的仿真結(jié)果可以看出,閃爍現(xiàn)象下,本文方法能夠很好地估計出旋翼目標(biāo)微動參數(shù),精度較高。對比表8和表9中改進(jìn)逆Radon變換與傳統(tǒng)逆Radon變換參數(shù)估計時間和兩種參數(shù)搜索方法的迭代次數(shù),也驗(yàn)證了改進(jìn)方法的時效性。

表9 不同算法迭代次數(shù)比較

4 結(jié) 論

旋翼目標(biāo)特征參數(shù)的提取對機(jī)型的穩(wěn)健識別意義重大。一般地,對閃爍現(xiàn)象下微動參數(shù)的提取方法研究不多,為此本文結(jié)合深度學(xué)習(xí)和變換域方法各自的優(yōu)勢,首先進(jìn)行時頻圖預(yù)處理,然后進(jìn)行參數(shù)估計,主要有以下結(jié)論:一是訓(xùn)練的去噪網(wǎng)絡(luò)和去閃爍網(wǎng)絡(luò)魯棒性較強(qiáng),從測試集結(jié)果來看,針對不同的時頻圖類型,都能夠顯著消除噪點(diǎn)、零頻帶和閃爍帶,為后續(xù)參數(shù)估計打下了基礎(chǔ);二是改進(jìn)逆Radon變換方法能在較高的精度要求下估計出旋翼目標(biāo)的微動參數(shù),時效性強(qiáng)。

實(shí)際情況中,目標(biāo)所處的噪聲環(huán)境十分復(fù)雜,DnCNN網(wǎng)絡(luò)的去噪能力和去閃爍能力會受到限制,不利于后續(xù)的參數(shù)估計。如何構(gòu)建一種具有去除實(shí)際噪聲的網(wǎng)絡(luò)模型是下一步研究的重點(diǎn)。