Impacts of climate change and rising atmospheric CO2 on future projected reference evapotranspiration in Emilia-Romagna (Italy)

The continuous increase of atmospheric CO2 content mainly due to anthropogenic CO2 emissions is causing a rise in temperature on earth, altering the hydrological and meteorological processes and affecting crop physiology. Evapotranspiration is an important component of the hydrological cycle. Thus, understanding the change in evapotranspiration due to global warming is essential for better water resources planning and management and agricultural production. In this study, the effect of climate change with a focus on the combined effect of temperature and elevated CO2 concentrations on reference evapotranspiration (ETo) was evaluated using the Penman–Monteith equation. A EURO-CORDEX regional climate model (RCM) ensemble was used to estimate ETo in five locations in the Emilia-Romagna region (Northern Italy) during the period 2021–2050. Then, its projected changes in response to different CO2 concentrations (i.e., 372 ppm and 550 ppm) under two Representative Concentration Pathways (RCP) scenarios (i.e., RCP4.5 and RCP8.5) were analyzed. Simulation results with both scenarios, without increasing CO2 levels (372 ppm), showed that the annual and summertime ETo for all locations increased by an average of 4 to 5.4% with regard to the reference period 1981–2005, for an increase of air temperature by 1 to 1.5 °C. When the effect of elevated CO2 levels (550 ppm) was also considered in combination with projected changes in temperature, changes in both annual and summer ETo demand for all locations varied from − 1.1 to 2.2% during the 2021–2050 period with regard to the reference period 1981–2005. This shows that higher CO2 levels moderated the increase in ETo that accompanies an increase in air temperature.


Introduction
The increase of anthropogenic emissions of carbon dioxide (CO 2 ) and other radiatively active gases, due to activities such as fossil fuel combustion and deforestation, is one of the most prominent causes of climate change (Hamza et al. 2020). The relentless rise of CO 2 concentrations will most probably influence water fluxes and resources at both local and global scales (Haddeland et al. 2014). This effect will be intensified by the projected increases in global surface air temperature (IPCC 2013) leading to an uneven distribution of water resources and causing several problems in water availability (Estrela et al. 2012;Boehlert et al. 2016;Akbarpour and Niksokhan 2018). This is especially true in the Mediterranean region as it has shown large climate shifts in the past and it has been identified as one of the most prominent "hot spots" in future climate change projections, with an expected greater than average warming, mostly in summer, and an increase of heat waves and dry spells (Olesen and Bindi 2002;Giorgi and Lionello 2008). Similar conclusions were drawn on a smaller scale in Northern Italy and Emilia Romagna region. The future projections constructed at a local scale through statistical downscaling techniques in the framework of different emission scenarios indicate that significant changes in temperature are very likely to occur during all seasons. A shift of the probability distribution functions to warmer values is projected for both minimum and maximum air temperature. The shift will be more 1 3 noticeable in summer, when the changes of mean values are expected to be about 2.5 °C in the 2021-2050 period (with respect to 1961-1990) (Tomozeiu et al. 2007(Tomozeiu et al. , 2014. These changes are projected to become more pronounced toward the end of the century. To better understand the effect of climate change on the hydrological cycle, evapotranspiration is widely used as it is the only link between the energy balance and water balance (Dong et al. 2020). Evapotranspiration represents the simultaneous processes of transfer of water to the atmosphere by transpiration and evaporation in a soil-plant system (Allen et al. 1998). In particular, reference evapotranspiration (ETo) is an important hydrological variable for irrigation water management and hydrological modeling (Immerzeel and Droogers 2008;Pereira et al. 2020). Since the agricultural sector accounts for 70% of freshwater withdrawals (Kummu et al. 2012) and around 90% of global freshwater consumption (Wada et al. 2012), improved understanding of changes in evapotranspiration is essential for dependable projections of global freshwater availability under conditions of climate change. An increase in air temperature is known to cause an increase in crop water requirements because it tends to stimulate plant transpiration (Allen et al. 2003), hence decreasing water use efficiency (WUE). This is mainly due to a decrease in leaf photosynthesis and an increase in stomatal conductance to water vapor (Kirschbaum 2004). Additionally, higher temperatures, especially when combined with lower precipitation, would raise the irrigation requirements of crops (Döll 2002;Moreno et al. 2005), alter the length of the growing season (Menzel and Fabian 1999;Linderholm 2006), and negatively affect crop yield (Lobell and Field 2007). On the other hand, the increase in CO 2 concentration has opposite physiological effects on plants, e.g. faster photosynthesis rate, greater leaf area, increase in biomass and yield, and decrease in stomatal conductance and evapotranspiration rate (Allen 1990;Ainsworth and Long 2004;van der Kooi et al. 2016), which could offset the consequences of air temperature increase. However, the magnitude of this effect is uncertain and still controversial within the scientific community (Allen 2000;Long et al. 2005;Tubiello et al. 2007) because evapotranspiration is a complex and nonlinear phenomenon depending on several other interacting factors such as humidity, wind speed, radiation, and type and growth stage of crops (Polley 2002). The impact of the interaction of these different factors on evapotranspiration was assessed using a modification of the standardized (FAO-PM) Penman-Monteith equation (Allen et al. 1998) by Snyder et al. (2011) and Moratiel et al. (2011). They adjusted the stomatal resistance factor (the reciprocal of stomatal conductance) to account for an increase in CO 2 concentration (Lovelli et al. 2010). Other approaches could also be used in association with the FAO-PM equation (Ben Hamouda and Ventura 2020). However, they are either more complex or insufficiently tested. The main issue with the Penman-Monteith equation is data availability and the difficulty to access to these data, especially in Italy where weather monitoring networks are managed by regional and national services, without a common data sharing policy and data distribution platform (Pavan et al. 2013;Pelosi et al. 2020). This prompted researchers to look for alternative ways for obtaining the input data: estimation (de Carvalho et al. 2013;Córdova et al. 2015), weather forecast (Cai et al. 2007;Silva et al. 2010;Lorite et al. 2015), satellite imagery (Montero et al. 2018), remote sensing (Teixeira 2010), a combination of remote sensing and machine learning , artificial neural networks (ANNs) (Kumar et al. 2002), fuzzy and neuro-fuzzy systems, and genetic algorithms (Shiri et al. 2012;Kisi and Cengiz 2013). Another alternative choice is using reanalysis data. Climate reanalysis combines past observations with models to generate consistent time series of multiple climate variables, providing comprehensive snapshots of conditions at regular intervals over long time periods often years or decades (Parker 2016). The main advantages of reanalysis products are their spatial and temporal resolution consistency over three or more decades, the large quantity of variables available, and the continuous improvements of their model resolution and biases (Schubert et al. 2008;Dee et al. 2016). This is why reanalysis products are among the most-used datasets in the study of weather and climate (Fuka et al. 2014;Essou et al. 2017;Bhattacharya et al. 2019;Uniyal et al. 2019;Luo et al. 2020). However, few studies compared ETo computed using reanalysis data with observed ETo data (Yao et al. 2014;Martins et al. 2017;Tian et al. 2018).
The aim of this study is to assess the effect of future climate, with a focus on evapotranspiration (ETo), in Emilia-Romagna (Northern Italy) over the 2021-2050 period and under different emission scenarios. We use two Representative Concentration Pathway scenarios (RCPs) and three Regional Climate Models (RCMs). ETo is first evaluated using reanalysis data as an alternative to observations. Then, the bias correction of the different RCMs is performed. Finally, changes in ETo under different RCPs with or without the effects of elevated CO 2 levels are evaluated.

Study area
Emilia-Romagna is an administrative region located in Northern Italy, between 43° 80′ and 45° 10′ N and 9° 20′ and 12° 75′ E. Agriculture is an important sector in the region: within the total area of 22 452 km 2 , 10 642 km 2 are devoted to agriculture (Regione Emilia-Romagna 2010). Climate in the region is mainly influenced by a variable geomorphology, represented by the Po River in the North, the Adriatic Sea in the East, and the Apennines Mountains to the South. The Köppen-Geiger climate classification indicates a temperate fully humid climate with hot summers (Cfa) in the North, central, and northeastern plain areas, and a fully humid climate with warm summers (Cfb) in the mountain and highlands areas, especially in the South and South-West (Kottek et al. 2006). The locations of the five agrometeorological stations selected for this study and their details are shown in Fig. 1 and Table 1. The choice of these stations is based mainly on their data availability and their representativeness of the regional territory.

Reanalysis data
Since it is difficult to find a continuous observed dataset long enough to make climatological estimates of ETo, it was decided to evaluate it starting from hourly ERA5-Land reanalysis data downloaded from the Copernicus Climate Data Store (Copernicus CDS) (https:// cds. clima te. coper nicus. eu). ERA5-Land is a reanalysis dataset developed by the Copernicus Climate Change Service (C3S) at the European Centre for Medium-Range Weather Forecasts (ECMWF), providing a consistent view of the evolution of land variables over several decades at 9 km resolution, which is an enhanced resolution compared to ERA5 (31 km). ERA5-Land has been produced by replaying the land component of the ECMWF ERA5 climate reanalysis. The core of ERA5-Land is the Tiled ECMWF Scheme for Surface Exchanges over Land incorporating land surface hydrology (H-TESSEL). Reanalysis data and specifically ERA5 and ERA5-Land give satisfactory results when compared with observations (Ben Hamouda et al. 2019;Mahto and Mishra 2019;Tarek et al. 2020;Pelosi et al. 2020).
ERA5-Land data are produced with an hourly time-step from 1981 to the present. The 25 years between 1981 and 2005 will be used as a reference period. Extracted data for that period include mean 2 m air temperature (K), mean 2 m dew point temperature (K), mean downward surface  solar radiation (J m −2 ), and mean 10 m u and v components of wind (m s −1 ). The retrieved 10 m wind speed components were adjusted for a 2 m height using a logarithmic wind speed profile (Allen et al. 1998), then wind speed was calculated as WS = √ u 2 + v 2 . Maximum and minimum air temperature were transformed to °C. All data were aggregated to a daily timescale.

Observed data
Observed data for the different meteorological variables needed for the calculation of ETo for the reference period 1981-2005 were not available or incomplete. Since data were only available for the period 2008-2018 for the five agrometeorological stations, they were hence used only for validation purposes. Data for the stations of Martorano, San Pancrazio, San Pietro Capofiume, and Volano were made available by the "Regional Agency for Prevention, Environment and Energy of Emilia-Romagna" (ARPAE) via their Dext3r application (https:// simc. arpae. it/ dext3r). Data from ARPAE are checked for plausibility and quality on a daily basis. Quality checks include a spatial comparison between all available data, consistency between all meteorological variables at an hourly frequency, including also precipitation data, which are not considered in the present work. During the quality control process, all possible sources of information are taken into account, including geographical characteristics of stations' locations, such as the distance from the coast, altitude, and the relative altitude with respect to surrounding locations, and synchronous radar data. These operational activities are executed at ARPAE since 2001 and cover all data used in the present study. The time series of the chosen stations are characterized by good data consistency, covering more than 80% of the validation period, and each station's location is representative of a relatively large surrounding area. No data filling procedures have been used to fill in the missing data. In all stations, all parameters relevant for the present study have been measured during the period considered.
Data include daily mean 2 m air temperature (°C), daily mean relative humidity (%), hourly mean global solar irradiance (W m −2 ), and daily mean 10 m wind speed (m s −1 ). The August-Roche-Magnus (Lawrence 2005) approximation was used in order to compute the dew point temperature from the mean relative humidity and air temperature. The retrieved 10 m wind speed was adjusted for 2 m height. The same meteorological variables were also available at an hourly basis from the agrometeorological station of the Department of Agricultural and Food Sciences of the University of Bologna (Distal), installed at the experimental farm of Cadriano (Bologna), whose data were also quality-controlled, with uncertain data compared to data from other instruments in the same station or from other stations close to it (Matzneller et al. 2009). The only difference with ARPAE data is that mean wind speed data were measured at 2 m height. Finally, all hourly data were aggregated to a daily timescale.

EURO-CORDEX climate projections
To study the future changes in ETo, we used data made available by EURO-CORDEX Regional Climate Models (RCMs) with a daily time step. In particular, data from three regional models were used, RCA4 (Christensen et al. 2007), RACMO22E (Van Meijgaard et al. 2008), and HIRHAM5 (Samuelsson et al. 2015), with a horizontal resolution of 0.11° (12 × 12 km). EURO-CORDEX is the European branch of the international CORDEX initiative, which is a program sponsored by the World Climate Research Program (WRCP). It is the largest ensemble of RCM simulations covering all of Europe (Jacob et al. 2014). The choice of the RCMs was mainly based on the availability of the data necessary for the calculation of ETo. Daily data were extracted from Copernicus CDS for the entire available period 1981-2005, considered as the reference period, and the future period 2021-2050, and for two Representative Concentration Pathways RCP4.5 and RCP8.5. Representative Concentration Pathways (RCPs) were defined by the Intergovernmental Panel on Climate Change (IPCC) for its 5 th Assessment Report (AR5) to provide time-dependent projections of atmospheric greenhouse gas (GHG) concentrations based on socio-economic scenarios of how global society grows and develops (van Vuuren et al. 2011). RCP4.5 was developed using the MiniCAM model ) and represents a future where climate policies limit and achieve stabilization of greenhouse gas radiative forcing to 4.5 W m −2 by 2100, while RCP8.5 was developed using the MESSAGE model  and is called a business-as-usual scenario, where high emissions of greenhouse gases continue as no policy changes are taken to reduce emissions. Data include 2 m mean air temperature (K), 2 m mean relative humidity (%), surface solar downward irradiance (W m −2 ), and 10 m wind speed (m s −1 ). Dew point temperature was calculated using the August-Roche-Magnus approximation and the 10 m wind speed was adjusted for 2 m height.

The Penman-Monteith model
Reference evapotranspiration (ETo) is the ET estimate from a standardized hypothetical crop having fixed height, albedo, and canopy resistance (Allen et al. 1998). In this study, daily ETo (mm day −1 ) was computed using the FAO-Penman Monteith equation (FAO-PM) where R n net radiation at the surface (MJ m −2 day −1 ). G soil heat flux density (MJ m 2 day −1 ). T 2m mean daily air temperature at 2 m height (°C). u 2 wind speed at 2 m height (m s −1 ); e s saturation vapor pressure (kPa); e a actual vapor pressure (kPa). (e s − e a ) saturation vapor pressure deficit (kPa). ∆ slope of the vapor-pressure curve (kPa °C −1 ). γ psychometric constant (kPa °C −1 ).
To modify the CO 2 concentration in the FAO-PM equation, it is possible to estimate a new canopy resistance r c value , which is represented in the equation above by the 0.34 in the denominator: where r a is the aerodynamic resistance to sensible and latent heat transfer (s m −1 ). The stomatal resistance r s of an actively transpiring C 3 grass leaf surface has a value of about r s = 100 s m −1 . Using a leaf area index (LAI = 24 × crop height (h) = 24 × 0.12 = 2.88) and, assuming only half of the LAI contributes to transpiration, the canopy resistance for 0.12 m tall C 3 species grass r c is calculated as: The canopy resistance value of 70 sm −1 in the FAO-PM equation is assumed to apply to a CO 2 concentration of 372 ppm. Thus, to estimate possible impacts of higher CO 2 concentrations on ETo, r c is modified using the following form of the FAO-PM equation: According to Long et al. (2004), the stomatal conductance of many C 3 plants grown in Free-Air CO 2 Enrichment (FACE) experiments decreased about 20% when the CO 2 concentration was increased from 372 to about 550 ppm for about 200 independent measurements. Assuming this is true for the stomatal conductance of 0.12 m tall C 3 species grass with a stomatal resistance of 100 s m −1 , the stomatal conductance for C 3 grass should decrease from about 10 to 8 mm s −1 , which corresponds to r s = 125 s m −1 . Using T 2m +273 u 2 e s − e a + γ 1 + r c ∕r a the same approach used to calculate r c in Eq. (3), the r c for 550 ppm is calculated as: Thus, increasing CO 2 concentration from 372 to 550 ppm should increase canopy resistance of 0.12 m tall C 3 grass from 70 to 87 s m −1 .
Reference evapotranspiration was calculated using daily data and then aggregated to annual and summer (June to August) values since the latter is the most important season for agriculture and irrigation purposes in the Northern Hemisphere.

Performance indicators
The performance of the reanalysis and RCMs data was evaluated at each agrometeorological station using the root mean square error (RMSE) and bias as statistical indicators, which are calculated as follows: where A is any of the model simulated data, B is any of the reference data, i is the day number, and n the number of examined days. A small RMSE indicates that the difference between the simulated and reference data is small, while the closer the bias is to zero, the better the simulation will be. A positive bias means that there is an overestimation, and a negative bias means that there is an underestimation.

RCM bias-correction
RCM simulations are subject to possible biases that are mostly due to parameterization of sub-grid processes, limited representation of local features, incorrect boundary conditions, and differences between spatial resolutions of the simulations and observations (Benestad 2010;Ehret et al. 2012). To minimize these biases, the quantile mapping (QM) technique was employed to bias-correct air temperature, wind speed, solar radiation, and dew point temperature data extracted from the three EURO-CORDEX RCMs, for both the historical and future period, using ERA5-Land data as reference. ETo calculated using raw weather variables was also bias-corrected as this direct bias-correction gave slightly lower RMSE and bias than ETo calculated using the bias-corrected weather variables. The QM approach has become increasingly popular and widely used because of its efficiency and low computational cost (Fang et al. 2015;Sun et al. 2019;Pastén-Zapata et al. 2020;Torma et al. 2020). For this purpose, the full calibration period 1981-2005 is considered. This means that evaluation and validation periods are identical, and no independent crossvalidation exercise is carried out. The latter was covered by previous works (e.g., Ivanov and Kotlarski 2017;Gutiérrez et al. 2019) which revealed that QM, in general, performs well in historical cross-validation setting with independent calibration and validation periods (Feigenwinter et al. 2018).

Taylor diagram
In addition to computing the mean bias and RMSE, the degree of correspondence between the raw and bias-corrected RCM simulations and reanalysis outputs can be quantified and illustrated schematically in the form of a Taylor diagram (Taylor 2001). These diagrams are especially useful in evaluating multiple aspects of complex models or in gauging the relative skill of many different models (Yin et al. 2016;Tegegne and Melesse 2020;Torma et al. 2020). The performance of the different models is evaluated by comparing their correlation (R), centered root-mean-square difference (cRMSD), and the amplitude of their variations represented by their standard deviations (StDev o and StDev s ): where X o is the observed value, X s is the simulated value, Xo is the average of the observed value, Xs is the average of the simulated value, and N is the number of values. For each model, the three aforementioned statistics are plotted: (1) the Pearson correlation coefficient (R), which is related to the azimuthal angle; (2) the centered RMS difference (cRMSD) in the simulated field, which is proportional to the distance from the point on the x-axis (black contours) that is identified as "observed" (represented by a black circle); and (3) the ratio of standard deviation derived from the RCM simulations against the observations, which is proportional to the radial distance from the origin (red contour). The colored symbols represent the different RCM simulations where the color refers to the RCM and the different symbols represent the raw and bias-corrected daily ETo means. Simulations that agree well with observations will lie nearest the reference point on the x-axis. Finally, the different steps followed by the methodology in this study are summarized in the flow chart of Fig. 2.

Validation of the ERA5-Land reanalysis with observed data
The performance of the reanalysis and the observed data of four climatic variables, i.e., mean air temperature (T mean ), solar radiation (R s ), mean dew point temperature (T dew ), and mean wind speed (WS), along with ETo, calculated using these variables in the validation period 2008-2018 (4018 days), were evaluated at each agrometeorological station using the root mean square error RMSE and the bias as statistical indicators. Results are shown in Table 2. For all stations, T mean , T dew , and WS RMSE and bias are less than 20% of the observed values and the coefficient of determination R 2 is higher than 0.5. For R s , RMSE and bias are slightly higher but less than 30% with a minimum R 2 of 0.8. Looking at ETo performance in Table 2, it is clear that ERA5-Land slightly underestimates ETo in most locations, with bias values ranging between − 0.8 and 0.03 mm day −1 , a minimum RMSE of 0.2 mm day −1 and a maximum RMSE of 1.2 mm day −1 . This could be due mainly to the effect of the solar radiation underestimation especially because it was demonstrated that solar radiation errors have a high impact on the estimated evapotranspiration (Perera et al. 2014;Pelosi et al. 2016).
For all stations, both RMSE and bias are lower than 30% of the ETo from observations, while the minimum R 2 is around 0.8. Martins et al. (2017) considered an ETo RMSE lower than 0.8 mm/day satisfactory when they compared ETo from reanalysis with observed ETo in the Iberian Peninsula. However, the variability of the performance of ERA5-Land data could be explained by the presence of missing data in the observed time series in all stations, especially for Martorano and Volano, that could negatively affect RMSE and bias values (Kidson and Trenberth 1988). In the latter location, it could be further explained by the fact that ERA5-land is produced by running the scheme of surface exchanges over land without taking into account the possible vicinity of the sea. Observational practice indicates that surface parameters in coastal areas may be influenced by the sea up to 5-10 km inland; as a consequence, it is possible that ERA5-land may not correctly represent local climate variability in this area (Pelosi et al. 2020). The extension of the analysis to the regional level could help to understand if that effect is relevant at that level or if it is mostly a very local issue. In general, it is known that reanalysis datasets may feature a certain degree of deviations from the "true" value (Maraun et al. 2010). Overall, these results are considered acceptable and since, in our case, observed data necessary for the calculation of the FAO-PM equation for the study period 1981-2005 are not available or have large gaps, it was decided to use ERA5-Land data as a surrogate of ground observations and for the bias-correction of the different RCMs.

Validation and bias correction of RCMs
After this first step, reanalysis data were used to evaluate data extracted from the three RCM models for our five reference stations. As stated before, those RCMs were selected Fig. 2 Flowchart of the proposed methodology used to assess future reference evapotranspiration, using the FAO-PM equation and historical/future data among others because they have daily data of all the variables necessary for the calculation of the FAO-PM equation. Figure 3 reports mean annual ETo calculated from the reanalysis data and the three RCMs for the reference period 1981-2005, for the five agrometeorological stations considered, with and without bias-correction (described in par. 2.4.2). All RCMs overestimate ETo in all locations when compared to ETo from ERA5-Land. For all stations, the mean annual ETo RMSE is higher and the positive bias is more pronounced in the HIRHAM5 model simulations, The same observations apply also to the ETo means over the summer season (June to August). Similarly, all RCMs overestimate ETo in all five stations (Fig. 4). However, the positive bias is higher for the RCA4 model instead. It is also interesting to note that RACMO22E with and without bias correction is among the best performing RCMs at all locations and in both time frames. The effect of the distance from the sea on the performance indices is also noticed as San Pancrazio and Cadriano still have the highest RMSE and bias values and Volano the lowest ones. In the summer season, bias varies between 17 and 69% and RMSE values range between 14 and 43 mm/season. Applying the QM bias correction technique reduced the biases to a range between − 0.07 and 0.06% and the  RMSE to a range between 6 and 7 mm/season, which is negligible.
A better evaluation of model performances may be obtained using a Taylor diagram that summarizes the relative skill with which the different RCMs simulate ETo in the five locations in comparison with the reanalysis outputs over the period 1981-2005. Looking at Fig. 5, note that all RCM simulations whether raw or bias-corrected, exhibit high correlation coefficients (above 0.97). On the other hand, a difference is seen when looking at RMSD and the standard deviation. Those metrics are lower for the bias-corrected RCMs, hence their closeness to the reference point. In contrast to the raw RCMs symbols, those representing the biascorrected RCMs are packed on the red arc meaning that they have the same standard deviation than the reanalysis data, which indicates that their variations are of the right amplitude. The previous observation about RACMO22E being the best performing RCM with and without bias correction for all locations is also confirmed. It is evident that the application of the bias correction substantially improved the ETo simulations by the different RCMs when compared to the ERA5-Land.

ETo projections
The same approach applied to the reference period was also applied to bias-correct RCMs data for the period 2021-2050 following RCP4.5 and RCP8.5 scenarios. Then, an ensemble for each scenario from the three bias-corrected RCMs was built to study the projected ETo changes with respect to the reference period. The combination of information provided by an ensemble approach from different RCMs simulations instead of only one specific RCM is a way to reduce further uncertainties associated with climate models projections and to produce more robust results (Christensen et al. 2010;Kendon et al. 2010;Giménez and García-Galiano 2018).

Future change in temperature
Before analyzing what may happen in different scenarios to ETo, we want to look at what is projected for mean air temperature in the different scenarios. For the RCP4.5 scenario (Table 3), annual T mean will range between 14.2 ± 0.5 °C in San Pancrazio and 15.8 ± 0.4 °C in Volano, while, in summer, it will range between a minimum of 24.4 ± 0.9 °C in San Pancrazio and a maximum of 25.3 ± 0.9 °C in San Pietro Capofiume. This implies a general increase of annual T mean by about 1 °C in 2021-2050 with respect to the reference period 1981-2005, going up to 1.3 °C in summer. For RCP8.5 scenario, the annual average value of T mean is projected to fluctuate between 14.5 ± 0.5 °C in San Pancrazio and 16.1 ± 0.4 °C in Volano. In summer, T mean will range between 24.7 ± 0.7 °C in San Pancrazio and 25.5 ± 0.8 °C in San Pietro Capofiume. Compared to the reference period, the increase of the average annual value will be by around 1.2 °C in all locations, reaching 1.5 °C in summer. Similar results were found when assessing air temperature changes for the period 2021-2050 with respect to the reference period 1982-2011 in the Po river basin (Mercogliano et al. 2014;Vezzoli et al. 2016). D'Oria et al. (2018 analyzed the climate change effects on temperature over the Taro, Parma, and Enza River basins, in the Emilia Romagna region, using an ensemble of 13 Regional Climate Models and the same emission scenarios (RCP4.5 and RCP8.5). At annual scale, they expected increments up to + 0.75 °C in the 2016 − 2035 period and + 1.5 °C in the 2046-2065 period under the RCP4.5, and higher, up to + 1 °C in 2016-2035 and + 2 °C in 2046-2065 with the RCP 8.5, using 1985-2005 as a reference period. These results are also close to the national projections that anticipate an increase in T mean ranging between 1.25 and 1.75 °C in 2021-2050 for RCP4.5 and between 1.5 and 2.0 °C for  (Desiato et al. 2005). Overall, the magnitudes of air temperature increase for RCP4.5 and RCP8.5 are close, due to the small difference between the two emission pathways before 2050. The same observation was made by Wang and Chen (2014) when assessing the magnitudes of temperature increase in China for the same RCP scenarios.

Future changes in solar radiation, wind speed, and dew point temperature
Certainly, air temperature has an important influence on evapotranspiration. However, other weather variables, involved in the FAO-PM equation, have also an effect on ETo. Thus, it is useful to assess the future change in solar radiation, wind speed, and dew point temperature for the two considered RCP scenarios. For both RCP4.5 and RCP8.5, annual solar radiation will increase by around 1.3% in all locations during the future period 2021-2050 (Table 4), while summer solar radiation will increase by around 0.8 to 1.3%. However, according to Gutiérrez et al. (2020), the amplitude of the changes in surface solar radiation likely depends on the RCM and on its aerosol forcing choice, suggesting that it is important to be cautious when using CORDEX projections for solar radiation and arguing for the inclusion of aerosol forcing evolution in the next generation of CORDEX simulations. Besides, such a small increase will likely not affect much ETo.
For wind speed (Table 5), annual values show an average increase of around 14% in all locations and for both RCP scenarios, with regards to the reference period 1981-2005. On the other hand, summer values will stay constant. It is known that wind speed has a positive effect on ETo, and its increase contributes to an increase of the reference evapotranspiration (Allen et al. 1998).
Finally, an increase of around 9.2 to 12.3% in annual T dew is projected for all locations during 2021-2050 in comparison with the reference period (Table 6). Similarly, summer T dew will increase by 5.6 to 7.1%. This means an increase in the quantity of water vapor, which has the effect of decreasing evapotranspiration, offsetting, at least in part, the effect of the increase in T mean and WS. The impact on ETo projections of the variation of all these quantities, to which an increase of the CO 2 level and consequently of the stomatal resistance is added, is studied in the following paragraphs with a particular focus on the change in air temperature as it is used often as a major argument against the effect of a change in CO 2 concentrations.

Future reference evapotranspiration scenarios
For the estimation of the future situation of the ETo in 2021-2050, with respect to the reference period 1981-2005, four scenarios were considered: the first scenario considers RCP4.5 without changing the CO 2 concentration in the FAO-PM equation, which was about 372 ppm when it was first developed; the second scenario considers RCP8.5 without changing the future CO 2 concentration in the FAO-PM Eq. (372 ppm); the third scenario uses RCP4.5 with an increase of the canopy resistance to account for a CO 2 concentration of 550 ppm; and the fourth scenario considers RCP8.5 taking in consideration a future CO 2 concentration of 550 ppm. They are summarized in Table 7.
Reference evapotranspiration calculated by the FAO-PM equation using future data and the four scenarios of Table 7 is reported in Table 8. The average annual ETo   in 1961-1990 and between 778 and 833 mm in 2011-2070. The highest annual and summer ETo values are observed in Volano for both the historical and future period, while the lowest values are observed in San Pancrazio. Based on the data available, this could be partially explained by the former location having the highest air temperature, wind speed, and solar radiation. (Tables 3, 4, and 5), while the latter location having the lowest of those values. Similar conclusions were found by Goyal (2004) and Liu and Zhang (2013) when assessing ETo in China and India respectively.
For scenario ET 1 , an increase of annual T mean by around 1 °C in all locations in 2021-2050 will contribute to increasing annual ETo by an average of 4% ranging between 3.5% in San Pancrazio and 4.4% in Volano. An increase of summertime T mean by around 1.3 °C will be followed by an increase in ETo by an average of 4.5% with a minimum change of 3.5% in San Pancrazio and a maximum of 5% in Volano.
For scenario ET 2 , an increase of annual T mean by around 1.2 °C in 2021-2050 will be met by an increase in ETo by an average of 5.1% with the minimum is San Pancrazio (4.8%) and the maximum in Volano (5.3%). The same pattern is projected in summer when an increase of around 1.5 °C in T mean will help to increase ETo by 5.3% reaching its minimum in San Pancrazio (4.9%) and its maximum in Volano (5.7%). Overall, considering both scenarios, the absolute magnitude of variation from the reference annual and summertime ETo for all locations was between 3.5 and 5.7% for an increase of T mean by 1 to 1.5 °C, which represent an average increase of 4 to 5.4% with regard to the reference period 1981-2005. According to Allen (2000), under wellwatered conditions, evapotranspiration will increase about 4 to 5% per 1 °C rise in temperature. In the Mediterranean region, Saadi et al. (2015) predicted that an increase of air temperature of 1.57 ± 0.27 °C would increase annual reference evapotranspiration by 6.7%. These observations are also proportional to findings by Goyal (2004) in India where they showed that evapotranspiration will increase by 14.8% with an increase in temperature by 20%. They are also proportional to findings by Ramírez and Finnerty (1996) that indicated that a 3 °C increase in temperature induces a 14% increase in evapotranspiration rates in Colorado.
When increasing atmospheric CO 2 concentration, which means increasing the canopy resistance component of the FAO-PM equation to account for a 550 ppm CO 2 concentration, the result is a counterbalance of the effect of the change in the different weather variables, with a quite small difference between actual and future ETo, both for yearly and summer quantities. In particular, for scenario ET 3 , this means that annual ETo will increase by only a maximum of 0.8% in San Pancrazio with regards to the reference period and it will even decrease by a maximum of 1.1% in Volano. In summertime, ETo is anticipated to increase between 0.4% in Volano and 2% in San Pietro Capofiume.
For scenario ET 4 , the increase of CO 2 concentration is projected to have the same effect as scenario ET 3 . Annual ETo will increase by a maximum of 1.1% in San Pancrazio while it will decrease by 0.1% in Volano in comparison with the reference period. On the other hand, in summer, the increase is predicted to be between 1.2% in Martorano and 2.2% in San Pietro Capofiume. It is interesting to notice that, unlike scenario ET 1 and ET 2 , the windiest location (Volano), had the least increase in both annual and summer ETo with regards to the historical period when compared to San Pancrazio, the least windy station. This is because an increase in CO 2 concentration increases the factor 0.34 that accompanies wind speed in the denominator of Eq. (1) causing smaller ETo increments compared to the situation without CO 2 variation . Similar observations were reported by Snyder (2017) when he found that, with an increase in CO 2 concentration, the annual ETo increased slightly where there were mean wind speeds less than 1.7 m·s −1 , and it decreased for wind speeds greater than 1.7 m·s −1 . However, in our case, the change in ETo does not follow a clear wind speed negative gradient, as there are also several other factors that might affect at a certain extent the magnitude of the variation, such as dew point temperature and solar radiation.
Overall, when the effect of an increase in CO 2 concentration from 372 to 550 ppm was considered in combination with the projected increase in temperature by 1 to 1.5 °C, changes in both annual and summer ETo demand varied from − 1.1 to 2.2% during the 2021-2050 period with regard to the reference period 1981-2005. Comparable results were obtained in Colorado by Islam et al. (2012) whose simulations showed that the effect of an increase in CO 2 levels up to 450 ppm to offset the effect of about 1 °C rise in temperature, whereas a 2 °C rise in temperature was offset by a doubling of the CO 2 concentration (660 ppm). Similarly, simulations by Priya et al. (2014) concluded that the effect of a 2.5 °C rise in temperature in Varanasi (India) was offset by an increase of the CO 2 concentration from 330 to 660 ppm.
Accounting for the projected increase in CO 2 concentration in the FAO-PM equation should give a more realistic estimation of WUE and irrigation requirements, especially for C 3 crops such as wheat, barley, potato, and sunflower, for which Emilia-Romagna is famous, as C 4 species are less responsive to increased atmospheric CO 2 . This means that farmers should be able to achieve higher crop yields per unit land area with similar or less amounts of water. However, elevated carbon dioxide concentrations may have detrimental effects on agriculture, as they can affect crop quality due to reductions in nutrients, including many that are important for overall health, such as iron, zinc, and protein (Erda et al. 2005;Smith et al. 2018). Besides, they can increase the damage caused by certain pests and boost weeds prevalence, because many weedy species are C 3 plants, which are favored in an elevated carbon dioxide environment (Ramírez and Finnerty 1996;Fuhrer 2003). This is enhanced even further by increasing competitiveness and reducing the herbicide sensitivity of some weed species (Matzrafi 2019). Besides, temperature increases along with CO 2 will accelerate plants' maturity and hinder pollination (Boote et al. 2005). Furthermore, if temperatures continue to rise, the overall WUE could actually decrease due to increased water requirements in warmer climates such as that of Emilia-Romagna and to possible seed yield reductions caused by higher temperatures (Allen 2000). In fact, CO 2 fertilizing effects should decrease once optimal temperatures are exceeded for a range of processes, especially plant water use Backlund et al. 2008;Lovelli et al. 2010;Turral et al. 2011). However, even without taking into consideration the increments in CO 2 levels, the studied locations should not imply dramatic increases in water demand with respect to the 1981-2005 period as ETo increase does not exceed 6% with the projected increase in air temperature. Using water balance and a calibrated Hargreaves and Samani equation (Hargreaves and Samani 1982), a similar conclusion was drawn by Tomei et al. (2010) for the future period 2021-2050 with respect to the 1991-2008 period. They explain that the result is mainly due to a possible increase in spring precipitation that could mitigate the impact of the further increase in temperature and evapotranspiration.

Summary and conclusion
The aim of this paper was to assess the effect of climate change on crop reference evapotranspiration in five locations in the Emilia-Romagna region during the near future period 2021-2050, with a focus on the effects of higher CO 2 levels and agrometeorological parameters on ETo. Since climate change studies usually need between 20 and 30 years of observed data and due to the unavailability or incompleteness of such dataset in our case, the high-resolution ERA5-Land reanalysis dataset (9 km) was evaluated to serve as an alternative. The comparison made against data observations for the period 2008-2018 gave satisfying results, allowing the reanalysis dataset to be used as a long-term reference dataset. This dataset was used to correct for systematic biases in the daily weather variables extracted from three high-resolution EURO-CORDEX RCMs (12 km), then used for the calculation of ETo. The three-member RCM ensembles assuming a medium (RCP4.5) and a high range emission scenario (RCP8.5) were assessed for the future period 2021-2050 with respect to the reference period 1981-2005. According to the bias-corrected ensemble simulations and both RCP scenarios, the mean air temperature in the five study locations is anticipated to increase by 1 to 1.5 °C annually and during summertime with respect to the 1981-2005 period. This increase is coincident with an increase in the other weather variables involved in the FAO-PM equation, mainly wind speed and dew point temperature. Focusing on summer season, which is important for irrigation, and maintaining the preset CO 2 concentration of 372 ppm, the FAO-PM equation predicts the increase in future summer ETo to range between 3.5 and 5.7% in all stations. By taking into account the increase of CO 2 concentration up to 550 ppm in the ETo calculations, the new range of summer ETo variation will be between 0.4 and 2.2%, thus reducing the effect of an increase in the different weather variables, especially air temperature.
Overall, the results of this study favor the adoption of reanalysis data to estimate ETo when full observation datasets are not available. Besides, it shows that taking into consideration the future projected CO 2 concentrations is important to have a more accurate idea about future crop water demands and to be able to have more efficient water resources management strategies in Emilia-Romagna. As a future step, this work could be extended to the end of the century scenario, when important increases in temperature in all seasons are projected to occur over the Emilia-Romagna region (Tomozeiu et al. 2007(Tomozeiu et al. , 2018. Consent to participate All the authors read and approved the final version of the manuscript. The authors agree to publish the article.

Consent for publication
The work described has not been published before and also is not under consideration for publication elsewhere.

Conflict of interest The authors declare no competing interests.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.