Next Article in Journal
Electrospinning Nanoparticles-Based Materials Interfaces for Sensor Applications
Next Article in Special Issue
Capturing Maize Stand Heterogeneity Across Yield-Stability Zones Using Unmanned Aerial Vehicles (UAV)
Previous Article in Journal
A Framework for the Joint Placement of Edge Service Infrastructure and User Plane Functions for 5G
Previous Article in Special Issue
Maize Silage Kernel Fragment Estimation Using Deep Learning-Based Object Recognition in Non-Separated Kernel/Stover RGB Images
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Integrating Geophysical and Multispectral Data to Delineate Homogeneous Management Zones within a Vineyard in Northern Italy

1
Department of Agricultural and Environmental Sciences–Production, Landscape, Agroenergy, Università degli Studi di Milano, via Celoria 2, 20133 Milan, Italy
2
Department of Civil and Environmental Engineering, Politecnico di Milano, P.zza Leonardo da Vinci 32, 20133 Milan, Italy
*
Author to whom correspondence should be addressed.
Submission received: 15 July 2019 / Revised: 8 September 2019 / Accepted: 12 September 2019 / Published: 14 September 2019
(This article belongs to the Special Issue Sensors and Systems for Smart Agriculture)

Abstract

:
Soil electrical conductivity (EC) maps obtained through proximal soil sensing (i.e., geophysical data) are usually considered to delineate homogeneous site-specific management zones (SSMZ), used in Precision Agriculture to improve crop production. The recent literature recommends the integration of geophysical soil monitoring data with crop information acquired through multispectral (VIS-NIR) imagery. In non-flat areas, where topography can influence the soil water conditions and consequently the crop water status and the crop yield, considering topography data together with soil and crop data may improve the SSMZ delineation. The objective of this study was the fusion of EC and VIS-NIR data to delineate SSMZs in a rain-fed vineyard located in Northern Italy (Franciacorta), and the assessment of the obtained SSMZ map through the comparison with data acquired by a thermal infrared (TIR) survey carried out during a hot and dry period of the 2017 agricultural season. Data integration is performed by applying multivariate statistical methods (i.e., Principal Component Analysis). The results show that the combined use of soil, topography and crop information improves the SSMZ delineation. Indeed, the correspondence between the SSMZ map and the CWSI map derived from TIR imagery was enhanced by including the NDVI information.

1. Introduction

In agriculture an effective and efficient management of inputs (i.e., water and nutrients) is fundamental to make the crop production sustainable, for both the environment and economics. Knowledge about the plant-soil system is fundamental to achieve this goal. Since the nineties, the Precision Agriculture (PA) approach has required a detailed description of the variability at field scale of soil and plant properties, in order to apply water and nutrients with variable rates, according to the actual irrigation and nutrient requirements, which can be extremely different within a field because of the spatial variability of soil and plant properties [1]. Following this approach, not only the water and nutrient use efficiencies, but also the quantity and the quality of crop yield may improve. Site-specific management of water and nutrients in PA requires the delineation in the field of sub-regions with similar soil and crop characteristics affecting crop yield (Site Specific Management Zone, SSMZ) [2]. Intensive and relatively time-saving measurements of soil electrical conductivity (EC) through geophysical proximal soil sensors are among the most frequently used approaches in PA to delineate SSMZs [3,4,5,6,7]. Statistical procedures [8] are used to classify the EC maps derived by interpolation of the geophysical data, resulting in few sub-field zones to be managed separately. EC maps are considered to delineate SSMZs because EC is influenced by a combination of soil physical-chemical properties affecting crop yield, including soluble salts, clay content and mineralogy, soil water content, bulk density, organic matter, and pH. Finally, a strong relation has been shown to exist between EC and the total available water-holding capacity (AWHC) of soils [9,10,11,12,13,14]. The information on spatial variability of AWHC at field scale obtained from EC measurements may be extremely important to optimize irrigation [9,10,11], allowing the construction of an irrigation prescription map to be used to prescribe variable water amounts (i.e., variable rate irrigation, VRI) according to soil variability.
Ancillary data acquired by spectral sensors in the visible (VIS), near-infrared (NIR) and thermal infrared (TIR) regions shall be combined with EC data to characterize the soil spatial variability [15] and to improve the delineation of SSMZs [16,17]. Indeed, NIR reflectance from soil has been correlated with many soil properties, including total C, total N, water content and texture [18,19]. VIS-NIR and TIR images of the bare soil, acquired by sensors mounted on Unmanned Aerial Vehicles (UAVs), have been used to evaluate the spatial distribution of soil water content [20,21]. The Normalized Difference Vegetation Index (NDVI), calculated as a combination of VIS and NIR data, has been related with soil organic carbon [22,23]. An approach to improve the delineation of SSMZ by considering multi-sensor data describing both crop and soil variability is illustrated in [24]. Multi-sensor data included data acquired through geophysical sensors to describe soil variability and spectral reflectance data derived from satellite imagery to describe crop variability through the construction of vegetation indices maps. The most recent literature stresses the importance to consider multi-sensor data to optimally delineate SSMZs [25,26,27,28,29].
In non-flat areas, as often is the case with vineyards, topography can influence soil water conditions and consequently crop yield and field zonation; in this situation SSMZ delineation may be improved by integrating geophysical soil monitoring data with topography data [12]. The authors of [12] additionally emphasized how spectral vegetation indices (mostly NDVI), constructed from reflectance data acquired through proximal sensing, UAV, aircraft or satellite imagery, are useful in zoning vineyards when the crop behavior is different from year to year because of important interactions between soil and climate influencing the vine vigor.
The main objective of this work is the implementation of data fusion procedures to delineate a SSMZ map in a vineyard of 1.5 ha located in Franciacorta (BS), since this information is crucial for the elaboration of irrigation prescription maps to be used for the design and/or the management of variable-rate irrigation systems. In particular, this study shows how the different types of data which can be involved in the SSMZ delineation are acquired and analyzed, and finally integrated through a data-fusion approach. To assess the reliability of the SSMZ maps obtained from different types of data, each of them was compared with data acquired by a thermal infrared (TIR) survey carried out during a hot and dry period of the 2017 agricultural season, in a phase of the crop phenology during which the vine is normally most sensitive to water stress. Data collected during the TIR survey were used to produce a crop water stress index (CWSI) map, which was demonstrated to be well correlated with the final crop yield [30,31,32,33].

2. Materials and Methods

2.1. Study Area

The experimental site is a rain-fed vineyard of 15,000 m2 located in Franciacorta (Erbusco, 575,813 E, 5,050,828 N, Northern Italy), a rolling hills area south-east of Lake Iseo (Figure 1). Soils at this site are sandy-loam in texture, according to the regional 1:250.000 soil map (http://www.geoportale.regione.lombardia.it). They belong to the land system of ‘intermediate moraine deposits’. In this area, the typic Paleudalf coarse loamy and poorly gravelly soils (CZO1), are associated with more skeletal soils, very deep, with moderately faster permeability and drainage (VBO1) [34].
In the Franciacorta region, the climate is continental, with the lake providing a mitigating effect in both summer and winter. Considering the average monthly values of the main agrometeorological variables registered at the Erbusco station (part of the Lombardy regional monitoring network), located 2 km far from the experimental site, for the period 2008–2018, it can be observed that the minimum and maximum monthly rainfall occur respectively in July (70 mm) and October (85 mm), while the minimum and maximum daily air temperatures vary, respectively, from 7 °C in October to 20 °C in July, and from 16 °C in October to 32 °C in July. Figure 2 shows the behavior of the main agrometeorological variables recorded at the Erbusco station during the experimental period June-August 2017.

2.2. Experimental Surveys

Different types of data were collected in the vineyard to describe all the factors affecting crop yield, related to the hydrological condition of the soil, as well as to the crop vigor and water status. The soil properties were detected through an electro-magnetic induction (EMI) sensor, while the topography of the vineyard and the crop properties were investigated through multispectral and thermal sensors mounted on UAV. Various combinations of these data (i.e., data fusion) were analyzed, to assess their effectiveness in improving the delineation of the SSMZs aimed at optimizing the crop yield (see Section 2.3 for more details).

2.2.1. Soil Survey through EMI Sensors

The soil variability was detected through an EMI survey on 14th June 2017, when the soil water content might be considered close to the field capacity (FC), few days after a three-day period of rainfall that resulted in 26 mm of rain (Figure 2).
The geophysical data were collected with the multi-frequency EMI sensor Profiler EMP-400 (GSSI Inc., Nashua, NH, USA). The EMI sensor worked with up to three different frequencies from 1 to 16 kHz (15 kHz, plus at most two other frequencies), corresponding to decreasing Depths of Exploration (i.e., DoE). Two frequencies were selected for the survey, 15 kHz (DoE about 1.5 m) and 10 kHz (DoE about 2.5 m), to explore the soil in contact with the vineyard’s root system, which is usually 2–3 m deep. The data were acquired along parallel rows with an interdistance of 10 m, while vineyard rows have an interdistance of 2 m. In the portions of the fields characterized by gravelly soils the EMI measurements showed not to be valid, since the EC values were found to be negative. This extreme soil texture was more present in the upper part of the soil profiles (investigated with the sensor operating at higher frequencies), and showed to have a lower weight as the depth increases.

2.2.2. Vegetation and Topography Survey through UAV Multispectral and Thermal Imagery

Vegetation survey was conducted by means of an aerial campaign with sensors mounted on an UAV. The survey took place on 19th July 2017, under sunny and clear blue-sky conditions. The daily average air temperature was 26 °C, with a maximum value of 31 °C during the central hours of the day. The survey was conducted during the veraison phenological stage, in which the crop is more sensitive to crop water stress.
The UAV employed for the survey was the HexaKopter (MikroKopter, Moormerland, Germany). It is a multirotor equipped with six brushless motors, it weighs about 1.2 kg, including batteries, and its maximum transportable payload is equal to 0.5 kg. It can be remotely controlled and programmed for automatic navigation through the free and open source Mission Planner software [35]. Its maximum transmission range is about 200 m and the flight duration is limited to 10 min.
The UAV was equipped with three different sensors, in order to collect imagery in different portions of the electromagnetic spectrum: VIS (450–720 nm), NIR (800–1000 nm) and TIR (7000–14,000 nm). The Survey 2 camera (MAPIR, San Diego, CA, USA) was used for VIS acquisitions, while a modified SJ400 camera (SJCAM, Shangxue Technology Park, Putian, Shenzhen, China) was used to collect NIR imagery. Both instruments are low cost and light-weight cameras, with a CMOS sensor of maximum size 16 Mpx. TIR data were collected by the thermal camera OPTRIS PI400 (Optris GmbH, Berlin, Germany), with a spectral response in the range 7.5–13 μm. The thermal camera acquires data in radiometric video sequences format (.RAVI). For the photogrammetric processing, single frames with resolution equal to 382 × 288 px are subsequently extracted from the video. Technical specifications of the three sensors used for the vegetation survey are reported in Table 1.
According to the UAV payload, two flights were required to collect images with the three sensors. During the first flight, the UAV mounted the Survey2 and the SJ4000 cameras simultaneously, in order to acquire a multispectral dataset. Considering sensors characteristics and study area, flight planning included six strips at an altitude of 60 m above ground level (AGL), with forward and side overlaps equal to 80% and 65%, respectively. Two blocks of data (VIS and NIR) were acquired, each amounting 164 images with ground resolution, namely Ground Sample Distance (GSD), equal to 0.017 m. The adopted plan for the multispectral flight is reported in Figure 3.
During the second flight, the UAV mounted the OPTRIS PI400 to collect data of vegetation temperature. The video sequences were acquired with nadiral orientation at a constant speed of 2.7 m/s and at the altitude fixed to 55 m AGL. The derived images had a GSD of about 0.150 m and forward and side overlaps equal to 80% and 40%, respectively.
The georeferencing and the accuracy of the photogrammetric products were achieved by means of some targets used as Ground Control Points (GCPs), whose center coordinates were measured through a Global Navigation Satellite System (GNSS) receiver. A Leica Viva GS14 GNSS receiver (Leica Geosystems, Heerbrugg, Switzerland) in Network Real Time Kinematic (NRTK) mode was used in this study, with horizontal and vertical accuracies of 2–3 cm and 5 cm, respectively. Different types of targets were used for multispectral and thermal surveys: 16 black and white plastic square panels (30 cm × 30 cm) were employed for the multispectral survey, while 16 polystyrene square panels (60 cm × 60 cm), covered with aluminum foil and marked with a copper cross to enhance the central point were used for the thermal survey. In order to have an optimal distribution of GCPs, the targets were placed both all around the perimeter of the vineyard, on the ground, and inside the investigated area, on the top of the vineyard poles to ensure their visibility. Moreover, some targets with known reflectance and thermal characteristics were imaged, to perform radiometric calibration of VIS-NIR data and atmospheric correction of TIR images. According to [36], four square polystyrene panels (60 cm × 60 cm), covered with plastic, where used for the atmospheric correction of the thermal images. The panels, two white panels (i.e., cold target) and two black panels (i.e., hot target) were placed outside the investigated area in two different positions.

2.3. Methodological Approach for Delineating SSMZs through Data Fusion

Different types of data—geophysical data acquired through EMI sensors, and topographic and crop data acquired through the multispectral VIS-NIR sensors mounted on the UAV—were variously combined to delineate SSMZs. A data fusion approach was considered, by applying multivariate statistical methods (i.e., Principal Component Analysis, PCA) to integrate the different types of data. Precisely, the PCA was applied to the maps elaborated from geophysical and VIS-NIR data as explained in Section 3.1.
Moreover, the CWSI map was elaborated from the imagery acquired through the TIR sensor mounted on UAV. CWSI was calculated as expressed in the following formula:
CWSI = T s T wet T dry T wet
where TS is the crop surface temperature, Twet is the lower boundary of crop temperature corresponding to the water status of a leaf with stomata fully open and a maximum transpiration rate, Tdry is the upper boundary of crop temperature corresponding to the water status of a non-transpiring leaf with stomata completely closed.
The CWSI map was used to assess the effectiveness of the SSMZ delineation obtained from different combination of data even though crop yield maps are usually considered for this purpose. In this study, in absence of this type of information, the CWSI map was used as a proxy of the crop yield map. As a matter of fact, the crop water status (described through the CWSI) is assumed to be the main environmental factor affecting crop yield in this rain-fed vineyard. Areas with a low value in the CWSI map (i.e., good crop water status) were expected to correspond to areas with a high soil water content (i.e., high EC values and/or high NDVI values and/or low topographic slope values).
Particularly, the effectiveness of the data fusion approach to enhance the delineation of SSMZs was assessed by applying the methodology hereinafter explained and illustrated in Figure 4. Two separated areas were defined within the vineyard, because of the occurrence of not valid EMI measurements for gravelly soils (Section 2.2.1). In the first area (called ‘a’), characterized with valid EMI measurements, geophysical and VIS-NIR data were available; in the second area (called ‘b’), characterized with not valid EMI measurements, only VIS-NIR data were available. For each area, maps produced from different combinations of data were fused by applying PCA. Consequently, the SSMZs were elaborated (for each area, ‘a’ and ‘b’) from the integrated maps produced through PCA (i.e., maps of the Principal Componens, PCs) by applying Cluster Analysis (CA) through the Management Zone Analyst (MZA) software [37]. MZA implements an unsupervised fuzzy classification method and determines the optimal number of SSMZs through the minimization of both the indices Normalized Classification Entropy index (NCE) and Fuzziness Performance Index (FPI); the NCE measures the degree of disorganization among zones (the larger the NCE, the higher is the amount of disorganization), the FPI measures the degree of separation between zones (the larger the FPI, the stronger is the membership sharing between zones). Specifically, the CA was applied considering only the PC maps representing most of the variability of the input maps.
For the area ‘a’, the following three cases were analyzed: (1) only geophysical data were considered: the SSMZs were delineated based on the EC maps relative to different soil depths; (2) geophysical data were considered together with the topographic data obtained from VIS-NIR imagery: elevation and slope maps as well as EC maps referred to different soil depths were used to delineate SSMZs; (3) the complete dataset including also crop data was considered: the SSMZs were delineated by integrating the NDVI map with all the previously illustrated maps. NDVI map was elaborated from VIS-NIR imagery, as the index is defined as the normalized difference between NIR and Red bands:
NDVI = NIR Red NIR + Red
For the area ‘b’, the following two cases were analyzed: (1) the topographic data obtained from VIS-NIR imagery were considered: elevation and slope maps were used to delineate SSMZs; (2) the complete dataset, including topographic and crop data, was considered: the SSMZs were delineated by integrating the NDVI map with elevation and slope maps.
For each case, the SSMZ map was validated through a comparison with the CWSI map. The accuracy of the correspondence between SSMZ and CWSI map was analyzed considering the distributions of CWSI values within each SSMZ.

3. Results and Discussion

3.1. Soil, Vegetation and Topography Mapping

3.1.1. EC Maps

The EC measurements obtained for each frequency used with the EMI sensor (15 kHz and 10 kHz) were interpolated on a grid with 2 m pixel size. Two EC maps were obtained (Figure 5), each one relative to a different DoE corresponding approximatively to 1.5 m and 2.5 m, respectively.
The red color area (area ‘b’) in each map of Figure 5 corresponds to gravelly soils, for which the EMI measurements were not valid, because of the very low EC values characterizing those soils. In these zones, negative EC values were obtained from the EMI survey. The total extent of these areas decreases with the increasing DoE.

3.1.2. Topography and Slope Maps

The VIS and NIR imagery blocks were processed through standard photogrammetric workflow [38] with the Agisoft Photoscan Professional software v. 1.2.6 [39]. Finally, the Digital Surface Model (DSM) was produced with a spatial resolution equal to 0.05 m, representing the height model for both vegetation and soil.
In order to reconstruct the soil topography, vegetation pixels were detected and removed from the DSM. Vegetation detection was performed on the DSM by using an algorithm developed by the authors, which assumes that pixels with higher elevation values correspond to vegetation [40]. The algorithm classifies as vegetation pixels all the pixels having a height value greater than a user-defined threshold within a moving window. In a second step, vegetation pixels were subtracted from the DSM, thus producing a model representing the height of the terrain, namely the Digital Terrain Model (DTM) of the study area. Figure 6 shows the DSM and the DTM of the vineyard, obtained after photogrammetric processing and vegetation detection and removal, respectively. The final DTM reported in Figure 6b, obtained after the application of a moving average smoothing filter, has a spatial resolution of 2 m.
Slope and contour line maps were derived from the DTM, by using the Raster Terrain Analysis functions in QGIS v 3.2 [41]. A fixed interval of 1 m was set for the generation of the contour lines, while the slope map was computed as the gradient of the terrain model, having the same spatial resolution of the DTM (i.e., 2 m). Final results are shown in Figure 7.

3.1.3. Vegetation Indices Maps (NDVI and CWSI)

Multispectral VIS-NIR and TIR imagery were used to compute the vegetation indices NDVI and CWSI, commonly adopted to describe vegetation vigor and crop water status, respectively. A VIS-NIR orthomosaic was generated in Digital Number (DN) with a spatial resolution of 0.05 m, then converted in reflectance values, through the radiometric calibration obtained with an empirical line correction approach [42]. Images of the radiometric targets were used to compute the linear regression coefficients of the DN values against the reflectance values from the target surface. The orthomosaic corrected through the radiometric calibration was used to obtain the NDVI map (Equation (2)). Moreover, soil was masked out, to avoid the inclusion of soil pixels in the vegetation maps. The soil mask was created by extracting pixels that were not considered to be vegetation, as described in Section 3.1.2. Figure 8 shows NDVI maps before and after the soil pixel removal. From the graphs showing the frequency distribution of the NDVI maps (Figure 8b,d), it is evident that most of the pixels with NDVI values lower than 0.7-corresponding to soil and inter-row grass—could be removed and only vegetation pixels (NDVI values greater than 0.7) were retained.
The CWSI map was obtained from the TIR orthomosaic. The TIR orthomosaic (with spatial resolution equal to 0.15 m) was generated through a specific procedure for TIR images. This procedure, including single frames extraction, format conversion and photogrammetric processing, is described in detail in [43]. Atmospheric correction of the obtained TIR orthomosaic was performed by using the thermal images of the cold and hot targets, representing respectively the minimum (Tmin) and the maximum (Tmax) temperature values within the investigated area. The temperature of the targets was recorded at the ground level as well. These temperature values, acquired at flight height and at ground level, were used to derive an atmospheric model [36] successively applied to the TIR orthomosaic to obtain the crop surface temperature (called crop surface TIR orthomosaic hereinafter).
The CWSI map was calculated using Equation (1). The values Twet and Tdry were initially determined by an empirical approach reported in many studies [44,45,46,47]. The value Tdry was calculated by using the current Tair plus 5 °K [48,49,50], while the value Twet was calculated as the mean of the coolest 5% vegetated pixels in the crop surface TIR orthomosaic. As for the NDVI map, soil was masked out from the crop surface temperature map. Figure 9 illustrates the CWSI maps before and after soil pixels removal.
Both maps show a zone with values greater than 1, due to the presence of pixels with surface temperature higher than Tdry. This could be due to the fact that, in this study, Tair (31°C) was the average hourly temperature registered during the central hours of the day (from 1:00 p.m. to 2:00 p.m., solar time) at the Erbusco agro-meteorological station, placed 2 km away from the experimental site and positioned over a standard grass surface as indicated by WMO (World Meteorological Organization). In order to estimate a more reliable value for Tdry, the same approach used for Twet was adopted: Tdry was calculated based on the temperature histogram [51,52,53] as the mean of the hottest 5% vegetated pixels in the crop surface TIR orthomosaic. The derived CWSI map, successively used in this study, is shown in Figure 10.

3.2. SSMZ Mapping

Firstly, the SSMZ map was elaborated from EC maps only (Section 3.2.1). Afterwards, topography and crop information were integrated with the EC maps through PCA, to improve the SSMZ map (Section 3.2.2 and Section 3.2.3). The Pearson’s correlation coefficients were computed separately for areas ‘a’ and ‘b’, respectively among the EC, DTM, Slope and NDVI values calculated at the grid nodes used to interpolate the EC data (2 m pixel size), with valid EMI measurements (Table 2), and among the DTM, Slope and NDVI values calculated at the nodes of the same grid, with not valid EMI measurements (Table 3). All the coefficient values are statistically significant with level 0.001 (p-value < 5 × 10−4). Moreover, the Moran Index was calculated to describe the spatial autocorrelation of the variables and the spatial cross-correlation between the variables. The values, calculated separately for areas ‘a’ and ‘b’, considering the nodes respectively inside and outside the vineyard’s area with valid EMI measurements, are reported respectively in Table 4 and Table 5. The Moran Index values were always statistically significant with level 0.001 (p-value < 5 × 10−4). The univariate Moran Index calculated for the different variables was always positive and greater than 0.60, showing a high spatial autocorrelation of all the variables. The bivariate Moran Index (between variables) describes the correlation based on the relationships between each point and the neighboring ones. According to [54], the correlation described through the Pearson’s coefficients can be decomposed in two components, related to a direct correlation without distance effect and an indirect correlation based on the distance effect. Consequently, the difference between the Pearson’s coefficient and the Moran Index quantifies the direct correlation component. For the two datasets described in Table 4 and Table 5, this difference was always less than 0.03 in absolute value, equal to the 25% of the correlation coefficient at most, except for the correlation between the variables EC-15 kHz and EC-10 kHz (Table 4) and between the variables DTM and Slope (Table 5). These results highlighted the relevant contribution of the spatial pattern in cross-correlation.

3.2.1. Delineation of SSMZs from EC Maps

The two EC maps relative to frequencies 15 kHz and 10 kHz, calculated within area ‘a’ considering only the valid EMI measurements (i.e., not negative EC values), were analyzed through the PCA. The SSMZ map (Figure 11) was obtained by applying CA to the first PC, explaining about 96% of the variability of both the EC maps. Three SSMZs were delineated within area ‘a’; another SSMZ (red color) was defined, corresponding to area ‘b’ characterized with not valid EMI measurements (i.e., negative EC values) occurred for gravelly soils (Figure 5). SSMZs from 1 to 4 (Figure 11) correspond to decreasing EC values.
SSMZ and CWSI maps were compared. High EC values (i.e., high soil water contents) were expected to correspond with low CWSI values (i.e., good crop water status). Instead, areas with high CWSI values, denoting crop water stress, were included in the SSMZ 1 (characterized by high EC values), while, vice-versa, areas with low CWSI values were present in the SSMZ 4 (very low EC values). The analysis showed that for the study vineyard, physical-chemical soil properties described by the EC values where not sufficient to explain the crop water status. As matter of fact, the spatial patter of the SSMZs is quite different from that one of the zones in the CWSI map correspondent to low index values (from 0 to 0.5) and high index values (from 0.5 to 1).

3.2.2. Data Fusion: Delineation of SSMZs from Slope and EC Maps

The combined effect of soil properties (i.e., EC values) and field topography (i.e., elevation and slope) was investigated to improve the delineation of SSMZs, looking for a better correspondence with the spatial distribution of the zones in CWSI map with low and high index values. EC, elevation and slope maps were analyzed through PCA. The SSMZ map was obtained by applying CA to the first and second PCs explaining most of the variability of all the considered maps. Particularly, PCA and CA were applied separately to the area ‘a’ with valid EMI measurements, corresponding to SSMZs from 1 to 3 in Figure 11, as well as to the gravelly soil area ‘b’, corresponding to SSMZ 4 in Figure 11.
In the former area, four SSMZs (numbered from 1 to 4 in Figure 12a), were delineated considering EC, elevation and slope maps. In the latter area, three SSMZs (numbered from A to C in Figure 12a) were recognized taking into account only elevation and slope maps. The resulting SSMZ map is shown in Figure 12a. This map, even though improved with respect to that obtained from EC maps only (Figure 11), could not completely explain the spatial variability detected in the CWSI map. Indeed, as illustrated also in Figure 13 showing the distributions of the CWSI values within each SSMZ, the SSMZs 1 and 2 (high EC values) corresponded to low CWSI values (as expected), as well as for the SSMZ A and part of the SSMZ B; on the other hands, SSMZs 3 and 4 included areas with both low and high CWSI values (not expected due to the low EC values), as for the cases of SSMZ C and part of B. This behavior highlighted how the zonation shown in Figure 12a did not consider all the factors affecting the crop water status.

3.2.3. Data Fusion: Delineation of SSMZs from Soil Maps (Slope and EC Maps) and NDVI Map

Finally, also the variability of crop vigor (described by the NDVI map) was taken into account to produce a more reliable SSMZ map, with better correspondence to the zones in the CWSI map characterized by high (from 0.5 to 1) and low values (from 0 to 0.5) of the index. EC, elevation, slope and NDVI maps were analyzed through PCA and CA. Following the same approach considered in the Section 3.2.2, the SSMZ map was obtained by applying PCA and CA firstly to data available within the area ‘a’ with valid EMI measurements, and afterwards to data available within the area ‘b’ characterized by gravelly soils.
For area ‘a’, CA was applied to the first three principal components PCa1, PCa2 and PCa3 (Table 6), obtained from the EC, elevation, slope and NDVI maps: PCa1 represented mainly the physical-chemical soil properties (correlation coefficients with EC greater than 0.90), and partly the DTM (correlation coefficient equal to −0.48); PCa2 represented the topography (correlation coefficients with DTM and Slope equal to 0.72 and −0.45, respectively); PCa3 represented both the topography (correlation coefficient with Slope equal to −0.61) and the crop vigor (correlation coefficient with NDVI equal to 0.59) which are negatively correlated (see Table 2). For area ‘b’, CA was applied to the first two principal components PCb1 and PCb2 (Table 7), obtained from the elevation, slope and NDVI maps: PCb1 represented the topography (correlation coefficients with DTM and Slope equal to −0.89 and 0.88, respectively); PCb2 represented the crop vigor (correlation coefficient with NDVI equal to 0.99).
Within areas ‘a’ and ‘b’, respectively, five SSMZs (numbered from 1 to 5) were delineated from EC, elevation, slope and NDVI maps, and three SSMZs (numbered from A to C) were delineated considering only elevation, slope and NDVI maps. The resulting SSMZ map is shown in Figure 12b. The SSMZs 1–3, A and C corresponded to low CWSI values, while SSMZs 4, 5, and B mostly corresponded to high CWSI values, except for the three small areas highlighted with the red circles in Figure 12b. As matter of fact, the SSMZ delineation in Figure 12b was improved compared to the one shown in Figure 12a, as illustrated in Figure 14: (i) the mean and the standard deviation of the CWSI values within SSMZ 1 decreased, while SSMZs 5 (corresponding to SSMZ 4 in Figure 12a) mostly included high CWSI values, with an increased mean value compared to SSMZ 4 in Figure 12a; (ii) also the mean of the CWSI values within SSMZs 3 increased with respect to the values for SSMZ 1 in Figure 12a; (iii) the SSMZ B was characterized by the CWSI values with the highest mean; (iv) the SSMZ C mostly included low CWSI values (the mean and the standard deviation of the CWSI values within this SSMZ decreased compared to the values for SSMZ C in Figure 12a). Finally, the integration of the NDVI data allowed the delineation of SSMZs each one corresponding respectively to low or high CWSI values (except for the three small areas highlighted with red circles in Figure 12b).
The spatial distributions of the PCs (Figure 15 and Figure 16) explain which factors prevailed in the SSMZ delineation through CA. The SSMZs 1-3 were mainly determined by soil properties (EC data, described by PCa1), while the SSMZ 4, as well as the SSMZs A and B, were mainly determined by topography (DTM and Slope data, described by PCa2 for the case of SSMZ 4, and by PCb1 for the case of SSMZs A and B). The SSMZs 5 and C were mainly determined by crop vigor (NDVI data, described by PCa3 for the case of SSMZ 5, and by PCb2 for the case of SSMZ C).
Moreover, Table 8 shows the correlation between CWSI and the variables (elevation, slope and NDVI) used to integrate the EMI measurements in order to improve the reliability of the SSMZ map obtained from only EC data. Correlation with NDVI was the highest (p-values less than 5 × 10−4), highlighting how NDVI data were able to explain the CWSI spatial variability within the whole field area. As matter of fact, for NDVI the difference between the Pearson’s coefficient and the Moran Index, quantifying the direct correlation component independent from the spatial variability, is almost −0.30, while this component is almost zero for the other variables.

3.3. Discussion

Obtaining a reliable SSMZ map is of practical relevance for farmers, since this map is an important tool to actuate variable rate practices in PA, in term of both designing and managing application systems (e.g., for the water and nutrient management). As matter of fact, the delineated SSMZs are zones where the factors influencing the crop yield (i.e., soil, topography, and micro-climate) result in affecting the crop water status and vigor in a different way. These factors need to be adequately described through thematic maps, to allow the delineation of SSMZs through their combination.
This work proposed a fusion approach integrating different thematic maps (EC, DTM, slope and NDVI maps) to compute a SSMZ map, whose effectiveness was assessed considering a CWSI map. The approach was applied in a rainfed vineyard to obtain a SSMZ map useful for the design and the management of a variable rate irrigation system. By actuating a variable rate irrigation accordingly to the SSMZ map, farmers would achieve a twofold result: first, to obtain a higher and more uniform production and second, to optimize the water use. Interesting general discussion points that emerge from the results of this study are the following:
(1)
the SSMZ map can vary greatly its spatial configuration depending on the information layers used for its production, it is therefore necessary to conduct more research aimed at understanding which information it may be appropriate to include, based not only on the prevailing factors affecting the crop yield, but also on the purpose for which the SSMZ map is being developed (e.g., nutrient management, water management);
(2)
the addition of the topographic information to the soil data included in the EC maps leads the SSMZ map to have a spatial distribution more similar to that shown by the CWSI map; it can be deduced that in a vineyard in slope conditions, the topographic information together with the soil distribution information are able to explain part of the variability illustrated in the vegetation maps;
(3)
in the specific case of this study (i.e., rain-fed vineyard in severe water stress conditions), NDVI data were able to explain the CWSI spatial variability within the whole field area; indeed, the NDVI map showed a strong correlation with the crop water status (CWSI);
(4)
while information on soil properties and topography are not very variable over time, crop data may vary from year to year if soils and topography are not the only factors conditioning their distribution; it would therefore be necessary to repeat the study for several years to verify the ‘stability over time’ of the spatial distribution of crop data;
(5)
if the SSMZ map is to be used to design a ‘rigid’ variable-rate irrigation system (i.e., drip irrigation system subdivided into different irrigation sectors), considering that the geometry of the irrigation system (and consequently the irrigation amounts distributed) cannot be changed from year to year, it is even more important to verify the ‘stability over time’ of the spatial distribution of crop data, and if therefore it makes sense to consider them in the delineation of the SSMZs.

4. Conclusions

Recent literature suggests that integrating EC information with elevation and slope maps, as well as with crop indices describing the crop vigor, can improve the delineation of SSMZs in vineyards. This was demonstrated in this study, focusing on the fusion of EC maps obtained by an EMI geophysical survey and VIS-NIR data collected through UAV-mounted cameras, to optimally delineate SSMZs in a rain-fed vineyard of 15,000 m2 located in Northern Italy. In the study, in absence of a spatially distributed crop yield map, the crop water status detected in a severe dry period during the grape’s veraison phenological stage was assumed to summarize the effect of the principal environmental factors acting on the crop production; consequently, the CWSI map was used as a proxy of the crop yield map itself.
The obtained results stressed how the crop water status in the study vineyard was actually affected significantly not only by the physical-chemical soil properties (described by the EC maps), but also by the elevation and slope of terrain. Moreover, the NDVI map allowed us to include in the analysis time-dependent factors influencing the production (i.e., interaction among soil, topography, micro-climate and vegetation). In fact, for the study vineyard, a good correspondence between the spatial pattern of SSMZ and CWSI maps was achieved only by integrating the NDVI data to the other types of data. Consequently, at least for the study case, a reliable SSMZ map to be used to design and manage irrigation within the vineyard showed to require the availability of both ‘stable-over-time information’, related to soil properties and topography, and ‘time-dependent information’, related to the crop development. In this study, crop information was available only for the 2017 season, but a good practice to obtain reliable SSMZ maps would require the acquisition of crop data during different seasons.

Author Contributions

Conceptualization, B.O., G.S. and A.F.; Data curation, B.O., G.R. and A.M.; Formal analysis, B.O., G.S., G.R. and A.M.; Investigation, B.O., G.S. and G.R.; Methodology, B.O., G.S., G.R. and A.F.; Supervision, B.O., G.S. and A.F.; Visualization, B.O., G.R. and A.M.; Writing–original draft, B.O. and G.R.; Writing–review & editing, B.O., G.S., G.R., A.M. and A.F.

Funding

This research received no external funding.

Acknowledgments

The authors want to acknowledge Ing. Daniele Passoni who planned UAV missions and carried out the flights.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Mulla, D.; Khosla, R. Historical evolution and recent advances in precision farming. In Soil Specific Farming; CRC Press: Boca Raton, FL, USA, 2016; Available online: https://www. taylorfrancis.com/books/97814822 45349 (accessed on 15 November 2018).
  2. Nawar, S.; Corstanje, R.; Halcro, G.; Mulla, D.; Mouazen, A.M. Delineation of soil management zones for variable-rate fertilization: A review. Adv. Agron. 2017, 143, 175–245. [Google Scholar]
  3. Corwin, D.L.; Lesch, S.M.; Shouse, P.J.; Soppe, R.; Ayars, J.E. Identifying soil properties that influence cotton yield using soil sampling directed by apparent soil electrical conductivity. Agron. J. 2003, 95, 352–364. [Google Scholar] [CrossRef]
  4. Morari, F.; Castrignanò, A.; Pagliarin, C. Application of multivariate geostatistics in delineating management zones within a gravelly vineyard using geo-electrical sensors. Comput. Electron. Agric. 2009, 68, 97–107. [Google Scholar] [CrossRef]
  5. Moral, F.J.; Terrón, J.M.; Da Silva, J.M. Delineation of management zones using mobile measurements of soil apparent electrical conductivity and multivariate geostatistical techniques. Soil Tillage Res. 2010, 106, 335–343. [Google Scholar] [CrossRef]
  6. Van Meirvenne, M.; Islam, M.M.; De Smedt, P.; Meerschman, E.; Van De Vijver, E.; Saey, T. Key variables for the identification of soil management classes in the aeolian landscapes of north–west Europe. Geoderma 2013, 199, 99–105. [Google Scholar] [CrossRef]
  7. Neupane, J.; Guo, W. Agronomic Basis and Strategies for Precision Water Management: A Review. Agronomy 2019, 9, 87. [Google Scholar] [CrossRef]
  8. Pascucci, S.; Carfora, M.; Palombo, A.; Pignatti, S.; Casa, R.; Pepe, M.; Castaldi, F. A Comparison between Standard and Functional Clustering Methodologies: Application to Agricultural Fields for Yield Pattern Assessment. Remote Sens. 2018, 10, 585. [Google Scholar] [CrossRef]
  9. Hedley, C.B.; Yule, I.J. A method for spatial prediction of daily soil water status for precise irrigation scheduling. Agric. Water Manage. 2009, 96, 1737–1745. [Google Scholar] [CrossRef]
  10. Hedley, C.B.; Bradbury, S.; Ekanayake, J.; Yule, I.J.; Carrick, S. Spatial irrigation scheduling for variable rate irrigation. In Proceedings of the New Zealand Grassland Association; New Zealand Grassland Association: Lincoln, UK, 2010; Volume 72, pp. 97–102. [Google Scholar]
  11. Hedley, C.B.; Roudier, P.; Yule, I.J.; Ekanayake, J.; Bradbury, S. Soil water status and water table depth modelling using electromagnetic surveys for precision irrigation scheduling. Geoderma 2013, 199, 22–29. [Google Scholar] [CrossRef]
  12. Priori, S.; Martini, E.; Andrenelli, M.C.; Magini, S.; Agnelli, A.E.; Bucelli, P.; Biagi, M.; Pellegrini, S.; Costantini, E.A.C. Improving wine quality through harvest zoning and combined use of remote and soil proximal sensing. Soil Sci. Soc. Am. J. 2013, 77, 1338–1348. [Google Scholar] [CrossRef]
  13. Fortes, R.; Millán, S.; Prieto, M.H.; Campillo, C. A methodology based on apparent electrical conductivity and guided soil samples to improve irrigation zoning. Precis. Agric. 2015, 16, 441–454. [Google Scholar] [CrossRef]
  14. Ortuani, B.; Chiaradia, E.A.; Priori, S.; L’Abate, G.; Canone, D.; Comunian, A.; Giudici, M.; Mele, M.; Facchi, A. Mapping Soil Water Capacity Through EMI Survey to Delineate Site-Specific Management Units Within an Irrigated Field. Soil Sci. 2016, 181, 252–263. [Google Scholar] [CrossRef]
  15. Corwin, D.L. Past, present, and future trends in soil electrical conductivity measurements using geophysical methods. In Handbook of Agricultural Geophysics; CRC Press, Taylor & Francis Group: Boca Raton, FL, USA, 2008; pp. 17–44. [Google Scholar]
  16. López-Lozano, R.; Casterad, M.A.; Herrero, J. Site-specific management units in a commercial maize plot delineated using very high resolution remote sensing and soil properties mapping. Comput. Electron. Agric. 2010, 73, 219–229. [Google Scholar] [CrossRef] [Green Version]
  17. Scudiero, E.; Teatini, P.; Corwin, D.L.; Deiana, R.; Berti, A.; Morari, F. Delineation of site-specific management units in a saline region at the Venice Lagoon margin, Italy, using soil reflectance and apparent electrical conductivity. Comput. Electron. Agric. 2013, 99, 54–64. [Google Scholar] [CrossRef]
  18. Chang, C.W.; Laird, D.A.; Mausbach, M.J.; Hurburgh, C.R. Near-infrared reflectance spectroscopy–principal components regression analyses of soil properties. Soil Sci. Soc. Am. J. 2001, 65, 480–490. [Google Scholar] [CrossRef]
  19. Rossel, R.V.; Walvoort, D.J.J.; McBratney, A.B.; Janik, L.J.; Skjemstad, J.O. Visible, near infrared, mid infrared or combined diffuse reflectance spectroscopy for simultaneous assessment of various soil properties. Geoderma 2006, 131, 59–75. [Google Scholar] [CrossRef]
  20. Maltese, A.; Minacapilli, M.; Cammalleri, C.; Ciraolo, G.; D’Asaro, F. A thermal inertia model for soil water content retrieval using thermal and multispectral images. In Remote Sensing for Agriculture, Ecosystems, and Hydrology XII; International Society for Optics and Photonics: Bellingham, WA, USA, 2010; Volume 7824, p. 78241G. [Google Scholar]
  21. Minacapilli, M.; Cammalleri, C.; Ciraolo, G.; D’Asaro, F.; Iovino, M.; Maltese, A. Thermal inertia modeling for soil surface water content estimation: A laboratory experiment. Soil Sci. Soc. Am. J. 2012, 76, 92–100. [Google Scholar] [CrossRef]
  22. Gomez, C.; Rossel, R.A.V.; McBratney, A.B. Soil organic carbon prediction by hyperspectral remote sensing and field vis-NIR spectroscopy: An Australian case study. Geoderma 2008, 146, 403–411. [Google Scholar] [CrossRef]
  23. Zhang, W.; Wang, K.; Chen, H.; He, X.; Zhang, J. Ancillary information improves kriging on soil organic carbon data for a typical karst peak cluster depression landscape. J. Sci. Food Agric. 2012, 92, 1094–1102. [Google Scholar] [CrossRef]
  24. De Benedetto, D.; Castrignanò, A.; Rinaldi, M.; Ruggieri, S.; Santoro, F.; Figorito, B.; Gualano, S.; Diacono, M.; Tamborrino, R. An approach for delineating homogeneous zones by using multi-sensor data. Geoderma 2013, 199, 117–127. [Google Scholar] [CrossRef]
  25. Shaddad, S.M.; Madrau, S.; Castrignanò, A.; Mouazen, A.M. Data fusion techniques for delineation of site-specific management zones in a field in UK. Precis. Agric. 2016, 17, 200–217. [Google Scholar] [CrossRef]
  26. Castrignanò, A.; Buttafuoco, G.; Quarto, R.; Vitti, C.; Langella, G.; Terribile, F.; Venezia, A. A combined approach of sensor data fusion and multivariate geostatistics for delineation of homogeneous zones in an agricultural field. Sensors 2017, 17, 2794. [Google Scholar] [CrossRef] [PubMed]
  27. Scudiero, E.; Teatini, P.; Manoli, G.; Braga, F.; Skaggs, T.; Morari, F. Workflow to establish time-specific zones in precision agriculture by spatiotemporal integration of plant and soil sensing data. Agronomy 2018, 8, 253. [Google Scholar] [CrossRef]
  28. Anastasiou, E.; Castrignanò, A.; Arvanitis, K.; Fountas, S. A multi-source data fusion approach to assess spatial-temporal variability and delineate homogeneous zones: A use case in a table grape vineyard in Greece. Sci. Total Environ. 2019, 684, 155–163. [Google Scholar] [CrossRef] [PubMed]
  29. Martínez-Casasnovas, J.; Arnó, J. Use of farmer knowledge in the delineation of potential management zones in precision agriculture: a case study in maize (Zea mays L.). Agriculture 2018, 8, 84. [Google Scholar] [CrossRef]
  30. Maes, W.H.; Steppe, K. Estimating evapotranspiration and drought stress with ground-based thermal remote sensing in agriculture: A review. J. Exp. Bot. 2012, 63, 4671–4712. [Google Scholar] [CrossRef] [PubMed]
  31. Orta, A.H.; Baser, I.; Sehirali, S.; Erdem, T.; Erdem, Y. Use of infrared thermometry for developing baseline equations and scheduling irrigation in wheat. Cereal Res. Commun. 2004, 32, 363–370. [Google Scholar]
  32. Ajayi, A.E.; Olufayo, A.A. Evaluation of two temperature stress indices to estimate grain sorghum yield and evapotranspiration. Agron. J. 2004, 96, 1282–1287. [Google Scholar] [CrossRef]
  33. Erdem, Y.; Erdem, T.; ORTA, A.H.; Okursoy, H. Irrigation scheduling for watermelon with crop water stress index (CWSI). J. Cent. Eur. Agric. 2005, 6, 449–460. [Google Scholar]
  34. USDA. Soil Taxonomy, A Basic System of Soil Classification for Making and Interpreting Soil Surveys; USDA: Pittsburgh, PA, USA, 1999.
  35. Mission Planner Home. Available online: http://ardupilot.org/planner/index.html (accessed on 10 July 2017).
  36. Labbé, S.; Lebourgeois, V.; Jolivot, A.; Marti, R. Thermal infra-red remote sensing for water stress estimation in agriculture. Options Méditerranéennes Série B. Etudes et Recherches 2012, 67, 175–184. [Google Scholar]
  37. Fridgen, J.J.; Kitchen, N.R.; Sudduth, K.A.; Drummond, S.T.; Wiebold, W.J.; Fraisse, C.W. Management zone analyst (MZA). Agron. J. 2004, 96, 100–108. [Google Scholar] [CrossRef]
  38. Ronchetti, G.; Pagliari, D.; Sona, G. DTM generation through UAV survey with a fisheye camera on a vineyard. Int. Arch. Photogramm. Remote Sens. Spatial Inf. Sci. 2018, 42, 983–989. [Google Scholar] [CrossRef]
  39. Agisoft PhotoScan Professional, V.1.2.6. Available online: http://www.agisoft.com/ (accessed on 16 February 2018).
  40. Geophysical Research Abstracts Vol. 21, EGU2019-13801, EGU General Assembly 2019. Available online: https://meetingorganizer.copernicus.org/EGU2019/EGU2019-13801.pdf (accessed on 13 September 2019).
  41. Open Source Geospatial Foundation. QGIS Geographic Information System. Available online: http://qgis.osgeo.org (accessed on 28 May 2019).
  42. Smith, G.M.; Milton, E.J. The use of the empirical line method to calibrate remotely sensed data to reflectance. Int. J. Remote Sens. 1999, 20, 2653–2662. [Google Scholar] [CrossRef]
  43. Tucci, G.; Parisi, E.I.; Castelli, G.; Errico, A.; Corongiu, M.; Sona, G.; Viviani, E.; Bresci, E.; Preti, F. Multi-Sensor UAV Application for Thermal Analysis on a Dry-Stone Terraced Vineyard in Rural Tuscany Landscape. ISPRS Int. J. Geoinf. 2019, 8, 87. [Google Scholar] [CrossRef]
  44. Meron, M.; Tsipris, J.; Charitt, D. Remote mapping of crop water status to assess spatial variability of crop stress. In Proceedings of the Fourth European Conference on Precision Agriculture, Berlin, Germany, 15–19 June 2003; pp. 405–410. [Google Scholar]
  45. Cohen, Y.; Alchanatis, V.; Meron, M.; Saranga, Y.; Tsipris, J. Estimation of leaf water potential by thermal imagery and spatial analysis. J. Exp. Bot. 2005, 56, 1843–1852. [Google Scholar] [CrossRef]
  46. Gonzalez-Dugo, V.; Zarco-Tejada, P.; Nicolás, E.; Nortes, P.A.; Alarcón, J.J.; Intrigliolo, D.S.; Fereres, E.J.P.A. Using high resolution UAV thermal imagery to assess the variability in the water status of five fruit tree species within a commercial orchard. Precis. Agric. 2013, 14, 660–678. [Google Scholar] [CrossRef]
  47. Gerhards, M.; Schlerf, M.; Rascher, U.; Udelhoven, T.; Juszczak, R.; Alberti, G.; Inoue, Y. Analysis of airborne optical and thermal imagery for detection of water stress symptoms. Remote Sens. 2018, 10, 1139. [Google Scholar] [CrossRef]
  48. Möller, M.; Alchanatis, V.; Cohen, Y.; Meron, M.; Tsipris, J.; Naor, A.; Ostrovsky, V.; Sprintsin., M.; Cohen, S. Use of thermal and visible imagery for estimating crop water status of irrigated grapevine. J. Exp. Bot. 2006, 58, 827–838. [Google Scholar] [Green Version]
  49. Ben-Gal, A.; Agam, N.; Alchanatis, V.; Cohen, Y.; Yermiyahu, U.; Zipori, I.; Presnov, E.V.; Sprintsin, M.; Dag, A. Evaluating water stress in irrigated olives: correlation of soil water status, tree water status, and thermal imagery. Irrig. Sci. 2009, 27, 367–376. [Google Scholar] [CrossRef]
  50. Alchanatis, V.; Cohen, Y.; Cohen, S.; Moller, M.; Sprinstin, M.; Meron, M.; Tsipris, J.; Saranga, Y.; Sela, E. Evaluation of different approaches for estimating and mapping crop water status in cotton with thermal imaging. Precis. Agric. 2010, 11, 27–41. [Google Scholar] [CrossRef]
  51. Park, S.; Ryu, D.; Fuentes, S.; Chung, H.; Hernández-Montes, E.; O’Connell, M. Adaptive estimation of crop water stress in nectarine and peach orchards using high-resolution imagery from an unmanned aerial vehicle (UAV). Remote Sens. 2017, 9, 828. [Google Scholar] [CrossRef]
  52. Veysi, S.; Naseri, A.A.; Hamzeh, S.; Bartholomeus, H. A satellite based crop water stress index for irrigation scheduling in sugarcane fields. Agric. Water Manage. 2017, 189, 70–86. [Google Scholar] [CrossRef]
  53. Bian, J.; Zhang, Z.; Chen, J.; Chen, H.; Cui, C.; Li, X.; Chen, S.; Fu, Q. Simplified Evaluation of Cotton Water Stress Using High Resolution Unmanned Aerial Vehicle Thermal Imagery. Remote Sens. 2019, 11, 267. [Google Scholar] [CrossRef]
  54. Chen, Y. A new methodology of spatial cross-correlation analysis. PLoS ONE 2015, 10, e0126158. [Google Scholar] [CrossRef] [PubMed]
Figure 1. The experimental site; Coordinate Reference System (CRS): WGS84/UTM zone 32 N. Map data: ©OpenStreetMap contributors.
Figure 1. The experimental site; Coordinate Reference System (CRS): WGS84/UTM zone 32 N. Map data: ©OpenStreetMap contributors.
Sensors 19 03974 g001
Figure 2. Precipitation and temperature daily data collected at the agrometeorological station of Erbusco, during the experimental period from June to August 2017.
Figure 2. Precipitation and temperature daily data collected at the agrometeorological station of Erbusco, during the experimental period from June to August 2017.
Sensors 19 03974 g002
Figure 3. Flight track-lines for multispectral images acquisition and Ground Control Points (GCPs) distribution. Map data: ©Google Satellite.
Figure 3. Flight track-lines for multispectral images acquisition and Ground Control Points (GCPs) distribution. Map data: ©Google Satellite.
Sensors 19 03974 g003
Figure 4. Scheme of the methodological approach adopted in this study.
Figure 4. Scheme of the methodological approach adopted in this study.
Sensors 19 03974 g004
Figure 5. The obtained EC maps (mS/m): (a) frequency 15 kHz, corresponding to a DoE of 1.5 m; (b) frequency 10 kHz, corresponding to a DoE of 2.5 m. Red color area (area ‘b’) corresponds to gravelly soils (EMI measurement not valid).
Figure 5. The obtained EC maps (mS/m): (a) frequency 15 kHz, corresponding to a DoE of 1.5 m; (b) frequency 10 kHz, corresponding to a DoE of 2.5 m. Red color area (area ‘b’) corresponds to gravelly soils (EMI measurement not valid).
Sensors 19 03974 g005
Figure 6. DSM (a) and DTM (b) produced through photogrammetric processing of multispectral (VIS-NIR) dataset.
Figure 6. DSM (a) and DTM (b) produced through photogrammetric processing of multispectral (VIS-NIR) dataset.
Sensors 19 03974 g006
Figure 7. Slope map and contour lines derived from the DTM.
Figure 7. Slope map and contour lines derived from the DTM.
Sensors 19 03974 g007
Figure 8. NDVI map before (a) and after (c) soil masking, together with their respective frequency distribution (b,d).
Figure 8. NDVI map before (a) and after (c) soil masking, together with their respective frequency distribution (b,d).
Sensors 19 03974 g008
Figure 9. CWSI map before (a) and after (c) soil masking. The frequency distribution of the crop surface temperatures, with the illustration of Twet and Tdry values calculated according to the empirical approach described in Section 3.1.3, is reported for each case (b,d).
Figure 9. CWSI map before (a) and after (c) soil masking. The frequency distribution of the crop surface temperatures, with the illustration of Twet and Tdry values calculated according to the empirical approach described in Section 3.1.3, is reported for each case (b,d).
Sensors 19 03974 g009
Figure 10. CWSI map after soil masking (a), derived considering the Twet and Tdry values calculated as the mean of the coolest 5% and the hottest 5% vegetated pixels in the crop surface TIR orthomosaic, respectively. The frequency distribution of the crop surface temperatures is also reported (b).
Figure 10. CWSI map after soil masking (a), derived considering the Twet and Tdry values calculated as the mean of the coolest 5% and the hottest 5% vegetated pixels in the crop surface TIR orthomosaic, respectively. The frequency distribution of the crop surface temperatures is also reported (b).
Sensors 19 03974 g010
Figure 11. The SSMZ map obtained from the EC maps relative to frequencies 15 kHz and 10 kHz. SSMZ from 1 to 4 corresponds to decreasing EC values; in particular, SSMZ 4 corresponds to negative EC values, due to gravelly soils.
Figure 11. The SSMZ map obtained from the EC maps relative to frequencies 15 kHz and 10 kHz. SSMZ from 1 to 4 corresponds to decreasing EC values; in particular, SSMZ 4 corresponds to negative EC values, due to gravelly soils.
Sensors 19 03974 g011
Figure 12. The SSMZ map obtained from: (a) EC, elevation and slope maps; (b) EC, elevation, slope and NDVI.
Figure 12. The SSMZ map obtained from: (a) EC, elevation and slope maps; (b) EC, elevation, slope and NDVI.
Sensors 19 03974 g012
Figure 13. Distribution of the CWSI values within each SSMZ shown in Figure 12a: (a) SSMZ 1; (b) SSMZ 2; (c) SSMZ 3; (d) SSMZ 4; (e) SSMZ A; (f) SSMZ B; (g) SSMZ C.
Figure 13. Distribution of the CWSI values within each SSMZ shown in Figure 12a: (a) SSMZ 1; (b) SSMZ 2; (c) SSMZ 3; (d) SSMZ 4; (e) SSMZ A; (f) SSMZ B; (g) SSMZ C.
Sensors 19 03974 g013
Figure 14. Distribution of the CWSI values within each SSMZ shown in Figure 12a: (a) SSMZ 1; (b) SSMZ 2; (c) SSMZ 3; (d) SSMZ 4; (e) SSMZ 5; (f) SSMZ A; (g) SSMZ B; (h) SSMZ C.
Figure 14. Distribution of the CWSI values within each SSMZ shown in Figure 12a: (a) SSMZ 1; (b) SSMZ 2; (c) SSMZ 3; (d) SSMZ 4; (e) SSMZ 5; (f) SSMZ A; (g) SSMZ B; (h) SSMZ C.
Sensors 19 03974 g014aSensors 19 03974 g014b
Figure 15. Results of the PCA applied in the area ‘a’: spatial distribution of PCa1 (a), PCa2 (b) and PCa3 (c).
Figure 15. Results of the PCA applied in the area ‘a’: spatial distribution of PCa1 (a), PCa2 (b) and PCa3 (c).
Sensors 19 03974 g015
Figure 16. Results of the PCA applied in the area ‘b’: spatial distribution of PCb1 (a), PCb2 (b).
Figure 16. Results of the PCA applied in the area ‘b’: spatial distribution of PCb1 (a), PCb2 (b).
Sensors 19 03974 g016
Table 1. Technical specifications of the three cameras used for the vegetation survey.
Table 1. Technical specifications of the three cameras used for the vegetation survey.
Survey 2SJ4000OPTRIS PI400
AcquisitionVISNIRTIR
Focal length (mm)4.354.358
Sensor size (mm)4.86 × 3.644.86 × 3.649.55 × 7.2
Sensor size (px)4032 × 30244032 × 3024382 × 288
Pixel size (μm)1.21.225
Field of View (FOV)82°82°62° × 49°
Output formatJPEG imageJPEG imageRAVI video
Weight (g)6464380
Table 2. Pearson’s correlation coefficients among the variables used to delineate SSMZ, estimated considering the grid nodes with valid EMI measurements (area ‘a’).
Table 2. Pearson’s correlation coefficients among the variables used to delineate SSMZ, estimated considering the grid nodes with valid EMI measurements (area ‘a’).
EC-15 kHzEC-10 kHzDTMSlopeNDVI
EC-15 kHz10.93 ***−0.20 ***0.23 ***−0.22 ***
EC-10 kHz 1−0.29 ***0.26 ***−0.21 ***
DTM 1−0.40 ***−0.19 ***
Slope 1−0.07 **
NDVI 1
*** p-value < 5 × 10−5; ** p-value < 5 × 10−4.
Table 3. Pearson’s correlation coefficients among the variables used to delineate SSMZ, estimated considering the grid nodes with not valid EMI measurements (area ‘b’).
Table 3. Pearson’s correlation coefficients among the variables used to delineate SSMZ, estimated considering the grid nodes with not valid EMI measurements (area ‘b’).
EC-15 kHzEC-10 kHzDTMSlopeNDVI
EC-15 kHz-----
EC-10 kHz ----
DTM 1−0.40 ***−0.12 **
Slope 1−0.08 *
NDVI 1
*** p-value < 5 × 10−5; ** p-value < 5 × 10−4; * p-value < 5 × 10−3.
Table 4. Moran Index among the variables used to delineate SSMZ, estimated (using GeoDa software, by Luc Anselin) considering the grid nodes with valid EMI measurements (area ‘a’).
Table 4. Moran Index among the variables used to delineate SSMZ, estimated (using GeoDa software, by Luc Anselin) considering the grid nodes with valid EMI measurements (area ‘a’).
EC-15 kHzEC-10 kHzDTMSlopeNDVI
EC-15 kHz0.82 **0.79 **−0.20 **0.23 **−0.24 **
EC-10 kHz 0.81 **−0.28 **0.27 **−0.22 **
DTM 0.99 **−0.40 **−0.19 **
Slope 0.63 **−0.07 **
NDVI 0.76 **
*** p-value < 5 × 10−5; ** p-value < 5 × 10−4.
Table 5. Moran Index among the variables used to delineate SSMZ, estimated (using GeoDa software, by Luc Anselin) considering the grid nodes with not valid EMI measurements (area ‘b’).
Table 5. Moran Index among the variables used to delineate SSMZ, estimated (using GeoDa software, by Luc Anselin) considering the grid nodes with not valid EMI measurements (area ‘b’).
EC-15 kHzEC-10 kHzDTMSlopeNDVI
EC-15 kHz-----
EC-10 kHz ----
DTM 0.96 **−0.59 **−0.15 **
Slope 0.69 **−0.09 **
NDVI 0.66 **
*** p-value < 5 × 10−5; ** p-value < 5 × 10−4; * p-value < 5 × 10−3.
Table 6. Results of PCA applied in the area ‘a’: variance of the principal components considered in CA and correlation coefficients with the variables used to delineate SSMZ.
Table 6. Results of PCA applied in the area ‘a’: variance of the principal components considered in CA and correlation coefficients with the variables used to delineate SSMZ.
VarianceCumulative VarianceEC-15 kHzEC-10 kHzDTMSlopeNDVI
PCa12.2645%0.900.93−0.480.52−0.27
PCa21.2971%0.250.180.72−0.45−0.69
PCa30.8788%0.280.260.03−0.610.59
Table 7. Results of PCA applied in the area ‘b’: variance of the principal components considered in CA and correlation coefficients with the variables used to delineate SSMZ.
Table 7. Results of PCA applied in the area ‘b’: variance of the principal components considered in CA and correlation coefficients with the variables used to delineate SSMZ.
VarianceCumulative VarianceDTMSlopeNDVI
PCb11.5652%−0.890.880.06
PCb21.0387%−0.14−0.200.99
Table 8. Correlation between CWSI and variables DTM, Slope and NDVI.
Table 8. Correlation between CWSI and variables DTM, Slope and NDVI.
DTMSlopeNDVI
Pearson’s coefficient0.29 ***0.02 *−0.71 ***
Moran Index0.29 **0.02 **−0.52 **
*** p-value < 5 × 10−5; ** p-value < 5 × 10−4; * p-value > 0.1.

Share and Cite

MDPI and ACS Style

Ortuani, B.; Sona, G.; Ronchetti, G.; Mayer, A.; Facchi, A. Integrating Geophysical and Multispectral Data to Delineate Homogeneous Management Zones within a Vineyard in Northern Italy. Sensors 2019, 19, 3974. https://0-doi-org.brum.beds.ac.uk/10.3390/s19183974

AMA Style

Ortuani B, Sona G, Ronchetti G, Mayer A, Facchi A. Integrating Geophysical and Multispectral Data to Delineate Homogeneous Management Zones within a Vineyard in Northern Italy. Sensors. 2019; 19(18):3974. https://0-doi-org.brum.beds.ac.uk/10.3390/s19183974

Chicago/Turabian Style

Ortuani, Bianca, Giovanna Sona, Giulia Ronchetti, Alice Mayer, and Arianna Facchi. 2019. "Integrating Geophysical and Multispectral Data to Delineate Homogeneous Management Zones within a Vineyard in Northern Italy" Sensors 19, no. 18: 3974. https://0-doi-org.brum.beds.ac.uk/10.3390/s19183974

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop