An Automated Approach for Extracting Forest Inventory Data from Individual Trees Using a Handheld Mobile Laser Scanner

Many dendrometric parameters have been estimated by light detection and ranging (LiDAR) technology over the last two decades. Handheld mobile laser scanning (HMLS), in particular, has come into prominence as a cost-effective data collection method for forest inventories. However, most pilot studies were performed in domesticated landscapes, where the environ mental settings were far from those presented by (near)natural forest ecosystems. Besides, these studies consisted of numerous data processing steps, which were challenging when employed by manual means. Here we present an automated approach for deriving key inven tory data using the HMLS method in natural forest areas. To this end, many algorithms (e.g., cylinder/circle/ellipse fitting) and machine learning models (e.g., random forest classifier) were used in the data processing stage for estimation of the tree diameter at breast height (DBH) and the number of trees. The estimates were then compared against the reference data obtained by field measurements from six forest sample plots. The results showed that correlations be tween the estimated and reference DBHs were very strong at the plot level ( r =0.83–0.99, p <0.05). The average RMSE for tree DBHs was 1.8 cm at the forest landscape level. As for tree detection, 92.5% of 292 trunks were correctly classified on point cloud data. In general, estimation accuracy was sufficient for operational forest inventory needs. However, they could markedly decrease in »hard plots« located at rocky terrains with dense undergrowth and ir regular trunks. We concluded that area-based forest inventories might hugely benefit from the HMLS method, particularly in »easy plots«. By improving the algorithmic performances, the accuracy levels can be further increased by future research.


Introduction
Forest management planning requires accurate and updated information for characterizing the current state of forest ecosystems. Periodical forest inventories are the primary data sources for this information flow (Kangas andMaltamo 2006, Ozkan andDemirel 2018). In inventory surveying, the number of trees and the diameter at breast height (DBH) are two essential parameters since they form the basis for both forest stand density and timber volume calculations (Wan et al. 2019). Unless having these data, neither the sustainable harvest rates nor the revenue of forest enterprises can be determined accurately (Bettinger et al. 2017, Bulut et al. 2016, Vatandaşlar and Zeybek 2020.
Conventional data collection methods (i.e., field measurements) are usually expensive, time-consuming, and labor-intensive in forest inventory surveying (Trotter et al. 1997). Therefore, remote sensing is widely used individually or combined with field measurements (Forsman et al. 2016, Gómez et al. 2019, Ozdemir and Karnieli 2011, Ucar et al. 2018. Laser scanners, in particular, are more common with the rapid development of light detection and ranging (LiDAR) technology worldwide (Balenović et al. 2019, Gómez et al. 2019, Hyyppa et al. 2008, Oveland et al. 2018 Croat. j. for. eng. 42(2021)3 2017, Wang et al. 2013). A comprehensive review, focused on different platforms of LiDAR (e.g., spaceborne, airborne, terrestrial, etc.), can be seen in (Van Leeuwen and Nieuwenhuis 2010). Thanks to this technology, forest managers can quickly collect forestrelated data on the ground at a reduced cost (Cabo et al. 2018, Del Perugia et al. 2019, Ryding et al. 2015. Static LiDAR platforms have some limitations for operational forest inventory, such as the occlusion effect, low speed of data acquisition, and platform weight. Therefore, mobile laser scanners (MLS) have become preferable amongst forestry professionals recently (Bauwens et al. 2016, Hyyppa et al. 2020, Gollob et al. 2020. Given the carrier platforms, MLS can be grouped as: (i) vehicles, (ii), backpacks and (iii) handheld systems. The handheld mobile laser scanner (HMLS) systems, in particular, seem to be promising in forest inventories because the surveying time can be largely reduced thanks to their high mobility and lightweight (Balenović et al. (2021)). Indeed, one operator can collect data from trees at centimeter-level accuracy by simply walking through the forest with an HMLS at hand. Moreover, the dependency on GNSS has been eliminated in the latest systems, thanks to the Simultaneous Localization and Mapping (SLAM) algorithms. Thus, degraded GNSS signal under the canopy is no longer a problem in the forest. With the SLAM-based HMLS, point clouds can be automatically produced in the local coordinate system. Then, 3D data may be georeferenced to an absolute position. Readers are encouraged to see the compherensive study by Balenović et al. (2021) for more-detailed information on HMLS systems.
Using a SLAM-based HMLS, Gollob et al. (2020) successfully detected individual trees and estimated their DBHs in various forest conditions. They also used a density-based algorithm for the automatic detection of tree trunks. The researchers compared the HMLSderived estimates with terrestrial laser scanner data and manual field measurements. In another study, Hyyppa et al. (2020) compared the accuracy of five LiDAR-based systems for operational forest inventory. They found that the HMLS and drone-mounted laser scanner (under-canopy) were superior to other systems for DBH and stem curve estimations. Similarly, Jurjević et al. (2020) estimated tree heights using HMLS and a drone-mounted laser scanner (over-canopy) in deciduous forest stands. The researchers compared the estimation results with traditional ground-based measurements. Their results showed a high level of agreement among the data sets. They stated that HMLS was a useful and reliable method for tree height measurement, especially for easy forest conditions in leaf-off season. More recently, Balenović et al. (2021) examined HMLS studies conducted in the forest inventory field. In this study, state-of-the-art applications of HMLS were reviewed, and the pros&cons of the method were documented in detail. They concluded that further research was needed to test the HMLS method potential for practical use in forest inventories.
Despite promising results from previous research (Bauwens et al. 2016, Gollob et al. 2020, Hyyppa et al. 2020, there is still a knowledge gap regarding the practical use of the HMLS, particularly in complex forest conditions. For instance, the extraction of individual tree parameters (post-processing) may still be challenging in dense and heterogenic plots located on poor sites. Researchers have examined different approaches (e.g., circle/cylinder fitting, etc.) to reach the best estimation accuracy for individual tree parameters derived by HMLS data (Koreň et al. 2017, Liu et al. 2018, Wu et al. 2018, Zhong et al. 2017, Zhou et al. 2019. Unfortunately, the accuracy requirements have not been fully met since trees generally show irregular morphological elements, such as noncircular trunks and complex branches. Therefore, these elements are poorly modeled in HMLS-based forest inventories. Indeed, more robust approaches are needed for automatically classifying and extracting tree trunks in »hard forest plots«. In this regard, the use of machine learning classifiers with untested morphological models (e.g., ellipse fitting for DBH estimation of a noncircular tree trunk) may solve the problem. This approach can minimize the estimation errors encountered in the post-processing stage. Thus, it needs to be tested in unmanaged natural forests presenting a non-homogeneous structure (Gadow et al. 2012).
The aim of this study is to develop a new and automated approach for extracting the number of trees and tree DBH data from 3D point clouds captured by a GeoSLAM ZEB-REVO HMLS device. Data accuracy was assessed based on comparing the extracted information with field measurements. To this end, 292 trees were scanned on six forest sample plots in an unmanaged natural forest. Several algorithms and shape fitting methods (i.e., circle, ellipse, or cylinder) were examined for improving data accuracy. The proposed approach is expected to promote the automation process in digital forest inventories, as well as to enhance the algorithmic performance in the data processing. Thus, forest inventory surveying is likely to be more cost-effective, especially in the mountainous countries with (near)natural forests. This will allow field personnel to save more time and money.

Study Area and Data Acquisition
The study area comprised six forest sample plots located at Artvin Province in the northeastern part of Turkey (Fig. 1). Artvin is the most mountainous city in the country. It is also in Caucasus Biodiversity Hotspot; thus, it hosts one of the highest rates of ecosystem diversity (Manvelidze et al. 2009). The main tree species in the forest are Picea orientalis, Pinus sylvestris, and Abies nordmanniana, along with Quercus dschorochensis, Fagus orientalis, and Castanea sativa. The average annual rainfall was 753 mm, while the average air temperature was 12.4 o C in the period 1949-2018 (SMS 2019). Forest sample plots were designated in different shapes and sizes. Detailed information on each plot is given in Table  1. Field surveys were carried out on October 1 and 2, 2018. First, essential data for forest inventory (the number of trees, DBH, slope degree, etc.) were measured and recorded into inventory sheets according to the conventional ground measurement methods described in the national guidelines (GDF 2017). At this stage, the lower tree DBH threshold was taken as 8 cm following GDF (2017). In the next stage, forest sampling plots were scanned using the ZEB-REVO HMLS device. Scanning time was 2-3 mins for a typical forest plot with 400 m 2 in size. The survey path was planned as a closed-loop, as suggested by Del Perugia et al. (2019). During surveying, the surveyor started to walk from the plot center to the border and turned back several times. The starting point was also the endpoint of the path. In this way, the same area was scanned multiple times resulting in an increased scanning alignment. Finally, scanned LAS data were recorded on a personal laptop as the zip file for data processing.

Instrumentation
The point cloud data was acquired with the ZEB-REVO HMLS device (Cadge 2016) (https://geoslam.com) ( Fig. 2). The device consists of a 2D time-of-flight laser scanner, as well as an integrated inertial measurement unit (IMU) on its rotary engine. Thanks to the lightweight (less than 1 kg), it can be easily used by one person at hand during labor-intensive forestry surveys. Furthermore, SLAM technology allows for generating a map of its surroundings and bear itself properly within the map. Thus, the map comes with a relative coordinate system. As a result, almost all dendrometric parameters (e.g., height, DBH, timber volume, etc.) can be modeled via 3D reconstructions. Technical specifications for the ZEB-REVO HMLS are presented in Table 2.

Methods
The methodology was mainly based on improving the estimation accuracy of dendrometric parameters collected by the ZEB-REVO HMLS using various algorithms. Then, the data were compared against the reference, conventionally measured on-the-ground. It should be noted that tree positions are unnecessary for most forest inventory studies since aggregated data are assessed at the plot-, stand-, or landscape levels. Therefore, only the plot centers were located using a handheld GPS, but individual trees were not registered into the real coordinate system. The general workflow can be seen in Fig. 3.

Pre-processing Stage
Raw data collected from the forest sample plots were pre-processed using GeoSLAM Hub desktop software (ZEB-REVO 2019). The SLAM algorithm was used to combine 2D coordinates with IMU data to generate 3D point clouds (ZEB-REVO 2019). The algorithm uses a method similar to the traverse technique (di Filippo et al. 2018) to convert raw data into point clouds. Unclassified point clouds were then clipped using the borders of forest sample plots based on coordinates of plot centers.

Data Processing and Modeling
As for ground/non-ground classification, the cloth simulation filtering (CSF) algorithm (Zhang et al. 2016) was implemented using the RCSF package in the R programming language (Roussel and Qi 2018). CSF uses a physical movement of the cloth model on the 3D point cloud. Optimal parameters were empirically tested for filtering. Accordingly, 0.2, 0.1, and 2.0 values were assigned for the class threshold, resolution, and rigidness parameters, respectively. These parameters may differ depending on landscapes because of the nature of the  filtering algorithm. After filtering, points were classified as ground and non-ground in a LAS file.
As for normalization, an accurate digital terrain model (DTM) is essential. Tree heights are calculated based on DTM. To this end, DTM was created using the knnidw algorithm, which is a spatial interpolation method in the lidR package (Roussel and Auty 2019). It interpolates the ground points and creates a regularized DTM. Then, non-ground points were subtracted from the DTM using Eq. (1).
Where: Z AGL distance from above ground points Z NG above ground points Z G ground points projected on DTM.
Finally, point cloud data were cross-sectioned between 1.0 and 2.0 meters above ground level (AGL). Thus, data density was reduced, and DBH estimations could be determined. However, some outliers might still exist in the cross-sections. Therefore, the statistical outlier removal (SOR) algorithm (Rusu 2010) was ap-plied for noise reduction. Irregularities are another limitation along the tree trunks. In order to solve this problem, the moving least squares (MLS) algorithm and VoxelGrid up-sampling were implemented (Alexa et al. 2003). Graphical representation of the cross-sectioning stage (a.k.a. data slicing) can be seen in Fig. 4.

Machine Learning Algorithm and Tree Extraction
After slicing the data, tree trunks need to be detected. This is a difficult task, particularly in complex forest plots with dense noise. The binary classification model, trained with Random Forest (RF) machine learning algorithm, can automatically detect trees within the point cloud data. In this way, tree DBH information can be extracted from classified data. It is a fundamental step for the automatic DBH estimation process.
To this end, the machine learning technique was performed based on geometric structures and the adjacency of point clouds. Thus, the points representing tree trunks were automatically extracted. Trunks were classified and labeled using the RF algorithm (Breiman 2001, Ni et al. 2017. For this, Caret package was used In the next step, surface normal and curvature, which are the two essential geometrical features, were estimated in the point cloud. The estimation of surface normal depends on the calculation of the approximate plane at each point position. The vector of the surface normal was then calculated from the eigenvectors and eigenvalues of the covariance matrix as described in (Montgomery et al. 2006). These elements allowed for clustering trunks, leaves, and branches from each other. Moreover, five different eigenvalue-based characteristics were used. These were omnivariance, planarity, linearity, surface variance, and anisotropy change (Ni et al. 2017). Regarding model accuracy, the k-fold cross-validation method was used for dividing data sets into k-subgroups. For this, the value of k was set to 10.
Finally, the Euclidean clustering extraction (ECE) algorithm (Rusu 2010) was applied for improving the fitting accuracy. Cluster analysis was used to classify and group a set of similar points that meet the threshold value in terms of distance. This approach, as a whole, was first applied to HMLS data for the forest inventory purpose.

Geometrical Shape Fitting
Given each cluster number, shape fitting is applied to the data that are separated as »tree« and »non-tree« trunks. Here, non-tree trunk points were excluded from the analysis in order to estimate tree DBHs robustly. For better DBH estimates, three morphological models, i.e., cylinder, circle and ellipse, were tested on each tree trunk. As for cylinder fitting, the analyses were performed in 1-m-thick sections between 1 m and 2 m AGL. For circle and ellipse fittings, 5-cm-thick sections were selected between 1.28 m and 1.33 m AGL (see Fig. 4).
The random sample consensus (RANSAC) technique was used for cylinder fitting (Nurunnabi et al. 2017, Schnabel et al. 2007. It was an iterative algorithm for estimating the parameters of a cylinder model. The maximum number of iterations was determined by Eq. (3): Where: N maximum number of iteration p probability value in statistical test (taken as 0.99) n minimum number of points for cylinder fitting w probability for choosing an inlier.
As for circle fitting, Izhak Bucher's method ( Bucher 1991) was used. The algorithm was developed in R programming language. The formula can be seen in Eq. Where: x c and y c coordinates for circle center r circle radius x and y coordinates for trunk points.
When Eq. is expanded using the parameters (a 0 , a 1 , a 2 ), Eq. is obtained as follows: x 2 + y 2 + a 0 x + a 1 y + a 2 = 0 Where: Circle parameters (x c , y c , r) were obtained by matrix solutions applied for trunk points.
Finally, the Conicfit package was used to ellipse fitting (Gama and Chernov 2015). The package uses Eq. (6) for fitting ellipse shape on trunk points: Where: k = (A,B,C,D,E,F) T the vector of parameters to be estimated. The ellipse fitting method is also known as the orthogonal distance regression (ODR) (Al-Sharadqah and Chernov 2012, Gander et al. 1994). The ODR algorithm solves the minimization problem using Eq. (7).
Where: d i orthogonal distance from n data points (x i , y i ) to ellipse P(x, y; k) = 0.
For decision-making between the ellipse and circle, the circle fitting method was first performed on each trunk data. Then, the distribution of distances to the center was calculated. Direction-dependent densities from the Y-axis were tested using a histogram. If the histogram is bimodal (see Fig. 5b), the first condition is met for the ellipse fitting. As for the second condition, the median absolute deviation (MAD) values are used. MAD is calculated based on the distances using Eq. (8): Where: c constant-coefficient used as 1.4826 x distance from center coordinates median(x) centroid.
The threshold value for MAD was set to 3 in this study. However, it may vary depending on tree species in different forest settings. If the MAD value is higher than the threshold, then the second condition is met, and ellipse fitting is permanently chosen, as shown in Fig. 5c. Pseudo-codes used in the decisionmaking process can be seen in Algorithm 1. All the equations and algorithms presented in this subsection were written in the R programming language by the authors.
Algorithm 1 Decision-making procedure for using a circle or an ellipse fitting method

Accuracy Assessment and Correlation Analysis
Estimation accuracy was assessed both for machine learning classification and single-tree detection. As for machine learning, k-fold cross-validation (Lantz 2015) was performed. For this purpose, point cloud data was partitioned into two groups: training (70%) and test (30%) samples.
As for tree detection accuracy, a confusion matrix was generated. True positive, false positive, and false negative values were documented by visual interpretation for true detected, wrong detected, and undetected trees, respectively. Accordingly, recall, precision ratios, and F-scores were calculated using Eq (9)-(11): Finally, scanned (estimation) and measured (reference) DBH data were subjected to Pearson's correlation analysis. In this step, only TP trees were used against the reference data.
All these processes were implemented in R and C++ programming language (Gama and Chernov 2015, Kuhn et al. 2019, Roussel and Auty 2019, Rusu 2010, Team 2019) using a modest laptop with Intel Core i5-3210M (2.5 GHz) processor.

Individual Tree Detection and DBH Estimation
As for the number of trees estimations, it was found that F-scores from the six sampling plots differed between 0.78 and 1.0 with a mean of 0.93 (Table 3). The best estimation was obtained in Plot 2, which was a clean forest site with a canopy cover of about 50%. All trees were correctly detected on this plot. The plot had no undergrowth such as shrubs, seedlings or lying deadwoods on the forest floor. These kinds of elements often create noise on the point cloud, which may cause misclassifications. The worst estimation, on the other hand, was in Plot 4, in which almost half of the trees were wrongly detected (Fig. 6). This plot was on a poor site hosting many broadleaved trees on rocky terrain. The irregularities on low-quality trunks led to misclassified results. Similarly, Carr and Slyder (2018) found a mean F-score of 0.93 for a mixed temperate forest in Fig. 5 Decision-making between circle and ellipse fitting methods Pennsylvania, USA. They stated that woody debris, dense shrub layer, and rough topography increased tree misclassifications, which was also the case in our study. Chen et al. (2019), on the other hand, correctly detected 93.3% of the tree trunks using a personal laser scanning device combined with SLAM technology in Beijing, China. It was seen that there was a good agreement for the individual tree detection (the number of trees) results obtained by the present study with those presented by other researchers. In another study by Gollob et al. (2020), tree stems in 20-m-radius circular plots were successfully mapped with a detection rate of 96%. The researchers stated that true detection rate was affected by plot size and the lower DBH threshold. The number of true detected trees generally increased with increasing threshold. In accordance with the national guidelines (GDF 2017), we only considered trees whose DBH≥8cm in the present study. It might be expected that if a higher threshold was taken, our F-scores could increase, too. In other words, the accuracy of HMLSbased tree detection is higher in mature forest stands than those of newly-developed or thin stands. This is also true for DBH estimates derived by HMLS data. Balenović et al. (2021) attributed this phenomena to: (i) high noise, (ii) low density of point clouds, (iii) low ranging accuracy, and (iv) high beam divergence provided by HMLS devices.
As for DBH, estimations were highly correlated with ground data in both cylinder and circle/ellipse fitting methods. Pearson's r coefficients are given in Table  4. They were between 0.91 and 0.98 for cylinder fitting, and 0.83 and 0.99 for circle/ellipse fitting. Except for Plot 1, no statistically significant differences were found between cylinder and circle/ellipse fitting methods (p<0.05). In Plot 1, however, there were six forked trees (twins), which resulted in poor estimations, especially for the case of circle fitting. When necessary conditions were met, the ellipse fitting method provided more accurate estimates than circle fitting. The cylinder fitting method, on the other hand, took the cross-sections up to 2 meters, and thus, they could separately identify the forked trees. That was why its accuracy was higher than   Gollob et al. (2020), for example, tried five different methods and stated that a natural cubic spline minimized the estimation errors for ZEB-REVO HORIZON HMLS data. For terrestrial laser scanner data, the best method was reported as ellipse fitting in the same study. They stated that small trees (DBH<10 cm) were generally overestimated, while large trees (>10 cm) were underestimated by HMLS regardless of the fitting method.
In the present study, the best DBH estimations were made for Plot 3, where regular and high-quality spruce trunks existed. Their smooth barks may also help to obtain more accurate results. In contrast, the plots with mature Scots pine trunks, such as Plot 4 and Plot 5, yielded poorer estimates because of their thick and rough barks. Their bark structures frequently generated noises along the trunks on the point cloud data. Nevertheless, the results were generally in line with the related literature. Xi et al. (2016), for example, achieved an r coefficient of 0.98 between the estimated and field data, which was close to ours. Slight differences may stem from different tree species or distinct topographic settings.
The RMSEs in DBH estimations showed a relatively high variation from 0.7 to 4.9 cm with a mean of 1.8 cm in the circle/ellipse fitting method (Table 4). In the cylinder fitting method, contrastingly, they were more stable. However, its RMSEs were generally higher than those of circle/ellipse fitting. That was why the average RMSE for circle/ellipse fitting was 0.1 cm lower than the cylinder fitting method. Chen et al. (2019) and Xi et al. (2016) conducted similar LiDAR-based studies in China and Finland. They achieved the RMSEs of 1.58 and 0.90 cm, respectively. The lower values can be attributed to different instrumentations used in these studies. Xi et al. (2016), for example, used a TLS instrument (Leica HDS6100), while Chen et al. (2019) used personal laser scanning with a real-time viewer (ZEB-REVO RT). On the other hand, Balenović et al. (2021) recently reviewed the HMLS studies in the forestry literature. In this review, RMSEs were higher in the studies, including smaller trees (<10 cm) into the analyses. Accordingly, RMSEs ranged from 2.3 cm to 3.1 cm in the studies by Ryding et al. (2015), Oveland et al. (2018), and Gollob et al. (2020). In the present study, the lower DBH threshold was taken as 8 cm, as done in the Turkish forest management system.
Importantly, an assessment for bias is needed for reliable DBH estimations. Pearson's r coefficients can be very high, even for heavily biased data. Therefore, bias values for DBH estimations are also presented in Table  4. They were 0.97 cm and −0.87 cm for cylinder and circle/ellipse fitting methods, respectively. The average values close to zero showed that DBH estimations were almost unbiased. At the plot level, however, relatively biased data were found in Plot 1 for both cylinder (2.63 cm) and circle/ellipse (−9.72 cm) methods. Hyyppa et al. (2020) estimated tree DBHs using HMLS in a forest type (pine+spruce+birch) similar to ours. They reported bias values as −0.39 cm and −0.44 cm for easy (sparse) and hard (obstructed) forest plots, respectively. The differences are attributable to the HMLS devices used by two studies. We used an older model of ZEB-REVO, which had a maximum outdoor range of ~15 m in practice. Hyyppa et al. (2020), on the other hand, used the newer model -ZEB-REVO HORIZON. The maximum range for ZEB-REVO HORIZON is reported as ~100 m in its user manual (https://geoslam.com).

Performance of Machine Learning Model
An RF algorithm was used for segmentation, training, and validation steps of the entire classification process. It was found that the surface normal at the Z-axis (nz), planarity, and omni-variance were three important features for the training of the model (Fig. 7). A total of 20,000 data points were used since additional points did not significantly improve the model accuracy according to Cohen's kappa index. As seen in Fig. 8, the overall accuracy and kappa index rates were 96.06% and 88.8%, respectively. This classification accuracy is satisfactory for most forestry applications, as suggested by Kangas and Maltamo (2006) and Vatandaşlar and Zeybek (2020).

Further Improvements in Algorithmic Performance
For further improvements in algorithmic performance, the MLS method was applied to reduce point  is thought to be bark roughness and the precision of laser footprint. 85% of the data points, on the other hand, were within a limit of 0−1.0 cm distance from the surface (Fig. 9c). To overcome the point density problem, the up-sampling interval in MLS was taken as 1.0 cm for some trees. Moreover, neighborhood and standard deviation values for the SOR application were set to 3 and 1, respectively. As a result, an agreement of 99.9% was obtained between the field and enhanced data at the single-tree level. The first step in the segmentation process is to cross-section the point cloud data. It is also known as data slicing, as graphically shown in Fig. 4. In fact, there is no consensus for an optimal slice thickness in the literature. Liu et al. (2018), for example, compared ten different values (from 1 to 10 cm) for circle fitting for oak, pine, and spruce species. They found 6-cmthick slices were the best for minimizing the errors in DBH estimates. In their case, the average RMSE for tree DBHs was found to be 2.3 cm. In our case, 5-cmthick slices yielded the best result with an RMSE range from 0.7 cm to 4.9 cm, depending on the sample plot characteristics.
As for the cylinder fitting, a thickness value of 1.0 m was applied between 1.0 m and 2.0 m AGL in the present study. Similarly as in the circle fitting, there is no optimal thickness value for the cylinder fitting method in the relevant studies. In a study conducted in China, Chen et al. (2019) proposed a thickness value of 40 cm between 1.1 m and 1.5 m AGL for the noncontact cylinders. Since there is considerable uncertainty in the forestry literature, future studies should focus on parameterization of the different tree species in this manner.

Conclusions
In this study, the ZEB-REVO HMLS was used to retrieve essential inventory data from the natural forest plots under harsh topographic conditions. A new approach was developed in the R programming language for fully-automating the data extraction process, which was manually done by commercial software before. Moreover, many algorithms were tested to minimize the classification error encountered during the tree detection stage.
The results showed that both the number of trees and tree DBH parameters could be estimated at acceptable accuracy levels. The data accuracy was generally more than 90%, which was enough for the plotlevel forest inventories. It was also seen that the ellipse fitting method, in general, provided a more convenient solution for DBH estimations than the circle and cyl-inder fitting methods. However, the data accuracy significantly decreased in some »hard forest plots«. These plots had noisy elements, such as undergrowth vegetation, dense branches, forked trees, and microhabitats along tree trunks.
In conclusion, HMLS offers high operational potential for forest inventory surveying. Therefore, automated approaches -such as those presented in this study -will likely be more prevalent in the forestry sector. Thus, less field time and labor force are likely to be needed for forest professionals in the near future. However, the potential of HMLS should be further investigated in other forest types, presenting different stand structures and environmental settings. Likewise, the feasibility of HMLS should be tested for many applied fields of forest engineering, including forest road maintenance, silvicultural interventions, and logging operations.