I/O Module

Run

class elm_diagnostics.io.run.Run(path, name=None, streams=None, chunks=None, chunk_mode='auto', chunk_target_mb=64, analysis_year=None, analysis_year_min=None, analysis_year_max=None, analysis_year_tolerance=0, strict_combine=False)[source]

Bases: object

Atomic unit of analysis: one ELM case’s history-file streams.

Parameters:

path (str | Path) – Directory containing *.elm.h*.nc files, or a glob pattern.
name (str | None) – Display name. Defaults to the case name extracted from filenames.
streams (dict[str, str] | None) – Explicit stream mapping, e.g. {"h0": "*.elm.h0.*.nc"}. If None, streams are auto-discovered.
chunks (dict | None) – Passed to xr.open_mfdataset for dask-backed lazy loading.
analysis_year (int | None) – Requested analysis year for early file narrowing before open.
analysis_year_min (int | None) – Inclusive lower bound for year-aware file narrowing.
analysis_year_max (int | None) – Inclusive upper bound for year-aware file narrowing.
analysis_year_tolerance (int) – Year-window half-width when narrowing files (0 means exact year).
strict_combine (bool) – If True, open streams using stricter multi-file combine options. Defaults to False.
chunk_mode (Literal['off', 'auto', 'manual'])
chunk_target_mb (int)

bounds_dataset(tape=None)[source]

Return the stream dataset used for flux-integration time bounds.

Balance modules need time_bounds (and the time coordinate) for cumulative integration. This returns the cached full stream for the tape that carries the bounds, opened once via _open_stream. With the compat="override" combine strategy, that open is already cheap even for many-file, many-variable streams, so no separate bounds-only open is needed.

Return type:: Dataset
Parameters:: tape (str | None)

property streams: dict[str, Dataset]: All streams as open datasets, keyed by tape name.

property cadence: dict[str, str | Timedelta]: Cadence per stream (inferred from time_bounds).

property tape_priority: list[str]: Tapes ordered by cadence (finest first).

get(varname, tape=None)[source]

Retrieve a variable, searching tapes in priority order.

If the variable is not found, attempts to derive it from available components (e.g., compute QFLX_EVAP_TOT from QSOIL + QVEGE + QVEGT).

Variables are cached with an LRU policy (max 15 variables) to balance performance and memory usage on large datasets.

The per-tape variable-name index (single file header) is used to route to the tape that owns varname without opening the other tapes’ datasets. The owning stream is then opened once (cached via _open_stream) and the variable sliced out; the open is cheap thanks to the compat="override" combine strategy.

Parameters:

varname (str) – History field name (e.g. "GPP", "SOILLIQ").
tape (str | None) – Specific tape to search. If None, searches all tapes in cadence-priority order (finest first).

Return type:

DataArray

Raises:

KeyError – If the variable is not found in any tape and cannot be derived.

has(varname)[source]

Check whether a variable exists in any tape or can be derived.

Uses the per-tape header index (single file read) rather than opening full streams, and checks derivability via component availability (DERIVABLE_REQUIREMENTS) instead of executing the derivation. Both keep this cheap on high-variable-count datasets.

Return type:: bool
Parameters:: varname (str)

close()[source]

Close all open datasets and clear caches.

Return type:: None

Comparison

class elm_diagnostics.io.run.Comparison(base, experiment, align='intersect')[source]

Bases: object

Pair of runs for side-by-side diagnostics.

Parameters:

base (Run) – Reference / control run.
experiment (Run) – Experiment / perturbation run.
align (Literal['intersect', 'union']) – How to align time axes. 'intersect' keeps only overlapping times; 'union' fills missing times with NaN.

get(varname, tape=None)[source]

Retrieve a variable from both runs, time-aligned.

Alignment preserves dask chunks for lazy evaluation. Computation is deferred until the data is actually used in plot generation.

Return type:

tuple[DataArray, DataArray]

Parameters:

varname (str)
tape (str | None)

Derived Variables

Compute derived ELM variables from available history output.

This module provides functions to compute commonly-needed variables that may not be in the default h0 output, but can be calculated from component variables.

Based on ELM source code analysis (April 2026): - Total ET from components: QSOIL + QVEGE + QVEGT - Total storage from vertical profiles: sum(SOILLIQ) + sum(SOILICE) + …

elm_diagnostics.io.derived.compute_total_et(run)[source]

Compute total evapotranspiration if QFLX_EVAP_TOT is missing.

Based on ELM source (SoilFluxesMod.F90, VegetationDataType.F90):: QFLX_EVAP_TOT = QSOIL + QVEGE + QVEGT
Where:: QSOIL = Ground evaporation (soil/snow evap + sublimation - dew) QVEGE = Canopy evaporation (evap from leaves and stems) QVEGT = Canopy transpiration (stomatal)

Parameters:: run (Run) – Run object containing the necessary component variables.
Returns:: Total evapotranspiration in mm/s (or units of components).
Return type:: DataArray
Raises:: ValueError – If required component variables are not available.

elm_diagnostics.io.derived.aggregate_vertical_storage(run, varname, vertical_dim='levgrnd')[source]

Aggregate a vertical profile storage variable to column total.

For variables like SOILLIQ(time, levgrnd, …) or SOILICE(time, levgrnd, …), sum over the vertical dimension to get total column storage.

Parameters:

run (Run) – Run object containing the variable.
varname (str) – Variable name (e.g., “SOILLIQ”, “SOILICE”).
vertical_dim (str) – Name of the vertical dimension. Default is “levgrnd”. Will auto-detect from [“levgrnd”, “levsoi”, “levdcmp”] if not specified.

Returns:

Column-total storage, summed over vertical levels.

Return type:

DataArray

elm_diagnostics.io.derived.compute_total_soil_water(run)[source]

Compute total soil water (liquid + ice) column storage.

Sums SOILLIQ and SOILICE over vertical levels.

Parameters:: run (Run) – Run object.
Returns:: Total soil water in kg/m² (or mm, depending on units).
Return type:: DataArray

elm_diagnostics.io.derived.compute_total_precip(run)[source]

Compute total precipitation rate from liquid and frozen components.

Based on ELM forcing conventions:: PRECT = RAIN + SNOW
Where:: RAIN = liquid precipitation rate (forc_rain) SNOW = frozen precipitation rate (forc_snow)

Parameters:: run (Run) – Run object containing RAIN and SNOW history variables.
Returns:: Total precipitation rate in the same units as RAIN (typically mm/s).
Return type:: DataArray
Raises:: ValueError – If either RAIN or SNOW is not available in the run.

elm_diagnostics.io.derived.can_derive(varname)[source]

Check if a variable can be derived from components.

Parameters:: varname (str) – Variable name to check.
Returns:: True if variable can be derived.
Return type:: bool

elm_diagnostics.io.derived.derive_variable(run, varname)[source]

Derive a variable from available components.

Parameters:

run (Run) – Run object.
varname (str) – Variable name to derive.

Returns:

Derived variable.

Return type:

DataArray

Raises:

ValueError – If variable cannot be derived or required components are missing.

Sub-gridcell Support

Sub-gridcell hierarchy detection and helpers.

elm_diagnostics.io.subgrid.detect_subgrid_dims(ds)[source]

Return the set of sub-gridcell dimensions present in a dataset.

If the set is empty, the dataset uses gridcell-averaged output (dov2xy = .true.).

Return type:: set[str]
Parameters:: ds (Dataset)

elm_diagnostics.io.subgrid.has_subgrid(ds)[source]

Check whether a dataset has sub-gridcell dimensions.

Return type:: bool
Parameters:: ds (Dataset)

elm_diagnostics.io.subgrid.validate_by_keyword(ds, by)[source]

Validate that the by keyword is compatible with the dataset.

Raises:

ValueError – If by is requested but the dataset is gridcell-averaged, or if the requested level isn’t present.

Return type:

None

Parameters:

ds (Dataset)
by (Literal['column', 'pft', 'landunit'] | None)

elm_diagnostics.io.subgrid.get_subgrid_level(da)[source]

Determine the sub-gridcell level of a DataArray, if any.

Return type:: Optional[Literal['column', 'pft', 'landunit']]
Parameters:: da (DataArray)

Units

Unit handling for ELM history variables using pint / pint-xarray.

elm_diagnostics.io.units.get_registry()[source]

Return the shared pint unit registry.

Return type:: UnitRegistry

elm_diagnostics.io.units.normalize_unit_string(raw)[source]

Convert an ELM unit string to a pint-parseable form.

Return type:: str
Parameters:: raw (str)

elm_diagnostics.io.units.parse_units(raw)[source]

Parse an ELM unit string into a pint Unit.

Return type:: Unit
Parameters:: raw (str)

elm_diagnostics.io.units.classify_variable(da)[source]

Classify a variable as flux, state, or intensive.

Uses (in priority order): 1. Known variable name lists. 2. cell_methods attribute (contains "time: mean" → flux-like). 3. Unit string inspection (contains /s → flux).

Return type:: Literal['flux', 'state', 'intensive']
Parameters:: da (DataArray)

elm_diagnostics.io.units.convert_flux_to_cumulative_units(units_str)[source]

Determine the target cumulative units for a flux.

Returns (target_unit_string, seconds_multiplier).

Return type:: tuple[str, float]
Parameters:: units_str (str)

Examples

>>> convert_flux_to_cumulative_units("mm/s")
('mm', 1.0)
>>> convert_flux_to_cumulative_units("gC/m^2/s")
('g/m**2', 1.0)
>>> convert_flux_to_cumulative_units("W/m^2")
('J/m**2', 1.0)

elm_diagnostics.io.units.convert_water_to_mm(da)[source]

Convert water storage variable to mm units.

Water mass per unit area (kg/m²) and water depth (mm) are numerically equivalent for liquid water: 1 kg/m² = 1 mm H2O (assuming density = 1000 kg/m³).

This function standardizes the units attribute to “mm” without changing values. Variables already in mm are returned unchanged.

Parameters:: da (DataArray) – Water storage variable with units in kg/m², mm, or variants.
Returns:: Variable with units standardized to “mm”. Values are unchanged.
Return type:: DataArray
Raises:: ValueError – If units cannot be converted to mm (e.g., temperature, pressure).

Examples

>>> soilliq = xr.DataArray([100.0, 150.0], attrs={"units": "kg/m2"})
>>> result = convert_water_to_mm(soilliq)
>>> result.attrs["units"]
'mm'
>>> result.values
array([100., 150.])