Skip to contents

Orchestrates the full VCA pipeline: slope thresholding, cost-distance analysis, flood surface modelling, and morphological cleanup to identify unconfined valley bottoms.

Usage

fl_valley_confine(
  dem,
  streams,
  area_field,
  slope = NULL,
  slope_threshold = 9,
  max_width = 2000,
  cost_threshold = 2500,
  flood_factor = 6,
  precip = NULL,
  waterbodies = NULL,
  channel_buffer = NULL,
  size_threshold = 5000,
  hole_threshold = 2500,
  field = NULL
)

Arguments

dem

A SpatRaster of elevation.

streams

An sf linestring object or a SpatRaster of rasterized streams. If sf, it is rasterized using area_field. If a SpatRaster, its cell values are used as-is and must already be upstream contributing area in hectares.

area_field

Character. Column of streams holding upstream contributing area in hectares, rasterized onto the DEM grid by fl_stream_rasterize(). Required when streams is sf — there is no default. The rasterized values become the drainage-area term of the bankfull regression in fl_flood_surface(), which accepts any positive numeric column without complaint: channel width in place of area returns a smaller floodplain with no error and no warning. Not used when streams is already a SpatRaster. That branch cannot inspect the values it is handed; it warns only when the layer's name gives it away, so a raster burned from any other wrong column carries the same defect one call earlier, undetected.

slope

A SpatRaster of percent slope. If NULL, derived from dem.

slope_threshold

Numeric. Maximum percent slope for valley floor. Default 9.

max_width

Numeric. Maximum valley width in map units (metres). Default 2000.

cost_threshold

Numeric. Maximum accumulated cost distance. Default 2500.

flood_factor

Numeric. Multiplier on bankfull depth. Default 6.

precip

A SpatRaster or numeric scalar of mean annual precipitation in millimetres, converted to cm/yr internally for the bankfull regression. Default NULL, which drops the precipitation term.

waterbodies

An sf polygon object of lakes and/or wetlands, or NULL (default). Waterbody polygons are rasterized onto the valley grid and added to the output after morphological cleanup. No buffer is applied — a lake or wetland in the valley is part of the flood system as-is. No spatial filtering is applied — all polygons are rasterized. Pre-filter to valley-bottom features before calling if headwater waterbodies are not wanted.

channel_buffer

Logical. Buffer streams by their channel_width attribute and add to the valley output. Default TRUE when streams is an sf object with a channel_width column, FALSE otherwise. The stream channel is floodplain but can be sub-pixel at coarse DEM resolution.

size_threshold

Numeric. Minimum valley patch area (m²). Default 5000.

hole_threshold

Numeric. Maximum hole area to fill (m²). Default 2500.

field

Deprecated. The former name of area_field, whose "channel_width" default was wrong for the flood model (#47). Supplying it warns and forwards to area_field; removal is tracked in flooded#53.

Value

A SpatRaster with binary values: 1 = unconfined valley, 0 = confined / hillslope, NA = outside analysis extent.

Details

The algorithm combines four criteria via intersection (AND):

  1. Slope mask — cells with slope <= slope_threshold

  2. Distance mask — cells within max_width / 2 of a stream

  3. Cost distance mask — cells with accumulated cost < cost_threshold

  4. Flood mask — cells identified as flooded by bankfull regression

The combined mask then undergoes morphological cleanup:

  • Closing filter (3x3) to bridge small gaps

  • Fill small holes (< hole_threshold)

  • Remove small patches (< size_threshold)

  • Majority filter (3x3) to smooth edges

After cleanup, optional features are added via logical OR:

  • Channel buffer — streams buffered by channel_width (DEM correction)

  • Waterbodies — user-supplied lake/wetland polygons rasterized as-is

Adapted from the USDA Valley Confinement Algorithm Toolbox (BlueGeo implementation by Devin Cairns, MIT license) and bcfishpass lateral habitat assembly (Simon Norris, Apache 2.0).

Performance

Several internal operations (focal filters, distance calculations, raster math) support multi-threading via terra::terraOptions(). Set threads before calling this function to speed up processing on large rasters:

terra::terraOptions(threads = 12)

On an Apple M4 Max (16 cores), 12 threads reduced runtime from ~3.5 minutes to ~1 minute for a 27M-cell raster (~2,700 km² at 10 m).

See also

fl_stream_rasterize() for how area_field is burned onto the grid, and fl_flood_surface() for the regression that consumes it. fl_mask(), fl_cost_distance(), fl_flood_model(), fl_patch_rm(), fl_valley_poly()

Examples

dem <- terra::rast(system.file("testdata/dem.tif", package = "flooded"))
streams <- sf::st_read(
  system.file("testdata/streams.gpkg", package = "flooded"),
  quiet = TRUE
)
precip_r <- fl_stream_rasterize(streams, dem, field = "map_upstream")

# Basic VCA (channel buffer auto-detected from streams$channel_width)
valleys <- fl_valley_confine(
  dem, streams,
  area_field = "upstream_area_ha",
  precip = precip_r
)
terra::plot(valleys, col = c("grey90", "darkgreen"), main = "Unconfined valleys")


# With waterbodies — fills lake/wetland donut holes
waterbodies <- sf::st_read(
  system.file("testdata/waterbodies.gpkg", package = "flooded"),
  quiet = TRUE
)
valleys_wb <- fl_valley_confine(
  dem, streams,
  area_field = "upstream_area_ha",
  precip = precip_r,
  waterbodies = waterbodies
)
terra::plot(valleys_wb, col = c("grey90", "darkgreen"), main = "With waterbodies")