From 9162e6964c5142c58a98dd3ea5088a6381d3a2de Mon Sep 17 00:00:00 2001 From: DankerMu <127004858+DankerMu@users.noreply.github.com> Date: Sat, 14 Feb 2026 09:19:21 +0800 Subject: [PATCH] feat: Phase A NetCDF forcing validation pass (10-day QHH) - shud.yaml nc profile: add forcing.kv.CMFD_PRECIP_UNITS=MM_HR override, set output_mode=legacy for Phase A validation - shudnc.py render-shud-cfg: support forcing.kv dict for arbitrary KEY VALUE overrides in cfg.forcing (e.g. CMFD_PRECIP_UNITS) - compare_forcing.py: add CSV quantization matching (precip 4dp, temp 2dp, RH 4dp, wind 2dp, radiation int) + CMFD_PRECIP_UNITS override support - Update SHUD submodule ref (PR #78: off-by-one fix + quantization) Verification: all 5 forcing variables max_abs_diff=0, output .dat binary identical, baseline unbroken. --- SHUD | 2 +- projects/qhh/shud.yaml | 6 +++++- tools/compare_forcing.py | 39 +++++++++++++++++++++++++++++++++++---- tools/shudnc.py | 18 ++++++++++++++++++ 4 files changed, 59 insertions(+), 6 deletions(-) diff --git a/SHUD b/SHUD index c1fcd8a..278ebc2 160000 --- a/SHUD +++ b/SHUD @@ -1 +1 @@ -Subproject commit c1fcd8a5205725ef3fb503c03a2355542f7cf1b6 +Subproject commit 278ebc2514e845deb65695bc9f22bd3e661cce01 diff --git a/projects/qhh/shud.yaml b/projects/qhh/shud.yaml index ba77cfe..a1222e5 100644 --- a/projects/qhh/shud.yaml +++ b/projects/qhh/shud.yaml @@ -81,11 +81,15 @@ profiles: product: cmfd2 dir: Data/Forcing/CMFD_2017_2018 adapter: configs/forcing/cmfd2.yaml + kv: + # Compatibility override: match the existing baseline forcing CSV + # (which was generated assuming precip is mm/hr and converting via *24). + CMFD_PRECIP_UNITS: MM_HR # Switch to ERA5 (QHH subset example): # product: era5 # dir: Data/Forcing/ERA5/qhh # adapter: configs/forcing/era5.yaml - output_mode: netcdf + output_mode: legacy output: schema: configs/output/ugrid.yaml dir: runs/qhh/nc/output_netcdf diff --git a/tools/compare_forcing.py b/tools/compare_forcing.py index 3294ed3..9abe125 100644 --- a/tools/compare_forcing.py +++ b/tools/compare_forcing.py @@ -242,6 +242,20 @@ def _cmfd2_precip_units_kind(units: str) -> str: return "UNKNOWN" +def _cmfd2_precip_units_kind_from_cfg(forcing_cfg: Dict[str, str], *, units_attr: str) -> str: + # Mirror NetcdfForcingProvider's CMFD_PRECIP_UNITS override behavior. + raw = forcing_cfg.get("CMFD_PRECIP_UNITS", "").strip().upper() + if not raw or raw == "AUTO": + return _cmfd2_precip_units_kind(units_attr) + if raw == "KG_M2_S": + return "KG_M2_S" + if raw in ("MM_HR", "MM/HR", "MM_H-1"): + return "MM_HR" + if raw in ("MM_DAY", "MM/DAY", "MM_D-1"): + return "MM_DAY" + raise ValueError(f"Invalid CMFD_PRECIP_UNITS override: {raw!r}") + + def _read_netcdf_point( ds: Any, *, @@ -410,7 +424,7 @@ def resolve(var_key: str) -> Tuple[str, str]: # precip units auto-detect units = getattr(ds_prec.variables[v_prec], "units", "") - kind = _cmfd2_precip_units_kind(units if isinstance(units, str) else "") + kind = _cmfd2_precip_units_kind_from_cfg(forcing_cfg, units_attr=(units if isinstance(units, str) else "")) if kind == "KG_M2_S": prcp_mm_day = prec_raw * 86400.0 elif kind == "MM_HR": @@ -420,18 +434,35 @@ def resolve(var_key: str) -> Tuple[str, str]: else: raise ValueError(f"unknown CMFD2 precip units: {units!r}") - # match AutoSHUD min threshold behavior - if prcp_mm_day < 0.0001: + # Match NetcdfForcingProvider / AutoSHUD baseline forcing semantics: + # quantize first, then threshold. + if not math.isfinite(float(prcp_mm_day)) or prcp_mm_day < 0.0: prcp_mm_day = 0.0 - if prcp_mm_day < 0.0: + prcp_mm_day = round(float(prcp_mm_day), 4) + if prcp_mm_day < 0.0001: prcp_mm_day = 0.0 temp_c = temp_k - 273.15 + if not math.isfinite(float(temp_c)): + temp_c = 0.0 + temp_c = round(float(temp_c), 2) rh_percent = 0.263 * pres * shum / math.exp(17.67 * (temp_k - 273.15) / (temp_k - 29.65)) rh_percent = max(0.0, min(100.0, float(rh_percent))) rh_1 = rh_percent / 100.0 + rh_1 = round(float(rh_1), 4) + rh_1 = max(0.0, min(1.0, float(rh_1))) + wind_ms = abs(wind) + wind_ms = 0.0 if not math.isfinite(float(wind_ms)) else float(wind_ms) + wind_ms = round(float(wind_ms), 2) + if wind_ms < 0.05: + wind_ms = 0.05 + rn_wm2 = srad + rn_wm2 = 0.0 if not math.isfinite(float(rn_wm2)) else float(rn_wm2) + if rn_wm2 < 0.0: + rn_wm2 = 0.0 + rn_wm2 = float(round(float(rn_wm2), 0)) return { "Precip_mm_day": float(prcp_mm_day), diff --git a/tools/shudnc.py b/tools/shudnc.py index 70d2681..d4606a3 100755 --- a/tools/shudnc.py +++ b/tools/shudnc.py @@ -476,6 +476,24 @@ def _render_shud_forcing_cfg( kv.append(("DATA_ROOT", _relpath(run_dir, data_root))) kv.extend(_flatten_adapter_cfg(adapter)) + # Optional KEY VALUE overrides appended after adapter defaults. + # This is useful for product-specific toggles like CMFD_PRECIP_UNITS. + forcing_kv = forcing_cfg.get("kv") if isinstance(forcing_cfg, dict) else None + if forcing_kv is not None: + if not isinstance(forcing_kv, dict): + raise ConfigError(f"Expected mapping for profiles.{profile}.shud.forcing.kv, got: {forcing_kv!r}") + for k, v in forcing_kv.items(): + if not isinstance(k, str) or not k.strip(): + raise ConfigError(f"Invalid forcing.kv key for profile {profile}: {k!r}") + if v is None: + raise ConfigError(f"Invalid forcing.kv value for key {k!r} (must not be null)") + sv = str(v).strip() + if not sv: + raise ConfigError(f"Invalid forcing.kv value for key {k!r} (empty string)") + if any(ch.isspace() for ch in sv): + raise ConfigError(f"Invalid forcing.kv value for key {k!r}: contains whitespace: {sv!r}") + kv.append((k.strip().upper(), sv)) + prjname = str(_get(cfg, "project.name")) out_path = run_dir / "input" / prjname / f"{prjname}.cfg.forcing" header = [