Access Booms

import xarray as xr

data_href1 = 'https://s3.waw4-1.cloudferro.com/EarthCODE/OSCAssets/booms/seascape-t1-owt.zarr/'
data_href2 = 'https://s3.waw4-1.cloudferro.com/EarthCODE/OSCAssets/booms/seascape-t2-biophysical.zarr/'
data_href3 = 'https://s3.waw4-1.cloudferro.com/EarthCODE/OSCAssets/booms/seascape-t3-4d-atlantic.zarr/'
data_href4 = 'https://s3.waw4-1.cloudferro.com/EarthCODE/OSCAssets/booms/seascape-t3-4d-southern.zarr/'

ds = xr.open_zarr(data_href1)
ds

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.colors as mcolors
import cartopy.crs as ccrs
import cartopy.feature as cfeature


def plot_med_water_class(da, label="smoothed mode dominant water class"):

    # 14 discrete OWT classes
    classes = np.arange(1, 15)

    # A categorical colour map
    cmap = plt.get_cmap("tab20", 14)

    # Boundaries put each integer in the centre of a colour bin
    boundaries = np.arange(0.5, 15.5, 1)
    norm = mcolors.BoundaryNorm(
        boundaries,
        ncolors=cmap.N
    )

    fig = plt.figure(figsize=(13, 7))

    ax = plt.axes(
        projection=ccrs.PlateCarree()
    )

    im = ax.pcolormesh(
        da.lon,
        da.lat,
        da,
        cmap=cmap,
        norm=norm,
        shading="nearest",
        transform=ccrs.PlateCarree()
    )

    ax.set_extent(
        [lon_min, lon_max, lat_min, lat_max],
        crs=ccrs.PlateCarree()
    )

    ax.add_feature(
        cfeature.LAND,
        facecolor="lightgrey",
        zorder=2
    )

    ax.add_feature(
        cfeature.COASTLINE,
        linewidth=0.6,
        zorder=3
    )

    ax.add_feature(
        cfeature.BORDERS,
        linewidth=0.3,
        zorder=3
    )

    cbar = plt.colorbar(
        im,
        ax=ax,
        ticks=classes,
        boundaries=boundaries,
        pad=0.03,
        shrink=0.85
    )

    cbar.set_label(label)

    date = da.time.values if "time" in da.coords else ds.time.isel(time=time_index).values

    ax.set_title(
        f"Bio-optical Seascape - Mediterranean\n"
        f"{np.datetime_as_string(date, unit='D')}"
    )

    plt.tight_layout()

    return fig, ax

ds = xr.open_zarr(data_href1)

# subset data
da = ds['smoothed_mode_dominant_water_class'].isel(time=0)

# Mediterranean bounding box
lon_min, lon_max = -6, 37
lat_min, lat_max = 30, 46.5


# Crop before plotting
da = da.where(
    (da.lon >= lon_min) &
    (da.lon <= lon_max) &
    (da.lat >= lat_min) &
    (da.lat <= lat_max),
    drop=True
)

plot_med_water_class(da)

ds = xr.open_zarr(data_href2)

# subset data
da = ds['smoothed_dominant_water_class'].isel(time=0)

# Mediterranean bounding box
lon_min, lon_max = -6, 37
lat_min, lat_max = 30, 46.5


# Crop before plotting
da = da.where(
    (da.lon >= lon_min) &
    (da.lon <= lon_max) &
    (da.lat >= lat_min) &
    (da.lat <= lat_max),
    drop=True
)

plot_med_water_class(da, label="smoothed dominant water class")

ds = xr.open_zarr(data_href3)

# subset data
da = ds['mode_seascape_4D'].isel(time=0, depth=0).squeeze()
da = da.rename({'longitude': 'lon', 'latitude': 'lat'})
da_plot = da.transpose("lat", "lon")

plot_med_water_class(da_plot, label="mode seascape 4D depth 0")

ds = xr.open_zarr(data_href4)

da = ds.mode_seascape_4D.isel(time=0, depth=0).transpose("latitude", "longitude")

ax = plt.axes(projection=ccrs.SouthPolarStereo())
ax.set_extent([-180, 180, -90, -40], ccrs.PlateCarree())
im = ax.pcolormesh(da.longitude, da.latitude, da, transform=ccrs.PlateCarree(), shading="auto")
ax.coastlines()
plt.colorbar(im, ax=ax)
plt.show()