![\"Cubes](\"https://raw.githubusercontent.com/EO-College/cubes-and-clouds/main/icons/cnc_3icons_process_circle.svg\"\n",){kind=icon}
"
]
},
{
"cell_type": "markdown",
"id": "43f2c4b3-88f7-4ebf-8248-a72728d14012",
"metadata": {},
"source": [
"# Reduce Operators with Pangeo\n",
"\n",
"## Reduce Operators\n",
"\n",
"When computing statistics over time or indices based on multiple bands, it is possible to use reduce operators.\n",
"\n",
"In Pangeo and `Xarray` we can use the `groupby` method, which applies a reducer to a data cube dimension by collapsing all the values along the specified dimension into an output value computed by the reducer."
]
},
{
"cell_type": "markdown",
"id": "5691278c-4d66-4ee2-8592-75431cbff3c7",
"metadata": {},
"source": [
"Reduce the temporal dimension to a single value, the mean for instance:"
]
},
{
"cell_type": "code",
"execution_count": 1,
"id": "8bb92f58-0ccc-425b-8f00-cba0cc2a60bb",
"metadata": {},
"outputs": [],
"source": [
"# STAC Catalogue Libraries\n",
"import pystac_client\n",
"import stackstac\n",
"import xarray as xr"
]
},
{
"cell_type": "code",
"execution_count": 2,
"id": "eacfff33-cb98-4665-90c7-072fcfdc6d88",
"metadata": {},
"outputs": [],
"source": [
"# West, South, East, North\n",
"spatial_extent = [11.259613, 46.461019, 11.406212, 46.522237]\n",
"temporal_extent = [\"2021-05-28T00:00:00Z\",\"2021-06-30T00:00:00Z\"]"
]
},
{
"cell_type": "code",
"execution_count": 3,
"id": "85c71d09-f0b6-4498-86df-75c3bdbd28aa",
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"/srv/conda/envs/notebook/lib/python3.11/site-packages/stackstac/prepare.py:408: UserWarning: The argument 'infer_datetime_format' is deprecated and will be removed in a future version. A strict version of it is now the default, see https://pandas.pydata.org/pdeps/0004-consistent-to-datetime-parsing.html. You can safely remove this argument.\n",
" times = pd.to_datetime(\n"
]
},
{
"data": {
"text/html": [
"<xarray.DataArray 'stackstac-b671372a6175b16d2e4eaf48f5606039' (time: 26,\n",
" band: 2,\n",
" y: 714, x: 1146)>\n",
"dask.array<fetch_raster_window, shape=(26, 2, 714, 1146), dtype=float64, chunksize=(1, 1, 714, 1024), chunktype=numpy.ndarray>\n",
"Coordinates: (12/53)\n",
" * time (time) datetime64[ns] 2021-05-28...\n",
" id (time) <U24 'S2B_32TPS_20210528_...\n",
" * band (band) <U3 'red' 'nir'\n",
" * x (x) float64 6.733e+05 ... 6.848e+05\n",
" * y (y) float64 5.155e+06 ... 5.148e+06\n",
" s2:unclassified_percentage (time) float64 0.2507 ... 2.393\n",
" ... ...\n",
" proj:transform object {0, 600000, 10, 5200020, ...\n",
" title (band) <U20 'Red (band 4) - 10m'...\n",
" common_name (band) <U3 'red' 'nir'\n",
" center_wavelength (band) float64 0.665 0.842\n",
" full_width_half_max (band) float64 0.038 0.145\n",
" epsg int64 32632\n",
"Attributes:\n",
" spec: RasterSpec(epsg=32632, bounds=(673310.0, 5147750.0, 684770.0...\n",
" crs: epsg:32632\n",
" transform: | 10.00, 0.00, 673310.00|\\n| 0.00,-10.00, 5154890.00|\\n| 0.0...\n",
" resolution: 10.0<xarray.DataArray 'stackstac-b671372a6175b16d2e4eaf48f5606039' (band: 2,\n",
" y: 714, x: 1146)>\n",
"dask.array<mean_agg-aggregate, shape=(2, 714, 1146), dtype=float64, chunksize=(1, 714, 1024), chunktype=numpy.ndarray>\n",
"Coordinates: (12/21)\n",
" * band (band) <U3 'red' 'nir'\n",
" * x (x) float64 6.733e+05 ... 6.848e+05\n",
" * y (y) float64 5.155e+06 ... 5.148e+06\n",
" s2:product_type <U7 'S2MSI2A'\n",
" s2:datatake_type <U8 'INS-NOBS'\n",
" proj:epsg int64 32632\n",
" ... ...\n",
" proj:transform object {0, 600000, 10, 5200020, ...\n",
" title (band) <U20 'Red (band 4) - 10m'...\n",
" common_name (band) <U3 'red' 'nir'\n",
" center_wavelength (band) float64 0.665 0.842\n",
" full_width_half_max (band) float64 0.038 0.145\n",
" epsg int64 32632<xarray.DataArray 'stackstac-b671372a6175b16d2e4eaf48f5606039' (time: 26,\n",
" y: 714, x: 1146)>\n",
"dask.array<mean_agg-aggregate, shape=(26, 714, 1146), dtype=float64, chunksize=(1, 714, 1024), chunktype=numpy.ndarray>\n",
"Coordinates: (12/48)\n",
" * time (time) datetime64[ns] 2021-05-28...\n",
" id (time) <U24 'S2B_32TPS_20210528_...\n",
" * x (x) float64 6.733e+05 ... 6.848e+05\n",
" * y (y) float64 5.155e+06 ... 5.148e+06\n",
" s2:unclassified_percentage (time) float64 0.2507 ... 2.393\n",
" s2:product_type <U7 'S2MSI2A'\n",
" ... ...\n",
" s2:degraded_msi_data_percentage (time) object 0.0112 0 ... 0.0127 0\n",
" earthsearch:payload_id (time) <U74 'roda-sentinel2/work...\n",
" gsd int64 10\n",
" proj:shape object {10980}\n",
" proj:transform object {0, 600000, 10, 5200020, ...\n",
" epsg int64 32632"
],
"text/plain": [
"<xarray.DataArray 'stackstac-b671372a6175b16d2e4eaf48f5606039' (time: 26,\n",
" y: 714, x: 1146)>\n",
"dask.array<truediv, shape=(26, 714, 1146), dtype=float64, chunksize=(1, 714, 1024), chunktype=numpy.ndarray>\n",
"Coordinates: (12/48)\n",
" * time (time) datetime64[ns] 2021-05-28...\n",
" id (time) <U24 'S2B_32TPS_20210528_...\n",
" * x (x) float64 6.733e+05 ... 6.848e+05\n",
" * y (y) float64 5.155e+06 ... 5.148e+06\n",
" s2:unclassified_percentage (time) float64 0.2507 ... 2.393\n",
" s2:product_type <U7 'S2MSI2A'\n",
" ... ...\n",
" s2:degraded_msi_data_percentage (time) object 0.0112 0 ... 0.0127 0\n",
" earthsearch:payload_id (time) <U74 'roda-sentinel2/work...\n",
" gsd int64 10\n",
" proj:shape object {10980}\n",
" proj:transform object {0, 600000, 10, 5200020, ...\n",
" epsg int64 32632"
],
"text/plain": [
"<xarray.DataArray 'stackstac-7710f08aba1196bc4f45cd73b40041b2' (time: 65,\n",
" band: 2)>\n",
"dask.array<nanmedian, shape=(65, 2), dtype=float64, chunksize=(1, 1), chunktype=numpy.ndarray>\n",
"Coordinates: (12/51)\n",
" * time (time) datetime64[ns] 2021-01-11...\n",
" id (time) <U24 'S2B_32TPS_20210111_...\n",
" * band (band) <U3 'red' 'nir'\n",
" s2:unclassified_percentage (time) object 0.001827 ... 5.000642\n",
" s2:product_type <U7 'S2MSI2A'\n",
" s2:sequence (time) <U1 '1' '3' '0' ... '0' '0'\n",
" ... ...\n",
" proj:transform object {0, 600000, 10, 5200020, ...\n",
" title (band) <U20 'Red (band 4) - 10m'...\n",
" common_name (band) <U3 'red' 'nir'\n",
" center_wavelength (band) float64 0.665 0.842\n",
" full_width_half_max (band) float64 0.038 0.145\n",
" epsg int64 32632<xarray.DataArray 'stackstac-7710f08aba1196bc4f45cd73b40041b2' (time: 65)>\n",
"dask.array<truediv, shape=(65,), dtype=float64, chunksize=(1,), chunktype=numpy.ndarray>\n",
"Coordinates: (12/46)\n",
" * time (time) datetime64[ns] 2021-01-11...\n",
" id (time) <U24 'S2B_32TPS_20210111_...\n",
" s2:unclassified_percentage (time) object 0.001827 ... 5.000642\n",
" s2:product_type <U7 'S2MSI2A'\n",
" s2:sequence (time) <U1 '1' '3' '0' ... '0' '0'\n",
" s2:generation_time (time) <U27 '2023-05-31T05:43:03...\n",
" ... ...\n",
" s2:degraded_msi_data_percentage (time) object 0.0223 0.025 ... 0 0\n",
" earthsearch:payload_id (time) <U74 'roda-sentinel2/work...\n",
" gsd int64 10\n",
" proj:shape object {10980}\n",
" proj:transform object {0, 600000, 10, 5200020, ...\n",
" epsg int64 32632"
],
"text/plain": [
"<xarray.DataArray 'stackstac-7710f08aba1196bc4f45cd73b40041b2' (time: 65)>\n",
"array([1. , 1. , 0.72450176, 1. , 0.69950249,\n",
" 1. , 0.55233069, 1. , 0.6202058 , 0.63488844,\n",
" 1. , 0.69789227, 1. , 0.69595376, 1. ,\n",
" 0.70536829, 1. , 0.71320346, 0.63970998, 1. ,\n",
" 0.6643002 , 1. , 0.66222645, 0.77107365, 1. ,\n",
" 0.70346552, 0.81553398, 1. , 0.80752841, 0.83191411,\n",
" 0.81575657, 0.82283105, 0.84372402, 1. , 0.83625219,\n",
" 0.81262327, 1. , 0.8295082 , 1. , 0.84407484,\n",
" 1. , 0.8255814 , 1. , 0.84149419, 1. ,\n",
" 0.82274742, 1. , 0.70364437, 1. , 0.82560706,\n",
" 0.77223851, 0.80659341, 1. , 0.80912863, 0.76804734,\n",
" 0.7275 , 0.77817963, 1. , 0.99874135, 0.86932447,\n",
" 0.82356729, 1. , 1. , 0.79028436, 0.78220286])\n",
"Coordinates: (12/46)\n",
" * time (time) datetime64[ns] 2021-01-11...\n",
" id (time) <U24 'S2B_32TPS_20210111_...\n",
" s2:unclassified_percentage (time) object 0.001827 ... 5.000642\n",
" s2:product_type <U7 'S2MSI2A'\n",
" s2:sequence (time) <U1 '1' '3' '0' ... '0' '0'\n",
" s2:generation_time (time) <U27 '2023-05-31T05:43:03...\n",
" ... ...\n",
" s2:degraded_msi_data_percentage (time) object 0.0223 0.025 ... 0 0\n",
" earthsearch:payload_id (time) <U74 'roda-sentinel2/work...\n",
" gsd int64 10\n",
" proj:shape object {10980}\n",
" proj:transform object {0, 600000, 10, 5200020, ...\n",
" epsg int64 32632"
],
"text/plain": [
"