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-6p5Afb/python-xarray-0.16.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-6p5Afb/python-xarray-0.16.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)
   2098                            else suppress())
   2099                     with ctx:
-> 2100                         self.figure.draw(renderer)
   2101                     bbox_artists = kwargs.pop("bbox_extra_artists", None)
   2102                     bbox_inches = self.figure.get_tightbbox(renderer,

/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)
   1733
   1734             self.patch.draw(renderer)
-> 1735             mimage._draw_list_compositing_images(
   1736                 renderer, self, artists, self.suppressComposite)
   1737

/usr/lib/python3/dist-packages/matplotlib/image.py in _draw_list_compositing_images(renderer, parent, artists, suppress_composite)
    135     if not_composite or not has_images:
    136         for a in artists:
--> 137             a.draw(renderer)
    138     else:
    139         # 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, **kwargs)
    477         self._done_img_factory = True
    478
--> 479         return matplotlib.axes.Axes.draw(self, renderer=renderer, **kwargs)
    480
    481     def _update_title_position(self, renderer):

/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)
   2628             renderer.stop_rasterizing()
   2629
-> 2630         mimage._draw_list_compositing_images(renderer, self, artists)
   2631
   2632         renderer.close_group('axes')

/usr/lib/python3/dist-packages/matplotlib/image.py in _draw_list_compositing_images(renderer, parent, artists, suppress_composite)
    135     if not_composite or not has_images:
    136         for a in artists:
--> 137             a.draw(renderer)
    138     else:
    139         # 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)
    153         except ValueError:
    154             warnings.warn('Unable to determine extent. Defaulting to global.')
--> 155         geoms = self._feature.intersecting_geometries(extent)
    156
    157         # Combine all the keyword args in priority order.

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

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

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

/usr/lib/python3/dist-packages/cartopy/io/shapereader.py in natural_earth(resolution, category, name)
    293     format_dict = {'config': config, 'category': category,
    294                    'name': name, 'resolution': resolution}
--> 295     return ne_downloader.path(format_dict)
    296
    297

/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)
    344         target_dir = os.path.dirname(target_path)
    345         if not os.path.isdir(target_dir):
--> 346             os.makedirs(target_dir)
    347
    348         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