基于時頻脊-Radon變換的海面小目標檢測方法

伍僖杰 丁 昊, 2 劉寧波 關 鍵

(1. 海軍航空大學， 山東煙臺 264001; 2. 軍事科學院國防科技創新研究院， 北京 100071)

1 引言

隨著雷達探測環境的復雜化，目標形態的多樣化，海面慢速小目標速度低、散射截面積(Radar Cross Section, RCS)小的特性[1]日益凸顯，導致雷達對弱小目標的檢測愈發成為一個重難點問題。為了更有效地檢測海面弱小目標，研究者們提出了很多方法。文獻[2]從9個方面出發，對這些方法做了分類綜述。其中，從檢測機理層面可以分為能量檢測方法和特征檢測方法。能量檢測方法主要依據海雜波的局部幅度或功率水平信息構造似然比，根據門限因子形成檢測門限并對目標存在與否做出判決[3-4]。該方法對目標的信雜比(Signal-to-Clutter Ratio, SCR)有較高要求，且在海尖峰密集場景下易引起大量虛警，導致實際檢測性能難以達到預期。

特征檢測方法通過充分挖掘海雜波與目標回波特征差異性，將二者從高重疊的觀測空間轉換到低重疊的特征空間，在特征空間內實現目標單元與雜波單元的甄別。與基于嚴格統計模型的能量檢測方法不同，特征提取和融合往往是直觀、經驗、啟發式的[5]，可以根據實際海面環境和雷達系統參數，從雷達回波的幅度、多普勒譜、時頻圖、極化信息等多個維度提取[6]，由于特征信息不再受低階矩或統計模型限制，因此具有更大的自由度，在改善弱小目標檢測性能方面具有更好的潛力。

目前已提出的特征檢測方法有很多種，如時/頻域分形特征檢測方法[7-9]、多普勒譜非廣延熵特征檢測方法[10]、三特征和多特征融合檢測方法[11-15]等。近年來，基于時頻分析的特征檢測方法同樣受到國內外許多專家學者的關注[16-20]，該方法能夠融合應用雷達回波的時域和頻域信息，且可以實現目標能量的相干積累。常用的時頻分析方法包括短時傅里葉變換(Short-Time Fourier Transform, STFT)、分數階傅里葉變換(Fractional Fourier Transform, FRFT)、Wigner-Ville分布(Wigner-Ville Distribution, WVD)、平滑偽Wigner-Ville分布(Smoothed Pesudo-Wigner-Ville Distribution, SPWVD)、小波變換(Wavelet Transform, WT)等，通過參數或非參數化方法在二維平面上提取與目標對應的特征來實現海雜波背景下的目標檢測[6]。例如，文獻[5]提出的特征檢測方法融合應用了歸一化SPWVD中的3個時頻特征，分別為時頻脊累積量(Ridge Integration, RI)、連通區域個數(Number of connected Regions, NR)和最大連通區域尺寸(Maximal Size of connected regions, MS)，經實測數據驗證，將上述特征應用到海面漂浮小目標檢測中，是提高檢測性能的一種有效途徑。

當雷達駐留時間在秒量級或百毫秒量級時，上述特征檢測器大多可在低SCR條件下達到較高的檢測概率。然而，在掃描觀測模式下，目標駐留時間通常難以達到上述量級，隨著積累時間的降低，特征檢測方法性能下滑嚴重。同樣以文獻[5]提出的方法為例，在虛警率為10-3情況下，當積累時間由1.024 s降為0.128 s時，HH極化數據的檢測概率由0.821降為0.629，而VV極化時則由0.789降為0.599。實際上，雷達采用的相干脈沖數通常較少(大多情況下為64個或更少)，可用積累時間基本在幾十毫秒量級，在該時間尺度上，易導致特征不穩定，進而影響海雜波與目標可分性。為提升該情形下雷達的檢測性能，本文提出一種基于時頻脊-Radon變換(Ridges-Radon Transform, RRT)的幀平滑雙特征檢測方法，在常規雜波抑制和時頻域特征提取的基礎上，通過多幀掃描歷史數據和當前幀數據的綜合應用，對特征做幀平滑處理以增強可分性，該方法基本思路為：首先采用分塊白化濾波方法對海雜波進行抑制，在獲得待檢測單元及雜波單元時頻分布的基礎上，通過時頻脊-Radon變換方法提取出時頻脊峰值和頻寬兩個特征，并沿幀時間維進行特征平滑處理，在此基礎上進行特征融合，在雙特征平面上利用凸包算法進行目標檢測。最后，分別采用2級和4級海況實測數據對本文方法的檢測性能進行了分析。

2 分塊白化雜波抑制方法

假設對海雷達在一個波束內發射一系列相干脈沖，并在每個距離單元接收到長度為N的回波時間序列，則海雜波中的目標檢測可以歸結為以下的二元假設檢驗問題：

(1)

式中，x(n)、xr(n)分別代表檢測單元、參考單元信號，s(n)、c(n)分別代表目標信號、海雜波信號，R為參考單元數目。

復合高斯模型是對海雜波進行建模的常用模型，該模型將海雜波分解為兩種分量：散斑分量(Speckle Component, SC)和紋理分量(Texture Component, TC)，其中SC為零均值復高斯過程，TC為非負隨機過程。當雷達的積累時間較短時，由于紋理分量可近似看作常數，復合高斯模型可以簡化為球不變隨機向量(Spherically Invariant Random Vector)模型，自適應檢測器基于該模型利用待檢測單元周圍的參考單元數據估計雜波的協方差矩陣(Sample Co-variance Matrix, SCM)，完成雜波白化。然而，在高分辨海雜波情形下，參考單元數往往小于雷達的相干脈沖數，雜波抑制可能出現較大誤差，解決該問題的一種有效方法是分塊白化[21]。具體流程如下：

首先將待測單元序列x(n)和參考單元序列xr(n)截成長度為L的不重疊向量ui、up,i，如下式：

為保證檢測器性能，參考單元數R與截斷長度L需滿足關系：R≥2L。待檢測單元的雜波協方差矩陣可用參考單元的歸一化樣本協方差矩陣(Normalized SCM, NSCM)估計，即：

(3)

(4)

類似地，對所有距離單元上的數據進行分塊白化雜波抑制，則海面小目標檢測問題可以轉化為雜波白化條件下的目標檢測問題：

(5)

3 雙特征融合檢測方法

3.1 時頻脊-Radon變換提取特征

利用Cohen類時頻分布方法對回波數據進行分析和特征提取。雖然Wigner-Ville分布具有良好的時頻分辨率，但由于其具有非線性特性，在實際應用時會出現嚴重的交叉項干擾，因此本文采用平滑偽Wigner-Ville分布作為時頻分析工具進行處理。分塊白化后回波信號的SPWVD離散形式為：

(6)

定義時頻脊序列Ridges(n)為DSPWVD(t,f)沿時間維的極大值集合，表示為：

Ridges(n)=max(DSPWVD(t,f)|t=t(n))，n=1,2,…,N

(7)

(8)

取平均信雜比取3 dB，速度取1.5 m/s。圖1顯示了白化前后目標單元與雜波單元的典型時頻分布圖像，其中相干脈沖數為64，頻率平滑窗和時間平滑窗分別采用長度5點和31點的Hamming窗，其中時頻脊序列用綠點標識。

由圖1發現，在當前目標運動狀態下，目標單元回波在分塊白化濾波后，時頻脊線向頻率軸兩端擴散程度較輕，基本呈線型；海雜波單元的時頻脊線經分塊白化濾波后，在整個平面的分布無明顯規律。經大量實測數據統計分析可知，白化濾波后目標單元和雜波單元在時頻脊分布特性上差異較大，因此目標檢測問題可以轉化為時頻平面上檢測時頻脊線性程度的問題，本文采用時頻脊-Radon變換方法提取特征以對其線性程度進行定量描述，該變換表示為：

圖1 塊白化前后雜波單元與目標單元的時頻分布圖Fig.1 The SPWVD of the sea clutter cell and target cell before and after block-whitening

RRidges(ρ,θ)=?DRidgesDRidges(t,f)δ(ρ-tcosθ-fsinθ)dtdf

(9)

其中，DRidges(t,f)是從SPWVD中提取的時頻脊圖像，表示為：

(10)

式中，location(Ridges(n))表示時頻脊序列Ridges(n)在時頻平面中的(t,f)坐標位置集合。由于目標單元的時頻脊具有良好的聚集性，近似呈水平線，故在實際檢測中可以通過限制Radon的角度θ以達到減少計算量的目的。

圖2給出了目標單元與海雜波單元的時頻脊-Radon變換圖像，其中，θ∈[80°,100°]。可以看出，目標單元的RRT峰值明顯高于海雜波單元，這主要與目標單元時頻聚集性較強有關。相比之下，海雜波單元的峰值較小。經大量實測數據驗證，這種規律普遍存在。因此，本文將RRT峰值作為檢測器的特征之一，定義為：

圖2 雜波單元與目標單元的時頻脊-Radon變換圖像Fig.2 Images of sea clutter and target cell after RRT

ξMAX=maxRRidges(ρ,θ)}

(11)

式中，ξMAX代表峰值特征，RRidges(ρ,θ)代表時頻脊-Radon變換后的特征空間。

顯然，峰值特征體現了時頻脊能量的集中性，它意味著時頻脊序列Ridges(n)最多有ξMAX個脊點在同一條直線上。但該特征并不足以用于精確描述時頻脊的線性程度，因為它的產生只利用了部分脊點，而未考慮其他脊點的分布特性。下面保持峰值點對應的θMAX角不變，以徑向坐標ρMAX為中心，向兩端擴散，記錄每個徑向坐標上存在的脊點數，可以得到脊點數隨徑向坐標ρ變化的散布圖，如圖3?？梢钥闯?，目標單元的時頻脊散布性較雜波單元要弱，這從另一角度描述了時頻脊能量的集中性。因此為了補足單峰值特征的缺陷，本文定義頻寬特征ξBW用于衡量時頻脊序列在徑向坐標軸ρ上的散布性，如式(12)所示：

圖3 雜波單元與目標單元的脊點散布圖Fig.3 Ridge points distribution of sea clutter and target cell

(12)

式中，COEFF是為了避免個別脊點分布的偶然性而引入的控制量，稱作頻寬系數；P為時頻脊點的總數目(數值上與每一幀的信號長度相等)；ρMAX和θMAX分別是峰值特征ξMAX在Radon空間中的坐標。

顯然，與峰值特征ξMAX不同，頻寬特征ξBW是從宏觀角度對時頻脊的集中性進行描述，它考察了所有脊點在徑向坐標軸上的分布特性，這是只分析某一直線上脊點數的峰值特征無法做到的，而且與峰值特征ξMAX的表現相反，由于時頻脊能量聚集性的差異，目標單元的ξBW特征值較雜波單元小。

3.2 幀間特征平滑方法

圖4 待檢測單元的時頻特征直方圖對比Fig.4 Time-frequency features’ histogram comparison of detection cell

為使特征檢測器在短脈沖數條件下仍能具備較好的目標區分能力，這里通過多幀掃描歷史數據和當前幀數據的綜合應用，對RRT特征做幀平滑處理以增強可分性，即保留多幀掃描歷史數據得到的特征值，并與當前幀特征值進行平滑處理，最終得到的特征提取結果表示為：

(13)

圖5 幀平滑后待檢測單元的時頻特征直方圖對比Fig.5 Time-frequency features’ histogram comparison of detection cell after frame smoothing

3.3 凸包算法檢測分類

構造一個雙特征平面，橫軸代表RRT峰值特征，縱軸代表RRT頻寬特征。由前文分析可知，當待測單元回波中包含目標時，由于在時頻平面上良好的能量聚集性，其峰值特征較大，頻寬特征較小，理應分布在特征平面的右下方；相對的，雜波單元的峰值特征較小，頻寬特征較大，應分布在特征平面的左上方。基于此，本文根據設定的虛警率利用凸包學習算法[5,23]對特征點進行檢測分類，流程如下：

1) 初始化：設置海雜波數據的特征訓練集合為H，個數為W，計算虛警數為L=W·pF，其中pF為設定的虛警率。令l=0。

2) 計算當前數據點的凸包CH(H)，其中CH(H)的頂點為{v1,v2,…,vr}。然后統計落入凸包CH(H)中的特征點數量，設為nall。

3) 設置一個循環變量q，令其從1到r，計算凸包CH(H-{vq})，即從集合H中去除頂點vq，然后計算處于新凸包CH(H-{vq})中的特征點數量，設為nq。

4) 比較集合{nall-n1,…,nall-nr}，去掉該集合中最大數值所對應的頂點vi。

5) 令H-{vi}=H,l+1=l。

6) 若l

綜合以上流程，本文提出一種海面小目標的雙特征檢測器，檢測器原理框圖如圖6所示。在檢測時，先對當前環境下的海雜波數據進行仿真分析，一方面，通過分析形成給定虛警概率條件下的判決空間，另一方面，確定出頻寬系數、平滑因子的適宜取值。然后，分別對待檢測單元數據做分塊白化濾波、時頻脊-Radon變換和特征提取、幀間特征平滑處理，形成待檢測單元的雙特征空間，若落入判決空間內，判定為雜波單元，假設成立；否則為目標單元，假設成立。需要說明的是，當海雜波特征因海況或雷達視角的改變而發生明顯變化時，需要重新采集海雜波數據以確定新的凸包判決區域。本文提出的檢測器雖然歷經的步驟較多，但主要的時間消耗依然在特征提取階段，而由于雷達的相參積累時間較短，極大削減了算法運行所需的時間，因此該檢測器可以滿足雷達實時檢測的要求，具備一定的實用性。

圖6 雙特征融合檢測器框圖Fig.6 Block diagram of double feature fusion detector

4 實測數據分析結果

4.1 低海況數據分析

該部分采用低海況數據對本文檢測器性能進行分析，數據來源于某X波段對海雷達試驗數據，相參體制，脈沖重復頻率(Pulse Repetition Frequency, PRF)4 kHz，帶寬10 MHz，HH極化，海況約為2級。合作目標為漁船，勻速運動且運動速度較慢，圖7是回波數據的距離-時間圖像。經估算，目標SCR約為5～6 dB。在后續分析時，假定每幀數據均包含64個脈沖，對應的相干處理時間為16 ms。

圖7 某型X波段雷達雷達接收數據的距離-時間圖像Fig.7 The range-time image of data received by an X-band radar

4.1.1檢測器性能的影響因子分析

首先抽取該數據集中的部分海雜波數據添加仿真目標(以1.5 m/s的速度勻速運動)，分析虛警率為10-3時頻寬系數對本文檢測器性能的影響，如圖8(a)。顯然，頻寬系數在0.75～0.85范圍內時，檢測器性能較好，圖8(b)繼續對該區間內的檢測器性能作了詳細分析，發現在該海況條件下較為合理的頻寬系數區間為[0.79, 0.81]，故本文取COEFF=0.8。

事實上，經大量數據集驗證，該頻寬系數的取值具有較強的普適性。經理論分析可知，頻寬特征是從宏觀上對時頻脊的集中性進行描述，其準確性依賴于所使用的脊點數，當COEFF值過小時，頻寬特征的形成只使用了少部分脊點，而未考慮其他可能有效的脊點，必然導致其準確性不足，使目標單元與海雜波單元之間的可分性不夠，檢測器性能下滑，如圖8(a)中COEFF∈[0.2,0.6]區間；反之，若COEFF值過大，目標單元的時頻譜中也可能存在個別脊點嚴重偏離直線(θMAX,ρMAX)，同樣會引起檢測器性能下滑，如圖8(a)中COEFF∈[0.9,1]。因此后文不再對頻寬系數作分析，均取定COEFF=0.8。

圖8 雙特征檢測器性能隨頻寬系數的變化曲線Fig.8 The performance curve of the double feature detector varies with the bandwidth coefficient

進一步對平滑因子的取值影響作分析。分別取仿真目標的速度為1.5 m/s、0 m/s，圖9給出了平滑因子對檢測性能的影響。顯然，隨著M值增大，檢測器對海上目標的檢測性能呈提升趨勢，且在[5,30]區間內改觀尤為明顯。當M∈[30,50]時，檢測器性能的提升程度非常小，說明幀間特征平滑方法對檢測器性能的改善是有限的。考慮到在實際場景中目標還有可能出現跨距離單元走動現象，進而影響能量累積和特征可分性，因此在兼顧檢測性能和處理復雜度的情況下，取平滑因子M=30。

圖9 平滑因子對檢測器性能的影響Fig.9 The effect of smoothing factor on detector performance

4.1.2檢測器的檢測結果及對比分析

確定好頻寬系數、平滑因子的取值之后，利用前文所述方法提取出雷達數據的RRT峰值和頻寬特征，并在雙特征平面上利用凸包算法進行檢測分類。圖10分別展示了本文提出的檢測器、時頻三特征檢測器[5]、頻域CFAR檢測器、分形檢測器在虛警率為10-3時的檢測結果圖，可以看出，本文方法得到的檢測結果較為連續、穩定，且虛警較少，在性能上優于已有方法。

圖10 四類檢測器的檢測結果Fig.10 The detection results of four detectors

進一步對檢測性能進行量化分析，結果如圖11所示，表1展示了其中三種虛警率條件下的檢測概率，虛警率最高者用粗體標識。顯然，四類檢測器性能均隨著虛警率增加逐步提升；同時經過對比發現，本文提出的檢測器性能在虛警率pF∈[10-3,0]時要優于其他三類檢測器，且由表1可知，性能差距超過20%，這是本文綜合運用雷達歷史掃描數據的做法體現出的優勢；時頻三特征檢測器與頻域CFAR檢測器性能在不同虛警率時各有優劣，而分形檢測器在低虛警率檢測概率接近于零，無法滿足目標檢測的需要。

表1 四類檢測器對低海況數據的檢測概率

圖11 四類檢測器的性能對比Fig.11 Performance comparison of four detectors

4.2 高海況數據分析

本節繼續采用高海況數據對雙特征檢測器性能進行分析，數據來源于CSIR數據庫，該數據庫由Fynmeet雷達于2006年在南非西南海岸采集，該雷達為VV極化，相參體制。本文采用文件名為TFA17-007的數據進行分析，對應的雷達工作參數和環境參數如表2所示[24-26]，海況約為4級。圖12給出了該組數據的距離-時間圖像，其中目標在第24距離單元勻速運動。經估算，目標SCR約為3 dB。在后續分析時，同樣假定每幀數據均包含64個脈沖，對應的相干處理時間為12.8 ms，虛警率為10-3。

圖12 TFA17-007數據的距離-時間圖像Fig.12 The range-time image of TFA17-007 dataset

表2 雷達性能及環境參數

與4.1.1節類似，首先抽取檢測當時的海雜波

數據，添加勻速運動(1.5 m/s)的仿真目標進行分析，結果如表3、圖13所示。顯然，此時由于海況等級較高、海浪起伏周期長的影響，平滑因子對檢測器性能的改善程度并不如圖9明顯，且存在平滑因子、信雜比增大時檢測概率反而下降的異?，F象。進一步分析可知，平滑因子M=50、60、70時，本文檢測器的性能在信雜比scr∈[0,6]區間內相近，由此推斷，當M>50時，平滑因子對檢測器性能的改善是微弱的。而當0

表3 本文檢測器的檢測概率

圖13 平滑因子對檢測器性能的影響Fig.13 The effect of smoothing factor on detector performance

在圖14中，依舊通過量化分析給出了四類檢測器的ROC曲線，且用表4展示三種虛警率條件下的檢測概率。可以看出，在高海況條件下，時頻三特征檢測器性能不佳，說明時頻特征容易受到海況條件的影響，且在低虛警率范圍內該檢測器與本文檢測器的檢測概率均無明顯變化，這是因為部分屬于目標單元的特征點落入凸包判決空間的位置較深，只有當虛警率較大時，凸包學習算法形成的判決空間較小，才能將它們正確分類。而頻域CAFR檢測器、分形檢測器分別因為檢測機理及判決方式的不同，其檢測性能均隨虛警率增加穩步上升，未出現此種現象，但在虛警率較低時，兩者檢測概率不足10%，無法用于檢測目標。進一步分析表4可知，由于平滑因子的作用，本文檢測器相對于其他三類檢測器的性能改善程度在pF=10-3時將近60%，能夠較好滿足目標檢測的需求。

圖14 四類檢測器的性能對比Fig.14 Performance comparison of four detectors

表4 四類檢測器的檢測概率

5 結論

針對相干脈沖數較少條件下特征檢測器性能損失嚴重的問題，本文提出一種基于時頻脊-Radon變換的幀平滑雙特征檢測方法。該方法采用分塊白化濾波實現初步的海雜波抑制和SCR提升，以SPWVD為時頻分析手段，對時頻脊做Radon變換并提取出時頻脊峰值和頻寬兩個特征。為使特征檢測器在短脈沖數條件下仍能具備較好的目標區分能力，本文通過多幀掃描歷史數據和當前幀數據的綜合應用，對RRT特征做幀平滑處理以增強可分性。最后，在雙特征平面上利用凸包算法形成判決空間并實現目標檢測。經2級、4級海況實測數據驗證，當相干脈沖數僅有64個的條件下，本文檢測器的檢測性能明顯優于頻域CFAR檢測器以及已有的兩類特征檢測器(分別是時頻三特征檢測器和分形檢測器)。后續還需進一步深入研究不同環境因素和不同目標運動狀態下雙特征檢測器參數優化選取問題。