xarray Rainbow5 backend¶

In this example, we read Rainbow5 data files using the wradlib rainbow xarray backend.

[1]:

import glob
import gzip
import io
import wradlib as wrl
import warnings

warnings.filterwarnings("ignore")
import matplotlib.pyplot as pl
import numpy as np
import xarray as xr

try:
    get_ipython().run_line_magic("matplotlib inline")
except:
    pl.ion()

/home/runner/micromamba-root/envs/wradlib-notebooks/lib/python3.11/site-packages/tqdm/auto.py:22: TqdmWarning: IProgress not found. Please update jupyter and ipywidgets. See https://ipywidgets.readthedocs.io/en/stable/user_install.html
  from .autonotebook import tqdm as notebook_tqdm

Load Rainbow5 Volume Data¶

[2]:

fpath = "rainbow/2013051000000600dBZ.vol"
f = wrl.util.get_wradlib_data_file(fpath)
vol = wrl.io.open_rainbow_dataset(f, reindex_angle=False)

Downloading file 'rainbow/2013051000000600dBZ.vol' from 'https://github.com/wradlib/wradlib-data/raw/pooch/data/rainbow/2013051000000600dBZ.vol' to '/home/runner/work/wradlib-notebooks/wradlib-notebooks/wradlib-data'.

Inspect RadarVolume¶

[3]:

display(vol)

<wradlib.RadarVolume>
Dimension(s): (sweep: 14)
Elevation(s): (0.6, 1.4, 2.4, 3.5, 4.8, 6.3, 8.0, 9.9, 12.2, 14.8, 17.9, 21.3, 25.4, 30.0)

Goereferencing¶

[6]:

swp = vol[0].copy().pipe(wrl.georef.georeference_dataset)

Plotting¶

[7]:

swp.DBZH.plot.pcolormesh(x="x", y="y")
pl.gca().set_aspect("equal")

../../_images/notebooks_fileio_wradlib_rainbow_backend_14_0.png

[8]:

fig = pl.figure(figsize=(10, 10))
swp.DBZH.wradlib.plot_ppi(proj="cg", fig=fig)

[8]:

<matplotlib.collections.QuadMesh at 0x7f05f2c52f50>

../../_images/notebooks_fileio_wradlib_rainbow_backend_15_1.png

[9]:

import cartopy
import cartopy.crs as ccrs
import cartopy.feature as cfeature

map_trans = ccrs.AzimuthalEquidistant(
    central_latitude=swp.latitude.values, central_longitude=swp.longitude.values
)

[10]:

map_proj = ccrs.AzimuthalEquidistant(
    central_latitude=swp.latitude.values, central_longitude=swp.longitude.values
)
pm = swp.DBZH.wradlib.plot_ppi(proj=map_proj)
ax = pl.gca()
ax.gridlines(crs=map_proj)
print(ax)

< GeoAxes: +proj=aeqd +ellps=WGS84 +lon_0=6.379967 +lat_0=50.856633 +x_0=0.0 +y_0=0.0 +no_defs +type=crs >

../../_images/notebooks_fileio_wradlib_rainbow_backend_17_1.png

[11]:

map_proj = ccrs.Mercator(central_longitude=swp.longitude.values)
fig = pl.figure(figsize=(10, 8))
ax = fig.add_subplot(111, projection=map_proj)
pm = swp.DBZH.wradlib.plot_ppi(ax=ax)
ax.gridlines(draw_labels=True)

[11]:

<cartopy.mpl.gridliner.Gridliner at 0x7f05f2bc7290>

../../_images/notebooks_fileio_wradlib_rainbow_backend_18_1.png

[12]:

import cartopy.feature as cfeature


def plot_borders(ax):
    borders = cfeature.NaturalEarthFeature(
        category="physical", name="coastline", scale="10m", facecolor="none"
    )
    ax.add_feature(borders, edgecolor="black", lw=2, zorder=4)


map_proj = ccrs.Mercator(central_longitude=swp.longitude.values)
fig = pl.figure(figsize=(10, 8))
ax = fig.add_subplot(111, projection=map_proj)

DBZH = swp.DBZH
pm = DBZH.where(DBZH > 0).wradlib.plot_ppi(ax=ax)
plot_borders(ax)
ax.gridlines(draw_labels=True)

[12]:

<cartopy.mpl.gridliner.Gridliner at 0x7f05f2be48d0>

../../_images/notebooks_fileio_wradlib_rainbow_backend_19_1.png

[13]:

import matplotlib.path as mpath

theta = np.linspace(0, 2 * np.pi, 100)
center, radius = [0.5, 0.5], 0.5
verts = np.vstack([np.sin(theta), np.cos(theta)]).T
circle = mpath.Path(verts * radius + center)

map_proj = ccrs.AzimuthalEquidistant(
    central_latitude=swp.latitude.values,
    central_longitude=swp.longitude.values,
)
fig = pl.figure(figsize=(10, 8))
ax = fig.add_subplot(111, projection=map_proj)
ax.set_boundary(circle, transform=ax.transAxes)

pm = swp.DBZH.wradlib.plot_ppi(proj=map_proj, ax=ax)
ax = pl.gca()
ax.gridlines(crs=map_proj)

[13]:

<cartopy.mpl.gridliner.Gridliner at 0x7f05ea9e7dd0>

../../_images/notebooks_fileio_wradlib_rainbow_backend_20_1.png

[14]:

fig = pl.figure(figsize=(10, 8))
proj = ccrs.AzimuthalEquidistant(
    central_latitude=swp.latitude.values, central_longitude=swp.longitude.values
)
ax = fig.add_subplot(111, projection=proj)
pm = swp.DBZH.wradlib.plot_ppi(ax=ax)
ax.gridlines()

[14]:

<cartopy.mpl.gridliner.Gridliner at 0x7f05ea94fb10>

../../_images/notebooks_fileio_wradlib_rainbow_backend_21_1.png

[15]:

swp.DBZH.wradlib.plot_ppi()

[15]:

<matplotlib.collections.QuadMesh at 0x7f05e867fc50>

../../_images/notebooks_fileio_wradlib_rainbow_backend_22_1.png

Create simple plot¶

Using xarray features a simple plot can be created like this. Note the sortby('rtime') method, which sorts the radials by time.

[17]:

swp.DBZH.sortby("rtime").plot(x="range", y="rtime", add_labels=False)

[17]:

<matplotlib.collections.QuadMesh at 0x7f05e8669190>

../../_images/notebooks_fileio_wradlib_rainbow_backend_26_1.png

[18]:

fig = pl.figure(figsize=(5, 5))
pm = swp.DBZH.wradlib.plot_ppi(proj={"latmin": 3e3}, fig=fig)

../../_images/notebooks_fileio_wradlib_rainbow_backend_27_0.png

Export to ODIM and CfRadial2¶

[21]:

vol[0].DBZH.sortby("rtime").plot(y="rtime")

[21]:

<matplotlib.collections.QuadMesh at 0x7f05e852af50>

../../_images/notebooks_fileio_wradlib_rainbow_backend_32_1.png

[22]:

vol.to_odim("rainbow_as_odim.h5")
vol.to_cfradial2("rainbow_as_cfradial2.nc")

Check equality¶

We have to sort accordingly and drop the time variables when checking equality.

[25]:

xr.testing.assert_allclose(vol.root, vola.root)
xr.testing.assert_allclose(
    vol[0].drop(["rtime", "time"]), vola[0].drop(["rtime", "time"])
)
xr.testing.assert_allclose(vol.root, volb.root)
xr.testing.assert_allclose(
    vol[0].drop(["rtime", "time"]), volb[0].drop(["rtime", "time"])
)
xr.testing.assert_allclose(vola.root, volb.root)
xr.testing.assert_allclose(
    vola[0].drop(["rtime", "time"]), volb[0].drop(["rtime", "time"])
)

xarray Rainbow5 backend¶

Load Rainbow5 Volume Data¶

Inspect RadarVolume¶

Inspect root group¶

Inspect sweep group(s)¶

Goereferencing¶

Plotting¶

Inspect radar moments¶

Create simple plot¶

Mask some values¶

Export to ODIM and CfRadial2¶

Import again¶

Check equality¶

More Rainbow5 loading mechanisms¶

Use `xr.open_dataset` to retrieve explicit group¶

use wradlib compatibility layer¶

use xradar-backend¶

xarray Rainbow5 backend¶

Load Rainbow5 Volume Data¶

Inspect RadarVolume¶

Inspect root group¶

Inspect sweep group(s)¶

Goereferencing¶

Plotting¶

Inspect radar moments¶

Create simple plot¶

Mask some values¶

Export to ODIM and CfRadial2¶

Import again¶

Check equality¶

More Rainbow5 loading mechanisms¶

Use xr.open_dataset to retrieve explicit group¶

use wradlib compatibility layer¶

use xradar-backend¶

Use `xr.open_dataset` to retrieve explicit group¶