J Cuxart
email: joan.cuxart@uib.cat

A A Boone

Evapotranspiration over land from a Boundary-Layer Meteorology perspective

Keywords: Evapotranspiration, Eddy-covariance systems, Lysimeters, Heterogeneous surfaces, Land-surface models

The precise determination of Evapotranspiration (ET) is challenging because it is a quantity which is difficult to measure and to parametrize.

Direct estimates include the determination of the change of mass of a volume of soil and vegetation that evapotranspirates using lysimeters, or direct measurements of turbulent water vapour fluxes by eddy-covariance systems.

Parametrized estimates that make use of the Monin-Obukhov similarity theory use vertical gradient measurements of temperature and moisture at one point, and line or area averages by means of scintillometers operating at high frequency.

A relation for ET from well-watered surfaces was initially developed by Penman and later expanded for vegetated surfaces and for heterogeneous croplands. A popular simplified expression was obtained by Priestley and Taylor.

The current challenge is to find expressions for ET in non-saturated conditions, which are common in arid and semi-arid climates, and for heterogeneous terrain.

In numerical models, the estimated actual ET over land is obtained as the result of the explicit representation of the different involved sub-processes taking place in the soil and the canopy, using so-called land-surface models. Usually these mechanisms are described in a simplified manner and rely on a number of adjustable parameters. The improvement of such descriptions relies in the availability of experimental measurements to make the physical models more complete and robust.

Introduction

Evapotranspiration as a key process in the Earth system

The evapotranspiration rate (ET), is the flux of water vapour from the surface to the atmosphere, traditionally expressed in kg s -1 m -2 or an equivalent unit such as mm day -1 or mm year -1 , in which 1 mm represents 1 kg m -2 . We will denote it through the variable E as the sum of the evaporation from the soil (E g ) and the transpiration of plants (E tr ),

E = E g + E tr
(1)

Water evaporation requires energy, denoted the latent heat flux λE, which is the product of the enthalpy of vaporization of water λ and E, expressed in W m -2 , and ET can be found in the literature referred to in both ways.

From an energy point of view, ET is the most important non-radiative process transmitting heat from the surface to the atmosphere globally, and it is larger than direct heating of the air by the surface sensible heat flux [START_REF] Bosilovich | Global energy and water budgets in MERRA[END_REF]. Furthermore it is the branch of the water cycle transporting water upwards feeding clouds and at the origin of precipitation. The process links the energy and water budgets, and it is quantitatively very significant in both from a global perspective [START_REF] Trenberth | Estimates of the global water budget and its annual cycle using observational and model data[END_REF].

There is intense research on ET in the fields of plant physiology [START_REF] Jones | Plants and microclimate: a quantitative approach to environmental plant physiology[END_REF], agronomy [START_REF] Pereira | Crop evapotranspi-ration estimation with FAO56: Past and future[END_REF], soil physics [START_REF] Or | Advances in soil evaporation physics -A review[END_REF], global and regional climate [START_REF] Hagemann | Evaluation of water and energy budgets in regional climate models applied over Europe[END_REF][START_REF] Katul | Evapotranspiration: a process driving mass transport and energy exchange in the soil-plant-atmosphere-climate system[END_REF], remote sensing [START_REF] Zhang | A review of remote sensing based actual evapotranspiration estimation[END_REF] and hydrology [START_REF] Liu | A worldwide evaluation of basin-scale evapotranspiration estimates against the water balance method[END_REF], and for these subjects we refer the reader to the referenced works and citations within. Furthermore, there are other important processes related to the phase changes at the Earth's surface, such as interception of precipitation and the corresponding evaporation [START_REF] Van Dijk | Rainfall interception and the coupled surface water and energy balance[END_REF], sublimation from solid to vapour [START_REF] Bliss | Sublimation and surface energy budget of Taylor Glacier, Antarctica[END_REF], or the reverse process to ET, condensation, taking place mostly in clear and calm nights [START_REF] Jacobs | Towards closing the surface energy budget of a mid-latitude grassland[END_REF][START_REF] Cuxart | Evaluation of the surface energy budget equation with experimental data and the ECMWF model in the Ebro Valley[END_REF].

A quantification of the importance of ET at the global scale is estimated by [START_REF] Miralles | Magnitude and variability of land evaporation and its components at the global scale[END_REF] using satellite information. The biome types for which ET is maximum are low latitude tropical forests (1182 mm year -1 ) and savannas (806 mm year -1 ). Mid-latitude biome types have larger ET in greener areas, such as croplands (542 mm year -1 ), temperate forests (512 mm year -1 ) and grasslands (462 mm year -1 ), than in semi-arid biomes such as shrublands (315 mm year -1 ). High latitude boreal forests evaporate 372 mm year -1 , while deserts (122 mm year -1 ) and permanent snow (57 mm year -1 ) are the areas with less ET. Their maps show that E tr is the dominant process in vegetated areas, while E g prevails in desert areas and is significant in semi-arid climates.

Snow sublimation is important in the northern hemisphere and interception is especially relevant in the forested areas, being maximal in the tropics. [START_REF] Trenberth | Earth's global energy budget[END_REF] estimate a global λE of 80 W m 2 , which is much larger than the average sensible heat flux (H) of 17 W m 2 , and the λE average values over land (39 to 52 W m 2 depending on the source) are also larger than H (25 to 27 W m 2 ). In the current global evolution of temperature and CO 2 concentration in the atmosphere, Sikma et al (2019) indicate that, for vegetated terrain, increases solely in CO 2 concentration would enhance H and decrease λE, whereas increases solely of temperature would have the opposite effect, and the combined increases of both factors would depend on their relative importance.

Considering the atmospheric boundary layer (ABL), ET is constrained by the net radiation (R n ), and is often a major term in the surface energy budget (SEB). When the terrain is dry, the incoming solar energy is employed in warming the soil and the ABL, whereas when it becomes wet the energy used in the phase change results in a lesser warming in both media. The atmosphere is fed by ET with water vapour, thereby increasing the buoyancy of the plumes that may result in ABL clouds which, in turn, enhance the exchange of mass and energy with the free atmosphere. In addition, the amount of water vapour in the ABL and above modulates the terrestrial longwave radiation captured and re-emitted, influencing the surface radiation and energy budgets [START_REF] Edwards | Radiative processes in the stable boundary layer: Part I. Radiative aspects[END_REF]. Dry air from above the ABL is entrained across the convective boundary layer (CBL) capping inversion, and its transport by turbulence to the surface layer increases [START_REF] De Bruin | A model for the Priestley-Taylor parameter?[END_REF][START_REF] Van Heerwaarden | Interactions between dry-air entrainment, surface evaporation and convective boundary-layer development[END_REF].

Physical mechanisms and interactions

The main physical processes that are involved in determining ET at one location are usually considered to be: a) R n , that provides most of the energy for the phase change; b) the turbulence intensity, that favours the renewal of the air close to the surface; c) the water demand of the atmosphere; d) the amount of water in the soil-vegetation system and its ability to reach the interface with the atmosphere; e) the lateral transport of atmospheric moisture related to the surrounding surface heterogeneities; f) large-scale advection of heat and moisture and g) the entrainment of dry air across the ABL top, with the latter three modifying the local atmospheric conditions. The coupling of the water and energy balances at the land surface induces feedbacks that must be considered for a proper representation of the mechanisms in play [START_REF] Brubaker | Analysis of feedback mechanisms in land-atmosphere interaction[END_REF]1996;[START_REF] Kim | Feedbacks in the land-surface and mixedlayer energy budgets[END_REF]. For a review of the history and the state of the art on determining ET some fifteen years ago, the reader is referred to [START_REF] Shuttleworth | Putting the vap into evaporation[END_REF].

As we depart from the classical idealized picture of an homogeneous terrain covered by grass on a clear-sky day, more processes have to be taken into account. This complex system involves a large number of interactions, as described in [START_REF] Ek | Influence of soil moisture on boundary layer cloud development[END_REF] see their Figure 1), and it is difficult to inspect analytically. A substantial part of the progress in understanding ET dynamics is based in the conceptual model of [START_REF] Raupach | Equilibrium evaporation and the convective boundary layer[END_REF], that has inspired large-eddy simulation (LES) studies, especially as they become coupled with land-surface models (LSMs). This coupled framework provides a more complete description of the soil-vegetation processes and their interac- tion to the atmosphere (as in [START_REF] Lohou | Surface energy balance and buoyancy response to shallow cumulus shading[END_REF]. The parametrization of the process depends on the spatio-temporal scales under study, which is discussed in Section 5. This is complemented with new strategies of data acquisition and analysis, often in synergy with numerical studies. The different ABL processes affecting ET over land are shown schematically in Fig. 1.

The effect of entrainment can be formally analyzed using the framework set by [START_REF] Raupach | Equilibrium evaporation and the convective boundary layer[END_REF], who summarized the effect of the external environment through a generic conductance g e . This allowed the formal treatment of lateral advection and top entrainment in what is defined as a thermodynam-

ically partly-open system, which provides a solution for the steady state that depends on the ratio of environmental and surface conductances. A CBL over a well-watered grassland may reach equilibrium evaporation, which happens when the evaporative fraction (taken as λE/(H+λE)) remains approximately constant. This occurs when the atmospheric water demand increases with temperature at a similar rate at which humidity does as the result of surface evaporation and entrainment at the top of the CBL. van Heerwaarden et al.

(2009) determined an expression for the evaporative fraction in the CBL using Raupach's framework by combining three feedback mechanisms: "heating" as ET increases with increasing temperature of the atmosphere, "moistening" as ET is reduced as humidity builds up in the CBL, and "drying" as entrainment brings drier air into the CBL.

As the spatial complexity of the vegetation increases, more detailed analyses are needed. Continuous forest canopies act as porous elevated surfaces with complex turbulent motions above and across the canopy top, and the surface can be decoupled from the above-canopy dynamics if the leaf-area index is large [START_REF] Belcher | The wind in the willows: flows in forest canopies in complex terrain[END_REF]. In the case of patches of forest and open terrain, the terrain heterogeneities induce circulations [START_REF] Garcia-Carreras | Impact of mesoscale vegetation heterogeneities on the dynamical and thermodynamic properties of the planetary boundary layer[END_REF] that may have significant effects on the local values of ET [START_REF] Simó | Observed atmospheric and surface variability on heterogeneous terrain at the hectometre scale and related advective transports[END_REF].

An adequate physical description of the effects of surface heterogeneity is still lacking, although it has been postulated that patches on scales of around a few hundreds of metres may explain part of the organized transport in the lower ABL [START_REF] De Roo | The influence of idealized surface heterogeneity on virtual turbulent flux measurements[END_REF]Mauder, 2018, Mauder et al., 2020). The scale of the surface heterogeneities has an impact on the the entrainment rate at the top For sake of completeness, the main concepts relating soil physics to evaporation are summarized here [START_REF] Or | Advances in soil evaporation physics -A review[END_REF]. Soil evaporation, E g , is usually considered as a two-stage process for a porous material provided that the necessary supply of energy for the phase change is available, and that the atmosphere in contact with the wet surface is not saturated. As long as there is a film region (with pore spaces occupied by air or water) between the surface and the underlying saturated region forming a continuously connected hydraulic network, water will flow from the saturated region to the surface by capillarity. This stage is termed 'constant rate period', as it is controlled by the atmospheric evaporative demand, and it is close to the evaporative rate of a free water surface under the same conditions. The second stage starts as a dry region is formed in the upper part of the soil, driving the film region downward, as the capillary rise is not sufficient to reach the surface. The top of the film region is termed the 'secondary drying front', and evaporation proceeds to the surface by vapour diffusion through the pore spaces in the soil, with a significant drop of E g . The soil evaporation rate depends on the physical properties of the soil, essentially its porosity and texture, as well as the adsorptive capacities of the soil materials [START_REF] Verhoef | Adsorption of water vapor by bare soil in an olive grove in southern Spain[END_REF] that may reduce the amount of vapour reaching the surface. A summary of equations describing the transport of water vapour in the soil can be found in García-González et al ( 2012), and a discussion on the parameterization of these processes is given in Subsection 5.3.

The scope of this paper

In this paper, the current status of the determination of ET over land is presented, with a focus on its experimental measurement, parametrization and representation in numerical models. We restrict ourselves to estimations for a The discussions are intended for a precise location and the effects of terrain heterogeneities will be limited to those related to the immediate surroundings at the scale of the hectometre, since typically lateral transport related to surface variabilities is already relatively small for the kilometre scale [START_REF] Cuxart | Estimation of the advection effects induced by surface heterogeneities in the surface energy budget[END_REF]. Sustained synoptic or mesoscale advections are typically of longer spatial and temporal scales. The estimation of areal values, such as drainage basins, satellite pixels or model grid elements, and their relation to values at a single station or of a network of stations, will not addressed here. However, the use of the approaches presented here within the framework of numerical modelling and the associated large range of space scales will be discussed in Section 5.

It is worth also mentioning the relation between the actual ET and a theoretical concept, the potential evapotranspiration (E p ), which is widely used in applications. It is assumed to represent the amount of water vapour coming from a surface without any restriction on water availability. This concept is a good approximation to ET when the soil-vegetation system is close to waterholding field capacity, which is the common case in mid-latitudes -except in dry summers-, and in the tropics in the wet seasons. Note that when the soil is not saturated, E p may not be a suitable representation of ET, that should be determined either by direct or approximate measurements, by land-surface models (LSMs) or by a specific theory on the actual ET.

In section 2 the current methodologies for the experimental estimation of ET are revisited, going from the most direct determinations to those that are highly parametrized. Section 3 gives a review of the E p formulations which are currently most used, since it is a concept of wide application as long as the soil-vegetation system is charged with water. In Section 4, the challenges of determining the actual ET for non-saturated terrain and for heterogeneous surfaces are discussed, with a succinct description of some recent propositions.

Section 5 summarises LSM-based representations of ET with some information on the current limitations of this approach. Finally, Section 6 contains a description of the perspectives related to ET in ABL-related research.

Experimental determination of evapotranspiration

The experimental determination of ET is not straightforward. When it is considered as a flux of matter from the surface to the atmosphere, eddy-covariance systems are considered to be the reference method. This equipment is expensive and its use is subject to certain constraints [START_REF] Hicks | Measurements of fluxes over landcapabilities, origins, and remaining challenges[END_REF].

A common surrogate method uses vertical gradients of temperature, humidity and wind to estimate the turbulent fluxes using empirical functions depending on the thermal stability. From a surface point of view, ET may be considered as the loss of water mass from a volume of soil and vegetation, that can be determined using lysimeters, which monitor the weight changes. The comparison between both methods must take into account that they may be sampling the ET process at different spatial and temporal scales depending on their configuration.

Direct methods

An eddy-covariance (EC) system measures the wind, the sonic temperature and the water vapour concentration at high frequency (typically 10 or 20 Hz). The wind vector is obtained by a three-dimensional sonic anemometer measuring the time that an ultrasound wave takes between transducers in both directions. In addition, the sound speed can be derived and from it the so-called 'sonic temperature', which is very similar to the virtual temperature.

High-frequency series of water vapour are obtained using an open-path optical system, which emits infrared radiation that is absorbed depending on the water vapour density (Lee et al., 2004;[START_REF] Foken | Micrometeorology[END_REF].

Averages and fluctuations of these atmospheric variables can be obtained, from which turbulence covariances (including vertical fluxes) are computed. If the flow is not stationary, the very meaning of the average is compromised, and if the surroundings are heterogeneous, the conditions required to apply the Reynolds axioms of averaging may not be fulfilled [START_REF] Mauder | Surface Energy Balance Closure over Land: A Review[END_REF].

Therefore the turbulent flux values as provided by an EC system would be valid as a reference only in the conditions described. Water drops on the transducers, from rain or dew, alter the measurements and may give values which do not comply with the quality control criteria, therefore often water vapour fluxes are not provided during these events, even if they seem to be significant [START_REF] Hirschi | A site-level comparison of lysimeter and eddy covariance flux measurements of evapotranspiration[END_REF]. Dust may have a similar effect and regular cleaning and maintenance are a must.

It is assumed that the EC system is located in the surface layer, where the turbulent fluxes are considered to vary little with height, making the values of the fluxes relatively insensitive to the actual position of the system above the surface. This hypothesis is usually not fulfilled in stably stratified nights, when the fluxes may change significantly in the first few metres over the surface [START_REF] Mahrt | Stably stratified atmospheric boundary layers[END_REF]. Besides, nocturnal strong stable stratification may cause evaporation or condensation staying close to the ground below the level of the EC system, questioning the representativeness of the measured fluxes.

Lysimeters provide the changes in mass (assumed to be water changes) of a control volume, allowing to estimate ET if rainfall, seepage and changes in the soil water storage are monitored [START_REF] López-Urrea | Testing evapotranspiration equations using lysimeter observations in a semiarid climate[END_REF]. They are considered a reference system in agronomy which is expensive to setup. With surfaces typically ranging between 0.5 and 10 m 2 and depths between 1 and 2 m, lysimeters can measure small changes of weight using several load cells with high precision [START_REF] Seneviratne | Swiss prealpine Rietholzbach research catchment and lysimeter: 32 year time series and 2003 drought event[END_REF][START_REF] Groh | Determining dew and hoar frost formation for a low mountain range and alpine grassland site by weighable lysimeter[END_REF]. A recent strategy to estimate rainfall is to make high frequency measurements (e.g. 1 minute), assuming that there is no ET when precipitation occurs [START_REF] Schrader | Estimating Precipitation and Actual Evapotranspiration from Precision Lysimeter Measurements[END_REF].

In principle, the depth of the lysimeter should be scaled with the level to which moisture changes at the temporal scale of interest. Therefore, for daily monitoring, the depth could be less than that used for yearly monitoring, as long as most of the root system is included [START_REF] Diaz-Espejo | Illustration of micro-scale advection using grid-pattern mini-lysimeters[END_REF][START_REF] Heusinkveld | An automated microlysimeter to study dew formation and evaporation in arid and semiarid regions[END_REF]. On the other hand, lysimeters may be isolated from their surroundings and neglect lateral runoff, which may be substantial and difficult to estimate, especially in complex terrain [START_REF] Girona | Evapotranspiration and soil water dynamics of peach trees under water deficits[END_REF]. However, it is possible to maintain the lysimeter with similar water table depth to the immediate outside soil by means of a bi-directional pumping system and suction candles, so that the lysimeter water dynamics are according to field dynamics [START_REF] Groh | How to control the lysimeter bottom boundary to investigate the effect of climate change on soil processes?[END_REF][START_REF] Pütz | TERENO-SOILCan: a lysimeter-network in Germany observing soil processes and plant diversity influenced by climate change[END_REF]. The so-called 'oasis effect', related to differences of the lysimeter surface with its immediate surroundings, is also a source of uncertainty that must be considered [START_REF] Gebler | Actual evapotranspiration and precipitation measured by lysimeters: a comparison with eddy covariance and tipping bucket[END_REF].

Comparing EC systems with lysimeters is an interesting exercise, since both approaches are usually taken as reference methods. When the systems are very close to each other, in the same terrain plot, the differences between them depend on the variability of the surface and the surroundings, in particular the footprint effects on the EC system, or to the presence of persistent hectometre-scale motions not well accounted for in the computation of the turbulent water vapour flux [START_REF] Mauder | Surface Energy Balance Closure over Land: A Review[END_REF]. Furthermore, lysimeters provide estimations of ET and condensation in clear and calm nights, whereas EC systems often have difficulties due to the formation of droplets on the transducers [START_REF] Groh | Quantification and prediction of nighttime evapotranspiration for two distinct grassland ecosystems[END_REF]. [START_REF] Gebler | Actual evapotranspiration and precipitation measured by lysimeters: a comparison with eddy covariance and tipping bucket[END_REF] found, for a well-watered grass site, that lysimeter-based ET monthly-averages are very close to E p estimations and that the annual averaged value obtained from the EC system is 5% smaller than the lysimeter value. [START_REF] Hirschi | A site-level comparison of lysimeter and eddy covariance flux measurements of evapotranspiration[END_REF] showed acceptable correspondance between both methods for measurements in nearby locations for monthly averages, whereas large dispersion was found when comparing hourly values.

Plant physiological measurements are intended to quantify transpiration from a plant. For example, gas chambers measure the amount of water vapour released at the leaf level [START_REF] Pérez-Priego | Analysing uncertainties in the calculation of fluxes using wholeplant chambers: random and systematic errors[END_REF]. Plants in transparent isolated containers with controlled air inflow and outflow along with monitoring the other components of the water balance, allow a quantification of the amount of water transpired by the plant physiological activity [START_REF] Escalona | Responses of leaf night transpiration to drought stress in Vitis vinifera L[END_REF]. These methods can be used to calibrate sap flow measurements determining its speed in the xylem by heating at one point and detecting the heated flow at another point [START_REF] Dragoni | Transpiration of apple trees in a humid climate using heat pulse sap flow gauges calibrated with whole-canopy gas exchange chambers[END_REF], which depend on the plant type. Extensively employed in plant physiology studies, these techniques may be of interest for ABL studies in the case of wide dense crop fields, for which most of ET is transpiration. They are more difficult to apply for natural vegetation, since it is usually a mixture of several plant species. A recent attempt to combine these methods with standard micrometeorological measurements is the CloudRoots field experiment (Vilà-Guerau de [START_REF] Vilà-Guerau De Arellano | CloudRoots: Integration of advanced instrumental techniques and process modelling of sub-hourly and sub-kilometre land-atmosphere interactions[END_REF], that shows that information at the leaf level is necessary to obtain accurate parameters for the mechanistic representation of photosynthesis and stomatal aperture, and that sun-induced fluorescence data can be used to estimate the spatial variability of ET.

Methods based on the Similarity Theory

Estimating ET when EC systems or lysimeters are not available is a necessity. Traditionally, Monin-Obukhov Similarity Theory (MOST, [START_REF] Monin | Basic laws of turbulent mixing in the surface layer of the atmosphere[END_REF]Foken 2006) is used, which relates turbulent fluxes and vertical gradients through an eddy diffusivity coefficient including a correction by a function which depends on the thermal stability. These empirical formulations for the surface layer were firstly provided by [START_REF] Businger | Fluxprofile relationships in the atmospheric surface layer[END_REF], and they are currently used in the revised form of [START_REF] Högström | Non-dimensional wind and temperature profiles in the atmospheric surface layer: A re-evaluation[END_REF].

Monin-Obukhov theory assumes that surface conditions are horizontally homogeneous and stationary. It defines dimensionless gradients that are a function of a stability parameter, such as the Obukhov length, L, or the Richardson number, Ri. In the particular case of computing E using L, this relation can be written as

E = ρ a w q = -ρκzu * ∂q ∂z 1 φ e (z/L) (2)
where u * is the friction velocity, ρ a the density of air and φ e (z/L) the fluxgradient universal function. The overline stands for the average values for the selected interval to compute the statistics of interest. It can be different between moisture and other variables (momentum or temperature) in unstably stratified conditions (L < 0), but they are usually taken as identical in stable conditions (L > 0). The classical expressions, which are identical to those for the temperature, according to [START_REF] Högström | Non-dimensional wind and temperature profiles in the atmospheric surface layer: A re-evaluation[END_REF] are

φ e (z/L) = 0.95 1 -11.6 z L -1 2 -2 < z/L < 0 (3) φ e (z/L) = 0.95 + 7.8 z L 0 < z/L < 1 (4)
which have a range of validity restricted to moderate values of stability and unstable stratification (Foken, 2006).

An alternative way to express the fluxes is using the inverse of the eddy diffusivity (called conductance in this framework, see [START_REF] Jones | Plants and microclimate: a quantitative approach to environmental plant physiology[END_REF] as an analog to the electrical resistance, usually taking differences, as

E = -ρ a q(z 2 ) -q(z 1 ) r ae (5)
where r ae is called the aerodynamic resistance for water vapour (in s m -1 ),

with z 1 and z 2 the levels of measurement of q. The integration of Eq. 2 yields

q(z 2 ) -q(z 1 ) = -w q 1 κu * [ln( z 2 z 1 ) -Φ e ( z 2 L ) + Φ e ( z 1 L )] (6)
in which the terms in large brackets are r ae , while Φ e are the integrated forms of Eq. 3 and Eq. 4 as given by Paulson ( 1970)

Φ e (z/L) = 2 ln 1 + x 2 2 z/L ≤ 0 (7) Φ e (z/L) = -5 z L z/L ≥ 0 (8) in which x = (1 -16z/L) 1/4
. Similar equations can be written for the momentum and the sensible heat fluxes [START_REF] Moene | Transport in the atmospherevegetation-soil continuum[END_REF]. This approach is simple to implement and makes use of an iterative method with an arbitrary initial value of L until convergence is reached, since L is a function of the three aforementioned turbulent fluxes.

In terms of data availability, the surface layer variables at the screen level and estimations for the ground level are taken. The wind is considered zero at the height of the roughness length z 0 , whereas the land-surface temperature and a related moisture value are used for the gradients of temperature and humidity. The surface temperature determined from radiation measurements has uncertainties of a few degrees [START_REF] Simó | Impact of the surface-atmosphere variables on the relation between air and Land Surface Temperatures[END_REF] and strong hypotheses are imposed to obtain a value of humidity at the surface, especially when the soil is not saturated. Therefore, the computation of the turbulent fluxes is the surface layer is likely less accurate using estimated data at the surface than it would be actually using two levels of measurement.

The similarity theory has been recently adapted for flows over forests with homogeneous canopies, since the presence of large coherent structures above the canopy top represent a large part of the turbulent transport [START_REF] Finnigan | Turbulence structure above a vegetation canopy[END_REF]. This "roughness sublayer" has a length scale for the turbulent mixing approximately constant (not increasing with height as in the surface layer) and wind profiles which are different than those over the ground. The interaction with the canopy is made through the prescription of a penetration depth scale that depends on the leaf area index and the Stanton number that quantifies the difference between scalar and momentum transport efficiency and the subsequent modification of the MOST expressions [START_REF] Harman | The role of roughness sublayer dynamics within surface exchange schemes[END_REF].

Scintillometers provide an alternative experimental estimation of ET also relying on MOST. A combination of microwave and optical signals travelling

between an emitter and a receiver located in the surface layer allows the determination of the humidity and temperature structure parameters [START_REF] Meijninger | Determination of area-averaged water vapour fluxes with large aperture and radio wave scintillometers over a heterogeneous surface-Flevoland field experiment[END_REF][START_REF] Ward | A critical revision of the estimation of the latent heat flux from two-wavelength scintillometry[END_REF]. These quantities, combined with MOST (Kooijmans and Hartogensis, 2016), provide estimations of the temperature and moisture fluxes over heterogeneous terrain. The 3D spatial distribution of water vapour can be estimated using Raman Lidars, obtaining ET from the profiles combined with MOST [START_REF] Eichinger | Spatial and temporal variation in evapotranspiration using Raman lidar[END_REF][START_REF] Wulfmeyer | Can water vapour Raman lidar resolve profiles of turbulent variables in the convective boundary layer?[END_REF].

The elevated cost of these remote-sensing instruments so far is not conducive for a widespread use of them.

Evapotranspiration as the residual of the surface energy budget

Estimating a term of the SEB as the residual of the sum of the others is a common technique. In numerical models, the SEB closure is imposed as a constraint and therefore it is a legitimate approach to determine a term of the budget as the residual, as in [START_REF] Noilhan | A simple parameterization of land surface processes for meteorological models[END_REF], where the ground flux G is equivalent to the net flux into or out of the combined soil-vegetationsurface scheme.

Experimental estimations of λE by remote sensing from satellite or aircraft widely use this approach. R n at the surface is determined with an atmospheric radiative transfer model, together with estimations of surface temperature and emissivity [START_REF] Roerink | S-SEBI: A simple remote sensing algorithm to estimate the surface energy balance[END_REF]. The ground flux is usually given as a function of R n and some surface parameters, while H is estimated using the radiatively-derived surface temperature, the air temperature from some external source and a computation of the aerodynamic resistance. The large uncertainties involved in these estimations make residual estimates of λE not very precise.

To reduce the uncertainties, experimental validation campaigns are needed for calibration, using λE experimental estimations [START_REF] Jia | Validation of remotely sensed evapotranspiration over the Hai River Basin, China[END_REF] or surface temperature over homogeneous terrain [START_REF] Coll | Ground measurements for the validation of land surface temperatures derived from AATSR and MODIS data[END_REF].

From the point of view of the ABL experimentalist, the lack of closure of the SEB is an unavoidable fact which has to be dealt with [START_REF] Mauder | Surface Energy Balance Closure over Land: A Review[END_REF].

Usually the sum of the turbulent latent and sensible heat fluxes is smaller than the so called 'available energy' R n -G. The experimental estimation with MOST using temperature and humidity at the screen level and at the surface is the most widely used technique because of its simplicity and low cost. However, the strong hypotheses imposed on the surface and the range of validity of the empirical formulas usually reduce the applicability of this approach. Experimental estimations of ET in heteroge- neous terrain, such as those using scintillometers and Raman Lidars, rely fully on MOST and therefore share the same limitations.

3 Estimations of evapotranspiration for well-watered surfaces

The basic equations

As mentioned in the Introduction, the potential evapotranspiration E p (or λE p in terms of energy) is a good approximation to the actual ET when there is no limitation on the water availability for the evaporating surface. This is the case when the soil has enough moisture to allow water to reach the surface or when it is covered by vegetation and the root system has access to sub-surface water, even if it the upper soil is not saturated, as long as stomata are open.

A number of expressions have been derived to estimate E p . The simplest ones correlate E p with solar irradiation, or with temperature as its surrogate. [START_REF] Thornthwaite | An approach toward a rational classification of climate[END_REF] computed an annual heat index using the mean monthly temperatures and then he estimated the monthly E p using the average temperature of the month weighting by the number of daylight hours. A member of this class of approaches, the equation of Hargreaves [START_REF] Hargreaves | History and evaluation of Hargreaves evapotranspiration equation[END_REF] is recommended by the Food and Agricultural Organization when only basic meteorological data are available and a daily E p is needed. Its computation requires the daily average of air temperature T d and the daily temperature range, the global daily solar radiation, expressed in terms of the radiation at the top of the atmosphere R o , and an adjustable coefficient c H , while λ is taken at 20 • :

λE p = c H R o (T d + 17.8) T max -T min (9)
If the focus is put on the representation of the relevant processes, the traditional approach has made use of the SEB concept, from which expressions such as those of Penman-Monteith, Priestley-Taylor or Makkink-DeBruin are derived.

The basic idea was given by [START_REF] Penman | Natural evaporation from open water, bare soil and grass[END_REF], who wrote a four-term SEB equation consisting in R n , G, H and λE p closing the budget, taking R n -G as the available energy and expressing H and λE p in resistance form as in Eq. 5 using the surface as the bottom level. To isolate λE p from this equation he needed to i) assume a saturated surface so that the vapour pressure there is taken as saturated, i.e. e s = e sat (T s ), ii) impose a local linear variation of e sat with temperature allowing to estimate e sat (T s ) from e sat (T a ), an approximation that holds for small vertical temperature gradients, iii) take the aerodynamic resistances equal for heat and moisture, iv) expect that a value for the available energy was at hand. The formula is usually written as

λE p = s (R n -G) s + γ + ρ C p r ae [e sat (T a ) -e a ] (s + γ) (10)
in which s is the slope of the saturated vapour pressure curve at the air temperature T a , C p is the specific heat capacity of air at constant pressure,

γ = Rv R d
Cp λ p is the psychrometric constant, ρ is the air density, e a is the vapour pressure, r ae is the aerodynamic resistance for heat and moisture, p the air pressure and R d , R v the specific gas constants for dry air and water vapour.

The first term is customarily called the "radiation term" as it represents the effect of the radiation available for evapotranspiration. The second term is noted "aerodynamic" as it depends on the state of the atmospheric surface layer, namely the turbulence intensity and the atmospheric water vapour deficit at screen level (Fig. 3).

Recently

The treatment of vegetation

As indicated in [START_REF] Moene | Transport in the atmospherevegetation-soil continuum[END_REF], the expansion for fully-vegetated surfaces was proposed independently by [START_REF] Monteith | Evaporation and environment[END_REF] and by [START_REF] Rijtema | An analysis of actual evapotranspiration[END_REF]. The big-leaf approach represents the actual vegetation as one leaf with one stomatal cavity, for which the water phase change takes place in saturated conditions, by adding an extra canopy resistance r c to represent vapour flowing to the leaf surface which is in series with the aerodynamic resistance. The resulting equation for transpiration is identical to Penman's except that γ is substituted by γ * = γ(1 + rc rae ).

A widely used expression for the canopy resistance [START_REF] Jarvis | The interpretation of the variations in leaf water potential and stomatal conductance found in canopies in the field[END_REF][START_REF] Stewart | Modelling surface conductance of pine forest[END_REF]) takes it proportional to a minimum value depending on the vegetation type and to four factors that are, respectively, functions of the incident solar radiation, the soil moisture, the atmospheric water vapour deficit and the temperature of the air, and inversely proportional to the leaf area index. A more elaborated approach makes use of plant physiology concepts, known as the 'A-g s ' model [START_REF] Farquhar | A biochemical model of photosynthetic CO[END_REF][START_REF] Goudriaan | A simple and fast numerical method for the computation of daily totals of crop photosynthesis[END_REF][START_REF] Ball | A model predicting stomatal conductance and its contribution to the control of photosynthesis under different environmental conditions[END_REF][START_REF] Jacobs | Stom-atal behaviour and photosynthetic rate of unstressed grapevines in semiarid conditions[END_REF]. There the net assimilation of carbon A, taken proportional to difference of carbon concentration between the air in contact with the leaf and the intercellular space, is related to the stomatal conductance g s , the inverse of the stomatal resistance, and this value is then related to r c . Another The canopy resistance is considered minimal when water availability and solar radiation are plentiful and the atmospheric demand is not very high.

Typical r c values range from 30 s m -1 for short and well-watered vegetation to more than 100 s m -1 for tall vegetation. When bare and vegetated surfaces coexist, [START_REF] Shuttleworth | Evaporation from sparse crops: an energy combination theory[END_REF] proposed a two-source scheme (or "combination" model) in which the total E p is the sum of the potential evaporation from the soil E gp and the potential transpiration from the vegetation E trp , each obtained using the Penman-Monteith equation, as

E p = C s E gp + C c E trp (11)
in which the multiplying coefficients are combinations of the aerodynamic, soil and canopy resistances. This approach can be expanded to more combinations of sources and layers as described in [START_REF] Zhang | A review of remote sensing based actual evapotranspiration estimation[END_REF].

Variations of the Penman-Monteith equation

In the case of a non-turbulent saturated surface layer over a well-watered surface, only the radiation term of the Penman equation would remain, an expression which is commonly called "equilibrium evaporation". Nevertheless, turbulent mixing usually removes the saturated air from near the surface and brings drier air from above, allowing the evaporation rate to increase as the aerodynamic term stays active. [START_REF] Priestley | On the assessment of surface heat flux and evaporation using large-scale parameters[END_REF] took the aerodynamic term as proportional to the radiation term for a well-watered surface and expressed λE p as proportional to the equilibrium evaporation,

λE p = α s (R n -G) s + γ (12)
and they proposed an estimated value of 1.26 for α, in the absence of significant advective contributions, implying that the aerodynamic term is about onefourth of the radiation term.

In Another simplified expression conceptually derived from Penman's equation is found in the Makkink-de Bruin approach for the daily evaporation. In 1957, Makkink showed that for a well-watered surface in the Netherlands, the daily λE p scaled with the incoming shortwave solar radiation as

λE p = a s s + γ S d + b (13)
where S d the downwelling global shortwave radiation at the surface, a=0-9

and b= 30 W m -2 . This expression is consistent with the Priestley-Taylor equation in the sense that G is approximately zero averaged over a daily cycle and that in wet surfaces R n is roughly proportional to the global radiation. De [START_REF] De Bruin | From Penman to Makkink[END_REF] showed that the formula could be simplified by taking a=0.65 and b=0, as later proven by de Bruin and Lablans (1998) when comparing it to Penman's.

An alternative expression proposes that b=20 W m -2 (noted β) as an expression of the dry entrainment at the top of the ABL (De Bruin and Holtslag 1982), with a=1 and using the solar radiation S d as a surrogate of R n as (de [START_REF] De Bruin | A thermodynamically based model for actual evapotranspiration of an extensive grass field close to FAO Reference, suitable for remote sensing application[END_REF]:

λE p = s s + γ (1 -A)S d -C s S d S dta + β (14)
in which A indicates the albedo of the surface, S dta the incident shortwave radiation at the top of the atmosphere and C s is an adjustable parameter taken as 110 W m -2 . This formula is currently applied to estimate daily λE p from satellite-derived fields [START_REF] Trigo | Validation of reference evapotranspiration from Meteosat Second Generation (MSG) observations[END_REF]. 4 Estimations of evapotranspiration in general conditions

Non-saturated surfaces and heterogeneous terrain

When the surface is not well-watered, E p may not be a suitable approximation to the actual value of ET, since it will be limited by the availability of water in the soil. Furthermore, if the terrain is heterogeneous, ET will vary between adjacent parcels and lateral transport of humidity will take place by turbulent motions or semi-permanent small-scale circulations between areas at different temperature, a process which is not included in the derivation of the E p equation.

The strategy to extend the E p Penman-type equations for non-saturated surfaces usually faces the challenge of prescribing the surface moisture for the partially-dry upper soil (the film and the dry regions) and incorporate the response of plant transpiration in these conditions [START_REF] Verhoef | Modeling plant transpiration under limited soil water: Comparison of different plant and soil hydraulic parameterizations and preliminary implications for their use in land surface models[END_REF]. [START_REF] Trugman | Soil moisture stress as a major driver of carbon cycle uncertainty[END_REF] show that the empirical modifications of surface conductivity due to changes in soil moisture cause large uncertainties which are the major cause of inter-model variability in the determination of the carbon cycle, indicating that progress is still needed.

A pragmatic approach in applications is to use the Priestley-Taylor equation with the α coefficient adjusted to the type of terrain. [START_REF] Flint | Use of the Priestley-Taylor evaporation equation for soil water limited conditions in a small forest clearcut[END_REF] related α to soil moisture to determine the actual ET for a specific set of atmospheric conditions. Jiang and Islam (1999) adjusted α to a linear variation between 1.26 and 0 for well-watered and dry conditions respectively for remote sensing applications.

Over heterogeneous terrain, the soil moisture distribution is uneven depending to a large extent on the soil composition and vegetation cover of the different parcels. [START_REF] Simó | Observed atmospheric and surface variability on heterogeneous terrain at the hectometre scale and related advective transports[END_REF] have shown that the variability at the hectometre scale can be very large in terms of the vertical gradient of temperature and moisture in the surface layer as well as of the upper soil water contents.

Furthermore they found that in dry conditions, the lateral advection of water vapour can be comparable to the ET at one point, changing the characteristics of the atmosphere in that location compared to an homogeneous site and correspondingly ET (Fig. 4). Eichinger et al. (1996) explored the estimated variability of ET over agricultural fields with a Raman lidar and they found large variability at the decametre and hectometre scales.

Fig. 4 The effect of terrain heterogeneity on ET: humidity advection from parcel "1" changes the air humidity over parcel "2" and subsequently its evaporation rate and the upper soil moisture.

Finally, for estimating water needs for agricultural applications, the concept of a reference evaporation, E 0 , is often used. It is defined as E p for an ideal surface of well-watered (non-water stressed) grass. It is then multiplied by a non-dimensional so-called crop coefficient, commonly referred to as K c , which represents the ratio of the true evapotranspiration to the reference value [START_REF] Allen | Crop evapotranspiration: guidelines for computing crop water requirements, Irrigation and Drainage Paper 56[END_REF]. It takes into account the physiological plant properties and the effect of water stress. It is empirically adjusted and can range from zero to a value slightly above unity (to account for crops which are transpiring at values above the reference grass rate). It can also be split into two parts, the dual crop coefficient, which account separately for evaporation (using a water stress coefficient) and transpiration (using a basal crop coefficient which accounts for plant physiology and is multiplied by a plant stress coefficient).

Such coefficients are used quite extensively at the field or plot scale, but satellite data can be used to adjust the seasonality of K c for spatially distributed estimates. Note that such methods are adapted for estimating daily ET, and thus are not used directly for studying the link between the surface vapour flux and ABL dynamics, thus we will not further develop this concept in the current study.

Alternative approaches to determine evapotranspiration

There exist methods to deal with the ET determination problem using approaches not following the Penman rationale and its implicit limitations. Here three of them are presented: i) the advection-aridity model, which obtains ET using E p and the pan evaporation E pa , without imposing either MOST or SEB closure, and implicitly incorporating advection; ii) the maximum production entropy theory that, given R n , estimates ET using only the temperature and humidity of the surface and, iii) the surface-flux equilibrium expression for daily ET, that needs solely the daily-averaged screen values of air temperature and humidity.

Advection-aridity model. Also known as the 'complementary principle', it provides estimations of ET without imposing the closure of the SEB and, in some applications, it can even go around MOST. This concept was first given by [START_REF] Bouchet | Evapotranspiration réelle et potentielle, signification climatique[END_REF], and then it was developed further by [START_REF] Brutsaert | An advection-aridity approach to estimate actual regional evapotranspiration[END_REF], being formulated recently in a comprehensive theoretical manner by [START_REF] Brutsaert | A generalized complementary principle with physical constraints for land-surface evaporation[END_REF].

The theory results from the relations between three types of ET in the same environment, considered locally uniform: E, E p and E pa , the latter of which is called the 'apparent E p ' which comes from a small surface of water which can be measured with an evaporation pan. When the surface has full availability of water, the three values are identical, otherwise they may be ordered as

E ≤ E p ≤ E pa .
As the surface dries, λE p -λE = ∆Q e , where ∆Q e is the flux of energy not used as evapotranspiration. As [START_REF] Parlange | An advection-aridity evaporation model[END_REF] describe, the energy not employed to evaporate water in a hypothetically well-watered surface is used in a real non-saturated surface for other processes, including sensible heat flux, biological processes or lateral transport in heterogeneous environments.

This idea allows the method to circumvent the need to impose the constraint of SEB closure, and this is the reason why the method is also called the "advection-aridity evaporation model".

This excess of energy can cause E pa to become larger than E p , because the environmental temperature may increase more for a non-saturated surface and therefore the evaporative demand would be larger. This increase above potential conditions (E pa -E p ) was considered equal to ∆Q e by [START_REF] Bouchet | Evapotranspiration réelle et potentielle, signification climatique[END_REF] or proportional to it by [START_REF] Brutsaert | Hydrologic cycle explains the evaporation paradox[END_REF], without a clear justification in any of both cases. [START_REF] Brutsaert | A generalized complementary principle with physical constraints for land-surface evaporation[END_REF], inspired by [START_REF] Han | A nonlinear function approach for the normalized complementary relationship evaporation model[END_REF], proposed a polynomial function for which the coefficients are determined using realistic conditions for moist, dry and drying situations, applicable for arid environments, finding

E = E p E pa 2 (2E pa -E p ) (15)
for which E p can be found using a Penman-type expression and E pa using the value from an evaporation pan or Penman for a surface of water. Brutsaert Maximum production of entropy. This approach is based on nonequilibrium thermodynamics [START_REF] Wang | A model of surface heat fluxes based on the theory of maximum entropy production[END_REF]2011) and provides estimations for the energy fluxes in the surface layer for conditions going from dry to saturated. A so-called "dissipation function", D, is built, which is proportional to the sum of the square of the heat fluxes, each of which is divided by its respective "thermal inertia" coefficient:

D(λE, H, G) ≡ 2 G 2 I s + 2 H 2 I a + 2 (λE) 2 I e ( 16
)
where I s is taken from a simple relation between the diurnal ranges of variation of land-surface temperature and G, I a is derived for H making use of the Monin-Obukhov relations and I e comes from I a = σI e , assuming that both turbulent fluxes go through similar physical processes. σ is a postulated dimensionless parameter, which is formally similar to the inverse of the Bowen ratio: it varies between zero for a dry soil and s γ for a saturated soil. Its general expression depends only on the temperature and moisture of the surface as

σ = λ 2 C p R v q s T 2 s ( 17
)
Minimizing D with the constraint that the SEB closes, expressions of H and λE solely in terms of R n and some physical parameters of the surface are obtained. The general expression for λE is

λE = B(σ)H (18) B(σ) = 6 1 + 11 36 σ -1 (19)
The extension for a fully vegetated surface implies assuming G as negligible providing modified expressions. This model provides a unique solution for λE, H and G given T s , q s and R n , but not the other way around. It is applicable for the entire range of soil wetness, from dry to saturated conditions. The challenge remains in terms of the determination of q s , which the authors suggest to take as the specific humidity of an infinitesimally-thin layer close to the evaporating surface. This theory, proposed for application in the field of hydrology, is formally not restricted to the daily scale, although its validation has been made with daily data.

A surface-flux equilibrium expression for daily evapotranspiration. Following the concept that the diurnal evolution of SL temperature and humidity contains information on the magnitude of the surface heat fluxes if advection is limited [START_REF] Betts | FIFE atmospheric boundary layer budget methods[END_REF], it is hypothesized that the surface moistening and surface heating terms in the near-surface relative humidity budget approximately balance [START_REF] Gentine | Evaporation esti-mates using weather station data and boundary layer theory[END_REF]. Assuming the conceptual framework of [START_REF] Raupach | Equilibrium evaporation and the convective boundary layer[END_REF] for a box, here limited to exchanges across its vertical top, ET can be expressed as [START_REF] Mccoll | Surface flux equilibrium theory explains an empirical estimate of water-limited daily evapotranspiration[END_REF][START_REF] Mccoll | Emergent simplicity of continental evapotranspiration[END_REF])

λE = R n -G 1 + B , B ≈ R v C p T 2 a λ 2 q a ( 20
)
In this expression, (1+B) -1 represents the evaporative fraction, and the Bowen ratio, B, is the result of the balance between the surface moistening and the surface heating in the budget. It is intended for locations without significant lateral transport of moisture, therefore it is best suited for inland continental areas. It depends only on T and q at the screen level, and it implies that ET increases with the relative humidity. This is only correct for time scales of days or longer, since more ET will reflect on more moisture in the atmosphere. Mc Coll et al. (2020) show that the formula provides good results for continental sites using climatological data at the daily scale. A schematic diagram of the methods described in the last two sections is given in Fig. 5.

5 Estimating evapotranspiration using land-surface models

The land-surface modeling approach

Land-surface models (LSMs) were originally implemented in numerical weather prediction (NWP) models to provide interactive lower boundary conditions for atmospheric radiation and turbulence schemes. Their main task there was the computation of the fluxes of heat, mass and momentum between the land and the atmosphere on a timescale compatible with the aforementioned atmospheric processes (i.e. the diurnal cycle must be adequately resolved). The fluxes from the soil and the vegetation are computed through the solution of one or more surface energy budget equations, and arguably ET is the most important flux since it couples the energy and water budgets of the land surface.

In the past two decades, LSMs have evolved considerably to include more processes which are strongly linked with ET in order to meet the growing demands of both the research and the user communities [START_REF] Van Den Hurk | Acceleration of Land Surface Model Development over a Decade of Glass[END_REF], such as photosynthesis and the associated Carbon fluxes, vegetation phenology (biomass evolution, net primary production), and both superficial and root zone soil moisture prediction (for example, to estimate drought severity or crop irrigation needs).

The estimation of ET from an LSM must be approached across scales, where processes are resolved in such models from the local (or the so-called agricultural 'parcel' scale, which is as small as on the order of 10s of m: the scale at which ET measurements are generally made using, for example, eddycovariance techniques, see Section 2), to the mesoscale (from the parcel up to several km), regional (up to 10s of km) to global (areas up to 100s of km).

Of course, the aforementioned discrete scale definitions are somewhat artificial and result either from historical reasons and/or for convenience (adapted for a particular application). The LSM parametrization treats unresolved scaledependent processes as a function of some grid-average state variable though a combination of conceptual models, empirical relationships, and theory. As such, they must account for scale change either statistically, through scaled input parameters, or a combination of both. Therefore, a combination of theory, modeling, analysis and observations (data from field campaigns, re-analysis data-sets and satellite-based data) at multiple scales have been used to develop LSMs.

Generally speaking, the research community has a fundamental knowledge of how to model ET using rather comprehensive approaches. However, progress in LSM development is made seemingly slowly at times in fully coupled or spatially distributed models since such numerical modeling must strike a sometimes delicate balance between computational efficiency, the ability to accurately define input parameters and their related uncertainties, and the complexity of the physical parametrization. Such factors weigh heavily for both operational applications (where tight schedules for getting products to end users must be respected) and climate research (where simulations can take months of CPU time). And of course, certain applications require ET estimates over large domains using dense computational grids (operational hydrological forecasting, e.g.s [START_REF] Habets | The SAFRAN-ISBA-MODCOU hydrometeorological model applied over France[END_REF][START_REF] Snow | A High-Resolution National-Scale Hydrologic Forecast System from a Global Ensemble Land Surface Model[END_REF] for significantly long integration periods (e.g., high resolution land surface reanalysis products: [START_REF] Carrera | The Canadian Land Data Assimilation System (CaLDAS): Description and Synthetic Evaluation Study[END_REF][START_REF] Bonan | An ensemble square root filter for the joint assimilation of surface soil moisture and leaf area index within the Land Data Assimilation System LDAS-Monde: application over the Euro-Mediterranean region[END_REF].

In order to strike a balance, researchers use the strategy of developing a parametrization which only considers what are deemed to be the main processes (or even only those for which we have reasonable estimates of input parameters). Such models are oftentimes developed and bench-marked based on more complex process-based schemes. As our understanding of processes, improved physiographic data (describing both soil and vegetation properties, and land use), computational resources and the availability of observational or satellite data sets all increase, more realistic and complete physically-based ET representations can be developed.

Evapotranspiration in Land-Surface models

Evapotranspiration is computed in LSMs from the surface energy balance equation, and generally (neglecting cold-season processes such sublimation from snow or ice, and intercepted water) ET from the land surface is parametrized as consisting in ground surface evaporation, E g , and transpiration from vegetation, E tr :

E = ρ a r ae [q c -q a ] = E g (T g ) + E tr (T v ) (21)
where T g and T v represent the ground surface and vegetation temperatures, respectively. The specific humidity of the canopy air space and overlying atmosphere at some reference level above the surface are represented by q c and q a , respectively. q c is estimated by substituting the appropriate expressions for E g and E tr into Eq. 21 (which will be presented in the following sub-sections).

Finally, the aerodynamic resistance, r ae , is parametrized using either MOST directly (Eq. 6, as described in Section 2.2) or a bulk Richardson number approach empirically fitted to the aforementioned functions [START_REF] Louis | A parametric model of vertical eddy fluxes in the atmosphere[END_REF].

Early simple LSMs used the composite soil-vegetation single surface energy balance approach (e.g.s Bhumralkar, 1975;[START_REF] Noilhan | A simple parameterization of land surface processes for meteorological models[END_REF] for which a single energy budget was solved thus the surface temperature can be represented as a single temperature (T s = T v = T g ) and r ae defines the aerodynamic resistance between the composite surface and the overlying surface layer. There has been widespread use of such composite schemes to compute ET fluxes in operational NWP for years and they are still in use multiple operational centers for regional to global (e.g.s Belair et al., 2003;[START_REF] Ek | Implementation of Noah land-surface model advances in the NCEP operational mesoscale Eta model[END_REF] and kilometric (e.g. [START_REF] Seity | The AROME-France Convective-Scale Operational Model[END_REF] scale applications, owing in large part, to their relatively small number of required input (physiographic) parameters and easy implementation into operational data assimilation schemes.

In the late 1980s and early 1990s, the so-called big-leaf approach became more prevalent in which two surface energy budgets (also referred to as a two-source approach) are considered, one for the ground surface (often representing evaporation from a relatively thin surface soil layer) and the other for a bulk vegetation canopy. [START_REF] Sellers | The design of a Simple Biosphere model (SiB) for use within general circulation models[END_REF] proposed one of the first highly comprehensive schemes for use in a general circulation models (GCMs) which is still fairly representative of many LSMs in use today. If a two-source energy budget approach is used, then MOST is applied above a certain within-canopy displacement height to model the resistance between the so-called canopy air space and the overlying surface layer (r ae ). Furthermore, functions which account for the effects of vegetation density and height on the below canopy turbulence (eddy diffusivity) are used to compute the aerodynamic resistance between the ground and the vegetation and the canopy air space (Fig. 6).

Many LSMs use approaches inspired by [START_REF] Deardorff | Efficient Prediction of Ground Surface Temperature and Moisture, with Inclusion of a Layer of Vegetation[END_REF] with refinements (e.g. [START_REF] Choudhury | A four?layer model for the heat budget of homogeneous land surfaces[END_REF][START_REF] Raupach | Simplified expressions for vegetation roughness length and zero-plane displacement as functions of canopy height and area index[END_REF]. In recent years, some

GCMs have begun adopting multi-layer vegetation canopy models to compute an integrated value of E tr using explicit within canopy turbulent diffusivity computations as opposed to assumed profile functions [START_REF] Naudts | A vertically discretised canopy description for ORCHIDEE (SVN r2290) and the modifications to the energy, water and carbon fluxes[END_REF] including more recent modifications to better represent the roughness sub-layer [START_REF] Bonan | Modeling canopy-induced turbulence in the Earth system: a unified parameterization of turbulent exchange within plant canopies and the roughness sublayer (CLM-ml v0)[END_REF].

Note that one could also include evaporation terms from water intercepted by the vegetation canopy. Generally speaking, evaporation of canopyintercepted water is a relatively small component of ET, except for certain land cover types such as dense tropical forests where it can be a significant part of the total evapotranspiration from vegetation. It is notoriously difficult to quantify, however some model-based estimates can be found in the literature. For example, Choudhury and DiGirolamo (1998) used a detailed physically-based model which used meteorological measurements as input to derive a global estimate of 20% of the total ET. Since that time, multi-model based LSM estimates have been made which give similar values (e.g. 17% from [START_REF] Dirmeyer | GSWP-2: Multimodel Analysis and Implications for Our Perception of the Land Surface[END_REF]. Finally, intercepted canopy water in irrigated fields can have a significant impact on ET by covering leaf pores or by humidification of the canopy air space, but it's estimation remains complex [START_REF] Jasim | A new approach to estimate canopy evaporation and canopy interception capacity from evapotranspiration and sap flow measurements during and following wetting[END_REF].

Thus there remains considerable uncertainty in how to accurately model this process in LSMs so that most LSMs tend to use rather simple approaches [START_REF] Wang | Evaluation of canopy interception schemes in land surface models[END_REF].

Finally, some LSMs include evaporation of ponded water. But this process is either not represented or parametrized in very different ways among LSMs (perhaps using a different temperature or even methodology, or included in separate wetland or river storage modules): the evaporation from water surfaces such as lakes, marshes or rivers are therefore beyond the scope of the current discussion.

Ground evaporation in Land-Surface models

Ground evaporation is governed by molecular diffusion from the water surface in the soil pore space in some near surface soil layer to the reference level of the ground surface defined as the humidity roughness height, z 0q . As water is lost from this layer, additional moisture can be supplied from below via capillarity. Then, water vapour is exchanged with the overlying atmosphere through laminar and turbulent processes. The process in LSMs is generally represented using one of the three methods described below.

The so-called alpha approach expresses E g as

E g = ρ a r ae,g [α q sat (T g ) -q c ] (22)
where r ae,g represents the resistance between the ground surface and the canopy air space (s m -1 ), and q sat represents the saturation specific humidity at the surface (at temperature T g ).

Fig. 6 The resistance pathways for the two-source energy budget model is shown on the left, with explicit resistances for explicit turbulent exchanges between the canopy air and the ground and the vegetation. The corresponding surface temperatures (for computing the evapotranspiration components) are shown. The single or composite energy budget scheme is shown on the right, for which both Etv and Eg are computed using the same aerodynamic resistance, rae, and surface temperature, Ts.

The relative humidity of the air is considered at the reference level (generally assumed to be at z 0h for moisture and heat) just above the surface. This quantity is usually unknown, and its modeling would require a significant amount of added complexity, therefore a proxy is used based on soil moisture in the surface soil layer. Generally, expressions of the form are commonly used:

α = min 1, a gα θ g b gα θ g + c gα (23) α = 1 2 1 -cos θ g θ f c π 2 (α = 1 if θ g ≥ θ f c ) ( 24
)
where θ g represents the volumetric soil water content in the upper soil layer (m 3 m -3 ) which as a depth (thickness) defined as z gs (m). Note that the values of the parameters a gα , b gα and c gα vary among different studies and should likely depend on soil properties (such as texture, organic content...) or they could be calibrated based on more detailed soil models or observations. The so-called field capacity volumetric water content, θ f c , is generally based on soil texture properties in many LSMs, although for local scale applications it is generally best to define it based on soil moisture observations, if available.

Note that it should be obvious that the two expressions for α generally have very different forms despite both depending on the same surface soil moisture.

The different forms for α used among LSMs are known to be a source of potentially significant differences in simulated E g (e.g. [START_REF] Mahfouf | Comparative study of various formulations of evaporation from bare soil using in situ data[END_REF].

Note that the thickness of the surface layer, z gs , has yet to be defined: the basic idea is that a diurnal cycle of the surface soil moisture can be represented if this layer is thin enough to represent the dry down during the day and the capillary recharge of the soil moisture from below during periods of low radiative forcing (generally at night). Thus, the chosen value of the surface layer thickness (which is model dependent) can have a very large impact on the E g via the selected α function since it is directly proportional to the surface water holding capacity. Finally, the application of the alpha equation to daytime conditions has been shown to work well, however, numerically speaking, certain modifications must be done at night to prevent excessive condensation mainly over very dry soils. In addition, there has been relatively

little study of the ability to accurately represent condensation in the LSM community to date.

The so-called beta method is given as

E g = ρ a r ae,g β g [h g q sat (T g ) -q c ] (25)
Compared to the alpha method, here the relative humidity of the upper soil h g is included as

h g = exp g ψ R v T g ( 26
)
where ψ represents the surface soil matric potential (m). The matric potential is related to volumetric soil water content in LSMs generally using the pedotransfer functions based upon [START_REF] Brooks | Hydraulic properties of porous media[END_REF] or [START_REF] Van Genuchten | A closed-form equation for prediction the hydraulic conductivity of unsaturated soils[END_REF]. For this method, evaporation is assumed to come from a wet layer which becomes deeper in the soil as it dries. Theoretically speaking, this method requires relatively thin soil layers in order to properly model the strong vertical gradients of soil temperature and matric potential in the uppermost (near surface) layers. Some authors set this relative humidity factor h g to unity and so β g then simply depends on the soil moisture (e.g. Deardorf, 1978). Currently, LSMs generally model the beta function as

β g = r ae,g r ae,g + r g (27)
Note that if a single composite energy budget is represented, then r ae,g is replaced by r ae in Eq. 27, and q c is replaced by q a and T g is replaced by T s in Eq. 22. r g represents a soil resistance term, that should depend on soil texture or other structural properties of the soil, but most currently used formulations are based on that of [START_REF] Sellers | Relations between surface conductance and spectral vegetation indices at intermediate (100 m 2 to 15 km 2 ) length scales[END_REF] which were calibrated for a single soil type and climate:

r g = exp b gβ -c gβ θ g θ sat ( 28
)
where b gβ and c gβ are calibrated empirical parameters. [START_REF] Merlin | An Analytical Model of Evaporation Efficiency for Unsaturated Soil Surfaceswith an Arbitrary Thickness[END_REF] suggest that the soil resistance parametrization has some conceptual problems which make its use in LSMs problematic: its use is more questionable when capillarity forces are strong or under conditions when gravitational drainage is dominating the drying of the near surface soil layers. They argue for an approach using Eq. 25 with h g = 1 and where β g is represented by an analytical equation for evaporation efficiency which depends on both near surface soil moisture and the potential evaporation. Despite very promising results for both energy and water driven regimes, use in an operational LSM would be more complex since estimating the potential evaporation at each time step within an LSM is not straightforward.

Finally, the so-called gamma method, which is based on the minimum of atmospheric demand, E d , and soil moisture supply, reads

E g = min E d , ρ a r ae,g [q sat (T g ) -q c ] (29)
In order to determine the supply limited value, an expression of the bare soil evaporation must be obtained and several formulations exist based on vertically integration of the expression for vertical water diffusive flux within the soil [START_REF] Mahrt | A 2-layer model of soil hydrology[END_REF][START_REF] Wetzel | Concerning the relationship between evapotranspiration and soil moisture[END_REF]. One advantage of this method is that it does not require a relative humidity (or proxy) term, however it has been found to be highly sensitive to the method to compute the effective depth over which the diurnal cycle of moisture is damped, the parametrization of the rate of moisture diffusion within the soil is difficult to define, and it is highly sensitive to soil properties. This method is not commonly used among LSMs currently: it is assumed this is related to the aforementioned issues, but also for more pragmatic reasons: the α and β g formulations are also more easily implemented numerically. It should be evident by now that all of the methods (Eq.s 23-24 and Eq. 29) can be summarized as depending on surface soil moisture over some near-surface layer thickness, while having very different mathematical forms (thus potentially very different predicted E g ) and theoretical underpinnings.

As a final note, some LSMs have introduced parametrizations for litter, but the approach can be very different from one to another depending on their complexity. The simplest approach is to modify or add an additional ground resistance [START_REF] Sakaguchi | Effects of soil wetness, plant litter,and under-canopy atmospheric stability on ground evaporation in the Community Land Model (CLM3.5)[END_REF], while the alternative is to model the litter using an explicit single or multi-layer model (e.g.s [START_REF] Wilson | The effect of soil surface litter residue on energy and carbon fluxes in a deciduous forest[END_REF][START_REF] Napoly | The Interactions between Soil-Biosphere-Atmosphere (ISBA) land surface model Multi-Energy Balance (MEB) option in SURFEX-Part 2: Model evaluation for local scale forest sites[END_REF]. Including a litter layer implies that E g would come from the litter directly, with a possibly soil contribution (if litter covers a fractional surface or capillary rise from the soil into the litter is assumed, but this is generally thought to be small and is neglected in simpler models). The thermal properties of litter are also generally quite different from soil, thereby having an impact on T g . And once again, the thickness of the source layer, the litter layer in this case, and the hydrological properties of the litter will have a large impact on the water storage capacity and therefore the simulated E g .

Transpiration in Land-Surface models

Transpiration refers to the loss of water from vegetation as water vapor. LSMs represent transpiration with varying degrees of complexity, usually based on the intended application. The basic processes which are represented can be summarized as follows: water (containing dissolved nutrients) is absorbed by roots and transported in liquid form via xylem into the plant cells for photosynthesis (for which the plant uses light energy and atmospheric CO 2 to produce organic material). At the leaf surface, this water is transformed from liquid to water vapour where the associated latent heat release acts to keep the vegetation relatively cool. This vapor leaves the plant through the stomata.

The subsequent increase in the matric potential gradient between the leaf and the soil then draws water into the plant from the soil. If there is insufficient soil moisture (or water stored in the plant), the stomata close and transpiration ceases.

Transpiration is generally represented using a beta formulation as

E tr = ρ a r ae,v β v [q sat (T g ) -q c ] (30)
where

β v = r ae,v r ae,v + r c (31)
in which r ae,v represents the resistance between the canopy air and the vegetation. As mentioned in the previous subsection, most LSMs use parametrizations inspired by [START_REF] Deardorff | Efficient Prediction of Ground Surface Temperature and Moisture, with Inclusion of a Layer of Vegetation[END_REF], [START_REF] Choudhury | A four?layer model for the heat budget of homogeneous land surfaces[END_REF] and [START_REF] Raupach | Simplified expressions for vegetation roughness length and zero-plane displacement as functions of canopy height and area index[END_REF]. If a single composite energy budget is represented, then r ae,v is replaced by r ae in Eq. 31, and q c is replaced by q a and T v is replaced by T s in Eq. 30 (Fig. 6). Note that in reality, the temperature of the transpiring surface is cooler than that given by the energy balance, and this can be modeled by replacing T g in Eq. 30 by an equivalent wet-bulb temperature, or a similar effect can be obtained by modification of the β v [START_REF] Milly | Potential evaporation and soil moisture in general circulation models[END_REF] or by multiplying the q sat (T g ) term by a reduction factor (akin to using an α type approach). The plant physical processes are encompassed in the r c term, which represents the stomatal resistance integrated over the entire vegetation canopy.

Stomata close under sub-optimal conditions corresponding to insufficient solar radiation, a large humidity deficit between the atmosphere and the surface, as the air temperature falls outside of a species dependent specific optimal range (heat stress), and soil water stress. Through the 1990s, most LSMs used the relatively simple stomatal resistance model of [START_REF] Noilhan | A simple parameterization of land surface processes for meteorological models[END_REF] which is based on the work of [START_REF] Jarvis | The interpretation of the variations in leaf water potential and stomatal conductance found in canopies in the field[END_REF]. It includes a representation of the aforementioned factors:

r c = R s,min LAI F 1 (SW ↓) F 2 (θ r ) F 3 (δe a ) F 4 (T a ) ( 32
)
where the inclusion of the leaf area index (LAI) indicates an integration of the resistance over the entire canopy. R s,min is the minimum or unconstrained stomatal resistance (s m -1 ), and many models use lookup tables (based on land cover type) to determine it. However, calibration is often required to obtain good results as r c and thus E tr are found to be highly sensitive to this parameter. F 1 models the effect of the photosynthetically active radiation (and is a function of the downwelling shortwave radiation, SW ↓). The remaining factors parametrize the dependence on root zone soil moisture, F 2 , the atmospheric vapour pressure deficit, F 3 , and temperature stress, F 4 .

In more recent years, r c formulations have been extended to represent photosynthesis. The stomatal conductance to water vapour can be expressed as

g c = 1.6 A n C s -C i (33)
where C s -C i (kg CO 2 m -3 ) corresponds to the gradient between the outside and the leaf intercellular CO 2 concentrations, respectively, and the factor 1.6

represents the proportionality factor of water vapour to CO 2 . The net assimilation rate, A n , (kg CO 2 m -2 s -1 ), represents the net flow of CO 2 through the plant stomata and includes environmental factors (based on available radiation, air temperature, vapor pressure deficit, etc.). There are many models which have been proposed to represent A n , and most methods used by LSMs are based on the work of [START_REF] Farquhar | A biochemical model of photosynthetic CO[END_REF], [START_REF] Collatz | Coupled photosynthesisstomatal conductance model for leaves of c4 plants[END_REF] and [START_REF] Jacobs | Stom-atal behaviour and photosynthetic rate of unstressed grapevines in semiarid conditions[END_REF]. Once A n is known, g c is then integrated from the leaf to the canopy level using as assumed canopy leaf vertical distribution and density via the LAI, and this final value is then proportional to r -1 c .

Such parametrizations permit the LSM to simulate the diurnal cycle of carbon and water vapour fluxes using CO 2 responsive representations of photosynthesis (e.g.s [START_REF] Calvet | An interactive vegetation SVAT model tested against data from six contrasting sites[END_REF][START_REF] Yang | The community Noah land surface model with multi-parameterization options (Noah-MP): 2. Evaluation over global river basins[END_REF]Bousetta et al., 2013;[START_REF] Van Den Hoof | Improved evaporative flux partitioning and carbon flux in the land surface model JULES: Impact on the simulation of land surface processes in temperate Europe[END_REF]. The relatively rapid impact of changes of net radiation and atmospheric conditions on r c and therefore ET has been explored using fully coupled surface-atmosphere models in different studies which are mentioned in Section 3.2.

The photosynthesis schemes can also be coupled with models that simulate plant phenology, which thereby can impact E tr additionally by potentially changing surface characteristics such as z 0 , LAI, and surface albedo (e.g.s [START_REF] Bonan | An ensemble square root filter for the joint assimilation of surface soil moisture and leaf area index within the Land Data Assimilation System LDAS-Monde: application over the Euro-Mediterranean region[END_REF], thus potentially improving ET estimates.

These complex coupled models add a significant number of parameters to the LSM. However, the modeling of such processes is key to making projections of feed-backs between the plant and the atmosphere for different future greenhouse gas scenarios, among other applications. For either of the r c approaches (including photosynthesis or not), the parametrization of soil moisture stress is of critical importance for computing E tr . Some LSMs using the approaches above also now include advanced moisture stress parametrizations which depend upon different types of drought responses and are distinguished for different types of vegetation (e.g. [START_REF] Calvet | Modelling forest transpiration and CO 2 fluxes -response to soil moisture stress[END_REF]). [START_REF] Verhoef | Modeling plant transpiration under limited soil water: Comparison of different plant and soil hydraulic parameterizations and preliminary implications for their use in land surface models[END_REF] compare several approaches for modeling the effect of soil moisture availability on the transpiration. They note that current methods used in most LSMs, mainly based on volumetric water content, are dated and more realistic simulations can be obtained based on using soil matric potential, and eventually chemical signaling and plant hydraulics models (including more sophisticated representations of rooting systems) in the future (since such methods require additional input parameters and observational data for a myriad of climates, soils and plant types).

The extraction of soil water by rooting systems is obviously of critical importance to transpiration, notably when defining the overall amount of water available to plants over an extended period. In terms of land-atmosphere interactions, [START_REF] Desborough | The Impact of Root Weighting on the Response of Transpiration to Moisture Stress in Land Surface Schemes[END_REF] found that the treatment of root distribution had a considerable impact on simulated E tr and discussed implications for fully coupled LSM-GCM models. LSMs have adopted essentially two approaches: using a single bulk root zone water content in the stress function or using a weighting based on the assumed vertical root zone profile which generally remains fixed in time. But several key processes are neglected using such simplified models.

For example, [START_REF] Meunier | Measuring and modeling hydraulic lift of Lolium multiflorum using stable water isotopes[END_REF] showed that roots were responsible for the transfer of moisture from deeper soil layers upwards toward the surface during the night, and this was responsible for upwards of 20% of the transpiration the following day. Also, obviously rooting systems are dynamic, and [START_REF] Zhu | Estimating the contribution of groundwater to the root zone of winter wheat, using root density distribution functions[END_REF] showed that by modeling the temporal evolution of the root density for a wheat crop substantially improved the model estimates of transpiration.

Perspectives

The study of ET brings together basic and applied scientific disciplines, indicating the need to work in an interdisciplinary framework to help progress understanding and to improve applications [START_REF] Cuxart | Current challenges in Evapotranspiration determination[END_REF][START_REF] Verhoef | Report on the first Determining Evapotranspiration workshop[END_REF]. Advances rely on the availability of comprehensive and trustable data, adequate understanding of the basics of the process for proper analysis of the data and consequently the development of new and improved parametrizations for modeling and applications.

Reference data are now provided by lysimeters and EC-systems which compare relatively well for monthly averages, but show significant discrepancies at the sub-daily scale. The study of fast response lysimeter data at a scale which resolves the diurnal cycle could pave a way forward in ABL research, since they can provide information when EC systems are not reliable, especially at night or in presence of rain and dew.

Indirect determinations of ET make use of MOST, which is valid for homogeneous terrain and for a range of thermal stabilities in the surface layer.

Research on similarity theory over heterogeneous terrain and for very stable and unstable cases is needed in order to increase the usability of these methods. In the case of using the surface energy budget, there is a debate on how to distribute the imbalance between H and λE, as some use the Bowen ratio and others propose to attribute it all to ET. Lysimeter data may help elucidating this aspect.

The classical Penman-Monteith and Priestley-Taylor approaches are still the backbones of ET parametrization for well-watered vegetated terrain. In the former, complex mechanisms in the atmosphere (entrainment) and in the soil and vegetation are all packed in one single parameter, the surface or canopy resistance. Priestley-Taylor even avoids the use of this part and substitutes it by a parameter which is determined experimentally and that is related to the entrainment rate and the canopy resistance. At the daily scale, the Makkinkde Bruin expression, which depends only on the incident solar radiation, is an alternative.

The aforementioned methods have a difficult time handling non-saturated surfaces, in which the prescription of a value for the moisture of the surface is a key research issue. As the canopy characteristics get more complicated, with varying vertical and horizontal vegetation cover, new challenges appear, which are now mostly addressed using LES coupled with detailed LSMs, in which the details of the canopy can be explicitly introduced, progressing towards firmly establishing the basis of new expressions for these complex environments.

The treatment of heterogeneity is one of the major challenges ahead for the ABL community [START_REF] Bou-Zeid | The persistent challenge of surface heterogeneity in boundary-layer meteorology[END_REF]. In the last decades, the Raupach and the advection-aridity frameworks have provided ways to conceptually address the analysis of the problem. The issue can also be explored using LES modelling, estimating the impact of prescribed surface heterogeneity, eventually combining it with experimental data from dense networks [START_REF] Mengelkamp | Evaporation over a heterogeneous land surface[END_REF]. The development of parametrizations for ET in heterogeneous terrain is an expected research activity in the coming years.

LSMs have evolved considerably since the 1990s, and development has greatly benefited in the past three decades from international collaborative efforts, as those hosted by the Global Energy and Water Exchanges (GEWEX) project, which is supported by the World Climate Research Programme (WCRP) within the World Meteorological Organization (WMO). In terms of improvements, processes which modulate the transpiration component of ET have been added and are the focus of current and future development, such as photosynthesis and plant phenology, the use of explicit single or multi-layer vegetation canopy schemes (along with more advanced representations of with-in canopy turbulence and radiative transfer) along with explicit forest litter. In terms of soil evaporation, parametrizations which are the focus of current research include multi-layer soil schemes which consider organic material, more sophisticated representations of soil water uptake by roots, more complex infiltration parametrizations, and agricultural practices, to name a few. But despite this progress, there are still well known issues, for example LSMs generally produce very different ET fluxes given the same input atmospheric forcing and basic physiographic information (e.g. [START_REF] Best | The plumbing of land surface models[END_REF].

With the continuing emphasis on representing more sophisticated processes arising from ever-increasing applications (drought monitoring, climate models which are able to model the full coupling of the Carbon cycle with the atmosphere, operational hydrological forecast systems, surface reanalysis products, seasonal weather forecasting for which the land holds a significant amount of predictability, etc.), a myriad of additional feedback mechanisms and empirical parameters are also being introduced. In order to keep making advances in terms of model physics, LSMs will require more experimental data and over more contrasting regions of the globe. Satellite-based products can assist in this regard, but it can still not replace data from high quality field sites which are critical for process understanding. But statistical methods are being increasingly used to optimally combine the aforementioned data with LSMs to provided improved estimates of ET, such as machine learning for improving the model input parameters [START_REF] Chaney | Deriving global parameter estimates for the Noah land surface model using FLUXNET and machine learning[END_REF] and data assimilation for providing improved flux estimates [START_REF] Bonan | An ensemble square root filter for the joint assimilation of surface soil moisture and leaf area index within the Land Data Assimilation System LDAS-Monde: application over the Euro-Mediterranean region[END_REF]. Thus, the future improvement in LSM physics and parameter estimation which modulate ET will likely continue to progress in tandem, and will continue to benefit from improved methods to monitor ET using the observational techniques discussed herein.

Finally, there is a growing effort in the scientific community to promote efforts to evaluate LSM schemes in fully coupled mode [START_REF] Van Den Hurk | Acceleration of Land Surface Model Development over a Decade of Glass[END_REF], despite the additional complexity compared to offline testing. It is anticipated that such efforts will lead to further improved understanding and modelling of coupled LSM ET processes. This is of significance since there is a need to improve future GCM projections of water resource evolution.

Fig. 1

 1 Fig. 1 A schematic of the ABL mechanisms related to ET. Orange arrows reflect the components of the solar radiation, dark red arrows the longwave radiation, yellow arrows the net radiation. The radiation budget at the surface is in the orange circle, which can be transient due to cloud passing. Blue arrows are evapotranspiration, dark brown sensible heat flux and light brown ground flux. Black lines with arrows represent circulations. The green circle summarizes the feedbacks as explained in van Heerwaarden et al. (2009).

  of the CBL and, consequently, on the evolution of the ABL(van Heerwaarden el al., 2014). The vegetation density, typically characterized by the leaf area index, relates to the amount of transpiration and modifies the incident solar radiation as it penetrates the canopy[START_REF] Pedruzo-Bagazgoitia | Direct and diffuse radiation in the shallow cumulus-vegetation system: Enhanced and decreased evapotranspiration regimes[END_REF].The presence of stratified clouds affects ET depending on their optical thickness. Dense clouds essentially prevent solar radiation from reaching the surface, in contrast to thin clouds, which may have larger ET than in a cloudless day because of an increase of the amount of diffuse light and a decrease of direct light, resulting in a reduced canopy resistance[START_REF] Pedruzo-Bagazgoitia | Direct and diffuse radiation in the shallow cumulus-vegetation system: Enhanced and decreased evapotranspiration regimes[END_REF]. Passages of cumulus clouds decrease both the surface latent and sensible heat fluxes, but λE/R n is larger under the cloud because the atmospheric water demand is not changing while the ground flux G is directed upwards causing more evaporation from the upper soil[START_REF] Lohou | Surface energy balance and buoyancy response to shallow cumulus shading[END_REF]. The effect of a cloud passage impacts the leaf response since the opening/closing of the stomata has a mechanical response time on the order of several minutes(Sikma and Vilà-Guerau de Arellano, 2019) leading to a decoupling between E tr and the incoming solar radiation. The treatment of such transient conditions is a challenge for current methods used to determine surface turbulent fluxes since they assume stationarity. Conditional sampling strategies have been developed over forests with dense canopies by[START_REF] Kivalov | Observing the whole-canopy shortterm dynamic response to natural step changes in incident light: characteristics of tropical and temperate forests[END_REF]: they found that ET takes between 5 to 10 min to respond to a change of solar irradiation, in contrast to the air and leaf temperature that change more rapidly.

  period of time of the order of an hour or less, taking into account the diurnal cycle, refering to daily values when appropriate. The explicit treatment of sublimation from ice or of condensation of water over surface elements are not part of this work. The evaporation from rain or irrigation water intercepted by the vegetation will be commented upon in Subsection 5.2, within the framework of LSMs.

Fig. 2

 2 Fig. 2 Schematic of the experimental methods to determine ET.

  McColl (2020) has revised some elements of the derivation, in particular substituting the linearized form of the Clausius-Clapeyron equation by the expression proposed by Vallis et al. (2019), leading to improvements when comparing to EC values from the FLUXNET database (fluxnet.ornl.gov).

  development was the use of two big leaves in one layer, one sunlit the other one shaded, introduced by[START_REF] Sinclair | Development and evaluation of simplified models for simulating canopy photosynthesis and transpiration[END_REF] for which[START_REF] Wang | A two-leaf model for canopy conductance, photosynthesis and partitioning of available energy I:: Model description and comparison with a multi-layered model[END_REF] incorporated an A-gs model and a radiative transfer model.The transpiration of a vegetated surface depends on its canopy resistance, which is determined by the physiology, structure and root system of the plants composing the canopy, each of which has its own particular stomatal resistance.The stomata open when the plant photosynthesizes with incoming solar radiation and close when the loss of water to the atmosphere compromises the integrity of the plant cells. As indicated above, stomata take several minutes to react to light changes[START_REF] Woods | Stomatal response to changing light by four tree species of varying shade tolerance[END_REF]) and cumulus cloud passages alternating direct and diffuse light may decouple ET from actual changes in photosynthetically-active radiation illumination, as explored by Pedruzo-Bagazgoitia et al. (2017) using a scheme with two big leaves.

  1983, de Bruin showed that when dry entrainment at the top of the CBL equilibrates with the surface evaporation, α varies with the canopy resistance r c and has a diurnal cycle. For well-watered surfaces (r c =0) the estimated values were between 1.2 and 1.4, in good agreement with Priestley and Taylor, whereas α diminished with increasing r c , taking values close to 1 for wellwatered grass (r c =60 m s -1 ) and down to 0.6 for surfaces with r c =250 m s -1 .In accordance with those findings,[START_REF] Cristea | Revised coefficients for Priestley-Taylor and Makkink-Hansen equations for estimating daily reference evapotranspiration[END_REF] summarized the values for α obtained experimentally for a number of sites, showing the lowest values for dry ventilated sites and the highest for wet locations. The values of α for wet surfaces are in good agreement with those obtained with the equation for this parameter derived by[START_REF] Van Heerwaarden | Interactions between dry-air entrainment, surface evaporation and convective boundary-layer development[END_REF].

Fig. 3

 3 Fig. 3 The two conceptual terms in Penman-Monteith equation. Left: the radiation term (E rad) as a function of the available energy (Rn -G); right: the aerodynamic term (E aero) as a function of the atmospheric conditions, the aerodynamic and the canopy resistances.

( 2015 )

 2015 used the Priestley-Taylor equation for E p and Penman's formula over open water for E pa , thereby providing a formula for the actual E. The temporal resolution of the formula is then just a matter of the corresponding resolution of the input data.

Fig. 5

 5 Fig. 5 Schematics of the existing approaches to estimate ET. The methods with an asterisk represent those based on daily values.

  The imbalance is usually closed by distributing it between λE and H following some proportionality, typically the Bowen ratio (B=H/λE) or the evaporative fraction. The validity of this assumption is unclear, as[START_REF] Wang | An extremum principle of evaporation[END_REF] showed that ET is maximized in balancing the surface energy budget.[START_REF] Gebler | Actual evapotranspiration and precipitation measured by lysimeters: a comparison with eddy covariance and tipping bucket[END_REF] chose to attribute all the energy imbalance to λE (in what they call the worst case scenario), finding increases in the monthly averages of ET up to 25% in the summer and up to 37% in the winter, resulting in values very similar to those provided by

	a lysimeter.

On the other hand, if H is available from a sonic anemometer, but we do not have a high-frequency determination of λE (a very common situation), then λE can be estimated directly as the residual if we proceed as

[START_REF] Gebler | Actual evapotranspiration and precipitation measured by lysimeters: a comparison with eddy covariance and tipping bucket[END_REF]

or corrections can be made using the Bowen ratio computed with the vertical gradients of temperature and humidity.

2.4 Section summary

Figure

2

summarizes the methods described in this section. The direct determinations of ET by EC systems or lysimeters each have their own limitations, but it is commonly accepted that, for monthly averages, EC systems tend to underestimate ET while lysimeters (if they are properly set) provide more accurate values at the scale of their plot. Flux chambers give direct estimations at ground and leaf level, but at the price of modifying the measurement conditions, which implies a need to define accurate sampling strategies. The plant physiology methods provide estimations of transpiration and need to be calibrated against EC or lysimeters while monitoring the other components of the water balance. The SEB residual method usually requires calibration and can be considered a direct determination of ET if λE is made equal to the energy imbalance with H provided by an EC-system, with corrections if necessary.

Acknowledgements The authors want to acknowledge the input of an anonymous reviewer that greatly contributed to the final form of the paper. Comments on the revised version were provided by Jordi Vilà-Guerau de Arellano, Anne Verhoef, Jannis Groh and Alexander Graf. Discussions with colleagues under the auspices of the GEWEX Evapotranspiration Working Group initiated at the 8th GEWEX Open Science Conference in Canmore, Canada in 2018 and then at the GEWEX Hydroclimatology Panel Workshop 'Determining Evapotranspiration', Sydney, Australia in 2019, have been a source of inspiration for this work. Funding has been provided by the Spanish research project of FEDER/Ministerio de Ciencia e Innovación and Agencia Estatal de Investigación, RTI2018-098693-B-C31.