In [1]:
%matplotlib inline
from matplotlib import pylab
import matplotlib.pyplot as plt
import netCDF4 as NC
import numpy as np
In [2]:
fT = NC.Dataset('../../../nemo-forcing/open_boundaries/source_files/bc4SA_extractedFromBCCBaseline.nc','r')
In [3]:
print fT.file_format
NETCDF3_CLASSIC
In [4]:
for dim in fT.dimensions.values():
print dim
<type 'netCDF4.Dimension'>: name = 'xi_rho', size = 4 <type 'netCDF4.Dimension'>: name = 'xi_u', size = 3 <type 'netCDF4.Dimension'>: name = 'xi_v', size = 4 <type 'netCDF4.Dimension'>: name = 'eta_rho', size = 15 <type 'netCDF4.Dimension'>: name = 'eta_u', size = 15 <type 'netCDF4.Dimension'>: name = 'eta_v', size = 14 <type 'netCDF4.Dimension'>: name = 's_rho', size = 30 <type 'netCDF4.Dimension'>: name = 'ocean_time', size = 336
In [5]:
for var in fT.variables.values():
print var
<type 'netCDF4.Variable'>
float64 ocean_time(ocean_time)
long_name: time since model initialization
units: seconds since 1990-01-01 00:00:00
calendar: gregorian
field: time, scalar, series
unlimited dimensions:
current shape = (336,)
<type 'netCDF4.Variable'>
float64 temp(ocean_time, s_rho, eta_rho, xi_rho)
long_name: potential temperature
units: Celsius
field: temperature, scalar, series
time: ocean_time
unlimited dimensions:
current shape = (336, 30, 15, 4)
<type 'netCDF4.Variable'>
float64 salt(ocean_time, s_rho, eta_rho, xi_rho)
long_name: salinity
units: PSU
field: salinity, scalar, series
time: ocean_time
unlimited dimensions:
current shape = (336, 30, 15, 4)
<type 'netCDF4.Variable'>
float64 zeta(ocean_time, eta_rho, xi_rho)
long_name: sea surface height
units: meter
field: SSH, scalar, series
time: ocean_time
unlimited dimensions:
current shape = (336, 15, 4)
<type 'netCDF4.Variable'>
float64 u(ocean_time, s_rho, eta_u, xi_u)
long_name: u-velocity
units: meter second-1
field: u-velocity, scalar, series
time: ocean_time
unlimited dimensions:
current shape = (336, 30, 15, 3)
<type 'netCDF4.Variable'>
float64 v(ocean_time, s_rho, eta_v, xi_v)
long_name: v-velocity
units: meter second-1
field: v-velocity, scalar, series
time: ocean_time
unlimited dimensions:
current shape = (336, 30, 14, 4)
<type 'netCDF4.Variable'>
float64 zr(ocean_time, s_rho, eta_rho, xi_rho)
long_name: s-level depth at rho points
units: m
field: depth, scalar, series
time: ocean_time
unlimited dimensions:
current shape = (336, 30, 15, 4)
<type 'netCDF4.Variable'>
float64 lon_rho(eta_rho, xi_rho)
long_name: Longitude at rho-points
units: Degrees
field: longitude, scalar, series
unlimited dimensions:
current shape = (15, 4)
<type 'netCDF4.Variable'>
float64 lat_rho(eta_rho, xi_rho)
long_name: Latitude at rho-points
units: Degrees
field: latitude, scalar, series
unlimited dimensions:
current shape = (15, 4)
In [6]:
latR = fT.variables['lat_rho'][:]
lonR = fT.variables['lon_rho'][:]
In [7]:
import sys
sys.path.append('../../bathymetry')
import bathy_tools
In [142]:
fB = NC.Dataset('../../../nemo-forcing/grid/bathy_meter_SalishSea.nc','r')
D = fB.variables['Bathymetry'][:]
lat = fB.variables['nav_lat'][:]
lon = fB.variables['nav_lon'][:]
In [10]:
import datetime
test = fT.variables['ocean_time'][:]
start_time = datetime.datetime(1990,1,1,0,0,0)
time = np.array([start_time + datetime.timedelta(seconds=i) for i in test])
In [11]:
t=15
fig = bathy_tools.plot_colourmesh(
fB, 'Salish Sea Bathymetry, '+str(time[t]),
axis_limits=(-126.3, -122.15, 47.05, 51))
#plt.plot(lonR,latR,'kx')
plt.contourf(lonR,latR,fT.variables['salt'][t,29,:,:],np.arange(28,33,0.1),cmap='RdPu')
plt.xlim((-125.5,-124.0))
plt.ylim((48.0,49.0))
Out[11]:
(48.0, 49.0)
In [12]:
zr = fT.variables['zr'][:]
plt.pcolor(zr[0,:,:14,2])
plt.colorbar()
Out[12]:
<matplotlib.colorbar.Colorbar instance at 0x54a9fc8>
In [145]:
#plot the depth of the bottom cell from south (x=0) to north (x=14) at time=0
for tt in range(0,29):
bottom_depth = zr[tt,0,:14,2]
depth1 = zr[tt,10,:14,2]
depth2 = zr[tt,20,:14,2]
depth3 = zr[tt,28,:14,2]
plt.plot(range(0,14),bottom_depth,'b')
plt.plot(range(0,14),depth1,'r')
plt.plot(range(0,14),depth2,'g')
plt.plot(range(0,14),depth3,'m')
plt.ylabel('depth of cell [m]')
plt.xlabel('eta_rho [number of grid cells]')
In [14]:
salt = fT.variables['salt'][:]
plt.contourf(np.transpose(salt[:,:,4,2]))
plt.colorbar()
Out[14]:
<matplotlib.colorbar.Colorbar instance at 0x63fe050>
In [77]:
plt.plot(time,fT.variables['salt'][:,29,7,1],'r')
plt.plot(time,fT.variables['salt'][:,15,7,1],'b')
plt.plot(time,fT.variables['salt'][:,4,7,1],'g')
plt.plot(time,fT.variables['salt'][:,29,14,1],'r:')
plt.plot(time,fT.variables['salt'][:,15,14,1],'b:')
plt.plot(time,fT.variables['salt'][:,4,14,1],'g:')
Out[77]:
[<matplotlib.lines.Line2D at 0x34850510>]
In [15]:
#get the min and max for the colorbar
def get_cbar_lims(salt):
cbarmin = round(np.min(salt))
cbarmax = round(np.max(salt[salt<9e35]))+1
return cbarmin, cbarmax
In [111]:
#define a function to find the index of the model cell closest to the required depth
def find_vert_cell_for_req_depth(req_depth,time,eta_rho_ind,xi_rho_ind):
all_depths = zr[time,:,eta_rho_ind,xi_rho_ind]
if np.logical_and(np.min(all_depths)<0,np.min(all_depths)<req_depth):
bottom_depth = np.min(all_depths)
#a range of tolerances is needed because depths do not increase evenly
tols = np.arange(0,-10,-0.01)
tol_count = 0
vert_cell_index = np.array([])
while np.size(vert_cell_index)==0:
tol = tols[tol_count]
vert_cell_index = np.where(np.logical_and(all_depths>(req_depth+tol),all_depths<(req_depth-tol)))
tol_count = tol_count+1
vert_cell_index = vert_cell_index[0].astype(int)
else:
vert_cell_index = np.ma.zeros((1,))
vert_cell_index[0] = np.ma.masked
return vert_cell_index[0]
In [127]:
a, b = fT.variables['lat_rho'].shape
req_depths = np.arange(0,-300,-20)
depth_ind = np.ma.zeros((a,b,len(req_depths)))
time = 0
for dep in range(0,len(req_depths)):
for i in range(0,a):
for j in range(0,b):
depth_ind[i,j,dep] = find_vert_cell_for_req_depth(req_depths[dep],time,i,j)
depth_ind = depth_ind.astype(int)
print(depth_ind[:,:,0])
print a,b
[[-- -- 29 29] [-- 29 29 29] [29 29 29 29] [29 29 29 29] [29 29 29 29] [29 29 29 29] [29 29 29 29] [29 29 29 29] [29 29 29 29] [29 29 29 29] [29 29 29 29] [29 29 29 29] [29 29 29 --] [29 29 29 --] [29 29 -- --]] 15 4
In [97]:
#define a function to plot all salinity data at a given depth, with 24 subplots per year
def salinity_contours_all_years(depth,salt,time):
cbarmin, cbarmax = get_cbar_lims(salt[:,depth,:,:])
numplots = 24
for t in range(0,len(time)/numplots):
plt.figure(figsize=(20,9))
for k in range(1,numplots+1):
fig = plt.subplot(1,numplots,k)
plt.contourf(salt[k+t*numplots-1,depth,:,:],np.arange(cbarmin,cbarmax,0.5))
plt.title(time[k+t*numplots-1].strftime('%m/%y'))
fig.axes.get_xaxis().set_ticks([])
fig.axes.get_yaxis().set_ticks([])
plt.colorbar()
In [ ]:
In [ ]:
In [116]:
salinity_contours_all_years(29,salt,time)
In [98]:
salinity_contours_all_years(20,salt,time)
In [117]:
salinity_contours_all_years(0,salt,time)
In [118]:
salinity_contours_all_years(0,temp,time)
In [256]:
salinity_contours_all_years(15,temp,time)
*Diane Masson's POM model results*
Salinity at s-rho = 29 (surface)
- varies between 27 and 32.5 psu throughout the year
- is higher in summer (30-32.5psu) and lower in winter (27-30psu)
Salinity at s-rho = 20
- varies between 28 psu and 33.5 psu throughout the year.
- is higher in summer (30.5-33.5psu) and lower in winter (28.5-31psu)
Salinity at s-rho = 0 (bottom)
- varies between 30.5 psu and 34.5 psu throughout the year.
- is higher in summer (??-??psu) and lower in winter (??-??psu)
Salinity at all depths
tends to be higher towards the south side of Juan de Fuca Strait and lower towards the north side of Juan de Fuca Strait.
high and low salinity water tends to appear at south end of the boundary and move into Juan de Fuca Strait over a period of weeks to months.
there are some anomalous salinity plots:
- 04/96 where salinity decreases in surface cells (s-rho > 20)
- 01/02 where salinity decreases at all depths
- maybe 05/05 where salinity slightly decreases in surface cells (s-rho > 20)
Potential temperature at s-rho = 0 (bottom)
- varies between 6 and 11 degrees throughout the year
In [ ]: