Skip to content

Commit fa2d382

Browse files
AntSimiAntSimi
committed
Add short example about loopers and eddies
1 parent fc9ab7a commit fa2d382

File tree

4 files changed

+264
-52
lines changed

4 files changed

+264
-52
lines changed

‎examples/12_external_data/pet_drifter_loopers.py‎

Lines changed: 67 additions & 49 deletions
Original file line numberDiff line numberDiff line change
@@ -1,11 +1,9 @@
11
"""
2-
[draft] Loopers vs eddies from altimetry
2+
Colocate looper with eddy from altimetry
33
========================================
44
5-
All loopers datas used in this example are a subset from the dataset describe in this article
5+
All loopers data used in this example are a subset from the dataset describe in this article
66
[Lumpkin, R. : Global characteristics of coherent vortices from surface drifter trajectories](https://doi.org/10.1002/2015JC011435)
7-
8-
97
"""
108

119
import re
@@ -71,65 +69,85 @@ def update_axes(ax, mappable=None):
7169
)
7270

7371
# %%
74-
#
72+
# Global view
73+
# ===========
7574
ax = start_axes("All drifter available in med from Lumpkin dataset")
7675
loopers_med.plot(ax, lw=0.5, color="r", ref=-10)
7776
update_axes(ax)
7877

7978
# %%
80-
# Only long period drifter
81-
long_drifter = loopers_med.extract_with_length((400, -1))
82-
ax = start_axes("Long period drifter")
83-
long_drifter.plot(ax, lw=0.5, color="r", ref=-10)
84-
update_axes(ax)
85-
print(np.unique(long_drifter.tracks))
86-
87-
# %%
88-
drifter_1 = long_drifter.extract_ids((3588,))
89-
fig = plt.figure(figsize=(8, 4))
79+
# One segment of drifter
80+
# ======================
81+
#
82+
# Get a drifter segment (index used have no correspondance with original dataset).
83+
looper = loopers_med.extract_ids((3588,))
84+
fig = plt.figure(figsize=(16, 6))
9085
ax = fig.add_subplot(111, aspect="equal")
86+
looper.plot(ax, lw=0.5, label="Original position of drifter")
87+
looper_filtered = looper.copy()
88+
looper_filtered.position_filter(1, 13)
89+
s = looper_filtered.scatter(
90+
ax,
91+
"time",
92+
cmap=plt.get_cmap("Spectral_r", 20),
93+
label="Filtered position of drifter",
94+
)
95+
plt.colorbar(s).set_label("time (days from 1/1/1950)")
96+
ax.legend()
9197
ax.grid()
92-
drifter_1.plot(ax, lw=0.5)
93-
# %%
94-
drifter_1.close_tracks(
95-
anticyclonic_eddies, method="close_center", delta=0.5, nb_obs_min=25
96-
).tracks
9798

9899
# %%
99-
# Animation which show drifter with colocated eddies
100-
ed = anticyclonic_eddies.extract_ids((4738, 4836, 4910,))
101-
x, y, t = drifter_1.lon, drifter_1.lat, drifter_1.time
102-
103-
104-
def update(frame):
105-
# We display last 5 days of loopers trajectory
106-
m = (t < frame) * (t > (frame - 5))
107-
a.func_animation(frame)
108-
line.set_data(x[m], y[m])
109-
110-
111-
a = Anim(ed, intern=True, figsize=(8, 8), cmap="magma_r", nb_step=10, dpi=60)
112-
line = a.ax.plot([], [], "r", lw=4, zorder=100)[0]
113-
kwargs = dict(frames=np.arange(*a.period, 1), interval=100)
114-
a.fig.suptitle("")
115-
_ = VideoAnimation(a.fig, update, **kwargs)
100+
# Try to find a detected eddies with adt at same place. We used filtered track to simulate an eddy center
101+
match = looper_filtered.close_tracks(
102+
anticyclonic_eddies, method="close_center", delta=0.1, nb_obs_min=50
103+
)
104+
fig = plt.figure(figsize=(16, 6))
105+
ax = fig.add_subplot(111, aspect="equal")
106+
looper.plot(ax, lw=0.5, label="Original position of drifter")
107+
looper_filtered.plot(ax, lw=1.5, label="Filtered position of drifter")
108+
match.plot(ax, lw=1.5, label="Matched eddy")
109+
ax.legend()
110+
ax.grid()
116111

117112
# %%
118-
fig = plt.figure(figsize=(12, 6))
113+
# Display radius of this 2 datasets.
114+
fig = plt.figure(figsize=(20, 8))
119115
ax = fig.add_subplot(111)
120-
ax.plot(drifter_1.time, drifter_1.radius_s / 1e3, label="loopers")
121-
drifter_1.median_filter(7, "time", "radius_s", inplace=True)
122-
drifter_1.loess_filter(15, "time", "radius_s", inplace=True)
116+
ax.plot(looper.time, looper.radius_s / 1e3, label="loopers")
117+
looper_radius = looper.copy()
118+
looper_radius.median_filter(1, "time", "radius_s", inplace=True)
119+
looper_radius.loess_filter(13, "time", "radius_s", inplace=True)
120+
ax.plot(
121+
looper_radius.time,
122+
looper_radius.radius_s / 1e3,
123+
label="loopers (filtered half window 13 days)",
124+
)
125+
ax.plot(match.time, match.radius_s / 1e3, label="altimetry")
126+
match_radius = match.copy()
127+
match_radius.median_filter(1, "time", "radius_s", inplace=True)
128+
match_radius.loess_filter(13, "time", "radius_s", inplace=True)
123129
ax.plot(
124-
drifter_1.time,
125-
drifter_1.radius_s / 1e3,
126-
label="loopers (filtered half window 15 days)",
130+
match_radius.time,
131+
match_radius.radius_s / 1e3,
132+
label="altimetry (filtered half window 13 days)",
127133
)
128-
ax.plot(ed.time, ed.radius_s / 1e3, label="altimetry")
129-
ed.median_filter(7, "time", "radius_s", inplace=True)
130-
ed.loess_filter(15, "time", "radius_s", inplace=True)
131-
ax.plot(ed.time, ed.radius_s / 1e3, label="altimetry (filtered half window 15 days)")
132134
ax.set_ylabel("radius(km)"), ax.set_ylim(0, 100)
133135
ax.legend()
136+
ax.set_title("Radius from loopers and altimeter")
134137
ax.grid()
135-
_ = ax.set_title("Radius from loopers and altimeter")
138+
139+
140+
# %%
141+
# Animation which show drifter with colocated eddy
142+
def update(frame):
143+
# We display last 5 days of loopers trajectory
144+
m = (looper.time < frame) * (looper.time > (frame - 5))
145+
anim.func_animation(frame)
146+
line.set_data(looper.lon[m], looper.lat[m])
147+
148+
149+
anim = Anim(match, intern=True, figsize=(8, 8), cmap="magma_r", nb_step=10, dpi=75)
150+
# mappable to show drifter in red
151+
line = anim.ax.plot([], [], "r", lw=4, zorder=100)[0]
152+
anim.fig.suptitle("")
153+
_ = VideoAnimation(anim.fig, update, frames=np.arange(*anim.period, 1), interval=125)

‎notebooks/python_module/12_external_data/pet_drifter_loopers.ipynb‎

Lines changed: 191 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,191 @@
1+
{
2+
"cells": [
3+
{
4+
"cell_type": "code",
5+
"execution_count": null,
6+
"metadata": {
7+
"collapsed": false
8+
},
9+
"outputs": [],
10+
"source": [
11+
"%matplotlib inline"
12+
]
13+
},
14+
{
15+
"cell_type": "markdown",
16+
"metadata": {},
17+
"source": [
18+
"\nColocate looper with eddy from altimetry\n========================================\n\nAll loopers data used in this example are a subset from the dataset describe in this article\n[Lumpkin, R. : Global characteristics of coherent vortices from surface drifter trajectories](https://doi.org/10.1002/2015JC011435)\n"
19+
]
20+
},
21+
{
22+
"cell_type": "code",
23+
"execution_count": null,
24+
"metadata": {
25+
"collapsed": false
26+
},
27+
"outputs": [],
28+
"source": [
29+
"import re\n\nimport numpy as np\nimport py_eddy_tracker_sample\nfrom matplotlib import pyplot as plt\nfrom matplotlib.animation import FuncAnimation\n\nfrom py_eddy_tracker import data\nfrom py_eddy_tracker.appli.gui import Anim\nfrom py_eddy_tracker.observations.tracking import TrackEddiesObservations"
30+
]
31+
},
32+
{
33+
"cell_type": "code",
34+
"execution_count": null,
35+
"metadata": {
36+
"collapsed": false
37+
},
38+
"outputs": [],
39+
"source": [
40+
"class VideoAnimation(FuncAnimation):\n def _repr_html_(self, *args, **kwargs):\n \"\"\"To get video in html and have a player\"\"\"\n content = self.to_html5_video()\n return re.sub(\n r'width=\"[0-9]*\"\\sheight=\"[0-9]*\"', 'width=\"100%\" height=\"100%\"', content\n )\n\n def save(self, *args, **kwargs):\n if args[0].endswith(\"gif\"):\n # In this case gif is use to create thumbnail which are not use but consume same time than video\n # So we create an empty file, to save time\n with open(args[0], \"w\") as _:\n pass\n return\n return super().save(*args, **kwargs)\n\n\ndef start_axes(title):\n fig = plt.figure(figsize=(13, 5))\n ax = fig.add_axes([0.03, 0.03, 0.90, 0.94], aspect=\"equal\")\n ax.set_xlim(-6, 36.5), ax.set_ylim(30, 46)\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]))"
41+
]
42+
},
43+
{
44+
"cell_type": "markdown",
45+
"metadata": {},
46+
"source": [
47+
"Load eddies dataset\n\n"
48+
]
49+
},
50+
{
51+
"cell_type": "code",
52+
"execution_count": null,
53+
"metadata": {
54+
"collapsed": false
55+
},
56+
"outputs": [],
57+
"source": [
58+
"cyclonic_eddies = TrackEddiesObservations.load_file(\n py_eddy_tracker_sample.get_demo_path(\"eddies_med_adt_allsat_dt2018/Cyclonic.zarr\")\n)\nanticyclonic_eddies = TrackEddiesObservations.load_file(\n py_eddy_tracker_sample.get_demo_path(\n \"eddies_med_adt_allsat_dt2018/Anticyclonic.zarr\"\n )\n)"
59+
]
60+
},
61+
{
62+
"cell_type": "markdown",
63+
"metadata": {},
64+
"source": [
65+
"Load loopers dataset\n\n"
66+
]
67+
},
68+
{
69+
"cell_type": "code",
70+
"execution_count": null,
71+
"metadata": {
72+
"collapsed": false
73+
},
74+
"outputs": [],
75+
"source": [
76+
"loopers_med = TrackEddiesObservations.load_file(\n data.get_demo_path(\"loopers_lumpkin_med.nc\")\n)"
77+
]
78+
},
79+
{
80+
"cell_type": "markdown",
81+
"metadata": {},
82+
"source": [
83+
"Global view\n===========\n\n"
84+
]
85+
},
86+
{
87+
"cell_type": "code",
88+
"execution_count": null,
89+
"metadata": {
90+
"collapsed": false
91+
},
92+
"outputs": [],
93+
"source": [
94+
"ax = start_axes(\"All drifter available in med from Lumpkin dataset\")\nloopers_med.plot(ax, lw=0.5, color=\"r\", ref=-10)\nupdate_axes(ax)"
95+
]
96+
},
97+
{
98+
"cell_type": "markdown",
99+
"metadata": {},
100+
"source": [
101+
"One segment of drifter\n======================\n\nGet a drifter segment (index used have no correspondance with original dataset).\n\n"
102+
]
103+
},
104+
{
105+
"cell_type": "code",
106+
"execution_count": null,
107+
"metadata": {
108+
"collapsed": false
109+
},
110+
"outputs": [],
111+
"source": [
112+
"looper = loopers_med.extract_ids((3588,))\nfig = plt.figure(figsize=(16, 6))\nax = fig.add_subplot(111, aspect=\"equal\")\nlooper.plot(ax, lw=0.5, label=\"Original position of drifter\")\nlooper_filtered = looper.copy()\nlooper_filtered.position_filter(1, 13)\ns = looper_filtered.scatter(\n ax,\n \"time\",\n cmap=plt.get_cmap(\"Spectral_r\", 20),\n label=\"Filtered position of drifter\",\n)\nplt.colorbar(s).set_label(\"time (days from 1/1/1950)\")\nax.legend()\nax.grid()"
113+
]
114+
},
115+
{
116+
"cell_type": "markdown",
117+
"metadata": {},
118+
"source": [
119+
"Try to find a detected eddies with adt at same place. We used filtered track to simulate an eddy center\n\n"
120+
]
121+
},
122+
{
123+
"cell_type": "code",
124+
"execution_count": null,
125+
"metadata": {
126+
"collapsed": false
127+
},
128+
"outputs": [],
129+
"source": [
130+
"match = looper_filtered.close_tracks(\n anticyclonic_eddies, method=\"close_center\", delta=0.1, nb_obs_min=50\n)\nfig = plt.figure(figsize=(16, 6))\nax = fig.add_subplot(111, aspect=\"equal\")\nlooper.plot(ax, lw=0.5, label=\"Original position of drifter\")\nlooper_filtered.plot(ax, lw=1.5, label=\"Filtered position of drifter\")\nmatch.plot(ax, lw=1.5, label=\"Matched eddy\")\nax.legend()\nax.grid()"
131+
]
132+
},
133+
{
134+
"cell_type": "markdown",
135+
"metadata": {},
136+
"source": [
137+
"Display radius of this 2 datasets.\n\n"
138+
]
139+
},
140+
{
141+
"cell_type": "code",
142+
"execution_count": null,
143+
"metadata": {
144+
"collapsed": false
145+
},
146+
"outputs": [],
147+
"source": [
148+
"fig = plt.figure(figsize=(20, 8))\nax = fig.add_subplot(111)\nax.plot(looper.time, looper.radius_s / 1e3, label=\"loopers\")\nlooper_radius = looper.copy()\nlooper_radius.median_filter(1, \"time\", \"radius_s\", inplace=True)\nlooper_radius.loess_filter(13, \"time\", \"radius_s\", inplace=True)\nax.plot(\n looper_radius.time,\n looper_radius.radius_s / 1e3,\n label=\"loopers (filtered half window 13 days)\",\n)\nax.plot(match.time, match.radius_s / 1e3, label=\"altimetry\")\nmatch_radius = match.copy()\nmatch_radius.median_filter(1, \"time\", \"radius_s\", inplace=True)\nmatch_radius.loess_filter(13, \"time\", \"radius_s\", inplace=True)\nax.plot(\n match_radius.time,\n match_radius.radius_s / 1e3,\n label=\"altimetry (filtered half window 13 days)\",\n)\nax.set_ylabel(\"radius(km)\"), ax.set_ylim(0, 100)\nax.legend()\nax.set_title(\"Radius from loopers and altimeter\")\nax.grid()"
149+
]
150+
},
151+
{
152+
"cell_type": "markdown",
153+
"metadata": {},
154+
"source": [
155+
"Animation which show drifter with colocated eddy\n\n"
156+
]
157+
},
158+
{
159+
"cell_type": "code",
160+
"execution_count": null,
161+
"metadata": {
162+
"collapsed": false
163+
},
164+
"outputs": [],
165+
"source": [
166+
"def update(frame):\n # We display last 5 days of loopers trajectory\n m = (looper.time < frame) * (looper.time > (frame - 5))\n anim.func_animation(frame)\n line.set_data(looper.lon[m], looper.lat[m])\n\n\nanim = Anim(match, intern=True, figsize=(8, 8), cmap=\"magma_r\", nb_step=10, dpi=75)\n# mappable to show drifter in red\nline = anim.ax.plot([], [], \"r\", lw=4, zorder=100)[0]\nanim.fig.suptitle(\"\")\n_ = VideoAnimation(anim.fig, update, frames=np.arange(*anim.period, 1), interval=125)"
167+
]
168+
}
169+
],
170+
"metadata": {
171+
"kernelspec": {
172+
"display_name": "Python 3",
173+
"language": "python",
174+
"name": "python3"
175+
},
176+
"language_info": {
177+
"codemirror_mode": {
178+
"name": "ipython",
179+
"version": 3
180+
},
181+
"file_extension": ".py",
182+
"mimetype": "text/x-python",
183+
"name": "python",
184+
"nbconvert_exporter": "python",
185+
"pygments_lexer": "ipython3",
186+
"version": "3.7.7"
187+
}
188+
},
189+
"nbformat": 4,
190+
"nbformat_minor": 0
191+
}

‎src/py_eddy_tracker/observations/network.py‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -1327,7 +1327,7 @@ def extract_light_with_mask(self, mask):
13271327
)
13281328
new.sign_type = self.sign_type
13291329
if nb_obs == 0:
1330-
logger.warning("Empty dataset will be created")
1330+
logger.info("Empty dataset will be created")
13311331
else:
13321332
logger.info(
13331333
f"{nb_obs} observations will be extracted ({nb_obs / self.shape[0]:.3%})"
@@ -1353,7 +1353,7 @@ def extract_with_mask(self, mask):
13531353
new = self.__class__.new_like(self, nb_obs)
13541354
new.sign_type = self.sign_type
13551355
if nb_obs == 0:
1356-
logger.warning("Empty dataset will be created")
1356+
logger.info("Empty dataset will be created")
13571357
else:
13581358
logger.debug(
13591359
f"{nb_obs} observations will be extracted ({nb_obs / self.shape[0]:.3%})"

‎src/py_eddy_tracker/observations/tracking.py‎

Lines changed: 4 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -502,7 +502,7 @@ def extract_with_mask(
502502
new = self.__class__.new_like(self, nb_obs)
503503
new.sign_type = self.sign_type
504504
if nb_obs == 0:
505-
logger.warning("Empty dataset will be created")
505+
logger.info("Empty dataset will be created")
506506
else:
507507
for field in self.obs.dtype.descr:
508508
logger.debug("Copy of field %s ...", field)
@@ -568,6 +568,9 @@ def close_tracks(self, other, nb_obs_min=10, **kwargs):
568568
It could be a costly operation for huge dataset
569569
"""
570570
p0, p1 = self.period
571+
p0_other, p1_other = other.period
572+
if p1_other < p0 or p1 < p0_other:
573+
return other.__class__.new_like(other, 0)
571574
indexs = list()
572575
for i_self, i_other, t0, t1 in self.align_on(other, bins=arange(p0, p1 + 2)):
573576
i, j, s = self.match(other, i_self=i_self, i_other=i_other, **kwargs)

0 commit comments

Comments
 (0)