aTrunk — An ALS-Based Trunk Detection Algorithm

This paper presents a rapid multi-return ALS-based (Airborne Laser Scanning) tree trunk detection approach. The multi-core Divide & Conquer algorithm uses a CBH (Crown Base Height) estimation and 3D-clustering approach to isolate points associated with single trunks. For each trunk, a principal-component-based linear model is fitted, while a deterministic modification of LO-RANSAC is used to identify an optimal model. The algorithm returns a vector-based model for each identified trunk while parameters like the ground position, zenith orientation, azimuth orientation and length of the trunk are provided. The algorithm performed well for a study area of 109 trees (about 2/3 Norway Spruce and 1/3 European Beech), with a point density of 7.6 points per m, while a detection rate of about 75% and an overall accuracy of 84% were reached. Compared to crown-based tree detection methods, the aTrunk approach has the advantages of a high reliability (5% commission error) and its high tree positioning accuracy (0.59m average difference and 0.78m RMSE). The usage of overlapping segments with parametrizable size allows a seamless detection of the tree trunks.


Relevance
For a quantification of environmentally sustainable forest management an exact knowledge of various parameters is needed, e.g., the available number and distribution of trees as well as their species, the timber volume or the LAI (Leaf Area Index) [1].Remote sensing based single tree characterisation procedures have been proven suitable to enhance forest inventories regarding, for example, ecology, wildlife or biodiversity (e.g., [2,3]).Moreover, growth simulations require single tree models (with features such as height, crown size, tree species, spatial distribution) as input variables [4].In addition, single tree models have the potential to characterise features of biodiversity (for example, the identification of habitat trees [5]), the tree species [6][7][8] or the spatial distribution of trees.

State of the Art
Various approaches have been developed focusing on single tree detection using remote sensing data, while initially spectroscopy images and today mostly ALS (Airborne Laser Scanning) have been used as input data [9].ALS data has been proven suitable to characterise forest stands on a single tree level ( [9][10][11]).An airborne laser scanner is an active remote sensing system, which provides 3D (three-dimensional) information.Full-waveform and discrete-return ALS systems are distinguished, which are processed in the form of point clouds with additional intensity information [12,13].Single tree detection relies on pattern recognition methods, which take advantage of the tree phenotype and, if ALS is used, optionally the neighbourhood data of the laser returns [3].Jakubowski et al. [9] mention a general trend towards the usage of ALS alone and a more complex analysis with an increased accuracy of results.For a 2D (two-dimensional) identification and delineation, image-based approaches have been developed, for example watershed or valley-following segmentation techniques [14,15], the marked point process technique [16] or object-based image analysis [9].An adequate identification of suppressed trees relies on the analysis of full-waveform or laser-point data ( [11,12,17]).A 3D delineation is possible using voxel-based or vector-based approaches, for example multi-layer imaging techniques [18], variants of k-means clustering [19][20][21], adaptive 3D clustering [22], multi-scale clustering [23] or graph-based methods [2,12].Generally the 3D approaches result in improved detection rates in vertically inhomogeneous forest stands and more complex and realistic crown shapes [9] compared to raster-based approaches.

Related Work
The mentioned 3D tree delineation methods are focusing on a crown segmentation which mostly relies on a detection of tree tops based on local maxima (cf.[24]).Alternative tree identification approaches focusing on a tree trunk detection are rare.
Reitberger et al. [12] present a single tree detection approach, which includes an identification and modelling of tree trunks (ALS point density of 10-25 points per m 2 ).This method distinguishes trees using a normalized cut segmentation, which is capable of using already known tree positions to increase the accuracy.The initial tree positions are derived by a watershed segmentation method related to Vincent and Soille [25].Within the segments the CBH (Crown Base Height) is estimated by analysing the vertical histogram.The trunk points are separated within a segment between the CBH and a ground cover threshold.The trunks are identified by a 2D hierarchical clustering, and subsequently modelled by a RANSAC-based (cf.[26]) line fitting.Reitberger et al. [12] defines several rules to the trunk models-minimum number of tree points, minimum length, maximum angle-which shall keep the commission errors low.
The method of Abd Rahman et al. [27] uses high density ALS data (70 points per m 2 ) to identify single trees and model their trunks.Crown segments are identified using a watershed segmentation approach, while the vertical histogram of the points within a segment is evaluated to estimate the crown base height (reference height).The trunk of a segment is identified using a three-dimensional top-down trunk-growing algorithm.The 3D distance, the vertical histogram as well as the trunk diameter are evaluated to decide which points are assigned to the trunk.
The method of Lu et al. [28] takes advantage of the ALS intensity values to isolate points associated with trunks from leaf-off deciduous trees (10 points per m 2 ).The extracted points are used for a three-dimensional bottom-up trunk-growing process, while the assignment is controlled by the least 3D distance and a horizontal distance threshold.To avoid false positives, the trunk length and the maximum height of the lowest trunk point are controlled by threshold values.Furthermore the identified trunks are used for a 3D crown delineation.
Crown-based tree detection approaches usually reach an RMSE in positioning of 0.7 m up to above 3 m [11,24,29], which relies-besides the detection technique-on the tree height class (cf.[11]).In combination with the derived trunk positions Reitberger et al. [12] notices an increased accuracy in positioning by up to 25 %.

Objective
The aTrunk approach shall focus on a rapid ALS-based trunk detection without a dependency on a previous crown identification and shall work seamlessly on the entire observation area.The results shall be suitable to support large scale operational forest surveys, which makes a limited computation effort and data reduction crucial (caused by cost limitations [3]).The unique linear representation of a trunk promises a high reliability of the models, while potentially suppressed trees or dead trunks could be detected too.Applications of the approach could be to use the models as initial positions for crown segmentation approaches (e.g., [12,16]) or to confirm alternatively derived tree positions.An analysis of the trunk characteristics, regarding trunk length, trunk orientation (zenith, azimuth), point distribution in growth direction or spatial tree distribution could give additional information for characterising forest stands and biodiversity.

Study Area
The study area-illustrated in Figure 1

ALS Data
The multi-return ALS data used in this study was collected by the state forest service of Rhineland-Palatinate for the whole federal state, which specifies an absolute horizontal accuracy of 0.30 m and a vertical accuracy for open terrains of 0.15 m while a point density of at least 4 points per m 2 is guaranteed.The data was available in form of ASCII-files with a spatial expanse of 1 km × 1 km for each dataset.For the evaluation of the aTrunk approach, a subset of the study area of about 125 m × 116 m (with a extent of [368,112.5, 5,519,462.5, 368,237.5, 5,519,578.0] in EPSG:25832 system) and an overall point density of about 7.6 points per m 2 was chosen.Figure 1b illustrates the 3D structure of the study area.In addition, a full 1 km 2 dataset with the extent [368,000, 5,519,000, 369,000, 5,520,000] and an average point density of 7.7 points per m 2 was analysed.

Validation Data
The validation data collection aimed especially at the position of the trunks because this information can be easily used to estimate the detection rate of the algorithm and the accuracy in modelling the trunk positions.Moreover, ground measurements of the diameter at breast height (DBH) were taken for each trunk, which resulted in an average DBH of 0.255 m with a standard deviation of 0.049 m.
For a measurement of the trunk positions, a differential GNSS (Global Navigation Satellite System) of the type Topcon HiPer V (cf.Topcon Corporation [31]) was used.Because of a poor signal quality, the total accuracies of the measurements were quite low, with location differences clearly above 0.5 m.
To gain accurate information about the topology of the trunks, eight terrestrial laser scans (TLS) were taken with a Laser Scanner Photon 120 of the manufacturer FARO R , which measures at a wavelength of 758 nm and reaches a ranging accuracy of ±2 mm at a distance of 25 m (cf.FARO Europe GmbH [32]).Each of the all-around scans had a scan size of 8044 × 3446 pixels with the scanner-specific parameters: 1/5 resolution and 3× quality.The positions of the scans were chosen in such a manner that the study area was mapped completely.The alignment of the single scans to each other was prepared by placing reference spheres which can be identified in the post-processing software.
These TLS datasets were used to estimate the trunk position and diameter by applying a slicing approach (similar to [33]).The slices were 0.05 m thick, beginning at a level of 1.3 m up to the crowning height with a vertical distance of 1 m.To compare these models to the aTrunk positions, the positions at ground level were estimated by fitting a linear model to the slice centres of each trunk.

Preprocessing
The aTrunk approach relies on a height-normalised ALS point cloud using a digital terrain model.The delivered ALS data of the study area had already been classified into ground points and non-ground points by the provider.A Delaunay-triangulation of the ground points was used to calculate the height above ground for the non-ground points using CloudCompare [34].

Assumptions on Trunk Representation
ALS point clouds of different forested areas were visually inspected to define rules for a knowledge-based trunk modelling.The 3D view of the vegetation layer (green highlighted dots in Figure 1b) shows the complex crown structure and allows to identify single trunks because of a linear structure of the corresponding points.The following characteristics of a trunk mapped by ALS are suggested, which form the basis of the trunk model concept outlined in Sections 2.5.3 and 2.5.4: The LiDAR points associated with a trunk are: • widely spatially separable from the crown portion and the ground covering vegetation.
• moderately surrounded by points associated with branches, foliage or other objects.
• arranged in a straight line, which is oriented along the growth direction of the trunk.The maximum deviation from this line depends on the length of the trunk, for example, caused by irregular growth or branching.• largely uniformly distributed in growth direction of the trunk, which is substantiated in the spatial resolution of the ALS data.

Divide & Conquer
The large amount of points makes a direct analysis of the point cloud almost unmanageable.For this reason the idea of introducing a splitting step is used to divide the point cloud into subsets that are easier to handle, which is sketched in Figure 2b.In addition, this procedure facilitates a local estimation of the CBH, which allows a separation of points associated with trunks from other points (see Section 2.5.2).
The input point cloud is stored as an object called sample S, which allows a splitting into two subsamples, while a split is always done in an xy-direction.The size of a sample needs to be selected in such a way that it is small enough to estimate the CBH accurately, but large enough to include at least one tree.Each sample is separated into two new samples until the extent of the sample falls below a predefined threshold (maxSampleSize).This successive splitting allows a multi-core implementation of the algorithm in which it is possible to add a new CPU core to the program after each split.
Nevertheless, the splitting of the point cloud results in some disadvantages, which have to be dealt with.It is possible that points associated with a single trunk are separated from each other at random.This would result in an omission of trees.Therefore an overlapping area along the cutting edge of width overlap is suggested.The implementation of an overlap inevitably results in an additional analysis effort and the possibility of a multiple detection of a single trunk.These repeatedly detected trunks are merged, as described in Section 2.5.5.

Separation of the Trunk Section
The separation of the trunk section is based on the assumption that the potential trunk points P T S lie between points which are associated with ground-covering vegetation and tree crowns (illustrated in Figures 2c and 3b).A first threshold plane at the user-specified height Z GCV detaches the low-growing vegetation.For a sample with a maximum elevation z max it is assumed that the crown base height The CBH is estimated for each sample by analysing the vertical histogram (similarly to [12]).To become independent from the tree height, the points are divided into N layers horizontal layers between Z GCV and z max .The histogram is normalised by dividing the number of points of each layer by the total number of points, and filtered by applying a moving average filter of width 3 (Figure 3a).The estimated CBH corresponds to the highest intersection point where the histogram exceeds the threshold value th.To reduce the dependency of th from the number of layers, it is defined by th := thCBH N layers , while the parameter thCBH is set by the user.If the estimated CBH is located outside the range R CBH -for example caused by a vertically uniform distribution of the points which occurs sometimes for dead trees-it is guessed by ρ d • z max , while the proportion ρ d is estimated by the user for the study area.

Clustering
The clustering step is used to identify points associated with trunks in the point set P T S .It is assumed that several points which are spatially close together will form a trunk.Isolated points without or with just a few neighbours are assumed to be noise or sparse vegetation.The major problem is to identify point groups without knowing the number of clusters a priori (in contrast to, for example, a k-Means-based clustering approach).The proposed cluster definition-as a multi-dimensional extension of the DBSCAN [36] approach-takes advantage of the pair-wise spatial neighbourhood of points and the point density to solve this problem.
Cluster Definition: A cluster C p 0 ,δ,cM inP ts ⊆ P of a point cloud P ⊆ R N with dimension N is developed around the initialisation point p 0 using the threshold distance δ ∈ R + (see Equation (1)).To limit the influence of outliers, each point of the cluster has to have at least cM inP ts neighbours with distance δ.
For the purpose of this paper, this cluster definition allows to identify points associated with a single trunk P Ct out of P T S with N = 3 (illustrated in Figure 2d).Due to the high vertical variability of the ALS points, the z-coordinates are scaled by the parameter zBuf f erScale (similar to [23]).

Trunk Model
To get the desired orientation information of the trunks, a linear model is fitted to the points associated with the trunk P Ct .To receive a 3D vector-based regression model, a principal component analysis (PCA) is used.
Basic PCA-Model: A PCA of a point cloud P ⊆ R N with N ∈ N dimensions provides N pair-wise orthogonal vectors − − → P C i (with i ∈ 1 . . .N ) which are called principal component (PC) vectors.The PCs define a projected coordinate system, in which the PC scores (scores (P ) ⊆ R N ) of the points P correspond to the projected coordinates of these points in the new coordinate system (cf.Wold et al. [37]).
The first principal component t − − → P C 1 of the point cloud P Ct is oriented in the direction of the highest variance.For a nearly perfect linear alignment of the points associated with the trunk, the t − − → P C 1 vector should be oriented in the direction of the trunk.As a data preprocessing, a mean centring of the P Ct points is applied.So the trunk model corresponds to the t − − → P C 1 vector which is translated to the original centre point of the defining P Ct points.The model residuals correspond to the PC scores of the second and third principal component, so these are calculated by Equation (2) ( cf. p. 9985).Figure 4 illustrates the expected principal component model for a trunk.Unfortunately, the assumption of perfectly linearly arranged points with no kind of outliers cannot be maintained.Therefore an attempt was made to find one PC model which fits most of the points by excluding outliers.
Idea of Best Model Selection: One way to solve the problem of fitting a regression model to highly outlier-affected data is the RANSAC approach by Fischler and Bolles [26], which was also used by Reitberger et al. [12].The RANSAC approach is based on the assumption that outlying points can be identified by fitting multiple models which each rely on the minimum number of necessary points.
A basic assumption is that the points are randomly chosen and all models are independent from each other.In the case of a PC model, exactly two points are needed for each RANSAC sample.After this initialisation, all those points which have a residuum to the model below a specific threshold are assumed to support the model, while the others are assumed to be outliers.The model with the smallest proportion of outliers is accepted as the best model.
In this study the idea of the LO-RANSAC approach [38] was used, which optimises the outlier identification by running an additional model fitting for each RANSAC sample, using the model-supporting points.This technique relies on k ≈ log that the number of points associated with a trunk is insufficient in order to find enough independent points for the common LO-RANSAC approach.This is caused by the low number of points associated with a trunk (usually clearly below twenty) and an extremely high proportion of outliers (often caused by ramifications in the upper trunk section or by low-growing vegetation).To derive a deterministic model, every point pair combination of P T S is defined as an initial pair, whether or not the models are independent.Model Quality: The quality analysis of a trunk model is divided into two steps.In a first step, the validity of the model is tested, where already one fail results in rejecting the model.In a second step, the model with most model-supporting points is assumed to be the optimum.In case of the same point number, the model with the smallest mean squared error (MSE) is selected.For the validity check, some properties of each model are evaluated according to the modelling parameters summarised in Table 1 (see p. 9988).
A model is assumed to be valid if ...
• it contains enough points to ensure an accurate adaptation but unlikely false detections.⇒ |t P | ≥ n min • it contains only some points, because it is assumed that a high number of neighboured points is probably caused by leaves or branches.The value ρ corresponds to the point density of the sample.
• the ratio between the z-range (height) and xy-range (width) is comprehensible.⇒ range(t P z ) max(range(t P x ),range(t P y )) ≥ hwRel • the zenith angle of the trunk is reasonable.
• the model has a favourable ratio between model-supporting points and outliers.
• the points associated with the trunk are largely uniformly distributed in t Due to the assumption of residuals increasing with the length of the trunk, an adaptable threshold is needed to identify outliers.Therefore a length-dependent quality criterion-called M EP L (Maximum Error Per Length)-is proposed, which shall privilege large trunks in residual weighting.Equation (9) (p.9987) illustrates the calculation of this criterion, while the still unknown length of the trunk is approximated by the z-values of the points.
Figure 5 illustrates the modified LO-RANSAC concept for the PC-based trunk modelling approach, while Figure 2e sketches the modelling results within a sample.

Merge Duplicated Trunks
The split of the samples-done in step 2.5.1-with an intended overlap area can result in a multiple detection of single trunks in the overlapping separation section.This undesired effect is compensated for by an identification and subsequent merging of duplicated trunks whose centres are close together.As it cannot be ruled out that the models differ from each other (caused by a different point basis), all raw points are joined, while duplicates are eliminated.After the merging, a new trunk model is fitted just as explained in Section 2.5.4.All the derived trunk models are combined again in a result dataset as sketched in Figure 2f.

Principal Component Model
The PC model is based on a set of points associated with the trunk t P , whose selection was explained in 2.5.4.In addition, the points assumed to be outliers t O , which surround the trunk, are provided by the model.The three principal components t − − → P C 1 , t  3)) describes the angle deviation of t − − → P C 1 (the trunk) to a vertical line.
The azimuth angle t α ∈ [0, 360[ (Equation ( 4)) describes the deviation of t − − → P C 1 from the northern direction which indicates the leaning direction of the trunk.

Position
The centre point t c (Equation ( 5)) of the model corresponds to the spatial centre of the points t P associated with the trunk which were used to fit the model.
An advantage of the principal component line fitting approach is that the coordinates of every point of the trunk can be calculated as a linear combination of the first principal component t − − → P C 1 and the centre point t c .Therefore a point p k on this regression line can be calculated by Equation ( 6), in which the scalar k corresponds to the desired height along the growth direction of the trunk relative to the centre point c t : This feature is used to estimate the position of the trunk t p , while the z-component of the t p coordinate should be zero.Therefore the parameter k is estimated by Equation (7), in which t ζ corresponds to the zenith angle of the trunk and h to the height of the centre point t c .
The trunk top position t top corresponds to the modelled centre of the trunk at the assumed CBH height.It can be calculated corresponding to the ground position by using Equation ( 6), while the parameter k is set to the assumed height of the trunk (see Section 2.6.4).

Trunk Height
The height of a trunk t h is implicitly given by the expected crown base height Z CBH (see Section 2.5.2).The length of a trunk t l is calculated by Equation (7), in which the variable h is set to Z CBH .

Quality Criteria
The quality of a trunk model is evaluated by different criteria.Apart from the classical MSE, the MEPL (mentioned in 2.5.4) or the uniform distribution criterion (t χ 2 ) of the points associated with the trunk can be used.The t M SE is based on the residuals (t P ) of the regression line (Equation ( 8)), while t M EP L is defined by Equation ( 9).The uniform distribution criterion of the inliers is calculated by applying a Chi-Square Test on the pair-wise distance according to their sorted t − − → P C 1 scores.

Methods of Evaluation
To carry out an evaluation of the aTrunk approach, it was applied to the ALS point cloud of the study area.The measured GNSS positions and derived TLS trunk positions served as reference datasets.To minimise systematic errors both the GNSS and the TLS positions were adjusted to the detected positions using an affine point set registration.A detected trunk position was assumed to correspond to its reference position when its distance was below 4 m.In addition, the selection was limited to the extent of the reference datasets.It should be noted that both reference datasets did not include all trunks because some trees were not measured or not detected by the slicing approach.
A local maxima-based watershed segmentation-similar to the approach used by Hyyppä et al. [3])-served as benchmark, which identifies the trees using assumed water basins and dams.A raster model with 0.5 m × 0.5 m pixels, using the maximum point elevation for each pixel, was generated and smoothed by a 5 m × 5 m Gaussian-filter.The trees were identified by the watershed-segmentation-tool of SAGA (cf.[39]).For each segment the coordinates of the pixel with maximum value were considered as a tree position when the value was above 15 m.
To get information about the potential accuracy improvement of the aTrunk approach to the watershed approach, the matching positions (distance below 4 m) were identified.A combined dataset was generated using both the aTrunk and watershed positions.In case of corresponding positions the aTrunk positions were preferred.An accuracy assessment of the aTrunk, watershed, matching and combined positions was performed using the TLS positions as reference.Due to the large location residuals, the GNSS reference dataset was rejected for further evaluation.
A sensitivity analysis of the the aTrunk approach was performed in two steps.Firstly the parameters were optimised by applying the algorithm repeatedly to the study area with a large number of randomly selected parameter sets.The TLS positions served as reference to quantify the accuracy, while the positions were assigned automatically (with a maximum distance of 3 m).A parameter combination (illustrated in Table 1) with a high detection rate maintained for the study area.For the second step of the sensitivity analysis, the chosen set served as initial values, while each parameter was varied in its reasonable range of values.

Sensitivity Analysis
For the evaluation of the results of the sensitivity analysis, the parameters were grouped by their qualitative effect on the results (see Table 2).Selected parameters are illustrated in Figures 6-8 while other parameters can be found in the supplement.Group 1: Group 1 controls the computation effort and has-in a reasonable range of values-a low influence on the accuracy.The parameter overlap should be larger than the assumed average tree distance to ensure identifying all points of a trunk.
Group 2: Group 2 defines the demands of the user to the trunk model accuracy while stricter demands result in reduced commission errors or improved positioning accuracies but also in reduced detection rates.Figure 6a illustrates that the detection rate decreases the higher n min is chosen.Although the effect of n min on the commission errors is low for the study area, a positive effect is expected for more dense point clouds and larger n min values.The parameter M EP L (Figure 6b) leads to reduced commission errors for values below 0.1 m while the detection rate decreases clearly for values below 0.05 m caused by an overestimation of the number of outliers.The parameter maxZenith ensures finding almost vertical trunk models which results in a slightly optimised RMSE.A value above 5 • is recommended to reduce omissions.The parameter χ 2 , which ensures a certain uniform distribution of the point in growing direction, has shown to be unsuitable to reduce the commission errors.
Group 3: The values of Group 3 conform to forest characteristics, which are not known exactly a priori.Nevertheless the parameter mergeBuf f er is widely stable for about the assumed minimum distance between two trunks.The parameters ρ min , ρ max and ρ d -which support the automated CBH estimation-are estimated by expert knowledge.Assuming a correct automated CBH estimation, these parameters should have no influence on the results.The parameter Z GCV reduces the detection rate for higher values, because the number of available points in the lower trunk section is decreased.For the study area the minimum trunk height minZRange has a low influence on the commission error, while values above 5 m reduce the detection rate.Nevertheless it is assumed that for vertically inhomogeneous forest stands stricter demands on minZRange will minimise the commission errors.
Group 4: The parameters of Group 4 are necessary due to the clustering concept and are affected by the point density, vertical point distribution and noise.It is recommended to optimise the parameters δ (Figure 7a), zBuf f erScale (Figure 7b) and cM inP ts for each dataset.
Group 5: Group 5 controls the CBH estimation accuracy.The parameter maxSampleSize (Figure 8b) ensures a local estimation of the CBH, which results in an increased detection rate for small values, but also in increased commission errors.These incorrect detections are caused by trees, which are divided into small point groups.Such points are identified as a trunk, because the low horizontal extent results in a vertical arrangement similar to a trunk.A maxSampleSize value of about the average crown diameter is proposed.The threshold value thCBH has a great impact on the detection rate (Figure 8a), because for values chosen too large the trunk models reach into the crown and for small values the available number of points is reduced, which both result in inaccurate trunk models and omissions.For the analysis of thCBH the parameter ρ d was set to zero, because the threshold range defined by ρ min and ρ max as well as ρ d conceal wrong CBH estimations.

Evaluation
Table 3 illustrates the accuracy assessment using the TLS-derived trunk positions as reference.The assignment of the detected positions to the reference positions is visualised in Figure 9.The detection rate, precision and overall accuracy is calculated according to Yu et al. [13].The detection rate of the aTrunk approach is clearly outnumbered by the detection rate of the watershed approach.This is mainly due to the effect that not all trunks are recorded by the ALS data.In addition the homogeneous stand structure benefits the crown identification.The precision of the aTrunk approach is clearly better compared to the watershed results, because the point structure of trunks is less ambiguous compared to crown shapes which are partially merged.Nevertheless, other crown segmentation approaches could reach better results for this dataset.
The actual benefit of a trunk detection compared to a crown detection is the high accuracy in positioning, which significantly increases the positioning accuracy of the combined approach compared to the watershed approach.The improvement of the mean positioning error of about 32% supports the results of Reitberger et al. [12] which mention an improvement of up to 25%.The remaining residuals are mainly caused by branches and the trunk diameter.In addition uncertainties of the ALS positions (horizontal accuracy of about 5 to 15 cm assumed by May and Toth [40]) and of the TLS reference positions are remaining.Nevertheless the inaccurate GNSS measurements imply that the trunk detection can be superior to ground-based GNSS positioning.
As the trunk detection is independent from the crown detection, the matching positions can be used as particularly reliable tree positions.The remaining commission error is partially caused by incomplete reference positions because of covered trees in the TLS data.The zenith is right-skewed distributed with a median of about 3.5 • , while most of the trunks clearly differ from a vertical orientation.The circular histogram of the trunk azimuth implies a preferred east-west or west-east orientation of the trunks.Under the assumption of low systematic measurement errors, this effect could be caused by the stand characteristics (e.g., soil conditions, topography or dominant wind direction).

Model Performance
The objective of this paper was the development of a rapid algorithm which is capable of detecting single tree trunks independent from a crown identification.It has been shown that the number of false detections of the aTrunk approach is low.On a Dell Precision T1500 with an Intel Core i5-750 processor and 8 GB RAM, the algorithm reached a computation time of about 13 s (0.2 s data loading and 12.7 s analysis) for the study area and about 7 min (13 s data loading and 461 s analysis) for the larger 1 km 2 sample.
A sensitivity analysis of the effects of data quality (e.g., point density, recording angle), stand structure (e.g., tree species, suppressed trees, inhomogeneous vertical structure, tree number, canopy density), seasonal conditions and terrain (e.g., morphology, steep or varying terrain) has to be performed in a further study including an evaluation of the orientation information, to provide a tool which is capable of supporting operational forest surveys.The additional information about the trunk orientation could help to improve the characterisation of forest stands.

Model Concept
The presented approach is related to the trunk detection module of Reitberger et al. [12], so the overall concepts-ground covering vegetation cut-off, CBH estimation, clustering and line fitting-were adopted.The RANSAC concept has been replaced by a deterministic modification of LO-RANSAC using a length-dependent quality criterion, and the 2D hierarchical clustering has been replaced by a 3D variation of DBSCAN.Especially, the non-overlapping crown segmentation has been exchanged with partially overlapping rectangular segments allowing a seamless tree trunk detection and a multi-core computation.
Current drawbacks are caused by the simplified CBH estimation, insufficient rules to decide if a point group represents a trunk, the simplified removal of ground covering vegetation, and not explicitly considering the characteristics of dominated trees.The CBH estimation could be optimised by more advanced methods (e.g., [41,42]).The sampling technique with overlapping areas leads to an additional computation effort.An alternative method which joins the points associated with trunks to a single dataset may reduce the computation effort under the assumption that the clustering can be done efficiently.Visual inspections had shown that trees which slightly cross the border of a segment can result in incorrect detections, because the isolated points may be aligned in a line.The rules which are used to decide if a point group represents a trunk need to be evaluated and refined to get optimal results.Independently from the performed sensitivity analysis, the modelling concept implies that the approach depends on a high point density in the trunk section and clearly separable trunks, caused by the assumption of a linear arrangement of points associated with the trunks.So the approach is limited to forest stands with spatially distinguishable trunks where the canopy is permeable for the ALS beams.In addition low-branched trunks (like beech compared to spruce) and leaf off conditions are preferred, while a higher point density could optimise the trunk detection accuracy.

Conclusions
This paper presented a trunk detection approach which is independent from a crown identification and works seamlessly on the entire observation area avoiding non-overlapping segments.It derives reliable models (below 5% commission error), which are accurate in positioning (0.78 m RMSE).The models provide-next to the position-information about the trunk length, trunk orientation and point distribution (e.g., in growing direction).With a general focus on the computation effort, the algorithm is designed to be a potential tool for supporting operational forest surveys.
It has also been confirmed that the combination of the trunk positions with a crown-based tree detection is suitable to significantly improve the accuracy in tree positioning (32% for the study site).The independent detections result in a mutual confirmation of the tree positions which leads to decreased commission errors.In addition a preferred leaning direction of trunks was observed, which should be evaluated in further studies.
The algorithm has to be further evaluated regarding data quality, stand structure, seasonal conditions and terrain.It should also be evaluated whether the trunk detection is suitable to reduce the dependency of the tree size on the accuracy of detected tree locations which results in more accurate positions for taller trees [11].

Figure 1 .
Figure 1.Study area.(a) Study area near city of Hermeskeil, Germany with trunk positions predicted by a TLS-based (Terrestrial Laser Scan) trunk detection approach and background WMS-Service [30].Coordinates are specified in the EPSG:25832 system; (b) 3D ALS view of the study area.Yellow-coloured points correspond to ground and green points to vegetation.The brightness values correspond to the elevation.Red lines mark the finally modelled trunks.

Figure 2 Figure 2 .
Figure 2 illustrates the major steps of the presented aTrunk approach.The current implementation was realised in Python [35].As base data, the algorithm uses a height-normalised ALS point cloud in the ASCII-format.

Figure 3 .
Figure 3. Crown Base Height (CBH) estimation and trunk section separation within a sample.(a) The CBH is estimated by identifying the intersection point between the smoothed vertical histogram and a threshold line.Here N layers = 20 horizontal layers were used; (b) Separation of potential trunk points (green dots) from crown points or low-growing vegetation (red dots).

Figure 4 .
Figure 4. Three-dimensional one component PC trunk modelling concept with P Ct points (black dots) and centre point (red dot).

Figure 5 .
Figure 5. LO-RANSAC-based line fitting approach.Two possible raw trunk models (Figure 5a,b) using two randomly chosen points (blue dots), inlier threshold (red dashed lines), identified inliers (green dots) and outliers (red dots).Grey dots correspond to excluded points.Final model (Figure 5c) using inlier points t P (black dots), centre point t c (orange dot) and trunk position t p (red cross).(a) Low quality model in terms of RANSAC; (b) High quality model in terms of RANSAC; (c) Final trunk model.

Figure 6 .Figure 7 .Figure 8 .
Figure 6.Sensitivity analysis of selected parameters controlling the model quality (Group 2).(a) Minimum number of points n min ; (b) Quality criterion M EP L.

Figure 9 .
Figure 9. Evaluation of the detected aTrunk positions using the TLS positions as reference.The left image illustrates the assignment of the positions using absolute coordinates.The right image illustrates the relative positional differences between the reference positions and the corresponding detected positions.

Figures 10 and 11
Figures 10 and 11 illustrate the modelling results of the study area and the larger 1 km 2 sample using the same parametrisation.The homogeneously tall trees of the study area result in long detected trunks, with a mean value of 16.5 m and a standard deviation (SD) of 2.5 m.In comparison, the 1 km 2 sample consists of different areas with varying tree heights.This results in a high standard deviation of about 2.7 m.The MSE distribution of the points associated with the trunk is right-skewed for both datasets.A median of 0.15 m is reached for the study area and of about 0.06 m for the 1 km 2 sample.

Figure 10 .
Figure 10.Modelling results of the study area (86 detected trunks).

Table 1 .
Model Parameters

Table 3 .
Accuracy assessment of the different approaches.