2

How can I create the smallest-possible circle which covers a polygon geometry using the sf (preferred) or sp packages in R?

One reprex using GADM data:

# Load packages
library(sf)

# Load Nigeria national boundary polygon
nigeria_raw <- readRDS(url("https://biogeo.ucdavis.edu/data/gadm3.6/Rsf/gadm36_NGA_0_sf.rds"))

# Create geometry of the smallest-possible circle
# which covers all points/entire polygon of nigeria_raw

I've not been able to find a solution online.

A solution using sf would be preferred.

Smallest-circle problem: https://en.wikipedia.org/wiki/Smallest-circle_problem

2
  • 1
    Maybe stackoverflow.com/questions/60234988/…
    – dieghernan
    Mar 1, 2021 at 20:30
  • Interesting, though I think that answer uses the centroid to find the distance to the polygon’s edge. I don’t think the smallest-possible circle is necessarily centred at the centroid.
    – Will M
    Mar 2, 2021 at 21:08

4 Answers 4

3

The library spatstats have a smallest-circle fitting function that handles sp objects if maptools is loaded.

library(sp)
library(maptools)
library(spatstat)

pol <- rbind(c(0,0), c(1,0), c(3,2), c(2,4), c(1,4), c(0,0))
pol <- Polygon(pol)
pol <- Polygons(list(pol), ID = "1")
pol <- SpatialPolygons(list(pol))
pol <- as.owin(pol)

plot(pol)
plot(boundingcircle(pol), add=TRUE, border="red", lwd=1)
plot(boundingcentre(pol), pch="+", cex=2, col="blue", add=TRUE)

Created on 2021-02-26 by the reprex package (v1.0.0)

2
  • Works great. Thanks!
    – Will M
    Feb 27, 2021 at 0:39
  • For anyone else using GADM data, I used this to convert from sf to owin (for spatstat): owin <- gadm %>% st_transform(3857) %>% as_Spatial() %>% as.owin()
    – Will M
    Feb 27, 2021 at 0:42
2

lwgeom::st_minimum_bounding_circle() yields the bounding circle in a convenient vector format (class sfc)

1

Building on Will M's answer, I show how to draw the smallest possible circle in sf around a polygon. Polygon centroid is not helpful. You need to figure which two points of the polygon are farthest from each other in order to draw the circle in the correct place and size.

# Load packages
library(sf)
library(tidyverse)
# Load Nigeria national boundary polygon
nigeria_raw <- readRDS(url("https://biogeo.ucdavis.edu/data/gadm3.6/Rsf/gadm36_NGA_0_sf.rds"))
# Format GADM data
nigeria <- nigeria_raw %>%
  transmute(country = as.character(NAME_0),
            geometry)
# Simplify country boundary for plotting
nigeria_simp <- nigeria %>%
  st_simplify(., dTolerance = 0.1)
#Extract point data from polygons (here is where our answers start to differ)
points <-st_as_text(nigeria_simp$geometry[[1]]) %>% 
  gsub("MULTIPOLYGON", "",.) %>%
  gsub("\\(|\\)","",.) %>%
  strsplit(.,",") %>%
  data.frame(.) 
colnames(points)<- c("latLong")
points$latLong<-trimws(points$latLong)
#Split point coordinates into x and y.
points<-reshape2::colsplit(points$latLong, "(\\s)", names=c("x","y"))
# Make a placeholder for calculating distance between all points..
v <- numeric(nrow(points)^2)
#loop through all points to calculate A^2 + B^2 = C^2 and store in V.
for (i in 1:nrow(points)){ 
for (z in 1:nrow(points)){
v[z+i*nrow(points)-nrow(points)] <- ((points$x[i]-points$x[z])^2+(points$y[i]-points$y[z])^2)^(1/2)
  }
}
#find the index of the two points that are farthest from one another.
#There are other ways to do this (such as in a matrix).
point_1_index<- floor(which.max(v)/nrow(points))+1
point_2_index<- which.max(v) %% nrow(points)
#Find the diameter of the circle 
diameter_of_circle<-((points$x[point_1_index] - points$x[point_2_index])^2 +(points$y[point_1_index]-points$y[point_2_index])^2)^(1/2)
#find the center of the circle. 
circle_y<-(points$y[point_1_index] + points$y[point_2_index]) /2
circle_x<-(points$x[point_1_index] + points$x[point_2_index]) /2
#make an sf point from this circle center
circle_center<-st_sfc(st_point(c(circle_x,circle_y)))
center_circle_sf= st_sf(a = 1, geom = circle_center)
st_crs(center_circle_sf) = 4326
st_geometry(center_circle_sf)
# Create the circle  using center and diameter within st_buffer
nigeria_circ <- center_circle_sf %>%
  st_buffer(., diameter_of_circle/2) 
# Plot
ggplot2::ggplot() +
  geom_sf(data = nigeria_circ, fill = NA) +
  geom_sf(data = nigeria_simp)

Smallest circle over complex polygon of the country of Nigeria

1
  • This is awesome. Thanks! :)
    – Will M
    Mar 10, 2021 at 19:42
0

This is one flawed solution I've come up with. It creates one small circle surrounding the geometry but not the smallest-possible.

This solution calculates the centroid of a polygon and then finds the largest distance from that centroid to the polygon boundary:

# Load packages
library(sf)
library(tidyverse)

# Load Nigeria national boundary polygon
nigeria_raw <- readRDS(url("https://biogeo.ucdavis.edu/data/gadm3.6/Rsf/gadm36_NGA_0_sf.rds"))

# Format GADM data
nigeria <- nigeria_raw %>%
  transmute(country = as.character(NAME_0),
            geometry
  )

# Simplify country boundary for plotting
nigeria_simp <- nigeria %>%
  st_simplify(., dTolerance = 0.1)

###### Create smallest-possible circle via centroid
# Convert GADM to line (to later find distance)
nigeria_line <- nigeria %>%
  st_cast(to = "MULTILINESTRING")

# Find centroid
## If a point is outside of country, use another function
nigeria_cent <- st_centroid(nigeria_line)

# Calculate largest distance 
nigeria_cent$distance <- st_distance(nigeria_line$geometry, nigeria_cent$geometry,
                                  by_element = TRUE,
                                  which = "Hausdorff"
)

# Create circle with buffer
nigeria_circ <- nigeria_cent %>%
  st_buffer(., .$distance)

# Plot
# Notice that the circle is small but not the smallest it could be
ggplot() +
  geom_sf(data = nigeria_circ, fill = NA) +
  geom_sf(data = nigeria_simp)

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.