Extraction of Leaf Biophysical Attributes Based on a Computer Graphic-based Algorithm Using Terrestrial Laser Scanning Data

Leaf attribute estimation is crucial for understanding photosynthesis, respiration, transpiration, and carbon and nutrient cycling in vegetation and evaluating the biological parameters of plants or forests. Terrestrial laser scanning (TLS) has the capability to provide detailed characterisations of individual trees at both the branch and leaf scales and to extract accurate structural parameters of stems and crowns. In this paper, we developed a computer graphic-based 3D point cloud segmentation approach for accurately and efficiently detecting tree leaves and their morphological features (i.e., leaf area and leaf angle distributions (leaf azimuthal angle and leaf inclination angle)) from single leaves. To this end, we adopted a sphere neighbourhood model with an adaptive radius to extract the central area points of individual leaves with different morphological structures and complex spatial distributions; meanwhile, four auxiliary criteria were defined to ensure the accuracy of the extracted central area points of individual leaf surfaces. Then, the density-based spatial clustering of applications with noise (DBSCAN) algorithm was used to cluster the central area points of leaves and to obtain the centre point corresponding to each leaf surface. We also achieved segmentation of individual leaf blades using an advanced 3D watershed algorithm based on the extracted centre point of each leaf surface and two morphology-related parameters. Finally, the leaf attributes (leaf area and leaf angle distributions) were calculated and assessed by analysing the segmented single-leaf point cloud. To validate the final results, the actual leaf area, leaf inclination and azimuthal angle data of designated leaves on the experimental trees were manually measured during field activities. In addition, a sensitivity analysis investigated the effect of the parameters in our segmentation algorithm. The results demonstrated that the segmentation accuracy of Ehretia macrophylla (94.0%) was higher than that of crape myrtle (90.6%) and Fatsia japonica (88.8%). The segmentation accuracy of Fatsia japonica was the lowest of the three experimental trees. In addition, the single-leaf area estimation accuracy for Ehretia macrophylla (95.39%) was still the highest among the three experimental trees, and the single-leaf area estimation accuracy for crape myrtle (91.92%) was lower than that for Ehretia macrophylla (95.39%) and Fatsia japonica (92.48%). Third, the method proposed in this paper provided accurate leaf inclination and azimuthal angles for the three experimental trees (Ehretia macrophylla: leaf inclination angle: R2 = 0.908, RMSE = 6.806◦ and leaf azimuth angle: R2 = 0.981, RMSE = 7.680◦; crape myrtle: leaf inclination angle: R2 = 0.901, RMSE = 8.365◦ and leaf azimuth angle: R2 = 0.938, RMSE = 7.573◦; Fatsia japonica: leaf inclination angle: R2 = 0.849, RMSE = 6.158◦ and leaf azimuth angle: R2 = 0.947, RMSE = 3.946◦). The results indicate that the proposed method is effective and operational for providing accurate, detailed Remote Sens. 2019, 11, 15; doi:10.3390/rs11010015 www.mdpi.com/journal/remotesensing Remote Sens. 2019, 11, 15 2 of 23 information on single leaves and vegetation structure from scanned data. This capability facilitates improvements in applications such as the estimation of leaf area, leaf angle distribution and biomass.


Introduction
As the main organs of vegetation, leaves play a crucial role in vegetation growth [1].Leaf attributes are critical to describing the interactions between the land surface and the atmosphere, particularly in close relation to many biological and physical processes, such as photosynthesis, respiration, transpiration, and carbon and nutrient cycling [2].Therefore, the estimation of leaf structural and biophysical parameters is important for vegetation growth monitoring [3].The leaf attributes extracted in this paper include leaf area, and leaf angle distributions [4].Leaf angle distributions can be described by the leaf inclination angle and leaf azimuthal angle, which are defined as the angle between the leaf surface normal and the zenith [5] and the clockwise angle between the north direction and the projection of the principal axis of foliage on the horizontal plane [6], respectively.
The two major approaches for estimating leaf area involve the use of either direct or indirect methods.Methods have been devised to facilitate the direct measurement of leaf area and included image analysis [7] and the use of a conventional planimeter or an electronic leaf area metre [8,9].Direct methods for measuring the leaf area yield accurate results but are extremely laborious and time-consuming.Therefore, indirect methods have been developed and are typically preferred for this measurement.Indirect methods, by which the leaf area is inferred via measurements of other variables, such as gap probability [10] or light transmission through canopies, are efficient, non-destructive, and amendable to automation [11]; most indirect methods rely on optical instruments, such as the LAI-2000 Plant Canopy Analyser [12], Architecture of Canopies TRAC [13], Hemispherical Photography [14], and HemiView [15] instruments.Numerous attempts have been made to indirectly estimate LAI using aerial or satellite remote sensing imagery [16,17].Most of these attempts use approaches based on light reflectance from the canopy [18], and virtually all provide canopy-level estimates of LAI.These methods have had limited success, but estimates are often not precise enough to support site-specific forest management decisions.However, existing methods for estimating leaf angle distributions also involve direct or indirect approaches.A simple device consisting of a ruler, magnetic compass, and protractor [19] and a mechanical instrument consisting of high-precision potentiometers with three protruding arms [20] are commonly used to directly measure leaf inclination and azimuthal angles.Although direct methods can produce highly accurate results, they are extremely laborious and time-consuming [21].An indirect method to estimate leaf inclination and azimuthal angles based on digital photography was introduced [22] and has shown the potential to overcome the disadvantages of direct methods.However, this photographic method is difficult to automate, and substantial user interaction is required to identify suitable leaves upon which to measure leaf inclination angles [23]; in addition, this method must be carried out on towers, poles, ladders, unmanned aerial vehicles, and other conventional platforms [24].
Terrestrial laser scanning (TLS), which is used in measuring vegetation structure information, has advantages, such as a favourable directivity, high angular resolution, and strong anti-interference capability.In recent decades, laser scanning measurements techniques have enabled centimetre-level-accurate information to be acquired for individual trees.TLS is also a popular tool in forest ecology.The use of TLS has been intensively studied, e.g., for the estimation of leaf area index [25], crown structure [26], tree height [27] and diameter at breast height (DBH) [28], leaf area distributions and leaf angle distributions.This research has led to the development of accurate and effective leaf area extraction methods [29] and a leaf angle distribution estimation technique [30].The individual tree 3D modelling methods [31,32] from laser scanning data were also applied for the detection of various tree attributes.LiDAR-based leaf attribute retrieval has become an appealing concept due to its ability to capture the structural information of canopies as 3D point cloud data.The point cloud data acquired using this technique record geometrical information for each leaf, which has greatly improved the accuracy of leaf attribute estimation and has been used in numerous studies to retrieve leaf areas [33] and leaf angle distributions [34] from various vegetation types.A number of statistical models have been used to estimate leaf area, including the Leaf Area Constant Model [35], a statistical model based on linear measurements such as leaflet length and width [36], and allometric models that use non-destructive measurements of leaf length and/or width for estimating the leaf area [37].The regression-based method is used to estimate the leaf area of an individual tree on the basis of a regression model with LiDAR-derived tree dimensions [38].Computer graphics and different space partitioning methods are applied for point cloud processing.For example, the leaf area and leaf edge are accurately estimated by calculating the number of pixels of the RGB image [39], a voxel-based approach to retrieve the leaf area distribution for individual trees [40,41], and a voxel-based gap size algorithm to retrieve forest canopy clumping index [42].For leaf area deduction, the projected outer surface or projected tree row surface (PTRS) is linearly related to the leaf area [43].In addition, the non-destructive estimation of leaf area based on an artificial neural network approach has also been studied [44,45].Furthermore, a number of methods have been developed to estimate the leaf angle distributions.Multi-angular spectral data have been applied for identifying vegetation leaf angle distributions either from the structure sensitive index [46] or by model inversion [47].The information on leaf inclination angles from measurements of radiation penetration of the canopy can be retrieved at different view angles [48].A new method, the photographic method, was developed based on the analysis of levelled digital camera images of canopies consisting of flat leaves, which allows a rapid, non-contact and accurate estimation of leaf angle distributions [49].Leaf normals can be estimated using an algorithm that fits a plane to neighbouring LiDAR leaf intersection points [50].Furthermore, average leaf normals can be calculated by manually identifying leaves within LiDAR point cloud data and using planar fits [51].The method triangulates laser-leaf intersection points recorded by the LiDAR scan to calculate normal vectors [34].The non-photosynthetic and photosynthetic parts of the tree be separated to more accurately estimate leaf attributes and tree crown characteristics via approaches such as the geometric method using the 3D coordinates of each point [52] and the use of a series of features for each point [53].With these methods, progress has been made in the estimation of individual leaf areas and leaf angle distributions, but the description of tree crown details using the LiDAR technique still needs further improvement.
Although many studies have been performed using tree leaf attribute estimation based on TLS data, initial TLS data points are extremely numerous and unorganised [54].The current research faces the following problems upon obtaining the forestry parameters from scanning data.(1) The scanning data always experience a significant occlusion effect, which results in a large error in leaf attribute estimation.(2) Because the point cloud data contain a large number of noise points and because of the limitations of the performance of the instrument itself, it is complicated to extract parameters from tree scanning point cloud data.(3) For individual leaf-scale information, it is still necessary to design a computer graphic-based algorithm to accurately extract individual leaves from massive point clouds.
Considering leaves with different morphological structures and complex spatial distributions, we aimed to devise a computer graphic-based 3D point cloud segmentation approach to automatically detect single leaves by processing the whole point cloud of the tree crown to obtain individual leaf-scale information, including leaf area and leaf angle distributions.The specific objectives of this study are as follows.(1) Based on computer graphics methods, a sphere neighbourhood model with an adaptive radius and four auxiliary defined criteria was adopted, enabling the centre points of every individual leaf surface for three different tree species to be extracted efficiently.(2) Based on the extracted centre point of each leaf combined with a 3D watershed algorithm, the segmentation of each leaf was achieved.
(3) After the segmentation of each leaf, the area of each leaf was calculated by triangulation, and leaf inclination and azimuthal angles were calculated by computing the vector angle of each leaf surface.A qualitative comparison verification was carried out between the manual data and the data obtained using our methods for investigating the usability of leaf segmentation and verifying the effectiveness of the algorithm.

Study Site and Data Collection
The experimental trees were selected on the campus of Nanjing Forestry University (32 • 08'N, 118 • 81'E, WGS-84) and included many well-isolated individual Ehretia macrophylla (Ehretia macrophylla Wall), crape myrtle (Lagerstroemia indica L) and Fatsia japonica (Fatsia japonica (Thunb.)Decne.et Planch.)trees.All object trees with broad leaves and canopies with varying densities were scanned using a TLS (Leica C10) instrument, and one side-lateral scan was used to obtain the data.
The data were collected with a Leica C10 TLS system on May 2018.The Leica C10 instrument is a 532-nm phase-based scanner with a 360 • × 270 • upward field-of-view and laser rate of 50,000 points per second.The range measurement accuracy is ±1.5 mm at a distance of 3 m.The circular laser beam diameter and beam divergence at the scanner exit are 3 mm and 0.22 mrad, respectively, yielding a minimum distance between consecutive beams of approximately 0.4 mm at a distance of 3 m from the instrument.For the Ehretia macrophylla, crape myrtle and Fatsia japonica trees, the distances of the experimental trees and TLS were 5 m, 5 m and 3 m respectively, depending on tree height, in order to get the phenotypic characteristics of the whole tree in one scan within the viewing window.The height of the scanner was 1.57 m, 1.5 m and 0.8 m respectively.In this study, the normal scanning precision of the Leica C10 instrument was employed.
Three experimental trees, i.e., one Ehretia macrophylla, one crape myrtle and one Fatsia japonica tree, were selected for testing.Their structures represented by TLS point clouds are illustrated in Figure 1.The structural features of these trees, specifically tree height, crown base height and canopy projection area are listed in Figure 1.The LiDAR point density values were 16,946.3points (pts) • m −2 , 16,188.76 pts • m −2 and 51,483.96pts • m −2 for the Ehretia macrophylla, crape myrtle and Fatsia japonica trees, respectively; the heights of the trees were 3.5 m, 3.42 m and 1.6 m respectively; the crown base heights of the trees were 0.97 m, 1.26 m and 0.1 m, respectively; the projected areas of the crowns were 1.274 m 2 , 1.975 m 2 and 3.126 m 2 , respectively.The diameter at breast height (DBH) statistics for the trees are also listed in Figure 1, where DBH is defined as the sum of the diameter of each branch chest height.These comparisons indicate that the basic structural parameters of the sample trees were similar among all three tree species.Overall, there were no extreme cases of falling leaves and noise point increases due to the occurrence of strong winds for the selected trees.Thus, the selected experimental trees were appropriate for leaf attribute estimation.
To obtain the true leaf area and leaf angle distributions and demonstrate the validity of our method, the single-leaf area and leaf inclination and azimuthal angle of all target trees were measured using an LI-3000C portable area metre and angle measurement device.To validate the final results, we sampled the three experimental trees with different leaf numbers and manually measured the leaf inclination and azimuthal angles of all their leaves using an angle measurement device.For each experimental tree, the number of sampled leaves accounted for 40% of the total number of leaves in the crown, and the randomly sampled leaves were evenly distributed in the crown.Taking Ehretia macrophylla as an example, the crown was divided into upper, middle and lower portions.Each portion was divided into the eastern, western, northern and southern parts.Then, the tree crown was divided into 12 parts, with the sampled leaves distributed evenly in each part.The leaf angle distribution box plots (see Figure 2) show that the analytic dataset accounted for 40% of the total number of leaves in the tree crown, and the leaf inclination angle and leaf azimuth angle values tended to be relatively stable.To obtain the true leaf area and leaf angle distributions and demonstrate the validity of our method, the single-leaf area and leaf inclination and azimuthal angle of all target trees were measured using an LI-3000C portable area metre and angle measurement device.To validate the final results, we sampled the three experimental trees with different leaf numbers and manually measured the leaf inclination and azimuthal angles of all their leaves using an angle measurement device.For each experimental tree, the number of sampled leaves accounted for 40% of the total number of leaves in the crown, and the randomly sampled leaves were evenly distributed in the crown.Taking Ehretia macrophylla as an example, the crown was divided into upper, middle and lower portions.Each portion was divided into the eastern, western, northern and southern parts.Then, the tree crown was divided into 12 parts, with the sampled leaves distributed evenly in each part.The leaf angle distribution box plots (see Figure 2) show that the analytic dataset accounted for 40% of the total number of leaves in the tree crown, and the leaf inclination angle and leaf azimuth angle values tended to be relatively stable.The plots show that when the number of leaves analysed reaches 40% of the total number of leaves in the tree crown, the leaf inclination angle and leaf azimuth angle values tend to be relatively stable.Therefore, in the verification below, the analytic dataset accounted for 40% of the total number of leaves in the experimental tree.

Wood-leaf Separation
Wood-leaf separation, which aims to classify LiDAR points into wood and leaf components, is an essential prerequisite for achieving leaf separation and deriving individual tree leaf characteristics.An unknown degree of wood components likely causes errors in leaf area estimates [55].Considering the importance of wood-leaf separation results, in the current study, a series of features for each point [53] was used to separate the leaves of the experimental trees.These features  (d-f): the distributions of the true leaf azimuth angles of three experimental trees with proportion of sampled leaves number increasing.The plots show that when the number of leaves analysed reaches 40% of the total number of leaves in the tree crown, the leaf inclination angle and leaf azimuth angle values tend to be relatively stable.Therefore, in the verification below, the analytic dataset accounted for 40% of the total number of leaves in the experimental tree.

Wood-leaf Separation
Wood-leaf separation, which aims to classify LiDAR points into wood and leaf components, is an essential prerequisite for achieving leaf separation and deriving individual tree leaf characteristics.An unknown degree of wood components likely causes errors in leaf area estimates [55].Considering the importance of wood-leaf separation results, in the current study, a series of features for each point [53] was used to separate the leaves of the experimental trees.These features were the normal vector, the structure tensor and the distribution of the point normal vector.Figure 3 shows the wood-leaf separation results, including the wood points and the leaf points for each of the three experimental trees.

Individual Leaf Segmentation
The leaf blade is the main part of the tree crown and is mostly a green flat body with a large surface area, which is conducive to gas exchange and the absorption of light energy; this is of great relevance for studying leaf structure parameters.The sizes and shapes of leaves vary with different tree species.However, leaf morphology is relatively stable among leaves in the same plant and can be used as a basis for identifying plant species and individual leaf segmentation.The leaf blade width of the same plant has a relatively stable numerical range.The leaves of the three experimental trees studied in this paper have plane features; therefore, the corresponding scanned points of each leaf also have plane features.
Based on the above features, the proposed algorithm of individual leaf segmentation consists of three main stages, which are shown in Figure 4.The first stage is the extraction of the point clouds of the leaf central with a sphere neighbourhood model with an adaptive radius; four auxiliary criteria are defined to ensure the accuracy of the extracted central area points of each leaf surface.In the second stage, the density-based spatial clustering of applications with noise algorithm (DBSCAN) is used to cluster the central area points of the leaves, thereby obtaining the centre point corresponding to each leaf surface.In the third stage, individual leaf segmentation is realized by an advanced leaf blade segmentation algorithm combined with a 3D watershed algorithm.

Extracting the Central Area Points of Individual Leaves
For leaves of different morphological structures and complex spatial distributions, it is especially important to determine the spatial plane where each leaf blade is located.Through analysis, the leaves of the three experimental trees can be approximately treated as a surface, and each leaf has its own centre point, which is usually at the centre of the leaf.However, each leaf has a different zenith angle and azimuth angle, making it difficult to directly express the plane of each leaf.However, central points must be able to extend in a certain spatial direction across the entire leaf blade to form a plane, which can always be found in 3D space, that is, the plane where the point cloud of the leaf blade is located.In the sphere neighbourhood model of the adaptive radius proposed in this paper, when the leaf blade point cloud uniformly fills the great circle plane in the whole sphere model, the centre of the sphere is the central area point of the leaf.Randomly distributed scanned points lie on a plane in 3D space on each leaf surface.The ternary equation  (x, y, z) = 0 is applicable to this situation, and any point  ( ,  ,  ), ( = 1,2,3, … ,  ) on the leaf blade satisfies  ( ,  ,  ) ∈  , ( = 1,2,3, … ,  ). represents the total number of scanned points in each tree crown. represents the total number of leaves in each tree crown, and the value of  is unkown. represents each leaf surface in 3D space.However, each leaf surface is irregular and cannot be easily expressed by a mathematical equation.To present the surface equation for each leaf,  can be represented by the central area points of each leaf surface.The extraction of points in the centre area of each leaf is also the basis for the subsequent individual leaf segmentation.Therefore, this paper primarily uses the sphere neighbourhood model to extract the central area points of each leaf, as shown in Figure 5(a).
From the above analysis, each leaf blade is approximately in the plane distribution, and for the central area points of each leaf, the neighbourhood points must be uniformly distributed on a certain great circle of the sphere neighbourhood model.This paper defines four auxiliary criteria, and a

Extracting the Central Area Points of Individual Leaves
For leaves of different morphological structures and complex spatial distributions, it is especially important to determine the spatial plane where each leaf blade is located.Through analysis, the leaves of the three experimental trees can be approximately treated as a surface, and each leaf has its own centre point, which is usually at the centre of the leaf.However, each leaf has a different zenith angle and azimuth angle, making it difficult to directly express the plane of each leaf.However, central points must be able to extend in a certain spatial direction across the entire leaf blade to form a plane, which can always be found in 3D space, that is, the plane where the point cloud of the leaf blade is located.In the sphere neighbourhood model of the adaptive radius proposed in this paper, when the leaf blade point cloud uniformly fills the great circle plane in the whole sphere model, the centre of the sphere is the central area point of the leaf.Randomly distributed scanned points lie on a plane in 3D space on each leaf surface.The ternary equation S ζ (x, y, z) = 0 is applicable to this situation, and any point p i (x i , y i , z i ), (i = 1, 2, 3, . . ., n 1 ) on the leaf blade satisfies p i (x i , y i , z i ) ∈ S ζ , (ζ = 1, 2, 3, . . ., n 2 ).n 1 represents the total number of scanned points in each tree crown.n 2 represents the total number of leaves in each tree crown, and the value of n 2 is unkown.S ζ represents each leaf surface in 3D space.However, each leaf surface is irregular and cannot be easily expressed by a mathematical equation.To present the surface equation for each leaf, S ζ can be represented by the central area points of each leaf surface.The extraction of points in the centre area of each leaf is also the basis for the subsequent individual leaf segmentation.Therefore, this paper primarily uses the sphere neighbourhood model to extract the central area points of each leaf, as shown in Figure 5a.From the above analysis, each leaf blade is approximately in the plane distribution, and for the central area points of each leaf, the neighbourhood points must be uniformly distributed on a certain great circle of the sphere neighbourhood model.This paper defines four auxiliary criteria, and a point p i (x i , y i , z i ) that satisfies these auxiliary criteria is a leaf central area point.The point-to-point verification method is used to extract the centre area points of each leaf.
In this study, the extraction of leaf central area points by the point-to-point through the sphere neighbourhood model is the key step in individual leaf segmentation.First, the following variables must be defined: the radius of each sphere neighbourhood model is r 1 , and the centre of each sphere is p i , as shown in Figure 5(a4).The determination of the sphere neighbourhood model radius r 1 is a critical step in the pointwise classification because this variable can affect the classification accuracy.In this paper, the radius of the searching ball was set to r 1 = W tree /4 based on our sensitivity analysis to balance the classification accuracy and computational efficiency.W tree represents the average leaf width in the whole tree crown of the current process tree.For a point p i (x i , y i , z i ) in the point cloud P, P ⊂ R 3 , the neighbourhood points of p i within the radius r 1 were defined as p i,j x i,j , y i,j , z i,j (j = 1, 2, 3, . . ., n 3 ) and satisfied the condition p i,j − p i ≤ r 1 .n 3 is the number of neighbourhood points p i,j for the point p i (i.e., the number of points in the sphere neighbourhood model).The value of n 3 is should greater than thresholds0.The thresholds0 for the three experimental trees (Ehretia macrophylla, crape myrtle and Fatsia japonica) is set at 15, 10 and 30, respectively, according to the apparent characteristics of plants, leaf size and leaf area density.
where p i denotes the mean of the neighbourhood points p i,j for the point p i .Criterion 1: Find the point p i that satisfies Equation (2).Equation ( 2) represents the distance between the current point p i and point p i , which is less than the thresholds1.Equation ( 2) ensures that the current point p i is close to its centre of the neighbourhood points.
Criterion 2: Find the point p i that satisfies Equation (3).Equation ( 3) expresses that the current point p i must be close to the fitting plane S i .S i is generated using the least squares method from the current point p i and its neighbourhood points p i,j .The fitting plane is defined as S i : Ax i + By i + Cz i + D = 0.The distance between the current point p i and the fitting plane S i is less than thresholds2.
where the value of thresholds2 is half of the thresholds1 value.Criterion 3: Ensure that the distribution composed of the current point p i and its neighbourhood points has spatial planar features.Equation (4) guarantees that the distance between all the neighbourhood points p i,j of the current point p i and the fitting plane S i are less than the thresholds2, which ensures that all the neighbourhood points p i,j are close to the fitting plane S i , i.e., the mean value of the distance from each neighbourhood point p i,j of the current point p i to the fitting S i is less than the thresholds2.
where n 3 is the number of neighbourhood points p i,j for point p i .Criterion 4: S c i represents the great circle obtained by the intersection of the fitting plane S i and the sphere neighbourhood model of the current point p i , as shown in Figure 5(b1,b2).The great circle plane S c i is divided into the t block area as shown in Figure 5(b1,b2).The projection points of p i and p i,j on the great circle plane S c i are p i and p i,j , respectively.It is reasonable that the points p i,j should be uniformly distributed around point p i on the great circle plane S c i .The variation ratio is used to measure the degree by which the point p i and its neighbourhood points p i,j are distributed on the great circle plane S c i .(Num p , i,j − Num v )/Num p , i,j denotes the variation ratio of discrete points in each v(v = 1, 2, . . ., t) block area.Num p , i,j denotes the total number of projection points on the great circle plane, and Num v is the number of the projection points in the vth block area Equation ( 5) guarantees that the point p i and its neighbourhood points p i,j are evenly distributed in each block area.This criterion eliminates candidate points at the leaf edge that overlap with other leaves, which can be mistakenly classified as The central area points p carea ξ (ξ = 1, 2, 3, . . ., n 4 ) satisfying Equations ( 1)-( 5) are extracted as black dots shown in Figure 5(a1).n 4 is the total number of central area points of each tree.

Clustering of Leaf Central Area Points
In real situations, the central area points of the leaves in each tree crown are divided into a large number of clusters.Therefore, we adopted the cluster algorithm proposed by Ester et al. [56], DBSCAN, to determine the central points of leaf segmentation.DBSCAN requires two input parameters containing the minimum number of points (MinPts) needed to form a cluster and the maximum radius of the neighbourhood from the core point (MaxR) (the maximum distance between clusters, here).The determination of the two input parameters is a critical step in the leaf central area points clustering process since these parameters can affect the accuracy of the clustering result.In this paper, the parameter MinPts for the three experimental trees (Ehretia macrophylla, crape myrtle and Fatsia japonica) was set at 15, 10 and 30, respectively, according to the scanned point density and area of each leaf surface, which can be adopted to obtain the appropriate parameter MaxR and achieve the classification result.The extracted data of the central area points of each leaf have a uniform distribution; then, MaxR can be calculated from MinPts and the size of the point cloud using the following equations [57]: , respectively.MaxR was obtained by Equation ( 6), calculated as half of the average leaf width for each tree crown.Then, the central area points were clustered based on the two input parameters MaxR and MinPts.After the clustering analysis using the DBSCAN algorithm, the central area points of each leaf were segmented into different clusters, as shown by the coloured dots in Figure 4d.After the central area points p carea ξ were divided into n 2 clusters, p carea ξ,k (ξ = 1, 2, 3, . . ., n 4 , k = 1, 2, 3, . . ., n 2 ) was obtained to represent the central area points of each leaf blade surface, where k indicates that the ξth centre area point belongs to the kth cluster.The centre point for central area points p carea ξ,k of each leaf blade is denoted as c k (k = 1, 2, 3 . . .n 2 ) and is labelled by red stars in Figure 4d.c k represents the centre point of each leaf in the tree crown.These points serve as the seed point for the 3D watershed algorithm (Section 2.3.3).Each cluster of central area points p carea ξ,k corresponds to one centre point of each leaf c k .

Individual Leaf Segmentation
The 3D watershed algorithm [58,59] is used to segment the remaining leaf point cloud of the tree crown, that is, the non-central area points p l (l = 1, 2, 3 . . .n 5 ), where p l ∪ p carea ξ = p i .n 5 is the total number of non-central area points of each tree.In our method, the watershed algorithm starting points are the centre points of each leaf c k , which were determined using our methods.The two morphology-related parameters, namely, the Euclidean distance between the centre point and the non-central area point and the cosine of the angle between the vector composed of the centre point and the non-central area point and the normal vector of the leaf blade, were used to segment individual leaves, which reduced the leaf oversegmentation problem and improved the precision of segmentation in the leaf edge detection.The ultimate goal of Equation ( 8) is to find the minimal value for each non-central area point p l , thereby completing the classification of each p l value with the corresponding centre point c k .To balance the two terms, the parameters a 1 and a 2 should be adjusted.The equivalent relationship of the parameters a 1 and a 2 is a 1 + a 2 = 1.In the experiment, the values of a 1 and a 2 were set to 0.7 and 0.3, respectively.Figure 6 illustrates Equation ( 8).In the experiment, the values of  and  were set to 0.7 and 0.3, respectively.Figure 6 illustrates Equation 8.

Leaf Attributes Calculation Based on Classified Leaf Points
After single-leaf segmentation via our method, Delaunay triangulation was adopted to deduce the area of each leaf.For segmented single-leaf scanning data, the normal vector of each leaf surface was obtained by fitting the leaf blade plane with the least squares method.For the known single-leaf point cloud data, the leaf inclination angle was used to calculate the included angle between the vector of the leaf surface normal and the zenith angle.According to the definition of the leaf azimuth angle, it is important to extract the principal axis of foliage.In this paper, the principal axis of each leaf blade is defined as the line between the most distant two points in a single-leaf point cloud.Thus, the leaf azimuthal angle can be obtained by calculating the clockwise angle between the north direction and the projection of the principal axis of foliage on the horizontal plane.

Leaf Attributes Calculation Based on Classified Leaf Points
After single-leaf segmentation via our method, Delaunay triangulation was adopted to deduce the area of each leaf.For segmented single-leaf scanning data, the normal vector of each leaf surface was obtained by fitting the leaf blade plane with the least squares method.For the known single-leaf point cloud data, the leaf inclination angle was used to calculate the included angle between the vector of the leaf surface normal and the zenith angle.According to the definition of the leaf azimuth angle, it is important to extract the principal axis of foliage.In this paper, the principal axis of each leaf blade is defined as the line between the most distant two points in a single-leaf point cloud.Thus, the leaf azimuthal angle can be obtained by calculating the clockwise angle between the north direction and the projection of the principal axis of foliage on the horizontal plane.

Plant Leaf Classification
After applying our classification method to the preliminary scanned datasets of the experimental trees, we acquired an excellent consistency between the point clouds and real plant dimensions.These good results can be largely attributed to the proper selection of features, including the normal vector, the structure tensor and the distribution of the point normal vector.The quantitatively assessed classification accuracies of each target tree class are listed in Table 1.The tree height is a field measurement.The crown projection area is calculated by convex hull algorithms [60] of the single trees based on the scanned points.One result for leaf point number and the branch point number is obtained by the wood-leaf separation method in this study and the other using the manual measurement results.The overall accuracy is the ratio of the results produced by the classification method in this paper to the benchmark results derived using manual measurements.

Leaf Segmentation
Three experimental trees were segmented using our method.The detailed results are shown in Figure 7 and demonstrate the classification results for Ehretia macrophylla (Figure 7a), crape myrtle (Figure 7b) and Fatsia japonica (Figure 7c).The accuracy of leaf segmentation is 94.0%, 90.6% and 88.8% for the Ehretia macrophylla, crape myrtle and Fatsia japonica trees, respectively.The accuracy of leaf segmentation in the upper tree crown (97.7%) was found for the Ehretia macrophylla trees than for the middle tree crown (90.6%) and lower tree crown (93.1%).The accuracy of leaf segmentation in the upper tree crown (92.2%) than for the lower tree crown (88.8%) was found for the crape myrtle trees.The central area points of each leaf for the experimental trees were accurately extracted by our method, and the extraction results are shown as black dots in Figure 7(a1,b1,c1).Then, the extracted central area points of the experimental trees were clustered using the DBSCAN algorithm, and the central area points of individual leaves were segmented into different clusters, as shown by the coloured dots in Figure 7(a2,b2,c2).The best segmentation effect was obtained for Ehretia macrophylla, which has relatively flat leaves in the overall canopy and an occlusion effect that is not obvious.The second best results were obtained for crape myrtle, which has a small leaf size and high leaf density; this tree type also achieved a good segmentation effect.Finally, Fatsia japonica, with the most complicated blade shape structure, also had a good segmentation effect.

Leaf Attribute Calculation
After segmentation via our method, Ehretia macrophylla was divided into upper, middle and lower sections with three different horizontal layers, crape myrtle was divided into upper and lower sections with 2 different horizontal layers, and Fatsia japonica was not stratified.The segmented leaf points in the crowns were used as preliminary data and were randomly chosen as the analytic dataset.For each target tree, the leaves in the analytic dataset were evenly distributed in four directions, east, west, north and south, and the number of leaves in the analytic dataset accounted for 40% of the total number of leaves in the crown of the target tree.In this study, Delaunay triangulation was adopted to estimate each designated leaf area, and vector angles were computed to obtain leaf azimuthal angle and leaf inclination angle values.The specific leaf parameters obtained (e.g., leaf length, leaf width, leaf area, leaf azimuthal angle and leaf inclination angle) using our methods for different tree species are listed in Table 2.A comparison of the number of leaves in each crown obtained via manual measurements and using our method is presented in Table 2, which indicates that our method obtained satisfactory estimates of the number of individual tree leaves.Table 2 shows that the leaf areas determined by the Delaunay triangulation algorithm are very similar to the actual leaf areas obtained by manual measurements.In addition, the values of the angle distribution calculated directly on the point clouds were effectively highly correlated with the measurements of the actual angles taken in the field.Figure 8 depicts the measured vs. estimated individual leaf attribute values, including leaf area, leaf azimuthal angle, and leaf inclination angle, obtained using our method and manual measurements.The leaf attribute estimations of the correctly detected trees using our method fall close to the 1:1 line on, the scatterplots of field-measured and method-estimated leaf attributes for cross-validation (Figure 8).The widths of the confidence intervals for Ehretia macrophylla are relatively narrow, indicating the high quality of Ehretia macrophylla's leaf segmentation and attribute estimation result.In general, there were no significant mean differences between the field-measured and method-estimated values for leaf area, inclination and azimuthal angles (see Table 3).This result indicated that the individual leaf area estimation was in good agreement with values obtained from actual measurements and that the individual leaf area in the upper tree crown was larger than that in the middle and lower tree crown (Figure 8(a1,b1)).The leaf inclination angle estimations were in good agreement with values obtained from actual measurements for the three experimental trees (R 2 = 0.908, RMSE = 6.806 • , R 2 = 0.901, RMSE = 8.365 • and R 2 = 0.849, RMSE = 6.158 • , respectively; Figure 8(a2,b2,c2)).A higher probability of large inclination angles in the upper crown was found for the Ehretia macrophylla and the crape myrtle trees than for the middle tree crown and lower tree crown (Figure 8(a2,b2)).One explanation for this result could be that upper crown leaves grow more vertically than lower crown leaves to increase the amount of solar radiation penetrating through the plant canopy, ultimately increasing the possibility of light penetration and benefiting the growth of lower foliage elements [23].Similarly, lower leaf elements tend to grow more horizontally than upper leaf elements to optimize the absorbed solar radiation.The leaf azimuthal angles of the estimations were in good agreement with values obtained from actual measurements of the three experimental trees for Ehretia macrophylla, crape myrtle and Fatsia japonica.(R 2 = 0.981, RMSE = 7.680 • , R 2 = 0.938, RMSE = 7.573 • and R 2 = 0.947, RMSE = 3.946 • , respectively; Figure 8(c1-c3)).The variations in the leaf azimuthal angle distributions with tree height were very small.In the Ehretia macrophylla tree, leaf azimuthal angles were randomly distributed from 3.3 • to 353.6 • across the whole crown (Figure 8(a3)).The azimuthal angle distribution of the crape myrtle tree from the upper to lower layers was very similar (Figure 8(b3)).In addition, the Fatsia japonica azimuthal angle probability ranged from 50.9 • to 266.3 • , which might be attributed to the light source orientation of vegetation leaves or the incomplete data collection of the one side-lateral TLS (Figure 8(c3)).These findings demonstrate that the segmentation method allows the robust determination of the distribution of orientation values.Comparisons of the estimated leaf attributes between our algorithm and the manual method for the three target trees: (a1,b1,c1): scatter plots of the reference values obtained using LI-3000C versus the LiDAR-estimated individual tree leaf areas obtained using our method; (a2,b2,c2): scatter plots of the reference values obtained by manual measurements and individual tree leaf inclination angles obtained using our method; (a3,b3,c3): scatter plots of the reference values obtained by manual measurements and the individual tree leaf azimuth angles obtained using our method.

Discussion
For the optimization of plant cultivation, the monitoring of crop growth and productivity, successive leaf attribute measurements are necessary [61,62].It is well established that leaf area is pivotal for light interception, photosynthesis and transpiration and leaf angle distribution is an important canopy structure parameter that has a great influence on the transmission of radiation within vegetation canopies and the distribution of incident photosynthetically active radiation [63,64].Thus these features are key parameters in several agronomic and physiological studies [65][66][67].TLS has emerged as among the most promising remote sensing technologies, providing detailed, spatially explicit, 3D information on plant structure, for operational applications in a wide range of disciplines related to the extraction of plant morphological features.These applications include, as reported in previous literature, the estimation of forest biophysical parameters such as gap fraction [68] and leaf area index [25,69].These methods mostly aim to extract the vertical and horizontal distributions of biophysical variables or to determine an empirical correlation formula between laser scanning and in situ data sets.The literature [70,71] has reported that the non-photosynthetic and photosynthetic parts of trees should be separated to achieve more accurate biomass estimations.The specific analysis of the photosynthetic part of trees by TLS improves leaf area index estimations.Gaussian mixture models [72] are successfully used for separating leaf and non-leaf classes of individual trees, and a support vector machine classification technique [53] has been presented to extract the leaves of trees.On this basis, the leaf separation of the canopy was successfully achieved in the present study, and parameter analysis at the single-leaf scale was completed.The following points must be considered to apply this method to leaf segmentation and leaf attribute measurements.

Tree Species Effects
Tree species showed a negligible effect on the final classification accuracy using the proposed method.For the experimental trees, the final accuracy for Ehretia macrophylla was 94.0% (see Table 2), the final accuracy for crape myrtle was 90.6% (see Table 2) and the final accuracy for Fatsia japonica was 88.8% (see Table 2).The Fatsia japonica tree, with the most complicated blade shape structure, also had good segmentation results.The study results show that our method is not suitable for coniferous trees because it is difficult to determine the centre points of coniferous leaves by using sphere neighbourhood model, but the phenotypic traits of broad leaves with wide surfaces are suitable for our algorithms.Although different species have varying leaf inclination angles and leaf density distributions, resulting in different degrees of occlusion and the need for different search radii and thresholds to be set, our method can still provide a valuable and reasonable approach for assessing the real leaf area.However, for TLS, there are more occlusions and a greater lack of scanned data in the upper part of the tree crown in tall trees, which leads to low accuracy in the segmentation of individual leaves by the algorithm.We can consider the combination of airborne laser scanning (ALS) [73] and TLS scan modes to obtain complete scanned data and use the method proposed in this study to calculate each leaf attribute in the whole tree crown.Thus, our proposed method can accurately and efficiently detect individual tree leaves and their morphological attributes from scanned point cloud data.

Sensitivity Analysis of Search Radius (r 1 )
The search radius r 1 affected both the preliminary and final individual leaf segmentation results.To better understand the effect of the search radius r 1 on the central area point extraction accuracy, a sensitivity analysis was conducted using the data from the three experimental trees.The sensitivity analysis showed that the correctly segmented leaf number (CSLN) curves for all experimental tree crowns changed as the search radius increased.As shown in Figure 9, the CSLN increased as the search radius r 1 increased for all experimental trees.The highest individual leaf segmentation accuracy was obtained when the search radius r 1 reached 3.6 cm, 0.9 cm and 4.2 cm for Ehretia macrophylla, crape myrtle and Fatsia japonica, respectively.In this case, the values of search radius r 1 (3.6 cm, 0.9 cm and 4.2 cm) were approximately equal to one quarter of the average leaf width in each tree crown.The average widths of the leaves in the tree crown for Ehretia macrophylla, crape myrtle and Fatsia japonica were 12.8 cm, 3.6 cm and 17.2 cm, respectively; in addition, the CSLN for the experimental tree crowns provided increasingly lower segmentation accuracy when the search radius r 1 was greater than 3.6 cm, 0.9 cm and 4.2 cm, respectively.In the case in which the scanning precision remain unchanged, the shorter the search radius, the fewer the points in the sphere neighbourhood model.When the number of points in the sphere neighbourhood model was greater than the thresholds0 and the four auxiliary criteria (discussed in Section 2.3.1)defined were satisfied, the point p i satisfied the conditions to be considered the central area point of the leaf.Therefore, when the radius of the sphere neighbourhood model decreases, the number of points in the great circle plane will be decreased and the number of calculated central area points will be gradually reduced, resulting in a downtrend in the recognition of leaves using our program.The value of the input parameter MaxR also affected the final individual leaf segmentation when leaf central area points were clustered using the DBSCAN algorithm.As shown in Figure 9, when the value of MaxR was larger, the number of leaf blades segmented was smaller; when the value of MaxR was smaller, the number of leaf blades segmented was larger.Overall, when the value of MaxR was approximately half of the average leaf width in the tree crown, the results showed the highest segmentation accuracy (see Figure 9).The intersection of two lines for each figure in Figure 9 represents the optimal detected leaf number result, which is closest to the true value.

Recommendations
The occlusion of forest trees is determined by various factors such as leaf area index, the location of the scanner, scanning parameter setting [74], brand of the scanner [75] and tree Our results suggest that for the three experimental trees, the CSLN changed as the search radius r 1 increased, and the highest segmentation accuracy was obtained when r 1 was one-quarter of the average width of the leaves of each tree crown.In addition, the best selected value of the input parameter MaxR was approximately half of the average leaf width for the crowns of the different tree species.

Recommendations
The occlusion of forest trees is determined by various factors such as leaf area index, the location of the scanner, scanning parameter setting [74], brand of the scanner [75] and tree topological structure.These factors deserve further exploration.The TLS setup could affect the classification accuracy, and the TLS-measured data were scanned from one perspective which might cause the representation accuracy of the leaves to be relatively low and not homogeneous, with the shadow effect on the other side of the tree.Considering that the occlusion inside the tree crown of the three experimental trees in this study is not serious, only one-sided scanning was sufficient to obtain nearly overall information of the tree crown and acquire the corresponding high-density point cloud.Furthermore, the spatial resolution of the point clouds decreased as the acquisition distance between the TLS instrument and the experimental tree increased, causing leaves or branches to be misclassified and tree stems or branches surrounded by shrubs or grasses to often be misclassified as random leaf points.As a result, the wood-leaf separation and leaf point occlusion were greater, which reduced the leaf segmentation accuracy.Finally, the noise points caused by wind had a great influence on the leaf segmentation and the calculation of the area of individual leaves.Thus, it was very important to remove noise points for the direct estimation of individual leaf areas.To increase the accuracy of the segmentation of individual leaves, we recommend using multiple scan locations to better capture the 3D structural characteristics of experiment trees.In addition to the above factors, the quality of the point cloud data (including accuracy and completeness) was another important factor that increased the algorithm accuracy of single leaf blade segmentation.The low segmentation precision rate for the leaf blade in the TLS-based data was primarily caused by the incompleteness of the wood-leaf separation.The quality of the TLS-based data was affected by three factors: outliers, woody phytoelement points, and occlusion.Although we removed outliers and woody phytoelement points, the influence of these three factors could not be completely eliminated.The quality of the TLS-based data still has room to improve.Therefore, the need to achieve better denoising and the separation of woody points from foliage points represents a great challenge to accurate leaf segmentation in future studies.
Our method's final individual leaf segmentation accuracy increased as the occlusion of the tree crown decreased.The results showed that the leaf recognition ratio decreased for TLS-based data as the leaf numbers and point cloud density increased (see Table 2).This phenomenon could be explained by the fact that the degree of overlapping in the leaf distribution increased as the leaf number increased, and a high degree of leaf overlapping increased the difficulty of leaf segmentation and decreased the leaf recognition ratio.The individual leaf segmentation accuracy was as low as 88.80% for Fatsia japonica, which was characterized by a high point cloud density (51,483.96pts • m −2 ).Among the three experimental trees, Fatsia japonica showed the lowest final leaf segmentation accuracy.This result was likely due to the more severe occlusion within the tree crown, leading to the incomplete wood-leaf separation, and the complicated leaf blade shape structure inherent in Fatsia japonica.In the case of the crape myrtle, which was characterized by a medium point cloud density (16,188.76pts • m −2 ), the individual leaf segmentation was as high as 93.62%.The leaf segmentation accuracy of the crape myrtle in two different horizontal layers was also negatively correlated with point density (see Table 2).Finally, for Ehretia macrophylla, which was characterized by a low point cloud density (16,946.3pts • m −2 ), the final individual leaf segmentation reached up to 94.0%.For the upper, middle and lower layers of the Ehretia macrophylla tree, which represented three different horizontal layers, the internal occlusion and point density of the upper tree crown were lower than those of the other layers, so the individual leaf segmentation accuracy was the highest; the middle tree crown had the highest number of leaves, and the individual leaf segmentation accuracy of the middle layer was slightly lower than the segmentation accuracy of the other layers of the tree crown (see Table 2).

Conclusions
This study developed a computer graphic-based 3D point cloud segmentation approach for accurately and efficiently detecting tree leaves and their morphological features, including leaf area, leaf azimuthal angle and leaf inclination angle, from individual leaves.First, based on the proposed computer graphics methods, the sphere neighbourhood model with an adaptive radius was used to extract the central area points of each leaf blade.In addition, four auxiliary criteria were defined to ensure the accuracy of the extracted central area points of each leaf surface.Second, DBSCAN was used to cluster the central area points of individual leaves.Third, based on the extracted centre point of each leaf combined with the advanced 3D watershed algorithm, individual leaf segmentation was achieved.Finally, leaf morphological features (leaf area and leaf angle distributions) were calculated from the single-leaf point cloud data of the segmented leaves to derive single-leaf-scale information and investigate the usability of leaf segmentation.

Figure 1 .
Figure 1.Experimental tree selection and scanning data collection.Terrestrial laser scanning (TLS) was used to scan the different tree species, namely (a) Ehretia macrophylla, (b) crape myrtle and (c) Fatsia japonica; (d), (e) and (f) show the scanned point clouds and structural features of the experimental trees.

Figure 1 .
Figure 1.Experimental tree selection and scanning data collection.Terrestrial laser scanning (TLS) was used to scan the different tree species, namely (a) Ehretia macrophylla, (b) crape myrtle and (c) Fatsia japonica; (d-f) show the scanned point clouds and structural features of the experimental trees.Remote Sens. 2018, 10, x FOR PEER REVIEW 6 of 25

Figure 2 .
Figure 2. Box plot showing the distributions of the true leaf inclination angles and leaf azimuth angles of the three experimental trees.(a), (b), (c): the distributions of the true leaf inclination angles of three experimental trees with proportion of sampled leaves number increasing.(d), (e), (f): the distributions of the true leaf azimuth angles of three experimental trees with proportion of sampled leaves number increasing.The plots show that when the number of leaves analysed reaches 40% of the total number of leaves in the tree crown, the leaf inclination angle and leaf azimuth angle values tend to be relatively stable.Therefore, in the verification below, the analytic dataset accounted for 40% of the total number of leaves in the experimental tree.

Figure 2 .
Figure 2. Box plot showing the distributions of the true leaf inclination angles and leaf azimuth angles of the three experimental trees.(a-c): the distributions of the true leaf inclination angles of three experimental trees with proportion of sampled leaves number increasing.(d-f): the distributions of the true leaf azimuth angles of three experimental trees with proportion of sampled leaves number increasing.The plots show that when the number of leaves analysed reaches 40% of the total number of leaves in the tree crown, the leaf inclination angle and leaf azimuth angle values tend to be relatively stable.Therefore, in the verification below, the analytic dataset accounted for 40% of the total number of leaves in the experimental tree.

18 Figure 3 .
Figure 3. Results of wood-leaf separation with a series of features for each point, i.e., the normal vector, the structure tensor and the distribution of the point normal vector, for Ehretia macrophylla (first row), crape myrtle (second row) and Fatsia japonica (third row).(a) Initial point clouds, (b) wood points and (c) leaf points.

Figure 4 .
Figure 4. Flowchart of the proposed algorithm of individual leaf segmentation.The entire experimental method consists of three stages.(a) The first stage is the extraction of the central area points of each leaf using a sphere neighbourhood model with an adaptive radius and four auxiliary criteria.The black dots extracted from the first stage comprise the central area points of each leaf.(b-d) In the second stage, the central area points of each leaf are clustered using the density-based spatial clustering of applications with noise algorithm (DBSCAN) algorithm, and the centre point of each leaf is obtained, as labelled by the red star.(d-e) In the third stage, individual leaf separation is achieved using a spatial watershed algorithm based on the centre point of each extracted leaf surface and the two morphology-related parameters.

Figure 4 .
Figure 4. Flowchart of the proposed algorithm of individual leaf segmentation.The entire experimental method consists of three stages.(a) The first stage is the extraction of the central area points of each leaf using a sphere neighbourhood model with an adaptive radius and four auxiliary criteria.The black dots extracted from the first stage comprise the central area points of each leaf.(b-d) In the second stage, the central area points of each leaf are clustered using the density-based spatial clustering of applications with noise algorithm (DBSCAN) algorithm, and the centre point of each leaf is obtained, as labelled by the red star.(d-e) In the third stage, individual leaf separation is achieved using a spatial watershed algorithm based on the centre point of each extracted leaf surface and the two morphology-related parameters.

18 Figure 5 .
Figure 5. Schematic illustrating the extraction of the centre point of each leaf.(a) The radius of the sphere is r width/4.Black dots indicate the central area points extracted from each leaf.(a2), (a3)and (a4) show that criteria 2 and 3 are satisfied, i.e., the distance of point  and its neighbourhood points  , to the fitting plane  is less than the thresholds.Criterion 4 guarantees that the neighbourhood points  , , are uniformly distributed around point  , on the great circle plane  .(b1)shows a case that satisfies criterion 4; (b2) shows a case that does not satisfy criterion 4.

Figure 5 .
Figure 5. Schematic illustrating the extraction of the centre point of each leaf.(a) The radius of the sphere is r 1 = width/4.Black dots indicate the central area points extracted from each leaf.(a2-a4) show that criteria 2 and 3 are satisfied, i.e., the distance of point p i and its neighbourhood points p i,j to the fitting plane S i is less than the thresholds.Criterion 4 guarantees that the neighbourhood points p i,j are uniformly distributed around point p i on the great circle plane S c i .(b1) shows a case that satisfies criterion 4; (b2) shows a case that does not satisfy criterion 4.

5 )Figure 5 (
Figure 5(b1) shows an example of satisfied criterion 4, and Figure 5(b2) shows an example of unsatisfied criterion 4.The central area points p carea

T = Max x p carea ξ −
Min n 4 denotes the number of central area points of each tree crown.The value of n is 3 in this study, representing the dimensionality of the points.Γ is the gamma function, and T is the volume of the experimental space formed by m points.p carea ξ denotes the data set of the central area points, and x p carea ξ , y p carea ξ , z p carea ξ are the x, y and z values of p carea ξ (l = 1, 2, 3 . . .n 5 , k = 1, 2, 3 . . .n 2 ) represents the point after classifying, as shown by coloured dots shown in Figure 4f.The first term on the right-hand side of equation p l − c k 2 2 represents the Euclidean distance between the centre point c k and the non-centre point p l .The second term on the right-hand side of the equation cos (∠( → (c k , p l ), → c k )) calculates the cosine of the included angle between the vector composed of the centre point c k and the non-central area point p l and the normal vector of the leaf blade surface → c k .
Remote Sens. 2018, 10, x FOR PEER REVIEW 12 of 25 segmentation in the leaf edge detection.The ultimate goal of Equation 8 is to find the minimal value for each non-central area point  , thereby completing the classification of each  value with the corresponding centre point  . , = argmin  ‖ −  ‖ +  cos ∠(( ,  ) ⃗ ,  ⃗) (8) where  , ( = 1,2,3 …  ,  = 1,2,3 …  ) represents the point after classifying, as shown by coloured dots shown in Figure 4f.The first term on the right-hand side of equation ‖ −  ‖ represents the Euclidean distance between the centre point  and the non-centre point  .The second term on the right-hand side of the equation cos ∠(( ,  ) ⃗ ,  ⃗) calculates the cosine of the included angle between the vector composed of the centre point  and the non-central area point  and the normal vector of the leaf blade surface  ⃗.To balance the two terms, the parameters  and  should be adjusted.The equivalent relationship of the parameters  and  is  +  = 1.

Figure 6 .
Figure 6.Figure 6 depicts the role of Equation 8 in the segmentation of an individual leaf, where the distance (i.e., the Euclidean distance between the centre point and the non-central area point) and the included angle (i.e., the included angle between the vector composed of the centre point and the non-central area point and the normal vector of the leaf blade) are used for leaf segmentation.(a) The distance is the main contributing factor to completing individual leaf segmentation, which means that the non-central area point  is closer to the real leaf centre point  .(b) The included angle between the vector (  ,  ) ⃗ and the normal vector  ⃗ of each leaf is the main contributing factor to completing individual leaf segmentation.The included angle of the two vectors noted above composed of points belonging to the same leaf surface is nearly 90 degrees, which decreases the value of the second term on the right-hand side of Equation 8.

Figure 6 .
Figure 6. Figure 6 depicts the role of Equation (8) in the segmentation of an individual leaf, where the distance (i.e., the Euclidean distance between the centre point and the non-central area point) and the included angle (i.e., the included angle between the vector composed of the centre point and the non-central area point and the normal vector of the leaf blade) are used for leaf segmentation.(a) The distance is the main contributing factor to completing individual leaf segmentation, which means that the non-central area point p l is closer to the real leaf centre point c k .(b) The included angle between the vector

Figure 7 .
Figure 7. Scatter plots of point cloud canopy data obtained using our proposed method to automatically segment and randomly colour leaves.The obtained individual foliage points were prepared for leaf parameter estimates.(a) Ehretia macrophylla; (b) crape myrtle; and (c) Fatsia japonica.(a1,b1,c1): The black dots indicate the central area points of each leaf extracted by the sphere neighbourhood model, and the green dots represent the entire leaf points of each leaf.(a2,b2,c2): The central area points were clustered using density-based spatial clustering with the noise algorithm and visualized with random colours.(a3,b3,c3): The leaf segmentation results, with partial close-up images presented in (a4,b4,c4), respectively.

Figure 8 .
Figure 8. Comparisons of the estimated leaf attributes between our algorithm and the manual method for the three target trees: (a1,b1,c1): scatter plots of the reference values obtained using LI-3000C versus the LiDAR-estimated individual tree leaf areas obtained using our method; (a2,b2,c2): scatter plots of the reference values obtained by manual measurements and individual tree leaf inclination angles obtained using our method; (a3,b3,c3): scatter plots of the reference values obtained by manual measurements and the individual tree leaf azimuth angles obtained using our method.

Figure 9 .
Figure 9. Sensitivity analysis of the effect of the search radius on the segmentation accuracy of individual leaves for the three experimental tree plots: (a) Ehretia macrophylla; (b) crape myrtle; and (c) Fatsia japonica.The intersection of two lines in each figure represents the optimal results for the detected leaf number, indicating the value that was closest to the true value.The experiments show that an r value of one quarter of the average leaf width of each tree crown, and a  value of half of the average leaf width of each tree crown obtain the optimal individual leaf segmentation precision.

Figure 9 .
Figure 9. Sensitivity analysis of the effect of the search radius on the segmentation accuracy of individual leaves for the three experimental tree plots: (a) Ehretia macrophylla; (b) crape myrtle; and (c) Fatsia japonica.The intersection of two lines in each figure represents the optimal results for the detected leaf number, indicating the value that was closest to the true value.The experiments show that an r 1 value of one quarter of the average leaf width of each tree crown, and a MaxR value of half of the average leaf width of each tree crown obtain the optimal individual leaf segmentation precision.

Table 1 .
Overall accuracy assessment of the classification of different tree species.

Table 2 .
Statistics of leaf attributes estimate for the three individual trees using our method and manual measurement.

Table 3 .
Statistics for the overall results of the leaf morphological features (leaf area, leaf azimuthal angle and leaf inclination angle).