Accessing glider data via the Glider DAC API with Python

IOOS provides an API for getting information on all the glider deployments available in the Glider DAC.

The raw JSON can be accessed at https://data.ioos.us/gliders/providers/api/deployment and it is quite simple to parse it with Python.

First, lets check how many glider deployments exist in the Glider DAC.

import requests


url = "http://data.ioos.us/gliders/providers/api/deployment"

response = requests.get(url)

res = response.json()

print(f"Found {res['num_results']} deployments!")
Found 994 deployments!

And here is the JSON of the last deployment found in the list.

deployments = res["results"]

deployment = deployments[-1]

deployment
{'archive_safe': True,
 'attribution': '',
 'checksum': '5e7eaf93157096305d5a0b7e78b38b2a',
 'completed': True,
 'compliance_check_passed': False,
 'created': 1620316396148,
 'dap': 'http://gliders.ioos.us/thredds/dodsC/deployments/castelao/sp025-20210505T2104/sp025-20210505T2104.nc3.nc',
 'delayed_mode': False,
 'deployment_date': 1620248640000,
 'deployment_dir': 'castelao/sp025-20210505T2104',
 'erddap': 'http://gliders.ioos.us/erddap/tabledap/sp025-20210505T2104.html',
 'estimated_deploy_date': None,
 'estimated_deploy_location': None,
 'glider_name': 'sp025',
 'id': '609410ecd51c2400b8d73cfd',
 'iso': 'http://gliders.ioos.us/erddap/tabledap/sp025-20210505T2104.iso19115',
 'latest_file': 'sp025-20210818T092500.nc',
 'latest_file_mtime': 1629300519457,
 'name': 'sp025-20210505T2104',
 'operator': 'Scripps Institution of Oceanography',
 'sos': 'http://gliders.ioos.us/thredds/sos/deployments/castelao/sp025-20210505T2104/sp025-20210505T2104.nc3.nc?service=SOS&request=GetCapabilities&AcceptVersions=1.0.0',
 'thredds': 'http://gliders.ioos.us/thredds/catalog/deployments/castelao/sp025-20210505T2104/catalog.html?dataset=deployments/castelao/sp025-20210505T2104/sp025-20210505T2104.nc3.nc',
 'updated': 1629858501105,
 'username': 'castelao',
 'wmo_id': '4801900'}

The metadata is very rich and informative. A quick way to get to the data is to read dap endpoint with iris.

import iris


# Get this specific glider because it looks cool ;-)
for deployment in deployments:
    if deployment["name"] == "sp064-20161214T1913":
        url = deployment["dap"]

# See https://github.com/Unidata/netcdf-c/issues/1299 for the explanation of `#fillmismatch`.
cubes = iris.load_raw("".join([url, "#fillmismatch"]))

print(cubes)
0: Profile ID / (unknown)              (-- : 1; -- : 931)
1: sea_water_density / (kg m-3)        (-- : 1; -- : 931; -- : 367)
2: precise_lon Variable Quality Flag / (no_unit) (-- : 1; -- : 931; -- : 367)
3: time status_flag / (no_unit)        (-- : 1; -- : 931; -- : 367)
4: Trajectory Name / (unknown)         (-- : 1; -- : 64)
5: precise_lat Variable Quality Flag / (no_unit) (-- : 1; -- : 931; -- : 367)
6: longitude / (degrees)               (-- : 1; -- : 931; -- : 367)
7: time / (seconds since 1970-01-01T00:00:00Z) (-- : 1; -- : 931; -- : 367)
8: sea_water_practical_salinity / (1)  (-- : 1; -- : 931; -- : 367)
9: sea_water_temperature / (Celsius)   (-- : 1; -- : 931; -- : 367)
10: sea_water_electrical_conductivity / (S m-1) (-- : 1; -- : 931; -- : 367)
11: latitude status_flag / (no_unit)    (-- : 1; -- : 931; -- : 367)
12: CTD Metadata / (1)                  (-- : 1; -- : 931; -- : 367)
13: longitude status_flag / (no_unit)   (-- : 1; -- : 931; -- : 367)
14: sea_water_pressure / (dbar)         (-- : 1; -- : 931; -- : 367)
15: eastward_sea_water_velocity / (m s-1) (-- : 1; -- : 931)
16: northward_sea_water_velocity / (m s-1) (-- : 1; -- : 931)
17: WMO ID / (unknown)                  (-- : 1; -- : 64)
18: latitude / (degrees)                (-- : 1; -- : 931; -- : 367)
19: Platform Metadata / (1)             (-- : 1; -- : 931; -- : 367)
20: northward_sea_water_velocity status_flag / (no_unit) (-- : 1; -- : 931; -- : 367)
21: time status_flag / (no_unit)        (-- : 1; -- : 931; -- : 367)
22: time status_flag / (no_unit)        (-- : 1; -- : 931; -- : 367)
23: time status_flag / (no_unit)        (-- : 1; -- : 931; -- : 367)
24: eastward_sea_water_velocity status_flag / (no_unit) (-- : 1; -- : 931; -- : 367)
/home/filipe/miniconda3/envs/IOOS/lib/python3.9/site-packages/iris/fileformats/cf.py:285: UserWarning: Missing CF-netCDF ancillary data variable 'profile_time_qc', referenced by netCDF variable 'time'
  warnings.warn(message % (name, nc_var_name))
/home/filipe/miniconda3/envs/IOOS/lib/python3.9/site-packages/iris/fileformats/cf.py:285: UserWarning: Missing CF-netCDF ancillary data variable 'profile_lat_qc', referenced by netCDF variable 'latitude'
  warnings.warn(message % (name, nc_var_name))
/home/filipe/miniconda3/envs/IOOS/lib/python3.9/site-packages/iris/fileformats/cf.py:285: UserWarning: Missing CF-netCDF ancillary data variable 'profile_lon_qc', referenced by netCDF variable 'longitude'
  warnings.warn(message % (name, nc_var_name))
/home/filipe/miniconda3/envs/IOOS/lib/python3.9/site-packages/iris/fileformats/cf.py:285: UserWarning: Missing CF-netCDF ancillary data variable 'lat_qc', referenced by netCDF variable 'precise_lat'
  warnings.warn(message % (name, nc_var_name))
/home/filipe/miniconda3/envs/IOOS/lib/python3.9/site-packages/iris/fileformats/cf.py:285: UserWarning: Missing CF-netCDF ancillary data variable 'lon_qc', referenced by netCDF variable 'precise_lon'
  warnings.warn(message % (name, nc_var_name))
/home/filipe/miniconda3/envs/IOOS/lib/python3.9/site-packages/iris/fileformats/cf.py:1156: UserWarning: Ignoring variable 'lon_uv_qc' referenced by variable 'lon_uv': Dimensions ('trajectory', 'profile', 'obs') do not span ('trajectory', 'profile')
  warnings.warn(msg)
/home/filipe/miniconda3/envs/IOOS/lib/python3.9/site-packages/iris/fileformats/cf.py:1156: UserWarning: Ignoring variable 'lat_uv_qc' referenced by variable 'lat_uv': Dimensions ('trajectory', 'profile', 'obs') do not span ('trajectory', 'profile')
  warnings.warn(msg)
/home/filipe/miniconda3/envs/IOOS/lib/python3.9/site-packages/iris/fileformats/cf.py:1156: UserWarning: Ignoring variable 'time_uv_qc' referenced by variable 'time_uv': Dimensions ('trajectory', 'profile', 'obs') do not span ('trajectory', 'profile')
  warnings.warn(msg)
/home/filipe/miniconda3/envs/IOOS/lib/python3.9/site-packages/iris/fileformats/cf.py:1156: UserWarning: Ignoring variable 'u_qc' referenced by variable 'u': Dimensions ('trajectory', 'profile', 'obs') do not span ('trajectory', 'profile')
  warnings.warn(msg)
/home/filipe/miniconda3/envs/IOOS/lib/python3.9/site-packages/iris/fileformats/cf.py:1156: UserWarning: Ignoring variable 'v_qc' referenced by variable 'v': Dimensions ('trajectory', 'profile', 'obs') do not span ('trajectory', 'profile')
  warnings.warn(msg)

In order to plot, for example sea water temperature data, one must clean the data first for missing values

import matplotlib.pyplot as plt
import numpy as np
import numpy.ma as ma
import seawater as sw
from mpl_toolkits.axes_grid1.inset_locator import inset_axes


def distance(x, y, units="km"):
    if ma.isMaskedArray(x):
        x = x.filled(fill_value=np.NaN)
    if ma.isMaskedArray(y):
        y = y.filled(fill_value=np.NaN)
    dist, pha = sw.dist(x, y, units=units)
    return np.r_[0, np.cumsum(dist)]


def apply_range(cube_coord):
    if isinstance(cube_coord, iris.cube.Cube):
        data = cube_coord.data.squeeze()
    elif isinstance(cube_coord, (iris.coords.AuxCoord, iris.coords.Coord)):
        data = cube_coord.points.squeeze()

    actual_range = cube_coord.attributes.get("actual_range")
    if actual_range is not None:
        vmin, vmax = actual_range
        data = ma.masked_outside(data, vmin, vmax)
    return data


def plot_glider(cube, cmap=plt.cm.viridis, figsize=(9, 3.75), track_inset=False):

    data = apply_range(cube)
    x = apply_range(cube.coord(axis="X"))
    y = apply_range(cube.coord(axis="Y"))
    z = apply_range(cube.coord(axis="Z"))
    t = cube.coord(axis="T")
    t = t.units.num2date(t.points.squeeze())

    fig, ax = plt.subplots(figsize=figsize)
    dist = distance(x, y)
    z = ma.abs(z)
    dist, _ = np.broadcast_arrays(dist[..., np.newaxis], z.filled(fill_value=np.NaN))
    dist, z = map(ma.masked_invalid, (dist, z))
    cs = ax.pcolor(dist, z, data, cmap=cmap, snap=True)
    kw = dict(orientation="horizontal", extend="both", shrink=0.65)
    cbar = fig.colorbar(cs, **kw)

    if track_inset:
        axin = inset_axes(
            ax,
            width=2,
            height=2,
            loc=4,
            bbox_to_anchor=(1.15, 0.35),
            bbox_transform=ax.figure.transFigure,
        )
        axin.plot(x, y, "k.")
        start, end = (x[0], y[0]), (x[-1], y[-1])
        kw = dict(marker="o", linestyle="none")
        axin.plot(*start, color="g", **kw)
        axin.plot(*end, color="r", **kw)
        axin.axis("off")

    ax.invert_yaxis()
    ax.invert_xaxis()
    ax.set_xlabel("Distance (km)")
    ax.set_ylabel("Depth (m)")
    return fig, ax, cbar

The functions above apply the actual_range metadata to the data, mask the invalid/bad values, and prepare the parameters for plotting.

The figure below shows the temperature slice (left), and glider track (right) with start and end points marked with green and red respectively.

Note: This glider was deployed off the west of the U.S.

%matplotlib inline

temp = cubes.extract_cube("sea_water_temperature")

fig, ax, cbar = plot_glider(
    temp, cmap=plt.cm.viridis, figsize=(9, 4.25), track_inset=True
)
/tmp/ipykernel_84279/1378783559.py:44: MatplotlibDeprecationWarning: shading='flat' when X and Y have the same dimensions as C is deprecated since 3.3.  Either specify the corners of the quadrilaterals with X and Y, or pass shading='auto', 'nearest' or 'gouraud', or set rcParams['pcolor.shading'].  This will become an error two minor releases later.
  cs = ax.pcolor(dist, z, data, cmap=cmap, snap=True)
../../../_images/2016-12-20-searching_glider_deployments_9_1.png

There are many things the user can do with the API. Here is another example that finds all glider deployments within a boundary box.

bbox = [[-125.72, 32.60], [-117.57, 36.93]]

The cell below defines two helper functions to parse the geometry from the JSON and convert the trajectory to a shapely LineString to prepare the data for GIS operations later.

from shapely.geometry import LineString


def parse_geometry(geometry):
    """
    Filters out potentially bad coordinate pairs as returned from
    GliderDAC. Returns a safe geometry object.

    :param dict geometry: A GeoJSON Geometry object

    """
    coords = []
    for lon, lat in geometry["coordinates"]:
        if lon is None or lat is None:
            continue
        coords.append([lon, lat])
    return {"coordinates": coords}


def fetch_trajectory(deployment):
    """
    Downloads the track as GeoJSON from GliderDAC

    :param dict deployment: The deployment object as returned from GliderDAC

    """
    track_url = "http://data.ioos.us/gliders/status/api/track/{}".format
    response = requests.get(track_url(deployment["deployment_dir"]))
    if response.status_code != 200:
        raise IOError(
            "Failed to get Glider Track for %s" % deployment["deployment_dir"]
        )
    geometry = parse_geometry(response.json())
    coords = LineString(geometry["coordinates"])
    return coords

Now it is easy to check which tracks lie inside the box.

res = response.json()["results"]
len(res[-100:])
100
from shapely.geometry import box


search_box = box(bbox[0][0], bbox[0][1], bbox[1][0], bbox[1][1])

inside = dict()
# Getting only the 10 deployments.
for deployment in response.json()["results"][-10:]:
    try:
        coords = fetch_trajectory(deployment)
    except IOError:
        continue
    if search_box.intersects(coords):
        inside.update({deployment["name"]: coords})

Finally, we can create an interactive map displaying the tracks found in the bounding box.

def plot_track(coords, name, color="orange"):
    x, y = coords.xy
    locations = list(zip(y.tolist(), x.tolist()))
    kw = {"fill": True, "radius": 10, "stroke": False}
    folium.CircleMarker(locations[0], color="green", **kw).add_to(m)
    folium.CircleMarker(locations[-1], color="red", **kw).add_to(m)

    folium.PolyLine(
        locations=locations, color=color, weight=8, opacity=0.2, popup=name
    ).add_to(m)
import folium

tiles = (
    "http://services.arcgisonline.com/arcgis/rest/services/"
    "World_Topo_Map/MapServer/MapServer/tile/{z}/{y}/{x}"
)

location = [search_box.centroid.y, search_box.centroid.x]

m = folium.Map(location=location, zoom_start=5, tiles=tiles, attr="ESRI")


for name, coords in inside.items():
    plot_track(coords, name, color="orange")


m
Make this Notebook Trusted to load map: File -> Trust Notebook