李 海,呼延澤,毛志杰
(1.中國(guó)民航大學(xué)電子信息與自動(dòng)化學(xué)院,天津 300300;2.國(guó)防科技大學(xué)信息通信學(xué)院,陜西 西安 710106)
低空風(fēng)切變是一種具有突發(fā)性、短持續(xù)、強(qiáng)破壞特點(diǎn)的危害性氣象現(xiàn)象,大多發(fā)生在低于600 m的低空區(qū)域。當(dāng)飛機(jī)飛行到該區(qū)域時(shí),極有可能被風(fēng)切變危害而引發(fā)慘重的空難事件。因此,準(zhǔn)確估計(jì)風(fēng)場(chǎng)速度,進(jìn)而探測(cè)到低空風(fēng)切變,并及時(shí)給飛機(jī)提供告警信息非常重要。
機(jī)載氣象雷達(dá)能夠幫助飛機(jī)實(shí)時(shí)監(jiān)視前方的氣象狀況,有助于飛行安全,其在工作時(shí)較多采用厘米波。而毫米波的波長(zhǎng)更短,波束更窄。機(jī)載氣象雷達(dá)發(fā)射毫米波,可以縮小天線,并且更加準(zhǔn)確地探測(cè)目標(biāo)。目前,機(jī)載氣象雷達(dá)普遍采用脈沖體制,而毫米波線性調(diào)頻連續(xù)波(linear frequency modulated continuous wave,LFMCW)體制信號(hào)處理系統(tǒng)結(jié)構(gòu)更為簡(jiǎn)單,并且重量輕、對(duì)發(fā)射峰值功率要求低。將LFMCW技術(shù)用于機(jī)載氣象雷達(dá),能夠使其向小型化發(fā)展,進(jìn)而提高載機(jī)的靈活性、機(jī)動(dòng)性。
檢測(cè)風(fēng)場(chǎng)時(shí),機(jī)載氣象雷達(dá)會(huì)嚴(yán)重被地雜波信號(hào)干擾,難以識(shí)別到目標(biāo)。此時(shí)常采用空時(shí)自適應(yīng)處理(space-time adaptive processing,STAP)技術(shù)壓制雜波。在此基礎(chǔ)上,已有多種風(fēng)速估計(jì)方法,但它們都未解決速度模糊問(wèn)題。毫米波雷達(dá)系統(tǒng)因其具有短波長(zhǎng)的特點(diǎn),很容易出現(xiàn)速度估計(jì)模糊,這給低空風(fēng)切變的檢測(cè)帶來(lái)阻礙。因此,必須考慮解模糊,目前主流的解模糊方法通常是以脈沖重復(fù)頻率組為基礎(chǔ)進(jìn)行的。文獻(xiàn)[8]提出了余差查表法,該算法利用目標(biāo)在各重頻上的余數(shù)之差進(jìn)行解模糊,其對(duì)存儲(chǔ)空間要求較大,并且在查表過(guò)程中會(huì)浪費(fèi)大量時(shí)間在不必要的模糊值上,實(shí)時(shí)性不高。文獻(xiàn)[9]提出剩余定理(孫子定理)法,通過(guò)對(duì)不同重頻下得到的同余方程組進(jìn)行求解來(lái)實(shí)現(xiàn)解模糊,該算法對(duì)測(cè)量精度要求較高,一旦存在測(cè)量誤差,解模糊的質(zhì)量難以保證。文獻(xiàn)[10]提出一維集算法,通過(guò)將所有可能的速度值從小到大排列,再尋找多普勒頻率組的最小均方值進(jìn)行解模糊,該方法存在計(jì)算量大的問(wèn)題。
速度模糊問(wèn)題,究其根本,是目標(biāo)的多普勒頻率大于了調(diào)頻周期重復(fù)頻率的一半,此時(shí)根據(jù)香農(nóng) 奈奎斯特采樣定理,會(huì)產(chǎn)生多普勒混疊??梢?jiàn),雷達(dá)系統(tǒng)可測(cè)的最大不模糊速度范圍是受香農(nóng) 奈奎斯特采樣定理限制的,近年來(lái)出現(xiàn)的壓縮感知(compressive sensing,CS)理論很好地突破了這一限制。因此,本文考慮用CS理論來(lái)解決速度模糊問(wèn)題。首先利用一個(gè)相干處理時(shí)間內(nèi)回波信號(hào)非均勻欠采樣特性在角度 多普勒域構(gòu)造無(wú)模糊冗余字典,利用CS技術(shù)提取雜波譜主要成分,估計(jì)雜波能量支撐域,構(gòu)造出雜波抑制矩陣。之后利用加權(quán)最小范數(shù)優(yōu)化模型實(shí)現(xiàn)雜波抑制與低空風(fēng)切變的速度解模糊。此方法只需對(duì)待測(cè)距離門數(shù)據(jù)進(jìn)行處理,無(wú)需像傳統(tǒng)STAP方法那樣估計(jì)雜波協(xié)方差矩陣,因此計(jì)算量大大減小。此外,該方法在角度 多普勒域進(jìn)行細(xì)分來(lái)構(gòu)造字典,因而對(duì)運(yùn)動(dòng)目標(biāo)有更高的角度和多普勒分辨率。
假設(shè)機(jī)載氣象雷達(dá)結(jié)構(gòu)如圖1所示。載機(jī)平臺(tái)距地面距離為,速度為,天線陣元數(shù)為,相鄰陣元距離=/2,為信號(hào)波長(zhǎng)。和分別對(duì)應(yīng)強(qiáng)地雜波的方位角和俯仰角,和φ分別對(duì)應(yīng)風(fēng)場(chǎng)的方位角和俯仰角。設(shè)定雷達(dá)在相干處理間隔內(nèi)發(fā)射個(gè)周期的調(diào)頻信號(hào)。
圖1 機(jī)載雷達(dá)示意圖Fig.1 Schematic diagram of airborne radar
為了解決毫米波段引起的測(cè)速模糊問(wèn)題,信號(hào)模型中采用兩重頻的脈組參差方式,兩重脈沖重復(fù)頻率(pulse repetition frequency,PRF)的脈組參差方式采樣時(shí)序圖如圖2所示。
圖2 兩重PRF的脈組參差方式Fig.2 Mode of pulse group stagger in double PRF
將一個(gè)完整的相干處理時(shí)間間隔(coherent processing interval,CPI)均分為兩個(gè)子CPI,設(shè)使用的兩重頻分別為f 、f ,則對(duì)應(yīng)的采樣間隔分別為T =1/f 、T =1/f 。假設(shè)機(jī)載氣象雷達(dá)在其工作范圍內(nèi)有個(gè)距離單元,第個(gè)距離單元的回波數(shù)據(jù)()=[(),()],其中,z ()為第(=0,1)個(gè)子CPI內(nèi)第個(gè)距離單元的空時(shí)二維快拍數(shù)據(jù):
式中:c ()是第(=0,1)個(gè)子CPI內(nèi)地雜波;n 為高斯白噪聲;s ()是低空風(fēng)切變信號(hào)。
令第(=0,1)個(gè)子CPI內(nèi)第(=1,2,…,/2)個(gè)周期的發(fā)射信號(hào)為
式中:ξ為回波的幅度;τ為回波信號(hào)時(shí)延。
設(shè)定第個(gè)距離單元混頻的參考信號(hào)為
式中:=2(-1)Δ/c表示時(shí)延量,Δ=c/2為距離單元寬度;c為光速。
根據(jù)LFMCW雷達(dá)信號(hào)處理原理,將式(3)和式(4)中的回波信號(hào)與參考信號(hào)進(jìn)行混頻、低通濾波得到差頻信號(hào),之后再沿著快時(shí)間維進(jìn)行快速傅里葉變換(fast Fourier transform,F(xiàn)FT),求得頻譜如下:
式中:=T +-τ為有效時(shí)長(zhǎng);Sa()=sin()/,R為該點(diǎn)目標(biāo)與雷達(dá)的斜距。根據(jù)差頻信號(hào)頻譜特點(diǎn)可得在第(=0,1)個(gè)子CPI內(nèi)雷達(dá)收到該點(diǎn)目標(biāo)的空時(shí)采樣數(shù)據(jù),即
基于上述針對(duì)點(diǎn)目標(biāo)的回波信號(hào)分析,再結(jié)合低空風(fēng)切變的分布式特點(diǎn),可得到第(=0,1)個(gè)子CPI內(nèi)第個(gè)距離單元低空風(fēng)切變的回波數(shù)據(jù)s ()為
如果維度為的原始信號(hào)在一個(gè)正交變換基上的投影結(jié)果具有稀疏性,那么存在與不相關(guān)的一個(gè)測(cè)量矩陣,用其對(duì)信號(hào)進(jìn)行線性測(cè)量可得到維度為的測(cè)量值,最后根據(jù)測(cè)量值建立優(yōu)化問(wèn)題,恢復(fù)出原信號(hào)。CS的基本流程圖如圖3所示。
圖3 CS理論基本流程圖Fig.3 Basic flow chart of CS theory
式中:=為CS矩陣。信號(hào)重構(gòu)即為由觀測(cè)數(shù)據(jù)重構(gòu),或重構(gòu)與之等價(jià)的。重構(gòu)的過(guò)程可以轉(zhuǎn)化為如下的數(shù)學(xué)優(yōu)化問(wèn)題:
雜波抑制加權(quán)CS解速度模糊方法的核心思想是利用CS從非均勻欠采樣的數(shù)據(jù)中恢復(fù)出均勻完整的不模糊數(shù)據(jù)。首先需要利用多重PRF下一個(gè)CPI內(nèi)回波的非均勻欠采樣特性構(gòu)造出解模糊的冗余字典;然后將冗余字典代入最小范數(shù)優(yōu)化問(wèn)題中提取強(qiáng)雜波的若干稀疏系數(shù),估計(jì)雜波能量支撐區(qū)并據(jù)此構(gòu)造出雜波抑制矩陣;接著用該雜波抑制矩陣對(duì)恢復(fù)向量加權(quán),即構(gòu)造出加權(quán)的最小范數(shù)優(yōu)化問(wèn)題,實(shí)現(xiàn)雜波抑制與低空風(fēng)切變?cè)诮嵌?多普勒域中均勻完整的無(wú)模糊稀疏向量恢復(fù),得到準(zhǔn)確的多普勒信息;最后根據(jù)多普勒頻率與速度的關(guān)系求得準(zhǔn)確的速度值。下面針對(duì)方法的四個(gè)關(guān)鍵步驟分別進(jìn)行闡述。
2.2.1 角度-多普勒域的無(wú)模糊感知矩陣構(gòu)建
當(dāng)采用本文第1節(jié)所述的脈組參差方式進(jìn)行模糊問(wèn)題解決時(shí),設(shè)第(=0,1)重頻f 對(duì)應(yīng)的多普勒維采樣間隔為()=1/f ,將重頻f 發(fā)射的第1個(gè)調(diào)頻周期作為起始調(diào)頻周期,則可將一個(gè)完整CPI內(nèi)每個(gè)調(diào)頻周期所對(duì)應(yīng)的多普勒維采樣時(shí)刻表示為
式中:int(·)為向下取整運(yùn)算;mod(·)為求余運(yùn)算。由式(11)分析可知,采樣時(shí)刻,,…,t 之間是等間隔的,采樣時(shí)刻t ,t ,…,t 之間也是等間隔的,這是因?yàn)樵诟髯覥PI內(nèi)重頻是固定的,但是對(duì)于一個(gè)完整的CPI而言,由于采用了兩個(gè)不同的重頻,則全部采樣時(shí)刻,,…,t 是非等間隔的。因此,對(duì)于第個(gè)距離單元,一個(gè)總CPI的回波信號(hào)()顯然具有非均勻采樣的特性。
若雷達(dá)回波可能的最大多普勒頻率為f ,要想恢復(fù)出0~f 范圍內(nèi)的無(wú)模糊的多普勒信息,根據(jù)Nyquist采樣定理,要求PRF值f ≥f ,并且該P(yáng)RF值遠(yuǎn)大于重頻參差中的PRF值,即f ?f 。此時(shí)在~t 的時(shí)間范圍內(nèi),當(dāng)采用一重頻f 對(duì)信號(hào)進(jìn)行等間隔均勻采樣時(shí),全采樣數(shù)應(yīng)滿足
且有?。由此可見(jiàn),重頻脈組參差方法下得到的采樣點(diǎn)數(shù)顯然不能滿足全采樣的要求,即一個(gè)總CPI的回波信號(hào)()是欠采樣的。
雷達(dá)探測(cè)低空風(fēng)切變時(shí),回波信號(hào)在角度 多普勒域中只分布在很小范圍,即在角度 多普勒域具有高度稀疏性,這為CS的使用提供了條件。本論文首先將無(wú)模糊多普勒頻率范圍0~f 均勻劃分為份,劃分間隔為Δ≈f /,則可得到一系列無(wú)模糊頻率劃分位置,,…,f ,接著建立基于不模糊多普勒頻率的時(shí)域字典A :
可以看到,A 的行代表某一采樣時(shí)刻目標(biāo)所有可能的不模糊多普勒分布,列代表某一多普勒頻率處信號(hào)的有限個(gè)非均勻采樣時(shí)刻。
之后再將角度軸均勻劃分為N 份,建立基于不同角度的空域字典:
將時(shí)域字典A 與空域字典A 進(jìn)行Kronecker積構(gòu)造得到×YN 維的CS字典矩陣:
此時(shí)角度-多普勒域被劃分為均勻的網(wǎng)格點(diǎn),每個(gè)網(wǎng)格對(duì)應(yīng)著不同的空間角度及多普勒頻率。
由于重頻參差下的采樣點(diǎn)數(shù)遠(yuǎn)小于不模糊時(shí)均勻完整的全采樣點(diǎn)數(shù),所以由重頻參差下的回波數(shù)據(jù)重構(gòu)不模糊的目標(biāo)信息時(shí)所用感知矩陣是欠定的。根據(jù)文獻(xiàn)[14]的相關(guān)論述,感知矩陣具有較好的限制等距性質(zhì),因此原始信號(hào)能夠以很大的概率從重頻參差的雷達(dá)回波中正確的恢復(fù)出來(lái)。
2.2.2 雜波抑制矩陣構(gòu)建
得到感知矩陣之后,低空風(fēng)切變解模糊問(wèn)題可以表示為如下所示的優(yōu)化問(wèn)題:
然而回波信號(hào)中必定存在強(qiáng)地雜波,如果對(duì)上式直接進(jìn)行求解,通常只能提取出方向圖主瓣附近的較大雜波。因此需要抑制雜波,方可進(jìn)行低空風(fēng)切變的解模糊。
2.2.3 基于加權(quán)范數(shù)優(yōu)化模型的無(wú)模糊稀疏向量恢復(fù)
根據(jù)所構(gòu)建的雜波抑制矩陣可以將式(16)轉(zhuǎn)化為如下形式:
2.2.4 低空風(fēng)切變無(wú)模糊速度估計(jì)
可以計(jì)算出第個(gè)距離單元不模糊的風(fēng)場(chǎng)速度值v 。
對(duì)個(gè)距離單元依次使用上述方法進(jìn)行處理,最終得到整個(gè)風(fēng)場(chǎng)的不模糊速度。
強(qiáng)雜波下基于CS的低空風(fēng)切變速度解模糊方法的流程圖如圖4所示。
圖4 強(qiáng)雜波下基于CS的風(fēng)速解模糊流程圖Fig.4 Flow chart of wind speed ambiguity resolution based on CS under strong clutter
為驗(yàn)證所提方法的有效性,進(jìn)行了相關(guān)的仿真實(shí)驗(yàn),所用參數(shù)如表1所示。仿真時(shí)采用兩重頻的脈組參差方式,在一個(gè)CPI內(nèi)先以重頻7 000 Hz發(fā)射32個(gè)脈沖,再以8 000 Hz發(fā)射32個(gè)脈沖。
表1 系統(tǒng)仿真參數(shù)Table 1 System simulation parameters
圖5仿真了PRF為7 000 Hz時(shí)雷達(dá)回波信號(hào)的空時(shí)二維譜??梢?jiàn),地雜波的功率譜擴(kuò)展到大于PRF的區(qū)域,又折疊進(jìn)入可觀測(cè)的多普勒空間,風(fēng)場(chǎng)的空時(shí)二維譜呈帶狀。由于強(qiáng)地雜波回波的存在,風(fēng)場(chǎng)已經(jīng)難以檢測(cè)。
圖5 回波空時(shí)譜Fig.5 Space time spectrum of echo
圖6為重頻7 000 Hz時(shí),仿真得到的低空風(fēng)切變信號(hào)和地雜波信號(hào)的距離-多普勒譜。從圖6(a)中可以明顯看到,風(fēng)切變信號(hào)主要存在于第57號(hào)~第110號(hào)距離單元內(nèi)(對(duì)應(yīng)的距離范圍為8.5~16.5 km),其多普勒頻率有正有負(fù),且多普勒頻率隨距離單元的變化已不能呈現(xiàn)出低空風(fēng)切變具有的典型反S形分布的特點(diǎn),說(shuō)明低空風(fēng)切變信號(hào)出現(xiàn)了多普勒模糊。圖6(b)為地雜波距離 多普勒譜,可以發(fā)現(xiàn)其在第57號(hào)~第110號(hào)距離單元內(nèi)接近零頻,且在近距離處變化明顯,說(shuō)明近程雜波具有距離依賴性。
圖6 雷達(dá)回波信號(hào)距離 多普勒譜Fig.6 Range-Doppler spectrum of radar echo signal
圖7是采用式(17)估計(jì)的不同距離單元雜波能量分布軌跡圖,可以看到估計(jì)出的雜波軌跡與真實(shí)雜波軌跡基本一致,這將為后續(xù)的雜波抑制和不模糊低空風(fēng)切變檢測(cè)打下良好基礎(chǔ)。另外,從仿真結(jié)果來(lái)看,距離較近的第5號(hào)單元和距離較遠(yuǎn)的第66號(hào)單元雜波的橢圓軌跡有明顯的變化,這正是雜波在空時(shí)平面內(nèi)能量分布在近距離單元處具有距離依賴性的表現(xiàn)。
圖7 不同距離單元CS方法估計(jì)的雜波能量支撐區(qū)Fig.7 Clutter energy support region estimated by different range unit CS methods
圖8分別是采用雜波抑制加權(quán)CS估計(jì)得到的不同距離單元目標(biāo)的空時(shí)二維譜。觀察仿真結(jié)果,計(jì)算可得第68號(hào)距離單元對(duì)應(yīng)的速度為48.5 m/s,估計(jì)結(jié)果與該距離單元的真實(shí)風(fēng)速47.6 m/s非常接近,而第78號(hào)距離單元的風(fēng)速估計(jì)結(jié)果為15 m/s(真實(shí)風(fēng)速為17.1 m/s),同理可以計(jì)算出其他距離單元的風(fēng)速??梢?jiàn)采用本文方法以后,可以準(zhǔn)確的進(jìn)行雜波抑制,并恢復(fù)出不模糊的目標(biāo)的速度值。
圖8 不同距離單元解模糊后風(fēng)場(chǎng)空時(shí)譜Fig.8 Space time spectrum of wind field after defuzzification with different distance units
圖9分別是得到的風(fēng)場(chǎng)多普勒估計(jì)結(jié)果圖、風(fēng)速估計(jì)結(jié)果圖。首先,多普勒頻率的估計(jì)結(jié)果隨著距離變化呈現(xiàn)出反S形分布的特點(diǎn),說(shuō)明此時(shí)已經(jīng)恢復(fù)出了不模糊的多普勒信息,之后根據(jù)多普勒頻率與速度之間的關(guān)系可以計(jì)算出風(fēng)場(chǎng)速度。從圖中可以看到風(fēng)速估計(jì)結(jié)果與原始速度分布較吻合,說(shuō)明風(fēng)速估計(jì)結(jié)果較好。
圖9 雜波抑制加權(quán)CS解模糊后風(fēng)場(chǎng)估計(jì)結(jié)果Fig.9 Wind field estimation results after clutter suppression weighted CS difuzzification
圖10對(duì)比了所提解模糊方法與其他解模糊方法的估計(jì)結(jié)果,可以看到剩余定理法解模糊后的估計(jì)結(jié)果與原始風(fēng)速分布基本吻合,個(gè)別距離處的風(fēng)速估計(jì)結(jié)果存在偏差,一維集算法和本文所提方法的估計(jì)結(jié)果都可以較好地與原始風(fēng)速分布進(jìn)行吻合。為了更好地對(duì)比不同方法的估計(jì)精度,分別計(jì)算了對(duì)應(yīng)的均方根誤差,如表2所示,可以看到剩余定理法對(duì)應(yīng)的均方根誤差大于一維集算法與所提方法的均方根誤差,說(shuō)明一維集算法與所提方法具有較高的估計(jì)精度。
圖10 不同方法解模糊后風(fēng)速估計(jì)結(jié)果對(duì)比Fig.10 Comparison of wind speed estimation results after different difuzzification methods
表2 不同方法估計(jì)結(jié)果的均方根誤差比較Table 2 Comparison of root mean square error of different estimation methods
表3給出了不同方法的運(yùn)行時(shí)間對(duì)比。分析可知表中前兩種方法需要在兩重頻下分別使用組合空時(shí)主通道自適應(yīng)處理(combined space-time main channel adaptive processing,CMCAP)抑制雜波、匹配風(fēng)場(chǎng)信號(hào),此時(shí)不可避免要求解雜波協(xié)方差矩陣及其逆矩陣,之后還要求解出兩重頻下的模糊速度值,再來(lái)解模糊,導(dǎo)致算法具有較高計(jì)算量及復(fù)雜度。而本文所提方法通過(guò)構(gòu)造雜波抑制矩陣對(duì)恢復(fù)向量加權(quán),求解加權(quán)的優(yōu)化問(wèn)題,同時(shí)實(shí)現(xiàn)雜波抑制與不模糊的風(fēng)場(chǎng)速度估計(jì),計(jì)算量顯著降低,運(yùn)行時(shí)間少。綜合上述分析可以看出,所提方法具有較好的估計(jì)精度,且計(jì)算量小,相較常用方法具有更強(qiáng)的實(shí)用性。
表3 不同方法運(yùn)行時(shí)間對(duì)比Table 3 Comparison of different methods’running time
針對(duì)毫米波LFMCW雷達(dá)檢測(cè)低空風(fēng)切變時(shí)存在的強(qiáng)雜波干擾以及風(fēng)速估計(jì)模糊問(wèn)題,本文采用CS進(jìn)行低空風(fēng)切變速度解模糊。首先在整個(gè)CPI內(nèi)采用重頻脈組參差方式發(fā)射信號(hào),得到非均勻欠采樣的回波數(shù)據(jù),之后在角度 多普勒域構(gòu)造冗余字典,利用CS技術(shù)恢復(fù)出大的雜波譜成分,估計(jì)雜波能量的支撐區(qū)域,構(gòu)造雜波抑制矩陣對(duì)恢復(fù)向量加權(quán),求解加權(quán)的最小范數(shù)優(yōu)化問(wèn)題實(shí)現(xiàn)雜波抑制與低空風(fēng)切變的速度解模糊。仿真結(jié)果表明所提方法能實(shí)現(xiàn)較好的解模糊效果,且計(jì)算量小,具有較強(qiáng)實(shí)用性。