Italy - Apulian Tavoliere

Specific Project Objectives & Deliverables


Earth Observation (EO) Data Used:

Sentinel-2 carries an optical payload with visible, near infrared and shortwave infrared sensors comprising 13 spectral bands: 4 bands at 10 m, 6 bands at 20 m and 3 bands at 60 m spatial resolution with a swath width of 290 km.

The 13 spectral bands guarantee consistent time series, showing variability in land surface conditions and minimising any artefacts introduced by atmospheric variability.

The mission orbits at a mean altitude of approximately 800 km and, with the pair of satellites in operation, has a revisit time of five days at the equator (under cloud-free conditions) and 2–3 days at mid-latitudes. 

Sentinel-2 combines a large swath, frequent revisit, and systematic acquisition of all land surfaces at high-spatial resolution and with a large number of spectral bands, all of which makes a unique mission to serve Copernicus. In particular Sentinel-2 incorporates three new spectral bands in the red-edge regions, which are centered at 705, 740 and 783 nm.

The Sentinel-2 Payload Data Ground Segment offers data processing options with bottom-of-atmosphere (BOA) reflectance (Level-2a), through a software toolbox SNAP with Sen2cor tool made available to users, to derive and enhance cloud masks from the top-of-atmosphere reflectance (Level-1c) (2). In this report the Sentinel-2 atmospheric correction was performed based on the algorithm proposed in the Atmospheric/Topographic Correction for Satellite Imagery (3).

Vegetation indices are simple numerical indicators that reduce multispectral (two or more spectral bands) data to a single variable for predicting and assessing vegetation characteristics (1). Vegetation indices (VI) are among the oldest and most widely used tools to provide land-cover maps, land-change detection maps, and plants geophysical variables, such as chlorophyll content per unit leaf area (Ch), leaf area index (LAI) and leaf water content.

In this study the following VIs were calculated:

Normalized Difference Rededge Index  (NDRE) as the difference of the reflectance at near infra-red (NIR) and Rededge  spectral bands normalized by the sum of the reflectance at these spectral bands at 20-m spatial resolution:

Canopy Chlorophyll Content Index (CCCI) was calculated from the Normalized Difference Red Edge index (NDRE) because it is correlated with chlorophyll and leaf N (4):

CCCI was downscaled to 10-m resolution.

AquaCrop Model Simulation

For model application the dataset collected at an experimental field in Capitanata area was used.

The application of the AquaCrop model (Heng et al., 2009) to a given crop requires that a series of inputs is determined. The inputs are related to:

  • Meteorological data: Daily values of maximum, minimum temperature, rain and reference evapotranspiration;
  • Management data: sowing date of 1ft December 2014 was considered; quantity of nitrogen (150 Kg of N ha-1); roots growth data (1.2 m in depth).
  • Physiological parameters of wheat calibrated and validated in the last years were used
  • Soil data: Water content at the field capacity (FC) and wilting point (WP) were measured at four points of the test sites (table 1 and figure 2).


Table 1 - Water content at the field capacity (FC) and wilting point (WP) collected in four points at an experiemntal field in Capitanata.

Figure 2. Four points at an experimental field in Capitanata

At each point daily soil water content was monitored.

For recording soil water status (W), 4 capacitive probes of 0.1 m in length (Decagon devices, 10HS, USA) were installed horizontally in the soil at 25 cm from the soil surface. The probes were linked to a Grillo datalogger (Tecno.El, ITA), which recorded the daily data of soil water content from December 15, 2015 to May 20, 2016.

The soil water status data in the effective root depth were used for  the validation of the AquaCrop  model. In particular, the measured data were compared with the output data from the model. The model performances were evaluated through two statistical tests: 


  1.  the Relative Root Mean Square Error (RRMSE), calculated from the following equation:

        2.The coefficient of determination R² is defined as the squared value of the Pearson correlation coefficient. r² signifies the proportion of the variance in measured data explained by the model. It ranges from 0 to 1, with values close to 1 indicating a good agreement, and typically values greater than 0.5 are considered acceptable in watershed simulations (Moriasi et al., 2007).


Figure 3  shows  the comparison of the soil water content in the effective root depth measured by probes with the one simulated by AquaCrop Model at each point of the experimental field of Capitanata

Figure 3. Values of soil water content in the effective root depth measured and simulated by AquaCrop model at the four points during the 2016 season. The point a, b, c, d are listed in table 1


The values of the statistical tests listed in table 2  show that model simulates adequately the soil water content over daily scale and in more points of the experimental area. In particular RRMSE is between 6 and 13% (excellent and good performance, respectively), and R2 is > 0.86 (figure 4).


Table 2 - Relative Root Mean Square Error (RRMSE in %) and coefficient of determination (R2) determined on wheat grown in Capitanata for daily values of soil water content (W) simulated by AquaCrop  

Figure 4. Correlation between soil water content measured and simulated in each point. The point a, b, c and d are listed in the table 1

The Physiological parameters calibrated in the 2014 season have been used in the simulation model  to predict the final biomass and yield of wheat at each of the 4 points in the field (figure 5 and table 3)

Figure 5. Trends of wheat biomass simulated at each point of the experimental field in Capitanata

Table 3 - Biomass and yield of the four test sites


Figure 5 and Table 3 show that the model was able to simulate a lower production at point C, where the field capacity and wilting point were lower.



In this report, it was not possible to apply the methodology described in the previous reports, to derive reference LAI maps from ground-based measurements over the study area within the JECAM site of Capitanata, owing to a malfunction of Li-COOR sensor.

However, we wanted to investigate the benefit of the Sentinel-2 mission for the agriculture domain, in particular for the rainfed durum wheat grown in a farm of Capitanata located within the JECAM site. The objective was to assess and model the spatial relationship between the yield map and the multi-band images of Sentinel-2 by using multivariate geostatistical techniques. We had planned to analyse multi-temporal images of Sentinel-2, in order to select the one or ones which showed the strongest relationship with yield and then the highest potential for grain-yield prediction. Due to the persistence of wet meteorological conditions during the winter and spring seasons over the study area, only the images of February and May were free enough from clouds. At the time of reporting, only the multi-band image of February was jointly analysed with the yield map.

Materials and Methods

The research was carried out at the Menichella Experimental Farm of CREA-SCA, located in the countryside of Foggia, southern Italy (41° 27’ N, 15° 36’ E, 90 m a.s.l.), within the study area of the JECAM site (Fig. 1) during the wheat season 2015-2016. The trial was conducted on two 5-ha fields cropped with rainfed durum wheat (Triticum durum Desf. Cv “Claudio”), under continuous cultivation. The two fields were differently tilled: the field on the right was uniformly ploughed, whereas the field on the left was submitted to a tillage trial: the experimental design was strip plot with two treatments (ploughed, not ploughed), replicated three times.

The soil, that is very common in Apulian Tavoliere, is silty-clay Vertisol of alluvial origin, classified as Fine Mesic Typic Cromoxerert by Soil Taxonomy USDA (Soil Survey Staff, 1999). The meteorological conditions were mostly those typical of a Mediterranean environment, characterised by dry season between May and September and cold season from October-November to March-April. The field was harvested on July 12th 2016 and grain yield  was normalised at 13% moisture.

Multivariate geostatistical techniques were applied to fuse on-line yield measurements, collected by a John Deer combine machine, with Sentinel-2 sensing data (in particular RGB and NIR data), in order to assess the spatial relationship between crop and remote data and then evaluate the potential of Sentinel 2 data to predict wheat production. The set of geostatistical techniques, to perform sensor data fusion, has been widely described in the previous reports, here we give only a recall of the main procedures.

 An overview of the geostatistical data fusion approach


The main geostatistical procedures applied to fuse the multiple data sets are briefly described below. All geostatistical analyses were performed with the software package ISATIS (Geovariances, 2016). The approach is based on fitting a Linear model of coregionalization (LMC) to the set of both direct and cross- experimental variograms, which considers all the studied variables as the result of the same independent physical processes, acting over different spatial scales.

A difficulty in the practical application of a multivariate approach occurs when the variables are of widely differing sizes. A solution is to standardize the individual variables to give each an average of zero and a variance of unity. Variogram modelling is further complicated by the presence of outliers in highly skewed data distributions. In this case, it is better to perform a normalization of data through Gaussian anamorphosis modelling. Gaussian anamorphosis is a mathematical function, which transforms a variable Y with a Gaussian standardized distribution into a new variable Z with any distribution: Z = ɸ(Y). As this function needs to be known for any Gaussian value, a model is required. This is made by fitting a finite expansion of Hermite polynomials (Chile`s and Delfiner, 1999).

Adopting a Gaussian model, a LMC was fitted to all experimental variograms of the transformed data, and then multicollocated cokriging (Castrignanò et al., 2009) was applied as conditional expectation estimator. Finally, the estimates were back-transformed to the raw values of the variables through the anamorphosis functions previously calculated.


Multi-collocated  cokriging


Ordinary cokriging is one of the most basic cokriging methods, which assumes the local mean to be constant but of unknown value. A way of integrating secondary fine-resolution information in primary sparse variable modelling is collocated cokriging, where the contribution of a secondary variable to the cokriging estimate relies only on the cross- correlation between the two variables (Goovaerts, 1997). The approach is quite similar to ordinary cokriging with the only difference being in the neighbourhood search. Since using all secondary exhaustive information contained within the neighbourhood may lead to an intractable solution, due to too much information, the secondary variable is used only at the target location and also at all the locations where the primary variable is defined within the neighbourhood. This solution has generally produced more reliable and stable results (Castrignano` et al. 2009, 2012). The modified version, also referred to as ‘‘Multi-collocated cokriging’’ in the literature (Rivoirard, 2001), is less precise than full cokriging, and does not use all the auxiliary information contained within the neighbourhood, but it is much less computationally demanding. However, because the co-located secondary datum tends to screen the influence of more distant secondary data, there is actually little loss of information.

The types of spatial data (yield and spectral data) involved in the study are the result of the measurements from different sensors, therefore their integration arises the problem called “the change of support problem” in geostatistics. The term “support” means the size or volume associated with each data value, but it also includes the geometrical size, shape and spatial orientation of the region associated with the measurements (Olea, 1991). The change of support creates a new variable, related to the original one but with different statistical and spatial properties. The crucial problem is to define an approach which relates two or more variables of different supports in a way that permits valid inference. Many of the statistical solutions to this problem were developed by Matheron (1963) in block (co)kriging, whose equations are valid regardless of the support of data. However, taking into account different supports in the calculation of autocovariance and cross-covariance functions is crucial to obtain valid inference. Once a consistent model for the point covariance functions was estimated, using linear model of coregionalization described above, block cokriging was applied by replacing the cross-covariances in point cokriging by their block-averaged counterparts through variogram regularization (Myers, 1984). The regularized variograms were then calculated in different directions using a discretization of the blocks into equal cells, where the blocks were replaced by the union of the cell centers.



In Figure 6 the maps of CCCI and Rededge are shown: at a visual inspection it is evident an indirect spatial relationship between them. This property can be explained since more vigourous plants (greater CCCI) have spectral signatures whith rededge more shifted towards the visible range and then with less reflectivity. The same sort of inverse relationship can be observed also for CCCI and rededge derived from the Sentinel-2 image for the two fields for which the yield maps were available (Fig. 7).

Figure 6. CCCI (a) and rededge (b) maps for the JECAM site of Capitanata

Figure 7. CCCI (a) rededge (b) maps for the two fields within JECAM site


To make more objective the relationship between yield and remote sensing data, we performed a set of statistical and geostatistical procedures. Since most data showed relevant shifts from normal distribution, as showed in the box plots (Fig. 8), all data were previously transformed through Gaussian anamorphosis. A LMC was then fitted to the transformed data, which was then regularised to apply block kriging over the Sentinel-2 grid with mesh of 10 m. The regularised LMC included two basic spatial structures: a spherical model with range of 32 m and a spherical model with range of 300 m (Fig. 9).

Figure 8. Box plot of the variables used in sensor data fusion approach


Figure 9. Regularised LMC in multivariate geostatistical approach


In Fig. 10 there are reported the block-cokriged maps of blue, green , red and NIR bands and the map of yield. From the remote sensing images it is quite manifest a clear differentiation between the two fields, submitted to different tillage and with different agronomical history: the field, uniformly ploughed on the right, has less reflectivity in blue, green, and red and higher reflectivity in NIR, which is consistent with the previous results obtained for CCCI and rededge maps. Moreover, the field on the left shows a characteristic banding due to the different tillage treatments. These differences between the two fields may be due to differences in soil roughness and in crop growth and development. However, in the yield map this differentiation disappears and the pattern of yield variation seems to be affected by different factors.



Figure 10. Block-cogriged maps of Blue, Green, Red (a) and NIR, Yield (b) data



  3. Barnes, E.M., Clarke, T.R., Richards, S.E., Colaizzi, P.D., Haberland, J., Kostrzewski, M., Waller, P., Choi, C., Riley, E., Thompson, T., Lascano, R.J., Li, H., Moran, M.S., 2000. Coincident detection of crop water stress, nitrogen status and canopy density using ground based multispectral data. In: Robert, P.C., Rust, R.H., Larson, W.E. (Eds.), Proc. 5th Int. Conf. Precis Agric. Bloomington, MN, USA.
  4. Castrignanò, A., Costantini, E., Barbetti, R., Sollitto, D., 2009. Accounting for extensive topographic and pedologic secondary information to improve soil mapping. Catena 77: 28-38.
  5. Castrignanò, A., Wong, M.T., Stelluti, M., De Benedetto, D., Sollitto, D., 2012. Use of EMI, gamma-ray emission and GPS height as multi-sensor data for soil characterisation. Geoderma 175:78–89.
  6. Chilès, J.P., Delfiner, P. Geostatistics: Modelling Spatial Uncertainty. Wiley, New York. 695 pp.; 1999.
  7. Géovariances. Isatis Technical Ref., ver. 2016. Geovariances & Ecole Des Mines De Paris. Avon Cedex, France; 2009.
  8. Goovaerts, P., 1997. Geostatistics for natural resources evaluation. Oxford University Press on Demand.
  9. Heng, L.K., Hsiao, T.C., Evett, S.R., Howell, T.A., Steduto, P., 2009. Testing of FAO AquaCrop model for rainfed and irrigated maize. Agron. J., 101, 3, 488–498.
  10. Jamieson, P. D., Porter, J. R., Wilson, D. R., 1991. A test of the computer simulation model ARCWHEAT on wheat crops grown in New Zealand. Fields Crop Research, 27, 337-345
  11. Matheron, G., 1963.Principle of geostatistics. Economic Geology, 58(8): 1246–1266.
  12. Moriasi, D. N., Arnold, J. G., Liew, M. W. V., Bingner, R. L., Harmel, R. D., and Veith, T. L., 2007. Model evaluation guidelines for systematic quantification of accuracy in watershed simulations. Transactions of The ASABE 50, 885–900.
  13. Myers, D. E., 1984. Cokriging: New Developments. Verly G. (1984, ed.), Geostatistics for Natural Resource Characterization, D. Reidel Pub. Co., Dordrecht, 295-305.
  14. Olea, R. A., 1991. Geostatistical Glossary and Multilingual Dictionary. International Association for Mathematical Geology Studies in Mathematical Geology. Oxford University Press, New York. No. 3, 177p.
  15. Richter, R., Schlapfer, D., Muller, A., 2011. Operational atmospheric correction for imaging spectrometers accounting for the smile effect. IEEE Transactions on Geoscience and Remote Sensing, 49 (2011), pp. 1772–1780.
  16. Rivoirard, J., 2001. Which models for collocated cokriging? Math. Geol., 33: 117-131.
©2015 Joint Experiment for Crop Assessment and Monitoring