Skip to article frontmatterSkip to article content
Site not loading correctly?

This may be due to an incorrect BASE_URL configuration. See the MyST Documentation for reference.

Global Mean Surface Temperature Anomalies (GMSTA) from CMIP6 data


Overview

In this notebook we will compute the Global Mean Surface Temperature Anomalies (GMSTA) from CMIP6 data and compare it with observations. This notebook is heavily inspired by the GMST example in the CMIP6 cookbook and we thank the authors for their workflow.

  1. We will get the CMIP6 temperature data from the AWS open data program via the us-west-2 origin

  2. In order to do this, we will use an intake-ESM catalog (hosted on NCAR’s GDEX) that uses pelicanFS backed links instead of https or s3 links

  3. We will grab observational data hosted on NCAR’s GDEX, which is accessible via the NCAR origin

  4. Please refer to the first chapter of this cookbook to learn more about OSDF, pelican or pelicanFS

  5. This notebook demonstrates that you can seamlessly stream data from multiple OSDF origins in your workflow

Prerequisites

ConceptsImportanceNotes
Intro to Intake-ESMNecessaryUsed for searching CMIP6 data
Understanding of ZarrHelpfulFamiliarity with metadata structure
SeabornHelpfulUsed for plotting
PelicanFSNecessaryThe python package used to stream data in this notebook
OSDFHelpfulOSDF is used to stream data in this notebook
  • Time to learn: 20 mins

Imports

from matplotlib import pyplot as plt
import xarray as xr
import numpy as np
from dask.diagnostics import progress
from tqdm.autonotebook import tqdm
import intake
import fsspec
import seaborn as sns
import aiohttp
import dask
from dask.distributed import LocalCluster
import pelicanfs 
/tmp/ipykernel_4455/2849986847.py:5: TqdmWarning: IProgress not found. Please update jupyter and ipywidgets. See https://ipywidgets.readthedocs.io/en/stable/user_install.html
  from tqdm.autonotebook import tqdm

We will use an intake-ESM catalog hosted on NCAR’s Geoscience Data Exchange. This is nothing but the AWS cmip6 catalog modified to use OSDF

# Load catalog URL
gdex_url    =  'https://data.gdex.ucar.edu/'
cat_url     = gdex_url +  'd850001/catalogs/osdf/cmip6-aws/cmip6-osdf-zarr.json'
print(cat_url)
https://data.gdex.ucar.edu/d850001/catalogs/osdf/cmip6-aws/cmip6-osdf-zarr.json

Set up local dask cluster

Before we do any computation let us first set up a local cluster using dask

cluster = LocalCluster()          
client = cluster.get_client()
/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/distributed/node.py:188: UserWarning: Port 8787 is already in use.
Perhaps you already have a cluster running?
Hosting the HTTP server on port 46785 instead
  warnings.warn(
# Scale the cluster
n_workers = 3
cluster.scale(n_workers)
cluster
Loading...

Data Loading

Load CMIP6 data from AWS

col = intake.open_esm_datastore(cat_url)
col
Loading...
# there is currently a significant amount of data for these runs
expts = ['historical', 'ssp245', 'ssp370']

query = dict(
    experiment_id=expts,
    table_id='Amon',
    variable_id=['tas'],
    member_id = 'r1i1p1f1',
    #activity_id = 'CMIP',
)

col_subset = col.search(require_all_on=["source_id"], **query)
col_subset
Loading...
  • Let us inspect the zarr store paths to see if we are using the pelican protocol.

  • We see that zstore column has paths that start with ‘osdf:///’ instead of ‘https://’ which tells us that we are not using a simple ‘https’ GET request to fetch the data.

  • In order to know more about the pelican protocol, please refer to the first chapter of this cookbook.

col_subset.df
Loading...

Grab some Observational time series data for comparison with ensemble spread

  • The observational data we will use is the HadCRUT5 dataset from the UK Met Office

  • The data has been downloaded to NCAR’s Geoscience Data Exchange (GDEX) from https://www.metoffice.gov.uk/hadobs/hadcrut5/

  • We will use an OSDF to access this copy from the GDEX. Again the links will start with ‘osdf:///’

%%time
obs_url    = 'osdf:///ncar/gdex/d850001/HadCRUT.5.0.2.0.analysis.summary_series.global.monthly.nc'
#
obs_ds = xr.open_dataset(obs_url, engine='h5netcdf').tas_mean
obs_ds
CPU times: user 1.83 s, sys: 238 ms, total: 2.07 s
Wall time: 3.87 s
Loading...

Some helpful functions

def drop_all_bounds(ds):
    drop_vars = [vname for vname in ds.coords
                 if (('_bounds') in vname ) or ('_bnds') in vname]
    return ds.drop_vars(drop_vars)

def open_dset(df):
    assert len(df) == 1
    mapper = fsspec.get_mapper(df.zstore.values[0])
    #path = df.zstore.values[0][7:]+".zmetadata"
    ds = xr.open_zarr(mapper, consolidated=True)
    return drop_all_bounds(ds)

def open_delayed(df):
    return dask.delayed(open_dset)(df)

from collections import defaultdict
dsets = defaultdict(dict)

for group, df in col_subset.df.groupby(by=['source_id', 'experiment_id']):
    dsets[group[0]][group[1]] = open_delayed(df)
dsets_ = dask.compute(dict(dsets))[0]
#calculate global means
def get_lat_name(ds):
    for lat_name in ['lat', 'latitude']:
        if lat_name in ds.coords:
            return lat_name
    raise RuntimeError("Couldn't find a latitude coordinate")

def global_mean(ds):
    lat = ds[get_lat_name(ds)]
    weight = np.cos(np.deg2rad(lat))
    weight /= weight.mean()
    other_dims = set(ds.dims) - {'time'}
    return (ds * weight).mean(other_dims)
#calculate global means
def get_lat_name(ds):
    for lat_name in ['lat', 'latitude']:
        if lat_name in ds.coords:
            return lat_name
    raise RuntimeError("Couldn't find a latitude coordinate")

def global_mean(ds):
    lat = ds[get_lat_name(ds)]
    weight = np.cos(np.deg2rad(lat))
    weight /= weight.mean()
    other_dims = set(ds.dims) - {'time'}
    return (ds * weight).mean(other_dims)

GMST computation

expt_da = xr.DataArray(expts, dims='experiment_id', name='experiment_id',
                       coords={'experiment_id': expts})

dsets_aligned = {}

for k, v in tqdm(dsets_.items()):
    expt_dsets = v.values()
    if any([d is None for d in expt_dsets]):
        print(f"Missing experiment for {k}")
        continue

    for ds in expt_dsets:
        ds.coords['year'] = ds.time.dt.year

    # workaround for
    # https://github.com/pydata/xarray/issues/2237#issuecomment-620961663
    dsets_ann_mean = [v[expt].pipe(global_mean).swap_dims({'time': 'year'})
                             .drop_vars('time').coarsen(year=12).mean()
                      for expt in expts]

    # align everything with the 4xCO2 experiment
    dsets_aligned[k] = xr.concat(dsets_ann_mean, join='outer',dim=expt_da)
  0%|          | 0/27 [00:00<?, ?it/s]
  4%|▎         | 1/27 [00:02<00:56,  2.17s/it]
  7%|▋         | 2/27 [00:03<00:41,  1.64s/it]
 11%|█         | 3/27 [00:04<00:34,  1.44s/it]
 15%|█▍        | 4/27 [00:06<00:34,  1.51s/it]
 19%|█▊        | 5/27 [00:07<00:34,  1.56s/it]
 22%|██▏       | 6/27 [00:08<00:22,  1.08s/it]
 26%|██▌       | 7/27 [00:09<00:25,  1.26s/it]
 30%|██▉       | 8/27 [00:11<00:25,  1.36s/it]
 33%|███▎      | 9/27 [00:12<00:25,  1.44s/it]
 37%|███▋      | 10/27 [00:14<00:25,  1.52s/it]
 41%|████      | 11/27 [00:16<00:24,  1.55s/it]
 44%|████▍     | 12/27 [00:17<00:23,  1.57s/it]
 48%|████▊     | 13/27 [00:19<00:20,  1.47s/it]
 52%|█████▏    | 14/27 [00:20<00:19,  1.53s/it]
 56%|█████▌    | 15/27 [00:22<00:18,  1.56s/it]
 59%|█████▉    | 16/27 [00:23<00:17,  1.58s/it]
 63%|██████▎   | 17/27 [00:25<00:15,  1.59s/it]
 67%|██████▋   | 18/27 [00:27<00:14,  1.60s/it]
 70%|███████   | 19/27 [00:28<00:12,  1.61s/it]
 74%|███████▍  | 20/27 [00:30<00:10,  1.50s/it]
 78%|███████▊  | 21/27 [00:31<00:08,  1.41s/it]
 81%|████████▏ | 22/27 [00:33<00:07,  1.58s/it]
 85%|████████▌ | 23/27 [00:35<00:06,  1.64s/it]
 89%|████████▉ | 24/27 [00:36<00:04,  1.63s/it]
 93%|█████████▎| 25/27 [00:38<00:03,  1.62s/it]
 96%|█████████▋| 26/27 [00:39<00:01,  1.51s/it]
100%|██████████| 27/27 [00:41<00:00,  1.59s/it]
100%|██████████| 27/27 [00:41<00:00,  1.53s/it]

%%time
with progress.ProgressBar():
    dsets_aligned_ = dask.compute(dsets_aligned)[0]
2026-02-16 02:30:54,645 - distributed.protocol.pickle - ERROR - Failed to serialize 500, message='Internal Server Error', url='https://sdsc-cache.nationalresearchplatform.org:8443/aws-opendata/us-west-2/cmip6-pds/CMIP6/CMIP/CMCC/CMCC-ESM2/historical/r1i1p1f1/Amon/tas/gn/v20210114/tas/4.0.0'.
Traceback (most recent call last):
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/distributed/worker.py", line 2988, in _run_task_simple
    result = task(data)
             ^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/dask/_task_spec.py", line 768, in __call__
    return self.func(*new_argspec)
           ^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/dask/_task_spec.py", line 199, in _execute_subgraph
    res = execute_graph(final, keys=[outkey])
          ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/dask/_task_spec.py", line 1087, in execute_graph
    cache[key] = node(cache)
                 ^^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/dask/_task_spec.py", line 768, in __call__
    return self.func(*new_argspec)
           ^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/dask/_task_spec.py", line 199, in _execute_subgraph
    res = execute_graph(final, keys=[outkey])
          ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/dask/_task_spec.py", line 1087, in execute_graph
    cache[key] = node(cache)
                 ^^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/dask/_task_spec.py", line 768, in __call__
    return self.func(*new_argspec)
           ^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/dask/array/core.py", line 140, in getter
    c = np.asarray(c)
        ^^^^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/core/indexing.py", line 659, in __array__
    return np.asarray(self.get_duck_array(), dtype=dtype, copy=copy)
                      ^^^^^^^^^^^^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/core/indexing.py", line 664, in get_duck_array
    return self.array.get_duck_array()
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/core/indexing.py", line 897, in get_duck_array
    return self.array.get_duck_array()
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/coding/common.py", line 80, in get_duck_array
    return self.func(self.array.get_duck_array())
                     ^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/core/indexing.py", line 737, in get_duck_array
    array = self.array[self.key]
            ~~~~~~~~~~^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/backends/zarr.py", line 262, in __getitem__
    return indexing.explicit_indexing_adapter(
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/core/indexing.py", line 1129, in explicit_indexing_adapter
    result = raw_indexing_method(raw_key.tuple)
             ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/backends/zarr.py", line 225, in _getitem
    return self._array[key]
           ~~~~~~~~~~~^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/zarr/core.py", line 798, in __getitem__
    result = self.get_orthogonal_selection(pure_selection, fields=fields)
             ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/zarr/core.py", line 1080, in get_orthogonal_selection
    return self._get_selection(indexer=indexer, out=out, fields=fields)
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/zarr/core.py", line 1343, in _get_selection
    self._chunk_getitems(
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/zarr/core.py", line 2179, in _chunk_getitems
    cdatas = self.chunk_store.getitems(ckeys, contexts=contexts)
             ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/zarr/storage.py", line 1435, in getitems
    raise v
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/fsspec/asyn.py", line 244, in _run_coro
    return await asyncio.wait_for(coro, timeout=timeout), i
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/asyncio/tasks.py", line 520, in wait_for
    return await fut
           ^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/fsspec/implementations/http.py", line 247, in _cat_file
    self._raise_not_found_for_status(r, url)
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/fsspec/implementations/http.py", line 230, in _raise_not_found_for_status
    response.raise_for_status()
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/aiohttp/client_reqrep.py", line 636, in raise_for_status
    raise ClientResponseError(
aiohttp.client_exceptions.ClientResponseError: 500, message='Internal Server Error', url='https://sdsc-cache.nationalresearchplatform.org:8443/aws-opendata/us-west-2/cmip6-pds/CMIP6/CMIP/CMCC/CMCC-ESM2/historical/r1i1p1f1/Amon/tas/gn/v20210114/tas/4.0.0'

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/distributed/protocol/pickle.py", line 63, in dumps
    result = pickle.dumps(x, **dump_kwargs)
             ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
TypeError: can't pickle multidict._multidict.CIMultiDictProxy objects

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/distributed/protocol/pickle.py", line 68, in dumps
    pickler.dump(x)
TypeError: can't pickle multidict._multidict.CIMultiDictProxy objects

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/distributed/protocol/pickle.py", line 80, in dumps
    result = cloudpickle.dumps(x, **dump_kwargs)
             ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/cloudpickle/cloudpickle.py", line 1544, in dumps
    cp.dump(obj)
  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/cloudpickle/cloudpickle.py", line 1313, in dump
    return super().dump(obj)
           ^^^^^^^^^^^^^^^^^
TypeError: can't pickle multidict._multidict.CIMultiDictProxy objects
2026-02-16 02:30:54,941 - distributed.worker - ERROR - Compute Failed
Key:       ('mean_chunk-mean_agg-aggregate-360930ee1a35f4ea31044d90b91e4630', 4)
State:     executing
Task:  <Task ('mean_chunk-mean_agg-aggregate-360930ee1a35f4ea31044d90b91e4630', 4) _execute_subgraph(...)>
Exception: "ClientResponseError(RequestInfo(url=URL('https://sdsc-cache.nationalresearchplatform.org:8443/aws-opendata/us-west-2/cmip6-pds/CMIP6/CMIP/CMCC/CMCC-ESM2/historical/r1i1p1f1/Amon/tas/gn/v20210114/tas/4.0.0'), method='GET', headers=<CIMultiDictProxy('Host': 'sdsc-cache.nationalresearchplatform.org:8443', 'Accept': '*/*', 'Accept-Encoding': 'gzip, deflate, br, zstd', 'User-Agent': 'Python/3.12 aiohttp/3.13.3')>, real_url=URL('https://sdsc-cache.nationalresearchplatform.org:8443/aws-opendata/us-west-2/cmip6-pds/CMIP6/CMIP/CMCC/CMCC-ESM2/historical/r1i1p1f1/Amon/tas/gn/v20210114/tas/4.0.0')), (), status=500, message='Internal Server Error', headers=<CIMultiDictProxy('Connection': 'Close', 'Server': 'XrootD/v5.9.1', 'Access-Control-Allow-Origin': '*', 'Access-Control-Allow-Headers': 'Authorization,Want-Digest,Content-Type,User-Agent', 'Content-Length': '155')>)"
Traceback: '  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/dask/array/core.py", line 140, in getter\n    c = np.asarray(c)\n        ^^^^^^^^^^^^^\n  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/core/indexing.py", line 659, in __array__\n    return np.asarray(self.get_duck_array(), dtype=dtype, copy=copy)\n                      ^^^^^^^^^^^^^^^^^^^^^\n  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/core/indexing.py", line 664, in get_duck_array\n    return self.array.get_duck_array()\n           ^^^^^^^^^^^^^^^^^^^^^^^^^^^\n  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/core/indexing.py", line 897, in get_duck_array\n    return self.array.get_duck_array()\n           ^^^^^^^^^^^^^^^^^^^^^^^^^^^\n  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/coding/common.py", line 80, in get_duck_array\n    return self.func(self.array.get_duck_array())\n                     ^^^^^^^^^^^^^^^^^^^^^^^^^^^\n  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/core/indexing.py", line 737, in get_duck_array\n    array = self.array[self.key]\n            ~~~~~~~~~~^^^^^^^^^^\n  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/backends/zarr.py", line 262, in __getitem__\n    return indexing.explicit_indexing_adapter(\n           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^\n  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/core/indexing.py", line 1129, in explicit_indexing_adapter\n    result = raw_indexing_method(raw_key.tuple)\n             ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^\n  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/backends/zarr.py", line 225, in _getitem\n    return self._array[key]\n           ~~~~~~~~~~~^^^^^\n  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/zarr/core.py", line 798, in __getitem__\n    result = self.get_orthogonal_selection(pure_selection, fields=fields)\n             ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^\n  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/zarr/core.py", line 1080, in get_orthogonal_selection\n    return self._get_selection(indexer=indexer, out=out, fields=fields)\n           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^\n  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/zarr/core.py", line 1343, in _get_selection\n    self._chunk_getitems(\n  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/zarr/core.py", line 2179, in _chunk_getitems\n    cdatas = self.chunk_store.getitems(ckeys, contexts=contexts)\n             ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^\n  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/zarr/storage.py", line 1435, in getitems\n    raise v\n  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/fsspec/asyn.py", line 244, in _run_coro\n    return await asyncio.wait_for(coro, timeout=timeout), i\n           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^\n  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/asyncio/tasks.py", line 520, in wait_for\n    return await fut\n           ^^^^^^^^^\n  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/fsspec/implementations/http.py", line 247, in _cat_file\n    self._raise_not_found_for_status(r, url)\n  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/fsspec/implementations/http.py", line 230, in _raise_not_found_for_status\n    response.raise_for_status()\n  File "/home/runner/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/aiohttp/client_reqrep.py", line 636, in raise_for_status\n    raise ClientResponseError(\n'

CPU times: user 5.1 s, sys: 779 ms, total: 5.88 s
Wall time: 39.5 s
---------------------------------------------------------------------------
Exception                                 Traceback (most recent call last)
Cell In[14], line 1
----> 1 get_ipython().run_cell_magic('time', '', 'with progress.ProgressBar():\n    dsets_aligned_ = dask.compute(dsets_aligned)[0]\n')

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/IPython/core/interactiveshell.py:2572, in InteractiveShell.run_cell_magic(self, magic_name, line, cell)
   2570 with self.builtin_trap:
   2571     args = (magic_arg_s, cell)
-> 2572     result = fn(*args, **kwargs)
   2574 # The code below prevents the output from being displayed
   2575 # when using magics with decorator @output_can_be_silenced
   2576 # when the last Python token in the expression is a ';'.
   2577 if getattr(fn, magic.MAGIC_OUTPUT_CAN_BE_SILENCED, False):

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/IPython/core/magics/execution.py:1447, in ExecutionMagics.time(self, line, cell, local_ns)
   1445 if interrupt_occured:
   1446     if exit_on_interrupt and captured_exception:
-> 1447         raise captured_exception
   1448     return
   1449 return out

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/IPython/core/magics/execution.py:1411, in ExecutionMagics.time(self, line, cell, local_ns)
   1409 st = clock2()
   1410 try:
-> 1411     exec(code, glob, local_ns)
   1412     out = None
   1413     # multi-line %%time case

File <timed exec>:2

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/dask/base.py:685, in compute(traverse, optimize_graph, scheduler, get, *args, **kwargs)
    682     expr = expr.optimize()
    683     keys = list(flatten(expr.__dask_keys__()))
--> 685     results = schedule(expr, keys, **kwargs)
    687 return repack(results)

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/core/indexing.py:659, in __array__()
    655 def __array__(
    656     self, dtype: DTypeLike | None = None, /, *, copy: bool | None = None
    657 ) -> np.ndarray:
    658     if Version(np.__version__) >= Version("2.0.0"):
--> 659         return np.asarray(self.get_duck_array(), dtype=dtype, copy=copy)
    660     else:
    661         return np.asarray(self.get_duck_array(), dtype=dtype)

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/core/indexing.py:664, in get_duck_array()
    663 def get_duck_array(self):
--> 664     return self.array.get_duck_array()

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/core/indexing.py:897, in get_duck_array()
    896 def get_duck_array(self):
--> 897     return self.array.get_duck_array()

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/coding/common.py:80, in get_duck_array()
     79 def get_duck_array(self):
---> 80     return self.func(self.array.get_duck_array())

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/core/indexing.py:737, in get_duck_array()
    734 from xarray.backends.common import BackendArray
    736 if isinstance(self.array, BackendArray):
--> 737     array = self.array[self.key]
    738 else:
    739     array = apply_indexer(self.array, self.key)

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/backends/zarr.py:262, in __getitem__()
    260 elif isinstance(key, indexing.OuterIndexer):
    261     method = self._oindex
--> 262 return indexing.explicit_indexing_adapter(
    263     key, array.shape, indexing.IndexingSupport.VECTORIZED, method
    264 )

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/core/indexing.py:1129, in explicit_indexing_adapter()
   1107 """Support explicit indexing by delegating to a raw indexing method.
   1108 
   1109 Outer and/or vectorized indexers are supported by indexing a second time
   (...)   1126 Indexing result, in the form of a duck numpy-array.
   1127 """
   1128 raw_key, numpy_indices = decompose_indexer(key, shape, indexing_support)
-> 1129 result = raw_indexing_method(raw_key.tuple)
   1130 if numpy_indices.tuple:
   1131     # index the loaded duck array
   1132     indexable = as_indexable(result)

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/xarray/backends/zarr.py:225, in _getitem()
    224 def _getitem(self, key):
--> 225     return self._array[key]

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/zarr/core.py:798, in __getitem__()
    796     result = self.vindex[selection]
    797 elif is_pure_orthogonal_indexing(pure_selection, self.ndim):
--> 798     result = self.get_orthogonal_selection(pure_selection, fields=fields)
    799 else:
    800     result = self.get_basic_selection(pure_selection, fields=fields)

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/zarr/core.py:1080, in get_orthogonal_selection()
   1077 # setup indexer
   1078 indexer = OrthogonalIndexer(selection, self)
-> 1080 return self._get_selection(indexer=indexer, out=out, fields=fields)

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/zarr/core.py:1343, in _get_selection()
   1340 if math.prod(out_shape) > 0:
   1341     # allow storage to get multiple items at once
   1342     lchunk_coords, lchunk_selection, lout_selection = zip(*indexer)
-> 1343     self._chunk_getitems(
   1344         lchunk_coords,
   1345         lchunk_selection,
   1346         out,
   1347         lout_selection,
   1348         drop_axes=indexer.drop_axes,
   1349         fields=fields,
   1350     )
   1351 if out.shape:
   1352     return out

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/zarr/core.py:2179, in _chunk_getitems()
   2177     if not isinstance(self._meta_array, np.ndarray):
   2178         contexts = ConstantMap(ckeys, constant=Context(meta_array=self._meta_array))
-> 2179     cdatas = self.chunk_store.getitems(ckeys, contexts=contexts)
   2181 for ckey, chunk_select, out_select in zip(ckeys, lchunk_selection, lout_selection):
   2182     if ckey in cdatas:

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/zarr/storage.py:1435, in getitems()
   1432     continue
   1433 elif isinstance(v, Exception):
   1434     # Raise any other exception
-> 1435     raise v
   1436 else:
   1437     # The function calling this method may not recognize the transformed
   1438     # keys, so we send the values returned by self.map.getitems back into
   1439     # the original key space.
   1440     results[keys_transformed[k]] = v

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/fsspec/asyn.py:244, in _run_coro()
    242 async def _run_coro(coro, i):
    243     try:
--> 244         return await asyncio.wait_for(coro, timeout=timeout), i
    245     except Exception as e:
    246         if not return_exceptions:

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/asyncio/tasks.py:520, in wait_for()
    517         raise TimeoutError from exc
    519 async with timeouts.timeout(timeout):
--> 520     return await fut

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/fsspec/implementations/http.py:247, in _cat_file()
    245 async with session.get(self.encode_url(url), **kw) as r:
    246     out = await r.read()
--> 247     self._raise_not_found_for_status(r, url)
    248 return out

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/fsspec/implementations/http.py:230, in _raise_not_found_for_status()
    228 if response.status == 404:
    229     raise FileNotFoundError(url)
--> 230 response.raise_for_status()

File ~/micromamba/envs/osdf-cookbook/lib/python3.12/site-packages/aiohttp/client_reqrep.py:636, in raise_for_status()
    633 if not self._in_context:
    634     self.release()
--> 636 raise ClientResponseError(
    637     self.request_info,
    638     self.history,
    639     status=self.status,
    640     message=self.reason,
    641     headers=self.headers,
    642 )

Exception: ClientResponseError(RequestInfo(url=URL('https://sdsc-cache.nationalresearchplatform.org:8443/aws-opendata/us-west-2/cmip6-pds/CMIP6/CMIP/CMCC/CMCC-ESM2/historical/r1i1p1f1/Amon/tas/gn/v20210114/tas/4.0.0'), method='GET', headers=<CIMultiDictProxy('Host': 'sdsc-cache.nationalresearchplatform.org:8443', 'Accept': '*/*', 'Accept-Encoding': 'gzip, deflate, br, zstd', 'User-Agent': 'Python/3.12 aiohttp/3.13.3')>, real_url=URL('https://sdsc-cache.nationalresearchplatform.org:8443/aws-opendata/us-west-2/cmip6-pds/CMIP6/CMIP/CMCC/CMCC-ESM2/historical/r1i1p1f1/Amon/tas/gn/v20210114/tas/4.0.0')), (), status=500, message='Internal Server Error', headers=<CIMultiDictProxy('Connection': 'Close', 'Server': 'XrootD/v5.9.1', 'Access-Control-Allow-Origin': '*', 'Access-Control-Allow-Headers': 'Authorization,Want-Digest,Content-Type,User-Agent', 'Content-Length': '155')>)
source_ids = list(dsets_aligned_.keys())
source_da = xr.DataArray(source_ids, dims='source_id', name='source_id',
                         coords={'source_id': source_ids})

big_ds = xr.concat([ds.reset_coords(drop=True)
                    for ds in dsets_aligned_.values()],
                    dim=source_da)

big_ds
# Compute annual mean temperatures anomalies of observational data
obs_gmsta = obs_ds.resample(time='YS').mean(dim='time')
# obs_gmsta

Compute anomlaies and plot

  • We will compute the temperature anomalies w.r.t 1960-1990 baseline period

  • Convert xarray datasets to pandas dataframes

  • Use Seaborn to plot GMSTA

df_all = big_ds.to_dataframe().reset_index()
df_all.head()
# Define the baseline period
baseline_df = df_all[(df_all["year"] >= 1960) & (df_all["year"] <= 1990)]

# Compute the baseline mean
baseline_mean = baseline_df["tas"].mean()

# Compute anomalies
df_all["tas_anomaly"] = df_all["tas"] - baseline_mean
df_all
obs_df = obs_gmsta.to_dataframe(name='tas_anomaly').reset_index()
# Convert 'time' to 'year' (keeping only the year)
obs_df['year'] = obs_df['time'].dt.year

# Drop the original 'time' column since we extracted 'year'
obs_df = obs_df[['year', 'tas_anomaly']]
obs_df

Almost there! Let us now use seaborn to plot all the anomalies

g = sns.relplot(data=df_all, x="year", y="tas_anomaly",
                hue='experiment_id', kind="line", errorbar="sd", aspect=2, palette="Set2")  # Adjust the color palette)

# Get the current axis from the FacetGrid
ax = g.ax

# Overlay the observational data in red
sns.lineplot(data=obs_df, x="year", y="tas_anomaly",color="red", 
             linestyle="dashed", linewidth=2,label="Observations", ax=ax)

# Adjust the legend to include observations
ax.legend(title="Experiment ID + Observations")

# Show the plot
plt.show()

Summary

In this notebook, we used surface air temperature data from several CMIP6 models for the ‘historical’, ‘SSP245’ and ‘SSP370’ runs to compute Global Mean Surface Temperature Anomaly (GMSTA) relative to the 1960-1990 baseline period and compare it with anomalies computed from the HadCRUT monthly surface temperature dataset. We used a modified intake-ESM catalog and pelicanFS to ‘stream/download’ temperature data from two different OSDF origins. The CMIP6 model data was streamed from the AWS OpenData origin in the us-west-2 region and the observational data was streamed from NCAR’s OSDF origin.

Resources and references

  1. Original notebook in the Pangeo Gallery by Henri Drake and Ryan Abernathey

  2. CMIP6 cookbook by Ryan Abernathey, Henri Drake, Robert Ford and Max Grover

  3. Coupled Model Intercomparison Project 6 was accessed from https://registry.opendata.aws/cmip6 using a modified intake-ESM catalog hosted on NCAR’s GDEX

  4. We thank the UK Met Office Hadley Center for providing the observational data