diff --git a/CHANGELOG.md b/CHANGELOG.md index 5465d7e..6b362ba 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -7,6 +7,12 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ## [Unreleased] +### Added +- `ivert validate -V` accepts `mtl` (mean tide level) and `dtl` (diurnal tide level), passed to transformez as `vdatum:mtl` and `vdatum:dtl` (#235). transformez 1.0.1 added both references. + +### Changed +- IVERT now needs `transformez>=1.0.1` (#235), for the `vdatum:mtl` and `vdatum:dtl` references. transformez 1.0.1 also extracts the NOAA VDatum model grids (such as xGEOID20B, about 200 MB) into its cache, where 1.0.0 extracted them into whatever folder `ivert` was run from. + ### Fixed - Vertical datum shifts are read from the right place in their shift grid (#231). transformez stores each grid value at its pixel centre, but IVERT read the grid as if its outer values sat on its edges, stretching it by half a pixel at each side. Heights moved by up to about 2 mm on IVERT's 3-arcsecond grids; points near the grid's edges were off the most. - File sizes shown by `ivert database size` and `ivert database delete` no longer repeat a digit when they round up (#232). A 2047-byte file was listed as "11.0 KB" instead of "2.0 KB", and one just under 1 MB as "10231.0 KB" instead of "1.0 MB". diff --git a/docs/source/user_guide/validate.md b/docs/source/user_guide/validate.md index 6d077c0..81e8bec 100644 --- a/docs/source/user_guide/validate.md +++ b/docs/source/user_guide/validate.md @@ -58,7 +58,7 @@ Which CRS wins: the horizontal CRS comes from `-p`, else the DEM file. The verti ICESat-2 photons are stored either in WGS84 'ellipsoid' (EPSG:4979) or EGM2008 'geoid' (EPSG:3855) elevations, depending on the setting used (see 'ivert options' for more information). IVERT automatically transforms them to the DEM's vertical datum before comparing. If the DEM's CRS already encodes the vertical datum, `-V` is not needed. -The vertical datum shift is done by [transformez](https://github.com/continuous-dems/transformez). Tidal datums (`mllw`, `mlw`, `msl`, `mhw`, `mhhw`) are passed to it as its own NOAA VDatum references (`vdatum:mllw`, …) rather than as EPSG codes, and transformez decides how each surface is applied. IVERT checks that transformez can transform the datum you give before it starts, so an unsupported one stops the run straight away. +The vertical datum shift is done by [transformez](https://github.com/continuous-dems/transformez). Tidal datums (`mllw`, `mlw`, `mtl`, `dtl`, `msl`, `mhw`, `mhhw`) are passed to it as its own NOAA VDatum references (`vdatum:mllw`, …) rather than as EPSG codes, and transformez decides how each surface is applied. IVERT checks that transformez can transform the datum you give before it starts, so an unsupported one stops the run straight away. ### Photon filtering diff --git a/pyproject.toml b/pyproject.toml index 7f45003..39fe3bf 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -41,7 +41,7 @@ dependencies = [ "xarray>=2023.10.0", "fetchez>=0.9.1", "globato>=0.3.9", - "transformez>=1.0.0", + "transformez>=1.0.1", ] [dependency-groups] diff --git a/src/ivert/vdatum_lookup.py b/src/ivert/vdatum_lookup.py index e5882c1..9f39247 100644 --- a/src/ivert/vdatum_lookup.py +++ b/src/ivert/vdatum_lookup.py @@ -43,13 +43,11 @@ # Tidal datums go to transformez as its own IDs rather than EPSG codes: it builds # them from NOAA VDatum grids, and it is the authority on whether each surface is a # height or a depth. -# -# transformez's VDatum support also lists the mean tide level (mtl) and diurnal -# tide level (dtl) surfaces, but it has no vdatum:mtl or vdatum:dtl reference to -# reach them through yet. Once it does, add "mtl" and "dtl" here. _TIDAL_NAMES: dict[str, str] = { "mllw": "vdatum:mllw", "mlw": "vdatum:mlw", + "mtl": "vdatum:mtl", + "dtl": "vdatum:dtl", "msl": "vdatum:msl", "mhw": "vdatum:mhw", "mhhw": "vdatum:mhhw", @@ -96,18 +94,6 @@ def resolve_vdatum(name: str | int | None) -> str | None: return None if epsg is None else f"EPSG:{epsg}" -def _reference_errors() -> tuple[type[Exception], ...]: - """Return the exceptions transformez raises for a reference it cannot use.""" - # transformez derives these from the built-in ReferenceError, not ValueError. - from transformez.reference.parser import ( # noqa: PLC0415 - slow import - InvalidReferenceError, - ReferenceInputError, - UnsupportedReferenceError, - ) - - return InvalidReferenceError, ReferenceInputError, UnsupportedReferenceError - - def check_vdatum(reference: str) -> None: """Check that transformez can transform heights to or from a vertical reference. @@ -125,14 +111,12 @@ def check_vdatum(reference: str) -> None: resolve_reference, ) - try: - parsed = parse_reference(reference) - if parsed.vertical is None: - msg = f"{reference!r} has no vertical component." - raise ValueError(msg) - resolve_reference(parsed) - except _reference_errors() as exc: - raise ValueError(str(exc)) from exc + # transformez raises a ValueError for a reference it cannot parse or transform. + parsed = parse_reference(reference) + if parsed.vertical is None: + msg = f"{reference!r} has no vertical component." + raise ValueError(msg) + resolve_reference(parsed) def describe_vdatum(reference: str) -> str: @@ -141,10 +125,9 @@ def describe_vdatum(reference: str) -> str: parse_reference, ) - errors: tuple[type[Exception], ...] = (*_reference_errors(), ValueError) try: vertical = parse_reference(reference).vertical - except errors: + except ValueError: return "" return vertical.name if vertical is not None else "" diff --git a/tests/test_vdatum_lookup.py b/tests/test_vdatum_lookup.py new file mode 100644 index 0000000..79e2b1a --- /dev/null +++ b/tests/test_vdatum_lookup.py @@ -0,0 +1,124 @@ +"""Tests for translating vertical datum names into references transformez can use.""" + +import pytest + +from ivert import vdatum_lookup + +# --- resolve_vdatum ----------------------------------------------------------------- + + +@pytest.mark.parametrize( + ("name", "reference"), + [ + ("navd88", "EPSG:5703"), + ("NAVD88", "EPSG:5703"), + (" egm2008 ", "EPSG:3855"), + ("EGM96", "EPSG:5773"), + ("prvd02", "EPSG:6641"), + ("cgvd2013", "EPSG:6647"), + ("navd88 height (ft)", "EPSG:8228"), + ], +) +def test_geoid_datum_names_resolve_to_epsg_codes(name, reference): + """Names are case-insensitive and ignore surrounding spaces.""" + assert vdatum_lookup.resolve_vdatum(name) == reference + + +@pytest.mark.parametrize("name", ["ellipsoid", "wgs84", "WGS84"]) +def test_ellipsoid_names_resolve_to_wgs84_3d(name): + """Ellipsoidal heights mean WGS84, the photon database's frame (EPSG:4979).""" + assert vdatum_lookup.resolve_vdatum(name) == "EPSG:4979" + + +@pytest.mark.parametrize( + ("name", "reference"), + [ + ("mllw", "vdatum:mllw"), + ("MLW", "vdatum:mlw"), + ("mtl", "vdatum:mtl"), + ("dtl", "vdatum:dtl"), + ("msl", "vdatum:msl"), + ("mhw", "vdatum:mhw"), + ("mhhw", "vdatum:mhhw"), + ], +) +def test_tidal_datum_names_resolve_to_vdatum_references(name, reference): + """Tidal datums go to transformez as its own VDatum references, not EPSG codes.""" + assert vdatum_lookup.resolve_vdatum(name) == reference + + +@pytest.mark.parametrize( + ("code", "reference"), + [ + (5703, "EPSG:5703"), + ("5703", "EPSG:5703"), + ("EPSG:5703", "EPSG:5703"), + ("vdatum:mllw", "vdatum:mllw"), + ("epsg:3855", "epsg:3855"), + ], +) +def test_codes_and_qualified_references_pass_through(code, reference): + """A bare code gains its 'EPSG:'; a qualified reference is left as given.""" + assert vdatum_lookup.resolve_vdatum(code) == reference + + +@pytest.mark.parametrize("name", [None, "", "not a datum", "navd89"]) +def test_unknown_names_resolve_to_none(name): + """A name the table doesn't know is None, not a guess.""" + assert vdatum_lookup.resolve_vdatum(name) is None + + +# --- check_vdatum ------------------------------------------------------------------- + + +@pytest.mark.parametrize("reference", sorted(vdatum_lookup.list_vdatums())) +def test_every_listed_datum_is_one_transformez_can_use(reference): + """Each reference a name resolves to passes the check, so -V accepts it.""" + vdatum_lookup.check_vdatum(reference) + + +@pytest.mark.parametrize( + ("reference", "match"), + [ + ("EPSG:4326", "no vertical component"), + ("EPSG:5868", "no operation binding"), + ("vdatum:xyz", "Unknown custom reference"), + ("not a datum", "Unsupported coordinate reference"), + ], +) +def test_references_transformez_cannot_use_are_a_value_error(reference, match): + """Horizontal-only, untransformable and unknown references are all ValueErrors.""" + with pytest.raises(ValueError, match=match): + vdatum_lookup.check_vdatum(reference) + + +# --- describe_vdatum and list_vdatums ----------------------------------------------- + + +@pytest.mark.parametrize( + ("reference", "description"), + [ + ("EPSG:5703", "NAVD88 height"), + ("vdatum:mtl", "NOAA VDatum Mean Tide Level"), + ], +) +def test_references_are_described_by_transformez(reference, description): + """The description is transformez's name for the vertical reference.""" + assert vdatum_lookup.describe_vdatum(reference) == description + + +@pytest.mark.parametrize("reference", ["EPSG:4326", "vdatum:xyz", "not a datum"]) +def test_a_reference_with_no_vertical_name_is_described_as_empty(reference): + """No vertical part, or one transformez doesn't know, gives an empty string.""" + assert vdatum_lookup.describe_vdatum(reference) == "" + + +def test_the_list_maps_each_reference_to_the_names_for_it(): + """Every listed name resolves to the reference it's listed under.""" + listed = vdatum_lookup.list_vdatums() + assert "navd88" in listed["EPSG:5703"] + assert listed["vdatum:mtl"] == ["mtl"] + for reference, names in listed.items(): + assert names + for name in names: + assert vdatum_lookup.resolve_vdatum(name) == reference