r.dem.lod
Compute Level of Detection (LoD) for DEM difference maps
r.dem.lod [dem=name] [reference=name] [dod=name] output=name method=string [confidence=float] [window=integer] [point_density=name] [nmad=float] [stable_mask=name] [floor=float] [min_stable=integer] [sigma_extra=name [,name,...]] [output_sigma=name] [output_domain=name] [--overwrite] [--verbose] [--quiet] [--qq] [--ui]
Example:
r.dem.lod dod=name output=name method=local
grass.tools.Tools.r_dem_lod(dem=None, reference=None, dod=None, output, method="local", confidence=0.95, window=21, point_density=None, nmad=None, stable_mask=None, floor=0.0, min_stable=25, sigma_extra=None, output_sigma=None, output_domain=None, overwrite=None, verbose=None, quiet=None, superquiet=None)
Example:
tools = Tools()
tools.r_dem_lod(dod="name", output="name", method="local")
This grass.tools API is experimental in version 8.5 and expected to be stable in version 8.6.
grass.script.run_command("r.dem.lod", dem=None, reference=None, dod=None, output, method="local", confidence=0.95, window=21, point_density=None, nmad=None, stable_mask=None, floor=0.0, min_stable=25, sigma_extra=None, output_sigma=None, output_domain=None, overwrite=None, verbose=None, quiet=None, superquiet=None)
Example:
gs.run_command("r.dem.lod", dod="name", output="name", method="local")
Parameters
dem=name
Post-event DEM (co-registered); requires reference
reference=name
Reference DEM; requires dem
dod=name
Precomputed (e.g. debiased) DEM of Difference; alternative to dem+reference
output=name [required]
Output LoD raster
method=string [required]
LoD method: global (uniform) or local (spatially variable)
Allowed values: global, local
Default: local
confidence=float
Confidence level (e.g., 0.95 for 95% CI)
Default: 0.95
window=integer
Moving window size (cells) for local LoD, must be odd
Default: 21
point_density=name
Point cloud density raster (pts per square meter), used in local LoD
nmad=float
Pre-computed NMAD (m), skips stable pixel estimation if provided
stable_mask=name
Raster mask of stable pixels for sigma estimation (1=stable, null=unstable)
floor=float
Vertical uncertainty floor (meters, 1 sigma); method=global adds it in quadrature with the NMAD, method=local treats it as the flight-wide stable-residual NMAD and adds only its long-wavelength excess over the windowed dispersion, once (two-scale decomposition)
Default: 0.0
min_stable=integer
Minimum stable cells inside the window for the local dispersion to be defined (guards the sparse-coverage sigma_win=0 degeneracy)
Default: 25
sigma_extra=name [,name,...]
Additional 1-sigma uncertainty rasters added in quadrature
output_sigma=name
Output combined 1-sigma uncertainty raster (output = z * this)
output_domain=name
Output significance-domain raster (1 = LoD defined, NULL elsewhere; method=local)
--overwrite
Allow output files to overwrite existing files
--help
Print usage summary
--verbose
Verbose module output
--quiet
Quiet module output
--qq
Very quiet module output
--ui
Force launching GUI dialog
dem : str | np.ndarray, optional
Post-event DEM (co-registered); requires reference
Used as: input, raster, name
reference : str | np.ndarray, optional
Reference DEM; requires dem
Used as: input, raster, name
dod : str | np.ndarray, optional
Precomputed (e.g. debiased) DEM of Difference; alternative to dem+reference
Used as: input, raster, name
output : str | type(np.ndarray) | type(np.array) | type(gs.array.array), required
Output LoD raster
Used as: output, raster, name
method : str, required
LoD method: global (uniform) or local (spatially variable)
Allowed values: global, local
Default: local
confidence : float, optional
Confidence level (e.g., 0.95 for 95% CI)
Default: 0.95
window : int, optional
Moving window size (cells) for local LoD, must be odd
Default: 21
point_density : str | np.ndarray, optional
Point cloud density raster (pts per square meter), used in local LoD
Used as: input, raster, name
nmad : float, optional
Pre-computed NMAD (m), skips stable pixel estimation if provided
stable_mask : str | np.ndarray, optional
Raster mask of stable pixels for sigma estimation (1=stable, null=unstable)
Used as: input, raster, name
floor : float, optional
Vertical uncertainty floor (meters, 1 sigma); method=global adds it in quadrature with the NMAD, method=local treats it as the flight-wide stable-residual NMAD and adds only its long-wavelength excess over the windowed dispersion, once (two-scale decomposition)
Default: 0.0
min_stable : int, optional
Minimum stable cells inside the window for the local dispersion to be defined (guards the sparse-coverage sigma_win=0 degeneracy)
Default: 25
sigma_extra : str | list[str], optional
Additional 1-sigma uncertainty rasters added in quadrature
Used as: input, raster, name
output_sigma : str | type(np.ndarray) | type(np.array) | type(gs.array.array), optional
Output combined 1-sigma uncertainty raster (output = z * this)
Used as: output, raster, name
output_domain : str | type(np.ndarray) | type(np.array) | type(gs.array.array), optional
Output significance-domain raster (1 = LoD defined, NULL elsewhere; method=local)
Used as: output, raster, name
overwrite : bool, optional
Allow output files to overwrite existing files
Default: None
verbose : bool, optional
Verbose module output
Default: None
quiet : bool, optional
Quiet module output
Default: None
superquiet : bool, optional
Very quiet module output
Default: None
Returns:
result : grass.tools.support.ToolResult | np.ndarray | tuple[np.ndarray] | None
If the tool produces text as standard output, a ToolResult object will be returned. Otherwise, None will be returned. If an array type (e.g., np.ndarray) is used for one of the raster outputs, the result will be an array and will have the shape corresponding to the computational region. If an array type is used for more than one raster output, the result will be a tuple of arrays.
Raises:
grass.tools.ToolError: When the tool ended with an error.
dem : str, optional
Post-event DEM (co-registered); requires reference
Used as: input, raster, name
reference : str, optional
Reference DEM; requires dem
Used as: input, raster, name
dod : str, optional
Precomputed (e.g. debiased) DEM of Difference; alternative to dem+reference
Used as: input, raster, name
output : str, required
Output LoD raster
Used as: output, raster, name
method : str, required
LoD method: global (uniform) or local (spatially variable)
Allowed values: global, local
Default: local
confidence : float, optional
Confidence level (e.g., 0.95 for 95% CI)
Default: 0.95
window : int, optional
Moving window size (cells) for local LoD, must be odd
Default: 21
point_density : str, optional
Point cloud density raster (pts per square meter), used in local LoD
Used as: input, raster, name
nmad : float, optional
Pre-computed NMAD (m), skips stable pixel estimation if provided
stable_mask : str, optional
Raster mask of stable pixels for sigma estimation (1=stable, null=unstable)
Used as: input, raster, name
floor : float, optional
Vertical uncertainty floor (meters, 1 sigma); method=global adds it in quadrature with the NMAD, method=local treats it as the flight-wide stable-residual NMAD and adds only its long-wavelength excess over the windowed dispersion, once (two-scale decomposition)
Default: 0.0
min_stable : int, optional
Minimum stable cells inside the window for the local dispersion to be defined (guards the sparse-coverage sigma_win=0 degeneracy)
Default: 25
sigma_extra : str | list[str], optional
Additional 1-sigma uncertainty rasters added in quadrature
Used as: input, raster, name
output_sigma : str, optional
Output combined 1-sigma uncertainty raster (output = z * this)
Used as: output, raster, name
output_domain : str, optional
Output significance-domain raster (1 = LoD defined, NULL elsewhere; method=local)
Used as: output, raster, name
overwrite : bool, optional
Allow output files to overwrite existing files
Default: None
verbose : bool, optional
Verbose module output
Default: None
quiet : bool, optional
Quiet module output
Default: None
superquiet : bool, optional
Very quiet module output
Default: None
DESCRIPTION
r.dem.lod computes a Level of Detection (LoD) raster for DEM differencing, the elevation-change magnitude below which a measured difference cannot be distinguished from noise. It supports a global (spatially uniform) and a local (spatially variable) mode.
The uncertainty of the difference is estimated from the Normalized Median Absolute Deviation (NMAD) of the elevation residuals on stable terrain:
NMAD = 1.4826 * median(|dh - median(dh)|)
where dh is the residual and z below is the two-tailed normal critical
value for the requested confidence level. The residual comes either from
dem minus reference, or from a precomputed (typically bias-corrected)
difference supplied via dod. Both paths difference the surfaces before
estimation, so the NMAD estimates the dispersion of the difference directly
and no sqrt(2) epoch factor applies on either path.
method=global
A single NMAD is estimated over the whole study area (or the stable_mask region when supplied), producing a uniform LoD applied everywhere:
LoD = z * sqrt(NMAD^2 + floor^2)
A non-zero floor here must be an independent registration budget (for example a shift-parameter standard error), not the same population's NMAD, or the two terms double-count.
method=local
A spatially variable LoD is computed in a Gaussian-weighted moving window of size window (the weighting matches r.dem.stats, so the two tools' windowed dispersion estimators are numerically comparable):
sigma_win(x,y) = 1.4826 * median_W(|dh - median_W(dh)|)
s_long^2 = max(0, floor^2 - median_stable(sigma_win^2))
LoD(x,y) = z * sqrt(sigma_win^2 + s_long^2 + sum(sigma_extra_i^2))
floor is the flight-wide stable-residual NMAD (1 sigma). Because the
windowed NMAD already measures the short-wavelength part of that same
error budget, only the long-wavelength excess s_long is added, once (a
plain quadrature of sigma_win and the floor would count the
short-wavelength component twice and operate a nominal 95% limit at
roughly 99%). s_long is reported in the message output; note that
output_sigma therefore merges a per-cell random component with a
spatially correlated one, so it must not be aggregated to volume
uncertainty by sqrt(N) averaging. The local dispersion is defined only
where at least min_stable stable cells fall inside the window (with
fewer, the windowed median degenerates: a single stable cell yields
sigma_win = 0 exactly); sigma_extra accepts additional 1-sigma rasters (for
example the bias-model coefficient SE from r.dem.bias). NULL in any
sigma_extra or point_density raster propagates to the LoD: a cell
whose uncertainty is unknown is untestable, so the significance domain
shrinks rather than silently assuming zero extra uncertainty. The optional
output_sigma raster stores the combined 1-sigma surface (output equals
z times it), and output_domain marks the cells where the LoD is
defined.
With a stable_mask, the windowed dispersion exists only within the window's reach of stable cells, and the LoD is deliberately NOT extended beyond that reach: coefficient-style uncertainty cannot express model-form error under extrapolation, so an extension requires an out-of-sample error envelope (against independent validation data), which is outside this tool's scope. The tested share of observed cells is always reported.
When a point_density raster is supplied, the local uncertainty is penalised in sparsely sampled areas before the LoD is computed.
NOTES
A precomputed nmad value can be supplied to skip the stable-pixel estimation, and a stable_mask restricts the residual statistics to terrain assumed unchanged (roads, parking lots, bare ground).
The output LoD raster can be passed directly to r.dem.change as the lod input for significance thresholding. For full per-source uncertainty propagation (combining several error rasters in quadrature) use r.dem.errprop.
The tool requires the Python scipy package.
EXAMPLES
The commands below use the example scene built in the r.dem toolset manual, which is derived from the North Carolina sample dataset. Build it there first.
A single detection limit for the whole map, estimated from the stable residuals:
g.region raster=elev_lid792_1m
r.dem.lod dod=dod_debiased output=lod_global method=global \
stable_mask=stable_lod confidence=0.95
The debiased residual on the stable cells is the injected noise, so the
result is z(0.95) times its NMAD:
NMAD: 0.0890 m
LoD: 0.1743 m (uniform)
A spatially variable limit, keeping the combined 1-sigma surface for r.dem.errprop and the domain raster that marks where the limit is defined:
r.dem.lod dod=dod_debiased output=lod_local method=local window=21 \
stable_mask=stable_lod output_sigma=sigma_combined \
output_domain=lod_domain confidence=0.95
Because the survey noise varies across this scene, so does the limit: it runs from roughly 0.12 m on the smooth fields to over 0.30 m under canopy, against a single uniform value of 0.175 m. That is the case for method=local, and the reason the mask must include forest.
However, the limit is undefined wherever no stable cell falls inside the window, which on this scene means the interior of the change features. Fall back to the uniform limit there before thresholding:
r.mapcalc "lod_filled = if(isnull(lod_local), lod_global, lod_local)"
Add the bias-model coefficient SE from r.dem.bias to the quadrature when the regression path was used:
r.dem.lod dod=dod_regression output=lod_with_se method=local \
stable_mask=stable_terrain sigma_extra=bias_se \
output_sigma=sigma_with_se
Figure: The spatially variable Level of Detection against the uniform one, on
a shared scale. White marks cells with no stable cell inside the window, where
the local limit is undefined.
REFERENCES
- Wheaton et al. (2010), Earth Surface Processes and Landforms 35:136-156.
- Hoehle and Hoehle (2009), ISPRS J. Photogramm. Remote Sens. 64:398-406.
SEE ALSO
r.dem, r.dem.change, r.dem.errprop, r.dem.stats, r.neighbors, r.univar
AUTHORS
Corey T. White, Center for Geospatial Analytics, NC State University
SOURCE CODE
Available at: r.dem.lod source code
(history)
Latest change: Thursday Aug 20 19:32:37 2026 in commit d08d385