Transforms Module

Transform functions for processing DEM and score data.

Geographic Transforms

src.terrain.transforms.reproject_raster(src_crs='EPSG:4326', dst_crs='EPSG:32617', nodata_value=nan, num_threads=4)[source]

Generalized raster reprojection function

Parameters:

  • src_crs – Source coordinate reference system

  • dst_crs – Destination coordinate reference system

  • nodata_value – Value to use for areas outside original data

  • num_threads – Number of threads for parallel processing

Returns:

Function that transforms data and returns (data, transform, new_crs)

Used in Great Lakes Elevation Visualization to convert WGS84 to UTM.

src.terrain.transforms.flip_raster(axis='horizontal')[source]

Create a transform function that mirrors (flips) the DEM data. If axis=’horizontal’, it flips top ↔ bottom. (In terms of rows, row=0 becomes row=(height-1).)

If axis=’vertical’, you could do left ↔ right (np.fliplr).

src.terrain.transforms.scale_elevation(scale_factor=1.0, nodata_value=nan)[source]

Create a raster elevation scaling transform function.

Multiplies all elevation values by the scale factor. Useful for reducing or amplifying terrain height without changing horizontal scale.

Parameters:

  • scale_factor (float) – Multiplication factor for elevation values (default: 1.0)

  • nodata_value – Value to treat as no data (default: np.nan)

Returns:

A transform function that scales elevation data

Return type:

function

src.terrain.transforms.downsample_raster(zoom_factor=0.1, method='average', nodata_value=nan)[source]

Create a raster downsampling transform function with specified parameters.

Parameters:

  • zoom_factor – Scaling factor for downsampling (default: 0.1)

  • method – Downsampling method (default: “average”) - “average”: Area averaging - best for DEMs, no overshoot - “lanczos”: Lanczos resampling - sharp, minimal aliasing - “cubic”: Cubic spline interpolation - “bilinear”: Bilinear interpolation - safe fallback

  • nodata_value – Value to treat as no data (default: np.nan)

Returns:

A transform function that downsamples raster data

Return type:

function

Smoothing Transforms

src.terrain.transforms.feature_preserving_smooth(sigma_spatial=3.0, sigma_intensity=None, nodata_value=nan)[source]

Create a feature-preserving smoothing transform using bilateral filtering.

Removes high-frequency noise while preserving ridges, valleys, and drainage patterns. Uses bilateral filtering: spatial Gaussian weighted by intensity similarity.

Parameters:

  • sigma_spatial – Spatial smoothing extent in pixels (default: 3.0). Larger = more smoothing. Typical range: 1-10 pixels.

  • sigma_intensity – Intensity similarity threshold in elevation units. Larger = more smoothing across elevation differences. If None, auto-calculated as 5% of elevation range.

  • nodata_value – Value to treat as no data (default: np.nan)

Returns:

Transform function compatible with terrain.add_transform()

Return type:

function

Bilateral filter that preserves ridges and edges.

src.terrain.transforms.wavelet_denoise_dem(nodata_value=nan, wavelet='db4', levels=3, threshold_sigma=2.0, preserve_structure=True)[source]

Create a transform that removes noise while preserving terrain structure.

Uses wavelet decomposition to separate terrain features (ridges, valleys, drainage patterns) from high-frequency noise. This is smarter than median filtering because it understands that terrain has structure at certain spatial frequencies.

How it works: 1. Decompose DEM into frequency bands using wavelets 2. Estimate noise level from finest (highest-frequency) band 3. Apply soft thresholding to remove coefficients below noise threshold 4. Reconstruct DEM from cleaned coefficients

The result preserves terrain structure while removing sensor noise, SRTM artifacts, and other high-frequency disturbances.

Parameters:

  • nodata_value – Value to treat as no data (default: np.nan)

  • wavelet (str) – Wavelet type (default: “db4” - Daubechies 4). Options: “db4” (smooth), “haar” (sharp edges), “sym4” (symmetric). - “db4”: Best for natural terrain (smooth transitions) - “haar”: Best for urban/artificial structures - “sym4”: Good balance, symmetric filtering

  • levels (int) – Decomposition levels (default: 3). More levels = coarser structure preserved. Each level halves the resolution. - 2: Preserves finer detail, removes less noise - 3: Good balance (recommended) - 4: Aggressive smoothing, may blur small features

  • threshold_sigma (float) – Noise threshold multiplier (default: 2.0). Higher = more aggressive denoising. - 1.5: Light denoising, preserves more detail - 2.0: Standard denoising (recommended) - 3.0: Aggressive, may remove subtle features

  • preserve_structure (bool) – If True, only denoise highest-frequency band to maximize structure preservation. If False, denoise all bands.

Returns:

A transform function for use with terrain.add_transform()

Return type:

function

Example

>>> # Standard terrain denoising
>>> terrain.add_transform(wavelet_denoise_dem())

>>> # Aggressive denoising for very noisy DEM
>>> terrain.add_transform(wavelet_denoise_dem(threshold_sigma=3.0, levels=4))

>>> # Light denoising, preserve maximum detail
>>> terrain.add_transform(wavelet_denoise_dem(threshold_sigma=1.5, levels=2))

Frequency-aware denoising. Used in Combined Render: Full-Featured Example.

Example:

terrain.add_transform(wavelet_denoise_dem(
    wavelet='db4',
    levels=3,
    threshold_sigma=2.0
))
src.terrain.transforms.slope_adaptive_smooth(slope_threshold=2.0, smooth_sigma=5.0, transition_width=1.0, nodata_value=nan, elevation_scale=1.0, edge_threshold=None, edge_window=5, strength=1.0)[source]

Create a transform that smooths flat areas more aggressively than hilly areas.

This addresses the problem of buildings/structures appearing as bumps in flat regions. Flat areas get strong Gaussian smoothing to remove these artifacts, while slopes and hills are preserved with minimal smoothing.

How it works:

  1. Compute local slope at each pixel using gradient magnitude

  2. Create a smooth weight mask: 1.0 where flat, 0.0 where steep

  3. Apply Gaussian blur to entire DEM

  4. Blend: output = original * (1-weight) + smoothed * weight

The transition from “flat” to “steep” is smooth (using sigmoid) to avoid visible boundaries in the output.

Parameters:

  • slope_threshold (float) – Slope angle in degrees below which terrain is considered “flat” (default: 2.0 degrees). - 1.0°: Very aggressive, only smooths nearly horizontal areas - 2.0°: Good default, smooths typical flat areas with buildings - 5.0°: Smooths gentle slopes too

  • smooth_sigma (float) – Gaussian blur sigma in pixels (default: 5.0). Controls the strength of smoothing in flat areas. - 3.0: Light smoothing - 5.0: Moderate smoothing (recommended) - 10.0: Very strong smoothing, may blur valid terrain features

  • transition_width (float) – Width of transition zone in degrees (default: 1.0). Controls how quickly smoothing fades off above threshold. - 0.5: Sharp transition - 1.0: Smooth transition (recommended) - 2.0: Very gradual transition

  • nodata_value – Value to treat as no data (default: np.nan)

  • elevation_scale (float) – Scale factor that was applied to elevation data (default: 1.0). If elevation was scaled (e.g., by scale_elevation(0.0001)), pass that factor here so slope computation uses real-world elevation differences. The gradient is divided by this factor to recover true slopes.

  • edge_threshold (float | None) – Elevation difference threshold for edge preservation (default: None). If set, sharp elevation discontinuities (like lake boundaries) are preserved. Areas where local elevation range exceeds this threshold are not smoothed. - None: Disabled (original behavior) - 5.0: Preserve edges with >5m elevation change - 10.0: Only preserve very sharp edges (>10m change) Recommended: 3-10m depending on terrain features to preserve.

  • edge_window (int) – Window size for edge detection (default: 5). Larger windows detect edges over broader areas but may over-protect.

  • strength (float) – Overall smoothing strength multiplier (default: 1.0). Scales the maximum smoothing effect in flat areas. - 1.0: Full smoothing (original behavior) - 0.5: Half the smoothing effect - 0.25: Gentle smoothing - 0.0: No smoothing (transform has no effect)

Returns:

A transform function for use with terrain.add_transform()

Return type:

function

Example

>>> # Standard: smooth flat areas with >2° slopes preserved
>>> terrain.add_transform(slope_adaptive_smooth())

>>> # Aggressive: smooth anything below 5° slope
>>> terrain.add_transform(slope_adaptive_smooth(slope_threshold=5.0, smooth_sigma=8.0))

>>> # Conservative: only smooth very flat areas, light blur
>>> terrain.add_transform(slope_adaptive_smooth(slope_threshold=1.0, smooth_sigma=3.0))

>>> # Compensate for prior scale_elevation(0.0001)
>>> terrain.add_transform(slope_adaptive_smooth(elevation_scale=0.0001))

>>> # Preserve lake boundaries and other sharp edges
>>> terrain.add_transform(slope_adaptive_smooth(edge_threshold=5.0))

>>> # Gentle smoothing (25% of full effect)
>>> terrain.add_transform(slope_adaptive_smooth(strength=0.25))

Smooths flat areas (buildings) while preserving hills. Used in Combined Render: Full-Featured Example.

src.terrain.transforms.remove_bumps(kernel_size=3, structure='disk', strength=1.0)[source]

Remove local maxima (bumps) from DEM using morphological opening.

Morphological opening = erosion followed by dilation. This operation: - Removes small bright features (buildings, trees, noise) - Never creates new local maxima (mathematically guaranteed) - Preserves larger terrain features and overall shape - Leaves valleys and depressions untouched

This is the standard approach for “removing buildings from DEMs” in geospatial processing.

Parameters:

  • kernel_size (int) – Size of the structuring element (default: 3). Controls the maximum size of bumps to remove: - 1: Removes features up to ~2 pixels across (very subtle) - 3: Removes features up to ~6 pixels across - 5: Removes features up to ~10 pixels across For 30m DEMs, size=3 removes ~180m features

  • structure (str) – Shape of structuring element (default: “disk”). - “disk”: Circular, isotropic (recommended) - “square”: Faster but may create artifacts on diagonals

  • strength (float) – Blend factor between original and opened result (default: 1.0). - 0.0: No effect (returns original) - 0.5: Half the bump removal effect (subtle) - 1.0: Full bump removal (original behavior) Values between 0 and 1 provide fine-grained control.

Returns:

A transform function for use with terrain.add_transform()

Return type:

function

Example

>>> # Remove small bumps (buildings on 30m DEM)
>>> terrain.add_transform(remove_bumps(kernel_size=3))

>>> # Subtle bump reduction (50% strength)
>>> terrain.add_transform(remove_bumps(kernel_size=1, strength=0.5))

>>> # More aggressive bump removal
>>> terrain.add_transform(remove_bumps(kernel_size=5))

Morphological opening to remove local maxima.

src.terrain.transforms.despeckle_dem(nodata_value=nan, kernel_size=3)[source]

Create a transform that removes isolated elevation noise using median filtering (GPU-accelerated).

Unlike bilateral smoothing (–smooth) which preserves edges but can look patchy, median filtering uniformly removes local outliers/speckles across the entire DEM. This is better for removing sensor noise or small DEM artifacts.

Uses PyTorch GPU acceleration when available (5-10x speedup on CUDA).

For smarter frequency-aware denoising that preserves terrain structure, use wavelet_denoise_dem() instead.

Parameters:

  • nodata_value – Value to treat as no data (default: np.nan)

  • kernel_size (int) – Size of median filter kernel (default: 3 for 3x3). Must be odd integer ≥3. Larger = more smoothing. - 3: Removes single-pixel noise (recommended) - 5: Removes 2x2 artifacts - 7: Stronger smoothing, may affect small terrain features

Returns:

A transform function for use with terrain.add_transform()

Return type:

function

Example

>>> terrain.add_transform(despeckle_dem(kernel_size=3))

Median filter for isolated outliers.

Score Data Transforms

src.terrain.transforms.smooth_score_data(scores, sigma_spatial=3.0, sigma_intensity=None)[source]

Smooth score data using bilateral filtering.

Applies feature-preserving smoothing to reduce blocky pixelation from low-resolution source data (e.g., SNODAS ~925m) when displayed on high-resolution terrain (~30m DEM).

Uses bilateral filtering: smooths within similar-intensity regions while preserving edges between different score zones.

Parameters:

  • scores (ndarray) – 2D numpy array of score values (typically 0-1 range)

  • sigma_spatial (float) – Spatial smoothing extent in pixels (default: 3.0). Larger = more smoothing. Typical range: 1-10 pixels.

  • sigma_intensity (float) – Intensity similarity threshold in score units. Larger = more smoothing across score differences. If None, auto-calculated as 15% of score range (good for 0-1 data).

Returns:

Smoothed score array with same shape as input. NaN values are preserved. Output is clipped to [0, 1] range.

Return type:

ndarray

Example

>>> # Smooth blocky SNODAS-derived scores
>>> sledding_scores = load_score_grid("sledding_scores.npz")
>>> smoothed = smooth_score_data(sledding_scores, sigma_spatial=5.0)

Reduces blockiness in low-resolution SNODAS data.

src.terrain.transforms.despeckle_scores(scores, kernel_size=3)[source]

Remove isolated speckles from score data using median filtering (GPU-accelerated).

Unlike bilateral filtering which preserves edges, median filtering replaces each pixel with the median of its neighborhood. This effectively removes isolated outlier pixels (speckles) while preserving larger regions.

Uses PyTorch GPU acceleration when available (5-10x speedup on CUDA).

Use case: SNODAS snow data upsampled to high-res DEM often has isolated low-score pixels (speckles) in otherwise high-score regions due to resolution mismatch. These appear as visual noise in the rendered terrain.

Parameters:

  • scores (ndarray) – 2D array of score values (typically 0-1 range)

  • kernel_size (int) – Size of median filter kernel (default: 3 for 3x3). Larger kernels remove larger speckle clusters but may affect legitimate small features. Common values: 3, 5, 7.

Returns:

Despeckled score array with same shape as input. NaN values are preserved.

Return type:

ndarray

Example

>>> # Remove single-pixel speckles
>>> despeckled = despeckle_scores(scores, kernel_size=3)
>>> # Remove up to 2x2 speckle clusters
>>> despeckled = despeckle_scores(scores, kernel_size=5)
src.terrain.transforms.upscale_scores(scores, scale=4, method='auto', nodata_value=nan)[source]

Upscale score grid to reduce blockiness when applied to terrain.

Uses AI super-resolution (Real-ESRGAN) when available, falling back to bilateral upscaling for edge-preserving smoothness.

Parameters:

  • scores (ndarray) – Input score grid (2D numpy array)

  • scale (int) – Upscaling factor (default: 4, meaning 4x resolution)

  • method (str) – Upscaling method: - “auto”: Try Real-ESRGAN, fall back to bilateral - “esrgan”: Use Real-ESRGAN (requires optional realesrgan package) - “bilateral”: Use bilateral filter upscaling (no extra dependencies) - “bicubic”: Simple bicubic interpolation

  • nodata_value (float) – Value treated as no data (default: np.nan)

Returns:

Upscaled score grid with smoother gradients

Return type:

ndarray

Note

The “esrgan” method requires the optional realesrgan package:

pip install realesrgan

Or install terrain-maker with the upscale extra:

pip install terrain-maker[upscale]

Without it, “auto” will fall back to “bilateral” which produces good results without ML dependencies.

Example

>>> scores_hires = upscale_scores(sledding_scores, scale=4)
>>> # Now scores_hires is 4x the resolution with smoother transitions

>>> # Force bilateral method (no ML dependencies)
>>> scores_hires = upscale_scores(scores, scale=4, method="bilateral")

AI super-resolution for score data.

Caching

class src.terrain.cache.TransformCache(cache_dir=None, enabled=True)[source]

Bases: object

Cache for transform pipeline results with dependency tracking.

Tracks chains of transforms (reproject -> smooth -> water_detect) and computes cache keys that incorporate the full dependency chain, ensuring downstream caches are invalidated when upstream params change.

Parameters:

  • cache_dir (Path | None)

  • enabled (bool)

cache_dir

Directory where cache files are stored

enabled

Whether caching is enabled

dependencies

Graph of transform dependencies

transforms

Registered transforms with their parameters

Example:

from src.terrain.cache import TransformCache

cache = TransformCache(cache_dir='.cache')
result = cache.get_or_compute(
    'my_transform',
    compute_fn=lambda: expensive_operation(),
    params={'key': 'value'}
)
__init__(cache_dir=None, enabled=True)[source]

Initialize transform cache.

Parameters:

  • cache_dir (Path | None) – Directory for cache files. If None, uses .transform_cache/

  • enabled (bool) – Whether caching is enabled (default: True)

compute_transform_hash(upstream_hash, transform_name, params)[source]

Compute cache key from upstream hash and transform parameters.

The key incorporates: - Upstream cache key (propagating the full dependency chain) - Transform name - All transform parameters (sorted for determinism)

Parameters:

  • upstream_hash (str) – Hash of upstream data/transform

  • transform_name (str) – Name of this transform (e.g., “reproject”, “smooth”)

  • params (dict) – Transform parameters dict

Returns:

SHA256 hash string (64 chars)

Return type:

str

get_cache_path(cache_key, transform_name)[source]

Get path for cache file.

Parameters:

  • cache_key (str) – Cache key hash

  • transform_name (str) – Name of transform

Returns:

Path to cache .npz file

Return type:

Path

get_metadata_path(cache_key, transform_name)[source]

Get path for metadata file.

Parameters:

  • cache_key (str) – Cache key hash

  • transform_name (str) – Name of transform

Returns:

Path to metadata .json file

Return type:

Path

save_transform(cache_key, data, transform_name, metadata=None)[source]

Save transform result to cache.

Parameters:

  • cache_key (str) – Cache key hash

  • data (ndarray) – Transform result array

  • transform_name (str) – Name of transform

  • metadata (dict | None) – Optional additional metadata

Returns:

Tuple of (cache_path, metadata_path) or (None, None) if disabled

Return type:

Tuple[Path | None, Path | None]

load_transform(cache_key, transform_name)[source]

Load transform result from cache.

Parameters:

  • cache_key (str) – Cache key hash

  • transform_name (str) – Name of transform

Returns:

Cached array or None if cache miss/disabled

Return type:

ndarray | None

register_dependency(child, upstream)[source]

Register a dependency between transforms.

Parameters:

  • child (str) – Name of dependent transform

  • upstream (str) – Name of upstream transform it depends on

Return type:

None

register_transform(name, upstream, params)[source]

Register a transform with its parameters.

Parameters:

  • name (str) – Transform name

  • upstream (str) – Name of upstream dependency

  • params (dict) – Transform parameters

Return type:

None

get_dependency_chain(transform_name)[source]

Get full dependency chain for a transform.

Parameters:

transform_name (str) – Name of transform

Returns:

List of transform names from root to target

Return type:

list[str]

get_full_cache_key(transform_name, source_hash)[source]

Compute full cache key incorporating dependency chain.

Parameters:

  • transform_name (str) – Target transform name

  • source_hash (str) – Hash of original source data

Returns:

Cache key hash

Return type:

str

invalidate_downstream(transform_name)[source]

Invalidate all caches downstream of a transform.

Parameters:

transform_name (str) – Name of transform whose downstream should be invalidated

Returns:

Number of cache files deleted

Return type:

int