Skip to content

AORC Precipitation Catalog for HEC-RAS Rain-on-Grid Models

This notebook demonstrates a complete workflow for using NOAA's Analysis of Record for Calibration (AORC) gridded precipitation data with HEC-RAS 2D rain-on-grid models.

Example Project: BaldEagleCrkMulti2D

Workflow

  1. Extract example project and create a labeled working copy
  2. Get project bounds from geometry HDF (with buffer)
  3. Generate a storm catalog from historical AORC data
  4. Download precipitation, export hyetographs, and create HEC-RAS plans
  5. Execute all plans in parallel (2 cores, 3 workers)
  6. Extract and compare results

Data Export

All precipitation data is exported to the Precipitation/ subfolder: - storm_catalog.csv - Complete storm catalog with metadata - storm_YYYYMMDD.nc - NetCDF precipitation files - hyetographs/ - PNG plots of each storm's precipitation - storm_catalog_summary.png - Overview plot of all storms

AORC Dataset Overview

  • Coverage: CONUS (1979-present), Alaska (1981-present)
  • Resolution: ~800 meters, hourly timesteps
  • Format: Cloud-optimized Zarr on AWS S3
  • Access: Anonymous (no authentication required)
Python
# Install dependencies (uncomment if needed)
# !pip install ras-commander[precip]  # Includes xarray, zarr, s3fs, netCDF4, rioxarray
Python
# =============================================================================
# DEVELOPMENT MODE TOGGLE
# =============================================================================
USE_LOCAL_SOURCE = False  # <-- TOGGLE THIS

if USE_LOCAL_SOURCE:
    import sys
    from pathlib import Path
    local_path = str(Path.cwd().parent)
    if local_path not in sys.path:
        sys.path.insert(0, local_path)
    print(f"📁 LOCAL SOURCE MODE: Loading from {local_path}/ras_commander")
else:
    print("📦 PIP PACKAGE MODE: Loading installed ras-commander")

# Import ras-commander
from ras_commander import init_ras_project, RasExamples, RasPlan, RasUnsteady, RasCmdr, RasUtils
from ras_commander.precip import PrecipAorc
from ras_commander.hdf.HdfProject import HdfProject

# Additional imports
import os
import shutil
import glob
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import xarray as xr

# Verify which version loaded
import ras_commander
print(f"✓ Loaded: {ras_commander.__file__}")
Text Only
📦 PIP PACKAGE MODE: Loading installed ras-commander


✓ Loaded: C:\Users\bill\.config\superpowers\worktrees\ras-commander-pr-docs\codex-pr251-timeseries-alignment\ras_commander\__init__.py

Parameters

Configure these values to customize the notebook for your project.

Python
# =============================================================================
# PARAMETERS - Edit these to customize the notebook
# =============================================================================
from pathlib import Path

# Project Configuration
PROJECT_NAME = "BaldEagleCrkMulti2D"  # Example project to extract
RAS_VERSION = "7.0"               # HEC-RAS version (6.3, 6.5, 6.6, etc.)
SUFFIX = "901"                    # Notebook identifier for project folder

# AORC Settings
ONLINE = True                     # Enable network requests

print(f"Project: {PROJECT_NAME}, Version: {RAS_VERSION}")
Text Only
Project: BaldEagleCrkMulti2D, Version: 7.0

AORC Catalog Verification

After generating precipitation catalog:

  • [ ] Storm events span simulation period + warmup time
  • [ ] Total depths comparable to Atlas 14 estimates (within 20%)
  • [ ] Inter-event separation prevents storm merging (8+ hours typical)
  • [ ] Spatial coverage includes entire watershed (50% buffer minimum)

Data Quality Checks:

Python
# Example: Compare catalog storms to Atlas 14
# (See NOAA Atlas 14: https://hdsc.nws.noaa.gov/hdsc/pfds/)
for storm in storm_catalog.iterrows():
    total_depth = storm['total_depth_in']
    duration = storm['duration_hours']
    # Compare to Atlas 14 for this duration/return period
    # Document if >20% different (grid vs point)

References: - NOAA Atlas 14: Precipitation-Frequency Atlas - AORC Documentation

Step 1: Extract Example Project and Create Working Copy

We extract the example project into a notebook-specific work root and create a labeled copy for our AORC analysis. Note: Existing notebook 901 AORC folders are removed with retry cleanup to ensure repeatability without touching notebook 900 worker folders.

Python
# Configuration
YEAR = 2020  # Year to analyze
TEMPLATE_PLAN = "06"  # Template plan with gridded precipitation and infiltration enabled
NUM_CORES = 2  # Cores per HEC-RAS instance
MAX_WORKERS = 3  # Parallel worker processes

# Keep this notebook isolated from notebook 900 parallel worker folders.
notebook_work_root = Path.cwd() / "working" / "example_projects_901_aorc_catalog"
notebook_work_root.mkdir(parents=True, exist_ok=True)
print(f"Notebook work root: {notebook_work_root}")

# Extract base example project into the notebook-specific work root.
base_project = RasExamples.extract_project(
    PROJECT_NAME,
    output_path=notebook_work_root,
    suffix=SUFFIX,
)
print(f"Base project extracted to: {base_project}")

# Clean up any existing AORC folders (for repeatability)
print("\nCleaning up existing AORC folders...")
aorc_pattern = f"{PROJECT_NAME}_AORC_*_{SUFFIX}*"
existing_folders = list(base_project.parent.glob(aorc_pattern))
removed_folders = 0
for folder in existing_folders:
    if folder.is_dir():
        print(f"  Removing: {folder.name}")
        if RasUtils.remove_with_retry(folder, ras_object=None):
            removed_folders += 1
        else:
            raise PermissionError(f"Unable to remove stale notebook 901 folder: {folder}")
print(f"  Removed {removed_folders} existing folders")

# Create labeled working copy (include SUFFIX for concurrent testing)
working_folder = base_project.parent / f"{PROJECT_NAME}_AORC_{YEAR}_{SUFFIX}"
print(f"\nCreating working copy: {working_folder}")
shutil.copytree(base_project, working_folder)

# Initialize project
ras = init_ras_project(working_folder, RAS_VERSION)
print(f"\nProject: {ras.project_name}")
print(f"Location: {ras.project_folder}")
print(f"Plans: {len(ras.plan_df)}")
Text Only
Notebook work root: C:\Users\bill\.config\superpowers\worktrees\ras-commander-pr-docs\codex-pr251-timeseries-alignment\examples\working\example_projects_901_aorc_catalog


2026-07-08 23:17:36 - ras_commander.RasExamples - INFO - Successfully extracted project 'BaldEagleCrkMulti2D' to working\example_projects_901_aorc_catalog\BaldEagleCrkMulti2D_901


Base project extracted to: C:\Users\bill\.config\superpowers\worktrees\ras-commander-pr-docs\codex-pr251-timeseries-alignment\examples\working\example_projects_901_aorc_catalog\BaldEagleCrkMulti2D_901

Cleaning up existing AORC folders...
  Removing: BaldEagleCrkMulti2D_AORC_2020_901


  Removed 1 existing folders

Creating working copy: C:\Users\bill\.config\superpowers\worktrees\ras-commander-pr-docs\codex-pr251-timeseries-alignment\examples\working\example_projects_901_aorc_catalog\BaldEagleCrkMulti2D_AORC_2020_901


2026-07-08 23:17:40 - ras_commander.RasPrj - INFO - ras-commander v0.98.2 | An open-source project of CLB Engineering Corporation (https://clbengineering.com/) | Docs: https://rascommander.info | GitHub: https://github.com/gpt-cmdr/ras-commander


2026-07-08 23:17:40 - ras_commander.RasPrj - INFO - Project initialized: BaldEagleDamBrk


2026-07-08 23:17:40 - ras_commander.RasPrj - INFO - Using HEC-RAS version 7.0


2026-07-08 23:17:40 - ras_commander.RasPrj - INFO - 
═══════════════════════════════════════════════════════════════════════
ras-commander | HEC-RAS Automation Library
Docs: https://rascommander.info/
Repo: https://github.com/gpt-cmdr/ras-commander
LLM agents: https://rascommander.info/llms.txt
═══════════════════════════════════════════════════════════════════════

PROJECT DATAFRAMES (single source of truth — use these, not file globbing):
  ras.plan_df        Plans, HDF paths, geometry/flow associations
  ras.geom_df        Geometry files and HDF preprocessor paths
  ras.flow_df        Steady flow files
  ras.unsteady_df    Unsteady flow files and configurations
  ras.boundaries_df  Boundary conditions (type, name, location)
  ras.results_df     Lightweight HDF results summaries
  ras.rasmap_df      RASMapper layers, terrain, land cover paths

KEY APIS (static classes — call directly, never instantiate):
  Execution:    RasCmdr.compute_plan() / compute_parallel() / compute_test_mode()
  Plan Files:   RasPlan.clone_plan() / clone_geom() / set_geom()
  Unsteady:     RasUnsteady — IC/BC management, gate openings, precipitation
  Geometry:     GeomCrossSection, GeomBridge, GeomStorage, GeomLateral, GeomMesh
  HDF Results:  HdfResultsPlan.get_wse() / get_compute_messages()
                HdfResultsMesh.get_mesh_max_ws() / get_mesh_cells_timeseries()
                HdfMesh.get_mesh_cell_points()
  QA/QC:        RasCheck.run_check() / RasFixit (geometry repair)
  DSS:          RasDss.get_timeseries() / check_pathname()
  USGS:         UsgsGaugeSpatial, GaugeMatcher, RasUsgsBoundaryGeneration
  Precipitation: StormGenerator, Atlas14Storm, PrecipAorc, Atlas14Variance
  Terrain:      RasTerrain.create_terrain_hdf() / RasTerrainMod

MULTI-PROJECT: Pass ras_object= to all API calls when using local RasPrj instances.

EXAMPLES: 100+ notebooks in examples/ (100s=execution, 200s=geometry, 300s=unsteady,
  400s=HDF results, 500s=remote, 800s=QA/QC, 900s=data integration).
  Review relevant notebooks before assembling new workflows.

PLATFORM: Most HEC-RAS operations require Windows. Linux/Wine support for
  headless execution, data access, geometry modification, and preprocessing
  is available via RasProcess (HEC-RAS 6.6+). See ras_commander/RasProcess.py.
  Remote distributed execution: ras_commander/remote/ (PsExec, Docker, SSH, cloud).
═══════════════════════════════════════════════════════════════════════



Project: BaldEagleDamBrk
Location: C:\Users\bill\.config\superpowers\worktrees\ras-commander-pr-docs\codex-pr251-timeseries-alignment\examples\working\example_projects_901_aorc_catalog\BaldEagleCrkMulti2D_AORC_2020_901
Plans: 11

Step 2: Get Project Bounds from Geometry HDF

Use HdfProject.get_project_bounds_latlon() to extract proper bounds with buffering. This handles CRS transformation and ensures precipitation coverage.

Python
# Get geometry HDF path from template plan
# The 'Geom File' column contains just the number (e.g., "09"),
# so we need to add the "g" prefix for the geometry file extension
_matched = ras.plan_df[ras.plan_df['plan_number'] == TEMPLATE_PLAN]
if _matched.empty:
    available = ras.plan_df['plan_number'].tolist()
    raise ValueError(f"TEMPLATE_PLAN '{TEMPLATE_PLAN}' not found in project. "
                     f"Available plan numbers: {available}")
template_row = _matched.iloc[0]
geom_number = template_row['Geom File']
geom_hdf = ras.project_folder / f"{ras.project_name}.g{geom_number}.hdf"

print(f"Template plan {TEMPLATE_PLAN} uses geometry: g{geom_number}")
print(f"Geometry HDF: {geom_hdf}")
print(f"Exists: {geom_hdf.exists()}")

# Get bounds with 50% buffer (default)
# This properly handles CRS transformation and includes 2D mesh, 1D elements, and storage areas
bounds = HdfProject.get_project_bounds_latlon(
    geom_hdf,
    buffer_percent=50.0,  # 50% buffer ensures full coverage after reprojection
    include_1d=True,
    include_2d=True,
    include_storage=True
)

west, south, east, north = bounds
print(f"\nProject Bounds (WGS84 with 50% buffer):")
print(f"  West:  {west:.4f}")
print(f"  South: {south:.4f}")
print(f"  East:  {east:.4f}")
print(f"  North: {north:.4f}")
print(f"  Size:  {east-west:.4f} x {north-south:.4f} degrees")
Text Only
Template plan 06 uses geometry: g09
Geometry HDF: C:\Users\bill\.config\superpowers\worktrees\ras-commander-pr-docs\codex-pr251-timeseries-alignment\examples\working\example_projects_901_aorc_catalog\BaldEagleCrkMulti2D_AORC_2020_901\BaldEagleDamBrk.g09.hdf
Exists: True

Project Bounds (WGS84 with 50% buffer):
  West:  -77.8672
  South: 40.9044
  East:  -77.2192
  North: 41.2408
  Size:  0.6480 x 0.3364 degrees

Step 3: Generate Storm Catalog

Analyze AORC data to identify all significant precipitation events for the year.

Python
# Generate storm catalog
storm_catalog = PrecipAorc.get_storm_catalog(
    bounds=bounds,
    year=YEAR,
    inter_event_hours=8.0,      # USGS standard for storm separation
    min_depth_inches=0.75,      # Minimum significant precipitation
    buffer_hours=48             # Simulation warmup buffer
)

print(f"\nStorm Catalog: {len(storm_catalog)} events for {YEAR}")
print("="*90)
print(storm_catalog[['storm_id', 'start_time', 'end_time', 'total_depth_in', 
                     'peak_intensity_in_hr', 'duration_hours', 'rank']].to_string(index=False))
Text Only
2026-07-08 23:17:40 - ras_commander.precip.PrecipAorc - INFO - Generating AORC storm catalog for 2020: bounds W=-77.8672, S=40.9044, E=-77.2192, N=41.2408; inter_event=8.0h, min_depth=0.75in, buffer=48h


2026-07-08 23:17:48 - ras_commander.precip.PrecipAorc - INFO - Storm catalog complete: 12 storms; depth 0.76-2.72 in; largest 2020-12-24 12:00 (2.72 in)



Storm Catalog: 12 events for 2020
==========================================================================================
 storm_id          start_time            end_time  total_depth_in  peak_intensity_in_hr  duration_hours  rank
        1 2020-02-07 09:00:00 2020-02-07 17:00:00           1.023                 0.287               9     9
        2 2020-03-19 01:00:00 2020-03-19 10:00:00           0.763                 0.209              10    12
        3 2020-03-28 08:00:00 2020-03-29 11:00:00           1.530                 0.156              28     6
        4 2020-04-17 17:00:00 2020-04-18 08:00:00           0.777                 0.093              16    11
        5 2020-04-30 04:00:00 2020-05-01 00:00:00           2.115                 0.203              21     2
        6 2020-05-28 12:00:00 2020-05-30 01:00:00           1.188                 0.284              38     8
        7 2020-06-03 06:00:00 2020-06-03 22:00:00           0.846                 0.360              17    10
        8 2020-09-29 11:00:00 2020-09-30 08:00:00           1.278                 0.177              22     7
        9 2020-10-29 07:00:00 2020-10-30 18:00:00           1.735                 0.127              36     5
       10 2020-11-11 12:00:00 2020-11-11 22:00:00           1.989                 0.448              11     3
       11 2020-12-16 17:00:00 2020-12-17 09:00:00           1.834                 0.174              17     4
       12 2020-12-24 12:00:00 2020-12-25 12:00:00           2.720                 0.202              25     1

Step 4: Download AORC Data and Export Precipitation Records

Download precipitation data, export storm catalog CSV, and generate hyetograph plots for documentation.

Python
# Create precipitation folder structure
precip_folder = ras.project_folder / "Precipitation"
precip_folder.mkdir(exist_ok=True)
hyetograph_folder = precip_folder / "hyetographs"
hyetograph_folder.mkdir(exist_ok=True)

# Export storm catalog to CSV
catalog_csv = precip_folder / "storm_catalog.csv"
storm_catalog.to_csv(catalog_csv, index=False)
print(f"Storm catalog exported to: {catalog_csv}")

# Download AORC data for all storms
print(f"\nDownloading AORC precipitation data for {len(storm_catalog)} storms...")
print("="*70)

precip_files = {}
for idx, storm in storm_catalog.iterrows():
    storm_id = storm['storm_id']
    date_str = storm['start_time'].strftime('%Y%m%d')
    precip_file = precip_folder / f"storm_{date_str}.nc"

    if not precip_file.exists():
        print(f"  Storm {storm_id:2d}: {storm['start_time'].strftime('%b %d')} - Downloading...")
        PrecipAorc.download(
            bounds=bounds,
            start_time=storm['sim_start'],
            end_time=storm['sim_end'],
            output_path=precip_file,
            target_crs="EPSG:5070",
            resolution=2000.0
        )
    else:
        print(f"  Storm {storm_id:2d}: {storm['start_time'].strftime('%b %d')} - Already downloaded")

    precip_files[storm_id] = precip_file

print(f"\nDownloaded {len(precip_files)} precipitation files to {precip_folder}")
Text Only
2026-07-08 23:17:48 - ras_commander.precip.PrecipAorc - INFO - Downloading AORC APCP_surface: 2020-02-05 09:00:00 to 2020-02-09 17:00:00, bounds W=-77.8672, S=40.9044, E=-77.2192, N=41.2408


Storm catalog exported to: C:\Users\bill\.config\superpowers\worktrees\ras-commander-pr-docs\codex-pr251-timeseries-alignment\examples\working\example_projects_901_aorc_catalog\BaldEagleCrkMulti2D_AORC_2020_901\Precipitation\storm_catalog.csv

Downloading AORC precipitation data for 12 storms...
======================================================================
  Storm  1: Feb 07 - Downloading...


2026-07-08 23:17:52 - ras_commander.precip.PrecipAorc - INFO - AORC download complete: storm_20200207.nc (0.2 MB, 120 timesteps)


2026-07-08 23:17:52 - ras_commander.precip.PrecipAorc - INFO - Downloading AORC APCP_surface: 2020-03-17 01:00:00 to 2020-03-21 10:00:00, bounds W=-77.8672, S=40.9044, E=-77.2192, N=41.2408


  Storm  2: Mar 19 - Downloading...


2026-07-08 23:17:56 - ras_commander.precip.PrecipAorc - INFO - AORC download complete: storm_20200319.nc (0.2 MB, 120 timesteps)


2026-07-08 23:17:56 - ras_commander.precip.PrecipAorc - INFO - Downloading AORC APCP_surface: 2020-03-26 08:00:00 to 2020-03-31 11:00:00, bounds W=-77.8672, S=40.9044, E=-77.2192, N=41.2408


  Storm  3: Mar 28 - Downloading...


2026-07-08 23:18:00 - ras_commander.precip.PrecipAorc - INFO - AORC download complete: storm_20200328.nc (0.2 MB, 144 timesteps)


2026-07-08 23:18:00 - ras_commander.precip.PrecipAorc - INFO - Downloading AORC APCP_surface: 2020-04-15 17:00:00 to 2020-04-20 08:00:00, bounds W=-77.8672, S=40.9044, E=-77.2192, N=41.2408


  Storm  4: Apr 17 - Downloading...


2026-07-08 23:18:04 - ras_commander.precip.PrecipAorc - INFO - AORC download complete: storm_20200417.nc (0.2 MB, 144 timesteps)


2026-07-08 23:18:04 - ras_commander.precip.PrecipAorc - INFO - Downloading AORC APCP_surface: 2020-04-28 04:00:00 to 2020-05-03 00:00:00, bounds W=-77.8672, S=40.9044, E=-77.2192, N=41.2408


  Storm  5: Apr 30 - Downloading...


2026-07-08 23:18:09 - ras_commander.precip.PrecipAorc - INFO - AORC download complete: storm_20200430.nc (0.2 MB, 144 timesteps)


2026-07-08 23:18:09 - ras_commander.precip.PrecipAorc - INFO - Downloading AORC APCP_surface: 2020-05-26 12:00:00 to 2020-06-01 01:00:00, bounds W=-77.8672, S=40.9044, E=-77.2192, N=41.2408


  Storm  6: May 28 - Downloading...


2026-07-08 23:18:13 - ras_commander.precip.PrecipAorc - INFO - AORC download complete: storm_20200528.nc (0.3 MB, 168 timesteps)


2026-07-08 23:18:13 - ras_commander.precip.PrecipAorc - INFO - Downloading AORC APCP_surface: 2020-06-01 06:00:00 to 2020-06-05 22:00:00, bounds W=-77.8672, S=40.9044, E=-77.2192, N=41.2408


  Storm  7: Jun 03 - Downloading...


2026-07-08 23:18:17 - ras_commander.precip.PrecipAorc - INFO - AORC download complete: storm_20200603.nc (0.2 MB, 120 timesteps)


2026-07-08 23:18:17 - ras_commander.precip.PrecipAorc - INFO - Downloading AORC APCP_surface: 2020-09-27 11:00:00 to 2020-10-02 08:00:00, bounds W=-77.8672, S=40.9044, E=-77.2192, N=41.2408


  Storm  8: Sep 29 - Downloading...


2026-07-08 23:18:21 - ras_commander.precip.PrecipAorc - INFO - AORC download complete: storm_20200929.nc (0.2 MB, 144 timesteps)


2026-07-08 23:18:21 - ras_commander.precip.PrecipAorc - INFO - Downloading AORC APCP_surface: 2020-10-27 07:00:00 to 2020-11-01 18:00:00, bounds W=-77.8672, S=40.9044, E=-77.2192, N=41.2408


  Storm  9: Oct 29 - Downloading...


2026-07-08 23:18:26 - ras_commander.precip.PrecipAorc - INFO - AORC download complete: storm_20201029.nc (0.2 MB, 144 timesteps)


2026-07-08 23:18:26 - ras_commander.precip.PrecipAorc - INFO - Downloading AORC APCP_surface: 2020-11-09 12:00:00 to 2020-11-13 22:00:00, bounds W=-77.8672, S=40.9044, E=-77.2192, N=41.2408


  Storm 10: Nov 11 - Downloading...


2026-07-08 23:18:32 - ras_commander.precip.PrecipAorc - INFO - AORC download complete: storm_20201111.nc (0.2 MB, 120 timesteps)


2026-07-08 23:18:32 - ras_commander.precip.PrecipAorc - INFO - Downloading AORC APCP_surface: 2020-12-14 17:00:00 to 2020-12-19 09:00:00, bounds W=-77.8672, S=40.9044, E=-77.2192, N=41.2408


  Storm 11: Dec 16 - Downloading...


2026-07-08 23:18:37 - ras_commander.precip.PrecipAorc - INFO - AORC download complete: storm_20201216.nc (0.2 MB, 144 timesteps)


2026-07-08 23:18:37 - ras_commander.precip.PrecipAorc - INFO - Downloading AORC APCP_surface: 2020-12-22 12:00:00 to 2020-12-27 12:00:00, bounds W=-77.8672, S=40.9044, E=-77.2192, N=41.2408


  Storm 12: Dec 24 - Downloading...


2026-07-08 23:18:42 - ras_commander.precip.PrecipAorc - INFO - AORC download complete: storm_20201224.nc (0.2 MB, 144 timesteps)



Downloaded 12 precipitation files to C:\Users\bill\.config\superpowers\worktrees\ras-commander-pr-docs\codex-pr251-timeseries-alignment\examples\working\example_projects_901_aorc_catalog\BaldEagleCrkMulti2D_AORC_2020_901\Precipitation
Python
# Generate and save individual hyetograph plots for each storm
print(f"Generating hyetograph plots for {len(storm_catalog)} storms...")
print("="*70)

for idx, storm in storm_catalog.iterrows():
    storm_id = storm['storm_id']
    date_str = storm['start_time'].strftime('%Y%m%d')
    precip_file = precip_files[storm_id]

    # Load precipitation data
    ds = xr.open_dataset(precip_file)
    da = ds['APCP_surface']
    hourly_mean = da.mean(dim=['x', 'y']).values
    times = pd.to_datetime(da.time.values)

    # Create figure
    fig, ax = plt.subplots(figsize=(12, 5))
    ax.bar(times, hourly_mean, width=0.03, color='steelblue', alpha=0.8, edgecolor='darkblue', linewidth=0.5)

    # Add storm info
    ax.set_title(f"Storm {storm_id}: {storm['start_time'].strftime('%B %d, %Y')}\n"
                 f"Total: {storm['total_depth_in']:.2f} in | Peak: {storm['peak_intensity_in_hr']:.3f} in/hr | "
                 f"Duration: {storm['duration_hours']:.0f} hours", fontsize=12, fontweight='bold')
    ax.set_xlabel('Date/Time', fontsize=11)
    ax.set_ylabel('Precipitation Rate (mm/hr)', fontsize=11)
    ax.tick_params(axis='x', rotation=45)
    ax.grid(True, alpha=0.3)

    # Add statistics box
    stats_text = (f"Max Rate: {hourly_mean.max():.2f} mm/hr\n"
                  f"Total: {hourly_mean.sum():.1f} mm\n"
                  f"Timesteps: {len(hourly_mean)}")
    ax.text(0.98, 0.95, stats_text, transform=ax.transAxes, fontsize=9,
            verticalalignment='top', horizontalalignment='right',
            bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.8))

    plt.tight_layout()

    # Save plot
    plot_path = hyetograph_folder / f"storm_{date_str}_hyetograph.png"
    plt.savefig(plot_path, dpi=150, bbox_inches='tight')
    plt.close()

    ds.close()
    print(f"  Storm {storm_id:2d}: {storm['start_time'].strftime('%b %d')} - Saved to {plot_path.name}")

print(f"\nHyetographs saved to: {hyetograph_folder}")
Text Only
Generating hyetograph plots for 12 storms...
======================================================================


  Storm  1: Feb 07 - Saved to storm_20200207_hyetograph.png


  Storm  2: Mar 19 - Saved to storm_20200319_hyetograph.png


  Storm  3: Mar 28 - Saved to storm_20200328_hyetograph.png


  Storm  4: Apr 17 - Saved to storm_20200417_hyetograph.png


  Storm  5: Apr 30 - Saved to storm_20200430_hyetograph.png


  Storm  6: May 28 - Saved to storm_20200528_hyetograph.png


  Storm  7: Jun 03 - Saved to storm_20200603_hyetograph.png


  Storm  8: Sep 29 - Saved to storm_20200929_hyetograph.png


  Storm  9: Oct 29 - Saved to storm_20201029_hyetograph.png


  Storm 10: Nov 11 - Saved to storm_20201111_hyetograph.png


  Storm 11: Dec 16 - Saved to storm_20201216_hyetograph.png


  Storm 12: Dec 24 - Saved to storm_20201224_hyetograph.png

Hyetographs saved to: C:\Users\bill\.config\superpowers\worktrees\ras-commander-pr-docs\codex-pr251-timeseries-alignment\examples\working\example_projects_901_aorc_catalog\BaldEagleCrkMulti2D_AORC_2020_901\Precipitation\hyetographs
Python
# Create storm catalog summary plot
fig, axes = plt.subplots(2, 2, figsize=(14, 10))

# 1. Total depth by storm
ax1 = axes[0, 0]
bars1 = ax1.bar(storm_catalog['storm_id'], storm_catalog['total_depth_in'], color='steelblue', alpha=0.8)
ax1.set_xlabel('Storm ID')
ax1.set_ylabel('Total Depth (inches)')
ax1.set_title('Precipitation Totals')
ax1.grid(True, alpha=0.3, axis='y')
# Add value labels
for bar, val in zip(bars1, storm_catalog['total_depth_in']):
    ax1.text(bar.get_x() + bar.get_width()/2, bar.get_height() + 0.02, f'{val:.2f}',
             ha='center', va='bottom', fontsize=8)

# 2. Peak intensity
ax2 = axes[0, 1]
bars2 = ax2.bar(storm_catalog['storm_id'], storm_catalog['peak_intensity_in_hr'], color='darkorange', alpha=0.8)
ax2.set_xlabel('Storm ID')
ax2.set_ylabel('Peak Intensity (in/hr)')
ax2.set_title('Peak Intensity')
ax2.grid(True, alpha=0.3, axis='y')

# 3. Duration
ax3 = axes[1, 0]
bars3 = ax3.bar(storm_catalog['storm_id'], storm_catalog['duration_hours'], color='green', alpha=0.8)
ax3.set_xlabel('Storm ID')
ax3.set_ylabel('Duration (hours)')
ax3.set_title('Storm Duration')
ax3.grid(True, alpha=0.3, axis='y')

# 4. Timeline
ax4 = axes[1, 1]
for idx, storm in storm_catalog.iterrows():
    ax4.barh(storm['storm_id'], storm['duration_hours'], left=storm['start_time'].dayofyear,
             color='steelblue', alpha=0.7, height=0.6)
ax4.set_xlabel(f'Day of Year ({YEAR})')
ax4.set_ylabel('Storm ID')
ax4.set_title('Storm Timeline')
ax4.grid(True, alpha=0.3, axis='x')

plt.suptitle(f'Storm Catalog Summary - {YEAR}\n({len(storm_catalog)} storms, '
             f'{storm_catalog["total_depth_in"].sum():.1f} inches total)', 
             fontsize=14, fontweight='bold')
plt.tight_layout()

# Save summary plot
summary_path = precip_folder / "storm_catalog_summary.png"
plt.savefig(summary_path, dpi=150, bbox_inches='tight')
print(f"Summary plot saved to: {summary_path}")
plt.show()
Text Only
Summary plot saved to: C:\Users\bill\.config\superpowers\worktrees\ras-commander-pr-docs\codex-pr251-timeseries-alignment\examples\working\example_projects_901_aorc_catalog\BaldEagleCrkMulti2D_AORC_2020_901\Precipitation\storm_catalog_summary.png

png

Python
# Create storm plans with precipitation data written to HDF
print(f"Creating HEC-RAS plans for {len(storm_catalog)} storms...")
print("="*70)

# Use create_storm_plans which calls set_gridded_precipitation
# This writes precipitation data directly to the .u##.hdf file
results = PrecipAorc.create_storm_plans(
    storm_catalog=storm_catalog,
    bounds=bounds,
    template_plan=TEMPLATE_PLAN,
    precip_folder="Precipitation",
    ras_object=ras,
    download_data=False  # Already downloaded above
)

# Show results
print(f"\nPlan Creation Results:")
print(results[['storm_id', 'start_time', 'total_depth_in', 'plan_number', 'status']].to_string(index=False))

# Refresh plan list
ras.plan_df = ras.get_plan_entries()
print(f"\nTotal plans in project: {len(ras.plan_df)}")
Text Only
2026-07-08 23:18:49 - ras_commander.precip.PrecipAorc - INFO - Creating storm plans from template plan 06 (unsteady 03); processing 12 storms


Creating HEC-RAS plans for 12 storms...
======================================================================


2026-07-08 23:18:49 - ras_commander.RasUnsteady - INFO - Configured gridded precipitation in BaldEagleDamBrk.u04: source=.\Precipitation\storm_20200207.nc, interpolation=Bilinear, dataset=APCP_surface


2026-07-08 23:18:50 - ras_commander.RasUnsteady - INFO - Imported gridded precipitation into BaldEagleDamBrk.u04.hdf: 120 timesteps, 720 cells, range=0.0-40.7 mm


2026-07-08 23:18:50 - ras_commander.RasPlan - INFO - Updated simulation date in plan file: BaldEagleDamBrk.p07


2026-07-08 23:18:51 - ras_commander.precip.PrecipAorc - INFO - Created storm 1: 2020-02-07 (1.02 in) -> plan 07, unsteady 04, with HDF time series


2026-07-08 23:18:51 - ras_commander.RasUnsteady - INFO - Configured gridded precipitation in BaldEagleDamBrk.u05: source=.\Precipitation\storm_20200319.nc, interpolation=Bilinear, dataset=APCP_surface


2026-07-08 23:18:51 - ras_commander.RasUnsteady - INFO - Imported gridded precipitation into BaldEagleDamBrk.u05.hdf: 120 timesteps, 720 cells, range=0.0-45.9 mm


2026-07-08 23:18:52 - ras_commander.RasPlan - INFO - Updated simulation date in plan file: BaldEagleDamBrk.p08


2026-07-08 23:18:52 - ras_commander.precip.PrecipAorc - INFO - Created storm 2: 2020-03-19 (0.76 in) -> plan 08, unsteady 05, with HDF time series


2026-07-08 23:18:53 - ras_commander.RasUnsteady - INFO - Configured gridded precipitation in BaldEagleDamBrk.u06: source=.\Precipitation\storm_20200328.nc, interpolation=Bilinear, dataset=APCP_surface


2026-07-08 23:18:53 - ras_commander.RasUnsteady - INFO - Imported gridded precipitation into BaldEagleDamBrk.u06.hdf: 144 timesteps, 720 cells, range=0.0-62.8 mm


2026-07-08 23:18:54 - ras_commander.RasPlan - INFO - Updated simulation date in plan file: BaldEagleDamBrk.p09


2026-07-08 23:18:54 - ras_commander.precip.PrecipAorc - INFO - Created storm 3: 2020-03-28 (1.53 in) -> plan 09, unsteady 06, with HDF time series


2026-07-08 23:18:55 - ras_commander.RasUnsteady - INFO - Configured gridded precipitation in BaldEagleDamBrk.u14: source=.\Precipitation\storm_20200417.nc, interpolation=Bilinear, dataset=APCP_surface


2026-07-08 23:18:55 - ras_commander.RasUnsteady - INFO - Imported gridded precipitation into BaldEagleDamBrk.u14.hdf: 144 timesteps, 720 cells, range=0.0-30.2 mm


2026-07-08 23:18:56 - ras_commander.RasPlan - INFO - Updated simulation date in plan file: BaldEagleDamBrk.p10


2026-07-08 23:18:56 - ras_commander.precip.PrecipAorc - INFO - Created storm 4: 2020-04-17 (0.78 in) -> plan 10, unsteady 14, with HDF time series


2026-07-08 23:18:56 - ras_commander.RasUnsteady - INFO - Configured gridded precipitation in BaldEagleDamBrk.u15: source=.\Precipitation\storm_20200430.nc, interpolation=Bilinear, dataset=APCP_surface


2026-07-08 23:18:57 - ras_commander.RasUnsteady - INFO - Imported gridded precipitation into BaldEagleDamBrk.u15.hdf: 144 timesteps, 720 cells, range=0.0-69.2 mm


2026-07-08 23:18:57 - ras_commander.RasPlan - INFO - Updated simulation date in plan file: BaldEagleDamBrk.p11


2026-07-08 23:18:57 - ras_commander.precip.PrecipAorc - INFO - Created storm 5: 2020-04-30 (2.12 in) -> plan 11, unsteady 15, with HDF time series


2026-07-08 23:18:58 - ras_commander.RasUnsteady - INFO - Configured gridded precipitation in BaldEagleDamBrk.u16: source=.\Precipitation\storm_20200528.nc, interpolation=Bilinear, dataset=APCP_surface


2026-07-08 23:18:58 - ras_commander.RasUnsteady - INFO - Imported gridded precipitation into BaldEagleDamBrk.u16.hdf: 168 timesteps, 720 cells, range=0.0-61.2 mm


2026-07-08 23:18:59 - ras_commander.RasPlan - INFO - Updated simulation date in plan file: BaldEagleDamBrk.p12


2026-07-08 23:18:59 - ras_commander.precip.PrecipAorc - INFO - Created storm 6: 2020-05-28 (1.19 in) -> plan 12, unsteady 16, with HDF time series


2026-07-08 23:19:00 - ras_commander.RasUnsteady - INFO - Configured gridded precipitation in BaldEagleDamBrk.u17: source=.\Precipitation\storm_20200603.nc, interpolation=Bilinear, dataset=APCP_surface


2026-07-08 23:19:00 - ras_commander.RasUnsteady - INFO - Imported gridded precipitation into BaldEagleDamBrk.u17.hdf: 120 timesteps, 720 cells, range=0.0-47.1 mm


2026-07-08 23:19:01 - ras_commander.RasPlan - INFO - Updated simulation date in plan file: BaldEagleDamBrk.p14


2026-07-08 23:19:01 - ras_commander.precip.PrecipAorc - INFO - Created storm 7: 2020-06-03 (0.85 in) -> plan 14, unsteady 17, with HDF time series


2026-07-08 23:19:02 - ras_commander.RasUnsteady - INFO - Configured gridded precipitation in BaldEagleDamBrk.u18: source=.\Precipitation\storm_20200929.nc, interpolation=Bilinear, dataset=APCP_surface


2026-07-08 23:19:02 - ras_commander.RasUnsteady - INFO - Imported gridded precipitation into BaldEagleDamBrk.u18.hdf: 144 timesteps, 720 cells, range=0.0-39.5 mm


2026-07-08 23:19:03 - ras_commander.RasPlan - INFO - Updated simulation date in plan file: BaldEagleDamBrk.p16


2026-07-08 23:19:03 - ras_commander.precip.PrecipAorc - INFO - Created storm 8: 2020-09-29 (1.28 in) -> plan 16, unsteady 18, with HDF time series


2026-07-08 23:19:04 - ras_commander.RasUnsteady - INFO - Configured gridded precipitation in BaldEagleDamBrk.u19: source=.\Precipitation\storm_20201029.nc, interpolation=Bilinear, dataset=APCP_surface


2026-07-08 23:19:04 - ras_commander.RasUnsteady - INFO - Imported gridded precipitation into BaldEagleDamBrk.u19.hdf: 144 timesteps, 720 cells, range=0.0-60.8 mm


2026-07-08 23:19:05 - ras_commander.RasPlan - INFO - Updated simulation date in plan file: BaldEagleDamBrk.p20


2026-07-08 23:19:05 - ras_commander.precip.PrecipAorc - INFO - Created storm 9: 2020-10-29 (1.74 in) -> plan 20, unsteady 19, with HDF time series


2026-07-08 23:19:05 - ras_commander.RasUnsteady - INFO - Configured gridded precipitation in BaldEagleDamBrk.u20: source=.\Precipitation\storm_20201111.nc, interpolation=Bilinear, dataset=APCP_surface


2026-07-08 23:19:06 - ras_commander.RasUnsteady - INFO - Imported gridded precipitation into BaldEagleDamBrk.u20.hdf: 120 timesteps, 720 cells, range=0.0-77.9 mm


2026-07-08 23:19:06 - ras_commander.RasPlan - INFO - Updated simulation date in plan file: BaldEagleDamBrk.p21


2026-07-08 23:19:06 - ras_commander.precip.PrecipAorc - INFO - Created storm 10: 2020-11-11 (1.99 in) -> plan 21, unsteady 20, with HDF time series


2026-07-08 23:19:07 - ras_commander.RasUnsteady - INFO - Configured gridded precipitation in BaldEagleDamBrk.u21: source=.\Precipitation\storm_20201216.nc, interpolation=Bilinear, dataset=APCP_surface


2026-07-08 23:19:07 - ras_commander.RasUnsteady - INFO - Imported gridded precipitation into BaldEagleDamBrk.u21.hdf: 144 timesteps, 720 cells, range=0.0-63.2 mm


2026-07-08 23:19:08 - ras_commander.RasPlan - INFO - Updated simulation date in plan file: BaldEagleDamBrk.p22


2026-07-08 23:19:08 - ras_commander.precip.PrecipAorc - INFO - Created storm 11: 2020-12-16 (1.83 in) -> plan 22, unsteady 21, with HDF time series


2026-07-08 23:19:09 - ras_commander.RasUnsteady - INFO - Configured gridded precipitation in BaldEagleDamBrk.u22: source=.\Precipitation\storm_20201224.nc, interpolation=Bilinear, dataset=APCP_surface


2026-07-08 23:19:09 - ras_commander.RasUnsteady - INFO - Imported gridded precipitation into BaldEagleDamBrk.u22.hdf: 144 timesteps, 720 cells, range=0.0-89.5 mm


2026-07-08 23:19:10 - ras_commander.RasPlan - INFO - Updated simulation date in plan file: BaldEagleDamBrk.p23


2026-07-08 23:19:10 - ras_commander.precip.PrecipAorc - INFO - Created storm 12: 2020-12-24 (2.72 in) -> plan 23, unsteady 22, with HDF time series


2026-07-08 23:19:10 - ras_commander.precip.PrecipAorc - INFO - Storm plan creation complete: 12/12 successful



Plan Creation Results:
 storm_id          start_time  total_depth_in plan_number  status
        1 2020-02-07 09:00:00           1.023          07 success
        2 2020-03-19 01:00:00           0.763          08 success
        3 2020-03-28 08:00:00           1.530          09 success
        4 2020-04-17 17:00:00           0.777          10 success
        5 2020-04-30 04:00:00           2.115          11 success
        6 2020-05-28 12:00:00           1.188          12 success
        7 2020-06-03 06:00:00           0.846          14 success
        8 2020-09-29 11:00:00           1.278          16 success
        9 2020-10-29 07:00:00           1.735          20 success
       10 2020-11-11 12:00:00           1.989          21 success
       11 2020-12-16 17:00:00           1.834          22 success
       12 2020-12-24 12:00:00           2.720          23 success

Total plans in project: 23

Rectify Sayers Dam Gate Schedule

The template plan includes Sayers Dam, an SA/2D connection with a low-level Gate #1 and a spillway/weir profile. The source gridded-precipitation example uses a constant 2.0 ft opening for this gate at its native gate-control interval.

Because observed Sayers Dam gate-operation data are not included in this example project, this notebook preserves the source model assumption (2.0 ft constant opening) and confirms that each cloned storm simulation has enough gate-opening values to cover its simulation window. This avoids changing the example hydraulics while making the time-series boundary condition consistent with the cloned plan windows.

Python
# Confirm Sayers Dam / Gate #1 operation schedule covers every storm plan.
# The example project uses a constant 2.0 ft opening; preserve that source-model assumption.
GATE_BOUNDARY_INDEX = 0  # only gate-bearing boundary in these cloned unsteady files
GATE_NAME = "Gate #1"
SAYERS_DAM_GATE_OPENING_FT = 2.0


def _parse_ras_simulation_date(sim_date_text):
    """Parse HEC-RAS Simulation Date text into pandas timestamps."""
    start_date, start_time, end_date, end_time = [part.strip() for part in str(sim_date_text).split(",")]
    start = pd.to_datetime(f"{start_date} {start_time}", format="%d%b%Y %H%M")
    end = pd.to_datetime(f"{end_date} {end_time}", format="%d%b%Y %H%M")
    return start, end


def _ras_interval_hours(interval_text):
    """Convert common HEC-RAS interval text to hours."""
    text = str(interval_text).strip().upper()
    if text.endswith("HOUR"):
        return float(text.replace("HOUR", "") or 1)
    if text.endswith("MIN"):
        return float(text.replace("MIN", "")) / 60.0
    raise ValueError(f"Unsupported gate interval: {interval_text!r}")


rectified_gate_rows = []
for _, storm in results[results["status"].eq("success")].iterrows():
    plan_num = str(storm["plan_number"]).zfill(2)
    plan_match = ras.plan_df[ras.plan_df["plan_number"].astype(str).str.zfill(2).eq(plan_num)]
    if plan_match.empty:
        raise ValueError(f"Plan {plan_num} not found in ras.plan_df")
    plan_row = plan_match.iloc[0]
    unsteady_num = str(storm.get("unsteady_number") or plan_row.get("unsteady_number")).zfill(2)

    if "sim_start" in storm.index and "sim_end" in storm.index:
        sim_start = pd.to_datetime(storm["sim_start"])
        sim_end = pd.to_datetime(storm["sim_end"])
    else:
        sim_start, sim_end = _parse_ras_simulation_date(plan_row["Simulation Date"])

    gate_before = RasUnsteady.get_gate_openings(
        unsteady_num,
        boundary_index=GATE_BOUNDARY_INDEX,
        ras_object=ras,
    )
    interval = gate_before["interval"]
    interval_hours = _ras_interval_hours(interval)
    duration_hours = (sim_end - sim_start) / pd.Timedelta(hours=1)
    required_values = int(np.ceil(duration_hours / interval_hours)) + 1
    existing_values = list(gate_before["values"])
    values_written = gate_before["count"]
    action = "already covered"

    if len(existing_values) < required_values:
        gate_values = [SAYERS_DAM_GATE_OPENING_FT] * required_values
        RasUnsteady.set_gate_openings(
            unsteady_num,
            values=gate_values,
            boundary_index=GATE_BOUNDARY_INDEX,
            gate_name=GATE_NAME,
            interval=interval,
            ras_object=ras,
        )
        gate_after = RasUnsteady.get_gate_openings(
            unsteady_num,
            boundary_index=GATE_BOUNDARY_INDEX,
            ras_object=ras,
        )
        values_written = gate_after["count"]
        action = "extended"

    rectified_gate_rows.append({
        "storm_id": storm["storm_id"],
        "plan_number": plan_num,
        "unsteady_number": unsteady_num,
        "sim_start": sim_start,
        "sim_end": sim_end,
        "duration_hours": duration_hours,
        "gate_name": gate_before["gate_name"],
        "gate_opening_ft": SAYERS_DAM_GATE_OPENING_FT,
        "required_values": required_values,
        "values_written": values_written,
        "interval": interval,
        "action": action,
    })

# Refresh boundary metadata after checking cloned unsteady files.
ras.boundaries_df = ras.get_boundary_conditions()
gate_rectification_df = pd.DataFrame(rectified_gate_rows)
print("Sayers Dam gate schedules checked using the source model's constant 2.0 ft opening:")
display(gate_rectification_df[[
    "storm_id", "plan_number", "unsteady_number", "duration_hours",
    "required_values", "values_written", "gate_opening_ft", "interval", "action"
]])
Text Only
Sayers Dam gate schedules checked using the source model's constant 2.0 ft opening:
storm_id plan_number unsteady_number duration_hours required_values values_written gate_opening_ft interval action
0 1 07 04 104.0 19 100 2.0 6HOUR already covered
1 2 08 05 105.0 19 100 2.0 6HOUR already covered
2 3 09 06 123.0 22 100 2.0 6HOUR already covered
3 4 10 14 111.0 20 100 2.0 6HOUR already covered
4 5 11 15 116.0 21 100 2.0 6HOUR already covered
5 6 12 16 133.0 24 100 2.0 6HOUR already covered
6 7 14 17 112.0 20 100 2.0 6HOUR already covered
7 8 16 18 117.0 21 100 2.0 6HOUR already covered
8 9 20 19 131.0 23 100 2.0 6HOUR already covered
9 10 21 20 106.0 19 100 2.0 6HOUR already covered
10 11 22 21 112.0 20 100 2.0 6HOUR already covered
11 12 23 22 120.0 21 100 2.0 6HOUR already covered

Step 5: Execute Storm Plans in Parallel

Run all storm plans using parallel execution with 2 cores per instance and 3 workers.

Python
# Get list of storm plan numbers
storm_plan_numbers = results[results['status'] == 'success']['plan_number'].tolist()
print(f"Plans to execute: {storm_plan_numbers}")
print(f"Execution config: {NUM_CORES} cores x {MAX_WORKERS} workers")
print("="*70)
Text Only
Plans to execute: ['07', '08', '09', '10', '11', '12', '14', '16', '20', '21', '22', '23']
Execution config: 2 cores x 3 workers
======================================================================
Python
# Execute plans in parallel
import time

print(f"Starting parallel execution of {len(storm_plan_numbers)} plans...")
start_time = time.time()

execution_results = RasCmdr.compute_parallel(
    plan_number=storm_plan_numbers,
    max_workers=MAX_WORKERS,
    num_cores=NUM_CORES,
    ras_object=ras,
    overwrite_dest=True
)

elapsed = time.time() - start_time

# Results summary
success_count = sum(1 for success in execution_results.values() if success)
fail_count = len(execution_results) - success_count

print(f"\nExecution complete in {elapsed/60:.1f} minutes")
print(f"  Successful: {success_count}")
print(f"  Failed: {fail_count}")
print(f"\nResults copied to: {ras.project_folder.parent / (ras.project_folder.name + ' [Computed]')}")
Text Only
2026-07-08 23:19:10 - ras_commander.RasCmdr - INFO - Filtered plans to execute: 12 plan(s) (07, 08, 09, 10, 11 ... 21, 22, 23)


2026-07-08 23:19:10 - ras_commander.RasCmdr - INFO - Adjusted max_workers to 3 based on the number of plans to compute: 12


Starting parallel execution of 12 plans...


2026-07-08 23:19:13 - ras_commander.RasCmdr - INFO - Prepared 3 worker folder(s) for parallel execution


2026-07-09 00:31:29 - ras_commander.RasCmdr - INFO - Consolidating worker artifacts back to original project folder: BaldEagleCrkMulti2D_AORC_2020_901


2026-07-09 00:40:06 - ras_commander.RasCmdr - INFO - Consolidated 39 worker artifact(s) to BaldEagleCrkMulti2D_AORC_2020_901


2026-07-09 00:40:06 - ras_commander.RasCmdr - INFO - Execution results: 12/12 plan(s) successful



Execution complete in 80.9 minutes
  Successful: 12
  Failed: 0

Results copied to: C:\Users\bill\.config\superpowers\worktrees\ras-commander-pr-docs\codex-pr251-timeseries-alignment\examples\working\example_projects_901_aorc_catalog\BaldEagleCrkMulti2D_AORC_2020_901 [Computed]

Step 6: Extract and Compare Results

Re-initialize from the computed folder and extract results summary.

Python
# Results are already in original folder (v0.88.1+) - no re-initialization needed

# Extract results summary
import h5py
import pathlib
from ras_commander.results.ResultsParser import ResultsParser

print("\nStorm Execution Results")
print("="*80)

parallel_results_df = getattr(execution_results, 'results_df', pd.DataFrame())
if parallel_results_df is None or parallel_results_df.empty:
    parallel_results_df = getattr(ras, 'results_df', pd.DataFrame())


def _existing_plan_hdf_path(plan_num, plan_row, results_df):
    """Return the first existing plan HDF path from ras-commander metadata."""
    plan_num_text = str(plan_num).zfill(2)
    candidates = []

    if results_df is not None and not results_df.empty and 'plan_number' in results_df.columns:
        result_rows = results_df.copy()
        result_rows['_plan_number_text'] = result_rows['plan_number'].astype(str).str.zfill(2)
        result_rows = result_rows[result_rows['_plan_number_text'] == plan_num_text]
        if not result_rows.empty:
            result_row = result_rows.iloc[0]
            for column in ('hdf_path', 'HDF_Results_Path'):
                if column in result_row.index:
                    value = result_row.get(column)
                    if pd.notna(value) and str(value).strip():
                        candidates.append(pathlib.Path(value))

    if len(plan_row) > 0:
        row = plan_row.iloc[0]
        if 'HDF_Results_Path' in row.index:
            value = row.get('HDF_Results_Path')
            if pd.notna(value) and str(value).strip():
                candidates.append(pathlib.Path(value))
        plan_path = pathlib.Path(row['full_path'])
        candidates.append(pathlib.Path(f"{plan_path}.hdf"))

    candidates.append(ras.project_folder / f"{ras.project_name}.p{plan_num_text}.hdf")
    legacy_computed_folder = ras.project_folder.parent / f"{ras.project_folder.name} [Computed]"
    candidates.append(legacy_computed_folder / f"{ras.project_name}.p{plan_num_text}.hdf")

    seen = set()
    for candidate in candidates:
        candidate = pathlib.Path(candidate)
        key = str(candidate).lower()
        if key in seen:
            continue
        seen.add(key)
        if candidate.exists():
            return candidate

    return None


def _read_compute_status(hdf_path):
    """Read concise compute status from a plan HDF."""
    status = {"compute_time": None, "has_errors": False, "has_warnings": False, "message_tail": ""}
    if hdf_path is None:
        return status
    try:
        with h5py.File(hdf_path, 'r') as f:
            if 'Results/Summary/Compute Processes' in f:
                cp = f['Results/Summary/Compute Processes'][:]
                if len(cp) > 0:
                    status['compute_time'] = cp[0]['Compute Time'].decode('utf-8').strip()
            if 'Results/Summary/Compute Messages (text)' in f:
                raw = f['Results/Summary/Compute Messages (text)'][0]
                text = raw.decode('utf-8', errors='replace') if isinstance(raw, bytes) else str(raw)
                parsed = ResultsParser.parse_compute_messages(text)
                status['has_errors'] = parsed['has_errors']
                status['has_warnings'] = parsed['has_warnings']
                lines = [line.strip() for line in text.replace('\r', '\n').split('\n') if line.strip()]
                status['message_tail'] = ' | '.join(lines[-3:])
    except Exception as exc:
        status['has_errors'] = True
        status['message_tail'] = f"Unable to read compute status: {exc}"
    return status


storm_results = []
for idx, row in results[results['status'] == 'success'].iterrows():
    storm_id = row['storm_id']
    plan_num = str(row['plan_number']).zfill(2)
    exec_success = execution_results.get(plan_num, execution_results.get(row['plan_number'], False))

    plan_row = ras.plan_df[ras.plan_df['plan_number'].astype(str).str.zfill(2) == plan_num]
    if len(plan_row) == 0:
        continue

    hdf_path = _existing_plan_hdf_path(plan_num, plan_row, parallel_results_df)
    hdf_exists = hdf_path is not None
    compute_status = _read_compute_status(hdf_path)
    hydraulics_complete = bool(exec_success and hdf_exists and not compute_status['has_errors'])

    result = {
        'storm_id': storm_id,
        'plan_number': plan_num,
        'start_time': row['start_time'],
        'total_depth_in': row['total_depth_in'],
        'exec_success': exec_success,
        'hydraulics_complete': hydraulics_complete,
        'hdf_exists': hdf_exists,
        'has_errors': compute_status['has_errors'],
        'has_warnings': compute_status['has_warnings'],
        'compute_time': compute_status['compute_time'],
        'message_tail': compute_status['message_tail'],
        'hdf_path': str(hdf_path) if hdf_path else None,
        'hdf_size_mb': hdf_path.stat().st_size / 1e6 if hdf_exists else 0,
    }

    storm_results.append(result)
    if hydraulics_complete:
        status = 'HYDRAULICS COMPLETE'
    elif hdf_exists and compute_status['has_errors']:
        status = 'COMPUTE ERROR'
    elif exec_success:
        status = 'COMPUTE OK - HDF NOT FOUND'
    else:
        status = 'FAILED'
    print(f"Storm {storm_id:2d} ({row['start_time'].strftime('%b %d')}): Plan {plan_num} - {status} - {result['hdf_size_mb']:.1f} MB")
    if compute_status['has_errors'] and compute_status['message_tail']:
        print(f"  Last compute messages: {compute_status['message_tail']}")

storm_results_df = pd.DataFrame(storm_results)
print(f"\nHydraulic results complete: {int(storm_results_df['hydraulics_complete'].sum())} of {len(storm_results_df)} storms")
print(f"Plan HDF files found: {int(storm_results_df['hdf_exists'].sum())} of {len(storm_results_df)} storms")
Text Only
Storm Execution Results
================================================================================
Storm  1 (Feb 07): Plan 07 - HYDRAULICS COMPLETE - 11670.2 MB
Storm  2 (Mar 19): Plan 08 - HYDRAULICS COMPLETE - 11869.9 MB
Storm  3 (Mar 28): Plan 09 - HYDRAULICS COMPLETE - 14850.5 MB
Storm  4 (Apr 17): Plan 10 - HYDRAULICS COMPLETE - 11044.3 MB
Storm  5 (Apr 30): Plan 11 - HYDRAULICS COMPLETE - 12526.2 MB
Storm  6 (May 28): Plan 12 - HYDRAULICS COMPLETE - 13604.4 MB
Storm  7 (Jun 03): Plan 14 - HYDRAULICS COMPLETE - 11301.1 MB
Storm  8 (Sep 29): Plan 16 - HYDRAULICS COMPLETE - 11011.3 MB
Storm  9 (Oct 29): Plan 20 - HYDRAULICS COMPLETE - 14755.4 MB
Storm 10 (Nov 11): Plan 21 - HYDRAULICS COMPLETE - 10602.8 MB
Storm 11 (Dec 16): Plan 22 - HYDRAULICS COMPLETE - 11169.7 MB
Storm 12 (Dec 24): Plan 23 - HYDRAULICS COMPLETE - 14002.2 MB

Hydraulic results complete: 12 of 12 storms
Plan HDF files found: 12 of 12 storms
Python
# Display concise execution status for all generated storm plans.
_display_cols = [
    'plan_number', 'start_time', 'total_depth_in', 'hydraulics_complete',
    'has_errors', 'has_warnings', 'compute_time', 'hdf_size_mb'
]
storm_results_df[_display_cols]
plan_number start_time total_depth_in hydraulics_complete has_errors has_warnings compute_time hdf_size_mb
0 07 2020-02-07 09:00:00 1.023 True False False 11670.151944
1 08 2020-03-19 01:00:00 0.763 True False False 11869.876277
2 09 2020-03-28 08:00:00 1.530 True False False 14850.520406
3 10 2020-04-17 17:00:00 0.777 True False False 11044.261680
4 11 2020-04-30 04:00:00 2.115 True False False 12526.177605
5 12 2020-05-28 12:00:00 1.188 True False False 13604.426073
6 14 2020-06-03 06:00:00 0.846 True False False 11301.050709
7 16 2020-09-29 11:00:00 1.278 True False False 11011.269451
8 20 2020-10-29 07:00:00 1.735 True False False 14755.350674
9 21 2020-11-11 12:00:00 1.989 True False False 10602.846230
10 22 2020-12-16 17:00:00 1.834 True False False 11169.743784
11 23 2020-12-24 12:00:00 2.720 True False False 14002.194917
Python
# Summary visualization
if storm_results:
    fig, axes = plt.subplots(1, 2, figsize=(14, 5))

    # Precipitation totals
    ax1 = axes[0]
    storm_dates = [r['start_time'].strftime('%m/%d') for r in storm_results]
    precip_totals = [r['total_depth_in'] for r in storm_results]
    colors = ['green' if r['exec_success'] else 'red' for r in storm_results]
    from matplotlib.patches import Patch
    status_legend = [
        Patch(facecolor='green', alpha=0.8, label='Successful'),
        Patch(facecolor='red', alpha=0.8, label='Failed'),
    ]

    bars = ax1.bar(range(len(storm_results)), precip_totals, color=colors, alpha=0.8)
    ax1.set_xlabel('Storm')
    ax1.set_ylabel('Total Precipitation (inches)')
    ax1.set_title('Storm Precipitation Totals')
    ax1.set_xticks(range(len(storm_results)))
    ax1.set_xticklabels(storm_dates, rotation=45)
    ax1.legend(handles=status_legend, title='Execution status', fontsize=9)
    ax1.grid(True, alpha=0.3, axis='y')

    # HDF sizes
    ax2 = axes[1]
    hdf_sizes = [r.get('hdf_size_mb', 0) for r in storm_results]
    hdf_count = sum(1 for r in storm_results if r.get('hdf_exists'))
    ax2.set_xlabel('Storm')
    ax2.set_ylabel('HDF File Size (MB)')
    ax2.set_title(f'Simulation Output Files ({hdf_count}/{len(storm_results)} HDFs found)')
    ax2.set_xticks(range(len(storm_results)))
    ax2.set_xticklabels(storm_dates, rotation=45)
    ax2.grid(True, alpha=0.3, axis='y')

    if any(size > 0 for size in hdf_sizes):
        hdf_bars = ax2.bar(range(len(storm_results)), hdf_sizes, color=colors, alpha=0.8)
        y_offset = max(hdf_sizes) * 0.02
        for bar, result in zip(hdf_bars, storm_results):
            if result.get('hdf_exists'):
                ax2.text(bar.get_x() + bar.get_width()/2, bar.get_height() + y_offset,
                         f"{result['hdf_size_mb']:.1f}", ha='center', va='bottom', fontsize=8)
            else:
                ax2.text(bar.get_x() + bar.get_width()/2, y_offset, 'No HDF',
                         ha='center', va='bottom', fontsize=8, rotation=90, color='dimgray')
    else:
        ax2.bar(range(len(storm_results)), [0] * len(storm_results), color='lightgray', alpha=0.25)
        ax2.set_ylim(0, 1)
        ax2.set_yticks([])
        ax2.set_ylabel('HDF File Size (MB)\n(no HDF results found)')
        ax2.grid(False)
        ax2.spines['left'].set_visible(False)
        ax2.text(0.5, 0.55,
                 'No plan HDF result files were found\n'
                 'Storms were cataloged and plans were created, but file sizes are unavailable.\n'
                 'Run or verify HEC-RAS execution to populate this panel.',
                 transform=ax2.transAxes, ha='center', va='center', fontsize=10,
                 bbox=dict(boxstyle='round,pad=0.5', facecolor='white', edgecolor='gray', alpha=0.9))

    plt.suptitle(f'AORC Storm Simulation Results - {YEAR}', fontsize=14, fontweight='bold')
    plt.tight_layout()
    plt.show()

png

Data Export Summary

All precipitation data has been exported to the project's Precipitation/ folder:

Python
# List all exported files
print(f"Precipitation Data Export Summary")
print(f"Location: {precip_folder}")
print("="*70)

# Count files
nc_files = list(precip_folder.glob("storm_*.nc"))
png_files = list(hyetograph_folder.glob("*.png"))
csv_files = list(precip_folder.glob("*.csv"))
summary_files = list(precip_folder.glob("storm_catalog_summary.png"))

print(f"\nExported Files:")
print(f"  Storm Catalog CSV:    {len(csv_files)} file(s)")
print(f"  NetCDF Precip Files:  {len(nc_files)} file(s)")
print(f"  Hyetograph PNGs:      {len(png_files)} file(s)")
print(f"  Summary Plot:         {len(summary_files)} file(s)")

# Calculate total size
total_size = sum(f.stat().st_size for f in nc_files + png_files + csv_files + summary_files)
print(f"\nTotal Size: {total_size / 1e6:.1f} MB")

print(f"\nFolder Structure:")
print(f"  Precipitation/")
print(f"    storm_catalog.csv")
print(f"    storm_catalog_summary.png")
print(f"    storm_YYYYMMDD.nc  (x{len(nc_files)})")
print(f"    hyetographs/")
print(f"      storm_YYYYMMDD_hyetograph.png  (x{len(png_files)})")
Text Only
Precipitation Data Export Summary
Location: C:\Users\bill\.config\superpowers\worktrees\ras-commander-pr-docs\codex-pr251-timeseries-alignment\examples\working\example_projects_901_aorc_catalog\BaldEagleCrkMulti2D_AORC_2020_901\Precipitation
======================================================================

Exported Files:
  Storm Catalog CSV:    1 file(s)
  NetCDF Precip Files:  12 file(s)
  Hyetograph PNGs:      12 file(s)
  Summary Plot:         1 file(s)

Total Size: 3.6 MB

Folder Structure:
  Precipitation/
    storm_catalog.csv
    storm_catalog_summary.png
    storm_YYYYMMDD.nc  (x12)
    hyetographs/
      storm_YYYYMMDD_hyetograph.png  (x12)

Summary

This notebook demonstrated a complete AORC precipitation workflow:

  1. Project Setup - Extracted example and created labeled working copy (with cleanup)
  2. Bounds Calculation - Used HdfProject.get_project_bounds_latlon() with 50% buffer
  3. Storm Catalog - Generated catalog using USGS standard parameters
  4. Precipitation Export - Downloaded AORC data, exported CSV catalog and hyetograph plots
  5. Plan Creation - Created HEC-RAS plans with precipitation written to HDF
  6. Parallel Execution - Ran all plans with 2 cores and 3 workers and checked compute-message status

Key Functions

Function Description
HdfProject.get_project_bounds_latlon() Get buffered bounds from geometry HDF
PrecipAorc.get_storm_catalog() Generate catalog of precipitation events
PrecipAorc.download() Download AORC data to NetCDF
PrecipAorc.create_storm_plans() Create plans with precipitation in HDF
RasCmdr.compute_parallel() Execute plans in parallel

Optional: Process All Available AORC Years (1979-2024)

The following cell processes storm catalogs for all available years in the AORC dataset. Each year gets its own project folder with full precipitation data export.

Warning: This can take a very long time and generate many files!

Python
# OPTIONAL: Process all AORC years
# Set PROCESS_ALL_YEARS = True to run

PROCESS_ALL_YEARS = False  # Set to True to run

if PROCESS_ALL_YEARS:
    import time

    # AORC is available from 1979-present for CONUS
    ALL_YEARS = list(range(1979, 2025))  # 1979-2024

    print(f"Processing {len(ALL_YEARS)} years of AORC data")
    print(f"Configuration: {NUM_CORES} cores x {MAX_WORKERS} workers per year")
    print("="*80)

    # Clean up this notebook's existing AORC folders first
    print("\nCleaning up existing AORC folders...")
    existing_folders = list(base_project.parent.glob(f"{PROJECT_NAME}_AORC_*_{SUFFIX}*"))
    removed_folders = 0
    for folder in existing_folders:
        if folder.is_dir():
            print(f"  Removing: {folder.name}")
            if RasUtils.remove_with_retry(folder, ras_object=None):
                removed_folders += 1
            else:
                raise PermissionError(f"Unable to remove stale notebook 901 folder: {folder}")
    print(f"  Removed {removed_folders} existing folders")

    all_year_results = {}

    for year in ALL_YEARS:
        try:
            print(f"\n{'='*80}")
            print(f"PROCESSING YEAR: {year}")
            print(f"{'='*80}")

            # Create year-specific working folder
            year_folder = base_project.parent / f"{PROJECT_NAME}_AORC_{year}_{SUFFIX}"
            print(f"  Creating: {year_folder.name}")
            shutil.copytree(base_project, year_folder)

            # Initialize
            ras_year = init_ras_project(year_folder, RAS_VERSION)

            # Generate storm catalog
            print(f"  Generating storm catalog for {year}...")
            catalog = PrecipAorc.get_storm_catalog(
                bounds=bounds,
                year=year,
                inter_event_hours=8.0,
                min_depth_inches=0.75,
                buffer_hours=48
            )
            print(f"  Found {len(catalog)} storms")

            if len(catalog) == 0:
                print(f"  No storms found for {year}, skipping")
                all_year_results[year] = {'storms': 0, 'plans': 0, 'success': 0, 'elapsed_min': 0}
                continue

            # Create precipitation folder structure
            year_precip_folder = year_folder / "Precipitation"
            year_precip_folder.mkdir(exist_ok=True)
            year_hyetograph_folder = year_precip_folder / "hyetographs"
            year_hyetograph_folder.mkdir(exist_ok=True)

            # Export storm catalog to CSV
            catalog.to_csv(year_precip_folder / "storm_catalog.csv", index=False)
            print(f"  Exported storm catalog CSV")

            # Download precipitation and create plans
            print(f"  Downloading precipitation and creating plans...")
            year_plan_results = PrecipAorc.create_storm_plans(
                storm_catalog=catalog,
                bounds=bounds,
                template_plan=TEMPLATE_PLAN,
                precip_folder="Precipitation",
                ras_object=ras_year,
                download_data=True
            )

            # Generate hyetographs for each storm
            print(f"  Generating hyetograph plots...")
            for idx, storm in catalog.iterrows():
                storm_id = storm['storm_id']
                date_str = storm['start_time'].strftime('%Y%m%d')
                precip_file = year_precip_folder / f"storm_{date_str}.nc"

                if precip_file.exists():
                    ds = xr.open_dataset(precip_file)
                    da = ds['APCP_surface']
                    hourly_mean = da.mean(dim=['x', 'y']).values
                    times = pd.to_datetime(da.time.values)

                    fig, ax = plt.subplots(figsize=(12, 5))
                    ax.bar(times, hourly_mean, width=0.03, color='steelblue', alpha=0.8)
                    ax.set_title(f"Storm {storm_id}: {storm['start_time'].strftime('%B %d, %Y')}\n"
                                 f"Total: {storm['total_depth_in']:.2f} in | Peak: {storm['peak_intensity_in_hr']:.3f} in/hr",
                                 fontsize=12, fontweight='bold')
                    ax.set_xlabel('Date/Time')
                    ax.set_ylabel('Precipitation Rate (mm/hr)')
                    ax.tick_params(axis='x', rotation=45)
                    ax.grid(True, alpha=0.3)
                    plt.tight_layout()
                    plt.savefig(year_hyetograph_folder / f"storm_{date_str}_hyetograph.png", dpi=150, bbox_inches='tight')
                    plt.close()
                    ds.close()

            # Create summary plot
            fig, axes = plt.subplots(1, 3, figsize=(15, 4))
            axes[0].bar(catalog['storm_id'], catalog['total_depth_in'], color='steelblue', alpha=0.8)
            axes[0].set_xlabel('Storm ID')
            axes[0].set_ylabel('Total Depth (in)')
            axes[0].set_title('Precipitation Totals')
            axes[1].bar(catalog['storm_id'], catalog['peak_intensity_in_hr'], color='darkorange', alpha=0.8)
            axes[1].set_xlabel('Storm ID')
            axes[1].set_ylabel('Peak Intensity (in/hr)')
            axes[1].set_title('Peak Intensity')
            axes[2].bar(catalog['storm_id'], catalog['duration_hours'], color='green', alpha=0.8)
            axes[2].set_xlabel('Storm ID')
            axes[2].set_ylabel('Duration (hours)')
            axes[2].set_title('Storm Duration')
            plt.suptitle(f'Storm Catalog Summary - {year}', fontsize=14, fontweight='bold')
            plt.tight_layout()
            plt.savefig(year_precip_folder / "storm_catalog_summary.png", dpi=150, bbox_inches='tight')
            plt.close()

            # Get plan numbers and execute
            plan_numbers = year_plan_results[year_plan_results['status'] == 'success']['plan_number'].tolist()
            print(f"  Created {len(plan_numbers)} plans")

            if len(plan_numbers) == 0:
                print(f"  No successful plans for {year}, skipping execution")
                all_year_results[year] = {'storms': len(catalog), 'plans': 0, 'success': 0, 'elapsed_min': 0}
                continue

            # Execute in parallel
            print(f"  Executing {len(plan_numbers)} plans ({NUM_CORES} cores x {MAX_WORKERS} workers)...")
            start_time = time.time()

            exec_results = RasCmdr.compute_parallel(
                plan_number=plan_numbers,
                max_workers=MAX_WORKERS,
                num_cores=NUM_CORES,
                ras_object=ras_year,
                overwrite_dest=True
            )

            elapsed = time.time() - start_time
            success = sum(1 for v in exec_results.values() if v)

            all_year_results[year] = {
                'storms': len(catalog),
                'plans': len(plan_numbers),
                'success': success,
                'elapsed_min': elapsed / 60
            }

            print(f"  Completed: {success}/{len(plan_numbers)} in {elapsed/60:.1f} minutes")

        except Exception as e:
            print(f"  ERROR processing {year}: {e}")
            all_year_results[year] = {'error': str(e)}

    # Final summary
    print("\n" + "="*80)
    print("ALL YEARS SUMMARY")
    print("="*80)

    total_storms = 0
    total_success = 0
    for year, res in sorted(all_year_results.items()):
        if 'error' in res:
            print(f"{year}: ERROR - {res['error']}")
        else:
            total_storms += res['storms']
            total_success += res['success']
            print(f"{year}: {res['storms']} storms, {res['success']}/{res['plans']} success, {res['elapsed_min']:.1f} min")

    print(f"\nTOTAL: {total_storms} storms, {total_success} successful simulations")
else:
    print("Set PROCESS_ALL_YEARS = True to process all AORC years (1979-2024)")
    print("Warning: This will take many hours and create ~50 project folders!")
    print("\nEach year will include:")
    print("  - Precipitation/storm_catalog.csv")
    print("  - Precipitation/storm_catalog_summary.png")
    print("  - Precipitation/storm_YYYYMMDD.nc (per storm)")
    print("  - Precipitation/hyetographs/storm_YYYYMMDD_hyetograph.png (per storm)")
Text Only
Set PROCESS_ALL_YEARS = True to process all AORC years (1979-2024)
Warning: This will take many hours and create ~50 project folders!

Each year will include:
  - Precipitation/storm_catalog.csv
  - Precipitation/storm_catalog_summary.png
  - Precipitation/storm_YYYYMMDD.nc (per storm)
  - Precipitation/hyetographs/storm_YYYYMMDD_hyetograph.png (per storm)