采用斜距波數子帶劃分的SAR圖像自聚焦方法

馬彥恒， 褚麗娜， 楊曉亮， 侯曉澤

(1.陸軍工程大學 石家莊校區，河北 石家莊 050003；2.中國電子科技集團 第五十四研究所，河北 石家莊 050000)

小型旋翼無人機載合成孔徑雷達(small rotor unmanned aerial vehicle synthetic aperture radar, SRUAV-SAR)可全天時、全天候，靈活有效地獲取觀測區域的二維高分辨圖像[1-3]。SRUAV-SAR平臺易受大氣湍流等因素的影響，偏離理想運動軌跡[4]，不再滿足方位不變假設，導致SAR圖像的散焦和幾何失真[5-6]，傳統的頻域處理算法[7-10]不能達到良好的聚焦性能。時域反向投影(back-projection, BP)[11]算法沿斜距歷程對成像網格上的每個像素點進行積分重構，較好地解決了距離-方位耦合的問題。BP算法對軌跡、帶寬、積累角等成像條件沒有限制，且沒有方位不變假設，更適用于SRUAV-SAR成像。但是，BP算法計算量大，成像效率低，快速因子分解BP(fast factorized BP, FFBP)算法以遞歸的方式顯著地減少了計算負擔，并保持了BP算法[12]的準確性和適用性。FFBP算法可以通過并行處理器硬件實現，在機載SAR成像領域具有廣闊的應用前景[13]。

通常，可以將運動誤差對SAR圖像的影響分為非系統距離單元徙動(non-systematic range cell migration, NsRCM)和方位相位誤差(azimuthal phase errors, APE)[14]，前者導致距離向散焦，影響APE自聚焦的性能，后者導致方位向聚焦退化。傳統的運動誤差補償方法是采用機載慣性導航系統(inertial navigation system, INS)和全球定位系統(global positioning system, GPS)直接測量的運動偏差在SAR成像之前進行APE和NsRCM補償。但是，SRUAV-SAR平臺的載荷能力有限，難以搭載高精度的運動傳感器，并且現有運動傳感器的精度也達不到高波段SAR成像精度要求[14]。因此基于回波數據的自聚焦技術對于SRUAV-SAR成像非常重要。

相位梯度自聚焦(phase gradient autofocus, PGA)[15]充分利用了圖像域與相位歷程域之間的快速傅里葉變換(fast Fourier transform, FFT)關系，是一種高效的空間不變相位誤差估計自聚焦方法。文獻[16，17]將基于擬極坐標系(quasi-polar grid, QPG)的FFBP成像和PGA或基于加權的PGA(weighted PGA, WPGA)相結合，實現了基于回波數據的圖像自聚焦，但未考慮沿方位向運動誤差；文獻[18]將FFBP與自動校準結合用來處理機載SAR數據，但沒有考慮NsRCM補償；文獻[19]通過一致修正NsRCM和APE來實現補償，但未考慮運動誤差的空變性。因此，基于無人機二維空變運動誤差校正的SAR數據自聚焦仍是一個亟待解決的問題。

本文提出一種斜距波數子帶劃分的SAR圖像自聚焦方法，用于SRUAV-SAR數據無INS自聚焦成像處理。首先基于QPG-FFBP成像模型，分析運動誤差在波數域的空變特征；然后基于時域到波數域的映射關系，推導NsRCM和APE的關系式，用于多子帶NsRCM補償；最后仿真和實測數據結果驗證該方法的有效性。

1 擬極坐標系SAR成像模型

1.1 SAR平臺運動模型

直角坐標系中無人機載SAR的幾何關系如圖1所示。圖中，X軸表示方位向，Y軸表示距離向，Z軸表示高度向，O為整個合成孔徑天線相位中心(antenna phase center, APC)時刻。理想飛行軌跡用綠色實線表示，與X軸重合，Ai為ta時刻理想軌跡的APC位置；實際飛行軌跡用紅色虛線表示，為不規則的曲線，Ar為ta時刻實際軌跡的APC位置。P是地面波束覆蓋區中任意一點，Ri和Rr分別表示ta時刻SAR到P點的理想和實際瞬時斜距。則P點的斜距誤差為

圖1 機載SAR幾何關系

設理想軌跡OAi=[xi(ta),yi(ta),zi(ta)]T，實際軌跡OAr=OAi+AiAr，其中AiAr=[Δx(ta),Δy(ta),Δz(ta)]T，表示實際軌跡與理想軌跡的軌跡偏差；P點坐標OP=[x0,y0,-H]T，其中H表示實際軌跡距離地面的高度。代入式(1)得

其中

則

ΔR(ta;x0,y0)=

Ri(ta;x0,y0)

(4)

式中：Δx(ta)?Ri(ta;x0,y0)，Δy(ta)?Ri(ta;x0,y0)，Δz(ta)?Ri(ta;x0,y0)。忽略式(4)中的軌跡偏差的高次項，得

根據圖1中的幾何關系可得

(6)

式中：θ表示瞬時斜視角，β表示APC的入射角。將式(6)代入式(5)，運動誤差可近似表示為

由于三維軌跡偏差中θ和β的存在，使得不同位置目標點的斜距誤差ΔR(ta;x0,y0)具有空變特性，傳統的自聚焦方法不再適用。

1.2 擬極坐標系下SAR成像模型

Bao等[20-24]在擬極坐標系下分析了運動誤差引起的方位相位誤差和距離徙動間的關系，并在圖像譜域推導APE和NsRCM之間的關系，與FFBP算法結合可實現NsRCM補償，在無INS情況下基于回波數據進行運動誤差補償，顯著提高圖像聚焦的質量。但該方法將圖像波數域的相位誤差線性映射到方位時域，這種空間不變的線性近似并不總能成立。

在此只考慮適用于自聚焦的聚束SAR或子孔徑SAR回波數據。點目標P和直角坐標系、擬極坐標系之間的幾何關系如圖2所示。其中：X為方位向，Y為斜距向；Ac為子孔徑APC中心，坐標為[xc,yc]T；Ai為雷達在子孔徑內運動的位置，坐標為[xi(ta),yc]T；P為雷達照射區域內任意一點，直角坐標系坐標為[x0,y0]T，擬極坐標系坐標為[r0,cosθ0]T，r0表示斜距平面內子孔徑中心Ac到目標點P的斜距，cosθ0為子孔徑中心斜視角θ0的余弦。由圖2可知，參數間有如下關系

圖2 斜距平面擬極坐標系和直角坐標系的幾何關系

(8)

當無高精度INS用于運動誤差補償時，Ai即為理想運動軌跡，根據余弦定理點P的理想瞬時斜距為

設Δxi(ta)=xi(ta)-xc，則

點P的瞬時斜距可表示為

設SAR回波為線性調頻信號，經距離壓縮后可表示為

式中：fr=γ·tr表示距離頻率，與距離向快時間tr呈線性關系，γ為常數，表示頻率變化率；Br為信號帶寬；ta為方位向慢時間，tc表示當前子孔徑中心時刻，Ts為當前子孔徑波束照射時間。式(11，12)中相位誤差可表示為

點P的擬極坐標FFBP成像可表示為

將式(11)代入式(14)，整理得

由式(10)可得

將式(16)在r=r0,cosθ=cosθ0處進行一階泰勒級數展開，有

式中：Δφ(Kr,Kcosθ;r0,cosθ0)表示變量代換后的相位誤差；QKcosθ、QKr分別表示Kcosθ和Kr的支撐區域。

根據式(18)，具有運動誤差的圖像的波數譜可表示為

I(Kr,Kcosθ;r0,cosθ0)=

exp[-j(Krr0+Kcosθcosθ0)+

jΔφ(Kr,Kcosθ;r0,cosθ0)]

(19)

文獻[23]提出，當Kcosθ和Kr滿足線性關系時，QPG-FFBP成像算法具備空間不變特性。

2 基于波數子帶劃分的SAR圖像自聚焦方法

2.1 時域到波數域映射分析

根據式(10)，圖像的斜距波數可表示為

SRUAV-SAR系統參數通常滿足Δxi(ta)?r0，忽略高次項近似得

(22)

即斜距波數Kr與距離頻率fr呈線性關系，具有頻率相關性。在相同的方位時間，不同波數下的斜距誤差和相位誤差位于不同的方位波數。

由Kr/Kcosθ≈-1/Δxi(ta)得

Kcosθ≈-KrΔxi(ta)

(23)

式(23)體現了運動誤差的空變特性，即不同目標在同一方位時刻的斜距誤差可能位于不同的方位波數。

假設運動誤差為0，無人機沿理想軌跡飛行，則式(18)可表示為

在式(24)成立的條件下，存在兩個傅里葉變換對，一是斜距r和圖像的斜距波數Kr，二是子孔徑中心斜視角的余弦cosθ和圖像的方位波數Kcosθ。利用FFT快速傅里葉變換，將FFBP成像的結果轉換到方位波數域進行運動誤差校正后，再逆變換回圖像域得到精聚焦的圖像。

實際飛行過程中，運動誤差不可能為0，因此必須討論存在誤差的條件下，式(24)的成立條件，即QKcosθ和QKr。利用成像網格邊緣的像素點推導約束條件[1,24]，得

式中：La為系統合成孔徑長度；λ為SAR系統波長；rmin為QPG的最近斜距。對于正側視聚束SAR，QPG的網格范圍為|sinθ|≤sin(θBW/2)，其中θBW為天線的方位角度。最大值|sinθ|max=sin(θBW/2)，代入式(25)得

(26)

式(26)是式(24)的約束條件，當式(26)成立時QPG可以提供圖像域(Kr和Kcosθ)與相位歷程域(r和cosθ)之間的傅里葉變換關系。

2.2 NsRCM補償

將式(19)中的運動誤差引起的相位誤差Δφ(Kr,Kcosθ;r0,cosθ0)沿方位向和斜距進行向量分解得

式中：Δρ(xi(ta))表示系統運動誤差，是與系統的方位位置xi(ta)有關的函數；Δρr(xi(ta))和Δρcosθ(xi(ta))分別表示Δρ(xi(ta))沿斜距波數Kr和方位波數Kcosθ方向的分量。

將式(23)代入Δxi(ta)=xi(ta)-xc整理得

(28)

將式(28)代入式(27)得

將式(29)的兩部分分別在斜距波數的中心Kr=Kr0進行泰勒級數展開，得

當下我國不動產測繪行業競爭十分激烈，這一嚴峻的形勢使得一些資質不夠的測繪機構遭到淘汰，這不僅順應了不動產測繪市場化發展的趨勢，同時對于整個行業發展而言也起到至關重要的作用。相對比其他測繪工作，不動產測繪工作要求更加嚴格，同時由于測繪結果對于不動產統一登記管理影響較大，所以必須要求測繪機構具備相關的資質和條件，否則不僅會對自身發展造成影響，也會使整個行業的發展造成限制。因此，不動產測繪工作的市場化發展具有必然性特點，只有在激烈的市場競爭中不斷完善和發展，才能確保自身不會被淘汰，同時提高整個行業的發展水平。

(30)

式(30)中其余項都與(Kr-Kr0)有關，可歸為NsRCM引起的相位誤差，記作

其中

根據SAR回波數據的相干性要求，可建立方位相位誤差ΔφAPE(Kr0,Kcosθ)和非系統性距離徙動ΔRNsRCM(Kr,Kcosθ)間的關系式

式(34)用方位相位誤差函數ΔφAPE(Kr0,Kcosθ)及其導數表示非系統性距離徙動ΔRNsRCM(Kr,Kcosθ)，為從一個誤差估計實現兩個誤差的估計奠定理論基礎。此外，式(34)是綜合考慮方位誤差和斜距誤差進行推導所得，可以大大簡化SAR回波數據的運動誤差補償。

2.3 基于波數子帶劃分的SAR圖像自聚焦方法流程

基于斜距波數子帶劃分的QPG-FFBP成像自聚焦方法流程如圖3所示，共分為6個主要步驟。

圖3 算法流程圖

步驟1擬極坐標系成像。使用式(14)對經過距離壓縮后的原始回波數據s(fr,ta;r0,cosθ0)進行擬極坐標系FFBP成像，形成擬極坐標系下SAR圖像i(r,cosθ;r0,cosθ0)。FFBP成像算法作為經典時域成像算法之一，能夠對任意運動軌跡的回波數據處理生成高分辨率的SAR圖像。當沒有GPS和IMU等導航數據滿足運動誤差補償需要時，可借助SAR回波進行運動誤差補償。由式(18)可知，圖像域(Kr和Kcosθ)與相位歷程域(r和cosθ)之間的傅里葉變換關系，借助IFFT將數據轉換到二維波數域i(Kr,Kcosθ;r0,cosθ0)。

步驟2斜距波數子帶劃分。二維波數域SAR圖像i(Kr,Kcosθ;r0,cosθ0)在斜距波數方向被分為幾個子帶信號，以克服空變的斜距波數與方位波數映射關系。每個子帶信號都可近似為一個非空變的信號，同時確保具有殘余RCM誤差的目標能量可以落入同一個斜距單元。為了使用WPGA沿方位波數實現準確的相位誤差估計，子帶寬應足夠小，以確保殘余RCM小于一個斜距單元。子帶劃分后，子帶圖像中的相位誤差可近似為非空變信號

步驟4子帶NsRCM融合。將每個子帶估計的局部NsRCM誤差ΔRNsRCM,i(Kr,Kcosθ)進行融合，得到帶有空變特性的ΔRNsRCM(Kr,Kcosθ)。由式(32)可得子帶NsRCM引起的相位誤差ΔφNsRCM,i(Kr,Kcosθ)與NsRCM之間的關系為

NsRCM引起的總的相位誤差ΔφNsRCM(Kr,Kcosθ)可以通過多個子帶相位誤差的加權平均來計算

(37)

式中WF,i為加權因子。它可根據子帶相位誤差的方差進行計算

(38)

式中D[·]表示方差。

使用估計出的NsRCM相位誤差ΔφNsRCM(Kr,Kcosθ)對圖像整體進行NsRCM校正。

步驟5圖像精聚焦。使用WPGA對圖像進行APE估計，補償殘余APE，實現SAR圖像的精聚焦，得到擬極坐標系中的精聚焦圖像ifine(Kr,Kcosθ;r0,cosθ0)。由于已經補償了NsRCM誤差，因此無需再次進行降采樣。

步驟6坐標轉換。經過運動誤差補償和NsRCM修正后，將擬極坐標系圖像ifine(Kr,Kcosθ;r0,cosθ0)投影到直角坐標系，獲得聚焦良好的SAR圖像ifine(x,y;r0,cosθ0)。

3 實驗結果分析與驗證

3.1 點目標仿真數據處理結果

點目標仿真數據的SAR參數和目標參數設置分別如表1、2所示，滿足式(26)要求，即式(24)中的傅里葉變換對成立，符合PGA使用條件。仿真回波數據中添加任意形式軌跡偏差，包括距離和高度兩個維度的運動誤差，如圖4所示。

圖4 運動誤差結果對比

表1 系統參數

表2 9個目標點的坐標

原始回波數據經過擬極坐標系FFBP成像后(記作QPG)，選用WPGA自聚焦(記作WPGA)、文獻[19]提出的運動誤差補償方法(記作NsRCM+WPGA)作為參考方法，與本文提出的斜距波數子帶劃分自聚焦方法(記作Subband+NsRCM+WPGA)對比，3種方法的運動誤差估計結果如圖4所示。由于WPGA沒有NsRCM校正，故該方法沒有斜距向運動誤差估計結果。相較于參考方法，本文所提方法Subband+NsRCM+WPGA具有最佳的估計精度。

仿真點目標數據成像結果如圖5所示，圖中彩色曲線為擬極坐標系中非線性采樣導致的方位坐標彎曲，故相同方位坐標不同斜距的目標在圖中的位置并不在一條直線上。取中心點A和邊緣點B局部放大。

由于軌跡偏差的存在，圖5(a、b、c)存在嚴重散焦；WPGA自聚焦處理和NsRCM+WPGA處理雖然可以對場景中心點A進行聚焦，但場景邊緣點B仍存在散焦；Subband+NsRCM+WPGA處理能同時提升場景中心點A和邊緣點B的聚焦效果。圖5的中心點A和邊緣點B方位向和斜距向剖面圖如圖6所示，本文所提方法可明顯提升方位向和斜距向的分辨率和峰值旁瓣比，尤其是邊緣點B的方位向性能。

圖5 仿真點目標成像結果

圖6 點目標A、B性能

3.2 實測數據處理結果

2021年9月，在中國石家莊長安區收集了基于小型旋翼無人機的SAR實驗數據，系統工作參數見表1。選取采集數據的合成孔徑時間約為60 s。場景為滹沱河高爾夫球場，光學圖像如圖7(a)所示。由于機載GPS和INS精度較低，無法滿足運動誤差補償需要，故采用接收的回波來估計運動誤差。

由于系統帶寬為0.75 GHz，將圖像在斜距波數域分為0.3 GHz的子帶，確保每個子帶與相鄰子帶有一定重疊。所用方法與仿真數據處理方法相同，處理結果如圖7(b)所示。將圖像中場景A和場景B局部放大，選取特征點，自聚焦成像方法效果對比如圖8所示。圖8(c—f)是圖8(a)中點C的局部放大，可見由于未知的平臺運動軌跡偏差，只用擬極坐標系下FFBP成像，圖8(c)圖像存在散焦；通過WPGA自聚焦處理后，受運動誤差空變性的影響，圖8(d)仍存在散焦；NsRCM+WPGA處理后，同樣受運動誤差空變性影響，圖8(e)也存在散焦；采用Subband+NsRCM+WPGA方法自動聚焦，圖8(f)方位向和斜距向的主瓣能量得到增強，聚焦效果得到提升。圖8(g—j)是圖8(b)中點D的局部放大，相較于圖8(g—i)，圖8(j)的目標主瓣能量更強，方位向分辨更清晰，聚焦效果最好。分別計算各方法處理后C點和D點的峰值旁瓣比(peak sidelobe ratio, PSLR)和積分旁瓣比(integrated sidelobe ratio, ISLR)，如表3所示，進一步驗證了本文所提方法的有效性。

圖7 場景圖像

圖8 擬極坐標系場景局部圖像

表3 圖8中C點和D點性能參數表

4 結論

本文針對小型旋翼無人機平臺的運動特點，提出了一種基于斜距波數子帶劃分的SAR成像自聚焦方法，該方法能在不限制斜視角的情況下完成SAR成像；采用斜距波數多子帶劃分和融合來適應運動誤差的二維空變特性；利用APE和NsRCM誤差的相干性，通過兩次WPGA估計分別校正NsRCM和APE誤差。經仿真和實測數據處理驗證，該方法真實有效，適用于無高精度INS的小型旋翼無人機載SAR自聚焦成像。