雙基地角時變下的ISAR等效旋轉中心估計

2019-03-29 06:45:48郭寶鋒孫慧賢胡文華尹文龍陳關軍

航空學報 2019年3期

郭寶鋒，孫慧賢，胡文華,*，尹文龍，陳關軍

1. 陸軍工程大學石家莊校區電子與光學工程系，石家莊 050003 2. 中國人民解放軍32136部隊，張家口 075000

雙基地雷達是發射機和接收機分置、且基線長度與目標距離可比擬的雷達系統，收發分置的工作方式使雷達在對抗“四大威脅”方面具有突出的優勢[1]。雙基地逆合成孔徑雷達[2-3]是基于雙基地雷達平臺的逆合成孔徑雷達(Inverse Synthetic Aperture Radar, ISAR)系統，在具備“四抗”特性的同時，利用接收的目標非后向散射回波進行成像，能夠較單基地雷達獲取更加豐富的目標信息[4-6]。為此，雙(多)基地雷達系統日益成為現代雷達研究的熱點問題[7-11]。距離-多普勒(Range-Doppler, RD)算法[12-13]是雙基地ISAR成像的經典算法，由于其物理意義明確、操作簡便被廣泛應用于成像仿真及實測數據的處理中，但當觀測目標尺寸過大或成像累積轉角過大時，該算法存在越分辨單元徙動問題[14-15]，造成圖像散焦，離旋轉中心越遠，散射點散焦現象越嚴重。為了完成越分辨單元徙動的校正從而提高成像質量，需要估計目標的等效旋轉中心。

目標等效旋轉中心估計可由基于單特顯點的運動補償算法近似得到，但該算法得到的旋轉中心位置精度較低，影響后續圖像校正的聚焦效果，并且，不是所有的圖像都有理想的單特顯點單元[13]，尤其是實測數據中，回波信噪比一般很低，難以找到單一的特顯點，因此，該算法的應用受到很大的限制。針對平動補償不能提供等效旋轉中心的情況，清華大學的葉春茂等[16-17]針對單基地ISAR旋轉中心估計問題，提出了利用3個散射點在兩幅圖像中的位置差提取旋轉中心的算法，但由于雙基地ISAR成像時，雙基地角變化會引起圖像畸變[18]，使得該算法無法適用于雙基地ISAR，文獻[19]利用相似原理實現了等效旋轉中心估計，但未考慮雙基地角的變化，也使得算法應用范疇受限。文獻[20]考慮到雙基地ISAR的畸變特性，基于圖像對比度最大準則估計了圖像的等效旋轉中心位置，該估計算法采用搜索方式，每搜索一個距離點都需要構造補償相位項，并作快速傅里葉變換(Fast Fourier Transform, FFT)變換完成成像，運算量很大，并且等效旋轉中心估計精度受限于所構造補償相位的精確程度。

基于此，本文針對雙基地角時變下的ISAR等效旋轉中心估計問題，提出了一種等效旋轉中心估計算法。首先，建立了雙基地ISAR成像模型，分析了雙基地ISAR成像特性，并推導了圖像畸變角度的表達式；然后，提出了等效旋轉中心估計算法，該算法將運動補償后脈沖壓縮數據分成兩組成像，并進行畸變校正，通過圖像旋轉相關搜索等效旋轉中心位置；最后，進行了仿真實驗，驗證了算法有效性。

1 雙基地ISAR成像模型

圖1 雙基地ISAR成像模型Fig.1 Imaging model for bistatic ISAR

由于雷達本身具有精確的測距能力，RTOm、RROm、L等距離信息可通過雷達獲得，則雙基地角β0、βm及ζ0、ζm的值可由三角形關系得到

(1)

(2)

(3)

(4)

(5)

(6)

此時累積轉角θm為

θm=η0-ηm

(7)

假設收發雙站雷達理想同步，發射站雷達以脈沖重復周期TPRT發射線性調頻信號

(8)

其基頻信號為

(9)

設第p個散射點A的散射系數為σp，tm時刻，散射點A到收發雙站的距離和為Rpm，即Rpm=RTpm+RRpm，則接收站雷達接收到的回波信號為

(10)

經過相參本振下變頻至零中頻，并使用匹配濾波器對目標回波進行脈沖壓縮，得到散射點A的一維距離像為

(11)

雙基地雷達系統中，目標尺寸遠小于到收發雙站距離，即d?RTO,RRO，由于αm=θm+α0，則散射點A到發射站和接收站的距離可表示為

(12)

(13)

由式(12)和式(13)可得pm到收發雙站距離和為

Rpm=(RTOm+RROm)+2(xpsinθm+

(14)

令ROm=RTOm+RROm，為目標質心到收發雙站的距離和，則散射點A到收發雙站距離和與目標質心到收發雙站距離和之差為

(15)

經過理想的運動補償(包括包絡對齊和相位校正兩個步驟)后，回波的一維距離像式(11)可表示為

(16)

假設成像期間轉動角速度恒定為ω，且雙基地角固定不變恒為β，由于雷達視角變化很小，累積轉角可作如下近似：sinθm≈θm，cosθm≈1，且θm=ωtm。因此，式(16)可表示為

(17)

(18)

2 雙基地ISAR成像特性分析

式(17)和式(18)的原理推導是在雙基地角恒定及目標等效轉速均勻的情況下給出的，并且累積轉角作了近似，而在實際成像過程中，雙基地角恒定的假設很少滿足。下面推導雙基地角時變情況下慢時間多普勒表示，對式(16)中的指數相位項求導可得散射點精確的瞬時多普勒信息

(19)

為便于理解，首先分析單基地情況下，令式(19)中雙基地角為0，并使fd=0，可得瞬時多普勒為零的散射點的集合(即雷達觀測坐標系uOmv坐標系中的縱軸)與x軸(方位向)的夾角γm為

(20)

為便于分析，單基地ISAR瞬時成像結果如圖2(a)所示。

從式(20)可知，由于目標轉動，成像初始時刻在y軸(距離向)上的點現在偏移了θm角度，單基地ISAR完全是散射點模型旋轉后的圖像顯示。

對于雙基地ISAR，根據式(19)，雙基地角時變時，慢時間瞬時多普勒fd是方位坐標xp、距離坐標yp的復雜函數關系。成像時，瞬時多普勒為零的散射點的集合(即雷達觀測坐標系uOmv坐標系中的縱軸)與x軸的夾角γm為

(21)

雙基地ISAR瞬時成像結果如圖2(b)所示，設多普勒為零的散射點與距離向y軸的夾角為φm，則

圖2 單/雙基地ISAR成像特性對比Fig.2 Comparison of imaging characteristics of monostatic and bistatic ISAR

(22)

以另一種形式書寫式(19):

(23)

(24)

(25)

(26)

3 雙基地ISAR等效旋轉中心估計

根據雙基地ISAR的成像特性分析，得到了圖像的畸變角度，據此可以構造平移量，通過物理平移消除雙基地ISAR的圖像畸變，這樣，圖像與散射點模型之間只存在一個繞旋轉中心的旋轉，旋轉角度由成像視角差決定。基于這個原理，本文提出基于圖像旋轉相關的雙基地ISAR等效旋轉中心估計算法，首先將成像過程按累積轉角和雙基地角把一維距離像序列分成兩組，并分別成像，然后，對兩幅圖像進行畸變校正，則此時兩幅圖像是散射點模型旋轉不同角度的顯示，那么這兩幅圖像繞旋轉中心旋轉可達到最大程度的吻合，即通過圖像旋轉相關性最大準則對旋轉中心位置進行搜索，相關性最大時對應的距離單元就是等效旋轉中心位置。

3.1 圖像畸變校正

假設成像期間積累脈沖個數為M，累積轉角為θM，每個脈沖快時間采樣點數為N，快時間采樣頻率為fs，畸變角度為φdis_A，平均雙基地角為βA，雷達發射LFM信號的載波波長為λ，則每個距離采樣單元和每個多普勒單元代表的長度為

(27)

(28)

為了消除畸變，若以第k個距離單元為基準，則第n個距離單元在方位向需要移動的數值為

Δyn=(n-k)Δytanφdis_A

(29)

對應移動的多普勒單元個數為

(30)

式中:round(·)表示對數值四舍五入取整，ΔMn>0時向左移位，ΔMn<0時向右移位。對每個距離單元按式(30)進行移位操作即可完成雙基地ISAR的圖像畸變校正。式(30)的取整操作會存在量化誤差的問題，影響成像結果，這可通過對圖像進行4～8倍插值細化，減小量化誤差對成像的影響。

3.2 等效旋轉中心估計

(31)

設散射點A在二維圖像中的像素位置可表示為(Xp,Yp)，則坐標(xp′,yp′)與像素位置(Xp,Yp)滿足如下關系

xp′=(Xp-Xc)Δx

(32)

yp′=(Yp-Yc)Δy

(33)

聯合式(31)～式(33)，整理可得

(34)

現將回波運動補償后的數據分成兩組，為盡可能減小后續的估計誤差，分組后應使兩組數據的多普勒分辨率相同，由于成像期間雙基地角的變化，會導致兩組回波脈沖個數未必相同。假設兩組數據的脈沖個數分別為M1、M2，不妨設M1>M2，對兩組數據分別成像，并對脈沖積累個數多的圖像(M1對應的圖像)以其多普勒中心進行截短，截取其中的M2個多普勒單元，使兩幅圖像數據維數一致，得到ISAR圖像1和ISAR圖像2，此時兩幅圖像數據維數一致。設兩幅圖像的累積轉角分別為θM1、θM2，平均雙基地角分別為βA1、βA2，由于兩幅圖像的方位向分辨率相同，則θM1cos(βA1/2)=θM2cos(βA2/2)，那么兩幅圖像的方位分辨率可表示為

(35)

距離定標因子分別為

(36)

(37)

由于兩幅圖像是由同一組運動補償后的數據分組得到，因此，等效旋轉中心在同一距離單元上，并且，旋轉中心處多普勒值必為零，因此，兩幅圖像的等效旋轉中心位置相同，設為(Xc,Yc)，以該中心對兩幅圖像進行畸變校正，消除畸變帶來的方位偏移，此時兩幅圖像僅有一定的視角差，視角差為Δθ=(θM1+θM2)/2。

設(X1,Y1)、(X2,Y2)分別為散射點A在兩幅圖像中的像素位置，依照式(34)可得到散射點坐標與在圖像上像素位置對應關系

(38)

(39)

根據式(38)和式(39)，散射點A在兩幅ISAR圖像中的像素位置關系滿足

(40)

(41)

式中：Φ為兩幅圖像的坐標轉換矩陣。式(41)表明，第1幅ISAR像的散射點相對等效旋轉中心的位置可由第2幅圖像相對等效旋轉中心的位置經變換矩陣Φ得到。雙基地ISAR的等效旋轉中心在零多普勒線上，只需假定等效旋轉中心的對應距離單元位置，對第2幅圖通過Φ進行矩陣變換，用得到的圖像與第1幅圖作互相關，假設兩幅圖像分別記為I1和I2，其有效距離單元和多普勒單元個數分別為N和M2，則相關系數的計算公式為

(42)

當相關性最大時，該假定的旋轉中心位置即為目標的等效旋轉中心。需要說明的是，由于成像時相位校正一般基于強散射點的相位實現，使得等效旋轉中心在散射強度較大的點附近，因此，進行距離單元搜索時，可選擇圖像亮度高的位置及其附近區域進行，以減小運算量。

3.3 算法運算量分析

本文算法需要對每個距離單元進行處理，包括畸變角度計算、圖像畸變校正、圖像旋轉、圖像相關等操作，這里給出算法的運算量分析。假設成像積累總脈沖個數為M，分組后兩組數據的脈沖個數分別為M1、M2，令M0=min(M1,M2)，每個脈沖壓縮后截取的有效距離單元個數為N，為減小運算量，等效旋轉中心估計選擇的距離單元搜索范圍為N0，各關鍵步驟對應的運算量如表1所示。

從表1可以看出，圖像畸變校正的移位操作的運算量較大，圖像旋轉會產生較多的實數乘法和實數加法運算，這兩個環節占用算法主要的運算量，其次是圖像相關步驟，畸變角度計算的運算量則很少。此外，運算量大小與選擇的距離單元搜索范圍N0直接相關，在實際操作時可根據目標的大致尺寸確定N0的大小，以最大限度減小運算量。

表1 運算量統計Table 1 Calculation statistics

4 算法流程

雙基地角時變下的ISAR等效旋轉中心估計流程如圖3所示，具體步驟為

步驟1對雙基地ISAR回波進行脈沖壓縮和運動補償，得到相位校正后的一維距離像序列。

步驟2依據累積轉角和雙基地角把一維距離像序列分成兩組，使兩組數據的多普勒分辨率一致，并分別成像得到ISAR圖像1和ISAR圖像2，對多普勒單元數目多的圖像進行數據截短，使兩幅圖像的數據維數一致。

步驟3假定某一距離單元為等效旋轉中心所在位置。

步驟4根據式(26)計算兩幅圖像的畸變角度，并通過式(30)對兩幅圖像進行畸變校正。

步驟5對圖像2依據式(41)進行旋轉，得到新的圖像，并與圖像1作互相關。

步驟6假設下一個距離單元為等效旋轉中心位置，重復步驟4～步驟5，直至距離單元遍歷結束。

步驟7根據遍歷過程，尋找兩幅圖像相關性最大時所假定的等效旋轉中心位置，該位置就是估計得到的等效旋轉中心。

圖3 等效旋轉中心估計算法流程Fig.3 Flow chart of estimation algorithm for equivalent rotation center

5 仿真實驗

為驗證等效旋轉中心估計算法的有效性，進行了仿真實驗。仿真場景如圖4所示，其中，“◇”為發射站位置，“○”為接收站位置，雷達基線長度為50 km，目標運動速度為400 m/s，目標軌道加粗部分為選擇的雷達成像區域，該區域雙基地角是時變的。仿真散射點模型如圖5所示，仿真參數如表2所示。

首先將采集的1 024個回波進行匹配濾波完成脈沖壓縮，并進行運動補償。然后根據累積轉角與雙基地角關系，將1 024個脈沖壓縮結果分為兩組，使兩組數據的多普勒分辨率相同，經計算，采用的分組方法是：第1～525個脈沖為第1組、第526～1 024個脈沖為第2組，對兩組脈沖構成的一維距離像序列分別成像，得到ISAR圖像1和ISAR圖像2，如圖6所示，這兩幅圖像都有一定的畸變，畸變角度理論值約為26°，向左偏移，根據圖中的散射點坐標也可計算得到畸變角度為26.1°，與理論值吻合。

圖4 仿真場景Fig.4 Simulation scenario

圖5 散射點模型Fig.5 Scatter model

表2 成像仿真參數Table 2 Imaging simulation parameters

由于第2幅圖像方位選取了499個脈沖，小于第1幅圖像的525個脈沖，為了使后續的圖像旋轉具有相同的數據維數，將第1幅圖像數據截取，保留中間的499個多普勒單元。假定等效旋轉中心在某一距離單元(等效旋轉中心在多普勒單元的位置為零多普勒線上，無需估計，即像素位置Xc=250)，依據畸變角度理論值，對兩幅圖像進行畸變校正；對第2幅圖像進行旋轉，并與第1幅圖像作相關，存儲相關系數；然后，假定下一個距離單元為旋轉中心位置，再進行兩幅圖像的畸變校正、旋轉相關，直至遍歷整個距離單元。為減小運算量，本實驗對目標附近的距離單元進行遍歷，得到的歸一化相關系數曲線如圖7中實線所示，峰值點出現在第500個距離單元處，則等效旋轉中心在兩幅圖像上的像素位置為(250,500),為驗證本文算法的性能，將其與文獻[19]給出的算法進行對比，文獻[19]未考慮雙基地角變化及圖像畸變，采用該算法得到等效旋轉中心搜索曲線如圖7中虛線所示，峰值位置在524處，與本文算法相差24個距離單元。圖8(a)和8(b)分別是ISAR圖像1和ISAR圖像2依據旋轉中心畸變校正后的圖像，此時，圖像的畸變現象消除。

根據得到的等效旋轉中心，采用文獻[20]的越分辨單元徙動校正算法對運動補償后的1 024個脈沖進行徙動校正，并分別成像。圖9是未作校正的1 024個脈沖雙基地ISAR原始成像結果，圖10是越距離單元徙動校正后的ISAR像，圖11(a)和11(b)分別是采用文獻[19]和文獻[20]的算法估計等效旋轉中心并進行越分辨單元徙動校正得到的ISAR像(文獻[20]算法估計旋轉中心的距離單元位置為527)，圖11(c)是本文算法估計旋轉中心并進行越分辨單元徙動校正得到的ISAR像，對比圖10和圖11可以看出，圖11的成像質量均優于圖10，但圖11(a)和圖11(b)較圖11(c)中散射點出現變“胖”現象，這是由于文獻[19]和文獻[20]的算法估計旋轉中心存在誤差，用該估計值未能完全校正圖10的越分辨單元徙動，而本文算法適用于雙基地角時變的情形，等效旋轉中心估計準確，用該估計中心校正越分辨單元徙動使圖像質量達到了最佳。

圖6 分組成像結果(畸變未校正)Fig.6 Imaging results after grouping (without distortion correction)

圖7 旋轉中心搜索曲線Fig.7 Rotation center estimation curves

圖8 分組成像結果(畸變校正后)Fig.8 Imaging results after grouping (after distortion correction)

圖9 雙基地ISAR原始圖像Fig.9 Original image of bistatic ISAR

圖10 越距離單元徙動校正后ISAR圖像Fig.10 ISAR image after migration through range cell correction

圖11 越分辨單元徙動校正后的ISAR圖像Fig.11 ISAR images after migration through resolution cell correction

此外，為定量分析成像質量，計算3幅圖像的圖像對比度，并統計圖像上所有散射點的距離向和方位向3 dB主瓣平均寬度，如表3所示，數據對比表明了本文等效旋轉中心估計算法的準確性。

表3 圖像質量定量分析Table 3 Quantitative analysis of image quality

雙基地ISAR的圖像畸變校正利于后續的目標識別，根據圖像畸變角度，可通過不同距離單元的移位操作實現圖像的畸變校正，畸變校正后的結果如圖12所示，其中圖12(a)是RD直接成像結果(圖9)校正后的圖像，圖12(b)是圖11(a)圖像畸變校正后的結果，圖12(c)是圖11(b)圖像畸變校正后的結果，圖12(d)是圖11(c)圖像畸變校正后的結果，可以看出畸變校正后的圖像與仿真模型具有一致性(可通過定量計算散射點之間的夾角確定)，只存在視角的差異。

此外，為驗證算法運算量，表4給出了算法及關鍵步驟運算時間統計，其中，距離單元搜索范圍分別設置為400個和800個。仿真使用MATLAB軟件，電腦配置為Win7操作系統、i7-4510 CPU、2.6 GHz主頻、8 G內存。從表4可以看出，圖像畸變校正過程所用時間最長，這是由于移位操作所用運算量較大，其次，圖像旋轉需要較大運算量。對比表中上下兩組數據，當搜索距離單元個數不同時，算法運算總時間差異較大，若要進一步減小算法運算量，可通過減小距離單元搜索范圍實現。