diff --git a/README.md b/README.md index 0891002..791bfb2 100644 --- a/README.md +++ b/README.md @@ -1,6 +1,6 @@ # `sgeop`: Street Geometry Processing Toolkit -[![Continuous Integration](https://github.com/uscuni/sgeop/actions/workflows/testing.yml/badge.svg)](https://github.com/uscuni/sgeop/actions/workflows/testing.yml) +[![Continuous Integration](https://github.com/uscuni/sgeop/actions/workflows/testing.yml/badge.svg)](https://github.com/uscuni/sgeop/actions/workflows/testing.yml) [![codecov](https://codecov.io/gh/uscuni/sgeop/branch/main/graph/badge.svg?token=VNn0WR5JWT)](https://codecov.io/gh/uscuni/sgeop) Street geometry processing toolkit @@ -15,7 +15,6 @@ This package developed & and maintained by: * [Anastassia Vybornova](https://github.com/anastassiavybornova) * [James D. Gaboardi](https://github.com/jGaboardi) - ## Introduction ## Documentation diff --git a/codecov.yml b/codecov.yml new file mode 100644 index 0000000..a79ddc6 --- /dev/null +++ b/codecov.yml @@ -0,0 +1,22 @@ +codecov: + notify: + after_n_builds: 6 +coverage: + range: 50..95 + round: nearest + precision: 1 + status: + project: + default: + threshold: 2% + patch: + default: + threshold: 2% + target: 80% + ignore: + - "tests/*" +comment: + layout: "reach, diff, files" + behavior: once + after_n_builds: 6 + require_changes: true diff --git a/sgeop/geometry.py b/sgeop/geometry.py index 5461171..be94906 100644 --- a/sgeop/geometry.py +++ b/sgeop/geometry.py @@ -13,9 +13,7 @@ from .nodes import consolidate_nodes -def _is_within( - line: np.ndarray[shapely.Geometry], poly: shapely.Polygon, rtol: float = 1e-4 -) -> bool: +def _is_within(line: np.ndarray, poly: shapely.Polygon, rtol: float = 1e-4) -> bool: """Check if the line is within a polygon with a set relative tolerance. Parameters @@ -97,14 +95,14 @@ def angle_between_two_lines( def voronoi_skeleton( - lines, - poly=None, - snap_to=None, - max_segment_length=1, - buffer=None, - secondary_snap_to=None, - limit_distance=2, - consolidation_tolerance=None, + lines: list | np.ndarray | gpd.GeoSeries, + poly: None | shapely.Polygon = None, + snap_to: None | gpd.GeoSeries = None, + max_segment_length: int = 1, + buffer: None | float | int = None, + secondary_snap_to: None | gpd.GeoSeries = None, + limit_distance: None | int = 2, + consolidation_tolerance: None | float = None, ): """ Returns average geometry. @@ -112,18 +110,24 @@ def voronoi_skeleton( Parameters ---------- lines : array_like - LineStrings connected at endpoints - poly : shapely.geometry.Polygon - polygon enclosed by `lines` - snap_to : gpd.GeoSeries - series of geometries that shall be connected to the skeleton - distance : float - distance for interpolation - buffer : float - optional custom buffer distance for dealing with Voronoi infinity issues - consolidation_tolerance : float - tolerance passed to node consolidation within the resulting skeleton. If None, - no consolidation happens + LineStrings connected at endpoints. If ``poly`` is passed in, ``lines`` + must be a ``geopandas.GeoSeries``. + poly : None | shapely.Polygon = None + Polygon enclosed by ``lines``. + snap_to : None | gpd.GeoSeries = None + Series of geometries that shall be connected to the skeleton. + max_segment_length: int = 1 + Additional vertices will be added so that all line segments + are no longer than this value. Must be greater than 0. + buffer : None | float | int = None + Optional custom buffer distance for dealing with Voronoi infinity issues. + secondary_snap_to : None | gpd.GeoSeries = None + ... + limit_distance : None | int = 2 + ... + consolidation_tolerance : None | float = None + Tolerance passed to node consolidation within the resulting skeleton. + If ``None``, no consolidation happens. Returns ------- @@ -133,6 +137,8 @@ def voronoi_skeleton( if buffer is None: buffer = max_segment_length * 20 if not poly: + if not isinstance(lines, gpd.GeoSeries): + lines = gpd.GeoSeries(lines) poly = shapely.box(*lines.total_bounds) # get an additional line around the lines to avoid infinity issues with Voronoi extended_lines = list(lines) + [poly.buffer(buffer).boundary] @@ -196,9 +202,7 @@ def voronoi_skeleton( edgeline = shapely.intersection(edgeline, limit) # in edge cases, this can result in a MultiLineString with one sliver part - if edgeline.geom_type == "MultiLineString": - parts = shapely.get_parts(edgeline) - edgeline = parts[np.argmax(shapely.length(parts))] + edgeline = _remove_sliver(edgeline) # check if a, b lines share a node intersection = shapely_lines[b].intersection(shapely_lines[a]) @@ -258,15 +262,38 @@ def voronoi_skeleton( edgelines = edgelines[edgelines != None] # noqa: E711 edgelines = shapely.line_merge(edgelines[shapely.length(edgelines) > 0]) + edgelines = _as_parts(edgelines) + edgelines = _consolidate(edgelines, consolidation_tolerance) + + return edgelines, splitters + + +def _remove_sliver( + edgeline: shapely.LineString | shapely.MultiLineString, +) -> shapely.LineString: + """Remove sliver(s) if present.""" + if edgeline.geom_type == "MultiLineString": + parts = shapely.get_parts(edgeline) + edgeline = parts[np.argmax(shapely.length(parts))] + return edgeline + + +def _as_parts(edgelines: np.ndarray) -> np.ndarray: + """Return constituent LineStrings if MultiLineString present.""" if np.unique(shapely.get_type_id(edgelines)).shape[0] > 1: edgelines = shapely.get_parts(edgelines) + return edgelines + +def _consolidate( + edgelines: np.ndarray, consolidation_tolerance: float | int +) -> np.ndarray: + """Return ``edgelines`` from consolidated nodes, if criteria met.""" if consolidation_tolerance and edgelines.shape[0] > 0: edgelines = consolidate_nodes( edgelines, tolerance=consolidation_tolerance, preserve_ends=True ).geometry.to_numpy() - - return edgelines, splitters + return edgelines def snap_to_targets(edgelines, poly, snap_to, secondary_snap_to=None): diff --git a/sgeop/tests/conftest.py b/sgeop/tests/conftest.py new file mode 100644 index 0000000..8eafde2 --- /dev/null +++ b/sgeop/tests/conftest.py @@ -0,0 +1,66 @@ +import geopandas.testing +import numpy +import pandas +import pytest +import shapely + +line_collection = ( + list[shapely.LineString] + | tuple[shapely.LineString] + | numpy.ndarray + | pandas.Series + | geopandas.GeoSeries +) + +geometry_collection = ( + list[shapely.GeometryCollection] + | tuple[shapely.GeometryCollection] + | numpy.ndarray + | pandas.Series + | geopandas.GeoSeries +) + + +def polygonize( + collection: line_collection, as_geom: bool = True +) -> shapely.Polygon | geopandas.GeoSeries: + """Testing helper -- Create polygon from collection of lines.""" + if isinstance(collection, pandas.Series | geopandas.GeoSeries): + _poly = geopandas.GeoSeries(collection).polygonize() + if as_geom: + return _poly.squeeze() + else: + return _poly + else: + return shapely.polygonize(collection).buffer(0) + + +def is_geopandas(collection: geometry_collection) -> bool: + return isinstance(collection, geopandas.GeoSeries | geopandas.GeoDataFrame) + + +def geom_test( + collection1: geometry_collection, + collection2: geometry_collection, + tolerance: float = 1e-1, +) -> bool: + """Testing helper -- geometry verification.""" + + if not is_geopandas(collection1): + collection1 = geopandas.GeoSeries(collection1) + + if not is_geopandas(collection2): + collection2 = geopandas.GeoSeries(collection2) + + assert shapely.equals_exact( + collection1.geometry.normalize(), + collection2.geometry.normalize(), + tolerance=tolerance, + ).all() + + +def pytest_configure(config): # noqa: ARG001 + """PyTest session attributes, methods, etc.""" + + pytest.polygonize = polygonize + pytest.geom_test = geom_test diff --git a/sgeop/tests/test_geometry.py b/sgeop/tests/test_geometry.py index 43b855b..3765f8d 100644 --- a/sgeop/tests/test_geometry.py +++ b/sgeop/tests/test_geometry.py @@ -1,3 +1,6 @@ +import geopandas.testing +import numpy +import pandas import pytest import shapely @@ -73,7 +76,7 @@ def test_not_within(self): assert known == observed -class TestAngelBetween2Lines: +class TestAngleBetween2Lines: def setup_method(self): self.line1 = shapely.LineString(((0, 0), (1, 0))) self.line2 = shapely.LineString(((1, 0), (1, 1))) @@ -116,3 +119,109 @@ def test_not_adjacent(self): ): observed = sgeop.geometry.angle_between_two_lines(self.line1, self.line4) assert observed == known + + +voronoi_skeleton_params = pytest.mark.parametrize( + "lines_type,as_poly,buffer", + [ + (list, False, None), + (list, True, 0.001), + (numpy.array, False, 0.01), + (numpy.array, True, 0.1), + (pandas.Series, False, 1), + (pandas.Series, True, 2.0), + (geopandas.GeoSeries, False, 5), + (geopandas.GeoSeries, True, 10.314), + ], +) + + +class TestVoronoiSkeleton: + def setup_method(self): + self.square = [ + shapely.LineString(((0, 0), (1000, 0))), + shapely.LineString(((1000, 0), (1000, 1000))), + shapely.LineString(((0, 0), (0, 1000))), + shapely.LineString(((0, 1000), (1000, 1000))), + ] + self.known_square_skeleton_edges = numpy.array( + [ + shapely.LineString(((1000, 0), (998, 2), (500, 500))), + shapely.LineString(((0, 0), (2, 2), (500, 500))), + shapely.LineString(((1000, 1000), (998, 998), (500, 500))), + shapely.LineString(((0, 1000), (2, 998), (500, 500))), + ] + ) + self.known_square_skeleton_splits = [shapely.Point(0, 0)] + self.known_square_skeleton_splits_snap_to = [ + shapely.Point(1000, 0), + shapely.Point(0, 0), + shapely.Point(0, 1000), + shapely.Point(1000, 1000), + ] + + @voronoi_skeleton_params + def test_square(self, lines_type, as_poly, buffer): + known_edges = self.known_square_skeleton_edges + known_splits = self.known_square_skeleton_splits + + lines = lines_type(self.square) + poly = pytest.polygonize(lines) if as_poly else None + observed_edges, observed_splits = sgeop.geometry.voronoi_skeleton( + lines, + poly=poly, + buffer=buffer, + ) + + pytest.geom_test(observed_edges, known_edges) + pytest.geom_test(observed_splits, known_splits) + + @voronoi_skeleton_params + def test_square_snap_to(self, lines_type, as_poly, buffer): + known_edges = self.known_square_skeleton_edges + known_splits = self.known_square_skeleton_splits_snap_to + + lines = lines_type(self.square) + poly = pytest.polygonize(lines) if as_poly else None + observed_edges, observed_splits = sgeop.geometry.voronoi_skeleton( + lines, + poly=poly, + buffer=buffer, + snap_to=( + pytest.polygonize(geopandas.GeoSeries(lines), as_geom=False) + .extract_unique_points() + .explode() + ), + ) + + pytest.geom_test(observed_edges, known_edges) + pytest.geom_test(observed_splits, known_splits) + + +line_100_900 = shapely.LineString(((1000, 1000), (1000, 9000))) +line_100_120 = shapely.LineString(((1000, 1020), (1020, 1020))) +lines_100_900_100_120 = shapely.MultiLineString((line_100_900, line_100_120)) +line_110_900 = shapely.LineString(((1000, 9000), (1100, 9000))) + + +def test_remove_sliver(): + known = line_100_900 + observed = sgeop.geometry._remove_sliver(lines_100_900_100_120) + assert observed == known + + +def test_as_parts(): + known = numpy.array([line_100_900, line_100_120, line_110_900]) + observed = sgeop.geometry._as_parts( + numpy.array([lines_100_900_100_120, line_110_900]) + ) + numpy.testing.assert_array_equal(observed, known) + + +@pytest.mark.parametrize("tolerance", [0.1, 1, 10, 100, 1_000, 10_000, 100_000]) +def test_consolidate(tolerance): + known = numpy.array([line_100_900, line_100_120, line_110_900]) + observed = sgeop.geometry._consolidate( + numpy.array([line_100_900, line_100_120, line_110_900]), tolerance + ) + numpy.testing.assert_array_equal(observed, known) diff --git a/sgeop/tests/test_simplify.py b/sgeop/tests/test_simplify.py index c99d9cd..62bdb01 100644 --- a/sgeop/tests/test_simplify.py +++ b/sgeop/tests/test_simplify.py @@ -1,7 +1,7 @@ import pathlib import geopandas -import shapely +import pytest from pandas.testing import assert_series_equal import sgeop @@ -17,7 +17,5 @@ def test_simplify(): geopandas.read_parquet(test_data / f"{ac}_original.parquet") ) - assert shapely.equals_exact( - known.geometry.normalize(), observed.geometry.normalize(), tolerance=1e-1 - ).all() + pytest.geom_test(known, observed) assert_series_equal(known._status, observed._status)