Making a nice table

import xarray as xr
xr.set_options(display_style='html')
import intake
import cftime
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import numpy as np
import pandas as pd
import datetime
import seaborn as sns

def global_mean(ds):
    if 'lat' in ds:
        la = 'lat'
        lo = 'lon'
        weights = np.cos(np.deg2rad(ds.lat))
    elif 'latitude' in ds:
        weights = np.cos(np.deg2rad(ds.latitude))
        la = 'latitude'
        lo = 'longitude'
    else:
        weights = np.cos(np.deg2rad(ds.lat))
        
    return ds.weighted(weights).mean(la)

#cat_url = "https://storage.googleapis.com/cmip6/pangeo-cmip6.json"
#cat_url = "/mnt/craas1-ns9989k-geo4992/data/cmip6.json"
cat_url = '/mnt/craas2-ns9988k/data/catalogs/cmip6.json'
col = intake.open_esm_datastore(cat_url)
col
varl = ['tas','clt']

models = ['CESM2','CanESM5-1','EC-Earth3-Veg-LR','GISS-E2-2-H','E3SM-1-0']
cat = col.search(source_id = models, 
                 experiment_id=['historical'], 
                 table_id=['Amon',], 
                 variable_id=varl, 
                 member_id=['r1i1p1f1'],
                )
cat.df.head()

	variable_id	table_id	source_id	experiment_id	member_id	grid_label	time_range	activity_id	institution_id	version	path	dcpp_init_year
0	tas	Amon	CanESM5-1	historical	r1i1p1f1	gn	185001-201412	CMIP	CCCma	v20190429	/mnt/craas2-ns9988k-dl-ns9560k/ESGF/CMIP6/CMIP...	NaN
1	clt	Amon	E3SM-1-0	historical	r1i1p1f1	gr	185001-187412	CMIP	E3SM-Project	v20190913	/mnt/craas2-ns9988k-dl-ns9560k/ESGF/CMIP6/CMIP...	NaN
2	clt	Amon	E3SM-1-0	historical	r1i1p1f1	gr	187501-189912	CMIP	E3SM-Project	v20190913	/mnt/craas2-ns9988k-dl-ns9560k/ESGF/CMIP6/CMIP...	NaN
3	clt	Amon	E3SM-1-0	historical	r1i1p1f1	gr	190001-192412	CMIP	E3SM-Project	v20190913	/mnt/craas2-ns9988k-dl-ns9560k/ESGF/CMIP6/CMIP...	NaN
4	clt	Amon	E3SM-1-0	historical	r1i1p1f1	gr	192501-194912	CMIP	E3SM-Project	v20190913	/mnt/craas2-ns9988k-dl-ns9560k/ESGF/CMIP6/CMIP...	NaN

cat.df['source_id'].unique()

array(['CanESM5-1', 'E3SM-1-0', 'EC-Earth3-Veg-LR', 'GISS-E2-2-H',
       'CESM2'], dtype=object)

cat.esmcat.aggregation_control.groupby_attrs = ['activity_id','experiment_id', 'source_id','table_id','grid_label']
cat.esmcat.aggregation_control.groupby_attrs

['activity_id', 'experiment_id', 'source_id', 'table_id', 'grid_label']

Create dictionary from the list of datasets we found

• This step may take several minutes so be patient!

dset_dict = cat.to_dataset_dict(zarr_kwargs={'use_cftime':True})

--> The keys in the returned dictionary of datasets are constructed as follows:
	'activity_id.experiment_id.source_id.table_id.grid_label'

100.00% [5/5 00:36<00:00]

/opt/conda/envs/pangeo-notebook/lib/python3.11/site-packages/xarray/conventions.py:286: SerializationWarning: variable 'tas' has multiple fill values {1e+20, 1e+20} defined, decoding all values to NaN.
  var = coder.decode(var, name=name)
/opt/conda/envs/pangeo-notebook/lib/python3.11/site-packages/xarray/conventions.py:286: SerializationWarning: variable 'clt' has multiple fill values {1e+20, 1e+20} defined, decoding all values to NaN.
  var = coder.decode(var, name=name)

list(dset_dict.keys())

['CMIP.historical.CanESM5-1.Amon.gn',
 'CMIP.historical.GISS-E2-2-H.Amon.gn',
 'CMIP.historical.CESM2.Amon.gn',
 'CMIP.historical.E3SM-1-0.Amon.gr',
 'CMIP.historical.EC-Earth3-Veg-LR.Amon.gr']

simplify keys:

ds_dict = {}
for k in dset_dict.keys():
    mod = k.split('.')[2]
    ds_dict[mod] = dset_dict[k]

Compute global mean (replace by whatever you want to calculate)

ds_dict[mod]

<xarray.Dataset> Size: 811MB
Dimensions:    (time: 1980, bnds: 2, lat: 160, lon: 320, member_id: 1)
Coordinates:
  * time       (time) object 16kB 1850-01-16 12:00:00 ... 2014-12-16 12:00:00
    time_bnds  (time, bnds) object 32kB dask.array<chunksize=(1, 2), meta=np.ndarray>
  * lat        (lat) float64 1kB -89.14 -88.03 -86.91 ... 86.91 88.03 89.14
    lat_bnds   (lat, bnds) float64 3kB dask.array<chunksize=(160, 2), meta=np.ndarray>
  * lon        (lon) float64 3kB 0.0 1.125 2.25 3.375 ... 356.6 357.8 358.9
    lon_bnds   (lon, bnds) float64 5kB dask.array<chunksize=(320, 2), meta=np.ndarray>
  * member_id  (member_id) object 8B 'r1i1p1f1'
    height     float64 8B 2.0
Dimensions without coordinates: bnds
Data variables:
    clt        (member_id, time, lat, lon) float32 406MB dask.array<chunksize=(1, 1, 160, 320), meta=np.ndarray>
    tas        (member_id, time, lat, lon) float32 406MB dask.array<chunksize=(1, 1, 160, 320), meta=np.ndarray>
Attributes: (12/51)
    Conventions:                      CF-1.7 CMIP-6.2
    activity_id:                      CMIP
    branch_method:                    standard
    branch_time_in_child:             0.0
    contact:                          cmip6-data@ec-earth.org
    data_specs_version:               01.00.31
    ...                               ...
    intake_esm_attrs:grid_label:      gr
    intake_esm_attrs:activity_id:     CMIP
    intake_esm_attrs:institution_id:  EC-Earth-Consortium
    intake_esm_attrs:version:         v20200217
    intake_esm_attrs:_data_format_:   netcdf
    intake_esm_dataset_key:           CMIP.historical.EC-Earth3-Veg-LR.Amon.gr

xarray.Dataset

• Dimensions:
  • time: 1980
  • bnds: 2
  • lat: 160
  • lon: 320
  • member_id: 1
• Coordinates: (8)
  • time
    (time)
    object
    1850-01-16 12:00:00 ... 2014-12-...

    bounds :

    time_bnds

    axis :

    T

    long_name :

    time

    standard_name :

    time

    array([cftime.DatetimeProlepticGregorian(1850, 1, 16, 12, 0, 0, 0, has_year_zero=True),
           cftime.DatetimeProlepticGregorian(1850, 2, 15, 0, 0, 0, 0, has_year_zero=True),
           cftime.DatetimeProlepticGregorian(1850, 3, 16, 12, 0, 0, 0, has_year_zero=True),
           ...,
           cftime.DatetimeProlepticGregorian(2014, 10, 16, 12, 0, 0, 0, has_year_zero=True),
           cftime.DatetimeProlepticGregorian(2014, 11, 16, 0, 0, 0, 0, has_year_zero=True),
           cftime.DatetimeProlepticGregorian(2014, 12, 16, 12, 0, 0, 0, has_year_zero=True)],
          dtype=object)

  • time_bnds
    (time, bnds)
    object
    dask.array<chunksize=(1, 2), meta=np.ndarray>

    Array Chunk Bytes 30.94 kiB 16 B Shape (1980, 2) (1, 2) Dask graph 1980 chunks in 336 graph layers Data type object numpy.ndarray

  • lat
    (lat)
    float64
    -89.14 -88.03 ... 88.03 89.14

    bounds :

    lat_bnds

    units :

    degrees_north

    axis :

    Y

    long_name :

    Latitude

    standard_name :

    latitude

    array([-89.141519, -88.029429, -86.910771, -85.790629, -84.669924, -83.548947,
           -82.427818, -81.306595, -80.18531 , -79.063982, -77.942624, -76.821243,
           -75.699844, -74.578432, -73.457008, -72.335576, -71.214136, -70.09269 ,
           -68.97124 , -67.849784, -66.728326, -65.606864, -64.485399, -63.363932,
           -62.242462, -61.120991, -59.999518, -58.878044, -57.756569, -56.635092,
           -55.513614, -54.392135, -53.270655, -52.149175, -51.027694, -49.906212,
           -48.784729, -47.663246, -46.541763, -45.420279, -44.298794, -43.177309,
           -42.055824, -40.934338, -39.812852, -38.691366, -37.56988 , -36.448393,
           -35.326906, -34.205418, -33.083931, -31.962443, -30.840955, -29.719467,
           -28.597979, -27.476491, -26.355002, -25.233514, -24.112025, -22.990536,
           -21.869047, -20.747558, -19.626069, -18.50458 , -17.383091, -16.261601,
           -15.140112, -14.018622, -12.897133, -11.775643, -10.654153,  -9.532664,
            -8.411174,  -7.289684,  -6.168194,  -5.046704,  -3.925215,  -2.803725,
            -1.682235,  -0.560745,   0.560745,   1.682235,   2.803725,   3.925215,
             5.046704,   6.168194,   7.289684,   8.411174,   9.532664,  10.654153,
            11.775643,  12.897133,  14.018622,  15.140112,  16.261601,  17.383091,
            18.50458 ,  19.626069,  20.747558,  21.869047,  22.990536,  24.112025,
            25.233514,  26.355002,  27.476491,  28.597979,  29.719467,  30.840955,
            31.962443,  33.083931,  34.205418,  35.326906,  36.448393,  37.56988 ,
            38.691366,  39.812852,  40.934338,  42.055824,  43.177309,  44.298794,
            45.420279,  46.541763,  47.663246,  48.784729,  49.906212,  51.027694,
            52.149175,  53.270655,  54.392135,  55.513614,  56.635092,  57.756569,
            58.878044,  59.999518,  61.120991,  62.242462,  63.363932,  64.485399,
            65.606864,  66.728326,  67.849784,  68.97124 ,  70.09269 ,  71.214136,
            72.335576,  73.457008,  74.578432,  75.699844,  76.821243,  77.942624,
            79.063982,  80.18531 ,  81.306595,  82.427818,  83.548947,  84.669924,
            85.790629,  86.910771,  88.029429,  89.141519])

  • lat_bnds
    (lat, bnds)
    float64
    dask.array<chunksize=(160, 2), meta=np.ndarray>

    Array Chunk Bytes 2.50 kiB 2.50 kiB Shape (160, 2) (160, 2) Dask graph 1 chunks in 821 graph layers Data type float64 numpy.ndarray

  • lon
    (lon)
    float64
    0.0 1.125 2.25 ... 357.8 358.9

    bounds :

    lon_bnds

    units :

    degrees_east

    axis :

    X

    long_name :

    Longitude

    standard_name :

    longitude

    array([  0.   ,   1.125,   2.25 , ..., 356.625, 357.75 , 358.875])

  • lon_bnds
    (lon, bnds)
    float64
    dask.array<chunksize=(320, 2), meta=np.ndarray>

    Array Chunk Bytes 5.00 kiB 5.00 kiB Shape (320, 2) (320, 2) Dask graph 1 chunks in 821 graph layers Data type float64 numpy.ndarray

  • member_id
    (member_id)
    object
    'r1i1p1f1'

    array(['r1i1p1f1'], dtype=object)

  • height
    ()
    float64
    2.0

    units :

    m

    axis :

    Z

    positive :

    up

    long_name :

    height

    standard_name :

    height

    array(2.)

• Data variables: (2)
  • clt
    (member_id, time, lat, lon)
    float32
    dask.array<chunksize=(1, 1, 160, 320), meta=np.ndarray>

    standard_name :

    cloud_area_fraction

    long_name :

    Total Cloud Cover Percentage

    comment :

    Total cloud area fraction (reported as a percentage) for the whole atmospheric column, as seen from the surface or the top of the atmosphere. Includes both large-scale and convective cloud.

    units :

    %

    cell_methods :

    area: time: mean

    cell_measures :

    area: areacella

    Array Chunk Bytes 386.72 MiB 200.00 kiB Shape (1, 1980, 160, 320) (1, 1, 160, 320) Dask graph 1980 chunks in 496 graph layers Data type float32 numpy.ndarray

  • tas
    (member_id, time, lat, lon)
    float32
    dask.array<chunksize=(1, 1, 160, 320), meta=np.ndarray>

    standard_name :

    air_temperature

    long_name :

    Near-Surface Air Temperature

    comment :

    near-surface (usually, 2 meter) air temperature

    units :

    K

    cell_methods :

    area: time: mean

    cell_measures :

    area: areacella

    Array Chunk Bytes 386.72 MiB 200.00 kiB Shape (1, 1980, 160, 320) (1, 1, 160, 320) Dask graph 1980 chunks in 496 graph layers Data type float32 numpy.ndarray

• Indexes: (4)
  • time
    PandasIndex

    PandasIndex(CFTimeIndex([1850-01-16 12:00:00, 1850-02-15 00:00:00, 1850-03-16 12:00:00,
                 1850-04-16 00:00:00, 1850-05-16 12:00:00, 1850-06-16 00:00:00,
                 1850-07-16 12:00:00, 1850-08-16 12:00:00, 1850-09-16 00:00:00,
                 1850-10-16 12:00:00,
                 ...
                 2014-03-16 12:00:00, 2014-04-16 00:00:00, 2014-05-16 12:00:00,
                 2014-06-16 00:00:00, 2014-07-16 12:00:00, 2014-08-16 12:00:00,
                 2014-09-16 00:00:00, 2014-10-16 12:00:00, 2014-11-16 00:00:00,
                 2014-12-16 12:00:00],
                dtype='object',
                length=1980,
                calendar='proleptic_gregorian',
                freq=None))

  • lat
    PandasIndex

    PandasIndex(Index([-89.1415194264611, -88.0294288679515,  -86.910770814124,
           -85.7906288836365, -84.6699240844466, -83.5489469125421,
           -82.4278175240078, -81.3065945226689, -80.1853098724772,
           -79.0639824814086,
           ...
            79.0639824814086,  80.1853098724772,  81.3065945226689,
            82.4278175240078,  83.5489469125421,  84.6699240844466,
            85.7906288836365,   86.910770814124,  88.0294288679515,
            89.1415194264611],
          dtype='float64', name='lat', length=160))

  • lon
    PandasIndex

    PandasIndex(Index([    0.0,   1.125,    2.25,   3.375,     4.5,   5.625,    6.75,   7.875,
               9.0,  10.125,
           ...
            348.75, 349.875,   351.0, 352.125,  353.25, 354.375,   355.5, 356.625,
            357.75, 358.875],
          dtype='float64', name='lon', length=320))

  • member_id
    PandasIndex

    PandasIndex(Index(['r1i1p1f1'], dtype='object', name='member_id'))

• Attributes: (51)

  Conventions :

  CF-1.7 CMIP-6.2

  activity_id :

  CMIP

  branch_method :

  standard

  branch_time_in_child :

  0.0

  contact :

  cmip6-data@ec-earth.org

  data_specs_version :

  01.00.31

  experiment :

  all-forcing simulation of the recent past

  experiment_id :

  historical

  external_variables :

  areacella

  forcing_index :

  1

  frequency :

  mon

  further_info_url :

  https://furtherinfo.es-doc.org/CMIP6.EC-Earth-Consortium.EC-Earth3-Veg-LR.historical.none.r1i1p1f1

  grid :

  T159L62-ORCA1L75

  grid_label :

  gr

  initialization_index :

  1

  institution :

  AEMET, Spain; BSC, Spain; CNR-ISAC, Italy; DMI, Denmark; ENEA, Italy; FMI, Finland; Geomar, Germany; ICHEC, Ireland; ICTP, Italy; IDL, Portugal; IMAU, The Netherlands; IPMA, Portugal; KIT, Karlsruhe, Germany; KNMI, The Netherlands; Lund University, Sweden; Met Eireann, Ireland; NLeSC, The Netherlands; NTNU, Norway; Oxford University, UK; surfSARA, The Netherlands; SMHI, Sweden; Stockholm University, Sweden; Unite ASTR, Belgium; University College Dublin, Ireland; University of Bergen, Norway; University of Copenhagen, Denmark; University of Helsinki, Finland; University of Santiago de Compostela, Spain; Uppsala University, Sweden; Utrecht University, The Netherlands; Vrije Universiteit Amsterdam, the Netherlands; Wageningen University, The Netherlands. Mailing address: EC-Earth consortium, Rossby Center, Swedish Meteorological and Hydrological Institute/SMHI, SE-601 76 Norrkoping, Sweden

  institution_id :

  EC-Earth-Consortium

  mip_era :

  CMIP6

  nominal_resolution :

  250 km

  parent_activity_id :

  CMIP

  parent_experiment_id :

  piControl

  parent_mip_era :

  CMIP6

  parent_source_id :

  EC-Earth3-Veg-LR

  parent_time_units :

  days since 1850-01-01

  parent_variant_label :

  r1i1p1f1

  physics_index :

  1

  product :

  model-output

  realization_index :

  1

  realm :

  atmos

  source :

  EC-Earth3-Veg-LR (2019): aerosol: none atmos: IFS cy36r4 (TL159, linearly reduced Gaussian grid equivalent to 320 x 160 longitude/latitude; 62 levels; top level 5 hPa) atmosChem: none land: HTESSEL (land surface scheme built in IFS) and LPJ-GUESS v4 landIce: none ocean: NEMO3.6 (ORCA1 tripolar primarily 1 degree with meridional refinement down to 1/3 degree in the tropics; 362 x 292 longitude/latitude; 75 levels; top grid cell 0-1 m) ocnBgchem: none seaIce: LIM3

  source_id :

  EC-Earth3-Veg-LR

  source_type :

  AOGCM

  sub_experiment :

  none

  sub_experiment_id :

  none

  table_id :

  Amon

  table_info :

  Creation Date:(24 July 2019) MD5:70649eeb16bc90c431e35b583fac7375

  title :

  EC-Earth3-Veg-LR output prepared for CMIP6

  variant_label :

  r1i1p1f1

  license :

  CMIP6 model data produced by EC-Earth-Consortium is licensed under a Creative Commons Attribution-ShareAlike 4.0 International License (https://creativecommons.org/licenses). Consult https://pcmdi.llnl.gov/CMIP6/TermsOfUse for terms of use governing CMIP6 output, including citation requirements and proper acknowledgment. Further information about this data, including some limitations, can be found via the further_info_url (recorded as a global attribute in this file) and at http://www.ec-earth.org. The data producers and data providers make no warranty, either express or implied, including, but not limited to, warranties of merchantability and fitness for a particular purpose. All liabilities arising from the supply of the information (including any liability arising in negligence) are excluded to the fullest extent permitted by law.

  cmor_version :

  3.5.0

  branch_time_in_parent :

  164359.0

  intake_esm_attrs:table_id :

  Amon

  intake_esm_attrs:source_id :

  EC-Earth3-Veg-LR

  intake_esm_attrs:experiment_id :

  historical

  intake_esm_attrs:member_id :

  r1i1p1f1

  intake_esm_attrs:grid_label :

  gr

  intake_esm_attrs:activity_id :

  CMIP

  intake_esm_attrs:institution_id :

  EC-Earth-Consortium

  intake_esm_attrs:version :

  v20200217

  intake_esm_attrs:_data_format_ :

  netcdf

  intake_esm_dataset_key :

  CMIP.historical.EC-Earth3-Veg-LR.Amon.gr

ds_vals_mod_var = {}
for mod in ds_dict.keys():
    print(mod)
    ds_vals_mod_var[mod]=dict()
    _ds = ds_dict[mod]
    _ds_m = global_mean(_ds.sel(time=slice('1990','2015'))).mean()
    for v in varl:
        if v in _ds_m:
            val = float(_ds_m[v].values)
            ds_vals_mod_var[mod][v] = val
        else:
            ds_vals_mod_var[mod][v] = np.nan
            

CanESM5-1
GISS-E2-2-H
CESM2
E3SM-1-0
EC-Earth3-Veg-LR

/opt/conda/envs/pangeo-notebook/lib/python3.11/site-packages/dask/array/numpy_compat.py:53: RuntimeWarning: invalid value encountered in divide
  x = np.divide(x1, x2, out)

ds_vals_mod_var_meas = {}
for mod in ds_dict.keys():
    ds_vals_mod_var_meas[mod]=dict()
    _ds = ds_dict[mod].copy()
    meas = {}
    
    meas['mean'] = global_mean(_ds.sel(time=slice('1990','2015')).mean('lon').mean('time'))
    meas['median']= global_mean(_ds.sel(time=slice('1990','2015')).mean('lon').median('time'))
    meas['std'] = global_mean(_ds.sel(time=slice('1990','2015')).mean('lon').std('time'))
    for v in varl:
        ds_vals_mod_var_meas[mod][v] =dict()
        for m in meas:
            if v in meas[m]:
                try:
                    ds_vals_mod_var_meas[mod][v][m] = float(meas[m][v].squeeze().values)
                except:
                    ds_vals_mod_var_meas[mod][v][m] = np.nan
                    
            else:
                ds_vals_mod_var_meas[mod][v][m] = np.nan
            

Extract units

_ds = ds_dict[models[0]]
u_dic = {v:_ds[v].units for v in varl}
units = pd.DataFrame.from_dict(u_dic,orient='index',columns=['units'] )#index=['units'])

units

	units
tas	K
clt	%

import pandas as pd

df_nice =pd.DataFrame.from_dict(ds_vals_mod_var)
df_nice = pd.concat([units,df_nice], axis=1)
#pd.options.display.float_format = '{:20,.2f}'.format
pd.set_option('display.float_format','{:20,.2f}'.format)
display(df_nice)

	units	CanESM5-1	GISS-E2-2-H	CESM2	E3SM-1-0	EC-Earth3-Veg-LR
tas	K	287.65	NaN	287.94	287.47	287.54
clt	%	NaN	NaN	69.74	67.78	64.10

ls_meas = list()
for m in ds_vals_mod_var_meas:
    _df = pd.DataFrame.from_dict(ds_vals_mod_var_meas[m])
    _df['model']=m
    ls_meas.append(_df)

df_out = pd.concat(ls_meas).reset_index().sort_values('index').set_index(['index','model'])

pd.set_option('display.float_format','{:20,.2f}'.format)

df_out

		tas	clt
index	model
mean	CanESM5-1	287.65	NaN
	GISS-E2-2-H	NaN	NaN
	CESM2	287.94	69.74
	E3SM-1-0	287.47	67.78
	EC-Earth3-Veg-LR	287.54	64.10
median	CanESM5-1	287.47	NaN
	GISS-E2-2-H	NaN	NaN
	CESM2	287.92	70.04
	E3SM-1-0	287.35	68.22
	EC-Earth3-Veg-LR	287.29	64.48
std	CanESM5-1	3.70	NaN
	GISS-E2-2-H	NaN	NaN
	CESM2	3.47	5.73
	E3SM-1-0	3.61	5.48
	EC-Earth3-Veg-LR	3.65	6.35

df_out.T

index	mean					median					std
model	CanESM5-1	GISS-E2-2-H	CESM2	E3SM-1-0	EC-Earth3-Veg-LR	CanESM5-1	GISS-E2-2-H	CESM2	E3SM-1-0	EC-Earth3-Veg-LR	CanESM5-1	GISS-E2-2-H	CESM2	E3SM-1-0	EC-Earth3-Veg-LR
tas	287.65	NaN	287.94	287.47	287.54	287.47	NaN	287.92	287.35	287.29	3.70	NaN	3.47	3.61	3.65
clt	NaN	NaN	69.74	67.78	64.10	NaN	NaN	70.04	68.22	64.48	NaN	NaN	5.73	5.48	6.35