ENHANCED DETECTION OF 3 D INDIVIDUAL TREES IN FORESTED AREAS USING AIRBORNE FULL-WAVEFORM LIDAR DATA BY COMBINING NORMALIZED CUTS WITH SPATIAL DENSITY CLUSTERING

A detailed understanding of the spatial distribution of forest understory is important but difficult. LiDAR remote sensing has been developing as a promising additional instrument to the conventional field work towards automated forest inventory. Unfortunately, understory (up to 50% of the top-tree height) in mixed and multilayered forests is often ignored due to a difficult observation scenario and limitation of the tree detection algorithm. Currently, the full-waveform (FWF) LiDAR with high penetration ability against overstory crowns can give us new hope to resolve the forest understory. Former approach based on 3D segmentation confirmed that the tree detection rates in both middle and lower forest layers are still low. Therefore, detecting sub-dominant and suppressed trees cannot be regarded as fully solved. In this work, we aim to improve the performance of the FWF laser scanner for the mapping of forest understory. The paper is to develop an enhanced methodology for detecting 3D individual trees by partitioning point clouds of airborne LiDAR. After extracting 3D coordinates of the laser beam echoes, the pulse intensity and width by waveform decomposition, the newly developed approach resolves 3D single trees are by an integrated approach, which delineates tree crowns by applying normalized cuts segmentation to the graph structure of local dense modes in point clouds constructed by mean shift clustering. In the context of our strategy, the mean shift clusters approximate primitives of (sub) single trees in LiDAR data and allow to define more significant features to reflect geometric and reflectional characteristics towards the single tree level. The developed methodology can be regarded as an object-based point cloud analysis approach for tree detection and is applied to datasets captured with the Riegl LMS-Q560 laser scanner at a point density of 25 points/m in the Bavarian Forest National Park, Germany, respectively under leaf-on and leaf-off conditions. The experiments lead to a detection rate of up to 67% for trees in the middle height layer and up to 53% for trees in the lower forest layer. It corresponds to an overall improvement in the detection rate of nearly 25% for forest understory compared to that obtained by the former method by extracting individual trees using normalized cuts segmentation solely. * Corresponding author.


INTRODUCTION
Laser scanning or LiDAR has been widely used in mapping the Earth's surface and especially in forest applications.Techniques for tree extraction from LiDAR data have been investigated for mapping forests at both plot and tree levels to identify important structural and biophysical parameters (Heurich, 2008;Korpela et al., 2010;Yao et al., 2012).Recent advances in LIDAR technology have generated new full waveform scanners that can trigger and record more backscattered pulses within the travel path of one laser ray, providing a higher spatial point density and additional information about the reflectional characteristics and vertical structure of trees (Stilla et al., 2007;Reitberger et al., 2008;Yao et al., 2010).
Figure 1 Mixed and multilayered forest Tree crowns are typically derived with the watershed algorithm, or by a region growing (Solberg et al., 2006) on the crown height model (CHM).New methods for single tree detection tackle conceptually the segmentation problem with a 3D approach instead of using only the CHM.In combination with full waveform data Reitberger et al. (2009) successfully demonstrated that the detection rate of single trees could be significantly improved in overall terms, especially in heterogeneous multilayered forest types (Figure 1) where groups of trees grow closely to each other.The fusion of 3D techniques with full waveform data seems to push the single tree approach to a new level of accuracy and completeness.Consequently, the estimation of tree shape parameters is enhanced using the 3D shape of segmented trees.Moreover, the analysis of the internal tree reflectional characteristics gains more insight into structure information which are significant for instance for tree species classification (Yao et al., 2012).
So far, little success has been achieved to identify individual trees in forest understory using LiDAR information.Former approach in Reitberger et al. (2009) also confirmed that the tree detection rates in both middle and lower forest layers are still low, although an improvement of up to 20% in the lower forest layers was found.Therefore, detecting sub-dominant and suppressed trees cannot be regarded as fully solved.In recent years there are certain authors who deal with the similar topics as forest understory monitoring by enhancing the tree detection method.Morsdorf et al. ( 2010) used airborne discrete-return LiDAR height and intensity information to identify individual vegetation strata on 5 m×5 m pixels in various forest conditions ISPRS Annals of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume II-5/W2, 2013ISPRS Workshop Laser Scanning 2013, 11 -13 November 2013, Antalya, Turkey and had little success in detecting the presence of the understory vegetation strata.Korpela et al. (2012) used airborne LiDAR data to study the understory trees by designing a conceptual compensation model for the transmission losses of laser pulses through overstory canopies.However, it was still an area-based detection and assessment method of the understory.Ferraz et al. (2012) has applied mean shift clustering to airborne LiDAR data of a multi-layered forest to extract single trees, assuming the spatial pattern of forest and boundaries of forest stratums is known in advance.They achieved a detection rate of 12.8% for the suppressed trees.
Full waveform LiDAR systems can overcome drawbacks of conventional laser scanners by detecting significantly more reflections in the understory forest strata, and providing the intensity and width of pulses as reflectional parameters.The objective of this paper is (i) to develop an enhanced approach that detects single trees for multilayered forests with an integrated 3D segmentation, (ii) to enable the new approach to utilizes the geometric and reflectional features derived for local dense modes object level) of point clouds (iii) to show how the detection and location of single trees across datasets of different properties are achieved using the developed approach.
The paper is divided into five sections.Section 2 focuses on the detection of single trees by combining normalized cuts with mean shift clustering.Section 3 shows the results which have been obtained from full waveform data acquired in the Bavarian Forest National Park.Finally, the results are discussed with conclusions in sections 4 and 5.

Decomposition of full waveform data
As usual, a single waveform is decomposed by fitting a series of Gaussian pulses to the waveform which contains N R reflections (Figure 2).
Figure 2 3D points and attributes derived from a waveform The vector ( , , , , ) as the 3D coordinates of the reflection.
Additionally, the points i X are given the width of the return pulse with i σ as the standard deviation and A i as the amplitude of the reflection i (Reitberger et al., 2009;Jutzi and Stilla, 2005).Note that basically each reflection can be detected by the waveform decomposition.
The sensor data are calibrated by referencing W i and I i to the pulse width e W and the intensity e I of the emitted Gaussian pulse and correcting the intensity with respect to the travel length s i of the laser beam and a nominal distance s 0 .
The correction assumes a target size larger or equal to the footprint (Wagner et al., 2006).The points from a waveform are subdivided into four point classes depending on the order of reflections within a waveform.

Mean shift clustering
Mean shift (MS) is a versatile tool for feature-space clustering.
MS has been successfully applied to image segmentation tasks by exploiting the spectral-spatial feature space (Comaniciu and Meer, 2002).As the feature-based analysis depends on the quality of selected features, the derivation of feature set play a fundamental role in design of a segmentation algorithm.Since we want to avoid the bias caused by deriving geometric features such as height textures, planarity and curvature caused by neighborhood selection, the 3D geographic space of forest stands spanned by i x ) , , ( coordinates of point clouds is chosen to explicitly represent the feature space.ALS point clouds convey a multimodal distribution in which each given mode defined as a local maximum in density correspond to a crown apex or a part of crown (Ferraz et al., 2012).MS vector is defined as where x is the center of the kernel (window), and h is a bandwidth parameter for the kernel.Given the function .   S I ={ I S }, the feature Z S is introduced as overall mean pulse intensity for the entire MS cluster, which uses the information provided for each point from waveform decomposition.The defined feature considered here is not depending on the pulse reflection types. S W ={ W S }, where W S is the overall mean pulse width for the entire MS cluster.However, the definition of this feature is limited to the single und first pulses, which could lead to a distinct broadening effect of pulses.

Normalized cuts segmentation
Within each sample plot the segmentation technique using normalized cuts (Shi and Malik, 2000) is used to extract point clouds associated to single trees based on a graph structure (Figure 5).This makes it possible to detect also smaller trees in the understory which cannot be indicated by local maxima in the CHM.This segmentation uses the positions (x i , y i , z i ) of the laser reflections within each MS cluster and optionally the pulse width W i and the intensity I i from the waveform decomposition.Additionally, stem positions or local maximums of CHM can be used as prior knowledge.The normalized cut segmentation applied to the local dense modes of point clouds is based on a graph structure G.The two disjoint segments A and B of the graph are found by minimizing the cost function:   as the total sum of weights between the segments A and B and , ( ,   as the sum of the weights of all edges ending in the segment A. The weights w ij specify the similarity between the MS clusters and are a function of the spatial distribution and features of MS clusters. A minimum solution for Eq.( 4) is found by means of a corresponding generalized eigenvalue problem.Consequently, the weighting matrix consists of following similarity functions between MS cluster i and j which is located within a cylinder of radius r XY around the MS cluster i: ) if ( ) ( , ) 0 otherwise e e e e D r w i j (5) The components P(i,j) and Z(i,j) weight the quadratic Euclidian distances between the MS clusters, where P(i,j) measures the horizontal distance S hp and the Z(i,j) is vertical distance S vp .The component F(i,j)  Figure 6 shows a sample area containing several coniferous trees of different sizes, which are grouped by normalized cuts based on the similarity between/ within MS clusters.The tree tops derived from the local maximums reasonably correspond to the reference trees (black lines), while a small understory tree is also correctly separated.Although 3D segmentation approach is not dependent on full waveform LiDAR data, the application to conventional LiDAR data just providing 3D point coordinates could undermine the performance, since the pulse intensity and width of MS clusters play an important role in constructing the similarity matrix controlling the graph cuts.

Control parameters
The MS clustering has only one control parameter -kernel bandwidth (h v , h h ) needed to be set in advance h = h v = 2.4m), while the normalized cut segmentation is controlled by several control parameters whose values can be optimized in experiments.Firstly, the most important parameter Ncut thres , which controls the subdivision (=cut) of a graph G, was set equal to 0.18.Moreover, a graph G is no longer subdivided if the number of MS clusters of the graph undershoots the limit of 2 clusters.Another important factor ADJ_RADIUS that define the radius of adjacencies of nodes in the graph structure is set to 9.7m.Finally, we used the empirical values σ f = 0.5, σ xy = 3.15m, σ z = 11.0m and σ G = 3.5m to control the influence of the impact factors F(i,j), P(i,j), Z(i,j) and G(i,j).The value for σ z is larger than σ xy assuming that the tree height is larger than the tree crown diameter.

Material
Experiments were conducted in the Bavarian Forest National Park which is located in south-eastern Germany along the border to the Czech Republic (49 o 3' 19" N, 13 o 12' 9" E). 2 sample plots with an area size between 1000 m 2 and 3600 m 2 were selected in the mixed mountain forests.The plots comprise forest in the regeneration phase and the late pole phase(Table 1).The test sites have suffered from tree disease due to bark beetle attack.Reference data for all trees with stem diameter larger than 7 cm have been collected in May 2006 and 2007 for 93 Norway spruces (Picea abies), 190 European beeches (Fagus sylvatica), 50 Fir (Abies), 24 Sycamore maples (Acer pseudoplatanus), 9 Norway maple (Acer platanoides) and 2 Tilia.Several tree parameters like the DBH, total tree height, stem position and tree species were measured and determined with the help of GPS, tacheometry and the 'Vertex' III system.A DTM with a grid size of 1 m and an absolute accuracy of 25 cm was available (Heurich et al., 2008)

Calibration
The calibration of the Riegl full-waveform system was performed by special calibration flights over an airfield.Several tracks were flown at different flying heights (200 m and 400 m) along and across the airfield.The mean intensity I i , corrected with respect to the emitted intensity I e , and the mean run length s i were calculated in four homogeneous areas (122 m 2 -133 m 2 ) for each track.

Flight 2006 Flight 2007
Calibrated parameter k 1.902 1.736 Table 3 Estimation of calibration parameter k According to Eq. ( 2), the best coefficient k was estimated from all possible observation equations which can be formulated for two tracks i and j flown at different heights.Table 3 shows the results obtained for the two flights of data sets.

Results
The enhanced procedure for 3D single tree detection was applied to sample plots in a batch procedure without any manual interaction.Tables 4 and 5 summarize the percentage of detected trees for the entire sample plots under different foliage conditions.The plots are subdivided into 3 forest layers with respect to the mean height h top of the 100 highest trees per ha.The lower layer contains all trees below 50 % of h top , the intermediate layer refers to all trees between 50 % and 80 % of h top , and, finally, the upper layer contains the rest of the trees.The tree detection results were evaluated by matching with single trees in reference data using two criterions: i). the distance of detected trees should be smaller than 60% of the mean tree spacing of the plot; ii) the height difference between detected and reference trees should be smaller than 15% of h top .
If a reference tree is assigned to more than one tree position, the tree position with the minimum distance to the reference is selected.Detected trees that are liked to one tree position are "detected trees" and detected trees without any link to a reference tree position are treated as "false positives".The overall detection rate of up to ca. 65 % can be achieved for the leaf-on case using the newly developed approach, while roughly a third of the detections correspond to false alarm, which indicates a moderate reliability.It can be seen that most of the trees are detected in the upper layer.In comparison, the detection rates in the intermediate and lower layer are smaller.However, the enhanced tree detection method based on 3D segmentation combining mean shift with normalized cuts has fully exploited advantages of full waveform data to detect supressed trees in the both middle and lower layers.It can be seen that the tree detection rate in the both forest layers have been improved by 25% compared to that obtained by the former method presented in Reitberger et al. (2009).The improvement seems to be uncorrelated to the foliage condition, since constant improvements between 10% and 25% are observed for all three forest layers.Additionally, if we compare dataset I (leaf-off) to dataset II (leaf-on) it can be addressed that the foliage condition could affect the improvement of the detection rate, but not significantly.As unexpected, the overall detection rate is worse by ca.1% in leaf-off situation, while the false positive has increased by 2%.It was caused by the fact that dense tree crowns in leaf-on condition and dead woods emerging in the meantime could benefit the detection of dominated trees in the overstory.Furthermore, we show the distribution of the detected trees from two methods over different stem diameters in Figure 7. Therefore, the capability of new 3D tree segmentation with respect to detecting supressed understory trees in a multilayered forest can be highlighted here.

Mean position error Leaf-off
Leaf-on Former method 1.53m 1.48m New method 1.61m 1.71m Table 7. Accuracy of the tree position determination Finally, Table 7 shows the absolute positional accuracy of the trees detected by two methods.The mean positioning error of detected trees using the new method gets worse by ca.5-10%, which corresponds to 10 cm in the leaf-off case.Namely, the overall accuracy of tree position determination is still better in the case of former tree detection method.

DISCUSSION
According to Rutzinger et al. (2008), the developed methodology can be regarded as an object-based point cloud analysis approach for tree detection.Conceptually, the presented approach to detect individual trees from airborne LIDAR data of forest areas goes one step further than 3D tree segmentation at the point level.This paper presents an enhanced scheme for detecting single trees using high density full waveform LiDAR data in a multilayered forest.In this study, based on the graph structure of MS clusters derived automatically in advance the detection of trees using waveformlaser measurements is performed in the context of global grouping.The computational costs needed for the bipartition of the affinity matrix used in normalized cuts can be essentially reduced due to the small number of generated mean shift clusters representing the graph nodes.The 3D segmentation for single tree detection by combining normalized cuts with mean shift could improve the detection rate in the lower forest layer by averagely 25% compared to the former 3D segmentation presented in Reitberger et al. (2009).Moreover, the leaf-on case seems to provide slightly better results for single tree detection in our experiment than leaf-off case.It could be caused by mixed reasons.Mostly important, under leaf-on condition the crown of overstory deciduous trees exhibits a more abundant type, which leads to stronger spatial and biophysical evidences for differing from coniferous trees making the detection easier.
When viewing at the results of tree detection in forest understory we observed an improvement of up to 25% in the detection rate for datasets under different foliage conditions.It is interesting to see that the improvement is almost independent on foliage condition for all the three forest layers, while the lower forest layer showed the largest improvement that is 10% better than the higher forest layer.It can also be retrieved from Figures 7 that low trees (with small stem diameter) in multilayered forests have a greater potential to be better resolved by full-waveform laser scanners when using a more advanced data handling technique.Overstory trees still is a key factor in preventing the laser pulse from reaching or transmitting through supressed trees in the understory, the tree detection rate for the middle and lower forest layers in the leafoff case has increased by 5% compared to leaf-on situation, and it does not much depend on the tree detection algorithm used.
However, false alarm rate is still not very low, and it should be a kind of over segmentation problem.Points of tree crown and ISPRS Annals of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume II-5/W2, 2013 ISPRS Workshop Laser Scanning 2013, 11 -13 November 2013, Antalya, Turkey branch which actually belong to adjacent trees undermine local dense structure and weaken the feature functions used in the graph-based segmentation, even leading to error detections.Such false points could happen to the tree segments owing to incompleteness of reference data acquired by filed work and dense forest stands.Tao et al. (2007) showed that the integration of mean shift and normalized cuts can help to improve the image segmentation results and reduce the oversegmentation effect as well.For our case, it is stated that the combined utilization of two mutually enhanced techniques detected more single trees especially in forest understory without producing more false positives.It can be attributed to the reason that the mean shift clusters approximate primitives of (sub) single trees in LiDAR data and allow to define more significant features to reflect geometric and reflectional characteristics towards single tree level.The similar improvement caused by synergic effect could also be witnessed in the application of urban object detection using LiDAR data (Yao et al., 2009).

CONCLUSIONS
The study presents an enhanced scheme for detecting single trees in mixed and multi-layered forest from full-waveform LiDAR data based on an integrated 3D segmentation method.
The developed methodology can be regarded as an object-based point cloud analysis approach for tree detection and an extension to the former work based on normalized cuts solely.
The results attained in a heterogeneous forest of different foliage conditions show that the overall detection rate of single trees could be improved by 15%, being less dependent on foliage condition, while the detection rate of depressed trees in the forest understory has increased by ca.25%.Future research could be focussed on assessing the influence of tree species on the detection results, since the mean pulse intensity and tree canopy geometry are clearly different for coniferous and deciduous trees according to previous study.Furthermore, extended features of MS clusters need to be examined, whether providing direct clues to reduce the false alarm rate.

Figure 3
Figure 3 Cylindrical shaped kernel for density estimation with horizontal Gaussian profile

Figure 4
Figure 4 Mean shift clusters of an exemplary sample plot 2.2.3 Feature derivation for mean shift clusters Deriving significant features describing each tree/sub-tree individually is a key step in detecting single trees using graphbased segmentation.This part is specifically to deal with feature derivation for constructing the affinity matrix in normalized cuts to measure the similarity function between each connected graph node.The mean shift spatial clustering in LiDAR point clouds provides for each segmented (sub) tree cluster corresponding laser points with associated properties.The four-group feature set S = {S hp , S vpr , S I, S W } can be defined to reflect geometric and reflectional properties.The features will be defined as follows, respectively:  S hp ={ X S , Y S }, which records the horizontal (x,y) position of MS clusters, respectively  S vp ={ Z S }, which records the vertical position of MS clusters

Figure 5
Figure 5 Schematic representation of graph-based segmentation using normalized cuts based on MS clusters describes the quadratic Euclidian distance between mean LiDAR intensities S I and width S W , while G(i,j) records the maximum horizontal distance of two MS clusters to the local maximum of CHM.The position and height of a segmented tree is defined as the centroid and maximum height of the tree segment from the point cloud.The explicit definitions of the four functions are same to those in Reitberger et al. (2009).(a) (b) Figure 6 Single tree segmentation: from MS clusters (a) to 3D singles trees (b) based on global grouping by normalized cuts.MS clusters Similarity W ISPRS Annals of the Photogrammetry, Remote Sensing and Spatial Information Sciences, Volume II-5/W2, 2013 ISPRS Workshop Laser Scanning 2013, 11 -13 November 2013, Antalya, Turkey

Figure 7
Figure 7 Distribution of stem diameter vs. tree detection improvement, top: leaf-off season bottom: leaf-on season

2.2 Singe tree detection 2.2.1 Local maximal
stem detection method to further detect sub-dominated trees which are not represented by local maximums, when sufficient stem points are available.

Table 1 .
. Characteristics of sample plots Full waveform data have been collected by Milan Flug GmbH with the Riegl LMS-Q560 scanner in May 2006 after snowmelt but prior to foliation and in May 2007 after foliation with an average point density of 25 points/m 2 (Table2).The vertical sampling distance was 15 cm, the pulse width at half maximum reached 4 ns and the laser wavelength was 1550 nm.The flying altitude of 400 m resulted in a footprint size of 20 cm.