Skip to content
Closed
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
Correction of global indexation, Add comments
Merge branch 'master' of https://github.com/AntSimi/py-eddy-tracker

# Conflicts:
#	src/py_eddy_tracker/observations/observation.py
  • Loading branch information
CoriPegliasco committed Jan 13, 2021
commit f0070e9c6999691e32e2ccce4cc065e95e4fbc68
4 changes: 2 additions & 2 deletions src/py_eddy_tracker/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -197,7 +197,7 @@ def parse_args(self, *args, **kwargs):
nc_attr=dict(
units="degrees_east",
axis="X",
comment="Longitude center of the fitted circle",
comment="Longitude center of the fit circle",
long_name="Eddy Center Longitude",
standard_name="longitude",
),
Expand All @@ -214,7 +214,7 @@ def parse_args(self, *args, **kwargs):
axis="Y",
long_name="Eddy Center Latitude",
standard_name="latitude",
comment="Latitude center of the fitted circle",
comment="Latitude center of the fit circle",
),
),
lon_max=dict(
Expand Down
92 changes: 53 additions & 39 deletions src/py_eddy_tracker/observations/tracking.py
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
# -*- coding: utf-8 -*-
"""
Class to manage observations gathered in track
Class to manage observations gathered in trajectories
"""
import logging
from datetime import datetime, timedelta
Expand Down Expand Up @@ -124,7 +124,7 @@ def __repr__(self):
return content

def add_distance(self):
"""Add a field of distance (m) between to consecutive observation, 0 for the last observation of each track"""
"""Add a field of distance (m) between two consecutive observations, 0 for the last observation of each track"""
if "distance_next" in self.observations.dtype.descr:
return self
new = self.add_fields(("distance_next",))
Expand Down Expand Up @@ -182,7 +182,7 @@ def filled_by_interpolation(self, mask):
)

def extract_longer_eddies(self, nb_min, nb_obs, compress_id=True):
"""Select eddies which are longer than nb_min"""
"""Select the trajectories longer than nb_min"""
mask = nb_obs >= nb_min
nb_obs_select = mask.sum()
logger.info("Selection of %d observations", nb_obs_select)
Expand Down Expand Up @@ -227,9 +227,9 @@ def set_global_attr_netcdf(self, h_nc):

def extract_with_period(self, period, **kwargs):
"""
Extract with a period
Extract within a time period

:param (int,int) period: two date to define period, must be specify from 1/1/1950
:param (int,int) period: two dates to define the period, must be specify from 1/1/1950
:param dict kwargs: look at :py:meth:`extract_with_mask`
:return: Return all eddy tracks which are in bounds
:rtype: TrackEddiesObservations
Expand All @@ -252,9 +252,9 @@ def extract_with_period(self, period, **kwargs):

def get_azimuth(self, equatorward=False):
"""
Return azimuth for each tracks.
Return azimuth for each track.

Azimuth is compute with first and last observation
Azimuth is computed with first and last observation

:param bool equatorward: If True, Poleward are positive and equatorward negative
:rtype: array
Expand Down Expand Up @@ -286,7 +286,7 @@ def compute_index(self):
"""
if self.__first_index_of_track is None:
s = self.tracks.max() + 1
# Doesn't work => core dump with numba, maybe he wait i8 instead of u4
# Doesn't work => core dump with numba, maybe he wants i8 instead of u4
# self.__first_index_of_track = -ones(s, self.tracks.dtype)
# self.__obs_by_track = zeros(s, self.observation_number.dtype)
self.__first_index_of_track = -ones(s, "i8")
Expand Down Expand Up @@ -334,12 +334,12 @@ def nb_obs_by_track(self):

@property
def lifetime(self):
"""Return for each observation lifetime"""
"""Return lifetime for each observation"""
return self.nb_obs_by_track.repeat(self.nb_obs_by_track)

@property
def age(self):
"""Return for each observation age in %, will be [0:100]"""
"""Return age in % for each observation, will be [0:100]"""
return self.n.astype("f4") / (self.lifetime - 1) * 100.0

def extract_ids(self, tracks):
Expand All @@ -348,10 +348,10 @@ def extract_ids(self, tracks):

def extract_toward_direction(self, west=True, delta_lon=None):
"""
Get eddy which go in same direction
Get trajectories going in the same direction

:param bool west: Only eastward eddy if True return westward
:param None,float delta_lon: Only eddy with more than delta_lon span in longitude
:param bool west: Only eastward eddies if True return westward
:param None,float delta_lon: Only eddies with more than delta_lon span in longitude
:return: Only eastern eddy
:rtype: __class__

Expand Down Expand Up @@ -398,10 +398,10 @@ def extract_in_direction(self, direction, value=0):

def extract_with_length(self, bounds):
"""
Return all observations in [b0:b1]
Return the observations within trajectories lasting between [b0:b1]

:param (int,int) bounds: length min and max of selected eddies, if use of -1 this bound is not used
:return: Return all eddy tracks which have length between bounds
:param (int,int) bounds: length min and max of the desired trajectories, if -1 this bound is not used
:return: Return all trajectories having length between bounds
:rtype: TrackEddiesObservations

.. minigallery:: py_eddy_tracker.TrackEddiesObservations.extract_with_length
Expand Down Expand Up @@ -461,11 +461,11 @@ def extract_with_mask(
Extract a subset of observations

:param array(bool) mask: mask to select observations
:param bool full_path: extract full path if only one part is selected
:param bool remove_incomplete: delete path which are not fully selected
:param bool compress_id: resample track number to use a little range
:param bool reject_virtual: if track are only virtual in selection we remove track
:return: same object with selected observations
:param bool full_path: extract the full trajectory if only one part is selected
:param bool remove_incomplete: delete trajectory if not fully selected
:param bool compress_id: resample trajectory number to use a smaller range
:param bool reject_virtual: if only virtual are selected, the trajectory is removed
:return: same object with the selected observations
:rtype: self.__class__
"""
if full_path and remove_incomplete:
Expand Down Expand Up @@ -510,7 +510,9 @@ def re_reference_index(index, ref):

def shape_polygon(self, intern=False):
"""
Get polygons which enclosed each track
Get the polygon enclosing each trajectory.

The polygon merges the non-overlapping bounds of the specified contours

:param bool intern: If True use speed contour instead of effective contour
:rtype: list(array, array)
Expand All @@ -520,9 +522,9 @@ def shape_polygon(self, intern=False):

def display_shape(self, ax, ref=None, intern=False, **kwargs):
"""
This function will draw shape of each track
This function will draw the shape of each trajectory

:param matplotlib.axes.Axes ax: ax where drawed
:param matplotlib.axes.Axes ax: ax to draw
:param float,int ref: if defined all coordinates will be wrapped with ref like west boundary
:param bool intern: If True use speed contour instead of effective contour
:param dict kwargs: keyword arguments for Axes.plot
Expand All @@ -547,10 +549,10 @@ def display_shape(self, ax, ref=None, intern=False, **kwargs):

def close_tracks(self, other, nb_obs_min=10, **kwargs):
"""
Get close from another atlas.
Get close trajectories from another atlas.

:param self other: Atlas to compare
:param int nb_obs_min: Minimal number of overlap for one track
:param int nb_obs_min: Minimal number of overlap for one trajectory
:param dict kwargs: keyword arguments for match function
:return: return other atlas reduce to common track with self

Expand All @@ -577,10 +579,10 @@ def format_label(self, label):

def plot(self, ax, ref=None, **kwargs):
"""
This function will draw path of each track
This function will draw path of each trajectory

:param matplotlib.axes.Axes ax: ax where drawed
:param float,int ref: if defined all coordinates will be wrapped with ref like west boundary
:param matplotlib.axes.Axes ax: ax to draw
:param float,int ref: if defined, all coordinates will be wrapped with ref like west boundary
:param dict kwargs: keyword arguments for Axes.plot
:return: matplotlib mappable
"""
Expand All @@ -595,7 +597,7 @@ def plot(self, ax, ref=None, **kwargs):
return ax.plot(x, y, **kwargs)

def split_network(self, intern=True, **kwargs):
"""Divide each group in track"""
"""Return each group (network) divided in segments"""
track_s, track_e, track_ref = build_index(self.tracks)
ids = empty(
len(self),
Expand All @@ -610,26 +612,32 @@ def split_network(self, intern=True, **kwargs):
],
)
ids["group"], ids["time"] = self.tracks, self.time
# To store id track
# Initialisation
# To store the id of the segments, the backward and forward cost associations
ids["track"], ids["previous_cost"], ids["next_cost"] = 0, 0, 0
# To store the indexes of the backward and forward observations associated
ids["previous_obs"], ids["next_obs"] = -1, -1
# At the end, ids["previous_obs"] == -1 means the start of a non-split segment
# and ids["next_obs"] == -1 means the end of a non-merged segment

xname, yname = self.intern(intern)
for i_s, i_e in zip(track_s, track_e):
if i_s == i_e or self.tracks[i_s] == self.NOGROUP:
continue
sl = slice(i_s, i_e)
local_ids = ids[sl]
# built segments with local indices
self.set_tracks(self[xname][sl], self[yname][sl], local_ids, **kwargs)
m = local_ids["previous_obs"] == -1
# shift the local indices to the total indexation for the used observations
m = local_ids["previous_obs"] != -1
local_ids["previous_obs"][m] += i_s
m = local_ids["next_obs"] == -1
m = local_ids["next_obs"] != -1
local_ids["next_obs"][m] += i_s
return ids

def set_tracks(self, x, y, ids, window):
"""
Will split one group in tracks
Will split one group (network) in segments

:param array x: coordinates of group
:param array y: coordinates of group
Expand All @@ -641,19 +649,22 @@ def set_tracks(self, x, y, ids, window):
nb = x.shape[0]
used = zeros(nb, dtype="bool")
track_id = 1
# build all polygon (need to check if wrap is needed)
# build all polygons (need to check if wrap is needed)
polygons = [Polygon(create_vertice_from_2darray(x, y, i)) for i in range(nb)]
for i in range(nb):
# If observation already in one track, we go to the next one
# If the observation is already in one track, we go to the next one
if used[i]:
continue
# Search a possible continuation (forward)
self.follow_obs(i, track_id, used, ids, polygons, *time_index, window)
track_id += 1
# Search a possible ancestor
# Search a possible ancestor (backward)
self.previous_obs(i, ids, polygons, *time_index, window)

@classmethod
def follow_obs(cls, i_next, track_id, used, ids, *args):
"""Associate the observations to the segments"""

while i_next != -1:
# Flag
used[i_next] = True
Expand All @@ -677,6 +688,8 @@ def follow_obs(cls, i_next, track_id, used, ids, *args):

@staticmethod
def previous_obs(i_current, ids, polygons, time_s, time_e, time_ref, window):
"""Backward association of observations to the segments"""

time_cur = ids["time"][i_current]
t0, t1 = time_cur - 1 - time_ref, max(time_cur - window - time_ref, 0)
for t_step in range(t0, t1 - 1, -1):
Expand All @@ -700,6 +713,7 @@ def previous_obs(i_current, ids, polygons, time_s, time_e, time_ref, window):

@staticmethod
def next_obs(i_current, ids, polygons, time_s, time_e, time_ref, window):
"""Forward association of observations to the segments"""
time_max = time_e.shape[0] - 1
time_cur = ids["time"][i_current]
t0, t1 = time_cur + 1 - time_ref, min(time_cur + window - time_ref, time_max)
Expand All @@ -710,9 +724,9 @@ def next_obs(i_current, ids, polygons, time_s, time_e, time_ref, window):
# No observation at the time step
if i0 == i1:
continue
# Intersection / union, to be able to separte in case of multiple inside
# Intersection / union, to be able to separate in case of multiple obs inside the current polygon in next_step
c = polygon_overlap(polygons[i_current], polygons[i0:i1], minimal_area=True)
# We remove low overlap
# We remove low overlap TO DO : add the value as a parameter (default 0.1)?
c[c < 0.1] = 0
# We get index of maximal overlap
i = c.argmax()
Expand Down
6 changes: 3 additions & 3 deletions src/py_eddy_tracker/poly.py
Original file line number Diff line number Diff line change
Expand Up @@ -452,7 +452,7 @@ def polygon_overlap(p0, p1, minimal_area=False):

:param list(Polygon) p0: List of polygon to compare with p1 list
:param list(Polygon) p1: List of polygon to compare with p0 list
:param bool minimal_area: If True, function will compute intersection/little polygon, else intersection/union
:param bool minimal_area: If True, function will compute intersection/smaller polygon, else intersection/union
:return: Result of cost function
:rtype: array
"""
Expand All @@ -462,10 +462,10 @@ def polygon_overlap(p0, p1, minimal_area=False):
p_ = p1[i]
# Area of intersection
intersection = (p0 & p_).area()
# we divide intersection with the little one result from 0 to 1
# we divide the intersection by the smaller area, result from 0 to 1
if minimal_area:
cost[i] = intersection / min(p0.area(), p_.area())
# we divide intersection with polygon merging result from 0 to 1
# we divide the intersection by the merged polygons area, result from 0 to 1
else:
cost[i] = intersection / (p0 + p_).area()
return cost
Expand Down