基于激光雷達的斜坡改進風切變預警算法

張千千,沈 淳,高 航,李健兵

(國防科技大學 電子科學學院 復雜電磁環(huán)境效應國家重點實驗室,湖南 長沙 410073)

1 引 言

低空風切變指的是在距地面600 m范圍內,風矢量(風速、風向)在水平或垂直距離上發(fā)生明顯改變的現(xiàn)象[1],具有時間短、類型多變、破壞力強等特點,且常伴隨著有極端天氣產(chǎn)生,對飛行器的起飛和降落危害極大[2]。根據(jù)國際民航組織(International Civil Organization,ICAO)的規(guī)定:風切變強度分為輕度、中度、重度、嚴重4個等級[2],可用風切變因子來表示,風切變因子是指空間兩點距離(水平距離或垂直距離)為30 m時風矢量差的大小[2]。低空風切變包括風的水平切變、風的垂直切變和垂直氣流的切變[2-3],主要由雷暴天氣、鋒面運動、逆溫現(xiàn)象和地形、地標建筑物影響等[4]造成。按照風場對飛機飛行姿態(tài)的影響,在航空氣象學中把低空風切變分為順風切變、逆風切變、側風切變和垂直氣流切變[5]。由于大部分飛行事故都發(fā)生在飛機起飛和進近階段[4],因此,機場的風切變預警測技術尤為重要[4]。

現(xiàn)有風場探測設備主要包括地表風速計、風廓線雷達(WPR)和多普勒天氣雷達(TDWR)等,在一定條件下為飛機起降提供了護航服務,但這些設備存在一定的缺陷。地表風速計探測高度有限,僅能探測水平風切變,且探測結果和儀器數(shù)量有關;風廓線雷達只能探測頂空風場情況,不能進行下滑道掃描,難以預警飛機進近走廊的低空風切變;多普勒天氣雷達探測范圍大,但其距離分辨率不夠精細,不能滿足小尺度低空風切變的探測需求[6-7]。相干激光測風雷達作為一種新型的探測裝置,體積小、重量輕、抗干擾能力強、分辨率高,能夠彌補其他探測手段的不足,尤其是在晴空條件下成為最有效地風場探測手段,對保障飛機起降安全、預警低空風切變有重要指導作用[8]。因此,利用激光測風雷達,準確地對低空風場預測逐漸成為保證飛行的重要工具。

自二十世紀70年代,著名氣象學家FUJITA開啟風切變的研究[9]。1983年,WOODFIELD和WOODS提出S因子算法[10],但這種僅通過求最大最小值之間風速差的算法極易產(chǎn)生較大誤差。2011年,香港天文臺陳柏瑋等人提出香港機場沿用至今的單斜坡低空風切變自動預警算法[11],但是該算法受到斜坡長度的限制僅能探測600 m尺度的風切變。2012年,CHEN提出F因子算法[12],但直接計算由激光雷達得到的逆風梯度數(shù)據(jù)會導致結果快速波動且無法計算垂直分量的風切變。2014年,他們又提出利用渦流耗散率[13]預警低空風切變,但因為渦流耗散率的閾值在國際中沒有明確的規(guī)定,該算法還需要在實際應用中繼續(xù)探索。2016年,JIANG等人提出單雙斜坡結合的算法[14]來預警風切變,不足的是該算法的目的是預警大尺度風切變。2019年,MA等人提出基于主成分分析(principal component analysis,PCA)和相位差校正的風切變預警算法[15],但是該算法的多普勒頻移和風切變閾值還需要進一步通過經(jīng)驗確定。鑒于現(xiàn)有算法存在不足,本文中提出了一種基于斜坡改進的激光雷達風切變預警算法,將重組逆風廓線劃分成不同的斜坡區(qū)間,再通過確定每個斜坡上的能量變化判斷是否有風切變發(fā)生,最終確定風切變發(fā)生的位置和時刻,實現(xiàn)風切變的穩(wěn)健預警。

2 原 理

2.1 激光雷達探測場景

激光雷達探測風切變主要采用下滑道掃描模式,該方式由香港天文臺陳柏瑋等人提出[10,16-18],利用同時改變激光雷達仰角和方位角的方式掃描飛機的進近走廊,直接探測低空風切變。傳統(tǒng)的PPI掃描完成一個周期的時間約為6 min,而下滑道掃描一個周期得到徑向數(shù)據(jù)僅需要2 min,可為飛機起降提供實時預警服務,保障飛行安全。在下滑道掃描過程中,激光光束與飛機起降跑道的夾角不小于30°。飛機起飛和降落時跑道與下滑道的夾角不同,起飛離場時夾角為6°,降落進場時夾角為3°。參考文獻[19]提供了香港國際機場的地理位置。圖1是激光雷達下滑道掃描示意圖[10],其中圖1(a)為降落場景,圖1(b)為香港機場下滑道與激光雷達雷達對應位置關系的三維圖。

圖1 激光雷達下滑道掃描示意圖

假設飛機的下滑道是被一條一定尺度的管道包裹,下滑道處在管道的中央,那么重組逆風廓線的徑向數(shù)據(jù)就是落在管道中的數(shù)據(jù)[20]。為了有效提取管道內數(shù)據(jù),以雷達為原點建立坐標系,利用笛卡爾坐標系將徑向風速進行轉換。以飛機降落為例,此時跑道與下滑道夾角為3°,將管道在笛卡爾坐標系三個平面投影,并根據(jù)激光波束與跑道夾角小于30°的約束條件,可得到x軸,y軸,z軸的坐標范圍:

(1)

式中,d是激光雷達與跑道的垂直距離;L是激光雷達到跑道的垂足與著陸點之間的距離;α是激光光束與跑道方向的夾角;h為下滑道上徑向數(shù)據(jù)點到y(tǒng)軸的垂足距離h=tan3°(y-l)。

2.2 算法原理

F因子算法是航空安全中使用較多的風切變預警方法,但其計算結果在某些情況下波動過于劇烈,長期存在虛警率過高問題;單斜坡算法可以實現(xiàn)風切變的快速檢測,但受到斜坡長度的限制對持續(xù)時間較長的風切變檢測結果不理想。基于現(xiàn)有算法存在的問題,本文提出一種將F因子與斜坡算法結合的改進斜坡算法。

該區(qū)域主要分布在流域的西部、西南部和東北部，海拔相對較高，雖然降水量較大，但降雨時空分布極為不均，且由于地勢的原因，土壤蓄水能力差、地表水資源利用率低、農田水利設施落后，是流域內最容易受旱的地區(qū)，也是農業(yè)干旱脆弱性最嚴重的地區(qū)。經(jīng)過對所搜集資料的綜合分析，該區(qū)域可分為東北沂蒙山丘區(qū)、豫西山丘區(qū)和淮南山丘區(qū)等3個二級區(qū)。

圖2為整個算法的流程圖,具體步驟如下:

圖2 斜坡改進算法的流程圖

①確定重組逆風廓線并進行信號處理。首先篩選數(shù)據(jù)構造逆風廓線;再根據(jù)角度得到下滑道上的數(shù)據(jù)并進行質量控制得到重組逆風廓線,隨后延長重組逆風廓線的尾部并進行一次滑動平均,最終得到?jīng)]有數(shù)據(jù)缺失、尾部延長且保留特征的重組逆風廓線。

②根據(jù)風速差求每個斜坡上的能量變化。首先根據(jù)經(jīng)典斜坡算法將重組逆風廓線劃分成長度不同的斜坡,再利用每個斜坡上最大最小值的風速差確定斜坡的能量變化,并與閾值比較。

③發(fā)出預警或預警解除。主要是判斷每個斜坡的能量是否超過閾值,超過閾值,發(fā)出風切變預警,此時該斜坡對應的位置即為風切變發(fā)生的位置;未超過閾值,則警報解除。如果計算結果遠超閾值,還應將該斜坡提取出來,并進一步確定風切變的位置。

本文提出算法的第一個創(chuàng)新點是區(qū)別于傳統(tǒng)方法通過計算相鄰兩點之間的能量變化進行預警,在得到每個斜坡上的速度變化之后,求得斜坡上的速度差,再計算每個斜坡上能量的變化情況。當計算結果超過閾值±0.05之后,發(fā)出風切變警告。此時,計算結果不再劇烈波動,虛警率也顯著下降。其中,每個斜坡的能量變化率表示為:

(2)

式中,ΔVm表示斜坡上最大風速和最小風速之差;D為斜坡長度;g表示重力加速度,而Vapp表示飛機的進近速度。

經(jīng)典F因子算法是求取重組逆風廓線上相鄰兩點之間的能量變化,并與閾值判斷,從而給出風切變預警信息,如圖3中Fc所示。而本文提到的斜坡改進算法中,將F因子與斜坡算法結合后,求取每個斜坡上能量的變化,然后與閾值判斷再給出風切變預警信息,如圖3中FI所示。

圖3 斜坡改進算法與經(jīng)典F因子算法區(qū)別

斜坡算法的本質是將重組逆風廓線劃分成不同的速度變化區(qū)間,改進斜坡算法在得到重組逆風廓線上不同斜坡之后,再比較每個斜坡上的能量變化,避免了相鄰兩風速快速變化造成的F因子算法快速波動,降低了F因子算法的虛警率,也將風切變發(fā)生的位置直接定位到某個斜坡。當某個斜坡上速度發(fā)生突變,必然會導致該斜坡的能量發(fā)生突變,若超出閾值很多,此時將該斜坡提取出來,可以準確判斷風切變發(fā)生的具體位置。

本文提出算法的第二個創(chuàng)新點是用最大最小值之差代替單斜坡算法中斜坡首尾的風速差。在單斜坡算法中,有兩種情況需要對斜坡進行延長:一是斜坡持續(xù)增速,但跨過斜坡后降速;二是斜坡持續(xù)降速,但跨過斜坡后增速[11],這兩種情況下,最長可延長斜坡至原始斜坡長度的1/2。第二種情況如圖4中的局部重組逆風廓線所示。這兩種情況在實際中很普遍,此時若還用單斜坡算法中斜坡首尾兩端的風速差作為斜坡速度的變化,顯然小于最大最小值的風速差,進而可能小于單斜坡算法中7.7 m/s的閾值,產(chǎn)生漏報。當斜坡需要收縮[11]的時候,最大最小值風速差與斜坡首尾風速差一致。利用最大最小值之差作為斜坡上速度的變化,當最大值出現(xiàn)在最小值之后,用最大值減最小值得到正的速度變化,表明該斜坡上風速增加,且風向是靠近激光雷達;當最大值出現(xiàn)在最小值之前,用最小值減最大值得到負的速度變化,表明該斜坡上風速減小,且風向是遠離激光雷達。圖4展現(xiàn)了原始斜坡長度取400 m,延長至600 m時,這一個斜坡上最大最小值風速差與單斜坡算法中首尾風速差的區(qū)別。

圖4 斜坡改進算法與單斜坡算法區(qū)別

尾部延長的意義是防止最后一部分數(shù)據(jù)不滿足最小原始斜坡長度,不能進行斜坡長度的確定造成信息缺失。尾部延長有兩種形式,尾部對稱性擴展法和尾部直接延長法[21]。尾部對稱性擴展法是選取重組逆風廓線最后一部分數(shù)據(jù)做對稱延長,一般選取的長度為5～10個數(shù)據(jù)點,當最后一段數(shù)據(jù)呈線性變化時,會人為造成最后一個點的突變;尾部直接延長法是將重組逆風廓線尾部最后一個元素直接擴展,擴展的長度至少為一個原始斜坡長度對應的數(shù)據(jù)點個數(shù),該方法保證在對重組逆風廓線最后一段數(shù)據(jù)檢測時,不會造成尾部缺失,故本文采用尾部直接延長法。

3 方法驗證、對比與分析

2020年開始,香港國際機場應用激光雷達對飛機起降下滑道進行了全天候掃描,場景如圖1所示。結合飛行員報告風切變發(fā)生的高度及風切變形成的原因,本文選取三組比較有代表性的數(shù)據(jù)對斜坡改進算法進行驗證。由于香港國際機場產(chǎn)生風切變的主要原因多為海風鋒,所以F因子閾值取±0.05[22-23],而非國際民航組織規(guī)定的±0.105,本文中對比的F因子、斜坡改進算法閾值亦取±0.05。

3.1 案例一:海風鋒情況

2020年2月24日,北京時間6:15、6:28、6:31在07R跑道有三趟飛機在進場降落報告有風切變發(fā)生,報警高度分別為300 m、200 m和100 m,塔臺分析是由海風鋒形成的風切變。論文提取了2020年2月24日6時32′49″的激光雷達掃描圖,利用質量控制和信號處理方法得到重組逆風廓線,如圖5所示,其中5(a)是激光雷達掃描多普勒速度分布圖,正值代表風場方向靠近激光雷達,負值代表風場方向遠離激光雷達,07R/25L、07L/25R表示不同方向的兩條跑道(虛線部分),跑道上黑色實線表示3°仰角的下滑道,飛機在07R/25L跑道降落。

圖5 2020年2月24日6:32:49時風場信息

圖5(b)是提取的重組逆風廓線及對應的下滑道高度信息,橫坐標為下滑道上距離降落點的距離,單位為km;左邊縱坐標表示重組逆風廓線對應的風速(實線部分),數(shù)值為正表示風向靠近激光雷達,數(shù)值為負表示風向遠離激光雷達,單位為m/s;右邊縱坐標表示下滑道距地面的高度(虛線部分),單位為m,下滑道的仰角為3°。

使用改進算法對提取的重組逆風廓線進行判斷,斜坡原始長度依次取400 m、800 m、1600 m、3200 m,檢測結果如圖6所示。由圖6可知,取四種不同的原始斜坡長度,均能檢測到風切變,只有實線所示的400 m原始斜坡檢測出4處風切變,對應的橫坐標依次為x1=3.7 km、x2=4.3 km、x3=4.9 km、x4=5.4 km,對應的下降高度分別為h1=193.91 m、h2=225.35 m、h3=256.80 m和h4=283.00 m,與飛行員報告的高度范圍一致。虛線所示的800 m原始斜坡檢測出2處風切變,對應的橫坐標為x1=3.8 km、x2=5 m,對應高度h1=199.15 m、h2=262.04 m,與飛行員報告一致。點線和點劃線均只檢測出一處風切變,橫坐標相同為x=3.1 km,對應的高度為162.46 m。此次實驗不同長度的原始斜坡均能檢測到風切變,但400m原始斜坡檢測結果覆蓋更廣,也與飛行員報告的高度更貼近。

圖6 斜坡改進算法檢測結果

再將斜坡改進算法與F因子、單斜坡算法進行對比,結果如圖7所示。圖中橫坐標表示下滑道上距離降落點的距離,單位是km;左邊縱坐標表示檢測結果,右邊縱坐標表示下滑道上風速,后文圖10、圖13坐標意義相同。圖7(a)是F因子算法結果,圖7(b)是原始斜坡取400 m的改進算法結果,圖7(c)是原始斜坡取400 m的單斜坡算法,圖中灰色陰影部分是風切變發(fā)生的范圍,點線是信號處理后的下滑道重組逆風廓線。根據(jù)圖7展示的結果,F因子波動較大,出現(xiàn)14處風切變預警,風切變的范圍在3.9～6.6 km,對應高度為204.39～345.89 m。斜坡改進算法有4處風切變,對應高度為193.91～283.00 m。單斜坡算法沒有檢測到風切變。由于飛行員報告的高度是100～300 m,因此本文提出斜坡改進算法較F因子算法給出的范圍更精確。

圖7 各種算法檢測結果對比

3.2 案例二:陣風鋒、對流、雷雨大風情況

2020年8月2日,北京時間9:07和9:23北京飛香港的KA901航班、悉尼飛香港的CX162航班在07L/25R跑道降落時報告有風切變,報警高度為500 m和200 m,符合低空風切變的高度。塔臺給出的風切變成因是:陣風鋒、對流和雷雨大風。利用質量控制及信號處理方法得到重組逆風廓線,如圖8所示,其中8(a)是激光雷達掃描多普勒速度分布圖,圖8(b)是提取的重組逆風廓線及對應的下滑道高度信息。

圖8 2020年8月2日9:22:57時風場信息

首先開展不同原始長度的斜坡改進算法對比,對比結果如圖9所示。原始斜坡取400 m、800 m、1600 m、3200 m檢測到風切變次數(shù)依次是:6、2、1、0,由重組逆風廓線可知,在距離降落點2～4 km范圍,風速變化雖然沒有很大,但風向發(fā)生了改變;在距離降落點4～7 km范圍,風速增加到最大,風向也由正向變?yōu)樨撓蛟僮優(yōu)檎?在7～8 km,風速從最大逐漸減小到0,再反向增加。按照判斷標準,這三個階段均是有風切變發(fā)生的,而原始斜坡取800 m、1600 m、3200 m的改進算法分別有不同層度的漏報,400 m的效果最好。

圖9 斜坡改進算法檢測結果

再進行不同算法直接的對比,圖10給出了對比的結果,圖10(a)是F因子算法結果,圖10(b)是原始斜坡取40 0 m的改進算法結果,圖10(c)是原始斜坡取400 m的單斜坡算法,圖中灰色陰影部分是風切變發(fā)生的范圍,綠線是信號處理后的下滑道重組逆風廓線。F因子從距離降落點0.3 km至距離降落點8.5 km均報告有風切變發(fā)生,對應的高度為15.72～445.47 m,檢測結果起伏劇烈,波動較大,涉及范圍較廣,夸大了風切變的程度。斜坡改進算法給出了三個位置有風切的預警報告,對應的高度范圍是:125.78～167.70 m、277.76 m、377.34～408.78 m。

圖10 各種算法檢測結果對比

相對于F因子的檢測結果,斜坡改進算法給出風切變發(fā)生的范圍更精確,呈現(xiàn)的結果也更清晰。單斜坡算法這次檢測出一處發(fā)生了風切變,對應的橫坐標為x=7.8 km,與斜坡改進算法中的最后一處一致。通過前面不同長度的改進算法對比分析,單斜坡算法產(chǎn)生了兩處漏報。

3.3 案例三:無風切變情況

為了進一步對比F因子算法和改進斜坡算法的區(qū)別,我們找到一組幾天內沒有明顯發(fā)生風切變、也沒有飛行員報告的數(shù)據(jù)進行驗證。

圖11是北京時間2020年4月1日08:08:22時,07L/25R跑道的激光雷達掃描多普勒速度分布圖,圖中最大負向風速約5 m/s,最大正向風速約6 m/s,提取下滑道重組逆風廓線后如圖11(b)所示。圖11(a)中風速看似起伏變換,但整個逆風廓線最大最小風速差為7.32 m/s,小于7.7 m/s的標準,確定明顯風切變發(fā)生。圖12是取不同長度的斜坡改進算法對比結果圖,圖13是不同算法之間的對比結果圖,其中圖13(a)是F因子算法結果,圖13(b)是原始斜坡取400 m的改進算法結果,圖13(c)是原始斜坡取400 m的單斜坡算法,圖中灰色陰影部分是風切變發(fā)生的范圍,虛線是信號處理后的下滑道重組逆風廓線。

圖11 2020年4月1日8:8:22時風場信息

圖12 斜坡改進算法檢測結果

如圖13所示,在確定沒有風切變發(fā)生的背景下,F因子檢測出4個范圍(6處)有風切變發(fā)生,虛警率過高,而不同長度的斜坡改進算法均沒有產(chǎn)生報警,且檢測結果遠小于報警閾值。案例三反應了F因子的不足。

將單斜坡算法的結果與斜坡改進算法的結果放在同一張圖上進行比較,如圖14所示,其中圖14(a)～(c)依次為案例一、案例二、案例三,圖中左邊縱坐標為斜坡改進算法,右邊縱坐標表示單斜坡算法,圖中所取的原始斜坡長度均為400 m。單斜坡算法和改進斜坡算法大部分走勢和形狀一致,但是兩種算法也有細微區(qū)別,圖中點線的圓圈部分展示了它們的差異。這是因為,在大部分時刻,改進斜坡算法中每個斜坡的最大最小風速與斜坡算法中首尾風速是一致的,但是在一些情況下,最大或最小風速就不出現(xiàn)在斜坡的首端或尾端,從而斜坡上風速差有細微區(qū)別,檢測結果在局部發(fā)生改變,使得斜坡的能量變化超過閾值發(fā)出預警信息,這也使得改進斜坡算法更優(yōu)于斜坡算法且彌補了斜坡算法產(chǎn)生漏報的不足。根據(jù)大量數(shù)據(jù)驗算的結果,圖14展示的有差異的情況很常見。

圖14 斜坡改進算法與單斜坡算法區(qū)別

3.4 數(shù)據(jù)統(tǒng)計分析

根據(jù)香港國際機場飛行員返回的報告,2020年全年機場有風切變報告的天數(shù)為74天。論文將所有發(fā)生風切變對應時刻的激光雷達掃描數(shù)據(jù)提取出來,進行風切變檢測,并對結果進行了統(tǒng)計,表1上半部分是有風切變發(fā)生時算法檢測結果的統(tǒng)計。由表1可知:F因子算法全部檢測出來,檢測率分別為100 %;單斜坡算法僅檢測到17次風切變,檢測率22.97 %,漏警率77.03 %;而對于斜坡改進算法,原始斜坡取400 m、800 m、1600 m和3200 m時,依次檢測到73、45、13次風切變,檢測率依次為98.65 %、60.81 %、17.57 %,漏警率依次為1.35 %、39.19 %、82.43 %。斜坡改進算法在取原始斜坡取400 m時檢測結果很好。

表1 各種算法檢測結果統(tǒng)計

再對沒有風切變報告的292天進行算法檢測。由于香港機場的激光測風雷達每3 min完成一次掃描,一小時能得到約20次下滑道數(shù)據(jù),每日隨機抽取連續(xù)2小時內共計40次下滑道數(shù)據(jù),因此共得到11680次無風切變報告的下滑道數(shù)據(jù)進行統(tǒng)計,檢測結果如表1下半部分所示。由表1可知,單斜坡算法和斜坡改進算法均沒有檢測到風切變,而F因子檢測出5022次風切變,此時虛警率高達43 %!

結合有風切變和無風切變的情況統(tǒng)計,F因子虛警率高,單斜坡算法漏警率高,而斜坡改進算法原始斜坡取400 m時,檢測率接近100 %,虛警率為0,有效彌補了傳統(tǒng)風切變檢測算法的不足。

斜坡改進算法將F因子和單斜坡算法結合起來,通過計算每個斜坡上飛機遭遇的能量變化來檢測風切變。不同原始長度的斜坡對比結果表明:斜坡長度取400 m結果最好。這是因此,當斜坡長度較長時,若風速變化不太劇烈,則每個斜坡的能量變化較小;而原始斜坡長度取400 m時,一個斜坡上最多僅有6個相鄰數(shù)據(jù),保證了計算的數(shù)量和質量。斜坡改進算法與經(jīng)典算法的對比結果表明:由于F因子是依次比較相鄰兩個數(shù)據(jù)點,使得結果波動劇烈,虛警過高;單斜坡算法是比較斜坡上首尾兩端的速度差,當大于7.7 m/s時才會發(fā)出報警,加大了漏警的概率。而斜坡算法恰好彌補了這兩種經(jīng)典算法的不足,且給出的風切變發(fā)生范圍更細致;利用斜坡上最大最小值作為風速差代替首尾兩端的風速差,使得風切變的計算結果也更精確。因此在后面的飛行保障中,可以采用原始斜坡取400 m的斜坡改進算法對風切變進行預警。

4 結 論

本文利用激光測風雷達在香港國際機場開展了風切變觀測試驗,提出了基于斜坡改進的激光雷達風切變預警算法,將F因子與單斜坡算法結合,計算每個斜坡上的能量變化,彌補F因子虛警較高和單斜坡算法漏警較高的不足;用最大最小值代替斜坡首尾的風速差,提高了檢測率和精度。通過對海風鋒和陣風鋒、對流、雷雨大風等情況下的風切變進行識別,并結合飛行員的報告對結果進行了分析,證實了斜坡改進算法的準確性和魯棒性,有利于風切變檢測,能為飛機起降提供技術保障。

基于激光雷達的斜坡改進算法,能有效彌補F因子虛警過高、單斜坡算法容易產(chǎn)生漏報等不足,正確預警飛機下滑道上是否有風切變發(fā)生及風切變發(fā)生的位置與范圍,為飛行安全提供有效保障。但是該算法對于風速增加緩慢的數(shù)據(jù)檢測不是很理想,后面將圍繞此類數(shù)據(jù)對算法進行優(yōu)化。