Skip to content

Commit ff5ba45

Browse files
committed
Various optimisation
1 parent 8ad4b97 commit ff5ba45

File tree

3 files changed

+211
-196
lines changed

3 files changed

+211
-196
lines changed

make_eddy_track_AVISO.py

Lines changed: 10 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -799,6 +799,7 @@ def pcol_2dxy(self, x, y):
799799
#date_str, date_end = 19980107, 19991110 #
800800
#date_str, date_end = 20081107, 20100110 #
801801
#date_str, date_end = 19921014, 20120718 #
802+
date_str, date_end = 20020101, 20020110 #
802803
date_str, date_end = 20020101, 20020630 #
803804

804805
# Choose type of diagnostic: either q-parameter ('Q') or sea level anomaly ('sla')
@@ -811,12 +812,14 @@ def pcol_2dxy(self, x, y):
811812
#savedir = '/marula/emason/aviso_eddy_tracking/pablo_exp/'
812813
#savedir = '/marula/emason/aviso_eddy_tracking/new_AVISO_test/'
813814
#savedir = '/marula/emason/aviso_eddy_tracking/new_AVISO_SUBSAMP-3days/'
814-
savedir = '/marula/emason/aviso_eddy_tracking/junk/'
815+
#savedir = '/marula/emason/aviso_eddy_tracking/junk/'
816+
savedir = '/marula/emason/aviso_eddy_tracking/junk2/'
815817
#savedir = '/marula/emason/aviso_eddy_tracking/new_AVISO_test/BlackSea/'
816818
#savedir = '/marula/emason/aviso_eddy_tracking/Corridor_V3_Dec2014/'
817819
#savedir = '/marula/emason/aviso_eddy_tracking/CLS_test_1_acd66ba338a9/'
818820
#savedir = '/marula/emason/aviso_eddy_tracking/CLS_test_2/'
819821
#savedir = '/marula/emason/aviso_eddy_tracking/CLS_test_3/'
822+
#savedir = '/marula/emason/aviso_eddy_tracking/CLS_test_4/'
820823
#savedir = '/path/to/save/your/outputs/'
821824

822825

@@ -917,6 +920,11 @@ def pcol_2dxy(self, x, y):
917920
lonmax = -25.
918921
latmin = 30.
919922
latmax = 38.
923+
924+
#lonmin = .15 # big test
925+
#lonmax = 359.15
926+
#latmin = -70.
927+
#latmax = 70.
920928

921929

922930
elif the_domain in 'MedSea':
@@ -1159,6 +1167,7 @@ def pcol_2dxy(self, x, y):
11591167

11601168
# Loop through the AVISO files...
11611169
start = True
1170+
11621171
start_time = time.time()
11631172

11641173
print '\nStart tracking'

make_eddy_tracker_list_obj.py

Lines changed: 108 additions & 71 deletions
Original file line numberDiff line numberDiff line change
@@ -37,6 +37,7 @@
3737
#from make_eddy_track import *
3838
from py_eddy_tracker_classes import *
3939

40+
#import find_closest_point_on_leg as fcpl
4041

4142
def haversine_distance_vector(lon1, lat1, lon2, lat2):
4243
'''
@@ -135,32 +136,34 @@ def strcompare(str1, str2):
135136

136137

137138

138-
def _find_closest_point_on_leg(p1, p2, p0):
139-
"""find closest point to p0 on line segment connecting p1 and p2"""
139+
#def _find_closest_point_on_leg(p1, p2, p0):
140+
#"""find closest point to p0 on line segment connecting p1 and p2"""
140141

141-
# handle degenerate case
142-
if np.all(p2 == p1):
143-
d = p0 - p1
144-
d **= 2
145-
return d.sum()
142+
##print type(p1), type(p2), type(p0)
146143

147-
d21 = p2 - p1
148-
d01 = p0 - p1
144+
## handle degenerate case
145+
#if np.all(p2 == p1):
146+
#d = p0 - p1
147+
#d **= 2
148+
#return d.sum()
149+
150+
#d21 = p2 - p1
151+
#d01 = p0 - p1
149152

150-
# project on to line segment to find closest point
151-
proj = np.dot(d01, d21)
152-
proj /= np.dot(d21, d21)
153-
if proj < 0:
154-
proj = 0
155-
elif proj > 1:
156-
proj = 1
157-
pc = p1 + proj * d21
153+
## project on to line segment to find closest point
154+
#proj = np.dot(d01, d21)
155+
#proj /= np.dot(d21, d21)
156+
#if proj < 0:
157+
#proj = 0
158+
#elif proj > 1:
159+
#proj = 1
160+
#pc = p1 + proj * d21
158161

159-
# find squared distance
160-
d = pc - p0
161-
d **= 2
162+
## find squared distance
163+
#d = pc - p0
164+
#d **= 2
162165

163-
return d.sum()#, pc
166+
#return d.sum()#, pc
164167

165168

166169
#def _find_closest_point_on_leg(p1, p2, p0):
@@ -222,73 +225,107 @@ def _find_closest_point_on_leg(p1, p2, p0):
222225

223226
#return dmin#, xcmin, legmin)
224227

225-
def _find_closest_point_on_path(lc, point):
226-
"""
227-
lc: coordinates of vertices
228-
point: coordinates of test point
229-
"""
230-
#import matplotlib.mlab as mlab; print 'delete me'
231-
# Find index of closest vertex for this segment
232-
ds = np.sum((lc - point[None, :])**2, 1)
233-
imin = ds.argmin()
234-
235-
dmin = np.inf
236-
closed = np.alltrue([lc[0,0] == lc[-1,0], lc[0,1] == lc[-1,1]])
237-
#closedd = mlab.is_closed_polygon(lc)
228+
#def _find_closest_point_on_path(lc, point):
229+
#"""
230+
#lc: coordinates of vertices
231+
#point: coordinates of test point
232+
#"""
233+
##import matplotlib.mlab as mlab; print 'delete me'
234+
## Find index of closest vertex for this segment
235+
#ds = np.sum((lc - point[None, :])**2, 1)
236+
#imin = ds.argmin()
237+
238+
#dmin = np.inf
239+
#closed = np.alltrue([lc[0,0] == lc[-1,0], lc[0,1] == lc[-1,1]])
240+
##closedd = mlab.is_closed_polygon(lc)
238241

239-
#assert closed == closedd, 'should be same'
240-
#print '***************************'
242+
##assert closed == closedd, 'should be same'
243+
##print '***************************'
241244

242-
# Build list of legs before and after this vertex
243-
legs = []
244-
if imin > 0 or closed:
245-
legs.append(((imin-1) % len(lc), imin))
246-
if imin < len(lc) - 1 or closed:
247-
legs.append((imin, (imin+1) % len(lc)))
248-
249-
for leg in legs:
250-
#d, xc = _find_closest_point_on_leg(lc[leg[0]], lc[leg[1]], point)
251-
d = _find_closest_point_on_leg(lc[leg[0]], lc[leg[1]], point)
252-
if d < dmin:
253-
dmin = d
245+
## Build list of legs before and after this vertex
246+
#legs = []
247+
#if imin > 0 or closed:
248+
#legs.append(((imin-1) % len(lc), imin))
249+
#if imin < len(lc) - 1 or closed:
250+
#legs.append((imin, (imin+1) % len(lc)))
254251

255-
return dmin
252+
#for leg in legs:
253+
##d, xc = _find_closest_point_on_leg(lc[leg[0]], lc[leg[1]], point)
254+
##dd = _find_closest_point_on_leg(lc[leg[0]], lc[leg[1]], point)
255+
##print type(point), type(lc[leg[0]]), type(lc[leg[1]])
256+
#d = fcpl.find_closest_point_on_leg(lc[leg[0]], lc[leg[1]], point)
257+
258+
##assert dd == d, 'not the same'
259+
260+
#if d < dmin:
261+
#dmin = d
256262

263+
#return dmin
257264

258-
def find_nearest_contour(contour_obj, x, y, indices):
259-
"""
260-
Finds contour that is closest to a point.
261265

262-
Returns a tuple containing the contour, segment of minimum point.
266+
#def find_nearest_contour(contour_obj, x, y, indices):
267+
#"""
268+
#Finds contour that is closest to a point.
263269

264-
Call signature::
270+
#Returns a tuple containing the contour, segment of minimum point.
265271

266-
conmin, segmin = find_nearest_contour(contour_obj, x, y, indices)
267-
"""
268-
dmin = np.inf
269-
#conmin = None
270-
segmin = None
272+
#Call signature::
273+
274+
#conmin, segmin = find_nearest_contour(contour_obj, x, y, indices)
275+
#"""
276+
#dmin = np.inf
277+
##conmin = None
278+
#segmin = None
271279

272-
point = np.array([x, y])
280+
#point = np.array([x, y])
273281

274-
for icon in indices:
275-
con = contour_obj.collections[icon]
276-
paths = con.get_paths()
282+
#for icon in indices:
283+
#con = contour_obj.collections[icon]
284+
#paths = con.get_paths()
277285

278-
for segNum, linepath in enumerate(paths):
279-
lc = linepath.vertices
280-
d = _find_closest_point_on_path(lc, point)
286+
#for segNum, linepath in enumerate(paths):
287+
#lc = linepath.vertices
288+
##d = _find_closest_point_on_path(lc, point)
289+
#d = fcpl.find_closest_point_on_path(lc, point)
281290

282-
if d < dmin:
283-
dmin = d
284-
#conmin = icon
285-
segmin = segNum
291+
#if d < dmin:
292+
#dmin = d
293+
##conmin = icon
294+
#segmin = segNum
295+
296+
#return segmin
297+
298+
286299

287-
return segmin
300+
def find_nearest_contour(contcoll, x, y):
301+
"""
302+
Finds contour that is closest to a point.
288303
304+
Returns a tuple containing the contour & segment.
289305
306+
Call signature::
290307
308+
segmin = find_nearest_contour(contcoll, x, y)
291309
310+
"""
311+
dmin = 1e10
312+
segmin = None
313+
linepathmin = None
314+
315+
paths = contcoll.get_paths()
316+
for segNum, linepath in enumerate(paths):
317+
lc = linepath.vertices
318+
ds = lc[:,0] - x
319+
ds **= 2
320+
dss = lc[:,1] - y
321+
dss **= 2
322+
ds += dss
323+
d = ds.min()
324+
if d < dmin:
325+
dmin = d
326+
segmin = segNum
327+
linepathmin = linepath
328+
return (segmin, linepathmin)
292329

293330

294331

0 commit comments

Comments
 (0)