Dataclass instead of namedtuples

mikeio/dataset/_data_plot.py

@@ -563,7 +563,7 @@ def _plot_FM_map(self, ax: Axes, **kwargs: Any) -> Axes:
             node_coordinates=geometry.node_coordinates,
             element_table=geometry.element_table,
             element_coordinates=geometry.element_coordinates,
-            boundary_polylines=geometry.boundary_polylines,
+            boundary_polylines=geometry.boundary_polylines.lines,
             projection=geometry.projection,
             z=values,
             ax=ax,

mikeio/spatial/_FM_geometry.py

@@ -1,5 +1,4 @@
 from __future__ import annotations
-from collections import namedtuple
 from functools import cached_property
 from pathlib import Path
 from typing import (
@@ -9,6 +8,7 @@
     Sized,
     TYPE_CHECKING,
 )
+import warnings
 
 
 import numpy as np
@@ -23,7 +23,8 @@
 from ._FM_utils import (
     _get_node_centered_data,
     _plot_map,
-    BoundaryPolylines,
+    BoundaryPolygons,
+    Polygon,
     _set_xy_label_by_projection,  # TODO remove
     _to_polygons,  # TODO remove
 )
@@ -113,7 +114,7 @@ def _plot_FM_map(self, ax: Axes, **kwargs: Any) -> Axes:
             node_coordinates=g.node_coordinates,
             element_table=g.element_table,
             element_coordinates=g.element_coordinates,
-            boundary_polylines=g.boundary_polylines,
+            boundary_polylines=g.boundary_polygons.lines,
             plot_type=plot_type,
             projection=g.projection,
             z=None,
@@ -158,9 +159,9 @@ def outline(
 
         linwid = 1.2
         out_col = "0.4"
-        for exterior in self.g.boundary_polylines.exteriors:
+        for exterior in self.g.boundary_polygons.exteriors:
             ax.plot(*exterior.xy.T, color=out_col, linewidth=linwid)
-        for interior in self.g.boundary_polylines.interiors:
+        for interior in self.g.boundary_polygons.interiors:
             ax.plot(*interior.xy.T, color=out_col, linewidth=linwid)
         if title is not None:
             ax.set_title(title)
@@ -843,9 +844,16 @@ def get_element_area(self) -> np.ndarray:
         return np.abs(area)
 
     @cached_property
-    def boundary_polylines(self) -> BoundaryPolylines:
-        """Lists of closed polylines defining domain outline."""
-        return self._get_boundary_polylines()
+    def boundary_polylines(self) -> BoundaryPolygons:
+        warnings.warn(
+            "boundary_polylines is renamed to boundary_polygons", FutureWarning
+        )
+        return self._get_boundary_polygons()
+
+    @cached_property
+    def boundary_polygons(self) -> BoundaryPolygons:
+        """Lists of polygons defining domain outline."""
+        return self._get_boundary_polygons()
 
     def contains(self, points: np.ndarray) -> np.ndarray:
         """Test if a list of points are contained by mesh.
@@ -861,34 +869,18 @@ def contains(self, points: np.ndarray) -> np.ndarray:
             True for points inside, False otherwise
 
         """
-        import matplotlib.path as mp  # type: ignore
-
         points = np.atleast_2d(points)
 
-        exterior = self.boundary_polylines.exteriors[0]
-        cnts = mp.Path(exterior.xy).contains_points(points)
-
-        if self.boundary_polylines.n_exteriors > 1:
-            # in case of several dis-joint outer domains
-            for exterior in self.boundary_polylines.exteriors[1:]:
-                in_domain = mp.Path(exterior.xy).contains_points(points)
-                cnts = np.logical_or(cnts, in_domain)
-
-        # subtract any holes
-        for interior in self.boundary_polylines.interiors:
-            in_hole = mp.Path(interior.xy).contains_points(points)
-            cnts = np.logical_and(cnts, ~in_hole)
-
-        return cnts
+        return self.boundary_polylines.contains(points)
 
     def __contains__(self, pt: np.ndarray) -> bool:
         return self.contains(pt)[0]
 
-    def _get_boundary_polylines_uncategorized(self) -> list[list[np.int64]]:
-        """Construct closed polylines for all boundary faces."""
+    def _get_boundary_polygons_uncategorized(self) -> list[list[np.int64]]:
+        """Construct closed polygons for all boundary faces."""
         boundary_faces = self._get_boundary_faces()
         face_remains = boundary_faces.copy()
-        polylines = []
+        polygons = []
         while face_remains.shape[0] > 1:
             n0 = face_remains[:, 0]
             n1 = face_remains[:, 1]
@@ -907,34 +899,26 @@ def _get_boundary_polylines_uncategorized(self) -> list[list[np.int64]]:
                     break
 
             face_remains = np.delete(face_remains, index_to_delete, axis=0)
-            polylines.append(polyline)
-        return polylines
+            polygons.append(polyline)
+        return polygons
 
-    def _get_boundary_polylines(self) -> BoundaryPolylines:
+    def _get_boundary_polygons(self) -> BoundaryPolygons:
         """Get boundary polylines and categorize as inner or outer by assessing the signed area."""
-        polylines = self._get_boundary_polylines_uncategorized()
-
-        poly_lines_int = []
-        poly_lines_ext = []
-        Polyline = namedtuple("Polyline", ["n_nodes", "nodes", "xy", "area"])
-
-        for polyline in polylines:
-            xy = self.node_coordinates[polyline, :2]
-            area = (
-                np.dot(xy[:, 1], np.roll(xy[:, 0], 1))
-                - np.dot(xy[:, 0], np.roll(xy[:, 1], 1))
-            ) * 0.5
-            poly_line = np.asarray(polyline)
-            xy = self.node_coordinates[poly_line, 0:2]
-            poly = Polyline(len(polyline), poly_line, xy, area)
-            if area > 0:
-                poly_lines_ext.append(poly)
+        polygons = self._get_boundary_polygons_uncategorized()
+
+        interiors = []
+        exteriors = []
+
+        for polygon in polygons:
+            polygon_np = np.asarray(polygon)
+            xy = self.node_coordinates[polygon_np, 0:2]
+            poly = Polygon(xy)
+            if poly.area > 0:
+                exteriors.append(poly)
             else:
-                poly_lines_int.append(poly)
+                interiors.append(poly)
 
-        n_ext = len(poly_lines_ext)
-        n_int = len(poly_lines_int)
-        return BoundaryPolylines(n_ext, poly_lines_ext, n_int, poly_lines_int)
+        return BoundaryPolygons(exteriors=exteriors, interiors=interiors)
 
     def _get_boundary_faces(self) -> np.ndarray:
         """Construct list of faces."""

mikeio/spatial/_FM_geometry_layered.py

@@ -3,6 +3,7 @@
 from pathlib import Path
 
 from typing import Any, Iterable, Literal, Sequence
+import warnings
 
 from matplotlib.axes import Axes
 import numpy as np
@@ -12,7 +13,7 @@
 from ._FM_geometry import GeometryFM2D, _GeometryFM, _GeometryFMPlotter
 from ._geometry import GeometryPoint3D
 
-from ._FM_utils import _plot_vertical_profile, BoundaryPolylines
+from ._FM_utils import _plot_vertical_profile, BoundaryPolygons
 
 from ._utils import _relative_cumulative_distance
 
@@ -653,7 +654,14 @@ def __init__(
         self.plot = _GeometryFMPlotter(self)
 
     @property
-    def boundary_polylines(self) -> BoundaryPolylines:
+    def boundary_polylines(self) -> BoundaryPolygons:
+        warnings.warn(
+            "boundary_polylines is renamed to boundary_polygons", FutureWarning
+        )
+        return self.geometry2d.boundary_polylines
+
+    @property
+    def boundary_polygons(self) -> BoundaryPolygons:
         return self.geometry2d.boundary_polylines
 
     def contains(self, points: np.ndarray) -> np.ndarray:

mikeio/spatial/_FM_utils.py

@@ -1,5 +1,7 @@
-from __future__ import annotations
+from dataclasses import dataclass
 from typing import Any, Literal, Sequence
+
+from numpy.typing import NDArray
 from matplotlib.axes import Axes
 from matplotlib.cm import ScalarMappable
 from matplotlib.collections import PatchCollection
@@ -8,25 +10,73 @@
 from matplotlib.tri import Triangulation
 import numpy as np
 from scipy.sparse import csr_matrix
-from collections import namedtuple
 
 from ._utils import _relative_cumulative_distance
 
 
 MESH_COL = "0.95"
 MESH_COL_DARK = "0.6"
 
-BoundaryPolylines = namedtuple(
-    "BoundaryPolylines",
-    ["n_exteriors", "exteriors", "n_interiors", "interiors"],
-)
+
+@dataclass
+class Polygon:
+    xy: NDArray
+
+    @property
+    def area(self) -> float:
+        return (
+            np.dot(self.xy[:, 1], np.roll(self.xy[:, 0], 1))
+            - np.dot(self.xy[:, 0], np.roll(self.xy[:, 1], 1))
+        ) * 0.5
+
+
+@dataclass
+class BoundaryPolygons:
+    exteriors: list[Polygon]
+    interiors: list[Polygon]
+
+    @property
+    def lines(self) -> list[Polygon]:
+        return self.exteriors + self.interiors
+
+    def contains(self, points: np.ndarray) -> np.ndarray:
+        """Test if a list of points are contained by mesh.
+
+        Parameters
+        ----------
+        points : array-like n-by-2
+            x,y-coordinates of n points to be tested
+
+        Returns
+        -------
+        bool array
+            True for points inside, False otherwise
+
+        """
+        import matplotlib.path as mp  # type: ignore
+
+        exterior = self.exteriors[0]
+        cnts = mp.Path(exterior.xy).contains_points(points)
+
+        if len(self.exteriors) > 1:
+            # in case of several dis-joint outer domains
+            for exterior in self.exteriors[1:]:
+                in_domain = mp.Path(exterior.xy).contains_points(points)
+                cnts = np.logical_or(cnts, in_domain)
+
+        # subtract any holes
+        for interior in self.interiors:
+            in_hole = mp.Path(interior.xy).contains_points(points)
+            cnts = np.logical_and(cnts, ~in_hole)
+
+        return cnts
 
 
 def _plot_map(
     node_coordinates: np.ndarray,
     element_table: np.ndarray,
     element_coordinates: np.ndarray,
-    boundary_polylines: BoundaryPolylines,
+    boundary_polylines: list[Polygon],
     projection: str = "",
     z: np.ndarray | None = None,
     plot_type: Literal[
@@ -347,14 +397,14 @@ def __plot_mesh_only(ax: Axes, nc: np.ndarray, element_table: np.ndarray) -> Non
     ax.add_collection(fig_obj)
 
 
-def __plot_outline_only(ax: Axes, boundary_polylines: BoundaryPolylines) -> Axes:
+def __plot_outline_only(ax: Axes, boundary_polylines: list[Polygon]) -> Axes:
     """plot outline only (no data).
 
     Parameters
     ----------
     ax : matplotlib.axes.Axes
         axes object
-    boundary_polylines : BoundaryPolylines
+    boundary_polylines : list[PolyLine]
         boundary polylines
 
     Returns
@@ -583,23 +633,22 @@ def __add_non_tri_mesh(
     ax.add_collection(p)
 
 
-def __add_outline(ax: Axes, boundary_polylines: BoundaryPolylines) -> None:
+def __add_outline(ax: Axes, boundary_polylines: list[Polygon]) -> None:
     """add outline to axes.
 
     Parameters
     ----------
     ax : matplotlib.axes.Axes
         axes object
-    boundary_polylines: BoundaryPolylines
+    boundary_polylines: list[PolyLine]
         boundary polylines
 
     Returns
     -------
     None
 
     """
-    lines = boundary_polylines.exteriors + boundary_polylines.interiors
-    for line in lines:
+    for line in boundary_polylines:
         ax.plot(*line.xy.T, color="0.4", linewidth=1.2)
 
 

mikeio/spatial/_geometry.py

@@ -1,11 +1,54 @@
 from abc import ABC, abstractmethod
 
-from collections import namedtuple
-from typing import Any
+from dataclasses import dataclass
+from typing import Any, Sequence
 
 from mikecore.Projections import MapProjection
 
-BoundingBox = namedtuple("BoundingBox", ["left", "bottom", "right", "top"])
+
+@dataclass
+class BoundingBox:
+    """Bounding box for spatial data."""
+
+    left: float
+    bottom: float
+    right: float
+    top: float
+
+    def __post_init__(self) -> None:
+        if self.left > self.right:
+            raise ValueError(
+                f"Invalid x axis, left: {self.left} must be smaller than right: {self.right}"
+            )
+
+        if self.bottom > self.top:
+            raise ValueError(
+                f"Invalid y axis, bottom: {self.bottom} must be smaller than top: {self.top}"
+            )
+
+    def overlaps(self, other: "BoundingBox") -> bool:
+        """Check if two bounding boxes overlap."""
+        return not (
+            self.left > other.right
+            or self.bottom > other.top
+            or self.right < other.left
+            or self.top < other.bottom
+        )
+
+    @staticmethod
+    def parse(
+        values: "BoundingBox" | Sequence[float],
+    ) -> "BoundingBox":
+        match values:
+            case BoundingBox():
+                bbox = values
+            case left, bottom, right, top:
+                bbox = BoundingBox(left, bottom, right, top)
+            case _:
+                raise ValueError(
+                    "values must be a bounding box of coordinates e.g. (-10.0, 10.0 20.0, 30.0)"
+                )
+        return bbox
 
 
 class _Geometry(ABC):