python4oceanographers

Turning ripples into waves

Packaging all my iris workarounds

After copying-and-pasting so many of my iris code snippets I decided to package everything in a new module: TARDIS.

Tardis is nothing more than a hackish wrapper around iris that expands some of the iris methods to work with 2D-Coords and to facilitate some of its methods calls.

Virtually everything that tardis does could be accomplished using the raw iris API, but it will be a little bit easier with tardis ;-)

For example, how to load cubes based on a variable list of CF standard names?

In [3]:
from tardis import quick_load_cubes


name_list = ['sea_surface_height',
             'sea_surface_elevation',
             'sea_surface_height_above_geoid',
             'sea_surface_height_above_sea_level',
             'water_surface_height_above_reference_datum',
             'sea_surface_height_above_reference_ellipsoid']

url = ("http://tds.marine.rutgers.edu/thredds/dodsC/roms/espresso/2013_da/his_Best/"
       "ESPRESSO_Real-Time_v2_History_Best_Available_best.ncd")

cube = quick_load_cubes(url, name_list, callback=None, strict=True)
/home/filipe/.virtualenvs/blog/lib/python2.7/site-packages/iris/fileformats/cf.py:1032: UserWarning: Ignoring formula terms variable u'h' referenced by data variable u'v' via variable u's_rho': Dimensions (u'eta_rho', u'xi_rho') do not span (u'time', u's_rho', u'eta_v', u'xi_v')
  warnings.warn(msg)
/home/filipe/.virtualenvs/blog/lib/python2.7/site-packages/iris/fileformats/cf.py:1032: UserWarning: Ignoring formula terms variable u'zeta' referenced by data variable u'v' via variable u's_rho': Dimensions (u'time', u'eta_rho', u'xi_rho') do not span (u'time', u's_rho', u'eta_v', u'xi_v')
  warnings.warn(msg)
/home/filipe/.virtualenvs/blog/lib/python2.7/site-packages/iris/fileformats/cf.py:1032: UserWarning: Ignoring formula terms variable u'h' referenced by data variable u'u' via variable u's_rho': Dimensions (u'eta_rho', u'xi_rho') do not span (u'time', u's_rho', u'eta_u', u'xi_u')
  warnings.warn(msg)
/home/filipe/.virtualenvs/blog/lib/python2.7/site-packages/iris/fileformats/cf.py:1032: UserWarning: Ignoring formula terms variable u'zeta' referenced by data variable u'u' via variable u's_rho': Dimensions (u'time', u'eta_rho', u'xi_rho') do not span (u'time', u's_rho', u'eta_u', u'xi_u')
  warnings.warn(msg)
/home/filipe/.virtualenvs/blog/lib/python2.7/site-packages/iris/fileformats/_pyke_rules/compiled_krb/fc_rules_cf_fc.py:1395: UserWarning: Gracefully filling 'time' dimension coordinate masked points
  warnings.warn(msg.format(str(cf_coord_var.cf_name)))
/home/filipe/.virtualenvs/blog/lib/python2.7/site-packages/iris/fileformats/netcdf.py:443: UserWarning: Unable to find coordinate for variable u'zeta'
  '{!r}'.format(name))
/home/filipe/.virtualenvs/blog/lib/python2.7/site-packages/iris/fileformats/netcdf.py:443: UserWarning: Unable to find coordinate for variable u'h'
  '{!r}'.format(name))
/home/filipe/.virtualenvs/blog/lib/python2.7/site-packages/iris/fileformats/netcdf.py:560: UserWarning: Unable to construct Ocean s-coordinate, generic form 1 factory due to insufficient source coordinates.
  warnings.warn('{}'.format(e))

Aside: Note all those warnings above? That is one of the reasons why it is complicated to deal with all ROMS variables as a single dataset. Hopefully pysgrid will save us =)

Let's check this cube out.

In [4]:
print(cube)
sea_surface_height / (meter)        (time: 18636; -- : 82; -- : 130)
     Dimension coordinates:
          time                           x           -        -
     Auxiliary coordinates:
          forecast_reference_time        x           -        -
          latitude                       -           x        x
          longitude                      -           x        x
     Attributes:
          CPP_options: MyCPP, ADD_FSOBC, ADD_M2OBC, ANA_BSFLUX, ANA_BTFLUX, ASSUMED_SHAPE, AVERAGES,...
          Conventions: CF-1.4, _Coordinates
          DODS_EXTRA.Unlimited_Dimension: ocean_time
          _ChunkSize: [  1  82 130]
          _CoordSysBuilder: ucar.nc2.dataset.conv.CF1Convention
          ana_file: ROMS/Functionals/ana_btflux.h, /home/julia/ROMS/espresso/RealTime/Comp...
          avg_base: espresso_avg_3450
          bry_file: ../Data/espresso_bdry_new.nc
          cdm_data_type: GRID
          clm_file: ../Data/espresso_clm_new.nc
          code_dir: /home/julia/ROMS/espresso/svn1409
          compiler_command: /opt/pgisoft/openmpi/bin/mpif90
          compiler_flags:  -O3 -Mfree
          compiler_system: pgi
          cpu: x86_64
          featureType: GRID
          field: free-surface, scalar, series
          file: espresso_his_3450_0005.nc
          flt_file: espresso_flt_3450.nc
          format: netCDF-4/HDF5 file
          fpos_file: /home/om/roms/espresso/Data/espresso_floats_glider.in
          frc_file_01: /home/om/roms/espresso/Data/espresso_tide_c05_20060101.nc
          frc_file_02: ../Data/espresso_river.nc
          frc_file_03: ../Data/rain_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_04: ../Data/swrad_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_05: ../Data/Tair_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_06: ../Data/Pair_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_07: ../Data/Qair_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_08: ../Data/lwrad_down_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_09: ../Data/Uwind_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_10: ../Data/Vwind_ncepnam_3hourly_MAB_and_GoM.nc
          grd_file: /home/om/roms/espresso/Data/espresso_grid_c05.nc
          header_dir: /home/julia/ROMS/espresso/RealTime/Compile/fwd
          header_file: espresso.h
          his_base: espresso_his_3450
          history: ROMS/TOMS, Version 3.5, Tuesday - June 16, 2015 -  6:22:08 AM ;
FMRC Best...
          ini_file: /home/julia/ROMS/espresso/RealTime/Storage/run04/espresso_ini_3450.nc
          location: Proto fmrc:espresso_2013_da_his_best
          os: Linux
          rst_file: espresso_rst_3450.nc
          script_file: nl_ocean_espresso.in
          summary: Operational nowcast/forecast system version 2 for MARACOOS project (http://maracoos.org)....
          svn_rev: exported
          svn_url: https://www.myroms.org/svn/src/trunk
          tiling: 004x002
          time: ocean_time
          title: ROMS ESPRESSO Real-Time Operational IS4DVAR Forecast System Version 2 (NEW)...
          type: ROMS/TOMS history file

The module is called tardis because it allows you to "travel" in space, time, and (variables) dimension.

The workhorse function is proc_cube. This function takes a cube and a series of optional constrains to "conform" a cube into the desired time span, space bounding box, and/or automatic units conversion.

In [5]:
import iris
import pytz
from datetime import datetime, timedelta

from tardis import proc_cube


# Time.
stop = datetime(2014, 7, 7, 12)
stop = stop.replace(tzinfo=pytz.utc)
start = stop - timedelta(days=1)

# Space.
bbox = [-87.40, 24.25, -74.70, 36.70]

# Units.
units = iris.unit.Unit('meters')

# New cube.
cube = proc_cube(cube, bbox=bbox, time=(start, stop), units=units)

print(cube)
sea_surface_height / (meter)        (time: 24; -- : 75; -- : 38)
     Dimension coordinates:
          time                           x        -        -
     Auxiliary coordinates:
          forecast_reference_time        x        -        -
          latitude                       -        x        x
          longitude                      -        x        x
     Attributes:
          CPP_options: MyCPP, ADD_FSOBC, ADD_M2OBC, ANA_BSFLUX, ANA_BTFLUX, ASSUMED_SHAPE, AVERAGES,...
          Conventions: CF-1.4, _Coordinates
          DODS_EXTRA.Unlimited_Dimension: ocean_time
          _ChunkSize: [  1  82 130]
          _CoordSysBuilder: ucar.nc2.dataset.conv.CF1Convention
          ana_file: ROMS/Functionals/ana_btflux.h, /home/julia/ROMS/espresso/RealTime/Comp...
          avg_base: espresso_avg_3450
          bry_file: ../Data/espresso_bdry_new.nc
          cdm_data_type: GRID
          clm_file: ../Data/espresso_clm_new.nc
          code_dir: /home/julia/ROMS/espresso/svn1409
          compiler_command: /opt/pgisoft/openmpi/bin/mpif90
          compiler_flags:  -O3 -Mfree
          compiler_system: pgi
          cpu: x86_64
          featureType: GRID
          field: free-surface, scalar, series
          file: espresso_his_3450_0005.nc
          flt_file: espresso_flt_3450.nc
          format: netCDF-4/HDF5 file
          fpos_file: /home/om/roms/espresso/Data/espresso_floats_glider.in
          frc_file_01: /home/om/roms/espresso/Data/espresso_tide_c05_20060101.nc
          frc_file_02: ../Data/espresso_river.nc
          frc_file_03: ../Data/rain_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_04: ../Data/swrad_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_05: ../Data/Tair_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_06: ../Data/Pair_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_07: ../Data/Qair_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_08: ../Data/lwrad_down_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_09: ../Data/Uwind_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_10: ../Data/Vwind_ncepnam_3hourly_MAB_and_GoM.nc
          grd_file: /home/om/roms/espresso/Data/espresso_grid_c05.nc
          header_dir: /home/julia/ROMS/espresso/RealTime/Compile/fwd
          header_file: espresso.h
          his_base: espresso_his_3450
          history: ROMS/TOMS, Version 3.5, Tuesday - June 16, 2015 -  6:22:08 AM ;
FMRC Best...
          ini_file: /home/julia/ROMS/espresso/RealTime/Storage/run04/espresso_ini_3450.nc
          location: Proto fmrc:espresso_2013_da_his_best
          os: Linux
          rst_file: espresso_rst_3450.nc
          script_file: nl_ocean_espresso.in
          summary: Operational nowcast/forecast system version 2 for MARACOOS project (http://maracoos.org)....
          svn_rev: exported
          svn_url: https://www.myroms.org/svn/src/trunk
          tiling: 004x002
          time: ocean_time
          title: ROMS ESPRESSO Real-Time Operational IS4DVAR Forecast System Version 2 (NEW)...
          type: ROMS/TOMS history file

One of the advantages of using iris is the ability to save the data in a CF-compliant netCDF file.

In [6]:
iris.save(cube, "subset.nc")
/home/filipe/.virtualenvs/blog/lib/python2.7/site-packages/iris/fileformats/netcdf.py:1785: UserWarning: NetCDF default saving behaviour currently assigns the outermost dimensions to unlimited. This behaviour is to be deprecated, in favour of no automatic assignment. To switch to the new behaviour, set iris.FUTURE.netcdf_no_unlimited to True.
  warnings.warn(msg)

In [7]:
cube = iris.load_cube("subset.nc")
In [8]:
print(cube)
sea_surface_height / (meter)        (time: 24; -- : 75; -- : 38)
     Dimension coordinates:
          time                           x        -        -
     Auxiliary coordinates:
          forecast_reference_time        x        -        -
          latitude                       -        x        x
          longitude                      -        x        x
     Attributes:
          CPP_options: MyCPP, ADD_FSOBC, ADD_M2OBC, ANA_BSFLUX, ANA_BTFLUX, ASSUMED_SHAPE, AVERAGES,...
          Conventions: CF-1.5
          DODS_EXTRA.Unlimited_Dimension: ocean_time
          _ChunkSize: [  1  82 130]
          _CoordSysBuilder: ucar.nc2.dataset.conv.CF1Convention
          ana_file: ROMS/Functionals/ana_btflux.h, /home/julia/ROMS/espresso/RealTime/Comp...
          avg_base: espresso_avg_3450
          bry_file: ../Data/espresso_bdry_new.nc
          cdm_data_type: GRID
          clm_file: ../Data/espresso_clm_new.nc
          code_dir: /home/julia/ROMS/espresso/svn1409
          compiler_command: /opt/pgisoft/openmpi/bin/mpif90
          compiler_flags:  -O3 -Mfree
          compiler_system: pgi
          cpu: x86_64
          featureType: GRID
          field: free-surface, scalar, series
          file: espresso_his_3450_0005.nc
          flt_file: espresso_flt_3450.nc
          format: netCDF-4/HDF5 file
          fpos_file: /home/om/roms/espresso/Data/espresso_floats_glider.in
          frc_file_01: /home/om/roms/espresso/Data/espresso_tide_c05_20060101.nc
          frc_file_02: ../Data/espresso_river.nc
          frc_file_03: ../Data/rain_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_04: ../Data/swrad_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_05: ../Data/Tair_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_06: ../Data/Pair_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_07: ../Data/Qair_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_08: ../Data/lwrad_down_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_09: ../Data/Uwind_ncepnam_3hourly_MAB_and_GoM.nc
          frc_file_10: ../Data/Vwind_ncepnam_3hourly_MAB_and_GoM.nc
          grd_file: /home/om/roms/espresso/Data/espresso_grid_c05.nc
          header_dir: /home/julia/ROMS/espresso/RealTime/Compile/fwd
          header_file: espresso.h
          his_base: espresso_his_3450
          history: ROMS/TOMS, Version 3.5, Tuesday - June 16, 2015 -  6:22:08 AM ;
FMRC Best...
          ini_file: /home/julia/ROMS/espresso/RealTime/Storage/run04/espresso_ini_3450.nc
          location: Proto fmrc:espresso_2013_da_his_best
          os: Linux
          rst_file: espresso_rst_3450.nc
          script_file: nl_ocean_espresso.in
          summary: Operational nowcast/forecast system version 2 for MARACOOS project (http://maracoos.org)....
          svn_rev: exported
          svn_url: https://www.myroms.org/svn/src/trunk
          tiling: 004x002
          time: ocean_time
          title: ROMS ESPRESSO Real-Time Operational IS4DVAR Forecast System Version 2 (NEW)...
          type: ROMS/TOMS history file

or,

In [9]:
!ncdump -h subset.nc
netcdf subset {
dimensions:
	time = UNLIMITED ; // (24 currently)
	dim1 = 75 ;
	dim2 = 38 ;
variables:
	float zeta(time, dim1, dim2) ;
		zeta:_FillValue = 1.e+37f ;
		string zeta:standard_name = "sea_surface_height" ;
		string zeta:long_name = "free-surface" ;
		zeta:units = "meter" ;
		zeta:coordinates = "lat_rho lon_rho time_run" ;
	double time(time) ;
		time:axis = "T" ;
		time:units = "hours since 2013-05-18T00:00:00Z" ;
		string time:standard_name = "time" ;
		string time:long_name = "Forecast time for ForecastModelRunCollection" ;
		time:calendar = "gregorian" ;
		time:_CoordinateAxisType = "Time" ;
	double time_run(time) ;
		time_run:units = "hours since 2013-05-18T00:00:00Z" ;
		string time_run:standard_name = "forecast_reference_time" ;
		string time_run:long_name = "run times for coordinate = time" ;
		time_run:calendar = "gregorian" ;
		time_run:_CoordinateAxisType = "RunTime" ;
	double lat_rho(dim1, dim2) ;
		lat_rho:units = "degrees" ;
		string lat_rho:standard_name = "latitude" ;
		string lat_rho:long_name = "latitude of RHO-points" ;
		lat_rho:_ChunkSize = 82, 130 ;
		lat_rho:_CoordinateAxisType = "Lat" ;
		lat_rho:field = "lat_rho, scalar" ;
	double lon_rho(dim1, dim2) ;
		lon_rho:units = "degrees" ;
		string lon_rho:standard_name = "longitude" ;
		string lon_rho:long_name = "longitude of RHO-points" ;
		lon_rho:_ChunkSize = 82, 130 ;
		lon_rho:_CoordinateAxisType = "Lon" ;
		lon_rho:field = "lon_rho, scalar" ;

// global attributes:
		:CPP_options = "MyCPP, ADD_FSOBC, ADD_M2OBC, ANA_BSFLUX, ANA_BTFLUX, ASSUMED_SHAPE, AVERAGES, BULK_FLUXES, CURVGRID, DEFLATE, DJ_GRADPS, DOUBLE_PRECISION, EAST_FSCHAPMAN, EAST_M2FLATHER, EAST_M3NUDGING, EAST_M3RADIATION, EAST_TNUDGING, EAST_TRADIATION, EMINUSP, FLOATS, FORWARD_WRITE, GLS_MIXING, HDF5, KANTHA_CLAYSON, LONGWAVE_OUT, M3CLIMATOLOGY, M3CLM_NUDGING, MASKING, MIX_GEO_TS, MIX_S_UV, MPI, NONLINEAR, NONLIN_EOS, NORTHERN_WALL, N2S2_HORAVG, POWER_LAW, PROFILE, K_GSCHEME, !RST_SINGLE, SALINITY, SOLAR_SOURCE, SOLVE3D, SOUTH_FSCHAPMAN, SOUTH_M2FLATHER, SOUTH_M3NUDGING, SOUTH_M3RADIATION, SOUTH_TNUDGING, SOUTH_TRADIATION, SPLINES, SSH_TIDES, TCLIMATOLOGY, TCLM_NUDGING, TS_A4HADVECTION, TS_A4VADVECTION, TS_DIF2, TS_PSOURCE, UV_ADV, UV_COR, UV_U3HADVECTION, UV_C4VADVECTION, UV_QDRAG, UV_PSOURCE, UV_TIDES, UV_VIS2, VAR_RHO_2D, WEST_FSCHAPMAN, WEST_M2FLATHER, WEST_M3NUDGING, WEST_M3RADIATION, WEST_TNUDGING, WEST_TRADIATION" ;
		:DODS_EXTRA.Unlimited_Dimension = "ocean_time" ;
		:_ChunkSize = 1, 82, 130 ;
		:_CoordSysBuilder = "ucar.nc2.dataset.conv.CF1Convention" ;
		:ana_file = "ROMS/Functionals/ana_btflux.h, /home/julia/ROMS/espresso/RealTime/Compile/fwd/ana_nudgcoef.h" ;
		:avg_base = "espresso_avg_3450" ;
		:bry_file = "../Data/espresso_bdry_new.nc" ;
		:cdm_data_type = "GRID" ;
		:clm_file = "../Data/espresso_clm_new.nc" ;
		:code_dir = "/home/julia/ROMS/espresso/svn1409" ;
		:compiler_command = "/opt/pgisoft/openmpi/bin/mpif90" ;
		:compiler_flags = " -O3 -Mfree" ;
		:compiler_system = "pgi" ;
		:cpu = "x86_64" ;
		:featureType = "GRID" ;
		:field = "free-surface, scalar, series" ;
		:file = "espresso_his_3450_0005.nc" ;
		:flt_file = "espresso_flt_3450.nc" ;
		:format = "netCDF-4/HDF5 file" ;
		:fpos_file = "/home/om/roms/espresso/Data/espresso_floats_glider.in" ;
		:frc_file_01 = "/home/om/roms/espresso/Data/espresso_tide_c05_20060101.nc" ;
		:frc_file_02 = "../Data/espresso_river.nc" ;
		:frc_file_03 = "../Data/rain_ncepnam_3hourly_MAB_and_GoM.nc" ;
		:frc_file_04 = "../Data/swrad_ncepnam_3hourly_MAB_and_GoM.nc" ;
		:frc_file_05 = "../Data/Tair_ncepnam_3hourly_MAB_and_GoM.nc" ;
		:frc_file_06 = "../Data/Pair_ncepnam_3hourly_MAB_and_GoM.nc" ;
		:frc_file_07 = "../Data/Qair_ncepnam_3hourly_MAB_and_GoM.nc" ;
		:frc_file_08 = "../Data/lwrad_down_ncepnam_3hourly_MAB_and_GoM.nc" ;
		:frc_file_09 = "../Data/Uwind_ncepnam_3hourly_MAB_and_GoM.nc" ;
		:frc_file_10 = "../Data/Vwind_ncepnam_3hourly_MAB_and_GoM.nc" ;
		:grd_file = "/home/om/roms/espresso/Data/espresso_grid_c05.nc" ;
		:header_dir = "/home/julia/ROMS/espresso/RealTime/Compile/fwd" ;
		:header_file = "espresso.h" ;
		:his_base = "espresso_his_3450" ;
		:history = "ROMS/TOMS, Version 3.5, Tuesday - June 16, 2015 -  6:22:08 AM ;\nFMRC Best Dataset" ;
		:ini_file = "/home/julia/ROMS/espresso/RealTime/Storage/run04/espresso_ini_3450.nc" ;
		:location = "Proto fmrc:espresso_2013_da_his_best" ;
		:os = "Linux" ;
		:rst_file = "espresso_rst_3450.nc" ;
		:script_file = "nl_ocean_espresso.in" ;
		:summary = "Operational nowcast/forecast system version 2 for MARACOOS project (http://maracoos.org). Grid configuration and assimilation system the same as for the ~6 km resolution ROMS Espresso Experiment. These simulations are forced by NCEP NAM operational atmospheric forcing, boundary conditions are from HYCOM-NCODA 1/12 deg model; tides are from ADCIRC model. Assimilated data: SST (AVHRR and blended MW+IR), HF Radar, along track SSH from Jason2." ;
		:svn_rev = "exported" ;
		:svn_url = "https://www.myroms.org/svn/src/trunk" ;
		:tiling = "004x002" ;
		:time = "ocean_time" ;
		:title = "ROMS ESPRESSO Real-Time Operational IS4DVAR Forecast System Version 2 (NEW) 2013-present FMRC History (Best)" ;
		:type = "ROMS/TOMS history file" ;
		:Conventions = "CF-1.5" ;
}

The tardis module jewel is the get_nearest_water function. (Which I just realized should be renamed to get_valid_data!) Together with make_tree it is possible to easily find any valid data point near to a certain position.

Note that iris nearest methods only works with 1D-Coords grids, while tardis will work with virtually any kind of grid: Ugrid, ROMS-Curvilinear (2D Coords), etc.

Right now make_tree uses Scipy's KDTree, but I am re-writing it to use rtree instead. In addition the tree object will be folded into get_nearest_water get_valid_data and cached, so we can skip one function call.

In [10]:
from tardis import get_nearest_water, make_tree

tree, lon, lat = make_tree(cube)

obs = dict(lon=-76.67, lat=34.72)

kw = dict(k=10, max_dist=0.04, min_var=0.01)
series, dist, idx = get_nearest_water(cube, tree, obs['lon'], obs['lat'], **kw)

Here is the result.

In [11]:
%matplotlib inline

from matplotlib import style
style.use('ggplot')

import cartopy.crs as ccrs
import matplotlib.pyplot as plt
from cartopy.feature import NaturalEarthFeature
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER

states = NaturalEarthFeature(category='cultural', scale='50m', facecolor='none',
                             name='admin_1_states_provinces_shp')

def make_map(projection=ccrs.PlateCarree()):
    fig, ax = plt.subplots(figsize=(8, 6),
                           subplot_kw=dict(projection=projection))
    ax.coastlines(resolution='50m')
    gl = ax.gridlines(draw_labels=True)
    gl.xlabels_top = gl.ylabels_right = False
    gl.xformatter = LONGITUDE_FORMATTER
    gl.yformatter = LATITUDE_FORMATTER
    return fig, ax
In [12]:
import numpy.ma as ma


# Latest time stamp for contouring.
c = cube[-1, ...]
c.data = ma.masked_invalid(c.data)

lon = c.coord(axis='X').points
lat = c.coord(axis='Y').points

extent = (lon.min(), lon.max(),
          lat.min(), lat.max())
In [13]:
from oceans import cm

from matplotlib.ticker import MaxNLocator
from mpl_toolkits.axes_grid1.inset_locator import mark_inset
from mpl_toolkits.axes_grid1.inset_locator import zoomed_inset_axes

fig, ax = make_map()
ax.set_extent(extent)
cs = ax.pcolormesh(lon, lat, c.data, cmap=cm.rscolmap, alpha=0.5)
ax.plot(obs['lon'], obs['lat'], 'o', color='darkgreen', label='Beaufort, NC', alpha=0.5)
ax.plot(lon[idx], lat[idx], 'o', color='crimson', label='ESPRESSO', alpha=0.5)
ax.set_title(c.attributes.get('title', None))
ax.add_feature(states, edgecolor='gray')
leg = ax.legend(numpoints=1, fancybox=True, framealpha=0, loc="upper left")

axins = zoomed_inset_axes(ax, 20, loc=5,
                          bbox_to_anchor=(1.1, 0.75),
                          bbox_transform=ax.figure.transFigure)
axins.pcolormesh(lon, lat, c.data, cmap=cm.rscolmap, alpha=0.5)
axins.plot(obs['lon'], obs['lat'], 'o', color='darkgreen', label='Beaufort, NC', alpha=0.5)
axins.plot(lon[idx], lat[idx], 'o', color='crimson', label='ESPRESSO', alpha=0.5)
axins.axis([-76.7, -76.64,
             34.67,  34.73])
mark_inset(ax, axins, loc1=2, loc2=4, ec="0.5")
# Remove ticks
_ = axins.axes.get_xaxis().set_ticks([])
_ = axins.axes.get_yaxis().set_ticks([])

That way we can easily find the nearest ROMS grid point to the Beaufort, NC station.

Let's take a look at the time-series at that point.

In [14]:
import iris.quickplot as qplt

fig, ax = plt.subplots(figsize=(9, 2.75))
l, = qplt.plot(series)
ax.grid(True)
In [15]:
HTML(html)
Out[15]:

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

Creative Commons License
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/.

Comments