# Clear environment
rm(list = ls(globalenv()))
# Load package
library(lidR)Read/Plot/Query
Relevant Resources
Overview
Welcome to this lidar processing tutorial using R and the lidR package! In this tutorial, you will learn how to read, visualize, and query lidar data. We’ll explore basic information about a lidar file including the header and data frame as well as visualize point clouds using different colour schemes based on attributes with plot(). We’ll use the select argument in readLAS() to load specific attributes and the filter argument to only load points of interest or apply transformations on-the-fly.
Let’s get started with processing lidar data efficiently using lidR and R! Happy learning!
Environment
We start each page by loading the necessary packages, clearing our current environment, and specifying that some warnings be turned off to make our outputs clearer. We will do this for each section in the tutorial.
Basic Usage
Load and Inspect lidar Data
Download these two .laz files to follow along with the next few pages:
Load the lidar point cloud data from a LAS file using the readLAS() function. The data is stored in the las object. We can inspect the header information and attributes of the las object by calling the object.
# Load the sample point cloud
las <- readLAS(files = "data/fm_norm.laz")
# Inspect header information and print a summary
las
#> class : LAS (v1.2 format 1)
#> memory : 14.2 Mb
#> extent : 254034.6, 254284.6, 5235452, 5235702 (xmin, xmax, ymin, ymax)
#> coord. ref. : NAD83(CSRS) / MTM zone 7
#> area : 63292 m²
#> points : 266.5 thousand points
#> type : airborne
#> density : 4.21 points/m²
#> density : 3.23 pulses/m²
# Check the file size of the loaded lidar data
format(object.size(las), "Mb")
#> [1] "20.3 Mb"Visualizing lidar Data
We can visualize the lidar data using the plot() function. We have several options to control the colours in the plot, such as selecting specific attributes from the data to be mapped.
plot() background colour
Set the background of plots to white using plot(las, bg = "white"). To keep the website code clean I’ve omitted this from examples.
plot(las)
Point Classification
Lidar systems and pre-processing steps add extra information to point clouds, beyond just the X, Y, and Z coordinates of each point.
For instance, points are typically assigned Classification values to differentiate between surfaces. This is often already done by the data provider or through user applied algorithms to classify for example noise and ground points. The meaning of these ID values depends on the LAS file version and source.
Note
LAS Standard Point classification is now standardized by the American Society for Photogrammetry and Remote Sensing (ASPRS). For a complete list of classes, refer to the LAS 1.4 specification (R15).
Our example ALS dataset, collected in 2016 over Forêt Montmorency in Quebec uses an older LAS 1.2 standard. Its classification codes are interpreted as follows:
| Class Code | Meaning | Interpretation | Number of points |
|---|---|---|---|
| 2 | Ground | The bare earth surface. | |
| 1 | Unclassified | These points have yet to be assigned a classification | |
plot(las, color = "Classification")
Lidar point clouds capture Return Intensity information; representing the strength of the returning signal. This attribute varies widely depending on sensor/acquisition characteristics and can be useful in biodiversity mapping.
plot(las, color = "Intensity")
Laser scanning systems record the scan angle of each pulse and attribute this information to the point cloud as ScanAngleRank.
plot(las, color = "ScanAngleRank")
Filtering Dataset
Often we do not need the entire set of points or attributes (columns) loaded in memory for our analysis. We can use Classification as well as other attributes to pre-filter point clouds before further processing.
There are two key ways in lidR to achieve this.
Filtering points on-the-fly (efficient)
First, we load subsets of the lidar points based on certain criteria using the filter argument in directly in readLAS().
The filter argument in readLAS() can be useful for tasks such as filtering specific classifications, isolating ground, removing noise, etc…
# Load point cloud keeping only Classification 2 (ground returns)
las <- readLAS(files = "data/fm_class.laz", filter = "-keep_class 2")Visualize the point cloud of ground returns using the plot() function.
plot(las)
-keep_class 2.
Each lidar pulse can record multiple discrete returns (points).
Here we use a filter during readLAS to subset only first returns.
# Load only the first return points
las <- readLAS(files = "data/fm_norm.laz", filter = "-keep_first")# Inspect the loaded points
las
#> class : LAS (v1.2 format 1)
#> memory : 10.9 Mb
#> extent : 254034.6, 254284.6, 5235452, 5235702 (xmin, xmax, ymin, ymax)
#> coord. ref. : NAD83(CSRS) / MTM zone 7
#> area : 63256 m²
#> points : 204.5 thousand points
#> type : airborne
#> density : 3.23 points/m²
#> density : 3.23 pulses/m²
# Check the memory size after loading only the filtered points
format(object.size(las), "Mb")
#> [1] "15.6 Mb"Visualize the point cloud of first returns using the plot() function.
plot(las)
-keep_first.
The readLAS() function also allows us to select specific attributes (columns) to be loaded into memory. This is useful to reduce memory requirements when working with large lidar datasets.
# Load only the xyz coordinates (X, Y, Z) and ignore other attributes
las <- readLAS(files = "data/fm_norm.laz", select = "xyz")# Inspect the loaded attributes
las@data
#> X Y Z
#> <num> <num> <num>
#> 1: 254034.6 5235698 3.84
#> 2: 254034.8 5235695 3.21
#> 3: 254034.6 5235695 3.11
#> 4: 254034.7 5235694 3.14
#> 5: 254034.6 5235692 4.73
#> ---
#> 266468: 254284.6 5235455 1.95
#> 266469: 254284.5 5235455 2.03
#> 266470: 254284.6 5235454 2.04
#> 266471: 254284.5 5235453 1.71
#> 266472: 254284.4 5235452 1.08
# Check the memory size (much smaller)
format(object.size(las), "Mb")
#> [1] "6.1 Mb"
See all
readLAS pre-built filters
run readLAS(filter = "-help") for a full list of these filters.
Filtering Points using filter_poi()
Alternatively LASobjects loaded into memory with readLAS() can be filtered immediately using the filter_poi() function. These filters can be custom made by combining boolean operators (==, !=, >, <, &, |, etc…) with point cloud attributes to formulate logical statements. Statements are applied to points to assign TRUE (kept) or FALSE (filtered out) values.
# Load the lidar file with all the all attributes
las <- readLAS(files = "data/fm_class.laz")# Filter points with Classification == 2
class_2 <- filter_poi(las = las, Classification == 2L)
# Combine queries to filter points with Classification 2 and ReturnNumber == 1
first_returns <- filter_poi(las = las, Classification == 2L & ReturnNumber == 1L)plot(class_2)
filter_poi().
plot(first_returns)
filter_poi().
Exercises and Questions
Using:
las <- readLAS(files = "data/fm_norm.laz")E1.
Using the plot() function plot the point cloud with a different attribute that has not been done yet
Try adding axis = TRUE, legend = TRUE to your plot argument plot(las, axis = TRUE, legend = TRUE)
E2.
Create a filtered las object of returns that have an Intensity greater that 50, and plot it.
E3.
Read in the las file with only xyz and intensity only. Hint go to the lidRbook section to find out how to do this
Conclusion
This tutorial section explored the basic usage of the lidR package in R for processing and analyzing lidar data. We covered loading lidar data, inspecting and visualizng the data, selceting specific attributes, and filtering points of interest.