diff --git a/momepy/elements.py b/momepy/elements.py index c0abe3da..6eeced5f 100644 --- a/momepy/elements.py +++ b/momepy/elements.py @@ -1,4 +1,7 @@ +import math import warnings +from collections import Counter, deque +from enum import Enum import geopandas as gpd import numpy as np @@ -153,6 +156,9 @@ def enclosed_tessellation( threshold: float = 0.05, simplify: bool = True, n_jobs: int = -1, + inner_barriers: GeoSeries | GeoDataFrame = None, + cell_size: float = 1, + neighbor_mode: str = "moore", **kwargs, ) -> GeoDataFrame: """Generate enclosed tessellation @@ -163,6 +169,11 @@ def enclosed_tessellation( street network, railway). Original morphological tessellation is used under the hood to partition each enclosure. + When ``inner_barriers`` are provided, the tessellation is derived using a Cellular + Automata implementation that recognizes dangling barriers (such as dead-end streets) + as valid limits of cell growth. This is more computationally intensive but handles + complex barrier configurations more accurately. + Tessellation requires data of relatively high level of precision and there are three particular patterns causing issues: @@ -206,9 +217,27 @@ def enclosed_tessellation( n_jobs : int, optional The number of jobs to run in parallel. -1 means using all available cores. By default -1 + inner_barriers: GeoSeries | GeoDataFrame, optional + Barriers that should be included in the tessellation process. When provided, + tessellation will be derived using a Cellular Automata implementation that + recognizes dangling barriers (such as dead-end streets or cul-de-sacs) as valid + limits of cell growth. This is more computationally intensive than the default + Voronoi-based approach but can handle inner barriers. By default None. + cell_size : float, optional + Grid cell size for the Cellular Automata implementation when ``inner_barriers`` + is not None. Smaller values provide higher resolution but increase computational + cost. This parameter controls the spatial granularity of the tessellation grid. + When ``inner_barriers`` is None, this parameter is ignored. By default 1.0 + neighbor_mode : str, optional + Choice of neighbor connectivity for the Cellular Automata implementation when + ``inner_barriers`` is not None. Options are 'moore' (8-connected, including + diagonal neighbors) or 'neumann' (4-connected, only orthogonal neighbors). + When ``inner_barriers`` is None, this parameter is ignored. By default 'moore'. **kwargs - Additional keyword arguments pased to libpysal.cg.voronoi_frames, such as - ``grid_size``. + Additional keyword arguments passed to :func:`libpysal.cg.voronoi_frames` when + ``inner_barriers`` is None, such as ``grid_size``. These arguments are ignored + when ``inner_barriers`` is provided and the Cellular Automata implementation is + used. Warnings -------- @@ -299,7 +328,18 @@ def enclosed_tessellation( # generate tessellation in parallel new = Parallel(n_jobs=n_jobs)( - delayed(_tess)(*t, threshold, shrink, segment, index_name, simplify, kwargs) + delayed(_tess)( + *t, + threshold, + shrink, + segment, + index_name, + simplify, + inner_barriers, + cell_size, + neighbor_mode, + kwargs, + ) for t in tuples ) @@ -332,8 +372,51 @@ def enclosed_tessellation( return pd.concat([new_df, singles.drop(columns="position"), clean_blocks]) -def _tess(ix, poly, blg, threshold, shrink, segment, enclosure_id, to_simplify, kwargs): - """Generate tessellation for a single enclosure. Helper for enclosed_tessellation""" +def _tess( + ix, + poly, + blg, + threshold, + shrink, + segment, + enclosure_id, + to_simplify, + inner_barriers, + cell_size, + neighbor_mode, + kwargs, +): + """Generate tessellation for a single enclosure. Helper for enclosed_tessellation. + + Parameters + ---------- + ix : int + Enclosure index. + poly : shapely.Geometry + Enclosure geometry. + blg : GeoSeries | GeoDataFrame + Buildings within the enclosure. + threshold : float + Threshold for building inclusion. + shrink : float + Shrink distance for tessellation. + segment : float + Segmentation distance. + enclosure_id : str + Column name for enclosure ID. + inner_barriers : GeoSeries | GeoDataFrame + Inner barriers for tessellation. + cell_size : float + Grid cell size when inner_barriers is not None. + neighbor_mode : str + Choice of neighbor connectivity ('moore' or 'neumann') + when inner_barriers is not None. + + Returns + ------- + GeoDataFrame + Tessellation for the enclosure. + """ # check if threshold is set and filter buildings based on the threshold if threshold: blg = blg[ @@ -342,19 +425,35 @@ def _tess(ix, poly, blg, threshold, shrink, segment, enclosure_id, to_simplify, ] if len(blg) > 1: - tess = voronoi_frames( - blg, - clip=poly, - shrink=shrink, - segment=segment, - return_input=False, - as_gdf=True, - **kwargs, - ) + if inner_barriers is None or inner_barriers.empty: + tess = voronoi_frames( + blg, + clip=poly, + shrink=shrink, + segment=segment, + return_input=False, + as_gdf=True, + **kwargs, + ) + tolerance = segment / 2 + + else: + tess = _voronoi_by_ca( + seed_geoms=blg, + barrier_geoms=poly, + cell_size=cell_size, + neighbor_mode=neighbor_mode, + barriers_for_inner=inner_barriers, + ) + # torelance set as the square root of the isosceles right triangle with 2 + # cells_size edges + tolerance = ((2 * cell_size) ** 2 / 2) ** 0.5 + if to_simplify: tess.geometry = shapely.coverage_simplify( - tess.geometry, tolerance=segment / 2, simplify_boundary=False + tess.geometry, tolerance=tolerance, simplify_boundary=False ) + tess[enclosure_id] = ix return tess @@ -372,6 +471,413 @@ def _tess(ix, poly, blg, threshold, shrink, segment, enclosure_id, to_simplify, ) +def _voronoi_by_ca( + seed_geoms: GeoSeries | GeoDataFrame, + barrier_geoms: GeoSeries | GeoDataFrame, + cell_size: float = 1.0, + neighbor_mode: str = "moore", + barriers_for_inner: GeoSeries | GeoDataFrame = None, +) -> GeoDataFrame: + """ + Generate an aggregated Voronoi tessellation as a GeoDataFrame via cellular automata. + + This unified function performs the following: + + - Ensures that the CRS of seed and barrier geometries are aligned. + - Combines inner barriers with the enclosure. + - Computes grid bounds covering both seed and barrier geometries. + - Marks barrier cells using a prepared geometry for fast intersection. + - Seeds the grid with seed geometries. + - Propagates seed values via a BFS expansion that respects barriers. + - Uses a voting mechanism to finalize boundary cell assignments. + - Converts grid cells into a GeoDataFrame and dissolves adjacent cells + with the same seed id. + - Clips the output cells by the barrier boundaries. + + Parameters + ---------- + seed_geoms : GeoSeries | GeoDataFrame + GeoDataFrame containing seed features. + barrier_geoms : GeoSeries | GeoDataFrame + GeoDataFrame containing barrier features or a shapely Polygon. + cell_size : float, optional + Grid cell size. By default 1.0 + neighbor_mode : str, optional + Choice of neighbor connectivity ('moore' or 'neumann'). By default 'moore'. + barriers_for_inner : GeoSeries | GeoDataFrame, optional + GeoDataFrame containing inner barriers to be included. By default None + + Returns + ------- + GeoDataFrame + A GeoDataFrame representing the aggregated Voronoi tessellation, clipped by + barriers. + """ + geom_type = barrier_geoms.geom_type + if geom_type not in {"Polygon", "MultiPolygon"}: + raise ValueError("Enclosure must be a Polygon or MultiPolygon") + + # Get inner barriers as intersection or containment of the barrier_geoms + inner_barriers = _get_inner_barriers(barrier_geoms, barriers_for_inner) + + # Handle barrier_geoms if it is a Polygon or MultiPolygon + if geom_type == "Polygon": + # Take buffer of polygon and extract its exterior boundary (10 cells) + barrier_geoms_buffered = GeoSeries( + [barrier_geoms.buffer(10 * cell_size).exterior], crs=seed_geoms.crs + ) + barrier_geoms = GeoSeries([barrier_geoms], crs=seed_geoms.crs) + + elif geom_type == "MultiPolygon": + # Process each polygon: take buffer then exterior boundary + # (to ensure there's no gap between enclosures) + parts = list(shapely.get_parts(barrier_geoms)) + exteriors = shapely.get_exterior_ring(parts) + buffered_exteriors = shapely.buffer(exteriors, 10 * cell_size) + barrier_geoms_buffered = GeoSeries( + list(buffered_exteriors), crs=seed_geoms.crs + ) + barrier_geoms = GeoSeries(parts, crs=seed_geoms.crs) + + outer_union = barrier_geoms_buffered.union_all() + inner_union = ( + inner_barriers.union_all() + if inner_barriers is not None and not inner_barriers.empty + else None + ) + + geoms_to_union = [ + geom + for geom in (outer_union, inner_union) + if geom is not None and not geom.is_empty + ] + prep_barrier = shapely.union_all(geoms_to_union) if geoms_to_union else None + + # Compute grid bounds + origin, grid_width, grid_height = _get_grid_bounds( + seed_geoms, barrier_geoms_buffered, cell_size + ) + + # Initialize grid states with UNKNOWN values. + states = np.full((grid_height, grid_width), CellState.UNKNOWN.value, dtype=int) + + xs, ys = np.meshgrid(np.arange(grid_width), np.arange(grid_height)) + xs_flat = xs.ravel() + ys_flat = ys.ravel() + cell_polys_array = shapely.box( + origin[0] + xs_flat * cell_size, + origin[1] + ys_flat * cell_size, + origin[0] + (xs_flat + 1) * cell_size, + origin[1] + (ys_flat + 1) * cell_size, + ) + + # Identify barrier cells in the grid + barrier_mask = np.zeros((grid_height, grid_width), dtype=bool) + if prep_barrier is not None and not prep_barrier.is_empty: + barrier_mask = shapely.intersects(cell_polys_array, prep_barrier).reshape( + grid_height, grid_width + ) + states[barrier_mask] = CellState.BARRIER.value + + # Seed the grid with seed geometries. + for site_id, geom in enumerate(seed_geoms.geometry): + if not geom.is_empty: + cells = _geom_to_cells(geom, origin, cell_size, grid_width, grid_height) + valid_cells = [ + (x, y) for x, y in cells if states[y, x] == CellState.UNKNOWN.value + ] + if valid_cells: + indices = np.array(valid_cells) + states[indices[:, 1], indices[:, 0]] = site_id + + # Initialize the BFS queue with all seeded cells’ neighbors. + # ... existing code ... + + queue: deque[tuple[int, int]] = deque() + seed_indices = np.argwhere(states >= 0) + for y, x in seed_indices: + _enqueue_neighbors(x, y, states, grid_width, grid_height, neighbor_mode, queue) + + # Process BFS to propagate seed values. + while queue: + # Dequeue the current cell and skip if it is not a frontier cell. + x_current, y_current = queue.popleft() + if states[y_current, x_current] != CellState.FRONTIER.value: + continue + + # Get neighbor cells that were already assigned a seed id or still unknown + # (state >= 0). + # Note that boundary or barrier cells are skipped (state < 0). + neighbor_seeds = [ + states[ny, nx] + for nx, ny in _get_neighbors( + x_current, y_current, grid_width, grid_height, mode=neighbor_mode + ) + if states[ny, nx] >= 0 + ] + if not neighbor_seeds: + continue + + # Assign as a boundary if multiple seed ids are found. + if len(set(neighbor_seeds)) > 1: + states[y_current, x_current] = CellState.BOUNDARY.value + # If not, equeue neighbor cells for further propagation. + else: + assigned_seed = set(neighbor_seeds).pop() + states[y_current, x_current] = assigned_seed + _enqueue_neighbors( + x_current, + y_current, + states, + grid_width, + grid_height, + neighbor_mode, + queue, + ) + + # Post-process barrier and boundary cells using a voting mechanism. + states = _assign_adjacent_seed_cells(states, neighbor_mode) + + # Create grid cell polygons and build a GeoDataFrame. + grid_polys = GeoSeries(cell_polys_array, crs=seed_geoms.crs) + grid_gdf = GeoDataFrame( + {"site_id": states.flatten()}, geometry=grid_polys, crs=seed_geoms.crs + ) + + # Include only cells with valid seed assignments and dissolve contiguous regions. + grid_gdf = ( + grid_gdf[grid_gdf["site_id"] >= 0] + .dissolve(by="site_id") + .reset_index() + .drop("site_id", axis=1) + ) + + # Clip by barriers + if barrier_geoms is not None and (not barrier_geoms.empty): + barrier_union = barrier_geoms.union_all() + polygonized = [] + if barrier_union is not None and not barrier_union.is_empty: + polygonized = [ + geom + for geom in shapely.get_parts(shapely.polygonize([barrier_union])) + if geom.geom_type in {"Polygon", "MultiPolygon"} + ] + if polygonized: + barrier_union = shapely.union_all(polygonized) + if barrier_union is not None and not barrier_union.is_empty: + grid_gdf["geometry"] = grid_gdf["geometry"].intersection(barrier_union) + + # cleanup possible collections + def sanitize_geometry(geom): + parts = ( + [ + part + for part in shapely.get_parts(geom) + if part.geom_type in {"Polygon", "MultiPolygon"} + ] + if geom is not None + else [] + ) + return shapely.union_all(parts) if parts else None + + grid_gdf["geometry"] = grid_gdf["geometry"].apply(sanitize_geometry) + grid_gdf = grid_gdf.dropna(subset=["geometry"]) + return grid_gdf + + +def _get_inner_barriers(enclosure, barriers): + """Get inner barriers that intersect or are contained within an enclosure. + + Parameters + ---------- + enclosure : GeoSeries | GeoDataFrame + The enclosure geometry. + barriers : GeoDataFrame + The barriers GeoDataFrame. + + Returns + ------- + shapely.Polygon + A single Polygon combining the enclosure and any intersecting barriers. + """ + # Clip those segments to stay within the enclosure + inner_barriers = gpd.clip(barriers, enclosure) + + # Only keep the geometry which is within the enclosure + inner_barriers = inner_barriers[inner_barriers.intersects(enclosure) | + inner_barriers.within(enclosure)] + + return inner_barriers + + +class CellState(Enum): + """ + Enumeration of cell states for grid processing for improving the readability, + instead of integers. + + Attributes + ---------- + UNKNOWN : int + Cell has not been processed. + BOUNDARY : int + Cell is at a junction between different seed regions. + BARRIER : int + Cell originally designated as a barrier. + FRONTIER : int + Cell queued for BFS expansion. + """ + + UNKNOWN = -1 + BOUNDARY = -2 + BARRIER = -3 + FRONTIER = -4 + + +def _get_neighbors( + x: int, y: int, max_x: int, max_y: int, mode: str = "moore" +) -> list[tuple[int, int]]: + """ + Retrieve valid neighboring cell indices based on connectivity. + + Parameters + ---------- + x : int + The current cell x index. + y : int + The current cell y index. + max_x : int + The maximum x dimension of the grid. + max_y : int + The maximum y dimension of the grid. + mode : str, optional + The connectivity mode, "moore" for 8-connected or "neumann" for 4-connected + neighbors. By default "moore". + + Returns + ------- + list[tuple[int, int]] + A list of (x, y) tuples for valid neighbor indices. + """ + neighbor_dirs = { + "moore": [(-1, -1), (-1, 0), (-1, 1), (0, -1), (0, 1), (1, -1), (1, 0), (1, 1)], + "neumann": [(0, -1), (-1, 0), (1, 0), (0, 1)], + } + directions = neighbor_dirs.get(mode) + if directions is None: + raise ValueError("Invalid neighbor_mode: choose 'moore' or 'neumann'") + return [ + (x + dx, y + dy) + for dx, dy in directions + if 0 <= x + dx < max_x and 0 <= y + dy < max_y + ] + + +def _get_grid_bounds( + seed_geoms: GeoSeries | GeoDataFrame, + barrier_geoms: GeoSeries | GeoDataFrame, + cell_size: float, +) -> tuple[tuple[float, float], int, int]: + seed_bounds = seed_geoms.total_bounds # [xmin, ymin, xmax, ymax] + barrier_bounds = barrier_geoms.total_bounds + + xmin = min(seed_bounds[0], barrier_bounds[0]) + ymin = min(seed_bounds[1], barrier_bounds[1]) + xmax = max(seed_bounds[2], barrier_bounds[2]) + ymax = max(seed_bounds[3], barrier_bounds[3]) + # expand bounds by 1 cell in each direction + new_xmin = xmin - cell_size + new_ymin = ymin - cell_size + new_xmax = xmax + cell_size + new_ymax = ymax + cell_size + grid_width = math.ceil((new_xmax - new_xmin) / cell_size) + grid_height = math.ceil((new_ymax - new_ymin) / cell_size) + return (new_xmin, new_ymin), grid_width, grid_height + + +def _geom_to_cells( + geom: shapely.geometry, + origin: tuple[float, float], + cell_size: float, + grid_width: int, + grid_height: int, +) -> list[tuple[int, int]]: + """ + Determine grid cell indices that intersect the given geometry. + """ + if isinstance(geom, shapely.geometry.Point): + sx = int((geom.x - origin[0]) // cell_size) + sy = int((geom.y - origin[1]) // cell_size) + return [(sx, sy)] if 0 <= sx < grid_width and 0 <= sy < grid_height else [] + + else: + minx, miny, maxx, maxy = geom.bounds + start_x = max(0, int((minx - origin[0]) // cell_size)) + start_y = max(0, int((miny - origin[1]) // cell_size)) + end_x = min(grid_width, int(math.ceil((maxx - origin[0]) / cell_size))) + end_y = min(grid_height, int(math.ceil((maxy - origin[1]) / cell_size))) + + x_range = np.arange(start_x, end_x) + y_range = np.arange(start_y, end_y) + xx, yy = np.meshgrid(x_range, y_range) + xx_flat = xx.ravel() + yy_flat = yy.ravel() + candidate_polys = shapely.box( + origin[0] + xx_flat * cell_size, + origin[1] + yy_flat * cell_size, + origin[0] + (xx_flat + 1) * cell_size, + origin[1] + (yy_flat + 1) * cell_size, + ) + mask = shapely.intersects(candidate_polys, geom) + return list(zip(xx_flat[mask], yy_flat[mask], strict=True)) + + +def _enqueue_neighbors( + x: int, + y: int, + states: np.ndarray, + grid_width: int, + grid_height: int, + neighbor_mode: str, + queue: deque, +) -> None: + """ + Enqueue valid neighboring cells for BFS expansion. + """ + for nx, ny in _get_neighbors(x, y, grid_width, grid_height, mode=neighbor_mode): + if states[ny, nx] == CellState.UNKNOWN.value: + states[ny, nx] = CellState.FRONTIER.value + queue.append((nx, ny)) + + +def _assign_adjacent_seed_cells( + states: np.ndarray, neighbor_mode: str = "moore" +) -> np.ndarray: + """ + Reassign border and barrier cells to the proximate seed areas using a voting + mechanism. + """ + new_states = states.copy() + indices = np.argwhere( + np.isin(states, [CellState.BARRIER.value, CellState.BOUNDARY.value]) + ) + grid_height, grid_width = states.shape + + for y, x in indices: + neighbor_seeds = [ + states[ny, nx] + for nx, ny in _get_neighbors( + x, y, grid_width, grid_height, mode=neighbor_mode + ) + if states[ny, nx] >= 0 + ] + if neighbor_seeds: + cnt = Counter(neighbor_seeds) + # In case of ties, choose the smaller seed id. + chosen_seed = min(cnt.items(), key=lambda item: (-item[1], item[0]))[0] + new_states[y, x] = chosen_seed + return new_states + + def verify_tessellation(tessellation, geometry): """Check whether result matches buildings and contains only Polygons. diff --git a/momepy/tests/test_elements.py b/momepy/tests/test_elements.py index 44d20dda..fcf66c57 100644 --- a/momepy/tests/test_elements.py +++ b/momepy/tests/test_elements.py @@ -10,9 +10,10 @@ from packaging.version import Version from pandas.testing import assert_index_equal from shapely import LineString, affinity -from shapely.geometry import MultiPoint, Polygon, box +from shapely.geometry import MultiPoint, MultiPolygon, Polygon, box import momepy as mm +import momepy.elements as elements SHPLY_GE_210 = Version(shapely.__version__) >= Version("2.1.0") LPS_G_4_13_0 = Version(libpysal.__version__) > Version("4.13.0") @@ -175,6 +176,180 @@ def test_enclosed_tessellation(self): assert (self.df_buildings.index.isin(tessellation_custom_index.index)).all() assert tessellation_custom_index.enclosure_index.isin(custom_index.index).all() + tessellation_inner_barrier = mm.enclosed_tessellation( + self.df_buildings, + self.enclosures.geometry, + simplify=False, + inner_barriers=self.df_streets, + cell_size=5, + neighbor_mode="neumann", + n_jobs=1, + ) + assert set(tessellation_inner_barrier.geom_type.unique()) <= {"Polygon", "MultiPolygon"} + assert tessellation_inner_barrier.crs == self.df_buildings.crs + assert len(tessellation_inner_barrier) == len(tessellation) + assert tessellation_inner_barrier.enclosure_index.isin( + self.enclosures.index + ).all() + + def test_enclosed_tessellation_inner_barrier_cellular(self): + # Define sample building geometries (including Points to cover Point handling). + blg_polygons = [ + shapely.geometry.Polygon([(15, 32), (35, 32), (35, 38), (15, 38)]), + shapely.geometry.Polygon([(15, 22), (35, 22), (35, 28), (15, 28)]), + shapely.geometry.Polygon([(15, 92), (35, 92), (35, 98), (15, 98)]), + shapely.geometry.Polygon([(15, 82), (35, 82), (35, 88), (15, 88)]), + shapely.geometry.Polygon([(45, 62), (65, 62), (65, 68), (45, 68)]), + shapely.geometry.Polygon([(45, 52), (65, 52), (65, 58), (45, 58)]), + shapely.geometry.Point(25, 50), # Point building to test Point handling + shapely.geometry.Point(55, 80), # Another Point building + ] + buildings = gpd.GeoDataFrame( + { + "building_id": list(range(1, len(blg_polygons) + 1)), + "geometry": blg_polygons, + }, + crs="EPSG:3857", + ) + + # Define sample barrier geometries. + barrier_geoms = [ + shapely.geometry.LineString([(0, 0), (80, 0)]), + shapely.geometry.LineString([(80, 0), (80, 120)]), + shapely.geometry.LineString([(80, 120), (0, 120)]), + shapely.geometry.LineString([(0, 120), (0, 0)]), + shapely.geometry.LineString([(40, 0), (40, 110)]), + shapely.geometry.LineString([(10, 30), (40, 30)]), + shapely.geometry.LineString([(10, 90), (40, 90)]), + shapely.geometry.LineString([(40, 60), (70, 60)]), + ] + + inner_barriers = gpd.GeoDataFrame( + { + "name": [ + "Bottom Edge", + "Right Edge", + "Top Edge", + "Left Edge", + "Main Vertical", + "Left Cul-de-Sac (Bottom)", + "Left Cul-de-Sac (Top)", + "Right Cul-de-Sac (Middle)", + ], + "geometry": barrier_geoms, + }, + crs="EPSG:3857", + ) + + # Create enclosure from the barrier boundaries + enclosures = gpd.GeoDataFrame( + {"geometry": [box(0, 0, 80, 120)]}, crs="EPSG:3857" + ) + + tess = mm.enclosed_tessellation( + buildings, + enclosures, + simplify=False, + inner_barriers=inner_barriers, + cell_size=1, + neighbor_mode="neumann", + threshold=None, + n_jobs=1, + ) + + assert set(tess.geom_type.unique()) <= {"Polygon", "MultiPolygon"} + assert tess.enclosure_index.unique().tolist() == [0] + assert set(buildings.index).issubset(tess.index) + + point_barriers = gpd.GeoDataFrame( + {"geometry": [shapely.Point(20, 40), shapely.Point(80, 40)]}, + crs="EPSG:3857", + ) + tess_point = mm.enclosed_tessellation( + buildings, + enclosures, + simplify=False, + inner_barriers=point_barriers, + cell_size=1, + neighbor_mode="moore", + threshold=None, + n_jobs=1, + ) + assert set(tess_point.geom_type.unique()) <= {"Polygon", "MultiPolygon"} + assert set(buildings.index).issubset(tess_point.index) + + empty_barriers = gpd.GeoDataFrame(geometry=[], crs="EPSG:3857") + tess_empty = mm.enclosed_tessellation( + buildings, + enclosures, + simplify=False, + inner_barriers=empty_barriers, + cell_size=1, + neighbor_mode="moore", + threshold=None, + n_jobs=1, + ) + assert set(tess_empty.geom_type.unique()) <= {"Polygon", "MultiPolygon"} + assert set(buildings.index).issubset(tess_empty.index) + + def test_enclosed_tessellation_invalid_enclosure_geometry(self): + crs = self.df_buildings.crs + buildings = gpd.GeoDataFrame( + {"uID": [1, 2]}, + geometry=[box(10, 10, 30, 30), box(70, 10, 90, 30)], + crs=crs, + ) + invalid_enclosures = gpd.GeoDataFrame( + {"geometry": [LineString([(0, 20), (100, 20)])]}, crs=crs + ) + with pytest.raises(ValueError, match="Enclosure must be a Polygon"): + mm.enclosed_tessellation( + buildings, + invalid_enclosures, + simplify=False, + inner_barriers=self.df_streets.head(1), + cell_size=1, + threshold=None, + n_jobs=1, + ) + + def test_enclosed_tessellation_invalid_neighbor_mode(self): + crs = self.df_buildings.crs + buildings = gpd.GeoDataFrame( + {"uID": [1, 2]}, + geometry=[box(10, 10, 30, 30), box(70, 10, 90, 30)], + crs=crs, + ) + enclosures = gpd.GeoDataFrame( + { + "geometry": [ + shapely.MultiPolygon([box(0, 0, 40, 80), box(60, 0, 100, 80)]) + ] + }, + crs=crs, + ) + with pytest.raises(ValueError, match="Invalid neighbor_mode"): + mm.enclosed_tessellation( + buildings, + enclosures, + simplify=False, + inner_barriers=self.df_streets.head(1), + cell_size=1, + neighbor_mode="invalid", + threshold=None, + n_jobs=1, + ) + + def test_morphological_and_enclosed_tessellation_require_shapely(self, monkeypatch): + monkeypatch.setattr(elements, "SHPLY_GE_210", False) + with pytest.raises(ImportError): + mm.morphological_tessellation(self.df_buildings.head(1)) + with pytest.raises(ImportError): + mm.enclosed_tessellation( + self.df_buildings.head(1), + self.enclosures.geometry, + ) + def test_verify_tessellation(self): df = self.df_buildings b = df.total_bounds