Stop dft_rast_classify() modifying the caller's raster (#89)

CLAUDE.md

@@ -142,7 +142,7 @@ Shape: 14651 x 11552 at 10 m (EPSG:32609), 169M cells, 97.7% `NA` outside the fl
 
 **Why:** the first `dft_transition_artifact()` held one in-memory 169M-cell raster per class per intermediate and was killed for memory here at eight classes, with 125 unit assertions green (#44). Only a run at scale finds that class of bug.
 
-**How to apply:** run any new or changed raster-pipeline function on this pair with an RSS sampler (`ps -o rss= -p $PID` every 2 s) before opening the PR, and record the numbers in the PR body. Write terra intermediates with `filename =` (see `code-check-spatial.md`).
+**How to apply:** run any new or changed raster-pipeline function on this pair with an RSS sampler (`ps -o rss= -p $PID` every 2 s) before opening the PR, and record the numbers in the PR body. For a call shorter than about a minute, use `/usr/bin/time -l R -f script.R` instead. It reports the kernel's true peak, while a 2 s sampler lands on arbitrary instants: #89 read 0.26 to 2.31 GiB for a run that peaked at 4.21 GiB. Write terra intermediates with `filename =` (see `code-check-spatial.md`).
 
 ### SRED cross-reference
 

DESCRIPTION

@@ -1,6 +1,6 @@
 Package: drift
 Title: Detecting Riparian and Inland Floodplain Transitions
-Version: 0.19.0
+Version: 0.19.1
 Date: 2026-09-29
 Authors@R: c(
     person("Allan", "Irvine", , "al@newgraphenvironment.com", role = c("aut", "cre"),

NEWS.md

@@ -1,3 +1,8 @@
+# drift 0.19.1
+
+- **`dft_rast_classify()` no longer modifies the raster it is given (#89).** It called the in-place `terra::set.cats()` on the caller's own object, so after classifying, the original was a factor named `class_name`. This affected file-backed, in-memory and list-element inputs alike; only a `remap =` that matched a class was spared, because that path already builds a new raster. The fix swaps the order of the two setters. `coltab<-` runs first and makes a copy, and the levels are then set on that copy. The output is unchanged: `identical()` `cats()` and `coltab()` against 0.19.0. The obvious alternative, `terra::deepcopy()` first, adds a second full copy. On BULK in memory (169M cells) the reorder peaks at 4.21 GiB, the same as 0.19.0, while `deepcopy()` peaks at 5.46 GiB. A test pins that the caller is unmodified, and it goes red if terra ever makes `coltab<-` work in place.
+- **Found on the way:** given a raster that is already a factor, such as the published floodplain rasters, `dft_rast_classify()` returns empty levels and no colours (#91). This is not new in 0.19.1.
+
 # drift 0.19.0
 
 - **Accuracy and error-adjusted area (#81).** drift reported mapped area only. Mapped area is biased wherever the map is wrong, and badly biased for change maps, because a single wrong label on either date manufactures a transition. Four new functions implement the standard remedy (Olofsson et al. 2014): a stratified random sample of points, reference labels, and estimators that weight the labels back to the population.

R/dft_rast_classify.R

@@ -44,14 +44,18 @@ dft_rast_classify <- function(x,
   present_codes <- terra::unique(x)[, 1]
   ct <- class_table[class_table$code %in% present_codes, ]
 
-  # Set factor levels (use set.cats to avoid namespace issues with levels<-)
-  lvl_df <- data.frame(id = ct$code, class_name = ct$class_name)
-  terra::set.cats(x, layer = 1, value = lvl_df)
-
-  # Set color table
+  # The order is load-bearing (#89). `x` is the caller's own raster, and
+  # set.cats() works in place, so it must only ever reach a copy. `coltab<-`
+  # deep-copies before setting (the fact strip_copy() in dft_rast_break_class.R
+  # also relies on), so colours first and levels on that copy: one copy, where
+  # an explicit deepcopy() before `coltab<-` would make two. Should terra make
+  # `coltab<-` in place, the caller-unmodified test goes red.
   coltab_df <- data.frame(value = ct$code, col = ct$color)
   terra::coltab(x) <- coltab_df
 
+  lvl_df <- data.frame(id = ct$code, class_name = ct$class_name)
+  terra::set.cats(x, layer = 1, value = lvl_df)
+
   x
 }
 

planning/archive/2026-09-issue-89-classify-mutation/README.md

@@ -0,0 +1,24 @@
+## Outcome
+
+`dft_rast_classify()` modified the raster it was given: it called the in-place `terra::set.cats()` on the caller's own object, so the original came back as a factor named `class_name`. That held for file-backed, in-memory and list-element inputs, and for a `remap =` that matched no class. The fix swaps the two setters. `coltab<-` deep-copies as its first statement at every supported terra (checked at the 1.8-10 floor and at 1.9.50), so the levels are now set on that copy. The output is `identical()` to 0.19.0, and no second copy is made. A new test pins the caller as unmodified and goes red on the old order with 6 failures. Code-check ran three rounds, all clean. The BULK run found a separate, pre-existing defect: a raster that is already a factor gets empty levels, because `terra::unique()` returns labels. It is filed as drift#91.
+
+## Measurement
+
+BULK `classified_2017.tif` (169,248,352 cells), terra 1.9.50. Kernel max RSS from `/usr/bin/time -l`, two reps, identical to 0.01 GiB:
+
+| variant | input | max RSS (GiB) | caller untouched |
+|---|---|---|---|
+| 0.19.0 order | in memory | 4.21 | no |
+| reorder (shipped) | in memory | 4.21 | yes |
+| `deepcopy()` then 0.19.0 order (the issue's proposal) | in memory | 5.46 | yes |
+| 0.19.0 / reorder | file-backed | 1.05 / 1.05 | no / yes |
+
+The reorder costs nothing. `deepcopy()` would add 1.25 GiB, one full double copy, per year classified in memory. That number decided the design.
+
+**Wrong turn, kept.** The first measurement used a 2 s `ps` sampler and read 0.26 GiB for the in-memory run. I blamed sampling the `Rscript` PID instead of R, then retracted that: `Rscript` execs in place and the PID is R's own. The real cause was a 2 s sampler on a call whose peak lasts under a second; its readings for one variant ranged from 0.26 to 2.31 GiB. The same run also briefly suggested the fix failed at scale. That came from a names-only check against an input that was already a factor, and a `cats()`/`coltab()` snapshot resolved it. The resulting CLAUDE.md note says to use `/usr/bin/time -l` for sub-minute scale checks.
+
+## Evidence
+
+The scripts were one-off and ran in the session scratchpad; the numbers above and in `findings.md` are the record.
+
+Closed by: PR for #89 (branch `89-dft-rast-classify-mutates-the-caller-s-r`)

planning/archive/2026-09-issue-89-classify-mutation/findings.md

@@ -0,0 +1,110 @@
+# Findings — dft_rast_classify() mutates the caller's raster in place (set.cats on the input) (#89)
+
+## Issue context
+
+## Problem
+
+`dft_rast_classify()` calls `terra::set.cats(x, ...)` on the raster it was given (`R/dft_rast_classify.R:49`). `set.cats()` works in place, so the caller's object turns into a factor named `class_name`. When `remap = NULL`, `x` is the caller's own SpatRaster, not a copy.
+
+```r
+r <- terra::rast(system.file("extdata", "example_2017.tif", package = "drift"))
+names(r); terra::is.factor(r)        # "data" FALSE
+cl <- dft_rast_classify(list("2017" = r), source = "io-lulc")
+names(r); terra::is.factor(r)        # "class_name" TRUE   <- the input changed
+```
+
+Found in #81. A test built `c(r17, r23)` from the raw tiles after an earlier test had classified `r17`, and the layer names were no longer the file's. Any caller that classifies a raster and then reuses the original for anything that reads names or levels is affected, for example stacking it or passing it back into a function that validates it.
+
+CLAUDE.md's spatial conventions already name the trap: `set.cats()` "mutates whatever raster it is given, so use it on a copy you own".
+
+## Fix
+
+Take a copy before `set.cats()` (`x <- terra::deepcopy(x)`, or `levels<-`, which copies first), and add a test asserting that the caller's raster keeps its names and stays a non-factor. Check the other `set.cats()` / `set.names()` / `set.values()` call sites in `R/` for the same shape.
+
+Not fixed in #81: it is outside that issue, and #81's tests now avoid relying on it.
+
+## Plan-mode exploration (2026-09-29)
+
+### Probe and call-site sweep
+
+`dft_rast_classify()` calls `terra::set.cats(x, ...)` (`R/dft_rast_classify.R:49`) on the raster it
+was handed. `set.cats()` works in place, so with `remap = NULL` the caller's object becomes a factor
+named `class_name`. Found in #81; NEWS 0.19.0 already names it as "found on the way".
+
+**Probed (terra 1.9.50, in memory, source tree):**
+
+| input | caller after the call |
+|---|---|
+| file-backed single raster | `class_name`, factor — **mutated** |
+| in-memory raster (`r * 1L`) | **mutated** |
+| element of a named list | **mutated** |
+| `remap =` supplied | untouched (`terra::classify()` returns a new raster) |
+
+**Fix found by probe, zero extra copy:** swap the two setters — `coltab<-` first (it deep-copies
+unconditionally, the same fact `strip_copy()` in `R/dft_rast_break_class.R:350` already relies on
+and documents), then `set.cats()` on that copy. Output is identical to today's
+(`identical()` on `cats()` and `coltab()`), and the caller comes back `data` / non-factor.
+
+Chosen over the issue's `x <- terra::deepcopy(x)`: that adds a second full copy, because `coltab<-`
+still copies afterwards. For an in-memory BULK year (what `dft_stac_fetch()` returns when it fits,
+~1.5 GB) that is +1.5 GB transient per year. The reorder's reliance on `coltab<-` copying is guarded
+by the new test, which goes red if terra ever makes it in place.
+
+**Other call sites checked — all act on rasters the function itself created, none need a change:**
+`dft_rast_transition.R` (141/158/177: `r_trans`/`r_removed` from arithmetic / `ifel`),
+`dft_rast_consensus.R:92` (`r_out <- terra::rast(x[[1]])`, an empty template),
+`dft_rast_break_category.R:153` (`app()` output), `dft_rast_break_class.R` (207: stack from
+`rast(list)`, already guarded by the caller-unmutated test at `test-dft_rast_break_class.R:329`;
+267/302: `out[["transition"]]` subset; 362: `strip_copy()` after `coltab<-`),
+`dft_transition_artifact.R:308` (explicit `deepcopy`), `dft_stac_cube.R:287` (`names<-` on its own
+`stk`).
+
+## BULK scale check (2026-09-29, Phase 3)
+
+`bulk_co_ff04/classified_2017.tif` (14651 x 11552, 169,248,352 cells), terra 1.9.50, 64 GB machine.
+Each variant a fresh `R -f` process under `/usr/bin/time -l` (kernel max RSS), two reps each, identical
+to 0.01 GiB. "mem" is `r * 1L` (in memory, the shape `dft_stac_fetch()` returns when it fits). Script:
+scratchpad `bulk_classify2.R`; `main` is `git show origin/main:R/dft_rast_classify.R` sourced into
+the namespace; `deepcopy` is main's order with `x <- terra::deepcopy(x)` before it (the issue's fix).
+
+| variant | input | max RSS (GiB) | classify (s) | caller untouched |
+|---|---|---|---|---|
+| main | file-backed | 1.05 | 0.8 | no |
+| branch | file-backed | 1.05 | 0.8 | yes |
+| main | in memory | 4.21 | 0.8 | no |
+| branch | in memory | 4.21 | 0.9 | yes |
+| deepcopy | in memory | 5.46 | 1.2 | yes |
+
+The reorder costs nothing; `deepcopy()` costs +1.25 GiB, one 169M-cell double copy, per year classified.
+Branch and main outputs are `identical()` in `cats()` and `coltab()` in both modes.
+
+"Caller untouched" for file-backed is from a `cats()`/`coltab()`/`names()`/`is.factor()` snapshot
+compared before and after. The published raster is ALREADY a factor (RAT with `class_name`, palette),
+so a names-only check reads it as mutated whatever happens — the first run of the script did exactly
+that and briefly looked like the fix failing at scale.
+
+### Wrong turns (kept as evidence)
+
+- **Diagnosed, then retracted: "`Rscript` forks `R`, so the sampler watched the wrapper."** A 2 s
+  `ps` sampler gave 0.26 GiB for a 169M-cell in-memory raster, and I blamed the PID. Wrong: measured
+  afterwards, `Rscript` execs in place. The backgrounded PID's `comm` is `.../bin/exec/R`, it has no
+  children, and it equals R's own `Sys.getpid()`. Switching to `R -f` changed nothing that mattered.
+- **The actual cause: a 2 s sampler on a ~1 s call.** Samples land at arbitrary instants, and the
+  classify peak lasts well under the interval. Readings for the same variant ranged from 0.26 to
+  2.31 GiB, and the deepcopy variant read 0.25 GiB, below the size of its own in-memory input.
+  `/usr/bin/time -l` reports the kernel's true peak and gave identical numbers across two reps. The
+  CLAUDE.md "RSS every 2 s" recipe fits the minutes-long pipeline functions it was written for, not
+  a sub-second call.
+
+### Found on the way: factor input gives empty levels (drift#91)
+
+On the published (already-factor) raster the output had no levels and no colours on main and on the
+branch. `terra::unique(x)[, 1]` returns labels on a factor, so no code matches `class_table$code`.
+Pre-existing, out of scope for #89; filed as drift#91. Round-1 review found it independently.
+
+## Errors Encountered
+
+| Error | Resolution |
+|-------|------------|
+| BULK RSS 0.26 GiB for an in-memory 169M-cell raster | Not the PID (`Rscript` execs in place, verified). A 2 s sampler misses a sub-second peak; use `/usr/bin/time -l` |
+| `cmd && V=... && R ... & PID=$!` lost `$PID` | `&` backgrounds the whole `&&` list; split setup and the backgrounded command |

planning/archive/2026-09-issue-89-classify-mutation/progress.md

@@ -0,0 +1,14 @@
+# Progress — dft_rast_classify() mutates the caller's raster in place (set.cats on the input) (#89)
+
+## Session 2026-09-29
+
+- Plan-mode exploration — phases approved by user
+- Created branch `89-dft-rast-classify-mutates-the-caller-s-r` off main
+- Scaffolded PWF baseline from issue #89 with approved phases
+- Next: start Phase 1
+- Phase 1: caller-unmodified test added; red on current code for all three inputs (file-backed, in-memory, list element), 6 failures
+- Phase 2: `coltab<-` moved before `set.cats()`; classify tests 32 pass; the Phase 1 red run against the original ordering is the restore-the-bug check
+- Full suite `[ FAIL 0 | WARN 0 | SKIP 15 | PASS 1318 ]` (49 s); changed lines lint clean (two pre-existing indentation lints in the remap tests, untouched); `document()` no change
+- Code-check: three rounds, all clean; round 1 flagged the stale #89 workaround comment in `test-dft_accuracy_sample.R` (removed). Committed `fe93a82`
+- Phase 3: BULK scale check. In-memory peak 4.21 GiB on main and on the branch, 5.46 GiB with the `deepcopy()` alternative; file-backed 1.05 GiB either way. Outputs are identical. Measurement wrong turn recorded in findings: I blamed the `Rscript` PID, then retracted it after verifying that `Rscript` execs in place. The real cause was a 2 s sampler on a ~1 s call; `/usr/bin/time -l` is the instrument
+- Filed drift#91: factor input returns empty levels (pre-existing, found at scale)

planning/archive/2026-09-issue-89-classify-mutation/review-round1.md

@@ -0,0 +1,39 @@
+# Review round 1 — #89 (`dft_rast_classify()` mutates the caller's raster)
+
+## Clean
+
+No issues found in the diff.
+
+## What was verified (terra 1.9.50, probes run from the scratchpad, repo untouched)
+
+- **`coltab<-` copies unconditionally.** Printed the `coltab<-` SpatRaster method: its first
+  statement is `x@pntr <- x@pntr$deepcopy()`, before any branch on `value` or `layer`, so it copies
+  for a data.frame value (this call) as well as for `NULL` (the `strip_copy()` call). `set.cats()`
+  therefore only ever reaches the copy.
+- **New ordering matches the old one.** Ran the new function against a copy of the old
+  set.cats-then-coltab body over nine inputs: file-backed; in-memory (`* 1L`); all-NA (empty
+  `ct`); no code in the class table (`* 0L + 99L`, empty `ct`); a single code; a caller already a
+  factor; a caller already carrying a coltab; an already-classified raster; and a two-layer stack.
+  In all nine `cats()`, `coltab()`, `names()`, `is.factor()`, `has.colors()` and the values were
+  `identical()` between the two orderings. `set.cats()` after `coltab<-` keeps the colour table.
+  The empty-`ct` case neither errors nor behaves differently in either order.
+- **The caller is untouched in all nine cases**, including the multi-layer stack. `names`,
+  `is.factor`, `has.colors`, `cats` and `coltab` were `identical()` before and after the call. The
+  remap path is also unaffected: `terra::classify()` already returns a new raster, and the remapped
+  output still carries both levels and colours.
+- **Restore the bug:** evaluated the new `test_that()` block against the old ordering and it fails
+  with 6 failures, 13 successes, which matches the author's number. The guard fires. Unlike a
+  `levels<-`/`coltab<-` strip (checklist, "copy before they strip"), this test *can* fail, because
+  the old code reached the caller through the in-place `set.cats()`.
+- Every in-package call site uses the return value. Nothing depended on the in-place mutation.
+
+## Observations outside this diff (not findings)
+
+- `tests/testthat/test-dft_accuracy_sample.R:137`: the comment "fresh tiles:
+  dft_rast_classify() mutates its input in place (#89)" is stale once this merges. The test still
+  passes, so this is housekeeping only.
+- A pre-existing issue, the same in both orderings: when the input is **already a factor**,
+  `terra::unique(x)[, 1]` returns the labels, not the codes. `ct` then comes back empty, and the
+  output has no colours and no usable levels. This happens when an already-classified raster is
+  re-classified, or when the caller's raster is a factor with other labels. The diff does not
+  cause it; it may be worth its own issue.