Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
R, metadata = transformation.dB_transform(R, metadata, threshold=0.1, zerovalue=-15.0)
# Nicely print the metadata
pprint(metadata)
###############################################################################
# Compute the nowcast
# -------------------
#
# The extrapolation nowcast is based on the estimation of the motion field,
# which is here performed using a local tacking approach (Lucas-Kanade).
# The most recent radar rainfall field is then simply advected along this motion
# field in oder to produce an extrapolation forecast.
# Estimate the motion field with Lucas-Kanade
oflow_method = motion.get_method("LK")
V = oflow_method(R[-3:, :, :])
# Extrapolate the last radar observation
extrapolate = nowcasts.get_method("extrapolation")
R[~np.isfinite(R)] = metadata["zerovalue"]
R_f = extrapolate(R[-1, :, :], V, n_leadtimes)
# Back-transform to rain rate
R_f = transformation.dB_transform(R_f, threshold=-10.0, inverse=True)[0]
# Plot the motion field
plot_precip_field(R_, geodata=metadata)
quiver(V, geodata=metadata, step=50)
###############################################################################
# Verify with FSS
# Plot the motion field
plot_precip_field(R_, geodata=metadata, title="LK")
quiver(V1, geodata=metadata, step=25)
###############################################################################
# Variational echo tracking (VET)
# -------------------------------
#
# This module implements the VET algorithm presented
# by Laroche and Zawadzki (1995) and used in the McGill Algorithm for
# Prediction by Lagrangian Extrapolation (MAPLE) described in
# Germann and Zawadzki (2002).
# The approach essentially consists of a global optimization routine that seeks
# at minimizing a cost function between the displaced and the reference image.
oflow_method = motion.get_method("VET")
V2 = oflow_method(R[-3:, :, :])
# Plot the motion field
plot_precip_field(R_, geodata=metadata, title="VET")
quiver(V2, geodata=metadata, step=25)
###############################################################################
# Dynamic and adaptive radar tracking of storms (DARTS)
# -----------------------------------------------------
#
# DARTS uses a spectral approach to optical flow that is based on the discrete
# Fourier transform (DFT) of a temporal sequence of radar fields.
# The level of truncation of the DFT coefficients controls the degree of
# smoothness of the estimated motion field, allowing for an efficient
# motion estimation. DARTS requires a longer sequence of radar fields for
# estimating the motion, here we are going to use all the available 10 fields.
#
# Also note that the radar image includes NaNs in areas of missing data.
# These are used by the optical flow algorithm to define the radar mask.
#
# Sparse Lucas-Kanade
# -------------------
#
# By setting the optional argument 'dense=False' in 'x,y,u,v = LK_optflow(.....)',
# the LK algorithm returns the motion vectors detected by the Lucas-Kanade scheme
# without interpolating them on the grid.
# This allows us to better identify the presence of wrongly detected
# stationary motion in areas where precipitation is leaving the domain (look
# for the red dots within the blue circle in the figure below).
# get Lucas-Kanade optical flow method
LK_optflow = motion.get_method("LK")
# Mask invalid values
R = np.ma.masked_invalid(R)
R.data[R.mask] = np.nan
# Use default settings (i.e., no buffering of the radar mask)
fd_kwargs1 = {"buffer_mask":0}
xy, uv = LK_optflow(R, dense=False, fd_kwargs=fd_kwargs1)
plt.imshow(ref_dbr, cmap=plt.get_cmap("Greys"))
plt.imshow(mask, cmap=colors.ListedColormap(["black"]), alpha=0.5)
plt.quiver(
xy[:, 0],
xy[:, 1],
uv[:, 0],
uv[:, 1],
color="red",
# Plot the motion field
plot_precip_field(R_, geodata=metadata, title="VET")
quiver(V2, geodata=metadata, step=25)
###############################################################################
# Dynamic and adaptive radar tracking of storms (DARTS)
# -----------------------------------------------------
#
# DARTS uses a spectral approach to optical flow that is based on the discrete
# Fourier transform (DFT) of a temporal sequence of radar fields.
# The level of truncation of the DFT coefficients controls the degree of
# smoothness of the estimated motion field, allowing for an efficient
# motion estimation. DARTS requires a longer sequence of radar fields for
# estimating the motion, here we are going to use all the available 10 fields.
oflow_method = motion.get_method("DARTS")
R[~np.isfinite(R)] = metadata["zerovalue"]
V3 = oflow_method(R) # needs longer training sequence
# Plot the motion field
plot_precip_field(R_, geodata=metadata, title="DARTS")
quiver(V3, geodata=metadata, step=25)
R, metadata = transformation.dB_transform(R, metadata, threshold=0.1, zerovalue=-15.0)
# Nicely print the metadata
pprint(metadata)
###############################################################################
# Lucas-Kanade (LK)
# -----------------
#
# The Lucas-Kanade optical flow method implemented in pysteps is a local
# tracking approach that relies on the OpenCV package.
# Local features are tracked in a sequence of two or more radar images. The
# scheme includes a final interpolation step in order to produce a smooth
# field of motion vectors.
oflow_method = motion.get_method("LK")
V1 = oflow_method(R[-3:, :, :])
# Plot the motion field
plot_precip_field(R_, geodata=metadata, title="LK")
quiver(V1, geodata=metadata, step=25)
###############################################################################
# Variational echo tracking (VET)
# -------------------------------
#
# This module implements the VET algorithm presented
# by Laroche and Zawadzki (1995) and used in the McGill Algorithm for
# Prediction by Lagrangian Extrapolation (MAPLE) described in
# Germann and Zawadzki (2002).
# The approach essentially consists of a global optimization routine that seeks
# at minimizing a cost function between the displaced and the reference image.
## threshold the data
R[R<0.1] = 0.0
## copy the original data
R_ = R.copy()
## set NaN equal to zero
R_[~np.isfinite(R_)] = 0.0
## transform to dBR
R_, _ = stp.utils.dB_transform(R_)
# Compute motion field
oflow_method = stp.motion.get_method("lucaskanade")
UV = oflow_method(R_)
# Perform the advection of the radar field
n_lead_times = 12
adv_method = stp.extrapolation.get_method("semilagrangian")
R_fct = adv_method(R_[-1,:,:], UV, n_lead_times, verbose=True)
## transform forecast values back to mm/h
R_fct, _ = stp.utils.dB_transform(R_fct, inverse=True)
# Plot the nowcast...
stp.plt.animate(R, R_fct=R_fct, UV=UV, nloops=5)