Problem
fl_cost_distance() treats every zero-friction cell as a stream seed, not just stream cells.
# R/fl_cost_distance.R:51-52
cost <- terra::ifel(!is.na(streams), 0, friction)
out <- terra::costDist(cost, target = 0)
costDist(target = 0) finds all cells equal to 0. Stream cells are set to 0 deliberately — but so
is any cell whose friction (percent slope) is already exactly 0. The roxygen says otherwise:
Stream cells are seed points with cost zero. […] streams identifies seed cells (cost = 0).
Demonstration
Synthetic 50x50 grid, friction 10% everywhere, one stream cell, plus a flat (0%) patch 380 m away:
cost at flat-patch centre : 0 (should be non-zero)
cost one cell beyond the patch : 50 (should be non-zero)
cost at equal distance, no patch : 1450
cells with cost 0 : 122 (stream cells: 1)
The bundled test DEM hides this — terrain()-derived slope has zero cells exactly equal to 0
(min 1.42e-14), so no existing test can catch it. Production DEMs do produce exact zeros:
integer-metre DEMs, hydro-flattened lake surfaces, and void-filled plateaus all quantize to flat.
Why it matters more now
Under fl_valley_confine() alone, a spurious flat seed lowers cost slightly in ground that was
probably valley anyway. Under per-group attribution (#40) it is worse: a flat patch becomes a free
zero-cost source for whichever group's corridor contains it, so that group's cost mask can spread
across ground its own streams never reach. Cost then stops discriminating between groups exactly
where floodplains are — on flat ground — and membership degenerates toward the distance buffer.
Fix (as landed in 0.4.1, PR #45)
Floor the friction before seeding:
fr <- terra::ifel(friction == 0, 1e-6, friction)
cost <- terra::ifel(!is.na(streams), 0, fr)
Corrected during implementation: == 0, not the <= 0 originally proposed here.
terra::costDist() rejects a negative cost surface outright (negative friction values not allowed). Flooring <= 0 would have silently disabled that guard, converting input terra
correctly refuses into plausible-looking output — the same failure direction this fix exists to
close. == 0 eliminates exactly the set that costDist(target = 0) would mistake for seeds and
nothing more. A test now asserts the guard still fires.
Also settled by measurement: costDist target matching is exact equality (a 1e-14 cell reads
7.07e-14, not 0), so == 0 is the correct comparison rather than a tolerance.
Note the original "should be ~1400" figures above were wrong in premise — a flat patch legitimately
is cheaper to traverse, so the correct post-fix cost is lower than the no-patch control, not equal
to it. The bug was that flat ground was a source, not that it was cheap to cross.
Blast radius — measured on both shipped DEMs
The fix can only raise a cost that was spuriously zero, so it strictly removes cells from the cost
mask and never adds.
| DEM shipped by this package |
exact-zero slope cells |
cost mask (< 2500) |
valley cells |
bundled dem.tif / slope.tif, 10 m |
0 of 45,726 (min 1.42e-14) |
unchanged |
53,635 -> 53,635 |
pars_dem.tif (MRDEM-30, 30 m, 20.9 Mcell) |
80 of 10.7 M |
-2,289 cells (214 ha), 0 added |
521,028 -> 521,028 |
So MRDEM-30 — the package default source — does produce exact zeros at watershed scale; a small
clip returning none is not evidence about the source. The cost mask moves. The delineation does not,
on either dataset, because the slope, distance and flood criteria plus morphological cleanup absorb
every affected cell — so the shipped vignette artifacts are still current, verified against
pars_valleys.tif directly. That is a property of these two datasets, not a guarantee: where cost
is the binding criterion the delineation will shrink.
Filed separately rather than folded into #40 because it changes VCA output on any DEM that contains
exact zeros — that is a result-changing decision that deserves its own diff and its own regression
check, not a side effect of an attribution feature. Bundled data is unaffected, so no existing test
output moves.
Also fix the roxygen at R/fl_cost_distance.R:9-12 and :19, which currently describes behaviour
the function does not have.
Found during the plan review for #40.
Problem
fl_cost_distance()treats every zero-friction cell as a stream seed, not just stream cells.costDist(target = 0)finds all cells equal to 0. Stream cells are set to 0 deliberately — but sois any cell whose friction (percent slope) is already exactly 0. The roxygen says otherwise:
Demonstration
Synthetic 50x50 grid, friction 10% everywhere, one stream cell, plus a flat (0%) patch 380 m away:
The bundled test DEM hides this —
terrain()-derived slope has zero cells exactly equal to 0(min 1.42e-14), so no existing test can catch it. Production DEMs do produce exact zeros:
integer-metre DEMs, hydro-flattened lake surfaces, and void-filled plateaus all quantize to flat.
Why it matters more now
Under
fl_valley_confine()alone, a spurious flat seed lowers cost slightly in ground that wasprobably valley anyway. Under per-group attribution (#40) it is worse: a flat patch becomes a free
zero-cost source for whichever group's corridor contains it, so that group's cost mask can spread
across ground its own streams never reach. Cost then stops discriminating between groups exactly
where floodplains are — on flat ground — and membership degenerates toward the distance buffer.
Fix (as landed in 0.4.1, PR #45)
Floor the friction before seeding:
Corrected during implementation:
== 0, not the<= 0originally proposed here.terra::costDist()rejects a negative cost surface outright (negative friction values not allowed). Flooring<= 0would have silently disabled that guard, converting input terracorrectly refuses into plausible-looking output — the same failure direction this fix exists to
close.
== 0eliminates exactly the set thatcostDist(target = 0)would mistake for seeds andnothing more. A test now asserts the guard still fires.
Also settled by measurement:
costDisttarget matching is exact equality (a1e-14cell reads7.07e-14, not 0), so== 0is the correct comparison rather than a tolerance.Note the original "should be ~1400" figures above were wrong in premise — a flat patch legitimately
is cheaper to traverse, so the correct post-fix cost is lower than the no-patch control, not equal
to it. The bug was that flat ground was a source, not that it was cheap to cross.
Blast radius — measured on both shipped DEMs
The fix can only raise a cost that was spuriously zero, so it strictly removes cells from the cost
mask and never adds.
< 2500)dem.tif/slope.tif, 10 mpars_dem.tif(MRDEM-30, 30 m, 20.9 Mcell)So MRDEM-30 — the package default source — does produce exact zeros at watershed scale; a small
clip returning none is not evidence about the source. The cost mask moves. The delineation does not,
on either dataset, because the slope, distance and flood criteria plus morphological cleanup absorb
every affected cell — so the shipped vignette artifacts are still current, verified against
pars_valleys.tifdirectly. That is a property of these two datasets, not a guarantee: where costis the binding criterion the delineation will shrink.
Filed separately rather than folded into #40 because it changes VCA output on any DEM that contains
exact zeros — that is a result-changing decision that deserves its own diff and its own regression
check, not a side effect of an attribution feature. Bundled data is unaffected, so no existing test
output moves.
Also fix the roxygen at
R/fl_cost_distance.R:9-12and:19, which currently describes behaviourthe function does not have.
Found during the plan review for #40.