import netCDF4 as nc
import numpy as np
from matplotlib import pyplot as plt
import cmocean as cm
import datetime as dt
import glob
from salishsea_tools import viz_tools, visualisations
%matplotlib inline
Look at the difference in biological fields between Feb 6 2015 and the same date in other years to try to get an idea of why resetting to 2015 biological fields on that date in 2017 worked so well.
f2015=nc.Dataset('/results2/SalishSea/hindcast/11jan15/SalishSea_1h_20150111_20150111_ptrc_T.nc')
f2016=nc.Dataset('/results2/SalishSea/hindcast/11jan16/SalishSea_1h_20160111_20160111_ptrc_T.nc')
f2017=nc.Dataset('/results2/SalishSea/hindcast/11jan17/SalishSea_1h_20170111_20170111_ptrc_T.nc')
f2018=nc.Dataset('/results2/SalishSea/hindcast/11jan18/SalishSea_1h_20180111_20180111_ptrc_T.nc')
f2015.variables.keys()
odict_keys(['nav_lat', 'nav_lon', 'bounds_lon', 'bounds_lat', 'area', 'deptht', 'deptht_bounds', 'nitrate', 'time_centered', 'time_centered_bounds', 'time_counter', 'time_counter_bounds', 'ammonium', 'silicon', 'diatoms', 'flagellates', 'ciliates', 'microzooplankton', 'dissolved_organic_nitrogen', 'particulate_organic_nitrogen', 'biogenic_silicon', 'Fraser_tracer', 'mesozooplankton'])
cmap = plt.get_cmap(cm.cm.balance) #('nipy_spectral')
cmap.set_bad('burlywood')
with nc.Dataset('/ocean/eolson/MEOPAR/NEMO-forcing/grid/mesh_mask201702_noLPE.nc') as fm:
tmask=np.copy(fm.variables['tmask'])
for var in ('nitrate', 'ammonium', 'silicon', 'diatoms', 'flagellates', 'ciliates', 'microzooplankton', 'dissolved_organic_nitrogen', 'particulate_organic_nitrogen', 'biogenic_silicon', 'mesozooplankton'):
fig,ax = plt.subplots(1,5,figsize=(24,8))
for ii in range(0,5):
viz_tools.set_aspect(ax[ii])
lev=0
iax=ax[0]
diff=f2017.variables[var][0,lev,200:800,50:]-f2015.variables[var][0,lev,200:800,50:]
cl=np.mean(f2017.variables[var][0,lev,200:800,50:])+np.std(f2017.variables[var][0,lev,200:800,50:])
mesh=iax.pcolormesh(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,diff),vmin=-cl,vmax=cl,cmap=cmap)
iax.set_title(var+' 2017-2015 \nz='+str(f2015.variables['deptht'][lev])+' m', fontsize=14)
cbar = fig.colorbar(mesh,ax=iax)
cbar.ax.tick_params(labelsize=14)
iax.get_xaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.set_xlim((50,348))
lev=9
iax=ax[1]
diff=f2017.variables[var][0,lev,200:800,50:]-f2015.variables[var][0,lev,200:800,50:]
cl=np.mean(f2017.variables[var][0,lev,200:800,50:])+np.std(f2017.variables[var][0,lev,200:800,50:])
mesh=iax.pcolormesh(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,diff),vmin=-cl,vmax=cl,cmap=cmap)
iax.set_title(var+' 2017-2015 \nz='+str(f2015.variables['deptht'][lev])+' m', fontsize=14)
cbar = fig.colorbar(mesh,ax=iax)
cbar.ax.tick_params(labelsize=14)
iax.get_xaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.set_xlim((50,348))
lev=20
iax=ax[2]
diff=f2017.variables[var][0,lev,200:800,50:]-f2015.variables[var][0,lev,200:800,50:]
cl=np.mean(f2017.variables[var][0,lev,200:800,50:])+np.std(f2017.variables[var][0,lev,200:800,50:])
mesh=iax.pcolormesh(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,diff),vmin=-cl,vmax=cl,cmap=cmap)
iax.set_title(var+' 2017-2015 \nz='+str(f2015.variables['deptht'][lev])+' m', fontsize=14)
cbar = fig.colorbar(mesh,ax=iax)
cbar.ax.tick_params(labelsize=14)
iax.get_xaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.set_xlim((50,348))
lev=25
iax=ax[3]
diff=f2017.variables[var][0,lev,200:800,50:]-f2015.variables[var][0,lev,200:800,50:]
cl=np.mean(f2017.variables[var][0,lev,200:800,50:])+np.std(f2017.variables[var][0,lev,200:800,50:])
mesh=iax.pcolormesh(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,diff),vmin=-cl,vmax=cl,cmap=cmap)
iax.set_title(var+' 2017-2015 \nz='+str(f2015.variables['deptht'][lev])+' m', fontsize=14)
cbar = fig.colorbar(mesh,ax=iax)
cbar.ax.tick_params(labelsize=14)
iax.get_xaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.set_xlim((50,348))
lev=30
iax=ax[4]
diff=f2017.variables[var][0,lev,200:800,50:]-f2015.variables[var][0,lev,200:800,50:]
cl=np.mean(f2017.variables[var][0,lev,200:800,50:])+np.std(f2017.variables[var][0,lev,200:800,50:])
#cl=np.max(np.abs(diff))
mesh=iax.pcolormesh(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,diff),vmin=-cl,vmax=cl,cmap=cmap)
iax.set_title(var+' 2017-2015 \nz='+str(f2015.variables['deptht'][lev])+' m', fontsize=14)
cbar = fig.colorbar(mesh,ax=iax)
cbar.ax.tick_params(labelsize=14)
iax.get_xaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.set_xlim((50,348))
for var in ('nitrate', 'ammonium', 'silicon', 'diatoms', 'flagellates', 'ciliates', 'microzooplankton', 'dissolved_organic_nitrogen', 'particulate_organic_nitrogen', 'biogenic_silicon', 'mesozooplankton'):
fig,ax = plt.subplots(1,5,figsize=(24,8))
for ii in range(0,5):
viz_tools.set_aspect(ax[ii])
lev=0
iax=ax[0]
diff=f2016.variables[var][0,lev,200:800,50:]-f2015.variables[var][0,lev,200:800,50:]
cl=np.mean(f2016.variables[var][0,lev,200:800,50:])+np.std(f2016.variables[var][0,lev,200:800,50:])
mesh=iax.pcolormesh(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,diff),vmin=-cl,vmax=cl,cmap=cmap)
iax.set_title(var+' 2016-2015 \nz='+str(f2015.variables['deptht'][lev])+' m', fontsize=14)
cbar = fig.colorbar(mesh,ax=iax)
cbar.ax.tick_params(labelsize=14)
iax.get_xaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.set_xlim((50,348))
lev=9
iax=ax[1]
diff=f2016.variables[var][0,lev,200:800,50:]-f2016.variables[var][0,lev,200:800,50:]
cl=np.mean(f2016.variables[var][0,lev,200:800,50:])+np.std(f2016.variables[var][0,lev,200:800,50:])
mesh=iax.pcolormesh(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,diff),vmin=-cl,vmax=cl,cmap=cmap)
iax.set_title(var+' 2016-2015 \nz='+str(f2015.variables['deptht'][lev])+' m', fontsize=14)
cbar = fig.colorbar(mesh,ax=iax)
cbar.ax.tick_params(labelsize=14)
iax.get_xaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.set_xlim((50,348))
lev=20
iax=ax[2]
diff=f2016.variables[var][0,lev,200:800,50:]-f2015.variables[var][0,lev,200:800,50:]
cl=np.mean(f2016.variables[var][0,lev,200:800,50:])+np.std(f2016.variables[var][0,lev,200:800,50:])
mesh=iax.pcolormesh(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,diff),vmin=-cl,vmax=cl,cmap=cmap)
iax.set_title(var+' 2016-2015 \nz='+str(f2015.variables['deptht'][lev])+' m', fontsize=14)
cbar = fig.colorbar(mesh,ax=iax)
cbar.ax.tick_params(labelsize=14)
iax.get_xaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.set_xlim((50,348))
lev=25
iax=ax[3]
diff=f2016.variables[var][0,lev,200:800,50:]-f2015.variables[var][0,lev,200:800,50:]
cl=np.mean(f2016.variables[var][0,lev,200:800,50:])+np.std(f2016.variables[var][0,lev,200:800,50:])
mesh=iax.pcolormesh(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,diff),vmin=-cl,vmax=cl,cmap=cmap)
iax.set_title(var+' 2016-2015 \nz='+str(f2015.variables['deptht'][lev])+' m', fontsize=14)
cbar = fig.colorbar(mesh,ax=iax)
cbar.ax.tick_params(labelsize=14)
iax.get_xaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.set_xlim((50,348))
lev=30
iax=ax[4]
diff=f2016.variables[var][0,lev,200:800,50:]-f2015.variables[var][0,lev,200:800,50:]
cl=np.mean(f2016.variables[var][0,lev,200:800,50:])+np.std(f2016.variables[var][0,lev,200:800,50:])
#cl=np.max(np.abs(diff))
mesh=iax.pcolormesh(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,diff),vmin=-cl,vmax=cl,cmap=cmap)
iax.set_title(var+' 2016-2015 \nz='+str(f2015.variables['deptht'][lev])+' m', fontsize=14)
cbar = fig.colorbar(mesh,ax=iax)
cbar.ax.tick_params(labelsize=14)
iax.get_xaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.set_xlim((50,348))
for var in ('nitrate', 'ammonium', 'silicon', 'diatoms', 'flagellates', 'ciliates', 'microzooplankton', 'dissolved_organic_nitrogen', 'particulate_organic_nitrogen', 'biogenic_silicon', 'mesozooplankton'):
fig,ax = plt.subplots(1,5,figsize=(24,8))
for ii in range(0,5):
viz_tools.set_aspect(ax[ii])
lev=0
iax=ax[0]
diff=f2017.variables[var][0,lev,200:800,50:]-f2016.variables[var][0,lev,200:800,50:]
cl=np.mean(f2017.variables[var][0,lev,200:800,50:])+np.std(f2017.variables[var][0,lev,200:800,50:])
mesh=iax.pcolormesh(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,diff),vmin=-cl,vmax=cl,cmap=cmap)
iax.set_title(var+' 2017-2016 \nz='+str(f2016.variables['deptht'][lev])+' m', fontsize=14)
cbar = fig.colorbar(mesh,ax=iax)
cbar.ax.tick_params(labelsize=14)
iax.get_xaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.set_xlim((50,348))
lev=9
iax=ax[1]
diff=f2017.variables[var][0,lev,200:800,50:]-f2016.variables[var][0,lev,200:800,50:]
cl=np.mean(f2017.variables[var][0,lev,200:800,50:])+np.std(f2017.variables[var][0,lev,200:800,50:])
mesh=iax.pcolormesh(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,diff),vmin=-cl,vmax=cl,cmap=cmap)
iax.set_title(var+' 2017-2016 \nz='+str(f2016.variables['deptht'][lev])+' m', fontsize=14)
cbar = fig.colorbar(mesh,ax=iax)
cbar.ax.tick_params(labelsize=14)
iax.get_xaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.set_xlim((50,348))
lev=20
iax=ax[2]
diff=f2017.variables[var][0,lev,200:800,50:]-f2016.variables[var][0,lev,200:800,50:]
cl=np.mean(f2017.variables[var][0,lev,200:800,50:])+np.std(f2017.variables[var][0,lev,200:800,50:])
mesh=iax.pcolormesh(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,diff),vmin=-cl,vmax=cl,cmap=cmap)
iax.set_title(var+' 2017-2016 \nz='+str(f2016.variables['deptht'][lev])+' m', fontsize=14)
cbar = fig.colorbar(mesh,ax=iax)
cbar.ax.tick_params(labelsize=14)
iax.get_xaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.set_xlim((50,348))
lev=25
iax=ax[3]
diff=f2017.variables[var][0,lev,200:800,50:]-f2016.variables[var][0,lev,200:800,50:]
cl=np.mean(f2017.variables[var][0,lev,200:800,50:])+np.std(f2017.variables[var][0,lev,200:800,50:])
mesh=iax.pcolormesh(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,diff),vmin=-cl,vmax=cl,cmap=cmap)
iax.set_title(var+' 2017-2016 \nz='+str(f2016.variables['deptht'][lev])+' m', fontsize=14)
cbar = fig.colorbar(mesh,ax=iax)
cbar.ax.tick_params(labelsize=14)
iax.get_xaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.set_xlim((50,348))
lev=30
iax=ax[4]
diff=f2017.variables[var][0,lev,200:800,50:]-f2016.variables[var][0,lev,200:800,50:]
cl=np.mean(f2017.variables[var][0,lev,200:800,50:])+np.std(f2017.variables[var][0,lev,200:800,50:])
#cl=np.max(np.abs(diff))
mesh=iax.pcolormesh(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,diff),vmin=-cl,vmax=cl,cmap=cmap)
iax.set_title(var+' 2017-2016 \nz='+str(f2016.variables['deptht'][lev])+' m', fontsize=14)
cbar = fig.colorbar(mesh,ax=iax)
cbar.ax.tick_params(labelsize=14)
iax.get_xaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.set_xlim((50,348))
for var in ('diatoms', 'flagellates', 'ciliates'):
fig,ax=plt.subplots(1,4,figsize=(10,5))
lev=0
iax=ax[0]
mesh=iax.pcolormesh(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,f2015.variables[var][0,lev,200:800,50:]),vmin=0,vmax=.01,cmap=cm.cm.thermal)
iax.set_title(var+' 2015 \nz='+str(f2015.variables['deptht'][lev])+' m', fontsize=14)
cbar = fig.colorbar(mesh,ax=iax)
cbar.ax.tick_params(labelsize=14)
iax.get_xaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.set_xlim((50,348))
iax=ax[1]
mesh=iax.pcolormesh(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,f2016.variables[var][0,lev,200:800,50:]),vmin=0,vmax=.01,cmap=cm.cm.thermal)
iax.set_title(var+' 2016 \nz='+str(f2016.variables['deptht'][lev])+' m', fontsize=14)
cbar = fig.colorbar(mesh,ax=iax)
cbar.ax.tick_params(labelsize=14)
iax.get_xaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.set_xlim((50,348))
iax=ax[2]
mesh=iax.pcolormesh(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,f2017.variables[var][0,lev,200:800,50:]),vmin=0,vmax=.01,cmap=cm.cm.thermal)
iax.set_title(var+' 2017 \nz='+str(f2017.variables['deptht'][lev])+' m', fontsize=14)
cbar = fig.colorbar(mesh,ax=iax)
cbar.ax.tick_params(labelsize=14)
iax.get_xaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.set_xlim((50,348))
iax=ax[3]
mesh=iax.pcolormesh(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,f2018.variables[var][0,lev,200:800,50:]),vmin=0,vmax=.01,cmap=cm.cm.thermal)
iax.set_title(var+' 2018 \nz='+str(f2017.variables['deptht'][lev])+' m', fontsize=14)
cbar = fig.colorbar(mesh,ax=iax)
cbar.ax.tick_params(labelsize=14)
iax.get_xaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.get_yaxis().set_visible(False)
iax.set_xlim((50,348))
print('2015 mean '+var+':',np.mean(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,f2015.variables[var][0,lev,200:800,50:])))
print('2016 mean '+var+':',np.mean(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,f2016.variables[var][0,lev,200:800,50:])))
print('2017 mean '+var+':',np.mean(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,f2017.variables[var][0,lev,200:800,50:])))
print('2018 mean '+var+':',np.mean(np.ma.masked_where(tmask[0,lev,200:800,50:]==0,f2018.variables[var][0,lev,200:800,50:])))
2015 mean diatoms: 0.38869720412796394 2016 mean diatoms: 0.009574445690582155 2017 mean diatoms: 0.003947068404121194 2018 mean diatoms: 0.010426087807798173 2015 mean flagellates: 0.40299248287124967 2016 mean flagellates: 0.5706577760075081 2017 mean flagellates: 0.4635242299481408 2018 mean flagellates: 0.46379328672619224 2015 mean ciliates: 0.08286303943507582 2016 mean ciliates: 0.2248189199196145 2017 mean ciliates: 0.17116243714666735 2018 mean ciliates: 0.17624804206420452
fb=nc.Dataset('/data/eolson/MEOPAR/NEMO-forcing-new/grid/bathymetry_201702.nc')
fm=nc.Dataset('/ocean/eolson/MEOPAR/NEMO-forcing/grid/mesh_mask201702_noLPE.nc')
for var in ('nitrate', 'ammonium', 'silicon', 'diatoms', 'flagellates', 'ciliates', 'microzooplankton', 'dissolved_organic_nitrogen', 'particulate_organic_nitrogen', 'biogenic_silicon', 'mesozooplankton'):
fig,ax = plt.subplots(1,1,figsize=(18,6))
diff=f2017.variables[var][0,:,:,:]-f2015.variables[var][0,:,:,:]
cl=np.mean(np.abs(diff))+np.std(np.abs(diff))
m=visualisations.contour_thalweg(ax,diff,fb,fm,method='pcolormesh',cmap=cm.cm.balance,mesh_args={'vmin':-cl,'vmax':cl})
ax.set_title(var+ '2017-2015')