基於馬氏距離重加權的非均勻雜波抑制方法與流程
2023-07-02 09:37:21
本發明屬於雷達信號處理技術領域，特別涉及基於馬氏距離重加權的非均勻雜波抑制方法，用於非均勻雜波環境下的雜波抑制及慢速目標檢測。
背景技術：
機載預警機雷達由於常處在下視工作狀態，其工作性能必然會受到地面雜波的強烈影響。同時，由於機載雷達平臺與地面的相對運動，使得機載雷達接收的雜波呈現出很強的空時耦合性，給弱目標、慢速目標的檢測帶來很大困難。在上世紀七十年代，brennan等人提出了空時自適應處理(space-timeadaptiveprocessor，stap)方法，通過在空域和時域的二維處理來進行雜波抑制達到弱目標和慢速目標檢測的目的。經過這些年的研究和發展，stap方法成為一種有效的雜波抑制手段，並且可以有效提高雷達對弱小、慢速目標檢測的能力。
通常，stap方法在理想情況下，即環境均勻情況下，可以取得較好的性能。stap方法在實施時需要利用訓練樣本來估計待檢測單元的雜波協方差矩陣並計算自適應權矢量。理想狀況下訓練樣本滿足獨立同分布且不包含目標信號時，stap方法可以取得較好的雜波抑制性能。然而，在真實情況中機載雷達面臨的實際雜波環境一般是非均勻的，尤其是在地面慢速運動目標檢測時，由於機載雷達波束照射範圍內運動目標比較多，訓練樣本中往往包含有目標信號，從而導致訓練樣本不符合獨立同分布狀態，此時樣本的均勻假設條件將不再成立。如果利用這些包含目標信號的非均勻的樣本來估計待檢測單元的雜波協方差矩陣並計算自適應權矢量，進而利用自適應權矢量進行濾波，將導致待檢測樣本中目標信號自相消，使得目標檢測概率下降。
為此，有學者提出了非均勻檢測器(nonhomogeneitydetector，nhd)，nhd能有效的檢測並剔除被汙染的訓練樣本，使其不參與權矢量的運算，因而能改善stap的性能。現有的nhd有：(1)廣義內積檢測器(generalizedinnerproduct，gip)，gip對於訓練樣本中存在明顯目標信號時檢測效果較好，然而由於廣義內積法只提取協方差矩陣中的大特徵值，故在多目標信號的情況下，其檢測性能將會下降；自適應功率剩餘(adaptivepowerresidue，apr)檢測器，同樣該檢測器在訓練樣本中的目標信號數量較少、信噪比較低時可以獲取較好的檢測效果，當訓練樣本中的目標信號數量較多、信噪比較高時，訓練樣本所構成的幹擾加噪聲協方差矩陣會受到目標信號嚴重擾動，apr檢測器性能下降；對角加載的廣義內積檢測器(ioadedgeneralizedinnerproduct，lgip)，該檢測器通過對由訓練樣本估計出的協方差矩陣對角加載來降低目標信號對檢測器的擾動，存在如何選擇合適的加載系量以提高檢測性能的難題；基於長球波函數的廣義內積檢測器(prolatespheroidalwavefunctionsgip，pswf-gip)，該檢測器通過長球波函數的基向量直接估計待檢測樣本的雜波協方差矩陣，避免了訓練樣本中目標信號的影響，在理論上可以取得較好的性能，但在計算長球波函數的基向量時要求雷達構型參數及陣列流型精確已知，而實際情況中的陣列誤差以及其他因素會影響陣列流型及雷達構型參數，使得檢測性能下降。
技術實現要素：
有鑑於此，本發明提供了基於馬氏距離重加權的非均勻雜波抑制方法，能夠在非均勻雜波情況下有效抑制雜波，從而提高檢測性能。
為達到上述目的，本發明採用如下技術方案：
一種基於馬氏距離重加權的非均勻雜波抑制方法，包括以下步驟：
步驟1，獲取機載雷達天線陣列接收到的原始回波數據，天線陣列為包括n個陣元且陣元間距為半波長的均勻線陣，原始回波數據包括目標分量、雜波分量以及噪聲分量；
步驟2，對原始回波數據進行下變頻處理和脈衝壓縮處理，得到空時快拍數據；
利用降維矩陣對空時快拍數據進行降維處理，得到降維後的空時快拍數據，降維後的空時快拍數據包含l個距離單元的數據；
獲取目標空時導向矢量，利用降維矩陣對目標空時導向矢量進行降維處理，得到降維後的目標空時導向矢量；
步驟3，初始化：令待檢測單元序號l＝1，即令l個距離單元中的第l個距離單元為當前待檢測單元；
步驟4，假設第l個距離單元存在對應的kl個採樣樣本，對第l個距離單元的kl個採樣樣本進行求平均，得到第l個距離單元的基準樣本，並計算第l個距離單元的基準樣本對應的基準協方差矩陣；
利用基準協方差矩陣，計算得到kl個採樣樣本中每個採樣樣本與基準樣本之間的馬氏距離，進而利用各採樣樣本與基準樣本之間的馬氏距離，構造每個採樣樣本對應的剔除係數；
利用各採樣樣本對應的剔除係數以及預設的選取門限，確定kl個採樣樣本中的幹擾樣本，並剔除kl個採樣樣本中的幹擾樣本，得到kl′個數據樣本；其中，幹擾樣本為剔除係數大於選取門限的採樣樣本；
步驟5，利用基準協方差矩陣，計算得到kl′個數據樣本與基準樣本間的馬氏距離；
判斷kl′是否大於或等於6n：若kl′大於或等於6n，則選取kl′個數據樣本中與基準樣本間馬氏距離最小的6n個數據樣本作為第l個距離單元的訓練樣本；若kl′小於6n，則選取kl′個數據樣本作為第l個距離單元的訓練樣本；
利用每個訓練樣本與基準樣本間馬氏距離，構造該訓練樣本對應的重加權係數；
利用各訓練樣本及其對應的重加權係數進行重加權，構造得到第l個距離單元對應的採樣協方差矩陣；
步驟6，利用降維後的目標空時導向矢量以及第l個距離單元對應的採樣協方差矩陣，構造得到第l個距離單元對應的自適應權矢量；
利用第l個距離單元對應的自適應權矢量對第l個距離單元的數據進行雜波抑制，得到第l個距離單元對應的輸出功率；
步驟7，令待檢測單元序號l加1，重複執行步驟4至步驟6，直至l＝l，即完成對全部l個距離單元的雜波抑制。
基於本發明上述方案，利用非均勻雜波數據的特性，通過採樣數據與基準數據間的馬氏距離確定各採樣數據對應的剔除係數，利用剔除係數及預設的選取門限剔除採樣樣本中的幹擾樣本得到數據樣本，進而利用所得到的數據樣本與基準樣本之間的馬氏距離，選取對應的訓練樣本，再進一步利用訓練樣本與基準樣本之間的馬氏距離構造對應的重加權係數進行重加權，構造得到對應的採樣協方差矩陣，最後利用所得的協方差矩陣構造自適應權矢量，利用該自適應權矢量進行雜波抑制，得到對應的輸出功率。因此，相比於傳統廣義內積法，本發明方法在採樣樣本數不同的情況下可以穩健的獲得較好的雜波抑制性能。同時，本發明方法對於複雜雜波背景中的雜波非均勻抑制具有普適性和穩健性，能夠有效保證stap在非均勻雜波環境中的雜波抑制性能。
附圖說明
為了更清楚地說明本發明實施例或現有技術中的技術方案，下面將對實施例或現有技術描述中所需要使用的附圖作簡單地介紹，顯而易見地，下面描述中的附圖僅僅是本發明的一些實施例，對於本領域普通技術人員來講，在不付出創造性勞動的前提下，還可以根據這些附圖獲得其他的附圖。
圖1為本發明實施例提供的一種基於馬氏距離重加權的非均勻雜波抑制方法的流程示意圖；
圖2為採樣樣本取2×3n時本發明方法和gip算法的雜波剩餘曲線對比圖，圖中橫坐標為都卜勒通道號，縱坐標為功率大小，功率單位為db；
圖3為採樣樣本取3×3n時本發明方法和gip算法的雜波剩餘曲線對比圖，圖中橫坐標為都卜勒通道號，縱坐標為功率大小，功率單位為db；
圖4為採樣樣本取4×3n時本發明方法和gip算法的雜波剩餘曲線對比圖，圖中橫坐標為都卜勒通道號，縱坐標為功率大小，功率單位為db；
圖5a為本發明方法在採樣樣本分別取2×3n、3×3n、4×3n情況下的雜波剩餘曲線圖；
圖5b為gip算法在採樣樣本分別取2×3n、3×3n、4×3n情況下的雜波剩餘曲線圖。
具體實施方式
下面將結合本發明實施例中的附圖，對本發明實施例中的技術方案進行清楚、完整地描述，顯然，所描述的實施例僅僅是本發明一部分實施例，而不是全部的實施例。基於本發明中的實施例，本領域普通技術人員在沒有做出創造性勞動前提下所獲得的所有其他實施例，都屬於本發明保護的範圍。
圖1所示為本發明實施例提供的一種基於馬氏距離重加權的非均勻雜波抑制方法的流程示意圖。
如圖1所示，本發明實施例提供的基於馬氏距離重加權的非均勻雜波抑制方法，包括以下步驟：
步驟1，獲取機載雷達天線陣列接收到的原始回波數據。
其中，機載雷達天線陣列為包括n個陣元且陣元間距為半波長的均勻線陣，原始回波數據包括目標分量、雜波分量以及噪聲分量。
為敘述方便，在本實施例中，假設機載雷達在一個相干處理間隔(coherentprocessinginterval，cpi)內發射m個脈衝，其脈衝重複頻率為fr，雷達工作波長為λ，載機高度為h，載機速度為v。
步驟2，對原始回波數據進行下變頻處理和脈衝壓縮處理，得到空時快拍數據，進而利用降維矩陣對空時快拍數據進行降維處理，得到降維後的空時快拍數據，以及獲取目標空時導向矢量，利用降維矩陣對目標空時導向矢量進行降維處理，得到降維後的目標空時導向矢量。
其中，降維後的空時快拍數據包含l個距離單元的數據。假設其中的第k個距離單元的空時快拍數據可表示為則該距離單元對應的二元假設檢驗可以表示為：
h0：xk＝ck+nk
h1：xk＝αss-t+xk/h0』
其中，h0假設表示該距離單元的回波數據不包含目標信號，h1假設表示該距離單元的回波數據包含目標信號；αss-t，ck，nk分別表示目標分量、雜波分量和噪聲分量，其中α代表目標的復幅度，ss-t表示目標空時導向矢量，式中代表克羅內克積，ss表示時間導向矢量，st表示空間導向矢量。
其中，時間導向矢量ss和空間導向矢量st又可分別表示為：
式中，表示歸一化都卜勒頻率，表示歸一化空間頻率，fd表示目標都卜勒頻率，θ表示方位角，表示俯仰角，d表示陣元間距。
具體的，步驟2中，利用降維矩陣對空時快拍數據進行降維處理，得到降維後的空時快拍數據，可以包括：
根據第三預設公式，利用降維矩陣對空時快拍數據進行降維處理，得到降維後的空時快拍數據及降維後的目標空時導向矢量。
其中，第三預設公式包括：
式中，x表示空時快拍數據，t表示降維矩陣，表示降維後的空時快拍數據，xl表示第l個距離單元的數據，l＝1，2，…，l，上標h表示共軛轉置操作，的維數為3n×l，x的維數為mn×l，t的維數為mn×3n，m表示機載雷達在一個相干處理間隔內發射的脈衝個數。其中，降維矩陣t可表示為其中f表示m×3維的傅立葉矩陣加權，in表示n×n維的單位矩陣，表示克羅內克積。
具體的，步驟2中，利用降維矩陣對空時快拍數據進行降維處理，得到降維後的空時快拍數據，具體可以包括：
根據第四預設公式，利用降維矩陣對目標空時導向矢量進行降維處理，得到降維後的目標空時導向矢量。
其中，第四預設公式包括：
式中，ss-t表示目標空時導向矢量，表示降維後的目標空時導向矢量。
步驟3，初始化：令待檢測單元序號l＝1，即令l個距離單元中的第l個距離單元為當前待檢測單元。
步驟4，假設第l個距離單元存在對應的kl個採樣樣本，對第l個距離單元的kl個採樣樣本進行求平均，得到第l個距離單元的基準樣本，並計算第l個距離單元的基準樣本對應的基準協方差矩陣；利用基準協方差矩陣計算得到kl個採樣樣本中每個採樣樣本與基準樣本之間的馬氏距離，進而利用各採樣樣本與基準樣本之間的馬氏距離，構造每個採樣樣本對應的剔除係數；利用各採樣樣本對應的剔除係數以及預設的選取門限，確定kl個採樣樣本中的幹擾樣本，並剔除kl個採樣樣本中的幹擾樣本，得到kl′個數據樣本。
其中，幹擾樣本為剔除係數大於選取門限的採樣樣本。需要說明的是，所述的選取門限具體可設置為0.3，即當某一採樣樣本對應的馬氏距離佔總體比重大於0.3時，即視其為幹擾數據將其剔除。當然，本領域技術人員可以理解，以上所述僅為本發明實施例提供的一種可選的實現方式，在實際應用中可根據具體需求設置選取門限，本發明實施例對此不作具體限定。
具體的，步驟4中，對第l個距離單元的kl個採樣樣本進行求平均，得到第l個距離單元的基準樣本，並計算第l個距離單元的基準樣本對應的基準協方差矩陣，具體可以包括以下步驟：
對第l個距離單元的kl個採樣樣本進行求平均，得到第l個距離單元的基準樣本
利用第l個距離單元的基準樣本μ，計算第l個距離單元的基準樣本對應的基準協方差矩陣ru＝μμt+σn2i。
其中，i＝1，2，…，kl，σn表示對角加載係數，σn取底噪水平，i表示單位矩陣，上標t表示轉置操作。
進一步的，步驟4中，利用基準協方差矩陣，計算得到kl個採樣樣本中每個採樣樣本與基準樣本之間的馬氏距離，進而利用各採樣樣本與基準樣本之間的馬氏距離，構造每個採樣樣本對應的剔除係數，具體可以包括以下步驟：
利用基準協方差矩陣，計算得到kl個採樣樣本中每個採樣樣本與基準樣本之間的馬氏距離；
計算各採樣樣本與基準樣本間馬氏距離的倒數，進而對各採樣樣本與基準樣本間馬氏距離的倒數進行歸一化處理；
將各採樣樣本對應的歸一化處理後得到的數值作為該採樣樣本的剔除係數。
其中，第i個採樣樣本與基準樣本之間的馬氏距離為第i個採樣樣本的剔除係數為其中上標-1表示求逆操作。
步驟5，利用基準協方差矩陣，計算得到kl′個數據樣本與基準樣本間的馬氏距離；判斷kl′是否大於或等於6n：若kl′大於或等於6n，則選取kl′個數據樣本中與基準樣本間馬氏距離最小的6n個數據樣本作為第l個距離單元的訓練樣本；若kl′小於6n，則選取kl′個數據樣本作為第l個距離單元的訓練樣本；利用每個訓練樣本與基準樣本間馬氏距離，構造該訓練樣本對應的重加權係數；利用各訓練樣本及其對應的重加權係數進行重加權，構造得到第l個距離單元對應的採樣協方差矩陣。
本發明實施例的一種優選實現方式中，步驟5中，利用每個訓練樣本與基準樣本間馬氏距離，構造該訓練樣本對應的重加權係數，具體可以包括：
計算得到每個訓練樣本與基準樣本間馬氏距離的倒數，進而對各訓練樣本與基準樣本間馬氏距離的倒數進行歸一化處理；
將各訓練樣本對應的歸一化處理後的得到的數值作為該訓練樣本的重加權係數。
其中，第l個距離單元的第k個訓練樣本的剔除係數為dlk表示第k個訓練樣本與基準樣本之間的馬氏距離，k＝1，2，…，k，k表示第l個距離單元的訓練樣本個數。
可以理解，k的具體取值取決於kl′。具體來說，當kl′≥6n時，k＝6n；當kl′＜6n時，k＝kl′。
具體的，步驟5中，利用各訓練樣本及其對應的重加權係數進行重加權，構造得到第l個距離單元對應的採樣協方差矩陣，具體包括：
利用各訓練樣本及其對應的重加權係數，根據第一預設公式進行重加權，構造得到第l個距離單元對應的採樣協方差矩陣。
其中，第一預設公式包括：
式中，表示第l個距離單元對應的採樣協方差矩陣，表示第l個距離單元的第k個訓練樣本，表示第l個距離單元的第k個訓練樣本對應的重加權係數，上標h表示共軛轉置操作，∑表示求和符號，k＝1，2，…，k，k表示第l個距離單元的訓練樣本個數。
步驟6，利用降維後的目標空時導向矢量以及第l個距離單元對應的採樣協方差矩陣，構造得到第l個距離單元對應的自適應權矢量；利用第l個距離單元對應的自適應權矢量對第l個距離單元的數據進行雜波抑制，得到第l個距離單元對應的輸出功率。
具體的，步驟6中，利用降維後的目標空時導向矢量以及第l個距離單元對應的採樣協方差矩陣，構造得到第l個距離單元對應的自適應權矢量，具體可以包括：
利用降維後的目標空時導向矢量以及第l個距離單元對應的採樣協方差矩陣，根據第二預設公式，構造得到第l個距離單元對應的自適應權矢量。
其中，第二預設公式包括：
式中，w表示第l個距離單元對應的自適應權矢量，表示第l個距離單元對應的採樣協方差矩陣，表示目標空時導向矢量，上標-1表示求逆操作，上標h表示共軛轉置操作。
進一步的，步驟6中，利用第l個距離單元對應的自適應權矢量對第l個距離單元的數據進行雜波抑制，得到第l個距離單元對應的輸出功率，具體包括：
根據第五預設公式，利用第l個距離單元對應的自適應權矢量對第l個距離單元的數據進行雜波抑制，得到第l個距離單元對應的輸出功率。
其中，第五預設公式包括：
yl＝whxl，
式中，w表示第l個距離單元對應的自適應權矢量，yl表示第l個距離單元對應的輸出功率。
步驟7，令待檢測單元序號l加1，重複執行步驟4至步驟6，直至l＝l，即完成對全部l個距離單元的雜波抑制。
至此，即完成了對全部l個距離單元的雜波抑制，本發明實施例提供的基於馬氏距離重加權的非均勻雜波抑制方法即結束。
基於本發明上述方案，利用非均勻雜波數據的特性，通過採樣數據與基準數據間的馬氏距離確定各採樣數據對應的剔除係數，利用剔除係數及預設的選取門限剔除採樣樣本中的幹擾樣本得到數據樣本，進而利用所得到的數據樣本與基準樣本之間的馬氏距離，選取對應的訓練樣本，再進一步利用訓練樣本與基準樣本之間的馬氏距離構造對應的重加權係數進行重加權，構造得到對應的採樣協方差矩陣，最後利用所得的協方差矩陣構造自適應權矢量，利用該自適應權矢量進行雜波抑制，得到對應的輸出功率。因此，相比於傳統廣義內積法，本發明方法在採樣樣本數不同的情況下可以穩健的獲得較好的雜波抑制性能。同時，本發明方法對於複雜雜波背景中的雜波非均勻抑制具有普適性和穩健性，能夠有效保證stap在非均勻雜波環境中的雜波抑制性能。
以下通過仿真實驗驗證本發明方法的正確性及上述有益效果：
1、實測數據：
本發明通過採用機載雷達的實測數據驗證算法的可行性，本發明採用的是mcarm數據，mcarm數據是在1996年由room實驗室為研究stap及其相關技術而錄取的一些列機載雷達數據。本仿真實驗具體採用的是第5次飛行中錄取序號為575的數據來進行實驗分析，具體參數如表1所示。
表1
2、仿真內容及結果分析
(1)仿真內容
鑑於mcarm數據的前150個距離單元數據受汙染，故本仿真實驗從mcarm數據的第151個距離門開始採樣，選取第151到第384個距離門的數據作為實驗數據。利用本發明方法和傳統的gip算法，分別在採樣樣本數為2×3n、3×3n、4×3n時進行雜波抑制，並對本發明方法和傳統的gip算法的雜波抑制性能進行比較。其中，圖3為採樣樣本取2×3n時本發明方法和gip算法的雜波剩餘曲線對比圖，圖4為採樣樣本取3×3n時本發明方法和gip算法的雜波剩餘曲線對比圖，圖5為採樣樣本取4×3n時本發明方法和gip算法的雜波剩餘曲線對比圖。其中，所述的雜波剩餘曲線指的是距離向取平均後的功率曲線。
(2)仿真結果分析
從圖2中可以看出，本發明方法和gip算法在採樣數為2×3n時的雜波抑制性能相差不大，雜波功率曲線總體比較吻合。
從圖3中可以看出，本發明方法和gip算法在採樣數為3×3n時，本發明中提及的算法的雜波抑制性能優於gip的雜波抑制性能。
從圖4中可以看出，本發明方法和gip算法在採樣數為4×3n時，本發明中提及的算法的雜波抑制性能優於gip的雜波抑制性能。
從圖5a、圖5b中可以看出，本發明方法和gip算法在採樣數分別為2×3n、3×3n、4×3n時，gip算法隨著採樣數的增多抑制雜波的性能逐漸下降，但是本發明提及算法在採樣數增多時可以保持一個相對穩健的雜波抑制性能，總體上來說雜波抑制性能優於傳統的gip算法。
本領域普通技術人員可以理解：實現上述方法實施例的全部或部分步驟可以通過程序指令相關的硬體來完成，前述的程序可以存儲於一計算機可讀取存儲介質中，該程序在執行時，執行包括上述方法實施例的步驟；而前述的存儲介質包括：rom、ram、磁碟或者光碟等各種可以存儲程序代碼的介質。
以上，僅為本發明的具體實施方式，但本發明的保護範圍並不局限於此，任何熟悉本技術領域的技術人員在本發明揭露的技術範圍內，可輕易想到變化或替換，都應涵蓋在本發明的保護範圍之內。因此，本發明的保護範圍應以權利要求的保護範圍為準。