Merge pull request AntSimi#66 from CoriPegliasco/master

fancy figures in ACC detection

Author: AntSimi

examples/02_eddy_identification/pet_eddy_detection_ACC.py

@@ -1,35 +1,58 @@
 """
-Eddy detection : Antartic circum polar
-======================================
+Eddy detection : Antartic Circumpolar Current
+=============================================
 
-Script will detect eddies on adt field, and compute u,v with method add_uv(which could use, only if equator is avoid)
+This script detect eddies on the ADT field, and compute u,v with the method add_uv (use it only if the Equator is avoided)
 
-Two ones with filtering adt and another without
+Two detections are provided : with a filtered ADT and without filtering
 
 """
 from datetime import datetime
 
 from matplotlib import pyplot as plt
+from matplotlib import style
 
 from py_eddy_tracker import data
 from py_eddy_tracker.dataset.grid import RegularGridDataset
 
+pos_cb = [0.1, 0.52, 0.83, 0.015]
+pos_cb2 = [0.1, 0.07, 0.4, 0.015]
+
 
 def quad_axes(title):
-    fig = plt.figure(figsize=(13, 8.5))
-    fig.suptitle(title, weight="bold")
+    style.use("default")
+    fig = plt.figure(figsize=(13, 10))
+    fig.suptitle(title, weight="bold", fontsize=14)
     axes = list()
-    for position in (
-        [0.05, 0.53, 0.44, 0.44],
-        [0.53, 0.53, 0.44, 0.44],
-        [0.05, 0.03, 0.44, 0.44],
-        [0.53, 0.03, 0.44, 0.44],
-    ):
+
+    ax_pos = dict(
+        topleft=[0.1, 0.54, 0.4, 0.38],
+        topright=[0.53, 0.54, 0.4, 0.38],
+        botleft=[0.1, 0.09, 0.4, 0.38],
+        botright=[0.53, 0.09, 0.4, 0.38],
+    )
+
+    for key, position in ax_pos.items():
         ax = fig.add_axes(position)
         ax.set_xlim(5, 45), ax.set_ylim(-60, -37)
         ax.set_aspect("equal"), ax.grid(True)
         axes.append(ax)
-    return axes
+        if "right" in key:
+            ax.set_yticklabels("")
+    return fig, axes
+
+
+def set_fancy_labels(fig, ticklabelsize=14, labelsize=14, labelweight="semibold"):
+    for ax in fig.get_axes():
+        ax.grid()
+        ax.grid(which="major", linestyle="-", linewidth="0.5", color="black")
+        if ax.get_ylabel() != "":
+            ax.set_ylabel(ax.get_ylabel(), fontsize=labelsize, fontweight=labelweight)
+        if ax.get_xlabel() != "":
+            ax.set_xlabel(ax.get_xlabel(), fontsize=labelsize, fontweight=labelweight)
+        if ax.get_title() != "":
+            ax.set_title(ax.get_title(), fontsize=labelsize, fontweight=labelweight)
+        ax.tick_params(labelsize=ticklabelsize)
 
 
 # %%
@@ -69,20 +92,25 @@ def quad_axes(title):
 # %%
 # Figures
 # -------
-axs = quad_axes("General properties field")
-m = g_raw.display(axs[0], "adt", vmin=-1, vmax=1, cmap="RdBu_r")
-axs[0].set_title("ADT(m)")
-m = g.display(axs[1], "adt_low", vmin=-1, vmax=1, cmap="RdBu_r")
-axs[1].set_title("ADT (m) large scale with cut at 700 km")
-m = g.display(axs[2], "adt", vmin=-1, vmax=1, cmap="RdBu_r")
-axs[2].set_title("ADT (m) high scale with cut at 700 km")
-cb = plt.colorbar(
-    m, cax=axs[0].figure.add_axes([0.03, 0.51, 0.94, 0.01]), orientation="horizontal"
-)
-cb.set_label("ADT(m)", labelpad=-2)
+kw_adt = dict(vmin=-1.5, vmax=1.5, cmap=plt.get_cmap("RdBu_r", 30))
+fig, axs = quad_axes("General properties field")
+m = g_raw.display(axs[0], "adt", **kw_adt)
+axs[0].set_title("Total ADT (m)")
+m = g.display(axs[1], "adt_low", **kw_adt)
+axs[1].set_title("ADT (m) large scale, cutoff at 700 km")
+m2 = g.display(axs[2], "adt", cmap=plt.get_cmap("RdBu_r", 20), vmin=-0.5, vmax=0.5)
+axs[2].set_title("ADT (m) high-pass filtered, a cutoff at 700 km")
+cb = plt.colorbar(m, cax=axs[0].figure.add_axes(pos_cb), orientation="horizontal")
+cb.set_label("ADT (m)", labelpad=0)
+cb2 = plt.colorbar(m2, cax=axs[2].figure.add_axes(pos_cb2), orientation="horizontal")
+cb2.set_label("ADT (m)", labelpad=0)
+set_fancy_labels(fig)
+
+# %%
+# The large-scale North-South gradient is removed by the filtering step.
 
 # %%
-axs = quad_axes("")
+fig, axs = quad_axes("")
 axs[0].set_title("Without filter")
 axs[0].set_ylabel("Contours used in eddies")
 axs[1].set_title("With filter")
@@ -91,13 +119,18 @@ def quad_axes(title):
 g.contours.display(axs[1], lw=0.5, only_used=True)
 g_raw.contours.display(axs[2], lw=0.5, only_unused=True)
 g.contours.display(axs[3], lw=0.5, only_unused=True)
+set_fancy_labels(fig)
+
+# %%
+# Removing the large-scale North-South gradient reveals closed contours in the
+# South-Western corner of the ewample region.
 
 # %%
 kw = dict(ref=-10, linewidth=0.75)
 kw_a = dict(color="r", label="Anticyclonic ({nb_obs} eddies)")
 kw_c = dict(color="b", label="Cyclonic ({nb_obs} eddies)")
 kw_filled = dict(vmin=0, vmax=100, cmap="Spectral_r", lut=20, intern=True, factor=100)
-axs = quad_axes("Comparison between two detection")
+fig, axs = quad_axes("Comparison between two detections")
 # Match with intern/inner contour
 i_a, j_a, s_a = a_.match(a, intern=True, cmin=0.15)
 i_c, j_c, s_c = c_.match(c, intern=True, cmin=0.15)
@@ -107,10 +140,8 @@ def quad_axes(title):
 c_.index(i_c).filled(axs[0], s_c, **kw_filled)
 m = c.index(j_c).filled(axs[1], s_c, **kw_filled)
 
-cb = plt.colorbar(
-    m, cax=axs[0].figure.add_axes([0.03, 0.51, 0.94, 0.01]), orientation="horizontal"
-)
-cb.set_label("Similarity index", labelpad=-5)
+cb = plt.colorbar(m, cax=axs[0].figure.add_axes(pos_cb), orientation="horizontal")
+cb.set_label("Similarity index (%)", labelpad=-5)
 a_.display(axs[0], **kw, **kw_a), c_.display(axs[0], **kw, **kw_c)
 a.display(axs[1], **kw, **kw_a), c.display(axs[1], **kw, **kw_c)
 
@@ -124,6 +155,12 @@ def quad_axes(title):
 
 for ax in axs:
     ax.legend()
+
+set_fancy_labels(fig)
+
+# %%
+# Very similar eddies have Similarity Indexes >= 40%
+
 # %%
 # Criteria for rejecting a contour :
 #  0. Accepted (green)
@@ -140,10 +177,10 @@ def quad_axes(title):
 
 for i, (label, field, factor, stop) in enumerate(
     (
-        ("speed radius (km)", "radius_s", 0.001, 120),
-        ("outter radius (km)", "radius_e", 0.001, 120),
-        ("amplitude (cm)", "amplitude", 100, 25),
-        ("speed max (cm/s)", "speed_average", 100, 25),
+        ("Speed radius (km)", "radius_s", 0.001, 120),
+        ("Effective radius (km)", "radius_e", 0.001, 120),
+        ("Amplitude (cm)", "amplitude", 100, 25),
+        ("Speed max (cm/s)", "speed_average", 100, 25),
     )
 ):
     ax = fig.add_subplot(2, 2, i + 1, title=label)
@@ -166,3 +203,5 @@ def quad_axes(title):
     ax.plot((0, 1000), (0, 1000), "g")
     ax.set_xlim(0, stop), ax.set_ylim(0, stop)
     ax.legend()
+
+set_fancy_labels(fig)