How finely resolved should our data be? 30 min, 1 hour, 4 hours?

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

References

0. Imports

In [15]:
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

1. Ariane set-up (original)

Comparing Dec 14-15 at 1h and 4h interval outputs

Condition 1: delta_t × frequency × nb_output < tunit × ntfic × lmt (if initial time is 0.5)
Condition 2: delta_t × frequency × nb_output < tunit × ntfic × (lmt + 0.5 - max(fl)) (if initial time is greater than 0.5)

a) 1 Hour

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

b) 4 Hour

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

2. Ariane set-up (half hour outputs)

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.

a) 0.5 Hour input // 0.5 Hour output

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

b) 1 Hour input // 0.5 Hour output

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

c) 4 Hour input // 0.5 Hour output

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

3. Ariane set-up (thalweg)

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

4. Data

In [16]:
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']
In [17]:
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']
In [18]:
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.

In [19]:
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']
In [20]:
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']
In [21]:
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

In [22]:
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']
In [23]:
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']
In [24]:
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']

5. Comparison (original)

Both 1 hour and 4 hour scenarios have 3 particles at different depths (1 metre apart)

In [25]:
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()
Out[25]:
<matplotlib.legend.Legend at 0x7effdd8f3e10>

6. Comparison: Trajectories at 30 minute intervals when their inputs were at 30 min, 1 hour, and 4 hour intervals

In [26]:
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()
Out[26]:
<matplotlib.legend.Legend at 0x7effdd8f9c50>

6. Comparison - Thalweg: Trajectories/Depths at 30 minute intervals when their inputs were at 30 min, 1 hour, and 4 hour intervals

In [32]:
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')