Polar stereographic projection#
Antarctica has the bad idea to be located at the South Pole, which makes the mapping of this region a bit different than what we’ve learned so far. However, as you will see, it’s not very difficult and the change is minimal.
First, the imports. They are the same as always, but I removed the figure size defaults:
# Import the tools we are going to need today:
import matplotlib.pyplot as plt # plotting library
import numpy as np # numerical library
import xarray as xr # netCDF library
import cartopy # Map projections libary
import cartopy.crs as ccrs # Projections list
Read and select the regional data#
Reading the data works as always:
ds = xr.open_dataset('../data/ERA5_LowRes_Invariant.nc')
To select the data for a specific region, we will use xray’s sel
function as we learned in the exercises. For now, we will select the data below 60° south, but you are free to make the domain bigger if you find it useful for your analyses:
ds = ds.sel(latitude=slice(-60, -90))
Read the variable:
z = ds.z / 9.81
Plot the data#
Until now, nothing special in comparison to other regions of the world. It’s for the plotting that things become a bit more complicated. Let’s do it the “usual” way:
# Try a regional plot as usual
ax = plt.axes(projection=ccrs.PlateCarree())
z.plot(ax=ax, transform=ccrs.PlateCarree(), vmin=0, cmap='terrain')
ax.coastlines(); # Add coastlines to the plot
OK. Not very useful. We will now use a better projection:
# The new projection:
ax = plt.axes(projection=ccrs.SouthPolarStereo())
# Limit the map to -60 degrees latitude and below:
ax.set_extent([-180, 180, -90, -60], ccrs.PlateCarree())
# Plot the data as usual
z.plot(ax=ax, transform=ccrs.PlateCarree(), vmin=0, cmap='terrain')
# Add details to the plot:
ax.coastlines();
ax.gridlines();
Looks better now!
Change some details#
The plot above looks fine for me. If you want you can change some details, starting with its size:
# Prepare the figure with the wanted size:
fig = plt.figure(figsize=(9, 7))
# The rest doesn't change:
ax = plt.axes(projection=ccrs.SouthPolarStereo())
ax.set_extent([-180, 180, -90, -60], ccrs.PlateCarree())
z.plot(ax=ax, transform=ccrs.PlateCarree(), vmin=0, cmap='terrain')
ax.coastlines();
ax.gridlines();
If you want, you can also make a circular plot instead of the quadratic one. First, you’ll have to run these few lines (only once for the notebook!):
import matplotlib.path as mpath
theta = np.linspace(0, 2*np.pi, 100)
map_circle = mpath.Path(np.vstack([np.sin(theta), np.cos(theta)]).T * 0.5 + [0.5, 0.5])
And then add one line to the plot commands:
# This did not change:
fig = plt.figure(figsize=(9, 7))
ax = plt.axes(projection=ccrs.SouthPolarStereo())
ax.set_extent([-180, 180, -90, -60], ccrs.PlateCarree())
# Add the following line:
ax.set_boundary(map_circle, transform=ax.transAxes)
# This did not change either:
z.plot(ax=ax, transform=ccrs.PlateCarree(), vmin=0, cmap='terrain')
ax.coastlines();
ax.gridlines();
OK, you should be good now!
Tired of writing so many lines?#
Note that it is possible to simplify your plotting commands by writing a function:
def prepare_plot(figsize=(9, 7)):
"""This function returns prepared axes for the polar plot.
Usage:
fig, ax = prepare_plot()
"""
fig = plt.figure(figsize=figsize)
ax = plt.axes(projection=ccrs.SouthPolarStereo())
ax.set_extent([-180, 180, -90, -60], ccrs.PlateCarree())
ax.set_boundary(map_circle, transform=ax.transAxes)
ax.coastlines(); ax.gridlines();
return fig, ax
Now, making a plot has become even easier:
fig, ax = prepare_plot()
z.plot(ax=ax, transform=ccrs.PlateCarree(), vmin=0, cmap='terrain');
Need to save the plot for your presentation?#
The easiest way is to use “right-click -> save as” on the image in the notebook.
Also, you can save the plot as pdf or png quite easily (examples below). But this might look quite different as the picture on screen sometimes…
fig, ax = prepare_plot()
z.plot(ax=ax, transform=ccrs.PlateCarree(), vmin=0, cmap='terrain');
plt.savefig('topo.pdf')
fig, ax = prepare_plot()
z.plot(ax=ax, transform=ccrs.PlateCarree(), vmin=0, cmap='terrain');
plt.savefig('topo.png')