# Daily CO2 flux estimates over Europe from continuous atmospheric measurements: 1, inverse methodology

P. Peylin, P. J. Rayner, Philippe Bousquet, C. Carouge, F. Hourdin, P.

Heinrich, Philippe Ciais

► To cite this version:

P. Peylin, P. J. Rayner, Philippe Bousquet, C. Carouge, F. Hourdin, et al. Daily CO2 flux estimates over Europe from continuous atmospheric measurements: 1, inverse methodology. Atmospheric Chemistry and Physics, 2005, 5 (12), pp.3173-3186. 10.5194/acp-5-3173-2005. hal-00328413

HAL Id: hal-00328413

https://hal.science/hal-00328413v1

Submitted on 18 Jun 2008

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Daily CO₂ flux estimates over Europe from continuous atmospheric measurements: 1, inverse methodology

P. Peylin<fnref n="1" />, P. J. Rayner<fnref n="2" />,<fnref n="3" />, P. Bousquet<fnref n="3" />, C. Carouge<fnref n="3" />, F. Hourdin<fnref n="4" />, P. Heinrich<fnref n="5" />, P. Ciais<fnref n="3" />, and AEROCARB contributors*

A. Adolphsen, F. Apadula, R. Graul, L. Haszpra, M. Ramonet, R. Santaguida, M. Schmidt

Received: 21 December 2004 – Published in Atmos. Chem. Phys. Discuss.: 18 March 2005

Revised: 16 September 2005 – Accepted: 1 November 2005 – Published: 25 November 2005

**Abstract.** This paper presents an inverse method for inferring trace gas fluxes at high temporal (daily) and spatial (model grid) resolution from continuous atmospheric concentration measurements. The method is designed for regional applications and for use in intensive campaigns. We apply the method to a one month inversion of fluxes over Europe. We show that the information added by the measurements depends critically on the smoothness constraint assumed among the source components. We show that the initial condition affects the inversion for 20 days, provided one has enough observing sites to constrain regional fluxes. We show that the impact of the far-field fluxes grows throughout the inversion and hence a reasonable global flux field is a prerequisite for a regional inversion.

## 1 Introduction

The task of determining the space-time structure of carbon fluxes to the atmosphere is one step in any attempt to monitor and possibly manage the carbon cycle. The task has generally been performed via two complementary approaches. In the so-called bottom-up methods, small-scale flux estimates are aggregated together to form regional totals. If the measurements are not spatially dense, the approach uses some kind of statistical or physical model to fill gaps. Spatially dense measurements are never of carbon fluxes directly (e.g. land cover) so a different kind of model (e.g. Potter et al.,

2003) is used to relate these measurements to carbon fluxes on some time and space scale.

The other approach, the so-called top-down or inverse approach, infers the space-time flux distribution from concentration signatures in the atmosphere. The approach faces many of the usual difficulties of inverse problems, principally a lack of concentration observations and reliance on uncertain atmospheric transport. A detailed explanation of the underlying principles is given in Enting (2002) and especially of the “matrix approach” used in this study.

The two approaches have different characteristics in almost all respects. Most importantly they give rise to different kinds of uncertainty. For example, if the bottom-up approach must use some kind of extrapolation, then any error in a point measurement will be propagated by the extrapolator and bias large-scale estimates. This is even clearer if some kind of physical model is used to relate the measurements to the carbon fluxes we seek; errors in the model are not random and so will not disappear as we move to larger and larger regions.

The case for the inverse approach is quite different. Here the sparsity of concentration data limits the resolution of the inverse procedure, so that estimates are more certain at large scales. In fact the ability of the method to resolve small scales is limited and most studies employ some form of regularization method. Regularization imposes some extra constraint on the solution to limit its sensitivity to individual data points and hence limit error amplification. Bayesian methods, in which prior estimates of fluxes are inserted as extra data into the problem, have been the norm (e.g. Enting et al., 1995; Rayner et al., 1999; Bousquet et al., 2000; Gurney et al., 2002) but some studies have employed other forms

<footnote>*Experimentalists that contributed to the atmospheric measurements:</footnote>

<footnote n="1">¹Laboratoire de Biogéochimie isotopique, CNRS-UPMC-INRA, Paris, France</footnote>

<footnote n="2">²CSIRO Atmospheric Research, Aspendale, Victoria, Australia</footnote>

<footnote n="3">³Laboratoire des Sciences du Climat et de l'Environnement, CEA, Gif-sur-Yvette, France</footnote>

<footnote n="4">⁴Laboratoire de Météorologie Dynamique, CNRS, Paris, France</footnote>

<footnote n="5">⁵Laboratoire DASE/LDG, CEA, Bruyeres-le-Chatel, France</footnote>

<footnote>*Correspondence to:* P. Peylin
(peylin@lsce.saclay.cea.fr)</footnote>
of regularization or have studied the impact of a choice of regularization methods, e.g. Fan et al. (1999); Baker (2001).

A crude alternative to regularization is to solve for fluxes only in a greatly reduced solution space. This is normally accomplished by solving only for broad-scale integrals of fluxes (i.e. large regions), with smaller scale structure being prescribed as part of the problem set-up. This was the approach chosen in the important study of Fan et al. (1998). Trampert and Snieder (1996) and Kaminski et al. (2001) demonstrated the pitfalls in this approach and Gurney et al. (2002) and Baker (2001) demonstrated the dependence of the results in Fan et al. (1998) on the choice.

Behind all these problems with regularization methods lies the lack of atmospheric concentration data. The most common source of data is a network of approximately 120 stations globally at which flask samples are taken at frequencies from weekly to monthly. The data is usually used in the consolidated form provided by GLOBALVIEW-CO₂ (2002), smoothed at roughly monthly frequencies although reported as weekly values. The network has a strong bias to the marine boundary layer or the free troposphere, reflecting both the historical motivation for the measurements and the difficulties of interpreting samples with strong continental influence. This dilutes the capacity of the measurements to constrain continental fluxes. Also the temporal smoothing of the data removes much of the signal tying concentration signatures to their antecedent fluxes. This information is contained in the tracer labels of individual air masses and varies on synoptic rather than monthly frequencies. It is this more detailed tracer information we intend to use in this study.

In a recent series of studies, Law et al. (2002, 2003, 2004) have pointed out the potential power of using data measured at higher time frequencies. The resolving power of the data is improved not merely because monthly mean concentration is better sampled but more importantly because the synoptic variation in flow acts as a differential sampling tool; a fixed sampling location constrains fluxes from several surrounding regions depending on flow characteristics. We note that this requires that variations in fluxes be predictable on timescales longer than the synoptic.

Law et al. (2002), in their observation system simulation experiment, showed that known fluxes could be recovered with reasonable biases and uncertainties provided firstly that the transport was well-known and secondly that the solution space for fluxes was detailed enough to represent the true fluxes reasonably. Law et al. (2003) showed, perhaps surprisingly, that the requirement for perfect transport could be considerably relaxed provided the data were smoothed to synoptic frequencies.

The main recommendation from these two studies was future deployment of instruments capable of making and reporting these high frequency measurements. There is, however, a considerable body of existing data from such instruments stretching back several decades. Some attempt has already been made to use this data to recover flux information.

The most common technique (e.g. Biraud et al., 2000, 2002) has been to use a so-called ratio or fingerprint approach in which an ancillary tracer is used to provide information on the synoptic transport and the ratio of signals at an observing site is used to solve for the unknown flux distribution. The requirements for such a method are stringent. Firstly the ancillary tracer must have a well-known distribution and secondly the distributions of the known and target fluxes must be very similar.

In this work we make a first attempt to use some of this high-frequency data directly in a “matrix” inversion. The work should be regarded as an exploration of some of the problems which will arise in the use of such data. The usual approach in such a preliminary study is to use synthetic data, that is data generated from a model and usually perturbed by some random noise. It is likely that the major problems in this kind of inversion will revolve around the veracity of the transport model and the diagnostics of quality we use for the inversion. Both of these will only become apparent with the use of real data and this is the approach taken here. The data density we use is not high enough yet to make confident statements about the flux distributions we infer.

There are computational, operational and conceptual issues which arise as we begin to use high-frequency concentration data. First we must address the problem of the spatial and temporal domain and resolution of our solution. The number of parameters (flux components) for which we can solve is limited by computational constraints. Often we wish to extract as much detail as possible over a target region. The problem arises of how to treat the rest of the world. Restricting our domain at all requires us to specify lateral boundary fluxes. In their study of the Mauna Loa site, Vukićević and Hess (2000) noted the high sensitivity of concentrations measured at the observing sites to fluxes at the lateral boundaries, as also pointed out in Chevillard et al. (2002). An inversion which neglects much of the world must deal with this problem in some way. Here we use the alternative approach of reducing the resolution of the solution far from the target region but retaining the global domain.

Often the high-frequency measurements we use in an inversion are gathered during intensive campaigns. Unlike the more usual concentration measurements from monitoring stations, these campaigns are usually too short to allow a spin-up of the concentration field. Thus the structure of the initial concentration field is important and we must include it in our problem somehow (see for instance Gerbig et al., 2003).

At a more conceptual level is the relationship between resolution (spatial and temporal) of the fluxes and the measured concentrations. Part of the power of the high-frequency measurements relies on ascribing high-frequency concentration variations to the impact of somewhat “well-known” high-frequency variations in transport convolved with slowly varying fluxes. More strictly, we require that time variations in fluxes be well-known on the time-scale of the concentration
measurements. The difference between “well-known” and “constant” in the above sentence is important. For example we do not need to solve for hourly fluxes just because there is a large diurnal cycle, but we must take the cycle into account somehow if we wish to use diurnally varying concentrations. However the use of source patterns (usually called basis functions) with large extent in the time domain runs the same risk of aggregation errors (Peylin et al., 1999), as noted for the space domain by Kaminski et al. (2001). Here we are faced with a serious trade-off. The higher the time resolution of the fluxes for which we solve, the less information will be available from the concentration field to constrain each component. However the more we average our flux components in advance, the more likely we are to suffer from aggregation errors in time. The problem has been studied recently by Law et al. (2004). They used an identical twin experiment with fluxes with high time variations and monthly mean basis functions. They concluded that one should account for the diurnal cycle if one wished to use concentrations at resolution higher than one day but that other sub-monthly variations were not a major problem. Here we use daily averaged concentrations and so hope to avoid the problems of the diurnal cycle. In order to explore the impact of varying time resolution of sources, we solve for fluxes at daily resolution and produce averages after the fact.

We cannot completely avoid the impact of the diurnal cycle. The covariance between diurnal cycles of vertical transport and CO₂ fluxes gives rise to nonzero average concentrations. This is the diurnal equivalent (Denning et al., 1996) of the well-known rectifier effect much studied on seasonal timescales (Keeling et al., 1989; Denning et al., 1995). We will not include such an effect in this first theoretical study (which is perhaps justified for November where the diurnal cycle is quite small) and the associated errors will be considered as part of the data uncertainties.

Our experimental design, then, requires us to calculate the relationship of fluxes from every pixel in the model domain and every day to daily averaged concentrations at six stations. In the normal forward model framework, this would require several hundred thousand model runs, which is obviously infeasible. The use of an adjoint model, however, requires only one run per datum rather than per unknown. This is feasible in this case. The construction of adjoint transport models can be a long task and Hourdin and Issartel (2000) have shown that an other approach, using the “retro-transport” tracer equation, is equivalent to the adjoint derivation (for an air-mass weighted scalar product) and can be constructed by changing the sign of some expressions in the forward model. In this study, we use such a “retro-transport” approach with the LMDZ global transport model (Hourdin et al., 2005a) to calculate the required response functions.

The outline of the paper is as follows. In Sect. 2, we describe the various special characteristics of our formulation, including the details of the use of LMDZ, treatment of the initial condition, and the imposition of spatial correlations

among increments in fluxes. Section 3 describes the fit to the data, the contributions of the various components of the space of unknowns to the fit and the reduction of uncertainty and change in fluxes produced in such an inversion. Section 4 discusses the relevance of these findings to various applications of such inversions.

## 2 Inverse method

### 2.1 Inverse framework and atmospheric data

We describe in this section the overall inverse set up that is used to assimilate continuous data in a limited space-time domain. Our domain is roughly Europe for November 1998, and we use 6 continuous sites. The inverse process involves a preliminary independent step, in which we perform a global inversion for one year using monthly data at 64 surface sites but no continuous measurements. In this step, we estimate the monthly magnitudes of surface fluxes for large regions. These fluxes are hence compatible with the global atmospheric record given our transport model. We then perform the inversion of the continuous records themselves, using the previous flux estimates as prior fluxes. We will mainly correct these values and their spatial and temporal distribution for the regions that directly influence our six continuous measurement sites, i.e. Europe and the eastern North Atlantic. For these regions we greatly refine the spatial and temporal resolutions of the inversion compared to the preliminary step. Both inverse steps are based on a “matrix” approach (Enting et al., 1995) and use the Bayesian formalism to regularize the inverse problem.

The preliminary monthly inversion is based on the classical approach (Gurney et al., 2002, e.g.) in which we divide the land and ocean into a few large regions. The set-up follows directly from Bousquet et al. (2000) and Peylin et al. (1999) with 18 land and 14 ocean regions and prior spatial and temporal patterns from the TURC model (Ruimy et al., 1996) for the land biosphere and from Takahashi et al. (2002) for the ocean (Fig. 1). The response at each station to each separate flux is computed by direct integration of the LMDZ zoomed transport model (see Sect. 2.2). Note that fossil emissions are set to fixed monthly values (Andres et al., 1996; Marland et al., 2001). We use atmospheric data from 64 sites (excluding the 6 continuous European sites) taken from GLOBALVIEW-CO₂ (2002). Sites are shown in Fig. 1. The data uncertainties are derived following Bousquet et al. (2000). We will not discuss the results of this “classical” inversion.

It is important to avoid double-counting of information, so for the major inverse step we use only the six continuous records for November 1998 obtained within the AERO-CARB European project (name, longitude, latitude, altitude in meters): Schauinsland (SCH, 7.9°, 47.9°, 1205), Mace Head (MHD, -9.8°, 53.3°, 25), Hegyhatsal tower (HUN,
<figure><img src="image_2.png" /><figcaption>**Fig. 1.** Upper map: contoured regions used for the preliminary monthly inversion and grey regions where fluxes are solved at the model resolution in the second step. The position of the stations are shown: circles for all stations of the preliminary step and stars for the continuous European sites. Lower map: zoomed mesh of the LMDZ model.</figcaption></figure>
16.6°, 46.9°, 82), Westerland (WES, 8.0°, 55.0°, 8), Monte Cimone (CMN, 10.7°, 44.2°, 2165), and Plateau Rosa (PRS, 7.7°, 45.9°, 3480). These sites are shown in Fig. 1.

In this test of methodology, we neglect diurnal variations in both concentration and flux. This overlooks two important complexities. First it assumes that concentrations from all times are equally well simulated in transport models. This choice might be critical in the case of low elevation continental towers (like HUN) where the day-to-day variability of the mean concentration might be primarily driven by the variability of stability of the nocturnal boundary layer rather than by advective effects. Such variability is still poorly represented in current transport models at this resolution, including LMDZ (see Geels et al., 2005¹), which could potentially bias the flux estimates. We also calculate the sensitivities of these diurnally averaged concentrations with respect to diurnally averaged fluxes. Hence we cannot include the impact of the diurnal rectifier inside the inversion framework. We could account for this by precalculating the diurnal rectifier and removing its signal from the data prior to inversion (a so-called presubtraction) but this relies on a very good simulation of the diurnal cycle by TURC of which we were not confident. Extensions of both the sampling and flux distributions will be considered in a future study.

The uncertainty on the daily average is defined as the standard deviation of the hourly concentration averages. This choice is rather crude, but still captures the general feature that model data mismatch is likely to be large for sites and days with large hourly concentration variations. For example, the magnitude of the diurnal rectifier is likely to be larger for larger diurnal cycles. A more rigorous approach will be investigated in a following paper.

As discussed in the introduction, we need to adjust the spatial and temporal resolutions of the flux field in order to make proper use of the daily variations of CO₂ concentrations (Law et al., 2003). There is no consensus about which resolution should be used for a particular data frequency, so we choose the conservative option of high spatial resolution for regions with strong impact on our observing sites. We divide the western European and north-eastern Atlantic regions (Fig. 1, dark regions) to the pixel level of our transport model and solve for daily fluxes in each pixel. With the LMDZ zoom grid that we use (centred over Europe, Fig. 1), this generates of order 200 000 unknown flux components. Such a choice imposes a strong constraint on the numerical formulation of the inverse problem, especially for the inclusion of the prior spatial correlations between the “pixel fluxes” to regularise the problem (see Sect. 2.4). For the rest of the

globe, we keep the same regions and monthly time resolution as for the first inversion step. To compute the solution (flux values X and uncertainties Pª), we thus use a matrix formulation (Tarantola, 1987, Eqs. 4.5 and 4.6) that only needs the inversion of a matrix of the size of the observation vector (180),

$$X = X^b + P^b H^T (H P^b H^T + R^0)^{-1} (Y^0 - H X^b) \quad (1)$$

$$P^a = P^b - P^b H^T (H P^b H^T + R^0)^{-1} H P^b \quad (2)$$

with H the model response functions, Y⁰ the observations, Xᵇ the prior unknowns (fluxes and initial conditions), R⁰, Pᵇ and Pª the error covariance matrices on the observations, the prior unknowns and the posterior unknowns, respectively.

## 2.2 Atmospheric transport

We use the general circulation model of Laboratoire de Météorologie Dynamique, LMDZ (Hourdin and Armengaud, 1999), that can be zoomed over a particular region of the globe. For this study, the grid is stretched as illustrated in Fig. 1, with a mean horizontal resolution over Europe of 50×50 km and 19 sigma-pressure layers up to 3 hPa. This corresponds to a vertical resolution of about 150, 180, 350, and 460 m for the first four levels, respectively, and to a resolution of about 2 km a the tropopause. In order to stay as close as possible to the observed synoptic events, we relax the simulated winds and temperature toward the analysed fields of ECMWF with a time constant of 2.5 h (“nudging” mode). The advection of tracers is calculated with the finite-volume, second-order scheme, proposed by Van-Leer (1977) as described in Hourdin and Armengaud (1999). Deep convection is parameterized according to the scheme of Tiedke (1989) and the turbulent mixing in the planetary boundary layer is based on a local second-order closure formalism (Hourdin and Armengaud, 1999). The high horizontal resolution and the “nudging” capabilities of LMDZ are crucial to simulate the large spatio-temporal variability of the CO₂ concentrations observed at continental sites (Geels et al., 2005¹). Note, finally, that with a 50-km grid, high altitude stations (PRS, CMN, SCH) are located a few levels above ground in the model. This choice is important especially if we consider the diurnal cycle (Geels et al., 2005¹).

As mentioned in the introduction, we use the “retro-transport” approach implemented in LMDZ to compute the response functions at all sites (as in Cosme et al., 2005). Let us denote a measurement performed at a given station by,

$$M = \int_{\Omega \times \tau} \rho \mu c \, dx dt \quad (3)$$

with ρ being the density of the air, μ(x, t) the distribution of the measurement (uniform over a day and within the station mesh in our case), c the tracer concentration per unit mass of air, Ω the spatial domain, and τ the time domain. It was shown by Hourdin et al. (2005a) that M could be explicitly expressed as a function of the initial condition at t₀ (c|t₀), and

<footnote n="1">¹Geels, C., Bousquet, P., Ciais, P., Gloor, M., Peylin, P., Vermeulen, A. T., Dargaville, R., Brandt, J., Christensen, J. H., Frohn, L. M., Heimann, M., Karstens, U., Rödenbeck, C., and Rivier, L.: Comparing Atmospheric Transport Models for Regional Inversions over Europe. Part 1: Mapping the CO₂ Atmospheric Signals, Telus, submitted, 2005.</footnote>
of the tracer emission s, using the “retro-tracer” concentration c*:

$$ M = \int_{\Omega} \rho c^* c_{|t_0} dx + \int_{\Omega \times \tau} \rho c^* s dxdt \\ = \text{Initial condition term} + \text{Source term} \qquad (4) $$

This equation corresponds to Eq. (16) of Hourdin et al. (2005a) with a null inflow term as is the case with a global model. The “retro-tracer” c* is the solution of the adjoint of the tracer transport equation with respect to the air-mass weighted scalar product and its distribution can be simply computed by reversing the sign of the different advection terms and keeping the sign of the unresolved diffusion terms<fnref n="2" />. These properties rely on the fundamental time symmetry of fluid transport. Note, however, that in the numerical world the symmetry is not necessarily preserved with the spatial and temporal discretisation of the transport equations. The introduction of slope limiters in the Van-Leer (1977) advection scheme used in LMDZ is one example. The “retro-transport” in this case differs from the adjoint of the direct numerical model. The choice between the two approaches to represent the transport characteristics is not straightforward as it depends on the use that will be made of the results. In the case of LMDZ, Hourdin et al. (2005b) have shown a very good agreement between the forward and the backward calculations, in the context of the European Transport Experiment (see their Figs. 2 and 3).

The backward simulation directly provides the sensitivity of a given measurement to all past surface sources (up to one month in this study) and to the concentration of all 3-D grid boxes at $t_0$ (the initial condition). Using the same matrix notation as for the expression of the solution (Eq. 1), we can simply write:

$$ M = \mathbf{H}_0 \mathbf{C}_{|t_0} + \mathbf{H}_s \mathbf{S} \qquad (5) $$

with $\mathbf{H}_0$ and $\mathbf{H}_s$ being the sensitivity of the model concentration to the initial conditions $\mathbf{C}_{|t_0}$ and to the sources $\mathbf{S}$, respectively. Note that the vector of unknowns $\mathbf{X}$ defined in Eq. (1) gathers the two components, $[\mathbf{C}_{|t_0}, \mathbf{S}]$.

In practice, we perform a forward meteorological simulation nudged by ECMWF winds and store on a 3-h time step all large scale mass fluxes as well as convective detrainments and entrainments. We then run the transport model backward in time using the previously stored quantities and inject a pulse of “retro-tracer” for each day in November and each site. Each pulse is run back separately from the others. Because we did not consider any selection criteria for the observations (we averaged all hourly data across the full 24 h), we injected the “retro-tracer” uniformly in time each day ($\mu$ is uniform). Note that a daytime sampling selection criteria could be simply implemented with the injection of the “retro-tracer” during daytime only. The backward simulations for

all “retro-tracers” are performed back to the first of November at 00:00 h. We then average the simulated “retro-tracer” concentrations for each day at the surface pixels to get the sensitivity of the measurement to all daily sources. Again, we can use a weighting time-average to account for a time distribution of the CO₂ sources, like the diurnal cycle of the biospheric fluxes. In this study, we assume constant daily sources and thus neglect the diurnal variation in biospheric fluxes in November. Finally, we also take the “retro-tracer” 3-D fields at the first November – 00:00 h ($t_0$) to obtain the sensitivities to the initial conditions.

## 2.3 Treatment of initial conditions

We see from Eqs. (4) and (5) that there is a contribution to the observed concentrations by the initial concentration field as well as the flux field. Normally we choose a study period so long that we can afford a spinup period. The diffusive nature of atmospheric transport means that, provided this period is longer than about one year, the spatial structure of the initial condition will have dissipated to a more or less uniform background value. Hence we can specify all we need for the impact of the initial condition on the inversion by setting this value. Normally we include it as an extra unknown in the solution. With our short study period, we cannot afford such a spin-up so we must solve for the initial condition directly. At first sight this is daunting since the full three-dimensional field of LMDZ contains more than $5 \times 10^5$ grid boxes. However, almost all combinations of these values are unobservable by our chosen network, since it takes observations only over one month and a limited domain. For example, it is unlikely a perturbation in the southern polar stratosphere will be observed at all. It is computationally extravagant to attempt to solve for all these unobserved components.

The task then is to define a subspace of the possible initial conditions which is important for the inverse problem and solve only for elements of that subspace. Fortunately this subspace is easy to define. Let us rewrite the contribution of the initial concentration field to the model station concentration: $\mathbf{H}_0 \mathbf{C}_{|t_0}$ (Eq. 5). Using the usual Singular Value Decomposition, we can reform the Jacobian (or mapping) matrix $\mathbf{H}_0$,

$$ \mathbf{H}_0 = \mathbf{U} \cdot \mathbf{\Lambda} \cdot \mathbf{V}^T \qquad (6) $$

Where $\mathbf{U}$ and $\mathbf{V}$ are orthonormal and $\mathbf{\Lambda}$ is diagonal with rank of the number of observations (180). Numbers in parenthesis indicate the size of the different matrices after dropping the null parts for $\mathbf{\Lambda}$ (and the corresponding part of $\mathbf{V}$). $\mathbf{U}$ spans the space of observations and, importantly $\mathbf{V}$ spans the observable subspace of possible initial conditions. If we now project $\mathbf{C}_{|t_0}$ onto this subspace $\mathbf{V}$, we get the new unknowns $\mathbf{C}_{|t_0}^*$,

$$ \mathbf{C}_{|t_0}^* = \mathbf{V}^T \cdot \mathbf{C}_{|t_0} \qquad (7) $$

<footnote n="2">²See Hourdin et al. (2005b) for a complete description of the treatment of the different subgrid-scale processes.</footnote>
and a simpler expression for the Jacobian $\mathbf{H}_0^*$,

$$ \mathbf{H}_0^* = \mathbf{U} \cdot \Lambda \qquad (8) $$
$$ (180, 180) \quad (180, 180) \quad (180, 180) $$

Most important is the reduction in the dimension of the unknowns from more than $5 \times 10^5$ to 180. The new unknowns are the multipliers of $\mathbf{H}_0^*$, and we can simply rewrite the model station concentration as,

$$ M = \mathbf{H}_0^* C_{l0}^* + \mathbf{H}_s \mathbf{S} \qquad (9) $$

The overall vector of unknowns $X$ thus becomes $[C_{l0}^*, S]$.

It remains to find the initial uncertainty for $C_{l0}^*$. By the usual rule for error propagation, the variance-covariance matrix, $\text{Cov}(C_{l0}^*)$, on $C_{l0}^*$ is calculated by

$$ \text{Cov}(C_{l0}^*) = \mathbf{V}^T \text{Cov}(C_{l0}) \mathbf{V} \qquad (10) $$

In the case where the variance-covariance matrix on the initial unknowns is uniform and diagonal, $\text{Cov}(C_{l0})=\sigma^2\mathbf{I}$, (i.e. $\sigma$ being the prior uncertainty), we see that the orthonormality of $\mathbf{V}$ yields

$$ \text{Cov}(C_{l0}^*) = \sigma^2 \mathbf{I} \qquad (11) $$

where the identity matrix here has the dimension of the number of observations. Any other form of $\text{Cov}(C_{l0})$ will yield a much more complicated structure for $\text{Cov}(C_{l0}^*)$. We use the diagonal simplification throughout this paper. Finally, note that $\text{Cov}(C_{l0}^*)$ is one part of the covariance matrix of all prior unknowns, $\mathbf{P}^b$.

In the above simplification we neglect any correlations between the errors on the initial concentration field. However, the direct simulation used to define the initial field at 00:00 h of November 1st will contain positively and negatively correlated errors. Among the main reason for these errors is the coarse resolution of the fluxes from the first inversion step. As an example, consider a region of 500 km × 500 km containing all levels in the troposphere (10 grid cells vertically). Assume a uniform uncertainty $\sigma$ for the concentrations in this region and an uncertainty correlation $r$. Further assume that our observing set-up can see only the average of these pixels: one value of $C_{l0}^*$ noted $x_0$. The eigen-vector of $\mathbf{H}_0$ corresponding to this set-up will be uniform with values $1/\sqrt{N}$, where $N$ is the number of pixels in the region. Applying the rule for error propagation it can be shown (with a little algebra) that the prior variance of the coefficient $x_0$ in our transformed basis is

$$ \text{Cov}(x_0) = \frac{N\sigma^2 + (N^2 - N)r\sigma^2}{N} = [(N - 1)r + 1]\sigma^2 \qquad (12) $$

With $r=0.5$ and $N=1000$ we see that positive correlation increases the uncertainty on the coefficient in the reduced subspace by a factor 22.4. We assume uncorrelated prior uncertainties on the values of the initial condition so we must

inflate the uncertainty to account for this correlation, otherwise the uncertainty in our reduced subspace will be artificially low. On the other hand, negative correlation would decrease the uncertainty. Given these considerations and to keep a conservative choice of uncertainty, we have prescribed an initial uncertainty on each element in the initial condition of $10 \text{ ppm} \times 20=200 \text{ ppm}$. This appears large, but recall that the initial condition is a nuisance variable in this problem and that a conservative choice is preferable to a wrong answer.

## 2.4 Prior flux covariance

Increasing the number of independent regions allows us to recover more information from the atmospheric data and to be less sensitive to the a priori spatial structure of the fluxes (Kaminski et al., 2001). However, the inverse problem becomes highly under-constrained and hence the reduction of the estimated error on each pixel remains small. One can replace the hard constraint of fixed patterns within large regions by a soft constraint of correlated fluxes defined at the resolution of the transport model.

In this study, the parameter space (unknowns) comprises the fluxes for the European land model pixels (2853), the Northeast Atlantic ocean model pixels (4377), and 26 large regions (rest of the world) and the initial concentration field. Recall that the time resolution is monthly for the large regions and daily for the pixels. We define a prior variance/covariance matrix $\mathbf{P}^b$ on all those unknowns. For the large regions, we use the estimated variances from the first inversion step, but neglect covariances. The case of the initial conditions has been treated in the previous section. For model-pixels' fluxes, we spread a summed uncertainty of $1/\sqrt{12 \text{ GtC/month}} (\approx 0.29)$ over the land pixels and $0.5/\sqrt{12 \text{ GtC/month}} (\approx 0.14)$ over the ocean pixels. This is equivalent to an uncertainty of 1 GtC/yr and 0.5 GtC/yr, respectively, for the traditional approach of solving for independent monthly fluxes. We apply the following recipe to distribute this uncertainty in space (pixels) and time (days)

1. We assume no covariances and apportion the total variance uniformly in time and according to the area of the pixels.

2. We define spatial correlations within (but not between) land and ocean regions and convert them into covariances by scaling with the variances from step 1 in order to get the matrix $\mathbf{P}^b$. We use no temporal correlations.

3. The positive covariances defined in step 2 inflate the variance of the area-integrated flux, so we rescale $\mathbf{P}^b$ to return to the original a priori estimates (0.29 and 0.14 GtC/month).

We are now left with the problem of defining the error correlations. We should in principle account for the different processes and driving factors that control the carbon sources and sinks at larger scales than model resolution (i.e. land cover
type, climatic variables, fertilization etc). If the inversion requires a flux correction along the path of the “retro-plume”, it is likely that the correction should affect the neighbouring pixels, or at least the pixels with the same climate and the same cover type. Defining such correlation is difficult. Sensitivity analysis of models that produce these fluxes can yield covariance matrices suitable for this task (Kaminski et al., 2002, e.g.) but these are under construction for most of the key processes. As a first approach in this study, we use the distance as the main criterion and define exponentially decaying error correlations $cor_{i,j}$ between model pixels $i$ and $j$, according to: $cor_{i,j} = \exp(-d_{i,j}/L)$ with $L$ being a correlation length and $d_{i,j}$ the distance between the two pixels. Rödenbeck et al. (2003) used this approach at monthly resolution. They determined correlation lengths of 1275 km and 1912 km over land and ocean, respectively, from an analysis of autocorrelation and cross-correlation functions of monthly CO₂ fluxes calculated by a set of land and ocean biogeochemical models. On a daily time step, the value of $L$ is probably smaller, as daily flux patterns are usually less homogeneous than monthly flux patterns. In this study, we investigate two different cases with $L$ set to 500 km or 2000 km for the land pixels. For the ocean pixels, we doubled the correlation length to account for the lower spatial variability as in Rödenbeck et al. (2003).

The size of $P^b$ (217 106×217 106) makes it impossible to store in memory. However, given the form we use for the inverse solution (Eqs. 1, 2), we only need to perform matrix multiplications and to invert a small matrix that has the size of the number of observations and not the number of sources. Such an approach is thus extensible to many more stations and even to hourly observations. We can also project both sides of Eq. (2) onto a subspace of interest, e.g. groups of pixels (line by line projection) and hence avoid the need for calculations with the flux covariance as a whole. Finally note that in the particular case of no temporal correlation (this paper) the block form of $P^b$ simplifies the above multiplications.

## 3 Results

We performed a series of inversions to investigate the sensitivity of the estimated fluxes and initial conditions to the critical parameters described above. The purpose of this paper is to validate the inverse method for a particular month using some statistical diagnostics and to identify the critical aspects of the inversion, like the sensitivity to the correlation lengths. In the following, we nevertheless mainly refer to a standard case with correlations of 500 km and 1000 km for land and ocean pixels, respectively.

### 3.1 Model-data fit

Figure 2 displays the a priori and the estimated model concentrations together with the daily data and their uncertainties for the six sites in November. As a general feature, optimised model concentrations capture most of the daily variations at all sites, except for day 24 at SCH and day 26 at CMN. The reduced chi-square diagnostic ($\chi^2$; i.e., twice the cost function at its minimum divided by the number of observations) is much lower than one, 0.09, which attest for the quality of the fit. Given the data uncertainties that we choose and the prior flux uncertainties of the standard configuration, there is enough freedom to fit all daily observations. Note that the model-data fit is significantly better at the beginning of the month than at the end (the daily reduced $\chi^2$ values increase from 0.01 to 0.2 during the month). This results from the optimisation of the initial conditions, since there are more degrees of freedom available to fit data near the start of the period. We also note that the fit degrades in the case with higher correlation lengths (2000 km and 4000 km for land and ocean, respectively) with a reduced $\chi^2$ of 0.6. At the other extreme, an inversion with no correlation produces a perfect model-data fit. Finally, one should mention that the optimized fluxes do not degrade the fit to the station used in the first inversion so that our two-step procedure has not introduced an inconsistency.

More important is the good agreement between the prior modelled concentrations (dashed lines) and the data, especially if we consider the phase of the different synoptic events. At high-altitude stations like CMN, PRS, and even SCH, the events are already well captured in the prior. This reflects the ability of the model to satisfactorily reproduce the major features of the transport over Europe during November 1998. At MHD and WES coastal sites, the phase of the synoptic events, that usually last for a few days, is also well captured in the prior fit but the model amplitude is too small. For these stations, the inversion will correct the nearby fluxes assuming that the transport and more importantly the vertical mixing close to the surface is correctly modelled. At the Hegyhatsal tower (HUN), the prior model misses the large synoptic event at the end of the month. Such a deficiency is also present at WES. Great care should thus be taken in the interpretation of the inverse estimates around those sites at the end of the month. Overall, we stress the importance of the prior model-data fit. Strong a priori inconsistency should always be screened before the inversion and checked for potential transport errors. In the case of LMDZ, we are confident that the transport model is able to reproduce most of the daily concentration variations, at least for the phase of the events, a prerequisite in this approach.

### 3.2 Contributions to the model concentration output

We now analyse the contribution from the different components, i.e. the sum of European land pixels, the sum of North
<figure><img src="image_3.png" /><figcaption>**Fig. 2.** Model prior (dashed lines) and posterior (solid lines) concentration together with the data and their uncertainties, for the six continuous sites used in this study. Data with no error bar were not used in the inversion (missing values).</figcaption></figure>

Atlantic pixels, the sum of all “big regions”, and the initial conditions to the total modelled concentration. Figure 3 displays the prior and posterior contributions only at the Hegyhatsal tower HUN and the coastal site MHD. As expected, it appears that

1. the contribution from the European land pixel dominates,

2. the contribution from the distant regions increases over time,

3. the contribution of the initial conditions decreases through time toward a non-zero mean value.

Note that by “contribution” we implicitly refer to the day to day variations induced by each component. We do not consider here the mean value as it is directly adjusted with an overall offset in the inverse procedure.

In November, the European land fluxes mostly control the synoptic variations of the concentration at the six European sites, the contribution from the North Atlantic being much smaller. At the high-altitude stations (CMN, PRS, SCH) or

the western MHD site, the contribution from the rest of the northern hemisphere (predominantly North America) is significant and even dominates at the end of the month. At the beginning of the month the air masses from North America have still not reached the European continent and their contribution increases through time, whereas the influence of the initial conditions is rather large in early November. The large a priori model-data misfits at MHD between the 13 and the 20 of November and at HUN between the 24 and the 27 of November (too-low model concentrations, Fig. 2) are both adjusted through an increase of the European land pixel fluxes.

The case of the initial conditions (IC) is very specific to this study. As a first outcome, the IC appear to be rather important and to produce day to day concentration variations up to 2 ppm (Fig. 3). Although their contribution significantly decreases with time, they still induce variations larger than 0.5 ppm after 20 days at most sites. Such variations are comparable for the coastal or high-altitude sites to the variations induced by the other components. It is thus quite important to account for an initial concentration field when doing an
<figure><img src="image_4.png" /><figcaption>**Fig. 3.** Contribution of the different components, i.e. the sum of the European pixel based land fluxes ("Pixel-land"), the sum of the North Atlantic pixel based fluxes ("Pixel-ocean"), the sum of the extra-European large regions ("big regions"), and the initial conditions, to the model concentration at 2 stations, Mace Head (MHD) and the Hegyhatsal tower (HUN). Plain lines are for the posterior values and dotted lines for the prior.</figcaption></figure>

inversion over a short time period (i.e. a month). This result is important when considering analysing continuous data acquired during short campaigns. The optimisation of the initial conditions brings also another aspect. After inversion, we obtain some concentration changes from the prior values on the order of 0.5 ppm for the first 20 days and much less afterward. Although not critical, these changes highlight a potential source of uncertainty in an inversion where the IC would not be optimised. In the case of an inversion with a much smaller prior uncertainty on the initial field (20 ppm instead of 200 ppm, see Sect. 2.3), the concentration changes from the prior become negligible. It is thus important to define the initial concentration field as accurately as possible, especially if we can not solve for it. However, one should

note that these later conclusions probably depend on the size and geometry of the network and we might anticipate that with more stations (especially "upwind" of Europe) we might produce a better initial condition and hence require smaller changes from the prior.

### 3.3 Retrieved fluxes and uncertainties

This paper, with its short study period and low data density, is intended as an exploration of the methodology for using continuous concentration observations. Thus we consider here the information added by the observations to our prior estimate rather than the final values of the fluxes themselves. A subsequent paper, with more stations and a longer study
period, will examine the European carbon cycle in the light of this new information. The information added by the observations is in two forms, an increment to the prior flux distribution and a reduction in its uncertainty.

Figure 4, upper panel, displays the flux changes (posterior minus prior) for two correlation lengths, 500 km and 2000 km over land pixels. The largest changes occur around the Baltic area, close to the WES station. Figure 2 shows a large underestimation of the concentration at WES from the prior flux distribution. Hence, the inversion increases the fluxes south and east of that station (the upwind direction during November) by more than 70 gC/m2/month. Whether these very large changes (more than 100% of the a priori values) are realistic or not is a difficult question to answer. There are significant anthropogenic emissions in this area, so part of the change may be a refinement of these estimates. It is also possible that transport errors, especially ventilation of the boundary layer, may result in incorrect amplitudes for synoptic events. The impact of such transport errors, which affect both the prior estimate of concentration and the response to any increment in flux, is hard to determine. In southern Europe around the Mediterranean, the fluxes are decreased from the prior by about 20 gC/m2/month, a value close to 80% of the initial field and also questionable. Finally, one should note that the flux changes over the ocean are much lower than over land, hardly appearing on the scale of Fig. 4.

If we now consider the difference between the two correlation lengths, we notice a large change in the spatial extent of the flux increment. With the 2000 km correlation length, the flux increment impacts the whole of Europe more uniformly than with the 500 km length scale. This parameter is critical in determining the spread of the information contained in the continuous measurements. However, given the lack of theoretical foundation for the correlation length, it is difficult to choose a particular value. The 2000 km length was chosen as an extreme case to emphasize the role of the flux error covariance for an inversion with only a few sites. The difference between these cases suggests that the use of large regions in previous inversion studies might have overestimated the power of isolated concentration measurements to constrain regional fluxes.

The lower panel of Fig. 4 shows the percentage reduction of uncertainty in estimated flux (monthly average). This quantity indirectly measures the information added by the network regarding fluxes at each pixel. Note that the mean values directly depend on the prior flux errors so that we will mainly discuss the differences between pixels. The spatial patterns of uncertainty reduction are rather independent of the pattern of the flux increment. The increment depends both on the connection of a given pixel with an observing site as well as the mismatch of prior simulated concentrations. The mismatch does not affect the uncertainty reduction. In the 500 km case, the error reduction is maximum around each station (except for the high-altitude sites, PRS and CMN, not

much influenced by regional pixels in November). The pattern also reflects the major patterns of flow. For example, we see a north-south component centred in Germany and Eastern France, and a West-East component around MHD that is visible on the ocean pixels. Over the central part of Europe, the mean error reduction is of order 15–20%, with maximum values at the stations up to 30% (values partly masked on Fig. 4 by the station symbols). When we increase the correlation length to 2000 km, the mean error reduction significantly increases both in magnitude and extent. As with the flux changes, this behaviour reflects the propagation of the information further from the main air flow paths.

Figure 5 displays the time evolution of the error reduction for some geographic regions, the Western, the Mediterranean, the Balkan, and the North plus Central European regions (see their contour in Fig. 1). We clearly see that our network of six sites strongly constrains Western Europe during November but also the Balkan area at the end of the month. The large day to day variations (up to 20%) reflect changes both in the track of the retro-plumes and in the amplitude of the convective mixing. The peak of reduction (49%) for the Balkan region around the 23 of November is somewhat anti-correlated with the error reduction over Western Europe, reflecting a period with easterly winds at HUN, CMN, PRS and SCH stations. Other anti-correlations, e.g. between the North Atlantic and the North+Central Europe are also apparent. With the increased correlation length (2000 km), the results show similar patterns (not shown). All curves are shifted towards larger error reductions (10–15% value) compared to the 500 km correlation length and the differences between the land regions decrease. If we increase the correlation length up to 10 000 km, we tend towards the case of a “big region” approach (“European pixels” being one region) and the error reduction is then identical for each sub-region, around 60%. This test suggests that the uncertainties might be underestimated using large regions.

Overall the error reduction diagnostic is crucial; it allows us not only to identify the regions that are relatively well constrained, but also to discuss the pertinence of the estimated fluxes. Before interpreting a large flux correction, one should always verify that it occurred over a well constrained area. Finally such a diagnostic helps to define the inverse set-up, defining quantities like the correlation lengths. In the 2000 km case, the spatial pattern of the error reduction does not reflect structures that can be related to the major air flow path around the different sites, like in the 500 km case. The pattern now reflects predominantly the distance from the stations, a feature that is probably oversimplified.

## 4 Summary and conclusions

This study makes use for the first time of continuous CO₂ measurements over Europe for one month (six sites in November 1998) in a full 3-D inverse approach specifically
<figure><img src="image_5.png" /><figcaption>**Fig. 4.** Maps of the posterior minus prior monthly mean fluxes and of the monthly mean errors reduction (mean across all daily error reductions, defined as posterior minus prior divided by the prior) expressed in percentage for two cases of correlation lengths, 500 km and 2000 km for land pixels (note that for ocean pixels the length is doubled). The stars represent the location of the six stations.</figcaption></figure>

designed to infer daily surface fluxes. The major features of this new inverse approach are

- to use both continuous data gathered on a daily time step and monthly data,

- to optimise not only the surface fluxes at the model spatial resolution over Europe but also the initial concentration field,

- to use the “retro-plume” approach (Hourdin et al., 2005a) within the LMDZ transport model to efficiently compute the response functions.

This study is primarily methodological, although it makes use of real data rather than the identical twin experiments of previous studies. The major outcomes can be summarized as

- the LMDZ model is able to satisfactorily represent the
<figure><img src="image_6.png" /><figcaption>**Fig. 5.** Time evolution of the error reduction for different regions, Western Europe, Mediterannean Europe, Balkan Europe, North+Central Europe, and North Atlantic (see their boundaries in Fig. 1). Results correspond to the “500 km case” of the land error correlation.</figcaption></figure>

phase of the major synoptic events in the concentration record for November over Europe. This is a necessary condition for using these synoptic events as information for inversions.

- The major adjustment to the flux and initial concentration field within the inversion is made to the pixels on the back trajectory from each station.

- The choice of correlation length or, more generally, prior covariance structure, strongly affects the patterns of both uncertainty reduction and flux increment.

- For campaign-style inversion studies, shorter than about 20 days, it is necessary to use an accurate field of the initial concentration and/or to include an optimization of that field in the inversion. We provide an algorithm for doing this feasibly.

Although preliminary, this work outlines the potential of continuous measurements in an inverse approach. There are, however, many simplifications in the current work. We chose to use complete diurnal averages, even though transport models have great difficulty simulating nocturnal boundary layers. We have made no special effort to deal with diurnal recertification. We have not attempted a comparison of our inversion results with independent bottom-up calculations, either from models or flux measurements. These will be addressed in developments of this work, once a sufficiently dense and extended dataset is available. A study for the year 2001 is in preparation.

*Acknowledgements.* The Commissariat à l'Energie Atomique partly funded this work, especially the visit of P. Rayner at LSCE (Paris) and it also contributed to the computing resources.

Edited by: M. Heimann

## References

Andres, R. J., Marland, G., Fung, I., and Matthews, E.: A 1°×1° distribution of carbon dioxide emissions from fossil fuel consumption and cement manufacture, 1950–1990, Global Biogeochem. Cycles, 10, 419–429, 1996.

Baker, D. F.: Sources and sinks of atmospheric CO₂ estimated from batch least-squares inversions of CO₂ concentration measurements, Ph.D. thesis, Princeton University, 2001.

Biraud, S., Ciais, P., Ramonet, M., Simmonds, P., Kazan, V., Monfray, P., O'Doherty, S., Spain, T. G., and Jennings, S. G.: European greenhouse gas emissions estimated from continuous atmospheric measurements and radon 222 at Mace Head, Ireland, J. Geophys. Res., 105, 1351–1366, 2000.

Biraud, S., Ciais, P., Ramonet, M., Simmonds, P., Kazan, V., Monfray, P., O'Doherty, S., Spain, G., and Jennings, S. G.: Quantification of carbon dioxide, methane, nitrous oxide and chloroform emissions over Ireland from atmospheric observations at Mace Head, Tellus, 54B, 41–60, 2002.

Bousquet, P., Peylin, P., Ciais, P., Quéré, C. L., Friedlingstein, P., and Tans, P. P.: Regional changes in carbon dioxide fluxes of land and oceans since 1980, Science, 290, 1342–1346, 2000.

Chevillard, A., Karstens, U., Ciais, P., Lafont, S., and Heimann, M.: Simulation of atmospheric CO₂ over Europe and western Siberia using the regional scale model REMO, Tellus, 54B, 872–894, 2002.

Cosme, E., Hourdin, F., Genthon, C., and Martinerie, P.: The origin of dimethylsulfide (DMS), non-sea-salt sulfate, and methane-sulfonic acid (MSA) in Eastern Antartica, J. Geophys. Res., 110(D3), D03302, doi:10.1029/2004JD004881, 2005.

Denning, A. S., Fung, I. Y., and Randall, D. A.: Gradient of Atmospheric CO₂ due to Seasonal Exchange with Land Biota, Nature, 376, 240–243, 1995.

Denning, A. S., Collatz, G. J., Zhang, C., Randall, D. A., Berry, J. A., Sellers, P. J., Colello, G. D., and Dazlich, D. A.: Simulations of terrestrial carbon metabolism and atmospheric CO₂ in a general circulation model, 1, Surface carbon fluxes, Tellus, 48B, 521–542, 1996.

Enting, I. G.: Inverse Problems in Atmospheric Constituent Transport, Cambridge University Press, 2002.

Enting, I. G., Trudinger, C. M., and Francey, R. J.: A synthesis inversion of the concentration and $\delta^{13}$C of atmospheric CO₂, Tellus, 47B, 35–52, 1995.

Fan, S., Gloor, M., Mahlman, J., Pacala, S., Sarmiento, J., Takahashi, T., and Tans, P.: A Large Terrestrial Carbon Sink in North America Implied by Atmospheric and Oceanic CO₂ Data and Models, Science, 282, 442–446, 1998.

Fan, S., Sarmiento, J. L., Gloor, M., and Pacala, S. W.: On the use of regularization techniques in the inverse modelling of atmospheric carbon dioxide, J. Geophys. Res., 104, 21 503–21 512, 1999.

Gerbig, C., Lin, J. C., Wofsy, S. C., Daube, B. C., Andrews, A. E., Stephens, B. B., Bakwin, P. S., and Grainger, C. A.: Toward constraining regional-scale fluxes of CO₂ with atmospheric ob-
servations over a continent: 2. Analysis of COBRA data us-
ing a receptor-oriented framework, J. Geophys. Res., 108, 4757-
ACH6, 2003.

GLOBALVIEW-CO₂: Cooperative Atmospheric Data Integration Project - Carbon Dioxide, CD-ROM, NOAA CMDL, Boulder, Colorado, also available on Internet via anonymous FTP to ftp: //ftp.cmdl.noaa.gov, Path: ccg/co2/GLOBALVIEW, 2002.

Gurney, K. R., Law, R. M., Denning, A. S., Rayner, P. J., Baker, D., Bousquet, P., Bruhwiler, L., Chen, Y.-H., Ciais, P., Fan, S., Fung, I. Y., Gloor, M., Heimann, M., Higuchi, K., John, J., Maki, T., Maksyutov, S., Masarie, K., Peylin, P., Prather, M., Pak, B. C., Randerson, J., Sarmiento, J., Taguchi, S., Takahashi, T., and Yuen, C.-W.: Towards robust regional estimates of CO₂ sources and sinks using atmospheric transport models, Nature, 415, 626–630, 2002.

Hourdin, F. and Armengaud, A.: Test of a hierarchy of finite-
volume schemes for transport of trace species in an atmospheric
general circulation model, Mon. Wea. Rev., 127, 822–837, 1999.

Hourdin, F. and Issartel, J.: Sub-surface nuclear tests monitor-
ing through the CTBT xenon network, Geophys. Res. Lett., 27,
2245–2248, 2000.

Hourdin, F., Talagrand, O., and Idelkadi, A.: Eulerian Backtracking of Atmospheric Tracers: I Adjoint derivation, parametrization of subgrid-scale transport, Q. J. R. Meteorological Soc., in press, 2005a.

Hourdin, F., Talagrand, O., and Idelkadi, A.: Eulerian Backtracking of Atmospheric Tracers: II numerical aspects, Q. J. R. Meteorological Soc., in press, 2005b.

Kaminski, T., Rayner, P. J., Heimann, M., and Enting, I. G.: On Aggregation Errors in Atmospheric Transport Inversions, J. Geophys. Res., 106, 4703–4715, 2001.

Kaminski, T., Knorr, W., Rayner, P., and Heimann, M.: Assimilat-
ing Atmospheric data into a Terrestrial Biosphere Model: A case
study of the seasonal cycle, Glob. Biogeochem. Cyc., 16, 1066,
doi:10.1029/2001GB001463, 2002.

Keeling, C. D., Piper, S. C., and Heimann, M.: A Three-
dimensional Model of Atmospheric CO₂ transport based on ob-
served winds: 4. Mean annual gradients and interannual varia-
tions, in: Aspects of Climate Variability in the Pacific and the
Western Americas, Geophysical Monograph 55, edited by: Pe-
terson, D. H., Washington, (USA), AGU, pp. 305–363, 1989.

Law, R. M., Rayner, P. J., Steele, L. P., and Enting, I. G.: Using high temporal frequency data for CO₂ inversions, Global Biogeochem. Cycles, 16, 1053, doi:10.1029/2001GB001593, 2002.

Law, R. M., Rayner, P. J., Steele, L. P., and Enting, I. G.: Data and modelling requirements for CO₂ inversions using high frequency data, Tellus, 55B, 512–521, doi:10.1034/j.1600-0560.2003.0029.x, 2003.

Law, R. M., Rayner, P. J., and Wang, Y. P.: Inversion of diurnally-varying synthetic CO₂: network optimisation for an Australian test case, Global Biogeochem. Cycles, 18, GB1044, doi:10.1029/2003GB002136, 2004.

Marland, G., Boden, T. A., and Andres, R. J.: Global, Regional, and National CO₂ Emissions, in: Trends: A Compendium of Data on Global Change, Carbon Dioxide Information Analysis Center, Oak Ridge National Laboratory, U.S. Department of Energy, Oak Ridge, USA, 2001.

Peylin, P., Bousquet, P., Ciais, P., and Monfray, P.: Differences of CO₂ Flux Estimates Based on a Time-Independent Versus a Time-Dependent Inversion Method, in: Inverse methods in global biogeochemical cycles, edited by: Kashibata, P., Heimann, M., Rayner, P., et al., American Geophysical Union, vol. 114, pp. 295–309, 1999.

Potter, C., Klooster, S., Myneni, R., Genovese, V., Tan, P.-N., and Kumar, V.: Continental-scale comparisons of terrestrial carbon sinks estimated from satellite data and ecosystem modeling 1982–1998, Glob. Planetary Change, 39, 201–213, 2003.

Rayner, P. J., Enting, I. G., Francey, R. J., and Langenfelds, R.: Re-
constructing The Recent Carbon Cycle From Atmospheric CO₂,
δ¹³C And O₂/N₂ Observations, Tellus, 51B, 213–232, 1999.

Rödenbeck, C., Houweling, S., Gloor, M., and Heimann, M.: Time-
dependent atmospheric CO₂ inversions based on interannually
varying tracer transport, Tellus B, 55, 488–497, 2003.

Ruimy, A., Dedieu, G., and Saugier, B.: TURC: A diagnostic model of continental gross primary productivity and net primary productivity, Global Biogeochem. Cycles, 10, 269–285, 1996.

Takahashi, T., Sutherland, S. C., Sweeney, C., Poisson, A., Metzl, N., Tilbrook, B., Bates, N., Wanninkhof, R. H., Feely, R. A., Sabine, C., Olafsson, J., and Noriji, Y.: Global sea-air CO₂ flux based on climatological surface ocean pCO₂, and seasonal biological and temperature effects, Deep Sea Res. II, 49, 2002.

Tarantola, A.: Inverse problem theory, Elsevier, Amsterdam, The Netherlands, 1987.

Tiedke, M.: A comprehensive mass flux scheme for cumulus parametrization in large scale models, Monthly Weather Review, 117, 1779–1800, 1989.

Trampert, J. and Snieder, R.: Model Estimations Biased by Trun-
cated Expansions: Possible Artifacts in Seismic Tomography,
Science, 271, 1257–1260, 1996.

Van-Leer, B.: Towards the ultimate conservative difference scheme: IV a new approach to numerical convection, J. Computational Physics, 23, 276–299, 1977.

Vukićević, T. and Hess, P.: Analysis of tropospheric transport in the Pacific Basin using the adjoint technique, J. Geophys. Res., 105D, 7213–7230, 2000.