Anika Seppelt

FORTLS is used for processing of point cloud data derived from terrestrial-based technologies such as Terrestrial Laser Scanning (TLS) or Simultaneous Localization and Mapping (SLAM). Point cloud data must be provided as .las or .laz files. The first obligatory step is the normalization of the point cloud applying the function normalize. The obtained normalized point clouds serve as input data for the tree detection functions tree.detection.single.scan, tree.detection.multi.scan and tree.detection.several.plots. The function tree.detection.single.scan detects trees from normalized TLS single-scan data and tree.detection.multi.scan from normalized TLS multi-scan (or SLAM) data. If data from more than one plot are to be analyzed automatically, the function tree.detection.several.plots should be used, which includes both the normalization and the tree detection functions and executes these functions on each input plot sequentially.


The aim of the normalization process is to obtain the coordinates relative to the plot’s center and the ground level. In this process, the functions readLAS, clip_circle, classify_ground, grid_terrain and normalize_height from the lidR package are used internally (Roussel et al., 20201). The following steps are executed:

The following figure shows the normalized point cloud which is used as example data below. The function plot from the lidR package is used to generate the figure.

multiLAS <- lidR::readLAS(paste(, "Galicia_multi_scan.laz", sep = "/"))
lidR::plot(multiLAS, color = "RGB")


The arguments of normalize

The normalize function is applied as follows:

pcd.single <- normalize(las = "Galicia_single_scan.laz", 
                 normalized = NULL,
        = 0, = 0,
                 max.dist = 10, min.height = NULL, max.height = NULL, 
                 algorithm.dtm = "knnidw", res.dtm = 0.2,
                 csf = list(cloth_resolution = 0.5),
                 RGB = TRUE,
                 scan.approach = "single",
                 id = NULL, file = "single.txt",
        =, save.result = FALSE, dir.result = NULL)

The name of the .las or .laz file containing the point cloud data is introduced in las argument and must include the .las/.laz extension. Optionally, the plot identification number (id) and the file name (file) can be defined. Both are set to NULL by default, which assigns 1 to the plot identification number and 1.txt (same name as the identification number) to the reduced point cloud saved in the working directory specified in dir.result.

The directory of the input .las/.laz files and the output file can be specified in and dir.result respectively. If not specified, the current working directory is used. The output .txt files containing the reduced point clouds will be saved if not otherwise specified in save.result (save.result = TRUE by default).

If the point cloud in the input file was already normalized, the argument normalized can be set to normalized = TRUE (default setting normalized = NULL). As a result, one part of the internal normalization process is skipped. Furthermore the scanning approach applied for data collection must be specified in scan.approach with "single" (set by default) indicating the TLS single-scan approach and "multi" indicating the TLS multi-scan and SLAM point clouds approaches.

Defining the point cloud

The planimetric coordinates \(x\) and \(y\) of the center are by default = 0 and = 0. If this does not coincide with the point cloud data, the coordinates of the plot center must be specified by and

Furthermore the size of the point cloud can be reduced by the arguments max.dist, min.height and max.height. If the maximum horizontal distance in meter to the plot center (max.dist) is set, points that are further away are discarded. Similarly, the minimum and maximum height in meters (min.height, max.height respectively) defines which points are discarded, that are those below the minimum height and those above the maximum height relative to the ground level. The default value for all three arguments is NULL. Hence, no points are discarded from the point cloud after normalization.

Adjusting the algorithms applied in normalize function

In order to generate the DTM, two different algorithms can be applied specified by algorithm.dtm. Spatial interpolation based on a k-nearest neighbor approach with inverse-distance weighting (knnidw) is selected by default. The second method is the Delaunay triangulation (tin). The resolution of the DTM (res.dtm) is set to 0.2 m by default but can be adjusted manually.

To adjust the CSF algorithm, a list with parameters (e.g. the cloth resolution which is set to 0.5 by default) can be introduced in csf.

When the point clouds are colorized, the RGB values can be used to improve the normalization and tree detection process (RGB). The colors serve to distinguish leaf from ground and stem points by the Green Leaf Algorithm (GLA, Louhaichi et al., 20014). If the GLA algorithm should be applied to remove some points from the point cloud (i.e. leaf points), it must be indicated by RGB = TRUE.

The output data frame

The normalize function generates the data frame shown below. Each row corresponds to one point of the point cloud of the input data. The columns id, file and point indicate the plot identification number, the file name and the point number respectively. The following columns contain the normalized Cartesian, cylindrical and spherical coordinates x (distance on x axis in m), y (distance on y axis in m), z (height relative to ground level in m), rho (horizontal distance in m), phi (angle in rad), r (radial distance in m) and theta (polar angle in rad). The column slope displays the slope of the terrain given in rad. If the GLA algorithm was used, the column GLA shows the results of that algorithm. Furthermore, a selection probability is assigned to each point by applying the PCP algorithm (prob) and the column shows the selected plots (indicated with 1) and discarded points (indicated with 0).

id file point x y z rho phi r theta slope R G B GLA prob prob.selec
2697 1 single.txt 1 -1.994 -8.925 0.035 9.145035 4.492581 9.145102 0.0038272 0.6136267 29696 24832 27136 -0.0673077 0.0990750 0
2699 1 single.txt 2 -1.982 -8.923 0.035 9.140473 4.493815 9.140540 0.0038291 0.6136267 33536 29184 32000 -0.0578512 0.0989762 0
2700 1 single.txt 3 -1.975 -8.918 0.036 9.134076 4.494444 9.134147 0.0039413 0.6136267 38144 33536 37376 -0.0592460 0.0988378 0
2701 1 single.txt 4 -1.977 -8.925 0.080 9.141343 4.494396 9.141693 0.0087512 0.6136267 29184 25600 29184 -0.0654206 0.0990012 0
2702 1 single.txt 5 -1.963 -8.911 0.110 9.124653 4.495563 9.125316 0.0120547 0.6136267 32000 27648 29184 -0.0505495 0.0986468 0
2703 1 single.txt 6 -1.833 -8.987 0.033 9.172026 4.511187 9.172085 0.0035979 0.6136267 34560 30976 33024 -0.0434783 0.0996605 0

The .txt file saved to the directory indicated by dir.result (if save.result = TRUE) contains a similar data frame to that shown above. However, the data frame will only include the reduced point cloud, i.e. only the selected points ( = 1). The data frame is saved without row names as .txt file by using the vroom_write function of the vroom package.

Tree detection

The tree detection functions include algorithms to detect as many trees as possible in the point clouds. Additionally, the diameter at 1.3 m above ground level (diameter at breast height, \(dbh\)) is estimated and the coordinates of the tree’s center are calculated for each detected tree. Depending on the TLS approach, different tree detection functions should be used.

Data from TLS single-scan approach

When the single-scan approach was used to collect the data, the function tree.detection.single.scan can be applied as follows:

tls.resolution = list(point.dist = 6.34, tls.dist = 10)

tree.list.single.tls <- tree.detection.single.scan(data = pcd.single,
                           dbh.min = 4, dbh.max = 200, h.min = 1.3,
                           ncr.threshold = 0.1,
                           tls.resolution = tls.resolution,
                  = NULL,
                           plot.attributes = NULL,
                           breaks = 1.3, stem.section = c(1,5),
                           save.result = FALSE, dir.result = NULL)
#> Statistical filtering of the whole point cloud
#> Detecting tree stem axes
#> Computing section: 1.3 m

The normalized and reduced point cloud, i.e. the output of the normalize function, is the input data frame for this function (data). The different arguments that can be specified are explained below.

Defining the range of diameters and heights of possible trees

With dbh.min and dbh.max, the range of possible tree diameters can be specified. Hence, only cluster of points with a bigger diameter than dbh.min and a smaller diameter than dbh.max will be considered as possible trees. Additionally, min.height defines the minimum height of a possible tree or point cluster to be considered as a tree. If not manually specified, the values are set to dbh.min = 4, dbh.max = 200 (values in cm) and h.min = 1.3 (value in m).

Resolution of the TLS

The resolution of the TLS scan (tls.resolution) can be defined either by the aperture angle or the distance between to consecutive points. The aperture angle is determined by the horizontal and vertical aperture angles (horizontal.angle and vertical.angle). When choosing the angle to define the TLS resolution, both elements must be part of the list required in tls.resolution = list(horizontal.angle, vertical.angle). The second option to determine the resolution considers the distance of two consecutive points (point.dist) at a certain distance from the TLS device (tls.dist) also given in a list as it is shown in the example above.

Including further information about the plots

In plot.attributes a data frame with attributes at plot level (e.g. strata) can be inserted for additional information. This data frame must contain a column named id coinciding with that used in the id argument of the function normalize. If there are strata, the column specifying the strata must be named stratum (numeric) for other functions (e.g., estimation.plot.size or metrics.variables). If this argument is not specified, it will be set to NULL by default and the function will not add possible plot attributes.

Algorithm to distinguish stem points and foliage points

In order to distinguish stem points from points belonging to thin branches or foliage, the local surface variation, also known as normal change rate (NCR) is calculated for each point. This is a quantitative measure of the curvature feature (Pauly et al., 20025). For each point, the NCR index is estimated in a local neighborhood with a radius of 5 cm. This radius is considered as suitable for the stem separation in forests (Ma et al. 20156; Xia et al., 20157). Higher NCR values indicate more curved surfaces e.g. branches and foliage. Therefore, a threshold (ncr.threshold) is established, which can be modified manually. By defalut it is set to 0.1 according to Zhang et al. (20198), meaning that points with a higher NCR value than that threshold are discarded.

Algorithms for identification of trees

In order to improve the detection of trees, the point cloud is reduced by removing parts of it with no trees. The argument stem.section serves to identify the part of the point cloud, i.e. a range of the coordinate \(z\), which contains less bushes, branches or other disruptive points. Hence, a range of the coordinate \(z\) and therefore a belt-like area is selected, either by defining the range manually or by an internal algorithm. This belt-like area includes predominantly the stems of the trees. Within this horizontal area, point clusters with higher density are chosen, which are supposedly the stems of the trees. Applying a circular buffer around the stems, vertical cylinders are created, which contain the stems. In the following algorithms only these vertical cylindrical parts of the point cloud are used to detect trees.

After the cylinders have been selected from the point cloud, breaks defines the height (in m) of horizontal slices on which the tree detection algorithms are applied. If not otherwise specified, slices are taken every 0.3 m starting at a height of 0.4 m until reaching the maximum height. The slices have a extension of 0.1 m (height of slice +/- 5 cm). On each slice the following algorithms are applied:

  • Removal of branches and foliage: The NCR values are calculated for each point and only those points are kept as stem points with a NCR value lower than the predefined threshold (ncr.threshold)
  • Clustering of the points: The clustering process is applied on the horizontal projection of the point’s Cartesian coordinates i.e., only their \(x\) and \(y\) coordinates are considered. The Density-Based Spatial Clustering of Applications based on the Noise (DBSCAN) method (Ester et al., 19969) implemented by the dbscan function of the dbscan package is used to perform the clustering. The radius of the epsilon neighborhood (eps) is defined as the minimum distance between two consecutive points at the furthest distance from the plot center in the respective horizontal slice
  • Removing of points belonging to remaining branches and foliage: Two characteristics are used to distinguish stem and branch clusters. First, the density of stem sections is higher than of branch and foliage sections. And second, stem points should have a predominant vertical distribution. Hence, when the point cloud is vertical projected and dissected into a grid, cells over stems have a higher point density than cells over leafs or branches. Points in cells with a lower density than the median density are removed from the point cloud
  • Calculation of the center of the potential tree section: The center is considered as the point in which the variance of the distances between all cluster points and the potential center point reaches its lowest value
  • Classification of clusters: In order to decide whether the clusters belong to a tree section, different geometric features are considered. First, the calculated center of the tree section must be located behind the stem point clusters relative to the TLS position. Second, no points should be located behind the tree surface relative to the TLS position. And third, the clusters should form an arc, which implicates that extreme points should lay further away from the TLS than points in the middle.

As explained above, these algorithms are applied on all horizontal slices (defined by breaks). Thus, tree sections are identified at different heights. Those sections that belong to the same tree are joined by applying the DBSCAN algorithm on the horizontal projection of the different sections. Thereafter, tree attributes can be estimated.

Estimation of tree attributes

The diameter of the detected trees (\(dbh\)) is obtained at 1.3 m as the double of radius. If the tree is not detected in the section at 1.3 m, the \(dbh\) is estimated by fitting a linear taper function with radius as response variable and the section heights as explanatory variables. Thus, this function allows to estimate the radius at 1.3 m and to calculate \(dbh\).

The argument defines the top stem diameter (in cm), which is used for the calculation of the commercial stem volume. If this argument is not specified, the commercial stem volume ( is not calculated.

The output data frame

id file tree x y phi phi.left phi.right h.dist dbh h v n.pts n.pts.est partial.occlusion
1 single.txt 1 -3.74725 -2.02375 3.636775 3.629975 3.676496 4.258808 35.06046 19.371 1.0023443 8235 244 4479.957 241.5228 1
1 single.txt 2 -0.61650 -5.20150 4.594416 4.603049 4.620664 5.237908 24.59479 19.725 0.5016175 1894 52 3142.674 169.4274 1
1 single.txt 3 -3.22650 -4.64500 4.105284 4.087843 4.117148 5.655646 25.72770 19.284 0.5374891 4440 116 3287.435 177.2317 0
1 single.txt 4 -6.92300 1.41650 2.939770 2.917678 2.960854 7.066428 35.82367 15.523 0.8538078 5364 268 4577.478 246.7803 0
1 single.txt 5 -2.28850 9.64400 1.803785 1.797801 1.810814 9.911810 20.74592 2.822 0.0884351 1266 180 2650.872 142.9134 0

The data frame shown above is the output of the tree.detect.single.scan function. Each row represents a detected tree (consecutively numbered in the column tree). The columns id and file display the plot identification number and the file name respectively equal to the columns in the normalize output table. The coordinates of the detected trees are given as Cartesian coordinates of the tree’s center (x and y, in m) and azimuthal angles of the center (phi in rad), the left border (phi.left in rad), the right border (phi.right in rad) and the horizontal distance from the tree’s center to the plot’s center (h.dist in m).

Furthermore, the tree attributes dbh (diameter at breast height in cm), h (total height in m), v (tree stem volume in m\(^3\)) are estimated. If was defined as argument, the volume of the stem from the ground to the height of the diameter given in is estimated (commercial stem volume, in m\(^3\)).

For each tree, the number of points of the normal section slice (1.3 m +/- 0.05 m) of the original point cloud and the reduced point cloud (n.pts and respectively) are calculated and also estimated (n.pts.est and respectively). The column partial.occlusion describes whether the the detected tree is partially occluded (1) or not (0).

The data frame is saved as .csv file without row names using the write.csv function from the utils package.

Data from TLS multi-scan approach

When multiple scans were performed in the same sampling plot (multi-scan approach) or SLAM devices were used, the function tree.detection.multi.scan can be applied as follows below. Additionally, the function normalize must be adjusted by specifying scan.approach = "multi".

pcd.multi <- normalize(las = "Galicia_multi_scan.laz", 
        = 0, = 0, RGB = TRUE,
                 scan.approach = "multi", file = "multi.txt",
        =, save.result = FALSE)

tree.list.multi.tls <- tree.detection.multi.scan(data = pcd.multi[pcd.multi$prob.selec == 1, ],
                          dbh.min = 4, dbh.max = 200, h.min = 1.3,
                          ncr.threshold = 0.1,
                          tls.precision = 0.03,
                          breaks = 1.3, stem.section = c(1,5),
                 = NULL,
                          plot.attributes = NULL,
                          save.result = FALSE, dir.result = NULL)
#> Statistical filtering of the whole point cloud
#> Detecting tree stem axes
#> Computing section: 1.3 m

The function tree.detection.multi.scan comes along with the same arguments as the function tree.detection.single.scan, which are described in “Data from TLS single-scan approach”. However, instead of specifying the resolution, the precision of the TLS (in m) can be defined in tls.precision. If not defined, the default value is 0.03 m. The procedure remains the same and the output data frame contains the all the columns explained above:

id file tree x y phi h.dist dbh h v n.pts n.pts.est partial.occlusion
1 multi.txt 1 0.31575 3.07225 1.468381 3.088433 32.07958 19.832 0.8576855 488 488 400.9930 400.9930 1
1 multi.txt 2 -3.72875 -2.02875 3.639882 4.244927 33.62085 20.656 0.9784851 75 75 420.2588 420.2588 0
1 multi.txt 3 4.69400 -0.37100 6.204312 4.708639 40.16431 22.135 1.4897246 58 58 502.0516 502.0516 0
1 multi.txt 4 -0.64700 -5.23150 4.589340 5.271357 32.43667 20.515 0.9049718 47 47 405.4565 405.4565 0
1 multi.txt 5 -3.22100 -4.65900 4.107492 5.664020 27.34435 16.210 0.5174679 234 234 341.8028 341.8028 1
1 multi.txt 6 -6.88775 1.40425 2.940472 7.029439 30.16053 21.270 0.8092723 494 494 377.0049 377.0049 1

The following figure shows the same point cloud ("Galicia_multi_scan.laz") as above. The trees that were detected by the tree.detection.multi.scan function are labeled with a red belt at 1.3 m.