|
15 | 15 | "cell_type": "markdown",
|
16 | 16 | "metadata": {},
|
17 | 17 | "source": [
|
18 |
| - "\n# Eddy detection : Antartic circum polar\n\nScript will detect eddies on adt field, and compute u,v with method add_uv(which could use, only if equator is avoid)\n\nTwo ones with filtering adt and another without\n" |
| 18 | + "\n# Eddy detection : Antartic Circumpolar Current\n\nThis script detect eddies on the ADT field, and compute u,v with the method add_uv (use it only if the Equator is avoided)\n\nTwo detections are provided : with a filtered ADT and without filtering\n" |
19 | 19 | ]
|
20 | 20 | },
|
21 | 21 | {
|
|
26 | 26 | },
|
27 | 27 | "outputs": [],
|
28 | 28 | "source": [
|
29 |
| - "from datetime import datetime\n\nfrom matplotlib import pyplot as plt\n\nfrom py_eddy_tracker import data\nfrom py_eddy_tracker.dataset.grid import RegularGridDataset\n\n\ndef quad_axes(title):\n fig = plt.figure(figsize=(13, 8.5))\n fig.suptitle(title, weight=\"bold\")\n axes = list()\n for position in (\n [0.05, 0.53, 0.44, 0.44],\n [0.53, 0.53, 0.44, 0.44],\n [0.05, 0.03, 0.44, 0.44],\n [0.53, 0.03, 0.44, 0.44],\n ):\n ax = fig.add_axes(position)\n ax.set_xlim(5, 45), ax.set_ylim(-60, -37)\n ax.set_aspect(\"equal\"), ax.grid(True)\n axes.append(ax)\n return axes" |
| 29 | + "from datetime import datetime\n\nfrom matplotlib import pyplot as plt\nfrom matplotlib import style\n\nfrom py_eddy_tracker import data\nfrom py_eddy_tracker.dataset.grid import RegularGridDataset\n\npos_cb = [0.1, 0.52, 0.83, 0.015]\npos_cb2 = [0.1, 0.07, 0.4, 0.015]\n\n\ndef quad_axes(title):\n style.use(\"default\")\n fig = plt.figure(figsize=(13, 10))\n fig.suptitle(title, weight=\"bold\", fontsize=14)\n axes = list()\n\n ax_pos = dict(\n topleft=[0.1, 0.54, 0.4, 0.38],\n topright=[0.53, 0.54, 0.4, 0.38],\n botleft=[0.1, 0.09, 0.4, 0.38],\n botright=[0.53, 0.09, 0.4, 0.38],\n )\n\n for key, position in ax_pos.items():\n ax = fig.add_axes(position)\n ax.set_xlim(5, 45), ax.set_ylim(-60, -37)\n ax.set_aspect(\"equal\"), ax.grid(True)\n axes.append(ax)\n if \"right\" in key:\n ax.set_yticklabels(\"\")\n return fig, axes\n\n\ndef set_fancy_labels(fig, ticklabelsize=14, labelsize=14, labelweight=\"semibold\"):\n for ax in fig.get_axes():\n ax.grid()\n ax.grid(which=\"major\", linestyle=\"-\", linewidth=\"0.5\", color=\"black\")\n if ax.get_ylabel() != \"\":\n ax.set_ylabel(ax.get_ylabel(), fontsize=labelsize, fontweight=labelweight)\n if ax.get_xlabel() != \"\":\n ax.set_xlabel(ax.get_xlabel(), fontsize=labelsize, fontweight=labelweight)\n if ax.get_title() != \"\":\n ax.set_title(ax.get_title(), fontsize=labelsize, fontweight=labelweight)\n ax.tick_params(labelsize=ticklabelsize)" |
30 | 30 | ]
|
31 | 31 | },
|
32 | 32 | {
|
|
80 | 80 | },
|
81 | 81 | "outputs": [],
|
82 | 82 | "source": [
|
83 |
| - "axs = quad_axes(\"General properties field\")\nm = g_raw.display(axs[0], \"adt\", vmin=-1, vmax=1, cmap=\"RdBu_r\")\naxs[0].set_title(\"ADT(m)\")\nm = g.display(axs[1], \"adt_low\", vmin=-1, vmax=1, cmap=\"RdBu_r\")\naxs[1].set_title(\"ADT (m) large scale with cut at 700 km\")\nm = g.display(axs[2], \"adt\", vmin=-1, vmax=1, cmap=\"RdBu_r\")\naxs[2].set_title(\"ADT (m) high scale with cut at 700 km\")\ncb = plt.colorbar(\n m, cax=axs[0].figure.add_axes([0.03, 0.51, 0.94, 0.01]), orientation=\"horizontal\"\n)\ncb.set_label(\"ADT(m)\", labelpad=-2)" |
| 83 | + "kw_adt = dict(vmin=-1.5, vmax=1.5, cmap=plt.get_cmap(\"RdBu_r\", 30))\nfig, axs = quad_axes(\"General properties field\")\nm = g_raw.display(axs[0], \"adt\", **kw_adt)\naxs[0].set_title(\"Total ADT (m)\")\nm = g.display(axs[1], \"adt_low\", **kw_adt)\naxs[1].set_title(\"ADT (m) large scale, cutoff at 700 km\")\nm2 = g.display(axs[2], \"adt\", cmap=plt.get_cmap(\"RdBu_r\", 20), vmin=-0.5, vmax=0.5)\naxs[2].set_title(\"ADT (m) high-pass filtered, a cutoff at 700 km\")\ncb = plt.colorbar(m, cax=axs[0].figure.add_axes(pos_cb), orientation=\"horizontal\")\ncb.set_label(\"ADT (m)\", labelpad=0)\ncb2 = plt.colorbar(m2, cax=axs[2].figure.add_axes(pos_cb2), orientation=\"horizontal\")\ncb2.set_label(\"ADT (m)\", labelpad=0)\nset_fancy_labels(fig)" |
| 84 | + ] |
| 85 | + }, |
| 86 | + { |
| 87 | + "cell_type": "markdown", |
| 88 | + "metadata": {}, |
| 89 | + "source": [ |
| 90 | + "The large-scale North-South gradient is removed by the filtering step.\n\n" |
84 | 91 | ]
|
85 | 92 | },
|
86 | 93 | {
|
|
91 | 98 | },
|
92 | 99 | "outputs": [],
|
93 | 100 | "source": [
|
94 |
| - "axs = quad_axes(\"\")\naxs[0].set_title(\"Without filter\")\naxs[0].set_ylabel(\"Contours used in eddies\")\naxs[1].set_title(\"With filter\")\naxs[2].set_ylabel(\"Closed contours but not used\")\ng_raw.contours.display(axs[0], lw=0.5, only_used=True)\ng.contours.display(axs[1], lw=0.5, only_used=True)\ng_raw.contours.display(axs[2], lw=0.5, only_unused=True)\ng.contours.display(axs[3], lw=0.5, only_unused=True)" |
| 101 | + "fig, axs = quad_axes(\"\")\naxs[0].set_title(\"Without filter\")\naxs[0].set_ylabel(\"Contours used in eddies\")\naxs[1].set_title(\"With filter\")\naxs[2].set_ylabel(\"Closed contours but not used\")\ng_raw.contours.display(axs[0], lw=0.5, only_used=True)\ng.contours.display(axs[1], lw=0.5, only_used=True)\ng_raw.contours.display(axs[2], lw=0.5, only_unused=True)\ng.contours.display(axs[3], lw=0.5, only_unused=True)\nset_fancy_labels(fig)" |
| 102 | + ] |
| 103 | + }, |
| 104 | + { |
| 105 | + "cell_type": "markdown", |
| 106 | + "metadata": {}, |
| 107 | + "source": [ |
| 108 | + "Removing the large-scale North-South gradient reveals closed contours in the\nSouth-Western corner of the ewample region.\n\n" |
95 | 109 | ]
|
96 | 110 | },
|
97 | 111 | {
|
|
102 | 116 | },
|
103 | 117 | "outputs": [],
|
104 | 118 | "source": [
|
105 |
| - "kw = dict(ref=-10, linewidth=0.75)\nkw_a = dict(color=\"r\", label=\"Anticyclonic ({nb_obs} eddies)\")\nkw_c = dict(color=\"b\", label=\"Cyclonic ({nb_obs} eddies)\")\nkw_filled = dict(vmin=0, vmax=100, cmap=\"Spectral_r\", lut=20, intern=True, factor=100)\naxs = quad_axes(\"Comparison between two detection\")\n# Match with intern/inner contour\ni_a, j_a, s_a = a_.match(a, intern=True, cmin=0.15)\ni_c, j_c, s_c = c_.match(c, intern=True, cmin=0.15)\n\na_.index(i_a).filled(axs[0], s_a, **kw_filled)\na.index(j_a).filled(axs[1], s_a, **kw_filled)\nc_.index(i_c).filled(axs[0], s_c, **kw_filled)\nm = c.index(j_c).filled(axs[1], s_c, **kw_filled)\n\ncb = plt.colorbar(\n m, cax=axs[0].figure.add_axes([0.03, 0.51, 0.94, 0.01]), orientation=\"horizontal\"\n)\ncb.set_label(\"Similarity index\", labelpad=-5)\na_.display(axs[0], **kw, **kw_a), c_.display(axs[0], **kw, **kw_c)\na.display(axs[1], **kw, **kw_a), c.display(axs[1], **kw, **kw_c)\n\naxs[0].set_title(\"Without filter\")\naxs[0].set_ylabel(\"Detection\")\naxs[1].set_title(\"With filter\")\naxs[2].set_ylabel(\"Contours' rejection criteria\")\n\ng_raw.contours.display(axs[2], lw=0.5, only_unused=True, display_criterion=True)\ng.contours.display(axs[3], lw=0.5, only_unused=True, display_criterion=True)\n\nfor ax in axs:\n ax.legend()" |
| 119 | + "kw = dict(ref=-10, linewidth=0.75)\nkw_a = dict(color=\"r\", label=\"Anticyclonic ({nb_obs} eddies)\")\nkw_c = dict(color=\"b\", label=\"Cyclonic ({nb_obs} eddies)\")\nkw_filled = dict(vmin=0, vmax=100, cmap=\"Spectral_r\", lut=20, intern=True, factor=100)\nfig, axs = quad_axes(\"Comparison between two detections\")\n# Match with intern/inner contour\ni_a, j_a, s_a = a_.match(a, intern=True, cmin=0.15)\ni_c, j_c, s_c = c_.match(c, intern=True, cmin=0.15)\n\na_.index(i_a).filled(axs[0], s_a, **kw_filled)\na.index(j_a).filled(axs[1], s_a, **kw_filled)\nc_.index(i_c).filled(axs[0], s_c, **kw_filled)\nm = c.index(j_c).filled(axs[1], s_c, **kw_filled)\n\ncb = plt.colorbar(m, cax=axs[0].figure.add_axes(pos_cb), orientation=\"horizontal\")\ncb.set_label(\"Similarity index (%)\", labelpad=-5)\na_.display(axs[0], **kw, **kw_a), c_.display(axs[0], **kw, **kw_c)\na.display(axs[1], **kw, **kw_a), c.display(axs[1], **kw, **kw_c)\n\naxs[0].set_title(\"Without filter\")\naxs[0].set_ylabel(\"Detection\")\naxs[1].set_title(\"With filter\")\naxs[2].set_ylabel(\"Contours' rejection criteria\")\n\ng_raw.contours.display(axs[2], lw=0.5, only_unused=True, display_criterion=True)\ng.contours.display(axs[3], lw=0.5, only_unused=True, display_criterion=True)\n\nfor ax in axs:\n ax.legend()\n\nset_fancy_labels(fig)" |
| 120 | + ] |
| 121 | + }, |
| 122 | + { |
| 123 | + "cell_type": "markdown", |
| 124 | + "metadata": {}, |
| 125 | + "source": [ |
| 126 | + "Very similar eddies have Similarity Indexes >= 40%\n\n" |
106 | 127 | ]
|
107 | 128 | },
|
108 | 129 | {
|
|
120 | 141 | },
|
121 | 142 | "outputs": [],
|
122 | 143 | "source": [
|
123 |
| - "i_a, j_a = i_a[s_a >= 0.4], j_a[s_a >= 0.4]\ni_c, j_c = i_c[s_c >= 0.4], j_c[s_c >= 0.4]\nfig = plt.figure(figsize=(12, 12))\nfig.suptitle(f\"Scatter plot (A : {i_a.shape[0]}, C : {i_c.shape[0]} matches)\")\n\nfor i, (label, field, factor, stop) in enumerate(\n (\n (\"speed radius (km)\", \"radius_s\", 0.001, 120),\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(2, 2, i + 1, title=label)\n ax.set_xlabel(\"Without filter\")\n ax.set_ylabel(\"With filter\")\n\n ax.plot(\n a_[field][i_a] * factor,\n a[field][j_a] * factor,\n \"r.\",\n label=\"Anticyclonic\",\n )\n ax.plot(\n c_[field][i_c] * factor,\n c[field][j_c] * factor,\n \"b.\",\n label=\"Cyclonic\",\n )\n ax.set_aspect(\"equal\"), ax.grid()\n ax.plot((0, 1000), (0, 1000), \"g\")\n ax.set_xlim(0, stop), ax.set_ylim(0, stop)\n ax.legend()" |
| 144 | + "i_a, j_a = i_a[s_a >= 0.4], j_a[s_a >= 0.4]\ni_c, j_c = i_c[s_c >= 0.4], j_c[s_c >= 0.4]\nfig = plt.figure(figsize=(12, 12))\nfig.suptitle(f\"Scatter plot (A : {i_a.shape[0]}, C : {i_c.shape[0]} matches)\")\n\nfor i, (label, field, factor, stop) in enumerate(\n (\n (\"Speed radius (km)\", \"radius_s\", 0.001, 120),\n (\"Effective 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(2, 2, i + 1, title=label)\n ax.set_xlabel(\"Without filter\")\n ax.set_ylabel(\"With filter\")\n\n ax.plot(\n a_[field][i_a] * factor,\n a[field][j_a] * factor,\n \"r.\",\n label=\"Anticyclonic\",\n )\n ax.plot(\n c_[field][i_c] * factor,\n c[field][j_c] * factor,\n \"b.\",\n label=\"Cyclonic\",\n )\n ax.set_aspect(\"equal\"), ax.grid()\n ax.plot((0, 1000), (0, 1000), \"g\")\n ax.set_xlim(0, stop), ax.set_ylim(0, stop)\n ax.legend()\n\nset_fancy_labels(fig)" |
124 | 145 | ]
|
125 | 146 | }
|
126 | 147 | ],
|
|
0 commit comments