API Reference

This page documents the complete public API exported from the mnished package. See also:

Buckets

The primary class for running a watershed simulation.

class mnished.Buckets(T_monthly_normals=None)[source]

Bases: object

Incorporates a list of reservoirs into a linear hierarchy that sends water either downwards or out to the surface. Reservoirs are ordered from top (nearest Earth’s surface) to bottom (deepest groundwater); this order controls the direction of infiltration between layers.

MNiShed is designed as a daily-timestep model. This is a deliberate design choice: the physical parameterisations — degree-day snowmelt, Thornthwaite ET, and linear reservoir drainage — are climatological approximations that are well-founded at daily resolution but lose physical meaning at finer scales. The daily timestep is enforced in initialize().

__init__(T_monthly_normals=None)[source]

Initialize the watershed model.

If using the ThornthwaiteChang2019 ET method, pass T_monthly_normals here so that the thermal index I and exponent a are computed once from climatological normals and remain fixed throughout the simulation.

Parameters:

T_monthly_normals (array-like of length 12, optional) – Long-term mean monthly temperatures (°C) used to compute the Thornthwaite thermal index I and exponent a per Chang et al. (2019), https://doi.org/10.1002/ird.2309. Required when evapotranspiration_method is ‘ThornthwaiteChang2019’.

export_Hlist()[source]

Return the current water depths in all subsurface reservoirs.

Useful for checkpointing reservoir state between a spin-up run and the main simulation, or for restarting a paused run.

Returns:

Water depth in each reservoir, ordered from shallowest (index 0) to deepest.

Return type:

list of float

initialize(config_file=None, enforce_water_balance=None)[source]

Set up the model from a YAML configuration file.

Reads the configuration file, loads the input time series, builds the reservoir stack, instantiates snowpack if temperature data are present, optionally computes the water-year ET multiplier, and runs any requested spin-up cycles. Part of the CSDMS Basic Model Interface.

Parameters:

  • config_file (str, optional) – Path to the YAML configuration file. If None, all required values must be set on the object directly before calling update().

  • enforce_water_balance (str or None, optional) –

    How to scale ET to close the water balance. When None (default), the value is read from general: enforce_water_balance in the YAML config, which itself defaults to 'water-year' if absent. Accepted values:

    • 'water-year' — Scale ET by a per-water-year multiplier so that P - Q - ET = 0 over each water year.

    • 'global' — Scale ET by a single multiplier computed from sum(P - Q_obs) / sum(ET_raw) over the full record. No per-year overfitting; does not add hidden degrees of freedom to AIC comparisons.

    • 'none' — Use raw ET without correction. Appropriate only when supplying trusted measured ET (e.g. eddy covariance). Using 'none' with ThornthwaiteChang2019 will raise a warning because Thornthwaite ET carries large systematic biases.

compute_water_year()[source]

Assign a water-year label to each row in self.hydrodata.

Adds a ‘Water Year’ column. A water year begins on water_year_start_month and is labelled by the calendar year in which it ends. For example, with a start month of October (USGS convention), 1 Oct 2020 – 30 Sep 2021 is water year 2021.

When water_year_start_month is 1 (January), the water year equals the calendar year and no offset is applied.

compute_global_ET_multiplier(start=None, end=None)[source]

Compute a single ET scale factor to close the long-term water balance.

Computes scale = sum(P - Q_obs) / sum(ET_raw) over all days where P, ET, and Q_obs are all finite, and stores it as self.global_et_multiplier. Unlike compute_ET_multiplier(), which fits one multiplier per water year, this uses a single ratio and does not overfit interannual variability. Appropriate when enforce_water_balance is set to ‘global’.

Parameters:

  • start (str or datetime-like, optional) – Start of the window used to compute the multiplier (inclusive). If None, uses the full record.

  • end (str or datetime-like, optional) – End of the window used to compute the multiplier (inclusive). If None, uses the full record.

compute_ET_multiplier()[source]

Compute per-water-year ET scaling factors to enforce water balance.

For each water year, computes the ratio of required ET (P - Q) to measured or computed ET, and stores this as ‘ET multiplier’ in self.hydrodata_WY_means. This multiplier is later applied in compute_ET() to scale ET so that P - Q - ET = 0 over each water year.

compute_ET()[source]

Build the ET time series used in the model.

Obtains raw daily ET from the input data file or the Thornthwaite–Chang 2019 equation (see evapotranspiration_Chang2019()). Five modes:

et_scale (default 1.0) is applied as a universal final multiplier in every mode. Four modes determine the base ET before et_scale:

  1. et_water_stress=True (or et_reservoir_draw + ‘none’): base = raw ET. et_scale is then the sole WB-closure mechanism.

  2. enforce_water_balance=’water-year’ (stress off): base = raw ET × per-water-year multiplier so P - Q - ET ≈ 0 each year (before et_scale).

  3. enforce_water_balance=’global’ (stress off): base = raw ET × single full-record multiplier from sum(P - Q_obs) / sum(ET_raw).

  4. enforce_water_balance=’none’ (stress off): base = raw ET.

In modes 2–4, et_scale=1.0 (the default) leaves existing behaviour unchanged. Set et_scale ≠ 1 to add a calibrated offset from the water-balance multiplier (e.g. to capture decade-specific land-cover or climate sensitivity); this relaxes exact WB closure.

The result is stored as ‘ET for model [mm/day]’ in self.hydrodata.

update(time_step=None)[source]

Advance the model by one time step.

Routes precipitation minus ET through the snowpack (if present) and then through each subsurface reservoir in order from shallowest to deepest. Stores modeled specific discharge, snowpack SWE, and total subsurface storage in self.hydrodata for the current time step. Part of the CSDMS Basic Model Interface.

Parameters:

time_step (int, optional) – Index into self.hydrodata for the time step to update. If None, uses and then increments the internal counter self._timestep_i.

evapotranspiration_Chang2019(Tmax=None, Tmin=None, photoperiod=None, k=0.69)[source]

Modified daily Thornthwaite ET₀ equation.

Chang et al. (2019), Eq. 1–4. https://doi.org/10.1002/ird.2309

Parameters:

  • Tmax (array-like) – Daily maximum temperature (°C).

  • Tmin (array-like) – Daily minimum temperature (°C).

  • photoperiod (array-like) – Photoperiod N (hours), computed from latitude and Julian day per Allen et al. (1998), Eqs. 2–4 of Chang et al. (2019).

  • k (float) – Calibration coefficient in the T_ef formula. Default 0.69, recommended by Pereira & Pruitt (2004) for daily ET₀ (https://doi.org/10.1016/j.agrformet.2004.01.005). Use 0.72 for monthly ET₀ per Camargo et al. (1999).

Returns:

ET0 – Daily reference evapotranspiration (mm day⁻¹).

Return type:

array-like

run(start=None, end=None, store_depths=False)[source]

Advance the model through time steps in self.hydrodata.

Resets the internal time counter before iterating, so run() is safe to call more than once (e.g. spin-up then main run). Captures storage at the start of the run for check_mass_balance(). Part of the CSDMS Basic Model Interface.

Parameters:

  • start (str or datetime-like, optional) – First date to simulate (inclusive). If None, starts from the beginning of self.hydrodata.

  • end (str or datetime-like, optional) – Last date to simulate (inclusive). If None, runs to the end of self.hydrodata.

Notes

When start or end is provided, only the matching rows are stepped through; reservoir states at the boundaries carry in/out naturally, making windowed runs chainable (e.g. spin-up on pre-decade data followed by a scored run on the decade itself).

finalize()[source]

Report model skill and display output plots.

Calls compute_NSE(verbose=True) to print the Nash–Sutcliffe Efficiency to stdout, then calls plot() to display a time-series comparison of observed and modeled specific discharge. Part of the CSDMS Basic Model Interface.

plot()[source]

Display a time-series comparison of precipitation and specific discharge.

Produces a dual-axis figure: precipitation as a bar chart on the left axis and both observed and modeled specific discharge as line plots on the right axis.

check_mass_balance(time_step=None)[source]

Compute the mass-balance discrepancy at a given time step.

Compares cumulative inputs (P - ET) from the start of the record through time_step with cumulative outputs (discharge) plus current storage (snowpack + subsurface reservoirs) and any carried-over deficit. Returns the excess mass still in the model; a value near zero indicates good mass conservation.

Parameters:

time_step (int, optional) – Row index in self.hydrodata at which to evaluate the balance. Defaults to the last row.

Returns:

excess_mass_in_model – Excess water remaining in the model budget (mm). Should be ~0 for a mass-conserving run.

Return type:

float

compute_NSE(return_nse=True, verbose=False)[source]

Compute the Nash–Sutcliffe Efficiency of the discharge simulation.

Compares modeled and observed specific discharge for all time steps where both values are non-missing. Stores the result as self.NSE.

Parameters:

  • return_nse (bool, optional) – If True (default), return the NSE value.

  • verbose (bool, optional) – If True, print the NSE value to stdout.

Returns:

NSE – Nash–Sutcliffe Efficiency coefficient. Returns None if return_nse is False. A value of 1 indicates perfect agreement; values below 0 indicate the model performs worse than the observed-mean predictor.

Return type:

float or None

Reservoir

A single power-law (or linear) reservoir.

class mnished.Reservoir(recession_coeff, f_to_discharge=1.0, Hmax=inf, pdm_H0=None, H0=0.0, f_tile=0.0, tau_tile=None, junction_type='fraction', leakance_R=None, H_threshold=0.0, multipath_threshold=None, multipath_timescale=None)[source]

Bases: object

Generic reservoir. Accepts new water (recharge), and sends it to other reservoirs and/or out of the system (discharge) at a rate that is proportional to the amount of water held in the reservoir.

__init__(recession_coeff, f_to_discharge=1.0, Hmax=inf, pdm_H0=None, H0=0.0, f_tile=0.0, tau_tile=None, junction_type='fraction', leakance_R=None, H_threshold=0.0, multipath_threshold=None, multipath_timescale=None)[source]

Initialize a reservoir.

Parameters:

  • recession_coeff (float) – Recession coefficient. For a linear reservoir (recession exponent b = 1) this equals the true e-folding timescale, in days. For b > 1 it is not a timescale — its units are day·mm^(b-1) and the actual mean residence time depends on storage level; use mean_residence_time() to obtain a physically comparable timescale at a given reference discharge.

  • f_to_discharge (float, optional) – Fraction of water lost each time step that exits as river discharge. The remainder (1 - f_to_discharge) infiltrates to the next-deeper reservoir. Default 1.0 (all to discharge). Used only when junction_type is ‘fraction’ or ‘threshold’.

  • Hmax (float, optional) – Maximum water depth the reservoir can hold. Default np.inf.

  • pdm_H0 (float or None, optional) – Characteristic storage depth [mm] for the probability-distributed model (PDM) of saturation-excess overland flow. Storage capacity is assumed exponentially distributed across the catchment with mean pdm_H0; the saturated fraction when reservoir depth is H is f_sat = 1 - exp(-H / pdm_H0). That fraction of each positive recharge pulse is shed immediately as overland flow. Mutually exclusive with a finite Hmax. Default None (PDM off).

  • H0 (float, optional) – Initial water depth at the start of the simulation. Default 0.

  • f_tile (float, optional) –

    Fraction of subsurface drainage (H_to_next) that is diverted to a fast tile-drain sub-reservoir instead of passing to the next-deeper reservoir. The tile sub-reservoir drains directly to stream with e-folding time tau_tile. This is a fractional-bypass tile representation: a constant fraction of recession outflow is routed through a downstream fast linear reservoir, regardless of current storage. Default 0 (no tile drainage). Requires tau_tile when > 0.

    For a threshold-activated parallel drainage representation, in which a second outflow path from the same storage activates only above a water-table depth, use multipath_threshold / multipath_timescale instead.

  • tau_tile (float or None, optional) – E-folding residence time [days] of the tile-drain sub-reservoir (used with f_tile). Typical values: 3–21 days for agricultural tile systems. Required when f_tile > 0; ignored when f_tile == 0. Default None.

  • multipath_threshold (float or None, optional) –

    Storage depth [mm] above which a parallel fast drainage path from this reservoir activates. The parallel path drains directly to stream with e-folding timescale multipath_timescale and is added to the discharge alongside the primary recession. The outflow is

    Q_multipath = max(0, H - multipath_threshold) / multipath_timescale.

    Below the threshold this term is exactly zero, so the reservoir behaves identically to the no-multipath case. Above it, a second linear drain operates in parallel with the primary recession, giving a two-timescale response: slow matrix drainage at low storage and faster combined drainage at high storage.

    Physical analogue: a tile-drain system installed at a fixed depth that activates only once the water table rises above the drain elevation. Contrast with f_tile/tau_tile, which is a constant-fraction fast-routing path through a separate downstream reservoir regardless of storage.

    Both multipath_threshold and multipath_timescale must be provided together, or both left as None (default = multipath disabled).

  • multipath_timescale (float or None, optional) – E-folding timescale [days] of the parallel multipath drain. Typical values: 3–21 days for agricultural tile systems represented this way. Required together with multipath_threshold; ignored when it is None. Default None.

  • junction_type (str, optional) –

    How drainage is partitioned between river discharge and the next-deeper reservoir. Options:

    'fraction' (default)

    Fixed split: f_to_discharge fraction exits as stream discharge; remainder infiltrates. Current default behaviour.

    'leakance'

    Head-difference driven flow to the next reservoir, physically representing Darcy flow through a confining unit (e.g. a shale layer). Q_leak = max(H_this - H_next, 0) / leakance_R, capped at the total drainage dH. The remainder goes to stream. Requires leakance_R.

    'threshold'

    Dead-storage threshold: the recession law is applied only to max(H - H_threshold, 0); water below H_threshold never drains. The above-threshold drainage splits by f_to_discharge as in the 'fraction' case. Models a stream-aquifer connection that activates only when the water table exceeds the streambed elevation.

  • leakance_R (float or None, optional) – Leakance resistance [days]. Q_leak = ΔH / leakance_R. Required when junction_type is 'leakance'. Default None.

  • H_threshold (float, optional) – Dead-storage threshold depth [mm]. Recession is applied only to max(H - H_threshold, 0). Used when junction_type is 'threshold'. Default 0.0 (no threshold; full storage drains).

Raises:

ValueError – If recession_coeff <= 0, f_to_discharge < 0 or > 1, Hmax < 0, pdm_H0 <= 0, f_tile < 0 or > 1, f_tile > 0 with no tau_tile, junction_type is unrecognised, junction_type is ‘leakance’ with no leakance_R, multipath_threshold < 0, multipath_timescale <= 0, or only one of the two multipath parameters is provided.

property has_multipath

True when both multipath parameters are configured (drain active).

recharge(H)[source]

Add or remove water from the reservoir.

Recharge H can be positive (net water input, e.g. P > ET) or negative (net deficit, e.g. ET > P). Sets H_excess if the reservoir overflows Hmax, or H_deficit if more water is removed than the reservoir holds.

Parameters:

H (float) – Depth of water added (positive) or removed (negative).

Raises:

ValueError – If Hwater is already negative before recharge is applied.

discharge(dt, H_next=None)[source]

Discharge water from the reservoir over one time step.

Computes water lost by the recession law, partitions it between river discharge (H_exfiltrated) and infiltration to the next-deeper reservoir (H_to_next) according to junction_type, and adds overflow from recharge() (H_excess) to H_discharge.

Parameters:

  • dt (float) – Time step duration (same units as recession_coeff; typically days).

  • H_next (float or None, optional) – Current water depth of the next-deeper reservoir [mm]. Used only when junction_type is 'leakance' to compute the head-difference driven flux. Ignored for other junction types.

mean_residence_time(Q_ref)[source]

Mean residence time [days] at a reference steady-state discharge.

For a nonlinear reservoir governed by \(Q = (H/\tau)\cdot(H/H_{\mathrm{ref}})^{b-1} = H^{b}/\tau_{\mathrm{eff}}\), with effective constant \(\tau_{\mathrm{eff}} = \tau\,H_{\mathrm{ref}}^{\,b-1}\), the steady-state storage at discharge Q_ref is \(H_{ss} = (Q_{\mathrm{ref}} \cdot \tau_{\mathrm{eff}})^{1/b}\), giving

\[\mathrm{MRT} = \frac{H_{ss}}{Q_{\mathrm{ref}}} = \frac{\tau_{\mathrm{eff}}^{1/b}}{Q_{\mathrm{ref}}^{\,1 - 1/b}}\]

With the default \(H_{\mathrm{ref}} = 1\) this is \(\tau^{1/b} / Q_{\mathrm{ref}}^{\,1-1/b}\).

For a linear reservoir (b = 1) this reduces exactly to \(\tau\). For b > 1, MRT is smaller than \(\tau\) whenever \(Q_{\mathrm{ref}} > 1\) mm/day, reflecting the faster drainage at realistic operating storage depths.

Unlike recession_coeff, which equals a timescale only when b=1 (linear reservoir), MRT is a physically comparable timescale across reservoirs with different recession exponents.

Parameters:

Q_ref (float) – Representative steady-state discharge from this reservoir [mm/day]. Use the long-term mean flux attributed to this layer (e.g. mean annual discharge partitioned by reservoir).

Returns:

Mean residence time [days].

Return type:

float

Raises:

ValueError – If Q_ref <= 0.

Snowpack

Snowpack accumulation and positive-degree-day melt.

class mnished.Snowpack(melt_factor=None)[source]

Bases: object

Snowpack reservoir driven by temperature.

Accumulates precipitation as snow when mean temperature is at or below 0 °C. Melts at a positive-degree-day rate when temperature is above 0 °C. All melt is routed to the top subsurface reservoir as infiltration; direct discharge to the river is not modeled.

Should precede the subsurface reservoir list in a watershed model. The melt factor is a positive-degree-day factor [mm/°C/day].

__init__(melt_factor=None)[source]

Initialize an empty snowpack.

Parameters:

melt_factor (float, optional) – Positive-degree-day melt factor (mm SWE °C⁻¹ day⁻¹). Can be set or updated later via set_melt_factor().

set_melt_factor(melt_factor)[source]

Set or update the positive-degree-day melt factor.

Parameters:

melt_factor (float) – Melt rate per positive degree-day (mm SWE °C⁻¹ day⁻¹).

set_temperature(T)[source]

Set the mean air temperature for the current time step.

Parameters:

T (float) – Mean air temperature (°C).

recharge(H)[source]

Apply net water input or deficit to the snowpack.

If T <= 0, positive H accumulates as snow (SWE). If T > 0, positive H bypasses the snowpack and is passed directly to the top subsurface reservoir via H_infiltrated. Negative H (ET > P) is removed from the snowpack as sublimation; any remainder that exceeds available SWE becomes H_deficit.

Parameters:

H (float) – Net water depth for this time step (mm). Positive = input (P - ET > 0); negative = deficit (ET - P > 0).

melt(dt, P=0.0)[source]

Compute positive-degree-day and rain-on-snow melt; update state.

Both terms are routed to H_infiltrated (→ top soil reservoir). If total available energy exceeds the SWE present, the leftover is returned as equivalent degree-days so the caller can credit it toward frozen-soil thawing (FGI reduction) rather than losing it.

Rain-on-snow sensible heat: water arriving at T_mean > 0 °C carries (c_p / L_f) · T · P mm SWE of latent-heat capacity. Spring snowpacks are near-isothermal at 0 °C, so cold-content corrections are negligible and the latent-heat term dominates.

References

McCabe et al. (2007) doi:10.1175/BAMS-88-3-319 Würzer et al. (2016) doi:10.1175/JHM-D-15-0181.1

Parameters:

  • dt (float) – Timestep [days].

  • P (float, optional) – Raw liquid precipitation [mm/day]. Used to compute rain-on-snow sensible-heat melt. Default 0 (PDD only).

Returns:

excess_dd – Leftover melt energy after the snowpack is exhausted, expressed as degree-day equivalent [°C·day] = leftover mm SWE / melt_factor. Zero when SWE is not fully depleted.

The melt factor (mm SWE °C⁻¹ day⁻¹) serves as the bridge between the PDD snowmelt representation and the frozen ground index (FGI): dividing excess melt depth (mm SWE) by melt_factor recovers the equivalent thermal forcing in °C·day, which is the currency the FGI uses. See Buckets._update_fgi().

Return type:

float

run_and_score

Run the model with a given parameter set and return goodness-of-fit metrics.

mnished.run_and_score(cfg, recession_coeff=None, f_to_discharge=None, Hmax=None, pdm_H0=None, f_tile=None, tau_tile=None, multipath_threshold=None, multipath_timescale=None, multipath_calibrated=0, leakance_R=None, leakance_R_calibrated=0, H_threshold=None, H_threshold_calibrated=0, melt_factor=None, fdd_threshold=None, snow_insulation_k=None, et_scale=None, et_alpha=None, wp_soil=None, wp_soil_sigma=None, recession_exponents=None, recession_exponents_calibrated=0, direct_runoff_fraction=None, baseflow_Q=None, modules=None, initial_states=None, post_spinup_states=None, post_spinup_k=0, start=None, end=None, spin_up_cycles=None, metric='KGE', routing_N=2, routing_K=None, enforce_water_balance='water-year')[source]

Run MNiShed and return a CalibResult named tuple.

Parameters:

  • cfg (str) – Path to a YAML configuration file. Should have spin_up_cycles: 0 so that spin-up is performed here with the calibrated parameters.

  • recession_coeff (list of float, optional) – Recession coefficients [days], one per reservoir (shallowest first). Overrides the values in cfg.

  • f_to_discharge (list of float, optional) – Exfiltration fractions to stream, one per reservoir except the deepest (which is always 1.0). Overrides the values in cfg. Used only for reservoirs whose junction_type is ‘fraction’ or ‘threshold’.

  • leakance_R (list of float or None, optional) – Leakance resistance [days] for each reservoir, one per reservoir. None entries leave the reservoir at its config-defined junction type. Non-None entries set junction_type to ‘leakance’ and assign the resistance. Q_leak = max(H_this - H_next, 0) / R. Default None (all fraction).

  • leakance_R_calibrated (int, optional) – Number of entries in leakance_R that are free calibration parameters (for AIC k-counting). Default 0.

  • H_threshold (list of float or None, optional) – Dead-storage threshold depth [mm] for each reservoir, one per reservoir. None entries leave the reservoir unchanged. Non-None entries set junction_type to ‘threshold’: only max(H - H_threshold, 0) drains. Models a stream-aquifer connection that activates above a threshold head. Default None (no threshold).

  • H_threshold_calibrated (int, optional) – Number of entries in H_threshold that are free calibration parameters (for AIC k-counting). Default 0.

  • multipath_threshold (list of float or None, optional) – Storage depth [mm] above which a parallel fast drain activates, per reservoir. None entries leave the reservoir’s multipath configuration unchanged from cfg. Non-None entries enable a threshold-activated parallel drain that adds max(0, H - thr)/multipath_timescale to discharge above the threshold. Distinct from f_tile/tau_tile (which is a constant-fraction bypass through a downstream sub-reservoir); see mnished.Reservoir for the contrast. Requires the matching multipath_timescale entry to be non-None. Default None.

  • multipath_timescale (list of float or None, optional) – E-folding timescale [days] of the parallel multipath drain, per reservoir. Required paired with multipath_threshold. Default None.

  • multipath_calibrated (int, optional) – Number of free parameters contributed by multipath_threshold + multipath_timescale combined (for AIC k-counting). Default 0.

  • Hmax (list of float, optional) – Maximum effective water depths [mm], one per reservoir. Overrides the values in cfg.

  • melt_factor (float, optional) – Degree-day snowmelt factor [mm SWE per degC per day]. Overrides the value in cfg.

  • direct_runoff_fraction (float or None, optional) – Fraction (0–1) of positive daily recharge that bypasses the reservoir cascade and exits directly as runoff. Conceptually inspired by Hortonian (infiltration-excess) overland flow, but at a daily timestep rainfall intensity is unavailable, so the fraction is not a rigorous physical representation – except in extreme events where intense rainfall dominates the daily total. In practice it acts as a calibrated fast-bypass fraction. None (default) leaves the value from the YAML config (itself defaulting to 0; off by default).

  • fdd_threshold (float or None, optional) – Frozen ground index threshold [°C·day]. The frozen ground index accumulates freezing degree-days and decays during warming (Molnau & Bissell 1983, https://westernsnowconference.org/bibliography/1983Molnau.pdf). When the index exceeds fdd_threshold, infiltration from the top reservoir to deeper layers is set to zero for that timestep (all drainage becomes direct runoff). None (default) disables the effect.

  • snow_insulation_k (float or None, optional) – Snow insulation decay constant [mm⁻¹ SWE]. Scales the effective air temperature reaching the soil as T_eff = T · exp(-k · SWE), reducing FGI accumulation under a deep snowpack. Applied to both freezing and thawing temperature forcing; excess degree-days from meltwater (excess_dd) are not scaled because meltwater delivers heat directly to the soil surface. None (default) leaves the value from cfg (itself defaulting to 0.0, i.e. no insulation). Literature values: LISFLOOD uses ~0.057 mm⁻¹; GSSHA default is 0.5 (units ambiguous in original source).

  • initial_states (dict, optional) –

    Starting reservoir water depths, snowpack SWE, and frozen ground index, as returned by a previous call’s CalibResult.final_states. Format:

    {'reservoirs': [H_shallow, H_deep, ...],
 'snowpack':    H_snow_SWE,
 'fgi':         frozen_ground_index}

    When provided, these override the H0 values from cfg. Use with spin_up_cycles=0 when chaining consecutive decades so that water storage and frozen-ground state are physically continuous. When None (default), reservoirs are initialised to their analytical steady-state depths before spin-up (see _steady_state_depths()).

  • post_spinup_states (dict, optional) – Reservoir water depths injected after spin-up completes and before the scored run begins. Same format as initial_states but only the 'reservoirs' key is used; individual entries may be None to leave that reservoir at its spin-up end state. Intended for calibrating decade-specific initial storage (e.g. log__H0_deep) when the spin-up equilibrium is poorly constrained by sparse pre-decade forcing. Only applied in decade mode (start is not None); ignored in full-record mode.

  • post_spinup_k (int, optional) – Number of free parameters contributing to post_spinup_states (for AIC counting). Default 0.

  • start (str or datetime-like, optional) – Start of the scoring window (inclusive). Score, AIC, BFI, and FDC are all computed within this window. Spin-up still uses the full record.

  • end (str or datetime-like, optional) – End of the scoring window (inclusive). Same as start.

  • spin_up_cycles (int or None, optional) – Number of times to loop the full record before the scored run. None (default) auto-computes as ceil(tau_max / record_length), where tau_max is the longest reservoir e-folding time after parameter overrides and record_length is the number of days in the input record. Because initial conditions are set to analytical steady-state depths, one e-folding time is sufficient to resolve the seasonal and inter-annual climate memory. Set to 0 when providing initial_states for chained decade runs.

  • metric ({'KGE', 'NSE', 'logKGE', 'KGE_logKGE', 'KGE_logKGE_logFDC',) – ‘KGE_logKGE_logFDC_BFI’}, optional Goodness-of-fit metric. Default is ‘KGE’. 'KGE_logKGE' returns 0.5*KGE + 0.5*logKGE, balancing peak and low-flow performance (Yilmaz et al. 2008). 'KGE_logKGE_logFDC_BFI' adds a BFI bias-ratio score (1 - abs(BFI_mod/BFI_obs - 1)) as a fourth equal-weight component.

  • routing_N (int, optional) – Number of identical linear reservoirs in the Nash cascade used for channel routing (shape parameter of the gamma IUH). Default is 2. Set routing_K to enable routing; routing_N is counted as a free parameter in k only when it is explicitly calibrated (the caller must add 1 to the k count via run_and_score if routing_N is varied).

  • routing_K (float or None, optional) – Storage time constant [days] of each Nash-cascade reservoir (scale parameter). Mean travel time through the cascade is routing_N * routing_K. When None (default), no routing is applied and the model output is compared directly to observed discharge. When provided, routing_K is counted as one free parameter.

  • modules (dict or None, optional) – Enable/disable optional process modules. Keys: 'snowpack', 'frozen_ground', 'rain_on_snow', 'direct_runoff'. Values: bool. Overrides the modules block in the YAML config. Absent keys leave the YAML/default value unchanged.

  • enforce_water_balance (str, optional) – ‘water-year’ (default): scale ET by a per-water-year multiplier so that P - Q - ET = 0 over each water year. ‘global’: scale ET by a single multiplier computed from sum(P - Q_obs) / sum(ET_raw) over the full record — no per-year overfitting, does not add hidden degrees of freedom to AIC. ‘none’: use raw ET without correction (only for measured ET). Legacy boolean True/’False are silently mapped to ‘water-year’/’none’. Overrides the enforce_water_balance key in the YAML config.

Returns:

Named tuple with fields score, aic, bfi_obs, bfi_mod, fdc_obs, fdc_mod, final_states, and buckets. See CalibResult for full field descriptions. All scalar fields are np.nan if the scoring window contains no valid overlapping data.

Return type:

CalibResult

Notes

When start is provided, the function operates in decade mode:

  • enforce_water_balance='global' closes the water balance over the full record (computed once at initialization), not re-fitted to the decade window. Re-fitting would assume ΔS = 0 over the decade, which is incorrect for transitional wet/dry decades.

  • Spin-up runs only the pre-decade data (record start through the day before start), so reservoir states at the beginning of the decade reflect pre-decade climatology rather than an arbitrary end-of-record state.

  • The scored run covers [start, end] only, continuing directly from the spin-up end state.

When start is None (default), the original full-record behaviour is preserved: spin-up and scoring both run the complete hydrodata record.

CalibResult

Named tuple returned by run_and_score().

class mnished.CalibResult(score, aic, bfi_obs, bfi_mod, kge_logfdc, fdc_obs, fdc_mod, final_states, buckets)

Named tuple returned by run_and_score().

score

Goodness-of-fit score (higher is better). Metric is one of KGE, NSE, or logKGE as requested; see run_and_score().

Type:

float

aic

Akaike Information Criterion computed on log-transformed flows (lower is better). Useful for comparing models with different numbers of free parameters.

Type:

float

bfi_obs

Baseflow index of the observed discharge, computed with the Eckhardt (2005) recursive digital filter.

Type:

float

bfi_mod

Baseflow index of the modelled discharge.

Type:

float

fdc_obs

Observed flow duration curve: discharge at 99 evenly-spaced exceedance probabilities (0.5–99.5 %), indexed by exceedance %.

Type:

pd.Series

fdc_mod

Modelled flow duration curve (same format as fdc_obs).

Type:

pd.Series

final_states

End-of-run reservoir states suitable for passing as initial_states to the next call:

{'reservoirs': [H_shallow, H_deep, ...],  # [mm]
 'snowpack':    H_snow_SWE,               # [mm SWE]
 'fgi':         frozen_ground_index}       # [°C·day]
Type:

dict

buckets

The Buckets object after the final run, including the full hydrodata DataFrame with modelled discharge.

Type:

Buckets

All scalar fields are ``np.nan`` if the scoring window contains no
valid overlapping data.

HydrographSeparation

Spectral and time-domain decomposition of a discharge record into reservoir-timescale components.

class mnished.HydrographSeparation(Q, dt=1.0, n_reservoirs=None, max_reservoirs=4, precip=None, recession_min_duration=10, recession_precip_threshold=0.5, recession_antecedent_days=3, tau_deep=None, tau_fast=None)[source]

Bases: object

Separate observed discharge into reservoir components and estimate timescales and initial storage depths.

Parameters:

  • Q (array-like) – Observed specific discharge [mm/day].

  • dt (float) – Timestep [days]. Default 1.0.

  • n_reservoirs (int or None) – Number of reservoirs. If None, selected by AIC on the spectral fit. The slowest reservoir (karst/deep) is always fitted by recession analysis regardless of this setting.

  • max_reservoirs (int) – Maximum number of reservoirs to test when n_reservoirs is None. Default 4.

  • precip (array-like or None) – Daily precipitation [mm/day], same length as Q. Used to mask recession periods; if None, only discharge decline is used.

  • recession_min_duration (int) – Minimum qualifying recession length [days]. Default 10.

  • recession_precip_threshold (float) – Maximum daily precipitation within a recession segment [mm/day]. Default 0.5.

  • recession_antecedent_days (int) – Number of days immediately before a recession that must also be dry. Default 3.

  • tau_deep (float or None) – Externally supplied timescale [days] for the deep (karst) reservoir. When provided, recession fitting is skipped for τ_karst and this value is used directly. Recommended when τ_karst >> T_record, since a decay too slow to observe cannot be reliably estimated from the discharge record alone.

  • tau_fast (array-like or None) – Externally supplied timescales [days] for the fast reservoirs (all except the deepest), e.g. from a prior calibration run. When provided, spectral fitting is skipped. Note: spectral corner periods are 2πτ, not τ — a 29-day reservoir has its -3 dB corner at 182 days, and a 466-day soil reservoir at ~8 years, which may lie near or beyond the record length and be unreliable from spectral fitting alone.

  • fit()) (Attributes (populated after)

  • ------------------------------------

  • tau (ndarray) – Reservoir timescales [days], ordered slowest to fastest.

  • H0 (ndarray) – Initial storage depths [mm], same order as tau.

  • n_reservoirs_fitted (int) – Number of reservoirs selected (spectral fast + 1 karst).

  • aic_scores (list of float) – AIC values for k = 1 … max_reservoirs-1 spectral fits.

  • tau_karst (float) – Deep-reservoir timescale [days] from recession fitting.

fit()[source]

Run the full separation pipeline. Returns self.

get_initial_conditions()[source]

Initial reservoir storage depths, ordered fastest to slowest (shallow first, karst last), matching the MNiShed reservoir index convention.

Returns:

{'H0': [H0_shallow, H0_soil, ..., H0_karst]} in mm.

Return type:

dict

get_parameter_priors()[source]

Suggested initial values and bounds (in params.yml units) for the MNiShed log__recession_coeff_* calibration parameters.

Bounds are ±0.5 in log10 around the estimated timescale — wide enough to allow the optimizer freedom, tight enough to stay physical.

Returns:

Keyed by MNiShed parameter name. Each value is a dict with 'initial', 'lower', 'upper' in log10(days).

Return type:

dict

summary()[source]

Print a brief summary of fitted timescales and initial conditions.