Skip to content

Commit 67d78e4

Browse files
erikvansebillereint-fischer
authored andcommitted
Adding a tutorial to show gsw use in Kernels
1 parent 58f02bd commit 67d78e4

File tree

3 files changed

+157
-0
lines changed

3 files changed

+157
-0
lines changed
Lines changed: 155 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,155 @@
1+
{
2+
"cells": [
3+
{
4+
"cell_type": "markdown",
5+
"id": "0",
6+
"metadata": {},
7+
"source": [
8+
"# Using the `gsw` toolbox to compute density"
9+
]
10+
},
11+
{
12+
"cell_type": "markdown",
13+
"id": "1",
14+
"metadata": {},
15+
"source": [
16+
"This tutorial shows how to use the [`gsw` toolbox](https://teos-10.github.io/GSW-Python/) (the Gibbs SeaWater Oceanographic Toolbox of TEOS-10) within Parcels to compute density from temperature and salinity fields. The `gsw` toolbox can be installed via `conda install gsw`."
17+
]
18+
},
19+
{
20+
"cell_type": "markdown",
21+
"id": "2",
22+
"metadata": {},
23+
"source": [
24+
"First, load the necessary libraries and the data:"
25+
]
26+
},
27+
{
28+
"cell_type": "code",
29+
"execution_count": null,
30+
"id": "3",
31+
"metadata": {},
32+
"outputs": [],
33+
"source": [
34+
"import numpy as np\n",
35+
"import xarray as xr\n",
36+
"\n",
37+
"import parcels\n",
38+
"\n",
39+
"# Load the CopernicusMarine data in the Agulhas region from the example_datasets\n",
40+
"example_dataset_folder = parcels.download_example_dataset(\n",
41+
" \"CopernicusMarine_data_for_Argo_tutorial\"\n",
42+
")\n",
43+
"\n",
44+
"ds = xr.open_mfdataset(f\"{example_dataset_folder}/*.nc\", combine=\"by_coords\")\n",
45+
"\n",
46+
"# TODO check how we can get good performance without loading full dataset in memory\n",
47+
"ds.load() # load the dataset into memory\n",
48+
"\n",
49+
"fieldset = parcels.FieldSet.from_copernicusmarine(ds)"
50+
]
51+
},
52+
{
53+
"cell_type": "markdown",
54+
"id": "4",
55+
"metadata": {},
56+
"source": [
57+
"Now, define a custom Particle class that includes temperature, salinity, and density as variables, and create a ParticleSet with one particle at a known location:"
58+
]
59+
},
60+
{
61+
"cell_type": "code",
62+
"execution_count": null,
63+
"id": "5",
64+
"metadata": {},
65+
"outputs": [],
66+
"source": [
67+
"GSWParticle = parcels.Particle.add_variable(\n",
68+
" [\n",
69+
" parcels.Variable(\"temp\", dtype=np.float32, initial=np.nan),\n",
70+
" parcels.Variable(\"salt\", dtype=np.float32, initial=np.nan),\n",
71+
" parcels.Variable(\"density\", dtype=np.float32, initial=np.nan),\n",
72+
" ]\n",
73+
")\n",
74+
"\n",
75+
"# Initiate one Argo float in the Agulhas Current\n",
76+
"pset = parcels.ParticleSet(\n",
77+
" fieldset=fieldset,\n",
78+
" pclass=GSWParticle,\n",
79+
" lon=[32],\n",
80+
" lat=[-31],\n",
81+
" depth=[200],\n",
82+
")"
83+
]
84+
},
85+
{
86+
"cell_type": "markdown",
87+
"id": "6",
88+
"metadata": {},
89+
"source": [
90+
"Now (as the core part of this tutorial) define a custom kernel that uses the `gsw` toolbox to compute density from temperature and salinity:"
91+
]
92+
},
93+
{
94+
"cell_type": "code",
95+
"execution_count": null,
96+
"id": "7",
97+
"metadata": {},
98+
"outputs": [],
99+
"source": [
100+
"def ParcelsGSW(particles, fieldset):\n",
101+
" import gsw\n",
102+
"\n",
103+
" particles.temp = fieldset.thetao[particles]\n",
104+
" particles.salt = fieldset.so[particles]\n",
105+
" pressure = gsw.p_from_z(-particles.depth, particles.lat)\n",
106+
" particles.density = gsw.density.rho(particles.salt, particles.temp, pressure)"
107+
]
108+
},
109+
{
110+
"cell_type": "markdown",
111+
"id": "8",
112+
"metadata": {},
113+
"source": [
114+
"Finally, run the `ParcelsGSW` Kernel for one timestep and check (for Continuous Integration purposes) that the computed density is as expected:"
115+
]
116+
},
117+
{
118+
"cell_type": "code",
119+
"execution_count": null,
120+
"id": "9",
121+
"metadata": {},
122+
"outputs": [],
123+
"source": [
124+
"pset.execute(ParcelsGSW, runtime=np.timedelta64(1, \"s\"), dt=np.timedelta64(1, \"s\"))\n",
125+
"\n",
126+
"np.testing.assert_allclose(pset.density, [1026.8281], rtol=1e-5)\n",
127+
"\n",
128+
"print(\n",
129+
" f\"Temperature: {pset.temp[0]:.2f}, Salinity: {pset.salt[0]:.2f}, Density: {pset.density[0]:.2f}\"\n",
130+
")"
131+
]
132+
}
133+
],
134+
"metadata": {
135+
"kernelspec": {
136+
"display_name": "parcels",
137+
"language": "python",
138+
"name": "python3"
139+
},
140+
"language_info": {
141+
"codemirror_mode": {
142+
"name": "ipython",
143+
"version": 3
144+
},
145+
"file_extension": ".py",
146+
"mimetype": "text/x-python",
147+
"name": "python",
148+
"nbconvert_exporter": "python",
149+
"pygments_lexer": "ipython3",
150+
"version": "3.13.5"
151+
}
152+
},
153+
"nbformat": 4,
154+
"nbformat_minor": 5
155+
}

parcels/_tutorial.py

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -47,6 +47,7 @@
4747
],
4848
"CopernicusMarine_data_for_Argo_tutorial": [
4949
"cmems_mod_glo_phy-cur_anfc_0.083deg_P1D-m_uo-vo_31.00E-33.00E_33.00S-30.00S_0.49-2225.08m_2024-01-01-2024-02-01.nc",
50+
"cmems_mod_glo_phy-so_anfc_0.083deg_P1D-m_so_31.00E-33.00E_33.00S-30.00S_0.49-2225.08m_2024-01-01-2024-02-01.nc",
5051
"cmems_mod_glo_phy-thetao_anfc_0.083deg_P1D-m_thetao_31.00E-33.00E_33.00S-30.00S_0.49-2225.08m_2024-01-01-2024-02-01.nc",
5152
],
5253
"DecayingMovingEddy_data": [

pixi.toml

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -65,6 +65,7 @@ tests-notebooks = "pytest --nbval-lax docs/examples"
6565
jupyter = "*"
6666
trajan = "*"
6767
matplotlib-base = ">=2.0.2"
68+
gsw = "*"
6869

6970
[feature.docs.dependencies]
7071
numpydoc = "!=1.9.0"

0 commit comments

Comments
 (0)