張朝星，趙興剛，韓曉軍

(1.解放軍32137部隊，河北 張家口 075000；2.解放軍66136部隊，北京 100042；3.解放軍31675部隊，河北 張家口 075000)

0 引 言

時頻分析[1]通過將一維時域信號轉換為二維時頻域圖像，從而提取有效的圖像特征來判斷有無目標，是實現(xiàn)海雜波背景下目標檢測[2]的一種有效方法。短時傅里葉變換(STFT)作為一種常用的時頻分布[3]，其基本思想是用一個確定寬度的窗函數(shù)沿時間軸逐段對信號進行傅里葉分析，得到信號的譜圖。STFT具有許多優(yōu)點，然而也有其缺點和不足，根據(jù)不確定性原理[4]，其時間分辨率與頻率分辨率的乘積恒定，不能同時提高，并且窗的寬度不能根據(jù)不同頻率實現(xiàn)自適應調整，時頻分辨力較低，導致使用STFT方法對海雜波背景下的弱小目標進行檢測時效果較差。

曲線擬合選擇一條合適的曲線來最佳地擬合觀測數(shù)據(jù)，在歐式空間中，常見的曲線類型有二維直線、拋物線、橢圓等。在功率譜的信息幾何理論[5]中，所有功率譜的集合可以看做一個黎曼流形[6]，則對回波數(shù)據(jù)進行STFT處理后所得到的不同時間點的功率譜就對應了流形上的若干點。類比于歐式空間中的曲線擬合，本文提出了一種使用流形上的距離和測地線[7]對不同時間點的功率譜進行曲線擬合以提高其時頻分辨力的方法。由于STFT是一種線性變換，我們選擇歐式空間中的直線在流形上的推廣(即測地線)對不同點的功率譜進行擬合。在具體的實現(xiàn)過程中，歐式空間使得給定數(shù)據(jù)與所選擇曲線對應點的歐式距離總的平方和最小(最小二乘法)[8]。轉換到流形上以后，實現(xiàn)方法不變，但歐式距離要替換為兩功率譜間的距離。我們選擇M-K運輸問題[9]中的運輸距離作為兩功率譜之間的差異度量，得到了相應的黎曼度規(guī)，對STFT后的譜圖實現(xiàn)了測地線擬合。擬合后的譜圖對于單一點的功率譜進行了重構，使其更加平滑和準確。其次擬合所用測地線對不同時間點的功率譜實現(xiàn)了時域上的平滑，時間分辨率也得到提高，進而整體提高了譜圖的時頻分辨力，STFT的檢測性能得到改善。

1 基于信息幾何的功率譜曲線擬合

1.1 運輸距離及測地線

在上節(jié)中已提到，具體實現(xiàn)流形上的擬合時，要將最小二乘法中的歐式距離替換為2個不同時間點功率譜之間的距離，且該距離除了滿足對稱性和三角不等式等距離的一般性質外，還必須要滿足弱連續(xù)性[10]，以使得功率譜能夠收斂，即隨著2個功率譜的“靠近”或“遠離”，該距離能相應地減小或增大。因此，我們選用了M-K運輸問題中的運輸距離來度量兩功率譜之間的差異，對于2個歸一化的密度函數(shù)f0和f1，其運輸距離定義如下[9]：

(1)

式中：X=[-π,π]；ψ(x)為一個f0和f1之間的保測度映射，即：|detψ(x)|f1(ψ(x))=f0(x)。

在一般情況下，ψ(x)可通過下式計算得到[9]：

F0(θ)=F1(ψ(θ))

(2)

(3)

當f1是f0的一個小的增量時，可由運輸距離得到一個黎曼度規(guī)，并稱之為運輸度規(guī)[11]：

(4)

(5)

定義了流形上的度規(guī)，就可得到f0和f1之間由梯度流決定的測地線[12]：

Fτ((1-τ)θ+τψ(θ))=F0(θ)

(6)

式中：ψ(θ)可以通過式(2)計算；Fτ為fτ的積累函數(shù)。

1.2 流形上的曲線擬合

有了2個功率譜之間的距離和測地線，我們就可以在流形上對功率譜進行曲線擬合。假設給定一個功率譜序列：G={gτi(θ):θ∈[-π,π],i=0,1,…,n}，其中τi是一個指代時間的歸一化序列，且τ0=0,τn=1。這些功率譜通?？梢酝ㄟ^對觀測數(shù)據(jù)進行短時傅里葉變換(STFT)得到，τi(i=0,1,…,n)可看做STFT中對應時間窗的中點。令d代表兩功率譜之間的距離，則曲線擬合的過程就是尋找一條測地線fτ,τ∈[0,1]，使得下式取得最小值：

(7)

利用運輸距離及其相應的測地線來實現(xiàn)該過程，則整個測地線擬合的目標就是要確定一條測地線fτ,τ∈[0,1]，使得下式取得最小值：

(8)

任意一條測地線都可以完全由兩“點”確定，在這里就是指f0和f1。另外，它也可以根據(jù)式(2)和式(6)中的Ψ確定。f0,f1和Ψ是根據(jù)STFT處理后得到的功率譜G來確定的，根據(jù)式(3),有：

(9)

根據(jù)式(6)，由ψ(θ)可得f0和f1之間測地線fτ的具體表達式為：

F0(θ)=F1(ψ(θ))=Fτi((1-τi)θ+τiψ(θ))

(10)

(11)

將式(11)代入式(9)，則目標函數(shù)可寫為：

(12)

(13)

(14)

所以Jg(f0,f1)可以通過以下的求和式進行近似：

(15)

所以，功率譜的測地線擬合就轉化為了以式(15)為目標函數(shù)的最優(yōu)化問題，其線性約束為：

(16)

圖1 功率譜測地線擬合示意圖

2 仿真實驗

本節(jié)采用加拿大McMaster大學提供的IPIX雷達#320號數(shù)據(jù)[14]對本文所提方法的檢測性能進行驗證，該數(shù)據(jù)文件包含14個距離單元，每個單元由131 072個采樣樣本構成，其中第7個單元為主目標單元，第6、8、9單元為次目標單元，其余為純雜波單元。

首先將第1個距離單元的采樣數(shù)據(jù)作為雜波背景，通過STFT分析海雜波本身的譜特性。取觀測數(shù)據(jù)樣本長度為256，窗函數(shù)寬度為32，時間窗重疊區(qū)域寬度為16，脈沖重復頻率fr=1 000 Hz，可得第1個距離單元的海雜波二維譜圖，如圖2所示。

圖2 第1個距離單元的海雜波二維譜圖

從圖2中可以看出，海雜波譜是非對稱的，由于海浪的移動，使得譜的主瓣偏離零頻，雜波能量主要集中在[-200 Hz,0 Hz]這個區(qū)間內，下面將仿真目標信號加到第1個單元的純海雜波數(shù)據(jù)中，設雷達在該距離單元觀測的復包絡信號為：

x(n)=s(n)+v(n),n=1,…,N

(17)

式中：s(n)=aej(2πfdn/fr)，a為信號幅度，fd為目標多普勒頻率。

令信噪比分別取0 dB和4 dB,目標多普勒頻率分別取-100 Hz(主雜波譜區(qū))和100 Hz(非主雜波譜區(qū)),分別給出相應的STFT變換后的譜圖，如圖3所示。

圖3 海雜波加仿真信號后STFT三維譜圖

從圖3(a)和圖3(c)中可以看出，當目標頻率處于非主雜波區(qū)(fd=100 Hz)時,可以非常輕易地從STFT變換后的時頻圖像中判斷目標的存在。但當目標頻率處于主雜波區(qū)(fd=-100 Hz)時，如圖3(b)和圖3(d)所示，信噪比較高時，還可以從譜圖中對目標進行辨別；但在低信噪比(SNR=0 dB)時，目標被雜波淹沒，受譜圖分辨力的限制，很難判斷目標的存在。

然后就可以用Matlab中的最優(yōu)化工具箱對該問題進行求解。對應圖2中的4種情況，擬合得到的4條測地線的2個端點如圖4所示。

圖4 擬合得到的測地線端點圖

圖5 擬合處理后的譜圖

從圖4中單個時間點的功率譜已經可以初步判斷各個情況下目標的存在，但當目標處于主雜波區(qū)時，還不夠明顯。根據(jù)式(11)，有了2個端點就可以確定一條測地線，進而就得到了測地線上各個時刻的功率譜，如圖5所示。

從圖5中可以看出，當目標處于主雜波區(qū)(fd=-100 Hz)并且信噪比較低時(SNR=0 dB)，與圖3(d)相比，圖5(d)中可以明顯地看到一條凸出的“山脊”，其所在頻率位置與目標多普勒頻率相近，進而就可以較容易地判斷目標的存在。

以上結果都是將仿真目標信號加到純雜波數(shù)據(jù)上進行處理的，下面直接對實際有目標的第7個距離單元數(shù)據(jù)進行STFT變換，并進一步進行測地線擬合，給出處理結果如圖6所示。

圖6 第7個單元(主目標單元)STFT+測地線擬合處理后的譜圖

IPIX實測數(shù)據(jù)是由雷達對漂浮于海面直徑約1 m的金屬球持續(xù)照射得到的，該金屬球會隨著海浪的波動進行漂移，但其運動速度與海浪運動速度相比相對較小，因而體現(xiàn)在功率譜上即其能量主要集中在零頻附近，且其運動是無規(guī)律的，目標RCS和速度都會隨著時間推移而發(fā)生變化，因而其功率譜也會隨時間變化。從圖5中可以非常明顯地看到，在零頻附近有一條凸出的脊線，而且這條脊線并不像加仿真目標信號時形狀比較規(guī)則，而是在時域上有一定起伏和變化，這與實際目標狀況是相符的，因而就可以判斷目標的存在。

3 結束語

本文根據(jù)功率譜的信息幾何理論，提出利用流形上的距離和測地線對STFT變換后的不同時間點的功率譜實現(xiàn)曲線擬合，以提高其時頻分辨力；并通過推導,將該擬合過程轉換為可用數(shù)值方法解決的具有線性約束的最優(yōu)化問題，使其容易求解；最后在仿真實驗中，分別利用IPIX實測純海雜波數(shù)據(jù)加仿真目標信號和本身帶目標的海雜波數(shù)據(jù)2種方式對所提方法的性能進行了測試，驗證了本文方法的有效性。