.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/plot_optical_flow.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note Click :ref:`here ` to download the full example code .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_examples_plot_optical_flow.py: Optical flow ============ This tutorial offers a short overview of the optical flow routines available in pysteps and it will cover how to compute and plot the motion field from a sequence of radar images. .. GENERATED FROM PYTHON SOURCE LINES 9-19 .. code-block:: default from datetime import datetime from pprint import pprint import matplotlib.pyplot as plt import numpy as np from pysteps import io, motion, rcparams from pysteps.utils import conversion, transformation from pysteps.visualization import plot_precip_field, quiver .. GENERATED FROM PYTHON SOURCE LINES 20-26 Read the radar input images --------------------------- First, we will import the sequence of radar composites. You need the pysteps-data archive downloaded and the pystepsrc file configured with the data_source paths pointing to data folders. .. GENERATED FROM PYTHON SOURCE LINES 26-31 .. code-block:: default # Selected case date = datetime.strptime("201505151630", "%Y%m%d%H%M") data_source = rcparams.data_sources["mch"] .. GENERATED FROM PYTHON SOURCE LINES 32-34 Load the data from the archive ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 34-54 .. code-block:: default root_path = data_source["root_path"] path_fmt = data_source["path_fmt"] fn_pattern = data_source["fn_pattern"] fn_ext = data_source["fn_ext"] importer_name = data_source["importer"] importer_kwargs = data_source["importer_kwargs"] timestep = data_source["timestep"] # Find the input files from the archive fns = io.archive.find_by_date( date, root_path, path_fmt, fn_pattern, fn_ext, timestep=5, num_prev_files=9 ) # Read the radar composites importer = io.get_method(importer_name, "importer") R, quality, metadata = io.read_timeseries(fns, importer, **importer_kwargs) del quality # Not used .. GENERATED FROM PYTHON SOURCE LINES 55-57 Preprocess the data ~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 57-70 .. code-block:: default # Convert to mm/h R, metadata = conversion.to_rainrate(R, metadata) # Store the reference frame R_ = R[-1, :, :].copy() # Log-transform the data [dBR] R, metadata = transformation.dB_transform(R, metadata, threshold=0.1, zerovalue=-15.0) # Nicely print the metadata pprint(metadata) .. rst-class:: sphx-glr-script-out Out: .. code-block:: none {'accutime': 5, 'cartesian_unit': 'm', 'institution': 'MeteoSwiss', 'product': 'AQC', 'projection': '+proj=somerc +lon_0=7.43958333333333 +lat_0=46.9524055555556 ' '+k_0=1 +x_0=600000 +y_0=200000 +ellps=bessel ' '+towgs84=674.374,15.056,405.346,0,0,0,0 +units=m +no_defs', 'threshold': -10.0, 'timestamps': array([datetime.datetime(2015, 5, 15, 15, 45), datetime.datetime(2015, 5, 15, 15, 50), datetime.datetime(2015, 5, 15, 15, 55), datetime.datetime(2015, 5, 15, 16, 0), datetime.datetime(2015, 5, 15, 16, 5), datetime.datetime(2015, 5, 15, 16, 10), datetime.datetime(2015, 5, 15, 16, 15), datetime.datetime(2015, 5, 15, 16, 20), datetime.datetime(2015, 5, 15, 16, 25), datetime.datetime(2015, 5, 15, 16, 30)], dtype=object), 'transform': 'dB', 'unit': 'mm/h', 'x1': 255000.0, 'x2': 965000.0, 'xpixelsize': 1000.0, 'y1': -160000.0, 'y2': 480000.0, 'yorigin': 'upper', 'ypixelsize': 1000.0, 'zerovalue': -15.0, 'zr_a': 316.0, 'zr_b': 1.5} .. GENERATED FROM PYTHON SOURCE LINES 71-79 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. .. GENERATED FROM PYTHON SOURCE LINES 79-88 .. code-block:: default oflow_method = motion.get_method("LK") V1 = oflow_method(R[-3:, :, :]) # Plot the motion field on top of the reference frame plot_precip_field(R_, geodata=metadata, title="LK") quiver(V1, geodata=metadata, step=25) plt.show() .. image-sg:: /auto_examples/images/sphx_glr_plot_optical_flow_001.png :alt: LK :srcset: /auto_examples/images/sphx_glr_plot_optical_flow_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 89-98 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. .. GENERATED FROM PYTHON SOURCE LINES 98-107 .. code-block:: default 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) plt.show() .. image-sg:: /auto_examples/images/sphx_glr_plot_optical_flow_002.png :alt: VET :srcset: /auto_examples/images/sphx_glr_plot_optical_flow_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out Out: .. code-block:: none Running VET algorithm original image shape: (3, 640, 710) padded image shape: (3, 640, 710) padded template_image image shape: (3, 640, 710) Number of sectors: 2,2 Sector Shape: (320, 355) Minimizing residuals 3102581.058171551 smoothness_penalty 0.0 original image shape: (3, 640, 710) padded image shape: (3, 640, 712) padded template_image image shape: (3, 640, 712) Number of sectors: 4,4 Sector Shape: (160, 178) Minimizing residuals 2507515.978176187 smoothness_penalty 0.4972321040151195 original image shape: (3, 640, 710) padded image shape: (3, 640, 720) padded template_image image shape: (3, 640, 720) Number of sectors: 16,16 Sector Shape: (40, 45) Minimizing residuals 2261599.3328208188 smoothness_penalty 40.0274102790147 original image shape: (3, 640, 710) padded image shape: (3, 640, 736) padded template_image image shape: (3, 640, 736) Number of sectors: 32,32 Sector Shape: (20, 23) Minimizing residuals 2283873.2019218802 smoothness_penalty 186.8423811170226 .. GENERATED FROM PYTHON SOURCE LINES 108-117 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. .. GENERATED FROM PYTHON SOURCE LINES 117-127 .. code-block:: default 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) plt.show() .. image-sg:: /auto_examples/images/sphx_glr_plot_optical_flow_003.png :alt: DARTS :srcset: /auto_examples/images/sphx_glr_plot_optical_flow_003.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out Out: .. code-block:: none Computing the motion field with the DARTS method. ----- DARTS ----- Computing the FFT of the reflectivity fields...Done in 0.29 seconds. Constructing the y-vector...Done in 0.32 seconds. Constructing the H-matrix...Done in 1.40 seconds. Solving the linear systems...Done in 0.16 seconds. --- 2.2079529762268066 seconds --- .. GENERATED FROM PYTHON SOURCE LINES 128-134 Anisotropic diffusion method (Proesmans et al 1994) --------------------------------------------------- This module implements the anisotropic diffusion method presented in Proesmans et al. (1994), a robust optical flow technique which employs the notion of inconsitency during the solution of the optical flow equations. .. GENERATED FROM PYTHON SOURCE LINES 134-145 .. code-block:: default oflow_method = motion.get_method("proesmans") R[~np.isfinite(R)] = metadata["zerovalue"] V4 = oflow_method(R[-2:, :, :]) # Plot the motion field plot_precip_field(R_, geodata=metadata, title="Proesmans") quiver(V4, geodata=metadata, step=25) plt.show() # sphinx_gallery_thumbnail_number = 1 .. image-sg:: /auto_examples/images/sphx_glr_plot_optical_flow_004.png :alt: Proesmans :srcset: /auto_examples/images/sphx_glr_plot_optical_flow_004.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** ( 1 minutes 6.525 seconds) .. _sphx_glr_download_auto_examples_plot_optical_flow.py: .. only :: html .. container:: sphx-glr-footer :class: sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_optical_flow.py ` .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_optical_flow.ipynb ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_