Skip to content

Instantly share code, notes, and snippets.

@mdsumner
Last active March 11, 2026 19:58
Show Gist options
  • Select an option

  • Save mdsumner/c6ad59afb600b10c5ff602693c40e65e to your computer and use it in GitHub Desktop.

Select an option

Save mdsumner/c6ad59afb600b10c5ff602693c40e65e to your computer and use it in GitHub Desktop.
# 1. Linear percent clip
# low/high are the natural UX knobs; 2%/98% is a solid default,
# tighter for low-contrast scenes.
s2_stretch_linear <- function(x, low = 0.02, high = 0.98) {
v <- terra::values(x, na.rm = TRUE)
lo <- quantile(v, low, na.rm = TRUE)
hi <- quantile(v, high, na.rm = TRUE)
terra::clamp(x, lo, hi, values = TRUE) |>
(\(r) (r - lo) / (hi - lo))()
}
# 2. Log stretch (per-band minimum offset, guards log(0))
# Use only on integer overview values, not native-res reflectance — see gotcha below.
s2_stretch_log <- function(x, offset = NULL) {
if (is.null(offset)) {
mn <- min(terra::global(x, "min", na.rm = TRUE)[[1]])
offset <- if (mn <= 0) abs(mn) + 1 else 0
}
log1p(x + offset) |> s2_stretch_linear()
}
# 3. Log-relative: subtract global minimum before log.
# + 1 guards against log(0) when a pixel equals the minimum exactly.
s2_stretch_log_relative <- function(x, low = 0.02, high = 0.98) {
mn <- min(terra::global(x, "min", na.rm = TRUE)[[1]])
log(x - mn + 1) |> s2_stretch_linear(low, high)
}
# Joint variant: single global minimum AND a single quantile clip across all
# bands combined — preserves inter-band luminance relationships for true-colour RGB.
# In practice the visual difference from per-band is subtle for well-exposed scenes;
# the real fix for colour balance is the shared minimum subtraction.
s2_stretch_log_relative_joint <- function(x, low = 0.02, high = 0.98) {
mn <- min(terra::global(x, "min", na.rm = TRUE)[[1]])
logged <- log(x - mn + 1)
v <- as.vector(terra::values(logged, na.rm = TRUE)) # flatten all bands
lo <- quantile(v, low, na.rm = TRUE)
hi <- quantile(v, high, na.rm = TRUE)
terra::clamp(logged, lo, hi, values = TRUE) |>
(\(r) (r - lo) / (hi - lo))()
}
# 4. Square root — gentler than log, good middle ground.
s2_stretch_sqrt <- function(x) {
mn <- min(terra::global(x, "min", na.rm = TRUE)[[1]])
offset <- if (mn < 0) abs(mn) else 0
sqrt(x + offset) |> s2_stretch_linear()
}
# 5. Histogram equalisation.
# Good when interesting detail spans simultaneously deep shadow and bright areas.
s2_stretch_histeq <- function(x, n = 256) {
v <- terra::values(x, na.rm = TRUE)
breaks <- quantile(v, probs = seq(0, 1, length.out = n + 1), na.rm = TRUE)
out <- x
terra::values(out) <- findInterval(terra::values(x), unique(breaks)) /
length(unique(breaks))
out
}
#Dispatcher
#joint defaults to TRUE for multi-band input so that RGB composites share a common stretch across all three channels.
s2_stretch <- function(x,
method = c("linear", "log", "sqrt", "histeq", "log_relative"),
low = 0.02, high = 0.98,
joint = terra::nlyr(x) > 1, ...) {
method <- match.arg(method)
switch(method,
linear = if (joint) {
v <- terra::values(x, na.rm = TRUE)
lo <- quantile(v, low, na.rm = TRUE)
hi <- quantile(v, high, na.rm = TRUE)
terra::clamp(x, lo, hi, values = TRUE) |>
(\(r) (r - lo) / (hi - lo))()
} else s2_stretch_linear(x, low, high),
log_relative = if (joint) s2_stretch_log_relative_joint(x, low, high)
else s2_stretch_log_relative(x, low, high),
log = s2_stretch_log(x, ...),
sqrt = s2_stretch_sqrt(x),
histeq = s2_stretch_histeq(x, ...)
)
}
# terra::plotRGB() requires [0, 255]; terra::stretch() is used solely for that
# byte-scaling side effect, applied to already-stretched [0, 1] data.
pltrgb <- function(x, ...) {
terra::plotRGB(terra::stretch(x), ...)
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment