近場單發(fā)多收合成孔徑雷達成像的頻域算法

朱榮強, 周劍雄, 付　強

(國防科技大學電子科學學院ATR實驗室, 湖南 長沙 410073)

0　引　言

合成孔徑雷達(synthetic aperture radar,SAR)由于其獨特的優(yōu)勢,已被廣泛用于遙感、監(jiān)測等領域。但是,傳統(tǒng)的SAR只能獲得目標的二維散射系數(shù)分布,因此丟失了目標的部分信息,并且成像結果存在一些缺點。近年來,將多站陣列與孔徑合成技術相結合的新成像體制成為了研究熱點[1-7]。該成像結構利用與孔徑合成方向相正交的陣列天線獲得三維成像能力,但是與SAR成像相比,由于多站陣列中收發(fā)陣元分置,使得其成像處理變得更為復雜。

在SAR成像中,距離遷徙算法(range migration algorithm,RMA)是一種常用的精確成像算法。該方法通過將回波數(shù)據(jù)變換至波數(shù)域進行補償,有效提高了成像效率并且保證了成像精度[8-9]。但是,RMA不能直接用于多站陣列成像。為了將RMA拓展用于多站陣列成像,文獻[10]根據(jù)等效相位中心原理將多發(fā)多收(multiple-input multiple-output,MIMO)陣列回波數(shù)據(jù)轉換為單站格式,并利用參考信號對轉換后的數(shù)據(jù)進行補償;文獻[11-12]利用多平面波近似,將RMA拓展至單發(fā)多收(single-input multiple-output,SIMO)陣列成像,并將波數(shù)空間進行分割以降低近似誤差;文獻[13]利用二階泰勒展開,推導了SIMO-SAR的近似波數(shù)譜。這些方法利用近似簡化成像處理的同時降低了成像精度。文獻[14]將RMA拓展用于近場MIMO陣列,但是要求收發(fā)陣元的采樣間隔必須滿足奈奎斯特采樣定理,無法適用于一般的MIMO陣列。文獻[15]提出了一種適用于SIMO陣列的頻域成像算法,該方法利用基于快速傅里葉變換(fast Fourier transform, FFT)的子成像結構有效提高了成像效率并且避免了插值,但是該方法不能直接用于SIMO-SAR成像。

在多站陣列成像處理中,后向投影算法(backprojection algorithm,BPA)是最為常用的算法之一。該方法逐點相干累加成像,因此具有很高的成像精度,并且可以適用于多種成像體制。但是,用于三維成像時,BPA的成像效率較低。

本文提出了一種適用于SIMO-SAR的頻域成像算法。該方法首先將回波數(shù)據(jù)變換至波數(shù)域,然后對波數(shù)域數(shù)據(jù)進行擴維以實現(xiàn)補償和插值,最后利用降維累加得到了散射系數(shù)分布的三維傅里葉變換。該方法在波數(shù)域精確補償了波前彎曲,因此可以用于近場三維高精度成像,并且該方法采用頻域成像結構顯著提高了成像效率。

1　成像算法

SIMO-SAR的成像幾何模型如圖1所示。

圖1　SIMO-SAR成像幾何模型

以目標中心為原點建立直角坐標系O-XYZ,X軸、Y軸和Z軸分別表示方位向、距離向和高度向。SIMO陣列由一個發(fā)射陣元和NR個均勻分布的接收陣元組成,陣列平行于X軸并沿高度向進行掃描。將接收陣元的位置記為(u,y0,h),發(fā)射陣元的位置記為(xt,y0,h),目標上任意一點的位置記為(x,y,z),對應散射系數(shù)為σ(x,y,z)。假設發(fā)射信號為步進頻信號,則回波信號可以表示為

s(u,k,h)=?σ(x,y,z)e-jkRte-jkRrdxdydz

(1)

式中,k=2πf/c表示波數(shù),f表示頻率,c表示光速;Rt和Rr分別表示發(fā)射陣元和接收陣元與任意目標間的距離,具體關系式為

(2)

(3)

由時域相關法的成像原理可得散射系數(shù)分布的估計為

?s(u,k,h)ejkRtejkRrdudkdh

(4)

式中，指數(shù)項表示球面波傳播,利用球面波分解原理可以將其分解為

(5)

(6)

式中,kv和ku的具體關系式為

(7)

(8)

將式(5)和式(6)代入式(4)可得

(9)

式中,S′(ku,k,h)表示回波信號沿u方向的傅里葉變換。式(9)中的第3項利用球面波分解可得

(10)

其中

(11)

利用式(10),式(9)可以簡化為

ej(ku+kv)xejkyyejkzzdkudkvdkdkz

(12)

式中，S(ku,k,kz)表示信號S′(ku,k,h)沿h方向的傅里葉變換。為了能夠補償指數(shù)項e-jkyy0,以及便于實現(xiàn)插值,將信號S(ku,k,kz)沿kv向進行擴維,擴維前后信號的關系為

(13)

ej(ku+kv)xejkyyejkzzdkudkvdkdkz

(14)

ejkyyejkzzdkudkvdkydkz

(15)

(16)

降維累加運算將分布在ku-kv域內(nèi)的數(shù)據(jù)累加至kx域,降維累加后的信號與散射系數(shù)分布之間的關系為

?F(kx,ky,kz)ejkxxejkyyejkzzdkxdkydkz

(17)

從式(17)可以看出,對降維累加后的數(shù)據(jù)利用三維逆傅里葉變換就可以得到散射系數(shù)分布。

2　算法實現(xiàn)及運算量分析

算法的流程圖如圖2所示。

具體成像步驟如下:

步驟1利用二維FFT得到回波信號沿u和h方向的傅里葉變換S(ku,k,kz)。為了避免波數(shù)譜產(chǎn)生混疊,回波數(shù)據(jù)在方位向和高度向的采樣間隔Δu和Δh必須滿足

(18)

(19)

式中,λmin表示最大工作頻率對應的波長;θx和θz分別表示天線在方位向和高度向的波束角。為了便于分析運算量,近似認為所有FFT均是基二FFT。該步驟的運算量近似為

C1≈5NuNkNzlog2(NuNz)FLOP

(20)

式中,Nu和Nz分別表示回波數(shù)據(jù)沿方位向和高度向傅里葉變換的點數(shù);Nk表示回波數(shù)據(jù)的頻率采樣點數(shù)，運算量的單位為浮點運算(floating-point operation,FLOP)。

步驟2將波數(shù)域數(shù)據(jù)沿kv向進行擴維,并補償指數(shù)項e-jkyy0。擴維后的信號在kv向的分布區(qū)間由式(7)決定。分布間隔應與波數(shù)域數(shù)據(jù)在ku向的分布間隔相等,以便于實現(xiàn)降維累加。該步驟的運算量為

C2=6NuNkvNkNzFLOP

(21)

式中,Nkv表示擴維后的信號在kv向的長度。

步驟3對擴維后的信號沿ky向進行插值,得到在ky域內(nèi)均勻分布的數(shù)據(jù)。采用線性插值,該步驟的運算量近似為

C3≈6NuNkvNyNzFLOP

(22)

式中,Ny表示插值后的數(shù)據(jù)在ky向的長度。

步驟4將分布在ku-kv-ky-kz域的數(shù)據(jù)降維累加至kx-ky-kz域。該步驟的運算量近似為

C4≈2(NuNkv-Nx)NyNzFLOP

(23)

式中,Nx表示降維累加后的信號在kx域中的分布長度。

步驟5利用三維IFFT得到成像結果。該步驟的運算量近似為

C5≈5NxNyNzlog2(NxNyNz) FLOP

(24)

由分析可得該算法的運算量約為

Cp≈5NuNkNz[log2(NuNz)+1.2Nkv]+8NuNkvNyNz+

5NxNyNz[log2(NxNyNz)-0.4]FLOP

(25)

在相同成像條件下,BPA的運算量為

CBPA≈8NRNhNxNyNzFLOP

(26)

式中,Nh表示回波數(shù)據(jù)在高度向的采樣點數(shù)。

比較式(25)和式(26)可以看出,本文的方法與BPA相比顯著降低了運算量。

3　實驗驗證

3.1　仿真數(shù)據(jù)驗證

采用仿真實驗比較近似RMA[11-12]、BPA和本文方法的成像效果。仿真目標由9個理想點目標組成,位置如圖3所示。SIMO陣列由1個發(fā)射陣元和101個接收陣元組成,接收陣元間隔6 mm,發(fā)射陣元位于接收陣列的中心。SIMO陣列沿高度向的掃描間隔為3.5 mm,采樣點數(shù)為161。仿真發(fā)射信號為步進頻信號,工作頻率為30～40 GHz,頻率步進間隔為100 MHz。目標距離陣列300 mm。3種算法的仿真成像結果如圖4所示。從圖4(a)和圖4(b)可以看出,采用近似RMA成像時,非原點處目標的成像結果出現(xiàn)了散焦和位置偏移,這是因為,近似RMA在成像處理中采用了平面波近似,只能對原點處的目標實現(xiàn)精確重構。從圖4(c)～圖4(f)可以看出,采用BPA和本文方法能很好地實現(xiàn)對目標的重構,并且本文方法成像結果的旁瓣水平略低于BPA成像結果的旁瓣水平。

圖3　仿真目標設置

采用成像結果的幾何均方根誤差和熵衡量算法的成像質(zhì)量。幾何均方根誤差可以反映成像結果的位置重構精度,幾何均方根誤差大于一個成像網(wǎng)格點間隔時,表明圖像存在位置重構誤差;熵可以反映圖像的聚焦效果,熵越小則圖像的聚焦越好。兩個指標的計算公式分別為

(27)

(28)

(29)

由式(27)和式(28)計算得到的幾何均方根誤差和熵如表1所示。

圖4　仿真成像結果

參數(shù)近似RMABPA本文方法幾何均方根誤差/mm6.920.460.46熵11.2911.1810.67

從表1可以看出,近似RMA的幾何均方根誤差大于一個像素間隔,表明近似RMA對目標的位置重構存在誤差,這會導致成像結果出現(xiàn)變形,影響后續(xù)圖像解譯。BPA和本文方法的幾何均方根誤差相同,均小于一個像素間隔,表明本文方法能精確地重構出目標的幾何分布。并且本文方法成像結果的熵要小于BPA成像結果的熵。

3.2　實測數(shù)據(jù)驗證

采用實測數(shù)據(jù)驗證本文方法對實際目標的成像有效性。目標由一把剪刀和一個金屬小球組成,如圖5所示。目標距離陣列約150 mm,其余參數(shù)與仿真實驗中的參數(shù)設置相同。采用本文方法和BPA的成像結果如圖6所示。為了便于比較兩種算法的性能,兩種算法的成像網(wǎng)格點均為256(方位向)×64(距離向)×256(高度向),并且成像處理中不采用加窗處理。

圖5　目標圖像

圖6　成像結果

比較兩個算法的成像結果可以看出,兩個算法均實現(xiàn)了對目標的精確重構,并且成像效果相同。表2給出了由式(25)和式(26)計算得到的兩個算法的近似理論運算量和實際成像時間。兩種算法均通過Matlab實現(xiàn),并且在相同環(huán)境下運行。從表2可以看出,本文方法的運算量要顯著小于BPA。

表2　理論運算量與成像時間

注1)：GFLOP=109FLOP。

4　結　論

本文針對SIMO-SAR提出了一種頻域成像算法。該方法通過將波數(shù)域數(shù)據(jù)進行擴維,在波數(shù)域實現(xiàn)了對波前彎曲的精確補償,因此成像精度較高,并且可用于近場三維成像。另外,從理論上分析了該方法的運算量要遠小于BPA。實驗結果也表明了該方法可以獲得與BPA相同的成像效果,但是顯著降低了運算量。本文的研究對于近場多站陣列成像具有一定的借鑒意義。

參考文獻：

[1] ZHUGE X,YAROVOY A G. A sparse aperture MIMO-SAR-based UWB imaging system for concealed weapon detection[J]. IEEE Trans.on Geoscience Remote Sensing,2011,49(1):509-518.

[2] GUMBMANN F, SCHMIDT L P. Millimeter-wave imaging with optimized sparse periodic array for short-range applications[J]. IEEE Trans.on Geoscience Remote Sensing, 2011, 49(10): 3629-3638.

[3] WEI M, GILLES M. Initial ARTINO radar experiments[C]∥Proc.of the 8th European Conference on Synthetic Aperture Radar, 2010: 1-4.

[4] KRIEGER G. MIMO-SAR: opportunities and pitfalls[J]. IEEE Trans.on Geoscience Remote Sensing, 2014, 52(5): 2628-2645.

[5] 溫鑫, 黃培康, 年豐,等. 主動式毫米波近距離圓柱掃描三維成像系統(tǒng)[J]. 系統(tǒng)工程與電子技術, 2014, 36(6): 1044-1049.

WEN X,HUANG P K, NIAN F, et al. Active millimeter-wave near-field cylindrical scanning three-dimensional imaging system[J].Systems Engineering and Electronics,2014,36(6):1044-1049.

[6] 周偉,劉永祥,黎湘,等.MIMO-SAR技術發(fā)展概況及應用淺析[J].雷達學報,2014,3(1): 10-18.

ZHOU W, LIU Y X, LI X, et al. Brief analysis on the development and application of multi-input multi-output synthetic aperture radar[J]. Journal of Radars, 2014, 3(1): 10-18.

[7] 董會旭, 張永順,馮存前,等.基于線陣的MIMO-ISAR二維成像方法[J]. 電子與信息學報, 2015, 37(2): 309-314.

DONG H X, ZHANG Y S, FENG C Q, et al. Two-dimensional imaging using MIMO radar and ISAR technique based on linear array[J]. Journal of Electronics and Information Technology, 2015, 37(2): 309-314.

[8] SOUMEKH M. Synthetic aperture radar signal processing with Matlab algorithms[M]. New York: Wiley, 1999.

[9] 譚維賢,洪文,王彥平,等.基于波數(shù)域積分的人體表面微波三維成像算法研究[J].電子與信息學報,2009,31(11):2541-2545.

TAN W X, HONG W, WANG Y P, et al. Three-dimensional microwave imaging algorithm for the surface of the human body based on wavenumber domain integration[J]. Journal of Electronics and Information Technology,2009,31(11):2541-2545.

[10] WILLIAM F M, JAMES D K, JANUSZ J M, et al. Development of a high-throughput microwave imaging system for concealed weapons detection[C]∥Proc.of the IEEE International Symposium on Phased Array System and Technology, 2016: 1-6.

[13] DU L, WANG Y P, HONG W, et al. A three-dimensional range migration algorithm for downward-looking 3D-SAR with single-transmitting and multiple-receiving linear array antennas[J]. EURASIP Journal on Advances in Signal Processing,2010,2010(1): 1-15.

[14] ZHUGE X, YAROVOY A G. Three-dimensional near-field MIMO array imaging using range migration techniques[J]. IEEE Trans.on Image Processing, 2012, 21(6): 3026-3033.

[15] ZHU R Q, ZHOU J X, TANG L, et al. Frequency domain imaging algorithm for single-input-multiple-output array[J]. IEEE Geoscience Remote Sensing Letters, 2016, 13(14): 1747-1751.