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

GRIB Data Example

GRIB format is commonly used to disseminate 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')
---------------------------------------------------------------------------
ModuleNotFoundError                       Traceback (most recent call last)
/build/python-xarray-aMeoUt/python-xarray-0.19.0/xarray/tutorial.py in open_dataset(name, cache, cache_dir, engine, **kws)
    113     try:
--> 114         import pooch
    115     except ImportError as e:

ModuleNotFoundError: No module named 'pooch'

The above exception was the direct cause of the following exception:

ImportError                               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-aMeoUt/python-xarray-0.19.0/xarray/tutorial.py in load_dataset(*args, **kwargs)
    223     open_dataset
    224     """
--> 225     with open_dataset(*args, **kwargs) as ds:
    226         return ds.load()
    227

/build/python-xarray-aMeoUt/python-xarray-0.19.0/xarray/tutorial.py in open_dataset(name, cache, cache_dir, engine, **kws)
    114         import pooch
    115     except ImportError as e:
--> 116         raise ImportError(
    117             "tutorial.open_dataset depends on pooch to download and manage datasets."
    118             " To proceed please install pooch."

ImportError: tutorial.open_dataset depends on pooch to download and manage datasets. To proceed please install pooch.

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, pad_inches, bbox_extra_artists, backend, **kwargs)
   2191                            else suppress())
   2192                     with ctx:
-> 2193                         self.figure.draw(renderer)
   2194
   2195                     bbox_inches = self.figure.get_tightbbox(

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

/usr/lib/python3/dist-packages/matplotlib/figure.py in draw(self, renderer)
   1861
   1862             self.patch.draw(renderer)
-> 1863             mimage._draw_list_compositing_images(
   1864                 renderer, self, artists, self.suppressComposite)
   1865

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

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

/usr/lib/python3/dist-packages/cartopy/mpl/geoaxes.py in draw(self, renderer, **kwargs)
    515         self._done_img_factory = True
    516
--> 517         return matplotlib.axes.Axes.draw(self, renderer=renderer, **kwargs)
    518
    519     def _update_title_position(self, renderer):

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

/usr/lib/python3/dist-packages/matplotlib/cbook/deprecation.py in wrapper(*inner_args, **inner_kwargs)
    409                          else deprecation_addendum,
    410                 **kwargs)
--> 411         return func(*inner_args, **inner_kwargs)
    412
    413     return wrapper

/usr/lib/python3/dist-packages/matplotlib/axes/_base.py in draw(self, renderer, inframe)
   2745             renderer.stop_rasterizing()
   2746
-> 2747         mimage._draw_list_compositing_images(renderer, self, artists)
   2748
   2749         renderer.close_group('axes')

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

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

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

/usr/lib/python3/dist-packages/cartopy/feature/__init__.py in intersecting_geometries(self, extent)
    295         """
    296         self.scaler.scale_from_extent(extent)
--> 297         return super().intersecting_geometries(extent)
    298
    299     def with_scale(self, new_scale):

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

/usr/lib/python3/dist-packages/cartopy/feature/__init__.py in geometries(self)
    277         key = (self.name, self.category, self.scale)
    278         if key not in _NATURAL_EARTH_GEOM_CACHE:
--> 279             path = shapereader.natural_earth(resolution=self.scale,
    280                                              category=self.category,
    281                                              name=self.name)

/usr/lib/python3/dist-packages/cartopy/io/shapereader.py in natural_earth(resolution, category, name)
    280     format_dict = {'config': config, 'category': category,
    281                    'name': name, 'resolution': resolution}
--> 282     return ne_downloader.path(format_dict)
    283
    284

/usr/lib/python3/dist-packages/cartopy/io/__init__.py in path(self, format_dict)
    201         else:
    202             # we need to download the file
--> 203             result_path = self.acquire_resource(target_path, format_dict)
    204
    205         return result_path

/usr/lib/python3/dist-packages/cartopy/io/shapereader.py in acquire_resource(self, target_path, format_dict)
    331         target_dir = os.path.dirname(target_path)
    332         if not os.path.isdir(target_dir):
--> 333             os.makedirs(target_dir)
    334
    335         url = self.url(format_dict)

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

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

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

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

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

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

/usr/lib/python3.9/os.py in makedirs(name, mode, exist_ok)
    223             return
    224     try:
--> 225         mkdir(name, mode)
    226     except OSError:
    227         # 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