What can we learn from European continuous atmospheric CO2 measurements to quantify regional fluxes – Part 1: Potential of the 2001 network

An inverse model using atmospheric CO 2 observations from a European network of stations to reconstruct daily CO2 fluxes and their uncertainties over Europe at 40 km resolution has been developed within a Bayesian framework. In this first part, a pseudo-data experiment is performed to assess the potential of continuous measurements over Europe using a network of 10 stations of the AEROCARB project such as in 2001 ( http://www.aerocarb.cnrs-gif.fr/ ). Under the assumptions of a small observation noise and a perfect atmospheric transport model, the reconstruction of daily CO 2 fluxes and in particular of their synoptic variability is best over Western Europe where the network is the densest. At least a 10 days temporal and a 1000 km spatial averaging of the inverted daily/40 km fluxes is required in order to obtain a good agreement between the estimated and the “true” fluxes in terms of correlation and variability. The performance of the inversion system rapidly degrades when fluxes are sought for a smaller temporal or spatial averaging.


Introduction
The problem of determining the space-time structure of surface CO 2 fluxes has gained considerable prominence along with the rising interest in anthropogenic climate change.The two classes of methods developed by scientists are distin-Correspondence to: C. Carouge (claire.carouge@lsce.ipsl.fr)guished as "bottom-up" and "top-down" methods.In the bottom-up approach, process models based on observations at very small scales are scaled up to regions of interest.In the top-down, or inverse, method, the integrated atmospheric signatures of the fluxes contained in atmospheric concentration gradients are disentangled to recover the spatial and temporal distribution of fluxes.Two advantages of the top-down method are that it integrates some of the small-scale heterogeneity that may not be of direct interest (either for policy or scientific applications) and that it does not require a direct knowledge of the processes giving rise to the fluxes.
The top-down method has hitherto been limited by a severe lack of data and biases in atmospheric transport models (Gurney et al., 2002;Stephens et al., 2007;Geels et al., 2007).Law et al. (2002) have proposed to use the records increasingly available from in-situ air sampling instruments.The use of continuous data creates stringent demands on atmospheric transport models (Geels et al., 2007;Gerbig et al., 2003), probably requiring higher spatial resolution, ability to reproduce diurnal PBL dynamics and synoptic shifts in transport.Law et al. (2002) also noted that the so-called aggregation error (Kaminski et al., 2001), due to the a priori spatial aggregation of surface fluxes to be optimized by inversions, was likely to be more serious with the use of such continuous data.To tackle this issue, Kaminski et al. (2001) suggested solving fluxes at model resolution in inversions, prescribing prior error correlations for fluxes in order to limit biases caused by lack of spatial resolution in the data, and the generation of non-physical solutions.
The combined requirements of a relatively high resolution for modeling atmospheric transport and an equally high resolution for determining sources in the inversion pose a significant computational challenge.In order to calculate the Jacobian matrix for the inversion problem, it is necessary to know the sensitivity of each concentration measurement to fluxes at all preceding times and places.In the conventional approach where a forward model run is attached to each source, this requires potentially millions of forward atmospheric transport simulations.An alternative is provided by the retrotransport of Hourdin et al. (2006a, b) who used the reversibility of transport for passive constituents to calculate the sensitivity of one measurement to all influencing sources, by emitting a pulse of passive tracer at the observing point and tracking its dispersion backwards in time.Thus only one tracer run is required per observation to map the sensitivity to all the sources.Another important although unrelated property of the model they used (LMDZ, originally described in Sadourny and Laval 1984) is its capability of a zoomed grid over a particular region.In this study, using LMDZ zoomed over Europe gives us the technical possibility of simulating the high-resolution transport necessary for matching the continuous European CO 2 stations, while retaining a globally coherent picture of CO 2 elsewhere.Peylin et al. (2005) used output from the same LMDZ model zoomed over Europe in their methodological study.That paper solved for daily fluxes over a region including Europe and parts of the North Atlantic during the test period of November 1998.Our work builds upon the results of Peylin et al. (2005) and goes beyond them in three important directions: 1) we use continuous daily data for an entire year which enables us to analyze the impact of very different meteorological and flux conditions, 2) we use data from ten rather than six continuous stations, and 3) we focus on so-called pseudo-data generated with the correct "true" value of the flux fields, and try to retrieve the fluxes under diverse, more or less optimistic, assumptions about the data.Pseudo-data experiment is a common tool for exploring the information content of potential observations (Gloor et al., 2000;Rayner et al., 2001;Law et al., 2002Law et al., , 2003)).
The main focus of this paper is the information content available from the existing network of continuous CO 2 measurements over Europe of the AEROCARB project in 2001 (http://www.aerocarb.cnrs-gif.fr/).In particular we wish to quantify the "optimal" spatial and temporal scales at which fluxes may be determined reliably.A companion paper (Carouge et al., 2010, CA08) considers extensions to the network and sensitivity tests of the inversion system to different error scenarios.The outline of the paper is as follows: in Sect.2, the inverse methodology is described with special reference to new developments from Peylin et al. (2005).Then, Sect. 3 presents and discusses the fluxes optimised in the pseudo-data experiment, where one year of pseudo-data with noise is inverted.

Grid-based regional inversion
For this study we describe below the inverse setup that is used to assimilate daily atmospheric CO 2 pseudo-data over Europe and retrieve daily fluxes at the model resolution (40 km).The use of a pseudo-data modeling framework allows investigation, in a "controlled environment", of the potential of the 2001 European network to infer regional CO 2 fluxes.In this section, we detail in turn the overall inverse approach (Sect.2.1), the transport model (Sect.2.2) and the Jacobian matrix calculation (Sect.2.3), the pseudo-data generation (Sect.2.4), the prior error covariance matrix (Sect.2.5), and a few critical technical choices (Sect.2.6).

Overall approach: the use of pseudo data
Following the methodology of Peylin et al. (2005), daily CO 2 fluxes are optimized over Europe at the LMDZ model spatial resolution of ∼40 km, using information contained in daily pseudo-observations.A Bayesian synthesis inverse method is applied, based on a matrix formulation, to invert one year of fluxes.
Pseudo-data allow testing the accuracy of the solution and the impact of different inverse setups by comparing the inverted fluxes to the "true" fluxes.The extent to which the results of those experiments can inform on a real-data case strongly depends on the realism of the inverse setup.In this study, we consider an ideal case compared to Peylin et al., 2005.We only optimize for daily land ecosystem fluxes over Europe and for air-sea fluxes over the eastern North Atlantic where small flux adjustments are allowed.We consider fluxes from other regions are null in building the pseudo-data and performing the inversion.Note that we chose to not consider fossil fuel emissions in our setup.Thus, we do not need to perform a first global inversion with all stations and all fluxes, as Peylin et al., 2005 did using real observations.Rather, we only have to define two sets of fluxes over Europe: the target or true fluxes and the first guess prior fluxes to be optimized.A rigorous approach is to perturb the true flux distribution according to a given error covariance matrix, P b , in order to define the prior fluxes (Chevallier et al., 2007).In that case, the inverse problem is statistically consistent (P b being used in the optimization process).However, considering the poor knowledge of land-ecosystem fluxes and their error covariance (Chevallier et al., 2006), we want to start with differences between prior and the true fluxes that are a plausible representation of the differences between any prior and the unknown true flux in the real world.This is the reason why we chose one ecosystem model (TURC) to generate the daily prior flux maps, and another independently developed ecosystem model (ORCHIDEE) to produce the true flux distributions (see Sect. 2.4).
The pseudo-data are calculated by applying the LMDZ transport model to the daily true fluxes.Note that the year 2001 was consistently chosen for atmospheric transport and for forcing the land ecosystem model calculating the true fluxes.The pseudo-data are perturbed in order to account for data and model errors.We choose a white noise of relatively small amplitude (±0.3 ppm).Such a small noise is representative of data errors only.The atmospheric transport is considered perfect as illustrative of an ideal case (CA08 explore a case with a large white noise).
The inversion of daily flux maps during one year is divided into a series of consecutive optimizations overlapping in time in order to cope with the numerical size of the problem (see Sect. 2.6).The quality of the results will be analyzed by comparing the optimized fluxes to the true fluxes using statistical diagnostics.

Global atmospheric transport model zoomed over Europe
We use the LMDZ transport model (Sadourny and Laval, 1984) and the same grid as Peylin et al. (2005) with a zoom centered over Europe leading to a maximum resolution of 40×40 km and 19 sigma-pressure layers up to 3 hPa (10 layers in the troposphere).In LMDZ, the advection of tracers is calculated based on the finite-volume, second-order scheme proposed by Van Leer (1977) as described by Hourdin and Armengaud (1999).Deep convection is parameterized according to the scheme of Tiedtke (1989) and the turbulent mixing in the planetary boundary layer is based on a local second-order closure formalism (Hourdin and Armengaud, 1999).Finally, model winds are relaxed towards analyzed fields of the European Center for Medium Range Weather Forecasting (ECMWF) for the year 2001 in order to remain as close as possible to the observed synoptic events (with a time constant of 2.5 h).The model has been widely used for climate studies (IPCC, 2007) and for direct and inverse modeling of CO 2 (Peylin et al., 2005) and of other atmospheric trace gases (Hauglustaine et al., 2004;Bousquet et al., 2005Bousquet et al., , 2006)).

Jacobian transport matrix calculation
Within a synthesis inverse approach, one needs to define the sensitivity of concentrations (at each site and each moment in time) to the surface flux of each source region (each model grid cell during one time step).The sensitivity of one measurement to all the sources is often called the influence function, or concentration footprint.We used the retro-transport formulation implemented in LMDZ to compute these sensitivities at ten atmospheric stations with daily observations (as in Peylin et al., 2005).The approach relies on the time symmetry of modeled transport which defines a "retro-tracer" (transport backward in time) equivalent to adjoint transport of a tracer without developing the adjoint model (Hourdin et al., 2006a;Issartel and Baverel, 2003).Time symmetry occurs because a tracer and a retro-tracer would symmetrically spread out along with time in LMDZ.If some passive tracer is emitted at one point, after enough time, it would be spread all over the world.A retro-tracer represents the locations of the air parcels composing a measurement at times before the measurement time.These locations would spread across the world if the simulation time goes back long enough, and thus should do a retro-tracer.So a tracer and a retro-tracer are symmetrically transported.
The distribution of the retro-tracer is computed by a single simulation backward in time, i.e. reversing the sign of the different advection and convection mass fluxes but keeping the sign of the unresolved diffusion terms.On account of discretization issues, the retro-tracer is not the exact solution of the adjoint of the transport equation but, in the case of LMDZ, Hourdin et al. (2006b) have shown a fair agreement between the forward and backward calculations in the analysis of the European Transport Experiment (ETEX).We also realized a comparison between forward and backward calculations for three months at the European stations used in this study.It shows daily differences smaller than 0.3 ppm at all stations.With this approach, computation of sensitivity to all fluxes requires only one backward simulation per observation, with a pulse of retro-tracer emitted backward at each time step for each station.

Pseudo-data and their error
A network comprising the 10 continuous surface stations that were operating in Europe in 2001 is used here in the context of the AEROCARB project (Table 1, Fig. 1, and http://www.aerocarb.cnrs-gif.fr/).Daily pseudo-data at these ten stations are generated with LMDZ for the whole year with daily Net Ecosystem Exchange (NEE) over Europe of the ORCHIDEE model (Krinner et al., 2005).For consistency, ORCHIDEE was forced by meteorological fields from ECMWF for the year 2001 (Uppala et al., 2005).OR-CHIDEE is a state of the art mechanistic model that computes the turbulent fluxes of CO 2 , H 2 O and energy on a half hourly basis, and the dynamics of ecosystem C and water pools (phenology, allocation, growth, mortality, soil   organic matter decomposition) on a daily basis.The modeled daily NEE compares relatively well with the measured NEE at specific eddy-covariance flux towers, with a mean standard deviation of the model-observed differences close to 2 gC/m 2 /day (Chevallier et al., 2006).Atmospheric transport models are known to have difficulties in simulating the dayto-day variability of the nocturnal planetary boundary layer (PBL) height (Geels et al., 2007).We thus selected the model concentrations for daytime only (11:00 to 16:00 local time) to generate the pseudo-data and the associated influence functions.
The choice of assigning a white noise of standard deviation 0.3 ppm for the pseudo-data at each station reflects an optimistic setup, with small measurement uncertainties and no transport error (i.e. the model is able to represent perfectly well the measurements at each site).The value of 0.3 ppm is consistent with instrumental noise in measuring CO 2 at continuous stations.In a companion paper (CA08), we investigate the impact of a larger and likely more realistic noise on the results.This Gaussian noise defines the error statistics of the pseudo-data in the inversion, no error correlations being assumed between different stations (i.e. the error covariance matrix, R 0 , is diagonal).

European land ecosystem fluxes and errors
The prior daily NEE flux maps are calculated from the TURC model (Lafont et al., 2002).TURC is a diagnostic model driven by 10-daily satellite vegetation index observations from VEGETATION-SPOT4.The model run with climate forcing data and vegetation index values corresponding to the period April 1998 to April 1999 (hereafter referred to as year 1998).This arbitrary choice of a different NEE model induces large differences between daily prior and true fluxes.The a priori vs. true differences follow approxi-mately a Gaussian distribution with a standard deviation of 1.3 gC m −2 day −1 .We doubled this value to define the uncertainty on prior daily NEE (3 gC m −2 day −1 on each gridcell), considering 1) that the tail of the distribution is larger than the ideal Gaussian case, 2) that the a priori vs. true NEE differences significantly vary through time, and 3) that in the real world we do not know the a priori flux error characteristics precisely.Note that Chevallier et al., 2006 found a similar value (∼2 gC m −2 day −1 ) by comparing ORCHIDEE and measured NEE (eddy-covariance data) at 34 eddy covariance sites worldwide.
When the fluxes are optimized in each grid cell, prior flux error covariances are crucial in propagating the information given by the station network.In a real data inversion, the prior flux errors have multiple causes that strongly depend on the underlying model.Some of these causes are likely to induce a large-scale spatial error correlation, for example, structural biases of the underlying vegetation model (e.g.parameter values or vegetation classification) or large-scale biases in the forcing data (Jung et al., 2007).On the other hand, meteorological events like frontal systems will certainly decorrelate daily flux errors between nearby regions if the response of the NEE model to changing meteorology is not perfect.Moreover, for a given day, errors in the meteorological forcing (e.g.heterogeneous cloud cover) are likely to cause random errors in prior NEE.Overall the sign and magnitude of the prior flux error correlations are difficult to assess with real data.Most global inversion studies have prescribed spatially correlated errors following an exponential decay with the distance between pixels when solving for weekly or monthly fluxes (Rödenbeck et al., 2003, Peylin et al., 2005).However when solving for daily fluxes, there seems to be no clear evidence of such long-range spatial error correlations, although temporal error correlations seem to be significant (Chevallier et al., 2006).
In our pseudo-data experiment, the error structure of the prior fluxes can be computed from the differences between TURC (prior) and ORCHIDEE (true).The autocorrelation in time of the TURC vs. ORCHIDEE differences shows for each pixel an exponential decrease with variations in the decreasing speed.We choose to define a long decay time because the correlations are reduced during the process to calculate the covariance matrix (see Sect. 2.5.3).The longest autocorrelated pixels have R 2 dropping down to 0.3 after 10 days.We thus defined exponentially decreasing temporal correlations for prior NEE errors, with a decay time of 10 days.For spatial correlations of NEE errors, Chevallier et al., 2006 show no significant correlations but, as for time correlations, we want to maximize them.So we use an exponentially decreasing function with distance with an e-folding length of 1000 km, representing the typical size of synoptic meteorological systems.We do not consider crosscorrelations between space and time to avoid computational complications (see below).In a companion paper (CA08) we test the impact of different NEE error correlation structures.

Eastern North Atlantic air-sea fluxes and errors
Only the north-east Atlantic fluxes are optimized in this study in addition to European daily NEE (Fig. 1).Over this ocean region, the prior air-sea fluxes are set to zero as for the generation of pseudo-data with a total regionally averaged uncertainty of 0.05 GtC year −1 (13 10 6 km 2 ).To calculate the grid point uncertainty, we first distribute the total uncertainty over all grid points assuming a constant uncertainty per day and per squared meter.Then, we build the full covariance matrix with spatial and temporal correlations.We then uniformly rescale the variances so that the total uncertainty of the full matrix is equal to 0.05 GtC year −1 .This is equivalent to an error of 0.5 gC m −2 day −1 over each grid point when taking into account the flux error correlations over the north-east Atlantic.Such a small regional error follows the hypothesis that fluxes outside of Europe are well constrained, so that only small adjustments of the "upwind" flux over the North Atlantic region are allowed.Prior air-sea flux error covariances between ocean grid points are set using an exponential decrease with a length scale of 1500 km in space and 10 days in time.

Calculation of the error covariance matrix
The calculation of a full covariance matrix, based on true minus prior fluxes, although possible, is rather difficult and uneasy to implement given the size of the inverse problem (more than 2 10 6 parameters (see companion paper CA08 for several tests on this).Thus, in a first approximation, we choose to disjointedly define spatial and temporal covariance matrices and to add them, thus neglecting cross-covariances between space and time.In other words, there are no correlations between flux errors unless the fluxes have either a common location or a common time.In this case, the resulting spatial and temporal correlations are necessarily reduced compared to the original "space-only" or "time-only" correlation matrix, in order to remain mathematically consistent.Technically the prior flux covariance matrix (P b ) is defined in four steps: Step 1.We calculate the total daily flux variance for each grid cell using a standard deviation of 3 gC m −2 day −1 for Europe and 0.5 gC m −2 day −1 for the eastern North Atlantic region.
Step 2. We define spatial and temporal correlations between land or ocean grid points in two separate matrices (matrix S' for spatial correlations and matrix T' for temporal correlations).
Step 3. We convert correlations into covariances matrices, S and T, by multiplying S' and T' by the variances from step 1.
Step 4. We compute the total a priori flux covariance matrix as P b = 1 2 [S+T].Note that spatial and temporal correlations in P b are divided by two compared to S and T, reaching a maximum value of 0.5 only.The variances are unchanged.With this approach, the total European a priori flux uncertainty is 0.15 GtC/y and 0.05 GtC/y over eastern North Atlantic.Additional choices for the structure of P b are tested in a companion paper (CA08).

Sequential inverse procedure
Solving for daily fluxes over each model grid point yields to ∼2 700 000 unknown fluxes each year.The number of daily observations is comparatively very small (3650 for 10 stations each year).Therefore, we choose a "data-oriented" expression (see Tarantola, 1987, page 70) to calculate the estimated fluxes X a : with X b the prior fluxes, Y 0 the observations, H the model response functions, R 0 and P b the observation and prior error covariance matrices.In this expression, the matrix to invert has the dimension of the observation space (3650×3650).The critical step is to compute the product HP b H T because P b and H are matrices of very large dimensions, 2 700 000 2 and (3650×2 700 000) elements respectively.This product can be decomposed using only few lines of H at a time and the result stored.The size of the resulting file could be reduced by taking advantage of the structure of these matrices.P b is symmetric but also very sparse as we do not use cross-correlations between space and time; and the size of H can be roughly divided by 2 as half of the matrix is filled by zeros.To further reduce the size of the inverse problem, we choose to use a sequential approach with consecutive inversions.Peylin et al., 2005 showed that the influence of the initial conditions is critical for the first 15 to 20 days.We thus choose windows of three months for the inversion of daily fluxes with a two months overlap between consecutive windows.Each three-month inversion sequence is shifted one month ahead from the previous one, and only estimated fluxes for the month two are kept to build estimated fluxes for all the year except for the beginning and the end of the year.
This sequential approach requires us to carry the influence of all fluxes from previous sequences.Rigorously, for a sequence starting at month "m" and ending at month "m + 3" this means transporting forward the initial conditions (i.e.resulting from all fluxes prior to month "m") in order to compute their contribution at each station for the period "m"-"m + 3".This process would involve too many forward transport simulations.To reduce computing time, we consider that, if we use a long enough spin-up time, the influence of initial conditions at each station can be approximated by an overall constant offset.Thus, we solve for this offset at each sequence to account for past sources.In addition, we replace the prior fluxes for month "m" by the estimates of the previous sequence.We also linearly decrease www.atmos-chem-phys.net/10/3107/2010/Atmos.Chem.Phys., 10, 3107-3117, 2010 the prior flux variance for month "m" from the standard value at the beginning of the month (3 gC m −2 day −1 over Europe and 0.5 gC m −2 day −1 over northeastern Atlantic) to 10% of that value at the end.This forces the system to start month "m + 1" with fluxes close to those optimized at the end of month "m" in the previous sequence.This prevents unrealistic flux variations between successive months.We checked for two following sequences that this simplified approach provides similar results as compared to the rigorous treatment of the influence of past sources.Note finally that for the first sequence, initial conditions are not solved, as in Peylin et al., 2005, given that the impact is limited to the first 20 days.

Daily Fluxes at the transport model grid scale
We first evaluate the potential of the chosen network to quantify CO 2 fluxes at model grid scale in Western Europe (Fig. 1, region in blue).This region has the highest density of observations.As the problem is still under-constrained, with only a few observations for more than 2000 unknowns each day, we expect little constraint on individual grid-point fluxes.Differences between the estimated (or prior) and the true daily fluxes for all the grid points of Western Europe are summarized using correlation (R) and normalized standard deviation (NSD) statistics.NSD is calculated as the ratio between the yearly standard deviation of the fluxes and the yearly standard deviation of the true fluxes.We prefer these two statistics to the conventional RMS diagnostic in order to separate mismatches in the phase (R) and amplitude (NSD) between retrieved and true fluxes.In the ideal case where R=1 and NSD=1, there is a perfect match between optimized fluxes and true fluxes.Averaged values of R and NSD over Western Europe are reported in Table 2.As expected, R and NSD at the grid scale level are very low.The inversion even degrades the correlations, with a mean a posteriori correlation, R APO =0.43±0.21(error from standard deviation of correlation over all grid points) compared to a mean a priori correlation, R APR =0.62±0.12.However, the inversion improves the NSD at the grid scale level.
In order to separate synoptic from seasonal variations, the prior, estimated and true fluxes are de-seasonalized in each grid point.A smooth curve comprising 4 harmonics and a 2nd-order polynomial is first subtracted from each daily NEE time series (Thoning et al., 1989).Then, residuals are smoothed in time with a low pass filter at 80 days to create deseasonalized fluxes.For the deseasonalized prior fluxes, R APR drops to almost zero and NSD APR remains large (1.15) showing that the correlations between prior and true fluxes over each grid-point are dominated by the seasonal cycle of NEE.These differences between synoptic variations in fluxes are not surprising because two distinct NEE models forced with climate data from different years were used for estimating the prior and true fluxes.For the deseasonalized optimized fluxes, the correlation at the grid scale level remain very small on average (R APO =0.11) while the NSD further increases from the prior.This result indicates that unrealistically large day-to-day NEE variations are introduced by the inversion in order to match the pseudo-data.Despite the spatial and temporal flux error correlations, these daily variations at the grid scale level are clearly too large and not even in phase with true fluxes.Our choice for the prior covariance matrix P b , with rather low spatial and temporal correlations, may limit the constraint on the amplitude of the estimated fluxes (see CA08).The constraint delivered by the network of ten stations, despite small observational errors, for NEE on each grid-point is thus quite poor on a daily basis.
As an example, we display the different flux time series for a specific pixel (called SP) located in Germany (5 • 10 E, 47 • 36 N), in the middle of a ring of 5 stations (SCH, CBW, SAC, PUY, PRS).The SP pixel starts out with a pretty good R APR (0.58) and a too large NSD APR of 1.37 (Table 2).At this location, TURC seasonal cycle differs from that of ORCHIDEE especially during the spring uptake.The TURC peak to peak amplitude is higher than in ORCHIDEE (Fig. 2a).R and NSD statistics for that particular pixel reflect those obtained on average over Western Europe (Table 2).The low NSD APO (0.72) is mainly due to a smaller seasonal cycle in the optimized fluxes than in the true fluxes at SP.An analysis of the residuals (Fig. 2b) shows that the exaggerated NSD and small correlation with the true fluxes result from noisy day-to-day variations in optimized fluxes while the true fluxes show only a few synoptic events.Even for a favorably located grid-point, the inversion cannot retrieve the day-to-day variations of NEE.In summary, the 2001 European network does not contain enough information to reliably estimate daily CO 2 fluxes at the grid-scale level, even under the optimistic assumption of fairly small data errors (i.e.all sites are perfectly modeled).The introduction of a priori spatial and temporal error correlations does not compensate for the under-constrained nature of this inverse problem.
In the following, we investigate whether aggregation in space and time of the retrieved fluxes can turn the inversion results from useless to useful.Seven of our ten stations are located in the region "Western Europe" of Fig. 1.Therefore, from now on, we focus our analysis on Western Europe and on the inversion accuracy for retrieving residual deseasonalized CO 2 fluxes rather than seasonal daily fluxes.

Effects of aggregation of inversion results in space and time
For each pixel within Western Europe (Fig. 1), we calculate R and NSD as a function of spatial and temporal flux aggregation levels (Fig. 3).Starting from each grid-point (40 km), fluxes of the neighboring land grid-points are progressively added to form flux aggregates and these aggregated fluxes are then deseasonalized.The process is repeated until all Western Europe becomes a big aggregate.At each step, R and NSD statistics are calculated for aggregated fluxes.The impact of temporal aggregation is estimated by smoothing in time the different spatially aggregated fluxes.We use a boxcar average of width ranging between 1 and 17 days.Longer temporal aggregations are not considered as they would be meaningless for already deseasonalized aggregated fluxes.After interpolation, the aggregated inversion fluxes show a regular evolution of R and NSD (Fig. 3), as a function of  temporal aggregation in the range 1 to 17 days, and of spatial aggregation, in the range of 40 km (grid point) to ≈1200 km (Western Europe).In addition, we estimated the statistical significance of the correlations and variance differences at all aggregation scales, based on Gaussian law and F-variance tests, respectively (see Saporta, 1990, p. 136 and 329).For time series of 365 points (one year of daily fluxes), we calculated a 95% confidence interval for correlations of ±0.1  which means that prior and estimated correlations differences higher than 0.2 are statistically significant.
Comparison between prior deseasonalized fluxes and true deseasonalized fluxes shows only a small improvement with spatial or temporal aggregation alone (Fig. 3a).Overall, correlations only increase from 0.05 to a maximum of 0.35 and NSD values increase from 1.25 to 1.45 degrading with spatial aggregation.The use of two different years of meteorological forcing for TURC (1998) and ORCHIDEE (2001) probably explains why synoptic events do not match even at a rather large scale (after aggregation).However, the statistics of optimized deseasonalized fluxes are clearly improved under combined time and space aggregations.R APO increases from 0.15 up to 0.75 and NSD varies between 0.9 and 1.2 for temporal aggregation scales longer than a few days and spatial aggregation lengths larger than a few hundred kilometers.More precisely, for too short time aggregations (<4 days), the isolines for NSD are almost parallel to the axis of spatial aggregation (Fig. 3).This indicates that temporal aggregation is the limiting factor at all spatial scales in this case.For longer time aggregation (>4 days), temporal aggregation is still the limiting factor but only for spatial aggregation lengths larger than 500 km.Below 500 km, the isolines for NSD become parallel to the temporal axis, indicating that spatial aggregation is the limiting factor for long time aggregations, especially after 10 days.This shows that daily grid-point fluxes need to be aggregated both in space and in time to produce realistic and accurate flux estimates.These improvements can be quantified by the differences between the statistics before and after inversion (Fig. 3c).Improvements from the prior remain small at low aggregation levels (0.2/0.15 for R/NSD around 500 km and 4 days) and only become substantial around 10 days and 1000 km (0.35/0.45 for correlation/NSD).At low temporal aggregation (<4 days), NSD is even slightly degraded after the optimization, indicating that the inversion introduces large short-term flux variations to match the pseudo-data.The 95% confidence interval in R shows that estimated correlations are different from the prior ones for temporal aggregation longer than 3 days and spatial aggregation larger than 200 km (Fig. 3c).At a 95% level, the variance of prior residuals is statistically different from the variance of the true residual for spatial aggregations higher than 500 km and all temporal aggregations.On the other hand, the estimated variances are not statistically different from the true variances for temporal aggregations longer than 3 days and all spatial aggregations (not shown), indicating thus a statistically significant improvement above 3 days and 200 km.
If one regards the inversion as successful to retrieve true fluxes if R >0.7 and NSD≈1, the network of 10 continuous stations over Europe (assumed to be perfectly modeled) allows us to retrieve the true fluxes over Western Europe at a spatial resolution greater than 1000 km and a temporal resolution greater than 10 days in agreement with the statistical significance analysis (Fig. 3).
Figure 4 displays the different fluxes over Western Europe smoothed at 10 days.Although the amplitude of the seasonal cycle is still underestimated compared to the truth, the 10day smoothed aggregated residual fluxes show a good agreement with the true deseasonalized fluxes (R=0.63,NSD=1.0,Table 2), a marked improvement from the prior (R=0.3,NSD=1.5).Thus, even with a prior estimate inconsistent with the truth (different land model and different years of driving meteorology), the current network of stations allows us to correct the prior and to successfully retrieve NEE on a 10-day average basis, for the relatively well-constrained and large Western Europe region.Only a few flux variations such as those observed in May are not retrieved well.Note that with real data, the prior fluxes may be as far from the unknown truth as in this particular case, which reinforces the general character of this finding.

Flux improvement over other European regions
Minimum requirements in terms of spatial and temporal aggregation to get satisfying flux retrieval depend on the region under consideration.For other large regions in Europe (Fig. 1), aggregation generally brings smaller improvements than for the Western Europe region, as fewer measurement sites constrain these regions.In some regions, the fluxes cannot be retrieved satisfactorily with the criteria defined above (R >0.7 and NSD=1).This is the case for the Mediterranean Europe region (Fig. 5a), with a maximum R APO of 0.5 (12-day aggregation over the whole region) and a large NSD for all aggregations.NSD is even significantly degraded from the prior, indicating again large unrealistic dayto-day variations in the inverted fluxes (Fig. 5a).Other regions show a relatively good agreement.Central Europe shows smaller R APO than Western Europe but minimum requirements can be fulfilled with a 15-day averaging of the residual fluxes over the whole region (Fig. 5b).For Scandinavia (not shown), the results are only slightly degraded as compared to Western Europe.However, the relatively good agreement between estimated and true fluxes in Scandinavia is due to the initial agreement between prior and true fluxes and only partially due to additional information delivered by atmospheric data.

Conclusions
We have built an inverse model to infer daily CO 2 fluxes over the European continent using continuous daily concentration observations and prior information on surface fluxes.This model is a classic Bayesian inversion.We used prior flux error correlations in space and time.Both correlations are exponentially decreasing and characterized by a correlation length.In a first attempt, we used a simplified prior covariance matrix by neglecting cross-correlations between space and time.This simplification mathematically implies a reduction of the correlations when combined together.We thus used longer correlation lengths, both in time and space, to counterbalance.We avoid the problem of dealing with huge matrices (H and P) by using overlapping time windows inversions of three months, with two months overlapping.
We have shown that, in the idealized perfect transport case, where the transport model operator can accurately represent each site (±0.3 ppm of white noise added to daily CO 2 pseudo-data), it is not possible to reconstruct European fluxes each day at the model resolution.Given the overwhelming number of unknown fluxes compared to the amount of data, it is clear that aggregation is necessary.When aggregating the inversion fluxes in space (∼1000 km) and time (∼10 days), the CO 2 budget of Western Europe, the best-observed region, can be retrieved.Other European regions where the network is less dense, show a limited ability of the inversion to retrieve the true fluxes, even at regional scale, highlighting the need of dense regional atmospheric observation networks.In this experiment, we only consider the prior flux errors and measurement errors.But, in real data inversions, there are several other sources of uncertainties: the transport error, the presence of uncertain fossil fuel emissions.Thus, the retrieval fluxes from real data inversion would be degraded compared to this experiment.
Extending the area where CO 2 fluxes could properly be retrieved in Europe requires enhancing the development of the atmospheric continuous network, especially in Eastern Europe and around the Mediterranean basin.The use of tall towers and small aircraft can also bring additional information from meso to regional scales (www.carboeurope.org).The potential of combining continuous surface measurements with upcoming satellite measurements is also a promising perspective, but it will require updated inversion tools, capable to handle larger volumes of data and associated error covariances.
Further improvements can come from more accurate prior flux scenarios, including knowledge of prior flux error covariance, both for NEE and for fossil fuel emissions.For the former, coupling between atmospheric models and vegetation models, developed for instance with mesoscale models (Lauvaux et al., 2008) or the use of Carbon Data Assimilation Systems (CCDAS) (Rayner et al., 2005)  for improving ecosystem modeling and atmospheric data fusion.For the latter, efforts have been recently made to improve fossil fuel emission scenarios and to provide hourly fossil fuel CO 2 emission maps (IER: http://carboeurope.ier.uni-stuttgart.de/,EDGAR: Van Aardenne, 2005).Differences between these scenarios have a significant seasonal impact on the concentration at most European stations (Peylin et al., 2008).
Maximizing the benefits of these new atmospheric and emission constraints requires progress on the quality of modeled transport, and integration of relevant scenarios for error correlations in inversions, both in the flux and in the observation domains.In the companion paper of this work (CA08), we propose a first analysis of the impact on inversion results, of transport model error and of scenarios for flux error correlations.
C. Carouge et al.: Potential of the 2001 network Figure 1.

Fig. 1 .
Fig. 1.Map of the 2001 European continuous stations represented by black triangles.After inversion, fluxes were aggregated over five different regions: "Western Europe" in blue, "MediterraneanEurope" in orange, "Balkans" in light green, "Central Europe" in red and "Scandinavia" in green.

Fig. 3 .
Fig. 3. Evolution of correlation and NSD with spatial and temporal aggregation over Western Europe for prior residual fluxes (a), posterior residual fluxes (b) and the difference between these panels (c).

Fig. 5 .
Fig. 5. Evolution of correlation and NSD with spatial and temporal aggregation over Mediterranean region for prior (a) and posterior (b) residual fluxes, and over Central Europe for prior (c) and posterior (d) residual fluxes.
are promising ways C. Carouge et al.: Potential of the 2001 network

Table 1 .
List of European continuous stations with their position.

Table 2 .
Correlation and NSD for raw and residual fluxes for SP pixel and aggregated flux over Western Europe and their mean over all pixels of Western Europe and their standard deviation.SP pixel (5• 10 E, 47 • 36 N) is chosen in the middle of a ring of 5 stations as representative of a well-constrained pixel.