Evaluation of the surface energy balance system (SEBS) applied to ASTER imagery with flux-measurements at the SPARC 2004 site (Barrax, Spain)

Accurate quantiﬁcation of the amount and spatial variation of evapotranspiration is important in a wide range of disciplines. Remote sensing based surface energy balance models have been developed to estimate turbulent surface energy ﬂuxes at di ﬀ erent scales. The objective of this study is to evaluate the Surface Energy Balance System 5 (SEBS) model on a landscape scale, using tower-based ﬂux measurements at di ﬀ erent land cover units during an overpass of the ASTER sensor over the SPARC 2004 experimental site in Barrax (Spain). Additionally, the e ﬀ ect of replacement of empirical roughness functions in the model by ﬁeld estimates or literature values is investigated. with ﬂux measurements over land cover compared to cases di land are in the measurement footprint. Furthermore SEBS underestimates sensible heat ﬂux, which is in one source models.


Introduction
Accurate quantification of the amount of evapotranspiration and its spatial distribution 15 is important in research in fields of hydrology, agronomy and meteorology. This information aids in precision irrigation, determining crop water stress and water use of vulnerable ecosystems, and predicting weather and climate change.
As surface processes contributing to evapotranspiration are complex and relevant physical parameters are difficult to measure, accurate mapping of evapotranspiration 20 remains an important challenge.
Evapotranspiration links the water balance to the surface energy balance. It is controlled by the availability of moisture and available energy at the surface. Models that simulate the interaction between the land surface and the atmosphere are known as Soil-Vegetation-Atmosphere Transfer schemes (SVAT). Conventional SVATs are based 25 on point measurements and give only reliable results at the local scale. These mod-els cannot be upscaled to larger areas because of the heterogeneity of land surfaces and the dynamic nature of heat transfer processes. Recently, surface energy balance models have been developed that use remote sensing data, e.g. TSEB (Norman et al., 1995), SEBAL (Bastiaanssen et al., 1998), and SEBS (Su, 2002). Nowadays a number of sensors is available with different spectral and spatial resolutions, airborne and 5 satellite based. Medium resolution satellite images, e.g. images from the Advanced Spaceborne Thermal Emission and Reflection Radiometer (ASTER) and Landsat TM, are inexpensive and can be used to derive relevant physical parameters for surface energy balance models at the landscape scale. Field measurements of these parameters are expensive and have a limited spatial support. The empirical derivation of surface 10 aerodynamic properties from remote sensing images, however, can give large errors. In this study the sensitivity of the Surface Energy Balance System (SEBS) model to all input parameters is investigated. Furthermore, the effect of replacing empirical calculations of surface aerodynamic parameters with values from field measurements or literature is tested using two scenarios: 15 1. Empirical Scenario: empirical equations are used to derive canopy height, h c , surface roughness for momentum transport, z 0M , and zero-plane displacement height, d 0 ; 2. Field Scenario: values from field measurements or literature are used to replace the empirical calculations of the Empirical Scenario. 20 Outcomes of remote sensing based energy balance models are usually evaluated with a limited number of tower-based flux measurements (Bastiaanssen et al., 1998;French et al., 2003;Su, 2002), supported by footprint modelling (Horst and Weil, 1992). The validity of these models on the landscape scale is uncertain, in particular for a-typical and heterogeneous land cover types, because flux measurements are typically situated 25 at homogeneous sites.
In this study the two scenarios of the SEBS model are evaluated with distributed flux measurements at the landscape scale. The focus is on the sensible heat flux as this is 1167 the most critical part of surface energy balance models. The main research questions are: 1. How do SEBS estimated sensible heat fluxes relate to flux measurements for homogeneous and heterogeneous land cover in the footprint?
2. What is the sensitivity of SEBS derived sensible heat fluxes to errors in input 5 data?
3. Do the results of SEBS improve significantly by using field measurements or literature values of surface aerodynamic properties instead of empirical formulations?
In this study the turbulent sensible heat fluxes are modelled with SEBS, using ASTER in combination with field measurements collected during the SPARC2004 field campaign 10 of July 2004 in Barrax, Spain (Su et al., 2008). First, the SEBS model will be introduced followed by a description of the data sets used in this research. Next, results of the Empirical Scenario will be evaluated as well as the sensitivity of the SEBS model. Finally, the differences between the Field Scenario and the Empirical Scenario will be discussed.

The Surface Energy Balance System (SEBS)
The Surface Energy Balance System (SEBS) is a single-source model, which estimates atmospheric turbulent fluxes and surface evaporative fraction from remote sensing data. Single-source models make no distinction between the energy balance, temperature and vapour regimes of the vegetation canopy and the soil surface (Friedl, 20 2002). SEBS uses an excess resistance term that accounts for the fact that the roughness lengths for heat and momentum are different for canopy and soil surface. Dual source models use two sets of resistances across which individual, local, single-source models are applied: a bare soil scheme and a vegetation resistance scheme. Within the canopy interaction between soil and vegetation components is included. Finally, 25 a single aerodynamic resistance connects the combined canopy with the atmosphere. A single-source model uses only one resistance and assumes that all surfaces can be represented by one effective temperature and humidity value. The physical detail of dual-source models requires more ancillary data and calibration compared to single-source models. Therefore, single-source schemes are more widely utilized for 5 operational monitoring and forecasting (Timmermans et al., 2005a).
In the current setup SEBS requires three sets of input data: (1) Products derived from remote sensing data: albedo, emissivity, temperature and the Normalized Difference Vegetation Index (NDVI) to derive local surface roughness parameters; (2) Meteorological parameters collected at a reference height (air pressure, temperature, 10 relative humidity, wind speed); (3) Radiation data (downward solar radiation, downward longwave radiation).
The SEBS algorithm (Su, 2002) was implemented with the PCRaster Python Library (Karssenberg et al., 2007). The model consists of three modules: (1) Derivation of energy balance terms; (2) Submodel to derive roughness length for heat transfer (Su 15 et al., 2001); (3) Submodel to derive stability parameters. Using these three modules, the energy balance for limiting cases (i.e. completely wet or dry pixels) can be resolved. Consequently, the energy balance terms, relative evaporation, evaporative fraction and evapotranspiration flux can be derived for all pixels.
The submodel to derive roughness length for heat transfer can use either field esti-20 mates and literature values or, when this data is not available, empirical relationships with NDVI for surface aerodynamic properties. The empirical relation between the roughness length of momentum transfer, z 0M [m], and NDVI used in this implementation of SEBS is (Su, 2001): where "max" is the NDVI for maximum vegetation cover in the image. The height of the 1169 canopy, h c [m], is derived with (Brutsaert, 1982): The displacement height, d 0 [m], is calculated with (Brutsaert, 1982): In this paper two scenarios will be evaluated and compared: the Empirical Scenario 5 uses Eqs. (1), (2) and (3), while the Field Scenario uses values from field measurements or literature for canopy height, h c , surface roughness for momentum transport, z 0M , and zero-plane displacement height, d 0 . Data used in the Field Scenario will be described in the next section. For a detailed description of SEBS, the reader is referred to Su (2002).

Data description
The data set (Su et al., 2008)  took place during two weeks in mid-summer when natural surfaces are under waterstress. This agricultural area, which is partly irrigated, comprises of land covers ranging from completely bare soil to fully vegetated parcels with canopy heights from several centimetres up to two meters. Corn, sun flower and other crops are irrigated by pivots. The area is a plateau and is situated at an average of 700 m above mean sea level. in 5 thermal infrared bands -are used from an overpass on 18 July at 11:00:29 UTC (Fig. 1). After atmospheric correction, broadband surface albedo is derived from 6 shortwave channels following Liang (2001) and vegetation cover from NDVI using 2 VNIR bands and a method described by Carlson and Ripley (1997). Surface temperature and sur-5 face emissivity are retrieved from a temperature-emissivity separation (TES) algorithm (Gillespie et al., 1999) using all five atmospherically corrected TIR bands.

Ground data
Other input data for the SEBS model consisted of meteorological data and radiation measurements (air temperature, relative humidity, wind speed, air pressure, incoming 10 solar radiation). Furthermore, validation data has been acquired from flux towers at six locations in the study area. These consist of measurements of incoming and outgoing shortwave and longwave radiation as well as soil and sensible heat fluxes at six sites. Sensible heat flux measurements were done by using Large Aperture Scintillometers (LAS), an eddy correlation system, a sonic anemometer and fast response thermo-15 couples (Su et al., 2008). The measurements were performed over typical land-cover units, comprising a forest nursery (F), a wheat stubble field (W), vineyard (V-LAS and V-EC), and a sunflower field (S). An additional measurement was carried out at the edge of a corn field (C) adjacent to the vineyard site, covering either site depending on wind-direction. The locations of the flux measurements are indicated in Fig. 1. 20 A land cover map ( Fig. 1) is derived from a combination of the SPARC 2004 landuse database and a supervised classification of the 15 m resolution ASTER imagery, which was afterwards resampled to the 90 m gridsize. This land cover map was used to provide the aerodynamic surface properties as a replacement for Eqs. (1), (2) and (3) in the Field Scenario. Assessment of canopy height, h c , was done using field 25 observations, whereas surface roughness for momentum transport, z 0M , and zeroplane displacement height, d 0 , were assigned using tabulated values from existing literature (Brutsaert, 1982) (Table 1). 1171 For a successful interpretation of the turbulent fluxes, it is necessary to know the composition and extension of the up-wind source area. In the footprint methodology, see for example Horst and Weil (1992), Schmid (2002) and Schuepp et al. (1990), the flux measured at a particular height is considered as originating from an array of point sources where the relative weight depends on the location relative to the measuring 5 instrument. By spatial integration of a source weight function, it is possible to determine the footprint, i.e. the contribution of each point source as a percentage of the total flux. Here, these relative contributions are calculated using the approximate analytical model described by Hsieh et al. (2000), using the distance along the main wind direction, the measuring height, the friction velocity and surface roughness at the tower 10 sites. The comparison of model results with tower observations is done by applying the relative contributions as a weighting function to the pixels concerned.
A footprint consists of 5 to 12 ASTER pixels, depending on the vertical location and the type of the measuring device. For comparison with the field measurements of fluxes the weighted average flux within the footprint is calculated using: where x i is the modelled flux at pixel i with relative contribution w i and N, the number of pixels in the footprint. The sum of the weights equals 1. The footprint weighted variance of the modelled fluxes is calculated as: 20 For the measurements by the Large Aperture Scintillometer (LAS) the situation is slightly different. The LAS consists of a transmitter and a receiver between which the scintillation, which is a measure for the turbulence, is measured along the pathlength between the transmitter and receiver. As such, the "location" of the instrument, in terms of determining the source area, is not a point but a line, and the source area is determined by assuming an array of point measurements along the pathlength. Wang et al. (1978) show that scintillations produced by turbulence near the centre of the path-5 length contribute more to the signal measured than scintillations near the transmitter and receiver, according to a spatial weighting function. Therefore the relative contribution determined by the approximate analytical model is combined with this spatial weighting function as described in Meijninger (2003).

Empirical Scenario
First, the model is run without the use of lookup tables for z 0M , h c and d 0 , and a land cover map, but solely using remote sensing and meteorological data. The output of the SEBS model consists of the spatial distribution of net radiation (R n ), soil heat flux (G 0 ), sensible heat flux (H), latent heat flux (λE ) and evaporative fraction (Λ) at the moment of satellite overpass (18 July 2004, 11:00:29 UTC) (Fig. 2). Although the land cover map was not used as an input, the results have a clear relationship with land cover (Fig. 1). This can be explained by the sharp contrast between irrigated fields and their surroundings for the model inputs of surface temperature and NDVI. 20 Irrigated fields have a high latent heat flux and a low sensible heat flux. Inside the pivot irrigation fields the sensible heat flux is at times negative due to the oasis effect (Brutsaert, 1982;Stull, 1988). Because of the high vegetation cover at the pivots, the soil heat flux is low at these sites. Net radiation is high in these areas, because the albedo of fully vegetated areas generally is much lower compared to bare soils. Table 2 25 shows the average results of the Empirical Scenario for each land cover type.

1173
Footprint modelling revealed that only the measurements over the sunflower field (S) and the forest nursery (F) cover most of these land cover types and are therefore representative for the land cover type. The source area of the instrument situated at the edge of the corn (C) is influenced by vineyard, while the measurement over wheat stubble (W) is slightly influenced by corn. The eddy correlation measurement in the 5 vineyard (V-EC) is a little influenced by wheat stubble. The LAS over vineyard (V-LAS) measures a mix of different land covers and is also influenced by wheat stubble. This means that the comparison of flux measurements with modelled values is for most measurement locations related to a number of different land cover types.
Comparison of the measured and modelled fluxes (weighted for the footprint) is 10 shown in Fig. 3. Figure 3 shows that SEBS estimates at locations with a homogeneous footprint, i.e. S, V-EC, and F, show the best correspondence with the ground measurements. The trends also compare well for these land cover types: the well irrigated sun flower pivot (S) has the lowest H, while the dryer vineyard (V-EC) and forest nursery (F) have a much higher H. 15 On the other hand, the results show that for the classes with a mixed footprint, i.e. C, W and V-LAS the results are worse. For these locations SEBS underestimates the sensible heat flux, compared to the ground measurements. These errors are probably introduced by the aggregation of heterogeneous surfaces to ASTER pixels. When all weighting pixels fall inside one homogeneous land cover type, the results are better. 20 This is the case for the eddy correlation measurements of sensible heat flux over the vineyard (V-EC), the scintillometer measurement over the sun flower field (S) and the sonic anemometer measurement over the forest nursery, which show a difference with the modelled H of 43, 36 and 65 W/m 2 respectively. The scintillometer measurement over the vineyard (V-LAS), however, shows a larger difference (91 W/m 2 ), because the 25 footprint also partly covers bare soil and wheat stubble and pixels in the ASTER image are mixed. The sonic anemometer that measured the corn field gave very high values of sensible heat flux for this irrigated field. The placement on the northern edge of the corn pivot in combination with a northern wind (5 • ) biases the measurement by including other land covers, mainly vineyard, in the footprint. Obviously, SEBS results are less comparable to flux measurements when the measurement footprint covers different land cover types. This is most likely caused by side effects in the footprint, when sharp transitions between the land cover types exist, as is the case with pivots. Sharp transitions have a large effect on roughness, which is 5 difficult to parameterize in surface energy balance models. Furthermore, SEBS and comparable models do not take into account wind direction, which influences z 0M . Moreover, care should be taken in these typical transition zones where the Monin-Obhukov similarity theory might be violated (Brutsaert, 1982).
At many locations SEBS underestimates the sensible heat flux; a phenomenon seen 10 more often at high sensible heat flux rates when dealing with one source models (Huntingford et al., 2000;Kustas et al., 1996).

Model sensitivity
The results shown in the previous section are computed using several empirical relations included in SEBS. We hypothesize that SEBS results can be improved by replac- 15 ing some of these functions by values derived from field measurement or literature. First, however, we need to determine the most sensitive parameters that should be replaced in the Field Scenario. Su (2002) derived the sensitivity of the sensible heat flux analytically. He estimated that the sensitivity is in the order of 20 W/m 2 when the input variables are within 50% of their actual values, which is around 20% relative to 20 the mean sensible heat flux. We perform sensitivity analysis on all input maps (surface temperature, emissivity, NDVI, albedo, DEM) and field measurements (reference temperature, wind speed, air pressure, relative humidity, height of boundary layer, incoming shortwave radiation) in a non-analytical way.

25
1175 Sensitivity (S i ) of a SEBS input is defined here as: where H 0 , H + and H − are the sensible heat flux predicted by SEBS when the input equals its reference value i 0 , 1.25·i 0 and 0.75·i 0 , respectively, with reference values used for all other inputs. For air temperature, however, a deviation of 1% was used, 5 since a 25% deviation exceeds its physical limits. For similar reasons absolute deviations of 0.1, 0.5 and 1.0 K were used for surface temperature. The sensitivity has been analysed for all validation locations separately (V-EC, V-LAS, S, F, W and C) and the average effect on each land cover class. Table 3 shows the parameters for which H modelled by SEBS is relatively insensitive 10 (S i <10%). The low sensitivity to changes in surface elevation (DEM), emissivity ( ), NDVI, albedo (α), relative humidity (H f ), incoming solar radiation (K ↓ s ) and height of the planetary boundary layer (Z pbl ) was expected, because the derivation of sensible heat flux requires only meteorological parameters at reference height and surface temperature. This means that the calculation of H in SEBS is independent of other surface 15 energy balance terms in contrast with most other models (Su, 2002). In some cases, however, other input parameters can have an effect on the sensible heat flux estimation. This is the case when both the calculation of H exceeds the model restriction that the minimal H is equal to the wet-limit sensible heat flux and the maximum of H equals the dry limit sensible heat flux. This is determined using net radiation and soil heat flux 20 calculations. The values of H outside the dry-and wet-limit occur when the iteration in the submodel for the derivation of stability parameters (Su, 2001) does not converge. This happens with some land cover types when varying the K ↓ s (Fig. 4). Figure 5 shows the parameters for which H modelled by SEBS is sensitive (S i >10%): windspeed (u), air temperature (T a ) and air pressure (P ). These parameters, that are 25 measured at the meteorological station, are directly used in the calculation of sensible heat flux.
Most parameters show a comparable sensitivity at the different locations. Exceptions to this are the shortwave incoming radiation (Fig. 4) and surface temperature (Fig. 6). Sensible heat flux at the sunflower pivot (S) is especially sensitive to a small error in surface temperature. Figure 7 shows the sensitivity of H to the roughness height of momentum (z 0M ), zero 5 plane displacement height (d 0 ) and canopy height (h c ). These parameters are empirically derived from their relationship with NDVI (Eqs. 1, 2 and 3) and are correlated. Table 4 shows that the differences between estimated values using NDVI and field estimations/literature values can differ more than the deviations used in the sensitivity analysis and the 50% limit stated by Su (2002). Furthermore, there is low correlation 10 (r 2 =0.4) between modelled and measured values of h c , z 0M and d 0 . For example, the modelled z 0M of grassland is higher than vineyard, which means that homogeneous grassland is rougher than heterogeneous vineyard. The relation of z 0M with the arrangement of surface objects, their spacing and physical height is not well represented by the relation between NDVI and z 0M (Eq. 1). This was also found by Tasumi et al.

15
(2000), Hasager and Jensen (1999), Jasinski and Crago (1999) and Timmermans et al. (2005b). The high z 0M value modelled for bare soil is probably related to errors in the classification. Because of these shortcomings of the Empirical Scenario, the Field Scenario will be evaluated in the next section.

25
A lowering of the sensible heat flux for most land cover types is expected, because of an overestimation of the roughness parameters by the empirical functions. Only for vineyard these parameters were underestimated (Table 4). 1177 Table 5 shows that the Field Scenario indeed results in a lower sensible heat flux compared to the Empirical Scenario for the land cover classes bare soil, wheat stubble, grassland and water, as was expected from the overestimation of these parameters by the empirical functions (Table 4). The expected increase of the average sensible heat flux for vineyard is also observed.

5
For forest nursery, crops and corn fields, however, an increase in sensible heat flux was found (Table 5), which was not expected from the decrease in z 0M , h c and d 0 (Table 4). This is most likely caused by misclassifications in the land cover map around the pivots. The round pivots, which in reality have sharp borders between wet and dry surfaces, are bordered by a large number of mixed pixels at 90 m resolution, causing 10 misclassifications at the edges of the pivots. Also the forest nursery is surrounded by pivots, which bias the results with mixed pixels. The results for sun flower pivots are comparable in both scenarios. Compared to the flux measurements (Fig. 8), the model results improve for forest nursery (F). Also vineyard (V-EC) results improve at the location of the eddy correlation 15 system. Surprisingly, the sensible heatfux for V-EC decreased in spite of an increase in aerodynamic roughness values. This is probably caused by the complex relation between H and z 0M , h c and d 0 . The scintillometer measurement (V-LAS), which was intended to measure a heterogeneous footprint, gives similar results for both scenarios.
The sensible heat flux estimates over the wheat stubble (W) get worse due to the 20 fact that pixels from a corn pivot are included in the footprint of the wheat stubble measurement, which gives side effects as explained before. Because the wheat stubble is very dry, the higher sensible heat flux measured by the eddy correlation system seems to be more realistic for this land cover than both values modelled with SEBS.
The results for the sunflower field (S) and the corn field (C) are the same for both 25 scenarios. The empirical algorithms for the derivation of surface aerodynamic properties seem to perform better for homogeneous covers like corn and sunflower, instead of the heterogeneous surfaces they are developed for. The roughness parameters that have been derived empirically are comparable with values from field measurements and literature for these homogeneous covers (see Table 4). However, the difference between modelled and measured sensible heat flux over the corn field is large for both scenarios. This can be related to the measurement position and wind direction at satellite overpass, as stated before. Intuitively, the modelled values seem to be more appropriate because at the irrigated field the evapotranspiration 5 is large, resulting in a high latent heat flux and a low, or even negative sensible heat flux.
The large improvement for the forest nursery by using the Field Scenario can be explained by the fact that the footprint mainly covers the forest nursery plot and field estimated roughness values are representative. The same can be concluded for one 10 of the vineyard sites (V-EC). However, it should be noted that the surface roughness values as given for the forest nursery and particularly the vineyard may vary considerably depending on wind direction. This is due to the structure of these orchards, i.e. row orientation. 15 In this paper SEBS has been evaluated at the landscape scale, using distributed field measurements of sensible heat flux H. SEBS is capable of estimating H in the same order of magnitude as the field measurements. Furthermore, the contrasts between well irrigated pivots and other land cover types are observed in both the field measurements and the SEBS results. Standard deviations in field measurements of H are 20 similar to standard deviations of H modelled by SEBS.

Discussion and conclusions
SEBS estimated sensible heat fluxes relate well to measured fluxes when the footprint of the measurements covers only one land cover type. When different land cover types are included in the footprint, errors are introduced by land surface variables aggregated to the size of an ASTER pixel. At the current modelling scale of 90 m this 25 becomes an important issue when dealing with e.g. pivot irrigation, causing mixed pixels on its fringe. Using higher resolution imagery (e.g. Airborne Hyperspectral Sensor, 1179 AHS) or disaggregating surface temperature (e.g. Kustas et al., 2003) to a scale in accordance with the objects in the footprint of the flux measurements can give better validation results. However, disaggregating surface temperature should be done very accurately as the sensitivity analysis showed that SEBS is most sensitive to the surface temperature. Especially in well irrigated fields, sensible heat flux estimates by SEBS 5 can deviate up to 70% with 0.5 K difference in surface temperature.
Although sensitivity of SEBS derived sensible heat flux to errors in surface aerodynamic parameters is smaller compared to surface temperature, the errors in the estimation of these parameters from remote sensing images using empirical relations can be larger and exceed the 50% limit of input accuracy for many land cover types. In 10 the previous paragraphs this has been investigated using two scenarios of input data. In the Empirical Scenario solely remote sensing data is used for derivation of z 0M , h c and d 0 . In the Field Scenario lookup tables with field observations and literature values of these parameters were used in combination with a land cover map. The average sensible heat fluxes per land cover unit have been calculated for both scenarios. 15 The Field Scenario, however, only resulted in a small improvement, compared to the Empirical Scenario, where the field flux measurements are placed within a homogeneous footprint. The Field Scenario can even worsen the result for heterogeneous footprints, by creating sharp borders related to the discrete borders in the land cover map. 20 It can be concluded that the use of higher resolution remote sensing data can better represent the shape of the footprint of the flux measurements and that the effect of wind direction on surface roughness for momentum transport should be incorporated in SEBS in order to relate SEBS results to flux measurements, independent of the location of the measurements. This should cope with heterogeneity within the footprint 25 and varying roughness for row crops.  hydrol-earth-syst-sci.net/6/85/2002/. 1167, 1169, 1170, 1175, 1176, 1177: An evaluation of two models for estimation of the roughness height for heat transfer between the land surface and the atmosphere, J. Appl. Meteorol., 40, 1933-1951, doi:10.1175/1520-0450(2001 (Fig. 4) and surface temperature (Fig. 6)). u = wind speed, T a = air temperature and P = air pressure.