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-tests/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/wradlib/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 0x7f08c48f49d0>

../../_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 0x7f08c494ca90>

../../_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 0x7f08c162f390>

../../_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 0x7f08c4981890>

../../_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 0x7f08c494f650>

../../_images/notebooks_fileio_wradlib_rainbow_backend_21_1.png

[15]:

swp.DBZH.wradlib.plot_ppi()

[15]:

<matplotlib.collections.QuadMesh at 0x7f08c155de90>

../../_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 0x7f08c157c590>

../../_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 0x7f08c11f9bd0>

../../_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

Contents

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

Contents

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#