Mediterranean Sea — surface temperature map¶

Demonstrates the regional-reanalysis map pattern: pull one summer day of the Mediterranean physics reanalysis potential temperature (thetao), clipped to the surface layer, and render the basin-wide SST field.

Dataset cmems_mod_med_phy-temp_my_4.2km_P1D-m is a multi-year reanalysis on a ~4 km grid (stable historical coverage, so a fixed date works). It is 4-D, so we clip minimum_depth / maximum_depth to the top few metres to fetch just the surface field.

Reads credentials from COPERNICUSMARINE_SERVICE_USERNAME / COPERNICUSMARINE_SERVICE_PASSWORD.

Setup¶

Imports up front: pyramids provides NetCDF (reading the downloaded cube), earthlens provides the unified EarthLens entry point and the CMEMS Catalog, and xarray / numpy / matplotlib handle the reduction and plot.

In [ ]:

import os
from pathlib import Path

import matplotlib.pyplot as plt
import numpy as np
import xarray as xr
from pyramids.netcdf import NetCDF

from earthlens import EarthLens
from earthlens.cmems import Catalog

import os from pathlib import Path import matplotlib.pyplot as plt import numpy as np import xarray as xr from pyramids.netcdf import NetCDF from earthlens import EarthLens from earthlens.cmems import Catalog

Request parameters¶

Pin the dataset id, the mid-summer day (strong basin-wide SST gradient), and an output directory for the downloaded NetCDF.

In [ ]:

OUT_DIR = Path('data/cmems-medsst')
OUT_DIR.mkdir(parents=True, exist_ok=True)

DATASET_ID = 'cmems_mod_med_phy-temp_my_4.2km_P1D-m'
DAY = '2020-08-15'  # mid-summer, strong basin-wide SST gradient

OUT_DIR = Path('data/cmems-medsst') OUT_DIR.mkdir(parents=True, exist_ok=True) DATASET_ID = 'cmems_mod_med_phy-temp_my_4.2km_P1D-m' DAY = '2020-08-15' # mid-summer, strong basin-wide SST gradient

Inspect the dataset metadata¶

Look up the dataset in the CMEMS Catalog to confirm its domain, cadence, and the units of the thetao variable before downloading.

In [ ]:

ds_meta = Catalog().get_dataset(DATASET_ID)
print(f'{DATASET_ID}: domain={ds_meta.domain}, cadence={ds_meta.cadence}')
print('thetao units:', ds_meta.variables['thetao'].units)

ds_meta = Catalog().get_dataset(DATASET_ID) print(f'{DATASET_ID}: domain={ds_meta.domain}, cadence={ds_meta.cadence}') print('thetao units:', ds_meta.variables['thetao'].units)

Download one day, surface layer, whole basin¶

Build the EarthLens request first — dataset, day window, basin bounding box, and the minimum_depth / maximum_depth clip that keeps only the top model level.

In [ ]:

el = EarthLens(
    data_source='cmems',
    start=DAY,
    end=DAY,
    cadence='daily',
    dataset=DATASET_ID,
    variables=['thetao'],
    aoi=[-6.0, 30.0, 36.5, 46.0],
    path=str(OUT_DIR),
    minimum_depth=0.0,
    maximum_depth=2.0,
    service_username=os.environ.get('COPERNICUSMARINE_SERVICE_USERNAME'),
    service_password=os.environ.get('COPERNICUSMARINE_SERVICE_PASSWORD'),
)

el = EarthLens( data_source='cmems', start=DAY, end=DAY, cadence='daily', dataset=DATASET_ID, variables=['thetao'], aoi=[-6.0, 30.0, 36.5, 46.0], path=str(OUT_DIR), minimum_depth=0.0, maximum_depth=2.0, service_username=os.environ.get('COPERNICUSMARINE_SERVICE_USERNAME'), service_password=os.environ.get('COPERNICUSMARINE_SERVICE_PASSWORD'), )

Run the download. download() returns the list of written file paths; this request writes a single NetCDF cube.

In [ ]:

paths = el.download()
print(paths)

paths = el.download() print(paths)

Open and reduce to a 2-D surface field¶

Read the cube with pyramids' NetCDF, hand it to xarray (decoding CF metadata), and close the file handle.

In [ ]:

nc = NetCDF.read_file(str(paths[0]), read_only=True)
ds = xr.decode_cf(nc.to_xarray())
nc.close()
print('dims:', dict(ds.sizes))

nc = NetCDF.read_file(str(paths[0]), read_only=True) ds = xr.decode_cf(nc.to_xarray()) nc.close() print('dims:', dict(ds.sizes))

Select the top model level on the single day to collapse the 4-D cube to a 2-D (latitude, longitude) field, then summarise the valid (non-NaN) surface temperatures.

In [ ]:

# Top model level, single day -> 2-D (latitude, longitude).
sst = ds['thetao'].isel(time=0, depth=0)
valid = sst.values[~np.isnan(sst.values)]
print(
    f'surface thetao min/mean/max: {valid.min():.1f} / {valid.mean():.1f} / {valid.max():.1f} degrees_C'
)

# Top model level, single day -> 2-D (latitude, longitude). sst = ds['thetao'].isel(time=0, depth=0) valid = sst.values[~np.isnan(sst.values)] print( f'surface thetao min/mean/max: {valid.min():.1f} / {valid.mean():.1f} / {valid.max():.1f} degrees_C' )

Map the surface temperature¶

Warm eastern-basin / Levantine waters versus the cooler Gulf of Lions and Aegean — the classic August Mediterranean SST signature.

In [ ]:

fig, ax = plt.subplots(figsize=(10, 4))
mesh = ax.pcolormesh(
    sst['longitude'].values,
    sst['latitude'].values,
    sst.values,
    cmap='inferno',
    shading='auto',
)
fig.colorbar(mesh, ax=ax, label='SST (degrees_C)')
ax.set_title(f'Mediterranean surface temperature ({DAY})')
ax.set_xlabel('lon')
ax.set_ylabel('lat')
ax.set_aspect('equal')
fig.tight_layout()

fig, ax = plt.subplots(figsize=(10, 4)) mesh = ax.pcolormesh( sst['longitude'].values, sst['latitude'].values, sst.values, cmap='inferno', shading='auto', ) fig.colorbar(mesh, ax=ax, label='SST (degrees_C)') ax.set_title(f'Mediterranean surface temperature ({DAY})') ax.set_xlabel('lon') ax.set_ylabel('lat') ax.set_aspect('equal') fig.tight_layout()