feat: IDSSE-256: GridProj

python/idsse_common/idsse/common/aws_utils.py

@@ -90,7 +90,6 @@ def aws_cp(self, path: str, dest: str) -> bool:
             finally:
                 pass
 
-
     def check_for(self, issue: datetime, valid: datetime) -> Optional[Tuple[datetime, str]]:
         """Checks if an object passed issue/valid exists
 
@@ -99,8 +98,8 @@ def check_for(self, issue: datetime, valid: datetime) -> Optional[Tuple[datetime
             valid (datetime): The valid date/time used to format the path to the object's location
 
         Returns:
-            Optional[Tuple[datetime, str]]: A tuple of the valid date/time (indicated by object's 
-                                            location) and location (path) of a object, or None 
+            Optional[Tuple[datetime, str]]: A tuple of the valid date/time (indicated by object's
+                                            location) and location (path) of a object, or None
                                             if object does not exist
         """
         lead = TimeDelta(valid-issue)
@@ -161,7 +160,7 @@ def get_valids(self,
 
         Returns:
             Sequence[Tuple[datetime, str]]: A sequence of tuples with valid date/time (indicated by
-                                            object's location) and the object's location (path). 
+                                            object's location) and the object's location (path).
                                             Empty Sequence if no valids found for given time range.
         """
         if valid_start and valid_start == valid_end:

python/idsse_common/idsse/common/grid_proj.py

@@ -0,0 +1,82 @@
+"""Module that wraps pyproj (cartographic/geodetic) library and transforms objects into other data forms"""
+# -------------------------------------------------------------------------------
+# Created on Mon Jul 31 2023
+#
+# Copyright (c) 2023 Colorado State University. All rights reserved. (1)
+#
+# Contributors:
+#     Mackenzie Grimes (1)
+#
+# -------------------------------------------------------------------------------
+# pylint: disable=invalid-name
+
+from typing import Self, Tuple, Any, Optional
+
+from pyproj import CRS, Transformer
+from pyproj.enums import TransformDirection
+
+
+class GridProj:
+    """Wrapper for pyproj instance with methods to convert to and from geographic coordinates, pixels, etc."""
+    def __init__(self,
+                 crs: CRS,
+                 latitude: float,
+                 longitude: float,
+                 width: float,
+                 height: float,
+                 dx: float,
+                 dy: Optional[float] = None):
+        # pylint: disable=too-many-arguments,unpacking-non-sequence
+        self._trans = Transformer.from_crs(crs.geodetic_crs, crs)
+        self._x_offset, self._y_offset = self._trans.transform(longitude, latitude)
+        self._w = width
+        self._h = height
+        self._dx = dx
+        self._dy = dy if dy else dx
+
+    @classmethod
+    def from_proj_grid_spec(cls, proj_string: str, grid_string: str) -> Self:
+        """ Create GridProj instance from projection grid specs
+
+        Args:
+            proj_string (str): pyproj projection string
+            grid_string (str): Grid string
+
+        Returns:
+            Self: New instance of GridProj which can then be converted to any format
+        """
+        crs = CRS.from_proj4(proj_string)
+        grid_args = {
+            key[1:]: float(value)
+            for key, value in (
+                pair.split('=') for pair in grid_string.split(' ')
+            )
+        }
+        return GridProj(crs,
+                        grid_args['lat_ll'], grid_args['lon_ll'],
+                        int(grid_args['w']), int(grid_args['h']),
+                        grid_args['dx'], grid_args['dy'])
+
+    def map_proj_to_pixel(self, x, y) -> Tuple[float, float]:
+        """Map projection to a pixel"""
+        return self.map_geo_to_pixel(*self._trans.transform(x, y))  # pylint: disable=not-an-iterable
+
+    def map_pixel_to_proj(self, x: float, y: float) -> Tuple[Any, Any]:
+        """Map pixel to a projection"""
+        return self._trans.transform(*self.map_pixel_to_geo(x, y), direction=TransformDirection.INVERSE)
+
+    def map_proj_to_geo(self, x, y) -> Tuple[float, float]:
+        """Map projection to geographic coordinates"""
+        return self._trans.transform(x, y)
+
+    def map_pixel_to_geo(self, x: float, y: float) -> Tuple[float, float]:
+        """Map pixel to geographic coordinates"""
+        return x * self._dx + self._x_offset, y * self._dy + self._y_offset
+
+    def map_geo_to_proj(self, x: float, y: float) -> Tuple[Any, Any]:
+        """Map geographic coordinates to a projection"""
+        return self._trans.transform(x, y, direction=TransformDirection.INVERSE)
+
+    def map_geo_to_pixel(self, x: float, y: float) -> Tuple[float, float]:
+        """Map geographic coordinates to a pixel"""
+        return (x - self._x_offset) / self._dx, (y - self._y_offset) / self._dy

python/idsse_common/idsse/common/path_builder.py

@@ -1,6 +1,6 @@
 """Module for building (and parsing) paths that are dependent on issue, valid, and lead.
 
-In this weather forecasting data context, 
+In this weather forecasting data context,
 - issue: the datetime when a given forecast was generated (or at least when generation began)
 - valid: the datetime when a given forecast "expires" or is superceded by newer forecasts
 - lead: the time duration (often in hours) between issue and valid, a sort of forecast "lifespan"
@@ -49,7 +49,7 @@ def from_dir_filename(cls, dir_fmt: str, filename_fmt: str) -> Self:
         """Construct a PathBuilder object from specified directory and filename formats.
         Format string args must follow Python format specifications to define expected fields:
         https://docs.python.org/3/library/string.html#format-specification-mini-language
-        
+
         For example, filepaths could be expected to contain "issue" and "lead" fields using format:
         'blend.{issue.year:04d}{issue.month:02d}{issue.day:02d}/{issue.hour:02d}/{lead.hour:03d}/'
 

python/idsse_common/idsse/common/publish_confirm.py

@@ -11,7 +11,7 @@
 import pika
 from pika.exchange_type import ExchangeType
 
-from idsse.common.log_util import get_default_log_config, set_corr_id_context_var
+from idsse.common.log_util import get_default_log_config
 
 logger = logging.getLogger(__name__)
 

python/idsse_common/test/test_config.py

@@ -116,9 +116,8 @@ def __init__(self, config: dict) -> None:
             self.b_key = None
             super().__init__(config, '', ignore_missing=True)
 
-    config = WithoutKeyConfig(
-        [{'a_key': 'value for a'}, {'b_key': 'value for b'}])
-
+    config = WithoutKeyConfig([{'a_key': 'value for a'}, {'b_key': 'value for b'}])
+
     assert config.a_key == 'value for a'
     assert config.next.b_key == 'value for b'
 

python/idsse_common/test/test_grid_proj.py

@@ -0,0 +1,128 @@
+"""Test suite for grid_proj.py"""
+# --------------------------------------------------------------------------------
+# Created on Wed Aug 2 2023
+#
+# Copyright (c) 2023 Colorado State University. All rights reserved. (1)
+#
+# Contributors:
+#     Mackenzie Grimes (1)
+#
+# --------------------------------------------------------------------------------
+# pylint: disable=missing-function-docstring,redefined-outer-name,invalid-name,protected-access
+
+from pytest import fixture, approx
+
+from idsse.common.grid_proj import GridProj
+
+# example data
+EXAMPLE_PROJ_SPEC = '+proj=lcc +lat_0=25.0 +lon_0=-95.0 +lat_1=25.0 +r=6371200'
+PROJ_SPEC_WITH_OFFSET = (
+    '+proj=lcc +lat_0=25.0 +lon_0=-95.0 +lat_1=25.0 +x_0=3271.152832031251 '
+    '+y_0=263.7934687500001 +r=6371.2, +units=km'
+)
+EXAMPLE_GRID_SPEC = "+dx=2.539703 +dy=2.539703 +w=2345 +h=1597 +lat_ll=19.229 +lon_ll=-126.2766"
+EXAMPLE_PROJECTION = (-126.2766, 19.229)
+EXAMPLE_PIXEL = (0, 1)
+EXAMPLE_GEO = (-3275807.35073, -260554.63043)
+
+
+# fixtures
+@fixture
+def grid_proj() -> GridProj:
+    return GridProj.from_proj_grid_spec(EXAMPLE_PROJ_SPEC, EXAMPLE_GRID_SPEC)
+
+
+# test class methods
+def test_from_proj_grid_spec(grid_proj: GridProj):
+    assert isinstance(grid_proj, GridProj)
+
+    assert grid_proj._dx == approx(2.539703)
+    assert grid_proj._dy == grid_proj._dx
+    assert grid_proj._x_offset == approx(-3275807.3507)
+    assert grid_proj._y_offset == approx(-260554.6304)
+    assert grid_proj._h == 1597
+    assert grid_proj._w == 2345
+
+    transformer = grid_proj._trans
+    assert transformer.source_crs is not None
+    assert transformer.source_crs.type_name == 'Geographic 2D CRS'
+    assert transformer.target_crs is not None
+    assert transformer.target_crs.type_name == 'Projected CRS'
+
+
+def test_from_proj_grid_spec_with_offset():
+    proj_with_offset = GridProj.from_proj_grid_spec(PROJ_SPEC_WITH_OFFSET, EXAMPLE_GRID_SPEC)
+    assert proj_with_offset._x_offset == approx(-3272536.1979)
+    assert proj_with_offset._y_offset == approx(-260290.8369)
+
+
+# transformation methods testing
+def test_map_proj_to_pixel(grid_proj: GridProj):
+    pixel_x, pixel_y = grid_proj.map_proj_to_pixel(*EXAMPLE_PROJECTION)
+    assert (pixel_x, pixel_y) == (0, 0)
+
+    pixel_x, pixel_y = grid_proj.map_proj_to_pixel(-126.2666, 19.2333)
+    assert (pixel_x, pixel_y) == (approx(448.1133842), approx(88.24123))
+
+
+def test_map_proj_to_geo(grid_proj: GridProj):
+    geo_x, geo_y = grid_proj.map_proj_to_geo(*EXAMPLE_PROJECTION)
+    assert (geo_x, geo_y) == (approx(EXAMPLE_GEO[0]), approx(EXAMPLE_GEO[1]))
+
+    geo_x, geo_y = grid_proj.map_proj_to_geo(-126.2666, 19.2333)
+    assert (geo_x, geo_y) == (approx(-3274669.2758), approx(-260330.52389))
+
+
+def test_map_pixel_to_geo(grid_proj: GridProj):
+    geo_x, geo_y = grid_proj.map_pixel_to_geo(*EXAMPLE_PIXEL)
+    assert (geo_x, geo_y) == (approx(-3275807.35073), approx(-260552.090732))
+
+    geo_x, geo_y = grid_proj.map_pixel_to_geo(448.1133842, 8824124)
+    assert (geo_x, geo_y) == (approx(-3274669.2758), approx(22150099.56473))
+
+
+def test_map_pixel_to_proj(grid_proj: GridProj):
+    proj_x, proj_y = grid_proj.map_pixel_to_proj(0, 0)
+    assert (proj_x, proj_y) == (approx(EXAMPLE_PROJECTION[0]), approx(EXAMPLE_PROJECTION[1]))
+
+    proj_x, proj_y = grid_proj.map_pixel_to_proj(2000, 1500)
+    assert (proj_x, proj_y) == (approx(-126.238035), approx(19.272773))
+
+
+def test_map_geo_to_proj(grid_proj: GridProj):
+    proj_x, proj_y = grid_proj.map_geo_to_proj(*EXAMPLE_GEO)
+    assert (proj_x, proj_y) == (approx(EXAMPLE_PROJECTION[0]), approx(EXAMPLE_PROJECTION[1]))
+
+    proj_x, proj_y = grid_proj.map_geo_to_proj(-3274669.2758, 22150099.56473)
+    assert (proj_x, proj_y) == (approx(-110.868173), approx(62.982259))
+
+
+def test_geo_to_pixel(grid_proj: GridProj):
+    pixel_x, pixel_y = grid_proj.map_geo_to_pixel(*EXAMPLE_GEO)
+    assert (round(pixel_x, 5), round(pixel_y, 5)) == (0, 0)
+
+    pixel_x, pixel_y = grid_proj.map_geo_to_pixel(-3274669.2758, 22150099.56473)
+    assert (pixel_x, pixel_y) == (approx(448.1133948), approx(8824124))
+
+
+def test_transformations_between_all_formats(grid_proj: GridProj):
+    # start with pixel, convert to projection
+    initial_pixel = (2000, 1500)
+    proj_x, proj_y = grid_proj.map_pixel_to_proj(*initial_pixel)
+    initial_proj = (proj_x, proj_y)
+
+    # convert projection to geographic coordinates
+    geo_x, geo_y = grid_proj.map_proj_to_geo(proj_x, proj_y)
+    initial_geo = geo_x, geo_y
+
+    # convert geographic coordinates back to pixel, full circle, and data should be unchanged
+    pixel_x, pixel_y = grid_proj.map_geo_to_pixel(geo_x, geo_y)
+    assert (pixel_x, pixel_y) == (approx(initial_pixel[0]), approx(initial_pixel[1]))
+
+    # convert pixel back to geographic coordinates
+    geo_x, geo_y = grid_proj.map_pixel_to_geo(pixel_x, pixel_y)
+    assert (geo_x, geo_y) == (approx(initial_geo[0]), approx(initial_geo[1]))
+
+    # convert geographic coordinates back to projection
+    proj_x, proj_y = grid_proj.map_geo_to_proj(geo_x, geo_y)
+    assert (proj_x, proj_y) == (approx(initial_proj[0]), approx(initial_proj[1]))