CN104834823A - 基于卫星-地基co2数据联合同化的碳源汇估测方法 - Google Patents

Abstract

本发明公开了基于卫星-地基CO₂数据联合同化的碳源汇估测方法，所述方法结合CO₂站点浓度、CO₂卫星柱浓度，采用集合卡尔曼平滑方法求解出最优化的CO₂浓度和通量时空分布场。本发明引进柱浓度同化方案，构建了基于卫星柱浓度及地基站点观测数据为基础的卫星-地基CO₂联合同化方法，改进了大气反演模型只能同化站点观测CO₂浓度的缺陷，实现了多种来源CO₂数据的联合同化，为精确的区域碳源/汇估算提供有效手段。

Description

基于卫星-地基CO2数据联合同化的碳源汇估测方法

技术领域

本发明涉及遥感技术、大气科学、生态学等学科领域，特别是涉及基于卫星-地基CO₂数据联合同化的碳源汇估测方法。

背景技术

人类活动引起的大气CO₂浓度升高导致全球气候变化，给人类社会的生存和发展带来了巨大威胁。陆地生态系统碳汇与大气CO₂间存在着紧密关系，它不仅被动地受大气CO₂的影响和控制，也对大气CO₂存在反馈作用，在气候系统中具有举足轻重的地位。因此，在气候与环境问题日益突出的今天，加深对陆地生态系统碳源汇的动态变化机制的理解，加强区域陆地碳源汇分布格局的认识，具有重要的科学价值。

目前，国际上已发展了一系列直接(如：箱式通量观测法及通量塔观测法)或间接(如：“自下而上”的陆面生态过程模型模拟、土地利用和土地覆被变化监测法，“自上而下”的大气CO₂反演法)的碳源汇观测与估算方法。其中，“自上而下”的大气反演法是碳源汇计算的重要方法，在估算区域和全球生态系统碳收支和预测未来气候变化和人类活动改变条件下生态系统碳循环过程提供了有效的手段。然而，传统的大气CO₂反演法是利用分散在世界各地的近地面CO₂浓度观测站(仅约200个点)数据作为观测值，结合数据同化方法进行的碳源汇推算。由于观测站点的数量不足、分布不匀及观测指标的不统一，使全球二氧化碳浓度/源汇估算结果具有相当大的不确定性。因此，地基碳观测数据有限的现状，对准确估算全球碳源汇分布以及大气CO₂浓度总体状况提出了严重的挑战。

大气CO₂反演法的模拟精度对观测数据的依赖性强，其观测数据的数量和质量决定了碳源汇模拟结果的时空分辨率和精度。由于目前可用的地基CO₂观测站点数量不够、数据分布不均及观测指标不统一等限制，大气反演法往往只能用有限的站点去估测几年内(10-20年)全球几十个(22-100)区的平均碳通量的分布情况，其模拟结果的精度、时间和空间分辨率都比较低。

鉴于以上问题，开发一种能够提高时间和空间分辨率以及估算精度的碳源汇估测方法具有重要意义。

发明内容

本发明的目的是提供一种基于卫星-地基CO₂数据联合同化的碳源汇估测方法，使其对碳源汇的估测结果具有较高的精度和时空分辨率，从而克服现有方法对碳源汇估测结果精度低、时空分辨率低的不足。

为实现上述目的，本发明采用如下技术方案：

基于卫星-地基CO₂数据联合同化的碳源汇估测方法，包括以下步骤：根据先验CO₂通量数据驱动大气传输模型，获得CO₂模拟浓度时空数据集；采用集合卡尔曼平滑方法同化CO₂观测数据，对所述CO₂观测数据与CO₂模拟浓度之差最小化求解，获得优化的CO₂浓度及CO₂通量的时空分布场；其中，所述CO₂观测数据包括CO₂站点浓度数据及CO₂卫星柱浓度数据：当所述CO₂观测数据为CO₂站点浓度时，从所述CO₂模拟浓度时空数据集中提取时空信息匹配的CO₂模拟浓度，对该CO₂模拟浓度与所述CO₂站点浓度之差最小化求解；当所述CO₂观测数据为CO₂卫星柱浓度时，从所述CO₂模拟浓度时空数据集中提取时空信息匹配的CO₂模拟浓度，先将该CO₂模拟浓度转化成CO₂模拟柱浓度，再对CO₂卫星柱浓度与转化后的CO₂模拟柱浓度之差最小化求解。

进一步地，所述CO₂站点浓度数据来自美国海洋与大气总署的地球系统研究实验室及澳大利亚联邦科学与工业研究组织的瓶装采样的flask数据；所述CO₂卫星柱浓度数据来自GOSAT-ACOS3.3中的L2数据产品。

进一步地，所述进行同化过程的CO₂卫星柱浓度数据选取GOSAT-ACOS V3.3中标记为Good，且满足以下条件的数据：卫星柱浓度数据与其对应的转化后的CO₂模拟柱浓度相差小于3ppm。

进一步地，所述同化过程是通过大气反演模型进行的。

进一步地，具体包括以下步骤：S1：根据碳通量调整系数λ_t的背景场利用方程1-1计算先验CO₂通量数据，即CO₂通量预报场；

x_(t)＝f(λ_t)＝[(F_bio,t+F_ocn,t]×λ_t+F_ff,t+F_fire,t         (1-1)

式中，x_(t)代表CO₂通量，F_bio,t、F_ocn,t、F_ff,t、F_fire,t分别代表陆表碳通量、海洋碳通量、人为碳排放以及野火燃烧碳排放，λ_t为碳通量调整系数，下标t代表时间，当t＝0时，λ_t为初始值1.0；S2：根据所述先验CO₂通量数据，驱动大气传输模式运行，生成CO₂模拟浓度时空数据集；S3：采用集合卡尔曼平滑方法对时空信息匹配的CO₂观测浓度与CO₂模拟浓度之差最小化求解，获得优化的碳通量调整系数即λ_t的分析场，再根据利用方程1-1计算优化的CO₂通量，即CO₂通量分析场；S4：根据优化的CO₂通量再次驱动大气传输模式，获得优化后的CO₂模拟浓度时空数据集；S5：将优化的碳通量调整系数作为下一时刻的碳通量调整系数的背景场，返回步骤S1。

进一步地，所述陆表碳通量F_bio,t是通过陆地生态系统模型获得的。

进一步地，所述将CO₂模拟浓度转化成CO₂模拟柱浓度是通过公式1-2进行的：

$X_{\mathrm{CO}2}^m=h^T{X}_a+h^{T}A(X_h-X_a)---(1-2)$

式中，是转化后的CO₂模拟柱浓度；h^T是卫星提供的气压加权函数；X_a是卫星提供的先验廓线；A是卫星提供的平均核函数；X_h是与卫星提供的先验浓度廓线X_a时空尺度相匹配的、大气传输模式模拟的模拟浓度廓线。

进一步地，所述模拟浓度廓线X_h通过以下方式实现与所述先验浓度廓线X_a时空尺度相匹配：大气传输模式首先根据卫星柱浓度数据的时间、空间信息对CO₂模拟浓度进行采样，获得时空信息相一致的CO₂模拟浓度，再根据先验浓度廓线X_a的分层高度对CO₂模拟浓度进行插值，以生成与先验浓度廓线X_a时空尺度相匹配的模拟浓度廓线X_h。

进一步地，所述同化过程通过最小化代价函数J获取要求解的CO₂通量x及其分析误差P^a的最优化估测：

x＝x₀+k(y⁰-H(x₀))             (1-6)

P^a＝(I-kH)P                   (1-8)

其中，k是一个增益矩阵，表示如下方程：

k＝(PH^T)/(HPH^T+R)             (1-7)

所述代价函数J：

$J=\frac{1}{2}{(y^0-H\left(x\right))}^T{R}^{-1}(y^0-H\left(x\right)))+\frac{1}{2}{(x-x_0)}^T{P}^{-1}(x-x_0)---(1-5)$

式中，y⁰是CO₂观测浓度，R＝var(ζ)是观测浓度误差协方差矩阵；x₀是先验CO₂通量，P是先验CO₂通量的误差协方差矩阵；上标T代表矩阵转置；

H为如下所示的观测算子矩阵：

其中H_ij代表像元点i、地基站点浓度j的地基站点观测算子；W_ij代表的是像元点i、卫星柱浓度j的卫星浓度观测算子。

由于采用了以上技术方案，本发明至少具有以下优点：

本发明建立了基于卫星柱浓度数据及地基站点观测数据为基础的卫星-地基联合同化方法，这种方法改进现有的大气反演模型只能同化站点观测数据的缺陷，在原来的地基站点数据同化方法上引进柱浓度数据同化方案，实现地基站点观测数据和卫星柱浓度数据的联合同化，为精确的估算区域碳源/汇提供有效手段。

本发明是大气CO₂反演法的一个崭新的尝试和发展，改变大气反演法只能用站点观测数据去估测几年内(10-20年)全球几十个(22-100)区的平均碳通量的现状，从而大大提高模拟结果的精度、时间和空间分辨率。

附图说明

上述仅是本发明技术方案的概述，为了能够更清楚了解本发明的技术手段，以下结合附图与具体实施方式对本发明作进一步的详细说明。

图1是本发明的总体技术路线图。

图2是卫星-地基CO₂数据联合同化技术路线图。

图3是CO₂浓度观测与模拟值间的误差分布图。

图4是CO₂浓度观测与模拟值间的相关性分析图。

图5A、5B分别是2010年间在站点WSA(united States；49.93°N，60.02°W，30m)的CO₂浓度观测与模拟值间的对比图、分析图。

图6是2010年间CONTRAIL的CO₂浓度观测与模拟值间在三个高度：(a)485-525hPa；(b)375-425hPa；(c)225-275hPa的对比分析图。

图7是2010年全球CO₂碳通量的组成图。

图8是2010年中国陆地CO₂通量的分布图。

图9是2010年“卫星-地基”联合同化实验、“地基”同化实验的中国陆地碳通量标准差对比图。

图10是2010年“卫星-地基”联合同化实验、“地基”同化实验在站点Minamitorishima(MNM，Japan；24.29°N，153.98°E，8m)的浓度对比分析图。

图11是“卫星-地基”联合同化实验、“地基”同化实验在中国陆地生态系统碳汇分布差异图。

具体实施方式

本发明提供了基于卫星-地基CO₂数据联合同化的碳源汇估测方法，该方法通过引入卫星柱浓度数据，构建了可对卫星柱浓度与站点观测数据进行同时同化的全球CO₂数据同化系统，改进了现有的大气反演模型，大大提高碳源汇估测结果精度以及时间和空间分辨率。

请参阅图1所示，本发明所述的碳源汇估测方法的总的技术路线如下：首先，以陆地生态系统模型(Remote sensing-based ecosystem-atmospheresimulation scheme,EASS)(Chen B.et al.,2007)的陆表碳通量模拟结果，结合海洋碳通量、人为碳排放以及野火燃烧碳排放等CO₂通量数据为初始化数据，根据方程(1-1)获得先验CO₂通量：

x_(t)＝f(λ_t)＝[(F_bio,t+F_ocn,t]×λ_t+F_ff,t+F_fire,t          (1-1)

式中，x_(t)代表CO₂通量，F_bio,t、F_ocn,t、F_ff,t、F_fire,t分别代表陆表碳通量、海洋碳通量、人为碳排放以及野火燃烧碳排放，λ_t为碳通量调整系数，下标t代表时间。

利用先验CO₂通量，结合气象数据等驱动大气传输模型，模拟出CO₂浓度的时空分布图。在此基础上，采用集合卡尔曼平滑(Ensemble KalmanSmooth,EnKS)方法，对比、分析CO₂观测浓度(包括站点浓度和卫星柱浓度)和CO₂模拟浓度，基于贝叶斯原理求解出最优化的CO₂浓度和通量分布场，实现高精度的CO₂估算。其中，所涉及的核心方法及关键技术主要包括：(1)卫星-地基CO₂数据联合同化技术和(2)集合卡尔曼同化技术，具体阐述如下：

1 卫星-地基CO₂数据联合同化技术

卫星监测数据具有覆盖面广、时空分辨率高的特点。因此，引入卫星监测的CO₂柱浓度数据可改善大气CO₂反演模型观测数据不足的现状。卫星-地基CO₂数据联合同化方法的关键技术是如何把大气传输模型模拟的垂直分层的CO₂浓度转化为与卫星柱浓度相匹配的时空尺度。我们利用Rodgers等(Rodgers C.D.and Connor B.J.2003.Intercomparison of remotesounding instruments.Journal of Geophysical Research,108,4116)提出的方法将大气传输模式模拟的分层CO₂浓度转化为柱平均浓度，公式如下：

$X_{\mathrm{CO}2}^m=h^T{X}_a+h^{T}A(X_h-X_a)---(1-2)$

其中，是转化后的模拟平均柱浓度；h^T是卫星提供的气压加权函数；X_a是卫星提供的先验浓度廓线；A是卫星提供的平均核函数；X_h是与卫星提供的先验浓度廓线X_a时空尺度相匹配的、大气传输模式模拟的浓度廓线。

将分层模拟浓度转化为柱浓度后，利用集合卡尔曼平滑方法(EnKS)对观测浓度(包含站点浓度、卫星柱浓度)进行统一同化。

请配合参阅图2，整个技术路线分2步走：

(1)根据观测点的经纬坐标(x，y)、时间(t)、高度(z)信息提取相应的模拟浓度，首先判断当前的观测数据类型：如果是站点数据，系统直接对模拟浓度进行采样提取相应浓度数据；如果不是站点数据，则要根据公式1-2对大气传输模式模拟的分层浓度进行转化，计算出与观测数据相对应的柱平均浓度。

(2)再通过最小化观测浓度与模拟浓度间的差值来同化CO₂的浓度与通量。

2集合卡尔曼同化方法

利用集合卡尔曼平滑方法对卫星和站点观测数据进行同化，实现优化区域CO₂浓度和通量时空数据的目标。

其主要步骤是：(1)将陆表碳通量、海洋碳通量、人为碳排放以及野火燃烧碳排放四个CO₂输入数据结合碳通量调整系数，计算出先验CO₂通量数据；(2)用大气传输模式(TM5)将上述先验CO₂通量数据扩散到大气中的不同位置和高度，获得大气CO₂模拟浓度时空数据集(x,y,z,t)；(3)根据观测数据时空信息对模拟浓度时空数据集进行采样和尺度转化，获取与CO₂观测浓度时空及尺度信息相匹配的点状模拟浓度；(4)计算出最小模拟和观测的大气CO₂浓度之差，得到优化的碳通量调整线性系数，据此调整先验碳通量场，获得优化的后验碳通量场，并进一步利用集合卡尔曼平滑EnKS方法确定下一时间步长的碳通量调整线性系数；(5)用后验碳通量场重新计算大气CO₂浓度分布。

在列出集合卡尔曼平滑方法的具体步骤之前，首先，需要列出所述CO₂浓度与CO₂通量间的线性关系：

y⁰＝G(f(λ))+A(y)+ξ             (1-3)

式中，y⁰是CO₂观测浓度；

ξ是CO₂观测浓度误差，代表的是一个平均值为0、协方差矩阵为R＝var(ξ)的向量；

G代表的是观测算子矩阵，其作用是通过CO₂通量值计算出浓度的变化量；

f(λ)是CO₂通量的线性方程(1-1)；λ是线性比例系数，是陆表碳通量、海洋碳通量的调整因子，集合卡尔曼平滑方法通过对λ进行优化，降低模拟浓度与观测浓度的误差，达到优化CO₂通量的目的；

A是与CO₂浓度初始场y相关的单位矩阵；

把公式(1-3)的矩阵G、A合并为H＝(G,A)，把CO₂通量x＝f(λ)、浓度初始场y合并为z＝(x^Ty^T)^T，则公式(1-3)可转化为y⁰＝H(z)+ξ。

用CO₂通量的线性方程(1-1)x＝f(λ)对公式z＝(x^Ty^T)^T进行简化，则公式(1-3)转化为：

y⁰＝H(x)+ξ             (1-4)

式中，H为观测算子矩阵，用于将通量数据转为CO₂浓度，并根据观测数据的时空信息对模拟浓度进行采样；观测算子矩阵H如下所示：

H_ij代表的是像元点i、地基站点浓度j的观测算子；W_ij代表的是像元点i、卫星柱浓度j的观测算子；

不同于地基站点观测算子H_ij，对于卫星浓度观测算子W_ij，首先需要利用H_ij将CO₂通量数据转为CO₂模拟浓度，再结合公式(1-2)将CO₂模拟浓度转化为CO₂模拟柱浓度；

公式(1-4)的代价函数J：

$J=\frac{1}{2}{(y^0-H\left(x\right))}^T{R}^{-1}(y^0-H\left(x\right)))+\frac{1}{2}{(x-x_0)}^T{P}^{-1}(x-x_0)---(1-5)$

式中，y⁰是CO₂观测浓度，R＝var(ζ)是观测浓度误差协方差矩阵；x₀是先验CO₂通量，P是先验CO₂通量的误差协方差矩阵；x是状态变量，代表同化模型中要求解的CO₂通量；上标T代表矩阵转置；通过最小化代价函数J，可以获取CO₂通量x的最优化估测：

x＝x₀+k(y⁰-H(x₀))                 (1-6)

式中，k是一个增益矩阵，表示如下方程：

k＝(PH^T)/(HPH^T+R)                 (1-7)

CO₂通量x的分析误差P^a表示如下：

P^a＝(I-kH)P                     (1-8)

公式1-6、1-8所示即CO₂通量和误差的最优化估测。

下面结合附图和实施例对本发明作进一步详细描述：

步骤1 数据准备

收集和整理驱动大气传输模式运行的先验通量数据(包括陆表碳通量、海洋碳通量、人为碳排放以及野火燃烧碳排放)、气象数据及同化所需的CO₂观测数据。

其中，站点浓度观测数据主要来自美国海洋与大气总署(NOAA，National Oceanic and Atmospheric Administration)的地球系统研究实验室(ESRL，Earth System Research Laboratory)及澳大利亚联邦科学与工业研究组织(CSIRO，Commonwealth Scientific and Industrial ResearchOrganization，Australia)的瓶装采样的flask数据。表2-1列出本研究采用的地基CO₂浓度观测数据。

卫星观测数据，选用从2009.6到2010.12(GOSAT只从2009.6才有连续的观测数据)的GOSAT-ACOS3.3(Version3.3 of the Atmospheric CO₂Observations from Space)中的L2数据产品作为卫星-地基联合同化系统的卫星观测数据。

气象数据是收集和整理2009-2010年间欧洲中期天气预报中心(ECMWF，European Centre for Medium-Range Weather Forecasts)的风、湿、温、压等气象数据，并处理为全球1×1°、3小时的标准数据，为大气CO₂浓度和通量模拟提供气象驱动数据。

步骤2 卫星-地基联合同化方案

以大气反演模型CarbonTracker(Peters，et al.，2005)为例，在此模型的工作基础上，构建大气反演同化系统，采用集合卡尔曼平滑方法发展卫星-地基联合同化方法。这里对卫星-地基CO₂观测浓度同化过程进行阐述：

2.1 CO₂观测浓度同化方法的设计

首先给出碳源/汇同化的过程控制方程：

式中各参数含义如下：

x_t：状态变量，M_t-1:预报算子，y⁰：CO₂浓度观测值，H_t：观测算子，ξ：观测误差，η：预报模型误差；下标t代表时间。

公式2-1所示的两个过程方程分别用来描述t时刻状态变量CO₂通量的预报过程和CO₂浓度值的模拟过程。根据公式2-1，可给出卫星-地基观测浓度同化方法的代价函数J：

$J=\frac{1}{2}{(y^0-H\left(x\right))}^T{R}^{-1}(y^0-H\left(x\right)))+\frac{1}{2}{(x-x_0)}^T{P}^{-1}(x-x_0)---(2-2)$

式中，上标T代表矩阵转置；y⁰为卫星和地基浓度观测值，R是y⁰的误差协方差矩阵(一般根据先验知识，给定观测误差)，R＝var(ζ)；x₀是先验CO₂通量，P是x₀的背景误差协方差矩阵(给定状态变量背景误差协方差矩阵)；x是状态变量，代表要求解的CO₂通量，而x是关于线性系数λ的函数方程：

$x=f\left(\mathrm{\λ}\right)=\underset{r=1}{\overset{N_{\mathrm{eco}}}{\mathrm{\Σ}}}{\mathrm{\λ}}_r^{\mathrm{eco}}{F}_{\mathrm{bio}}(x,y,t)+\underset{r=1}{\overset{N_{\mathrm{oce}}}{\mathrm{\Σ}}}{\mathrm{\λ}}_r^{\mathrm{oce}}{F}_{\mathrm{oce}}(x,y,t)+F_{\mathrm{ff}}(x,y,t)+F_{\mathrm{fire}}(x,y,t)$

因此，λ才是公式2-2最终要求解的状态变量，式中F_bio和F_oce分别代表陆地和海洋的先验通量，其时间分辨率为3小时，空间分辨率为1×1°；F_ff和F_fire分别代表给定的人为碳排放源和野火燃烧碳排放源，其时间分辨率为月，空间分辨率为1×1°；H是观测算子，这里的观测算子是大气传输模式TM5，用以模拟CO₂浓度值并根据观测数据的空间、时间信息来对模拟浓度进行采样，为同化提供数据准备。与一般的数据同化方法一样，在同一个同化过程中，卫星-地基CO₂观测浓度联合同化方法分两个步骤来求解最优CO₂通量：(1)状态分析过程和(2)状态预报过程。以下对这两个步骤进行介绍。

2.1.1 状态分析

求解的状态变量λ和它的协方差Q的计算公式如下：

${\mathrm{\λ}}_t^a={\mathrm{\λ}}_t^b+K(y^o-H\left(f\left({\mathrm{\λ}}_t^b\right)\right))---(2-3)$

$Q_t^a=(I-\mathrm{KH})Q_t^b---(2-4)$

式中，下标t代表时间，上标b代表背景值，而上标a代表分析值；为t时刻CO₂通量线性系数的分析场，为t时刻CO₂通量线性系数的背景场，注意：我们在公式(2-3)中用 $x=f\left(\mathrm{\λ}\right)=\underset{r=1}{\overset{N_{\mathrm{eco}}}{\mathrm{\Σ}}}{\mathrm{\λ}}_r^{\mathrm{eco}}{F}_{\mathrm{bio}}(x,y,t)+\underset{r=1}{\overset{N_{\mathrm{oce}}}{\mathrm{\Σ}}}{\mathrm{\λ}}_r^{\mathrm{oce}}{F}_{\mathrm{oce}}(x,y,t)+F_{\mathrm{ff}}(x,y,t)+F_{\mathrm{fire}}(x,y,t)$ 来求得t时刻的先验CO₂通量，并通过观测算子H来转化为模拟浓度，因此在公式(2-3)中的模拟浓度为 是t时刻状态变量λ的背景误差协方差矩阵，是t时刻状态变量λ的分析误差协方差矩阵；观测算子H(.)将状态变量投影到观测空间；I表示单位矩阵。

对于地基站点观测数据，表示传输模式TM5预报的CO₂浓度(即CO₂模拟浓度)。而对于GOSAT卫星观测数据，表示传输模式TM5预报的三维CO₂浓度模拟场，需要进行进一步尺度转换才能获得相应的CO₂模拟柱浓度由于GOSAT的平均柱浓度X_CO2是各层CO₂浓度的加权平均，因此在通量反演过程中，必须使用相同的加权平均算子来计算根据Rodgers等(Rodgers C.D.and Connor B.J.2003.Intercomparison of remote sounding instruments.Journal ofGeophysical Research,108,4116)，推导出模拟CO₂的柱平均浓度计算公式如下：

$y_t^b=X_{\mathrm{CO}2}^m=h^T{X}_a+h^{T}A\left(S\right[H\left(f\left({\mathrm{\λ}}_t^b\right)\right)]-X_a)---(2-5)$

其中，和均代表为转化后的CO₂模拟柱浓度；h^T是卫星提供的气压加权函数；X_a是卫星提供的先验浓度廓线；A是卫星提供的平均核函数；S(.)是空间插值算子，它将TM5模拟得到的三维CO₂浓度场插值到GOSAT卫星观测的星下点，得到该点的CO₂垂直廓线。由于同化效果与模式系统和观测资料的质量都密切相关，同化中只挑选了ACOS V3.3中标记为Good，且与TM5模拟所得的平均柱浓度相差小于3ppm的卫星柱浓度数据进入同化过程。为了简单考虑，假设所有的观测误差是不相关的，即卫星观测误差方差矩阵R是值为9(3²)的对角矩阵。

式2-3中，K为卡尔曼增益系数，其计算公式如下：

$K=\left(Q_t^b{H}^T\right)/({\mathrm{HQ}}_t^b{H}^T+R)---(2-6)$

为了使方程求解方便、简洁，在同化系统中对状态向量λ进行细分，把它定义为其均值偏差λ'之和，即：

$\mathrm{\λ}=\stackrel{\‾}{\mathrm{\λ}}+{\mathrm{\λ}}^{\′}---(2-7)$

则状态向量偏差λ'(维度N)每一列的矩阵定义为X：

$\begin{array}{c}X=\frac{1}{\sqrt{N-1}}{({\mathrm{\λ}}_1^{\′},{\mathrm{\λ}}_2^{\′},...............{\mathrm{\λ}}_N^{\′})}^T\\ =\frac{1}{\sqrt{N-1}}{({\mathrm{\λ}}_1-\stackrel{\‾}{\mathrm{\λ}},{\mathrm{\λ}}_2-\stackrel{\‾}{\mathrm{\λ}},...............{\mathrm{\λ}}_N-\stackrel{\‾}{\mathrm{\λ}})}^T\end{array}---(2-8)$

根据公式2-8，方差矩阵可以写成：

Q^a＝XX^T                 (2-9)

为了减少计算量，在计算卡尔曼增益系数K时，对公式中的HQH^T,QH^T进行了简化处理，HQH^T,QH^T可定义为：

${\mathrm{HQH}}^T\≈\frac{1}{N-1}(H\left(f\left({\mathrm{\λ}}_1^{\′}\right)\right),H\left(f\left({\mathrm{\λ}}_2^{\′}\right)\right),...H\left(f\left({\mathrm{\λ}}_N^{\′}\right)\right))\×(H\left(f\left({\mathrm{\λ}}_1^{\′}\right)\right),H\left(f\left({\mathrm{\λ}}_2^{\′}\right)\right),...H\left(f\left({\mathrm{\λ}}_N^{\′}\right)\right){)}^T---(2-10)$

${\mathrm{QH}}^T\≈\frac{1}{N-1}({\mathrm{\λ}}_1^{\′},{\mathrm{\λ}}_2^{\′},......{\mathrm{\λ}}_N^{\′})\×(H\left(f\left({\mathrm{\λ}}_1^{\′}\right)\right),H\left(f\left({\mathrm{\λ}}_2^{\′}\right)\right),...H\left(f\left({\mathrm{\λ}}_N^{\′}\right)\right){)}^T---(2-11)$

对于单独一个CO₂观测值来说，公式2-10变成了2个矩阵的点乘，使HQH^T则变成了一个[1×1]一维的标量，而QH^T则变成了一个[n×1]一维的向量。利用公式2-10、2-11计算出来的HQH^T,QH^T很容易就计算卡尔曼增益系数K，实现公式2-6(卡尔曼增益系数)的求解。

卡尔曼增益系数K是用来更新平均状态向量结合公式2-3，实现平均状态变量分析值的更新。平均状态变量更新的同时，对状态向量的偏差λ'也进行了相应的更新，更新公式如下：

${\mathrm{\λ}}_i^{\′a}={\mathrm{\λ}}_i^{\′b}-\tilde{k}H\left(f\left({\mathrm{\λ}}_i^{\′b}\right)\right)---(2-12)$

公式2-12中，向量与卡尔曼增益系数K之间存在相关性，可用公式表示如下：

$\tilde{k}=K*\mathrm{\α}---(2-13)$

$\mathrm{\α}={(1+\sqrt{\frac{R}{{\mathrm{HQ}}^b{H}^T+R}})}^{-1}---(2-14)$

这种对状态向量和状态向量偏差λ'分别更新的方式，可有效防止系统性低估Q^a量。且α的计算过程十分简单，只需估算出向量R和HQH^T，就可计算出α。

对应于状态向量λ的更新，系统也对取样点处的CO₂浓度模拟值(即：TM5预报的CO₂浓度)进行更新。一般来说，CO₂浓度模拟值更新的最直接方式就是利用更新后的λ的计算出一系列新的CO₂通量，通过重新运行观测算子TM5估算出新的CO₂浓度模拟值，达到取样点处的CO₂浓度模拟值更新的目的。然而，这种通过运行TM5进行浓度更新的方式计算代价十分巨大，相当于CO₂浓度观测值每更新一次状态向量时TM5就要重新运行一次。因此，为了提高运行效率，对模拟浓度数据采取了类似于状态向量和状态向量偏差λ'的更新方式，CO₂浓度模拟值的更新也通过卡尔曼增益系数K来实现。取样点m处的CO₂浓度模拟值的更新包括其均值的更新：

${\stackrel{\‾}{y}}_t^a=H{\left(f\left({\stackrel{\‾}{\mathrm{\λ}}}_t^a\right)\right)}_m=H{\left(f\left({\stackrel{\‾}{\mathrm{\λ}}}_t^b\right)\right)}_m+H_{m}K(y_t^0-H\left(f\left({\stackrel{\‾}{\mathrm{\λ}}}_t^b\right)\right))---(2-15)$

以及其偏差的更新：

$y_t^{\′a}=H{\left(f\left({\mathrm{\λ}}_i^{\′a}\right)\right))}_m=H{\left(f\left({\mathrm{\λ}}_i^{\′b}\right)\right)}_m-H_m\tilde{k}\left(H(f\left({\mathrm{\λ}}_i^b\right)\right)---(2-16)$

在完成状态向量和状态向量偏差λ'和浓度偏差y'更新后，系统将状态向量作为背景值推动同化系统向下一个时刻运行。

2.1.2 状态预报

卫星-地基观测浓度联合同化的另一个重要过程就是状态变量的预报过程(用M来表示)。该过程描述了状态变量在时间上的更新，它以前一时刻的状态变量为背景量，通过预报算子M估测出下一时刻的状态变量及方差矩阵的预报值：

${\mathrm{\λ}}_{t+1}^b=M\left({\mathrm{\λ}}_t^a\right)---(2-17)$

$Q_{t+1}^b=M{Q}_t^a{M}^T---(2-18)$

状态变量λ是CO₂通量的线性系数，即CO₂通量是关于λ的线性方程，由于现在的科学界无法寻找到一个科学的机理模型来描述CO₂通量的预报过程，进而无法找到状态变量λ的预报过程。我们按照现有国际的惯例，简化状态变量λ的预报过程，用一个单位矩阵I来代替预报算子(即前一时刻的状态变量分析值是下一时刻的状态变量背景值)，用以预报出下一时刻的状态。

步骤3 模型的运行

在模型运行之前需要对各种参数进行设定。

3.1 模型输入与输出设定

在模型运行之前，对模型输入及输出项进行设定。模型的驱动数据主要包括CO₂站点浓度及卫星观测数据、4个先验通量和气象数据，数据的来源及时空分辨率见表3-1。其中，风、湿、温、压等气象数据以及陆地碳通量、海洋碳通量、人为碳排放以及野火燃烧碳排放等CO₂通量数据主要用来驱动大气传输模式TM5的运行，为数据同化提供模拟浓度；CO₂站点及卫星浓度数据主要是数据同化过程提供浓度观测值，通过最小观测与模拟浓度之差，最终实现CO₂浓度和通量的优化。

表3-1 输入数据信息表

*Takahashi et al.[2009]^a:Takahashi,T.,S.C.Sutherland,C.Sweeney,R.A.Feely,D.W.Chipman,B.Hales,G.

Friederich,F Chavez,C.Sabine,A.Watson,D.C.E.Bakker,U.Schuster,N.Metzl,H.Yoshikawa-Inoue,M.Ishii,T.

Midorikawa,Y.Nojiri,A.Kortzinger,T.Steinhoff,M.Hoppema,J.Olaffson,T.S.Anarson,B.Tilbrook,T.

Johannessen,A.Olsen,R.Bellerby,C.S.Wong,B.Delille,N.R.Bates,and H.J.W.de Baar(2009),Climatologicalmean and decadal change in surface ocean pCO,and net sea-air CO2flux over the global oceans,Deep-Sea Resarch II,56,554--577.

EC：Environment Canada，加拿大环境部

NCAR：National Center for Atmospheric Research，美国国家大气研究中心

LBNL：Lawrence Berkeley National Laboratory，劳伦斯伯克利国家实验室

NOAA：National Oceanic and Atmospheric Administration，美国国家海洋和大气管理局

ESRL：Earth System Research Laboratory，地球系统研究实验室

CDIAC：Carbon Dioxide Information and Analysis Center，二氧化碳信息分析中心

GFEDv2：Global Fire Emissions Database version 2，全球火灾排放数据库版(v2)

CSIRO：Commonwealth Scientific and Industrial Research Organization，Australia,澳大利亚联邦科学与工业研究组织

模型输入项设定了模式运行的起始时间、模型运行频率及滞后窗口见表3-2。

表3-2 模型输入参数设定

为了同时获取全球生态系统和中国陆地生态系统时空分布特征，本例中嵌套设定了3个不同空间分辨率研究区，用以综合模拟分析全球生态系统碳源汇时空特征，重点分析中国陆地生态系统碳源汇的变化趋势，确定中国对全球生态系统碳汇的贡献。本例中以全球(6×4°)为背景区，在亚洲设定一个一级嵌套区(过渡性质的嵌套区，用来过渡中国嵌套区)，其空间分辨率为3×2°，在亚洲嵌套区的基础上再设立一个二级嵌套区(1×1°)，用来研究中国的陆地碳汇分布情况。

此外，模型也对模型输出项进行了设定，最终输出结果为周步长的、全球(6×4°)、亚洲(3×2°)、中国(1×1°)CO₂浓度及源汇分布图，具体设置见表3-3。

表3-3 模型输出设定

步骤4 模型验证

主要从两个方面对模型进行验证，一是浓度的误差分析，包括模拟(同化后)与观测浓度的对比分析和模拟(独立验证)与观测浓度的对比分析，二是对全球碳源汇的模拟与分析。此外，本发明还对大气中不同高度的CO₂浓度进行了单独验证。

4.1 浓度的误差分析

4.1.1 模拟(同化后)与观测浓度的对比

用同化后的10,317个CO₂浓度模拟值与观测值进行对比，分析浓度/通量的模拟效果，如图3和图4所示，浓度模拟值与观测值之间的误差为-0.02±1.83ppm(R＝0.88，P<0.05)，小于2个ppm，满足碳同化模型的要求，说明卫星-地基CO₂联合同化系统能够比较合理模拟出现实中浓度/通量的时空分布。

4.1.2 模拟(独立验证)与观测浓度的对比

除了对同化后的浓度进行验证，还利用站点及飞机观测数据进行了独立验证。所谓的“独立验证”是指，不同于图3和图4，观测数据都没有被卫星-地基CO₂联合同化系统同化过。

图5显示了CO₂浓度观测与模拟值在站点WSA间(“独立验证”)的对比分析图。总体来说，站点WSA的浓度模拟值与观测值之间的误差为-0.47±1.87ppm，小于2个ppm，说明卫星-地基CO₂联合同化系统模拟的浓度/通量基本能够描述出现实世界的浓度/通量时空分布特征。

此外，用CONTRAIL飞机观测数据对大气中不同高度的CO₂浓度进行了单独验证。图6结果表明，485-525hPa、375-425hPa及225-275hPa三层CO₂的模拟结果良好，与观测浓度的相关系数分别为R＝0.89(P<0.05)、R＝0.86(P<0.05)及R＝0.84(P<0.05)。同时也发现，随着大气压的降低(与海拔相反，其气压的降低表明海拔在升高)，CO₂浓度的模拟精度在逐渐下降，其原因有两个：一是越靠近地表，CO₂浓度观测控制点越多，有效的提高了CO₂的模拟精度；二是大气传输模型的垂直扩散过程存在很多不确定性，大大影响了CO₂的垂直模拟结果。尽管CO₂的模拟还存在一系列问题，其模拟的精度受大气传输、扩散过程影响，垂直模拟数字化过程还有待进一步完善，但不管是站点的模拟(如图5中的WSA站点)还是垂直分层对比(图6的CONTRAIL数据的垂直模拟)，其模拟值与观测值十分接近，进一步说明卫星-地基CO₂联合同化系统能够比较合理地模拟出现实中CO₂浓度/通量的时空分布。

步骤5 碳源汇的结果分析

5.1 全球碳源汇的模拟与分析

表5-1 全球按TransCom 23个区的先验、后验碳通量及不确定统计结果(2010年平均)

注：flux Diff(差异)＝Post.Flux(后验通量)-Prior Flux(先验通量)

uncertainty Reduced rate％(通量不确定性减少率％)＝(Prior Flux Uncertainties(先验通量不确定性)-Post.Flux Uncertainties后验通量不确定性))/Prior Flux Uncertainties(先验通量不确定性)*100％；

未处理区(Non-optimized)指模型系统数据缺失、无法优化的区域，其所占全球比例很小(1％)。

表5-1给出了全球TransCom 23个区的先验、后验碳通量及不确定性统计结果。由于CO₂观测的控制及同化作用，全球(Global Total)34.45％的先验碳通量的不确定性被降低。

图7给出了2010年年平均全球陆地、海洋、野火燃烧及化石燃料碳通量分布的柱状图。由图可见，2010年全球净碳收支量为4.72Pg C/yr，其中全球陆地收支量为-3.81Pg C/yr，全球海洋收支量为-2.88Pg C/yr，野火燃烧释放的碳源为2.37Pg C/yr，化石燃烧释放的碳源为9.03Pg C/yr。2010年全球净碳收支量估测结果(4.72Pg C/yr，其中陆地碳汇为-1.44Pg C/yr，海洋碳分量为-2.88Pg C/yr)与CT2011(5.23Pg C/yr，其中陆地碳汇为-1.41PgC/yr，海洋碳分量为-2.49Pg C/yr，http://www.esrl.noaa.gov/gmd/ccgg/carbontracker)及CTE(5.37Pg C/yr，其中陆地碳汇为-1.39Pg C/yr，洋碳分量为-2.17Pg C/yr，http://carbontracker.eu/)的估测结果十分相近。

5.2 中国区碳源汇的模拟与分析

对中国区的陆地碳源汇估测结果进行了更加详细的分析。

5.2.1 中国陆地碳源汇分析

从年平均碳源汇分布分析结果表明(如图8和表5-3所示)：中国陆地生态系统是一个大气碳汇(用深灰色表示)，2010年间吸收了-0.29Pg C。我国森林生态系统吸收碳汇量为-0.13Pg C/yr，占中国陆地碳汇总量的44.63％。其中，-0.03Pg C/yr由针叶林吸收，-0.03Pg C/yr由阔叶林吸收，-0.05Pg C/yr由混交林吸收，-0.02Pg C/yr由其它林种吸收。中国农田生态系统的CO₂吸收量为-0.06Pg C/yr，占中国陆地碳汇总量的20.60％。其中，绝大部分农田碳汇分布在中国的华北、华中、东南。2010年，估测的农田碳汇的单位面积碳吸收能力为37g C/m²/yr。这种碳吸收能力与我国农地的作业方式及农产品栽培技术相关。因此，卫星-地基联合同化系统可能因为无法捕捉到农作物收割及农产品消耗这种由于人类活动所引起的“侧向传输过程”而高估了农田碳汇量。2010年，估测的草地/灌木林吸收的CO₂量为-0.08Pg C/yr，占中国陆地碳汇总量的27.23％。草地/灌木碳汇的单位面积碳吸收量为30g C/m²/yr，其中，最强的草地/灌木碳吸收量出现在内蒙古东部。除森林、草地/灌木及农田生态系统外，其它类型的陆地生态系统(如东北地区的半苔原Semitundra)的CO₂吸收量为-0.02Pg C/yr，占中国陆地碳汇总量的6.91％，其单位面积碳吸收量为2g C/m²/yr。

表5-3 中国陆地生态系统不同植被的碳汇统计表(2010年平均)

5.2.2 卫星观测对中国陆地碳汇的影响

为了测试卫星观测数据对中国陆表碳浓度/通量的影响，设置了一个卫星观测数据对比实验。图9给了“卫星-地基”和“地基”的两个实验结果的标准差分布图。从图中可以看出，“卫星-地基”实验中由于卫星观测数据的加入，导致中国陆表通量的标准差分布大大不同：相对于不加卫星观测数据的“地基”实验，“卫星-地基”实验的标准差大大降低，且这种误差的减少出现明显的季节性变化，最大误减小出现在在夏季8月份。这是由于额外观测数据GOSAT的加入，进一步提高了陆地碳汇的估测精度，而这种精度的提高在植被生长季的剧烈CO₂变化过程中被放大，使之表现为明显的不确定性降低。

图10给定了“卫星-地基”实验与“地基”实验在MNM站点间CO₂浓度观测和模拟值的对比分析图。从图10(a)及(b)中可以看出，观测数据设置不同，两个实验在站点MNM的浓度观测与模拟值之间的相关系数也不同。“卫星-地基”实验的浓度观测与模拟值间的相关系数明显高于不加卫星观测数据的实验，这与图9给定的结果相一致。然则图10(c)及(d)的直方图分析，“卫星-地基”实验的误差分析(0.61ppm)却大于“地基”实验(0.24ppm)，这是由于在统计过程中，正负误差的相互加减令平均误差不能完全描述出观测与模拟间的差异情况；而均方根误差比平均误差更加科学地描述出模拟值与观测值之间的偏差。在图10中(c)及(d)中，“地基”实验的均方根误差RMS(±0.79)大于加卫星观测数据的实验(±0.78)，这也和图10中(a)、(b)的结果相一致：额外观测数据(GOSAT)提高了中国碳汇的估测精度。

表5-4 卫星观测数据对中国陆地碳汇的影响(2010年平均)

通过对比“卫星-地基”联合同化实验、“地基”实验的CO₂通量结果来分析碳卫星观测数据GOSAT对中国陆地碳汇的影响。表5-4分别列出了“卫星-地基”实验、“地基”实验的中国陆地碳汇先验、后验值、两者之差及它们的标准差(G-uncertainties)、标准差降低百分比(＝(先验误差-后验误差)/先验误差×100％)。“卫星-地基”实验、“地基”实验的中国(后验)陆地碳汇估测结果分别为-0.29和-0.21Pg C/yr，两者相差0.08Pg C/yr，这就表明碳卫星观测数据GOSAT的加入提高了2010年中国陆地生态系统碳汇强度。

图11显示了2010年中国陆地生态系统在“卫星-地基”实验和“地基”实验间碳汇分布差异图。从图中可以看出，碳卫星观测数据GOSAT的加入令“卫星-地基”实验的陆地碳汇分布发生了重大改变。注意，在图中深灰色表示“卫星-地基”实验陆地碳汇的增加，白色表示“卫星-地基”实验陆地碳汇的减少。陆地碳汇的增加区为中国的华南、东南、内蒙和西南的大部分地区，而东北、华北及华中地区则表现为明显的陆地碳汇的减少。西北及青藏高原区的碳汇增加不是十分明显。与“地基”实验相比，“卫星-地基”实验的不确定性(G-uncertainties)也比“地基”实验降低了17％(表5-3，17％＝43％-26％)，这表明：GOSAT观测数据的加入将会进一步的优化陆地碳汇的估测结果，提高模拟的精度。因此，碳卫星观测数能够提高系统对中国陆地生态系统CO₂通量的估测精度和强度。然则，从前人的实验结果也得知，碳卫星观测本身也存在局限性，有些GOSAT观测数据本身的不确定性就比较大。

5.2.4 其它结果的对比

表5-5列出了基于碳同化方法的中国陆地碳汇估测结果的对比表。从表中可以看出，本研究估测的2010年的陆地碳汇与CT2011(-0.27PgC/yr)、Jiang et al.，(Jiang F,Wang H W,Chen J M,et al.2013.Nestedatmospheric inversion for the terrestrial carbon sources and sinks in China.Biogeosciences,10,5311-5324.)结果(-0.28±0.18Pg C/yr，2002-2008平均)十分相近，但略高于Zhang et al.,(Zhang H F,Chen B Z,van der Laan-LuijkxI T,et al.2014b.Estimating Asian terrestrial carbon fluxes from CONTRAILaircraft and surface CO₂ observations for the period 2006-2010.AtmosphericChemistry and Physics,(11)：5807-5824)的估测结果(-0.22±0.36Pg C/yr)，略低于Piao et al.，(Piao S,Fang J,Ciais P,et al.2009.The carbon balance ofterrestrial ecosystems in China.Nature,458,1009-1013)的估测结果(-0.28±0.18Pg C/yr，2002-2008平均)。

表5-5 中国陆地碳汇估测结果的对比表(Pg C/yr)

^aCT2011的结果下载于ftp://aftp.cmdl.noaa.gov/products/carbontracker/co2/fluxes/monthly/

综上所述，由于采用了以上技术方案，本发明成功构建了基于卫星柱浓度及地基站点观测数据为基础的卫星-地基CO₂联合同化方法，改进大气反演模型只能同化站点观测的缺陷，在原来的地基站点同化方法上引进柱浓度同化方案，实现多种来源CO₂数据的联合同化，为精确的区域碳源/汇估算提供有效手段。

以上所述，仅是本发明的较佳实施例而已，并非对本发明作任何形式上的限制，本领域技术人员利用上述揭示的技术内容做出些许简单修改、等同变化或修饰，均落在本发明的保护范围内。

Claims (9)

1.基于卫星-地基CO₂数据联合同化的碳源汇估测方法，其特征在于，包括以下步骤：

根据先验CO₂通量数据驱动大气传输模型，获得CO₂模拟浓度时空数据集；

采用集合卡尔曼平滑方法同化CO₂观测数据，对所述CO₂观测数据与CO₂模拟浓度之差最小化求解，获得优化的CO₂浓度及CO₂通量的时空分布场；

其中，所述CO₂观测数据包括CO₂站点浓度数据及CO₂卫星柱浓度数据：

当所述CO₂观测数据为CO₂站点浓度时，从所述CO₂模拟浓度时空数据集中提取时空信息匹配的CO₂模拟浓度，对该CO₂模拟浓度与所述CO₂站点浓度之差最小化求解；

当所述CO₂观测数据为CO₂卫星柱浓度时，从所述CO₂模拟浓度时空数据集中提取时空信息匹配的CO₂模拟浓度，先将该CO₂模拟浓度转化成CO₂模拟柱浓度，再对CO₂卫星柱浓度与转化后的CO₂模拟柱浓度之差最小化求解。

2.根据权利要求1所述的基于卫星-地基CO₂数据联合同化的碳源汇估测方法，其特征在于，所述CO₂站点浓度数据来自美国海洋与大气总署的地球系统研究实验室及澳大利亚联邦科学与工业研究组织的瓶装采样的flask数据；所述CO₂卫星柱浓度数据来自GOSAT-ACOS3.3中的L2数据产品。

3.根据权利要求1所述的基于卫星-地基CO₂数据联合同化的碳源汇估测方法，其特征在于，所述进行同化过程的CO₂卫星柱浓度数据选取GOSAT-ACOS V3.3中标记为Good，且满足以下条件的数据：卫星柱浓度数据与其对应的转化后的CO₂模拟柱浓度相差小于3ppm。

4.根据权利要求1所述的基于卫星-地基CO₂数据联合同化的碳源汇估测方法，其特征在于，所述同化过程是通过大气反演模型进行的。

5.根据权利要求1所述的基于卫星-地基CO₂数据联合同化的碳源汇估测方法，其特征在于，具体包括以下步骤：

S1：根据碳通量调整系数λ_t的背景场利用方程1-1计算先验CO₂通量数据，即CO₂通量预报场；

x_(t)＝f(λ_t)＝[(F_bio,t+F_ocn,t]×λ_t+F_ff,t+F_fire,t    (1-1)

式中，x_(t)代表CO₂通量，F_bio,t、F_ocn,t、F_ff,t、F_fire,t分别代表陆表碳通量、海洋碳通量、人为碳排放以及野火燃烧碳排放，λ_t为碳通量调整系数，下标t代表时间，当t＝0时，λ_t为初始值1.0；

S2：根据所述先验CO₂通量数据，驱动大气传输模式运行，生成CO₂模拟浓度时空数据集；

S3：采用集合卡尔曼平滑方法对时空信息匹配的CO₂观测浓度与CO₂模拟浓度之差最小化求解，获得优化的碳通量调整系数λ_t ^a，即λ_t的分析场，再根据λ_t ^a利用方程1-1计算优化的CO₂通量，即CO₂通量分析场；

S4：根据优化的CO₂通量再次驱动大气传输模式，获得优化后的CO₂模拟浓度时空数据集；

S5：将优化的碳通量调整系数λ_t ^a作为下一时刻的碳通量调整系数的背景场，返回步骤S1。

6.根据权利要求5所述的基于卫星-地基CO₂数据联合同化的碳源汇估测方法，其特征在于，所述陆表碳通量F_bio,t是通过陆地生态系统模型获得的。

7.根据权利要求1所述的基于卫星-地基CO₂数据联合同化的碳源汇估测方法，其特征在于，所述将CO₂模拟浓度转化成CO₂模拟柱浓度是通过公式1-2进行的：

$X_{\mathrm{CO}2}^m=h^T{X}_a+h^{T}A(X_h-X_a)---(1-2)$

式中，是转化后的CO₂模拟柱浓度；h^T是卫星提供的气压加权函数；X_a是卫星提供的先验廓线；A是卫星提供的平均核函数；X_h是与卫星提供的先验浓度廓线X_a时空尺度相匹配的、大气传输模式模拟的模拟浓度廓线。

8.根据权利要求7所述的基于卫星-地基CO₂数据联合同化的碳源汇估测方法，其特征在于，所述模拟浓度廓线X_h通过以下方式实现与所述先验浓度廓线X_a时空尺度相匹配：

大气传输模式首先根据卫星柱浓度数据的时间、空间信息对CO₂模拟浓度进行采样，获得时空信息相一致的CO₂模拟浓度，再根据先验浓度廓线X_a的分层高度对CO₂模拟浓度进行插值，以生成与先验浓度廓线X_a时空尺度相匹配的模拟浓度廓线X_h。

9.根据权利要求1所述的基于卫星-地基CO₂数据联合同化的碳源汇估测方法，其特征在于，所述同化过程通过最小化代价函数J获取要求解的CO₂通量x及其分析误差P^a的最优化估测：

x＝x₀+k(y⁰-H(x₀))        (1-6)

P^a＝(I-kH)P            (1-8)

其中，k是一个增益矩阵，表示如下方程：

k＝(PH^T)/(HPH^T+R)        (1-7)

所述代价函数J：

$J=\frac{1}{2}{(y^0-H\left(x\right))}^T{R}^{-1}(y^0-H\left(x\right)))+\frac{1}{2}{(x-x_0)}^T{P}^{-1}(x-x_0)---(1-5)$

式中，y⁰是CO₂观测浓度，R＝var(ζ)是观测浓度误差协方差矩阵；x₀是先验CO₂通量，P是先验CO₂通量的误差协方差矩阵；上标T代表矩阵转置；

H为如下所示的观测算子矩阵：

其中H_ij代表像元点i、地基站点浓度j的地基站点观测算子；W_ij代表的是像元点i、卫星柱浓度j的卫星浓度观测算子。