Creating GIS Polygons and Points from a Point Cloud

Introduction

In this tutorial, we’ll convert a point cloud to a shapefile by generating building footprints and replacing trees with a single point. This will simplify the point cloud, representing the same area in a more understandable visual format.

A point cloud’s classification component describes which parts of the dataset are trees, buildings, and other features, so that is what we will use to create the points and polygons. Key transformers include:

• PointCloudFilter, used to extract points that have the desired classification (e.g. buildings are class 6; trees are class 5)
• PointCloudToPointCoercer, used to convert LiDAR points into GIS points
• HullReplacer, used to replace points with polygons surrounding groups of points, which is helpful in this case when a high density of points is used to represent the same feature

We will also add important data validation and refinement steps to the workspace since point clouds often contain misclassification errors.
 

Step-by-Step Instructions

In this scenario, we have a LiDAR scan of an area containing buildings and trees. We want to convert it into a GIS file that contains polygons representing buildings and points representing trees. We will create an FME workspace that transforms the LAS point cloud into a shapefile. To QA the resulting dataset, we will overlay it on an orthophoto.

Part 1: Creating Building Footprint Polygons

First, let’s create polygons to represent the buildings in the LiDAR scan.

1. Create a New Workspace
Open FME Workbench and create a blank workspace.

2. Add a LAS reader
The first step is to add the LiDAR source file into the workspace.

Click Add Reader and set the following parameters:

• Format: ASPRS Lidar Data Exchange Format (LAS)
• Dataset: <Tutorial download>/data/granville_lidar.las
• Coordinate System: UTM83-10

Click OK to add the reader to the workspace.

3. Add a PointCloudFilter
According to the ASPRS specification, buildings are classification 6, so we will use a PointCloudFilter transformer to filter out points with a classification of 6.

Add a PointCloudFilter after the reader feature type. Open the parameters and add an expression:

• Expression: @Component(classification)==6
• Output Port: buildings
• Output Unfiltered Points: No

4. Add a PointCloudToPointCoercer
Now we will convert the LiDAR points to GIS points.

Add a PointCloudToPointCoercer and connect it after the “buildings” output port. Open the parameters and set the following:
 

• Output Geometry: Single Multipoint

If you run the workspace up until here and inspect the data, you’ll see that all building points have been converted to a single multipoint feature.
 

5. Add a HullReplacer
Now we want to replace the point geometry with a polygon that surrounds the groups of points. For this, we will use a HullReplacer.

Add a HullReplacer and connect it after the “Coerced” output port. Open the parameters and set:

• Hull Type: Concave )(
• Alpha Value: 2.5

6. Add a Deaggregator
Now that we have building outlines, we would like each building to be referenced as an individual feature rather than all being referenced as a single multi-polygon feature. For this, we will use a Deaggregator, which decomposes an aggregate feature into its components.

Add a Deaggregator and connect it after the “Hull” output port. The parameters can remain set to the defaults.

The workspace should look like this:



The next steps are for refining the dataset, including:
 

• Filtering out polygons that likely came from misclassified points in the original point cloud
• Merging polygons that have holes in them due to gaps in the original point cloud

7. Add an AreaCalculator
Since LiDAR scans usually contain at least a few misclassified points, we will attempt to remove some by removing areas with a low point density (i.e. they are too small to be buildings).

Add an AreaCalculator and connect it after the Deaggreagator. Leave the default parameters.

If you connect an Inspector and run the workspace, you can inspect the results in Visual Preview. Note the values in the _area column.



Based on this, we can say that polygons must be larger than 50 units to be considered buildings.

8. Add a Tester
The second part of our QA step: now we want to filter out polygons with an _area of 50 or less.

Add a Tester after the AreaCalculator, and add the following test clause:

• Left Value: _area
• Operator: >
• Right Value: 50

In other words, the _area must be greater than 50.



9. Add a Bufferer
Next, we will fill in the gaps between building polygons. We will use a Bufferer to create solid polygons.

Add a Bufferer and connect it after the “Passed” output port. In the parameters, set:
 

• Buffer Distance: 4
• End Cap Style: Square
• Corner Style: Miter

 

10. Add Another Bufferer
Since the above Bufferer also extended the edges of the polygons, we would now like to bring them back.

Connect another Bufferer and set the following parameters:

• Buffer Distance: -2.5
• End Cap Style: Square
• Corner Style: Miter

The negative buffer will bring the edges of the polygons in again. We use -2.5 instead of -4 so that we have a safe overestimation of the boundaries, ensuring all building parts are contained within the polygons.



11. Add a Generalizer
Finally, to smooth the edges of the polygons, we will use a Generalizer.

Add a Generalizer. In the parameters, set:

• Algorithm: McMaster (Smooth)

12. Add a Shapefile Writer
The last part of the workspace is to write the building polygons to an Esri Shapefile.

Add a writer and set the following parameters:
 

• Format: Esri Shapefile
• Dataset: <Tutorial download>/data/output
• Shapefile Definition: Automatic…

13. Run the Workspace and Inspect the Results
Run the workspace and inspect the output. An orthophoto of this area is provided in the tutorial downloads if you would like to use it for comparison. You may also view the results in Visual Preview over a background map.

The image below shows some of the final buildings over top of the orthophoto. Note that although the polygons are very close to representing the actual buildings, there are still some remaining errors due to the nature of LiDAR datasets. These errors will vary depending on the scenario and can be mitigated with further iterations of buffering and testing.
 

Part 2: Creating Tree Points

In the second part of this scenario, we will create a new workspace that generates GIS points representing the trees in the LiDAR scan.

1. Create a New Workspace
Open FME Workbench and create a blank workspace.

2. Add a LAS reader
The first step is to add the LiDAR source file to the workspace.

Click Add Reader and set the following parameters:

• Format: ASPRS Lidar Data Exchange Format (LAS)
• Dataset: <Tutorial download>/data/granville_lidar.las
• Coordinate System: UTM83-10

Click OK to add the reader to the workspace.

3. Add a PointCloudFilter
According to the ASPRS specification, High Vegetation is classification 5, so we will use a PointCloudFilter transformer to filter out points with a classification of 5.

Add a PointCloudFilter after the reader feature type. Open the parameters and add an expression:

• Expression: @Component(classification)==5
• Output Port: trees
• Output Unfiltered Points: No

4. Add a PointCloudToPointCoercer
Now we will convert the LiDAR points to GIS points.

Add a PointCloudToPointCoercer and connect it after the “trees” output port. Open the parameters and set the following:
 

• Output Geometry: Single Multipoint

If you run the workspace up until here and inspect the data, you’ll see that all tree points have been converted to a single multipoint feature.



5. Add a HullReplacer
In a LiDAR scan, one tree is represented by many points, so we’ll first want to create polygons to represent the trees. A polygon might represent one tree or a cluster of trees.

Add a HullReplacer and connect it after the “Coerced” output port. Open the parameters and set:

• Hull Type: Concave )(
• Alpha Value: 2.5

6. Add a Shapefile Reader
We will use the buildings shapefile we wrote in Part 1 to ensure we will not accidentally place trees inside a building polygon.

Add a reader and set the following parameters:

• Format: Esri Shapefile
• Dataset: <Tutorial download>/data/granville_gis.shp (you can either use the building shapefile written in Part 1 or the data attached to the tutorial)

7. Add a Clipper
Add a Clipper. Connect the building Shapefile to the “Clipper” input port and the HullReplacer’s “Hull” output to “Candidate”.

The workspace should look as follows:



8. Add a Deaggregator
Now we would like each polygon to be individual features rather than a multi-polygon feature. For this, we will use a Deaggregator, which decomposes an aggregate feature into its components.

Add a Deaggregator and connect it after the “Outside” output port. The parameters can remain set to the defaults.

The output currently looks like this when inspected in Visual Preview:
 

9. Add an AreaCalculator
Some tree polygons represent a single tree, while other large polygons represent a cluster of trees. We want to create one point per tree, which means the large polygons need further processing.

Add an AreaCalculator and connect it after the Deaggregator. Leave the parameters as the defaults.

If you connect an Inspector and run the workspace, you can inspect the results in Visual Preview. Note the values in the _area column.



Based on this, we will apply further processing to trees with an area of over 100 units.

10. Add a Tester
Now we want to filter out polygons with an _area less than 100.

Add a Tester after the AreaCalculator, and add the following test clause:

• Left Value: _area
• Operator: <
• Right Value: 100

In other words, the _area must be less than 100.



11. Add a Chopper
To split the large tree polygons uniformly, we’ll use a Chopper. Add one and connect it after the “Failed” output port from the Tester. Set the following parameters:

• Mode: By Length
• Approximate Length: 15

This represents the approximate height and width of the newly chopped areas.



Up to here, the output data looks like this:
 

12. Add a CenterPointReplacer
Now we will replace the polygons with points to represent the trees.

Add a CenterPointReplacer. Connect it after the “Passed” tree polygons, as well as the “Chopped”, and “Untouched” output ports.

13. Add a Shapefile Writer
Finally, we will write the tree points to an Esri Shapefile.

Add a writer after the CenterPointReplacer and set the following parameters:

• Format: Esri Shapefile
• Dataset: <Tutorial download>/data/output
• Shapefile Definition: Automatic…

14. Run the Workspace and Inspect the Results
Run the workspace and inspect the output Shapefile. An orthophoto of this area is provided in the tutorial downloads if you would like to use it for comparison. You may also view the results in Visual Preview over a background map.

The image below shows some of the final points over top of the orthophoto. As with the buildings in Part 1, many of the output points are quite good, but inaccuracies should be expected due to the nature of LiDAR datasets, and further iterations may be required.
 

Data Attribution

The data used here originates from data made available by the City of Vancouver, British Columbia. It contains information licensed under the Open Government License - Vancouver.