diff --git a/.gitignore b/.gitignore index 8798d74..48ffc2c 100644 --- a/.gitignore +++ b/.gitignore @@ -7,3 +7,6 @@ docs /doc/ /Meta/ Rplots.pdf + +# Calibration report cache for data-raw/make_camera_formats.R (~25 MB, refetchable) +data-raw/.cache/ diff --git a/CLAUDE.md b/CLAUDE.md index 31ec7a3..bcff48f 100644 --- a/CLAUDE.md +++ b/CLAUDE.md @@ -14,8 +14,11 @@ metadata from the BC Data Catalogue, and georeference the images onto their esti - One exported function per file: `R/fly_footprint.R` → `tests/testthat/test-fly_footprint.R` - `inst/testdata/` — Upper Bulkley River floodplain near Houston, BC (20 photos, dual scale). All 1968 film, -`Film - BW`, focal 153 — there is no digital frame in it, and `data-raw/make_testdata.R` sources a film-only AOI, -so digital coverage cannot come from there +`Film - BW`, focal 153 — so the *sampled* photos cannot exercise the digital path. The **AOI is not +film-only**, though: it holds 181 `Digital - Colour` frames from two cameras, 24 of which now ship as +`inst/testdata/photo_centroids_digital.gpkg` (fly#32). An earlier note here said digital coverage could not +come from there; that came from a `BBOX(SHAPE, ...)` CQL query which returns 0 features for this AOI even for +a positive control — a broken probe, not an absence. Query bboxes through `bcdata::filter(BBOX(...))` - `inst/testdata/dem.tif` — MRDEM-30 clip (NRCan 30 m bare-earth DTM), buffered 5.4 km past the centroids. It is one CRS at one resolution, so it **cannot** exercise the reprojection, coarse-grid, truncating-extent or wide-spread branches of the terrain code — see `inst/notes/terrain-correction.md` diff --git a/DESCRIPTION b/DESCRIPTION index 831cf9d..ba7934b 100644 --- a/DESCRIPTION +++ b/DESCRIPTION @@ -1,7 +1,7 @@ Package: fly Title: Historic Airphoto Footprints, Selection and Georeferencing for British Columbia -Version: 0.5.1 +Version: 0.6.0 Date: 2026-08-29 Authors@R: c( person("Allan", "Irvine", , "al@newgraphenvironment.com", role = c("aut", "cre"), diff --git a/NEWS.md b/NEWS.md index 66cde2a..5107200 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,5 +1,18 @@ # fly (development version) +## 0.6.0 (2026-08-30) + +- `fly_footprint()` now sizes digital frames, closing the gap #30 made honest but left open ([#32](https://github.com/NewGraphEnvironment/fly/issues/32)). Province-wide that is 223,667 of 1,670,471 frames — the package was quietly film-only for anything after ~2010 +- Sensor dimensions are read from the camera calibration reports the catalogue itself links to through `camera_calibration_url`, and shipped as `inst/extdata/camera_formats.csv` (built by `data-raw/make_camera_formats.R`). 14 calibrations covering 169,688 frames, plus focal-length fallback rows for frames carrying no calibration +- **The catalogue's `SCALE` is not the true image scale for a digital frame, and is no longer used for one.** Measured against terrain on 40 UltraCam Eagle frames it gives 34% of true width: it is a derived nominal figure, implying a pixel pitch of ~12.5 um for every camera regardless of model against real pitches of 3.9-12 um. A digital frame is sized as `pixel count x ground_sample_distance` instead, which needs neither `scale` nor a DEM. `ground_sample_distance` is centimetres +- **Footprints are no longer always square.** Digital sensors run from 1.10:1 (Leica DMC II) to 1.80:1 (Intergraph DMC), so a square footprint was up to 76% too deep. Non-square footprints are rotated onto the flight line via `fly_bearing()`. Film stays square and its output is unchanged +- New `width_source` column names the calibration file or fallback rule behind every digital footprint, and `footprint_terrain` gains `"gsd_scaled"` — `nominal_scale` is documented as "sized from the reported scale", which is the one thing this route never does +- Frames whose calibration could not be corroborated are refused rather than inferred, listed with the reason in `inst/extdata/camera_formats_excluded.csv`. Two are medium-format bodies about half the width of everything else in the record, so inferring one from focal length would have been ~1.95x too wide +- `fly_georef()` excludes rotated footprints with a warning: its corner mapping applies its own bearing rotation, calibrated for axis-aligned squares, and would count the rotation twice +- The shipped numbers are parsed from the reports, never typed, and gated on four checks before the table is written — `px x pitch` against the stated image size, report focal against the catalogue's, plausibility bounds, and an implied ground elevation that must be a real BC elevation. The last two caught a camera whose catalogue metadata contradicts its own report, which is withheld +- `format_size` is unchanged and still takes precedence, so a caller who knows their camera can override the shipped table +- New `inst/testdata/photo_centroids_digital.gpkg`: 24 real digital frames over the AOI that already ships, from two cameras 0.46 apart in aspect ratio + ## 0.5.1 (2026-08-29) - Fix `fly_footprint()` silently dropping `footprint_basis`, `footprint_terrain`, `height_agl` and `dem_coverage` whenever its input carried the `tbl_df` class ([#35](https://github.com/NewGraphEnvironment/fly/issues/35)). `bcdata::collect()` returns exactly that class, so every caller querying `WHSE_IMAGERY_AND_BASE_MAPS.AIMG_PHOTO_CENTROIDS_SP` — the documented source for this package — lost the whole reporting surface 0.4.0 and 0.5.0 added, and the documented "filter on `footprint_basis`" and "filter on `dem_coverage`" workflows were unreachable from it diff --git a/R/fly_bearing.R b/R/fly_bearing.R index 7995224..c207605 100644 --- a/R/fly_bearing.R +++ b/R/fly_bearing.R @@ -41,6 +41,10 @@ fly_bearing <- function(photos_sf) { ord <- order(photos_sf$film_roll, photos_sf$frame_number) + # `isTRUE()` on the roll comparisons below, rather than a bare `==`: an NA in + # `film_roll` makes the comparison NA, which aborts an `if` with "missing value where + # TRUE/FALSE needed". A frame with no roll simply has no neighbour to take a bearing + # from, which is what NA already means here. bearing <- rep(NA_real_, nrow(photos_sf)) rolls <- photos_sf$film_roll[ord] @@ -48,12 +52,12 @@ fly_bearing <- function(photos_sf) { y <- coords[ord, 2] for (i in seq_along(ord)) { - if (i < length(ord) && rolls[i] == rolls[i + 1]) { + if (i < length(ord) && isTRUE(rolls[i] == rolls[i + 1])) { # Forward bearing to next frame on same roll dx <- x[i + 1] - x[i] dy <- y[i + 1] - y[i] bearing[ord[i]] <- (atan2(dx, dy) * 180 / pi) %% 360 - } else if (i > 1 && rolls[i] == rolls[i - 1]) { + } else if (i > 1 && isTRUE(rolls[i] == rolls[i - 1])) { # Last frame on roll: use bearing from previous dx <- x[i] - x[i - 1] dy <- y[i] - y[i - 1] diff --git a/R/fly_camera_format.R b/R/fly_camera_format.R new file mode 100644 index 0000000..fb2bb86 --- /dev/null +++ b/R/fly_camera_format.R @@ -0,0 +1,150 @@ +# Recording-format dimensions for the digital cameras in the BC air photo catalogue. +# +# `AIMG_PHOTO_CENTROIDS_SP` carries no sensor size, which is why #30 refused to size +# digital frames rather than invent one. The number is recoverable from the calibration +# report each frame links to through `camera_calibration_url`, and +# `data-raw/make_camera_formats.R` parses it out of those reports into +# `inst/extdata/camera_formats.csv`. See fly#32. +# +# Two key types share one table so one lookup reads both: +# +# `calib_file` — keyed on the calibration file the frame names. Exact: the report +# gives the array size and pixel pitch, and the millimetres are +# checked against them. +# `focal_length` — keyed on the catalogue's `focal_length`, for the ~20% of digital +# frames carrying no calibration URL. Inferred, and carrying +# `width_spread_pct` so the room for error travels with the number. + +fly_camera_cache <- new.env(parent = emptyenv()) + +fly_camera_read <- function(file) { + if (is.null(fly_camera_cache[[file]])) { + path <- system.file("extdata", file, package = "fly") + if (!nzchar(path)) { + stop("`", file, "` is missing from the installed package.", call. = FALSE) + } + # `key` must stay character: the fallback keys are focal lengths, and read.csv would + # type them numeric, so `"80"` and `80` would stop matching between the two halves + # of the same table. + fly_camera_cache[[file]] <- utils::read.csv( + path, stringsAsFactors = FALSE, colClasses = c(key = "character") + ) + } + fly_camera_cache[[file]] +} + +# The shipped format table. +fly_camera_table <- function() fly_camera_read("camera_formats.csv") + +# Calibrations the catalogue offers that are deliberately not shipped, each with the +# reason. Kept beside the table rather than dropped, so "we have not looked at this" +# and "we looked and could not use it" stay distinguishable. +fly_camera_excluded <- function() fly_camera_read("camera_formats_excluded.csv") + + +# `GROUND_SAMPLE_DISTANCE` is recorded in CENTIMETRES. +# +# Worth a named function rather than a bare `/ 100`, because getting it wrong is a +# factor of 100 in every digital footprint and the field name says nothing about units. +# Confirmed against the catalogue's own arithmetic: an UltraCam Eagle frame at GSD 30 +# gives 20010 px x 0.30 m = 6003 m across, which agrees with sizing the same frame from +# `104.052 mm x (height above ground / focal length)`. In metres it would be 100x. +fly_gsd_m <- function(gsd) gsd / 100 + + +# Resolve each row to a recording format. +# +# Keyed on `camera_calibration_url`, not on `media` or `focal_length`. `media` is a +# single value (`Digital - Colour`) across all 14 cameras in the record, and focal +# length is ambiguous — catalogue focal 92 spans an 87.1 mm DMC II and a 100.3 mm +# DMC III, a 15% difference. The calibration file is exact, and it is the only key that +# separates the two cameras the catalogue files under serial 20814295: an UltraCam +# Eagle through 2017 and a different body in 2018, whose own report numbers it 22814295. +# +# Where no calibration URL is present — about a fifth of digital frames — the catalogue +# focal length is the only remaining discriminator, so it is used and the row is marked +# inferred. That direction is defensible for WIDTH, which spreads 1-3% at a given focal, +# and not for PIXEL COUNT, which spreads 32-83%; the fallback rows carry no pixel counts +# for exactly that reason, which keeps them off the `px * GSD` route by construction. +# +# Returns one row per input row, all-NA where nothing resolved. +fly_camera_format <- function(centroids_sf) { + n <- nrow(centroids_sf) + none <- data.frame( + width_mm = rep(NA_real_, n), height_mm = rep(NA_real_, n), + px_cross = rep(NA_real_, n), px_along = rep(NA_real_, n), + camera = rep(NA_character_, n), width_source = rep(NA_character_, n), + # `resolved` and `inferred` are separate on purpose. A row that resolved to nothing + # is also `inferred = FALSE`, so `!inferred` — the natural filter for "trustworthy" + # — would sweep up every unresolved frame as well. + resolved = rep(FALSE, n), inferred = rep(FALSE, n), stringsAsFactors = FALSE + ) + if (n == 0 || !"media" %in% names(centroids_sf)) { + return(none) + } + + media <- as.character(centroids_sf$media) + # Film is sized from `negative_size`; this table describes sensors only. Restricting + # to digital also stops a fallback row keyed on focal length from quietly resolving a + # film frame that happens to share the focal length. + digital <- !is.na(media) & !(media %in% fly_film_media()) + if (!any(digital)) { + return(none) + } + + out <- none + tbl <- fly_camera_table() + calib <- tbl[tbl$key_type == "calib_file", ] + fb <- tbl[tbl$key_type == "focal_length", ] + + take <- function(rows, src, from, inferred) { + out$width_mm[rows] <<- from$width_mm + out$height_mm[rows] <<- from$height_mm + out$px_cross[rows] <<- from$px_cross + out$px_along[rows] <<- from$px_along + out$camera[rows] <<- from$camera + out$width_source[rows] <<- src + out$resolved[rows] <<- TRUE + out$inferred[rows] <<- inferred + } + + matched <- rep(FALSE, n) + if ("camera_calibration_url" %in% names(centroids_sf)) { + u <- as.character(centroids_sf$camera_calibration_url) + has_url <- digital & !is.na(u) & nzchar(u) + key <- rep(NA_character_, n) + # `basename(character(0))` is character(0), so guard rather than assign into a + # zero-length subscript. + if (any(has_url)) { + key[has_url] <- sub("\\.zip$", "", basename(u[has_url])) + } + m <- match(key, calib$key) + matched <- !is.na(m) + if (any(matched)) { + take(matched, key[matched], calib[m[matched], ], FALSE) + } + + # A frame whose calibration was deliberately withheld must not fall through to + # focal-length inference. The two withheld medium-format cameras are about 53 mm + # wide against the ~104 mm large-format bodies that share their focal neighbourhood, + # so inferring one would be ~1.95x too wide and 3.8x too much ground area. Today + # they happen to return NA because no fallback row exists at their focal lengths — + # that is safety by coincidence of the current table, and this makes it structural. + ex <- fly_camera_excluded() + refused <- !is.na(key) & key %in% ex$key[ex$key_type == "calib_file"] + out$width_source[refused] <- paste0("withheld:", key[refused]) + matched <- matched | refused + } + + if ("focal_length" %in% names(centroids_sf) && nrow(fb)) { + # Match numerically rather than on a formatted string: `as.character(100)` and + # `as.character(100L)` agree, but a double that prints as "1e+02" would not. + m <- match(as.numeric(centroids_sf$focal_length), as.numeric(fb$key)) + use <- digital & !matched & !is.na(m) + if (any(use)) { + take(use, paste0("focal_length=", fb$key[m[use]]), fb[m[use], ], TRUE) + } + } + + out +} diff --git a/R/fly_footprint.R b/R/fly_footprint.R index e828878..ec88dea 100644 --- a/R/fly_footprint.R +++ b/R/fly_footprint.R @@ -100,30 +100,78 @@ fly_dem_sample <- function(dem, rects) { list(elev = elev, covered = covered) } -# Build axis-aligned squares of `half_side` metres about each coordinate pair. -# A half-side that is NA or non-finite yields an empty polygon — the #30 -# contract for a frame whose recording format could not be resolved, and the -# only safe answer for a frame whose metadata divides by zero. An infinite -# rectangle is not merely meaningless: it reaches the sampler, which then tries -# to size a raster to it. -fly_rectangles <- function(coords, half_side) { +# Build rectangles of `half_cross` by `half_along` metres about each coordinate pair. +# +# `half_cross` spans the across-flight axis and `half_along` the along-flight one. Film +# is square, so the two are equal and the distinction costs nothing; a digital sensor is +# not — the Leica DMC III is 100.3 x 56.9 mm — and drawing it square would be 76% too +# deep. Defaulting `half_along` to `half_cross` keeps every film caller unchanged. +# +# A half-dimension that is NA, non-finite **or zero** yields an empty polygon. Zero +# matters as much as the others and is easy to miss: `is.finite(0)` is TRUE, so a zero +# would build a rectangle with five identical vertices, which `st_is_empty()` reports as +# FALSE and `fly_warn_unsized()` therefore never mentions, while it silently covers +# nothing downstream. That is strictly worse than the empty geometry #30 chose, and it +# is reachable — `ground_sample_distance` is 0 on every frame of some digital rolls. +# +# `bearing` rotates the rectangle to the flight line, and is applied only where the two +# half-dimensions differ. A square is unchanged by rotation up to vertex order, so +# leaving film alone keeps its output identical rather than merely equivalent. +# +# Vertex order is preserved as BL, BR, TR, TL, BL in the rectangle's own frame, because +# `fly_georef()` maps image corners onto footprint corners positionally and that order +# is its contract. +fly_rectangles <- function(coords, half_cross, half_along = half_cross, bearing = NULL) { sf::st_sfc(lapply(seq_len(nrow(coords)), function(i) { - w <- half_side[i] - if (!is.finite(w)) { + hc <- half_cross[i] + ha <- half_along[i] + if (!is.finite(hc) || !is.finite(ha) || hc <= 0 || ha <= 0) { return(sf::st_polygon()) } cx <- coords[i, 1] cy <- coords[i, 2] - sf::st_polygon(list(matrix(c( - cx - w, cy - w, - cx + w, cy - w, - cx + w, cy + w, - cx - w, cy + w, - cx - w, cy - w - ), ncol = 2, byrow = TRUE))) + xy <- matrix( + c(-hc, -ha, hc, -ha, hc, ha, -hc, ha, -hc, -ha), + ncol = 2, byrow = TRUE + ) + + b <- if (is.null(bearing)) NA_real_ else bearing[i] + if (!isTRUE(all.equal(hc, ha)) && is.finite(b)) { + # `bearing` is degrees clockwise from north, so the along-track axis is the local + # +y. Rotating (x, y) by b clockwise sends (0, 1) to (sin b, cos b) — the heading + # itself — which is what puts the long axis across the flight line rather than + # along it. + rad <- b * pi / 180 + rot <- matrix(c(cos(rad), sin(rad), -sin(rad), cos(rad)), nrow = 2) + xy <- xy %*% rot + } + + sf::st_polygon(list(cbind(xy[, 1] + cx, xy[, 2] + cy))) }), crs = 3005) } +# Which footprints are squares. +# +# This is the predicate `fly_georef()` needs, and "is it axis-aligned" is not — that is +# a proxy for "was it rotated", and a flight line running exactly north or east produces +# a rotated footprint whose edges are still axis-parallel, so the proxy misses it. +# +# Squareness is the property that actually matters. `georef_one()` maps image corners +# onto footprint corners positionally and then shifts that mapping by a +# 90-degree-quantized bearing, a scheme calibrated against north-up 9x9 negatives. With +# a square footprint a wrong-by-90 shift is harmless; with a 1.76:1 rectangle it maps a +# landscape image onto a portrait quad. And only non-square footprints are rotated, so +# this catches the double-rotation case too. +fly_is_square <- function(footprints) { + g <- sf::st_geometry(sf::st_transform(footprints, 3005)) + vapply(seq_along(g), function(i) { + if (sf::st_is_empty(g[i])) return(TRUE) + xy <- sf::st_coordinates(g[[i]])[, 1:2, drop = FALSE] + d <- sqrt(rowSums(diff(xy)^2)) + isTRUE(all.equal(max(d), min(d))) + }, logical(1)) +} + #' Estimate photo footprint polygons from centroids and scale #' #' Creates rectangular polygons representing the estimated ground coverage @@ -136,8 +184,10 @@ fly_rectangles <- function(coords, half_side) { #' 9" x 9"). Applies to film frames, and to every frame when there is no #' `media` column. It never sizes a digital frame — see `format_size`. #' @param format_size Named numeric vector of recording-format widths in inches, -#' keyed by `media` value, merged over the shipped film defaults. Supply this -#' to size frames whose format `fly` does not know — see Details. +#' keyed by `media` value, merged over the shipped film defaults. Frames it names are +#' sized from the reported `scale`, as film is, and it takes precedence over the +#' shipped camera table — it is the escape hatch for a camera `fly` does not know. +#' See Details. #' @param dem Optional elevation raster used to size each frame from its true #' height above ground rather than the reported scale. A `terra::SpatRaster`, #' a file path, or a `/vsicurl/` URL. Requires `flying_height` and @@ -177,15 +227,46 @@ fly_rectangles <- function(coords, half_side) { #' #' \describe{ #' \item{the `media` value}{format resolved from the format table} +#' \item{`"inferred_format"`}{digital frame with no calibration, sized from a +#' format inferred from its `focal_length`} #' \item{`"assumed_default"`}{no `media` column; `negative_size` applied} #' \item{`"unknown_format"`}{`media` present but unknown; empty geometry} #' } #' -#' Shipped defaults cover film only. Digital frames resolve to -#' `"unknown_format"` rather than an invented number, because the sensor width -#' they would need is not in the centroid metadata — and neither is the pixel -#' count that would let `ground_sample_distance` stand in for it. Supply -#' `format_size` if you know the camera: +#' @section Digital frames: +#' +#' A digital frame has no negative, and the catalogue mixes film and digital in one +#' layer — 223,667 of 1,670,471 frames province-wide are `Digital - Colour`. Sensor +#' dimensions are not in the centroid metadata, but they are recoverable from the +#' calibration report each frame links to through `camera_calibration_url`, and `fly` +#' ships them (`inst/extdata/camera_formats.csv`, built by +#' `data-raw/make_camera_formats.R`). +#' +#' Digital frames are sized as `pixel count x ground_sample_distance`, which needs +#' neither `scale` nor a DEM. **`scale` is never used for a digital frame `fly` sized +#' itself.** That field is not the true image scale for digital: measured against +#' terrain on 40 UltraCam Eagle frames it gives 34% of true width, because it is a +#' derived nominal figure — the pixel pitch it implies is about 12.5 um for every +#' camera regardless of model, against real pitches of 3.9 to 12 um. +#' `ground_sample_distance` is in centimetres. +#' +#' Where a frame carries no `camera_calibration_url` — about a fifth of digital frames — +#' the format is inferred from `focal_length` and `footprint_basis` records +#' `"inferred_format"`. Sensor width spreads only 1-3% at a given focal length, but +#' pixel count spreads 32-83%, so an inferred frame can only be sized through a DEM +#' (`width x height above ground / focal length`) and never from its GSD. +#' +#' `width_source` names the calibration file or fallback rule per row, so every +#' footprint traces back to a source. Calibrations that could not be corroborated are +#' listed in `inst/extdata/camera_formats_excluded.csv` with the reason, and frames +#' naming one are refused rather than inferred. +#' +#' **Digital footprints are not square** — sensors run from 1.10:1 (Leica DMC II) to +#' 1.80:1 (Intergraph DMC) — so they are rotated onto the flight line using +#' [fly_bearing()]. Where no bearing can be computed the rectangle stays axis-aligned +#' and `width_source` says so. Film stays square and is unaffected. +#' +#' Supply `format_size` to size a frame `fly` cannot, or to override it: #' #' ```r #' fly_footprint(photos, format_size = c("Digital - Colour" = 3.54)) @@ -232,6 +313,9 @@ fly_rectangles <- function(coords, half_side) { #' \describe{ #' \item{`"nominal_scale"`}{sized from the reported scale (no `dem`, or a #' fallback — see below)} +#' \item{`"gsd_scaled"`}{digital frame sized from its pixel count and ground +#' sample distance; used neither the reported scale nor a DEM, so `height_agl` +#' and `dem_coverage` are `NA`} #' \item{`"dem_agl"`}{sized from height above ground} #' \item{`"no_dem_coverage"`}{`dem` supplied but does not cover the frame} #' \item{`NA`}{no footprint to place — see `footprint_basis`} @@ -349,7 +433,13 @@ fly_footprint <- function(centroids_sf, negative_size = 9, format_size = NULL, input_crs <- sf::st_crs(centroids_sf) pts_3005 <- sf::st_transform(centroids_sf, 3005) coords <- sf::st_coordinates(pts_3005) - scale_num <- as.numeric(stringr::str_remove(centroids_sf$scale, "1:")) + # An unparseable `scale` is an expected input, not an exception: it yields NA, the + # frame gets no footprint, and that is reported by name below. Base R's + # "NAs introduced by coercion" would arrive alongside that as a second, vaguer warning + # pointing at no column in particular. + scale_num <- suppressWarnings( + as.numeric(stringr::str_remove(centroids_sf$scale, "1:")) + ) n <- nrow(centroids_sf) film <- fly_film_media() @@ -371,7 +461,24 @@ fly_footprint <- function(centroids_sf, negative_size = 9, format_size = NULL, basis <- as.character(ifelse(is.na(width_in), "unknown_format", media)) } - unresolved <- is.na(width_in) + # Sensor dimensions for digital frames, from the shipped camera table. Consulted only + # where `formats` did not already resolve the row, so a caller's own `format_size` + # still wins \u2014 it is the documented escape hatch for a camera `fly` does not know. + fmt <- fly_camera_format(centroids_sf) + from_table <- is.na(width_in) & !is.na(fmt$width_mm) + + width_source <- rep(NA_character_, n) + width_source[from_table] <- fmt$width_source[from_table] + # A frame naming a withheld calibration is not `from_table` — it resolved to nothing — + # so without this the refusal is computed and then dropped, and the frame is + # indistinguishable from one whose media was simply unknown. + withheld <- is.na(width_in) & !is.na(fmt$width_source) & + startsWith(fmt$width_source, "withheld:") + width_source[withheld] <- fmt$width_source[withheld] + basis[from_table] <- ifelse(fmt$inferred[from_table], "inferred_format", + as.character(centroids_sf$media)[from_table]) + + unresolved <- is.na(width_in) & !from_table if (any(unresolved)) { unknown <- sort(unique(as.character(centroids_sf$media)[unresolved])) warning( @@ -383,19 +490,85 @@ fly_footprint <- function(centroids_sf, negative_size = 9, format_size = NULL, ) } - half_side <- width_in * scale_num * 0.0254 / 2 + # Film, and anything `format_size` names: ground width is the format width times the + # reported scale, and the negative is square. + half_cross <- width_in * scale_num * 0.0254 / 2 + half_along <- half_cross + + # Digital: ground width is the pixel count times the ground sample distance. This uses + # neither `scale` nor a DEM. + # + # `scale` is deliberately not used. Measured against terrain on 40 UltraCam Eagle + # frames, sizing a digital frame from the catalogue's `SCALE` gives 34% of its true + # width, because that field is a derived nominal figure rather than the image scale \u2014 + # the pixel pitch it implies is ~12.5 um for every camera regardless of model, against + # real pitches of 3.9 to 12 um. A third-size footprint would still draw, still overlap + # its neighbours and still yield a coverage percentage, which is the failure #30 was + # written to prevent. + gsd_m <- rep(NA_real_, n) + if ("ground_sample_distance" %in% names(centroids_sf)) { + gsd_m <- fly_gsd_m(as.numeric(centroids_sf$ground_sample_distance)) + } + by_gsd <- from_table & !is.na(fmt$px_cross) & !is.na(fmt$px_along) & + !is.na(gsd_m) & gsd_m > 0 + half_cross[by_gsd] <- fmt$px_cross[by_gsd] * gsd_m[by_gsd] / 2 + half_along[by_gsd] <- fmt$px_along[by_gsd] * gsd_m[by_gsd] / 2 + + # Flight-line azimuth, for rotating a non-square footprint onto the flight line. + # + # Decided from the FORMAT's aspect ratio, not from the half-dimensions. The + # half-dimensions are NA for every camera-table row until a sizing route fills them, + # and the DEM route fills them *after* this point — so keying on them would leave + # `non_square` FALSE for exactly the frames the DEM exists to size, drawing them + # axis-aligned while `fly_bearing()` had a perfectly good azimuth for them. It would + # also make the answer depend on what else was in the batch. The aspect ratio is known + # before any route runs, which is what makes it the right thing to key on. + # + # This is the same NA-by-construction fact that has now bitten three separate + # conditions in this function; a value that only some routes populate is not a safe + # thing to branch on. + fmt_aspect_cross <- ifelse(is.na(width_in), fmt$width_mm, width_in * 25.4) + fmt_aspect_along <- ifelse(is.na(width_in), fmt$height_mm, width_in * 25.4) + non_square <- !is.na(fmt_aspect_cross) & !is.na(fmt_aspect_along) & + abs(fmt_aspect_cross - fmt_aspect_along) > + sqrt(.Machine$double.eps) * pmax(fmt_aspect_cross, fmt_aspect_along) + + # `fly_bearing()` stops rather than returning NA when its columns are absent, and + # `fly_footprint()` requires only `scale`, so the guard belongs here. + bearing <- rep(NA_real_, n) + if (any(non_square) && all(c("film_roll", "frame_number") %in% names(centroids_sf))) { + bearing <- fly_bearing(centroids_sf)$bearing + } + if (any(non_square & !is.finite(bearing))) { + width_source[non_square & !is.finite(bearing)] <- paste0( + width_source[non_square & !is.finite(bearing)], "; axis_aligned_no_bearing" + ) + } - # Keyed on half_side, not width_in: an unparseable `scale` also leaves a frame - # with no footprint, and a frame with no footprint has had no terrain - # treatment to report. - terrain <- as.character(ifelse(is.na(half_side), NA_character_, "nominal_scale")) + # Keyed on the half-dimensions, not on width_in: an unparseable `scale` also leaves a + # frame with no footprint, and a frame with no footprint has had no terrain treatment + # to report. + # + # A frame sized from its ground sample distance did not come from the reported scale + # and did not come from a DEM, so it gets its own value rather than borrowing + # `nominal_scale` \u2014 which is documented as "sized from the reported scale" and would + # be a false claim about the one route that deliberately avoids it. + terrain <- as.character(ifelse(is.na(half_cross), NA_character_, "nominal_scale")) + terrain[by_gsd] <- "gsd_scaled" height_agl <- rep(NA_real_, n) dem_coverage <- rep(NA_real_, n) if (!is.null(dem)) { focal_m <- centroids_sf$focal_length / 1000 + + # Format width in metres, whichever way the row was resolved: `negative_size` is + # inches, the camera table is millimetres. + fmt_cross_m <- ifelse(is.na(width_in), fmt$width_mm / 1000, width_in * 0.0254) + fmt_along_m <- ifelse(is.na(width_in), fmt$height_mm / 1000, width_in * 0.0254) + resize <- function(e) { - width_in * ((centroids_sf$flying_height - e) / focal_m) * 0.0254 / 2 + k <- (centroids_sf$flying_height - e) / focal_m + list(cross = fmt_cross_m * k / 2, along = fmt_along_m * k / 2) } # Two passes. The first averages the DEM over the nominal-scale rectangle, @@ -407,14 +580,41 @@ fly_footprint <- function(centroids_sf, negative_size = 9, format_size = NULL, # by a median 14%; the second pass moves it by at most a further 0.53%, and # a third by 0.03%. It is worth one more extract, not the emphasis a bare # "iterates until it converges" would imply. - first <- fly_dem_sample(dem, fly_rectangles(coords, half_side)) - second <- fly_dem_sample(dem, fly_rectangles(coords, resize(first$elev))) + # The sampling windows are built exactly as the returned footprint is — same two + # half-dimensions, same rotation. Sampling an unrotated rectangle and returning a + # rotated one would make `dem_coverage` and the mean elevation describe a shape the + # caller never receives, which is the returned-versus-measured mismatch this + # function already guards against elsewhere. + # Which frames the DEM route may size. + # + # NOT `!is.na(half_cross)`. That is the nominal-scale half-side, and it is NA for + # every row resolved from the camera table — `width_in` is NA there by definition — + # so keying on it makes the DEM route unreachable for exactly the frames it exists + # to serve. An inferred-format frame carries no pixel count and so has no GSD route + # at all; the DEM is its only route, and it would have come back empty with a DEM + # supplied and no warning. + dem_eligible <- !by_gsd & (!is.na(half_cross) | from_table) + + # A camera-table frame has no nominal rectangle to sample the first pass over, so + # seed one from the whole flying height — terrain at sea level, which is the largest + # plausible window and therefore certain to contain the true footprint. The second + # pass then averages over the rectangle the first produced, as it does for film. + seed_cross <- half_cross + seed_along <- half_along + need_seed <- dem_eligible & is.na(seed_cross) + if (any(need_seed)) { + at_sea_level <- resize(0) + seed_cross[need_seed] <- at_sea_level$cross[need_seed] + seed_along[need_seed] <- at_sea_level$along[need_seed] + } + + first <- fly_dem_sample(dem, fly_rectangles(coords, seed_cross, seed_along, bearing)) + r1 <- resize(first$elev) + second <- fly_dem_sample(dem, fly_rectangles(coords, r1$cross, r1$along, bearing)) # Keep the first pass wherever the second could not improve on it, so a # frame is never lost to the resize alone. elev <- ifelse(is.na(second$elev), first$elev, second$elev) - - sized <- !is.na(half_side) agl <- centroids_sf$flying_height - elev candidate <- resize(elev) @@ -424,9 +624,17 @@ fly_footprint <- function(centroids_sf, negative_size = 9, format_size = NULL, # an empty geometry here is indistinguishable from one whose recording # format was never resolved, so the frame would vanish under a warning # pointing at `format_size` rather than at its own metadata. - corrected <- sized & is.finite(candidate) & candidate > 0 - uncovered <- sized & !corrected & is.na(elev) - unusable <- sized & !corrected & !uncovered + # + # `!by_gsd` is what keeps the DEM from overwriting a frame already sized from its + # ground sample distance. Without it, passing `dem` would silently switch a digital + # frame onto a different route — and every downstream consumer forwards `dem`, so + # that is the ordinary path rather than an edge case. The two routes agree to about + # 1%, but agreeing is not the same as being interchangeable: the GSD route is the + # measurement and the DEM route is the estimate. + corrected <- dem_eligible & is.finite(candidate$cross) & candidate$cross > 0 & + is.finite(candidate$along) & candidate$along > 0 + uncovered <- dem_eligible & !corrected & is.na(elev) + unusable <- dem_eligible & !corrected & !uncovered # Coverage has to describe the footprint that is actually returned. A # corrected frame ships the second pass's rectangle, so it takes the second @@ -441,7 +649,7 @@ fly_footprint <- function(centroids_sf, negative_size = 9, format_size = NULL, # a frame we cannot correct is still a frame. if (any(uncovered)) { warning( - sum(uncovered), " of ", sum(sized), " frames fall outside the DEM's ", + sum(uncovered), " of ", sum(dem_eligible), " frames fall outside the DEM's ", "coverage and were sized from nominal scale instead. See ", "`footprint_terrain`.", call. = FALSE @@ -449,7 +657,7 @@ fly_footprint <- function(centroids_sf, negative_size = 9, format_size = NULL, } if (any(unusable)) { warning( - sum(unusable), " of ", sum(sized), " frames have `flying_height` or ", + sum(unusable), " of ", sum(dem_eligible), " frames have `flying_height` or ", "`focal_length` values that give no usable height above ground \u2014 ", "missing, zero, or terrain at or above the aircraft. Check that ", "`flying_height` is metres above sea level. Sized from nominal scale ", @@ -476,18 +684,58 @@ fly_footprint <- function(centroids_sf, negative_size = 9, format_size = NULL, ) } - half_side[corrected] <- candidate[corrected] + half_cross[corrected] <- candidate$cross[corrected] + half_along[corrected] <- candidate$along[corrected] height_agl[corrected] <- agl[corrected] # Reported for every frame that had a footprint to sample, not only the # corrected ones: `no_dem_coverage` is a measured zero, and leaving it NA # makes the documented "filter on dem_coverage" workflow impossible. - dem_coverage[sized] <- covered[sized] - terrain[sized] <- "nominal_scale" + # + # Excluding `by_gsd`, though. A DEM was sampled under those frames, but the number + # describes a window that had no bearing on the geometry returned — the footprint + # came from the pixel count and the ground sample distance. Reporting it would be a + # coverage figure for a shape the caller never receives. + dem_coverage[dem_eligible] <- covered[dem_eligible] + terrain[dem_eligible] <- "nominal_scale" terrain[corrected] <- "dem_agl" terrain[uncovered] <- "no_dem_coverage" terrain[unusable] <- "nominal_scale" } + # A frame with no rectangle has had no terrain treatment to report, whichever route + # failed to produce one. Keeping the invariant "footprint_terrain is NA exactly where + # the geometry is empty" is what lets a caller read the column at all. + no_geom <- is.na(half_cross) | is.na(half_along) | half_cross <= 0 | half_along <= 0 + terrain[no_geom] <- NA_character_ + height_agl[no_geom] <- NA_real_ + dem_coverage[no_geom] <- NA_real_ + + # A frame whose format resolved but which could not be sized is the quiet case: its + # `footprint_basis` names a real format and `width_source` names a calibration, so + # nothing about the row says the geometry is empty. It is not covered by the + # unknown-format warning above, and with no `dem` it is covered by none of the terrain + # warnings either. + # Split by cause: a film frame reaches this state through an unparseable `scale`, and + # telling its owner to supply a ground sample distance points at the wrong column. + unsized_digital <- from_table & no_geom + unsized_film <- !is.na(width_in) & no_geom + if (any(unsized_digital)) { + warning( + sum(unsized_digital), " of ", n, " frames have a known recording format but no ", + "way to size it, so they have no footprint. A digital frame needs either a ", + "`ground_sample_distance` and a calibrated pixel count, or `dem` together with ", + "`flying_height` and `focal_length`. See `footprint_basis` and `width_source`.", + call. = FALSE + ) + } + if (any(unsized_film)) { + warning( + sum(unsized_film), " of ", n, " frames have a known recording format but no ", + "usable `scale`, so they have no footprint. Expected a value like \"1:12000\".", + call. = FALSE + ) + } + # Assign the reporting columns onto the attribute frame rather than passing # them to `st_sf()` as trailing arguments. `st_sf()` builds its attribute # frame with `else if (inherits(x[[1]], c("tbl_df", "tbl"))) x[[1]]`, so a @@ -498,10 +746,14 @@ fly_footprint <- function(centroids_sf, negative_size = 9, format_size = NULL, attrs <- sf::st_drop_geometry(pts_3005) attrs$footprint_basis <- basis attrs$footprint_terrain <- terrain + attrs$width_source <- width_source attrs$height_agl <- height_agl attrs$dem_coverage <- dem_coverage - result <- sf::st_sf(attrs, geometry = fly_rectangles(coords, half_side)) + result <- sf::st_sf( + attrs, + geometry = fly_rectangles(coords, half_cross, half_along, bearing) + ) sf::st_transform(result, input_crs) } diff --git a/R/fly_georef.R b/R/fly_georef.R index dd757bc..843217f 100644 --- a/R/fly_georef.R +++ b/R/fly_georef.R @@ -119,6 +119,29 @@ fly_georef <- function(fetch_result, photos_sf, footprints <- fly_footprint(photos_sf, dem = dem) |> sf::st_transform(3005) fly_warn_unsized(footprints, "georeferencing") + # `georef_one()` maps image corners onto footprint corners positionally and shifts + # that mapping by a 90-degree-quantized bearing — a scheme calibrated against north-up + # 9x9 negatives. A non-square footprint breaks it twice over: the footprint is already + # rotated onto its flight line so the bearing would be counted a second time, and a + # wrong-by-90 mapping that is harmless on a square maps a landscape image onto a + # portrait quad on a 1.76:1 rectangle. + # + # Skipped rather than guessed at. The right corner mapping for a pre-rotated + # rectangle depends on camera mounting relative to flight direction, which is what + # `bearing_to_rotation()` was empirically fitted to and cannot be re-derived without + # imagery to check against. Digital frames had no footprint at all before fly#32, so + # this is the same coverage as before rather than a regression — but explicit now + # instead of silent. See fly#38. + rotated <- !fly_is_square(footprints) + if (any(rotated)) { + warning( + sum(rotated), " of ", nrow(footprints), " frames have a non-square footprint ", + "and are excluded from georeferencing: the corner mapping is calibrated for ", + "square, axis-aligned footprints. See fly#38.", + call. = FALSE + ) + } + # Match fetch results to photos by airp_id ids <- fetch_result$airp_id @@ -163,6 +186,8 @@ fly_georef <- function(fetch_result, photos_sf, # Find matching footprint fp_idx <- which(photos_sf[["airp_id"]] == results$airp_id[i]) if (length(fp_idx) == 0) next + # Warned about once, above, rather than per frame. + if (rotated[fp_idx[1]]) next fp <- footprints[fp_idx[1], ] # No footprint means no ground control to warp onto; leave success = FALSE diff --git a/data-raw/make_camera_formats.R b/data-raw/make_camera_formats.R new file mode 100644 index 0000000..79404bd --- /dev/null +++ b/data-raw/make_camera_formats.R @@ -0,0 +1,606 @@ +#!/usr/bin/env Rscript +# +# make_camera_formats.R +# +# Build inst/extdata/camera_formats.csv — the recording-format dimensions for every +# digital camera in the BC air photo catalogue, read out of the camera calibration +# reports the catalogue itself links to. +# +# Why this exists: `AIMG_PHOTO_CENTROIDS_SP` carries no sensor size, so #30 refused to +# size digital frames rather than invent one. `CAMERA_CALIBRATION_URL` is populated for +# ~80% of them and its basename identifies the calibration exactly, so the number is +# recoverable — see fly#32. +# +# The numbers are PARSED, never typed. A transcription error is invisible to any test +# that reads the CSV the error produced, so the human reviews a parse rather than a +# keyboard entry, and the QA section below constrains every field with at least two +# checks drawn from different sources. +# +# Network: BC WFS (openmaps.gov.bc.ca) for the frame attributes, and the same host for +# the calibration zips (~25 MB, cached under data-raw/.cache/). +# +# Run from fly repo root: Rscript data-raw/make_camera_formats.R + +pkgload::load_all(quiet = TRUE) # source tree, never the installed package + +library(dplyr) + +CACHE <- "data-raw/.cache" +OUT <- "inst/extdata/camera_formats.csv" +WFS <- "https://openmaps.gov.bc.ca/geo/pub/wfs" +LAYER <- "pub:WHSE_IMAGERY_AND_BASE_MAPS.AIMG_PHOTO_CENTROIDS_SP" +RETRIEVED <- "2026-08-30" # the snapshot date these rows describe + +fs::dir_create(fs::path(CACHE, "zip")) +fs::dir_create(fs::path(CACHE, "pdf")) +fs::dir_create("inst/extdata") + + +# --- WFS helpers ------------------------------------------------------------------ + +# One page of the layer as a CSV data frame. `count` is capped server-side at 10000, +# which is not reported as an error — reconcile against resultType=hits, never assume a +# short page means the end. +wfs_page <- function(cql, props, count = 10000, start = 0) { + q <- c( + service = "WFS", version = "2.0.0", request = "GetFeature", + typeNames = LAYER, CQL_FILTER = cql, propertyName = props, + sortBy = "AIRP_ID", count = count, startIndex = start, + outputFormat = "csv" + ) + url <- paste0(WFS, "?", paste0(names(q), "=", + vapply(q, utils::URLencode, character(1), + reserved = TRUE), + collapse = "&")) + f <- tempfile(fileext = ".csv") + on.exit(unlink(f), add = TRUE) + utils::download.file(url, f, mode = "wb", quiet = TRUE) + utils::read.csv(f, stringsAsFactors = FALSE, colClasses = c(SCALE = "character")) +} + +wfs_hits <- function(cql) { + q <- c(service = "WFS", version = "2.0.0", request = "GetFeature", + typeNames = LAYER, CQL_FILTER = cql, resultType = "hits") + url <- paste0(WFS, "?", paste0(names(q), "=", + vapply(q, utils::URLencode, character(1), + reserved = TRUE), + collapse = "&")) + as.integer(sub('.*numberMatched="([0-9]+)".*', "\\1", + paste(readLines(url, warn = FALSE), collapse = " "))) +} + +wfs_all <- function(cql, props) { + n <- wfs_hits(cql) + message(" ", format(n, big.mark = ","), " frames to page") + out <- lapply(seq(0, n, by = 10000), function(s) wfs_page(cql, props, 10000, s)) + d <- dplyr::bind_rows(out) + # Reconcile: a silently short page is the failure mode this guards. + stopifnot(nrow(d) == n, !anyDuplicated(d$FID)) + d +} + + +# --- 1. Discover every digital frame and its calibration -------------------------- + +message("Querying the catalogue for digital frames ...") +DIGITAL <- "MEDIA LIKE 'Digital%'" +props <- paste("PHOTO_YEAR", "FOCAL_LENGTH", "GROUND_SAMPLE_DISTANCE", "SCALE", + "FLYING_HEIGHT", "CAMERA_CALIBRATION_URL", sep = ",") +frames <- wfs_all(DIGITAL, props) + +frames$calib_url <- ifelse( + is.na(frames$CAMERA_CALIBRATION_URL) | !nzchar(frames$CAMERA_CALIBRATION_URL), + NA_character_, frames$CAMERA_CALIBRATION_URL +) +frames$key <- ifelse(is.na(frames$calib_url), NA_character_, + sub("\\.zip$", "", basename(frames$calib_url))) + +urls <- sort(unique(stats::na.omit(frames$calib_url))) +# The cache, and `frames$key`, are keyed on the basename rather than the full URL. Two +# calibrations in different directories sharing a basename would silently merge into +# one camera, so the assumption is asserted rather than left implicit. +stopifnot(!anyDuplicated(basename(urls))) +message(" ", nrow(frames), " digital frames, ", + sum(!is.na(frames$key)), " with a calibration, ", + length(urls), " distinct calibration files") + + +# --- 2. Fetch and unpack the calibration reports ---------------------------------- + +message("Fetching calibration reports (cached) ...") +for (u in urls) { + z <- fs::path(CACHE, "zip", basename(u)) + # Guard on non-empty, not existence: `download.file` truncates its target before it + # runs, so a failed fetch leaves a 0-byte file that an existence check blesses forever. + if (!fs::file_exists(z) || fs::file_size(z) == 0) { + tmp <- paste0(z, ".part") + utils::download.file(u, tmp, mode = "wb", quiet = TRUE) + if (fs::file_size(tmp) == 0) stop("empty download: ", u, call. = FALSE) + fs::file_move(tmp, z) + } + # Same guard-on-existence trap as the zip above, one level down: `unzip()` WARNS + # rather than errors on a bad archive, so a gate on `dir_exists()` caches an empty + # directory forever. The key then lands in the excluded file carrying the specific + # claim "present only as a scanned image", which is false and permanent. Extract to a + # scratch directory and only publish it once it actually holds a PDF. + d <- fs::path(CACHE, "pdf", sub("\\.zip$", "", basename(u))) + if (!fs::dir_exists(d) || !length(fs::dir_ls(d, regexp = "\\.pdf$", recurse = TRUE))) { + part <- paste0(d, ".part") + if (fs::dir_exists(part)) fs::dir_delete(part) + fs::dir_create(part) + utils::unzip(z, exdir = part) + if (!length(fs::dir_ls(part, regexp = "\\.pdf$", recurse = TRUE))) { + fs::dir_delete(part) + stop("no PDF extracted from ", basename(z), call. = FALSE) + } + if (fs::dir_exists(d)) fs::dir_delete(d) + fs::file_move(part, d) + } +} + + +# --- 3. Parsers ------------------------------------------------------------------- + +pdf_lines <- function(path) unlist(lapply(pdftools::pdf_text(path), function(p) strsplit(p, "\n")[[1]])) + +num <- function(x) as.numeric(gsub("[^0-9.]", "", x)) + +# Poppler renders a zero as a capital O in some of the older Vexcel reports — the 2013 +# UltraCam Eagle certificate reads `2001Opixel` where it means 20010. +# +# Note what happens without this: `num()` strips non-digits, so an unrepaired `2001O` +# becomes 2001 rather than failing. A silently-wrong pixel count is far worse than no +# row, so the substitution is made explicitly and check B — `px * pitch == the stated +# millimetres`, two numbers that did not pass through this — is what proves it right. +unglyph <- function(s) gsub("(?<=[0-9])[Oo]", "0", s, perl = TRUE) + +# Vexcel UltraCam calibration report. +# +# The panchromatic and multispectral blocks are formatted identically and the +# multispectral one sits directly below, so anchoring on the heading is load-bearing — +# taking the wrong block is the one parse error the numeric QA cannot see. +parse_vexcel <- function(lines) { + anchor <- grep("Large Format Panchromatic Output Image", lines) + if (!length(anchor)) return(NULL) + blk <- lines[anchor[1]:min(length(lines), anchor[1] + 25)] + # stop before the multispectral block if it follows within the window + stop_at <- grep("Multispectral", blk) + if (length(stop_at)) blk <- blk[seq_len(stop_at[1] - 1)] + + lt <- unglyph(grep("long track", blk, value = TRUE)[1]) + ct <- unglyph(grep("cross track", blk, value = TRUE)[1]) + ps <- grep("Pixel Size", blk, value = TRUE)[1] + fl <- grep("Focal Length|ck\\s*=", blk, value = TRUE)[1] + if (any(is.na(c(lt, ct, ps)))) return(NULL) + + g <- function(s, pat) regmatches(s, regexpr(pat, s, perl = TRUE)) + list( + height_mm = num(g(lt, "[0-9.]+\\s*mm")), + px_along = num(g(lt, "[0-9,]+\\s*pixel")), + width_mm = num(g(ct, "[0-9.]+\\s*mm")), + px_cross = num(g(ct, "[0-9,]+\\s*pixel")), + # Take the first number on the Pixel Size line rather than matching a unit. These + # reports write the micron sign three different ways and the 2013 one loses it + # altogether — `Pixel Size 5.200 m*5.200 m`, which reads as metres. Anchoring on + # the label and letting check B (px * pitch == the stated millimetres) and the + # plausibility bound settle the unit is safer than trusting the glyph. + pitch_um = num(g(ps, "[0-9.]+")), + focal_mm = if (!is.na(fl)) num(g(fl, "[0-9.]+\\s*mm")) else NA_real_, + # The report states millimetres independently of pixels and pitch, so check B is a + # real constraint on this row rather than a restatement of its own arithmetic. + stated_mm = TRUE + ) +} + +# QSI boresight report — a specifications table, used where the Vexcel appendix is +# abridged and carries no image-format block. +parse_qsi <- function(lines) { + txt <- paste(lines, collapse = " | ") + if (!grepl("Array Size", txt)) return(NULL) + gr <- function(pat) { + m <- regmatches(txt, regexpr(pat, txt, perl = TRUE)) + if (!length(m)) NA_character_ else m + } + arr <- gr("Array Size\\s*\\|?\\s*\\|?\\s*([0-9,]+)\\s*x\\s*([0-9,]+)") + ps <- gr("CCD Pixel Size[^0-9]*([0-9.]+)") + fl <- gr("Focal Length[^0-9]*([0-9.]+)") + if (is.na(arr)) return(NULL) + dims <- as.numeric(gsub(",", "", regmatches(arr, gregexpr("[0-9,]{4,}", arr))[[1]])) + pitch <- num(sub(".*?([0-9.]+)\\s*$", "\\1", ps)) + list( + px_cross = max(dims), px_along = min(dims), pitch_um = pitch, + width_mm = max(dims) * pitch / 1000, + height_mm = min(dims) * pitch / 1000, + focal_mm = num(sub("Focal Length[^0-9]*", "", fl)), + # This report gives array size and pitch but no image size, so the millimetres are + # DERIVED. Check B is vacuous here and the tests must not pretend otherwise. + stated_mm = FALSE + ) +} + +# Leica DMC II / DMC III calibration certificate. States pixel count, pixel size and +# image size independently, which is what makes check B a real constraint here. +# +# The micron sign in these reports is U+F06D — a Private Use Area codepoint from a +# Symbol font, not `µ` (U+00B5). So the pixel-size label reads as `Pixel Size [m]`, +# which a literal `\[m\]` misses and a human reading the extracted text sees as `[m]` +# and takes for METRES. Match the unit permissively and record it as microns. +UNIT_MICRON <- "\\[.{0,3}m\\]" + +parse_leica <- function(lines) { + rc <- grep("Number of rows/columns \\[pixels\\]", lines, value = TRUE)[1] + ps <- grep(paste0("Pixel Size\\s*", UNIT_MICRON), lines, value = TRUE)[1] + im <- grep("Image Size \\[mm\\]", lines, value = TRUE)[1] + fl <- grep("Focal Length \\[mm\\]", lines, value = TRUE)[1] + if (any(is.na(c(rc, ps, im)))) return(NULL) + pair <- function(s) as.numeric(regmatches(s, gregexpr("[0-9]+\\.?[0-9]*", s))[[1]]) + d <- pair(sub(".*\\[pixels\\]", "", rc)) + p <- pair(sub(paste0(".*", UNIT_MICRON), "", ps)) + i <- pair(sub(".*\\[mm\\]", "", im)) + list( + px_cross = max(d), px_along = min(d), pitch_um = p[1], + width_mm = max(i), height_mm = min(i), + focal_mm = if (!is.na(fl)) pair(sub(".*\\[mm\\]", "", fl))[1] else NA_real_, + stated_mm = TRUE + ) +} + +# Intergraph DMC — reports a "virtual" image, the resampled composite the frames are +# actually delivered as. The high-resolution (panchromatic) block is the first one. +parse_dmc <- function(lines) { + a <- grep("Virtual Focal Length", lines) + if (!length(a)) return(NULL) + blk <- lines[a[1]:min(length(lines), a[1] + 4)] + fl <- grep("Virtual Focal Length", blk, value = TRUE)[1] + sz <- grep("Virtual Sensor Size", blk, value = TRUE)[1] + ps <- grep("Virtual Pixel Size", blk, value = TRUE)[1] + if (any(is.na(c(fl, sz, ps)))) return(NULL) + d <- as.numeric(regmatches(sz, gregexpr("[0-9]+", sz))[[1]]) + pitch <- num(sub(paste0(".*", UNIT_MICRON), "", ps)) + list( + px_cross = max(d), px_along = min(d), pitch_um = pitch, + width_mm = max(d) * pitch / 1000, height_mm = min(d) * pitch / 1000, + focal_mm = num(sub(".*\\[m\\]", "", fl)) * 1000, + stated_mm = FALSE # derived, as parse_qsi + ) +} + +# The serial the report gives itself, which is not always the one the catalogue's URL +# basename implies: the 2018 UltraCam report is `UC-EpII-1-22814295-f80` where the +# catalogue files it under 20814295. Recorded so that mismatch stays visible. +report_serial <- function(lines) { + txt <- paste(lines, collapse = " ") + m <- regmatches(txt, regexpr("UC-[A-Za-z0-9]+-[0-9]+-[0-9]+(-f[0-9]+)?", txt)) + if (length(m)) return(m[1]) + m <- regmatches(txt, regexpr("(DMC ?I{0,3}|DMC0?[0-9]+) ?-? ?[0-9]{4,8}", txt)) + if (length(m)) m[1] else NA_character_ +} + +camera_name <- function(lines) { + txt <- paste(lines, collapse = " ") + pats <- c("UltraCam Eagle Prime II", "UltraCam Eagle-M3", "UltraCam Eagle M3", + "UltraCam Falcon M2", "UltraCamXp", "UltraCam Xp", "UltraCam Eagle", + "UltraCamEagle", "UltraCam X", "DMC III", "DMC II", "DMC") + for (p in pats) if (grepl(p, txt, fixed = TRUE)) return(p) + NA_character_ +} + +# --- 4. Parse every report --------------------------------------------------------- + +NUMERIC_FIELDS <- c("width_mm", "height_mm", "px_cross", "px_along", "pitch_um", "focal_mm") + +parse_one <- function(k) { + dir <- fs::path(CACHE, "pdf", k) + pdfs <- fs::dir_ls(dir, regexp = "\\.pdf$", recurse = TRUE, type = "file") + # A boresight report restates specs loosely; prefer a real calibration certificate, + # but fall back to one rather than losing the row. + pdfs <- c(pdfs[!grepl("oresight", pdfs)], pdfs[grepl("oresight", pdfs)]) + for (p in pdfs) { + lines <- tryCatch(pdf_lines(p), error = function(e) character(0)) + if (!length(lines)) next # image-only PDF: 0 extractable characters + for (fn in list(parse_vexcel, parse_leica, parse_dmc, parse_qsi)) { + r <- tryCatch(fn(lines), error = function(e) NULL) + if (is.null(r)) next + if (!all(vapply(r[NUMERIC_FIELDS], function(x) length(x) == 1 && is.finite(x), + logical(1)))) next + return(data.frame( + key = k, camera = camera_name(lines), report_serial = report_serial(lines), + source_pdf = fs::path_file(p), width_mm = r$width_mm, height_mm = r$height_mm, + px_cross = r$px_cross, px_along = r$px_along, pitch_um = r$pitch_um, + focal_mm = r$focal_mm, stated_mm = r$stated_mm, stringsAsFactors = FALSE + )) + } + } + NULL +} + +message("Parsing reports ...") +keys <- sort(unique(stats::na.omit(frames$key))) +parsed <- dplyr::bind_rows(lapply(keys, parse_one)) +unparsed <- setdiff(keys, parsed$key) +message(" parsed ", nrow(parsed), " of ", length(keys), " calibration files") + + +# --- 5. QA ------------------------------------------------------------------------ +# +# These gate the write. A row that cannot be corroborated is withheld rather than +# shipped: a footprint we know fails its own consistency check is exactly the mystery +# error this section exists to prevent. + +fail <- character(0) +note <- function(...) message(" ", ...) + +## B — px * pitch reproduces the image size the report states, both axes. +## Only meaningful where the report states millimetres independently; where the parser +## derived them the check is vacuous and is skipped rather than counted as a pass. +b <- parsed[parsed$stated_mm, ] +b_cross <- abs(b$px_cross * b$pitch_um / 1000 - b$width_mm) / b$width_mm +b_along <- abs(b$px_along * b$pitch_um / 1000 - b$height_mm) / b$height_mm +note("B px*pitch == stated mm: ", nrow(b), " of ", nrow(parsed), " rows constrainable, ", + "max rel. err ", format(max(c(b_cross, b_along)), digits = 3)) +if (any(c(b_cross, b_along) > 1e-6)) { + fail <- c(fail, paste("B failed:", paste(b$key[b_cross > 1e-6 | b_along > 1e-6], + collapse = ", "))) +} + +## C — the report's focal length against the catalogue's, which is an independent field. +cat_focal <- tapply(frames$FOCAL_LENGTH, frames$key, function(x) median(x, na.rm = TRUE)) +parsed$focal_catalogue <- as.numeric(cat_focal[parsed$key]) +# `focal_catalogue` is NA when a key's frames all carry a missing FOCAL_LENGTH, and +# `any(NA)` is NA rather than FALSE — which errors in an `if`. Absence of a catalogue +# focal is not a disagreement, so it resolves to FALSE. +parsed$focal_disagrees <- !is.na(parsed$focal_catalogue) & + abs(parsed$focal_mm - parsed$focal_catalogue) > 1 +note("C report focal vs catalogue: ", sum(!parsed$focal_disagrees), " of ", nrow(parsed), + " agree within 1 mm", + if (any(parsed$focal_disagrees)) + paste0(" (disagree: ", paste(parsed$key[parsed$focal_disagrees], collapse = ", "), ")") + else "") + +## D — plausibility bounds. Catches a unit slip (m / mm / um), which is the error class +## that survives B by being internally self-consistent. +parsed$aspect <- parsed$width_mm / parsed$height_mm +## The width bound is 30-200 mm, not the 80-170 mm the large-format cameras occupy. +## A tighter bound would reject a medium-format sensor as implausible — the two the +## catalogue actually holds are ~53 mm wide — and a guard that refuses valid data is +## worse than the unit slip it is trying to catch. 30-200 mm still separates a +## micrometre read as a millimetre (~5 mm) or a metre read as one (~0.0001 mm), which +## is what D exists for. +## `focal_mm` is bounded here too. It was previously the only shipped number with no +## gating check at all: D did not cover it, C reports a disagreement without failing the +## build (there is a known legitimate one), and F is skipped for any key whose frames +## carry no GSD. A focal that fell out of a missed unit anchor would have shipped. +bad_d <- with(parsed, width_mm < 30 | width_mm > 200 | aspect < 1 | aspect > 2 | + pitch_um < 3 | pitch_um > 13 | focal_mm < 40 | focal_mm > 300) +note("D plausibility bounds: ", sum(!bad_d), " of ", nrow(parsed), " within range") +if (any(bad_d)) fail <- c(fail, paste("D failed:", paste(parsed$key[bad_d], collapse = ", "))) + +## F — implied ground elevation. Invert for the ground the aircraft was over: +## terrain = flying_height - (GSD / pitch) * focal. Nothing here comes from the table +## except pitch and focal; FLYING_HEIGHT and GROUND_SAMPLE_DISTANCE are the catalogue's +## own. The result must be a plausible BC elevation. +## +## GROUND_SAMPLE_DISTANCE is in CENTIMETRES — see gsd_m() in R/fly_footprint.R. Getting +## that wrong is a factor of 100, so it is asserted rather than assumed. +f <- merge(frames[!is.na(frames$key), ], parsed[, c("key", "pitch_um", "focal_mm")], by = "key") +f <- f[!is.na(f$GROUND_SAMPLE_DISTANCE) & f$GROUND_SAMPLE_DISTANCE > 0 & + !is.na(f$FLYING_HEIGHT), ] +f$terrain <- f$FLYING_HEIGHT - (f$GROUND_SAMPLE_DISTANCE / 100) / (f$pitch_um * 1e-6) * + (f$focal_mm / 1000) +terr <- tapply(f$terrain, f$key, median) +parsed$terrain_implied <- round(as.numeric(terr[parsed$key])) +parsed$terrain_ok <- is.na(parsed$terrain_implied) | + (parsed$terrain_implied > -50 & parsed$terrain_implied < 2800) +note("F implied ground elevation plausible: ", sum(parsed$terrain_ok), " of ", nrow(parsed), + if (any(!parsed$terrain_ok)) + paste0(" (implausible: ", paste(parsed$key[!parsed$terrain_ok], collapse = ", "), ")") + else "") + +if (length(fail)) stop("QA failed, CSV not written:\n ", paste(fail, collapse = "\n "), + call. = FALSE) + +## Withhold anything F rejects. A frame with no footprint is honest; a frame with a +## footprint contradicted by its own metadata is not. +parsed$ship <- parsed$terrain_ok + + +# --- 6. Focal-length fallback rows ------------------------------------------------- +# +# ~20% of digital frames carry no calibration URL. Focal length is the only other +# discriminator, and it is a weak one — so the rows are DERIVED by rule here rather than +# chosen by hand, and each records how many calibrated cameras stand behind it and how +# far apart they are. Consumers mark frames sized this way as inferred. + +shipped <- parsed[parsed$ship, ] +by_key <- as.data.frame(table(frames$key), stringsAsFactors = FALSE) +names(by_key) <- c("key", "frames") +shipped <- merge(shipped, by_key, by = "key", all.x = TRUE) + +# Which catalogue focal lengths need a fallback: those on frames with no calibration. +# Extrapolation is OPT-IN, per catalogue focal length, with a reason. +# +# The first version of this was a refusal list, which encoded the one instance that had +# been measured rather than the property it stood for — so a focal length nobody had +# looked at yet got a row by default. Worse, the obvious generalisation does not work: +# extrapolation *distance* does not predict the error. The refused focal-83 row reached +# across a 3.6% focal gap and was 93% wrong, while focal 120 reaches across a larger +# 5.5% gap and is right. Only camera identity separates them, and the generator cannot +# see that. +# +# So the safe default is no row, and each exception is named here with the argument for +# it. A focal length absent from this list and absent from the calibrated set simply +# gets no fallback, and is recorded in the excluded file. +FALLBACK_EXTRAPOLATE <- c( + "120" = paste("no calibrated camera at focal 120; the 11,984 frames are 2011 and the", + "only large-format 120 mm body in the record is the Z/I DMC, whose own", + "report gives a virtual focal of exactly 120 mm - so its 165.888 mm", + "format is used. Marked extrapolated: this rests on camera identity,", + "not on the catalogue focal being close to 127") +) + +need <- sort(unique(frames$FOCAL_LENGTH[is.na(frames$key) & !is.na(frames$FOCAL_LENGTH)])) + +shipped <- parsed[parsed$ship, ] +by_key <- as.data.frame(table(frames$key), stringsAsFactors = FALSE) +names(by_key) <- c("key", "frames") +shipped <- merge(shipped, by_key, by = "key", all.x = TRUE) + +# Ground truth: what widths do calibrated frames at each catalogue focal actually have? +gt <- merge(frames[!is.na(frames$key), c("key", "FOCAL_LENGTH")], + shipped[, c("key", "width_mm", "height_mm")], by = "key") + +has_exact <- vapply(need, function(fl) any(gt$FOCAL_LENGTH == fl), logical(1)) +allowed <- as.character(need) %in% names(FALLBACK_EXTRAPOLATE) +build <- need[has_exact | allowed] +no_fallback <- need[!has_exact & !allowed] + +# A declared exception that never fires is a stale decision nobody is told about. +unused <- setdiff(names(FALLBACK_EXTRAPOLATE), as.character(need[!has_exact])) +if (length(unused)) { + stop("FALLBACK_EXTRAPOLATE declares focal ", paste(unused, collapse = ", "), + " but nothing needs extrapolating there any more.", call. = FALSE) +} + +fallback <- dplyr::bind_rows(lapply(build, function(fl) { + at <- gt[gt$FOCAL_LENGTH == fl, ] + exact <- nrow(at) > 0 + if (!exact) { + near <- gt$FOCAL_LENGTH[which.min(abs(gt$FOCAL_LENGTH - fl))] + at <- gt[gt$FOCAL_LENGTH == near, ] + } + if (!nrow(at)) return(NULL) + # Modal width by frame count. Selected by INDEX, not by reparsing the name back out + # of `table()`: `as.numeric(as.character(100.3392))` does not round-trip exactly, and + # comparing the reparsed value with `==` would silently match nothing and leave the + # paired height NA. + tw <- table(at$width_mm) + w <- at$width_mm[at$width_mm %in% as.numeric(names(tw)[which.max(tw)])][1] + if (is.na(w)) w <- at$width_mm[which.max(tabulate(match(at$width_mm, at$width_mm)))] + h <- at$height_mm[match(w, at$width_mm)] + data.frame( + key = as.character(fl), key_type = "focal_length", + camera = if (exact) "inferred from focal length" else "extrapolated from nearest focal", + report_serial = NA_character_, source_pdf = NA_character_, + width_mm = w, height_mm = h, + # Pixel counts are NOT carried: at a given focal they spread 32-83% where width + # spreads 1-3%, so a fallback frame can be sized by width x AGL/focal but never by + # px x GSD. Leaving them NA is what stops that route being taken. + px_cross = NA_real_, px_along = NA_real_, + pitch_um = NA_real_, focal_mm = fl, mm_stated = NA, + n_cameras = length(unique(at$width_mm)), + # Spread is dispersion among the SOURCE cameras, so a single-source row is + # structurally 0 however wrong the inference. Reporting 0 on an extrapolated row + # would give the least-supported rows in the file the most confident label, so an + # extrapolated row reports NA instead — unknown, which is the truth. + width_spread_pct = if (exact) round(100 * (max(at$width_mm) / min(at$width_mm) - 1), 1) else NA_real_, + extrapolated = !exact, + note = if (exact) "modal width among calibrated cameras at this catalogue focal" + else FALLBACK_EXTRAPOLATE[[as.character(fl)]], + stringsAsFactors = FALSE + ) +})) + + +# --- 7. Write ---------------------------------------------------------------------- + +calib_rows <- data.frame( + key = shipped$key, key_type = "calib_file", camera = shipped$camera, + report_serial = shipped$report_serial, source_pdf = shipped$source_pdf, + width_mm = shipped$width_mm, height_mm = shipped$height_mm, + px_cross = shipped$px_cross, px_along = shipped$px_along, + pitch_um = shipped$pitch_um, focal_mm = shipped$focal_mm, + # Whether the report stated the millimetres itself. Where it did not, `width_mm` is + # px * pitch and check B on that row compares the arithmetic with itself — the tests + # must skip it rather than count a vacuous pass. + mm_stated = shipped$stated_mm, + n_cameras = 1L, width_spread_pct = 0, extrapolated = FALSE, + note = ifelse(shipped$focal_disagrees, + paste0("catalogue records focal ", shipped$focal_catalogue, + "; the report says ", shipped$focal_mm, " - report preferred"), + NA_character_), + stringsAsFactors = FALSE +) + +out <- rbind(calib_rows, fallback) +out$retrieved <- RETRIEVED +out <- out[order(out$key_type, out$key), ] +utils::write.csv(out, OUT, row.names = FALSE, na = "") + +# Every discovered key must land in exactly one of the two files, each with a reason. +# A guard that lets an entry sit in neither is how drift becomes invisible. +# Each block is guarded on a non-empty key vector. `paste0()` returns length 1 for a +# zero-length argument, so an unguarded `data.frame(key = character(0), reason = ...)` +# aborts with "arguments imply differing number of rows: 0, 1" — and the second block +# fires in the HEALTHY case, the moment check F rejects nothing. That would kill the run +# after every download and all the QA, with camera_formats.csv already on disk and the +# excluded file and manifest not: exactly the inconsistent state the manifest guard +# exists to detect. +# `excl()` is what makes each block safe on an empty key vector. `paste0()` returns +# length 1 for a zero-length argument, so `data.frame(key = character(0), reason = ...)` +# aborts at CONSTRUCTION with "arguments imply differing number of rows: 0, 1" — +# subsetting afterwards is too late. The withheld block fires in the HEALTHY case, the +# moment check F rejects nothing, which would kill the run after every download and all +# the QA with camera_formats.csv already on disk and the excluded file and manifest not: +# exactly the inconsistent state the manifest guard exists to detect. +excl <- function(key, key_type, reason) { + if (!length(key)) { + return(data.frame(key = character(0), key_type = character(0), + reason = character(0), stringsAsFactors = FALSE)) + } + data.frame(key = key, key_type = key_type, reason = reason, stringsAsFactors = FALSE) +} + +excluded <- rbind( + # These carry their calibration only as a scanned image, so no text pass can reach it. + # An independent visual reading (fly#32, check E) recovered the specs and they are + # recorded here so the knowledge is not lost — but they are NOT shipped, because a row + # this generator cannot reproduce would break the guarantee that re-running it + # reproduces the table. Two of the three are medium-format bodies roughly half the + # width of everything else in the record, which is also why no fallback may reach for + # their focal lengths. + excl(unparsed, "calib_file", paste0( + "calibration present only as a scanned image, not machine-readable", + ifelse(unparsed == "10210206_2015", + "; visually read as UltraCam Eagle 20010x13080 @ 5.2um = 104.052x68.016mm f100.5", ""), + ifelse(unparsed == "11937933_2009", + "; visually read as AIC Pro (P65+) 8984x6732 @ 6.0um = 53.904x40.392mm f60.68", ""), + ifelse(unparsed == "12335326_2017", + "; visually read as PhaseOne IXU-RS-1000 11608x8708 @ 4.6um = 53.4x40.1mm f51.56", "") + )), + excl(parsed$key[!parsed$ship], "calib_file", + paste0("catalogue GSD/FLYING_HEIGHT contradict the report: implied ground ", + "elevation ", parsed$terrain_implied[!parsed$ship], " m")), + # Focal lengths deliberately left without a fallback row, in the same file so every + # refusal carries its reason rather than living only in a code comment. + excl(as.character(no_fallback), "focal_length", + paste0("no fallback row: no calibrated camera at focal ", no_fallback, + ", and extrapolating from a different focal is not warranted without ", + "knowing the camera - focal distance does not predict the error")) +) +excluded$frames <- as.integer(by_key$frames[match(excluded$key, by_key$key)]) +excluded$retrieved <- RETRIEVED +utils::write.csv(excluded[order(excluded$key), ], + "inst/extdata/camera_formats_excluded.csv", row.names = FALSE, na = "") + +# The keys the catalogue actually offered at RETRIEVED, written independently of the +# two dispositions above. Committing it lets the test suite re-check offline that every +# discovered key is still dispositioned — including after a hand-edit to either file, +# which is the drift a self-referential check cannot see. +utils::write.csv(data.frame(key = keys, retrieved = RETRIEVED), + "inst/extdata/camera_formats_manifest.csv", row.names = FALSE) + +stopifnot(setequal( + keys, + c(out$key[out$key_type == "calib_file"], + excluded$key[excluded$key_type == "calib_file"]) +)) +stopifnot(setequal(as.character(need), + c(out$key[out$key_type == "focal_length"], as.character(no_fallback)))) + +message("\nWrote ", OUT, ": ", sum(out$key_type == "calib_file"), " calibration rows + ", + sum(out$key_type == "focal_length"), " fallback rows") +message("Wrote inst/extdata/camera_formats_excluded.csv: ", nrow(excluded), " keys (", + format(sum(excluded$frames, na.rm = TRUE), big.mark = ","), " frames)") +message("Coverage: ", + format(sum(shipped$frames, na.rm = TRUE), big.mark = ","), " of ", + format(nrow(frames), big.mark = ","), " digital frames resolvable by calibration") diff --git a/data-raw/make_testdata.R b/data-raw/make_testdata.R index 0bb575c..37ae237 100644 --- a/data-raw/make_testdata.R +++ b/data-raw/make_testdata.R @@ -177,3 +177,49 @@ message("dem.tif: ", paste(dim(dem_clip)[1:2], collapse = "x"), " cells at ", " m, ", round(file.size(dem_path) / 1024), " KB") message("\nDone. Test data in: ", outdir) + + +# --- Digital centroids: real frames over the same AOI (#32) ---------------- +# +# The 20 photos above are the 1968 film sample, so nothing in them can exercise the +# digital path. The AOI itself is not film-only though — it holds 181 `Digital - Colour` +# frames from two cameras with very different sensor shapes, which is what lets a test +# tell a correctly-shaped footprint from a square one: +# +# 121201_2011 Leica DMC II 230 87.1 x 79.2 mm aspect 1.10 +# 20814295_2018 UltraCam Eagle M3 105.8 x 68.0 mm aspect 1.56 +# +# Fetched through bcdata rather than a hand-built WFS URL: a `BBOX(SHAPE, ...)` CQL +# filter returns 0 features for this AOI — including for a control that certainly has +# frames — so the obvious query reports a false absence. + +message("Fetching digital centroids (network step) ...") +dig <- bcdata::bcdc_query_geodata("WHSE_IMAGERY_AND_BASE_MAPS.AIMG_PHOTO_CENTROIDS_SP") |> + bcdata::filter(bcdata::BBOX(local(test_bbox), crs = "EPSG:4326")) |> + dplyr::filter(MEDIA == "Digital - Colour") |> + bcdata::collect() + +stopifnot(nrow(dig) > 0) +names(dig) <- tolower(names(dig)) + +# Keep whole flight lines rather than a random sample: `fly_bearing()` needs consecutive +# frames on a roll, and a scattered sample leaves every frame bearing-less. +set.seed(42) +keep <- dig |> + dplyr::group_by(.data$film_roll) |> + dplyr::arrange(.data$frame_number, .by_group = TRUE) |> + dplyr::slice_head(n = 6) |> + dplyr::ungroup() + +keep <- keep |> + dplyr::select(dplyr::any_of(c( + "airp_id", "photo_year", "photo_date", "scale", "film_roll", "frame_number", + "media", "photo_tag", "nts_tile", "focal_length", "flying_height", + "ground_sample_distance", "thumbnail_image_url", "flight_log_url", + "camera_calibration_url", "patb_georef_url", "geometry" + ))) + +st_write(keep, file.path(outdir, "photo_centroids_digital.gpkg"), + delete_dsn = TRUE, quiet = TRUE) +message("photo_centroids_digital.gpkg: ", nrow(keep), " frames, cameras ", + paste(sort(unique(basename(keep$camera_calibration_url))), collapse = ", ")) diff --git a/inst/extdata/camera_formats.csv b/inst/extdata/camera_formats.csv new file mode 100644 index 0000000..3ec53ae --- /dev/null +++ b/inst/extdata/camera_formats.csv @@ -0,0 +1,20 @@ +"key","key_type","camera","report_serial","source_pdf","width_mm","height_mm","px_cross","px_along","pitch_um","focal_mm","mm_stated","n_cameras","width_spread_pct","extrapolated","note","retrieved" +"121201_2011","calib_file","DMC II",,"Manuf_calib_DMCii_027.pdf",87.0912,79.2064,15552,14144,5.6,92.0145,TRUE,1,0,FALSE,,"2026-08-30" +"20114172_2019","calib_file","UltraCam Falcon M2","UC-Fp-1-20114172-f70","QSI_Calibration_UC-Fp-1-20114172-f70_20190926.pdf",103.86,67.86,17310,11310,6,70.5,FALSE,1,0,FALSE,,"2026-08-30" +"20814295_2013","calib_file","UltraCam Eagle","UC-Eagle-1-20814295-f80","2013_20814295.pdf",104.052,68.016,20010,13080,5.2,79.8,TRUE,1,0,FALSE,,"2026-08-30" +"20814295_2014","calib_file","UltraCam Eagle","UC-Eagle-1-20814295-f80","UC-E-1-20814295-f80-Rev04.01_V01.pdf",104.052,68.016,20010,13080,5.2,79.8,TRUE,1,0,FALSE,,"2026-08-30" +"20814295_2016","calib_file","UltraCam Eagle","UC-Eagle-1-20814295-f80","UC-E-1-20814295-f80-Rev07_V07_Short.pdf",104.052,68.016,20010,13080,5.2,79.8,TRUE,1,0,FALSE,,"2026-08-30" +"20814295_2017","calib_file","UltraCam Eagle","UC-E-1-20814295-f80","2017_092EFLK102I_Camera_Calibration_Report.pdf",104.052,68.016,20010,13080,5.2,79.8,TRUE,1,0,FALSE,,"2026-08-30" +"20814295_2018","calib_file","UltraCam Eagle M3","UC-EpII-1-22814295-f80","2018_calib_UC-Eagle_1_22814295-f80.pdf",105.84,68.016,26460,17004,4,79.8,TRUE,1,0,FALSE,,"2026-08-30" +"20910461_2016","calib_file","UltraCamXp","UC-SXp-1-20910461","20910461_16.pdf",103.86,67.86,17310,11310,6,100.5,TRUE,1,0,FALSE,,"2026-08-30" +"40112365_2015","calib_file","UltraCamXp","UC-SXp-1-40112365","40112365_15.pdf",103.86,67.86,17310,11310,6,100.5,TRUE,1,0,FALSE,,"2026-08-30" +"50311261_2014","calib_file","UltraCam Eagle","UC-Eagle-1-50311261-f100","50311261_14.pdf",104.052,68.016,20010,13080,5.2,100.5,TRUE,1,0,FALSE,,"2026-08-30" +"70912643_2015","calib_file","UltraCamXp","UC-SX-1-70912643","70912643_15.pdf",103.896,67.824,14430,9420,7.2,100.5,TRUE,1,0,FALSE,,"2026-08-30" +"dmc100039_2006","calib_file","DMC","DMC01 - 0039","dmc100039_06.pdf",165.888,92.16,13824,7680,12,120,FALSE,1,0,FALSE,"catalogue records focal 127; the report says 120 - report preferred","2026-08-30" +"dmc327542_2017","calib_file","DMC III","DMC III 27542","CalibProtocol_DMCIII_27542.pdf",100.3392,56.9088,25728,14592,3.9,92,TRUE,1,0,FALSE,,"2026-08-30" +"dmc327550_2018","calib_file","DMC III","DMC III 27550","Calibration Certificate_27550.pdf",100.3392,56.9088,25728,14592,3.9,92,TRUE,1,0,FALSE,,"2026-08-30" +"100","focal_length","inferred from focal length",,,103.86,67.86,,,,100,,3,0.2,FALSE,"modal width among calibrated cameras at this catalogue focal","2026-08-30" +"120","focal_length","extrapolated from nearest focal",,,165.888,92.16,,,,120,,1,,TRUE,"no calibrated camera at focal 120; the 11,984 frames are 2011 and the only large-format 120 mm body in the record is the Z/I DMC, whose own report gives a virtual focal of exactly 120 mm - so its 165.888 mm format is used. Marked extrapolated: this rests on camera identity, not on the catalogue focal being close to 127","2026-08-30" +"127","focal_length","inferred from focal length",,,165.888,92.16,,,,127,,1,0,FALSE,"modal width among calibrated cameras at this catalogue focal","2026-08-30" +"80","focal_length","inferred from focal length",,,104.052,68.016,,,,80,,3,5.5,FALSE,"modal width among calibrated cameras at this catalogue focal","2026-08-30" +"92","focal_length","inferred from focal length",,,100.3392,56.9088,,,,92,,2,15.2,FALSE,"modal width among calibrated cameras at this catalogue focal","2026-08-30" diff --git a/inst/extdata/camera_formats_excluded.csv b/inst/extdata/camera_formats_excluded.csv new file mode 100644 index 0000000..69ee8d8 --- /dev/null +++ b/inst/extdata/camera_formats_excluded.csv @@ -0,0 +1,6 @@ +"key","key_type","reason","frames","retrieved" +"10210206_2015","calib_file","calibration present only as a scanned image, not machine-readable; visually read as UltraCam Eagle 20010x13080 @ 5.2um = 104.052x68.016mm f100.5",495,"2026-08-30" +"11937933_2009","calib_file","calibration present only as a scanned image, not machine-readable; visually read as AIC Pro (P65+) 8984x6732 @ 6.0um = 53.904x40.392mm f60.68",443,"2026-08-30" +"12335326_2017","calib_file","calibration present only as a scanned image, not machine-readable; visually read as PhaseOne IXU-RS-1000 11608x8708 @ 4.6um = 53.4x40.1mm f51.56",7140,"2026-08-30" +"72914123_2019","calib_file","catalogue GSD/FLYING_HEIGHT contradict the report: implied ground elevation -418 m",1545,"2026-08-30" +"83","focal_length","no fallback row: no calibrated camera at focal 83, and extrapolating from a different focal is not warranted without knowing the camera - focal distance does not predict the error",,"2026-08-30" diff --git a/inst/extdata/camera_formats_manifest.csv b/inst/extdata/camera_formats_manifest.csv new file mode 100644 index 0000000..a21b473 --- /dev/null +++ b/inst/extdata/camera_formats_manifest.csv @@ -0,0 +1,19 @@ +"key","retrieved" +"10210206_2015","2026-08-30" +"11937933_2009","2026-08-30" +"121201_2011","2026-08-30" +"12335326_2017","2026-08-30" +"20114172_2019","2026-08-30" +"20814295_2013","2026-08-30" +"20814295_2014","2026-08-30" +"20814295_2016","2026-08-30" +"20814295_2017","2026-08-30" +"20814295_2018","2026-08-30" +"20910461_2016","2026-08-30" +"40112365_2015","2026-08-30" +"50311261_2014","2026-08-30" +"70912643_2015","2026-08-30" +"72914123_2019","2026-08-30" +"dmc100039_2006","2026-08-30" +"dmc327542_2017","2026-08-30" +"dmc327550_2018","2026-08-30" diff --git a/inst/testdata/photo_centroids_digital.gpkg b/inst/testdata/photo_centroids_digital.gpkg new file mode 100644 index 0000000..08d4d38 Binary files /dev/null and b/inst/testdata/photo_centroids_digital.gpkg differ diff --git a/man/fly_footprint.Rd b/man/fly_footprint.Rd index 9d236b5..397f50d 100644 --- a/man/fly_footprint.Rd +++ b/man/fly_footprint.Rd @@ -16,8 +16,10 @@ recording format per frame when present.} \code{media} column. It never sizes a digital frame — see \code{format_size}.} \item{format_size}{Named numeric vector of recording-format widths in inches, -keyed by \code{media} value, merged over the shipped film defaults. Supply this -to size frames whose format \code{fly} does not know — see Details.} +keyed by \code{media} value, merged over the shipped film defaults. Frames it names are +sized from the reported \code{scale}, as film is, and it takes precedence over the +shipped camera table — it is the escape hatch for a camera \code{fly} does not know. +See Details.} \item{dem}{Optional elevation raster used to size each frame from its true height above ground rather than the reported scale. A \code{terra::SpatRaster}, @@ -64,15 +66,47 @@ records the outcome: \describe{ \item{the \code{media} value}{format resolved from the format table} +\item{\code{"inferred_format"}}{digital frame with no calibration, sized from a +format inferred from its \code{focal_length}} \item{\code{"assumed_default"}}{no \code{media} column; \code{negative_size} applied} \item{\code{"unknown_format"}}{\code{media} present but unknown; empty geometry} } - -Shipped defaults cover film only. Digital frames resolve to -\code{"unknown_format"} rather than an invented number, because the sensor width -they would need is not in the centroid metadata — and neither is the pixel -count that would let \code{ground_sample_distance} stand in for it. Supply -\code{format_size} if you know the camera: +} +\section{Digital frames}{ + + +A digital frame has no negative, and the catalogue mixes film and digital in one +layer — 223,667 of 1,670,471 frames province-wide are \code{Digital - Colour}. Sensor +dimensions are not in the centroid metadata, but they are recoverable from the +calibration report each frame links to through \code{camera_calibration_url}, and \code{fly} +ships them (\code{inst/extdata/camera_formats.csv}, built by +\code{data-raw/make_camera_formats.R}). + +Digital frames are sized as \verb{pixel count x ground_sample_distance}, which needs +neither \code{scale} nor a DEM. \strong{\code{scale} is never used for a digital frame \code{fly} sized +itself.} That field is not the true image scale for digital: measured against +terrain on 40 UltraCam Eagle frames it gives 34\% of true width, because it is a +derived nominal figure — the pixel pitch it implies is about 12.5 um for every +camera regardless of model, against real pitches of 3.9 to 12 um. +\code{ground_sample_distance} is in centimetres. + +Where a frame carries no \code{camera_calibration_url} — about a fifth of digital frames — +the format is inferred from \code{focal_length} and \code{footprint_basis} records +\code{"inferred_format"}. Sensor width spreads only 1-3\% at a given focal length, but +pixel count spreads 32-83\%, so an inferred frame can only be sized through a DEM +(\verb{width x height above ground / focal length}) and never from its GSD. + +\code{width_source} names the calibration file or fallback rule per row, so every +footprint traces back to a source. Calibrations that could not be corroborated are +listed in \code{inst/extdata/camera_formats_excluded.csv} with the reason, and frames +naming one are refused rather than inferred. + +\strong{Digital footprints are not square} — sensors run from 1.10:1 (Leica DMC II) to +1.80:1 (Intergraph DMC) — so they are rotated onto the flight line using +\code{\link[=fly_bearing]{fly_bearing()}}. Where no bearing can be computed the rectangle stays axis-aligned +and \code{width_source} says so. Film stays square and is unaffected. + +Supply \code{format_size} to size a frame \code{fly} cannot, or to override it: \if{html}{\out{
}}\preformatted{fly_footprint(photos, format_size = c("Digital - Colour" = 3.54)) }\if{html}{\out{
}} @@ -86,6 +120,7 @@ the bundled test data. Note the field is \code{SCALE}, not \code{PHOTO_SCALE} latter returns all \code{NULL}, which reads as missing data rather than a wrong field name. } + \section{Terrain}{ @@ -118,6 +153,9 @@ at most 0.5\% against the correction's own 14\% — and a third moves it by \describe{ \item{\code{"nominal_scale"}}{sized from the reported scale (no \code{dem}, or a fallback — see below)} +\item{\code{"gsd_scaled"}}{digital frame sized from its pixel count and ground +sample distance; used neither the reported scale nor a DEM, so \code{height_agl} +and \code{dem_coverage} are \code{NA}} \item{\code{"dem_agl"}}{sized from height above ground} \item{\code{"no_dem_coverage"}}{\code{dem} supplied but does not cover the frame} \item{\code{NA}}{no footprint to place — see \code{footprint_basis}} diff --git a/planning/archive/2026-08-issue-32-digital-sensor-formats/README.md b/planning/archive/2026-08-issue-32-digital-sensor-formats/README.md new file mode 100644 index 0000000..125f977 --- /dev/null +++ b/planning/archive/2026-08-issue-32-digital-sensor-formats/README.md @@ -0,0 +1,43 @@ +# Issue #32 — sensor widths for digital frames + +`fly_footprint()` marked every digital frame `unknown_format` with an empty geometry. +#30 made that gap honest without closing it: 223,667 of 1,670,471 frames province-wide, +so the package was quietly film-only for anything after ~2010. + +Digital frames are now sized from their camera's sensor, read out of the calibration +reports the catalogue links to through `camera_calibration_url` and shipped as +`inst/extdata/camera_formats.csv` — 14 calibrations covering 169,688 frames, plus +focal-length fallback rows for frames carrying no calibration. + +Two findings during exploration reshaped the work, both confirmed with the user before +implementation: + +- **The catalogue's `SCALE` is not the true image scale for a digital frame.** Sizing + from it gives 34% of true width; it is a derived nominal figure implying a ~12.5 um + pixel pitch for every camera against real pitches of 3.9-12 um. Digital frames are + sized as `pixel count x ground_sample_distance` instead, needing neither `scale` nor + a DEM. +- **Sensors are not square** (1.10:1 to 1.80:1, widths 87.1-165.9 mm), so footprints + became rectangles rotated onto the flight line via `fly_bearing()`. Film unchanged. + +The numbers are parsed, never typed, and gated on five checks. Two caught real problems +before shipping: a camera whose catalogue metadata contradicts its own report (withheld), +and — via an independent second reading — that the focal-83 frames are a 53.9 mm +medium-format body, so the fallback that would have inferred a 104 mm UltraCam for them +is refused. That was a 93% error. + +**What this issue is worth remembering for.** Three defects were found only by restoring +them: the rotation matrix was backwards, the orientation test was decoration (`%% 180` +discarded the sign it existed to catch), and the 2013 UltraCam report renders `20010` as +`2001O` and drops the micron sign, so `Pixel Size 5.200 m` reads as metres. + +And code-check rounds 2 and 3 each found a defect **inside the previous round's fix**, +both with the same root cause: branching on `half_cross`, which is NA for every +camera-table row by construction. Round 2 — the DEM route unreachable for exactly the +frames it exists to serve. Round 3 — rotation never applied to those same frames. Asking +for the mechanism rather than more instances is what ended it: rotation is now decided +from the format's aspect ratio, known before any sizing route runs, and a structural +invariant sweep over 12 input shapes guards the reporting columns that had drifted. + +Closed by PR (see `Fixes #32`). Suite 1162 pass / 0 fail / 0 warn. +Follow-up: #38 (georeferencing digital frames). Related: #30, #9. diff --git a/planning/archive/2026-08-issue-32-digital-sensor-formats/findings.md b/planning/archive/2026-08-issue-32-digital-sensor-formats/findings.md new file mode 100644 index 0000000..5c59582 --- /dev/null +++ b/planning/archive/2026-08-issue-32-digital-sensor-formats/findings.md @@ -0,0 +1,129 @@ +# Findings — sensor widths for digital frames (#32) + +All numbers below measured 2026-08-29/30 against the live BC Data Catalogue WFS +(`WHSE_IMAGERY_AND_BASE_MAPS.AIMG_PHOTO_CENTROIDS_SP`) and the calibration reports it +links to. + +## Population + +| | frames | share | +|---|---|---| +| whole catalogue | 1,670,471 | | +| `Digital - Colour` (the only digital media value) | 223,667 | 13.4% | +| digital with a `CAMERA_CALIBRATION_URL` | 179,311 | 80.2% of digital | +| digital resolvable to a sensor size | 171,233 | 76.6% of digital | + +18 distinct calibration files, 14 distinct serials. Not the two cameras the issue +assumed, but a tractable set rather than a long tail. + +## `SCALE` is not the true image scale for digital frames + +The finding that most changes the work. Measured on 40 UltraCam Eagle frames against +MRDEM-30 terrain: + +| route | width | +|---|---| +| `pixel_count x GSD` | 6003 m | +| `sensor_width x (flying_height - terrain) / focal` | 6070 m | +| `sensor_width x SCALE` (today's arithmetic) | **2081 m** | + +The first two agree to ~1% and reproduce the catalogue's own `GROUND_SAMPLE_DISTANCE` +to a 1.011 median ratio — an independent check, since GSD was not used to derive them. +The implied pixel pitch from `GSD x 10000 / scale_denominator` is ~12.5 um for almost +every camera regardless of model, against real pitches of 3.9-12 um, which is the tell +that `SCALE` is a derived nominal figure rather than a measurement. + +Consequence: shipping a sensor width while keeping the scale path would draw digital +footprints at a third of true size — still drawing, still overlapping, still producing +coverage percentages. Worse than #30's refusal. + +## Sensors are not square + +| camera | w x h mm | aspect | +|---|---|---| +| Leica DMC II 230 | 87.0912 x 79.2064 | 1.10 | +| UltraCam (all models) | ~103.9-105.8 x ~67.9-68.0 | 1.53-1.56 | +| Leica DMC III | 100.3392 x 56.9088 | 1.76 | +| Intergraph DMC | 165.888 x 92.160 | 1.80 | + +`fly_rectangles()` builds squares. A square DMC III footprint is 76% too deep. Width +alone also spans 87.1-165.9 mm, so a single invented "digital" default would have been +wrong by up to 90%. + +Orientation matters once footprints are non-square: the wide dimension is cross-track +(the Vexcel reports label it `cross track` explicitly), so the rectangle must be rotated +to the flight line. `fly_bearing()` already computes that. + +## Keying + +`media` is one value across all 14 cameras. `focal_length` is ambiguous — measured +against ground truth on frames that do have a calibration: + +| catalogue focal | distinct widths | spread | +|---|---|---| +| 70, 79, 90, 127 | 1 each | exact | +| 100 | 4 | 1.9% | +| 80 | 3 | 5.5% | +| 92 | 2 (87.09 vs 100.34) | 15.2% | + +Pixel count is far worse at the same keys — 32% spread at focal 80, 83% at focal 100 — +so a focal-keyed frame can be sized by width x AGL/focal but **not** by px x GSD. + +Serial 20814295 is genuinely two cameras: UltraCam Eagle (20010 px @ 5.2 um) through +2017, Eagle Prime II (26460 px @ 4.0 um) in 2018 — and the 2018 report numbers itself +22814295, so the catalogue's URL basename is not a reliable serial. Only the full +calibration-file key distinguishes them. + +## Unresolvable reports + +| file | why | frames | +|---|---|---| +| `11937933_2009` (Rollei P65) | image-only PDF, 0 extractable characters | 443 | +| `12335326_2017` (PhaseOne) | image-only PDF, 0 extractable characters | 7,140 | +| `10210206_2015` | aerial-triangulation report, not a calibration certificate | 495 | + +## QA results (pre-implementation run) + +| check | result | +|---|---| +| B: `px x pitch == stated image size`, both axes | **15/15**, max rel. err 1.8e-16 | +| C: report focal vs catalogue focal | 14/15 — `dmc100039_2006` reports 120, catalogue 127 | +| F: implied ground elevation plausible | 13/14 — `72914123_2019` implies -414 m | +| G: two sizing routes agree | 1.1% median deviation over 40 frames | + +Every report states pixel count, pixel size **and** image size independently, so B is a +real constraint on all rows rather than a restatement. + +F's sensitivity, measured rather than assumed: doubling every pitch makes only 8 of 14 +cameras implausible. It is a gross-error net; B is what gives precision. + +### Two rows the QA flagged + +- `dmc100039_2006` (4,052 frames) — report gives virtual focal 120 mm, catalogue records + 127. Ship with both recorded; GSD is 0 for all these frames so only the DEM route + applies and focal enters linearly. +- `72914123_2019` (1,545 frames) — catalogue `GSD` 12 cm, `SCALE` 9600 and + `FLYING_HEIGHT` ~2600 m are mutually inconsistent under the report's own 4.0 um / + 100.5 mm: implied AGL 3015 m exceeds the aircraft's height. GSD 8 cm or pitch 6.0 um + would reconcile it; neither is supported. **Withheld.** + +## Fixture + +The bundled Houston AOI holds **181 digital frames** over the DEM that already ships: +76 `121201_2011` (Leica DMC II, aspect 1.10) and 105 `20814295_2018` (UltraCam Eagle +Prime II, aspect 1.56). Two very different aspect ratios, so the fixture can reach the +non-square branch. + +This corrects `CLAUDE.md`, which says the AOI has no digital frame and that "digital +coverage cannot come from there". True of the 20 *sampled* photos (`photo_year == 1968`), +not of the AOI. + +## Errors Encountered + +| Error | Resolution | +|-------|------------| +| `BBOX(SHAPE,...)` CQL filter returned 0 for a **positive control** (an AOI known to hold 1,584 frames) | Not a data finding — a broken probe. Query bboxes through `bcdata::filter(BBOX(...))`, and always run a positive control before reporting an absence | +| `sed 1d f1 f2 f3` inside `find -exec ... +` strips only the first file's header | 23 stray header rows entered a 223k-row analysis silently. Loop per file: `for f in ...; do sed 1d "$f"; done` | +| `sed -n '/X/,$d' file` prints nothing | `-n` suppresses auto-print, so a delete-to-end script emits an empty file. Drop `-n` | +| WFS caps `GetFeature` at 10,000 features with no error | Page with `startIndex`, and reconcile the total against `resultType=hits` | +| MRDEM `/vsicurl` extraction over points scattered across BC exceeded a 10-minute timeout | For QA that only needs plausibility, invert for implied terrain from `FLYING_HEIGHT` and `GSD` instead — no DEM needed | diff --git a/planning/archive/2026-08-issue-32-digital-sensor-formats/progress.md b/planning/archive/2026-08-issue-32-digital-sensor-formats/progress.md new file mode 100644 index 0000000..b5e9541 --- /dev/null +++ b/planning/archive/2026-08-issue-32-digital-sensor-formats/progress.md @@ -0,0 +1,32 @@ +# Progress — sensor widths for digital frames (#32) + +## Session 2026-08-30 + +- Plan-mode exploration against the live catalogue; phases approved by user +- Established the provincial digital population (223,667 frames, 18 calibration files, + 14 serials) — the issue's "unknown at filing" +- Found that `SCALE` is not the true image scale for digital frames, and that sensors + are not square; both confirmed with the user and folded into the plan +- Ran the numeric QA before writing any code: B 15/15, C 14/15, F 13/14, G 1.1% median. + Two rows flagged and dispositioned +- Created branch `32-establish-sensor-widths-for-digital-fram` off main +- Scaffolded PWF baseline with approved phases +- Phase 1: generator + `camera_formats.csv` (14 calibration + 5 fallback rows) +- Phase 2: QA B/C/D/E/F all run. E (independent double entry) agreed exactly on all 14 + shipped rows, and additionally recovered three scanned reports plus the AIC Pro's + 53.9 mm width — which is why the focal-83 fallback is now refused +- Phase 3-4: non-square rotated footprints, resolver, `width_source`, `gsd_scaled` +- Phase 5: real digital testdata, check G, docs, NEWS, `CLAUDE.md` correction +- Plan review (10 findings) and code-check round 1 (12 findings) both folded in +- Suite 426 pass / 0 fail; lint +3 over baseline, all confirmed installed-vs-source + artifacts for the new internal functions +- Code-check rounds 2 and 3 both found a defect inside the previous round's fix, and + both had the SAME root cause: branching on `half_cross`, which is NA for every + camera-table row by construction. Round 2 - the DEM route unreachable for those rows; + round 3 - rotation never applied to them. Rotation is now decided from the format's + aspect ratio, which is known before any sizing route runs +- Added a structural invariant sweep (`test-fly_footprint_invariants.R`) over 12 input + shapes, since three separate conditions had drifted from the same fact +- Suite 1162 pass / 0 fail / 0 warn; lint baseline+3, all confirmed artifacts +- Filed #38 for georeferencing digital frames, deferred rather than guessed +- Next: PR diff --git a/planning/archive/2026-08-issue-32-digital-sensor-formats/review-32.md b/planning/archive/2026-08-issue-32-digital-sensor-formats/review-32.md new file mode 100644 index 0000000..5ebf0f7 --- /dev/null +++ b/planning/archive/2026-08-issue-32-digital-sensor-formats/review-32.md @@ -0,0 +1,73 @@ +# Plan review — #32 (Plan agent, 2026-08-30) + +Spawned after the PWF baseline, per `planning.md`. Read `fly_footprint.R` in full plus +every downstream consumer. Findings verified against the code before acting; the +verification result is recorded beside each. + +## Blockers + +**1. Sizing precedence is the reverse of what the code does.** `fly_footprint.R:417-427,479` +overwrites `half_side[corrected]` unconditionally whenever a DEM yields a usable AGL. A +digital frame sized by `px x GSD` arrives with `sized == TRUE`, so passing `dem` silently +switches its route — and every downstream consumer forwards `dem`. **Confirmed.** Fix: +gate `corrected` on `!gsd_sized`, and test that a digital frame is identical with and +without `dem`. + +**2. `footprint_terrain` makes a false claim for `px x GSD` frames.** Line 391 labels any +non-NA half-side `"nominal_scale"`, documented as "sized from the reported scale" — the +exact claim the work promises never to make for a digital frame. **Confirmed.** Needs a +fifth value, plus honest `height_agl` (NA) and `dem_coverage` (NA) for that route. + +**3. Rotation must be applied inside the two-pass DEM sampler, not just to the output.** +Lines 410-411 build both sampling windows through `fly_rectangles()`. Rotating only the +returned geometry makes `dem_coverage` describe an unrotated rectangle — the +returned/measured mismatch `test-fly_footprint.R:614-649` exists to catch. **Confirmed.** + +**4. `GROUND_SAMPLE_DISTANCE` units.** Reviewer flagged that findings.md is ambiguous +between metres and centimetres and that Phase 4 writes `px x GSD` with no conversion. +**Verified: the field is CENTIMETRES.** UltraCam Eagle, GSD 30, 20010 px -> 20010 x 0.30 m += 6003 m, which matches `104.052 mm x (AGL/focal)` = 6070 m. In metres it would be 100x. +`mixed_media_fixture()` sets `0.10`, a metres assumption, and the shipped GeoPackage +column is `int` and all-NA — so nothing currently exercises it. Convert explicitly and +assert the unit. + +## Gaps + +**5. `GSD == 0` builds a degenerate polygon, not an empty one.** `fly_rectangles()` guards +`is.finite(w)`; zero is finite, so a five-identical-vertex POLYGON passes `st_is_empty()` +and is invisible to `fly_warn_unsized()`. **Confirmed, and reachable:** `dmc100039_2006` +has GSD 0 on all 4,052 frames. Add `w > 0`. + +**6. The focal-length fallback silently resolves `mixed_media_fixture()` row 4** (focal +100, no calibration URL), breaking four existing assertions at `test-fly_footprint.R:52, +61, 70, 85`. **Confirmed.** Gate the fallback on digital `media`, and update those tests +deliberately rather than discovering them red. + +**7. `fly_georef()` applies its own bearing rotation to footprint corners** +(`fly_georef.R:199-244, 300-304`), calibrated against axis-aligned squares. A pre-rotated +non-square footprint would be double-corrected, and a 90-degree error maps landscape onto +portrait. **Confirmed.** Needs an explicit decision, not silence. + +**8. `fly_bearing()` stops rather than returning NA** when `film_roll`/`frame_number` are +absent (`fly_bearing.R:31-34`), and takes the cross-track jump as the bearing at each +flight-line turn. **Confirmed.** Guard in `fly_footprint()`; document the turn case. + +## Acceptance + +**9. The proposed rotation tests cannot fail the way that matters.** Area is invariant +under *every* rotation and a 90-degree bbox swap tests the mechanism, not the angle — so a +long-axis/short-axis transposition and a sign-flipped azimuth both pass. **Agreed.** +Discriminating assertion: the long axis must be perpendicular to `bearing`, and +consecutive same-roll frames must show along-track overlap (~60%) not side-lap (~30%). + +**10. "`SCALE` is never used for a digital frame" contradicts keeping `format_size`.** +`test-fly_footprint.R:89-98` pins `format_size` sizing a digital frame from `SCALE`. +**Confirmed.** Restate as "never used for a frame `fly` sized itself". Also: `format_size` +is a named numeric vector and cannot carry two dimensions per key without a type change. + +## Non-finding (recorded so it is not re-litigated) + +Rotation does **not** break the bounded-allocation invariant. `fly_dem_grid()` sizes to one +footprint's extent; a rotated `w x h` bbox peaks at `(w+h)^2/2`, about 2.2x unrotated for +aspect 1.8 — bounded and per-frame, nowhere near the union-grid defect. Centre-to-corner +distance is invariant under rotation, so the existing DEM buffer rule already covers it. diff --git a/planning/archive/2026-08-issue-32-digital-sensor-formats/review-round1.md b/planning/archive/2026-08-issue-32-digital-sensor-formats/review-round1.md new file mode 100644 index 0000000..b1f5667 --- /dev/null +++ b/planning/archive/2026-08-issue-32-digital-sensor-formats/review-round1.md @@ -0,0 +1,195 @@ +# Review round 1 — fly#32 camera format table + +Reviewed: staged diff plus the unstaged `fly_camera_format()` resolver in +`R/fly_camera_format.R` (working tree moved during the review; findings below were +re-verified against the tree state at the end, and the ones the author fixed mid-review +— the failing manifest drift guard, the generator/CSV divergence over focal 83 — are +**not** listed). + +Verification used: live BC WFS query of all 223,667 `MEDIA LIKE 'Digital%'` frames +(`numberMatched` reconciled, 0 duplicate FIDs), plus R probes for each claimed R +semantic. The parsed numbers themselves check out against published camera +specifications — UltraCam Eagle 20010x13080 @ 5.2 um, Eagle M3 26460x17004 @ 4.0, +UltraCamXp 17310x11310 @ 6.0, UltraCamX 14430x9420 @ 7.2, DMC II 230 15552x14144 @ 5.6, +DMC III 25728x14592 @ 3.9, Z/I DMC 13824x7680 @ 12 f120 — and the cross-track/along-track +assignment is correct in every parser. No transcription error found in the shipped table. + +## Findings + +- **[bug]** `data-raw/make_camera_formats.R:483-500` — both exclusion-block + `data.frame()` calls abort the script when their key vector is empty. `paste0()` + returns length 1 for a zero-length argument (and `ifelse()` on a `logical(0)` test + returns `logical(0)`, which `paste0` also absorbs), so `reason` is length 1 while + `key` is length 0. Verified both: + ``` + unparsed empty -> ERROR: arguments imply differing number of rows: 0, 1 + nothing withheld -> ERROR: arguments imply differing number of rows: 0, 1 + ``` + The second one fires in the **healthy** case: the moment check F rejects nothing, + `parsed$key[!parsed$ship]` is `character(0)` and the run dies at line 496 — after + ~25 MB of downloads, all parsing and all QA, with `camera_formats.csv` already + written at line 471 but the excluded file and manifest not. That leaves the three + artifacts inconsistent on disk, which is the state the manifest guard exists to + detect. Guard each block with `if (length(...))`, as line 504 already does for + `refused`. + +- **[bug]** `data-raw/make_camera_formats.R:118-121` — the extract is guarded on + directory *existence*, and `fs::dir_create(d)` runs before `utils::unzip()`. `unzip` + does not error on a bad archive; verified: + ``` + warning: error 1 in extracting from zip file + unzip result: no error dir still exists: TRUE files: 0 + ``` + So an interrupted or corrupt extract leaves an empty directory that every later run + skips. `parse_one()` then finds no PDFs, the key lands in `unparsed`, and it is + written to `camera_formats_excluded.csv` with the reason **"calibration present only + as a scanned image, not machine-readable"** — a specific factual claim that is now + false and permanent until someone clears the cache by hand. This is the same + guard-on-existence mistake the author explicitly fixed one block earlier for the zip + (line 111, "Guard on non-empty, not existence"); the fix was not carried to the + extract. Extract to a `.part` directory and rename, or gate on the directory + containing at least one `.pdf`. + +- **[bug]** `R/fly_camera_format.R:117-131` — the resolver never consults + `camera_formats_excluded.csv`, so a frame whose calibration key was *deliberately + refused* falls straight through to focal-length inference. Verified against real + catalogue focal lengths: + ``` + key focal width_mm width_source inferred + 72914123_2019 100 103.860 focal_length=100 TRUE (1,545 frames) + 10210206_2015 100 103.860 focal_length=100 TRUE ( 495 frames) + 11937933_2009 83 NA FALSE ( 443 frames) + 12335326_2017 53 NA FALSE (7,140 frames) + ``` + The two dangerous ones (AIC Pro 53.904 mm, PhaseOne 53.4 mm) come back NA **only + because no fallback row happens to exist at focal 83 or 53** — safety by coincidence + of the current table contents, not by construction. `FALLBACK_REFUSE` enforces this + in the generator; nothing enforces it in the consumer. Any future table that gains a + row at those focals silently sizes 7,583 frames at ~1.95x their true width (3.8x + area). Match `key` against the excluded set and refuse, rather than relying on the + fallback set staying sparse. + +- **[bug]** `data-raw/make_camera_formats.R:401-408` — `FALLBACK_REFUSE` encodes the + one instance that was measured, not the property it stands for. Two consequences, + both live: + + 1. **Focal 53 is not listed.** `12335326_2017` is a 53.4 mm PhaseOne with 7,140 + frames. It has no no-URL frames today, so `need` excludes 53 and no row is built. + The moment BC publishes one focal-53 frame without a calibration URL, the + generator builds a row extrapolated from the nearest calibrated focal (70, the + UltraCam Falcon at **103.86 mm**) — 1.95x too wide, 3.8x too much ground area, + shipped with `width_spread_pct = 0`. + 2. **Extrapolation distance does not predict the error, so no distance bound fixes + this.** The refused focal-83 row extrapolated from focal 80 — a 3.6% focal gap + producing a **93%** width error. The shipped row `"120"` extrapolates from focal + 127 — a *larger* 5.5% gap — and is fine. Nothing in the generator's data separates + them; only camera identity does. Row `"120"` covers 11,984 frames (5.4% of the + whole digital catalogue) on exactly the warrant just refused as unsound. + + Also note the refusal is silent if it stops applying: if the catalogue ever populates + a URL for the focal-83 frames, `refused` becomes empty, the excluded row disappears, + and nothing reports that a declared safety refusal no longer fires. Assert that every + `FALLBACK_REFUSE` name was actually matched. + +- **[fragile]** `data-raw/make_camera_formats.R:428,439` — `width_spread_pct` and + `n_cameras` are a proxy that does not guard what they are read as. They measure + dispersion among the *source* cameras, so a single-source row is structurally `0` + however bad the inference is. The two extrapolated rows carry the most confident + labels in the file: + ``` + "120" extrapolated from nearest focal 165.888 n_cameras 1 width_spread_pct 0 (11,984 frames) + "127" inferred from focal length 165.888 n_cameras 1 width_spread_pct 0 ( 845 frames) + ``` + 12,829 frames are labelled zero-spread on a width derived from one camera at a + different catalogue focal. The refused focal-83 row would have shipped + `width_spread_pct = 5.5` while being 93% wrong. An extrapolated row should not be + able to report 0. + + (Row `"120"` is *probably* right — a 120 mm digital frame is very likely the Z/I DMC, + which is 165.888 mm. But the generator did not reason that way: it reached it via + "nearest catalogue focal is 127", using the same catalogue focal field that the calib + row for that very camera declares wrong — `"catalogue records focal 127; the report + says 120 - report preferred"`. Right by luck through a field the script distrusts.) + +- **[fragile]** `data-raw/make_camera_formats.R:246` + QA sections C/F — `focal_mm` is + the only shipped number with no gating bound. `parse_dmc` computes + `num(sub(".*\\[m\\]", "", fl)) * 1000`; the unit anchor is the strict literal `[m]`, + which does not match `[mm]`. When the sub misses, it is a no-op and `num()` — which + strips non-digits and concatenates rather than failing — returns a finite number from + the whole line, which then passes the `is.finite` gate at line 287 and is multiplied + by 1000. Nothing catches it: check D bounds width, aspect and pitch but not focal; + check C computes `focal_disagrees` and only prints it, never appending to `fail`; and + check F is explicitly skipped where unmeasurable — + `terrain_ok <- is.na(terrain_implied) | ...` — so a key whose frames carry no + `GROUND_SAMPLE_DISTANCE`/`FLYING_HEIGHT` has **no** gating check on its focal length + at all. Add a plausibility band on `focal_mm` the way D does for width. + +- **[fragile]** `data-raw/make_camera_formats.R:342-345` and + `tests/testthat/test-camera_formats.R:41` — check D's lower bound of 80 mm fails + toward abort on legitimate input. The catalogue demonstrably contains medium-format + bodies at 53.4 mm and 53.9 mm (7,583 frames, per the author's own exclusion notes). + If either calibration ever becomes machine-readable and parses correctly, `bad_d` + fires and `stop("QA failed, CSV not written")` discards the entire run including the + 14 good rows — and the test asserting `all(tbl$width_mm > 80)` fails on a correct + table. The bound is tuned to large-format aerial cameras and the record is not + exclusively large-format. + +- **[fragile]** `data-raw/make_camera_formats.R:335` — `if (any(parsed$focal_disagrees))` + errors with "missing value where TRUE/FALSE needed" when a key's catalogue focal is + NA and no other row disagrees (`any()` returns NA, not FALSE). `focal_catalogue` is + NA whenever a key's frames all carry a missing `FOCAL_LENGTH`. Same NA then flows + into `ifelse(shipped$focal_disagrees, ...)` at line 461, producing an NA note. + +- **[fragile]** `data-raw/make_camera_formats.R:426-427` — `w <- as.numeric(names(tw)[1])` + round-trips a double through `as.character`, then `at$width_mm == w` compares it for + exact equality. That round trip is not lossless for several widths this pipeline + produces; verified: + ``` + 25728 x 3.9 -> "100.3392" roundtrip_ok = FALSE + 14592 x 3.9 -> "56.9088" roundtrip_ok = FALSE + 14144 x 5.6 -> "79.2064" roundtrip_ok = FALSE + ``` + It does not bite today only because the DMC III's 100.3392 is a *stated* value read + from the report, not the `px * pitch` product. A derived width of the same shape + yields `h = NA` and a fallback row with a missing height. Compare on the index + (`which.max(table(...))` against the original vector) rather than on a reparsed value. + +- **[fragile]** `data-raw/make_camera_formats.R:108,117` — the zip and pdf caches are + keyed on `basename(u)`, but the identity of a calibration is the full URL. Two URLs + in different directories sharing a basename would collapse: the second finds a + non-empty cached zip, skips the download, and silently inherits the first camera's + dimensions. Verified 18 distinct URLs / 18 distinct basenames today, so this is + latent, not live — but `frames$key` is built from the basename too, so the two frames + sets would already have merged before the cache was consulted. + +- **[fragile]** `R/fly_camera_format.R:77,104` — `inferred = FALSE` is emitted both for + an exact calibration match and for a row that resolved to nothing (see the probe in + the third finding: `11937933_2009` and `12335326_2017` come back `width_mm = NA, + inferred = FALSE`). A consumer filtering `!inferred` to get trustworthy rows also + picks up every unresolved frame. `width_mm` is NA there so a careful caller is safe, + but `!inferred` is the natural filter to write. + +- **[fragile]** `data-raw/make_camera_formats.R:511` — `match(excluded$key, by_key$key)` + pools two key namespaces. `by_key` holds calibration keys only, so the focal-keyed + refusal row gets `frames = NA` where the answer is 47 (no-URL frames at focal 83). + The summary at line 531-532 then understates the excluded frame count. Reporting + only, no effect on shipped numbers. + +## Checked and clean + +- `nzchar(NA)` ordering is correct at `make_camera_formats.R:92` (`is.na(...) | !nzchar(...)`) + and `fly_camera_format.R:110` (`!is.na(u) & nzchar(u)`). +- Empty CSV fields read back as `""` not NA for character columns, so the + `!is.na(x) & nzchar(x)` assertions in the tests are the right form and do bite. +- `fly_camera_cache[[file]]` on a missing name returns NULL (does not error), so the + cache short-circuit is sound; the cache key covers everything the read depends on. +- `stopifnot(nrow(d) == n, !anyDuplicated(d$FID))` is not vacuous — the WFS CSV does + carry `FID` regardless of `propertyName`, confirmed over all 223,667 rows, 0 duplicates. +- No roxygen blocks in the new R file, NAMESPACE unchanged at 9 exports — no + `@export` rebinding. +- `on.exit()` at `wfs_page:56` is inside a function, so it fires. +- `download.file` truncation is correctly handled via the `.part` + non-empty check. +- Check B is correctly restricted to `mm_stated` rows in both the generator and the + test, and the test asserts its premise (`expect_gt(nrow(b), 10)`) rather than passing + on an empty set. +- All six tests in `test-camera_formats.R` pass at the current tree state. diff --git a/planning/archive/2026-08-issue-32-digital-sensor-formats/review-round2.md b/planning/archive/2026-08-issue-32-digital-sensor-formats/review-round2.md new file mode 100644 index 0000000..7ddd6dc --- /dev/null +++ b/planning/archive/2026-08-issue-32-digital-sensor-formats/review-round2.md @@ -0,0 +1,180 @@ +# Review — round 2 (fly#32, branch `32-establish-sensor-widths-for-digital-fram`) + +Scope: staged diff (~1900 lines). Every claim below was reproduced by running the +code, not by reading it. Suite state at time of review: `FAIL 0 | PASS 426`. + +## Findings + +- **[bug]** `R/fly_footprint.R:564` (`sized <- !is.na(half_cross)`) with `:581` + (`corrected <- sized & !by_gsd & ...`) — **no digital frame can ever be + DEM-corrected, so every table-resolved frame that is not on the GSD route gets an + empty geometry even when a DEM is supplied.** + + `half_cross` at line 479 is `width_in * scale_num * 0.0254 / 2`, and `width_in` is + `NA` for every row resolved from the camera table (that is what `from_table` at + :458 means). So `sized` is FALSE for all of them. `corrected` requires `sized`, and + `corrected` additionally excludes `by_gsd`. The two conditions are disjoint: + `sized` is only ever TRUE for a digital row *because* `by_gsd` set `half_cross` at + :497 — and `!by_gsd` then removes exactly those. Compounding it, the first-pass + sampling rectangle at :556 is built from the same NA `half_cross`, so it is an + empty polygon and `first$elev` is NA regardless. + + Reproduced (digital fixture with `camera_calibration_url` stripped, so rows take + the focal-92 fallback; synthetic 100 m DEM covering all frames): + + ``` + footprint_basis: inferred_format 24 + footprint_terrain: 24 + empty geometries: 24 of 24 + areas (km2): 0 0 0 0 0 0 + ``` + + This contradicts the documentation added in this same diff + (`R/fly_footprint.R:252-253`): *"pixel count spreads 32-83%, so an inferred frame + can only be sized through a DEM (`width x height above ground / focal length`) and + never from its GSD."* The route the docs name as the **only** route for an inferred + frame is unreachable. `fmt_cross_m`/`fmt_along_m` (:534-535) were written precisely + to serve it and are dead for these rows. + + Same root cause hits calib-matched frames whose `ground_sample_distance` is + missing or zero — a population the code itself calls out as real at :115 + (*"`ground_sample_distance` is 0 on every frame of some digital rolls"*), and which + `digital_fixture()` row 6 exists to represent. Reproduced with the real digital + fixture and `ground_sample_distance <- 0`: 24/24 empty. + +- **[bug]** `R/fly_footprint.R:465` (`unresolved <- is.na(width_in) & !from_table`) + and `:462` (`basis[from_table] <- ...`) — **the frames in the finding above are + silent, and their reporting columns claim they were resolved.** + + Because `from_table` is TRUE, they are excluded from the `unresolved` warning; and + because `sized` is FALSE, none of the three DEM warnings (`uncovered`, `unusable`, + `partial`) can fire either. A direct `fly_footprint()` call emits **no warning at + all** while returning empty geometry for every digital frame. Meanwhile: + + | column | value | what it means per the roxygen | + |---|---|---| + | `footprint_basis` | `"inferred_format"` / `"Digital - Colour"` | "sized from a format inferred from its focal_length" / "format resolved from the format table" | + | `width_source` | `"focal_length=92"` / `"121201_2011"` | "names the calibration file … so every footprint traces back to a source" | + | `footprint_terrain` | `NA` | "no footprint to place — see `footprint_basis`" | + + `footprint_basis` and `footprint_terrain` directly contradict each other, and a + caller filtering on `footprint_basis != "unknown_format"` — the workflow the + `@examples` block teaches at :371 — keeps a set that is entirely empty geometry. + Before #32 these rows were `"unknown_format"` with a warning naming `format_size`, + so this is a regression in signalling even where the geometry outcome is unchanged. + + `fly_warn_unsized()` does not rescue this: it is only called from `fly_coverage()` + / `fly_georef()`, never from `fly_footprint()` itself. + +- **[bug]** `R/fly_footprint.R:507-508` — `fly_footprint()` now aborts on an `NA` in + `film_roll`, with an opaque base-R message. + + The guard is `all(c("film_roll", "frame_number") %in% names(centroids_sf))`, which + tests for *presence*, not usability. `fly_bearing()` then does + `if (i < length(ord) && rolls[i] == rolls[i + 1])` (`R/fly_bearing.R:53`); with an + NA roll that condition is `NA` and `if` errors. Reproduced on the bundled digital + fixture with one `film_roll` set to NA: + + ``` + ERROR: missing value where TRUE/FALSE needed + ``` + + This is latent in `fly_bearing()` but *newly reachable from the package's primary + function*: before this diff `fly_footprint()` never called it. `flight_log_url` is + already NA throughout the bundled digital fixture, so NA string fields in this + layer are not hypothetical. Film-only input is unaffected (`non_square` is all + FALSE, so `fly_bearing()` is not called), which is why the suite is green. + + Either guard on usable values before calling, or make `fly_bearing()` NA-safe + (`identical(rolls[i], rolls[i+1])` / `%in%`). + +- **[fragile]** `R/fly_camera_format.R:134-136` vs `R/fly_footprint.R:460-461` — the + `withheld:` diagnostic is computed and then thrown away. + + `fly_camera_format()` sets `width_source = "withheld:"` for a refused + calibration, but those rows have `resolved = FALSE` / `width_mm = NA`, so + `from_table` is FALSE and `width_source[from_table] <- ...` never copies it. + Reproduced with the two withheld medium-format keys: + + ``` + resolver width_source: withheld:11937933_2009 withheld:12335326_2017 + footprint width_source: NA NA + warning: "... no known recording format (Digital - Colour) ... + see the `format_size` argument" + ``` + + The roxygen at `:257-258` states *"frames naming one are refused rather than + inferred"* and `:255` that `width_source` traces every footprint to a source — + neither is observable from the returned object. The user is told the camera is + unknown when `fly` in fact holds a visually-read spec it deliberately declined to + ship, which is a materially different situation and points them at the wrong fix. + +- **[fragile]** `R/fly_footprint.R:160-169` (`fly_axis_aligned`) — the georef + double-rotation guard tests a proxy, not the property. + + `fly_georef()` skips a frame when its ring is not axis-aligned, standing in for + "was rotated onto its flight line". A bearing at an exact multiple of 90° produces + a footprint that *was* rotated and *is* axis-aligned, so it is not skipped and + `bearing_to_rotation()` applies the shift a second time — the exact + landscape-onto-portrait failure the comment at `R/fly_georef.R:122-127` describes. + The tolerance (`sqrt(eps) * max(abs(xy))` ≈ 1.5 cm at BC Albers magnitudes) makes + this narrow, but the property is knowable exactly: `fly_footprint()` already knows + which rows it rotated. Carrying that as a column (or reusing the + `axis_aligned_no_bearing` marker already in `width_source`) removes the guess. + +- **[fragile]** `R/fly_footprint.R:496` — `by_gsd` gates on `px_cross` only, then + `:525` sets `terrain[by_gsd] <- "gsd_scaled"` unconditionally. + + A table row carrying `px_cross` but `NA` `px_along` would leave `half_along` NA at + :498 → empty polygon at :128, while `footprint_terrain` claims `"gsd_scaled"`, i.e. + a sized frame. Not reachable from today's shipped CSV (every `calib_file` row has + both), and `test-camera_formats.R:65-66` only pins the *fallback* rows as NA — + nothing asserts calib rows carry both. One `& !is.na(fmt$px_along)` closes it. + +## Test gaps (assertions that cannot fail on the above) + +- `tests/testthat/test-fly_camera_format.R:44-59` and `:62-78` are the only DEM + tests on digital input, and both restrict themselves to `sized` rows / rows 1-4 — + which are exactly the `by_gsd` rows. Nothing anywhere asserts that + `digital_fixture()` row 5 (inferred format, the row the docs say *requires* a DEM) + gets a footprint when a DEM is supplied. Adding + `expect_false(sf::st_is_empty(sf::st_geometry(with_dem)[5]))` fails today. +- `tests/testthat/test-fly_camera_format.R:53` computes `sized` from the **no-DEM** + run and then indexes the DEM run by it, so a DEM run that sizes strictly fewer + frames than the flat run cannot be detected. +- `tests/testthat/test-fly_camera_format.R:145-156` covers the "columns absent" path + only. There is no case with `film_roll` present and NA, which is the crash. +- `width_source` was not added to `fly_reported_cols()` (`tests/testthat/setup.R:159`), + so the four-shape class sweep — including the DEM-path one at + `test-fly_footprint.R:692-716`, which is where #35 actually lived — never compares + its values. The replacement at `test-fly_camera_format.R:176-189` sweeps two shapes + on the flat path only, and its loop compares `out$plain` against itself for the + `plain` element. + +## Verified correct (checked, no action) + +- Rotation matrix `R/fly_footprint.R:145`: `xy %*% rot` sends `(0,1)` to + `(sin b, cos b)`. Confirmed numerically against a returned ring at bearing 270.5° — + the along-track edge lies along the heading and the long cross-track edge is + perpendicular. Vertex order BL/BR/TR/TL/BL is preserved through the rotation. +- `bearing` row alignment: `fly_bearing()` writes via `bearing[ord[i]]`, i.e. back + into input row positions, and `coords` comes from `st_transform(centroids_sf)` with + no reordering. Aligned. +- The `!by_gsd` gate at `:581`/`:645` works as intended: with a DEM supplied, GSD-sized + frames keep byte-identical areas and `dem_coverage`/`height_agl` stay NA. +- Film output is bit-identical to before, including in a mixed film+digital frame + (geometry `all.equal` TRUE, areas TRUE), and `fly_axis_aligned()` survives the + 3005 → 4326 → 3005 round trip. +- `take()`'s `<<-` in `R/fly_camera_format.R:100-109`: `rows` is a length-`n` logical + and `from` has `sum(rows)` rows in both call sites; `matched` / `refused` / `digital` + are all length `n` and consistently derived from `key`. No length or index + disagreement found. +- Zero-row input: `fly_footprint(x[0, ])` returns character-typed reporting columns + and does not error; `paste0("withheld:", character(0))` assigns into a zero-length + subscript harmlessly. +- Round-1 fixes spot-checked and correct: the `excl()` empty-key guard, the unzip + guard on PDF presence rather than directory existence, the excluded-table lookup in + the resolver, `FALLBACK_EXTRAPOLATE` as an opt-in property with a stale-entry check, + the `!is.na()`-before-`any()` on `focal_disagrees`, the modal-width index selection + (with an `is.na(w)` backstop for the `%in%` reparse), and the `focal_mm` bound in + check D. diff --git a/planning/archive/2026-08-issue-32-digital-sensor-formats/review-round3.md b/planning/archive/2026-08-issue-32-digital-sensor-formats/review-round3.md new file mode 100644 index 0000000..226d7be --- /dev/null +++ b/planning/archive/2026-08-issue-32-digital-sensor-formats/review-round3.md @@ -0,0 +1,144 @@ +# Review round 3 — `R/fly_footprint.R` + +Scoped to `R/fly_footprint.R`, with `R/fly_camera_format.R`, `R/fly_bearing.R`, +`tests/testthat/test-fly_camera_format.R` and `tests/testthat/test-camera_formats.R` +read for context. Every claim below was reproduced against the working tree. + +## Findings + +- **[bug]** `R/fly_footprint.R:516-525` (with `:666-667` and `:724`) — **a footprint sized + by the DEM route is never rotated onto the flight line, and nothing records that it + wasn't.** This is the mechanism behind round 2's blocker showing up one line later: + `dem_eligible` was corrected to see camera-table rows, but `non_square` — computed at + line 516, *before* the DEM block — was not. + + `non_square` is derived from `half_cross`/`half_along` as they stand at line 516. For + every camera-table row those are `NA` (`width_in` is NA there by definition — the same + fact that made round 2's `sized <- !is.na(half_cross)` unreachable). So `non_square` is + `FALSE` for every DEM-sized digital row, which has two consequences: + + 1. `bearing` is left all-`NA` (line 518 never fires), so the final + `fly_rectangles(coords, half_cross, half_along, bearing)` at line 724 draws an + axis-aligned rectangle. + 2. The `; axis_aligned_no_bearing` annotation at lines 521-525 is also keyed on the + stale `non_square`, so `width_source` stays silent. The documented contract — + *"Where no bearing can be computed the rectangle stays axis-aligned and + `width_source` says so"* (roxygen line 266) — is violated in the one direction that + matters: the bearing **was** computable and the rectangle is axis-aligned anyway. + + Measured, on the suite's own fixture (`digital_fixture()` with + `camera_calibration_url <- NA`, `ground_sample_distance <- NA`, `dem = dem.tif` — + i.e. verbatim the setup of the test at `test-fly_camera_format.R:270`): + + ``` + basis: inferred_format x6 + terrain: dem_agl x6 + width_source: focal_length=92 focal_length=92 focal_length=80 ... <- no annotation + fly_bearing(p)$bearing: 90.49 90.45 90.42 90.42 NA NA <- computable + row 1: cross=2595.1 along=1471.8 along-edge azimuth = 360.00 <- should be ~90.5 + ``` + + Row 1 is a DMC III (100.3 x 56.9 mm, 1.76:1) on a flight line running due east. The + long axis is drawn **along** the heading instead of across it — the footprint is + transposed 90 degrees. Area is unchanged, so `st_area` cannot see it, but the ground + actually covered is wrong, and `fly_coverage()`, `fly_overlap()` and photo selection all + read that ground. + + The DEM sampling windows use the same all-`NA` bearing, so `dem_coverage` and + `height_agl` are internally consistent — consistent with the wrong shape. + + Second, quieter symptom of the same mechanism: whether a row is rotated depends on + **what else is in the same call**, because `bearing` is assigned as a whole vector + under `if (any(non_square))`. Same fixture, same row 1, differing only in whether + unrelated rows carry a GSD: + + ``` + rows 1-4 keep GSD (by_gsd -> non_square TRUE) -> row 1 azimuth = 90.49 (rotated) + rows 1-4 GSD set to NA (all on the DEM route) -> row 1 azimuth = 0.00 (not rotated) + ``` + + A frame's geometry changing because a *different* frame in the batch had a + `ground_sample_distance` is the "two conditions that only happen to agree" shape. + + Why the suite can't see it: the two tests that exercise the DEM route for digital + (`test-fly_camera_format.R:270` "an inferred-format frame IS sized when a dem is + supplied" and `:291` "a calibrated frame with no usable gsd is sized by the dem + instead") assert `footprint_basis`, non-emptiness, `footprint_terrain` and + `height_agl > 0` — nothing about orientation. The only rotation assertions + (`:119` and `:148`) both run **without** a DEM, so they only ever reach the + `by_gsd` path where `non_square` happens to be computable in time. No assertion in the + package currently reaches a rotated DEM-sized footprint. + + Fix shape: `non_square` (and therefore `bearing`, and therefore the + `axis_aligned_no_bearing` annotation) has to be decided from the dimensions the row + will actually ship — i.e. the format's aspect, which is known from `fmt$width_mm` / + `fmt$height_mm` before any sizing route runs — not from the half-dimensions as they + stand mid-function. Note that `bearing` must be available *before* line 590, because + the seed and second-pass sampling rectangles use it too. + + Regression test that would have caught it: repeat the `:119` azimuth assertion + ("along-track edge must point along the heading, mod 360") on the DEM route, and assert + `width_source` gains `axis_aligned_no_bearing` **only** where `fly_bearing()` returns + NA. + +- **[fragile]** `R/fly_footprint.R:697-706` — the new `resolved_unsized` warning fires on a + **film** frame whose only problem is an unparseable `scale`, and its text points + entirely at digital metadata. Reproduced on the bundled film centroids with + `scale[2] <- "unknown"`: + + ``` + WARN: 1 of 3 frames have a known recording format but could not be sized, so they + have no footprint. A digital frame needs either a `ground_sample_distance` and a + calibrated pixel count, or `dem` together with `flying_height` and `focal_length`. + ``` + + `basis` is `Film - BW`, `width_source` is `NA`, and neither `ground_sample_distance` + nor `dem` would help — the fix is the `scale` string. Same misfire on the + no-`media`-column path (`assumed_default`). This case produced no `fly_footprint()` + warning at all before this diff, so it is new behaviour, and it sends the reader to the + wrong column. Splitting the message on `from_table` (digital advice) versus + `!is.na(width_in) & is.na(scale_num)` (say `scale` could not be parsed) costs one + branch. + +## Checked and clean + +- **`half_cross` / `half_along` mutual consistency.** `half_along` is only ever set + independently on the two paths that also set `half_cross` (`by_gsd` at 508-509, + `corrected` at 666-667). They cannot disagree in NA-ness, so `no_geom` at 687 and + `fly_rectangles`' own guard classify identically for every value reachable from + catalogue data. +- **Invariant "footprint_terrain is NA exactly where the geometry is empty".** Holds on + every path I could reach. The one structural gap is that `no_geom` tests `is.na` while + `fly_rectangles` tests `is.finite`, so an `Inf` half-dimension would ship an empty + geometry with a non-NA `footprint_terrain` — but I could not reach `Inf` from any + realistic input: the DEM route filters on `is.finite(candidate$...)` before assigning, + and `by_gsd` requires a finite positive GSD. Noted, not reported. +- **Sea-level seed is safe for film.** `need_seed <- dem_eligible & is.na(seed_cross)` is + `FALSE` for every film row with a parseable scale, so film seeds from its nominal + rectangle exactly as before. A film row with an unparseable scale is not `dem_eligible` + at all (`!is.na(half_cross) | from_table` is FALSE), so it never reaches `resize(0)`. + Film output is byte-identical: `test-camera_formats.R:97` passes, and I re-confirmed + all 20 bundled footprints are square to 4.1e-12 relative. +- **Seed non-finiteness.** `resize(0)` returns NA/Inf for an NA or zero `focal_length`; + `fly_rectangles` turns that into an empty polygon, `fly_dem_sample` returns NA elev, and + the row lands in `uncovered` (rather than `unusable`, which is a slightly misleading + warning but the final state — empty geometry, NA terrain — is correct). +- **`fly_is_square()` tolerance.** Not a false-positive risk: real film footprints + round-tripped through EPSG:4326 show a maximum relative edge spread of 4.1e-12 against + `all.equal`'s 1.5e-8 — four orders of margin. All 20 classify square from both a 3005 + and a 4326 input. +- **`fly_bearing()` row order.** Assigns via `bearing[ord[i]]`, so the returned vector is + in input order and aligns with `coords` / `half_cross`. Confirmed by the 90.4-90.5 + values landing on rows 1-4. +- **`dem_coverage` / `terrain` exclusion of `by_gsd`.** Correct on both branches: the + `gsd_scaled` rows are sampled (wastefully) by the seed pass but excluded from + `dem_eligible`, so neither `terrain[dem_eligible]` nor `dem_coverage[dem_eligible]` + touches them. +- **Test assertions.** No assertion in either test file is structurally incapable of + failing. `test-camera_formats.R:18` (check B) correctly excludes the `!mm_stated` rows + that would compare arithmetic against itself, and asserts the premise. The + `expect_gt(..., 5 * X / 5)` at `test-fly_camera_format.R:26` is an odd expression but + does discriminate the scale route from the GSD route (8.80e6 vs 3.04e6). The + drift-guard-fires test at `test-camera_formats.R:106` re-implements the setdiff rather + than calling the guard, so it duplicates rather than proves — but it does have a real + failure mode, so it is not decoration. diff --git a/planning/archive/2026-08-issue-32-digital-sensor-formats/task_plan.md b/planning/archive/2026-08-issue-32-digital-sensor-formats/task_plan.md new file mode 100644 index 0000000..4552b56 --- /dev/null +++ b/planning/archive/2026-08-issue-32-digital-sensor-formats/task_plan.md @@ -0,0 +1,65 @@ +# Task: Establish sensor widths for digital frames and ship them as format defaults (#32) + +## Problem + +`fly_footprint()` sizes each frame from its `media` value and returns an empty geometry +for anything it cannot resolve (#30). That made the gap honest but did not close it: +every `Digital - Colour` frame still has no footprint — **223,667 of 1,670,471 frames +(13.4%)** provincially. The package is quietly film-only for anything after ~2010. + +What was missing is one number per camera — the sensor width — which is not in +`AIMG_PHOTO_CENTROIDS_SP`. Exploration established it is recoverable from +`CAMERA_CALIBRATION_URL`, and turned up two things that change the shape of the work: +the catalogue's `SCALE` is not the true image scale for digital frames (sizing from it +gives 34% of true width), and digital sensors are not square (aspect 1.10 to 1.80). + +Full plan: `~/.claude/plans/federated-mixing-wilkes.md` + +## Phase 1: Camera format table + +- [x] `data-raw/make_camera_formats.R` — fetch the 18 calibration zips, extract per vendor +- [x] `inst/extdata/camera_formats.csv` — calibration rows + focal-length fallback rows +- [x] Record source PDF and retrieval date per row (snapshot, not contract) +- [x] Tests: key uniqueness, key_type domain, source_pdf present on calibration rows +- [x] Drift guard: CSV key set matches what the generator discovers; prove the alarm fires + +## Phase 2: QA on the transcribed numbers + +- [x] B — `px x pitch == stated image size`, both axes, all rows +- [x] C — report focal vs catalogue `FOCAL_LENGTH` within 1 mm +- [x] D — plausibility bounds (width 80-170 mm, aspect 1.0-2.0, pitch 3-13 um) +- [x] E — independent second extraction diffed against the CSV (double entry) +- [x] F — implied ground elevation `H - (GSD/pitch)*focal` is a plausible BC elevation +- [x] Disposition: ship `dmc100039_2006` with its focal note; withhold `72914123_2019` + +## Phase 3: Resolver and non-square footprints + +- [x] `fly_camera_format()` internal resolver, keyed on `camera_calibration_url` +- [x] `fly_rectangles()` gains a second half-dimension and optional rotation +- [x] Rotate non-square footprints to the flight line via `fly_bearing()` +- [x] Tests: film output unchanged; rotation swaps bbox and preserves area +- [x] Restore the square-only version and confirm the non-square test goes red + +## Phase 4: Wire into `fly_footprint()` + +- [x] Sizing precedence: film scale / `px x GSD` / `width x AGL/focal` / unknown +- [x] `SCALE` never used for a digital frame +- [x] New `width_source` column; `footprint_basis` gains inferred values +- [x] `format_size` accepts width or width x height; roxygen example corrected +- [x] Warning text names refused vs inferred counts +- [x] Tests: `centroid_shapes()` sweep, zero-row, all-unresolved, downstream consumers + +## Phase 5: Fixtures and verification + +- [x] `inst/testdata/photo_centroids_digital.gpkg` — real frames, two aspect ratios +- [x] G — the two sizing routes agree within 2% (the issue's own acceptance criterion) +- [x] H — warn when implied AGL exceeds `flying_height` +- [x] Digital fixture in `setup.R` covering every resolver branch +- [x] Vignette, `NEWS.md`, correct the `CLAUDE.md` claim about digital frames in the AOI + +## Validation + +- [x] Tests pass +- [x] `/code-check` clean on each commit +- [x] PWF checkboxes match landed work +- [x] `/planning-archive` on completion diff --git a/tests/testthat/setup.R b/tests/testthat/setup.R index d286b46..316113a 100644 --- a/tests/testthat/setup.R +++ b/tests/testthat/setup.R @@ -4,21 +4,30 @@ testdata_path <- function(...) { } -# Synthesized mixed film/digital fixture. +# Synthesized mixed film/digital fixture, where the digital frames CANNOT be resolved. # -# The bundled centroids are 100% `Film - BW`: they come from a 1968 AOI near -# Houston (data-raw/make_testdata.R), which has no digital coverage at all. A -# digital frame therefore has to be constructed rather than sampled. +# The bundled centroids are 100% `Film - BW`, so a digital frame has to be constructed +# rather than sampled. # -# Focal lengths and GSD mirror the catalogue's two tells for a sensor: 92/100 mm -# against 153/305 mm for film, and a populated GROUND_SAMPLE_DISTANCE. +# Focal length 53 is deliberate and is not an arbitrary number: it is the catalogue +# focal of the PhaseOne rolls, whose calibration report is an image-only PDF and which +# `camera_formats_excluded.csv` therefore records as unresolvable. So it is a real +# unresolvable case rather than one invented to dodge the format table — and because +# no frame at focal 53 lacks a calibration, there is no focal-length fallback row for +# it either. +# +# That premise is asserted in the tests that depend on it (`test-fly_footprint.R`), +# so a future table that starts shipping focal 53 fails naming the real cause rather +# than failing on the behaviour under test. +# +# GROUND_SAMPLE_DISTANCE is CENTIMETRES in the catalogue — see `fly_gsd_m()`. mixed_media_fixture <- function() { sf::st_sf( airp_id = 1:4, scale = c("1:12000", "1:12000", "1:15000", "1:15000"), media = c("Film - BW", "Film - Colour", "Digital - Colour", "Digital - Colour"), - focal_length = c(153, 305, 92, 100), - ground_sample_distance = c(NA, NA, 0.10, 0.10), + focal_length = c(153, 305, 53, 53), + ground_sample_distance = c(NA, NA, 10, 10), geometry = sf::st_sfc( sf::st_point(c(-126.60, 54.40)), sf::st_point(c(-126.58, 54.40)), @@ -30,6 +39,49 @@ mixed_media_fixture <- function() { } +# Digital frames that DO resolve, one per key type and with two very different sensor +# shapes. +# +# Rows 1-2 are a Leica DMC II (87.1 x 79.2 mm, aspect 1.10) and rows 3-4 an UltraCam +# Eagle (105.8 x 68.0 mm, aspect 1.56), both keyed by their real calibration file. Two +# aspect ratios that far apart are what lets a test tell a correctly-shaped footprint +# from a square one — a fixture carrying a single shape cannot. +# +# Row 5 has no calibration URL and focal 100, so it takes the focal-length fallback and +# is marked inferred. Row 6 has a calibration but GSD 0, which is the state that would +# build a degenerate five-identical-vertex polygon if the zero guard were removed. +# +# `film_roll` and `frame_number` are present so `fly_bearing()` can resolve a flight +# line; frames 1-4 run west to east. +digital_fixture <- function() { + url <- function(k) paste0("https://openmaps.gov.bc.ca/thumbs/calib_report_zips/", k, ".zip") + sf::st_sf( + airp_id = 1:6, + scale = c("1:20000", "1:20000", "1:20000", "1:20000", "1:20000", "1:20000"), + media = rep("Digital - Colour", 6), + film_roll = c("a", "a", "a", "a", "b", "c"), + frame_number = c(1, 2, 3, 4, 1, 1), + focal_length = c(92, 92, 80, 80, 100, 120), + flying_height = rep(3000, 6), + ground_sample_distance = c(20, 20, 15, 15, 25, 0), + camera_calibration_url = c( + url("121201_2011"), url("121201_2011"), + url("20814295_2018"), url("20814295_2018"), + NA, url("dmc100039_2006") + ), + geometry = sf::st_sfc( + sf::st_point(c(-126.62, 54.40)), + sf::st_point(c(-126.58, 54.40)), + sf::st_point(c(-126.54, 54.40)), + sf::st_point(c(-126.50, 54.40)), + sf::st_point(c(-126.46, 54.40)), + sf::st_point(c(-126.42, 54.40)), + crs = 4326 + ) + ) +} + + # Skip a terrain test when terra is unavailable. # # `terra` is in Suggests, not Imports — the DEM path is optional. A test that @@ -107,3 +159,29 @@ mixed_media_shapes <- function() { fly_reported_cols <- function() { c("footprint_basis", "footprint_terrain", "height_agl", "dem_coverage") } + + +# Every input shape `fly_footprint()` can be handed, for the invariant sweep in +# test-fly_footprint_invariants.R. Lives here rather than in the test file so its +# dependencies — the other fixtures and `testdata_path()` — are defined alongside it. +footprint_cases <- function() { + dem <- testdata_path("dem.tif") + film <- sf::st_read(testdata_path("photo_centroids.gpkg"), quiet = TRUE) + digital <- sf::st_read(testdata_path("photo_centroids_digital.gpkg"), quiet = TRUE) + no_gsd <- digital + no_gsd$ground_sample_distance <- NA + list( + "film, no dem" = list(x = film, dem = NULL), + "film, dem" = list(x = film, dem = dem), + "digital, no dem" = list(x = digital, dem = NULL), + "digital, dem" = list(x = digital, dem = dem), + "digital, no gsd, dem" = list(x = no_gsd, dem = dem), + "digital, no gsd, no dem" = list(x = no_gsd, dem = NULL), + "unknown format" = list(x = mixed_media_fixture(), dem = NULL), + "mixed resolvability" = list(x = digital_fixture(), dem = NULL), + "mixed resolvability + dem" = list(x = digital_fixture(), dem = dem), + "terrain fixture" = list(x = terrain_fixture(), dem = dem), + "no media column" = list(x = film[, setdiff(names(film), "media")], dem = NULL), + "empty input" = list(x = digital_fixture()[0, ], dem = NULL) + ) +} diff --git a/tests/testthat/test-camera_formats.R b/tests/testthat/test-camera_formats.R new file mode 100644 index 0000000..5ef007d --- /dev/null +++ b/tests/testthat/test-camera_formats.R @@ -0,0 +1,181 @@ +test_that("the shipped camera format table is structurally sound", { + tbl <- fly_camera_table() + + expect_gt(nrow(tbl), 0) + expect_equal(anyDuplicated(tbl$key), 0L) + expect_true(all(tbl$key_type %in% c("calib_file", "focal_length"))) + expect_true(all(nzchar(tbl$key))) + + # Every calibration row must name the PDF it came from. Provenance is the point of + # the artifact — a row whose source cannot be reopened cannot be re-checked. + calib <- tbl[tbl$key_type == "calib_file", ] + expect_true(all(!is.na(calib$source_pdf) & nzchar(calib$source_pdf))) + expect_true(all(!is.na(calib$camera) & nzchar(calib$camera))) + expect_true(all(!is.na(tbl$retrieved))) +}) + + +test_that("check B: pixel count times pitch reproduces the stated image size", { + # The strongest guard on the transcription, and only on the rows where it is real: + # `mm_stated` is FALSE where the report gave array size and pitch but no millimetres, + # so `width_mm` there IS `px * pitch` and comparing them tests arithmetic against + # itself. Skipping those is the difference between a check and a vacuous pass. + tbl <- fly_camera_table() + b <- tbl[!is.na(tbl$mm_stated) & tbl$mm_stated, ] + + # Assert the premise, so a future table that loses its constrainable rows fails here + # naming the real cause rather than passing on an empty set. + expect_gt(nrow(b), 10) + + expect_lt(max(abs(b$px_cross * b$pitch_um / 1000 - b$width_mm) / b$width_mm), 1e-6) + expect_lt(max(abs(b$px_along * b$pitch_um / 1000 - b$height_mm) / b$height_mm), 1e-6) +}) + + +test_that("check D: every shipped format is physically plausible", { + # Catches a unit slip (m / mm / um), which is the error class that survives check B by + # being internally self-consistent. The reports write the micron sign three different + # ways and one drops it altogether, so this is not hypothetical. + tbl <- fly_camera_table() + + # 30-200 mm, matching the generator. Deliberately wider than the 87-166 mm the + # shipped large-format cameras occupy: the catalogue also holds medium-format bodies + # about 53 mm wide, and a bound that rejects valid data is worse than the unit slip it + # guards against. This range still separates a micrometre read as a millimetre (~5 mm) + # from a real sensor. + expect_true(all(tbl$width_mm > 30 & tbl$width_mm < 200)) + expect_true(all(tbl$height_mm > 20 & tbl$height_mm < 150)) + aspect <- tbl$width_mm / tbl$height_mm + expect_true(all(aspect >= 1 & aspect <= 2)) + + pitch <- tbl$pitch_um[!is.na(tbl$pitch_um)] + expect_true(all(pitch > 3 & pitch < 13)) +}) + + +test_that("fallback rows carry no pixel count, so they cannot take the GSD route", { + # At a given catalogue focal length, sensor WIDTH spreads 1-3% but PIXEL COUNT spreads + # 32-83%. A focal-keyed frame can therefore be sized by width x AGL/focal and never by + # px x GSD. Withholding the pixel counts is what enforces that, rather than a comment + # asking the resolver to remember. + tbl <- fly_camera_table() + fb <- tbl[tbl$key_type == "focal_length", ] + + expect_gt(nrow(fb), 0) + expect_true(all(is.na(fb$px_cross))) + expect_true(all(is.na(fb$px_along))) + expect_true(all(!is.na(fb$width_mm) & !is.na(fb$height_mm))) + # Each records how much room for error it carries. Spread is required only on rows + # inferred from cameras at the SAME focal length — an extrapolated row has no + # meaningful spread and says so with NA, asserted separately below. + expect_true(all(!is.na(fb$n_cameras))) + expect_true(all(!is.na(fb$width_spread_pct[!fb$extrapolated]))) + expect_true(all(!is.na(fb$note) & nzchar(fb$note))) +}) + + +test_that("every calibration the catalogue offers is dispositioned, with a reason", { + # The drift guard. `camera_formats_manifest.csv` records the keys the catalogue + # actually offered when the table was built, written independently of the two + # dispositions — so a key hand-edited out of either file, or a new camera added to + # the manifest without being dispositioned, fails here rather than silently becoming + # an unresolvable frame. + manifest <- utils::read.csv( + system.file("extdata", "camera_formats_manifest.csv", package = "fly"), + stringsAsFactors = FALSE + ) + excluded <- fly_camera_excluded() + # The excluded file also records focal lengths deliberately left without a fallback + # row, which are not calibrations and so are not in the manifest. Compare like with + # like rather than pooling them — a pooled comparison passes for exactly the drift + # the guard exists to catch. + excluded_calib <- excluded$key[excluded$key_type == "calib_file"] + shipped <- fly_camera_table() + shipped <- shipped$key[shipped$key_type == "calib_file"] + + expect_gt(nrow(manifest), 0) + expect_setequal(manifest$key, c(shipped, excluded_calib)) + expect_equal(length(intersect(shipped, excluded_calib)), 0L) + + # An excluded entry without a reason is a backlog note pretending to be a decision. + expect_true(all(!is.na(excluded$reason) & nzchar(excluded$reason))) + expect_true(all(excluded$key_type %in% c("calib_file", "focal_length"))) +}) + + +test_that("the drift guard fires on an undeclared calibration", { + # A guard nobody has seen fail is decoration. Feed it a key that is in the manifest + # and in neither disposition, and confirm it is reported. + manifest <- utils::read.csv( + system.file("extdata", "camera_formats_manifest.csv", package = "fly"), + stringsAsFactors = FALSE + ) + shipped <- fly_camera_table() + shipped <- shipped$key[shipped$key_type == "calib_file"] + ex <- fly_camera_excluded() + excluded <- ex$key[ex$key_type == "calib_file"] + + undeclared <- c(manifest$key, "ultracam999_2031") + expect_false(setequal(undeclared, c(shipped, excluded))) + expect_equal(setdiff(undeclared, c(shipped, excluded)), "ultracam999_2031") +}) + + +test_that("an extrapolated fallback row does not claim zero spread", { + # `width_spread_pct` measures dispersion among the SOURCE cameras, so a single-source + # row is structurally 0 however wrong the inference. Reporting 0 would give the + # least-supported rows in the file the most confident label. An extrapolated row + # reports NA — unknown, which is the truth. + tbl <- fly_camera_table() + ex <- tbl[!is.na(tbl$extrapolated) & tbl$extrapolated, ] + + expect_gt(nrow(ex), 0) # premise: something is extrapolated + expect_true(all(is.na(ex$width_spread_pct))) + # And it must say why, at length — an extrapolation is the one row type whose + # warrant cannot be read off the numbers. + expect_true(all(nchar(ex$note) > 40)) +}) + + +test_that("a withheld calibration cannot fall through to focal-length inference", { + # The two withheld medium-format cameras are ~53 mm wide against the ~104 mm bodies + # sharing their focal neighbourhood, so inferring one would be 1.95x too wide. Today + # they return NA because no fallback row happens to exist at focal 53 or 83 — safety + # by coincidence. This pairs each withheld key with focal 100, which DOES have a + # fallback row, so only the refusal itself can stop the fall-through. + ex <- fly_camera_excluded() + withheld <- ex$key[ex$key_type == "calib_file"] + expect_gt(length(withheld), 0) + + tbl <- fly_camera_table() + expect_true("100" %in% tbl$key[tbl$key_type == "focal_length"]) # premise + + photos <- sf::st_sf( + media = rep("Digital - Colour", length(withheld)), + focal_length = rep(100, length(withheld)), + camera_calibration_url = paste0( + "https://openmaps.gov.bc.ca/thumbs/calib_report_zips/", withheld, ".zip" + ), + geometry = sf::st_sfc( + lapply(seq_along(withheld), function(i) sf::st_point(c(-126, 54))), crs = 4326 + ) + ) + got <- fly_camera_format(photos) + + expect_true(all(!got$resolved)) + expect_true(all(is.na(got$width_mm))) + expect_true(all(grepl("^withheld:", got$width_source))) +}) + + +test_that("resolved and inferred are distinct, so !inferred is a safe filter", { + # A row that resolved to nothing is also `inferred = FALSE`, so without `resolved` the + # natural "give me the trustworthy rows" filter sweeps up every unresolved frame. + photos <- digital_fixture() + got <- fly_camera_format(photos) + + expect_true(all(got$resolved[1:5])) + expect_false(got$inferred[1]) # exact calibration + expect_true(got$inferred[5]) # focal-length fallback + expect_true(all(!is.na(got$width_mm[got$resolved]))) +}) diff --git a/tests/testthat/test-fly_camera_format.R b/tests/testthat/test-fly_camera_format.R new file mode 100644 index 0000000..0c35429 --- /dev/null +++ b/tests/testthat/test-fly_camera_format.R @@ -0,0 +1,389 @@ +# Sizing digital frames from the camera format table (#32). + +test_that("digital frames resolve from their calibration url", { + fp <- suppressWarnings(fly_footprint(digital_fixture())) + + # Rows 1-4 carry a calibration; row 5 falls back to focal length; row 6 has a + # calibration but no ground sample distance. + expect_equal(fp$width_source[1:2], rep("121201_2011", 2)) + expect_equal(fp$width_source[3:4], rep("20814295_2018", 2)) + expect_equal(fp$footprint_basis[1:4], rep("Digital - Colour", 4)) + expect_equal(fp$footprint_basis[5], "inferred_format") + expect_match(fp$width_source[5], "^focal_length=100") +}) + + +test_that("a digital footprint is the pixel count times the ground sample distance", { + # The measurement this whole issue rests on. GROUND_SAMPLE_DISTANCE is centimetres, so + # a Leica DMC II at GSD 20 covers 15552 px x 0.20 m = 3110.4 m across-track and + # 14144 px x 0.20 m = 2828.8 m along-track. + fp <- suppressWarnings(fly_footprint(digital_fixture())) + g <- sf::st_transform(fp[1, ], 3005) + + expect_equal(as.numeric(sf::st_area(g)), 15552 * 0.20 * 14144 * 0.20, tolerance = 1e-6) + # Sized from the sensor, not from `scale`: the scale route would give + # 87.0912 mm x 20000 = 1742 m, and the area would be off by a factor of 3.2. + expect_gt(as.numeric(sf::st_area(g)), 5 * (87.0912e-3 * 20000)^2 / 5) +}) + + +test_that("scale is never used to size a frame fly resolved itself", { + # Two frames identical but for `scale`. If `scale` reached the digital route at all, + # their footprints would differ. + a <- digital_fixture() + b <- a + b$scale <- "1:99999" + fa <- suppressWarnings(fly_footprint(a)) + fb <- suppressWarnings(fly_footprint(b)) + + area <- function(x) as.numeric(sf::st_area(sf::st_transform(x[1:5, ], 3005))) + expect_equal(area(fa), area(fb)) +}) + + +test_that("a resolved digital frame keeps its footprint when a dem is supplied", { + # The precedence rule. Every downstream consumer forwards `dem`, so if the DEM route + # overwrote a frame already sized from its ground sample distance, passing `dem` + # would silently change the answer on the ordinary path rather than in an edge case. + skip_if_no_terra() + photos <- digital_fixture() + no_dem <- suppressWarnings(fly_footprint(photos)) + with_dem <- suppressWarnings(fly_footprint(photos, dem = testdata_path("dem.tif"))) + + sized <- !sf::st_is_empty(sf::st_geometry(no_dem)) + expect_gt(sum(sized), 0) + expect_equal( + as.numeric(sf::st_area(sf::st_transform(with_dem[sized, ], 3005))), + as.numeric(sf::st_area(sf::st_transform(no_dem[sized, ], 3005))) + ) +}) + + +test_that("a gsd-sized frame does not claim a terrain treatment it did not get", { + # `nominal_scale` is documented as "sized from the reported scale", which is the one + # thing this route deliberately never does. Claiming it would be a false affirmative + # about the frame's provenance. + fp <- suppressWarnings(fly_footprint(digital_fixture())) + expect_equal(fp$footprint_terrain[1:4], rep("gsd_scaled", 4)) + + # And with a DEM supplied, the sampled coverage describes a window that had no + # bearing on the returned geometry, so it is not reported for those frames. + skip_if_no_terra() + with_dem <- suppressWarnings( + fly_footprint(digital_fixture(), dem = testdata_path("dem.tif")) + ) + expect_equal(with_dem$footprint_terrain[1:4], rep("gsd_scaled", 4)) + expect_true(all(is.na(with_dem$dem_coverage[1:4]))) + expect_true(all(is.na(with_dem$height_agl[1:4]))) +}) + + +test_that("a zero ground sample distance yields an empty geometry, not a degenerate one", { + # `is.finite(0)` is TRUE, so without an explicit guard a zero builds a rectangle with + # five identical vertices: `st_is_empty()` reports FALSE, `fly_warn_unsized()` never + # mentions it, and it covers nothing downstream while looking like a real footprint. + # Reachable — GSD is 0 on every frame of the dmc100039 rolls. + photos <- digital_fixture() + expect_equal(photos$ground_sample_distance[6], 0) + + # And it is reported: a frame whose format resolved but which could not be sized is + # otherwise silent — `footprint_basis` names a real format and `width_source` names a + # calibration, so nothing about the row says the geometry is empty. + expect_warning(fp <- fly_footprint(photos), "no way to size it") + expect_true(sf::st_is_empty(sf::st_geometry(fp)[6])) +}) + + +test_that("film is sized and shaped exactly as before", { + # The regression net. Explicit assertions rather than a snapshot, which skips on CRAN + # and therefore in CI — the run where it would matter. + centroids <- sf::st_read(testdata_path("photo_centroids.gpkg"), quiet = TRUE) + fp <- fly_footprint(centroids) + + expect_true(all(fp$footprint_basis == "Film - BW")) + expect_true(all(is.na(fp$width_source))) + + side <- 9 * as.numeric(sub("1:", "", centroids$scale)) * 0.0254 + expect_equal(as.numeric(sf::st_area(sf::st_transform(fp, 3005))), side^2, + tolerance = 1e-6) + + # Square, so unrotated: every footprint's bounding box is as wide as it is tall. + bb <- vapply(sf::st_geometry(sf::st_transform(fp, 3005)), function(g) { + b <- sf::st_bbox(g) + unname((b["xmax"] - b["xmin"]) / (b["ymax"] - b["ymin"])) + }, numeric(1)) + expect_equal(bb, rep(1, nrow(fp)), tolerance = 1e-8) +}) + + +test_that("a non-square footprint puts its long axis ACROSS the flight line", { + # The assertion that actually discriminates. Area is invariant under every rotation + # and a 90-degree bbox swap only tests the mechanism, so both a transposed long axis + # and a sign-flipped azimuth would pass those. This measures the angle itself. + # + # Frames 1-4 of the fixture run west to east, so bearing is ~90 degrees and the long + # (cross-track) axis must point NORTH-SOUTH. + photos <- digital_fixture() + fp <- sf::st_transform(suppressWarnings(fly_footprint(photos)), 3005) + bearing <- fly_bearing(photos)$bearing + expect_equal(bearing[1], 90, tolerance = 1) + + ring <- sf::st_coordinates(sf::st_geometry(fp)[[1]])[1:4, 1:2] + # Edge 1->2 is the cross-track edge; edge 2->3 the along-track one. + edge <- function(a, b) sqrt(sum((ring[b, ] - ring[a, ])^2)) + cross_len <- edge(1, 2) + along_len <- edge(2, 3) + expect_gt(cross_len, along_len) # DMC II is 87.1 x 79.2 mm + + # The along-track edge must point ALONG the heading, and the modulus here is + # load-bearing: taken mod 180 this assertion passes under a sign-flipped rotation + # (269.5 becomes 89.5), which is the single error it exists to catch. Verified by + # restoring that defect and watching this go red — mod 180 it stayed green. + v <- ring[3, ] - ring[2, ] + az <- unname((atan2(v[1], v[2]) * 180 / pi) %% 360) + expect_equal(az, 90, tolerance = 1) +}) + + +test_that("rotation is skipped, and recorded, when no bearing can be computed", { + photos <- digital_fixture() + photos$film_roll <- NULL + photos$frame_number <- NULL + fp <- suppressWarnings(fly_footprint(photos)) + + expect_match(fp$width_source[1], "axis_aligned_no_bearing") + bb <- sf::st_bbox(sf::st_geometry(sf::st_transform(fp, 3005))[[1]]) + # Axis-aligned: the bounding box is the sensor's own shape, 15552 x 14144 px at 0.20 m + expect_equal(unname(bb["xmax"] - bb["xmin"]), 15552 * 0.20, tolerance = 1e-6) + expect_equal(unname(bb["ymax"] - bb["ymin"]), 14144 * 0.20, tolerance = 1e-6) +}) + + +test_that("the unknown-format contract still holds for a format nothing resolves", { + # Premise first, so a future table that starts shipping focal 53 fails here naming the + # real cause rather than failing on the behaviour under test. + photos <- mixed_media_fixture() + expect_equal(unique(photos$focal_length[photos$media == "Digital - Colour"]), 53) + + tbl <- fly_camera_table() + expect_false("53" %in% tbl$key[tbl$key_type == "focal_length"]) + ex <- fly_camera_excluded() + expect_true(any(grepl("12335326", ex$key))) # the focal-53 camera, unresolvable + + fp <- suppressWarnings(fly_footprint(photos)) + expect_equal(fp$footprint_basis[3:4], rep("unknown_format", 2)) + expect_true(all(sf::st_is_empty(sf::st_geometry(fp)[3:4]))) +}) + + +test_that("width_source survives every class shape a caller can supply", { + # #35 shipped through two releases because `st_sf()` silently drops trailing columns + # when its first argument is a tibble, which is what `bcdata::collect()` returns. + # The sweep is only worth running on a fixture where the new column actually varies. + mm <- digital_fixture() + shapes <- list(plain = mm, tbl = sf::st_as_sf(dplyr::as_tibble(mm))) + stopifnot(!inherits(shapes$plain, "tbl_df"), inherits(shapes$tbl, "tbl_df")) + + out <- lapply(shapes, function(x) suppressWarnings(fly_footprint(x))) + for (nm in names(out)) { + expect_true("width_source" %in% names(out[[nm]]), info = nm) + expect_equal(out[[nm]]$width_source, out$plain$width_source, info = nm) + } +}) + + +test_that("reporting columns keep their types on empty input", { + # `ifelse(logical(0), ...)` returns `logical(0)`, so a query that matched no frames + # would report a character column as logical and fail to bind to a populated result. + fp <- fly_footprint(digital_fixture()[0, ]) + + expect_equal(nrow(fp), 0L) + expect_type(fp$width_source, "character") + expect_type(fp$footprint_basis, "character") + expect_type(fp$footprint_terrain, "character") +}) + + +test_that("check G: the two independent sizing routes agree on real frames", { + # The issue's own acceptance criterion. `px x GSD` and `width x AGL/focal` share no + # arithmetic — the first uses the pixel count and the catalogue's ground sample + # distance, the second the millimetre width, `flying_height` and a DEM — so a + # mis-keyed camera breaks the agreement. + # + # What this can and cannot catch, measured rather than assumed. GROUND_SAMPLE_DISTANCE + # is stored as INTEGER centimetres, so at GSD 12 one unit is 8% and the quantization + # alone bounds the agreement at about +/-4%; at GSD 30 it is 1.7%. Measured on the + # bundled frames: DMC II (GSD 30) agrees to 1.4%, UltraCam Eagle M3 (GSD 12) to 7.5%. + # So this guards against a camera keyed to the wrong row — the widths in the table run + # from 87.1 to 165.9 mm, up to 90% apart — and NOT against a small width error. + # Check B is what guards precision. + skip_if_no_terra() + p <- sf::st_read(testdata_path("photo_centroids_digital.gpkg"), quiet = TRUE) + fmt <- fly_camera_format(p) + expect_true(all(fmt$resolved)) # premise + + route_gsd <- fmt$px_cross * fly_gsd_m(p$ground_sample_distance) + + # Mean elevation under the footprint actually returned, not a centroid reading — the + # two differ by up to 140 m on a frame this wide. + fp <- fly_footprint(p) + dem <- terra::rast(testdata_path("dem.tif")) + samp <- fly_dem_sample(dem, sf::st_geometry(sf::st_transform(fp, 3005))) + route_dem <- (fmt$width_mm / 1000) * + (p$flying_height - samp$elev) / (p$focal_length / 1000) + + ok <- !is.na(route_dem) & samp$covered > 0.95 + expect_gt(sum(ok), 15) # premise: enough covered frames + + ratio <- route_dem[ok] / route_gsd[ok] + expect_lt(max(abs(ratio - 1)), 0.15) + + # The camera whose GSD is least quantized must agree much more closely. Without this + # the 15% bound above would pass a genuinely mis-sized frame. + fine <- ok & p$ground_sample_distance >= 25 + expect_gt(sum(fine), 5) + fine_ratio <- route_dem[fine] / route_gsd[fine] + expect_lt(abs(median(fine_ratio) - 1), 0.03) +}) + + +test_that("real digital frames get non-square footprints in two distinct shapes", { + # A fixture carrying one sensor shape cannot tell a correctly-shaped footprint from a + # square one. These are real catalogue frames from two cameras 0.46 apart in aspect. + p <- sf::st_read(testdata_path("photo_centroids_digital.gpkg"), quiet = TRUE) + fmt <- fly_camera_format(p) + + aspect <- round(fmt$width_mm / fmt$height_mm, 2) + expect_setequal(unique(aspect), c(1.10, 1.56)) + + fp <- fly_footprint(p) + expect_false(any(sf::st_is_empty(sf::st_geometry(fp)))) + expect_true(all(fp$footprint_terrain == "gsd_scaled")) + + # Non-square, which is the property `fly_georef()` keys on — a footprint rotated onto + # a flight line running due east is still axis-aligned, so axis-alignment is a proxy + # that misses it. + expect_false(any(fly_is_square(fp))) +}) + + +test_that("an inferred-format frame IS sized when a dem is supplied", { + # The documentation says a frame carrying no calibration can only be sized through a + # DEM, because pixel count spreads 32-83% at a given focal length. That route has to + # actually be reachable: keying eligibility on the nominal-scale half-side made it + # unreachable for every camera-table row at once — `width_in` is NA there by + # definition — so these came back empty with a DEM supplied, and silently, because + # they are excluded from the unknown-format warning too. + skip_if_no_terra() + p <- digital_fixture() + p$camera_calibration_url <- NA # force the focal-length fallback + p$ground_sample_distance <- NA # remove the GSD route entirely + + fp <- suppressWarnings(fly_footprint(p, dem = testdata_path("dem.tif"))) + + expect_equal(fp$footprint_basis, rep("inferred_format", nrow(p))) + expect_false(any(sf::st_is_empty(sf::st_geometry(fp)))) + expect_equal(fp$footprint_terrain, rep("dem_agl", nrow(p))) + expect_true(all(fp$height_agl > 0)) +}) + + +test_that("a calibrated frame with no usable gsd is sized by the dem instead", { + # `ground_sample_distance` is 0 on every frame of some digital rolls, so this is a + # real population rather than a constructed one. + skip_if_no_terra() + p <- digital_fixture() + p$ground_sample_distance <- 0 + + fp <- suppressWarnings(fly_footprint(p, dem = testdata_path("dem.tif"))) + expect_false(any(sf::st_is_empty(sf::st_geometry(fp)))) + expect_equal(fp$footprint_terrain, rep("dem_agl", nrow(p))) + + # Without a DEM the same frames have no route at all, and must say so. + expect_warning(fly_footprint(p), "no way to size it") +}) + + +test_that("an NA film_roll does not abort fly_footprint", { + # `fly_bearing()` compared rolls with a bare `==`, so an NA gave NA to an `if` and + # aborted. Latent until non-square footprints made `fly_footprint()` call it. + p <- digital_fixture() + p$film_roll[2] <- NA + + expect_no_error(fp <- suppressWarnings(fly_footprint(p))) + expect_equal(nrow(fp), nrow(p)) +}) + + +test_that("a frame naming a withheld calibration says so", { + # The refusal is computed in the resolver; without carrying it through, the frame is + # indistinguishable from one whose media was simply unknown. + ex <- fly_camera_excluded() + withheld <- ex$key[ex$key_type == "calib_file"][1] + + p <- digital_fixture()[1, ] + p$camera_calibration_url <- paste0( + "https://openmaps.gov.bc.ca/thumbs/calib_report_zips/", withheld, ".zip" + ) + fp <- suppressWarnings(fly_footprint(p)) + + expect_equal(fp$width_source, paste0("withheld:", withheld)) + expect_true(sf::st_is_empty(sf::st_geometry(fp))) + expect_true(is.na(fp$footprint_terrain)) +}) + + +test_that("a DEM-sized footprint is rotated onto the flight line too", { + # The rotation tests above run without a DEM, so they only reach the GSD route. The + # DEM route fills the half-dimensions *after* the point where rotation is decided, so + # keying that decision on them left DEM-sized frames axis-aligned while `fly_bearing()` + # had a good azimuth for them — a 1.76:1 footprint transposed 90 degrees, with the area + # unchanged and the ground covered wrong. + skip_if_no_terra() + p <- digital_fixture() + p$ground_sample_distance <- NA # force the DEM route + bearing <- fly_bearing(p)$bearing + + fp <- sf::st_transform( + suppressWarnings(fly_footprint(p, dem = testdata_path("dem.tif"))), 3005 + ) + expect_equal(fp$footprint_terrain[1:4], rep("dem_agl", 4)) + + along_azimuth <- function(i) { + ring <- sf::st_coordinates(sf::st_geometry(fp)[[i]])[1:4, 1:2] + v <- ring[3, ] - ring[2, ] + unname((atan2(v[1], v[2]) * 180 / pi) %% 360) + } + for (i in 1:4) expect_equal(along_azimuth(i), bearing[i], tolerance = 0.01) +}) + + +test_that("rotation does not depend on what else is in the batch", { + # Rotation was decided from a vector only some sizing routes populate, so whether a + # given row came out rotated depended on whether *other* rows carried a GSD. + skip_if_no_terra() + dem <- testdata_path("dem.tif") + mixed <- digital_fixture() # rows 1-5 carry a GSD + no_gsd <- digital_fixture() + no_gsd$ground_sample_distance <- NA + + ring1 <- function(p) { + fp <- sf::st_transform(suppressWarnings(fly_footprint(p, dem = dem)), 3005) + r <- sf::st_coordinates(sf::st_geometry(fp)[[1]])[1:4, 1:2] + v <- r[3, ] - r[2, ] + unname((atan2(v[1], v[2]) * 180 / pi) %% 360) + } + expect_equal(ring1(mixed), ring1(no_gsd), tolerance = 0.01) +}) + + +test_that("an unsizable film frame is told about scale, not about GSD", { + # A film frame reaches "resolved but unsized" through an unparseable `scale`. Telling + # its owner to supply a ground sample distance points at the wrong column entirely. + photos <- mixed_media_fixture()[1:2, ] + photos$scale <- c("not a scale", "1:12000") + + expect_warning(fp <- fly_footprint(photos), "no usable `scale`") + expect_true(sf::st_is_empty(sf::st_geometry(fp)[1])) + expect_false(sf::st_is_empty(sf::st_geometry(fp)[2])) +}) diff --git a/tests/testthat/test-fly_footprint_invariants.R b/tests/testthat/test-fly_footprint_invariants.R new file mode 100644 index 0000000..500b112 --- /dev/null +++ b/tests/testthat/test-fly_footprint_invariants.R @@ -0,0 +1,71 @@ +# Structural invariants of the reporting columns, swept across every input shape a +# caller can present. +# +# A global invariant beats more examples here: it cannot be gamed by fixture choice, and +# the reporting surface is exactly where this package has been bitten — #30 claimed a +# terrain treatment for a frame with no geometry, #35 dropped the columns entirely on a +# tibble, and #32 shipped a basis and a width_source alongside an empty geometry with +# nothing to say so. + +test_that("footprint_terrain is NA exactly where the geometry is empty", { + skip_if_no_terra() + cases <- footprint_cases() + for (nm in names(cases)) { + cs <- cases[[nm]] + fp <- suppressWarnings(fly_footprint(cs$x, dem = cs$dem)) + empty <- sf::st_is_empty(sf::st_geometry(fp)) + + # Both directions. A terrain value on an empty geometry claims a treatment for a + # frame that was never placed; an NA on a real footprint hides how it was sized. + expect_identical(is.na(fp$footprint_terrain), empty, info = nm) + expect_true(all(is.na(fp$height_agl[empty])), info = nm) + expect_true(all(is.na(fp$dem_coverage[empty])), info = nm) + } +}) + + +test_that("every frame gets a basis, and the reporting columns keep their types", { + skip_if_no_terra() + cases <- footprint_cases() + for (nm in names(cases)) { + cs <- cases[[nm]] + fp <- suppressWarnings(fly_footprint(cs$x, dem = cs$dem)) + + expect_true(all(!is.na(fp$footprint_basis)), info = nm) + # `ifelse(logical(0), ...)` returns `logical(0)`, so an empty result would report a + # character column as logical and fail to bind to a populated one. + expect_type(fp$footprint_basis, "character") + expect_type(fp$footprint_terrain, "character") + expect_type(fp$width_source, "character") + expect_type(fp$height_agl, "double") + expect_type(fp$dem_coverage, "double") + expect_identical(nrow(fp), nrow(cs$x), info = nm) + } +}) + + +test_that("a sized footprint is a closed rectangle of positive area", { + skip_if_no_terra() + cases <- footprint_cases() + for (nm in names(cases)) { + cs <- cases[[nm]] + fp <- suppressWarnings(fly_footprint(cs$x, dem = cs$dem)) + g <- sf::st_geometry(sf::st_transform(fp, 3005)) + sized <- !sf::st_is_empty(g) + if (!any(sized)) next + + # Catches the degenerate five-identical-vertex polygon a zero half-dimension builds: + # `st_is_empty()` reports FALSE for it, so it would pass every emptiness check while + # covering nothing. + expect_true(all(as.numeric(sf::st_area(g[sized])) > 0), info = nm) + for (i in which(sized)) { + xy <- sf::st_coordinates(g[[i]])[, 1:2, drop = FALSE] + expect_identical(nrow(xy), 5L, info = nm) + expect_equal(xy[1, ], xy[5, ], info = nm) + # Opposite edges equal: still a rectangle after any rotation. + d <- sqrt(rowSums(diff(xy)^2)) + expect_equal(d[1], d[3], info = nm) + expect_equal(d[2], d[4], info = nm) + } + } +}) diff --git a/vignettes/airphoto-selection.Rmd b/vignettes/airphoto-selection.Rmd index 4043bee..4d4e96f 100644 --- a/vignettes/airphoto-selection.Rmd +++ b/vignettes/airphoto-selection.Rmd @@ -51,12 +51,33 @@ table(footprints$footprint_basis) sized <- footprints[footprints$footprint_basis != "unknown_format", ] ``` -Digital frames resolve to `"unknown_format"` because the sensor width they need -is not in the centroid metadata. Supply it yourself if you know the camera: +A digital frame is sized from its camera instead. Sensor dimensions are not in +the centroid metadata, but they are recoverable from the calibration report each +frame links to through `camera_calibration_url`, and `fly` ships them for the +cameras in the provincial record. The footprint is then +`pixel count x ground_sample_distance`, which uses neither the reported `scale` +nor a DEM: + +```r +digital <- st_read(system.file("testdata/photo_centroids_digital.gpkg", package = "fly")) +dfp <- fly_footprint(digital) +table(dfp$width_source) +``` + +Two things differ from the film case. **`scale` is not used**: for a digital +frame the catalogue's `SCALE` is a derived nominal figure rather than the true +image scale, and sizing from it gives about a third of the real width. +**Footprints are not square** — sensors run from 1.10:1 to 1.80:1 — so they are +rotated onto the flight line using `fly_bearing()`. + +`width_source` names the calibration file behind each footprint, or the +fallback rule where a frame carries no calibration; those are marked +`"inferred_format"` in `footprint_basis`. `format_size` still overrides +everything if you know your camera: `fly_footprint(photos, format_size = c("Digital - Colour" = 3.54))`. -The photos used throughout this vignette are 1968 film, so every footprint below -is sized from the 9-inch negative. +The photos used throughout the rest of this vignette are 1968 film, so every +footprint below is sized from the 9-inch negative. By default footprints are sized from the reported scale, which assumes flat ground at whatever elevation that scale was computed for. Supplying a DEM