Global waves — significant wave height¶

Demonstrates the 3-hourly wave-field pattern: pull a few days of the global wave reanalysis (MFWAM) significant wave height VHM0 over the North Atlantic, then show both a snapshot map and a point time-series through a winter storm.

Dataset cmems_mod_glo_wav_my_0.2deg_PT3H-i is a multi-year reanalysis (stable historical coverage) on a 0.2° grid with a 3-hourly step, so a fixed date works fine here.

Reads credentials from COPERNICUSMARINE_SERVICE_USERNAME / COPERNICUSMARINE_SERVICE_PASSWORD.

Setup¶

Consolidate the imports up front: xarray and pyramids' NetCDF for reading the downloaded field, plus the unified EarthLens entry point and the CMEMS Catalog. We also pick an output directory for the download.

In [ ]:

import os
from pathlib import Path

import xarray as xr
from pyramids.netcdf import NetCDF

from earthlens import EarthLens
from earthlens.cmems import Catalog

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

import os from pathlib import Path import xarray as xr from pyramids.netcdf import NetCDF from earthlens import EarthLens from earthlens.cmems import Catalog OUT_DIR = Path('data/cmems-wave') OUT_DIR.mkdir(parents=True, exist_ok=True)

Inspect the dataset in the catalog¶

Before downloading, look the dataset up in the CMEMS Catalog to confirm its domain, cadence, and the units of the VHM0 variable. We also fix the point (NE Atlantic, west of Brittany) used later for the time-series.

In [ ]:

DATASET_ID = 'cmems_mod_glo_wav_my_0.2deg_PT3H-i'
POINT_LAT, POINT_LON = 48.0, -16.0  # NE Atlantic, west of Brittany

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

DATASET_ID = 'cmems_mod_glo_wav_my_0.2deg_PT3H-i' POINT_LAT, POINT_LON = 48.0, -16.0 # NE Atlantic, west of Brittany ds_meta = Catalog().get_dataset(DATASET_ID) print(f'{DATASET_ID}: domain={ds_meta.domain}, cadence={ds_meta.cadence}') print('VHM0 units:', ds_meta.variables['VHM0'].units)

Download three days over the North Atlantic¶

3-hourly × 3 days = 24 time steps of a 2-D field over a ~25°×15° box. We build the request first — source, date window, cadence, dataset, variable, bounding box, output path, and the Copernicus Marine credentials.

In [ ]:

el = EarthLens(
    data_source='cmems',
    start='2014-02-08',
    end='2014-02-10',
    cadence='hourly',
    dataset=DATASET_ID,
    variables=['VHM0'],
    aoi=[-30.0, 40.0, -5.0, 55.0],
    path=str(OUT_DIR),
    service_username=os.environ.get('COPERNICUSMARINE_SERVICE_USERNAME'),
    service_password=os.environ.get('COPERNICUSMARINE_SERVICE_PASSWORD'),
)

el = EarthLens( data_source='cmems', start='2014-02-08', end='2014-02-10', cadence='hourly', dataset=DATASET_ID, variables=['VHM0'], aoi=[-30.0, 40.0, -5.0, 55.0], path=str(OUT_DIR), service_username=os.environ.get('COPERNICUSMARINE_SERVICE_USERNAME'), service_password=os.environ.get('COPERNICUSMARINE_SERVICE_PASSWORD'), )

With the request built, download() fetches the subset to the output directory and returns the list of written file paths.

In [ ]:

paths = el.download()
print(paths)

paths = el.download() print(paths)

Open the field and pull a point series¶

Read the downloaded NetCDF with pyramids' NetCDF, hand it to xarray via decode_cf, and report the cube dimensions.

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))

Significant wave height at the chosen point¶

Select the VHM0 field and pull the nearest grid point to (48N, 16W), then report the peak significant wave height as the storm passes.

In [ ]:

vhm0 = ds['VHM0']
point = vhm0.sel(latitude=POINT_LAT, longitude=POINT_LON, method='nearest')
peak = point.max()
print(f'peak Hs at ({POINT_LAT}N, {abs(POINT_LON)}W): {float(peak):.1f} m')

vhm0 = ds['VHM0'] point = vhm0.sel(latitude=POINT_LAT, longitude=POINT_LON, method='nearest') peak = point.max() print(f'peak Hs at ({POINT_LAT}N, {abs(POINT_LON)}W): {float(peak):.1f} m')

Snapshot map + point time-series¶

Left: the significant-wave-height field at the stormiest step. Right: the 3-hourly Hs series at the chosen point as the storm passes. First find the peak time step and slice out that snapshot.

In [ ]:

import matplotlib.pyplot as plt

peak_step = int(point.argmax('time'))
snapshot = vhm0.isel(time=peak_step)

import matplotlib.pyplot as plt peak_step = int(point.argmax('time')) snapshot = vhm0.isel(time=peak_step)

Draw the two panels¶

A pcolormesh map of the snapshot (with the point marked by a red star) alongside the 3-hourly Hs time-series at that point.

In [ ]:

fig, (axm, axt) = plt.subplots(1, 2, figsize=(12, 4))
mesh = axm.pcolormesh(
    snapshot['longitude'].values,
    snapshot['latitude'].values,
    snapshot.values,
    cmap='viridis',
    shading='auto',
)
fig.colorbar(mesh, ax=axm, label='Hs (m)')
axm.plot([POINT_LON], [POINT_LAT], 'r*', markersize=12)
axm.set_title(f'VHM0 snapshot {str(snapshot["time"].values)[:13]}')
axm.set_xlabel('lon')
axm.set_ylabel('lat')

axt.plot(point['time'].values, point.values, marker='o', ms=3)
axt.set_title('Significant wave height at 48N, 16W')
axt.set_ylabel('Hs (m)')
axt.set_xlabel('time')
axt.grid(alpha=0.3)
fig.autofmt_xdate()
fig.tight_layout()

fig, (axm, axt) = plt.subplots(1, 2, figsize=(12, 4)) mesh = axm.pcolormesh( snapshot['longitude'].values, snapshot['latitude'].values, snapshot.values, cmap='viridis', shading='auto', ) fig.colorbar(mesh, ax=axm, label='Hs (m)') axm.plot([POINT_LON], [POINT_LAT], 'r*', markersize=12) axm.set_title(f'VHM0 snapshot {str(snapshot["time"].values)[:13]}') axm.set_xlabel('lon') axm.set_ylabel('lat') axt.plot(point['time'].values, point.values, marker='o', ms=3) axt.set_title('Significant wave height at 48N, 16W') axt.set_ylabel('Hs (m)') axt.set_xlabel('time') axt.grid(alpha=0.3) fig.autofmt_xdate() fig.tight_layout()