# Joint assimilation of eddy covariance flux measurements and FAPAR products over temperate forests within a process-oriented biosphere model

C. Bacour, P. Peylin, N. Macbean, P. Rayner, F. Delage, F. Chevallier, M. Weiss, J. Demarty, D. Santaren, F. Baret, et al.

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

## ► To cite this version:

C. Bacour, P. Peylin, N. Macbean, P. Rayner, F. Delage, et al. Joint assimilation of eddy covariance flux measurements and FAPAR products over temperate forests within a process-oriented biosphere model. Journal of Geophysical Research: Biogeosciences, 2015, 120 (9), pp.1839 - 1857. 10.1002/2015JG002966. [hal-01805695](https://hal.science/hal-01805695)

HAL Id: hal-01805695

https://hal.science/hal-01805695

Submitted on 28 May 2020

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published 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_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
RESEARCH ARTICLE

10.1002/2015JG002966

**Key Points:**

* Compatibility of in situ NEE and LE data and various FAPAR products through model-data fusion

* FAPAR mainly constrains phenology; when assimilated alone, it may degrade the modeled fluxes

* Combining the two data streams is preferable for improving a process-based vegetation model

**Correspondence to:**

C. Bacour,
cedric.bacour@noveltis.fr

**Citation:**

Bacour, C., et al. (2015), Joint assimilation of eddy covariance flux measurements and FAPAR products over temperate forests within a process-oriented biosphere model, *J. Geophys. Res. Biogeosci.*, 120, 1839–1857, doi:10.1002/2015JG002966.

Received 20 FEB 2015

Accepted 10 AUG 2015

Accepted article online 14 AUG 2015

Published online 30 SEP 2015

Joint assimilation of eddy covariance flux measurements
and FAPAR products over temperate forests within
a process-oriented biosphere model

C. Bacour¹,², P. Peylin², N. MacBean², P. J. Rayner²,³, F. Delage²,⁴, F. Chevallier², M. Weiss⁵,
J. Demarty⁵,⁶, D. Santaren⁷,⁸, F. Baret⁵, D. Berveiller⁹, E. Dufrêne⁹, and P. Prunet¹

¹NOVELTIS, Labège, France, ²Laboratoire des Sciences du Climat et de l'Environnement, CEA/CNRS/UVSQ, Gif-sur-Yvette, France, ³Now at School of Earth Sciences, University of Melbourne, Melbourne, Victoria, Australia, ⁴Now at Centre for Australian Weather and Climate Research, Bureau of Meteorology, Docklands, Melbourne, Victoria, Australia, ⁵UMR1114 EMMAH, Environnement Méditerranéen et Modélisation des Agro-Hydrosystèmes, INRA, Avignon, France, ⁶Now at UMR HydroSciences Montpellier, Université de Montpellier2, Montpellier, France, ⁷Environmental Physics, Institute of Biogeochemistry and Pollutant Dynamics, ETH Zurich, Zurich, Switzerland, ⁸Now at Laboratoire des Sciences du Climat et de l'Environnement, CEA/CNRS/UVSQ, Gif-sur-Yvette, France, ⁹Laboratoire Ecologie Systématique et Evolution, CNRS, Université Paris-Sud AgroParisTech, Orsay, France

**Abstract** We investigate the benefits of assimilating in situ and satellite data of the fraction of photosynthetically active radiation (FAPAR) relative to eddy covariance flux measurements for the optimization of parameters of the ORCHIDEE (Organizing Carbon and Hydrology in Dynamic Ecosystem) biosphere model. We focus on model parameters related to carbon fixation, respiration, and phenology. The study relies on two sites—Fontainebleau (deciduous broadleaf forest) and Puechabon (Mediterranean broadleaf evergreen forest)—where measurements of net carbon exchange (NEE) and latent heat (LE) fluxes are available at the same time as FAPAR products derived from ground measurements or derived from spaceborne observations at high (SPOT (Satellite Pour l'Observation de la Terre)) and medium (MERIS (MEdium Resolution Imaging Spectrometer)) spatial resolutions. We compare the different FAPAR products, analyze their consistency with the in situ fluxes, and then evaluate the potential benefits of jointly assimilating flux and FAPAR data. The assimilation of FAPAR data leads to a degradation of the model-data agreement with respect to NEE at the two sites. It is caused by the change in leaf area required to fit the magnitude of the various FAPAR products. Assimilating daily NEE and LE fluxes, however, has a marginal impact on the simulated FAPAR. The results suggest that the main advantage of including FAPAR data is the ability to constrain the timing of leaf onset and senescence for deciduous ecosystems, which is best achieved by normalizing FAPAR time series. The joint assimilation of flux and FAPAR data leads to a model-data improvement across all variables similar to when each data stream is used independently, corresponding, however, to different and likely improved parameter values.

**1. Introduction**

The terrestrial biosphere plays a key role in the control of the exchange of energy and matter (in particular carbon and water) between the land surface and the atmosphere [Pielke *et al.*, 1998]. The use of land surface models (LSMs) that describe these main governing processes is of growing importance for improving our understanding of the fate of the terrestrial ecosystems to environmental changes [Pitman, 2003; Sitch *et al.*, 2008]. LSMs rely on generic hypotheses and fixed parameterizations that were derived from a limited number of observations, from the scale of individual plant organs to the scale of the plant community, and under specific environmental conditions. Therefore, large uncertainties remain in their ability to reliably represent the spatial and temporal variations of the ecosystem characteristics and the carbon cycle under current or future climate conditions [Field *et al.*, 1995; Friedlingstein *et al.*, 2006; Wullschleger *et al.*, 2014]. Data assimilation techniques are increasingly used to reduce these uncertainties by improving the model parameters [Wang *et al.*, 2001; Kaminski *et al.*, 2013] while also highlighting possible model deficiencies [Verbeeck *et al.*, 2011; Kuppel *et al.*, 2012; Keenan *et al.*, 2013].

In this context, in situ eddy covariance flux measurements have mainly been used to constrain the model parameters controlling the processes of carbon and water exchange [Wang *et al.*, 2001; Braswell *et al.*, 2005; Knorr and Kattge, 2005; Santaren *et al.*, 2007; Moore *et al.*, 2008; Williams *et al.*, 2009; Groenendijk *et al.*, 2011; Kuppel *et al.*, 2014]. Eddy flux data alone may not be sufficient to disentangle different concurrent
processes, as for instance the partitioning of net carbon flux into gross photosynthesis and ecosystem respiration [Sacks et al., 2006], and the assimilation of additional biometric ground-based measurements has already been proven to reduce the model uncertainty for simulating ecosystem carbon exchange [Richardson et al., 2010; Ricciuto et al., 2011]. The number of instrumented sites is however limited, though increasing, and some ecosystems remain poorly monitored [Baldocchi, 2008; Williams et al., 2009].

On the other hand, satellite-derived biophysical products cover the whole globe over a time period unrivaled by local measurements of terrestrial ecosystems. The history of the use of remotely sensed products of vegetation photosynthetic activity (FAPAR) in semidiagnostic ecosystem models to estimate the gross primary productivity is rich [McCallum et al., 2009; Seixas et al., 2009; Jung et al., 2011; Cheng et al., 2014]. Most of these studies rely on the concept of light use efficiency [Monteith, 1977]. The use of satellite products in assimilation to constrain the parameters of process-based models is, however, more limited whilst expanding. Assimilated satellite-derived observations range from raw reflectance data [Quaife et al., 2008] and vegetation indices derived from reflectance data [Migliavacca et al., 2009] to higher added value products monitoring the phenology of vegetation activity more directly [Stöckli et al., 2008; Knorr et al., 2010]. Finally, atmospheric CO₂ mole fraction measurements have also been used to provide large-scale constraints on the model parameters related to net ecosystem exchange (NEE) [Rayner et al., 2005].

There have been relatively few attempts to assimilate simultaneously satellite-derived products of vegetation activity and in situ measurements, despite the potential increased constraint on the model processes brought by adding different sources of information. Kaminski et al. [2012] assimilated monthly MERIS (MEdium Resolution Imaging Spectrometer) FAPAR products at very coarse spatial resolution (10° by 8°) together with monthly mean atmospheric CO₂ mole fraction data; They showed an added value in combining the two data streams for improving the hydrological (evapotranspiration) and carbon (net primary productivity) fluxes simulated by the BETHY-TM2 land-atmosphere model at the global scale but only moderate benefits of using FAPAR data to constrain the net ecosystem exchange. Kato et al. [2013] jointly assimilated 10 day SeaWiFs FAPAR data at medium spatial resolution (1.5 km) and in situ latent heat (LE) measurements. Zobitz et al. [2014] combined in situ NEE measurements and spatially averaged MODIS (MODerate resolution Imaging Spectroradiometer) FAPAR 8 day composite products at 1 km at different time scales to optimize an ecosystem process model for a subalpine coniferous forest. The two latter studies point out how assimilating one data stream may degrade the model simulations to the other observable (in particular, the assimilation of FAPAR deteriorates the model fit to the measured flux) and how combining the two may increase the model agreement with the various data sets. At the same time, their findings regarding the time scales at which the combination of these data streams can be beneficial are different. Such multi-data stream optimization studies are hampered by inconsistencies between the measurements, the equations of the LSMs (including the spatial resolution), and the observation operators that simulate the measurements from these LSMs. The diversity of horizontal resolutions used in the models and in the observations also contributes to these inconsistencies, and one can note that only medium (if not coarse) spatial resolution satellite products have previously been used, leaving the potential of remote sensing products at higher spatial resolution unexplored by joint data assimilation studies. Moreover, there is no clear evidence that satellite products are compatible with local-scale measurements of ecosystem functioning, given the different spatial scales involved: in situ flux measurements are representative of the exchange of mass and energy between ecosystem and atmosphere within only a few hundred meters (mostly) around the flux towers, while satellite observations range from about 10 m to often a few kilometers.

In this study, we evaluate the benefit of simultaneously assimilating FAPAR data, satellite-derived products at high and medium spatial resolution and ground-based measurements, along with flux tower measurements of both carbon and water fluxes using the state-of-the-art mechanistic terrestrial biosphere model ORCHIDEE (Organizing Carbon and Hydrology in Dynamic Ecosystems [Krinner et al., 2005]). The approach relies on a variational assimilation system that has been designed for ORCHIDEE (ORCHIDAS), which has already been used at site level to constrain the model using in situ flux measurements [Santaren et al., 2007; Verbeeck et al., 2011; Kuppel et al., 2012; Peng et al., 2013; Santaren et al., 2014]. Our study focuses on two flux tower sites corresponding to deciduous broadleaf (Fontainebleau) and Mediterranean broadleaf evergreen (Puechabon) forests, for which NEE and latent heat (LE) flux measurements are available at a half-hourly time step. Time series of FAPAR estimates are both measured in situ and derived from satellite observations, at the scale of the flux tower footprint from SPOT (Satellite Pour l'Observation de la Terre) observations (at 40 m spatial resolution) and at medium-scale spatial resolution from MERIS observations (at 1 km).
The study focuses on assessing the complementarity of flux and FAPAR data to constrain the mean seasonal cycle of NEE and LE simulated by ORCHIDEE, with no specific investigation of the year to year variations [see Santaren et al., 2014]. For the two sites, however, we use flux measurements from additional years in order to evaluate the model performance after the calibration over 1 year. In order to evaluate the potential of a simultaneous assimilation of flux and FAPAR data, we will address the following questions:

1. What is the impact of assimilating in situ flux data on the simulated FAPAR and vice versa?

2. How important are the differences among the different FAPAR products in terms of their ability to constrain the parameters of the ORCHIDEE model?

3. What is the impact of jointly assimilating flux and FAPAR data in terms of improvement of the model-data misfit?

4. To what extent are the FAPAR products compatible with the in situ flux data measurements and with the ORCHIDEE model?

The paper is structured as follows. Section 2 describes the data used (2.1), the ORCHIDEE model (2.2), and the inversion method and the assimilation configurations considered (2.3). Section 3 presents the results of different assimilation scenarios conducted to answer the above questions and the corresponding impacts in terms of parameter retrieval and uncertainty. The results and their implications are discussed in section 4. Section 5 concludes the paper.

## 2. Materials and Methods

### 2.1. Data

#### 2.1.1. In Situ Flux Measurements

This study focuses on the flux tower sites of Fontainebleau and Puechabon. The Fontainebleau forest site (48.4763°N; 2.7801°E) is located southeast of Paris, France. Deciduous broadleaf trees (oak (*Quercus petraea* and *Quercus robur*), beech (*Fagus sylvatica*), and hornbeam (*Carpinus betulus*)) are the dominant species in the vicinity of the flux tower. The Puechabon site, located in the south of France (43.7414°N; 3.5958°E), is a Mediterranean evergreen forest that is mainly composed of green oak trees (*Quercus ilex*).

Fluxes and meteorological variables were measured at the sites using the standardized CARBO-EUROFLUX protocol [Aubinet et al., 2000] and delivered on a half-hourly basis. The data were processed (correction, gap-filling, and partitioning) using standard methodologies [Reichstein et al., 2005].

In the assimilation studies conducted further (sections 3.1 to 3.4.3), we use only 1 year of NEE and LE data given that a single year of FAPAR products was available at each site; we use measurements performed in 2006 at Fontainebleau and 2004 at Puechabon. In addition, we also use in situ NEE and LE fluxes available for additional years at both sites so as to cross-validate the optimized parameters and assess the improvement of the model optimized using only 1 year of data over additional years (section 3.5): at Fontainebleau, measurements over 2006–2010 are considered and over 2000–2006 at Puechabon.

Although the fluxes are measured at a half-hourly time step, we only use their daily means in the assimilation procedure in order to be consistent with the temporal constraint brought by the FAPAR data (section 2.1.2) and given our focus on the seasonal cycle. Days with more than 20% of missing half-hourly observations were removed from the analysis. The flux time series (both observed and simulated) are also further smoothed using a centered 15 day moving average window in order to remove high-frequency variations in the data that are not the focus of this study and that may complicate the optimization of the mean seasonal cycle [see Santaren et al., 2007]. We also use gross primary productivity (GPP) estimates made from the flux measurements using the flux partitioning method of Reichstein et al. [2005] as a diagnostic to assess the constraint brought by each data stream on the ecosystem productivity.

#### 2.1.2. FAPAR Products

In situ measurements of FAPAR are available for both sites. For Fontainebleau, the temporal monitoring of FAPAR is inferred from the measurements of incoming, reflected, and transmitted PAR (photosynthetic active radiation), the latter being performed by seven sensors distributed under the canopy to sample the local spatial heterogeneity. For Puechabon, FAPAR is derived from 14 sensors on the ground. Note that these in situ FAPAR measurements are not fully representative of the tower footprint which may extend several kilometers
under low-turbulence conditions. Satellite FAPAR estimates were derived from SPOT and MERIS 10 × 10 km² reflectance images encompassing the tower sites. They were derived from a Neural Network estimation algorithm trained with radiative transfer model simulations of FAPAR [Bacour et al., 2006; Baret et al., 2007; Weiss et al., 2007] using the PROSPECT + SAIL reflectance model of vegetation canopy [Jacquemoud et al., 2009] coupled to the SMAC atmosphere radiative transfer model [Rahman and Dedieu, 1994]. FAPAR is simulated by PROSPECT + SAIL considering sun geometries at 10:00 local solar time, which corresponds to a good approximation of the daily integrated value during clear-sky days [Baret et al., 2007]. For each instrument, the Neural Network algorithm is fed with the observed reflectances in each spectral band and with the corresponding observation geometries at 10:00. Finally, FAPAR is estimated on a pixel basis for each date of acquisition.

FAPAR products for SPOT are derived at 40 m resolution from the aggregation of original SPOT pixels at 20 m in order to smooth geometrical errors induced by geolocation and coregistration. The final values correspond to the mean of the FAPAR values associated with the pixel of the flux tower and those of the eight surrounding pixels. For MERIS products, the mean is performed over the 1 km pixel encompassing the flux tower and four neighboring pixels having the same land cover composition (based on a classification performed at high spatial resolution using SPOT images). Compared to SPOT, MERIS observations offer a denser temporal monitoring of the sites as MERIS products are provided every 8 days, whereas only four acquisition dates for SPOT are available for Fontainebleau (and none over Puechabon). To overcome the scarcity of SPOT temporal monitoring at high spatial resolution, SPOT FAPAR products at Fontainebleau were extrapolated in time based on MERIS data (and hence delivered at an 8 day frequency). For this site, SPOT and MERIS FAPAR products are therefore not fully independent, but, as they show strong discrepancies, they are used further to illustrate the model response following the assimilation of such contrasting data streams. In order to remove the high-frequency noise and not to focus on synoptic variations, the FAPAR products are smoothed using a moving average window of ±20 days around the centered date, with a Gaussian function truncated at ±20 days. Weekly mean values were used in the assimilation for all FAPAR products.

We evaluate further the use of normalized FAPAR time series in the assimilation for the Fontainebleau site (deciduous forest). The normalization is performed by scaling the data between the 5th and 95th percentiles of their distribution (to avoid spuriously high and low data).

## 2.2. ORCHIDEE Vegetation Model

### 2.2.1. Model Description

ORCHIDEE is a state-of-the-art mechanistic vegetation model that simulates the exchanges of carbon dioxide, water, and heat fluxes within the soil-vegetation-atmosphere continuum and the evolution of water and carbon pools [Krinner et al., 2005]. The model describes the exchanges of water, carbon, and energy, between biosphere and atmosphere, as well as the soil water budget, at a half-hourly time step, and the slow components of the terrestrial carbon cycle (including carbon allocation in plant reservoirs, soil carbon dynamics, and litter decomposition) on a daily basis. It is part of one of the Earth System Models that were used to assess future climate changes for the fifth assessment report of Intergovernmental Panel on Climate Change [Ciais et al., 2013]. As in most global vegetation models, vegetation is described for main plant functional types (PFTs), with 12 different types of vegetation plus bare soil. All processes in ORCHIDEE follow the same governing equations, except the phenology, but many parameter values are PFT dependent. The leaf phenology models follow Botta et al. [2000]. The start of the growing season for deciduous PFTs follows a typical Growing Degree Day (GDD) model and is triggered when the daily calculated GDD exceeds a calculated threshold: For deciduous trees, that threshold is an increased function of the number of chilling days in order to fulfill chilling physiological requirements, while for natural C3 grass, the GDD threshold depends on the long-term mean annual air temperature. For deciduous trees which are driven by sensitivity to cold temperatures, the senescence begins when the monthly air surface temperature goes below a threshold temperature. In addition, the turnover of the leaves is continually affected by the aging of the leaves during the growing season for all tree PFTs. Details of the senescence processes are provided in MacBean et al. [2015]. FAPAR is computed from the simulated leaf area index (LAI) using the classical Beer-Lambert law for the extinction of the direct illumination within the canopy:

$$FAPAR(t) = 1 - \exp(-0.5 \times LAI(t)) \quad (1)$$
We have chosen to use a constant extinction coefficient fixed to 0.5 (that corresponds to a spherical distribution of leaves with an illumination angle at nadir) as ORCHIDEE equations related to photosynthesis have been established and calibrated based on this assumption.

### 2.2.2. Model Setup

For the following simulations, we consider that the vegetation at Puechabon is mainly temperate evergreen broadleaf forest (TempEBF) with a small portion of bare soil (10%). For Fontainebleau, it consists mainly of temperate deciduous broadleaf forest (TempDBF hereafter), and we account for a small fraction of C3 grass to represent both the vegetation understory in the vicinity of the flux tower and the impact of neighboring grass plots in satellite observations (10%). ORCHIDEE is forced by local measurements of the meteorological fields (incoming shortwave and longwave radiations, near surface air temperature, specific humidity, atmospheric pressure, precipitation, and wind speed). Spin-up runs were performed by cycling ECMWF (European Centre for Medium-Range Weather Forecasts) meteorological forcing for 2001–2006 over a 3000 year period, in order to bring the different soil carbon reservoirs to realistic values. Note that the spin-up includes a transient simulation anterior to the year of measurements at each site in order to account for the increase of atmospheric CO₂ concentrations. Although this step is mandatory as not enough data exist to calibrate the different soil carbon pools (following the CENTURY model [Parton et al., 1987]) at these two sites, the drawback is that the spin-up puts the model in a steady state (the mean annual NEE is close to zero) which is not the case for the two sites considered. The assimilation is expected to correct for this bias.

## 2.3. Assimilation Methodology

### 2.3.1. System Description

The ORCHIDAS variational data assimilation system allows NEE and LE eddy covariance measurements and FAPAR estimates (in situ or remotely sensed) to be assimilated simultaneously in order to optimize the parameters of ORCHIDEE. The latest description of the data assimilation system is provided in Kuppel et al. [2012]. It relies on a Bayesian framework with the hypothesis of Gaussian errors, which leads to the minimization of the following misfit function $J(\mathbf{x})$, that compares (1) the observations $\mathbf{y}$ and corresponding model outputs $\mathbf{H}(\mathbf{x})$ and (2) a priori (background) $\mathbf{x}_b$ and optimized parameter $\mathbf{x}$ values. Both terms are weighted by the error covariance matrices of the observations $\mathbf{R}$ and of the parameters ($\mathbf{B}$) [Tarantola, 2005]:

$$J(\mathbf{x}) = (\mathbf{y} - \mathbf{H}(\mathbf{x}))^T \mathbf{R}^{-1} (\mathbf{y} - \mathbf{H}(\mathbf{x})) + (\mathbf{x} - \mathbf{x}_b)^T \mathbf{B}^{-1} (\mathbf{x} - \mathbf{x}_b) \quad (2)$$

The determination of the optimal set of parameters that minimizes $J(\mathbf{x})$ is performed with a variational approach, based on the L-BFGS-B algorithm [Byrd et al., 1995], which is specifically dedicated to solving large nonlinear optimization problems subject to simple bounds on the parameters. This L-BFGS-B quasi-Newton optimization algorithm explores each parameter space simultaneously along the gradient of the misfit function and uses an approximation of the Hessian (second derivative) of $J$ that is updated at each iteration. The latter is based on the estimation of the gradient of the misfit function, which requires the derivative of each model output with respect to the parameters. This is provided by the tangent linear model of ORCHIDEE that was automatically generated by the numerical TAF tool (Transformation of Algorithms in FORTRAN) [Giering et al., 2005] for all parameters except $K_{\text{pheno\_crit}}$ and $T_{\text{senes}}$ (see Table 1). As these parameters are involved in threshold functions, the tangent linear model may be null for certain parameter values. To overcome this, the gradient of the misfit function with respect to these two parameters is determined using a finite difference approach [Santaren et al., 2007].

As compared to global search methods, e.g., Monte Carlo approaches, of which computational cost may become too prohibitive with an increasing number of parameters, a gradient-based optimization algorithm such as L-BFGS-B may converge toward a local minimum rather than toward the global minimum of the misfit function [Santaren et al., 2014]. However, we have run several tests using different sets of first guess parameters and have checked that the use of a single initial parameter set did not affect the interpretation of the results.

Only diagonal elements are accounted for in the prior error variance-covariance matrix $\mathbf{B}$. The standard deviation of the uncertainty for each parameter is set to 40% of an allowed physical range of variation (defined from expert knowledge and literature survey) provided to L-BFGS-B for the optimization (Table 1). To improve the minimization efficiency of the algorithm, the minimization is preconditioned by providing L-BFGS-B with
<table><caption>**Table 1.** ORCHIDEE Parameters to be Optimized, Their Prior Value, and Their Authorized Variation Interval for the Temperate Deciduous Broadleaf Forest, C3 Grass, and Temperate Evergreen Broadleaf Forests, Vegetation Types</caption><thead><tr><th rowspan="3">Parameter Name</th><th rowspan="3">Description (Unit)</th><th colspan="3">Min/Prior/Max</th></tr><tr><th colspan="2">Fontainebleau</th><th>Puechabon</th></tr><tr><th>TempDBF</th><th>C3grass</th><th>TempEBF</th></tr></thead><tbody><tr><td colspan="5"><strong>Photosynthesis</strong></td></tr><tr><td>V<sub>cmax</sub></td><td>maximum carboxylation rate (µmol m<sup>-2</sup> s<sup>-1</sup>)</td><td>30/55/80</td><td>38/70/102</td><td>25/45/65</td></tr><tr><td>G<sub>s,slope</sub></td><td>Ball-Berry slope</td><td>6/9/12</td><td>6/9/12</td><td>6/9/12</td></tr><tr><td>T<sub>opt</sub></td><td>optimal photosynthesis temperature (°C)</td><td>18/26/34</td><td>19.25/27.25/35.25</td><td>24/32/40</td></tr><tr><td>SLA</td><td>specific leaf area (m<sup>2</sup> g<sup>-1</sup>)</td><td>0.013/0.026/0.05</td><td>0.013/0.026/0.05</td><td>0.01/0/02/0.04</td></tr><tr><td colspan="5"><strong>Energy Balance</strong></td></tr><tr><td>K<sub>albedo,veg</sub></td><td>multiplicative factor of the vegetation albedo</td><td colspan="3">0.8/1/1.2<sup>a</sup></td></tr><tr><td colspan="5"><strong>Soil water availability</strong></td></tr><tr><td>Hum<sub>cste</sub></td><td>root profile (m<sup>-1</sup>)</td><td>0.2/0.8/3</td><td>1/4/10</td><td>0.2/0.8/3</td></tr><tr><td colspan="5"><strong>Respiration</strong></td></tr><tr><td>Q<sub>10</sub></td><td>temperature dependency of heterotrophic respiration</td><td colspan="3">1/1.99372/3<sup>a</sup></td></tr><tr><td>K<sub>soilC</sub></td><td>multiplicative factor of the initial soil carbon pools</td><td colspan="3">0.5/1/2<sup>a</sup></td></tr><tr><td>GR<sub>frac</sub></td><td>fraction of biomass available for growth respiration</td><td>0.2/0.28/0.36</td><td>0.2/0.28/0.36</td><td>0.2/0.28/0.36</td></tr><tr><td colspan="5"><strong>Phenology</strong></td></tr><tr><td>LAI<sub>MAX</sub></td><td>maximum LAI value</td><td>3/5/8</td><td>1.5/2.5/3.5</td><td>3/5/8</td></tr><tr><td>LAI<sub>init</sub></td><td>value of LAI at t = 0</td><td>0/0/8</td><td>0/2.27/3.5</td><td>1/4.71/8</td></tr><tr><td>L<sub>agecrit</sub></td><td>average critical age of leaves (days)</td><td>90/180/240</td><td>60/120/180</td><td>490/730/970</td></tr><tr><td>τ<sub>leaf,init</sub></td><td>number of days of use of the carbohydrate reserves (days)</td><td>5/10/30</td><td>5/10/30</td><td>5/10/30</td></tr><tr><td>K<sub>pheno,crit</sub></td><td>multiplicative parameter of the threshold that determines the start of the growing season</td><td>0.65/1/1.65</td><td>0.65/1/1.65</td><td>0.65/1/1.65</td></tr><tr><td>K<sub>LAI,happy</sub></td><td>LAI threshold to stop using carbohydrate reserves</td><td>0.35/0.5/0.7</td><td>0.35/0.5/0.7</td><td>0.35/0.5/0.7</td></tr><tr><td>T<sub>senes</sub></td><td>temperature threshold for senescence (°C)</td><td>9/12/16</td><td>-4.375/-1.375/1.625</td><td>-</td></tr></tbody></table>

<sup>a</sup>No PFT dependency.

scaled parameters **x'**, using the a priori values and uncertainties rather than **x** [Chevallier et al., 2005] as this homogenizes the range of variation of the optimized parameters:

$$ \mathbf{x}' = \mathbf{B}^{-1/2}(\mathbf{x} - \mathbf{x}_b) \quad (3) $$

The posterior error covariance matrix of the optimized parameters (A) is computed under the hypotheses of model linearity in the vicinity of the solution. It is given by

$$ \mathbf{A} = [\mathbf{H}_{\infty}^{T}\mathbf{R}^{-1}\mathbf{H}_{\infty} + \mathbf{B}^{-1}]^{-1} \quad (4) $$

where **H**<sub>∞</sub> is the Jacobian matrix associated to the gradients of the model outputs with respect to the parameters at the solution. For a given parameter, the comparison of the posterior errors between scenarios will allow an evaluation of the performance of the assimilation and hence will highlight the contribution that each type of information brings to our knowledge on the system.

### 2.3.2. Model-Observation Uncertainties

The error matrix **R** should account for both the error in the measurements and the error in the model structure (i.e., error in the representation of the processes) [Tarantola, 2005]. For flux data, the measurement error is known to vary in time (depending on the magnitude of the fluxes) and can be estimated as the residual of the gap-filling algorithm [Richardson et al., 2008]. For NEE and LE fluxes, the measurement error is usually small as compared to the model error and has a correlation structure that is negligible on a daily timescale [Lasslop et al., 2008]. Model errors are rather difficult to assess and may be much larger than the measurement error itself. Kuppel et al. [2013] showed that the model error in ORCHIDEE dominates the error budget: for NEE for instance, it is on the order of 1.5–1.7 g C/m²/d when the measurement error is between 0.2 and 0.8 g C/m²/d [Richardson et al., 2008]. For this reason, we choose to define the model-observation uncertainty as the root-mean-square error (RMSE) between the various data and the ORCHIDEE a priori simulations and choose to keep **R** diagonal as have most assimilation studies (similar to Santaren et al. [2014] and Kuppel et al. [2014]).

We also checked that each data stream had similar weight in the prior cost function; a smaller uncertainty (RMSE divided by 2) has been attributed to FAPAR principally because it has fewer available observations (weekly versus daily for each of the two fluxes).
### 2.3.3. Optimized Parameters

Among the many parameters of the ORCHIDEE model, a preliminary sensitivity analysis allowed us to select a subset of 16 parameters to be optimized; they are the main drivers of the net CO₂ and latent heat fluxes and FAPAR dynamics. As our focus is mainly on the constraints brought by FAPAR and flux measurements on the carbon cycle (GPP), only a few parameters related to water and energy exchanges were considered. The list is provided in Table 1, together with the associated range of variation and prior value. The parameters, and the processes they are involved in, are detailed in Kuppel et al. [2012] and Santaren et al. [2014]. Similarly to Santaren et al. [2007], emphasis has been put on processes causing rapid variations of the observables rather than on those driving the long-term changes in carbon and water budgets (e.g., tree growth and soil carbon turnover). The prior values $x_b$ of the parameters correspond to the standard values of ORCHIDEE. Most of the parameters (apart from $Q_{10}$, $K_{soilC}$, and $K_{albedo,veg}$) are PFT dependent, which in the case of Fontainebleau leads us to consider two sets of parameters for forest (TempDBF) and grasses (C3grass), while only the parameters for TempEBF are accounted for at Puechabon. For some parameters, we have chosen to keep relatively large intervals of variation (though still within a physical range) to allow more flexibility to the optimization in order to focus on the compatibility and complementarity between flux and FAPAR data.

### 2.3.4. Assimilation Scenarios

In the following, we test different assimilation scenarios in order to understand the respective constraint brought by each source of data and to investigate their compatibility. For each scenario, a specific set of parameters is considered for optimization. Four different subsets are used:

**P1:** only the parameters of ORCHIDEE that control photosynthesis (and hence gross primary productivity—GPP) and phenology (all parameters in Table 1);

**P2:** all phenological parameters plus a few parameters related to photosynthesis;

**P3:** all phenological parameters only; and

**P4:** all phenological parameters minus $LAI_{init}$ and $LAI_{MAX}$.

In total, five different assimilation scenarios are evaluated:

**A1:** assimilation of flux data alone, with the set of parameters P1;

**A2:** assimilation of original FAPAR products alone, with the sets of parameters P2;

**A3:** assimilation of original FAPAR products alone, with the sets of parameters P3;

**A4:** assimilation of normalized FAPAR products, with the set of parameters P3 minus $LAI_{init}$ and $LAI_{MAX}$ (P4); and

**A5:** assimilation of both flux and FAPAR data (FAPAR being normalized for the deciduous site), with the set of parameters P1.

Scenarios A1 to A3 are designed to evaluate the compatibility of the two data streams together and their compatibility with ORCHIDEE, by analyzing the impact of the assimilation of one data stream (flux or FAPAR) on the model-data mismatch with the other data stream. Comparison between scenarios A2, A3, and A4 will provide insight on the amount of information that can be inferred from FAPAR data to improve the ORCHIDEE model. Finally, scenario A5 will allow us to assess the potential benefit of a joint assimilation of flux and FAPAR data (in terms of model improvement).

For assimilation studies with FAPAR data only (scenarios A2, A3, and A4), we use only a subset of the parameters in the optimization, as the information content of FAPAR products is more limited and allows constraining only the (i) intensity (related to GPP) and (ii) timing of the vegetation cycle through a smaller number of control parameters (parameter sets P2, P3, and P4). The assimilation test with normalized FAPAR data (scenario A4) is designed to illustrate that FAPAR constrains mainly the phenology and the seasonality of the growing season. The scenario using normalized FAPAR (see section 2.1.2) is tested only for the deciduous site (Fontainebleau) that shows a strong seasonality. Whenever the observed FAPAR are normalized, the model simulations are processed using the same procedure.

## 3. Results

### 3.1. Assimilation of In Situ Flux Measurements

#### 3.1.1. Prior Fit to Observations

Figure 1 compares the agreement between the NEE and LE measurements with the prior ORCHIDEE simulations and with the optimized model after assimilation of the flux data for Fontainebleau and Puechabon
<figure><img src="image_3.png" /><figcaption>**Figure 1.** For (a) Fontainebleau and (b) Puechabon, comparison of NEE and LE observations (black) with corresponding ORCHIDEE model simulations before (blue line) and after assimilation (red line) for scenario A1. The 1 sigma uncertainty on the flux observations used for the computation of the misfit function (equation (2)) is shown in grey. The RMSE measures the fit of the model prior and posterior simulations with the corresponding observations. Simulated FAPAR temporal variations before and after assimilation are also shown in the right column and are compared to the in situ data.</figcaption></figure>

(scenario A1). For Fontainebleau, the prior ORCHIDEE simulation shows a relatively good agreement with the observations: the seasonal patterns are well captured for both fluxes, as well as the synoptic (10 to 20 day) variations. The prior model, however, overestimates the ecosystem respiration in winter by about 1 g C/m²/d for both sites. This feature follows from an incorrect initialization (overestimation) of the carbon pools after the spin-up runs, where the soil carbon pool may represent the maximum soil carbon content for such forest, without considering past land use and disturbances. At the same time, the model simulates a lower amount of carbon uptake during the growing season than that measured for this relatively young forest. The model-data agreement is lower for Puechabon, in particular for the LE, which ORCHIDEE strongly overestimates in summer (by about 400%) when the local measurements show a strong hydric stress.

### 3.1.2. Improvement After Optimization

The improvement of the model-data fit resulting from the assimilation of the flux data is particularly evident for Fontainebleau, with a reduction of RMSE from 1.3 to 0.57 g C/m²/d for NEE and from 6.44 to 4.85 W/m² for LE. The excess of respiration in winter for the NEE seasonal cycle has been corrected, mainly as a result of the tuning of the initial carbon pools as shown by the decrease of the K<sub>soilC</sub> parameter (see section 3.4.1). The amplitude of the synoptic events is also well captured.

The improvement of the model-data fit for NEE and LE is smaller for Puechabon in summer, with an inappropriate model response to summer heat and dryness. The observations show a strong decrease in evaporation around day 200 (end of June) that the model is unable to represent after optimization (see discussion in section 4.1).

In contrast, the assimilation of the in situ flux data have a small impact on the simulated FAPAR time series; the optimized simulations remain very close to the prior values, with nearly no change in RMSE with respect to the in situ FAPAR data (see values in Figure 1). A small change in the simulated phenology can, however, be observed for Fontainebleau, with a later start of the growing season (~10 days), as explained by the increase of K<sub>pheno_crit</sub> (see section 3.4.1). However, this results in an even stronger bias in the timing of the start of the growing season depicted by FAPAR measurements compared to the prior model.

## 3.2. Assimilation of FAPAR Products

### 3.2.1. Compatibility of the A Priori Model With the FAPAR Products

Figure 2 shows strong discrepancies in the magnitude of FAPAR between the various products and the prior ORCHIDEE simulation at both sites. The satellite products are systematically lower than the model and than
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Comparison of the prior ORCHIDEE FAPAR simulation for the Fontainebleau and Puechabon tower pixels, with the various FAPAR products used in this study. For Fontainebleau, the vertical lines indicate the actual dates of the SPOT observations (the two last acquisitions are separated by 1 day only).</figcaption></figure>

the in situ measurements, in particular for Puechabon where the magnitude of the MERIS FAPAR is up to half that of the in situ data. This discrepancy might be due to errors resulting from the satellite FAPAR processing chain (radiative transfer modeling, including incorrect representation of the optical properties of Mediterranean oak species; atmospheric correction) and also partly from differences in spatial footprint (given the heterogeneity of the landscape around the tower). The large bias seen in Puechabon has been reported in other studies for evergreen broadleaf forests between other satellite products [Weiss et al., 2007]. Contrary to this, the magnitudes of the in situ FAPAR products are very close to those of the prior model simulations.

Differences in seasonality are also noticeable. For the deciduous forest site of Fontainebleau, the satellite products present smooth temporal variations, whereas both the simulated FAPAR time series and the in situ measurements show more abrupt changes at leaf onset and senescence. In addition, satellite-derived and in situ FAPAR indicate an earlier start of the growing season for Fontainebleau compared to the ORCHIDEE simulation; however, the agreement is stronger between the model and the in situ FAPAR measurements than between any of the FAPAR products (model or observations).

<figure><img src="image_5.png" /><figcaption>**Figure 3.** For (a) Fontainebleau and (b) Puechabon, results of the assimilations of FAPAR data (scenarios A2 and A3) with respect to the model-data agreement for NEE, LE, and the various FAPAR products considered. The observed data are shown in black and their 1 sigma uncertainty in grey.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Ratio of the RMSE of the fit after optimization to the prior RMSE between model and data (FAPAR, NEE, and LE), for the data assimilations conducted with original and normalized FAPAR time series (scenarios A3-P3 and A4-P4), considering the in situ, MERIS, and SPOT products. Values < 1 (respectively, > 1) correspond to an improvement (degradation) of the model-data agreement.</figcaption></figure>

products, as measured by the decrease of the RMSEs (around 40% in most cases). However, some features can still not be represented after optimization. For Fontainebleau, the model underestimates the measured in situ FAPAR values in winter. ORCHIDEE also remains unable to reproduce the very high in situ FAPAR values in the summer measured at the two sites, the change in magnitude being capped by the imposed maximum vegetation fraction (90% at Puechabon). The discrepancy may also arise due to the assumptions used to simulate FAPAR from LAI in the model (discussed in section 4.2). For Puechabon, the simulated FAPAR values at the beginning of the year are adjusted to the observations (either in situ or satellite). However, the model cannot reproduce the low satellite FAPAR values for the rest of the year, which is possibly due to a too high fraction of temperate evergreen broadleaved trees in the model for the satellite footprint. For the deciduous forest, assimilating in situ FAPAR leads to an earlier onset (opposite to the results obtained with the flux data, Figure 1) for the two assimilation scenarios considered (scenarios A2 and A3); the shift in senescence, however, depends on the set of parameters optimized: it is delayed with the extended set of parameters P2, where it is advanced using only the phenological parameters P3 (in agreement with the assimilation of flux data).

Concerning the agreement of the model to the flux data (NEE and LE), assimilating FAPAR produces contrasting results depending on the product used (in situ or satellite), the fluxes (NEE or LE), the site considered, and the set of parameters optimized. The change in the simulated fluxes is mainly proportional to the difference in magnitude between the prior ORCHIDEE FAPAR simulations and the observations. For the LE, the assimilation of satellite-derived FAPAR improves the fit to the flux data at both sites, while using in situ FAPAR degrades it. This is because the necessary decrease of leaf area (section 3.4.2) for fitting SPOT or MERIS FAPAR observations induces a decrease in the modeled LE. The impact on the NEE may be more or less detrimental depending on the FAPAR data set used. The use of in situ FAPAR does not result in much change of NEE at either site. For Fontainebleau, the use of additional model parameters related to photosynthesis (scenario A2 with set of parameters P2) leads to an overestimation of the carbon uptake in summer and hence increases the model-data mismatch. Only marginal changes are observed in winter after the assimilation of in situ FAPAR because the NEE is principally driven by soil and litter respiration not impacted by FAPAR data. There is only a small improvement for Puechabon early in the season. Fitting satellite FAPAR results in an even stronger negative impact on the simulated NEE than in situ data. Especially, the assimilation of MERIS FAPAR with low values during summer time largely degrades the fit to the NEE. The apparent inconsistency between FAPAR products in summer and the magnitude of NEE measured on the site will be further discussed in sections 4.2 and 4.3.

### 3.2.3. Impact of Normalizing the FAPAR Data on the Assimilation

The previous results have highlighted that correcting the magnitude of the simulated FAPAR to fit that of the observations may drastically degrade the fit to the NEE. We therefore tested a data assimilation scenario where we normalize FAPAR (scenario A4) to constrain the timing of the growing season (phenology) only. Note that in this set up, only five phenological parameters that control the seasonal cycle of the vegetation

Overall, the prior comparison indicates strong discrepancies in mean absolute values of the different FAPAR estimates, while their temporal variations are consistent.

### 3.2.2. Improvement After Optimization

Figure 3 shows the impact of assimilating the different FAPAR data sets on the simulated fluxes and FAPAR for scenarios A3 (phenology parameters, P3) and A2 (extended parameters, P2). In all cases, apart from scenario A2 applied only to in situ FAPAR measurements, the seven phenological parameters of Table 1 have been optimized (section 2.3.3).

The assimilation significantly corrects the discrepancies in the magnitude between the prior model and the various FAPAR
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Ratio between the posterior and prior RMSE of fit, between the model simulations and different observed variables, considering assimilations performed with (first column) only flux data (scenario A1), (left bars in the second to fourth columns) FAPAR data only (normalized data for Fontainebleau scenario A4-P4, original data for Puechabon scenario A3-P3) and (right bars in the second to fourth columns) the combination of the two data streams (scenario A5).</figcaption></figure>

are optimized (see Table 1: phenology parameters P3 minus LAI<sub>init</sub> and LAI<sub>MAX</sub>). The assimilations are conducted for Fontainebleau only as the evergreen site of Puechabon exhibits a too low seasonality.

Figure 4 illustrates how the normalization of FAPAR impacts the assimilation, compared to the previous results with unnormalized (original) data. It shows the ratio of the RMSE between model and data after and prior to the assimilation. Values less than (greater than) 1.0 show an improvement (degradation) in the model with respect to the data. As expected, the improvement of the model-data agreement with respect to original FAPAR is lower with the normalization, as the correction of the magnitude of the modeled FAPAR is not sought. Nevertheless, the normalization still improves the modeled FAPAR (with the exception of MERIS products), while it mostly reduces the degradation of the model-data fit for NEE. However, the results still tend to indicate that it may not be possible to improve the modeled NEE as compared to the prior simulations when assimilating FAPAR alone (either using original or normalized data). For LE data, the results are more variable between the various FAPAR products as different errors may compensate each other (with in situ FAPAR, the normalization improves the fit to LE, while it slightly degrades it for the satellite FAPAR).

### 3.3. Joint Assimilation of In Situ Flux Measurements and FAPAR Products

Figure 5 synthesizes the improvement/degradation in model-data fit (posterior to prior RMSE ratio) for several observations (NEE, LE, and FAPAR) when considering successively the assimilation of in situ flux data alone (scenario A1), FAPAR products alone (scenarios A3 and A4), or combining both flux and FAPAR data (scenario A5). In order to make the best possible use of the FAPAR products, we have chosen to normalize the time series for Fontainebleau (scenario A4, see section 4.2) and assimilate the original data for Puechabon (scenario A3) given the small seasonal variations. Note that we also use the GPP as a diagnostic (see section 2.1.1).

The joint assimilation of NEE and LE in situ flux measurements and FAPAR products (scenario A5) seen in Figure 5 reconciles the two sources of information and the model, thus dealing with the inconsistencies described above when only one data stream is assimilated. The optimized simulations improve the fit to both the flux and FAPAR data compared to the a priori model simulation at both sites. It results in a similar model-data agreement as that obtained when each data stream is assimilated independently. The assimilation of both data streams together prevents the degradation seen for the variable not included in the individual data
<figure><img src="image_8.png" /><figcaption>**Figure 6.** Values of some optimized parameters related to photosynthesis and respiration (see Table 1) for the (a) TempDBF (Fontainebleau) and (b) TempEBF (Puechabon) PFTs and for the data main assimilation scenarios considered (scenarios A1, A2, and A5). When assimilating FAPAR alone, original data are used; when assimilating jointly flux and FAPAR data, normalized FAPAR are used for Fontainebleau while original products are considered for Puechabon. For each parameter, the box corresponds to its range of variation and the horizontal line to its prior value; the posterior uncertainty is provided for each case.</figcaption></figure>

stream optimizations (scenarios A1 to A4). This is particularly noticeable for the assimilation of satellite FAPAR products, considering that in five out of the seven cases assimilating FAPAR data alone causes an increase in the RMSE for NEE data of up to 80% (Figure 3).

For NEE, there is an improvement of ≥50% or more at Fontainebleau and Puechabon when flux and FAPAR products (normalized at Fontainebleau) are assimilated together, which is about the same as that when just fluxes are used in the optimization. At Fontainebleau, the highest improvements are obtained with the satellite products, whereas at Puechabon the use of in situ FAPAR data (unnormalized) leads to the best results. For LE, the 10–35% reduction of the prior RMSE obtained at both sites when both data streams are included is lower than when fluxes only are assimilated, but it corresponds always to an improvement with respect to the assimilation of only FAPAR. For FAPAR, the reduction in RMSE obtained by the joint assimilation is even higher than when FAPAR data alone are used in the assimilation (around 40%), which is partly explained by the optimization of additional parameters related to photosynthesis in the joint assimilation. For GPP, the degradation of the model-data fit following the assimilation of FAPAR data alone is similar or higher to that observed for NEE. In contrast, the joint assimilation generally leads to an improvement of the model-data agreement for GPP, just as for NEE and LE flux data. For Fontainebleau, the GPP improvement is, however, lower than for NEE, which reveals that there are compensating errors between GPP and ecosystem respiration.

<figure><img src="image_9.png" /><figcaption>**Figure 7.** Same as Figure 6 for some optimized parameters related to phenology (see Table 1). The results of the estimation of only the phenological parameters with FAPAR data are also presented for Puechabon (scenario A3) and Fontainebleau with normalized data (scenario A4).</figcaption></figure>
### 3.4. Optimized Parameter Values and Uncertainties

For the main data assimilation scenarios described (scenarios A1 to A5) above, the optimized values and the associated uncertainty are presented for both sites for a selection of parameters controlling photosynthesis and respiration in Figure 6 and phenology in Figure 7. For Fontainebleau, only the parameters for the TempDBF PFT are presented, although the grass PFT parameters were also optimized.

#### 3.4.1. Assimilation of Flux Data: Impact on Parameters

For both sites, we note significant changes of the parameters that control the ecosystem respiration: $K_{soilC}$ (multiplier of soil initial C pools) which is decreased and $Q_{10}$ (dependency of heterotrophic respiration to soil temperature) which is increased. The changes lead to a decrease of the heterotrophic respiration at both sites, which is the principal cause of the model-data improvement in winter. The optimization of these two parameters explains more than 60% of the decrease of the RMSE for NEE for the two sites (results not shown). This points out the importance of the model spin-up that largely determines the amplitude of ecosystem respiration and hence the fit to NEE during winter [Carvalhais et al., 2010; Kuppel et al., 2012]. The errors from these two respiration parameters are highly correlated (as depicted by the a posteriori error correlation inferred from equation (4)) for Fontainebleau (0.78) but merely for Puechabon (0.25).

At Fontainebleau, the increase of $V_{cmax}$ and SLA mainly explains the increase of the carbon uptake in summer, while $L_{AI_{MAX}}$ is reduced. The error correlations between the estimated parameters are usually lower than 0.35 except between few photosynthesis parameters and/or respiration parameters: $V_{cmax}$ with $G_{s,slope}$ (-0.72), $V_{cmax}$ with $T_{opt}$ (0.51), and $V_{cmax}$ with $GR_{frac}$ (0.59), or between the phenological parameters $L_{agecrit}$ and $\tau_{leaf,init}$ (-0.54).

At Puechabon, the decrease of the modal data mismatch for LE is achieved by decreasing $V_{cmax}$, $G_{s,slope}$, and $H_{um_{cste}}$ principally. Like for Fontainebleau, the estimate for $V_{cmax}$ is strongly correlated with the estimate for $G_{s,slope}$ (-0.73), $T_{opt}$ (0.44), and $GR_{frac}$ (0.55), and also with SLA (-0.71) and $K_{soilC}$ (0.41).

The assimilation of flux data results in a shortening of the growing season for the deciduous forest at Fontainebleau, even though this is not clearly visible in the simulated FAPAR time series in Figure 3. This is mainly achieved by a decrease of $L_{age,crit}$ of about 30 days (corresponding to an earlier loss of the leaf photosynthetic efficiency) combined with an increase of $K_{pheno,crit}$ (later leaf onset). At Puechabon, the evergreen phenology model is nearly only sensitive to the leaf age ($L_{age,crit}$).

Overall, the absolute values of the error correlations between parameters are small with an average (±standard deviation) of 0.05 ± 0.10 for Fontainebleau and 0.09 ± 0.15 for Puechabon.

#### 3.4.2. Assimilation of FAPAR Data: Impact on Parameters

In spite of their very different temporal profiles, in situ FAPAR data and MERIS products induce an earlier senescence at Fontainebleau, mostly explained by the reduction of $L_{age,crit}$. An increase of that parameter is obtained with SPOT products without any noticeable change on the entering into senescence (Figure 3). The fit to both in situ and SPOT FAPAR data lead to an earlier start of the growing season with the decrease of $K_{pheno,crit}$ (conversely to the optimization with flux data). Note that the uncertainty attached to $K_{pheno,crit}$ is smaller when assimilating FAPAR data (and even smaller with the in situ ones) than when assimilating flux data. The large changes in $L_{AI_{MAX}}$ (increase for in situ FAPAR and decrease for the satellite FAPAR) are the main cause of the degradation of the modeled NEE.

For Puechabon, the optimization of the mean FAPAR magnitude to match in situ or MERIS data is achieved by changing the initial LAI ($L_{AI_{init}}$) parameter (not shown in the figures). When assimilating in situ data, it increases from 4.7 to reach the maximum value allowed of 8. Conversely, when assimilating MERIS products, it decreases to 1.43.

#### 3.4.3. Joint Assimilation of Flux and FAPAR Data: Impact on Parameters

Figures 6 and 7 reveal that the values of the retrieved parameters may change significantly when flux data are assimilated alone (scenario A1) or combined with FAPAR data (scenario A5). It depends on the initial bias between the model and each FAPAR product and the correlations between model parameters. This is the case for some photosynthesis ($V_{cmax}$ at both sites; $H_{um_{cste}}$ or SLA at Puechabon) or phenology ($K_{LAI,happy}$ and $K_{pheno,crit}$ at Fontainebleau; $L_{age,crit}$ at both sites) parameters, which show high variability in posterior values, depending on the assimilation scenario. For most of the parameters, however, the optimized values are closer to the ones retrieved when assimilating flux data alone, indicating the stronger weight of NEE and LE in the
<table><caption>**Table 2.** Root-Mean-Square Errors of Fit Between Daily Measured NEE and LE Fluxes at Fontainebleau and Puechabon and (1) the A Priori Model and (2) the Optimized Models Using Either Flux Data Alone (Scenario A1) or Flux and FAPAR Data (Scenario A5), for Years Outside the Calibration Period</caption><thead><tr><th rowspan="2"></th><th colspan="2">Prior</th><th colspan="2">Flux</th><th colspan="2">Flux + FAPAR (In Situ)</th><th colspan="2">Flux + FAPAR (SPOT)</th><th colspan="2">Flux + FAPAR (MERIS)</th></tr><tr><th>NEE</th><th>LE</th><th>NEE</th><th>LE</th><th>NEE</th><th>LE</th><th>NEE</th><th>LE</th><th>NEE</th><th>LE</th></tr></thead><tbody><tr><td colspan="11" style="text-align:center; font-style:italic;">Fontainebleau</td></tr><tr><td>2006<sup>a</sup></td><td>1.30</td><td>6.43</td><td>0.66</td><td>5.83</td><td>0.94</td><td>6.38</td><td>0.74</td><td>6.02</td><td>0.61</td><td>6.06</td></tr><tr><td>2007</td><td>1.85</td><td>11.43</td><td>1.37</td><td>15.26</td><td>1.05</td><td>14.49</td><td>1.30</td><td>14.86</td><td>1.35</td><td>15.39</td></tr><tr><td>2008</td><td>1.64</td><td>9.13</td><td>1.06</td><td>13.46</td><td>1.01</td><td>12.52</td><td>1.03</td><td>12.58</td><td>1.06</td><td>13.42</td></tr><tr><td>2009</td><td>1.14</td><td>13.23</td><td>0.84</td><td>17.25</td><td>0.77</td><td>17.14</td><td>0.75</td><td>17.01</td><td>0.76</td><td>17.57</td></tr><tr><td>2010</td><td>1.50</td><td>10.51</td><td>0.75</td><td>14.84</td><td>1.02</td><td>14.46</td><td>0.80</td><td>14.39</td><td>0.80</td><td>15.19</td></tr><tr><td colspan="11" style="text-align:center; font-style:italic;">Puechabon</td></tr><tr><td>2000</td><td>0.65</td><td>30.70</td><td>0.44</td><td>18.14</td><td>0.49</td><td>21.25</td><td></td><td></td><td>0.63</td><td>30.40</td></tr><tr><td>2001</td><td>1.11</td><td>28.26</td><td>0.60</td><td>18.62</td><td>0.66</td><td>20.93</td><td></td><td></td><td>1.08</td><td>27.65</td></tr><tr><td>2002</td><td>1.09</td><td>23.81</td><td>0.56</td><td>17.45</td><td>0.60</td><td>18.09</td><td></td><td></td><td>1.05</td><td>23.22</td></tr><tr><td>2003</td><td>1.06</td><td>29.99</td><td>0.63</td><td>20.35</td><td>0.76</td><td>20.88</td><td></td><td></td><td>1.02</td><td>29.40</td></tr><tr><td>2004<sup>a</sup></td><td>1.01</td><td>24.82</td><td>0.53</td><td>14.80</td><td>0.66</td><td>16.22</td><td></td><td></td><td>1.03</td><td>23.91</td></tr><tr><td>2005</td><td>0.94</td><td>32.30</td><td>0.59</td><td>20.91</td><td>0.55</td><td>26.64</td><td></td><td></td><td>0.89</td><td>31.41</td></tr><tr><td>2006</td><td>0.96</td><td>31.01</td><td>0.56</td><td>23.81</td><td>0.61</td><td>24.98</td><td></td><td></td><td>0.93</td><td>31.06</td></tr></tbody></table>

assimilation process (this is particularly the case for the parameters related to respiration—$Q_{10}$ and $K_{soilc}$). The estimated parameter uncertainties are usually lower for the joint assimilation as expected from the expression of A (equation (4)) using an increased number of observations that are treated as independent data.

The assimilation of FAPAR products tends to diminish the correlation error between the phenological parameters as compared to scenario A1 (relying only of flux data) except between $LAI_{MAX}$ and $K_{LAI,happy}$ (values between -0.98 with in situ FAPAR and -0.58 with MERIS at Fontainebleau). In general, however, adding FAPAR in the assimilation procedure slightly increases the correlation errors between the other parameters. For Fontainebleau (Puechabon), the mean and standard deviation of the absolute values of the error correlations go from 0.05 ± 0.10 (0.09 ± 0.15) when assimilating fluxes only to 0.06 ± 0.13 (0.12 ± 0.20) when adding in situ, SPOT, or MERIS FAPAR data.

### 3.5. Temporal Validation of the Optimized Model

Using in situ flux measurements available at the two sites for additional years, we evaluate (i) if the model improvement obtained over the calibration year still holds for other years (different meteorological forcing) and (ii) how it may differ depending on the data streams used for calibrating the model. The analysis is performed for scenarios A1 (flux only) and A5 (flux + FAPAR, with normalized FAPAR products at Fontainebleau). Table 2 displays the RMSE for the calibration and validation years.

For both sites, the model-data fit is improved for the validation years for all scenarios. At Fontainebleau, the agreement of the optimized model to the flux data for the validation years is usually reinforced by adding FAPAR data in the assimilation procedure. The RMSEs are lower with flux + FAPAR than with flux only between 2007 and 2010 when using SPOT or in situ FAPAR (for the latter, however, the benefit in 2010 is obtain only for LE). A similar finding is also obtained with MERIS FAPAR products but in a lower extent, the gain for NEE being achieved from 2006 to 2009 (compared to the assimilation of fluxes only), and only for 2008 for LE. The benefit of using FAPAR is less pronounced at Puechabon although there is always an improvement of the optimized model with respect to the standard ORCHIDEE parameterization. Finally, for the two sites, the higher model-data agreement is usually achieved with the in situ FAPAR measurements rather than with the satellite FAPAR products.

## 4. Discussion

### 4.1. Constraints and Insights Brought by Each Data Stream

The assimilation of each individual data stream (either flux or FAPAR data) leads to a significant improvement of the model-data agreement with respect to the variables that are assimilated while also possibly highlighting some model deficiencies.

<footnote><sup>a</sup>Years used for the model calibration.</footnote>
For the assimilation of flux data, the strong reduction of the model-data misfit with respect to NEE and LE confirms previous studies that show the benefit of assimilating local flux measurements into ORCHIDEE for deciduous ecosystems [Kuppel et al., 2012; Santaren et al., 2014]. The poorer performance obtained for the Mediterranean forest points out some model structural deficiencies that require not only a change in the model parameters but also in the model structure to better simulate the behavior of semiarid ecosystems, as for instance using a refined soil hydrology scheme to better calculate the dynamic of the soil water stress and incorporating the effect of drought stress on the photosynthetic capacity not related to stomatal closure [Keenan et al., 2009]. The assimilation compensates for the bad representation of soil hydrology by overfitting other parameters, with therefore suboptimal values. Also, at Puechabon, both FAPAR products show bowl shaped temporal variations, with higher values in winter. This could reflect leaf adaptation to dry summer conditions and a reduction of the photosynthetic activity, or simply a change of the extinction coefficient with the illumination angle (with higher optical paths for oblique illuminations), neither of which are represented in ORCHIDEE.

The large discrepancies between the various FAPAR products lead to different optimized model trajectories after assimilation, with most of the time a large negative impact on the NEE and LE fluxes when FAPAR data are assimilated alone.

## 4.2. Limitations of Using FAPAR Data

Although remotely sensed products of vegetation activity have been widely used in semidiagnostic ecosystem models [McCallum et al., 2009; Seixas et al., 2009; Jung et al., 2011], our study highlights some challenges to be overcome before they can be fully exploited to optimize process-based model parameters. In particular, the accuracy of the products should be improved, especially the mean value during the peak of the growing season.

In addition, our ability to exploit any observation for monitoring ecosystem status depends on the ability of the observation operator we build to reproduce and account for the various processes and vegetation components contributing to the observables. FAPAR estimates are impacted by nongreen components of the canopy (as for instance litter, tree trunks, and branches) when simulated FAPAR usually relies on green LAI (see equation (1)). Also, we used a constant value for the extinction coefficient (0.5), whereas in reality it varies in time with the sun angle, as well as with the angle distribution and clumping of leaves that modify the canopy capacity to intercept the incoming solar radiation.

Note that we could foresee assimilating satellite-derived LAI estimates rather than FAPAR as using LAI does not necessitate the implementation of an additional observation operator (with its own assumptions and error characteristics), and one could then consider that LAI brings a more direct constraint on the land surface model as it directly relates to model outputs (i.e., leaf C mass). Nevertheless, the discrepancies between different satellite-derived LAI estimates are also considerable [Garrigues et al., 2008]; the associated uncertainty is usually greater than that of FAPAR [Claverie et al., 2013] as the reflectance saturates with large LAI values. Also, the nonlinear relationship between reflectance and LAI [Liang, 2000] results in a higher increase of LAI retrieval errors with the pixel heterogeneity than for FAPAR [Weiss and Baret, 1999]. Note finally that the differences between the FAPAR products presented here are of the same order than in other studies and are mainly attributed the differences in the radiative transfer models used in the processing chain [Meroni et al., 2013; D'Odorico et al., 2014]. Given the current level of uncertainty in the magnitude of these remote sensing products, FAPAR or LAI should be used to constrain only the timing of the seasonal cycle of the vegetation (phenology) of LSMs when assimilated alone.

For forests, the understory below the canopy and the vertical mixing of vegetation are still not represented in nearly all models (including ORCHIDEE), whereas this plays a role for the dynamics of the in situ and satellite FAPAR observations. Likewise, the model is not able to account for the spatial variability in the phenological stages that occur in a real landscape. Because of different environmental conditions (different species for a given PFT; different soil characteristics impacting water and nutrient availability; and microclimate), the dates of leaf onset and senescence of the various plant individuals are distributed in time, giving rise to a smoother temporal profile (as seen in the satellite products) than the sudden changes simulated by most land surface models. Then, rather than optimizing discrete phenological parameters when using satellite observations, a solution could be to account for a distribution of parameter values in the phenology equations, as proposed by Knorr et al. [2010].
The footprint of high spatial resolution satellite observations (HR) is more compatible with that of the flux tower. Hence, HR data are likely more easily combined with in situ flux measurements than medium spatial resolution (MR) data. HR satellite products have, however, their own limitations (low revisit frequency, lower signal to noise ratio, etc.) that make MR observations still valuable to evaluate/optimize ecosystem models [Knorr et al., 2010; Kato et al., 2013]. Also, while a few vegetation types (one or two typically) suffice to characterize HR pixels, an increased number of different PFTs are required at lower resolution, depending on the scene heterogeneity. Hence, the optimization with MR data is likely more complex with a larger set of parameters to optimize simultaneously.

### 4.3. Compatibility of Flux With FAPAR and With ORCHIDEE

The impact on the model-data agreement with respect to the other data stream not included in the assimilation is a first assessment of the consistency between the different observations and the ORCHIDEE model. The assimilation of flux data had almost no impact on the corresponding simulated FAPAR, which indicates that there are enough degrees of freedom with the optimization of photosynthetic rate and respiration parameters to fit the NEE without drastically changing the maximum LAI and the other phenological parameters and, thus, FAPAR.

In contrast, the assimilation of original (unnormalized) FAPAR time series often exerts a strong negative impact on the modeled NEE and LE. The increase of the model misfit with respect to the in situ flux measurements is proportional to the differences between the a priori ORCHIDEE FAPAR values and the in situ/satellite products. In particular, during the growing season, the fit of the simulated FAPAR to the observed magnitude of the satellite products has a drastic influence on the rate of carbon uptake, due to the decrease in LAI that reduces the carbon fixation capacity of the ecosystems (Figure 3). Using in situ FAPAR still leads to a negative impact on the fluxes (due to an increase of LAI), although much lower than with the satellite products. This reveals that the optimization with in situ FAPAR as a sole constraint may lead to a set of parameters that is not optimal to simulate carbon fluxes.

Note also that the better compatibility between the model and in situ FAPAR data relates to the shape of the time series. In particular for the deciduous forest, the increase rate of FAPAR after leaf onset and senescence are quite abrupt in the in situ observations typical of the behavior of an even-age forest stand, as represented by ORCHIDEE. This feature is not captured by the satellite products that have smoother temporal variations, the origin of which is twofold: FAPAR processing (smoothing using a moving window) and spatial heterogeneity.

Opposite changes in the start of leaf onset are obtained whether one or the other data stream is used: while a later leaf unfolding is estimated when flux data are assimilated, an earlier one is retrieved with all FAPAR products but MERIS. This may result from the fact that measured FAPAR and NEE have different drivers: the temporality of NEE is governed by the combination of the temporal profiles of carbon respiration and carbon uptake and is therefore different to that of GPP, which is likely more correlated to FAPAR and LAI [Connolly et al., 2009]. This suggests that both data streams bring different information on the ecosystem functioning and thus may be complementary to optimize model parameters (provided there are not significant biases associated with each observation and incompatibilities between the model and the observations).

### 4.4. Normalization of FAPAR Data

Given the issues described above, this study proposes a better use of FAPAR products by only constraining the timing of the growing season through the optimization of a limited set of phenological parameters (consistent with Knorr et al. [2010]). For deciduous canopies with pronounced seasonal cycle, we found that the assimilation conducted with normalized FAPAR time series still improves the model-data agreement with respect to FAPAR and mostly reduces the degradation of the model-data fit for NEE and LE. In addition, the normalization of the FAPAR time series decreases the differences between the different satellite products and reduces the strong uncertainties attached to the magnitude of the FAPAR data.

### 4.5. Joint Assimilation of Fluxes and FAPAR

The joint assimilation of flux measurements with any of the FAPAR products (considering either original or normalized FAPAR values) succeeded in reconciling the two sources of information by finding a set of
parameters that improved the model-data fit with respect to each data stream. Note that we performed an additional test at Fontainebleau using original FAPAR data in the joint assimilation which resulted in an improvement similar to that obtained when using normalized FAPAR values (results not shown). This indicates that different sets of parameters, within physical ranges, are able to provide a reasonable fit to the observed NEE and to different LAI values (constrained through FAPAR). Optimizing the model with only one data stream may thus potentially lead to non-optimal parameters.

The fit resulting from the joint assimilation (in terms of RMSE) is a compromise between the cases where the two data streams are considered independently and always results in an improvement of the model. Satellite FAPAR products helped to better characterize some critical parameters controlling the phenology of ORCHIDEE (turnover of the leaves, time of leaf onset), parameters that are less constrained by NEE/LE data alone. The differences in the estimated parameter values (in particular the phenology parameters) when the FAPAR and flux data are combined in the assimilation, as compared to when they are assimilated independently, reveals some levels of complementarity between the data streams. This result highlights the importance of combining several data streams. It may strongly impact the model predictions as different parameter sets will lead to different carbon balance trajectories in the future. Given the correlations between the model parameters, a joint assimilation seems preferable to a stepwise approach which may penalize the retrieval of parameters unless these correlations are properly accounted for.

The benefit of combining FAPAR products with flux measurements was confirmed by the temporal validation experiments where the optimized model was evaluated with respect to ancillary flux data outside the calibration period. The additional use of FAPAR helps better constraining the vegetation phenology and hence better partitioning the carbon (heterotrophic respiration and photosynthesis) and energy fluxes (sensible and latent heat).

## 5. Conclusion

In the context of increasing availability of in situ flux data and satellite observations of vegetation activity, we have investigated the benefit of combining the two data streams in a data assimilation framework for optimizing the parameters of the ORCHIDEE process-based vegetation model. In doing so, we have highlighted the main difficulties that need to be worked out before being able to fully exploit such remote sensing products.

The upcoming new generation of satellite observing systems, including the ESA Sentinel-2 mission, will monitor the Earth's surface at both spatial and temporal high resolutions and should therefore permit significant improvements of our knowledge of the vegetation dynamics. However, the uncertainty of derived products such as FAPAR or LAI may limit the potential benefits of future satellite observations. Reducing the observation uncertainty becomes even more crucial for evergreen vegetation with low seasonality as the magnitude of the satellite products becomes the only information that can be used in an assimilation framework. We need to prepare advanced land surface data assimilation systems capable of assimilating these new satellite data, in combination with in situ local measurements, and other data streams (such as atmospheric CO2 concentration data, biomass data, and surface temperature). The complementarity of FAPAR and flux data demonstrated in this study raises hopes that current process-based models are capable of the task, but we must see whether the result holds for more ecosystems and a wider array of data.

## References

* Aubinet, M., et al. (2000), Estimates of the annual net carbon and water exchange of forests: The EUROFLUX methodology, *Adv. Ecol. Res.*, **30**, 113–175.

* Bacour, C., F. Baret, D. Béal, M. Weiss, and K. Pavageau (2006), Neural network estimation of LAI, fAPAR, fCover and LAIxCab, from top of canopy MERIS reflectance data: Principles and validation, *Remote Sens. Environ.*, **105**(4), 313–325.

* Baldocchi, D. (2008), 'Breathing' of the terrestrial biosphere: Lessons learned from a global network of carbon dioxide flux measurement systems, *Aust. J. Bot.*, **56**, 1–26, doi:10.1071/BT07151.

* Baret, F., et al. (2007), LAI, fAPAR and fCover CYCLOPES global products derived from VEGETATION—Part 1: Principles of the algorithm, *Remote Sens. Environ.*, **110**(3), 275–286, doi:10.1016/j.rse.2007.02.018.

* Botta, A., N. Viovy, P. Ciais, P. Friedlingstein, and P. Monfray (2000), A global prognostic scheme of leaf onset using satellite data, *Global Change Biol.*, **6**, 709–725.

* Braswell, B. H., W. J. Sacks, E. Linder, and D. S. Schimel (2005), Estimating diurnal to annual ecosystem parameters by synthesis of a carbon flux model with eddy covariance net ecosystem exchange observations, *Global Change Biol.*, **11**, 335–355, doi:10.1111/j.1365-2486.2005.00897.x.

## Acknowledgments

The authors are very grateful to Serge Rambal and Jean Marc Ourcival for providing the in situ flux and FAPAR measurements at Puechabon and to Jean-Yves Pontallier and Guerric le Maire for the data at Fontainebleau. The satellite data presented in this article are available from the authors on request. In situ measurements can be obtained upon request to the persons listed above. This study was partly financed by the European Space Agency under the CAMELIA project (20050/06/NL/HE). P.J. Rayner is in receipt of an Australian Professorial fellowship (DP1096309). The authors thank the computer team at LSCE for the computing support and resources. Finally, the author would like to thank two anonymous reviewers and the Editor (A. Desai) for their comments, questions, and suggestions, which significantly contributed to improve this paper.
Byrd, R. H., P. Lu, J. Nocedal, and C. Zhu (1995), A limited memory algorithm for bound constrained optimization, *SIAM J. Sci. Stat. Comput.*, **16**(5), 1190–1208.

Carvalhais, N., M. Reichstein, P. Ciais, G. J. Collatz, M. D. Mahecha, L. Montagnani, D. Papale, S. Rambal, and J. Seixas (2010), Identification of vegetation and soil carbon pools out of equilibrium in a process model via eddy covariance and biometric constraints, *Global Change Biol.*, **16**(10), 2813–2829.

Cheng, Y.-B., Q. Zhang, A. I. Lyapustin, Y. Wang, and E. M. Middleton (2014), Impacts of light use efficiency and fPAR parameterization on gross primary production modeling, *Agric. Forest Meteorol.*, **189**, 187–197.

Chevallier, F., M. Fisher, P. Peylin, S. Serrar, P. Bousquet, F.-M. Bréon, A. Chédin, and P. Ciais (2005), Inferring CO₂ sources and sinks from satellite observations: Method and application to TOVS data, *J. Geophys. Res.*, **110**, D24309, doi:10.1029/2005JD006390.

Ciais, P., et al. (2013), Carbon and other biogeochemical cycles, in *Climate Change 2013: The Physical Science Basis. Contribution of Working Group I to the Fifth Assessment Report of the Intergovernmental Panel on Climate Change*, edited by T. F. Stocker et al., Cambridge Univ. Press, Cambridge, U. K., and New York.

Claverie, M., E. F. Vermote, M. Weiss, F. Baret, O. Hagolle, and V. Demarez (2013), Validation of coarse spatial resolution LAI and FAPAR time series over cropland in southwest France, *Remote Sens. Environ.*, **139**, 216–230.

Connolly, J., N. T. Roulet, J. W. Seaquist, N. M. Holden, P. M. Lafleur, E. R. Humphreys, B. W. Heumann, and S. M. Ward (2009), Using MODIS derived fPAR with ground based flux tower measurements to derive the light use efficiency for two Canadian peatlands, *Biogeosciences*, **6**, 225–234.

D'Odorico, P., A. Gonsamo, B. Pinty, N. Gobron, N. Coops, E. Mendez, and M. E. Schaepman (2014), Intercomparison of fraction of absorbed photosynthetically active radiation products derived from satellite data over Europe, *Remote Sens. Environ.*, **142**, 141–154.

Field, C. B., R. B. Jackson, and H. A. Mooney (1995), Stomatal responses to increased CO₂: Implications from the plant to the global scale, *Plant, Cell Environ.*, **18**, 1214–1225, doi:10.1111/j.1365-3040.1995.tb00630.x.

Friedlingstein, P., et al. (2006), Climate carbon cycle feedback analysis: Results from the C4MIP model intercomparison, *J. Clim.*, **19**, 3337–3353.

Garrigues, S., et al. (2008), Validation and intercomparison of global leaf area index products derived from remote sensing data, *J. Geophys. Res.*, **113**, G02028, doi:10.1029/2007JG000635.

Giering, R., T. Kaminski, and T. Slawig (2005), Generating efficient derivative code with TAF, *Future Gener. Comput. Syst.*, **21**, 1345–1355, doi:10.1016/j.future.2004.11.003.

Groenendijk, M., et al. (2011), Assessing parameter variability in a photosynthesis model within and between plant functional types using global Fluxnet eddy covariance data, *Agric. Forest. Meteorol.*, **151**, 22–38, doi:10.1016/j.agrformet.2010.08.013.

Jacquemoud, S., W. Verhoef, F. Baret, C. Bacour, P. J. Zarco-Tejada, G. P. Asner, C. Francois, and S. L. Ustin (2009), PROSPECT plus SAIL models: A review of use for vegetation characterization, *Remote Sens. Environ.*, **113**, S6–66, doi:10.1016/j.rse.2008.01.026.

Jung, M., et al. (2011), Global patterns of land–atmosphere fluxes of carbon dioxide, latent heat, and sensible heat derived from eddy covariance, satellite, and meteorological observations, *J. Geophys. Res.*, **116**, G00J07, doi:10.1029/2010JG001566.

Kaminski, T., W. Knorr, M. Scholze, N. Gobron, B. Pinty, R. Giering, and P.-P. Mathieu (2012), Consistent assimilation of MERIS FAPAR and atmospheric CO₂ into a terrestrial vegetation model and interactive mission benefit analysis, *Biogeosciences*, **9**(8), 3173–3184.

Kaminski, T., et al. (2013), The BETHY/JSBACH Carbon Cycle Data Assimilation System: Experiences and challenges, *J. Geophys. Res. Biogeosci.*, **118**, 1414–1426, doi:10.1002/jgrg.20118.

Kato, T., W. Knorr, M. Scholze, E. Veenendaal, T. Kaminski, J. Kattge, and N. Gobron (2013), Simultaneous assimilation of satellite and eddy covariance data for improving terrestrial water and carbon simulations at a semi-arid woodland site in Botswana, *Biogeosciences*, **10**, 789–802, doi:10.5194/bg-10-789-2013.

Keenan, T., R. Garcia, A. D. Friend, S. Zaehle, C. Gracia, and S. Sabate (2009), Improved understanding of drought controls on seasonal variation in Mediterranean forest canopy CO₂ and water fluxes through combined in situ measurements and ecosystem modelling, *Biogeosciences*, **6**(8), 1423–1444.

Keenan, T. F., E. Davidson, J. W. Munger, and A. D. Richardson (2013), Rate my data: Quantifying the value of ecological data for models of terrestrial carbon cycle, *Ecol. Appl.*, **23**(1), 273–286, doi:10.1890/12-0747.1.

Knorr, W., and J. Kattge (2005), Inversion of terrestrial ecosystem model parameter values against eddy covariance measurements by Monte Carlo sampling, *Global Change Biol.*, **11**, 1333–1351, doi:10.1111/j.1365-2486.2005.00977.x.

Knorr, W., T. Kaminski, M. Scholze, N. Gobron, B. Pinty, R. Giering, and P.-P. Mathieu (2010), Carbon cycle data assimilation with a generic phenology model, *J. Geophys. Res.*, **115**, G04017, doi:10.1029/2009JG001119.

Krinner, G., N. Viovy, N. de Noblet-Ducoudré, J. Ogee, J. Polcher, P. Friedlingstein, P. Ciais, S. Sitch, and I. C. Prentice (2005), A dynamic global vegetation model for studies of the coupled atmosphere-biosphere system, *Global Biogeochem. Cycles*, **19**, GB1015, doi:10.1029/2003GB002199.

Kuppel, S., P. Peylin, F. Chevallier, C. Bacour, F. Maignan, and A. D. Richardson (2012), Constraining a global ecosystem model with multi-site eddy-covariance data, *Biogeosciences*, **9**, 3757–3776.

Kuppel, S., F. Chevallier, and P. Peylin (2013), Quantifying the model structural error in carbon cycle data assimilation systems, *Geosci. Model Dev.*, **6**, 45–55, doi:10.5194/gmd-6-45-2013.

Kuppel, S., P. Peylin, F. Maignan, F. Chevallier, G. Kiely, L. Montagnani, and A. Cescatti (2014), Model–data fusion across ecosystems: From multi-site optimizations to global simulations, *Geosci. Model Dev. Discuss.*, **7**, 2961–3011, doi:10.5194/gmdd-7-2961-2014.

Lasslop, G., M. Reichstein, J. Kattge, and D. Papale (2008), Influences of observation errors in eddy flux data on inverse model parameter estimation, *Biogeosciences*, **5**, 1311–1324, doi:10.5194/bg-5-1311-2008.

Liang, S. (2000), Numerical experiments on the spatial scaling of land surface albedo and leaf area index, *Remote Sens. Rev.*, **19**, 225–242.

MacBean, N., F. Maignan, P. Peylin, C. Bacour, F.-M. Bréon, and P. Ciais (2015), Using satellite data to improve the leaf phenology of a global Terrestrial Biosphere Model, *Biogeosci. Discuss.*, **12**, 13,311–13,373, doi:10.5194/bgd-12-13311-2015.

McCallum, I., W. Wagner, C. Schmullius, A. Shvidenko, M. Obersteiner, S. Fritz, and S. Nilsson (2009), Satellite-based terrestrial production efficiency modeling, *Carbon Balance Manage.*, **4**, 8, doi:10.1186/1750-0680-4-8.

Meroni, M., C. Atzberger, C. Vancutsem, N. Gobron, F. Baret, R. Lacaze, H. Eerens, and O. Leo (2013), Evaluation of agreement between space remote sensing SPOT-VEGETATION FAPAR time series, *IEEE Trans. Geosci. Remote Sens.*, **51**, 1951–1962.

Migliavacca, M., M. Meroni, L. Busetto, R. Colombo, T. Zenone, G. Matteucci, G. Manca, and G. Seufert (2009), Modeling gross primary production of agro-forestry ecosystems by assimilation of satellite-derived information in a process-based model, *Sensors*, **9**, 922–942.

Monteith, J. L. (1977), Climate and the efficiency of crop production in Britain, *Philos. Trans. R. Soc. London, Ser. B*, **281**(980), 277–294.

Moore, D. J. P., J. Hu, W. J. Sacks, D. S. Schimel, and R. K. Monson (2008), Estimating transpiration and the sensitivity of carbon uptake to water availability in a subalpine forest using a simple ecosystem process model informed by measured net CO₂ and H₂O fluxes, *Agric. For. Meteorol.*, **148**, 1467–1477.
Parton, W. J., D. S. Schimel, C. V. Cole, and D. S. Ojima (1987), Analysis of factors controlling soil organic matter levels in Great Plains grasslands, *Soil Sci. Soc. Am. J.*, **51**(5), 1173-1179.

Peng, S. S., S. L. Piao, Z. H. Shen, P. Ciais, Z. Z. Sun, S. P. Chen, C. Bacour, P. Peylin, and A. P. Chen (2013), Precipitation amount, seasonality and frequency regulate carbon cycling of a semi-arid grassland ecosystem in Inner Mongolia, China: A modeling analysis, *Agric. For. Meteorol.*, **178**, 46-55, doi:10.1016/j.agrformet.2013.02.002.

Pielke, R. A., Sr., R. Avissar, M. Raupach, A. J. Dolman, X. Zeng, and A. S. Denning (1998), Interactions between the atmosphere and terrestrial ecosystems: Influence on weather and climate, *Global Change Biol.*, **4**, 461-475, doi:10.1046/j.1365-2486.1998.t01-1-00176.x.

Pitman, A. J. (2003), The evolution of, and revolution in, land surface schemes designed for climate models, *Int. J. Climatol.*, **23**, 479-510.

Quaife, T., P. Lewis, M. De Kauwe, M. Williams, B. E. Law, M. Disney, and P. Bowyer (2008), Assimilating canopy reflectance data into an ecosystem model with an Ensemble Kalman Filter, *Remote Sens. Environ.*, **112**(4), 1347-1364, doi:10.1016/j.rse.2007.05.020.

Rahman, H., and G. Dedieu (1994), SMAC: A Simplified Method for the Atmospheric Correction of satellite measurements in the solar spectrum, *Int. J. Remote Sens.*, **16**(1), 123-143.

Rayner, P. J., M. Scholze, W. Knorr, T. Kaminski, R. Giering, and H. Widmann (2005), Two decades of terrestrial carbon fluxes from a carbon cycle data assimilation system (CCDAS), *Global Biogeochem. Cycles*, **19**, GB2026, doi:10.1029/2004GB002254.

Reichstein, M., et al. (2005), On the separation of net ecosystem exchange into assimilation and ecosystem respiration: Review and improved algorithm, *Global Change Biol.*, **11**(9), 1424-1439, doi:10.1111/j.1365-2486.2005.001002.x.

Ricciuto, D. M., A. W. King, D. Dragoni, and W. M. Post (2011), Parameter and prediction uncertainty in an optimized terrestrial carbon cycle model: Effects of constraining variables and data record length, *J. Geophys. Res.*, **116**, G01033, doi:10.1029/2010JG001400.

Richardson, A. D., et al. (2008), Statistical properties of random CO₂ flux measurement uncertainty inferred from model residuals, *Agric. For. Meteorol.*, **148**, 38-50, doi:10.1016/j.agrformet.2007.09.001.

Richardson, A. D., et al. (2010), Estimating parameters of a forest ecosystem C model with measurements of stocks and fluxes as joint constraints, *Oecologia*, **164**(1), 25-40, doi:10.1007/s00442-010-1628-y.

Sacks, W. J., D. S. Schimel, R. K. Monson, and B. H. Braswell (2006), Model-data synthesis of diurnal and seasonal CO₂ fluxes, at Niwot Ridge, Colorado, *Global Change Biol.*, **12**, 240-259, doi:10.1111/j.1365-2486.2005.01059.x.

Santaren, D., P. Peylin, N. Viovy, and P. Ciais (2007), Optimizing a process-based ecosystem model with eddy-covariance flux measurements: A pine forest in southern France, *Global Biogeochem. Cycles*, **21**, GB2013, doi:10.1029/2006GB002834.

Santaren, D., P. Peylin, C. Bacour, P. Ciais, and B. Longdoz (2014), Ecosystem model optimization using in-situ flux observations: Benefit of Monte-Carlo vs. variational schemes and analyzes of the year-to-year model performances, *Biogeosci. Discuss.*, **10**(100), 18,009-18,064.

Seixas, J., N. Carvalhais, C. Nunes, and A. Benali (2009), Comparative analysis of MODIS-FAPAR and MERIS-MGVI datasets: Potential impacts on ecosystem modeling, *Remote Sens. Environ.*, **113**(12), 2547-2559, doi:10.1016/j.rse.2009.07.018.

Sitch, S., et al. (2008), Evaluation of the terrestrial carbon cycle, future plant geography and climate-carbon cycle feedbacks using 5 Dynamic Global Vegetation Models (DGVMs), *Global Change Biol.*, **14**, 2015-2039, doi:10.1111/j.1365-2486.2008.01626.x.

Stöckli, R., T. Rutishauser, D. Dragoni, J. O'Keefe, P. E. Thornton, M. Jolly, L. Lu, and A. S. Denning (2008), Remote sensing data assimilation for a prognostic phenology model, *J. Geophys. Res.*, **113**, G04021, doi:10.1029/2008JG000781.

Tarantola, A. (2005), *Inverse Problem Theory and Methods for Model Parameter Estimations*, 342 pp., Society for Industrial and Applied Mathematic, Philadelphia, Pa.

Verbeeck, H., P. Peylin, C. Bacour, D. Bonal, K. Steppe, and P. Ciais (2011), Seasonal patterns of CO₂ fluxes in Amazon forests: Fusion of eddy covariance data and the ORCHIDEE model, *J. Geophys. Res.*, **116**, G02018, doi:10.1029/2010JG001544.

Wang, Y.-P., R. Leuning, H. A. Cleugh, and P. A. Coppin (2001), Parameter estimation in surface exchange models using nonlinear inversion: How many parameters can we estimate and which measurements are most useful?, *Global Change Biol.*, **7**, 495-510.

Weiss, M., and F. Baret (1999), Evaluation of canopy biophysical variable retrieval performances from the accumulation of large swath satellite data, *Remote Sens. Environ.*, **70**, 293-306.

Weiss, M., F. Baret, S. Garrigues, and R. Lacaze (2007), LAI and fAPAR CYCLOPES global products derived from VEGETATION. Part 2: Validation and comparison with MODIS collection 4 products, *Remote Sens. Environ.*, **110**(3), 317-331, doi:10.1016/j.rse.2007.03.001.

Williams, M., A. D. Richardson, M. Reichstein, P. C. Stoy, P. Peylin, H. Verbeeck, N. Carvalhais, M. Jung, D. Y. Hollinger, and J. Kattge (2009), Improving land surface models with FLUXNET data, *Biogeosci. Discuss.*, **6**(2), 2785-2835.

Wullschleger, S. D., H. E. Epstein, E. O. Box, E. S. Euskirchen, S. Goswami, C. M. Iversen, J. Kattge, R. J. Norby, P. van Bodegom, and X. Xu (2014), Plant functional types in Earth system models: Past experiences and future directions for application of dynamic vegetation models in high-latitude ecosystems, *Ann. Bot.*, **114**(1), 1-16.

Zobitz, J. M., D. J. P. Moore, T. Quaife, B. H. Braswell, A. Bergeson, J. A. Anthony, and R. K. Monson (2014), Joint data assimilation of satellite reflectance and net ecosystem exchange data constrains ecosystem carbon fluxes at a high-elevation subalpine forest, *Agric. For. Meteorol.*, **195-196**, 73-88.