Continuous Daily Evapotranspiration with Optical Spaceborne Observations at Sub-Kilometre Spatial Resolution

: Evapotranspiration (ET) is a key parameter in the description of the energy and water ﬂuxes over land. Continuous and spatially detailed ET simulations are thus required for a number of scientiﬁc and management-related purposes. These conditions are determined by the modelling approach and the composition of the forcing dataset. This study aimed at simulating daily ET in a diversity of climate and land cover conditions at a spatial resolution of ∼ 1 km and higher. The modelling approach was based on the algorithm driving the ET product developed and set in operations in the framework of the Satellite Application Facility on Land Surface Analysis programme (LSA-SAF). The implemented algorithm allowed the ingestion of biophysical parameters derived from SPOT-V and PROBA-V observations developed by the Copernicus Global Land Programme, as well as other model parameters at a similar spatial resolution. The model was tested at an ∼ 1 km spatial resolution in over 40 sites located in different climate and land cover contexts. The implementation at ∼ 300 m was tested in the upper Biebrza basin, in Poland. The simulations correlated well with the validation dataset (r2 > 0.75 in 80% of sites) and exhibited root mean squared values lower than 1 mm/day in 80% of the cases. The results also pointed to the need for reﬁning the accuracy of soil moisture data sources, especially in dry areas. The results showed the ability of the modelling approach and the SPOT-V/PROBA-V missions to support the generation of long ET time series. They also opened the gate to incorporate Sentinel-3 in ET continuous modelling.


Introduction
The timing and magnitude of energy and water fluxes at the Earth's surface result from the interplay of a number of processes and conditions in the soil-plant-atmosphere continuum. The accurate description of those processes is a major concern towards understanding the mechanisms and drivers of terrestrial ecosystems.
The release of water vapour from the Earth's surface into the atmosphere is called evapotranspiration (ET). ET comprises the volume of water evaporated from water bodies and reservoirs, the upper layers of the soil and intercepted water in different layers of the vegetation canopy, as well as the water transpired by photosynthetic organisms. Needless to say, the accurate and timely estimation of this component of the water balance is relevant in a large number of domains like agriculture monitoring and yield prediction [1], drought assessment [2][3][4], water catchments' management [5,6], urban planning [7,8], etc.
ET plays a major role in the global climatic system as it is a key component of both the carbon and water cycle [9]. The heating effect of solar energy on Earth is tightly connected to ET as a substantial fraction of the solar radiation triggers the evaporative process (latent heat flux). The resulting ET rates also depend on local conditions (like water vapour demand, soil moisture, land cover, etc.).
The increasing interest in understanding and estimating ET at different spatial and temporal resolutions and in a continuous way has driven the development of algorithms and the search for optimal ways of exploiting available datasets to force those models. In this respect, a great deal of the developments and data exploration efforts are related to the use of spaceborne data in ET models [10] as the only means to visit every point on Earth on a regular basis virtually. These developments have been boosted in the last few years by the launch of satellite missions offering increased sampling frequency, spatial detail and timely availability in the generated datasets. More spatial detail in the observations implies better accounting for landscape heterogeneity, and a higher frequency of observations aids the construction of continuous time series of relevant land surface parameters.
This study aimed to generate continuous series of daily ET at an ∼1 km spatial resolution or higher in different ecosystems across the globe by combining remote sensing (RS) data, meteorological data and land cover-related parameters in a surface energy balance setting. Much attention was paid to the notion of continuity whereby: (a) daily ET series should cover several years and; (b) the number of missing values should be reduced to a minimum. These features rely on the availability of long and consistent time series of the model forcing variables and on algorithm-related characteristics coping with non-ideal conditions for RS observations; e.g., clouds.
As for continuity in RS data, the compatibility between the SPOT-Végétation (S-Vgt) and the PROBA-V mission offers interesting opportunities. The S-Vgt mission provided a wealth of spaceborne observations in the course of its operational life-between 1998 and 2014-at an ∼1 km spatial resolution. An important criterion in the design of the PROBA-V mission-operational since 2013-was the continuity of the S-Vgt observations and products. This aspect has been already explored and exploited in several studies [11][12][13]. In addition, the PROBA-V mission conducts observations with enhanced spatial detail (∼300 m and ∼100 m) while having equivalent radiometric characteristics across the three spatial resolutions.
The ET model in this study was based on the algorithm driving the operational ET product developed in the framework of the Land Surface Analysis Satellite Application Facility (LSA-SAF) [14][15][16]. The LSA-SAF ET algorithm has been in operation for about 10 years now and has proven its robustness in a number of applications [6,[17][18][19][20] and the ability to operate in all-weather conditions [21].
This article presents an implementation of the principles at the basis of the LSA-SAF ET algorithm at an ∼1 km spatial resolution in various locations representing different ecosystems around the globe where in situ measured validation data were available. The study also includes ET simulations over a particular region in northern Poland where the spatial grid of PROBA-V 300 m products was the target spatial resolution. Finally, the study concludes with some considerations towards the opportunities given by Sentinel-3 to be ingested in the ET algorithm after the operational period of PROBA-V.

The ET Algorithm
The fundamentals of the algorithm can be expressed with the well known Equations (1) and (2) of the surface energy balance.
where: R net : net radiation; LE: latent heat flux; H: sensible heat flux; G ground heat flux; α: albedo; S ↓ and L ↓: downwelling short and long wave solar radiation; ε: emissivity; σ: the Stefan-Boltzmann constant; T sk : skin temperature. Various ET algorithms depart from these basic expressions, but differ in the approach followed to estimate the components, as well as in the forcing and the discretization of the temporal and spatial domains. The methods applied in this study followed greatly the algorithm powering the ET operational product designed in the framework of EUMETSAT's programme Land Surface Analysis Satellite Application Facility (LSA-SAF) [14]. Thorough descriptions of this algorithm can be found in Ghilain et al. [15] and Ghilain et al. [16] and in the documentation of the LSA-SAF ET product (http://lsa-saf.eumetsat.int). The general lines of the procedure of ET estimation are described in the pseudocode presented in Figure 1.
Require: S ↓, L ↓, α, ε, SM , veg, weather, lcover procedure ET calculation for c in raster cells do for t in tiles do set up initialization Figure 1. General procedure for ET estimation where: LE t : latent heat flux; H t : sensible heat flux; G t ground heat flux; α: albedo; S ↓ and L ↓: downwelling short and long wave solar radiation; ε: emissivity; T sk : skin temperature; SM: soil moisture; weather: meteorological forcing; lcover: parameters related to land cover.
The procedure indicated in Figure 1 was implemented at a time step of 30 min. The half-hourly estimates were then aggregated to obtain the target product: daily ET. The basic spatial unit was a cell c belonging to a gridded representation of the area of interest.
The procedure was fed by a series of variables and parameters containing information on incoming solar radiation and energy reflected by the surface (S ↓, L ↓, ε, α), vegetation type and dynamics (grouped as veg), soil moisture (SM), meteorological conditions (various variables grouped under weather) and the land cover classes (lcover) and their occupation fraction (ζ) in each cell.
One cell may contain up to 4 tiles t, which corresponded to homogeneous land cover types. The computation of the energy balance was performed at the tile level. It occurred in an iterative way whereby the difference between R net and the aggregated energy fluxes was not allowed to exceed a threshold of 0.1; i.e., the energy balance closure was ensured. The occupation fraction of each tile ζ was then used to weigh the contribution of each tile to the ET cell value. The estimation of LE and H was preceded by the computation of the response of vegetation to environmental determinants (canopy resistance) and the aerodynamic resistance. The former was driven by vegetation-related parameters (veg), soil moisture (SM) and the atmospheric water pressure deficit, which is related to the meteorological conditions (weather). The latter was determined by the atmosphere stability, derived from weather and lcover. The soil heat flux (G) was computed as a function of R net and veg.
This algorithm was designed such that the temporal resolution of the forcing would account for the variability in major drivers of ET throughout the day; notably, solar radiation and meteorological conditions. It is noteworthy to point out that this algorithm computes the T sk term of Equation (2) in an iterative way. An alternative approach would be to retrieve land surface temperature (LST) from spaceborne data sources. Such an approach, however, would entail disadvantages related to the dependence of LST observations on cloudiness; it would then produce many gaps in the sub-daily ET estimates. The iterative estimation of T sk in an energy balance setting has proven effective to generate estimates equivalent to RS LST [21].

Study Sites
The set of study sites encompassed a large range of conditions in terms of climate and vegetation cover. The sites were selected such that in situ measured validation data were available; specifically, energy fluxes measured at eddy covariance (EC) towers. Most of the simulations were done on sites belonging to the FLUXNET network [22]. The FLUXNET data used in this study followed the guidelines of the CC-BY-4.0 data usage license (Attribution 4.0 International (CC BY 4.0); https: //creativecommons.org/licenses/by/4.0/). The location of the study sites is shown in the map in Figure 2, which also shows the distribution of the Köppen-Geiger climatic classes as obtained by Beck et al. [23]. More details on these sites are shown in Table 1. The simulations on the sites listed in Table 1 were performed at a spatial resolution of 1 km (details on the forcing variables are given in further sections of this article). In addition to that, the algorithm was implemented at a 300 m resolution across the upper part of the Biebrza basin in Poland. The Biebrza basin is one of the largest and most valued protected areas in Europe. By virtue of its ecological value, it is listed in the Natura2000 network (https://ec.europa.eu/environment/nature/natura2000/) and belongs to RAMSAR convention for wetlands (https://www.ramsar.org).
The study area in the Biebrza basin comprised over 35,000 hectares, indicated in the map of Figure 3. The map also shows the land cover classes in the study area according to the CORINE land cover map [24]. The simulations in this area were contrasted against in situ measurement conducted at a flux tower located in Rogozynek (53 • 42'8.19", 23 • 24'46.9").

Meteorological Variables
The meteorological forcing for the simulations on the sites listed in Table 1 was obtained from in situ measurements contained in the FLUXNET2015 dataset (https://fluxnet.fluxdata.org/data/ fluxnet2015-dataset/) and ERA5, a dataset of reanalysed meteorological and atmospheric observation generated at the European Centre for Medium-Range Weather Forecasts (ECMWF) [25] (ERA stands for ECMWF Re-Analysis). Half-hourly data on incoming solar radiation (S ↓, L ↓) and wind speed were taken from FLUXNET2015. Other meteorological variables were taken from the hourly records of ERA5 as the in situ air temperature and relative humidity data series had gaps. The ERA5 meteorological data were then interpolated to match the half-hourly solar radiation time series. Data on soil moisture and soil temperature also were taken from ERA5 for 4 soil layers.
ERA5 was the main source of meteorological forcing for the simulation over the upper Biebrza basin as well. However, the radiation terms for the simulations over the Biebrza basin were taken from the LSA-SAF downwelling short and long solar radiation products [26,27]. These products are based on observations by the SEVIRI sensor onboard the geostationary Meteosat Second Generation Satellite (MSG) in the electromagnetic regions 0.3-4.0 µm and 4-100 µm, respectively. These products are generated with a 0.5 h temporal resolution.

Polar-Orbit Satellites
The amount of vegetation foliage in the canopy exerts a great impact on the composition of the ET flux. It determines, amongst other, the transpiration flux resulting from photosynthesis, the amount of rainfall water intercepted by the canopy, the soil water uptake, the resistance to wind currents, the fraction of R net stored in the upper soil layers, etc.
In order to account for vegetation dynamics and albedo, the ET model ingested estimates derived from observations by polar-orbit satellites. These products are able to capture the seasonal patterns of vegetation and are normally available at moderate/high spatial resolution. The model variables taken from polar-orbit satellites were leaf area index (LAI), albedo (α) and fraction of vegetation cover. The LAI values were extracted from the dedicated product generated in the framework of the Copernicus Global Land Programme (CGLP) (https://land.copernicus.eu/) [28]. This product consists of a continuous archive of LAI estimates built on the basis of coupling the operation period of SPOT-V and PROBA-V [29]; which makes more than 20 years of global vegetation-related information. Albedo [30] and fraction of vegetation cover are part of the CGLP portfolio as well. The former was required in the computation of R net , and the latter served as a basis for weighing the contribution of each tile within the cells.
The products used for the simulation over the sites in Table 1 had a spatial resolution of ∼1 km. The simulations over the upper part of the Biebrza basin were conducted over the grid of the PROBA-V imagery at an ∼300 m spatial resolution [31]. All needed RS products but albedo were available in the CGLP portfolio at an ∼300 m spatial resolution for the period starting in late 2013 (the start of operations of PROBA-V). The albedo values used at the ∼300 m resolution were obtained from the resampling of the albedo MODIS MCD43A2 product [32], which combines observations by the MODIS sensor onboard the TERRA and AQUA satellites.

Land Cover
A number of important model parameters (minimum stomatal resistance, root density in different soil layers and vegetation height, amongst others) was associated with the land cover. In this study, the Ecoclimap I database [33] was used to extract information on the type of ecosystem at each studied location in Table 1. The Ecoclimap database covers the globe at an ∼1 km spatial resolution. Every point on Earth is associated with an ecosystem type, which in turn is composed of tiles representing functional vegetation types, water bodies, bare soil, urban environments, impermeable surfaces or snow. The Ecoclimap I database also contains data on the occupation fraction of the different tiles in each ecosystem. This spatial information was taken from Ecoclimap I, and the actual parameter values associated with the tiles were the same as those used in the operational ET LSA-SAF product [34].
The information on land cover over the studied site in the upper Biebrza basin was obtained from a more recent version of the Ecoclimap database, the Ecoclimap Second-Generation (ECOCSG) [35]. Novel elements of ECOCSG with respect to older versions of Ecoclimap are: (a) higher spatial detail (∼300 m); and (b) each cell in the ECOCSG map corresponds to only one tile. The dataset was built on the basis of the ESA-CCI map and contains 33 land cover classes, out of which 9 are urban.
The criterion of the Ecoclimap database was followed for deriving the value of emissivity ε as a function of the fraction of vegetation cover ( f cov): ε = f cov * 0.97 + (1 − f cov) * 0.94. Table 2 presents an overview of the data sources described in the previous paragraphs and indicates the simulation domain (flux tower locations or upper Biebrza basin) where each data source was used. Table 2. Overview of data sources to force the ET model at flux tower locations and across the study area in the upper Biebrza basin, Poland. The SPOT-Végétation (S-Vgt) and PROBA-V products were generated in the framework of the Copernicus Global Land Service at a 1 km and a 300 spatial resolution, respectively. LSA-SAF, Land Surface Analysis programme (LSA-SAF); ECOCSG, Ecoclimap Second-Generation.

Validation
The estimated daily ET was validated with in situ measurements at the sites listed in Table 1 and compiled in the FLUXNET2015 dataset. The LE data actually used for the validation were selected based on the quality flag provided in FLUXNET2015 and were taken from series where no correction of energy balance closure was performed.
The validation data for the simulation over the upper Biebrza basin were obtained from measurements at a flux tower located within the study area, in the village of Rogozynek. This dataset proved useful in previous research studies concerning evaluation of surface energy balance and soil-vegetation-atmosphere transfer models [36] and finding a correlation between spectral indices and greenhouse gas fluxes [37].
The simulated data series were compared with the measurement at the flux towers on the basis of commonly used statistical parameters: root mean squared error, bias, coefficient of determination and the Nash-Sutcliffe efficiency coefficient (NSE) [38]. Table 3 presents statistical parameters that indicate the ability of the ET model to simulate the daily ET derived from the measurement at flux towers. Table 3. Bias (mm/day), root mean squared error (RMSE) (mm/day), coefficient of determination (r 2 ) and the Nash-Sutcliffe efficiency index obtained from contrasting the simulated daily ET against in situ measurements at flux tower stations as compiled in the FLUXNET 2015 datasets.  In most of the studied sites, the simulated daily ET exhibited high correlation with the validation data. This translated into r 2 values greater than 0.75 in most of the sites in all plant functional types. The r 2 values were even higher than 0.9 in various sites. In the best cases, high correlation between simulation and validation data series was accompanied by high NSE values, which indicated similarity between modelled and measured values; i.e., low bias and low RMSE. Such a condition was present in sites associated with all considered plant functional types.

Simulations at Flux Towers' Locations
The plots of Figure 4 depict the results of the simulations in terms of the obtained r 2 and RMSE per site in association with the corresponding plant functional type and climatic zone. These plots show that all plant functional types were represented in the r 2 > 0.5 and RMSE < 1 graph quadrant. Most of the simulations falling outside that area exhibited RMSE values higher than one and belonged to climatic zones associated with tropical regions or dry climate. The plots in Figure 5 depict the variability in simulated and validation daily ET throughout the year for sites in the considered plant functional types. Although this analysis could vary from site to site, the examples shown in Figure 5 indicate the nature of the similarities and dissimilarities between simulations and validation data expressed in the statistical parameters of Table 3. In general, the plots show the ability of the model to capture the most important vegetation phenology events; namely, the rise and decline of the growing season. Likewise, the plots show the periods of overor under-estimation of the actual ET. For instance, the simulations for Brasschaat Belgium (BE-Bra) tended to be higher than the measured values; the model tended to slightly overestimate ET in Tchizalamou Congo (CG-Tch) during the period of lowest ET rates; a short period where the highest ET rates in Puechabon France (FR-Pue) were underestimated, etc. Nevertheless, the plots show an important degree of overlap between the areas representing the variability of the simulated and validation datasets. The output of daily ET simulations over the sites in Table 1 are made available to the reader as supplementary material of this article. Figure 5. Variability of daily ET as measured at the flux towers (light blue area) and simulated (shaded area). Variability was computed as the average ± the standard deviation.

Simulations in the Upper Biebrza Basin
The simulations conducted for the years 2015-2016 at the location of the flux tower in Rogozynek, in the upper part of the Biebrza basin, exhibited a positive bias of 0.1685 mm/day with a root mean squared error equal to 0.6930 mm/day (n = 311). The coefficient of determination (r 2 ) was equal to 0.8705, and the value of the Nash-Sutcliffe Efficiency index was 0.6598. The time series of simulated and validation daily ET are presented in Figure 6. The statistical parameters presented in the previous paragraph and the features of the time series shown in Figure 6 indicated overall good agreement between the output of the simulations and the in situ measured data. The simulated time series correlated well with the validation data, and the time series plot allowed the identification of the particular periods of the growing season where both time series differed the most.
The modelled daily ET at the location of the flux tower in Rogozynek was extracted from simulations conducted across a larger area at a 300 m spatial resolution; specifically, over the grid of PROBA-V 300 m imagery. Examples of those simulations for the growing season 2015 are presented in Figure 7. Although no means were available to validate simulations outside the footprint of the flux tower, the daily ET maps of Figure 7 depict spatial patterns that could be expected in the considered dates and for the landscape arrangement of the study area (see Figure 3). The highest ET rates were located in the river channel and adjacent wetlands. The maps also allow distinguishing the ET patterns in the forest area in the northwest side of the study area, as well as a more scattered pattern in other areas where mosaics of grassland and crops were located.
As described in Section 2.1, the ET algorithm in this study did not ingest LST as part of the forcing dataset. It rather estimated the T sk term of the energy balance through an iterative process. One should expect, however, that the estimate of T sk would be similar in magnitude and shape to LST obtained from spaceborne observations. This could serve as an additional test for the coherence of the simulations. Figure 8 presents scatterplots that contrast, for four land cover classes, the T sk estimated internally by the model against LST data retrieved from the MODIS MOD11A1 LST daily product at a 1 km resolution [39]. The comparison of T sk with MODIS LST does not imply that the latter should be considered as reference for validation. Nevertheless, finding high correlation between these two independent and different in nature datasets is a positive indicator of coherence in the results.

Discussion
The main motivation for implementing the LSA-SAF ET model at a higher spatial resolution setting (∼300 m and ∼1 km) was rooted in the need for continuous estimates of land surfaces variables at useful spatial resolutions for monitoring management units (forest stands, agriculture regions, wetlands, protected areas, etc.). The results presented in the previous section revealed the potential of the approach followed here and confirmed the value of a number of principles in the LSA-SAF ET model that proved useful and valid in the generation of daily ET estimates at a sub-kilometre spatial resolution. At least two aspects can be underlined in that respect: (a) The fruitful combination of geostationary and polar-orbit satellites towards the generation of ET estimates: The approach followed in this study strived to utilise spaceborne information optimally such that ET drivers exhibiting high variability during the day (namely, short and long wave incoming radiation) would be obtained from RS data sources based on geostationary satellites. Other ET determinants that were closely related to land cover and exhibited seasonal temporal patterns (like LAI, albedo, fraction of vegetation cover) may be better derived from high spatial resolution data sources, i.e., polar-orbit satellites like SPOT-V, PROBA-V and, more recently, Sentinel-2 and Sentinel-3.
(b) The continuity element in daily ET estimates was emphasized in the Introduction Section of this article. The non-dependency of the LSA-SAF ET algorithm on RS-based LST was a crucial factor in reducing the gaps in the generated time series due to the high sensitivity of spaceborne LST to cloudiness. The comparison with MODIS LST (Figure 8) indicated that this approach yielded good results.
Moving towards more spatial detail entails the challenge of building forcing datasets compatible with the target spatial resolution. The advent of spatial missions with enhanced performance has contributed to make this possible. Nevertheless, currently available data of soil moisture and/or soil characteristics might still be too coarse, with respect to the target spatial resolution, and may result in poor simulations of ET. This might partially explain the somewhat low statistical scores in some sites in climatic zones where soil moisture was the main driver of ET. The discrepancy might only be present in particular periods of the year, as shown in Figure 9, where the simulated daily ET rates were underestimated during the dry season in that part of the world. This was, in turn, related to underestimating the actual available volume of soil water for transpiration. It also illustrates the sensitivity of the model to proper information on soil moisture and soil physical properties determining the actual availability of soil water. This study considered the ERA5 dataset as the source of soil moisture. Future implementations should test the performance of other soil moisture data sources. In this respect, the approach adopted in the most recent version of the LSA-SAF ET product is worth considering. This development was based on the principle of thermal inertia to derive estimates of soil moisture and relied on LST data derived from the geostationary MSG satellite [40]. This approach proved successful in arid and semi-arid areas.
The daily ET simulations over the upper Biebrza basin were a first experience in the use of ECOCSG as forcing of the ET model presented here. The simulations of this study adopted the parameterization of the LSA-SAF ET model and applied it at a higher spatial resolution. Future simulations based on ECOCSG, and implemented in other regions, might benefit from revising and fine tuning the parameters as this dataset included new vegetation types that were not present in previous versions of Ecoclimap.

A Glance into the Future
The results of this study, as well as the findings in previous research [13], showed the suitability of products derived from SPOT-V and PROBA-V observations to integrate the forcing dataset of the ET model implemented here. The PROBA-V operational phase allowed taking a step towards more spatial detail in daily ET simulations. However, in view of the upcoming end of the operational life of PROBA-V, it is pertinent to state future implementations of this algorithm in terms of products derived from Sentinel-3 observations.
The ET algorithm used in this study ingested RS products derived from observations of optical sensors. Therefore, the Sentinel-3 sensor of interest is the Ocean and Land Color Instrument (OLCI), which senses the Earth in 21 spectral bands in the region 400-1020 nm. This region overlaps with three of the four bands of PROBA-V. The full resolution products are delivered at a 300 m spatial resolution.
Including OLCI observations in ET simulations, as those performed in this study, would require the translation of its spectral information into broad-band albedo, LAI and fraction of vegetation cover. These parameters are the basic ingredients of the LSA-SAF ET algorithm, as well as of other similar tools simulating different processes of the land surface. Although currently delivered OLCI data contain the necessary spectral information for deriving such parameters, the interest of the land surface modelling community in OLCI could be stimulated, making it part of the list of standard delivered OLCI products. Previous efforts in the derivation of land biophysical values with MERIS could be the starting point [41,42].
Currently, the OLCI Global and Terrestrial Chlorophyll Indices (OGVI and OTCI, respectively) are available as L2 products. These products are intended to give continuity to the analogous MGVI and MTCI based on ENVISAT/MERIS [43,44]. The ET LSA-SAF algorithm is not designed to ingest any of these products; moreover, none of them is a full equivalent of LAI. However, MGVI-the precursor of OGVI-represented FAPAR [45], which is known to be related to LAI [46]. The maps and plots in Figure 10 show the existence of common spatial patterns and a substantial degree of correlation between OGVI and PROBA-V LAI in our study area in Biebrza, Poland. On the basis of the relationship PROBA-V LAI and OLCI OGVI, one could derive empirical models to estimate LAI and implement the LSA-SAF ET for local or regional studies. Estimates of the fraction of vegetation cover and albedo could be approached in a similar way for fast and local appraisals. As suggested above, the full exploitation of OLCI observations in land surface models would be eased by developing such biophysical parameters in an operational framework.

Conclusions
This study aimed at conducting daily ET simulations at an ∼1 km spatial resolution in a variety of climate and land cover conditions. The approach relied on the principles behind the LSA-SAF ET algorithm and the adequacy of biophysical variables derived from SPOT-V and PROBA-V to integrate the forcing dataset.
The approach delivered satisfactory results when contrasted against measurements done at eddy covariance stations. The obtained daily ET estimates correlated well with the validation time series, and the comparison of the internally computed T sk in the upper Biebrza basin with spaceborne LST exhibited great similarity. These were all indicators that featured the approach presented here as a promising tool to support applications where continuous daily ET at sub-kilometre spatial resolution may be required (e.g., watershed management, water accounting, regional drought monitoring, etc.).
The results also pointed to avenues of future research envisaging further improvement of the model and the exploitation of OLCI Sentinel-3 observations. The adoption of ECOCSG in the simulations in the upper Biebrza basin delivered good results; future tests in other climate zones and landscape configurations should be performed to verify its suitability in other conditions. In parallel, the best methods to translate current OLCI products into biophysical parameters to feed the ET algorithm (and other land surface models) should be investigated. Moreover, the step towards higher resolution might entail the upgrade of some land cover-related parameters, which may help with reducing the error of the simulations.