引言

随着同化技术的发展，常规观测资料同化已经不能够满足数值预报对初始场质量的要求，不断提高对非常规观测(如卫星和雷达)资料的同化使用能力是目前国内外众多科研和业务单位的主要研究方向。FY-3A成功发射后，在数值天气预报中的同化应用情况立即受到广泛关注，并利用欧洲数值预报中心的工作平台进行了“FY-3A卫星资料在欧洲数值预报中心(ECMWF)同化应用”的研究(李俊等，2012)。

卫星辐射资料直接进入模式需要通过辐射传输模式，目前数值模式中接入的快速辐射传输模式主要有CRTM和RTTOV两种。对于卫星红外通道的资料同化，均需要预先剔除其中的云区资料(董佩明等，2008)；而对于微波通道，在不考虑散射作用的前提下，则可以将云区的观测假设为晴空观测来使用。Xu等(2006)使用CRTM，结合NCEP的GSI三维变分同化模块对ATOVS卫星资料做了一系列直接同化试验，结果表明：卫星资料的直接同化能够对暴雨的24 h预报提供较好的正效果。希爽等(2010)利用在WRF三维变分同化系统中同化了AMSU-A微波资料，发现同化试验对台风预报的正效果要大于负效果。黄兵等(2006)利用GRAPES 3Dvar针对中尺度强暴雨过程进行了ATOVS资料的直接变分同化和数值模拟试验，结果表明：同化对对流层中上层温度场和对流层低层湿度场的影响较大，而且对对流层中下层风场也有一定改进。为了进一步加大对卫星资料的使用率，提高预报效果，RTTOV和CRTM都加入了水物质辐射及散射的模块，这使得在水物质的分布和天气系统配置合理并符合云物理基本特征的前提下，由于考虑了水物质的辐射效应而极大改善了卫星观测模拟效果(董佩明等，2008；2009；刘君等，2013)。

与卫星自大气层顶对大气进行观测相反，雷达是从大气底层对大气进行观测的，并能提供区域内高分辨率的大气三维风场和云水信息，因此合理恰当地使用雷达观测资料可以为卫星探测提供有益的补充。徐枝芳等(2002, 2003)在中尺度系统的模拟中同时同化卫星和雷达两种资料，改善了模式对MCS的预报能力。Milan等(2008)利用改善过的物理初始化方案在中尺度模式中同化雷达降雨预报产品和卫星的云顶高产品，得到了较好的临近预报效果。Sokol等(2009)将雷达反射率和卫星亮温转换为降雨率在COSMO NWP模式中进行同化，对1~3 h的降水预报起到了较好的效果。李媛等(2011)使用GRAPES 3Dvar同化系统对多普勒径向速度以及基本反射率因子的直接同化可以有效改进中尺度降水数值预报的效果，并且指出反射率直接同化的主要影响体现在湿度场上。此外，在WRF模式上雷达资料的直接同化模块也已经被开发使用(Xiao et al，2006；2007；2008)。

由于微波辐射能在实际大气中的衰减包含了气态水(即水汽)的吸收作用和液(固)态水(如云滴、雨滴、云冰等)的散射作用。因此若只单一考虑水汽的吸收作用，将辐射能的衰减全部认为是被水汽所吸收，往往会对大气中的水汽含量估计过高。解决这一问题的根本在于，卫星同化的观测算子是否能够较好地模拟大气对辐射能的散射作用。自RTTOV 8以来，微波的散射模块被开发出来并逐渐完善，但却从未被模式的同化系统正式引入。

相比RTTOV的晴空辐射模块(以下简称RTTOV_DIRECT)，微波散射模块(以下简称RTTOV_SCATT)的输入变量中增加了对大气水凝物含量的需求，为了使RTTOV_SCATT能够在WRF三维变分同化系统(以下简称WRF-3Dvar)中有效运行，本文在WRF-3Dvar的接口模块中设计编写了各种云水参数的计算程序，以支撑RTTOV_SCATT的计算，并在个例中加以检验。

1 大气水凝物含量的计算方案

由于RTTOV-SCATT的主要作用是计算水凝物对辐射的散射作用，因此WRF-3Dvar要调用RTTOV-SCATT，首先需要计算出大气中每层的水凝物含量，包括：云水混合比、云冰混合比、雨水流量，以及云量。这些物理量的计算方案设计如下，并以此编写FORTRAN程序加入到WRF-3Dvar的原代码中：

(1) 云水和云冰混合比：将通过WRF模式模拟得出的28层云水、云冰物理量从模式空间插值到RTTOV所要求的51层气压层上，其中超出模式上下边界的高度层上的值定为0，其余高度层上的值通过三线性插值计算得出。

(2) 雨水流量：可以定义为在某一气压层单位时间通过单位面积的雨水质量，单位为kg·m^-2·s^-1，根据定义，其计算公式可设计为以下形式：

$rain = rho \times {V_{\rm{T}}} \times (\frac{{{q_m}}}{{{q_m} + 1}})$

式中：rho为空气在该层的总密度，单位为:kg·m^-3；V_T为雨滴的末速度，单位:m·s^-1；q_rn为该层的雨水混合比，单位:kg·kg^-1。根据Sun等(1998)的计算公式求出V_T：

$\begin{array}{l} {V_{\rm{T}}} = 5.40 \times a \times q_{rn}^{0.125}\\ a = {(\frac{{{p_0}}}{{\bar p}})^{0.4}} \end{array}$

订正因子a由基本状态气压p和地面气压p₀求出。

(3) 云量：根据Smagorinsky(1960)的经验诊断方法，当格点内的平均相对湿度大于该层的“阈值”时，则该层有云存在，且云量可以用式计算得出：

$f(i) = {[\frac{{RH(i) - R{H_c}(i)}}{{1 - R{H_c}(i)}}]^n}$

f(i)∈[0, 1]，i代表气压层，RH_c(i)为i层的相对湿度阈值，RH(i)为i层的实际相对湿度，当RH(i) < RH_c(i)时，f(i)=0，n是经验常数，通常取2。参照Koch et al.(1997)确定的f-RH关系，对RH_c分高层( < 450 hPa)、中层(800~450 hPa)和底层(1000~800 hPa)分别取80%、60%和40%。

由于相对湿度不是模式背景场的物理量，因此需要用模式的已有物理量进行转换计算：

$\begin{array}{l} RH = \frac{q}{{{q_s}}}\\ {q_s} = 0.622\frac{E}{{p - 0.378E}} \end{array}$

式中：q和q_s分别为比湿和饱和比湿，E为饱和水汽压，根据Rogers等(1989)， $E = 6.112{{\rm{e}}^{\frac{{17.67{t_c}}}{{{t_c} + 243.5}}}}$ ，t_c为摄氏温度，这样相对湿度就转换为温度与气压的函数，而温度和气压都是模式中的常规变量，可直接调用。

2 个例的选取和资料预处理 2.1 个例的选取

文中选取的两个研究个例分别是2010年8月19日四川地区的暴雨个例(简称8·19个例)，研究区域选取为(26°~35°N、97°~109°E)和2011年7月24日北京的暴雨个例(简称7·24个例)，研究区域选取为(39°~42°N、115°~118°E)。

2.2 卫星观测资料介绍

卫星观测资料来自FY-3A卫星的微波湿度计(简称MWHS)，其中8·19个例的卫星资料取2010年8月19日09：17和10：58(北京时，下同)(如图 1a和1b)两次观测，将这两个时次资料合并(如图 1d)，以弥补单次观测的不全面。7·24个例的卫星资料取2011年7月24日21：21。并且在同化前需剔除其中以地面探测为主的通道1和通道2，并对其余3个通道做极值检验，剔除观测资料中亮温值大于350 K和小于120 K的点。

卫星资料的偏差订正通过WRF-3Dvar自带的VarBC变分法在资料同化之前进行偏差订正。由于两个时次最边缘的扫描点正是此次个例降水的中心位置，剔除之后对降水的模拟会产生很大的负效果，因此权衡利弊后选择保留，并对这两个时次观测中出现重合的边缘点，进行求平均，以此将边缘误差的影响尽量减小。

2.3 雷达资料及同化方法简介

本文两个个例所使用的雷达体扫基数据分别来自成都的单偏振多普勒天气雷达(简称成都雷达)和北京的双偏振多普勒天气雷达(简称北京雷达)，两种基数据的主要参数如表 1所示。表 2给出了卫星与雷达观测的时刻对照。

成都雷达的观测区域为(29.4°~32.0°N、102.6°~105.6°E)；北京雷达的观测区域为(39°~41.4°N、114.5°~118°E)。

雷达反射率因子的直接同化，通过将反射率因子Z与雨水混合比q_r的关系：z=43.1+17.5log(ρq_r) (Sun等，1997)作为雷达反射率因子的观测算子，经过WRF-3DVar进行变分同化计算，得出初始场中最优的雨水混合比。

3 试验方案的设计和结果分析 3.1 试验方案的设计

表 3显示了两个个例的同化时刻和时间窗的设定。对比试验的方案设计如表 4所示，其中：(1) 控制试验为不同化任何观测的预报模拟试验；(2) 对比试验1为使用RTTOV_DIRECT在晴空假设条件下对初始场进行MWHS资料同化，并做未来时刻的预报模拟(目前国内外最普遍的方法)；(3) 对比试验2为在同化计算中使用RTTOV_SCATT代替RTTOV_DIRECT做观测算子，对初始场进行MWHS资料同化，并做未来时刻的预报模拟；(4) 对比试验3为对初始场单独同化雷达反射率因子；(5) 对比试验4为将试验3中同化雷达资料后的分析场作为新的初始场，并在此基础上重复试验2。

通过将4组对比试验分别与控制试验比较可以对同化前后的效果进行检验，将试验1和试验2比较可以得出两种不同观测算子的优劣，将试验4与试验2比较可以对雷达资料同化前后，使用RTTOV_SCATT算子的效果进行检验。

3.2 同化结果分析 3.2.1 8·19个例同化模拟结果分析 3.2.1.1 水汽场的对比分析

本文以模式的预报场做为背景场，因此背景场中包含了模式通过模拟计算得出的大气垂直速度和各种水物质分布。图 2是同化时刻水汽混合比分布图。从图中可以看出，受高原地形影响，在高层400 hPa，四川西部的水汽含量要明显高于东部，达到3.5 g·kg^-1以上(图 2a)；而在低层800 hPa，只有海拔较低的东部地区水汽较为丰沛，其最大中心值达15 g·kg^-1以上，而西部高原地区由于地形高度要高于800 hPa的高度，因此没有数值(图 2b)。

根据微波湿度计的通道特性，本文主要对比分析了三个水汽通道的主要探测高度上(800、600和400 hPa)截取同化之后的水汽混合比增量，并用这三个高度分别代表大气的低层、中层和高层。

如图 3所示，其中图 3a、3b和3c分别是对比试验1相对控制试验在不同高度的水汽混合比增量图。从图中可以看出，WRF-3Dvar同化模块调用晴空观测算子RTTOV-DIRECT同化FY-3A微波湿度资料后，对初始场水汽含量的改善作用并不明显，仅有相对较少区域的水汽含量得到订正，且幅度在±0.1 g·kg^-1之间，增幅较小。由此可见，在本个例中以模式预报场为背景场，配合RTTOV晴空观测算子进行同化计算的结果并不理想，不能较好发挥卫星观测对初始场的订正和改善作用。

图 3a′、3b′和3c′是不同高度上的试验1相对控制试验模拟到预报时刻时的水汽场混合比增量图。与图 3a、3b和3c相比，预报场中水汽增量的变化略为明显。在中高层增量数值达到±1 g·kg^-1左右；在低层，被订正的格点数明显增多。虽然预报场中产生了较多的增量，但通常并不认为它们是同化对初始场改善作用的体现，不具备任何指导意义，因为这些增量格点的出现很大程度上是产生自模式变量场对初始场中新增扰动的适应调整过程。由此可见，目前卫星同化研究所普遍使用的同化方法(试验1方案)在某些情况下并不能十分有效地支持对个例的研究，这就需要通过其他方法增大对观测资料中大气信息的提取程度，如微波湿度探测亮温中的大气水凝物的分布状态。

图 4是使用散射算子RTTOV_SCATT进行同化的对比试验2相对控制试验在不同高度的水汽混合比增量图。从图 4a、4b和4c中可以看出，相对背景场，同化后的分析场有了一定的订正。在高层，四川中部存在一小负值区，而东部和高原地区上空有大面积的正增量，且增幅达到1.2 g·kg^-1以上。在中低层，四川东部地区水汽场有较大范围的负增量，其中负极大值出现在大气中层，达到-1 g·kg^-1左右。相比试验1中晴空算子的同化结果，大气中的水凝物信息在试验2中得到有效利用，辐射能在传输中的消耗就由晴空模式中单一的水汽吸收变成了吸收和散射共同作用的结果，更加符合大气的真实情况，减少了对水汽含量的不当高估。因此在辐射的强衰减区域，通过散射算子得到的同化结果中的水汽含量要少于通过晴空算子得到的结果。

从图 4a′、4b′和4c′中可以看出，对于此次个例中降水总量较大的四川中东部区域，在大气的高层，水汽含量明显减小，数值达-1.5 g·kg^-1以上，在大气的低层，水汽含量有明显的正增量，数值达1 g·kg^-1以上，而初始时刻却恰恰相反，表现为高层正增量、低层负增量。由此可见，试验2的模拟结果相比控制试验来说，高层偏干，低层偏湿，并且水汽在垂直方向上有向下传输的趋势。

图 5是同化时刻对比试验3相对控制试验的雨水混合比增量图。从图中可以看出，经过同化，只有在有雷达观测资料的区域内，初始场中的雨水变量才会被同化更新。分析发现：经试验3同化之后各层的雨水增量均为正增量，高层的大值区出现在远离观测中心的北部，数值较小；中层的大值区出现在观测的偏东区域，数值最高为0.08 g·kg^-1；在低层，大值区主要出现在观测区域的中心附近，且数值高达0.1 g·kg^-1。

试验4是以试验3输出的经过反射率因子订正的分析场作背景场，并进行卫星资料的同化计算。图 6是试验4相对控制试验的水汽增量分布图，从图 6a、6b和6c中可以看出，已经在试验3中被雷达反射率订正过的区域，再用散射算子进行微波亮温资料同化时，该区域内各高度层上的水汽混合比的减小幅度要明显大于其他区域。

与试验2单独使用卫星资料进行同化的模拟结果不同，由于同化单部雷达探测的有效半径为150 km，无法覆盖全部的模拟区域，因此试验4的预报场在800 hPa高度层上只有雷达探测范围内的地区有水汽的正增量，其余地区均为负增量(见图 6a′、6b′和6c′)，并不能很好地反映出试验4对四川东部降水集中区的低层水汽增量预报的改善作用。

3.2.1.2 降水场的对比分析

图 7a是20日08时24 h累计降水量。从图中可以看出，实况降水的最大量值中心有两个，分别位于(30°~31°N、103°~104°E)和(27°~28°N、104°~105°E)，其中最大降水量值达55 mm以上。

从4组试验预报降水量图的总体分布来看，4组降水落区的预报相差不大，但强中心位置均较实况都略有偏西。从累计降水量值来看，最大累计降水达到55 mm的只有控制试验，试验1和试验2的最大降水均在45 mm左右，但试验4只有35 mm。考虑到模拟的降水累计时间跨度比实况少两个小时，所以试验1和试验2模拟的降水量是最合理的。与其他三个试验相比，试验4中累计降水量大值区强度偏弱，而在其他地区，降水量基本一致。另外，结合图 6c′可以发现，图 7e中降水减少的区域与图 6c′水汽增加的区域几乎完全吻合，由此可见同化雷达反射率可以通过改变水汽间接影响到降水的预报，而试验1和试验2之间区别不大。

3.2.2 7·24个例同化模拟结果分析 3.2.2.1 水汽场的对比分析

相对8·19个例，7·24个例研究区域较小，在观测资料的覆盖程度上大大优于前者，这也在理论上保证了良好的同化效果。

对于使用晴空算子的试验1方案，在同化时刻(图 8a、8b和8c)，水汽增量的最大值出现在600 hPa，增量达3.2 g·kg^-1，且范围覆盖了北京中部和西北部的大部分地区。而在图 8a′、8b′和8c′中可以看出，水汽含量在预报场中形成了以北京中部水汽零增量带为界的水汽正增量区和水汽负增量区。在零增量带的西偏南区域，大气水汽增量在高层和低层均为正增量区，中层为负增量区，这表明了水汽在垂直方向上存在有辐散的增量；但在东半区域，水汽在高层和低层表现为负增量，中层为正增量，表征着水汽辐合增量的存在。

图 9是对比试验2相对控制试验的水汽混合比增量分布图。同化时刻(见图 9a、9b和9c)，试验2的水汽增量分布与试验1的结果比较一致，但增量强度略小于试验1。同样，11 h的预报场中水汽混合比增量(见图 9a′、9b′和9c′)也与试验1的结果比较一致，但正增量值要小于试验1。由此可见，在试验1中快速增长的水汽增量范围、增量数值在试验2中都得到了有效抑制，但是散射算子只能订正由晴空算子计算得到的变量数值，而对于变量的散度、涡度等特性只能与晴空算子保持一致。

图 10是试验4相对控制试验的水汽混合比增量分布图。从图 10a、10b和10c中可以看出，雷达资料的加入对水汽负增量的增强作用明显。与图 9相比，在同化时刻，试验4的水汽混合比增量在数值上明显偏低，正增量的区域明显缩小，负增量增多。这种现象同样表现在模拟结果中(图 10a′、10b′和10c′所示)。从试验4与试验3的对比可以看出，雷达资料的引入对水汽含量负增量的产生有较大影响，这主要表现为试验3中的正增量区域在试验4中表现为负增量或0值区，且正的峰值减小，负的峰值增大。

3.2.2.2 降水场的对比分析

从降水实况来看，在(40.2°~40.5°N、116.7°~117°E)存在12 h的累计降水中心，降水量高达120 mm。根据4组试验对降水的预报(图 11b、11c、11d和11e)，从落区来看，试验4的中心大值区与实况位置最为接近，试验1和试验2的位置略有偏北；从降水量的预报来看，试验1、试验2和试验4的最大降水量均达到了100 mm以上。就预报的总体效果而言，试验4的效果要明显好于其他三组试验，而试验1和试验2的预报效果相差不大，与8·19个例的所得结论一致。

4 结论

通过对四川2010年8·19暴雨过程和北京2011年7·24暴雨过程进行同化研究对比，结果表明：

(1) 在微波亮温的资料同化中，改用RTTOV_SCATT将大气中水凝物的信息纳入同化的计算中，可以有效地解决使用以往的方法所导致的对大气中水汽含量的高估现象。可见，相比传统的微波资料晴空假设方法，新观测算子的使用具有可行性和有效性。

(2) 使用RTTOV_SCATT算子引入水凝物对辐射的散射作用后，水汽增量范围和数值能够得到有效抑制，但对于变量的散度、涡度等特性只能与RTTOV_DIRECT的结果保持一致。

(3) 加入雷达反射率因子同化后，同化区域内的大气水凝物含量明显增加，从而在该区域内的卫星资料同化结果中，大气中水汽含量的减小量要明显大于其他区域。

(4) RTTOV_SCATT的使用对降水预报的同化改善并不明显。而在雷达同化的区域内，随着雷达观测资料覆盖范围的增大，对比试验4的同化方案对降水预报有显著的正效果。

本文通过引入简单的水汽凝结物计算方法，使得观测算子在直接同化FY-3微波湿度资料时能粗略反映出大气中的水凝物信息，从而改进了试验的模拟结果。若在业务中将此方法用于对暴雨的预报，可以结合多个极轨卫星组成的观测网，来弥补单个卫星在时间分辨率上的不足，以及使用多部多普勒天气雷达资料，来弥补单站雷达空间覆盖率的不足，使预报具有时效性。