Skip to content

Commit

Permalink
Merge branch 'feature/zs_numpy' of https://github.com/GEOSYS/earthdai…
Browse files Browse the repository at this point in the history
…ly-python-client into feature/zs_numpy
  • Loading branch information
nicolasK committed Aug 2, 2024
2 parents 0d96dcf + a92dfd1 commit eb3ecca
Show file tree
Hide file tree
Showing 3 changed files with 18 additions and 13 deletions.
13 changes: 9 additions & 4 deletions earthdaily/earthdatastore/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -1062,14 +1062,19 @@ def datacube(
{
"clear_percent": (
"time",
np.multiply(np.divide(xr_datacube["clear_pixels"].data,xr_datacube.attrs["usable_pixels"]), 100
np.multiply(
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()

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:
Expand Down
11 changes: 5 additions & 6 deletions earthdaily/earthdatastore/cube_utils/_zonal.py
Original file line number Diff line number Diff line change
Expand Up @@ -48,27 +48,25 @@ def _rasterize(gdf, dataset, all_touched=False):
return feats, yx_pos



def zonal_stats_numpy(
dataset,
gdf,
operations=dict(mean=np.nanmean),
all_touched=False,
preload_time=False,
batch_time=1
batch_time=1,
):
tqdm_bar = tqdm.tqdm(total=len(dataset.data_vars) * int(np.ceil(dataset.time.size/batch_time)))
dataset = dataset.rio.clip_box(*gdf.to_crs(dataset.rio.crs).total_bounds)#.load()

feats, yx_pos = _rasterize(gdf, dataset, all_touched=all_touched)
vals = []



positions = [np.asarray(yx_pos[i + 1]) for i in np.arange(gdf.shape[0])]
for t_ in range(0,dataset.time.size,batch_time):
ts = np.arange(t_,np.min((t_+batch_time,dataset.time.size)))
dataset_time = dataset.isel(time=ts)
if batch_time>1:
if batch_time > 1:
dataset_time = dataset_time.load()
for t in ts:
vals.append({})
Expand Down Expand Up @@ -109,6 +107,7 @@ def zonal_stats_numpy(
tqdm_bar.close()
return da.transpose("feature", "time", "stats")


def zonal_stats(
dataset,
gdf,
Expand Down Expand Up @@ -151,7 +150,7 @@ def zonal_stats(
DESCRIPTION.
"""

if method == "geocube":
from geocube.api.core import make_geocube
from geocube.rasterize import rasterize_image
Expand Down
7 changes: 4 additions & 3 deletions earthdaily/earthdatastore/mask/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -143,14 +143,15 @@ def compute_clear_coverage(
"clear_percent": (
"time",
np.multiply(
self._obj['n_pixels_as_labels'] / self._obj.attrs["usable_pixels"],
self._obj["n_pixels_as_labels"]
/ self._obj.attrs["usable_pixels"],
100,
).astype(np.int8),
)
}
)
self._obj['clear_pixels'] = self._obj['clear_pixels'].load()
self._obj['clear_percent'] = self._obj['clear_percent'].load()
self._obj["clear_pixels"] = self._obj["clear_pixels"].load()
self._obj["clear_percent"] = self._obj["clear_percent"].load()

return self._obj

Expand Down

0 comments on commit eb3ecca

Please sign in to comment.