Vegetation dynamics and soil water balance in a water-limited Mediterranean ecosystem on Sardinia, Italy

Mediterranean ecosystems are commonly heterogeneous savanna-like ecosystems, with contrasting plant functional types (PFTs, e.g., grass and woody vegetation) competing for the water use. Mediterranean ecosystems are also commonly characterized by strong inter-annual rainfall variability, which inﬂuences the distributions of PFTs 5 that vary spatially and temporally. With the objective to investigate interactions between vegetation dynamics, soil water budget and land-surface ﬂuxes in a water-limited ecosystem, an extensive ﬁeld campaign in a Mediterranean setting was performed. Also a vegetation dynamic model (VDM) is coupled to a 3-component (bare soil, grass and woody vegetation) Land surface model (LSM). The case study is in Orroli, situated 10 in the mid-west of Sardegna within the Flumendosa river basin. The landscape is a mixture of Mediterranean patchy vegetation types: trees, including wild olives and cork oaks, di ﬀ erent shrubs and herbaceous species. Land surface ﬂuxes, soil moisture and vegetation growth were monitored during the May 2003–June 2006 period. Interest-ingly, hydrometeorological conditions of the monitored years strongly di ﬀ er, with dry 15 and wet years in turn, such that a wide range of hydrometeorological conditions can be analyzed. The coupled VDM-LSM model is successfully tested for the case study, demonstrating high model performance for the wide range of eco-hydrologic conditions. The use of the VDM in the LSM is demonstrated to be essential when study-ing the climate-soil-vegetation interactions of these water-limited ecosystems. Results 20 demonstrate also that vegetation dynamics are strongly inﬂuenced by the inter-annual variability of atmospheric forcing, with grass leaf area index changing signiﬁcantly each spring season according to seasonal rainfall amount.


EGU
An attractive compromise is the coupled VDM and LSM of Montaldo et al. (2005), which starting from Cayrol et al. (2000) and Nouvellon et al. (2000) VDMs developed a parsimonious and robust coupled model for grass dynamics only. In the model the VDM provides the grass leaf area index (LAI) evolution through time, and the LSM uses this to compute the land surface fluxes and update the soil water contents. They suc-5 cessfully tested the coupled model for two case studies of water-limited grass fields in California (USA) and North Carolina (USA). Even in such "simple" ecosystems characterized by only one PFT, they demonstrated the significant role of vegetation dynamics on soil water balance modeling in water-limited conditions, and the importance of including the VDM for correctly predicting land surface fluxes and soil water balance. 10 Here we further develop the coupled model of Montaldo et al. (2005) for including 3 cover types (bare soil, grass and woody vegetation) typical of more complex heterogeneous ecosystems, and test the model for a sufficient long data set including years characterized by different hydro-meteorological conditions.
The case study site is within the Flumendosa river basin on Sardinia, which is one of 15 the regions of Italy most affected by water deficits. There is therefore an urgent need to exploit advanced observation and simulation technologies to provide a better understanding of the water balance regime for the entire island and for its major catchments. In this sense, the dam system of the Flumendosa river constitutes the water supply for much of southern Sardinia, including the island's largest city, Cagliari (about 350 000 20 inhabitants in the urban area). The case study site is a natural patchy mixture of grass and woody vegetation, typical of Mediterranean ecosystems. During May 2003-June 2006 a micrometeorological tower has been installed and an extended field campaign has been conducted (Detto et al., 2006). This paper addresses the following objectives: The landscape is a patchy mixture of Mediterranean vegetation types: trees, mainly wild olive (Olea sylvestris) of height approximately 3.5-4.5 m, and a few cork oaks (Quercus suber ) of height approximately 6-7 m, shrubs (Asparagus acutifolius and Rubus ulmifolius), creepers of the wild olive trees (Crataegus azarolus and Smilax aspera), and C3 herbaceous (grass) species (Asphodelus microcarpus, Ferula comunis, 15 Bellium bellidioides, Scolymus hispanicum, Sonchus arvensis, Vicia sativa, Euphorbia characias, Dancus cerota, Bellis perennis; monocotyledons: Avena fatua, Hordeum murinum) that are present in live form only during wet seasons and reach heights of approximately 0.5 m. The soil thickness varies from 15-40 cm, bounded from below by a rocky layer of basalt. This impervious layer leads to tree and shrub rooting systems 20 that expand horizontally. The root zone depth is coincident with the soil depth for these thin soils. The climate at this site is typically Mediterranean-maritime, with a mean historical  annual precipitation of 690 mm (raingage data from the nearby village of Nurri), and mean historical monthly precipitations ranging between 103 mm in Decem- 25 ber and 12 mm in July (Fig. 1a) EGU annual value of 13.9 • C, mean monthly values ranging between 23.1 • C in July and 6.1 • C in January (Fig. 1b).

Field measurements
An extended field campaign was carried out from May 2003 to June 2006, during which micrometeorological, soil moisture (θ), and vegetation dynamics measurements were 5 conducted.

Micrometeorological tower
A 10 m tower was instrumented to measure land-atmosphere fluxes of energy, water, and carbon in addition to key state variables. The tower is surrounded by wild olive trees, grass and bare soils. It includes a Campbell Scientific CSAT-3 sonic anemometer 10 and a Licor-7500 CO 2 /H 2 O infrared gas analyzer at 10 m above the ground to measure velocity, temperature and gas concentrations at 10 Hz for the estimation of latent heat (LE), and sensible heat fluxes (H) through standard eddy-correlation methods (e.g., Brutsaert, 1982;Garratt, 1992). Half hourly statistics were computed and recorded by a 23X data logger (Campbell Scientific Inc., Logan, Utah). The effect of the gentle 15 slope of the plateau was removed by an axis rotation (Detto et al., 2006) and the Webb-Pearman-Leuning adjustment (Webb et al., 1980) was applied. Three infrared transducers, IRTS-P (Apogee Instrument, accuracy of 0.3 • C) were used to measure the surface temperature (T s ) of the different PFTs. One IRTS-P observed the skin temperature of a tree (wild olive) canopy at 3.5 m height above the 20 ground and with a canopy view zenith angle of ∼70 • , another observed either bare soil or grass (depending on the season) at 1.6 m above the ground with a canopy view zenith angle of ∼50 • , and the third sensor was placed at a greater height (10 m above the ground, view zenith angle of ∼40 • ) to observe a composite mixture of trees and soil or trees and grasses (depending on the season). Introduction

EGU
The incoming and outgoing shortwave and longwave radiation components were measured by a CNR-1 (Kipp & Zonen) integral radiometer positioned at 10 m with a hemispherical field of view. Soil heat flux (G) was measured at two different locations close to the tower, one in an open patch (4 m from the tower) and one under a tree canopy of wild olive (5.5 m from the tower), with thermopile plates, HFT3 (REBS), 5 buried at 8 cm below the soil surface. Two thermocouples (per plate) were buried at 2 and 6 cm, and one frequency domain reflectometer probe (FDR Campbell CS615) per plate was buried horizontally at 5 cm, as needed to estimate changes in the stored energy above the plates (see HFT3 instruction manual edited by Campbell Sci.).
Precipitation was measured by an ARG100 (Environmental Measurements Limited) 10 tipping bucket raingauge. Recorded precipitation time series are shown in Fig. 1 (monthly) and Fig. 2 (daily). Data gaps (13.5% of the total half hour values) exist mainly due to power supply failures and maintenance operations. Rainfall and meteorological observations during the data gaps are filled with data of nearby stations located close to the town of Nurri (∼5 km from the tower).

15
The two-dimensional footprint model of Detto et al. (2006), previously tested for this site, was used for interpreting eddy-correlation measurements in the context of the contributing land cover area. The footprint (source area) of eddy-correlation flux measurements changes in size and direction through time with the wind speed and stability of the flow. This variation can be exploited to sample various mixtures of the relative 20 fractions of the different surface types.
The observed years were characterized by strongly contrasted hydro-meteorological conditions and offer a wide range of conditions understanding the ecosystem behavior. A comparison between monthly mean historical  and observed (2003)(2004)(2005)(2006) precipitation and temperature is provided in Fig. 1 10 LAI was measured indirectly through a ceptometer (Accupar model PAR-80, Decagon devices inc., Washington USA), which measures the Photosynthetically Active Radiation (PAR) in the 400-700 nm waveband, and estimates the LAI from these readings (details in the instruction manual edited by Decagon devices inc.). LAI measurements were performed during the entire observation period, especially during the 15 grass growth season (Fig. 3). LAI of the woody vegetation (Olea sylvestris that is the predominant woody vegetation type and which includes the creepers) changes moderately throughout the year (Fig. 3), whereas, the green leaf area of the herbaceous species increases rapidly with winter and spring precipitation and vanishes for the dry summer (Fig. 3). 20 Finally, specific leaf areas (LAI divided by dry biomass) of predominant grass (=0.01 m 2 gDM −1 ) and woody vegetation (=0.005 m 2 gDM −1 ) species were measured directly (weighing the dry biomass). These values are needed to connect the biomass estimates of the vegetation dynamics model (discussed below) with the traditional LAI values, as reported in Fig. 3 (Fig. 4). The two images depict the contrast between the spring (bottom) and summer (top) land cover present 5 surrounding the field site. The top image characterizes the land cover when the soil moisture conditions are very dry (θ ≈0.08) and green herbaceous cover is absent such as is typical in the Sardinian Summer. The bottom image depicts the land cover conditions after a long wet period (θ ≈0.4, Fig. 2) particularly propitious for plant growing, so that the bare soil was nearly absent while the flourishing grasses reached their 10 maximum growth in those days. The 6S model of Vermote et al. (1997) was used to correct the images for atmospheric effects. Details are provided in Detto et al. (2006). A supervised classification scheme based on the parallepiped algorithm (Richards, 1999) allows for distinguishing "woody-vegetation" (WV) from "non-woody-vegetation" (NWV, i.e., bare soil or grass 15 according to the time period) from the images. The widely used normalized difference vegetation index (NDVI) (e.g. Gamon et al., 1995;Carlson and Ripley, 1997;Scanlon et al., 2002) was computed from the surface reflectance values averaged over ranges of wavelengths in the visible red and NIR regions of the spectrum. Following Detto et al. (2006), in each map pixel the fraction of woody vegetation cover is estimated as 20 NDVI i j /NDVI max , where NDVI i j is the NDVI value of a particular grid cell and NDVI max is the spatial maximum of the particular NDVI map. The NDVI/NDVI max map of the field around the tower (the tower is in the center of the map) for the DOY=220, 2003 is computed (Fig. 5a). Note that NDVI/NDVI max values of WV pixels are greater than 0.6, so that the color bar of the Fig. 5a is modified for a better contrast of the WV pixels. 25 The combined use of the footprint model (see Sect. 2.2.1) and the high-resolution satellite images allows us to interpret the eddy-correlation observed surface flux and distinguish the source area of each PFT and bare soil to the measured flux, using the EGU methodology of Detto et al. (2006). Figure 5b reports the frequency distribution of the fraction of WV cover (f f p,W V ) in the footprint of the micrometeorological observations for the entire dataset.
We note that f f p,W V is mainly in the range of 0.1-0.22 with the peak of the distribution close to 0.15.

5
In this section we describe the land surface model (LSM) and the vegetation dynamics model (VDM). The essence of this modeling coupling is that the VDM provides the leaf area index (LAI) evolution through time for each PFT, which are then used by the LSM for computations of the energy partitioning between soil and vegetation. Model parameters are defined in Table 1.

The land surface model
The LSM predicts dynamics of water and energy fluxes at the land surface on a halfhour time step. It is derived from the LSM of Montaldo and Albertson (2001) including three components in the land surface: bare soil and two vegetated components (e.g., grass and WV). The states of surface temperature and moisture are estimated through 15 the force-restore method (Noilhan and Planton, 1989;Montaldo and Albertson, 2001). The root zone supplies the bare-soil and vegetation with soil moisture for evapotranspiration, and controls the infiltration and runoff mechanisms. The base of the root zone represents the lower boundary of the LSM. Note that the equations for surface temperature and three components (H, G and the net radiation, R n ) of the energy balance are 20 the same as Noilhan and Planton (1989), with the only difference that in the proposed LSM they are applied separately for each land cover component.
In the unsaturated soil the Clapp and Hornberger (1978) relationships are used to describe the non-linear dependencies of volumetric soil moisture (θ) and hydraulic conductivity (k) on the matric potential (ψ). The soil water balance equation of the where θ rz is the soil moisture of the root zone, d rz is the root zone depth, I bs is the infiltration rate on bare soil, I wv and I gr are the throughfall rates infiltrating into the soil covered by WV and grass respectively, q D the rate of drainage out of the bottom of 5 the root zone, which is estimated using the unit head gradient assumption (Albertson and Kiely, 2001), E bs is the rate of bare soil evaporation, T wv and T gr are the rates of transpiration of WV and grass respectively, f v,W V is the fraction of green WV area per unit of ground area, f v,gr is the fraction of green grass vegetation area per unit of ground area, and f bs (=1−f vt,W V −f vt,gr ) is the fraction of bare soil, where f vt,W V and f vt,gr are 10 the total WV and grass vegetation area (including dead vegetation) respectively. The total evapotranspiration, E T , is equal to f bs E bs +f v,W V T wv +f v,gr T gr . As in the original Noilhan and Planton (1989) model, the throughfall rate is modeled through a balance equation of the intercepted water by the canopy reservoir (its capacity is a function of the LAI), which produces throughfall when the reservoir is saturated. 15 In the original model of Montaldo and Albertson (2001) the infiltration rate was computed through a saturation excess mechanism, which is not suitable for this case study, typically characterized by Hortonian overland flow due to the thin soil and the semi-arid conditions (e.g., Chow et al., 1988). Hence, the infiltration model was updated for including the infiltration excess mechanism. According to this mechanism, the actual 20 infiltration rate of the x th land cover type, I x , is taken as the minimum of the rainfall rate (or throughfall in the case of vegetated components) and an infiltration capacity, I * , based on the Philip's infiltration equation where t k is the time since the onset of infiltration, S s is the sorptivity, and A is a con- EGU soil properties and the root zone moisture content at the start of the storm event. For eliminating t k , the Milly (1986) approach (based on the time compression approximation) is used, so that I * only depends on the cumulative infiltration (in addition to S s and A).
Transpiration rates are estimated by the Penman-Monteith equation (e.g., Brutsaert, 5 1982, p. 224), in which the aerodynamic resistance and the canopy resistance are estimated for each PFT distinctly. The canopy resistances that account for environmental stresses are estimated following Montaldo et al. (2005), and are described in Appendix A. The aerodynamic resistances are estimated as function of wind velocity through the transfer coefficient for water vapor, C E (Garratt, 1999, equation 3.57), according to the Monin-Obukhov similarity theory. C E and the heat transfer coefficient (used in H estimates) account for atmosphere stability (Garratt, 1999, equation 3.47), with the flux profile functions for stable and unstable conditions estimated through equations (3.35), (3.36), (3.39) and (3.40) of Garratt (1999). Finally, the actual rate of bare soil evaporation is determined by where E p is the potential evaporation estimated by the Penman equation (e.g., Brutsaert, 1982, equations 10.15, 10.16 and 10.19), θ g is the surface soil moisture, and α(θ g ) is a rate-limiting function, estimated by the polynomial function of Parlange et al. (1999).

The vegetation dynamic model
The VDM computes change in biomass over time from the difference between the rates of biomass production (photosynthesis) and loss, such as occur through respiration and senescence (e.g., Larcher, 1995;Cayrol et al., 2000). The VDM distinguishes WV and grass components, and is adapted from Montaldo et al. (2005), which derived 25 a VDM for grass species starting from the Nouvellon et al. (2000) model. Since we are modeling semi-arid regions, we assume that water availability is the major factor 230 Introduction  (Larcher, 1995;Mouillot et al., 2001).
In the VDM of WV, four separate biomass states (compartments) are tracked: green leaves (B g ), stem (B s ), living root (B r ), and standing dead (B d ). The biomass [g DM m −2 ] components are simulated by ordinary differential equations integrated nu-5 merically at a daily time step (Nouvellon et al., 2000;Cayrol et al., 2000;Arora and Boer, 2005;Montaldo et al., 2005): where P g is the gross photosynthesis, a a , a s and a r are allocation (partitioning) coefficients to leaves, stem and root compartments (a a + a s + a r = 1), R g , R s and R r are the respiration rates from leaves, stem and root biomass, respectively, S g , S s and S r are the senescence rates of leaves, stem and root biomass, respectively, and L a is the 15 litter fall.
The key term of the VDM, P g , is computed using the approach of Montaldo et al. (2005) (Table 2). Starting from a simplified form of Fick's law applied to gas exchange in plants (Larcher, 1995;Lambers et al., 1998) and the Nouvellon et al. (2000) model, Montaldo et al. (2005) derived a simplified expression that estimates P g by mainly the photosinthetically active radiation, PAR, and other routinely monitored variables (wind velocity, and air humidity and temperature, EGU terms of the VDM are computed as in Montaldo et al. (2005), while the photosynthesis allocation to leaves, stem and roots is estimated by adapting the new approach of Arora and Boer (2005) (see this reference for a deeper discussion of the allocation coefficients behavior). The equations for the photosynthesis allocation estimates and the other terms of (4)-(7) are described in Table 2 and model parameters are defined 5 in Table 1. The VDM of grass distinguishes only three biomass compartments (green leaves, roots and standing dead) and the biomass components are simulated using (4), (6) and (7), respectively.
Leaf area index values are estimated from the biomass through linear relationships 10 (Hanson et al., 1988;Nouvellon et al., 2000;Arora, 2003;Montaldo et al., 2005): where LAI and LAI d are the green and dead leaf area index of the x th land cover type, respectively. The total leaf area index LAI t is then estimated by

Coupling the land surface model and the vegetation dynamic model
The LSM is then coupled with the VDM. VDM provides LAI values of WV and grass daily by (8), which are then used by the LSM for computing the evapotranspiration estimate (e.g., equation A1), energy flux and the soil water content in the root-zone by 20 (1) at half-hour time step. Leaf area index values are also used for updating the total fraction of vegetation cover, f vt , of the generic PFT values through (Montaldo et al., 2005): with f v the fraction of vegetation of the generic PFT. not available (Fig. 3). An analysis of the influence of key environmental factors on the vegetation dynamics interannual variability is finally provided.

The coupled VDM-LSM
The VDM-LSM coupled model was calibrated for the case study, comparing observed 10 and simulated time series of the energy balance terms, θ, and LAI through a trialand-error procedure. Note that for comparing micrometeorological observations and model predictions of LE and H we used the time varying footprint of the tower for estimating the fraction of WV cover (see Sect. 2.2.1), while for the soil moisture budget we used the fraction of land covers of the field monitored by the soil moisture probes 15 (f vt,W V =0.25 and f vt,gr =0.6). The Table 1 reports the calibrated parameter values. Note that all parameters are held constant throughout the study period.
Observed surface temperature, a key indicator of the energy balance, is well simulated for the three land cover components (Fig. 6). The accuracy of the coupled model for predicting energy balance terms is demonstrated by the results shown in Figs. 7, 20 8 and 9. Net radiation (R n ), sensible heat flux (H) and soil heat flux (G) dynamics are all well estimated on the whole (Fig. 7). The scatter plots of Fig. 8 are on a daily time scale.
Daily E T rates and cumulative E T are shown in Fig. 9, and a scatter plot of E T time series are in Fig. 8c. Cumulative E T , which is important for soil water balance EGU missing observed E T values (28% of the study period) are replaced by modeled values at that time (30% of total observed E T ). Root zone soil moisture dynamics are also sufficiently well simulated with rmse of 0.076 for the calibration period and 0.053 for the validation period (Fig. 2). Note that main inaccuracy of soil moisture modeling are during Summer rain events because 5 measured soil moisture peaks do not agree with the rain gage input.
LAI is very well simulated for both grass and WV as can be seen in Figure 3. In particular, the grass model calibrated for the 2003 and 2004 years is able to simulate well the 2006 growth and the following decrease of LAI. The dynamics of the WV LAI confirms the strong tolerance of the WV species to prolonged droughts.

Inter-annual variability of vegetation dynamics
The different hydro-meteorological conditions of the observed years ( Fig. 1 and Sect. 2.2.1) affect significantly grass vegetation dynamics, as can be well depicted by Fig. 10, where LAI grass time series during the early part (DOY 30-180, which coincides with the grass growth season in Sardinia) of each year are compared (Fig. 10a). 15 Soil moisture dynamics (Fig. 10b) and E p , which reflects the relevant atmospheric forcing, including solar radiation, air humidity and wind velocity (Fig. 10c), are also compared. In this case study the dynamics of grass LAI responded readily to meteorological forcing due the limited soil depth and the absence of available groundwater, which is typical in Sardinian basins. In the year 2003 after a average January precipitation 20 (Fig. 1a), an extremely high precipitation occurred in February but then precipitation strongly decreased during the key months for the LAI growth in Sardinia (March, April and May) so that low LAI values were observed during the high E p period, and finally a very dry Summer occurred. During the Spring 2004 the best hydrologic conditions occurred for this site-high soil moisture until the end of the Spring season, when high 25 values of E p are observed-so that extremely high LAI values were predicted (solid line in Fig. 10a). In 2005 the environmental conditions were not so favorable and less grass growth were predicted. Indeed grass growth was limited before (DOY 100-120) due to 234 Introduction Finally an interesting correlation between the grass LAI dynamics and precipitation during the two typical growth months, April and May, was found. In Fig. 11  produced favorable soil moisture conditions for grass growth but was not sufficient due to the low E p (Fig. 10c).

Conclusions
The monitored 3 hydrologic years in Orroli (Sardinia) were characterized by strong inter-annual variability of hydro-meteorological conditions, such as is typical of Mediter-20 ranean semi-arid ecosystems. The inter-annual variability of atmospheric forcing significantly impacts soil moisture and vegetation dynamics, in particular during the Spring and early Summer seasons, which are key seasons for Sardinian water resources planning and management.
The yearly variability of hydro-meteorological conditions offered a wide range of con- EGU observation and was able to accurately predict vegetation dynamics, soil water balance and land surface fluxes. In particular the evapotranspiration, a key term in these ecosystem, is very well predicted (to within 99% of the total observed evapotranspiration of the study period). The typical woody vegetation species of Sardinia, representative of the broader 5 Mediterranean water-limited region, confirm a strong tolerance to prolonged drought, such as occurred in the Summer of 2003. Even with the extreme dry conditions the WV species didn't wilt and LAI was still high (>3), showing moderate changes throughout the year. Instead, the dynamics of grass LAI responded readily to meteorological forcing due the limited soil depth and the absence of available groundwater, which is typical in Sardinian basins. This allowed to find an interesting correlation between the precipitation, and the grass LAI dynamics during the Spring season, the growth season in Sardinia. The correlation was found to be high when the values of precipitation and LAI are aggregated at 15-day time intervals, and there is a sufficient time lag (15-days) 15 between the forcing (precipitation) and the answer (LAI). The relationship between LAI and precipitation is not linear showing a threshold value of LAI for the highest precipitation values. These results highlight the high correlation between grass dynamics and precipitation forcing in these ecosystem, and seem encouraging the climate change warning that is pervading the Mediterranean European community. Indeed, if climate 20 change effects are producing warmer conditions together with decreases of precipitation and increases of drought durations in many semi-arid regions (Feddema, 1999;Lelieveld et al.;2002;Moonen et al., 2002;Ragab and Prudhomme, 2002;Ventura et al., 2002) the impact on vegetation (i.e., land cover types and distribution) will be significant. Introduction

Appendix A
The canopy resistance (r c ) is used both in the LSM for the transpiration estimates using the Penman-Monteith equation and in the VDM for the photosynthesis estimate (see Table 2). Following Montaldo et al. (2005) it is estimated with a typical Jarvis (1976) approach: where f 1 , f 2 and f 3 are stress functions of soil moisture, air temperature (T a ), and vapor pressure deficit (VPD). The soil moisture effect is treated differently for grass and WV due to the particular resistance to water stress of WV species modeled in this Sardinian ecosystem (see Detto et al. (2006) for details on this function) where θ l i m and θ wp depend on the type of vegetation (e.g., Larcher, 1995;Eagleson, 2002). The effect of temperature on the stomata is treated by (Nouvellon et al., 2000;Larcher, 1995)

EGU
where T a,min , T a,opt and T a,max are characteristics of the plant types (Larcher, 1995;Eagleson, 2002). Finally, the effect of VPD on stomata opening is treated by (Jarvis, 1976) Acknowledgements. This work was supported by the Ministero dell'Università e della   [-] PAR extinction coefficient 0.5 0.5 E2002 ξ a [-] Parameter controlling allocation to leaves 0.6 0.55 cal. ξ s [-] Parameter controlling allocation to stem 0. ξ a + ξ s + ξ r = 1; λ = e −k e LAI For grass a a = ξ a +Ωλ 1+Ω[1+λ−f 1 (θ)] a r = ξ r +Ω(1−f 1 (θ)) 1+Ω[1+λ−f 1 (θ)] ξ a + ξ r = 1 Respiration R g = m a f 4 (T ) B g + g a a a P g C1986; R s = m s f 4 (T ) B g + g s a s P g N2000; R r = m r f 4 (T ) B r + g r a r P g M2005  The map of NDVI/NDVI max of WV of the field around the tower (the tower is in the center of the map) determined from the Quickbird image for DOY=220, 2003; note that NDVI/NDVIMAX values of WV pixels are greater than 0.6, so that the color bar is modified for a better contrast of the WV pixels.