Assign color to Voronoi

Question

I would like to assign each Thiessen polygon the same color as the point used to generate the tessellation. We can use the meuse dataset.

library(sf) # 'simple features' representations of spatial objects
#- read
data("meuse", package = "sp")
#- to sf
meuse.sf <- st_as_sf(meuse, coords = c("x","y"))
#- crs
st_crs(meuse.sf) <- 28992

#- Voronoi tesselation
pnts <- st_union(meuse.sf)
voronoi_grid <- st_voronoi(pnts)

#-- plot
plot(voronoi_grid, col = NA, #meuse.sf[,"zinc"],
      xlim=c(178605,181390),
      ylim=c(329714,333611))
points(st_coordinates(meuse.sf)[,1], st_coordinates(meuse.sf)[,2],
       pch=21,
       bg=sp::bpy.colors(length(meuse.sf$zinc))[rank(meuse.sf$zinc)],
       cex=0.4*meuse.sf$zinc/max(meuse.sf$zinc))
#grid()
title("Tesselated Voronoi Surface")

How do I transfer the meuse.sf['zinc'] values to the Voronoi and visualize the polygons colored like the points please?

Answer 1

You can extract polygons from the GEOMETRYCOLLECTION that's returned by st_voronoi(), convert the resulting set to an sf object and spatially join it to (a subset of) meuse.sf:

library(sf) 
data("meuse", package = "sp")
meuse.sf <- st_as_sf(meuse, coords = c("x","y"))
st_crs(meuse.sf) <- 28992

pnts <- st_union(meuse.sf)
voronoi_grid <- st_voronoi(pnts)

# check returned object
voronoi_grid
#> Geometry set for 1 feature 
#> Geometry type: GEOMETRYCOLLECTION
#> Dimension:     XY
#> Bounding box:  xmin: 174708 ymin: 325817 xmax: 185287 ymax: 337508
#> Projected CRS: Amersfoort / RD New
#> GEOMETRYCOLLECTION (POLYGON ((174708 332989.6, ...

# extract POLYGONs, turn resulting geomerty set to sf and 
# spatially join to subset of meuse.sf:
voronoi_sf <- st_sf(geometry = st_collection_extract(voronoi_grid, "POLYGON")) |>
  st_join(meuse.sf[,"zinc"])
voronoi_sf
#> Simple feature collection with 155 features and 1 field
#> Geometry type: POLYGON
#> Dimension:     XY
#> Bounding box:  xmin: 174708 ymin: 325817 xmax: 185287 ymax: 337508
#> Projected CRS: Amersfoort / RD New
#> First 10 features:
#>    zinc                       geometry
#> 1   553 POLYGON ((174708 332989.6, ...
#> 2   332 POLYGON ((178621.5 330091.1...
#> 3   907 POLYGON ((174708 335077.5, ...
#> 4   577 POLYGON ((174708 333882.1, ...
#> 5   560 POLYGON ((178796.1 330570.9...
#> 6   601 POLYGON ((177424 325817, 17...
#> 7   783 POLYGON ((174708 325817, 17...
#> 8   342 POLYGON ((178775.3 330076.5...
#> 9  1383 POLYGON ((178463.4 331081.4...
#> 10  420 POLYGON ((178770.5 330462.1...

plot(voronoi_sf[,"zinc"],
     xlim=c(178605,181390),
     ylim=c(329714,333611))

^{Created on 2023-09-05 with reprex v2.0.2}