@@ -504,12 +504,8 @@ def get_uavg(Eddy, CS, collind, centlon_e, centlat_e, poly_eff,
504
504
any_inner_contours = False
505
505
506
506
citer = np .nditer (CS .cvalues , flags = ['c_index' ])
507
- #print '************************************************'
507
+
508
508
while not citer .finished :
509
-
510
- #print '--------------' * 10
511
- #print '--------------' * 10
512
- #print '\nciter.index ::', citer.index
513
509
## Get contour around centlon_e, centlat_e at level [collind:][iuavg]
514
510
#segii, poly_i = eddy_tracker.find_nearest_contour(
515
511
#CS.collections[citer.index], centlon_e, centlat_e)
@@ -518,7 +514,6 @@ def get_uavg(Eddy, CS, collind, centlon_e, centlat_e, poly_eff,
518
514
519
515
# Leave loop if no contours at level citer.index
520
516
if Eddy .swirl .level_slice is None :
521
- #print '------------------------------------------'
522
517
citer .iternext ()
523
518
continue
524
519
@@ -561,13 +556,12 @@ def get_uavg(Eddy, CS, collind, centlon_e, centlat_e, poly_eff,
561
556
citer .iternext ()
562
557
563
558
564
- if any_inner_contours :
559
+ if any_inner_contours : # set speed based contour parameters
565
560
cx , cy = Eddy .M (theseglon , theseglat )
566
- # Get speed based eddy radius (eddy_radius_s)
567
561
centx_s , centy_s , eddy_radius_s , junk = fit_circle (cx , cy )
568
562
centlon_s , centlat_s = Eddy .M .projtran (centx_s , centy_s , inverse = True )
569
563
570
- else :
564
+ else : # use the effective contour
571
565
centlon_s , centlat_s = centlon_e , centlat_e
572
566
eddy_radius_s = eddy_radius_e
573
567
inner_seglon , inner_seglat = theseglon , theseglat
@@ -576,6 +570,7 @@ def get_uavg(Eddy, CS, collind, centlon_e, centlat_e, poly_eff,
576
570
if not save_all_uavg :
577
571
return (uavg , centlon_s , centlat_s , eddy_radius_s ,
578
572
theseglon , theseglat , inner_seglon , inner_seglat )
573
+
579
574
else :
580
575
return (uavg , centlon_s , centlat_s , eddy_radius_s ,
581
576
theseglon , theseglat , inner_seglon , inner_seglat , all_uavg )
@@ -775,14 +770,19 @@ def collection_loop(CS, grd, rtime, A_list_obj, C_list_obj,
775
770
776
771
elif 'SLA' in Eddy .DIAGNOSTIC_TYPE :
777
772
773
+ args = (Eddy , CS , collind ,
774
+ centlon_e , centlat_e ,
775
+ cont , grd , eddy_radius_e )
776
+
778
777
if not Eddy .TRACK_EXTRA_VARIABLES :
779
- uavg , centlon_s , centlat_s , eddy_radius_s ,\
780
- contlon_s , contlat_s , inner_contlon , inner_contlat = get_uavg ( Eddy , CS , collind ,
781
- centlon_e , centlat_e , cont , grd , eddy_radius_e )
778
+ ( uavg , centlon_s , centlat_s ,
779
+ eddy_radius_s , contlon_s , contlat_s ,
780
+ inner_contlon , inner_contlat ) = get_uavg ( * args )
782
781
else :
783
- uavg , centlon_s , centlat_s , eddy_radius_s ,\
784
- contlon_s , contlat_s , inner_contlon , inner_contlat , uavg_profile = get_uavg (Eddy , CS , collind ,
785
- centlon_e , centlat_e , cont , grd , eddy_radius_e , save_all_uavg = True )
782
+ (uavg , centlon_s , centlat_s ,
783
+ eddy_radius_s , contlon_s , contlat_s ,
784
+ inner_contlon , inner_contlat ,
785
+ uavg_profile ) = get_uavg (* args , save_all_uavg = True )
786
786
787
787
788
788
centlon_lmi , centlat_lmi , junk , junk = fit_circle (inner_contlon ,
@@ -800,8 +800,6 @@ def collection_loop(CS, grd, rtime, A_list_obj, C_list_obj,
800
800
properties .rtime = rtime
801
801
properties .teke = teke
802
802
803
-
804
-
805
803
# Update Q eddy properties
806
804
if 'Q' in Eddy .DIAGNOSTIC_TYPE :
807
805
# We pass eddy_radius_e as a dummy for eddy_radius_s
@@ -913,9 +911,9 @@ def track_eddies(Eddy, first_record):
913
911
debug_dist = False
914
912
915
913
# We will need these in m for ellipse below
916
- old_x , old_y = Eddy .M (np .asarray (Eddy .old_lon ), np .asarray (Eddy .old_lat ))
917
- new_x , new_y = Eddy .M (np .asarray (Eddy .new_lon_tmp ),
918
- np .asarray (Eddy .new_lat_tmp ))
914
+ old_x , old_y = Eddy .M (np .array (Eddy .old_lon ), np .array (Eddy .old_lat ))
915
+ new_x , new_y = Eddy .M (np .array (Eddy .new_lon_tmp ),
916
+ np .array (Eddy .new_lat_tmp ))
919
917
920
918
X_old = np .array ([Eddy .old_lon , Eddy .old_lat ]).T
921
919
X_new = np .array ([Eddy .new_lon_tmp , Eddy .new_lat_tmp ]).T
@@ -1055,7 +1053,7 @@ def track_eddies(Eddy, first_record):
1055
1053
# corresponding new_eddy_inds set to False
1056
1054
new_eddy_inds [np .nonzero (Eddy .new_lon_tmp ==
1057
1055
Eddy .new_lon_tmp [new_ind ])] = False
1058
- dist_mat [:,new_ind ] = 1e9 # km
1056
+ dist_mat [:, new_ind ] = 1e9 # km
1059
1057
1060
1058
1061
1059
if Eddy .TRACK_EXTRA_VARIABLES :
@@ -1093,8 +1091,8 @@ def track_eddies(Eddy, first_record):
1093
1091
1094
1092
# Choice of using effective or speed-based...
1095
1093
delta_area = np .r_ [delta_area ,
1096
- np .abs (np .diff ([np .pi * (Eddy .old_radii_e [old_ind ]** 2 ),
1097
- np .pi * (new_rd_e [i ]** 2 )]))]
1094
+ np .abs (np .diff ([np .pi * (Eddy .old_radii_e [old_ind ]** 2 ),
1095
+ np .pi * (new_rd_e [i ]** 2 )]))]
1098
1096
delta_amp = np .r_ [delta_amp ,
1099
1097
np .abs (np .diff ([Eddy .old_amp [old_ind ], new_am [i ]]))]
1100
1098
@@ -1147,7 +1145,7 @@ def track_eddies(Eddy, first_record):
1147
1145
1148
1146
# Use backup_ind to reinsert distances into dist_mat for the unused eddy/eddies
1149
1147
for i , bind in enumerate (backup_ind [dx_unused ]):
1150
- dist_mat [:,bind ] = dist_mat_copy [:,bind ]
1148
+ dist_mat [:, bind ] = dist_mat_copy [:, bind ]
1151
1149
1152
1150
if debug_dist :
1153
1151
print 'backup_ind[dx_unused].shape' ,backup_ind [dx_unused ].shape
@@ -1229,10 +1227,14 @@ def accounting(Eddy, old_ind, centlon, centlat,
1229
1227
new_eddy : flag indicating a new eddy
1230
1228
first_record : flag indicating that we're on the first record
1231
1229
"""
1232
- if first_record : # is True then all eddies are new...
1230
+ if first_record : # is True then all eddies are new
1231
+ new_eddy = True
1233
1232
if Eddy .VERBOSE :
1234
1233
print '------ writing first record'
1235
- new_eddy = True
1234
+
1235
+ kwargs = {'temp' :cent_temp , 'salt' :cent_salt ,
1236
+ 'contour_e' :contour_e , 'contour_s' :contour_s ,
1237
+ 'uavg_profile' :uavg_profile , 'shape_error' :shape_error }
1236
1238
1237
1239
if not new_eddy : # it's an old (i.e., active) eddy
1238
1240
@@ -1257,19 +1259,7 @@ def accounting(Eddy, old_ind, centlon, centlat,
1257
1259
args = (old_ind , centlon , centlat , rtime , uavg , teke ,
1258
1260
eddy_radius_s , eddy_radius_e , amplitude )
1259
1261
1260
- if 'ROMS' in Eddy .DATATYPE :
1261
-
1262
- Eddy .update_track (* args ,
1263
- temp = cent_temp , salt = cent_salt ,
1264
- contour_e = contour_e , contour_s = contour_s ,
1265
- uavg_profile = uavg_profile , shape_error = shape_error )
1266
-
1267
- elif 'AVISO' in Eddy .DATATYPE :
1268
-
1269
- Eddy .update_track (* args ,
1270
- temp = None , salt = None ,
1271
- contour_e = contour_e , contour_s = contour_s ,
1272
- uavg_profile = uavg_profile , shape_error = shape_error )
1262
+ Eddy .update_track (* args , ** kwargs )
1273
1263
1274
1264
else : # it's a new eddy
1275
1265
@@ -1298,19 +1288,8 @@ def accounting(Eddy, old_ind, centlon, centlat,
1298
1288
1299
1289
args = (centlon , centlat , rtime , uavg , teke ,
1300
1290
eddy_radius_s , eddy_radius_e , amplitude )
1301
- kwargs = {'contour_e' :contour_e , 'contour_s' :contour_s ,
1302
- 'uavg_profile' :uavg_profile , 'shape_error' :shape_error }
1303
1291
1304
- if 'ROMS' in Eddy .DATATYPE :
1305
-
1306
- kwargs ['temp' ] = cent_temp
1307
- kwargs ['salt' ] = cent_salt
1308
-
1309
- Eddy .add_new_track (* args , ** kwargs )
1310
-
1311
- elif 'AVISO' in Eddy .DATATYPE :
1312
-
1313
- Eddy .add_new_track (* args , ** kwargs )
1292
+ Eddy .add_new_track (* args , ** kwargs )
1314
1293
1315
1294
Eddy .index += 1
1316
1295
@@ -1363,14 +1342,14 @@ def hann2d_fast(var, ni, nj):
1363
1342
#print 'jsz, isz',jsz, isz
1364
1343
var_ext = np .ma .zeros ((nj , ni )) # add 1-more line parallell to
1365
1344
var_ext [1 :- 1 , 1 :- 1 ] = var # each of 4-sides
1366
- var_ext [1 :- 1 , 0 ] = var [:,0 ] # duplicate W-side
1367
- var_ext [1 :- 1 ,- 1 ] = var [:,- 1 ] # duplicate E-side
1368
- var_ext [0 ,1 :- 1 ] = var [0 ,: ] # duplicate N-side
1369
- var_ext [- 1 ,1 :- 1 ] = var [- 1 ,: ] # duplicate S-side
1370
- var_ext [0 ,0 ] = np .nan # NW-corner
1371
- var_ext [0 ,- 1 ] = np .nan # NE-corner
1372
- var_ext [- 1 ,0 ] = np .nan # SW-corner
1373
- var_ext [- 1 ,- 1 ] = np .nan # SE-corner
1345
+ var_ext [1 :- 1 , 0 ] = var [:, 0 ] # duplicate W-side
1346
+ var_ext [1 :- 1 , - 1 ] = var [:, - 1 ] # duplicate E-side
1347
+ var_ext [0 , 1 :- 1 ] = var [0 ] # duplicate N-side
1348
+ var_ext [- 1 , 1 :- 1 ] = var [- 1 ] # duplicate S-side
1349
+ var_ext [0 , 0 ] = np .nan # NW-corner
1350
+ var_ext [0 , - 1 ] = np .nan # NE-corner
1351
+ var_ext [- 1 , 0 ] = np .nan # SW-corner
1352
+ var_ext [- 1 , - 1 ] = np .nan # SE-corner
1374
1353
1375
1354
# npts is used to count number of valid neighbors
1376
1355
npts = ne .evaluate ('var_ext * 0. + 1.' )
@@ -1383,7 +1362,7 @@ def hann2d_fast(var, ni, nj):
1383
1362
cc = np .ma .zeros ((var .shape ))
1384
1363
varS = np .ma .zeros ((var .shape ))
1385
1364
1386
- cc = npts [1 :nj - 1 , 1 :ni - 1 ] * (npts [0 :nj - 2 , 1 :ni - 1 ] + npts [2 :nj , 1 :ni - 1 ] +
1365
+ cc = npts [1 :nj - 1 , 1 :ni - 1 ] * (npts [0 :nj - 2 , 1 :ni - 1 ] + npts [2 :nj , 1 :ni - 1 ] +
1387
1366
npts [1 :nj - 1 , 0 :ni - 2 ] + npts [1 :nj - 1 , 2 :ni ])
1388
1367
1389
1368
varS = (var_ext [0 :nj - 2 , 1 :ni - 1 ] + var_ext [2 :nj , 1 :ni - 1 ] +
@@ -1405,7 +1384,7 @@ def hann2d_fast(var, ni, nj):
1405
1384
1406
1385
def get_circle (x0 , y0 , r , npts ):
1407
1386
"""
1408
- Return points on a circle, with specified (x0,y0) center and radius
1387
+ Return points on a circle, with specified (x0, y0) center and radius
1409
1388
(and optional number of points too!).
1410
1389
1411
1390
Input : 1 - x0, scalar, center X of circle
0 commit comments