oxigdal-terrain 0.1.3

Advanced terrain analysis and DEM processing for OxiGDAL - Derivatives, hydrology, viewshed, and geomorphometry
Documentation

OxiGDAL Terrain

Crates.io Documentation License Rust

Advanced terrain analysis and DEM (Digital Elevation Model) processing library for OxiGDAL. Provides comprehensive geospatial analysis capabilities for terrain derivatives, hydrological modeling, visibility analysis, and geomorphometric classification—all in 100% Pure Rust.

Features

  • Terrain Derivatives: Slope (Horn, Zevenbergen-Thorne), aspect, curvature (profile, plan, total), hillshade (traditional, multidirectional, combined), TPI, TRI, roughness
  • Hydrological Analysis: Flow direction (D8, D-Infinity), flow accumulation, sink filling, watershed delineation, stream network extraction with Strahler ordering
  • Visibility Analysis: Viewshed computation (binary and cumulative), line-of-sight analysis
  • Geomorphometry: Landform classification (Weiss, Iwahashi-Pike), convergence index, positive/negative openness
  • Parallel Processing: Optional parallel computation with Rayon for large datasets
  • Pure Rust: No C/Fortran dependencies—100% safe, idiomatic Rust implementation
  • SciRS2 Integration: Seamless integration with SciRS2-Core for scientific computing
  • Comprehensive Error Handling: Rich error types with no unwrap panic policy
  • Feature-Gated: Modular design with optional feature flags for selective compilation

Installation

Add to your Cargo.toml:

[dependencies]
oxigdal-terrain = "0.1.3"

# With all features enabled
oxigdal-terrain = { version = "0.1.3", features = ["all_features", "parallel"] }

Feature Flags

  • default: Enables std, derivatives, hydrology
  • std: Standard library support (enabled by default)
  • derivatives: Terrain derivatives (slope, aspect, curvature, etc.)
  • hydrology: Hydrological analysis (flow, watershed, streams)
  • visibility: Visibility analysis (viewshed, line-of-sight)
  • geomorphometry: Geomorphometric features (landforms, convergence, openness)
  • parallel: Parallel processing support with Rayon
  • all_features: Enables all analysis features (excludes parallel)

Quick Start

Basic Slope Calculation

use oxigdal_terrain::derivatives::{slope_horn, SlopeUnits};
use scirs2_core::prelude::*;

fn main() -> Result<(), Box<dyn std::error::Error>> {
    // Create a simple DEM (100x100 pixels, 10m cell size)
    let dem = Array2::from_elem((100, 100), 100.0_f32);

    // Calculate slope in degrees using Horn's method
    let slope = slope_horn(&dem, 10.0, SlopeUnits::Degrees, None)?;

    println!("Slope grid shape: {:?}", slope.dim());
    Ok(())
}

Hydrological Analysis

use oxigdal_terrain::hydrology::*;
use scirs2_core::prelude::*;

fn main() -> Result<(), Box<dyn std::error::Error>> {
    let dem = Array2::from_elem((100, 100), 100.0_f32);

    // Fill sinks first (important for D8 flow analysis)
    let filled = fill_sinks(&dem, None)?;

    // Calculate flow direction (D8 algorithm)
    let flow_dir = flow_direction(&filled, 10.0, FlowAlgorithm::D8, None)?;

    // Calculate flow accumulation
    let flow_accum = flow_accumulation(&flow_dir, None)?;

    // Extract streams (cells with flow accumulation > 100)
    let streams = extract_streams(&flow_accum, 100.0, None)?;

    println!("Number of stream cells: {:?}",
             streams.iter().filter(|&&v| v > 0.0).count());

    Ok(())
}

Visibility Analysis

use oxigdal_terrain::visibility::*;
use scirs2_core::prelude::*;

fn main() -> Result<(), Box<dyn std::error::Error>> {
    let dem = Array2::from_elem((100, 100), 100.0_f32);

    // Binary viewshed from observer at (50, 50) with observer height 2m
    let viewshed = viewshed_binary(
        &dem,
        10.0,  // cell size in meters
        50, 50, // observer position (x, y)
        2.0,    // observer height above terrain
        50.0,   // target height above terrain
        None,   // radius limit (None = whole DEM)
    )?;

    println!("Visible cells: {}",
             viewshed.iter().filter(|&&v| v > 0.0).count());

    Ok(())
}

Geomorphometric Classification

use oxigdal_terrain::geomorphometry::*;
use scirs2_core::prelude::*;

fn main() -> Result<(), Box<dyn std::error::Error>> {
    let dem = Array2::from_elem((100, 100), 100.0_f32);

    // Classify landforms using Weiss method
    let landforms = classify_weiss(
        &dem,
        10.0,   // cell size
        21,     // analysis radius in cells
        None,   // nodata value
    )?;

    println!("Landform classification complete");

    Ok(())
}

API Overview

Module Description Key Functions
derivatives Terrain surface characteristics slope, aspect, curvature, hillshade, tpi, tri, roughness
hydrology Hydrological flow and watersheds flow_direction, flow_accumulation, fill_sinks, watershed_from_point, extract_streams
visibility Viewshed and line-of-sight viewshed_binary, viewshed_cumulative, line_of_sight
geomorphometry Landform and surface properties classify_weiss, classify_iwahashi_pike, convergence_index, openness

Derivatives Module

  • slope() - Calculate slope (auto-select algorithm)
  • slope_horn() - Horn's method for slope calculation
  • slope_zevenbergen_thorne() - Zevenbergen-Thorne method
  • aspect() - Calculate aspect direction with flat handling options
  • curvature() - Calculate total curvature
  • profile_curvature() - Slope of slope along steepest descent
  • plan_curvature() - Slope of aspect perpendicular to steepest descent
  • tangential_curvature() - Surface curvature perpendicular to slope direction
  • hillshade() - Shaded relief for visualization (auto-select algorithm)
  • hillshade_traditional() - Traditional hillshade rendering
  • hillshade_multidirectional() - 16-directional hillshade
  • hillshade_combined() - Combined hillshade with slope adjustment
  • tpi() - Topographic Position Index (elevation relative to neighbors)
  • tpi_parallel() - Parallel TPI computation
  • tri() - Terrain Ruggedness Index
  • tri_parallel() - Parallel TRI computation
  • tri_riley() - Riley's TRI variant
  • roughness() - Surface roughness (auto-select method)
  • roughness_range() - Range of elevations in neighborhood
  • roughness_stddev() - Standard deviation of elevations
  • vector_ruggedness_measure() - VRM for roughness quantification

Hydrology Module

  • flow_direction() - Flow direction calculation (D8 or D-Infinity)
  • flow_direction_d8() - D8 (8-directional) flow algorithm
  • flow_direction_dinf() - D-Infinity continuous flow algorithm
  • flow_accumulation() - Cumulative flow contribution
  • fill_sinks() - Hydrologically correct DEM preprocessing
  • extract_streams() - Extract stream network from flow accumulation
  • strahler_order() - Assign Strahler stream order
  • watershed_from_point() - Delineate watershed from pour point

Visibility Module

  • viewshed_binary() - Boolean visibility map
  • viewshed_cumulative() - Weighted visibility accumulation
  • line_of_sight() - Line-of-sight visibility between two points

Geomorphometry Module

  • classify_weiss() - Landform classification using Weiss (2001) method
  • classify_iwahashi_pike() - Iwahashi & Pike (2007) classification
  • convergence_index() - Surface convergence/divergence index
  • positive_openness() - Exposed terrain openness
  • negative_openness() - Enclosed terrain openness

Algorithms

Slope Calculation

  • Horn (1981): 3x3 neighborhood with smooth gradient estimation
  • Zevenbergen & Thorne (1987): Improved edge handling for irregular DEMs

Flow Direction

  • D8 (O'Callaghan & Mark, 1984): 8-directional steepest descent
  • D-Infinity (Tarboton, 1997): Continuous flow direction on triangulated surface

Landform Classification

  • Weiss (2001): Profile/plan curvature-based classification
  • Iwahashi & Pike (2007): Slope/curvature-based landform types

Performance

Benchmarks on Apple Silicon (M3 Max) with synthetic DEM:

Operation 100x100 500x500 1000x1000
Slope (Horn) 0.15ms 3.2ms 13ms
Aspect (Horn) 0.18ms 3.5ms 14ms
Curvature (Profile) 0.16ms 3.1ms 12ms
Hillshade (Traditional) 0.22ms 4.2ms 17ms
TPI (radius=1) 0.14ms 2.8ms 11ms
TRI (Standard) 0.13ms 2.6ms 10ms

Note: Performance scales approximately O(n) with grid size. Parallel variants (with parallel feature) show 3-4x speedup on large datasets (1000x1000+).

No Unwrap Policy

This library strictly follows the "no unwrap" policy. All fallible operations return Result<T, TerrainError> with descriptive error variants:

  • InvalidDimensions - DEM dimensions out of bounds
  • InvalidCellSize - Non-positive cell size
  • InvalidObserverPosition - Position outside DEM bounds
  • InvalidAzimuth / InvalidAltitude - Invalid angle parameters
  • FlowDirectionError - Flow computation issues
  • WatershedError - Watershed delineation problems
  • ViewshedError - Visibility computation errors
  • ComputationError - Generic computation failures
  • InsufficientMemory - Memory allocation issues

Example error handling:

use oxigdal_terrain::{derivatives::slope_horn, SlopeUnits, TerrainError};
use scirs2_core::prelude::*;

fn analyze_terrain(dem: &Array2<f32>) -> Result<Array2<f32>, TerrainError> {
    slope_horn(dem, 10.0, SlopeUnits::Degrees, None)
}

fn main() {
    match analyze_terrain(&dem) {
        Ok(slope) => println!("Success: {:?}", slope.dim()),
        Err(e) => eprintln!("Error: {}", e),
    }
}

Pure Rust

This library is 100% Pure Rust with no C/Fortran dependencies. All numerical operations use:

  • SciRS2-Core: Scientific computing primitives and linear algebra
  • Standard Rust: Safe, idiomatic Rust without external C bindings

Examples

The examples directory contains complete, runnable examples:

  • terrain_derivatives.rs - Calculate all terrain derivatives
  • hydrological_analysis.rs - Watershed delineation workflow
  • visibility_analysis.rs - Viewshed computation
  • landform_classification.rs - Geomorphometric classification

Run examples with:

cargo run --example terrain_derivatives --features all_features
cargo run --example hydrological_analysis --features hydrology
cargo run --example visibility_analysis --features visibility
cargo run --example landform_classification --features geomorphometry

Documentation

Integration with OxiGDAL

oxigdal-terrain is a specialized module within the OxiGDAL ecosystem:

use oxigdal::raster::Band;
use oxigdal_terrain::derivatives::slope_horn;

// Load raster from GDAL source
let band = Band::from_file("dem.tif")?;
let dem = band.to_array2::<f32>()?;

// Apply terrain analysis
let slope = slope_horn(&dem, 10.0, SlopeUnits::Degrees, None)?;

Contributing

Contributions welcome! Please ensure:

  1. No unwrap() calls in production code
  2. Comprehensive error handling with TerrainError
  3. Benchmark critical operations with Criterion
  4. Write doc tests for public APIs
  5. Follow COOLJAPAN coding standards

Testing

Run the full test suite:

cargo test --all-features
cargo test --doc --all-features
cargo test --test terrain_test --all-features

Run benchmarks:

cargo bench --bench terrain_bench --all-features

Code Statistics

  • 2,933 lines of Rust code
  • 25 source files organized by feature
  • ~130 lines of documentation
  • 100% Pure Rust implementation

License

This project is licensed under Apache-2.0. See LICENSE for details.

Related COOLJAPAN Projects

  • OxiGDAL - Core geospatial library
  • OxiBLAS - Pure Rust BLAS
  • OxiFFT - Pure Rust FFT
  • SciRS2 - Scientific computing ecosystem
  • NumRS2 - Numerical computing (NumPy-like)
  • Oxicode - Serialization framework

Part of the COOLJAPAN Pure Rust Ecosystem