One huge commit that included two new notebooks

docs/_config.yml

@@ -8,6 +8,7 @@ logo: logo.png
 # Force re-execution of notebooks on each build.
 # See https://jupyterbook.org/content/execute.html
 execute:
+  timeout: 300
   execute_notebooks: 'auto'
   stderr_output: "remove"
   exclude_patterns:

docs/_toc.yml

@@ -31,6 +31,8 @@ parts:
     sections:
     - file: notebook/global_ECCO
     - file: sciserver_notebooks/IGP
+    - file: sciserver_notebooks/Fjord
+    - file: sciserver_notebooks/LLC4320
     - file: notebook/AVISO
     - file: notebook/topology_tutorial
 - caption: Contributing

docs/contribute/contr_doc.md

@@ -61,6 +61,15 @@ od = ospy.open_oceandataset.from_catalog('NameOfDataset')
 ds = od._ds
 ```
 
+Note that since you are using an `Oceanography` image, most packages are already downloaded. We still need two dependencies just to convert the notebooks.
+
+```
+pip install -U jupyter-book
+pip install jupytext
+```
+
+And that's it, you don't have to follow the steps of preparing environment.
+
 Now, follow these steps:
 
 1. [Fork and clone](use_git.md) this repo adjacent to the seaduck directory

docs/get_started.md

@@ -16,7 +16,7 @@ If you have questions about any function, you can always go to [API references](
 Running Lagrangian particle on regional dataset
 ^^^
 
-```{image} https://github.com/MaceKuailv/seaduck_sciserver_notebook/blob/master/IGP_files/IGP_27_1.png?raw=true
+```{image} https://github.com/MaceKuailv/seaduck_sciserver_notebook/blob/master/stable_images/IGP_32_0.png?raw=true
 :height: 200
 ```
 
@@ -34,14 +34,17 @@ ECCO dataset for now, but I want to have a 4320 example notebook.
 :::
 
 :::{grid-item-card}
-:link: notebook/AVISO
+:link: sciserver_notebooks/Fjord
 :link-type: doc
 :class-header: bg-light
 
-Sketch the ACC
+Interpolation near a fjord in Eastern Greenland
 ^^^
 
-Could have a prettier example using ETOPO.
+```{image} https://github.com/MaceKuailv/seaduck_sciserver_notebook/blob/master/stable_images/Fjord_29_0.png?raw=true
+:height: 200
+```
+
 :::
 
 ::::

docs/network_of_object.md

@@ -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.                           | [horizontal stream function conservation](idealize_test/hor_stream) |
+| [`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.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)   |

docs/sciserver_notebooks/Fjord.md

@@ -0,0 +1,223 @@
+---
+jupytext:
+  text_representation:
+    extension: .md
+    format_name: myst
+    format_version: 0.13
+    jupytext_version: 1.14.6
+kernelspec:
+  display_name: Oceanography
+  language: python
+  name: oceanography
+---
+
+# Demonstrate Position object with Fjord
+
+Author: Wenrui Jiang, 14 June 2023
+> **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/Fjord.ipynb. The `eulerian.Position` object is really what connect a point and the numerical model. Its `interpolate` method really is the core of this package. We are going to use a rather interesting example to demonstrate the functionalities of `eulerian.Position`.
+
+```{code-cell} ipython3
+import oceanspy as ospy
+import numpy as np
+import matplotlib.pyplot as plt
+import matplotlib as mpl
+import seaduck as sd
+from functools import partial
+import xarray as xr
+import cmocean
+
+mpl.rcParams["figure.dpi"] = 300
+```
+
+Here on sciserver, we have an interesting dataset simulating the interaction between background circulation and Kangerdlugssuaq Fjord. More information can be found below:
+
+```{code-cell} ipython3
+---
+jupyter:
+  outputs_hidden: true
+---
+fjord = ospy.open_oceandataset.from_catalog("KangerFjord")
+```
+
+```
+Opening KangerFjord.
+A realistic numerical model constructed to simulate the oceanic conditions
+and circulation in a large southeast Greenland fjord (Kangerdlugssuaq) and
+the adjacent shelf sea region during winter 2007–2008.
+Citation:
+ * Fraser et al., 2018 - JGR.
+```
+
++++
+
+Let's first look at what the dataset looks like. We are going to use ETOPO dataset to give you an idea where this domain is located.
+
+```{code-cell} ipython3
+etopo = ospy.open_oceandataset.from_catalog("ETOPO")
+etopo = etopo._ds.sel(X=slice(-40.6, -13.1), Y=slice(62.4, 70.5)).where(
+    etopo._ds.Depth > 0
+)
+```
+
+```{code-cell} ipython3
+:tags: [hide-input]
+
+vmax = 0.2
+cmap = cmocean.cm.balance
+```
+
+```{code-cell} ipython3
+:tags: [hide-input]
+
+c = plt.pcolormesh(
+    fjord._ds.XC, fjord._ds.YC, fjord._ds.Eta[263], vmax=vmax, vmin=-vmax, cmap=cmap
+)
+plt.pcolormesh(
+    etopo.X,
+    etopo.Y,
+    etopo.Depth,
+    cmap="bone",
+    vmin=-3000,
+    vmax=4500,
+    zorder=20,
+)
+plt.gca().set_facecolor("lightsteelblue")
+plt.title("Model domain and the sea surface height")
+plt.xlabel("Longitude")
+plt.ylabel("Latitude")
+plt.colorbar(c, label="m")
+plt.show()
+```
+![png](https://github.com/MaceKuailv/seaduck_sciserver_notebook/blob/master/Fjord_files/Fjord_8_0.png?raw=true)
+
++++ {"tags": ["mdformat-skip"]}
+
+## Doing Interpolation
+We are going to use the sea surface height field $\eta$ as an example.
+
+First, we are going to convert the `xarray.Dataset` to `seaduck.OceData`.
+
+```{code-cell} ipython3
+tub = sd.OceData(fjord._ds)
+```
+
+For now, the interpolation points are simply the gridpoints of the dataset.
+
+```{code-cell} ipython3
+:tags: [hide-input]
+
+x = tub.XC.ravel()
+y = tub.YC.ravel()
+z = np.ones_like(x) * -5.0
+t = np.ones_like(x) * sd.utils.convert_time("2008-01-01")
+```
+
+Let's create the `eulerian.Position` object
+
+```{code-cell} ipython3
+p = sd.Position().from_latlon(x, y, z, t, data=tub)
+```
+
+Two kernels are defined here, both of them are the default kernel used by the package. However, we are going to "hack" one of them to demonstrate the "cascade" capacity of interpolation.
+
+```{code-cell} ipython3
+kernel = sd.KnW()
+kernel_to_be_hacked = sd.KnW()
+```
+
++++ {"tags": ["mdformat-skip"]}
+
+First, we do the normal interpolation on $\eta$ in the normal way, and plot it
+
+```{code-cell} ipython3
+eta = p.interpolate("Eta", kernel)
+```
+
+```{code-cell} ipython3
+:tags: [hide-input]
+
+plt.scatter(x, y, c=eta, s=0.7, vmax=vmax, vmin=-vmax, cmap=cmap)
+plt.title("Interpolated sea surface height")
+plt.colorbar(label="m")
+plt.xlabel("Longitude")
+plt.ylabel("Latitude")
+plt.show()
+```
+![png](https://github.com/MaceKuailv/seaduck_sciserver_notebook/blob/master/Fjord_files/Fjord_19_0.png?raw=true)
+
+We are going to "hack" the code to make it return the size of the kernels used.
+
+The details of the hack is not very important. But if you are interested you can read the inline comments below.
+
+```{code-cell} ipython3
+:tags: [hide-input]
+
+# The partials of this function is going to replace all the weight functions.
+# For legal weight functions, the sum of contribution (weight) should be one,
+# but here every neighboring point has contribution of one
+def dummy_weight(rx, ry, ans):
+    n = len(rx)
+    return np.ones((n, ans))
+
+
+# This is the step that overwrites the existing weight function. Highly not recommended.
+kernel_to_be_hacked.hfuncs = [partial(dummy_weight, ans=len(k)) for k in kernel.kernels]
+# Create a dataset of ones with the same shape as Eta.
+tub["ones"] = xr.ones_like(tub["Eta"])
+```
+
+```{code-cell} ipython3
+how_many = p.interpolate("ones", kernel_to_be_hacked)
+```
+
+Now, we can look at what the interpolation package does under the hood.
+
+```{code-cell} ipython3
+:tags: [hide-input]
+
+plt.scatter(x, y, c=how_many, s=0.7, cmap="Set1", vmax=9.5, vmin=0.5)
+plt.title("The number of grid points used for interpolation")
+plt.colorbar(label="Point count")
+plt.xlabel("Longitude")
+plt.ylabel("Latitude")
+plt.show()
+```
+![png](https://github.com/MaceKuailv/seaduck_sciserver_notebook/blob/master/Fjord_files/Fjord_24_0.png?raw=true)
+
+## Filling between
+
+As you have seen earlier, the grid of this dataset is very uneven.
+
+Well, the strength of the `eulerian.Position` object is filling in information in between, so let's do that.
+
+```{code-cell} ipython3
+:tags: [hide-input]
+
+Nlon = 600
+Nlat = 100
+Ndep = 1
+ax, ay, az, at = sd.utils.easy_3d_cube(
+    (-34.5, -28.5, Nlon), (66.5, 67.0, Nlat), (-5.0, -5.0, Ndep), "2008-01-01"
+)
+```
+
+Create `eulerian.Position` in between, and plot it!
+
+```{code-cell} ipython3
+addition = sd.Position().from_latlon(ax, ay, az, at, data=tub)
+more_eta = addition.interpolate("Eta", kernel)
+```
+
+```{code-cell} ipython3
+:tags: [hide-input]
+
+plt.scatter(x, y, c=eta, s=0.7, vmax=vmax, vmin=-vmax, cmap=cmap)
+plt.scatter(ax, ay, c=more_eta, s=0.7, vmax=vmax, vmin=-vmax, cmap=cmap)
+plt.title("Interpolated sea surface height (Reinforced)")
+plt.xlabel("Longitude")
+plt.ylabel("Latitude")
+plt.show()
+```
+![png](https://github.com/MaceKuailv/seaduck_sciserver_notebook/blob/master/Fjord_files/Fjord_29_0.png?raw=true)
+
+I'd say the filling-in is done pretty well!