USEAST -- glider

In [1]:

    

import warnings
import iris

url = ('http://tds.marine.rutgers.edu/thredds/dodsC/projects/wilkin/'
       'mab-glider/Gridded/20130911T000000_20130920T000000_gp2013_modena.nc')

with warnings.catch_warnings():
    warnings.simplefilter('ignore')
    glider = iris.load(url)


x = glider.extract_strict('Longitude').data
y = glider.extract_strict('Latitude').data
z = glider.extract_strict('Temperature').coord('depth').points
t = glider.extract_strict('Temperature').coord('time')
t = t.units.num2date(t.points)
s = glider.extract_strict('Temperature').data


    

In [2]:

    

url = ('http://ecowatch.ncddc.noaa.gov/thredds/dodsC/'
       'ncom_us_east_agg/US_East_Apr_05_2013_to_Current_best.ncd')

with warnings.catch_warnings():
    warnings.simplefilter('ignore')
    cube = iris.load_cube(url, 'sea_water_temperature')


    

In [3]:

    

cube


    

    
Out[3]:





<iris 'Cube' of sea_water_temperature / (degC) (time: 12597; depth: 40; latitude: 664; longitude: 544)>

In [4]:

    

from oceans import wrap_lon360
from ioos_tools.tardis import time_slice

bbox = wrap_lon360(min(x))[0], wrap_lon360(max(x))[0], min(y), max(y)

dx = dy = 0.08
cube = cube.intersection(longitude=(bbox[0]-dx, bbox[1]+dx),
                         latitude=(bbox[2]-dy, bbox[3]+dy))
cube = time_slice(cube, start=t[0], stop=t[-1])


    

In [5]:

    

cube


    

    
Out[5]:





<iris 'Cube' of sea_water_temperature / (degC) (time: 69; depth: 40; latitude: 17; longitude: 18)>

In [6]:

    

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


def distance(x, y, units='km'):
    dist, pha = sw.dist(x, y, units=units)
    return np.r_[0, np.cumsum(dist)]


def plot_glider(x, y, z, t, data, cmap=plt.cm.viridis,
                figsize=(15, 3), track_inset=False):

    fig, ax = plt.subplots(figsize=figsize)
    dist = distance(x, y, units='km')
    z = np.abs(z)
    dist, z = np.broadcast_arrays(dist[..., np.newaxis], z)
    cs = ax.pcolor(dist, z, data, cmap=cmap, snap=True)
    kw = dict(orientation='vertical', extend='both', shrink=0.65)
    fig.colorbar(cs, **kw)

    if track_inset:
        axin = inset_axes(ax, width="25%", height="30%", loc=4)
        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.set_xlabel('Distance (km)')
    ax.set_ylabel('Depth (m)')
    title = '{0:%Y-%m-%d %H:%M:%S} -- {1:%Y-%m-%d %H:%M:%S}'.format
    ax.set_title(title(t[0], t[-1]).format(t[0], t[-1]))
    return fig, ax


    

In [7]:

    

%matplotlib inline
import matplotlib.pyplot as plt

import cartopy.crs as ccrs
import cartopy.feature as cfeature
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER

LAND = cfeature.NaturalEarthFeature('physical', 'land', '10m',
                                    edgecolor='face',
                                    facecolor=cfeature.COLORS['land'])


def plot_cube(cube, projection=ccrs.PlateCarree()):
    fig, ax = plt.subplots(figsize=(9, 13),
                           subplot_kw=dict(projection=projection))
    gl = ax.gridlines(draw_labels=True)
    gl.xlabels_top = gl.ylabels_right = False
    gl.xformatter = LONGITUDE_FORMATTER
    gl.yformatter = LATITUDE_FORMATTER

    lon = cube.coord(axis='X').points
    lat = cube.coord(axis='Y').points
    extent = [lon.min(), lon.max(),
              lat.min(), lat.max()]
    ax.set_extent(extent)
    ax.add_feature(LAND, zorder=0)
    ax.coastlines('10m', zorder=1)
    return fig, ax


    

In [8]:

    

from gridgeo import GridGeo
from oceans import wrap_lon180

fig, ax = plot_cube(cube)

# Grid outline.
grid = GridGeo(cube)
xb, yb = grid.outline.boundary.xy
polygon = zip(xb, yb)

lon = wrap_lon180(cube.coord(axis='X').points)
lat = cube.coord(axis='Y').points
lon, lat = np.meshgrid(lon, lat)
ax.plot(lon, lat, 'k.')

# Glider track.
ax.plot(wrap_lon180(x), y, 'r')
ax.plot(wrap_lon180(x[0]), y[0], 'go')
ax.plot(wrap_lon180(x[-1]), y[-1], 'ro')


    

    
Out[8]:





[<matplotlib.lines.Line2D at 0x7fe5649e25f8>]






    

In [9]:

    

import cmocean

fig, ax = plot_cube(cube)

# Glider track.
ax.plot(wrap_lon180(x), y, 'r')
ax.plot(wrap_lon180(x[0]), y[0], 'go')
ax.plot(wrap_lon180(x[-1]), y[-1], 'ro')


# Variable at surface and last time step.
# (There is no _FillValue/missing_value!)
temp = ma.masked_greater_equal(cube[-1, 0, ...].data, 9.99999993e+36)
kw = dict(orientation='vertical', extend='both', shrink=0.45)
cs = ax.pcolormesh(lon, lat, temp,
                   cmap=cmocean.cm.thermal, zorder=0)
cbar = fig.colorbar(cs, **kw)


    

In [10]:

    

from iris.analysis.trajectory import interpolate

track = [('longitude', wrap_lon180(x)),
         ('latitude', y),
         ('time', t)]

vgliver = interpolate(cube, track, method='linear')


    

    




/home/filipe/miniconda3/envs/IOOS/lib/python3.6/site-packages/iris/fileformats/netcdf.py:394: RuntimeWarning: invalid value encountered in greater
  data = variable[keys]

In [11]:

    

vertical = cube.coord('depth').points

fig, ax = plt.subplots(figsize=(15, 3))
ax.pcolormesh(t, vertical, np.ma.masked_invalid(vgliver.data),
              cmap=cmocean.cm.haline)
ax.invert_yaxis()
ax.set_ylim(80, 0)


    

    
Out[11]:





(80, 0)






    

In [12]:

    

fig, ax = plot_glider(x, y, z, t, s,
                      cmap=cmocean.cm.haline, track_inset=True)
ax.set_ylim(80, 0)


    

    
Out[12]:





(80, 0)