基于激光雷达的斜坡改进风切变预警算法

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

(国防科技大学 电子科学学院 复杂电磁环境效应国家重点实验室,湖南 长沙 410073)

1 引 言

低空风切变指的是在距地面600 m范围内,风矢量(风速、风向)在水平或垂直距离上发生明显改变的现象[1],具有时间短、类型多变、破坏力强等特点,且常伴随着有极端天气产生,对飞行器的起飞和降落危害极大[2]。根据国际民航组织(International Civil Organization,ICAO)的规定:风切变强度分为轻度、中度、重度、严重4个等级[2],可用风切变因子来表示,风切变因子是指空间两点距离(水平距离或垂直距离)为30 m时风矢量差的大小[2]。低空风切变包括风的水平切变、风的垂直切变和垂直气流的切变[2-3],主要由雷暴天气、锋面运动、逆温现象和地形、地标建筑物影响等[4]造成。按照风场对飞机飞行姿态的影响,在航空气象学中把低空风切变分为顺风切变、逆风切变、侧风切变和垂直气流切变[5]。由于大部分飞行事故都发生在飞机起飞和进近阶段[4],因此,机场的风切变预警测技术尤为重要[4]。

现有风场探测设备主要包括地表风速计、风廓线雷达(WPR)和多普勒天气雷达(TDWR)等,在一定条件下为飞机起降提供了护航服务,但这些设备存在一定的缺陷。地表风速计探测高度有限,仅能探测水平风切变,且探测结果和仪器数量有关;风廓线雷达只能探测顶空风场情况,不能进行下滑道扫描,难以预警飞机进近走廊的低空风切变;多普勒天气雷达探测范围大,但其距离分辨率不够精细,不能满足小尺度低空风切变的探测需求[6-7]。相干激光测风雷达作为一种新型的探测装置,体积小、重量轻、抗干扰能力强、分辨率高,能够弥补其他探测手段的不足,尤其是在晴空条件下成为最有效地风场探测手段,对保障飞机起降安全、预警低空风切变有重要指导作用[8]。因此,利用激光测风雷达,准确地对低空风场预测逐渐成为保证飞行的重要工具。

自二十世纪70年代,著名气象学家FUJITA开启风切变的研究[9]。1983年,WOODFIELD和WOODS提出S因子算法[10],但这种仅通过求最大最小值之间风速差的算法极易产生较大误差。2011年,香港天文台陈柏玮等人提出香港机场沿用至今的单斜坡低空风切变自动预警算法[11],但是该算法受到斜坡长度的限制仅能探测600 m尺度的风切变。2012年,CHEN提出F因子算法[12],但直接计算由激光雷达得到的逆风梯度数据会导致结果快速波动且无法计算垂直分量的风切变。2014年,他们又提出利用涡流耗散率[13]预警低空风切变,但因为涡流耗散率的阈值在国际中没有明确的规定,该算法还需要在实际应用中继续探索。2016年,JIANG等人提出单双斜坡结合的算法[14]来预警风切变,不足的是该算法的目的是预警大尺度风切变。2019年,MA等人提出基于主成分分析(principal component analysis,PCA)和相位差校正的风切变预警算法[15],但是该算法的多普勒频移和风切变阈值还需要进一步通过经验确定。鉴于现有算法存在不足,本文中提出了一种基于斜坡改进的激光雷达风切变预警算法,将重组逆风廓线划分成不同的斜坡区间,再通过确定每个斜坡上的能量变化判断是否有风切变发生,最终确定风切变发生的位置和时刻,实现风切变的稳健预警。

2 原 理

2.1 激光雷达探测场景

激光雷达探测风切变主要采用下滑道扫描模式,该方式由香港天文台陈柏玮等人提出[10,16-18],利用同时改变激光雷达仰角和方位角的方式扫描飞机的进近走廊,直接探测低空风切变。传统的PPI扫描完成一个周期的时间约为6 min,而下滑道扫描一个周期得到径向数据仅需要2 min,可为飞机起降提供实时预警服务,保障飞行安全。在下滑道扫描过程中,激光光束与飞机起降跑道的夹角不小于30°。飞机起飞和降落时跑道与下滑道的夹角不同,起飞离场时夹角为6°,降落进场时夹角为3°。参考文献[19]提供了香港国际机场的地理位置。图1是激光雷达下滑道扫描示意图[10],其中图1(a)为降落场景,图1(b)为香港机场下滑道与激光雷达雷达对应位置关系的三维图。

图1 激光雷达下滑道扫描示意图

假设飞机的下滑道是被一条一定尺度的管道包裹,下滑道处在管道的中央,那么重组逆风廓线的径向数据就是落在管道中的数据[20]。为了有效提取管道内数据,以雷达为原点建立坐标系,利用笛卡尔坐标系将径向风速进行转换。以飞机降落为例,此时跑道与下滑道夹角为3°,将管道在笛卡尔坐标系三个平面投影,并根据激光波束与跑道夹角小于30°的约束条件,可得到x轴,y轴,z轴的坐标范围:

(1)

式中,d是激光雷达与跑道的垂直距离;L是激光雷达到跑道的垂足与着陆点之间的距离;α是激光光束与跑道方向的夹角;h为下滑道上径向数据点到y轴的垂足距离h=tan3°(y-l)。

2.2 算法原理

F因子算法是航空安全中使用较多的风切变预警方法,但其计算结果在某些情况下波动过于剧烈,长期存在虚警率过高问题;单斜坡算法可以实现风切变的快速检测,但受到斜坡长度的限制对持续时间较长的风切变检测结果不理想。基于现有算法存在的问题,本文提出一种将F因子与斜坡算法结合的改进斜坡算法。

该区域主要分布在流域的西部、西南部和东北部，海拔相对较高，虽然降水量较大，但降雨时空分布极为不均，且由于地势的原因，土壤蓄水能力差、地表水资源利用率低、农田水利设施落后，是流域内最容易受旱的地区，也是农业干旱脆弱性最严重的地区。经过对所搜集资料的综合分析，该区域可分为东北沂蒙山丘区、豫西山丘区和淮南山丘区等3个二级区。

图2为整个算法的流程图,具体步骤如下:

图2 斜坡改进算法的流程图

①确定重组逆风廓线并进行信号处理。首先筛选数据构造逆风廓线;再根据角度得到下滑道上的数据并进行质量控制得到重组逆风廓线,随后延长重组逆风廓线的尾部并进行一次滑动平均,最终得到没有数据缺失、尾部延长且保留特征的重组逆风廓线。

②根据风速差求每个斜坡上的能量变化。首先根据经典斜坡算法将重组逆风廓线划分成长度不同的斜坡,再利用每个斜坡上最大最小值的风速差确定斜坡的能量变化,并与阈值比较。

③发出预警或预警解除。主要是判断每个斜坡的能量是否超过阈值,超过阈值,发出风切变预警,此时该斜坡对应的位置即为风切变发生的位置;未超过阈值,则警报解除。如果计算结果远超阈值,还应将该斜坡提取出来,并进一步确定风切变的位置。

本文提出算法的第一个创新点是区别于传统方法通过计算相邻两点之间的能量变化进行预警,在得到每个斜坡上的速度变化之后,求得斜坡上的速度差,再计算每个斜坡上能量的变化情况。当计算结果超过阈值±0.05之后,发出风切变警告。此时,计算结果不再剧烈波动,虚警率也显著下降。其中,每个斜坡的能量变化率表示为:

(2)

式中,ΔVm表示斜坡上最大风速和最小风速之差;D为斜坡长度;g表示重力加速度,而Vapp表示飞机的进近速度。

经典F因子算法是求取重组逆风廓线上相邻两点之间的能量变化,并与阈值判断,从而给出风切变预警信息,如图3中Fc所示。而本文提到的斜坡改进算法中,将F因子与斜坡算法结合后,求取每个斜坡上能量的变化,然后与阈值判断再给出风切变预警信息,如图3中FI所示。

图3 斜坡改进算法与经典F因子算法区别

斜坡算法的本质是将重组逆风廓线划分成不同的速度变化区间,改进斜坡算法在得到重组逆风廓线上不同斜坡之后,再比较每个斜坡上的能量变化,避免了相邻两风速快速变化造成的F因子算法快速波动,降低了F因子算法的虚警率,也将风切变发生的位置直接定位到某个斜坡。当某个斜坡上速度发生突变,必然会导致该斜坡的能量发生突变,若超出阈值很多,此时将该斜坡提取出来,可以准确判断风切变发生的具体位置。

本文提出算法的第二个创新点是用最大最小值之差代替单斜坡算法中斜坡首尾的风速差。在单斜坡算法中,有两种情况需要对斜坡进行延长:一是斜坡持续增速,但跨过斜坡后降速;二是斜坡持续降速,但跨过斜坡后增速[11],这两种情况下,最长可延长斜坡至原始斜坡长度的1/2。第二种情况如图4中的局部重组逆风廓线所示。这两种情况在实际中很普遍,此时若还用单斜坡算法中斜坡首尾两端的风速差作为斜坡速度的变化,显然小于最大最小值的风速差,进而可能小于单斜坡算法中7.7 m/s的阈值,产生漏报。当斜坡需要收缩[11]的时候,最大最小值风速差与斜坡首尾风速差一致。利用最大最小值之差作为斜坡上速度的变化,当最大值出现在最小值之后,用最大值减最小值得到正的速度变化,表明该斜坡上风速增加,且风向是靠近激光雷达;当最大值出现在最小值之前,用最小值减最大值得到负的速度变化,表明该斜坡上风速减小,且风向是远离激光雷达。图4展现了原始斜坡长度取400 m,延长至600 m时,这一个斜坡上最大最小值风速差与单斜坡算法中首尾风速差的区别。

图4 斜坡改进算法与单斜坡算法区别

尾部延长的意义是防止最后一部分数据不满足最小原始斜坡长度,不能进行斜坡长度的确定造成信息缺失。尾部延长有两种形式,尾部对称性扩展法和尾部直接延长法[21]。尾部对称性扩展法是选取重组逆风廓线最后一部分数据做对称延长,一般选取的长度为5～10个数据点,当最后一段数据呈线性变化时,会人为造成最后一个点的突变;尾部直接延长法是将重组逆风廓线尾部最后一个元素直接扩展,扩展的长度至少为一个原始斜坡长度对应的数据点个数,该方法保证在对重组逆风廓线最后一段数据检测时,不会造成尾部缺失,故本文采用尾部直接延长法。

3 方法验证、对比与分析

2020年开始,香港国际机场应用激光雷达对飞机起降下滑道进行了全天候扫描,场景如图1所示。结合飞行员报告风切变发生的高度及风切变形成的原因,本文选取三组比较有代表性的数据对斜坡改进算法进行验证。由于香港国际机场产生风切变的主要原因多为海风锋,所以F因子阈值取±0.05[22-23],而非国际民航组织规定的±0.105,本文中对比的F因子、斜坡改进算法阈值亦取±0.05。

3.1 案例一:海风锋情况

2020年2月24日,北京时间6:15、6:28、6:31在07R跑道有三趟飞机在进场降落报告有风切变发生,报警高度分别为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;左边纵坐标表示重组逆风廓线对应的风速(实线部分),数值为正表示风向靠近激光雷达,数值为负表示风向远离激光雷达,单位为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的单斜坡算法,图中灰色阴影部分是风切变发生的范围,点线是信号处理后的下滑道重组逆风廓线。根据图7展示的结果,F因子波动较大,出现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检测到风切变次数依次是:6、2、1、0,由重组逆风廓线可知,在距离降落点2～4 km范围,风速变化虽然没有很大,但风向发生了改变;在距离降落点4～7 km范围,风速增加到最大,风向也由正向变为负向再变为正向;在7～8 km,风速从最大逐渐减小到0,再反向增加。按照判断标准,这三个阶段均是有风切变发生的,而原始斜坡取800 m、1600 m、3200 m的改进算法分别有不同层度的漏报,400 m的效果最好。

图9 斜坡改进算法检测结果

再进行不同算法直接的对比,图10给出了对比的结果,图10(a)是F因子算法结果,图10(b)是原始斜坡取40 0 m的改进算法结果,图10(c)是原始斜坡取400 m的单斜坡算法,图中灰色阴影部分是风切变发生的范围,绿线是信号处理后的下滑道重组逆风廓线。F因子从距离降落点0.3 km至距离降落点8.5 km均报告有风切变发生,对应的高度为15.72～445.47 m,检测结果起伏剧烈,波动较大,涉及范围较广,夸大了风切变的程度。斜坡改进算法给出了三个位置有风切的预警报告,对应的高度范围是:125.78～167.70 m、277.76 m、377.34～408.78 m。

图10 各种算法检测结果对比

相对于F因子的检测结果,斜坡改进算法给出风切变发生的范围更精确,呈现的结果也更清晰。单斜坡算法这次检测出一处发生了风切变,对应的横坐标为x=7.8 km,与斜坡改进算法中的最后一处一致。通过前面不同长度的改进算法对比分析,单斜坡算法产生了两处漏报。

3.3 案例三:无风切变情况

为了进一步对比F因子算法和改进斜坡算法的区别,我们找到一组几天内没有明显发生风切变、也没有飞行员报告的数据进行验证。

图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的标准,确定明显风切变发生。图12是取不同长度的斜坡改进算法对比结果图,图13是不同算法之间的对比结果图,其中图13(a)是F因子算法结果,图13(b)是原始斜坡取400 m的改进算法结果,图13(c)是原始斜坡取400 m的单斜坡算法,图中灰色阴影部分是风切变发生的范围,虚线是信号处理后的下滑道重组逆风廓线。

图11 2020年4月1日8:8:22时风场信息

图12 斜坡改进算法检测结果

如图13所示,在确定没有风切变发生的背景下,F因子检测出4个范围(6处)有风切变发生,虚警率过高,而不同长度的斜坡改进算法均没有产生报警,且检测结果远小于报警阈值。案例三反应了F因子的不足。

将单斜坡算法的结果与斜坡改进算法的结果放在同一张图上进行比较,如图14所示,其中图14(a)～(c)依次为案例一、案例二、案例三,图中左边纵坐标为斜坡改进算法,右边纵坐标表示单斜坡算法,图中所取的原始斜坡长度均为400 m。单斜坡算法和改进斜坡算法大部分走势和形状一致,但是两种算法也有细微区别,图中点线的圆圈部分展示了它们的差异。这是因为,在大部分时刻,改进斜坡算法中每个斜坡的最大最小风速与斜坡算法中首尾风速是一致的,但是在一些情况下,最大或最小风速就不出现在斜坡的首端或尾端,从而斜坡上风速差有细微区别,检测结果在局部发生改变,使得斜坡的能量变化超过阈值发出预警信息,这也使得改进斜坡算法更优于斜坡算法且弥补了斜坡算法产生漏报的不足。根据大量数据验算的结果,图14展示的有差异的情况很常见。

图14 斜坡改进算法与单斜坡算法区别

3.4 数据统计分析

根据香港国际机场飞行员返回的报告,2020年全年机场有风切变报告的天数为74天。论文将所有发生风切变对应时刻的激光雷达扫描数据提取出来,进行风切变检测,并对结果进行了统计,表1上半部分是有风切变发生时算法检测结果的统计。由表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 各种算法检测结果统计

再对没有风切变报告的292天进行算法检测。由于香港机场的激光测风雷达每3 min完成一次扫描,一小时能得到约20次下滑道数据,每日随机抽取连续2小时内共计40次下滑道数据,因此共得到11680次无风切变报告的下滑道数据进行统计,检测结果如表1下半部分所示。由表1可知,单斜坡算法和斜坡改进算法均没有检测到风切变,而F因子检测出5022次风切变,此时虚警率高达43 %!

结合有风切变和无风切变的情况统计,F因子虚警率高,单斜坡算法漏警率高,而斜坡改进算法原始斜坡取400 m时,检测率接近100 %,虚警率为0,有效弥补了传统风切变检测算法的不足。

斜坡改进算法将F因子和单斜坡算法结合起来,通过计算每个斜坡上飞机遭遇的能量变化来检测风切变。不同原始长度的斜坡对比结果表明:斜坡长度取400 m结果最好。这是因此,当斜坡长度较长时,若风速变化不太剧烈,则每个斜坡的能量变化较小;而原始斜坡长度取400 m时,一个斜坡上最多仅有6个相邻数据,保证了计算的数量和质量。斜坡改进算法与经典算法的对比结果表明:由于F因子是依次比较相邻两个数据点,使得结果波动剧烈,虚警过高;单斜坡算法是比较斜坡上首尾两端的速度差,当大于7.7 m/s时才会发出报警,加大了漏警的概率。而斜坡算法恰好弥补了这两种经典算法的不足,且给出的风切变发生范围更细致;利用斜坡上最大最小值作为风速差代替首尾两端的风速差,使得风切变的计算结果也更精确。因此在后面的飞行保障中,可以采用原始斜坡取400 m的斜坡改进算法对风切变进行预警。

4 结 论

本文利用激光测风雷达在香港国际机场开展了风切变观测试验,提出了基于斜坡改进的激光雷达风切变预警算法,将F因子与单斜坡算法结合,计算每个斜坡上的能量变化,弥补F因子虚警较高和单斜坡算法漏警较高的不足;用最大最小值代替斜坡首尾的风速差,提高了检测率和精度。通过对海风锋和阵风锋、对流、雷雨大风等情况下的风切变进行识别,并结合飞行员的报告对结果进行了分析,证实了斜坡改进算法的准确性和鲁棒性,有利于风切变检测,能为飞机起降提供技术保障。

基于激光雷达的斜坡改进算法,能有效弥补F因子虚警过高、单斜坡算法容易产生漏报等不足,正确预警飞机下滑道上是否有风切变发生及风切变发生的位置与范围,为飞行安全提供有效保障。但是该算法对于风速增加缓慢的数据检测不是很理想,后面将围绕此类数据对算法进行优化。