Add visibility_frequency wrapper (#1145)

Author: brendancol

xrspatial/tests/test_visibility.py

@@ -245,3 +245,32 @@ def test_preserves_coords_and_dims(self):
         np.testing.assert_array_equal(result.coords['y'].values,
                                       raster.coords['y'].values)
         assert result.attrs.get('crs') == 'EPSG:4326'
+
+
+from xrspatial.visibility import visibility_frequency
+
+
+class TestVisibilityFrequency:
+    def test_flat_terrain_all_ones(self):
+        data = np.zeros((10, 10), dtype=float)
+        raster = _make_raster(data)
+        observers = [
+            {'x': 2.0, 'y': 2.0, 'observer_elev': 10},
+            {'x': 7.0, 'y': 7.0, 'observer_elev': 10},
+        ]
+        result = visibility_frequency(raster, observers)
+        assert result.dtype == np.float64
+        np.testing.assert_allclose(result.values, 1.0)
+
+    def test_equals_cumulative_divided_by_n(self):
+        data = np.random.RandomState(7).rand(15, 15).astype(float) * 100
+        raster = _make_raster(data)
+        observers = [
+            {'x': 3.0, 'y': 3.0, 'observer_elev': 50},
+            {'x': 10.0, 'y': 10.0, 'observer_elev': 50},
+            {'x': 7.0, 'y': 2.0, 'observer_elev': 50},
+        ]
+        freq = visibility_frequency(raster, observers)
+        cum = cumulative_viewshed(raster, observers)
+        expected = cum.values.astype(np.float64) / 3.0
+        np.testing.assert_allclose(freq.values, expected)

xrspatial/visibility.py

@@ -234,3 +234,23 @@ def cumulative_viewshed(
 
     return xarray.DataArray(count, coords=raster.coords,
                             dims=raster.dims, attrs=raster.attrs)
+
+
+def visibility_frequency(
+    raster: xarray.DataArray,
+    observers: list,
+    target_elev: float = 0,
+    max_distance: float = None,
+) -> xarray.DataArray:
+    """Fraction of observers that can see each cell.
+
+    Parameters are the same as :func:`cumulative_viewshed`.
+
+    Returns
+    -------
+    xarray.DataArray
+        Float64 raster with values in [0, 1].
+    """
+    cum = cumulative_viewshed(raster, observers, target_elev, max_distance)
+    freq = cum.astype(np.float64) / len(observers)
+    return freq