引言

由于洋面上缺少足够有效的观测资料，台风初始化对于改进全球数值预报业务模式的涡旋质量仍然是至关重要的技术难题。尽管拥有大量丰富的卫星探测资料可供选择，但是因为受到云和降水的“污染”，这些资料的利用目前仍然处于技术探索阶段。为了能使T639全球数值天气业务模式对台风产生一个客观的、合理的分析和预报产品，国家气象中心应用了一套完整可行的台风初始化数值方案，主要有初始涡旋形成、涡旋重定位和涡旋调整三部分技术组成(瞿安祥等，2009a；2009b；麻素红等，2012)，其技术路线是：(1) 在台风刚形成的第一个预报时刻，在背景场中嵌入一个与观测接近的bogus涡旋环流(如果背景场台风涡旋环流非常弱或不存在的情况下)，然后进行资料同化分析-模式积分预报流程。(2) 经过第一个预报时刻的涡旋处理，第二个预报时刻的背景场(即第一个时刻的6 h预报场)中往往就会存在结构比较完整的台风涡旋环流系统，只是这个涡旋与当前时刻的观测不匹配，主要表现在强度弱(强)、中心位置不准确上，利用涡旋重定位、涡旋动力学公式对台风涡旋环流进行位置平移和强度调整，使之与预报员分析的观测接近，然后再进行资料同化分析-模式积分预报流程。(3) 在后续的循环滚动预报时刻，重复利用与第二个时刻同样的方法对背景场中的台风涡旋进行初始化，直至台风消亡。从2009年国家气象中心将T639全球数值模式系统投入业务运行至今，该台风初始化方案在提高模式的台风预报能力方面发挥了重大作用，也产生了可观的社会和经济效益(许映龙等，2011；曹晓岗等，2013；陈国民等，2013；2014；2015；张守峰等，2015)。但与此同时，几年来的业务运行也暴露出了该台风初始化方案在实际预报中存在的一些技术问题。由前面所述可知，该方案的关键在于台风生成所在第一个预报时刻的涡旋形成技术，初始涡旋的好坏直接影响到后续滚动时刻的台风涡旋质量乃至预报效果，尤其在缺乏有效台风观测资料的洋面上更是如此。但是，目前业务中应用的嵌入式bogus涡旋技术主要依据统计的经验模型建立，其环流结构与T639全球数值模式的动力学、物理学属性并不匹配，并且直接嵌入背景场也会引发与周围环境场协调性问题，并且这种不适应性会影响到后续循环同化分析时刻的涡旋自身属性及其与环境场之间的平衡性。

当台风发生且背景场没有涡旋存在的情况下，美国国家环境预报中心(NCEP)在其全球数值预报模式系统中采取的方案是：根据实际台风强度、尺寸构造一些伪风场资料进入资料同化系统，从而强迫同化分析出初始场中的涡旋环流系统(Kleist，2011)。日本气象厅(JMA)在其全球数值模式系统中也采取了类似的方案，通过构造海平面气压、风场两类伪台风观测要素来同化分析出台风涡旋结构。同样，英国气象局针对台风有效观测缺失的情况下，也采用了构造伪台风观测资料的方案(Heming，2009)。实际上，许多科学研究也表明，利用相关的统计模型构造一些代表台风特征的伪观测资料(Zou et al，2000；Xiao et al，2006)，依靠成熟的变分/集合卡尔曼资料同化技术(Chen et al，2007；Wu et al，2010)，将构造的伪观测资料与背景场有效的“同化融合”在一起，就可以大幅度提升分析场中初始涡旋的描述质量，从而可以显著提高台风的路径和强度预报水平。

事实上，目前国家气象中心T639全球变分资料同化系统的水平格点分辨率已经达到30 km左右，垂直层次有60层，已经具备分析表达出现实的台风涡旋环流结构(包括与台风相互作用的环境场)的能力。据与此，本文设计发展了一套升级版的台风初始化业务方案，即依据预报员提供的台风实时特征参数，利用背景环流信息，构建匹配的台风模型，然后从中抽取一些气压、风场要素作为台风的伪观测资料进入变分同化系统，分析同化出初始预报时刻的台风涡旋环流，而后续预报时刻的涡旋重定位、涡旋强度调整技术保持不变；基于国家气象中心的T639全球数值预报系统的3Dvar变分同化平台，本文完成了该技术方案的具体实现过程；然后应用新升级的台风初始化方案，对2011—2012年生成于西北太平洋的27个台风进行了数值试验，分析了该方案在实现台风涡旋环流结构方面的效果和能力，并对试验结果进行详细的比较，来展望新方案在台风数值预报业务上的应用前景。

1 台风观测资料构造 1.1 台风模型构建

台风发生时(本文中提到的台风是指底层中心附近最大平均风速大于10.8 m·s^-1的发生在西北太平洋热带气旋的统称)^*，预报员会根据实时搜集到的各种观测信息(卫星云图、雷达探测等)对台风的位置和强度进行诊断分析，并得到一系列可信度较高的估计参数值，包括中心定位、中心最低海平面气压、7级风圈半径(台风等级达到热带风暴才产生)、最大风速等；本文所采用的台风模型就是在这些相对质量较高的分析参数的基础上，综合考虑台风实际所处周围环境背景场构建而成。

^*中央气象台定义西北太平洋台风为6个等级，包括热带低压(底层中心最大平均风速10.8~17.1 m·s^-1)，热带风暴(底层中心最大平均风速17.2~24.4 m·s^-1)，强热带风暴(底层中心最大平均风速24.5~32.6 m·s^-1)，台风(底层中心最大平均风速32.7~41.4 m·s^-1)，强台风(底层中心最大平均风速41.5~50.9 m·s^-1)，超强台风(底层中心最大平均风速≥51 m·s^-1)。

首先应用绝对角动量守恒定理，确定计算实际台风环流的水平半径范围R_tc，计算公式如下：

$ {R_{{\rm{tc}}}}\left({1 + \frac{{2{V_{15}}}}{{f{R_{15}}}}} \right){R_{15}}, {V_{15}} = 15 $

(1)

式中，R₁₅为15 m·s^-1的台风切向风圈所在的半径(其值根据预报员分析提供的台风7级风半径所得)，f为科氏力参数，由台风中心所处纬度位置确定。在确定台风范围的基础上，应用Fujita(1952)经验公式可以计算模拟出轴对称结构的台风海平面气压分布模型场(Zou et al，2000)：

$ P\left(r \right) = {P_\infty } - \left({{P_\infty } - {P_{\rm{c}}}} \right){\left[ {1 + \frac{1}{2}{{\left({\frac{r}{{{R_{\rm{p}}}}}} \right)}^2}} \right]^{-\frac{1}{2}}}, $

其中

$ {P_\infty } = {P_{\rm{c}}} + \left({{P_{\rm{e}}} - {P_{\rm{c}}}} \right)/\left\{ {1 - {{\left[ {1 + \frac{1}{2}{{\left({\frac{{{R_{{\rm{tc}}}}}}{{{R_{\rm{p}}}}}} \right)}^2}} \right]}^{-\frac{1}{2}}}} \right\} $

(2)

式中，P(r)表示距离台风中心径向距离为r的海平面气压分布值；P_c为已知的预报员分析提供的台风中心最低海平面气压值，P_e为已知的台风边界区域上的环境海平面气压值(可通过背景场计算得出)，P_∞为距离台风中心无限远处的海平面气压值，R_p为台风区域内海平面气压径向梯度最大值所处半径。P_∞和R_p两个未知参数可经由台风15 m·s^-1切向风圈半径处所满足的梯度风平衡动力约束关系来确定(ρ为空气密度)，即：

$ {P_\infty } - {P_{\rm{c}}} = G\frac{{2R_{\rm{p}}^2}}{{{R_{15}}}}{\left[ {1 + \frac{1}{2}{{\left({\frac{{{R_{15}}}}{{{R_{\rm{p}}}}}} \right)}^2}} \right]^{\frac{3}{2}}}, $

其中

$ G = \rho \left({\frac{{V_{15}^2}}{{{R_{15}}}} + f{V_{15}}} \right) $

(3)

联合前面的公式，就可以计算出台风R_tc半径范围内任意点的海平面气压场分布值(注意，其值是呈轴对称分布的)。

接下来，可以通过如下的公式计算标准垂直气压层P的台风位势高度场分布：

$ \begin{array}{l} \Delta H\left({r, P} \right) = H\left({r, P} \right) - {H_{\rm{e}}}\left(P \right)\\ = \frac{{{R_{\rm{d}}}}}{g}\int_{P\left(r \right)}^P {T{\rm{dln}}P + \frac{{{R_{\rm{d}}}}}{g}\int_{{P_{\rm{e}}}}^P {{T_e}{\rm{dln}}P} } \\ = \frac{{{R_{\rm{d}}}}}{g}\int_{P\left(r \right)}^P {\left({T - {T_{\rm{e}}}} \right)} {\rm{dln}}P + \frac{{{R_{\rm{d}}}}}{g}\int_{{P_{\rm{e}}}}^{P\left(r \right)} {{T_{\rm{e}}}{\rm{dln}}P} \end{array} $

(4)

式中，ΔH(r, P)为待求解的垂直气压层P上、距离台风中心径向距离r的位势高度场距平值；P(r)为距离台风中心径向距离为r的海平面气压值，T_e为台风边界环境场的温度值，R_d为干空气气体常数，g为重力加速度。计算过程中采用了如下的台风模型假设(Frank，1977)：(1) 台风涡旋的暖心在250 hPa层次上达到最大值，并且在100 hPa层次上消失；(2) 台风位势高度场距平值在100 hPa层次上按比例衰减至环境场消失。

得到标准垂直层位势高度场分布后，可通过梯度风平衡关系式计算出相应气压层的切向风场V(r)分布，公式如下(g为重力加速度)：

$ \frac{{V\left(r \right)}}{r} + fV\left(r \right) = g\frac{{\partial H\left(r \right)}}{\partial } $

(5)

其中，为了考虑到台风边界层摩擦效应，将低层的切向风场进行两方面的经验调整：风速大小按一定比例进行垂直层次消减；风向台风内部进行一定角度的扭转，加强辐合效应。

1.2 台风观测资料构造

在构建完台风三维涡旋模型的基础上，本文设计抽取了海平面气压场和风场这两大要素作为类探空常规观测数据进入变分同化系统。要素的空间位置设计成围绕台风中心而进行的，即以台风中心为圆心，在半径r为100、200、300、400、600 km的圆圈上间隔选取4、6、8、8、12个点作为伪观测数据位置，并且相邻圈上的点在径向上错落分开，数据径向边界取值范围由前面计算的台风环流的水平半径R_tc确定。垂直方向上选取了地面、1000、925、850、800、700、600、500、400、300 hPa共10层气压层作为台风观测数据的分布层次。同时，台风中心位置的气压值也作为伪观测资料的一部分进入了同化系统。

为了能使伪观测资料在分析同化时对形成台风涡旋三维结构起积极促进的作用，本文对构造好的资料设计了初步的质量控制检查方案。针对海平面气压数据，选取的标准是其值小于同一位置的背景场值，否则剔除该数据；对于垂直空间各点的风场数据，则进行两方面条件的控制，一方面条件是该点的位势高度值要小于背景场位势高度的值，另外一方面条件是该点相对于台风中心的切向风速要大于背景场的切向风速值。经过这样的质量控制之后，就会确保：如果背景场没有涡旋环流的情况，进入分析同化的台风伪观测资料就会促进台风涡旋环流结构产生，如果背景场存在浅涡旋环流的情况下，台风伪观测资料也会增强背景场中的涡旋环流分量。

在资料误差确定方面，由于整个伪观测资料是基于预报员分析的台风特征参数基础上构建而成，这就似乎造成越靠近台风中心的伪观测数据具有相对较高的质量，而越远离台风中心，其相对质量较低；因此针对不同位置的海平面气压伪观测数据，本文给定了不同的观测误差，其值最大为1.6 hPa，最小为0.8 hPa，并且依据远离台风中心位置的距离而变化；而风场伪观测数据的误差确定，则参考采用云迹风这一类常规观测资料的误差给定值。同时，为了保证这些伪观测资料在变分同化时不被质量控制步骤剔除，需要对这些资料的总量控制参数做出适当放宽调整。

2 同化效果分析

基于国家气象中心T639模式的全球资料同化分析系统, 本文应用实现了新设计的基于同化伪观测资料的初始台风涡旋形成技术方案，并对该技术方案进行了台风个例的同化试验，来检验该技术方案在形成台风涡旋结构方面的实际效果。

以1203号台风玛娃(Mawar)为例，该台风于业务对应时间2012年6月01日06时(UTC，下同)在菲律宾以东大约500 km的西北太平洋洋面上生成。预报员分析估计的其特征参数为：热带风暴等级，中心位置：(14.9°N、125.5°E)；中心海平面气压：998 hPa；7级风半径：180 km。按照台风初始化方案的流程，系统首先判断出当前时刻的背景场中存在一个非常弱的涡旋系统(图 1)，搜寻所得中心海平面气压只有1004 hPa，比起预报员实际提供的观测数据998 hPa相差了6 hPa，同时，海平面气压场的闭合线表现得非常稀疏(图 1a)；从底层风场分析来看(图 1b和1c)，涡旋环流的形势也不是很明显，风速量级都很小。如果不对背景场中这样的涡旋进行初始化处理，势必会影响到分析场中台风涡旋的质量，也会影响到T639模式积分输出场中台风涡旋预报和其周围形势场预报的正确性，从而进一步影响到下一个预报时刻背景场中台风涡旋的质量。

应用本文新设计的台风观测资料构造技术，我们根据预报员分析的台风特征参数，计算出其水平径向范围约为420 km，则构造的台风伪观测数据分布(类似探空观测站)为27个站点(类似探空观测站)：1(台风中心点)+4(100 km径圈上)+6(200 km径圈上)+8(300 km径圈上)+ 8(400 km径圈上)；然后将这些伪数据作为一类有效台风资料投入变分同化系统。图 2显示的是变分同化系统吸收伪观测资料后在分析场中形成的台风涡旋环流情况。从海平面气压场来看(对比图 1a和图 2a)，比起原背景场的浅涡旋系统，同化吸收伪观测资料之后，分析场中台风涡旋的气压梯度明显加深，实际数据显示，此时的强度已经达到了998 hPa，与预报员的分析值十分接近。对比洋面风场(对比图 1b，1c和图 2b，2c)可以看出，背景场中的涡旋环流特征非常弱，但是经过同化分析后，台风的气旋性环流结构明显清晰，特别在台风涡旋中心区，同等值的风场闭合线得到了有效的扩大，相应得，风速的量值也明显加强。实际的数据显示来看，台风附近的风场分析增量普遍在6 m·s^-1左右，最高能达到12 m·s^-1以上。类似地，对比经过台风中心的垂直风场的纬向剖面图(对比图 1d和图 2d)也可以看出涡旋环流系统的台风结构特征明显清晰，强度和范围得到了有效的增强和加深。

从上面个例分析可以看出，本文新设计的初始台风涡旋形成技术方案十分有效，它通过构造包含台风特征的伪观测资料进入变分同化系统，从而迭代分析出一个三维结构比较完整的台风涡旋环流系统；这个台风涡旋不但保留了背景场中隐含的模式动力学和物理学属性，而且吸收了伪观测资料的台风涡旋特征，具有很不错的代表性；并且相比于业务应用的嵌入人造bogus涡旋而言，它大大减少了人为的主观经验技术对台风三维结构的影响；实际上，业务采用的bogus技术形成的涡旋环流结构是完全对称的，嵌入背景场时并没有保留背景场任何环流信息(图略)，因此并不会形成如图 2所示的台风涡旋非对称环流结构信息，而台风的非对程结构信息，特别是风场信息对某些台风的路径和强度预报至关重要(张玲等，2014)。同时，从上面分析也可以看出，通过同化分析产生的涡旋不仅与和预报员分析的观测数据相匹配，而且与周围的环境场也比较协调(得益于变分同化中背景场误差协方差相关属性及其动力学约束)，大大优于直接嵌入人造bogus涡旋而引起的边界和环境场适应性问题(尽管在嵌入bogus涡旋的时候可以选取各种各样的边界耦合平滑方案，但并不能从根本上解决与周围环境场的协调性问题)。

实际上，在循环滚动预报试验(见下文)中本文也发现，无论背景场是否存在弱涡旋环流，通过同化本文设计的台风伪观测资料，总可以在分析场中形成不错质量的台风涡旋环流，并且这样形成的涡旋环流会明显改善台风的预报路径和预报性能。

3 平行试验

应用新设计的技术方案，本文将台风初始化过程融入到了国家气象中心T639全球资料分析同化-预报循环系统中，并建立了相应的台风预报业务流程，如图 3所示，从图中可以看出，当有某个台风发生时，相关的台风初始化作业流程就会启动，针对不同编报时刻执行涡旋形成、涡旋重定位或涡旋调整过程；而如果多个台风并发时，就会执行不同编号台风应用不同的初始化过程。

利用新建立的数值预报系统，本文对2011—2012年发生在西北太平洋的27个台风(编号1108~1118，1203~1218) 进行了平行数值试验。试验的标准严格按照实时业务预报的需求，资料也取自实时业务运行的数据。每个台风样本试验选取的时段从预报员初始编报开始，一直到消亡停编为止，在此生命史期间，每天共进行四个时刻的循环滚动预报，分别对应于业务上的00、06、12、18 UTC的资料分析时刻，每个时刻输出120 h时效的台风路径预报产品。

3.1 试验样本

本文试验选取的27个西北太平洋台风样本(如图 4所示，来源于中央气象台实时编号数据)中，有10个登陆我国，7个热带气旋虽未在我国登陆，但对我国近海海域有影响。这些登陆或影响沿海陆地的台风一方面给我国南方大部地区带来了丰沛的降水，缓解了江南、华南等地的干旱和酷热天气，使农业生产获益，另外一方面，热带风暴带来的狂风、暴雨及大潮也给部分地区造成了严重损失。因此本文选取的这27个台风样本具有较典型的代表意义。

3.2 结果分析

为了对比检验新设计的台风涡旋形成方案对这27个台风的路径预报影响，本文也应用了目前业务中使用的嵌入人造bogus台风涡旋形成方案进行了平行试验，即试验分两个方案进行，一个为控制试验，即初始涡旋形成应用业务中的嵌入bogus涡旋技术；另一个为新方案预报路径试验，应用前面所述的基于伪观测资料同化形成初始台风涡旋技术；在背景场台风涡旋结构已经存在并达到一定强度后，两个试验都采用同样的涡旋重定位和强度调整技术。试验结束后本文将两种方案的预报结果进行对比，来检验新方案在相比于业务方案的优劣好坏。

图 5显示了1206、1208这两个西行登陆台风在各自生命史中，分别应用业务方案和新方案预报的台风路径和实际观测路径的对比情况。首先看1206号台风的预报路径，从图 5a中可以看出，在该台风生成初期的第一个预报时刻，应用业务方案的控制试验产生了一条明显北偏的预报路径，并且在后期预报登陆穿过了台湾岛，这与台风的实际移动观测路径产生了较大的差异；实际上，如果用这个预报结果给预报员做指导性参考产品，就会给预报员造成错误导向。但是应用新台风涡旋形成技术(图 5b)后，尽管预报的路径也稍许北偏，但与控制试验相比已大幅度改善，而且台风预报的主体移动方向已与后期的实际观测路线相差不多。再看1208号台风的预报情况，从图 5c中可以看出，相比与台风实际观测路径，控制试验在台风生命史的第一个预报时刻产生了路径移动明显转折北上的错误趋势，并且这种错误趋势直接影响到了后续的三个预报时刻。但是反观应用新方案(图 5d)，则会在很大程度上纠正改善这种错误趋势，只是在第一个预报时刻的后期路径稍微转折北上，但是幅度已减缓了很多，而且后续的两三个预报时刻已逐渐调整趋向于台风实际移动路线方向。

接下来分析另外两种类型台风：北上行台风1207号和西北行登陆台风1209号的预报路径对比情况，如图 6所示。从图 6a中可以看出，应用业务方案后，1207号台风在其生成初期的前三个预报时刻都产生了明显左偏的预报路径，但是新方案则大大纠正了这种左偏的趋势(图 6b)，无论是第一个预报时刻，还是后续的二、三次的预报时刻，都产生了与实况比较一致的路径移向。再看1209号台风的情况，业务方案对这个台风的预报算是比较失败的(图 6c)，整个生命史的前3 d(约十几个预报时刻)产生的预报路径明显东北移向，与实际北偏西的路径大相径庭。而应用新方案后，只是在前三个预报时刻(不到1 d)产生了相对比较右偏的路径，而且偏移的幅度已大大减轻。尽管整个生命史的预报情况还不尽如意，但是相比于业务方案的异常预报路径，已经大幅调整改进了许多。

纵观前面的几个台风分析，应用基于资料同化的涡旋形成技术方案后，T639模式对台风的预报能力得到了一定幅度的提升，相比于业务的嵌入bogus涡旋技术方案，路径预报的改善主要集中在台风生命史的前几个预报时刻，持续时间因台风个例而异，从一个到十几个预报时刻不等。实际上，这种类似结论对试验中选取的大部分台风都有效，限于篇幅，这里没有做一一展示。这就说明：新设计的利用资料同化技术形成的台风涡旋大大改善台风初期时刻产生的异常偏折的预报路径，并且间接提高了后续滚动时刻产生的预报路径质量。

实际的检验数据也显示了这点，本文对选取的27个台风前后共进行672次连续滚动的预报试验，图 7为应用两种方案后，全球数值模式T639在120 h时效内预报的台风路径平均误差的对比分析。从图 7中可以看出，应用新台风涡旋形成方案后，台风平均路径预报误差普遍降低，在36~120 h的预报中，路径误差有了几千米到几十千米的下降；计算表明，在48~120 h的预报中，平均路径误差普遍降低了3%~15%，其中84~120 h的下降百分比都在10%以上。这说明，在间隔为6 h的120 h预报时效内，基于新方案预报的台风路径平均误差整体上都要小于基于业务方案的路径误差。这个统计结果还是比较让人欢欣鼓舞的。

同时，本文试验研究中也发现，新方案也有一些劣势表现。如图 7显示，在短时效的36 h之内，平均预报路径误差有了小幅增长，造成这方面的主要原因是：部分台风样本在形成的初始预报时刻，往往在背景场中存在一个中心偏离——强度非常弱的涡旋环流，并且这种涡旋环流经过伪观测资料的有效同化后其中心位置并不能得到有效纠正，与实际观测仍然存在一定偏差，而且这种偏差会在模式积分时延续一段时间，从而短时效的路径误差有所增长。另外，本文在试验中也发现，新方案对那些初始时刻比较强，范围比较小的台风表现的也差强人意，因为这时利用伪观测资料同化产生的涡旋往往范围较大，强度过弱(这与变分同化技术本身特性相关，分析同化产生的涡旋一定是背景场和伪观测资料折中的结果)，与实际台风环流相差较远，从而造成预报失败。

4 结论和讨论

基于T639全球资料分析同化-预报循环系统，本文发展了一套升级版的台风初始化业务方案，与旧方案的主要区别在于升级方案利用同化台风伪观测资料来形成初始涡旋而不是直接嵌入人造bogus涡旋；这在很大程度上减少了人为主观因素对台风初始结构的影响，而更多地是依靠资料同化的寻优迭代过程来协调约束产生三维空间的涡旋环流；该过程巧妙利用背景场误差协方差相关属性，不但在分析场中形成一个动力属性相对比较平衡的台风涡旋环流结构，而且改善了周围环境场与台风的相互协调性，使得涡旋环流和环境场更加溶为一体，大大减少了人为主观因素对台风涡旋结构及其环境场的影响。

应用新方案，本文对生成于西北太平洋2011—2012年的27个台风进行了数值试验，初步分析表明，新设计的方案在形成台风初始涡旋环流结构方面有着不错的表现，通过同化分析形成的台风涡旋环流系统不但保留了背景场中隐含的模式动力学和物理学属性，而且吸收了伪观测资料的台风涡旋特征，具有很不错的代表性。批量的试验表明，相比于业务使用的嵌入bogus涡旋形成方案而言，新方案台风路径平均预报误差有了不错幅度的降低，并且这种降低对早期阶段的预报效果改善更加明显。基于这个试验结果，国家气象中心已于2014年5月将升级的台风初始化方案投入了业务运行。经过近一年实际应用显示，新系统在向预报员提供指导性的台风路径数值预报产品方面起到了良好效果。