Main Content

segmentGroundSMRF

Segment ground from lidar data using a SMRF algorithm

Since R2021a

Description

groundPtsIdx = segmentGroundSMRF(ptCloud) segments the input point cloud ptCloud into ground and nonground points using a simple morphological filter (SMRF) algorithm. For more information on the SMRF algorithm, see Simple Morphological Filter.

example

groundPtsIdx = segmentGroundSMRF(ptCloud,gridResolution) additionally specifies the dimension of the grid elements.

[groundPtsIdx,nonGroundPtCloud,groundPtCloud] = segmentGroundSMRF(___) returns the ground points and nonground points as individual pointCloud objects, in addition to the ground point indices, using any of the input argument combinations from previous syntaxes.

[___] = segmentGroundSMRF(___,Name=Value) specifies options using one or more name-value arguments. For example, ElevationThreshold=0.4 sets the elevation threshold for identifying nonground points to 0.4.

example

Examples

collapse all

Segment the ground in an unorganized aerial point cloud.

Create a lasFileReader object to access the data of aerialLidarData2.las.

fileName = fullfile(toolboxdir("lidar"),"lidardata","las", ...
           "aerialLidarData2.las");
lasReader = lasFileReader(fileName);

Read the point cloud data from the LAS file using the readPointCloud function.

ptCloud = readPointCloud(lasReader);

Segment the ground data from the point cloud.

[groundPtsIdx,nonGroundPtCloud,groundPtCloud] = segmentGroundSMRF(ptCloud);

Visualize the ground and nonground points.

figure
pcshowpair(groundPtCloud,nonGroundPtCloud)

Figure contains an axes object. The axes object contains 2 objects of type scatter.

Segment and remove the ground from an organized point cloud.

Load the point cloud data into the workspace. This point cloud was captured in a highway scenario.

ld = load("drivingLidarPoints.mat");

Display the input point cloud.

pcshow(ld.ptCloud)
xlim([-40 40])
ylim([-50 50])

Figure contains an axes object. The axes object contains an object of type scatter.

Segment the ground data from the point cloud.

[~,nonGroundPtCloud,groundPtCloud] = segmentGroundSMRF( ...
  ld.ptCloud,MaxWindowRadius=5,ElevationThreshold=0.1,ElevationScale=0.25);

Visualize the nonground points.

figure
pcshow(nonGroundPtCloud)
xlim([-40 40])
ylim([-50 50])

Figure contains an axes object. The axes object contains an object of type scatter.

Input Arguments

collapse all

Point cloud data, specified as a pointCloud object.

Dimension of each grid element, specified as a positive scalar. The function samples the input point cloud into grids along the xy-direction using the grid size to create a minimum surface map. Decreasing the value of grid resolution may return nonground points as ground.

Data Types: single | double | int8 | int16 | int32 | int64 | uint8 | uint16 | uint32 | uint64

Name-Value Arguments

Specify optional pairs of arguments as Name1=Value1,...,NameN=ValueN, where Name is the argument name and Value is the corresponding value. Name-value arguments must appear after other arguments, but the order of the pairs does not matter.

Before R2021a, use commas to separate each name and value, and enclose Name in quotes.

Example: segmentGroundSMRF(ptCloud,ElevationThreshold=0.4) sets the elevation threshold for identifying nonground points to 0.4.

Maximum radius for the disc-shaped structuring element in the morphological opening operation, specified as a positive integer. You can segment large buildings as ground by specifying a smaller radius. Increasing this value can increase the computation time of the function.

Note

The default value works effectively for aerial lidar data. For better performance on terrestrial data, set MaxWindowRadius to a smaller value, such as 8.

Data Types: single | double | int8 | int16 | int32 | int64 | uint8 | uint16 | uint32 | uint64

Slope threshold for identifying grid elements as ground or nonground in a minimum elevation surface map, specified as a nonnegative scalar. The function classifies a grid element as ground if its slope is less than SlopeThreshold. Increase this value to classify steeper slopes as ground.

Data Types: single | double | int8 | int16 | int32 | int64 | uint8 | uint16 | uint32 | uint64

Elevation threshold for identifying points as ground or nonground in the estimated elevation model, specified as a nonnegative scalar. The function classifies a point as ground if the elevation difference between the point and estimated ground surface is less than ElevationThreshold. To include more points from bumpy ground, increase this value.

Note

The default value works effectively for aerial lidar data. For best results on terrestrial data, set ElevationThreshold to a smaller value, such as 0.1.

Data Types: single | double | int8 | int16 | int32 | int64 | uint8 | uint16 | uint32 | uint64

Elevation threshold scaling factor, specified as a nonnegative scalar. You can identify ground points on steep slopes by increasing this value.

Data Types: single | double | int8 | int16 | int32 | int64 | uint8 | uint16 | uint32 | uint64

Output Arguments

collapse all

Binary map of the segmented point cloud, returned as a logical matrix or a logical vector.

  • For an organized point cloud of the form M-by-N-by-3, groundPtsIdx is returned as an M-by-N logical matrix.

  • For an unorganized point cloud of the form M-by-3, groundPtsIdx is returned as an M element logical vector.

Elements of this output that correspond to ground points in the point cloud are true and nonground points are false.

Point cloud of nonground points, returned as a pointCloud object.

Point cloud of ground points, returned as a pointCloud object.

Algorithms

A simple morphological filter (SMRF) algorithm [1] segments point cloud data into ground and nonground points. The algorithm consists of three stages:

  1. Create a minimum elevation surface map from the point cloud data.

  2. Segment the surface map into ground and nonground grid elements.

  3. Segment the original point cloud data.

Create Minimum Elevation Surface Map

  1. Divide the point cloud data into grids along the xy- plane (bird's-eye view). Specify the grid size using gridResolution.

  2. Find the lowest elevation (Zmin) value for each grid element (pixel).

  3. Combine all the Zmin values into a 2-D matrix (raster image) to create a minimum elevation surface map.

Segment Surface Map

  1. Apply a morphological opening operation on the minimum surface map. This applies an erosion filter followed by a dilation filter. For more information on morphological opening, see Types of Morphological Operations.

  2. The shape of the structuring element and its window radius define the search neighbourhood for the morphological opening. Use a disc-shaped structuring element, and start with a window radius of 1 pixel. For more information, see Structuring Elements.

  3. Calculate the slope between the minimum surface and opened surface maps at each grid element. If the difference is greater than the elevation threshold, classify the pixel as nonground.

  4. Execute steps 1 through 3 iteratively. Increase the window radius by 1 pixel in each iteration until it reaches the maximum radius specified by MaxWindowRadius.

  5. The end result of the iteration process is a binary mask in which each pixel of the point cloud is classified as either ground or nonground.

Segment Original Point Cloud

  1. Apply the binary mask to the original minimum surface map to remove nonground grids.

  2. Fill the unfilled grids using image interpolation techniques to create an estimated elevation model.

  3. Calculate the elevation difference between each point in the original point cloud and the corresponding point in the estimated elevation model. If the difference is greater than ElevationThreshold, classify the pixel as nonground.

  4. Multiply the slope of the elevation model at the each point by ElevationScale, and add the result to the ElevationThreshold value to identify ground points on steep slopes.

References

[1] Pingel, Thomas J., Keith C. Clarke, and William A. McBride. “An Improved Simple Morphological Filter for the Terrain Classification of Airborne LIDAR Data.” ISPRS Journal of Photogrammetry and Remote Sensing 77 (March 2013): 21–30.https://doi.org/10.1016/j.isprsjprs.2012.12.002.

Extended Capabilities

GPU Code Generation
Generate CUDA® code for NVIDIA® GPUs using GPU Coder™.

Version History

Introduced in R2021a