Skip to content

Commit 1b2a824

Browse files
AntSimiAntSimi
committed
Add example from Evan Mason about normalized lifetime
1 parent 8d13115 commit 1b2a824

File tree

4 files changed

+199
-1
lines changed

4 files changed

+199
-1
lines changed

‎examples/10_tracking_diagnostics/pet_normalised_lifetime.py‎

Lines changed: 78 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,78 @@
1+
"""
2+
Normalised Eddy Lifetimes
3+
=========================
4+
5+
Example from Evan Mason
6+
"""
7+
from matplotlib import pyplot as plt
8+
from numba import njit
9+
from numpy import interp, linspace, zeros
10+
from py_eddy_tracker_sample import get_demo_path
11+
12+
from py_eddy_tracker.observations.tracking import TrackEddiesObservations
13+
14+
15+
# %%
16+
@njit(cache=True)
17+
def sum_profile(x_new, y, out):
18+
"""Will sum all interpolated given array"""
19+
out += interp(x_new, linspace(0, 1, y.size), y)
20+
21+
22+
class MyObs(TrackEddiesObservations):
23+
def eddy_norm_lifetime(self, name, nb, factor=1):
24+
"""
25+
:param str,array name: Array or field name
26+
:param int nb: size of output array
27+
"""
28+
y = self.parse_varname(name)
29+
x = linspace(0, 1, nb)
30+
out = zeros(nb, dtype=y.dtype)
31+
nb_track = 0
32+
for i, b0, b1 in self.iter_on("track"):
33+
y_ = y[i]
34+
size_ = y_.size
35+
if size_ == 0:
36+
continue
37+
sum_profile(x, y_, out)
38+
nb_track += 1
39+
return x, out / nb_track * factor
40+
41+
42+
# %%
43+
# Load atlas
44+
# ----------
45+
kw = dict(include_vars=("speed_radius", "amplitude", "track"))
46+
a = MyObs.load_file(
47+
get_demo_path("eddies_med_adt_allsat_dt2018/Anticyclonic.zarr"), **kw
48+
)
49+
c = MyObs.load_file(get_demo_path("eddies_med_adt_allsat_dt2018/Cyclonic.zarr"), **kw)
50+
51+
nb_max_a = a.nb_obs_by_track.max()
52+
nb_max_c = c.nb_obs_by_track.max()
53+
54+
# %%
55+
# Compute normalize lifetime
56+
# --------------------------
57+
58+
# Radius
59+
AC_radius = a.eddy_norm_lifetime("speed_radius", nb=nb_max_a, factor=1e-3)
60+
CC_radius = c.eddy_norm_lifetime("speed_radius", nb=nb_max_c, factor=1e-3)
61+
# Amplitude
62+
AC_amplitude = a.eddy_norm_lifetime("amplitude", nb=nb_max_a, factor=1e2)
63+
CC_amplitude = c.eddy_norm_lifetime("amplitude", nb=nb_max_c, factor=1e2)
64+
65+
# %%
66+
# Figure
67+
# ------
68+
fig, axs = plt.subplots(nrows=2, figsize=(8, 6))
69+
70+
axs[0].set_title("Normalised Mean Radius")
71+
axs[0].plot(*AC_radius), axs[0].plot(*CC_radius)
72+
axs[0].set_ylabel("Radius (km)"), axs[0].grid()
73+
axs[0].set_xlim(0, 1), axs[0].set_ylim(0, None)
74+
75+
axs[1].set_title("Normalised Mean Amplitude")
76+
axs[1].plot(*AC_amplitude, label="AC"), axs[1].plot(*CC_amplitude, label="CC")
77+
axs[1].set_ylabel("Amplitude (cm)"), axs[1].grid(), axs[1].legend()
78+
_ = axs[1].set_xlim(0, 1), axs[1].set_ylim(0, None)

‎notebooks/python_module/10_tracking_diagnostics/pet_normalised_lifetime.ipynb‎

Lines changed: 119 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,119 @@
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# Normalised Eddy Lifetimes\n\nExample from Evan Mason\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 numba import njit\nfrom numpy import interp, linspace, zeros\nfrom py_eddy_tracker_sample import get_demo_path\n\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+
"@njit(cache=True)\ndef sum_profile(x_new, y, out):\n \"\"\"Will sum all interpolated given array\"\"\"\n out += interp(x_new, linspace(0, 1, y.size), y)\n\n\nclass MyObs(TrackEddiesObservations):\n def eddy_norm_lifetime(self, name, nb, factor=1):\n \"\"\"\n :param str,array name: Array or field name\n :param int nb: size of output array\n \"\"\"\n y = self.parse_varname(name)\n x = linspace(0, 1, nb)\n out = zeros(nb, dtype=y.dtype)\n nb_track = 0\n for i, b0, b1 in self.iter_on(\"track\"):\n y_ = y[i]\n size_ = y_.size\n if size_ == 0:\n continue\n sum_profile(x, y_, out)\n nb_track += 1\n return x, out / nb_track * factor"
41+
]
42+
},
43+
{
44+
"cell_type": "markdown",
45+
"metadata": {},
46+
"source": [
47+
"## Load atlas\n\n"
48+
]
49+
},
50+
{
51+
"cell_type": "code",
52+
"execution_count": null,
53+
"metadata": {
54+
"collapsed": false
55+
},
56+
"outputs": [],
57+
"source": [
58+
"kw = dict(include_vars=(\"speed_radius\", \"amplitude\", \"track\"))\na = MyObs.load_file(\n get_demo_path(\"eddies_med_adt_allsat_dt2018/Anticyclonic.zarr\"), **kw\n)\nc = MyObs.load_file(get_demo_path(\"eddies_med_adt_allsat_dt2018/Cyclonic.zarr\"), **kw)\n\nnb_max_a = a.nb_obs_by_track.max()\nnb_max_c = c.nb_obs_by_track.max()"
59+
]
60+
},
61+
{
62+
"cell_type": "markdown",
63+
"metadata": {},
64+
"source": [
65+
"## Compute normalize lifetime\n\n"
66+
]
67+
},
68+
{
69+
"cell_type": "code",
70+
"execution_count": null,
71+
"metadata": {
72+
"collapsed": false
73+
},
74+
"outputs": [],
75+
"source": [
76+
"# Radius\nAC_radius = a.eddy_norm_lifetime(\"speed_radius\", nb=nb_max_a, factor=1e-3)\nCC_radius = c.eddy_norm_lifetime(\"speed_radius\", nb=nb_max_c, factor=1e-3)\n# Amplitude\nAC_amplitude = a.eddy_norm_lifetime(\"amplitude\", nb=nb_max_a, factor=1e2)\nCC_amplitude = c.eddy_norm_lifetime(\"amplitude\", nb=nb_max_c, factor=1e2)"
77+
]
78+
},
79+
{
80+
"cell_type": "markdown",
81+
"metadata": {},
82+
"source": [
83+
"## Figure\n\n"
84+
]
85+
},
86+
{
87+
"cell_type": "code",
88+
"execution_count": null,
89+
"metadata": {
90+
"collapsed": false
91+
},
92+
"outputs": [],
93+
"source": [
94+
"fig, axs = plt.subplots(nrows=2, figsize=(8, 6))\n\naxs[0].set_title(\"Normalised Mean Radius\")\naxs[0].plot(*AC_radius), axs[0].plot(*CC_radius)\naxs[0].set_ylabel(\"Radius (km)\"), axs[0].grid()\naxs[0].set_xlim(0, 1), axs[0].set_ylim(0, None)\n\naxs[1].set_title(\"Normalised Mean Amplitude\")\naxs[1].plot(*AC_amplitude, label=\"AC\"), axs[1].plot(*CC_amplitude, label=\"CC\")\naxs[1].set_ylabel(\"Amplitude (cm)\"), axs[1].grid(), axs[1].legend()\n_ = axs[1].set_xlim(0, 1), axs[1].set_ylim(0, None)"
95+
]
96+
}
97+
],
98+
"metadata": {
99+
"kernelspec": {
100+
"display_name": "Python 3",
101+
"language": "python",
102+
"name": "python3"
103+
},
104+
"language_info": {
105+
"codemirror_mode": {
106+
"name": "ipython",
107+
"version": 3
108+
},
109+
"file_extension": ".py",
110+
"mimetype": "text/x-python",
111+
"name": "python",
112+
"nbconvert_exporter": "python",
113+
"pygments_lexer": "ipython3",
114+
"version": "3.7.7"
115+
}
116+
},
117+
"nbformat": 4,
118+
"nbformat_minor": 0
119+
}

‎share/nrt_global_allsat_phy_l4_20190223_20190226.nc‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1 @@
1+
../src/py_eddy_tracker/data/nrt_global_allsat_phy_l4_20190223_20190226.nc

‎src/py_eddy_tracker/observations/network.py‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -575,7 +575,7 @@ def numbering_segment(self, start=0):
575575
New numbering of segment
576576
"""
577577
for i, _, _ in self.iter_on("track"):
578-
new_numbering(self.segment[i])
578+
new_numbering(self.segment[i], start)
579579

580580
def numbering_network(self, start=1):
581581
"""

0 commit comments

Comments
 (0)