Skip to content

Commit 3b38417

Browse files
committed
Add method to get first and last
1 parent 50979e0 commit 3b38417

File tree

6 files changed

+310
-33
lines changed

6 files changed

+310
-33
lines changed
Lines changed: 81 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,81 @@
1+
"""
2+
Birth and death
3+
===============
4+
5+
Following figures are based on https://doi.org/10.1016/j.pocean.2011.01.002
6+
7+
"""
8+
from matplotlib import pyplot as plt
9+
from py_eddy_tracker.observations.tracking import TrackEddiesObservations
10+
import py_eddy_tracker_sample
11+
12+
13+
# %%
14+
def start_axes(title):
15+
fig = plt.figure(figsize=(13, 5))
16+
ax = fig.add_axes([0.03, 0.03, 0.90, 0.94])
17+
ax.set_xlim(-6, 36.5), ax.set_ylim(30, 46)
18+
ax.set_aspect("equal")
19+
ax.set_title(title)
20+
return ax
21+
22+
23+
def update_axes(ax, mappable=None):
24+
ax.grid()
25+
if mappable:
26+
plt.colorbar(mappable, cax=ax.figure.add_axes([0.95, 0.05, 0.01, 0.9]))
27+
28+
29+
# %%
30+
# Load an experimental med atlas over a period of 26 years (1993-2019)
31+
kwargs_load = dict(
32+
include_vars=(
33+
"longitude",
34+
"latitude",
35+
"observation_number",
36+
"track",
37+
"time",
38+
"speed_contour_longitude",
39+
"speed_contour_latitude",
40+
)
41+
)
42+
a = TrackEddiesObservations.load_file(
43+
py_eddy_tracker_sample.get_path("eddies_med_adt_allsat_dt2018/Anticyclonic.zarr")
44+
)
45+
c = TrackEddiesObservations.load_file(
46+
py_eddy_tracker_sample.get_path("eddies_med_adt_allsat_dt2018/Cyclonic.zarr")
47+
)
48+
49+
# %%
50+
t0, t1 = a.period
51+
step = 0.125
52+
bins = ((-10, 37, step), (30, 46, step))
53+
kwargs = dict(cmap="terrain_r", factor=100 / (t1 - t0), name="count", vmin=0, vmax=1)
54+
55+
# %%
56+
# Cyclonic
57+
# --------
58+
ax = start_axes("Birth cyclonic frenquency (%)")
59+
g_c_first = c.first_obs().grid_count(bins, intern=True)
60+
m = g_c_first.display(ax, **kwargs)
61+
update_axes(ax, m)
62+
63+
# %%
64+
ax = start_axes("Death cyclonic frenquency (%)")
65+
g_c_last = c.last_obs().grid_count(bins, intern=True)
66+
m = g_c_last.display(ax, **kwargs)
67+
update_axes(ax, m)
68+
69+
# %%
70+
# Anticyclonic
71+
# ------------
72+
ax = start_axes("Birth anticyclonic frequency (%)")
73+
g_a_first = a.first_obs().grid_count(bins, intern=True)
74+
m = g_a_first.display(ax, **kwargs)
75+
update_axes(ax, m)
76+
77+
# %%
78+
ax = start_axes("Death anticyclonic frequency (%)")
79+
g_a_last = a.last_obs().grid_count(bins, intern=True)
80+
m = g_a_last.display(ax, **kwargs)
81+
update_axes(ax, m)

examples/10_tracking_diagnostics/pet_center_count.py

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,6 @@
11
"""
22
Count center
3-
======================
3+
============
44
55
"""
66
from matplotlib import pyplot as plt
Lines changed: 152 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,152 @@
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+
"\n# Birth and death\n\nFollowing figures are based on https://doi.org/10.1016/j.pocean.2011.01.002\n"
19+
]
20+
},
21+
{
22+
"cell_type": "code",
23+
"execution_count": null,
24+
"metadata": {
25+
"collapsed": false
26+
},
27+
"outputs": [],
28+
"source": [
29+
"from matplotlib import pyplot as plt\nfrom py_eddy_tracker.observations.tracking import TrackEddiesObservations\nimport py_eddy_tracker_sample"
30+
]
31+
},
32+
{
33+
"cell_type": "code",
34+
"execution_count": null,
35+
"metadata": {
36+
"collapsed": false
37+
},
38+
"outputs": [],
39+
"source": [
40+
"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]))"
41+
]
42+
},
43+
{
44+
"cell_type": "markdown",
45+
"metadata": {},
46+
"source": [
47+
"Load an experimental med atlas over a period of 26 years (1993-2019)\n\n"
48+
]
49+
},
50+
{
51+
"cell_type": "code",
52+
"execution_count": null,
53+
"metadata": {
54+
"collapsed": false
55+
},
56+
"outputs": [],
57+
"source": [
58+
"kwargs_load = dict(\n include_vars=(\n \"longitude\",\n \"latitude\",\n \"observation_number\",\n \"track\",\n \"time\",\n \"speed_contour_longitude\",\n \"speed_contour_latitude\",\n )\n)\na = TrackEddiesObservations.load_file(\n py_eddy_tracker_sample.get_path(\"eddies_med_adt_allsat_dt2018/Anticyclonic.zarr\")\n)\nc = TrackEddiesObservations.load_file(\n py_eddy_tracker_sample.get_path(\"eddies_med_adt_allsat_dt2018/Cyclonic.zarr\")\n)"
59+
]
60+
},
61+
{
62+
"cell_type": "code",
63+
"execution_count": null,
64+
"metadata": {
65+
"collapsed": false
66+
},
67+
"outputs": [],
68+
"source": [
69+
"t0, t1 = a.period\nstep = 0.125\nbins = ((-10, 37, step), (30, 46, step))\nkwargs = dict(cmap=\"terrain_r\", factor=100 / (t1 - t0), name=\"count\", vmin=0, vmax=1)"
70+
]
71+
},
72+
{
73+
"cell_type": "markdown",
74+
"metadata": {},
75+
"source": [
76+
"## Cyclonic\n\n"
77+
]
78+
},
79+
{
80+
"cell_type": "code",
81+
"execution_count": null,
82+
"metadata": {
83+
"collapsed": false
84+
},
85+
"outputs": [],
86+
"source": [
87+
"ax = start_axes(\"Birth cyclonic frenquency (%)\")\ng_c_first = c.first_obs().grid_count(bins, intern=True)\nm = g_c_first.display(ax, **kwargs)\nupdate_axes(ax, m)"
88+
]
89+
},
90+
{
91+
"cell_type": "code",
92+
"execution_count": null,
93+
"metadata": {
94+
"collapsed": false
95+
},
96+
"outputs": [],
97+
"source": [
98+
"ax = start_axes(\"Death cyclonic frenquency (%)\")\ng_c_last = c.last_obs().grid_count(bins, intern=True)\nm = g_c_last.display(ax, **kwargs)\nupdate_axes(ax, m)"
99+
]
100+
},
101+
{
102+
"cell_type": "markdown",
103+
"metadata": {},
104+
"source": [
105+
"## Anticyclonic\n\n"
106+
]
107+
},
108+
{
109+
"cell_type": "code",
110+
"execution_count": null,
111+
"metadata": {
112+
"collapsed": false
113+
},
114+
"outputs": [],
115+
"source": [
116+
"ax = start_axes(\"Birth anticyclonic frequency (%)\")\ng_a_first = a.first_obs().grid_count(bins, intern=True)\nm = g_a_first.display(ax, **kwargs)\nupdate_axes(ax, m)"
117+
]
118+
},
119+
{
120+
"cell_type": "code",
121+
"execution_count": null,
122+
"metadata": {
123+
"collapsed": false
124+
},
125+
"outputs": [],
126+
"source": [
127+
"ax = start_axes(\"Death anticyclonic frequency (%)\")\ng_a_last = a.last_obs().grid_count(bins, intern=True)\nm = g_a_last.display(ax, **kwargs)\nupdate_axes(ax, m)"
128+
]
129+
}
130+
],
131+
"metadata": {
132+
"kernelspec": {
133+
"display_name": "Python 3",
134+
"language": "python",
135+
"name": "python3"
136+
},
137+
"language_info": {
138+
"codemirror_mode": {
139+
"name": "ipython",
140+
"version": 3
141+
},
142+
"file_extension": ".py",
143+
"mimetype": "text/x-python",
144+
"name": "python",
145+
"nbconvert_exporter": "python",
146+
"pygments_lexer": "ipython3",
147+
"version": "3.7.7"
148+
}
149+
},
150+
"nbformat": 4,
151+
"nbformat_minor": 0
152+
}

notebooks/python_module/10_tracking_diagnostics/pet_center_count.ipynb

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -15,7 +15,7 @@
1515
"cell_type": "markdown",
1616
"metadata": {},
1717
"source": [
18-
"\nCount center\n======================\n"
18+
"\n# Count center\n"
1919
]
2020
},
2121
{
@@ -62,7 +62,7 @@
6262
},
6363
"outputs": [],
6464
"source": [
65-
"fig = plt.figure(figsize=(12, 18.5))\nax_a = fig.add_axes([0.03, 0.75, 0.90, 0.25])\nax_a.set_title(\"Anticyclonic center frequency\")\nax_c = fig.add_axes([0.03, 0.5, 0.90, 0.25])\nax_c.set_title(\"Cyclonic center frequency\")\nax_all = fig.add_axes([0.03, 0.25, 0.90, 0.25])\nax_all.set_title(\"All eddies center frequency\")\nax_ratio = fig.add_axes([0.03, 0.0, 0.90, 0.25])\nax_ratio.set_title(\"Ratio cyclonic / Anticyclonic\")\n\n# Count pixel used for each center\ng_a = a.grid_count(bins, intern=True, center=True)\ng_a.display(ax_a, **kwargs_pcolormesh)\ng_c = c.grid_count(bins, intern=True, center=True)\ng_c.display(ax_c, **kwargs_pcolormesh)\n# Compute a ratio Cyclonic / Anticyclonic\nratio = g_c.vars[\"count\"] / g_a.vars[\"count\"]\n\n# Mask manipulation to be able to sum the 2 grids\nm_c = g_c.vars[\"count\"].mask\nm = m_c & g_a.vars[\"count\"].mask\ng_c.vars[\"count\"][m_c] = 0\ng_c.vars[\"count\"] += g_a.vars[\"count\"]\ng_c.vars[\"count\"].mask = m\n\nm = g_c.display(ax_all, **kwargs_pcolormesh)\ncb = plt.colorbar(m, cax=fig.add_axes([0.94, 0.27, 0.01, 0.7]))\ncb.set_label(\"Eddies by 1\u00b0^2 by day\")\n\ng_c.vars[\"count\"] = ratio\nm = g_c.display(ax_ratio, name=\"count\", vmin=0.1, vmax=10, norm=LogNorm())\nplt.colorbar(m, cax=fig.add_axes([0.94, 0.02, 0.01, 0.2]))\n\nfor ax in (ax_a, ax_c, ax_all, ax_ratio):\n ax.set_aspect(\"equal\")\n ax.set_xlim(-6, 36.5), ax.set_ylim(30, 46)\n ax.grid()"
65+
"fig = plt.figure(figsize=(12, 18.5))\nax_a = fig.add_axes([0.03, 0.75, 0.90, 0.25])\nax_a.set_title(\"Anticyclonic center frequency\")\nax_c = fig.add_axes([0.03, 0.5, 0.90, 0.25])\nax_c.set_title(\"Cyclonic center frequency\")\nax_all = fig.add_axes([0.03, 0.25, 0.90, 0.25])\nax_all.set_title(\"All eddies center frequency\")\nax_ratio = fig.add_axes([0.03, 0.0, 0.90, 0.25])\nax_ratio.set_title(\"Ratio cyclonic / Anticyclonic\")\n\n# Count pixel used for each center\ng_a = a.grid_count(bins, intern=True, center=True)\ng_a.display(ax_a, **kwargs_pcolormesh)\ng_c = c.grid_count(bins, intern=True, center=True)\ng_c.display(ax_c, **kwargs_pcolormesh)\n# Compute a ratio Cyclonic / Anticyclonic\nratio = g_c.vars[\"count\"] / g_a.vars[\"count\"]\n\n# Mask manipulation to be able to sum the 2 grids\nm_c = g_c.vars[\"count\"].mask\nm = m_c & g_a.vars[\"count\"].mask\ng_c.vars[\"count\"][m_c] = 0\ng_c.vars[\"count\"] += g_a.vars[\"count\"]\ng_c.vars[\"count\"].mask = m\n\nm = g_c.display(ax_all, **kwargs_pcolormesh)\ncb = plt.colorbar(m, cax=fig.add_axes([0.94, 0.27, 0.01, 0.7]))\ncb.set_label(\"Eddies by 1\u00b0^2 by day\")\n\ng_c.vars[\"count\"] = ratio\nm = g_c.display(ax_ratio, name=\"count\", vmin=0.1, vmax=10, norm=LogNorm(), cmap='coolwarm_r')\nplt.colorbar(m, cax=fig.add_axes([0.94, 0.02, 0.01, 0.2]))\n\nfor ax in (ax_a, ax_c, ax_all, ax_ratio):\n ax.set_aspect(\"equal\")\n ax.set_xlim(-6, 36.5), ax.set_ylim(30, 46)\n ax.grid()"
6666
]
6767
}
6868
],

src/py_eddy_tracker/observations/observation.py

Lines changed: 60 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1353,6 +1353,43 @@ def set_global_attr_netcdf(self, h_nc):
13531353
for key, item in self.global_attr.items():
13541354
h_nc.setncattr(key, item)
13551355

1356+
def extract_with_area(self, area, **kwargs):
1357+
"""
1358+
Extract with a bounding box
1359+
1360+
:param dict area: 4 coordinates in a dictionary to specify bounding box (lower left corner and upper right corner)
1361+
:param dict kwargs: look at :py:meth:`extract_with_mask`
1362+
:return: Return all eddy tracks which are in bounds
1363+
:rtype: EddiesObservations
1364+
1365+
.. minigallery:: py_eddy_tracker.EddiesObservations.extract_with_area
1366+
"""
1367+
mask = (self.latitude > area["llcrnrlat"]) * (self.latitude < area["urcrnrlat"])
1368+
lon0 = area["llcrnrlon"]
1369+
lon = (self.longitude - lon0) % 360 + lon0
1370+
mask *= (lon > lon0) * (lon < area["urcrnrlon"])
1371+
return self.extract_with_mask(mask, **kwargs)
1372+
1373+
def extract_with_mask(self, mask):
1374+
"""
1375+
Extract a subset of observations
1376+
1377+
:param array(bool) mask: mask to select observations
1378+
:return: same object with selected observations
1379+
:rtype: self
1380+
"""
1381+
nb_obs = mask.sum()
1382+
new = self.__class__.new_like(self, nb_obs)
1383+
new.sign_type = self.sign_type
1384+
if nb_obs == 0:
1385+
logger.warning("Empty dataset will be created")
1386+
else:
1387+
for field in self.obs.dtype.descr:
1388+
logger.debug("Copy of field %s ...", field)
1389+
var = field[0]
1390+
new.obs[var] = self.obs[var][mask]
1391+
return new
1392+
13561393
def scatter(self, ax, name, ref=None, factor=1, **kwargs):
13571394
"""
13581395
:param matplotlib.axes.Axes ax: matplotlib axes use to draw
@@ -1449,6 +1486,29 @@ def display(
14491486
lon_e, lat_e = wrap_longitude(lon_e, lat_e, ref, cut=True)
14501487
ax.plot(lon_e, lat_e, linestyle="-.", **kwargs_e)
14511488

1489+
def first_obs(self):
1490+
"""
1491+
Get first obs of each tracks.
1492+
1493+
:rtype: __class__
1494+
1495+
.. minigallery:: py_eddy_tracker.EddiesObservations.first_obs
1496+
"""
1497+
return self.extract_with_mask(self["n"] == 0)
1498+
1499+
def last_obs(self):
1500+
"""
1501+
Get Last obs of each tracks.
1502+
1503+
:rtype: __class__
1504+
1505+
.. minigallery:: py_eddy_tracker.EddiesObservations.last_obs
1506+
"""
1507+
m = zeros(len(self), dtype="bool")
1508+
m[-1] = True
1509+
m[:-1][self["n"][1:] == 0] = True
1510+
return self.extract_with_mask(m)
1511+
14521512
def grid_count(self, bins, intern=False, center=False):
14531513
"""
14541514
Compute count of eddies in each bin (use of all pixel in each contour)

0 commit comments

Comments
 (0)