ENH: add orientation and shared_walls functional versions (#553)

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 *

momepy/functional/_distribution.py

@@ -0,0 +1,91 @@
+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 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
+    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",
+    )
+
+
+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

momepy/functional/tests/test_distribution.py

@@ -0,0 +1,41 @@
+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)
+
+        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,
+            "sum": 5310.17039728293,
+            "min": 0,
+            "max": 106.20917523555639,
+        }
+        r = mm.shared_walls(self.df_buildings)
+        assert_result(r, expected, self.df_buildings)