In R you can do this with the dev-version of quadmesh
, and the sf
package.
library(raster)
## create a raster
r <- raster(volcano)
## convert to indexed quad form, the rgl mesh3d format
qm <- quadmesh::quadmesh(r)
## split quad index into triangle index
tri <- quadmesh::triangulate_quads(qm$ib)
## convert to simple features
ps <- lapply(split(rbind(tri, tri[1, ]), rep(seq_len(ncol(tri)), each = 4)),
function(idx) sf::st_polygon(list(t(qm$vb[1:3, idx]))))
p_xyz <- sf::st_sf(geometry = sf::st_sfc(ps, crs = projection(r)))
p_xyz
Simple feature collection with 10614 features and 0 fields
geometry type: POLYGON
dimension: XYZ
bbox: xmin: 1.490116e-08 ymin: 1.490116e-08 xmax: 1 ymax: 1
epsg (SRID): NA
proj4string: NA
First 10 features:
geometry
1 POLYGON Z ((1.490116e-08 1 ...
2 POLYGON Z ((0.01639344 1 99...
3 POLYGON Z ((0.01639344 1 99...
4 POLYGON Z ((0.03278689 1 10...
5 POLYGON Z ((0.03278689 1 10...
6 POLYGON Z ((0.04918033 1 10...
7 POLYGON Z ((0.04918033 1 10...
8 POLYGON Z ((0.06557377 1 10...
9 POLYGON Z ((0.06557377 1 10...
10 POLYGON Z ((0.08196721 1 10...
Then convert to text form with
sf::st_as_text(sf::st_geometry(p_xyz))
For triangulate_quads
the development version is required, e.g.
devtools::install_github("hypertidy/quadmesh")
https://github.com/hypertidy/quadmesh