From 4e89058f8e9cdf0e7543d3d63ad4104b32cbc037 Mon Sep 17 00:00:00 2001 From: akrherz Date: Fri, 20 Dec 2024 12:03:30 -0600 Subject: [PATCH 1/5] =?UTF-8?q?=F0=9F=9A=80=20Introduce=20pyiem.grid.nav?= =?UTF-8?q?=20with=20a=20fancy=20interface?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- CHANGELOG.md | 1 + MANIFEST.in | 1 + src/pyiem/grid/nav.py | 102 ++++++++++++++++++ src/pyiem/grid/nav.pyi | 14 +++ src/pyiem/models/gridnav.py | 162 ++++++++++++++++++++++++++++ tests/grid/test_grid_nav.py | 13 +++ tests/models/test_models_gridnav.py | 38 +++++++ 7 files changed, 331 insertions(+) create mode 100644 src/pyiem/grid/nav.py create mode 100644 src/pyiem/grid/nav.pyi create mode 100644 src/pyiem/models/gridnav.py create mode 100644 tests/grid/test_grid_nav.py create mode 100644 tests/models/test_models_gridnav.py diff --git a/CHANGELOG.md b/CHANGELOG.md index 7a548c00f..e025d97f4 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -17,6 +17,7 @@ All notable changes to this library are documented in this file. ### New Features - Add `allowed_as_list` option to `iemapp()` helper to stop lists. +- Add `pyiem.grid.nav` with IEM grid information in a fancy form. - Add `MapPlot.imshow` with some optimized panel plotting. - Add maximum risk threshold within SPC outlook message (#969). - Add `pyiem.era5land` with IEM grid reference information. diff --git a/MANIFEST.in b/MANIFEST.in index 733eb7666..3c53b4f4b 100644 --- a/MANIFEST.in +++ b/MANIFEST.in @@ -1,3 +1,4 @@ recursive-include pyiem *.py +recursive-include pyiem *.pyi recursive-include pyiem *.txt recursive-include pyiem *.npy diff --git a/src/pyiem/grid/nav.py b/src/pyiem/grid/nav.py new file mode 100644 index 000000000..fe26219f8 --- /dev/null +++ b/src/pyiem/grid/nav.py @@ -0,0 +1,102 @@ +"""Build Grid Navigation Metadata from the NetCDF file.""" + +from pyiem.models.gridnav import CartesianGridNavigation + +_GRID_CONFIGS = { + "IEMRE": { + "left_edge": -126.0625, + "bottom_edge": 22.9375, + "dx": 0.125, + "dy": 0.125, + "nx": 488, + "ny": 216, + }, + "IEMRE_CHINA": { + "left_edge": 69.9375, + "bottom_edge": 14.9375, + "dx": 0.125, + "dy": 0.125, + "nx": 560, + "ny": 320, + }, + "IEMRE_EUROPE": { + "left_edge": -10.0625, + "bottom_edge": 34.9375, + "dx": 0.125, + "dy": 0.125, + "nx": 400, + "ny": 280, + }, + # Lamely hardcoded for now + "ERA5LAND": { + "left_edge": -126.05, + "bottom_edge": 22.95, + "dx": 0.1, + "dy": 0.1, + "nx": 610, + "ny": 270, + }, + "ERA5LAND_CHINA": { + "left_edge": 69.95, + "bottom_edge": 14.95, + "dx": 0.1, + "dy": 0.1, + "nx": 700, + "ny": 400, + }, + "ERA5LAND_EUROPE": { + "left_edge": -10.05, + "bottom_edge": 34.95, + "dx": 0.1, + "dy": 0.1, + "nx": 500, + "ny": 350, + }, + "STAGE4": { + "crs": ( + "+proj=stere +a=6371200 +b=6371200 +lat_0=90 " + "+lon_0=-105, +lat_ts=60" + ), + "left_edge": -1_904_912.924, + "bottom_edge": -7_619_986.180, + "dx": 4_762.5, + "dy": 4_762.5, + "nx": 1121, + "ny": 881, + }, + "STAGE4_PRE2002": { + "crs": ( + "+proj=stere +a=6371200 +b=6371200 +lat_0=90 " + "+lon_0=-105, +lat_ts=60" + ), + "left_edge": -2_097_827.439, + "bottom_edge": -7_622_315.608, + "dx": 4_763.0, + "dy": 4_763.0, + "nx": 1160, + "ny": 880, + }, + "MRMS_IEMRE": { # Specific to the IEM and not in general + "left_edge": -126.0, + "bottom_edge": 23.0, + "dx": 0.01, + "dy": 0.01, + "nx": 6100, + "ny": 2700, + }, + "PRISM": { + "left_edge": -125.0 - (1 / 24.0) / 2.0, + "bottom_edge": 24.083333 - (1 / 24.0) / 2.0, + "dx": 1 / 24.0, + "dy": 1 / 24.0, + "nx": 1405, + "ny": 621, + }, +} + + +def __getattr__(name: str): + """Build stuff on the fly.""" + if name in _GRID_CONFIGS: + return CartesianGridNavigation(**_GRID_CONFIGS[name]) + raise AttributeError(f"module '{__name__}' has no attribute '{name}'") diff --git a/src/pyiem/grid/nav.pyi b/src/pyiem/grid/nav.pyi new file mode 100644 index 000000000..6975e0334 --- /dev/null +++ b/src/pyiem/grid/nav.pyi @@ -0,0 +1,14 @@ +"""Stub file.""" + +from pyiem.models.gridnav import CartesianGridNavigation + +class IEMRE(CartesianGridNavigation): ... +class IEMRE_CHINA(CartesianGridNavigation): ... +class IEMRE_EUROPE(CartesianGridNavigation): ... +class ERA5LAND(CartesianGridNavigation): ... +class ERA5LAND_CHINA(CartesianGridNavigation): ... +class ERA5LAND_EUROPE(CartesianGridNavigation): ... +class STAGE4(CartesianGridNavigation): ... +class STAGE4_PRE2002(CartesianGridNavigation): ... +class MRMS_IEMRE(CartesianGridNavigation): ... +class PRISM(CartesianGridNavigation): ... diff --git a/src/pyiem/models/gridnav.py b/src/pyiem/models/gridnav.py new file mode 100644 index 000000000..7561cc982 --- /dev/null +++ b/src/pyiem/models/gridnav.py @@ -0,0 +1,162 @@ +"""Grid Navigation Metadata.""" + +from typing import Optional, Union + +import numpy as np +from pydantic import BaseModel, ConfigDict, Field, model_validator +from pyproj import CRS, Proj +from rasterio.transform import from_origin + + +class CartesianGridNavigation(BaseModel): + """Navigation for cartesian grid with (0,0) in lower left. + + The `left_edge` and `bottom_edge` are the only required fields. The + rest are optional, but you need to have enough information to define + the grid, ie provide `dx` and `dy` or `nx` and `ny`. + """ + + model_config = ConfigDict(arbitrary_types_allowed=True) + + crs: Union[str, CRS] = Field( + default="EPSG:4326", + description="The coordinate reference system of the grid", + ) + left_edge: float = Field( + default=..., + description="The left edge of the grid in projection units", + ) + bottom_edge: float = Field( + default=..., + description="The bottom edge of the grid in projection units", + ) + top_edge: float = Field( + default=None, + description="The top edge of the grid in projection units", + ) + right_edge: float = Field( + default=None, + description="The right edge of the grid in projection units", + ) + dx: float = Field( + default=None, + description="The grid cell width in projection units", + gt=0, + ) + dy: float = Field( + default=None, + description="The grid cell height in projection units", + gt=0, + ) + nx: int = Field( + default=None, + description="The number of grid cells in the x direction", + gt=0, + ) + ny: int = Field( + default=None, + description="The number of grid cells in the y direction", + gt=0, + ) + + @property + def x_points(self) -> np.ndarray: + """These are the centers of the cells in the x direction.""" + return np.arange(self.nx) * self.dx + self.left + + @property + def y_points(self) -> np.ndarray: + """These are the centers of the cells in the y direction.""" + return np.arange(self.ny) * self.dy + self.bottom + + @property + def x_edges(self) -> np.ndarray: + """These are the edges of the x cells (n=NX + 1).""" + return np.arange(self.nx + 1) * self.dx + self.left_edge + + @property + def y_edges(self) -> np.ndarray: + """These are the edges of the y cells (n=NY + 1).""" + return np.arange(self.ny + 1) * self.dy + self.bottom_edge + + @property + def left(self) -> float: + """The centroid of the left most grid cell.""" + return self.left_edge + (self.dx / 2.0) + + @property + def right(self) -> float: + """The centroid of the right most grid cell.""" + return self.right_edge - (self.dx / 2.0) + + @property + def bottom(self) -> float: + """The centroid of the bottom most grid cell.""" + return self.bottom_edge + (self.dy / 2.0) + + @property + def top(self) -> float: + """The centroid of the top most grid cell.""" + return self.top_edge - (self.dy / 2.0) + + @property + def affine(self): + """Return the affine transformation.""" + return from_origin(self.left_edge, self.bottom_edge, self.dx, self.dy) + + @property + def affine_image(self): + """Return the transformation associated with upper left origin.""" + return from_origin(self.left_edge, self.top_edge, self.dx, -self.dy) + + @model_validator(mode="before") + def complete_definition(cls, values): + """Use information that was provided to compute other fields.""" + # We have required fields left_edge, bottom_edge + # Require that either dx/dy is provided or nx/ny is provided + if values.get("top_edge") is None: + values["top_edge"] = values["bottom_edge"] + ( + values["ny"] * values["dy"] + ) + if values.get("right_edge") is None: + values["right_edge"] = values["left_edge"] + ( + values["nx"] * values["dx"] + ) + if values.get("dx") is None: + values["dx"] = (values["right_edge"] - values["left_edge"]) / ( + values["nx"] + ) + if values.get("dy") is None: + values["dy"] = (values["top_edge"] - values["bottom_edge"]) / ( + values["ny"] + ) + # Be a bit more careful here that our grid generates a near integer + for key, spacing, edges in [ + ("nx", "dx", ["left_edge", "right_edge"]), + ("ny", "dy", ["bottom_edge", "top_edge"]), + ]: + if values.get(key) is not None: + continue + nn = (values[edges[1]] - values[edges[0]]) / values[spacing] + if abs(nn - int(nn)) > 0.01: + msg = f"Computed {key} is not approximately an integer" + raise ValueError(msg) + values[key] = int(nn) + + return values + + def find_ij( + self, lon: float, lat: float + ) -> tuple[Optional[int], Optional[int]]: + """Find the grid cell that contains the given lon/lat (EPSG: 4326).""" + x, y = Proj(self.crs)(lon, lat) + if ( + x < self.left_edge + or x >= self.right_edge + or y < self.bottom_edge + or y >= self.top_edge + ): + return None, None + i = int((x - self.left_edge) / self.dx) + j = int((y - self.bottom_edge) / self.dy) + return i, j diff --git a/tests/grid/test_grid_nav.py b/tests/grid/test_grid_nav.py new file mode 100644 index 000000000..8a6a2c157 --- /dev/null +++ b/tests/grid/test_grid_nav.py @@ -0,0 +1,13 @@ +"""Test pyiem.grid.nav""" + +from pyiem.grid import nav + + +def test_api(): + """Test basic things.""" + assert nav.IEMRE.bottom == 23.0 + + +def test_prism_calc(): + """Test that PRISM works out to what we expect.""" + assert (nav.PRISM.right - -66.50) < 0.01 diff --git a/tests/models/test_models_gridnav.py b/tests/models/test_models_gridnav.py new file mode 100644 index 000000000..74376f6e8 --- /dev/null +++ b/tests/models/test_models_gridnav.py @@ -0,0 +1,38 @@ +"""Test the gridnav module.""" + +import pytest + +from pyiem.models.gridnav import CartesianGridNavigation + + +@pytest.fixture +def cgn() -> CartesianGridNavigation: + """Return a basic CartesianGridNavigation.""" + return CartesianGridNavigation( + left_edge=0, + bottom_edge=0, + dx=1, + dy=1, + nx=10, + ny=10, + ) + + +def test_api(cgn): + """Test basic things.""" + assert cgn.bottom == 0.5 + assert len(cgn.x_points) == 10 + assert len(cgn.y_points) == 10 + assert len(cgn.x_edges) == 11 + assert len(cgn.y_edges) == 11 + assert cgn.right == cgn.x_points[-1] + + +def test_find_ij(cgn): + """See if we can get the right cell.""" + i, j = cgn.find_ij(0.5, 0.5) + assert i == 0 + assert j == 0 + i, j = cgn.find_ij(0.5, 10.5) + assert i is None + assert j is None From 0893b089be81c5f80029557c27167ac25bb554a6 Mon Sep 17 00:00:00 2001 From: akrherz Date: Fri, 20 Dec 2024 12:19:21 -0600 Subject: [PATCH 2/5] =?UTF-8?q?=E2=9C=8F=EF=B8=8F=20Address=20lint?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- src/pyiem/models/gridnav.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/pyiem/models/gridnav.py b/src/pyiem/models/gridnav.py index 7561cc982..05f006de4 100644 --- a/src/pyiem/models/gridnav.py +++ b/src/pyiem/models/gridnav.py @@ -107,7 +107,7 @@ def affine(self): @property def affine_image(self): """Return the transformation associated with upper left origin.""" - return from_origin(self.left_edge, self.top_edge, self.dx, -self.dy) + return from_origin(self.left_edge, self.top_edge, self.dx, 0 - self.dy) @model_validator(mode="before") def complete_definition(cls, values): @@ -149,7 +149,7 @@ def find_ij( self, lon: float, lat: float ) -> tuple[Optional[int], Optional[int]]: """Find the grid cell that contains the given lon/lat (EPSG: 4326).""" - x, y = Proj(self.crs)(lon, lat) + x, y = Proj(self.crs)(lon, lat) # skipcq if ( x < self.left_edge or x >= self.right_edge From 652d26e770c8e018c8cd94607e8bcf139382cf19 Mon Sep 17 00:00:00 2001 From: akrherz Date: Fri, 20 Dec 2024 12:58:52 -0600 Subject: [PATCH 3/5] =?UTF-8?q?=E2=99=BB=EF=B8=8F=20Address=20coverage?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- tests/grid/test_grid_nav.py | 3 +++ tests/models/test_models_gridnav.py | 27 +++++++++++++++++++++++++++ 2 files changed, 30 insertions(+) diff --git a/tests/grid/test_grid_nav.py b/tests/grid/test_grid_nav.py index 8a6a2c157..a795869b4 100644 --- a/tests/grid/test_grid_nav.py +++ b/tests/grid/test_grid_nav.py @@ -6,6 +6,9 @@ def test_api(): """Test basic things.""" assert nav.IEMRE.bottom == 23.0 + assert nav.IEMRE_CHINA.top == 54.875 + assert nav.IEMRE.affine + assert nav.IEMRE_EUROPE.affine_image def test_prism_calc(): diff --git a/tests/models/test_models_gridnav.py b/tests/models/test_models_gridnav.py index 74376f6e8..9008968e8 100644 --- a/tests/models/test_models_gridnav.py +++ b/tests/models/test_models_gridnav.py @@ -18,6 +18,33 @@ def cgn() -> CartesianGridNavigation: ) +def test_non_even_spacing(): + """Test that this errors.""" + with pytest.raises(ValueError): + CartesianGridNavigation( + left_edge=0, + bottom_edge=0, + top_edge=10, + right_edge=10, + dx=2.2, + dy=0.7, + ) + + +def test_computing_nxny(): + """Test that ny and ny can be computed.""" + cgn = CartesianGridNavigation( + left_edge=0, + bottom_edge=0, + dx=1, + dy=1, + right_edge=10, + top_edge=10, + ) + assert cgn.nx == 10 + assert cgn.ny == 10 + + def test_api(cgn): """Test basic things.""" assert cgn.bottom == 0.5 From b08a536067c9fdfe4c546e37b7de50cc1505f386 Mon Sep 17 00:00:00 2001 From: akrherz Date: Fri, 20 Dec 2024 13:13:29 -0600 Subject: [PATCH 4/5] =?UTF-8?q?=F0=9F=90=9B=20Fix=20intermediate=20affine?= =?UTF-8?q?=20bug?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- src/pyiem/models/gridnav.py | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/src/pyiem/models/gridnav.py b/src/pyiem/models/gridnav.py index 05f006de4..1e17de52f 100644 --- a/src/pyiem/models/gridnav.py +++ b/src/pyiem/models/gridnav.py @@ -3,9 +3,9 @@ from typing import Optional, Union import numpy as np +from affine import Affine from pydantic import BaseModel, ConfigDict, Field, model_validator from pyproj import CRS, Proj -from rasterio.transform import from_origin class CartesianGridNavigation(BaseModel): @@ -102,12 +102,14 @@ def top(self) -> float: @property def affine(self): """Return the affine transformation.""" - return from_origin(self.left_edge, self.bottom_edge, self.dx, self.dy) + return Affine(self.dx, 0, self.left_edge, 0, self.dy, self.bottom_edge) @property def affine_image(self): """Return the transformation associated with upper left origin.""" - return from_origin(self.left_edge, self.top_edge, self.dx, 0 - self.dy) + return Affine( + self.dx, 0, self.left_edge, 0, 0 - self.dy, self.top_edge + ) @model_validator(mode="before") def complete_definition(cls, values): From 74efd5957f060adba80c77b198e96562d2b13f7a Mon Sep 17 00:00:00 2001 From: akrherz Date: Fri, 20 Dec 2024 13:43:30 -0600 Subject: [PATCH 5/5] =?UTF-8?q?=E2=9C=A8=20Add=20pyiem.grid.nav.get=5Fnav?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- src/pyiem/grid/nav.py | 7 +++++++ src/pyiem/grid/nav.pyi | 2 ++ tests/grid/test_grid_nav.py | 6 ++++++ tests/models/test_models_gridnav.py | 20 +++++++++++++++++--- 4 files changed, 32 insertions(+), 3 deletions(-) diff --git a/src/pyiem/grid/nav.py b/src/pyiem/grid/nav.py index fe26219f8..cc729416b 100644 --- a/src/pyiem/grid/nav.py +++ b/src/pyiem/grid/nav.py @@ -95,6 +95,13 @@ } +def get_nav(name: str, dom: str) -> CartesianGridNavigation: + """Helper to remove some boilerplate for fetching gridnav.""" + extra = f"_{dom.upper()}" if dom != "" else "" + key = f"{name.upper()}{extra}" + return CartesianGridNavigation(**_GRID_CONFIGS[key]) + + def __getattr__(name: str): """Build stuff on the fly.""" if name in _GRID_CONFIGS: diff --git a/src/pyiem/grid/nav.pyi b/src/pyiem/grid/nav.pyi index 6975e0334..c34ea8a0e 100644 --- a/src/pyiem/grid/nav.pyi +++ b/src/pyiem/grid/nav.pyi @@ -12,3 +12,5 @@ class STAGE4(CartesianGridNavigation): ... class STAGE4_PRE2002(CartesianGridNavigation): ... class MRMS_IEMRE(CartesianGridNavigation): ... class PRISM(CartesianGridNavigation): ... + +def get_nav(name: str, region: str) -> CartesianGridNavigation: ... diff --git a/tests/grid/test_grid_nav.py b/tests/grid/test_grid_nav.py index a795869b4..bcdbee05c 100644 --- a/tests/grid/test_grid_nav.py +++ b/tests/grid/test_grid_nav.py @@ -14,3 +14,9 @@ def test_api(): def test_prism_calc(): """Test that PRISM works out to what we expect.""" assert (nav.PRISM.right - -66.50) < 0.01 + + +def test_get_nav(): + """Test that we can get a helper.""" + assert nav.get_nav("iemRE", "") + assert nav.get_nav("era5LAND", "china") diff --git a/tests/models/test_models_gridnav.py b/tests/models/test_models_gridnav.py index 9008968e8..709ada204 100644 --- a/tests/models/test_models_gridnav.py +++ b/tests/models/test_models_gridnav.py @@ -31,9 +31,23 @@ def test_non_even_spacing(): ) +def test_computing_dxdy(): + """Test that we can get a dx and dy.""" + _cgn = CartesianGridNavigation( + left_edge=0, + bottom_edge=0, + nx=10, + ny=10, + right_edge=10, + top_edge=10, + ) + assert _cgn.dx == 1 + assert _cgn.dy == 1 + + def test_computing_nxny(): """Test that ny and ny can be computed.""" - cgn = CartesianGridNavigation( + _cgn = CartesianGridNavigation( left_edge=0, bottom_edge=0, dx=1, @@ -41,8 +55,8 @@ def test_computing_nxny(): right_edge=10, top_edge=10, ) - assert cgn.nx == 10 - assert cgn.ny == 10 + assert _cgn.nx == 10 + assert _cgn.ny == 10 def test_api(cgn):