Diagnostic plots: vertical sections

Created by Ivan Lima on Thu Jul 20 2017 14:36:38 -0400

%matplotlib notebook
import xarray as xr
import numpy as np
import cartopy.crs as ccrs
import cartopy.feature as cfeature
import matplotlib.pyplot as plt
import os
from matplotlib import ticker, colors
from mpl_utils import spd, center_cmap
from cesm_utils import *
from cases import *
from datetime import datetime
print('Last updated on %s'%datetime.now().ctime())

Last updated on Tue Apr 30 09:26:48 2019

Case description

case = 'tr3he.GIAF.06'
print(case,':',case_info[case])

tr3he.GIAF.06 : 186-year run with interannualy varying forcing (CORE 1948-2009). Physics is
spun-up for 126 years (2 IAF cycles: 1824-1949). Model is spun-up with pre-bomb tritium fluxes and
historic tritium deposition starts in the third IAF cycle (1950-2009). Time & latitude dependent
tritium concentration in precipitation (Cp) is estimated from GNIP data (Michaela's blended fits).
Tritium concentration in water vapor (Cv) is set to 0.7 Cp. Relative humidity has a max value of
0.95. Added tritium river fluxes. Time & latitude dependent tritium concentration in rivers is
estimated from GNIP data (Michaela's blended fits). Tritium concentration in precipitation over
ocean and land are smoothed over a 1 degree latitudinal grid. Run includes vapor flux correction.

Read data

datadir = '/home/ivan/Data/Postproc/Tritium-3H'
ds_sec_j12  = xr.open_dataset(os.path.join(datadir,'%s_sec_j12.nc'%case),decode_times=False)
ds_sec_j225 = xr.open_dataset(os.path.join(datadir,'%s_sec_j225.nc'%case),decode_times=False)
ds_sec_j12 = fixtime(ds_sec_j12,1824)
ds_sec_j225 = fixtime(ds_sec_j225,1824)

Monthly means from last year of simulation

def secfig(title=''):
    fig, ax = plt.subplots(figsize=(9,5))
    ax.set_title('%s (%d annual mean)'%(title,date['time.year']))
    #ax.set_axis_bgcolor('#b0b0b0')
    ax.set(xlabel='latitude',ylabel='depth (m)')
    return fig, ax

Compute excess $^3\!$He

ds_sec_j12['EXCESS_HELIUM3'] = ds_sec_j12.HELIUM3 - ds_sec_j12.INERT_HELIUM3
ds_sec_j225['EXCESS_HELIUM3'] = ds_sec_j225.HELIUM3 - ds_sec_j225.INERT_HELIUM3

Position of sections

fig, ax = plt.subplots(subplot_kw={'projection':ccrs.Robinson(central_longitude=210)})#,figsize=(9.5,7.125))
fig.subplots_adjust(top=0.95)
ax.set_global()
ax.add_feature(cfeature.LAND,facecolor='#b0b0b0')
l1, = ax.plot(ds_sec_j12.TLONG, ds_sec_j12.TLAT,'C0--',transform=ccrs.Geodetic())
l2, = ax.plot(ds_sec_j225.TLONG, ds_sec_j225.TLAT,'C1--',transform=ccrs.Geodetic())

Monthly means from last year of simulation

North-south section along 26.5$^{\circ}$W

cmap = plt.cm.plasma

varinfo = {
    'HELIUM3' : {
        'label' : r'pmol m$^{-3}$',
        'title' : r'$^3\!$He',
    },
    'INERT_HELIUM3' : {
        'label' : r'pmol m$^{-3}$',
        'title' : r'inert $^3\!$He',
    },
    'EXCESS_HELIUM3' : {
        'label' : r'pmol m$^{-3}$',
        'title' : r'excess $^3\!$He',
    },
    'HELIUM4' : {
        'label' : r'$\mu$mol m$^{-3}$',
        'title' : r'$^4\!$He',
    },
    'NEON' : {
        'label' : r'$\mu$mol m$^{-3}$',
        'title' : r'Ne',
    },    
    'ARGON' : {
        'label' : r'mmol m$^{-3}$',
        'title' : r'Ar',
    },    
    'TRITIUM' : {
        'label' : r'pmol m$^{-3}$',
        'title' : r'$^3\!$H',
    },
    'TRITIUM_DECAY' : {
        'label' : r'pmol m$^{-3}$ y$^{-1}$',
        'title' : r'$^3\!$H decay',
        'cfac'  : spd * 365, # pmol/m^3/sec -> pmol/m^3/year
    },
}

for vname in sorted(varinfo):
    data, date = ds_sec_j12[vname][-12:,...].mean(dim='time'), ds_sec_j12.time[-1:,...]
    data = np.ma.masked_where(np.isnan(data.values),data.values) * varinfo[vname].get('cfac',1)
    data[data<0] = 0
    fig, ax = secfig(varinfo[vname]['title'])
    pm = ax.pcolormesh(ds_sec_j12.TLAT,-ds_sec_j12.z_t/100,data,cmap=cmap)
    cs = ax.contour(ds_sec_j12.TLAT,-ds_sec_j12.z_t/100,data,12,colors='k',linewidths=0.5)
    cb = fig.colorbar(pm,ax=ax)
    cb.set_label(varinfo[vname]['label'])

North-south section along 145.9$^{\circ}$W

for vname in sorted(varinfo):
    data, dates = ds_sec_j225[vname][-12:,...].mean(dim='time'), ds_sec_j225.time[-1,...]
    data = np.ma.masked_where(np.isnan(data.values),data.values)
    data[data<0] = 0
    fig, ax = secfig(varinfo[vname]['title'])
    pm = ax.pcolormesh(ds_sec_j225.TLAT,-ds_sec_j225.z_t/100,data,cmap=cmap)
    cs = ax.contour(ds_sec_j225.TLAT,-ds_sec_j225.z_t/100,data,12,colors='k',linewidths=0.5)

    cb = fig.colorbar(pm,ax=ax)
    cb.set_label(varinfo[vname]['label'])