Skip to content

scenario_daily

apply_dry_day_rule(anomaly: xr.Dataset, target: xr.Dataset, dry_day_rule: str) -> xr.Dataset

Stop the eps offset from manufacturing rain on days the model reports as dry.

The multiplicative anomaly (T + eps)/(R + eps) is strictly positive even when the model reports no rain at all, so a rainless model day still receives E_ref(month) * a(d) > 0 of the ERA5 climatology. Worse, eps dominates the numerator for such a day, so every rainless day in a cell-month gets the identical anomaly 1/(R_m + eps) -- a flat positive floor rather than a dry spell. Wherever that floor clears the 0.1 mm/day wet-day threshold the pipeline reports a wet day that neither ERA5 nor the GCM has, which is what makes precipitation_days step at the boundary.

preserve zeroes the anomaly on those days and rescales the surviving days of the same cell-month so the month's summed anomaly is unchanged. Because that sum is preserved per GCM cell and interpolate_to_target_latlon is linear, the monthly -- and therefore the annual -- total on the target grid is untouched to floating point. Only the distribution across days moves. That is deliberate: this is a shape fix, exactly orthogonal to the Jensen de-bias, which is a level fix that leaves the shape alone. The two commute, because the de-bias scales a whole month uniformly and this rescale is scale-invariant.

A cell-month the model reports dry on every day has nothing to renormalise onto. Those are left exactly as they are rather than zeroed. Zeroing them is the variant that was measured and rejected: it loses up to 1.2% of the population-weighted annual total across the 1-3% of cell-months that are all-dry. Keeping them is what makes this rule total-preserving, and it is the whole difference between the two.

Source code in src/climate_data/generate/scenario_daily.py
def apply_dry_day_rule(
    anomaly: xr.Dataset,
    target: xr.Dataset,
    dry_day_rule: str,
) -> xr.Dataset:
    """Stop the ``eps`` offset from manufacturing rain on days the model reports as dry.

    The multiplicative anomaly ``(T + eps)/(R + eps)`` is strictly positive even when the
    model reports no rain at all, so a rainless model day still receives
    ``E_ref(month) * a(d) > 0`` of the ERA5 climatology. Worse, ``eps`` dominates the
    numerator for such a day, so *every* rainless day in a cell-month gets the identical
    anomaly ``1/(R_m + eps)`` -- a flat positive floor rather than a dry spell. Wherever that
    floor clears the 0.1 mm/day wet-day threshold the pipeline reports a wet day that neither
    ERA5 nor the GCM has, which is what makes ``precipitation_days`` step at the boundary.

    ``preserve`` zeroes the anomaly on those days and rescales the surviving days of the same
    cell-month so the month's *summed* anomaly is unchanged. Because that sum is preserved per
    GCM cell and ``interpolate_to_target_latlon`` is linear, the monthly -- and therefore the
    annual -- total on the target grid is untouched to floating point. Only the distribution
    across days moves. That is deliberate: this is a shape fix, exactly orthogonal to the
    Jensen de-bias, which is a level fix that leaves the shape alone. The two commute, because
    the de-bias scales a whole month uniformly and this rescale is scale-invariant.

    A cell-month the model reports dry on *every* day has nothing to renormalise onto. Those
    are left exactly as they are rather than zeroed. Zeroing them is the variant that was
    measured and rejected: it loses up to 1.2% of the population-weighted annual total across
    the 1-3% of cell-months that are all-dry. Keeping them is what makes this rule
    total-preserving, and it is the whole difference between the two.
    """
    if dry_day_rule == "none":
        return anomaly
    if dry_day_rule != "preserve":
        msg = f"Unknown dry-day rule: {dry_day_rule}"
        raise ValueError(msg)

    wet = target > cdc.DRY_DAY_THRESHOLD_MM
    kept = anomaly.where(wet, 0.0)

    month_total = anomaly.groupby("date.month").sum("date")
    kept_total = kept.groupby("date.month").sum("date")
    has_wet_day = kept_total > 0
    # 1.0 on an all-dry cell-month, so the restore below survives the rescale untouched.
    rescale = (month_total / kept_total.where(has_wet_day)).fillna(1.0)

    has_wet_day_daily = has_wet_day.sel(month=anomaly["date"].dt.month).drop_vars(
        "month"
    )
    kept = kept.where(has_wet_day_daily, anomaly)

    rescaled = kept.groupby("date.month") * rescale
    return rescaled.drop_vars("month")

apply_value_bounds(ds: xr.Dataset, bounds: tuple[float, float]) -> xr.Dataset

Clip every data variable to bounds (inclusive); NaN stays NaN.

Source code in src/climate_data/generate/scenario_daily.py
def apply_value_bounds(ds: xr.Dataset, bounds: tuple[float, float]) -> xr.Dataset:
    """Clip every data variable to ``bounds`` (inclusive); NaN stays NaN."""
    lo, hi = bounds
    if not lo < hi:
        msg = f"value bounds must satisfy lo < hi, got {bounds!r}."
        raise ValueError(msg)
    return ds.clip(lo, hi)

check_bounded_variable_scheme(target_variable: str, anomaly_scheme: str) -> None

Refuse a stabilised or yearly anomaly scheme for a value-bounded variable.

A bounded variable (relative humidity) has its raw inputs clipped to cdc.VALUE_BOUNDS before the ratio is formed, so the floor already does the job of the +1 / tapered eps in the eps-bearing schemes -- applying both would stabilise twice and shift the ratio away from the construction that was evaluated. The yearly family drops the monthly anchor that evaluation was made on. Called from the worker and, via variables_for_anomaly_scheme, from the launcher, so --target-variable all under the default scheme skips the variable with a message instead of queueing doomed jobs.

Source code in src/climate_data/generate/scenario_daily.py
def check_bounded_variable_scheme(target_variable: str, anomaly_scheme: str) -> None:
    """Refuse a stabilised or yearly anomaly scheme for a value-bounded variable.

    A bounded variable (relative humidity) has its raw inputs clipped to
    ``cdc.VALUE_BOUNDS`` before the ratio is formed, so the floor already does the job of
    the +1 / tapered eps in the eps-bearing schemes -- applying both would stabilise twice
    and shift the ratio away from the construction that was evaluated. The yearly family
    drops the monthly anchor that evaluation was made on. Called from the worker and, via
    ``variables_for_anomaly_scheme``, from the launcher, so ``--target-variable all`` under
    the default scheme skips the variable with a message instead of queueing doomed jobs.
    """
    bounds = cdc.VALUE_BOUNDS.get(target_variable)
    if bounds is None or anomaly_scheme in cdc.BOUNDED_VARIABLE_SCHEMES:
        return
    msg = (
        f"{target_variable!r} is value-bounded (inputs clipped to {bounds['input']}, output "
        f"to {bounds['output']}) and runs under {list(cdc.BOUNDED_VARIABLE_SCHEMES)} only; "
        f"got anomaly_scheme={anomaly_scheme!r}. Pass --anomaly-scheme "
        f"{cdc.BOUNDED_VARIABLE_SCHEMES[0]} for this variable."
    )
    raise ValueError(msg)

check_debias_variable(target_variable: str, debias_method: str, dry_day_rule: str = 'none') -> None

Refuse a correction for a variable it has not been validated against.

Called from the launchers as well as from the worker, so that --target-variable all fails in a second rather than after submitting thousands of doomed jobs.

Source code in src/climate_data/generate/scenario_daily.py
def check_debias_variable(
    target_variable: str, debias_method: str, dry_day_rule: str = "none"
) -> None:
    """Refuse a correction for a variable it has not been validated against.

    Called from the launchers as well as from the worker, so that ``--target-variable all``
    fails in a second rather than after submitting thousands of doomed jobs.
    """
    if debias_method != "none" and target_variable not in cdc.DEBIAS_VARIABLES:
        msg = (
            f"debias_method={debias_method!r} is not validated for {target_variable!r}. "
            f"Allowed: {list(cdc.DEBIAS_VARIABLES)}. Name the variable explicitly rather "
            "than using 'all'."
        )
        raise ValueError(msg)
    if dry_day_rule != "none" and target_variable not in cdc.DRY_DAY_VARIABLES:
        msg = (
            f"dry_day_rule={dry_day_rule!r} is not validated for {target_variable!r}. "
            f"Allowed: {list(cdc.DRY_DAY_VARIABLES)}. Name the variable explicitly rather "
            "than using 'all'."
        )
        raise ValueError(msg)

check_scheme_compatibility(anomaly_scheme: str, anomaly_type: str, debias_method: str, dry_day_rule: str) -> None

Reject combinations of the two correction axes that cannot both apply.

debias_method and dry_day_rule correct the (T + eps) / (R + eps) construction. Both monthly (constant eps) and monthly-taper (tapered eps) use it, so both accept them. The yearly and monthly-ratio families have no eps for those corrections to act on, so asking for either against them is a mistake rather than a no-op -- and silently ignoring the request would produce a file whose attrs claim a correction that was never applied.

Source code in src/climate_data/generate/scenario_daily.py
def check_scheme_compatibility(
    anomaly_scheme: str,
    anomaly_type: str,
    debias_method: str,
    dry_day_rule: str,
) -> None:
    """Reject combinations of the two correction axes that cannot both apply.

    `debias_method` and `dry_day_rule` correct the `(T + eps) / (R + eps)` construction.
    Both `monthly` (constant eps) and `monthly-taper` (tapered eps) use it, so both accept
    them. The yearly and monthly-ratio families have no eps for those corrections to act
    on, so asking for either against them is a mistake rather than a no-op -- and silently
    ignoring the request would produce a file whose attrs claim a correction that was
    never applied.
    """
    has_eps = anomaly_scheme in cdc.EPS_BEARING_SCHEMES
    if not has_eps and (debias_method != "none" or dry_day_rule != "none"):
        msg = (
            f"debias_method={debias_method!r} and dry_day_rule={dry_day_rule!r} cannot "
            f"be combined with anomaly_scheme={anomaly_scheme!r}: they correct the eps "
            f"stabiliser, which this scheme does not use. Schemes that carry an eps: "
            f"{list(cdc.EPS_BEARING_SCHEMES)}."
        )
        raise ValueError(msg)
    if anomaly_type != "multiplicative":
        msg = f"Anomaly scheme '{anomaly_scheme}' only applies to multiplicative variables."
        raise ValueError(msg)

compute_anomaly(reference: xr.Dataset, target: xr.Dataset, anomaly_type: str, *, debias_method: str, dry_day_rule: str, anomaly_scheme: str = cdc.ANOMALY_SCHEME_MONTHLY, eps_floor: float = cdc.DEFAULT_EPS_FLOOR, anomaly_cap: float | None = cdc.DEFAULT_ANOMALY_CAP) -> xr.Dataset

The forecast anomaly, optionally capped.

A thin wrapper so the ceiling applies to every scheme without threading it through each one: the scheme dispatch below has several return points, and duplicating the clip at each is how one of them ends up missing it.

Applied on the GCM grid, BEFORE interpolate_to_target_latlon. Capping afterwards would leave a blown-up cell already smeared across its neighbours by the regrid.

The cap bounds each cell's multiplier at the granularity the scheme anchors at -- annual for the yearly family, per calendar month for the monthly family -- and rescales the daily series to meet it.

It cannot be an elementwise clip. anomaly carries the target's daily date dimension, because the denominator is a reference-window mean with no date dim and the division broadcasts. Clipping elementwise bounds every day at anomaly_cap times the cell's reference mean, and against an annual denominator that ceiling sits inside the ordinary distribution of daily rainfall: a cell averaging 2 mm/day is cut at 40 mm/day, an unremarkable tropical wet day. The daily-over-annual ratio also carries the seasonal cycle, so such a clip bites hardest where seasonality is strongest rather than where the pathology is.

Measured on yearly-delta ssp126 before this change, an elementwise clip at 20 altered 20.5% of land pixels in 2083 and 15.5% in 2050; 82% of those had an annual anomaly at or below 2 -- cells projecting essentially no change -- and the global total moved -19.1% and -1.1% respectively.

Rescaling instead leaves any cell at or below the ceiling bit-identical to uncapped, brings a cell above it exactly to the ceiling, and preserves within-period shape so a wet day stays proportionally a wet day.

The cap still only removes precipitation -- it cannot conserve, and whatever it removes is gone from the total. That is intended: the values it removes are ones no cell's observed climatology supports. But it means a cap moves the level, so it must be chosen on evidence rather than set defensively.

Source code in src/climate_data/generate/scenario_daily.py
def compute_anomaly(
    reference: xr.Dataset,
    target: xr.Dataset,
    anomaly_type: str,
    *,
    debias_method: str,
    dry_day_rule: str,
    anomaly_scheme: str = cdc.ANOMALY_SCHEME_MONTHLY,
    eps_floor: float = cdc.DEFAULT_EPS_FLOOR,
    anomaly_cap: float | None = cdc.DEFAULT_ANOMALY_CAP,
) -> xr.Dataset:
    """The forecast anomaly, optionally capped.

    A thin wrapper so the ceiling applies to every scheme without threading it through
    each one: the scheme dispatch below has several return points, and duplicating the
    clip at each is how one of them ends up missing it.

    Applied on the GCM grid, BEFORE `interpolate_to_target_latlon`. Capping afterwards
    would leave a blown-up cell already smeared across its neighbours by the regrid.

    The cap bounds each cell's multiplier at the granularity the scheme anchors at --
    annual for the yearly family, per calendar month for the monthly family -- and
    rescales the daily series to meet it.

    It cannot be an elementwise clip. `anomaly` carries the target's daily date
    dimension, because the denominator is a reference-window mean with no date dim and
    the division broadcasts. Clipping elementwise bounds every *day* at `anomaly_cap`
    times the cell's reference mean, and against an annual denominator that ceiling sits
    inside the ordinary distribution of daily rainfall: a cell averaging 2 mm/day is cut
    at 40 mm/day, an unremarkable tropical wet day. The daily-over-annual ratio also
    carries the seasonal cycle, so such a clip bites hardest where seasonality is
    strongest rather than where the pathology is.

    Measured on `yearly-delta` ssp126 before this change, an elementwise clip at 20
    altered 20.5% of land pixels in 2083 and 15.5% in 2050; 82% of those had an annual
    anomaly at or below 2 -- cells projecting essentially no change -- and the global
    total moved -19.1% and -1.1% respectively.

    Rescaling instead leaves any cell at or below the ceiling bit-identical to uncapped,
    brings a cell above it exactly to the ceiling, and preserves within-period shape so
    a wet day stays proportionally a wet day.

    The cap still only removes precipitation -- it cannot conserve, and whatever it
    removes is gone from the total. That is intended: the values it removes are ones no
    cell's observed climatology supports. But it means a cap moves the level, so it must
    be chosen on evidence rather than set defensively.
    """
    anomaly = _compute_anomaly_uncapped(
        reference,
        target,
        anomaly_type,
        debias_method=debias_method,
        dry_day_rule=dry_day_rule,
        anomaly_scheme=anomaly_scheme,
        eps_floor=eps_floor,
    )
    if anomaly_cap is None:
        return anomaly
    if anomaly_type != "multiplicative":
        msg = (
            f"anomaly_cap={anomaly_cap!r} was requested for an additive anomaly. The cap "
            "bounds a multiplier; an additive anomaly is a difference in the variable's "
            "own units, where a ceiling of this kind is meaningless."
        )
        raise ValueError(msg)
    if anomaly_cap <= 0:
        msg = f"anomaly_cap must be positive, got {anomaly_cap!r}."
        raise ValueError(msg)
    if anomaly_scheme in cdc.YEARLY_ANOMALY_SCHEMES:
        # anchored to the annual level, so bound the annual multiplier
        over = anomaly.mean("date")
        factor = (anomaly_cap / over).where(over > anomaly_cap, 1.0)
        return anomaly * factor

    # the monthly family anchors each calendar month to ERA5 separately, so the ceiling
    # belongs per month. Bounding the annual mean instead would let one blown-up month
    # drag the rescale factor down and crush the other eleven, which are fine.
    over = anomaly.groupby("date.month").mean("date")
    factor = (anomaly_cap / over).where(over > anomaly_cap, 1.0)
    return anomaly.groupby("date.month") * factor

jensen_debias_factor(reference: xr.Dataset, reference_monthly: xr.Dataset, debias_method: str, eps: xr.Dataset | float = 1.0) -> xr.Dataset

The factor to divide a multiplicative anomaly by, per month and per GCM cell.

The anomaly is (T + 1) / (R + 1) with R a monthly mean over only five reference years. 1/(R + 1) is convex, so by Jensen's inequality the anomaly averages above 1 even when the target year is drawn from the same distribution as the reference period -- a level bias on every forecast year. This returns an estimate of that inflation.

loo -- leave-one-out. For each held-out reference year, form the multiplier of that year against the mean of the other years and average over the folds. Between reference years there is no climate signal, so an unbiased estimator would return 1; the excess is the bias, measured from the data with no series expansion. The held-out denominator averages n-1 years while the pipeline averages n, and the bias goes as 1/n, so the excess is rescaled by (n-1)/n.

This is provably >= 1: with u_y = T_y + 1 and S = sum_y u_y the held-out denominator is (S - u_y)/(n-1), so each fold is (n-1)*u_y/(S - u_y), convex in u_y, and Jensen gives mean_y >= f(S/n) = 1 with equality iff every reference year is identical. So dividing by it can only shrink the anomaly, never inflate it -- which is what makes the effect on a threshold count such as precipitation_days sign-definite.

analytic -- the second-order expansion 1 + Var(Rbar)/(R + 1)^2. Cheaper to reason about but it is a truncated series, and the neglected terms matter exactly where the correction is largest (near-zero R, where eps dominates the denominator). Kept for comparison; loo is the estimator this was built for.

reference_monthly is the pipeline's own denominator, passed in rather than recomputed so the analytic form squares precisely the value the anomaly divides by.

Source code in src/climate_data/generate/scenario_daily.py
def jensen_debias_factor(
    reference: xr.Dataset,
    reference_monthly: xr.Dataset,
    debias_method: str,
    eps: xr.Dataset | float = 1.0,
) -> xr.Dataset:
    """The factor to divide a multiplicative anomaly by, per month and per GCM cell.

    The anomaly is ``(T + 1) / (R + 1)`` with ``R`` a monthly mean over only five reference
    years. ``1/(R + 1)`` is convex, so by Jensen's inequality the anomaly averages above 1 even
    when the target year is drawn from the same distribution as the reference period -- a level
    bias on every forecast year. This returns an estimate of that inflation.

    ``loo`` -- leave-one-out. For each held-out reference year, form the multiplier of that
    year against the mean of the *other* years and average over the folds. Between reference
    years there is no climate signal, so an unbiased estimator would return 1; the excess is
    the bias, measured from the data with no series expansion. The held-out denominator
    averages ``n-1`` years while the pipeline averages ``n``, and the bias goes as ``1/n``, so
    the excess is rescaled by ``(n-1)/n``.

    This is provably ``>= 1``: with ``u_y = T_y + 1`` and ``S = sum_y u_y`` the held-out
    denominator is ``(S - u_y)/(n-1)``, so each fold is ``(n-1)*u_y/(S - u_y)``, convex in
    ``u_y``, and Jensen gives ``mean_y >= f(S/n) = 1`` with equality iff every reference year
    is identical. So dividing by it can only shrink the anomaly, never inflate it -- which is
    what makes the effect on a threshold count such as ``precipitation_days`` sign-definite.

    ``analytic`` -- the second-order expansion ``1 + Var(Rbar)/(R + 1)^2``. Cheaper to reason
    about but it is a truncated series, and the neglected terms matter exactly where the
    correction is largest (near-zero ``R``, where ``eps`` dominates the denominator). Kept for
    comparison; ``loo`` is the estimator this was built for.

    ``reference_monthly`` is the pipeline's own denominator, passed in rather than recomputed so
    the analytic form squares precisely the value the anomaly divides by.
    """
    by_year = _monthly_means_by_reference_year(reference)
    n_years = by_year.sizes["reference_year"]

    if debias_method == "loo":
        mean_year = by_year.mean("reference_year")
        folds = []
        valid = []
        for i in range(n_years):
            held_out = by_year.isel(reference_year=i)
            others = (n_years * mean_year - held_out) / (n_years - 1)
            fold_denominator = others + eps
            folds.append(
                (held_out + eps) / fold_denominator.where(fold_denominator > 0)
            )
            valid.append((fold_denominator > 0).astype("int8"))
        raw = xr.concat(folds, dim="reference_year").mean("reference_year")
        factor = 1.0 + ((n_years - 1) / n_years) * (raw - 1.0)
        # Where any fold is undefined the estimate is undefined, so apply no correction.
        # Averaging the surviving folds is NOT a repair: the undefined one is precisely
        # the large fold that lifts the mean above 1, so dropping it collapses the factor
        # (0.2 for a one-wet-four-dry cell-month) and the de-bias INFLATES the anomaly --
        # the opposite of its purpose. Only reachable when eps can be zero, i.e. under the
        # taper; the constant eps = 1 floors every denominator at 1.
        factor = factor.where(sum(valid) == n_years, 1.0).fillna(1.0)
    elif debias_method == "analytic":
        variance = by_year.var("reference_year", ddof=1)
        factor = 1.0 + (variance / n_years) / (reference_monthly + eps) ** 2
    else:
        msg = f"Unknown debias method: {debias_method}"
        raise ValueError(msg)

    factor = factor.drop_vars("reference_year", errors="ignore")
    values = factor.to_dataarray()
    if not bool(np.isfinite(values).all()):
        msg = (
            f"Non-finite value in the {debias_method} de-bias factor. Interpolation would "
            "silently fill it from a neighbour rather than surface it, so refusing to proceed."
        )
        raise ValueError(msg)
    minimum = float(values.min())
    if minimum < 1.0 - _DEBIAS_FACTOR_TOLERANCE:
        # Both estimators are >= 1 by construction -- loo by Jensen on a convex fold,
        # analytic because it adds a non-negative term. A value below 1 means the
        # construction's precondition failed, and dividing by it would inflate rather
        # than shrink the anomaly. Refuse rather than ship an anti-correction.
        msg = (
            f"The {debias_method} de-bias factor reached {minimum!r}, below its "
            "guaranteed lower bound of 1. Dividing by it would inflate the anomaly."
        )
        raise ValueError(msg)
    return factor

load_and_shift_longitude_and_correct_time(member_path: str | Path, year: str) -> xr.Dataset

Put a member's year onto the real Gregorian calendar, day by day.

The conversion is by DATE, never by value. interp_calendar used to resample onto the target axis by linear interpolation, so whenever the source calendar's year length differed from the target's -- a noleap member in a leap year -- every target day became a blend of two source days. convert_calendar maps each source day onto its own date instead and leaves 29 February missing; the reindex holds the output axis to exactly this year's days whatever the source calendar, and interpolate_na fills the gap from the nearest real day.

Every variable comes through here, not just precipitation, but what the blending cost depends on the field. Temperature is smooth, so blending barely moves it and the threshold measures it feeds -- days_over_30C, the suitability maps -- were knocked either way and largely cancelled in the spatial mean. Precipitation is spiky and mostly zero against a 0.1 mm cut sitting on the floor, so smearing a wet day onto its dry neighbours could only push them UP over the line, never below it. That one-way ratchet is why precipitation_days took a coherent ~13.5 d per noleap member in all 19 leap years 2024-2096 while everything else came out as cancelling noise. (CLIMATE-35)

No align_on is passed: it only takes effect when a 360_day calendar is involved, and no member in the extracts uses one. If that changes, align_on needs a deliberate choice rather than xarray's default.

Source code in src/climate_data/generate/scenario_daily.py
def load_and_shift_longitude_and_correct_time(
    member_path: str | Path,
    year: str,
) -> xr.Dataset:
    """Put a member's year onto the real Gregorian calendar, day by day.

    The conversion is by DATE, never by value. `interp_calendar` used to resample onto the
    target axis by linear interpolation, so whenever the source calendar's year length
    differed from the target's -- a `noleap` member in a leap year -- every target day
    became a blend of two source days. `convert_calendar` maps each source day onto its own
    date instead and leaves 29 February missing; the `reindex` holds the output axis to
    exactly this year's days whatever the source calendar, and `interpolate_na` fills the
    gap from the nearest real day.

    Every variable comes through here, not just precipitation, but what the blending cost
    depends on the field. Temperature is smooth, so blending barely moves it and the
    threshold measures it feeds -- `days_over_30C`, the suitability maps -- were knocked
    either way and largely cancelled in the spatial mean. Precipitation is spiky and mostly
    zero against a 0.1 mm cut sitting on the floor, so smearing a wet day onto its dry
    neighbours could only push them UP over the line, never below it. That one-way ratchet
    is why `precipitation_days` took a coherent ~13.5 d per noleap member in all 19 leap
    years 2024-2096 while everything else came out as cancelling noise. (CLIMATE-35)

    No `align_on` is passed: it only takes effect when a `360_day` calendar is involved,
    and no member in the extracts uses one. If that changes, `align_on` needs a deliberate
    choice rather than xarray's default.
    """
    time_slice = slice(f"{year}-01-01", f"{year}-12-31")
    time_range = pd.date_range(f"{year}-01-01", f"{year}-12-31")
    ds = load_and_shift_longitude(member_path, time_slice)
    ds = (
        ds.assign_coords(time=ds.time.dt.floor("D"))
        .convert_calendar("standard")
        .reindex(time=time_range)
        .interpolate_na(dim="time", method="nearest", fill_value="extrapolate")
        .rename({"time": "date"})
    )
    return ds

variables_for_anomaly_scheme(target_variables: list[str], anomaly_scheme: str, anomaly_types: dict[str, str]) -> list[str]

Drop target variables the anomaly scheme cannot be applied to.

The yearly schemes are defined for multiplicative variables only -- compute_anomaly raises for anything else -- while --target-variable defaults to ALL. Without this filter the default invocation of a yearly scheme submits additive jobs that are certain to fail, and they fail only after being scheduled and retried.

Skipped variables are named rather than dropped quietly, and selecting nothing but additive variables is an error rather than an empty run.

Source code in src/climate_data/generate/scenario_daily.py
def variables_for_anomaly_scheme(
    target_variables: list[str],
    anomaly_scheme: str,
    anomaly_types: dict[str, str],
) -> list[str]:
    """Drop target variables the anomaly scheme cannot be applied to.

    The yearly schemes are defined for multiplicative variables only -- `compute_anomaly`
    raises for anything else -- while `--target-variable` defaults to ALL. Without this
    filter the default invocation of a yearly scheme submits additive jobs that are
    certain to fail, and they fail only after being scheduled and retried.

    Skipped variables are named rather than dropped quietly, and selecting nothing but
    additive variables is an error rather than an empty run.
    """
    keep = []
    skip = []
    skip_bounded = []
    for variable in target_variables:
        if (
            variable in cdc.VALUE_BOUNDS
            and anomaly_scheme not in cdc.BOUNDED_VARIABLE_SCHEMES
        ):
            # Value-bounded variables are floored explicitly; see
            # check_bounded_variable_scheme for why the other schemes are refused.
            skip_bounded.append(variable)
        elif (
            anomaly_scheme != cdc.ANOMALY_SCHEME_MONTHLY
            and anomaly_types[variable] != "multiplicative"
        ):
            skip.append(variable)
        else:
            keep.append(variable)

    if skip:
        print(
            f"Anomaly scheme '{anomaly_scheme}' applies to multiplicative variables only;"
            f" skipping {len(skip)}: {', '.join(skip)}."
        )
    if skip_bounded:
        print(
            f"Anomaly scheme '{anomaly_scheme}' does not apply to value-bounded variables,"
            f" which run under {', '.join(cdc.BOUNDED_VARIABLE_SCHEMES)} only;"
            f" skipping {len(skip_bounded)}: {', '.join(skip_bounded)}."
        )
    if not keep:
        msg = (
            f"No multiplicative variables selected, so anomaly scheme '{anomaly_scheme}'"
            f" has nothing to run. Selected: {', '.join(target_variables)}."
        )
        raise click.UsageError(msg)
    return keep