r-lidar / lidR

Airborne LiDAR data manipulation and visualisation for forestry application
https://CRAN.R-project.org/package=lidR
GNU General Public License v3.0
607 stars 132 forks source link

Optional extent argument in grid_canopy and grid_terrain #114

Closed jmmonnet closed 6 years ago

jmmonnet commented 6 years ago

Hi, when creating gridded data from las files, it is sometimes useful to be able to specify the extent in order to make analyses with other existing rasters easier (e.g. optical images). Could it be possible to add extent or (xmin, xmax, ymin, ymax) as optional arguments in grid_canopy and grid_terrain functions ? Thank you ! JM

Jean-Romain commented 6 years ago

Hi, thank you for the suggestion. 4 additional parameters sounds to much to me. But we can imagine to be able use a raster as input a return results only within the non empty cells of the given raster using an extent parameter. This solution sounds flexible and convenient to me and fulfill your requirements.

Jean-Romain commented 6 years ago

Warning: the following is expected to work but is currently in development (it should work but it is not documented and subject to modifications).


In grid_metric functions the parameter res can be a RasterLayer. In that case metrics will be computed only in non empty cells of the raster (all other parameters are not considered)

LASfile <- system.file("extdata", "Megaplot.laz", package="lidR")
las = readLAS(LASfile)
ext = extent(las)

# make a raster smaller than the data and with NAs
r = raster(ext-40, res = 20)
n = ncell(r)
x = runif(n)
x[sample(1:n, 10)] = NA
r[] <- x

# Display to check
plot(ext, asp = 1)
plot(r, add = T)

# Compute using a raster as reference
X = grid_metrics(las, mean(Z), res = r)
plot(X)

Also it works with a LAScatalog

ctg = catalog("~/Documents/Haliburton dataset/Landscape LiDAR/")
cores(ctg) <- 1

r = raster(ext + 500, res = 20)
n = ncell(r)
x = runif(n)
r[] <- x

plot(ctg)
plot(extent(r), add = T, col = "red")

ctg

X = grid_metrics(ctg, mean(Z), res = r)
plot(X)

To use this, install rlas v1.2.0 currently in the development branch of rlas and lidR v1.5.0in the development branch of lidR.

About grid_terrain , grid_canopy and grid_tincanopy. They are driven by other mechanisms. I will update this thread later.

Jean-Romain commented 6 years ago

Works also with grid_terrain

LASfile <- system.file("extdata", "Topography.laz", package="lidR")
las = readLAS(LASfile)

ext = extent(las)

r = raster(ext-100, res = 1)
n = ncell(r)
x = runif(n)
x[sample(1:n, 10)] = NA
r[] <- x

plot(ext)
plot(extent(r), add = T, col = "red")

dtm = grid_terrain(las, "delaunay", res = r)
ctg = catalog("~/Documents/Montmorency dataset/las/")
cores(ctg) <- 1

plot(ctg)
ext = extent(330000, 330800, 5230000, 5230800)

r = raster(ext, res = 1)
n = ncell(r)
x = runif(n)
x[sample(1:n, 10)] = NA
r[] <- x

dtm = grid_terrain(ctg, method = "delaunay", res = r)
jmmonnet commented 6 years ago

Thank you very much, I will install the devel. versions and tests it !

JM

Jean-Romain commented 6 years ago

I'm not done yet. Still missing grid_canopy which is more difficult to manage because it is written in C++

zmpeg commented 6 years ago

This looks awesome! Any chance it works for grid_tincanopy and catalog?

Jean-Romain commented 6 years ago

Sure it will. But I don't know when. It requires some free time. I guess I can do it relatively quickly for grid_tincanopy. As you can see I planned to do it for the release of 1.6.0.

Jean-Romain commented 6 years ago

6 months later (following the Debian mantra when it's ready), this suggestion is fully implemented in the draft of lidR v2.0.