可业务化应用的L波段探空系统高空风改进算法

陈磊，卞建春，刘毅，石鹏飞，韩燕，唐敬

（1.中国科学院大气物理研究所中层大气和全球环境探测重点实验室，北京 100029；2.中国科学院大学，北京 100049；3.解放军63655部队，新疆 乌鲁木齐 841700）

可业务化应用的L波段探空系统高空风改进算法

陈磊1，2，3，卞建春1，刘毅1，石鹏飞3，韩燕3，唐敬3

（1.中国科学院大气物理研究所中层大气和全球环境探测重点实验室，北京 100029；2.中国科学院大学，北京 100049；3.解放军63655部队，新疆 乌鲁木齐 841700）

为进一步挖掘L波段高空气象探测业务系统原始资料潜力、优化改进现行测风业务算法，本文基于新疆博斯腾湖东岸L波段系统机动站累积的2 a探空原始资料，提出了可业务化运行的高空风改进算法。该算法首先对雷达原始秒点坐标进行严格的质量控制，采用低通滤波、加权最小二乘法、线性补缺等方法去除探空仪摆动、雷达测量误差等对秒点坐标造成的扰动；然后采用逐秒点坐标滑动计算矢量平均风，通过与同球施放的GPS探空做比对分析得出，在全程使用50～60 s计算时间窗口或前50 min使用30～40 s时间窗口、50 min以后使用50～60 s时间窗口条件下，雷达矢量平均风廓线与GPS矢量平均风廓线吻合较好；规定高度层风和固定垂直分辨率高度层风采用查找表方法确定，其结果不仅能在风场结构上与现行业务算法一致，同时能呈现出明显的风层中小尺度扰动信息。

L波段探空系统；测风算法；业务流程

高空风探测作为大气探测的重要内容，在气象学中具有非常特殊的意义。一是因为高空风所代表的空气运动标志着天气过程的发生或演变，在天气预报中具有重要作用[1]；二是高空风在北半球中高纬度地区平流层有显著的冬夏季节转换特征，高空风资料在平流层—对流层耦合机制研究和区域气候特征统计方面具有重要意义[2-3]；三是高空风的探测有别于温度、气压、湿度等气象要素的介质中直接测量，而是利用气流示踪物进行间接探测，其计算与分析处理方式有其自身的特点[4]。近年来，平流层飞行器逐渐成为人类拓展大气活动空间的重要方式，寻找适宜的飞行时段与区域是未来的研究热点与难点[5]，而高空风探测业务系统是目前可提供大范围地理区域下平流层的垂直分辨率数据的重要手段，在临近空间领域具有不可替代的重要作用；另一方面，在全球尺度的中高层大气动力学数值模拟中，需要对重力波效应进行参数化以获得合理的大气环流结构，而无线电探空观测由于其较高的垂直分辨率、较为全面的测量参数和良好的测量实施经济性，因此在获取重力波特性参数气候特征方面具有独一无二的作用[6-7]。

L波段高空气象探测系统是我国高空气象探测业务的主力设备，其数据采集率、测量准确度、操作自动化程度和稳定性相比701雷达有显著提高[8-9]，目前已推广到我国境内120个探空站。高探测量要素中，温度、气压、湿度的测量精度主要取决于传感器性能和环境状况，但由于风是根据气流示踪物坐标间接计算，因此其探测精度不仅取决于仪器跟踪性能，还取决于测风算法。对现行高探业务中的L波段二次测风雷达而言，测角测距精度、雷达工作状态、探空仪信号质量、目标物离测站距离、测风算法等都会对风场数据质量产生影响；而国际上先进的GPS探空系统采用卫星导航体制测风，对示踪物的定位精度较之雷达要更为准确且稳定性好，测风精度仅与测风算法有关[10]。以上因素导致L波段雷达测风质量与国际先进水平仍有一定差距，测风准确度指标尚难以满足WMO高空风观测精度要求[11-13]。在国产GPS探空难以大范围业务化实现、现行L波段雷达体制仍将延续多年的情况下，如何充分发挥现有探测仪器效能并优化改进全国测站数年累积资料，值得深入探讨和尝试。姚雯[14]等通过在L波段雷达探空仪中增加GPS模块的方式进行比对试验得出，L波段雷达的测风精度在对秒量级气球定位数据进行平滑处理的情况下可以接近GPS测风的水平。王缅[15]等通过滑动平均法处理秒量级的瞬时风，其计算的高空风与GPS-RS92测风结果一致性较好。梁建平[4]等提出了完全采用矢量平均风算法来求取规定高度风的业务流程方案，认为计算层时间间隔为30～50 s的滑动矢量平均风可以进一步提高与GPS测风结果的一致性。本文基于新疆博斯腾湖东岸车载L波段雷达获取的原始数据，在总结前人工作的基础上提出了可业务化实施的测风算法流程，为进一步开发利用高垂直分辨率原始资料、优化改进现行测风业务算法提供重要参考。

1 雷达测风算法简述及业务算法存在的问题

雷达测风采用“量得风层”算法，即求取示踪物在固定时间间隔内相对测站的位移水平投影矢量差，进而得到风向风速[15-17]。如图1所示，矢量差V→H的长度除以时间间隔即为风速，V→H同正北方向的夹角即为风向。“量得风层”算法有几个需要注意的地方：一是该算法本质上是一种矢量平均方法，其获取的风向风速反映的是计算层厚度内流场的平均状况。理论上计算层的时间间隔越小，“量得风层”算法计算出的风向风速应该越接近“瞬时风”状况。但由于气球—绳子—探空仪系统存在钟摆振荡以及大气湍流运动，因此基于球载探空仪跟踪系统采用“量得风层”算法计算的“瞬时风”并不能表示某时刻大气流场的真实运动，即利用雷达测量秒点坐标计算“秒点风”没有实际意义；二是不同探测仪器在同一“量得风层”下，计算的风向风速是有差别的，这主要取决于仪器对示踪物的定位精度、定位相应时间等，这也是目前GPS测风精度好于L波段雷达的主要原因。

现行L波段探空系统利用整分钟的气球坐标计算“量得风层”矢量平均风，即前20 min采用1 min的计算层，20～40 min采用2 min的计算层，40 min以上采用4 min计算层，并采用整分钟前后5点坐标来减小整分钟时刻的坐标误差，一定程度上提高了数据利用率和整分钟坐标精度[18]，但仍然存在以下问题：一是对雷达原始数据的利用率过低，如果雷达探测时间为80 min，即有4800条原始资料，实际只使用了约400条原始资料，利用率为8%；现行业务算法在低空和高空采用不同的“量得风层”时间间隔，仍然是基于701雷达设备精度和数据获取率的设计，直接导致了风廓线过于平滑。由于L波段雷达探测精度比701雷达有显著提高，因此完全可以摒弃现有业务算法对计算层时间间隔的限制；现行业务算法仅计算了整分钟时刻的矢量平均风，其余时刻风向风速通过内插获取，这种处理不仅导致内插前的原始风廓线时间分辨率过低，而且在特定高度层风的求取上引入了标量计算误差。虽然现行业务算法在获取基本风场结构特征上是足够的，但已经不满足精细化预报、中小尺度风场扰动等研究对实测探空资料的要求。

图1 “量得风层”矢量平均风算法示意图

2 L波段雷达原始数据特点与质量控制

L波段雷达和GPS探空的探空仪坐标数据获取率是一致的，已有分析指出二者在仰角不太低、离测站不太远的情况下定位精度相近，说明在一定测量时段L波段雷达测风可以接近或达到GPS测风水平[9]。因此非常有必要分析L波段雷达原始秒点坐标数据特点，在计算“量得风层”前对之做相应的质量控制。本文对2013年5月—2015年5月连续两年的雷达原始秒点数据筛查得出：在雷达跟踪状态良好、探空仪信号正常的情况下，单次可获取4500条以上的原始秒点（图2）；斜距、方位角和仰角曲线能直观体现该时次雷达跟踪性能，除方位角和仰角在近地面层可能呈现剧烈变化外，其余在整体上应呈现平缓变化；气压整体上应呈现平滑的指数递减变化规律，温度整体上应符合对流层和平流层典型温度变化特点，探空高度由气压通过压高公式反算得出；方位角和仰角曲线在100 s以下的尺度上呈现出明显的扰动，造成扰动的原因有探空仪摆动、雷达跟踪定位误差和大气湍流等，该扰动可以通过低通滤波进行减弱或去除,如图2局部放大图中的红线所示；气压、压高和温度秒点曲线在100 s以下的尺度上无明显扰动，因此无需对三者进行平滑处理。

图2 L波段雷达原始数据状况及局部细节

探空仪定位秒点坐标中，斜距出现的问题较多，其异常大致可以归为三类情况：一是单独出现的野点群和小段数据异常，可以通过基于加权最小二乘法和二阶多项式模型的局部回归方法去除，在局部回归过程中，对绝对偏差大于6倍平均绝对偏差的秒点斜距的权重设置为0，确保既能剔除野点群，又能保留原始有用信息（图3a、图3b）；二是大段平滑变化的数据异常，则必须对该异常数据段进行人工剔除，然后进行线性补缺（图3c）；三是斜距连续的剧烈扰动，说明相关时段雷达测距状况差或探空仪信号差，则该时次剔除或相关时段去除并进行线性补缺（图3d）。

图3 斜距异常状况分类及处理效果

3 矢量平均风计算与业务流程

本文采用逐秒点滑动计算矢量平均风的方法，即采用第1秒和第n秒的坐标数据计算第N2秒的风向风速，采用第2秒和第n+1的坐标数据计算第n/2+1秒的风向风速，以此类推，获取秒点矢量平均风。接下来的问题是如何选择矢量平均风计算时间窗口，需要兼顾2个方面：一是时间窗口要足够大，能够有效过滤探空仪摆动、雷达测角测距误差效应；二是时间窗口要足够小，计算出的矢量平均风能够保留风层的真实扰动信息。虽然雷达秒点坐标经过低通滤波已经去除了部分探空仪摆动效应和雷达测量误差，但二者的残留仍是造成“量得风层”计算偏差的一个重要因素。例如，探空仪通常垂直悬挂于气球下方30 m，在极端情况下，探空仪摆动可以在计算时间窗口内的探空仪水平位移方向上造成60 m左右的偏差，如果计算时间窗口小于30 s，则会造成2 m/s以上的风速计算误差。因此为进一步规避探空仪摆动效应和雷达测量误差，建议计算时间窗口应在30～60 s之间。为验证不同计算时间窗口下L波段雷达测风的准确性，本文采用2015年12月下旬24个时次同球施放的长峰GPS探空仪做比对试验。比对方法为：对GPS探空仪和L波段雷达数字探空仪获取的原始坐标秒点经过质量控制和时间匹配后，采用相同的时间窗口计算秒点矢量平均风。以GPS探空仪计算结果为真值，计算L波段探空仪的矢量平均风绝对误差，将所有比对时次的矢量平均风绝对误差按秒点进行分组，每组100个秒点，共计48组，并对每个分组做绝对误差分布箱型图，结果如图4所示。

从图4中可以看出，所有比对时次L波段雷达风速绝对误差在30～60 s的计算时间窗口条件下整体较小，且误差分布在探空前期和后期变化不大，50%的风速误差均在1 m/s以下，且随着计算时间窗口的加大，所有比对时次75%的风速误差均可控制在1 m/s以下；所有比对时次L波段雷达风向绝对误差在探空前期和后期变化较大，施放后50 min以内风向误差整体较小，75%的风向误差大都在5°以内，50 min以后风向误差明显变大，但随着计算时间窗口的加大，75%的风向误差大部分可以控制在5°以内。综上，在50～60 s计算时间窗口下，L波段雷达测风精度可以接近或是达到GPS测风精度，因此建议在业务中使用全程50～60 s时间窗口或是前50 min使用30～40 s时间窗口，50 min以后采用50～60 s时间窗口计算矢量平均风。

根据矢量平均风以及对应秒点的气压高度，采用查找表方法确定固定垂直分辨率高度层上的风向风速。为避免标量计算误差，将与固定垂直分辨率高度格点距离最近的秒点高度对应的风向风速作为高度格点上的风向风速，改进后的算法流程如图5所示。现行业务算法计算的固定垂直分辨率高度风和本文测风算法计算的固定垂直分辨率高度风对比如图6所示，从中可以看出，改进的算法可以呈现明显的中小尺度风速扰动情况，可以提供细节更为丰富的高空风场资料。

图4 L波段雷达风向（e，f，g，h）风速（a，b，c，d）绝对误差箱型图

图5 规定高度层风向风速计算流程

图6 改进算法（a）和现行业务算法（b）规定高度层风对比

4 结论

现行业务算法过低的探空仪坐标数据使用率直接导致了风廓线的过于平滑，并在后期固定分辨率高度风求取中引入标量计算误差，严重影响了现行L波段雷达探测效能的发挥和高空风数据的应用拓展。基于现行业务系统的探测能力，可以通过改进算法来挖掘现有装备探测潜力以期接近或达到GPS测风水平。研究结果表明：

（1）在L波段雷达和数字探空仪工作状况良好的情况下，获取的原始坐标秒点曲线整体平滑，由探空仪摆动、雷达测量误差、大气湍流等效应造成的扰动可以通过低通滤波进行有效去除；方位角和仰角秒点数据质量整体较好，斜距秒点问题较多，大部分的斜距问题可以通过加权最小二乘法或是人工剔除再线性补缺的方法进行质量控制。

（2）对质量控制后的秒点坐标采用逐点滑动计算矢量平均风的方法，获得时间分辨率为1 s的矢量平均风向风速。通过与同球施放的GPS探空做比对分析得出，在50～60 s计算时间窗口下，L波段雷达风廓线与GPS风廓线吻合较好，建议在业务中使用全程50～60 s时间窗口或前50 min使用30～40 s时间窗口，50 min以后采用50～60 s时间窗口计算矢量平均风。规定高度风和固定垂直分辨率高度风采用查找表方法确定，规避了标量运算。

（3）改进的算法在体现整层风场结构特征上与现行业务算法具有良好的一致性，同时可以呈现明显的风层中小尺度风层扰动，形成细节更为丰富的高空风场资料。

[1] 孙学金，王晓蕾，李浩，等.大气探测学[M] .北京：气象出版社，2010：172-173.

[2] 陶梦初，何金海，刘毅.平流层准零风层统计特征及准两年周期振荡对其影响分析[J] .气候与环境研究，2012，17（1）：92-102.

[3] 张爱英，任国玉，郭军，等.近30年我国高空风速变化趋势分析.高原气象，2009，28（3）：680-686.

[4] 梁建平，李宇中，黎洁波，等.L波段高空气象探测系统测风算法改进探讨[J] .气象科技，2014，42（5）：753-758.

[5] 张玉李.平流层冬—夏季节转换及多尺度动力过程的影响[D] .中国科学院研究生院（大气物理研究所），2015.

[6] 卞建春，陈洪滨，吕达仁.用垂直高分辨率探空资料分析北京上空下平流层重力波的统计特性[J] .中国科学D辑，地球科学，2004，34（8）：748-756.

[7] Wang L,Geller M A.Morphology of gravity-wave energy as observed from 4 years（1998-2001）of high vertical resolution U.S.radiosonde data[J] .J Geophys Res，2003，108（D16），4489，doi:10.1029/2002JD002786.

[8] 姚琪.新型L波段二次测风雷达系统测量精度分析[J] .现代雷达，2006，28（11）：22-24.

[9] 范跃，陈中钰，赵晓莉，等.L波段与59-701探空系统观测资料差异评估[J] .气象，2014，40（2）：238-246.

[10] 姚雯，马颖，黄炳勋，等.利用GPS定位资料分析L波段雷达测风性能[J] .应用气象学报，2009，20（2）：195-202.

[11] 李毅，邢毅，杨荣康.国产GPS探空系统探测能力分析[J] .气象科技，2012，40（4）：513-519.

[12] WMO.气象仪器和观测方法指南[M] .第6版.北京：气象出版社，2005：222-235.

[13] 王缅，李伟，陈永清，等.L波段探空系统高空风平滑计算方法探讨[J] .气象，2011，37（1）：85-91.

[14] 中国气象局监测网络司.常规高空气象探测规范[M] .北京：气象出版社，2002：10-12.

[15] 中国气象局监测网络司.L波段高空气象探测系统业务操作手册[M] .北京：气象出版社，2005：24-25.

[16] 张伟星，王晓蕾.WGS-84地心坐标系中高空风计算方法[J] .气象科学，2005，25（5）：484-489.

[17] 吴泓，李永，郑清华，盛丽.改进高空测风算法的试验[J] .气象科学，2011，31（1）：113-118.

An Improved Upper-air Wind Measurement Algorithm for L-band Radiosonde Sounding System

CHEN Lei1，2，3，BIAN Jianchun1，LIU Yi1，SHI Pengfei2，HAN Yan3，TANG Jing3
（1.Key Laboratory of Middle Atmosphere and Global Environment observation，Institute of Atmsopheric Physics，Chinese Academy of Sciences，Beijing 100029，China；2.University of Chinese Academy of Sciences，Beijing 100049，China；3.Unit of 63655 of the Chinese People’s Liberation Army，Urumqi 841700，China）

In order to further mining raw-data potential of L-band radiosonde sounding system and improve current upper-air wind measurement operational algorithm,based on the two years of sounding data obtained by the L-band radiosonde sounding system in Xinjiang Bosteng lake area, an improved upper-air wind measurement algorithm could be used in routine operation is presented. Firstly,itrun a strict quality control procedure on original second-ordered coordinate data,using low-pass filter,weighted least squares and linear compensation methods to remove coordinate data disturbance caused by radiosounder swinging,radar measurement error and atmospheric turbulence. Then,it adopt the vector average wind calculation method in certain time-window second by second to obtain second-ordered wind speed/direction.Through comparative analysis with same-balloon GPS sounding,the radar sounding wind-profile is in better agreement with GPS when using 50～60 s time-window whole process,or using the 30～40 s time-window before the 50min and 50～60 s time-window after the 50min;fixed vertical resolution altitude wind-profile is obtained from the lookout-table method using second-ordered wind-profile,the result is not only in agreement with current operational algorithm on wind-field structure,but also can provide obvious wind-field fluctuation information.

L-band sounding system；wind measurement algorithm；operational procedure

P412

B

1002-0799（2017）01-0022-06

10.12057/j.issn.1002-0799.2017.01.003

2016-06-16

国家863项目（2015AA7054034）。

陈磊（1984-），男，工程师，主要从事大气动力学方面的研究。E-mail：Legend31v7@163.com

陈磊，卞建春，刘毅，等.可业务化应用的L波段探空系统高空风改进算法[J] .沙漠与绿洲气象，2017，11（1）：22-27.