Add method to list merging/spliting event from network

Author: AntSimi

examples/16_network/pet_relative.py

@@ -15,6 +15,11 @@
 # ---------
 # Load data where observations are put in same network but no segmentation
 e = TrackEddiesObservations.load_file(data.get_path("c568803.nc"))
+# FIXME : Must be rewrote
+e.lon[:] = (e.lon + 180) % 360 - 180
+e.contour_lon_e[:] = ((e.contour_lon_e.T - e.lon + 180) % 360 - 180 + e.lon).T
+e.contour_lon_s[:] = ((e.contour_lon_s.T - e.lon + 180) % 360 - 180 + e.lon).T
+##############
 n = NetworkObservations.from_split_network(e, e.split_network(intern=False, window=5))
 
 # %%
@@ -113,18 +118,46 @@
 # %%
 fig = plt.figure(figsize=(15, 5))
 ax = fig.add_axes([0.04, 0.06, 0.90, 0.88])
-close_to_i = n.relative(i, order=1)
-ax.set_title(f"Close segments ({close_to_i.infos()})")
-_ = close_to_i.display_timeline(ax)
+close_to_i1 = n.relative(i, order=1)
+ax.set_title(f"Close segments ({close_to_i1.infos()})")
+_ = close_to_i1.display_timeline(ax)
 # %%
 fig = plt.figure(figsize=(15, 5))
 ax = fig.add_axes([0.04, 0.06, 0.90, 0.88])
-close_to_i = n.relative(i, order=2)
-ax.set_title(f"Close segments ({close_to_i.infos()})")
-_ = close_to_i.display_timeline(ax)
+close_to_i2 = n.relative(i, order=2)
+ax.set_title(f"Close segments ({close_to_i2.infos()})")
+_ = close_to_i2.display_timeline(ax)
 # %%
 fig = plt.figure(figsize=(15, 5))
 ax = fig.add_axes([0.04, 0.06, 0.90, 0.88])
-close_to_i = n.relative(i, order=3)
-ax.set_title(f"Close segments ({close_to_i.infos()})")
-_ = close_to_i.display_timeline(ax)
+close_to_i3 = n.relative(i, order=3)
+ax.set_title(f"Close segments ({close_to_i3.infos()})")
+_ = close_to_i3.display_timeline(ax)
+
+# %%
+# Display track on map
+# --------------------
+fig = plt.figure(figsize=(15, 8))
+ax = fig.add_axes([0.04, 0.06, 0.90, 0.88], projection="full_axes")
+close_to_i2.plot(ax)
+_ = ax.set_xlim(-10, 20), ax.set_ylim(-37, -21), ax.grid()
+
+# %%
+# Get merging event
+# -----------------
+fig = plt.figure(figsize=(15, 8))
+ax = fig.add_axes([0.04, 0.06, 0.90, 0.88], projection="full_axes")
+merging = close_to_i2.merging_event()
+merging.display(ax)
+ax.set_xlim(-10, 20), ax.set_ylim(-37, -21), ax.grid()
+merging
+
+# %%
+# Get spliting event
+# ------------------
+fig = plt.figure(figsize=(15, 8))
+ax = fig.add_axes([0.04, 0.06, 0.90, 0.88], projection="full_axes")
+spliting = close_to_i2.spliting_event()
+spliting.display(ax)
+ax.set_xlim(-10, 20), ax.set_ylim(-37, -21), ax.grid()
+spliting

src/py_eddy_tracker/observations/network.py

@@ -8,7 +8,7 @@
 from numba import njit
 from numpy import arange, array, bincount, empty, ones, uint32, unique
 
-from ..generic import build_index
+from ..generic import build_index, wrap_longitude
 from ..poly import bbox_intersection, vertice_overlap
 from .observation import EddiesObservations
 from .tracking import TrackEddiesObservations
@@ -232,16 +232,84 @@ def scatter_timeline(self, ax, name, factor=1, event=True, **kwargs):
         return mappables
 
     def insert_virtual(self):
+        # TODO
         pass
 
     def merging_event(self):
-        pass
+        indices = list()
+        for i, b0, b1 in self.iter_on("segment"):
+            nb = i.stop - i.start
+            if nb == 0:
+                continue
+            i_n = self.next_obs[i.stop - 1]
+            if i_n != -1:
+                indices.append(i.stop - 1)
+        indices = list(set(indices))
+        nb = len(indices)
+        new = EddiesObservations(
+            nb,
+            track_extra_variables=self.track_extra_variables,
+            track_array_variables=self.track_array_variables,
+            array_variables=self.array_variables,
+            only_variables=self.only_variables,
+            raw_data=self.raw_data,
+        )
+
+        for k in new.obs.dtype.names:
+            new[k][:] = self[k][indices]
+        new.sign_type = self.sign_type
+        return new
 
     def spliting_event(self):
-        pass
+        indices = list()
+        for i, b0, b1 in self.iter_on("segment"):
+            nb = i.stop - i.start
+            if nb == 0:
+                continue
+            i_p = self.previous_obs[i.start]
+            if i_p != -1:
+                indices.append(i.start)
+        indices = list(set(indices))
+        nb = len(indices)
+        new = EddiesObservations(
+            nb,
+            track_extra_variables=self.track_extra_variables,
+            track_array_variables=self.track_array_variables,
+            array_variables=self.array_variables,
+            only_variables=self.only_variables,
+            raw_data=self.raw_data,
+        )
+
+        for k in new.obs.dtype.names:
+            new[k][:] = self[k][indices]
+        new.sign_type = self.sign_type
+        return new
 
     def fully_connected(self):
         self.only_one_network()
+        # TODO
+
+    def plot(self, ax, ref=None, **kwargs):
+        """
+        This function will draw path of each trajectory
+
+        :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: a list of matplotlib mappables
+        """
+        mappables = list()
+        if "label" in kwargs:
+            kwargs["label"] = self.format_label(kwargs["label"])
+        for i, b0, b1 in self.iter_on("segment"):
+            nb = i.stop - i.start
+            if nb == 0:
+                continue
+            x, y = self.lon[i], self.lat[i]
+            if ref is not None:
+                x, y = wrap_longitude(x, y, ref, cut=True)
+            mappables.append(ax.plot(x, y, **kwargs)[0])
+        return mappables
 
     def remove_dead_branch(self, nobs=3):
         """"""