API Reference

This page provides detailed documentation for the NWM Coastal Python API.

Configuration Classes

CoastalCalibConfig

CoastalCalibConfig dataclass

CoastalCalibConfig(
    simulation,
    boundary,
    paths,
    model_config,
    monitoring=MonitoringConfig(),
    download=DownloadConfig(),
)

Complete coastal calibration workflow configuration.

Supports both SCHISM and SFINCS models via the polymorphic :attr:model_config field. The concrete type is selected by the model key in the YAML file and resolved through :data:MODEL_REGISTRY.

model property

model

Model identifier string (convenience accessor).

from_yaml classmethod

from_yaml(config_path)

Load configuration from YAML file with optional inheritance.

Supports variable interpolation using ${section.key} syntax. Variables are resolved from other config values, e.g.:

• ${user} -> value of $USER environment variable
• ${simulation.coastal_domain} -> value of simulation.coastal_domain
• ${model} -> the model type string ("schism" or "sfincs")

PARAMETER	DESCRIPTION
config_path	Path to YAML configuration file. TYPE: Path or str

RETURNS	DESCRIPTION
CoastalCalibConfig	Loaded configuration.

RAISES	DESCRIPTION
FileNotFoundError	If the configuration file does not exist.
YAMLError	If the YAML file is malformed.

Source code in src/coastal_calibration/config/schema.py

@classmethod
def from_yaml(cls, config_path: Path | str) -> CoastalCalibConfig:
    """Load configuration from YAML file with optional inheritance.

    Supports variable interpolation using ${section.key} syntax.
    Variables are resolved from other config values, e.g.:

    - ``${user}`` -> value of ``$USER`` environment variable
    - ``${simulation.coastal_domain}`` -> value of ``simulation.coastal_domain``
    - ``${model}`` -> the model type string (``"schism"`` or ``"sfincs"``)

    Parameters
    ----------
    config_path : Path or str
        Path to YAML configuration file.

    Returns
    -------
    CoastalCalibConfig
        Loaded configuration.

    Raises
    ------
    FileNotFoundError
        If the configuration file does not exist.
    yaml.YAMLError
        If the YAML file is malformed.
    """
    config_path = Path(config_path)
    if not config_path.exists():
        raise FileNotFoundError(f"Configuration file not found: {config_path}")

    try:
        data = yaml.safe_load(config_path.read_text())
    except yaml.YAMLError as e:
        raise yaml.YAMLError(f"Invalid YAML in {config_path}: {e}") from e

    if data is None:
        raise ValueError(f"Configuration file is empty: {config_path}")

    if "_base" in data:
        base_path = Path(data.pop("_base"))
        if not base_path.is_absolute():
            base_path = config_path.parent / base_path
        base_config = cls.from_yaml(base_path)
        data = _deep_merge(base_config.to_dict(), data)

    # Ensure model key has a default before interpolation
    data.setdefault("model", "schism")

    # Interpolate variables after merging
    data = _interpolate_config(data)

    return cls.from_dict(data)

from_dict classmethod

from_dict(data)

Create config from a plain dictionary.

PARAMETER	DESCRIPTION
data	Configuration dictionary with the same structure as the YAML file (see :meth:to_dict for the expected keys). The dict is read but not mutated. TYPE: dict

RETURNS	DESCRIPTION
CoastalCalibConfig

Source code in src/coastal_calibration/config/schema.py

@classmethod
def from_dict(cls, data: dict[str, Any]) -> CoastalCalibConfig:
    """Create config from a plain dictionary.

    Parameters
    ----------
    data : dict
        Configuration dictionary with the same structure as the YAML
        file (see :meth:`to_dict` for the expected keys). The dict
        is read but not mutated.

    Returns
    -------
    CoastalCalibConfig
    """
    if "model" not in data:
        raise ValueError("'model' is required (e.g., model: schism or model: sfincs)")
    model_type: str = data["model"]

    # Read but do not mutate the caller's dict.
    model_config_data = data.get("model_config") or {}

    sim_data = data.get("simulation", {})
    if "start_date" in sim_data:
        sim_data["start_date"] = pd.to_datetime(sim_data["start_date"]).to_pydatetime()
    simulation = SimulationConfig(**sim_data)

    boundary_data = data.get("boundary", {})
    if boundary_data.get("stofs_file"):
        boundary_data["stofs_file"] = Path(boundary_data["stofs_file"])
    boundary = BoundaryConfig(**boundary_data)

    paths_data = data.get("paths", {})
    paths = PathConfig(**paths_data)

    monitoring_data = data.get("monitoring", {})
    if monitoring_data.get("log_file"):
        monitoring_data["log_file"] = Path(monitoring_data["log_file"])
    monitoring = MonitoringConfig(**monitoring_data)

    download_data = data.get("download", {})
    download = DownloadConfig(**download_data)

    if model_type not in MODEL_REGISTRY:
        msg = (
            f"Unknown model type: {model_type!r}. Supported models: {', '.join(MODEL_REGISTRY)}"
        )
        raise ValueError(msg)

    model_cls = MODEL_REGISTRY[model_type]
    model_config = model_cls(**model_config_data)

    return cls(
        simulation=simulation,
        boundary=boundary,
        paths=paths,
        model_config=model_config,  # pyright: ignore[reportArgumentType]
        monitoring=monitoring,
        download=download,
    )

to_yaml

to_yaml(path)

Write configuration to YAML file.

PARAMETER	DESCRIPTION
path	Path to YAML output file. Parent directories will be created if they don't exist. TYPE: Path or str

Source code in src/coastal_calibration/config/schema.py

def to_yaml(self, path: Path | str) -> None:
    """Write configuration to YAML file.

    Parameters
    ----------
    path : Path or str
        Path to YAML output file. Parent directories will be created
        if they don't exist.
    """
    path = Path(path)
    path.parent.mkdir(parents=True, exist_ok=True)
    path.write_text(yaml.dump(self.to_dict(), default_flow_style=False, sort_keys=False))

to_dict

to_dict()

Convert config to dictionary.

Source code in src/coastal_calibration/config/schema.py

def to_dict(self) -> dict[str, Any]:
    """Convert config to dictionary."""
    return {
        "model": self.model,
        "simulation": {
            "start_date": self.simulation.start_date.isoformat(),
            "duration_hours": self.simulation.duration_hours,
            "coastal_domain": self.simulation.coastal_domain,
            "meteo_source": self.simulation.meteo_source,
            "timestep_seconds": self.simulation.timestep_seconds,
        },
        "boundary": {
            "source": self.boundary.source,
            "stofs_file": (str(self.boundary.stofs_file) if self.boundary.stofs_file else None),
            "tidal_model": self.boundary.tidal_model,
        },
        "paths": {
            "work_dir": str(self.paths.work_dir),
            "raw_download_dir": (
                str(self.paths.raw_download_dir) if self.paths.raw_download_dir else None
            ),
            "hot_start_file": (
                str(self.paths.hot_start_file) if self.paths.hot_start_file else None
            ),
            **({"parm_dir": str(self.paths.parm_dir)} if self.paths.parm_dir else {}),
            **({"nwm_dir": str(self.paths.nwm_dir)} if self.paths.nwm_dir else {}),
            **(
                {"tidal_atlas_dir": str(self.paths.tidal_atlas_dir)}
                if self.paths.tidal_atlas_dir
                else {}
            ),
        },
        "model_config": self.model_config.to_dict(),
        "monitoring": {
            "log_level": self.monitoring.log_level,
            "log_file": (str(self.monitoring.log_file) if self.monitoring.log_file else None),
            "enable_progress_tracking": self.monitoring.enable_progress_tracking,
            "enable_timing": self.monitoring.enable_timing,
        },
        "download": {
            "enabled": self.download.enabled,
            "timeout": self.download.timeout,
            "raise_on_error": self.download.raise_on_error,
            "limit_per_host": self.download.limit_per_host,
        },
    }

validate

validate()

Validate configuration and return list of errors.

Source code in src/coastal_calibration/config/schema.py

def validate(self) -> list[str]:
    """Validate configuration and return list of errors."""
    from coastal_calibration.data.downloader import validate_date_ranges

    errors: list[str] = []

    if self.simulation.duration_hours <= 0:
        errors.append("simulation.duration_hours must be positive")

    # Model-specific validation
    errors.extend(self.model_config.validate(self))

    # Shared boundary validation
    errors.extend(self._validate_boundary_source())

    # Date range validation
    if self.download.enabled:
        sim = self.simulation
        start_time = sim.start_date
        end_time = start_time + timedelta(hours=sim.duration_hours)
        date_errors = validate_date_ranges(
            start_time,
            end_time,
            sim.meteo_source,
            self.boundary.source,
            sim.coastal_domain,
        )
        errors.extend(date_errors)

    return errors

SimulationConfig

SimulationConfig dataclass

SimulationConfig(
    start_date,
    duration_hours,
    coastal_domain,
    meteo_source,
    timestep_seconds=200,
)

Simulation time and domain configuration.

start_date is normalized to naive UTC in __post_init__ so the rest of the pipeline can compare and serialize it without crossing the naive/aware boundary. Tz-aware values are converted to UTC then stripped; tz-naive values are passed through (assumed UTC by the project's data contract — NWM/STOFS are published on UTC days).

start_pdy property

start_pdy

Return start date as YYYYMMDD string.

start_cyc property

start_cyc

Return start cycle (hour) as HH string.

inland_domain property

inland_domain

Inland domain directory name for this coastal domain.

nwm_domain property

nwm_domain

NWM domain identifier for this coastal domain.

geo_grid property

geo_grid

Geogrid filename for this coastal domain.

meteo_crs property

meteo_crs

PROJ string for this domain's NWM meteorological forcing grid.

nwm_domain_for classmethod

nwm_domain_for(coastal_domain)

NWM domain identifier for an arbitrary coastal domain name.

Same mapping as :attr:nwm_domain, but callable without a configured simulation and lenient about names outside :data:CoastalDomain ("conus" is passed straight through, as the downloader accepts it). Used for naming the download cache, where the NWM domain is the right key: atlgulf and pacific pull byte-identical CONUS forcing and should share one copy.

Source code in src/coastal_calibration/config/schema.py

@classmethod
def nwm_domain_for(cls, coastal_domain: str) -> str:
    """NWM domain identifier for an arbitrary coastal domain name.

    Same mapping as :attr:`nwm_domain`, but callable without a
    configured simulation and lenient about names outside
    :data:`CoastalDomain` (``"conus"`` is passed straight through, as
    the downloader accepts it).  Used for naming the download cache,
    where the NWM domain is the right key: ``atlgulf`` and ``pacific``
    pull byte-identical CONUS forcing and should share one copy.
    """
    return cls._NWM_DOMAIN.get(coastal_domain, coastal_domain)

BoundaryConfig

BoundaryConfig dataclass

BoundaryConfig(
    source="harmonic",
    stofs_file=None,
    tidal_model="TPXO10-atlas-v2-nc",
)

Boundary condition configuration.

PARAMETER	DESCRIPTION
source	Boundary forcing source. harmonic predicts tides locally via pyTMD against the atlas at :attr:PathConfig.tidal_atlas_dir; stofs regrids the NOAA STOFS product (and falls back to harmonic past the STOFS 180 h window when the simulation runs longer). "tpxo" is accepted as a deprecated alias for "harmonic" and is normalized at construction time. TYPE: ('harmonic', 'stofs') DEFAULT: "harmonic"
stofs_file	STOFS NetCDF (only used when source == "stofs"). TYPE: Path DEFAULT: None
tidal_model	pyTMD model identifier (see pyTMD.io.load_database()). Defaults to TPXO10-atlas-v2 in netcdf form. Set to e.g. "FES2014" or "GOT4.10" to predict against another atlas without code changes — only the files under :attr:PathConfig.tidal_atlas_dir change. TYPE: str DEFAULT: 'TPXO10-atlas-v2-nc'

PathConfig

PathConfig dataclass

PathConfig(
    work_dir,
    raw_download_dir=None,
    hot_start_file=None,
    parm_dir=None,
    nwm_dir=None,
    tidal_atlas_dir=None,
)

Path configuration for data and executables.

Only work_dir is required. All other fields are optional and only needed by specific workflow stages.

parm_nwm property

parm_nwm

Parameter files directory (requires parm_dir).

download_dir property

download_dir

Effective download directory (fallback to work_dir/downloads).

meteo_subdir classmethod

meteo_subdir(meteo_source, coastal_domain)

Relative meteo path, meteo/<source>/<nwm domain>.

Every NWM domain names its hourly forcing YYYYMMDDHH.LDASIN_DOMAIN1, so files from different domains collide unless each domain gets its own directory: a cached Hawaii file would otherwise be served for a PRVI run covering the same hour.

Source code in src/coastal_calibration/config/schema.py

@classmethod
def meteo_subdir(cls, meteo_source: str, coastal_domain: str) -> Path:
    """Relative meteo path, ``meteo/<source>/<nwm domain>``.

    Every NWM domain names its hourly forcing ``YYYYMMDDHH.LDASIN_DOMAIN1``,
    so files from different domains collide unless each domain gets its
    own directory: a cached Hawaii file would otherwise be served for a
    PRVI run covering the same hour.
    """
    return (
        Path(cls.METEO_SUBDIR) / meteo_source / SimulationConfig.nwm_domain_for(coastal_domain)
    )

meteo_dir

meteo_dir(meteo_source, coastal_domain)

Directory for meteorological data.

Source code in src/coastal_calibration/config/schema.py

def meteo_dir(self, meteo_source: str, coastal_domain: str) -> Path:
    """Directory for meteorological data."""
    return self.download_dir / self.meteo_subdir(meteo_source, coastal_domain)

streamflow_subdir classmethod

streamflow_subdir(coastal_domain)

Relative streamflow path, hydro/nwm/<nwm domain>.

Keyed the same way as :meth:meteo_subdir, so a domain reads the same everywhere under the download directory even though NWM spells it differently in its own URLs (puertorico there, prvi here).

Only nwm_ana streamflow is downloaded; nwm_retro is read straight from the S3 Zarr store, so it has no directory here.

Source code in src/coastal_calibration/config/schema.py

@classmethod
def streamflow_subdir(cls, coastal_domain: str) -> Path:
    """Relative streamflow path, ``hydro/nwm/<nwm domain>``.

    Keyed the same way as :meth:`meteo_subdir`, so a domain reads the
    same everywhere under the download directory even though NWM
    spells it differently in its own URLs (``puertorico`` there,
    ``prvi`` here).

    Only ``nwm_ana`` streamflow is downloaded; ``nwm_retro`` is read
    straight from the S3 Zarr store, so it has no directory here.
    """
    return Path(cls.HYDRO_SUBDIR) / "nwm" / SimulationConfig.nwm_domain_for(coastal_domain)

streamflow_dir

streamflow_dir(coastal_domain='conus')

Directory for downloaded nwm_ana streamflow data.

Source code in src/coastal_calibration/config/schema.py

def streamflow_dir(self, coastal_domain: str = "conus") -> Path:
    """Directory for downloaded ``nwm_ana`` streamflow data."""
    return self.download_dir / self.streamflow_subdir(coastal_domain)

coastal_dir

coastal_dir(coastal_source)

Directory for coastal boundary data.

Source code in src/coastal_calibration/config/schema.py

def coastal_dir(self, coastal_source: str) -> Path:
    """Directory for coastal boundary data."""
    return self.download_dir / self.COASTAL_SUBDIR / coastal_source

geogrid_file

geogrid_file(sim)

Geogrid file path for the given domain (requires parm_dir).

Source code in src/coastal_calibration/config/schema.py

def geogrid_file(self, sim: SimulationConfig) -> Path:
    """Geogrid file path for the given domain (requires ``parm_dir``)."""
    return self.parm_nwm / sim.inland_domain / sim.geo_grid

ModelConfig

ModelConfig

Bases: ABC

Abstract base class for model-specific configuration.

Each concrete subclass owns its compute parameters, environment variable construction, stage ordering, validation, and SLURM script generation. This keeps model-specific concerns out of the shared configuration and makes adding new models straightforward: create a new subclass, implement the abstract methods, and register it in :data:MODEL_REGISTRY.

ATTRIBUTE	DESCRIPTION
omp_num_threads	Number of OpenMP threads per process. TYPE: int
runtime_env	Extra environment variables for the model run subprocess. Merged last so they can override any auto-detected value. Only used by model run stages (schism_run, sfincs_run). TYPE: dict[str, str]

model_name abstractmethod property

model_name

Return the model identifier string (e.g. 'schism', 'sfincs').

stage_order abstractmethod property

stage_order

Ordered list of stage names for this model's pipeline.

build_environment abstractmethod

build_environment(env, config)

Add model-specific environment variables to env (mutating).

Called by :meth:WorkflowStage.build_environment after shared variables (OpenMP pinning, HDF5 file locking) have been populated.

Source code in src/coastal_calibration/config/schema.py

@abstractmethod
def build_environment(self, env: dict[str, str], config: CoastalCalibConfig) -> dict[str, str]:
    """Add model-specific environment variables to *env* (mutating).

    Called by :meth:`WorkflowStage.build_environment` after shared
    variables (OpenMP pinning, HDF5 file locking) have been populated.
    """

validate abstractmethod

validate(config)

Return model-specific validation errors.

Source code in src/coastal_calibration/config/schema.py

@abstractmethod
def validate(self, config: CoastalCalibConfig) -> list[str]:
    """Return model-specific validation errors."""

create_stages abstractmethod

create_stages(config, monitor)

Construct and return the {name: stage} dictionary.

Source code in src/coastal_calibration/config/schema.py

@abstractmethod
def create_stages(self, config: CoastalCalibConfig, monitor: Any) -> dict[str, Any]:
    """Construct and return the ``{name: stage}`` dictionary."""

to_dict abstractmethod

to_dict()

Serialize model-specific fields to a dictionary.

Source code in src/coastal_calibration/config/schema.py

@abstractmethod
def to_dict(self) -> dict[str, Any]:
    """Serialize model-specific fields to a dictionary."""

SchismModelConfig

SchismModelConfig dataclass

SchismModelConfig(
    prebuilt_dir=None,
    geogrid_file=None,
    nodes=1,
    ntasks_per_node=0,
    exclusive=True,
    nscribes=0,
    omp_num_threads=2,
    oversubscribe=False,
    schism_exe=None,
    include_noaa_gages=False,
    discharge_file=None,
    create_water_level_animation=False,
    animation_fps=10,
    animation_time_stride=1,
    obs_points_csv=None,
    output_freq_hours=1.0,
    single_output_file=False,
    run_param_overrides=dict(),
    runtime_env=dict(),
)

Bases: ModelConfig

SCHISM model configuration.

Contains compute parameters (MPI layout, SCHISM binary), the path to a prebuilt model directory, and the geogrid file used for atmospheric forcing regridding.

PARAMETER	DESCRIPTION
prebuilt_dir	Path to the directory containing the pre-built SCHISM model files (hgrid.gr3, vgrid.in, param.nml, etc.). TYPE: Path DEFAULT: None
geogrid_file	Path to the WRF geogrid file (e.g. geo_em_HI.nc) used by the atmospheric forcing regridding stage. TYPE: Path DEFAULT: None
nodes	Number of SLURM nodes. Defaults to 1; set higher for multi-node HPC jobs. TYPE: int DEFAULT: 1
ntasks_per_node	MPI tasks per node. When <= 0 (the default), this is auto-set to get_cpu_count() // omp_num_threads so a single-node run fills the available physical cores (see :func:~coastal_calibration.utils.get_cpu_count). TYPE: int DEFAULT: 0
exclusive	Request exclusive node access. TYPE: bool DEFAULT: True
nscribes	Number of SCHISM scribe processes. When <= 0 (the default), this is auto-detected from the prebuilt's param.nml: a count of uncommented iof_*(N) = 1 flags plus one for iout_sta when include_noaa_gages is enabled. SCHISM aborts at init when nscribes is below the actual number of output variables. TYPE: int DEFAULT: 0
omp_num_threads	OpenMP threads per MPI rank. Defaults to 2 (typical SCHISM hybrid layout); combined with the auto-detected ntasks_per_node this fills one node's physical cores. TYPE: int DEFAULT: 2
oversubscribe	Pass --oversubscribe to mpiexec. Only honored under OpenMPI; silently ignored under MPICH (see :func:~coastal_calibration.utils.build_mpi_cmd). Defaults to False; set True when intentionally launching more MPI ranks than physical cores. TYPE: bool DEFAULT: False
schism_exe	Path to a compiled SCHISM executable. When set, the schism_run stage uses this binary instead of discovering pschism on PATH. Normally not needed -- SCHISM is compiled automatically when activating a pixi environment with the schism feature. Set this to a system-compiled binary on WCOSS2 or other clusters where the model is built against system MPI/HDF5/NetCDF. TYPE: Path DEFAULT: None
include_noaa_gages	When True, automatically query NOAA CO-OPS for water level stations within the model domain (computed from the concave hull of open boundary nodes in hgrid.gr3), write a station.in file, set iout_sta = 1 in param.nml, and generate sim-vs-obs comparison plots after the run. Requires the plot optional dependencies. TYPE: bool DEFAULT: False
discharge_file	Path to a nwmReaches.csv file mapping NWM reach feature IDs to SCHISM source/sink elements. When None (default), the discharge stage is skipped and no river forcing is generated. TYPE: Path DEFAULT: None
create_water_level_animation	When True, the schism_plot stage loads the 2-D elevation field from outputs/out2d_*.nc and renders an MP4 animation to figs/water_level.mp4 using :func:coastal_calibration.plotting.animate_water_level. Requires an ffmpeg binary on PATH. Independent of include_noaa_gages. Defaults to False. TYPE: bool DEFAULT: False
animation_fps	Frames per second for the animation output. TYPE: int DEFAULT: 10
animation_time_stride	Keep every animation_time_stride-th frame from the model time series; useful for long runs. TYPE: int DEFAULT: 1
obs_points_csv	Path to a CSV with columns id, lon, lat specifying extra observation points for water-level extraction after the run. The schism_plot stage interpolates water-surface elevation at each point (and at any NOAA CO-OPS gauges when include_noaa_gages is enabled) and writes the combined time series to obs_water_level.parquet in the work directory. TYPE: Path DEFAULT: None
output_freq_hours	How often SCHISM writes field outputs, in hours. Translated into the nspool parameter of param.nml. Defaults to 1.0 (hourly output, matching the previous hardcoded behavior). TYPE: float DEFAULT: 1.0
single_output_file	When True, set ihfskip to the full simulation length so SCHISM keeps appending to a single output file instead of rotating to a new file every nspool steps. Useful on shared filesystems where each file rotation costs an MPI barrier and metadata round-trips. Defaults to False (matching the previous ihfskip = nspool behavior). TYPE: bool DEFAULT: False
run_param_overrides	Arbitrary key/value pairs written into param.nml after the template values, output_freq_hours, and single_output_file have been applied. Use this to override any other namelist parameter (e.g. {"dt": 100, "iwbl": 1}). Mirrors the SFINCS run_param_overrides option. Validation catches mismatches between ihfskip, nhot_write, and nspool_sta before SCHISM is launched. TYPE: dict DEFAULT: dict()

total_tasks property

total_tasks

Total number of MPI tasks (nodes * ntasks_per_node).

coastal_parm property

coastal_parm

Directory containing prebuilt SCHISM model files.

geogrid_path property

geogrid_path

WRF geogrid file used for atmospheric regridding.

schism_mesh property

schism_mesh

SCHISM ESMF mesh file path.

resolved_discharge_file property

resolved_discharge_file

Resolve the NWM-reaches discharge CSV, or None to skip discharge.

Resolution order:

1. discharge_file explicitly set → use it if it exists, else None. An explicit configuration is treated as exclusive: it does not silently fall back to the prebuilt-directory convention.
2. discharge_file unset → look for nwmReaches.csv next to the prebuilt model (matching the Pacific/Hawaii/PRVI/AtlGulf convention). Use it if present.
3. Otherwise → None (river forcing is skipped and SCHISM is configured with if_source = 0).

A missing optional file degrades gracefully — the caller skips discharge rather than aborting.

elevation_correction_csv property

elevation_correction_csv

Return elevation_correction.csv next to the prebuilt model.

Returns None when prebuilt_dir is unset or when the correction file is absent. Callers can pass the result straight into readers that accept an optional path without re-doing the exists() check.

SfincsModelConfig

SfincsModelConfig dataclass

SfincsModelConfig(
    prebuilt_dir,
    model_root=None,
    discharge_locations_file=None,
    merge_discharge=False,
    include_precip=False,
    include_wind=False,
    include_pressure=False,
    meteo_res=None,
    forcing_to_mesh_offset_m=0.0,
    vdatum_mesh_to_msl_m=0.0,
    sfincs_exe=None,
    omp_num_threads=0,
    run_param_overrides=dict(),
    floodmap_dem=None,
    floodmap_hmin=0.05,
    floodmap_enabled=True,
    floodmap_land_only=True,
    create_water_level_animation=False,
    animation_fps=10,
    animation_time_stride=1,
    obs_points_csv=None,
    runtime_env=dict(),
)

Bases: ModelConfig

SFINCS model configuration.

SFINCS runs on a single node using OpenMP (all available cores). There is no MPI or multi-node support.

PARAMETER	DESCRIPTION
prebuilt_dir	Path to the directory containing the pre-built model files (sfincs.inp, sfincs.nc, region.geojson, etc.). TYPE: Path
model_root	Output directory for the built model. Defaults to {work_dir}/sfincs_model. TYPE: Path DEFAULT: None
discharge_locations_file	Path to a SFINCS .src or GeoJSON with discharge source point locations. TYPE: Path DEFAULT: None
merge_discharge	Whether to merge with pre-existing discharge source points. TYPE: bool DEFAULT: False
include_precip	When True, add precipitation forcing from the meteorological data catalog entry (derived from simulation.meteo_source). TYPE: bool DEFAULT: False
include_wind	When True, add spatially-varying wind forcing (wind10_u, wind10_v) from the meteorological data catalog entry. TYPE: bool DEFAULT: False
include_pressure	When True, add spatially-varying atmospheric pressure forcing (press_msl) and enable barometric correction (baro=1). TYPE: bool DEFAULT: False
meteo_res	Output resolution (m) for gridded meteorological forcing (precipitation, wind, pressure). When None (default) the resolution is determined from the SFINCS quadtree grid — it equals the base cell size (coarsest level) so that the meteo grid is never finer than needed. Setting an explicit value (e.g. 2000) overrides the automatic calculation. .. note:: Without this parameter the HydroMT reproject call retains the source-data resolution (≈ 1 km for NWM), and the LCC → UTM reprojection can inflate the output to the full CONUS extent, producing multi-GB files and very slow simulations. TYPE: float DEFAULT: None
forcing_to_mesh_offset_m	Vertical offset in meters added to the boundary-condition water levels before they enter SFINCS. Tidal-only sources (harmonic prediction) provide oscillations centered on zero (MSL) but carry no information about where MSL sits on the mesh's vertical datum. This parameter anchors the forcing signal to the correct geodetic height on the mesh. Set it to the elevation of MSL in the mesh datum obtained from VDatum (e.g. 0.171 for a NAVD88 mesh on the Texas Gulf coast, where MSL is 0.171 m above NAVD88). For sources that already report water levels in the mesh datum (e.g. STOFS on a NAVD88 mesh) set this to 0.0. Defaults to 0.0. TYPE: float DEFAULT: 0.0
vdatum_mesh_to_msl_m	Vertical offset in meters added to the simulated water level before comparison with NOAA CO-OPS observations (which are in MSL). The model output inherits the mesh vertical datum, so this converts it to MSL (e.g. 0.171 for a NAVD88 mesh on the Texas Gulf coast). Defaults to 0.0. TYPE: float DEFAULT: 0.0
sfincs_exe	Path to a compiled SFINCS executable. When set, the sfincs_run stage uses this binary instead of discovering sfincs on PATH. Normally not needed -- SFINCS is compiled automatically when activating a pixi environment with the sfincs feature. TYPE: Path DEFAULT: None
omp_num_threads	Number of OpenMP threads. Defaults to the number of physical CPU cores on the current machine (see :func:~coastal_calibration.utils.get_cpu_count). On HPC nodes this auto-detects correctly; on a local laptop it avoids over-subscribing the system. TYPE: int DEFAULT: 0
run_param_overrides	Arbitrary key/value pairs written to sfincs.inp just before the model is written to disk. Use this to override physics parameters that HydroMT-SFINCS sets by default (e.g. advection: 0, nuvisc: 0.01). Keys must be valid sfincs.inp parameter names. Mirrors the SCHISM run_param_overrides option. TYPE: dict DEFAULT: dict()
create_water_level_animation	When True, the sfincs_plot stage loads the time-dependent water level field from sfincs_map.nc and renders an MP4 animation to figs/water_level.mp4 using :func:coastal_calibration.plotting.animate_water_level. Requires an ffmpeg binary on PATH. Defaults to False. TYPE: bool DEFAULT: False
animation_fps	Frames per second for the animation output. TYPE: int DEFAULT: 10
animation_time_stride	Keep every animation_time_stride-th frame from the model time series; useful for long runs. TYPE: int DEFAULT: 1
obs_points_csv	Path to a CSV with columns id, lon, lat specifying extra observation points for water-level extraction after the run. The sfincs_plot stage interpolates water-surface elevation at each point (and at any NOAA CO-OPS gauges found in obs_station_map.json) and writes the combined time series to obs_water_level.parquet alongside sfincs_map.nc. TYPE: Path DEFAULT: None

MonitoringConfig

MonitoringConfig dataclass

MonitoringConfig(
    log_level="INFO",
    log_file=None,
    enable_progress_tracking=True,
    enable_timing=True,
)

Workflow monitoring configuration.

DownloadConfig

DownloadConfig dataclass

DownloadConfig(
    enabled=True,
    timeout=600,
    raise_on_error=True,
    limit_per_host=4,
)

Data download configuration.

SFINCS Creation Configuration

SfincsCreateConfig

SfincsCreateConfig dataclass

SfincsCreateConfig(
    aoi,
    output_dir,
    download_dir=None,
    grid=GridConfig(),
    elevation=ElevationConfig(),
    mask=MaskConfig(),
    subgrid=SubgridConfig(),
    data_catalog=DataCatalogConfig(),
    monitoring=MonitoringConfig(),
    river_discharge=None,
    aoi_simplify_neck_m=0.0,
    add_noaa_gages=False,
    observation_points=list(),
    observation_locations_file=None,
    merge_observations=False,
    obs_snap_depth_threshold=-2.0,
    obs_snap_search_radius_m=1000.0,
)

Root configuration for SFINCS model creation workflow.

Loaded from YAML via :meth:from_yaml. All paths are resolved to absolute paths during construction.

stage_order property

stage_order

Ordered list of creation stages to execute.

Roughness is embedded in the quadtree subgrid tables, so there is no separate roughness stage. The create_discharge stage is included only when :attr:river_discharge is configured.

from_yaml classmethod

from_yaml(config_path)

Load configuration from a YAML file.

PARAMETER	DESCRIPTION
config_path	Path to YAML configuration file. TYPE: Path or str

RETURNS	DESCRIPTION
SfincsCreateConfig	Loaded configuration.

RAISES	DESCRIPTION
FileNotFoundError	If the configuration file does not exist.
YAMLError	If the YAML file is malformed.

Source code in src/coastal_calibration/config/create_schema.py

@classmethod
def from_yaml(cls, config_path: Path | str) -> SfincsCreateConfig:
    """Load configuration from a YAML file.

    Parameters
    ----------
    config_path : Path or str
        Path to YAML configuration file.

    Returns
    -------
    SfincsCreateConfig
        Loaded configuration.

    Raises
    ------
    FileNotFoundError
        If the configuration file does not exist.
    yaml.YAMLError
        If the YAML file is malformed.
    """
    config_path = Path(config_path)
    if not config_path.exists():
        raise FileNotFoundError(f"Configuration file not found: {config_path}")

    try:
        data = yaml.safe_load(config_path.read_text())
    except yaml.YAMLError as e:
        raise yaml.YAMLError(f"Invalid YAML in {config_path}: {e}") from e

    if data is None:
        raise ValueError(f"Configuration file is empty: {config_path}")

    cls._resolve_relative_paths(data, config_path.parent)
    return cls.from_dict(data)

from_dict classmethod

from_dict(data)

Create config from a plain dictionary.

PARAMETER	DESCRIPTION
data	Configuration dictionary with the same structure as the YAML file (see :meth:to_dict for the expected keys). TYPE: dict

RETURNS	DESCRIPTION
SfincsCreateConfig

Source code in src/coastal_calibration/config/create_schema.py

@classmethod
def from_dict(cls, data: dict[str, Any]) -> SfincsCreateConfig:
    """Create config from a plain dictionary.

    Parameters
    ----------
    data : dict
        Configuration dictionary with the same structure as the YAML
        file (see :meth:`to_dict` for the expected keys).

    Returns
    -------
    SfincsCreateConfig
    """
    aoi = data.get("aoi")
    if aoi is None:
        raise ValueError("'aoi' is required (path to AOI polygon file)")
    output_dir = data.get("output_dir")
    if output_dir is None:
        raise ValueError("'output_dir' is required (model output directory)")

    grid_data = data.get("grid", {})
    refinement_raw = grid_data.pop("refinement", None)
    if refinement_raw is not None:
        grid_data["refinement"] = [RefinementLevel(**r) for r in refinement_raw]
    grid = GridConfig(**grid_data)

    elev_data = data.get("elevation", {})
    datasets_raw = elev_data.pop("datasets", None)
    if datasets_raw is not None:
        elev_data["datasets"] = [ElevationDataset(**d) for d in datasets_raw]
    elevation = ElevationConfig(**elev_data)

    mask_data = data.get("mask", {})
    mask = MaskConfig(**mask_data)

    subgrid_data = data.get("subgrid", {})
    if subgrid_data.get("reclass_table"):
        subgrid_data["reclass_table"] = Path(subgrid_data["reclass_table"])
    subgrid_cfg = SubgridConfig(**subgrid_data)

    catalog_data = data.get("data_catalog", {})
    data_catalog = DataCatalogConfig(**catalog_data)

    monitoring_data = data.get("monitoring", {})
    if monitoring_data.get("log_file"):
        monitoring_data["log_file"] = Path(monitoring_data["log_file"])
    monitoring = MonitoringConfig(**monitoring_data)

    river_discharge: RiverDischargeConfig | None = None
    nwm_data = data.get("river_discharge")
    if nwm_data is not None:
        nwm_data["flowlines"] = Path(nwm_data["flowlines"])
        river_discharge = RiverDischargeConfig(**nwm_data)

    download_dir_raw = data.get("download_dir")
    download_dir = Path(download_dir_raw) if download_dir_raw else None

    add_noaa_gages = data.get("add_noaa_gages", False)
    observation_points = data.get("observation_points", [])
    obs_file_raw = data.get("observation_locations_file")
    observation_locations_file = Path(obs_file_raw) if obs_file_raw else None
    merge_observations = data.get("merge_observations", False)
    aoi_simplify_neck_m = float(data.get("aoi_simplify_neck_m", 0.0))
    obs_snap_depth_threshold = float(data.get("obs_snap_depth_threshold", -2.0))
    obs_snap_search_radius_m = float(data.get("obs_snap_search_radius_m", 1000.0))

    return cls(
        aoi=Path(aoi),
        output_dir=Path(output_dir),
        download_dir=download_dir,
        grid=grid,
        elevation=elevation,
        mask=mask,
        subgrid=subgrid_cfg,
        data_catalog=data_catalog,
        monitoring=monitoring,
        river_discharge=river_discharge,
        aoi_simplify_neck_m=aoi_simplify_neck_m,
        add_noaa_gages=add_noaa_gages,
        observation_points=observation_points,
        observation_locations_file=observation_locations_file,
        merge_observations=merge_observations,
        obs_snap_depth_threshold=obs_snap_depth_threshold,
        obs_snap_search_radius_m=obs_snap_search_radius_m,
    )

to_yaml

to_yaml(path)

Write configuration to a YAML file.

PARAMETER	DESCRIPTION
path	Path to YAML output file. Parent directories are created automatically. TYPE: Path or str

Source code in src/coastal_calibration/config/create_schema.py

def to_yaml(self, path: Path | str) -> None:
    """Write configuration to a YAML file.

    Parameters
    ----------
    path : Path or str
        Path to YAML output file.  Parent directories are created
        automatically.
    """
    path = Path(path)
    path.parent.mkdir(parents=True, exist_ok=True)
    path.write_text(yaml.dump(self.to_dict(), default_flow_style=False, sort_keys=False))

to_dict

to_dict()

Convert configuration to a plain dictionary.

Source code in src/coastal_calibration/config/create_schema.py

def to_dict(self) -> dict[str, Any]:
    """Convert configuration to a plain dictionary."""
    return {
        "aoi": str(self.aoi),
        "output_dir": str(self.output_dir),
        **({"download_dir": str(self.download_dir)} if self.download_dir else {}),
        "grid": {
            "resolution": self.grid.resolution,
            "crs": self.grid.crs,
            "rotated": self.grid.rotated,
            "refinement": [
                {
                    "polygon": str(r.polygon),
                    "level": r.level,
                    **({"buffer_m": r.buffer_m} if r.buffer_m != 0.0 else {}),
                }
                for r in self.grid.refinement
            ],
        },
        "elevation": {
            "datasets": [
                {
                    "name": d.name,
                    "zmin": d.zmin,
                    **({"source": d.source} if d.source else {}),
                    **({"noaa_dataset": d.noaa_dataset} if d.noaa_dataset else {}),
                    **({"coastal_domain": d.coastal_domain} if d.coastal_domain else {}),
                    **({"offset": d.offset} if d.offset else {}),
                }
                for d in self.elevation.datasets
            ],
            "buffer_cells": self.elevation.buffer_cells,
        },
        "mask": {
            "zmin": self.mask.zmin,
            "boundary_zmax": self.mask.boundary_zmax,
            "reset_bounds": self.mask.reset_bounds,
            "keep_largest_only": self.mask.keep_largest_only,
        },
        "subgrid": {
            "nr_subgrid_pixels": self.subgrid.nr_subgrid_pixels,
            "lulc_dataset": self.subgrid.lulc_dataset,
            **({"lulc_source": self.subgrid.lulc_source} if self.subgrid.lulc_source else {}),
            "reclass_table": (
                str(self.subgrid.reclass_table) if self.subgrid.reclass_table else None
            ),
            "manning_land": self.subgrid.manning_land,
            "manning_sea": self.subgrid.manning_sea,
        },
        "data_catalog": {
            "data_libs": self.data_catalog.data_libs,
        },
        "monitoring": {
            "log_level": self.monitoring.log_level,
            "log_file": (str(self.monitoring.log_file) if self.monitoring.log_file else None),
            "enable_progress_tracking": self.monitoring.enable_progress_tracking,
            "enable_timing": self.monitoring.enable_timing,
        },
        "river_discharge": (
            {
                "flowlines": str(self.river_discharge.flowlines),
                "nwm_id_column": self.river_discharge.nwm_id_column,
                "max_snap_distance_m": self.river_discharge.max_snap_distance_m,
            }
            if self.river_discharge is not None
            else None
        ),
        "aoi_simplify_neck_m": self.aoi_simplify_neck_m,
        "add_noaa_gages": self.add_noaa_gages,
        "observation_points": self.observation_points,
        "observation_locations_file": (
            str(self.observation_locations_file) if self.observation_locations_file else None
        ),
        "merge_observations": self.merge_observations,
        "obs_snap_depth_threshold": self.obs_snap_depth_threshold,
        "obs_snap_search_radius_m": self.obs_snap_search_radius_m,
    }

validate

validate()

Validate configuration and return a list of error messages.

RETURNS	DESCRIPTION
list of str	Validation errors (empty when the config is valid).

Source code in src/coastal_calibration/config/create_schema.py

def validate(self) -> list[str]:
    """Validate configuration and return a list of error messages.

    Returns
    -------
    list of str
        Validation errors (empty when the config is valid).
    """
    errors: list[str] = []

    if not self.aoi.exists():
        errors.append(f"AOI file not found: {self.aoi}")

    if self.grid.resolution <= 0:
        errors.append(f"grid.resolution must be positive, got {self.grid.resolution}")

    errors.extend(self._validate_elevation())
    errors.extend(self._validate_subgrid())

    for ref in self.grid.refinement:
        if not ref.polygon.exists():
            errors.append(f"refinement polygon not found: {ref.polygon}")
        if ref.level < 1:
            errors.append(f"refinement level must be >= 1, got {ref.level}")

    if self.river_discharge is not None:
        nd = self.river_discharge
        if not nd.flowlines.exists():
            errors.append(f"river_discharge.flowlines not found: {nd.flowlines}")

    return errors

GridConfig

GridConfig dataclass

GridConfig(
    resolution=50.0,
    crs="utm",
    rotated=True,
    refinement=list(),
)

Grid generation configuration.

ElevationConfig

ElevationConfig dataclass

ElevationConfig(
    datasets=(
        lambda: [
            ElevationDataset(
                name="copdem_30m",
                zmin=0.001,
                source="copdem_30m",
            ),
            ElevationDataset(
                name="gebco_15arcs",
                zmin=-20000,
                source="gebco_15arcs",
            ),
        ]
    )(),
    buffer_cells=1,
)

Elevation and bathymetry configuration.

MaskConfig

MaskConfig dataclass

MaskConfig(
    zmin=-5.0,
    boundary_zmax=-5.0,
    reset_bounds=True,
    keep_largest_only=False,
)

Active-cell mask and boundary configuration.

SubgridConfig

SubgridConfig dataclass

SubgridConfig(
    nr_subgrid_pixels=5,
    lulc_dataset="esa_worldcover",
    lulc_source="esa_worldcover",
    reclass_table=None,
    manning_land=0.04,
    manning_sea=0.02,
)

Subgrid table configuration.

Roughness parameters are included here because for quadtree grids the Manning coefficients are embedded directly in the subgrid tables.

RiverDischargeConfig

RiverDischargeConfig dataclass

RiverDischargeConfig(
    flowlines, nwm_id_column, max_snap_distance_m=2000.0
)

River discharge source point configuration.

Derives discharge source points from user-provided flowline geometries (e.g. exported from the QGIS plugin). Each flowline's downstream endpoint (closest to the AOI boundary) is registered as a SFINCS discharge source location.

Workflow Runners

CoastalCalibRunner

CoastalCalibRunner

CoastalCalibRunner(config)

Main workflow runner for coastal model calibration.

This class orchestrates the entire calibration workflow, managing stage execution and progress monitoring.

Supports both SCHISM (model="schism", default) and SFINCS (model="sfincs") pipelines. The model type is selected via config.model.

Initialize the workflow runner.

PARAMETER	DESCRIPTION
config	Coastal calibration configuration. TYPE: CoastalCalibConfig

Source code in src/coastal_calibration/runner.py

def __init__(self, config: CoastalCalibConfig) -> None:
    """Initialize the workflow runner.

    Parameters
    ----------
    config : CoastalCalibConfig
        Coastal calibration configuration.
    """
    self.config = config

    # Ensure log directory exists early so file logging can start.
    config.paths.work_dir.mkdir(parents=True, exist_ok=True)

    # Set up file logging *before* creating the monitor so that
    # every message (including third-party) is captured on disk.
    if not config.monitoring.log_file:
        log_path = generate_log_path(config.paths.work_dir)
        configure_logger(file=str(log_path), file_level="DEBUG")

    # Silence noisy third-party loggers (HydroMT, xarray, ...)
    silence_third_party_loggers()

    # Recover from `srun --pty bash` constraining the step to one CPU:
    # if SLURM allocated more than the inherited mask allows, expand
    # so OpenMP threads can actually spread.  No-op outside SLURM.
    expand_cpu_affinity_if_constrained()

    self.monitor = WorkflowMonitor(config.monitoring)
    self._stages: dict[str, WorkflowStage] = {}
    self._results: dict[str, Any] = {}

validate

validate()

Validate configuration and prerequisites.

RETURNS	DESCRIPTION
list of str	List of validation error messages (empty if valid).

Source code in src/coastal_calibration/runner.py

def validate(self) -> list[str]:
    """Validate configuration and prerequisites.

    Returns
    -------
    list of str
        List of validation error messages (empty if valid).
    """
    errors = []

    config_errors = self.config.validate()
    errors.extend(config_errors)

    self._init_stages()
    for name, stage in self._stages.items():
        stage_errors = stage.validate()
        errors.extend(f"[{name}] {error}" for error in stage_errors)

    return errors

run

run(start_from=None, stop_after=None, dry_run=False)

Execute the calibration workflow.

PARAMETER	DESCRIPTION
start_from	Stage name to start from (skip earlier stages). TYPE: str DEFAULT: None
stop_after	Stage name to stop after (skip later stages). TYPE: str DEFAULT: None
dry_run	If True, validate but don't execute. TYPE: bool DEFAULT: False

RETURNS	DESCRIPTION
WorkflowResult	Result with execution details.

Source code in src/coastal_calibration/runner.py

def run(
    self,
    start_from: str | None = None,
    stop_after: str | None = None,
    dry_run: bool = False,
) -> WorkflowResult:
    """Execute the calibration workflow.

    Parameters
    ----------
    start_from : str, optional
        Stage name to start from (skip earlier stages).
    stop_after : str, optional
        Stage name to stop after (skip later stages).
    dry_run : bool, default False
        If True, validate but don't execute.

    Returns
    -------
    WorkflowResult
        Result with execution details.
    """
    start_time = utc_now()
    stages_completed: list[str] = []
    stages_failed: list[str] = []
    outputs: dict[str, Any] = {}
    errors: list[str] = []

    validation_errors = self.validate()
    # When resuming mid-pipeline, verify that earlier stages completed.
    if not validation_errors and start_from:
        validation_errors = self._check_prerequisites(start_from)
    if validation_errors:
        return WorkflowResult(
            success=False,
            job_id=None,
            start_time=start_time,
            end_time=utc_now(),
            stages_completed=[],
            stages_failed=[],
            outputs={},
            errors=validation_errors,
        )

    if dry_run:
        self.monitor.info("Dry run mode - validation passed, no execution")
        return WorkflowResult(
            success=True,
            job_id=None,
            start_time=start_time,
            end_time=utc_now(),
            stages_completed=[],
            stages_failed=[],
            outputs={"dry_run": True},
            errors=[],
        )

    self.monitor.register_stages(self.STAGE_ORDER)
    self.monitor.start_workflow()
    self.monitor.info("-" * 40)

    # Clean generated files from previous runs when starting fresh.
    # When resuming (start_from is set), preserve existing outputs.
    if not start_from:
        from coastal_calibration.schism.prep import clean_run_directory

        clean_run_directory(self.config.paths.work_dir)

    stages_to_run = self._get_stages_to_run(start_from, stop_after)

    current_stage = ""
    try:
        for current_stage in stages_to_run:
            stage = self._stages[current_stage]

            with self.monitor.stage_context(current_stage, stage.description):
                result = stage.run()
                self._results[current_stage] = result
                outputs[current_stage] = result
                stages_completed.append(current_stage)
                self._save_stage_status(current_stage)

        self.monitor.end_workflow(success=True)
        success = True

    except Exception as e:
        self.monitor.error(f"Workflow failed: {e}")
        self.monitor.end_workflow(success=False)
        errors.append(str(e))
        stages_failed.append(current_stage)
        success = False

    result = WorkflowResult(
        success=success,
        job_id=None,
        start_time=start_time,
        end_time=utc_now(),
        stages_completed=stages_completed,
        stages_failed=stages_failed,
        outputs=outputs,
        errors=errors,
    )

    result_file = self.config.paths.work_dir / "workflow_result.json"
    result.save(result_file)
    self.monitor.save_progress(self.config.paths.work_dir / "workflow_progress.json")

    return result

SfincsCreator

SfincsCreator

SfincsCreator(config)

Runner that orchestrates the SFINCS model creation pipeline.

Mirrors :class:~coastal_calibration.runner.CoastalCalibRunner but operates on a :class:SfincsCreateConfig and delegates to :class:~coastal_calibration.sfincs.create.CreateStage instances.

Source code in src/coastal_calibration/sfincs/create.py

def __init__(self, config: SfincsCreateConfig) -> None:
    self.config = config
    self._pending_status: dict[str, str] = {}

    # Ensure the output directory exists early so file logging can start.
    config.output_dir.mkdir(parents=True, exist_ok=True)

    # Set up file logging before creating the monitor.
    if not config.monitoring.log_file:
        log_path = generate_log_path(config.output_dir, prefix="sfincs-create")
        configure_logger(file=str(log_path), file_level="DEBUG")

    silence_third_party_loggers()

    self.monitor = WorkflowMonitor(config.monitoring)
    self._stages: dict[str, CreateStage] = {}
    self._results: dict[str, Any] = {}

run

run(start_from=None, stop_after=None, dry_run=False)

Execute the SFINCS model creation workflow.

PARAMETER	DESCRIPTION
start_from	Stage name to start from (skip earlier stages). TYPE: str DEFAULT: None
stop_after	Stage name to stop after (skip later stages). TYPE: str DEFAULT: None
dry_run	If True, validate but don't execute. TYPE: bool DEFAULT: False

RETURNS	DESCRIPTION
WorkflowResult	Result with execution details.

Source code in src/coastal_calibration/sfincs/create.py

def run(
    self,
    start_from: str | None = None,
    stop_after: str | None = None,
    dry_run: bool = False,
) -> WorkflowResult:
    """Execute the SFINCS model creation workflow.

    Parameters
    ----------
    start_from : str, optional
        Stage name to start from (skip earlier stages).
    stop_after : str, optional
        Stage name to stop after (skip later stages).
    dry_run : bool, default False
        If True, validate but don't execute.

    Returns
    -------
    WorkflowResult
        Result with execution details.
    """
    start_time = utc_now()
    stages_completed: list[str] = []
    stages_failed: list[str] = []
    outputs: dict[str, Any] = {}
    errors: list[str] = []

    validation_errors = self.validate()
    if not validation_errors and start_from:
        validation_errors = self._check_prerequisites(start_from)
    if validation_errors:
        return WorkflowResult(
            success=False,
            job_id=None,
            start_time=start_time,
            end_time=utc_now(),
            stages_completed=[],
            stages_failed=[],
            outputs={},
            errors=validation_errors,
        )

    if dry_run:
        self.monitor.info("Dry run mode - validation passed, no execution")
        return WorkflowResult(
            success=True,
            job_id=None,
            start_time=start_time,
            end_time=utc_now(),
            stages_completed=[],
            stages_failed=[],
            outputs={"dry_run": True},
            errors=[],
        )

    self.monitor.register_stages(self.stage_order)
    self.monitor.start_workflow()
    self.monitor.info("-" * 40)

    stages_to_run = self._get_stages_to_run(start_from, stop_after)

    # When resuming from a later stage, load the existing model so
    # that stages which reference ``self.sfincs`` can find it.
    if start_from and "create_grid" not in stages_to_run:
        if "create_fetch_data" not in stages_to_run:
            self._register_fetched_catalogs()
        _load_existing_model(self.config)

    current_stage = ""
    try:
        with suppress_hydromt_output():
            for current_stage in stages_to_run:
                stage = self._stages[current_stage]

                with self.monitor.stage_context(current_stage, stage.description):
                    result = stage.run()
                    self._results[current_stage] = result
                    outputs[current_stage] = result
                    stages_completed.append(current_stage)
                    self._record_stage(current_stage)
                    if current_stage == _WRITE_STAGE:
                        self._commit_status()

        self.monitor.end_workflow(success=True)
        success = True

    except Exception as e:
        self.monitor.error(f"Workflow failed: {e}")
        self.monitor.end_workflow(success=False)
        errors.append(str(e))
        stages_failed.append(current_stage)
        success = False

    result = WorkflowResult(
        success=success,
        job_id=None,
        start_time=start_time,
        end_time=utc_now(),
        stages_completed=stages_completed,
        stages_failed=stages_failed,
        outputs=outputs,
        errors=errors,
    )

    self._save_result(result, start_from)
    self.monitor.save_progress(self.config.output_dir / "create_progress.json")

    return result

WorkflowResult

WorkflowResult dataclass

WorkflowResult(
    success,
    job_id,
    start_time,
    end_time,
    stages_completed,
    stages_failed,
    outputs,
    errors,
)

Result of a workflow execution.

duration_seconds property

duration_seconds

Get workflow duration in seconds.

to_dict

to_dict()

Convert to dictionary.

Source code in src/coastal_calibration/runner.py

def to_dict(self) -> dict[str, Any]:
    """Convert to dictionary."""
    return {
        "success": self.success,
        "job_id": self.job_id,
        "start_time": self.start_time.isoformat(),
        "end_time": self.end_time.isoformat() if self.end_time else None,
        "duration_seconds": self.duration_seconds,
        "stages_completed": self.stages_completed,
        "stages_failed": self.stages_failed,
        "outputs": self.outputs,
        "errors": self.errors,
    }

save

save(path)

Save result to JSON file.

PARAMETER	DESCRIPTION
path	Path to output JSON file. Parent directories will be created if they don't exist. TYPE: Path or str

Source code in src/coastal_calibration/runner.py

def save(self, path: Path | str) -> None:
    """Save result to JSON file.

    Parameters
    ----------
    path : Path or str
        Path to output JSON file. Parent directories will be created
        if they don't exist.
    """
    path = Path(path)
    path.parent.mkdir(parents=True, exist_ok=True)
    path.write_text(json.dumps(self.to_dict(), indent=2))

Plotting

SfincsGridInfo

SfincsGridInfo dataclass

SfincsGridInfo(
    grid_type,
    crs,
    base_resolution,
    levels,
    n_faces=None,
    n_edges=None,
    shape=None,
    _verts=None,
    _level_per_face=None,
    _mask=None,
    _grid_extent=None,
)

Summary of a SFINCS model grid.

Use :meth:from_model_root to construct from a SFINCS model directory. The instance carries enough pre-computed state to drive :func:plot_mesh without re-loading the model.

Examples:

>>> info = SfincsGridInfo.from_model_root("run/sfincs_model")
>>> print(info)
SfincsGridInfo(quadtree, EPSG:32619)
  Faces:     293,850
  Edges:     596,123
  Level 1:    7,090 cells (512 m)
  ...

from_model_root classmethod

from_model_root(model_root)

Load grid metadata from a SFINCS model directory.

PARAMETER	DESCRIPTION
model_root	Path to the SFINCS model directory (must contain sfincs.inp and, for quadtree models, sfincs.nc). TYPE: Path | str

Source code in src/coastal_calibration/sfincs/plotting.py

@classmethod
def from_model_root(
    cls,
    model_root: Path | str,
) -> SfincsGridInfo:
    """Load grid metadata from a SFINCS model directory.

    Parameters
    ----------
    model_root:
        Path to the SFINCS model directory (must contain ``sfincs.inp``
        and, for quadtree models, ``sfincs.nc``).
    """
    from coastal_calibration.logging import suppress_hydromt_output
    from coastal_calibration.sfincs._hydromt_compat import apply_all_patches

    apply_all_patches()

    with suppress_hydromt_output():
        from hydromt_sfincs import SfincsModel

        sf = SfincsModel(root=str(model_root), mode="r+")
        sf.read()
    return cls._from_loaded_model(sf)

plot_mesh

plot_mesh

plot_mesh(
    info,
    *,
    ax=None,
    title=None,
    basemap=True,
    basemap_source=None,
    basemap_zoom=11,
    figsize=(11, 7),
)

Plot the SFINCS mesh colored by refinement level.

PARAMETER	DESCRIPTION
info	Grid metadata from :meth:SfincsGridInfo.from_model_root. TYPE: SfincsGridInfo
ax	Existing axes to plot into. A new figure is created when None. TYPE: Axes | None DEFAULT: None
title	Plot title. Defaults to a description derived from info. TYPE: str | None DEFAULT: None
basemap	If True (default), overlay satellite imagery via contextily. TYPE: bool DEFAULT: True
basemap_source	Tile provider passed to contextily.add_basemap. Defaults to cx.providers.Esri.WorldImagery. TYPE: Any | None DEFAULT: None
basemap_zoom	Zoom level for the basemap tiles. TYPE: int DEFAULT: 11
figsize	Figure size when ax is None. TYPE: tuple[float, float] DEFAULT: (11, 7)

RETURNS	DESCRIPTION
(Figure, Axes)

Source code in src/coastal_calibration/sfincs/plotting.py

def plot_mesh(
    info: SfincsGridInfo,
    *,
    ax: Axes | None = None,
    title: str | None = None,
    basemap: bool = True,
    basemap_source: Any | None = None,
    basemap_zoom: int = 11,
    figsize: tuple[float, float] = (11, 7),
) -> tuple[Figure, Axes]:
    """Plot the SFINCS mesh colored by refinement level.

    Parameters
    ----------
    info:
        Grid metadata from :meth:`SfincsGridInfo.from_model_root`.
    ax:
        Existing axes to plot into.  A new figure is created when *None*.
    title:
        Plot title.  Defaults to a description derived from *info*.
    basemap:
        If *True* (default), overlay satellite imagery via *contextily*.
    basemap_source:
        Tile provider passed to ``contextily.add_basemap``.  Defaults to
        ``cx.providers.Esri.WorldImagery``.
    basemap_zoom:
        Zoom level for the basemap tiles.
    figsize:
        Figure size when *ax* is *None*.

    Returns
    -------
    (Figure, Axes)
    """
    import matplotlib.pyplot as plt

    if ax is None:
        fig, ax = plt.subplots(figsize=figsize)
    else:
        fig = ax.get_figure()
        if fig is None:  # pragma: no cover
            msg = "ax must be attached to a Figure"
            raise ValueError(msg)

    if info.grid_type == "quadtree":
        _plot_quadtree(info, ax)
    else:
        _plot_regular(info, ax)

    if title is None:
        title = f"SFINCS {info.grid_type} mesh ({info.base_resolution:.0f} m)"
    ax.set_title(title)

    if basemap:
        _add_basemap(ax, info.crs, basemap_source, basemap_zoom)

    return fig, ax  # pyright: ignore[reportReturnType]

plot_floodmap

plot_floodmap

plot_floodmap(
    floodmap_path,
    *,
    ax=None,
    title=None,
    basemap=True,
    basemap_source=None,
    basemap_zoom=12,
    max_display_px=2000,
    vmax_percentile=98,
    figsize=(11, 7),
    color_map="viridis_r",
)

Plot a flood-depth COG with an optional satellite basemap.

Reads at an overview level that keeps the longest axis under max_display_px pixels, masks dry / NaN pixels, and renders with a reverse viridis color map.

PARAMETER	DESCRIPTION
floodmap_path	Path to the flood-depth GeoTIFF (e.g. floodmap_hmax.tif). TYPE: Path or str
ax	Existing axes to plot into. A new figure is created when None. TYPE: Axes DEFAULT: None
title	Plot title. Defaults to "Flood depth (hmax)". TYPE: str DEFAULT: None
basemap	If True (default), overlay satellite imagery via contextily. TYPE: bool DEFAULT: True
basemap_source	Tile provider passed to contextily.add_basemap. TYPE: optional DEFAULT: None
basemap_zoom	Zoom level for the basemap tiles, by default 12. TYPE: int DEFAULT: 12
max_display_px	Target maximum dimension (in pixels) for the rendered raster. Controls which overview level is read, by default 2000. TYPE: int DEFAULT: 2000
vmax_percentile	Upper percentile for the color-map range. TYPE: float DEFAULT: 98
figsize	Figure size when ax is None, by default (11, 7). TYPE: tuple DEFAULT: (11, 7)
color_map	Name of the Matplotlib colormap to use for plotting the flood depth, by default "viridis_r". TYPE: str DEFAULT: 'viridis_r'

RETURNS	DESCRIPTION
(Figure, Axes)

Source code in src/coastal_calibration/sfincs/plotting.py

def plot_floodmap(
    floodmap_path: Path | str,
    *,
    ax: Axes | None = None,
    title: str | None = None,
    basemap: bool = True,
    basemap_source: Any | None = None,
    basemap_zoom: int = 12,
    max_display_px: int = 2000,
    vmax_percentile: float = 98,
    figsize: tuple[float, float] = (11, 7),
    color_map: str = "viridis_r",
) -> tuple[Figure | SubFigure, Axes]:
    """Plot a flood-depth COG with an optional satellite basemap.

    Reads at an overview level that keeps the longest axis under
    *max_display_px* pixels, masks dry / NaN pixels, and renders
    with a reverse viridis color map.

    Parameters
    ----------
    floodmap_path : Path or str
        Path to the flood-depth GeoTIFF (e.g. ``floodmap_hmax.tif``).
    ax : Axes, optional
        Existing axes to plot into. A new figure is created when *None*.
    title : str, optional
        Plot title. Defaults to ``"Flood depth (hmax)"``.
    basemap : bool, optional
        If *True* (default), overlay satellite imagery via *contextily*.
    basemap_source : optional
        Tile provider passed to ``contextily.add_basemap``.
    basemap_zoom : int, optional
        Zoom level for the basemap tiles, by default 12.
    max_display_px : int, optional
        Target maximum dimension (in pixels) for the rendered raster.
        Controls which overview level is read, by default 2000.
    vmax_percentile : float, optional
        Upper percentile for the color-map range.
    figsize : tuple, optional
        Figure size when *ax* is *None*, by default (11, 7).
    color_map : str, optional
        Name of the Matplotlib colormap to use for plotting the flood depth,
        by default "viridis_r".

    Returns
    -------
    (Figure, Axes)
    """
    import matplotlib.pyplot as plt
    import rasterio

    floodmap_path = Path(floodmap_path)
    if not floodmap_path.exists():
        raise FileNotFoundError(
            f"Flood map not found: {floodmap_path}; "
            "ensure floodmap_dem is set and sfincs_map.nc contains zsmax."
        )

    if color_map not in plt.colormaps:
        raise ValueError(
            f"Invalid color_map: {color_map}. Must be a valid Matplotlib colormap name."
        )

    # ── Read metadata at full resolution ─────────────────────────
    with rasterio.open(floodmap_path) as src:
        bounds = src.bounds
        raster_crs = src.crs

        overviews = src.overviews(1)
        ovr_idx = next(
            (
                i
                for i, f in enumerate(overviews)
                if max(src.height, src.width) / f <= max_display_px
            ),
            len(overviews) - 1,
        )

    # ── Read at overview level ───────────────────────────────────
    with rasterio.open(floodmap_path, overview_level=ovr_idx) as src:
        hmax = src.read(1)

    hmax_masked = np.ma.masked_where(~np.isfinite(hmax) | (hmax <= 0), hmax)
    if hmax_masked.count() == 0:
        raise ValueError("No valid flood depth values found in the raster.")

    if ax is None:
        fig, ax = plt.subplots(figsize=figsize)
    else:
        fig = ax.get_figure()
        if fig is None:  # pragma: no cover
            msg = "ax must be attached to a Figure"
            raise ValueError(msg)

    extent = (bounds.left, bounds.right, bounds.bottom, bounds.top)
    cmap = plt.colormaps[color_map].copy()
    cmap.set_bad(alpha=0)

    im = ax.imshow(
        hmax_masked,
        extent=extent,
        origin="upper",
        cmap=cmap,
        vmin=0,
        vmax=np.percentile(hmax_masked.compressed(), vmax_percentile),
        interpolation="nearest",
        zorder=2,
    )
    fig.colorbar(im, ax=ax, label="Flood depth (m)", shrink=0.6, pad=0.02, extend="both")

    if title is None:
        title = "Flood depth (hmax)"
    ax.set_title(title)

    if basemap:
        import contextily as cx

        if basemap_source is None:
            basemap_source = cx.providers.Esri.WorldImagery  # pyright: ignore[reportAttributeAccessIssue]
        cx.add_basemap(ax, crs=raster_crs, source=basemap_source, zoom=basemap_zoom)  # pyright: ignore[reportArgumentType]

    return fig, ax

plot_station_comparison

plot_station_comparison

plot_station_comparison(
    runs, station_ids, figs_dir, *, obs_ds=None
)

Create station comparison figures for one or more simulated runs.

PARAMETER	DESCRIPTION
runs	Maps each run's label to a (times, elevation) pair. Labels become the legend entries. Each run's elevation must have the same n_stations columns aligned with station_ids; per-run time axes may differ. TYPE: Mapping[str, tuple[ArrayLike, NDArray]]
station_ids	NOAA station IDs, one per column of every run's elevation array. TYPE: list[str]
figs_dir	Output directory for figures (created if needed). TYPE: Path or str
obs_ds	Observed water levels with a water_level variable indexed by station and time. When None, the plots are pure model-vs-model comparisons. TYPE: Dataset DEFAULT: None

RETURNS	DESCRIPTION
list[Path]	Paths to the saved PNG figures, in pagination order (4 stations per 2x2 figure).

Notes

A station is plotable if any run or the observations has at least one finite value at that station. Stations with no data anywhere are silently skipped.

Runs are drawn on top of observations. Colors are drawn from the tab10 colormap (cycling at 10 runs); markers cycle through a 10-shape palette independently.

Source code in src/coastal_calibration/plotting/stations.py

def plot_station_comparison(
    runs: Mapping[str, tuple[ArrayLike, NDArray[np.float64]]],
    station_ids: list[str],
    figs_dir: Path | str,
    *,
    obs_ds: Any | None = None,
) -> list[Path]:
    """Create station comparison figures for one or more simulated runs.

    Parameters
    ----------
    runs : Mapping[str, tuple[ArrayLike, NDArray]]
        Maps each run's label to a ``(times, elevation)`` pair. Labels
        become the legend entries. Each run's *elevation* must have the
        same ``n_stations`` columns aligned with *station_ids*; per-run
        time axes may differ.
    station_ids : list[str]
        NOAA station IDs, one per column of every run's elevation array.
    figs_dir : Path or str
        Output directory for figures (created if needed).
    obs_ds : xr.Dataset, optional
        Observed water levels with a ``water_level`` variable indexed by
        ``station`` and ``time``. When *None*, the plots are pure
        model-vs-model comparisons.

    Returns
    -------
    list[pathlib.Path]
        Paths to the saved PNG figures, in pagination order (4 stations
        per 2x2 figure).

    Notes
    -----
    A station is plotable if any run or the observations has at least one
    finite value at that station. Stations with no data anywhere are
    silently skipped.

    Runs are drawn on top of observations. Colors are drawn from the
    ``tab10`` colormap (cycling at 10 runs); markers cycle through a
    10-shape palette independently.
    """
    if len(runs) == 0:
        msg = "plot_station_comparison: `runs` is empty."
        raise ValueError(msg)

    runs_dict = dict(runs)
    n_stations = len(station_ids)
    bad_shapes: list[str] = []
    for label, (_, elev) in runs_dict.items():
        arr = np.asarray(elev)
        if arr.ndim != 2 or arr.shape[1] != n_stations:
            bad_shapes.append(f"{label!r}: shape {arr.shape}")
    if bad_shapes:
        msg = (
            f"plot_station_comparison: each run's elevation must be 2-D with "
            f"{n_stations} columns (one per station in `station_ids`). "
            f"Offending runs: {', '.join(bad_shapes)}"
        )
        raise ValueError(msg)

    import sys

    import matplotlib

    if "ipykernel" not in sys.modules:
        matplotlib.use("Agg")

    import matplotlib.pyplot as plt

    figs_dir = Path(figs_dir)
    figs_dir.mkdir(parents=True, exist_ok=True)

    stations = _plotable_stations(station_ids, runs_dict, obs_ds)
    if not stations:
        return []

    palette = plt.colormaps["tab10"]
    labels = list(runs_dict.keys())
    run_colors = {label: palette(i % 10) for i, label in enumerate(labels)}
    run_markers = {label: _RUN_MARKERS[i % len(_RUN_MARKERS)] for i, label in enumerate(labels)}

    n_plotable = len(stations)
    n_figures = math.ceil(n_plotable / _STATIONS_PER_FIGURE)

    saved: list[Path] = []
    for fig_idx in range(n_figures):
        start = fig_idx * _STATIONS_PER_FIGURE
        end = min(start + _STATIONS_PER_FIGURE, n_plotable)
        batch = stations[start:end]
        batch_size = len(batch)

        nrows = 2 if batch_size > 2 else 1
        ncols = 2 if batch_size > 1 else 1

        fig, axes = plt.subplots(
            nrows,
            ncols,
            figsize=(16, 5 * nrows),
            squeeze=False,
        )
        axes_flat = axes.ravel()

        for i, (sid, col_idx) in enumerate(batch):
            _draw_station_panel(
                axes_flat[i], sid, col_idx, runs_dict, obs_ds, run_colors, run_markers
            )

        for j in range(batch_size, nrows * ncols):
            axes_flat[j].remove()

        fig.tight_layout()
        fig_path = figs_dir / f"stations_comparison_{fig_idx + 1:03d}.png"
        fig.savefig(fig_path, dpi=300, bbox_inches="tight")
        plt.close(fig)
        saved.append(fig_path)

    return saved

plot_water_level

plot_water_level

plot_water_level(
    ds,
    time=0,
    *,
    variable=None,
    ax=None,
    cmap="viridis",
    vmin=None,
    vmax=None,
    colorbar=True,
    title=None,
    figsize=(10, 7),
    shading_regular="auto",
    mask_dry=True,
    dry_threshold=0.05,
    basemap=False,
    basemap_source=None,
    basemap_zoom=None,
    crs=None,
)

Render a single water-level (or water-depth) frame from ds.

PARAMETER	DESCRIPTION
ds	Dataset produced by one of the load_* readers in :mod:coastal_calibration.schism.outputs / :mod:coastal_calibration.sfincs.outputs. Must carry a mesh_type attribute of "regular", "ugrid-triangle-or-quad", or "ugrid-quadtree". TYPE: Dataset
time	Time selector. Integers are positional (isel); datetime-like values (numpy.datetime64, pandas.Timestamp, ISO strings) are passed to sel(..., method="nearest"). TYPE: int or datetime - like DEFAULT: ``0``
variable	Data variable to plot. Common choices: "zs" / "elevation" — water-surface elevation (the default, auto-detected when omitted) "h" — water depth TYPE: str DEFAULT: None
ax	Axes to draw into. A new figure + axes are created when None. TYPE: Axes DEFAULT: None
cmap	Matplotlib colormap name. TYPE: str DEFAULT: ``"viridis"``
vmin	Colormap limits. When None, both are computed from the 1st/99th percentiles of the full time series after dry-cell masking so outlier bed-elevation cells do not dominate the scale. TYPE: float DEFAULT: None
vmax	Colormap limits. When None, both are computed from the 1st/99th percentiles of the full time series after dry-cell masking so outlier bed-elevation cells do not dominate the scale. TYPE: float DEFAULT: None
colorbar	If True, attach a colorbar to the axes. TYPE: bool DEFAULT: ``True``
title	Plot title. Defaults to "<variable> @ <time>" when the dataset has a datetime time axis. TYPE: str DEFAULT: None
figsize	Figure size when ax is None. TYPE: (float, float) DEFAULT: ``(10, 7)``
shading_regular	shading passed through to pcolormesh on the regular-grid path. Ignored for unstructured meshes. TYPE: str DEFAULT: ``"auto"``
mask_dry	Mask cells classified as dry. The mask is sourced (in priority order): dryFlagNode == 0 if the SCHISM dry flag is in the dataset. h > dry_threshold otherwise (water depth above threshold). When neither variable is present the dataset is plotted unmasked. TYPE: bool DEFAULT: ``True``
dry_threshold	Water-depth threshold (m) used by the fallback mask. 5 cm is a common SFINCS/HydroMT convention for "wet enough to plot". TYPE: float DEFAULT: ``0.05``
basemap	If True, overlay a satellite basemap via :mod:contextily. The data CRS is taken from crs, falling back to ds.attrs["crs"]. Raises if no CRS is available. TYPE: bool DEFAULT: ``False``
basemap_source	Tile provider passed to :func:contextily.add_basemap. Defaults to cx.providers.Esri.WorldImagery (satellite imagery). TYPE: optional DEFAULT: None
basemap_zoom	Zoom level for the basemap tiles; None lets contextily pick. TYPE: int DEFAULT: None
crs	Override the dataset CRS for basemap reprojection (e.g. "EPSG:32614"). When None, ds.attrs["crs"] is used. TYPE: str DEFAULT: None

RETURNS	DESCRIPTION
(ax, collection)	The axes used for drawing and the primitive collection (QuadMesh for regular grids, TriMesh / PolyCollection for unstructured). Returning the collection allows :func:coastal_calibration.plotting.animate.animate_water_level to call collection.set_array(...) per frame.

RAISES	DESCRIPTION
KeyError	If ds lacks the mesh_type attribute or variable is not a data variable.
ValueError	If no water-level variable can be auto-detected or the mesh type is unknown.

Notes

When ax is None a new figure is created and the caller owns it — remember to plt.close(ax.get_figure()) when looping over many frames outside a notebook, otherwise matplotlib accumulates figures and memory grows unbounded.

Source code in src/coastal_calibration/plotting/spatial.py

def plot_water_level(
    ds: xr.Dataset,
    time: int | Any = 0,
    *,
    variable: str | None = None,
    ax: Axes | None = None,
    cmap: str = "viridis",
    vmin: float | None = None,
    vmax: float | None = None,
    colorbar: bool = True,
    title: str | None = None,
    figsize: tuple[float, float] = (10, 7),
    shading_regular: PColorShading = "auto",
    mask_dry: bool = True,
    dry_threshold: float = 0.05,
    basemap: bool = False,
    basemap_source: Any | None = None,
    basemap_zoom: int | None = None,
    crs: str | None = None,
) -> tuple[Axes, Collection]:
    """Render a single water-level (or water-depth) frame from *ds*.

    Parameters
    ----------
    ds : xarray.Dataset
        Dataset produced by one of the ``load_*`` readers in
        :mod:`coastal_calibration.schism.outputs` /
        :mod:`coastal_calibration.sfincs.outputs`. Must carry a
        ``mesh_type`` attribute of ``"regular"``,
        ``"ugrid-triangle-or-quad"``, or ``"ugrid-quadtree"``.
    time : int or datetime-like, default ``0``
        Time selector. Integers are positional (``isel``); datetime-like
        values (``numpy.datetime64``, ``pandas.Timestamp``, ISO strings)
        are passed to ``sel(..., method="nearest")``.
    variable : str, optional
        Data variable to plot. Common choices:

        - ``"zs"`` / ``"elevation"`` — water-surface elevation (the default,
          auto-detected when omitted)
        - ``"h"`` — water depth
    ax : matplotlib.axes.Axes, optional
        Axes to draw into. A new figure + axes are created when *None*.
    cmap : str, default ``"viridis"``
        Matplotlib colormap name.
    vmin, vmax : float, optional
        Colormap limits. When *None*, both are computed from the 1st/99th
        percentiles of the full time series **after dry-cell masking** so
        outlier bed-elevation cells do not dominate the scale.
    colorbar : bool, default ``True``
        If *True*, attach a colorbar to the axes.
    title : str, optional
        Plot title. Defaults to ``"<variable> @ <time>"`` when the dataset
        has a datetime time axis.
    figsize : (float, float), default ``(10, 7)``
        Figure size when *ax* is *None*.
    shading_regular : str, default ``"auto"``
        ``shading`` passed through to ``pcolormesh`` on the regular-grid
        path. Ignored for unstructured meshes.
    mask_dry : bool, default ``True``
        Mask cells classified as dry. The mask is sourced (in priority order):

        1. ``dryFlagNode == 0`` if the SCHISM dry flag is in the dataset.
        2. ``h > dry_threshold`` otherwise (water depth above threshold).

        When neither variable is present the dataset is plotted unmasked.
    dry_threshold : float, default ``0.05``
        Water-depth threshold (m) used by the fallback mask. 5 cm is a
        common SFINCS/HydroMT convention for "wet enough to plot".
    basemap : bool, default ``False``
        If *True*, overlay a satellite basemap via :mod:`contextily`. The
        data CRS is taken from *crs*, falling back to ``ds.attrs["crs"]``.
        Raises if no CRS is available.
    basemap_source : optional
        Tile provider passed to :func:`contextily.add_basemap`. Defaults
        to ``cx.providers.Esri.WorldImagery`` (satellite imagery).
    basemap_zoom : int, optional
        Zoom level for the basemap tiles; *None* lets contextily pick.
    crs : str, optional
        Override the dataset CRS for basemap reprojection (e.g.
        ``"EPSG:32614"``). When *None*, ``ds.attrs["crs"]`` is used.

    Returns
    -------
    (ax, collection)
        The axes used for drawing and the primitive collection
        (``QuadMesh`` for regular grids, ``TriMesh`` / ``PolyCollection``
        for unstructured). Returning the collection allows
        :func:`coastal_calibration.plotting.animate.animate_water_level`
        to call ``collection.set_array(...)`` per frame.

    Raises
    ------
    KeyError
        If *ds* lacks the ``mesh_type`` attribute or *variable* is not a
        data variable.
    ValueError
        If no water-level variable can be auto-detected or the mesh type
        is unknown.

    Notes
    -----
    When *ax* is *None* a new figure is created and the caller owns it —
    remember to ``plt.close(ax.get_figure())`` when looping over many
    frames outside a notebook, otherwise matplotlib accumulates figures
    and memory grows unbounded.
    """
    import matplotlib.pyplot as plt

    mesh_type = ds.attrs.get("mesh_type")
    if mesh_type is None:
        msg = (
            "Dataset is missing the 'mesh_type' attribute. "
            "Use coastal_calibration.{schism,sfincs}.outputs loaders to produce it."
        )
        raise KeyError(msg)

    if variable is None:
        variable = _auto_variable(ds)
    if variable not in ds.data_vars:
        available = sorted(str(name) for name in ds.data_vars)
        msg = f"Variable {variable!r} not in dataset (have: {available})"
        raise KeyError(msg)

    if mask_dry:
        ds = _apply_wet_mask(ds, variable, dry_threshold)

    vmin, vmax = _resolve_limits(ds, variable, vmin, vmax)

    if ax is None:
        _, ax = plt.subplots(figsize=figsize)

    coll = _dispatch_plot(
        ds,
        variable,
        time,
        ax,
        mesh_type=mesh_type,
        cmap=cmap,
        vmin=vmin,
        vmax=vmax,
        shading_regular=shading_regular,
    )

    # Equal aspect keeps lon/lat or projected meters looking correct.
    # ``adjustable="box"`` resizes the axes in figure-space instead of
    # clipping the data limits — important when a contextily basemap
    # comes in afterward, since contextily sets its own xlim/ylim and
    # ``adjustable="datalim"`` would otherwise emit a warning about
    # ignoring them.
    ax.set_aspect("equal", adjustable="box")

    ax.set_title(title if title is not None else _default_title(ds, variable, time))
    if mesh_type == _REGULAR:
        ax.set_xlabel("x")
        ax.set_ylabel("y")
    else:
        ax.set_xlabel("longitude")
        ax.set_ylabel("latitude")

    if colorbar:
        fig = ax.get_figure()
        if fig is not None:
            units = ds[variable].attrs.get("units", "m")
            extend = _colorbar_extend(np.asarray(ds[variable].to_numpy()), vmin, vmax)
            fig.colorbar(
                coll,
                ax=ax,
                label=f"{variable} ({units})",
                shrink=0.8,
                extend=extend,
            )

    if basemap:
        crs_label = crs if crs is not None else ds.attrs.get("crs")
        if not crs_label:
            msg = (
                "basemap=True but no CRS available. Pass `crs='EPSG:...'` or "
                "ensure the loader populated ds.attrs['crs']."
            )
            raise ValueError(msg)
        _add_basemap(ax, crs_label, basemap_source, basemap_zoom)

    return ax, coll

animate_water_level

animate_water_level

animate_water_level(
    ds,
    outfile,
    *,
    variable=None,
    fps=10,
    time_stride=1,
    dpi=150,
    writer="auto",
    cmap="viridis",
    vmin=None,
    vmax=None,
    figsize=(10, 7),
    title_prefix=None,
    mask_dry=True,
    dry_threshold=0.05,
)

Render a time-animation of the water-level field to a movie file.

PARAMETER	DESCRIPTION
ds	Canonical dataset from a load_* reader in :mod:coastal_calibration.schism.outputs / :mod:coastal_calibration.sfincs.outputs. Must carry a mesh_type attribute that :func:plot_water_level recognizes. TYPE: Dataset
outfile	Destination path. The suffix selects the writer: .mp4 / .mov / .avi use FFMpegWriter (requires an ffmpeg binary on PATH); .gif uses the pure-Python PillowWriter. TYPE: str or Path
variable	Variable to animate. Defaults to the result of :func:_auto_variable. TYPE: str DEFAULT: None
fps	Frames per second in the output. TYPE: int DEFAULT: 10
time_stride	Keep every time_stride-th frame from the dataset's time axis. Useful for previews or reducing file size on long runs. TYPE: int DEFAULT: 1
dpi	Output resolution. TYPE: int DEFAULT: 150
writer	Writer selector; "auto" infers from the output suffix. TYPE: ('auto', 'ffmpeg', 'pillow') DEFAULT: "auto"
cmap	Matplotlib colormap name. TYPE: str DEFAULT: ``"viridis"``
vmin	Colormap limits — shared across all frames. When None, each is filled in from the 1st/99th percentile of the full time series (see :func:plot_water_level). TYPE: float DEFAULT: None
vmax	Colormap limits — shared across all frames. When None, each is filled in from the 1st/99th percentile of the full time series (see :func:plot_water_level). TYPE: float DEFAULT: None
figsize	Figure size. TYPE: tuple[float, float] DEFAULT: ``(10, 7)``
title_prefix	Prefix prepended to the auto-generated per-frame title. TYPE: str DEFAULT: None
mask_dry	Mask dry cells in every frame using dryFlagNode == 0 (preferred) or h > dry_threshold as a fallback. See :func:plot_water_level for details. TYPE: bool DEFAULT: ``True``
dry_threshold	Water-depth threshold (m) for the fallback mask. TYPE: float DEFAULT: ``0.05``

RETURNS	DESCRIPTION
Path	The resolved output path.

RAISES	DESCRIPTION
RuntimeError	If an .mp4/.mov/.avi is requested but ffmpeg is not on PATH.
ValueError	If the output suffix is unrecognized and no explicit writer was given.

Notes

The renderer is built from :func:plot_water_level, so any future additions to the frame layout (basemaps, projections, annotations) automatically flow through to animations.

Source code in src/coastal_calibration/plotting/animate.py

def animate_water_level(
    ds: xr.Dataset,
    outfile: str | Path,
    *,
    variable: str | None = None,
    fps: int = 10,
    time_stride: int = 1,
    dpi: int = 150,
    writer: Literal["auto", "ffmpeg", "pillow"] = "auto",
    cmap: str = "viridis",
    vmin: float | None = None,
    vmax: float | None = None,
    figsize: tuple[float, float] = (10, 7),
    title_prefix: str | None = None,
    mask_dry: bool = True,
    dry_threshold: float = 0.05,
) -> Path:
    """Render a time-animation of the water-level field to a movie file.

    Parameters
    ----------
    ds : xarray.Dataset
        Canonical dataset from a ``load_*`` reader in
        :mod:`coastal_calibration.schism.outputs` /
        :mod:`coastal_calibration.sfincs.outputs`. Must carry a
        ``mesh_type`` attribute that :func:`plot_water_level` recognizes.
    outfile : str or pathlib.Path
        Destination path. The suffix selects the writer:
        ``.mp4`` / ``.mov`` / ``.avi`` use ``FFMpegWriter`` (requires an
        ``ffmpeg`` binary on PATH); ``.gif`` uses the pure-Python
        ``PillowWriter``.
    variable : str, optional
        Variable to animate. Defaults to the result of :func:`_auto_variable`.
    fps : int, default 10
        Frames per second in the output.
    time_stride : int, default 1
        Keep every ``time_stride``-th frame from the dataset's time axis.
        Useful for previews or reducing file size on long runs.
    dpi : int, default 150
        Output resolution.
    writer : {"auto", "ffmpeg", "pillow"}, default ``"auto"``
        Writer selector; ``"auto"`` infers from the output suffix.
    cmap : str, default ``"viridis"``
        Matplotlib colormap name.
    vmin, vmax : float, optional
        Colormap limits — shared across all frames. When *None*, each is
        filled in from the 1st/99th percentile of the full time series
        (see :func:`plot_water_level`).
    figsize : tuple[float, float], default ``(10, 7)``
        Figure size.
    title_prefix : str, optional
        Prefix prepended to the auto-generated per-frame title.
    mask_dry : bool, default ``True``
        Mask dry cells in every frame using ``dryFlagNode == 0`` (preferred)
        or ``h > dry_threshold`` as a fallback. See
        :func:`plot_water_level` for details.
    dry_threshold : float, default ``0.05``
        Water-depth threshold (m) for the fallback mask.

    Returns
    -------
    pathlib.Path
        The resolved output path.

    Raises
    ------
    RuntimeError
        If an ``.mp4``/``.mov``/``.avi`` is requested but ``ffmpeg`` is not
        on PATH.
    ValueError
        If the output suffix is unrecognized and no explicit *writer* was
        given.

    Notes
    -----
    The renderer is built from :func:`plot_water_level`, so any future
    additions to the frame layout (basemaps, projections, annotations)
    automatically flow through to animations.
    """
    import matplotlib.pyplot as plt
    from matplotlib.animation import FuncAnimation

    outfile = Path(outfile).expanduser().resolve()
    outfile.parent.mkdir(parents=True, exist_ok=True)

    # Validate the writer selection up front so a bad suffix raises before
    # we build the figure / animation.
    movie_writer = _pick_writer(outfile, fps=fps, writer=writer)

    if variable is None:
        variable = _auto_variable(ds)

    # Compute shared color limits on the *masked* dataset so dry-cell
    # outliers do not stretch the colormap. The original ds is untouched
    # by _apply_wet_mask (it returns a shallow copy).
    ds_for_limits = _apply_wet_mask(ds, variable, dry_threshold) if mask_dry else ds
    vmin_r, vmax_r = _resolve_limits(ds_for_limits, variable, vmin, vmax)

    stride = int(time_stride)
    if stride < 1:
        msg = f"time_stride must be >= 1 (got {time_stride})"
        raise ValueError(msg)
    time_indices = np.arange(0, ds.sizes["time"], stride, dtype=np.int64)
    if time_indices.size == 0:
        msg = "No frames to animate: time dimension is empty."
        raise ValueError(msg)

    fig, ax = plt.subplots(figsize=figsize)
    # Build the first frame; `plot_water_level` attaches a colorbar and
    # applies the same dry-cell mask we will reuse in every frame.
    _, coll = plot_water_level(
        ds,
        time=int(time_indices[0]),
        variable=variable,
        ax=ax,
        cmap=cmap,
        vmin=vmin_r,
        vmax=vmax_r,
        colorbar=True,
        title=_frame_title(ds, variable, int(time_indices[0]), title_prefix),
        mask_dry=mask_dry,
        dry_threshold=dry_threshold,
    )

    base_update = _make_updater(
        coll,
        ds,
        variable,
        time_indices,
        mask_dry=mask_dry,
        dry_threshold=dry_threshold,
    )

    def update(i: int) -> tuple[Any, ...]:
        artists = base_update(i)
        ax.set_title(_frame_title(ds, variable, int(time_indices[i]), title_prefix))
        return artists

    anim = FuncAnimation(
        fig,
        update,
        frames=len(time_indices),
        interval=max(1, int(1000 / fps)),
        blit=False,
    )
    anim.save(str(outfile), writer=movie_writer, dpi=dpi)
    plt.close(fig)

    _log.info("wrote %s (%d frames, %d fps)", outfile, len(time_indices), fps)
    return outfile

Output Readers

load_schism_elevation

load_schism_elevation

load_schism_elevation(
    run_dir, *, time_slice=None, correction_file=None
)

Load SCHISM 2-D elevation across all output blocks.

PARAMETER	DESCRIPTION
run_dir	SCHISM run directory. Either the project root (the reader will look under outputs/) or the outputs/ directory itself. TYPE: str or Path
time_slice	Slice applied along the concatenated time dimension before the dataset is returned. Useful for previews of long runs. TYPE: slice DEFAULT: None
correction_file	Path to elevation_correction.csv. When provided, the conversion_factor is interpolated from the boundary nodes listed in the CSV onto every mesh node (nearest-neighbour) and added to the elevations to convert from the model's mesh datum back to MSL — the inverse of the boundary correction applied by correct_elevation. This preserves the spatial variation of the geoid-MSL offset across the domain instead of applying a single boundary mean (which under-corrects nodes whose local offset deviates from the average). Use this when the model was forced with STOFS via :func:make_stofs_boundary so that simulated values are MSL-referenced (matching NOAA CO-OPS observations). TYPE: str or Path DEFAULT: None

RETURNS	DESCRIPTION
Dataset	Dataset with variables elevation(time, node), node_x(node), node_y(node), and face_nodes(face, face_node). The time coordinate is an absolute datetime64[ns] axis reconstructed from the SCHISM base_date attribute. Attributes include mesh_type = "ugrid-triangle-or-quad" and source_run.

RAISES	DESCRIPTION
NotADirectoryError	If run_dir is not a directory.
FileNotFoundError	If no out2d_*.nc files are found.
KeyError	If any required mesh variable is missing from the first block.

Notes

The reader does not require dask. Blocks are concatenated manually via :func:xarray.concat after opening each file with the default engine. Mesh variables (node coordinates and face connectivity) are read from the first block only — SCHISM keeps these constant across a run.

Source code in src/coastal_calibration/schism/outputs.py

def load_schism_elevation(
    run_dir: str | Path,
    *,
    time_slice: slice | None = None,
    correction_file: str | Path | None = None,
) -> xr.Dataset:
    """Load SCHISM 2-D elevation across all output blocks.

    Parameters
    ----------
    run_dir : str or pathlib.Path
        SCHISM run directory. Either the project root (the reader will look
        under ``outputs/``) or the ``outputs/`` directory itself.
    time_slice : slice, optional
        Slice applied along the concatenated time dimension before the
        dataset is returned. Useful for previews of long runs.
    correction_file : str or pathlib.Path, optional
        Path to ``elevation_correction.csv``.  When provided, the
        ``conversion_factor`` is interpolated from the boundary nodes
        listed in the CSV onto every mesh node (nearest-neighbour) and
        added to the elevations to convert from the model's mesh datum
        back to MSL — the inverse of the boundary correction applied by
        ``correct_elevation``.  This preserves the spatial variation of
        the geoid-MSL offset across the domain instead of applying a
        single boundary mean (which under-corrects nodes whose local
        offset deviates from the average).  Use this when the model was
        forced with STOFS via :func:`make_stofs_boundary` so that
        simulated values are MSL-referenced (matching NOAA CO-OPS
        observations).

    Returns
    -------
    xarray.Dataset
        Dataset with variables ``elevation(time, node)``, ``node_x(node)``,
        ``node_y(node)``, and ``face_nodes(face, face_node)``. The time
        coordinate is an absolute ``datetime64[ns]`` axis reconstructed from
        the SCHISM ``base_date`` attribute. Attributes include
        ``mesh_type = "ugrid-triangle-or-quad"`` and ``source_run``.

    Raises
    ------
    NotADirectoryError
        If *run_dir* is not a directory.
    FileNotFoundError
        If no ``out2d_*.nc`` files are found.
    KeyError
        If any required mesh variable is missing from the first block.

    Notes
    -----
    The reader does not require ``dask``. Blocks are concatenated manually
    via :func:`xarray.concat` after opening each file with the default engine.
    Mesh variables (node coordinates and face connectivity) are read from
    the first block only — SCHISM keeps these constant across a run.
    """
    run_dir = Path(run_dir)
    outputs_dir = _resolve_outputs_dir(run_dir)

    files = sorted(outputs_dir.glob(_OUT2D_GLOB))
    if not files:
        # Unreachable because _resolve_outputs_dir checks, but kept for clarity.
        raise FileNotFoundError(f"No {_OUT2D_GLOB} under {outputs_dir}")

    # Open the first block to discover which optional time-varying variables
    # are actually present (e.g. dryFlagNode), then open every block and pull
    # only what we need.
    with xr.open_dataset(files[0], decode_times=False) as ds0_probe:
        missing = [v for v in _REQUIRED_DATA_VARS if v not in ds0_probe.variables]
        if missing:
            raise KeyError(
                f"SCHISM output file {files[0].name} is missing required "
                f"variables: {', '.join(missing)}"
            )
        present_optional = tuple(v for v in _OPTIONAL_TIME_VARS if v in ds0_probe.variables)

    time_vars = ("elevation", "time", *present_optional)
    per_block: list[xr.Dataset] = []
    for path in files:
        with xr.open_dataset(path, decode_times=False) as ds:
            per_block.append(ds[list(time_vars)].load())

    raw = xr.concat(per_block, dim="time", data_vars="minimal").sortby("time")

    node_x, node_y, face_nodes, depth, base_date, crs_label = _load_mesh_geometry(files[0])

    # SCHISM writes seconds-since-base_date as the raw time values.
    seconds = raw["time"].to_numpy().astype("float64")
    times = base_date + pd.to_timedelta(seconds, unit="s")

    elevation = raw["elevation"].to_numpy()

    # Optionally convert from mesh datum back to MSL by applying the
    # inverse of the boundary forcing correction. See
    # :func:`_apply_elevation_correction` for the per-node lookup.
    if correction_file is not None:
        corr_path = Path(correction_file)
        if corr_path.exists():
            elevation = _apply_elevation_correction(elevation, node_x, node_y, corr_path)

    # Water depth = water surface + bathymetric depth (depth is positive when
    # the bed is below the datum). For dry nodes, SCHISM reports
    # ``elevation = -depth + h_min``, so h ≈ 0 there.
    h = elevation + depth[np.newaxis, :]

    data_vars: dict[str, tuple[tuple[str, ...], NDArray[Any]]] = {
        "elevation": (("time", "node"), elevation),
        "h": (("time", "node"), h),
        "depth": (("node",), depth),
        "node_x": (("node",), node_x),
        "node_y": (("node",), node_y),
        "face_nodes": (("face", "face_node"), face_nodes),
    }
    if "dryFlagNode" in present_optional:
        data_vars["dryFlagNode"] = (
            ("time", "node"),
            raw["dryFlagNode"].to_numpy().astype(np.int8, copy=False),
        )

    out = xr.Dataset(
        data_vars=data_vars,
        coords={
            "time": times,
            "node": np.arange(node_x.size, dtype=np.int64),
            "face": np.arange(face_nodes.shape[0], dtype=np.int64),
        },
        attrs={
            "mesh_type": "ugrid-triangle-or-quad",
            "source_run": str(run_dir),
        },
    )
    if crs_label is not None:
        out.attrs["crs"] = crs_label

    _set_dataset_attrs(out, has_dry_flag="dryFlagNode" in out.data_vars)

    if time_slice is not None:
        out = out.isel(time=time_slice)
    return out

load_sfincs_water_level

load_sfincs_water_level

load_sfincs_water_level(run_dir, *, time_slice=None)

Load SFINCS time-dependent water-level + depth from sfincs_map.nc.

PARAMETER	DESCRIPTION
run_dir	SFINCS run directory containing sfincs_map.nc, or the path to the file itself. TYPE: str or Path
time_slice	Slice applied along the time dimension before the dataset is returned. Useful for previews of long runs. TYPE: slice DEFAULT: None

RETURNS	DESCRIPTION
Dataset	Canonical dataset containing both zs(time, ...) (water-surface elevation) and h(time, ...) (water depth, derived as zs - zb), the static bed elevation zb(...), and (when present) the static cell mask msk(...). The mesh_type attribute is "regular" or "ugrid-quadtree". Any 2-D geographic coordinates xc(y, x) / yc(y, x) found in the source file are preserved as auxiliary coordinates.

RAISES	DESCRIPTION
NotADirectoryError	If run_dir is neither a directory nor a sfincs_map.nc file.
FileNotFoundError	If the map output file cannot be located, or if a structured (n, m) map has no sfincs.nc alongside it to supply the mesh topology.
KeyError	If a required variable (zs or zb) is missing, or if sfincs.nc lacks the n / m / mesh geometry variables.
ValueError	If zs does not have the expected dims for the detected layout, or if sfincs.nc indexes cells outside the map grid.

Notes

The reader intentionally opens sfincs_map.nc directly with xarray rather than going through :class:hydromt_sfincs.SfincsModel. This keeps post-processing independent of a full HydroMT model setup and works in lightweight environments where only xarray + netcdf4 are available. A HydroMT-based path can be added in a follow-up if needed for quadtree outputs.

Source code in src/coastal_calibration/sfincs/outputs.py

def load_sfincs_water_level(
    run_dir: str | Path,
    *,
    time_slice: slice | None = None,
) -> xr.Dataset:
    """Load SFINCS time-dependent water-level + depth from ``sfincs_map.nc``.

    Parameters
    ----------
    run_dir : str or pathlib.Path
        SFINCS run directory containing ``sfincs_map.nc``, or the path to
        the file itself.
    time_slice : slice, optional
        Slice applied along the time dimension before the dataset is
        returned. Useful for previews of long runs.

    Returns
    -------
    xarray.Dataset
        Canonical dataset containing both ``zs(time, ...)`` (water-surface
        elevation) and ``h(time, ...)`` (water depth, derived as
        ``zs - zb``), the static bed elevation ``zb(...)``, and (when
        present) the static cell mask ``msk(...)``. The ``mesh_type``
        attribute is ``"regular"`` or ``"ugrid-quadtree"``.

        Any 2-D geographic coordinates ``xc(y, x)`` / ``yc(y, x)`` found in
        the source file are preserved as auxiliary coordinates.

    Raises
    ------
    NotADirectoryError
        If *run_dir* is neither a directory nor a ``sfincs_map.nc`` file.
    FileNotFoundError
        If the map output file cannot be located, or if a structured
        ``(n, m)`` map has no ``sfincs.nc`` alongside it to supply the
        mesh topology.
    KeyError
        If a required variable (``zs`` or ``zb``) is missing, or if
        ``sfincs.nc`` lacks the ``n`` / ``m`` / mesh geometry variables.
    ValueError
        If ``zs`` does not have the expected dims for the detected layout,
        or if ``sfincs.nc`` indexes cells outside the map grid.

    Notes
    -----
    The reader intentionally opens ``sfincs_map.nc`` directly with *xarray*
    rather than going through :class:`hydromt_sfincs.SfincsModel`. This
    keeps post-processing independent of a full HydroMT model setup and
    works in lightweight environments where only ``xarray`` + ``netcdf4``
    are available. A HydroMT-based path can be added in a follow-up if
    needed for quadtree outputs.
    """
    run_dir = Path(run_dir)
    map_file = _resolve_map_file(run_dir)

    with xr.open_dataset(map_file) as ds:
        layout = _detect_layout(ds)
        if layout == "structured-nm":
            out = _build_structured_nm(ds, map_file.parent / _GRID_FILENAME)
        elif layout == "ugrid-quadtree":
            out = _build_quadtree(ds)
        elif layout == "regular":
            out = _build_regular(ds)
        else:
            msg = (
                f"Unrecognized SFINCS map layout in {map_file.name}. "
                f"Dims: {sorted(str(d) for d in ds.dims)}"
            )
            raise ValueError(msg)

    out.attrs["source_run"] = str(run_dir)
    crs_label = _detect_crs(map_file)
    if crs_label is not None:
        out.attrs["crs"] = crs_label
    _annotate_attrs(out)

    if time_slice is not None:
        out = out.isel(time=time_slice)
    return out

Observation Points

load_obs_points

load_obs_points

load_obs_points(path)

Read a user-supplied observation-points CSV.

PARAMETER	DESCRIPTION
path	Path to a CSV with columns id, lon, lat. TYPE: str or Path

RETURNS	DESCRIPTION
DataFrame	A frame with exactly those three columns. id is coerced to str; lon / lat to float64.

RAISES	DESCRIPTION
FileNotFoundError	If path does not exist.
ValueError	If any required column is missing, if IDs are not unique, or if any coordinate is non-finite.

Source code in src/coastal_calibration/observations.py

def load_obs_points(path: str | Path) -> pd.DataFrame:
    """Read a user-supplied observation-points CSV.

    Parameters
    ----------
    path : str or pathlib.Path
        Path to a CSV with columns ``id``, ``lon``, ``lat``.

    Returns
    -------
    pandas.DataFrame
        A frame with exactly those three columns. ``id`` is coerced to
        ``str``; ``lon`` / ``lat`` to ``float64``.

    Raises
    ------
    FileNotFoundError
        If *path* does not exist.
    ValueError
        If any required column is missing, if IDs are not unique, or if
        any coordinate is non-finite.
    """
    path = Path(path)
    if not path.is_file():
        raise FileNotFoundError(f"Observation points CSV not found: {path}")

    df = pd.read_csv(path)
    missing = [c for c in _REQUIRED_COLS if c not in df.columns]
    if missing:
        msg = (
            f"{path.name} is missing required column(s): {', '.join(missing)}. "
            f"Expected columns: {_REQUIRED_COLS}."
        )
        raise ValueError(msg)

    df = df[list(_REQUIRED_COLS)].copy()
    df["id"] = df["id"].astype(str)
    df["lon"] = pd.to_numeric(df["lon"], errors="coerce").astype("float64")
    df["lat"] = pd.to_numeric(df["lat"], errors="coerce").astype("float64")

    if df["lon"].isna().any() or df["lat"].isna().any():
        bad = df[df["lon"].isna() | df["lat"].isna()]["id"].tolist()
        msg = f"Non-numeric lon/lat in {path.name} for id(s): {bad}"
        raise ValueError(msg)

    if df["id"].duplicated().any():
        dupes = df.loc[df["id"].duplicated(keep=False), "id"].unique().tolist()
        msg = f"Duplicate ids in {path.name}: {sorted(dupes)}"
        raise ValueError(msg)

    return df.reset_index(drop=True)

validate_points_in_domain

validate_points_in_domain

validate_points_in_domain(points, ds)

Raise if any point in points falls outside the model WGS84 bbox.

PARAMETER	DESCRIPTION
points	Must have lon and lat columns in WGS84. TYPE: DataFrame
ds	Canonical dataset from a load_* reader. TYPE: Dataset

RETURNS	DESCRIPTION
Polygon	The WGS84 bounding polygon used for the check. Callers typically discard this; it is returned for inspection / logging.

RAISES	DESCRIPTION
ValueError	If one or more points are outside the domain. The error message lists offending IDs (up to 10).

Source code in src/coastal_calibration/observations.py

def validate_points_in_domain(points: pd.DataFrame, ds: xr.Dataset) -> shapely.Polygon:
    """Raise if any point in *points* falls outside the model WGS84 bbox.

    Parameters
    ----------
    points : pandas.DataFrame
        Must have ``lon`` and ``lat`` columns in WGS84.
    ds : xarray.Dataset
        Canonical dataset from a ``load_*`` reader.

    Returns
    -------
    shapely.Polygon
        The WGS84 bounding polygon used for the check. Callers typically
        discard this; it is returned for inspection / logging.

    Raises
    ------
    ValueError
        If one or more points are outside the domain. The error message
        lists offending IDs (up to 10).
    """
    import shapely

    bbox = _domain_bbox_wgs84(ds)
    pts = shapely.points(points["lon"].to_numpy(), points["lat"].to_numpy())
    # Use `covers` rather than `contains` so points exactly on the bbox
    # edge pass — `contains` excludes the boundary.
    inside = shapely.covers(bbox, pts)
    if not np.asarray(inside).all():
        bad_ids = points.loc[~np.asarray(inside), "id"].tolist()
        head = bad_ids[:10]
        more = "" if len(bad_ids) <= 10 else f" (and {len(bad_ids) - 10} more)"
        minx, miny, maxx, maxy = bbox.bounds
        msg = (
            f"{len(bad_ids)} observation point(s) fall outside the model "
            f"domain WGS84 bbox "
            f"[{minx:.3f}, {miny:.3f}, {maxx:.3f}, {maxy:.3f}]: "
            f"{head}{more}"
        )
        raise ValueError(msg)
    return bbox

extract_water_level_series

extract_water_level_series

extract_water_level_series(ds, points, *, variable=None)

Extract per-point water-level time series by nearest-cell lookup.

PARAMETER	DESCRIPTION
ds	Canonical dataset from a load_* reader. TYPE: Dataset
points	Must carry columns id, lon, lat in WGS84. Typically the concatenation of user-supplied points and NOAA gauges. TYPE: DataFrame
variable	Data variable to extract. When None (default), auto-detects zs then elevation. TYPE: str DEFAULT: None

RETURNS	DESCRIPTION
DataFrame	Indexed by the dataset's time axis, with one column per input point (named by its id). Values are water-surface elevation (meters, MSL) interpolated by nearest-cell lookup in the model mesh.

Notes

The unstructured-mesh lookup is brute-force — one argmin scan of all mesh nodes/faces per observation point, implemented in Python. This is intentional (avoids a SciPy runtime dep and keeps peak memory bounded), but it scales as O(n_cells * n_points). Up to a few hundred obs points against a multi-million-node mesh stays well under a second; past that, pre-build a spatial index at the caller instead of calling this function inside a loop.

Source code in src/coastal_calibration/observations.py

def extract_water_level_series(
    ds: xr.Dataset,
    points: pd.DataFrame,
    *,
    variable: str | None = None,
) -> pd.DataFrame:
    """Extract per-point water-level time series by nearest-cell lookup.

    Parameters
    ----------
    ds : xarray.Dataset
        Canonical dataset from a ``load_*`` reader.
    points : pandas.DataFrame
        Must carry columns ``id``, ``lon``, ``lat`` in WGS84. Typically
        the concatenation of user-supplied points and NOAA gauges.
    variable : str, optional
        Data variable to extract. When *None* (default), auto-detects
        ``zs`` then ``elevation``.

    Returns
    -------
    pandas.DataFrame
        Indexed by the dataset's time axis, with one column per input
        point (named by its ``id``). Values are water-surface elevation
        (meters, MSL) interpolated by nearest-cell lookup in the model
        mesh.

    Notes
    -----
    The unstructured-mesh lookup is brute-force — one ``argmin`` scan
    of all mesh nodes/faces per observation point, implemented in
    Python. This is intentional (avoids a SciPy runtime dep and keeps
    peak memory bounded), but it scales as ``O(n_cells * n_points)``.
    Up to a few hundred obs points against a multi-million-node mesh
    stays well under a second; past that, pre-build a spatial index at
    the caller instead of calling this function inside a loop.
    """
    if variable is None:
        variable = _auto_variable(ds)
    if variable not in ds.data_vars:
        available = sorted(str(v) for v in ds.data_vars)
        msg = f"Variable {variable!r} not in dataset (have: {available})"
        raise KeyError(msg)
    if points.empty:
        return pd.DataFrame(
            index=pd.DatetimeIndex(ds["time"].to_numpy(), name="time"),
        )

    qx, qy = _query_points_in_model_crs(points, ds)
    locator = _nearest_cell_indices(ds, qx, qy)
    values = _pull_time_series(ds, variable, locator)

    time = pd.DatetimeIndex(ds["time"].to_numpy(), name="time")
    ids = points["id"].tolist()
    return pd.DataFrame(values, index=time, columns=ids)

Flood Depth Map

create_flood_depth_map

create_flood_depth_map

create_flood_depth_map(
    model_root,
    dem_path,
    output_path=None,
    *,
    index_path=None,
    create_index=True,
    hmin=0.05,
    land_only=True,
    dem_offset=0.0,
    reproj_method="nearest",
    nrmax=2000,
    model=None,
    log=None,
)

Create a downscaled flood depth map from SFINCS output.

Reads the maximum water surface elevation (zsmax) from the SFINCS map output, optionally builds an index COG that maps DEM pixels to SFINCS grid cells, then downscales onto a high-resolution DEM to produce a Cloud Optimized GeoTIFF of maximum flood depth.

PARAMETER	DESCRIPTION
model_root	Path to the SFINCS model directory (must contain sfincs_map.nc). TYPE: Path or str
dem_path	Path to a high-resolution DEM GeoTIFF covering the model domain. TYPE: Path or str
output_path	Output flood depth COG path. Defaults to <model_root>/floodmap_hmax.tif. TYPE: Path or str DEFAULT: None
index_path	Path for the index COG (DEM pixels -> SFINCS cell mapping). Defaults to <model_root>/floodmap_index.tif. TYPE: Path or str DEFAULT: None
create_index	If True (default), (re)generate the index COG via :func:make_index_cog. The index speeds up the downscaling significantly for large DEMs. TYPE: bool DEFAULT: True
hmin	Minimum flood depth (m) to classify a pixel as flooded. TYPE: float DEFAULT: 0.05
land_only	Drop pixels the model shows as permanently wet, so the map is inundation rather than inundation plus the sea. TYPE: bool DEFAULT: True
dem_offset	Vertical offset (m) added to dem_path to put it on the model's datum. Non-zero when the dataset this DEM came from was merged with an offset; see :attr:~coastal_calibration.config.create_schema.ElevationDataset.offset. TYPE: float DEFAULT: 0.0
reproj_method	Reprojection method ("nearest" or "bilinear"). TYPE: str DEFAULT: 'nearest'
nrmax	Maximum cells per processing block (controls peak memory). TYPE: int DEFAULT: 2000
model	An already-loaded :class:SfincsModel instance. When provided, model_root is still used to resolve default output paths but the model is not re-read from disk. TYPE: SfincsModel DEFAULT: None
log	Logging callback accepting a single message; falls back to the module logger when None. TYPE: callable DEFAULT: None

RETURNS	DESCRIPTION
Path	Path to the generated flood depth COG.

RAISES	DESCRIPTION
FileNotFoundError	If the DEM or sfincs_map.nc cannot be found.
FloodmapInputError	If zsmax is not present in the SFINCS map output.

Source code in src/coastal_calibration/sfincs/floodmap.py

def create_flood_depth_map(
    model_root: Path | str,
    dem_path: Path | str,
    output_path: Path | str | None = None,
    *,
    index_path: Path | str | None = None,
    create_index: bool = True,
    hmin: float = 0.05,
    land_only: bool = True,
    dem_offset: float = 0.0,
    reproj_method: str = "nearest",
    nrmax: int = 2000,
    model: SfincsModel | None = None,
    log: Any = None,
) -> Path:
    """Create a downscaled flood depth map from SFINCS output.

    Reads the maximum water surface elevation (``zsmax``) from the SFINCS
    map output, optionally builds an index COG that maps DEM pixels to
    SFINCS grid cells, then downscales onto a high-resolution DEM to
    produce a Cloud Optimized GeoTIFF of maximum flood depth.

    Parameters
    ----------
    model_root : Path or str
        Path to the SFINCS model directory (must contain ``sfincs_map.nc``).
    dem_path : Path or str
        Path to a high-resolution DEM GeoTIFF covering the model domain.
    output_path : Path or str, optional
        Output flood depth COG path.  Defaults to
        ``<model_root>/floodmap_hmax.tif``.
    index_path : Path or str, optional
        Path for the index COG (DEM pixels -> SFINCS cell mapping).
        Defaults to ``<model_root>/floodmap_index.tif``.
    create_index : bool
        If True (default), (re)generate the index COG via
        :func:`make_index_cog`.  The index speeds up the downscaling
        significantly for large DEMs.
    hmin : float
        Minimum flood depth (m) to classify a pixel as flooded.
    land_only : bool
        Drop pixels the model shows as permanently wet, so the map is
        inundation rather than inundation plus the sea.
    dem_offset : float
        Vertical offset (m) added to *dem_path* to put it on the model's
        datum. Non-zero when the dataset this DEM came from was merged with
        an ``offset``; see
        :attr:`~coastal_calibration.config.create_schema.ElevationDataset.offset`.
    reproj_method : str
        Reprojection method (``"nearest"`` or ``"bilinear"``).
    nrmax : int
        Maximum cells per processing block (controls peak memory).
    model : SfincsModel, optional
        An already-loaded :class:`SfincsModel` instance.  When provided,
        ``model_root`` is still used to resolve default output paths but
        the model is **not** re-read from disk.
    log : callable, optional
        Logging callback accepting a single message; falls back to the module
        logger when *None*.

    Returns
    -------
    Path
        Path to the generated flood depth COG.

    Raises
    ------
    FileNotFoundError
        If the DEM or ``sfincs_map.nc`` cannot be found.
    FloodmapInputError
        If ``zsmax`` is not present in the SFINCS map output.
    """
    # registers the ``DataArray.raster`` accessor used below
    import hydromt  # noqa: F401  # pyright: ignore[reportUnusedImport]

    # -- Ensure patches are applied before any hydromt-sfincs call --
    from coastal_calibration.logging import suppress_hydromt_output
    from coastal_calibration.sfincs._hydromt_compat import apply_all_patches

    apply_all_patches()

    # Import *after* patches so local references pick up the fixed versions.
    with suppress_hydromt_output():
        from hydromt_sfincs.workflows.downscaling import make_index_cog

    model_root = Path(model_root)
    dem_path = Path(dem_path)

    if not dem_path.exists():
        raise FileNotFoundError(f"DEM not found: {dem_path}")

    map_file = model_root / "sfincs_map.nc"
    if not map_file.exists():
        raise FileNotFoundError(f"SFINCS map output not found: {map_file}")

    output_path = model_root / "floodmap_hmax.tif" if output_path is None else Path(output_path)
    index_path = model_root / "floodmap_index.tif" if index_path is None else Path(index_path)

    def _info(msg: str) -> None:
        if log is not None:
            log(msg)
        else:
            _log.info(msg)

    # -- Load model and read output ---------------------------------
    if model is None:
        with suppress_hydromt_output():
            from hydromt_sfincs import SfincsModel as _Sfincs

            # Use "r+" (same as the pipeline) so all components are writable
            # and the quadtree grid loads correctly.
            model = _Sfincs(root=str(model_root), mode="r+")
            model.read()

    with suppress_hydromt_output():
        model.output.read()

    if "zsmax" not in model.output.data:
        msg = (
            "Variable 'zsmax' not found in SFINCS map output. "
            "Ensure SFINCS was configured to write zsmax (storzsmax = 1 in sfincs.inp)."
        )
        raise FloodmapInputError(msg)
    zsmax = model.output.data["zsmax"]
    _info(f"Loaded zsmax from {map_file}")

    # -- (Re)create index COG ---------------------------------------
    # Always regenerate so the index stays consistent with the
    # current grid and the patched ``get_indices_at_points``.
    if create_index:
        _info(f"Creating index COG: {index_path}")
        index_path.parent.mkdir(parents=True, exist_ok=True)
        output_path.parent.mkdir(parents=True, exist_ok=True)
        with suppress_hydromt_output():
            make_index_cog(
                model=model,
                indices_fn=str(index_path),
                topobathy_fn=str(dem_path),
                nrmax=nrmax,
            )
        _assert_index_hits_the_model(index_path, dem_path)
        _ensure_overviews(index_path, _info)
        _info(f"Index COG created ({index_path.stat().st_size / 1e6:.1f} MB)")

    # -- Downscale --------------------------------------------------
    # Read the DEM and index at full resolution.  The upstream
    # ``downscale_floodmap`` defaults to ``overview_level=0`` when
    # reading rasters from disk, which silently halves the resolution
    # and can mismatch with the index.  By loading the rasters into
    # memory first and passing them as DataArrays we bypass that bug.
    _info("Downscaling flood depth map")
    output_path.parent.mkdir(parents=True, exist_ok=True)

    _write_floodmap_cog(
        zsmax=_blank_dry_cells(zsmax, model.output.data, _info),
        dem_path=dem_path,
        index_path=index_path if index_path.exists() else None,
        output_path=output_path,
        hmin=hmin,
        reproj_method=reproj_method,
        nrmax=nrmax,
        dem_offset=dem_offset,
        baseline=_baseline_water_surface(model.output.data) if land_only else None,
    )

    _ensure_overviews(output_path, _info)

    size_mb = output_path.stat().st_size / 1e6
    _info(f"Flood depth map written: {output_path} ({size_mb:.1f} MB)")

    return output_path

Downloader

validate_date_ranges

validate_date_ranges

validate_date_ranges(
    start_time,
    end_time,
    meteo_source,
    coastal_source,
    domain,
)

Validate that requested dates are within available ranges.

Source code in src/coastal_calibration/data/downloader.py

def validate_date_ranges(
    start_time: datetime,
    end_time: datetime,
    meteo_source: str,
    coastal_source: str,
    domain: str,
) -> list[str]:
    """Validate that requested dates are within available ranges."""
    errors: list[str] = []

    meteo_range = get_date_range(meteo_source, domain)
    if meteo_range:
        error = meteo_range.validate(start_time, end_time)
        if error:
            errors.append(error)

    if coastal_source != "harmonic":
        coastal_range = get_date_range(coastal_source, domain)
        if coastal_range:
            error = coastal_range.validate(start_time, end_time)
            if error:
                errors.append(error)

    return errors

NOAA CO-OPS API

COOPSAPIClient

COOPSAPIClient

COOPSAPIClient(timeout=120)

Client for interacting with NOAA CO-OPS API.

Initialize COOPS API client.

PARAMETER	DESCRIPTION
timeout	Request timeout in seconds, by default 120 TYPE: int DEFAULT: 120

RAISES	DESCRIPTION
ImportError	If plot optional dependencies are not installed.

Source code in src/coastal_calibration/data/coops_api.py

def __init__(self, timeout: int = 120) -> None:
    """Initialize COOPS API client.

    Parameters
    ----------
    timeout : int, optional
        Request timeout in seconds, by default 120

    Raises
    ------
    ImportError
        If plot optional dependencies are not installed.
    """
    _check_plot_deps()
    self.timeout = timeout
    self._stations_metadata = self._get_stations_metadata()

stations_metadata property

stations_metadata

Get metadata for all water level stations as a GeoDataFrame.

RETURNS	DESCRIPTION
GeoDataFrame	GeoDataFrame with station metadata and Point geometries.

validate_parameters

validate_parameters(
    product, datum, units, time_zone, interval
)

Validate API parameters.

PARAMETER	DESCRIPTION
product	Data product type TYPE: str
datum	Vertical datum TYPE: str
units	Unit system TYPE: str
time_zone	Time zone TYPE: str
interval	Time interval for predictions TYPE: str | int | None

RAISES	DESCRIPTION
ValueError	If any parameter is invalid

Source code in src/coastal_calibration/data/coops_api.py

def validate_parameters(
    self,
    product: str,
    datum: str,
    units: str,
    time_zone: str,
    interval: str | int | None,
) -> None:
    """Validate API parameters.

    Parameters
    ----------
    product : str
        Data product type
    datum : str
        Vertical datum
    units : str
        Unit system
    time_zone : str
        Time zone
    interval : str | int | None
        Time interval for predictions

    Raises
    ------
    ValueError
        If any parameter is invalid
    """
    if product not in self.valid_products:
        raise ValueError(
            f"Invalid product '{product}'. Must be one of: {', '.join(self.valid_products)}"
        )

    if datum.upper() not in self.valid_datums:
        raise ValueError(
            f"Invalid datum '{datum}'. Must be one of: {', '.join(self.valid_datums)}"
        )

    if units not in self.valid_units:
        raise ValueError(
            f"Invalid units '{units}'. Must be one of: {', '.join(self.valid_units)}"
        )

    if time_zone not in self.valid_timezones:
        raise ValueError(
            f"Invalid time_zone '{time_zone}'. Must be one of: {', '.join(self.valid_timezones)}"
        )

    if (
        product == "predictions"
        and interval is not None
        and str(interval) not in self.valid_intervals
    ):
        raise ValueError(
            f"Invalid interval '{interval}' for predictions. "
            f"Must be one of: {', '.join(self.valid_intervals)}"
        )

build_url

build_url(
    station_id,
    begin_date,
    end_date,
    product,
    datum,
    units,
    time_zone,
    interval,
)

Build API request URL for a station.

PARAMETER	DESCRIPTION
station_id	Station ID TYPE: str
begin_date	Start date TYPE: str
end_date	End date TYPE: str
product	Data product TYPE: str
datum	Vertical datum TYPE: str
units	Unit system TYPE: str
time_zone	Time zone TYPE: str
interval	Time interval for predictions TYPE: str | int | None

RETURNS	DESCRIPTION
str	Complete API request URL

Source code in src/coastal_calibration/data/coops_api.py

def build_url(
    self,
    station_id: str,
    begin_date: str,
    end_date: str,
    product: str,
    datum: str,
    units: str,
    time_zone: str,
    interval: str | int | None,
) -> str:
    """Build API request URL for a station.

    Parameters
    ----------
    station_id : str
        Station ID
    begin_date : str
        Start date
    end_date : str
        End date
    product : str
        Data product
    datum : str
        Vertical datum
    units : str
        Unit system
    time_zone : str
        Time zone
    interval : str | int | None, optional
        Time interval for predictions

    Returns
    -------
    str
        Complete API request URL
    """
    params = {
        "begin_date": begin_date,
        "end_date": end_date,
        "station": station_id,
        "product": product,
        "datum": datum,
        "units": units,
        "time_zone": time_zone,
        "format": "json",
        "application": "coastal_calibration_coops",
    }

    if product == "predictions" and interval is not None:
        params["interval"] = str(interval)

    query_parts = [f"{k}={v}" for k, v in params.items()]
    return f"{self.base_url}?{'&'.join(query_parts)}"

fetch_data

fetch_data(urls)

Fetch data from API for multiple URLs.

PARAMETER	DESCRIPTION
urls	List of API request URLs TYPE: list[str]

RETURNS	DESCRIPTION
list[dict | None]	List of JSON responses (None for failed requests)

Source code in src/coastal_calibration/data/coops_api.py

def fetch_data(self, urls: list[str]) -> list[dict[str, Any] | None]:
    """Fetch data from API for multiple URLs.

    Parameters
    ----------
    urls : list[str]
        List of API request URLs

    Returns
    -------
    list[dict | None]
        List of JSON responses (None for failed requests)
    """
    logger.info("  Fetching data from %d station(s)", len(urls))

    return fetch(
        urls,
        "json",
        request_method="get",
        timeout=self.timeout,
        raise_status=False,
    )

get_datums

get_datums(station_ids: str) -> StationDatum

get_datums(station_ids: list[str]) -> list[StationDatum]

get_datums(station_ids)

Retrieve datum information for one or more stations.

PARAMETER	DESCRIPTION
station_ids	Single station ID or list of station IDs TYPE: str | list[str]

RETURNS	DESCRIPTION
StationDatum | list[StationDatum]	Single StationDatum object if input is str, list of StationDatum if input is list

RAISES	DESCRIPTION
COOPSUnavailableError	If every request failed at the transport level (API outage).
ValueError	If responses came back but none had valid datum data.

Source code in src/coastal_calibration/data/coops_api.py

def get_datums(self, station_ids: str | list[str]) -> StationDatum | list[StationDatum]:
    """Retrieve datum information for one or more stations.

    Parameters
    ----------
    station_ids : str | list[str]
        Single station ID or list of station IDs

    Returns
    -------
    StationDatum | list[StationDatum]
        Single StationDatum object if input is str,
        list of StationDatum if input is list

    Raises
    ------
    COOPSUnavailableError
        If every request failed at the transport level (API outage).
    ValueError
        If responses came back but none had valid datum data.
    """
    import numpy as np

    single_input = isinstance(station_ids, str)
    ids: list[str] = [station_ids] if isinstance(station_ids, str) else list(station_ids)

    datum_base_url = "https://api.tidesandcurrents.noaa.gov/mdapi/prod/webapi/stations"
    urls = [f"{datum_base_url}/{sid}/datums.json" for sid in ids]
    logger.info("  Fetching datum information for %d station(s)", len(ids))
    responses = self.fetch_data(urls)
    datum_objects = []
    n_failed = 0
    for station_id, response in zip(ids, responses, strict=False):
        if response is None:
            n_failed += 1
            logger.warning("  No datum data returned for station %s", station_id)
            continue

        if "error" in response:
            logger.warning(
                "  Datum API error for station %s: %s",
                station_id,
                response["error"].get("message", "Unknown error"),
            )
            continue

        raw_datums = response.get("datums") or []
        datum_values = [
            DatumValue(
                name=datum_dict.get("name", ""),
                description=datum_dict.get("description", ""),
                value=np.float64(datum_dict.get("value", np.nan)),
            )
            for datum_dict in raw_datums
        ]

        station_datum = StationDatum(
            station_id=station_id,
            accepted=response.get("accepted", ""),
            superseded=response.get("superseded", ""),
            epoch=response.get("epoch", ""),
            units=response.get("units", ""),
            orthometric_datum=response.get("OrthometricDatum", ""),
            datums=datum_values,
            lat=np.float64(response.get("LAT", np.nan)),
            lat_date=response.get("LATdate", ""),
            lat_time=response.get("LATtime", ""),
            hat=np.float64(response.get("HAT", np.nan)),
            hat_date=response.get("HATdate", ""),
            hat_time=response.get("HATtime", ""),
            min_value=np.float64(response.get("min", np.nan)),
            min_date=response.get("mindate", ""),
            min_time=response.get("mintime", ""),
            max_value=np.float64(response.get("max", np.nan)),
            max_date=response.get("maxdate", ""),
            max_time=response.get("maxtime", ""),
            datum_analysis_period=response.get("DatumAnalysisPeriod") or [],
            ngs_link=response.get("NGSLink", ""),
            ctrl_station=response.get("ctrlStation", ""),
        )

        datum_objects.append(station_datum)

    if not datum_objects:
        # A None response means the fetch failed at the transport layer
        # (non-2xx / connection error); a station that merely lacks datums
        # still returns a 200 dict. So if every request failed and there was
        # at least one, the API is unreachable rather than the stations being
        # data-less. Empty input falls through to the ValueError below.
        if ids and n_failed == len(ids):
            raise COOPSUnavailableError(
                f"NOAA CO-OPS API unreachable: all {n_failed} datum request(s) failed. "
                "Likely a server error or API throttling rather than a station data "
                "problem, but the HTTP status is not surfaced here. Retry later."
            )
        raise ValueError("No valid datum data returned for any station")

    if single_input:
        return datum_objects[0]
    return datum_objects

query_coops_byids

query_coops_byids

query_coops_byids(
    station_ids,
    begin_date,
    end_date,
    *,
    product="water_level",
    datum="MLLW",
    units="metric",
    time_zone="gmt",
    interval=None,
)

Fetch water level data from NOAA CO-OPS API for multiple stations.

PARAMETER	DESCRIPTION
station_ids	List of station IDs to retrieve data for. TYPE: list[str]
begin_date	Start date in format: yyyyMMdd, yyyyMMdd HH:mm, MM/dd/yyyy, or MM/dd/yyyy HH:mm TYPE: str
end_date	End date in same format as begin_date. TYPE: str
product	Data product to retrieve, by default water_level. TYPE: ('water_level', 'hourly_height', 'high_low', 'predictions') DEFAULT: "water_level"
datum	Vertical datum for water levels, by default "MLLW". TYPE: str DEFAULT: 'MLLW'
units	Units for data, by default "metric". TYPE: ('metric', 'english') DEFAULT: "metric"
time_zone	Time zone for returned data, by default "gmt". TYPE: ('gmt', 'lst', 'lst_ldt') DEFAULT: "gmt"
interval	Time interval for predictions product only, by default None. TYPE: str | int | None DEFAULT: None

RETURNS	DESCRIPTION
Dataset	Dataset containing water level data with dimensions (time, station).

RAISES	DESCRIPTION
ValueError	If invalid parameters are provided or if API returns errors.

Source code in src/coastal_calibration/data/coops_api.py

def query_coops_byids(
    station_ids: list[str],
    begin_date: str,
    end_date: str,
    *,
    product: Literal[
        "water_level",
        "hourly_height",
        "high_low",
        "predictions",
    ] = "water_level",
    datum: str = "MLLW",
    units: Literal["metric", "english"] = "metric",
    time_zone: Literal["gmt", "lst", "lst_ldt"] = "gmt",
    interval: str | int | None = None,
) -> xr.Dataset:
    """Fetch water level data from NOAA CO-OPS API for multiple stations.

    Parameters
    ----------
    station_ids : list[str]
        List of station IDs to retrieve data for.
    begin_date : str
        Start date in format: yyyyMMdd, yyyyMMdd HH:mm, MM/dd/yyyy, or MM/dd/yyyy HH:mm
    end_date : str
        End date in same format as begin_date.
    product : {"water_level", "hourly_height", "high_low", "predictions"}, optional
        Data product to retrieve, by default ``water_level``.
    datum : str, optional
        Vertical datum for water levels, by default "MLLW".
    units : {"metric", "english"}, optional
        Units for data, by default "metric".
    time_zone : {"gmt", "lst", "lst_ldt"}, optional
        Time zone for returned data, by default "gmt".
    interval : str | int | None, optional
        Time interval for predictions product only, by default None.

    Returns
    -------
    xr.Dataset
        Dataset containing water level data with dimensions (time, station).

    Raises
    ------
    ValueError
        If invalid parameters are provided or if API returns errors.
    """
    client = COOPSAPIClient()
    client.validate_parameters(product, datum, units, time_zone, interval)
    begin_dt = client.parse_date(begin_date)
    end_dt = client.parse_date(end_date)

    if end_dt <= begin_dt:
        raise ValueError("end_date must be after begin_date")

    begin_str = begin_dt.strftime("%Y%m%d %H:%M")
    end_str = end_dt.strftime("%Y%m%d %H:%M")

    logger.info(
        "  Requesting %s data for %d station(s) from %s to %s",
        product,
        len(station_ids),
        begin_str,
        end_str,
    )

    urls = [
        client.build_url(
            station_id=station_id,
            begin_date=begin_str,
            end_date=end_str,
            product=product,
            datum=datum,
            units=units,
            time_zone=time_zone,
            interval=interval,
        )
        for station_id in station_ids
    ]

    return _process_responses(
        responses=client.fetch_data(urls),
        station_ids=station_ids,
        product=product,
        datum=datum,
        units=units,
        time_zone=time_zone,
    )

query_coops_bygeometry

query_coops_bygeometry

query_coops_bygeometry(
    geometry,
    begin_date,
    end_date,
    *,
    product="water_level",
    datum="MLLW",
    units="metric",
    time_zone="gmt",
    interval=None,
)

Fetch water level data from NOAA CO-OPS API for stations within a geometry.

PARAMETER	DESCRIPTION
geometry	Geometry to select stations within (Point, Polygon, etc.) TYPE: BaseGeometry
begin_date	Start date in format: yyyyMMdd, yyyyMMdd HH:mm, MM/dd/yyyy, or MM/dd/yyyy HH:mm TYPE: str
end_date	End date in same format as begin_date. TYPE: str
product	Data product to retrieve, by default water_level. TYPE: ('water_level', 'hourly_height', 'high_low', 'predictions') DEFAULT: "water_level"
datum	Vertical datum for water levels, by default "MLLW". TYPE: str DEFAULT: 'MLLW'
units	Units for data, by default "metric". TYPE: ('metric', 'english') DEFAULT: "metric"
time_zone	Time zone for returned data, by default "gmt". TYPE: ('gmt', 'lst', 'lst_ldt') DEFAULT: "gmt"
interval	Time interval for predictions product only, by default None. TYPE: str | int | None DEFAULT: None

RETURNS	DESCRIPTION
Dataset	Dataset containing water level data for stations within the geometry.

Source code in src/coastal_calibration/data/coops_api.py

def query_coops_bygeometry(
    geometry: BaseGeometry,
    begin_date: str,
    end_date: str,
    *,
    product: Literal[
        "water_level",
        "hourly_height",
        "high_low",
        "predictions",
    ] = "water_level",
    datum: str = "MLLW",
    units: Literal["metric", "english"] = "metric",
    time_zone: Literal["gmt", "lst", "lst_ldt"] = "gmt",
    interval: str | int | None = None,
) -> xr.Dataset:
    """Fetch water level data from NOAA CO-OPS API for stations within a geometry.

    Parameters
    ----------
    geometry : shapely.geometry.base.BaseGeometry
        Geometry to select stations within (Point, Polygon, etc.)
    begin_date : str
        Start date in format: yyyyMMdd, yyyyMMdd HH:mm, MM/dd/yyyy, or MM/dd/yyyy HH:mm
    end_date : str
        End date in same format as begin_date.
    product : {"water_level", "hourly_height", "high_low", "predictions"}, optional
        Data product to retrieve, by default ``water_level``.
    datum : str, optional
        Vertical datum for water levels, by default "MLLW".
    units : {"metric", "english"}, optional
        Units for data, by default "metric".
    time_zone : {"gmt", "lst", "lst_ldt"}, optional
        Time zone for returned data, by default "gmt".
    interval : str | int | None, optional
        Time interval for predictions product only, by default None.

    Returns
    -------
    xr.Dataset
        Dataset containing water level data for stations within the geometry.
    """
    import numpy as np
    import shapely

    client = COOPSAPIClient()
    if not all(shapely.is_valid(np.atleast_1d(geometry))):  # pyright: ignore[reportCallIssue, reportArgumentType]
        raise ValueError("Invalid geometry provided.")

    stations_gdf = client.stations_metadata
    selected_stations = stations_gdf[stations_gdf.intersects(geometry)]

    if selected_stations.empty:
        raise ValueError("No stations found within the specified geometry and buffer.")

    station_ids = selected_stations["station_id"].tolist()
    return query_coops_byids(
        station_ids,
        begin_date,
        end_date,
        product=product,
        datum=datum,
        units=units,
        time_zone=time_zone,
        interval=interval,
    )

Type Aliases

# Model type
ModelType = Literal["schism", "sfincs"]

# Meteorological data source
MeteoSource = Literal["nwm_retro", "nwm_ana"]

# Coastal domain identifier
CoastalDomain = Literal["prvi", "hawaii", "atlgulf", "pacific"]

# Boundary condition source
BoundarySource = Literal["tpxo", "stofs"]

# Logging level
LogLevel = Literal["DEBUG", "INFO", "WARNING", "ERROR", "CRITICAL"]

Constants

Default Paths

DEFAULT_NFS_MOUNT = Path("/ngen-test")

Default Path Templates

DEFAULT_WORK_DIR_TEMPLATE = (
    "/ngen-test/coastal/${user}/"
    "${model}_${simulation.coastal_domain}_${boundary.source}_${simulation.meteo_source}/"
    "${model}_${simulation.start_date}"
)

DEFAULT_RAW_DOWNLOAD_DIR_TEMPLATE = (
    "/ngen-test/coastal/${user}/"
    "${model}_${simulation.coastal_domain}_${boundary.source}_${simulation.meteo_source}/"
    "raw_data"
)

Model Registry

MODEL_REGISTRY: dict[str, type[ModelConfig]] = {
    "schism": SchismModelConfig,
    "sfincs": SfincsModelConfig,
}