Measure coastal frames' sea surface: an elevation, not an error (#65)

CLAUDE.md

@@ -54,6 +54,12 @@ the note and NEWS quote. Ships `inst/extdata/flying_height_sweep.csv` and
 `flying_height_population.csv`, which `test-fly_footprint_height.R` reads so the three constants
 are checked against the data. Use a PSOCK cluster for remote DEM reads — `mclapply()` forks
 abort on GDAL's curl handles on macOS while the wrapper exits 0
+- `data-raw/dem_measure-coastal_water.R` — fly#65: whether a coastal frame's sea is an error.
+Ray-casts the true footprint of sampled coastal and inland frames onto MRDEM and scores the
+package's land-and-sea mean, the land-only mean and per-side sizing on area and on the land
+edge. Ships `inst/extdata/dem_coastal_frames.csv`, `_population.csv` and `_sites.csv`, which
+`test-fly_footprint_coastal.R` recomputes. `FLY_COASTAL_SMOKE=1` runs a handful of frames
+into a separate cache and writes nothing
 - `data-raw/height_calibrate-lower_tail_rolls.R` — settles the lower tail of `flying_height`
 per roll (fly#60) with adjacent-frame spacing and the hand transcription of the province's
 logbook scans in `data-raw/flying_height_logbooks.csv`, which is an input and is never
@@ -201,10 +207,37 @@ rather than a property of this code.
   digital, 0 of 2,975 random controls. The issue's own "18 of 416" does not reproduce; a
   DEM cropped to that run's AOI is the likely cause and is recorded as inference. So the
   frequency is a property of how the caller crops, not of the catalogue. A near-zero ocean
-  surface (not nodata, and with no exact zeros) is a coverage-1 failure the sweep cannot
-  generate and `dem_coverage` cannot see — out of scope, stated as a bound. Read
+  surface (not nodata, and with no exact zeros) is a coverage-1 case the sweep cannot
+  generate — recorded here as a failure, and fly#65 found it is not one (below). Read
   `inst/notes/terrain-correction.md`
 
+- **A coastal frame's sea is an elevation, and averaging it in stays** (fly#65, no code
+  change) — MRDEM carries near-shore sea at ~0.14 m, which fly#58 called a mean "dragged
+  toward sea level". Over water the imaged surface *is* the sea. Against a **ray-cast of the
+  true footprint** on bare earth, the package's land-and-sea mean (W) beats the land-only
+  mean (L) on area, and L beats W on where the land edge falls. The pre-registered test
+  asked, pooled, whether the sea makes W worse than inland, and it passed, so nothing
+  changed. Two qualifications are recorded, not acted on. Matched for relief (the coastal
+  sample is flatter), the land-edge pass disappears, so the sea does worsen W's land edge,
+  though not its area. And a canopy on the land can reverse the area verdict: first-order
+  only, a DTM effect everywhere, fly#80.
+
+  **Three things are load-bearing.**
+  - The instrument changed before any frame was measured. Adjacent-frame spacing measures
+    when the shutter fired, not what the photo covered, so it is secondary. The ray-cast's
+    flat and step synthetic controls converge.
+  - "Which is right" has two answers, one per consumer: area for
+    `fly_coverage()`/`fly_overlap()`, the land edge for `fly_filter()`/`fly_select()`
+    against a land AOI.
+  - The prose was the defect class. Three review rounds each found published causes,
+    denominators or comparisons written from a story rather than read off a producer line,
+    including inside the previous round's fix. So the note's tables are now rebuilt row by
+    row by the test and counted.
+
+  LidarBC was mostly nodata over sea in the two Howe Sound tiles probed, so the same frame is
+  sized ≈L there and, under 0.95 coverage, warns about interior nodata. Do not "fix" this by
+  masking near-zero cells: in the Roberts Bank window 14.8% of delta land reads under 1 m. Read `inst/notes/terrain-correction.md`
+
 - **Terrain error is a datum offset, not slope** (v0.5.0, #9) — `FLYING_HEIGHT` is metres **above sea level**,
 and reported scale is referenced to an elevation above the ground the photos cover, so it understates footprint
 area by a median 13.8% *always in the same direction*. `fly_footprint(dem =)` sizes each frame from

NEWS.md

@@ -1,5 +1,9 @@
 # fly (development version)
 
+- **A coastal frame's sea surface is an elevation, not an error, and `fly_footprint()` is unchanged** ([#65](https://github.com/NewGraphEnvironment/fly/issues/65)). v0.14.0's notes said MRDEM's near-zero surface over near-shore sea drags a coastal frame's mean "toward sea level". MRDEM carries that sea at about 0.14 m, not as nodata, and over water it is the surface the photo images. Against a ray-cast of the true footprint on MRDEM's bare earth, over 1,242 coastal frames, averaging land and sea (what the package does) is closer on **area** than the land-only mean the issue proposed: 0.34% against 0.67% of width, median. The land-only mean places the **land edge** better, misplacing 3.67% of the frame's land against 4.43%. The rule fixed before the run asked whether the package's 95th-percentile error on coastal frames exceeds inland. Pooled, it does not: 2.16% against 3.18% on area, and 15.05% against 15.87% on the land edge. So no remedy was triggered. Coastal frames are 95,222 of 1,437,147 DEM-eligible film frames (6.63%).
+- **Two qualifications, recorded rather than acted on.** First, the land-edge pass is an artefact of relief: the coastal sample is flatter, and matched for relief the inland land-edge 95th is 15.06% against 15.05% coastal. Within four of five relief bins coastal is worse by 2 to 9 points, so the sea does make the package's land edge worse. It does not make its area worse, even matched. Second, a canopy can reverse the area verdict. To first order, a uniform 30 m canopy on the land ties the two means and 60 m favours the land-only one. That is a property of sizing from a bare-earth DEM anywhere, filed as [#80](https://github.com/NewGraphEnvironment/fly/issues/80).
+- **LidarBC was mostly nodata over sea in the two tiles probed**, so the same coastal frame is sized nearer the land-only mean there. Where that leaves it under 95% covered, it warns about nodata inside the DEM's extent. `data-raw/dem_measure-coastal_water.R` prints every figure above, and ships `inst/extdata/dem_coastal_frames.csv`, `dem_coastal_population.csv` and `dem_coastal_sites.csv`, which the test suite rebuilds the note's tables from. `inst/notes/terrain-correction.md` explains why the instrument first planned, adjacent-frame spacing, was dropped before any frame was measured.
+
 ## 0.19.0 (2026-09-28)
 
 - **Breaking: grayscale georeferenced output is now 2 bands, Gray + Alpha, where it was 1 band with `NoData = 0`** ([#56](https://github.com/NewGraphEnvironment/fly/issues/56)). `fly_georef()` warps every frame with `-dstalpha`, so RGB stays 4 bands and, from a grayscale or RGB source, no output carries a NoData value. A consumer that reads band 1 alone, or treats 0 as fill, must now read the last band: 0 is fill, and on these outputs the opaque value is the datatype's maximum (255 on Byte, 65535 on UInt16). Under `-dstnodata 0`, GDAL rewrote genuine black inside the frame as 1 so it would not read as fill (or deleted it where `srcnodata = "0"` was also given). Measured on the output of the 182 grayscale frames the mask was calibrated on, it hit **161 frames** warped axis-aligned (every frame with no flight bearing): typically a few dozen pixels, up to 3.6% of a frame on roll bcb94081 (1250 x 1250 thumbnails on square footprints; an anisotropic scan resamples even when axis-aligned). On a footprint rotated 30 degrees it hit 41 frames, at most 0.26%, because resampling mixes a lone 0 with its neighbours. The Windows CI runner also gave a masked grayscale frame 2 bands where ubuntu and macOS gave 1 ([#68](https://github.com/NewGraphEnvironment/fly/issues/68)). Band count is now one number, measured on macOS with sf's GDAL 3.8.5 and asserted by the test suite on all three CI platforms. This supersedes 0.11.0's "nothing downstream changes band count". `stac_airphoto_bc` needs a COG-writer change and a republish of grayscale frames ([stac_airphoto_bc#36](https://github.com/NewGraphEnvironment/stac_airphoto_bc/issues/36)).

R/fly_footprint.R

@@ -808,6 +808,15 @@ fly_is_square <- function(footprints) {
 #' Falling back to nominal scale is **not** the better answer at any coverage: the
 #' DEM route beats it in the median even below 20% covered.
 #'
+#' **A coastal frame is sized from land and sea together.** MRDEM carries
+#' near-shore sea at about 0.14 m, not nodata, so the mean takes in the sea
+#' surface the photo images and `dem_coverage` stays near 1. Against a ray-cast
+#' of the true footprint on bare earth, that is closer on area than the
+#' land-only mean, and the land-only mean places the land edge better. A canopy
+#' can reverse the first, and matched for relief the sea does make the land edge
+#' worse than it is inland. LidarBC was mostly nodata over sea where probed, so
+#' there the same frame is sized nearer the land-only mean. See `inst/notes/terrain-correction.md` (fly#65).
+#'
 #' Buffer past the **corner** of the widest footprint, not its half-side: the
 #' far point of a square is `half_side * sqrt(2)`, which at 1:31680 is 5.1 km
 #' rather than 3.6 km. Allow more again for the correction itself, which