jlm_100_parts.R

library(tidyverse)
library(sf)
library(esri2sf)
jlm <- esri2sf("https://gisviewer.jerusalem.muni.il/arcgis/rest/services/BaseLayers/MapServer/41")
coords <- st_union(jlm) %>% 
  st_transform(2039) %>% 
  st_sample(10000) %>% 
  st_coordinates()
clus <- coords %>% 
  kmeans(100)
clus$centers %>% 
  as_tibble() %>% 
  st_as_sf(crs = 2039, coords = c("X","Y")) %>% 
  st_union() %>%
  st_voronoi(st_geometry(st_union(jlm) %>% st_transform(2039) )) %>%
  st_collection_extract() %>% 
  st_intersection(st_union(jlm) %>% st_transform(2039)) %>% 
  st_area()
  mapview::mapview()