陳是扦,彭志科,邢冠培
(1.上海交通大學(xué) 機(jī)械系統(tǒng)與振動(dòng)國(guó)家重點(diǎn)實(shí)驗(yàn)室,上海 200240; 2. 上海航天電子技術(shù)研究所,上海 201109)
雷達(dá)通過(guò)發(fā)射、接收目標(biāo)反射的電磁波進(jìn)行探測(cè)。當(dāng)目標(biāo)以一定速度運(yùn)動(dòng)時(shí),其反射信號(hào)將偏離雷達(dá)發(fā)射信號(hào)的頻率,該現(xiàn)象稱為多普勒效應(yīng)。如果目標(biāo)除了整體平移運(yùn)動(dòng)外,還進(jìn)行局部振動(dòng)、轉(zhuǎn)動(dòng)等微運(yùn)動(dòng),則會(huì)使雷達(dá)回波信號(hào)產(chǎn)生頻率調(diào)制,并在目標(biāo)整體平移運(yùn)動(dòng)引起的多普勒頻率周圍產(chǎn)生邊帶,這種現(xiàn)象稱為微多普勒效應(yīng)[1-3]。微多普勒調(diào)制頻率可反映目標(biāo)的微動(dòng)特性和結(jié)構(gòu)特征,被廣泛應(yīng)用于雷達(dá)目標(biāo)分類識(shí)別中[4-5]。旋轉(zhuǎn)運(yùn)動(dòng)是航天領(lǐng)域中最為常見(jiàn)的微運(yùn)動(dòng),如衛(wèi)星天線轉(zhuǎn)動(dòng)、彈道中段彈頭自旋運(yùn)動(dòng)等,因此本文將主要關(guān)注旋轉(zhuǎn)目標(biāo)的微多普勒頻率。
目標(biāo)運(yùn)動(dòng)使微多普勒信號(hào)頻率具有時(shí)變特征[6]。目前,國(guó)內(nèi)外已有大量利用時(shí)頻分析方法處理這類非平穩(wěn)信號(hào)的報(bào)道,如:文獻(xiàn)[7]利用經(jīng)驗(yàn)?zāi)J椒纸夥椒ǚ蛛x出目標(biāo)旋轉(zhuǎn)部件產(chǎn)生的微多普勒特征分量,以消除微多普勒特征對(duì)目標(biāo)主體部分ISAR(逆合成孔徑雷達(dá))成像的影響,但經(jīng)驗(yàn)?zāi)J椒纸夥椒ㄈ菀桩a(chǎn)生邊界效應(yīng)、模式混疊等問(wèn)題,難以處理航天領(lǐng)域各種復(fù)雜的雷達(dá)信號(hào);文獻(xiàn)[8-9]分別利用Fourier-Bessel級(jí)數(shù)和Chirplet變換提取雷達(dá)回波信號(hào)中的微多普勒特征,利用Fourier-Bessel級(jí)數(shù)分離雷達(dá)信號(hào)中的非平穩(wěn)成分與平穩(wěn)成分,而非平穩(wěn)成分往往與雷達(dá)微多普勒特征相對(duì)應(yīng),不能進(jìn)一步分離不同散射點(diǎn)的微多普勒特征。上述文獻(xiàn)主要介紹了微多普勒信號(hào)分解,未見(jiàn)信號(hào)頻率特征提取的報(bào)道。文獻(xiàn)[10-11]通過(guò)檢測(cè)回波信號(hào)時(shí)頻分布的局部峰值點(diǎn),估計(jì)目標(biāo)的微多普勒頻率,但由于目標(biāo)不同散射點(diǎn)的微多普勒頻率分量通常在時(shí)頻分布上相互重疊,該方法難以得到感興趣的頻率分量,且信號(hào)時(shí)頻分布的分辨率有限,因此較難估計(jì)高精度的微多普勒頻率;文獻(xiàn)[12]提出一種參數(shù)化解調(diào)方法,用于估計(jì)多項(xiàng)式頻率調(diào)制信號(hào)的瞬時(shí)頻率,而該方法能有效處理瞬時(shí)頻率交叉、重疊的信號(hào),但其采用的多項(xiàng)式模型在逼近復(fù)雜瞬時(shí)頻率曲線時(shí)存在龍格現(xiàn)象,無(wú)法估計(jì)高度振蕩的微多普勒頻率。
鑒于旋轉(zhuǎn)目標(biāo)微多普勒信號(hào)具有正弦頻率調(diào)制特征,本文提出了一種基于正弦模型的參數(shù)化解調(diào)方法。該方法與文獻(xiàn)[12]中的多項(xiàng)式模型解調(diào)方法不同,將頻率估計(jì)問(wèn)題轉(zhuǎn)化為模型參數(shù)優(yōu)化問(wèn)題,利用粒子群算法優(yōu)化參數(shù),對(duì)旋轉(zhuǎn)目標(biāo)呈周期振蕩的微多普勒頻率進(jìn)行有效估計(jì)。
將旋轉(zhuǎn)目標(biāo)簡(jiǎn)化為繞固定旋轉(zhuǎn)中心運(yùn)動(dòng)的散射點(diǎn)[13],雷達(dá)與旋轉(zhuǎn)散射點(diǎn)的幾何關(guān)系如圖1所示。圖中:散射點(diǎn)P繞回轉(zhuǎn)中心O以角頻率fr逆時(shí)針旋轉(zhuǎn);l0為雷達(dá)到回轉(zhuǎn)中心的距離;l1為回轉(zhuǎn)半徑;LP為雷達(dá)到散射點(diǎn)的瞬時(shí)距離;φ0、φt分別為起始時(shí)刻和t時(shí)刻散射點(diǎn)P的旋轉(zhuǎn)角,φt=φ0+2πfrt。
圖1 雷達(dá)與旋轉(zhuǎn)散射點(diǎn)的幾何關(guān)系Fig.1 Geometrical relationship between radar and rotating scattering point
由圖1中的幾何關(guān)系得到雷達(dá)與散射點(diǎn)之間距離,即
(1)
假設(shè)散射體位于雷達(dá)遠(yuǎn)場(chǎng),即(l1/l0)2→0,則式(1)可簡(jiǎn)化為
LP(t)≈l0+l1cos(φ0+2πfrt)
(2)
因此,點(diǎn)P的雷達(dá)回波基帶信號(hào)可表示為
(3)
(4)
(5)
式中:K為散射點(diǎn)數(shù)目;σi為第i個(gè)散射點(diǎn)的散射系數(shù);Li(t)為雷達(dá)與第i個(gè)散射點(diǎn)之間的瞬時(shí)距離。
假設(shè)所有微動(dòng)散射點(diǎn)的旋轉(zhuǎn)頻率為fr,則參照式(2)可得到瞬時(shí)距離函數(shù),其表達(dá)式為
Li(t)≈l0+li1cos(φi0+2πfrt)
(6)
式中:li1為第i個(gè)散射點(diǎn)的旋轉(zhuǎn)半徑;φi0為第i個(gè)散射點(diǎn)的起始旋轉(zhuǎn)角。參照式(4)可得到第i個(gè)散射點(diǎn)的微多普勒頻率,即
(7)
由式(7)可得:不同散射點(diǎn)的微多普勒頻率是一些幅值和相位不同、角頻率相同的正弦函數(shù),這些微多普勒頻率通常在時(shí)頻面上相互重疊,難以提取。
參數(shù)化解調(diào)方法可估計(jì)特定信號(hào)的模型參數(shù)。廣義調(diào)頻信號(hào)模型為
(8)
式中:A為信號(hào)幅值;ψ0為信號(hào)初始相位;fI(t)為信號(hào)瞬時(shí)頻率。若式(8)為雷達(dá)回波信號(hào),則fI(t)為微多普勒頻率,其表達(dá)式見(jiàn)式(4)。
參數(shù)化解調(diào)算子的表達(dá)式為
Φ(t;P,fc)=
(9)
式中:fc為載波頻率;g(t;P)為特定的參數(shù)化解調(diào)函數(shù)[12];P為解調(diào)函數(shù)的參數(shù)集合。將式(8)乘以式(9)得到解調(diào)信號(hào),即
(10)
當(dāng)且僅當(dāng)fI(t)=g(t;P)時(shí),式(10)為單頻信號(hào)sd(t)=Aexp[j(2πfct+ψ0)],即信號(hào)能量集中在頻率fc處,信號(hào)頻譜值在fc處達(dá)到最大。因此,可通過(guò)優(yōu)化信號(hào)模型參數(shù)P使sd(t)在載波頻率fc處的頻譜值達(dá)到最大,即
(11)
需要注意的是,原始信號(hào)由解調(diào)信號(hào)調(diào)制得到,即
s(t)=sd(t)×conj[Φ(t;P,fc)]=
(12)
式中:conj為復(fù)數(shù)共軛。
由于旋轉(zhuǎn)目標(biāo)的微多普勒頻率隨時(shí)間正弦變化,因此提出基于正弦模型的參數(shù)化解調(diào)方法,用于提取目標(biāo)微多普勒頻率。將式(7)整理為
(13)
式中:ai為多普勒頻率幅值,ai=4πfrli1/λ;φi0=φi0+π。參照式(9),構(gòu)造與式(13)相匹配的參數(shù)化解調(diào)算子,即
(14)
式中:fc為載波頻率,fc=0(參考式(9));Pi={ai,fr,φi0}為待估計(jì)的模型參數(shù)集。參照式(11)可估計(jì)微多普勒頻率參數(shù),即
(15)
式中:sr(t)為雷達(dá)回波信號(hào)(見(jiàn)式(5))。采用粒子群優(yōu)化算法對(duì)式(15)進(jìn)行優(yōu)化求解[14]。
回波信號(hào)包含K個(gè)信號(hào)分量,由于強(qiáng)信號(hào)分量對(duì)弱信號(hào)分量的估計(jì)有抑制作用,通常只能估計(jì)到能量最強(qiáng)(即σi最大)的信號(hào)分量參數(shù)。為得到所有K個(gè)散射點(diǎn)的微多普勒頻率參數(shù),采用迭代估計(jì)算法,將估計(jì)到的強(qiáng)信號(hào)分量從原始信號(hào)中剔除(消除強(qiáng)信號(hào)分量的干擾),從剩余信號(hào)中估計(jì)其他的弱信號(hào)分量參數(shù)。
算法流程具體步驟如下:
1) 令i=1;
7) 重復(fù)執(zhí)行步驟2)~6),直到i≥K。
上述步驟中,步驟6)將當(dāng)前提取到的強(qiáng)信號(hào)分量從回波信號(hào)中移除,因此該算法能在下次迭代中從剩余信號(hào)估計(jì)到其他能量較弱的信號(hào)分量參數(shù)。上述算法經(jīng)過(guò)K次迭代可提取所有散射點(diǎn)的微多普勒頻率參數(shù)。算法流程如圖2所示。
圖2 算法流程Fig.2 Flowchart of algorithm
雷達(dá)波長(zhǎng)λ=0.06 m;雷達(dá)與旋轉(zhuǎn)目標(biāo)中心的距離l0=500 m;目標(biāo)包含3個(gè)強(qiáng)散射點(diǎn),其旋轉(zhuǎn)半徑分別為l11=3 m,l21=2.5 m,l31=1 m;初始旋轉(zhuǎn)角分別為φ10=π/3 rad,φ20=5π/6 rad,φ30=4π/3 rad;散射系數(shù)σ1=3,σ2=2.8,σ3=2;目標(biāo)旋轉(zhuǎn)角頻率fr=2 Hz。假設(shè)信號(hào)被高斯白噪聲污染,信噪比為0 dB,則利用短時(shí)傅里葉變換得到含噪聲的雷達(dá)回波信號(hào)的時(shí)頻分布,仿真信號(hào)時(shí)頻分布如圖3所示。圖中可見(jiàn),3個(gè)散射點(diǎn)的微多普勒頻率曲線在時(shí)頻面上相互相交、重疊。
圖3 仿真信號(hào)時(shí)頻分布Fig.3 Time-frequency distribution of simulated signal
采用本文方法估計(jì)式(13)中的微多普勒模型參數(shù),包括微多普勒振幅ai、振蕩頻率fr(目標(biāo)旋轉(zhuǎn)角頻率)、初始相位φi0。進(jìn)一步估計(jì)散射點(diǎn)的旋轉(zhuǎn)半徑li1,li1=λai/(4πfr)。3個(gè)散射點(diǎn)的參數(shù)估計(jì)結(jié)果見(jiàn)表1。結(jié)果表明:本文方法在強(qiáng)噪聲干擾下能得到高精度的參數(shù)估計(jì)結(jié)果。
表1 仿真信號(hào)微多普勒參數(shù)估計(jì)結(jié)果
圖4 仿真信號(hào)微多普勒頻率提取結(jié)果對(duì)比Fig.4 Comparison of extracted micro-Doppler frequencies of simulated signal
分別用文獻(xiàn)[11]中的方法和本文方法進(jìn)行微多普勒頻率提取,對(duì)比結(jié)果如圖4所示。文獻(xiàn)[11]通過(guò)檢測(cè)回波信號(hào)時(shí)頻分布(見(jiàn)圖3)的局部峰值點(diǎn)提取微多普勒頻率曲線。由圖4可見(jiàn):不同曲線在交點(diǎn)處模糊,局部峰值檢測(cè)法在通過(guò)微多普勒頻率交叉點(diǎn)后容易識(shí)別到錯(cuò)誤信號(hào)分量,而本文方法能正確識(shí)別交叉信號(hào)分量。
為驗(yàn)證本文方法在實(shí)際應(yīng)用中的有效性,對(duì)雷達(dá)實(shí)測(cè)數(shù)據(jù)的微多普勒頻率進(jìn)行估計(jì)。實(shí)驗(yàn)裝置與文獻(xiàn)[15]中的類似,如圖5所示。圖中:2個(gè)金屬圓筒以180°的間隔布置在轉(zhuǎn)盤上,轉(zhuǎn)盤以一定速度旋轉(zhuǎn)運(yùn)動(dòng);將2個(gè)金屬圓筒視為2個(gè)旋轉(zhuǎn)目標(biāo),并用發(fā)射頻率為10.5 GHz(即波長(zhǎng)為0.028 5 m)的連續(xù)波雷達(dá)探測(cè)2個(gè)目標(biāo)的旋轉(zhuǎn)運(yùn)動(dòng);2個(gè)圓筒距離旋轉(zhuǎn)中心的距離都為0.1 m。
實(shí)測(cè)信號(hào)的時(shí)頻分布如圖6所示。由圖可見(jiàn):2個(gè)旋轉(zhuǎn)目標(biāo)的微多普勒頻率曲線隨時(shí)間正弦變化,兩曲線相位差約為180°。
用本文方法估計(jì)實(shí)驗(yàn)信號(hào)的微多普勒頻率參數(shù),結(jié)果見(jiàn)表2。實(shí)驗(yàn)中轉(zhuǎn)盤的實(shí)際轉(zhuǎn)速未知,無(wú)法提供各頻率參數(shù)的真實(shí)值,但由估計(jì)到的頻率參數(shù)計(jì)算得到的目標(biāo)旋轉(zhuǎn)半徑li1與實(shí)際半徑(0.1 m)的誤差不超過(guò)2%。實(shí)驗(yàn)結(jié)果表明:本文方法對(duì)于實(shí)測(cè)雷達(dá)數(shù)據(jù)仍具有較高的頻率估計(jì)精度。
本文方法與傳統(tǒng)時(shí)頻峰值檢測(cè)法得到的微多普勒頻率曲線對(duì)比如圖7所示。由圖可見(jiàn):傳統(tǒng)時(shí)頻峰值檢測(cè)法無(wú)法正確識(shí)別交叉的信號(hào)分量,且由于其受其他信號(hào)分量的干擾,提取到的頻率曲線出現(xiàn)了局部波動(dòng)和毛刺;本文方法采用與信號(hào)匹配的參數(shù)化模型,抗干擾能力強(qiáng),能得到精確的微多普勒頻率曲線。
本文對(duì)旋轉(zhuǎn)目標(biāo)的微多普勒頻率提取問(wèn)題進(jìn)行了研究,分析了基于等效散射點(diǎn)模型的旋轉(zhuǎn)目標(biāo)雷達(dá)回波信號(hào)解析形式。針對(duì)信號(hào)模型特征,提出了基于正弦模型的參數(shù)化解調(diào)方法,以提取旋轉(zhuǎn)目標(biāo)的微多普勒頻率。在算法具體實(shí)施中,采用迭代頻率提取方法,逐個(gè)獲取每個(gè)散射點(diǎn)的微多普勒頻率。在仿真實(shí)驗(yàn)中,將本文方法與傳統(tǒng)時(shí)頻峰值檢測(cè)法進(jìn)行對(duì)比。結(jié)果表明:本文方法能有效提取相互交叉的微多普勒頻率曲線,算法精度高,抗干擾能力強(qiáng),得到的微多普勒頻率可進(jìn)一步用于目標(biāo)運(yùn)動(dòng)參數(shù)識(shí)別提取,在識(shí)別、監(jiān)測(cè)彈道導(dǎo)彈飛行參數(shù)等方面具有重要應(yīng)用前景。后續(xù)將基于該方法開(kāi)展航天領(lǐng)域雷達(dá)目標(biāo)識(shí)別應(yīng)用的系統(tǒng)研究。