In [1]:
# Normal Python libraries
%load_ext autoreload
%autoreload 2
import numpy as np
import matplotlib
import netCDF4 as netCDF
import tracpy
import tracpy.plotting
from tracpy.tracpy_class import Tracpy
import datetime
matplotlib.rcParams.update({'font.size': 20})
%matplotlib qt
from matplotlib import pyplot as plt
In [34]:
# Data file
hisname = "runs/topoeddy/runew-04/" + "ocean_his.nc"
grdname = "runs/topoeddy/runew-04/" + "config/tes_grd.nc"
# Number of days to run the drifters.
ndays = 150
# Start date in date time formatting
date = datetime.datetime(1991, 3, 15, 0)
# Time between model outputs in seconds
tseas = 24*3600 # 24 hours between outputs, in seconds
# Time units
time_units = 'seconds since 1991-01-01'
# Sets a smaller limit than between model outputs for when to force interpolation if hasn't already occurred.
nsteps = 5
# Controls the sampling frequency of the drifter tracks.
N = 4
# Use ff = 1 for forward in time and ff = -1 for backward in time.
ff = 1
# Sub-grid scale diffusion
ah = 1. # m^2/s
av = 1. # m^2/s
# turbulence/diffusion flag
doturb = 0
# calculate volume transports for lagrangian stream function
dostream = 0
# simulation name, used for saving results into netcdf file
name = 'eddyshelf_test'
num_layers = 72
In [35]:
# Input starting locations as real space lon,lat locations
# DO NOT START AT ZERO - SOME BUG WITH VFLUX(2,0,1,1)
lon0, lat0 = np.meshgrid(np.linspace(20000, 300000, 50),
np.linspace(2000,40000,50)) # whole domain, 20 km
# for 3d flag, do3d=0 makes the run 2d and do3d=1 makes the run 3d
do3d = 1
# initial z-locations (array with shape of lon0)
if do3d:
z0 = -30 * np.ones_like(lon0)
zpar = 'fromZeta'
else:
z0 = 's'
zpar = num_layers-1
if dostream:
T0 = np.zeros_like(lon0).ravel()
else:
T0 = None
# Initialize Tracpy class
tp = Tracpy(hisname, grid_filename=grdname, name=name, tseas=tseas, ndays=ndays, nsteps=nsteps, usebasemap=False,
N=N, ff=ff, ah=ah, av=av, doturb=doturb, dostream=dostream, do3d=do3d, z0=z0, zpar=zpar, time_units=time_units,
usespherical=False)
# read in grid
#tp._readgrid()
# Eliminate points that are outside domain or in masked areas
#lon0, lat0 = tracpy.tools.check_points(lon0, lat0, tp.grid)
# initial z-locations (array with shape of lon0)
#if do3d:
# tp.z0 = -30 * np.ones_like(lon0)
# zpar = 'fromMSL'
#else:
# tp.z0 = 's'
# zpar = num_layers-1
In [36]:
# Note in timing that the grid was already read in
lonp, latp, zp, t, T0, U, V = tracpy.run.run(tp, date, lon0, lat0)
In [37]:
lonp, latp, zp, time = tracpy.inout.loadtracks('eddyshelf_test')
#U, V, lon0, lat0, T0 = tracpy.inout.loadtransport('eddyshelf_test')
tracpy.plotting.tracks(lonp, latp, tp.name, tp.grid, isll=False)
#step = 40
#plt.plot(lonp[1:-1:step,:].T, zp[1:-1:step,:].T);
In [38]:
tracpy.plotting.hist(lonp, latp, tp.name, grid=tp.grid, which='hexbin', bins=(50,50), isll=False)