notebook.community

Edit and run

MOM6 diagnostics for KPP single column test case

Results from this notebook:

  1. Basic diagnostics of KPP boundary layer and prognostic fields.

Caveats regarding this notebook:

  1. This notebook is written for the MOM6-examples/ocean_only/SCM_KPP_tests.

Hopes for the use of this notebook:

  1. To provide a starting point to document single column model tests;
  2. To illustrate a self-contained iPython notebook of use for MOM6 analysis.

This iPython notebook was originally developed at NOAA/GFDL, and it is provided freely to the MOM6 community. GFDL scientists developing MOM6 make extensive use of Python for diagnostics. We solicit modifications/fixes that are useful to the MOM6 community.





In [1]:



    


import matplotlib.pyplot as plt                  # For plotting
import netCDF4                                   # For reading data
pylab.rcParams['figure.figsize'] = (17.0, 10.0)  # Large figures
#path = '../WSwPSBF.B20cm/'
#path = '../WSwPSBF.B1m/'
path = '../WSwPSBF.B10m/'
#path = '../WSwPSBF.BCM4/'
#path = './'

#deltaz = '20cm'
#deltaz = '1m'
deltaz = '10m'
#deltaz = 'CM4'

#fname_deltaz = '_20cm'
#fname_deltaz = '_1m'
fname_deltaz = '_10m'
#fname_deltaz = '_CM4'

dpi=200



    











In [2]:



    


visc = netCDF4.Dataset(path+'visc.nc')
for v in visc.variables: print v,

# recall data layout is (t,z,y,x)

time = visc.variables['Time'][:]

# KPP boundary layer depth as function of time (metre)
h    = visc.variables['KPP_OBLdepth'][:,0,0]

# KPP tracer diffusivity as function of time and depth (m2/sec)
Kt   = visc.variables['KPP_Kheat'][:,:,0,0]

# KPP velocity diffusivity as function of time and depth (m2/sec)
Ku   = visc.variables['KPP_Kv'][:,:,0,0]

# surface (and penetrating) buoyancy flux, as used by [CVmix] KPP (m2/s3)
KPP_buoyFlux = visc.variables['KPP_buoyFlux'][:,:,0,0]

# Temperature tendency due to non-local transport of heat, as calculated by KPP (K/s)
KPP_dTdt = visc.variables['KPP_dTdt'][:,:,0,0]



    


















    






xh yh zi Time nv zl xq yq Kd_interface KPP_uStar KPP_buoyFlux KPP_QminusSW KPP_netSalt KPP_OBLdepth KPP_sigma KPP_BulkRi KPP_Ws KPP_N KPP_N2 KPP_Vt2 KPP_BulkUz2 KPP_BulkDrho KPP_Kd_in KPP_Kheat KPP_Ksalt KPP_Kv KPP_NLtransport_heat KPP_NLtransport_salt KPP_dTdt KPP_dSdt KPP_Tsurf KPP_Ssurf KPP_Usurf KPP_Vsurf MLD_003 average_T1 average_T2 average_DT Time_bounds

















In [3]:



    


prog  = netCDF4.Dataset(path+'prog.nc')
for v in prog.variables: print v,
   
# depth of cell interface     
zi = prog.variables['zi'][:]

# depth of cell center 
zl = prog.variables['zl'][:]

# zonal velocity as function of time and depth
u  = prog.variables['u'][:,:,0,0]

# zonal velocity as function of time and depth
v  = prog.variables['v'][:,:,0,0]

# temperature as function of time and depth 
temp = prog.variables['temp'][:,:,0,0]



    


















    






xq yh zl Time xh yq zi u v h e temp salt wd

















In [4]:



    


# KPP related fields 

plt.subplot(2,2,1)
data  = Kt
field = np.ma.masked_array(data, mask=[data==0.])
#plt.contourf(time, -zi, field.T, levels=numpy.arange(0,.01,.001))
plt.pcolormesh(time, -zi, field.T, vmin=0, vmax=.01)
plt.colorbar()
plt.ylim((-50,0))
plt.ylabel('z (m)',fontsize=14)
plt.title(r'$\kappa_{\Theta} (m^2/s)$ with $\Delta z$ ='+deltaz,fontsize=18)

plt.subplot(2,2,2)
data  = Ku
field = np.ma.masked_array(data, mask=[data==0.])
plt.pcolormesh(time, -zi, field.T, vmin=0, vmax=.01)
#plt.contourf(time, -zi, field.T, levels=numpy.arange(0,.01,.001))
plt.colorbar()
plt.ylim((-50,0))
plt.ylabel('z (m)',fontsize=14)
plt.title(r'$\kappa_u (m^2/s)$ with $\Delta z$ ='+deltaz,fontsize=18)

plt.subplot(2,2,3)
data  = KPP_buoyFlux
field = np.ma.masked_array(data, mask=[data==0.])
plt.pcolormesh(time, -zi, field.T)
plt.colorbar()
plt.ylim((-50,0))
plt.ylabel('z (m)',fontsize=14)
plt.title(r'Buoyancy flux $(m^2/s^3)$ with $\Delta z$ ='+deltaz,fontsize=18)
plt.xlabel('Time (days)',fontsize=14)

plt.subplot(2,2,4)
data  = KPP_dTdt
field = np.ma.masked_array(data, mask=[data==0.])
plt.pcolormesh(time, -zl, field.T)
plt.colorbar()
plt.ylim((-50,0))
plt.ylabel('z (m)',fontsize=14)
plt.title(r'Non-local $\Theta$ tendency $(K/s)$ with $\Delta z$ ='+deltaz ,fontsize=18)
plt.xlabel('Time (days)',fontsize=14);



    


















    






/usr/local/x64/python/2.7.3/lib/python2.7/site-packages/matplotlib/colorbar.py:808: RuntimeWarning: invalid value encountered in divide
  z = np.take(y, i0) + (xn-np.take(b,i0))*dy/db
/usr/local/x64/python/2.7.3/lib/python2.7/site-packages/numpy/ma/core.py:3785: UserWarning: Warning: converting a masked element to nan.
  warnings.warn("Warning: converting a masked element to nan.")










    
























In [5]:



    


# Boundary layer depth and SST 

figure(4)
fig, ax = plt.subplots(figsize=(16, 10), dpi=dpi)
CS = plt.plot(time, -h,'-',linewidth=3)
yticks = [-50, -40, -30, -20, -10, 0]
ax.set_yticks(yticks)
xticks = [0, 3, 6, 9, 12, 15]
ax.set_xticks(xticks)
plt.grid()
plt.xlabel('Time (days)',fontsize=14)
plt.ylabel('z (m)',fontsize=14)
plt.title(r'KPP boundary layer depth for WSwPSBF.B with $\Delta z$ ='+deltaz,fontsize=16)
ax.set_yticklabels(["$%.1f$" % y for y in yticks], fontsize=18)
ax.set_xticklabels(["$%.1f$" % x for x in xticks], fontsize=18);

figure(5)
fig, ax = plt.subplots(figsize=(16, 10), dpi=dpi)
CS = plt.plot(time, temp[:,0],'-',linewidth=3)
yticks = [14.5, 15.0, 15.5, 16.0, 16.5, 17.0]
ax.set_yticks(yticks)
xticks = [0, 3, 6, 9, 12, 15]
ax.set_xticks(xticks)
plt.grid()
plt.xlabel('Time (days)',fontsize=14)
plt.ylabel('z (m)',fontsize=14)
plt.title(r'SST for WSwPSBF.B with $\Delta z$ ='+deltaz,fontsize=16)
ax.set_yticklabels(["$%.1f$" % y for y in yticks], fontsize=18)
ax.set_xticklabels(["$%.1f$" % x for x in xticks], fontsize=18);



    


















    








<matplotlib.figure.Figure at 0x5406fd0>










    

















    
























In [6]:



    


# u,v,temp

plt.subplot(2,2,1)
data  = u
field = np.ma.masked_array(data, mask=[data==0.])
#plt.contourf(time, -zl, field.T, levels=numpy.arange(-.15,.15,.01))
plt.pcolormesh(time, -zl, field.T, vmin=-.15, vmax=.15)
plt.colorbar()
plt.ylim((-50,0))
plt.ylabel('z (m)',fontsize=14)
plt.title('u (m/s) with $\Delta z$ ='+deltaz,fontsize=16)

plt.subplot(2,2,2)
data  = v
field = np.ma.masked_array(data, mask=[data==0.])
#plt.contourf(time, -zl, field.T, levels=numpy.arange(-.15,.15,.01))
plt.pcolormesh(time, -zl, field.T, vmin=-.15, vmax=.15)
plt.colorbar()
plt.ylim((-50,0))
plt.ylabel('z (m)',fontsize=14)
plt.title('v (m/s) with $\Delta z$ ='+deltaz,fontsize=16)

plt.subplot(2,2,3)
data  = temp
field = np.ma.masked_array(data, mask=[data==0.])
#CS = plt.contourf(time, -zl, temp.T,levels=numpy.arange(14,17,.025))
plt.pcolormesh(time, -zl, field.T, vmin=14, vmax=17)
plt.colorbar()
C = plt.contour(time, -zl, temp.T, 8,linewidth=.5, colors='black')
plt.clabel(C, inline=1, fontsize=10)
plt.ylim((-50,0))
plt.ylabel('z (m)',fontsize=14)
plt.title(r'$\Theta$ ($\degree$C) with $\Delta z$ ='+deltaz,fontsize=16)
plt.xlabel('Time (days)',fontsize=14);

plt.subplot(2,2,4)
difftemp = temp[:,:] - temp[0,:]
data  = difftemp
field = np.ma.masked_array(data, mask=[data==0.])
#CS = plt.contourf(time, -zl, field.T,levels=numpy.arange(-.6,2,.02))
plt.pcolormesh(time, -zl, field.T, vmin=-.6, vmax=2)
plt.colorbar()
C = plt.contour(time, -zl, field.T, 8,linewidth=.5, colors='black')
plt.clabel(C, inline=1, fontsize=10)
plt.ylim((-50,0))
plt.ylabel('z (m)',fontsize=14)
plt.title(r'$\Theta-\Theta(t=0)$ ($\degree$C) with $\Delta z$ ='+deltaz,fontsize=16);
plt.xlabel('Time (days)',fontsize=14);



    


















    
























In [16]:



    


# save figures to a file


figure(1)
fig = plt.figure(figsize=(16,10), dpi=dpi)
plt.contourf(time, -zl, temp.T, levels=numpy.arange(14,17,.1) )
plt.colorbar()
C = plt.contour(time, -zl, temp.T, linewidth=.5, colors='black',levels=numpy.arange(14,17,.1))
plt.clabel(C, inline=1, fontsize=10)
plt.ylim((-50,0))
plt.ylabel('z (m)',fontsize=14)
plt.title(r'$\Theta$ ($\degree$C) with $\Delta z$ ='+deltaz,fontsize=16)
plt.xlabel('Time (days)',fontsize=14)
fname = 'temp'+fname_deltaz+'.png'
fig.savefig(fname,dpi=dpi);

figure(2)
fig = plt.figure(figsize=(16,10), dpi=dpi)
difftemp = temp[:,:] - temp[0,:]
data  = difftemp
field = np.ma.masked_array(data, mask=[data==0.])
#CS = plt.contourf(time, -zl, difftemp.T,levels=numpy.arange(-.6,2,.02))
plt.pcolormesh(time, -zl, field.T, vmin=-.6, vmax=2)
plt.colorbar()
C = plt.contour(time, -zl, difftemp.T, 8,linewidth=.5, colors='black')
plt.clabel(C, inline=1, fontsize=10)
plt.ylim((-50,0))
plt.ylabel('z (m)',fontsize=14)
plt.xlabel('Time (days)',fontsize=14)
plt.title(r'$\Theta-\Theta(t=0)$ ($\degree$C) with $\Delta z$ ='+deltaz,fontsize=16);
plt.xlabel('Time (days)',fontsize=14)
fname = 'dtemp'+fname_deltaz+'.png'
fig.savefig(fname,dpi=dpi);

figure(3)
fig = plt.figure(figsize=(16,10), dpi=dpi)
data  = Kt
field = np.ma.masked_array(data, mask=[data==0.])
#CS = plt.contourf(time, -zi, Kt.T, levels=numpy.arange(0,.01,.001))
plt.pcolormesh(time, -zi, field.T, vmin=0, vmax=0.01)
plt.colorbar()
#C = plt.contour(time, -zi, Kt.T, levels=numpy.arange(0,.01,.001), linewidth=.5, colors='black')
plt.ylim((-50,0))
plt.ylabel('z (m)',fontsize=14)
plt.xlabel('Time (days)',fontsize=14)
plt.title(r'$\kappa_{\Theta}$ ($m^2/s$) with $\Delta z$ ='+deltaz,fontsize=16);
plt.xlabel('Time (days)',fontsize=14)
fname = 'KPP_diffusivity'+fname_deltaz+'.png'
fig.savefig(fname,dpi=dpi);


figure(4)
fig, ax = plt.subplots(figsize=(16, 10), dpi=dpi)
CS = plt.plot(time, -h,'-',linewidth=3)
yticks = [-50, -40, -30, -20, -10, 0]
ax.set_yticks(yticks)
xticks = [0, 3, 6, 9, 12, 15]
ax.set_xticks(xticks)
plt.grid()
plt.xlabel('Time (days)',fontsize=14)
plt.ylabel('z (m)',fontsize=14)
plt.title(r'KPP boundary layer depth for WSwPSBF.B with $\Delta z$ ='+deltaz,fontsize=16)
ax.set_yticklabels(["$%.1f$" % y for y in yticks], fontsize=18)
ax.set_xticklabels(["$%.1f$" % x for x in xticks], fontsize=18)
fname = 'KPP_bldepth'+fname_deltaz+'.png'
fig.savefig(fname,dpi=dpi);

figure(5)
fig, ax = plt.subplots(figsize=(16, 10), dpi=dpi)
CS = plt.plot(time, temp[:,0],'-',linewidth=3)
yticks = [14.5, 15.0, 15.5, 16.0, 16.5, 17.0]
ax.set_yticks(yticks)
xticks = [0, 3, 6, 9, 12, 15]
ax.set_xticks(xticks)
plt.grid()
plt.xlabel('Time (days)',fontsize=14)
plt.ylabel('z (m)',fontsize=14)
plt.title(r'SST for WSwPSBF.B with $\Delta z$ ='+deltaz,fontsize=16)
ax.set_yticklabels(["$%.1f$" % y for y in yticks], fontsize=18)
ax.set_xticklabels(["$%.1f$" % x for x in xticks], fontsize=18)
fname = 'SST'+fname_deltaz+'.png'
fig.savefig(fname,dpi=dpi);



    


















    








<matplotlib.figure.Figure at 0x4d0e690>










    

















    

















    

















    

















    
























In [7]:



    






    











In [ ]:

Content source: aidanheerdegen/MOM6-examples

Similar notebooks:

notebook.community | gallery | about