Add DriverVectorCube.get_bounding_box_area

openeo_driver/_version.py

@@ -1 +1 @@
-__version__ = "0.21.1a1"
+__version__ = "0.21.2a1"

openeo_driver/datacube.py

@@ -18,6 +18,7 @@
 from openeo_driver.datastructs import SarBackscatterArgs, ResolutionMergeArgs, StacAsset
 from openeo_driver.errors import FeatureUnsupportedException, InternalException
 from openeo_driver.util.ioformats import IOFORMATS
+from openeo_driver.util.utm import area_in_square_meters
 from openeo_driver.utils import EvalEnv
 
 log = logging.getLogger(__name__)
@@ -320,6 +321,7 @@ def to_legacy_save_result(self) -> Union["AggregatePolygonResult", "JSONResult"]
         )
 
     def get_bounding_box(self) -> Tuple[float, float, float, float]:
+        # TODO: cache bounding box?
         return tuple(self._geometries.total_bounds)
 
     def get_bounding_box_geometry(self) -> shapely.geometry.Polygon:
@@ -328,6 +330,12 @@ def get_bounding_box_geometry(self) -> shapely.geometry.Polygon:
     def get_bounding_box_geojson(self) -> dict:
         return shapely.geometry.mapping(self.get_bounding_box_geometry())
 
+    def get_bounding_box_area(self) -> float:
+        """Bounding box area in square meters"""
+        return area_in_square_meters(
+            self.get_bounding_box_geometry(), crs=self.get_crs()
+        )
+
     def get_geometries(self) -> Sequence[shapely.geometry.base.BaseGeometry]:
         return self._geometries.geometry
 

openeo_driver/util/utm.py

@@ -55,7 +55,7 @@ def auto_utm_crs_for_geometry(geometry: BaseGeometry, crs: str) -> str:
     return 'epsg:' + str(epsg)
 
 
-def geometry_to_crs(geometry, crs_from, crs_to):
+def geometry_to_crs(geometry: BaseGeometry, crs_from, crs_to):
     # Skip if CRS definitions are exactly the same
     if crs_from == crs_to:
         return geometry
@@ -71,7 +71,7 @@ def project(x, y, z=0):
     return shapely.ops.transform(project, geometry)
 
 
-def area_in_square_meters(geometry, crs: Union[str, pyproj.CRS]):
+def area_in_square_meters(geometry: BaseGeometry, crs: Union[str, pyproj.CRS]):
     if isinstance(crs, str):
         crs = "+init=" + crs  # TODO: this is deprecated
 

tests/data/geojson/FeatureCollection06.json

@@ -0,0 +1,37 @@
+{
+  "type": "FeatureCollection",
+  "features": [
+
+    {
+      "type": "Feature",
+      "properties": {},
+      "geometry": {
+        "type": "Polygon",
+        "coordinates": [
+          [
+            [
+              3.67466926574707,
+              51.037804967049205
+            ],
+            [
+              3.70737075805664,
+              51.037804967049205
+            ],
+            [
+              3.70737075805664,
+              51.05793176907366
+            ],
+            [
+              3.67466926574707,
+              51.05793176907366
+            ],
+            [
+              3.67466926574707,
+              51.037804967049205
+            ]
+          ]
+        ]
+      }
+    }
+  ]
+}

tests/test_vectorcube.py

@@ -1,11 +1,13 @@
 import geopandas as gpd
+import numpy.testing
 import pyproj
 import pytest
 import xarray
 from shapely.geometry import Polygon, MultiPolygon
 
 from openeo_driver.datacube import DriverVectorCube
 from openeo_driver.testing import DictSubSet
+
 from .data import get_path
 
 
@@ -254,3 +256,12 @@ def test_get_bounding_box(self, gdf):
                 ((1.0, 1.0), (1.0, 4.0), (5.0, 4.0), (5.0, 1.0), (1.0, 1.0)),
             ),
         }
+
+    def test_get_bounding_box_area(self):
+        path = str(get_path("geojson/FeatureCollection06.json"))
+        vc = DriverVectorCube(gpd.read_file(path))
+        area = vc.get_bounding_box_area()
+        # TODO: area of FeatureCollection06 (square with side of approx 2.3km, based on length of Watersportbaan Gent)
+        #       is roughly 2.3 km * 2.3 km = 5.29 km2,
+        #       but current implementation gives result that is quite a bit larger than that
+        numpy.testing.assert_allclose(area, 8e6, rtol=0.1)