<figure><img src="image_1.png" /></figure>

HAL
open science

# Water Footprint of Cereals by Remote Sensing in Kairouan Plain (Tunisia)

Vetiya Dellaly, Aicha Chahbi Bellakanji, Hedia Chakroun, Sameh Saadi, Gilles Boulet, Mehrez Zribi, Zohra Lili Chabaane

► To cite this version:

Vetiya Dellaly, Aicha Chahbi Bellakanji, Hedia Chakroun, Sameh Saadi, Gilles Boulet, et al.. Water Footprint of Cereals by Remote Sensing in Kairouan Plain (Tunisia). Remote Sensing, 2024, 16 (3), pp.491. 10.3390/rs16030491. hal-04477495

HAL Id: hal-04477495

https://hal.inrae.fr/hal-04477495

Submitted on 26 Feb 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_1.png" /></figure>

Article

# Water Footprint of Cereals by Remote Sensing in Kairouan Plain (Tunisia)

Vetiya Dellaly <fnref n="1" />*, Aicha Chahbi Bellakanji <fnref n="2" /> ID, Hedia Chakroun <fnref n="1" /> ID, Sameh Saadi <fnref n="3" /> ID, Gilles Boulet <fnref n="4" /> ID,
Mehrez Zribi <fnref n="4" /> ID and Zohra Lili Chabaane <fnref n="2" /> ID

<figure><img src="image_2.png" /></figure>

check for
updates

**Citation:** Dellaly, V.; Chahbi Bellakanji, A.; Chakroun, H.; Saadi, S.; Boulet, G.; Zribi, M.; Lili Chabaane, Z. Water Footprint of Cereals by Remote Sensing in Kairouan Plain (Tunisia). *Remote Sens.* **2024**, *16*, 491. https://doi.org/10.3390/rs16030491

Academic Editor: Guido D'Urso

Received: 29 November 2023

Revised: 20 January 2024

Accepted: 22 January 2024

Published: 26 January 2024

<figure><img src="image_3.png" /></figure>

**Copyright:** © 2024 by the authors.
Licensee MDPI, Basel, Switzerland.
This article is an open access article
distributed under the terms and
conditions of the Creative Commons
Attribution (CC BY) license (https://
creativecommons.org/licenses/by/
4.0/).

**Abstract:** This article aims to estimate the water footprint (WF) of cereals—specifically, wheat and barley—in the Kairouan plain, located in central Tunisia. To achieve this objective, two components must be determined: actual evapotranspiration (ETa) and crop yield. The study covers three growing seasons from 2010 to 2013. The ETa estimation employed the S-SEBI (simplified surface energy balance index) model, utilizing Landsat 7 and 8 optical and thermal infrared spectral bands. For yield estimation, an empirical model based on the normalized difference vegetation index (NDVI) was applied. Results indicate the effectiveness of the S-SEBI model in estimating ETa, demonstrating an R<fnref n="2" /> of 0.82 and an RMSE of 0.45 mm/day. Concurrently, yields mapped over the area range between 6 and 77 qx/ha. Globally, cereals' average WF varied from 1.08 m<fnref n="3" />/kg to 1.22 m<fnref n="3" />/kg over the three study years, with the majority below 1 m<fnref n="3" />/kg. Notably in dry years, the importance of the blue WF is emphasized compared to years with average rainfall (WFb-2013 = 1.04 m<fnref n="3" />/kg, WFb-2012 = 0.61 m<fnref n="3" />/kg, WFb-2011 = 0.41 m<fnref n="3" />/kg). Moreover, based on an in-depth agronomic analysis combining yields and WF, four classes were defined, ranging from the most water efficient to the least, revealing that over 30% of cultivated areas during the study years (approximately 40% in 2011 and 2012 and 29% in 2013) exhibited low water efficiency, characterized by low yields and high WF. A unique index, the WFI, is proposed to assess the spatial variability of green and blue water. Spatial analysis using the WFI highlighted that in 2012, 40% of cereal plots with low yields but high water consumption were irrigated (81% blue water compared to 6% in 2011).

**Keywords:** crop evapotranspiration; S-SEBI; yield; water footprint; remote sensing

## 1. Introduction

The concepts of green, blue, and gray water have been widely debated over the last two decades [1–3], and they lead us to revisit discussions on crop production, irrigation, water use efficiency, water productivity (WP) and WF. Specifically, the term “green water” commonly denotes the use of rainwater in the soil or rainwater evaporation [4], while “blue water” encompasses groundwater and surface water resources brought to the crops in the irrigation process. On the other hand, “gray water” refers to water that becomes polluted during production processes, as outlined by [5]. Within the hydrological cycle, blue, green, and gray water are used in the WF evaluation of crops, defined as the volume of fresh water used in crop production measured across the entire supply chain. For the purpose of this study, our focus lies solely on the blue and green components. The assessment of green

<footnote n="1">¹ LR99ES19 Laboratory of Modelling in Hydraulics and Environment (LMHE), National Engineering School of Tunis (ENIT), University of Tunis El Manar, BP 37, Tunis 1002, Tunisia; hedia.chakroun@enit.utm.tn</footnote>

<footnote n="2">² LR17AGR01 InteGRatEd Management of Natural Resources: Remote Sensing, Spatial Analysis and Modeling (GREEN-TEAM), National Agronomic Institute of Tunisia, Carthage University, 43 Avenue Charles Nicolle, Tunis 1082, Tunisia; aicha.chahbi@inat.ucar.tn (A.C.B.); zohra.lilichabaane@inat.ucar.tn (Z.L.C.)</footnote>

<footnote n="3">³ Dynafor UMR 1201, Engineering School of Purpan, University of Toulouse, 75 voie du TOEC, BP57611, 31076 Toulouse, Cedex 3, France; sameh.saadi@purpan.fr</footnote>

<footnote n="4">⁴ CNES/CNRS/INRAE/IRD/UT3-Paul Sabatier, CESBIO, University of Toulouse, 18 Avenue Edouard Belin, 31401 Toulouse, France; gilles.boulet@ird.fr (G.B.); mehrez.zribi@ird.fr (M.Z.)</footnote>

<footnote>* Correspondence: vetiya.dellaly@enit.utm.tn</footnote>
<figure><img src="image_3.png" /></figure>

Article

# Water Footprint of Cereals by Remote Sensing in Kairouan Plain (Tunisia)

Vetiya Dellaly <fnref n="1" />*, Aicha Chahbi Bellakanji <fnref n="2" /> ID, Hedia Chakroun <fnref n="1" /> ID, Sameh Saadi <fnref n="3" /> ID, Gilles Boulet <fnref n="4" /> ID,
Mehrez Zribi <fnref n="4" /> ID and Zohra Lili Chabaane <fnref n="2" /> ID

<figure><img src="image_4.png" /></figure>

check for
updates

**Citation:** Dellaly, V.; Chahbi Bellakanji, A.; Chakroun, H.; Saadi, S.; Boulet, G.; Zribi, M.; Lili Chabaane, Z. Water Footprint of Cereals by Remote Sensing in Kairouan Plain (Tunisia). *Remote Sens.* **2024**, *16*, 491. https://doi.org/10.3390/rs16030491

Academic Editor: Guido D'Urso

Received: 29 November 2023

Revised: 20 January 2024

Accepted: 22 January 2024

Published: 26 January 2024

<figure><img src="image_5.png" /></figure>

**Copyright:** © 2024 by the authors.
Licensee MDPI, Basel, Switzerland.
This article is an open access article
distributed under the terms and
conditions of the Creative Commons
Attribution (CC BY) license (https://
creativecommons.org/licenses/by/
4.0/).

**Abstract:** This article aims to estimate the water footprint (WF) of cereals—specifically, wheat and barley—in the Kairouan plain, located in central Tunisia. To achieve this objective, two components must be determined: actual evapotranspiration (ETₐ) and crop yield. The study covers three growing seasons from 2010 to 2013. The ETₐ estimation employed the S-SEBI (simplified surface energy balance index) model, utilizing Landsat 7 and 8 optical and thermal infrared spectral bands. For yield estimation, an empirical model based on the normalized difference vegetation index (NDVI) was applied. Results indicate the effectiveness of the S-SEBI model in estimating ETₐ, demonstrating an R<fnref n="2" /> of 0.82 and an RMSE of 0.45 mm/day. Concurrently, yields mapped over the area range between 6 and 77 qx/ha. Globally, cereals' average WF varied from 1.08 m<fnref n="3" />/kg to 1.22 m<fnref n="3" />/kg over the three study years, with the majority below 1 m<fnref n="3" />/kg. Notably in dry years, the importance of the blue WF is emphasized compared to years with average rainfall (WFb-2013 = 1.04 m<fnref n="3" />/kg, WFb-2012 = 0.61 m<fnref n="3" />/kg, WFb-2011 = 0.41 m<fnref n="3" />/kg). Moreover, based on an in-depth agronomic analysis combining yields and WF, four classes were defined, ranging from the most water efficient to the least, revealing that over 30% of cultivated areas during the study years (approximately 40% in 2011 and 2012 and 29% in 2013) exhibited low water efficiency, characterized by low yields and high WF. A unique index, the WFI, is proposed to assess the spatial variability of green and blue water. Spatial analysis using the WFI highlighted that in 2012, 40% of cereal plots with low yields but high water consumption were irrigated (81% blue water compared to 6% in 2011).

**Keywords:** crop evapotranspiration; S-SEBI; yield; water footprint; remote sensing

## 1. Introduction

The concepts of green, blue, and gray water have been widely debated over the last two decades [1–3], and they lead us to revisit discussions on crop production, irrigation, water use efficiency, water productivity (WP) and WF. Specifically, the term “green water” commonly denotes the use of rainwater in the soil or rainwater evaporation [4], while “blue water” encompasses groundwater and surface water resources brought to the crops in the irrigation process. On the other hand, “gray water” refers to water that becomes polluted during production processes, as outlined by [5]. Within the hydrological cycle, blue, green, and gray water are used in the WF evaluation of crops, defined as the volume of fresh water used in crop production measured across the entire supply chain. For the purpose of this study, our focus lies solely on the blue and green components. The assessment of green

<footnote n="1">¹ LR99ES19 Laboratory of Modelling in Hydraulics and Environment (LMHE), National Engineering School of Tunis (ENIT), University of Tunis El Manar, BP 37, Tunis 1002, Tunisia; hedia.chakroun@enit.utm.tn</footnote>

<footnote n="2">² LR17AGR01 InteGRatEd Management of Natural Resources: Remote Sensing, Spatial Analysis and Modeling (GREEN-TEAM), National Agronomic Institute of Tunisia, Carthage University, 43 Avenue Charles Nicolle, Tunis 1082, Tunisia; aicha.chahbi@inat.ucar.tn (A.C.B.); zohra.lilichabaane@inat.ucar.tn (Z.L.C.)</footnote>

<footnote n="3">³ Dynafor UMR 1201, Engineering School of Purpan, University of Toulouse, 75 voie du TOEC, BP57611, 31076 Toulouse, Cedex 3, France; sameh.saadi@purpan.fr</footnote>

<footnote n="4">⁴ CNES/CNRS/INRAE/IRD/UT3-Paul Sabatier, CESBIO, University of Toulouse, 18 Avenue Edouard Belin, 31401 Toulouse, France; gilles.boulet@ird.fr (G.B.); mehrez.zribi@ird.fr (M.Z.)</footnote>

<footnote>* Correspondence: vetiya.dellaly@enit.utm.tn</footnote>
and blue water components is important, as it allows us to distinguish between rainfed and irrigated crops' needs satisfaction.

In Tunisia, a region facing water scarcity (with approximately 5.109 m³/year of water resources are entirely mobilized) [6], assessing the water footprint (WF) provides insights into water usage efficiency. This evaluation is crucial for optimal freshwater resource allocation amid prevailing water scarcity. Despite allocating nearly 78% of available blue water resources in the country to irrigation, a significant portion of food production relies on rainfed agriculture. Notably, irrigated crops, though having higher added value, lag in water resource planning considerations [6,7]. Holistic approaches on the economic value of water in Tunisia have been developed in [2,7,8]. For cereals, deemed an essential crop in Tunisia, the WF is expressed in cubic meters per kilogram, and its inverse, water productivity (WP), expressed as yield grain per unit of consumption by evapotranspiration (kilograms/cubic meter) of cereals, have been studied extensively [9–13]. These studies, based on existing methods, computed the WF by utilizing data from national statistics, reports and climatic databases. However, the incorporation of remote sensing data enhances both the spatial and temporal coverage to traditional in situ measurements. Notably, the spatial resolution sees enhancement, particularly when contrasted with statistics-driven or interpolated climatic station measurements, thereby facilitating the monitoring of spatial parameter variability.

On a regional scale, remote sensing has made it possible to develop tools for the global assessment of crops' water footprint (e.g., the FAO-WAPOR allows you to evaluate the WP, and among other things, the WF; available at http://www.fao.org/in-action/remote-sensing-for-water-productivity/en/, accessed on 12 April 2022). However, at local scales, accurate knowledge of spatiotemporal tracking of main crop water need remains insufficient.

A good estimate of water footprint requires a precise evaluation of the two main components, namely, crop water requirement (CWR) and crop yield for both rainfed and irrigated cases. The ETₐ allows the determination of CWR using either empirical methods or residual methods based on surface energy balance. In this study, for calculating ETₐ we have chosen to use two methods, as follows. (1) The most used empirical method is the FAO-56 [14,15] based on the reference ET₀, which depends on atmospheric demand to estimate CWR. To approximate ET₀, the standard method is the Penman–Monteith equation [16]. When there is no water limitation, the CWR equals the crop ETₐ, considering the crop characteristics and the soil factors. The FAO-56 method has been expressed in a single crop coefficient variant based on a unique crop coefficient (Kc) that combines both vegetation transpiration and soil evaporation under unrestricted water availability. An extended review of CWR estimation methods based on optical remote sensing during the last three decades was conducted by [17]. (2) The S-SEBI model (simplified surface energy balance index) [18], one of the physics-based models using residual methods, controls water exchange and partitions it into sensible and latent heat in the soil–vegetation–atmosphere continuum [19]. This model combines remote sensing variables of surface radiometric temperature, albedo, and vegetation cover with ground-based meteorological data to determine the different components of the energy balance, namely, the latent and the sensible heat fluxes. It is easy to use and does not require a lot of input data. Additionally, it has been validated in several semi-arid zones [20,21].

Regarding the crop yield component required for the evaluation of water footprint in crops, the estimation of agricultural production at national or regional scales is more crucial than ever in the current context of climate change and its impact on humans and the environment [22]. Various methods are employed for cereal yield estimation, each having distinct advantages and limitations. In this study, we utilized an approach widely applied at national and regional levels, based on the use of satellite imagery to capture information about the crop's spectral reflectance and vegetation indices throughout the agricultural season. These data were then correlated with ground-based measurements to estimate yield using regression models. The most frequently used index in these models is the normalized
and blue water components is important, as it allows us to distinguish between rainfed and irrigated crops' needs satisfaction.

In Tunisia, a region facing water scarcity (with approximately 5.109 m³/year of water resources are entirely mobilized) [6], assessing the water footprint (WF) provides insights into water usage efficiency. This evaluation is crucial for optimal freshwater resource allocation amid prevailing water scarcity. Despite allocating nearly 78% of available blue water resources in the country to irrigation, a significant portion of food production relies on rainfed agriculture. Notably, irrigated crops, though having higher added value, lag in water resource planning considerations [6,7]. Holistic approaches on the economic value of water in Tunisia have been developed in [2,7,8]. For cereals, deemed an essential crop in Tunisia, the WF is expressed in cubic meters per kilogram, and its inverse, water productivity (WP), expressed as yield grain per unit of consumption by evapotranspiration (kilograms/cubic meter) of cereals, have been studied extensively [9–13]. These studies, based on existing methods, computed the WF by utilizing data from national statistics, reports and climatic databases. However, the incorporation of remote sensing data enhances both the spatial and temporal coverage to traditional in situ measurements. Notably, the spatial resolution sees enhancement, particularly when contrasted with statistics-driven or interpolated climatic station measurements, thereby facilitating the monitoring of spatial parameter variability.

On a regional scale, remote sensing has made it possible to develop tools for the global assessment of crops' water footprint (e.g., the FAO-WAPOR allows you to evaluate the WP, and among other things, the WF; available at http://www.fao.org/in-action/remote-sensing-for-water-productivity/en/, accessed on 12 April 2022). However, at local scales, accurate knowledge of spatiotemporal tracking of main crop water need remains insufficient.

A good estimate of water footprint requires a precise evaluation of the two main components, namely, crop water requirement (CWR) and crop yield for both rainfed and irrigated cases. The ETₐ allows the determination of CWR using either empirical methods or residual methods based on surface energy balance. In this study, for calculating ETₐ we have chosen to use two methods, as follows. (1) The most used empirical method is the FAO-56 [14,15] based on the reference ET₀, which depends on atmospheric demand to estimate CWR. To approximate ET₀, the standard method is the Penman–Monteith equation [16]. When there is no water limitation, the CWR equals the crop ETₐ, considering the crop characteristics and the soil factors. The FAO-56 method has been expressed in a single crop coefficient variant based on a unique crop coefficient (Kc) that combines both vegetation transpiration and soil evaporation under unrestricted water availability. An extended review of CWR estimation methods based on optical remote sensing during the last three decades was conducted by [17]. (2) The S-SEBI model (simplified surface energy balance index) [18], one of the physics-based models using residual methods, controls water exchange and partitions it into sensible and latent heat in the soil–vegetation–atmosphere continuum [19]. This model combines remote sensing variables of surface radiometric temperature, albedo, and vegetation cover with ground-based meteorological data to determine the different components of the energy balance, namely, the latent and the sensible heat fluxes. It is easy to use and does not require a lot of input data. Additionally, it has been validated in several semi-arid zones [20,21].

Regarding the crop yield component required for the evaluation of water footprint in crops, the estimation of agricultural production at national or regional scales is more crucial than ever in the current context of climate change and its impact on humans and the environment [22]. Various methods are employed for cereal yield estimation, each having distinct advantages and limitations. In this study, we utilized an approach widely applied at national and regional levels, based on the use of satellite imagery to capture information about the crop's spectral reflectance and vegetation indices throughout the agricultural season. These data were then correlated with ground-based measurements to estimate yield using regression models. The most frequently used index in these models is the normalized
difference vegetation index (NDVI) [23–26]. Additionally, several techniques based on remote sensing data are commonly employed for cereal yield estimation. Vegetation indices, such as the NDVI or EVI (enhanced vegetation index), quantify vegetation vigor and health based on specific light-wavelength reflectance. These indices provide valuable information about crop biomass and photosynthetic activity, both linked to yield potential. The advantage of this method is its provision of accurate measurements and its efficiency in terms of time and labor, especially for large-scale assessments, as is the case with other methods.

This paper aims to explore the potential of remote sensing data for estimating and spatializing two key components: (1) the actual evapotranspiration (ETₐ) using the S-SEBI model, including its validation and comparison with the FAO-56 method, and (2) the green and blue components of the water footprint (WF) for both irrigated and rain-fed cereals in the semi-arid climate of the Kairouan plain in central Tunisia. Additionally, the assessment of irrigated and rainfed cereals, along with yield, involves the integration of crop growth modeling and satellite imagery.

## 2. Materials and Methods

### 2.1. *Study Area*

The plain of Kairouan is located in the center of Tunisia, in the eastern part of the Tunisian Atlas between 9°30'E and 10°15'E and 35°N and 35°45'N (Figure 1). This plain is characterized by the presence of irrigated perimeters and by arboriculture, particularly the olive tree and cereal growing [27]. It is characterized by a semi-arid climate with an average annual rainfall of 300 mm/year and a reference evapotranspiration of around 1600 mm/year [28]. In central Tunisia, the rainy season extends from September to April, while in summer there is practically no precipitation.

<figure><img src="image_4.png" /><figcaption>**Figure 1.** Location map of Kairouan plain region (Google Earth image).</figcaption></figure>
difference vegetation index (NDVI) [23–26]. Additionally, several techniques based on remote sensing data are commonly employed for cereal yield estimation. Vegetation indices, such as the NDVI or EVI (enhanced vegetation index), quantify vegetation vigor and health based on specific light-wavelength reflectance. These indices provide valuable information about crop biomass and photosynthetic activity, both linked to yield potential. The advantage of this method is its provision of accurate measurements and its efficiency in terms of time and labor, especially for large-scale assessments, as is the case with other methods.

This paper aims to explore the potential of remote sensing data for estimating and spatializing two key components: (1) the actual evapotranspiration (ETₐ) using the S-SEBI model, including its validation and comparison with the FAO-56 method, and (2) the green and blue components of the water footprint (WF) for both irrigated and rain-fed cereals in the semi-arid climate of the Kairouan plain in central Tunisia. Additionally, the assessment of irrigated and rainfed cereals, along with yield, involves the integration of crop growth modeling and satellite imagery.

## 2. Materials and Methods

### 2.1. Study Area

The plain of Kairouan is located in the center of Tunisia, in the eastern part of the Tunisian Atlas between 9°30'E and 10°15'E and 35°N and 35°45'N (Figure 1). This plain is characterized by the presence of irrigated perimeters and by arboriculture, particularly the olive tree and cereal growing [27]. It is characterized by a semi-arid climate with an average annual rainfall of 300 mm/year and a reference evapotranspiration of around 1600 mm/year [28]. In central Tunisia, the rainy season extends from September to April, while in summer there is practically no precipitation.

<figure><img src="image_6.png" /><figcaption>**Figure 1.** Location map of Kairouan plain region (Google Earth image).</figcaption></figure>
Figure 2 depicts the temporal variability in annual rainfall recorded between 2010 and 2013, considering agricultural years from September to August (depending on data availability). The agricultural year 2012–2013 marked the driest period, with total rainfall of 177.2 mm, while the wettest agricultural year was 2011–2012, recording 344 mm of annual rainfall. In the year 2010–2011, a total of 328.8 mm of rainfall was recorded, with the majority occurring towards the end of the agricultural season.

<figure><img src="image_5.png" /><figcaption>**Figure 2.** Rainfall variations over 3 agricultural years (2010–2011, 2011–2012, and 2012–2013) on the plain of Kairouan recorded by the Ben Salem station.</figcaption></figure>

## 2.2. Satellite and Ground Data

### 2.2.1. Satellite Data

A total of 13 cloud-free Landsat 7 (ETM+) and Landsat 8 (OLI + TRIS) images were procured from http://earthexplorer.usgs.gov, accessed on 6 May 2021, aligning with the available images spanning the three seasons of cereal crops from 2010 to 2013 in the Kairouan plain region. Atmospheric and radiometric corrections for these images were conducted using the semiautomated classification (SCP) plugin [29] based on dark object subtraction (DOS) [30]. For Landsat 7 (ETM+) images, the gdal_fillnodata algorithm [31], was employed to address the filling of no data pixels. This algorithm is grounded in inverse distance weighting interpolation and includes a smoothening process to mitigate artefacts. Subsequently, the calculation of biophysical variables, including albedo, NDVI, and emissivity, was carried out. Surface temperature was determined from at-satellite brightness temperature, as outlined in [32].

To monitor surface state parameters and comprehend vegetation dynamics on the Kairouan plain, 20 optical images from SPOT 4 and SPOT 5 spanning agricultural seasons from 2010–2011 to 2012–2013 were utilized for NDVI calculations. SPOT 4 and 5 constitute a multi-spectral, multi-band imaging mission offering high spatial resolution (20 m and 10 m, respectively) and a revisit time of approximately 21 days. Each image consists of four bands: green, red, near-infrared PIR, and mid-infrared MIR. All acquired images underwent orthorectification using the ground control point technique and were atmospherically corrected [33]. Among the available methods for atmospheric correction, the SMAC 6 s method (second simulation of the satellite signal in the solar spectrum) was employed [34–36]. This model represents an enhancement of the 5 s model developed by [37] and incorporates considerations for the bidirectional reflectance of observed surfaces and their altitudes.

Monitoring land cover is essential in the management of water resources, irrigation, and agronomy, which is the focus of our case study. Satellite image classification is a
Figure 2 depicts the temporal variability in annual rainfall recorded between 2010 and 2013, considering agricultural years from September to August (depending on data availability). The agricultural year 2012–2013 marked the driest period, with total rainfall of 177.2 mm, while the wettest agricultural year was 2011–2012, recording 344 mm of annual rainfall. In the year 2010–2011, a total of 328.8 mm of rainfall was recorded, with the majority occurring towards the end of the agricultural season.

<figure><img src="image_7.png" /><figcaption>**Figure 2.** Rainfall variations over 3 agricultural years (2010–2011, 2011–2012, and 2012–2013) on the plain of Kairouan recorded by the Ben Salem station.</figcaption></figure>

## 2.2. Satellite and Ground Data

### 2.2.1. Satellite Data

A total of 13 cloud-free Landsat 7 (ETM+) and Landsat 8 (OLI + TRIS) images were procured from http://earthexplorer.usgs.gov, accessed on 6 May 2021, aligning with the available images spanning the three seasons of cereal crops from 2010 to 2013 in the Kairouan plain region. Atmospheric and radiometric corrections for these images were conducted using the semiautomated classification (SCP) plugin [29] based on dark object subtraction (DOS) [30]. For Landsat 7 (ETM+) images, the gdal_fillnodata algorithm [31], was employed to address the filling of no data pixels. This algorithm is grounded in inverse distance weighting interpolation and includes a smoothening process to mitigate artefacts. Subsequently, the calculation of biophysical variables, including albedo, NDVI, and emissivity, was carried out. Surface temperature was determined from at-satellite brightness temperature, as outlined in [32].

To monitor surface state parameters and comprehend vegetation dynamics on the Kairouan plain, 20 optical images from SPOT 4 and SPOT 5 spanning agricultural seasons from 2010–2011 to 2012–2013 were utilized for NDVI calculations. SPOT 4 and 5 constitute a multi-spectral, multi-band imaging mission offering high spatial resolution (20 m and 10 m, respectively) and a revisit time of approximately 21 days. Each image consists of four bands: green, red, near-infrared PIR, and mid-infrared MIR. All acquired images underwent orthorectification using the ground control point technique and were atmospherically corrected [33]. Among the available methods for atmospheric correction, the SMAC 6 s method (second simulation of the satellite signal in the solar spectrum) was employed [34–36]. This model represents an enhancement of the 5 s model developed by [37] and incorporates considerations for the bidirectional reflectance of observed surfaces and their altitudes.

Monitoring land cover is essential in the management of water resources, irrigation, and agronomy, which is the focus of our case study. Satellite image classification is a
crucial technique employed to generate land-cover maps. Classification involves grouping pixels based on their spectral similarity or determining pixel contours to form spatial units interpretable in terms of classes or thematic categories [38]. Classifications can be based on either the pixel approach or the object approach. The pixel approach groups identical pixels based on their spectral signatures, while the object approach involves segmenting the image into several objects. In this study, we utilized land-cover maps for the three agricultural years under investigation, generated through pixel-based classification to distinguish cereals from other crops. Land-use mapping relied on a decision tree approach [26,28] using SPOT/HRG images for the three agricultural years. We defined eleven land-use classes, including water cover (wadi beds and dams), coastal salt flat areas known as "sebkhas," mountainous areas (reliefs), urban areas, non-irrigated olive trees, irrigated vegetables (summer and winter vegetables), cereals (wheat and barley), and irrigated olive trees. Annually, ground truth fields (approximately 100 fields) were collected in the study area. The first set was used to identify empirical NDVI limits between different vegetation classes, and the second set was used for validating our land-use classification approach [11,19,26]. The identification of NDVI thresholds for each class is carried out during the crop growth period. For instance, for summer crops, a threshold NDVI value of 0.3 is established for images acquired in summer (July, or August). In the case of wheat classes (irrigated or non-irrigated wheat and barley), we utilize empirical NDVI thresholds with images acquired between 15 March and 15 April for all three agricultural years (NDVI > 0.4). For instance, the land-cover map for the 2011–2012 agricultural year (Figure 3) achieves an overall accuracy coefficient of approximately 93.57% and a kappa accuracy coefficient of around 80.4% [26]. In the case of the cereal class, which is of primary interest in our study, 98% of the pixels are accurately classified. The method leverages various types of data, facilitating accurate image classification. Additionally, decision trees prove effective even with a relatively small calibration dataset, providing explicit and easily interpretable outcomes.

<figure><img src="image_6.png" /><figcaption>Figure 3. Land-use map for the 2011–2012 agricultural season.</figcaption></figure>

### 2.2.2. Ground Data

Standard meteorological measurements, including precipitation, air temperature, solar radiation, wind speed, and relative humidity, were recorded at half-hour intervals using an automated weather station situated in the study area (Figure 1). This station is referred to as the Ben Salem meteorological station (35°33′1.44″N; 9°55′18.11″E). Additionally, an eddy covariance (EC) flux station, known as the Ben Salem 1 flux station (located a few dozen meters away from the meteorological station), was used to measure the four components of the surface energy budget, along with soil properties of irrigated wheat, from
crucial technique employed to generate land-cover maps. Classification involves grouping pixels based on their spectral similarity or determining pixel contours to form spatial units interpretable in terms of classes or thematic categories [38]. Classifications can be based on either the pixel approach or the object approach. The pixel approach groups identical pixels based on their spectral signatures, while the object approach involves segmenting the image into several objects. In this study, we utilized land-cover maps for the three agricultural years under investigation, generated through pixel-based classification to distinguish cereals from other crops. Land-use mapping relied on a decision tree approach [26,28] using SPOT/HRG images for the three agricultural years. We defined eleven land-use classes, including water cover (wadi beds and dams), coastal salt flat areas known as "sebkhas," mountainous areas (reliefs), urban areas, non-irrigated olive trees, irrigated vegetables (summer and winter vegetables), cereals (wheat and barley), and irrigated olive trees. Annually, ground truth fields (approximately 100 fields) were collected in the study area. The first set was used to identify empirical NDVI limits between different vegetation classes, and the second set was used for validating our land-use classification approach [11,19,26]. The identification of NDVI thresholds for each class is carried out during the crop growth period. For instance, for summer crops, a threshold NDVI value of 0.3 is established for images acquired in summer (July, or August). In the case of wheat classes (irrigated or non-irrigated wheat and barley), we utilize empirical NDVI thresholds with images acquired between 15 March and 15 April for all three agricultural years (NDVI > 0.4). For instance, the land-cover map for the 2011–2012 agricultural year (Figure 3) achieves an overall accuracy coefficient of approximately 93.57% and a kappa accuracy coefficient of around 80.4% [26]. In the case of the cereal class, which is of primary interest in our study, 98% of the pixels are accurately classified. The method leverages various types of data, facilitating accurate image classification. Additionally, decision trees prove effective even with a relatively small calibration dataset, providing explicit and easily interpretable outcomes.

<figure><img src="image_8.png" /><figcaption>Figure 3. Land-use map for the 2011–2012 agricultural season.</figcaption></figure>

### 2.2.2. Ground Data

Standard meteorological measurements, including precipitation, air temperature, solar radiation, wind speed, and relative humidity, were recorded at half-hour intervals using an automated weather station situated in the study area (Figure 1). This station is referred to as the Ben Salem meteorological station (35°33′1.44″N; 9°55′18.11″E). Additionally, an eddy covariance (EC) flux station, known as the Ben Salem 1 flux station (located a few dozen meters away from the meteorological station), was used to measure the four components of the surface energy budget, along with soil properties of irrigated wheat, from
November 2012 to June 2013. The processing sequence of EC data to calculate turbulent fluxes from raw data is detailed in [33]. Furthermore, field observations and information from farmers regarding irrigated and rainfed crops were determined for 82 plots in the study area between 2010 and 2012. Table 1 provides a synthesis of the datasets and the preprocessing steps required for different models used in water footprint assessment.

<table><caption>**Table 1.** Summary of variables and models required in WF.</caption><thead><tr><th>Model</th><th>Earth Observation</th><th>Ground and Auxiliary Data</th><th>Process</th><th>WF Component</th></tr></thead><tbody><tr><td>S-SEBI</td><td>Landsat 7-ETM+<br>Landsat 8-OLI<br>B2 to B7: 30 m res<br>B10 to B11: 100 m res<br>-LST, Albedo, NDVI, ε<br>(2010–2013)</td><td>—</td><td>Daily:<br>-R<sub>n</sub><br>-Evaporative fraction Λ<br>-Fluxes: H, G<sub>0</sub>, λE</td><td>Yearly cumulative<br>ET<sub>g/b</sub></td></tr><tr><td>FAO-56<br>Single crop coefficient</td><td>—</td><td>-Air temperature<br>-Air humidity<br>-Solar radiation<br>-Rainfall<br>-Pressure</td><td>ET<sub>0</sub><br>(Penman–Monteith, 1965)</td><td>ET<sub>g/b</sub>–FAO56<br>[36]</td></tr><tr><td>Yield</td><td>SPOT 5<br>NDVI</td><td>-Training site<br>-Yield (82 plots)</td><td>-LAI (regression)<br>-Classification<br>(algorithm decision tree)</td><td>LAI, land use, yield<br>[29]</td></tr></tbody></table>

## 2.3. Water Footprint Background

The WF quantifies the amount of water needed to produce what we consume. It is essential for developing both multi-scale and holistic approaches for suitable water resource management. To determine crop WF, ref. [10] separates green and blue components to determine the ratio between crop water use (CWU) and the yield (Y) (Equations (1) and (2)). Green and blue evapotranspiration (ET<sub>g</sub> and ET<sub>b</sub>, respectively) are defined as the total water consumption of crops over the growing period from planting to harvest. The ET<sub>g</sub> represents green water evapotranspiration, calculated as the minimum of total cumulative actual evapotranspiration ET<sub>a</sub> and the cumulative effective rainfall R<sub>e</sub> during the growing season of the crops (Equation (3)). The R<sub>e</sub> refers to the part of the total precipitation that is retained by the soil so that it is potentially available for meeting the water needs of the crop [10]. The ET<sub>b</sub> represents blue water evapotranspiration; it is equal to the difference between the total crop evapotranspiration and the effective rainfall R<sub>e</sub>. When the effective rainfall exceeds the crop evapotranspiration, ET<sub>b</sub> is equal to zero (Equation (4)).

$$WF_g = \frac{CWU_g}{Y} = \frac{10\sum_{d=1}^{lgp} ET_g}{Y} \qquad (1)$$

$$WF_b = \frac{CWU_b}{Y} = \frac{10\sum_{d=1}^{lgp} ET_b}{Y} \qquad (2)$$

WF: water footprint (m³/kg).

CWU: crop water use (m³/ha).

Y: yield (kg).

ET: evapotranspiration determined for each period from day 1 to lgp (length of growing period) (mm/day).

$$ET_g = \text{Min}(ET_c, R_e) \qquad (3)$$

$$ET_b = \text{Max}(0, ET_c - R_e) \qquad (4)$$
November 2012 to June 2013. The processing sequence of EC data to calculate turbulent fluxes from raw data is detailed in [33]. Furthermore, field observations and information from farmers regarding irrigated and rainfed crops were determined for 82 plots in the study area between 2010 and 2012. Table 1 provides a synthesis of the datasets and the preprocessing steps required for different models used in water footprint assessment.

<table><caption>**Table 1.** Summary of variables and models required in WF.</caption><thead><tr><th>Model</th><th>Earth Observation</th><th>Ground and Auxiliary Data</th><th>Process</th><th>WF Component</th></tr></thead><tbody><tr><td>S-SEBI</td><td>Landsat 7-ETM+<br>Landsat 8-OLI<br>B2 to B7: 30 m res<br>B10 to B11: 100 m res<br>-LST, Albedo, NDVI, ε<br>(2010–2013)</td><td>—</td><td>Daily:<br>-R<sub>n</sub><br>-Evaporative fraction Λ<br>-Fluxes: H, G<sub>0</sub>, λE</td><td>Yearly cumulative<br>ET<sub>g/b</sub></td></tr><tr><td>FAO-56<br>Single crop coefficient</td><td>—</td><td>-Air temperature<br>-Air humidity<br>-Solar radiation<br>-Rainfall<br>-Pressure</td><td>ET<sub>0</sub><br>(Penman–Monteith, 1965)</td><td>ET<sub>g/b</sub>–FAO56<br>[36]</td></tr><tr><td>Yield</td><td>SPOT 5<br>NDVI</td><td>-Training site<br>-Yield (82 plots)</td><td>-LAI (regression)<br>-Classification<br>(algorithm decision tree)</td><td>LAI, land use, yield<br>[29]</td></tr></tbody></table>

## 2.3. Water Footprint Background

The WF quantifies the amount of water needed to produce what we consume. It is essential for developing both multi-scale and holistic approaches for suitable water resource management. To determine crop WF, ref. [10] separates green and blue components to determine the ratio between crop water use (CWU) and the yield (Y) (Equations (1) and (2)). Green and blue evapotranspiration (ETg and ETb, respectively) are defined as the total water consumption of crops over the growing period from planting to harvest. The ETg represents green water evapotranspiration, calculated as the minimum of total cumulative actual evapotranspiration ETa and the cumulative effective rainfall Re during the growing season of the crops (Equation (3)). The Re refers to the part of the total precipitation that is retained by the soil so that it is potentially available for meeting the water needs of the crop [10]. The ETb represents blue water evapotranspiration; it is equal to the difference between the total crop evapotranspiration and the effective rainfall Re. When the effective rainfall exceeds the crop evapotranspiration, ETb is equal to zero (Equation (4)).

$$WF_g = \frac{CWU_g}{Y} = \frac{10\sum_{d=1}^{lgp} ET_g}{Y} \qquad (1)$$

$$WF_b = \frac{CWU_b}{Y} = \frac{10\sum_{d=1}^{lgp} ET_b}{Y} \qquad (2)$$

WF: water footprint (m³/kg).

CWU: crop water use (m³/ha).

Y: yield (kg).

ET: evapotranspiration determined for each period from day 1 to lgp (length of growing period) (mm/day).

$$ET_g = \text{Min}(ET_c, R_e) \qquad (3)$$

$$ET_b = \text{Max}(0, ET_c - R_e) \qquad (4)$$
The effective rain $R_e$ was calculated by a model of interception proposed by [39] cited in [40] in Equation (5) as a function of the rainfall $R$, the vegetation leaf area index (LAI) and the vegetation fraction cover $F_c$. In the study region for cereals, the LAI has been calculated by an exponential-type relationship from the NDVI, which was obtained on the test plots of the Kairouan plain in [26] (Equation (6)). The vegetation fraction cover $F_c$ was determined in [33] by a linear relation from NDVI [41] derived from SPOT images (Equation (7)).

$$R_e = R - a \times (LAI) \times (1 - 1/(1 + (F_c \times R)/(a \times (LAI)))) \quad (5)$$

R: actual rainfall (mm/day).

a: empirical parameter representing the crop saturation per unit foliage area (~0.28 for most crops).

$$LAI = \left( \frac{-1}{0.77} \right) \times a \log \left( \frac{(0.75 - NDVI)}{0.6} \right) \quad (6)$$

$$F_c = 1.25 \times NDVI - 0.13 \quad (7)$$

We used ETa modeled by S-SEBI integrating earth observation for biophysical variable estimation. The cumulative is the summation of ETa-i, where (i) corresponds to each period separating two satellite image acquisitions throughout the growing season.

## 2.4. WF Components by Earth Observation

### 2.4.1. Evapotranspiration by S-SEBI Model

The S-SEBI model estimates the crop evapotranspiration using the basic energy balance (Equation (8)). The evapotranspiration corresponds to the latent heat flux term $\lambda E$ of the energy balance, as proposed in [18,42]. More precisely, E is evapotranspiration and $\lambda$ is the latent heat of vaporization, corresponding to the energy required to vaporize water [43]. The net radiation given by (Equation (9)) is defined as the result of all incoming and outgoing radiation [18], which can be measured directly or calculated from incoming radiation, and the remotely sensed biophysical variables of albedo, surface temperature, and the surface emissivity considered either constant over the whole study region (equal to 0.98 in [44] or as a function of vegetation indices [45]).

$$R_n = G_0 + H + \lambda E \quad (8)$$

$R_n$: net radiation (W × m⁻²).

$G_0$: soil heat flux (W × m⁻²).

$H$: sensible heat flux (W × m⁻²).

$\lambda E$: latent heat flux (W × m⁻²).

$$R_n = (1 - \alpha) \tau K_{\text{sun}}^{\downarrow} - \sigma \epsilon_0 T_0^4 + L^{\downarrow} \quad (9)$$

$$K_{\text{sun}}^{\downarrow} = S \times \cos Z \quad (10)$$

$$\cos(Z) = \sin(\lambda) \sin(\delta) + \cos(\lambda) \cos(h) \quad (11)$$

$$\delta = 23.45 \sin [360/365(384 + julien)] \quad (12)$$

$\alpha$: surface albedo [-].

$\tau$: the atmospheric transmittance (-) is a function of $\theta s$, where: $\theta s$ is the solar zenith angle described earlier (degree) after [46].

$K_{\text{sun}}^{\downarrow}$: exoatmospherical incoming solar radiation (W × m⁻²).

$S$: is the exa-atmospheric incident solar radiation (solar constant = 1395 (W.m⁻²)).

$\lambda$: latitude (degree).

$\delta$: declination of the Sun.

$h$: hour angle ($h$ = (12-h) × 15) (degree).

$\sigma$: Stefan-Boltzmann constant (5.67 × 10⁻⁸ W × m⁻² × K⁻⁴).

$\epsilon_0$: surface emissivity.
The effective rain $R_e$ was calculated by a model of interception proposed by [39] cited in [40] in Equation (5) as a function of the rainfall $R$, the vegetation leaf area index (LAI) and the vegetation fraction cover $F_c$. In the study region for cereals, the LAI has been calculated by an exponential-type relationship from the NDVI, which was obtained on the test plots of the Kairouan plain in [26] (Equation (6)). The vegetation fraction cover $F_c$ was determined in [33] by a linear relation from NDVI [41] derived from SPOT images (Equation (7)).

$$R_e = R - a \times (LAI) \times (1 - 1/(1 + (F_c \times R)/(a \times (LAI)))) \quad (5)$$

R: actual rainfall (mm/day).

a: empirical parameter representing the crop saturation per unit foliage area (~0.28 for most crops).

$$LAI = \left(\frac{-1}{0.77}\right) \times a \log \left(\frac{0.75 - NDVI}{0.6}\right) \quad (6)$$

$$F_c = 1.25 \times NDVI - 0.13 \quad (7)$$

We used ET$_{a}$ modeled by S-SEBI integrating earth observation for biophysical variable estimation. The cumulative is the summation of ET$_{a-i}$, where (i) corresponds to each period separating two satellite image acquisitions throughout the growing season.

## 2.4. WF Components by Earth Observation

### 2.4.1. Evapotranspiration by S-SEBI Model

The S-SEBI model estimates the crop evapotranspiration using the basic energy balance (Equation (8)). The evapotranspiration corresponds to the latent heat flux term $\lambda E$ of the energy balance, as proposed in [18,42]. More precisely, E is evapotranspiration and $\lambda$ is the latent heat of vaporization, corresponding to the energy required to vaporize water [43]. The net radiation given by (Equation (9)) is defined as the result of all incoming and outgoing radiation [18], which can be measured directly or calculated from incoming radiation, and the remotely sensed biophysical variables of albedo, surface temperature, and the surface emissivity considered either constant over the whole study region (equal to 0.98 in [44] or as a function of vegetation indices [45]).

$$R_n = G_0 + H + \lambda E \quad (8)$$

$R_n$: net radiation (W × m⁻²).

$G_0$: soil heat flux (W × m⁻²).

$H$: sensible heat flux (W × m⁻²).

$\lambda E$: latent heat flux (W × m⁻²).

$$R_n = (1 - \alpha) \tau K_{\text{sun}}^{\downarrow} - \sigma \epsilon_0 T_0^4 + L^{\downarrow} \quad (9)$$

$$K_{\text{sun}}^{\downarrow} = S \times \cos Z \quad (10)$$

$$\cos(Z) = \sin(\lambda) \sin(\delta) + \cos(\lambda) \cos(h) \quad (11)$$

$$\delta = 23.45 \sin [360/365(384 + julien)] \quad (12)$$

$\alpha$: surface albedo [-].

$\tau$: the atmospheric transmittance (-) is a function of $\theta s$, where: $\theta s$ is the solar zenith angle described earlier (degree) after [46].

$K_{\text{sun}}^{\downarrow}$: exoatmospherical incoming solar radiation (W × m⁻²).

$S$: is the exa-atmospheric incident solar radiation (solar constant = 1395 (W.m⁻²)).

$\lambda$: latitude (degree).

$\delta$: declination of the Sun.

$h$: hour angle ($h$ = (12-h) × 15) (degree).

$\sigma$: Stefan-Boltzmann constant (5.67 × 10⁻⁸ W × m⁻² × K⁻⁴).

$\epsilon_0$: surface emissivity.
$T_0$: surface temperature (°K).

$L^↓$: long-wave incoming radiation. It is measured as per [18] and is equal to 350.75 (W × m⁻²).

As mentioned in [47], the soil flux $G_0$ is generally considered a residual term of the energy balance equation or assumed to be negligible on daily timescales [48]. It varies with vegetation cover, and can range from 0.10 for full cover [49] to a value between 0.10 and 0.50 over sparse canopies and bare soils [50,51]. Amongst the various empirical formulations of $G_0$, we chose the formula proposed by [52] in semi-arid context and found that it allows a good fit of modeled ETa versus measured ones compared to the $G_0$ proposed in [18].

$$G_0 = R_n \times 0.5 \times e^{(-2.13 \times \text{MSAVI})} \qquad (13)$$

$$MSAVI = 0.5(2 \times \rho_R + 1) - 0.5 \times \sqrt{(2 \times \rho_R + 1)^2 - 8(\rho_R - \rho_{NIR})} \qquad (14)$$

$\rho_{NIR}$: reflectance (planetary albedo) of 650–670 nm band (red) (-).

$\rho_R$: reflectance (planetary albedo) of 841–876 nm band (near infrared) (-).

The S-SEBI model is developed for the estimation of turbulent atmospheric fluxes and evaporation fraction using remote sensing data. It has been used in other studies, [18,21,47] in a semi-arid and heterogeneous context. The model requires the estimation of biophysical variables by remote sensing products, namely, surface albedo, emissivity, and surface temperature. S-SEBI determines the sensible heat flux H based on the assumption that the evaporative fraction $\Lambda$ equals the ratio between the latent heat flux and the summation of the sensible and latent heat flux. $\Lambda$ is determined from the scatter plot of the albedo and the surface temperature for a given area where different soil moisture and water availability conditions prevail, as in mixed irrigated and rainfed agriculture region. Dry and wet edges are determined from the scatterplot of each image, where the water status of each pixel is determined relative to these edges [18] (Figure 4). These two limits, corresponding to the maximum sensible heat flux $H_{max}$ ($\alpha$) for the dry edge and to the maximum latent heat flux $\lambda E_{max}$ ($\alpha$) for the wet edge, allow a simple way to compute the evaporative fraction $\Lambda$ (Equation (15)). Finally, $\lambda E$ is determined as $(R_n - G_0)$ multiplied by the evaporative fraction ($\Lambda$).

$$\Lambda = \frac{\lambda E}{\lambda E + H} = \frac{a_{\lambda E} + b_{\lambda E} \times \alpha}{a_H + a_{\lambda E} + (b_H + b_{\lambda E}) \times \alpha} \qquad (15)$$

<figure><img src="image_7.png" /><figcaption>**Figure 4.** Example of a scatter with the dry edge (red) and the wet edge (green) from albedo and $T_s$ of LANDSAT images.</figcaption></figure>
$T_0$: surface temperature (°K).

$L^↓$: long-wave incoming radiation. It is measured as per [18] and is equal to 350.75 (W × m⁻²).

As mentioned in [47], the soil flux $G_0$ is generally considered a residual term of the energy balance equation or assumed to be negligible on daily timescales [48]. It varies with vegetation cover, and can range from 0.10 for full cover [49] to a value between 0.10 and 0.50 over sparse canopies and bare soils [50,51]. Amongst the various empirical formulations of $G_0$, we chose the formula proposed by [52] in semi-arid context and found that it allows a good fit of modeled ETa versus measured ones compared to the $G_0$ proposed in [18].

$$G_0 = R_n \times 0.5 \times e^{(-2.13 \times \text{MSAVI})} \quad (13)$$

$$MSAVI = 0.5(2 \times \rho_R + 1) - 0.5 \times \sqrt{(2 \times \rho_R + 1)^2 - 8(\rho_R - \rho_{NIR})} \quad (14)$$

$\rho_{NIR}$: reflectance (planetary albedo) of 650–670 nm band (red) (-).

$\rho_R$: reflectance (planetary albedo) of 841–876 nm band (near infrared) (-).

The S-SEBI model is developed for the estimation of turbulent atmospheric fluxes and evaporation fraction using remote sensing data. It has been used in other studies, [18,21,47] in a semi-arid and heterogeneous context. The model requires the estimation of biophysical variables by remote sensing products, namely, surface albedo, emissivity, and surface temperature. S-SEBI determines the sensible heat flux H based on the assumption that the evaporative fraction $\Lambda$ equals the ratio between the latent heat flux and the summation of the sensible and latent heat flux. $\Lambda$ is determined from the scatter plot of the albedo and the surface temperature for a given area where different soil moisture and water availability conditions prevail, as in mixed irrigated and rainfed agriculture region. Dry and wet edges are determined from the scatterplot of each image, where the water status of each pixel is determined relative to these edges [18] (Figure 4). These two limits, corresponding to the maximum sensible heat flux $H_{max}$ ($\alpha$) for the dry edge and to the maximum latent heat flux $\lambda E_{max}$ ($\alpha$) for the wet edge, allow a simple way to compute the evaporative fraction $\Lambda$ (Equation (15)). Finally, $\lambda E$ is determined as $(R_n - G_0)$ multiplied by the evaporative fraction ($\Lambda$).

$$\Lambda = \frac{\lambda E}{\lambda E + H} = \frac{a_{\lambda E} + b_{\lambda E} \times \alpha}{a_H + a_{\lambda E} + (b_H + b_{\lambda E}) \times \alpha} \quad (15)$$

<figure><img src="image_9.png" /><figcaption>Figure 4. Example of a scatter with the dry edge (red) and the wet edge (green) from albedo and $T_s$ of LANDSAT images.</figcaption></figure>
In the S-SEBI algorithm, for a given pixel with (a, T_s), the evaporative fraction Λ is based on the maximum and minimum surface temperature scatter of all pixels in the study area, corresponding to the dry edge with a maximum sensible heat flux (H_max), and to the wet edge with a maximum latent flux (λE_max), respectively [47].

The sensible heat flux and the latent heat flux are given by the Equations (16) and (17), respectively, according to [18].

$$H = (1 - \Lambda)(R_n - G_0) \qquad (16)$$

$$\lambda E = \Lambda \times (R_n - G_0) \qquad (17)$$

Actual evapotranspiration ET_a is therefore determined by the S-SEBI model (ET_a-S-SEBI) for each downloaded image based on the latent heat flux λE. This allows an empirical evaluation of the crop coefficient (K_c), a function of ET_a-S-SEBI and ET_0 (Equation (18)). ET_0 is the reference evapotranspiration calculated from the climatic date. Then, to obtain the cumulative ET_a-cum over the whole season, we used the calculated ET_a for two image acquisition dates, and K_c is interpolated for this period (Equation (19)).

$$K_c = ET_{a-S-SEBI} / ET_0 \qquad (18)$$

$$ET_{a-cum} = \sum (K_{c-int} \times ET_0) \qquad (19)$$

K_{c-int}: interpolated K_c between two image acquisition dates.

### 2.4.2. Earth Observation-Based Yield Evaluation

The cereal crop yield maps used in this study reveal the methodology developed by [26]. They are determined using SPOT-HRG optical data at high spatial resolution. The method is based on a statistical analysis between the most frequently used vegetation index, NDVI, and yields measured on test plots for two agricultural years (2010–2011 and 2011–2012) on the plain of Kairouan. These measurements were divided into two sets: 50% for model training and 50% for validation [26]. A relationship was established between the NDVI and grain and straw yields (irrigated and non-irrigated cereals) for the various dates on which satellite images were available. Grain yields can be estimated from mid-March onwards, when cereals reach the ripening phase. Prior to this date, maximum canopy development depends on sowing date, type of irrigation, fertilization management, etc. This can create temporary variability unrelated to grain filling. After mid-April, the NDVI begins to fall, which can also lead to errors in yield estimation.

An estimate of cereal grain yields is possible from the NDVI for mid-March using Equation (17). A good correlation coefficient (R² = 0.65) for the average grain weight is obtained. The validation of the estimates obtained by remote sensing shows that this approach is robust, with an error equal to 8.5 qx/ha for grain yields. To estimate water footprint for the cereal crop and for the three agricultural years, the same proposed empirical relationships was applied for 2012–2013. Grain yield maps are proposed for the class of cereals that includes wheat and barley at the site studied and for all agricultural years using the Equation (20).

$$Y_{\text{grain}} = 4.34 \exp (3.38 \times \text{NDVI}) \qquad (20)$$

## 3. Results

### 3.1. ET_a by S-SEBI and FAO56 Models

The validation of the S-SEBI estimate of ET_a was conducted by comparing the model's instantaneous fluxes with the measured fluxes. This validation was limited to the images acquired for the year 2012–2013 due to the absence of measured data for other years. Table 2 presents a statistical comparison between the heat flux, evaporative fraction and ET_a measured and estimated by S-SEBI. The results indicated a satisfactory agreement with the field measurement data for the primary surface energy balance fluxes R_n and λE, and for the daily ET_a values, with RMSE of about 110.96 W/m², 33.35 W/m² and 0.45 mm/day, respectively, and correlation coefficients of R² = 0.54, 0.85 and 0.82, respectively.
In the S-SEBI algorithm, for a given pixel with (a, T_s), the evaporative fraction Λ is based on the maximum and minimum surface temperature scatter of all pixels in the study area, corresponding to the dry edge with a maximum sensible heat flux (H_max), and to the wet edge with a maximum latent flux (λE_max), respectively [47].

The sensible heat flux and the latent heat flux are given by the Equations (16) and (17), respectively, according to [18].

$$H = (1 - \Lambda)(R_n - G_0) \qquad (16)$$

$$\lambda E = \Lambda \times (R_n - G_0) \qquad (17)$$

Actual evapotranspiration ET_a is therefore determined by the S-SEBI model (ET_a-S-SEBI) for each downloaded image based on the latent heat flux λE. This allows an empirical evaluation of the crop coefficient (K_c), a function of ET_a-S-SEBI and ET_0 (Equation (18)). ET_0 is the reference evapotranspiration calculated from the climatic date. Then, to obtain the cumulative ET_a-cum over the whole season, we used the calculated ET_a for two image acquisition dates, and K_c is interpolated for this period (Equation (19)).

$$K_c = ET_{a-S-SEBI} / ET_0 \qquad (18)$$

$$ET_{a-cum} = \sum (K_{c-int} \times ET_0) \qquad (19)$$

K_{c-int}: interpolated K_c between two image acquisition dates.

### 2.4.2. Earth Observation-Based Yield Evaluation

The cereal crop yield maps used in this study reveal the methodology developed by [26]. They are determined using SPOT-HRG optical data at high spatial resolution. The method is based on a statistical analysis between the most frequently used vegetation index, NDVI, and yields measured on test plots for two agricultural years (2010–2011 and 2011–2012) on the plain of Kairouan. These measurements were divided into two sets: 50% for model training and 50% for validation [26]. A relationship was established between the NDVI and grain and straw yields (irrigated and non-irrigated cereals) for the various dates on which satellite images were available. Grain yields can be estimated from mid-March onwards, when cereals reach the ripening phase. Prior to this date, maximum canopy development depends on sowing date, type of irrigation, fertilization management, etc. This can create temporary variability unrelated to grain filling. After mid-April, the NDVI begins to fall, which can also lead to errors in yield estimation.

An estimate of cereal grain yields is possible from the NDVI for mid-March using Equation (17). A good correlation coefficient (R² = 0.65) for the average grain weight is obtained. The validation of the estimates obtained by remote sensing shows that this approach is robust, with an error equal to 8.5 qx/ha for grain yields. To estimate water footprint for the cereal crop and for the three agricultural years, the same proposed empirical relationships was applied for 2012–2013. Grain yield maps are proposed for the class of cereals that includes wheat and barley at the site studied and for all agricultural years using the Equation (20).

$$Y_{\text{grain}} = 4.34 \exp (3.38 \times \text{NDVI}) \qquad (20)$$

## 3. Results

### 3.1. ET_a by S-SEBI and FAO56 Models

The validation of the S-SEBI estimate of ET_a was conducted by comparing the model's instantaneous fluxes with the measured fluxes. This validation was limited to the images acquired for the year 2012–2013 due to the absence of measured data for other years. Table 2 presents a statistical comparison between the heat flux, evaporative fraction and ET_a measured and estimated by S-SEBI. The results indicated a satisfactory agreement with the field measurement data for the primary surface energy balance fluxes R_n and λE, and for the daily ET_a values, with RMSE of about 110.96 W/m², 33.35 W/m² and 0.45 mm/day, respectively, and correlation coefficients of R² = 0.54, 0.85 and 0.82, respectively.
<table><caption>**Table 2.** Comparison between S-SEBI model estimated and measured instantaneous fluxes, evaporative fraction and ET<sub>a</sub> for 2012–2013 (RMSE: root mean squared error, R<sup>2</sup>: correlation coefficient).</caption><thead><tr><th></th><th>R<sub>n</sub> (W × m<sup>-2</sup>)</th><th>G<sub>0</sub> (W × m<sup>-2</sup>)</th><th>H (W × m<sup>-2</sup>)</th><th>λE (W × m<sup>-2</sup>)</th><th>Λ (-)</th><th>ET<sub>a</sub> (mm/Day)</th></tr></thead><tbody><tr><td>RMSE</td><td>110.96</td><td>68.28</td><td>100.87</td><td>33.35</td><td>0.18</td><td>0.45</td></tr><tr><td>R<sup>2</sup></td><td>0.54</td><td>0.09</td><td>0.05</td><td>0.85</td><td>0.51</td><td>0.82</td></tr></tbody></table>

Figure 5 shows the ET<sub>a</sub> average simulated by S-SEBI, measured by the flux tower and calculated using the FAO-56 single crop coefficient method. Both the FAO-56 and S-SEBI estimates of ET<sub>a</sub> capture the dynamic nature of the growing stages and their impact on the ET process. However, S-SEBI values demonstrate a closer alignment with the measurements compared to FAO-56. Figure 6 exhibits the yearly spatial distribution of modeled ET<sub>a</sub> by S-SEBI from 2010 to 2013, showcasing values spanning 0 to 577 mm/year. Regions characterized by heightened ET<sub>a</sub> are emphasized in red, indicating alignment with irrigated zones. These visual representations serve as useful tools for water managers, enabling the effective monitoring of irrigated areas and supporting informed decisions in water resource management.

<figure><img src="image_8.png" /><figcaption>**Figure 5.** Comparison of ET<sub>a</sub> of cereals by FAO-56 and S-SEBI based on Landsat 7 (2012–2013).</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Figure 6.** *Cont.*</figcaption></figure>
<table><caption>**Table 2.** Comparison between S-SEBI model estimated and measured instantaneous fluxes, evaporative fraction and ET<sub>a</sub> for 2012–2013 (RMSE: root mean squared error, R<sup>2</sup>: correlation coefficient).</caption><thead><tr><th></th><th>R<sub>n</sub> (W × m<sup>-2</sup>)</th><th>G<sub>0</sub> (W × m<sup>-2</sup>)</th><th>H (W × m<sup>-2</sup>)</th><th>λE (W × m<sup>-2</sup>)</th><th>Λ (-)</th><th>ET<sub>a</sub> (mm/Day)</th></tr></thead><tbody><tr><td>RMSE</td><td>110.96</td><td>68.28</td><td>100.87</td><td>33.35</td><td>0.18</td><td>0.45</td></tr><tr><td>R<sup>2</sup></td><td>0.54</td><td>0.09</td><td>0.05</td><td>0.85</td><td>0.51</td><td>0.82</td></tr></tbody></table>

Figure 5 shows the ET<sub>a</sub> average simulated by S-SEBI, measured by the flux tower and calculated using the FAO-56 single crop coefficient method. Both the FAO-56 and S-SEBI estimates of ET<sub>a</sub> capture the dynamic nature of the growing stages and their impact on the ET process. However, S-SEBI values demonstrate a closer alignment with the measurements compared to FAO-56. Figure 6 exhibits the yearly spatial distribution of modeled ET<sub>a</sub> by S-SEBI from 2010 to 2013, showcasing values spanning 0 to 577 mm/year. Regions characterized by heightened ET<sub>a</sub> are emphasized in red, indicating alignment with irrigated zones. These visual representations serve as useful tools for water managers, enabling the effective monitoring of irrigated areas and supporting informed decisions in water resource management.

<figure><img src="image_10.png" /><figcaption>**Figure 5.** Comparison of ET<sub>a</sub> of cereals by FAO-56 and S-SEBI based on Landsat 7 (2012–2013).</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>**Figure 6.** *Cont.*</figcaption></figure>
<figure><img src="image_10.png" /><figcaption>**Figure 6.** ET<sub>a</sub> of cereals by S-SEBI model for the 3 agricultural years 2010–2011, 2011–2012, and 2012–2013 on the plain of Kairouan.</figcaption></figure>

3.2. Yield Mapping

Equation (20) was utilized to invert the NDVI and generate maps of grain cereal
yield for three consecutive crop years. The spatial modeling of the model covered the
entire Kairouan plain, utilizing NDVI images acquired on 17 March 2011 for the 2010–2011
agricultural year, on 31 March 2012 for the 2011–2012 agricultural year, and on 19 March
2013 for the 2012–2013 agricultural year. Cereal masks derived from land-cover maps
established for each of the three agricultural years were applied to the NDVI images. This
ensured the inversion process exclusively considered cereal plots across the entire study
area (see Figure 7). It is important to note that these maps do not incorporate growth factors
such as irrigation, fertilization, crop rotation, or cereal varieties (e.g., Karim, Razzak, Khiar,
and Om Rabii). The approach relies solely on the NDVI index during the grain maturation
phase, occurring between mid-March and mid-April.
<figure><img src="image_12.png" /><figcaption>**Figure 6.** ET<sub>a</sub> of cereals by S-SEBI model for the 3 agricultural years 2010–2011, 2011–2012, and 2012–2013 on the plain of Kairouan.</figcaption></figure>

## 3.2. Yield Mapping

Equation (20) was utilized to invert the NDVI and generate maps of grain cereal yield for three consecutive crop years. The spatial modeling of the model covered the entire Kairouan plain, utilizing NDVI images acquired on 17 March 2011 for the 2010–2011 agricultural year, on 31 March 2012 for the 2011–2012 agricultural year, and on 19 March 2013 for the 2012–2013 agricultural year. Cereal masks derived from land-cover maps established for each of the three agricultural years were applied to the NDVI images. This ensured the inversion process exclusively considered cereal plots across the entire study area (see Figure 7). It is important to note that these maps do not incorporate growth factors such as irrigation, fertilization, crop rotation, or cereal varieties (e.g., Karim, Razzak, Khiar, and Om Rabii). The approach relies solely on the NDVI index during the grain maturation phase, occurring between mid-March and mid-April.
<figure><img src="image_11.png" /><figcaption>**Figure 7.** Yield map for the 3 agricultural years 2010–2011, 2011–2012, and 2012–2013 on the plain of Kairouan.</figcaption></figure>

### 3.3. Water Footprint Mapping

#### 3.3.1. Cereal WF Global Evaluation

The WF of cereals was determined for three agricultural years (2010–2013) based on water consumption assessed by ETₐ and yields, both derived from remote sensing data. The
<figure><img src="image_13.png" /><figcaption>**Figure 7.** Yield map for the 3 agricultural years 2010–2011, 2011–2012, and 2012–2013 on the plain of Kairouan.</figcaption></figure>

### 3.3. Water Footprint Mapping

#### 3.3.1. Cereal WF Global Evaluation

The WF of cereals was determined for three agricultural years (2010–2013) based on water consumption assessed by ETₐ and yields, both derived from remote sensing data. The
comparison between rainfall, cultivated area, production, volume of green water, volume of blue water, green WF, and blue WF is presented in Table 3. The table reveals that across a 444 km² study area, cereal crops occupied 16%, 35%, and 6% of the total area in 2011, 2012, and 2013, respectively. This variation in the proportion of cultivated area primarily stems from differences in water availability, whether from rainfall or irrigation, during the growing season. These disparities in water availability inevitably impact yields, volume of green/blue water, and consequently the water footprint (WF).

<table><caption>**Table 3.** Water footprint of cereals in Kairouan plain based on remote sensing modeling (2010–2013).</caption><thead><tr><th>Year</th><th>Rain (mm/year)</th><th>Study Area (km²)</th><th>Cereal Area (km²)</th><th>% Area Cereals</th><th>Production (Tons)</th><th>Green Water (10<sup>6</sup> m³/year)</th><th>Blue Water (10<sup>6</sup> m³/year)</th><th>Average Green WF (m³/kg)</th><th>Average Blue WF (m³/kg)</th><th>Total Average WF (m³/kg)</th></tr></thead><tbody><tr><td>2010–2011</td><td>328.8</td><td>444</td><td>72.48</td><td>16</td><td>20,996</td><td>7.6</td><td>4.60</td><td>0.67 (δ = 1.01)</td><td>0.41 (δ = 1.09)</td><td>1.08</td></tr><tr><td>2011–2012</td><td>344</td><td>444</td><td>156.68</td><td>35</td><td>40,304</td><td>13.82</td><td>18.41</td><td>0.47 (δ = 0.63)</td><td>0.61 (δ = 1.30)</td><td>1.08</td></tr><tr><td>2012–2013</td><td>177.2</td><td>444</td><td>25.48</td><td>6</td><td>8054</td><td>1.02</td><td>6.60</td><td>0.17 (δ = 0.38)</td><td>1.04 (δ = 1.72)</td><td>1.22</td></tr></tbody></table>

Building upon the references cited in the introduction regarding global and local methods for evaluating cereals’ water footprint worldwide and in Tunisia [6,8,9,11,12], we established threshold values for water footprint (WF) as follows: a low WF, indicating high water productivity, corresponds to values less than 0.66 m³/kg, a medium WF falls between 0.66 and 1 m³/kg, and a high WF exceeds 1 m³/kg.

Table 3 shows that the total average WF was the same (1.08 m³/kg) for the 2010–2011 and 2011–2012 growing seasons, with very close amounts of total rainfall (328.8 mm and 344 mm, respectively). However, the highest total WF was found for the driest growing season of 2012–2013 (total rainfall of about 177 mm), with total WF exceeding 1 m³/kg indicating low water productivity, which can be explained by the multiple effect of the cereal yield drop and actual evapotranspiration increase during the 2012–2013 growing season. These findings align with results for irrigated cereals, for which WF values were generally below 1 m³/kg in the northern region of Tunisia [53] and approximately 1 m³/kg for all types of cereals [54], which is in agreement with our results, since we are dealing with cereals cropped in central Tunisia. We found WF values of about 1 m³/kg for the growing seasons with usual rainfall amounts.

### 3.3.2. Comparing and Assessing Spatiotemporal Variations in WF Cereals

To assess the green and blue WF of cereals along with their spatial distribution, we computed a water footprint index (WFI) using Equation (21). A positive value of the WFI signifies that the WF is predominantly attributed to the green component, while a negative WFI value indicates a prevalence of blue water in the footprint.

$$WFI = \frac{WF_g - WF_b}{WF_g + WF_b} \qquad (21)$$

Figure 8 depicts the WFI of cereals for the three agricultural seasons. In the seasons 2010–2011 and 2011–2012, the WFI predominantly showcases the influence of the green water component in the first season, while the second season is characterized by the prevalence of the blue water component. Conversely, during the 2012–2013 season, there is an almost complete absence of the green water component in cereal production within the studied area.

Taking into account the consensus that attributes blue water to irrigation and green water to the consumption of rainfed crops, we assessed the WFI by comparing it to field data on rainfed and irrigated plots between 2010 and 2012. As detailed in Table 4, 69% of rainfed plots exhibited a positive WFI, signifying a prevalence of green water, while only 52% of irrigated plots had a negative WFI, indicating a prevalence of blue water.
comparison between rainfall, cultivated area, production, volume of green water, volume of blue water, green WF, and blue WF is presented in Table 3. The table reveals that across a 444 km² study area, cereal crops occupied 16%, 35%, and 6% of the total area in 2011, 2012, and 2013, respectively. This variation in the proportion of cultivated area primarily stems from differences in water availability, whether from rainfall or irrigation, during the growing season. These disparities in water availability inevitably impact yields, volume of green/blue water, and consequently the water footprint (WF).

<table><caption>**Table 3.** Water footprint of cereals in Kairouan plain based on remote sensing modeling (2010–2013).</caption><thead><tr><th>Year</th><th>Rain (mm/year)</th><th>Study Area (km²)</th><th>Cereal Area (km²)</th><th>% Area Cereals</th><th>Production (Tons)</th><th>Green Water (10<sup>6</sup> m³/year)</th><th>Blue Water (10<sup>6</sup> m³/year)</th><th>Average Green WF (m³/kg)</th><th>Average Blue WF (m³/kg)</th><th>Total Average WF (m³/kg)</th></tr></thead><tbody><tr><td>2010–2011</td><td>328.8</td><td>444</td><td>72.48</td><td>16</td><td>20,996</td><td>7.6</td><td>4.60</td><td>0.67 (δ = 1.01)</td><td>0.41 (δ = 1.09)</td><td>1.08</td></tr><tr><td>2011–2012</td><td>344</td><td>444</td><td>156.68</td><td>35</td><td>40,304</td><td>13.82</td><td>18.41</td><td>0.47 (δ = 0.63)</td><td>0.61 (δ = 1.30)</td><td>1.08</td></tr><tr><td>2012–2013</td><td>177.2</td><td>444</td><td>25.48</td><td>6</td><td>8054</td><td>1.02</td><td>6.60</td><td>0.17 (δ = 0.38)</td><td>1.04 (δ = 1.72)</td><td>1.22</td></tr></tbody></table>

Building upon the references cited in the introduction regarding global and local methods for evaluating cereals’ water footprint worldwide and in Tunisia [6,8,9,11,12], we established threshold values for water footprint (WF) as follows: a low WF, indicating high water productivity, corresponds to values less than 0.66 m³/kg, a medium WF falls between 0.66 and 1 m³/kg, and a high WF exceeds 1 m³/kg.

Table 3 shows that the total average WF was the same (1.08 m³/kg) for the 2010–2011 and 2011–2012 growing seasons, with very close amounts of total rainfall (328.8 mm and 344 mm, respectively). However, the highest total WF was found for the driest growing season of 2012–2013 (total rainfall of about 177 mm), with total WF exceeding 1 m³/kg indicating low water productivity, which can be explained by the multiple effect of the cereal yield drop and actual evapotranspiration increase during the 2012–2013 growing season. These findings align with results for irrigated cereals, for which WF values were generally below 1 m³/kg in the northern region of Tunisia [53] and approximately 1 m³/kg for all types of cereals [54], which is in agreement with our results, since we are dealing with cereals cropped in central Tunisia. We found WF values of about 1 m³/kg for the growing seasons with usual rainfall amounts.

### 3.3.2. Comparing and Assessing Spatiotemporal Variations in WF Cereals

To assess the green and blue WF of cereals along with their spatial distribution, we computed a water footprint index (WFI) using Equation (21). A positive value of the WFI signifies that the WF is predominantly attributed to the green component, while a negative WFI value indicates a prevalence of blue water in the footprint.

$$WFI = \frac{WF_g - WF_b}{WF_g + WF_b} \qquad (21)$$

Figure 8 depicts the WFI of cereals for the three agricultural seasons. In the seasons 2010–2011 and 2011–2012, the WFI predominantly showcases the influence of the green water component in the first season, while the second season is characterized by the prevalence of the blue water component. Conversely, during the 2012–2013 season, there is an almost complete absence of the green water component in cereal production within the studied area.

Taking into account the consensus that attributes blue water to irrigation and green water to the consumption of rainfed crops, we assessed the WFI by comparing it to field data on rainfed and irrigated plots between 2010 and 2012. As detailed in Table 4, 69% of rainfed plots exhibited a positive WFI, signifying a prevalence of green water, while only 52% of irrigated plots had a negative WFI, indicating a prevalence of blue water.
Notably, in the 2010–2011 period, all irrigated plots displayed a positive WFI, contrary to the expected trend.

<figure><img src="image_12.png" /><figcaption>**Figure 8.** WFI of cereals for the three agricultural years 2010–2011, 2011–2012, and 2012–2013 on the plain of Kairouan.</figcaption></figure>
Notably, in the 2010–2011 period, all irrigated plots displayed a positive WFI, contrary to the expected trend.

<figure><img src="image_14.png" /><figcaption>**Figure 8.** WFI of cereals for the three agricultural years 2010–2011, 2011–2012, and 2012–2013 on the plain of Kairouan.</figcaption></figure>
<table><caption>**Table 4.** Comparison of the WF index with plot field data for 2010–2011 and 2011–2012.</caption><thead><tr><th rowspan="2"></th><th colspan="2">2010–2011</th><th colspan="2">2011–2012</th></tr><tr><th>Irrigated Field</th><th>Rainfed Field</th><th>Irrigated Field</th><th>Rainfed Field</th></tr></thead><tbody><tr><td>Field data (number of plots)</td><td>18</td><td>8</td><td>36</td><td>18</td></tr><tr><td>WFI &lt; 0</td><td>0</td><td>0</td><td>28</td><td>8</td></tr><tr><td>Mean WFI</td><td>-</td><td>-</td><td>-0.197</td><td>-0.108</td></tr><tr><td>Standard deviation</td><td>-</td><td>-</td><td>0.126</td><td>0.144</td></tr><tr><td>WFI &gt; 0</td><td>18</td><td>8</td><td>8</td><td>10</td></tr><tr><td>Mean WFI</td><td>0.318</td><td>0.544</td><td>0.118</td><td>0.199</td></tr><tr><td>Standard deviation</td><td>0.071</td><td>0.086</td><td>0.118</td><td>0.317</td></tr></tbody></table>

### 3.3.3. WF-Yield Relationship and Water Management

The scatterplot of the WF plotted against the yields provides a valuable tool for exploring the spatiotemporal relationship between these two components, offering insights into water use efficiency. Similar approaches have been employed in related studies by [55–57] for evaluating spatiotemporal patterns of crop water productivity in various climatic contexts. This method aids in understanding agronomic practices and water efficiency, contributing to improved water management strategies.

Grain yield and WF thresholds were established at 25 qx/ha (the average yield in our study area) following [26] and 1 m³/kg following [54,58]. The scatterplot, with the X-axis representing grain yield and the Y-axis representing WF, was constructed using the mean values of yield and WF in our study area (Figure 9). This approach aims to offer a clear visualization of the distribution of yield and WF concerning these average values. Figure 9 illustrates the scatterplots for the three years of experimentation, revealing four distinct classes ranging from the most water efficient (C1) to the least water efficient (C4). Our objective is to evaluate WF, encompassing water consumption, in terms of yields. While there are methods that compare WF with evapotranspiration to determine water efficiency, a choice not explored in our case, the combination of yields and WF serves as the starting point for further agronomic analysis. This analysis aims to identify the practices and conditions that promote high cereal yields. Notably, it emerges from these analyses that more than 30% of cropped areas during the study years (approximately 40% in 2011 and 2012 and 29% in 2013, as depicted in Figure 10) exhibit low water efficiency with high WF and low yields. We also determined for each class of the combined maps (WF and yield) the distribution of green and blue water based on the WFI, as illustrated in Figure 11. Dark green and light green represent high and low green WF, respectively, while dark blue and light blue represent high and low blue WF, respectively. From this analysis, it becomes apparent that 40% of the cereal plots in 2012, exhibiting low production while consuming a substantial amount of water, are identified as irrigated according to the WFI mapping (81% blue water compared to 6% in 2011). Therefore, managers should prioritize these plots for interventions aimed at improving water efficiency and agronomic practices. To facilitate this, a spatial distribution of plots with low yields but high water consumption was created (Figure 12).

<figure><img src="image_13.png" /><figcaption>**Figure 9.** Scatters of yields and total WF of cereals in the study area.</figcaption></figure>
<table><caption>**Table 4.** Comparison of the WF index with plot field data for 2010–2011 and 2011–2012.</caption><thead><tr><th rowspan="2"></th><th colspan="2">2010–2011</th><th colspan="2">2011–2012</th></tr><tr><th>Irrigated Field</th><th>Rainfed Field</th><th>Irrigated Field</th><th>Rainfed Field</th></tr></thead><tbody><tr><td>Field data (number of plots)</td><td>18</td><td>8</td><td>36</td><td>18</td></tr><tr><td>WFI &lt; 0</td><td>0</td><td>0</td><td>28</td><td>8</td></tr><tr><td>Mean WFI</td><td>-</td><td>-</td><td>-0.197</td><td>-0.108</td></tr><tr><td>Standard deviation</td><td>-</td><td>-</td><td>0.126</td><td>0.144</td></tr><tr><td>WFI &gt; 0</td><td>18</td><td>8</td><td>8</td><td>10</td></tr><tr><td>Mean WFI</td><td>0.318</td><td>0.544</td><td>0.118</td><td>0.199</td></tr><tr><td>Standard deviation</td><td>0.071</td><td>0.086</td><td>0.118</td><td>0.317</td></tr></tbody></table>

### 3.3.3. WF-Yield Relationship and Water Management

The scatterplot of the WF plotted against the yields provides a valuable tool for exploring the spatiotemporal relationship between these two components, offering insights into water use efficiency. Similar approaches have been employed in related studies by [55–57] for evaluating spatiotemporal patterns of crop water productivity in various climatic contexts. This method aids in understanding agronomic practices and water efficiency, contributing to improved water management strategies.

Grain yield and WF thresholds were established at 25 qx/ha (the average yield in our study area) following [26] and 1 m³/kg following [54,58]. The scatterplot, with the X-axis representing grain yield and the Y-axis representing WF, was constructed using the mean values of yield and WF in our study area (Figure 9). This approach aims to offer a clear visualization of the distribution of yield and WF concerning these average values. Figure 9 illustrates the scatterplots for the three years of experimentation, revealing four distinct classes ranging from the most water efficient (C1) to the least water efficient (C4). Our objective is to evaluate WF, encompassing water consumption, in terms of yields. While there are methods that compare WF with evapotranspiration to determine water efficiency, a choice not explored in our case, the combination of yields and WF serves as the starting point for further agronomic analysis. This analysis aims to identify the practices and conditions that promote high cereal yields. Notably, it emerges from these analyses that more than 30% of cropped areas during the study years (approximately 40% in 2011 and 2012 and 29% in 2013, as depicted in Figure 10) exhibit low water efficiency with high WF and low yields. We also determined for each class of the combined maps (WF and yield) the distribution of green and blue water based on the WFI, as illustrated in Figure 11. Dark green and light green represent high and low green WF, respectively, while dark blue and light blue represent high and low blue WF, respectively. From this analysis, it becomes apparent that 40% of the cereal plots in 2012, exhibiting low production while consuming a substantial amount of water, are identified as irrigated according to the WFI mapping (81% blue water compared to 6% in 2011). Therefore, managers should prioritize these plots for interventions aimed at improving water efficiency and agronomic practices. To facilitate this, a spatial distribution of plots with low yields but high water consumption was created (Figure 12).

<figure><img src="image_15.png" /><figcaption>**Figure 9.** Scatters of yields and total WF of cereals in the study area.</figcaption></figure>
<figure><img src="image_14.png" /></figure>

<figure><img src="image_15.png" /></figure>

<figure><img src="image_16.png" /><figcaption>Figure 10. Combination of yields and total WF.</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>Figure 11. Green and blue WF of each cereal class, dark green and light green indicate high and low green WF, while dark blue and light blue indicate high and low blue WF.</figcaption></figure>

<figure><img src="image_18.png" /><figcaption>Figure 12. *Cont.*</figcaption></figure>
<figure><img src="image_16.png" /><figcaption>Figure 10. Combination of yields and total WF.</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>Figure 11. Green and blue WF of each cereal class, dark green and light green indicate high and low green WF, while dark blue and light blue indicate high and low blue WF.</figcaption></figure>

<figure><img src="image_18.png" /><figcaption>Figure 12. *Cont.*</figcaption></figure>
<figure><img src="image_19.png" /><figcaption>**Figure 12.** Class 4 (C4) presenting low yields and high water consumption for the 3 agricultural years 2010–2011, 2011–2012, and 2012–2013 on the plain of Kairouan.</figcaption></figure>

**4. Discussion**

Comparing the average ETₐ simulated by S-SEBI, measured by the flux tower and calculated by the FAO-56 single crop coefficient method, reveals notable differences. The FAO-56 method tends to overestimate ETₐ, a discrepancy attributed to its assumption of optimal conditions throughout the entire crop growing period, a scenario not always realized. In contrast, the ETₐ estimated by the S-SEBI model, utilizing Landsat 7–8 data, showed good agreement with the data measured by the flux tower station, with an R² = 0.82 and an RMSE = 45 mm/day in Kairouan plain region. The model’s estimates are deemed acceptable, consistent with findings in similar semi-arid contexts [21,47,52]. The three maps presented in Figure 5 highlight noticeable variations in the density of cereals from one year to another, underscoring the impact of drought on agricultural practices in our study area. These observations are consistent with prior literature [26,58]. The differences can be attributed to the predominant irrigation practices in the plain of Kairouan, where a substantial portion of the land is irrigated. In this region, the majority of plots, often small, are not exclusively designated for cereal cultivation. Additionally, crop rotation practices are implemented, especially if the onset of the vegetative cycle (in November or December)
<figure><img src="image_19.png" /><figcaption>**Figure 12.** Class 4 (C4) presenting low yields and high water consumption for the 3 agricultural years 2010–2011, 2011–2012, and 2012–2013 on the plain of Kairouan.</figcaption></figure>

**4. Discussion**

Comparing the average ETₐ simulated by S-SEBI, measured by the flux tower and calculated by the FAO-56 single crop coefficient method, reveals notable differences. The FAO-56 method tends to overestimate ETₐ, a discrepancy attributed to its assumption of optimal conditions throughout the entire crop growing period, a scenario not always realized. In contrast, the ETₐ estimated by the S-SEBI model, utilizing Landsat 7–8 data, showed good agreement with the data measured by the flux tower station, with an R² = 0.82 and an RMSE = 45 mm/day in Kairouan plain region. The model's estimates are deemed acceptable, consistent with findings in similar semi-arid contexts [21,47,52]. The three maps presented in Figure 5 highlight noticeable variations in the density of cereals from one year to another, underscoring the impact of drought on agricultural practices in our study area. These observations are consistent with prior literature [26,58]. The differences can be attributed to the predominant irrigation practices in the plain of Kairouan, where a substantial portion of the land is irrigated. In this region, the majority of plots, often small, are not exclusively designated for cereal cultivation. Additionally, crop rotation practices are implemented, especially if the onset of the vegetative cycle (in November or December)
lacks significant precipitation. In such cases, climatic factors do not favor cereal cultivation in the plain. Given their long growth cycle and the need for substantial irrigation, cereals are not the preferred choice for farmers under these circumstances.

The regression model spatialized the grain cereal yield across the entire Kairouan plain using NDVI images for the three agricultural years 2010–2011, 2011–2012, and 2012–2013. The agricultural year 2011/2012, characterized by heavy rains during the cereal vegetative cycle, witnessed a larger number of fields cultivated with cereals [59,60]. A comparison between the 2010–2011 and 2011–2012 seasons reveals that plots abandoned in the first year became productive in the subsequent year due to the abundant rainfall, even though their yields were relatively low. The increase in cereal production during the 2011–2012 season was also influenced by a government policy aimed at boosting cereal production by reducing the prices of seeds [28]. However, the year 2012–2013 marked the driest period, with the fewest cultivated plots and the lowest yields.

The objective of this article is to estimate the two components of cereal WF (ETₐ and yield) through remote sensing for the three crop years. The overall results of the cereal WF (as shown in Table 3) indicate that in the agricultural year 2010–2011, some plots in the field were abandoned due to delayed rainfall (as depicted in Figure 2), despite the total rainfall being relatively high (328.8 mm). Additionally, regarding irrigation, surveys conducted among farmers suggest that the Ministry of Agriculture permitted a limited release of water, resulting in the abandonment of certain plots. In terms of the impact on the WF, it is observed that the average green water (0.67 m³/kg) exceeds the average blue water (0.41 m³/kg) in cereal production during this year. In the year 2011–2012, which was an average year with 344 mm of rainfall, irrigation was employed, as 202 mm of this rainfall occurred before sowing (20 September) and after the cereal harvest (25 May). This explains why blue water is higher than the average green water (0.61 m³/kg compared to 0.47 m³/kg). The driest year in this study, 2012–2013, witnessed the lowest average green water (0.17 m³/kg) and the highest average blue water (1.04 m³/kg) used in cereal cultivation, resulting in the highest total average WF (1.22 m³/kg), primarily supplied by blue water. The total water footprint (WF) values in our study generally align with the literature, indicating that WF is around 1 m³/kg for various cereal types [53,54]. However, they deviate from the average values reported in [2,61] for the period 1996–2005, showing lower WF values in our study. It is crucial to note that our study did not differentiate between wheat and barley, unlike these aforementioned studies. Given that the majority of plots are cultivated with wheat, we compared our values only to those found for wheat. The average WF for the three years (2010–2013) in our study was 1.13 m³/kg, while [2,61] reported values of 2.5 m³/kg and 4.2 m³/kg, respectively, for the period 1996–2005. Therefore, compared to our study, which is based on modeling ETₐ and yields using remote sensing data, there appears to be an overestimation in the average WF for cereals in the global statistical data. One plausible explanation for this discrepancy is that the estimation of water consumption in global methods heavily relies on the FAO-56 method, which tends to overestimate ET compared to the measured ETₐ, as illustrated in Figure 5.

The water footprint index (WFI) (Equation (21)), was calculated to provide a more accurate assessment of the green and blue waters of cereals and their spatial distribution, with WFI values ranging between -1 (mostly blue WF) and 1 (mostly green WF). The WFI analysis for the two seasons 2010–2011 and 2011–2012 reveals a predominance of the green water component in the first season, while the second season is characterized by the prevalence of blue water. This observation aligns with the reports from farmers regarding low irrigation during the 2010–2011 period. In contrast, during the 2012–2013 season, there is an almost total absence of the green water component in cereal production in the studied area. Therefore, the WFI emerges as a valuable indicator for assessing the interannual variation of cereal WF concerning total precipitation and its occurrence, along with the impact on irrigation water management. The WFI was validated for the years 2010–2011 and 2011–2012, where field data were available. Around 70% of rainfed plots
lacks significant precipitation. In such cases, climatic factors do not favor cereal cultivation in the plain. Given their long growth cycle and the need for substantial irrigation, cereals are not the preferred choice for farmers under these circumstances.

The regression model spatialized the grain cereal yield across the entire Kairouan plain using NDVI images for the three agricultural years 2010–2011, 2011–2012, and 2012–2013. The agricultural year 2011/2012, characterized by heavy rains during the cereal vegetative cycle, witnessed a larger number of fields cultivated with cereals [59,60]. A comparison between the 2010–2011 and 2011–2012 seasons reveals that plots abandoned in the first year became productive in the subsequent year due to the abundant rainfall, even though their yields were relatively low. The increase in cereal production during the 2011–2012 season was also influenced by a government policy aimed at boosting cereal production by reducing the prices of seeds [28]. However, the year 2012–2013 marked the driest period, with the fewest cultivated plots and the lowest yields.

The objective of this article is to estimate the two components of cereal WF (ETₐ and yield) through remote sensing for the three crop years. The overall results of the cereal WF (as shown in Table 3) indicate that in the agricultural year 2010–2011, some plots in the field were abandoned due to delayed rainfall (as depicted in Figure 2), despite the total rainfall being relatively high (328.8 mm). Additionally, regarding irrigation, surveys conducted among farmers suggest that the Ministry of Agriculture permitted a limited release of water, resulting in the abandonment of certain plots. In terms of the impact on the WF, it is observed that the average green water (0.67 m³/kg) exceeds the average blue water (0.41 m³/kg) in cereal production during this year. In the year 2011–2012, which was an average year with 344 mm of rainfall, irrigation was employed, as 202 mm of this rainfall occurred before sowing (20 September) and after the cereal harvest (25 May). This explains why blue water is higher than the average green water (0.61 m³/kg compared to 0.47 m³/kg). The driest year in this study, 2012–2013, witnessed the lowest average green water (0.17 m³/kg) and the highest average blue water (1.04 m³/kg) used in cereal cultivation, resulting in the highest total average WF (1.22 m³/kg), primarily supplied by blue water. The total water footprint (WF) values in our study generally align with the literature, indicating that WF is around 1 m³/kg for various cereal types [53,54]. However, they deviate from the average values reported in [2,61] for the period 1996–2005, showing lower WF values in our study. It is crucial to note that our study did not differentiate between wheat and barley, unlike these aforementioned studies. Given that the majority of plots are cultivated with wheat, we compared our values only to those found for wheat. The average WF for the three years (2010–2013) in our study was 1.13 m³/kg, while [2,61] reported values of 2.5 m³/kg and 4.2 m³/kg, respectively, for the period 1996–2005. Therefore, compared to our study, which is based on modeling ETₐ and yields using remote sensing data, there appears to be an overestimation in the average WF for cereals in the global statistical data. One plausible explanation for this discrepancy is that the estimation of water consumption in global methods heavily relies on the FAO-56 method, which tends to overestimate ET compared to the measured ETₐ, as illustrated in Figure 5.

The water footprint index (WFI) (Equation (21)), was calculated to provide a more accurate assessment of the green and blue waters of cereals and their spatial distribution, with WFI values ranging between -1 (mostly blue WF) and 1 (mostly green WF). The WFI analysis for the two seasons 2010–2011 and 2011–2012 reveals a predominance of the green water component in the first season, while the second season is characterized by the prevalence of blue water. This observation aligns with the reports from farmers regarding low irrigation during the 2010–2011 period. In contrast, during the 2012–2013 season, there is an almost total absence of the green water component in cereal production in the studied area. Therefore, the WFI emerges as a valuable indicator for assessing the interannual variation of cereal WF concerning total precipitation and its occurrence, along with the impact on irrigation water management. The WFI was validated for the years 2010–2011 and 2011–2012, where field data were available. Around 70% of rainfed plots
exhibited a positive WFI, while only half of irrigated plots had a negative WFI. This can be explained by the significant amount of rainfall towards the end of the 2011 agricultural season, as the WFI considers all precipitation regardless of when it occurs. As previously explained, the factor of irrigation versus the occurrence of rains in relation to the growing season plays a crucial role in interpreting the interannual variation of the WF. However, other factors such as cultural practices or soil characteristics can also contribute to the significant heterogeneity observed in WF, as noted in other studies [55–57]. The use of scatterplots (Figure 9, representing yields and WF) facilitated the illustration of this complexity. For instance, in the 2013 scatterplot, an abundance of pixels with low yield but varying WF from 1 to 6 can be observed, indicating the influence of diverse soils and cultural practices. These specific areas identified through this spatial study should undergo closer examination to validate this hypothesis [26,59].

**5. Conclusions**

To the best of our knowledge, this study represents the first attempt to assess both components of the water footprint (WF), evapotranspiration and grain yield, using remote sensing techniques for the Kairouan plain in Tunisia over three agricultural years (2010–2013). The study area is confronted with challenges stemming from the variability in precipitation and irrigation across different years, as well as spatial heterogeneity in the agricultural landscape and variable plot sizes. Therefore, to overcome the challenge of accurate knowledge of crop water consumption over large areas, we estimated the following.

For actual evapotranspiration from the S-SEBI model using Landsat 7–8, the results demonstrate a commendable agreement with data measured by the flux tower station, with an R² = 0.82 and an RMSE = 45 mm/day.

On the other hand, yields were mapped over the area based on an empirical model using SPOT NDVI, and modeled yield ranged between 6 and 77 qx/ha.

Globally, the cereals average WF varied between 1.08 m³/kg and 1.22 m³/kg over the three studied years. These results are generally in agreement with the literature indicating that WF is around 1 m³/kg for all cereal types [58]. The results confirm that in a dry year, the blue WF is more important than a year of average rainfall (WFb-2013 = 1.04 m³/kg, WFb-2012 = 0.61 m³/kg, WFb-2011 = 0.41 m³/kg).

For a more accurate interpretation of WF, the WFI, varying between -1 (WF mostly blue) and 1 (WF mostly green), was developed and mapped. This index makes it possible to distinguish between plots with various cereal WF.

In comparison to field data, it was observed that almost 70% of rainfed plots exhibited a positive WFI, while only half of the irrigated plots had a negative WFI. This discrepancy can be attributed to the calculation of the water footprint based on total precipitation, regardless of when it occurs.

A comprehensive agronomic analysis, integrating both yields and WF, identified four classes ranging from the most water-efficient to the least. This analysis revealed that more than 30% of cultivated areas during the study years (approximately 40% in 2011 and 2012 and 29% in 2013) exhibited low water efficiency, characterized by low yields and high WF. Further exploration within each class determined the distribution of green and blue water based on the WFI.

The study also brought to light that in 2012, 40% of cereal plots with low yields but high water consumption were irrigated (81% blue water compared to 6% in 2011). Consequently, it is recommended that managers prioritize these plots to enhance water efficiency and agronomic practices. Future studies should encompass a more extended period and incorporate additional data, as these techniques, with their inherent capabilities, are expected to offer further clarity and improvements in the assessment of the water footprint.

**Author Contributions:** Conceptualization, V.D., H.C. and A.C.B.; methodology, V.D., A.C.B., H.C., S.S. and G.B.; software, V.D., A.C.B., S.S. and H.C.; validation, G.B., S.S., A.C.B., V.D. and M.Z.; formal analysis, V.D., A.C.B., H.C., S.S. and G.B.; investigation, V.D., A.C.B., H.C., S.S. and G.B.; resources, V.D., A.C.B., H.C., S.S. and G.B.; data curation, V.D., H.C., A.C.B., S.S., G.B., M.Z. and
exhibited a positive WFI, while only half of irrigated plots had a negative WFI. This can be explained by the significant amount of rainfall towards the end of the 2011 agricultural season, as the WFI considers all precipitation regardless of when it occurs. As previously explained, the factor of irrigation versus the occurrence of rains in relation to the growing season plays a crucial role in interpreting the interannual variation of the WF. However, other factors such as cultural practices or soil characteristics can also contribute to the significant heterogeneity observed in WF, as noted in other studies [55–57]. The use of scatterplots (Figure 9, representing yields and WF) facilitated the illustration of this complexity. For instance, in the 2013 scatterplot, an abundance of pixels with low yield but varying WF from 1 to 6 can be observed, indicating the influence of diverse soils and cultural practices. These specific areas identified through this spatial study should undergo closer examination to validate this hypothesis [26,59].

**5. Conclusions**

To the best of our knowledge, this study represents the first attempt to assess both components of the water footprint (WF), evapotranspiration and grain yield, using remote sensing techniques for the Kairouan plain in Tunisia over three agricultural years (2010–2013). The study area is confronted with challenges stemming from the variability in precipitation and irrigation across different years, as well as spatial heterogeneity in the agricultural landscape and variable plot sizes. Therefore, to overcome the challenge of accurate knowledge of crop water consumption over large areas, we estimated the following.

For actual evapotranspiration from the S-SEBI model using Landsat 7–8, the results demonstrate a commendable agreement with data measured by the flux tower station, with an R² = 0.82 and an RMSE = 45 mm/day.

On the other hand, yields were mapped over the area based on an empirical model using SPOT NDVI, and modeled yield ranged between 6 and 77 qx/ha.

Globally, the cereals average WF varied between 1.08 m³/kg and 1.22 m³/kg over the three studied years. These results are generally in agreement with the literature indicating that WF is around 1 m³/kg for all cereal types [58]. The results confirm that in a dry year, the blue WF is more important than a year of average rainfall (WFb-2013 = 1.04 m³/kg, WFb-2012 = 0.61 m³/kg, WFb-2011 = 0.41 m³/kg).

For a more accurate interpretation of WF, the WFI, varying between -1 (WF mostly blue) and 1 (WF mostly green), was developed and mapped. This index makes it possible to distinguish between plots with various cereal WF.

In comparison to field data, it was observed that almost 70% of rainfed plots exhibited a positive WFI, while only half of the irrigated plots had a negative WFI. This discrepancy can be attributed to the calculation of the water footprint based on total precipitation, regardless of when it occurs.

A comprehensive agronomic analysis, integrating both yields and WF, identified four classes ranging from the most water-efficient to the least. This analysis revealed that more than 30% of cultivated areas during the study years (approximately 40% in 2011 and 2012 and 29% in 2013) exhibited low water efficiency, characterized by low yields and high WF. Further exploration within each class determined the distribution of green and blue water based on the WFI.

The study also brought to light that in 2012, 40% of cereal plots with low yields but high water consumption were irrigated (81% blue water compared to 6% in 2011). Consequently, it is recommended that managers prioritize these plots to enhance water efficiency and agronomic practices. Future studies should encompass a more extended period and incorporate additional data, as these techniques, with their inherent capabilities, are expected to offer further clarity and improvements in the assessment of the water footprint.

**Author Contributions:** Conceptualization, V.D., H.C. and A.C.B.; methodology, V.D., A.C.B., H.C., S.S. and G.B.; software, V.D., A.C.B., S.S. and H.C.; validation, G.B., S.S., A.C.B., V.D. and M.Z.; formal analysis, V.D., A.C.B., H.C., S.S. and G.B.; investigation, V.D., A.C.B., H.C., S.S. and G.B.; resources, V.D., A.C.B., H.C., S.S. and G.B.; data curation, V.D., H.C., A.C.B., S.S., G.B., M.Z. and
Z.L.C.; writing—original draft preparation, V.D. and H.C.; writing—review and editing, A.C.B., S.S., G.B., M.Z. and Z.L.C.; visualization, V.D., H.C. and A.C.B.; supervision, H.C.; project administration, H.C.; funding acquisition, H.C. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by the ERANETMED CHAAMS project (ERANET3-062).

**Data Availability Statement:** Data are contained within the article.

**Conflicts of Interest:** The authors declare no conflicts of interest.

## References

1. Falkenmark, M. Shift in Thinking to Address the 21st Century Hunger Gap. *Water Resour. Manag.* **2007**, *21*, 3–18. [CrossRef]

2. Mekonnen, M.M.; Hoekstra, A.Y. The Green, Blue and Grey Water Footprint of Crops and Derived Crop Products. *Hydrol. Earth Syst. Sci. Discuss.* **2011**, *8*, 763–809. [CrossRef]

3. Hoekstra, A.Y.; Chapagain, A.K.; Van Oel, P.R. Progress in Water Footprint Assessment. In *Progress in Water Footprint Assessment*; MDPI AG: Basel, Switzerland, 2019. [CrossRef]

4. Schyns, J.F.; Hoekstra, A.Y.; Booij, M.J. Review and Classification of Indicators of Green Water. *Hydrol. Earth Syst. Sci.* **2015**, *19*, 4581–4608. [CrossRef]

5. Hoekstra, A.Y. (Ed.) Virtual water trade. In *Proceedings of the International Expert Meeting on Virtual Water Trade, Delft, The Netherlands, 12–13 December 2002*; Value of Water Research Report Series No. 12; UNESCO-IHE: Delft, The Netherlands, 2003.

6. Besbes, M.; Chahed, J.; Hamdane, A. *National Water Security: Case Study of an Arid Country: Tunisia*; Springer: Berlin/Heidelberg, Germany, 2019; ISBN 9783319754994.

7. Besbes, M.; Chahed, J.; Hamdane, A. *A Sécurité Hydrique de la Tunisie*; L'Harmattan: Paris, France, 2014.

8. Chouchane, H.; Krol, M.S.; Hoekstra, A.Y. Virtual Water Trade Patterns in Relation to Environmental and Socioeconomic Factors: A Case Study for Tunisia. *Sci. Total Environ.* **2018**, *613–614*, 287–297. [CrossRef]

9. Chapagain, A.K.; Hoekstra, A.Y. *Water Footprints of Nations*; Value of Water Research Report Series No. 16; UNESCOIHE: Delft, The Netherlands, 2004; Volume 1–2.

10. Aldaya, M.M.; Chapagain, A.K.; Hoekstra, A.Y.; Mekonnen, M.M. *The Water Footprint Assessment Manual*; Setting the Global Standard; Routledge: London, UK, 2011; Volume 31, ISBN 9781849712798.

11. FAO. *Evaluation de L'approvisionnement Alimentaire dans un Contexte de Pénurie d'eau dans le Région NENA*; FAO: Rome, Italy, 2018; 165p, ISBN 978-92-5-130826-4.

12. Latrech; Lasram, A.; M'Nassri, S.; Nciri, R.; Gharbi, W.; Masmoudi, M.M.; Mechlia, N.B. Combining Remote Sensing Data and AquaCrop for Assessement of Yield and Water Productivity of Durum Wheat in Semi-Arid Conditions in Tunisia. Available online: https://swdcc2022.com/wp-content/uploads/2022/08/Session-FAO (accessed on 28 November 2023).

13. Zwart, S.J.; Bastiaanssen, W.G.M.; de Fraiture, C.; Molden, D.J. WATPRO: A Remote Sensing Based Model for Mapping Water Productivity of Wheat. *Agric. Water Manag.* **2010**, *97*, 1628–1636. [CrossRef]

14. Allen, R.G.; Pereira, L.S.; Raes, D.; Smith, M. *Crop Evapotranspiration-Guidelines for Computing Crop Water Requirements-FAO Irrigation and Drainage Paper 56*; Fao: Rome, Italy, 1998; Volume 300, D05109.

15. Allen, R.G.; Pereira, L.S.; Smith, M.; Raes, D.; Wright, J.L. *FAO-56 Dual Crop Coefficient Method for Estimating Evaporation from Soil and Application Extensions*; Fao: Rome, Italy, 2005.

16. Monteith, J.L. Evaporation and environment. *Symp. Soc. Exp. Biol.* **1965**, *19*, 205–234. [PubMed]

17. Calera, A.; Campos, I.; Osann, A.; D'urso, G.; Menenti, M. Remote Sensing for Crop Water Management: From ET Modelling to Services for the End Users. *Sensors* **2017**, *17*, 1104. [CrossRef]

18. Roerink, G.J.; Su, Z.; Menenti, M. S-SEBI: A simple remote sensing algorithm to estimate the surface energy balance. *Phys. Chem. Earth Part B Hydrol. Oceans Atmos.* **2000**, *25*, 147–157. [CrossRef]

19. Li, Z.-L.; Tang, R.; Wan, Z.; Bi, Y.; Zhou, C.; Tang, B.; Yan, G.; Zhang, X. A review of current methodologies for regional evapotranspiration estimation from remotely sensed data. *Sensors* **2009**, *9*, 3801–3853. [CrossRef]

20. Galleguillos, M.; Jacob, F.; Prévot, L.; French, A.; Lagacherie, P. Comparison of two temperature differencing methods to estimate daily evapotranspiration over a Mediterranean vineyard watershed from ASTER data. *Remote Sens. Environ.* **2011**, *115*, 1326–1340. [CrossRef]

21. Acharya, B.; Sharma, V. Comparison of satellite driven surface energy balance models in estimating crop evapotranspiration in semi-arid to arid inter-mountain region. *Remote Sens.* **2021**, *13*, 1822. [CrossRef]

22. Van Ittersum, M.K.; Ewert, F.; Heckelei, T.; Wery, J.; Olsson, J.A.; Andersen, E.; Bezlepkina, I.; Brouwer, F.; Donatelli, M.; Flichman, G.; et al. Integrated assessment of agricultural systems—A component-based framework for the European Union (SEAMLESS). *Agric. Syst.* **2008**, *96*, 150–165. [CrossRef]

23. Balaghi, R.; Tychon, B.; Eerens, H.; Jlibene, M. Empirical Regression Models Using NDVI, Rainfall and Temperature Data for the Early Prediction of Wheat Grain Yields in Morocco. *Int. J. Appl. Earth Obs. Geoinf.* **2008**, *10*, 438–452. [CrossRef]

24. Becker-Reshef, I.; Vermote, E.; Lindeman, M.; Justice, C. A generalized regression-based model for forecasting winter wheat yields in Kansas and Ukraine using MODIS Data. *Remote Sens. Environ.* **2010**, *114*, 1312–1323. [CrossRef]
Z.L.C.; writing—original draft preparation, V.D. and H.C.; writing—review and editing, A.C.B., S.S., G.B., M.Z. and Z.L.C.; visualization, V.D., H.C. and A.C.B.; supervision, H.C.; project administration, H.C.; funding acquisition, H.C. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by the ERANETMED CHAAMS project (ERANET3-062).

**Data Availability Statement:** Data are contained within the article.

**Conflicts of Interest:** The authors declare no conflicts of interest.

## References

1. Falkenmark, M. Shift in Thinking to Address the 21st Century Hunger Gap. *Water Resour. Manag.* **2007**, *21*, 3–18. [CrossRef]

2. Mekonnen, M.M.; Hoekstra, A.Y. The Green, Blue and Grey Water Footprint of Crops and Derived Crop Products. *Hydrol. Earth Syst. Sci. Discuss.* **2011**, *8*, 763–809. [CrossRef]

3. Hoekstra, A.Y.; Chapagain, A.K.; Van Oel, P.R. Progress in Water Footprint Assessment. In *Progress in Water Footprint Assessment*; MDPI AG: Basel, Switzerland, 2019. [CrossRef]

4. Schyns, J.F.; Hoekstra, A.Y.; Booij, M.J. Review and Classification of Indicators of Green Water. *Hydrol. Earth Syst. Sci.* **2015**, *19*, 4581–4608. [CrossRef]

5. Hoekstra, A.Y. (Ed.) Virtual water trade. In *Proceedings of the International Expert Meeting on Virtual Water Trade, Delft, The Netherlands, 12–13 December 2002*; Value of Water Research Report Series No. 12; UNESCO-IHE: Delft, The Netherlands, 2003.

6. Besbes, M.; Chahed, J.; Hamdane, A. *National Water Security: Case Study of an Arid Country: Tunisia*; Springer: Berlin/Heidelberg, Germany, 2019; ISBN 9783319754994.

7. Besbes, M.; Chahed, J.; Hamdane, A. *A Sécurité Hydrique de la Tunisie*; L'Harmattan: Paris, France, 2014.

8. Chouchane, H.; Krol, M.S.; Hoekstra, A.Y. Virtual Water Trade Patterns in Relation to Environmental and Socioeconomic Factors: A Case Study for Tunisia. *Sci. Total Environ.* **2018**, *613–614*, 287–297. [CrossRef]

9. Chapagain, A.K.; Hoekstra, A.Y. *Water Footprints of Nations*; Value of Water Research Report Series No. 16; UNESCOIHE: Delft, The Netherlands, 2004; Volume 1–2.

10. Aldaya, M.M.; Chapagain, A.K.; Hoekstra, A.Y.; Mekonnen, M.M. *The Water Footprint Assessment Manual*; Setting the Global Standard; Routledge: London, UK, 2011; Volume 31, ISBN 9781849712798.

11. FAO. *Evaluation de L'approvisionnement Alimentaire dans un Contexte de Pénurie d'eau dans le Région NENA*; FAO: Rome, Italy, 2018; 165p, ISBN 978-92-5-130826-4.

12. Latrech; Lasram, A.; M'Nassri, S.; Nciri, R.; Gharbi, W.; Masmoudi, M.M.; Mechlia, N.B. Combining Remote Sensing Data and AquaCrop for Assessement of Yield and Water Productivity of Durum Wheat in Semi-Arid Conditions in Tunisia. Available online: https://swdde2022.com/wp-content/uploads/2022/08/Session-FAO (accessed on 28 November 2023).

13. Zwart, S.J.; Bastiaanssen, W.G.M.; de Fraiture, C.; Molden, D.J. WATPRO: A Remote Sensing Based Model for Mapping Water Productivity of Wheat. *Agric. Water Manag.* **2010**, *97*, 1628–1636. [CrossRef]

14. Allen, R.G.; Pereira, L.S.; Raes, D.; Smith, M. *Crop Evapotranspiration-Guidelines for Computing Crop Water Requirements-FAO Irrigation and Drainage Paper 56*; Fao: Rome, Italy, 1998; Volume 300, D05109.

15. Allen, R.G.; Pereira, L.S.; Smith, M.; Raes, D.; Wright, J.L. *FAO-56 Dual Crop Coefficient Method for Estimating Evaporation from Soil and Application Extensions*; Fao: Rome, Italy, 2005.

16. Monteith, J.L. Evaporation and environment. *Symp. Soc. Exp. Biol.* **1965**, *19*, 205–234. [PubMed]

17. Calera, A.; Campos, I.; Osann, A.; D'urso, G.; Menenti, M. Remote Sensing for Crop Water Management: From ET Modelling to Services for the End Users. *Sensors* **2017**, *17*, 1104. [CrossRef]

18. Roerink, G.J.; Su, Z.; Menenti, M. S-SEBI: A simple remote sensing algorithm to estimate the surface energy balance. *Phys. Chem. Earth Part B Hydrol. Oceans Atmos.* **2000**, *25*, 147–157. [CrossRef]

19. Li, Z.-L.; Tang, R.; Wan, Z.; Bi, Y.; Zhou, C.; Tang, B.; Yan, G.; Zhang, X. A review of current methodologies for regional evapotranspiration estimation from remotely sensed data. *Sensors* **2009**, *9*, 3801–3853. [CrossRef]

20. Galleguillos, M.; Jacob, F.; Prévot, L.; French, A.; Lagacherie, P. Comparison of two temperature differencing methods to estimate daily evapotranspiration over a Mediterranean vineyard watershed from ASTER data. *Remote Sens. Environ.* **2011**, *115*, 1326–1340. [CrossRef]

21. Acharya, B.; Sharma, V. Comparison of satellite driven surface energy balance models in estimating crop evapotranspiration in semi-arid to arid inter-mountain region. *Remote Sens.* **2021**, *13*, 1822. [CrossRef]

22. Van Ittersum, M.K.; Ewert, F.; Heckelei, T.; Wery, J.; Olsson, J.A.; Andersen, E.; Bezlepkina, I.; Brouwer, F.; Donatelli, M.; Flichman, G.; et al. Integrated assessment of agricultural systems—A component-based framework for the European Union (SEAMLESS). *Agric. Syst.* **2008**, *96*, 150–165. [CrossRef]

23. Balaghi, R.; Tychon, B.; Eerens, H.; Jlibene, M. Empirical Regression Models Using NDVI, Rainfall and Temperature Data for the Early Prediction of Wheat Grain Yields in Morocco. *Int. J. Appl. Earth Obs. Geoinf.* **2008**, *10*, 438–452. [CrossRef]

24. Becker-Reshef, I.; Vermote, E.; Lindeman, M.; Justice, C. A generalized regression-based model for forecasting winter wheat yields in Kansas and Ukraine using MODIS Data. *Remote Sens. Environ.* **2010**, *114*, 1312–1323. [CrossRef]
25. Franch, B.; Vermote, E.F.; Becker-Reshef, I.; Claverie, M.; Huang, J.; Zhang, J.; Justice, C.; Sobrino, J.A. Improving the Timeliness of Winter Wheat Production Forecast in the United States of America, Ukraine and China Using MODIS Data and NCAR Growing Degree Day Information. *Remote Sens. Environ.* **2015**, *161*, 131–148. [CrossRef]

26. Chahbi, A.; Zribi, M.; Lili-Chabaane, Z.; Duchemin, B.; Shabou, M.; Mougenot, B.; Boulet, G. Estimation of the dynamics and yields of cereals in a semi-arid area using remote sensing and the SAFY growth model. *Int. J. Remote Sens.* **2014**, *35*, 1004–1028. [CrossRef]

27. Nougaret, G.; Jenhaoui, Z.; Leduc, C. Ressources en Eau Dans Le Kairouanais: Évolutions des Disponibilités et des Usages DEPUIS 2000 Ans. 2019. Available online: http://www5.funceme.br/arid/wp-content/uploads/2020/05/guide_terrain_Merguellil_04-10-2019-GN-CL.pdf (accessed on 28 November 2023).

28. Zribi, M.; Chahbi, A.; Shabou, M.; Lili-Chabaane, Z.; Duchemin, B.; Baghdadi, N.; Amri, R.; Chehbouni, A. Soil Surface Moisture Estimation over a Semi-Arid Region Using ENVISAT ASAR Radar Data for Soil Evaporation Evaluation. *Hydrol. Earth Syst. Sci.* **2011**, *15*, 345–358. [CrossRef]

29. Luca, C.; Michele, M.; Silvia, M. *Investigating the Relationship between Land Cover and Vulnerability* T0 Climate Change in dar es Salaam; Working Paper April 2013; Sapienza University: Rome, Italy, 2013.

30. Chavez, P.S. An improved dark-object subtraction technique for atmospheric scattering correction of multispectral data. *Remote Sens. Environ.* **1988**, *24*, 459–479. [CrossRef]

31. Warmerdam, F. The Geospatial Data Abstraction Library. *Open Source Approaches Spat. Data Handl.* **2008**, *2*, 87–104. [CrossRef]

32. Weng, Q.; Lu, D.S.; Schubring, J. Estimation of land surface temperature—Vegetation abundance relationship for urban heat island studies. *Remote Sens. Environ.* **2004**, *89*, 467–483. [CrossRef]

33. Saadi, S.; Simonneaux, V.; Boulet, G.; Raimbault, B.; Mougenot, B.; Fanise, P.; Ayari, H.; Lili-Chabaane, Z. Monitoring Irrigation Consumption Using High Resolution NDVI Image Time Series: Calibration and Validation in the Kairouan Plain (Tunisia). *Remote Sens.* **2015**, *7*, 13005–13028. [CrossRef]

34. Rahman, H.; Dedieu, G. SMAC: A simplified method for the atmospheric correction of satellite measurements in the solar spectrum. *Int. J. Remote Sens.* **1994**, *15*, 123–143. [CrossRef]

35. Berthelot, B.; Dedieu, G. Correction of Atmospheric Effects for VEGETATION Data. In *Physical Measurements and Signatures in Remote Sensing*; Guyot, G., Phulpin, T., Eds.; Courchevel: Balkema, France, 1997; pp. 19–25.

36. Kotchenova, S.Y.; Vermote, E.F.; Matarrese, R.; Klemm, F.J., Jr. Validation of a vector version of the 6S radiative transfer code for atmospheric correction of satellite data. Part I: Path radiance. *Appl. Opt.* **2006**, *45*, 6762–6774. [CrossRef] [PubMed]

37. Tanré, D.; Deroo, C.; Duhaut, P.; Herman, M.; Morcrette, J.J.; Perbos, J.; Deschamps, P.Y. Description of a computer code to simulate the satellite signal in the solar spectrum: The 5S code. *Int. J. Remote Sens.* **1990**, *11*, 659–668. [CrossRef]

38. Caloz, R.; Collet, C. *Précis de Télédétection, Volume 3 Traitements Numériques D'images de Télédétection*; Presses de l'Université du Québec/AUPELF: Sainte-Foy, QC, Canada, 2001; 386p.

39. Braden, H. Ein Energiehaushalts- und Verdunstungsmodell for Wasser und Stoffhaushaltsuntersuchungen landwirtschaftlich genutzer Einzugsgebiete. *M. Dtsch. Bodenkd. Gesellschaft* **1985**, *42*, 294–299.

40. Moran, M.S.; Clarke, T.R.; Inoue, Y.; Vidal, A. Estimating crop water deficit using the relation between surface-air temperature and spectral vegetation index. *Remote Sens. Environ.* **1994**, *49*, 246–263. [CrossRef]

41. Tucker, J. Red and Photographic Infrared Linear Combinations for Monitoring Vegetation. *Remote Sens. Environ.* **1979**, *150*, 127–150. [CrossRef]

42. Norman, J.M.; Kustas, W.P.; Humes, K.S. Source Approach for Estimating Soil and Vegetation Energy Fluxes in Observations of Directional Radiometric Surface Temperature. *Agric. For. Meteorol.* **1995**, *77*, 263–293. [CrossRef]

43. Chemin, Y. *Evapotranspiration of Crops by Remote Sensing Using the Energy Balance Based Algorithms*; International Water Management Institute: Pathumithani, Thailand, 2003.

44. Choudhury, B.J.; Ahmed, N.U.; Idso, S.B.; Reginato, R.J.; Daughtry, C.S.T. Relations between evaporation coefficients and vegetation indices studied by model simulations. *Remote Sens. Environ.* **1994**, *50*, 1–17. [CrossRef]

45. Valor, E.; Caselles, V. Mapping land surface emissivity from NDVI: Application to European, African, and South American areas. *Remote Sens. Environ.* **1996**, *57*, 167–184. [CrossRef]

46. Soudani, K.; Christophe, F.; Et Valérie, L.D. Analyse comparative d'IKONOS, Données SPOT et ETM+ pour l'estimation de l'indice de surface foliaire des conifères et des feuillus tempérés peuplements forestiers. *Télédétection L'environnement* **2006**, *102*, 161–175.

47. Chakroun, H.; Zemni, N.; Benhmid, A.; Dellaly, V.; Slama, F.; Bouksila, F.; Berndtsson, R. Evapotranspiration in Semi-Arid Climate: Remote Sensing vs. Soil Water Simulation. *Sensors* **2023**, *23*, 2823. [CrossRef]

48. Kustas, W.P.; Daughtry, C.S. Estimation of the soil heat flux/net radiation ratio from spectral data. *Agric. For. Meteorol.* **1990**, *49*, 205–223. [CrossRef]

49. Ogée, J.; Lamaud, E.; Brunet, Y.; Berbigier, P.; Bonnefond, J.-M. A long-term study of soil heat flux under a forest canopy. *Agric. For. Meteorol.* **2001**, *106*, 173–186. [CrossRef]

50. Choudhury, B.J. Relationships between vegetation indices, radiation absorption, and net photosynthesis evaluated by a 638 sensitivity analysis. *Remote Sens. Environ.* **1987**, *22*, 209–233. [CrossRef]

51. Santanello, J.A., Jr.; Friedl, M.A. Diurnal covariation in soil heat flux and net radiation. *J. Appl. Meteorol.* **2003**, *42*, 851–862. [CrossRef]
25. Franch, B.; Vermote, E.F.; Becker-Reshef, I.; Claverie, M.; Huang, J.; Zhang, J.; Justice, C.; Sobrino, J.A. Improving the Timeliness of Winter Wheat Production Forecast in the United States of America, Ukraine and China Using MODIS Data and NCAR Growing Degree Day Information. *Remote Sens. Environ.* **2015**, *161*, 131–148. [CrossRef]

26. Chahbi, A.; Zribi, M.; Lili-Chabaane, Z.; Duchemin, B.; Shabou, M.; Mougenot, B.; Boulet, G. Estimation of the dynamics and yields of cereals in a semi-arid area using remote sensing and the SAFY growth model. *Int. J. Remote Sens.* **2014**, *35*, 1004–1028. [CrossRef]

27. Nougaret, G.; Jenhaoui, Z.; Leduc, C. Ressources en Eau Dans Le Kairouanais: Évolutions des Disponibilités et des Usages DEPUIS 2000 Ans. 2019. Available online: http://www5.funceme.br/arid/wp-content/uploads/2020/05/guide_terrain_Merguellil_04-10-2019-GN-CL.pdf (accessed on 28 November 2023).

28. Zribi, M.; Chahbi, A.; Shabou, M.; Lili-Chabaane, Z.; Duchemin, B.; Baghdadi, N.; Amri, R.; Chehbouni, A. Soil Surface Moisture Estimation over a Semi-Arid Region Using ENVISAT ASAR Radar Data for Soil Evaporation Evaluation. *Hydrol. Earth Syst. Sci.* **2011**, *15*, 345–358. [CrossRef]

29. Luca, C.; Michele, M.; Silvia, M. *Investigating the Relationship between Land Cover and Vulnerability* T0 Climate Change in dar es Salaam; Working Paper April 2013; Sapienza University: Rome, Italy, 2013.

30. Chavez, P.S. An improved dark-object subtraction technique for atmospheric scattering correction of multispectral data. *Remote Sens. Environ.* **1988**, *24*, 459–479. [CrossRef]

31. Warmerdam, F. The Geospatial Data Abstraction Library. *Open Source Approaches Spat. Data Handl.* **2008**, *2*, 87–104. [CrossRef]

32. Weng, Q.; Lu, D.S.; Schubring, J. Estimation of land surface temperature—Vegetation abundance relationship for urban heat island studies. *Remote Sens. Environ.* **2004**, *89*, 467–483. [CrossRef]

33. Saadi, S.; Simonneaux, V.; Boulet, G.; Raimbault, B.; Mougenot, B.; Fanise, P.; Ayari, H.; Lili-Chabaane, Z. Monitoring Irrigation Consumption Using High Resolution NDVI Image Time Series: Calibration and Validation in the Kairouan Plain (Tunisia). *Remote Sens.* **2015**, *7*, 13005–13028. [CrossRef]

34. Rahman, H.; Dedieu, G. SMAC: A simplified method for the atmospheric correction of satellite measurements in the solar spectrum. *Int. J. Remote Sens.* **1994**, *15*, 123–143. [CrossRef]

35. Berthelot, B.; Dedieu, G. Correction of Atmospheric Effects for VEGETATION Data. In *Physical Measurements and Signatures in Remote Sensing*; Guyot, G., Phulpin, T., Eds.; Courchevel: Balkema, France, 1997; pp. 19–25.

36. Kotchenova, S.Y.; Vermote, E.F.; Matarrese, R.; Klemm, F.J., Jr. Validation of a vector version of the 6S radiative transfer code for atmospheric correction of satellite data. Part I: Path radiance. *Appl. Opt.* **2006**, *45*, 6762–6774. [CrossRef] [PubMed]

37. Tanré, D.; Deroo, C.; Duhaut, P.; Herman, M.; Morcrette, J.J.; Perbos, J.; Deschamps, P.Y. Description of a computer code to simulate the satellite signal in the solar spectrum: The 5S code. *Int. J. Remote Sens.* **1990**, *11*, 659–668. [CrossRef]

38. Caloz, R.; Collet, C. *Précis de Télédétection, Volume 3 Traitements Numériques D'images de Télédétection*; Presses de l'Université du Québec/AUPELF: Sainte-Foy, QC, Canada, 2001; 386p.

39. Braden, H. Ein Energiehaushalts- und Verdunstungsmodell for Wasser und Stoffhaushaltsuntersuchungen landwirtschaftlich genutzer Einzugsgebiete. *M. Dtsch. Bodenkd. Gesellschaft* **1985**, *42*, 294–299.

40. Moran, M.S.; Clarke, T.R.; Inoue, Y.; Vidal, A. Estimating crop water deficit using the relation between surface-air temperature and spectral vegetation index. *Remote Sens. Environ.* **1994**, *49*, 246–263. [CrossRef]

41. Tucker, J. Red and Photographic Infrared Linear Combinations for Monitoring Vegetation. *Remote Sens. Environ.* **1979**, *150*, 127–150. [CrossRef]

42. Norman, J.M.; Kustas, W.P.; Humes, K.S. Source Approach for Estimating Soil and Vegetation Energy Fluxes in Observations of Directional Radiometric Surface Temperature. *Agric. For. Meteorol.* **1995**, *77*, 263–293. [CrossRef]

43. Chemin, Y. *Evapotranspiration of Crops by Remote Sensing Using the Energy Balance Based Algorithms*; International Water Management Institute: Pathumithani, Thailand, 2003.

44. Choudhury, B.J.; Ahmed, N.U.; Idso, S.B.; Reginato, R.J.; Daughtry, C.S.T. Relations between evaporation coefficients and vegetation indices studied by model simulations. *Remote Sens. Environ.* **1994**, *50*, 1–17. [CrossRef]

45. Valor, E.; Caselles, V. Mapping land surface emissivity from NDVI: Application to European, African, and South American areas. *Remote Sens. Environ.* **1996**, *57*, 167–184. [CrossRef]

46. Soudani, K.; Christophe, F.; Et Valérie, L.D. Analyse comparative d'IKONOS, Données SPOT et ETM+ pour l'estimation de l'indice de surface foliaire des conifères et des feuillus tempérés peuplements forestiers. *Télédétection L'environnement* **2006**, *102*, 161–175.

47. Chakroun, H.; Zemni, N.; Benhmid, A.; Dellaly, V.; Slama, F.; Bouksila, F.; Berndtsson, R. Evapotranspiration in Semi-Arid Climate: Remote Sensing vs. Soil Water Simulation. *Sensors* **2023**, *23*, 2823. [CrossRef]

48. Kustas, W.P.; Daughtry, C.S. Estimation of the soil heat flux/net radiation ratio from spectral data. *Agric. For. Meteorol.* **1990**, *49*, 205–223. [CrossRef]

49. Ogée, J.; Lamaud, E.; Brunet, Y.; Berbigier, P.; Bonnefond, J.-M. A long-term study of soil heat flux under a forest canopy. *Agric. For. Meteorol.* **2001**, *106*, 173–186. [CrossRef]

50. Choudhury, B.J. Relationships between vegetation indices, radiation absorption, and net photosynthesis evaluated by a 638 sensitivity analysis. *Remote Sens. Environ.* **1987**, *22*, 209–233. [CrossRef]

51. Santanello, J.A., Jr.; Friedl, M.A. Diurnal covariation in soil heat flux and net radiation. *J. Appl. Meteorol.* **2003**, *42*, 851–862. [CrossRef]
52. Sobrino, J.A.; Gómez, M.; Jiménez-Muñoz, J.C.; Olioso, A. Application of a Simple Algorithm to Estimate Daily Evapotranspiration from NOAA-AVHRR Images for the Iberian Peninsula. *Remote Sens. Environ.* **2007**, *110*, 139-148. [CrossRef]

53. Abdallah, A. *Eau Virtuelle et Sécurité Alimentaire en Tunisie: Du Constat à l'Appui au Développement (EVSAT-CAD)*; Projet de Recherche Developpement: Mograne, Tunisia, 2015.

54. Chahed, J.; Hamdane, A.; Besbes, M. A comprehensive water balance of Tunisia: Blue water, green water and virtual water. *Water Int.* **2008**, *33*, 415-424. [CrossRef]

55. Ghorbanpour, A.K.; Kisekka, I.; Afshar, A.; Hessels, T.; Taraghi, M.; Hessari, B.; Touri-an, M.J.; Duan, Z. Crop Water Productivity Mapping and Benchmarking Using Re-mote Sensing and Google Earth Engine Cloud Computing. *Remote Sens.* **2022**, *14*, 4934. [CrossRef]

56. Bastiaanssen, W.G.M.; Steduto, P. The Water Productivity Score (WPS) at Global and Regional Level: Methodology and First Results from Remote Sensing Measurements of Wheat, Rice and Maize. *Sci. Total Environ.* **2017**, *575*, 595-611. [CrossRef]

57. Blatchford, M.L.; Karimi, P.; Bastiaanssen, W.; Nouri, H. From Global Goals to Local Gains—A Framework for Crop Water Productivity. *ISPRS Int. J. Geo-Inf.* **2018**, *7*, 414. [CrossRef]

58. Chahed, Y.; Hassan, F.A. 2012 Grain and Feed Update Tunisia. In *Global Agricultural Information Network*; TS1204; USDA Foreign Agriculture Service: Washington, DC, USA, 2012.

59. Khlif, M.; Escorihuela, M.J.; Bellakanji, A.C.; Paolini, G.; Kassouk, Z.; Chabaane, Z.L. Multi-Year Cereal Crop Classification Model Using Sentinel 2 and 3 Landsat 7-8 Data. *Agriculture* **2023**, *13*, 1633. [CrossRef]

60. Bousbih, S.; Zribi, M.; Lili-Chabaane, Z.; Baghdadi, N.; El Hajj, M.; Gao, Q.; Mougenot, B. Potential of Sentinel-1 Radar Data for the Assessment of Soil and Cereal Cover Parameters. *Sensors* **2017**, *17*, 2617. [CrossRef]

61. Chouchane, H. Economic Allocation of Water to Crops in International Context: A National and Global Perspective. Ph.D. Thesis, University of Twente, Enschede, The Netherlands, 2019.

**Disclaimer/Publisher's Note:** The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.
52. Sobrino, J.A.; Gómez, M.; Jiménez-Muñoz, J.C.; Olioso, A. Application of a Simple Algorithm to Estimate Daily Evapotranspiration from NOAA-AVHRR Images for the Iberian Peninsula. *Remote Sens. Environ.* **2007**, *110*, 139-148. [CrossRef]

53. Abdallah, A. *Eau Virtuelle et Sécurité Alimentaire en Tunisie: Du Constat à l'Appui au Développement (EVSAT-CAD)*; Projet de Recherche Developpement: Mograne, Tunisia, 2015.

54. Chahed, J.; Hamdane, A.; Besbes, M. A comprehensive water balance of Tunisia: Blue water, green water and virtual water. *Water Int.* **2008**, *33*, 415-424. [CrossRef]

55. Ghorbanpour, A.K.; Kisekka, I.; Afshar, A.; Hessels, T.; Taraghi, M.; Hessari, B.; Touri-an, M.J.; Duan, Z. Crop Water Productivity Mapping and Benchmarking Using Re-mote Sensing and Google Earth Engine Cloud Computing. *Remote Sens.* **2022**, *14*, 4934. [CrossRef]

56. Bastiaanssen, W.G.M.; Steduto, P. The Water Productivity Score (WPS) at Global and Regional Level: Methodology and First Results from Remote Sensing Measurements of Wheat, Rice and Maize. *Sci. Total Environ.* **2017**, *575*, 595-611. [CrossRef]

57. Blatchford, M.L.; Karimi, P.; Bastiaanssen, W.; Nouri, H. From Global Goals to Local Gains—A Framework for Crop Water Productivity. *ISPRS Int. J. Geo-Inf.* **2018**, *7*, 414. [CrossRef]

58. Chahed, Y.; Hassan, F.A. 2012 Grain and Feed Update Tunisia. In *Global Agricultural Information Network*; TS1204; USDA Foreign Agriculture Service: Washington, DC, USA, 2012.

59. Khlif, M.; Escorihuela, M.J.; Bellakanji, A.C.; Paolini, G.; Kassouk, Z.; Chabaane, Z.L. Multi-Year Cereal Crop Classification Model Using Sentinel 2 and 3 Landsat 7-8 Data. *Agriculture* **2023**, *13*, 1633. [CrossRef]

60. Bousbih, S.; Zribi, M.; Lili-Chabaane, Z.; Baghdadi, N.; El Hajj, M.; Gao, Q.; Mougenot, B. Potential of Sentinel-1 Radar Data for the Assessment of Soil and Cereal Cover Parameters. *Sensors* **2017**, *17*, 2617. [CrossRef]

61. Chouchane, H. Economic Allocation of Water to Crops in International Context: A National and Global Perspective. Ph.D. Thesis, University of Twente, Enschede, The Netherlands, 2019.

**Disclaimer/Publisher's Note:** The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.