3  Ground Classification

Classification of ground points is an important step in processing point cloud data. Distinguishing between ground and non-ground points allows creation of a continuous model of terrain elevation (see Chapter 4)). Many algorithms have been reported in the literature and lidR currently provides three of them: Progressive Morphological Filter (PMF), Cloth Simulation Function (CSF) and Multiscale Curvature Classification (MCC) usable with the function classify_ground().

3.1 Progressive Morphological Filter

The implementation of PMF algorithm in lidR is based on the method described in Zhang et al. (2003) with some technical modifications. The original method is raster-based, while lidR performs point-based morphological operations because lidR is a point cloud oriented software. The main step of the methods are summarised in the figure below:

The pmf() function requires defining the following input parameters: ws (window size or sequence of window sizes), and th (threshold size or sequence of threshold heights). More experienced users may experiment with these parameters to achieve best classification accuracy, however lidR contains util_makeZhangParam() function that includes the default parameter values described in Zhang et al. (2003).

LASfile <- system.file("extdata", "Topography.laz", package="lidR")
las <- readLAS(LASfile, select = "xyzrn")
las <- classify_ground(las, algorithm = pmf(ws = 5, th = 3))

We can now visualize the result:

plot(las, color = "Classification", size = 3, bg = "white") 

To better illustrate the classification results we can generate and plot a cross section of the point cloud (see Section 2.3.5)).

We can see that although the classification worked, there are multiple points above terrain that are classified 2 (i.e. “ground” according to ASPRS specifications). This clearly indicates that additional filtering steps are needed and that both ws and th parameters should be adjusted. Below we use multiple values for the two parameters instead of a single value in the example above.

ws <- seq(3, 12, 3)
th <- seq(0.1, 1.5, length.out = length(ws))
las <- classify_ground(las, algorithm = pmf(ws = ws, th = th))

After this adjustment the classification result changed, and points in the canopy are no longer classified as “ground”.

3.2 Cloth Simulation Function

Cloth simulation filtering (CSF) uses the Zhang et al 2016 algorithm and consists of simulating a piece of cloth draped over a reversed point cloud. In this method the point cloud is turned upside down and then a cloth is dropped on the inverted surface. Ground points are determined by analyzing the interactions between the nodes of the cloth and the inverted surface. The cloth simulation itself is based on a grid that consists of particles with mass and interconnections that together determine the three-dimensional position and shape of the cloth.

The csf() functions use the default values proposed by Zhang et al 2016 and can be used without providing any arguments.

las <- classify_ground(las, algorithm = csf())

Similar to the previous examples, classification results can be assessed using a cross section:

While the default parameters of csf() are designed to be universal and provide accurate classification results, according to the original paper, it’s apparent that the algorithm did not work properly in our example because a significant portion of points located in the ground were not classified. In such cases the algorithm parameters need to be tuned to improve the result. For this particular data set a set of parameters that resulted in an improved classification result were formulated as follows:

mycsf <- csf(sloop_smooth = TRUE, class_threshold = 1, cloth_resolution = 1, time_step = 1)
las <- classify_ground(las, mycsf)

We can also subset only the ground points to display the results in 3D

gnd <- filter_ground(las)
plot(gnd, size = 3, bg = "white") 

3.3 Multiscale Curvature Classification (MCC)

Multiscale Curvature Classification (MCC) uses the Evans and Hudak 2016 algorithm originally implemented in the mcc-lidar software.

las <- classify_ground(las, mcc(1.5,0.3))

3.4 Edge Artifacts

No matter which algorithm is used in lidR or other software, ground classification will be weaker at the edges of point clouds as the algorithm must analyze the local neighbourhood (which is missing on edges). To find ground points, an algorithm need to analyze the local neighborhood or local context that is missing at edge areas. When processing point clouds it’s important to always consider a buffer around the region of interest to avoid edge artifacts. lidR has tools to manage buffered tiles and this advanced use of the package will be covered in Chapter 14).

3.5 How to Choose a Method and Its Parameters?

Identifying the optimal algorithm parameters is not a trivial task and often requires several trial runs. lidR proposes several algorithms and may introduce even more in future versions; however, the main goal is to provide a means to compare outputs. We don’t know which one is better or which parameters best suit a given terrain. It’s likely that parameters need to be dynamically adjusted to the local context, as parameters that work well in one file may yield inadequate results in another.

Note

lidR is an research and development package, not a production package. The goal of lidR is to provide numerous and handy ways to manipulate, test, compare point cloud processing methods. If available, we recommend using classifications provided by the data provider. The classify_ground() function is useful for small to medium-sized unclassified regions of interest because it is feasible to visually assess classification results. For large acquisitions where visual assessment is no longer feasible, we do not recommend performing ground classification without first studying its accuracy.

3.6 Other Methods

Ground segmentation is not limited to the 3 methods described above. Many more have been presented and described in the literature. In Chapter 19 we will learn how to create a plugin algorithm and test it seamlessly in R with lidR using a syntax like:

las <- classify_ground(las, algorithm = my_new_method(param1 = 0.05))