v0.2.4

CHANGELOG.md

@@ -4,7 +4,17 @@ All notable changes to this project will be documented in this file.
 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.2.3] - 2024-23-07
+
+## [0.2.4] - 2024-08-13
+
+### Changed
+
+- `zonal_stats` now uses homemade `numpy` method. Using the accessor 
+`smart_load` is set by default to `True`. stats methods are given to 
+the  `reducers` argument.
+
+
+## [0.2.3] - 2024-07-23
 
 ### Added
 

earthdaily/__init__.py

@@ -7,7 +7,7 @@
 # to hide warnings from rioxarray or nano seconds conversion
 # warnings.filterwarnings("ignore")
 
-__version__ = "0.2.3"
+__version__ = "0.2.4"
 
 
 def EarthDataStore(

earthdaily/accessor/__init__.py

@@ -132,17 +132,36 @@ def sel_nearest_dates(
     def zonal_stats(
         self,
         geometry,
-        operations: list = ["mean"],
-        raise_missing_geometry: bool = False,
+        reducers: list = ["mean"],
+        label: str = None,
+        smart_load=True,
+        **kwargs,
     ):
+        """
+        Zonal stats from dtacube
+
+        Parameters
+        ----------
+        geometry : str,gpd.GeoDataFrame
+            A geometry (wkt, geopandas...)
+        stats : list, optional
+            The default is ["mean"].
+        raise_missing_geometry : bool, optional
+            DESCRIPTION. The default is False.
+        **kwargs : dict
+            Any kwargs for xvec.zonal_stats.
+
+        Returns
+        -------
+        xr.dataset
+            With new dimension "zonal_statistics" and "geometry".
+
+        """
         from ..earthdatastore.cube_utils import zonal_stats, GeometryManager
 
         geometry = GeometryManager(geometry).to_geopandas()
         return zonal_stats(
-            self._obj,
-            geometry,
-            operations=operations,
-            raise_missing_geometry=raise_missing_geometry,
+            self._obj, geometry, reducers=reducers, smart_load=smart_load, **kwargs
         )
 
     def lee_filter(self, window_size: int):

earthdaily/earthdatastore/__init__.py

@@ -60,6 +60,7 @@ def _parallel_search(func):
     def _search(*args, **kwargs):
         from joblib import Parallel, delayed
 
+        t0 = time.time()
         kwargs.setdefault("batch_days", "auto")
         batch_days = kwargs.get("batch_days", None)
         datetime = kwargs.get("datetime", None)
@@ -79,6 +80,7 @@ def _search(*args, **kwargs):
                     raise _no_item_msg
         if not need_parallel:
             items = func(*args, **kwargs)
+        logging.info(f"Search/load items : {np.round(time.time()-t0,3)}s.")
         return items
 
     return _search
@@ -795,6 +797,7 @@ def datacube(
         cross_calibration_collection: (None | str) = None,
         properties: (bool | str | list) = False,
         groupby_date: str = "mean",
+        cloud_search_kwargs={},
         **kwargs,
     ) -> xr.Dataset:
         """
@@ -995,7 +998,9 @@ def datacube(
                     "ag_cloud_mask": {"agriculture-cloud-mask": "ag_cloud_mask"},
                     "cloud_mask": {"cloud-mask": "cloud_mask"},
                 }
-                acm_items = self.find_cloud_mask_items(items, cloudmask=mask_with)
+                acm_items = self.find_cloud_mask_items(
+                    items, cloudmask=mask_with, **cloud_search_kwargs
+                )
                 acm_datacube = datacube(
                     acm_items,
                     intersects=intersects,
@@ -1055,21 +1060,26 @@ def datacube(
             n_pixels_as_labels = xy - null_pixels
 
             xr_datacube = xr_datacube.assign_coords(
-                {"clear_pixels": ("time", n_pixels_as_labels.load().values)}
+                {"clear_pixels": ("time", n_pixels_as_labels.data)}
             )
 
             xr_datacube = xr_datacube.assign_coords(
                 {
                     "clear_percent": (
                         "time",
                         np.multiply(
-                            xr_datacube["clear_pixels"].values
-                            / xr_datacube.attrs["usable_pixels"],
+                            np.divide(
+                                xr_datacube["clear_pixels"].data,
+                                xr_datacube.attrs["usable_pixels"],
+                            ),
                             100,
                         ).astype(np.int8),
                     )
                 }
             )
+
+            xr_datacube["clear_pixels"] = xr_datacube["clear_pixels"].load()
+            xr_datacube["clear_percent"] = xr_datacube["clear_percent"].load()
         if mask_with:
             xr_datacube = xr_datacube.drop(mask_with)
         if clear_cover:
@@ -1253,7 +1263,9 @@ def search(
             raise _no_item_msg
         return items_collection
 
-    def find_cloud_mask_items(self, items_collection, cloudmask="ag_cloud_mask"):
+    def find_cloud_mask_items(
+        self, items_collection, cloudmask="ag_cloud_mask", **kwargs
+    ):
         """
         Search the catalog for the ag_cloud_mask items matching the given items_collection.
         The ag_cloud_mask items are searched in the `ag_cloud_mask_collection_id` collection using the
@@ -1293,9 +1305,9 @@ def ag_cloud_mask_from_items(items):
         for items_start_idx in range(0, len(ids_), step):
             items = self.search(
                 collections=collections,
-                # intersects=self.intersects,
                 ids=ids_[items_start_idx : items_start_idx + step],
                 limit=step,
+                **kwargs,
             )
             items_list.extend(list(items))
         return ItemCollection(items_list)

earthdaily/earthdatastore/cube_utils/__init__.py

@@ -7,7 +7,7 @@
 from rasterio.enums import Resampling
 from shapely.geometry import box
 from .geometry_manager import GeometryManager
-from ._zonal import zonal_stats, zonal_stats_numpy
+from ._zonal import zonal_stats
 from .harmonizer import Harmonizer
 from .asset_mapper import AssetMapper
 import rioxarray