
import xarray as xr
import numpy as np
import hvplot.xarray
import holoviews as hv
import os, datetime, warnings
from cartopy import crs as ccrs
from cartopy import feature as cfeature
print('Last updated on {}'.format(datetime.datetime.now().ctime()))

Last updated on Mon Nov  8 17:15:03 2021


datadir = '/home/ivan/Data/Postproc/1DPOP_CESM2'
case = 'cesm2.B1850.f09_g17.1dpop.005'
proj = ccrs.Robinson(central_longitude=210)

vname = 'TEMP'
ds = xr.open_dataset(os.path.join(datadir, '{}.{}.annual.nc'.format(case,vname)), decode_times=False)
ds_surf = ds.isel(time=-1, z_t=0, z_w=0)

sst = ds_surf.TEMP.where((ds.TLAT>=-60) & (ds.TLAT<=60))
map1 = sst.hvplot.quadmesh(x='TLONG', y='TLAT', geo=True, coastline=True, projection=proj, project=True,
                           cmap='kbc', rasterize=True, global_extent=True, title='Near global', colorbar=False)

sst = ds_surf.TEMP.where(ds.REGION_MASK==6)
sst = sst.where(ds.TLAT>0)
map2 = sst.hvplot.quadmesh(x='TLONG', y='TLAT', geo=True, coastline=True, projection=proj, project=True,
                           cmap='kbc', rasterize=True, global_extent=True, title='North Atlantic', colorbar=False)

sst = ds_surf.TEMP.where(ds.REGION_MASK==2)
sst = sst.where((ds.TLAT>=-5) & (ds.TLAT<=5))
map3 = sst.hvplot.quadmesh(x='TLONG', y='TLAT', geo=True, coastline=True, projection=proj, project=True,
                           cmap='kbc', rasterize=True, global_extent=True, title='Tropical Pacific', colorbar=False)

sst = ds_surf.TEMP.where(ds.REGION_MASK==8)
map4 = sst.hvplot.quadmesh(x='TLONG', y='TLAT', geo=True, coastline=True, projection=proj, project=True,
                           cmap='kbc', rasterize=True, global_extent=True, title='Labrador Sea', colorbar=False)

sst = ds_surf.TEMP.where(ds.REGION_MASK==10)
map5 = sst.hvplot.quadmesh(x='TLONG', y='TLAT', geo=True, coastline=True, projection=proj, project=True,
                           cmap='kbc', rasterize=True, global_extent=True, title='Arctic', colorbar=False)

sst = ds_surf.TEMP.where(ds.REGION_MASK==1)
map6  = sst.hvplot.quadmesh(x='TLONG', y='TLAT', geo=True, coastline=True, projection=proj, project=True,
                           cmap='kbc', rasterize=True, global_extent=True, title='Southern Ocean', colorbar=False)

(map1 + map2 + map3 + map4 + map5 + map6).cols(2)

Map showing the regions for which diagnostics were extracted¶