rplotvoronoi

st_voronoi and basic plot with R


I would like to execute an st_voronoi() with a basic plot in R

text.txt=

name,long,lat,water_level,elevation,depth
EM_01,18.553392,-34.07027,14.4,20.358,63.0
EM_27,18.574777,-34.068709,16.196,19.966,48.0
EM_29,18.613985,-34.053271,18.766,25.477,39.0
EM_20,18.654089,-34.045177,20.102,36.502,45.0

with R code:

library(sf) # 'simple features' representations of spatial objects

file = 'text.txt'
#-- import
cfaq <- read.csv(file, header = 1, sep = ',', dec = '.') 
#- set as a spatial feature with xy coords an existing projection
cfaq.sf <- st_as_sf(cfaq, coords=c("long", "lat"), crs = 4326) #wgs83
#- transform to local crs
cfaq.sf <- st_transform(cfaq.sf, crs = 32734) #utm 34s
# Voronoi tesselation
voronoi_grid <- st_voronoi(cfaq.sf)

# basic plot with points on top of the Voronoi 
plot(voronoi_grid, col = "lightblue", border = "black", lwd = 1.5)
plot(cfaq.sf, col = "red", pch = 16, cex = 2, add = TRUE)
#text(cfaq.sf$long, cfaq.sf$lat, labels = cfaq.sf$name, pos = 3, cex = 0.8)
# Add a legend
#legend("bottomright", legend = "Points", col = "red", pch = 16, bty = "n", pt.cex = 1.5)
# Add a title
title(main = "Voronoi Diagram with Points")

result:

Error in if (xsize * ysize * n > prod(total_size)) {: missing value where TRUE/FALSE needed
Traceback:

1. plot(voronoi_grid, col = "lightblue", border = "black", lwd = 1.5)
2. plot(voronoi_grid, col = "lightblue", border = "black", lwd = 1.5)
3. plot.sf(voronoi_grid, col = "lightblue", border = "black", lwd = 1.5)
4. .get_layout(st_bbox(x), min(max.plot, length(cols)), par("din"), 
 .     key.pos, key.width)
5. vapply(1:n, function(x) size(x, n, asp), 0)
6. FUN(X[[i]], ...)
7. size(x, n, asp)

how do I execute a st_voronoi() and then plot the result over the original points with the basic R plot() function please?


Solution

  • Ref. examples guide you to call st_voronoi() on a MULTIPOINT and extract polygons with st_collection_extract(), something like this:

    library(sf)
    #> Linking to GEOS 3.9.3, GDAL 3.5.2, PROJ 8.2.1; sf_use_s2() is TRUE
    cfaq.sf <- 
      read.csv(text = "name,long,lat,water_level,elevation,depth
                       EM_01,18.553392,-34.07027,14.4,20.358,63.0
                       EM_27,18.574777,-34.068709,16.196,19.966,48.0
                       EM_29,18.613985,-34.053271,18.766,25.477,39.0
                       EM_20,18.654089,-34.045177,20.102,36.502,45.0") |>
      st_as_sf(coords=c("long", "lat"), crs = 4326) |>
      st_transform(cfaq.sf, crs = 32734)
    
    voronoi_grid <- 
      st_union(cfaq.sf) |> 
      st_voronoi() |> 
      st_collection_extract()
    
    plot(voronoi_grid, col = "lightblue", border = "black", lwd = 1.5)
    plot(cfaq.sf, col = "red", pch = 16, cex = 2, add = TRUE)
    #> Warning in plot.sf(cfaq.sf, col = "red", pch = 16, cex = 2, add = TRUE):
    #> ignoring all but the first attribute
    title(main = "Voronoi Diagram with Points")
    

    Created on 2023-08-02 with reprex v2.0.2