2016, 27 (3): 257-272
2. 南京大学大气科学学院，南京 210046;
3. 美国国家大气研究中心，科罗拉多州博尔德 CO 80307;
4. 中国气象科学研究院，北京 100081;
5. 北京市气象局，北京 100089
2. School of Atmospheric Sciences, Nanjing University, Nanjing 210046;
3. National Center for Atmospheric Research, Boulder, CO 80307;
4. Chinese Academy of Meteorological Sciences, Beijing 100081;
5. Beijing Meteorological Service, Beijing 100089
目前，在强对流天气临近预报预警业务中应用的客观方法和自动化系统，主要以天气雷达资料为基础，通过计算机自动识别对流风暴 (包括可能产生的强天气类型) 并进行跟踪和外推预报[1-3]。虽然部分系统也利用了数值模式结果，但并没有过多考虑中低层精细的大气动力、热力特征对风暴生消发展的影响，因此，很难对快速发展和局地新生的对流风暴做出客观准确的临近预报和预警[4]。
近30年来，大量的中尺度观测分析和对流尺度数值模拟研究已证实，对流风暴的新生和发展演变通常与低层三维动力和热力条件及其相互作用存在密切联系。一方面，低层三维风场指示了近风暴前环境的辐合上升、垂直风切变、风暴相关螺旋度以及风暴尾部入流等低层动力特征，一定程度上也决定了风暴及其所产生的强天气类型[5-6]。另一方面，近地面冷池是大多数对流风暴的显著特征，低层温度场直接反映了风暴下部近地面冷池的强弱和深度，从而决定了风暴出流边界 (阵风锋) 的发展演变[7-8]。而二者相互作用，成为判断特定类型对流风暴发生发展、传播移动和生命史特征的重要判据之一[9-11]。研究也表明，边界层辐合线 (包括地面冷锋、阵风锋、海风锋、干线等) 与对流风暴的新生、演变及生命史具有密切关系[12-14]。因此，如何获得低层大气精细的动力和热力特征，成为提升强对流天气特别是快速发展和局地新生对流风暴临近预报水平的关键。
基于上述问题，近20年来，大量学者开展了雷达资料四维变分同化技术 (4DVar) 研究，致力于通过4DVar的先进同化技术，将高时空分辨率雷达观测信息引入数值模式，实现对精细的大气动力和热力特征的有效模拟，从而提升强对流天气的短时临近预报水平[15]。Sun等[16-17]基于雷达资料4DVar技术和一个简化的三维数值云模式及其伴随模式，发展了一套变分多普勒天气雷达分析系统 (VDRAS)，能够反演出对流尺度的三维动力和热力场，并得到示范应用[18-19]。随后，一系列研究和测试对VDRAS进行改进[20-23]。2008年同化北京和天津两部多普勒天气雷达资料的VDRAS系统在北京奥运会期间的强天气临近预报中得到成功的实时示范应用[23]。2012年以来，在VDRAS基础上通过拓展研究，逐渐形成了一个基于雷达资料快速更新四维变分同化 (RR4DVar) 技术和三维数值云模式的对流尺度数值模拟系统——快速更新雷达四维变分分析系统，并对不同类型对流风暴个例的动力和热力特征进行模拟分析[11, 24-28]。RR4DVar技术的优势就在于通过简化的三维数值云模式和它的伴随模式，可在短时间内 (10~20 min) 实现对逐6 min更新的多部雷达观测资料 (包括径向速度和反射率因子) 的快速更新4DVar同化和模拟分析。
但对于快速更新雷达四维变分分析系统，还需要开展一系列敏感性试验研究和评估，进一步确定系统各模块特别是雷达资料RR4DVar的合理配置方案，为后续改进及其在强对流天气临近预报预警中的应用奠定基础，也为目前国内正在开展的雷达资料同化和对流尺度数值预报研究提供借鉴。
1 系统介绍 1.1 资料及其预处理雷达资料预处理和质量控制是实现雷达资料有效同化并改进对流尺度数值模拟和预报的关键和前提[29]。目前，系统使用的雷达资料包括京津冀地区4部S波段 (北京、天津、石家庄、秦皇岛) 和2部C波段 (张北、承德) 多普勒天气雷达基数据 (径向速度和反射率因子)。对雷达资料的预处理和质量控制分3步。第1步是基于北京自动临近预报系统 (BJ-ANC)[13]的雷达资料一级预处理和质量控制，这里不再赘述。第2步是在变分分析系统中的雷达资料二级预处理和质量控制，主要包括如下几方面：①径向速度退模糊。S波段雷达和C波段雷达在业务扫描模式 (VCP21) 下最低几个仰角的最大不模糊速度 (Nyquist速度) 大约分别为±27.1 m·s-1和±13.55 m·s-1，因此，在第1步质量控制基础上，基于Lim等[30]的方法，根据不同雷达每个仰角的Nyquist速度表，寻找可能的速度模糊区域并进行退模糊处理。②阈值界定。根据地物回波和超折射地物回波的径向速度接近于0以及反射率因子廓线随高度快速递减的特征，设定在雷达最低两个仰角的径向速度绝对值小于0.25 m·s-1及反射率因子垂直递减率大于25 dBZ的资料为无效值；另外，认为小于12 dBZ的回波为非降水回波，不被同化。③资料平滑和填充。对雷达资料进行插值、平滑、填充等一系列处理。第3步是雷达资料的观测误差统计。对于本系统，基于局地标准偏差计算方案，在资料插值过程中，通过对每个格点与九点平均值之间的标准偏差，估计雷达观测误差。
除雷达资料外，也使用京津冀5 min自动气象站 (AWS) 观测及北京快速更新循环数值预报系统 (BJ-RUC) 3 km分辨率结果，用于计算4DVar的中尺度背景场、初猜场和边界条件。
1.2 云模式及RR4DVar 1.2.1 云模式及其微物理方案简化本系统使用了一个简化的三维数值云模式表征对流尺度大气运动，方程采用滞弹性近似，包含3个动量方程及热力方程、雨水方程和总水方程。模式预报变量有三维风场、液态水位温、雨水混合比和总水混合比，其中液态水位温作为热力变量，而水汽混合比、云水混合比、扰动温度和扰动气压则由预报变量诊断得到。模式采用Kessler微物理方案，包含了云水凝结和蒸发、雨水蒸发、云雨转换、云增长、雨滴沉降等物理过程的参数化[16]。
为适应聚焦低层的快速实时运算，本文对微物理方案进行了简化，重点考虑与低层热力变化特别是冷池及出流边界 (阵风锋) 发生发展、辐合上升等密切相关的雨水蒸发冷却和雨滴沉降过程。热力变量控制方程被简化为
|
(1) |
式 (1) 中，θ是扰动位温，在不同模式层由每一层的总位温减去平均值得到，κ是热力学扩散系数，Lv=2.5×106 J/kg是蒸发潜热，
|
(2) |
在式 (2) 中，β=0.0486 s-1，qv和qr分别是水汽混合比和雨水混合比，qvs是饱和水汽混合比，可由模式的温度场诊断得到。
1.2.2 RR4DVar设计简言之，4DVar同化的目的是通过求解模式方程，找到与实际大气观测结果在时空分辨率上尽可能接近的模式变量场，一般通过对一个预先定义的代价函数进行递归最小化得到。这个代价函数则描述了模式变量场与观测和背景场之间的偏差程度。本系统中雷达资料4DVar同化的代价函数定义如下：
|
(3) |
式 (3) 中，Jo为观测偏差项，代表了模式变量与雷达观测之间的偏差。根据Sun等研究[16, 18]，在假定误差协方差矩阵是对角阵即观测误差不相关时，Jo被写为
|
(4) |
式 (4) 中，vio和Zio代表了第i个雷达观测的径向速度和反射率因子，vi和Zi代表了与它们相对应的云模式变量，F表示从模式格点到PPI雷达资料点的变量转换函数，ηv和ηz代表径向速度和反射率因子的权重系数，σ和τ则分别代表同化窗的空间域和时间域。
式 (3) 中的Jb为背景场偏差项。在雷达资料4DVar中背景场非常重要，主要用来提供模式侧边界条件、对雷达观测空缺区域进行填充，并为式 (3) 中代价函数J的最小化提供一个初猜场，因此, 有利于快速降低模式噪声并提高分析精度。Jb可以表示为
|
(5) |
式 (5) 为矢量关系式，其中, x0表示同化开始前的模式变量场，xb表示背景场，B表示背景场误差协方差矩阵，所以式 (5) 代表了模式变量场与背景场之间的偏差。在早期的VDRAS中，B被简化为一个对角阵[23]。目前，一个多通道递归滤波客观分析方案[31]被引入，并在满足高斯分布的各向同性相关函数约束下，用来估计B矩阵。在这个滤波方案中，设定不同的预报背景场误差水平不相关长度，其中微物理变量为8 km，速度场和温度场为32 km。
目前，系统采用一个冷启动和多个热启动的快速更新循环方式运行，因此，对于冷启动和热启动的背景场计算方案有所差异。在冷启动时，背景场为
|
(6) |
而在热启动时，背景场为
|
(7) |
式 (6) 和式 (7) 中，xmb为中尺度分析背景场 (初猜场)，通过基于Barnes插值方案和最小二乘拟合的三维客观分析方案，计算流程如下：①基于每部雷达2.4°仰角观测 (假定雷达扫描模式为VCP21)，得到雷达VAD风廓线；②计算25 km间隔BJ-RUC模式探空廓线；③基于Barnes插值方案将雷达VAD风廓线和BJ-RUC模式探空进行融合，并插值到本系统模式水平网格点上，形成地面以上中尺度分析场；④基于Barnes插值方案将AWS观测插值到模式水平网格点上，形成地面中尺度分析场；⑤基于局地线性最小二乘拟合法，将地面和地面以上中尺度分析场进行融合，形成最终的三维中尺度分析背景场 (初猜场)xmb。
式 (7) 中，xfb表示上一循环的云模式预报背景场，λ为权重系数 (范围从0到1)。可见，热启动背景场实际由上一循环预报背景场xfb与中尺度分析背景场xmb加权得到。在以前的VDRAS中，因考虑到中尺度观测资料稀疏，在热启动时全部使用xfb作为背景场。目前，为确保大范围高分辨率模拟结果 (特别是低层) 不至于和中尺度分析背景偏差太大，并在地面AWS观测资料较为丰富、BJ-RUC模式探空较为可靠情况下[32-33]，在热启动时也对中尺度分析背景场xmb给予一定权重。
式 (3) 中，Jp为约束项，主要用来对模式风场 (u, v, w) 和扰动位温θ从时间和空间上实施强迫约束和平滑，以确保在没有雷达资料的区域模式解不会偏离观测和中尺度背景场太远，从而确保代价函数能够加速收敛[18]。Jp项中的一组权重系数一般通过测试和误差计算得到，但存在很大主观性。在目前的系统设计中，仅在冷启动时需要Jp项，而在随后的热启动循环中则使用上述递归滤波方案和背景场计算方案即可满足中尺度约束和代价函数收敛条件，并确保快速更新循环稳定运行。
至此，主要问题归结为在云模式方程约束下求解式 (3) 中代价函数J的最小值问题。系统使用了一个有限元准牛顿迭代算法用来求解代价函数的最小值。该算法需要计算代价函数及其梯度。代价函数通过对云模式正向积分得到，代价函数梯度通过对云模式的伴随模式反向积分得到。
1.3 系统运行流程简言之，雷达四维变分分析系统是在VDRAS基础上发展而来，核心是雷达资料RR4DVar同化、一个简化三维数值云模式及其伴随模式。系统以快速更新循环方式运行，流程如图 1所示，不同长度同化窗主要用于敏感性试验。
|
|
| 图 1. 系统的快速更新循环流程 Fig 1. Rapid-refresh cycling processes of the system | |
2 敏感性试验 2.1 个例说明
敏感性试验选取2009年7月22日发生在京津冀地区的一次在弱天气尺度系统强迫下的强风暴过程。2009年7月20日—8月20日北京城市气象研究所开展了2009年环北京夏季暴雨观测试验，选取该个例的优势是可以利用这次观测试验期间的加密资料对结果进行检验。图 2a给出了进行检验所用的加密观测站点分布，包括风廓线、微波辐射仪和铁塔。图 2b是这次风暴过程40 dBZ以上强回波演变过程，可以看出，从14:29(北京时，下同) 开始，风暴在从西北向东南下山传播的过程中组织性逐渐增强，范围逐渐增大，18:29达到最强，然后从19:29开始逐渐分散并最终减弱，风暴持续6 h以上。对比图 2a和图 2b也可发现，观测试验仪器的分布基本位于风暴传播路径上，因此, 观测数据对试验结果的检验对比非常有利。
|
|
| 图 2. 系统基本配置及观测信息 (a) 系统运算范围 (长虚线框)、结果分析范围 (短虚线框) 及观测站点位置 (加号表示天气雷达，空心三角表示风廓线，实心三角表示微波辐射仪，叉号表示观测塔; 图中粗实线为200 m地形等高线)，(b)40 dBZ以上组合反射率因子拼图指示2009年7月22日对流风暴在不同时间的强回波形态及演变过程 Fig 2. Basic system configuration and observational information (a) domain of system running (long dash line box), domain of result analysis (short dash line box), and location of observation sites (the plus denotes CINRAD radar site, the open triangle denotes wind profiler, the closed triangle denotes microwave radiometer, the multiplication sign denotes observation tower, thick black line indicates 200 m elevation), (b) patterns and evolution of strong storm echoes at different time on 22 Jul 2009 using mosaic of composite reflectivity no less than 40 dBZ | |
2.2 试验方案设计
除特别说明，在所有试验中模式均从2009年7月22日14:29积分开始，20:29积分结束，基本覆盖这次强风暴演变过程。如引言所述，本文所关注的是与风暴生消发展密切相关的低层三维动力和热力场，且最终目标是建立实时运行的系统，这就要求确保每个RR4DVar循环在其墙钟时间内计算完成。为节约计算时间，默认为低层配置。在所有试验中，模拟范围为540 km×540 km (图 2a中长虚线)，基本覆盖风暴发展区域。模式水平分辨率为3 km，垂直分辨率为375 m。除了全对流层试验，其他试验垂直均为15层，其中，在3 km以下同化雷达资料，而3 km以上为海绵边界层，由背景场、4DVar中的云模式积分动力调整及重力波阻尼算法控制，用来抑制重力波对模式上边界的影响和反射效应，同时确保有无雷达资料同化垂直层间的平滑过渡。模式最低层为垂直分辨率一半, 即187.5 m。根据测试，在代价函数最小化迭代达45次时，代价函数梯度下降效率显著降低，设定迭代次数为50后输出结果，并利用云模式进行6 min预报，用于计算下一热启动循环背景场。试验共4组：控制试验、全对流层对比试验、背景场试验、雷达资料同化试验 (表 1)。根据以往的调试经验，在此配置1组默认最佳参数作为控制试验，将其作为对比基准。
|
表 1 敏感性试验设置 Table 1 Configuration of sensitivity experiments |
2.3 试验结果
选取2009年7月22日强风暴个例模拟的3个时次15:29，18:29和20:17，分别代表风暴的发展、成熟和减弱阶段，进行不同试验结果对比分析，主要聚焦于代表低层动力和热力特征变化的风场和扰动温度场。
2.3.1 控制试验从控制试验EXP_CTRL最低层187.5 m的风场和扰动温度场看 (图 3)，在风暴整个发展演变过程中，前方平原地区低层环境风场为较强偏南风，带来充沛的低层暖湿空气，有利于这一地区长生命史、有组织对流系统的发展[34-35]。在风暴发展增强阶段 (图 3a)，降水导致的近地面冷池已形成，中心最低扰动温度为-10.3℃，风暴前沿风场也出现出流特征。风暴在下山过程中强度明显增大，这与地形强迫下的低层偏西南暖湿气流在山前抬升关系密切[34-35]。风暴在下山后逐渐发展成熟 (图 3b)，随着降水的增加和持续，雨水蒸发冷却导致的风暴后部冷空气不断下沉，使得冷池强度和范围均显著增加，中心最低扰动温度达到-15.4℃，风暴出流加强。风暴减弱阶段 (图 3c)，冷池强度减弱 (中心最低扰动温度升至-11.1℃)，但范围扩大，出流风场也向外不断扩展。从定性分析看，EXP_CTRL试验得到的低层动力和热力场与该风暴的演变密切相关，冷池导致的低层出流与近风暴低层环境入流在风暴前沿形成强的低层暖湿空气辐合上升，有利于风暴前新对流单体的生成和整个风暴系统的发展、传播及生命史延续。
|
|
| 图 3. EXP_CTRL试验模拟的187.5 m风场 (矢量) 和扰动温度场 (填色)(白色等值线为40 dBZ以上雷达反射率因子观测值，间隔为10 dBZ) (a)15：29，(b)18：29，(c)20：17 Fig 3. Simulated winds (vectors) and perturbation temperature (the shaded) of 187.5 m level at 1529 BT (a), 1829 BT (b) and 2017 BT (c) from EXP_CTRL experiment (white contours indicate observations of radar reflectivity no less than 40 dBZ with 10 dBZ interval) | |
2.3.2 全对流层试验
为检查简化Kessler微物理方案和低层配置对低层动力和热力场的影响，执行全对流层对比试验EXP_DEEP (表 1)。从最低层187.5 m风场 (图略) 看，除了一些小尺度扰动外，水平风场并未有明显差异，平原地区偏西南气流及风暴前沿出流大风特征均与EXP_CTRL试验结果类似。从扰动温度场看 (图略)，低层温度场形势与EXP_CTRL试验非常接近，冷池结构也得到很好模拟，但冷池强度稍有偏差，在风暴成熟阶段冷池中心最大扰动温度偏高0.9℃，这主要是因为EXP_CTRL试验与EXP_DEEP试验的模式微物理方案设置存在一定差异造成的。试验表明，如聚焦模式积分时间很短的3 km以下低层动力场和热力场模拟分析 (而非长时间预报) 时，低层配置和简化微物理方案完全满足这一需求。后文观测检验也证明了低层设置的合理性。该试验对比结果也与Sun等[22]试验结论类似。
2.3.3 背景场试验背景场试验主要考察地面AWS观测、BJ-RUC模式结果及中尺度风场动力约束的影响 (表 1)。从NO_AWS试验的187.5 m模拟结果 (图 4) 看，如无地面观测，低层风偏差明显增加，风暴前平原地区低层风向转为偏西风，风暴消散减弱阶段的出流风场也异常偏强。低层环境温度偏低，冷暖结构不明显，特别是冷池强度显著减弱, 形态也发生明显变化。可见，足够密集的AWS资料对雷达资料4DVar同化及近地面热、动力分析非常重要，特别是这样的低层RR4DVar同化分析，融合稠密的高频次AWS观测能有效弥补天气雷达对近地面观测的不足，从而获得更为可信的、高时空分辨率的低层动力和热力三维分析场。
|
|
| 图 4. NO_AWS试验模拟的187.5 m风场 (矢量) 和扰动温度场 (填色)(其他说明同图 3) Fig 4. Simulated winds (vectors) and perturbation temperature (the shaded) of 187.5 m level from NO_AWS experiment (others same as in Fig. 3) | |
从CBG_NORUC试验模拟的187.5 m风场和扰动温度场 (图略) 看，尽管使用了与冷启动时间最为接近的探空资料，在风暴发展阶段 (15:29)，其与EXP_CTRL试验结果差异也非常明显，低层环境风场存在偏差，风暴前沿出流偏差明显，冷池则显著偏弱。在风暴成熟及减弱阶段，随着热启动循环不断调整，模拟结果差异明显减小。总之，在冷启动计算中尺度背景场 (初猜场) 时，如采用时空分辨率很稀疏的探空资料，虽然也融合了雷达VAD及地面AWS资料，但不论从时间上还是空间上，对背景场的分析偏差均比使用较为可靠的逐小时更新的BJ-RUC模式探空结果的偏差要大很多，从而影响到对随后几个热启动循环的数值模拟。在热启动时，因为使用上一循环云模式的6 min预报作为初猜场，同时融合地面AWS观测、雷达VAD风廓线和BJ-RUC模式风场得到背景场，云模式6 min预报又是在上一循环同化6部雷达资料后进行的，因此, 在几个热启动循环后，冷启动背景场对模拟结果的影响逐渐减小。CBG_RUCOUT试验用于检验在冷启动时单独使用BJ-RUC模式结果计算背景场的效果 (图略)。在距离冷启动时刻较近的风暴发展阶段的15:29，与EXP_CTRL试验相比，低层环境风场和温度场结构差异明显，冷池也有所偏强。可见，在冷启动时单独使用中尺度模式资料而不考虑地面观测和雷达VAD风廓线进行背景场计算，会对随后几个循环的雷达资料4DVar同化和低层热、动力结果模拟分析产生明显负面影响。另外，研究也表明，相比直接插值中尺度模式结果计算背景场讲，抽取模式探空计算背景场能够导致雷达资料4DVar同化迭代更好地收敛[22]。与CBG_NORUC试验类似，由于热启动循环的不断调整，几个循环后的低层风场、温度场与EXP_CTRL试验结果差异不明显，受冷启动背景场的影响也明显减小。
随后4个试验是测试热启动中尺度背景场调整的影响。在这些试验里，为突出验证中尺度背景场的动力约束作用，设定不同λ仅对风场起约束作用。从WBG_MAWIND0结果 (图 5) 可以看出，低层风场小尺度扰动较多，风速偏强，风向偏差较大。随着热启动循环次数的增加，低层风和扰动温度的误差逐渐增大，特别是在东南平原地区和渤海湾附近，对流风暴成熟和减弱阶段的风暴前环境低层风向与EXP_CTRL结果完全相反，风速则异常偏大；扰动温度场也表明冷池范围有所扩大，但强度明显偏弱。即在没有中尺度背景风场动力约束下，RR4DVar热启动循环误差累积导致低层动力和热力分析场出现明显偏差。从WBG_MAWIND20结果 (图略) 看，引入一定比例的中尺度风场进行动力约束后，模拟得到显著改善。从WBG_MAWIND60结果 (图略) 看，比较接近于EXP_CTRL试验结果，主要差异是风场过度平滑，风暴出流稍微偏小，冷池强度也出现一定偏差。可见给予中尺度分析风场更多的权重后，与风暴发展演变相关的一些小尺度信息可能会被平滑和过滤掉。从WBG_MAWIND100结果 (图略) 看，风暴前出流强度和冷池范围明显偏小，在风暴成熟阶段的冷池强度显著减弱。这是由于中尺度分析得到的水平风场不能有效反映对流尺度信息，因此，融合来自云模式上一循环预报风场非常必要。更多的试验和分析也表明λ值0.3~0.5最适合，这样，恰当的中尺度动力约束既可以平滑掉不必要的模式小尺度累积噪声，又可以保留与风暴发展演变密切相关的出流、冷池等重要的对流尺度信息。
|
|
| 图 5. WBG_MAWIND0试验模拟的187.5 m风场 (矢量) 和扰动温度场 (填色) (其他说明同图 3) Fig 5. Simulated winds (vectors) and perturbation temperature (the shaded) of 187.5 m level from WBG_MAWIND0 experiment (others same as in Fig. 3) | |
2.3.4 雷达资料同化试验
在NO_CLOSERADAR2C试验中，不同化张北和承德两部C波段雷达资料。结果显示，在上述两个雷达站附近的低层风场和扰动温度场出现偏差，但在远离这两个雷达站的地方，低层动力场、热力场结构特征与EXP_CTRL试验结果非常相似 (图略)。风暴发展初期，因风暴距离两部C波段雷达相对较近，缺失这两部雷达资料导致风暴前出流减小，冷池强度减弱。随后，风暴向东南传播和移动，逐渐离开这两部雷达的有效探测范围，对与风暴密切相关的冷池和出流等的影响逐渐减小。说明至少对该系统而言，有效同化C波段雷达资料对对流尺度动力和热力场的模拟具有正面效果。从风暴移动路径来看，北京和天津雷达对风暴的整个发展演变观测最为有效，通过NO_CLOSERADAR2S试验分析不同化上述两部S波段雷达资料的影响。由图 6可知，低层风场和扰动温度场明显变差，冷池结构未合理模拟出来，冷池强度显著减弱；风暴附近低层风存在明显偏差，风暴前出流显著减弱。对比发现，至少对RR4DVar来说，S波段雷达资料同化效果明显优于C波段雷达，这与两种型号雷达观测性能差异存在直接关系。一般地，S波段雷达信号衰减明显低于C波段雷达；S波段雷达的10 cm波长可以探测由于扰动梯度形成的Bragg散射，可获得更多的有效观测数据；S波段和C波段雷达在最低几个仰角的最大不模糊速度一般分别为±27.1 m·s-1和±13.55 m·s-1，S波段雷达多普勒优势明显，导致能够被有效同化的径向速度资料明显多于C波段雷达。而对于雷达资料变分同化 (3DVar/4DVar) 技术而言，其优势更多依赖于雷达径向速度观测，因此，缺失这两部S波段雷达资料，对低层动力场、热力场分析的影响必然更明显。即对于对流尺度数值模拟来说，有效同化雷达资料 (特别是径向速度) 非常重要。从NO_ALLRADAR结果 (图略) 可以看出，如不同化雷达资料模拟结果差异更大, 15:29的冷池未被分析出来，其他时次冷池明显偏弱；低层风速明显偏小，不能模拟出风暴出流。这说明不进行雷达资料RR4DVar同化，能定性模拟出基本的大气中尺度特征，但很难模拟出与对流风暴发展演变密切相关的出流、冷池等对流尺度结构。在北京“7.21”特大暴雨数值模拟中也能看到雷达资料RR4DVar的优势[25]。
|
|
| 图 6. NO_CLOSERADAR2S试验模拟的187.5 m风场 (矢量) 和扰动温度场 (填色)(其他说明同图 3) Fig 6. Simulated winds (vectors) and perturbation temperature (the shaded) of 187.5 m level from NO_CLOSEERADAR2S experiment (others same as in Fig. 3) | |
在RADAR_3V和RADAR_6V试验中，分别设定同化每部雷达3组和6组体扫资料 (图 1)。需要说明的是，同化窗时段内未将第1组雷达资料6 min体扫时间计算在内。另外，为得到上述3个时间点的模拟结果，两个试验模式积分时间从14:17开始。从RADAR_3V试验 (图略) 定性看，与EXP_CTRL比较接近，但730 s同化窗导致分析时间间隔增加为18 min，这对于短生命史或快速演变的风暴来说，有时不能得到精确的动力、热力演变特征。从RADAR_6V试验 (图略) 看，与EXP_CTRL试验存在明显差异，冷池范围扩大 (冷池前沿已经伸展到风暴之前)，但中心强度降低，特别是在风暴发展成熟阶段；低层风偏大，出流异常偏强，低层辐合明显加强。当4DVar同化窗加长时，相当于在1个同化循环中引入过多的连续观测信息。就该试验而言，相当于30 min内每部雷达逐6 min观测的风暴短时演变信息在1个循环内被同化，导致冷池范围异常扩展但中心强度下降，也就是说，30 min内雷达连续观测信息与短时间内风暴快速发展演变的对流尺度特征不协调。另外，随着4DVar同化时间窗长度的增加，模式误差累积和观测相关性也可能导致分析质量变差。所以在1个RR4DVar同化窗内，同化每部雷达2~3组观测设置合理。
3 试验结果对比检验针对2009年7月22日强风暴个例，对14:29—20:29各个试验的模拟结果进行检验，包括风速、风向和温度廓线的站点的检验，以及风速、风向和温度序列的站点对比分析。
3.1 风和温度廓线检验首先，利用位于北京延庆、北京海淀、北京市观象台和天津塘沽的风廓线雷达逐6 min间隔观测进行敏感性试验的风场统计检验 (包括风速和风向)。根据试验设置，仅考虑低层2.8125 km以下结果的检验。在2.8125 km以下，上述4部风廓线的垂直观测间隔为50~240 m，比模式垂直层密集，因此，首先利用3点滑动平均法，对风廓线观测实施平滑。然后利用样条函数和距离权重垂直插值方法，将风廓线资料垂直插值到模式层上。再利用水平距离权重插值方法，将距离风廓线站点最近的4个模式格点资料插值到风廓线站点上。在时间上，取6 min间隔风廓线观测与试验模拟结果最接近时次作为二者对比的时次。这样，得到模拟与观测的廓线资料并进行相应检验，包括偏差 (BS) 和均方根误差 (RMSE)。另外，因海拔高度原因，北京延庆风廓线检验从模式第2层即562.5 m开始。从风速检验 (图 7) 和风向检验 (图 8) 来看，与定性分析基本一致。EXP_CTRL试验较为理想，最低层风速偏差和均方根误差分别在-0.5 m·s-1和1.3 m·s-1以内，风向偏差和均方根误差分别在8°和14°以内。EXP_DEEP试验和RADAR_3V试验与EXP_CTRL试验类似。没有雷达观测被同化的NO_ALLRADAR试验结果偏差最大，最低层风速偏差和均方根误差分别超过-2.6 m/s和4.2 m/s，风向偏差和均方根误差分别超过33°和68°。不同化S波段雷达资料的NO_CLOSERADAR2S试验结果偏差次之，最低层风速偏差和均方根误差分别超过-2.3 m·s-1和3.6 m·s-1，风向偏差和均方根误差分别超过30°和60°。不使用地面AWS观测的NO_AWS结果在最低几层也存在明显偏差，尤其最低层风向偏差相比其他试验最大，偏差和均方根误差分别超过36°和78°。不恰当的中尺度模式背景风场约束、不同化两部C波段雷达资料、过长的雷达资料同化窗的试验结果偏差也相对较为明显，而其他试验结果偏差较小。从风场检验也可看出，与大部分试验结果不同的是，缺失雷达观测的NO_ALLRADAR试验和NO_CLOSERADAR2S试验的风速偏差随着高度逐渐减小，这表明在RR4DVar同化策略下，雷达资料对低层动力场分析最为重要，也说明了因为天气雷达探测特征，决定了有效同化最低几个仰角面上雷达资料的效果最为明显。
|
|
| 图 7. 试验结果与风廓线雷达观测的风速偏差 (BS，实线) 和均方根误差 (RMSE，虚线) Fig 7. Bias (BS, solid lines) and root mean square error (RMSE, dashed lines) of wind speed profiles between results from experiments and observations from wind profiler radars | |
|
|
| 图 8. 试验结果与风廓线雷达观测的风向偏差 (BS，实线) 和均方根误差 (RMSE，虚线) Fig 8. Bias (BS, solid lines) and root mean square error (RMSE, dashed lines) of wind direction profiles between results from experiments and observations from wind profiler radars | |
其次，使用北京顺义和天津塘沽的地基微波辐射仪观测，进行温度廓线统计检验。对温度观测和模拟廓线的时空处理和检验方法与风廓线站点类似。因为微波辐射仪探测时间间隔为1 min，所以在模拟结果时间点上，将模拟与观测进行统计检验即可。从温度检验 (图 9) 看，与风速风向检验类似，在不同高度上，EXP_CTRL，EXP_DEEP试验和RADAR_3V试验结果偏差和均方根误差相对最小，最低层偏差和均方根误差分别在-0.5℃和-1℃以内。没有雷达观测被同化的NO_ALLRADAR试验结果偏差最明显，其中最低层的偏差和均方根误差分别超过-2.4℃和4.7℃。不同化两部S波段雷达资料的NO_CLOSERADAR2S试验结果偏差也较明显，最低层偏差和均方根误差分别超过-2.2℃和4.3℃。不使用地面AWS观测的NO_AWS试验结果在最低几层的偏差也非常大，尤其在最低层，位于顺义站的偏差和均方根误差分别达到2.8℃和5.4℃。
|
|
| 图 9. 试验结果与微波辐射仪观测的温度偏差 (BS，实线) 和均方根误差 (RMSE，虚线) Fig 9. Bias (BS, solid lines) and root mean square error (RMSE, dashed lines) of temperature profiles between results from experiments and observations from microwave radiometers | |
从NO_AWS试验的风速、风向和温度检验对比看，随着高度增加，偏差和均方根误差逐渐减小，地面AWS资料对模拟结果的影响可达到2000 m以上，在1000 m以下则较为明显，特别是对最低层的影响最为显著。通过对比本文EXP_CTRL试验的检验结果与以往的VDRAS检验结果[22-23]看，目前发展的雷达四维变分分析系统对低层热、动力场的分析效果已进一步改善，风速、风向和温度廓线的偏差和均方根误差均有不同程度减小。
需要指出，风廓线雷达本身也存在一定的观测误差，虽然风廓线雷达具有特定的观测数据质量控制算法，但由于大气不均匀、垂直上升及强降水等导致的观测误差依然存在[36]。微波辐射技术对大气温度探测也存在一定偏差。另外, 廓线资料插值也会带来误差。这导致该检验具有一定不确定性, 但作为对各敏感性试验结果相对好坏的定量评估来说，上述检验有显著意义。
3.2 风和温度序列对比利用位于北京中国科学院大气物理研究所 (以下简称大气所) 和位于天津市气象局的边界层观测铁塔数据，针对全部试验，进行模式最低层187.5 m风速、风向和温度时间序列对比。利用水平距离权重插值方法，将距离观测塔站点最近的4个模式格点的风速、风向和温度资料插值到观测塔站点上。观测选用最接近模式最低层高度的边界层观测塔180 m高度数据。其中，大气所观测塔风速使用东南和西北两个观测臂获得的最大风速值，风速取10 min内最大风速和10 min内平均风速两个结果，风向则取10 min内最大风速对应风向及10 min内风向频率最高象限风向平均值。而对于温度，大气所观测塔数据取30 min平均值。因为仅能得到天津市气象局观测塔数据1 h平均值，所以仅进行温度对比。
首先，分析大气所观测塔的风速、风向随时间的变化 (图略)。17:00左右，也就是当对流风暴经过观测塔站点时，因为风暴伴随的局地短时大风的影响，观测的180 m风速和风向均出现剧烈变化，风速剧增了8 m·s-1以上，风向发生了100°以上的突变。从风速、风向观测的10 min平均值和10 min最大值看，不论是风速还是风向，二者趋势基本一致，也就是说，10 min平均风速、风向能够代表风速、风向的实际变化趋势。在上述试验中，雷达资料同化窗的长度决定了大部分试验结果为12 min内平均值，而RADAR_3V和RADAR_6V试验结果更是反映了18 min和36 min内平均值。因此，对于风速、风向，利用观测塔180 m高度的10 min平均值与模式最低层187.5 m分析结果进行对比更加合理。可以看出，不论是风速还是风向，EXP_CTRL, EXP_DEEP和RADAR_3V试验结果最接近10 min平均观测值，17:00左右当风暴经过观测塔时出现的风速剧增和风向突变特征在上述4个试验中也得到很好反映。与前述风廓线检验类似，NO_ALLRADAR, NO_CLOSERADAR2S，NO_AWS 3个试验结果与观测之间偏差最大，风速明显偏小，风向变化与观测也明显不一致。从冷启动背景场试验CBG_NORUC和CBG_RUCOUT结果看，冷启动中尺度背景场的影响只在前面几个循环比较明显，随着模式调整，后续循环的模拟结果与EXP_CTRL试验结果非常接近。从中尺度模式背景风场约束试验看，不合适的中尺度动力场权重约束结果也与EXP_CTRL试验结果及观测偏差较大。过长的雷达资料同化窗试验RADAR_6V结果也存在较为明显的偏差。不同化两部C波段雷达资料的NO_CLOSERADAR2C试验的风速虽然偏差相对较小，但风向与EXP_CTRL试验相比，却存在一定偏差。
其次，分析两个观测塔站点温度随时间变化 (图略)。由观测可知，17:00左右在大气所站点、19:30左右在天津市气象局站点，观测的180 m高度的温度均出现陡降，两个站点分别下降了7℃和10℃，这是因为此时对流风暴经过站点，风暴形成的冷池导致温度剧烈下降。从试验结果看，EXP_CTRL, EXP_DEEP和RADAR_3V试验结果与两个站点观测最为接近，而NO_ALLRADAR，NO_CLOSE-RADAR2S，NO_AWS 3个试验结果与观测相比出现显著偏高，特别是在大气所站点，温度陡降特征不明显，表明这3个试验没有很好分析出风暴冷池特征。其余试验或多或少也存在比EXP_CTRL试验结果更大的偏差。总之，针对所有试验来看，温度与风速、风向的对比结果类似，这与前面定性分析和廓线检验结论基本一致。
最后需要强调，不同化任何雷达资料的NO_ALLRADAR试验，实际上反映的是中尺度背景场的动力降尺度结果。不论是从定性对比分析还是定量检验来看，NO_ALLRADAR试验结果偏差更大。可以说，一方面，实现对雷达观测的RR4DVar同化，对于边界层对流尺度动力和热力特征模拟分析具有明显正面效果；另一方面，说明另外两部S波段雷达 (石家庄和秦皇岛) 资料对RR4DVar同化和分析结果也具有比较明显的正面影响。
4 小结本文在VDRAS基础上，基于雷达资料RR4DVar同化方案和三维数值云模式，研发建立了一个雷达四维变分分析系统，可快速模拟得出与强对流生消发展密切相关的区域低层大气三维动力、热力和微物理场精细结构。主要结论如下：
1) 采用雷达资料RR4DVar同化，能够模拟出环境风场结构、低层近风暴环境入流、低层小尺度辐合和暖区等风暴环境特征，以及风暴形成的冷池、出流等与风暴本身发展传播和单体新生密切相关的对流尺度特征。与观测塔数据的对比表明，合理配置的模拟结果能够精确反映出风暴经过观测站时的风速剧增、风向突变及温度陡降等。
2) 如果仅聚焦与局地突发强对流密切相关的低层三维动力和热力结构模拟分析，简化微物理方案的低层配置与全微物理方案的全对流层配置的模拟结果相近。
3) 地面AWS资料对RR4DVar同化和模式最低几层热、动力特征模拟分析非常重要，AWS观测能够有效弥补天气雷达对近地面观测的不足，从而获得更为可信的低层对流尺度动力和热力三维分析场。
4) 冷启动中尺度背景场的改变对冷启动及随后几个热启动循环模拟结果产生明显影响。在每个RR4DVar同化模拟循环中，利用合理的中尺度动力约束，不但可平滑掉不必要的模式小尺度累积噪声，还可保留与对流风暴发展演变密切相关的对流尺度特征。
5) 有效的雷达资料RR4DVar同化，是模拟出与对流风暴发展演变密切相关的对流尺度特征的关键所在。就目前设置而言，S波段雷达比C波段雷达的资料同化效果更明显，且同化每个雷达2~3组观测可获得最为合理的模拟分析结果。
目前，该系统已成为聚焦京津冀0~12 h短时临近预报的快速更新多尺度分析和预报系统 (RMAPS) 的重要组成部分。在后续研究中，将对京津冀多个典型强对流天气个例进行模拟，对获得的低层三维动力场和热力场进行分析和检验，进一步评估系统在对流尺度数值模拟及强对流天气临近预报预警中应用的可靠性。另外，也将基于雷达资料RR4DVar初始化方案，开展对流尺度数值临近预报研究。
致谢 美国国家大气研究中心 (NCAR) 资深科学家James Wilson对模拟结果分析提供了良好建议。谨此致谢。| [1] | 陈明轩, 俞小鼎, 谭晓光, 等. 对流天气临近预报技术的发展与研究进展. 应用气象学报, 2004, 15, (6): 754–766. |
| [2] | 俞小鼎, 周小刚, 王秀明. 雷暴与强对流临近天气预报技术进展. 气象学报, 2012, 70, (3): 311–337. DOI:10.11676/qxxb2012.030 |
| [3] | 郑永光, 周康辉, 盛杰, 等. 强对流天气监测预报预警技术进展. 应用气象学报, 2015, 26, (6): 641–657. DOI:10.11898/1001-7313.20150601 |
| [4] | Wilson J W, Feng Y, Chen M, et al. Nowcasting challenges during the Beijing Olympics:Successes, failures, and implications for future nowcasting systems. Wea Forecasting, 2010, 25: 1691–1714. DOI:10.1175/2010WAF2222417.1 |
| [5] | Houze R A Jr, Smull B F, Dodge P. Mesoscale organization of springtime rainstorms in Oklahoma. Mon Wea Rev, 1990, 118: 613–654. DOI:10.1175/1520-0493(1990)118<0613:MOOSRI>2.0.CO;2 |
| [6] | Weisman M L, Rotunno R. The use of vertical wind shear versus helicity in interpreting supercell dynamics. J Atmos Sci, 2000, 57: 1452–1472. DOI:10.1175/1520-0469(2000)057<1452:TUOVWS>2.0.CO;2 |
| [7] | Droegemeier K K, Wilhelmson R B. Numerical simulation of thunderstorm outflow dynamics.Part Ⅰ:Outflow sensitivity experiments and turbulence dynamics. J Atmos Sci, 1987, 44: 1180–1210. DOI:10.1175/1520-0469(1987)044<1180:NSOTOD>2.0.CO;2 |
| [8] | Mueller C K, Carbone R E. Dynamics of a thunderstorm outflow. J Atmos Sci, 1987, 44: 1879–1898. DOI:10.1175/1520-0469(1987)044<1879:DOATO>2.0.CO;2 |
| [9] | Rotunno R, Klemp J B, Weisman M L. A theory for strong, long-lived squall lines. J Atmos Sci, 1988, 45: 463–485. DOI:10.1175/1520-0469(1988)045<0463:ATFSLL>2.0.CO;2 |
| [10] | Coniglio M C, Corfidi S F, Kain J S. Views on applying RKW theory:An illustration using the 8 May 2009 derecho-producing convective system. Mon Wea Rev, 2012, 140: 1023–1043. DOI:10.1175/MWR-D-11-00026.1 |
| [11] | 陈明轩, 王迎春. 低层垂直风切变和冷池相互作用影响华北地区一次飑线过程发展维持的数值模拟. 气象学报, 2012, 70, (3): 371–386. DOI:10.11676/qxxb2012.033 |
| [12] | Wilson J W, Megenhardt D L. Thunderstorm initiation, organization and lifetime associated with Florida boundary layer convergence lines. Mon Wea Rev, 1997, 125: 1507–1525. DOI:10.1175/1520-0493(1997)125<1507:TIOALA>2.0.CO;2 |
| [13] | 陈明轩, 高峰, 孔荣, 等. 自动临近预报系统及其在北京奥运期间的应用. 应用气象学报, 2010, 21, (4): 395–404. DOI:10.11898/1001-7313.20100402 |
| [14] | 王彦, 于莉莉, 李艳伟, 等. 边界层辐合线对强对流系统形成和发展的作用. 应用气象学报, 2011, 22, (6): 724–731. DOI:10.11898/1001-7313.20110610 |
| [15] | Sun J, Xue M, Wilson J. Use of NWP for nowcasting convective precipitation:Recent progress and chuenges. Bull Amer Meteor Soc, 2014, 95: 409–426. DOI:10.1175/BAMS-D-11-00263.1 |
| [16] | Sun J, Crook N A. Dynamical and microphysical retrieval from Doppler radar observations using a cloud model and its adjoint:Ⅰ.model development and simulated data experiments. J Atmos Sci, 1997, 54: 1642–1661. DOI:10.1175/1520-0469(1997)054<1642:DAMRFD>2.0.CO;2 |
| [17] | Sun J, Crook N A. Dynamical and microphysical retrieval from Doppler radar observations using a cloud model and its adjoint:Ⅱ.Retrieval experiments of an observed Florida convective storm. J Atmos Sci, 1998, 55: 835–852. DOI:10.1175/1520-0469(1998)055<0835:DAMRFD>2.0.CO;2 |
| [18] | Sun J, Crook N A. Real-time low-level wind and temperature analysis using single WSR-88D data. Wea Forecasting, 2001, 16: 117–132. DOI:10.1175/1520-0434(2001)016<0117:RTLLWA>2.0.CO;2 |
| [19] | Crook N A, Sun J. Analysis and forecasting of the low-level wind during the Sydney 2000 Forecast Demonstration Project. Wea Forecasting, 2004, 19: 151–167. DOI:10.1175/1520-0434(2004)019<0151:AAFOTL>2.0.CO;2 |
| [20] | Chen M X, Sun J, Wang Y C. A Frequent-updating High-resolution Analysis System Based on Radar Data for the 2008 Summer Olympics.Preprint. The 33rd International Conference on Radar Meteorology, Cairns, Australia, 2007. |
| [21] | Sun J, Zhang Y. Analysis and prediction of a squall line observed during IHOP using multiple WSR-88D observations. Mon Wea Rev, 2008, 136: 2364–2388. DOI:10.1175/2007MWR2205.1 |
| [22] | Sun J, Chen M X, Wang Y C. A frequent-updating analysis system based on radar, surface, and mesoscale model data for the Beijing 2008 Forecast Demonstration Project. Wea Forecasting, 2010, 25: 1715–1735. DOI:10.1175/2010WAF2222336.1 |
| [23] | 陈明轩, 王迎春, 高峰, 等. 基于雷达资料4DVar的低层热动力反演系统及其在北京奥运期间的初步应用分析. 气象学报, 2011, 69, (1): 64–78. DOI:10.11676/qxxb2011.006 |
| [24] | 陈明轩, 王迎春, 肖现, 等. 基于雷达资料四维变分同化和三维云模式对一次超级单体风暴发展维持热动力机制的模拟分析. 大气科学, 2012, 36, (5): 929–944. DOI:10.3878/j.issn.1006-9895.2012.11132 |
| [25] | 陈明轩, 王迎春, 肖现, 等. 北京"7.21"暴雨雨团的发生和传播机理. 气象学报, 2013, 71, (4): 569–592. DOI:10.11676/qxxb2013.053 |
| [26] | 孙继松, 何娜, 郭锐, 等. 多单体雷暴的形变与列车效应传播机制. 大气科学, 2013, 37, (1): 137–148. DOI:10.3878/j.issn.1006-9895.2012.12015 |
| [27] | 肖现, 王迎春, 陈明轩, 等. 基于雷达资料四维变分同化技术对北京地区一次下山突发性增强风暴热动力机制的模拟分析. 气象学报, 2013, 71, (5): 797–816. DOI:10.11676/qxxb2013.077 |
| [28] | 肖现, 陈明轩, 高峰, 等. 弱天气系统强迫下北京地区对流下山演变的热动力机制. 大气科学, 2015, 39, (1): 100–124. DOI:10.3878/j.issn.1006-9895.1403.13318 |
| [29] | Sun J. Convective-scale assimilation of radar data:Progress and challenges. Q J R Meteor Soc, 2005, 131: 3439–3463. DOI:10.1256/qj.05.149 |
| [30] | Lim E, Sun J. A velocity dealiasing technique using rapidly updated analysis from a four-dimensional variational doppler radar data assimilation system. J Atmos Oceanic Technol, 2010, 27: 1140–1152. DOI:10.1175/2010JTECHA1300.1 |
| [31] | Hayden C M, Purser R J. Recursive filter objective analysis of meteorological fields:Applications to NESDIS operational processing. J Appl Meteor, 1995, 34: 3–15. DOI:10.1175/1520-0450-34.1.3 |
| [32] | 陈敏, 范水勇, 郑祚芳, 等. 基于BJ-RUC系统的临近探空及其对强对流发生潜势预报的指示性能初探. 气象学报, 2011, 69, (1): 181–194. DOI:10.11676/qxxb2011.016 |
| [33] | 刘梦娟, 陈敏. BJ-RUC系统对北京夏季边界层的预报性能评估. 应用气象学报, 2014, 25, (2): 212–221. DOI:10.11898/1001-7313.20140211 |
| [34] | Chen M X, Wang Y C, Gao F, et al. Diurnal variations in convective storm activity over contiguous North China during the warm-season based on radar mosaic climatology. J Geophys Res, 2012, 117: D20115. DOI:10.1029/2012JD018158 |
| [35] | Chen M X, Wang Y C, Gao F, et al. Diurnal evolution and distribution of warm-season convective storms in different prevailing wind regimes over contiguous North China. J Geophys Res Atmos, 2014, 119: 2742–2763. DOI:10.1002/2013JD021145 |
| [36] | 邓闯, 阮征, 魏鸣, 等. 风廓线雷达测风精度评估. 应用气象学报, 2012, 23, (5): 523–533. DOI:10.11898/1001-7313.20120502 |