## ----include = FALSE----------------------------------------------------------
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  dev = "svglite",
  fig.ext = "svg",
  fig.width = 7,
  fig.height = 5
)
# Save par settings and restore on exit (CRAN policy)
oldpar <- par(no.readonly = TRUE)
knitr::knit_hooks$set(document = function(x) { par(oldpar); x })

## ----setup--------------------------------------------------------------------
library(areaOfEffect)
library(sf)

## ----dataframe-example, eval=FALSE--------------------------------------------
# # Your occurrence data
# observations <- data.frame(
#   species = c("Oak", "Beech", "Pine", "Spruce"),
#   lon = c(14.5, 15.2, 16.8, 20.0),
#   lat = c(47.5, 48.1, 47.2, 48.5)
# )
# 
# # Classify relative to Austria
# result <- aoe(observations, "Austria")
# result$aoe_class
# #> [1] "core" "core" "halo"

## ----sf-example, eval=FALSE---------------------------------------------------
# result <- aoe(pts_sf, "AT")

## ----austria-visual, fig.cap="Austria (black) with its area of effect (dashed blue). The halo has equal area to the core."----
# Get Austria and transform to equal-area projection
austria <- get_country("AT")
austria_ea <- st_transform(austria, "ESRI:54009")

# Create a point inside Austria
dummy_pt <- st_centroid(austria_ea)

# Run aoe() to get geometries (uses buffer method by default)
result <- aoe(dummy_pt, austria_ea)
geoms <- aoe_geometry(result, "both")

# Extract geometries
austria_geom <- geoms[geoms$type == "original", ]
aoe_geom <- geoms[geoms$type == "aoe", ]

# Plot
par(mar = c(1, 1, 1, 1), bty = "n")
plot(st_geometry(aoe_geom), border = "steelblue", lty = 2, lwd = 1.5)
plot(st_geometry(austria_geom), border = "black", lwd = 2, add = TRUE)
legend("topright",
       legend = c("Austria (core)", "Area of Effect"),
       col = c("black", "steelblue"),
       lty = c(1, 2),
       lwd = c(2, 1.5),
       inset = 0.02)

## ----support------------------------------------------------------------------
support <- st_as_sf(
  data.frame(id = 1),
  geometry = st_sfc(st_polygon(list(
    cbind(c(0, 100, 100, 0, 0), c(0, 0, 100, 100, 0))
  ))),
  crs = 32631
)

## ----points-------------------------------------------------------------------
pts <- st_as_sf(
  data.frame(
    id = 1:5,
    value = c(10, 20, 15, 25, 30)
  ),
  geometry = st_sfc(
    st_point(c(50, 50)),   # center
    st_point(c(10, 10)),   # near corner
    st_point(c(95, 50)),   # near edge
    st_point(c(120, 50)),  # outside, in halo
    st_point(c(250, 250))  # far outside
  ),
  crs = 32631
)

## ----aoe----------------------------------------------------------------------
result <- aoe(pts, support)
print(result)

## ----class--------------------------------------------------------------------
result$aoe_class

## ----multiple-----------------------------------------------------------------
# Two adjacent regions
supports <- st_as_sf(
  data.frame(region = c("A", "B")),
  geometry = st_sfc(
    st_polygon(list(cbind(c(0, 50, 50, 0, 0), c(0, 0, 100, 100, 0)))),
    st_polygon(list(cbind(c(50, 100, 100, 50, 50), c(0, 0, 100, 100, 0))))
  ),
  crs = 32631
)

# Points that may fall in overlapping AoEs
pts_multi <- st_as_sf(
  data.frame(id = 1:3),
  geometry = st_sfc(
    st_point(c(25, 50)),   # inside A
    st_point(c(50, 50)),   # on boundary
    st_point(c(75, 50))    # inside B
  ),
  crs = 32631
)

result_multi <- aoe(pts_multi, supports)
print(result_multi)

## ----mask-example, fig.cap="AoE with land mask. The AoE is clipped to the land boundary."----
# Create a coastal support
support_coast <- st_as_sf(
  data.frame(id = 1),
  geometry = st_sfc(st_polygon(list(
    cbind(c(40, 80, 80, 40, 40), c(20, 20, 60, 60, 20))
  ))),
  crs = 32631
)

# Create land mask (irregular coastline)
land_mask <- st_as_sf(
  data.frame(id = 1),
  geometry = st_sfc(st_polygon(list(cbind(
    c(0, 100, 100, 70, 50, 30, 0, 0),
    c(0, 0, 50, 60, 55, 70, 60, 0)
  )))),
  crs = 32631
)

# Create some points
pts_coast <- st_as_sf(
  data.frame(id = 1:4),
  geometry = st_sfc(
    st_point(c(60, 40)),  # core
    st_point(c(50, 30)),  # core
    st_point(c(30, 40)),  # halo (on land)
    st_point(c(90, 70))   # would be halo but in sea
  ),
  crs = 32631
)

# Apply with mask
result_coast <- aoe(pts_coast, support_coast, mask = land_mask)

# Get geometries for visualization
aoe_masked <- aoe_geometry(result_coast, "aoe")
support_geom <- aoe_geometry(result_coast, "original")

par(mar = c(1, 1, 1, 1), bty = "n")
plot(st_geometry(land_mask), col = NA, border = "steelblue", lwd = 2,
     xlim = c(-10, 110), ylim = c(-10, 90))
plot(st_geometry(aoe_masked), col = rgb(0.5, 0.5, 0.5, 0.3),
     border = "steelblue", lty = 2, add = TRUE)
plot(st_geometry(support_geom), border = "black", lwd = 2, add = TRUE)

# Add points with colors
cols <- ifelse(result_coast$aoe_class == "core", "forestgreen", "darkorange")
plot(st_geometry(result_coast), col = cols, pch = 16, cex = 1.5, add = TRUE)

# Show pruned point
plot(st_geometry(pts_coast)[4], col = "gray60", pch = 4, cex = 1.2, add = TRUE)

text(85, 75, "SEA", col = "steelblue", font = 2, cex = 1.2)

legend("topleft",
       legend = c("Support", "AoE (masked)", "Coastline", "Core", "Halo", "Pruned"),
       col = c("black", "steelblue", "steelblue", "forestgreen", "darkorange", "gray60"),
       lty = c(1, 2, 1, NA, NA, NA),
       lwd = c(2, 1, 2, NA, NA, NA),
       pch = c(NA, NA, NA, 16, 16, 4),
       pt.cex = c(NA, NA, NA, 1.5, 1.5, 1.2),
       inset = 0.02)

## ----portugal-mask, fig.cap="Portugal with land-masked AoE. The halo extends into Spain but not into the Atlantic."----
# Create a point inside Portugal (approximate center of mainland)
dummy <- st_as_sf(
  data.frame(id = 1),
  geometry = st_sfc(st_point(c(-8, 39.5))),
  crs = 4326
)

# Without mask
result_no_mask <- aoe(dummy, "PT")
aoe_no_mask <- aoe_geometry(result_no_mask, "aoe")

# With mask + area=1 for equal land area
result_masked <- aoe(dummy, "PT", mask = "land", area = 1)
aoe_masked <- aoe_geometry(result_masked, "aoe")

# Get support geometry
support_geom <- aoe_geometry(result_masked, "original")

# Transform to equal area for plotting
crs_ea <- st_crs("+proj=laea +lat_0=39.5 +lon_0=-8 +datum=WGS84")
aoe_no_mask_ea <- st_transform(aoe_no_mask, crs_ea)
aoe_masked_ea <- st_transform(aoe_masked, crs_ea)
support_ea <- st_transform(support_geom, crs_ea)

# Plot - expand xlim for legend, crop bottom margin
bbox <- st_bbox(aoe_no_mask_ea)
x_range <- bbox[3] - bbox[1]
y_range <- bbox[4] - bbox[2]
par(mar = c(1, 1, 1, 1), bty = "n")
plot(st_geometry(aoe_no_mask_ea), border = "gray50", lty = 2, lwd = 1.5,
     xlim = c(bbox[1], bbox[3]),
     ylim = c(bbox[2] + y_range * 0.25, bbox[4]),
     axes = FALSE, xaxt = "n", yaxt = "n")
plot(st_geometry(aoe_masked_ea), col = rgb(0.3, 0.5, 0.7, 0.3),
     border = "steelblue", lty = 2, lwd = 1.5, add = TRUE)
plot(st_geometry(support_ea), border = "black", lwd = 2, add = TRUE)

legend("topright",
       legend = c("Portugal", "AoE (unmasked)", "AoE (land only)"),
       col = c("black", "gray50", "steelblue"),
       lty = c(1, 2, 2),
       lwd = c(2, 1.5, 1.5),
       bty = "n",
       inset = 0.05)

## ----scale--------------------------------------------------------------------
# Default: equal core/halo areas (scale = sqrt(2) - 1)
result_default <- aoe(pts, support)

# Scale = 1: larger halo (3:1 area ratio)
result_large <- aoe(pts, support, scale = 1)

## ----area-param, eval=FALSE---------------------------------------------------
# # Halo area = original area (same as scale = sqrt(2) - 1 without mask)
# result <- aoe(pts, support, area = 1)
# 
# # Halo area = half of original
# result <- aoe(pts, support, area = 0.5)

## ----area-masked, eval=FALSE--------------------------------------------------
# # Target area = 1 means halo = original, even after coastline clipping
# result <- aoe(pts, support, area = 1, mask = "land")

## ----expand-example-----------------------------------------------------------
# Create sparse data
set.seed(42)
pts_sparse <- st_as_sf(
  data.frame(id = 1:15),
  geometry = st_sfc(c(
    lapply(1:5, function(i) st_point(c(runif(1, 20, 80), runif(1, 20, 80)))),
    lapply(1:10, function(i) st_point(c(runif(1, -50, 150), runif(1, -50, 150))))
  )),
  crs = 32631
)

# Expand until at least 10 points are captured
result_expand <- aoe_expand(pts_sparse, support, min_points = 10)

## ----expand-caps, eval=FALSE--------------------------------------------------
# # Strict caps
# result <- aoe_expand(pts, support,
#                      min_points = 50,
#                      max_area = 1.5,    # halo ≤ 1.5× original
#                      max_dist = 500)    # max 500m expansion

## ----expand-info--------------------------------------------------------------
info <- attr(result_expand, "expansion_info")
info

## ----sample-example-----------------------------------------------------------
# Create imbalanced data (many core, few halo)
set.seed(42)
pts_imbal <- st_as_sf(
  data.frame(id = 1:60),
  geometry = st_sfc(c(
    lapply(1:50, function(i) st_point(c(runif(1, 10, 90), runif(1, 10, 90)))),
    lapply(1:10, function(i) st_point(c(runif(1, 110, 140), runif(1, 10, 90))))
  )),
  crs = 32631
)

result_imbal <- aoe(pts_imbal, support, scale = 1)

# Default: balance core/halo (downsamples core to match halo)
set.seed(123)
balanced <- aoe_sample(result_imbal)
table(balanced$aoe_class)

## ----sample-custom------------------------------------------------------------
# Fixed n with 70/30 split
set.seed(123)
sampled <- aoe_sample(result_imbal, n = 20, ratio = c(core = 0.7, halo = 0.3))
table(sampled$aoe_class)

## ----sample-support, eval=FALSE-----------------------------------------------
# sampled <- aoe_sample(result_multi, by = "support")

## ----border-example, fig.cap="Border classification. Points are classified by side (blue vs orange) and distance (core vs halo) from the border line."----
# Create a diagonal border line
border_line <- st_as_sf(
  data.frame(id = 1),
  geometry = st_sfc(st_linestring(matrix(
    c(0, 0,
      100, 100), ncol = 2, byrow = TRUE
  ))),
  crs = 32631
)

# Create points on both sides
set.seed(42)
pts_border <- st_as_sf(
  data.frame(id = 1:30),
  geometry = st_sfc(c(
    # Points on side 1 (above the line)
    lapply(1:15, function(i) st_point(c(runif(1, 10, 90), runif(1, 10, 90) + 20))),
    # Points on side 2 (below the line)
    lapply(1:15, function(i) st_point(c(runif(1, 10, 90), runif(1, 10, 90) - 20)))
  )),
  crs = 32631
)

# Classify by distance from border
result_border <- aoe_border(
  pts_border, border_line,
  width = 30,
  side_names = c("north", "south")
)

# Built-in plot method
plot(result_border)

## ----border-area, eval=FALSE--------------------------------------------------
# # Each side's core zone has area 5000 (in CRS units²)
# result <- aoe_border(pts, border, area = 5000)

## ----border-sample------------------------------------------------------------
# Balance by side (equal north/south)
set.seed(123)
balanced_side <- aoe_sample(result_border, ratio = c(north = 0.5, south = 0.5))
table(balanced_side$side)

# Balance by distance class
set.seed(123)
balanced_class <- aoe_sample(result_border, by = "class")
table(balanced_class$aoe_class)

## ----summary------------------------------------------------------------------
aoe_summary(result)