v0.0.10

v0.0.10

Author: nkarasiak

CHANGELOG.md

@@ -5,8 +5,20 @@ The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/),
 and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
 
 
-## [0.0.10] - Unreleased
+## [0.0.10] - 2024-03-05
 
+### Added
+
+-  `ed.whittaker` adapted from pywapor github.
+- `ed.zonal_stats` using new `geocube` zonal_stats engine.
+
+### Changed
+
+- `zonal_stats` now uses `geocube` as default processing method.
+
+### Removed
+
+- `ed.plot_index` deprecated, use `ed.plot_band` instead.
 
 ## [0.0.9] - 2024-02-29
 
@@ -19,7 +31,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 ### Added
 
-- better management of `col_wrap` in  `ed` xarray accessor.
+- better management of `col_wrap` in `ed` xarray accessor.
 
 ### Fixed
 

earthdaily/accessor/__init__.py

@@ -334,3 +334,9 @@ def whittaker(
             max_value=max_value,
             max_iter=max_iter,
         )
+
+    def zonal_stats(self, geometry, operations: list = ["mean"]):
+        from ..earthdatastore.cube_utils import zonal_stats, GeometryManager
+
+        geometry = GeometryManager(geometry).to_geopandas()
+        return zonal_stats(self._obj, geometry, operations=operations)

earthdaily/earthdatastore/cube_utils/_zonal.py

@@ -107,9 +107,24 @@ def zonal_stats(
     gdf,
     operations=["mean"],
     all_touched=False,
-    method="optimized",
+    method="geocube",
     verbose=False,
 ):
+    if method == "geocube":
+        from geocube.api.core import make_geocube
+
+        gdf["feature"] = list(gdf.index)
+        out_grid = make_geocube(
+            gdf,
+            measurements=["feature"],
+            like=dataset,  # ensure the data are on the same grid
+        )
+        cube = dataset.groupby(out_grid.feature)
+        zonal_stats = xr.concat(
+            [getattr(cube, operation)() for operation in operations], dim="stats"
+        )
+        zonal_stats["stats"] = operations
+        return zonal_stats
     tqdm_bar = tqdm.tqdm(total=gdf.shape[0])
 
     if dataset.rio.crs != gdf.crs:

examples/field_evolution.py

@@ -55,11 +55,7 @@
 # Compute zonal stats for the pivot
 # ----------------------------------------------------
 
-zonal_stats = earthdatastore.cube_utils.zonal_stats(
-    pivot_cube, pivot, operations=["mean", "max", "min"],
-    method="standard"
-)
-zonal_stats = zonal_stats.load()
+zonal_stats = pivot_cube.ed.zonal_stats(pivot, ['mean','max','min'])
 
 zonal_stats.isel(feature=0).to_array(dim="band").plot.line(
     x="time", col="band", hue="stats"

setup.py

@@ -37,7 +37,8 @@
         "dask",
         "spyndex",
         "dask-image",
-        "numba"
+        "numba",
+        "geocube"
     ],
     include_package_data=True,
     package_data={"":['*.geojson','*.json']},