2.1. Study Area and Vineyard Description
The study area is located in the Apulia region (
Figure 1) in the southeastern part of the Italian peninsula, where over 70% of the total area is occupied by agriculture and vineyards cover about 26% of the crop area [
40]. The regional topography is mainly flat or slightly sloping, with the exception of the Gargano area, situated in the northwest of the region. The climate of the study area is classified as Mediterranean semi-arid, characterized by moderately cold and rainy winters and dry summer seasons. Mean annual precipitation is about 550 mm, concentrated in the fall and winter; long-term mean air temperature is 15.4 °C, while the minimum and maximum yearly temperatures are 3.5 and 29.5 °C, respectively. The climate, rainfall in particular, exhibits a marked inter-annual variability that makes water availability a permanent threat to the economic development and ecosystem conservation of the region [
41]. The soil is typical of the Capitanata region, with a well-differentiated profile, characterized by Quaternary marine sediments of illitic/smectitic composition that originate clay soils with vertic features (cambisols, entisols, vertisols) [
42]. Due to its geo-climatic conditions, the region suffers from overall water exploitation, which exposes the water supply system to severe water scarcity events.
Figure 1.
Study area location with overview of LAI plot measurements. Natural color composition of Landsat 8 image acquired on 26 August 2014.
Figure 1.
Study area location with overview of LAI plot measurements. Natural color composition of Landsat 8 image acquired on 26 August 2014.
The investigations have been carried within the Land Reclamation and Irrigation Consortium Capitanata, one of the largest agricultural irrigation districts of Southern Italy (
Figure 1). This region mainly contains vineyards and orchards over flat terrain. In particular, field data have been collected in a 13 ha multiplot vineyard (41°17ʹ29ʹʹ North latitude, 16°07ʹ22ʹʹ West longitude) trained on a “tendone” system, with commercial table grapes (cultivars Italia, Victoria, and RedGlobe), uncovered and drip irrigated (
Figure 2). The row and vine spacing in the vineyard was generally 2.3 × 2.3 m, about 2 m high with a plant density of 1890 plants·ha
−1. Agronomic management generally maintains the vegetation under non-water stress conditions and soil surface free of weeds and cover crop.
Figure 2.
(a) Overall view of the vineyard study field displayed on a RapidEye image acquired on 2 August 2014 (infrared false color composition 532). Yellow borders show the experimental plots. (b) Canopy growth of the vineyard inter-row acquired on 24 July 2014.
Figure 2.
(a) Overall view of the vineyard study field displayed on a RapidEye image acquired on 2 August 2014 (infrared false color composition 532). Yellow borders show the experimental plots. (b) Canopy growth of the vineyard inter-row acquired on 24 July 2014.
2.3. Methods
The application of the procedure is achieved by the integration of EO techniques and meteorological data of the study area. The conceptualization has been developed by D’Urso [
45] and implemented into the SIMODIS model (SImulation and Management of On-Demand Irrigation System). The estimation of crop evapotranspiration under standard conditions (ET
p—mm/d)—disease-free, adequate fertilization, and soil water availability—has been carried out using the FAO’s Penman–Monteith method, which requires standard meteorological data such as solar radiation (S), air temperature (T
a), air humidity (RH), wind speed (U), and biophysical parameters that are crop-specific, such as albedo (
r), LAI, and crop height (
hc) [
20]. The value of ET
p is a function of the crop parameters and meteorological data, and the procedure, known as the one-step approach, is computed as:
Since crop parameters (mainly
r and LAI) are very difficult to obtain in the open field without specific instruments, ET
p is calculated in the routine irrigation management by multiplying ET
0 by K
c, starting from an ET
0 assuming an ideal crop with standard parameters (
hc = 0.12 m,
r = 0.23, LAI = 2.88) and a K
c extracted from the FAO-56 table, following a single K
c approach [
20].
Figure 3.
Field estimates of LAI collected on 23–24 July 2014 with a LICOR LAI-2000 Plant Canopy Analyzer. The filled circles show the average LAI with the standard error.
Figure 3.
Field estimates of LAI collected on 23–24 July 2014 with a LICOR LAI-2000 Plant Canopy Analyzer. The filled circles show the average LAI with the standard error.
It should be noted that ET
p calculated in this way represents the maximum value, (
i.e., minimum canopy resistance), so the resulting K
c values, if compared with values reported in the literature, should be used with caution. The comparison between the results of the proposed methodology and field measurements—including fluxes of latent heat by means of micrometeorological instrumentations—requires fully controlled conditions (especially for irrigation), which were not available at the current study site. However, we assume that the model can be applied to the studied canopy, having been validated in the course of field experiments on different fully irrigated crops [
29,
39]. Actual evapotranspiration from a well-watered crop will generally approach ET
p during the active growing stage, but it falls below it toward the end of the growing season as the plants begin to dry out. Diversely from our ET
0 estimates, actual evapotranspiration depends on the conditions of soil water availability, and is only possible to compute through complex measurements in the open field or by filling out a detailed hydrological balance, taking into account all the parameters for the simulation of water transport/flow in the Soil, Plant Atmosphere system, using numerical schematization. In this case, the value of ET
p determined with the approach proposed here is used as an upper limiting condition [
46]. This kind of analysis has been carried out for some perennial tree crops in Sicily [
47], where soil water and surface energy balances with visible, near, and thermal infrared data have been carried out. These previous studies support the assumption of validity, provided that the canopy parameters,
i.e., LAI, are adequately described. Finally, crop water requirements can be accounted for in a simple way by subtracting the net precipitation (
Pn) from ET
p.
In order to describe the amount of intercepted water from the plant surface,
Pn is calculated as a function of the actual precipitation (
P), LAI, and fractional cover. The semi-empirical model of interception [
48] is described by the following relationship:
where
P (cm·d
−1) and
a (cm·d
−1) are empirical parameters representing the crop saturation per unit foliage area (~0.28 for most crops), and
fcover is the fractional vegetation cover derived from LAI using a polynomial empirical expression where coefficients are determined from field measurements and are valid for a wide range of crops (LAI ≤ 5 m
2·m
−2).
K
c values are extremely variable, even within the same type of crop, depending on many factors,
i.e., date and seeding density, intake of nutrients, nature of the soil, and agronomic practices. K
c is basically the ratio of the ET
p to the ET
0, and can be computed as follows [
45,
49]:
where
r*, LAI
* and
hc* represent the parameters of the crop vegetation present at the time of the satellite overpass, which can be related to satellite observations, while
K↓ represent the global incoming short-wave radiation flux density. The value of K
c depends both on the crop parameters (except albedo, which is also influenced by the degree of humidity of the soil), which can be considered relatively stable for a period of 1–2 weeks, and the weather conditions of the study area, which instead vary in time. A sensitivity analysis of the value of K
c with respect to these parameters showed that the value of K
c is more closely related to vegetation parameters
r* and LAI
*, and that the sensitivity of the change in the crop height is higher in autumn than in summer [
45]. This is due to the fact that, at the daily scale, the aerodynamic component of evapotranspiration is of much less importance than the radiative. A flow-chart of the processing chain of this procedure can be found in Vuolo
et al. [
39]. The validation has been carried out over different crops (including vineyards), thereby fulfilling the “standard conditions” for ET
p by means of eddy-covariance independent measurements [
29,
39]. Since LAI represents the main parameter, it can be taken as the reference variable for the calibration of this method to local conditions, as described later.
Recent studies conducted in semi-arid areas confirm this trend [
50]; in fact, a percentage change of 50% in
hc* corresponds to a variation of 5% of K
c. Therefore, a constant value of
hc = 0.4 m has been set in this study, valid for the climatic conditions of the Mediterranean area, reducing the calculation to the estimation of
r and LAI.
Data Processing for Deriving EO-Based Crop Development Maps
A total of 14 Landsat 8 images and 10 RapidEye images were used for both growing seasons. The pre-processing phase of satellite images was achieved in two steps. Firstly, Landsat 8 scenes were orthorectified by using additional Ground Control Points and a Digital Elevation Model, while RapidEye scenes were orthorectified using only the image’s rational functions. Secondly, the images have been atmospherically corrected using the ATCOR-2 software suite in ERDAS imagine, which incorporates the MODerate resolution atmospheric TRANsmission (MODTRAN version 4) model and uses look-up tables with pre-calculated model simulations for different satellite sensor types and a range of atmospheric conditions [
51,
52]. Furthermore, supervised classification of Landsat 8 images was performed to identify at each acquisition the plots of agricultural land cover with vegetation and areas with bare soil and/or urban areas.
In order to implement the FAO-56 model, for each satellite image a set of maps of K
c were derived on the basis of LAI and
r maps. The surface albedo, needed to derive the net radiant flux, is an approximation of the hemispherical and spectrally integrated surface albedo, and was extracted as part of the value-added products available in ATCOR-2 after the atmospheric correction [
53]. LAI is a biophysical surface parameter defined as the total one-sided area of photosynthetic tissue per unit of ground area [
54]. The estimate of the LAI was performed based on the simplified Clevers’ Leaf Area Index by Reflectance model, CLAIR [
55], based on the Weighted Difference Vegetation Index (WDVI). The WDVI is a radiometric index calculated, for each pixel of the image, from the values of reflectance r
sRED and r
sNIR, respectively, in the bands of red (0.63–0.69 μm) and the near infrared (0.76–0.90 μm). The ratio of near-infrared and red soil reflectance is also known in the literature as the “soil line slope,” usually between 0.9 and 1.3. The effect of weighting the red band with the slope of the soil line is the maximization of the vegetation signal in the near-infrared band and the minimization of the effect of soil brightness. The LAI is related with WDVI of the observed vegetation, through the following expression, deduced from a simplified analysis of the radiative behavior of different types of crops [
50,
55,
56]:
where α* is an extinction coefficient (increase of WDVI for an unitary increase of LAI—an empirical shape parameter, mainly depending on canopy architecture and computed from field measurement and considered to be 0.39 in this study, corresponding to the minimum error between the observed and estimated LAI), and WDVI∞ is the asymptotical value of WDVI for LAI→∞, and is calculated at each acquisition at the target “vegetation cover (not woody) very dense and uniform” or maximum vegetation (agricultural) cover (usually between 0.55 and 0.75).
In situ LAI measurements have been used to produce LAI maps from Landsat 8 and RapidEye images by applying the semi-empirical CLAIR model based on an inverse relationship between the WDVI and LAI [
55,
57].
The CLAIR model needs to be calibrated by estimating the extinction coefficient α*. Equation (4) must be inverted, fixing the LAI, WDVI, and WDVI∞ values. WDVI∞ was calculated (r
sNIR/r
sRED = 1.35, WDVI∞ = 0.55) so as to create a correspondence of pixels with maximum vegetation cover. The resulted α* values were in the range of 0.34–0.70. The final value of α* = 0.39 corresponds to the minimum error between observed and estimated LAI, leading to an average error of 25% in the estimation of the LAI. To test the model prediction accuracies, we used the
root-mean-square error (
RMSE) and the coefficient of determination (R
2) between measured and predicted LAI. According to the FAO-56 model, the ET
0 computed from the agrometeorological variables using the Penman–Monteith equation is used to calculate the value of K
c (according to the analytical method), where K
c is the ratio between ET
p and ET
0 [
20]. K
c incorporates and synthesizes all the effects on the evapotranspiration related to morpho-physiological characteristics of the different crops—phenological stage, degree of soil cover, soil and climate conditions—that make them different from the reference crop.
The mean value of the ratio between ETp and ET0 over a period of about seven days around the date of satellite’s passage is applied twice, once with standard values of LAI, r, and crop height for a reference surface, and once by using the LAI, r, and crop height (this latter is fixed for tree crops) from satellite images. To speed up this processing in a GIS environment, we derive a polynomial expression of degree 4, whose coefficients are derived by calculating the Kc intervals of LAI and albedo targets (i.e., LAI in range 0–5 and step 0.2; albedo in range 0.05 to 0.5 and step 0.01). The coefficients mostly depend on meteorological data, and they are dimensionless as LAI and Kc. This procedure is not strictly needed, since we could directly use the FAO-56 method on a daily basis, but it has the advantage of showing the direct link between Kc, LAI, and albedo. Once we derive albedo and LAI from the satellite images and calculate the coefficients for a given time interval, in a GIS environment we derive the Kc map on a pixel-by-pixel basis. One important advantage of deriving canopy parameters or crop coefficients from spectral measurements is that their values do not depend on other variables such as planting date and density, but on the effective cover.