name = '2015-09-07-iris_outline'
title = 'Iris "hidden" feature'
%matplotlib inline
import os
from datetime import datetime
from IPython.core.display import HTML
with open('creative_commons.txt', 'r') as f:
html = f.read()
hour = datetime.utcnow().strftime('%H:%M')
comments="true"
date = '-'.join(name.split('-')[:3])
slug = '-'.join(name.split('-')[3:])
metadata = dict(title=title,
date=date,
hour=hour,
comments=comments,
slug=slug,
name=name)
markdown = """Title: {title}
date: {date} {hour}
comments: {comments}
slug: {slug}
{{% notebook {name}.ipynb cells[2:] %}}
""".format(**metadata)
content = os.path.abspath(os.path.join(os.getcwd(),
os.pardir,
os.pardir,
'{}.md'.format(name)))
with open('{}'.format(content), 'w') as f:
f.writelines(markdown)
html = '''
<small>
<p> This post was written as an IPython notebook.
It is available for <a href='https://ocefpaf.github.com/python4oceanographers/downloads/notebooks/%s.ipynb'>download</a>
or as a static <a href='https://nbviewer.ipython.org/url/ocefpaf.github.com/python4oceanographers/downloads/notebooks/%s.ipynb'>html</a>.</p>
<p></p>
%s''' % (name, name, html)
I had to create a figure showing the data used to create the World Ocean Atlas climatology (WOA). WOA has that information in each file already, so that part was easy. Let's check the metadata:
import iris
url = ('http://data.nodc.noaa.gov/thredds/dodsC/woa/WOA09/NetCDFdata/'
'silicate_annual_5deg.nc')
cube = iris.load_cube(url, 'Number of Observations')
print(cube)
Number of Observations / (1) (time: 1; depth: 33; latitude: 36; longitude: 72) Dimension coordinates: time x - - - depth - x - - latitude - - x - longitude - - - x Attributes: CVS_ID: 1.0 Conventions: CF-1.4 Metadata_Convention: Unidata Dataset Discovery v1.0 acknowledgment: This work was made possible by a grant from the NOAA Climate and Global... cdm_data_type: Grid contributor_name: Ocean Climate Laboratory creator_email: NODC.Services@noaa.gov creator_name: NODC creator_url: http://www.nodc.noaa.gov date_created: 2011-1-11T18:07:10Z date_issued: 2011-1-11T18:07:10Z date_modified: 2011-1-11T18:07:10Z description: The number of observations of Silicate in each one-degree square of the... geospatial_lat_max: 90.0 geospatial_lat_min: -90.0 geospatial_lat_resolution: 5.0 geospatial_lat_units: degrees_east geospatial_lon_max: 360.0 geospatial_lon_min: 0.0 geospatial_lon_resolution: 5.0 geospatial_lon_units: degrees_north id: 7872719a-ac0b-4f86-9f71-477eeb2ed035 institution: National Oceanographic Data Center(NODC), NOAA keywords: OCEANS > OCEAN CHEMISTRY > SILICATE keywords_vocabulary: Global Change Master Directory license: This data is free for anyone to use, reuse and redistribute without any... naming_authority: gov.noaa.nodc processing_level: synthesized product project: World Ocean Atlas - 2009 publisher_email: NODC.Services@noaa.gov publisher_name: National Oceanographic Data Center(NODC) publisher_url: http://www.nodc.noaa.gov references: Garcia, H. E., R. A. Locarnini, T. P. Boyer, and J. I. Antonov, 2010. World... source: Please see the acknowledgment and References standard_name_vocabulary: CF summary: World Ocean Atlas 2009 (WOA09) is a set of objectively analyzed (1 degree... time_coverage_duration: P1Y time_coverage_end: 2008-02-16 time_coverage_resolution: P1Y time_coverage_start: 1933-07-24 title: World Ocean Atlas 09: Silicate - annual
What I did not know is that iris has a grid outline plotting method!
With the outline
method and pcolormesh
it is easy to see the data cells
and its values.
import numpy as np
import numpy.ma as ma
import iris.plot as iplt
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER
surface = cube[0, 0, ...]
surface.data = ma.log10(ma.masked_equal(surface.data, 0))
fig, ax = plt.subplots(figsize=(15, 9),
subplot_kw=dict(projection=ccrs.PlateCarree()))
ax.coastlines()
iplt.outline(surface)
cs = iplt.pcolormesh(surface, cmap=plt.cm.rainbow)
pt = iplt.points(surface, alpha=0.5, marker='.')
gl = ax.gridlines(draw_labels=True)
gl.xlabels_top = gl.ylabels_right = False
gl.xformatter = LONGITUDE_FORMATTER
gl.yformatter = LATITUDE_FORMATTER
cbar = fig.colorbar(cs, orientation='vertical', extend='both',
shrink=0.7)
tick_levels = [0, 1, 2, 3, 4]
cbar.set_ticks(tick_levels)
cbar.set_label('Number of obs [log10]')
Now I have to explain to our students that, due to the lack of vision of our navy and oil companies, so few observations from Brazilians surveys make into the climatology :-(
HTML(html)
This post was written as an IPython notebook. It is available for download or as a static html.