聯(lián)合時(shí)頻分析和譜估計(jì)的機(jī)動(dòng)目標(biāo)ISAR成像

周子鉑， 王鑫奎， 蔡萬勇， 劉建衛(wèi)， 張朝偉

(1． 空軍預(yù)警學(xué)院， 湖北武漢 430019; 2. 中國人民解放軍63696部隊(duì)， 江蘇江陰 214400)

0 引言

逆合成孔徑雷達(dá)(Inverse Synthetic Aperture Radar, ISAR)成像作為微波成像領(lǐng)域的代表性技術(shù)，通過積累回波可以獲得目標(biāo)的微波電磁圖像。具體而言，ISAR通過發(fā)射大帶寬的信號獲得距離向的高分辨，方位向的高分辨則由目標(biāo)與雷達(dá)之間的相對轉(zhuǎn)動(dòng)實(shí)現(xiàn)[1-3]。對于勻速轉(zhuǎn)動(dòng)目標(biāo)，傳統(tǒng)的距離-多普勒(Range-Doppler, RD)算法在完成平動(dòng)補(bǔ)償之后就能得到比較清晰的目標(biāo)圖像。然而對于機(jī)動(dòng)目標(biāo)而言，同一距離單元內(nèi)多個(gè)散射點(diǎn)回波的頻譜易發(fā)生混疊，致使傳統(tǒng)的RD成像方法獲得的ISAR圖像非常模糊。

此前，學(xué)者們已經(jīng)提出了很多方法來解決機(jī)動(dòng)目標(biāo)ISAR成像的方位向散焦問題。機(jī)動(dòng)目標(biāo)旋轉(zhuǎn)一般包含有加速度甚至加加速度，目標(biāo)的旋轉(zhuǎn)加速度產(chǎn)生二次相位項(xiàng)，高機(jī)動(dòng)目標(biāo)的旋轉(zhuǎn)加加速度將產(chǎn)生三次相位項(xiàng)。因此，散射點(diǎn)回波的多普勒頻率都是時(shí)變的[4-5]。文獻(xiàn)[6]提出一種利用參數(shù)自相關(guān)方程和廣義擴(kuò)展傅里葉變換的參數(shù)估計(jì)方法來重構(gòu)機(jī)動(dòng)目標(biāo)的距離-瞬時(shí)多普勒圖像。鑒于機(jī)動(dòng)目標(biāo)方位向回波為多分量線性調(diào)頻信號，自相關(guān)方程產(chǎn)生的交叉項(xiàng)致使ISAR圖像質(zhì)量惡化。由于目標(biāo)回波多普勒頻率的時(shí)變特性，能反映信號頻率隨時(shí)間變化的時(shí)頻分析工具通常被用來獲取高分辨的瞬時(shí)ISAR圖像，一些有效的時(shí)頻分析方法被提出。在文獻(xiàn)[7]中，回波信號被建模為多分量的三次相位信號，綜合利用距離-瞬時(shí)多普勒算法和距離-瞬時(shí)調(diào)頻率算法獲得高分辨的瞬時(shí)ISAR圖像。文獻(xiàn)[8]提出了基于距離-瞬時(shí)多普勒導(dǎo)數(shù)的機(jī)動(dòng)目標(biāo)ISAR成像算法，該算法通過多普勒頻率的任意階導(dǎo)數(shù)實(shí)現(xiàn)方位分辨，可視為傳統(tǒng)RID算法和RIC(Range-Instantaneous-Chirp-Rate)算法的擴(kuò)展。然而，上述算法中復(fù)雜的迭代參數(shù)估計(jì)步驟增加了算法負(fù)擔(dān)并且會(huì)導(dǎo)致較嚴(yán)重的誤差傳導(dǎo)。

對于弱機(jī)動(dòng)目標(biāo)，上述算法一般可以產(chǎn)生高分辨的ISAR圖像，但是對于強(qiáng)機(jī)動(dòng)目標(biāo)，受相干積累時(shí)間內(nèi)方位向高階相位的影響，上述方法成像分辨率會(huì)顯著惡化。在相干積累時(shí)間一定的條件下，超分辨成像可有效提高ISAR成像的分辨率。鑒于多重信號分類(Multiple Signal Classification, MUSIC)譜估計(jì)算法在波達(dá)方向(Direction of Arrival, DOA)應(yīng)用方面的高分辨特點(diǎn)，Odendaal等在文獻(xiàn)[9]中首先提出利用Schmidt提出的空間譜估計(jì)方法來重構(gòu)大轉(zhuǎn)角條件下的二維超分辨雷達(dá)圖像。南京電子技術(shù)研究所的凌牧等在文獻(xiàn)[10]利用基于AR-CAPON聯(lián)合譜估計(jì)的方法實(shí)現(xiàn)ISAR圖像的超分辨，該方法對實(shí)測ISAR數(shù)據(jù)建立自回歸(Auto Regressive, AR)模型，然后利用Burg算法求解AR模型進(jìn)行數(shù)據(jù)外推，提高CAPON濾波器的長度上限，再利用CAPON算法進(jìn)行超分辨ISAR成像。上述所述均為譜估計(jì)方法應(yīng)用于普通目標(biāo)，由于機(jī)動(dòng)目標(biāo)在成像時(shí)間內(nèi)的復(fù)雜運(yùn)動(dòng)，將譜估計(jì)方法應(yīng)用于經(jīng)由時(shí)頻工具處理之后的回波，不僅可有效避免目標(biāo)機(jī)動(dòng)對ISAR成像的影響，而且可克服小轉(zhuǎn)角對方位分辨率的限制，進(jìn)而獲得機(jī)動(dòng)目標(biāo)的高分辨ISAR圖像。

由于機(jī)動(dòng)目標(biāo)的方位向信號可被建模為多分量的線性調(diào)頻信號，本文提出利用短時(shí)傅里葉變換(Short Time Fourier Transform, STFT)獲得單個(gè)距離單元特定瞬間的頻譜，進(jìn)而減弱長成像時(shí)間內(nèi)多普勒頻率的時(shí)變特性。由于STFT的分辨率受限于“不確定性”原則，本文提出在對回波經(jīng)過STFT處理之后，再利用MUSIC譜估計(jì)方法重構(gòu)機(jī)動(dòng)目標(biāo)的高分辨ISAR圖像。本方法不僅可以有效避免其他時(shí)頻分析方法的交叉項(xiàng)影響，而且可以克服“不確定”原則對STFT分辨率的限制，進(jìn)而可獲得高分辨的ISAR圖像。

1 ISAR成像模型

機(jī)動(dòng)目標(biāo)的ISAR成像模型如圖1所示，XAY為雷達(dá)坐標(biāo)系，A為目標(biāo)的旋轉(zhuǎn)中心，即雷達(dá)所在的位置。R0表示目標(biāo)參考點(diǎn)O到雷達(dá)的距離;w,α分別表示目標(biāo)旋轉(zhuǎn)的初速度和加速度。θ表示目標(biāo)的瞬時(shí)旋轉(zhuǎn)角度：

(1)

此處，tm代表慢時(shí)間，目標(biāo)上的任一散射點(diǎn)Qi到雷達(dá)的距離可表示為

R=R0+x0sinθ+y0cosθ

(2)

圖1 ISAR成像模型

成像時(shí)間內(nèi)的目標(biāo)旋轉(zhuǎn)角度通常較小(3°～5°)，因此，有如下近似sinθ≈θ,cosθ≈1。進(jìn)而散射點(diǎn)到雷達(dá)的距離可以表示為

(3)

Qi散射點(diǎn)的多普勒頻率可以表示為

(4)

此處，λ表示雷達(dá)發(fā)射信號的波長。重構(gòu)Qi散射點(diǎn)的相位如下:

(5)

此處，φ0代表初始相位。P表示一個(gè)距離單元內(nèi)包含的散射點(diǎn)個(gè)數(shù)，在完成平動(dòng)補(bǔ)償之后，第l個(gè)距離單元的回波信號可以表示為多分量的線性調(diào)頻信號[11]:

(6)

根據(jù)公式(5)和(6)可知，機(jī)動(dòng)目標(biāo)的慢時(shí)間域回波可以表示為多分量的線性調(diào)頻信號，經(jīng)過傅里葉變換之后，各散射點(diǎn)回波的寬帶頻譜會(huì)產(chǎn)生重疊。因此，傳統(tǒng)RD算法獲取的機(jī)動(dòng)目標(biāo)ISAR圖像將是方位向散焦的。下面介紹本文提取的獲取機(jī)動(dòng)目標(biāo)的高分辨ISAR圖像的方法。

2 機(jī)動(dòng)目標(biāo)的高分辨ISAR成像

2.1 基于STFT方法的距離-瞬時(shí)多普勒成像

STFT通過截取信號的一小段并應(yīng)用傅里葉變換獲得此一小段信號的頻譜，具體實(shí)現(xiàn)如圖2所示。

圖2 STFT截取信號

(7)

(8)

(9)

(10)

2.2 基于STFT方法的距離-瞬時(shí)多普勒成像

鑒于“不確定”原則對STFT時(shí)間分辨率和頻率分辨率的限制，在窗函數(shù)截取信號使信號長度變短之后，頻率分辨率變差，進(jìn)而ISAR成像結(jié)果變得模糊。因此，此部分主要介紹空間譜估計(jì)方法對STFT結(jié)果的超分辨處理，克服“不確定”原則的限制。

如式(8)和(9)所示，在截取窗的短時(shí)間內(nèi)，時(shí)間tm可以近似為常數(shù)，每個(gè)散射點(diǎn)的多普勒頻率僅僅與散射點(diǎn)的橫向坐標(biāo)xi相關(guān)。可得第l個(gè)距離單元的回波可表示為

(11)

已知高分辨空間譜估計(jì)方法常用于波達(dá)方向(Direction of Arrival, DOA)估計(jì)，MUSIC估計(jì)方法可以估計(jì)得到多個(gè)信源的準(zhǔn)確方位，類比可以將MUSIC估計(jì)方法應(yīng)用于單個(gè)距離單元內(nèi)散射點(diǎn)橫向位置估計(jì)，即多普勒頻率域的超分辨。首先，由式(11)可將第l個(gè)距離單元的回波以矩陣形式表示。由于慢時(shí)間維選取了N個(gè)脈沖用來成像，因此第l個(gè)距離單元的維度為N×1，如式(12)所示。

(12)

X=A·S+N

(13)

(14)

(15)

(16)

注意到，MUSIC譜估計(jì)方法用到了信號的統(tǒng)計(jì)特性，這需要多次數(shù)據(jù)采集來實(shí)現(xiàn)。不同于陣列天線可以在不同的時(shí)間獲取多次快拍以獲得信號的統(tǒng)計(jì)特性，ISAR成像僅有一次快拍，所以經(jīng)典的MUSIC譜估計(jì)方法不能直接應(yīng)用在此處。滑窗技術(shù)被用于解決此問題，以獲得多次快拍，并消除回波之間的相干性，具體如圖3所示。其中，滑窗的長度為(d×1,d

(17)

第k個(gè)子矢量Xk的協(xié)方差矩陣可寫作

(18)

此處H代表復(fù)共軛轉(zhuǎn)置。完整回波矢量的協(xié)方差矩陣為

(19)

(20)

圖3 回波滑窗技術(shù)

MUSIC譜估計(jì)算法實(shí)現(xiàn)高分辨，是基于數(shù)據(jù)自相關(guān)矩陣信號特征子空間和噪聲子空間的正交關(guān)系。協(xié)方差矩陣RXX可以分解得到N個(gè)特征值，從大到小依次為(λ1>λ2>λ3>…>λN)，對應(yīng)特征值的特征向量設(shè)為v。信號子空間由P個(gè)大的特征值對應(yīng)的P個(gè)特征向量張成Vs=[v1,v2,v3,…,vP]，噪聲子空間則由剩余的N-P個(gè)小特征值對應(yīng)的特征向量張成Vn=[vP+1,vP+2,vP+3,…,vN]。因此，協(xié)方差矩陣RXX的分解如式(21)所示：

(21)

此處∑表示P個(gè)大特征值組成的信號特征值，因此∑為(P×P)維的對角矩陣。σ表示N-P個(gè)小特征值，噪聲特征值可近似視為相等。最后，MUSIC的譜估計(jì)結(jié)果可由式(22)求解得到。在特定的范圍內(nèi)搜索散射點(diǎn)的橫向坐標(biāo)x，峰值出現(xiàn)的位置即為散射點(diǎn)所在的位置，相比于經(jīng)典的傅里葉變換而言，基于MUSIC算法的空間譜估計(jì)算法能獲得更好的頻率分辨率。

(22)

2.3 具體成像步驟

步驟1 對所有回波進(jìn)行距離壓縮和平動(dòng)補(bǔ)償，得到方位向回波為

(23)

步驟2 對所有距離單元的方位向回波S(f,tm)進(jìn)行STFT處理，并通過最小熵準(zhǔn)則來選擇聚焦效果最佳的多普勒頻譜S(f,fm)。

上述為基于聯(lián)合時(shí)頻分析和MUSIC超分辨的機(jī)動(dòng)目標(biāo)ISAR成像算法，具體的算法執(zhí)行流程如圖4所示。

圖4 算法執(zhí)行步驟

3 ISAR成像實(shí)驗(yàn)與分析

本部分通過ISAR成像實(shí)驗(yàn)驗(yàn)證上述所提出方法的有效性，包含兩組實(shí)驗(yàn)：散射點(diǎn)仿真模型和波音B727實(shí)測數(shù)據(jù)成像實(shí)驗(yàn)。雷達(dá)參數(shù)和仿真目標(biāo)參數(shù)如表1所示。所有實(shí)驗(yàn)的成像模型如圖1所示，并假設(shè)平動(dòng)補(bǔ)償已經(jīng)完成。仿真散射點(diǎn)模型如圖5所示。

表1 目標(biāo)和雷達(dá)的仿真參數(shù)

圖5 仿真模型

為形象準(zhǔn)確分析目標(biāo)機(jī)動(dòng)對于回波多普勒頻率的影響，此處先對第132個(gè)距離單元進(jìn)行時(shí)頻分析，分別采用STFT和魏格納-威利分布(WVD)對其進(jìn)行處理，具體兩種算法獲得的時(shí)頻分布如圖6所示。圖6(a)表示STFT算法獲得的時(shí)頻分布圖，圖中顯示包含有4條亮線，即表示包含有4個(gè)散射點(diǎn)，但是4條線在頻率軸方向都比較粗，表明STFT獲得的時(shí)頻分布圖的頻率分辨率比較低，即為前面講到的STFT分辨率受“不確定”原則限制。圖6(b)表示W(wǎng)VD算法獲得的時(shí)頻分布圖，根據(jù)WVD的定義可知WVD分布圖包含很多交叉項(xiàng)，并且圖中無法區(qū)分出明確的4條線，但是最外圍的兩條自項(xiàng)產(chǎn)生的時(shí)頻分布線比STFT結(jié)果較細(xì)，即WVD可以獲得較高的頻率分辨率。然而，兩個(gè)分布圖都表示多普勒頻率隨慢時(shí)間逐漸變化，表明機(jī)動(dòng)目標(biāo)回波多普勒頻率的時(shí)變特性。

(a) STFT計(jì)算結(jié)果

(b) WVD計(jì)算結(jié)果圖6 第132個(gè)距離單元的時(shí)頻分析結(jié)果

3.1 散射點(diǎn)仿真模型實(shí)驗(yàn)

仿真散射點(diǎn)模型的機(jī)動(dòng)目標(biāo)成像實(shí)驗(yàn)的回波模型如圖1所示，此時(shí)目標(biāo)的旋轉(zhuǎn)參數(shù)分別為：初始角速度為0.04 rad/s，角加速度為0.03 rad/s2。圖7顯示了4種方法的仿真模型機(jī)動(dòng)目標(biāo)成像結(jié)果。圖7(a)為傳統(tǒng)RD方法獲得的機(jī)動(dòng)目標(biāo)ISAR成像結(jié)果，從圖中可以看到由于目標(biāo)的機(jī)動(dòng)飛行，所成ISAR圖像在方位維會(huì)產(chǎn)生嚴(yán)重的散焦現(xiàn)象。圖7(b)和圖7(c)為時(shí)頻分析工具WVD和STFT算法所成的ISAR圖像，在不同的窗長和時(shí)間點(diǎn)上對回波信號進(jìn)行WVD或者STFT處理，然后根據(jù)最小熵準(zhǔn)則選擇聚焦效果最佳的一維方位像，最終組成所得到的ISAR圖像。觀察圖像可見，一方面，由于WVD中眾多交叉項(xiàng)的存在，在所成的ISAR圖像中會(huì)產(chǎn)生眾多對應(yīng)的“假點(diǎn)”；另一方面，雖然STFT算法所獲得的圖像沒有交叉項(xiàng)產(chǎn)生的“假點(diǎn)”，但是由于“不確定”原則的限制，STFT算法的方位向分辨率較差，最終導(dǎo)致所成ISAR圖像也比較模糊，如圖中白色圓圈部分所示。而且，根據(jù)圖像可進(jìn)一步觀察到，散射點(diǎn)越遠(yuǎn)離目標(biāo)的方位軸心，散射點(diǎn)在方位維的模糊現(xiàn)象越嚴(yán)重。圖7(d)所示為本文所提新方法獲得的ISAR圖像，相比于上述三種成像方法，新方法獲得的圖像更加清晰，具體而言，白色圓圈內(nèi)的模糊現(xiàn)象得到很好的解決。

(a) 經(jīng)典RD算法成像結(jié)果

(b) WVD算法成像結(jié)果

(c) STFT算法成像結(jié)果

(d) 所提出方法成像結(jié)果圖7 仿真模型ISAR成像結(jié)果

為進(jìn)一步表明所提出的聯(lián)合時(shí)頻分析和譜估計(jì)成像方法可獲得更好的聚焦效果，選取上述圖7中的第132個(gè)距離單元，其歸一化曲線圖如圖8所示，觀察圖像可知，此距離單元具體包含有4個(gè)散射點(diǎn)。圖8(a)和圖8(b)分別表示FFT和WVD的多普勒頻率聚焦效果。因?yàn)镕FT呈現(xiàn)的是回波的整個(gè)的多普勒頻譜，由于目標(biāo)的機(jī)動(dòng)操作，各散射點(diǎn)的多普勒頻率為時(shí)變的，即為寬帶頻譜，多個(gè)散射點(diǎn)的多個(gè)頻譜發(fā)生混疊，最終多普勒頻率聚焦效果被嚴(yán)重影響，導(dǎo)致4個(gè)散射點(diǎn)無法被完全分辨。另一方面，WVD的交叉項(xiàng)導(dǎo)致原本4個(gè)散射點(diǎn)的4個(gè)譜峰連接在一起，并且多個(gè)交叉項(xiàng)的疊加會(huì)產(chǎn)生偽峰，從而導(dǎo)致各散射點(diǎn)無法被完全區(qū)分，甚至產(chǎn)生“假點(diǎn)”。圖8(c)和圖8(d)呈現(xiàn)STFT算法和所提出算法的更好的散射點(diǎn)分辨效果。為具體說明MUSIC譜估計(jì)方法對于STFT算法的估計(jì)效果，圖8(c)和圖8(d)的第二個(gè)峰值(圖8(c)中X=121，圖8(d)中X=58)的寬度數(shù)據(jù)如表2所示。由表2數(shù)據(jù)可知，所提出方法的峰值下降速度比STFT算法的峰值下降速度更快，可以說明圖8(d)中的峰值更窄，即所提出方法可以獲得更好的方位分辨率。

(a) FFT算法聚焦效果

(b) WVD算法聚焦效果

(c) STFT算法聚焦效果

(d) 所提出方法聚焦效果圖8 第132個(gè)距離單元的方位向聚焦效果

表2 仿真數(shù)據(jù)第二個(gè)譜峰的寬度

3.2 波音B727實(shí)測數(shù)據(jù)實(shí)驗(yàn)

此部分，美國海軍實(shí)驗(yàn)室的波音B727實(shí)測數(shù)據(jù)[13]被用來進(jìn)一步實(shí)驗(yàn)驗(yàn)證所提成像方法的有效性。其中，回波的載頻和帶寬分別為9 GHz和150 MHz，脈沖重復(fù)頻率為20 kHz，此處選擇256個(gè)脈沖用于ISAR成像。目標(biāo)回波平動(dòng)補(bǔ)償之后的一維距離像如圖9所示。

圖9 波音B727實(shí)測數(shù)據(jù)一維距離像

圖10給出了4種方法對于波音B727飛機(jī)的實(shí)測數(shù)據(jù)成像結(jié)果。同樣，圖10(a)為傳統(tǒng)RD算法獲得的ISAR圖像，顯然由于目標(biāo)的機(jī)動(dòng)飛行，ISAR圖像在方位向出現(xiàn)嚴(yán)重的散焦現(xiàn)象。然后，為克服散焦現(xiàn)象，WVD和STFT算法分別用來實(shí)現(xiàn)機(jī)動(dòng)目標(biāo)ISAR成像，如圖10(b)和圖10(c)所示。與上述散射點(diǎn)模型實(shí)驗(yàn)類似，WVD算法中生成的交叉項(xiàng)會(huì)產(chǎn)生眾多“假點(diǎn)”，從而導(dǎo)致所成圖像幾乎無法識別出目標(biāo)的輪廓。相比于WVD算法，STFT算法可以獲得目標(biāo)的圖像，但是，鑒于“不確定”原則對分辨率的限制，且STFT中的滑窗降低了信號的長度，導(dǎo)致STFT所成圖像在方位維分辨率較低，即各散射點(diǎn)依舊會(huì)在方位維產(chǎn)生模糊，如圖10(c)中呈現(xiàn)的結(jié)果所示。由于MUSIC空間譜估計(jì)方法可以獲得更好的頻譜分辨率，所以MUSIC對STFT算法處理之后的回波進(jìn)行譜估計(jì)可以克服“不確定”原則對頻率分辨率的限制，獲得更佳的ISAR成像結(jié)果，本文所提新方法的成像結(jié)果如圖10(d)所示，并且所得圖像比上述三種方法更佳清晰。

(a) 經(jīng)典RD算法成像結(jié)果

(b) WVD算法成像結(jié)果

(c) STFT算法成像結(jié)果

(d) 所提方法成像結(jié)果圖10 波音B727實(shí)測數(shù)據(jù)成像結(jié)果

同樣，為進(jìn)一步表明多提出方法對頻率分辨率的改善，圖10中的第31個(gè)距離單元的聚焦效果被選取出顯示于圖11。圖11(a)為STFT算法對方位向的聚焦效果，圖11(b)為本文所提方法對方位向聚焦的效果，由圖可見，該距離單元主要包含兩個(gè)散射點(diǎn)，此處選擇第一個(gè)峰值(圖11(a)中X=119，圖11(b)中X=56)用比較頻譜寬度，具體譜峰寬度數(shù)據(jù)如表3所示。可見，所提成像方法的峰值的下降速度比STFT算法的峰值下降速度快，即圖11(b)中的峰值寬度更窄，即所提方法能獲得更好的方位向聚焦效果。

上述實(shí)驗(yàn)結(jié)果已經(jīng)證實(shí)基于聯(lián)合時(shí)頻分析和譜估計(jì)的成像方法的有效性，更進(jìn)一步的機(jī)動(dòng)目標(biāo)ISAR超分辨成像算法的研究將在未來的工作呈現(xiàn)。

(a) STFT算法聚焦效果

(b) 所提方法聚焦效果圖11 波音B727實(shí)測數(shù)據(jù)第31個(gè)距離單元的聚焦效果

表3 實(shí)測數(shù)據(jù)第一個(gè)譜峰的寬度

4 結(jié)束語

本文介紹了一種聯(lián)合時(shí)頻分析和譜估計(jì)的機(jī)動(dòng)目標(biāo)ISAR成像算法，該方法在整個(gè)成像時(shí)間內(nèi)選取一段短時(shí)間對目標(biāo)進(jìn)行方位向壓縮，進(jìn)而克服長時(shí)間內(nèi)的多普勒頻率的時(shí)變特性，然后應(yīng)用MUSIC譜估計(jì)方法對上述方位向壓縮結(jié)果進(jìn)行超分辨處理，最終在該段短時(shí)間內(nèi)實(shí)現(xiàn)機(jī)動(dòng)目標(biāo)的高分辨成像。相比于單一的時(shí)頻分析方法，所提方法可以克服“不確定”原則對成像分辨率的限制，獲得高分辨成像結(jié)果。最后，通過仿真模型和波音B727實(shí)測數(shù)據(jù)證實(shí)了所提方法的有效性。