Multi-Parametric Approach to Management Zone Delineation in a Hazelnut Grove in Italy

: The increase in high-density hazelnut ( Corylus avellana ) areas drives the interest in practices of precision management. This work addressed soil apparent electrical conductivity (ECa), RGB aerial (UAV) images, proximal sensing, and ﬁeld scouting in delineating and validating management zones (MZs) in a 2.96 ha hazelnut grove in Italy. ECa data were ﬁtted to a semi-variogram, interpolated (simple kriging), and clustered, resulting in two MZs that were subjected to soil analysis. RGB imagery was used to extract tree canopies from the soil background and determine two vegetation indices (VIs) of general crop status: the Visible Atmospherically Resistant Index ( VARI ) and the Normalized Green-Red Difference Index ( NGRDI ). Then, plant growth parameters were manually assessed (tree height, crown size, etc.) and a proximal VI, the Canopy Index (CI), was determined with the MECS-VINE ® vertical multisensor. MZ1 was characterized by lower ECa values than MZ2. This was associated with a lower clay content (9% vs. 21% in MZ1 vs. MZ2) and organic matter content (1.03% vs. 1.51% in MZ1 vs. MZ2), indicating lower soil fertility in MZ1 vs. MZ2. Additionally, hazelnut trees had signiﬁcantly smaller canopies (1.42 vs. 1.94 m 2 tree − 1 ) and slightly lower values of VARI , NGRDI , and CI in MZ1 vs. MZ2. In conclusion, our approach used ECa to identify homogeneous ﬁeld areas, which showed differences in soil properties inﬂuencing tree growth. This is the premise for differential hazelnut management in view of better efﬁciency and sustainability in the use of crop inputs.


Introduction
Italy is the second largest hazelnut producer in the world, and the Viterbo province (Central Italy) ranks first in hazelnut production with 48,400 tons [1,2]. Currently, hazelnut groves cover 43% of this province's agricultural area; they are mainly organised as monocultures with an increasing trend towards high-density large plantations [2]. Therefore, Italy plays a relevant role in the supply of high-quality hazelnuts, the demand for which is soaring under the needs of the confectionery industry [2]. This drive toward intensive production systems is coupled with increased adoption of irrigation and mechanical harvesting [3].
From this perspective, hazelnut intensification could benefit from precise, site-specific management in terms of crop input efficiency, as well as environmental quality. However, proof of hazelnut suitability for site-specific management is still lacking in the literature. Therefore, the adoption of these practices is based on uncertain grounds at present. At the same time, progress in hazelnut cultivation is expected to come from advanced tools and strategies, such as wireless sensor networks, aerial and terrestrial unmanned vehicles, image analysis, and big data management [4][5][6].
Soil apparent electrical conductivity (ECa) via on-the-go sensors is considered a useful and easily applied system for the rapid prospection of soil properties affecting crop productivity, while saving in labour-intensive and expensive soil analysis. Based on this, ECa has an active role in site-specific soil management [7,8]. Soil mapping with ECa and crop sensing with spectral vegetation indices (VIs) have mainly addressed field crops, with fewer applications in orchards and tree groves. However, differences in soil properties can affect tree growth, fruit yield, and quality as well. Few studies have addressed orchards/groves, and even fewer under Mediterranean conditions. Käthner and Zude-Sasse [9] reported significant correlations between ECa, fruiting, and final yield in a plum orchard. This supports the quality of ECa-based field zoning as a master variable in a common fruit tree species, which was, nevertheless, very different from hazelnut trees.
The validation of ECa-based homogeneous management zones (MZs) requires validation by means of other georeferenced information layers, including soil parameters, plant growth traits, and canopy parameters from remote and proximal instrumental sensing [10]. In tree groves, analysis of the data obtained from remote sensing is hindered by the presence of the soil background, requiring it to be separated from the canopy. In fact, the canopy projection area is important for determining tree growth, water requirements during the growing season, and yield potential [11]. This issue is not given sufficient attention, as demonstrated by a recent review on tree crop management [12], where no study addressing this issue was reported for hazelnuts.
The Excess Green (ExG) VI, extracted from unmanned aerial vehicle (UAV) RGB images, has been proven to perform well in separating plants from non-plant portions in both low and tall species [13,14]. Moreover, canopy reflectance using RGB images, i.e., in the visible spectrum, has been successfully used to generate canopy Vis, including the Visible Atmospheric Resistant Index (VARI) and the Normalized Red-Green Difference Index (NGRDI), which are used to determine the aerial biomass and nutrient status of various crop species [15,16]. It has been, therefore, proposed that canopy reflectance in the visible spectrum could simultaneously perform different tasks: canopy extraction from the surrounding area and proper canopy sensing for plant growth status.
In addition to aerial platforms, several ground-based systems designed for tree crop sensing have appeared on the market. A specific multi-parametric sensor of this kind, the MECS-VINE ® (TEAM Group, Piacenza, Italy), was used to characterize plant growth and environmental conditions in vineyards and olive groves, where displayed a valid pattern of important canopy features and was a good predictor of yield potential [17,18]. It has been envisaged that such instruments could contribute to a reduction of the uncertainties in management options and also in less investigated crops, such as hazelnuts.
Under this premise, this study addressed the use ECa as a master variable for establishing MZs in an intensive hazelnut grove. The aims were: (1) to characterize the spatial variability of soil ECa and establish homogeneous zones potentially suited for differential management, and (2) to use proximal RGB imagery, field scouting, and the MECS-VINE ® -derived Canopy Index to validate the previously established MZs. The final object was probing the reliability of ECa-based MZs under multiple viewpoints in a representative case study placed in a vacated area. This was considered the premise for the subsequent step of differential crop input management according to specific plant needs and yield potential.

Materials and Methods
The study was carried out on a 2.96 ha surface subjected to uniform management involving rainfed cropping. The surface was part of a 7.3 ha hazelnut grove located near Viterbo in the Lazio Region, Italy (42.377034 • N, 12.018766 • E, 192 m asl) ( Figure 1). The grove was five years old with a 5 × 5 m plantation pattern, using the Giffoni cultivar for 90% of the plants intermixed with Tonda Gentile Romana and Nocchione cultivars as pollinators, constituting the residual 10% of plants. In the soil map of the Lazio Region,  According to Metzger et al. [19], the study area belongs to the Mediterranean North environmental zone. This zone features: (i) average minima and maxima temperatures of 8.2 and 18.1 °C throughout the year, with limited freezing in the winter time; (ii) an average precipitation of 478 and 218 mm in October-April and May-September, respectively.

Soil Parameters
A survey was conducted on 26 March 2022 with the proximal ECa sensor Top Soil Mapper (Geoprospectors ® GmbH, Traiskirchen, Austria), which measures ECa at four soil depths, such as 0-0.25 m (R1), 0-0.4 m (R2), 0-0.6 m (R3), and 0-0.8 m (R4), and provides a relative water total content (RWTC) index that defines the spatial distribution of soil moisture through an algorithm calculating volumetric water content based on electrical conductivity. Data were georeferenced using the GNSS receiver STONEX ® S8 PLUS (Stonex Srl, Monza, Italy) and a real-time kinematic (RTK) correction system with centimetre-level accuracy. ECa values above or below ±2.5 standard deviation (SD) were considered outliers and removed from the original dataset.
The ECa empirical semi-variograms were modelled in ArcGIS geo-statistical module (ArcGIS 10.3 ® ESRI, Redlands, CA, USA) to determine the spatial dependence (SD) in ECa data by means of optimized semi-variogram parameters: nugget (Co), sill (Co + C), and range. Spatial dependence (SD) was calculated as the nugget-to-sill ratio for determining the spatial structure in the ECa dataset; the degree of spatial dependence (DSD) (weak, medium, strong) was interpreted according to Cambardella et al. [20]. Prediction accuracy According to Metzger et al. [19], the study area belongs to the Mediterranean North environmental zone. This zone features: (i) average minima and maxima temperatures of 8.2 and 18.1 • C throughout the year, with limited freezing in the winter time; (ii) an average precipitation of 478 and 218 mm in October-April and May-September, respectively.

Soil Parameters
A survey was conducted on 26 March 2022 with the proximal ECa sensor Top Soil Mapper (Geoprospectors ® GmbH, Traiskirchen, Austria), which measures ECa at four soil depths, such as 0-0.25 m (R1), 0-0.4 m (R2), 0-0.6 m (R3), and 0-0.8 m (R4), and provides a relative water total content (RWTC) index that defines the spatial distribution of soil moisture through an algorithm calculating volumetric water content based on electrical conductivity. Data were georeferenced using the GNSS receiver STONEX ® S8 PLUS (Stonex Srl, Monza, Italy) and a real-time kinematic (RTK) correction system with centimetre-level accuracy. ECa values above or below ±2.5 standard deviation (SD) were considered outliers and removed from the original dataset.
The ECa empirical semi-variograms were modelled in ArcGIS geo-statistical module (ArcGIS 10.3 ® ESRI, Redlands, CA, USA) to determine the spatial dependence (SD) in ECa data by means of optimized semi-variogram parameters: nugget (Co), sill (Co + C), and range. Spatial dependence (SD) was calculated as the nugget-to-sill ratio for determining the spatial structure in the ECa dataset; the degree of spatial dependence (DSD) (weak, medium, strong) was interpreted according to Cambardella et al. [20]. Prediction accuracy was determined in terms of root mean square standardised (RMSS) error. Based on this, the best-fit semi-variogram model was determined for ECa data (two soil depths, plus ECa-derived soil moisture parameters), as it was considered a master parameter for the identification of homogeneous zones, and the ordinary kriging method was used for smoothing the ECa data at a 1 m grid size using the best fitting semi-variogram model.
Thereafter, homogenous sub-field areas were developed through the k-Means Clustering method based on ECa interpolated data. The pixels were grouped into no more than two different management zones (MZs) to reduce the fragmentation of the small field surface [21]. A composite soil sample was taken in each MZ to describe the average soil properties (soil texture, pH, organic carbon, total nitrogen, C:N ratio), reflecting ECa values.

Remote Sensing
On 14 July 2022, ultra-high resolution RGB imagery was collected by the UAV DJI Mavic Pro equipped with a DJI FC220 camera with a 12 Megapixels resolution. The flight was performed at 30 m above ground. The digital camera captured three wavebands, centred at 450 nm (blue), 550 nm (green), and 750 nm (red). The spatial resolution of the image was 0.01 m.
To discriminate the hazelnut canopy projection area from the vegetated ground cover, the Excess Green (ExG) VI was calculated [13,22], as follows (Equation (1)): where R* RED , GREEN, BLUE were normalized (0-1) RGB values. The ExG tonal raster image was binarized, creating a good contrast between the plant canopy and the background, thus enabling an easy separation [23]. In this process, the disturbance of lateral tree crowns is minimized by overlaying the processed image obtained from the RGB camera with the polygonal canopy buffer centred on the trees.
The tonalised regions were then assumed to represent the plant canopy projection areas of hazelnut trees [24]. A point vector layer with a grid corresponding to the planting layout was overlapped on the raster layer, so that each point corresponded to a tree position centroid. The centroid positions of trees placed in an irregular pattern were corrected by manually forcing each point of the grid to correspond to the actual position of each tree. A buffer layer with a 1.5 m radius was then created around each centroid so that each circle circumscribed the canopy of each tree. Thus, a multi-polygon vector layer was obtained, the centroids of which represented tree positions and the pixels inside the polygons defined the tree canopy cover. Finally, each tree in the plot was identified as an individual object, while the pixels outside were identified as the background. The number of pixels that fell within each buffer estimated the tree size.
Based on the same RGB drone images, two additional VIs of general crop status were calculated: the Visible Atmospherically Resistant Index (VARI) (Equation (2)), and the Normalized Green-Red Difference Index (NGRDI) (Equation (3)). In both cases, data were normalised to a 0-1 range [25,26].

Proximal Sensing
Proximal, on-the-go sensing of tree canopy parameters was performed with the Micro Environment and Canopy Sensor MECS-VINE ® (TEAM Group, Piacenza, Italy). This is a multiparametric sensor integrating a GPS receiver with sensors for canopy status, air temperature and relative humidity, surface temperature, and distance from the canopy wall [17,18]. The canopy status derives from the processing of the information provided by two RGB optical matrix imaging sensors operating on the right and left sides of the travelling tractor. The system's internal algorithm, Canopyct, linearly combines RGB bands, image brightness, and hue information. The collected data are filtered, and the resulting Canopy Index (CI) consists of a pure number (range, 1-1000) indicating canopy density and healthiness [17]. This multi-sensor performs a segmented evaluation of the canopy up to 16 sectors differentiated in height so as to cover the entire tree profile.
The survey was performed on 10 November 2022 mounting the multi-sensor behind a tractor travelling at a forward speed of 5.5 km h −1 . The optical sensors were oriented to intercept, with their cone of vision, the two vegetative walls of the trees (Figure 2). The data returned by the MECS-VINE ® consisted of 16 + 1 raster layers, a 1 × 1 m grid size, representing the distribution maps of the values assumed by the CI on 16 different height sectors, plus a layer with the cumulative value. Based on the thematic raster map returned by the device and the associated legend divided into 10 CI intervals, it was possible to associate each pixel's CI to the corresponding class value in the 16 + 1 raster maps.
Proximal, on-the-go sensing of tree canopy parameters was performed with the Micro Environment and Canopy Sensor MECS-VINE ® (TEAM Group, Piacenza, Italy). This is a multiparametric sensor integrating a GPS receiver with sensors for canopy status, air temperature and relative humidity, surface temperature, and distance from the canopy wall [17,18]. The canopy status derives from the processing of the information provided by two RGB optical matrix imaging sensors operating on the right and left sides of the travelling tractor. The system's internal algorithm, Canopyct, linearly combines RGB bands, image brightness, and hue information. The collected data are filtered, and the resulting Canopy Index (CI) consists of a pure number (range, 1-1000) indicating canopy density and healthiness [17]. This multi-sensor performs a segmented evaluation of the canopy up to 16 sectors differentiated in height so as to cover the entire tree profile.
The survey was performed on 10 November 2022 mounting the multi-sensor behind a tractor travelling at a forward speed of 5.5 km h −1 . The optical sensors were oriented to intercept, with their cone of vision, the two vegetative walls of the trees (Figure 2). The data returned by the MECS-VINE ® consisted of 16 + 1 raster layers, a 1 × 1 m grid size, representing the distribution maps of the values assumed by the CI on 16 different height sectors, plus a layer with the cumulative value. Based on the thematic raster map returned by the device and the associated legend divided into 10 CI intervals, it was possible to associate each pixel's CI to the corresponding class value in the 16 + 1 raster maps.

Field Scouting
On 11 November 2022, a field survey was carried out on 84 randomly selected plants, identifying each plant with the corresponding code of the multi-polygon canopy map to geo-reference the parameters of the surveyed plants. For each tree, the following traits were measured: plant height, two orthogonal crown diameters, average stump diameter, and the average diameter of the three most vigorous shoots.

Data Analysis
Data from remote sensing, proximal sensing, and field scouting were considered auxiliary parameters useful for validating the differences between the two MZs established by ECa as the master parameter.
Based on the WGS84/UTM zone 33N as the coordinate reference system, a 1 × 1 m grid was created and used to intersect the raster images of the following layers: Vis (RGB imagery), CI (MECS-VINE ® ), and ECa-derived MZs. The multi-polygon canopy layer was used to intersect the raster data, which led to splitting each previous raster layer into two layers: one for the entire field and the other for the background. The multi-polygon canopy layer was then input into the QGIS zonal statistics module to extract the data (min, max, mean) of the two VIs and 16 MECS-VINE ® raster layers, to associate each tree with

Field Scouting
On 11 November 2022, a field survey was carried out on 84 randomly selected plants, identifying each plant with the corresponding code of the multi-polygon canopy map to geo-reference the parameters of the surveyed plants. For each tree, the following traits were measured: plant height, two orthogonal crown diameters, average stump diameter, and the average diameter of the three most vigorous shoots.

Data Analysis
Data from remote sensing, proximal sensing, and field scouting were considered auxiliary parameters useful for validating the differences between the two MZs established by ECa as the master parameter.
Based on the WGS84/UTM zone 33N as the coordinate reference system, a 1 × 1 m grid was created and used to intersect the raster images of the following layers: Vis (RGB imagery), CI (MECS-VINE ® ), and ECa-derived MZs. The multi-polygon canopy layer was used to intersect the raster data, which led to splitting each previous raster layer into two layers: one for the entire field and the other for the background. The multi-polygon canopy layer was then input into the QGIS zonal statistics module to extract the data (min, max, mean) of the two VIs and 16 MECS-VINE ® raster layers, to associate each tree with its values. The data were stored in a point and a multi-polygon layer where the trees were represented by their centroids.
Data from UAV imagery (n = 582 and 428 in MZ1 and MZ2, respectively), MECS-VINE ® sensing (n = 474 and 378 in the respective MZs), and field scouting (n = 41 and 43 in the respective MZs) were subjected to the two-tailed t-test at p ≤ 0.05 and p ≤ 0.01.

ECa Data Interpolation
By evaluating the parameters (nugget, sill, range) of the principal semi-variogram models applied to ECa data, it was found that the exponential model was always the best Sustainability 2023, 15, 10106 6 of 12 performer, namely in terms of the high sill and range values. Therefore, the exponential model was always chosen to represent ECa properties (Table 1). Additionally, this model always exhibited zero nuggets, i.e., all of the sill was due to the variable component, meaning that the spatial dependence was always very high, and the resulting DSD was always strong. The three ECa properties and their spatial variability maps are presented in terms of semi-variogram parameters (Table 1) and kriged maps (Figure 3). The ECa values at two soil depths and the ECa-based RWTC expressing soil moisture status modestly varied in terms of their semi-variogram parameters. Additionally, the exponential semivariogram could minimize the error in data fitting, as shown by the low residual values (RMSS) originating from the respective equations. in the respective MZs) were subjected to the two-tailed t-test at p ≤ 0.05 and p ≤ 0.01.

ECa Data Interpolation
By evaluating the parameters (nugget, sill, range) of the principal semi-variogram models applied to ECa data, it was found that the exponential model was always the best performer, namely in terms of the high sill and range values. Therefore, the exponential model was always chosen to represent ECa properties (Table 1). Additionally, this model always exhibited zero nuggets, i.e., all of the sill was due to the variable component, meaning that the spatial dependence was always very high, and the resulting DSD was always strong. The three ECa properties and their spatial variability maps are presented in terms of semi-variogram parameters (Table 1) and kriged maps (Figure 3). The ECa values at two soil depths and the ECa-based RWTC expressing soil moisture status modestly varied in terms of their semi-variogram parameters. Additionally, the exponential semi-variogram could minimize the error in data fitting, as shown by the low residual values (RMSS) originating from the respective equations.

Management Zone Delineation
Upon performing the ordinary kriging of the ECa data with the exponential model, clustering sub-divided the 2.96 ha field into two MZs: MZ1 (1.65 ha) and MZ2 (1.31 ha) (Figure 4), featuring lower/higher ECa values, respectively ( Table 2). Upon performing the ordinary kriging of the ECa data with the exponential mode clustering sub-divided the 2.96 ha field into two MZs: MZ1 (1.65 ha) and MZ2 (1.31 ha (Figure 4), featuring lower/higher ECa values, respectively (Table 2).  The soil samples taken in the two separate MZs exhibited compositional differences which were consistent with their ECa values. MZ1 vs. MZ2 had a coarser texture (sand loam vs. sandy clay loam), neutral vs. slightly acidic pH, lower organic C and total N contents, and similar C-to-N ratios ( Table 2). It was evidenced that MZ2 featured a mor fertile soil thanks to its greater ability to maintain soil moisture (higher clay content) unde rainfed cultivation, and more favourable biological properties (higher organic C content

Canopy Projection Area
In separating the tree canopy projection from the soil background, a low discriminat ing threshold (0.1) was used in the ExG assessment, thanks to the clear contrast provided by the VI between plants and soil. This, in addition to inserting the canopy delimitatio buffer (Figure 5), enabled the removal of any residual background images linked to th presence of grass.
The image made it possible to discriminate all of the 1010 trees in the field and iden tify an average canopy area of 1.65 m 2 obtained from the sum of the pixels correspondin to the horizontal canopy projection excluding the empty pixels. The soil samples taken in the two separate MZs exhibited compositional differences, which were consistent with their ECa values. MZ1 vs. MZ2 had a coarser texture (sandy loam vs. sandy clay loam), neutral vs. slightly acidic pH, lower organic C and total N contents, and similar C-to-N ratios ( Table 2). It was evidenced that MZ2 featured a more fertile soil thanks to its greater ability to maintain soil moisture (higher clay content) under rainfed cultivation, and more favourable biological properties (higher organic C content).

Canopy Projection Area
In separating the tree canopy projection from the soil background, a low discriminating threshold (0.1) was used in the ExG assessment, thanks to the clear contrast provided by the VI between plants and soil. This, in addition to inserting the canopy delimitation buffer (Figure 5), enabled the removal of any residual background images linked to the presence of grass.

Remote Canopy Vegetation Indices and Field Scouting
The two VIs assessed from UAV at peak time during the growing season outlined small differences between MZ2 and MZ1 (1.3% and 8.6% for the respective VARI and NGRDI) (Table 3). However, these differences were proven to be highly significant by the statistical analysis, indicating a very limited variation for both VIs around mean values in both MZs (the coefficient of variation was always below 1.0% and 4.0% for the respective VARI and NGRDI). The image made it possible to discriminate all of the 1010 trees in the field and identify an average canopy area of 1.65 m 2 obtained from the sum of the pixels corresponding to the horizontal canopy projection excluding the empty pixels.

Remote Canopy Vegetation Indices and Field Scouting
The two VIs assessed from UAV at peak time during the growing season outlined small differences between MZ2 and MZ1 (1.3% and 8.6% for the respective VARI and NGRDI) (Table 3). However, these differences were proven to be highly significant by the statistical analysis, indicating a very limited variation for both VIs around mean values in both MZs (the coefficient of variation was always below 1.0% and 4.0% for the respective VARI and NGRDI). Compared to this, the differences in the five growth traits between MZ2 and MZ1 were more relevant (Table 3): 26%, 123%, 34%, 23%, and 23% for tree height, crown area, crown diameter, stump diameter, and most vigorous shoots diameter, respectively. It should be noted that all traits originated from field scouting, except the crown area, which was obtained from remote sensing and subsequent pixel extraction.

MECS-VINE ® Proximal Sensing
The CIs shown by the two MZs described a slope from the base to the top sector ( Figure 6). Regarding the average of the 16 sectors (not shown in Figure 6), MZ2 staged a significantly higher CI value with respect to MZ1 (309 vs. 290). However, three canopy layers could be identified, where the two MZs exhibited different behaviour in terms of CI data. The two bottom canopy sectors (1 and 2) outlined insignificant/modest CI differences between the two MZs. CI differences between MZs increased to approximately 60-80 points in the lower-mid canopy portion (sectors 3 to 7), which featured intrinsically high CI values. Lastly, CI differences faded in the top canopy portion (sectors 8 to 16), where intrinsically low CI values were generally shown. Regarding the average of the 16 sectors (not shown in Figure 6), MZ2 staged a significantly higher CI value with respect to MZ1 (309 vs. 290). However, three canopy layers could be identified, where the two MZs exhibited different behaviour in terms of CI data. The two bottom canopy sectors (1 and 2) outlined insignificant/modest CI differences between the two MZs. CI differences between MZs increased to approximately 60-80 points in the lower-mid canopy portion (sectors 3 to 7), which featured intrinsically high CI values.

Discussion
Lastly, CI differences faded in the top canopy portion (sectors 8 to 16), where intrinsically low CI values were generally shown.

Discussion
Hazelnuts are becoming increasingly more strategic for industries processing these nuts into products due to increasing market demand. To keep pace with this, new hazelnut orchards have been designed with higher planting densities with single-trunk training systems to support mechanization; dripping or sub-surface irrigation systems have been increasingly adopted [11]. Precision management technologies represent agronomic strategies, the effects of which are expected to positively influence hazelnut production efficiency as well as environmental sustainability [27]. This also concerns water management, owing to the fact that water is often a limiting factor in the Mediterranean region [28]. Therefore, it comes as no surprise that the principles of site-specific management are gaining the interest of a growing community focusing on this crop.
This work contributes to filling the gap of knowledge in the specific sector, demonstrating how spatial variability can be determined in hazelnut groves. The result was the establishment of appropriate MZs as a preliminary step for the differential application of crop inputs. In fact, the hazelnut plant, despite its specific growth habit and general suitability for varying ambient conditions, was shown to be influenced by different soil characteristics, which, in turn, may be demonstrated by a rapid prospection as that provided by ECa. This saved a large number of soil samples, which alternatively would have needed to be collected and analysed in view of proper interpolation and subsequent MZ establishment [29]. Compared to this, the cheaper and quicker ECa mapping was shown to be a reliable method for establishing different MZs in the surveyed case study.
As an additional requisite, delineation was aimed at obtaining a reasonably low number of MZs with respect to the given area: in this experiment, a three-hectare surface was sub-divided into two MZs of similar surfaces, and their respective distributions were almost continuous in terms of space ( Figure 4). This is the premise for reducing field patchiness in view of simpler crop management according to a site-specific approach [8].
Another delicate issue addressed in this work is tree canopy separation from the vegetable undercover in the soil background by means of aerial (UAV) images. This work did not try to characterize the canopy volume as other works have [11,27], but succeeded in assessing the size of the canopy area projected onto the ground, which represents a very important parameter for characterizing hazelnut growth and validating the reliability of the ECa-based MZs. It is worth noting that this result was obtained in a very simple way using low-cost RGB cameras in lieu of the more expensive multispectral sensors used to generate an NDVI map [16]. There are doubts pertaining to the risk that increasing the canopy size affects the quality of image segregation with respect to the background and the surrounding plants. However, we believe that, as proposed in this study, by overlaying the image obtained from the RGB camera with the polygonal buffer of the canopy centred on the trees, the disturbance of the lateral tree crowns should be limited.
The methodology adopted was based on two steps: first, cutting out round buffers around plant centroids to reduce background RGB reflectance around core areas for tree canopy detection; second, discriminating the actual tree canopy within each buffer by means of the ExG index, for which purpose a high discriminating threshold was not required. In practice, this double step may be considered sounder than the simple ExG index. In fact, the double step was proven to be successful in discriminating between vegetation layers. The non-negligible side effect was that common aerial images in the visible spectrum, i.e., from an unprofessional camera, were proven to be suitable for this specific task [23].
The two MZs based on ECa values exhibited not only different soil characteristics (Table 2), on average, but also different tree attributes both at ground level as well as from the aerial survey (Table 3). It was perceived that crown size (i.e., crown area and diameter) was the combined attribute most sensitive to soil fertility differences between the two MZs.
It must be noted that the variation in the crown area between MZ2 and MZ1 (+123%) was much larger than in the crown diameter (+34%). The crown area directly influences crop canopy and, therefore, light interception. Therefore, higher soil fertility originated from higher clay and organic C content [30] in MZ2 vs. MZ1, reflected in the plant trait that is the main driver for productivity, and this trait can be smartly sensed using an RGB camera mounted on a drone.
Optical sensing of leaf walls from ground level with the MECS-VINE ® multi-layer sensor performed equally well in discriminating between the two MZs ( Figure 6). The trend of declining CI values from the bottom to the top leaf layers was consistent with a previous application of the specific sensor on a vineyard, the arboreal plantation for which this system was initially devised [17]. However, the test in an olive grove evidenced a time trend of CI values across the growing season, which was more relevant than the CI variation in the vertical distribution of tree vegetation [18]. It is perceived, therefore, that new instrumentation, which is being made increasingly available for use on various crops, requires specific calibration to obtain a univocal interpretation of the data retrieved. Despite this potential limitation, it may not pass unnoticed that ground-based, optically sensed traits also revealed differential characteristics in the crown structure between the two MZs in our experiment.

Conclusions
Hazelnut groves are experiencing a process of crop intensification driven by rising market demand. This is made possible by intrinsic plant characteristics, namely the dry nut product, enabling relatively safe handling compared to fresh fruit products. From this perspective, techniques pertaining to the specific management of field portions featuring intrinsic soil differences are of increasing interest. However, the development of site-specific management is still at an early stage, and doubts persist regarding the means to address the spatial variability in hazelnut groves.
This work evidences the intrinsic soil and plant differences between two areas in a relatively small area (2.96 ha), and epitomizes the tools and approaches that may be adopted to assess them. In this respect, soil apparent electrical conductivity (ECa) was shown to be a good proxy for clay and organic matter content, which are important soil fertility parameters, especially in the case of rainfed grove management. Once delineated, two different management zones, tree canopy sizes, and two spectral vegetation indices were obtained from aerial RGB images, and tree growth parameters from field scouting and the vertical canopy structure from proximal sensing were consistent with ECa values in the two management zones.
Overall, it was demonstrated that more fertile soil conditions (higher clay and organic matter content) in one of the two management zones resulted in a larger tree crown area, which is the premise for accrued light interception and final yield. It remains necessary to address, in future studies, whether in the less fertile zone, a stronger support with external inputs (fertilisers, irrigation, etc.) could bridge the gap with the more fertile zone.