diff --git a/eolearn/core/utils/testing.py b/eolearn/core/utils/testing.py index ae1cf60d..340267cb 100644 --- a/eolearn/core/utils/testing.py +++ b/eolearn/core/utils/testing.py @@ -23,7 +23,7 @@ from ..constants import FeatureType from ..eodata import EOPatch -from ..types import Feature +from ..types import FeaturesSpecification from ..utils.parsing import FeatureParser DEFAULT_BBOX = BBox((0, 0, 100, 100), crs=CRS("EPSG:32633")) @@ -46,7 +46,7 @@ def __post_init__(self) -> None: def generate_eopatch( - features: list[Feature] | None = None, + features: FeaturesSpecification | None = None, bbox: BBox = DEFAULT_BBOX, timestamps: list[dt.datetime] | None = None, seed: int = 42, diff --git a/eolearn/geometry/morphology.py b/eolearn/geometry/morphology.py index 8939091d..dbf039fe 100644 --- a/eolearn/geometry/morphology.py +++ b/eolearn/geometry/morphology.py @@ -10,11 +10,10 @@ import itertools as it from enum import Enum -from typing import Callable +from functools import partial +import cv2 import numpy as np -import skimage.filters.rank -import skimage.morphology from eolearn.core import EOPatch, EOTask, MapFeatureTask from eolearn.core.types import FeaturesSpecification, SingleFeatureSpec @@ -44,7 +43,7 @@ def __init__( parsed_mask_feature = parse_renamed_feature(mask_feature, allowed_feature_types=lambda fty: fty.is_array()) self.mask_type, self.mask_name, self.new_mask_name = parsed_mask_feature - self.disk = skimage.morphology.disk(disk_radius) + self.disk = cv2.getStructuringElement(cv2.MORPH_ELLIPSE, (disk_radius, disk_radius)) self.erode_labels = erode_labels self.no_data_label = no_data_label @@ -57,7 +56,7 @@ def execute(self, eopatch: EOPatch) -> EOPatch: erode_labels = set(erode_labels) - {self.no_data_label} other_labels = set(all_labels) - set(erode_labels) - {self.no_data_label} - eroded_masks = [skimage.morphology.binary_erosion(feature_array == label, self.disk) for label in erode_labels] + eroded_masks = [cv2.erode((feature_array == label).astype(np.uint8), self.disk) for label in erode_labels] other_masks = [feature_array == label for label in other_labels] merged_mask = np.logical_or.reduce(eroded_masks + other_masks, axis=0) @@ -75,20 +74,18 @@ class MorphologicalOperations(Enum): CLOSING = "closing" DILATION = "dilation" EROSION = "erosion" - MEDIAN = "median" @classmethod - def get_operation(cls, morph_type: MorphologicalOperations) -> Callable: + def get_operation(cls, morph_type: MorphologicalOperations) -> int: """Maps morphological operation type to function :param morph_type: Morphological operation type """ return { - cls.OPENING: skimage.morphology.opening, - cls.CLOSING: skimage.morphology.closing, - cls.DILATION: skimage.morphology.dilation, - cls.EROSION: skimage.morphology.erosion, - cls.MEDIAN: skimage.filters.rank.median, + cls.OPENING: cv2.MORPH_OPEN, + cls.CLOSING: cv2.MORPH_CLOSE, + cls.DILATION: cv2.MORPH_DILATE, + cls.EROSION: cv2.MORPH_ERODE, }[morph_type] @@ -103,15 +100,7 @@ def get_disk(radius: int) -> np.ndarray: :param radius: Radius of disk :return: The structuring element where elements of the neighborhood are 1 and 0 otherwise. """ - return skimage.morphology.disk(radius) - - @staticmethod - def get_diamond(radius: int) -> np.ndarray: - """ - :param radius: Radius of diamond - :return: The structuring element where elements of the neighborhood are 1 and 0 otherwise. - """ - return skimage.morphology.diamond(radius) + return cv2.getStructuringElement(cv2.MORPH_ELLIPSE, (radius, radius)) @staticmethod def get_rectangle(width: int, height: int) -> np.ndarray: @@ -120,7 +109,7 @@ def get_rectangle(width: int, height: int) -> np.ndarray: :param height: Height of rectangle :return: A structuring element consisting only of ones, i.e. every pixel belongs to the neighborhood. """ - return skimage.morphology.rectangle(width, height) + return cv2.getStructuringElement(cv2.MORPH_RECT, (height, width)) @staticmethod def get_square(width: int) -> np.ndarray: @@ -128,7 +117,7 @@ def get_square(width: int) -> np.ndarray: :param width: Size of square :return: A structuring element consisting only of ones, i.e. every pixel belongs to the neighborhood. """ - return skimage.morphology.square(width) + return cv2.getStructuringElement(cv2.MORPH_RECT, (width, width)) class MorphologicalFilterTask(MapFeatureTask): @@ -139,7 +128,7 @@ def __init__( input_features: FeaturesSpecification, output_features: FeaturesSpecification | None = None, *, - morph_operation: MorphologicalOperations | Callable, + morph_operation: MorphologicalOperations, struct_elem: np.ndarray | None = None, ): """ @@ -153,21 +142,19 @@ def __init__( output_features = input_features super().__init__(input_features, output_features) - if isinstance(morph_operation, MorphologicalOperations): - self.morph_operation = MorphologicalOperations.get_operation(morph_operation) - else: - self.morph_operation = morph_operation + self.morph_operation = MorphologicalOperations.get_operation(morph_operation) self.struct_elem = struct_elem def map_method(self, feature: np.ndarray) -> np.ndarray: """Applies the morphological operation to a raster feature.""" feature = feature.copy() + morph_func = partial(cv2.morphologyEx, kernel=self.struct_elem, op=self.morph_operation) if feature.ndim == 3: for channel in range(feature.shape[2]): - feature[..., channel] = self.morph_operation(feature[..., channel], self.struct_elem) + feature[..., channel] = morph_func(feature[..., channel]) elif feature.ndim == 4: for time, channel in it.product(range(feature.shape[0]), range(feature.shape[3])): - feature[time, ..., channel] = self.morph_operation(feature[time, ..., channel], self.struct_elem) + feature[time, ..., channel] = morph_func(feature[time, ..., channel]) else: raise ValueError(f"Invalid number of dimensions: {feature.ndim}") diff --git a/tests/geometry/test_morphology.py b/tests/geometry/test_morphology.py index 74852c4e..646e80ec 100644 --- a/tests/geometry/test_morphology.py +++ b/tests/geometry/test_morphology.py @@ -10,9 +10,8 @@ import pytest from numpy.testing import assert_array_equal -from sentinelhub import CRS, BBox - from eolearn.core import EOPatch, FeatureType +from eolearn.core.utils.testing import PatchGeneratorConfig, generate_eopatch from eolearn.geometry import ErosionTask, MorphologicalFilterTask, MorphologicalOperations, MorphologicalStructFactory CLASSES = (0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10) @@ -21,6 +20,15 @@ # ruff: noqa: NPY002 +@pytest.fixture(name="patch") +def patch_fixture() -> EOPatch: + config = PatchGeneratorConfig(max_integer_value=10, raster_shape=(50, 100), depth_range=(3, 4)) + patch = generate_eopatch([MASK_FEATURE, MASK_TIMELESS_FEATURE], config=config) + for feat in [MASK_FEATURE, MASK_TIMELESS_FEATURE]: + patch[feat] = patch[feat].astype(np.uint8) + return patch + + @pytest.mark.parametrize("invalid_input", [None, 0, "a"]) def test_erosion_value_error(invalid_input): with pytest.raises(ValueError): @@ -28,58 +36,60 @@ def test_erosion_value_error(invalid_input): def test_erosion_full(test_eopatch): - mask_before = test_eopatch.mask_timeless["LULC"].copy() - - erosion_task = ErosionTask((FeatureType.MASK_TIMELESS, "LULC", "LULC_ERODED"), 1) + erosion_task = ErosionTask((FeatureType.MASK_TIMELESS, "LULC", "LULC_ERODED"), disk_radius=3) eopatch = erosion_task.execute(test_eopatch) - mask_after = eopatch.mask_timeless["LULC_ERODED"].copy() + mask_after = eopatch.mask_timeless["LULC_ERODED"] - assert not np.all(mask_before == mask_after) - - for label in CLASSES: - if label == 0: - assert np.sum(mask_after == label) >= np.sum(mask_before == label), "Error in the erosion process" - else: - assert np.sum(mask_after == label) <= np.sum(mask_before == label), "Error in the erosion process" + assert_array_equal(np.unique(mask_after, return_counts=True)[1], [1942, 6950, 1069, 87, 52]) def test_erosion_partial(test_eopatch): - mask_before = test_eopatch.mask_timeless["LULC"].copy() - # skip forest and artificial surface specific_labels = [0, 1, 3, 4] erosion_task = ErosionTask( - mask_feature=(FeatureType.MASK_TIMELESS, "LULC", "LULC_ERODED"), disk_radius=1, erode_labels=specific_labels + mask_feature=(FeatureType.MASK_TIMELESS, "LULC", "LULC_ERODED"), disk_radius=3, erode_labels=specific_labels ) eopatch = erosion_task.execute(test_eopatch) - mask_after = eopatch.mask_timeless["LULC_ERODED"].copy() + mask_after = eopatch.mask_timeless["LULC_ERODED"] - assert not np.all(mask_before == mask_after) + assert_array_equal(np.unique(mask_after, return_counts=True)[1], [1145, 7601, 1069, 87, 198]) - for label in CLASSES: - if label == 0: - assert np.sum(mask_after == label) >= np.sum(mask_before == label), "Error in the erosion process" - elif label in specific_labels: - assert np.sum(mask_after == label) <= np.sum(mask_before == label), "Error in the erosion process" - else: - assert_array_equal(mask_after == label, mask_before == label, err_msg="Error in the erosion process") - -@pytest.mark.parametrize("morph_operation", MorphologicalOperations) @pytest.mark.parametrize( - "struct_element", [None, MorphologicalStructFactory.get_disk(5), MorphologicalStructFactory.get_rectangle(5, 6)] + ("morph_operation", "struct_element", "mask_counts", "mask_timeless_counts"), + [ + ( + MorphologicalOperations.DILATION, + None, + [6, 34, 172, 768, 2491, 7405, 19212, 44912], + [1, 2, 16, 104, 466, 1490, 3870, 9051], + ), + (MorphologicalOperations.EROSION, MorphologicalStructFactory.get_disk(11), [74957, 42, 1], [14994, 6]), + (MorphologicalOperations.OPENING, MorphologicalStructFactory.get_disk(11), [73899, 1051, 50], [14837, 163]), + (MorphologicalOperations.CLOSING, MorphologicalStructFactory.get_disk(11), [770, 74230], [425, 14575]), + ( + MorphologicalOperations.OPENING, + MorphologicalStructFactory.get_rectangle(5, 6), + [48468, 24223, 2125, 169, 15], + [10146, 4425, 417, 3, 9], + ), + ( + MorphologicalOperations.DILATION, + MorphologicalStructFactory.get_rectangle(5, 6), + [2, 19, 198, 3929, 70852], + [32, 743, 14225], + ), + ], ) -def test_morphological_filter(morph_operation, struct_element): - eopatch = EOPatch(bbox=BBox((0, 0, 1, 1), CRS(3857)), timestamps=["2015-7-7"] * 10) - eopatch[MASK_FEATURE] = np.random.randint(20, size=(10, 100, 100, 3), dtype=np.uint8) - eopatch[MASK_TIMELESS_FEATURE] = np.random.randint(20, 50, size=(100, 100, 5), dtype=np.uint8) - +def test_morphological_filter(patch, morph_operation, struct_element, mask_counts, mask_timeless_counts): task = MorphologicalFilterTask( [MASK_FEATURE, MASK_TIMELESS_FEATURE], morph_operation=morph_operation, struct_elem=struct_element ) - task.execute(eopatch) + task.execute(patch) - assert eopatch[MASK_FEATURE].shape == (10, 100, 100, 3) - assert eopatch[MASK_TIMELESS_FEATURE].shape == (100, 100, 5) + assert patch[MASK_FEATURE].shape == (5, 50, 100, 3) + assert patch[MASK_TIMELESS_FEATURE].shape == (50, 100, 3) + assert_array_equal(np.unique(patch[MASK_FEATURE], return_counts=True)[1], mask_counts) + assert_array_equal(np.unique(patch[MASK_TIMELESS_FEATURE], return_counts=True)[1], mask_timeless_counts)