diff --git a/examples/02_eddy_identification/pet_contour_circle.py b/examples/02_eddy_identification/pet_contour_circle.py index 5cbba26b..45f1554d 100644 --- a/examples/02_eddy_identification/pet_contour_circle.py +++ b/examples/02_eddy_identification/pet_contour_circle.py @@ -13,7 +13,7 @@ a = EddiesObservations.load_file(data.get_path("Anticyclonic_20190223.nc")) # %% -# Plot +# Plot the speed and effective (dashed) contours fig = plt.figure(figsize=(15, 8)) ax = fig.add_axes((0.05, 0.05, 0.9, 0.9)) ax.set_aspect("equal") @@ -21,7 +21,7 @@ ax.set_ylim(-50, -25) a.display(ax, label="Anticyclonic contour", color="r", lw=1) -# Replace contour by circle +# Replace contours by circles using center and radius (effective is dashed) a.circle_contour() a.display(ax, label="Anticyclonic circle", color="g", lw=1) ax.legend(loc="upper right") diff --git a/examples/02_eddy_identification/pet_display_id.py b/examples/02_eddy_identification/pet_display_id.py index a2c2f84d..5d061c24 100644 --- a/examples/02_eddy_identification/pet_display_id.py +++ b/examples/02_eddy_identification/pet_display_id.py @@ -14,7 +14,7 @@ c = EddiesObservations.load_file(data.get_path("Cyclonic_20190223.nc")) # %% -# filled contour with amplitude field +# Fill effective contour with amplitude fig = plt.figure(figsize=(15, 8)) ax = fig.add_axes([0.03, 0.03, 0.90, 0.94]) ax.set_aspect("equal") @@ -24,11 +24,11 @@ a.display(ax, **kwargs), c.display(ax, **kwargs) a.filled(ax, "amplitude", cmap="magma_r", vmin=0, vmax=.5) m = c.filled(ax, "amplitude", cmap="magma_r", vmin=0, vmax=.5) -colorbar = plt.colorbar(m, cax=ax.figure.add_axes([0.95, 0.05, 0.01, 0.9])) +colorbar = plt.colorbar(m, cax=ax.figure.add_axes([0.95, 0.03, 0.02, 0.94])) colorbar.set_label('Amplitude (m)') # %% -# draw contour +# Draw speed contours fig = plt.figure(figsize=(15, 8)) ax = fig.add_axes([0.03, 0.03, 0.94, 0.94]) ax.set_aspect("equal") diff --git a/examples/02_eddy_identification/pet_eddy_detection.py b/examples/02_eddy_identification/pet_eddy_detection.py index 130b35ca..ecfe76b6 100644 --- a/examples/02_eddy_identification/pet_eddy_detection.py +++ b/examples/02_eddy_identification/pet_eddy_detection.py @@ -20,7 +20,7 @@ def start_axes(title): ax = fig.add_axes([0.03, 0.03, 0.90, 0.94]) ax.set_xlim(-6, 36.5), ax.set_ylim(30, 46) ax.set_aspect("equal") - ax.set_title(title) + ax.set_title(title, weight='bold') return ax @@ -31,23 +31,23 @@ def update_axes(ax, mappable=None): # %% -# Load Input grid, ADT will be used to detect eddies +# Load Input grid, ADT is used to detect eddies g = RegularGridDataset( data.get_path("dt_med_allsat_phy_l4_20160515_20190101.nc"), "longitude", "latitude" ) ax = start_axes("ADT (m)") -m = g.display(ax, "adt", vmin=-0.15, vmax=0.15) +m = g.display(ax, "adt", vmin=-0.15, vmax=0.15, cmap="RdBu_r") update_axes(ax, m) # %% -# Get u/v +# Get geostrophic speed u,v # ------- -# U/V are deduced from ADT, this algortihm are not usable around equator (~+- 2°) +# U/V are deduced from ADT, this algortihm is not ok near the equator (~+- 2°) g.add_uv("adt") ax = start_axes("U/V deduce from ADT (m)") ax.set_xlim(2.5, 9), ax.set_ylim(37.5, 40) -m = g.display(ax, "adt", vmin=-0.15, vmax=0.15) +m = g.display(ax, "adt", vmin=-0.15, vmax=0.15, cmap="RdBu_r") u, v = g.grid("u").T, g.grid("v").T ax.quiver(g.x_c, g.y_c, u, v, scale=10) update_axes(ax, m) @@ -55,94 +55,92 @@ def update_axes(ax, mappable=None): # %% # Pre-processings # --------------- -# Apply high filter to remove long scale to highlight mesoscale +# Apply a high-pass filter to remove the large scale and highlight the mesoscale g.bessel_high_filter("adt", 500) ax = start_axes("ADT (m) filtered (500km)") -m = g.display(ax, "adt", vmin=-0.15, vmax=0.15) +m = g.display(ax, "adt", vmin=-0.15, vmax=0.15, cmap="RdBu_r") update_axes(ax, m) # %% # Identification # -------------- -# run identification with slice of 2 mm +# Run the identification step with slices of 2 mm date = datetime(2016, 5, 15) a, c = g.eddy_identification("adt", "u", "v", date, 0.002, shape_error=55) # %% -# All closed contour found in this input grid (Display only 1 contour every 4) -ax = start_axes("ADT closed contour (only 1 / 4 levels)") +# Display of all closed contours found in the grid (only 1 contour every 4) +ax = start_axes("ADT closed contours (only 1 / 4 levels)") g.contours.display(ax, step=4) update_axes(ax) # %% -# Contours include in eddies -ax = start_axes("ADT contour used as eddies") +# Contours included in eddies +ax = start_axes("ADT contours used as eddies") g.contours.display(ax, only_used=True) update_axes(ax) # %% -# Post analyse +# Post analysis # ------------ -# Contours reject from several origin (shape error to high, several extremum in contour, ...) -ax = start_axes("ADT contour reject") +# Contours can be rejected for several reasons (shape error to high, several extremum in contour, ...) +ax = start_axes("ADT rejected contours") g.contours.display(ax, only_unused=True) update_axes(ax) # %% -# Contours reject criterion -# +# Creteria for rejecting a contour # 0. - Accepted (green) -# 1. - Reject for shape error (red) -# 2. - Masked value in contour (blue) -# 3. - Under or over pixel limit bound (black) +# 1. - Rejection for shape error (red) +# 2. - Masked value within contour (blue) +# 3. - Under or over the pixel limit bounds (black) # 4. - Amplitude criterion (yellow) -ax = start_axes("Contour reject criterion") -g.contours.display(ax, only_unused=True, lw=0.25, display_criterion=True) +ax = start_axes("Contours' rejection criteria") +g.contours.display(ax, only_unused=True, lw=0.5, display_criterion=True) update_axes(ax) # %% -# Display shape error of each tested contour, the limit of shape error is set to 55 % +# Display the shape error of each tested contour, the limit of shape error is set to 55 % ax = start_axes("Contour shape error") m = g.contours.display( - ax, lw=0.5, field="shape_error", bins=arange(20, 90.1, 5), cmap="seismic" + ax, lw=0.5, field="shape_error", bins=arange(20, 90.1, 5), cmap="PRGn_r" ) update_axes(ax, m) # %% -# Contours closed which contains several eddies -ax = start_axes("ADT contour reject but which contain eddies") +# Some closed contours contains several eddies (aka, more than one extremum) +ax = start_axes("ADT rejected contours containing eddies") g.contours.label_contour_unused_which_contain_eddies(a) g.contours.label_contour_unused_which_contain_eddies(c) g.contours.display( - ax, only_contain_eddies=True, color="k", lw=1, label="Could be interaction contour" + ax, only_contain_eddies=True, color="k", lw=1, label="Could be a contour of interaction" ) -a.display(ax, color="r", linewidth=0.5, label="Anticyclonic", ref=-10) -c.display(ax, color="b", linewidth=0.5, label="Cyclonic", ref=-10) +a.display(ax, color="r", linewidth=0.75, label="Anticyclonic", ref=-10) +c.display(ax, color="b", linewidth=0.75, label="Cyclonic", ref=-10) ax.legend() update_axes(ax) # %% # Output # ------ -# Display detected eddies, dashed lines represent effective contour -# and solid lines represent contour of maximum of speed. See figure 1 of https://doi.org/10.1175/JTECH-D-14-00019.1 +# When displaying the detected eddies, dashed lines are for effective contour, solide lines for the contour of the maximum mean speed. See figure 1 of https://doi.org/10.1175/JTECH-D-14-00019.1 -ax = start_axes("Eddies detected") -a.display(ax, color="r", linewidth=0.5, label="Anticyclonic", ref=-10) -c.display(ax, color="b", linewidth=0.5, label="Cyclonic", ref=-10) +ax = start_axes("Detected Eddies") +a.display(ax, color="r", linewidth=0.75, label="Anticyclonic", ref=-10) +c.display(ax, color="b", linewidth=0.75, label="Cyclonic", ref=-10) ax.legend() update_axes(ax) # %% -# Display speed radius of eddies detected -ax = start_axes("Eddies speed radius (km)") -a.scatter(ax, "radius_s", vmin=10, vmax=50, s=80, ref=-10, cmap="jet", factor=0.001) -m = c.scatter(ax, "radius_s", vmin=10, vmax=50, s=80, ref=-10, cmap="jet", factor=0.001) +# Display the speed radius of the detected eddies +ax = start_axes("Speed Radius (km)") +a.scatter(ax, "radius_s", vmin=10, vmax=50, s=80, ref=-10, cmap="magma_r", factor=0.001) +m = c.scatter(ax, "radius_s", vmin=10, vmax=50, s=80, ref=-10, cmap="magma_r", factor=0.001) update_axes(ax, m) # %% -# Display speed radius of eddies detected with filled eddy contours -ax = start_axes("Eddies speed radius (km)") +# Filling the effective radius contours with the effective radius values +ax = start_axes("Effective Radius (km)") kwargs = dict(vmin=10, vmax=80, cmap="magma_r", factor=0.001, lut=14, ref=-10) a.filled(ax, "radius_e", **kwargs) m = c.filled( diff --git a/examples/02_eddy_identification/pet_eddy_detection_gulf_stream.py b/examples/02_eddy_identification/pet_eddy_detection_gulf_stream.py index 40bd379a..b2dd562c 100644 --- a/examples/02_eddy_identification/pet_eddy_detection_gulf_stream.py +++ b/examples/02_eddy_identification/pet_eddy_detection_gulf_stream.py @@ -21,18 +21,18 @@ def start_axes(title): ax = fig.add_axes([0.03, 0.03, 0.90, 0.94]) ax.set_xlim(279, 304), ax.set_ylim(29, 44) ax.set_aspect("equal") - ax.set_title(title) + ax.set_title(title, weight="bold") return ax def update_axes(ax, mappable=None): ax.grid() if mappable: - plt.colorbar(mappable, cax=ax.figure.add_axes([0.95, 0.05, 0.01, 0.9])) + plt.colorbar(mappable, cax=ax.figure.add_axes([0.94, 0.05, 0.01, 0.9])) # %% -# Load Input grid, ADT will be used to detect eddies +# Load Input grid, ADT is used to detect eddies margin = 30 g = RegularGridDataset( data.get_path("nrt_global_allsat_phy_l4_20190223_20190226.nc"), @@ -46,108 +46,106 @@ def update_axes(ax, mappable=None): ) ax = start_axes("ADT (m)") -m = g.display(ax, "adt", vmin=-0.15, vmax=1) +m = g.display(ax, "adt", vmin=-1, vmax=1, cmap="RdBu_r") # Draw line on the gulf stream front great_current = Contours(g.x_c, g.y_c, g.grid("adt"), levels=(0.35,), keep_unclose=True) great_current.display(ax, color="k") update_axes(ax, m) # %% -# Get u/v +# Get geostrophic speed u,v # ------- -# U/V are deduced from ADT, this algortihm are not usable around equator (~+- 2°) +# U/V are deduced from ADT, this algortihm is not ok near the equator (~+- 2°) g.add_uv("adt") # %% # Pre-processings # --------------- -# Apply high filter to remove long scale to highlight mesoscale +# Apply a high-pass filter to remove the large scale and highlight the mesoscale g.bessel_high_filter("adt", 700) ax = start_axes("ADT (m) filtered (700km)") -m = g.display(ax, "adt", vmin=-0.25, vmax=0.25) +m = g.display(ax, "adt", vmin=-0.4, vmax=0.4, cmap="RdBu_r") great_current.display(ax, color="k") update_axes(ax, m) # %% # Identification # -------------- -# run identification with slice of 2 mm +# Run the identification step with slices of 2 mm date = datetime(2016, 5, 15) a, c = g.eddy_identification("adt", "u", "v", date, 0.002, shape_error=55) # %% -# All closed contour found in this input grid (Display only 1 contour every 5) -ax = start_axes("ADT closed contour (only 1 / 5 levels)") +# Display of all closed contours found in the grid (only 1 contour every 5) +ax = start_axes("ADT closed contours (only 1 / 5 levels)") g.contours.display(ax, step=5, lw=1) great_current.display(ax, color="k") update_axes(ax) # %% -# Contours include in eddies -ax = start_axes("ADT contour used as eddies") +# Contours included in eddies +ax = start_axes("ADT contours used as eddies") g.contours.display(ax, only_used=True, lw=0.25) great_current.display(ax, color="k") update_axes(ax) # %% -# Post analyse +# Post analysis # ------------ -# Contours reject from several origin (shape error to high, several extremum in contour, ...) -ax = start_axes("ADT contour reject") +# Contours can be rejected for several reasons (shape error to high, several extremum in contour, ...) +ax = start_axes("ADT rejected contours") g.contours.display(ax, only_unused=True, lw=0.25) great_current.display(ax, color="k") update_axes(ax) # %% -# Contours reject criterion -# +# Criteria for rejecting a contour # 0. - Accepted (green) -# 1. - Reject for shape error (red) -# 2. - Masked value in contour (blue) -# 3. - Under or over pixel limit bound (black) +# 1. - Rejection for shape error (red) +# 2. - Masked value within contour (blue) +# 3. - Under or over the pixel limit bounds (black) # 4. - Amplitude criterion (yellow) -ax = start_axes("Contour reject criterion") -g.contours.display(ax, only_unused=True, lw=0.25, display_criterion=True) +ax = start_axes("Contours' rejection criteria") +g.contours.display(ax, only_unused=True, lw=0.5, display_criterion=True) update_axes(ax) # %% -# Display shape error of each tested contour, the limit of shape error is set to 55 % +# Display the shape error of each tested contour, the limit of shape error is set to 55 % ax = start_axes("Contour shape error") m = g.contours.display( - ax, lw=0.5, field="shape_error", bins=arange(20, 90.1, 5), cmap="seismic" + ax, lw=0.5, field="shape_error", bins=arange(20, 90.1, 5), cmap="PRGn_r" ) update_axes(ax, m) # %% -# Contours closed which contains several eddies -ax = start_axes("ADT contour reject but which contain eddies") +# Some closed contours contains several eddies (aka, more than one extremum) +ax = start_axes("ADT rejected contours containing eddies") g.contours.label_contour_unused_which_contain_eddies(a) g.contours.label_contour_unused_which_contain_eddies(c) g.contours.display( - ax, only_contain_eddies=True, color="k", lw=1, label="Could be interaction contour" + ax, only_contain_eddies=True, color="k", lw=1, label="Could be a contour of interaction" ) -a.display(ax, color="r", linewidth=0.5, label="Anticyclonic", ref=-10) -c.display(ax, color="b", linewidth=0.5, label="Cyclonic", ref=-10) +a.display(ax, color="r", linewidth=0.75, label="Anticyclonic", ref=-10) +c.display(ax, color="b", linewidth=0.75, label="Cyclonic", ref=-10) ax.legend() update_axes(ax) # %% # Output # ------ -# Display detected eddies, dashed lines represent effective contour -# and solid lines represent contour of maximum of speed. See figure 1 of https://doi.org/10.1175/JTECH-D-14-00019.1 +# When displaying the detected eddies, dashed lines are for effective contour, solide lines for the contour of the maximum mean speed. See figure 1 of https://doi.org/10.1175/JTECH-D-14-00019.1 ax = start_axes("Eddies detected") -a.display(ax, color="r", linewidth=0.5, label="Anticyclonic", ref=-10) -c.display(ax, color="b", linewidth=0.5, label="Cyclonic", ref=-10) +a.display(ax, color="r", linewidth=0.75, label="Anticyclonic", ref=-10) +c.display(ax, color="b", linewidth=0.75, label="Cyclonic", ref=-10) ax.legend() great_current.display(ax, color="k") update_axes(ax) # %% -# Display speed radius of eddies detected -ax = start_axes("Eddies speed radius (km)") +# Display the effective radius of the detected eddies +ax = start_axes("Effective radius (km)") a.filled(ax, "radius_e", vmin=10, vmax=150, cmap="magma_r", factor=0.001, lut=14) m = c.filled(ax, "radius_e", vmin=10, vmax=150, cmap="magma_r", factor=0.001, lut=14) great_current.display(ax, color="k") diff --git a/examples/02_eddy_identification/pet_filter_and_detection.py b/examples/02_eddy_identification/pet_filter_and_detection.py index ce064dc4..03c163bf 100644 --- a/examples/02_eddy_identification/pet_filter_and_detection.py +++ b/examples/02_eddy_identification/pet_filter_and_detection.py @@ -16,65 +16,66 @@ def start_axes(title): ax = fig.add_axes([0.03, 0.03, 0.90, 0.94]) ax.set_xlim(-6, 36.5), ax.set_ylim(30, 46) ax.set_aspect("equal") - ax.set_title(title) + ax.set_title(title, weight="bold") return ax def update_axes(ax, mappable=None): ax.grid() if mappable: - plt.colorbar(mappable, cax=ax.figure.add_axes([0.95, 0.05, 0.01, 0.9])) + plt.colorbar(mappable, cax=ax.figure.add_axes([0.94, 0.05, 0.01, 0.9])) # %% -# Load Input grid, ADT will be used to detect eddies +# Load Input grid, ADT is used to detect eddies. +# Add a new filed to store the high-pass filtered ADT g = RegularGridDataset( data.get_path("dt_med_allsat_phy_l4_20160515_20190101.nc"), "longitude", "latitude", ) g.add_uv("adt") g.copy("adt", "adt_high") -wavelength = 400 +wavelength = 800 g.bessel_high_filter("adt_high", wavelength) date = datetime(2016, 5, 15) # %% -# Run algorithm of detection -a_f, c_f = g.eddy_identification("adt_high", "u", "v", date, 0.002) -merge_f = a_f.merge(c_f) -a_r, c_r = g.eddy_identification("adt", "u", "v", date, 0.002) -merge_r = a_r.merge(c_r) +# Run the detection for the total grid and the filtered grid +a_filtered, c_filtered = g.eddy_identification("adt_high", "u", "v", date, 0.002) +merge_f = a_filtered.merge(c_filtered) +a_tot, c_tot = g.eddy_identification("adt", "u", "v", date, 0.002) +merge_t = a_tot.merge(c_tot) # %% -# Display detection +# Display the two detections ax = start_axes("Eddies detected over ADT") m = g.display(ax, "adt", vmin=-0.15, vmax=0.15) -merge_f.display(ax, lw=0.5, label="Eddy from filtered grid", ref=-10, color="k") -merge_r.display(ax, lw=0.5, label="Eddy from raw grid", ref=-10, color="r") +merge_f.display(ax, lw=0.75, label="Eddies in the filtered grid", ref=-10, color="k") +merge_t.display(ax, lw=0.75, label="Eddies without filter", ref=-10, color="r") ax.legend() update_axes(ax, m) # %% -# Parameters distribution +# Amplitude and Speed Radius distributions # ----------------------- fig = plt.figure(figsize=(12, 5)) -ax_a = fig.add_subplot(121, xlabel="amplitdue(cm)") -ax_r = fig.add_subplot(122, xlabel="speed radius (km)") +ax_a = fig.add_subplot(121, xlabel="Amplitude (cm)") +ax_r = fig.add_subplot(122, xlabel="Speed Radius (km)") ax_a.hist( merge_f["amplitude"] * 100, bins=arange(0.0005, 100, 1), - label="Eddy from filtered grid", + label="Eddies in the filtered grid", histtype="step", ) ax_a.hist( - merge_r["amplitude"] * 100, + merge_t["amplitude"] * 100, bins=arange(0.0005, 100, 1), - label="Eddy from raw grid", + label="Eddies without filter", histtype="step", ) ax_a.set_xlim(0, 10) ax_r.hist(merge_f["radius_s"] / 1000.0, bins=arange(0, 300, 5), histtype="step") -ax_r.hist(merge_r["radius_s"] / 1000.0, bins=arange(0, 300, 5), histtype="step") +ax_r.hist(merge_t["radius_s"] / 1000.0, bins=arange(0, 300, 5), histtype="step") ax_r.set_xlim(0, 100) ax_a.legend() @@ -82,16 +83,16 @@ def update_axes(ax, mappable=None): # Match detection and compare # --------------------------- -i_, j_, c = merge_f.match(merge_r, cmin=0.1) +i_, j_, c = merge_f.match(merge_t, cmin=0.1) # %% -# where is lonely eddies -kwargs_f = dict(lw=1.5, label="Lonely eddy from filtered grid", ref=-10, color="k") -kwargs_r = dict(lw=1.5, label="Lonely eddy from raw grid", ref=-10, color="r") +# Where are the lonely eddies? +kwargs_f = dict(lw=1.5, label="Lonely eddies in the filtered grid", ref=-10, color="k") +kwargs_t = dict(lw=1.5, label="Lonely eddies without filter", ref=-10, color="r") ax = start_axes("Eddies with no match, over filtered ADT") mappable = g.display(ax, "adt_high", vmin=-0.15, vmax=0.15) merge_f.index(i_, reverse=True).display(ax, **kwargs_f) -merge_r.index(j_, reverse=True).display(ax, **kwargs_r) +merge_t.index(j_, reverse=True).display(ax, **kwargs_t) ax.legend() update_axes(ax, mappable) @@ -101,26 +102,26 @@ def update_axes(ax, mappable=None): u, v = g.grid("u").T, g.grid("v").T ax.quiver(g.x_c, g.y_c, u, v, scale=10, pivot="mid", color="gray") merge_f.index(i_, reverse=True).display(ax, **kwargs_f) -merge_r.index(j_, reverse=True).display(ax, **kwargs_r) +merge_t.index(j_, reverse=True).display(ax, **kwargs_t) ax.legend() update_axes(ax, mappable) # %% fig = plt.figure(figsize=(12, 12)) -fig.suptitle(f"Scatter plot ({i_.shape[0]} matches)") +fig.suptitle(f"Scatter plot ({i_.shape[0]} matches)", weight="bold") for i, (label, field, factor, stop) in enumerate( ( - ("speed radius (km)", "radius_s", 0.001, 80), - ("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, 80), + ("Effective radius (km)", "radius_e", 0.001, 120), + ("Amplitude (cm)", "amplitude", 100, 25), + ("Maximum Speed (cm/s)", "speed_average", 100, 25), ) ): ax = fig.add_subplot( - 2, 2, i + 1, xlabel="filtered grid", ylabel="raw grid", title=label + 2, 2, i + 1, xlabel="Filtered grid", ylabel="Without filter", title=label ) - ax.plot(merge_f[field][i_] * factor, merge_r[field][j_] * factor, ".") + ax.plot(merge_f[field][i_] * factor, merge_t[field][j_] * factor, ".") ax.set_aspect("equal"), ax.grid() ax.plot((0, 1000), (0, 1000), "r") ax.set_xlim(0, stop), ax.set_ylim(0, stop) diff --git a/notebooks/python_module/02_eddy_identification/pet_contour_circle.ipynb b/notebooks/python_module/02_eddy_identification/pet_contour_circle.ipynb index 05f7d14e..23c2e6a3 100644 --- a/notebooks/python_module/02_eddy_identification/pet_contour_circle.ipynb +++ b/notebooks/python_module/02_eddy_identification/pet_contour_circle.ipynb @@ -51,7 +51,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "Plot\n\n" + "Plot the speed and effective (dashed) contours\n\n" ] }, { @@ -62,7 +62,7 @@ }, "outputs": [], "source": [ - "fig = plt.figure(figsize=(15, 8))\nax = fig.add_axes((0.05, 0.05, 0.9, 0.9))\nax.set_aspect(\"equal\")\nax.set_xlim(10, 70)\nax.set_ylim(-50, -25)\na.display(ax, label=\"Anticyclonic contour\", color=\"r\", lw=1)\n\n# Replace contour by circle\na.circle_contour()\na.display(ax, label=\"Anticyclonic circle\", color=\"g\", lw=1)\nax.legend(loc=\"upper right\")" + "fig = plt.figure(figsize=(15, 8))\nax = fig.add_axes((0.05, 0.05, 0.9, 0.9))\nax.set_aspect(\"equal\")\nax.set_xlim(10, 70)\nax.set_ylim(-50, -25)\na.display(ax, label=\"Anticyclonic contour\", color=\"r\", lw=1)\n\n# Replace contours by circles using center and radius (effective is dashed)\na.circle_contour()\na.display(ax, label=\"Anticyclonic circle\", color=\"g\", lw=1)\nax.legend(loc=\"upper right\")" ] } ], @@ -82,7 +82,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.7.7" + "version": "3.7.6" } }, "nbformat": 4, diff --git a/notebooks/python_module/02_eddy_identification/pet_display_id.ipynb b/notebooks/python_module/02_eddy_identification/pet_display_id.ipynb index 9273910a..fc22f34b 100644 --- a/notebooks/python_module/02_eddy_identification/pet_display_id.ipynb +++ b/notebooks/python_module/02_eddy_identification/pet_display_id.ipynb @@ -51,7 +51,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "filled contour with amplitude field\n\n" + "Fill effective contour with amplitude \n\n" ] }, { @@ -62,14 +62,14 @@ }, "outputs": [], "source": [ - "fig = plt.figure(figsize=(15, 8))\nax = fig.add_axes([0.03, 0.03, 0.90, 0.94])\nax.set_aspect(\"equal\")\nax.set_xlim(0, 140)\nax.set_ylim(-80, 0)\nkwargs = dict(extern_only=True, color=\"k\", lw=1)\na.display(ax, **kwargs), c.display(ax, **kwargs)\na.filled(ax, \"amplitude\", cmap=\"magma_r\", vmin=0, vmax=.5)\nm = c.filled(ax, \"amplitude\", cmap=\"magma_r\", vmin=0, vmax=.5)\ncolorbar = plt.colorbar(m, cax=ax.figure.add_axes([0.95, 0.05, 0.01, 0.9]))\ncolorbar.set_label('Amplitude (m)')" + "fig = plt.figure(figsize=(15, 8))\nax = fig.add_axes([0.03, 0.03, 0.90, 0.94])\nax.set_aspect(\"equal\")\nax.set_xlim(0, 140)\nax.set_ylim(-80, 0)\nkwargs = dict(extern_only=True, color=\"k\", lw=1)\na.display(ax, **kwargs), c.display(ax, **kwargs)\na.filled(ax, \"amplitude\", cmap=\"magma_r\", vmin=0, vmax=.5)\nm = c.filled(ax, \"amplitude\", cmap=\"magma_r\", vmin=0, vmax=.5)\ncolorbar = plt.colorbar(m, cax=ax.figure.add_axes([0.95, 0.03, 0.02, 0.94]))\ncolorbar.set_label('Amplitude (m)')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "draw contour\n\n" + "Draw speed contours\n\n" ] }, { @@ -129,7 +129,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.7.7" + "version": "3.7.6" } }, "nbformat": 4, diff --git a/notebooks/python_module/02_eddy_identification/pet_eddy_detection.ipynb b/notebooks/python_module/02_eddy_identification/pet_eddy_detection.ipynb index 9e59b80f..18412d39 100644 --- a/notebooks/python_module/02_eddy_identification/pet_eddy_detection.ipynb +++ b/notebooks/python_module/02_eddy_identification/pet_eddy_detection.ipynb @@ -37,14 +37,14 @@ }, "outputs": [], "source": [ - "def start_axes(title):\n fig = plt.figure(figsize=(13, 5))\n ax = fig.add_axes([0.03, 0.03, 0.90, 0.94])\n ax.set_xlim(-6, 36.5), ax.set_ylim(30, 46)\n ax.set_aspect(\"equal\")\n ax.set_title(title)\n return ax\n\n\ndef update_axes(ax, mappable=None):\n ax.grid()\n if mappable:\n plt.colorbar(mappable, cax=ax.figure.add_axes([0.94, 0.05, 0.01, 0.9]))" + "def start_axes(title):\n fig = plt.figure(figsize=(13, 5))\n ax = fig.add_axes([0.03, 0.03, 0.90, 0.94])\n ax.set_xlim(-6, 36.5), ax.set_ylim(30, 46)\n ax.set_aspect(\"equal\")\n ax.set_title(title, weight='bold')\n return ax\n\n\ndef update_axes(ax, mappable=None):\n ax.grid()\n if mappable:\n plt.colorbar(mappable, cax=ax.figure.add_axes([0.94, 0.05, 0.01, 0.9]))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "Load Input grid, ADT will be used to detect eddies\n\n" + "Load Input grid, ADT is used to detect eddies\n\n" ] }, { @@ -55,14 +55,14 @@ }, "outputs": [], "source": [ - "g = RegularGridDataset(\n data.get_path(\"dt_med_allsat_phy_l4_20160515_20190101.nc\"), \"longitude\", \"latitude\"\n)\n\nax = start_axes(\"ADT (m)\")\nm = g.display(ax, \"adt\", vmin=-0.15, vmax=0.15)\nupdate_axes(ax, m)" + "g = RegularGridDataset(\n data.get_path(\"dt_med_allsat_phy_l4_20160515_20190101.nc\"), \"longitude\", \"latitude\"\n)\n\nax = start_axes(\"ADT (m)\")\nm = g.display(ax, \"adt\", vmin=-0.15, vmax=0.15, cmap=\"RdBu_r\")\nupdate_axes(ax, m)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "## Get u/v\nU/V are deduced from ADT, this algortihm are not usable around equator (~+- 2\u00b0)\n\n" + "## Get geostrophic speed u,v\nU/V are deduced from ADT, this algortihm is not ok near the equator (~+- 2\u00b0)\n\n" ] }, { @@ -73,14 +73,14 @@ }, "outputs": [], "source": [ - "g.add_uv(\"adt\")\nax = start_axes(\"U/V deduce from ADT (m)\")\nax.set_xlim(2.5, 9), ax.set_ylim(37.5, 40)\nm = g.display(ax, \"adt\", vmin=-0.15, vmax=0.15)\nu, v = g.grid(\"u\").T, g.grid(\"v\").T\nax.quiver(g.x_c, g.y_c, u, v, scale=10)\nupdate_axes(ax, m)" + "g.add_uv(\"adt\")\nax = start_axes(\"U/V deduce from ADT (m)\")\nax.set_xlim(2.5, 9), ax.set_ylim(37.5, 40)\nm = g.display(ax, \"adt\", vmin=-0.15, vmax=0.15, cmap=\"RdBu_r\")\nu, v = g.grid(\"u\").T, g.grid(\"v\").T\nax.quiver(g.x_c, g.y_c, u, v, scale=10)\nupdate_axes(ax, m)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "## Pre-processings\nApply high filter to remove long scale to highlight mesoscale\n\n" + "## Pre-processings\nApply a high-pass filter to remove the large scale and highlight the mesoscale\n\n" ] }, { @@ -91,14 +91,14 @@ }, "outputs": [], "source": [ - "g.bessel_high_filter(\"adt\", 500)\nax = start_axes(\"ADT (m) filtered (500km)\")\nm = g.display(ax, \"adt\", vmin=-0.15, vmax=0.15)\nupdate_axes(ax, m)" + "g.bessel_high_filter(\"adt\", 500)\nax = start_axes(\"ADT (m) filtered (500km)\")\nm = g.display(ax, \"adt\", vmin=-0.15, vmax=0.15, cmap=\"RdBu_r\")\nupdate_axes(ax, m)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "## Identification\nrun identification with slice of 2 mm\n\n" + "## Identification\nRun the identification step with slices of 2 mm\n\n" ] }, { @@ -116,7 +116,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "All closed contour found in this input grid (Display only 1 contour every 4)\n\n" + "Display of all closed contours found in the grid (only 1 contour every 4)\n\n" ] }, { @@ -127,14 +127,14 @@ }, "outputs": [], "source": [ - "ax = start_axes(\"ADT closed contour (only 1 / 4 levels)\")\ng.contours.display(ax, step=4)\nupdate_axes(ax)" + "ax = start_axes(\"ADT closed contours (only 1 / 4 levels)\")\ng.contours.display(ax, step=4)\nupdate_axes(ax)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "Contours include in eddies\n\n" + "Contours included in eddies\n\n" ] }, { @@ -145,14 +145,14 @@ }, "outputs": [], "source": [ - "ax = start_axes(\"ADT contour used as eddies\")\ng.contours.display(ax, only_used=True)\nupdate_axes(ax)" + "ax = start_axes(\"ADT contours used as eddies\")\ng.contours.display(ax, only_used=True)\nupdate_axes(ax)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "## Post analyse\nContours reject from several origin (shape error to high, several extremum in contour, ...)\n\n" + "## Post analysis\nContours can be rejected for several reasons (shape error to high, several extremum in contour, ...)\n\n" ] }, { @@ -163,14 +163,14 @@ }, "outputs": [], "source": [ - "ax = start_axes(\"ADT contour reject\")\ng.contours.display(ax, only_unused=True)\nupdate_axes(ax)" + "ax = start_axes(\"ADT rejected contours\")\ng.contours.display(ax, only_unused=True)\nupdate_axes(ax)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "Contours reject criterion\n\n0. - Accepted (green)\n1. - Reject for shape error (red)\n2. - Masked value in contour (blue)\n3. - Under or over pixel limit bound (black)\n4. - Amplitude criterion (yellow)\n\n" + "Creteria for rejecting a contour\n0. - Accepted (green)\n1. - Rejection for shape error (red)\n2. - Masked value within contour (blue)\n3. - Under or over the pixel limit bounds (black)\n4. - Amplitude criterion (yellow)\n\n" ] }, { @@ -181,14 +181,14 @@ }, "outputs": [], "source": [ - "ax = start_axes(\"Contour reject criterion\")\ng.contours.display(ax, only_unused=True, lw=0.25, display_criterion=True)\nupdate_axes(ax)" + "ax = start_axes(\"Contours' rejection criteria\")\ng.contours.display(ax, only_unused=True, lw=0.5, display_criterion=True)\nupdate_axes(ax)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "Display shape error of each tested contour, the limit of shape error is set to 55 %\n\n" + "Display the shape error of each tested contour, the limit of shape error is set to 55 %\n\n" ] }, { @@ -199,14 +199,14 @@ }, "outputs": [], "source": [ - "ax = start_axes(\"Contour shape error\")\nm = g.contours.display(\n ax, lw=0.5, field=\"shape_error\", bins=arange(20, 90.1, 5), cmap=\"seismic\"\n)\nupdate_axes(ax, m)" + "ax = start_axes(\"Contour shape error\")\nm = g.contours.display(\n ax, lw=0.5, field=\"shape_error\", bins=arange(20, 90.1, 5), cmap=\"PRGn_r\"\n)\nupdate_axes(ax, m)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "Contours closed which contains several eddies\n\n" + "Some closed contours contains several eddies (aka, more than one extremum)\n\n" ] }, { @@ -217,14 +217,14 @@ }, "outputs": [], "source": [ - "ax = start_axes(\"ADT contour reject but which contain eddies\")\ng.contours.label_contour_unused_which_contain_eddies(a)\ng.contours.label_contour_unused_which_contain_eddies(c)\ng.contours.display(\n ax, only_contain_eddies=True, color=\"k\", lw=1, label=\"Could be interaction contour\"\n)\na.display(ax, color=\"r\", linewidth=0.5, label=\"Anticyclonic\", ref=-10)\nc.display(ax, color=\"b\", linewidth=0.5, label=\"Cyclonic\", ref=-10)\nax.legend()\nupdate_axes(ax)" + "ax = start_axes(\"ADT rejected contours containing eddies\")\ng.contours.label_contour_unused_which_contain_eddies(a)\ng.contours.label_contour_unused_which_contain_eddies(c)\ng.contours.display(\n ax, only_contain_eddies=True, color=\"k\", lw=1, label=\"Could be a contour of interaction\"\n)\na.display(ax, color=\"r\", linewidth=0.75, label=\"Anticyclonic\", ref=-10)\nc.display(ax, color=\"b\", linewidth=0.75, label=\"Cyclonic\", ref=-10)\nax.legend()\nupdate_axes(ax)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "## Output\nDisplay detected eddies, dashed lines represent effective contour\nand solid lines represent contour of maximum of speed. See figure 1 of https://doi.org/10.1175/JTECH-D-14-00019.1\n\n" + "## Output\nWhen displaying the detected eddies, dashed lines are for effective contour, solide lines for the contour of the maximum mean speed. See figure 1 of https://doi.org/10.1175/JTECH-D-14-00019.1\n\n" ] }, { @@ -235,14 +235,14 @@ }, "outputs": [], "source": [ - "ax = start_axes(\"Eddies detected\")\na.display(ax, color=\"r\", linewidth=0.5, label=\"Anticyclonic\", ref=-10)\nc.display(ax, color=\"b\", linewidth=0.5, label=\"Cyclonic\", ref=-10)\nax.legend()\nupdate_axes(ax)" + "ax = start_axes(\"Detected Eddies\")\na.display(ax, color=\"r\", linewidth=0.75, label=\"Anticyclonic\", ref=-10)\nc.display(ax, color=\"b\", linewidth=0.75, label=\"Cyclonic\", ref=-10)\nax.legend()\nupdate_axes(ax)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "Display speed radius of eddies detected\n\n" + "Display the speed radius of the detected eddies\n\n" ] }, { @@ -253,14 +253,14 @@ }, "outputs": [], "source": [ - "ax = start_axes(\"Eddies speed radius (km)\")\na.scatter(ax, \"radius_s\", vmin=10, vmax=50, s=80, ref=-10, cmap=\"jet\", factor=0.001)\nm = c.scatter(ax, \"radius_s\", vmin=10, vmax=50, s=80, ref=-10, cmap=\"jet\", factor=0.001)\nupdate_axes(ax, m)" + "ax = start_axes(\"Speed Radius (km)\")\na.scatter(ax, \"radius_s\", vmin=10, vmax=50, s=80, ref=-10, cmap=\"magma_r\", factor=0.001)\nm = c.scatter(ax, \"radius_s\", vmin=10, vmax=50, s=80, ref=-10, cmap=\"magma_r\", factor=0.001)\nupdate_axes(ax, m)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "Display speed radius of eddies detected with filled eddy contours\n\n" + "Filling the effective radius contours with the effective radius values\n\n" ] }, { @@ -271,7 +271,7 @@ }, "outputs": [], "source": [ - "ax = start_axes(\"Eddies speed radius (km)\")\nkwargs = dict(vmin=10, vmax=80, cmap=\"magma_r\", factor=0.001, lut=14, ref=-10)\na.filled(ax, \"radius_e\", **kwargs)\nm = c.filled(\n ax, \"radius_e\", vmin=10, vmax=80, cmap=\"magma_r\", factor=0.001, lut=14, ref=-10\n)\nupdate_axes(ax, m)" + "ax = start_axes(\"Effective Radius (km)\")\nkwargs = dict(vmin=10, vmax=80, cmap=\"magma_r\", factor=0.001, lut=14, ref=-10)\na.filled(ax, \"radius_e\", **kwargs)\nm = c.filled(\n ax, \"radius_e\", vmin=10, vmax=80, cmap=\"magma_r\", factor=0.001, lut=14, ref=-10\n)\nupdate_axes(ax, m)" ] } ], @@ -291,7 +291,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.7.7" + "version": "3.7.6" } }, "nbformat": 4, diff --git a/notebooks/python_module/02_eddy_identification/pet_eddy_detection_gulf_stream.ipynb b/notebooks/python_module/02_eddy_identification/pet_eddy_detection_gulf_stream.ipynb index 1e4f484e..644171bd 100644 --- a/notebooks/python_module/02_eddy_identification/pet_eddy_detection_gulf_stream.ipynb +++ b/notebooks/python_module/02_eddy_identification/pet_eddy_detection_gulf_stream.ipynb @@ -37,14 +37,14 @@ }, "outputs": [], "source": [ - "def start_axes(title):\n fig = plt.figure(figsize=(13, 8))\n ax = fig.add_axes([0.03, 0.03, 0.90, 0.94])\n ax.set_xlim(279, 304), ax.set_ylim(29, 44)\n ax.set_aspect(\"equal\")\n ax.set_title(title)\n return ax\n\n\ndef update_axes(ax, mappable=None):\n ax.grid()\n if mappable:\n plt.colorbar(mappable, cax=ax.figure.add_axes([0.95, 0.05, 0.01, 0.9]))" + "def start_axes(title):\n fig = plt.figure(figsize=(13, 8))\n ax = fig.add_axes([0.03, 0.03, 0.90, 0.94])\n ax.set_xlim(279, 304), ax.set_ylim(29, 44)\n ax.set_aspect(\"equal\")\n ax.set_title(title, weight=\"bold\")\n return ax\n\n\ndef update_axes(ax, mappable=None):\n ax.grid()\n if mappable:\n plt.colorbar(mappable, cax=ax.figure.add_axes([0.94, 0.05, 0.01, 0.9]))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "Load Input grid, ADT will be used to detect eddies\n\n" + "Load Input grid, ADT is used to detect eddies\n\n" ] }, { @@ -55,14 +55,14 @@ }, "outputs": [], "source": [ - "margin = 30\ng = RegularGridDataset(\n data.get_path(\"nrt_global_allsat_phy_l4_20190223_20190226.nc\"),\n \"longitude\",\n \"latitude\",\n # Manual area subset\n indexs=dict(\n longitude=slice(1116 - margin, 1216 + margin),\n latitude=slice(476 - margin, 536 + margin),\n ),\n)\n\nax = start_axes(\"ADT (m)\")\nm = g.display(ax, \"adt\", vmin=-0.15, vmax=1)\n# Draw line on the gulf stream front\ngreat_current = Contours(g.x_c, g.y_c, g.grid(\"adt\"), levels=(0.35,), keep_unclose=True)\ngreat_current.display(ax, color=\"k\")\nupdate_axes(ax, m)" + "margin = 30\ng = RegularGridDataset(\n data.get_path(\"nrt_global_allsat_phy_l4_20190223_20190226.nc\"),\n \"longitude\",\n \"latitude\",\n # Manual area subset\n indexs=dict(\n longitude=slice(1116 - margin, 1216 + margin),\n latitude=slice(476 - margin, 536 + margin),\n ),\n)\n\nax = start_axes(\"ADT (m)\")\nm = g.display(ax, \"adt\", vmin=-1, vmax=1, cmap=\"RdBu_r\")\n# Draw line on the gulf stream front\ngreat_current = Contours(g.x_c, g.y_c, g.grid(\"adt\"), levels=(0.35,), keep_unclose=True)\ngreat_current.display(ax, color=\"k\")\nupdate_axes(ax, m)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "## Get u/v\nU/V are deduced from ADT, this algortihm are not usable around equator (~+- 2\u00b0)\n\n" + "## Get geostrophic speed u,v\nU/V are deduced from ADT, this algortihm is not ok near the equator (~+- 2\u00b0)\n\n" ] }, { @@ -80,7 +80,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "## Pre-processings\nApply high filter to remove long scale to highlight mesoscale\n\n" + "## Pre-processings\nApply a high-pass filter to remove the large scale and highlight the mesoscale\n\n" ] }, { @@ -91,14 +91,14 @@ }, "outputs": [], "source": [ - "g.bessel_high_filter(\"adt\", 700)\nax = start_axes(\"ADT (m) filtered (700km)\")\nm = g.display(ax, \"adt\", vmin=-0.25, vmax=0.25)\ngreat_current.display(ax, color=\"k\")\nupdate_axes(ax, m)" + "g.bessel_high_filter(\"adt\", 700)\nax = start_axes(\"ADT (m) filtered (700km)\")\nm = g.display(ax, \"adt\", vmin=-0.4, vmax=0.4, cmap=\"RdBu_r\")\ngreat_current.display(ax, color=\"k\")\nupdate_axes(ax, m)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "## Identification\nrun identification with slice of 2 mm\n\n" + "## Identification\nRun the identification step with slices of 2 mm\n\n" ] }, { @@ -116,7 +116,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "All closed contour found in this input grid (Display only 1 contour every 5)\n\n" + "Display of all closed contours found in the grid (only 1 contour every 5)\n\n" ] }, { @@ -127,14 +127,14 @@ }, "outputs": [], "source": [ - "ax = start_axes(\"ADT closed contour (only 1 / 5 levels)\")\ng.contours.display(ax, step=5, lw=1)\ngreat_current.display(ax, color=\"k\")\nupdate_axes(ax)" + "ax = start_axes(\"ADT closed contours (only 1 / 5 levels)\")\ng.contours.display(ax, step=5, lw=1)\ngreat_current.display(ax, color=\"k\")\nupdate_axes(ax)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "Contours include in eddies\n\n" + "Contours included in eddies\n\n" ] }, { @@ -145,14 +145,14 @@ }, "outputs": [], "source": [ - "ax = start_axes(\"ADT contour used as eddies\")\ng.contours.display(ax, only_used=True, lw=0.25)\ngreat_current.display(ax, color=\"k\")\nupdate_axes(ax)" + "ax = start_axes(\"ADT contours used as eddies\")\ng.contours.display(ax, only_used=True, lw=0.25)\ngreat_current.display(ax, color=\"k\")\nupdate_axes(ax)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "## Post analyse\nContours reject from several origin (shape error to high, several extremum in contour, ...)\n\n" + "## Post analysis\nContours can be rejected for several reasons (shape error to high, several extremum in contour, ...)\n\n" ] }, { @@ -163,14 +163,14 @@ }, "outputs": [], "source": [ - "ax = start_axes(\"ADT contour reject\")\ng.contours.display(ax, only_unused=True, lw=0.25)\ngreat_current.display(ax, color=\"k\")\nupdate_axes(ax)" + "ax = start_axes(\"ADT rejected contours\")\ng.contours.display(ax, only_unused=True, lw=0.25)\ngreat_current.display(ax, color=\"k\")\nupdate_axes(ax)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "Contours reject criterion\n\n0. - Accepted (green)\n1. - Reject for shape error (red)\n2. - Masked value in contour (blue)\n3. - Under or over pixel limit bound (black)\n4. - Amplitude criterion (yellow)\n\n" + "Criteria for rejecting a contour\n0. - Accepted (green)\n1. - Rejection for shape error (red)\n2. - Masked value within contour (blue)\n3. - Under or over the pixel limit bounds (black)\n4. - Amplitude criterion (yellow)\n\n" ] }, { @@ -181,14 +181,14 @@ }, "outputs": [], "source": [ - "ax = start_axes(\"Contour reject criterion\")\ng.contours.display(ax, only_unused=True, lw=0.25, display_criterion=True)\nupdate_axes(ax)" + "ax = start_axes(\"Contours' rejection criteria\")\ng.contours.display(ax, only_unused=True, lw=0.5, display_criterion=True)\nupdate_axes(ax)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "Display shape error of each tested contour, the limit of shape error is set to 55 %\n\n" + "Display the shape error of each tested contour, the limit of shape error is set to 55 %\n\n" ] }, { @@ -199,14 +199,14 @@ }, "outputs": [], "source": [ - "ax = start_axes(\"Contour shape error\")\nm = g.contours.display(\n ax, lw=0.5, field=\"shape_error\", bins=arange(20, 90.1, 5), cmap=\"seismic\"\n)\nupdate_axes(ax, m)" + "ax = start_axes(\"Contour shape error\")\nm = g.contours.display(\n ax, lw=0.5, field=\"shape_error\", bins=arange(20, 90.1, 5), cmap=\"PRGn_r\"\n)\nupdate_axes(ax, m)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "Contours closed which contains several eddies\n\n" + "Some closed contours contains several eddies (aka, more than one extremum)\n\n" ] }, { @@ -217,14 +217,14 @@ }, "outputs": [], "source": [ - "ax = start_axes(\"ADT contour reject but which contain eddies\")\ng.contours.label_contour_unused_which_contain_eddies(a)\ng.contours.label_contour_unused_which_contain_eddies(c)\ng.contours.display(\n ax, only_contain_eddies=True, color=\"k\", lw=1, label=\"Could be interaction contour\"\n)\na.display(ax, color=\"r\", linewidth=0.5, label=\"Anticyclonic\", ref=-10)\nc.display(ax, color=\"b\", linewidth=0.5, label=\"Cyclonic\", ref=-10)\nax.legend()\nupdate_axes(ax)" + "ax = start_axes(\"ADT rejected contours containing eddies\")\ng.contours.label_contour_unused_which_contain_eddies(a)\ng.contours.label_contour_unused_which_contain_eddies(c)\ng.contours.display(\n ax, only_contain_eddies=True, color=\"k\", lw=1, label=\"Could be a contour of interaction\"\n)\na.display(ax, color=\"r\", linewidth=0.75, label=\"Anticyclonic\", ref=-10)\nc.display(ax, color=\"b\", linewidth=0.75, label=\"Cyclonic\", ref=-10)\nax.legend()\nupdate_axes(ax)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "## Output\nDisplay detected eddies, dashed lines represent effective contour\nand solid lines represent contour of maximum of speed. See figure 1 of https://doi.org/10.1175/JTECH-D-14-00019.1\n\n" + "## Output\nWhen displaying the detected eddies, dashed lines are for effective contour, solide lines for the contour of the maximum mean speed. See figure 1 of https://doi.org/10.1175/JTECH-D-14-00019.1\n\n" ] }, { @@ -235,14 +235,14 @@ }, "outputs": [], "source": [ - "ax = start_axes(\"Eddies detected\")\na.display(ax, color=\"r\", linewidth=0.5, label=\"Anticyclonic\", ref=-10)\nc.display(ax, color=\"b\", linewidth=0.5, label=\"Cyclonic\", ref=-10)\nax.legend()\ngreat_current.display(ax, color=\"k\")\nupdate_axes(ax)" + "ax = start_axes(\"Eddies detected\")\na.display(ax, color=\"r\", linewidth=0.75, label=\"Anticyclonic\", ref=-10)\nc.display(ax, color=\"b\", linewidth=0.75, label=\"Cyclonic\", ref=-10)\nax.legend()\ngreat_current.display(ax, color=\"k\")\nupdate_axes(ax)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "Display speed radius of eddies detected\n\n" + "Display the effective radius of the detected eddies\n\n" ] }, { @@ -253,7 +253,7 @@ }, "outputs": [], "source": [ - "ax = start_axes(\"Eddies speed radius (km)\")\na.filled(ax, \"radius_e\", vmin=10, vmax=150, cmap=\"magma_r\", factor=0.001, lut=14)\nm = c.filled(ax, \"radius_e\", vmin=10, vmax=150, cmap=\"magma_r\", factor=0.001, lut=14)\ngreat_current.display(ax, color=\"k\")\nupdate_axes(ax, m)" + "ax = start_axes(\"Effective radius (km)\")\na.filled(ax, \"radius_e\", vmin=10, vmax=150, cmap=\"magma_r\", factor=0.001, lut=14)\nm = c.filled(ax, \"radius_e\", vmin=10, vmax=150, cmap=\"magma_r\", factor=0.001, lut=14)\ngreat_current.display(ax, color=\"k\")\nupdate_axes(ax, m)" ] } ], @@ -273,7 +273,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.7.7" + "version": "3.7.6" } }, "nbformat": 4, diff --git a/notebooks/python_module/02_eddy_identification/pet_filter_and_detection.ipynb b/notebooks/python_module/02_eddy_identification/pet_filter_and_detection.ipynb index c894797c..5a1e173c 100644 --- a/notebooks/python_module/02_eddy_identification/pet_filter_and_detection.ipynb +++ b/notebooks/python_module/02_eddy_identification/pet_filter_and_detection.ipynb @@ -37,14 +37,14 @@ }, "outputs": [], "source": [ - "def start_axes(title):\n fig = plt.figure(figsize=(13, 5))\n ax = fig.add_axes([0.03, 0.03, 0.90, 0.94])\n ax.set_xlim(-6, 36.5), ax.set_ylim(30, 46)\n ax.set_aspect(\"equal\")\n ax.set_title(title)\n return ax\n\n\ndef update_axes(ax, mappable=None):\n ax.grid()\n if mappable:\n plt.colorbar(mappable, cax=ax.figure.add_axes([0.95, 0.05, 0.01, 0.9]))" + "def start_axes(title):\n fig = plt.figure(figsize=(13, 5))\n ax = fig.add_axes([0.03, 0.03, 0.90, 0.94])\n ax.set_xlim(-6, 36.5), ax.set_ylim(30, 46)\n ax.set_aspect(\"equal\")\n ax.set_title(title, weight=\"bold\")\n return ax\n\n\ndef update_axes(ax, mappable=None):\n ax.grid()\n if mappable:\n plt.colorbar(mappable, cax=ax.figure.add_axes([0.94, 0.05, 0.01, 0.9]))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "Load Input grid, ADT will be used to detect eddies\n\n" + "Load Input grid, ADT is used to detect eddies.\nAdd a new filed to store the high-pass filtered ADT\n\n" ] }, { @@ -55,14 +55,14 @@ }, "outputs": [], "source": [ - "g = RegularGridDataset(\n data.get_path(\"dt_med_allsat_phy_l4_20160515_20190101.nc\"), \"longitude\", \"latitude\",\n)\ng.add_uv(\"adt\")\ng.copy(\"adt\", \"adt_high\")\nwavelength = 400\ng.bessel_high_filter(\"adt_high\", wavelength)\ndate = datetime(2016, 5, 15)" + "g = RegularGridDataset(\n data.get_path(\"dt_med_allsat_phy_l4_20160515_20190101.nc\"), \"longitude\", \"latitude\",\n)\ng.add_uv(\"adt\")\ng.copy(\"adt\", \"adt_high\")\nwavelength = 800\ng.bessel_high_filter(\"adt_high\", wavelength)\ndate = datetime(2016, 5, 15)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "Run algorithm of detection\n\n" + "Run the detection for the total grid and the filtered grid\n\n" ] }, { @@ -73,14 +73,14 @@ }, "outputs": [], "source": [ - "a_f, c_f = g.eddy_identification(\"adt_high\", \"u\", \"v\", date, 0.002)\nmerge_f = a_f.merge(c_f)\na_r, c_r = g.eddy_identification(\"adt\", \"u\", \"v\", date, 0.002)\nmerge_r = a_r.merge(c_r)" + "a_filtered, c_filtered = g.eddy_identification(\"adt_high\", \"u\", \"v\", date, 0.002)\nmerge_f = a_filtered.merge(c_filtered)\na_tot, c_tot = g.eddy_identification(\"adt\", \"u\", \"v\", date, 0.002)\nmerge_t = a_tot.merge(c_tot)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "Display detection\n\n" + "Display the two detections\n\n" ] }, { @@ -91,14 +91,14 @@ }, "outputs": [], "source": [ - "ax = start_axes(\"Eddies detected over ADT\")\nm = g.display(ax, \"adt\", vmin=-0.15, vmax=0.15)\nmerge_f.display(ax, lw=0.5, label=\"Eddy from filtered grid\", ref=-10, color=\"k\")\nmerge_r.display(ax, lw=0.5, label=\"Eddy from raw grid\", ref=-10, color=\"r\")\nax.legend()\nupdate_axes(ax, m)" + "ax = start_axes(\"Eddies detected over ADT\")\nm = g.display(ax, \"adt\", vmin=-0.15, vmax=0.15)\nmerge_f.display(ax, lw=0.75, label=\"Eddies in the filtered grid\", ref=-10, color=\"k\")\nmerge_t.display(ax, lw=0.75, label=\"Eddies without filter\", ref=-10, color=\"r\")\nax.legend()\nupdate_axes(ax, m)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "## Parameters distribution\n\n" + "## Amplitude and Speed Radius distributions\n\n" ] }, { @@ -109,7 +109,7 @@ }, "outputs": [], "source": [ - "fig = plt.figure(figsize=(12, 5))\nax_a = fig.add_subplot(121, xlabel=\"amplitdue(cm)\")\nax_r = fig.add_subplot(122, xlabel=\"speed radius (km)\")\nax_a.hist(\n merge_f[\"amplitude\"] * 100,\n bins=arange(0.0005, 100, 1),\n label=\"Eddy from filtered grid\",\n histtype=\"step\",\n)\nax_a.hist(\n merge_r[\"amplitude\"] * 100,\n bins=arange(0.0005, 100, 1),\n label=\"Eddy from raw grid\",\n histtype=\"step\",\n)\nax_a.set_xlim(0, 10)\nax_r.hist(merge_f[\"radius_s\"] / 1000.0, bins=arange(0, 300, 5), histtype=\"step\")\nax_r.hist(merge_r[\"radius_s\"] / 1000.0, bins=arange(0, 300, 5), histtype=\"step\")\nax_r.set_xlim(0, 100)\nax_a.legend()" + "fig = plt.figure(figsize=(12, 5))\nax_a = fig.add_subplot(121, xlabel=\"Amplitude (cm)\")\nax_r = fig.add_subplot(122, xlabel=\"Speed Radius (km)\")\nax_a.hist(\n merge_f[\"amplitude\"] * 100,\n bins=arange(0.0005, 100, 1),\n label=\"Eddies in the filtered grid\",\n histtype=\"step\",\n)\nax_a.hist(\n merge_t[\"amplitude\"] * 100,\n bins=arange(0.0005, 100, 1),\n label=\"Eddies without filter\",\n histtype=\"step\",\n)\nax_a.set_xlim(0, 10)\nax_r.hist(merge_f[\"radius_s\"] / 1000.0, bins=arange(0, 300, 5), histtype=\"step\")\nax_r.hist(merge_t[\"radius_s\"] / 1000.0, bins=arange(0, 300, 5), histtype=\"step\")\nax_r.set_xlim(0, 100)\nax_a.legend()" ] }, { @@ -127,14 +127,14 @@ }, "outputs": [], "source": [ - "i_, j_, c = merge_f.match(merge_r, cmin=0.1)" + "i_, j_, c = merge_f.match(merge_t, cmin=0.1)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ - "where is lonely eddies\n\n" + "Where are the lonely eddies?\n\n" ] }, { @@ -145,7 +145,7 @@ }, "outputs": [], "source": [ - "kwargs_f = dict(lw=1.5, label=\"Lonely eddy from filtered grid\", ref=-10, color=\"k\")\nkwargs_r = dict(lw=1.5, label=\"Lonely eddy from raw grid\", ref=-10, color=\"r\")\nax = start_axes(\"Eddies with no match, over filtered ADT\")\nmappable = g.display(ax, \"adt_high\", vmin=-0.15, vmax=0.15)\nmerge_f.index(i_, reverse=True).display(ax, **kwargs_f)\nmerge_r.index(j_, reverse=True).display(ax, **kwargs_r)\nax.legend()\nupdate_axes(ax, mappable)\n\nax = start_axes(\"Eddies with no match, over filtered ADT (zoom)\")\nax.set_xlim(25, 36), ax.set_ylim(31, 35.25)\nmappable = g.display(ax, \"adt_high\", vmin=-0.15, vmax=0.15)\nu, v = g.grid(\"u\").T, g.grid(\"v\").T\nax.quiver(g.x_c, g.y_c, u, v, scale=10, pivot=\"mid\", color=\"gray\")\nmerge_f.index(i_, reverse=True).display(ax, **kwargs_f)\nmerge_r.index(j_, reverse=True).display(ax, **kwargs_r)\nax.legend()\nupdate_axes(ax, mappable)" + "kwargs_f = dict(lw=1.5, label=\"Lonely eddies in the filtered grid\", ref=-10, color=\"k\")\nkwargs_t = dict(lw=1.5, label=\"Lonely eddies without filter\", ref=-10, color=\"r\")\nax = start_axes(\"Eddies with no match, over filtered ADT\")\nmappable = g.display(ax, \"adt_high\", vmin=-0.15, vmax=0.15)\nmerge_f.index(i_, reverse=True).display(ax, **kwargs_f)\nmerge_t.index(j_, reverse=True).display(ax, **kwargs_t)\nax.legend()\nupdate_axes(ax, mappable)\n\nax = start_axes(\"Eddies with no match, over filtered ADT (zoom)\")\nax.set_xlim(25, 36), ax.set_ylim(31, 35.25)\nmappable = g.display(ax, \"adt_high\", vmin=-0.15, vmax=0.15)\nu, v = g.grid(\"u\").T, g.grid(\"v\").T\nax.quiver(g.x_c, g.y_c, u, v, scale=10, pivot=\"mid\", color=\"gray\")\nmerge_f.index(i_, reverse=True).display(ax, **kwargs_f)\nmerge_t.index(j_, reverse=True).display(ax, **kwargs_t)\nax.legend()\nupdate_axes(ax, mappable)" ] }, { @@ -156,7 +156,7 @@ }, "outputs": [], "source": [ - "fig = plt.figure(figsize=(12, 12))\nfig.suptitle(f\"Scatter plot ({i_.shape[0]} matches)\")\n\nfor i, (label, field, factor, stop) in enumerate(\n (\n (\"speed radius (km)\", \"radius_s\", 0.001, 80),\n (\"outter radius (km)\", \"radius_e\", 0.001, 120),\n (\"amplitude (cm)\", \"amplitude\", 100, 25),\n (\"speed max (cm/s)\", \"speed_average\", 100, 25),\n )\n):\n ax = fig.add_subplot(\n 2, 2, i + 1, xlabel=\"filtered grid\", ylabel=\"raw grid\", title=label\n )\n ax.plot(merge_f[field][i_] * factor, merge_r[field][j_] * factor, \".\")\n ax.set_aspect(\"equal\"), ax.grid()\n ax.plot((0, 1000), (0, 1000), \"r\")\n ax.set_xlim(0, stop), ax.set_ylim(0, stop)" + "fig = plt.figure(figsize=(12, 12))\nfig.suptitle(f\"Scatter plot ({i_.shape[0]} matches)\", weight=\"bold\")\n\nfor i, (label, field, factor, stop) in enumerate(\n (\n (\"Speed radius (km)\", \"radius_s\", 0.001, 80),\n (\"Effective radius (km)\", \"radius_e\", 0.001, 120),\n (\"Amplitude (cm)\", \"amplitude\", 100, 25),\n (\"Maximum Speed (cm/s)\", \"speed_average\", 100, 25),\n )\n):\n ax = fig.add_subplot(\n 2, 2, i + 1, xlabel=\"Filtered grid\", ylabel=\"Without filter\", title=label\n )\n ax.plot(merge_f[field][i_] * factor, merge_t[field][j_] * factor, \".\")\n ax.set_aspect(\"equal\"), ax.grid()\n ax.plot((0, 1000), (0, 1000), \"r\")\n ax.set_xlim(0, stop), ax.set_ylim(0, stop)" ] } ], @@ -176,7 +176,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.7.7" + "version": "3.7.6" } }, "nbformat": 4, diff --git a/src/py_eddy_tracker/dataset/grid.py b/src/py_eddy_tracker/dataset/grid.py index 4a4ebeb8..5d55cf8b 100644 --- a/src/py_eddy_tracker/dataset/grid.py +++ b/src/py_eddy_tracker/dataset/grid.py @@ -304,7 +304,7 @@ def __init__( @property def is_centered(self): - """Give information if pixel is describe with center position or + """Give True if pixel is described with its center's position or a corner :return: True if centered @@ -316,7 +316,7 @@ def is_centered(self): return self.centered def load_general_features(self): - """Load attrs to be store in object + """Load attrs to be stored in object """ logger.debug( "Load general feature from %(filename)s", dict(filename=self.filename) @@ -343,9 +343,9 @@ def load_general_features(self): self.global_attrs = {attr: getattr(h, attr) for attr in h.ncattrs()} def write(self, filename): - """Write dataset output with same format like input + """Write dataset output with same format as input - :param str filename: filename which will be used to save grid + :param str filename: filename used to save the grid """ with Dataset(filename, "w") as h_out: for dimension, size in self.dimensions.items(): @@ -449,7 +449,7 @@ def variables(self): def copy(self, grid_in, grid_out): """ - Duplicate a variable + Duplicate the variable from grid_in in grid_out :param grid_in: :param grid_out: @@ -468,17 +468,17 @@ def add_grid(self, varname, grid): """ Add a grid in handler - :param str varname: name of future grid + :param str varname: name of the future grid :param array grid: grid array """ self.vars[varname] = grid def grid(self, varname, indexs=None): - """give grid required + """Give the grid required :param str varname: Variable to get - :param dict,None indexs: If defined dict must have dimensions name like key - :return: array asked reduce with indexs + :param dict,None indexs: If defined dict must have dimensions name as key + :return: array asked, reduced by the indexes :rtype: array .. minigallery:: py_eddy_tracker.GridDataset.grid @@ -518,7 +518,7 @@ def grid(self, varname, indexs=None): return self.vars[varname] def grid_tiles(self, varname, slice_x, slice_y): - """give grid tiles required, without buffer system + """Give the grid tiles required, without buffer system """ coordinates_dims = list(self.x_dim) coordinates_dims.extend(list(self.y_dim)) @@ -550,20 +550,25 @@ def grid_tiles(self, varname, slice_x, slice_y): return data def high_filter(self, grid_name, w_cut, **kwargs): - """create a high filter with a low one + """Return the grid high-pass filtered, by substracting to the grid the low-pass filter (default: order=1) + :param grid_name: the name of the grid + :param int, w_cut: the half-power wavelength cutoff (km) + """ result = self._low_filter(grid_name, w_cut, **kwargs) self.vars[grid_name] -= result def low_filter(self, grid_name, w_cut, **kwargs): - """low filtering + """Return the grid low-pass filtered (default: order=1) + :param grid_name: the name of the grid + :param int, w_cut: the half-power wavelength cutoff (km) """ result = self._low_filter(grid_name, w_cut, **kwargs) self.vars[grid_name] -= self.vars[grid_name] - result @property def bounds(self): - """Give bound + """Give bounds """ return ( self.x_bounds.min(), @@ -588,21 +593,21 @@ def eddy_identification( **kwargs, ): """ - Compute eddy identification on specified grid + Compute eddy identification on the pecified grid - :param str grid_height: Grid name of height + :param str grid_height: Grid name of Sea Surface Height :param str uname: Grid name of u speed component :param str vname: Grid name of v speed component - :param datetime.datetime date: Date which will be store in object to date data + :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 allow for outter contour in % - :param int sampling: Sampling of contour and speed profile + :param float,int shape_error: Maximal error allowed for outter contour in % + :param int sampling: Number of points to store contours and speed profile :param (int,int),None pixel_limit: - Min and max of pixel which must be inside inner and outer contour to be considered like an eddy - :param float,None precision: Truncate value at the defined precision in m - :param str force_height_unit: Unit to used for height unit - :param str force_speed_unit: Unit to used for speed unit - :param dict kwargs: Argument give to amplitude + Min and max number of pixels inside the inner and the outer 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 :rtype: py_eddy_tracker.observations.observation.EddiesObservations @@ -635,7 +640,7 @@ def eddy_identification( if precision is not None: precision /= factor - # Get h grid + # Get ssh grid data = self.grid(grid_height).astype("f8") # In case of a reduce mask if len(data.mask.shape) == 0 and not data.mask: @@ -686,7 +691,7 @@ def eddy_identification( "contour_lat_s", "uavg_profile", ] - # Compute cyclonic and anticylonic research: + # Complete cyclonic and anticylonic research: a_and_c = list() for anticyclonic_search in [True, False]: eddies = list() @@ -731,7 +736,7 @@ def eddy_identification( contour.reject = 2 continue - # Test to know cyclone or anticyclone + # Test of the rotating sense: cyclone or anticyclone if has_value( data, i_x_in, i_y_in, cvalues, below=anticyclonic_search ): @@ -774,7 +779,7 @@ def eddy_identification( centlon_e = x[centi, centj] centlat_e = y[centi, centj] - # centlat_e and centlon_e must be index of maximum, we will loose some inner contour, if it's not + # centlat_e and centlon_e must be index of maximum, we will loose some inner contour if it's not ( max_average_speed, speed_contour, @@ -792,7 +797,7 @@ def eddy_identification( pixel_min=pixel_limit[0], ) - # FIXME : Instantiate new EddyObservation object (high cost need to be review) + # FIXME : Instantiate new EddyObservation object (high cost need to be reviewed) obs = EddiesObservations( size=1, track_extra_variables=track_extra_variables, @@ -831,7 +836,7 @@ def eddy_identification( centlon_s, centlat_s, eddy_radius_s, aerr_s = _fit_circle_path( create_vertice(*xy_s) ) - # Computed again to use resample contour + # Compute again to use resampled contour _, _, eddy_radius_e, aerr_e = _fit_circle_path( create_vertice(*xy_e) ) @@ -939,7 +944,7 @@ def get_uavg( if not poly_contain_poly(original_contour.vertices, level_contour.vertices): break # 3. Respect size range (for max speed) - # nb_pixel properties need call of pixels_in before with a grid of pixel + # nb_pixel properties need to call pixels_in before with a grid of pixel level_contour.pixels_in(self) # Interpolate uspd to seglon, seglat, then get mean level_average_speed = self.speed_coef_mean(level_contour) @@ -1002,7 +1007,7 @@ def get_amplitude( class UnRegularGridDataset(GridDataset): - """Class which manage unregular grid + """Class managing unregular grid """ __slots__ = ( @@ -1218,7 +1223,7 @@ def ystep(self): return self._y_step def compute_pixel_path(self, x0, y0, x1, y1): - """Give a series of index which describe the path between to position + """Give a series of indexes which describe the path between to position """ return compute_pixel_path( x0, @@ -1238,7 +1243,7 @@ def clean_land(self): pass def is_circular(self): - """Check if grid is circular + """Check if the grid is circular """ if self._is_circular is None: self._is_circular = ( @@ -1464,7 +1469,7 @@ def bessel_high_filter(self, grid_name, wave_length, order=1, lat_max=85, **kwar """ :param str grid_name: grid to filter, data will replace original one :param float wave_length: in km - :param int order: order to use, if > 1 negativ value will be present in kernel + :param int order: order to use, if > 1 negative values of the cardinal sinus are present in kernel :param float lat_max: absolute latitude above no filtering apply :param dict kwargs: look at :py:meth:`RegularGridDataset.convolve_filter_with_dynamic_kernel` @@ -1723,7 +1728,7 @@ def add_uv_lagerloef(self, grid_height, uname="u", vname="v", schema=15): def add_uv(self, grid_height, uname="u", vname="v", stencil_halfwidth=4): """Compute a u and v grid - :param str grid_height: grid name where funtion will apply stencil method + :param str grid_height: grid name where the funtion will apply stencil method :param str uname: future name of u :param str vname: future name of v :param int stencil_halfwidth: largest stencil could be apply @@ -1772,7 +1777,7 @@ def add_uv(self, grid_height, uname="u", vname="v", stencil_halfwidth=4): ) def speed_coef_mean(self, contour): - """some nan can be compute over contour if we are near border, + """Some nan can be computed over contour if we are near border, something to explore """ return mean_on_regular_contour( @@ -1819,11 +1824,11 @@ def contour(self, ax, name, factor=1, **kwargs): def regrid(self, other, grid_name, new_name=None): """ - Interpolate another grid on current grid position + Interpolate another grid at the current grid position :param RegularGridDataset other: - :param str grid_name: var name to interpolate - :param str new_name: name use to store, if None method will use current ont + :param str grid_name: variable name to interpolate + :param str new_name: name used to store, if None method will use current ont .. minigallery:: py_eddy_tracker.RegularGridDataset.regrid """ @@ -1845,7 +1850,7 @@ def interp(self, grid_name, lons, lats): """ Compute z over lons, lats - :param str grid_name: Grid which will be interp + :param str grid_name: Grid to be interpolated :param lons: new x :param lats: new y @@ -1861,7 +1866,7 @@ def interp(self, grid_name, lons, lats): @njit(cache=True, fastmath=True) def compute_pixel_path(x0, y0, x1, y1, x_ori, y_ori, x_step, y_step, nb_x): - """Give a series of index which describe the path between to position + """Give a serie of indexes describing the path between two position """ # index nx = x0.shape[0] diff --git a/src/py_eddy_tracker/observations/observation.py b/src/py_eddy_tracker/observations/observation.py index 5270d9da..4544b859 100644 --- a/src/py_eddy_tracker/observations/observation.py +++ b/src/py_eddy_tracker/observations/observation.py @@ -209,7 +209,7 @@ def hist(self, varname, x, bins, percent=False, mean=False, nb=False): :param array bins: :param bool percent: normalize by sum of all bins :param bool mean: compute mean by bins - :param nb mean: only count by bins + :param bool nb: only count by bins :return: value by bins :rtype: array """ @@ -230,7 +230,7 @@ def box_display(value): def __repr__(self): """ - return general informations on dataset as a string + Return general informations on dataset as a string :return: informations on datasets :rtype: str @@ -311,7 +311,7 @@ def add_rotation_type(self): def circle_contour(self): """ - Set contour like a circle with radius and center data + Set contour as a circle with radius and center data .. minigallery:: py_eddy_tracker.EddiesObservations.circle_contour """ @@ -351,7 +351,7 @@ def dtype(self): @property def elements(self): - """Return all variable name + """Return all the names of the variables """ elements = [i for i in self.ELEMENTS] if self.track_array_variables > 0: @@ -412,7 +412,7 @@ def reset(self): @property def obs(self): - """return an array observations + """Return an array observations """ return self.observations @@ -425,10 +425,10 @@ def __iter__(self): def iter_on(self, xname, bins=None): """ - Yield observation group for each bins + Yield observation group for each bin :param str varname: - :param array bins: bounds og each bins , + :param array bins: bounds of each bin , :return: Group observations :rtype: self.__class__ """ @@ -454,6 +454,9 @@ def iter_on(self, xname, bins=None): yield slice(index[0], index[-1] + 1), bins[i_], bins[i_ + 1] def align_on(self, other, var_name="time", **kwargs): + """ + Align the time indexes of two datasets + """ iter_self, iter_other = ( self.iter_on(var_name, **kwargs), other.iter_on(var_name, **kwargs), @@ -543,6 +546,23 @@ def zarr_dimension(filename): @classmethod def load_file(cls, filename, **kwargs): + """ + Load the netcdf or the zarr file + Load only latitude and longitude on the first 300 obs : + + ```python + kwargs_latlon_300 = dict( + include_vars=[ + "longitude", + "latitude", + ], + indexs=dict(obs=slice(0, 300)), + ) + small_dataset = TrackEddiesObservations.load_file( + nc_file.nc, + **kwargs_latlon_300) + ``` + """ filename_ = ( filename.filename if isinstance(filename, ExFileObject) else filename ) @@ -871,20 +891,23 @@ def intern(flag, public_label=False): return labels def match(self, other, method="overlap", intern=False, cmin=0, **kwargs): - """return index and score compute with area + """Return index and score computed with the chosen method on the effective contour :param EddiesObservations other: Observations to compare :param str method: - if method is "overlap" method will use contour to compute score, - if method is "circle" method will apply a formula of circle overlap - :param bool intern: if True, speed contour will be used - :param float cmin: 0 < cmin < 1, return only couple above cmin + "overlap": the score is computed with contours; + "circle": circles are computed and used for score + :param bool intern: if True, speed contour is used (default = effective contour) + :param float cmin: 0 < cmin < 1, return only couples with score >= cmin :param dict kwargs: look at :py:meth:`vertice_overlap` - :return: return index of couple in self and other and cost value + :return: return the indexes of the eddies in self coupled with eddies in + other and their associated score :rtype: (array(int), array(int), array(float)) .. minigallery:: py_eddy_tracker.EddiesObservations.match """ + # if method is "overlap" method will use contour to compute score, + # if method is "circle" method will apply a formula of circle overlap x_name, y_name = self.intern(intern) if method == "overlap": i, j = bbox_intersection( @@ -902,7 +925,7 @@ def match(self, other, method="overlap", intern=False, cmin=0, **kwargs): self.latitude, self.longitude, other.latitude, other.longitude, **kwargs ) - m = c > cmin + m = c >= cmin # ajout >= pour garder la cmin dans la sélection return i[m], j[m], c[m] @classmethod @@ -947,6 +970,18 @@ def mask_function(self, other, distance): @staticmethod def cost_function(records_in, records_out, distance): + """Return cost function between obs to associate + ```python + CF = sqrt ( [( amplitude_in - amplitude_out)/amplitude_in]^2 + + [( speed_radius_in - speed_radius_out)/speed_radius_in]^2 + + [ distance / 125]^2 ) + ``` + + :param records_in: starting observations + :param records_out: observations to associate + :param distance: computed between in and out + + """ cost = ( (records_in["amplitude"] - records_out["amplitude"]) / records_in["amplitude"] @@ -1052,6 +1087,7 @@ def solve_conflict(cost): @staticmethod def solve_simultaneous(cost): + """Write something""" mask = ~cost.mask # Count number of link by self obs and other obs self_links = mask.sum(axis=1) @@ -1423,9 +1459,12 @@ def set_global_attr_netcdf(self, h_nc): def extract_with_area(self, area, **kwargs): """ - Extract with a bounding box + Extract geographically with a bounding box :param dict area: 4 coordinates in a dictionary to specify bounding box (lower left corner and upper right corner) + ```python + area = dict(llcrnrlon=x0, llcrnrlat=y0, urcrnrlon=x1, urcrnrlat=y1) + ``` :param dict kwargs: look at :py:meth:`extract_with_mask` :return: Return all eddy tracks which are in bounds :rtype: EddiesObservations @@ -1446,6 +1485,12 @@ def extract_with_mask(self, mask): :return: same object with selected observations :rtype: self """ + # ça n'existe plus ça?? + # full_path=False, + # remove_incomplete=False, + # compress_id=False, + # reject_virtual=False + nb_obs = mask.sum() new = self.__class__.new_like(self, nb_obs) new.sign_type = self.sign_type @@ -1460,10 +1505,10 @@ def extract_with_mask(self, mask): def scatter(self, ax, name=None, ref=None, factor=1, **kwargs): """ - :param matplotlib.axes.Axes ax: matplotlib axes use to draw - :param str,None name: - var which will be use to fill contour, if None all element of collection will have same color - :param float,None ref: if define use like west bound + :param matplotlib.axes.Axes ax: matplotlib axe used to draw + :param str, None name: + variable used to fill the contour, if None all elements have the same color + :param float, None ref: if define use like west bound ? :param float factor: multiply value by :param dict kwargs: look at :py:meth:`matplotlib.axes.Axes.scatter` :return: scatter mappable @@ -1492,14 +1537,14 @@ def filled( **kwargs, ): """ - :param matplotlib.axes.Axes ax: matplotlib axes use to draw - :param str,array varname, None: var which will be use to fill contour, or an array of same size of obs - :param float,None ref: if define use like west bound - :param bool intern: if True draw speed contour instead of effective contour + :param matplotlib.axes.Axes ax: matplotlib axe used to draw + :param str,array varname, None: variable used to fill the contours, or an array of same size than obs + :param float,None ref: if define use like west bound? + :param bool intern: if True draw speed contours instead of effective contours :param str cmap: matplotlib colormap name - :param int,None lut: Number of division of colormaps - :param float,None vmin: - :param float,None vmax: + :param int,None lut: Number of colors in the colormap + :param float,None vmin: Min value of the colorbar + :param float,None vmax: Max value of the colorbar :param float factor: multiply value by :return: Collection drawed :rtype: matplotlib.collections.PolyCollection @@ -1535,10 +1580,10 @@ def filled( def bins_stat(self, xname, bins=None, yname=None, method=None): """ - :param str,array xname: var which will be use - :param array, None bins: bins to perform statistics,if None method will get min and max of variable - :param None,str yname: var which will be use to apply method - :param None,str method: If None method count in each bin, could also do mean, std + :param str,array xname: variable to compute stats on + :param array, None bins: bins to perform statistics, if None bins = arange(variable.min(), variable.max() + 2) + :param None,str yname: variable used to apply method + :param None,str method: If None method counts the number of observations in each bin, can be "mean", "std" :return: x array and y array :rtype: array,array @@ -1558,12 +1603,12 @@ def bins_stat(self, xname, bins=None, yname=None, method=None): def display( self, ax, ref=None, extern_only=False, intern_only=False, nobs=True, **kwargs ): - """ - :param matplotlib.axes.Axes ax: matplotlib axes use to draw - :param float,None ref: if define use like west bound - :param bool extern_only: if True draw effective contour only - :param bool intern_only: if True draw speed contour only - :param bool nobs: if True add number of eddies in label + """Plot the speed and effective (dashed) contour of the eddies + :param matplotlib.axes.Axes ax: matplotlib axe used to draw + :param float,None ref: western longitude reference used + :param bool extern_only: if True, draw only the effective contour + :param bool intern_only: if True, draw only the speed contour + :param bool nobs: if True, add the number of eddies in label :param dict kwargs: look at :py:meth:`matplotlib.axes.Axes.plot` .. minigallery:: py_eddy_tracker.EddiesObservations.display @@ -1591,7 +1636,7 @@ def display( def first_obs(self): """ - Get first obs of each tracks. + Get first obs of each trajectory. :rtype: __class__ @@ -1601,7 +1646,7 @@ def first_obs(self): def last_obs(self): """ - Get Last obs of each tracks. + Get Last obs of each trajectory. :rtype: __class__ @@ -1614,10 +1659,10 @@ def last_obs(self): def is_convex(self, intern=False): """ - Get flag of eddy convexity + Get flag of the eddy's convexity - :param bool intern: If true use speed contour instead of effective contour - :return: true if contour is convex + :param bool intern: If True use speed contour instead of effective contour + :return: True if the contour is convex :rtype: array[bool] """ xname, yname = self.intern(intern) @@ -1625,7 +1670,7 @@ def is_convex(self, intern=False): def inside(self, x, y, intern=False): """ - True for each postion inside an eddy + True for each point inside the effective contour of an eddy :param array x: longitude :param array y: latitude @@ -1638,12 +1683,12 @@ def inside(self, x, y, intern=False): def grid_count(self, bins, intern=False, center=False): """ - Compute count of eddies in each bin (use of all pixel in each contour) + Count the eddies in each bin (use all pixels in each contour) - :param (numpy.array,numpy.array) bins: bins to compute count + :param (numpy.array,numpy.array) bins: bins (grid) to count :param bool intern: if True use speed contour only :param bool center: if True use of center to count - :return: return grid of count + :return: return the grid of counts :rtype: py_eddy_tracker.dataset.grid.RegularGridDataset .. minigallery:: py_eddy_tracker.EddiesObservations.grid_count @@ -1662,7 +1707,7 @@ def grid_count(self, bins, intern=False, center=False): lat=(y_bins[1:] + y_bins[:-1]) / 2, ), variables_description=dict( - count=dict(long_name="Number of times pixel is in eddies") + count=dict(long_name="Number of times the pixel is within an eddy") ), centered=True, ) @@ -1690,9 +1735,9 @@ def grid_box_stat(self, bins, varname, method=50, data=None): """ Compute mean of eddies in each bin - :param (numpy.array,numpy.array) bins: bins to compute count - :param str varname: name of variable to compute mean - :param str,float method: method to apply at each set + :param (numpy.array,numpy.array) bins: bins (grid) to count + :param str varname: variable to apply the method + :param str,float method: method to apply. If float, use ? :param array data: Array used to compute stat if defined :return: return grid of method :rtype: py_eddy_tracker.dataset.grid.RegularGridDataset @@ -1733,12 +1778,12 @@ def grid_box_stat(self, bins, varname, method=50, data=None): def grid_stat(self, bins, varname, data=None): """ - Compute mean of eddies in each bin + Return the mean of the eddies' variable in each bin - :param (numpy.array,numpy.array) bins: bins to compute count - :param str varname: name of variable to compute mean + :param (numpy.array,numpy.array) bins: bins (grid) to compute the mean on + :param str varname: name of variable to compute the mean on :param array data: Array used to compute stat if defined - :return: return grid of mean + :return: return the gridde mean variable :rtype: py_eddy_tracker.dataset.grid.RegularGridDataset .. minigallery:: py_eddy_tracker.EddiesObservations.grid_stat @@ -1807,7 +1852,7 @@ def interp_grid( @property def period(self): """ - Give time coverage + Give the time coverage :return: first and last date :rtype: (int,int) @@ -1818,7 +1863,7 @@ def period(self): @njit(cache=True) def grid_count_(grid, i, j): """ - Add one for each indice + Add one to each index """ for i_, j_ in zip(i, j): grid[i_, j_] += 1 @@ -1856,12 +1901,12 @@ def grid_box_stat(x_c, y_c, grid, mask, x, y, value, circular=False, method=50): """ Compute method on each set (one set by box) - :param array_like x_c: longitude coordinate of grid - :param array_like y_c: latitude coordinate of grid - :param array_like grid: grid to store result of method - :param array[bool] mask: grid to store not used box - :param array_like x: longitude of positions - :param array_like y: latitude of positions + :param array_like x_c: grid longitude coordinates + :param array_like y_c: grid latitude coordinates + :param array_like grid: grid to store the result + :param array[bool] mask: grid to store unused boxes + :param array_like x: longitude of observations + :param array_like y: latitude of observations :param array_like value: value to group to apply method :param bool circular: True if grid is wrappable :param float method: percentile @@ -1907,10 +1952,10 @@ def grid_box_stat(x_c, y_c, grid, mask, x, y, value, circular=False, method=50): @njit(cache=True) def grid_stat(x_c, y_c, grid, x, y, result, circular=False, method="mean"): """ - Compute mean of grid for each contour + Compute the mean or the max of the grid for each contour - :param array_like x_c: longitude coordinate of grid - :param array_like y_c: latitude coordinate of grid + :param array_like x_c: the grid longitude coordinates + :param array_like y_c: the grid latitude coordinates :param array_like grid: grid value :param array_like x: longitude of contours :param array_like y: latitude of contours