The novel use of proximal photogrammetry and terrestrial LiDAR to quantify the structural complexity of orchard trees

Within the agrifood sector, the production of high yields is a driver for UK orchard husbandry. Currently, orchard tree management is typically a non-discriminatory method with all trees subjected to the same interventions. Previous studies indicate that structural complexity of individual orchard trees is an indicator for future yield, which can guide the management of individual trees. However, data on the structure of individual trees is often limited. This study investigated the suitability of using remote sensing methods to capture data that can be used to quantify tree structure. Descriptive metrics based on the mathematical assessment of self-affinity and dimensionality were applied to the remotely-sensed data to quantify tree structure, and were also analysed for suitability as a predictor of fruit yield. The findings suggest that while proximal photogrammetry is informative, terrestrial LiDAR data can be used to quantify structural complexity most effectively and this approach holds greater potential for informing orchard management.


Introduction
Fruit orchards are a high value element of the UK agrifood production network and make a significant contribution to the gross value-added sector of the UK economy (DEFRA et al. 2013;Gongal et al. 2016). The top fruit sector increased production between 1997 and 2016, with orchard output increasing 27% to 300 000t annually, and apple production accounting for 87% of this (DEFRA, 2017). As demand for production increases, legislative pressure, such as the Pesticide Adjustment to the Crop Environment (Walklate and 1 3 Cross 2014), and commercial pressures are driving the need for new innovation within the sector.
For effective management of orchard crops, precise knowledge of the size, shape and spatial distribution of the canopy allow decisions to be made about operational interventions (Gongal et al. 2016). In forestry, the structural assessment of tree crowns has been shown to be a reliable proxy of tree condition (Zhang et al. 2007;Murray et al. 2018). In agriculture, quantifying individual tree structure is important both for informing pruning regimes and for the precision application of irrigation, fertilisation and other inputs considered necessary to increase fruit yield (Lauri et al. 2004). Previous studies have shown that spatial measurements of canopy structure may be used to estimate the yield in singlevariety experiments (Sinoquet and Rivet 1997;Gongal et al. 2016).
Traditional arboricultural techniques for estimating tree structure condition are labourintensive and require the use of skilled personnel. The utilisation of modern remote sensing (RS) techniques has huge potential to reduce both the time required to make such measurements, as well as minimising the subjective elements of manual tree assessments. Studies have shown that terrestrial light detection and ranging (LiDAR) data may be used to measure canopy size, density and orientation (Llorens et al. 2011;Beyer et al. 2017), while recent studies in viticulture have demonstrated the potential of laser scanning to characterise the extent of vine vigour as an indicator of yield (Tagarakis et al. 2017). However, researchers have noted that the overlapping nature of tree crown structures can make LiDAR measurements of dense orchards problematic, due to the trees in the foreground obscuring the path of the laser and shadowing the trees behind (Olschofsky et al. 2016). A study by Arnó et al. (2017) investigated the use of tractor mounted (mobile) terrestrial laser scanning (TLS) in a sampling scheme along linear vine trellis. This study found that inter-row distances of less than 0.5 m from the TLS were unsuitable for scanning and analysis. However, the co-registration of LiDAR scans (the merging of several different LiDAR viewpoints) may be used to overcome issues of shadowing (Olschofsky et al. 2016) and the use of ground measurement targets allows the increase in point density and optimisation of the co-registered 3D point cloud (Van der Zande et al. 2006). Following the capture of scanner data, analysis using RS metrics typically rely upon the quantification of entire woodland or forest canopies, where canopy gap fraction, canopy openness or the calculation of leaf area index (LAI) are used as a proxy of tree condition Alonzo et al. 2015;Rautiainen and Stenberg 2015). LAI can be calculated on an individual tree level from aerial laser scanning (ALS) however, the method as described by Lin and West (2016) requires a combined ALS campaign with supporting TLS validation at high resolution, making this approach a prohibitively complex undertaking for many applications. Jonckheere et al. (2004) suggested that hemispherical imagery can be used to overcome erroneous data collection issues frequently associated with LAI analyses. Furthermore, Liang et al. (2014) describe that for individual tree crown investigations, off-the-shelf digital single-lens reflex cameras can provide equivalent data to laser scanning, including when fitted with a hemispherical lens and used to image the canopy from beneath. This technique allows the capture of subject detail not available from other perspectives (Talbot et al. 2017;Murray et al. 2018). Processing these images allows the quantification of tree structure through measuring the fractal dimension (D f ) of the image (Fig. 1 a, b), which as Mandelbrot (1982) describes, provides opportunities to quantify the complexities of patterns within natural structures, particularly in tree crowns which are formed by both the genetic predisposition and environmental influences of "chance, irregularities and non-uniformity". The D f approach has been used as a reliable measure of tree structure complexity in mature trees (Murray et al. 2018). Other measures of dimensionality, such as correlation dimension (D c ) are also considered to be effective at quantifying dimensionality with low numbers of data points (Zhang et al. 2007;Theiler 1990).
In this study, two techniques were compared for the measurement of canopy structure (photogrammetric and LiDAR) in a research orchard. The canopy structure measurements were related to historical fruiting yield and ground reference data for each tree. Data processing pipelines for both the 2D (photogrammetric) and 3D (LiDAR) datasets were developed to quantify measures of dimensionality on a per-tree basis.

Study site
The apple trees were grown in a research plot in Kent, southeast England (51°17′09.8″N 0°27′09.5″E), and replicated a production orchard layout. The local topography was predominantly level, with even exposure to solar radiation and localised wind forces on all sides. The trees were established on a 3 × 3 m grid, totalling 140 tree locations in seven rows of 20 stems, and were approximately 5-8 years of age. All trees had dwarf species provenance, with ground reference (GR) tree heights measured between 1.4 to 3.3 m. Fig. 1 The visualisation of a single tree measured using two different field techniques: a hemispherical imagery taken looking towards the zenith viewpoint, b following processing to remove unwanted image artefacts and conversion to a binary form where the occupied (white) parts of the image are analysed using a box-counting function, and, c a LiDAR point cloud acquired from terrestrial laser scanning, classified by range from the cloud centre Bordering the entire compartment was a shelterbelt of apple trees, one row deep, which were omitted from the analysis. All trees were staked to 2400 × 100 mm round posts, driven approximately 600 mm into the ground, secured with standard tree ties with individual identification tags for rapid identification. The trees were a mixture of fruiting scions (Malus domestica 'Braeburn', M. domestica 'Gala') grafted onto 12 different varieties of M. domestica root-stocks.

Ground reference (GR) measurements
Using traditional arboricultural techniques, direct, manual measurements of each tree were taken for GR data and to enable later RS data validation. The study site (Fig. 2) was visited during the winter of 2017/18 while the trees were in leaf-off condition, which allowed direct access to the tree structure. The total tree height from the ground to the crowns highest point, the cardinal crown spread (measured from the centre stem to the maximum extent of the four cardinal points), and stem girth at 150 mm above the graft line, were measured. The trees were classified into four condition categories; poor, moderate, good and vigorous, using the visual tree assessment (VTA) and traditional arboricultural assessment (TAA) methods (Mattheck and Breloer 1994;Murray et al. 2018). VTA and TAA consider whether the potential of the tree biology, biomechanics or the physiological function have been compromised by any biological (e.g. fungi, disease), or mechanical (e.g. structural cracks, wind stress) factors, using a systematic, visual field method. The study orchard is located at an agricultural research facility where seasonal formative pruning and records of fruiting yield (kg of fruit) were acquired for each tree. Due to measuring the tree structure in winter, these yield records from the previous season were used as the current yield output for each tree in the study. Prior to the fieldwork and scanning, the trees were last pruned in December 2015.

Proximal photogrammetry
Using the proximal photogrammetric techniques as described in Murray et al. (2018), hemispherical imagery was taken from beneath each tree crown, looking upwards towards the zenith viewpoint (Fig. 1a). The photography equipment used was a highresolution digital single-lens reflex (dSLR) camera (EOS 550D DS126271, Canon, Oita, Japan) with an 18 mm lens with hemispherical lens adaptor (Super Wide AF Fisheye Lens 0.20×, Opteka, Yokohama, Japan). The lens adaptor converted the standard focal length into a 3.6 mm circular lens. Images taken with the hemispherical lens permitted the photographing of the maximum amount of tree crown within each image as possible. In Murray et al. (2018), the camera was placed on a standard photographic tripod and was leveled ~ 500 mm from the ground. Due to the low canopy of the dwarf trees in this study, the body of the camera was leveled flat on the ground with the hemispherical lens ~ 200 mm above the ground surface. The distance between the camera and the lowest branches was ~ 500 mm above the lens. At each tree, images were captured at the midpoint of the south western crown axis in order to capture as much of the crown as possible in the image (Murray et al. 2018). All images were subject to post-processing corrections, removing unwanted image artefacts such as blurring or chromatic aberration (CA). CA is the mis-registration of the image red, green and blue (RGB) channels and interference with the image sensor which could disrupt the classification of image pixels (Schwalbe et al. 2009). CA errors were removed by converting the image to the G channel only for analysis. Image blurring was removed using an un-sharp masking algorithm, where a blurred negative copy of the image was removed, leaving a visually sharper image. Finally, barrel distortion of the image was removed using a distortion correction algorithm which transformed the barrel shaped image edge, to be similar to that of an image taken at normal focal length (de Vries 2012).
The photogrammetric D f analysis was undertaken as described in Moisy (2008), and Murray et al. (2018), and is described in Eq. (1) as: where N is the number of boxes required to cover the occupied elements of the binarised shape, R is the box unit size, and N(R) is the box frequency needed to fulfil the image region. Lim is the limit of R, as R approaches infinity (Bonnet et al. 2001;Moisy 2008). This approach is considered a 2D photogrammetric analysis, as the images were taken of a 3D structure and presented as a 2D image. Each image was assessed for self-affinity by the calculation of the logarithmic mean of the D f in each image (Murray et al. 2018). The D f is calculated using a box-counting function, which derives a local D f over the pixel occupancy of each binarised image (Moisy 2008). This generates a fractal region for each image and, as Murray et al. (2018) describe, this quantifies the structural complexity of the tree crown. (

Terrestrial laser scanning
Before scanning, five retro-reflectors (100 x 100 mm cylinder reflectors, individually leveled on tribrachs to increase accuracy) were placed on standard surveyor tripods within the orchard area, and used as tie-points for creating a multi-perspective, co-registered scan. A Riegl LMS Z210ii scanner (Riegl Laser Management Systems GmbH, Horn, Austria) was used to scan the trees, which were reconstructed as a 3D point cloud representative of the tree structures (Fig. 1c). The LiDAR scans were undertaken from two separate positions, orientated perpendicularly to the survey plot and the point clouds combined in Riegl RiScan-Pro software to increase the point cloud density (Fig. 2). All other analyses were completed using the Python 3 programming language (Rossum and Drake 2009). The 3D point clouds were represented by an n × 3 NumPy array, describing the X, Y and Z co-ordinates of each point in real space. The point cloud was filtered to contain only points within a volume bounding the orchard area. A ground model was generated by binning the point cloud into an evenly spaced grid in the X-Y plane and finding the minimum altitude for each bin. This was smoothed by convolution with a Gaussian kernel. Points representing trees were isolated by removing all points above the ground model. The ground-excluded dataset was segmented into individual trees using the DBSCAN algorithm (Pedregosa et al. 2011). The D c equation used for the analysis is defined in Eq.
(2) as: where r is a positive real number and C(r) is the fraction of pairs of points which have a Euclidean separation smaller than r. For practical implementation of this, a distance matrix was calculated and subsequently the gradient of the least squares regression for r and C(r) where r min < r < 0.5(r max ), as discussed in Theiler (1990).

Statistical analyses
Python 3 (Rossum and Drake 2009) was used to complete all statistical analysis using the SciPy.stats package for ANOVA (testing for class differences) and Pearson's correlation (testing for correlation between continuous variables). The plots; Figs. 3, 4 were created using MatPlotLib (Hunter 2007) and Seaborn (Waskom et al. 2016) packages.

Results
During field data collection, it was noted that ten trees had died previously and had been removed from the site. Of the remaining 130 trees, two were classified as poor during the GR survey, and were excluded from analysis due to the low sample number in this condition category. There were significant inter-category (moderate, good and vigorous) differences in both the photogrammetric D f [F( 2125 ) = 953.202, p < 0.001] and the LiDAR D c measures of dimension [F( 2125 ) = 765.708, p < 0.001] at the median (Fig. 3). Measured photogrammetric D f values were in the range 1.38 > D f > 1.64 and LiDAR D c values in the range 1.02 > D c > 2.13. Both datasets were rescaled for comparison, where the minimum and maximum values were zero and one respectively. In both techniques, normalised dimension incrementally improved from the lowest values in the moderate category, to the highest values in the vigorous category (Fig. 3). A post hoc pairwise comparison between log r the three categories was made using Tukey's Range Test. For the photogrammetric dataset, the larger difference in dimensionality was between 'moderate' and 'good' (0.053, p < 0.05) compared to 'good' and 'vigorous' (0.025, p < 0.05). The same trend was observed in the LiDAR dataset where the difference in dimensionality was also larger between 'moderate' and 'good' (0.252, p < 0.05) than between 'good' and 'vigorous' (0.160, p < 0.05). Trees in the moderate category had significantly lower cumulative yields (median = 4.5 kg per tree) when compared to good (median = 7.2 kg), and vigorous (median = 6.6 kg) trees [F( 2125 ) = 198.319, p < 0.001] (Fig. 4a). There was a weak but statistically significant positive relationship between LiDAR D c value and cumulative yield (R 2 = 0.205, p = 0.020) (Fig. 4b). There was no relationship between photogrammetric D f and cumulative yield (R 2 = − 0.055, p = 0.541) (Fig. 4b).

Discussion
This study describes the novel measurement of orchard tree structure, using RS methods with fractal (D f ) and correlation (D c ) dimensionality analyses. This research showed that RS techniques can be used to objectively assign trees to standard management categories based on the assessment of tree structure. This is similar to the observations of Murray et al. (2018) who reported structural similarity probabilities that the quantified tree structure were represenatative of the observed referrence standard (i.e. trees in the best Fig. 3 Comparison of dimensionality values between two remote sensing methods, proximal photogrammetry and terrestrial laser scanning, of dwarf orchard trees in different condition categories. All boxes, and both investigative methods (photogrammetry, LiDAR) follow the trend of lower dimension scores for the moderate categorised trees, increasing at the median through the categories up to the vigorous category. There is some overlap in the interquartile range boxes (50% data distribution) for the good to vigorous categories, also across both investigative methods structural condition) to the order of; good (~ 99%), moderate (~ 89%), poor (~ 29%), dead (< 1%), using this field technique. Proximal photogrammetry and LiDAR were used to develop innovative quantification of structural variability in orchard trees. This is a novel approach for this specific application in an orchard environment. The findings of this study demonstrate that these RS methods can successfully discriminate between the structural nuances of heavily managed orchard trees. This approach provides opportunities for the development of additional methods and analyses that, it is believed, can lead to the replacement of subjectively driven, manual, management intervention techniques as are currently used.
Tree structure is cited as being a useful indicator for identifying the fruiting potential of orchard crops (Llorens et al. 2011;Gongal et al. 2016). However, despite the use of repeatable, high-resolution and high-accuracy RS data collection methods, this research shows that within young, healthy orchard trees, yield values are highly variable and difficult to predict simply from structural assessment. The D f method, previously shown to be an objective measure of structural condition in mature trees (Murray et al. 2018), was not able to satisfactorily predict cumulative yield. Similarly, due to the low degree of model fit and high variability of the yield at higher D c values, the D c analysis also does not offer a satisfactory proxy of yield (Fig. 4b). It is believed that the many combinations of scion to rootstock in the orchard, and the dwarf provenance of the trees used in the study affected this outcome. A further factor affecting the analysis is the low number of trees in condition categories other than 'good'.
The orchard is heavily managed and subject to repeated interventions, principally structural pruning, therefore natural structural variety will have been affected. However, an interesting observation shows a drop off in yield for vigorous trees (Fig. 4a). This response Fig. 4 a Cumulative yield (kg) of dwarf orchard trees, arranged by ground reference categories. The boxes show an increase in yield at the median from moderate to good categories, but a reduction in yield for the vigorous category when compared to the good category at the median. b Normalised dimensionality of photogrammetric and LiDAR assessment methods, modelled as a predictor of fruiting yield. The solid line indicates the regression line for the photogrammetric method, while the dashed line is the LiDAR method regression is indicative of the findings of Stephan et al. (2008) who compared two orchard pruning systems. Stephan et al. (2008) showed that vegetative (non-fruiting) shoots greatly increased as a response to central leader and branch pruning, causing proliferation of latent buds close to the pruning wound. The additional vegetative shoots increased light interception around the fruiting shoots and lowered crop yield overall. The increased structural complexity reduced the available energy resources for fruit production. Therefore, the vigorous trees endure periods of lower yields until an equilibrium is reached. It is suggested that LiDAR with D c analysis for quantifying pruning impacts in yield, particularly where repeat iterations of quantifying yield over several years is possible which was outwith the scope of this study, is an important area for further research and could lead to the adoption of RS based management in commercial orchards.

Conclusion
The RS methods demonstrated in this research are an improvement over traditional manual measurement of tree structural attributes, as high levels of analytical detail have been achieved on an individual tree basis. This research showed that D c analysis applied to terrestrial LiDAR data provided the more reliable approach for quantifying tree structure as a proxy for orchard crop yield when compared to the photogrammetric D f method. However, due to the complexities of investigating dwarf tree structures, combined with multiple combinations of scion and rootstock varieties, both assessment methods and dimensional analyses are imperfect predictors of yield for these orchard conditions. Based on the findings of this initial research, it is recommended that further investigation is required for the use of LiDAR data for quantifying tree structure as a proxy for crop yield, particularly focussing on the impacts of pruning on cumulative yield.