## -----------------------------------------------------------------------------
#| label: setup
#| include: false

# Limit threads to avoid a CRAN NOTE
Sys.setenv(OMP_THREAD_LIMIT = 2)

options(scipen = 999)

## Packages needed to *build* the vignette (not to run the benchmarks)
library(ggplot2)
library(patchwork)

## ---------------------------------------------------------------------------
## Where the benchmark results live.
## Works both when the working directory is the package root (interactive use)
## and when it is vignettes/ (what happens at build time).
## ---------------------------------------------------------------------------
bench_dir <- if (dir.exists("vignettes")) "vignettes/data" else "data"
bench_csv <- file.path(bench_dir, "bench-results.csv")

if (!file.exists(bench_csv)) {
  stop(
    "Benchmark results not found at '", bench_csv, "'.\n",
    "Run the benchmark chunks of this vignette locally and execute the ",
    "'export-results' chunk at the bottom before building."
  )
}

bench_results <- utils::read.csv(bench_csv, stringsAsFactors = FALSE)

## ---------------------------------------------------------------------------
## Helpers to turn the stored numbers into prose
## ---------------------------------------------------------------------------

# Subset one operation
bench_op <- function(op) {
  out <- bench_results[bench_results$op == op, ]
  if (nrow(out) == 0) stop("No benchmark results stored for op = '", op, "'")
  out
}

# One value (median seconds or mem_alloc bytes) for one op / size / package
bench_get <- function(op, n, pkg, metric) {
  row <- bench_results[
    bench_results$op == op & bench_results$n == n & bench_results$pkg == pkg,
  ]
  if (nrow(row) != 1) {
    stop(
      "Expected exactly 1 row for op = '", op, "', n = ", n,
      ", pkg = '", pkg, "'; found ", nrow(row),
      ". The results CSV is probably stale."
    )
  }
  row[[metric]]
}

# sf / duckspatial ratio: > 1 means duckspatial wins
bench_ratio <- function(op, n, metric = c("median", "mem_alloc")) {
  metric <- match.arg(metric)
  bench_get(op, n, "sf", metric) / bench_get(op, n, "duckspatial", metric)
}

fmt_num <- function(x) {
  if (x >= 10) format(round(x)) else format(round(x, 1))
}

# "about 50 times faster than {sf}" / "about 1.2 times slower than {sf}"
bench_speed <- function(op, n) {
  r <- bench_ratio(op, n, "median")
  if (r >= 1) {
    paste0("about ", fmt_num(r), " times faster than {sf}")
  } else {
    paste0("about ", fmt_num(1 / r), " times slower than {sf}")
  }
}

# "allocated 18 times less memory" / "allocated 2 times more memory"
bench_mem <- function(op, n) {
  r <- bench_ratio(op, n, "mem_alloc")
  if (r >= 1) {
    paste0("allocated ", fmt_num(r), " times less memory")
  } else {
    paste0("allocated ", fmt_num(1 / r), " times more memory")
  }
}

## ---------------------------------------------------------------------------
## Plotting helpers (now fed by plain numeric columns, not bench objects)
## ---------------------------------------------------------------------------
ggplot_benchmark <- function(
  data,
  log = FALSE,
  show.legend = TRUE,
  ...
  ) {

    sizes <- sort(unique(data$n))
    data$n_lab <- factor(
      data$n,
      levels = sizes,
      labels = format(sizes, big.mark = ".", scientific = FALSE, trim = TRUE)
    )

    p <- ggplot(
      data = data,
      aes(
        x     = n_lab,
        y     = median,
        color = pkg,
        group = pkg
    )) +
      geom_line(
        linewidth   = 0.7,
        linetype    = "dashed",
        alpha       = 0.6,
        show.legend = show.legend
      ) +
      geom_point(
        aes(size = mem_alloc),
        alpha       = 0.9,
        show.legend = show.legend
      ) +
      scale_size_binned(
        name     = "Memory allocated",
        labels   = scales::label_bytes(),
        n.breaks = 4,
        range    = c(1, 8)
      ) +
      scale_color_brewer(palette = "Set1") +
      labs(
        color = "Package",
        x     = "Data size (rows)",
        y     = "Computation time",
        ...
      ) +
      theme_minimal(base_size = 13) +
      theme(
        panel.grid.minor = element_blank(),
        plot.title       = element_text(face = "bold")
      )

    time_lab <- scales::label_number(suffix = " s", accuracy = 0.01)

    if (log) {
      p + scale_y_log10(labels = time_lab)
    } else {
      p + scale_y_continuous(labels = time_lab)
    }
}

ggplot_assemble <- function(plot1, plot2, fun_name) {

  plot1 +
    plot2 +
    plot_annotation(
      title    = paste("Benchmark Comparison for", fun_name),
      subtitle = "Execution time vs. data size",
      theme    = theme(
        plot.title    = element_text(face = "bold", hjust = .5),
        plot.subtitle = element_text(hjust = .5)
      )
    )

}


## -----------------------------------------------------------------------------
#| eval: false
#| message: false
#| warning: false
#| code-fold: true
#| code-summary: "Set-up"

# # Load necessary packages
# library(duckspatial)
# library(bench)
# library(dplyr)
# library(sf)
# library(ggplot2)
# options(scipen = 999)
# 
# # Function to generate random points
# make_points <- function(n_points) {
#     points_df <- data.frame(
#       id = 1:n_points,
#       x = runif(n_points, min = -180, max = 180),
#       y = runif(n_points, min = -90, max = 90),
#       value = rnorm(n_points, mean = 100, sd = 15),
#       category = sample(c("A", "B", "C", "D"), n_points, replace = TRUE)
#   ) |>
#     sf::st_as_sf(coords = c("x", "y"), crs = 4326)
# }
# 
# # Turn a bench::mark() result into plain, CSV-friendly columns:
# #  - median    : seconds   (numeric)
# #  - mem_alloc : bytes     (numeric)
# bench_tidy <- function(bm, op, n) {
#   data.frame(
#     op        = op,
#     n         = n,
#     pkg       = bm$pkg,
#     median    = as.numeric(bm$median),
#     min       = as.numeric(bm$min),
#     mem_alloc = as.numeric(bm$mem_alloc),
#     n_itr     = bm$n_itr,
#     stringsAsFactors = FALSE
#   )
# }
# 
# # Generate datasets of different sizes
# withr::with_seed(27, {
#   points_sf_100k <- make_points(1e5)
#   points_sf_1mi  <- make_points(1e6)
#   points_sf_3mi  <- make_points(3e6)
# })
# 
# # Generate polygons
# # Create large polygon dataset (e.g., administrative regions, zones, etc.)
# n_polygons <- 10000
# polygons_list <- vector("list", n_polygons)
# 
# for(i in 1:n_polygons) {
#   # Random center point with buffer from edges
#   center_x <- runif(1, min = -170, max = 170)
#   center_y <- runif(1, min = -80, max = 80)
# 
#   # Create simple rectangular polygons to avoid geometry issues
#   width <- runif(1, min = 0.5, max = 3)
#   height <- runif(1, min = 0.5, max = 3)
# 
#   # Create rectangle coordinates (must be closed: first point = last point)
#   x_coords <- c(
#     center_x - width/2,
#     center_x + width/2,
#     center_x + width/2,
#     center_x - width/2,
#     center_x - width/2  # Close the polygon
#   )
# 
#   y_coords <- c(
#     center_y - height/2,
#     center_y - height/2,
#     center_y + height/2,
#     center_y + height/2,
#     center_y - height/2  # Close the polygon
#   )
# 
#   # Create polygon matrix
#   coords <- cbind(x_coords, y_coords)
# 
#   # Create polygon (wrapped in list as required by st_polygon)
#   polygons_list[[i]] <- st_polygon(list(coords))
# }
# 
# polygons_sf <- st_sf(
#   poly_id    = 1:n_polygons,
#   region     = sample(c("North", "South", "East", "West"), n_polygons, replace = TRUE),
#   population = sample(1000:1000000, n_polygons, replace = TRUE),
#   geometry   = st_sfc(polygons_list, crs = 4326)
# )


## -----------------------------------------------------------------------------
#| eval: false
#| message: false
#| code-fold: true
#| code-summary: "Benchmark code - ddbs_join"

# # Helper to run the benchmark
# run_join_benchmark <- function(points_sf) {
#   temp <- bench::mark(
#     iterations  = 3,
#     check       = FALSE,
#     duckspatial = ddbs_join(points_sf, polygons_sf, join = "within"),
#     sf          = st_join(points_sf, polygons_sf, join = st_within)
#   )
#   temp$pkg <- c("duckspatial", "sf")
#   bench_tidy(temp, op = "join", n = nrow(points_sf))
# }
# 
# # Run the benchmark
# df_bench_join <- lapply(
#   X   = list(points_sf_100k, points_sf_1mi, points_sf_3mi),
#   FUN = run_join_benchmark
# ) |>
#   dplyr::bind_rows()


## -----------------------------------------------------------------------------
#| label: fig-st-join
#| fig-cap: "Benchmark that compares spatial join for duckspatial and sf"
#| echo: false
#| warning: false
#| fig-width: 12
#| fig-height: 6

ggplot_assemble(
  plot1 = ggplot_benchmark(
    bench_op("join"), log = FALSE, show.legend = TRUE,
    subtitle = "A) Normal scale"
  ),
  plot2 = ggplot_benchmark(
    bench_op("join"), log = TRUE, show.legend = FALSE,
    subtitle = "B) Log scale"
  ),
  fun_name = "ddbs_join()"
)


## -----------------------------------------------------------------------------
#| eval: false
#| message: false
#| code-fold: true
#| code-summary: "Benchmark code - ddbs_filter"

# # Helper to run the benchmark
# run_filter_benchmark <- function(points_sf) {
#   temp <- bench::mark(
#     iterations  = 3,
#     check       = FALSE,
#     duckspatial = ddbs_filter(points_sf, polygons_sf),
#     sf          = st_filter(points_sf, polygons_sf)
#   )
#   temp$pkg <- c("duckspatial", "sf")
#   bench_tidy(temp, op = "filter", n = nrow(points_sf))
# }
# 
# # Run the benchmark
# df_bench_filter <- lapply(
#   X   = list(points_sf_100k, points_sf_1mi, points_sf_3mi),
#   FUN = run_filter_benchmark
# ) |>
#   dplyr::bind_rows()


## -----------------------------------------------------------------------------
#| label: fig-st-filter
#| fig-cap: "Benchmark that compares spatial filter for duckspatial and sf"
#| echo: false
#| warning: false
#| fig-width: 12
#| fig-height: 6

ggplot_assemble(
  plot1 = ggplot_benchmark(
    bench_op("filter"), log = FALSE, show.legend = TRUE,
    subtitle = "A) Normal scale"
  ),
  plot2 = ggplot_benchmark(
    bench_op("filter"), log = TRUE, show.legend = FALSE,
    subtitle = "B) Log scale"
  ),
  fun_name = "ddbs_filter()"
)


## -----------------------------------------------------------------------------
#| eval: false
#| message: false
#| code-fold: true
#| code-summary: "Benchmark code - ddbs_distance"

# # Helper to run the benchmark
# run_distance_benchmark <- function(n) {
# 
#   points_sf <- withr::with_seed(27, make_points(n))
# 
#   temp <- bench::mark(
#     iterations  = 1,
#     check       = FALSE,
#     duckspatial = ddbs_distance(points_sf, points_sf),
#     sf          = st_distance(points_sf, points_sf)
#   )
#   temp$pkg <- c("duckspatial", "sf")
#   bench_tidy(temp, op = "distance", n = n)
# }
# 
# df_bench_distance <- lapply(
#   X   = c(1000, 5000, 10000),
#   FUN = run_distance_benchmark
# ) |>
#   dplyr::bind_rows()


## -----------------------------------------------------------------------------
#| label: fig-st-distance
#| fig-cap: "Benchmark that compares spatial distance for duckspatial and sf"
#| echo: false
#| warning: false
#| fig-width: 12
#| fig-height: 6

ggplot_assemble(
  plot1 = ggplot_benchmark(
    bench_op("distance"), log = FALSE, show.legend = TRUE,
    subtitle = "A) Normal scale"
  ),
  plot2 = ggplot_benchmark(
    bench_op("distance"), log = TRUE, show.legend = FALSE,
    subtitle = "B) Log scale"
  ),
  fun_name = "ddbs_distance()"
)


## -----------------------------------------------------------------------------
#| eval: false
#| message: false
#| code-fold: true
#| code-summary: "Benchmark code - ddbs_union_agg"

# # Helper to run the benchmark
# run_union_benchmark <- function(points_sf) {
#   temp <- bench::mark(
#     iterations  = 3,
#     check       = FALSE,
#     duckspatial = ddbs_union_agg(points_sf, by = "category"),
#     sf          = points_sf |>
#       group_by(category) |>
#       summarise(geometry = st_union(geometry))
#   )
#   temp$pkg <- c("duckspatial", "sf")
#   bench_tidy(temp, op = "union_agg", n = nrow(points_sf))
# }
# 
# # Run the benchmark
# df_bench_union <- lapply(
#   X   = list(points_sf_100k, points_sf_1mi, points_sf_3mi),
#   FUN = run_union_benchmark
# ) |>
#   dplyr::bind_rows()


## -----------------------------------------------------------------------------
#| label: fig-st-dissolve
#| fig-cap: "Benchmark that compares dissolving geometries by group for duckspatial and sf"
#| echo: false
#| warning: false
#| fig-width: 12
#| fig-height: 6

ggplot_assemble(
  plot1 = ggplot_benchmark(
    bench_op("union_agg"), log = FALSE, show.legend = TRUE,
    subtitle = "A) Normal scale"
  ),
  plot2 = ggplot_benchmark(
    bench_op("union_agg"), log = TRUE, show.legend = FALSE,
    subtitle = "B) Log scale"
  ),
  fun_name = "ddbs_union_agg()"
)


## -----------------------------------------------------------------------------
#| eval: false
#| message: false
#| code-fold: true
#| code-summary: "Benchmark code - ddbs_intersects"

# # Helper to run the benchmark
# run_predicate_benchmark <- function(points_sf) {
#   temp <- bench::mark(
#     iterations  = 1,
#     check       = FALSE,
#     duckspatial = ddbs_intersects(points_sf, polygons_sf),
#     sf          = st_intersects(points_sf, polygons_sf)
#   )
#   temp$pkg <- c("duckspatial", "sf")
#   bench_tidy(temp, op = "intersects", n = nrow(points_sf))
# }
# 
# # Run the benchmark
# df_bench_predicate <- lapply(
#   X   = list(points_sf_100k, points_sf_1mi, points_sf_3mi),
#   FUN = run_predicate_benchmark
# ) |>
#   dplyr::bind_rows()


## -----------------------------------------------------------------------------
#| label: fig-st-intersects
#| fig-cap: "Benchmark that compares the geometry predicate ST_Intersects(x, y) for duckspatial and sf"
#| echo: false
#| warning: false
#| fig-width: 12
#| fig-height: 6

ggplot_assemble(
  plot1 = ggplot_benchmark(
    bench_op("intersects"), log = FALSE, show.legend = TRUE,
    subtitle = "A) Normal scale"
  ),
  plot2 = ggplot_benchmark(
    bench_op("intersects"), log = TRUE, show.legend = FALSE,
    subtitle = "B) Log scale"
  ),
  fun_name = "ddbs_intersects()"
)


## -----------------------------------------------------------------------------
#| label: export-results
#| eval: false
#| code-fold: true
#| code-summary: "Export the results to CSV"

# bench_export <- dplyr::bind_rows(
#   df_bench_join,
#   df_bench_filter,
#   df_bench_distance,
#   df_bench_union,
#   df_bench_predicate
# )
# 
# bench_export$date          <- as.character(Sys.Date())
# bench_export$v_duckspatial <- as.character(utils::packageVersion("duckspatial"))
# bench_export$v_duckdb      <- as.character(utils::packageVersion("duckdb"))
# bench_export$v_sf          <- as.character(utils::packageVersion("sf"))
# 
# dir.create(bench_dir, recursive = TRUE, showWarnings = FALSE)
# utils::write.csv(bench_export, bench_csv, row.names = FALSE)

