Consistency between In Situ , Model-Derived and High-Resolution-Image-Based Soil Temperature Endmembers : Towards a Robust Data-Based Model for Multi-Resolution Monitoring of Crop Evapotranspiration

Due to their image-based nature, “contextual” approaches are very attractive to estimate evapotranspiration (ET) from remotely-sensed land surface temperature (LST) data. Their application is however limited to highly heterogeneous areas where the soil and vegetation temperature endmembers (Tends) can be observed at the thermal sensor resolution. This paper aims to develop a simple theoretical approach to estimate Tends independently from LST images. Soil Tends are simulated by a soil energy balance model forced by meteorological data. Vegetation Tends are obtained from soil Tends and air temperature. Model-derived soil Tends are first evaluated with in situ measurements made over an irrigated area in Morocco. The root mean square difference (RMSD) between modeled and ground-based soil Tends is estimated as 2.4 ◦C. Model-derived soil Tends are next compared with the soil Tends retrieved from 90-m resolution ASTER (Advanced Remote Sens. 2015, 7 10445 Spaceborne Thermal Emission and Reflection Radiometer) data collected over two irrigated areas in Mexico and Spain. Such a comparison reveals a strong consistency between model-derived and high-resolution image-based soil Tends. A recent contextual ET model (SEB-1S) is then applied to 90-m resolution and to 1-km resolution (aggregated) ASTER data using the model-derived or image-based Tends as the input. The RMSD between 90-m resolution SEB-1S and in situ ET is estimated as 65 and 82 W·m−2, and the RMSD between 1-km resolution SEB-1S and aggregated SEB-1S ET is estimated as 78 and 56 W·m−2 for the image-based and model-derived Tends, respectively. In light of the above results, Tends should be estimated a priori when contextual models are applied to low resolution images. Moreover, the consistency over highly heterogeneous areas between model-derived and high-resolution image-based Tends provides a meaningful basis for developing mixed modeling observational approaches.


Introduction
Crop evapotranspiration (ET) is by far the main outward water flux over semi-arid irrigated areas.In the Mediterranean countries, for instance, ET can represent up to 80% of the consumptive uses of water [1].In those regions, accurate estimates of spatially-averaged ET are hence required for efficient crop irrigation and water resources management.Remote sensing is in this regard one of the most promising and cost-effective techniques for mapping and monitoring ET over broad areas.Especially multi-sensor/multi-resolution remote sensing data have a strong potential to provide ET data at multiple scales: the crop field scale where irrigation volumes are allocated and the regional scale where long-term management decisions are made.
Among the considerable variety of existing approaches to estimate ET from remote sensing data, the most widely-used approach is to force the FAO-56 method [2][3][4] with NDVI (Normalized Difference Vegetation Index) data [5].It consists of computing ET as a function of the potential ET derived from meteorological measurements, a crop-specific coefficient and the crop phenology estimated from NDVI-derived fractional green vegetation cover (f vg ).The FAO-56 method forced by NDVI time series has excellent operational application capabilities, notably due to a remarkable coherence between the complexity of the model parameterization and the availability of input (meteorological and NDVI) data.However, one main limitation of the FAO-56 method is that the coefficient(s) used (either in its single-coefficient or dual-coefficient form) is calibrated using field measurements, which, given the spatial and temporal heterogeneity in these coefficients, restricts its application to large areas [6].Moreover, the soil moisture constraint on evapotranspiration in the dual-coefficient form is based on a soil water budget model and a priori knowledge of irrigation volumes and dates.
Another remote sensing variable relevant to ET monitoring is the land surface temperature (LST).LST is a signature of both ET and soil water availability, especially under non-energy-limited conditions.It is currently observed by thermal sensors, such as ASTER (Advanced Spaceborne Thermal Emission and Reflection Radiometer), Landsat-8 and MODIS (Moderate-Resolution Imaging Spectroradiometer) at a spatial resolution of 90 m, 100 m and 1 km, respectively.The operational use of thermal remote sensing for hydrological and water resource management studies has focused on regional-scale applications at the kilometric spatial resolution [7,8].The reason is that the revisit time of high-spatial resolution sensors (16 days for ASTER and Landsat) is relatively long compared to the quick soil moisture changes that generally occur in irrigated semi-arid regions.The temporal frequency of 1-km resolution MODIS data (one or two times per day) is instead very relevant for monitoring ET on an operational basis.
LST-based ET models fall into two broad ranges of models based on (i) the residual method and (ii) the evaporative fraction (EF, defined as the ratio of latent heat to available energy) method.The first category of models evaluates ET as the residual term of the surface energy balance, which is the available energy at the surface minus the sensible heat flux, e.g., TSEB (Two-Source Energy Balance) [9] and SEBS (Surface Energy Balance System) [10].The second category evaluates ET as the available energy times EF.The latter models are often called "contextual" because the methodology involves determining the wet and dry boundaries of LST, which may or may not be present at the observation resolution within the observation scene.An overview of these contextual models is presented in [11,12].In practice, physical boundaries of LST are estimated in a two-step procedure.First, the temperature endmembers (Tends) corresponding to fully dry and wet conditions for both soil/vegetation components are located within the space defined by LST and α [13,14] and/or the space defined by LST and f vg [8,14,15].Then, the vertices of the obtained polygons are connected by straight lines [16] or curves [17] to interpolate the dry and wet boundaries over the full range of vegetation cover.Once the LST boundaries have been determined, for any data point in the LST − α or in the LST − f vg space, EF is derived as the ratio of the distance separating the point from the line identified as the dry edge to the distance separating the dry edge and the line identified as the wet edge.In general, both submodels are fully independent, meaning that the different Tends algorithms can be implemented using the same EF submodel.The contextual nature of EF-based models is hence attributed to the Tends submodel.
EF-based ET models are very attractive for operational applications.Due to their image-based nature, they are easily transferable to different areas.In particular, they do not require calibrating the, often site-specific, surface resistances of residual-based methods [18].The downside is the applicability of the image-based Tends submodel, which requires: (i) uniform atmospheric conditions; (ii) a relatively flat area; and above all; (iii) that the temperature endmembers be actually observed within the study domain at the thermal sensor resolution [15,[19][20][21][22][23].In a nutshell, image-based Tends submodels should be quite well adapted for semi-arid irrigated areas provided that the inherent high heterogeneity be resolved at the observation scale, meaning that high-spatial resolution data are used [24].In contrast, image-based Tends submodels are likely to be inadequate for less heterogeneous rainfed areas and when using medium to low (∼1 km) resolution thermal remote sensing data.Long and Singh offer an in-depth study of how the determination of the image-based Tends impact the EF retrievals in terms of performance and error propagation [21].
To remove some of the above limitations, attempts have been made to de-contextualize Tends submodels, notably by developing Tends algorithms independent of the (contextual) surface conditions within the study area and, hence, of the spatial resolution of input LST data [15,16].Tends are thus modeled by running a surface energy balance model forced by meteorological data in both dry and wet conditions and for soil and vegetation components separately.This involves implementing aerodynamic resistances for four extreme surfaces: water-saturated bare soil, dry bare soil, well-watered vegetation, fully-stressed vegetation.Note that the four modeled extreme conditions can be reduced to two-fully dry bare soil and water-stressed vegetation-by setting the cold edge to air temperature for all vegetation covers [15,18].Whereas the parameterization of bare soil aerodynamic resistance is relatively straightforward, the computation of canopy resistance is more complex, as it involves additional parameters to determine the zero plane displacement, such as vegetation height [15,16,25].
In this context, this study aims to develop a simple theoretical approach to estimate Tends independently from LST images.Soil Tends are simulated by a soil energy balance model, and two formulations of aerodynamic resistance are tested over bare soil.Vegetation Tends are obtained from soil Tends and air temperature.Model-derived soil Tends are first evaluated with the in situ measurements made in an irrigated area in Morocco.They are next compared to the soil Tends retrieved from 90-m resolution ASTER data collected over two irrigated areas in Spain and Mexico using an algorithm that combines both LST − α and LST − f vg spaces.Model-derived and image-based Tends are then used as input to a contextual (SEB-1S, [14]) ET model, and results in terms of 90-m resolution and 1-km resolution ET estimates are evaluated over the Mexican site, where ground-based ET measurements are available.
The three study sites in Morocco, Spain and Mexico and their associated experimental data are described in Section 3. In Section 4, the SEB-1S modeling approach is reiterated, and the Tends models forced by either meteorological data or remote sensing data are presented.In Section 5, Tends models are intercompared, and the ET estimates are evaluated at 90-m and 1-km resolutions.

Sites and Data Description
The study is based on data collected over three irrigated semi-arid areas: an irrigation perimeter situated in the Haouz plain in central Morocco (hereby referred to as "R3"), the Urgell area in Catalonia, in northeastern Spain, and the Yaqui valley (Sonora State) in northwestern Mexico.A visual representation of the three study sites is provided in Figure 1.

R3 Perimeter
The Haouz plain in central Morocco covers around 6000 km 2 .The climate is semi-arid, typically Mediterranean, with an average annual precipitation of about 250 mm.About 85% of available water is used for agriculture.The R3 irrigated perimeter, situated in the Haouz plain, includes mainly wheat, olive trees, maize and beet.The typical size of crops is 3 up to 4 ha.Wheat is generally sown between mid-November and mid-January, depending on climatic conditions and the start of the rainfall season.Harvest usually occurs about 5 up to 6 months later, in May or June.More details about the description of the R3 site can be found in [26,27].A field experiment was conducted in a maize parcel of the R3 perimeter on five days: 14 April, 22 April, 30 April, 8 May and 16 May 2014.Maize was seeded in early April, so that the crop field was mostly bare soil during the whole experiment.The experiment focused on measuring the minimum and maximum soil Tends (T s,wet and T s,dry ) to serve as a reference for the direct validation of Tends modeling approaches using meteorological data.Weather conditions were mostly sunny, with one cloudy day (16 May) and one post-rainy day (22 April).The different atmospheric conditions met throughout the campaign ensured testing the model in various conditions.The temperature measurements were taken over two small soil plots (2 m by 2 m).One soil plot was kept fully dry, whereas the soil moisture of the other soil plot was maintained equal to the soil moisture at saturation.The dry soil conditions were fulfilled by isolating the respective soil plot with a plastic sheeting, buried at a 15-cm depth to prevent water infiltrations underneath from capillary rising.The dry plot was also covered by plastic sheeting between two consecutive sampling days to protect the soil from possible rainfall events before the experiments took place.On each sampling day, the saturated soil plot was manually and continuously irrigated from about 10-20 min before the start of the temperature measurements and until the time of the last measurement.Much attention was given to maintain the soil plot as wet: irrigation water was poured in sufficient quantity to reach soil saturation, but as slow as possible to avoid flooding and runoff.Figure 2 shows an illustration of the two 2 m by 2 m soil plots where the experiment was undertaken.Note that the soil surface structure was not preserved during the re-filling of the soil pits, but the soil surface was made relatively flat to mimic the surrounding (bare soil) corn field.Continuous measurements of T s,wet and T s,dry were taken between 10:00 a.m. and 11:30 a.m.UTC (Coordinated Universal Time).The observation time is compatible with the overpass times of ASTER, Landsat-7/Landsat-8 and MODIS thermal sensors.The soil temperature for dry and wet soils were monitored using DS1921G-F5 Thermochron iButton sensors.These sensors are rugged, low-cost and self-sufficient systems that record and store time and temperature.They are equipped with a wide (from −40-85 • C) temperature range and fit for a moderately long (from 1-255 min) logging rate.Their accuracy spans from ± 1 • C in the −30-70 • C range to ± 1.3 • C outside this range, with a field resolution of 0.5 • C ( [28]).In order to increase the confidence in the experimental data and to have a better estimate of the uncertainty in measurements, several Thermochron iButtons were set up at each soil plot.In practice, iButtons were fixed on 7 cm by 7 cm acrylic glass plates, for a set of three iButtons per plate.Because more variability was found in the dry soil measurements, two plates (six iButtons) were used to measure T s,dry , whereas one plate was used to estimate T s,wet .The plates were buried very close to the soil surface, by making sure that iButtons were not directly exposed to solar radiation.Thermochron data were logged at a 5-min time interval.Entries registered by all Thermochron iButtons were averaged per soil plot, and per 30-min bin, cumulated at 10:30 a.m., 11:00 a.m. and 11:30 a.m.UTC.This resulted in a total of three sets (T s,wet and T s,dry ) of values per day.
In addition to soil temperature measurements, meteorological data, such as air temperature, relative air humidity, solar radiation and wind speed, were acquired every half hour by a permanent meteorological station located in an alfalfa crop field about 200 m from the dry/wet soil plots.The acquisition height of the meteorological measurements was 2 m.
Figure 3 illustrates the mean (the average of Thermochron iButton entries) soil temperature measurements for wet and dry conditions (T s,wet and T s,dry , respectively) obtained on 8 May.The reason why the wet and dry temperatures are very similar before 10:00 a.m. is that the sensors were set up right before and removed right after each field experiment.They were not measuring soil temperatures corresponding to wet and dry conditions, but they were stored in the same media (air) and, thus, logging similar temperatures.Soil Tends measurements perfectly match until 30 min before the start of the experiment, which reflects excellent confidence in the iButtons' embedded calibration.By computing the average on each soil plot of the standard deviation of instantaneous measurements during the experiment duration, the uncertainty in T s,wet and T s,dry is estimated as 0.90 • C and 1.1 • C for wet and dry conditions, respectively.

Urgell Site
The Urgell area presents a semi-arid continental Mediterranean climate, with a mean yearly air temperature of 16 • C, precipitation of 400 mm and a number of days with rain of 60.The area of interest includes both irrigated crops consisting of wheat, maize, alfalfa and various fruit (apple and pear) trees, as well as dryland crops, such as barley, olive trees, vineyards and almond trees [29].The Urgell experiment, conducted in 2011, focused on a 20 km by 20 km area, situated within the Urgell irrigation region.Meteorological data are continuously acquired over the area.The meteorological data were collected by the automatic weather station of Golmes (WC) as part of the XEMA network operated by the Meteorological Service of Catalonia (SMC).Air temperature, relative air humidity, solar radiation and wind speed were recorded every 30 min at a 2-m height by the Golmes (WC) station.

Yaqui Site
The Yaqui valley in northwestern Mexico presents a semi-arid climate, with an average annual potential evapotranspiration of 2233 mm and an annual precipitation of 290 mm [30].The main crops include wheat (about 50%), broccoli, beans, chili pepper, potatoes, chickpea, safflower, orange and corn.Irrigation takes up to 90% of the total water expenditure.
The experiment was performed from December 2007-May 2008 on an irrigated area situated in the Yaqui valley.Measurements were focused on a 4 km by 4 km area, while the study area in this paper is defined as a 16 km by 10 km area containing the 4 km by 4 km Yaqui experimental area.Seven micro-meteorological stations equipped with an eddy-covariance flux measurement system were set up in different fields [14].The six stations, providing data concurrently, with at least four ASTER overpasses, are used in the ET validation analysis.The latent heat flux was measured with KH20 fast response hygrometers at a frequency of 10 Hz and converted to the 30-min average.Meteorological data, including air temperature, relative air humidity, solar radiation and wind speed, were acquired every half hour, from 27 December 2007-17 May 2008 from a station located in the center of the region.The acquisition height was 10 m.Further details on the experiment can be found in [14,[30][31][32].

ASTER
ASTER is an imaging instrument on-board NASA's Terra satellite, launched in 1999 and having an equator crossing time at 11 a.m., with a 16-day revisit cycle.
The ASTER thermal sensor has five thermal bands, with a spatial resolution of 90 m.ASTER official surface kinetic temperature products (AST_08) were downloaded from the Earth Observing System Data Gateway and extracted over the 16  In order to test the contextual ET models at multiple resolutions, ASTER data are also used at a 1-km (MODIS thermal) resolution by aggregating the original 90-m resolution ASTER temperature using a simple average [31].

Formosat-2
The Formosat-2 satellite is an Earth observation satellite of the National Space Organization of Taiwan, launched in 2004.Its high (8 m) resolution images are provided daily at 9:30 a.m.
equator-crossing time.The reason behind using Formosat-2 derived albedo resides in the unusability of the ASTER shortwave infrared data on its seven overpass days over the Yaqui site [30].Its high temporal resolution (one day revisit cycle) enabled us to use data collected on the nearest date from each of the seven ASTER overpass dates in order to estimate both NDVI and α from the red and near-infrared reflectances, aggregated at ASTER thermal sensor resolution.

LST, Fractional Green Vegetation Cover, NDVI and Surface Albedo
In order to retrieve the LST, the "temperature and emissivity separation" algorithm [33,34] was applied.
Fractional green (photosynthetically active) vegetation cover and surface albedo are derived over the Urgell and Yaqui study areas from satellite surface reflectances in the red and near-infrared bands.For the Urgell area, reflectance data are obtained from 15-m resolution ASTER red and near-infrared data on its two overpass days.For the Yaqui area, reflectance data are obtained from Formosat-2 data, available on dates closest to the ASTER overpass dates.In both cases, surface reflectances are aggregated from the observation resolution to ASTER thermal resolution (90 m) beforehand.
The Gutman and Ignatov (1998) formula is applied in order to obtain the fractional green vegetation cover (f vg ): with NDVI vg corresponding to fully-covering green vegetation and NDVI s to bare soil or to bare soil partially covered by senescent (non-photosynthetically active) vegetation.In this study, NDVI vg and NDVI s are set to the maximum and minimum value of the NDVI observed during the agricultural season within the study domain.NDVI endmembers are estimated as 0.93 and 0.018, and 0.93 and 0.18 for the Urgell and Yaqui areas, respectively.The reason behind the small value registered for NDVI s for Urgell is the presence of open water (reservoir of irrigation water) in the region.NDVI is computed as the ratio of the difference between (re-sampled) near-infrared and red reflectances to their sum.Surface albedo (α) is estimated as a weighted sum of (re-sampled) red and near-infrared reflectances, using the coefficients provided by [35] and corroborated by [36] and in [30].

Models
The main goal of this study is to assess the impact of different modeling approaches of Tends on the ET estimated at multiple resolutions by a contextual model.Herein, the contextual model used for the estimation of ET is the recently-developed SEB-1S model.Like a number of contextual models, SEB-1S computes ET as EF times the available energy, with EF being estimated as the ratio of the distance separating a pixel in the LST − α space from the line defined as the dry edge to the distance separating the dry edge and the line identified as the wet edge.The dry and wet edges are defined from Tends.
Two different approaches are investigated to estimate Tends using (1) available remotely-sensed data (referred to as Tends_RS) and ( 2) a soil energy balance model forced by available meteorological data (referred to as Tends_EBsoil).The energy balance model offers the possibility to derive Tends estimates independently of surface conditions and of the spatial resolution of remote sensing data.

SEB-1S: An EF Model (EF_SEB-1S)
SEB-1S relies on the same LST − α space as the one used in the Simplified Surface Energy Balance Index (S-SEBI, [13]).It computes EF similarly to S-SEBI, as the ratio of the distance separating the point from the dry edge to the distance separating the wet and dry edges [37][38][39][40][41][42].However, SEB-1S provides a different interpretation of the wet edge, which has proven to significantly improve ET estimates over the Yaqui area [14].Figure 4 presents a graphic illustration of the computation methodology applied in SEB-1S.ET is defined as EF times the available energy (Rn − G), where Rn (W•m −2 ) is the surface net radiation and G (W•m −2 ) the ground heat flux.Further details describing the methodology can be found in [14].

Image-Based Tends Algorithm (Tends_RS)
The Tends_RS approach provides soil and vegetation Tends from a synergy between the LST − α space and the LST − f vg space, as in the original version of the SEB-1S Tends submodule.Note that previous studies [8,15,19,20,43] use only one space in the determination of the extreme temperatures.Figure 5 illustrates side-by-side the SEB-1S interpretation of the edges and vertices of the polygon in the LST − α space and the polygon in the LST − f vg space for data collected on 11 April, over the Yaqui site.
In the LST − α space, SEB-1S identifies four Tends and three α endmembers (α_ends), which correspond to bare soil (α s ), green vegetation (α vg ) and senescent vegetation (α vs ) albedo.The edges are interpreted as follows: [AB] as "bare soil", [BC] as "wet surface", [CD] as "full-cover vegetation" and [DA] as "dry surface".In the LST − f vg space, SEB-1S identifies the same four Tends and both f vg endmembers (0 and 1).The edges are interpreted as follows: [AB] as "mixed soil and senescent vegetation", [BC] as "wet surface", [CD] as "full-cover green vegetation" and [DA] as "dry surface".The dry edge [AD] has a negative slope in both S-SEBI [13] and SEB-1S [14] models due to a radiative and convective control in dry conditions, respectively.In SEB-1S, the wet edge BC identified in the LST − α space has a negative slope, also, due to a convective control in wet conditions [14], consistent with the common interpretation of the LST − f vg space.As in the original version of the SEB-1S Tends submodule, the minimum vegetation temperature T v,wet is set to the air temperature T a .The maximum soil temperature T s,dry is the same for both spaces and is set to the maximum surface temperature T max .The other Tends are estimated by combining the Tends retrieved from the LST − α and LST − f vg spaces separately.The minimum soil temperature for each space is computed by extrapolating the wet edge passing through the well-watered full-cover green vegetation index until the bare soil line.The maximum vegetation temperature is computed such that all data points are located in between the dry [AD] and the wet [BC] edges.This process is first undertaken for each space separately and then sets T v,dry to the average of the two retrieved values.Please refer to Appendix A for further details concerning the original version of the algorithm and the adjustments made to it.

Model-Derived Tends (Tends_EBsoil)
A different methodology is proposed to estimate soil and vegetation Tends independently of remotely-sensed data.A soil energy balance model forced by available meteorological data is used to derive soil Tends.Two simple assumptions are then made for deriving vegetation Tends from model-derived soil Tends.Alternative to the use of a complete surface energy balance model for both soil and vegetation components, the energy balance is modeled herein for bare soil conditions only (and not for full-cover vegetation).The idea is to develop a meteorological-based approach for Tends determination using a minimum number of a priori input parameters.

Soil Energy Balance Model
The maximum and minimum soil temperatures are modeled by setting the near-surface soil moisture to zero and to the soil moisture at saturation respectively, as an input to an energy balance model for bare soils (EBsoil).In practice, EBsoil operates an iterative loop on the soil temperature T s (initialized to T a ), until thermal equilibrium is established [9]: with H s (W•m −2 ) being the soil sensible heat flux and LE s (W•m −2 ) the soil latent heat flux.The latent heat flux is calculated as: with ρ (kg•m −3 ) being the air density, C p (J•kg −1 •K −1 ) the specific heat of air at constant pressure, γ (Pa•K −1 ) the psychrometric constant, e sat (T s ) (Pa) the saturated vapor pressure at soil temperature, e sat (T a ) (Pa) the saturated vapor pressure at air temperature, r ss (s•m −1 ) the soil evaporation resistance and r ah (s•m −1 ) the aerodynamic resistance to heat transfer.r ss is given by the formulation in [44], as a function of soil moisture: where the two best fit parameters A (unitless) and B (unitless) are considered as 8 and 5, respectively [45,46].SM represents the surface (0-5 cm) soil moisture, whereas SM f c the soil moisture at field capacity.Given the constraint of Equation ( 2), EBsoil computes T s,dry and T s,wet by setting SM to 0 and SM sat , respectively.

Aerodynamic Resistance Formulations for Bare Soil
The aerodynamic resistance that intervenes in the computation of the sensible heat flux, as well the latent heat flux (see Appendix B) plays an essential role in the prediction of soil Tends.Different parameterizations have been proposed based on the Monin-Obukhov (MO) similarity theory [47][48][49][50].These parameterizations fall within three different categories: those that follow the MO similarity theory, those that are based on empirical methods and, finally, the so-called semi-empirical ones [51].Different formulations were developed with the purpose of saving computation time, by introducing various assumptions and simplifications in regards to the stability effects, which perform differently when comparing estimates to measurements of r ah [51].
Under stable conditions, an exact solution of the aerodynamic resistance can be easily obtained due to the linearity of the functions of the stability parameters.However, under non-stable conditions, these functions are highly non-linear, and an iterative algorithm is used in order to obtain the exact solution.This method is known as the "MO similarity" method.
To assess the impact of r ah on soil Tends, two different formulations are implemented within EBsoil.The first formulation RI we consider in our implementation takes into account the Richardson number R i (unitless), which represents the importance of natural relative to forced convection and is expressed as in [48]: with η (unitless) being a coefficient equal to 0.75 in unstable conditions (T s > T a ) and to 2 in stable conditions (T s < T a ).We define a simple formulation, based on an "empirical" method, which represents the aerodynamic resistance that neglects natural convection: with Z 0m (m) the roughness length for momentum transfer over bare soil, Z r (m) the reference height at which the wind speed is measured, k (unitless) the von Kármán constant and u a (m•s −2 ) the wind speed.
The bulk Richardson number is a measure of the influence of atmospheric stability on the flux-gradient relationship in the surface layer: with g (m•s −2 ) the gravitational constant.This method is referred herein to as the "empirical" method.
The second formulation takes into account the MO length, and is computed as: with ψ h (unitless) the stability correction factor for heat transport and u * (m•s −1 ) the friction velocity.
For further details on the derivation of r ah from the MO theory, please refer to Appendix C.

Vegetation Tends
Two assumptions are made for deriving vegetation Tends from model-derived soil Tends.The minimum vegetation temperature T v,wet corresponding to well-watered unstressed vegetation is set to the air temperature [52][53][54].This hypothesis is common to many models, such as SEBAL (Surface Energy Balance Algorithm for Land) [54], S-SEBI [13] and SEB-1S [14].The work in [55] considers T v,wet as a fraction of the air temperature, using a correction factor determined as 0.993.This correction factor was introduced after observing a systematic bias between T s and T a of well-watered vegetation.
The maximum vegetation temperature T v,dry corresponding to fully water-stressed vegetation is estimated as: by assuming that the dry and wet edges in the LST − f vg space are parallel.This is based on the observation that the lines [AD] and [BC] in the LST − f vg space are practically parallel, especially when using high-resolution remote sensing data in highly heterogeneous conditions [14].From a physical point of view, it is assumed that the water-saturated bare soil to full-cover unstressed vegetation aerodynamic resistance difference is approximately the same as the bone dry bare soil to full-cover water stressed vegetation aerodynamic resistance difference.In other words, the impact of water stress on vegetation height is neglected, so that the decreasing rate in surface temperature with increasing vegetation cover (and the related increase in roughness length for momentum transfer from bare soil to full-cover vegetation) is practically similar for both dry and wet surface edges.Note that by setting the minimum vegetation temperature to T a , the roughness length for momentum transfer over full-cover unstressed vegetation is implicitly constrained.Equation 9 represents a similar constraint on the roughness length over full-cover water-stressed vegetation (for T v,dry estimation), which is fully consistent with the T v,min = T a strategy.
Deriving vegetation Tends from soil Tends and T a data provides the advantage of not requiring the vegetation height information, as well as the zero plane displacement, the roughness length for momentum transfer and the roughness length for heat transfer.The proposed approach only relies on the roughness length for momentum transfer (Z 0m ) over bare soils, which can be uniformly set to 0.001 m [56].

Validation Strategy
Two validation strategies are applied, in terms of both Tends and ET estimates.Figure 6 offers a schematic diagram of the models used, including corresponding inputs and outputs, alongside the performed validations.First, the soil Tends simulated by EBsoil (soil Tends_EBsoil) are compared to the soil Tends measured in situ (soil Tends_IS).Second, the soil Tends retrieved from remote sensing data at 90-m resolution (Tends_RS_90m) are compared to the soil Tends_EBsoil.Third, the soil Tends retrieved from remote sensing data at 1-km resolution (Tends_RS_1km) are compared to the soil Tends_EBsoil.
Concerning ET estimates, first of all, a validation at 90-m resolution with in situ measurements is presented.SEB-1S is then applied to 1-km resolution data, which were derived from the aggregation (simple average) of 90-m resolution data.This 1-km resolution modeled ET is compared against the ET simulated by SEB-1S at the original 90-m ASTER resolution and subsequently aggregated to 1 km.

Results and Discussion
In this section, the soil Tends simulated by the soil energy balance model (Tends_EBsoil) using the RI or the MO r ah formulation are evaluated against: (1) the in situ measurements (Tends_IS) in the R3 area; and (2) the soil Tends retrieved from 90-m resolution ASTER data (Tends_RS_90m) in the Yaqui and Urgell areas.The model-derived and image-based Tends are then used separately as input to SEB-1S when applied to 90-m resolution and to 1-km resolution (aggregated) ASTER data.The impact of Tends on contextual ET estimates is finally discussed in terms of surface conditions and observation resolution.

Direct Validation Using In Situ Measurements
When comparing the two sets of results, it is clear that the latter gives the best results in terms of root mean square error (RMSE), bias, slope of linear regression and correlation coefficient (see Table 1).The mean RMSE is equal to 2.4 • C for the MO formulation, compared to 3.6 • C for the RI formulation.The lowest mean bias is also obtained for the MO formulation, −1.6 • C, as opposed to −2.2 • C for the RI formulation.The mean correlation is increased from 0.95 (RI) up to 0.97 (MO formulation, respectively).
It is found (not shown) that T s,dry accounts for more bias than T s,wet , and its correlation with in situ measurements is slightly weaker than for T s,wet .This is possibly due to the T s,dry measurements, which present greater variability than the T s,wet measurements, as explained previously.
Table 1.Root mean square error (RMSE), bias, correlation coefficient (R) and slope of the linear regression between soil Tends simulated by the energy balance model for bare soils (EBsoil) with the RI and MO aerodynamic resistance and the soil Tends observed in situ at the soil plots of the R3 perimeter.As the calculation of Tends depends on r ah , which depends on wind speed, a sensitivity study of simulated r ah with respect to wind speed is also conducted (see Figure 8).It appears that a higher sensitivity is obtained for wind speeds inferior to 2 m•s −1 and that T s,dry is impacted slightly more than T s,wet .Both RI and MO formulations are similar in terms of T s,wet sensibility to wind speed.However, this is not the case for T s,dry , as the RI formulation is more sensitive to wind speeds inferior to 1 m•s −1 .As wind speed decreases from 2 m•s −1 , the resistance to heat transfer becomes very large and is no longer representative of the turbulent processes occurring.Results are consistent with those obtained in [57], stating that one of the parameters contributing to important uncertainties in the estimation of aerodynamic resistances is the wind speed.
As a summary, the soil Tends simulated using the MO formulation compare quite well to in situ measurements and seem to be relatively stable even in the conditions of low wind speeds.The uncertainty in modeled T s,dry is larger than in modeled T s,wet , consistent with a significantly larger range of variation of T s,dry .

Consistency with Remotely-Sensed Soil Tends
The proposed two r ah formulations are now implemented over the Urgell and Yaqui areas.The soil Tends issued from the physical aerodynamic resistance model are compared to the soil Tends derived by the SEB-1S Tends submodel from 90-m resolution ASTER data.Figure 9 plots model-derived (Tends_EBsoil) versus image-based (Tends_RS_90m) soil Tends for the Yaqui and Urgell areas separately.
Statistical results in terms of daily bias in T s,wet and T s,dry are provided in Table 2, while values of the slope of the linear regression, correlation coefficient and RMSD between modeled and remotely-sensed Tends can be found in Table 3, for each r ah formulation.When analyzing the results obtained over the Yaqui and Urgell areas combined, the MO r ah formulation gives the best overall statistical results, with a RMSD estimated as 4.1 • C, compared to 6.4 • C for the RI formulation.The consistency obtained between model-derived and image-based soil Tends reflects the robustness and precision of the image-based Tends_RS_90m over highly heterogeneous irrigated areas.Such a consistency is quite satisfying, as it validates both the physical approach based on the MO theory and the remote sensing Tends algorithm based on both LST − α and LST − f vg spaces.The approach is also implemented using 1-km resolution (aggregated) ASTER data over the two sites.Figure 10 shows side by side the 90-m resolution (grey dots) and the 1-km resolution (black dots) LST − α and LST − f vg spaces plotted for three ASTER overpass dates over the Yaqui site (30 December,11 April and 6 May, respectively).The polygons obtained using Tends_EBsoil (RI and MO formulations), Tends_RS_90m and Tends_RS_1km are overlaid on each space.The temporal variability in the boundaries of the two spaces can be explained by the presence or absence of vegetation and by the change in vegetation color during the senescence period, as noted in [14].Nonetheless, one observes that the contour of 90-m resolution data points is relatively well represented by the Tends_RS_90m polygons.The other three polygons (Tends_RS_1km and Tends_EBsoil for RI and MO formulations) are significantly different, especially in the LST − α space.The albedo endmembers for Tends_EBsoil and Tends_RS_1km are derived from 1-km resolution ASTER data.This is the reason why the three polygons are narrower and closer to the 1-km resolution data points.However, the MO r ah formulation provides dry and wet edges closer to the Tends_RS_90m polygon than the RI r ah formulation (consistent with a better performance in terms of soil Tends).Moreover, the MO r ah significantly outperforms the Tends_RS_1km algorithm, which systematically underestimates the dry edge and overestimates the wet edge.Consequently, one concludes that the major factor that contributes to the estimation of physics-based Tends is the aerodynamic resistance and that the MO formulation appears to be a better candidate than the RI formulation, especially for the dry edge determination.As the spatial resolution of observation decreases, so does the number of data points within both LST − α and LST − f vg spaces.A small number of data points represents a strong constraint on the image-based determination of Tends.In fact, a larger domain size would be needed at 1-km resolution to capture more variability in surface (especially dry and wet) conditions.However, the domain size is another strong constraint on the contextual Tends algorithms, because the application of image-based methods implicitly assumes that meteorological conditions are uniform within the study domain.In particular, Tends are by definition zero-dimensional, that is Tends are not scale-dependent.Therefore, the relevance of image-based methods over a given application domain is a compromise between an extent sufficiently large to meet the heterogeneous requirement (so that dry and wet conditions can really be observed) and an extent sufficiently small to meet the uniform atmospheric forcing requirement (so that Tends are uniform, too).In this paper, image-based algorithms were applied to ∼10 km-sized areas to limit the variability of atmospheric conditions within the study domain.
Table 3. Root mean square difference (RMSD), correlation coefficient (R) and slope of the linear regression between the soil Tends simulated by EBsoil for the RI and MO aerodynamic resistance retrieved from 1-m resolution data (RS_1km) and the soil Tends retrieved from 90-m resolution data, for the Yaqui and Urgell regions combined.

Application to ET Estimation at 90-m and 1-km Resolutions
One key advantage of model-derived Tends lies within the possibility to "de-contextualize" Tends submodels, meaning to extend the applicability of the so-called contextual ET models to remote sensing data at multiple spatial resolutions and to regions less heterogeneous than irrigated areas.In this study, the performance of the Tends algorithms is assessed by applying SEB-1S to both 90-m resolution and 1-km resolution data, using the model-derived or image-based Tends as input.The issue of the scale mismatch between remotely-sensed and ground-based ET estimates [12] is addressed by the following stepwise approach: 90-m resolution SEB-1S ET is evaluated with eddy covariance measurements, and the 1-km resolution SEB-1S ET is compared to the 90-m SEB-1S ET aggregated at 1-km resolution.
Figure 11 plots the ET simulated by SEB-1S at 90-m resolution versus the ground-based observations of the six stations of the Yaqui experiment.SEB-1S ET estimates are obtained using the available energy (Rn − G) observed by the flux stations and either the model-derived (RI or MO r ah formulation) or the image-based Tends as input.The differences between the three simulation results (see the distribution of points in Figure 11) indicate that ET is greatly influenced by EF, which strongly depends on Tends.Moreover, the r ah modeling approach for computing Tends has a significant impact on ET estimates.The RMSD between modeled and ground-based ET is 65, 82 and 120 W•m −2 for Tends_RS_90m, Tends_EBsoil with the MO r ah formulation and Tends_EBsoil with the RI r ah formulation, respectively.Tends_RS_90m-based ET estimates are in good agreement with in situ measurements, especially given the wide observation range of ET (0-700 W•m −2 ) [14].The comparison between SEB-1S and in situ ET is also undertaken using the available energy (Rn − G) modeled by SEB-1S.In this case, the RMSD between modeled and ground-based ET is 85, 106 and 130 W•m −2 for Tends_RS_90m, Tends_EBsoil with the MO r ah formulation and Tends_EBsoil with the RI r ah formulation, respectively.The RMSD (65-85 W•m −2 ) of ET obtained when using image-based Tends as input is still in the range reported by the published literature [58][59][60][61].In the case of the Yaqui experiment, RMSD values tend to be high, because most components of the energy budget are often large in semi-arid lands at low latitudes.In addition, a fairly strong bias is present on net radiation, which tends to worsen ET estimations [30].Note that the significant difference in RMSD when using image-based versus model-derived Tends as input to SEB-1S could be attributed to uncertainties in the vegetation Tends (simply derived from soil Tends and air temperature).Nonetheless, the MO r ah still provides acceptable ET estimates at 90-m resolution.
Figure 12 plots the ET simulated by SEB-1S at 1-km resolution (from aggregated ASTER data) versus the 1-km resolution ET chosen as a reference over the Yaqui area, on all seven ASTER overpass days.Given the good results obtained with Tends_RS_90m in terms of 90-m resolution ET estimates, the reference dataset for evaluating the 1-km resolution ET models is obtained by aggregating the 90-m resolution Tends_RS_90m-forced SEB-1S ET at 1-km resolution.The Tends used as input to SEB-1S are either modeled using Tends_EBsoil with the RI or MO r ah formulation, or retrieved from (1-km resolution) remote sensing data using Tends_RS_1km.Statistical results, including RMSD, bias, slope of the linear regression and correlation coefficient, between simulated, and reference ET is reported in Table 4, for each Tends algorithm and for each observation day separately.The MO formulation yet again gives the best results, with a mean RMSD of 56 W•m −2 compared to 77 W•m −2 for the RI formulation.The mean RMSD obtained when using Tends_RS_1km is close to that obtained when using the RI-derived Tends, which is 78 W•m −2 .However, the lowest mean bias (15 W•m −2 ) in ET is obtained for Tends_RS_1km compared to 23 W•m −2 and 27 W•m −2 for Tends_EBsoil with the MO and RI formulations, respectively.The mean correlation coefficient is increased from 0.95 for Tends_RS_1km up to 0.98 for both Tends_EBsoil algorithms.The ET simulated at 1-km resolution is plotted against the 90-m resolution simulated ET, aggregated at 1-km resolution, for data over the Yaqui region, on the seven ASTER overpass dates separately.The Tends used as input to the EF_SEB-1S model consist of (a) the Tends simulated by EBsoil with the RI and (b) MO aerodynamic resistance, (c) the Tends retrieved from 1-km resolution data, (d) the Tends simulated with the mixed modeling approach using αends_1km and (e) the Tends simulated with the mixed modeling approach using αends_90m.Table 4. Root mean square difference (RMSD), bias, correlation coefficient (R) and slope of the linear regression between the ET simulated at 1-km resolution and the 90-m resolution simulated ET aggregated at 1-km resolution, for data over the Yaqui region, on the seven ASTER overpass dates separately.The Tends used as input to EF_SEB-1S consist of Tends simulated by EBsoil with the RI and MO aerodynamic resistance and the Tends retrieved from 1-km resolution data (RS).At 1-km resolution, the Tends_EBsoil algorithm with MO formulation improves ET predictions.However, the presence of a significant bias can be observed on the ET estimates computed using the MO-derived Tends, on 6 May.In fact, Tends_EBsoil are computed independently of the LST dataset and, thus, can lead to the polygon not being centered on the cluster points in the LST − α and LST − f vg spaces.This can introduce a bias in model-derived Tends and, as a result, in EF/ET.In contrast, the polygon obtained from Tends_RS_1km is centered on the cluster points, because the computation of Tends is made using the observed data points.Overall, the ET estimates derived from Tends_RS_1km have a low bias.However, the mean slope of the linear regression between modeled and reference ET is estimated as 1.7.This is due to the fact that Tends_RS_1km overestimates the wet edge and underestimates the dry edge, which translates into the polygon boundaries being too close to the point cluster.The data points are thus closer to the boundaries at 1-km resolution than at 90-m resolution, influencing the range of variability and leading to slopes significantly larger than one.SEB-1S, as a contextual model, performs well, provided that the extreme conditions associated with fully dry and wet soil/vegetation components are encountered in the same image.Many studies have documented the sensitivity of contextual models to the wet and dry edges.For instance, Long et al. [19] indicated that SEBAL shows great sensitivity to extreme soil Tends, which implies that the extreme conditions play a significant role in the estimation of ET.The dry and wet edges impact the ET estimates in a similar manner: increasing Tends for the wet/dry extremes tend to increase ET based on energy conservation.Yet the distance between the observed dry and wet edges depends the range of EF, which is linked to the growth stage of vegetation and to the spatial representativeness of the study area [8].It also depends on the domain size.The impact of the domain size on ET estimates can reach 50-80 W•m −2 [19].When applied to moderate or low spatial resolution images, contextual models could fail to discriminate extreme soil wetness conditions and, hence, to detect EF, especially for relatively small study sites [20].Conversely, deriving the physical boundaries independently of remote sensing data could help improve ET estimates, especially by reducing the underestimation and overestimation of the dry and wet edge, respectively.The theoretical approach developed herein, based on a simple soil energy balance model and the MO r rah formulation, is well adapted for large-scale applications using moderate to low resolution data.

A Mixed Modeling Remote Sensing Approach for Improved Tends and ET estimates
As the Tends_EBsoil with the MO formulation and the Tends_RS_1km provide complementary results in terms of bias and slope of the linear regression between simulated and reference ET, a mixed modeling approach is proposed for T s,dry .The new formulation considers T s,dry as the maximum between the two estimates: Furthermore, a new constraint concerning the albedo endmembers is applied.It consists of using either the 1-km resolution α endmembers (αends_1km) or the 90-m resolution α endmembers (αends_90m) in the estimation of Tends_RS_1km.A comparison between the 1-km resolution mixed-modeled and reference ET is also presented in Figure 12, for the Yaqui area, on all seven ASTER overpass days, and using αends_1km (e) and αends_90m (f).Results in terms of RMSD, bias, correlation coefficient and slope of the linear regression between simulated and reference ET are reported in Table 5, for each of the two mixed Tends formulations (MX_αends_1km and MX_αends_90m, respectively) and for the two previous configurations (MO-derived Tends_EBsoil and Tends_RS_1km, respectively).The statistical results are presented for all seven ASTER overpass dates combined.Both mixed modeling formulations significantly improve the results, with RMSD values of 43 W•m −2 and 52 W•m −2 (for MX_αends_90m and MX_αends_1km, respectively), as compared to 59 W•m −2 and 79 W•m −2 (for the MO_1km and RS_1km formulation, respectively).Table 5. Root mean square difference (RMSD), bias, correlation coefficient (R) and slope of the linear regression between the ET simulated at 1-km resolution and the 90-m resolution simulated ET aggregated at 1-km resolution, for data over the Yaqui region, on all of the seven ASTER overpass dates combined.The Tends used as input to EF_SEB-1S consist of Tends simulated by EBsoil with the MO aerodynamic resistance, retrieved from 1-km resolution data (RS), simulated with the mixed modeling approach using αends_1km (MX_αends_1km) and with αends_90m (MX_αends_90m).As expected, when using the αends_90m as input to the Tends_RS_1km, the computed boundaries are much closer to the boundaries estimated at 90-m resolution.This results in the lowest bias obtained between the mixed-modeled ET at 1-km resolution (using MX_αends_90m derived Tends) and the 1-km resolution reference ET.Nevertheless, the improvement provided by the αends_90m configuration is relatively small, so that the MX_αends_1km strategy seems to be a good compromise in terms of ET accuracy and data availability.Especially the MX_αends_1km does not require high (90 m) resolution remote sensing data for calibrating αends.

Conclusions
This paper develops a simple theoretical approach to estimate Tends independently from LST images, for use as input to contextual ET models.The soil Tends simulated by a soil energy balance model using the RI and MO r ah formulations are evaluated with in situ measurements made over an irrigated area in Morocco.Results indicate that the MO formulation provides the best soil Tends estimates, with a mean RMSE of 2.4 • C, compared to 3.6 • C for the RI formulation.Model-derived soil Tends are also compared to the image-based soil Tends derived from the SEB-1S Tends submodel, which combines both LST − α and LST − f vg spaces as additional constraints on wet and dry edges.The approach is tested using 90-m resolution ASTER data over two semi-arid irrigated areas in Mexico and Spain.The MO formulation still provides the best results, with an RMSD between simulated and image-based soil Tends estimated as 4.1 • C and 6.4 • C for the MO and RI formulations, respectively.
The performance of Tends algorithms in terms of ET estimation is assessed by applying the SEB-1S EF submodel to both 90-m resolution and 1-km resolution (aggregated) ASTER data using the model-derived or image-based Tends as input.In the case of the theoretical Tends algorithm, vegetation Tends are estimated from model-derived soil Tends and air temperature.The RMSD between 90-m resolution SEB-1S and in situ ET is estimated as 65 and 82 W•m −2 , and the RMSD between 1-km resolution SEB-1S and aggregated 90-m resolution SEB-1S ET is estimated as 78 and 56 W•m −2 for the image-based and model-derived Tends (MO formulation), respectively.In light of the above results, one concludes that Tends should be estimated a priori when contextual models are applied to low resolution images.Nonetheless, image-based Tends provide accurate ET estimates when high-resolution images are used in highly heterogeneous (irrigated) landscapes.In these conditions, the consistency between model-derived and high-resolution image-based Tends provides a meaningful basis for developing mixed modeling observational approaches.
Regarding the estimation of ET using contextual models with image-and physics-based edge determination, this study differs from previous studies in four main aspects: (1) Tends are derived from a r ah formulation for bare soil only, although most theoretical Tends algorithms are based on r ah formulations for both soil and vegetation [10,19,43]; (2) to the authors' knowledge, this is the first time that soil Tends have been evaluated using in situ measurements; (3) the image-based edge determination relies on both the LST − α and LST − f vg spaces, while previous studies have been based on one space only [8,15,19,20,43]; and (4) the sensitivity of ET simulations to observation resolution is assessed in terms of image-based and physics-based edge determination.
The successful implementation of a soil energy balance model in contextual ET models increases the accuracy in ET estimates at multiple resolutions and has many other potential applications.First, the consistency between model-derived and image-based Tends could lead to improving the parameterization of aerodynamic resistance by remote sensing means.Secondly, the minimum and maximum soil temperatures are the boundary conditions for LST-based soil evaporative efficiency (SEE, being defined as the ratio of actual to potential evaporation) models [62,63].Improving the accuracy in soil Tends opens the path for monitoring SEE from optical/thermal remote sensing.Thirdly, LST-derived SEE is tightly coupled to the near-surface soil moisture retrieved from microwave radiometers, such as SMOS (Soil Moisture and Ocean Salinity), especially under semi-arid non-energy-limited conditions [64].The DISPATCH (DISaggregation based on a Physical And Theoretical scale CHange; [29]) algorithm improves the spatial resolution of SMOS-like data by converting LST-derived SEE into high resolution soil moisture fields.Any improvement in LST-based SEE models would contribute to further developments of DISPATCH and to associated refinements of soil moisture products.Last, but not least, remotely-sensed estimates of LST-derived SEE and microwave-derived soil moisture data would be key for better constraining the partitioning of LST-derived ET into soil evaporation and plant transpiration [65].Such information is crucial for better quantifying crop water needs and managing water resources over semi-arid areas.lower than 90 m.Modifications mainly consist of defining new image-based thresholds to estimate the slopes of the wet and dry edges more accurately using data that are not necessarily representative of extreme (wet/dry) soil and vegetation conditions.Below, we provide a summary of the methodology for the estimation of minimum vegetation temperatures (T v,wet,1 and T v,wet,2 ), minimum soil temperatures (T s,wet,1 and T s,wet,2 ), maximum soil temperatures (T s,dry,1 and T s,dry,2 ) and maximum vegetation temperatures (T v,dry,1 and T v,dry,2 ) for the LST − α and LST − f vg spaces, respectively.With respect to the LST − α space, no modification is made to the computation of the minimum vegetation temperature T v,wet,1 and the maximum soil temperature T s,dry,1 , respectively: • T v,wet,1 (if α = α vg ) is set to the air temperature T a ; • T s,dry,1 (if α = α s ) is set to the maximum temperature (T max ) observed within the study area.
Several modifications are made in the estimation of T s,wet,1 and T v,dry,1 : • The minimum soil temperature T s,wet,1 is defined as the intercept at α = α s of the line passing through the point (α vg , T a ) and the point with α < α th,1 , such that the slope of the line is maximum (meaning that all of the other data points with α < α th,1 are located above the wet surface edge), with α th,1 being the average between α vg and α s .In the original version of SEB-1S, this threshold was set to α vg .• The maximum vegetation temperature T v,dry,1 is defined as the intercept at α = α vs of the line passing through (α s , T s,dry,1 ) and the point with α > α avg , such that the slope of the line is maximum (meaning that all of the other data points with α > α avg are located below the dry surface edge) with α avg being the average of all α values within the study area.
In a similar manner, regarding the LST − f vg space, T v,wet,2 and T v,dry,2 remain unchanged: • T v,wet,2 (if f vg = 1) is set to the air temperature T a ; • T s,dry,2 (if f vg = 0) is set to T max .
Several modifications are made in the estimation of T s,wet,2 and T v,dry,2 : • The minimum soil temperature T s,wet,2 is computed as the intercept (at f vg = 0) of the line passing through the point (1, T a ) and the point with f vg < f vg,avg , such that the slope of the line is maximum (meaning that all of the other data points with f vg < f vg,avg are located above the wet surface edge) with f vg,avg being the average of all f vg values within the study area.In the original version of SEB-1S, the threshold value (f vg,avg ) was set to 0.5.It is now computed for each day separately.• The maximum vegetation temperature T v,dry,2 is defined as the intercept (at f vg = 1) of the line passing through the point (0, T s,dry,2 ) and the point with f vg > f vg,avg , such that the slope of the line is maximum (meaning that all of the other data points with f vg > f vg,avg are located below the dry surface edge).
Finally, as in the original version of SEB-1S, an estimation of the four Tends is given by averaging the two Tends sets obtained for each space separately: T s,dry = T s,dry,1 = T s,dry,2 = T max (11) T s,wet = (T s,wet,1 + T s,wet,2 )/2 T v,wet = T v,wet,1 = T v,wet,2 = T a T v,dry = (T v,dry,1 + T v,dry,2 )/2 B. Energy Balance Model for Bare Soil The net radiation at the soil surface is given by: Where α s (unitless) represents the soil albedo, Rg (W•m −2 ) the incident solar radiation at short wavelengths, s (unitless) the soil emissivity, set to 0.96, Ra (W•m −2 ) the incident thermal radiation at large wavelengths, σ (W•m −2 •K −4 ) the Stefan-Boltzmann constant and T s ( • C) the soil temperature.
The soil heat conduction flux can be approximated as [10,66]: The soil sensible heat flux is computed as: r ah (17) with ρ (k•g•m −3 ) being the air density and C p (J•kg −1 •K −1 ) the specific heat of air at constant pressure.

C. Aerodynamic Resistance Modeling (MO Formulation)
The MO length L mo is computed as follows: u * (ms −1 ) represents the friction velocity and is expressed as: ψ m represents the stability correction factor for momentum transport and is given by: The stability correction factor for heat transport is expressed as: where x is a function of the MO length L mo and of the reference height for wind speed observations: x = (1 − 16 Z r L mo ) 0.25 (22)

Figure 1 .
Figure 1.The three study sites included in this paper: (a) the R3 perimeter, (b) the Urgell area and (c) the Yaqui valley.

Figure 2 .
Figure 2. Illustration of the dry soil plot in the R3 perimeter.(a) The pit was dug at a 15-cm depth and (b) the soil plot was kept fully dry by covering it with a plastic sheeting in between experiments.

Figure 3 .
Figure 3. Mean and standard deviation of the in situ measurements performed on 8 May on the dry and wet soil plots of the R3 perimeter, Haouz plain.The two vertical lines stand for the beginning and the end of the experiment.
km by 10 km Yaqui and 20 km by 20 km Urgell areas.Two cloud-free ASTER images were acquired over the Urgell site at around 11:00 a.m.UTC on 16 August and 3 October 2011.Seven cloud-free ASTER images were obtained over the Yaqui area during the agricultural season of 2007-2008, at around 11:30 a.m.local solar time on 30 December, 23 February, 10 March, 11 April, 27 April, 6 May and 13 May.

Figure 4 .
Figure 4.For a given J point, EF is computed as the ratio of IJ to IK. Underlying grey points correspond to 90-m resolution data acquired on 27 April, over the Yaqui area.

Figure 5 .
Figure 5. (a) Interpretation of the edges and vertices of the LST − α and (b) the LST − f vg polygons.Underlying points correspond to 90-m resolution data acquired on 11 April, over the Yaqui site.

Figure 6 .
Figure 6.Schematic diagram presenting an overview of the models used (including inputs, outputs and the relationships between them) and the validation strategy of temperature endmembers (Tends) and ET estimates.

Figure
Figure 7a,b plot the soil Tends modeled using the RI and MO version of r ah , respectively.

Figure 7 .
Figure 7.The soil Tends simulated by EBsoil for the (a) RI and (b) MO aerodynamic resistance are plotted against the soil Tends observed in situ at the two soil plots in the R3 perimeter.

Figure 8 .
Figure 8. Absolute value of the simulated to observed maximum (top) and minimum (bottom) soil temperature difference as a function of wind speed for the RI and MO aerodynamic resistance separately.

Figure 9 .
Figure 9.The soil Tends simulated by EBsoil for the (a) RI and (b) MO aerodynamic resistance are plotted against the image-based Tends for the Yaqui and Urgell regions separately.

Figure 10 .
Figure 10.The LST − α and LST − f vg spaces corresponding to 1-km resolution data (black dots) and 90-m resolution data (grey dots) are overlaid with the polygons built using Tends_RS_1km (blue), Tends_RS_90m (light blue) and the Tends_EBsoil for the RI (magenta) and MO (red) aerodynamic resistance.Results are shown over the Yaqui area (30 December, 11 April and 13 May, respectively).

Figure 11 .
Figure 11.The ET simulated at 90-m resolution is plotted against station measurements for data over the Yaqui region.The Tends used as input to the EF_SEB-1S model consist of (a) the Tends retrieved from 90-m resolution data, (b) the Tends simulated by EBsoil with the RI and (c) MO aerodynamic resistance.

Figure 12 .
Figure 12.The ET simulated at 1-km resolution is plotted against the 90-m resolution simulated ET, aggregated at 1-km resolution, for data over the Yaqui region, on the seven ASTER overpass dates separately.The Tends used as input to the EF_SEB-1S model consist of (a) the Tends simulated by EBsoil with the RI and (b) MO aerodynamic resistance, (c) the Tends retrieved from 1-km resolution data, (d) the Tends simulated with the mixed modeling approach using αends_1km and (e) the Tends simulated with the mixed modeling approach using αends_90m.

Table 2 .
Difference between the soil Tends simulated by EBsoil with the RI and MO aerodynamic resistance and the soil Tends retrieved from 90-m resolution data for the minimum and maximum soil Tends, for the Yaqui and Urgell regions separately.
Threshold albedo, computed as the average between α vg and α s Stability correction factor for momentum transport r ah Aerodynamic resistance to heat transfer (s•m −1 ) r ah,RI Aerodynamic resistance to heat transfer based on the Richardson number (s•m −1 ) r ah,M O Aerodynamic resistance to heat transfer based on the Monin-Obukhov length (s•m −1 )