You can run this notebook in a live session 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-kTVyut/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-kTVyut/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