Skip to content

Instantly share code, notes, and snippets.

@tim-salabim
Last active March 10, 2020 09:14
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save tim-salabim/9de17804156d1aaf745fe6df8d42a64c to your computer and use it in GitHub Desktop.
Save tim-salabim/9de17804156d1aaf745fe6df8d42a64c to your computer and use it in GitHub Desktop.
use leafgl to visualise mesh3d
library(ggplot2)
library(leafgl)
library(leaflet)
library(sf)
sfx <- sf::st_transform(sf::st_cast(dplyr::filter(silicate::inlandwaters, Province == "Tasmania"), "POLYGON")[2, ], 4326)
topo <- ceramic::cc_elevation(sfx, zoom = 9)
mesh <- anglr::as.mesh3d(anglr::copy_down(anglr::DEL(sfx, max_area = .0001),
topo))
st_as_sf.mesh3d <- function(x, ...) {
idx = if (!is.null(x$it)) x$it else x$ib
idx = rbind(idx, idx[1, ])
nc = dim(idx)[2L]
idx = as.vector(idx)
mat = cbind(x$vb[1L, idx], x$vb[2L, idx])
grp = rep(seq_len(nc), each = 4L)
lst = split.data.frame(mat, grp)
sf::st_sf(
average_z = tapply(x$vb[3L, idx], grp, mean, na.rm = TRUE)
, geometry = sf::st_sfc(
lapply(
lst
, function(i) sf::st_polygon(list(i))
)
)
, ...
)
}
mesh_sf = st_as_sf(mesh, crs = 4326)
clr = makeColorMatrix(~average_z, mesh_sf, palette = "inferno")
leaflet() %>%
addProviderTiles("CartoDB.DarkMatter") %>%
addGlPolygons(mesh_sf, color = clr, popup = TRUE) %>%
setView(lng = 146.5, lat = -42.5, zoom = 7)
@tim-salabim
Copy link
Author

Screenshot from 2020-03-07 12-23-33

@tim-salabim
Copy link
Author

The above st_as_sf.mesh3d can probably severley sped up. This is only a POC!

@dcooley
Copy link

dcooley commented Mar 8, 2020

I had a go at implementing it in Rcpp in mapdeck a while ago. You could test mapdeck:::mesh_to_sf() ?

@tim-salabim
Copy link
Author

Yip, it's much faster. I updated the original gist with a much faster version of st_as_sf.mesh3d, Rcpp still eats it for breakfast!

library(microbenchmark)

microbenchmark(
  sf = st_as_sf(mesh)
  , md = mapdeck:::mesh_to_sf(mesh, c("vb", "it"))
  , times = 10
)

Unit: milliseconds
 expr        min         lq      mean     median        uq       max neval
   sf 5517.50783 5825.09536 6232.2307 6232.57181 6524.7249 6998.7204    10
   md   82.27326   94.87767  107.7158   99.13353  105.7012  166.2793    10

@tim-salabim
Copy link
Author

@mdsumner maybe the above at_as_sf.mesh3d should be included in anglr for round-tripping purposes? It isn't as performant as @dcooley 's Rcpp version, but it isn't half bad either.

@mdsumner
Copy link

mdsumner commented Mar 9, 2020

I'll be into this in a few days going riding 👍👍

@tim-salabim
Copy link
Author

Have fun!

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment