大气科学  2016, Vol. 40 Issue (6): 1127-1142   PDF    
EnSRF雷达资料同化在一次飑线过程中的应用研究
高士博1,2 , 闵锦忠1,2 , 黄丹莲1,2     
1 南京信息工程大学气象灾害预报预警与评估协同创新中心, 南京 210044
2 南京信息工程大学气象灾害教育部重点实验室, 南京 210044
摘要: 本文利用包含复杂冰相微物理过程的WRF(Weather Research and Forecasting)模式,针对2007年4月23日发生在我国华南地区的一次典型飑线天气过程,分别进行了确定性预报和集合预报试验,发现确定性预报能大致捕捉到飑线系统的发生发展过程,但对飑线后部的层云区模拟效果较差。集合预报能够有效地减少模式的不确定性,大部分集合成员对飑线的模拟效果优于确定性预报。进一步将集合预报得到的40个成员作为背景场,采用EnSRF(Ensemble Square Root Filter)同化多普勒天气雷达资料,并将分析得到的集合作为初始场进行集合预报,通过与未同化雷达资料的集合对比,考察了EnSRF同化多部雷达资料对飑线系统的影响。结果表明:EnSRF雷达资料同化增加了模式初始场的中小尺度信息,大部分集合成员的分析场能够较准确地再现飑线的热力场、动力场和微物理场的细致特征,并且模拟出飑线后部的层云结构。通过对EnSRF分析的集合进行模拟发现,大部分集合成员较未同化雷达资料时模拟效果有明显改善。同化后的集合预报ETS(Equitable Threat Score)评分最高,其次是未同化的集合预报,确定性预报的最低。
关键词 集合均方根滤波      资料同化      确定性预报      集合预报      飑线     
The Simulation of a Squall Line with Doppler Radar Data Assimilation Using the EnSRF Method
GAO Shibo1,2, MIN Jinzhong1,2, HUANG Danlian1,2     
1 Collaborative Innovation Center on Forecast and Evaluation of Meteorological Disasters, Nanjing University of Information Science & Technology, Nanjing 210044
2 Key Laboratory of Meteorological Disaster of Ministry of Education, Nanjing University of Information Science & Technology, Nanjing 210044
Abstract: An ensemble forecast and a deterministic forecast of a squall line that occurred in southern China on 23 April 2007 have been conducted using the Weather Research and Forecasting (WRF) model with microphysical schemes that include complex ice and snow processes.It is found that the deterministic forecast can capture the main characteristics of the squall line, but the simulated squall line is inaccurate, especially in the back stratus cloud region.The ensemble forecast technique can reduce the uncertainty in the model simulation and the majority of the members in the ensemble show a better performance than the deterministic forecast.The analysis members, which are obtained from radar data assimilation using the EnSRF (Ensemble Square Root Filter) method with outputs of the 40 members in the ensemble experiment as backgrounds, are used to provide initial conditions for the ensemble forecast.Differences in results among the ensemble members with and without radar data assimilation reflect the impact of EnSRF radar data assimilation on the simulation of the squall line.The analysis members with radar data assimilation provide more mesoscale and microscale information of the convective cells in the squall line system.Most members can capture the thermal-dynamical structure of the squall line system and successfully simulate the suqall line in the back stratus cloud region.Analysis of the simulations in the ensemble forecast with radar data assimilation indicates that most members perform better than that without radar data assimilation.The ETS (Equitable Threat Score) of the ensemble forecast with radar data assimilation is higher than that without radar data assimilation, and the ETS of the deterministic forecast is lower than that of the ensemble forecast.
Key words: EnSRF method      Data assimilation      Deterministic forecast      Ensemble forecast      Squall line     
1 引言

飑线是由多个活跃雷暴单体排列成线状的中尺度对流系统,常伴随短时强暴雨、大风、冰雹和龙卷等严重的灾害性天气现象。由于飑线水平尺度较小,持续时间短,发生发展迅速,在时空上很难被常规气象观测网所识别,而新一代多普勒天气雷达网的建立,为研究飑线类强对流天气的精细结构提供了有利条件。

多普勒天气雷达观测的主要是径向速度和回波强度,不能直接提供模式变量信息。目前多采用变分资料同化方法,其中3DVar (Three Dimensional Variational Data Assimilation)因计算量小而被广泛使用,如Gao et al.(1999)在风暴个例中应用ARPS (Advanced Regional Prediction System)模式和3DVar系统同化雷达观测资料,得到了风暴内部清晰的上升和下沉气流的结构。Hu et al.(2006a, 2006b)等利用ARPS的云分析和3DVar同化系统,分别对三次龙卷过程进行雷达资料同化,取得了较好的同化效果。Xiao and Sun (2007)研究发现,同化雷达径向风和反射率资料,能够为飑线预报提供精细的初始场,进而提高定量降水的预报效果。但3DVar的背景误差统计量一般认为是均匀、各向同性、随时间少变的,不能很好地反映预报误差的统计特征。为此,科学家致力于发展更高级的同化技术,4DVar (Four Dimensional Variational Data Assimilation)可以给出同化时间窗内模式背景场和雷达观测的最佳拟合,一定程度上克服了3DVar的缺点。Sun et al.(1991)首次采用4DVar开展多普勒雷达资料反演工作,随后Sun and Crook (1994)用此方法,对美国的一次阵风锋过程进行了研究,结果表明此方法对阵风锋风场结构具有较好的反演能力。孙虎林等(2011)利用4DVar方法同化高分辨率地面观测资料以及多普勒天气雷达资料,分析了一次强飑线天气过程的形成背景和成熟阶段特征。但目前4DVar同化需要模式反向积分以及设计比较准确的切线性和伴随模式,而数值预报模式的切线性和伴随模式一般都是高度复杂的,且初始的背景误差还必须先验地指定(Fisher and Andersson, 2001),因此极大限制了其在研究和业务领域的应用。

集合卡尔曼滤波(Ensemble Kalman Filter,简称EnKF)利用一组短期集合预报扰动估计背景误差协方差,进行最优化分析,而后实现对误差协方差的更新。与变分同化方法相比,EnKF能够方便地得到具有流依赖性质的背景误差协方差,并且算法简单,不需要编写切线性和伴随模式,因此广泛应用于各种天气尺度资料同化中(Hamil and Synder, 2000; Keppenne, 2000; Anderson, 2001; Mitchell et al., 2002)。近年来,EnKF也逐步应用到对流尺度雷达资料同化中,如Snyder and Zhang (2003)通过同化模拟雷达资料,初步验证了EnKF在风暴尺度中的应用效果。随后其他学者也在这方面进行了一系列的试验,进一步验证了EnKF的可利用性和有效性(Dowell et al., 2004; Zhang et al., 2004; Tong and Xue, 2005, 2008; 许小永等, 2006; Xue et al., 2006; Aksoy et al., 2009; 兰伟仁等, 2010a, 2010b; 闵锦忠等, 2011)。Whitaker and Hamill (2002)提出了集合均方根滤波(Ensemble Square Root Filter,简称EnSRF)方法,这种方法不仅继承了EnKF的优点,而且能够避免EnKF方法中由于扰动观测所造成的采样误差。近期,秦琰琰等(2012)应用EnSRF方法,针对发生在山东省的一次飑线天气过程,开展基于多普勒雷达资料的同化应用试验,取得了较好的效果。

值得注意的是,EnKF与集合预报有着非常密切的关系,它可以为集合预报提供一组具有统计意义最优的初值。众所周知,确定性预报只能代表大气真实状态的一种可能,而集合预报通过一定的数学方法,获得在合理误差范围内具有某种概率密度函数分布特征的集合成员,其中每个集合成员都能代表大气的真实情况。相比确定性预报,集合预报能够有效地减少数值预报的不确定性。近期,Harnish and Keil (2015)基于EnKF对一次强对流过程进行同化分析,并利用得到的集合成员进行集合预报,取得了较好的效果。

由此可见,EnKF在对流尺度上有较好的同化效果,但其在飑线上的应用在国内尚不多见,目前大部分研究使用的是3DVar或云分析同化方法。另外,在飑线系统的模拟上,中尺度集合预报技术相对于确定性预报技术具有一定优势。本文综合利用EnSRF同化方法和集合预报技术,对2007年4月23日发生在我国华南地区的一次典型飑线天气过程进行了同化和模拟,通过与未进行雷达资料同化的集合成员对比,考察了EnSRF雷达资料同化对WRF模式所模拟的飑线回波结构、热力场、动力场和微物理量场的影响。

2 资料和方法 2.1 资料

模式初始场为NCEP (National Centers for Environmental Prediction)提供的1°×1°再分析资料,侧边界条件每6小时更新一次,本文使用中国新一代多普勒天气雷达资料,包括厦门、福州、建阳、广州、深圳、汕头、阳江、韶关以及桂林等9部S波段天气雷达资料(见图 1)。9部雷达均采用VCP21(Volume Coverage Pattern,scan strategy 2,version 1)观测模式进行连续体扫,仰角在0.5°~19.5°之间,体扫时间间隔为6 min。雷达资料径向风分辨率为250 m,反射率因子分辨率为1 km。

2.2 集合均方根滤波方法

根据Evensen (1994)的设计,集合卡尔曼滤波(EnKF)同化循环的计算分为集合预报和分析更新两部分,具体公式如下:

$\mathit{\boldsymbol{X}}_i^{\rm{b}}(t) = \mathit{\boldsymbol{MX}}_i^{\rm{a}}(t - 1),$ (1)
$\mathit{\boldsymbol{X}}_i^{\rm{a}} = \mathit{\boldsymbol{X}}_i^{\rm{b}} + \mathit{\boldsymbol{K}}(y_j^{\rm{o}} - \mathit{\boldsymbol{HX}}_i^{\rm{b}}),$ (2)
$\mathit{\boldsymbol{K}} = {\mathit{\boldsymbol{P}}^{\rm{b}}}{\mathit{\boldsymbol{H}}^{\rm{T}}}{(\mathit{\boldsymbol{H}}{\mathit{\boldsymbol{P}}^{\rm{b}}}\mathit{\boldsymbol{H}} + \mathit{\boldsymbol{R}})^{ - 1}},$ (3)
$\mathit{\boldsymbol{H}}{\mathit{\boldsymbol{P}}^{\rm{b}}}{\mathit{\boldsymbol{H}}^{\rm{T}}} \approx [(\mathit{\boldsymbol{HX}}_i^{\rm{b}} - \overline {\mathit{\boldsymbol{H}}{\mathit{\boldsymbol{X}}^{\rm{b}}}} ){(\mathit{\boldsymbol{HX}}_i^{\rm{b}} - \overline {\mathit{\boldsymbol{H}}{\mathit{\boldsymbol{X}}^{\rm{b}}}} )^{\rm{T}}}],$ (4)

其中,公式(1)为集合预报方程,X表示模式预报的状态变量,M表示预报模式,上标a表示分析估计,上标b表示预报估计,下标i表示第i个成员,t表示当前时刻;公式(2)、(3)为更新方程,K为卡尔曼增益矩阵,H为线性观测算子,yjo表示在原始时刻观测上加入扰动后得到的观测,j表示第j个扰动观测,P表示预报误差协方差矩阵,R表示观测误差协方差矩阵,上标T表示矩阵的转置;公式(4)中的HPbHT表示观测变量的背景误差协方差,上划线表示集合成员的平均值。

根据Whitaker and Hamill (2002)的研究,EnKF在对观测进行扰动时,会引入额外的采样误差,导致分析误差协方差被低估,引起滤波发散。为此,他们提出了集合均方根滤波(EnSRF)方法,通过在EnSRF方法中引入了一个小参数α,令$\tilde K = \alpha K$来调整卡尔曼增益矩阵K。在单一观测条件下$\alpha = {[1 + \sqrt {\mathit{\boldsymbol{R}}{{(\mathit{\boldsymbol{H}}{P^{\rm{b}}}{\mathit{\boldsymbol{H}}^{\rm{T}}} + \mathit{\boldsymbol{R}})}^{ - 1}}} ]^{ - 1}}$,EnSRF的更新方程变为

$\overline {{\mathit{\boldsymbol{X}}^{\rm{a}}}} = \overline {{\mathit{\boldsymbol{X}}^{\rm{b}}}} + \mathit{\boldsymbol{K}}\left( {{y^{\rm{o}}} - \overline {\mathit{\boldsymbol{H}}{\mathit{\boldsymbol{X}}^{\rm{b}}}} } \right),$ (5)
$\mathit{\boldsymbol{X'}}_i^{\rm{a}} = \mathit{\boldsymbol{X'}}_i^{\rm{b}} - \alpha \mathit{\boldsymbol{KHX'}}_i^{\rm{b}},$ (6)

其中,X代表集合平均,${\mathit{\boldsymbol{X'}}_i}$代表集合的扰动。通过公式(6)不断更新,避免了将扰动观测带来的采样误差引入到EnKF同化中。

2.3 雷达资料质量控制及观测算子

在进行多普勒雷达资料同化之前,利用NCAR (National Center for Atmospheric Research)开发的SOLOⅡ软件对雷达资料进行质量控制,包括去除噪声点、地物回波、二次回波、速度退模糊(Oye et al., 1995)和去距离折叠。考虑到水平分辨率较高,只在水平方向上将雷达资料插值到模式网格上。为了减少计算成本,在插值过程中进行了雷达观测资料的稀疏化。

在EnSRF方法中,通过统计集合成员之间的协方差矩阵来建立观测变量和模式变量之间的相关,可以直接同化雷达径向风和反射率因子,径向风和反射率因子的观测算子同Tong and Xue (2005)

3 试验设计

应用非静力中尺度数值模式WRF3.5对飑线过程进行模拟试验,模式水平分辨率为3 km,网格数为450×400×50。试验以(24°N,114°E)为中心,垂直方向取50个不等距分层,模式顶层取50 hPa。模拟区域主要覆盖中国广东、广西、福建以及江西、湖南、浙江大部分地区,地形分辨率为30 s (见图 1)。微物理过程为WSM6冰相方案(Hong et al., 2004),长波辐射为RRTM方案,短波辐射为Dudhia方案,不使用积云对流参数化方案。

图 1 模拟区域与雷达位置[图中黑色矩形框代表模拟范围、倒三角代表雷达位置,大圆圈代表半径为230 km雷达观测的范围,其中9部雷达分别为厦门(XMRD)、福州(FZRD)、建阳(JYRD)、广州(GZRD)、深圳(SZRD)、汕头(STRD)、阳江(YJRD)、韶关(SGRD)以及桂林(GLRD)雷达] Figure 1 The model domain and radar locations. Black rectangle: area of simulation; circles: the coverage of 230 km radar observations; inverted triangles: locations of the 9 radar stations at Xiamen (XMRD), Fuzhou (FZRD), Jianyang (JYRD), Guangzhou (GZRD), Shenzhen (SZRD), Shantou (STRD), Yangjiang (YJRD), Shaoguan (SGRD), and Guilin (GLRD)

本文设计了三组试验,分别用来检验集合预报和雷达资料同化对飑线模拟的改进效果:

第一组为控制试验,只用NCEP 1°×1°再分析资料进行模拟,即确定性预报试验。模式预报时间为2007年4月23日12:00(协调世界时,下同)至4月24日02:00。

第二组为集合预报试验,首先对23日12:00的NCEP资料进行8小时“spin up”后,在23日20:00加入初始扰动产生40个集合成员。扰动方法为随机背景误差协方差扰动方法,该方法利用GFS (Globe Forecast System)资料统计的背景误差协方差产生随机扰动,经垂直EOF (Empirical Orthogonal Function)分析、递归滤波和平衡方程约束后变换回模式空间,所产生的各个变量的扰动之间具有一定的关系和平衡约束(Skamarock et al., 2008),因此对该扰动方法产生的集合成员进行集合预报,效果可能会优于确定性预报。其中风场的水平分量uv标准差约为2 m s-1,位温pt的标准差约为1 K,水汽混合比qv标准差约为0.5 g kg-1。初始场没有同化任何资料,模式预报时间为2007年4月23日20:00至4月24日02:00。

第三组为同化试验,利用ARPS EnSRF同化系统的WRFAPI功能来开展实际雷达资料的同化试验。WRFAPI模块通过增加WRF模式变量的接口,利用ARPS EnSRF同化系统将观测同化进WRF模式的初始场,该模块可同化WRF模式的uvwptphqvqcqrqiqs以及qg变量(uvw分别为扰动风场的纬向分量、经向分量、垂直分量,pt为扰动位温,ph是位势高度变量,qh为冰雹混合比,qv为水汽混合比,qc为云水混合比,qr为雨水混合比,qi为冰晶混合比,qs为雪混合比,qg为霰混合比)。分别在23日23:00和24日00:00,通过EnSRF同化方法将9部多普勒雷达资料同化进背景场。协方差局地化为五阶距离相关函数方案(Gaspari and Cohn, 1999),参考Putnam et al.(2014),水平局地化距离为10 km,垂直局地化距离为4 km,协方差膨胀参数为1.25。最后利用EnSRF同化得到的集合成员,进行2 h的集合预报。模式预报时间为2007年4月24日00:00至4月24日02:00。

4 飑线过程及环境背景

2007年4月23日,受中层的低压槽,地面准静止锋和夜间低空急流的影响,华南地区发生了一次典型的飑线过程。飑线系统在20:00进入广西和广东省交界处,随后以17 m s-1向东南方向移动,24日02:00飑线强度达到最大,水平尺度达800 km左右,基本覆盖了整个广东省,一直到24日03:00逐渐消散,维持了约7 h。经统计,广东省累计降水超过50 mm,部分县市雷雨大风为8~9级,强降水超过100 mm。局部地区发生冰雹灾害,最大风速可达24 m s-1

从23日12:00的850 hPa位势高度场可以看出(图 2a),在(25°N,105°E)附近存在一个低压中心,从12:00到18:00,该低压中心与来自中纬度的一个短波槽合并。等位温线的狭窄密集带表示温度梯度较大的区域,一般出现在锋面附近。从位温场的分布特征可以看出,在25°N附近存在一个准静止锋即梅雨锋,大致呈东西走势,从中国西部一直延伸到日本。18:00在广西省西北部(图 2b),低层风速的北风分量突然加强,来自北方的冷空气推动准静止锋西部南下。北方干冷空气与南方的暖湿空气的辐合不断加强,另外从可降水量可以看出,我国广西、广东地区超过60 kg m-2,并且有一只西南气流不断地将水汽向两广地区输送,这都为飑线的发生提供了有利的环流背景形势和水汽条件。

图 2 2007年4月(a)23日12:00(协调世界时,下同)和(b)23日18:00基于NCEP再分析资料的850 hPa环境场变量。阴影为可降水量(单位:kg m−2),黑实线为等高线(单位:dagpm),红虚线为等位温线(单位:K),矢量为风场(单位:m s−1 Figure 2 Environmental variables at 850 hPa extracted from the NCEP reanalysis data, including the column-integrated precipitable water (shaded, units:kg m−2), geopotential height (black solid lines, units: dagpm), potential temperature (red dashed lines, units: K), and wind vectors (blue vectors, units: m s−1) at (a) 1200 UTC and (b) 1800 UTC, 23 April 2007
5 试验结果 5.1 控制试验模拟能力评估

从实况雷达组合反射率因子(图 3)可见,23日22:00,飑线移近广东和广西交界处,形成一条东北-西南向呈线状的强对流回波区,东西方向的长度为300 km,南北方向的宽度为30 km,强回波中心呈弓形,最强回波超过50 dBZ,此时飑线系统正处于发展阶段。在系统移动的后方有大范围的层云形成,这是因为在对流中心内部存在一支倾斜上升气流,该气流从飑线系统前部低层进入,逐渐上升从飑线后部顶端流出,期间会把飑线前部线状对流带内小的降水粒子向飑线系统后部抛洒,从而形成较为宽广的层云带。层云内有两块呈东北-西南走向的次强回波,回波强度均超过30 dBZ(图 3a)。1小时后,强对流回波和次强回波之间出现一条弱回波过渡带,且回波向东移动,强度加强,福建北部相对较弱的回波向西南方向移动,同时在移动方向的前方不断有新的对流单体生成(图 3b)。24日00:00,飑线已经移动到广东省中部,广东和广西原来分裂的两个强对流区域合并,形成一条长度约为500 km的新线性对流系统,同时弱回波过渡区更加明显,宽度约为35 km。24日01:00,层云内的回波强度进一步加强,最强回波超过40 dBZ。24日02:00,飑线发展到全盛期,强对流回波区呈现显著的“V”型结构,这表明近地面风速较大,飑线前段对流单体的移动速度很快。另外,层云区的回波范围和强度继续增大。

图 3 雷达观测组合反射率因子(单位:dBZ)2007年4月:(a)23日22:00;(b)23日23:00;(c)24日00:00;(d)24日01:00;(e)24日02:00 Figure 3 Observed composite radar reflectivity (units: dBZ) in the squall line at (a) 2200 UTC 23, (b) 2300 UTC 23, (c) 0000 UTC 24, (d) 0100 UTC 24, and (e) 0200 UTC 24 April 2007

图 4为控制试验模拟的雷达组合反射率因子,与图 3对比可知,确定性预报基本上模拟出了飑线前部的线状强对流回波区,但模拟的范围明显偏大,强度偏强,“V”型结构也不明显,未模拟出飑线的层云区域,模拟的飑线位置偏南,移动速度偏慢。23日22:00和23:00,强对流回波区在南北方向加宽,宽度约为50 km,弓形形状不明显,福建省北部的对流区域也未模拟出来。从23日22:00到24日02:00,福建的西南部的回波开始逐渐增强,一直到24日02:00,最强回波强度超过50 dBZ

图 4图 3,但为确定性预报试验模拟结果 Figure 4 Same as Fig. 3, but for the deterministic forecast experiment
5.2 集合预报模拟结果

上述分析表明,确定性预报技术能大概模拟出强对流回波区的主要特征,但对飑线后部的层云区、弱回波过渡带以及飑线位置的模拟与观测有较大的偏差。研究表明:对于风暴(Aberson, 2001)、飓风(Stensrud and Weiss, 2002)及局地强降水等中尺度极端天气事件,中尺度集合预报能够提供更多的概率预报信息,相比确定性预报表现出一定的优势。因此本文通过随机背景误差协方差扰动方法产生40个集合成员,进一步探讨集合预报的模拟效果。集合成员是利用WRFDA的三维变分程序,对背景场添加随机的背景误差协方差而产生的。

参考Parker and Johnson (2000)Meng and Zhang (2012)的方法,在飑线中,如果模拟的线状对流中心(回波强度超过40 dBZ)与观测的距离不超过200 km,则认为模拟出了飑线的主要特征。经分析基本上所有成员都达到了上述标准,并且计算各个集合成员的组合回波与雷达实况组合回波的均方根误差,其中以第2、3、12、13、16、17、19、21、22、36以及37个集合的均方根误差最低(图略),效果最佳,相比于确定性预报,飑线在西部的回波增强,“V”型结构更明显,更接近于观测场。考虑到集合成员的数目较多,以上述均方根误差最低的11个集合成员为例。图 5为这11个集合成员24日00:00雷达组合反射率因子的分布。与观测(图 5a)对比发现,集合预报模拟的线状强对流回波区在南北方向尺度同样偏大。在确定性预报中,广东省西部的回波为分散的点状,而集合成员模拟的对流回波在此处强度增强,且回波为团状,对流中心的弓形形状更明显,与观测场更接近。不同集合成员的模拟结果也存在一定差异,其中以第12个和第37个集合成员的均方根误差最低,这说明飑线系统对初始场是比较敏感的。

图 5 集合预报试验模拟的2007年4月24日00:00雷达组合反射率因子(单位:dBZ):(a)观测(OBS);(b-l)分别是第2、3、12、13、16、17、19、21、22、36和37个集合成员 Figure 5 The composite radar reflectivity (units: dBZ) in the squall line from observations and simulated by the ensemble forecast experiment at 0000 UTC 24 April 2007: (a) Observations (OBS); (b-l) simulations of member 2, member 3, member 12, member 13, member 16, member 17, member 19, member 21, member 22, member 36, and member 37
5.3 同化结果

5.2节以雷达观测为标准,分别对比了确定性预报和集合预报试验的组合回波分布特征,结果表明集合预报试验能够更加细致地描述飑线的中尺度结构。下面我们分别在23日23:00和24日00:00,利用EnSRF方法同化9部多普勒雷达资料。为了和同化雷达资料前的模拟结果对比,这里同样以24日00:00(同化结束时刻)的11个集合成员为例。

由雷达观测实况(图 6a)可知,在24日00:00,飑线已经发展为一个典型的成熟系统,飑线的三个组成部分,包括前部强对流回波区、中部对流向层云的过渡区和后部的层云区都清晰可见。其中强对流回波区为狭长的线状结构,回波整体东西走向,南北方向尺度较小,“V”型结构明显,最强回波强度能达到50 dBZ以上。层云区的亮带结构明显,过渡区的回波强度较弱,仅25 dBZ左右。对比图 6图 5可以看出,集合预报试验和同化试验的回波在空间分布上存在一定差异。同化前飑线前部的强对流回波区在南北方向的范围偏大,强度略强,且没有模拟出层云回波(图 5)。经过同化后,各个集合基本上呈现出了飑线成熟阶段主要特征,其中强对流回波区的范围和强度有所减少,线状结构更加明显,飑线后部的层云回波和过渡区也有所体现。另外,同化改善了福建处的回波结构,回波强度大于20 dBZ的范围明显增大(图 6)。可见,雷达资料同化有效地改善了飑线系统的回波结构(尤其是层云区)和强度,使背景场较同化前更接近实况,因此在数值模式中同化雷达资料对飑线的结构影响很大。

图 6图 5,但为EnSRF同化试验结果 Figure 6 Same as Fig. 5, but for EnSRF final analysis results

为了进一步考察雷达资料同化对飑线系统结构的影响,以第13个集合成员为例,图 7分别对确定性预报、集合预报(同化前)和同化试验的反射率因子过(24°N,115°E),沿纬向作垂直剖面。由图 7可以看出,三组试验均给出了飑线前缘的强对流回波区。相比确定性预报和集合预报试验,同化试验线状对流后部的层云回波范围有所增加。强对流区后方的层云水平范围达到200 km左右,平均伸展高度约为10 km,最强回波接近45 dBZ,这是因为层云上方固态水凝物降落到约4 km时表面开始融化,导致粒子群的雷达反射率迅速增大而产生的回波大值区。

图 7 三组试验2007年4月24日00:00雷达反射率因子的垂直剖面(单位:dBZ):(a)确定性预报;(b)未同化雷达资料的第13个集合成员;(c)同化雷达资料的第13个集合成员的分析场 Figure 7 Vertical cross sections of radar reflectivity in the squall line at 0000 UTC 24 April 2007 (units: dBZ): (a) Deterministic forecast, (b) member 13 forecast without radar data assimilation, (c) member 13 analysis with radar data assimilation

冷池对飑线系统的发展和维持有着非常重要的作用。由图 8a-c可以发现,三组试验在地表的位温场均存在低值中心,即低层形成明显的冷池结构。冷池中心最低扰动位温能达到-7 K,这与飑线的强对流回波区相对应。飑线具有明显的冷池结构是因为降水导致冷空气不断下沉扩散,从而使得近地面层扰动位温下降。对比三组试验发现,确定性预报试验模拟的冷池范围最小,强度最弱,集合预报试验次之,同化试验的冷池范围最大,强度最强,具有明显的冷空气丘结构(王焱等,2008)。在冷池的西侧,同化试验的风速明显大于确定性预报和集合预报试验,气流从冷空气堆处向外辐散流出,并与前部暖区的气流相汇合。为了更清楚地对比同化前后扰动位温的差异,图 8d为同化试验与集合预报试验扰动位温的差值场,发现同化后地面冷池更明显,扰动位温明显降低,最多约降低3 K。

图 8 三组试验2007年4月24日00:00,地表扰动位温(阴影,单位:K)和水平速度(矢量,单位:m s-1):(a)确定性预报;(b)未同化雷达资料的第13个集合成员;(c)同化雷达资料的第13个集合成员的分析场;(d)第13个集合成员同化试验与集合预报试验地表面扰动位温的差 Figure 8 Surface potential temperature perturbations (shaded, units: K) and horizontal winds (vector, units: m s-1) in the squall line at 0000 UTC 24 April 2007: (a) Deterministic forecast, (b) member 13 forecast without radar data assimilation, (c) member 13 analysis with radar data assimilation, (d) the difference in surface potential temperature perturbation between simulations of member 13 with and without radar data assimilation

由于冷池本身是由高密度的冷空气组成,与周围相比是一高气压区(雷暴高压),因此强冷池的形成将会在近地面形成强雷暴高压。由图 9可以看出,三组试验在广东省均存在雷暴高压中心,且在其周围有很强的向四周辐散的风,前沿的辐散风速达到10 m s-1以上。三组试验中,同化试验的雷暴高压气压最高,范围最大,雷暴高压前缘辐散的西北气流和东南暖湿气流辐合也最明显。从同化前后地表气压的差值场(图 9d)可知,在飑线的前部,气压有所降低,最多约下降3 hPa,在飑线的后部,地面气压差值存在高值区,最多约升高4 hPa。

图 9图 8,但为地表面气压(阴影,单位:hPa)和水平速度 Figure 9 Same as Fig. 8, but for surface pressure (shaded, units: hPa) and horizontal winds

雷达资料同化对水凝物总量(包括云水混合比qc、雨水混合比qr、冰晶混合比qi、雪混合比qs以及雹霰混合比qg)分布也有一定的改善,由图 10a-c可知,三组试验在广东省附近均有一条呈带状的水凝物高值区,这与飑线的线状强对流回波区是相对应。由图 4c图 5e可知,确定性预报和集合预报试验给出的飑线线状强对流区域强度偏强,且宽度偏大,因此该处的水凝物的范围也偏大,强度偏强。同化试验的水凝物高值中心,与观测的线状强对流中心(红虚线)的位置最接近,其它两组试验位置偏南,同化试验在水凝物高值带后方,也同化出水凝物次高值区,这与观测的层云区是相吻合的。对比同化前后的差异(图 10d)发现,同化后飑线线状对流带内的水凝物有所下降,最多约下降2 g kg-1,飑线后部层云区位置的水凝物明显增加,最多约升高2 g kg-1

图 10图 8,但为2 km水凝物总量(阴影,单位:g kg-1)和水平速度 Figure 10 Same as Fig. 8, but for water mixing ratio (shaded, units: g kg-1) and horizontal winds (vector, units: m s-1) in the squall line at 2 km
5.4 预报结果

以上分析表明,经过雷达资料同化后,大部分集合成员的飑线回波结构明显地改善,且能够较准确地分析出飑线的热力场、动力场和微物理量场的细致特征。为了进一步检验其同化效果,对同化雷达资料的集合成员进行了两个小时的预报,模式预报时间为2007年4月24日00:00至4月24日02:00。考虑到集合成员的数目较多,仍以上述的11个集合成员为例。

图 11a可知,在24日02:00,大于5 dBZ的雷达回波水平尺度超过800 km,并且在飑线系统内部存在一条强对流回波带,东西方向的长度为500 km,南北宽度为25 km,回波中心强度达到50 dBZ,“V”形结构比较明显,这说明飑线系统此时非常强盛。相比于图 5a,强对流回波带后方的层云回波强度增强,覆盖范围增大。广东省和福建省原来两块分离的回波合并成一条线状回波,回波强度均超过30 dBZ。对比图 11b-l可以看出,同化前各个集合成员均能模拟出强对流回波带,也大致模拟出“V”型结构,但相比于观测,强对流回波带的强度偏强,南北方向的宽度偏大,约为60 km,结构较为松散。同时可以看到,强对流回波带向东南方向倾斜,一直延伸至海上,飑线后方的层云区也不明显,层云区大部分为零星的点状回波。飑线系统对初始场较为敏感,各个集合成员的模拟效果也存在一定差异,其中第13、36个集合成员的模拟效果较好,“V”型结构最明显,且飑线的线型对流中心位置有了一定程度地改善,但对后部层云结构的模拟效果仍欠佳。

图 11图 5,但为02:00 Figure 11 Same as Fig. 5, but for 0200 UTC

图 12为基于EnSRF雷达资料同化试验的分析场,上述11个集合成员所模拟的雷达组合反射率因子。由图可以看出,加入雷达资料后,对飑线的预报效果有一定改进,线状强对流回波区的强度有所减弱,且南北方向宽度减少,这都更接近于实况。值得注意的是,飑线线状强对流回波带后部的层云区也有所体现,虽然较观测有一定差距,但是效果明显优于未同化雷达资料的集合预报试验。

图 12 EnSRF雷达资料同化试验模拟的2007年4月24日02:00雷达组合反射率因子(单位:dBZ):(a)观测;(b-l)分别是第2、3、12、13、16、17、19、21、22、36和37个集合成员 Figure 12 The composite radar reflectivity (unit: dBZ) in the squall line at 0200 UTC 24 April 2007 from observations and from simulations based on the EnSRF final analysis field: (a) Observations, (b-l) results of member 2, member 3, member 12, member 13, member 16, member 17, member 19, member 21, member 22, member 36, member 37

同样以第13个集合成员为例,图 13为其反射率因子过(24°N,115°E)沿着纬向作的垂直剖面。由三组试验可以看出,在飑线发展强盛期,强对流回波区位于飑线的前缘,确定性预报试验模拟的线状强对流区域宽度偏大,为80 km左右,同化试验模拟的线状强对流区域尺度约为30 km左右,与观测组合回波显示的线状对流区域宽度更加吻合。对于强对流回波区后方的层云区,同化试验模拟的层云水平范围最大,为150 km左右,集合预报试验模拟的层云区水平范围次之,确定性预报试验则基本没有模拟出层云结构。

图 13 三组试验2007年4月24日02:00模拟的雷达反射率因子的垂直剖面(单位:dBZ):(a)确定性预报;(b)未同化雷达资料的第13个集合成员;(c)同化雷达资料的第13个集合成员的预报场 Figure 13 The simulated cross sections of reflectivity of the squall line at 02:00 UTC on April 24, 2007(unit: dBZ): (a) Deterministic; (b) member 13 forecast with no radar data; (c) member 13 forecast with radar data

为了进一步定量地评估同化及模拟效果,分别计算三组试验从23日22:00到24日02:00,40个集合成员关于组合反射率因子,阀值为20 dBZ的ETS (Equitable Threat Score)评分(图略)。可以发现,大部分集合成员,同化雷达资料后ETS评分最高,其次是未同化雷达资料时,确定性预报最小。图 14给出了上述11个集合成员雷达组合反射率因子(20 dBZ)的ETS评分,由图可以看出,未同化雷达资料时,集合成员的ETS评分均高于确定性预报试验,其中以第13个和第21个集合成员最明显,ETS评分提高了0.06以上。经过同化后,集合成员的ETS评分均高于同化前的集合成员,其中第3、13、17、19和37个集合成员的ETS评分比同化前提高0.1以上。第2个集合成员在三组试验的ETS评分差异不大。可见,雷达资料同化能够有效地提高模式对飑线系统的模拟水平。

图 14 2007年4月23日22:00至24日02:00,确定性预报试验、集合预报试验和同化试验的第2、3、12、13、16、17、19、21、22、36和37个集合成员的雷达组合反射率因子(20 dBZ)的ETS评分 Figure 14 The ETS scores for the composite reflectivity (20 dBZ) simulated by member 2, member 3, member 12, member 13, member 16, member 17, member 19, member 21, member 22, member 36 and member 37 for the period from 2200 UTC 23 April 2007 to 0200 UTC 24 April 2007
6 结论

本文采用中尺度WRF模式,针对2007年4月23日发生在华南地区的一次典型飑线过程,探讨了中尺度集合预报技术对飑线的模拟能力,并在此基础上,基于EnSRF方法同化雷达资料,分别对40个集合成员进行了同化及模拟试验,验证了EnSRF同化多部雷达资料的效果。主要结论如下:

(1)相比确定性预报,集合预报具有更好的模拟能力。确定性预报能大致捕捉到飑线系统的发生发展过程,但模拟的飑线强对流区域范围偏大,“V”型结构不明显,飑线系统后部的层云结构也未能模拟出来。集合预报通过对多个集合成员分别积分,能最大限度地模拟出飑线的可能状态。大部分集合成员均模拟出飑线的主要特征,其中有11个集合成员均方根误差最低,模拟效果最佳。

(2) EnSRF能有效地同化实际雷达资料,增加模式初始场的中小尺度信息,使各个集合成员的分析场在飑线强对流回波区的范围和强度上更加准确,且同化试验分析出了飑线后部的层云回波和过渡区。相比确定性预报和集合预报,同化后得到的冷池面积增大,强度变强,“V”型结构更明显,飑线后部雷暴高压强度增大,水凝物高值中心的位置更加准确,且在层云区存在次高值区。

(3)利用集合预报技术,通过对同化后的集合成员进行模拟发现,大部分集合成员较未同化雷达资料时模拟效果有明显改善,模拟的雷达回波更接近观测,强对流回波区的线状回波和“V”型结构明显,并且较好地模拟出飑线后部的层云区。同化后的集合预报ETS评分最高,其次是未同化的集合预报,确定性预报最低。

需要指出的是,本文只针对一次飑线过程进行了初步研究,未来还需要结合更多的飑线个例进行检验。另外,随着预报时间的增加,2小时后,同化试验的ETS评分会逐渐降低,这可能与飑线系统内部误差增长较快有关,然而目前如何提高预报时效是一个难点,也有待进一步的研究。

参考文献
[] Aberson S D. 2001. The ensemble of tropical cyclone track forecasting models in the North Atlantic Basin (1976-2000)[J]. Bull.Amer.Meteor.Soc., 82(9) : 1895–1904 DOI:10.1175/1520-0477(2001)082<0000:TEOTCT>2.3.CO; 2
[] Aksoy A, Dowell D C, Snyder C. 2009. A multicase comparative assessment of the ensemble Kalman filter for assimilation of radar observations.Part I:Storm-scale analyses[J]. Mon.Wea.Rev., 137(6) : 1805–1824 DOI:10.1175/2008MWR2691.1
[] Anderson J L. 2001. An ensemble adjustment Kalman filter for data assimilation[J]. Mon.Wea.Rev., 129(12) : 2884–2903 DOI:10.1175/1520-0493(2001)129<2884:AEAKFF>2.0.CO; 2
[] Stensrud D J, Weiss S J. 2002. Mesoscale model ensemble forecasts of the 3 May 1999 Tornado outbreak[J]. Wea.Forecasting, 17(3) : 526–543 DOI:10.1175/1520-0434(2002)017<0526:MMEFOT>2.0.CO; 2
[] Dowell D C, Zhang F Q, Wicker L J, et al. 2004. Wind and temperature retrievals in the 17 May 1981 Arcadia, Oklahoma, supercell:Ensemble Kalman filter experiments[J]. Mon.Wea.Rev., 132(8) : 1982–2005 DOI:10.1175/1520-0493(2004)132<1982:WATRIT>2.0.CO; 2
[] Evensen G. 1994. Sequential data assimilation with a nonlinear quasi-geostrophic model using Monte Carlo methods to forecast error statistics[J]. J.Geophys.Res.:Oceans (1978-2012), 99(C5) : 10143–10162 DOI:10.1029/94JC00572
[] Fisher M, Andersson E. 2001. Developments in 4D-Var and Kalman filtering[J]. ECMWFTechinical Memorandum, No.347, Berkshire:ECMWF : 1–36
[] Gao J D, Xue M, Shapiro A, et al. 1999. A variational method for the analysis for three-dimensional wind fields from two Doppler radars[J]. Mon.Wea.Rev., 127(9) : 2128–2142 DOI:10.1175/1520-0493(1999)127<2128:AVMFTA>2.0.CO; 2
[] Gaspari G, Cohn S E. 1999. Construction of correlation functions in two and three dimensions[J]. Quart.J.Roy.Meteoro.Soc., 125(554) : 723–757 DOI:10.1002/qj.49712555417
[] Hamill T M, Synder C. 2000. A hybrid ensemble Kalman filter-3D variational analysis scheme[J]. Mon.Wea.Rev., 128(8) : 2905–2919 DOI:10.1175/1520-0493(2000)128<2905:AHEKFV>2.0.CO; 2
[] Harnish F, Keil C. 2015. Initial conditions for convective-scale ensemble forecasting provided by ensemble data assimilation[J]. Mon.Wea.Rev., 143(5) : 1583–1600 DOI:10.1175/MWR-D-14-00209.1
[] Hong S Y, Dudhia J, Chen S H. 2004. A revised approach to ice microphysical processes for the bulk parameterization of clouds and precipitation[J]. Mon.Wea.Rev., 132(1) : 103–120 DOI:10.1175/1520-0493(2004)132<0103:ARATIM>2.0.CO; 2
[] Hu M, Xue M, Brewster K. 2006a. 3DVAR and cloud analysis with WSR-88D Level-Ⅱ data for the prediction of the Fort Worth, Texas, tornadic thunderstorms.Part I:Cloud analysis and its impact[J]. Mon.Wea.Rev., 134(2) : 675–698 DOI:10.1175/MWR3092.1
[] Hu M, Xue M, Gao J D, et al. 2006b. 3DVAR and cloud analysis with WSR-88D Level-Ⅱ data for the prediction of the Fort Worth, Texas, Tornadic thunderstorms.Part Ⅱ:Impact of radial velocity analysis via 3DVAR[J]. Mon.Wea.Rev., 134(2) : 699–721 DOI:10.1175/MWR3093.1
[] Keppenne C L. 2000. Data assimilation into a primitive-equation model with a parallel ensemble Kalman filter[J]. Mon.Wea.Rev., 128(6) : 1971–1981 DOI:10.1175/1520-0493(2000)128<1971:DAIAPE>2.0.CO; 2
[] 兰伟仁, 朱江, XueMing, 等. 2010a. 风暴尺度天气下利用集合卡尔曼滤波模拟多普勒雷达资料同化试验I.不考虑模式误差的情形[J]. 大气科学, 34(3) : 640–652. Lan Weiren, Zhu Jiang, Xue Ming, et al. 2010a. Storm-scale ensemble Kalman filter data assimilation experiments using simulated Doppler radar data.Part I:Perfect model tests[J]. Chinese Journal of Atmospheric Sciences (in Chinese), 34(3) : 640–652 DOI:10.3878/j.issn.1006-9895.2010.03.15
[] 兰伟仁, 朱江, XueMing, 等. 2010b. 风暴尺度天气下利用集合卡尔曼滤波模拟多普勒雷达资料同化试验Ⅱ.不考虑模式误差的情形[J]. 大气科学, 34(4) : 737–753. Lan Weiren, Zhu Jiang, Xue Ming, et al. 2010b. Storm-scale ensemble Kalman filter data assimilation experiments using simulated Doppler radar data.Part Ⅱ:Imperfect model tests[J]. Chinese Journal of Atmospheric Sciences (in Chinese), 34(4) : 737–753 DOI:10.3878/j.issn.1006-9895.2010.04.07
[] Meng Z Y, Zhang Y J. 2012. On the squall lines preceding landfalling tropical cyclones in China[J]. Mon.Wea.Rev., 140(2) : 445–470 DOI:10.1175/MWR-D-10-05080.1
[] 闵锦忠, 陈杰, 王世璋, 等. 2011. WRF-EnSRF同化系统的效果检验及其应用[J]. 气象科学, 31(2) : 135–144. Min Jinzhong, Chen Jie, Wang Shizhang, et al. 2011. The tests of WRF-EnSRF data assimilation system and its applications in storm scale events[J]. J.Meteor.Sci., 31(2) : 135–144 DOI:10.3969/j.issn.1009-0827.2011.02.003
[] Mitchell H L, Houtekamer P L, Pellerin G. 2002. Ensemble size, balance, and model-error representation in an ensemble Kalman filter[J]. Mon.Wea.Rev., 130(11) : 2791–2808 DOI:10.1175/1520-0493(2002)130<2791:ESBAME>2.0.CO; 2
[] Oye R, Mueller C, Smith S.1995.Software for radar translation, visualization, editing, and interpolation[C]//Proceedings of the 27th Conference on Radar Meteorology.Boston, USA:American Meteorological Society, 359-361.
[] Parker M D, Johnson R H. 2000. Organizational modes of midlatitude mesoscale convective systems[J]. Mon.Wea.Rev., 128(10) : 3413–3436 DOI:10.1175/1520-0493(2001)129<3413:OMOMMC>2.0.CO; 2
[] Putnam B J, Xue M, Jung Y, et al. 2014. The analysis and prediction of microphysical states and polarimetric radar variables in a mesoscale convective system using double-moment microphysics, multinetwork radar data, and the ensemble Kalman filter[J]. Mon.Wea.Rev, 142(1) : 141–162 DOI:10.1175/MWR-D-13-00042.1
[] 秦琰琰, 龚建东, 李泽椿. 2012. 集合卡尔曼滤波同化多普勒雷达资料的观测系统模拟试验[J]. 气象, 38(5) : 513–525. Qin Yanyan, Gong Jiandong, Li Zechun. 2012. Assimilation of Doppler radar observations with an ensemble Kalman filter:OSS experiments[J]. Meteor.Mon.(in Chinese), 38(5) : 513–525 DOI:10.7519/j.issn.1000-0526.2012.5.001
[] Skamarock W C, Klemp J B, Dudhia J, et al.2008.A description of the advanced research WRF version 3[R].NCAR Technical Note TN-475+STR, 113.
[] Snyder C, Zhang F Q. 2003. Assimilation of simulated Doppler radar observations with an ensemble Kalman filter[J]. Mon.Wea.Rev., 131(8) : 1663–1677 DOI:10.1175//2555.1
[] 孙虎林, 罗亚丽, 张人禾, 等. 2011. 2009年6月3~4日黄淮地区强飑线成熟阶段特征分析[J]. 大气科学, 35(1) : 105–120. Sun Hulin, Luo Yali, Zhang Renhe, et al. 2011. Analysis on the mature-stage features of the severe squall line occurring over the Yellow River and Huaihe River basins during 3-4 June 2009[J]. Chinese Journal of Atmospheric Sciences (in Chinese), 35(1) : 105–120 DOI:10.3878/j.issn.1006-9895.2011.01.09
[] Sun J Z, Crook A. 1994. Wind and thermodynamic retrieval from single-Doppler measurements of a gust front observed during Phoenix Ⅱ[J]. Mon.Wea.Rev., 122(6) : 1075–1091 DOI:10.1175/1520-0493(1994)122<1075:WATRFS>2.0.CO; 2
[] Sun J Z, Flicker D W, Lilly D K. 1991. Recovery of three-dimensional wind and temperature fields from simulated single-Doppler radar data[J]. J.Atmos.Sci., 48(6) : 876–892 DOI:10.1175/1520-0469(1991)048<0876:ROTDWA>2.0.CO; 2
[] Tong M Q, Xue M. 2005. Ensemble Kalman filter assimilation of Doppler radar data with a compressible nonhydrostatic model:OSS experiments[J]. Mon.Wea.Rev., 133(7) : 1789–1807 DOI:10.1175/MWR2898.1
[] Tong M Q, Xue M. 2008. Simultaneous estimation of microphysical parameters and atmospheric state with simulated radar data and ensemble square root Kalman filter.Part Ⅱ:Parameter estimation experiments[J]. Mon.Wea.Rev., 136(5) : 1649–1668 DOI:10.1175/2007MWR2071.1
[] 王焱, 潘益农, 潘玉洁. 2008. 一次飑线过程的数值模拟及诊断分析[J]. 南京大学学报(自然科学版), 44(6) : 583–597. Wang Yan, Pan Yinong, Pan Yujie. 2008. A case study on squall line:Numerical simulation and diagnosis analysis[J]. Nanjing Univ.(Nat.Sci.), 44(6) : 583–597 DOI:10.3321/j.issn:0469-5097.2008.06.002
[] Whitaker J S, Hamill T M. 2002. Ensemble data assimilation without perturbed observations[J]. Mon.Wea.Rev., 130(7) : 1913–1924 DOI:10.1175/1520-0493(2002)130<1913:EDAWPO>2.0.CO; 2
[] Xiao Q N, Sun J Z. 2007. Multiple-radar data assimilation and short-range quantitative precipitation forecasting of a squall line observed during IHOP_2002[J]. Mon.Wea.Rev., 135(10) : 3381–3404 DOI:10.1175/MWR3471.1
[] 许小永, 刘黎平, 郑国光. 2006. 集合卡尔曼滤波同化多普勒雷达资料的数值试验[J]. 大气科学, 30(4) : 712–728. Xu Xiaoyong, Liu Liping, Zheng Guoguang. 2006. Numerical experiment of assimilation of Doppler radar data with an ensemble Kalman filter[J]. Chinese Journal of Atmospheric Sciences (in Chinese), 30(4) : 712–728 DOI:10.3878/j.issn.1006-9895.2006.04.16
[] Xue M, Tong M J, Droegemeier K K. 2006. An OSSE framework based on the ensemble square-root Kalman filter for evaluating the impact of data from radar networks on thunderstorm analysis and forecasting[J]. J.Atmos.Oceanic Technol., 23(1) : 46–66 DOI:10.1175/JTECH1835.1
[] Zhang Fuqing, Snyder C, Sun Juanzhen. 2004. Impacts of initial estimate and observation availability on convective-scale data assimilation with an ensemble Kalman filter[J]. Mon.Wea.Rev., 132(5) : 1238–1253 DOI:10.1175/1520-0493(2004)132<1238:IOIEAO>2.0.CO; 2