Identification of a general light use efficiency model for gross primary production

Non-stationary and non-linear dynamic time series analysis tools are applied to multi-annual eddy covariance and micrometeorological data from 44 FLUXNET sites to derive a light use efficiency model for gross primary production on a daily basis. The extracted typical behaviour of the canopies in response to meteorological forcing leads to a model formulation allowing for a variable influence of the environmental drivers temperature and moisture availability modulating the light use efficiency. Thereby, the model is applicable to a broad range of vegetation types and climatic conditions. The proposed model explains large proportions of the variation of the gross carbon uptake at the study sites while the optimized set of six parameters is well defined. With the parameters showing explainable and meaningful relations to site-specific environmental conditions, the model has the potential to serve as basis for general regionalization strategies for large scale carbon flux predictions.


Introduction
The atmosphere and the terrestrial biosphere are tightly coupled through the exchange of energy and matter (Monteith and Unsworth, 2008). A central component of this 15 coupling is the assimilation and release of CO 2 by photosynthesis and respiration; these opposed fluxes modulate substantially the global carbon cycle (Schimel et al., 2001). The rising CO 2 concentration in the atmosphere and the associated changing climate factors have implications on the functioning of ecosystems and hence provoke a feedback on the carbon cycle in turn (Cox et al., 2000). 20 Consequently, quantifying the global carbon balance under current conditions and predicting its characteristics in a future environment with enhanced atmospheric CO 2 concentrations implies the quantification of CO 2 uptake and respiration rates of ecosystem as well as descriptions of their main drivers (Cramer et al., 2001). With plants trading water vapour with CO 2 , the CO 2 assimilation additionally has effects on the Some studies go even further to overcome the restrictions of lacking information content in the available data to constrain model processes and to make the quantification of gross primary productivity better applicable for lager scales and chose very parsimonious model structures without the implementation of explicit physiological processes occurring at cell and leaf scale. Such models treat canopies as functional units 5 aggregating and averaging processes over space and time. A very popular approach uses the concept of light use efficiency which represents the ratio of carbon biomass production per unit of absorbed light (Watson, 1947;Monteith, 1972;Monteith and Unsworth, 2008). Various studies have proven the light use efficiency to be quite constant over the day, a fact which makes the light use efficiency concept suitable for 10 daily-step models (Ruimy et al., 1995;Rosati and Dejong, 2003;Sims et al., 2005). The light use efficiency parameter has been implemented in ecosystem models as a constant (Landsberg and Waring, 1997;Veroustraete et al., 2002) or modified by restricting environmental factors such as temperature and vapour pressure deficit with predefined functions (Potter et al., 1993;McMurtrie et al., 1994;Prince et al., 1995;15 Veroustraete et al., 2002;Xiao et al., 2004;Yuan et al., 2007;Makela et al., 2008). The LUE approach has been used as a stand-alone application (Yuan et al., 2007;Makela et al., 2008) as well as integrated in ecosystem models (Coops et al., 2005), it has been driven with ground measurement data as well as combined with remote sensing data (Potter et al., 1993;Law and Waring, 1994;Prince et al., 1995). The MODIS-GPP 20 algorithm (Running et al., 1999;Zhao et al., 2005) is principally based on the light use efficiency approach, too. But despite numerous models proposed, many questions have remained unanswered and primary production modelling on landscape scale is still "an active area of research with issues remaining to be solved on the leaf, stand, and landscape level" (Hilker et al., 2008, p. 418). 25 Basic relationships like the light use efficiency are well suited as a starting point for model identification procedures in a top-down fashion. Such methodologies try to develop models specifically at the scale of interest beginning with the most robust functional relationships which are iteratively refined according to data analysis results. In 7676 Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | types. Finally, found parameters are tested for patterns which relate themselves to site specific characteristics serving -as a first step -the final aim to regionalize the model parameters.
2 Data and methods 2.1 Micro-meteorological data 5 44 forest and grassland FLUXNET sites in climate zones reaching from boreal to semiarid were chosen as data base for this study ( Table 1). The selection criterion was the existence of at least three measurement years and no lengthy measurement gaps. The selected sites are located in North America and Europe ( Fig. 1) comprising 12 coniferous forest sites, 18 deciduous, 5 mixed, 2 evergreen forests as well as 7 grasslands. Table 1 summarizes their characteristics. The data were downloaded from the web gateways of the regional FLUXNET sub-networks AmeriFlux and CarboEurope. The downloaded data including energy and carbon fluxes along with meteorological variables have measurement gaps which were filled in the following way: Short gaps up to three hours of meteorological variables were linearly interpolated. The average values 15 of the respective values at the time of day in a 14-day moving time window around the gap (Falge et al., 2001) served to fill gaps of medium length up to 4 days. Large gaps were replaced with the respective values averaged over the whole time series available. Missing data in the time series of the net CO 2 flux, F N , are replenished with the multidimensional semi-parametric spline interpolation scheme explained in Stauch and 20 Jarvis (2006). The methodology not only fills missing data but also extracts the deterministic component of the measured time series and compares well to other gap-filling techniques for eddy covariance net carbon fluxes (Moffat et al., 2007). The signal component of the net flux is finally split up into respiration and the gross flux component, F G , according to the semi-parametric method in Desai et al. (2008). The so-derived 25 flux F G is used in our further considerations.

MODIS LAI/FPAR product
The fraction of absorbed photosynthetically active radiation (FPAR) is downloaded as MODIS Land Products subset (ORNL DAAC, 2010) for each FLUXNET site. These MOD15A2 and MYD15A2 subsets provide a grid of 7×7 pixels with a size of 1 km 2 centred on the tower. MOD15A2 and MYD15A2 data retrieved from the satellites Terra 5 and Aqua are merged into one dataset according to Yang et al. (2006). Each of the 49 pixels with the same land class as the study site according to the MODIS product MOD12Q1 is taken into account. If neither Terra nor Aqua delivered a FPAR value with the main algorithm the mean of all available years (2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010) at this day of year is taken instead. The final FPAR time series is retrieved by a cubic smoothing spline 10 fitted through all data points; due to the noisiness of the data, each value is weighted according to its inverse difference to the multi-annual mean. The absorbed photosynthetically active radiation (APAR) is finally calculated as product of PAR measured at the FLUXNET sites and FPAR. 15 To derive functional descriptions for the seasonal evolution of the gross carbon flux F G a non-stationary regression framework such as the Captain Toolbox for Matlab is perfectly suited (Pedregal et al., 2007). Two of its powerful tools based on Kalman filtering and smoothing techniques are applied in this study: Dynamic linear regression, DLR, and state dependent parameter analysis, SDP, allowing for the evolution 20 of parameters to be estimated directly from time series data and hence identification of any non-stationary and/or state dependency of these parameters. In particular, the underlying regression type model in case of DLR is of the form where y is the dependent variable, x i are the regressors, c i are time dependent regression parameters and ζ is the model error series. i is the increment running from 1 to the number of regressors, N. DLR extracts the incremental temporal variations in c assuming that the parameters gradually vary with time in the form of a random walk. In this process, each sampling instant depends on the data in its direct vicinity 5 using a Gaussian weighting function. The "bandwidth" of this Gaussian window function is characterized by the noise variance ratio (NVR) which is optimized from the data via maximum likelihood prediction error decomposition as proposed by Young (1999). Kalman filtering and fixed interval smoothing algorithms are employed to estimate the model, based on the optimized NVR values and, finally, error bounds are calculated. 10 In contrast, the SDP analysis presumes that the regression parameters do not vary with time but with a state of the considered non-linear system. This means that each sample depends on the data in its vicinity in the sorted state space out of temporal order. The state variable, however, can again vary with time: 15 with u being variables representing system states. For evaluation purposes, the Nash-Sutcliff efficiency criterion, EC, is applied beside the typically used coefficient of determination, r 2 , and the bias as difference of the means of measured and modelled time series. EC is defined as the variance of residuals of predicted (P ) and observed (O) values normalized by the variance of the 20 observed values, and subtracted from one (Legates and McCabe Jr, 1999):

Data analysis and pre-processing methods
In contrary to r 2 , EC is therewith sensitive to additive and proportional differences between measured and modelled data.

7680
Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | 3 Model identification

Evaluation of the Jarvis-model
Analysing flux data at two temperate forests (Harvard Forest, UMBS) with SDP, Jarvis et al. (2004) identified the light use efficiency as expressed in the following equation: to follow a sigmoid relationship with T F , the time-delayed soil temperature (T S ): where T H is the inflection point of between its minimum and maximum level, max , and k T the rate of change. These three parameters together with the time lag parameter 10 for T F were calibrated against data of the two study sites. The model was validated at Harvard Forest over a 6-year period. The resulting parameters were shown to be well-defined. This promising light use efficiency approach serves as a starting point for identifying an generalized model scheme applicable to a broader spectrum of vegetation 15 and climate types. To do so, the original Jarvis-model is first applied to all study sites regardless of their vegetational and climatological characteristics whereas the four parameters were optimized with the non-linear least squares method. The Jarvis-model reproduces well the gross CO 2 uptake F G of boreal and temperate forests. The model performs particularly well at deciduous forests with strong seasonal dynamics (Figs. 2a Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | though, but this is often just a result of the nature of the optimization procedure to -up to a certain degree -counterbalance shortcomings of model formulations; specifically, α seems to compensate inappropriate model structures. F G of forests experiencing hot summers as well as most grasslands cannot be simulated by the Jarvis-model at all, because the assumed sigmoid temperature dependency of does not exist, letting the 5 conclusion to be drawn that a water availability proxy is lacking. Even at the fully humid study sites analysed by Jarvis et al. (2004), cross-correlations between the model residuals and a water availability measure were found. Consequently, not only the dependency of to T S has to be reconsidered, but also appropriate water availability measures have to be identified and their functional relationship to has to be derived.

Finding new model structures
In order to identify any state dependency of SDP was applied first with T S as only state variable (Eq. 2) confirming that at many sites a distinct decrease of with increasing T S (Fig. 4a) occurs. At those sites where the Jarvis-model is not able to properly reproduce the carbon flux dynamics, SDP also fails to find a clear temperature dependency 15 indicating some additional control on the -dynamics. SDP is used here to also analyse several water availability measures as potential further controls (in addition to T S ) on including the evaporative fraction (EF) being the fraction of the latent heat and the available energy, the vapour pressure deficit (VPD), the soil water content (SWC, Fig. 4b), and the antecedent precipitation index (API, Fig. 4c). Besides API, these vari-20 ables have demonstrated before to significantly affecting primary production and are frequently used in light use efficiency models as moisture availability indicator (Makela et al., 2008;Yuan et al., 2007;Potter et al., 1993;Prince et al., 1995;Heinsch et al., 2006). DLR is applied to estimate EF with the aim to better capture seasonal variations and reduce the impact of short term fluctuations by noise effects: with the sensible heat flux H. API is calculated by a weighted sum of daily precipitation values in a time window w before the current time step t (Linsley et al., 1982;Samaniego-Eguiguren, 2003): where d denotes the number of time steps before t and κ is a recession constant, 5 commonly ranging between 0.85 and 0.98 (Chow et al., 1964). Indeed, in several cases the SDP-analysis shows a distinct dependency on the water availability state variables. At most sites, however, neither the temperature nor a water availability proxy alone can explain F G satisfyingly. Therefore, an additive SDP-model is applied using T S and one of the mentioned water availability proxies (W ). The model performance varies between the water availability surrogates at the study sites, but no of them delivers the best results in every case. EF, however, appears to perform most consistently throughout the sites, whereas the use of API tends to lead to somewhat 15 higher uncertainties, i.e. the related parameter estimations are less uniquely identifiable compared to the other state variables. The nonparametric relationship between max and T S in the additive SDP model either has a sigmoidal form, or it can be described by a (sigmoidal) peak function (Fig. 5a), or no clear relationship can be identified at all (Fig. 5c). If a clear relationship between and the water availability state 20 variable exists -as in the majority of cases -it shows a threshold-like behaviour as demonstrated exemplarily in Fig. 5b and d.

Formulating the generalized model
To overcome the applicability restrictions of the basic model with the lessons learned in the SDP analysis, the sigmoid temperature function is changed to a logistic peak Introduction function, f T , which enables a decrease of with increasing temperatures after a sigmoidal shift from the minimum to the maximum level: with k T being the rate of change and T opt being the temperature at which the function is 1.

5
To allow for the effect of water availability fluctuations, a sigmoid function is used since SDP shows the tendency that at very low and very high values of the respective water availability proxies W (EF, SWC, API and VPD) there is no change of the influence on ; this behaviour gets even more obvious when taking additively both temperature and moisture in one SDP-model into account (Fig. 5). The function allowing for 10 the influence of W on max , f W , is therefore chosen to be a sigmoid function: with k W being the rate of change between the minimum and maximum levels of f W and I W being the inflection point. Both f T as well as f W are scaled between 0 and 1. To account for lag effects between the response of to temperature variations we 15 allow again for the lag effect using the lag-parameter α applied to T S (Eq. 6) as it has proven to be significant in similar light use efficiency model approaches as proposed by (Makela et al., 2006(Makela et al., , 2008 for sites in temperate and boreal climates. However, in cases of W being the main driver of as it is the case in semi-arid climates, the lag function is applied to W instead of T S . In applying α only to the main driver, the number 20 of free parameters is minimized, and the lag is anyway only apparent in a distinctive manner on a daily time step basis when the canopy has to regenerate and redevelop green tissue after a dormant period; and these periods are largely determined by the main driver as the temperature in temperate and boreal climates and a moisture proxy in semi-arid climates. Introduction The final model is formulated as follows:

Conclusions References
where APAR is the absorbed photosynthetically active radiation in MJ m −2 d −1 as product of FPAR and PAR, and p is a parameter between 0 and 1. If both temperature and humidity conditions are optimal max is reached. If no humidity dependency can 5 be detected, because there is always enough water available, and -variations can be explained by the temperature, p approaches 1 and the second term approaches zero, and vice versa. 1−p is consequently indirectly a measure for the strength of the water availability influence on a vegetation stand. APAR was chosen instead of S 0 as in the Jarvis et al. (2004) study for better comparability with other studies and because of the overwhelming evidences for the significance of the leaf area index or FPAR as scaling-factor for soil-vegetation-atmosphere-transfer processes (Watson, 1958;Monteith, 1977;Tucker and Sellers, 1986;Goetz and Prince, 1999;Gower et al., 1999;Lindroth et al., 2008), an explanatory power and intrinsic scaling factor which cannot be compensated by other environmental variables used in the light use efficiency model 15 approach.

Model calibration and evaluation
The final model formulation given by Eq. (12) comprises seven parameters including k W , max , p, T opt , k T , I W and α. To explore the sensitivity and variability for each of these parameters among different sites, a set of 750 000 Monte Carlo simulations 20 was executed at each location allowing the seven parameters to vary randomly within predefined (bio-physiologically meaningful) ranges. Using the sum of squared errors between measured and modelled F G as a performance criteria we identified the site specific parameters distributions of the 1000 best model runs. here dependent on the dominant control (temperature or water availability, see below). Analysing the 1000 best solutions of all sites reveals that the best parameters exploit a wide range within the assigned upper and lower boundaries. In contrast, the parameter k W shows some non-uniqueness in that very good model results can be achieved over some range of the parameter values. However, when compared between differ-5 ent sites, the possible range of optimal parameter values reduces to a relatively narrow range (without figure). Therefore, to avoid over-parameterization of the model we treated k W constant at the median of all the best k W -values and only the remaining six parameters, max , p, T opt , k T , I W and α are calibrated in the following using a non-linear least-squares optimization routine ("lsqnonlin"  In the following, results are described in more detail for the model with EF instead of SWC, API and VPD as water availability proxy, since (i) the calibrations with EF resulted in the best results with the lowest confidence intervals, (ii) overall, EF already performed best in the SDP-analysis, (iii) it was often successfully employed as water availability proxy in similar studies (Kustas et al., 1994;Barr et al., 2007;Yuan et al., 5 2007) and proved to be superior as explaining variable for in similar analysis (Garbulsky et al., 2010), and (iv) from a regionalization point of view there is the potential it can be retrieved by remote sensing (Crago, 1996). Results for the other moisture surrogates are shown exemplarily (Figs. 6b and 10). The resulting parameters of the optimization procedure, their confidence intervals, r 2 -values and model efficiencies EC are given in Table 2 for the EF-model.
The parameters of f T have in general greater confidence intervals than I W of f W because f T has two free parameters. However, fixation of one parameter of f T deteriorates the model performance too much. The higher p, thus the more the temperature dominates the variations of F G , and the smaller the confidence intervals of the parameters 15 k T and T opt of the corresponding temperature function tend to be. Even more pronounced is the effect vice versa: the smaller p, thus the higher the influence of EF, the smaller are the confidence intervals of I W (Fig. 8a). This fact can also be seen in the sum of error squares related to the parameter spaces of the Monte Carlo simulations described above: For high values of p, thus a high influence of the temperature on 20 F G , the T opt -values of the best parameter sets are located in a relatively narrow range, and vice versa, if p is low and EF dominates F G variations, I W is better defined. This observed characteristic is an advantage of the proposed model structure: The model is on the one hand flexible enough to simulate daily fluxes of sites with very different characteristics, but on the other hand, gives less weight to the less influencing vari-Introduction The model parameter max varied at forest sites between 0.78 g C m −2 MJ −1 at a needleleaf forest in Canada and 1.93 g C m −2 MJ −1 at a German needleleaf forest in a temperate climate with mild summers (Cfc). Deciduous forests form the highest average max with a mean of 1.25 g C m −2 MJ −1 followed by the mixed forests (1.18 g C m −2 MJ −1 ) and evergreen needleleaf forests (1.16 g C m −2 MJ −1 ). Evergreen 5 broadleaf forests have with 1.14 g C m −2 MJ −1 the lowest average max due to low values in boreal climates and those with dry summers. Regarding the climate classes with more than one forest site, Cfb-class reveals the highest average max , closely followed by Dfb; Csb and Dfc have the lowest max . At grasslands sites, max -values are surprisingly high: max -estimations reach 2.25 at Oensingen and 1.80 at Neustift and lead 10 to an average max of 1.50 g C m −2 MJ −1 . The highest max -values at Oensingen are attained in spring and autumn when temperatures are favourable but APAR is relatively low.
Optimized parameter values for p, indicating the influence of T S and EF, range between 0 and 1 and only the lowest values are omitted: F G at the Mediterranean Roc-15 carespampani and at Audubon, Arizona, is largely explained by EF with a value of p of 0.14 and 0.24, respectively, F G at Griffin, England, follows highly the course of temperature (p=0.98). Most forest sites, however, have a medium p between 0.4 and 0.8. The low values of p at Hainich and especially at Boreas and Howland can be explained by a an especially high correlation between and EF in the seasonal course. At Hainich, 20 particularly the distinct summer drought of (Reichstein et al., 2007 with a strong decrease of in late summer leads to an higher influence than at the other sites in this climate class. Values of p greater than 0.6 are mainly clustered in the forest and fully humid climate classes, whereas estimates of p at forests at summery sites as well as grasslands take values in the medium to lower range. 25 The temperature T opt at which the sites reach max is smaller than 14 • C for the all needleleaf forests but the warmest fully humid study site, a result of lower average temperatures and a high efficiency in spring and autumn. the leafs have emerged following the rise of temperature in spring and before the loss of the leafs or they occur even in summer when temperatures get not too high and there is no lack of water. The grasslands located at sites with hot and semi-arid conditions are dominated by medium to high T opt -values up to 24.5 • C in case humid and warm periods coincide. The alpine and northern Dfb grassland sites are characterized 5 by medium to low T opt -estimates corresponding to the mild average temperatures and highest efficiencies in spring and autumn. In addition to the information at which temperature max occurs, the parameter k T characterizes the steepness of the temperature sensitivity in relation to the temperature range and the vegetation period. Accordingly, deciduous forests, especially at colder sites, as well as grasslands at semi-arid sites have rather low k T -values corresponding to a sharper peak of the temperature function, whereas evergreen sites particularly with a relatively small annual temperature range feature medium to high k T -values leading to a flatter and wider peak.
The seem to react rather promptly to the triggering variables. Finally, the parameter p for the further water availability proxies SWC, API and VPD is presented in Fig. 10 predominating in B-and summery climates as well as grasslands with exception of the alpine sites. More often than for SWC and API, a higher explaining power is attributed to VPD, especially at warmer coniferous sites. The highest influence of a W substitute, however, is assigned to API at the desert grassland Audubon with a value of 0.15 assigned to p, thus a contribution of 85%. 5

Discussion
State dependent parameter estimation was used to determine typical non-parametric relationships between the light use efficiency and relevant state variables. The relationship between and T S was determined as a (sigmoidal) peak function. Limiting functions with respect to temperature have been used before: Makela et al. (2008) tested 10 a model at several coniferous study sites and used a site-specific piecewise function with a linearly increasing part and a constant value above a threshold temperature; this approach contrasts with this study, in which especially the coniferous forest sites show a relatively small efficiency amplitude during the year with the highest -efficiencies at lower temperatures followed by an efficiency-decrease at higher temperatures. Pot- 15 ter and Klooster (1999) and Yuan et al. (2007) who also modelled F G across a broad range of conditions applied a symmetric peak function, too. They estimated an optimum temperature varying with the geographical latitude and estimated the optimum by non-linear optimization merging the data from all study sites, respectively. In our study, however, it is shown that highest efficiencies occur at temperatures that vary signifi-20 cantly across study sites. Regarding moisture availability measures, various functional forms as well as different proxies have been applied previously: For example, Yuan et al. (2007) and Heinsch et al. (2006) chose linear relationships between and EF, respectively. Makela et al. (2008) represented the relationship between and VPD by an exponential function and between and SWC by Weibull-function or a sigmoidal function following Landsberg and Waring (1997). Here in most cases, SDP revealed a threshold-like response of to all proxies given was sensitive to them. But in 7690 accordance with other studies (e.g., Makela et al., 2008) is the fact that SDP-analysis showed similar functional forms of the responses to the driving variables across various boundary conditions and vegetation types. Moreover, they are appropriate means to model F G if the functions describing the -dependencies are flexible enough. Indeed, the weighted additive model formulation rather than a multiplicative approach 5 has proven to serve as a robust approach. Most light use efficiency models use a multiplication of down-regulating scalars (e.g., Potter et al., 1993;Running et al., 1999;Yuan et al., 2007) what can lead to the parameterization of maximum -values which are often not reached in reality, especially if a down-regulating variable has no strong influence at a given site, as for example the soil water content in a forest with deep roots can be. Furthermore, insensitive variables are prone to high calibration uncertainties of respective parameters; these are less important if their influence is relativized by a weighting factor as realized in the proposed model. Additionally, maximum -values vary if a modifying scalar is added or omitted in a multiplicative approach (Makela et al., 2008). Instead, the proposed additive model with a site-specific weighting of the variables' influence on leads to a maximum which is actually realized by the considered vegetation stands. In contrast to Yuan et al. (2007Yuan et al. ( , 2010 this approach is based on the assumption that maximum -values and other model parameters on a daily basis vary between forest stands and grasslands (Turner et al., 2003;Bradford et al., 2005;Schwalm et al., 2006;Kjelgaard et al., 2008;Stoy et al., 2008) and no glob-20 ally valid maximum which is reachable by all vegetation stand under ideal conditions exists and biochemical processes are universal across species. Predictions with this model approach obviously require sufficiently long measurement time series covering optimal periods for vegetation growth. The minimum required measurement period of three years presupposed in this study can be critical in this sense (Nouvellon et al., 25 2000). But this restriction will become less important when FLUXNET measurement time series get longer and are made available for such calibration studies.
With the data available the calibrated parameter max varied strongly between sites. Compared to max -values of several studies presented in the review of Goetz and 7691 tential light use efficiency in Yuan et al. (2010). Both relatively high max -values at Oensingen and also Neustift are certainly a consequence of the agricultural management with mowing (4-5 and 2-3 times a year at Oensingen and Neustift, respectively) and fertilizing Ammann et al., 2009;Schmitt et al., 2009). Furthermore, it is not unlikely that CO 2 -measurement and FPAR-retrievals are influenced 10 by the relatively small scale landscape mosaic the study site Oensingen is located in; fallow fields in vicinity of the study site or fields covered by senescent vegetation could lead to an underestimation of MODIS FPAR values and consequently an overestimation of (Ammann, 2010). Schwalm et al. (2006)  The high max at Tharandt was also observed in the study of Makela et al. (2008) and explained by thinning of the forest stand. In average highest max -values at deciduous broadleaf forest sites and lower values in mixed and deciduous needleleaf forests match with the pattern obtained by the light use efficiency model of Yuan et al. 10 (2007). In the MODIS GPP algorithm, a maximum of 1.01 g C m −2 MJ −1 is stored in the algorithm's look-up table for evergreen needleleaf forests, compared to a higher mean max of 1.16 g C m −2 MJ −1 found in this study; the respective values for deciduous broadleaf forests are 1.16 and 1.25 g C m −2 MJ −1 , and for mixed forests 1.12 and 1.18 g C m −2 MJ −1 . Thus the model values of this study are somewhat higher 15 than those applied in the MODIS GPP algorithm -and that against the background of a rather negative bias of the model (see above) thus actually even higher max -values in the measurements. The calibration of p shows that temperature has indeed a high influence on , especially in cooler ecosystems, as shown by numerous studies (Runyon et al., 1994;Chen 20 et al., 1999;Nouvellon et al., 2000;Turner et al., 2003;Schwalm et al., 2006;Yuan et al., 2007). SWC as modulating variable had the highest impact at summerdry sites and grasslands what is not surprising considering the short rooting depth and the low depth at which the SWC-measurements were made. VPD appears to influence not only in dry areas as values of p around 0.6 in boreal and temperate forests indicate; 25 indeed, the interrelation between VPD and via stomatal conductance has often been shown (Wang and Leuning, 1998;Goetz and Prince, 1999;Lagergren et al., 2005;Katul et al., 2003;McCaughey et al., 2006)  the performance of the other model configurations, but at sites with strong periodic water shortage they can explain -variations with the lowest values of p, thus the highest contribution of all f W -functions compared to the other W -variables. API is therefore considered as promising variable. Overall however, the optimization procedure assigns EF the most often the highest explaining capability on between the hydrological relevant 5 variables applied and model runs with EF lead to the best results. This is not surprising considering the correlation of and EF (Monteith and Greenwood, 1986;Schulz and Jarvis, 2004) and thus EF being an "integrator" of environmental conditions. The model with EF consequently leads to a somewhat better model performance and even allows the modelling of managed sites such as Oensingen and Neustift to a certain degree. This behaviour is supported by Stoy et al. (2009) who performed a orthonormal wavelet transformation analysis on measured CO 2 -fluxes and found a high importance of "endogenous" variables compared to purely meteorological variables and a strong coupling between λE and F G . In their correlation analysis EF (Garbulsky et al., 2010) is also determined to be the best explaining variable of . In their study, EF alone 15 explained best and their model deteriorated when adding another variable. This contradicts with our calibration which assigned EF a significant influence but never the only contribution to the variation of : values of p of nearly 1 occurred but only in one case a value of smaller than 0.2 was determined. Overall, both SDP-analysis and the model application show that both temperature and water availability influence the variation of 20 but can not explained by one variable alone. Throwing a glance at the distribution of the other model parameters in the climatevegetation matrix may lead to the assumption that the optimized parameter values have no bearing on site specific characteristics. However, in the majority of cases the parameter values can be related to the vegetation class (i.e. deciduous or evergreen), 25 the length of the vegetation period (higher or lower k T ), the season in which gets maximal, the seasonal fluctuation of LAI and the degree of its minimization in dormant periods, the start of the vegetation period in relation to the course of temperature, the temperature amplitude, or the degree of superposition of seasonal temperature and humidity course. Noticeable are for example the differences of such similar forest sites as Howland and Harvard; both are located in the same category in the climatevegetation matrix: mixed forests in a boreal climate with warm summers. While the mean annual temperature is somewhat lower at Harvard, maximum -values at Howland occur in spring and autumn whereas the maximum at Harvard occurs rather in 5 summer which explains the higher T opt indicating a higher fraction of deciduous trees than at Howland; the higher at Harvard is in line with these assumptions that Harvard needs more time to develop leafs and reaching max . The model performance as determined by r 2 , EC and the degree of parameter uncertainty militates in favour of the proposed model particularly with regard to the wide variety of ecosystem characteristics of the study sites. A further study could test the benefit of integrating other drivers such as the often discussed influence of leaf nitrogen concentrations (Sinclair and Horie, 1989;Dewar, 1996;Kergoat et al., 2008) as intrinsic variable, the saturating behaviour of for high PAR values (Ruimy et al., 1995;Turner et al., 2003;Lagergren et al., 2005;Hilker et al., 2008) or the ratio of diffuse to 15 total PAR with a proxy for cloudiness (Schwalm et al., 2006;Jenkins et al., 2007).