#!/usr/bin/env python3 """ Build the static OAT web interface for a NorESM aerosol PPE. Renders, for every perturbed parameter and every mapped AeroCom variable, the five-panel ADF-style difference figure low - default high - default %diff low %diff high high - low plus a per-parameter table of area-weighted global means, and emits an index.html with a parameter drop-down that swaps between them. Output is a self-contained directory of relative-path assets -- copy or rsync it anywhere static files can be served, no server process required. Input is the regridded cache written by prepare_oat_data.py. Parameter values (default / low / high) are read from the campaign's check_member_report.json so each page states what was actually perturbed. The interface is deliberately built around a `layout` of (experiment -> label) groupings so the planned Latin-hypercube interface can reuse the same page shell with an ensemble-member drop-down instead of a parameter one. Run with the ml-notebook conda environment, e.g.: /nird/scratch/ovewh/ml-notebook/bin/python oat_web/build_oat_site.py \\ --ppe-dir aerosol_ppe/oatv5_aerosolPPE_20260818 """ import argparse import json import sys from pathlib import Path import matplotlib matplotlib.use("Agg") import cartopy.crs as ccrs # noqa: E402 import matplotlib.pyplot as plt # noqa: E402 import numpy as np # noqa: E402 import xarray as xr # noqa: E402 from matplotlib.colors import TwoSlopeNorm # noqa: E402 sys.path.insert(0, str(Path(__file__).resolve().parent.parent)) from aerocom_oat_checker import VARIABLE_MAP # noqa: E402 UNITS = { "od550aer": "1", "abs550aer": "1", "angstrm550_865": "1", "ssa550": "1", "ccns.3": "cm-3", "cdnc_incl_ct": "cm-3", "clt": "fraction", "cllvi": "kg m-2", "clivi": "kg m-2", "pr": "kg m-2 s-1", "scre": "W m-2", "lcre": "W m-2", "fnet": "W m-2", } def _storage_name(var): return var.replace(".", "_") def area_weights(lat, lon): return np.cos(np.deg2rad(lat))[:, None] * np.ones((1, len(lon))) def global_mean(field, weights): ok = np.isfinite(field) if not ok.any(): return float("nan") return float(np.sum(field[ok] * weights[ok]) / np.sum(weights[ok])) def symmetric_limit(arrays, percentile=99.0): """Robust symmetric colour limit shared by a set of difference fields.""" vals = np.concatenate([np.abs(a[np.isfinite(a)]).ravel() for a in arrays if a is not None]) if vals.size == 0: return 1.0 lim = float(np.percentile(vals, percentile)) return lim if lim > 0 else 1.0 def relative_diff(pert, base, weights): """Percent difference, masked where the base field is too close to zero. Dividing by a near-zero base (cloud water in a dry column, CCN over a clean ocean) produces enormous percentages that carry no information and would otherwise set the colour scale for the whole panel. """ scale = abs(global_mean(base, weights)) floor = 0.01 * scale if scale > 0 else 0.0 with np.errstate(divide="ignore", invalid="ignore"): rel = (pert - base) / base * 100.0 return np.where(np.abs(base) > floor, rel, np.nan) def _add_map(fig, rect, field, cmap, norm, title, weights): ax = fig.add_axes(rect, projection=ccrs.PlateCarree()) ax.set_global() ax.coastlines(linewidth=0.4, color="#333333") mesh = None if field is not None and np.isfinite(field).any(): mesh = ax.pcolormesh( _LON, _LAT, field, cmap=cmap, norm=norm, transform=ccrs.PlateCarree(), shading="auto", ) ax.set_title(f"{title}\nglobal mean {global_mean(field, weights):.4g}", fontsize=9, pad=4) else: ax.set_title(title, fontsize=9, pad=4) ax.text(0.5, 0.5, "not available", ha="center", va="center", transform=ax.transAxes, fontsize=9, color="#999999") return mesh def render_figure(var, low, high, base, lat, lon, weights, out_path, param): global _LAT, _LON _LAT, _LON = lat, lon d_low = None if low is None else low - base d_high = None if high is None else high - base d_hl = None if (low is None or high is None) else high - low r_low = None if low is None else relative_diff(low, base, weights) r_high = None if high is None else relative_diff(high, base, weights) abs_lim = symmetric_limit([d_low, d_high, d_hl]) rel_lim = symmetric_limit([r_low, r_high]) abs_norm = TwoSlopeNorm(vmin=-abs_lim, vcenter=0.0, vmax=abs_lim) rel_norm = TwoSlopeNorm(vmin=-rel_lim, vcenter=0.0, vmax=rel_lim) cmap = plt.get_cmap("RdBu_r") fig = plt.figure(figsize=(11, 11)) unit = UNITS.get(var, "") fig.suptitle(f"{var} [{unit}] — {param}", fontsize=13, y=0.965) w, h = 0.40, 0.20 m_abs = _add_map(fig, [0.06, 0.70, w, h], d_low, cmap, abs_norm, "low − default", weights) m2 = _add_map(fig, [0.54, 0.70, w, h], d_high, cmap, abs_norm, "high − default", weights) m_rel = _add_map(fig, [0.06, 0.44, w, h], r_low, cmap, rel_norm, "relative difference, low [%]", weights) m4 = _add_map(fig, [0.54, 0.44, w, h], r_high, cmap, rel_norm, "relative difference, high [%]", weights) m5 = _add_map(fig, [0.30, 0.18, w, h], d_hl, cmap, abs_norm, "high − low", weights) m_abs = m_abs if m_abs is not None else (m2 if m2 is not None else m5) m_rel = m_rel if m_rel is not None else m4 if m_abs is not None: cax = fig.add_axes([0.10, 0.105, 0.33, 0.013]) cb = fig.colorbar(m_abs, cax=cax, orientation="horizontal", extend="both") cb.set_label(f"difference [{unit}]", fontsize=9) cb.ax.tick_params(labelsize=8) if m_rel is not None: cax = fig.add_axes([0.57, 0.105, 0.33, 0.013]) cb = fig.colorbar(m_rel, cax=cax, orientation="horizontal", extend="both") cb.set_label("relative difference [%]", fontsize=9) cb.ax.tick_params(labelsize=8) spec = VARIABLE_MAP[var] footer = f"[{spec['status']}] {var} ← {spec.get('noresm') or getattr(spec.get('compute'), '__name__', '?')}" if "note" in spec: footer += f" — {spec['note']}" fig.text(0.5, 0.055, footer, ha="center", fontsize=8, color="#555555", wrap=True) fig.savefig(out_path, dpi=110, bbox_inches="tight", facecolor="white") plt.close(fig) def load_param_values(report_json): """{param: {'component':.., 'default':.., 'low':.., 'high':..}} from the PPE report.""" with open(report_json) as f: report = json.load(f) info = {} for member in report["members"]: for p in member.get("perturbed_params", []): entry = info.setdefault(p["param"], {"component": p.get("component"), "default": p.get("default")}) value, default = p.get("value"), p.get("default") if default is None or value is None or value == default: continue entry["low" if value < default else "high"] = value return info def main(): parser = argparse.ArgumentParser(description=__doc__, formatter_class=argparse.RawDescriptionHelpFormatter) parser.add_argument("--ppe-dir", type=Path, required=True) parser.add_argument("--maps", type=Path, default=None, help="default: oat_web/data/_maps.nc") parser.add_argument("--report-json", type=Path, default=None) parser.add_argument("-o", "--output", type=Path, default=None, help="default: oat_web/site/") parser.add_argument("--only-param", nargs="+", default=None, help="render just these parameters") parser.add_argument("--no-plots", action="store_true", help="rebuild index.html only") args = parser.parse_args() here = Path(__file__).resolve().parent ppe_name = args.ppe_dir.name maps_path = args.maps or here / "data" / f"{ppe_name}_maps.nc" report_json = args.report_json or Path("paramfile") / ppe_name / "check_member_report.json" out_dir = args.output or here / "site" / ppe_name (out_dir / "plots").mkdir(parents=True, exist_ok=True) ds = xr.open_dataset(maps_path) lat, lon = ds["lat"].values, ds["lon"].values weights = area_weights(lat, lon) experiments = [str(e) for e in ds["experiment"].values] variables = [v for v in VARIABLE_MAP if _storage_name(v) in ds.data_vars] if "base" not in experiments: sys.exit("error: no 'base' experiment in the regridded cache") params = sorted({e[:-2] for e in experiments if e.endswith(("_L", "_H"))}) if args.only_param: params = [p for p in params if p in args.only_param] param_values = load_param_values(report_json) fields = {v: ds[_storage_name(v)].values for v in variables} idx = {e: i for i, e in enumerate(experiments)} def get(var, exp): return None if exp not in idx else fields[var][idx[exp]] manifest = {} for param in params: low_e, high_e = f"{param}_L", f"{param}_H" pdir = out_dir / "plots" / param pdir.mkdir(parents=True, exist_ok=True) rows, identical = [], True for var in variables: base = get(var, "base") low, high = get(var, low_e), get(var, high_e) if not args.no_plots: render_figure(var, low, high, base, lat, lon, weights, pdir / f"{var}.png", param) gm_base = global_mean(base, weights) gm_low = None if low is None else global_mean(low, weights) gm_high = None if high is None else global_mean(high, weights) for gm in (gm_low, gm_high): if gm is not None and gm != gm_base: identical = False rows.append({ "variable": var, "status": VARIABLE_MAP[var]["status"], "units": UNITS.get(var, ""), "default": gm_base, "low": gm_low, "high": gm_high, "rel_low": None if gm_low is None or gm_base == 0 else (gm_low - gm_base) / gm_base * 100, "rel_high": None if gm_high is None or gm_base == 0 else (gm_high - gm_base) / gm_base * 100, }) info = param_values.get(param, {}) manifest[param] = { "component": info.get("component"), "default": info.get("default"), "low": info.get("low"), "high": info.get("high"), "has_low": low_e in idx, "has_high": high_e in idx, "identical_to_base": identical, "rows": rows, } print(f"built {param} ({len(variables)} variables)", flush=True) payload = { "ppe": ppe_name, "variables": variables, "params": params, "manifest": manifest, "notes": {v: VARIABLE_MAP[v].get("note", "") for v in variables}, "sources": { v: (VARIABLE_MAP[v].get("noresm") or getattr(VARIABLE_MAP[v].get("compute"), "__name__", "?")) for v in variables }, } (out_dir / "index.html").write_text(render_html(payload)) (out_dir / "data.json").write_text(json.dumps(payload)) print(f"\nwrote {out_dir}/index.html ({len(params)} parameters x {len(variables)} variables)") def render_html(payload): data = json.dumps(payload) return HTML_TEMPLATE.replace("__DATA__", data) HTML_TEMPLATE = r""" NorESM aerosol PPE — OAT diagnostics

NorESM aerosol PPE — OAT diagnostics

""" if __name__ == "__main__": main()