Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Twnh testing sci server #68

Merged
merged 3 commits into from
Jul 6, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
5 changes: 3 additions & 2 deletions docs/_toc.yml
Original file line number Diff line number Diff line change
Expand Up @@ -34,9 +34,10 @@ parts:
- file: notebook/topology_tutorial
- file: sciserver_example
sections:
- file: sciserver_notebooks/IGP
- file: sciserver_notebooks/Fjord
- file: sciserver_notebooks/IGPwinter
- file: sciserver_notebooks/KangerFjord
- file: sciserver_notebooks/LLC4320
- file: sciserver_notebooks/ECCO_plot_stations
- caption: Contributing
chapters:
- file: guide_for_developer
Expand Down
4 changes: 2 additions & 2 deletions docs/get_started.md
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@ If you have questions about any function, you can always go to [API references](
:gutter: 3

:::{grid-item-card}
:link: sciserver_notebooks/IGP
:link: sciserver_notebooks/IGPwinter
:link-type: doc
:class-header: bg-light

Expand All @@ -34,7 +34,7 @@ ECCO dataset for now, but I want to have a 4320 example notebook.
:::

:::{grid-item-card}
:link: sciserver_notebooks/Fjord
:link: sciserver_notebooks/KangerFjord
:link-type: doc
:class-header: bg-light

Expand Down
4 changes: 2 additions & 2 deletions docs/network_of_object.md
Original file line number Diff line number Diff line change
Expand Up @@ -7,6 +7,6 @@ Another good way to learn about the package is to learn how different objects ar
| [`seaduck.OceData`](api_reference/apiref_OceData.rst) | Interface to data, translate between lat-lon and indices. | [AVISO](notebook/AVISO.ipynb) |
| [`seaduck.Topology`](api_reference/apiref_topology.rst) | Describe the how the grids are connected. Similar to `xgcm`. | [topology tutorial](notebook/topology_tutorial.ipynb) |
| [`seaduck.KnW`](api_reference/apiref_kernelNweight.rst) | Define what interpolation/derivative to perform. | [ECCO](notebook/global_ECCO.ipynb) |
| [`seaduck.Position`](api_reference/apiref_eulerian.rst) | Interpolate at Eulerian positions. | [Interpolate in a fjord](sciserver_notebooks/Fjord.md) |
| [`seaduck.Particle`](api_reference/apiref_lagrangian.rst) | Perform Lagrangian particle simulation. Sub class of `seaduck.Position` | [Regional simulation](sciserver_notebooks/IGP.md) |
| [`seaduck.Position`](api_reference/apiref_eulerian.rst) | Interpolate at Eulerian positions. | [Interpolate in a fjord](sciserver_notebooks/KangerFjord.md) |
| [`seaduck.Particle`](api_reference/apiref_lagrangian.rst) | Perform Lagrangian particle simulation. Sub class of `seaduck.Position` | [Regional simulation](sciserver_notebooks/IGPwinter.md) |
| [`seaduck.OceInterp`](api_reference/apiref_OceInterp.rst) | Uniform interface to Lagrangian and Eulerian operations. | [one minute guide](one_min_guide.ipynb), [ECCO](notebook/global_ECCO.ipynb) |
162 changes: 162 additions & 0 deletions docs/sciserver_notebooks/ECCO_plot_stations.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,162 @@
---
jupytext:
text_representation:
extension: .md
format_name: myst
format_version: 0.13
jupytext_version: 1.14.7
kernelspec:
display_name: Oceanography
language: python
name: oceanography
---

# `OceInterp` interpolation in the ECCO global state estimate.

Extract temperature/salinity profiles at specified longitudes, latitudes, and times. This process mimics hydrographic stations from in-situ measurements (like a ship CTD station, or Argo float profile).

This notebook uses [Oceanspy](https://oceanspy.readthedocs.io/en/latest/) and demonstrates the interface to the Poseidon-viewer on SciServer.

Author: Tom Haine & Wenrui Jiang, Jun '23
> **Warning**⚠️ : the notebook was last ran on **2023-07-06** with **seaduck 0.1.dev433+g99d4d3b**. You can find the executable version at https://github.com/MaceKuailv/seaduck_sciserver_notebook/blob/master/ECCO_plot_stations.ipynb.
```{code-cell} ipython3
import seaduck as sd
import oceanspy as ospy
from poseidon_viewer import get_shapes

import numpy as np

import matplotlib.pyplot as plt
plt.rcParams['figure.figsize'] = 12, 8
```

```{admonition} Access ECCO
The global MITgcm run is the ECCO state estimate and is publicly available on [SciServer](sciserver.org) (from the Oceanography container). The simulation output can be opened using the [OceanSpy](https://github.com/hainegroup/oceanspyhttps://github.com/hainegroup/oceanspy) package using the [`from_catalog`](https://oceanspy.readthedocs.io/en/latest/generated/oceanspy.open_oceandataset.from_catalog.html#oceanspy.open_oceandataset.from_catalog) method.

`ecco = ospy.open_oceandataset.from_catalog("ECCO")`

Click [here](https://dev-poseidon-ocean.pantheonsite.io/products/datasets/) for a full list of the dataset hosted and [here](https://oceanspy.readthedocs.io/en/latest/datasets.html#igpwinterhttps://oceanspy.readthedocs.io/en/latest/datasets.html#igpwinter) to find out more.
```

```{code-cell} ipython3
od = ospy.open_oceandataset.from_catalog('ECCO')

od._ds = od._ds.drop_vars({'k', 'k_u', 'k_p1', 'k_l'})
co_list = [var for var in od._ds.variables if "time" not in od._ds[var].dims]
od._ds = od._ds.set_coords(co_list)

od._ds = od._ds.drop_vars('time_bnds')
od._ds['Temp'] = od._ds['THETA']
od._ds['S'] = od._ds['SALT']
this_time = '1992-03-16T12:00:00.000000000' # Select time of interest. This is cut and pasted (with editing to convert format) from Poseidon-viewer, but in future should be returned in the JSON object.

varList = ['Temp', 'S'] # Select variables of interest.
Nvars = len(varList)
od._ds = od._ds.drop_vars([var for var in od._ds.data_vars if var not in varList])
```

Define coordinates for interpolation from the Poseidon-viewer (or enter them manually):

```{code-cell} ipython3
ps = [{"type":"Point","coordinates":[-57.710530333876235,33.062070736515295]},{"type":"Point","coordinates":[-89.55599488497559,25.453020491528747]},{"type":"Point","coordinates":[-74.70605786848537,13.952759739624767]},{"type":"Point","coordinates":[-59.06355926811892,21.999755420813273]},{"type":"Point","coordinates":[-41.69588507736159,22.43955659234939]},{"type":"Point","coordinates":[-45.113026771172656,14.563543828761837]},{"type":"Point","coordinates":[-26.081872200732885,6.414099524482438]},{"type":"Point","coordinates":[-17.396102656758963,-4.381322875209875]},{"type":"Point","coordinates":[-26.702603318403906,-7.125636489486197]},{"type":"Point","coordinates":[-32.51011235240231,-22.847802807885373]}]
lons, lats = ospy.utils.viewer_to_range(ps)
Nstations = len(lons)
```

### `seaduck.OceInterp` interpolates the temperature and salinity to the specified coordinates.
This process makes synthetic hydrographic profiles from ECCO. Compute the potential density anomaly from the T/S data too.

```{code-cell} ipython3
tub = sd.OceData(od._ds)
times = sd.utils.convert_time(this_time)
Ntimes = times.size

depths = tub.Z.ravel()
Ndepths = len(depths)

sd_lons,sd_times,sd_depths = np.meshgrid(lons,times,depths,indexing='ij')
sd_lats,_,_ = np.meshgrid(lats,times,depths,indexing='ij')
sd_lons = np.ravel(sd_lons)
sd_lats = np.ravel(sd_lats)
sd_times = np.ravel(sd_times)
sd_depths = np.ravel(sd_depths)

[sd_Temp,sd_S] = np.array(sd.OceInterp(tub,varList, sd_lons, sd_lats, sd_depths, sd_times))
sd_depths2 = sd_depths.reshape(Nstations,Ntimes,Ndepths).squeeze()
sd_S = sd_S.reshape(Nstations,Ntimes,Ndepths).squeeze()
sd_Temp = sd_Temp.reshape(Nstations,Ntimes,Ndepths).squeeze()

sd_sigma0 = ospy.utils.densjmd95(sd_S,sd_Temp,0) - 1000.0
```

### Plot station locations:

```{code-cell} ipython3
cut_od = od.subsample.cutout(XRange = [min(lons)-10,max(lons)+10],YRange = [min(lats)-5,max(lats)+5])
fig = plt.figure(figsize=(12, 5))
ax = cut_od.plot.horizontal_section(varName="Depth", cmap='Greys_r')
colors=['r', 'b', 'k', 'orange', 'darkgreen', 'm', 'yellow', 'c', 'indigo', 'magenta', 'green']
legestations_ds = np.arange(Nstations)
for i in range(Nstations):
plt.plot(lons[i], lats[i], 'o', color=colors[i], label='station ' + str(legestations_ds[i]))
plt.legend()
plt.show()
```
![png](https://github.com/MaceKuailv/seaduck_sciserver_notebook/blob/master/ECCO_plot_stations_files/ECCO_plot_stations_11_0.png?raw=true)

### Plot vertical T/S station profiles:

```{code-cell} ipython3
for i in range(Nstations):

plt.subplot(1,3,1)
plt.plot(sd_Temp[i,:],depths,'o', color=colors[i], label='station ' + str(legestations_ds[i]))
plt.xlabel('Temperature [oC]')
plt.ylabel('Depth [m]')
plt.grid()

plt.subplot(1,3,2)
plt.plot(sd_S[i,:],depths,'o', color=colors[i], label='station ' + str(legestations_ds[i]))
plt.xlabel('Salinity [g/kg]')
plt.grid()


plt.subplot(1,3,3)
plt.plot(sd_sigma0[i,:],depths,'o', color=colors[i], label='station ' + str(legestations_ds[i]))
plt.xlabel('Density anomaly [kg/m^3]')
plt.grid()

plt.subplot(1,3,1)
plt.xlabel('Temperature [oC]')
plt.ylabel('Depth [m]')
plt.grid()

plt.subplot(1,3,2)
plt.xlabel('Salinity [g/kg]')
plt.grid()

plt.subplot(1,3,3)
plt.xlabel('Density anomaly [kg/m^3]')
plt.grid()
plt.legend()

plt.show()
```
![png](https://github.com/MaceKuailv/seaduck_sciserver_notebook/blob/master/ECCO_plot_stations_files/ECCO_plot_stations_13_0.png?raw=true)

### Plot T/S diagram.

```{code-cell} ipython3
for i in range(Nstations):
plt.plot(sd_S[i,:],sd_Temp[i,:],'o', color=colors[i], label='station ' + str(legestations_ds[i]))
plt.xlabel('Salinity [g/kg]')
plt.ylabel('Temperature [oC]')
plt.legend()
plt.grid()
plt.show()
```
![png](https://github.com/MaceKuailv/seaduck_sciserver_notebook/blob/master/ECCO_plot_stations_files/ECCO_plot_stations_9_1.png?raw=true)

```{code-cell} ipython3

```
Original file line number Diff line number Diff line change
Expand Up @@ -4,17 +4,17 @@ jupytext:
extension: .md
format_name: myst
format_version: 0.13
jupytext_version: 1.14.6
jupytext_version: 1.14.7
kernelspec:
display_name: Oceanography
language: python
name: oceanography
---

# `Particle` in a East Greenland regional simulation
# `Particle` in an East Greenland regional simulation

Author: Wenrui Jiang, Tom Haine Feb '23
> **Warning**⚠️ : the notebook was last ran on **2023-06-15** with **seaduck 0.1.3.dev44+g8d60702**. You can find the executable version at https://github.com/MaceKuailv/seaduck_sciserver_notebook/blob/master/IGP.ipynb.
> **Warning**⚠️ : the notebook was last ran on **2023-07-06** with **seaduck 0.1.dev433+g99d4d3b**. You can find the executable version at https://github.com/MaceKuailv/seaduck_sciserver_notebook/blob/master/IGPwinter.ipynb.
```{code-cell} ipython3
:tags: [hide-input]

Expand All @@ -25,29 +25,13 @@ import oceanspy as ospy
import seaduck as sd
```

```````{admonition} Access IGP
The regional MITgcm run is the IGPwinter simulationT and is publicly available on SciServer. The simulation output could opened using the OceanSpy package using the [`from_catalog`](https://oceanspy.readthedocs.io/en/latest/generated/oceanspy.open_oceandataset.from_catalog.html#oceanspy.open_oceandataset.from_catalog) method.
```{admonition} Access IGPwinter
The regional MITgcm run is the IGPwinter simulation and is publicly available on [SciServer](sciserver.org) (from the Oceanography container). The simulation output can be opened using the [OceanSpy](https://github.com/hainegroup/oceanspyhttps://github.com/hainegroup/oceanspy) package using the [`from_catalog`](https://oceanspy.readthedocs.io/en/latest/generated/oceanspy.open_oceandataset.from_catalog.html#oceanspy.open_oceandataset.from_catalog) method.

`ecco = ospy.open_oceandataset.from_catalog("ECCO")._ds`

or the daily-mean data ('daily_ecco').

`ecco = ospy.open_oceandataset.from_catalog('daily_ecco')._ds`
`ecco = ospy.open_oceandataset.from_catalog("IGPwinter")`

Click [here](https://dev-poseidon-ocean.pantheonsite.io/products/datasets/) for a full list of the dataset hosted and [here](https://oceanspy.readthedocs.io/en/latest/datasets.html#igpwinterhttps://oceanspy.readthedocs.io/en/latest/datasets.html#igpwinter) to find out more.
``````{admonition} Access full ECCO dataset
The ECCO dataset is publicly available on SciServer. The simulation output could opened using the OceanSpy package using the [`from_catalog`](https://oceanspy.readthedocs.io/en/latest/generated/oceanspy.open_oceandataset.from_catalog.html#oceanspy.open_oceandataset.from_catalog) method.

Choose between the monthly-mean data ('ECCO')

`ecco = ospy.open_oceandataset.from_catalog("ECCO")._ds`

or the daily-mean data ('daily_ecco').

`ecco = ospy.open_oceandataset.from_catalog('daily_ecco')._ds`

Click [here](https://dev-poseidon-ocean.pantheonsite.io/products/datasets/) for a full list of the dataset hosted.
```````
```

```{code-cell} ipython3
od = ospy.open_oceandataset.from_catalog("IGPwinter")
Expand All @@ -63,17 +47,17 @@ Citation:

+++

We are going to artificially create an open boundary in depth. This is to demonstrate that there is no weird behavior associated with open vertical boundaries.
We're going to artificially create an open boundary in depth. This is to demonstrate that there's no weird behavior associated with open vertical boundaries.

```{code-cell} ipython3
ds = od._ds.isel(Z=slice(0, 50), Zl=slice(0, 50))
```

## Prepare the simulation

Since this notebook works with an open domain, it is going to be a little different a few other notebooks.
Since this notebook works with an open domain (regional model, not global), it's going to differ a bit from other `seaduck` notebooks.

We are interested in looking at the coastal current. First, we need to prepare a `seaduck.OceData` object...
We're interested in looking at the coastal current. First, we need to prepare a `seaduck.OceData` object...

```{code-cell} ipython3
oce = sd.OceData(ds)
Expand All @@ -99,7 +83,7 @@ start_time = "2018-01-02"
t = np.ones_like(x) * sd.utils.convert_time(start_time)
```

...and integrate forward in time for one month.
...and integrate forward in time for one month (the actual trajectory calculation is made below; this is just setting parameters).

```{code-cell} ipython3
end_time = "2018-02-01"
Expand All @@ -119,13 +103,13 @@ plt.ylabel("Latitude")
plt.title("Bathymetry of model domain and particle initial position")
plt.show()
```
![png](https://github.com/MaceKuailv/seaduck_sciserver_notebook/blob/master/IGP_files/IGP_16_0.png?raw=true)
![png](https://github.com/MaceKuailv/seaduck_sciserver_notebook/blob/master/IGPwinter_files/IGPwinter_16_0.png?raw=true)

**Fig.1** The initial position of the particles released (red line) on a horizontal map of the surrounding regions.

+++

We are going to follow a cold puff of fresh water near the Greenland shelf. Since we have `OceanSpy` installed here at SciServer, let me quickly show you a vertical section of this region. (It really takes no effort at all!)
We are going to follow a cold puff of fresh water near the Greenland shelf. Since we have `OceanSpy` installed here at SciServer, let me quickly show you how to make a vertical section of this region. (It really takes no effort at all!)

```{code-cell} ipython3
od_surv = od.subsample.survey_stations(
Expand All @@ -136,15 +120,15 @@ od_surv.plot.vertical_section(varName="Temp", contourName="Sigma0")
plt.ylim([-750, 0])
plt.show()
```
![png](https://github.com/MaceKuailv/seaduck_sciserver_notebook/blob/master/IGP_files/IGP_19_1.png?raw=true)
![png](https://github.com/MaceKuailv/seaduck_sciserver_notebook/blob/master/IGPwinter_files/IGPwinter_19_1.png?raw=true)

**Fig.2** Vertical section at 71N. This vertical section goes a little further east than the particle position to include the shelf break. Colors denotes the potential temperature, while the contours is the potential density anomaly calculated using OceanSpy (with the equation of state the model used).
**Fig.2** Vertical section at 71N. This vertical section goes a little further east than the initial particle positions to include the shelf break. Colors denotes the potential temperature, while the contours are the potential density anomaly calculated using OceanSpy (with the equation of state the model used).

+++

Since we are in an open domain (in all three dimensions), given long enough time, some particles will go out of the domain and jeopardise the entire simulation. We can define a callback function to stop that from happening.
Since we are in an open domain (in all three dimensions), given long enough time, some particles will leave the domain. This will jeopardise the entire simulation! We can define a callback function to stop that from happening.

If provided, the function will be called every time particles cross walls. It can be used to manipulate the particles. But here we just use it as a stop/continue criterion.
If provided, the function will be called every time particles cross walls. It can be used to manipulate the particles and catch the out-of-domain issue. But here we just use it as a stop/continue criterion.

```{code-cell} ipython3
def continue_criterion(pt):
Expand All @@ -162,7 +146,7 @@ oce["vtrans"] = oce["V"] * oce["drF"] * oce["dxG"]
oce["wtrans"] = oce["W"] * oce["rA"]
```

Finally, create the particle object. And that's all of the preparation we need.
Finally, create the particle object:

```{code-cell} ipython3
p = sd.Particle(
Expand All @@ -182,6 +166,8 @@ p = sd.Particle(

## Perform the particle trajectory simulation

Run the simulation! This can take some time, so grab a break...

```{code-cell} ipython3
stops, raw = p.to_list_of_time([t[0], tf])
```
Expand Down Expand Up @@ -210,9 +196,9 @@ plt.ylabel("Latitude")
plt.title("Particle trajectories overlaid on bathymetry map")
plt.show()
```
![png](https://github.com/MaceKuailv/seaduck_sciserver_notebook/blob/master/IGP_files/IGP_32_0.png?raw=true)
![png](https://github.com/MaceKuailv/seaduck_sciserver_notebook/blob/master/IGPwinter_files/IGPwinter_32_0.png?raw=true)

**Fig.3** The particle trajectories overlaid on **Fig.1**. The coloring scheme of bathymetry is also the same as **Fig.1**
**Fig.3** The particle trajectories overlaid on **Fig.1**. The color scheme for bathymetry is the same as **Fig.1**.

```{code-cell} ipython3

Expand Down
Loading