空中非合作目標旋轉部件位置估計方法*

高 磊，曾勇虎，汪連棟，王 偉

(1. 電子信息系統復雜電磁環境效應國家重點實驗室， 河南 洛陽 471003；2. 中山大學 電子與通信工程學院， 廣東 廣州 510006)

逆合成孔徑雷達(Inverse Synthetic Aperture Radar, ISAR)觀測的目標通常都是非合作目標。當使用ISAR觀測空中非合作目標時，如果存在旋轉部件，按照常規方式開展成像處理時，旋轉部件對信號的調制作用會導致成像質量的下降[1-4]，從而降低后續散射特征提取與識別結果準確性。

在存在旋轉部件的情況下，為了提高ISAR圖像質量，常見的處理思路是從高分辨距離像(High Resolution Range Profile, HRRP)數據中分離出旋轉部件對應的微動分量，文獻[5-9]分別按照不同的方式實現了這個目標。在實現旋轉部件對應的微動分量分離過程中，一個重要的前提是確定微動分量在HRRP中的位置。

文獻[1]使用經驗模態分解(Empirical-Mode Decomposition, EMD)得到的固有模式函數(Intrinsic-Mode Function, IMF)分量中的零點個數來估計旋轉部件位置，文獻[10]則根據徑向距離單元熵值的大小來確定旋轉部件位置。文獻[1]的做法需要對全部距離單元進行EMD運算，計算量較大；而文獻[10]中根據熵值大小估計位置的原則不夠穩健，如果存在散射幅度強的非旋轉部件，則也會被認為是旋轉部件。

基于上述情況，在分析旋轉部件在ISAR圖像中分布特點基礎上，本文提出基于方位向距離單元熵和局部徑向距離單元熵的旋轉部件位置估計方法。

1 包含旋轉部件的寬帶成像特點分析

光學區條件下，目標一般都滿足散射點模型假設，即雷達回波可等效為一系列散射點回波的疊加。下面建立包含旋轉部件目標的散射點回波模型。

實際ISAR在對空中目標進行成像時，目標相對于雷達的運動可以分解為兩部分：徑向運動和切向運動。徑向運動是指目標沿著雷達波照射方向的運動，切向運動是指目標垂直于雷達波照射方向的運動。徑向運動不利于方位向散射點的高分辨，故在成像處理中需要通過運動補償的方式，補償徑向運動，僅保留目標的切向運動，進而將一般運動目標轉化成轉臺運動目標[11]。為簡化分析，這里直接考慮轉臺目標的成像觀測模型。為描述目標相對雷達的運動情況，建立雷達坐標系T-UV、目標坐標系O-XY和局部坐標系o-ξη，如圖1所示。假設目標由主體部分和微動部分組成，不失一般性，可設微動部件運動為旋轉運動。在成像期間，目標主體以角速度ω0繞轉動中心O勻速轉動，旋轉部件在伴隨主體運動的同時，還以角速度ω1繞其自身轉動中心o做高速旋轉。

圖1 含旋轉部件目標的成像觀測模型Fig.1 Imaging observation geometry model for target with rotating part

設雷達發射寬帶線性調頻信號為：

(1)

(2)

注意到，第二個求和符號中，旋轉部件散射點位置隨慢時間的變化規律不再是線性變化，與主體散射點相比，多了一個余弦函數。

對比式(2)中主體散射點與旋轉部件散射點位置隨慢時間變化的規律可知：

1)主體散射點位置隨慢時間變化呈現出線性變化，在旋轉部件散射強度不是特別強的前提下，經過包絡對齊之后，主體散射點在不同慢時間HRRP中的位置被對準，進而，在處理多幀HRRP進行方位向分辨得到ISAR圖像時，能夠實現良好聚焦。

2)旋轉部件散射點位置隨慢時間變化呈現線性變化疊加正弦變化，經過包絡對齊之后，旋轉部件散射點在相鄰HRRP中的位置無法對準。設信號分辨率為Δr，上述旋轉部件散射點在不同慢時間HRRP中的位置會在?2rq/Δr」個距離單元范圍內振動。進而，在處理多幀HRRP進行方位向分辨得到ISAR圖像時，對于存在旋轉部件的距離單元，沿著方位向會呈現嚴重散焦，即在ISAR圖像中，存在一個或多個徑向距離單元沿方位向呈現條帶狀散射結構分布。

2 旋轉部件位置估計方法

結合前一部分對主體散射點和旋轉部件散射點在ISAR圖像中表現特點的分析結果，本部分在徑向距離單元熵的基礎上，定義方位向距離單元熵和局部徑向距離單元熵，進而提出旋轉部件位置估計方法。

文獻[10]定義的徑向距離單元熵如下：設ISAR圖像矩陣為g(k,n)，k(k=1,2,…,K)表示方位向距離單元序號，n(n=1,2,…,N)表示徑向距離單元序號，則徑向距離單元熵ε(n)為：

(3)

下面給出方位向距離單元熵和局部徑向距離單元熵的定義。

定義1方位向距離單元熵εC(k)為：

(4)

式(4)反映了ISAR圖像中對應每個方位向距離單元上的能量與整體能量的加權比值。

式(3)和式(4)的形式具有相似性，但其含義是不同的，圖2直觀地給出了這兩個公式的含義對比。在K×N維ISAR圖像矩陣中，式(3)對應的是矩陣中第n列元素在整體圖像能量中的加權比值，而式(4)對應的是矩陣中第k行元素在整體圖像能量中的加權比值。ISAR圖像中的行序號和列序號具有不同含義，行序號對應方位向距離單元序號，而列序號對應徑向距離單元序號，因此式(3)和式(4)實際上是從兩個角度說明ISAR圖像中的能量分布情況的。

圖2 兩種距離單元熵含義對比Fig.2 Comparison of two range cell entropy

定義2設Γ為ISAR圖像矩陣g(k,n)中按照一定規則確定的部分行序號集合(即部分方位向距離單元序號集合)，則局部徑向距離單元熵為：

(5)

旋轉部件在HRRP中的位置估計方法步驟如下：

Step1：對ISAR圖像數據按照式(4)計算方位向距離單元熵εC(k)；

Step2：找出方位向距離單元序號kmax，使得εC(kmax)=max(εC(k))；

備注1：kmax通常在?K/2」附近，如果存在極端情況，則考慮按照循環移位的方式在1,2,…,K范圍內確定kmax前后?K/W」區間作為目標區域。

備注2：Λ為預先設定值，如果成像結果中噪聲較強，則經驗取值為5～8；如果成像結果中噪聲很弱，則經驗取值為10～20。

上述估計方法的內在物理意義有兩點：第一，沿著方位向，計算每個距離單元上的能量與整體能量的加權比值，以最大值前后?K/W」區間作為目標區域；第二，在存在旋轉部件情況下，在目標主體區域之外的部分，旋轉部件對應的徑向距離單元能量應該遠高于其余的徑向距離單元。

通過上述處理，可以實現ISAR圖像中旋轉部件位置的自動估計，進而能夠支撐從HRRP中分離旋轉部件信號分量，以提升目標主體部分ISAR圖像質量。

3 數據驗證與性能分析

3.1 仿真數據驗證

通過仿真實驗，對上述方法進行驗證。仿真驗證設置如下：

1)目標模型：主體散射結構是由44個散射點組成的飛機狀散射點模型，如圖3所示(圖中坐標軸對應目標坐標系)。主體散射結構散射點強度設置為3～8之間不等，但將機身中間位置個別散射點的強度設置為12。在上述模型中，設置兩個旋轉部件，位于機翼兩側，如圖3中用五角星標出的位置。兩個旋轉部件的轉動參數相同：旋轉半徑為0.2 m，旋轉周期20 r/s。兩個旋轉部件散射點強度設置為11。

圖3 仿真目標散射點模型Fig.3 Scatter model of target for simulation

2)觀測位置：雷達坐標系與目標坐標系的坐標軸平行，目標坐標系原點在雷達坐標系中的坐標是(8000 m,8000 m)，目標坐標系在雷達坐標系中的X軸方向速度為-230 m/s，Y軸方向速度為-50 m/s。

3)雷達參數：載頻為10 GHz；帶寬為1 GHz；脈寬為50 μs；寬帶信號脈沖重復頻率為500 Hz。

在上述參數設置下，成像觀測時間為1 s，對500幀寬帶雷達回波進行成像處理，其ISAR圖像結果和徑向距離單元熵結果如圖4所示。

在圖4(a)中，由于旋轉部件的存在，使得在ISAR圖像中，沿著方位向，出現兩個明顯的散射條帶，這就是旋轉部件存在的距離單元范圍。對比圖3和圖4(a)可知，ISAR圖像中旋轉部件的位置分布與目標模型中旋轉部件的位置設定是一致的。具體來說，在ISAR圖像結果中，存在旋轉部件的徑向距離單元序號是300，340。而由圖4(b)可知，其徑向距離單元熵最大的距離單元序號為320，并不在圖4(a)ISAR圖像提示的旋轉部件距離單元范圍內。

(a) 仿真數據ISAR圖像(a) ISAR image of simulation data

(b) 徑向距離單元熵(b) Range cell entropy圖4 飛機模型仿真數據ISAR圖像與徑向距離單元熵Fig.4 ISAR image of simulation data for airplane model and its range cell entropy

按照第2節的方法對旋轉部件進行位置估計。相應的方位向距離單元熵和局部徑向距離單元熵結果如圖5所示。根據圖5(a)的結果確定目標主體部分的位置后，對ISAR圖像余下的部分進行計算，得到圖5(b)的結果。

從圖5(b)中的局部徑向距離單元熵計算結果可知，旋轉部件對應的位置分別是：序號300和340附近的距離單元，這與ISAR圖像結果以及仿真設置是相同的，從而驗證了所提方法的有效性。在仿真數據中未疊加噪聲，故在進行旋轉部件位置估計時，門限設定為熵均值的20倍，即Λ=20。

(a) 方位向距離單元熵(a) Cross range cell entropy

(b) 局部徑向距離單元熵(b) Local range cell entropy圖5 位置估計過程中得到的兩種熵計算結果Fig.5 Result of two kinds entropy in position estimation

3.2 實測數據驗證

進一步應用某飛機目標寬帶雷達實測數據進行驗證。通過觀察其ISAR圖像結果可知，存在顯著的旋轉部件。處理該目標實測數據的500幀回波數據得到的ISAR圖像以及徑向距離單元熵如圖6所示。圖6(a)中兩條沿方位向的不均勻條帶即為旋轉部件導致的方位向散焦。圖6(b)中徑向距離單元熵最大的位置對應的距離單元序號是512，從ISAR圖像中可以看出，這個距離單元序號周圍并不對應旋轉部件，即文獻[10]給出的旋轉部件距離單元確定規則不夠穩健。

(a) 實測數據ISAR圖像(a) ISAR image of real data

(b) 徑向距離單元熵(b) Range cell entropy圖6 某飛機實測數據ISAR圖像與徑向距離單元熵Fig.6 ISAR image of real data for an airplane and its range cell entropy

按照第2節的方法，進行位置估計，其方位向距離單元熵和局部徑向距離單元熵結果如圖7所示。從圖7(a)中確定目標主體部分所對應的方位向距離單元序號后，對ISAR圖像中余下部分進行局部徑向距離單元熵的計算，得到圖7(b)的結果。由圖7(b)可知，旋轉部件對應的距離單元序號分別對應490，610，這個結果與ISAR圖像中的觀測結果一致。在圖7(b)中，因為噪聲的存在，在進行旋轉部件位置檢測時，對應的門限是熵均值的5倍，即Λ=5。

(a) 方位向距離單元熵(a) Cross range cell entropy

(b) 局部徑向距離單元熵(b) Local range cell entropy圖7 位置估計過程中得到的兩種熵計算結果Fig.7 Result of two kinds entropy in position estimation

3.3 穩健性分析

通過上述實測數據處理，驗證了位置估計方法的有效性，進一步考慮在實測數據中疊加不同幅度的額外噪聲，分析該方法的穩健性。

所疊加的噪聲模型為：

(6)

式中，randn(1,N)表示長度為N的高斯隨機序列，Am為噪聲幅度。結合實測數據的情況，N=1024，將Am取兩個值1500和2500。

疊加不同幅度的額外噪聲之后的ISAR圖像與未疊加額外噪聲的ISAR圖像之間的對比如圖8所示。在Am=1500的情況下，目標主體散射結構出現一定模糊，尤其是主體部分，但仍能分辨出成像結果；而在Am=2500的情況下，已經幾乎難以分辨目標的散射結構。

按照上述方法對疊加不同額外噪聲情況下的ISAR圖像進行旋轉部件位置估計，所得到的局部徑向距離單元熵結果如圖9所示。從圖中可以看出，在上述兩種不同額外噪聲幅度下，局部徑向距離單元熵的幾個極大值對應的位置幾乎相同，因此估計出的旋轉部件位置是相同的，從而可以認為上述方法對存在噪聲的情況具有良好的穩健性，只要ISAR圖像中能勉強分辨出目標的散射結構，便可以完成位置估計。

(a) Am=0的ISAR圖像(a) ISAR image for Am=0

(b) Am=1500的ISAR圖像(b) ISAR image for Am=1500

(c) Am=2500的ISAR圖像(c) ISAR image for Am=2500圖8 疊加不同幅度額外噪聲時的實測ISAR圖像Fig.8 ISAR images of real data when additional noise with different amplitude existed

(a) Am=1500局部徑向距離單元熵(a) Local range cell entropy when Am=1500

(b) Am=2500局部徑向距離單元熵(b) Local range cell entropy when Am=2500圖9 疊加不同幅度額外噪聲時的局部徑向距離單元熵Fig.9 Local range cell entropy result when additional noise with different amplitude was added

4 結論

仿真數據和實測數據結果驗證了所提旋轉部件位置估計方法的有效性，通過在實測數據中疊加不同幅度的額外噪聲，進一步檢驗了所提出方法對噪聲的穩健性。數據處理結果顯示，只要能夠從數據中獲取較為清晰的散射結構，就可以應用本文方法實現ISAR圖像中旋轉部件位置的自動估計。通過估計旋轉部件位置，可以有效支撐HRRP中旋轉部件分量分離，進而提升主體部分成像質量。