概要设计¶
OSE¶
观测系统影响试验（OSEs）是直接采用同化的方法检验观测影响的一种方法。该方法有别于仿真观测评估（OSSE）、以及观测系统贡献评价（FSO）等在模拟或者理想状态下对观测资料的贡献评价，OSE针对的是观测资料在实际业务同化预报系统中的实际影响，考察的是被考察资料在所有业务同化观测资料中的真实贡献。
OSE试验流程概念图
在观测站网设计与评估的其它工具、手段中，OSSE采用的是模拟的观测资料，资料来源于数值模式试验中假定的“真值”；TO方法考虑的是资料的布局（空间分布、水平分布）以及天气或气候系统对关键区域的敏感性；FSO在综合评价各种观测资料的作用时，采用的是特定的资料同化方法，该方法可能与实际业务中的数值模式、资料同化方法并不一致。OSE的方法针对的就是实际业务模式、实际业务同化资料系统以及实际同化的观测资料，考察的是在真实的业务资料同化系统中，当去除或添加某一类或某一部分资料时，资料同化对数值预报各天气、气候系统的预报效果的影响。
在研究和业务工作中，OSE一直被广泛应用于对观测资料的同化效果的评估。当需要评估特定的观测资料的同化效果时，一般设计如下两个对比试验来进行考核：
- 设置一“控制试验”，该试验可以同化所有设定的观测资料，或者不包含需要考察的观测资料的其它资料；
- 设置“敏感性试验”，该试验在控制试验的基础上去除或者添加需要被考察的观测资料，该被考察的资料可以是某一类资料、某一时段资料、某一资料的时、空分辨率的改变等。
OSSE¶
OSSE是用一个性能良好的模式，模式自由运行一段时间模拟大气作为大气的真值，称为Nature run, 并且用此真值加上随机误差构造观测，然后将该构造的资料应用在另一个数值模式系统，做同化预报，结果和Nature run 比较，以此来评估新增观测资料的影响。其流程结构如下图所示：
OSSE试验流程概念图
OSSE试验一般由以下功能模块构成：
- 大气三维综合仿真观测构建模块：建立Nature Run参考大气，构建大气三维综合仿真观测系统，模拟实现天地空各类观测。
- 同化预报模块：同化模拟观测资料并基于同化结果进行数值天气预报。
- 试验结果评估模块：建观测系统配置组合影响评估模块，利用GRAPES模式后处理模块以及fortran、matlab、grads、ncl等编程画图软件计算预报结果的预报偏差、降水TS评分、预报与实况对比、地面要素检验等，再进行统计对比分析，对模拟观测的贡献及影响进行评估。
基于集合预报的FSO方法¶
FSO是评估预报误差对观测敏感性的一种分析方法。它不但可以评估观测系统中具体组成部分（如探空观测系统，地面观测系统）对改进预报技巧的有效性，而且能够定量给出各观测变量（如风场，温度）在当前预报中的相对影响，同时给出这种影响的空间分布。FSO应用的一个重要假设位观测作用可以线性分解与叠加的，即预报误差中的主要部分是可以用非线性模式的切线性模式或是一组集合预报非线性扰动的线性组合表示。
目前，FSO方法的实现上主要有两种方式：1）基于伴随模式的敏感性分析方案；2）基于集合预报的敏感性分析方案。两种方案的数学原理是一致的，主要不同之处实践过程中如何描述扰动的线性演变与同化分析过程中增益矩阵的处理。
FSO的数学原理简要描述如下（Baker，2000）。在一个循环资料同化预报系统中，当前分析时刻的分析场一般由前一个时刻分析场产生的短时预报场经过当前观测修正而得到。令 \(x_a\) 表示分析时刻 \(t=0\) 的分析场，而 \(x_b\) 表示前一个同化预报循环产生的短期预报背景场。则观测对分析的影响可以表示为:
(1)¶\[x_a = x_b + \mathrm{K} \delta y \; ;\; \delta y = y - H(x_b)\]
式 (1) 中 \(x_a\) 为分析场，\(x_b\) 为背景场，\(y\) 为观测值，\(\mathrm{K}\) 为 Kalman 增益矩阵，\(H\) 为观测算子，\(\delta y\) 为更新变量。从 (1) 可以推导出分析场对观测和背景场的敏感梯度可写为：
(2)¶\[\frac{\partial x_a}{\partial y} = \mathrm{K}^\intercal \; ; \; \frac{\partial x_a}{\partial x_b} = \mathbf{I} - \mathrm{H}^\intercal \mathrm{K}^\intercal \; ; \; \mathrm{H} = \frac{\partial x_a}{\partial x_b}\]
式 (2) 中 \(\mathrm{H}\) 为切线性观测算子,上标 \(\intercal\) 代表矩阵转置，\(\mathrm{K}^\intercal\) 代表了同化系统的伴随系统。
考虑一个预报，其预报方程可以写为:
(3)¶\[x^f = M(x^0)\]
式 (3) 中 \(x^f\) 为模式预报场，\(x^0\) 为模式初始场，\(M\) 为非线性数值预报模式。
定义在某时刻的一个模式预报误差函数:
(4)¶\[J = \frac{1}{2} ( x^f - x^t )^\intercal \mathrm{C} ( x^f - x^t )\]
式 (4) 中，\(\mathrm{C}\) 为权重矩阵，一般为对角矩阵。从 (4) 中可以推导出预报误差函数 \(J\) 对初始场（如果经过了同化分析，初始场就是同化系统中的分析场，\(x^0\) 即为 \(x_a\) ）的敏感性梯度：
(5)¶\[\frac{\partial J}{\partial x_a} = \frac{\partial x^f}{\partial x_a} \mathrm{C} ( x^f - x^t ) = \frac{\partial M(x_a)}{\partial x_a} \mathrm{C} ( x^f - x^t ) = \mathrm{M}^\intercal \mathrm{C} ( x^f - x^t )\]
式 (5) 中 \(x^t\) 为真值，\(\mathrm{M}\) 为非线性模式的切线性模式，\(\mathrm{M}^\intercal\) 为 \(\mathrm{M}\) 模式的伴随模式。联合式 (1)，(4), (5) 可推导出模式预报误差函数 \(J\) 对观测的敏感梯度：
(6)¶\[\frac{\partial J}{\partial y} = \frac{\partial x_a}{\partial y} \frac{\partial J}{\partial x_a} = \mathrm{K}^\intercal \frac{\partial J}{\partial x_a} = \mathrm{K}^\intercal \mathrm{M}^\intercal \mathrm{C} ( x^f - x^t )\]
因此，从式中可以看到，预报误差对观测的敏感性梯度需要计算模式的伴随和同化系统的伴随。在同化系统中，观测资料对预报误差的贡献（Observation impacts）一级近似为：
(7)¶\[\delta J = \langle \frac{\partial J}{\partial y}\, , \delta y \rangle = \langle \mathrm{K}^\intercal \frac{\partial J}{\partial x_a }\, , y - H(x_b) \rangle= \langle \frac{\partial J}{\partial x_a} \, , \mathrm{K} ( y - H(x_b)) \rangle = \langle \frac{\partial J}{\partial x_a} \, , \delta x_a \rangle\]
式 (7) 中，\(\langle \, , \rangle\) 为内积算符。
由式 (7) 可知，观测对于预报的影响可以解释为：\(\delta y = y- H(x_b)\) 在预报误差对于观测敏感性梯度上的投影；或是，分析增量在预报对于初值的敏感性的投影。从式 (7) 可以看出，当误差改变为负值时，预报误差减少，同化观测提高预报能力；当 \(\delta e\) 为正值时，预报误差增大，同化观测降低预报能力。具体实现上，在四维变分同化系统中，可以利用既有的伴随模式及隐式求解增益矩阵 来评估分析增量及每种观测的对于预报误差的贡献。
在实际的同化系统中，同化资料后主要考察在背景场的基础上同化资料对预报误差的贡献。由背景场作为初始场的预报误差和分析场作为初始场的预报误差分别为:
定义 \(e_a\) 和 \(e_b\) 的差为：
(10)¶\[\delta e = e_a - e_b\]
所以
(13)¶\[\delta e = \langle \frac{\partial e_a}{\partial x_a^f} , (x_a^f - x^t) \rangle - \langle \frac{\partial e_b}{\partial x_b^f} , (x_b^f - x^t) \rangle\]
(14)¶\[\delta e = \langle (x_a^f - x_b^f) , \frac{\partial e_a}{\partial x_a^f} + \frac{\partial e_b}{\partial x_b^f} \rangle = \langle (x_a - x_b) , \frac{\partial e_a}{\partial x_a} + \frac{\partial e_b}{\partial x_b} \rangle\]
式 (14) 利用式 (5) \(x_a^f-x_b^f\) 可以近似展开 (15) 可以推导出右边的表达式
(15)¶\[x_a^f - x_b^f = \mathrm{M}_b^\intercal (x_a -x_b) = \mathrm{M}_a^\intercal (x_a - x_b)\]
因此采用不同的近似，(14) 有不同的表达形式。Gelaro等（2007）给出了5种不同的表达形式：
(16)¶\[\delta e_1 = 2(x_a - x_b^\intercal \mathrm{M}_b^\intercal \mathrm{C} (x_a^f - x^t)\](17)¶\[\delta e_2 = (x_a - x_b)^\intercal [\mathrm{M}_b^\intercal \mathrm{C} (x_a^f - x^t) + \mathrm{M}_a^\intercal \mathrm{C} (x_b^f - x^t)]\](18)¶\[\delta e_3 = (x_a - x_b)^\intercal [\mathrm{M}_b^\intercal \mathrm{C} (x_b^f - x^t) + \mathrm{M}_a^\intercal \mathrm{C} (x_a^f - x^t)]\](19)¶\[\delta e_4 = (x_a - x_b)^\intercal [\mathrm{M}_a^\intercal \mathrm{C} (x_a^f - x^t) + \mathrm{M}_a^\intercal \mathrm{C} (x_b^f - x^t)]\](20)¶\[\delta e_5 = (x_a - x_b)^\intercal [\mathrm{M}_b^\intercal \mathrm{C} (x_b^f - x^t) + \mathrm{M}_b^\intercal \mathrm{C} (x_a^f - x^t)]\]
(21)¶\[\delta e_1 = 2 \delta y^\intercal \mathrm{K}^\intercal \mathrm{M}_b^\intercal \mathrm{C} (x_a^f - x^t)\]
基于集合的FSO，关键在于利用一组集合扰动的线性组合来近似表示扰动的时间演变，进而省略非线性模式的切线性与伴随模式。在集合预报系统中，已知分析扰动与预报扰动的前提下，式 (6) - (7) 中，增益矩阵及线性模式与增益矩阵的组合可以表示为 (Kalnay 等，2012):
(22)¶\[\mathrm{K} = \frac{1}{K-1} \mathrm{X}_0^a {\mathrm{X}_0^a}^\intercal \mathrm{H}^\intercal \mathrm{R}^{-1}\](23)¶\[\mathrm{M} \mathrm{K} = \frac{1}{K-1} \mathrm{M} \mathrm{X}_0^a {\mathrm{X}_0^a}^\intercal \mathrm{H}^\intercal \mathrm{R}^{-1} \approx \frac{1}{K-1} \mathrm{X}_{t|0}^f {\mathrm{Y}_0^a}^\intercal \mathrm{R}^{-1}\](24)¶\[\bigtriangleup e^2 \approx \frac{1}{K-1} \delta y_0^\intercal \mathrm{R}^{-1} \mathrm{Y}_0^a \mathrm{X}_{t|0}^{f^\intercal} \mathrm{C} (e_{t|0} + e_{t|-6})\]
其中：\(\mathrm{Y}_0^a = \mathrm{H} \mathrm{X}_0^a\) 为分析扰动在观测空间的表示。
本项目建立基于集合的FSO（Kalnay等，2012）。具体实现表示如下：我们将应用一个时间序列的集合预报的扰动，来近似表示非线性模式的线性模式进而获得预报误差对初值的敏感梯度；以及观测集合的扰动观测来简化线性模式 \(\mathrm{M}\) 与增益矩阵 \(\mathrm{K}\) 的求解，进而实现观测对预报贡献的定量评估。