xarray CfRadial2 backend¶

In this example, we read CfRadial2 data files using the xarray cfradial2 backend.

[1]:

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 CfRadial2 Volume Data¶

[2]:

fpath = "netcdf/cfrad.20080604_002217_000_SPOL_v36_SUR_cfradial2.nc"
f = wrl.util.get_wradlib_data_file(fpath)
vol = wrl.io.open_cfradial2_dataset(f)

Downloading file 'netcdf/cfrad.20080604_002217_000_SPOL_v36_SUR_cfradial2.nc' from 'https://github.com/wradlib/wradlib-data/raw/pooch/data/netcdf/cfrad.20080604_002217_000_SPOL_v36_SUR_cfradial2.nc' to '/home/runner/work/wradlib-notebooks/wradlib-notebooks/wradlib-data'.

Inspect RadarVolume¶

[3]:

display(vol)

<wradlib.RadarVolume>
Dimension(s): (sweep: 9)
Elevation(s): (0.5, 1.1, 1.8, 2.6, 3.6, 4.7, 6.5, 9.1, 12.8)

Goereferencing¶

[6]:

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

Plotting¶

[7]:

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

../../_images/notebooks_fileio_wradlib_cfradial2_backend_14_0.png

[8]:

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

[8]:

<matplotlib.collections.QuadMesh at 0x7f588b0fdc10>

../../_images/notebooks_fileio_wradlib_cfradial2_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.DBZ.wradlib.plot_ppi(proj=map_proj)
ax = pl.gca()
ax.gridlines(crs=map_proj)
print(ax)

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

../../_images/notebooks_fileio_wradlib_cfradial2_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.DBZ.wradlib.plot_ppi(ax=ax)
ax.gridlines(draw_labels=True)

[11]:

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

../../_images/notebooks_fileio_wradlib_cfradial2_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)

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

[12]:

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

../../_images/notebooks_fileio_wradlib_cfradial2_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.DBZ.wradlib.plot_ppi(proj=map_proj, ax=ax)
ax = pl.gca()
ax.gridlines(crs=map_proj)

[13]:

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

../../_images/notebooks_fileio_wradlib_cfradial2_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.DBZ.wradlib.plot_ppi(ax=ax)
ax.gridlines()

[14]:

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

../../_images/notebooks_fileio_wradlib_cfradial2_backend_21_1.png

[15]:

swp.DBZ.wradlib.plot_ppi()

[15]:

<matplotlib.collections.QuadMesh at 0x7f587bdc2ed0>

../../_images/notebooks_fileio_wradlib_cfradial2_backend_22_1.png

Create simple plot¶

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

[17]:

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

[17]:

<matplotlib.collections.QuadMesh at 0x7f587bc1ecd0>

../../_images/notebooks_fileio_wradlib_cfradial2_backend_26_1.png

[18]:

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

../../_images/notebooks_fileio_wradlib_cfradial2_backend_27_0.png

Mask some values¶

[19]:

swp["DBZ"] = swp["DBZ"].where(swp["DBZ"] >= 0)
swp["DBZ"].plot()

[19]:

<matplotlib.collections.QuadMesh at 0x7f587bfa5850>

../../_images/notebooks_fileio_wradlib_cfradial2_backend_29_1.png

Export to ODIM and CfRadial2¶

[20]:

vol.to_odim("cfradial2_as_odim.h5")
vol.to_cfradial2("cfradial2_as_cfradial2.nc")

Import again¶

[21]:

vola = wrl.io.open_odim_dataset(
    "cfradial2_as_odim.h5",
    decode_coords=True,
    backend_kwargs=dict(keep_azimuth=True, keep_elevation=True, reindex_angle=False),
)

[22]:

volb = wrl.io.open_cfradial2_dataset("cfradial2_as_cfradial2.nc")

Check equality¶

Some variables need to be dropped, since they are not exported to the other standards or differ slightly (eg. re-indexed ray times).

[23]:

drop = set(vol[0]) ^ set(vola[0]) | set({"rtime"})
xr.testing.assert_allclose(vol.root, vola.root)
xr.testing.assert_allclose(
    vol[0].drop_vars(drop), vola[0].drop_vars(drop, errors="ignore")
)
xr.testing.assert_allclose(vol.root, volb.root)
xr.testing.assert_equal(vol[0], volb[0])
xr.testing.assert_allclose(vola.root, volb.root)
xr.testing.assert_allclose(
    vola[0].drop_vars(drop, errors="ignore"), volb[0].drop_vars(drop, errors="ignore")
)

Read the Docs v: 1.19.0

Versions: 1.19.0; 1.18.0; 1.17.0; 1.16.1; 1.16.0; 1.15.0; 1.14.0; 1.13.0; 1.11.0; 1.10.0; 1.9.0; 1.8.0; 1.7.0; 1.6.0; 1.5.0; 1.4.0; 1.3.0; 1.2.0; 1.1.0; 1.0.0; main

Downloads

On Read the Docs: Project Home; Builds