物理滤波初始化四维变分在临近预报中的应用

姜文静 梁旭东

(中国气象科学研究院灾害天气国家重点实验室， 北京 100081)

引 言

降水短时临近预报指0～12 h的预报，其中0～2 h预报为临近预报。20世纪60—80年代基于雷达外推的临近预报技术得到长足发展[1-4]。其中基于TREC(tracking radar echo by correlation，相关方法跟踪回波运动)技术[5-8]的追踪方法得到广泛应用，并发展出FRONTIER(forecasting rain optimized using new techniques of interactively enhanced radar and satellite，基于雷达和卫星交互增强新技术优化的雨量预报)[9],TITAN(thunderstorm identification, tracking, analysis and nowcasting，雷暴识别、跟踪、分析和临近预报)[10]等多种业务化系统。基于雷达回波的雷暴预报时效往往不超过3 h，其中主要原因是外推方法无法解决雷暴演变问题[11]。后来发展出将雷达回波外推与数值预报结合的方法[12]，Woflson等[13]指出，通常情况下，随着预报时效的延长，在外推预报性能降低的同时，数值预报的性能逐渐提高。但数值预报在临近预报的预报性能较低，如何提高数值预报的临近预报性能成为值得研究的问题，通过资料同化改进模式初始场是解决该问题的途径之一[14-15]。因此提高数值模式在0～3 h的预报水平是应用难点，也是本文的研究重点。俞小鼎等[16]在对国内外雷暴和强对流临近天气预报技术进行全面回顾后指出，临近预报系统应建立在风暴尺度数值预报模式集合预报基础上，采用适当的初值产生方法、模式不确定性处理方法以及合适的资料同化技术。随着高分辨率数值模式的发展和计算机计算水平的提高，越来越多的学者将数值模式应用到临近预报中，其中包括利用雷达资料和其他高分辨率的观测资料进行资料同化[17-18]，得到更精确的初始场，提高临近预报水平。该方法逐渐发展起来，但有待完善。Abhilash等[19]研究多普勒天气雷达反射率因子和径向速度对热带风暴及相关降雨事件短距离预报的影响，证明同化雷达反射率因子和径向速度引起的分析增量会影响动力场和热力场。雷达数据同化对短时降水预报有相当大的影响，而且对强降水(阈值大于7 mm)预报能力的改善比弱降水(阈值为1 mm和3 mm)更明显。

Takano等[20]通过研究各种初始化方法(包括用降水过程进行非线性正规模式初始化和利用根据降水观测估计的潜热强迫进行非线性正规模式初始化)，指出由于非绝热加热较弱，降水过程的非线性正态初始化并没有显著改变初始场，数值天气预报模式需要调整适应(spin-up)时间(几小时)才能产生足够的降水。初始场与模式的不协调是导致模式在预报初始阶段存在明显调整适应过程的主要原因之一，通过提高初始场与模式的协调性可以缩短模式调整适应时间，从而提高模式在短时临近时段的预报能力。Lynch等[21]提出DFI(digital filter initialization，数字滤波初始化)提高初始场与模式协调性以改进短时预报能力，该方法被Polavarapu等[22]以及Wee等[23]应用到四维变分中。另一类方法是采用物理约束提高初始场与模式的协调性[24-26]。Liang等[27]在三维变分同化技术的基础上提出MC-3DVar(model constrained three-dimensional variational analysis，模式约束三维变分技术)并应用到台风模式初始化中[28]。Peng等[29](以下称PW17)将MC-3DVar技术延伸到四维变分中，应用于WRF-4DVar系统，提出物理滤波初始化四维变分同化方法。PW17的研究表明：物理滤波初始化四维变分同化能够明显缩短模式调整适应时间，提高降水临近预报能力。该研究也指出，降低计算量和优化参数设置是将该方法用于临近预报需要研究的问题。

物理滤波初始化四维变分同化方案在传递背景场误差和观测信息的同时，遵循模式中的物理模型和动力模型，更好地协调变量之间关系，进而得到较平衡的分析场，缩短模式积分的调整适应时间。但临近预报对时效性要求很高，而PW17采用3 h的同化时间总长度，需要大量计算资源。另外，临近预报需要尽量使用最近的观测资料进行快速分析和预报。因此本研究将探索如何在PW17研究的基础上缩短同化窗口长度，减小计算量，并提高数值模式的临近预报性能。

1 数值模式和资料同化系统

1.1 同化方法和模式设置

参照PW17的研究方案，采用WRF模式以及基于WRFDA建立的物理滤波初始化四维变分同化和3DVar系统。PW17中的同化窗为3 h，约束项的权重系数σ取为4×10-5，本研究为缩短同化时间，将同化时间总长度设置为12 min，分3个时间窗。

背景场选取1°×1°的NCEP全球预报系统GFS的6 h预报场，模拟范围的中心点位置为37°N，115°E，水平网格数为150×150，水平分辨率为9 km。模拟试验采用的物理过程参数化方案包括Thompson微物理方案、RRTM长波辐射方案、Dudhia短波辐射方案、改进的MM5 Monin-Obukhov近地面层方案、Thermal Diffusion 陆面方案和Yonsei University边界层方案。

根据Sun等[30]的研究，采用以纬向风、经向风(u，v)分量代替流函数和势函数的选项更适合于中小尺度系统的资料同化，本研究中背景误差协方差矩阵B使用WRFDA中u，v分量独立的CV7选项。 矩阵B需要根据模式区域设置进行计算，2018年8月1日00：00至8月31日12:00(世界时，下同)每12 h进行12 h和24 h预报，预报结果使用WRFDA系统中基于NMC(National Meteorological Center，美国国家气象中心)方法[31]的gen_be模块生成背景误差协方差矩阵B，得到模拟范围内模式变量在模式动力和模式物理上相协调的误差统计结构[32]。

1.2 个例及同化参数

1.2.1 个例介绍

本研究选取2018年8月华北地区的17个降水个例(表1)，采用物理滤波初始化四维变分同化方法研究对临近预报的改进。每个降水个例的初始时刻为当日00：00或12:00，并根据观测资料，选取在预报起始时刻即有降水的个例。试验中降水分析范围比模式范围略小(30°～44°N，108°～122°E)以减小边界影响。

1.2.2 影响半径和方差权重参数

背景误差协方差的不同影响半径和方差权重参数可明显影响同化结果。WRFDA 系统中用于调节矩阵B影响半径(尺度因子)和方差权重(方差因子)的参数分别为len_scaling和var_scaling。len_scaling控制同化过程中观测信息在背景场中传递的空间尺度范围，len_scaling越大则观测信息传递范围越广。var_scaling控制同化过程中对背景场或观测场的信任程度，var_scaling越大则观测资料的权重越大，得到的分析场越接近观测资料(远离背景场)。

表1 2018年8月华北地区17个降水个例Table 1 Seventeen selected precipitation cases in North China in Aug 2018

为选取同化参数，本节选用降水强度大、范围广、强降水中心主要集中在观测资料密集的京津冀地区的2018年8月14日00:00降水个例。设置5个试验，分别采用不同的参数设置。通过对比6 h累积降水的TS评分(表略)，预报性能最好的方差因子var_scaling为1.5，尺度因子len_scaling为0.05。虽然采用大量样本试验可能得到更优的方差因子和尺度因子选择，但本研究的重点在于物理滤波初始化四维变分同化方法的适用性，因此本文采用以上取值开展试验。

1.2.3 约束项权重系数

物理滤波初始化四维变分同化代价函数中罚项的权重系数σ作用是平衡代价函数中分析场中模式变量与观测场或背景场的距离和罚项。如果权重系数设置为0，则物理滤波初始化四维变分同化方法和传统4DVar完全等价；权重系数过小或过大，会导致罚项的滤波作用不能发挥，或观测资料信息不能被同化。因此，合适的权重系数在物理滤波初始化四维变分同化的试验中也较为重要。基于2018年8月14日00:00降水个例，分别选取权重系数为0，4×10-9，4×10-8，4×10-7，4×10-6，4×10-5，4×10-4，4×10-3，4×10-2，4×10-1和1进行物理滤波初始化四维变分同化以及未同化试验。在6 h累积降水的TS评分中，评分较高的是权重系数为4×10-5和4×10-6的试验。结合逐1 h各量级降水的预报TS评分，本研究最终选取的约束项权重系数为4×10-6。

2 批量试验

2.1 试验设计

针对表1中17个降水个例进行物理滤波初始化四维变分同化试验，观测资料来自国家气象信息中心，包括高空定时观测(00：00或12:00)，地面逐分钟气压、温度和湿度观测，天气雷达6 min内的径向速度。雷达径向速度采用Liang等[33]的质量控制方法进行处理，检验分析采用地面站逐小时降水观测。模式和同化系统采用与本文第1章试验相同的设置。

2.2 批量试验结果

对17个降水个例，均采用与第1章相同的同化设置进行物理滤波初始化四维变分同化并作6 h预报(PFI试验)，同时与未进行同化的预报试验(CTL试验)和三维变分试验(3DV试验)进行比较，总体ETS(equitable threat score,去除随机事件后的公平风险评分)评分计算结果见表2和表3。

在前3 h的降水预报中，PFI试验除0.1 mm阈值的降水评分略低外，其他量级的降水预报评分明显高于CTL试验，特别是每小时超过13.0 mm的较强降水。物理滤波初始化四维变分同化试验的预报评分明显提高，表明在预报初始时段，物理滤波初始化四维变分同化试验已能够预报出比较强的降水，同样的结果也体现在6 h累积降水预报的评分中(表3)。在大于13.0 mm量级降水预报中，PFI试验有较好性能。对于较小量级的降水(4.0 mm及以下量级)CTL试验的预报效果更好，原因是物理滤波初始化四维变分同化后试验的预报降水较强，存在增加空报的情况。采用三维变分同化在第1,2小时的预报中，由于资料的引入导致初始场不协调，同化没有提升预报效果；直到第3小时，经过模式积分协调后，预报效果开始提升(表2第3小时)，这导致前6 h的整体评分比CTL试验偏低(表3)。批量试验结果说明物理滤波初始化四维变分同化方法明显提高临近时段降水预报效果，下面将选取降水个例进行详细分析，探究物理滤波初始化四维变分同化试验改进降水预报的具体原因。

表2 2018年8月17个降水个例批量试验的逐小时累积降水ETS评分Table 2 ETS scores of one-hour accumulated precipitation forecasts of seventeen cases in Aug 2018

表3 2018年8月17个降水个例批量试验的6 h累积降水ETS评分Table 3 ETS scores of six-hour accumulated precipitation forecasts of seventeen cases in Aug 2018

3 个例分析

3.1 2018年8月11日降水个例实况及观测资料

选取降水预报改进效果明显的2018年8月11日个例进行详细分析。8月11—12日受高空槽及低层切变影响，冷空气东移南下，在副热带高压北侧，西北地区东部、华北和东北一带冷暖交汇并形成降水，华北、东北等地出现较强降水天气过程。2018年8月11日00:00—12日00:00内蒙古东部偏南地区、陕西北部、山西北部、河北中北部、北京等地的部分地区有大雨或暴雨，局地伴有短时强降水、雷暴大风等强对流天气，小时降水量为20～50 mm，局地达70 mm。 2018年8月11日12：00—18:00华北地区6 h累积降水实况如图1所示。

图1 2018年8月11日12:00—18:00 6 h累积降水实况Fig.1 Observed six-hour accumulated precipitation during 1200-1800 UTC on 11 Aug 2018

3.2 同化试验降水预报能力

在以2018年8月11日12:00为初始时刻的物理滤波初始化四维变分同化试验中，同化迭代21步，使用墙钟时间25 min，同化后进行6 h预报。从目标函数下降曲线(图略)可以看到，迭代约10步后，目标函数下降已不明显，如果需要进一步减少同化计算时间，可以同化迭代10步左右，所需计算时间可缩短到约13 min。

未同化观测资料的对照试验(CTL试验)和物理滤波初始化四维变分同化的试验(PFI试验)前6 h 累积降水与观测对比如图2所示。

对比图2中CTL与PFI试验的降水分布可见，在整个模拟区域中，PFI和CTL试验的前6 h累积降水分布区域主要差别在河北地区(图中红圈A区)和北京地区(图中红圈B区)。在A区CTL试验存在明显漏报，而在B区存在明显空报。PFI试验预报A区降水，B区降水减小，更接近实况。两组试验降水预报产生差异的主要范围是37°～ 41°N，114°～118°E(北部的降水差异较明显，但由于靠近模式边界，因此不考虑)，并着重对这一区域(图中黑实线方框)进行分析。

图2 2018年8月11日观测与未同化(CTL)、同化(PFI)试验6 h累积降水预报Fig.2 Six-hour accumulated precipitation in observation,CTL experiment(without data assimilation),and PFI experiment(with data assimilation) on 11 Aug 2018

在0～3 h预报时段数值模式处于调整适应阶段，模式预报能力相对较低，而该时段正是提高数值模式临近预报能力的重点。选择图2中实线方框区域，分析不同试验前3 h累积降水(如图3所示)。

由图3可以看到，12:00—13:00最显著的特征是存在一条西南—东北向雨带，13:00—14:00北侧发展出第2条雨带，14:00—15:00存在两条雨带，均在40°N以南，北侧的雨带增强，而且两条雨带相接的西南侧，降水进一步增强。CTL试验中，雨带特征不明显，表现为比较集中的降水雨团，位置在39°N到40°N之间，并向北延伸。第1小时雨区开始出现在40°N，116°E附近，和观测相比，CTL试验的雨区总体偏北，导致南部降水明显漏报和北侧降水明显空报。PFI试验中，第1小时即出现西南—东北向的雨带特征，预报出CTL试验中漏报的南部降水，雨带比观测略偏北。第2小时雨带向北发展，与观测较一致，但是南侧的雨带较观测减弱更快。与CTL试验相比，PFI试验较好地预报出雨带特征，同时40°N，116°E附近的空报也明显减少。到第3小时，雨带维持，模拟区域南部的降水更接近观测，北部降水量较CTL试验小，空报减少。尽管位于西南部的雨带明显漏报，与CTL试验相比，PFI试验预报的降水在空间分布更接近观测。整个区域中物理滤波初始化四维变分同化试验比未同化试验的降水出现更早，降水范围也更接近于实况分布，特别是西南部区域，CTL试验始终未能预报出。同化资料后，减少了北部地区降水的空报。以上对降水预报的改进同样反映在ETS评分中(表4)。从ETS评分看，PFI试验评分在各量级上均高于CTL试验。

图3 2018年8月11日12:00—13:00，13:00—14:00，14:00—15:00观测与未同化(CTL)、同化试验(PFI)的逐小时累积降水分布Fig.3 One-hour accumulated precipitation in observation,CTL(control experiment without data assimilation) and PFI(data assimilation experiment) during 1200-1300 UTC,1300-1400 UTC and 1400-1500 UTC on 11 Aug 2018

表4 2018年8月11日个例1～3 h逐小时累积降水预报ETS评分Table 4 ETS scores of one-hour accumulated precipitation forecasts from the 1st to the 3rd hour on 11 Aug 2018

3.3 资料同化对初始场的影响

沿模拟区域降水分布线(图3中红色实线)从西南到东北方向的剖面(图4)可以看到，物理滤波初始化四维变分同化后的分析场比背景场的低层水汽更充足，特别是在700～850 hPa之间，分析场的相对湿度更大，同化后初始场的水汽明显增多，因此这些区域有利于降水的产生。

对比不同高度风场的背景场和分析场(图5和图6)可知，物理滤波初始化四维变分同化后，850 hPa高度分析增量有明显的偏北风，38°N，116.5°E附近散度增量场有增强的切变和辐合，这有利于北方的冷空气南压，使区域北部引发降水的槽线加深并南移，在38°N附近与南方暖湿空气交汇形成更强辐合，从而有利于雨区向南压。西南侧的辐散增强有利于减弱该地区的降水。在700 hPa和500 hPa散度增量场上，38°N附近辐散增强，流场变化可增强该区域低层辐合和高层辐散，更有利于降水出现，这与CTL试验中漏报的降水区域对应。另一方面，冷空气南压也有利于减弱CTL试验中北部区域空报的降水，40°N，116°E附近散度增量场在850 hPa 辐合减弱，700 hPa及500 hPa辐散减弱均有利于该地区降水的减弱(减少空报)。

图4 2018年8月11日个例同化试验背景场和分析场相对湿度剖面Fig.4 Section of relative humidity in background and analysis fields on 11 Aug 2018

图5 2018年8月11日850 hPa，700 hPa 和500 hPa风场的背景场和分析场以及分析增量Fig.5 Winds in the background, analysis fields and analysis increments at 850 hPa,700 hPa and 500 hPa on 11 Aug 2018

续图5

图6 同图5，但是为散度场Fig.6 The same as in Fig.5,but for divergence

4 降水预报的快速响应机制

为分析同化资料提高短时临近降水预报时效的主要原因，对于2018年8月11日降水个例分别分析模式中的格点可分辨率降水和次网格降水(图略)。从6 h累积降水的格点降水看，CTL和PFI试验中的降水相对较小，主要降水是次网格(积云对流参数化)降水。由此可见，该试验中对于短时临近降水预报时效影响较大的是积云对流参数化降水，而不是云微物理过程计算的降水。本研究还同时采用相同初始场进行更高分辨率(水平分辨率3 km)的模拟，不采用积云对流参数化，发现在短时临近时效阶段云微物理过程不能预报出38°N附近的降水(图略)。

积云参数化方案为Kain-Fritsch (new Eta)方案，是在Kain-Fritsch方案基础上改进而成。该方案利用一个伴有水汽上升下沉的简单云模式，考虑云中上升气流卷入和下曳气流卷出及相对简单的微物理过程影响。本文用物理滤波初始化四维变分方法同化地面和探空资料时对热力场有影响，但尚不足以对降水产生明显影响，而同化雷达径向风资料改变动力场后，通过积云参数化方案对降水预报影响较明显，且积云参数化方案有诊断的特性，能在短时间内响应，因此在模式预报的初始时刻就能产生降水。从逐小时累积降水的演变(图略)可以看到,第1 小时所有降水均来自积云对流参数化，第2小时和第3小时开始有格点降水(来自云微物理过程)，但降水强度小于1.5 mm·h-1。

对于批量试验，将2018年8月17个降水个例的模拟降水按照格点可分辨降水和次网格降水分别进行平均得到平均降水分布(图略)。由6 h累积降水平均分布可见，CTL和PFI试验中的格点降水均小于次网格降水。由此可见，批量试验中对短时临近预报降水影响较大的依然是积云对流参数化，同时还可看到，同化观测资料对格点降水预报效果也有改进。

由17个个例逐小时平均降水分布(图7)可以看到，同化资料对改进格点可分辨降水(PFI-RAINNC)有效，但大范围降水主要来自积云对流参数化降水(PFI-RAINC)。从批量试验结果也可看到，积云对流参数化降水对短时临近时效的降水影响较大(特别是降水分布)。Donner[34]介绍数值天气预报模型中积云对流参数化的初始化过程，发现积云对流通过加热高层大气和干燥低层大气直接影响大尺度环流，能够加快数值预报模式中的降水过程。以上研究表明：积云对流参数化具有快速响应的特性，对提高临近预报能力有一定帮助。

图7 2018年8月17个降水个例批量试验观测与未同化(CTL)、同化(PFI)试验的第1～3小时的逐小时平均累积降水中的积云对流降水(RAINC)和格点可分辨降水(RAINNC)分布Fig.7 Averaged convective parameterization(RAINC) and averaged grid-resolvable(RAINNC) precipitation accumulated within one-hour from the 1st to the 3rd hour in observation,CTL(without data assimilation),and PFI experiment(with data assimilation) of seventeen cases in Aug 2018

续图7

5 结论与讨论

本研究采用物理滤波四维变分同化技术进行同化试验以提高模式降水临近预报的能力。由于四维变分同化需要较大计算量，而临近预报需要在尽量短的时间给出预报结果，因此本研究通过缩短同化时间窗以减少同化所需的计算量。研究采用12 min的同化窗，并调试同化相关参数，进行个例和批量试验，结果表明：

1) 采用物理滤波初始化四维变分同化方法能明显提高模式在临近时段的降水预报能力。在模式起报的第1 小时即可出现与观测接近的降水。而12 min同化窗的四维变分同化需要的计算量可在13～25 min完成，使用更多计算资源可缩短完成时间，这对实际预报是可接受的。

2) 同化地面观测、探空观测和雷达观测径向风速度改善初始场中的动力场，进而改进模式降水预报。雷达资料的同化对预报改进更明显。在临近预报阶段，积云对流参数化(次网格降水)对提高降水预报具有较明显的贡献，这也提示，由于积云对流参数化的诊断特性，能在较短时间响应初始场的改变，可能是提高数值模式在临近阶段降水预报能力的一条有效途径。

3) 采用物理滤波初始化四维变分同化技术能明显改进模式对于降水临近预报能力，特别是大量级降水的预报能力。物理滤波初始化四维变分同化技术能将批量试验6 h累积降水大于25.0 mm量级的ETS评分由0.125提高到0.190，大于60.0 mm量级的ETS评分由0.016提高到0.081。

值得讨论的是，如果同化试验以改进热力场为主，如同化雷达反射率因子，是否可以通过格点可分辨降水(云微物理过程)提高临近预报效果，如果同时改进热力场和动力场，能否提高数值模式在短时临近阶段的降水预报能力。由于目前比较有效的同化雷达反射率因子的技术是云分析技术，而本研究重点是讨论物理滤波初始化四维变分同化的应用，因此相关的研究有待以后开展。