process_otsu_rasters: Process Burned Area Rasters Using Otsu Thresholding or...

In OtsuFire: Fire Scars, Severity and Regeneration Mapping Using 'Otsu' Thresholding

Process Burned Area Rasters Using Otsu Thresholding or Percentile Clipping

Description

This function computes binary burned area masks from a severity index raster (RBR or dNBR) using Otsu's thresholding or percentile clipping. The segmentation can be applied to the full raster or stratified by: - CORINE land cover classes ('corine_raster_path') - WWF ecoregions ('ecoregion_shapefile_path') - or the intersection of both (CORINE ? Ecoregion)

## Key features: - Supports Otsu thresholding after pre-filtering by minimum index values ('otsu_thresholds') or using original values - Supports percentile-based thresholding via 'trim_percentiles' - Enforces a minimum threshold ('min_otsu_threshold_value') when Otsu yields low values - Allows stratification by CORINE, ecoregions, or their intersection - Enables CORINE reclassification using a custom 'reclass_matrix' - Automatically computes RBR or dNBR if 'nbr_pre_path' and 'nbr_post_path' are provided - Outputs: - Burned area binary rasters - Polygon shapefiles (dissolved by threshold label) - Histogram and inter-class variance plots per threshold - Threshold log files

Arguments

raster_path	Path to a single-band RBR or dNBR raster.
nbr_pre_path, nbr_post_path	Optional. Paths to pre- and post-fire NBR rasters for RBR/dNBR calculation.
output_dir	Directory to save output files.
year	Optional year label.
otsu_thresholds	Numeric vector of minimum values to filter raster before applying Otsu.
trim_percentiles	Optional data.frame with 'min' and 'max' columns to define percentile clipping thresholds.
use_original	Logical. Use raw values without filtering (ignored if 'trim_percentiles' is set).
corine_raster_path	Path to CORINE raster for stratified segmentation.
reclassify_corine	Logical. If TRUE, reclassifies CORINE using 'reclass_matrix'.
reclass_matrix	Matrix with original and new values for CORINE reclassification.
peninsula_shapefile	Shapefile used to crop CORINE before reclassification.
output_corine_raster_dir	Directory to save reclassified CORINE raster.
output_corine_vector_dir	Directory to save reclassified CORINE shapefile.
reproject	Logical. Reproject reclassified CORINE raster to EPSG:3035.
resolution	Numeric. Target resolution for reprojected CORINE raster.
corine_classes	Optional vector of reclassified CORINE classes to keep.
ecoregion_shapefile_path	Path to WWF ecoregions shapefile.
ecoregion_field	Column in ecoregion shapefile used for class labels.
ecoregion_classes	Optional. Vector of selected ecoregion class names.
segment_by_intersection	Logical. If TRUE, intersects CORINE and ecoregions.
min_otsu_threshold_value	Minimum acceptable threshold. Used if Otsu value is lower.
python_exe	Path to Python executable.
gdal_polygonize_script	Path to 'gdal_polygonize.py' script.
gdalwarp_path	Path to 'gdalwarp' executable (default = "gdalwarp").
n_rows, n_cols	Number of rows/columns for raster tiling.
tile_overlap	Overlap buffer (in map units) for tiles.
tile	Logical. Apply raster tiling before polygonization.
index_type	"RBR" or "dNBR". Used if 'nbr_pre_path'/'nbr_post_path' are provided.
output_format	Output vector file format: "shp" or "geojson".

Details

- Raster values are internally rescaled to [0, 255] before Otsu thresholding. - Histogram smoothing and variance curves are used for enhanced threshold detection. - Output shapefiles can be in ESRI Shapefile or GeoJSON format. - Intermediate tile shapefiles and rasters are cleaned after processing.

Value

Named list of 'sf' polygons grouped by segmentation. Output files saved to 'output_dir'.

Note

Examples require large external raster files (hosted on Zenodo) and depend on external software (Python, GDAL). Therefore, they are wrapped in dontrun to avoid errors during R CMD check and to ensure portability.

Examples

## Not run: 
# Example 1: Basic Otsu thresholds (global)
process_otsu_rasters(
  raster_path = "data/rbr_1985.tif",
  output_dir = "output/1985",
  otsu_thresholds = c(0, 50, 100),
  python_exe = "C:/Python/python.exe",
  gdal_polygonize_script = "C:/Python/Scripts/gdal_polygonize.py"
)

# Example 2: Global with percentile trimming
process_otsu_rasters(
  raster_path = "data/rbr_1985.tif",
  output_dir = "output/1985",
  trim_percentiles = data.frame(min = 0.01, max = 0.99),
  python_exe = "C:/Python/python.exe",
  gdal_polygonize_script = "C:/Python/Scripts/gdal_polygonize.py",
  output_format = "geojson"
)

# Example 3: Stratified by reclassified CORINE
process_otsu_rasters(
  raster_path = "data/rbr_1985.tif",
  output_dir = "output/1985",
  corine_raster_path = "data/corine_1990_etrs89.tif",
  reclassify_corine = TRUE,
  reclass_matrix = matrix(c(
    22, 1,  # grasslands ? 1
    26, 1,
    32, 1,
    27, 2,  # shrublands ? 2
    29, 2,
    33, 3,  # burnt areas ? 3
    23, 4,  # broadleaved forest ? 4
    25, 5,  # mixed forest ? 5
    24, 6   # coniferous ? 6
  ), ncol = 2, byrow = TRUE),
  corine_classes = c(1, 2, 4, 5, 6),
  peninsula_shapefile = "data/peninsula_clip.shp",
  output_corine_raster_dir = "output/corine",
  output_corine_vector_dir = "output/corine",
  otsu_thresholds = c(50),
  min_otsu_threshold_value = 150,
  python_exe = "C:/Python/python.exe",
  gdal_polygonize_script = "C:/Python/Scripts/gdal_polygonize.py"
)

# Example 4: Stratified by WWF ecoregions (default field: "EnS_name")
process_otsu_rasters(
  raster_path = "data/rbr_1985.tif",
  output_dir = "output/1985",
  ecoregion_shapefile_path = "data/ecoregions_olson.shp",
  otsu_thresholds = c(50),
  min_otsu_threshold_value = 150,
  python_exe = "C:/Python/python.exe",
  gdal_polygonize_script = "C:/Python/Scripts/gdal_polygonize.py"
)

# Example 5: Stratified by WWF ecoregions with custom field
process_otsu_rasters(
  raster_path = "data/rbr_1985.tif",
  output_dir = "output/1985",
  ecoregion_shapefile_path = "data/ecoregions_olson.shp",
  ecoregion_field = "EnZ_name",
  ecoregion_classes = c("Iberian Conifer Forests", "Mediterranean Shrublands"),
  otsu_thresholds = c(50),
  min_otsu_threshold_value = 150,
  python_exe = "C:/Python/python.exe",
  gdal_polygonize_script = "C:/Python/Scripts/gdal_polygonize.py"
)

# Example 6: Stratified by CORINE ? Ecoregion intersection
process_otsu_rasters(
  raster_path = "data/rbr_1985.tif",
  output_dir = "output/1985",
  corine_raster_path = "data/corine_1990_etrs89.tif",
  reclassify_corine = TRUE,
  reclass_matrix = matrix(c(
    22, 1, 26, 1, 32, 1,  # grasslands
    27, 2, 29, 2,         # shrublands
    23, 4, 25, 5, 24, 6   # forests
  ), ncol = 2, byrow = TRUE),
  corine_classes = c(1, 2, 4, 5, 6),
  peninsula_shapefile = "data/peninsula_clip.shp",
  output_corine_raster_dir = "output/corine",
  output_corine_vector_dir = "output/corine",
  ecoregion_shapefile_path = "data/ecoregions_olson.shp",
  ecoregion_field = "EnZ_name",
  segment_by_intersection = TRUE,
  otsu_thresholds = c(50),
  min_otsu_threshold_value = 150,
  python_exe = "C:/Python/python.exe",
  gdal_polygonize_script = "C:/Python/Scripts/gdal_polygonize.py"
)

## End(Not run)

OtsuFire documentation built on June 13, 2025, 3:01 a.m.