Skip to content

Commit

Permalink
Merge pull request #113 from British-Oceanographic-Data-Centre/develop
Browse files Browse the repository at this point in the history
Update notebooks and configutation gallery
  • Loading branch information
soutobias authored Dec 11, 2023
2 parents fcfbee1 + cabd424 commit 14d7792
Show file tree
Hide file tree
Showing 27 changed files with 1,366 additions and 191 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -59,14 +59,14 @@ nemo = coast.Gridded(fn_nemo_dat, fn_nemo_dom, config=fn_nemo_config)
altimetry = coast.Altimetry(fn_altimetry, config=fn_altimetry_config)
```

/usr/share/miniconda/envs/coast/lib/python3.10/site-packages/xarray/core/dataset.py:278: UserWarning: The specified chunks separate the stored chunks along dimension "time_counter" starting at index 2. This could degrade performance. Instead, consider rechunking after loading.
/usr/share/miniconda/envs/coast/lib/python3.10/site-packages/xarray/core/dataset.py:282: UserWarning: The specified chunks separate the stored chunks along dimension "time_counter" starting at index 2. This could degrade performance. Instead, consider rechunking after loading.


././config/example_altimetry.json
Altimetry object at 0x560d5d3e2980 initialised
Altimetry object at 0x562b801a3980 initialised


/usr/share/miniconda/envs/coast/lib/python3.10/site-packages/xarray/core/dataset.py:278: UserWarning: The specified chunks separate the stored chunks along dimension "time" starting at index 1000. This could degrade performance. Instead, consider rechunking after loading.
/usr/share/miniconda/envs/coast/lib/python3.10/site-packages/xarray/core/dataset.py:282: UserWarning: The specified chunks separate the stored chunks along dimension "time" starting at index 1000. This could degrade performance. Instead, consider rechunking after loading.


### Subsetting
Expand All @@ -82,7 +82,7 @@ ind = ind[::4]
altimetry = altimetry.isel(t_dim=ind)
```

Subsetting Altimetry object at 0x560d5d3e2980 indices in [-10, 10], [45, 60]
Subsetting Altimetry object at 0x562b801a3980 indices in [-10, 10], [45, 60]


### Model interpolation
Expand All @@ -99,7 +99,7 @@ altimetry.obs_operator(nemo, mod_var_name="ssh", time_interp="nearest")
# to see for yourself.
```

Interpolating Gridded object at 0x560d5d3e2980 "ssh" with time_interp "nearest"
Interpolating Gridded object at 0x562b801a3980 "ssh" with time_interp "nearest"



Expand All @@ -118,7 +118,7 @@ altimetry.obs_operator(nemo, mod_var_name="ssh", time_interp="nearest")
stats = altimetry.basic_stats("ocean_tide_standard_name", "interp_ssh")
```

Altimetry object at 0x560d5d3e2980 initialised
Altimetry object at 0x562b801a3980 initialised



Expand Down Expand Up @@ -146,7 +146,7 @@ crps = altimetry.crps(nemo, model_var_name="ssh", obs_var_name="ocean_tide_stand
#crps.dataset # uncomment to print data object summary
```

Altimetry object at 0x560d5d3e2980 initialised
Altimetry object at 0x562b801a3980 initialised


### Plotting data
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -70,7 +70,7 @@ nemo_data = coast.Gridded(fn_data=fn_nemo_dat,

```

/usr/share/miniconda/envs/coast/lib/python3.10/site-packages/xarray/core/dataset.py:278: UserWarning: The specified chunks separate the stored chunks along dimension "time_counter" starting at index 2. This could degrade performance. Instead, consider rechunking after loading.
/usr/share/miniconda/envs/coast/lib/python3.10/site-packages/xarray/core/dataset.py:282: UserWarning: The specified chunks separate the stored chunks along dimension "time_counter" starting at index 2. This could degrade performance. Instead, consider rechunking after loading.


Calculate the climatology for temperature and sea surface height (ssh) as an example:
Expand Down
4 changes: 2 additions & 2 deletions content/en/docs/Examples/Notebooks/General/eof_tutorial.md
Original file line number Diff line number Diff line change
Expand Up @@ -89,7 +89,7 @@ eof_data.EOF.sel(mode=[1,2,3,4]).plot.pcolormesh(col='mode',col_wrap=2,x='longit



<xarray.plot.facetgrid.FacetGrid at 0x7fb34d977520>
<xarray.plot.facetgrid.FacetGrid at 0x7f0104539b40>



Expand All @@ -107,7 +107,7 @@ eof_data.temporal_proj.sel(mode=[1,2,3,4]).plot(col='mode',col_wrap=2,x='time')



<xarray.plot.facetgrid.FacetGrid at 0x7fb2f86bc760>
<xarray.plot.facetgrid.FacetGrid at 0x7f00ebff6290>



Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -58,7 +58,7 @@ ofile = "example_export_output.nc" # The target filename for output
nemo = coast.Gridded(fn_nemo_dat, fn_nemo_dom, config=config)
```

/usr/share/miniconda/envs/coast/lib/python3.10/site-packages/xarray/core/dataset.py:278: UserWarning: The specified chunks separate the stored chunks along dimension "time_counter" starting at index 2. This could degrade performance. Instead, consider rechunking after loading.
/usr/share/miniconda/envs/coast/lib/python3.10/site-packages/xarray/core/dataset.py:282: UserWarning: The specified chunks separate the stored chunks along dimension "time_counter" starting at index 2. This could degrade performance. Instead, consider rechunking after loading.


### We can export the whole xr.DataSet to a netCDF file
Expand Down
47 changes: 28 additions & 19 deletions content/en/docs/Examples/Notebooks/General/mask_maker_tutorial.md

Large diffs are not rendered by default.

255 changes: 255 additions & 0 deletions content/en/docs/Examples/Notebooks/General/polar_plotting.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,255 @@
---
title: "Polar plotting"
linkTitle: "Polar plotting"
weight: 5

description: >
Polar plotting example.
---
# Polar Plotting Tutorial

This demonstration will show how to re-project the NEMO velocities for quiver plotting in polar coordinates.

NEMO velocities are usually calculated and saved in along grid i and j directions. This causes an issue when plotting velocities as vectors on a map where it is assumed that i and j velocities point eastwards and northwards but the grid is curvilinear.

There are additional isses when plotting quivers over the poles that we will cover.

The data we use here comes from a global model that has been cropped to the Arctic region. If you want to see the images, please go to the configuration gallery codes in `example_scripts/configuration_gallery/polar_plotting.py` or in the [documentation](https://british-oceanographic-data-centre.github.io/COAsT/docs/examples/configs_gallery/).

**IMPORTANT: you have to install cartopy in order to plot the maps**


```python
# !pip install cartopy
```


```python
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy.feature as cfeature
import coast
```

/usr/share/miniconda/envs/coast/lib/python3.10/site-packages/pydap/lib.py:5: DeprecationWarning: pkg_resources is deprecated as an API. See https://setuptools.pypa.io/en/latest/pkg_resources.html
/usr/share/miniconda/envs/coast/lib/python3.10/site-packages/pkg_resources/__init__.py:2871: DeprecationWarning: Deprecated call to `pkg_resources.declare_namespace('pydap')`.
Implementing implicit namespace packages (as specified in PEP 420) is preferred to `pkg_resources.declare_namespace`. See https://setuptools.pypa.io/en/latest/references/keywords.html#keyword-namespace-packages
/usr/share/miniconda/envs/coast/lib/python3.10/site-packages/pkg_resources/__init__.py:2871: DeprecationWarning: Deprecated call to `pkg_resources.declare_namespace('pydap.responses')`.
Implementing implicit namespace packages (as specified in PEP 420) is preferred to `pkg_resources.declare_namespace`. See https://setuptools.pypa.io/en/latest/references/keywords.html#keyword-namespace-packages
/usr/share/miniconda/envs/coast/lib/python3.10/site-packages/pkg_resources/__init__.py:2350: DeprecationWarning: Deprecated call to `pkg_resources.declare_namespace('pydap')`.
Implementing implicit namespace packages (as specified in PEP 420) is preferred to `pkg_resources.declare_namespace`. See https://setuptools.pypa.io/en/latest/references/keywords.html#keyword-namespace-packages
/usr/share/miniconda/envs/coast/lib/python3.10/site-packages/pkg_resources/__init__.py:2871: DeprecationWarning: Deprecated call to `pkg_resources.declare_namespace('pydap.handlers')`.
Implementing implicit namespace packages (as specified in PEP 420) is preferred to `pkg_resources.declare_namespace`. See https://setuptools.pypa.io/en/latest/references/keywords.html#keyword-namespace-packages
/usr/share/miniconda/envs/coast/lib/python3.10/site-packages/pkg_resources/__init__.py:2350: DeprecationWarning: Deprecated call to `pkg_resources.declare_namespace('pydap')`.
Implementing implicit namespace packages (as specified in PEP 420) is preferred to `pkg_resources.declare_namespace`. See https://setuptools.pypa.io/en/latest/references/keywords.html#keyword-namespace-packages
/usr/share/miniconda/envs/coast/lib/python3.10/site-packages/pkg_resources/__init__.py:2871: DeprecationWarning: Deprecated call to `pkg_resources.declare_namespace('pydap.tests')`.
Implementing implicit namespace packages (as specified in PEP 420) is preferred to `pkg_resources.declare_namespace`. See https://setuptools.pypa.io/en/latest/references/keywords.html#keyword-namespace-packages
/usr/share/miniconda/envs/coast/lib/python3.10/site-packages/pkg_resources/__init__.py:2350: DeprecationWarning: Deprecated call to `pkg_resources.declare_namespace('pydap')`.
Implementing implicit namespace packages (as specified in PEP 420) is preferred to `pkg_resources.declare_namespace`. See https://setuptools.pypa.io/en/latest/references/keywords.html#keyword-namespace-packages
/usr/share/miniconda/envs/coast/lib/python3.10/site-packages/pkg_resources/__init__.py:2871: DeprecationWarning: Deprecated call to `pkg_resources.declare_namespace('sphinxcontrib')`.
Implementing implicit namespace packages (as specified in PEP 420) is preferred to `pkg_resources.declare_namespace`. See https://setuptools.pypa.io/en/latest/references/keywords.html#keyword-namespace-packages


### Usage of coast._utils.plot_util.velocity_grid_to_geo()

Plotting velocities with curvilinear grid and/or on a polar stereographic projection.


```python
root = "./"
# Paths to a single or multiple data files.
dn_files = root + "example_files/"
fn_nemo_dat_t = dn_files + "coast_nemo_quiver_thetao.nc"
fn_nemo_dat_u = dn_files + "coast_nemo_quiver_uo.nc"
fn_nemo_dat_v = dn_files + "coast_nemo_quiver_vo.nc"
fn_nemo_config_t = root + "config/gc31_nemo_grid_t.json"
fn_nemo_config_u = root + "config/gc31_nemo_grid_u.json"
fn_nemo_config_v = root + "config/gc31_nemo_grid_v.json"

# Set path for domain file if required.
fn_nemo_dom = dn_files + "coast_nemo_quiver_dom.nc"

# Define output filepath (optional: None or str)
fn_out = './quiver_plot.png'
```


```python
# Read in multiyear data (This example uses NEMO data from a single file.)
nemo_data_t = coast.Gridded(fn_data=fn_nemo_dat_t,
fn_domain=fn_nemo_dom,
config=fn_nemo_config_t,
).dataset
nemo_data_u = coast.Gridded(fn_data=fn_nemo_dat_u,
fn_domain=fn_nemo_dom,
config=fn_nemo_config_u,
).dataset
nemo_data_v = coast.Gridded(fn_data=fn_nemo_dat_v,
fn_domain=fn_nemo_dom,
config=fn_nemo_config_v,
).dataset
```

Select surface u and v as an example:


```python
# Select specific data variables.
data_u = nemo_data_u[["u_velocity"]]
data_v = nemo_data_v[["v_velocity"]]

# Select one time step and surface currents
data_u = data_u.isel(t_dim=0, z_dim=0)
data_v = data_v.isel(t_dim=0, z_dim=0)

# Calculate speed
speed = ((data_u.to_array().values[0, :, :] ** 2 + data_v.to_array().values[0, :, :] ** 2) ** 0.5)
```

Calculate adjustment for the curvilinear grid

This function may take a while


```python
uv_velocities = [data_u.to_array().values[0, :, :], data_v.to_array().values[0, :, :]]
u_new, v_new = coast._utils.plot_util.velocity_grid_to_geo(
nemo_data_t.longitude.values, nemo_data_t.latitude.values,
uv_velocities,
polar_stereo_cartopy_bug_fix=False)
```

100%|██████████| 59/59 [06:59<00:00, 7.11s/it]



```python
# Apply the CartoPy stereographic polar correction
# NOTE: This could have been applied automatically with `polar_stereo=True` in
# coast._utils.plot_util.velocity_grid_to_geo()

u_pol, v_pol = coast._utils.plot_util.velocity_polar_bug_fix(u_new, v_new, nemo_data_t.latitude.values)
```


```python
# Set things up for plotting North Pole stereographic projection

# Data projection
data_crs = ccrs.PlateCarree()
# Plot projection
mrc = ccrs.NorthPolarStereo(central_longitude=0.0)
```

Below shows the u and v velocities when plotted with and without adjustment.

The plot shows three cases: 1:no correction, 2:the NEMO grid correction, and 3:the NEMO grid correction with polar stereographic plot correction. We also plot the final corrected u and v velocities as streamlines.
- In case 1, the lower latitude velocities aren't too bad but become irregular further north as the grid lines deviate form lat and lon.
- In case 2, the irregularity still persists even with the grid correction. This is the result of a CartoPy bug which also worsens at high latitudes.
- In case 3, both corrections have been applied and the velocities quivers now align with the route of strong current speed as would be expected.

`!WARNING: UNCOMMENT THE LINES BELOW TO PLOT THE MAPS`


```python
# Subplot axes settings
n_r = 2 # Number of subplot rows
n_c = 2 # Number of subplot columns
figsize = (10, 10) # Figure size
subplot_padding = 0.5 # Amount of vertical and horizontal padding between plots
fig_pad = (0.075, 0.075, 0.1, 0.1) # Figure padding (left, top, right, bottom)

# Labels and Titles
fig_title = "Velocity Plot" # Whole figure title

# Create plot and flatten axis array
fig, ax = plt.subplots(n_r, n_c, subplot_kw={"projection": mrc}, sharey=True, sharex=True, figsize=figsize)
cax = fig.add_axes([0.3, 0.96, 0.4, 0.01])


ax = ax.flatten()
for rr in range(n_r * n_c):
ax[rr].add_feature(cfeature.LAND, zorder=100)
ax[rr].gridlines()
ax[rr].set_extent([-180, 180, 70, 90], crs=data_crs)
coast._utils.plot_util.set_circle(ax[rr])


cs = ax[0].pcolormesh(nemo_data_t.longitude.values, nemo_data_t.latitude.values, speed, transform=data_crs, vmin=0, vmax=0.3)
ax[0].quiver(nemo_data_t.longitude.values, nemo_data_t.latitude.values,
data_u.to_array().values[0, :, :], data_v.to_array().values[0, :, :],
color='w', transform=data_crs, angles='xy', regrid_shape=40)

ax[1].pcolormesh(nemo_data_t.longitude.values, nemo_data_t.latitude.values, speed, transform=data_crs, vmin=0, vmax=0.3)
ax[1].quiver(nemo_data_t.longitude.values, nemo_data_t.latitude.values,
u_new, v_new,
color='w', transform=data_crs, angles='xy', regrid_shape=40)

ax[2].pcolormesh(nemo_data_t.longitude.values, nemo_data_t.latitude.values, speed, transform=data_crs, vmin=0, vmax=0.3)
ax[2].quiver(nemo_data_t.longitude.values, nemo_data_t.latitude.values,
u_pol, v_pol,
color='w', transform=data_crs, angles='xy', regrid_shape=40)

ax[3].pcolormesh(nemo_data_t.longitude.values, nemo_data_t.latitude.values, speed, transform=data_crs, vmin=0, vmax=0.3)
ax[3].streamplot(nemo_data_t.longitude.values, nemo_data_t.latitude.values,
u_pol, v_pol, transform=data_crs, linewidth=1, density=2, color='w', zorder=101)

ax[0].set_title('1: No Correction')
ax[1].set_title('2: NEMO Grid Correction')
ax[2].set_title('3: Grid Correction and Polar Correction')
ax[3].set_title('As left but with streamlines')

fig.colorbar(cs, cax=cax, orientation='horizontal')
cax.set_xlabel('U (m s$^{-1}$)')

#fig.tight_layout(w_pad=subplot_padding, h_pad=subplot_padding)
#fig.subplots_adjust(left=(fig_pad[0]), bottom=(fig_pad[1]), right=(1 - fig_pad[2]), top=(1 - fig_pad[3]))

plt.show()
# uncomment this line to save an output image
# fig.savefig(fn_out)
```
![png](/COAsT/polar_plotting_files/polar_plotting_14_1.png)

Below shows the temperature plots


```python
# Data to plot
data_bathy = nemo_data_t.bathymetry.values
data_temp = nemo_data_t.temperature.isel(t_dim=0, z_dim=0)
```

`!WARNING: UNCOMMENT THE LINES BELOW TO PLOT THE MAPS`


```python
figsize = (5, 5) # Figure size
fig = plt.figure(figsize=figsize)
ax1 = fig.add_axes([0.1, 0.1, 0.8, 0.75], projection=mrc)
cax = fig.add_axes([0.3, 0.96, 0.4, 0.01])

ax1.add_feature(cfeature.LAND, zorder=105)
ax1.gridlines()
ax1.set_extent([-180, 180, 70, 90], crs=data_crs)
coast._utils.plot_util.set_circle(ax1)

# We use to a function to re-project the data for plotting contours over the pole.
cs1 = coast._utils.plot_util.plot_polar_contour(
nemo_data_t.longitude.values, nemo_data_t.latitude.values, data_bathy, ax1, levels=6, colors="k", zorder=101
)
cs2 = ax1.pcolormesh(
nemo_data_t.longitude.values, nemo_data_t.latitude.values, data_temp, transform=data_crs, vmin=-2, vmax=8
)
cax.set_xlabel(r"SST ($^{\circ}$C)")
fig.colorbar(cs2, cax=cax, orientation="horizontal")

# fig.savefig(fn_out, dpi=120)

```
![png](/COAsT/polar_plotting_files/polar_plotting_17_3.png)


```python

```
Loading

0 comments on commit 14d7792

Please sign in to comment.