1 引　言

数值模式分辨率的提高以及同化循环周期的缩短对参加同化的观测数据提出更高要求（陈朝平等，2010），地面自动气象站（简称自动站）观测具备观测要素均为模式变量的优势，观测量可以直接同化，而且自动站观测资料兼备高时、空分辨率特征，满足数值预报的快速同化发展需求，使其成为快速循环同化系统的重要资料来源，其在数值模式中所起的作用也愈加显著。

自动站观测时间频率为5 min，数据量丰富，然而受限于数据质量（Shafer，et al，2000；陶士伟等，2007）和自动站同化技术的发展（徐枝芳等，2006），自动站观测资料在数值模式中的同化应用进展较为缓慢（徐枝芳等，2006）。近年自动站观测数据的质量在逐步提升，如以往因区域级自动站缺少气压和相对湿度的观测而导致同化吸收的湿度信息较少的情况已经得到明显改善。观测要素量的增加配合科学合理的观测资料质量控制方法使得自动站观测数据质量问题已经基本得到解决（任芝花等，2015）。目前中国2800多个国家级自动站均为全要素观测站，且数据经过严格质量控制，可以为数值模式同化提供良好的观测数据支撑。因此，优化地面观测同化方案是进一步提高自动站观测数据同化应用效果的主要途径。

以往中外学者对地面观测资料同化方案进行了很多研究和探索，但近年相关工作则开展较少。Ruggiero等（1996）针对站点高度所代表的实际地形高度与模式地形高度不匹配的情况设计了相应的地面资料同化方法，若站点高度与插值到该点的模式地形高度差超过100 m，该站点观测资料被剔除；若站点高度高于模式地形高度，该站点观测数据作为高空观测数据加以同化应用；若站点高度低于模式地形高度，则将观测值由站点高度订正到模式地形高度。上述方案考虑到地形落差对于地面资料同化的影响，但是观测数据剔除率高，且各要素的订正方案需要调整优化。郭永润等（2002）设计了不同的地面资料同化方法，即假定观测站所在位置的地形高度与模式地形高度无差异，遂选用相似理论建立2 m气温、2 m湿度和10 m风场的观测算子以及伴随模式，再利用极小化算法将地面气压观测值折算到模式面上，Barker等（2003，2004）将该方法引入到中尺度数值模式MM5（Mesoscale Model 5）中。采用该方法的优势是可以同化所有的地面观测数据，但没有考虑地形高度误差因素对同化造成的影响。其他科研工作者在以上两种同化方法的基础上进行过相关探索工作，如Devenyi等（2003））和Benjamin等（2004）采取计算气象要素局地递减率的方法，将地面观测要素值由实际地形高度订正到模式地形高度上，一定程度解决了地面资料同化中地形差异的问题，但是局地递减率的不确定性导致订正后的观测数值存在一定误差；Lee等（2005）以莫宁-奥布霍夫相似理论为基础设计了近地面综合观测算子。徐枝芳等（2007a）在郭永润方法的基础上将地形高度的代表性误差加入到观测误差中，一定程度上解决了郭永润方法中地形差异的问题，但是在观测误差中增加地形代表性误差只改变了资料同化时观测数据的权重，并没有完全解决地形高度差所引起的观测误差问题。邵长亮等（2015）采用集合均方根滤波同化10 m风场、2 m露点温度、2 m气温和地表气压，可以改善降水预报，但除去地表气压外，其他要素同化时未考虑高度差异。

欧洲中心主要采用最优插值方法同化地面2 m气温、2 m比湿、土壤温度、雪深和积雪温度等观测资料，综合考虑背景场和观测误差的特征来给定插值权重，而对于土壤湿度观测，采用简化的卡尔曼滤波方法同化该资料（Drusch，2007）。目前欧洲中心同化使用的地面自动站资料主要来源于全球电信系统（Global Telecommunications System，GTS），站点10000个左右；英国主要采用混合四维变分同化系统同化常规观测资料（Dixon，et al，2009），同化时间窗为4 h，并通过背景场和观测的偏差（OMB）统计分析设定黑名单，去除观测质量较差的观测站数据，再通过常规观测资料对预报误差敏感性分析即FSO（Forecast Sensitivity to Observation）试验来分析评估常规观测对于预报结果的影响。

目前WRF-DA同化系统中包含两种地面资料同化方案。方案1基于Ruggiero等（1996）的研究成果搭建而成，方案2则根据Guo等(2002)的研究成果建立，二者最大的区别是方案1中探究了地形高度差异对同化效果的影响，方案2中缺少该因素的考量。中国地形复杂多样，模式地形高度与实际地形高度的差甚至超过千米，而地形高度差异对于地面观测资料同化的影响非常大。因此，同化时应考虑地形高度差异的影响（Lazarus，et al，2002；Deng，et al，2007）；不同地面观测站所处位置的地表类型和地形高度差异较大，而不同地表结构和起伏的地形对边界层低层动量、水汽的垂直输送起到重要作用（余兴等，1997）。因此，在设计地面资料同化方案时应该探究地形高度对近地面层气象要素的影响（徐枝芳等，2007a，2007b）。

2 数值模式和观测资料

选取北京城市气象研究院研发的快速更新多尺度分析及预报系统短期预报子系统RMAPS-ST（Rapid-refresh Multi-scale Analysis and Prediction System-Short Term）进行数值模拟试验。模式采取两重嵌套设计，外区分辨率为9 km，覆盖中国范围，内区分辨率为3 km，覆盖华北地区，两重嵌套的垂直层数均为51层，采用ECMWF（European Centre for Medium-Range Weather Forecasting）细网格资料作为模式的驱动场。RMAPS-ST的同化模块是基于WRF-DA搭建而成，其中地面资料同化方案沿用了WRF-DA中的同化方案1。

文中着重介绍改进同化方案对于模式预报性能的影响，重点分析区域为模式的内区。数值模拟试验同化的观测数据类型包括经过严格质量控制的国家级地面自动站观测数据（SYNOP）、探空数据（SOUND）、飞机报数据（AMDAR）以及地基全球定位系统（GPS）天顶总延迟数据（GPSZTD）。4种观测数据在模式内区的分布如图1所示。不同的同化时次，各类数据的到报量差别较小，同化的数据量比较稳定。由各类观测数据的空间分布情况可以看出，探空数据分布比较均匀；全球定位系统天顶总延迟数据和飞机报数据则主要集中在模拟的中部和东部地区，西部比较稀疏；而地面观测空间分辨率高，自动站点基本全覆盖，但同样西北部地区的数据覆盖率较低。由图1可以看出模式区域内有约10%的地面自动站所处位置的观测地形高度高于模式地形高度（冷色散点），而90%的自动站观测高度均位于相应模式地形高度以下。就站点所处位置的观测地形与模式地形的高度差而言，约94%的站点位置高度差在100 m范围内，仅有不到1%的站点高度差大于500 m。整体而言该区域内地面自动站的观测地形高度与所处位置的模式地形高度差异并不大，可以同化吸收的有效数据量丰富。需注意的是，更新后的同化方案效果只作用于地面自动站观测资料。

3 WRF-DA同化系统中地面观测资料同化方案的改进 3.1 WRF-DA中地面观测资料同化方案

在WRF-DA资料同化系统中，地面观测资料同化方案根据实际地形高度与模式地形高度的差值判断地面观测站的数据是否符合同化要求，即地面观测站的高度为 $ {h}_{{\rm{o}}} $ ，对应点的模式最低层高度为 $ {h}_{{\rm{m}}} $ ，若 $ {h}_{{\rm{o}}} $ 与 $ {h}_{{\rm{m}}} $ 的差超过100 m，则该观测站数据被剔除；其余观测站的数据根据 $ {h}_{{\rm{o}}} $ 与 $ {h}_{{\rm{m}}} $ 的大小关系分为两类，并采用不同的方案进行数据处理。采取探空资料同化方案同化 $ {h}_{{\rm{o}}} $ 大于 $ {h}_{{\rm{m}}} $ 的地面观测数据，用地面资料同化方案同化 $ {h}_{{\rm{o}}} $ 小于 $ {h}_{{\rm{m}}} $ 的数据。地面资料同化方案旨在解决地形高度不匹配所导致的地面观测数据需要重新订正的问题，即利用不同方法将各个地面观测要素值由观测高度订正到模式最低层高度，得到订正后的“观测数据”用于计算各要素的更新向量。地面资料同化流程如图2所示，中国35%的国家级自动站观测数据因为不满足地形高度差小于100 m的条件而被剔除，7%的观测数据进入探空资料同化流程中，最终满足实际地形高度低于模式地形高度且地形高度差小于100 m条件的自动站观测数据仅占总数据的58%。订正后观测要素值的准确性直接影响地面资料的同化效果，不同要素的订正方法如图2所示，需要指出的是湿度没有普适的垂直变化规律，因此同化方案中不包含湿度订正。地面气压的同化方案可靠性高，因此文中主要优化地面风和地面温度的同化方案。

3.1.1 WRF-DA中地面风的同化方案

风的同化方法是基于莫宁-奥布霍夫（M-O）相似理论，将观测高度上的风订正到模式最低层高度上。由相似理论可以推知地面以上高度 $ {h}_{1} $ 和 $ {h}_{2} $ 处的风速 $ { U}_{1} $ 及 $ { U}_{2} $ ，公式如下（Barthelmie，2001）

$\left\{ \begin{array}{l} {U_1} = \dfrac{{{U_*}}}{k} \cdot \ln \left({\dfrac{{{h_1}}}{{{z_0}}}} \right)\\ {U_2} = \dfrac{{{U_*}}}{k} \cdot \ln \left({\dfrac{{{h_2}}}{{{z_0}}}} \right) \end{array}\right. $

(1)

式中， $ k $ 为常数， $ {z}_{0} $ 为地面粗糙度， $ {U}_{*} $ 为地面摩擦速度。根据式（1）可以得出不同高度 $ {h}_{1} $ 和 $ {h}_{2} $ 上的风速 $ {U}_{1}{\text{与}}{U}_{2} $ 有以下关系

$ {U}_{2}={U}_{1} \cdot \ln \left(\frac{{h}_{2}}{{z}_{0}}\right)\bigg/\ln \left(\frac{{h}_{1}}{{z}_{0}}\right) $

(2)

基于式（2）将地面风速观测值换算到模式最低层高度时，地面观测站的高度相当于0，数学公式中0没有对数值，因此需要将地面风速观测值（ $ {U}_{\rm{sfc}} $ ）和某一高度（ $ {h}_x $ ）的风速观测值（ $ {U}_{x} $ ）建立相关，从而使用式（2）对不同高度间的风速进行转换。经过大量观测数据的分析（Grell，et al，1994），地面以上40 m高度的风速（ $ {U}_{40} $ ）和地面风速（ $ {U}_{\rm{sfc}} $ ）之间可以建立如下关系

$ {U_{40}} = {U_{{\rm{sfc}}}} {\text •} \alpha $

(3)

式中， $\alpha =\left\{ \begin{array}{l} { 1.169 + 0.315{z_0}\quad\;\;\;{{z_0} {\text{≥}} 0.2} }\\{1.000 + 0.320{z_0}^{0.2}\quad {{z_0} {\text{＜}} 0.2} }\end{array} \right.$

通过式（2）和（3）可以将地面风速（ $ {U}_{\rm{sfc}} $ ）订正到模式最低层高度 $({ h}_{{\rm{m}}})$ 上，订正后的风速为 $ {U}_{\rm{m}} $ ，订正系数为C_correc，其中地形高度差 $ {Z}_{{\rm{a}}}={h}_{{\rm{m}}}-{h}_{{\rm{o}}} $ ， $ {h}_{{\rm{o}}} $ 为观测站所在的地形高度（或称地面观测站高度）

$ \left\{\begin{array}{l} {U_{\rm{m}}} = {U_{{\rm{sfc}}}} {\text •} C_{\rm{correc}}\\ C_{\rm{correc}} = \alpha {\text •} \beta \\ \beta = \dfrac{{\ln \left({{Z_{\rm{a}}}/{z_0}} \right)}}{{\ln (40/{z_0})}} \end{array}\right. $

(4)

式（4）中将地面风速订正到模式最低层高度时，计算 $ \beta $ 用 $ {Z}_{{\rm{a}}} $ 而不是 $ {h}_{{\rm{m}}} $ 。在式（4）的基础上进行调整并加入限定条件后建立地面风的原同化方案，原方案中增加的限定条件如下

$\left\{ \begin{array}{l} Ri_b {\text{＜}} 0\\ \dfrac{{{Z_{\rm{a}}}}}{L} {\text{＞}} 1.5 \end{array}\right. $

(5)

式中， $ Ri_b $ 为总体理查森数， $ L $ 为莫宁-奥布霍夫（M-O）长度。原同化方案限定只有观测站点所在位置的大气层结处于非常不稳定状态下（即同时满足式（5）中的2个条件），才依据式（4）将该站点的风速观测订正到模式面；否则设定站点观测风的订正系数（C_correc）为1，即相当于不对地面风进行高度差订正。除此以外，原方案计算 $ \beta $ 时将变量 $ {z}_{0} $ 设定为常数0.05和1.0，从而计算得到2个 $ \beta $ 值，再取二者平均代替式（4）中的 $ \beta $ 。

3.1.2 WRF-DA中地面温度的同化方案

温度的同化方法主要基于位温垂直变化率做订正。原同化方案中采用最接近高于模式地形100（ ${\rm{k}}_{100{\rm{hPa}}}$ ）和200 hPa（ ${\rm{k}}_{200{\rm{hPa}}}$ ）的模式层高度（ $ h $ ）及位温（ $ \theta $ ）来计算位温垂直变化率（θ_lapes）（Ruggiero，et al，1996）

$ \theta_{\rm{lapes}}=\frac{{\theta }_{{\rm{k}}_{100{\rm{hPa}}}}-{\theta }_{{\rm{k}}_{200{\rm{hPa}}}}}{{h}_{{\rm{k}}_{100{\rm{hPa}}}}-{h}_{{\rm{k}}_{200{\rm{hPa}}}}} $

(6)

利用位温垂直变化率计算得到位于观测高度和模式最低层高度上的位温，分别为 $ {\theta }_{\rm{os}} $ 和 $ {\theta }_{\rm{ms}} $

$ \left\{\begin{array}{l} {\theta _{\rm{os}}} = {\theta _{{{\rm{k}}_{200{\rm{hPa}}}}}} + \theta_{\rm{lapes}}{\text •}({h_{\rm{o}}} - {h_{{{\rm{k}}_{200{\rm{hPa}}}}}})\\ {\theta _{\rm{ms}}} = {\theta _{{{\rm{k}}_{200{\rm{hPa}}}}}} + \theta_{\rm{lapes}}{\text •}({h_{\rm{m}}} - {h_{{{\rm{k}}_{200{\rm{hPa}}}}}}) \end{array} \right.$

(7)

通过观测温度和气压可获得真实位温（ $ {\theta }_{\rm{obs}} $ ），并得到 $ {\theta }_{\rm{obs}} $ 和 $ {\theta }_{\rm{os}} $ 之间的差值（ $ {\Delta \theta }_{{\rm{o}}} $ ），将 $ {\Delta \theta }_{{\rm{o}}} $ 由观测高度线性插值到模式最低层高度后再赋值给 $ {\theta }_{\rm{ms}} $ ，即算出订正后的位温“观测”（ $ {\theta }_{\rm{sfc}} $ ）

$\left\{\begin{array}{l} {\theta }_{\rm{sfc}}={\theta }_{\rm{ms}}+{K{\text •}\Delta \theta }_{{\rm{o}}}\\ K=\dfrac{{h}_{{\rm{k}}_{200{\rm{hPa}}}}-{h}_{{\rm{m}}}}{{h}_{{\rm{k}}_{200{\rm{hPa}}}}-{h}_{{\rm{o}}}}\end{array}\right. $

(8)

将位温（ $ {\theta }_{\rm{sfc}} $ ）换算为温度（ $ {T}_{\rm{sfc}} $ ）， $ {T}_{\rm{sfc}} $ 即为订正后的温度值。

3.2 改进WRF-DA中地面资料同化方案

大气层结包括稳定、中性和不稳定状态，只有中性大气层结状态下的风速垂直变化情况（Stull，1988）符合式（1），大气处于稳定状态或不稳定状态下，原同化方案对于风速的订正不准确；同样，大气处于稳定状态或不稳定状态下，位温在近地面的垂直变化特征与在边界层高度外的垂直变化特征有一定的差异（Stull，1988），原同化方案通过高于模式地面200和100 hPa的模式层位温计算出的位温垂直变化率来订正地面温度，使用的模式层高度较高，因此原方案计算出的位温变化率不太适用于近地面层。

采用原地面风的同化方案进行试验，统计2017年7月1—5日的平均结果发现，00时（世界时，下同）进入同化系统的1400多个有效站点中仅有11个站点的大气层结状态满足原方案所定义的不稳定条件，进而风的订正系数大于1，其余99%以上站点观测风的订正系数均为1，由于原方案中对应点的 $ {h}_{{\rm{m}}} $ 一定高于 $ {h}_{{\rm{o}}} $ ，订正系数为1表示使用原方案后该站点所在位置的风速和风向并没有随高度的升高而发生变化。同样，选取地面温度的原同化方案进行试验，5 d进入同化系统的有效温度数据约8000个，基本上计算出的所有站点的位温垂直变化率均小于0.75 K/（100 m）。通过同时段的探空观测数据发现近地面温度随高度变化迅速，换算为位温后，其垂直变化率甚至可超过3.0 K/（100 m）。可见需要对风速和温度的原同化方案进行调整。

优化地面风和温度的原同化方案后得到更新方案，需要指出的是原方案和更新方案中采用的大气层结状态和位温垂直变化率均由数值模式背景场所决定。更新方案同样将地面观测数据订正到测站所在位置对应的模式最低层高度上，但采用不同方案所得到的风速及温度的订正值存在差异，为了检验改进后方案的效果，下文将同时次与模式最低层相同高度的秒探空观测数据作为真值来验证更新方案的可靠性。

3.2.1 改进WRF-DA中地面风的同化方案

风的原方案中采用的式（1）对中性层结大气下近地面风的垂直变化特征表述准确，而大气处于稳定或不稳定状态下，式（1）需要进行调整。另外，地面粗糙度（ $ {z}_{0} $ ）是一个随下垫面和季节变化的量，设置为常数会影响风速订正的效果（Barthelmie，et al，1996）。对风的同化方案进行以下调整和优化得到更新方案（或称新同化方案）：去除原方案中式（5）的判定条件，改为在式（1）中加入大气稳定度的判定项 $ \varPsi \left(\tau \right) $ ，将式（1）变更为式（9）（Stull，1988）

$ \left\{\begin{array}{l} {U}_{1}=\dfrac{{U}_{*}}{k} {\text •} \left(\ln \left(\dfrac{{h}_{1}}{{z}_{0}}\right)-\varPsi \left(\tau \right)\right)\\ {U}_{2}=\dfrac{{U}_{*}}{k} {\text •} \left(\ln \left(\dfrac{{h}_{2}}{{z}_{0}}\right)-\varPsi \left(\tau \right)\right) \end{array}\right. $

(9)

式中， $ \varPsi \left(\tau \right) $ 为 $ \tau $ 的函数，而 $ \tau ={Z}_{{\rm{a}}}/L $ ， $ \varPsi \left(\tau \right) $ 的计算公式如下

$\varPsi \left(\tau \right)=\Biggr {\{} \begin{array}{l} 2\ln \left(\displaystyle\frac{1+\tau }{2}\right)+\ln \left(\displaystyle\frac{1+{\tau }^{2}}{2}\right)-{\rm{arctan}}\left(\tau +\displaystyle\frac{\rm{\pi }}{2}\right)\quad\tau \;{\text{＜}}0 \\ -5\tau\quad \tau {\text{≥}} 0 \end{array} $

(10)

大气处于稳定（ $ \tau $ ＞0）或不稳定（ $ \tau $ ＜0）的条件下，近地面风的变化符合式（9），而大气处于中性层结（ $ \tau $ =0）下推算出的式（1）即为式（9）的特殊情况；其次，在计算订正系数时，重新引入变量 $ {z}_{0} $ 而不是将其设置为固定值。更新方案得到新的订正后的风速（ $ {U_{{\rm{m}}\_{\rm{new}}}} $ ）

$ \left\{\begin{array}{l} {U}_{{\rm{m}}\_{\rm{new}}}={U}_{\rm{sfc}} {\text •} C_{\rm{correc}} \\ C_{\rm{correc}}=\alpha {\text •} \beta \\ \beta =\dfrac{\ln \left({Z}_{{\rm{a}}}/{z}_{0}\right)-\varPsi \left(\tau \right)}{{\ln}(40/{z}_{0})-\varPsi \left(\tau \right)} \end{array}\right. $

(11)

3.2.2 改进WRF-DA中地面温度的同化方案

准确表述近地面位温垂直变化规律是地面温度同化的核心（Benjamin，et al，2004）。通过位温垂直廓线可以明确看到，在边界层低层和边界层以上两个不同高度范围内位温垂直变化规律差异明显（Stull，1988）。原方案计算位温垂直变化率所用的模式层高度 ${h}_{{\rm{k}}_{200{\rm{hPa}}}}$ 和 ${h}_{{\rm{k}}_{100{\rm{hPa}}}}$ 在每天的大多数时段均超过模式边界层高度，尤其是夜间边界层高度较低，因此原方案推导出的位温变化率不能准确地表达位温在近地层的垂直变化规律。地面温度订正时应选取更合理的层次来计算位温垂直变化率，综合考虑地形对边界层低层气象要素影响大的因素（余兴等，1997）及边界层高度的日变化特征因素，选取模式边界层高度和近地层高度（边界层高度的十分之一）（Stull，1988）及与其相同位置的模式位温来计算垂直变化率进而订正地面温度。通过 $ Ri_{\rm{b}} $ 计算模式边界层高度，公式如下（张碧辉等，2012）

$ Ri_{\rm{b}}=\frac{g\left(\Delta z{\Delta \theta }_{{\rm{o}}}\right)}{\Delta U{\Delta \theta }_{v{\rm{o}}}} $

(12)

式中， $ g $ 是重力加速度， $ \Delta z $ 是模式K层与模式最低层的垂直距离， ${\Delta \theta }_{{\rm{o}}}$ 、 $ \Delta U $ 和 ${\Delta \theta }_{{\rm{vo}}}$ 分别是位温、平均风速和虚位温在垂直方向上的差值。 $ Ri_{\rm{b}} $ ＞0时的 $ \Delta z $ 即为模式边界层高度（H_PBL）（张碧辉等，2012）。因此，可以计算出新的位温垂直变化率（ ${\theta }_{\rm{lapes\_new}}$ ）

${\theta }_{\rm{lapes\_new}} =\frac {{\theta _{{\rm{PBLH}}}} - {\theta _{10{\text{%}} {\rm{PBLH}}}}}{H_{\rm{PBL}} - 10{\text{%}} H_{\rm{PBL}}}$

(13)

式中， $ {\theta _{{\rm{PBLH}}}}$ 和 ${\theta _{10{\text{%}} {\rm{PBLH}}}}$ 分别为位于边界层高度和近地层高度处的模式位温。通过式（6）、（7）、（8）和（13）可推知新的位温订正值（ $ {\theta }_{\rm{sfc\_new}} $ ），并可将其换算为新的温度订正值

$\left\{ \begin{array}{l} {\theta }_{\rm{sfc\_new}}={\theta }_{\rm{ms}}+{K{\text •}\Delta \theta }_{{\rm{o}}}\\ K=\dfrac{H_{\rm{PBL}}-{h}_{{\rm{m}}}}{H_{\rm{PBL}}-{h}_{{\rm{o}}}} \end{array}\right.$

(14)

3.3 不同方案的订正效果分析

根据原方案计算出的风速订正系数基本上均为1。但是从式（11）可知，采用更新方案后，近地面风速在稳定条件下，随着高度升高而增大，但是地形高度差的大小对订正系数的影响较小，订正系数值在1.2左右；而在不稳定条件下，风速的订正值则随着地形高度差的增大而呈指数增长。更新方案采用随时间变化的边界层高度值来计算位温垂直变化率，稳定条件下，采用更新方案计算的各站点位温的垂直变化率小于1.0 K/（100 m），而不稳定条件下的位温垂直变化率较大，最大可接近2.0 K/（100 m），而原方案计算的位温垂直变化率则小于0.75 K/（100 m）。

更新方案对风速及温度的订正是否更加接近实际情况，需要通过与实际观测的对比才能进一步验证。

利用中国所有探空观测站每天00和12时，以及个别站点如54511站在汛期每天06时增加的秒级探空观测数据，采用不同订正方案将与模式最低层高度同高位置的探空观测数据作为真值，对比地面观测的订正值与真值的偏离程度。2018年6—8月00和12时共获取7701个有效的风速真值数据，如图3a₁、b₁、c₁所示，使用原方案，00和12时风速订正值的绝对偏差分别为1.0和1.5 m/s，采用更新方案后，相应的绝对偏差值降为0.9和1.3 m/s；图3a₂、b₂、c₂是温度的秒探空观测值与订正值的对比。因为温度的秒探空观测大量缺测，00和12时共获取5702个有效数据。采取更新方案后温度订正值的绝对偏差分别由原来的1.1和1.4℃下降为1.0和1.2℃。而2018年6—8月每天06时均有探空观测的站点较少，图3中06时对比为54511站的结果，分别获得65个有效风速真值数据和62个有效温度真值数据。使用更新方案后风速订正值的绝对偏差平均值由0.88 m/s降低到0.76 m/s，温度则由0.95℃降低到0.86℃。

因为受实际观测资料的局限，所以新方案在风速订正时通过模式背景场来判定大气层结状态，采用背景场与采用实际观测推算出的大气层结状态有差异。以06时的54511站为例，2018年6—8月有探空观测的65 d中，模式和观测推算大气层结状态一致的为22 d，不一致的为43 d。表1是模式背景场与观测在大气稳定度一致（情况A）和不一致时（情况B）不同方案对于温度和风速的订正效果。可见两种情况下，采用更新方案后统计温度订正值与同高度的观测值的偏差（Mean Error，ME）及均方根误差（RMSE）均降低；更新方案的风速订正值的偏差和均方根误差也明显优于原方案。由于实际探空观测的站点数少且观测频率低，因此同化地面风速观测时使用模式背景场来判定大气稳定度情况，且不论背景场与观测得到的稳定度状况是否一致，新方案的订正效果均有明显优势。

以上结果表明，更新方案对近地面风速的描述更准确，订正后的风速值更接近实际观测值，且夜间（12时）大气稳定状态下更新方案的风速订正效果更具优势；夜间边界层高度相对较低，原方案利用较高层计算出的位温垂直变化率和边界层低层内位温垂直变化率的差别相对较大，因此使用更新方案后夜间的温度订正值更准确，且无论白天还是夜间更新方案对于地面温度的订正误差更小。

4 试验及结果分析

设置不同的试验旨在分析更新方案对数值预报结果产生的影响。单点试验用于验证地面资料同化更新方案的同化可靠性，并通过个例分析和批量试验结果来验证更新方案对于模式地面要素预报性能的改进。

4.1 单点试验

单点试验即同化一个观测站点的单要素观测数据，旨在分析各要素分析增量的空间分布特征。试验选取不同的地面资料同化方案（表2），单次分别同化地面自动站53391（41.9°N，114.0°E）的温度、风和相对湿度观测数据，站点53391的模式和观测的地形高度差为30 m，同化模块使用控制变量7，得到各要素的水平和垂直响应如图4和5所示。

两组试验各要素的分析增量均满足各向同性的高斯分布特征，符合三维变分同化特征。可以看到，无论是中心点的强度还是影响范围，使用更新方案后的温度场的响应明显增强，说明采用不同的温度订正方案对同化效果的影响非常大且更新方案对背景场温度的调整力度更大；同样使用更新方案后风场的中心点强度和影响范围明显更大，这与订正方案中地面粗糙度的改变以及稳定度项的加入有关；更新方案中对湿度观测的同化方式未做任何调整，所以两组试验湿度的分析增量变化一致（图略）。

4.2 同化试验结果及分析

选取2017年7月5日的个例，采用不同的地面资料同化方案进行同化试验，对比分析两组试验的分析场和预报场结果，旨在清楚了解更新方案对预报性能的影响。同化试验的设计方案如表3所示，00时的模式采取冷启动方式。

4.2.1 分析增量结果

00时冷启动的模式分析场如图5所示，原方案试验的温度场分析增量在陆面上以负增量为主，在山东北部出现增量极大值中心。在此增量变化基础上，更新方案试验的温度增量进一步呈现模拟区域负的变化，极大值分别出现在山东和山西境内。不同试验温度场的分析增量的分布特征有显著差异，而两个试验风场的分析增量差异分布没有区域性规律，可能是00时大气层结相对稳定，边界层高度较低，不同方案推导出的温度垂直变化特征的差异明显，而不同方案计算得出风的订正系数本身差异小且地面风速观测值较低，多数站点订正后的风速值无明显差别，因此风场的分析增量差异并不如温度的增量差异明显。

4.2.2 地面要素预报结果

图6a 为原方案试验的6 h温度预报偏差（图6中仅给出偏差绝对值不小于0.1℃的结果，后文不再赘述），可见明显偏高。使用更新方案后同时刻参与对比分析的750个观测站中62%的站点温度预报改善明显（图6b），其中京津冀地区的降温效果尤为显著；23%的站点温度预报无明显差异即两组试验的温度预报差值在±0.1℃；15%的站点预报偏差幅度增大。对模拟区域平均而言，使用更新方案后温度的预报偏差由原来1.27℃降低为1.02℃，降幅达20%，均方根误差也由2.62℃降低为2.5℃。

同样，6 h风速预报的性能也有所提高，模拟区域平均的风速预报偏差从0.5 m/s降低到0.4 m/s。使用更新方案后，风速预报性能提高的站点（图6d）占比为36%，不同试验风场的分析增量差别较小，因此58%的站点风速预报值无明显差异，即预报差值小于0.1 m/s。

由模拟区域24 h预报结果的客观统计得知，两组试验分析时刻和预报时刻的温度偏差均为正。原方案试验对于白天温度预报的效果较差，相比之下更新方案试验中白天的温度预报效果有明显提升。由图7可知更新方案试验分析场的温度与实况更接近，即更新方案可以提高地面温度资料同化效果并改善模式的分析场，进而提高温度预报的准确度；通过两组试验的风速预报偏差和均方根误差也可看出，使用更新方案后前15 h风速预报正偏差有所降低。虽然更新方案中风速的订正效果会受到模式地面粗糙度精度和推算的模式大气稳定度的影响，但是更新方案对于近地面风速的描述更准确，因此使用后可以改善风速的预报效果。

采取相同参数配置进行循环同化试验，进一步分析更新方案对采用逐3 h循环同化方式运转的数值预报系统产生的影响。选取00时启动的模式3 h预报结果作为03时启动的模式同化所需的背景场进行循环同化试验，预报时效为24 h。03时启动的更新方案和原方案试验不仅同化方案有区别，且背景场也有差异。结果（图8）可见，经过循环同化后更新方案试验的温度和风速的预报偏差全面降低，其均方根误差也更小。综上所述，更新方案中温度和风的同化方法更合理，使用后能提高预报的准确性，而循环同化时，前一个启动时刻精准的预报场作为模式背景场输入到同化系统中再叠加更新方案的效果，可以显著提升循环同化后模式的预报性能。

4.3 批量试验结果客观分析

为避免个例分析的偶然性，选取2017年7月1—31日00时起预报的模式结果进行批量试验，两组批量试验的唯一区别为使用的地面资料同化方案不同，试验的客观统计结果如图9所示。W₁₀、Q₂和T₂分别表示10 m风速、2 m比湿和2 m气温。可见，所有预报时次使用更新方案后温度的均方根误差和偏差均有改善，前9 h的温度预报性能有明显提升；大部分时次的风速预报效果也有所提高，但是风速预报的效果没有温度好；不同方案中均没有针对湿度观测的订正，因此湿度的预报效果基本上无变化。批量试验得出的结果与个例分析一致，同样说明更新方案可以提高温度和风速预报的准确度。

5 结论和讨论

基于WRF-DA中包含地形差造成影响的地面资料同化方案1进行优化，主要改进了原同化方案中地面风和温度的订正方案。风速同化方案基于相似理论设计而成，将风速观测值由观测高度订正到模式最低层高度。更新方案在原风速订正方案中增加大气层结稳定度特征参数和模式地表粗糙度参数，将仅适用于中性层结的原风速订正方案改进为适用于任意大气层结状态下的近地面风速订正方案；温度同化方案基于位温垂直变化率特征，将站点的温度观测值由观测高度订正到模式高度，因此位温垂直变化率的准确与否直接影响订正效果，更新方案中最主要的更改内容是将原方案中采用边界层高度外的两个高层计算的位温垂直变化率更正为适用于边界层内的位温垂直变化率，进而对地面观测温度进行高度订正。通过个例分析和批量试验结果评估得出以下结论：

（1）使用更新方案将地面温度观测订正到模式高度后的订正值与实况更接近。同样采用更新方案后地面风速观测的订正值与实况也更为一致。说明使用更新方案对于地面要素的订正效果更好。

（2）2017年7月5日个例分析以及2017年7月整月的批量试验结果评估表明，冷启动时采用更新方案对模式分析场的温度和风速有明显调整，24 h内模式温度和风速的预报准确度提高，尤其是前9 h的模式预报效果有明显改善。

（3）选取不同方案进行循环同化试验发现，上一个启动时刻更新方案对于模式预报的改善作用可以通过背景场传递到本次同化时刻，再叠加更新方案对于模式预报的正面效果，因此选取更新方案进行循环同化后，模式温度和风速的预报性能显著提升。

观测资料在测量、传输和解算等过程中会产生随机误差或系统性误差，而三维变分同化系统中假设观测误差服从高斯概率分布来确定大气状态的最优线性无偏估计。由于无法获得真值，不能直接验证观测误差是否满足高斯分布特征，因此常用的验证方式为诊断观测与相应背景场的差值（观测背景差）是否满足高斯分布特征。文中使用的地面自动站观测数据为经过严格质量控制后的国家级自动站观测，在下一步工作中将开展中国加密自动站数据同化应用研究，因此需要通过面向同化应用的地面观测资料质量控制，使得观测背景差所对应的观测资料满足同化算法关于观测性质的假设，有助于地面加密观测在模式系统中的同化应用。

除此以外，湿度在垂直方向上的变化较为复杂，更新方案和原同化方案中均没有地面湿度观测的订正方案，即假设在一定地形差范围内的湿度不随垂直高度的改变而变化，因此地面自动站中的湿度信息同化时存在一定误差，在下一步的研究工作中拟通过大量的湿度探空观测值来拟合出具有普适意义的湿度垂直廓线，并将其应用到同化系统中对湿度观测进行垂直订正，改善模式对地面水汽的预报效果。