民机阵风载荷测量试飞颠簸潜势诊断技术研究

闫文辉，黄兴友，赵钰锦，张莹，杨涛

(1. 中国飞行试验研究院气象台，陕西 西安 710089；2. 南京信息工程大学气象灾害预报预警与评估协同创新中心，江苏 南京 210044)

1 引 言

阵风载荷指飞行器穿越阵风时，升力面产生的载荷增量，是大气扰动的客观表现[1]。气象学上把阵风载荷称为飞机颠簸，颠簸特别严重时会破坏飞机结构，形成安全隐患，甚至造成飞行事故。2017年 6月 18日，东航空客 330-200E 客机在俄罗斯境内遭遇极重度颠簸，导致28 人受伤。Sharman 等[2]估计美国民航飞机每年遭遇6 300 次中到重度颠簸，5 000 次重度到极重度颠簸。Williams等[3]估计每年飞机颠簸造成全美约2亿美元的经济损失。2019年中国飞行试验研究院试飞任务超5 000架次，试飞任务的不断增长致使飞行器遭遇飞机颠簸的风险也不断增加，严重影响试飞安全与效率。

飞机颠簸与大气湍流有着密切的联系，Kelvin-Helmholtz (K-H)不稳定是最为公认的大气湍流产生机制[4]。大气湍流区域的中值长度约60 km，中值厚度约1 km[5]。大气湍流区域存在着尺度不同的湍涡，飞机颠簸是由那些与飞机尺度相当的、混乱的湍涡造成的，这种湍涡的直径一般在15～150 m[6]。Meneguz 等[7]研究表明只有约 14%飞机颠簸是发生在对流天气附近，这意味着造成飞机颠簸的主要原因是切变重力波破碎和地形因素。Jaeger等[8]发现约66.7%的飞机颠簸发生在高空急流附近，因此高空急流是诱发飞机颠簸最主要的原因。

俞飞等[9]发现高空急流、高空槽、切变线以及高空脊是造成华北地区飞机颠簸的主要天气形势。这种基于天气形势分析的方法主观经验性太强，易造成很高的漏报率和空报率。因此，运用数值模式输出量计算飞机颠簸指数进行颠簸诊断和预报目前已经成为主流。长期研究建立了大量飞机颠簸单一指数，主要有Brown 指数、Dutton 指数、L-P指数、Ri 数以及 Ellrod 指数等[10-19]。徐佳男等[20]选取诊断效果较好的3种指数，统计了中国地区飞机颠簸的气候特征。但不同指数重点考虑的因素有所差异，因此在不同区域、不同季节及不同成因的飞机颠簸中，指数自身及指数之间的预报能力有所不同。Sharman 等[2]提出集成指数预报算法（Graphical Turbulence Guidance, GTG），将几种颠簸指数组合使用，其总体预报性能比使用单一指数诊断效果要好。Kim 等[21]将GTG 方法应用于东亚晴空颠簸的研究，并估计了当前晴空颠簸数值预报水平。沈强等[22]选取了9个预报效果较好的颠簸指数，根据预报得分计算权重，构建了飞机颠簸综合指数。

为拓展民机遭遇飞机颠簸时的飞机包线、提升飞机品质、完善民航试航取证大纲，中国飞行试验研究院以ARJ21 支线客机作为试验机，承接民机阵风载荷测量试飞型号研发项目。阵风载荷测量试飞在国内是首次进行的，具有技术新、找风难、风险高及工作量大等特点。为准确预报潜在飞机颠簸区域，解决“找风难”问题，同时避免其他航空危险天气对试飞安全的影响，本文选取对流层高层（7 500 m 以上）诱发飞机颠簸的常见天气系统（高空急流），利用欧洲中期天气预报中心（European Centre for Medium - Range Weather Forecast, ECMWF）提供的 ERA5 再分析资料，从飞机颠簸飞行员天气报告（Pilot Reports, PIREPs）入手，基于10个应用成熟的飞机颠簸指数，通过加权集成给出了颠簸潜势综合指数，同时对颠簸潜势综合指数的性能进行了评估。

2 资料介绍

2.1 ERA再分析资料

ERA5 再分析资料是 ECMWF 利用 4DVar 同化产生的第5代大气再分析全球气候数据，该数据提供1979年至今的逐小时大气、陆地和海洋的气候变量估计值，其中大气数据的空间分辩率精确到了0.25 °×0.25 °网格，等压面数据垂直方向包括1～1 000 hPa的 37 层等压面[23]。数据包主要包括位势高度、温度、经向风、纬向风、垂直风、相对湿度、比湿、散度及涡度等参量的网格数据。

大气再分析数据在飞机颠簸指数的诊断和预报中具体广泛的应用，为了进一步说明ERA5再分析资料的准确性，本文随机下载了2020年9月21日08 时（北京时间，下同）300 hPa 天气图(图1)，图1a 为 ERA5 再分析资料，图 1b 为 ECMWF 细网格预报资料中08 时的天气实况。对比可见，二者基本一致，欧亚中高纬度均成两槽一脊型，中纬度为纬向平直气流，中低纬度受反气旋控制，并且风场、等温线及等高线的分布、大小基本一致。

图1 2020年9月21日08时300 hPa天气图

2.2 飞行员天气报告（PIREPs）

在飞行过程中若遇到飞机颠簸、飞机积冰等现象时，机组人员会以话音通信方式将现象、时间、位置、高度及强度等信息传输给管制单位，管制人员再将该信息发送给气象人员，气象人员将其记录并通过专业技术分析增添造成原因、备注等信息，最后进行统计保存，这种资料称为飞行员天气报告（PIREPs）[24]。为了准确评估各单一颠簸指数的预报性能，并推算颠簸潜势综合指数，本文筛选了 2017年 10月—2020年 12月期间，中国范围内7 500 m 高度以上高空急流背景下的170 例飞机颠簸报告及170 例无颠簸样本，其中280 例（140 例发生颠簸、140 无颠簸）样本用于构建颠簸潜势综合指数，剩余的60 例（30 例发生颠簸、30 无颠簸）样本用于评估颠簸潜势综合指数。具体报告位置分布如图2所示，颠簸产本基本分布于我国中纬度地区，与温带急流和副热带急流位置对应较为良好。需要注意的是，无颠簸样本选取阎良机场附近9 000 m左右的高度上一个固定位置，试飞期间当飞机未发生颠簸时取样一次，所以只有一个点。

图2 飞机颠簸报告位置分布图

3 单一颠簸指数介绍

飞机颠簸和大气湍流活动强弱有密切的关系，两者皆发生在风场垂直切变区、水平切变区、辐合或辐散区、水平形变区以及强的水平温度梯度区[10]。由于不同指数侧重因素不同，其预报效果也不同，综合考虑这五个方面，参考前人研究[2,21,26]，选取反映上述特性的，预报效果较为理想的10个飞机颠簸单一指数来构建综合指数。

(1)L-P指数。

前苏联水文气象中心综合考虑了温度切变、水平风场的水平切变和垂直切变情况，构建了L-P指数来预报潜在飞机颠簸区[10]，L-P指数采用如下计算方式：

式中，SV为水平风场的垂直切变(单位：m/(s·100 m))；ST为水平温度梯度(单位：K/(100 km))；SH为水平风场的水平切变(单位：m/(s·100 km))；u和v为水平风分量(单位：m/s)；T为气温（单位：K）；z为垂直高度(单位：位势米)。当L＜0 时，预报没有颠簸，反之在此基础上根据L指数计算概率P：

P的大小反映了飞机颠簸强度。

(2) -Ri数。

Richardson 数（Ri）被定义为表征大气稳定度的 Brunt-Vӓisӓlӓ 频率的平方和表征破坏大气稳定度的垂直风切变的平方之比[11]，其表达式为：

式中，θv为虚位温(单位：K)；SV为水平风场的垂直切变(单位：s-1)；g为重力加速度(单位：m/s2)，由于Ri数越小，颠簸强度越大，为了方便计算，取-Ri数作为本文颠簸指数。对于低静力稳定度和高垂直风切变区域，例如高空急流或锋区，如果晴空湍流将发展，通常-Ri大于0.25。

(3) Ellrod指数。

Ellrod 等[12]对锋生函数进行简化，综合考虑了水平风场的垂直切变、总形变以及辐合辐散，提出了两种飞机颠簸指数TI1和TI2，其表达式如下：

式中，SV为水平风场的垂直切变(单位：s-1)；DEF为水平风场的总形变(单位：s-1)，包括拉伸形变和切变形变两项；Div为水平风场的散度(单位：s-1)。在业务应用过程中，王洪芳等[25]发现散度项的指示性不强，所以本文采用TI1指数，TI1高值区对应强的飞机颠簸区。

(4) MOSCAT指数。

MOSCAT 指数是由NECP 嵌套网格模式NGM（Nested Grid Model）基 于 MOS（Model Output Statistics）方法统计输出的晴空湍流指数，又称概率预报因子指数[13]，其表达式为：

式中，|V|为水平合成风速(单位：m/s)，DEF为水平风场的总形变(单位：s-1)。

(5) Brown指数。

Brown[14]利用热成风关系简化了Ri 倾向方程，提出了Brown指数：

式中，ζa=ζ+f为绝对涡度(单位：s-1)，ζ为相对涡度，f=2Ωsinφ为科氏参数。DSH为切变形变项(单位：s-1)，DST为拉伸形变项(单位：s-1)，SV为水平风场的垂直切变(单位：s-1)。由于大气运动为准水平，所以将涡度简化为垂直方向的分量。

(6) Dutton指数。

Dutton 指数[15]综合考虑水平风场的水平切变和垂直切变对飞机颠簸强度的贡献，得到一个经验指数，其表达式为：

式中，SH为水平风场的水平切变(单位：m/(s·100 km))；SV为水平风场的垂直切变(单位：m/(s·km))；10.5为经验常数。

(7) 曲率方法指数。

综合分析飞机颠簸报告和环流形势，发现高空槽脊的曲率与飞机颠簸强度相关，而曲率与涡度有关[16]，定义曲率指数为：

式中，ζ2单位为 s-2，u和v为水平风分量(单位：m/s)，此处涡度也取垂直方向的分量。

(8) 水平温度梯度指数。

水平温度梯度指数是一个从热成风关系中计算形变和垂直风切变的方法[17]，其表达式为：

式中，HTG单位为K/m，T为气温 (单位：K)。

(9) 水平散度指数。

水平风场的强辐合辐散区与风场的切变区域对应良好，并且辐合可以对锋生产生作用[2]。定义水平散度指数为：

式中，DIV单位为s-1，u和v为水平风分量。

(10) 风相关指数。

急流背景下的飞机颠簸与大的风速相关[11]，定义风相关指数如下：

式中，S 单位为 m/s，u和v为水平风分量(单位：m/s)。

4 颠簸潜势综合指数算法

由于不同指数重点考虑的因素有所差异，因而指数自身及指数之间的预报能力也有所差异，为确保颠簸指数能够有效识别多种复杂条件下的飞机颠簸，本文参考美国国家大气研究中心（National Center for Atmospheric Research,NCAR）Sharman 等[2]提供的 GTG 算法，给出了我国高空急流背景下的颠簸潜势综合指数。

4.1 数据归一化

由于每个颠簸指数的单位和数值不同，计算结果的量级也相差较大，难以进行比较和加权计算，因此需要进行归一化处理。本文根据各单一颠簸指数的阈值，基于分段线性函数将计算结果进行归一化处理，将飞机颠簸强度分为无、轻度、中度、重度和极重度5个量级，中值分别对应数值0.0、0.25、0.5、0.75 和1.0。单一颠簸指数的具体阈值选取参考了 Sharman 等[2]、Kim 等[21]及 Willims等[26]的工作，其分类方法和原理为根据大量飞机颠簸报告强度等级与指数计算值之间的对应关系，统计各颠簸指数计算值的中值为该指数在当前量级的强度阈值，具体阈值如表1所示。

表1 高空急流背景下10种颠簸指数强度阈值划分

4.2 确定颠簸指数权重

在完成各颠簸指数归一化处理之后，为了定量评估各指数的预报水准，引入基于检测概率的评估方法POD(Probability of Detection)[27]，包括发生颠簸的预报准确率PODY、未发生颠簸的预报准确率PODN 及总体预报准确率PODA，其表达式如下：

式中，YY、YN 分别表示发生了颠簸，对应的指数预报正确和错误，NN、NY 分别表示未发生颠簸，对应的指数预报正确和错误。PODY 值越大，漏报率越低；PODN值越大，空报率越低；PODA值越大，总体预报效果越好。

为了更加精准地评估各指数的预报水准，再引入了 TSS（True Skill Score）评分[28]来评估总体预报准确性，其表达式如下：

TSS 取值越大，表明该指数的诊断效果越好。根据预报准确率计算各颠簸指数预报得分φ[2]，计算公式如下：

式中，C 和p 为常数，用于调节fMOG的重要性，通常取C = 1，p = 0.25。fMOG为在数值计算过程中得到的颠簸强度为中度及以上的格点数与总网格点数的比值。

根据得到的预报得分，可计算出每个颠簸指数的权重Wn[2]，具体表达式如下：

4.3 综合颠簸指数

CAT(i, j, k)的取值范围为[0, 1]，阈值为 0、0.25、0.5、0.75 及1，分别对应颠簸强度为无、轻度、中度、重度和极重度5个量级的中值。

5 颠簸潜势综合指数性能评估

选取与飞机报告时间、高度最临近的ERA5再分析资料作为输入量，计算各颠簸报告点的单一指数归一化值，然后利用PODY、PODN、PODA、TSS 及fMOG对各单一指数进行综合评估，并计算预报得分φ，最终结果如表2 所示。可见，在单一指数中，L-P指数发生颠簸时的预报准确率PODY 最高，达 97.9%，漏报率最低，TI 指数和 Dutton 指数次之，PODY 分别达到了67.9%和67.1%。Dutton指数未发生颠簸时的预报准确率PODN 最高，达82.1%，空报率最低，Brown 指数和TI1 指数次之，PODN 分别达到了77.1%和69.3%。从总体预报准确率PODA 和TSS诊断评分来看，L-P指数因具有最高的PODA(76.1%)和TSS(0.521)，总体预报准确率最好，Dutton 指数及Brown 指数的总体预报效果次之，PODA分别为74.6%及71.8%，TSS评分分别为0.493 及0.436。同时，从中度及以上颠簸所占面积fMOG来看，Dutton 指数、Brown 指数及TI1 指数的fMOG分别为7.8%、7.9%及12.4%，接近徐佳男等[20]利用AMADR 资料统计出的中度及以上颠簸仅有 10%左右的结果，较为合理。L-P指数识别到的中度及以上颠簸所占面积fMOG为20.4%，远大于10%，识别面积过大。综合考虑预报准确率和中度以上颠簸所占面积，Dutton 指数具有最高的φ(1.043)，总体预报效果最好，Brown指数、L-P指数及TI1 指数次之。另外，水平散度指数 DIV、曲率方法指数ζ2及-Ri 指数是这 10个单一指数中预报效果最差的3个。

表2 280例建模样本的单一颠簸指数性能评估表

在高空急流背景下，Dutton指数、Brown指数、L-P指数及TI1 指数预报颠簸的效果好，分析物理机理可见，水平风场的水平切变、垂直切变以及水平形变的贡献最大，强的水平温度梯度和较大的风速也有一定的影响。而水平散度指数DIV、曲率方法指数ζ2及-Ri指数诊断错误的原因可能是，对于对流层高层大气，大气环流以纬向型为主，辐合辐散、气旋反气旋特征不明显，且大气稳定度较高。

然后根据第3 节中介绍的流程对各单一指数进行加权集成，形成最终的颠簸潜势综合指数，利用剩余的60例颠簸样本对颠簸潜势综合指数进行验证。为了更加深入地评价颠簸潜势综合指数的优越性，本文额外引入了模型性能评价体系中广泛 使 用 的 ROC 曲 线（Receiver Operating Characteristic Curve）及 AUC 值（Area Under Curve）[29]。ROC 曲线，它是以 1-PODN 为横轴，PODY 为纵轴绘制的曲线，它越靠近左上角，表明指数模型区分发生颠簸和未发生颠簸的性能越好。AUC 就是ROC 曲线和x轴（1-PODN 轴）之间的面积，反映了模型把发生颠簸的样本排在未发生颠簸样本前面的比例（如果AUC=1，说明模型100%的将所有发生颠簸的样本排在未发生颠簸样本前面）。

60 例验证样本的颠簸指数性能评价结果和ROC 曲线如表 3 及图 3 所示，图 3a 为 10个单一指数的 ROC 曲线，图 3b 为 9个综合指数的 ROC 曲线，其中CAT2 到CAT10 分别表示根据预报得分φ的高低，由2个到10个单一指数组成的综合指数。可见，10个单一指数的预报效果均好于随机预报模型，且综合指数较好地集成了各个单一指数的优点，预报水平始终优于各单一指数。CAT5 为这9个颠簸潜势综合指数中最好的综合指数，由预报效果最好的五个单一指数（Dutton 指数、Brown 指数、L-P指数、TI1 指数及风相关指数S）加权集成，其 权 重Wn分 别 为 0.237、0.220、0.205、0.186 及0.152，预报得分φ达 1.225，且 PODY 达 86.7%，PODN 达 93.3%，PODA 达 90%，TSS 达 0.80，中度及以上颠簸所占面积fMOG为9.2%，各项指标均较为优秀，总体效果最好。同时，发现60例验证样本统计得到的单一指数性能与280 例建模样本统计的单一指数性能略有差异，说明样本量的大小和样本间的差异性会影响综合指数建模的好坏，但总体上Dutton指数、Brown指数及L-P指数效果最好，水平散度指数DIV及曲率方法指数ζ2效果最差。

图3 60例验证样本的ROC曲线

表3 60例验证样本的单一颠簸指数及最优综合指数性能评估表

根据预报得分φ的大小逐步增加单一指数数量，验证颠簸潜势综合指数对单一指数数量的敏感性，结果如图4a所示。可见，当单一指数增加到5个的时候，综合指数效果最优，预报得分φ达1.225，之后随着单一指数数量的增加，对综合指数预报效果开始产生负作用。为了进一步验证颠簸潜势综合指数对参与建模样本数量的敏感性，将参与建模的样本从40个逐步增加到280个，计算得到的颠簸潜势综合指数最优组合的预报得分φ如图4b 所示。可见，当建模样本数量为40个时，综合指数的预报得分φ最大，但随机误差较大，综合指数预报准确性置信度较低。综合指数的预报得分φ总体波动性不强，但是也随样本量的增加，综合指数性能逐步提升。

图4 颠簸潜势综合指数预报得分φ对单一指数数量敏感性图(a)及对参与建模样本量敏感性图(b)

本文选取一例样本来展示诸单一指数和综合指数的预报差异。2020年1月6日12:52，在蒲城及以东地区8 500 m 高度上，试飞院参试飞机（C919）出现中到重度颠簸，持续20 分钟以上。选取与颠簸报告出现时间、高度最临近的ERA5再分析资料（2020年 1月 6日 13:00，300 hPa、350 hPa及400 hPa三层等压面）作为数据源，350 hPa高空天气图如图5a 所示，陕西中南部及以东地区出现高空急流，急流轴风速大于40 m/s，颠簸报告点位于急流轴北侧风速切变最明显和温度梯度最大的区域。10个单一颠簸指数和颠簸潜势综合指数预报结果如图5b～5l所示，图中填色表示颠簸的区域和强度。可见，除了风相关指数S 外，其余指数位置都诊断正确，飞机颠簸报告点落在预报区域内。同时，各指数诊断的颠簸区域面积、强度分布均有较大差异。

图5 2020年1月6日12:52高空急流背景下天气形势图及各颠簸指数诊断结果

6 结论与讨论

本文利用飞机颠簸报告资料和ERA5 再分析资料，基于10个应用成熟的单一颠簸指数，通过加权集成构建了我国高空急流背景下的颠簸潜势综合指数。

(1) Dutton 指数因具有最高的预报得分φ(1.043)，是这10个单一指数中总体预报效果最好的指数，Brown指数、L-P指数及TI1指数次之。另外，水平散度指数DIV、曲率方法指数ζ2及-Ri 指数是这10个单一指数中预报效果最差的3个指数。

(2) 随着单一指数数量的增加，综合指数的预报效果先增大后减小，当单一指数增加到5个的时候，即由Dutton指数、Brown指数、L-P指数、TI1指数及风相关指数S 加权集成的时候，综合指数总体效果最优，PODY 达 86.7%，PODN 达 93.3%，PODA 达 90%，TSS 达 0.80，中度及以上颠簸所占面积fMOG为9.2%，预报得分φ达1.225，漏报率和空报率都得到了很好的改善，各项指标均较为优秀，总体效果最好。

(3) 随着参与建模的飞机颠簸样本量增加，综合指数性能逐步提升，随机误差逐步减小，置信度逐步提升。

该颠簸潜势综合指数可用于民机阵风载荷测量试飞解决“找风难”问题。飞机颠簸潜势综合指数的准确性取决于ERA5再分析资料的精度，同时该综合指数的评估效果取决于PIREPs 资料的数量，数量越多，评估效果越好。PIREPs资料在反映实际飞机颠簸强度的时候，仅仅依靠机组人员的主观意识，客观性弱且机型差异等因素影响较大，这给本文算法阈值的选择带来一定干扰，还值得进一步改进和优化。