diff --git a/movement/roi/__init__.py b/movement/roi/__init__.py new file mode 100644 index 000000000..35657bc1f --- /dev/null +++ b/movement/roi/__init__.py @@ -0,0 +1,2 @@ +from movement.roi.line import LineOfInterest +from movement.roi.polygon import PolygonOfInterest diff --git a/movement/roi/base.py b/movement/roi/base.py new file mode 100644 index 000000000..d91b8aac3 --- /dev/null +++ b/movement/roi/base.py @@ -0,0 +1,160 @@ +"""Class for representing 1- or 2-dimensional regions of interest (RoIs).""" + +from __future__ import annotations + +from collections.abc import Sequence +from typing import Literal, TypeAlias + +import shapely +from shapely.coords import CoordinateSequence + +from movement.utils.logging import log_error + +LineLike: TypeAlias = shapely.LinearRing | shapely.LineString +PointLike: TypeAlias = tuple[float, float] +PointLikeList: TypeAlias = Sequence[PointLike] +RegionLike: TypeAlias = shapely.Polygon +SupportedGeometry: TypeAlias = LineLike | RegionLike + + +class BaseRegionOfInterest: + """Base class for representing regions of interest (RoIs). + + Regions of interest can be either 1 or 2 dimensional, and are represented + by appropriate ``shapely.Geometry`` objects depending on which. Note that + there are a number of discussions concerning subclassing ``shapely`` + objects; + + - https://github.com/shapely/shapely/issues/1233. + - https://stackoverflow.com/questions/10788976/how-do-i-properly-inherit-from-a-superclass-that-has-a-new-method + + To avoid the complexities of subclassing ourselves, we simply elect to wrap + the appropriate ``shapely`` object in the ``_shapely_geometry`` attribute, + accessible via the property ``region``. This also has the benefit of + allowing us to 'forbid' certain operations (that ``shapely`` would + otherwise interpret in a set-theoretic sense, giving confusing answers to + users). + + This class is not designed to be instantiated directly. It can be + instantiated, however its primary purpose is to reduce code duplication. + """ + + __default_name: str = "Un-named region" + + _name: str | None + _shapely_geometry: SupportedGeometry + + @property + def coords(self) -> CoordinateSequence: + """Coordinates of the points that define the region. + + These are the points passed to the constructor argument ``points``. + + Note that for Polygonal regions, these are the coordinates of the + exterior boundary, interior boundaries must be accessed via + ``self.region.interior.coords``. + """ + return ( + self.region.coords + if self.dimensions < 2 + else self.region.exterior.coords + ) + + @property + def dimensions(self) -> int: + """Dimensionality of the region.""" + return shapely.get_dimensions(self.region) + + @property + def is_closed(self) -> bool: + """Return True if the region is closed. + + A closed region is either: + - A polygon (2D RoI). + - A 1D LoI whose final point connects back to its first. + """ + return self.dimensions > 1 or ( + self.dimensions == 1 + and self.region.coords[0] == self.region.coords[-1] + ) + + @property + def name(self) -> str: + """Name of the instance.""" + return self._name if self._name else self.__default_name + + @property + def region(self) -> SupportedGeometry: + """``shapely.Geometry`` representation of the region.""" + return self._shapely_geometry + + def __init__( + self, + points: PointLikeList, + dimensions: Literal[1, 2] = 2, + closed: bool = False, + holes: Sequence[PointLikeList] | None = None, + name: str | None = None, + ) -> None: + """Initialise a region of interest. + + Parameters + ---------- + points : Sequence of (x, y) values + Sequence of (x, y) coordinate pairs that will form the region. + dimensions : Literal[1, 2], default 2 + The dimensionality of the region to construct. + '1' creates a sequence of joined line segments, + '2' creates a polygon whose boundary is defined by ``points``. + closed : bool, default False + Whether the line to be created should be closed. That is, whether + the final point should also link to the first point. + Ignored if ``dimensions`` is 2. + holes : sequence of sequences of (x, y) pairs, default None + A sequence of items, where each item will be interpreted like + ``points``. These items will be used to construct internal holes + within the region. See the ``holes`` argument to + ``shapely.Polygon`` for details. Ignored if ``dimensions`` is 1. + name : str, default None + Human-readable name to assign to the given region, for + user-friendliness. Default name given is 'Un-named region' if no + explicit name is provided. + + """ + self._name = name + if len(points) < dimensions + 1: + raise log_error( + ValueError, + f"Need at least {dimensions + 1} points to define a " + f"{dimensions}D region (got {len(points)}).", + ) + elif dimensions < 1 or dimensions > 2: + raise log_error( + ValueError, + "Only regions of interest of dimension 1 or 2 are supported " + f"(requested {dimensions})", + ) + elif dimensions == 1 and len(points) < 3 and closed: + raise log_error( + ValueError, + "Cannot create a loop from a single line segment.", + ) + if dimensions == 2: + self._shapely_geometry = shapely.Polygon(shell=points, holes=holes) + else: + self._shapely_geometry = ( + shapely.LinearRing(coordinates=points) + if closed + else shapely.LineString(coordinates=points) + ) + + def __repr__(self) -> str: # noqa: D105 + return str(self) + + def __str__(self) -> str: # noqa: D105 + display_type = "-gon" if self.dimensions > 1 else " line segment(s)" + n_points = len(self.coords) - 1 + return ( + f"{self.__class__.__name__} {self.name} " + f"({n_points}{display_type})\n" + ) + " -> ".join(f"({c[0]}, {c[1]})" for c in self.coords) diff --git a/movement/roi/line.py b/movement/roi/line.py new file mode 100644 index 000000000..5359c830e --- /dev/null +++ b/movement/roi/line.py @@ -0,0 +1,58 @@ +"""1-dimensional lines of interest.""" + +from movement.roi.base import BaseRegionOfInterest, PointLikeList + + +class LineOfInterest(BaseRegionOfInterest): + """Representation of boundaries or other lines of interest. + + This class can be used to represent boundaries or other internal divisions + of the area in which the experimental data was gathered. These might + include segments of a wall that are removed partway through a behavioural + study, or coloured marking on the floor of the experimental enclosure that + have some significance. Instances of this class also constitute the + boundary of two-dimensional regions (polygons) of interest. + + An instance of this class can be used to represent these "one dimensional + regions" (lines of interest, LoIs) in an analysis. The basic usage is to + construct an instance of this class by passing in a list of points, which + will then be joined (in sequence) by straight lines between consecutive + pairs of points, to form the LoI that is to be studied. + """ + + def __init__( + self, + points: PointLikeList, + loop: bool = False, + name: str | None = None, + ) -> None: + """Create a new line of interest (LoI). + + Parameters + ---------- + points : tuple of (x, y) pairs + The points (in sequence) that make up the line segment. At least + two points must be provided. + loop : bool, default False + If True, the final point in ``points`` will be connected by an + additional line segment to the first, creating a closed loop. + (See Notes). + name : str, optional + Name of the LoI that is to be created. A default name will be + inherited from the base class if not provided, and + defaults are inherited from. + + Notes + ----- + The constructor supports 'rings' or 'closed loops' via the ``loop`` + argument. However, if you want to define an enclosed region for your + analysis, we recommend you create a ``PolygonOfInterest`` and use + its ``boundary`` property instead. + + See Also + -------- + movement.roi.base.BaseRegionOfInterest + The base class that constructor arguments are passed to. + + """ + super().__init__(points, dimensions=1, closed=loop, name=name) diff --git a/movement/roi/polygon.py b/movement/roi/polygon.py new file mode 100644 index 000000000..b3c8e4c8e --- /dev/null +++ b/movement/roi/polygon.py @@ -0,0 +1,105 @@ +"""2-dimensional regions of interest.""" + +from __future__ import annotations + +from collections.abc import Sequence + +from movement.roi.base import BaseRegionOfInterest, PointLikeList +from movement.roi.line import LineOfInterest + + +class PolygonOfInterest(BaseRegionOfInterest): + """Representation of a two-dimensional region in the x-y plane. + + This class can be used to represent polygonal regions or subregions + of the area in which the experimental data was gathered. These might + include the arms of a maze, a nesting area, a food source, or other + similar areas of the experimental enclosure that have some significance. + + An instance of this class can be used to represent these regions of + interest (RoIs) in an analysis. The basic usage is to construct an + instance of this class by passing in a list of points, which will then be + joined (in sequence) by straight lines between consecutive pairs of points, + to form the exterior boundary of the RoI. Note that the exterior boundary + (accessible as via the ``.exterior`` property) is a (closed) + ``LineOfInterest``, and may be treated accordingly. + + The class also supports holes - subregions properly contained inside the + region that are not part of the region itself. These can be specified by + the ``holes`` argument, and define the interior boundaries of the region. + These interior boundaries are accessible via the ``.interior_boundaries`` + property, and the polygonal regions that make up the holes are accessible + via the ``holes`` property. + """ + + def __init__( + self, + exterior_boundary: PointLikeList, + holes: Sequence[PointLikeList] | None = None, + name: str | None = None, + ) -> None: + """Create a new region of interest (RoI). + + Parameters + ---------- + exterior_boundary : tuple of (x, y) pairs + The points (in sequence) that make up the boundary of the region. + At least three points must be provided. + holes : sequence of sequences of (x, y) pairs, default None + A sequence of items, where each item will be interpreted as the + ``exterior_boundary`` of an internal hole within the region. See + the ``holes`` argument to ``shapely.Polygon`` for details. + name : str, optional + Name of the RoI that is to be created. A default name will be + inherited from the base class if not provided. + + See Also + -------- + movement.roi.base.BaseRegionOfInterest : The base class that + constructor arguments are passed to, and defaults are inherited + from. + + """ + super().__init__( + points=exterior_boundary, dimensions=2, holes=holes, name=name + ) + + @property + def exterior_boundary(self) -> LineOfInterest: + """The exterior boundary of this RoI.""" + return LineOfInterest( + self.region.exterior.coords, + loop=True, + name=f"Exterior boundary of {self.name}", + ) + + @property + def holes(self) -> tuple[PolygonOfInterest, ...]: + """The interior holes of this RoI. + + Holes are regions properly contained within the exterior boundary of + the RoI that are not part of the RoI itself (like the centre of a + doughnut, for example). A region with no holes returns the empty tuple. + """ + return tuple( + PolygonOfInterest( + int_boundary.coords, name=f"Hole {i} of {self.name}" + ) + for i, int_boundary in enumerate(self.region.interiors) + ) + + @property + def interior_boundaries(self) -> tuple[LineOfInterest, ...]: + """The interior boundaries of this RoI. + + Interior boundaries are the boundaries of holes contained within the + polygon. A region with no holes returns the empty tuple. + """ + return tuple( + LineOfInterest( + int_boundary.coords, + loop=True, + name=f"Interior boundary {i} of {self.name}", + ) + for i, int_boundary in enumerate(self.region.interiors) + ) diff --git a/pyproject.toml b/pyproject.toml index 27348c291..911a4e4e4 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -19,6 +19,7 @@ dependencies = [ "attrs", "pooch", "tqdm", + "shapely", "sleap-io", "xarray[accel,viz]", "PyYAML", diff --git a/tests/test_unit/test_roi/conftest.py b/tests/test_unit/test_roi/conftest.py new file mode 100644 index 000000000..8ba15cb94 --- /dev/null +++ b/tests/test_unit/test_roi/conftest.py @@ -0,0 +1,21 @@ +import numpy as np +import pytest + + +@pytest.fixture() +def unit_square_pts() -> np.ndarray: + return np.array( + [ + [0.0, 0.0], + [1.0, 0.0], + [1.0, 1.0], + [0.0, 1.0], + ], + dtype=float, + ) + + +@pytest.fixture() +def unit_square_hole(unit_square_pts: np.ndarray) -> np.ndarray: + """Hole in the shape of a 0.5 side-length square centred on (0.5, 0.5).""" + return 0.25 + (unit_square_pts.copy() * 0.5) diff --git a/tests/test_unit/test_roi/test_instantiate.py b/tests/test_unit/test_roi/test_instantiate.py new file mode 100644 index 000000000..4d8f7cb10 --- /dev/null +++ b/tests/test_unit/test_roi/test_instantiate.py @@ -0,0 +1,108 @@ +import re +from typing import Any + +import numpy as np +import pytest +import shapely + +from movement.roi.base import BaseRegionOfInterest + + +@pytest.mark.parametrize( + ["input_pts", "kwargs_for_creation", "expected_results"], + [ + pytest.param( + "unit_square_pts", + {"dimensions": 2, "closed": False}, + {"is_closed": True, "dimensions": 2, "name": "Un-named region"}, + id="Polygon, closed is ignored", + ), + pytest.param( + "unit_square_pts", + {"dimensions": 1, "closed": False}, + {"is_closed": False, "dimensions": 1}, + id="Line segment(s)", + ), + pytest.param( + "unit_square_pts", + {"dimensions": 1, "closed": True}, + {"is_closed": True, "dimensions": 1}, + id="Looped lines", + ), + pytest.param( + "unit_square_pts", + {"dimensions": 2, "name": "elephant"}, + {"is_closed": True, "dimensions": 2, "name": "elephant"}, + id="Explicit name", + ), + pytest.param( + np.array([[0.0, 0.0], [1.0, 0.0]]), + {"dimensions": 2}, + ValueError("Need at least 3 points to define a 2D region (got 2)"), + id="Too few points (2D)", + ), + pytest.param( + np.array([[0.0, 0.0]]), + {"dimensions": 1}, + ValueError("Need at least 2 points to define a 1D region (got 1)"), + id="Too few points (1D)", + ), + pytest.param( + np.array([[0.0, 0.0], [1.0, 0.0]]), + {"dimensions": 1}, + {"is_closed": False}, + id="Borderline enough points (1D)", + ), + pytest.param( + np.array([[0.0, 0.0], [1.0, 0.0]]), + {"dimensions": 1, "closed": True}, + ValueError("Cannot create a loop from a single line segment."), + id="Cannot close single line segment.", + ), + pytest.param( + "unit_square_pts", + {"dimensions": 3, "closed": False}, + ValueError( + "Only regions of interest of dimension 1 or 2 " + "are supported (requested 3)" + ), + id="Bad dimensionality", + ), + ], +) +def test_creation( + input_pts, + kwargs_for_creation: dict[str, Any], + expected_results: dict[str, Any] | Exception, + request, +) -> None: + if isinstance(input_pts, str): + input_pts = request.getfixturevalue(input_pts) + + if isinstance(expected_results, Exception): + with pytest.raises( + type(expected_results), match=re.escape(str(expected_results)) + ): + BaseRegionOfInterest(input_pts, **kwargs_for_creation) + else: + roi = BaseRegionOfInterest(input_pts, **kwargs_for_creation) + + expected_dim = kwargs_for_creation.pop("dimensions", 2) + expected_closure = kwargs_for_creation.pop("closed", False) + if expected_dim == 2: + assert isinstance(roi.region, shapely.Polygon) + assert len(roi.coords) == len(input_pts) + 1 + string_should_contain = "-gon" + elif expected_closure: + assert isinstance(roi.region, shapely.LinearRing) + assert len(roi.coords) == len(input_pts) + 1 + string_should_contain = "line segment(s)" + else: + assert isinstance(roi.region, shapely.LineString) + assert len(roi.coords) == len(input_pts) + string_should_contain = "line segment(s)" + assert string_should_contain in roi.__str__() + assert string_should_contain in roi.__repr__() + + for attribute_name, expected_value in expected_results.items(): + assert getattr(roi, attribute_name) == expected_value diff --git a/tests/test_unit/test_roi/test_polygon_boundary.py b/tests/test_unit/test_roi/test_polygon_boundary.py new file mode 100644 index 000000000..0ecd1ee0f --- /dev/null +++ b/tests/test_unit/test_roi/test_polygon_boundary.py @@ -0,0 +1,64 @@ +import numpy as np +import pytest +import shapely + +from movement.roi.line import LineOfInterest +from movement.roi.polygon import PolygonOfInterest + + +@pytest.mark.parametrize( + ["exterior_boundary", "interior_boundaries"], + [ + pytest.param("unit_square_pts", tuple(), id="No holes"), + pytest.param( + "unit_square_pts", tuple(["unit_square_hole"]), id="One hole" + ), + pytest.param( + "unit_square_pts", + ( + np.array([[0.0, 0.0], [0.25, 0.0], [0.0, 0.25]]), + np.array([[0.75, 0.0], [1.0, 0.0], [1.0, 0.25]]), + ), + id="Corners shaved off", + ), + ], +) +def test_boundary(exterior_boundary, interior_boundaries, request) -> None: + if isinstance(exterior_boundary, str): + exterior_boundary = request.getfixturevalue(exterior_boundary) + interior_boundaries = tuple( + request.getfixturevalue(ib) if isinstance(ib, str) else ib + for ib in interior_boundaries + ) + tolerance = 1.0e-8 + + polygon = PolygonOfInterest( + exterior_boundary, holes=interior_boundaries, name="Holey" + ) + expected_exterior = shapely.LinearRing(exterior_boundary) + expected_interiors = tuple( + shapely.LinearRing(ib) for ib in interior_boundaries + ) + expected_holes = tuple(shapely.Polygon(ib) for ib in interior_boundaries) + + computed_exterior = polygon.exterior_boundary + computed_interiors = polygon.interior_boundaries + computed_holes = polygon.holes + + assert isinstance(computed_exterior, LineOfInterest) + assert expected_exterior.equals_exact(computed_exterior.region, tolerance) + assert isinstance(computed_interiors, tuple) + assert isinstance(computed_holes, tuple) + assert len(computed_interiors) == len(expected_interiors) + assert len(computed_holes) == len(expected_holes) + assert len(computed_holes) == len(computed_interiors) + for i, interior_line in enumerate(computed_interiors): + assert isinstance(interior_line, LineOfInterest) + + assert expected_interiors[i].equals_exact( + interior_line.region, tolerance + ) + for i, interior_hole in enumerate(computed_holes): + assert isinstance(interior_hole, PolygonOfInterest) + + assert expected_holes[i].equals_exact(interior_hole.region, tolerance)