Add fly_overlap, fly_select modes, and document footprint assumptions

NAMESPACE

@@ -3,6 +3,7 @@
 export(fly_coverage)
 export(fly_filter)
 export(fly_footprint)
+export(fly_overlap)
 export(fly_query_habitat)
 export(fly_query_lakes)
 export(fly_select)

R/fly_footprint.R

@@ -12,6 +12,24 @@
 #' per side. Rectangles are constructed in BC Albers (EPSG:3005) for accurate
 #' metric distances, then transformed back to the input CRS.
 #'
+#' The scale denominator is parsed from the `scale` column string (e.g.
+#' `"1:12000"` becomes `12000`).
+#'
+#' **9x9 assumption:** the default `negative_size = 9` (inches) reflects the
+#' standard 228 mm format used by BC aerial survey cameras (e.g. Wild RC-10,
+#' Zeiss RMK). The BC Air Photo Database records camera focal length per roll
+#' (Type 02 field 3.2.2) but this is not available in the simplified centroid
+#' data from the catalogue. If working with non-standard format photography,
+#' override `negative_size` accordingly.
+#'
+#' **Flat-terrain assumption:** footprints are estimated assuming flat ground
+#' beneath the aircraft. In reality terrain slope changes the actual ground
+#' coverage — downhill slopes increase the true footprint (ground falls away
+#' from the camera), while uphill slopes reduce it. In steep terrain typical
+#' of BC valleys, true footprints may differ meaningfully from these estimates.
+#' Coverage and overlap calculations downstream (e.g. [fly_coverage()],
+#' [fly_overlap()]) inherit this limitation.
+#'
 #' @examples
 #' centroids <- sf::st_read(system.file("testdata/photo_centroids.gpkg", package = "fly"))
 #' footprints <- fly_footprint(centroids)

R/fly_overlap.R

@@ -0,0 +1,81 @@
+#' Compute pairwise overlap between photo footprints
+#'
+#' For each pair of photos whose footprints intersect, computes the overlap
+#' area and the percentage of each photo's footprint that overlaps.
+#' Most useful on same-scale photos from the same flight.
+#'
+#' Overlap percentages are estimates based on flat-terrain footprints from
+#' [fly_footprint()]. See that function for details on terrain limitations.
+#'
+#' @param photos_sf An sf point object with a `scale` column.
+#' @return A tibble with columns `photo_a`, `photo_b`, `overlap_km2`,
+#'   `pct_of_a`, and `pct_of_b`. Only pairs with non-zero overlap are returned.
+#'
+#' @examples
+#' centroids <- sf::st_read(system.file("testdata/photo_centroids.gpkg", package = "fly"))
+#' aoi <- sf::st_read(system.file("testdata/aoi.gpkg", package = "fly"))
+#' photos_12k <- centroids[centroids$scale == "1:12000", ]
+#' selected <- fly_select(photos_12k, aoi, mode = "all")
+#' fly_overlap(selected)
+#'
+#' @export
+fly_overlap <- function(photos_sf) {
+  sf::sf_use_s2(FALSE)
+  on.exit(sf::sf_use_s2(TRUE))
+
+  footprints <- fly_footprint(photos_sf) |> sf::st_transform(3005)
+  n <- nrow(footprints)
+
+  if (n < 2) {
+    return(dplyr::tibble(
+      photo_a = integer(0), photo_b = integer(0),
+      overlap_km2 = numeric(0), pct_of_a = numeric(0), pct_of_b = numeric(0)
+    ))
+  }
+
+  fp_areas <- as.numeric(sf::st_area(footprints))
+  pairs <- sf::st_intersects(footprints)
+
+  ids <- if ("airp_id" %in% names(footprints)) {
+    footprints$airp_id
+  } else {
+    seq_len(n)
+  }
+
+  results <- list()
+  for (i in seq_len(n)) {
+    neighbors <- pairs[[i]]
+    neighbors <- neighbors[neighbors > i]
+    if (length(neighbors) == 0) next
+
+    for (j in neighbors) {
+      overlap_geom <- tryCatch(
+        sf::st_intersection(sf::st_geometry(footprints[i, ]),
+                            sf::st_geometry(footprints[j, ])) |>
+          sf::st_make_valid(),
+        error = function(e) NULL
+      )
+      if (is.null(overlap_geom) || length(overlap_geom) == 0) next
+
+      overlap_area <- as.numeric(sf::st_area(overlap_geom))
+      if (overlap_area <= 0) next
+
+      results <- c(results, list(dplyr::tibble(
+        photo_a = ids[i],
+        photo_b = ids[j],
+        overlap_km2 = round(overlap_area / 1e6, 3),
+        pct_of_a = round(overlap_area / fp_areas[i] * 100, 1),
+        pct_of_b = round(overlap_area / fp_areas[j] * 100, 1)
+      )))
+    }
+  }
+
+  if (length(results) == 0) {
+    return(dplyr::tibble(
+      photo_a = integer(0), photo_b = integer(0),
+      overlap_km2 = numeric(0), pct_of_a = numeric(0), pct_of_b = numeric(0)
+    ))
+  }
+
+  dplyr::bind_rows(results)
+}

R/fly_select.R

@@ -1,23 +1,60 @@
-#' Select minimum photo set to cover an AOI (greedy set cover)
+#' Select photos covering an AOI
 #'
-#' Iteratively picks the photo whose footprint covers the most uncovered
-#' area until the target coverage is reached.
+#' Two modes: `"minimal"` picks the fewest photos to reach target coverage
+#' (greedy set-cover); `"all"` returns every photo whose footprint intersects
+#' the AOI.
 #'
 #' @param photos_sf An sf point object with a `scale` column
 #'   (pre-filtered to target year/scale).
 #' @param aoi_sf An sf polygon to cover.
+#' @param mode Either `"minimal"` (fewest photos to reach target) or `"all"`
+#'   (every photo touching the AOI).
 #' @param target_coverage Stop when this fraction is reached (default 0.95).
-#' @return An sf object (subset of `photos_sf`) with added columns
-#'   `selection_order` and `cumulative_coverage_pct`.
+#'   Only used when `mode = "minimal"`.
+#' @return An sf object (subset of `photos_sf`). For `mode = "minimal"`,
+#'   includes `selection_order` and `cumulative_coverage_pct` columns.
 #'
 #' @examples
 #' centroids <- sf::st_read(system.file("testdata/photo_centroids.gpkg", package = "fly"))
 #' aoi <- sf::st_read(system.file("testdata/aoi.gpkg", package = "fly"))
-#' selected <- fly_select(centroids, aoi, target_coverage = 0.80)
-#' selected[, c("airp_id", "scale", "selection_order", "cumulative_coverage_pct")]
+#'
+#' # Fewest photos to reach 80% coverage
+#' fly_select(centroids, aoi, mode = "minimal", target_coverage = 0.80)
+#'
+#' # All photos touching the AOI
+#' fly_select(centroids, aoi, mode = "all")
 #'
 #' @export
-fly_select <- function(photos_sf, aoi_sf, target_coverage = 0.95) {
+fly_select <- function(photos_sf, aoi_sf, mode = "minimal",
+                       target_coverage = 0.95) {
+  mode <- match.arg(mode, c("minimal", "all"))
+
+  if (mode == "all") {
+    return(fly_select_all(photos_sf, aoi_sf))
+  }
+
+  fly_select_minimal(photos_sf, aoi_sf, target_coverage)
+}
+
+#' @noRd
+fly_select_all <- function(photos_sf, aoi_sf) {
+  sf::sf_use_s2(FALSE)
+  on.exit(sf::sf_use_s2(TRUE))
+
+  footprints <- fly_footprint(photos_sf)
+  aoi_union <- sf::st_transform(aoi_sf, sf::st_crs(footprints)) |>
+    sf::st_union() |>
+    sf::st_make_valid()
+
+  touches <- sf::st_intersects(footprints, aoi_union, sparse = FALSE)[, 1]
+  result <- photos_sf[touches, ]
+  message("Selected ", nrow(result), " of ", nrow(photos_sf),
+          " photos intersecting the AOI")
+  result
+}
+
+#' @noRd
+fly_select_minimal <- function(photos_sf, aoi_sf, target_coverage) {
   sf::sf_use_s2(FALSE)
   on.exit(sf::sf_use_s2(TRUE))
 
@@ -66,7 +103,7 @@ fly_select <- function(photos_sf, aoi_sf, target_coverage = 0.95) {
       error = function(e) aoi_albers
     )
 
-    pct <- as.numeric(sf::st_area(covered_in_aoi)) / aoi_area
+    pct <- sum(as.numeric(sf::st_area(covered_in_aoi))) / aoi_area
     coverage_pcts <- c(coverage_pcts, pct)
 
     if (length(selected_idx) %% 10 == 0 || pct >= target_coverage) {

tests/testthat/test-fly_overlap.R

@@ -0,0 +1,55 @@
+test_that("fly_overlap returns tibble with expected columns", {
+  centroids <- sf::st_read(testdata_path("photo_centroids.gpkg"), quiet = TRUE)
+  result <- fly_overlap(centroids)
+  expect_s3_class(result, "tbl_df")
+  expect_true(all(c("photo_a", "photo_b", "overlap_km2",
+                     "pct_of_a", "pct_of_b") %in% names(result)))
+})
+
+test_that("fly_overlap finds overlapping same-scale photos", {
+  centroids <- sf::st_read(testdata_path("photo_centroids.gpkg"), quiet = TRUE)
+  photos_12k <- centroids[centroids$scale == "1:12000", ]
+  result <- fly_overlap(photos_12k)
+  # adjacent flight-line photos should have some overlap
+
+  expect_gt(nrow(result), 0)
+  expect_true(all(result$overlap_km2 > 0))
+  expect_true(all(result$pct_of_a >= 0 & result$pct_of_a <= 100))
+  expect_true(all(result$pct_of_b >= 0 & result$pct_of_b <= 100))
+})
+
+test_that("fly_overlap uses airp_id when available", {
+  centroids <- sf::st_read(testdata_path("photo_centroids.gpkg"), quiet = TRUE)
+  photos_12k <- centroids[centroids$scale == "1:12000", ]
+  result <- fly_overlap(photos_12k)
+  if (nrow(result) > 0) {
+    expect_true(all(result$photo_a %in% photos_12k$airp_id))
+    expect_true(all(result$photo_b %in% photos_12k$airp_id))
+  }
+})
+
+test_that("fly_overlap returns empty tibble for single photo", {
+  centroids <- sf::st_read(testdata_path("photo_centroids.gpkg"), quiet = TRUE)
+  result <- fly_overlap(centroids[1, ])
+  expect_equal(nrow(result), 0)
+  expect_s3_class(result, "tbl_df")
+})
+
+test_that("fly_overlap pairs are unique (no duplicates)", {
+  centroids <- sf::st_read(testdata_path("photo_centroids.gpkg"), quiet = TRUE)
+  result <- fly_overlap(centroids)
+  if (nrow(result) > 0) {
+    pair_keys <- paste(result$photo_a, result$photo_b, sep = "-")
+    expect_equal(length(pair_keys), length(unique(pair_keys)))
+  }
+})
+
+test_that("fly_overlap larger scale has larger overlaps", {
+  centroids <- sf::st_read(testdata_path("photo_centroids.gpkg"), quiet = TRUE)
+  overlap_12k <- fly_overlap(centroids[centroids$scale == "1:12000", ])
+  overlap_31k <- fly_overlap(centroids[centroids$scale == "1:31680", ])
+  # 1:31680 footprints are ~7x larger so overlap area should be larger
+  if (nrow(overlap_12k) > 0 && nrow(overlap_31k) > 0) {
+    expect_gt(max(overlap_31k$overlap_km2), max(overlap_12k$overlap_km2))
+  }
+})