0.5 hour frequency input: SalishSea_30m_20061214_20061215_grid_T/U/V.nc
1 hour frequency input: SalishSea_1h_20061214_20061215_grid_T/U/V.nc
4 hour frequency input input: SalishSea_4h_20061211_20061217_grid_T/U/V.nc
import netCDF4 as NC
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.cm as cm
from salishsea_tools import (nc_tools,viz_tools)
from salishsea_tools import tidetools
%matplotlib inline
Comparing Dec 14-15 at 1h and 4h interval outputs
delta_t = 3600 (hour unit, output), freq = 1 (1 hour sampling, output), nb_out = 47 (2 days)
nmax = 3, tunit = 3600 (hour unit, input), ntfic = 1 (1 hour sampling, input), lmt = 48 (2 days)
fl = 0.5 (start at Dec 14 at 00:00)
Note: Since initial time is 0.5, nb_out is 47 instead of 48
delta_t = 3600 (hour unit, output), freq = 4 (4 hour sampling, output), nb_out = 12 (2 days)
nmax = 3, tunit = 3600 (hour unit, input), ntfic = 4 (4 hour sampling, input), lmt = 42 (6 days)
fl = 18.5 (start at Dec 14 at 00:00 instead of Dec 11 or Dec 14 at 02:00)
Note: Condition 1 holds
For this set-up, the model has been re-run for Dec 14-15 and the data is now at intervals of 1 hour, 4 hours, and 30 minutes. The goal is to run Ariane at an interval of 30 minutes for all three scenarios. The input data for the Ariane runs are stored in a different directory now.
delta_t = 1800 (half hour unit, output), frequency = 1 (every half hour, output), nb_out = 95 (one less than 96 - 2 days)
tunit = 1800 (half hour unit, input), ntfic = 1 (half hour sampling, input), lmt = 96 (2dysx24hrsx2halves)
nmax = 3, fl = 0.5, Condition 1 holds
delta_t = 1800 (half hour unit, output), frequency = 1 (every half hour, output), nb_out = 95 (one less than 96 - 2 days)
tunit = 1800 (half hour unit, input), ntfic = 2 (1 hour sampling, input), lmt = 48 (2dysx24hrsx1hr)
nmax = 3, fl = 0.5, Condition 1 holds
delta_t = 1800 (half hour unit, output), frequency = 1 (every half hour, output), nb_out = 95 (one less than 96 - 2 days)
tunit = 1800 (half hour unit, input), ntfic = 8 (4 hour sampling, input), lmt = 12 (2 dysx6intervalsperday)
nmax = 3, fl = 0.5, Condition 1 holds
The purpose of this set is to find out if the time resolution has much of an effect on depth. We want to have the same set-up as the "half hour outputs" set - we just want to change the location i.e. edit the initial positions file. This file will be the same as the file used in http://nbviewer.ipython.org/urls/bitbucket.org/salishsea/analysis/raw/tip/Idalia/Particles_ThalwegIslands.ipynb
grid = NC.Dataset('/ocean/imachuca/MEOPAR/NEMO-forcing/grid/bathy_meter_SalishSea2.nc','r')
bathy, X, Y = tidetools.get_bathy_data(grid)
lats = grid.variables['nav_lat']
lons = grid.variables['nav_lon']
bath = grid.variables['Bathymetry']
h1R = NC.Dataset('/ocean/imachuca/MEOPAR/Ariane/results/time_res/riverparticle/1hour/ariane_trajectories_qualitative.nc','r')
lon1r=h1R.variables['traj_lon']
lat1r=h1R.variables['traj_lat']
dep1r=h1R.variables['traj_depth']
x1r=h1R.variables['init_x']
y1r=h1R.variables['init_y']
t1r=h1R.variables['traj_time']
h4R = NC.Dataset('/ocean/imachuca/MEOPAR/Ariane/results/time_res/riverparticle/4hour/ariane_trajectories_qualitative.nc','r')
lon4r=h4R.variables['traj_lon']
lat4r=h4R.variables['traj_lat']
dep4r=h4R.variables['traj_depth']
x4r=h4R.variables['init_x']
y4r=h4R.variables['init_y']
t4r=h4R.variables['traj_time']
Data below is for a run similar to above except for initial locations.
half = NC.Dataset('/ocean/imachuca/MEOPAR/Ariane/results/time_res/riverparticle2/30min/ariane_trajectories_qualitative.nc','r')
lonh=half.variables['traj_lon']
lath=half.variables['traj_lat']
deph=half.variables['traj_depth']
xh=half.variables['init_x']
yh=half.variables['init_y']
th=half.variables['traj_time']
one = NC.Dataset('/ocean/imachuca/MEOPAR/Ariane/results/time_res/riverparticle2/1hour/ariane_trajectories_qualitative.nc','r')
lon1=one.variables['traj_lon']
lat1=one.variables['traj_lat']
dep1=one.variables['traj_depth']
x1=one.variables['init_x']
y1=one.variables['init_y']
t1=one.variables['traj_time']
four = NC.Dataset('/ocean/imachuca/MEOPAR/Ariane/results/time_res/riverparticle2/4hour/ariane_trajectories_qualitative.nc','r')
lon4=four.variables['traj_lon']
lat4=four.variables['traj_lat']
dep4=four.variables['traj_depth']
x4=four.variables['init_x']
y4=four.variables['init_y']
t4=four.variables['traj_time']
Data below is for Ariane runs with half hour outputs
half_thal = NC.Dataset('/ocean/imachuca/MEOPAR/Ariane/results/time_res/thalweg/30min/ariane_trajectories_qualitative.nc','r')
lonh_thal=half_thal.variables['traj_lon']
lath_thal=half_thal.variables['traj_lat']
deph_thal=half_thal.variables['traj_depth']
xh_thal=half_thal.variables['init_x']
yh_thal=half_thal.variables['init_y']
th_thal=half_thal.variables['traj_time']
one_thal = NC.Dataset('/ocean/imachuca/MEOPAR/Ariane/results/time_res/thalweg/1hour/ariane_trajectories_qualitative.nc','r')
lon1_thal=one_thal.variables['traj_lon']
lat1_thal=one_thal.variables['traj_lat']
dep1_thal=one_thal.variables['traj_depth']
x1_thal=one_thal.variables['init_x']
y1_thal=one_thal.variables['init_y']
t1_thal=one_thal.variables['traj_time']
four_thal = NC.Dataset('/ocean/imachuca/MEOPAR/Ariane/results/time_res/thalweg/4hour/ariane_trajectories_qualitative.nc','r')
lon4_thal=four_thal.variables['traj_lon']
lat4_thal=four_thal.variables['traj_lat']
dep4_thal=four_thal.variables['traj_depth']
x4_thal=four_thal.variables['init_x']
y4_thal=four_thal.variables['init_y']
t4_thal=four_thal.variables['traj_time']
Both 1 hour and 4 hour scenarios have 3 particles at different depths (1 metre apart)
fig, ax = plt.subplots(1,1,figsize=(10,6.6))
n = np.arange(3)
colors = cm.winter(np.linspace(0, 1, len(n)))
labss = ['0-1m','1-2m','2-3m']
for N,c,labs in zip(n,colors,labss):
ax.scatter(lon1r[1:,N],lat1r[1:,N],color=c, label=labs)
ax.scatter(lon1r[0,N],lat1r[0,N],color='0.30',marker='s')
viz_tools.plot_land_mask(ax,grid,coords='map')
ax.set_xlim([-123.6,-122.8])
ax.set_ylim([48.5,49.5])
ax.set_title('1 hour')
ax.set_xlabel('lon')
ax.set_ylabel('lat')
ax.legend()
fig, ax = plt.subplots(1,1,figsize=(10,6.6))
n = np.arange(3)
colors = cm.autumn(np.linspace(0, 1, len(n)))
for N,c,labs in zip(n,colors,labss):
ax.scatter(lon4r[1:,N],lat4r[1:,N],color=c, label=labs)
ax.scatter(lon4r[0,N],lat4r[0,N],color='0.30',marker='s')
viz_tools.plot_land_mask(ax,grid,coords='map')
ax.set_xlim([-123.6,-122.8])
ax.set_ylim([48.5,49.5])
ax.set_title('4 hour')
ax.set_xlabel('lon')
ax.set_ylabel('lat')
ax.legend()
# *******************************************
n = np.arange(3)
colors1 = cm.winter(np.linspace(0, 1, len(n)))
colors4 = cm.autumn(np.linspace(0, 1, len(n)))
fig, ax = plt.subplots(1,1,figsize=(10,6.6))
for N,c1,c4 in zip(n,colors1,colors4):
ax.scatter(lon1r[1:,N],lat1r[1:,N],color=c1, label='1 hr')
ax.scatter(lon4r[1:,N],lat4r[1:,N],color=c4, label='4 hr')
viz_tools.plot_land_mask(ax,grid,coords='map')
ax.set_xlim([-123.4,-123])
ax.set_ylim([48.9,49.25])
ax.set_title('1 and 4 hour')
ax.set_xlabel('lon')
ax.set_ylabel('lat')
ax.legend()
<matplotlib.legend.Legend at 0x7effdd8f3e10>
fig, ax = plt.subplots(1,1,figsize=(15,10))
viz_tools.plot_land_mask(ax,grid,coords='map')
ax.set_xlim([-123.5,-123])
ax.set_ylim([48.96,49.17])
ax.set_title('[Grid 1 depth, 2 days] output=30 min, input=30 min/1 hour/4 hour')
ax.set_xlabel('lon')
ax.set_ylabel('lat')
labss = ['30 min','1 hour','4 hour']
lonn = [lonh,lon1,lon4]
latt = [lath,lat1,lat4]
colors = ['DarkGreen','DarkBlue','DeepPink']
for lon,lat,labs,c in zip(lonn, latt, labss, colors):
ax.scatter(lon[1:,0],lat[1:,0],color=c, marker='.', label=labs)
ax.scatter(lon[0,0],lat[0,0],color='0.30',marker='s')
ax.legend()
# ******
fig, ax = plt.subplots(1,1,figsize=(15,10))
viz_tools.plot_land_mask(ax,grid,coords='map')
ax.set_xlim([-123.3,-123])
ax.set_ylim([48.99,49.15])
ax.set_title('[Grid 2 depth, 2 days] output=30 min, input=30 min/1 hour/4 hour')
ax.set_xlabel('lon')
ax.set_ylabel('lat')
labss = ['30 min','1 hour','4 hour']
lonn = [lonh,lon1,lon4]
latt = [lath,lat1,lat4]
colors = ['DarkGreen','DarkBlue','DeepPink']
for lon,lat,labs,c in zip(lonn, latt, labss, colors):
ax.scatter(lon[1:,1],lat[1:,1],color=c, marker='.', label=labs)
ax.scatter(lon[0,1],lat[0,1],color='0.30',marker='s')
ax.legend()
# ******
fig, ax = plt.subplots(1,1,figsize=(15,10))
viz_tools.plot_land_mask(ax,grid,coords='map')
ax.set_xlim([-123.2,-123])
ax.set_ylim([49.02,49.15])
ax.set_title('[Grid 3 depth, 2 days] output=30 min, input=30 min/1 hour/4 hour')
ax.set_xlabel('lon')
ax.set_ylabel('lat')
labss = ['30 min','1 hour','4 hour']
lonn = [lonh,lon1,lon4]
latt = [lath,lat1,lat4]
colors = ['DarkGreen','DarkBlue','DeepPink']
for lon,lat,labs,c in zip(lonn, latt, labss, colors):
ax.scatter(lon[1:,2],lat[1:,2],color=c, marker='.', label=labs)
ax.scatter(lon[0,2],lat[0,2],color='0.30',marker='s')
ax.legend()
<matplotlib.legend.Legend at 0x7effdd8f9c50>
fig, axs = plt.subplots(1, 3, figsize=(24, 10), sharey=True)
n=np.arange(3)
titless = ['1','2','3']
for ax, N, titles in zip(axs, n, titless):
ax.plot(th_thal[0,N]*2, deph_thal[0,N]*-1, '-rs', label = 'Initial Position')
ax.plot(th_thal[1:,N]*2, deph_thal[1:,N]*-1, '-', label = '30 min', color='DarkGreen')
ax.plot(t1_thal[1:,N]*2, dep1_thal[1:,N]*-1, '-', label = '1 hour', color='DarkBlue')
ax.plot(t4_thal[1:,N]*2, dep4_thal[1:,N]*-1, '-', label = '4 hour', color='DeepPink')
ax.set_ylim(200, 0)
ax.legend(loc='lower left')
ax.set_xlabel('time [days]')
ax.set_ylabel('depth [m]')
ax.set_title(titles)
fig, axs = plt.subplots(1, 3, figsize=(24,24), sharey=True)
lonn = [lonh_thal,lon1_thal,lon4_thal]
latt = [lath_thal,lat1_thal,lat4_thal]
labss = ['30 min','1 hour','4 hour']
colors = ['DarkGreen','DarkBlue','DeepPink']
for ax,titles,N in zip(axs,titless,n):
aspect = viz_tools.set_aspect(ax, coords='map', lats=lats)
cmap = plt.get_cmap('Blues')
cmap.set_bad('burlywood')
mesh = ax.pcolormesh(lons[:], lats[:], bathy[:], cmap=cmap)
ax.set_xlim([-123.6,-122.9])
ax.set_ylim([48.6,49.05])
for lon,lat,c,labs in zip(lonn,latt,colors,labss):
ax.scatter(lon[1:,N],lat[1:,N],color=c, marker='.', label=labs)
ax.scatter(lon[0,N],lat[0,N],color='0.30',marker='s')
ax.legend()
ax.set_title(titles)
ax.set_xlabel('lon')
ax.set_ylabel('lat')
#
fig, axs = plt.subplots(1, 3, figsize=(24, 10), sharey=True)
m=n+3
titless = ['4','5','6']
for ax, M, titles in zip(axs, m, titless):
ax.plot(th_thal[0,M]*2, deph_thal[0,M]*-1, '-rs', label = 'Initial Position')
ax.plot(th_thal[1:,M]*2, deph_thal[1:,M]*-1, '-', label = '30 min', color='DarkGreen')
ax.plot(t1_thal[1:,M]*2, dep1_thal[1:,M]*-1, '-', label = '1 hour', color='DarkBlue')
ax.plot(t4_thal[1:,M]*2, dep4_thal[1:,M]*-1, '-', label = '4 hour', color='DeepPink')
ax.set_ylim(200, 0)
ax.legend(loc='lower left')
ax.set_xlabel('time [days]')
ax.set_ylabel('depth [m]')
ax.set_title(titles)
fig, axs = plt.subplots(1, 3, figsize=(24,24), sharey=True)
for ax,titles,M in zip(axs,titless,m):
aspect = viz_tools.set_aspect(ax, coords='map', lats=lats)
cmap = plt.get_cmap('Blues')
cmap.set_bad('burlywood')
mesh = ax.pcolormesh(lons[:], lats[:], bathy[:], cmap=cmap)
ax.set_xlim([-123.6,-122.8])
ax.set_ylim([48.35,49.0])
for lon,lat,c,labs in zip(lonn,latt,colors,labss):
ax.scatter(lon[1:,M],lat[1:,M],color=c, marker='.', label=labs)
ax.scatter(lon[0,M],lat[0,M],color='0.30',marker='s')
ax.legend()
ax.set_title(titles)
ax.set_xlabel('lon')
ax.set_ylabel('lat')