Implementation and validation of a new irrigation scheme in the ISBA land surface model

With an increase in the number of natural processes represented, global land surface models (LSMs) have become more and more accurate in representing natural terrestrial ecosystems. However, they are still limited, especially in the representation of the impact of agriculture on land surface variables. This is particularly true for agro-hydrological processes related to a strong human control on freshwater. While most LSMs consider natural processes only, the development of human-related 15 processes, e.g. crop phenology and irrigation in LSMs, is key. In this study we present the implementation of a new irrigation scheme in the ISBA (Interaction between Soil, Biosphere, and Atmosphere) LSM. This highly flexible scheme is designed to account for various configurations and can be applied at different spatial scales. For each vegetation type within a model grid cell, three irrigation systems can be used at the same time. A limited number of parameters are used to control (1) 20 the amount of water used for irrigation, (2) irrigation triggering (based on the soil moisture stress) and (3) crop seasonality (emergence, harvesting). After a presentation of the simulations of the new scheme at a plot scale, an evaluation is proposed over Nebraska (USA). This region is chosen for its high irrigation density and because independent observations of irrigation practices can be used to verify the simulated irrigation amounts. The ISBA simulations with and without the irrigation scheme are 25 compared to different satellite-based observations. The comparison shows that the irrigation scheme improves the simulated vegetation variables such as leaf area index and gross primary productivity and other variables largely impacted by irrigation such as evapotranspiration and land surface temperature. https://doi.org/10.5194/gmd-2021-332 Preprint. Discussion started: 11 October 2021 c © Author(s) 2021. CC BY 4.0 License.

, thereby precluding inter-annual variability of vegetation density and the impact of irrigation on vegetation growth. Having a more complete irrigation description is needed to follow the irrigation seasonality, and to represent possible changes in crop phenology such as emergence and harvest dates. The impact of changing irrigation characteristics in a context of climate change could 85 thereby be evaluated, such as increasing irrigation efficiency (currently around 56%; FAO, 2014) and freshwater saving potential (Perry et al., 2017;Koech and Langat, 2018).
The objective of this work is to develop and evaluate a more detailed representation of irrigation practices into the ISBA LSM within the SURFEX (SURFace EXternalisée) modelling platform (Masson et al., 2013). In SURFEX, land cover is described by ECOCLIMAP-II (Faroux et al., 2013). 90 This study takes advantage of the ECOCLIMAP-SG (Calvet and Champeaux, 2020; Supplement S1) major update of ECOCLIMAP-II. While the SURFEX framework allows the coupling of terrestrial processes with atmospheric and hydrological models, only offline ISBA simulations are considered in this study.
The evaluation of the new irrigation scheme is made over the state of Nebraska (United States of 95 America, USA). This area presents a high density of irrigated fields (Fig. 1) and large freely available observational datasets for evaluation.
Section 2 presents the observational datasets, the current version of the ISBA LSM, the description of the new irrigation scheme, followed by a description of the validation protocol. Section 3 illustrates the impact of the new irrigation scheme when compared to a model run without irrigation. An 100 evaluation of the performance of the model is made over Nebraska. Section 4 discusses the added value and the limits of the newly implemented irrigation scheme. Finally, section 5 presents the conclusions and future research directions. https://doi.org/10.5194/gmd-2021-332 Preprint. Discussion started: 11 October 2021 c Author(s) 2021. CC BY 4.0 License.

Irrigation map
One of the main challenges of this study was to obtain an upgraded map of irrigation at the global scale, to be consistent with the resolution (300 m × 300 m) of the European Space Agency -Climate Change Initiative (ESA-CCI) land cover map used in ECOCLIMAP-SG. The 1 km × 1 km resolution global irrigation map proposed by Meier et al. (2018), based on a statistical approach and satellite data, was 110 used. A reason to choose this product was that its development process was based (amongst other) on the ESA-CCI land cover product (v1.6.1), the same as the one used to develop the ECOCLIMAP-SG vegetation map (Supplement S1).
In order to transfer the Meier irrigation map (1 km × 1 km) to ECOCLIMAP-SG (300 m × 300 m), a spatial rescaling of the Meier map was performed. A simple majority rule was used by assigning 115 to each 300 m × 300 m grid point of ECOCLIMAP-SG the irrigation status (irrigated or rainfed) of the main corresponding grid-cell of the Meier 1 km × 1 km map. An irrigation map at a spatial resolution of 300 m × 300 m was obtained, with a single vegetation type attributed to each grid cell together with the irrigation status. The main limitation of this map is that there is no information on the type of irrigation.
In this study, we considered that all irrigation is of "sprinkler" type as this is the most common 120 irrigation type in the USA (AQUASTAT and FAO, 2019), where the testbed area of this study is located. This entails that irrigation water is added to the precipitation forcing over the irrigated agricultural parcels.

Atmospheric forcing
The simulations presented in this study were not coupled with the atmosphere. They were forced by a 125 simulated atmospheric dataset of the European Centre for Medium-Range Weather Forecasts (ECMWF): the ERA-5 atmospheric reanalysis at 0.25° × 0.25° (Hersbach and Dee, 2016). This global dataset was successfully used to force the ISBA LSM in previous studies (e.g. Albergel et al., 2019, Bonan et al., 2020. Beck et al. (2019) showed that the ERA-5 precipitation dataset is reasonably https://doi.org/10.5194/gmd-2021-332 Preprint. Discussion started: 11 October 2021 c Author(s) 2021. CC BY 4.0 License. consistent with gauge-radar data over CONUS, except for mountainous areas. A subset of the ERA-5 130 forcing over Nebraska was used for the time period from 1979 to 2018. The following atmospheric variables were used to force the ISBA LSM and were taken from ERA-5 at an hourly time step: air temperature, wind speed, air specific humidity, atmospheric pressure, shortwave and longwave downwelling radiation and precipitation (liquid and solid).

Validation datasets 135
A large variety of data was used to evaluate the simulations over Nebraska. Six observation datasets were used (Table 1): the water used for irrigation, satellite-derived Leaf Area Index (LAI), gross primary production (GPP), evapotranspiration, land surface temperature (LST), and precipitation.
The water use records were provided by the US Geological Survey (USGS) through the National Water Information System (available at https://waterdata.usgs.gov/ne/nwis/wu, last access August 2021). Every 5 years from 1985 onward, the annual raw amount of water collected for irrigation was available by county, associated to the surface area of the irrigated vegetation. This allowed us to 145 compute the amount of water used for irrigation per unit surface area (in mm) over the specific studied zone in Nebraska (Fig. 1e). The USGS data we used cover the 1985-2019 time period. Because of inconsistencies in the record for 1995, this year was not taken into account.
The simulated LAI was compared with a satellite-derived LAI product at 0.01° × 0.01° spatial resolution derived from SPOT-VGT and PROBA-V satellite data (up to May 2014, andafter May 2014, 150 respectively) by the European Copernicus Global Land Service (CGLS). This LAI product is described in Baret et al. (2013). We used Version 2 of this product (GEOV2). It is available every 10 days for all simulation years. The LAI time series is available from 1999 onward. It does not cover the whole simulation time period (1979 to 2018). The simulated GPP is compared to an upscaled estimate of GPP available at 0.25° from 1980 to 2013, from the FLUXCOM project (Jung et al., 2017). The simulated 155 evapotranspiration is compared to the GLEAM satellite-driven model estimates of land  Martens et al., 2017) at a spatial resolution of 0.25° × 0.25°. Finally, the simulated LST at 12h00 (local solar time) is compared to the LST derived from geostationary meteorological satellites by CGLS at 12h00 (local solar time). This product has a spatial resolution of 0.05° × 0.05° and is available from 2009 to 2018 (Freitas et al., 160 2013).

The ISBA land surface model
The ISBA model (originally described in Noilhan and Planton, 1989) is the LSM developed by the research department of Météo-France (Centre National de Recherches Météorologiques, CNRM). It is embedded into the SURFEX modelling platform (Masson et al., 2013;Voldoire et al., 2017;Le Moigne 165 et al., 2018), and can provide initial land surface conditions to various atmospheric models (e.g. Fischer et al., 2005), or be forced by atmospheric conditions in offline (i.e. stand-alone) mode. In SURFEX, the evolution of land surface states (surface temperature, albedo, roughness…) and fluxes (evaporation, sensible heat, …) is simulated for four different tiles: natural and cultivated lands (e.g. deciduous and broadleaf forests, tropical, temperate and boreal grasslands, crops, deserts, …), 170 urban areas, oceans and inland waters (such as lakes). The ISBA LSM is used to simulate natural and cultivated lands.

ALADIN in
In this study, the version of ISBA including photosynthesis and temporal dynamical LAI evolution in response to environmental conditions was used (ISBA-A-gs; Calvet et al., 1998;Gibelin et al., 2006), together with the multi-layer soil hydrology scheme described in Decharme et al. (2019). 175 Moreover, an updated land cover description was used: ECOCLIMAP-SG (see Supplement S1). The simulations were based on the SURFEX v8.1 version, which is similar to the v8.0 version (Voldoire et al., 2017), but with new technical developments (Le Moigne et al., 2018). Since this study focusses on irrigation, only the tile of natural and cultivated lands was simulated with ISBA, representing the evolution of soil (temperature and water profiles), vegetation (leaf-level and canopy-level 180 photosynthesis, biomass, LAI and carbon fluxes), surface hydrology (runoff and drainage) and snow conditions. To represent the global-scale heterogeneity of continental natural surfaces, twenty different https://doi.org/10.5194/gmd-2021-332 Preprint. surface types (hereafter referred to as "nature types") can be used in ECOCLIMAP-SG to represent the evolution of landscapes with low vegetation, with wooded vegetation, and without vegetation.

Irrigation modelling concept 185
In this study, a pre-existing simple irrigation scheme (Calvet et al., 2008) within the ISBA LSM was upgraded to build a new version able to work at a global scale and to represent several types of irrigation practices. The new irrigation scheme is operated using the ECOCLIMAP-SG land cover classification within SURFEX. The best achievable spatial resolution of ECOCLIMAP-SG is 300 m × 300 m. The irrigation can be activated for ISBA versions able to simulate interactive vegetation biomass 190 and LAI.

Irrigation processes
Within ISBA, irrigation is represented by imposing an additional water flux forcing to the soil-plant system. Water is applied at a given time and over a certain period of time. A number of irrigation variables need to be simulated such as the irrigation amount, the irrigation rate, the irrigation start and 195 end times. A parsimonious approach is used in order to limit the number of parameters of the model. Table 2 lists the parameters and the values used by default in this study. Irrigation is triggered using thresholds of the simulated extractable soil moisture content. Moreover, specific crop phenology parameters such as emergence and harvest dates are used for irrigated crops. Three irrigation types are considered in Lawston et al. (2015): sprinkler irrigation, flood irrigation and drip irrigation. In the new 200 version of ISBA the same irrigation types are represented but a different modelling approach is used. In this study, the sprinkler irrigation type is used. The new irrigation algorithm is based on several steps described below.
First, it is determined whether fields within the grid cell can be irrigated, i.e. they are equipped for irrigation (e.g. water supply, valves, pipes…). This information is given by the irrigation map 205 described in section 2.1.1.
Secondly, it is checked that the vegetation growth stage is compatible with irrigation. For crops, irrigation can be triggered after the emergence and until a few days before the harvest (by default two weeks). In practice, two dates are prescribed: emergence and harvest. This is a simple way to represent specific crop phenology attributes of irrigated crops. Between these two dates, irrigation is possible. 210 Before the emergence and after the harvest, LAI is fixed at the model's minimum value (LAI = 0.3 m 2 m -2 ). This new irrigation scheme is able to support up to three plant growth seasons per year. The crop phenology parameters are not applied to wooded vegetation (trees and shrubs), and can be applied without irrigation. Irrigation can optionally be triggered without considering any specific crop phenology parameter. 215 The availability of resources (equipment or local water distribution) is taken into account trough a default minimum return time period between two irrigations. This default parameter value is a constant (7 days by default) but maps of this parameter could be used when available.
Irrigation can only be triggered when vegetation growths is limited by the extractable soil moisture availability. The plant water stress level is evaluated using a soil wetness index (SWI) along 220 the root profile. The root-zone SWI (SWI root_zone ) is a unitless weighted average SWI value based on the soil volumetric water content profile (Wc i , m 3 m -3 ), the field capacity volumetric water content profile (Wfc i , m 3 m -3 ) and the wilting point profile (Wwilt i , depending on clay and sand fraction, m 3 m -3 ), for each soil layer i. The root fraction inside each soil layer ( ) is used as a weighting factor: where n soil is the number of top soil layers containing roots. This value depends on the considered vegetation type. For example, n soil = 9 for crops, with a rooting depth of 1.5 m.
A SWI root_zone value close to one corresponds to a well-watered soil, while a value close to zero indicates extreme stress. In order to trigger irrigation, the SWI root_zone value is compared to predefined SWI thresholds given as input parameters. These SWI thresholds are evolving during the irrigation 230 season and default values used in this study are fixed to 0.7 for the first irrigation, 0.55 for the second irrigation, 0.4 for the third irrigation, and 0.25 afterwards (following Voirin-Morel, 2003 andCalvet et al., 2008).
If all of these conditions are satisfied, irrigation is triggered with a predefined quantity of water of 30 mm (by default), following Calvet et al. (2008). The irrigation water flux is evenly distributed 235 over a period of time of 8 hours (by default) and is applied on top of the vegetation canopy like precipitation for sprinkler irrigation. In this case, the irrigation water can be intercepted by vegetation and a fraction evaporates. In the case of drip or flood irrigation, the water flux is applied directly to the soil surface, with no leaf interception. In this study, only sprinkler irrigation is considered. Irrigation simulations are illustrated in Supplements S2 and S3 over southwestern France and over the Hampton 240 irrigated area in Nebraska (Fig. 1e), respectively. Observed monthly precipitation in Nebraska is presented for contrasting years in Supplement S4.
All the values of the new parameters presented above have been set within a default configuration. These values can be user-defined for each nature type and for each grid cell, including, when possible, seasonal variations. See Supplement S5 for configuration details and possibilities. 245

Aggregation of irrigated and rainfed vegetation
In contrast to previous versions of ISBA (Voirin-Morel, 2003;Calvet et al., 2008), there is no specific irrigated nature type in the new ECOCLIMAP-SG vegetation description. On the other hand, irrigation of all nature types is possible. The new irrigation scheme is able to represent the sub-grid heterogeneity of the irrigation fractional coverage. For each nature type, an irrigated and a non-irrigated fraction are 250 considered at the simulation resolution. In order to prevent an excessive increase in the number of simulated nature types (potentially 20 non-irrigated and 20 irrigated times 3 irrigation types, i.e. a total of 120 types), involving a large increase of complexity, memory and computing cost, some choices were made for the implementation:

1.
Selection of a limited number of irrigated nature types. The default implementation consists 255 in six irrigated nature types. Temperate deciduous and evergreen trees types (No 8 and 10 in Table S1.2, respectively) can be used to represent fruits trees or olive trees for example, respectively. Shrub type (No 15) can be used to represent, among others, vine plants, and types No 19, 20 and 21 may represent irrigated crops (e.g. wheat, soybean, and corn, respectively). 260

2.
Selection of the main irrigation method used for each grid cell and nature type, considering that in one grid cell there is only one dominant method for a given nature type (e.g. flooded rice in China or sprinkled corn in France). Finally, the system state variables (soil water content, surface and soil temperature, vegetation biomass, etc.) differ in irrigated and non-irrigated parts of the cell. This implies to (1) duplicate a nature 265 type if it is partially irrigated, (2) attribute for each grid cell the corresponding irrigated fraction, and (3) select the irrigation type for the irrigated fraction. Lastly, the two irrigated and non-irrigated nature types are treated separately but the same rooting depth and secondary parameters (see Table S1.1) are used.
In order to limit the computing time, vegetation types can (optionally) be gathered. In this case 270 vegetation "patches" are created. In ISBA, patch aggregation (Masson et al., 2013) is a method used to reduce the number of simulated nature types. It is based on the aggregation of the fractions of nature types, as shown in Fig. 2. The model primary parameters such as rooting depth, LAI and tree height are weighted using the fractional coverage of each nature type in the grid cell. The mean parameter values are calculated following different laws: dominant, arithmetic averaging, inverse averaging or inverse of 275 square logarithm averaging, depending on the considered parameters, as described in Noilhan and Lacarrère (1995) and Noilhan et al. (1997). In practice, the nature types to be aggregated (see the list in  Fig. 2). For secondary parameters (e.g. photosynthesis parameters in this study) a minimum number of patches is needed in order to avoid combining incompatible vegetation types (e.g. C3 crops and C4 crops).
In a first step (step 0 in Fig. 2) the differentiation between irrigated and rainfed nature types is done by computing the irrigated (and rainfed) fraction for each nature type and for each grid cell. 285 Arithmetic averaging is used to cross information from the nature type fractional coverage and from the global irrigation fraction map described in Section 2.1.1. The ECOCLIMAP-SG land cover classification uses this additional data layer to compute the fraction of irrigated vegetation at the spatial resolution of the model simulations. Each nature type considered as irrigated (by default 6) is duplicated (meaning that for each of them a new nature type is created with the same parameters). This ensures the 290 distinction of irrigated and rainfed soil water budget types. Then, as before, the nature types are aggregated by patch and the primary parameter values are computed (step 1 and step 2 in Fig. 2, respectively).
This change of the code structure based on the aggregation tool is a way to (1) maintain the continuity with previous versions of the code, (2) ensure flexibility for the number of irrigated nature 295 types to be considered and (3) simulate distinct irrigated and rainfed fractions of a nature type.

Model implementation and evaluation
The objective of the model evaluation is to demonstrate that the model is able to reproduce irrigation and that the irrigation scheme improves vegetation modelling and the associated surface fluxes as compared to observations. We chose to study an area where the use of irrigation is well documented:  (Table 3): "ISBA_ref" without any irrigation (the benchmark), "ISBA_pheno" with only crop phenology attributes (emergence and harvest dates) and the complete "ISBA_pheno_irr" simulation 315 with irrigation and crop phenology attributes. For the comparison (ISBA_ref vs. ISBA_pheno vs. ISBA_pheno_irr), we selected areas where the irrigation fractional coverage is larger than 50 % as determined from the irrigation map.
In order to assess the consistency of the simulated irrigation process with observations, the simulated number of yearly irrigation events on irrigated areas in Nebraska was compared with the 320 USGS irrigation water amount estimates (Table 1). Since the latter consisted of 5-yearly means, only values of the mean and standard deviation of the yearly irrigation number were compared. Irrigation water amount was converted to a number of irrigation events using the model default irrigation water amount of 30 mm per irrigation event. The comparison was made for the irrigated croplands (either C3 or C4 crops) as defined using the irrigation map (Section 2.1.1) within the studied irrigated area in 325 Nebraska (Fig. 1e).
The reference ISBA_ref LAI simulations were compared with those from ISBA_pheno and ISBA_pheno_irr experiments, and with the 0.01° × 0.01° LAI satellite observations over areas in Nebraska where the vegetation is considered as C3 or C4 irrigated crops by ECOCLIMAP-SG.
In order to compare the time series simulations with observations, the correlation coefficient (R) 330 and the root-mean-square difference (RMSD) scores were used. For water and carbon fluxes, they were calculated using daily values.

Results
The comparison between the model without irrigation (ISBA_ref experiment) before and after the changes in the code structure (section 2.3) did not permit the detection of any impact on the model 335 outputs (not shown). The results presented below are thus focused on the impacts of the crop phenology and irrigation implementation on the simulated land surface variables over Nebraska. In addition to these results, illustrations of the response to irrigation of simulated key land surface variables (SWI, LAI, GPP, evapotranspiration, LST) are shown over southwestern France and over the Hampton area in Nebraska in Supplements S2 and S3, respectively. In the case of Hampton, it can be noticed that the 340 simulated irrigation mainly occurs in July and August (Fig. S3.1).

Irrigation: water use
In Fig. 3, the mean yearly number of irrigation events for C3 and C4 crops for the ISBA_pheno_irr experiment is compared to the values derived from the observations from the USGS. The simulated irrigation numbers present a large interannual variability, with a minimum of 2 in 1993 and a maximum 345 close to 13 in 2002. It must be noticed that 1993 was one of the wettest year recorded at the Lincoln weather station (https://lincolnweather.unl.edu/records/annual.asp, last access August 2021) and that year 2012 was characterized by a prolonged drought. The mean simulated value of the yearly irrigation water amount used for irrigation (271±75 mm year -1 ) is almost identical to the observed one (264±65 mm year -1 ), with a difference of only +2.7%. This difference is small, although the model does not take 350 into account the availability of the water resource yet. Figure 4 illustrates the mean seasonal and interannual variability of LAI in the most densely irrigated part of Nebraska for areas with a fraction of irrigated crops larger than 50 % in Fig. 1e, from 1999 to 2018. Table 4 presents the peak LAI characteristics. In all ISBA LAI simulations, the start of the 355 growing season corresponds to a gradual increase in LAI from the initial value of LAI = 0.30 m 2 m -2 imposed to the model in winter. The observed LAI presents a smaller minimum LAI value of 0.15 m 2 m -2 , starts increasing in April and a value of 0.30 m 2 m -2 is reached at the end of April. Then, plant growth continues at about the same low rate till the end of May. The LAI growth rate increases in June and LAI reaches a mean peak value of 4.9 (±0.8) m 2 m -2 is observed on 31 July (Table 4). The observed 360 LAI then sharply decreases to reach its minimum value at about the end of September.

Irrigation: plant growth
The ISBA_ref LAI simulations do not mirror the observed late growing season and rapid senescence. The ISBA_ref vs. observations comparison shows that without irrigation the simulated LAI generally starts increasing in March. On average, a peak LAI value of 3.6 (±0.2) m 2 m -2 is simulated by ISBA_ref on 2 July, before slowly decreasing until the end of December. The ISBA_ref growing season 365 is much longer than observed. It starts two months before the observations and stops three months after the observations. The simulated LAI peaks one month before the observations. The simulated yearly LAI amplitude is 28 % smaller than observed.
The ISBA_pheno LAI simulation is much more consistent with the LAI observations. The growing season starts in mid-May and the senescence ends at the end of September. However, the 370 simulated peak LAI is still 30 % smaller than observed (LAI = 3.5 (±0.2) m 2 m -2 ). The peak LAI is reached on 26 August, much later that the ISBA_ref peak LAI, and about one month after the observed peak. The sharp decrease of LAI in September results from harvests at random dates in September.
Adding irrigation (ISBA_pheno_irr) does not change the general pattern of the LAI curve, but increases the LAI amplitude, with a mean peak LAI value of 3.7 (±0.1) m 2 m -2 on 28 August, larger (+8%) than 375 for ISBA_pheno but still below the observation (-24%). Fig. 4b, from 2002 to 2008. The ISBA_ref LAI presents a systematic drop in summer, which is not present in the observations nor simulated by the ISBA_pheno and ISBA_pheno_irr experiments. Without the regular seasonality imposed by crop phenology parameters, the model may simulate a re-growth of vegetation 380 in autumn (e.g. in 2003), that is not present in the observations. The ISBA_pheno and ISBA_pheno_irr simulations are more consistent with the observed seasonality.

Impact of crop phenology and irrigation on LAI at a regional scale
This section is focused on the impact of irrigation practices for the south Nebraska zone (as defined in

Impact on the GPP flux
As the vegetation productivity is linked to LAI, the seasonality pattern of GPP (Fig. 7) is comparable to the one of LAI (Fig. 5) but the observed GPP peak (9.2 ± 2.1 g [C].m -2 .day -1 ) occurs on mid-July while 405 the observed LAI peaks on 31 July. During the plant growing period, the smallest differences between all the simulations and the observations occur at about the same time as the observed GPP peak. For all simulations, a GPP plateau at a value of 9.0 ±1.8 g[C] m -2 day -1 is reached at the beginning of July and lasts until mid-July. Finally the simulated GPP decreases in September with a delay of about two weeks with respect to the observations. 410 Before July, the ISBA_ref photosynthetic activity is well in advance as compared to the observations, of about 20 days in May. This is consistent with the very large LAI values simulated by ISBA_ref in May: about 2 m 2 m -2 , while the mean LAI observation hardly exceeds 0.5 m 2 m -2 . The simulated GPP maximum (9.7 ± 2.0 g[C] m -2 day -1 ) is reached before the end of June. After a sharp decrease at the end of June, the ISBA_ref GPP decreases at a slower rate than the observations. From 415 mid-September to the end of October, the simulated GPP is much larger than the observed GPP.
The ISBA_ref flaws are much less pronounced in ISBA_pheno and ISBA_pheno_irr experiments. In the latter simulations, the increase of the GPP occurs at about the same time as in the observations. The GPP values are systematically larger with irrigation in July and August than for other https://doi.org/10.5194/gmd-2021-332 Preprint. Discussion started: 11 October 2021 c Author(s) 2021. CC BY 4.0 License. simulations. As for LAI, the GPP R and RMSD scores ( Fig. 7b and 7c, respectively) are better for 420 ISBA_pheno_irr than for ISBA_ref, with an improvement on 87 % and 81 % over the domain, respectively.

Impact on evapotranspiration
Investigating evapotranspiration is a way to assess the impact of irrigation on the hydrological system. On the contrary, from mid-July to mid-August, all ISBA simulations tend to underestimate evapotranspiration with respect to GLEAM, by up to 1.3 mm day -1 for ISBA_ref. Accounting for crop 435 phenology and irrigation into the model has a substantial impact on this variable and reduces the bias.
Over the whole irrigation period, the mean bias goes from -0.4 ± 0.4 mm day -1 (-13 ± 12 %) for ISBA_ref to -0.2 ± 0.3 mm day -1 (-7 ± 11 %) and -0.1 ± 0.3 mm day -1 (-2 ± 11 %) for ISBA_pheno and ISBA_pheno_irr, repectively. Evapotranspiration is overestimated again after the harvest, from mid-September to November by 0.38±0.18 mm day -1 (42±20 % compared to the observations). 440 The newly implemented processes have a small but positive impact on the bias before and after the irrigation period. During the growing season, from April to June, the overestimation decreases from 38 % in ISBA_ref to 33 % and 34 % for ISBA_pheno and ISBA_pheno_irr, respectively. From mid-September to November, the overestimation decreases from 42 % to 35 % and 36 %, respectively. The R and RMSD differences between ISBA_ref and ISBA_pheno_irr (Fig. 8)  is small (less than 0.1) and heterogeneous in time and space. Figure 9 shows that the R score is mainly improved in August and in September, before the harvest. The improvement of RMSD is more stable, and can be observed from May to October, the impact being more pronounced in July and August.

Impact on LST 450
In order to evaluate the impact of irrigation on the land surface energy budget, Figure 10 shows land surface temperature at 12h00 local time simulated by the three model configurations and derived from the CGLS product. Overall, ISBA tends to overestimate LST at noon, especially in April-May, up to 7 °C in Fig. 10a. The bias is reduced during the summer.
Due to the difficulty of observing the differences between the simulations, Figure 10b presents 455 differences of ISBA_pheno and ISBA_pheno_irr versus ISBA_ref. With crop phenology (with or without irrigation) the simulated LST is globally higher from April to June and from mid-September to November. The maximum difference with respect to ISBA_ref is +0.7±0.3 °C. It is observed for ISBA_pheno in September. During the summer (July and August) the new model versions tend to present lower LST values, with temperature differences close to -0.2 ± 0.1 °C in ISBA_pheno_irr. 460 Moreover, from May to mid-September the temperature in ISBA_pheno_irr is lower than in ISBA_pheno, and this difference can reach locally -0.9°C in summer. Figure 11 presents the monthly R and RMSD scores of ISBA_ref with respect to CGLS LST observations and the ISBA_pheno_irr score difference is shown. It shows that there is a seasonal dependence of these statistical values, with slightly better R and RMSD scores observed for ISBA_pheno_irr in July and August during the irrigation 465 period. However, the representation of irrigation tends to degrade RMSD before (April, May) and after (October, November) the irrigation period.

Discussion and perspectives
The results presented in Section 3 show that the new version of ISBA is able to produce a realistic yearly irrigation water amount (Fig. 3). It also markedly improves the LAI and GPP simulations (Fig. 4-470 limited impact on the evapotranspiration and on the LST simulations and is not able to significantly reduce the strong model biases that are observed for these variables before and after the irrigation time period (Fig. 8-10).

Could the new irrigation scheme be further improved? 475
The results of our numerical experiments over Nebraska show that considering crop phenology and irrigation improves the consistency of the simulations with LAI and GPP observations. The corresponding correlation and RMSD scores are improved. Two new developments can explain this behaviour: (1) the crop phenology parameters used to force emergence and harvest dates reduce the length of the growing season, delay spring growth and avoid a regrowth in the autumn, and (2) the 480 irrigation limits the water stress and enhances plant growth at summertime. Nevertheless they both have shortcomings and their performance could be limited by difficulties in simulating processes that are not directly related to irrigation.
Firstly, the same emergence and harvest dates are imposed for all years, while in reality crop phenology may present an inter-annual variability related to climate conditions. This is particularly the 485 case for Nebraska because the start of the growing season depends on the snow melt and soil thawing dates. These processes are represented in ISBA and crop phenology parameters could be related to snow melting and soil thawing, but this would require extensive developments to be implemented at a global scale. Moreover, the representation of the cold season processes is not perfect in ISBA (Decharme et al. 2019) and this could explain biases in soil temperature and LST simulations before and after the 490 irrigation time period. Figure 10 shows that LST values below the freezing level can be observed in April and that their model counterparts are about 7 °C warmer. The earlier thawing in model simulations is reflected in the much earlier leaf onset in LAI simulations. Figure 5 shows that while the observed LAI does not exceed 0.5 m 2 m -2 at the end of April, the ISBA_ref LAI reaches the same value about one month earlier. The unrealistically early leaf onset is consistent with the warm model bias at 495 the end of the cold season. This shows that improving the representation of the cold season by assimilating satellite-derived or in situ snow cover fraction observations could improve the simulation of the crop growing period in this area. However, the currently used empirical approach to establish the crop season provides robust results over the irrigated grid cells (Fig. 4).
Secondly, the irrigation itself is based on fixed parameter values such as the minimum period 500 between two consecutive irrigations (one week) and SWI levels triggering irrigation turns. The simulations over the Hampton grid cell show that the first irrigation can start at quite low levels of the SWI (Fig. S3.1), even below the second irrigation threshold of 0.55 defined in Section 2.3.1.
Suppressing the one week constraint of irrigation turns improves the simulation of the peak LAI, which otherwise is rather poorly simulated (Fig. 4). However, this change triggers unrealistic large irrigation 505 water amounts (not shown). A lack of irrigation water amount cannot explain the excessive soil water deficit. One could also challenge the quality of the ERA5 precipitation. Figure S4.1 and Fig. S4.2 show that ERA5 precipitation compares well with in situ observations and that the seasonal and inter-annual variability is fairly represented. A more plausible explanation could be that the initial soil water storage value between the end of the cold season and the first irrigation turn is withdrawn too quickly from the 510 soil by the model. This explanation would be consistent with the marked overestimation of evapotranspiration in spring, from April to June (Fig. 8), before the irrigation time period.
In order to investigate the evapotranspiration bias in spring, the evaporation components were plotted in Fig. S3.3 and Fig. S3.4 for the Hampton irrigated area in 2018. Figure S3.3 shows that total evapotranspiration of ISBA_ref and ISBA_pheno_irr are quite similar. This is consistent with the small 515 impact of crop phenology on total evapotranspiration showed in Fig. 8. On the other hand, soil evaporation and plant transpiration differ. In the ISBA_pheno_irr simulation, transpiration is reduced in spring by 30 % to more than 100 %, in comparison with ISBA_ref. The lower transpiration is compensated by larger soil evaporation values. As a result, total evapotranspiration does not change much and the bias is not reduced in ISBA_pheno_irr. Also, Fig. S3.4 shows that the new irrigation 520 module does not affect interception much. Therefore, the ISBA_pheno_irr evapotranspiration bias in spring could be caused by the large soil evaporation. The evaporation component could be overestimated because (1) the soil is too warm in relation to a poor representation of thawing or because (2) crop residues at the soil surface are not represented. Wortmann et al. (2012) showed that in this area, not harvesting crop residues tends to reduce soil evaporation and increase crop yield, limit water runoff, 525 soil erosion, and contributes to maintaining soil fertility. The ISBA model includes a representation of litter in forests (Napoly et al. 2017) that will be generalized to low vegetation in the next version of SURFEX. Using this new capability could improve our simulations.

Is the irrigation scheme flexible enough?
In this study, sprinkling irrigation was considered. The model is able to represent other irrigation 530 systems such as flooding irrigation but more developments are needed to limit the runoff to the irrigated plot and this options needs to be validated. The newly implemented irrigation processes, along with the new ECOCLIMAP-SG vegetation description let users choose which nature type should be irrigated.
Irrigation can be represented at various spatial scales, ranging from the field scale for agricultural studies to the global scale for climate studies. Model parameters can be specified using new datasets or 535 local characteristics. For example, in this article we used a unique date for starting and ending the crop growing season with a random variability, but more accurate dates could be prescribed (varying spatially and from one vegetation type to another, or using crop calendars). Moreover, the better spatial resolution of ECOCLIMAP-SG allows the use of high resolution atmospheric forcing. This provides new opportunities for assessing the impact of irrigation on local climate and water resource conditions. 540 This study is mainly focused on a zone in the south of Nebraska where the irrigation density is relatively high (Fig 3), and results could differ in other regions. Except for the fixed emergence and harvesting dates corresponding to regional crop phenology (from USDA and NASS, 2010), default values are used for all the other parameters (Section 2.4). Tests performed in southwestern France (Supplement S2) allowed ensuring that the model is able to work in contrasting climate conditions. 545 In this study, the ISBA simulations were not coupled to the atmosphere, nor to the CTRIP river routing system. Such coupled numerical experiments could be performed thanks to the SURFEX modelling platform. However, more developments are needed in order to ensure water conservation in the hydrological system. In particular, irrigation water amounts should be consistent with the available water resource in rivers, groundwater, and dams. 550 A new uncalibrated irrigation scheme was implemented within the ISBA land surface model in order to improve the representation of vegetation over agricultural areas. A case study over an irrigated area in the state of Nebraska (USA) was performed to validate the new scheme. Simple crop phenology rules represent emergence and harvesting and improve the seasonality of plant growth, while the additional 555 water supply from the irrigation mostly impacts the peak LAI value. The model is able to produce a realistic yearly irrigation water amount and markedly improves the LAI and GPP. It is shown that model performance can be limited by processes not directly related to irrigation, such as thawing or crop residues. The irrigation scheme has many possible configurations and the code is highly flexible.
With this capability, ancillary data on farming practices such as emergence and harvest dates, or the 560 amount of water per irrigation event, could be used.

Code availability
The ISBA land surface model is available as open source via the SURFEX modelling platform, available at https://www.umr-cnrm.fr/surfex/spip.php?article387. It is under a CECILL-C License 565 (French equivalent to the L-GPL licence). The version developed and use for the experiment in this study is available in the git-hub: https://github.com/ArseneD/OPEN_SURFEX_V81_IRR. It based on the SURFEX version 8.1 (ref f70f6457). For future use, it is strongly recommended to use the newest version of ISBA, from the version 9.0 (scheduled for release in 2020) from which the irrigation developed will be included by default.     . 4) and all vegetation types (see Fig. 5).