From fc4b3161caba0365a10df26509521a2d84eec5b9 Mon Sep 17 00:00:00 2001 From: Martin Fleischmann Date: Thu, 29 Feb 2024 16:13:04 +0100 Subject: [PATCH 1/4] ENH: add orientation and shared_walls functional versions --- momepy/__init__.py | 1 + momepy/functional/_distribution.py | 101 +++++++++++++++++++ momepy/functional/tests/test_distribution.py | 32 ++++++ 3 files changed, 134 insertions(+) create mode 100644 momepy/functional/_distribution.py create mode 100644 momepy/functional/tests/test_distribution.py diff --git a/momepy/__init__.py b/momepy/__init__.py index f9a41de5..60f51162 100644 --- a/momepy/__init__.py +++ b/momepy/__init__.py @@ -8,6 +8,7 @@ from .diversity import * from .elements import * from .functional._dimension import * +from .functional._distribution import * from .functional._shape import * from .graph import * from .intensity import * diff --git a/momepy/functional/_distribution.py b/momepy/functional/_distribution.py new file mode 100644 index 00000000..05f4132d --- /dev/null +++ b/momepy/functional/_distribution.py @@ -0,0 +1,101 @@ +import geopandas as gpd +import numpy as np +import shapely +from geopandas import GeoDataFrame, GeoSeries +from packaging.version import Version +from pandas import Series + +__all__ = ["orientation", "shared_walls"] + +GPD_GE_013 = Version(gpd.__version__) >= Version("0.13.0") + + +def orientation(geometry: GeoDataFrame | GeoSeries) -> Series: + """Calculate the orientation of object. + + The 'orientation' is defined as the deviation of orientation of the bounding + rectangle from cardinal directions. As such it is within a range 0 - 45. The + orientation of LineStrings is represented by the orientation of the line connecting + the first and the last point of the segment. + + Adapted from :cite:`schirmer2015`. + + Parameters + ---------- + geometry : GeoDataFrame | GeoSeries + A GeoDataFrame or GeoSeries containing polygons to analyse. + + Returns + ------- + Series + """ + geom_types = geometry.geom_type + poly_mask = geom_types.str.contains("Polygon") + line_mask = geom_types.str.contains("Line") + + result = np.full(len(geometry), np.nan, dtype=float) + + if poly_mask.any(): + bboxes = shapely.minimum_rotated_rectangle(geometry.geometry.loc[poly_mask]) + coords = shapely.get_coordinates(bboxes) + pt0 = coords[::5] + pt1 = coords[1::5] + angle = np.arctan2(pt1[:, 0] - pt0[:, 0], pt1[:, 1] - pt0[:, 1]) + result[poly_mask] = np.degrees(angle) + + if line_mask.any(): + first = shapely.get_point(geometry.geometry, 0) + last = shapely.get_point(geometry.geometry, -1) + pt0 = shapely.get_coordinates(first) + pt1 = shapely.get_coordinates(last) + angle = np.arctan2(pt1[:, 0] - pt0[:, 0], pt1[:, 1] - pt0[:, 1]) + result[line_mask] = np.degrees(angle) + + return Series( + np.abs((result + 45) % 90 - 45), + index=geometry.index, + dtype=float, + name="orientation", + ) + # buildings + # 3.23 s ± 77 ms per loop (mean ± std. dev. of 7 runs, 1 loop each) + # 932 ms ± 12.5 ms per loop (mean ± std. dev. of 7 runs, 1 loop each) + # streets + # 374 ms ± 5.56 ms per loop (mean ± std. dev. of 7 runs, 1 loop each) + # 62 ms ± 681 µs per loop (mean ± std. dev. of 7 runs, 10 loops each) + + +def shared_walls(geometry: GeoDataFrame | GeoSeries) -> Series: + """Calculate the length of shared walls of adjacent elements (typically buildings). + + Note that data needs to be topologically correct. Overlapping polygons will lead to + incorrect results. + + Adapted from :cite:`hamaina2012a`. + + Parameters + ---------- + geometry : GeoDataFrame | GeoSeries + A GeoDataFrame or GeoSeries containing polygons to analyse. + + Returns + ------- + Series + """ + if GPD_GE_013: + inp, res = geometry.sindex.query(geometry.geometry, predicate="touches") + else: + inp, res = geometry.sindex.query_bulk(geometry.geometry, predicate="touches") + left = geometry.geometry.take(inp).reset_index(drop=True) + right = geometry.geometry.take(res).reset_index(drop=True) + intersections = left.intersection(right).length + walls = intersections.groupby(inp).sum() + walls.index = geometry.index.take(walls.index) + + results = Series(0.0, index=geometry.index, name="shared_walls") + results.loc[walls.index] = walls + return results + + # 8.68 s ± 28.7 ms per loop (mean ± std. dev. of 7 runs, 1 loop each) + # 10.3 s ± 17.2 ms per loop (mean ± std. dev. of 7 runs, 1 loop each) + # but more correct diff --git a/momepy/functional/tests/test_distribution.py b/momepy/functional/tests/test_distribution.py new file mode 100644 index 00000000..b6225c65 --- /dev/null +++ b/momepy/functional/tests/test_distribution.py @@ -0,0 +1,32 @@ +import geopandas as gpd + +import momepy as mm + +from .test_shape import assert_result + + +class TestDistribution: + def setup_method(self): + test_file_path = mm.datasets.get_path("bubenec") + self.df_buildings = gpd.read_file(test_file_path, layer="buildings") + self.df_streets = gpd.read_file(test_file_path, layer="streets") + + def test_orientation(self): + expected = { + "mean": 20.983859394267952, + "sum": 3021.6757527745854, + "min": 7.968673890244247, + "max": 42.329365250279125, + } + r = mm.orientation(self.df_buildings) + assert_result(r, expected, self.df_buildings) + + def test_shared_walls(self): + expected = { + "mean": 36.87618331446485, + "sum": 5310.17039728293, + "min": 0, + "max": 106.20917523555639, + } + r = mm.shared_walls(self.df_buildings) + assert_result(r, expected, self.df_buildings) From 59c5127996daa1804cc257eb272e48910a43274a Mon Sep 17 00:00:00 2001 From: Martin Fleischmann Date: Thu, 29 Feb 2024 16:15:48 +0100 Subject: [PATCH 2/4] cleanup --- momepy/functional/_distribution.py | 10 ---------- 1 file changed, 10 deletions(-) diff --git a/momepy/functional/_distribution.py b/momepy/functional/_distribution.py index 05f4132d..ea4bc14f 100644 --- a/momepy/functional/_distribution.py +++ b/momepy/functional/_distribution.py @@ -57,12 +57,6 @@ def orientation(geometry: GeoDataFrame | GeoSeries) -> Series: dtype=float, name="orientation", ) - # buildings - # 3.23 s ± 77 ms per loop (mean ± std. dev. of 7 runs, 1 loop each) - # 932 ms ± 12.5 ms per loop (mean ± std. dev. of 7 runs, 1 loop each) - # streets - # 374 ms ± 5.56 ms per loop (mean ± std. dev. of 7 runs, 1 loop each) - # 62 ms ± 681 µs per loop (mean ± std. dev. of 7 runs, 10 loops each) def shared_walls(geometry: GeoDataFrame | GeoSeries) -> Series: @@ -95,7 +89,3 @@ def shared_walls(geometry: GeoDataFrame | GeoSeries) -> Series: results = Series(0.0, index=geometry.index, name="shared_walls") results.loc[walls.index] = walls return results - - # 8.68 s ± 28.7 ms per loop (mean ± std. dev. of 7 runs, 1 loop each) - # 10.3 s ± 17.2 ms per loop (mean ± std. dev. of 7 runs, 1 loop each) - # but more correct From a456d18ee4a335dabf5fafb9b0b0ccd0fde0a87d Mon Sep 17 00:00:00 2001 From: Martin Fleischmann Date: Thu, 29 Feb 2024 16:36:41 +0100 Subject: [PATCH 3/4] test coverage --- momepy/functional/tests/test_distribution.py | 9 +++++++++ 1 file changed, 9 insertions(+) diff --git a/momepy/functional/tests/test_distribution.py b/momepy/functional/tests/test_distribution.py index b6225c65..50657fa6 100644 --- a/momepy/functional/tests/test_distribution.py +++ b/momepy/functional/tests/test_distribution.py @@ -21,6 +21,15 @@ def test_orientation(self): r = mm.orientation(self.df_buildings) assert_result(r, expected, self.df_buildings) + expected = { + "mean": 21.176405050561755, + "sum": 741.1741767696615, + "min": 0.834911325974133, + "max": 44.83357900046826, + } + r = mm.orientation(self.df_streets) + assert_result(r, expected, self.df_streets) + def test_shared_walls(self): expected = { "mean": 36.87618331446485, From 36fc85a7968c1bbbab9bcb248905faa928f9b1b6 Mon Sep 17 00:00:00 2001 From: Martin Fleischmann Date: Thu, 29 Feb 2024 17:44:47 +0100 Subject: [PATCH 4/4] Update _distribution.py Co-authored-by: James Gaboardi --- momepy/functional/_distribution.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/momepy/functional/_distribution.py b/momepy/functional/_distribution.py index ea4bc14f..55732de5 100644 --- a/momepy/functional/_distribution.py +++ b/momepy/functional/_distribution.py @@ -11,7 +11,7 @@ def orientation(geometry: GeoDataFrame | GeoSeries) -> Series: - """Calculate the orientation of object. + """Calculate the orientation of objects. The 'orientation' is defined as the deviation of orientation of the bounding rectangle from cardinal directions. As such it is within a range 0 - 45. The