You can run this notebook in a live session Binder or view it on Github.

GRIB Data Example

GRIB format is commonly used to disemminate atmospheric model data. With Xarray and the cfgrib engine, GRIB data can easily be analyzed and visualized.

[1]:
import xarray as xr
import matplotlib.pyplot as plt

To read GRIB data, you can use xarray.load_dataset. The only extra code you need is to specify the engine as cfgrib.

[2]:
ds = xr.tutorial.load_dataset('era5-2mt-2019-03-uk.grib', engine='cfgrib')
---------------------------------------------------------------------------
FileNotFoundError                         Traceback (most recent call last)
<ipython-input-2-783584127f97> in <module>
----> 1 ds = xr.tutorial.load_dataset('era5-2mt-2019-03-uk.grib', engine='cfgrib')

/build/python-xarray-tvkIUx/python-xarray-0.15.0/xarray/tutorial.py in load_dataset(*args, **kwargs)
    111     open_dataset
    112     """
--> 113     with open_dataset(*args, **kwargs) as ds:
    114         return ds.load()
    115

/build/python-xarray-tvkIUx/python-xarray-0.15.0/xarray/tutorial.py in open_dataset(name, cache, cache_dir, github_url, branch, **kws)
     76         # May want to add an option to remove it.
     77         if not _os.path.isdir(longdir):
---> 78             _os.mkdir(longdir)
     79
     80         url = "/".join((github_url, "raw", branch, fullname))

FileNotFoundError: [Errno 2] No such file or directory: '/sbuild-nonexistent/.xarray_tutorial_data'

Let’s create a simple plot of 2-m air temperature in degrees Celsius:

[3]:
ds = ds - 273.15
ds.t2m[0].plot(cmap=plt.cm.coolwarm)
---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
<ipython-input-3-1c96aded89da> in <module>
----> 1 ds = ds - 273.15
      2 ds.t2m[0].plot(cmap=plt.cm.coolwarm)

NameError: name 'ds' is not defined

With CartoPy, we can create a more detailed plot, using built-in shapefiles to help provide geographic context:

[4]:
import cartopy.crs as ccrs
import cartopy
fig = plt.figure(figsize=(10,10))
ax = plt.axes(projection=ccrs.Robinson())
ax.coastlines(resolution='10m')
plot = ds.t2m[0].plot(cmap=plt.cm.coolwarm, transform=ccrs.PlateCarree(), cbar_kwargs={'shrink':0.6})
plt.title('ERA5 - 2m temperature British Isles March 2019')
---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
<ipython-input-4-246f01288f90> in <module>
      4 ax = plt.axes(projection=ccrs.Robinson())
      5 ax.coastlines(resolution='10m')
----> 6 plot = ds.t2m[0].plot(cmap=plt.cm.coolwarm, transform=ccrs.PlateCarree(), cbar_kwargs={'shrink':0.6})
      7 plt.title('ERA5 - 2m temperature British Isles March 2019')

NameError: name 'ds' is not defined
---------------------------------------------------------------------------
PermissionError                           Traceback (most recent call last)
/usr/lib/python3/dist-packages/IPython/core/formatters.py in __call__(self, obj)
    339                 pass
    340             else:
--> 341                 return printer(obj)
    342             # Finally look for special method names
    343             method = get_real_method(obj, self.print_method)

/usr/lib/python3/dist-packages/IPython/core/pylabtools.py in <lambda>(fig)
    246
    247     if 'png' in formats:
--> 248         png_formatter.for_type(Figure, lambda fig: print_figure(fig, 'png', **kwargs))
    249     if 'retina' in formats or 'png2x' in formats:
    250         png_formatter.for_type(Figure, lambda fig: retina_figure(fig, **kwargs))

/usr/lib/python3/dist-packages/IPython/core/pylabtools.py in print_figure(fig, fmt, bbox_inches, **kwargs)
    130         FigureCanvasBase(fig)
    131
--> 132     fig.canvas.print_figure(bytes_io, **kw)
    133     data = bytes_io.getvalue()
    134     if fmt == 'svg':

/usr/lib/python3/dist-packages/matplotlib/backend_bases.py in print_figure(self, filename, dpi, facecolor, edgecolor, orientation, format, bbox_inches, **kwargs)
   2054                     # The first save command (to a BytesIO) is just to estimate
   2055                     # the bounding box of the figure.
-> 2056                     result = print_method(
   2057                         io.BytesIO(),
   2058                         dpi=dpi,

/usr/lib/python3/dist-packages/matplotlib/backends/backend_agg.py in print_png(self, filename_or_obj, metadata, pil_kwargs, *args, **kwargs)
    525
    526         else:
--> 527             FigureCanvasAgg.draw(self)
    528             renderer = self.get_renderer()
    529             with cbook._setattr_cm(renderer, dpi=self.figure.dpi), \

/usr/lib/python3/dist-packages/matplotlib/backends/backend_agg.py in draw(self)
    386         self.renderer = self.get_renderer(cleared=True)
    387         with RendererAgg.lock:
--> 388             self.figure.draw(self.renderer)
    389             # A GUI class may be need to update a window using this draw, so
    390             # don't forget to call the superclass.

/usr/lib/python3/dist-packages/matplotlib/artist.py in draw_wrapper(artist, renderer, *args, **kwargs)
     36                 renderer.start_filter()
     37
---> 38             return draw(artist, renderer, *args, **kwargs)
     39         finally:
     40             if artist.get_agg_filter() is not None:

/usr/lib/python3/dist-packages/matplotlib/figure.py in draw(self, renderer)
   1706
   1707             self.patch.draw(renderer)
-> 1708             mimage._draw_list_compositing_images(
   1709                 renderer, self, artists, self.suppressComposite)
   1710

/usr/lib/python3/dist-packages/matplotlib/image.py in _draw_list_compositing_images(renderer, parent, artists, suppress_composite)
    133     if not_composite or not has_images:
    134         for a in artists:
--> 135             a.draw(renderer)
    136     else:
    137         # Composite any adjacent images together

/usr/lib/python3/dist-packages/matplotlib/artist.py in draw_wrapper(artist, renderer, *args, **kwargs)
     36                 renderer.start_filter()
     37
---> 38             return draw(artist, renderer, *args, **kwargs)
     39         finally:
     40             if artist.get_agg_filter() is not None:

/usr/lib/python3/dist-packages/cartopy/mpl/geoaxes.py in draw(self, renderer, inframe)
    385         self._done_img_factory = True
    386
--> 387         return matplotlib.axes.Axes.draw(self, renderer=renderer,
    388                                          inframe=inframe)
    389

/usr/lib/python3/dist-packages/matplotlib/artist.py in draw_wrapper(artist, renderer, *args, **kwargs)
     36                 renderer.start_filter()
     37
---> 38             return draw(artist, renderer, *args, **kwargs)
     39         finally:
     40             if artist.get_agg_filter() is not None:

/usr/lib/python3/dist-packages/matplotlib/axes/_base.py in draw(self, renderer, inframe)
   2645             renderer.stop_rasterizing()
   2646
-> 2647         mimage._draw_list_compositing_images(renderer, self, artists)
   2648
   2649         renderer.close_group('axes')

/usr/lib/python3/dist-packages/matplotlib/image.py in _draw_list_compositing_images(renderer, parent, artists, suppress_composite)
    133     if not_composite or not has_images:
    134         for a in artists:
--> 135             a.draw(renderer)
    136     else:
    137         # Composite any adjacent images together

/usr/lib/python3/dist-packages/matplotlib/artist.py in draw_wrapper(artist, renderer, *args, **kwargs)
     36                 renderer.start_filter()
     37
---> 38             return draw(artist, renderer, *args, **kwargs)
     39         finally:
     40             if artist.get_agg_filter() is not None:

/usr/lib/python3/dist-packages/cartopy/mpl/feature_artist.py in draw(self, renderer, *args, **kwargs)
    162         except ValueError:
    163             warnings.warn('Unable to determine extent. Defaulting to global.')
--> 164         geoms = self._feature.intersecting_geometries(extent)
    165
    166         # Combine all the keyword args in priority order.

/usr/lib/python3/dist-packages/cartopy/feature/__init__.py in intersecting_geometries(self, extent)
    301         """
    302         self.scaler.scale_from_extent(extent)
--> 303         return super(NaturalEarthFeature, self).intersecting_geometries(extent)
    304
    305     def with_scale(self, new_scale):

/usr/lib/python3/dist-packages/cartopy/feature/__init__.py in intersecting_geometries(self, extent)
    118             extent_geom = sgeom.box(extent[0], extent[2],
    119                                     extent[1], extent[3])
--> 120             return (geom for geom in self.geometries() if
    121                     geom is not None and extent_geom.intersects(geom))
    122         else:

/usr/lib/python3/dist-packages/cartopy/feature/__init__.py in geometries(self)
    283         key = (self.name, self.category, self.scale)
    284         if key not in _NATURAL_EARTH_GEOM_CACHE:
--> 285             path = shapereader.natural_earth(resolution=self.scale,
    286                                              category=self.category,
    287                                              name=self.name)

/usr/lib/python3/dist-packages/cartopy/io/shapereader.py in natural_earth(resolution, category, name)
    357     format_dict = {'config': config, 'category': category,
    358                    'name': name, 'resolution': resolution}
--> 359     return ne_downloader.path(format_dict)
    360
    361

/usr/lib/python3/dist-packages/cartopy/io/__init__.py in path(self, format_dict)
    220         else:
    221             # we need to download the file
--> 222             result_path = self.acquire_resource(target_path, format_dict)
    223
    224         return result_path

/usr/lib/python3/dist-packages/cartopy/io/shapereader.py in acquire_resource(self, target_path, format_dict)
    408         target_dir = os.path.dirname(target_path)
    409         if not os.path.isdir(target_dir):
--> 410             os.makedirs(target_dir)
    411
    412         url = self.url(format_dict)

/usr/lib/python3.8/os.py in makedirs(name, mode, exist_ok)
    211     if head and tail and not path.exists(head):
    212         try:
--> 213             makedirs(head, exist_ok=exist_ok)
    214         except FileExistsError:
    215             # Defeats race condition when another thread created the path

/usr/lib/python3.8/os.py in makedirs(name, mode, exist_ok)
    211     if head and tail and not path.exists(head):
    212         try:
--> 213             makedirs(head, exist_ok=exist_ok)
    214         except FileExistsError:
    215             # Defeats race condition when another thread created the path

/usr/lib/python3.8/os.py in makedirs(name, mode, exist_ok)
    211     if head and tail and not path.exists(head):
    212         try:
--> 213             makedirs(head, exist_ok=exist_ok)
    214         except FileExistsError:
    215             # Defeats race condition when another thread created the path

/usr/lib/python3.8/os.py in makedirs(name, mode, exist_ok)
    211     if head and tail and not path.exists(head):
    212         try:
--> 213             makedirs(head, exist_ok=exist_ok)
    214         except FileExistsError:
    215             # Defeats race condition when another thread created the path

/usr/lib/python3.8/os.py in makedirs(name, mode, exist_ok)
    211     if head and tail and not path.exists(head):
    212         try:
--> 213             makedirs(head, exist_ok=exist_ok)
    214         except FileExistsError:
    215             # Defeats race condition when another thread created the path

/usr/lib/python3.8/os.py in makedirs(name, mode, exist_ok)
    211     if head and tail and not path.exists(head):
    212         try:
--> 213             makedirs(head, exist_ok=exist_ok)
    214         except FileExistsError:
    215             # Defeats race condition when another thread created the path

/usr/lib/python3.8/os.py in makedirs(name, mode, exist_ok)
    221             return
    222     try:
--> 223         mkdir(name, mode)
    224     except OSError:
    225         # Cannot rely on checking for EEXIST, since the operating system

PermissionError: [Errno 13] Permission denied: '/sbuild-nonexistent'
<Figure size 720x720 with 1 Axes>

Finally, we can also pull out a time series for a given location easily:

[5]:
ds.t2m.sel(longitude=0,latitude=51.5).plot()
plt.title('ERA5 - London 2m temperature March 2019')
---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
<ipython-input-5-d6577a89d149> in <module>
----> 1 ds.t2m.sel(longitude=0,latitude=51.5).plot()
      2 plt.title('ERA5 - London 2m temperature March 2019')

NameError: name 'ds' is not defined