Article Combining Satellite Remote Sensing Data with the FAO-56 Dual Approach for Water Use Mapping In Irrigated Wheat Fields of a Semi-Arid Region

The aim of this study was to combine the FAO-56 dual approach and remotely-sensed data for mapping water use (ETc) in irrigated wheat crops of a semi-arid region. The method is based on the relationships established between Normalized Difference Vegetation Index (NDVI) and crop biophysical variables such as basal crop coefficient, cover fraction and soil evaporation. A time series of high spatial resolution SPOT and Landsat images acquired during the 2002/2003 agricultural season has been used to generate the profiles of NDVI in each pixel that have been related to crop biophysical parameters which were used in conjunction with FAO-56 dual source approach. The obtained results showed that the spatial distribution of seasonal ETc varied between 200 and 450 mm depending to sowing date and the development of the vegetation. The validation of spatial results showed that the ETc estimated by FAO-56 corresponded well with actual ET measured by eddy covariance system over test sites of wheat, especially when soil evaporation and plant water stress are not encountered.


Introduction
Estimates of land surface evapotranspiration (ET) using remote sensing data are essential in effective planning of irrigation water use in arid and semiarid regions.The Haouz plain that surrounds the city of Marrakech (Central Morocco) is classified among the regions in the country facing strong water shortages.This is mainly due the combined effect of persistent drought and the increase of water demands related to increases in irrigated surfaces, urbanization and tourism recreational projects.There is thus a crucial need to develop tools for better management of irrigation water use through accurate estimates of crop water requirement (ET c ) at a regional scale.
Over the last decades, several methods for estimating ET based on combination of crop modelling and remotely-sensed data have been developed.These methods ranged from simple ones such as the Vegetation Index/Temperature trapezoid [1], to complex land surface models [2].However, due to low frequency of the required high resolution satellite data combined the intermittent presence of cloud makes these approaches of little interest for operational purposes [3].
The most common and practical approach used for estimating crop water requirement is the FAO-56 method [4].For the operational monitoring of soil-plant water balance, the FAO-56 model is often preferred to complex soil-vegetation-atmosphere-transfer models.In the FAO-56 approach, ET c is estimated through the combination of a reference evapotranspiration (ET 0 ) and crop coefficients.There are two different FAO-56 approaches: single and dual crop coefficients.In the single crop coefficient approach, the plant transpiration and soil evaporation are combined into a single crop coefficient (K c ).The dual crop coefficient approach uses two coefficients to separate the respective contribution of plant transpiration (K cb ) and soil evaporation (K e ).
Remotely sensed spectral reflectance may provide an indirect estimate of crop coefficient or basal crop coefficients.Indeed, several authors have shown similarity between the seasonal patterns of vegetation indices and transpiration over annual crops [5][6][7][8][9][10][11][12].Consequently, K c can be estimated from spectral vegetation indices since both of them are related to leaf area index and fractional ground cover [13][14][15].In this context, several relationships between K c and vegetation indices have been established.However, there is no agreement on the nature and generality of these relationships [16].Some studies [6][7][8]11,14,17] have shown that these relationships are linear, but others have found non linearity relationships [9,10,12,16].Therefore, establishing a unique relationship between crop coefficient and spectral vegetation indices is an ongoing research topic.
In this study, we focused on the FAO-56 dual crop coefficient model to derive ET c maps from remote sensing data through the use of relationships between the Normalized Difference Vegetation Index (NDVI) and crop biophysical variables such as basal crop coefficient, cover fraction.The objective is to incorporate these relationships into the FAO-56 dual crop coefficient model to derive seasonal variation of evapotranspiration maps.This approach has been previously evaluated at a local scale using local field measurements [12].Here, our aim is to extend the approach to larger scale using NDVI data derived from SPOT and Landsat images during the 2002/2003 agricultural season.The approach is implemented over an irrigated district located in the Haouz plain, 40 km East of Marrakech (Central Morocco), as part of the SudMed project [18].

Study Area
The region of interest is an irrigated area located in the Haouz plain in the centre of the Tensift basin (Central Morocco), 40 km east of Marrakech city.The climate is of semi arid Mediterranean type with an average annual precipitation of about 250 mm of which 70 % falls during winter and spring.The soil type is clay to loam.The area covers about 2,800 ha and is mostly flat.The main land cover classes are cereals; mostly wheat, then barley (Figure 1).The distribution of dam irrigation water is organized by a regional public agency (ORMVAH) together with the local farmer associations.For wheat crops during the 2002-2003 seasons, three irrigation rounds were decided on, with an amount of around 30 mm in equivalent water depth each time.Field 1 was irrigated six times using ground water.More details about the site description and experimental data are given in [11,12].

Data Description
Meteorological measurements were recorded within the irrigated area using an automated weather station consisting of measurements of solar radiation, air temperature and humidity, wind speed, net radiation and rainfall.Daily averaged values of climatic data were calculated in order to compute the daily reference evapotranspiration ET 0 (mm/day), according to the FAO-56 Penman-Monteith parameterization [4].Three eddy-covariance systems were also installed over three fields of wheat to provide continuous measurements of evapotranspiration flux (Figure 1).These data were used to validate spatial ET estimates zooming on each of these fields.Detail on the climatic and flux measurements can be found in [12].
In addition of climatic and flux measurements, a time series of 20 high spatial resolution images acquired by Landsat and SPOT was collected during the growing season of wheat.Due to cloudiness or uncertainty in atmospheric corrections, only 10 images have been used in this study.These images were radiometrically calibrated and atmospherically corrected based on the reflectance of an invariant objects and transformed to NDVI maps [19].The NDVI was derived from red and near infrared reflectances bands after resampling of the time series of images at 20 m.The geometric correction of satellite images was performed using the Ground Control Points collected during the experiment.The radiometric correction was achieved in three steps.Firstly, we corrected a "reference" image from atmospheric effects using the SMAC correction algorithm [20] and standard values of atmospheric components.Secondly, the radiometry of each image were homogenised against this reference images thanks to a set of reliable invariant features [21].This normalisation was performed by applying linear relationships which were established on these invariants between the digital numbers of raw images and the reflectance values of the reference image.Thirdly, an additional linear correction has been applied between the satellite NDVI values and the NDVI ground measurements collected using the hand-held Cropscan MultiSpectral Radiometer [12].This inter-calibration ensures a maximal agreement between satellite and surface NDVI values.

Model Description
The model used in this study is the FAO-56 dual crop coefficient developed by [4].This model describes the relationship between the crop evapotranspiration under standard non-stressed conditions (ET c ) and reference evapotranspiration (ET 0 ) by separating crop coefficient (K c ) into the basal crop coefficient (K cb ) and soil water evaporation (K e ) coefficients: where ET 0 is calculated at daily time step by the FAO Penman-Monteith equation (Equation 6in FAO-56 paper).In order to integrate the remote sensing data into FAO model, the others parameters (K cb and K e ) for the equation 1 were derived from NDVI by equations 2 and 3, respectively.More details about these equations are available in [12].The local calibration performed in [12] was also retained here, the only difference being the source of the NDVI data (satellite or hand-held radiometer).Thus there is no contradiction to state that local calibration is necessary and to find that the same method driven by satellite data is accurate: where NDVI min and NDVI max are the minimum and the maximum values of NDVI associated with bare soil and dense vegetation, respectively.The values 0.14 and 0.93 are used in this study.f c is the vegetation fraction cover given by [12]: Equation 2 presents the general form of the equation for deriving K cb from vegetation indices as shown by [15,16].An empirical function between K cb and NDVI, consisting of two regressions, has been found by [9] for cotton crops.Other authors [6,14,17,22] have found linear K cb -VI relationships.
For Equation 3, it is a simple formalism to account for the fact that soil evaporation occurs predominantly from the exposed soil.Soil evaporation coefficient (K e ) is clearly correlated with (1-f c ).This equation showed that K e decreases with the development of crop (increase of fraction of ground cover, f c ) and becomes negligible when the crop is well developed and completely covers the soil (f c = 1).The relation between K e and (1-f c ) is also modulated with the value 0.25.The value 0.25 was determined according to the figure 29 in the FAO-56 book [4] based on the observed frequency of water supply (≈10 days) and the average value of ET 0 (4 mm/day) during the growing season.

Results and Discussion
Before applying the FAO dual approach in each pixel, we masked the images to keep only wheat fields.A binary mask (wheat) was built from the land cover map [11], and was superimposed on surface NDVI images to keep only wheat pixels.After masking the wheat planting, we used an unsupervised classification method (K-means, [23]) to regroup the pixels which have similar NDVI profiles.50 different classes were extracted (Figure 2).The principle of K-means is summarized in the Appendix.The inter-class variability is due to the differences in agricultural practices (sowing date, irrigation, fertilization…).Plant sowing dates were computed based on the time profiles of NDVI (Figure 3).The sowing date is determined based on the vegetation fraction cover (f c ) derived from NDVI through Equation 4. Sowing dates are determined from emergence dates (identified using a f c = 0.1 criteria) by subtracting 13 days for early sown crops and 20 days for late sown crops.The values (13 days and 20 days) were determined based on the observations in some cultivated fields.The time lag between sowing and emergence is larger for late sown crops because lower temperatures were experienced at this time of year.It can be seen in Figure 3 that there is a large spatial variability in sowing dates, with two distinct periods: early (before December 15) and late (after January 15).This was found consistent with field observations of sowing dates.After, we derived the map of the seasonal ET c (mm) of winter wheat by cumulating daily values of ET c between sowing and full senescence stages (Figure 4).Daily crop water requirement (ET c ) was calculated by multiplying daily K c of daily ET 0 (Equation 1).As mentioned above, K c was derived through the relationships between K cb and NDVI (Equation 2), and between K e and NDVI (Equation 3).In addition, a cubic interpolation was used to determine the values of NDVI between two satellite overpass.Figure 4 shows that the values of seasonal ET c varied between 200 and 450 mm, with an average value of 330 mm.The observed variability of ET c can be mainly explained by two major factors.The first one is the length of crop cycle, due to differences in sowing dates.This is well captured using spectral vegetation index.The cycle is very short for the late sowing classes and larger for early sowing ones.As it can be seen in Figure 2, the amplitude of NDVI (or K c ) is lower for the late sowing classes (<0.5) and higher for early sowing classes (>0.8).This result is in agreement with other studies [24,25].The second factor is cropping practices i.e., irrigation, application or not of nutrients, wheat variety.Globally, the obtained values of ET c are also in accordance with others studies [26,27], which found that the ET c of wheat cultivated in Morocco is around 480 mm depending to the wheat variety and sowing date.The obtained results reveal the potential of remote sensing data to estimate crop water requirements (ET c ) on an operational basis and consumption at a regional scale.To assess the performance of the FAO-56 model, we compared ET c calculated by this approach with the measured ET over three fields of wheat (Figure 5a).According to this figure, it is safe to state that the FAO-56 correctly simulates ET c using remote sensing data.The corresponding Root Mean Square Error (RMSE) was about 0.88 mm per day can be considered acceptable regarding to the average values of ET c (about 3.6 mm per day).The coefficient of determination (R 2 ) is around 0.5, the efficiency (E) is around 0.4 and the slope is close to 1.
However Figure 5a shows that some discrepancies between observed and simulated ET occurred under specific conditions.To identify the reasons of these discrepancies, we analysed soil moisture data collected in these fields [12].According to soil moisture data, the model seems to over-estimate ET when the plant is stressed and under-estimate ET when the soil fully evaporates.This is because the FAO model simulates the evapotranspiration under no limiting water supply and with an average value for soil evaporation (0.25) while the measurements are for actual conditions (presence of the stress).The approximation of K e is not satisfactory for spatial monitoring.It may vary greatly regarding the large heterogeneity in the water supply in terms of amount and periods of distribution.
In order to quantify the gain due to the driving of K c by NDVI, the comparison between (ET 0 , i.e., K c = 1) and actual ET measured at three fields of wheat using eddy covariance system has been made (Figure 5b).According to this comparison, measured ET appeared evidently overestimated and there is no clear relationship between the two variables.The coefficient of determination R 2 is equal to 0.02, showing that the two variables are not linearly related.The corresponding RMSE (1.9 mm) and the efficiency (E = −1.77)are also not acceptable.
In order to test the accuracy of our approach, we compared the ET c provided by this model under standard conditions with the measurements for wheat field 1.This field was monitored during the period of maturity (high vegetation recovering, thus reduced soil evaporation) and was irrigated six times to prevent plant water stress.In this case, there is a close agreement between simulated and measured ET c (Figure 6a); the coefficient of determination and the efficiency are larger than 0.85 and the RMSE is around 0.40 mm per day.Knowing that the K c values of the mature wheat under standard non-stressed conditions are close to 1 [4], it is of interest to check this and to quantify also the contribution of K c in calculating crop evapotarnspiration (ET c ).For this purpose, we plotted (ET 0 , i.e., K c = 1) against actual ET measured over field 1 under standard conditions (Figure 6b).This comparison revealed practically perfect agreement almost similarly to figure 6a (R 2 = 0.84, with E = 0.80 and RMSE = 0.48 mm) between the two variables, justifying that the contribution of K c was minimal since the crop is mature and moisture is unlimited.This is in agreement with other studies (e.g., [10,12]) when they found that the K c derived from NDVI was quasi-equal to 1 for the mature wheat developed under standards conditions.The difference is the source of NDVI data; the ground remotely-sensed measurements have been used for deriving K c instead of satellite remote sensing data used in this study.Based on the calculation of ET c , the amount of irrigation needed for the wheat was derived as the difference between ET c and precipitation.The spatial and temporal distributed of irrigation needs was presented in Figure 7. Temporally, the irrigation amount appears the largest in March and April when it is the peak period of winter wheat development and the reference evapotranspiration begins to be large.It is the smallest in December and May, when winter wheat was in the initial and senescence stages.The obtained cumulative amount of irrigation needs during the whole growing season varied between 50 and 220 mm, depending of sowing date.This value is centered on (about 135 mm) the amount which has been supplied during this agricultural season by the regional office which is in charge of the distribution of dam water for irrigation.Maps of irrigations could be used by farmers to assist in water management and irrigation scheduling potentially saving water and improving yield.

Conclusions and Perspectives
This study showed that how remote sensing data can be integrated in FAO dual approach for mapping water use (ET c ) of wheat crop in semi arid region.The method consists of linking the main crop biophysical parameters (basal crop coefficient, cover fraction) to the Normalized Difference Vegetation Index (NDVI).The results showed that remote sensing estimates of ET c compare very satisfactorily with ground measurements, especially when the soil evaporation and plant water stress are negligible.The RMSE between measured and estimated ET was about 0.40 mm per day.
For actual estimates of ET, developments should be done to get information about plant water stress.Additional information such as surface temperature (which is the most indicator of the stress) from TIR sensors can be very useful.In the near future we apply spatially an innovative approach tested at local scale by [28].It consists of assimilating daily ET from TIR data combined to an energy balance model into FAO approach.Finally, it should be noted that our approach was tested only for the wheat crop.The method should be tested over an extended spatially distributed dataset with different crops.

Figure 1 .
Figure 1.Land cover (coloured shapes) and location of three wheat fields (delimited with white rectangles) under study during 2002-2003 agricultural season in the Tensift AL Haouz, Marrakech, Morocco.The towers equipped with evapotranspiration measurement systems (black disks) are located in each field.

Figure 2 .
Figure 2. Profiles of NDVI corresponding to 50 classes of winter wheat identified by Kmeans classification.A cubic interpolation method was used to determine the values of NDVI between two satellite overpass.

Figure 3 .
Figure 3. Map of sowing date of winter wheat.The value 1 corresponds to November 7th 2002.Note that the sowing date was calculated by subtracting 13 days for early sown crops and 20 days for late sown crops from the date corresponding to c f = 0.1.

Figure 4 .
Figure 4. Map of seasonal crop water requirement ET c (in mm) of winter wheat obtained by applying the FAO-56 dual approach.

Figure 5 .
Figure5.Scatter plot between daily estimated evapotranspiration (ET c ) by the FAO-56 dual approach (Figure5a), reference ET 0 (Figure5b), and measured actual (ET) by Eddy covariance technique over three fields of wheat in the Tensift basin (center of Morocco).The graph includes the 1:1 line.

Figure 6 .
Figure 6.Same as Figure5, but for one field (field 1) of mature wheat under standard conditions (no stress) and the soil evaporation was negligible.

Figure 7 .
Figure 7. Spatial and temporal distribution of irrigation amount (mm) of winter wheat in the Tensift ALHaouz basin during 2002-2003 growing season.