GST時頻濾波SAR射頻干擾抑制算法

賀音光， 譚小敏， 楊娟娟， 龔紫翼， 閆 偉

(西安空間無線電技術研究所， 陜西西安 710100)

0 引言

合成孔徑雷達(Synthetic Aperture Radar，SAR)是一種主動微波成像遙感設備，在工作過程中容易受到射頻干擾[1-2](Radio Frequency Interference，RFI)的影響，使得接收回波的信干噪比(Signal to Interference plus Noise Ratio，SINR)降低，嚴重時會影響SAR成像質量及后續圖像解譯，因此在回波中對RFI信號進行有效抑制具有重要意義[3-4]。

SAR射頻干擾抑制方法的研究是一項重要的課題。在缺乏先驗信息和參數模型的情況下，非參數化方法可在時域或頻域實現回波數據中干擾的自適應抑制，與參數化方法相比有較大優勢。頻域陷波法[5]和子空間濾波法[6]是兩種經典的非參數化方法，它們處理流程簡單，不足是忽視了干擾的非平穩局部特性，而且頻域陷波法會完全損失與干擾信號處于同頻帶內的有用信號，從而造成后續成像質量的下降。

本文提出一種基于GST時頻濾波的SAR射頻干擾抑制算法，有效利用了時頻分布關注干擾的非平穩局部特性，在有效抑制干擾的同時避免有用信號的損失。該方法首先利用配對樣本T檢驗[7]進行干擾的方位檢測，接著分別處理包含干擾數據的實部和虛部：進行GST[8]，得到其時頻分布，在此基礎上，利用子空間濾波法去除每一個時間窗內數據中的干擾分量，逆GST后得到干擾抑制后的數據。最后，實驗結果驗證了所提方法的有效性。

1 SAR射頻干擾的特性及模型

SAR系統受到射頻信號干擾的接收回波x(t)有如下形式：

x(t)=s(t)+rfi(t)+noi(t)

(1)

式中，s(t)為目標回波，rfi(t)為射頻干擾，noi(t)為加性噪聲，近似為高斯白噪聲。

實測大場景SAR回波數據的統計分析表明：

1) 不含干擾的SAR目標回波信號近似服從于復高斯分布[9]；

2) RFI具有比較復雜的時頻特性：有窄帶性和非平穩性；

3) RFI從源到SAR接收通道主要是單程傳播，干擾功率一般大于SAR接收回波功率。

基于上述特點，目標回波和噪聲信號可以被統一視作寬帶高斯白噪聲信號，所以SAR射頻干擾抑制問題就轉化為在高斯白噪聲背景下的窄帶干擾抑制問題。因此，式(1)可表示為

x(t)=sn(t)+rfi(t)

(2)

式中，sn(t)=s(t)+n(t)。與x(t)相對應的離散表達式為

x(n)=s(n)+rfi(n)+noi(n)=

sn(n)+rfi(n)

(3)

式中，n=1,2,…，Nr，Nr為距離向采樣點數，上式中的元素均是式(1)和式(2)的離散形式。

對于窄帶RFI信號，具有類似線性調頻(LFM)信號的形式：

(4)

式中，Bk，fk和γk分別表示第k個干擾分量的幅度、中心頻率及調頻率。

2 GST時頻濾波干擾抑制原理

2.1 配對樣本T檢驗原理

目前合成孔徑雷達RFI抑制算法的處理對象多是快時間距離向回波數據。然而對全部接收回波進行干擾抑制，在理論上可以極大地使算法發揮作用，但會增加無意義運算，降低信號處理效率，所以有必要進行RFI檢測，標記出不存在干擾的方位。文獻[10]指出快時間回波數據之間具有相關性，并提出基于相關系數檢測的方法，然而實際純凈回波的相關程度較難判斷，無法按照常規設置門限進行檢測。實際上，這種相關性也體現在純凈回波信號期望值μ小的起伏程度；同時，純凈SAR回波的實部與虛部均近似服從高斯分布。通過以上分析，本文提出了一種基于配對樣本T檢驗的方位向干擾檢測方法。其步驟為：

1) 計算所有距離向快時間數據的功率并進行比較，考慮受RFI影響的回波具有相對大的功率，進一步地選取功率值較小的快時間信號xk(n)作為不含干擾的參考信號，其中k為方位向標號。

2) 將參考信號xk(n)分別與其他所有距離快時間信號xs(n)進行配對樣本T檢驗，其中s= 1,…,Na，且s≠k，Na為最大方位向標號。第一步，建立原假設和備擇假設，分別記為H0:μk=μs和H1:μk≠μs，確定顯著性水平參數α，為了防止包含弱干擾方位的漏檢，實際中取α=0.1可以滿足要求。第二步，按照以下公式計算統計量t′值：

(5)

該算法利用了SAR一維距離向時域回波的相關性及起伏不大的特性，通過計算統計量t′，并比較其顯著性以判斷干擾的存在，計算量小，實用性強。

2.2 RFI抑制流程

假設RFI檢測后被標記的某條距離向復信號為xi(n)，它由實部和虛部組成：

xi(n)=ai(n)ejφi(n)=pi(n)+jqi(n)

(6)

(7)

φi(n)=arctg(qi(n)/pi(n))

(8)

式中，n=0,1,…,Nr-1，ai(n)為幅度，φi(n)為相位，pi(n)為實部，qi(n)為虛部，且pi(n)和qi(n)均為實數。

由于RFI信號具有比SAR回波高的功率及具有未知的相位，會同時影響接收回波的幅度與相位，使得接收回波信號的實部和虛部電平被抬高；同時為有效減少處理后可能存在的殘余相位與相位畸變，分別處理實部和虛部可作為一種選擇。

GST是一種用于分析非平穩信號的有效工具，高頻處有高的時間分辨率，通過改變調節因子可以獲得高的頻率分辨率，高的時頻分辨率可令基于該變換域的干擾抑制處理精細化，即更大程度地保留與干擾信號譜圖不重疊的有用信號；同時，子空間濾波能有效地分離混合信號，降低與干擾信號譜圖重疊的有用信號損失。

2.2.1GST及窗函數調節因子

廣義S變換是S變換的推廣，通過引入窗函數調節因子緩解了窗口函數寬度減小而導致的頻率分辨率降低問題，從而可以獲得理想的二維分辨率。

由小波變換可以推出如式(9)的一維連續信號x(t)的S變換表達式，選用的小波變換基函數為式(10)的高斯窗函數，滿足積分為1。其中，ψ為小波基函數，τ為時移因子，f為頻率。

(9)

(10)

對于式(9),引入窗函數調節因子λ,p，得到GST：

e-j2πftdt

(11)

通過調節兩個因子可以調節二維時頻分辨率。當且僅當λ=1,p=1時，GST退化為S變換；當調節λ<1或p<1時，可以使得頻率分辨率提高而時間分辨率下降；當調節λ>1或p>1時，可以使得時間分辨率提高而頻率分辨率下降。一般來說，兩個因子取值不宜過大或過小，應取在1的附近。文獻[11]指出p值的改變對窗函數窗寬、幅值、窗脊的變化具有決定性作用，而λ對窗的形態改變遠不及p值，在實際應用過程中可以起到對窗口微調的作用。

為了對(λ,p)有效取值，利用反映信號在二維域聚集特性(二維分辨率)兩種度量[12]χ1和χ2對不同(λ,p)組合的GST時頻分辨率進行評估。假設待處理的SAR一維距離向數據實部為pi(n)，其經GST得到的時頻域分布為GSi(該符號代表GST的離散形式)，GSi(m,v)為矩陣中的元素，m和v分別為頻率和時間標號，則兩個度量可以表示為

(12)

(13)

χ1和χ2越大表示能量越集中。

為選取合適的因子，一般需要進行二維尋優，并根據兩個度量來評價因子選取的合適性，這需要利用每一組(λ,p)進行GST，這必然會降低處理效率。實際過程中，提出利用以下步驟選取(λ,p)：

1)λ取1，在p的選擇區間0.5～1.5中，計算度量χ1和χ2，選擇使得兩個度量最大時對應的p值(一般取p<1，以達到良好的頻率分辨率)；

2) 確定p后，微調λ，補償時間分辨率(一般如果p<1，λ選擇略大于1)。

2.2.2GST時頻濾波干擾抑制算法

針對被標記的一維距離向信號的實部pi(n)，調節(λ,p)(也可以統一預設)，將pi(n)進行GST得到GSi，列向量為各時間窗內的數據譜。GSi按列作IFFT，記gsv為其第v個時間窗內的數據向量，即gsv=[gsv(1),gsv(2),…,gsv(L)]Τ，L為數據向量長度，(·)Τ表示轉置。接下來，利用子空間濾波對每個時間窗內數據進行干擾抑制處理：

1)gsv延遲嵌套構造Hankel矩陣Gv：

(14)

H和J為矩陣的行數與列數，且滿足：J=L-H+1和(L-H+1)≥2H。

2) 構造Rv并特征分解：

Rv=GvGvΗ=UΛUΗ

(15)

式中，Λ=diag[λ1,λ2,…,λH]是H維對角陣，λ表示特征值且λ1≥λ2≥…≥λH，對應特征向量排成矩陣為U=[u1,u2,…,uH]。此時的干擾子空間由代表干擾信號的前r個大特征值對應的特征向量Ur張成，進一步地，Gr=UrUrΗGv表示Gv向干擾子空間投影后得到的干擾信號矩陣。大特征值個數r利用AIC準則估計：

(16)

式中，LLF(r)和PAIC(H,r)分別為

LLF(r)=(H-r)·

(17)

PAIC(H,r)=r(2H-r)

(18)

3) 從Gr中重建一維干擾信號序列gsjv(l)：

(19)

4) 依據式(20)將gsjv(l)從gsv(l)中減去，完成對一個時間窗內數據的干擾抑制處理：

gssv(l)=gsv(l)-gsjv(l),l=1,2,…,L

(20)

(21)

對所有被標記受干擾影響的SAR回波信號進行上述操作，即可完成SAR回波的干擾抑制。圖1是整個干擾抑制流程的框圖。

圖1 GST時頻濾波干擾抑制流程

3 實驗結果與分析

3.1 基于SAR一維距離向回波的實驗

本實驗中，引入兩個指標來定量評估干擾抑制方法：干擾抑制比(ISR)和信號失真度(SDR)[13]。ISR表示干擾抑制的程度。SDR表示干擾抑制后對有用信號的損失程度。實驗數據來自Radarsat-1，干擾為人為添加，部分參數如表1所示。

表1 實驗的參數設置

圖2(a)和2(b)分別為原始距離向回波頻譜和時頻譜。在數據中加入干擾后，兩種譜圖如圖2(c)和2(d)所示，圖中的干擾具有窄帶及非平穩的特性，其中的非平穩性體現在：其一，干擾信號的持續時間可能小于回波的持續時間；其二，干擾的頻率隨時間的變化。采用頻域陷波法后的結果如圖2(e)和2(f)所示，從頻譜圖可見存在干擾的頻帶被置零，在抑制干擾的同時嚴重損失了該頻帶的有用回波信號，在時頻譜圖中，損失了所有時間窗內干擾頻點處的信號。圖2(g)和2(h)是采用本文方法處理后的結果，選取(p,λ)=(0.5,1.2)，與頻域陷波法相比，該方法有效利用了干擾的窄帶及非平穩特性，即對混合信號GST瞬時譜濾波，所以在抑制干擾的同時，對不含干擾的瞬時譜和與干擾重疊的有用信號譜的損失很小。同時，表2中干擾抑制性能指標的比較表明，本文方法的ISR與頻域陷波法相當，但具有比頻域陷波法更低的SDR，從而表明其性能優于頻域陷波法。

圖2 距離向回波的干擾抑制分析

表2 距離向回波干擾抑制后的結果

3.2 基于SAR回波的成像實驗

為有效定量評估利用干擾抑制處理后的數據成像的結果，同樣利用文獻[13]中所采取的兩個指標：對比度D和信噪比γ。兩個指標值的大小與干擾抑制方法的有效性高低呈正相關。

如圖3(a)為受到干擾影響的SAR回波數據成像結果(橫向為距離向，縱向為方位向)，可見RFI在SAR圖像中的表現形式是沿著距離向密集的亮紋，該亮紋遮蓋了圖中的目標繼而影響了圖像的判讀，使得后續基于SAR圖像的各種應用變得困難。圖3(b)為采用頻域陷波法后的成像結果，RFI對成像的影響有所降低，干擾覆蓋場景的輪廓可以分辨，但是仍然丟失了圖像細節信息，原因是部分干擾的殘留和處理造成的有用信號損失。圖3(c)為利用本文所提方法得到的結果，場景輪廓很清晰，而且更多的細節信息被保留，成像質量明顯優于圖3(b)。同時表3中的圖像對比度和信噪比指標的定量對比也證明了本文方法的有效性。

(c) 本文所提方法處理的圖像

表3 干擾抑制后的成像結果

4 結束語

針對SAR射頻干擾抑制，本文提出一種基于GST時頻濾波的方法。該方法在處理前首先進行RFI方位檢測，提高了處理效率；同時在GST時頻域進行子空間濾波處理，與頻域陷波法相比，在抑制干擾的同時，盡可能地保留了有用回波信號，進而提高了SAR圖像質量。實驗的干擾抑制結果均表明所提方法優于頻域陷波法。