From 826415904cad7ea4a228cb9d4a4bcbbee32da949 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Fri, 10 Nov 2023 18:10:05 +0000 Subject: [PATCH] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- datacube.py | 44 +++++++++++++++++++++----------------------- 1 file changed, 21 insertions(+), 23 deletions(-) diff --git a/datacube.py b/datacube.py index 6e6031d2..c25840a2 100644 --- a/datacube.py +++ b/datacube.py @@ -22,20 +22,14 @@ import random from datetime import datetime, timedelta + import geopandas as gpd -import pystac -import pystac_client +import numpy as np import planetary_computer as pc +import pystac_client import stackstac import xarray as xr -import numpy as np -from shapely.geometry import box, shape, Point -import leafmap -from rasterio.errors import RasterioIOError -from rasterio.warp import transform_bounds -from getpass import getpass -import os -import json +from shapely.geometry import Point, box STAC_API = "https://planetarycomputer.microsoft.com/api/stac/v1" S2_BANDS = ["B02", "B03", "B04", "B05", "B06", "B07", "B08", "B8A", "B11", "B12", "SCL"] @@ -75,8 +69,8 @@ def get_week(year, month, day): date = datetime(year, month, day) start_of_week = date - timedelta(days=date.weekday()) end_of_week = start_of_week + timedelta(days=6) - start_date_str = start_of_week.strftime('%Y-%m-%d') - end_date_str = end_of_week.strftime('%Y-%m-%d') + start_date_str = start_of_week.strftime("%Y-%m-%d") + end_date_str = end_of_week.strftime("%Y-%m-%d") return f"{start_date_str}/{end_date_str}" @@ -114,7 +108,9 @@ def search_sentinel2(week, aoi): CENTROID = aoi.centroid BBOX = aoi.bounds - geom_CENTROID = Point(CENTROID.x, CENTROID.y) # Create point geom object from centroid + geom_CENTROID = Point( + CENTROID.x, CENTROID.y + ) # Create point geom object from centroid catalog = pystac_client.Client.open(STAC_API, modifier=pc.sign_inplace) @@ -243,9 +239,7 @@ def search_dem(BBOX, catalog, week, s2_items): Returns: - tuple: A tuple containing DEM items and DEM GeoDataFrame. """ - search = catalog.search( - collections=["cop-dem-glo-30"], bbox=BBOX - ) + search = catalog.search(collections=["cop-dem-glo-30"], bbox=BBOX) dem_items = search.item_collection() print(f"Found {len(dem_items)} items") @@ -283,7 +277,9 @@ def make_dataarrays(s2_items, s1_items, dem_items, BBOX, resolution, epsg): ) da_sen1: xr.DataArray = stackstac.stack( - items=s1_items[1:], # To only accept the same orbit state and date. Need better way to do this. + items=s1_items[ + 1: + ], # To only accept the same orbit state and date. Need better way to do this. assets=["vh", "vv"], # SAR polarizations epsg=epsg, bounds_latlon=BBOX, # W, S, E, N @@ -320,7 +316,7 @@ def make_dataarrays(s2_items, s1_items, dem_items, BBOX, resolution, epsg): fill_value=np.nan, ) - _, index = np.unique(da_dem['time'], return_index=True) # Remove redundant time + _, index = np.unique(da_dem["time"], return_index=True) # Remove redundant time da_dem = da_dem.isel(time=index) return da_sen2, da_sen1, da_dem @@ -338,7 +334,7 @@ def merge_datarrays(da_sen2, da_sen1, da_dem): Returns: - xr.DataArray: Merged xarray DataArray. """ - da_merge = xr.merge([da_sen2, da_sen1, da_dem], compat='override') + da_merge = xr.merge([da_sen2, da_sen1, da_dem], compat="override") print(da_merge) return da_merge @@ -367,16 +363,18 @@ def process(year1, year2, aoi, resolution, epsg): s1_items, s1_gdf = search_sentinel1(BBOX, catalog, week) # s1_items, s1_gdf = search_sentinel1_calc_max_area(BBOX, catalog, week) # WIP - dem_items, dem_gdf = search_dem(BBOX, catalog, week, s2_items) - da_sen2, da_sen1, da_dem = make_dataarrays(s2_items, s1_items, dem_items, BBOX, resolution, epsg) + da_sen2, da_sen1, da_dem = make_dataarrays( + s2_items, s1_items, dem_items, BBOX, resolution, epsg + ) da_merge = merge_datarrays(da_sen2, da_sen1, da_dem) return da_merge + # EXAMPLE -california_tile = gpd.read_file("ca.geojson") +california_tile = gpd.read_file("ca.geojson") sample = california_tile.sample(1) aoi = sample.iloc[0].geometry -process(2017, 2023, aoi, 10, 26910) # UTM Zone 10N and spatial resolution of 10 metres \ No newline at end of file +process(2017, 2023, aoi, 10, 26910) # UTM Zone 10N and spatial resolution of 10 metres