diff --git a/CHANGELOG.md b/CHANGELOG.md index 1630c43..c75cd0e 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -5,6 +5,54 @@ uses Semantic Versioning for public releases. ## Unreleased ### Fixed +- `inspect_variable` reported statistics over an unstated subset of a + variable. `sst` returned `mean: 80.5`; the same field with 145 of its 162 + faces masked returned `mean: 153.0` in the same shape, with nothing saying + the number describes 17 faces. A fully masked field returned + `min/max/mean: nan` — not valid JSON — after emitting three + `RuntimeWarning: All-NaN slice encountered` on stderr the caller never + sees. Statistics are now taken over a boolean index of the finite entries; + a partly masked variable also reports `n_finite` and `n_total`, and a field + with nothing finite reports `null` rather than `nan`. Integers and booleans + carry no missing value to skip and keep the plain reductions with no extra + keys, matching `summarize_array`. +- `mesh_coverage` never reached a result computed on HPC, so `MESH_NOT_GLOBAL` + could not fire on any remote reply whatever the mesh. `AllCodeStrategies` + ships one function's code and nothing else, so the worker has no + `uxarray_mcp` to import; the measurement is now nested inside + `remote_inspect_mesh` and `remote_calculate_area`, following the same + worker-side inlining already used for `profile_coverage` and + `source_coverage`. `tests/test_remote_mesh_coverage.py` compares the three + copies key for key on a global mesh, a regional patch and a polar-hole grid, + which is the only thing standing between them and drift. +- Every result computed on HPC failed its own response contract. + `_run_on_hpc` stamps `tool=func.__name__`, so a worker reply arrived as + `remote_calculate_area`; nothing declares a contract under that name, so + `attach_provenance` skipped the required `operation` field and + `validate_response("calculate_area", ...)` answered + `{missing_fields: ["operation"], verdict: "malformed_envelope"}` for a reply + whose science was fine — and an SDK validating `structuredContent` against + the published schema rejects such a reply outright. The four `remote_*` + functions now alias to the operations they answer, and `operation` names the + contract rather than the venue, which is already in `_provenance`. +- `export` claimed success without checking anything. The three export + functions contained no `stat`, `getsize` or `nbytes` call, so + `status: "complete"` meant `to_csv` returned without raising, and + `rows_written` was `len(frame)` read off the in-memory DataFrame. Measured + on a 9-face dataset: CSV dropped all 8 attributes, including every unit, + and wrote its one NaN as an empty field with no sentinel and no mention; + single-variable NetCDF export dropped `salinity` and the CF `crs` container + while leaving `sst` pointing at it, so the file references a coordinate + reference system it does not contain. Exports now return an + `export_fidelity` block measured on the written file — `bytes_written`, + `rows_written` against `rows_expected`, `variables_dropped`, + `attributes_dropped`, `missing_values` with `missing_written_as`, + `dangling_grid_mapping` — and raise `EXPORT_ATTRIBUTES_DROPPED`, + `EXPORT_MISSING_VALUES_UNMARKED`, `EXPORT_VARIABLES_DROPPED`, + `EXPORT_DANGLING_GRID_MAPPING`, `EXPORT_ROW_COUNT_MISMATCH` or + `EXPORT_EMPTY_FILE`. The file is still written and its path still returned; + a lossy export is the export that was asked for, and only the silence was + wrong. The datasets the exporters opened are now closed. - Nothing in an area result said how much of the sphere it was summed over. A 5-degree mesh spanning 0-40E/0-40N with `sphere_radius: 6371000.0` returned `total_area: 22936016715559.137 m^2` — 4.4967% of `4*pi*R^2` — with @@ -215,6 +263,18 @@ uses Semantic Versioning for public releases. toolregistry-server 0.5.0, uxarray 2026.8.1, holoviews 1.19.0 and matplotlib 3.9.0. +### Added +- `scripts/measure_payload.py` says where a reply's bytes go, which + `tests/test_payload_budget.py` can only pass or fail on. It shares the + budget test's fixtures so a figure printed here and a budget asserted there + describe the same payload, and reports per-key and per-category + breakdowns plus the tool catalog. Pooled over five replies: 8199 bytes, + 28.9% answer, 26.2% provenance, 24.8% checks, 12.6% advice, 6.1% status. + Token counts are printed only when `tiktoken` is installed, and the + encoding is named; they are never estimated from a bytes-per-token + constant, which moves with how much of a payload is JSON punctuation and + would be wrong for every caller who does not share the guess. + ## 0.3.1 — 2026-09-05 ### Fixed - `calculate_zonal_mean` and `azimuthal_mean` now count how many of their bins diff --git a/docs/api.rst b/docs/api.rst index 08b2ec8..956e558 100644 --- a/docs/api.rst +++ b/docs/api.rst @@ -20,6 +20,10 @@ These modules contain the pure computation logic, separate from MCP and I/O. :members: :undoc-members: +.. automodule:: uxarray_mcp.domain.export_fidelity + :members: + :undoc-members: + .. automodule:: uxarray_mcp.domain.variable :members: :undoc-members: diff --git a/docs/tools.md b/docs/tools.md index 1db1f0c..7f4cc4d 100644 --- a/docs/tools.md +++ b/docs/tools.md @@ -372,6 +372,41 @@ The default also reads from the `UXARRAY_MCP_VERDICT_POLICY` environment variable, and an unrecognized value is rejected before the computation runs rather than after it has already cost something. +## Export fidelity + +`export` used to be the one operation that reported `status: "complete"` +with nothing behind it — no `stat`, no read-back, no comparison against the +source. It now returns an **`export_fidelity`** block measured on the file +that was written: `bytes_written`, `rows_written` and `rows_expected`, +`n_variables_written`, `variables_dropped`, `attributes_dropped`, +`missing_values` with `missing_written_as`, and `dangling_grid_mapping`. +`summary.rows_written` is counted out of the file too; it used to be +`len(frame)`, the number of rows the writer intended. + +Measured on a 9-face dataset carrying `sst` (units K, a long_name, a CF +`grid_mapping: crs`), `salinity` (units psu), a scalar `crs` container and +global `title`/`Conventions`: + +| Export | Before | Now | +|---|---|---| +| CSV, whole dataset | 129 bytes, `complete`, no codes | `attributes_dropped: 8`, `missing_values: 1` written as an empty field, `warning` | +| NetCDF, `sst` only | 8264 bytes, `complete`, no codes | `variables_dropped: ["crs", "salinity"]`, `dangling_grid_mapping: ["crs"]`, `warning` | + +The NetCDF case is the one worth reading twice: `sst` keeps its +`grid_mapping: "crs"` attribute into a file that no longer contains `crs`, +so a CF reader is told where to find the coordinate reference system and +finds nothing. + +A lossy export is still written and still returns its path — the caller +asked for CSV and CSV is what a CSV can hold. What changed is that +`physically_interpretable` goes `false` when a unit or a CRS did not +survive: numbers in a file that no longer says what they measure are not +interpretable, and the export is the last point at which anyone can see it. +The codes are `EXPORT_EMPTY_FILE`, `EXPORT_ROW_COUNT_MISMATCH`, +`EXPORT_VARIABLES_DROPPED`, `EXPORT_ATTRIBUTES_DROPPED`, +`EXPORT_MISSING_VALUES_UNMARKED` and `EXPORT_DANGLING_GRID_MAPPING`. A +whole-dataset NetCDF copy raises none of them. + ## Response contracts (`contract/`) Two tools under the `contract/` namespace let a caller ask what shape a diff --git a/scripts/measure_payload.py b/scripts/measure_payload.py new file mode 100644 index 0000000..8516ff3 --- /dev/null +++ b/scripts/measure_payload.py @@ -0,0 +1,263 @@ +"""Measure what a reply costs the caller, and where the bytes go. + +Every result this server returns is carried in the conversation and re-sent +on each later turn, so its size is paid for repeatedly. ``tests/ +test_payload_budget.py`` turns that into a ratchet with a pass/fail number +per operation; this script is the other half, the one that says *why* a +number is what it is -- which key holds the bytes, and how much of the reply +is answer rather than envelope. + +It uses the same fixtures as the budget test on purpose, so a figure printed +here and a budget asserted there describe the same payload. + +Bytes are exact. Tokens are not: they depend on the tokenizer, and this +repository does not depend on one. With ``tiktoken`` installed the counts +are real and the encoding is named in the output; without it the token +columns are omitted rather than estimated from a bytes-per-token ratio that +would be wrong for every caller who does not share our guess. + +Usage:: + + uv run python scripts/measure_payload.py + uv run python scripts/measure_payload.py --json + uv run --with tiktoken python scripts/measure_payload.py +""" + +from __future__ import annotations + +import argparse +import collections +import json +import os +import tempfile +import warnings +from typing import Any, Callable + +#: Categories match ``tests/test_payload_budget.py``. ``preconditions``, +#: ``postconditions`` and ``scientific_status`` count as signal: they are the +#: checked answer, not decoration around it. +CATEGORY = { + "_provenance": "provenance", + "recommended_next_steps": "advice", + "preconditions": "checks", + "postconditions": "checks", + "scientific_status": "status", + "grid_info": "grid_info", +} + +SIGNAL_CATEGORIES = ("answer", "checks", "status") + +#: Earth's mean radius, so ``calculate_area`` measures a completed payload +#: rather than the refusal it returns without a radius. +EARTH_RADIUS_M = 6371000.0 + + +def _category(key: str) -> str: + return CATEGORY.get(key, "answer") + + +def _size(obj: Any) -> int: + return len(json.dumps(obj, default=str)) + + +def _load_tokenizer() -> tuple[Callable[[str], int] | None, str | None]: + """A real token counter, or nothing. + + Returning ``None`` is deliberate. A bytes-per-token constant looks like a + measurement and is not one; the ratio moves with how much of a payload is + JSON punctuation, and this server's replies are unusually punctuation- + heavy. + """ + try: + import tiktoken + except ImportError: + return None, None + encoding = tiktoken.get_encoding("cl100k_base") + return (lambda text: len(encoding.encode(text))), "cl100k_base" + + +def build_fixtures(tmp_dir: str) -> tuple[str, str]: + """A 162-face global mesh and one face-centred field on it.""" + import numpy as np + import uxarray as ux + import xarray as xr + + grid = ux.Grid.from_structured( + lon=np.arange(0, 360, 20.0), lat=np.arange(-80, 81, 20.0) + ) + grid_file = os.path.join(tmp_dir, "grid.nc") + data_file = os.path.join(tmp_dir, "data.nc") + grid.to_xarray().to_netcdf(grid_file) + rng = np.random.default_rng(11) + xr.Dataset( + {"temperature": (["n_face"], 250 + 30 * rng.random(int(grid.n_face)))} + ).to_netcdf(data_file) + return grid_file, data_file + + +def measure(grid_file: str, data_file: str) -> dict[str, dict[str, Any]]: + """Run one call per operation family and take each reply apart.""" + from uxarray_mcp.tools.frontdoor import run_analysis + + calls: dict[str, dict[str, Any]] = { + "inspect_mesh": {}, + "calculate_area": {"sphere_radius": EARTH_RADIUS_M}, + "inspect_variable": {"variable_name": "temperature", "data_path": data_file}, + "calculate_zonal_mean": { + "variable_name": "temperature", + "data_path": data_file, + }, + "validate_dataset": {"data_path": data_file}, + } + + measured: dict[str, dict[str, Any]] = {} + with warnings.catch_warnings(): + warnings.simplefilter("ignore") + for operation, kwargs in calls.items(): + result = run_analysis(operation=operation, grid_path=grid_file, **kwargs) + per_key = {key: _size({key: value}) for key, value in result.items()} + per_category: collections.Counter[str] = collections.Counter() + for key, size in per_key.items(): + per_category[_category(key)] += size + measured[operation] = { + "bytes": _size(result), + "per_key": per_key, + "per_category": dict(per_category), + "signal_bytes": sum(per_category[c] for c in SIGNAL_CATEGORIES), + "serialized": json.dumps(result, default=str), + } + return measured + + +def measure_catalog() -> dict[str, Any]: + """How much context the tool catalog itself occupies before any call.""" + from uxarray_mcp.app import make_registry + + schemas = { + schema.get("function", schema)["name"]: schema + for schema in make_registry().get_schemas() + } + by_name = {name: _size(schema) for name, schema in schemas.items()} + total = sum(by_name.values()) + ordered = sorted(by_name.values()) + run_analysis_schema = schemas["run_analysis"] + return { + "n_tools": len(schemas), + "bytes": total, + "mean_bytes": total // len(schemas), + "median_bytes": ordered[len(ordered) // 2], + "largest": sorted(by_name.items(), key=lambda kv: -kv[1])[:5], + "run_analysis_params": len( + run_analysis_schema.get("function", run_analysis_schema)["parameters"][ + "properties" + ] + ), + } + + +def _report( + measured: dict[str, dict[str, Any]], + catalog: dict[str, Any], + count_tokens: Callable[[str], int] | None, + encoding_name: str | None, +) -> None: + rule = "=" * 74 + print(rule) + if count_tokens is None: + print("tokens: not counted (install tiktoken to count them)") + else: + print(f"tokens: {encoding_name}") + header = f"{'operation':24s} {'bytes':>7s} {'signal%':>8s}" + if count_tokens is not None: + header += f" {'tokens':>7s} {'B/token':>8s}" + print(header) + + for operation, entry in measured.items(): + total = entry["bytes"] + line = f"{operation:24s} {total:7d} {100 * entry['signal_bytes'] / total:7.1f}%" + if count_tokens is not None: + tokens = count_tokens(entry["serialized"]) + entry["tokens"] = tokens + line += f" {tokens:7d} {total / tokens:8.2f}" + print(line) + + pooled: collections.Counter[str] = collections.Counter() + for entry in measured.values(): + pooled.update(entry["per_category"]) + grand = sum(pooled.values()) + print("-" * 74) + print(f"POOLED over {len(measured)} replies: {grand} bytes") + for category, size in pooled.most_common(): + print(f" {category:14s} {size:6d} B {100 * size / grand:5.1f}%") + + print(rule) + print("calculate_area key breakdown:") + for key, size in sorted( + measured["calculate_area"]["per_key"].items(), key=lambda kv: -kv[1] + ): + print(f" {key:26s} {size:5d} B") + + print(rule) + print( + f"TOOL CATALOG: {catalog['n_tools']} tools, {catalog['bytes']} B " + f"(mean {catalog['mean_bytes']} B, median {catalog['median_bytes']} B)" + ) + if count_tokens is not None: + print(f" run_analysis params: {catalog['run_analysis_params']}") + for name, size in catalog["largest"]: + print( + f" {name:26s} {size:6d} B " + f"{100 * size / catalog['bytes']:5.1f}% of catalog" + ) + three_tools = 3 * catalog["bytes"] // catalog["n_tools"] + print( + f"RETRIEVAL: 3 tools ~= {three_tools} B vs {catalog['bytes']} B full " + f"catalog = {100 * (1 - 3 / catalog['n_tools']):.1f}% reduction" + ) + + +def main(argv: list[str] | None = None) -> int: + parser = argparse.ArgumentParser(description=__doc__.splitlines()[0]) + parser.add_argument( + "--json", + action="store_true", + help="emit the measurements as JSON instead of a table", + ) + args = parser.parse_args(argv) + + os.environ.setdefault( + "UXARRAY_MCP_STATE_DIR", tempfile.mkdtemp(prefix="payload-state") + ) + tmp_dir = tempfile.mkdtemp(prefix="payload") + grid_file, data_file = build_fixtures(tmp_dir) + measured = measure(grid_file, data_file) + catalog = measure_catalog() + count_tokens, encoding_name = _load_tokenizer() + + if args.json: + if count_tokens is not None: + for entry in measured.values(): + entry["tokens"] = count_tokens(entry["serialized"]) + for entry in measured.values(): + # The full payload is reproducible from the fixtures and would + # dominate the output. + entry.pop("serialized", None) + print( + json.dumps( + { + "token_encoding": encoding_name, + "operations": measured, + "catalog": catalog, + }, + indent=2, + default=str, + ) + ) + return 0 + + _report(measured, catalog, count_tokens, encoding_name) + return 0 + + +if __name__ == "__main__": + raise SystemExit(main()) diff --git a/src/uxarray_mcp/domain/export_fidelity.py b/src/uxarray_mcp/domain/export_fidelity.py new file mode 100644 index 0000000..659fa46 --- /dev/null +++ b/src/uxarray_mcp/domain/export_fidelity.py @@ -0,0 +1,157 @@ +"""What reached the file, measured on the file rather than on the plan. + +``export`` was the one operation that claimed success without checking +anything. The front door stamped ``scientific_status: {"status": +"complete", "physically_interpretable": null, "warning_codes": []}`` on +every export, and ``advanced.py`` had no ``st_size``, ``getsize``, +``stat`` or ``nbytes`` call anywhere in the three export functions -- so +"complete" meant "``to_csv`` returned", not "the data is in the file". + +Measured on a 9-face grid carrying ``sst`` (units K, long_name, and a CF +``grid_mapping: crs``), ``salinity`` (units psu), a scalar ``crs`` +container, and global ``title``/``Conventions``: + +CSV export of the whole dataset produced 129 bytes and reported +``{"rows_written": 9}``, a number read off the in-memory DataFrame and +never off the file. All 8 attributes were gone -- 6 on the variables, 2 +global -- so the ``sst`` column is bare numbers with no unit anywhere in +the artifact. The one NaN was written as ``2,,2.0,0``: an empty field, +with no sentinel and nothing in the reply mentioning it. The scalar +``crs`` container became a column of nine zeros, which reads as data. + +NetCDF export of ``sst`` alone produced 8264 bytes containing exactly one +variable. Its ``grid_mapping: "crs"`` attribute survived and ``crs`` did +not, so the file references a coordinate reference system it does not +contain -- a CF reader following that attribute finds nothing. + +None of that makes an export useless. A CSV of unlabelled numbers is +still the numbers. What it cannot be is silent: the caller is the one who +decides whether the loss matters, and they can only decide it if the +reply says what was lost. +""" + +from __future__ import annotations + +import csv +from pathlib import Path +from typing import Any + + +def measure_written_csv(destination: Path) -> tuple[int, int | None]: + """Bytes on disk and data rows read back out of the file. + + The row count comes from parsing the file rather than from + ``len(frame)``, which is what the summary used to report. The two + agree right up until the write is the thing that failed, which is the + only moment the number was worth having. ``csv.reader`` rather than a + newline count, because a quoted field may legally contain one. + """ + size = destination.stat().st_size + if size == 0: + return 0, 0 + try: + with destination.open("r", newline="") as handle: + rows = sum(1 for _ in csv.reader(handle)) + except (OSError, UnicodeDecodeError): # pragma: no cover - defensive + return size, None + return size, max(rows - 1, 0) + + +def count_dropped_attributes(dataset: Any) -> int: + """Every attribute the CSV writer had no column for. + + Variable attributes and global attributes are summed rather than + reported separately: the caller's question is whether the artifact + still says what the numbers mean, and either kind going missing + answers it the same way. + """ + total = len(getattr(dataset, "attrs", {}) or {}) + for name in getattr(dataset, "variables", {}): + total += len(dataset[name].attrs or {}) + return total + + +def dangling_grid_mappings(dataset: Any, written: set[str]) -> list[str]: + """CF ``grid_mapping`` targets that a written variable names and the file lacks. + + A variable that keeps ``grid_mapping: "crs"`` into a file with no + ``crs`` is worse than one that never declared a CRS: it tells a + reader where to look and the place is empty. + """ + dangling: list[str] = [] + for name in written: + try: + target = dataset[name].attrs.get("grid_mapping") + except (KeyError, AttributeError): # pragma: no cover - defensive + continue + if isinstance(target, str): + # CF allows an extended form ("crs: lat lon"); the container + # name is the first token either way. + container = target.split(":")[0].strip() + if container and container not in written and container not in dangling: + dangling.append(container) + return dangling + + +def export_fidelity( + *, + format: str, + destination: Path, + rows_expected: int | None = None, + rows_written: int | None = None, + n_variables_written: int | None = None, + variables_dropped: list[str] | None = None, + attributes_dropped: int = 0, + missing_values: int = 0, + missing_written_as: str | None = None, + dangling_grid_mapping: list[str] | None = None, +) -> dict[str, Any]: + """Assemble the block, filling in the file measurements that are free.""" + block: dict[str, Any] = { + "format": format, + "bytes_written": destination.stat().st_size if destination.exists() else 0, + "variables_dropped": variables_dropped or [], + "attributes_dropped": int(attributes_dropped), + "missing_values": int(missing_values), + "dangling_grid_mapping": dangling_grid_mapping or [], + } + if n_variables_written is not None: + block["n_variables_written"] = int(n_variables_written) + if rows_written is not None: + block["rows_written"] = int(rows_written) + if rows_expected is not None: + block["rows_expected"] = int(rows_expected) + if missing_values and missing_written_as is not None: + block["missing_written_as"] = missing_written_as + return block + + +def export_fidelity_warning_codes(fidelity: dict[str, Any]) -> list[str]: + """Codes for the ways an export can be lossy or wrong. + + An empty block yields nothing. An export written by an older build + sends no measurement, and absent measurement is not a clean bill. + """ + if not fidelity: + return [] + + codes: list[str] = [] + if fidelity.get("bytes_written") == 0: + codes.append("EXPORT_EMPTY_FILE") + + written = fidelity.get("rows_written") + expected = fidelity.get("rows_expected") + if written is not None and expected is not None and written != expected: + codes.append("EXPORT_ROW_COUNT_MISMATCH") + + if fidelity.get("variables_dropped"): + codes.append("EXPORT_VARIABLES_DROPPED") + if fidelity.get("attributes_dropped"): + # Units are attributes. A column of numbers whose unit was dropped + # is not interpretable, however faithfully the numbers were copied. + codes.append("EXPORT_ATTRIBUTES_DROPPED") + if fidelity.get("missing_values"): + codes.append("EXPORT_MISSING_VALUES_UNMARKED") + if fidelity.get("dangling_grid_mapping"): + codes.append("EXPORT_DANGLING_GRID_MAPPING") + return codes diff --git a/src/uxarray_mcp/domain/variable.py b/src/uxarray_mcp/domain/variable.py index eb07b48..ff37e74 100644 --- a/src/uxarray_mcp/domain/variable.py +++ b/src/uxarray_mcp/domain/variable.py @@ -1,10 +1,61 @@ -"""Shared variable inspection logic.""" +"""Shared variable inspection logic. + +The statistics here are the same kind of claim ``summarize_array`` makes, and +they had the same defect. Measured on a 162-face global grid with a field +masked over 145 of its faces, this module returned +``{"min": 145.0, "max": 161.0, "mean": 153.0}``: the true statistics of the +seventeen faces that held a value, presented as the statistics of the +variable, with nothing saying nine tenths of it was absent. Land masks are +ordinary in this data, so that was most fields. An all-NaN field was worse -- +``{"min": nan, "max": nan, "mean": nan}``, three ``RuntimeWarning``s on +stderr, and ``nan`` is not a JSON number. + +``n_finite``/``n_total`` follow ``summarize_array``'s contract exactly: they +appear only when they differ, because a count that always equals the size +costs payload on every call and tells the caller nothing ``shape`` does not +already say. +""" from __future__ import annotations from typing import Any, Optional +def _numeric_statistics(values: Any) -> dict[str, Any]: + """Statistics over the entries that are there, saying when some were not. + + Masked entries are dropped with a boolean index rather than by calling + ``np.nanmin`` and friends. Those emit ``RuntimeWarning: All-NaN slice + encountered`` on a fully masked field and then return ``NaN`` anyway, so + the caller got a warning on stderr they never see and a number that is + not JSON. + + Integers and booleans carry no missing value to skip, so they take the + plain reductions and never grow the two extra keys. + """ + import numpy as np + + if not np.issubdtype(values.dtype, np.inexact): + return { + "min": float(values.min()), + "max": float(values.max()), + "mean": float(values.mean()), + } + + finite = np.isfinite(values) + n_finite = int(finite.sum()) + usable = values[finite] + stats: dict[str, Any] = { + "min": float(usable.min()) if n_finite else None, + "max": float(usable.max()) if n_finite else None, + "mean": float(usable.mean()) if n_finite else None, + } + if n_finite != values.size: + stats["n_finite"] = n_finite + stats["n_total"] = int(values.size) + return stats + + def compute_variable_info(uxds: Any, variable_name: Optional[str] = None) -> dict: """Extract variable metadata and statistics from a UXarray dataset. @@ -18,7 +69,14 @@ def compute_variable_info(uxds: Any, variable_name: Optional[str] = None) -> dic Returns ------- dict - Keys: variables (list of metadata dicts), grid_info + Keys: variables (list of metadata dicts), grid_info. + + Each variable's ``statistics`` is ``{min, max, mean}`` over the + entries that are actually present, plus ``n_finite`` and ``n_total`` + when those differ. ``min``/``max``/``mean`` are ``None`` when nothing + is finite -- an honest absence rather than a ``NaN`` that no strict + JSON decoder will accept. ``statistics`` itself is ``None`` for a + variable that has no numeric statistics at all. """ import numpy as np @@ -55,12 +113,7 @@ def compute_variable_info(uxds: Any, variable_name: Optional[str] = None) -> dic try: if np.issubdtype(var.dtype, np.number): - values = var.values - var_info["statistics"] = { - "min": float(np.nanmin(values)), - "max": float(np.nanmax(values)), - "mean": float(np.nanmean(values)), - } + var_info["statistics"] = _numeric_statistics(np.asarray(var.values)) else: var_info["statistics"] = None except Exception: diff --git a/src/uxarray_mcp/provenance.py b/src/uxarray_mcp/provenance.py index d9d037d..827c3de 100644 --- a/src/uxarray_mcp/provenance.py +++ b/src/uxarray_mcp/provenance.py @@ -87,24 +87,37 @@ def attach_provenance( # against the schema rejects the reply outright. Set it here, at the # single point every contracted result already passes through, so the # promise and the payload cannot drift apart again. - if _has_contract(tool): - result.setdefault("operation", tool) + # + # The name written is the contract's, not the function's. A worker runs + # ``remote_calculate_area`` and the caller receives ``calculate_area``; + # the operation is the same one either way, and a venue name in a + # contract field would be a second thing for a client to special-case. + contract_name = _contract_name(tool) + if contract_name is not None: + result.setdefault("operation", contract_name) return result -def _has_contract(tool: str) -> bool: - """Whether this operation declares a response contract. +def _contract_name(tool: str) -> str | None: + """The declared operation this tool answers as, or ``None``. Only contracted families gain an ``operation`` field: the contract is what makes the field a promise, and adding it to results that never - promised it would be noise. + promised it would be noise. Aliases resolve here, so a ``remote_*`` + worker function reports the operation the caller asked for. """ try: from .response_contract import _CONTRACTS, _normalize - return _normalize(tool) in _CONTRACTS + name = _normalize(tool) + return name if name in _CONTRACTS else None except Exception: # pragma: no cover - defensive - return False + return None + + +def _has_contract(tool: str) -> bool: + """Whether this operation declares a response contract.""" + return _contract_name(tool) is not None def attach_scientific_status( diff --git a/src/uxarray_mcp/remote/compute_functions.py b/src/uxarray_mcp/remote/compute_functions.py index 77f3704..7986781 100644 --- a/src/uxarray_mcp/remote/compute_functions.py +++ b/src/uxarray_mcp/remote/compute_functions.py @@ -164,7 +164,8 @@ def remote_inspect_mesh(file_path: str) -> Dict[str, Any]: Returns ------- dict - Mesh topology including n_face, n_node, n_edge, source + Mesh topology including n_face, n_node, n_edge, source, and the + same ``mesh_coverage`` block the local path attaches. Notes ----- @@ -175,6 +176,103 @@ def remote_inspect_mesh(file_path: str) -> Dict[str, Any]: import uxarray as ux + # The same measurement domain/mesh_coverage.py makes, nested here rather + # than imported: AllCodeStrategies ships this function's code and nothing + # else, so the worker has no uxarray_mcp to import from. Nested, not + # module-level, for the same reason. tests/test_remote_mesh_coverage.py + # asserts the two implementations agree on the same grid, which is the + # only thing standing between them and drift. + def _mesh_coverage(_grid, _steradians=None): + import math as _math + + import numpy as _np + + _MAX_FACES = 250_000 + _DP = 6 + _cov: dict = { + "sphere_fraction": None, + "closed": None, + "euler_characteristic": None, + "lon_extent": None, + "lat_extent": None, + } + if _steradians is None: + try: + _steradians = float(_np.asarray(_grid.face_areas).sum()) + except Exception: + _steradians = None + if _steradians is not None and _math.isfinite(_steradians): + _cov["sphere_fraction"] = round(float(_steradians) / (4.0 * _math.pi), _DP) + try: + _lon = _np.asarray(_grid.node_lon, dtype=float) + _lat = _np.asarray(_grid.node_lat, dtype=float) + if _lon.size and _lat.size: + _cov["lon_extent"] = [ + round(float(_lon.min()), _DP), + round(float(_lon.max()), _DP), + ] + _cov["lat_extent"] = [ + round(float(_lat.min()), _DP), + round(float(_lat.max()), _DP), + ] + except Exception: + pass + try: + _n_face = int(_grid.n_face) + except Exception: + return _cov + # Kept identical to the local limit, and it matters more here: the + # meshes that justify an HPC endpoint are the ones above it. + if _n_face > _MAX_FACES: + _cov["topology_skipped"] = ( + f"{_n_face} faces exceeds the {_MAX_FACES}-face limit for " + "counting edge incidences; closure was not checked." + ) + return _cov + try: + _cov["euler_characteristic"] = ( + int(_grid.n_node) - int(_grid.n_edge) + _n_face + ) + except Exception: + pass + try: + _conn = _np.asarray(_grid.face_node_connectivity) + _clon = _np.asarray(_grid.node_lon, dtype=float) % 360.0 + _clat = _np.asarray(_grid.node_lat, dtype=float) + except Exception: + return _cov + _seen: dict = {} + _ids = [] + for _x, _y in zip(_clon, _clat): + if abs(abs(_y) - 90.0) < 1e-9: + _key = f"pole{_y:+.1f}" + else: + _key = f"{round(_x, _DP) % 360:.6f}_{round(_y, _DP):.6f}" + _ids.append(_seen.setdefault(_key, len(_seen))) + _n_node = len(_ids) + _incidence: dict = {} + for _face in _conn: + _nodes: list = [] + for _raw in _face: + _index = int(_raw) + if not 0 <= _index < _n_node: + continue + _node = _ids[_index] + if not _nodes or _nodes[-1] != _node: + _nodes.append(_node) + if len(_nodes) > 1 and _nodes[0] == _nodes[-1]: + _nodes.pop() + if len(_nodes) < 3: + continue + for _index, _node in enumerate(_nodes): + _other = _nodes[(_index + 1) % len(_nodes)] + _edge = (min(_node, _other), max(_node, _other)) + _incidence[_edge] = _incidence.get(_edge, 0) + 1 + _cov["closed"] = bool(_incidence) and all( + _count == 2 for _count in _incidence.values() + ) + return _cov + if file_path.lower().startswith("healpix:"): grid = ux.Grid.from_healpix(int(file_path.split(":")[1])) elif os.path.splitext(file_path.lower())[1] in [".shp", ".geojson"]: @@ -187,6 +285,7 @@ def remote_inspect_mesh(file_path: str) -> Dict[str, Any]: "n_node": int(grid.n_node), "n_edge": int(grid.n_edge), "source": file_path, + "mesh_coverage": _mesh_coverage(grid), "_worker_runtime": { "hostname": __import__("socket").gethostname(), "python_version": __import__("platform").python_version(), @@ -287,7 +386,9 @@ def remote_calculate_area(file_path: str) -> Dict[str, Any]: Returns ------- dict - Area statistics including total_area, mean_area, min_area, max_area + Area statistics including total_area, mean_area, min_area, + max_area, and the ``mesh_coverage`` block that says what fraction + of the sphere the total was summed over. Notes ----- @@ -299,6 +400,99 @@ def remote_calculate_area(file_path: str) -> Dict[str, Any]: import numpy as np import uxarray as ux + # Nested copy of domain/mesh_coverage.py -- see remote_inspect_mesh for + # why it cannot be imported or shared. Without it a remote total is the + # bare number the local path stopped shipping in #33: nothing in the + # payload says whether it was summed over the planet or over a patch. + def _mesh_coverage(_grid, _steradians=None): + import math as _math + + import numpy as _np + + _MAX_FACES = 250_000 + _DP = 6 + _cov: dict = { + "sphere_fraction": None, + "closed": None, + "euler_characteristic": None, + "lon_extent": None, + "lat_extent": None, + } + if _steradians is None: + try: + _steradians = float(_np.asarray(_grid.face_areas).sum()) + except Exception: + _steradians = None + if _steradians is not None and _math.isfinite(_steradians): + _cov["sphere_fraction"] = round(float(_steradians) / (4.0 * _math.pi), _DP) + try: + _lon = _np.asarray(_grid.node_lon, dtype=float) + _lat = _np.asarray(_grid.node_lat, dtype=float) + if _lon.size and _lat.size: + _cov["lon_extent"] = [ + round(float(_lon.min()), _DP), + round(float(_lon.max()), _DP), + ] + _cov["lat_extent"] = [ + round(float(_lat.min()), _DP), + round(float(_lat.max()), _DP), + ] + except Exception: + pass + try: + _n_face = int(_grid.n_face) + except Exception: + return _cov + if _n_face > _MAX_FACES: + _cov["topology_skipped"] = ( + f"{_n_face} faces exceeds the {_MAX_FACES}-face limit for " + "counting edge incidences; closure was not checked." + ) + return _cov + try: + _cov["euler_characteristic"] = ( + int(_grid.n_node) - int(_grid.n_edge) + _n_face + ) + except Exception: + pass + try: + _conn = _np.asarray(_grid.face_node_connectivity) + _clon = _np.asarray(_grid.node_lon, dtype=float) % 360.0 + _clat = _np.asarray(_grid.node_lat, dtype=float) + except Exception: + return _cov + _seen: dict = {} + _ids = [] + for _x, _y in zip(_clon, _clat): + if abs(abs(_y) - 90.0) < 1e-9: + _key = f"pole{_y:+.1f}" + else: + _key = f"{round(_x, _DP) % 360:.6f}_{round(_y, _DP):.6f}" + _ids.append(_seen.setdefault(_key, len(_seen))) + _n_node = len(_ids) + _incidence: dict = {} + for _face in _conn: + _nodes: list = [] + for _raw in _face: + _index = int(_raw) + if not 0 <= _index < _n_node: + continue + _node = _ids[_index] + if not _nodes or _nodes[-1] != _node: + _nodes.append(_node) + if len(_nodes) > 1 and _nodes[0] == _nodes[-1]: + _nodes.pop() + if len(_nodes) < 3: + continue + for _index, _node in enumerate(_nodes): + _other = _nodes[(_index + 1) % len(_nodes)] + _edge = (min(_node, _other), max(_node, _other)) + _incidence[_edge] = _incidence.get(_edge, 0) + 1 + _cov["closed"] = bool(_incidence) and all( + _count == 2 for _count in _incidence.values() + ) + return _cov + if file_path.lower().startswith("healpix:"): grid = ux.Grid.from_healpix(int(file_path.split(":")[1])) elif os.path.splitext(file_path.lower())[1] in [".shp", ".geojson"]: @@ -314,14 +508,18 @@ def remote_calculate_area(file_path: str) -> Dict[str, Any]: area_attrs = getattr(areas, "attrs", {}) or {} units = area_attrs.get("units") values = areas.values if hasattr(areas, "values") else np.asarray(areas) + steradians = float(np.sum(values)) return { - "total_area": float(np.sum(values)), + "total_area": steradians, "mean_area": float(np.mean(values)), "min_area": float(np.min(values)), "max_area": float(np.max(values)), "area_units": units, "n_face": int(grid.n_face), + # Measured on the unit sphere, before any radius the caller applies + # locally -- the same place compute_area_stats attaches it. + "mesh_coverage": _mesh_coverage(grid, steradians), "_worker_runtime": { "hostname": __import__("socket").gethostname(), "python_version": __import__("platform").python_version(), diff --git a/src/uxarray_mcp/response_contract.py b/src/uxarray_mcp/response_contract.py index 1f47940..ee646f4 100644 --- a/src/uxarray_mcp/response_contract.py +++ b/src/uxarray_mcp/response_contract.py @@ -238,12 +238,24 @@ def _field( } #: Operations that share another operation's declared shape. +#: +#: The ``remote_*`` entries are the same operation run on a worker, and the +#: contract is about what the caller receives, not where it was computed. +#: Without them ``_has_contract`` was false for every HPC result, so +#: ``attach_provenance`` skipped the required ``operation`` field and +#: ``validate_response`` returned ``verdict: malformed_envelope`` for a reply +#: whose science was fine. An SDK validating ``structuredContent`` against +#: the published schema rejects such a reply outright. _ALIASES: dict[str, str] = { "area": "calculate_area", "mesh": "inspect_mesh", "zonal_mean": "calculate_zonal_mean", "verify": "verification", "check": "verification", + "remote_calculate_area": "calculate_area", + "remote_inspect_mesh": "inspect_mesh", + "remote_calculate_zonal_mean": "calculate_zonal_mean", + "remote_validate_dataset": "validate_dataset", } diff --git a/src/uxarray_mcp/tools/advanced.py b/src/uxarray_mcp/tools/advanced.py index d76cd19..ec63deb 100644 --- a/src/uxarray_mcp/tools/advanced.py +++ b/src/uxarray_mcp/tools/advanced.py @@ -14,6 +14,12 @@ from matplotlib.path import Path as MplPath from uxarray_mcp.domain.dims import FACE_DIMS +from uxarray_mcp.domain.export_fidelity import ( + count_dropped_attributes, + dangling_grid_mappings, + export_fidelity, + measure_written_csv, +) from uxarray_mcp.domain.mesh import load_dataset, load_grid from uxarray_mcp.domain.remap_coverage import ( compute_scattered_coverage, @@ -1664,6 +1670,10 @@ def export_to_netcdf( destination = Path(output_path) destination.parent.mkdir(parents=True, exist_ok=True) + dropped: list[str] = [] + dangling: list[str] = [] + n_written: int | None = None + if result_handle is not None: stored_result = get_result(result_handle) artifact_path = stored_result.get("artifact_path") @@ -1685,19 +1695,41 @@ def export_to_netcdf( written = copy_artifact(dataset["data_path"], output_path) summary = {"copied_source": dataset["data_path"]} else: - ds = xr.open_dataset(dataset["data_path"]) - if variable_name not in ds: - raise ValueError( - f"Variable '{variable_name}' not found in {dataset['data_path']}." - ) - ds[[variable_name]].to_netcdf(output_path) + with xr.open_dataset(dataset["data_path"]) as ds: + if variable_name not in ds: + raise ValueError( + f"Variable '{variable_name}' not found in " + f"{dataset['data_path']}." + ) + subset = ds[[variable_name]] + subset.to_netcdf(output_path) + # Exporting one variable is a request, not a mistake; what + # was missing is the reply saying which siblings stayed + # behind. A CF grid_mapping container is the case that + # bites: the written variable keeps pointing at it. + kept = set(subset.variables) + dropped = sorted(set(ds.variables) - kept) + dangling = dangling_grid_mappings(ds, kept) + n_written = len(kept) + summary = summarize_dataset(subset) written = str(destination) - summary = summarize_dataset(ds[[variable_name]]) else: raise ValueError("Provide either result_handle or dataset_handle.") tracker.succeed("NetCDF export complete.") - response: dict[str, Any] = {"output_path": written, "summary": summary} + response: dict[str, Any] = { + "output_path": written, + "summary": summary, + # NetCDF carries attributes and its own fill values, so the two + # losses CSV cannot avoid do not arise here. + "export_fidelity": export_fidelity( + format="netcdf", + destination=Path(written), + n_variables_written=n_written, + variables_dropped=dropped, + dangling_grid_mapping=dangling, + ), + } response = attach_provenance( response, tool="export_to_netcdf", @@ -1725,6 +1757,12 @@ def export_to_csv( destination = Path(output_path) destination.parent.mkdir(parents=True, exist_ok=True) + rows_expected: int | None = None + attributes_dropped = 0 + missing_values = 0 + dropped: list[str] = [] + n_written: int | None = None + if result_handle is not None: stored_result = get_result(result_handle) artifact_path = stored_result.get("artifact_path") @@ -1737,16 +1775,21 @@ def export_to_csv( writer = csv.DictWriter(handle, fieldnames=sorted(payload)) writer.writeheader() writer.writerow(payload) - summary = {"rows_written": 1} + rows_expected = 1 else: try: - data = xr.open_dataarray(artifact) - frame = data.to_dataframe(name=data.name or "value").reset_index() + with xr.open_dataarray(artifact) as data: + frame = data.to_dataframe(name=data.name or "value").reset_index() + attributes_dropped = len(data.attrs or {}) + n_written = 1 except ValueError: - dataset_artifact = xr.open_dataset(artifact) - frame = dataset_artifact.to_dataframe().reset_index() + with xr.open_dataset(artifact) as dataset_artifact: + frame = dataset_artifact.to_dataframe().reset_index() + attributes_dropped = count_dropped_attributes(dataset_artifact) + n_written = len(dataset_artifact.data_vars) frame.to_csv(destination, index=False) - summary = {"rows_written": int(len(frame))} + rows_expected = int(len(frame)) + missing_values = int(frame.isna().sum().sum()) elif dataset_handle is not None: if session_id is None: raise ValueError("session_id is required when exporting a dataset_handle.") @@ -1757,20 +1800,46 @@ def export_to_csv( ) if dataset.get("data_path") is None: raise ValueError("Dataset handle does not include a data file to export.") - ds = xr.open_dataset(dataset["data_path"]) - if variable_name is not None and variable_name not in ds: - raise ValueError( - f"Variable '{variable_name}' not found in {dataset['data_path']}." - ) - export_ds = ds if variable_name is None else ds[[variable_name]] - frame = export_ds.to_dataframe().reset_index() - frame.to_csv(destination, index=False) - summary = {"rows_written": int(len(frame))} + with xr.open_dataset(dataset["data_path"]) as ds: + if variable_name is not None and variable_name not in ds: + raise ValueError( + f"Variable '{variable_name}' not found in {dataset['data_path']}." + ) + export_ds = ds if variable_name is None else ds[[variable_name]] + frame = export_ds.to_dataframe().reset_index() + frame.to_csv(destination, index=False) + rows_expected = int(len(frame)) + missing_values = int(frame.isna().sum().sum()) + attributes_dropped = count_dropped_attributes(export_ds) + n_written = len(export_ds.data_vars) + dropped = sorted(set(ds.variables) - set(export_ds.variables)) else: raise ValueError("Provide either result_handle or dataset_handle.") + # Counted out of the file, not off the frame: the old summary reported + # len(frame), which is the number of rows we meant to write. + _, rows_written = measure_written_csv(destination) + summary = {"rows_written": rows_written} + tracker.succeed("CSV export complete.") - response: dict[str, Any] = {"output_path": str(destination), "summary": summary} + response: dict[str, Any] = { + "output_path": str(destination), + "summary": summary, + # CSV has nowhere to put an attribute and no missing-value + # convention, so both losses are certain and only the disclosure + # was ever in question. + "export_fidelity": export_fidelity( + format="csv", + destination=destination, + rows_expected=rows_expected, + rows_written=rows_written, + n_variables_written=n_written, + variables_dropped=dropped, + attributes_dropped=attributes_dropped, + missing_values=missing_values, + missing_written_as="empty field", + ), + } response = attach_provenance( response, tool="export_to_csv", diff --git a/src/uxarray_mcp/tools/frontdoor.py b/src/uxarray_mcp/tools/frontdoor.py index e0a4774..a42303d 100644 --- a/src/uxarray_mcp/tools/frontdoor.py +++ b/src/uxarray_mcp/tools/frontdoor.py @@ -11,6 +11,7 @@ from typing import Any from uxarray_mcp.domain.anomaly_coverage import anomaly_coverage_warning_codes +from uxarray_mcp.domain.export_fidelity import export_fidelity_warning_codes from uxarray_mcp.domain.mesh_coverage import mesh_coverage_warning_codes from uxarray_mcp.domain.profile_coverage import profile_coverage_warning_codes from uxarray_mcp.domain.subset_coverage import subset_coverage_warning_codes @@ -352,6 +353,23 @@ def _finalize_analysis_result( if operation == "anomaly" else evaluate_temporal_preconditions(operation, coverage) ) + elif operation == "export" and "export_fidelity" in result: + # Keyed off the block being present for the same reason as the + # coverage gates. Before this, export was the one operation that + # returned `status: complete` with no check of any kind behind it. + # + # A lossy export is not a failed one -- the caller asked for CSV + # and CSV is what it got -- so the file still gets written and the + # reply still carries its path. What changes is that the loss is + # named. `physically_interpretable` goes False rather than staying + # unjudged when a unit or a CRS did not make it: numbers in a file + # that no longer says what they measure are not interpretable, and + # the export is the last point at which anyone can see that. + codes = export_fidelity_warning_codes(result["export_fidelity"]) + if codes: + status = "warning" + physically_interpretable = False + warning_codes.extend(codes) elif operation == "calculate_area" and "area_basis" in result: # Keyed off the block being present for the same reason as the # coverage gates: a remote worker on an older build sends none, and diff --git a/src/uxarray_mcp/tools/inspection.py b/src/uxarray_mcp/tools/inspection.py index 44c5bd5..a9a1183 100644 --- a/src/uxarray_mcp/tools/inspection.py +++ b/src/uxarray_mcp/tools/inspection.py @@ -135,7 +135,12 @@ def _inspect_variable_local( - dtype: Data type string - location: "faces", "nodes", "edges", or "other" - attrs: Variable attributes dict (units, long_name, etc.) - - statistics: {min, max, mean} if numeric, None otherwise + - statistics: {min, max, mean} over the entries that are + present, if numeric, None otherwise. Gains n_finite and + n_total when the variable is partly masked, and reports + min/max/mean as null rather than NaN when none of it is + present. A land-masked field is the common case, and its + mean is the mean of the ocean, not of the field. - grid_info: Brief grid summary {n_face, n_node, n_edge} Example: diff --git a/tests/test_export_fidelity.py b/tests/test_export_fidelity.py new file mode 100644 index 0000000..7fc19ba --- /dev/null +++ b/tests/test_export_fidelity.py @@ -0,0 +1,228 @@ +"""An export must say what it lost on the way to the file. + +``export`` was the one operation returning ``status: "complete"`` with no +check of any kind behind it: the three functions in ``tools/advanced.py`` +contained no ``stat``, ``getsize`` or ``nbytes`` call, so "complete" meant +``to_csv`` returned without raising. + +Measured on a 9-face grid carrying ``sst`` (units K, a long_name, and a CF +``grid_mapping: crs``), ``salinity`` (units psu), a scalar ``crs`` +container and global ``title``/``Conventions``: + +- CSV of the whole dataset: 129 bytes, all 8 attributes gone, the single + NaN written as ``2,,2.0,0``, and ``rows_written: 9`` read off the + in-memory DataFrame rather than the file. +- NetCDF of ``sst`` alone: 8264 bytes, one variable, ``salinity`` and + ``crs`` dropped, and ``sst`` still carrying ``grid_mapping: "crs"`` + into a file with no ``crs`` in it. + +Both replies said ``complete`` with an empty ``warning_codes``. +""" + +from __future__ import annotations + +import numpy as np +import pytest +import xarray as xr + +from uxarray_mcp.domain.export_fidelity import export_fidelity_warning_codes +from uxarray_mcp.state import create_session, register_dataset +from uxarray_mcp.tools.frontdoor import run_analysis + +ux = pytest.importorskip("uxarray") + + +@pytest.fixture +def lossy_dataset(state_dir, tmp_path): + """A dataset carrying everything CSV cannot represent. + + Units, a long_name, global attributes, one missing value and a CF + grid-mapping container -- the five things an export can silently drop, + on one small grid so the counts in the assertions are checkable by + hand. + """ + grid_file = tmp_path / "grid.nc" + data_file = tmp_path / "data.nc" + grid = ux.Grid.from_structured( + lon=np.arange(0, 360, 120.0), lat=np.arange(-30, 31, 30.0) + ) + grid.to_xarray().to_netcdf(grid_file) + + n_face = int(grid.n_face) + sst = np.arange(n_face, dtype="float64") + sst[2] = np.nan + xr.Dataset( + { + "sst": ( + ["n_face"], + sst, + {"units": "K", "long_name": "sea surface temp", "grid_mapping": "crs"}, + ), + "salinity": ( + ["n_face"], + np.arange(n_face, dtype="float64"), + {"units": "psu"}, + ), + "crs": ( + (), + np.int32(0), + { + "grid_mapping_name": "latitude_longitude", + "earth_radius": 6371000.0, + }, + ), + }, + attrs={"title": "demo", "Conventions": "CF-1.8"}, + ).to_netcdf(data_file) + + session = create_session("export-fidelity")["session_id"] + handle = register_dataset( + session, grid_path=str(grid_file), data_path=str(data_file) + )["dataset_handle"] + return session, handle, n_face + + +def _export(session, handle, output_path, output_format, **kwargs): + return run_analysis( + operation="export", + output_path=str(output_path), + output_format=output_format, + session_id=session, + dataset_handle=handle, + **kwargs, + ) + + +class TestTheCsvSaysWhatItCouldNotCarry: + def test_the_dropped_attributes_are_counted(self, lossy_dataset, tmp_path): + session, handle, _ = lossy_dataset + result = _export(session, handle, tmp_path / "out.csv", "csv") + fidelity = result["export_fidelity"] + # 6 on the three variables, 2 global. + assert fidelity["attributes_dropped"] == 8 + assert ( + "EXPORT_ATTRIBUTES_DROPPED" in result["scientific_status"]["warning_codes"] + ) + + def test_the_missing_value_is_named_and_so_is_its_spelling( + self, lossy_dataset, tmp_path + ): + """An empty CSV field is indistinguishable from an unwritten one. + + Nothing in the file says which of the two it is, so the reply has + to. + """ + session, handle, _ = lossy_dataset + output = tmp_path / "out.csv" + result = _export(session, handle, output, "csv") + fidelity = result["export_fidelity"] + assert fidelity["missing_values"] == 1 + assert fidelity["missing_written_as"] == "empty field" + assert ",," in output.read_text() + + def test_the_row_count_comes_off_the_file(self, lossy_dataset, tmp_path): + session, handle, n_face = lossy_dataset + output = tmp_path / "out.csv" + result = _export(session, handle, output, "csv") + fidelity = result["export_fidelity"] + assert fidelity["rows_written"] == n_face + assert fidelity["rows_expected"] == n_face + assert result["summary"]["rows_written"] == n_face + # The header is not a row, and the file is the only witness. + assert len(output.read_text().strip().splitlines()) == n_face + 1 + assert fidelity["bytes_written"] == output.stat().st_size + + def test_an_export_of_one_variable_names_the_ones_left_behind( + self, lossy_dataset, tmp_path + ): + session, handle, _ = lossy_dataset + result = _export( + session, handle, tmp_path / "one.csv", "csv", variable_name="sst" + ) + fidelity = result["export_fidelity"] + assert fidelity["variables_dropped"] == ["crs", "salinity"] + assert ( + "EXPORT_VARIABLES_DROPPED" in result["scientific_status"]["warning_codes"] + ) + + +class TestTheNetcdfSaysWhatItLeftBehind: + def test_a_single_variable_export_reports_the_dangling_crs( + self, lossy_dataset, tmp_path + ): + """``sst`` keeps pointing at a container that is not in the file. + + Worse than declaring no CRS at all: the attribute tells a CF + reader where to look and the place is empty. + """ + session, handle, _ = lossy_dataset + output = tmp_path / "one.nc" + result = _export(session, handle, output, "netcdf", variable_name="sst") + fidelity = result["export_fidelity"] + assert fidelity["dangling_grid_mapping"] == ["crs"] + assert fidelity["variables_dropped"] == ["crs", "salinity"] + assert fidelity["n_variables_written"] == 1 + + with xr.open_dataset(output) as back: + assert back["sst"].attrs["grid_mapping"] == "crs" + assert "crs" not in back.variables + + def test_netcdf_keeps_the_attributes_csv_cannot(self, lossy_dataset, tmp_path): + """The two formats lose different things, and the block says which.""" + session, handle, _ = lossy_dataset + result = _export( + session, handle, tmp_path / "one.nc", "netcdf", variable_name="sst" + ) + assert result["export_fidelity"]["attributes_dropped"] == 0 + assert result["export_fidelity"]["missing_values"] == 0 + + def test_a_whole_dataset_copy_loses_nothing(self, lossy_dataset, tmp_path): + """A byte copy is the one export with nothing to disclose.""" + session, handle, _ = lossy_dataset + result = _export(session, handle, tmp_path / "all.nc", "netcdf") + fidelity = result["export_fidelity"] + assert export_fidelity_warning_codes(fidelity) == [] + assert fidelity["bytes_written"] > 0 + assert result["scientific_status"]["status"] == "complete" + + +class TestTheStatusFollowsTheMeasurement: + def test_a_lossy_export_is_no_longer_reported_as_complete( + self, lossy_dataset, tmp_path + ): + session, handle, _ = lossy_dataset + result = _export(session, handle, tmp_path / "out.csv", "csv") + status = result["scientific_status"] + assert status["status"] == "warning" + # Numbers in a file that no longer says what they measure are not + # interpretable, and the export is the last place anyone can see it. + assert status["physically_interpretable"] is False + assert sorted(status["warning_codes"]) == [ + "EXPORT_ATTRIBUTES_DROPPED", + "EXPORT_MISSING_VALUES_UNMARKED", + ] + + def test_the_file_is_still_written(self, lossy_dataset, tmp_path): + """Warned, not refused. A lossy export is still the export asked for.""" + session, handle, _ = lossy_dataset + output = tmp_path / "out.csv" + result = _export(session, handle, output, "csv") + assert output.exists() + assert result["output_path"] == str(output) + assert result["outcome"] == "complete" + + +class TestAnAbsentBlockIsNotACleanBill: + def test_no_measurement_raises_no_codes(self): + assert export_fidelity_warning_codes({}) == [] + + def test_an_empty_file_is_a_code_of_its_own(self): + assert export_fidelity_warning_codes({"bytes_written": 0}) == [ + "EXPORT_EMPTY_FILE" + ] + + def test_a_short_write_is_caught_by_comparing_the_two_counts(self): + codes = export_fidelity_warning_codes( + {"bytes_written": 42, "rows_written": 3, "rows_expected": 9} + ) + assert codes == ["EXPORT_ROW_COUNT_MISMATCH"] diff --git a/tests/test_remote_mesh_coverage.py b/tests/test_remote_mesh_coverage.py new file mode 100644 index 0000000..81a3385 --- /dev/null +++ b/tests/test_remote_mesh_coverage.py @@ -0,0 +1,127 @@ +"""The worker's copy of ``mesh_coverage`` must agree with the local one. + +``AllCodeStrategies`` ships one function's code and nothing else, so a +compute function cannot import ``uxarray_mcp``; the coverage measurement is +nested inside ``remote_inspect_mesh`` and ``remote_calculate_area`` as a +literal second and third copy of ``domain/mesh_coverage.py``. Three copies +of an algorithm drift, and the drift would be invisible: a remote reply that +disagreed with a local one about whether a mesh is global would look like a +finding about the mesh. + +These tests are the thing standing between them. They run the compute +functions in-process -- they are plain Python and need no endpoint -- and +compare their block to the local one key by key. + +The gap they close was measured on the shipped code before this change: a +remote ``calculate_area`` reply carried no ``mesh_coverage`` at all, so +``mesh_coverage_warning_codes`` saw an empty dict and ``MESH_NOT_GLOBAL`` +could not fire on any HPC result, whatever the mesh. +""" + +from __future__ import annotations + +import numpy as np +import pytest + +from uxarray_mcp.domain.mesh_coverage import ( + compute_mesh_coverage, + mesh_coverage_warning_codes, +) +from uxarray_mcp.remote.compute_functions import ( + remote_calculate_area, + remote_inspect_mesh, +) + +ux = pytest.importorskip("uxarray") + +#: Mirrors tests/test_mesh_coverage.py: a global mesh, a regional patch, and +#: a structured grid that stops half a cell short of each pole. +GRIDS = { + "global": (np.arange(0, 360, 20.0), np.arange(-80, 81, 20.0)), + "patch": (np.arange(0, 41, 5.0), np.arange(0, 41, 5.0)), + "polar_hole": (np.arange(0, 360, 20.0), np.arange(-89, 90, 1.0)), +} + + +@pytest.fixture(scope="module") +def grid_files(tmp_path_factory): + directory = tmp_path_factory.mktemp("remote-coverage") + written = {} + for name, (lon, lat) in GRIDS.items(): + path = directory / f"{name}.nc" + ux.Grid.from_structured(lon=lon, lat=lat).to_xarray().to_netcdf(path) + written[name] = str(path) + return written + + +def _local(path): + return compute_mesh_coverage(ux.open_grid(path)) + + +class TestTheThreeCopiesAgree: + @pytest.mark.parametrize("name", sorted(GRIDS)) + def test_remote_inspect_mesh_matches_the_local_block(self, grid_files, name): + remote = remote_inspect_mesh(grid_files[name])["mesh_coverage"] + assert remote == _local(grid_files[name]) + + @pytest.mark.parametrize("name", sorted(GRIDS)) + def test_remote_calculate_area_matches_the_local_block(self, grid_files, name): + remote = remote_calculate_area(grid_files[name])["mesh_coverage"] + assert remote == _local(grid_files[name]) + + def test_the_two_remote_copies_agree_with_each_other(self, grid_files): + path = grid_files["polar_hole"] + assert ( + remote_inspect_mesh(path)["mesh_coverage"] + == remote_calculate_area(path)["mesh_coverage"] + ) + + +class TestAWorkerResultCanNowBeWarnedAbout: + def test_a_remote_patch_raises_the_code(self, grid_files): + result = remote_calculate_area(grid_files["patch"]) + codes = mesh_coverage_warning_codes(result["mesh_coverage"]) + assert codes == ["MESH_NOT_GLOBAL"] + # The fraction is the one the local path measures on the same mesh. + assert result["mesh_coverage"]["sphere_fraction"] == pytest.approx( + 0.044967, abs=1e-6 + ) + + def test_a_remote_global_mesh_raises_nothing(self, grid_files): + result = remote_calculate_area(grid_files["global"]) + assert mesh_coverage_warning_codes(result["mesh_coverage"]) == [] + assert result["mesh_coverage"]["closed"] is True + + def test_the_remote_total_is_still_on_the_unit_sphere(self, grid_files): + # mesh_coverage is attached before any radius is applied locally, so + # sphere_fraction and total_area describe the same measurement. + result = remote_calculate_area(grid_files["global"]) + assert result["total_area"] == pytest.approx(4.0 * np.pi, rel=1e-5) + assert result["mesh_coverage"]["sphere_fraction"] == pytest.approx( + 1.0, abs=1e-5 + ) + + +class TestTheWorkerCopyIsSelfContained: + def test_neither_function_imports_uxarray_mcp(self): + import inspect + + # Comments name the package on purpose; an executable line that + # imports it would fail on a worker that has only uxarray. + for function in (remote_inspect_mesh, remote_calculate_area): + code = [ + line.split("#", 1)[0] + for line in inspect.getsource(function).splitlines() + ] + offenders = [line for line in code if "uxarray_mcp" in line] + assert offenders == [], (function.__name__, offenders) + + def test_the_size_guard_is_the_same_number_in_all_three_copies(self): + import inspect + + from uxarray_mcp.domain import mesh_coverage as local + + limit = str(local.TOPOLOGY_MAX_FACES) + for function in (remote_inspect_mesh, remote_calculate_area): + source = inspect.getsource(function).replace("_", "") + assert limit.replace("_", "") in source, function.__name__ diff --git a/tests/test_response_contract_envelope.py b/tests/test_response_contract_envelope.py index 93362e4..956550d 100644 --- a/tests/test_response_contract_envelope.py +++ b/tests/test_response_contract_envelope.py @@ -108,6 +108,66 @@ def test_the_envelope_is_not_reported_back_as_extra(self, completed_results): assert verdict["valid"] is True, operation +class TestAWorkerResultIsTheSameOperation: + """A reply computed on HPC answers the contract the caller asked about. + + ``_run_on_hpc`` stamps ``tool=func.__name__``, so every remote reply + carried ``remote_calculate_area``. Nothing declares a contract under + that name, so ``attach_provenance`` skipped ``operation`` and + ``validate_response("calculate_area", )`` answered + ``{'missing_fields': ['operation'], 'verdict': 'malformed_envelope'}`` + for a result whose science was fine. + """ + + def _remote_area(self, grid_file): + from uxarray_mcp.provenance import attach_provenance + from uxarray_mcp.remote.compute_functions import remote_calculate_area + + return attach_provenance( + remote_calculate_area(grid_file), + tool="remote_calculate_area", + inputs={"args": [grid_file]}, + venue="hpc:test", + ) + + def test_a_remote_area_reply_validates_as_calculate_area( + self, earth_radius_mesh_files + ): + grid_file, _ = earth_radius_mesh_files + result = self._remote_area(grid_file) + assert result["operation"] == "calculate_area" + + verdict = validate_response("calculate_area", result) + assert verdict["missing_fields"] == [] + assert verdict["valid"] is True + + def test_the_venue_stays_out_of_the_contract_field(self, earth_radius_mesh_files): + """``operation`` names what was computed, not where. + + A caller matching on ``operation`` should not have to know that the + same question answered on a worker comes back under a different + name; the venue is already in ``_provenance``. + """ + grid_file, _ = earth_radius_mesh_files + result = self._remote_area(grid_file) + assert result["_provenance"]["tool"] == "remote_calculate_area" + assert result["_provenance"]["execution_venue"] == "hpc:test" + assert "remote" not in result["operation"] + + def test_an_uncontracted_tool_gains_no_operation_field(self): + """Silence, not a guess. + + ``operation`` is a contract field. Stamping it on a family that + declares no shape would advertise a promise nothing checks. + """ + from uxarray_mcp.provenance import attach_provenance + + result = attach_provenance( + {"ok": True}, tool="not_a_contracted_tool", inputs={} + ) + assert "operation" not in result + + class TestSchemaAgreesWithTheFrontDoor: def test_envelope_blocks_reuse_the_front_door_definition(self): """One description per field, whichever schema a client reads. diff --git a/tests/test_variable_statistics.py b/tests/test_variable_statistics.py new file mode 100644 index 0000000..aabc631 --- /dev/null +++ b/tests/test_variable_statistics.py @@ -0,0 +1,139 @@ +"""A variable's statistics must say what they were taken over. + +Measured before the change, on a 162-face global grid (20-degree cells, +lat -80 to 80) carrying three float fields: + + sst {'min': 0.0, 'max': 161.0, 'mean': 80.5} + sst_masked {'min': 145.0, 'max': 161.0, 'mean': 153.0} + sst_allnan {'min': nan, 'max': nan, 'mean': nan} + +``sst_masked`` holds a value on 17 of 162 faces. 153.0 is the true mean of +those seventeen, and the payload said nothing to distinguish it from the +mean of the field -- which is 80.5, forty percent lower. Land masks are +ordinary in this data, so this was most fields. + +``sst_allnan`` was worse in two ways: three ``RuntimeWarning``s went to +stderr (``All-NaN slice encountered`` twice, ``Mean of empty slice`` once), +and ``nan`` is not a JSON number, so a strict decoder loses the response +rather than the field. + +``state.summarize_array`` already had this right; the contract here is +copied from it deliberately, so the two cannot drift. +""" + +from __future__ import annotations + +import warnings + +import numpy as np +import pytest +import xarray as xr + +from uxarray_mcp.domain.variable import compute_variable_info + +ux = pytest.importorskip("uxarray") + + +@pytest.fixture(scope="module") +def masked_dataset(): + """A global mesh with a whole field, a mostly-absent one, and an empty one.""" + grid = ux.Grid.from_structured( + lon=np.arange(0, 360, 20.0), lat=np.arange(-80, 81, 20.0) + ) + n = int(grid.n_face) + whole = np.arange(n, dtype=float) + masked = whole.copy() + masked[: n - 17] = np.nan + return ux.UxDataset( + xr.Dataset( + { + "sst": (("n_face",), whole), + "sst_masked": (("n_face",), masked), + "sst_allnan": (("n_face",), np.full(n, np.nan)), + "flags": (("n_face",), np.arange(n, dtype="int32")), + } + ), + uxgrid=grid, + ) + + +def _stats(dataset, name): + info = compute_variable_info(dataset, name) + return info["variables"][0]["statistics"] + + +class TestStatisticsSayWhatTheyCovered: + def test_a_masked_field_reports_how_much_of_it_was_there(self, masked_dataset): + stats = _stats(masked_dataset, "sst_masked") + assert stats["n_finite"] == 17 + assert stats["n_total"] == 162 + # The value itself is unchanged and still correct: 153.0 is the mean + # of the seventeen faces that hold one. Only the silence was wrong. + assert stats["mean"] == pytest.approx(153.0) + + def test_the_masked_mean_differs_from_the_whole_field_mean(self, masked_dataset): + whole = _stats(masked_dataset, "sst") + masked = _stats(masked_dataset, "sst_masked") + assert whole["mean"] == pytest.approx(80.5) + assert masked["mean"] == pytest.approx(153.0) + # Nearly double. A caller who reads one as the other is wrong about + # the field by more than any rounding this server does. + assert masked["mean"] > whole["mean"] * 1.5 + + def test_a_whole_field_is_not_charged_for_the_two_extra_keys(self, masked_dataset): + stats = _stats(masked_dataset, "sst") + assert "n_finite" not in stats + assert "n_total" not in stats + + def test_an_integer_field_is_not_charged_either(self, masked_dataset): + stats = _stats(masked_dataset, "flags") + assert "n_finite" not in stats + assert stats["min"] == 0.0 + assert stats["max"] == 161.0 + + +class TestAnAbsentFieldReportsAbsence: + def test_nothing_finite_reports_null_rather_than_nan(self, masked_dataset): + stats = _stats(masked_dataset, "sst_allnan") + assert stats["min"] is None + assert stats["max"] is None + assert stats["mean"] is None + assert stats["n_finite"] == 0 + assert stats["n_total"] == 162 + + def test_the_result_survives_a_strict_json_encoder(self, masked_dataset): + import json + + info = compute_variable_info(masked_dataset) + # allow_nan=False is what a decoder in any other language enforces. + # Before the change this raised: Out of range float values are not + # JSON compliant: nan. + json.dumps( + [v["statistics"] for v in info["variables"]], + allow_nan=False, + ) + + def test_an_empty_field_emits_no_runtime_warning(self, masked_dataset): + # Recorded rather than raised. ``simplefilter("error")`` turns the + # warning into an exception that ``compute_variable_info``'s own + # ``except Exception`` swallows into ``statistics: None``, so the + # test would pass against the old code for the wrong reason. + with warnings.catch_warnings(record=True) as caught: + warnings.simplefilter("always") + stats = _stats(masked_dataset, "sst_allnan") + runtime = [w for w in caught if issubclass(w.category, RuntimeWarning)] + assert runtime == [], [str(w.message) for w in runtime] + assert stats is not None + + +class TestTheStatisticsStayOptional: + def test_a_non_numeric_field_still_reports_none(self): + grid = ux.Grid.from_structured( + lon=np.arange(0, 360, 40.0), lat=np.arange(-60, 61, 40.0) + ) + n = int(grid.n_face) + dataset = ux.UxDataset( + xr.Dataset({"label": (("n_face",), np.array(["x"] * n))}), + uxgrid=grid, + ) + assert _stats(dataset, "label") is None