The gapfraction package for R is designed for processing
airborne laser scanning (ALS) light-detection-and-ranging (LiDAR) data
of forests. A version of the package was implemented in the final
chapter of my doctoral dissertation at University of British
Columbia[1]. The package is designed to be used with LiDAR data
pre-processed with LAStools[2] or USDA Fusion[3]. The main input
to functions in the gapfraction package are LAS format height-normalized
point clouds, typically LiDAR plots corresponding to field plots. The
compressed LAZ format is not yet supported, as I rely on the rLiDAR
readLAS
function for importing data[4].
The gapfraction package implements my new fast-pit-free canopy height model (CHM) algorithm based on Khosravipour et al. (2013)[5], two new LiDAR metrics of canopy gap fraction (Po) and angular canopy closure (ACC), several recent individual tree crown (ITC) detection methods, canopy distance and direction metrics, effective leaf area index (Le) and apparent clumping index (ACI) estimation methods, as well as four mathematical fisheye (hemispherical) lens models: equi-angular, equi-distant, stereographic, and orthographic. An alphabetical list of functions in the gapfraction package is provided below.
chm
- Simple canopy height modelchm.pf
- Fast-pit-free canopy height modeldd.canopy
- Euclidean distance and direction to nearest canopy
pixel from plot centerdd.crown
- Euclidean distance and direction to nearest tree crown
from plot centersfc.aci
- Above-height cover index of fractional canopy coverfc.bl
- Beer-Lambert-Law-modified intensity-return ratio of
fractional canopy coverfc.cv
- 2-D Cartesian Voronoi fractional canopy coverfc.fci
- First-echo cover index of fractional canopy coverfc.fr
- Canopy-to-first-return ration of fractional canopy coverfc.ir
- Intesity-return ratio of fractional canopy coverfc.p
- Canopy-to-total-pixel ratio of fractional canopy coverfc.r
- Canopy-to-total-return ratio of fractional canopy coverfc.sci
- Solberg's cover index of fractional canopy covergf.hv
- Hemipsherical Voronoi canopy gap fractiongf.hv.par
- Parallel hemispherical Voronoi canopy gap fraction
with SOCKSgf.laie.aci
- Point-density-normalized canopy gap fraction,
effective LAI, and ACIitc.mw
- Variable-window individual tree crown detectionitc.mw.h
- Hierarchical variable-window individual tree crown
detectionitc.wat
- Watershed segmentation individual tree crown detectionitc.wat.h
- Hierarchical watershed segmentation individual tree
crown detectionlai.e
- Ground-to-total-return ratio with a spherical leaf angle
distributionlai.n
- Contact frequency and fractional canopy cover-based
effective LAIradial.grid.hemi
- Modified radial.grid function supporting
hemispherical lens geometriessun.path
- Modified solar position plots of Thomas SteinerYou can write math expressions: Y = Xβ + ϵ
$$\begin{align*} \sum_{k=1}^n c x_k &=cx_1+cx_2+\cdots+cx_n\\ &=c(x_1+x_2+\cdots+x_n)\\ &=c\sum_{k=1}^nx_k \end{align*}$$
For LiDAR data without ground point classifications, height-normalized
point clouds can be produced either with two LAStools command line
functions, lasground
and lasheight
, or with three functions in
USDA Fusion, GroundFilter
, GridSurfaceCreate
, and CanopyModel
.
If the ground points are already classified then you only need to use
the lasheight
function of LAStools, while the process for Fusion
still requires the same three functions. Hence, I recommend that you
use, as the LAStools ground point classification algorithm is also
superior to that of Fusion, producing more accurate
height-normalized point clouds. This is because Fusion uses the
Kraus and Pfeifer (1998) algorithm[6], while LAStools implements an
optimized version of the Axelsson (1999) algorithm[7]. For more
information, read Maguya, Junttila, and Kauranne (2014)[8]. An example
application of lasground
and lasheight
, implemented in Command
Prompt or Bash, is provided below. Very simple!
lasground -i lidar.las -o lidar_g.las
lasheight -i lidar_g.las -o lidar_n.las -replace_z
In order to run these functions, you need to istall LAStools. For
Windows, don't forget to add the LAStools bin
directory to your
system path[9]. For a single LiDAR plot, this is simple to run without
leaving your R session. You can call these functions using the system
function included in base R, as shown below.
setwd('C:/lidar')
system(lasground -i lidar.las -o lidar_g.las)
system(lasheight -i lidar_g.las -o lidar_n.las -replace_z)
To loop through LAS files stored in a folder, the syntax follows something like this:
folder <- 'C:/lidar'
files <- list.files(folder, pattern="\\.las$", full.names=TRUE)
for (i in 1:length(files)) {
file <- files[i]
basenm <- basename(file)
filenm <- strsplit(basenm,'.',fixed=TRUE)[[1]][1]
ground <- paste(folder,filenm,'_ground.las',sep='')
htnorm <- paste(folder,filenm,'_norm.las',sep='')
system(paste('lasground -i ',file,' -o ',ground, sep=''))
while (!file.exists(ground)) { Sys.sleep(1) }
system(paste('lasheight -i ',ground,' -o ',htnorm,' -replace_z', sep=''))
while (!file.exists(htnorm)) { Sys.sleep(1) }
}
What this for
loop does is read in each LAS file path, extract the
name of the file without extension, create the filenames of the ground
and height-normalized outputs, execute lasground
, wait for the output,
execute lasheight
using the ground file as the input, wait for the
output, then proceed to the next iteration. The code should be simple to
parallelize using the foreach
package, with each fork running in a new
Command Prompt or Bash window.
After loading the gapfraction package with the library
function, the
example data can be loaded by calling data(las)
. The included data
consists of X, Y, Z coordinates in UTM 11N and meters, along with
8-bit unsigned interger values for intensity and return number. The data
represents a 100-meter diameter LiDAR plot in the foothills of western
Alberta, Canada, where I conducted my PhD research. Based on previous
research, I recommend using plots minimally of this size for comparison
to ground data (e.g., hemispherical photography) to capture edge
effects.
Once the data is loaded, you can proceed to call functions from the
gapfraction
package.
Comparison of canopy height model (CHM) algorithms.
chm(las)
## class : RasterLayer
## dimensions : 100, 100, 10000 (nrow, ncol, ncell)
## resolution : 0.98802, 0.987822 (x, y)
## extent : 472568.5, 472667.4, 5882398, 5882496 (xmin, xmax, ymin, ymax)
## coord. ref. : NA
## data source : in memory
## names : layer
## values : 0, 21.83 (min, max)
chm.pf(las)
## class : RasterLayer
## dimensions : 100, 100, 10000 (nrow, ncol, ncell)
## resolution : 0.98802, 0.987822 (x, y)
## extent : 472568.5, 472667.4, 5882398, 5882496 (xmin, xmax, ymin, ymax)
## coord. ref. : NA
## data source : in memory
## names : index_1
## values : 0, 20.41133 (min, max)
Creating a pit-free CHM and performing individual tree crown (ITC) detection with the standard variable-window and watershed algorithms.
mw <- itc.mw(chm, ht2rad=function(x) 0.15746*x)
wat <- itc.wat(chm, ht2rad=function(x) 0.15746*x)
Creating a stacked pit-free CHM and performing individual tree crown (ITC) detection with the hierarchical variable-window and watershed algorithms.
mw <- itc.mw.h(chm, ht2rad=function(x) 0.15746*x, silent=TRUE)
wat <- itc.wat.h(chm, ht2rad=function(x) 0.15746*x, silent=TRUE)
Then you can use the chunk option fig.cap = "Your figure caption."
in
knitr.
The html_vignette
template includes a basic CSS theme. To override
this theme you can specify your own CSS in the document metadata as
follows:
output:
rmarkdown::html_vignette:
css: mystyles.css
Also a quote using >
:
"He who gives up [code] safety for [code] speed deserves neither." (via)
[1] Erickson, A. (2016) Forecasting Brown Bear (Ursus Arctos) Habitat Through the Integration of Remote Sensing, a Process-based Tree Establishment Model, and a Forest Landscape Model. University of British Columbia.
[2] LAStools: http://rapidlasso.com/lastools/
[3] USDA Fusion: http://forsys.cfr.washington.edu/fusion/fusionlatest.html
[4] Silva, Crookston, Hudak, and Vierling (2015). rLiDAR: LiDAR Data Processing and Visualization. R package version 0.1. https://CRAN.R-project.org/package=rLiDAR
[5] Khosravipour et al. (2013) Development of an algorithm to generate a LiDAR pit-free canopy height model. http://www.riegl.com/uploads/tx_pxpriegldownloads/khosravipour_SilviLaser2013.pdf
[6] Kraus and Pfeifer (1998) Determination of terrain models in wooded areas with airborne laser scanner data. http://www.sciencedirect.com/science/article/pii/S0924271698000094
[7] Axelsson (1999) Processing of laser scanner data—algorithms and applications. http://www.sciencedirect.com/science/article/pii/S0924271699000088
[8] Maguya, Junttila, and Kauranne (2014) http://www.mdpi.com/2072-4292/6/7/6524
[9] http://www.computerhope.com/issues/ch000549.htm
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.