Irrigation Timing Retrieval at the Plot Scale Using Surface Soil Moisture Derived from Sentinel Time Series in Europe

The difficulty of calculating the daily water budget of irrigated fields is often due to the uncertainty surrounding irrigation amounts and timing. The automated detection of irrigation events has the potential to greatly simplify this process, and the combination of high-resolution SAR (Sentinel-1) and optical satellite observations (Sentinel-2) makes the detection of irrigation events feasible through the use of a surface soil moisture (SSM) product. The motivation behind this study is to utilize a large irrigation dataset (collected during the ESA Irrigation + project over five sites in three countries over three years) to analyze the performance of an established algorithm and to test potential improvements. The study’s main findings are (1) the scores decrease with SSM observation frequency; (2) scores decrease as irrigation frequency increases, which was supported by better scores in France (more sprinkler irrigation) than in Germany (more localized irrigation); (3) replacing the original SSM model with the force-restore model resulted in an improvement of about 6% in the F-score and narrowed the error on cumulative seasonal irrigation; (4) the Sentinel-1 configuration (incidence angle, trajectory) did not show a significant impact on the retrieval of irrigation, which supposes that the SSM is not affected by these changes. Other aspects did not allow a definitive conclusion on the irrigation retrieval algorithm: (1) the lower scores obtained with small NDVI compared to large NDVI were counter-intuitive but may have been due to the larger number of irrigation events during high vegetation periods; (2) merging different runs and interpolating all SSM data for one run produced comparable F-scores, but the estimated cumulative sum of irrigation was around −20% lower compared to the reference dataset in the best cases.


Introduction
Irrigation water use is a critical factor in managing and optimizing agricultural water resources amid the increasing global demand for irrigation water [1]. The rise in demand is particularly noticeable in response to the production of bioenergy [2], while the decrease in groundwater availability amplifies competition with other sectors, including urban water demand [3]. An irrigation, like a precipitation event, has a specific start and end time, during which a given amount of water is applied to the field, usually referred to as the amount or dose of irrigation. This quantity is expressed either as a volume (m 3 ) or a water depth (mm) at the field (plot) level. Quantifying irrigation has been mainly carried out at large scales using a method based on coarse-scale surface soil moisture products derived from either passive or active remote sensing sensors [4][5][6][7]. These techniques operate in the microwave domain where the dielectric properties of soil are primarily governed by its water content, which directly affects the emissivity and backscattering of the soil surface [8]. However, such methods operate at much lower resolutions than the field scale, and detecting irrigation events at the field scale is an emerging scientific area.
Land surface temperature provides an attractive proxy of root-zone soil moisture conditions as hydrologically stressed vegetation heats up in response to stomatal closure [9]. Olivera-Guerra et al. [10] derived crop water stress from Landsat land surface temperature data to feed a simple crop balance model for estimating irrigation water timing and amount. However, this method faces limitations, mainly attributed to the low revisit time of available thermal sensors into orbit. Nevertheless, this method is promising in the short term with the launch of new thermal sensors such as TRISHNA (Thermal InfraRed Imaging Satellite for High-resolution Natural resource Assessment, https://www.eoportal.org/satellitemissions/trishna, accessed the 2 February 2023) expected for 2024 or LSTM (Land Surface Temperature Monitoring, https://www.eoportal.org/satellite-missions/lstm accessed the 2 February 2023), expected for 2028.
In addition to the work mentioned above, most of the approaches developed recently rely on surface soil moisture (SSM) estimates derived from Sentinel-1. The dielectric properties of the soil phases (i.e., water, air, and solids) affect the emissivity and backscattering of microwaves from the soil surface [8], making it possible to retrieve soil moisture from microwave observations. The Sentinel-1 C-band satellites have penetration ranges from 1.87 to 3.75 cm and 0.75 to 1.5 cm for bare soil and cultivated soil, respectively, making it theoretically possible to estimate SSM at the nominal resolution of Sentinel-1 of 10 m. However, to decrease the well-known noise of radar observations due to interference from the coherent sum of scattered signals, also called "speckle noise," the results are generally provided at a lower resolution.
Bazzi et al. [11] took advantage of the difference in footprint between rainfall that waters a large area (>1 km) and an irrigation event whose impact is limited to the field by comparing the dynamic of the radar signal at both scales. This work has been updated [12] with a new study over summer crops, where the main findings indicated that the density of available Sentinel-1 affects the irrigation detection accuracy, especially for temperate areas. Le Page et al. [13] compared the relative change of surface soil moisture predicted by a simple soil water budget model, with no representation of irrigation and derived from Sentinel-1 [14], on maize plots to determine the timing of irrigation events. Ouaadi et al. [15] proposed an approach based on the assimilation of high-resolution soil moisture products [16] into a simple soil water budget model and demonstrated that both irrigation timing and water amounts could be retrieved with reasonable accuracy if the irrigation techniques (flooding or drip irrigation) were known. A recent study of Zappa et al. [4] analyzed the accuracy of irrigation retrieval from space as a function of the characteristics of soil moisture products, and indirectly the irrigation techniques demonstrated that (1) improved metrics are obtained with high spatial and temporal resolution products; (2) better results are obtained on fields where large amounts of water are applied such as for a traditional flooding irrigation method. Brombacher et al. [17] have proposed an alternative approach based on the comparison of evapotranspiration of natural areas to the evapotranspiration of irrigated fields. The difference of evapotranspiration is attributed to irrigation water amounts.
The detection of irrigation events in [13] relied on a relative comparison between two assessments of SSM. The first assessment was obtained from satellite observations, while the second assessment was obtained from a water budget model based on the FAO-56 model [18]. The comparison between the two assessments was conducted at each satellite observation. The method of irrigation timing retrieval was developed and tested on a small dataset of six maize plots over one season in southwest France. Consequently, the examination of the model's performance did not permit a detailed analysis of its effectiveness. The ESA Irrigation+ project dataset, which provides a wide variety of crops, soils, and irrigation methods in different European countries, is much larger and thus used 3 of 25 to evaluate the method with respect to crop type and irrigation method. This study aims at evaluating this method according to the kind of crop, and the type of irrigation. We also have organized the study in order to analyze impact of the configuration of the Sentinel-1 acquisitions (incidence angle, and orbit), and two classes of vegetation development by separating between small NDVI and high NDVI.
In addition to these evaluations, the force-restore approach [19] is tested as an alternative to the original bucket model in order to estimate SSM. Finally, two different strategies for using the observations of SSM are also tested.
The article is organized into four parts. The first part describes the study area and analyzes the forcing and validation datasets. The ESA Irrigation+ dataset is examined to demonstrate the diversity of crops, soils, irrigation methods, and schedules. The weather data, Sentinel-2 NDVI, and Sentinel-1 SSM datasets are then described. The second part outlines the methodology, briefly describing the model and providing the model equations in different annexes (Appendix A). Additionally, the statistics used are described. The third part presents the results, and the final part discusses conclusions and future works.

Dataset and Irrigation Detection Approach
The irrigation detection method relies on a simple water budget based on the FAO-56 model [20] that is forced by precipitation, the reference evapotranspiration ET0, the basal crop coefficient Kc b , the fractional cover derived from Sentinel-2 NDVI at a daily time step [21], and satellite-derived products of surface soil moisture. The performance of the approach is assessed over different plots in Europe described in the following section. A short description of the methodology is provided afterward, and the system of equations is fully described in Appendices A.2 and A.3.

Study Sites and Irrigation Data
The large dataset consists of a time series of irrigation events characterized by their timing and water amount for different plots in Europe ( Figure 1) The data were gathered in five different regions: Niedersachsen and Brandenburg in Germany, Lot and Tarn in France, and Budrio in Italy, within the frame of the ESA Irrigation+ project [5]. According to the Koppen-Geiger climate classification, Niedersachsen, Lot and Tarn are located in a temperate oceanic climate (Cfb); Brandenburg is located in a humid continental climate (Dfb); and Budrio in a humid subtropical climate (Cfa). There are 85 different plots with data acquired over three seasons from September 2017 to September 2019 depending on the plot for a total of 105 complete growing seasons (see Appendix A.5 for details). The size of the plots ranges from 0.3 to 10 ha. The soil texture has been retrieved from the OpenLandMap dataset using Google Earth Engine (GEE), a dataset based on machine learning predictions from a global compilation of soil profiles and samples. It is an enhancement of SoilGrids from Hengl et al. [22]. Tarn and Lot have clay loam soils (~30% clay and~30% sand); Budrio has silty clay loam soils. The two German regions have more sandy textures: sandy loam for Brandenburg and loam for Niedersachsen [23] (Figure 1). Table 1 describes the number of available plots/seasons by region and crop type. Most of the crops are cereals, consisting of maize, rye, wheat, and barley. The second most important crops in the dataset are tubers (potato and sugar beet). There are also some vegetables (mainly soybean), and crops used for fodder (maize, rapeseed, grass). Plots seeded with walnut and tobacco complete the dataset. However, as full descriptions (eg. Kc b , Fc, rooting depth) of some of those crops were unknown, they were given the configuration of a similar crop (for example, wheat for barley or maize for rye).  The irrigation characteristics of the five study areas are shown in Figure 2 in terms of amount, number per year, and frequency. Note that the amount of water per event is not known for Brandenburg, and was set at 15 mm/event. The average number of irrigation events per season is between 7 and 9 for the regions of Brandenburg, Tarn, and Lot, where sprinkler irrigation dominates. Budrio, with an average of 18 events per season, has a mix of sprinkler and drip irrigation. On average, Niedersachsen has a high number of irrigation events (24 on average, 80 maximum) with an average amount of 18 mm, but with a big variation of the amount between fields. As such, this area is dominated with drip irrigation but also has sprinkler irrigated fields.
Le Page et al. [13] showed that low amounts (less than 10 mm per event) make irrigation detection difficult. In this dataset, only a few events, located mainly in Nieder-  The irrigation characteristics of the five study areas are shown in Figure 2 in terms of amount, number per year, and frequency. Note that the amount of water per event is not known for Brandenburg, and was set at 15 mm/event. The average number of irrigation events per season is between 7 and 9 for the regions of Brandenburg, Tarn, and Lot, where sprinkler irrigation dominates. Budrio, with an average of 18 events per season, has a mix of sprinkler and drip irrigation. On average, Niedersachsen has a high number of irrigation events (24 on average, 80 maximum) with an average amount of 18 mm, but with a big variation of the amount between fields. As such, this area is dominated with drip irrigation but also has sprinkler irrigated fields. sachsen and Budrio, are below this threshold. Likewise, the detection was also shown to be difficult if the frequency of irrigation is higher than the time of the satellite overpass (6 days for the combination of Sentinel-1A and Sentinel-1B), such as in the Niedersachsen region where the average irrigation frequency is 2.8 days.   Figure 4 presents the season of irrigation for the main crop types. The season duration can extend for as long as 150 days in Brandenburg and be as short as one day at the French locations (one irrigation event only). The irrigation season begins around mid-April and ends around mid-September on average, except for in France where the irrigation season is very short (from June to mid-September) because maize, irrigated during the summer months, is the dominant crop in this region. Regarding the irrigation season of the main crops investigated in this study, wheat is generally irrigated from May to June, potato from May to August, and tomato from May to September.  Le Page et al. [13] showed that low amounts (less than 10 mm per event) make irrigation detection difficult. In this dataset, only a few events, located mainly in Niedersachsen and Budrio, are below this threshold. Likewise, the detection was also shown to be difficult if the frequency of irrigation is higher than the time of the satellite overpass (6 days for the combination of Sentinel-1A and Sentinel-1B), such as in the Niedersachsen region where the average irrigation frequency is 2.8 days. Figure 3 shows the distribution over time of the irrigation season per region, while Figure 4 presents the season of irrigation for the main crop types. The season duration can extend for as long as 150 days in Brandenburg and be as short as one day at the French locations (one irrigation event only). The irrigation season begins around mid-April and ends around mid-September on average, except for in France where the irrigation season is very short (from June to mid-September) because maize, irrigated during the summer months, is the dominant crop in this region. Regarding the irrigation season of the main crops investigated in this study, wheat is generally irrigated from May to June, potato from May to August, and tomato from May to September. sachsen and Budrio, are below this threshold. Likewise, the detection was also shown to be difficult if the frequency of irrigation is higher than the time of the satellite overpass (6 days for the combination of Sentinel-1A and Sentinel-1B), such as in the Niedersachsen region where the average irrigation frequency is 2.8 days.  Figure 3 shows the distribution over time of the irrigation season per region, while Figure 4 presents the season of irrigation for the main crop types. The season duration can extend for as long as 150 days in Brandenburg and be as short as one day at the French locations (one irrigation event only). The irrigation season begins around mid-April and ends around mid-September on average, except for in France where the irrigation season is very short (from June to mid-September) because maize, irrigated during the summer months, is the dominant crop in this region. Regarding the irrigation season of the main crops investigated in this study, wheat is generally irrigated from May to June, potato from May to August, and tomato from May to September.

Weather Data
The weather data needed to compute the water budget are the daily precipitation and reference evapotranspiration. When available, the rainfall measured by rain gauges is located close to the plot (Germany and Italy). In France, those measurements were not available outside of the irrigation season. For the latter or in the case of missing data, the amounts were provided by the AgERA5 dataset [24], at a spatial resolution of 0.1° × 0.1°. Likewise, while the FAO-56 water budget approach relies on the reference evapotranspiration (ET0, Allen et al. [20]), the potential evapotranspiration (PET) relying on MODIS (MODerate resolution Imaging Spectroradiometer) observations was used for ease of implementation. PET rates from the 8-day aggregated MOD16A2 product at 500 m spatial resolution were extracted through the GEE platform.

Vegetation INDEX, LAND COVER, and Crop Season
NDVI and Land cover are needed to estimate the basal crop coefficient Kcb for each crop and plot. Time series of NDVI from the Sentinel-2A and 2B satellites were retrieved and averaged over the plot using the COPERNICUS/S2_SR collection from GEE. This dataset is atmospherically corrected with the Sen2cor algorithm [25]. Although other algorithms can perform better [26], this processing provides a coherent time series of surface reflectance. The time series were further cloud-masked using the COPERNI-CUS/S2_CLOUD_PROBABILITY collection, which is a GEE precomputed dataset obtained with the s2cloudless algorithm. An automated cloud-detection algorithm is based on a gradient boosting algorithm which is fast and performs well [27]. As some clouds might still persist in the time series, the well-known Savitsky-Golay filter [28] has been applied in order to smooth the times series of NDVI. Figure 5 displays an example of NDVI time series before and after smoothing. The seasonal land cover type used for the relations between NDVI and fractional vegetation cover, and with NDVI and basal crop coefficient, was retrieved directly from the different dataset providers. The duration of the crop season was determined by an ad hoc peak and valley algorithm. The algorithm searches for local minimums and maximums on the smoothed time series of NDVI, and it deletes consecutive minimums or maximums. It then selects the local maximums above 0.6, and only keeps them if there is a duration of at least 90 days between two local maximums. The beginning and end of each season are found as the local minimums in the time series before and after each selected local maximum. An example result is given in Figure 5.

Weather Data
The weather data needed to compute the water budget are the daily precipitation and reference evapotranspiration. When available, the rainfall measured by rain gauges is located close to the plot (Germany and Italy). In France, those measurements were not available outside of the irrigation season. For the latter or in the case of missing data, the amounts were provided by the AgERA5 dataset [24], at a spatial resolution of 0.1 • × 0.1 • . Likewise, while the FAO-56 water budget approach relies on the reference evapotranspiration (ET0, Allen et al. [20]), the potential evapotranspiration (PET) relying on MODIS (MODerate resolution Imaging Spectroradiometer) observations was used for ease of implementation. PET rates from the 8-day aggregated MOD16A2 product at 500 m spatial resolution were extracted through the GEE platform.

Vegetation INDEX, LAND COVER, and Crop Season
NDVI and Land cover are needed to estimate the basal crop coefficient Kc b for each crop and plot. Time series of NDVI from the Sentinel-2A and 2B satellites were retrieved and averaged over the plot using the COPERNICUS/S2_SR collection from GEE. This dataset is atmospherically corrected with the Sen2cor algorithm [25]. Although other algorithms can perform better [26], this processing provides a coherent time series of surface reflectance. The time series were further cloud-masked using the COPERNI-CUS/S2_CLOUD_PROBABILITY collection, which is a GEE precomputed dataset obtained with the s2cloudless algorithm. An automated cloud-detection algorithm is based on a gradient boosting algorithm which is fast and performs well [27]. As some clouds might still persist in the time series, the well-known Savitsky-Golay filter [28] has been applied in order to smooth the times series of NDVI. Figure 5 displays an example of NDVI time series before and after smoothing. The seasonal land cover type used for the relations between NDVI and fractional vegetation cover, and with NDVI and basal crop coefficient, was retrieved directly from the different dataset providers. The duration of the crop season was determined by an ad hoc peak and valley algorithm. The algorithm searches for local minimums and maximums on the smoothed time series of NDVI, and it deletes consecutive minimums or maximums. It then selects the local maximums above 0.6, and only keeps them if there is a duration of at least 90 days between two local maximums. The beginning and end of each season are found as the local minimums in the time series before and after each selected local maximum. An example result is given in Figure 5.

Surface Soil Moisture "S2SM" Dataset
The S2SM dataset (Sentinel-1 and Sentinel-2 Soil Moisture) used in this study is derived from the approach proposed by El Hajj et al. [14]. It is based on the training of a neural network on a synthetic noisy dataset obtained from the coupling of the Integral Equation Model (IEM) and the Water Cloud Model (WCM, [29]) forced by Sentinel-2 NDVI used as the vegetation descriptor for the WCM. IEM + WCM is run for a wide range of incidence angles, soil roughness, soil moisture, and vegetation index (NDVI). This approach was implemented to produce a dataset of SSM from 2017 to 2020 for the different study areas of the project. It is shown that the best retrieval performance is obtained for NDVI lower than 0.7 and incidence angles close to 39°.
Sentinel-1 A and B acquisition geometries are different between northern (German) and southern (French and Italian) regions. Because of the overlap of different orbits on the site, the plots have about 2.5 acquisitions every 6 days in Germany, while there is around one acquisition per 6 days in France. Finally, the plots in France have an incidence angle closer to the optimal 39° of the S2SM product, while in Germany and Italy, the incidence angles range from 34° to 42° ( Figure 6). A second interesting issue to consider is that the SSM retrievals from S2SM might present a noticeable variability according to the different incidence and azimuth angles of

Surface Soil Moisture "S2SM" Dataset
The S2SM dataset (Sentinel-1 and Sentinel-2 Soil Moisture) used in this study is derived from the approach proposed by El Hajj et al. [14]. It is based on the training of a neural network on a synthetic noisy dataset obtained from the coupling of the Integral Equation Model (IEM) and the Water Cloud Model (WCM, [29]) forced by Sentinel-2 NDVI used as the vegetation descriptor for the WCM. IEM + WCM is run for a wide range of incidence angles, soil roughness, soil moisture, and vegetation index (NDVI). This approach was implemented to produce a dataset of SSM from 2017 to 2020 for the different study areas of the project. It is shown that the best retrieval performance is obtained for NDVI lower than 0.7 and incidence angles close to 39 • .
Sentinel-1 A and B acquisition geometries are different between northern (German) and southern (French and Italian) regions. Because of the overlap of different orbits on the site, the plots have about 2.5 acquisitions every 6 days in Germany, while there is around one acquisition per 6 days in France. Finally, the plots in France have an incidence angle closer to the optimal 39 • of the S2SM product, while in Germany and Italy, the incidence angles range from 34 • to 42 • ( Figure 6).

Surface Soil Moisture "S2SM" Dataset
The S2SM dataset (Sentinel-1 and Sentinel-2 Soil Moisture) used in this study is derived from the approach proposed by El Hajj et al. [14]. It is based on the training of a neural network on a synthetic noisy dataset obtained from the coupling of the Integral Equation Model (IEM) and the Water Cloud Model (WCM, [29]) forced by Sentinel-2 NDVI used as the vegetation descriptor for the WCM. IEM + WCM is run for a wide range of incidence angles, soil roughness, soil moisture, and vegetation index (NDVI). This approach was implemented to produce a dataset of SSM from 2017 to 2020 for the different study areas of the project. It is shown that the best retrieval performance is obtained for NDVI lower than 0.7 and incidence angles close to 39°.
Sentinel-1 A and B acquisition geometries are different between northern (German) and southern (French and Italian) regions. Because of the overlap of different orbits on the site, the plots have about 2.5 acquisitions every 6 days in Germany, while there is around one acquisition per 6 days in France. Finally, the plots in France have an incidence angle closer to the optimal 39° of the S2SM product, while in Germany and Italy, the incidence angles range from 34° to 42° ( Figure 6). A second interesting issue to consider is that the SSM retrievals from S2SM might present a noticeable variability according to the different incidence and azimuth angles of A second interesting issue to consider is that the SSM retrievals from S2SM might present a noticeable variability according to the different incidence and azimuth angles of the ascending pass (6 PM) and the descending pass (6 AM). Figure 7 shows this variability for a plot in Brandenburg where two orbits overlap (red down-triangles and yellow uptriangles). The decrease in SSM at the end of the growing phase is probably related to the lack of wetting events (days 75 to 110). When two wetting events occur around day 125, there is no meaningful increase in SSM. In this case, the well-developed vegetation might prevent the penetration of the radar signal down to the soil. Large fluctuations of SSM occur, and this is encouraging for the detection of wetting events even with well-developed vegetation. In addition, smaller fluctuations from one acquisition to another are noticeable. The peaks are not necessarily related to irrigation or irrigation events but might be related to the different orbits and the associated change in incidence or azimuth angles, and time of view. For example, morning dew might be identified as a higher SSM. The detection algorithm must deal with those small fluctuations. the ascending pass (6 PM) and the descending pass (6 AM). Figure 7 shows this variability for a plot in Brandenburg where two orbits overlap (red down-triangles and yellow up-triangles). The decrease in SSM at the end of the growing phase is probably related to the lack of wetting events (days 75 to 110). When two wetting events occur around day 125, there is no meaningful increase in SSM. In this case, the well-developed vegetation might prevent the penetration of the radar signal down to the soil. Large fluctuations of SSM occur, and this is encouraging for the detection of wetting events even with well-developed vegetation. In addition, smaller fluctuations from one acquisition to another are noticeable. The peaks are not necessarily related to irrigation or irrigation events but might be related to the different orbits and the associated change in incidence or azimuth angles, and time of view. For example, morning dew might be identified as a higher SSM. The detection algorithm must deal with those small fluctuations. The main graph shows the rainfall and irrigation events as black and blue bars, respectively. The NDVI is drawn with a dark-blue dashed line when lower than 0.75 and red dashes when higher than 0.75. The S2SM product is drawn with yellow arrows (ascending orbit) and red arrows (descending orbit); Upper panel: the incidence angle of acquisition is shown.

Model and Workflow
Three datasets are needed as inputs in order to obtain the estimation of irrigation as an output. The workflow to do this is summarized in Figure 8. The input time series (Sentinel-1 backscattering data, Sentinel-2 NDVI, and ERA5Land weather data) are shown in light gray boxes. The daily input variables of the models (Kcb, Fc, sowing and harvest date, rainfall, ET0, and S2SM) are shown in dark gray boxes. S2SM (1) has been previously described in [14]; the way to estimate Kcb and Fc (2) is described in [30]; the estimation of sowing and harvest dates (3) is described in Section 2.3. The computation of ET0 (4) is performed according to [20], and the water budget (5) is calculated following [21]. The equations of the water budget are detailed in Annex 2. The two alternative SSM water budgets are then computed according to [13] and [15]. The equations of the modified bucket approach are detailed in Annex 3. The equations of the force-restore model are detailed in Annex 4. As described in Section 3.2, S2SM is either used with the actual observations or is interpolated (7). Irrigation detection (8) and exact date choice (9) are The main graph shows the rainfall and irrigation events as black and blue bars, respectively. The NDVI is drawn with a dark-blue dashed line when lower than 0.75 and red dashes when higher than 0.75. The S2SM product is drawn with yellow arrows (ascending orbit) and red arrows (descending orbit); Upper panel: the incidence angle of acquisition is shown.

Model and Workflow
Three datasets are needed as inputs in order to obtain the estimation of irrigation as an output. The workflow to do this is summarized in Figure 8. The input time series (Sentinel-1 backscattering data, Sentinel-2 NDVI, and ERA5Land weather data) are shown in light gray boxes. The daily input variables of the models (Kc b , Fc, sowing and harvest date, rainfall, ET0, and S2SM) are shown in dark gray boxes. S2SM (1) has been previously described in [14]; the way to estimate Kc b and Fc (2) is described in [30]; the estimation of sowing and harvest dates (3) is described in Section 2.3. The computation of ET0 (4) is performed according to [20], and the water budget (5) is calculated following [21]. The equations of the water budget are detailed in Appendix A.2. The two alternative SSM water budgets are then computed according to [13] and [15]. The equations of the modified bucket approach are detailed in Appendix A.3. The equations of the force-restore model are detailed in Appendix A.4. As described in Section 3.2, S2SM is either used with the actual observations or is interpolated (7). Irrigation detection (8) and exact date choice (9) are performed according to [13]. Finally, if necessary, a merging (10) of the different runs is performed.
Remote Sens. 2023, 15, 1449 9 of 26 performed according to [13]. Finally, if necessary, a merging (10) of the different runs is performed. , and, on the other side, the variables that will allow the computing of a water budget for the estimation of the surface soil moisture (basal crop coefficient (Kcb), Fraction cover (Fc), sowing and harvest date, precipitation and reference evapotranspiration (ET0)). Steps 5 to 9 are the irrigation detection steps. The two diamonds (6 and 7) indicate sub-methods for dealing with soil moisture observation and soil moisture estimation, respectively. The final step (10) consists of merging results.
The root layer has been partitioned into an upper soil layer and root zone layer according to Raes [31]. If there is a significant increase in S2SM and a significant decrease in the predicted SSM (SSMmodel), this might be due to an irrigation event during the two soil moisture retrievals. The approach evaluates the variation in SSM of the model and satellite-derived products separately. Indeed, there are many uncertainties on both sides. On the model side, uncertainties can be explained by the input data, in particular the precipitation, but also the hydraulic characteristics of the soil. On the side of the satellite product, changes in soil roughness or vegetation density can significantly affect the accuracy of the product. Stated differently, more weight is given to the relative variations than to the absolute values of SSM. The model then searches for the best date of irrigation between the two S2SM observations. The method is fully described in Le Page et al.
(2020). The approach was applied for six sprinkler-irrigated maize plots in south-west France, where Sentinel-1 provided an observation every six days. The results were very encouraging with an average F-score of 0.69 (see below for the description of the F-score).
As an alternative to the soil bucket model of the FAO-56, a force-restore approach [32,33] is tested within this study, as it has been shown to better describe surface and root-zone coupling and, in particular, the capillary rise that is not considered in FAO-56. In addition, this approach resolves the daily cycle of soil moisture evolution allowing for a comparison of predicted and retrieved soil moisture at the time of the satellite overpass , and, on the other side, the variables that will allow the computing of a water budget for the estimation of the surface soil moisture (basal crop coefficient (Kcb), Fraction cover (Fc), sowing and harvest date, precipitation and reference evapotranspiration (ET0)). Steps 5 to 9 are the irrigation detection steps. The two diamonds (6 and 7) indicate sub-methods for dealing with soil moisture observation and soil moisture estimation, respectively. The final step (10) consists of merging results.
The root layer has been partitioned into an upper soil layer and root zone layer according to Raes [31]. If there is a significant increase in S2SM and a significant decrease in the predicted SSM (SSMmodel), this might be due to an irrigation event during the two soil moisture retrievals. The approach evaluates the variation in SSM of the model and satellitederived products separately. Indeed, there are many uncertainties on both sides. On the model side, uncertainties can be explained by the input data, in particular the precipitation, but also the hydraulic characteristics of the soil. On the side of the satellite product, changes in soil roughness or vegetation density can significantly affect the accuracy of the product. Stated differently, more weight is given to the relative variations than to the absolute values of SSM. The model then searches for the best date of irrigation between the two S2SM observations. The method is fully described in Le Page et al. (2020) [13]. The approach was applied for six sprinkler-irrigated maize plots in south-west France, where Sentinel-1 provided an observation every six days. The results were very encouraging with an average F-score of 0.69 (see below for the description of the F-score).
As an alternative to the soil bucket model of the FAO-56, a force-restore approach [32,33] is tested within this study, as it has been shown to better describe surface and root-zone coupling and, in particular, the capillary rise that is not considered in FAO-56. In addition, this approach resolves the daily cycle of soil moisture evolution allowing for a comparison of predicted and retrieved soil moisture at the time of the satellite overpass in contrast to the daily prediction of the FAO. An example of the trajectory of the two models is shown in Figure 9. The daily SSM estimate dries faster than the hourly SSM which takes into account capillary rise. Compared to the bucket model, which partitioned the root layer into an upper soil layer and root zone layer according to [15], this model should reach the following goals:

•
Obtain SSM for a thinner (10 cm) and fixed depth of soil; • Obtain an hourly estimate of SSM that can be used for the times of the ascending and descending orbits; • Take into account the capillary rise to the upper soil layer. in contrast to the daily prediction of the FAO. An example of the trajectory of the two models is shown in Figure 9. The daily SSM estimate dries faster than the hourly SSM which takes into account capillary rise. Compared to the bucket model, which partitioned the root layer into an upper soil layer and root zone layer according to [15], this model should reach the following goals: • Obtain SSM for a thinner (10 cm) and fixed depth of soil; • Obtain an hourly estimate of SSM that can be used for the times of the ascending and descending orbits; • Take into account the capillary rise to the upper soil layer. For each plot, the start and the end of the season are estimated as previously stated. An empirical linear relation between NDVI and Fc and NDVI and Kcb is used according to each field land cover in the same way as [13]. Because of the lack of in situ information, the maximum rooting depth is set to 60 cm independently of location or crop. For each single plot, the irrigation amount is set to the average known irrigation of the season. This is compliant with the original objective of this algorithm which was to retrieve irrigation events automatically at the plot scale for a farmer within the framework of a decision support system for irrigation planning. Wilting point and field capacity of the soil are computed with the Saxton equations [34] (Appendix A) according to the clay and sand content derived from the Openlandmap dataset (see Annex 1 for equations).

Implementation and Performance Assessment
In order to better exploit the wealth of Sentinel-1 observations, two strategies are deployed: 1. The first approach aims to retrieve irrigation timing separately for each Sentinel-1 viewing configuration (i.e., with different azimuth and incidence angles). This ensemble of irrigation timing retrievals are then merged by keeping one irrigation event within a delay of ± n days only. Different values of n (1, 2, and 3 days) were tested, which means that the time period containing only one irrigation event would be 3, 5, or 7 days. Note that as n increases, the number of possible detection of irrigation events mechanically lowers. Large n values thus degrade the detection of high irrigation rate. 2. The second approach consists in interpolating the S2SM estimates before irrigation detection. The time series of SSM are linearly interpolated to the daily time step regardless of the acquisition configuration. The linearly interpolated time series is then smoothed with the Savitsky-Golay algorithm (window of seven days and For each plot, the start and the end of the season are estimated as previously stated. An empirical linear relation between NDVI and Fc and NDVI and Kc b is used according to each field land cover in the same way as [13]. Because of the lack of in situ information, the maximum rooting depth is set to 60 cm independently of location or crop. For each single plot, the irrigation amount is set to the average known irrigation of the season. This is compliant with the original objective of this algorithm which was to retrieve irrigation events automatically at the plot scale for a farmer within the framework of a decision support system for irrigation planning. Wilting point and field capacity of the soil are computed with the Saxton equations [34] (Appendix A) according to the clay and sand content derived from the Openlandmap dataset (see Appendix A.1 for equations).

Implementation and Performance Assessment
In order to better exploit the wealth of Sentinel-1 observations, two strategies are deployed: 1.
The first approach aims to retrieve irrigation timing separately for each Sentinel-1 viewing configuration (i.e., with different azimuth and incidence angles). This ensemble of irrigation timing retrievals are then merged by keeping one irrigation event within a delay of ± n days only. Different values of n (1, 2, and 3 days) were tested, which means that the time period containing only one irrigation event would be 3, 5, or 7 days. Note that as n increases, the number of possible detection of irrigation events mechanically lowers. Large n values thus degrade the detection of high irrigation rate.

2.
The second approach consists in interpolating the S2SM estimates before irrigation detection. The time series of SSM are linearly interpolated to the daily time step regardless of the acquisition configuration. The linearly interpolated time series is then smoothed with the Savitsky-Golay algorithm (window of seven days and second degree polynomial) in order to remove small peaks of SSM. As the best time step was previously identified to be three days for irrigation timing detection [13], 3-day SSM are used for irrigation retrieval.
The results are analyzed as follows: • The incidence angle and the acquisition geometry: in the case of separate runs (approach 1), the impact of orbit and incidence angle can be analyzed. When two orbits are available such as over the German plots, four runs will be considered (AAC The soil water budget model: In both approaches, the force-restore model is compared to the original bucket scheme.

•
The merging of separate runs: The impact of the merging delay n, within which it is assumed that an irrigation event can occur, is analyzed with values ranging from 1 to 3 days.

Calculation of Scores
The performance of the retrieval of irrigation timing can be considered as a classification problem. If the simulated date of irrigation is equal to an observed date of irrigation, the retrieval can be qualified as perfect. However, this definition is extremely restrictive, so that in the calculation of the score, we have considered that an irrigation event is welldetected (true positive) if it is located within three days of an actual irrigation (observed event). This definition is still very strict because it is important to have good accuracy of irrigation timing in order to eventually replace actual irrigation events with those that are simulated irrigation events. The score itself is computed with the F-score (or F1-score, or Fmeasure). F-score is computed as the harmonic mean of precision and recall (Equation (1)). Precision expresses how many retrieved items are relevant. It is computed by the number of true positive results divided by the number of all positive results. Recall expresses how many relevant items have been retrieved. It is computed as the number of true positive results divided by the number of all samples that should have been identified as positive. An F-score of 1 means a perfect simulation, while an F-score of zero is obtained when either precision or recall are zero.
F-score does not necessarily give a good idea about seasonal performance. To do so, we simply compared the seasonal sum of irrigation amounts of each event (Equation (2)).

Seasonal Bias
where n is the number of irrigation events simulated, amount is the average observed amount of water per event (mm) during the season, m is the actual number of observed irrigation events, and amount is the observed amount of irrigation (mm) of each observed irrigation event.

Results with Separate Runs
This sub-section examines the results obtained with separate runs. The different satellite and field configurations are analyzed, and the two soil moisture models are compared. The final result is obtained by merging the different runs, which are analyzed in terms of merging intervals.

Performance according to Different Configurations
In this section, the performance of the irrigation method is analyzed as a function of the Sentinel-1 acquisition configurations (orbit and incidence angle), the irrigation frequency, the crop types, and the region. The crop stage is also separated into two parts: low NDVI and high NDVI. For the latter, a threshold of 0.75 [14] is applied to NDVI value to separate irrigation retrieval between small NDVI (smallNDVI) and large NDVI values (largeNDVI). In the smallNDVI analysis, irrigation events are retrieved when NDVI < 0.75, and observed irrigation events are used when NDVI > 0.75. Inversely, in the largeNDVI analysis, actual irrigation events are used in the model when NDVI < 0.75, while the retrieval method is applied when NDVI > 0.75. In this way, the computed soil moisture is not disturbed by previous time step estimates of irrigation. The statistics that are presented in Table 2 are obtained by averaging the statistics of the different runs according to different configurations. For example, if there are four different configurations of the Sentinel-1 acquisition, the general F-score is the average of the four F-scores obtained for each run, and the "ascending" and "descending" F-scores concern two runs each. In Table 3, the orbit is separated into ascending and descending modes, and the incidence angles in two classes: a range (37-41 • ) corresponding to the optimal angle of the S2SM (39 • ) and all other acquisitions (<37 • or >41 • ).  A slightly better performance (3 points of the F-score) is achieved using soil moisture retrieved from the descending orbits (around 6 AM) versus those obtained with the ascending (6 PM) orbit. F-scores are similar whatever the incidence angle, except that negative bias is higher for data outside of the range 37-41 • .
The scores decrease as the frequency of irrigation increases. It is particularly low (F-score around 20) for very frequent irrigation events (around 1 event/3 days). This low F-score is also associated with a significant negative bias. This is corroborated by the comparison between France and Germany. Indeed, the plots from Germany (in particular, Niedersachsen) are likely drip-irrigated with high irrigation frequency while the plots from France, dominated by Maize (36 among 54 plots), are more commonly sprinkler irrigated. Interestingly, the seasonal amounts of water for France for low and average irrigation frequency are close to the observed seasonal amounts. As potatoes are located solely in Germany and with frequent irrigation while the French plots are mainly cropped with maize characterized by low irrigation frequencies, the F-score by crops must be treated with caution. The better performance for maize compared to potato or wheat is probably due to a lower irrigation frequency.
The biases are always negative apart from plots with low-frequency irrigation (and, to a lesser extent, average) frequency (sprinkler). Interestingly enough, the best scores are obtained with higher values of NDVI, as the penetration of the radar signal is decreased by dense vegetation, and thus the accuracy of the SSM product is reduced [14]. Nevertheless, the number of irrigation events when the crops are well developed (i.e., with high NDVI values) largely outweighs the number of irrigation events outside this period.

Merge of Separate Runs
This is the final step using the separate runs. Table 3 presents the scores of merging with 1-day, 2-day, and 3-day intervals. The merging was only performed with the original bucket model. The results are presented according to the crop type, the country, and the irrigation frequency.
The best scores are obtained with a 3-day interval (1 day before and 1 day after the central day). The precision using different intervals hardly changed, but 1-day improved the recall, in particular for the wheat crop and the very frequent irrigation frequencies.
The mean F-score varies between 43 for low and very frequent irrigation to 56 for average frequency. Regarding seasonal amounts, the score is −16% for maize and wheat with an interval of 1 day, while it reaches −57% for the potato crop. Plots with low irrigation frequency exhibit an overestimation of 36%, which means there was some over-detection. Average and frequent irrigation events range from −17 to −24%, while there is a significant underestimation obtained for the very high irrigation frequencies (53%). Better F-scores were obtained in Germany, but with a higher negative bias. This is correlated with the lowest scores obtained over maize which is the dominant crop in the French dataset.

Impact of the Force-Restore Model
As stated before, the force-restore approach has been implemented to represent the daily cycle of soil moisture in order to improve the comparison of the model prediction with the satellite retrieval of SSM. Following Section 4.1, the mean of n runs according to acquisition configuration is considered. Table 4 summarizes the results as a function of irrigation frequency. Compared to the bucket approach, the F-score increases by about six points for low and average irrigation frequencies. There is also a small increase for frequent and very frequent irrigation event frequency. The seasonal sums are a bit better for frequent irrigation events. However, the scores are lower for the smaller NDVIs. Table 5 shows the differences in precision, recall, and F-score for the different configurations observed before. It must be noted that the improvement in F-score is obtained through a lower precision but a better recall. As a conclusion, a force-restore soil water budget improved the irrigation event retrieval with regards to the bucket model with a small loss of precision.

Results with Interpolated S2SM
In this section, SSM data were interpolated and smoothed regardless of acquisition configuration. Then, only one out of each three days of this time series was kept in order to be congruent with the finding of [13]. Table 6 presents the F-score and bias for different crop types, countries, and irrigation frequencies. The F-scores are very similar to the results of the merging approach but seem more stable from one region to another, and from one crop to another. The F-score values obtained in Germany are consistent with those found in [23], where the authors used a different Sentinel-1 soil moisture product over the same fields. A higher range of biases according to irrigation frequency is also observed with a large overestimation for low-frequency irrigation and an underestimation for high-frequency irrigation. This is corroborated by the positive biases on maize which is sprinkler irrigated. This approach seems to be rather sensitive to irrigation frequency. The counter-intuitive results already described regarding the NDVI values are also observed, as better scores are obtained with largeNDVI than with smallNDVI. Table 6. Scores (F-score, bias) of the interpolated approach according to different configurations orbit, incidence angle, crop, country, irrigation frequency), and the level of NDVI (all data, small NDVI values, high NDVI values). The force-restore model improves the performance in the same way shown above with a small decrease in precision and a good increase in recall which gives an improvement of the F-score of about 5 points. In all cases, the bias is reduced compared to the bucket model. Figure 10 shows a comparison of the three methods: interpolated SSM with the original model (interp-Original), interpolated SSM with force restore model (interp-FR), and merged results of different runs (Merge-FR). The merged approach gives better results for very frequent irrigation, while the results are generally better with the interpolated-FR approach for low and average irrigation frequency. This is true both for F-scores and seasonal sums of irrigation. As such, in France, the interp-FR supersedes the merge-FR method by 9.7 points on the F-score, while in Germany, the merge-FR supersedes the interp-FR by 5.6 points. Regarding seasonal sums of irrigation, interp-FR overestimated them by 43.3% in France, while Merge-FR underestimated them by −18%. However, the better performance for very frequent irrigation in Germany might also be due to a side-effect of the four Sentinel-1 acquisitions in this region.

Some Examples of Irrigation Timing Retrieval
In this last section, some examples of good and bad retrievals have been selected. Figure 11 shows an example of low irrigation frequency on winter rye in Brandenburg in 2017. Apart from the first irrigation, the four last events are correctly detected. Figure 12 shows an example of frequent irrigation on a potato crop in Brandenburg in 2018. About half of the irrigation events are correctly detected; a quarter of irrigation events are detected more than three days from the actual events; and about a quarter of the events are not detected. Figure 13

Some Examples of Irrigation Timing Retrieval
In this last section, some examples of good and bad retrievals have been selected. Figure 11 shows an example of low irrigation frequency on winter rye in Brandenburg in 2017. Apart from the first irrigation, the four last events are correctly detected. Figure 12 shows an example of frequent irrigation on a potato crop in Brandenburg in 2018. About half of the irrigation events are correctly detected; a quarter of irrigation events are detected more than three days from the actual events; and about a quarter of the events are not detected. Figure 13 illustrates the case of a drip-irrigated potato plot in Niedersachsen, in 2018, where irrigation events are very frequent with small amounts. In this difficult case, the algorithm underestimates the number of irrigation events. The high permanent level of SSM during the full development stages of vegetation could be an indicator of the irrigation method (localized versus sprinkler or flood).

Conclusions
A new dataset of irrigation events was used to analyze the retrieval of irrigation timing at the plot level using satellite observations. The Irrigation+ dataset was challenging due to the variety of crops, locations, soils, and irrigation methods, with most plots using localized irrigation, which is not well suited for the retrieval method designed for low-pace irrigation (flooding, sprinkler). Two strategies were tested: one that considered SSM observations separately according to the Sentinel-1 acquisition configura- Figure 13. Same as Figure 12 for a potato field in Niedersachsen in 2018 (high irrigation frequency and low amounts of irrigation).

Conclusions
A new dataset of irrigation events was used to analyze the retrieval of irrigation timing at the plot level using satellite observations. The Irrigation+ dataset was challenging due to the variety of crops, locations, soils, and irrigation methods, with most plots using localized irrigation, which is not well suited for the retrieval method designed for lowpace irrigation (flooding, sprinkler). Two strategies were tested: one that considered SSM observations separately according to the Sentinel-1 acquisition configuration and then merged the different runs, and another that merged and interpolated all observations and ran the model once with a constant revisit time. The impact of using a force-restore approach instead of the original bucket soil water budget model on both strategies was also analyzed.
Performance was evaluated using the F-score and the seasonal sum of irrigation. The first strategy allowed the analysis of the impact of the Sentinel-1 incidence angle and orbit on irrigation retrieval. The scores were low, especially with only a six-day frequency of observation, and there was no significant difference between the Sentinel-1 configurations (ascending/descending, close to 39 • /far from 39 • ). The lower scores obtained with small NDVI compared to large NDVI were almost systematic, which is counter-intuitive. However, as the score increased when taking all NDVI values, this might have been due to the larger number of irrigation events during high vegetation periods. Interestingly enough, it seems that this algorithm does not suffer a degradation of the score due to the known degradation of the SSM assessment associated with the development of vegetation. The scores decreased as irrigation frequency increased, which was substantiated by the fact that the scores were better in France (more sprinkler irrigation) than in Germany. Finally, the strategy of merging different runs versus the strategy of interpolating all SSM data for one run has produced very comparable results. However, compared to each separate run, the best scores were obtained for average and high frequencies of irrigation. The cumulative sum of irrigation was around −20% in those cases. The replacement of the bucket SSM model by the force-restore gave an improvement of about six points on the F-score, and also narrowed the error on cumulative seasonal irrigation.
The present analysis has revealed some important insights. Firstly, the frequency of irrigation has a significant impact on the detection of irrigation events. High-frequency localized irrigation events are difficult to detect due to several reasons, such as partial wetting of soil, wetting beneath the canopy, and small amounts leading to minimal variations of SSM that could fall within the error range of the SSM product. Moreover, for drip irrigation, detecting irrigation events is not relevant for irrigation scheduling, and a weekly assessment of root zone soil moisture using, for instance, an energy budget-based approach would be more appropriate. Alternatively, adjusting the thresholding parameters based on the type of irrigation method, as performed in [15], could be considered. Although the seasonal pattern of SSM may offer some indication of the irrigation method (flooding/sprinkler/local), it can only be analyzed retrospectively. Recently, Paolini et al. [35] employed a supervised machine learning approach, utilizing different earth observation products to classify the irrigation technique at the field scale. This could be a useful preliminary step before attempting to identify irrigation events.
Secondly, it is necessary to take into account the field context more comprehensively. In this work, the whole crop season is taken, when actually the irrigation season is generally shorter. Examining the irrigation calendars by crop type, as shown in Figures 4 and 5, would provide a probability estimation of irrigation, which would eliminate detection of out-of-season irrigation. This probability map, which could include information on crop type and region, would be used to qualify detected irrigation events. The MIRCA2000 dataset [36] could be used to derive this probability map. The approach proposed by Bazzi et al. [11,12] can also be considered a contextual approach. An irrigation event is declared when the high-resolution variation of SSM (the field) is not observed on the low-resolution variation (5 km).
Thirdly, it was observed that excessive saturation or attenuation of the radar signal adversely affects SSM estimation. Although the current algorithm accounts for this uncertainty, it would be beneficial to consider the confidence of irrigation detection based on vegetation development, such as using NDVI. Furthermore, it would be worthwhile to evaluate retrieval performance using other high-resolution soil moisture datasets that will likely be available in the future. Additionally, combining SAR-estimated SSM with opticalbased indices of surface water content, such as an NDWI based on medium-infrared and near-infrared [37], could be of interest. For future research, it is recommended to not only consider SSM but also the status of root zone soil moisture, particularly for low-frequency irrigation as it is likely that irrigation events will occur when the readily available root zone soil moisture is between 50 and 70%. Funding: This research was funded by European Space Agency (contract 4000129870/20/I-NB, project "Irrigation+") and Agence Nationale de la Recherche (ANR-19-CE01-0017, project "HILIAISE").

Data Availability Statement:
The irrigation dataset for Italy has been provided by CER (Canale Emiliano Romagnolo). The dataset for Germany has been provided by the Leibniz Centre for Agricultural Landscape Research (ZALF). The dataset for France has been provided separately by the Chambre d'Agriculture of the Tarn, and the Chambre d'Agriculture of the Lot under a special agreement which does not allow to distribute the data; however, we will provide all necessary connections to reach the data providers. The meteorological datasets were provided by the national rainfall network, managed by the Department of Civil and Environmental Protection (DPC) of Italy for the Italian pilot areas and by the German Weather Service (DWD) Climate Data Center (CDC).

Acknowledgments:
The authors acknowledge the support from the European Space Agency (ESA) under the IRRIGATION+ project (contract n. 4000129870/20/I-NB); for further details please visit https://esairrigationplus.org/, accessed on 1 January 2022. We acknowledge the Chambre d'Agriculture of Tarn, and Lot in France for providing their dataset of irrigation. We wish to thank CER (Canale Emiliano Romagnolo) for providing the irrigation amounts over the Italian pilot site. We also thank Claas Nendel and Beate Zimmerman for providing the Germany irrigation data.

Conflicts of Interest:
The authors declare no conflict of interest.  (2) Horizontal partition according to fraction cover and fraction wet First, the horizontal partition of the upper layer is carried on, comparing the fraction wet (fw) to the fraction cover derived from NDVI (Fc sat ). The soil wetted by both irrigation and precipitation (Few i ) and the soil wetted by precipitation (Few p ) are computed as follows: Rooting depth is actualized when maximum root depth is different from minimum root depth, in fact for annual crops. Rooting Depth is derived from a linear relation to Fraction Cover which assumes that the crop will reach its full rooting at maximum fraction cover.
where Zr j = Root Depth at day j, with Zrj > Ze+0.01. Fcsat j = Fraction Cover of the day j. Fcsat max = Maximum Fraction Cover of the time series. MaxZr(LC) = parameter for Max Root Depth of this land cover (LC). As roots may enter different soil horizons, average wilting point and field capacity are taken from the pre-processed soil map mentioned earlier. Depletions are then actualized in the deep (Dd) and root (Dr) bucket for each day j.
Readily Available Water at day j (RAW j ) is actualized with the equation from FAO-56, where EThb is described in part 5: RAW j = RAW j−1 + (TAW * p + 0.04 * 5 − EThb j−1 (4) Water balance, part 1 The water balance proceeds to input water from Irrigation (I) which either comes from actual or simulated irrigation and rainfall (R). The depletion of the evaporation bucket wetted with rain only (Dep), the bucket wetted with rain and irrigation (De), and the root bucket (Dr) are actualized: Stresses are calculated. They represent the stress after receiving water input, and after root actualization induced from the EO of the day. As remote sensing images are generally shot at mid-day, it appears to be a good compromise. Kr, the reduction coefficient for evaporation, is expressed for the fewi and fewp fractions as A weighting coefficient (W) for partitioning the energy available for evaporation in the two wetted, exposed fractions of the surface layer, depending on water availability is calculated:    [31], 40% of the transpiration comes from the first upper quarter of soil (Zr/4). The main bucket is divided into two layers, with the upper layer having a depth equal to one-quarter of the root layer. Dr top is the depletion of water in the upper layer. TAW top and TAW down correspond to the total available water content in the upper and lower layers, and are computed in the same way as TAW in the FAO-56 model, but are fractioned into the upper (25%) and lower (75%) percentages. The depletion of the bottom layer Dr down is computed as the difference in depletion between the full root layer and the top root layer. TAW down = θ f c − θ wp * Zr * 0.75 In a second step, Dr top is updated in accordance with known wetting events, and depletion is constrained between 0 and the maximum water content TAW top : The stress coefficient of the top layer (Ks top ) was computed in the exactly the same way as K s in the FAO-56 model, but with the water quantities of the top layer, where p adjust is the average fraction of TAW that can be depleted from the root zone before water stress occurs. The transpiration of the top layer (T top ) was then obtained with the basal crop coefficients Kc b and previously computed Ks top . In accordance with the findings of Raes [31], 40% of the transpiration is affected by the top layer. T top cannot be greater than the total amount of transpiration T computed with the regular FAO-56 model: The force-restore model receives the volumetric water content of the evaporation and root compartment, and precipitation, evaporation, and evapotranspiration daily fluxes from the FAO-56 model. The soil parameters and coefficients are computed according to [32]: θ sat , the saturated volumetric water content (m 3 m−3), θ wp , the wilting point volumetric water content (m 3 m−3), and θ f c the field capacity volumetric water content (m 3 m−3), the slope of the retention Curve (b), and the coefficients c2 re f and c1 sat at saturation. Then, as we are moving from a daily time step to a shorter time step, some unit conversions are needed. Fluxes of precipitation and evaporation and evapotranspiration are converted from mm.day−1 to mm.s−1. For each time step of one hour (tstep = τ 24 , with τ = 86400) in one day, the force restore model is applied. The total soil depth (d2) is set to 0.60 m, and the evaporation depth (d1 of 0.1 m) is implicitly included in the c1sat regression equation.
First, the hourly evapotranspiration (Etr) and evaporation (Eg) are obtained from the daily values, based on the cosines of the zenith angle (coszenith). The coszenith function is the same as the one described in [20] to compute the zenithal angle for a position (lat, lon) and hour at some Julian day. Etr = Etr daily * coszenith Eg = Eg daily * coszenith The variation of soil moisture is updated according to the force-restore approach from [33] derived from [19], where Pge is the hourly precipitation, which in this case supposes that the precipitation is evenly distributed over the whole day (Pg = P/24): x = θ 2 /θ sat θ eq = θ sat * x − a * x p * 1 − x 8p c2 = c2 re f * (w2/(θ sat − w2 + 1e −10 c1 = c1 sat * (θ sat /max θ g , θ wp