hexbin_metrics cell area wrong?

[hansoleorka]

Converting the results from hexbin_metrics seems to return the wrong cell area. Example under are supposed to provide 400 m^2 cells but will return 534 m^2 using the sp::HexPoints2SpatialPolygons to convert to polygons:

library(lidR)
library(sf)
LASfile <- system.file("extdata", "Megaplot.laz", package="lidR")
lidar = readLAS(LASfile,filter = "-keep_random_fraction 0.5")
col = grDevices::colorRampPalette(c("blue", "cyan2", "yellow", "red"))

# Maximum elevation with a resolution of 8 m
#source("R/main/hexbin_metrics_my.R")
hm = hexbin_metrics(lidar, ~max(Z), 20)
hexbin::plot(hm, colramp = col, main = "Max Z")
xy <- hexbin::hcell2xy(hm)
xys <- data.frame(x=xy$x,y=xy$y)
coordinates(xys) <- ~x+y
hex <- sp::HexPoints2SpatialPolygons(xys)
hex <- sf::st_as_sf(hex)
st_area(hex)

The reason seems to be that hexbin_metrics assumes the same cell height/width in both X and Y directions. Updating code to have different height and width provides the correct area:

resy <- round(sqrt((res^2)/(6 * sqrt(3)/4)), 2)
resx <- round(((sqrt(3)/2)*resy), 2)
ext <- raster::extent(las)
xmin <- round_any(ext@xmin, resx)
xmax <- round_any(ext@xmax, resx)
ymin <- round_any(ext@ymin, resy)
ymax <- round_any(ext@ymax, resy)
if (xmax < ext@xmax) 
  xmax <- xmax + resx
if (xmin > ext@xmin) 
  xmin <- xmin - resx
if (ymax < ext@ymax) 
  ymax <- ymax + resy
if (ymin > ext@ymin) 
  ymin <- ymin - resy
dx <- (xmax - xmin)
dy <- (ymax - ymin)
xbins <- (xmax - xmin)/(2 * resx)
layout <- hexbin::hexbin(las@data$X, las@data$Y, shape = dy/dx,xbins = xbins, xbnds = c(xmin, xmax), IDs = TRUE)

[hansoleorka]

Thanks for updating! But sad to hear that you would remove it from this package - hexagons have some nice properties for some applications.

[Jean-Romain]

Yes but

1. I need an extra dependency to make and store hexagons and I want to reduce the number of dependencies
2. The support of hexagons is somewhat limited in hexbin with no native way to store stacks of hexagons, bind hexagons, merge hexagons. This is why there is no LAScatalog version for this function
3. There is not way to write hexagon in a native format AFAIK. Yes we can convert to polygons but... no... it does not scale for broad areas.

I cannot provide to hexbin_metrics the same amount of features that grid_metrics has. This is the main reason I removing it. However v3.2.0 is not released yet and if you have ideas to share or arguments against mine I'm totally open to discussion.