雷达资料同化在可致龙卷的强暴雨 系统分析中的一些应用

傅承浩 郜吉东 王云恒, 唐杰 周长青 叶成志 庄照荣,

（1 湖南省气象台，长沙 410118；2 美国国家强风暴实验室，诺曼 73072，美国；3 美国俄克拉何马大学中尺度气象合作研究所，诺曼 OK 73072，美国；4 国家气象中心，中国气象局，北京 100081）

0 引言

因为雷达径向风和数值模式风场之间是线性关系，同化雷达径向风到对流尺度的数值模式相对来说比较容易，尤其是在变分同化的框架内。在过去30年中，变分资料同化的应用无论是三维变分还是四维变分主要目的是用来初始化数值预报模式[1-11]。但也有很多研究应用在雷达资料的反演过程中[12-20]。在雷达资料同化的研究中，NCAR的Sun等[12, 17]做了很多四维变分同化方面的工作。近年来，集合卡尔曼滤波也在雷达资料同化的应用中取得了很大的进展。本文主要是对作者近10年来在俄克拉荷马大学强风暴中心和NOAA强风暴实验室参与的三维变分研究以及在此基础上发展的集合变分同化系统和混合变分集合同化系统等作些简要的回顾和展望。

在雷达资料研究中，Gao等[18]较早提出了一个利用两个以上雷达径向风场的分析方案。这个方案允许较灵活地同时应用雷达资料和其他与风场有关的资料，例如探空资料、地面站资料和VAD廓线等。这个方案还包括了两个弱约束：质量连续性方程和平滑算子。无论是在理想试验还是在实时观测数据的试验中，这个方案都得到了较好的结果。但是这个方案也有缺陷，特别是在确定约束项的最优权重过程中存在一些困难。因此，发展了一个更标准化的三维变分同化系统[20]。在这个新系统中，递归滤波方法[21]被用于模拟背景场误差距阵。利用递归滤波可以很有效地把每个雷达资料同化到观测位置附近的模式格点上[22]。同时能降低变分系统的条件数，从而使得数值求解收敛得更快些。对每个弱约束项最优权重的确定是通过对各项相对误差的估计来实现。后来这个变分系统被应用到很多理想个例研究和实际观测个例的研究中，取得了较好的效果[23-26]。

以上的研究工作主要集中在对雷达径向风资料的同化上，但如何有效地利用雷达反射率资料是一个较为困难的课题。首先需要发展一个能把反射率和数值模式水凝物变量联系起来的较准确的正向算子。因为这个问题和数值模式复杂的微物理过程有关，因而有很多不确定性。其次是反射率和数值模式水凝物变量有很强的非线性关系，这使得在变分问题中数值解的收敛性存在困难。由于这些困难的存在，早期的反射率同化方案选择利用反射率做云分析。Albers等[27]发展了一个主要利用反射率的“局地分析和预报系统”。这个系统主要利用反射率来估计模式水凝物，包括雨水、雪和冰雹等，同时对云内水汽和温度做出调整。Zhang等[28]对这个系统做了些改进，使它更适合于对流尺度的暴雨系统。Hu等[23-24]、Schenkman等[26]利用此改进后的云分析系统同化反射率，并利用上述变分同化系统[20]分析雷达径向风到一个适用于对流尺度的区域预报系统[29-30]。结果显示该系统对短期强风暴，包含龙卷和冰雹等的超级单体有一定的预报能力。Hu等[23-24]的研究表明反射率的同化对暴雨的发展很重要。但是，云分析方案包含很多关于水凝物和反射率之间的经验关系，因此存在很多不确定性。这些不确定性使得云分析在对流尺度的天气数值预报中有一定的局限[31]。

Sun等[6]较早意识到同时定量同化径向风和反射率的重要性。在他们发展的适合于对流尺度的四维变分同化系统中，对几种反射率的同化作了系统的比较。在理想实验的框架里，他们发现最好的反射率同化方案是把反射率转化为雨水混合比，然后直接把雨水混合比同化到对流尺度数值模式中去。但这一方案的缺点是忽视了反射率中冰相（包括雪、冰雹等）的贡献，从而在他们的模式只有关于暖雨的微物理过程。他们也尝试直接同化反射率，发现这样做会增加问题的非线性，使得四维变分问题的求解更加困难。集合卡尔曼滤波是另外一个很有前途的同化雷达反射率和径向风的方法[32-39]。一般的结论是用此种方法同化雷达反射率和径向风资料能提高对流尺度的数值预报模式的预报水平。集合卡尔曼滤波同化雷达反射率资料主要的问题仍然是对流尺度模式的误差一般较大。因而和反射率有关的模式变量（如雨水、雪、冰雹等的混和比）的误差会通过协相关传播到该模式的其他变量（例如温度、水汽等），从而得到不合理的低层冷池分析结果[38]。

在Gao等[40]的研究中，与冰相微物理过程有关的对流尺度数值模式变量被引入进到三维变分的目标函数中。他们根据对流尺度数值模式的温度变量来决定雷达反射率的正向算子中水凝物（主要是雨、雪、冰雹等）的分类，这个研究使得在数值模式中对雷达反射率的同化更容易些。三维变分同化方法具有计算效率高、容易并行、编码简单等特性，所以比较适用于对流尺度的同化和预报系统。

这篇文章的主要目的是介绍一个基于三维变分同化的适用于对流尺度的雷达资料同化方法，特别是近几年来这个方法在实时分析和预报中的应用。

1 资料同化方法简介

1.1 目标函数的定义

正如上一节介绍的，Gao等发展了一个适用于对流尺度雷达资料同化的三维变分方案[18-20, 23-24, 41-44]。不同于大尺度的资料同化系统，这个方案强调利用质量连续性方程[18]、从适用于对流尺度的数值预报方程导出的气压诊断方程[44]作为弱约束来建立模式变量之间的联系。从而在同化中对没有观测的变量，例如垂直速度、温度等也做出一定程度的更新和改善。通过定义一个包括数值模式背景场项、雷达资料的观测项和弱约束项的目标函数，然后通过应用数学上的下降算法（例如共轭梯度法）使目标函数极小化来求得同化分析结果。同化分析结果可用于初始化对流尺度的数值模式。目前这个三维变分同化系统已经被用来初始化ARPS模式和WRF模式。

三维变分同化的目标函数可写为以下形式，

右边第一项定义分析场和模式背景场之间的距离。B是背景误差协方差距阵。第二项是观测项，代表分析变量和观测向量yo之间的距离。在最初的研究中，yo只包含雷达径向风，现在加入了雷达反射率和一些常规观测资料及风廓线资料等。关于雷达径向风和反射率的正向算子将在下文介绍。如果忽略观测资料之间的相关性，观测误差的协方差距阵R可以简化为对角阵，只包含了与各个观测相对应的仪器误差和代表性误差。第三项，JC(x)可以是包含质量连续性方程和数值模式方程在内的弱约束项。在目前的系统中，分析变量包括风场的3个分量、位温、气压、水汽混合比和模式水凝物（包括云水、云、冰、雨、雪、冰雹的混合比）等。

1.2 径向风正向算子

在方程（1）右边第二项的观测项可包含任何观测资料，侧重雷达资料的应用。根据Doviak等[45]，径向风算子可用以下公式表示：

式中，vr代表雷达径向风，r是雷达观测点与雷达位置的距离，h是雷达观测点离地球表面的高度，s是雷达射线投影到地面上的距离，是雷达观测点的径向角。如果在径向风算子中也考虑雷达的射线宽度，那么根据Rihan等[46]的研究，较复杂一点的径向风算子可表述为：

式中，是在某一高度角上的径向风，是数值模式的垂直格距，G描述了与雷达射线宽度有关的放大率函数[47]。

如果雷达径向风的观测误差为1或2 m/s[48-49]，方程（1）中的第二项可改写为

式中，H(vr,e)定义为把三维风场变量u，v，w投射到径向风方向的正向算子（式（2）或（3））和把径向风从模式格点投影到观测点上的三维线性插值算子。

1.3 反射率正向算子

早期的三维变分框架内，不包含对反射率的直接同化，反射率的分析是通过云分析来进行的[23-24]。Gao等[40]提出了一个在变分框架内直接同化反射率的方案。在这个方案中，首先发展了一个利用数值模式的温度背景场来对模式的水凝物进行分类的正向反射率算子。根据Lin等[50]、Gilmore等[51]和Doswell等[38]的研究，反射率因子的正向算子可表述为雨、雪、冰雹混合比的函数如下：

根据Smith等[52]的研究，雨水对反射率的贡献可表述为

ρ是空气密度。如果气温低于0℃，那么干雪对反射率的贡献可表述为，

如果气温高于0℃，那么湿雪对反射率的贡献可表述为，

根据Lin 等[50]和 Gilmore 等[51]的研究，冰雹对反射率的贡献为

用式（5）～（9）把反射率同化到数值模式中去不是很理想，存在着一定的缺陷。这是因为反射率是雨、雪、冰雹混合比三者的函数，直接用以上公式可能引起数值解的欠定性。例如，分析中可能导致在模式的低层产生雪的情况，但实际上模式低层由于温度较高却只能有雨。为解决这个问题，式（5）可被改进写成

式中，α随背景场温度的区间（Tb=－5℃和Tb=5℃）在0～1变化，Tb可从数值模式中得到。实际上用式（5）和（10）得到的反射率的因子大致相同。但在式（10）中，隐含着不同种类水凝物对反射率因子的贡献可以根据模式提供的温度场有所区分。从而当把反射率同化到数值模式时，使不同种类水凝物（雨、雪和冰雹的混合比）的分析结果更加合理。对反射率因子方程（10）取对数就可计算反射率，

这里，ZdB的单位是dBz。所以目标函数（1）中第二项关于反射率的观测项可写为，

类似于式（4），式（12）中的H(ZdB)定义为用方程（5）～（11）把模式变量，雨、雪和冰雹的混合比转化为反射率的非线性算子和把反射率从模式格点插值到雷达观测点上线性算子。

这种根据数值模式提供的温度场来确定水凝物对分射率贡献的正向算子使得对反射率的同化更精确，一定程度上去掉了解的不确定性。根据以前的对比试验，用温度场来确定水凝物的正向算子使得同化的结果更好，特别是对低层雨水的分析更合理。作者正在NOAA的“基于数值模式的预警”研究课题的资助下做进一步试验，希望将来能把此方法应用到实际业务中。

1.4 目标函数中的弱约束项Jc

在三维变分同化中，一些能把模式变量联系起来的物理约束（例如质量连续方程）和模式方程以及其导出的相关公式可以作为弱约束结合到目标函数中。在目前的系统中，Jc主要包括两个弱约束项即：质量连续性方程和从ARPS模式导出的气压诊断方程。这里主要介绍前者，后者可以使模式各个变量的分析更平衡[44]。

在这个三维变分同化中，3个风场分量之间的平衡关系是通过一个滞弹性近似的质量连续性方程来实现的，它的定义如下，

式中是给定层的平均空气密度。是个权重系数，它控制目标函数中这个约束项的相对重要性。一般取值是= 5.0×10-4。这个值也可通过对分析方案的调试得到。JDP定义气压诊断方程弱约束项[44]。Gao等[18,20]发现这个弱约束建立了风场3个分量之间的紧密联系，因而对垂直速度的分析非常有效。正如Gao等[18]在研究中指出的，在三维变分同化方案中把质量连续性方程作为弱约束而不是强约束[53-55]，避免了在垂直方向上积分连续性方程引起的误差累积。所以垂直速度的分析相对来说能更准确些。

1.5 变量替换和条件化

为了使数值求解更快，或者说使目标函数下降更快，定义控制变量v，使得C=，这样目标函数（式1）变为如下形式[56-57]：

式中H是H的线性化算子目标函数Jinc关于v的梯度和Hessian算子可从公式（15）导出，

式中，I是单位距阵。从式（17）中可看到，从式（1）到（15）的变换使相对应于新控制变量v的本征值不可能是1，从而降低了条件数，求解过程收敛会更快。

在式（15）中距阵C可分解为

式中，D是背景标准误差的对角化距阵，L是一个放大的系数。从理论上来说，F是一个对角元素是1，而非对角元素是背景协相关误差系数的矩阵；但是这个矩阵维数巨大，不可能用在数值求解过程中显示表示。Purser等[21]提出了一个利用对等空间递归滤波的方案来模拟F。一价的递归滤波定义如下

Xi是定义在格点i上的初始值，Yi是向前滤波后在格点i的值，Zi是向后滤波后在格点i的值，α是按式（20）定义的系数[58]，

式中，L是水平相关尺度，∆是水平格距，N是实施滤波的次数。这个递归滤波方法要求必须在正反两个方向各进行一次，以保证位相误差为0。这个递归滤波方法可分别应用在x、y、z这3个方向上。多尺度滤波可通过多次使用式（19）实现，每次应用不同的系数α。Xie等[59-60]在他们的三维变分研究中证明使用多尺度滤波比单尺度滤波效果会更好。

2 几个实时个例分析

在上述三维变分的框架上，发展了一个实时的可根据天气形势移动分析区域的分析系统。这个系统可用数值预报产品和实时收集的雷达资料进行分析[41]。这个系统有以下优点：1）可用多个雷达的资料，以及任何有关的常规数据；2）自动识别强天气可能发生的地点并提供有较高空间分辨率（1 km）和时间分辨率（5 min）的分析结果；也可由预报员根据当前天气形势确定分析区域；3）这个系统能够帮助预报员估计超级单体中层中尺度气旋的强弱从而在发布强天气预警时提供参考和帮助。这个系统曾经被应用到NOAA灾害性天气试验基地，并得到预报员的好评[61-62]。下面给出几个在美国产生强龙卷的灾害性天气试验个例和一个发生在湖南的强暴雨个例的分析结果。

2.Fig.年在俄克拉荷马州产生龙卷的一个超级单体

2013年5月20日，一个产生了强龙卷的超级单体袭击了美国俄克拉荷马州穆尔市（MOORE）。这个龙卷后来被美国气象部门归类于5级——龙卷强度的最高级别。这个龙卷共造成了23个人死亡，377个人受伤，并催毁了很多房屋（图1）。龙卷产生的瞬时最大风速达到了每小时340 km。龙卷的直径达2.1 km。根据事后收集的资料，龙卷在当地时间14：56从附近纽尔斯尔市开始形成，在地面穿行了大约27 km，大部分经过穆尔市并于当地时间15：35离开地面。这个产生了龙卷的超级单体被附近的多普勒雷达网监测到，其中包括有4个S波段的WSR-88D雷达，KTLX雷达离风暴系统最近，观测了系统的中低层的结构；KFDR，KINX和KVNX离风暴较远，主要观测到了风暴的高层结构。分析中的背景场是由美国环境预报中心（NCEP）的北美中尺度模式（NAM）提供的。在实时运行中，水平方向的分辨率设为1 km格距，模式区域共有200×200个格点。在垂直方向上有31层，采用逐步递增的不等距格距，平均格距是400 m。当天系统自动定位的分析中心位于（35.36°N，97.18°W）。分析主要集中在三维风场以及由此导出的垂直涡度和垂直速度上。

图1 2013年5月20日强龙卷袭击后的俄克拉荷马州穆尔市的住宅区 Fig.1 Selected images of tornado damage area in Moore, Oklahoma on 20 May 2013

2013年5月20日早晨，美国天气服务局暴雨预报中心发布了从密苏里州东南部到得克萨斯中北部的强天气预警。当地时间14：00左右，产生强龙卷的超级单体开始在雷达回波中显现。这个超级单体从出现到成熟发展得很快，在不到一个小时内（14：56），在纽尔斯尔市产生了第一个龙卷（图2）。当年，作者并没有在NOAA的强风暴试验基地每天实时运行分析系统，但在预报的强天气系统可能性提高时，分析系统会自动启动。系统在当地时间13：00开始实时接收雷达数据。系统自动定位了靠近穆尔市的分析区域。图2和图3显示了分析得到的位于地面高度以上3 km，从14：00—15：35的水平风场、垂直涡度和雷达反射率、以及相应的垂直剖面图。在这段时间内，超级单体经历了从生成、发展到成熟的阶段。14：00，在分析区域的西边界出现了几个雷达回波较弱的普通暴雨单体（图2a和图3a）。产生了穆尔市龙卷的暴雨单体迅速加强，仅仅25 min后，演变成超级单体。在对地面高度以上3 km上空的风场分析表明存在一个很强的中尺度气旋（图2b和图3b）。虽然分析也表明在这个超级单体的北部和南部各存在较强的超级单体，但它们都没产生龙卷（图2b和图2c）。穆尔市的龙卷在14：56下沉接连地面，在14：55的分析表明在雷暴中心存在超强的旋转和上升运动（图2e和图3e）。最强的垂直速度达到20 m/s，最强的垂直涡度达到0.015/s。在近一个小时期间，强大的气旋一直存在于超级单体内，在垂直方向上延伸几千米。当地时间15：20，这个超级单体分裂，强度减弱，龙卷于15：35时消散。为了验证分析是否合理，把从当地时间14：00—17：00分析中得到的垂直涡旋度的轨迹和气象局事后统计得到的龙卷轨迹一起画在图4中。从图中看出，三维变分分析的垂直涡旋度的轨迹前一段与龙卷造成损害的轨迹基本重合。三维变分的分析区域水平格距是1 km，这样的分辨率虽不可能分辨龙卷，但与龙卷紧密相关的中尺度气旋却能很清晰地被分析出来，并且伴随龙卷的发生一直非常强大。对这个例中类似的超级单体的预报是个非常困难的事情，因为这种超级单体在很短时间内快速加强并生成龙卷，对预报员提前准确发布强天气预警具有很大挑战。但快速三维变分同化的分析产品可以迅速地实时计算出气流旋转的强度，从而给预报员发布强天气预警提供及时帮助。

2.Fig.年在阿拉巴马州产生龙卷的多个超级单体

图2 分析得到的俄克拉荷马穆尔市地面以上3 km的多部雷达合成的反射率（彩色阴影区）、水平风场（风矢量）和垂直涡度（黑色等值线）（紫色实线表明图3垂直剖面图的位置） （a） 当地时间14：00，（b）14：20，（c）14：45，（d）14：50，（e）14：55，（f）15：00 Fig.2 The synthesized reflectivity (color shaded) and analyzed horizontal wind fields (vectors), and vertical vorticity (black contours) at 3 km AGL at (a) 14:00, (b) 14:20, (c) 14:45, (d) 14:50, (e) 14:55, and (f) 15:00 CST, 20 May 2013 near Moore, Oklahoma. Maroon line denotes location of cross-sections for vertical slices in Fig.3

图3 图例说明同图2，是图2中相应位置的垂直剖面图 Fig.3 Wind vectors and reflectivity (color shaded) through a vertical slice in Fig.2

2011年4月27日，在美国东南部几个州包括阿拉巴马、阿肯色、乔治亚、密西西比和田纳西等地暴发了剧烈的灾害性天气事件，是美国历史上死亡人数最多，损失最惨重的天气事件之一（图5）。总共有292个龙卷报告，大约造成了60亿美元的经济损失，348个人丧失生命。其中阿拉巴马州的损失最严重，4个龙卷被定性为5级（龙卷所能够造成破坏的最高级别）。仅在这个州就有238人死于这些龙卷造成的破坏。在该灾害事件发生的过程中，实时分析系统虽然在运行，但由于太多的超级单体，分析区域只覆盖了其中的一小部分。模式格点数为400×300，水平格距为1 km，垂直方向仍为不等距的31层。分析的中心位于（33.9°N，86.4°W），大约在伯明翰市西北80 km的地方。在这个分析区域，一共有7部WSR-88D的雷达资料被应用在分析中。分析区间是从当地时间17：00—18：30。图6显示了这段时间超级单体群从西南到东北方向移动。在分析区域中每个超级单体都可能成为一个移动的“小炸弹”，产生龙卷，袭击沿路的建筑及其他目标。选择对其中一个超级单体在18时的风场进行较为详细的分析（图6中的黑色方框）。这个超级单体产生了级别为4级的龙卷，袭击了阿拉巴马州的多拉市（Dora）。在这个时刻，龙卷已经产生。图7a显示离地面1.5 km高度上的风场环流叠加在反射率的钩状回波上。在中层5 km高度上形成直径为8～10 km的闭合环流（图7b），在较高层（10 km高度上）风场的形势是伴随着旋转的向下风方向的辐散。在垂直剖面上存在着大的上升运动区，与强反射率区域重合。最大上升速度达到20 m/s（图7d）。图8显示无论在低层、中层还是高层，暴雨区域的强垂直涡度区和强垂直速度区相重合。表明在这种产生了龙卷的超级单体中存在着极强的与旋转相伴随的上升运动。但在垂直方向上最大涡度中心和最强上升中心不重合，前者所在层次低于后者。在很多其他产生龙卷的超级单体的分析中也发现了类似现象。另外，在风暴的上层存在重力波从风暴中心向风暴的下风方向的传播。把从三维变分分析得到的每个格点垂直方向上最大涡度和最大垂直速度做成一个风暴移动轨迹显示在图9中。在分析中得到的对应于超级单体的垂直涡度中心和垂直速度中心的移动与美国天气服务局暴雨中心得到的龙卷的破坏规迹相重合（图10）。

图4 由三维变分分析得到的俄克拉荷马州穆尔市地面到3 km的当地时间14：00—17：00的垂直涡旋度（Helicity）（图中深蓝线是被龙卷破坏的轨迹图） Fig.4 3DVAR 0-3 km MSL updraft helicity accumulated over the period from 14:00 to 17:00 CST on 20 May 2013 in Central Oklahoma. The dark blue lines indicate tornado damage path as determined by the National Weather Service Forecast Office at Norman, Oklahoma

图5 2011年4月27日强龙卷袭击后的阿拉巴马州住宅区 Fig.5 Selected images of tornado damage area in Alabama on 27 April 2011

2.Fig.年一个在得克萨斯州产生龙卷的超级单体

2010年4月20日，一个产生了龙卷的超级单体袭击了得克萨斯州西北部区域，这一区域又被称为得州刀把区。这个超级单体生成于当地时间17时，终止于当地时间21时，生命期大约为4 h。在这个超级单体的生命期内有2个龙卷生成，美国气象局暴雨预报中心事后也有收集到关于冰雹和强阵风的报告（图11），幸运的是并无人员伤亡。这个风暴附近的3个WSR-88D的雷达KAMA、KLBB和KFDX实时观测到了这个风暴的演变。

这个超级单体是在得州和新墨西哥州边界形成的。在形成的最初两个小时内，它的风场旋转不太强。但到当地时间接近19时，风暴中的中尺度气旋突然加强。19：04，第一个较弱龙卷形成，之后风暴继续发展和加强。实时分析系统“跟踪了”这个风暴的演变。分析得到的水平风场、垂直涡度和反射率都表明风暴从19：30—20：00发展到顶峰（图12）。在这半小时内，地面以上3 km高度上的最大垂直涡度一直维持在0.01 s-1左右。当地时间19：45和19：50，由垂直涡度表征的中尺度气旋在垂直方向上延伸达到了10 km，最强的涡度中心和强反射率区域相重合（图13）。在风暴中心地面上空4 km处存在一个弱回波区和较强的低层入流区（图13a、13c）。第二个龙卷于19：47生成（NOAA暴雨中心的报告，图11a）。从19：35到20：00的反射率演变可以看到勾状云逐渐锢囚或闭合的一个过程（图14a～14f）。第二个龙卷刚好是在这个锢囚过程中形成的。在这个过程中，这个超级单体也降下了很多冰雹。

2.Fig.年发生在湖南的一场大暴雨过程

北京时间2010年5月6日00：00—06：00，湖南省中南部怀化、娄底等6个地区出现暴雨，其中7个县级市出现大暴雨，日最大降水出现在娄底新化为147.8 mm；据区域气象观测站资料，153乡镇自动气象站出现暴雨，77乡镇自动气象站出现大暴雨。本次降水特点为突发性强、短时雨强大、强降水区呈狭长带状分布。此次过程造成娄底、怀化、长沙等8个市（州）27个县（市、区）336个乡镇272.35万人受灾，因灾死亡11人，失踪1人，紧急转移人口6.01万人，倒塌房屋9944间，直接经济损失14.53亿元。

图6 阿拉巴马州在地面以上3 km的多部雷达合成的反射率（彩色阴影区）、水平风场（风矢量）和垂直涡度（黑色等值线）（图6e中的黑色小框是图7和图8的分析区域） （a）当地时间17：00，（b）17：15，（c）17：30，（d） 17：45，（e）18：00，（f）18：15 Fig.6 The synthesized reflectivity (color shaded) and analyzed horizontal wind fields (vectors), and vertical vorticity (black contours) at 3 km AGL, at (a) 17:00, (b) 17:15, (c) 17:30, (d) 17:45, (e) 18:00, and (f) 18:15 CST, 27 April 2011 near North Alabama and nearby area. The small box in Fig.6e is the domain shown in Fig.7 and Fig.8

鉴于此次暴雨过程从6日当地时间02：00开始加强，至08：00减弱，因此主要分析这个时段的雷达资料。该时段的降水情况如图15显示。对湖南境内的两部多普勒雷达资料进行了三维变分分析，图16显示的是2010年5月6日暴雨过程的分析区域和雷达位置，这两部雷达分别是位于长沙的CINRAD/SA雷达和位于常德的CINRAD/SB雷达。同化的水平分辨率是1 km，平均垂直分辨率为500 m。模式水平格点数为403×403，垂直层数为25。

图7 从图6中放大的阿拉巴马州多拉市在地面以上的多部雷达合成的反射率（彩色阴影区，单位：dBz）和水平风场（风矢量，单位：m/s） （a）地面以上1.5 km，（b）地面以上5 km，（c） 地面以上10 km ，（d）从y=140 km的垂直剖面图 Fig.7 The synthesized reflectivity (color shaded, unit: dBz) and horizontal wind fields (vectors, unit: m/s) which were zoomed in from the box in Fig.6e, at (a) 1.5 k m AGL, (b) 5 km AGL, (c) 10 km AGL, (d) vertical profile at y=140 km, on 27 April 2011 near Cordova, Alabama

从同化后雷达回波和风场的图中可以看出（图17），当地时间6日03时，在湖南中部地区出现从常德南部到益阳西南部和从湘西自治州东南部到怀化中东部两条回波带与娄底—湘潭中部—株洲西北部回波带汇合，呈“人”字型回波带（图17a）；04：00时回波北扩，强回波范围明显增大，＞45 dBz的回波连成线状，此时反射率因子达到最强，强回波仍然稳定在湘中一线；05：00—06：00时，回波持续北扩，同时在06：00衡阳北部地区又生成了新的东西向的回波带，但回波强度明显减弱（图17d）。07：00，湘中偏北位置的回波东移减弱，衡阳北部的回波仍然持续（图17e）；6日08：00后回波逐渐减弱（图17f）。

多普勒雷达径向速度无法直接反映风场信息，而同化雷达资料后的风场更为直观，并且包含更多中尺度信息：03：00出现了一个呈西北—东南向的弱辐合区，这预示着将有一个较强的中尺度对流系统发展（图17a）。04：00这个辐合区更加明显，同时，在湘中强回波带对应的区域，风的流线有明显的气旋性弯曲且有风速上的辐合，为强降水的发展提供了有利的动力条件（图17b）。

需要说明的是，目前中国天气雷达主要是从美国引进并国产化的新一代多普勒雷达，因此从类型上来说与美国主力雷达系统没有太大差异。并且WRF，ARPS等模式都能直接处理国内的雷达资料。但国内雷达数据质量还不够好。造成质量问题的原因很多，首先在客观上来说，中国地形较为复杂，有些雷达选址不当就会导致挡角影响明显，往往会出现底层几个仰角没有有效观测数据。比如湖南省怀化雷达的挡角问题就比较明显。其次，雷达的生产工艺及安装后的系统参数定标不当，也会影响雷达数据质量，比如杂波偏多，系统显示不清晰，有时对同一系统在两个相邻雷达的回波强度或径向速度有明显不同。由于国外雷达数据质量相对较好，因此国外的数值模式对国内的雷达数据往往水土不服，同化结果中，天气系统显示不够清晰，杂波偏多，影响了同化后的预报效果，也不利于天气系统的分析。因此，同化国内雷达数据时要特别加强雷达质量控制。这些也是今后需要研究的课题。另外做的研究例子还不是很多，将来也应做些关于中国龙卷的例子。

图8 图例说明同图7，图中彩色阴影区为垂直速度，黑色等值线为垂直涡度（粉色线表明重力波传播的位置） Fig.8 As same as in Fig.7, color shaded is for vertical velocity and contours are for vertical vorticity (pink curve indicates the propagation of gravity waves)

图9 阿拉巴马州中北部及其附近区域，由三维变分分析得到的垂直涡度（a）和垂直速度（b） Fig.9 The composite vertical vorticity and vertical velocity tracks derived from 3DVAR analyses at North and Central Alabama and nearby areas

图10 2011年4月27日美国东南部的几个州发生的龙卷、冰雹和大风实况报告 Fig.10 The storm report for tornadoes (red), hails (green) and damaging winds (blue) in southeastern U.S. on 27 April 2011, from NOAA/storm prediction center

3 结论与展望

在过去几年中，发展了一个基于三维变分的高分辨率（1 km）和高频率（每5 min）的实时分析系统。这个系统可根据天气实况确定分析区域，实时读取分析数据，并做出分析。在过去的一些研究和应用中作者对这个系统不断地完善。在这个三维变分同化中，对代价函数的定义只包含一个背景场、雷达径向风和一个弱约束项。其中的背景误差协方差距阵是由简单的递归滤波来实现。最近几年在变分同化的框架内增加了对反射率的同化，同时也引进了用集合扰动场来导出背景误差协方差距阵，和引进更多观测资料项等，以便使分析更精确稳定，同时有能力对对流尺度的数值预报模式提供初始场。

本文介绍了基于三维变分的实时分析系统。这个系统在NOAA的强灾害天气试验基地进行了准业务运行。它能根据实时天气形势确定分析区域，并利用实时雷达观测数据和数值预报产品快速做出分析。特别是能实时估计出强风暴中中尺度气旋环流的强度，从而给预报员发布强天气预警提供帮助。在试验时得到很多来试验基地访问的预报员的好评。文中典型超级单体和暴雨的分析个例表明这个系统能很好地模拟超级单体的发展和较准确地提供龙卷的预警。但这个分析系统也存在很多需要改进的地方。首先分析产品仍然滞后7～8 min，所以接收数据的时间需缩短，后处理的时间也需缩短。另外这个系统对于较弱的风暴分析不是很合理，还有待提高。

这个三维变分同化系统在其他对流尺度的分析和应用研究方面也发挥了很大作用（如Gao等[63]；Stensrud等[42]；Xue等[64]；Schenkman等[26]）。例如，在美国大气遥感工程研究中心（The Engineering Research Center for Collaborative Adaptive Sensing of the Atmosphere）建立的X波段雷达网的观测区域进行的高分辨率低层风场分析和超短期的预报[63, 65]。从2008年以来，这个系统也被应用到NOAA灾害性天气试验基地春季预报试验阶段来为高分辨率WRF预报模式的确定性和集合预报提供初值。这个试验的分析预报区域覆盖全美大陆地区，用到几乎所有WSR-88D雷达网提供的数据[64, 66]。这个春季试验主要是为了用计算机和数值模式提高对灾害性强风暴天气的预报能力，从而对NOAA暴雨预报中心的业务提供支撑。

图11 （a）2010年4月20日得克萨斯州刀把区的龙卷和冰雹实况报告；（b）最内黑框为分析区域（200 km×200 km）（外黑框用来确定所用雷达数目（500×500 km）；有3个WSD-88D雷达（代码分别是：KAMA、KFDX 和KLBB）被分析采用（美国气象学会提供）） Fig.11 (a) Reports of the tornadic storm event in Texas Panhandle on 20 April 2010; (b) The 3 radars of WSD-88D: the KAMA, KFDX and KLBB (in black box of 500×500 km) used for analysis, and the analysis area (internal 200 km×200 km) (courtesy of the American Meteorology Society)

图12 2010年4月20日分析得到的得克萨斯州刀把区域在地面以上3 km的多部雷达合成的反射率（彩色阴影区）、水平风场（风矢量）和垂直涡度（黑色等值线） （a）19：35，（b）19：40，（c）19：45，（d）19：50，（e）19：55，（f）20：00（以上均为当地时间，紫色实线表明图13中垂直剖面图的位置；图12a中小框为图14的分析区域（美国气象学会提供）） Fig.12 The synthesized reflectivity (color shaded, dBz), horizontal wind fields (vectors, m/s), and vertical vorticity (contour interval, 2×10-3 s-1) at 3 km AGL observed by 88D radars shown in Fig.6b, at (a) 19:35, (b) 19:40, (c) 19:45, (d) 19:50, (e) 19:55, and (f) 20:00 CST, 20 April 2010 near Umbarger, Texas. Maroon line denotes location of cross-sections in Fig.13. Black box in Fig.12a is the zoom-in area for Fig.14 (courtesy of the American Meteorology Society)

图13 图12中的两个时次的垂直剖面图 （a）和（b）为垂直涡度（黑色等值线）叠加于反射率； （c）和（d）为风矢量叠加于反射率（美国气象学会提供） Fig.13 Vertical vorticity (black contours in a and b), wind vectors (in c and d), and reflectivity (color shaded) along a vertical slice indicated in Fig.12, at 19:45 CST (left panel), and 19:50 CST (right panel), on 20 April 2010 near Umbarger, Texas. Between these two time levels, a tornado touched down (courtesy of the American Meteorology Society)

但三维变分系统在理论上不是最优的方法，存在着很多理论上的缺陷。例如不能引入观测资料的时间分布信息，背景误差协方差距阵不能随天气形势演变等。由于这些缺陷，不少对流尺度的资料同化和预报的研究工作中也探索了运用集合卡尔曼滤波和四维变分同化的方法。但是这些方法各有优缺点[67-68]。综合考虑各种因素，例如现阶段对流尺度，或对流尺度数值模式的准确性，预报又必须及时发布（时效性），仍然认为用改进的三维变分同化方法对高分辨率数值预报提供初值仍然是一个方便、高效和实用的选择。首先可以加入一些平衡约束来使分析变量之间保持平衡，同时也可把集合预报的信息引入到三维变分系统用来改进系统的误差协方差距阵，使之随流型的改变而变。Lorenc[67]首先提出了将变分同化和集合卡尔曼滤波相结合的方法。在这个方法中，变分同化系统可以是三维变分同化，也可以是四维变分同化。

在最近几年出现了大量的混和同化方面的研究，大部分研究集中在对大尺度和中尺度天气的资料同化及数值预报方面[9, 69-75]。对对流尺度的研究相对较少。Gao等[76]在他们的研究中证明混合三维变分和集合卡尔曼滤波方法在对流尺度的数值预报也很有用。把从集合预报中估计的背景误差协方差距阵用在三维变分同化中可以有效地提高雷达资料同化的精度，甚至在集合数量不太大的情况下也是如此。用理想数据的数值试验表明这个方法优于三维变分的同化方法。但在这个方法中必须同时使用三维变分和集合卡尔曼滤波。因为雷达资料分辨率高，资料量巨大，使得计算量巨大，特别是对高分辨率的对流尺度的数值模式。为减少计算量，同时利用变分和集合卡尔曼方法的优点，Gao等[77]提出了采用集合变分同化的方法把雷达资料同化到对流尺度的数值模式中。集合变分同化的概念和方法首先在法国气象局和欧州中期数值预报中心使用[78-79]。在这里Gao等[77]利用增加alpha控制变量的方法[67]来计算集合误差协方差距阵。如果集合变分资料同化的核心框架用的是四维变分同化，那么这个方法计算量也会很大。但是Gao等采用的是计算效率较高的三维变分同化方法，同时这种方法不需运行集合卡尔曼滤波，所以此方法在集合数量不是太大的情况下计算代价较小，容易并行化。这种方法的优点是能够应用从集合预报中导出随时间和流型而变的集合误差背景场协方差距阵，所以较适用于对流尺度的雷达资料同化。现在作者正在集中精力调试这些方案，希望将来能有效地应用到NOAA的灾害性天气试验基地的短期天气预报试验中，能对基于这些预报的强风暴预警提供支撑。

图14 图例说明同图12，但分析区域为图12a中的被放达的黑边小框区域（美国气象学会提供） Fig.14 As same as in Fig 12, but the analysis area is for the zoom-in area in Fig.12a (courtesy of the American Meteorology Society)

说明和致谢：本文有几位作者在过去几年中在俄克拉荷马大学作访问学者时参与了此软件的研究与开发，包括湖南省气象台傅承浩、唐杰、周长青，以及国家气象中心庄照荣等。在研究过程中曾得到俄克拉荷马大学超级计算中心Henry Neeman博士和他的团队在计算机资源及应用方面的大力帮助。在此一并表示感谢！

图15 2010年5月6日湖南省自动雨量站逐小时累积雨量图 Fig.15 Precipitation observed by the surface rain gauge

图16 2010年5月6日湖南暴雨过程的分析区域和雷达位置（内黑框为分析区域）Fig.16 Positions of two radars (9731, 9736) and the analysis domain for Hunan heavy rain event on 6 May 2010 (The internal domain is the analysis domain in Hunan Province)

图17 2010年5月6日同化得到的3 km高度上的雷达回波和风场叠加 （a）03：00，（b）04：00，（c）05：00，（d）06：00，（e）07：00，（f）08：00 Fig.17 The analyzed u-v wind vector and reflectivity (colored) at z=3 km AGL for a local heavy rain and severe convective weather event occurred at central south of Hunan Province in China: at 03:00 (a), 04:00 (b), 05:00 (c), 06:00 (d), 07:00 (e), 08:00 (f) BT on 6 May 2010

[1] Lewis J, Derber J. The use of adjoint equations to solve a variational adjustment problem with advective constraints. Tellus, 1985, 37A: 307-322.

[2] Le Dimet F, Talagrand O. Variational algorithms for analysis and assimilation of meteorological observations: theoretical aspects. Tellus, 1986, 38: 97-110.

[3] Talagrand O, Courtier P. Variational assimilation of meteorological observations with adjoint vorticity equation. I: theory. Quarterly Journal of the Royal Meteorological Society, 1987, 113: 1311-1328.

[4] Thacker C, Long R. Fitting dynamics to data. Journal of Geophysical Research, 1988, 93: 1227-1240.

[5] Daley R. Atmospheric data analysis. Cambridge: Cambridge University Press, 1991: 457.

[6] Sun J, Crook N A. Dynamical and microphysical retrieval from Doppler radar observations using a cloud model and its adjoint. Part I: model development and simulated data experiments. Journal of the Atmospheric Sciences, 1997, 54: 1642-1661.

[7] Rabier F, Järvinen H, et al. The ECMWF operational implementation of 4D variational assimilation. Part I: experimental results with simplified physics. Quarterly Journal of the Royal Meteorological Society, 2000, 126: 1143-1170.

[8] Barker D, Huang W, Guo Y, et al. A three-dimensional (3DVAR) data assimilation system for use with MM5: implementation and initial results. Monthly Weather Review, 2004, 132: 897-914.

[9] Barker D, Huang X, et al. The Weather research and forecasting model’s community variational/ensemble data assimilation system: WRFDA. Bulletin of the American Meteorological Society, 2012: 93, 831-843.

[10] Xiao Q, Kuo Y H, Sun J, et al. Assimilation of Doppler radar observations with a regional 3DVAR system: impact of Doppler velocities on forecasts of a heavy rainfall case. Journal of Applied Meteorology, 2005, 44: 768-788.

[11] Lewis J, Lakshmivarahan S, Dhall S. Dynamic data assimilation: a least squares approach. Cambridge University Press, 2006: 654.

[12] Sun J, Flicker D, Lilly D. Recovery of three-dimensional wind and temperature fields from simulated Doppler radar data. Journal of the Atmospheric Sciences, 1991, 48: 876-890.

[13] Xu Q, Qiu C. Simple adjoint methods for single-Doppler wind analysis with a strong constraint of mass conservation. Journal of Atmospheric and Oceanic Technology, 1994, 11: 289-298.

[14] Xu Q, Qiu C. Adjoint-method retrievals of low-altitude wind fields from single-Doppler reflectivity and radial-wind. Journal of Atmospheric and Oceanic Technology, 1995, 12: 1111-1119.

[15] Xu Q, Gu H, Yang S. Simple adjoint method for three-dimensional wind retrievals from single-Doppler radar. Quarterly Journal of the Royal Meteorological Society, 2001, 127: 1053-1067.

[16] Qiu C, Xu Q. Least Squares Retrieval of Microburst Winds from Single-Doppler radar data. Monthly Weather Review, 1996, 124: 1132-1144.

[17] Sun, J, Crook N A. Real-time low-level wind and temperature analysis using single WSR-88D data. Weather and Forecasting, 2001, 16: 117-132.

[18] Gao J, Xue M, et al. A variational method for the analysis of threedimensional wind fields from two Doppler radars. Monthly Weather Review, 1999, 127: 2128-2142.

[19] Gao J, Xue M, et al. New development of a 3DVAR system for a nonhydrostatic NWP model. Preprint, 15th Conf. Numerical Weather Prediction and 19th Conf. on Weather Analysis and Forecasting, San Antonio, TX, American Meteorological Society, 2002, 339-341.

[20] Gao J, Xue M, et al. A three-dimensional variational data assimilation method with recursive filter for single-Doppler radar, Journal of Atmospheric and Oceanic Technology, 2004, 21: 457-469.

[21] Purser R J, Wu W, et al. Numerical aspects of the application of recursive filters to variational statistical analysis. Part I: Spatially homogeneous and isotropic Gaussian covariances. Monthly Weather Review, 2003, 131: 1524-1535.

[22] Wu W, Purser R J, Parrish D. Three-dimensional variational analysis with spatially inhomogeneous covariances. Monthly Weather Review, 2002, 130: 2905-2916.

[23] Hu M, Xue M, Brewster K. 3DVAR and cloud analysis with WSR-88D Level-II data for the prediction of the fort worth tornadic thunderstorms. Part I: cloud analysis and its impact. Monthly Weather Review, 2006, 134: 675-698.

[24] Hu M, Xue M, Gao J, et al. 3DVAR and cloud analysis with WSR-88D Level-II data for the prediction of the fort worth tornadic thunderstorms. Part II: impact of radial velocity analysis via 3DVAR. Monthly Weather Review, 2006, 134: 699-721.

[25] Gao J, Stensrud D, et al. Three-dimensional analyses of several thunderstorms observed during VORTEX2 field operations. 34th Conference on Radar Meteorology, Willimsburg, VA, 2009, Online publication.

[26] Schenkman A, Xue M, et al.The analysis and prediction of the 8-9 May 2007 Oklahoma tornadic mesoscale convective system by assimilating WSR-88D and CASA radar data using 3DVAR. Monthly Weather Review, 2011, 139: 224-246.

[27] Albers S, McGinley J, et al. The local analysis and prediction system (LAPS): analysis of clouds, precipitation and temperature. Weather and Forecasting, 1996, 11: 273-287.

[28] Zhang J, Carr F, Brewster K. ADAS cloud analysis. Preprints, 12th Conf on Numerical Weather Prediction, Phoenix, AZ, American Meteorological, Society, 1998: 185-188.

[29] Xue M, Droegemeier K K, Wong V. The advanced regional prediction system (ARPS): a multiscale nonhydrostatic atmospheric simulation and prediction tool. Part I: model dynamics and verification. Meteorology and Atmospheric Physics, 2000, 75: 161-193.

[30] Xue M, Wang D, et al. The advanced regional prediction system (ARPS), storm scale numerical weather prediction and data assimilation. Meteorology and Atmospheric Physics, 2003, 82: 139-170.

[31] Gao J, Ge G, et al. The relative importance of assimilating radial velocity and reflectivity data to storm-scale analysis and forecast. 23nd Conf Wea Anal Forecasting/19th Conf Num Wea Pred, Omaha, NB, American Meteorology Society, 2009, Paper 16A.3.

[32] Snyder C, Zhang F. Assimilation of simulated Doppler radar observations with an ensemble Kalman filter. Monthly Weather Review, 2003, 131: 1663-1677.

[33] Zhang F, Snyder C, Sun J. Impacts of initial estimate and observations on the convective-scale data assimilation with an ensemble Kalman filter. Monthly Weather Review, 2004, 132: 1238-1253.

[34] Tong M, Xue M. Ensemble Kalman filter assimilation of Doppler radar data with a compressible nonhydrostatic model: OSS experiments. Monthly Weather Review, 2005, 133: 1789-1807.

[35] Meng Z, Zhang F. Tests of an ensemble Kalman filter for mesoscale and regional-scale data assimilation. Part II: imperfect model experiments. Monthly Weather Review, 2007, 135: 1403-1423.

[36] Meng Z, Zhang F. Limited-area ensemble-based data assimilation. Monthly Weather Review, 2011, 139: 2025-2045.

[37] Gao J, Xue M, An efficient dual-resolution approach for ensemble data assimilation and tests with assimilated Doppler radar data. Monthly Weather Review, 2008, 136: 945-963.

[38] Dowell, D, Wicker L, et al. Ensemble Kalman filter assimilation of radar observations of the 8 May 2003 Oklahoma City supercell: influences of reflectivity observations on storm-scale analyses. Monthly Weather Review, 2011, 132: 1982-2005.

[39] Yussouf N, Mansell E R, et al. The ensemble Kalman filter analyses and forecasts of the 8 May 2003 Oklahoma city tornadic supercell storm using single-and double-moment microphysics schemes. Monthly Weather Review, 2013, 141: 3388-3412.

[40] Gao J, Stensrud D. Assimilation of reflectivity data in a convectivescale, cycled 3DVAR framework with hydrometeor classification. Journal of the Atmospheric Sciences, 2012, 69: 1054-1065.

[41] Gao J, Smith T, et al. A realtime weather-adaptive 3DVAR analysis system for severe weather detections and warnings with automatic storm positioning capability. Weather and Forecasting, 2013, 28: 727-745.

[42] Stensrud D, Gao J. Importance of horizontally inhomogeneous environmental initial conditions to ensemble storm-scale radar data assimilation and very short range forecasts. Monthly Weather Review, 2010, 138: 1250-1272.

[43] Ge G, Gao J, et al. Impacts of beam broadening and earth curvature in radar data assimilation. Journal of Atmospheric and Oceanic Technology, 2010, 27: 617-636.

[44] Ge G, Gao J, Xue M. Diagnostic pressure equation as a weak constraint in a storm-scale three dimensional variational radar data assimilation system. Journal of Atmospheric and Oceanic Technology, 2012, 29: 1075-1092.

[45] Doviak, R, Zrnic D. Doppler Radar and Weather Observations. 2nd ed. Academic Press, 1993: 562.

[46] Rihan F A, Collier C G, et al. Assimilation of Doppler radial winds into a 3D-Var system: errors and impact of radial velocities on the variational analysis and model forecasts. Quarterly Journal of the Royal Meteorological Society, 2008, 134: 1701-1716.

[47] Wood V T, Brown R A. Effects of radar sampling on single-Doppler velocity signatures of mesocylcones and tornadoes. Weather and Forecasting, 1997, 12: 928-938.

[48] Doviak, R, Ray P, et al. Error estimation in wind fields derived from dual-Doppler radar measurement. Journal of Applied Meteorology, 1976, 15: 868-878.

[49] Miller L J, Sun J. Initialization and forecasting of thunderstorm: specification of radar measurement errors. Preprints, 31st Conf. on Radar Meteorology, Seattle, WA, American Meteorological Society, 2003: 146-149.

[50] Lin Y, Farley R D, Orville H D. Bulk parameterization of the snow field in a cloud model. Journal of Applied Meteorology, 1983, 22: 1065-1092.

[51] Gilmore M, Straka J, Rasmussen E. Precipitation and evolution sensitivity in simulated deep convective storms: comparisons between liquid-only and simple ice and liquid phase microphysics. Monthly Weather Review, 2004, 132: 1897-1916.

[52] Smith P, Myers C, Orville H. Radar reflectivity factor calculations in numerical cloud models using bulk parameterization of precipitation processes. Journal of Applied Meteorology, 1975, 14: 1156-1165.

[53] Sasaki Y. Some basic formalisms in numerical variational analysis. Monthly Weather Review, 1970, 98: 875-883.

[54] Sasaki Y, Numerical variational analysis formulated under the constraints as determined by longwave equations and a lowpass filter. Monthly Weather Review, 1970, 98: 884-898.

[55] Sasaki Y. Numerical variational analysis with weak constraint and application to surface analysis of severe storm gust. Monthly Weather Review, 1970, 98: 899-910.

[56] Courtier P, Thépaut J, et al. A strategy for operational implementation of 4D-Var, using an incremental approach. Quarterly Journal of the Royal Meteorological Society, 1994, 120: 1367-1388.

[57] Courtier P. Dual formulation of four dimensional variational assimilation. Quarterly Journal of the Royal Meteorological Society, 1997, 123: 2449-2461.

[58] Lorenc A. Iterative analysis using covariance functions and filters. Quarterly Journal of the Royal Meteorological Society, 1992, 118: 569-591.

[59] Xie Y, Koch S E, et al. A sequential variational analysis approach for mesoscale data assimilation. Preprints, 21st Conf. on Weather Analysis and Forecasting/17th Conf on Numerical Weather Prediction, Washington, DC, American Meteorological Society, 2005, 15B.7. http://ams.confex.com/ams/pdfpapers/93468.pdf.

[60] Xie Y, Koch S E, et al. A Space and time multiscale analysis system: a sequential variational analysis approach. Monthly Weather Review, 2011, 139: 1224-1240.

[61] Smith T, Gao J, et al. Performance of a real-time 3DVAR analysis system in the hazardous weather testbed. Weather and Forecasting, 2014, 29: 63-77.

[62] Calhoun K, Smith T, et al. Forecaster use and evaluation of realtime 3DVAR analyses during severe thunderstorm and tornado warning operations in the hazardous weather testbed. Weather and Forecasting, 2014, 29: 601-613.

[63] Gao J, Brewster K, et al, Real-time, low-level wind analysis including CASA and WSR-88D radar data using the ARPS 3DVAR. 25th Conf Severe Local Storms, American Meteorological Society, 2010: 7B.4 (online publication).

[64] Xue M, Kong F, et al. CAPS realtime storm scale ensemble and high resolution forecasts for the NOAA hazardous weather testbed 2010 spring experiment. 25th Conf Severe Local Storms, American Meteorological Society, 2010: 7B.3.

[65] Brewster K, Thomas K, et al. A nowcasting system using full physics numerical weather prediction initialized with CASA and NEXRAD radar data. Preprints, 25th Conf Severe Local Storms, Denver, CO, American Meteorological Society, Denver, CO, 2010: 9.4.

[66] Kain J S, Xue M, Coniglio M C, et al. Assessing advances in the assimilation of radar data within a collaborative forecastingresearch environment. Weather and Forecasting, 2010, 25: 1510-1521.

[67] Lorenc A. The potential of the ensemble Kalman filter for NWP: a comparison with 4DVar. Quarterly Journal of the Royal Meteorological Society, 2003, 129: 3183-3204.

[68] Kalnay E, Li H, et al. 4-D-Var or ensemble Kalman filter? Tellus, 2007, 59: 758-773.

[69] Buehner M. Ensemble-derived stationary and flow-dependent background-error covariances: evaluation in a quasi-operational NWP setting. Quarterly Journal of the Royal Meteorological Society, 2005,131: 1013-1043.

[70] Buehner M, Houtekamer P, et al. Intercomparison of variational data assimilation and the ensemble Kalman filter for global deterministic NWP. Part I: description and single-observation experiments. Monthly Weather Review, 2010, 138: 1550-1566.

[71] Buehner M, Houtekamer P, et al. Intercomparison of variational data assimilation and the ensemble Kalman filter for global deterministic NWP. Part II: one-month experiments with real observations. Monthly Weather Review, 2010, 138: 1567-1586.

[72] Wang X, Barker D M, et al. A Hybrid ETKF–3DVAR data assimilation scheme for the WRF model. Part I: observing system simulation experiment. Monthly Weather Review, 2008, 136: 5116-5131.

[73] Wang X, Barker D M, et al. A Hybrid ETKF–3DVAR data assimilation scheme for the WRF model. Part II: real observation experiments. Monthly Weather Review, 2008, 136: 5116-5131.

[74] Wang X, Parrish D, et al. GSI 3DVar-based ensemble-variational hybrid data assimilation for NCEP global forecast system: singleresolution experiments. Monthly Weather Review, 2013, 141: 4098-4117.

[75] Zhang F, Zhang M, Poterjoy J. E3DVar: coupling an ensemble Kalman filter with three-dimensional variational data assimilation in a limited-area weather prediction model and comparison to E4DVar. Monthly Weather Review, 2013, 141: 900-917.

[76] Gao J, Stensrud D. Some observing system simulation experiments with a hybrid 3DEnVAR system for storm scale radar data assimilation. Monthly Weather Review, 2014, 142: 3326-3346.

[77] Gao J, Fu C, et al. OSSEs for an Ensemble 3DVAR data assimilation system with radar observations of convective storms. Journal of the Atmospheric Sciences, 2016, 73: 2403-2426.

[78] Bonavita M, Isaksen L, Holm E. On the use of EDA background error variances in theECMWF 4D-Var. Quarterly Journal of the Royal Meteorological Society, 2012, 138: 1540-1559.

[79] Berre L, Pannekoucke O, Desroziers G, et al. A variational assimilation ensemble and the spatial filtering of its error covariances: increase of sample size by local spatial averaging. Proceedings of the ECMWF Workshop on flow-dependent aspects of data assimilation ECMWF, 2007: 151-168.