Raster Processing Performance Profiling

Measure stored-map and VRT-to-GeoTIFF configurations on your own HEC-RAS project. RasMap.store_all_maps() is the canonical map orchestration function; focused profiling logic remains in RasProcess and RasTerrain. This notebook defines a configuration matrix, graphs the resource evidence, and writes an HTML decision report.

Large-watershed mapping can consume several uncompressed Float32 surfaces per helper. Start with the estimator, retain the default memory_policy="enforce", and use a copied project because stored-map generation temporarily edits the .rasmap.

Python

import html as html_lib
import io
import json
import os
import subprocess
import sys
from pathlib import Path

import matplotlib.pyplot as plt
import pandas as pd

from ras_commander import (
    GeoTiffWriteOptions,
    RasMap,
    RasProcess,
    StoreMapPerformanceOptions,
    init_ras_project,
)
from ras_commander.terrain import RasTerrain


def semantic_raster_signature(signature):
    """Comparison fields that must not change across tuning candidates."""
    fields = (
        "pixel_sha256", "width", "height", "count", "dtypes",
        "crs", "transform", "nodata", "overviews_by_band",
    )
    return {name: signature.get(name) for name in fields}


def semantic_output_signatures(signatures):
    return sorted(
        json.dumps(semantic_raster_signature(item), sort_keys=True)
        for item in signatures.values()
    )

Configure a copied test project

Set RAS_PROFILE_PROJECT to a writable project copy with computed plan results. The Spring River BLE project is a useful large-watershed fixture; plan 02 is the measured Depth Max example. Outputs and reports go under working/ by default and are not committed. All heavy profiling switches default to False.

Python

PROJECT_FOLDER_TEXT = os.getenv("RAS_PROFILE_PROJECT", "")
RAS_VERSION = os.getenv("RAS_PROFILE_VERSION", "7.0")
PLAN_NUMBER = os.getenv("RAS_PROFILE_PLAN", "02")
PROFILE_ROOT = Path(os.getenv("RAS_PROFILE_OUTPUT", "working/raster_profiles"))
VRT_PATH_TEXT = os.getenv("RAS_PROFILE_VRT", "")
RUN_CONFIGURED_STORE_MAPS = False
RUN_STORE_MAP_PROFILES = False
NATIVE_TIFF_MANIFEST_TEXT = os.getenv(
    "RAS_PROFILE_TIFF_MANIFEST",
    "working/native_tiff_experiments/hr70-acd6ada0-tiff-parallel-v2/patch_manifest.json",
)
RUN_TIFF_PROFILES = False
RUN_NATIVE_TIFF_PROFILES = False

ras = None
if PROJECT_FOLDER_TEXT:
    ras = init_ras_project(Path(PROJECT_FOLDER_TEXT), RAS_VERSION)
else:
    print("Set RAS_PROFILE_PROJECT before running project-dependent cells.")

Run through the unified map API

Use one RasMap.store_all_maps() surface for historic native behavior, selected products, timesteps, or every computed plan. This opt-in cell demonstrates mode="selected"; a plain historic RasMap.store_all_maps(PLAN_NUMBER) call remains native.

Python

configured_map_summary = None
if RUN_CONFIGURED_STORE_MAPS:
    if ras is None:
        raise RuntimeError("Set RAS_PROFILE_PROJECT first.")
    configured_map_summary = RasMap.store_all_maps(
        PLAN_NUMBER,
        mode="selected",
        output_path=PROFILE_ROOT / "configured_maps",
        map_types=("wse", "depth", "velocity"),
        profile="Max",
        performance=StoreMapPerformanceOptions(max_workers=None),
        ras_object=ras,
        ras_version=RAS_VERSION,
        timeout=7200,
    )
configured_map_summary

Estimate memory and worker limits

The estimate uses all checked terrain TIFF dimensions from the .rasmap. Its conservative worker formula is max(floor, 4 × Float32 terrain surfaces + 512 MiB). Selection is bounded by job count, logical CPU count, available physical memory, Windows commit headroom, and the requested worker ceiling. Admission is checked again immediately before each helper launch.

Python

store_map_settings = {
    "serial": StoreMapPerformanceOptions(max_workers=1),
    "auto": StoreMapPerformanceOptions(max_workers=None),
    "auto_cache_256": StoreMapPerformanceOptions(
        max_workers=None,
        gdal_cachemax_mb=256,
    ),
    "two_helpers_cache_256": StoreMapPerformanceOptions(
        max_workers=2,
        gdal_cachemax_mb=256,
    ),
}

estimate_rows = []
if ras is not None:
    for label, settings in store_map_settings.items():
        estimate = RasProcess.estimate_store_map_resources(
            PLAN_NUMBER,
            map_types=("wse", "depth", "velocity"),
            performance=settings,
            ras_object=ras,
        )
        estimate_rows.append({"configuration": label, **estimate.to_dict()})

estimate_df = pd.DataFrame(estimate_rows)
estimate_df[[
    "configuration", "estimated_worker_private_mb", "effective_reserve_mb",
    "cpu_limit", "memory_limit", "selected_workers", "warnings",
]] if not estimate_df.empty else estimate_df

Profile StoreMaps configurations

Each run gets a separate output folder. Reports include elapsed time; peak RSS/private memory; minimum physical and commit headroom; interval CPU utilization; process-attributable throughput, IOPS, and mean request size; host-wide disk busy/latency and network rates; inferred phase summaries; simultaneous helper count; and blockwise pixel hashes plus geospatial metadata for equivalence checks. Run serial first as the baseline.

Python

store_profiles = {}
if RUN_STORE_MAP_PROFILES:
    if ras is None:
        raise RuntimeError("Set RAS_PROFILE_PROJECT first.")
    for label, settings in store_map_settings.items():
        run_root = PROFILE_ROOT / "store_maps" / label
        store_profiles[label] = RasProcess.profile_store_maps(
            PLAN_NUMBER,
            output_path=run_root / "maps",
            report_path=run_root / "profile.json",
            map_types=("wse", "depth", "velocity"),
            performance=settings,
            profile="Max",
            ras_object=ras,
            ras_version=RAS_VERSION,
            timeout=7200,
        )

Python

store_rows = []
for label, result in store_profiles.items():
    store_rows.append({
        "configuration": label,
        "elapsed_seconds": result.elapsed_seconds,
        "peak_private_mb": result.peak_tree_private_mb,
        "peak_rss_mb": result.peak_tree_rss_mb,
        "minimum_available_mb": result.minimum_available_memory_mb,
        "minimum_commit_headroom_mb": result.minimum_commit_headroom_mb,
        "selected_workers": result.resource_estimate.selected_workers,
        "observed_helpers": result.maximum_simultaneous_helpers,
        "average_machine_cpu_percent": result.performance_summary[
            "average_machine_cpu_percent"
        ],
        "peak_tree_cpu_percent": result.performance_summary[
            "peak_tree_cpu_percent"
        ],
        "average_write_mib_per_second": result.performance_summary[
            "average_write_mib_per_second"
        ],
        "average_write_iops": result.performance_summary[
            "average_write_iops"
        ],
        "mean_write_size_kib": (
            result.performance_summary["mean_write_size_bytes"] / 1024
            if result.performance_summary["mean_write_size_bytes"]
            else None
        ),
        "peak_host_disk_busy_percent": result.performance_summary[
            "peak_host_disk_busy_equivalent_percent"
        ],
        "peak_host_network_send_mib_per_second": result.performance_summary[
            "peak_host_network_send_mib_per_second"
        ],
        "write_operations": sum(result.write_operations_by_process.values()),
        "write_bytes": sum(result.write_bytes_by_process.values()),
        "semantic_signatures": semantic_output_signatures(
            result.output_signatures
        ),
    })

store_df = pd.DataFrame(store_rows)
if not store_df.empty:
    baseline_seconds = store_df.loc[
        store_df["configuration"] == "serial", "elapsed_seconds"
    ].iloc[0]
    store_df["speedup_vs_serial"] = baseline_seconds / store_df["elapsed_seconds"]
    store_df["semantics_match_serial"] = (
        store_df["semantic_signatures"]
        == store_df.iloc[0]["semantic_signatures"]
    )
    store_df = store_df.sort_values("elapsed_seconds")
store_df

Python

store_phase_rows = []
for label, result in store_profiles.items():
    for phase, metrics in result.phase_summary.items():
        store_phase_rows.append({
            "configuration": label,
            "phase": phase,
            **metrics,
        })

store_phase_df = pd.DataFrame(store_phase_rows)
store_phase_df

Python

STORE_TIMELINE_CONFIGURATION = next(iter(store_profiles), None)
store_timeline_df = pd.DataFrame()
if STORE_TIMELINE_CONFIGURATION:
    store_timeline_df = pd.DataFrame(
        sample.to_dict()
        for sample in store_profiles[STORE_TIMELINE_CONFIGURATION].samples
    )
    store_timeline_df = store_timeline_df[[
        "elapsed_seconds", "inferred_phase",
        "machine_cpu_percent", "tree_private_mb",
        "read_mib_per_second", "write_mib_per_second",
        "read_iops", "write_iops",
        "host_disk_busy_equivalent_percent",
        "host_network_receive_mib_per_second",
        "host_network_send_mib_per_second",
    ]]
store_timeline_df

Profile TIFF consolidation settings

GeoTiffWriteOptions exposes the GDAL parameters that materially affect writer throughput and memory: tile size, compression and level, predictor, GDAL thread count, GDAL cache, BigTIFF policy, and overview settings. The profiler verifies decoded pixels, dimensions, data types, CRS, transform, NoData, and required overview levels independently of compression and tiling. Include the legacy-compatible setting as the baseline. Codec availability is version-specific; for example, the tested HEC-RAS 7.0 GDAL build does not provide ZSTD. A failed candidate writes a failure JSON report; this matrix records the error and continues with the remaining candidates.

Python

tiff_settings = {
    "legacy_lzw": GeoTiffWriteOptions(),
    "lzw_512_cache_256": GeoTiffWriteOptions(
        compression="LZW", tile_size=512, gdal_cachemax_mb=256,
    ),
    "lzw_512_threads_cache": GeoTiffWriteOptions(
        compression="LZW", predictor=3,
        tile_size=512, gdal_num_threads="ALL_CPUS", gdal_cachemax_mb=256,
        overview_compression="LZW",
    ),
    "lzw_1024_threads_cache": GeoTiffWriteOptions(
        compression="LZW", predictor=3,
        tile_size=1024, gdal_num_threads="ALL_CPUS", gdal_cachemax_mb=512,
        overview_compression="LZW",
    ),
}

tiff_profiles = {}
tiff_failures = []
if RUN_TIFF_PROFILES:
    if not VRT_PATH_TEXT:
        raise RuntimeError("Set RAS_PROFILE_VRT to the source VRT.")
    for label, settings in tiff_settings.items():
        run_root = PROFILE_ROOT / "vrt_to_tiff" / label
        try:
            tiff_profiles[label] = RasTerrain.profile_vrt_to_tiff(
                VRT_PATH_TEXT,
                output_path=run_root / "terrain.tif",
                report_path=run_root / "profile.json",
                write_options=settings,
                hecras_version=RAS_VERSION,
            )
        except Exception as exc:
            tiff_failures.append({
                "configuration": label,
                "error_type": type(exc).__name__,
                "error": str(exc),
                "report_path": str(run_root / "profile.json"),
            })

tiff_failure_df = pd.DataFrame(tiff_failures)
tiff_failure_df

Python

tiff_rows = []
for label, result in tiff_profiles.items():
    counters = result.process_counters
    tiff_rows.append({
        "configuration": label,
        "elapsed_seconds": result.elapsed_seconds,
        "peak_private_mb": result.peak_tree_private_mb,
        "minimum_commit_headroom_mb": result.minimum_commit_headroom_mb,
        "average_machine_cpu_percent": result.performance_summary[
            "average_machine_cpu_percent"
        ],
        "average_write_mib_per_second": result.performance_summary[
            "average_write_mib_per_second"
        ],
        "average_write_iops": result.performance_summary[
            "average_write_iops"
        ],
        "mean_write_size_kib": (
            result.performance_summary["mean_write_size_bytes"] / 1024
            if result.performance_summary["mean_write_size_bytes"]
            else None
        ),
        "peak_host_disk_busy_percent": result.performance_summary[
            "peak_host_disk_busy_equivalent_percent"
        ],
        "write_operations": sum(
            int(values["write_operations"]) for values in counters.values()
        ),
        "write_bytes": sum(int(values["write_bytes"]) for values in counters.values()),
        "file_size_bytes": result.output_signature["file_size_bytes"],
        "semantic_signature": json.dumps(
            semantic_raster_signature(result.output_signature), sort_keys=True
        ),
    })

tiff_df = pd.DataFrame(tiff_rows)
if not tiff_df.empty:
    baseline = tiff_df.iloc[0]
    tiff_df["speedup_vs_legacy"] = baseline["elapsed_seconds"] / tiff_df["elapsed_seconds"]
    tiff_df["semantics_match_legacy"] = (
        tiff_df["semantic_signature"] == baseline["semantic_signature"]
    )
    tiff_df = tiff_df.sort_values("elapsed_seconds")
tiff_df

Python

tiff_phase_rows = []
for label, result in tiff_profiles.items():
    for phase, metrics in result.phase_summary.items():
        tiff_phase_rows.append({
            "configuration": label,
            "phase": phase,
            **metrics,
        })

tiff_phase_df = pd.DataFrame(tiff_phase_rows)
tiff_phase_df

Optional copied-TiffAssist research matrix

This optional matrix profiles the version/hash-pinned copied HEC-RAS 7.0 writer. It does not modify the installed runtime and is inert unless the benchmark launches the copied helper with its explicit opt-in environment. The matrix compares the installed original, batching with workers off, and 1/2/4/8 independent tile preparation and Deflate workers. Queue depth 2 bounds tile memory. Keep the switch false unless the copied assembly has already been built with prepare_copied_assembly.py, use a disposable project copy, and repeat promising cases before choosing a setting. These controls are experimental and are not part of the stable RasProcess API.

Spring River showed why profiling matters: the isolated writer scaled strongly, but the real map producer could not feed more than one compression worker continuously. More workers are not automatically faster.

Python

native_tiff_df = pd.DataFrame()
if RUN_NATIVE_TIFF_PROFILES:
    if not PROJECT_FOLDER_TEXT:
        raise RuntimeError("Set RAS_PROFILE_PROJECT first.")
    manifest_path = Path(NATIVE_TIFF_MANIFEST_TEXT)
    if not manifest_path.is_file():
        raise FileNotFoundError(
            "Build the copied assembly before profiling: "
            f"{manifest_path}"
        )
    native_root = PROFILE_ROOT / "native_tiff"
    command = [
        sys.executable,
        "scripts/benchmarks/native_tiff_experiments/benchmark_copied_tiffassist.py",
        "store-map",
        "--manifest", str(manifest_path),
        "--output-root", str(native_root),
        "--project-folder", PROJECT_FOLDER_TEXT,
        "--plan-number", PLAN_NUMBER,
        "--ras-version", RAS_VERSION,
        "--map-type", "Depth",
        "--profile", "Max",
        "--store-map-batch-bytes", "2097152",
        "--store-map-statistics-modes", "serial",
        "--store-map-pipeline-workers", "0,1,2,4,8",
        "--store-map-pipeline-queue-depth", "2",
        "--timeout", "7200",
    ]
    subprocess.run(command, check=True)
    native_tiff_df = pd.read_csv(native_root / "summary.csv")

native_columns = [
    "id", "seconds", "effective_logical_cpus",
    "peak_private_mib", "process_write_iops",
    "process_mean_write_bytes", "pipeline_workers",
    "pipeline_queue_depth", "pipeline_wall_seconds",
    "prepare_seconds", "deflate_seconds",
    "raw_commit_seconds", "maximum_owned_tiles",
    "equivalent", "result_hdf_preserved", "accepted",
]
native_tiff_df[native_columns] if not native_tiff_df.empty else native_tiff_df

Decision figures and HTML report

These figures put elapsed time beside memory, CPU, throughput, and IOPS so a fast configuration that exhausts headroom is not selected accidentally. The next cell writes a self-contained SVG-based HTML report when at least one profiling matrix has results.

Python

decision_figures = {}

if not store_df.empty:
    ordered = store_df.sort_values("elapsed_seconds")
    fig, axes = plt.subplots(1, 2, figsize=(12, 4.2))
    ordered.plot.bar(x="configuration", y="elapsed_seconds", ax=axes[0], legend=False, color="#0b6e69")
    ordered.plot.scatter(x="peak_private_mb", y="elapsed_seconds", s=90, ax=axes[1], color="#d97706")
    for _, row in ordered.iterrows():
        axes[1].annotate(row["configuration"], (row["peak_private_mb"], row["elapsed_seconds"]), xytext=(5, 4), textcoords="offset points", fontsize=8)
    axes[0].set_ylabel("wall seconds (lower is better)")
    axes[1].set_xlabel("peak private memory MiB")
    axes[1].set_ylabel("wall seconds")
    fig.suptitle("StoreMap time and memory tradeoff")
    fig.tight_layout()
    decision_figures["StoreMap time and memory tradeoff"] = fig

if not store_timeline_df.empty:
    fig, axes = plt.subplots(3, 1, figsize=(12, 8), sharex=True)
    timeline = store_timeline_df
    axes[0].plot(timeline["elapsed_seconds"], timeline["machine_cpu_percent"], color="#2563eb")
    axes[0].set_ylabel("machine CPU %")
    axes[1].plot(timeline["elapsed_seconds"], timeline["tree_private_mb"] / 1024, color="#9333ea")
    axes[1].set_ylabel("private GiB")
    axes[2].plot(timeline["elapsed_seconds"], timeline["read_iops"], label="read IOPS", color="#0b6e69")
    axes[2].plot(timeline["elapsed_seconds"], timeline["write_iops"], label="write IOPS", color="#d97706")
    axes[2].set_ylabel("process IOPS")
    axes[2].set_xlabel("elapsed seconds")
    axes[2].legend()
    fig.suptitle(f"Resource timeline: {STORE_TIMELINE_CONFIGURATION}")
    fig.tight_layout()
    decision_figures["StoreMap CPU, memory, and IOPS timeline"] = fig

if not tiff_df.empty:
    fig, axes = plt.subplots(1, 2, figsize=(12, 4.2))
    ordered = tiff_df.sort_values("elapsed_seconds")
    ordered.plot.bar(x="configuration", y="elapsed_seconds", ax=axes[0], legend=False, color="#0b6e69")
    ordered.plot.bar(x="configuration", y="peak_private_mb", ax=axes[1], legend=False, color="#9333ea")
    axes[0].set_ylabel("wall seconds")
    axes[1].set_ylabel("peak private MiB")
    fig.suptitle("VRT-to-TIFF time and memory")
    fig.tight_layout()
    decision_figures["VRT-to-TIFF time and memory"] = fig

if not native_tiff_df.empty:
    fig, axes = plt.subplots(1, 2, figsize=(12, 4.2))
    native_tiff_df.plot.bar(x="id", y="seconds", ax=axes[0], legend=False, color="#0b6e69")
    native_tiff_df.plot.scatter(x="process_write_iops", y="seconds", s=90, ax=axes[1], color="#d97706")
    axes[0].set_ylabel("wall seconds")
    axes[1].set_xlabel("process write IOPS")
    axes[1].set_ylabel("wall seconds")
    fig.suptitle("Copied TiffAssist latency and write pressure")
    fig.tight_layout()
    decision_figures["Copied TiffAssist latency and write pressure"] = fig

for figure in decision_figures.values():
    display(figure)

Python

def figure_as_svg(figure):
    buffer = io.StringIO()
    figure.savefig(buffer, format="svg", bbox_inches="tight")
    svg = buffer.getvalue()
    return svg[svg.index("<svg"):]

report_tables = {
    "StoreMap configurations": store_df.drop(columns=["semantic_signatures"], errors="ignore"),
    "StoreMap phase summary": store_phase_df,
    "VRT-to-TIFF configurations": tiff_df.drop(columns=["semantic_signature"], errors="ignore"),
    "Copied TiffAssist configurations": native_tiff_df,
}
report_tables = {name: table for name, table in report_tables.items() if not table.empty}
PROFILE_REPORT_HTML = PROFILE_ROOT / "raster_performance_profile.html"
if decision_figures or report_tables:
    sections = []
    for title, figure in decision_figures.items():
        sections.append(f"<section><h2>{html_lib.escape(title)}</h2>{figure_as_svg(figure)}</section>")
    for title, table in report_tables.items():
        sections.append(f"<section><h2>{html_lib.escape(title)}</h2>{table.to_html(index=False, escape=True)}</section>")
    document = f"""<!doctype html><html lang='en'><head><meta charset='utf-8'><meta name='viewport' content='width=device-width,initial-scale=1'><title>RAS raster performance profile</title><style>body{{max-width:1200px;margin:auto;padding:28px;font:16px/1.5 system-ui;color:#152329;background:#f5f7f6}}section{{background:white;border:1px solid #d5dfdc;border-radius:12px;padding:20px;margin:18px 0;overflow-x:auto}}table{{border-collapse:collapse;width:100%}}th,td{{border-bottom:1px solid #d5dfdc;padding:7px;text-align:right}}th:first-child,td:first-child{{text-align:left}}svg{{max-width:100%;height:auto}}</style></head><body><h1>RAS raster performance profile</h1><p>Choose the fastest repeatable configuration that preserves raster semantics and physical/commit memory headroom. Process I/O includes cached and SMB transfers; host counters are machine-wide.</p>{''.join(sections)}</body></html>"""
    PROFILE_REPORT_HTML.parent.mkdir(parents=True, exist_ok=True)
    PROFILE_REPORT_HTML.write_text(document, encoding="utf-8")
    print(f"Wrote {PROFILE_REPORT_HTML.resolve()}")
else:
    print("Run at least one profiling matrix to generate the HTML report.")

Interpret utilization and storage counters

tree_cpu_percent uses 100% for each fully occupied logical CPU; machine_cpu_percent divides that value by the machine's logical CPU count. Process I/O counters attribute file transfers to the Python/native/GDAL process tree, including cached and SMB I/O. Host disk and network counters show machine-wide load and can include unrelated activity. Phase labels are inferred from active descendants; when parallel helpers overlap, a phase row can contain mixed whole-tree work. A high process write IOPS rate with a small mean write size supports a small-write bottleneck; low CPU plus higher elapsed time on network storage supports latency or serialization; high CPU near the available-core ceiling supports a compute or compression bottleneck. Exact physical-device queue depth and request latency require ETW, Performance Monitor, or storage/server telemetry outside this portable profiler.

Select a setting for this machine

Reject configurations with non-equivalent pixels or unacceptable physical/commit headroom. Among the remaining rows, prefer the fastest repeatable result rather than assuming more workers, more GDAL threads, or a larger cache will help. Repeat important candidates at least three times because antivirus activity, filesystem cache state, compression ratio, and concurrent workloads can change timings. Do not use worker_memory_override_mb with enforced admission; an override requires an explicit memory_policy="warn" or "ignore" acknowledgement.