- add ellipse fit

- minor english modifs

Author: CoriPegliasco

examples/16_network/pet_ioannou_2017_case.py

@@ -188,3 +188,47 @@ def update_axes(ax, mappable=None):
 m = close_to_i3.scatter_timeline(ax, "shape_error_e", vmin=14, vmax=70, **kw)
 cb = update_axes(ax, m["scatter"])
 cb.set_label("Effective shape error")
+
+# %%
+# Rotation angle
+# --------------
+from py_eddy_tracker.generic import coordinates_to_local, local_to_coordinates
+from py_eddy_tracker.poly import fit_ellips
+theta_ = list()
+a_ = list()
+b_ = list()
+for obs in close_to_i3:
+    x, y = obs['contour_lon_s'], obs['contour_lat_s']
+    x0_, y0_ = x.mean(), y.mean()
+    x_, y_ = coordinates_to_local(x, y, x0_, y0_)
+    x0, y0, a, b, theta = fit_ellips(x_, y_)
+    theta_.append(theta)
+    a_.append(a)
+    b_.append(b)
+a_=array(a_)
+b_=array(b_)
+
+# %%
+# Theta
+ax = timeline_axes()
+m = close_to_i3.scatter_timeline(ax, theta_, vmin=-pi/2, vmax=pi/2, cmap='hsv')
+cb = update_axes(ax, m["scatter"])
+
+# %%
+# A
+ax = timeline_axes()
+m = close_to_i3.scatter_timeline(ax, a_ * 1e-3, vmin=0, vmax=80, cmap='Spectral_r')
+cb = update_axes(ax, m["scatter"])
+
+# %%
+# B
+ax = timeline_axes()
+m = close_to_i3.scatter_timeline(ax, b_ * 1e-3, vmin=0, vmax=80, cmap='Spectral_r')
+cb = update_axes(ax, m["scatter"])
+
+# %%
+# A/B
+ax = timeline_axes()
+m = close_to_i3.scatter_timeline(ax, a_/b_, vmin=1, vmax=2, cmap='Spectral_r')
+cb = update_axes(ax, m["scatter"])
+

src/py_eddy_tracker/dataset/grid.py

@@ -609,17 +609,17 @@ def eddy_identification(
         :param str vname: Grid name of v speed component
         :param datetime.datetime date: Date which will be stored in object to date data
         :param float,int step: Height between two layers in m
-        :param float,int shape_error: Maximal error allowed for outter contour in %
+        :param float,int shape_error: Maximal error allowed for outermost contour in %
         :param int sampling: Number of points to store contours and speed profile
         :param str sampling_method: Method to resample 'uniform' or 'visvalingam'
         :param (int,int),None pixel_limit:
-            Min and max number of pixels inside the inner and the outer contour to be considered as an eddy
+            Min and max number of pixels inside the inner and the outermost contour to be considered as an eddy
         :param float,None precision: Truncate values at the defined precision in m
         :param str force_height_unit: Unit used for height unit
         :param str force_speed_unit: Unit used for speed unit
         :param dict kwargs: Argument given to amplitude
 
-        :return: Return a list of 2 elements: Anticyclone and Cyclone
+        :return: Return a list of 2 elements: Anticyclones and Cyclones
         :rtype: py_eddy_tracker.observations.observation.EddiesObservations
 
         .. minigallery:: py_eddy_tracker.GridDataset.eddy_identification
@@ -729,7 +729,8 @@ def eddy_identification(
                 for contour in contour_paths:
                     if contour.used:
                         continue
-                    # FIXME : center could be not in contour and fit on raw sampling
+                    # FIXME : center could be outside the contour due to the fit
+                    # FIXME : warning : the fit is made on raw sampling
                     _, _, _, aerr = contour.fit_circle()
 
                     # Filter for shape
@@ -752,6 +753,7 @@ def eddy_identification(
                     ):
                         continue
 
+                    # Test the number of pixels within the outermost contour
                     # FIXME : Maybe limit max must be replace with a maximum of surface
                     if (
                         contour.nb_pixel < pixel_limit[0]
@@ -760,6 +762,9 @@ def eddy_identification(
                         contour.reject = 3
                         continue
 
+                    # Here the considered contour passed shape_error test, masked_pixels test,
+                    # values strictly above (AEs) or below (CEs) the contour, number_pixels test)
+
                     # Compute amplitude
                     reset_centroid, amp = self.get_amplitude(
                         contour,

src/py_eddy_tracker/eddy_feature.py

@@ -31,7 +31,7 @@
 class Amplitude(object):
     """
     Class to calculate *amplitude* and counts of *local maxima/minima*
-    within a closed region of a sea level anomaly field.
+    within a closed region of a sea surface height field.
     """
 
     EPSILON = 1e-8
@@ -66,8 +66,8 @@ def __init__(
         :param array data:
         :param float interval:
         :param int mle: maximum number of local maxima in contour
-        :param int nb_step_min: number of interval to consider like an eddy
-        :param int nb_step_to_be_mle: number of interval to be consider like another maxima
+        :param int nb_step_min: number of intervals to consider an eddy
+        :param int nb_step_to_be_mle: number of intervals to be considered as an another maxima
         """
 
         # Height of the contour
@@ -102,8 +102,9 @@ def __init__(
         self.nb_pixel = i_x.shape[0]
 
         # Only pixel in contour
+        # FIXME : change sla by ssh as the grid can be adt?
         self.sla = data[contour.pixels_index]
-        # Amplitude which will be provide
+        # Amplitude which will be provided
         self.amplitude = 0
         # Maximum local extrema accepted
         self.mle = mle
@@ -115,9 +116,10 @@ def within_amplitude_limits(self):
     def all_pixels_below_h0(self, level):
         """
         Check CSS11 criterion 1: The SSH values of all of the pixels
-        are below a given SSH threshold for cyclonic eddies.
+        are below (above) a given SSH threshold for cyclonic (anticyclonic)
+        eddies.
         """
-        # In some case pixel value must be very near of contour bounds
+        # In some cases pixel value may be very close to the contour bounds
         if self.sla.mask.any() or ((self.sla.data - self.h_0) > self.EPSILON).any():
             return False
         else:
@@ -293,10 +295,10 @@ class Contours(object):
 
     Attributes:
       contour:
-        A matplotlib contour object of high-pass filtered SLA
+        A matplotlib contour object of high-pass filtered SSH
 
       eddy:
-        A tracklist object holding the SLA data
+        A tracklist object holding the SSH data
 
       grd:
         A grid object
@@ -406,7 +408,7 @@ def __init__(self, x, y, z, levels, wrap_x=False, keep_unclose=False):
         fig = Figure()
         ax = fig.add_subplot(111)
         if wrap_x:
-            logger.debug("wrapping activate to compute contour")
+            logger.debug("wrapping activated to compute contour")
             x = concatenate((x, x[:1] + 360))
             z = ma.concatenate((z, z[:1]))
         logger.debug("X shape : %s", x.shape)
@@ -602,8 +604,8 @@ def display(
             Must be 'shape_error', 'x', 'y' or 'radius'.
             If define display_criterion is not use.
             bins argument must be define
-        :param array bins: bins use to colorize contour
-        :param str cmap: Name of cmap to use for field display
+        :param array bins: bins used to colorize contour
+        :param str cmap: Name of cmap for field display
         :param dict kwargs: look at :py:meth:`matplotlib.collections.LineCollection`
 
         .. minigallery:: py_eddy_tracker.Contours.display
@@ -688,7 +690,7 @@ def display(
             ax.autoscale_view()
 
     def label_contour_unused_which_contain_eddies(self, eddies):
-        """Select contour which contain several eddies"""
+        """Select contour containing several eddies"""
         if eddies.sign_type == 1:
             # anticyclonic
             sl = slice(None, -1)

src/py_eddy_tracker/poly.py

@@ -474,13 +474,13 @@ def polygon_overlap(p0, p1, minimal_area=False):
     return cost
 
 
-# FIXME: only one function are needed
+# FIXME: only one function is needed
 @njit(cache=True)
 def fit_circle(x, y):
     """
     From a polygon, function will fit a circle.
 
-    Must be call with local coordinates (in m, to get a radius in m).
+    Must be called with local coordinates (in m, to get a radius in m).
 
     :param array x: x of polygon
     :param array y: y of polygon
@@ -510,11 +510,11 @@ def fit_circle(x, y):
     radius **= 0.5
     x0 *= scale
     y0 *= scale
-    # radius of fitted circle
+    # radius of fit circle
     radius *= scale
-    # center X-position of fitted circle
+    # center X-position of fit circle
     x0 += x_mean
-    # center Y-position of fitted circle
+    # center Y-position of fit circle
     y0 += y_mean
 
     err = shape_error(x, y, x0, y0, radius)