Skip to content

Commit a82c258

Browse files
evanmasonevanmason
committed
Back at CLS
1 parent 56c37ec commit a82c258

File tree

4 files changed

+243
-99
lines changed

4 files changed

+243
-99
lines changed

‎eddy_tracker_configuration.yaml‎

Lines changed: 16 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -12,32 +12,36 @@ DOMAIN:
1212
THE_DOMAIN: 'Global'
1313
#THE_DOMAIN = 'BlackSea'
1414
#THE_DOMAIN = 'MedSea' # not yet implemented
15-
LONMIN: -38.
16-
LONMAX: -25.
17-
LATMIN: 30.
18-
LATMAX: 38.
19-
DATE_STR: 20020101
20-
DATE_END: 20020630
15+
LONMIN: -65. # -38.
16+
LONMAX: -5.5 # -25.
17+
LATMIN: 11.5 # 30.
18+
LATMAX: 38.5 # 38.
19+
DATE_STR: 19930101
20+
DATE_END: 20131231
2121

2222
# Specify the AVISO product
2323
AVISO:
2424
AVISO_DT14: Yes # Use new AVISO DT14 data (e.g., Capet et al, 2014)
2525
AVISO_FILES: 'dt_global_twosat_msla_h_????????_20140106.nc'
26-
AVISO_DT14_SUBSAMP: No # subsample the dqily AVISO DT14
27-
DAYS_BTWN_RECORDS: 3 # sampling rate (days)
26+
AVISO_DT14_SUBSAMP: Yes # subsample the dqily AVISO DT14
27+
DAYS_BTWN_RECORDS: 7 # sampling rate (days)
2828

2929
PATHS:
3030
# Path to AVISO data
31-
DATA_DIR: '/data/PVA/Externe/global/delayed-time/grids/msla/two-sat-merged/h/2002/'
31+
#DATA_DIR: '/data/PVA/Externe/global/delayed-time/grids/msla/two-sat-merged/h/2002/'
32+
DATA_DIR: '/marula/emason/data/altimetry/global/delayed-time/grids/msla/two-sat-merged/h/'
3233
# Path and filename of Chelton et al (1998) Rossby radius data
3334
# Obtain file from:
3435
# http://www-po.coas.oregonstate.edu/research/po/research/rossby_radius/
35-
RW_PATH: '/homelocal/emason/rossrad.dat'
36+
#RW_PATH: '/homelocal/emason/rossrad.dat'
37+
RW_PATH: '/home/emason/data_tmp/chelton_eddies/rossrad.dat'
3638

3739
# Path for saving of outputs
38-
SAVE_DIR: '/homelocal/emason/outputs/'
40+
#SAVE_DIR: '/homelocal/emason/outputs/'
41+
# SAVE_DIR: '/marula/emason/aviso_eddy_tracking/fd0cdb4b2bd9/'
42+
SAVE_DIR: '/marula/emason/aviso_eddy_tracking/fd0cdb4b2bd9/'
3943

40-
# Reference Julian day (Julian date at Jan 1, 1992)
44+
# Reference Julian day (Julian date at Jan 1, 1992)
4145
JDAY_REFERENCE: 2448623.
4246

4347
# Define contouring parameters

‎make_eddy_track_AVISO.py‎

Lines changed: 4 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -478,7 +478,7 @@ def __init__(self, AVISO_FILE, LONMIN, LONMAX, LATMIN, LATMAX,
478478
try: # new AVISO (2014)
479479
self._lon = self.read_nc(AVISO_FILE, 'lon')
480480
self._lat = self.read_nc(AVISO_FILE, 'lat')
481-
self.fillval = self.read_nc_att(AVISO_FILE, 'SLA', '_FillValue')
481+
self.fillval = self.read_nc_att(AVISO_FILE, 'sla', '_FillValue')
482482
base_date = self.read_nc_att(AVISO_FILE, 'time', 'units')
483483
self.base_date = dt.date2num(
484484
parser.parse(base_date.split(' ')[2:4][0]))
@@ -514,9 +514,9 @@ def get_AVISO_data(self, AVISO_FILE):
514514
if self.ZERO_CROSSING:
515515

516516
try: # new AVISO (2014)
517-
ssh1 = self.read_nc(AVISO_FILE, 'SLA',
517+
ssh1 = self.read_nc(AVISO_FILE, 'sla',
518518
indices='[:, self.jp0:self.jp1, :self.ip0]')
519-
ssh0 = self.read_nc(AVISO_FILE, 'SLA',
519+
ssh0 = self.read_nc(AVISO_FILE, 'sla',
520520
indices='[:, self.jp0:self.jp1, self.ip1:]')
521521
ssh0, ssh1 = ssh0.squeeze(), ssh1.squeeze()
522522
ssh0 *= 100. # m to cm
@@ -533,7 +533,7 @@ def get_AVISO_data(self, AVISO_FILE):
533533
else:
534534

535535
try: # new AVISO (2014)
536-
zeta = self.read_nc(AVISO_FILE, 'SLA',
536+
zeta = self.read_nc(AVISO_FILE, 'sla',
537537
indices='[:, self.jp0:self.jp1, self.ip0:self.ip1]')
538538
zeta = zeta.squeeze()
539539
zeta *= 100. # m to cm

‎make_eddy_tracker_list_obj.py‎

Lines changed: 52 additions & 46 deletions
Original file line numberDiff line numberDiff line change
@@ -514,7 +514,7 @@ def __init__(self, DATATYPE, TRACK_DURATION_MIN, TRACK_EXTRA_VARIABLES):
514514
self.index = 0 # counter
515515
self.ncind = 0 # index to write to nc files, will increase and increase...
516516
self.ch_index = 0 # index for Chelton style nc files
517-
517+
self.PAD = 2
518518
# Check for a correct configuration
519519
assert DATATYPE in ('ROMS', 'AVISO'), "Unknown string in 'DATATYPE' parameter"
520520

@@ -956,7 +956,7 @@ def write2netcdf(self, rtime):
956956
self._reduce_inactive_tracks()
957957

958958
# Update old_lon and old_lat...
959-
self.old_lon = list(self.new_lon[lasti:])
959+
self.old_lon = self.new_lon[lasti:]
960960
self.old_lat = self.new_lat[lasti:]
961961
self.old_radii_s = self.new_radii_s[lasti:]
962962
self.old_radii_e = self.new_radii_e[lasti:]
@@ -1003,7 +1003,8 @@ def get_distances(self, centlon, centlat):
10031003
clat = self.old_lat.copy()
10041004
clon.fill(centlon)
10051005
clat.fill(centlat)
1006-
return haversine_distance_vector(clon, clat, self.old_lon, self.old_lat)
1006+
return haversine_distance_vector(clon, clat,
1007+
self.old_lon, self.old_lat)
10071008

10081009

10091010
def insert_at_index(self, xarr, ind, x):
@@ -1028,90 +1029,93 @@ def insert_at_index(self, xarr, ind, x):
10281029
try:
10291030
self.new_lon[ind] = x
10301031
except:
1031-
self.new_lon.extend([0] * (ind - len(self.new_lon) + 1))
1032+
self.new_lon.extend([0] *
1033+
(ind - len(self.new_lon) + 1))
10321034
self.new_lon[ind] = x
10331035
elif strcompare('new_lat', xarr):
10341036
try:
10351037
self.new_lat[ind] = x
10361038
except:
1037-
self.new_lat.extend([0] * (ind - len(self.new_lat) + 1))
1039+
self.new_lat.extend([0] *
1040+
(ind - len(self.new_lat) + 1))
10381041
self.new_lat[ind] = x
10391042
elif strcompare('new_radii_s', xarr):
10401043
try:
10411044
self.new_radii_s[ind] = x
10421045
except:
1043-
self.new_radii_s.extend([0] * (ind - len(self.new_radii_s) + 1))
1046+
self.new_radii_s.extend([0] *
1047+
(ind - len(self.new_radii_s) + 1))
10441048
self.new_radii_s[ind] = x
10451049
elif strcompare('new_radii_e', xarr):
10461050
try:
10471051
self.new_radii_e[ind] = x
10481052
except:
1049-
self.new_radii_e.extend([0] * (ind - len(self.new_radii_e) + 1))
1053+
self.new_radii_e.extend([0] *
1054+
(ind - len(self.new_radii_e) + 1))
10501055
self.new_radii_e[ind] = x
10511056
elif strcompare('new_amp', xarr):
10521057
try:
10531058
self.new_amp[ind] = x
10541059
except:
1055-
self.new_amp.extend([0] * (ind - len(self.new_amp) + 1))
1060+
self.new_amp.extend([0] *
1061+
(ind - len(self.new_amp) + 1))
10561062
self.new_amp[ind] = x
10571063
elif strcompare('new_uavg', xarr):
10581064
try:
10591065
self.new_uavg[ind] = x
10601066
except:
1061-
self.new_uavg.extend([0] * (ind - len(self.new_uavg) + 1))
1067+
self.new_uavg.extend([0] *
1068+
(ind - len(self.new_uavg) + 1))
10621069
self.new_uavg[ind] = x
10631070
elif strcompare('new_teke', xarr):
10641071
try:
10651072
self.new_teke[ind] = x
10661073
except:
1067-
self.new_teke.extend([0] * (ind - len(self.new_teke) + 1))
1074+
self.new_teke.extend([0] *
1075+
(ind - len(self.new_teke) + 1))
10681076
self.new_teke[ind] = x
10691077
elif strcompare('new_temp', xarr):
10701078
try:
10711079
self.new_temp[ind] = x
10721080
except:
1073-
self.new_temp.extend([0] * (ind - len(self.new_temp) + 1))
1081+
self.new_temp.extend([0] *
1082+
(ind - len(self.new_temp) + 1))
10741083
self.new_temp[ind] = x
10751084
elif strcompare('new_salt', xarr):
10761085
try:
10771086
self.new_salt[ind] = x
10781087
except:
1079-
self.new_salt.extend([0] * (ind - len(self.new_salt) + 1))
1088+
self.new_salt.extend([0] *
1089+
(ind - len(self.new_salt) + 1))
10801090
self.new_salt[ind] = x
10811091
elif strcompare('new_shape_error', xarr):
10821092
try:
10831093
self.new_shape_error[ind] = x
10841094
except:
1085-
self.new_shape_error.extend([0] * (ind - len(self.new_shape_error) + 1))
1095+
self.new_shape_error.extend([0] *
1096+
(ind - len(self.new_shape_error) + 1))
10861097
self.new_shape_error[ind] = x
1087-
#elif strcompare('new_bounds', xarr):
1088-
#if ind < tmp.shape[0]:
1089-
#newsize = tmp.size
1090-
#else:
1091-
#newsize = ind + 1
1092-
#self.new_bounds = np.zeros((newsize, 4))
1093-
#self.new_bounds[:tmp.shape[0]] = tmp
1094-
#self.new_bounds[ind] = x
1095-
10961098

10971099
elif strcompare('new_contour_e', xarr):
10981100
try:
10991101
self.new_contour_e[ind] = x
11001102
except:
1101-
#self.new_contour_e += [0] * (ind - len(self.new_contour_e) + 1)
1102-
self.new_contour_e.append([0] * (ind - len(self.new_contour_e) + 1))
1103+
self.new_contour_e.append([0] *
1104+
(ind - len(self.new_contour_e) + 1))
11031105
self.new_contour_e[ind] = x
11041106
elif strcompare('new_contour_s', xarr):
11051107
try:
11061108
self.new_contour_s[ind] = x
11071109
except:
1108-
self.new_contour_s.append([0] * (ind - len(self.new_contour_s) + 1))
1110+
self.new_contour_s.append([0] *
1111+
(ind - len(self.new_contour_s) + 1))
11091112
self.new_contour_s[ind] = x
11101113
elif strcompare('new_uavg_profile', xarr):
11111114
try:
11121115
self.new_uavg_profile[ind] = x
11131116
except:
1114-
self.new_uavg_profile.append([0] * (ind - len(self.new_uavg_profile) + 1))
1117+
self.new_uavg_profile.append([0] *
1118+
(ind - len(self.new_uavg_profile) + 1))
11151119
self.new_uavg_profile[ind] = x
11161120
else:
11171121
raise Exception
@@ -1148,7 +1152,7 @@ def insert_at_index(self, xarr, ind, x):
11481152
#self.jmax = np.max([a_j, b_j, c_j, d_j]) + 5
11491153
#return self
11501154

1151-
def set_bounds(self, contlon, contlat, radius, i, j, grd):
1155+
def set_bounds(self, contlon, contlat, grd):
11521156
"""
11531157
Get indices to a bounding box around the eddy
11541158
WARNING won't work for a rotated grid
@@ -1164,18 +1168,19 @@ def set_bounds(self, contlon, contlat, radius, i, j, grd):
11641168
jarr = np.array([bl_j, tl_j, br_j, tr_j])
11651169
self.imin, self.imax = iarr.min(), iarr.max()
11661170
self.jmin, self.jmax = jarr.min(), jarr.max()
1167-
pad = 2
1171+
11681172
# For indexing the mins must not be less than zero
1169-
self.imin = np.maximum(self.imin - pad, 0)
1170-
self.jmin = np.maximum(self.jmin - pad, 0)
1171-
self.imax += pad + 1
1172-
self.jmax += pad + 1
1173+
self.imin = np.maximum(self.imin - self.PAD, 0)
1174+
self.jmin = np.maximum(self.jmin - self.PAD, 0)
1175+
self.imax += self.PAD + 1
1176+
self.jmax += self.PAD + 1
11731177
return self
11741178

11751179

11761180
def set_mask_eff(self, contour, grd):
11771181
"""
1178-
Set points within bounding box around eddy and calculate mask
1182+
Set points within bounding box around eddy and calculate
1183+
mask for effective contour
11791184
"""
11801185
self.points = np.array([grd.lon()[self.jmin:self.jmax,
11811186
self.imin:self.imax].ravel(),
@@ -1296,7 +1301,8 @@ def _get_rlongwave_spd(self, x, y):
12961301
class SearchEllipse (object):
12971302

12981303
def __init__(self, domain, DAYS_BTWN_RECORDS, rw_path, limits):
1299-
"""Class to construct a search ellipse/circle around a specified point.
1304+
"""
1305+
Class to construct a search ellipse/circle around a specified point.
13001306
13011307
"""
13021308
self.domain = domain
@@ -1333,13 +1339,13 @@ def _set_global_ellipse(self):
13331339
"""
13341340
self._set_east_ellipse()._set_west_ellipse()
13351341
e_verts = self.east_ellipse.get_verts()
1336-
e_size = e_verts[:,0].size
1342+
e_size = e_verts[:, 0].size
13371343
e_size *= 0.5
13381344
w_verts = self.west_ellipse.get_verts()
1339-
w_size = w_verts[:,0].size
1345+
w_size = w_verts[:, 0].size
13401346
w_size *= 0.5
1341-
ew_x = np.hstack((e_verts[:e_size,0], w_verts[w_size:,0]))
1342-
ew_y = np.hstack((e_verts[:e_size,1], w_verts[w_size:,1]))
1347+
ew_x = np.hstack((e_verts[:e_size, 0], w_verts[w_size:, 0]))
1348+
ew_y = np.hstack((e_verts[:e_size, 1], w_verts[w_size:, 1]))
13431349
self.ellipse_path = path.Path(np.array([ew_x, ew_y]).T)
13441350
return self#.ellipse_path
13451351

@@ -1350,8 +1356,8 @@ def _set_black_sea_ellipse(self):
13501356
self.black_sea_ellipse = patch.Ellipse((self.x, self.y),
13511357
2. * self.rw_d_mod, 2. * self.rw_d_mod)
13521358
verts = self.black_sea_ellipse.get_verts()
1353-
self.ellipse_path = path.Path(np.array([verts[:,0],
1354-
verts[:,1]]).T)
1359+
self.ellipse_path = path.Path(np.array([verts[:, 0],
1360+
verts[:, 1]]).T)
13551361
return self
13561362

13571363

@@ -1362,15 +1368,17 @@ def set_search_ellipse(self, x, y):
13621368
self.y = y
13631369

13641370
if self.domain in ('Global', 'ROMS'):
1365-
self.rw_d[:] = self.rwv.get_rwdistance(x, y, self.DAYS_BTWN_RECORDS)
1371+
self.rw_d[:] = self.rwv.get_rwdistance(x, y,
1372+
self.DAYS_BTWN_RECORDS)
13661373
self.rw_d_mod[:] = 1.75
13671374
self.rw_d_mod *= self.rw_d
13681375
self.rw_d_mod[:] = np.maximum(self.rw_d_mod, self.n_s_minor)
13691376
self.rw_d_mod *= 2.
13701377
self._set_global_ellipse()
13711378
elif 'BlackSea' in self.domain:
13721379
self.rw_d_mod[:] = 1.75
1373-
self.rw_d[:] = self.rwv.get_rwdistance(x, y, self.DAYS_BTWN_RECORDS)
1380+
self.rw_d[:] = self.rwv.get_rwdistance(x, y,
1381+
self.DAYS_BTWN_RECORDS)
13741382
self.rw_d_mod *= self.rw_d
13751383
self._set_black_sea_ellipse()
13761384
else: Exception
@@ -1386,10 +1394,8 @@ def view_search_ellipse(self, Eddy):
13861394
ax = plt.subplot(111)
13871395
ax.set_title('Rossby def. rad %s m' %self.rw_d[0])
13881396
Eddy.M.scatter(self.x, self.y, c='b')
1389-
#Eddy.M.plot(ee[:,0], ee[:,1], 'b')
1390-
#Eddy.M.plot(ww[:,0], ww[:,1], 'g')
1391-
Eddy.M.plot(self.ellipse_path.vertices[:,0],
1392-
self.ellipse_path.vertices[:,1], 'r')
1397+
Eddy.M.plot(self.ellipse_path.vertices[:, 0],
1398+
self.ellipse_path.vertices[:, 1], 'r')
13931399
Eddy.M.drawcoastlines()
13941400
plt.show()
13951401

0 commit comments

Comments
 (0)