Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

🚀 Introduce pyiem.grid.nav with a fancy interface #989

Merged
merged 5 commits into from
Dec 20, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions CHANGELOG.md
Original file line number Diff line number Diff line change
Expand Up @@ -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.
Expand Down
1 change: 1 addition & 0 deletions MANIFEST.in
Original file line number Diff line number Diff line change
@@ -1,3 +1,4 @@
recursive-include pyiem *.py
recursive-include pyiem *.pyi
recursive-include pyiem *.txt
recursive-include pyiem *.npy
109 changes: 109 additions & 0 deletions src/pyiem/grid/nav.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,109 @@
"""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 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:
return CartesianGridNavigation(**_GRID_CONFIGS[name])
raise AttributeError(f"module '{__name__}' has no attribute '{name}'")
16 changes: 16 additions & 0 deletions src/pyiem/grid/nav.pyi
Original file line number Diff line number Diff line change
@@ -0,0 +1,16 @@
"""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): ...

def get_nav(name: str, region: str) -> CartesianGridNavigation: ...
164 changes: 164 additions & 0 deletions src/pyiem/models/gridnav.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,164 @@
"""Grid Navigation Metadata."""

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


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 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 Affine(
self.dx, 0, self.left_edge, 0, 0 - self.dy, self.top_edge
)

@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) # skipcq
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
22 changes: 22 additions & 0 deletions tests/grid/test_grid_nav.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,22 @@
"""Test pyiem.grid.nav"""

from pyiem.grid import nav


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():
"""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")
Loading
Loading