李 海，呼延澤，毛志杰

（1.中國民航大學電子信息與自動化學院，天津 300300；2.國防科技大學信息通信學院，陜西 西安 710106）

0 引 言

低空風切變是一種具有突發(fā)性、短持續(xù)、強破壞特點的危害性氣象現象，大多發(fā)生在低于600 m的低空區(qū)域。當飛機飛行到該區(qū)域時，極有可能被風切變危害而引發(fā)慘重的空難事件。因此，準確估計風場速度，進而探測到低空風切變，并及時給飛機提供告警信息非常重要。

機載氣象雷達能夠幫助飛機實時監(jiān)視前方的氣象狀況，有助于飛行安全，其在工作時較多采用厘米波。而毫米波的波長更短，波束更窄。機載氣象雷達發(fā)射毫米波，可以縮小天線，并且更加準確地探測目標。目前，機載氣象雷達普遍采用脈沖體制，而毫米波線性調頻連續(xù)波（linear frequency modulated continuous wave，LFMCW）體制信號處理系統(tǒng)結構更為簡單，并且重量輕、對發(fā)射峰值功率要求低。將LFMCW技術用于機載氣象雷達，能夠使其向小型化發(fā)展，進而提高載機的靈活性、機動性。

檢測風場時，機載氣象雷達會嚴重被地雜波信號干擾，難以識別到目標。此時常采用空時自適應處理（space-time adaptive processing，STAP）技術壓制雜波。在此基礎上，已有多種風速估計方法，但它們都未解決速度模糊問題。毫米波雷達系統(tǒng)因其具有短波長的特點，很容易出現速度估計模糊，這給低空風切變的檢測帶來阻礙。因此，必須考慮解模糊，目前主流的解模糊方法通常是以脈沖重復頻率組為基礎進行的。文獻［8］提出了余差查表法，該算法利用目標在各重頻上的余數之差進行解模糊，其對存儲空間要求較大，并且在查表過程中會浪費大量時間在不必要的模糊值上，實時性不高。文獻［9］提出剩余定理（孫子定理）法，通過對不同重頻下得到的同余方程組進行求解來實現解模糊，該算法對測量精度要求較高，一旦存在測量誤差，解模糊的質量難以保證。文獻［10］提出一維集算法，通過將所有可能的速度值從小到大排列，再尋找多普勒頻率組的最小均方值進行解模糊，該方法存在計算量大的問題。

速度模糊問題，究其根本，是目標的多普勒頻率大于了調頻周期重復頻率的一半，此時根據香農 奈奎斯特采樣定理，會產生多普勒混疊?？梢?，雷達系統(tǒng)可測的最大不模糊速度范圍是受香農 奈奎斯特采樣定理限制的，近年來出現的壓縮感知（compressive sensing，CS）理論很好地突破了這一限制。因此，本文考慮用CS理論來解決速度模糊問題。首先利用一個相干處理時間內回波信號非均勻欠采樣特性在角度 多普勒域構造無模糊冗余字典，利用CS技術提取雜波譜主要成分，估計雜波能量支撐域，構造出雜波抑制矩陣。之后利用加權最小范數優(yōu)化模型實現雜波抑制與低空風切變的速度解模糊。此方法只需對待測距離門數據進行處理，無需像傳統(tǒng)STAP方法那樣估計雜波協方差矩陣，因此計算量大大減小。此外，該方法在角度 多普勒域進行細分來構造字典，因而對運動目標有更高的角度和多普勒分辨率。

1 信號模型

假設機載氣象雷達結構如圖1所示。載機平臺距地面距離為，速度為，天線陣元數為，相鄰陣元距離＝／2，為信號波長。和分別對應強地雜波的方位角和俯仰角，和φ分別對應風場的方位角和俯仰角。設定雷達在相干處理間隔內發(fā)射個周期的調頻信號。

圖1 機載雷達示意圖Fig.1 Schematic diagram of airborne radar

為了解決毫米波段引起的測速模糊問題，信號模型中采用兩重頻的脈組參差方式，兩重脈沖重復頻率（pulse repetition frequency，PRF）的脈組參差方式采樣時序圖如圖2所示。

圖2 兩重PRF的脈組參差方式Fig.2 Mode of pulse group stagger in double PRF

將一個完整的相干處理時間間隔（coherent processing interval，CPI）均分為兩個子CPI，設使用的兩重頻分別為f 、f ，則對應的采樣間隔分別為T ＝1／f 、T ＝1／f 。假設機載氣象雷達在其工作范圍內有個距離單元，第個距離單元的回波數據（）＝［（），（）］，其中，z （）為第（＝0，1）個子CPI內第個距離單元的空時二維快拍數據：

式中：c （）是第（＝0，1）個子CPI內地雜波；n 為高斯白噪聲；s （）是低空風切變信號。

令第（＝0，1）個子CPI內第（＝1，2，…，／2）個周期的發(fā)射信號為

式中：ξ為回波的幅度；τ為回波信號時延。

設定第個距離單元混頻的參考信號為

式中：＝2（-1）Δ／c表示時延量，Δ＝c／2為距離單元寬度；c為光速。

根據LFMCW雷達信號處理原理，將式（3）和式（4）中的回波信號與參考信號進行混頻、低通濾波得到差頻信號，之后再沿著快時間維進行快速傅里葉變換（fast Fourier transform，FFT），求得頻譜如下：

式中：＝T ＋-τ為有效時長；Sa（）＝sin（）／，R為該點目標與雷達的斜距。根據差頻信號頻譜特點可得在第（＝0，1）個子CPI內雷達收到該點目標的空時采樣數據，即

基于上述針對點目標的回波信號分析，再結合低空風切變的分布式特點，可得到第（＝0，1）個子CPI內第個距離單元低空風切變的回波數據s （）為

2 強雜波下基于CS的速度解模糊

2.1 CS基本原理

如果維度為的原始信號在一個正交變換基上的投影結果具有稀疏性，那么存在與不相關的一個測量矩陣，用其對信號進行線性測量可得到維度為的測量值，最后根據測量值建立優(yōu)化問題，恢復出原信號。CS的基本流程圖如圖3所示。

圖3 CS理論基本流程圖Fig.3 Basic flow chart of CS theory

式中：＝為CS矩陣。信號重構即為由觀測數據重構，或重構與之等價的。重構的過程可以轉化為如下的數學優(yōu)化問題：

2.2 雜波抑制加權CS解速度模糊

雜波抑制加權CS解速度模糊方法的核心思想是利用CS從非均勻欠采樣的數據中恢復出均勻完整的不模糊數據。首先需要利用多重PRF下一個CPI內回波的非均勻欠采樣特性構造出解模糊的冗余字典；然后將冗余字典代入最小范數優(yōu)化問題中提取強雜波的若干稀疏系數，估計雜波能量支撐區(qū)并據此構造出雜波抑制矩陣；接著用該雜波抑制矩陣對恢復向量加權，即構造出加權的最小范數優(yōu)化問題，實現雜波抑制與低空風切變在角度 多普勒域中均勻完整的無模糊稀疏向量恢復，得到準確的多普勒信息；最后根據多普勒頻率與速度的關系求得準確的速度值。下面針對方法的四個關鍵步驟分別進行闡述。

2.2.1 角度-多普勒域的無模糊感知矩陣構建

當采用本文第1節(jié)所述的脈組參差方式進行模糊問題解決時，設第（＝0，1）重頻f 對應的多普勒維采樣間隔為（）＝1／f ，將重頻f 發(fā)射的第1個調頻周期作為起始調頻周期，則可將一個完整CPI內每個調頻周期所對應的多普勒維采樣時刻表示為

式中：int（·）為向下取整運算；mod（·）為求余運算。由式（11）分析可知，采樣時刻，，…，t 之間是等間隔的，采樣時刻t ，t ，…，t 之間也是等間隔的，這是因為在各子CPI內重頻是固定的，但是對于一個完整的CPI而言，由于采用了兩個不同的重頻，則全部采樣時刻，，…，t 是非等間隔的。因此，對于第個距離單元，一個總CPI的回波信號（）顯然具有非均勻采樣的特性。

若雷達回波可能的最大多普勒頻率為f ，要想恢復出0～f 范圍內的無模糊的多普勒信息，根據Nyquist采樣定理，要求PRF值f ≥f ，并且該PRF值遠大于重頻參差中的PRF值，即f ?f 。此時在～t 的時間范圍內，當采用一重頻f 對信號進行等間隔均勻采樣時，全采樣數應滿足

且有?。由此可見，重頻脈組參差方法下得到的采樣點數顯然不能滿足全采樣的要求，即一個總CPI的回波信號（）是欠采樣的。

雷達探測低空風切變時，回波信號在角度 多普勒域中只分布在很小范圍，即在角度 多普勒域具有高度稀疏性，這為CS的使用提供了條件。本論文首先將無模糊多普勒頻率范圍0～f 均勻劃分為份，劃分間隔為Δ≈f ／，則可得到一系列無模糊頻率劃分位置，，…，f ，接著建立基于不模糊多普勒頻率的時域字典A ：

可以看到，A 的行代表某一采樣時刻目標所有可能的不模糊多普勒分布，列代表某一多普勒頻率處信號的有限個非均勻采樣時刻。

之后再將角度軸均勻劃分為N 份，建立基于不同角度的空域字典：

將時域字典A 與空域字典A 進行Kronecker積構造得到×YN 維的CS字典矩陣：

此時角度-多普勒域被劃分為均勻的網格點，每個網格對應著不同的空間角度及多普勒頻率。

由于重頻參差下的采樣點數遠小于不模糊時均勻完整的全采樣點數，所以由重頻參差下的回波數據重構不模糊的目標信息時所用感知矩陣是欠定的。根據文獻［14］的相關論述，感知矩陣具有較好的限制等距性質，因此原始信號能夠以很大的概率從重頻參差的雷達回波中正確的恢復出來。

2.2.2 雜波抑制矩陣構建

得到感知矩陣之后，低空風切變解模糊問題可以表示為如下所示的優(yōu)化問題：

然而回波信號中必定存在強地雜波，如果對上式直接進行求解，通常只能提取出方向圖主瓣附近的較大雜波。因此需要抑制雜波，方可進行低空風切變的解模糊。

2.2.3 基于加權范數優(yōu)化模型的無模糊稀疏向量恢復

根據所構建的雜波抑制矩陣可以將式（16）轉化為如下形式：

2.2.4 低空風切變無模糊速度估計

可以計算出第個距離單元不模糊的風場速度值v 。

對個距離單元依次使用上述方法進行處理，最終得到整個風場的不模糊速度。

3 算法流程

強雜波下基于CS的低空風切變速度解模糊方法的流程圖如圖4所示。

圖4 強雜波下基于CS的風速解模糊流程圖Fig.4 Flow chart of wind speed ambiguity resolution based on CS under strong clutter

4 仿真結果及分析

4.1 仿真條件

為驗證所提方法的有效性，進行了相關的仿真實驗，所用參數如表1所示。仿真時采用兩重頻的脈組參差方式，在一個CPI內先以重頻7 000 Hz發(fā)射32個脈沖，再以8 000 Hz發(fā)射32個脈沖。

表1 系統(tǒng)仿真參數Table 1 System simulation parameters

4.2 仿真結果分析

圖5仿真了PRF為7 000 Hz時雷達回波信號的空時二維譜?？梢姡仉s波的功率譜擴展到大于PRF的區(qū)域，又折疊進入可觀測的多普勒空間，風場的空時二維譜呈帶狀。由于強地雜波回波的存在，風場已經難以檢測。

圖5 回波空時譜Fig.5 Space time spectrum of echo

圖6為重頻7 000 Hz時，仿真得到的低空風切變信號和地雜波信號的距離-多普勒譜。從圖6（a）中可以明顯看到，風切變信號主要存在于第57號～第110號距離單元內（對應的距離范圍為8.5～16.5 km），其多普勒頻率有正有負，且多普勒頻率隨距離單元的變化已不能呈現出低空風切變具有的典型反S形分布的特點，說明低空風切變信號出現了多普勒模糊。圖6（b）為地雜波距離 多普勒譜，可以發(fā)現其在第57號～第110號距離單元內接近零頻，且在近距離處變化明顯，說明近程雜波具有距離依賴性。

圖6 雷達回波信號距離 多普勒譜Fig.6 Range-Doppler spectrum of radar echo signal

圖7是采用式（17）估計的不同距離單元雜波能量分布軌跡圖，可以看到估計出的雜波軌跡與真實雜波軌跡基本一致，這將為后續(xù)的雜波抑制和不模糊低空風切變檢測打下良好基礎。另外，從仿真結果來看，距離較近的第5號單元和距離較遠的第66號單元雜波的橢圓軌跡有明顯的變化，這正是雜波在空時平面內能量分布在近距離單元處具有距離依賴性的表現。

圖7 不同距離單元CS方法估計的雜波能量支撐區(qū)Fig.7 Clutter energy support region estimated by different range unit CS methods

圖8分別是采用雜波抑制加權CS估計得到的不同距離單元目標的空時二維譜。觀察仿真結果，計算可得第68號距離單元對應的速度為48.5 m／s，估計結果與該距離單元的真實風速47.6 m／s非常接近，而第78號距離單元的風速估計結果為15 m／s（真實風速為17.1 m／s），同理可以計算出其他距離單元的風速?？梢姴捎帽疚姆椒ㄒ院?，可以準確的進行雜波抑制，并恢復出不模糊的目標的速度值。

圖8 不同距離單元解模糊后風場空時譜Fig.8 Space time spectrum of wind field after defuzzification with different distance units

圖9分別是得到的風場多普勒估計結果圖、風速估計結果圖。首先，多普勒頻率的估計結果隨著距離變化呈現出反S形分布的特點，說明此時已經恢復出了不模糊的多普勒信息，之后根據多普勒頻率與速度之間的關系可以計算出風場速度。從圖中可以看到風速估計結果與原始速度分布較吻合，說明風速估計結果較好。

圖9 雜波抑制加權CS解模糊后風場估計結果Fig.9 Wind field estimation results after clutter suppression weighted CS difuzzification

圖10對比了所提解模糊方法與其他解模糊方法的估計結果，可以看到剩余定理法解模糊后的估計結果與原始風速分布基本吻合，個別距離處的風速估計結果存在偏差，一維集算法和本文所提方法的估計結果都可以較好地與原始風速分布進行吻合。為了更好地對比不同方法的估計精度，分別計算了對應的均方根誤差，如表2所示，可以看到剩余定理法對應的均方根誤差大于一維集算法與所提方法的均方根誤差，說明一維集算法與所提方法具有較高的估計精度。

圖10 不同方法解模糊后風速估計結果對比Fig.10 Comparison of wind speed estimation results after different difuzzification methods

表2 不同方法估計結果的均方根誤差比較Table 2 Comparison of root mean square error of different estimation methods

表3給出了不同方法的運行時間對比。分析可知表中前兩種方法需要在兩重頻下分別使用組合空時主通道自適應處理（combined space-time main channel adaptive processing，CMCAP）抑制雜波、匹配風場信號，此時不可避免要求解雜波協方差矩陣及其逆矩陣，之后還要求解出兩重頻下的模糊速度值，再來解模糊，導致算法具有較高計算量及復雜度。而本文所提方法通過構造雜波抑制矩陣對恢復向量加權，求解加權的優(yōu)化問題，同時實現雜波抑制與不模糊的風場速度估計，計算量顯著降低，運行時間少。綜合上述分析可以看出，所提方法具有較好的估計精度，且計算量小，相較常用方法具有更強的實用性。

表3 不同方法運行時間對比Table 3 Comparison of different methods’running time

5 結 論

針對毫米波LFMCW雷達檢測低空風切變時存在的強雜波干擾以及風速估計模糊問題，本文采用CS進行低空風切變速度解模糊。首先在整個CPI內采用重頻脈組參差方式發(fā)射信號，得到非均勻欠采樣的回波數據，之后在角度 多普勒域構造冗余字典，利用CS技術恢復出大的雜波譜成分，估計雜波能量的支撐區(qū)域，構造雜波抑制矩陣對恢復向量加權，求解加權的最小范數優(yōu)化問題實現雜波抑制與低空風切變的速度解模糊。仿真結果表明所提方法能實現較好的解模糊效果，且計算量小，具有較強實用性。