This post is a quick example on how to use the World Ocean Atlas database.

Recently I found myself re-creating several of the figures I use in my slides for the Introduction to Physical Oceanography course. My final goal is to create a repository with all the slides for the course using a Creative Commons License, so I need all the figures to be either re-done or downloaded from a source that provides a compatible license.

The two main functions we will use are woa_subset and plot_section from oceans and ctd respectively.

In [2]:

import os

import numpy as np

from pandas import read_hdf
from mpl_toolkits.basemap import Basemap
from mpl_toolkits.axes_grid1.inset_locator import inset_axes

from ctd import plot_section
from oceans.colormaps import cm
from oceans.datasets import woa_subset

Let's extract a section across the Atlantic at $25.5^{\circ}$ W. Note that we download the data just once, to avoid bugging the servers too much with our requests.

In [3]:

longitude = -25.5

fname = './data/cross_section_salinity.h5'

if os.path.isfile(fname):
    dataset = read_hdf(fname, 'salinity')
else:
    dataset = woa_subset(var='salinity', clim_type='annual',
                            resolution='1deg', levels=slice(0, 40),
                            llcrnrlat=-70.5, urcrnrlat=70.5,
                            llcrnrlon=longitude,
                            urcrnrlon=longitude)
    dataset = dataset['OA Climatology']['annual'].squeeze().T
    dataset.to_hdf(fname, 'salinity')

dataset = dataset.clip(27.5, 37.5)
levels = np.arange(dataset.min().min(), dataset.max().max(), 0.02)
dataset.lat = dataset.columns.values.astype(float)
dataset.lon = len(dataset.lat) * [longitude]

Now, using the plot_section and Basemap's bluemarble, we plot the section. Can you spot the major Water masses there?

In [4]:

fig, ax, cbar = plot_section(dataset, cmap=cm.avhrr, marker=None,
                             levels=levels, figsize=(10, 4))
X = ax.get_xticks()
offset = 1.01
ax.set_xlim(X.min() * offset, X.max() * offset)
new_labels = np.linspace(dataset.lat.min(), dataset.lat.max(), len(X))
new_labels = [u'%i\u00B0' % num for num in new_labels]
ax.set_xticklabels(new_labels)
ax.set_xlabel(u'WOA09 Salinity Section at Longitude %3.1f\u00B0' % longitude)
ax.set_ylabel('Depth [m]')

axin = inset_axes(ax, width="40%", height="40%", loc=4)
inmap = Basemap(projection='ortho', lon_0=longitude, lat_0=0,
                ax=axin, anchor='NE')
inmap.bluemarble()
_ = inmap.plot(dataset.lon, dataset.lat, 'r', latlon=True)

In [5]:

HTML(html)

Out[5]:

This post was written as an IPython notebook. It is available for download or as a static html.

python4oceanographers by Filipe Fernandes is licensed under a Creative Commons Attribution-ShareAlike 4.0 International License.
Based on a work at https://ocefpaf.github.io/.

python4oceanographers

Turning ripples into waves

Downloading and plotting World Ocean Atlas 2009 - Part 01

Comments