基于偏振分析的瞬態干擾抑制方法

劉 偉 樸勝春 祝捍皓

(哈爾濱工程大學水聲技術重點實驗室 哈爾濱 150001)

1 引言

矢量水聽器可以同步共點測量聲場標量和矢量信息，為后置信號處理提供更多的信息。它的諸多優點使其一經出現，便受到了國內外學者的普遍關注，已廣泛應用于海洋開發及研究中[1,2]。通常情況下，實驗海區并非安靜海域，漁民作業、過往船只和天氣等因素都會導致接收信號中不可避免地存在一些瞬態干擾。這類干擾具有持續時間短、強度高和頻帶寬的特點，直接影響信號處理的性能，因此必須對瞬態干擾進行抑制。

目前在雷達信號處理中已有很多方法可以抑制瞬態干擾[3-8]，文獻[3]利用復數據經驗模式分解實現瞬態干擾的檢測、抑制和數據重構。文獻[4-7]都是先找到瞬態干擾位置，將瞬態干擾去除，再利用文獻[8]中所描述的線性預測方法恢復信號。文獻[4]利用分形檢測，提取信號的相似性特征，找到瞬態干擾位置；文獻[5]利用高分辨算法TLS-ESPRIT檢測干擾發生位置；文獻[6]應用信號子空間檢測, 實時判決信號的奇異值，找到干擾位置；文獻[7]利用線性預測AR模型遞歸最小二乘(RLS)求解零點后聚類估計瞬態干擾的位置。這些方法均是針對雷達信號特點而提出的，算法相對較復雜。

本文利用矢量水聽器的特點，提出一種新的瞬態干擾消除方法。首先對接收的矢量信號進行偏振特性分析，獲得偏振參數；再根據信號和瞬態干擾在偏振狀態上區別，設計極化濾波器，對信號進行濾波，去除瞬態干擾；最后通過信號的偏振模型恢復瞬態干擾處數據。

2 矢量信號的偏振參數提取

偏振的概念最早起源于光學，現已廣泛應用于地震學和電磁學等領域[9-12]。偏振分析又稱極化分析，是矢量信號的共有屬性，用于描述質點運動軌跡[12]。水中傳播的聲波為縱波，當入射信號中僅存在直達聲時，接收信號的空間質點運動軌跡是直線；由于多途結構的存在，入射聲波為多個同頻平面波疊加，合成的振速矢量在此情況下是旋轉矢量，因此，在流體和氣體聲場中，偏振的概念也成立，并可用偏振分析將橢圓的屬性表示為定量值[2]。

將橢圓軌跡上樣點看作一批空間點的集合，通過統計的方法獲得橢圓的屬性。偏振分析常用的方法為特征分解、奇異值分解和時頻分析方法，本文采用特征分解方法。設t時刻接收矢量信號為S(t)=[vx(t)vy(t)vz(t)]T+N(t),T＞t＞0,N(t)是均值為0且功率為σ2的高斯白噪聲。由此獲得信號的協方差矩陣[13]

若信號為線性偏振，橢圓率可能無限大，故常用直線性描述信號偏振狀態[13]，即

時間窗長度的選擇對偏振參數的獲取產生一定影響。若時間窗過長，一個時間窗內存在幾種偏振狀態，偏振分析不能反映質點運動軌跡的真實情況；若時間窗長度過小，信號的統計特性變差，協方差矩陣易受噪聲影響。時間窗長度的選擇可采用“三次零點法”。將信號看作“+”，“0”，“-”組成的序列，一個時間窗內“0”最多出現 3次，“+”或“-”出現次數不受限制。若發射信號為單頻信號，3個“0”符號即組成一個信號周期，文章中采用信號周期長度作為時間窗長度。

3 基于偏振分析的瞬態干擾抑制方法

3.1 信號模型

若信號中不存在噪聲和瞬態干擾時，多途環境下目標信號可以化簡為

其 中S1(t)=ΣAkcosφkcos(ωt-ωτk),S2(t)=ΣAksinφkcos(ωt-ωτk),θ0為目標方位角，Ak,τk和φk為每個途徑的幅度、時延、俯仰角及多途距離。通過推導可知，多途環境下，接收信號質點運動為具有固定偏振參數的橢圓。

實際接收信號中不僅存在目標信號，還混有噪聲和瞬態干擾，即

式中i(t)為瞬態噪聲。無瞬態干擾時質點運動軌跡為具有小抖動的橢圓，橢圓的主軸方向在小角度內搖擺。當瞬態干擾存在時，由于其具有高強度特點，目標信號被淹沒，干擾占主導地位。由于瞬態干擾的偏振狀態是雜亂無章的，質點偏離固定軌跡向隨機方向運動。

3.2 瞬態干擾抑制

接收信號在未存在瞬態干擾和存在時，具有不同偏振狀態。Flinn設計了極化濾波器，利用調制函數將具有某種特定偏振狀態的信號保留，并減弱具有其他狀態的信號。信號分離公式為[13]

式中D1(t)為與橢圓反橢圓率有關的加權系數，D2(t)為與偏振方向有關的加權系數，S(t)為需要降噪的信號。D1(t),D2(t)共同作用，構成調制函數，實現濾波功能。Alkaz提出了將錐形高斯函數作為權值[12]

式中λ0(t)為期望的直線性，σλ為反橢圓率標準偏差，Δθ(t)是計算得到的偏振方向與期望的偏振方向之間的角度差，σθ為角度標準偏差。

理論上，最佳的期望直線性和角度為無瞬態干擾時信號的偏振參數，但由于不能確定瞬態干擾位置，實際上并未如此設置。考慮到瞬態干擾持續時間非常短的特點，受瞬態干擾影響的信號長度遠小于無干擾信號長度，采用平均值作為公式的參數，可將瞬態干擾的影響降低。將期望的直線性和角度差設置為整段信號偏振參數的平均值，標準偏差為整段信號偏振參數的均方差。

從式(9)和式(10)中可看出，加權系數對偏振參數的變化具有一定的“容忍”區間，若偏振參數在此區間內，加權系數比較穩定，在0.9-1.0之間；若超出此區間，輸出結果便會迅速減小，甚至會接近于0。瞬態干擾的存在，導致偏振參數不規則變化，且超過“容忍”區間，加權函數值急劇下降。利用式(8)進行極化濾波后，瞬態干擾處信號幅度會遠小于無干擾處信號幅度，甚至接近于 0，實現瞬態干擾抑制。

從上述分析可知，當調制函數值急劇下降，瞬態干擾到達；調制函數值急劇上升時，已無瞬態干擾影響。對調制函數值設定門限，認為低于門限時存在瞬態干擾。噪聲的存在，導致兩個加權函數值在0.9-1.0間抖動。綜合考慮兩個加權函數，門限值應小于0.8，具體數值應根據不同情況分析。本文中以D1(t)D2(t)≤0.707作為干擾存在的判斷依據，此點為半功率點，即信號功率下降為原來一半。

窄帶信號具有明顯的偏振特性，寬帶信號則不盡然。由于傳輸介質的頻散特性，不同頻率聲波傳播速度不同，這種頻散現象在低頻段尤為明顯。若分析信號包含頻散特性明顯的頻率區域，不同頻率成分到達時間有一定差別。此時，極化狀態不穩定，隨時間變化，需要將寬帶信號看作多個窄帶信號，對每段窄帶信號單獨處理，再合成寬帶結果。對于頻散特性不明顯的接收信號，需要考慮兩種不同情況：信號包絡隨時間劇烈變化和信號包絡緩變(或基本不變)。對于前一種信號，同樣需要將寬帶信號看作多個窄帶信號；后一種信號包絡基本不變，且頻散影響較小，比如線性調頻信號，其橢圓軌跡比較固定，可以直接采用本文所提出方法進行瞬態干擾抑制。

3.3 數據恢復

當信號呈現某種特定偏振狀態時，矢量信號可以由偏振參數表達

式中Δ為坐標點與原點連線方向和主軸方向夾角，A和B分別為偏振橢圓長短軸長度，θ為偏振方位角，φ為偏振俯仰角，γ為橢圓短軸方向與y軸夾角。上述參數中除Δ外均可由偏振分析獲得，則獲得Δ值便可對瞬態干擾處信號進行恢復。Δ值的估計方法有很多種，包括BP神經網絡，AR模型預測和曲線擬合等。曲線多項式擬合方法適用于規律較簡單的曲線，對于不滿足多項式規律的曲線誤差較大；BP神經網絡為非線性預測方法，適用于已知樣本點比較隨機的情況；AR模型可以從已知樣本點尋找規律，對某一時刻信號幅度進行預測。由于質點沿著固定的軌跡運動，夾角Δ值滿足某一規律，為了更好地獲得參數規律，本文中選擇AR模型對參數進行預測。若已知{Δ(m-P),Δ(m-P+1),…,Δ(m-1)}，則可以利用它們來預測m時刻值Δ(m)

式中αk,k=1,…,P為AR模型參數。改進的協方差方法具有較好的估計性能，文中利用此方法確定AR模型參數值[14]。

4 處理結果分析

4.1 仿真結果分析

仿真環境：發射頻率為400 Hz信號，目標方位角為50°，俯仰角為20°。無噪聲和瞬態干擾的3維矢量信號可由幾何關系得到，將3個振速分量分別記為vx,vy,vz，它們的功率為-7.3892 dB,-5.8655 dB和-12.3293 dB，噪聲功率為-20 dB，在0.0501-0.053 s時間段內存在瞬態干擾。

為了將瞬態干擾去除，采用滑動窗方法分析信號的偏振特性。根據文中方法設計加權函數，構建極化濾波器，極化濾波器的加權函數見圖1。從圖1中可以看出，第1個加權函數的權值在0.05 s附近很小，在0-0.2范圍內變化；理論上權值應該在0.05 s點上迅速下降，但從圖中觀察到，權值是逐漸下降的，由于采用滑動窗處理方法，相當于對輸出結果進行了平滑處理；而其他時間段內權值在0.9-1.0范圍內變化。

第 2個加權函數權值在 0.5s附近也較小，在0-0.9范圍內變化；其他時間段內權值在0.9-1.0范圍內變化。

按照上述方法對接收信號進行極化濾波，再對0.05 s附近信號進行數據恢復。原始信號與濾波后信號時域圖見圖2，對濾波前后信號進行時頻分析，分析結果見圖3。

通過觀察圖2可以知道：瞬態干擾基本被抑制；由于過渡帶的存在，濾波后輸出信號與原信號稍微有些不同，考慮到過渡帶時間遠小于信號脈沖時間，可以基本忽略其對信號的負面影響；其他時間段內信號均被保留，信號所攜帶的信息基本沒有損失。圖3中，濾波后瞬態干擾完全被抑制，信號的線譜成分更加明顯。

4.2 實測數據結果分析

實測數據的中心頻率為600 Hz的單頻矩形脈沖信號，截取存在瞬態干擾的一段信號進行分析，截取信號長度為0.05s，采樣頻率為10 kHz。

對信號加滑動窗，分析不同時刻偏振狀態，提取極化參數，利用偏振狀態變化情況對瞬態干擾進行抑制。聯合利用3個極化參數，選取合適的標準差，設計極化濾波器，加權函數見圖4。

圖1 根據偏振參數設計的極化濾波器加權函數

圖2 原始信號與濾波后信號時域對比圖

圖3 原始信號與處理后信號LOFAR圖

圖4 實測信號對應的極化濾波器加權函數

圖5 實測信號與濾波后信號對比圖

第1個加權函數，在0.012-0.014 s時間段內，權值在0.7-0.9范圍內變化；在0.031-0.033 s時間段內，權值的變化范圍為0.2-0.5；其余時刻的權值在0.9-1.0范圍內變化。第2個加權函數，在0.012-0.014 s時間段內，權值的變化范圍是0.7-0.9；在0.031-0.033 s時間段內，權值的變化范圍是 0.5-0.9；其余時刻權值在0.9-1.0范圍內變化，兩個權函數變化規律基本相同。

利用圖4所示權函數對實測信號進行加權，實現濾波功能，對0.012-0.014 s和0.031-0.033 s時間段內信號進行數據恢復。原始信號與濾波后信號時域圖見圖5。

在圖5中，瞬態干擾被很好地抑制，由于過渡帶的存在，濾波后瞬態干擾旁邊的信號與原始信號稍微有些不同。由于過渡帶時間非常短暫，且濾波引起的變化很小，可以認為該方法能夠抑制瞬態干擾，保留未受干擾污染的信號，并保證信號攜帶的信息基本沒有損失。圖6中，0.032 s附近的脈沖強度較大，0.13 s附近的瞬態干擾強度較弱，對比濾波前后信號的時頻分析結果，可以看出通過濾波及數據恢復已將兩處瞬態干擾去除。

5 結束語

當接收的矢量信號中混有瞬態干擾時，信號處理結果會受其影響，處理性能下降，必須對瞬態干擾進行抑制。本文提出了一種以偏振分析為基礎的瞬態干擾消除方法。此方法通過極化濾波可直接去除瞬態干擾，利用調制函數能夠比較直觀的判斷瞬態干擾的存在位置，不需要復雜計算，便于實現。仿真結果和實驗數據均證明了該方法的有效性。

圖6 實測信號與處理后信號LOFAR圖

[1]楊德森, 洪連進. 矢量水聽器原理及應用引論[M]. 第1版, 北京: 科學出版社, 2009: 41-45.

Yang De-sen and Hong Lian-jin. Vector Sensor Principle and Application[M]. 1st Edition, Beijing: Science Press, 2009:41-45.

[3]周忠根, 水鵬朗. 基于復數據經驗模式分解的天波超視距雷達瞬態干擾抑制[J]. 電子與信息學報, 2011, 33(12):2831-2836.

Zhou Zhong-gen and Shui Peng-lang.Transient interference suppression based on complex empirical mode decomposition in over-the-horizon radar[J].Journal of Electronics&Information Technology, 2011, 33(12): 2831-2836.

[4]王階, 劉濤, 陳曉旭, 等. 天波超視距雷達流星余跡瞬態干擾抑制方法[J]. 電子測量與儀器學報, 2009, 23(10): 67-72.

Wang Jie, Liu Tao, Chen Xiao-xu,et al..Method of meteor trail impulsive interference suppression in OTHR[J].Journal of Electronic Measurement and Instrument, 2009, 23(10):67-72.

[5]甄紅飛, 龔耀寰, 劉濤. 基于TLS-ESPRIT 算法的OTHR瞬態干擾抑制[J]. 雷達科學與技術, 2011, 9(4): 341-346.

Zhen Hong-fei, Gong Yao-huan, and Liu Tao.Research on impulsive interference suppression in OTHR based on TLS-ESPRIT algorithm[J].Radar Science and Technology,2011, 9(4): 341-346.

[6]陳希信, 黃銀河. 基于矩陣奇異值分解的高頻雷達瞬態干擾抑制[J]. 電子與信息學報, 2005, 27(12): 1879-1882.

Chen Xi-xin and Huang Yin-he.A SVD-based approach of suppressing transient interference in high-frequency radar[J].Journal of Electronics&Information Technology, 2005,27(12): 1879-1882.

[7]Liu Tao, Chen Xiao-xu, and Wang Jie. Subspace impulsive interference suppression in OTHR[J].Progressin Electromagnetics Research, 2009, C(7): 167-181.

[8]Turley M. Impulsive noise rejection in HF radar using a linear prediction technique[C]. Radar Conference, United States, May 5-8, 2003, 9: 358-362.

[9]戴幻堯, 李永禎, 王雪松, 等. 基于和差波束極化特性的目標極化散射矩陣測量方法研究[J]. 電子與信息學報, 2010, 32(4):913-918.

Dai Huan-yao, Li Yong-zhen, Wang Xue-song,et al..A new target scattering matrix measurement algorithm based on polarization characteristics of sum-and-difference beam[J].Journal of Electronics&Information Technology, 2010, 32(4):913-918.

[10]劉勇, 李永禎, 戴幻堯, 等. 基于極化二元陣雷達的空域虛擬極化濾波算法[J]. 電子與信息學報, 2010, 32(11): 2746-2750.

Liu Yong, Li Yong-zhen, Dai Huan-yao,et al..The spatial virtual polarization filtering algorithm based on dual polarized antenna radar[J].Journal of Electronics&Information Technology, 2010, 32(11): 2746-2750.

[11]王克讓, 朱曉華, 何勁. 基于矢量傳感器MIMO雷達的DOD DOA和極化聯合估計算法[J]. 電子與信息學報, 2012, 34(1):160-165.

Wang Ke-rang, Zhu Xiao-hua, and He Jin.Joint DOD DOA and polarization estimation for MIMO radar with electromagnetic vector sensors[J].Journal of Electronics&Information Technology, 2012, 34(1): 160-165.

[12]Alkaz V G, Ornofrash N I, and Perelberg A I. Polarization Analysis of Seismic Waves[M]. 1st Edition, Russian, Shtiintca Press, 1977: 100-105.

[13]Flinn E A. Signal analysis using rectilinearity and direction of particle motion[J].Proceedings of IEEE, 1965, 53(12):1874-1876.

[14]胡廣書. 數字信號處理理論、算法與實現[M]. 第 2版, 北京:清華大學出版社, 2006: 527-574.

Hu Guang-shu. Digital Signal Processing Theory, Algorithm and Realization[M]. 2nd Edition, Beijing: Tsinghua University Press, 2006: 527-574.