## ----setup, include = TRUE, echo = FALSE--------------------------------------
knitr::opts_knit$set(
  echo = TRUE,
  message = FALSE,
  warning = FALSE,
  out.width = "900px"
)


## -----------------------------------------------------------------------------
library(chopin)
library(dplyr)
library(sf)
library(terra)
library(anticlust)
library(igraph)
library(future)
library(future.mirai)

# old par
lastpar <- par(mfrow = c(1, 1))
options(sf_use_s2 = FALSE)

## ----nc-gen-points------------------------------------------------------------
ncpoly <- system.file("shape/nc.shp", package = "sf")
ncsf <- sf::read_sf(ncpoly)
ncsf <- sf::st_transform(ncsf, "EPSG:5070")
plot(sf::st_geometry(ncsf))

# sampling clustered point
library(spatstat.random)
set.seed(202404)
ncpoints <-
  sf::st_sample(
    x = ncsf,
    type = "Thomas",
    mu = 20,
    scale = 1e4,
    kappa = 1.25e-9
  )
ncpoints <- ncpoints[seq_len(2e3L), ]

ncpoints <- sf::st_as_sf(ncpoints)
ncpoints <- sf::st_set_crs(ncpoints, "EPSG:5070")
ncpoints$pid <- sprintf("PID-%05d", seq(1, nrow(ncpoints)))
plot(sf::st_geometry(ncpoints))

# convert to terra SpatVector
ncpoints_tr <- terra::vect(ncpoints)


## ----grid-plain---------------------------------------------------------------
compregions <-
  chopin::par_pad_grid(
    ncpoints_tr,
    mode = "grid",
    nx = 8L,
    ny = 5L,
    padding = 1e4L
  )

# a list object
class(compregions)
# length of 2
names(compregions)

par(mfrow = c(2, 1))
plot(compregions$original, main = "Original grids")
plot(compregions$padded, main = "Padded grids")

## ----gen-compregions----------------------------------------------------------
compregions <-
  chopin::par_pad_grid(
    ncpoints_tr,
    mode = "grid",
    nx = 8L,
    ny = 5L,
    padding = 1e4L
  )

## ----compare-compregions, fig.width = 8, fig.height = 8-----------------------
names(compregions)

oldpar <- par()
par(mfrow = c(2, 1))
terra::plot(compregions$original, main = "Original grids")
terra::plot(compregions$padded, main = "Padded grids")

par(mfrow = c(1, 1))
terra::plot(compregions$original, main = "Original grids")
terra::plot(ncpoints_tr, add = TRUE, col = "red", cex = 0.4)



## ----split-points-grid--------------------------------------------------------
grid_quantiles <-
  chopin::par_pad_grid(
    input = ncpoints_tr,
    mode = "grid_quantile",
    quantiles = seq(0, 1, length.out = 5),
    padding = 1e4L
  )

names(grid_quantiles)
par(mfrow = c(2, 1))
terra::plot(grid_quantiles$original, main = "Original grids")
terra::plot(grid_quantiles$padded, main = "Padded grids")

par(mfrow = c(1, 1))
terra::plot(grid_quantiles$original, main = "Original grids")
terra::plot(ncpoints_tr, add = TRUE, col = "red", cex = 0.4)


## ----split-points-grid-adv----------------------------------------------------
grid_advanced1 <-
  chopin::par_pad_grid(
    input = ncpoints_tr,
    mode = "grid_advanced",
    nx = 15L,
    ny = 8L,
    padding = 1e4L,
    grid_min_features = 25L,
    merge_max = 5L
  )

par(mfrow = c(2, 1))
terra::plot(grid_advanced1$original, main = "Original grids")
terra::plot(grid_advanced1$padded, main = "Padded grids")

par(mfrow = c(1, 1))
terra::plot(grid_advanced1$original, main = "Merged grids (merge_max = 5)")
terra::plot(ncpoints_tr, add = TRUE, col = "red", cex = 0.4)


## ----split-points-count-------------------------------------------------------

ncpoints_tr$n <- 1
n_points <-
  terra::zonal(
    ncpoints_tr,
    grid_advanced1$original,
    fun = "sum"
  )[["n"]]

grid_advanced1g <- grid_advanced1$original
grid_advanced1g$n_points <- n_points

terra::plot(grid_advanced1g, "n_points", main = "Number of points in each grid")


## ----split-points-grid-adv2---------------------------------------------------
grid_advanced2 <-
  chopin::par_pad_grid(
    input = ncpoints_tr,
    mode = "grid_advanced",
    nx = 15L,
    ny = 8L,
    padding = 1e4L,
    grid_min_features = 30L,
    merge_max = 4L
  )

par(mfrow = c(2, 1))
terra::plot(grid_advanced2$original, main = "Original grids")
terra::plot(grid_advanced2$padded, main = "Padded grids")

par(mfrow = c(1, 1))
terra::plot(grid_advanced2$original, main = "Merged grids (merge_max = 8)")
terra::plot(ncpoints_tr, add = TRUE, col = "red", cex = 0.4)


## ----split-points-grid-adv3---------------------------------------------------
grid_advanced3 <-
  chopin::par_pad_grid(
    input = ncpoints_tr,
    mode = "grid_advanced",
    nx = 15L,
    ny = 8L,
    padding = 1e4L,
    grid_min_features = 25L,
    merge_max = 3L
  )

par(mfrow = c(2, 1))
terra::plot(grid_advanced3$original, main = "Original grids")
terra::plot(grid_advanced3$padded, main = "Padded grids")

par(mfrow = c(1, 1))
terra::plot(grid_advanced3$original, main = "Merged grids (merge_max = 3)")
terra::plot(ncpoints_tr, add = TRUE, col = "red", cex = 0.4)


## ----par-group-balanced-------------------------------------------------------
# ngroups should be the exact divisor of the number of points!
group_bal_grid <-
  chopin::par_pad_balanced(
    points_in = ncpoints_tr,
    ngroups = 10L,
    padding = 1e4
  )
group_bal_grid$original$CGRIDID <- as.factor(group_bal_grid$original$CGRIDID)

par(mfrow = c(2, 1))
terra::plot(group_bal_grid$original, "CGRIDID",
            legend = FALSE,
            main = "Assigned points (ngroups = 10)")
terra::plot(group_bal_grid$padded, main = "Padded grids")

# revert to the original par
par(lastpar)