From 542d74a94254d76a338ebf562bb5add6e0365448 Mon Sep 17 00:00:00 2001
From: Dominique Brunet <dominique.brunet@canada.ca>
Date: Sat, 10 Jun 2023 06:57:15 +0000
Subject: [PATCH 01/32] Initial commit for Farneback optical flow

---
 pysteps/motion/farneback.py | 216 ++++++++++++++++++++++++++++++++++++
 pysteps/motion/interface.py |  23 ++--
 2 files changed, 229 insertions(+), 10 deletions(-)
 create mode 100644 pysteps/motion/farneback.py

diff --git a/pysteps/motion/farneback.py b/pysteps/motion/farneback.py
new file mode 100644
index 000000000..24b80e608
--- /dev/null
+++ b/pysteps/motion/farneback.py
@@ -0,0 +1,216 @@
+# -*- coding: utf-8 -*-
+"""
+pysteps.motion.farneback
+========================
+
+The Farneback dense optical flow module.
+
+This module implements the interface to the local `Farneback`_ routine
+available in OpenCV_.
+
+.. _OpenCV: https://opencv.org/
+
+.. _`Farneback`:\
+    https://docs.opencv.org/3.4/dc/d6b/group__video__track.html#ga5d10ebbd59fe09c5f650289ec0ece5af
+
+.. autosummary::
+    :toctree: ../generated/
+
+    farneback
+"""
+
+import numpy as np
+from numpy.ma.core import MaskedArray
+import time
+
+from pysteps.decorators import check_input_frames
+from pysteps.exceptions import MissingOptionalDependency
+from pysteps.utils.images import morph_opening
+
+try:
+    import cv2
+    cv2_imported = True
+except ImportError:
+    cv2_imported = False
+
+@check_input_frames(2)
+def farneback(
+    input_images,
+    pyr_scale=0.5,
+    levels=3,
+    winsize=15,
+    iterations=3,
+    poly_n=5,
+    poly_sigma=1.1,
+    flags=0,
+    size_opening=3,
+    verbose=False,
+):
+    """Run the Farneback optical flow routine.
+
+    .. _OpenCV: https://opencv.org/
+
+    .. _`Farneback`:\
+        https://docs.opencv.org/3.4/dc/d6b/group__video__track.html#ga5d10ebbd59fe09c5f650289ec0ece5af
+
+    .. _MaskedArray:\
+        https://docs.scipy.org/doc/numpy/reference/maskedarray.baseclass.html#numpy.ma.MaskedArray
+
+    .. _ndarray:\
+    https://docs.scipy.org/doc/numpy/reference/generated/numpy.ndarray.html
+
+    Parameters
+    ----------
+    input_images: ndarray_ or MaskedArray_
+        Array of shape (T, m, n) containing a sequence of *T* two-dimensional
+        input images of shape (m, n). The indexing order in **input_images** is
+        assumed to be (time, latitude, longitude).
+
+        *T* = 2 is the minimum required number of images.
+        With *T* > 2, all the resulting sparse vectors are pooled together for
+        the final interpolation on a regular grid.
+
+        In case of ndarray_, invalid values (Nans or infs) are masked,
+        otherwise the mask of the MaskedArray_ is used. Such mask defines a
+        region where features are not detected for the tracking algorithm.
+
+    pyr_scale : float, optional
+        Parameter specifying the image scale (<1) to build pyramids for each
+        image; pyr_scale=0.5 means a classical pyramid, where each next layer
+        is twice smaller than the previous one.  This and parameter documented 
+        below are taken directly from the original documentation. 
+        (See https://docs.opencv.org).
+        
+    levels : int, optional
+        number of pyramid layers including the initial image; levels=1 means
+        that no extra layers are created and only the original images are used.
+        
+    winsize : int, optional
+        Averaging window size; larger values increase the algorithm robustness
+        to image noise and give more 
+        Small windows (e.g. 10) lead to unrealistic motion.
+        
+    iterations : int, optional
+        Number of iterations the algorithm does at each pyramid level.
+        
+    poly_n : int
+        Size of the pixel neighborhood used to find polynomial expansion in
+        each pixel; larger values mean that the image will be approximated with
+        smoother surfaces, yielding more robust algorithm and more blurred
+        motion field, typically poly_n = 5 or 7.
+        
+    poly_sigma : float
+        Standard deviation of the Gaussian that is used to smooth derivatives
+        used as a basis for the polynomial expansion; for poly_n=5, you can set
+        poly_sigma=1.1, for poly_n=7, a good value would be poly_sigma=1.5.
+        
+    flags : int, optional
+        Operation flags that can be a combination of the following:
+
+        OPTFLOW_USE_INITIAL_FLOW uses the input 'flow' as an initial flow
+        approximation.
+
+        OPTFLOW_FARNEBACK_GAUSSIAN uses the Gaussian winsize x winsize filter
+        instead of a box filter of the same size for optical flow estimation;
+        usually, this option gives z more accurate flow that with a box filter,
+        at the cost of lower speed; normally, winsize for a Gaussian window
+        should be set to larger value to achieve the same level of robustness.
+        
+    size_opening : int, optional
+        Non-OpenCV parameter:
+        The structuring element size for the filtering of isolated pixels [px].
+
+    verbose: bool, optional
+        If set to True, print some information about the program.
+
+    Returns
+    -------
+    out : ndarray_, shape (2,m,n)
+        Return the advection field having shape
+        (2, m, n), where out[0, :, :] contains the x-components of the motion
+        vectors and out[1, :, :] contains the y-components.
+        The velocities are in units of pixels / timestep, where timestep is the
+        time difference between the two input images.
+        Return a zero motion field of shape (2, m, n) when no motion is
+        detected.
+        
+    References
+    ----------
+    Farnebäck, G.: Two-frame motion estimation based on polynomial expansion, 
+    In Image Analysis, pages 363–370. Springer, 2003.
+
+    """
+
+    if len(input_images.shape) != 3:
+        raise ValueError("input_images has %i dimensions, but a "
+                         "three-dimensional array is expected" 
+                         % len(input_images.shape))
+
+    input_images = input_images.copy()
+
+    if verbose:
+        print("Computing the motion field with the Farneback method.")
+        t0 = time.time()
+
+    if not cv2_imported:
+        raise MissingOptionalDependency(
+            "opencv package is required for the Farneback method "
+            "optical flow method but it is not installed")
+
+    nr_fields = input_images.shape[0]
+    domain_size = (input_images.shape[1], input_images.shape[2])
+
+    for n in range(nr_fields - 1):
+        # extract consecutive images
+        prvs_img = input_images[n, :, :].copy()
+        next_img = input_images[n + 1, :, :].copy()
+
+        # Check if a MaskedArray is used. If not, mask the ndarray
+        if not isinstance(prvs_img, MaskedArray):
+            prvs_img = np.ma.masked_invalid(prvs_img)
+        np.ma.set_fill_value(prvs_img, prvs_img.min())
+        
+        if not isinstance(next_img, MaskedArray):
+            next_img = np.ma.masked_invalid(next_img)
+        np.ma.set_fill_value(next_img, next_img.min())
+
+        # scale between 0 and 255
+        im_min = prvs_img.min()
+        im_max = prvs_img.max()
+        if (im_max - im_min) > 1e-8:
+            prvs_img = (prvs_img.filled() - im_min) / (im_max - im_min) * 255
+        else:
+            prvs_img = prvs_img.filled() - im_min
+    
+        im_min = next_img.min()
+        im_max = next_img.max()
+        if (im_max - im_min) > 1e-8:
+            next_img = (next_img.filled() - im_min) / (im_max - im_min) * 255
+        else:
+            next_img = next_img.filled() - im_min
+    
+        # convert to 8-bit
+        prvs_img = np.ndarray.astype(prvs_img, "uint8")
+        next_img = np.ndarray.astype(next_img, "uint8")
+
+        # remove small noise with a morphological operator (opening)
+        if size_opening > 0:
+            prvs_img = morph_opening(prvs_img, prvs_img.min(), size_opening)
+            next_img = morph_opening(next_img, next_img.min(), size_opening)
+
+        flow = cv2.calcOpticalFlowFarneback(prvs_img, next_img, None, 
+                                            pyr_scale, levels, winsize,
+                                            iterations, poly_n, poly_sigma, 
+                                            flags)
+        
+        fa,fb = np.dsplit(flow,2)
+        u = fa.reshape(domain_size)
+        v = fb.reshape(domain_size)
+        
+        UV = np.stack([u,v])
+
+    if verbose:
+        print("--- %s seconds ---" % (time.time() - t0))
+
+    return UV
+
diff --git a/pysteps/motion/interface.py b/pysteps/motion/interface.py
index eac81aa46..518f8a31a 100644
--- a/pysteps/motion/interface.py
+++ b/pysteps/motion/interface.py
@@ -1,4 +1,3 @@
-# -*- coding: utf-8 -*-
 """
 pysteps.motion.interface
 ========================
@@ -30,23 +29,24 @@
 from pysteps.motion.lucaskanade import dense_lucaskanade
 from pysteps.motion.proesmans import proesmans
 from pysteps.motion.vet import vet
+from pysteps_mods.motion.farneback import farneback
 
 _methods = dict()
-_methods["constant"] = constant
-_methods["lk"] = dense_lucaskanade
-_methods["lucaskanade"] = dense_lucaskanade
-_methods["darts"] = DARTS
-_methods["proesmans"] = proesmans
-_methods["vet"] = vet
+_methods['constant'] = constant
+_methods['lk'] = dense_lucaskanade
+_methods['lucaskanade'] = dense_lucaskanade
+_methods['darts'] = DARTS
+_methods['proesmans'] = proesmans
+_methods['vet'] = vet
+_methods['farneback'] = farneback
 _methods[None] = lambda precip, *args, **kw: np.zeros(
     (2, precip.shape[1], precip.shape[2])
 )
 
 
 def get_method(name):
-    """
-    Return a callable function for the optical flow method corresponding to
-    the given name. The available options are:\n
+    """Return a callable function for the optical flow method corresponding to
+    the given name. The available options are:
 
     +--------------------------------------------------------------------------+
     | Python-based implementations                                             |
@@ -72,6 +72,9 @@ def get_method(name):
     |                   | Laroche and Zawadzki (1995) and                      |
     |                   | Germann and Zawadzki (2002)                          |
     +-------------------+------------------------------------------------------+
+    |  farneback        | OpenCV implementation of the Farneback method which  |
+    |                   | is dense (vector at every pixel) by default.         |
+    +-------------------+------------------------------------------------------+
 
     +--------------------------------------------------------------------------+
     | Methods implemented in C (these require separate compilation and linkage)|

From 5eb402ea15c4f0eec2bebac5993346bd627a43e6 Mon Sep 17 00:00:00 2001
From: Dominique Brunet <dbr000@hpcr7-in4.science.gc.ca>
Date: Fri, 20 Mar 2026 22:14:20 +0000
Subject: [PATCH 02/32] Fix import of farneback in interface.py

---
 pysteps/motion/interface.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/pysteps/motion/interface.py b/pysteps/motion/interface.py
index 3cfb5ec19..221f8952c 100644
--- a/pysteps/motion/interface.py
+++ b/pysteps/motion/interface.py
@@ -30,7 +30,7 @@
 from pysteps.motion.lucaskanade import dense_lucaskanade
 from pysteps.motion.proesmans import proesmans
 from pysteps.motion.vet import vet
-from pysteps_mods.motion.farneback import farneback
+from pysteps.motion.farneback import farneback
 
 _methods = dict()
 _methods['constant'] = constant

From 91deaf2ac0aa73195feeba5fa08610c2055f1d17 Mon Sep 17 00:00:00 2001
From: Dominique Brunet <dominique.brunet@ec.gc.ca>
Date: Sat, 21 Mar 2026 03:45:03 +0000
Subject: [PATCH 03/32] Black format style for farneback.py and interface.py

---
 pysteps/motion/farneback.py | 43 +++++++++++++++++++++++--------------
 pysteps/motion/interface.py | 14 ++++++------
 2 files changed, 34 insertions(+), 23 deletions(-)

diff --git a/pysteps/motion/farneback.py b/pysteps/motion/farneback.py
index 24b80e608..536813753 100644
--- a/pysteps/motion/farneback.py
+++ b/pysteps/motion/farneback.py
@@ -29,10 +29,12 @@
 
 try:
     import cv2
+
     cv2_imported = True
 except ImportError:
     cv2_imported = False
 
+
 @check_input_frames(2)
 def farneback(
     input_images,
@@ -142,9 +144,10 @@ def farneback(
     """
 
     if len(input_images.shape) != 3:
-        raise ValueError("input_images has %i dimensions, but a "
-                         "three-dimensional array is expected" 
-                         % len(input_images.shape))
+        raise ValueError(
+            "input_images has %i dimensions, but a "
+            "three-dimensional array is expected" % len(input_images.shape)
+        )
 
     input_images = input_images.copy()
 
@@ -155,7 +158,8 @@ def farneback(
     if not cv2_imported:
         raise MissingOptionalDependency(
             "opencv package is required for the Farneback method "
-            "optical flow method but it is not installed")
+            "optical flow method but it is not installed"
+        )
 
     nr_fields = input_images.shape[0]
     domain_size = (input_images.shape[1], input_images.shape[2])
@@ -169,7 +173,7 @@ def farneback(
         if not isinstance(prvs_img, MaskedArray):
             prvs_img = np.ma.masked_invalid(prvs_img)
         np.ma.set_fill_value(prvs_img, prvs_img.min())
-        
+
         if not isinstance(next_img, MaskedArray):
             next_img = np.ma.masked_invalid(next_img)
         np.ma.set_fill_value(next_img, next_img.min())
@@ -181,14 +185,14 @@ def farneback(
             prvs_img = (prvs_img.filled() - im_min) / (im_max - im_min) * 255
         else:
             prvs_img = prvs_img.filled() - im_min
-    
+
         im_min = next_img.min()
         im_max = next_img.max()
         if (im_max - im_min) > 1e-8:
             next_img = (next_img.filled() - im_min) / (im_max - im_min) * 255
         else:
             next_img = next_img.filled() - im_min
-    
+
         # convert to 8-bit
         prvs_img = np.ndarray.astype(prvs_img, "uint8")
         next_img = np.ndarray.astype(next_img, "uint8")
@@ -198,19 +202,26 @@ def farneback(
             prvs_img = morph_opening(prvs_img, prvs_img.min(), size_opening)
             next_img = morph_opening(next_img, next_img.min(), size_opening)
 
-        flow = cv2.calcOpticalFlowFarneback(prvs_img, next_img, None, 
-                                            pyr_scale, levels, winsize,
-                                            iterations, poly_n, poly_sigma, 
-                                            flags)
-        
-        fa,fb = np.dsplit(flow,2)
+        flow = cv2.calcOpticalFlowFarneback(
+            prvs_img,
+            next_img,
+            None,
+            pyr_scale,
+            levels,
+            winsize,
+            iterations,
+            poly_n,
+            poly_sigma,
+            flags,
+        )
+
+        fa, fb = np.dsplit(flow, 2)
         u = fa.reshape(domain_size)
         v = fb.reshape(domain_size)
-        
-        UV = np.stack([u,v])
+
+        UV = np.stack([u, v])
 
     if verbose:
         print("--- %s seconds ---" % (time.time() - t0))
 
     return UV
-
diff --git a/pysteps/motion/interface.py b/pysteps/motion/interface.py
index 221f8952c..4f6ceb77f 100644
--- a/pysteps/motion/interface.py
+++ b/pysteps/motion/interface.py
@@ -33,13 +33,13 @@
 from pysteps.motion.farneback import farneback
 
 _methods = dict()
-_methods['constant'] = constant
-_methods['lk'] = dense_lucaskanade
-_methods['lucaskanade'] = dense_lucaskanade
-_methods['darts'] = DARTS
-_methods['proesmans'] = proesmans
-_methods['vet'] = vet
-_methods['farneback'] = farneback
+_methods["constant"] = constant
+_methods["lk"] = dense_lucaskanade
+_methods["lucaskanade"] = dense_lucaskanade
+_methods["darts"] = DARTS
+_methods["proesmans"] = proesmans
+_methods["vet"] = vet
+_methods["farneback"] = farneback
 _methods[None] = lambda precip, *args, **kw: np.zeros(
     (2, precip.shape[1], precip.shape[2])
 )

From 6c8d8993df9d4dac2f9c46b939babfdea1b65fc1 Mon Sep 17 00:00:00 2001
From: Dominique Brunet <dominique.brunet@ec.gc.ca>
Date: Sat, 21 Mar 2026 18:57:08 +0000
Subject: [PATCH 04/32] Black formatting of farneback.py

---
 pysteps/motion/farneback.py | 20 ++++++++++++++++++--
 1 file changed, 18 insertions(+), 2 deletions(-)

diff --git a/pysteps/motion/farneback.py b/pysteps/motion/farneback.py
index 536813753..dbf4fc7df 100644
--- a/pysteps/motion/farneback.py
+++ b/pysteps/motion/farneback.py
@@ -46,6 +46,7 @@ def farneback(
     poly_sigma=1.1,
     flags=0,
     size_opening=3,
+    sigma=60.0,
     verbose=False,
 ):
     """Run the Farneback optical flow routine.
@@ -122,6 +123,12 @@ def farneback(
         Non-OpenCV parameter:
         The structuring element size for the filtering of isolated pixels [px].
 
+    sigma : float, optional
+        Non-OpenCV parameter:
+        The smoothing bandwidth of the motion field. The motion field amplitude
+        is adjusted by multiplying by the ratio of average magnitude before and
+        after smoothing to avoid damping of the motion field.
+
     verbose: bool, optional
         If set to True, print some information about the program.
 
@@ -218,8 +225,17 @@ def farneback(
         fa, fb = np.dsplit(flow, 2)
         u = fa.reshape(domain_size)
         v = fb.reshape(domain_size)
-
-        UV = np.stack([u, v])
+        if sigma > 0:
+            uv2 = u * u + v * v  # squared magnitude of motion field
+            us = sndi.gaussian_filter(u, sigma, mode="nearest")
+            vs = sndi.gaussian_filter(v, signa, mode="nearest")
+            uvs2 = us * us + vs * vs  # squared magnitude of smoothed motion field
+            mult = np.sqrt(np.nanmean(uv2) / np.nanmean(uvs2))
+        else:
+            mult = 1.0
+        if verbose:
+            print("mult factor of smoothed motion field=", mult)
+        UV = np.stack([u * mult, v * mult])
 
     if verbose:
         print("--- %s seconds ---" % (time.time() - t0))

From c3bb273fa51e37a0d59260be6e30662391586b3f Mon Sep 17 00:00:00 2001
From: Dominique Brunet <dominique.brunet@ec.gc.ca>
Date: Sat, 21 Mar 2026 19:00:33 +0000
Subject: [PATCH 05/32] Add import of scipy.image in farneback.py

---
 pysteps/motion/farneback.py | 1 +
 1 file changed, 1 insertion(+)

diff --git a/pysteps/motion/farneback.py b/pysteps/motion/farneback.py
index dbf4fc7df..2748f4685 100644
--- a/pysteps/motion/farneback.py
+++ b/pysteps/motion/farneback.py
@@ -21,6 +21,7 @@
 
 import numpy as np
 from numpy.ma.core import MaskedArray
+import scipy.ndimage as sndi
 import time
 
 from pysteps.decorators import check_input_frames

From c3f8e2d98d43773a6f28cea9ba854556810ad777 Mon Sep 17 00:00:00 2001
From: Dominique Brunet <dominique.brunet@ec.gc.ca>
Date: Mon, 23 Mar 2026 13:30:14 +0000
Subject: [PATCH 06/32] Fix bug in farneback.py (typo)

---
 pysteps/motion/farneback.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/pysteps/motion/farneback.py b/pysteps/motion/farneback.py
index 2748f4685..0bc690450 100644
--- a/pysteps/motion/farneback.py
+++ b/pysteps/motion/farneback.py
@@ -229,7 +229,7 @@ def farneback(
         if sigma > 0:
             uv2 = u * u + v * v  # squared magnitude of motion field
             us = sndi.gaussian_filter(u, sigma, mode="nearest")
-            vs = sndi.gaussian_filter(v, signa, mode="nearest")
+            vs = sndi.gaussian_filter(v, sigma, mode="nearest")
             uvs2 = us * us + vs * vs  # squared magnitude of smoothed motion field
             mult = np.sqrt(np.nanmean(uv2) / np.nanmean(uvs2))
         else:

From 65a0719d22c5ef4fcb6648b3fc84217e10719aec Mon Sep 17 00:00:00 2001
From: Dominique Brunet <dominique.brunet@ec.gc.ca>
Date: Mon, 23 Mar 2026 13:31:47 +0000
Subject: [PATCH 07/32] Black formatting for plot_optical_flow.py

---
 examples/plot_optical_flow.py | 19 ++++++++++++++++++-
 1 file changed, 18 insertions(+), 1 deletion(-)

diff --git a/examples/plot_optical_flow.py b/examples/plot_optical_flow.py
index 2ed26d864..a45fbc717 100644
--- a/examples/plot_optical_flow.py
+++ b/examples/plot_optical_flow.py
@@ -130,7 +130,7 @@
 #
 # 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.
+# inconsistency during the solution of the optical flow equations.
 
 oflow_method = motion.get_method("proesmans")
 R[~np.isfinite(R)] = metadata["zerovalue"]
@@ -141,4 +141,21 @@
 quiver(V4, geodata=metadata, step=25)
 plt.show()
 
+################################################################################
+# Farnebäck smoothed method
+# -------------------------
+#
+# This module implements the pyramidal decomposition method for motion estimation
+# of Farnebäck as implemented in OpenCV, with an option for smoothing and
+# renormalization of the motion fields proposed by Driedger.
+
+oflow_method = motion.get_method("farneback")
+R[~np.isfinite(R)] = metadata["zerovalue"]
+V5 = oflow_method(R[-3:, :, :], verbose=True)
+
+# Plot the motion field
+plot_precip_field(R_, geodata=metadata, title="Farneback")
+quiver(V5, geodata=metadata, step=25)
+plt.show()
+
 # sphinx_gallery_thumbnail_number = 1

From 798f88e4709b200df8d1d3dece817b2b10f613e6 Mon Sep 17 00:00:00 2001
From: Dominique Brunet <dominique.brunet@ec.gc.ca>
Date: Mon, 23 Mar 2026 13:35:57 +0000
Subject: [PATCH 08/32] Change number of input images to 2 for Farneback method
 in plot_optical_flow example

---
 examples/plot_optical_flow.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/examples/plot_optical_flow.py b/examples/plot_optical_flow.py
index a45fbc717..2a0474a48 100644
--- a/examples/plot_optical_flow.py
+++ b/examples/plot_optical_flow.py
@@ -151,7 +151,7 @@
 
 oflow_method = motion.get_method("farneback")
 R[~np.isfinite(R)] = metadata["zerovalue"]
-V5 = oflow_method(R[-3:, :, :], verbose=True)
+V5 = oflow_method(R[-2:, :, :], verbose=True)
 
 # Plot the motion field
 plot_precip_field(R_, geodata=metadata, title="Farneback")

From 0c57950f51c473f17795ce5ad6ec74598e8b6304 Mon Sep 17 00:00:00 2001
From: Dominique Brunet <dominique.brunet@ec.gc.ca>
Date: Mon, 23 Mar 2026 15:39:40 +0000
Subject: [PATCH 09/32] Black formatting

---
 pysteps/motion/farneback.py | 39 +++++++++++++++++++++++--------------
 1 file changed, 24 insertions(+), 15 deletions(-)

diff --git a/pysteps/motion/farneback.py b/pysteps/motion/farneback.py
index 0bc690450..8a48421a2 100644
--- a/pysteps/motion/farneback.py
+++ b/pysteps/motion/farneback.py
@@ -170,9 +170,11 @@ def farneback(
         )
 
     nr_fields = input_images.shape[0]
+    nr_pairs = nr_fields - 1
     domain_size = (input_images.shape[1], input_images.shape[2])
-
-    for n in range(nr_fields - 1):
+    u_sum = np.zeros(domain_size)
+    v_sum = np.zeros(domain_size)
+    for n in range(nr_pairs):
         # extract consecutive images
         prvs_img = input_images[n, :, :].copy()
         next_img = input_images[n + 1, :, :].copy()
@@ -224,19 +226,26 @@ def farneback(
         )
 
         fa, fb = np.dsplit(flow, 2)
-        u = fa.reshape(domain_size)
-        v = fb.reshape(domain_size)
-        if sigma > 0:
-            uv2 = u * u + v * v  # squared magnitude of motion field
-            us = sndi.gaussian_filter(u, sigma, mode="nearest")
-            vs = sndi.gaussian_filter(v, sigma, mode="nearest")
-            uvs2 = us * us + vs * vs  # squared magnitude of smoothed motion field
-            mult = np.sqrt(np.nanmean(uv2) / np.nanmean(uvs2))
-        else:
-            mult = 1.0
-        if verbose:
-            print("mult factor of smoothed motion field=", mult)
-        UV = np.stack([u * mult, v * mult])
+        u_sum += fa.reshape(domain_size)
+        v_sum += fb.reshape(domain_size)
+
+    # Compute the average motion field
+    u = u_sum / nr_pairs
+    v = v_sum / nr_pairs
+
+    # Smoothing
+    if sigma > 0:
+        uv2 = u * u + v * v  # squared magnitude of motion field
+        us = sndi.gaussian_filter(u, sigma, mode="nearest")
+        vs = sndi.gaussian_filter(v, sigma, mode="nearest")
+        uvs2 = us * us + vs * vs  # squared magnitude of smoothed motion field
+        mult = np.sqrt(np.nanmean(uv2) / np.nanmean(uvs2))
+    else:
+        mult = 1.0
+    if verbose:
+        print("mult factor of smoothed motion field=", mult)
+
+    UV = np.stack([us * mult, vs * mult])
 
     if verbose:
         print("--- %s seconds ---" % (time.time() - t0))

From 62d677734647468d5de58c3894823b18066bcfbe Mon Sep 17 00:00:00 2001
From: Dominique Brunet <dominique.brunet@ec.gc.ca>
Date: Mon, 23 Mar 2026 16:25:57 +0000
Subject: [PATCH 10/32] Black formatting

---
 pysteps/motion/farneback.py | 8 +++++++-
 1 file changed, 7 insertions(+), 1 deletion(-)

diff --git a/pysteps/motion/farneback.py b/pysteps/motion/farneback.py
index 8a48421a2..d7619ba01 100644
--- a/pysteps/motion/farneback.py
+++ b/pysteps/motion/farneback.py
@@ -239,7 +239,13 @@ def farneback(
         us = sndi.gaussian_filter(u, sigma, mode="nearest")
         vs = sndi.gaussian_filter(v, sigma, mode="nearest")
         uvs2 = us * us + vs * vs  # squared magnitude of smoothed motion field
-        mult = np.sqrt(np.nanmean(uv2) / np.nanmean(uvs2))
+
+        mean_uv2 = np.nanmean(uv2)
+        mean_uvs2 = np.nanmean(uvs2)
+        if mean_uvs2 > 0:
+            mult = np.sqrt(mean_uv2 / mean_uvs2)
+        else:
+            mult = 1.0
     else:
         mult = 1.0
     if verbose:

From ac607eee26c5668973de90785f7c1de796351880 Mon Sep 17 00:00:00 2001
From: Dominique Brunet <dominique.brunet@ec.gc.ca>
Date: Mon, 23 Mar 2026 16:28:09 +0000
Subject: [PATCH 11/32] Black formatting

---
 examples/optical_flow_methods_convergence.py | 18 ++++++++++++++++++
 1 file changed, 18 insertions(+)

diff --git a/examples/optical_flow_methods_convergence.py b/examples/optical_flow_methods_convergence.py
index c5948ac4c..1d513f386 100644
--- a/examples/optical_flow_methods_convergence.py
+++ b/examples/optical_flow_methods_convergence.py
@@ -374,4 +374,22 @@ def plot_optflow_method_convergence(input_precip, optflow_method_name, motion_ty
 # ~~~~~~~~~~~~~~~~~
 plot_optflow_method_convergence(reference_field, "DARTS", "rotor")
 
+################################################################################
+# Farneback
+# -----
+#
+# Constant motion x-direction
+# ~~~~~~~~~~~~~~~~~~~~~~~~~~~
+plot_optflow_method_convergence(reference_field, "farneback", "linear_x")
+
+################################################################################
+# Constant motion y-direction
+# ~~~~~~~~~~~~~~~~~~~~~~~~~~~
+plot_optflow_method_convergence(reference_field, "farneback", "linear_y")
+
+################################################################################
+# Rotational motion
+# ~~~~~~~~~~~~~~~~~
+plot_optflow_method_convergence(reference_field, "farneback", "rotor")
+
 # sphinx_gallery_thumbnail_number = 5

From 1fa70636e0718c54696bf716b867d9ed598ef166 Mon Sep 17 00:00:00 2001
From: dominique-brunet-eccc
 <136009011+dominique-brunet-eccc@users.noreply.github.com>
Date: Mon, 23 Mar 2026 15:21:22 -0400
Subject: [PATCH 12/32] Update farneback.py, small refactor

---
 pysteps/motion/farneback.py | 6 ++----
 1 file changed, 2 insertions(+), 4 deletions(-)

diff --git a/pysteps/motion/farneback.py b/pysteps/motion/farneback.py
index d7619ba01..46917b802 100644
--- a/pysteps/motion/farneback.py
+++ b/pysteps/motion/farneback.py
@@ -71,8 +71,7 @@ def farneback(
         assumed to be (time, latitude, longitude).
 
         *T* = 2 is the minimum required number of images.
-        With *T* > 2, all the resulting sparse vectors are pooled together for
-        the final interpolation on a regular grid.
+        With *T* > 2, all the resulting motion vectors are averaged together.
 
         In case of ndarray_, invalid values (Nans or infs) are masked,
         otherwise the mask of the MaskedArray_ is used. Such mask defines a
@@ -169,8 +168,7 @@ def farneback(
             "optical flow method but it is not installed"
         )
 
-    nr_fields = input_images.shape[0]
-    nr_pairs = nr_fields - 1
+    nr_pairs = input_images.shape[0] - 1
     domain_size = (input_images.shape[1], input_images.shape[2])
     u_sum = np.zeros(domain_size)
     v_sum = np.zeros(domain_size)

From 508c19c7c84e9213826fd300b9cf1111a58e6cbf Mon Sep 17 00:00:00 2001
From: dominique-brunet-eccc
 <136009011+dominique-brunet-eccc@users.noreply.github.com>
Date: Mon, 23 Mar 2026 15:23:16 -0400
Subject: [PATCH 13/32] Update interface.py, shorten description for farneback

---
 pysteps/motion/interface.py | 3 +--
 1 file changed, 1 insertion(+), 2 deletions(-)

diff --git a/pysteps/motion/interface.py b/pysteps/motion/interface.py
index 4f6ceb77f..8400460c5 100644
--- a/pysteps/motion/interface.py
+++ b/pysteps/motion/interface.py
@@ -73,8 +73,7 @@ def get_method(name):
     |                   | Laroche and Zawadzki (1995) and                      |
     |                   | Germann and Zawadzki (2002)                          |
     +-------------------+------------------------------------------------------+
-    |  farneback        | OpenCV implementation of the Farneback method which  |
-    |                   | is dense (vector at every pixel) by default.         |
+    |  farneback        | OpenCV implementation of the Farneback (2003) method.|
     +-------------------+------------------------------------------------------+
 
     +--------------------------------------------------------------------------+

From a0808c84827decbdc83aaaa1ddb0937d3a4798b5 Mon Sep 17 00:00:00 2001
From: dominique-brunet-eccc
 <136009011+dominique-brunet-eccc@users.noreply.github.com>
Date: Mon, 23 Mar 2026 15:24:38 -0400
Subject: [PATCH 14/32] Update plot_optical_flow.py

---
 examples/plot_optical_flow.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/examples/plot_optical_flow.py b/examples/plot_optical_flow.py
index 2a0474a48..2e5d92da5 100644
--- a/examples/plot_optical_flow.py
+++ b/examples/plot_optical_flow.py
@@ -147,7 +147,7 @@
 #
 # This module implements the pyramidal decomposition method for motion estimation
 # of Farnebäck as implemented in OpenCV, with an option for smoothing and
-# renormalization of the motion fields proposed by Driedger.
+# renormalization of the motion fields proposed by N. Driedger.
 
 oflow_method = motion.get_method("farneback")
 R[~np.isfinite(R)] = metadata["zerovalue"]

From 4a0268935f0f43d4b38de319b6200896a74e3287 Mon Sep 17 00:00:00 2001
From: Dominique Brunet <dominique.brunet@ec.gc.ca>
Date: Mon, 23 Mar 2026 19:32:05 +0000
Subject: [PATCH 15/32] Minor edits of interface.py to avoid unnecessary
 changes for merge

---
 pysteps/motion/interface.py | 6 ++++--
 1 file changed, 4 insertions(+), 2 deletions(-)

diff --git a/pysteps/motion/interface.py b/pysteps/motion/interface.py
index 4f6ceb77f..95e42dbad 100644
--- a/pysteps/motion/interface.py
+++ b/pysteps/motion/interface.py
@@ -1,3 +1,4 @@
+# -*- coding: utf-8 -*-
 """
 pysteps.motion.interface
 ========================
@@ -46,8 +47,9 @@
 
 
 def get_method(name):
-    """Return a callable function for the optical flow method corresponding to
-    the given name. The available options are:
+    """
+    Return a callable function for the optical flow method corresponding to
+    the given name. The available options are:\n
 
     +--------------------------------------------------------------------------+
     | Python-based implementations                                             |

From eb24430d5fd6696f489a694f9c8f72904fbbdc21 Mon Sep 17 00:00:00 2001
From: dominique-brunet-eccc
 <136009011+dominique-brunet-eccc@users.noreply.github.com>
Date: Tue, 24 Mar 2026 10:17:32 -0400
Subject: [PATCH 16/32] Update pysteps/motion/farneback.py

Fix typo.

Co-authored-by: Copilot <175728472+Copilot@users.noreply.github.com>
---
 pysteps/motion/farneback.py | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/pysteps/motion/farneback.py b/pysteps/motion/farneback.py
index 46917b802..8e7fb750d 100644
--- a/pysteps/motion/farneback.py
+++ b/pysteps/motion/farneback.py
@@ -115,9 +115,9 @@ def farneback(
 
         OPTFLOW_FARNEBACK_GAUSSIAN uses the Gaussian winsize x winsize filter
         instead of a box filter of the same size for optical flow estimation;
-        usually, this option gives z more accurate flow that with a box filter,
+        usually, this option gives a more accurate flow than with a box filter,
         at the cost of lower speed; normally, winsize for a Gaussian window
-        should be set to larger value to achieve the same level of robustness.
+        should be set to a larger value to achieve the same level of robustness.
         
     size_opening : int, optional
         Non-OpenCV parameter:

From e46d8a0b299c5d84f35d2b09723b43831023c68c Mon Sep 17 00:00:00 2001
From: dominique-brunet-eccc
 <136009011+dominique-brunet-eccc@users.noreply.github.com>
Date: Tue, 24 Mar 2026 10:18:52 -0400
Subject: [PATCH 17/32] Update examples/optical_flow_methods_convergence.py

Fix styling.

Co-authored-by: Copilot <175728472+Copilot@users.noreply.github.com>
---
 examples/optical_flow_methods_convergence.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/examples/optical_flow_methods_convergence.py b/examples/optical_flow_methods_convergence.py
index 1d513f386..9208af594 100644
--- a/examples/optical_flow_methods_convergence.py
+++ b/examples/optical_flow_methods_convergence.py
@@ -376,7 +376,7 @@ def plot_optflow_method_convergence(input_precip, optflow_method_name, motion_ty
 
 ################################################################################
 # Farneback
-# -----
+# ---------
 #
 # Constant motion x-direction
 # ~~~~~~~~~~~~~~~~~~~~~~~~~~~

From cbdaa8dd0ce6d0b1f0cea2d1d60699818d00cf90 Mon Sep 17 00:00:00 2001
From: dominique-brunet-eccc
 <136009011+dominique-brunet-eccc@users.noreply.github.com>
Date: Tue, 24 Mar 2026 10:21:14 -0400
Subject: [PATCH 18/32] Update farneback.py to initialize us, vs when sigma <=
 0

---
 pysteps/motion/farneback.py | 2 ++
 1 file changed, 2 insertions(+)

diff --git a/pysteps/motion/farneback.py b/pysteps/motion/farneback.py
index 8e7fb750d..ada00cda4 100644
--- a/pysteps/motion/farneback.py
+++ b/pysteps/motion/farneback.py
@@ -246,6 +246,8 @@ def farneback(
             mult = 1.0
     else:
         mult = 1.0
+        us = u
+        vs = v
     if verbose:
         print("mult factor of smoothed motion field=", mult)
 

From 380a49ecf7b16b273f4a4563025308943cff974e Mon Sep 17 00:00:00 2001
From: dominique-brunet-eccc
 <136009011+dominique-brunet-eccc@users.noreply.github.com>
Date: Tue, 24 Mar 2026 10:23:45 -0400
Subject: [PATCH 19/32] Update pysteps/motion/farneback.py all caps for
 optional dependency

Co-authored-by: Copilot <175728472+Copilot@users.noreply.github.com>
---
 pysteps/motion/farneback.py | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/pysteps/motion/farneback.py b/pysteps/motion/farneback.py
index ada00cda4..12303714f 100644
--- a/pysteps/motion/farneback.py
+++ b/pysteps/motion/farneback.py
@@ -31,9 +31,9 @@
 try:
     import cv2
 
-    cv2_imported = True
+    CV2_IMPORTED = True
 except ImportError:
-    cv2_imported = False
+    CV2_IMPORTED = False
 
 
 @check_input_frames(2)

From b8569198599565a50783a9c82d3c159cf988d021 Mon Sep 17 00:00:00 2001
From: dominique-brunet-eccc
 <136009011+dominique-brunet-eccc@users.noreply.github.com>
Date: Tue, 24 Mar 2026 10:24:27 -0400
Subject: [PATCH 20/32] Update pysteps/motion/farneback.py fix error message
 wording

Co-authored-by: Copilot <175728472+Copilot@users.noreply.github.com>
---
 pysteps/motion/farneback.py | 3 +--
 1 file changed, 1 insertion(+), 2 deletions(-)

diff --git a/pysteps/motion/farneback.py b/pysteps/motion/farneback.py
index 12303714f..d9693ad2f 100644
--- a/pysteps/motion/farneback.py
+++ b/pysteps/motion/farneback.py
@@ -164,8 +164,7 @@ def farneback(
 
     if not cv2_imported:
         raise MissingOptionalDependency(
-            "opencv package is required for the Farneback method "
-            "optical flow method but it is not installed"
+            "OpenCV (cv2) is required for the Farneback optical flow method, but it is not installed"
         )
 
     nr_pairs = input_images.shape[0] - 1

From 33a0580b59e1d24d59bcf4b0005687c90a09637a Mon Sep 17 00:00:00 2001
From: dominique-brunet-eccc
 <136009011+dominique-brunet-eccc@users.noreply.github.com>
Date: Tue, 24 Mar 2026 10:25:39 -0400
Subject: [PATCH 21/32] Fix pysteps/motion/farneback.py grammar in description
 of inputs parameters.

Co-authored-by: Copilot <175728472+Copilot@users.noreply.github.com>
---
 pysteps/motion/farneback.py | 21 ++++++++++-----------
 1 file changed, 10 insertions(+), 11 deletions(-)

diff --git a/pysteps/motion/farneback.py b/pysteps/motion/farneback.py
index d9693ad2f..4131dc8f3 100644
--- a/pysteps/motion/farneback.py
+++ b/pysteps/motion/farneback.py
@@ -78,21 +78,20 @@ def farneback(
         region where features are not detected for the tracking algorithm.
 
     pyr_scale : float, optional
-        Parameter specifying the image scale (<1) to build pyramids for each
-        image; pyr_scale=0.5 means a classical pyramid, where each next layer
-        is twice smaller than the previous one.  This and parameter documented 
-        below are taken directly from the original documentation. 
-        (See https://docs.opencv.org).
-        
+        Parameter specifying the image scale (<1) used to build pyramids for
+        each image; pyr_scale=0.5 means a classical pyramid, where each next
+        layer is twice smaller than the previous one. This and the following
+        parameter descriptions are adapted from the original OpenCV
+        documentation (see https://docs.opencv.org).
+
     levels : int, optional
-        number of pyramid layers including the initial image; levels=1 means
+        Number of pyramid layers including the initial image; levels=1 means
         that no extra layers are created and only the original images are used.
-        
+
     winsize : int, optional
         Averaging window size; larger values increase the algorithm robustness
-        to image noise and give more 
-        Small windows (e.g. 10) lead to unrealistic motion.
-        
+        to image noise and give more stable motion estimates. Small windows
+        (e.g. 10) lead to unrealistic motion.
     iterations : int, optional
         Number of iterations the algorithm does at each pyramid level.
         

From 881a2af55b74d3ebccc64988edf4e4ce55dff914 Mon Sep 17 00:00:00 2001
From: dominique-brunet-eccc
 <136009011+dominique-brunet-eccc@users.noreply.github.com>
Date: Tue, 24 Mar 2026 10:42:52 -0400
Subject: [PATCH 22/32] Update test_motion.py to include test for masked arrays
 for both Farneback  and LK

---
 pysteps/tests/test_motion.py | 6 ++++--
 1 file changed, 4 insertions(+), 2 deletions(-)

diff --git a/pysteps/tests/test_motion.py b/pysteps/tests/test_motion.py
index 8d198960a..147a8e2eb 100644
--- a/pysteps/tests/test_motion.py
+++ b/pysteps/tests/test_motion.py
@@ -393,16 +393,18 @@ def test_vet_cost_function():
     assert (returned_values[0] - 1548250.87627097) < 0.001
 
 
-def test_lk_masked_array():
+@pytest.mark.parametrize("method", ["LK", "farneback"])
+def test_motion_masked_array(method):
     """
     Passing a ndarray with NaNs or a masked array should produce the same results.
+    Tests for both LK and Farneback motion estimation methods.
     """
     pytest.importorskip("cv2")
 
     __, precip_obs = _create_observations(
         reference_field.copy(), "linear_y", num_times=2
     )
-    motion_method = motion.get_method("LK")
+    motion_method = motion.get_method(method)
 
     # ndarray with nans
     np.ma.set_fill_value(precip_obs, -15)

From 9c42f51363a951c686f97696ab7ae47401c6cc6a Mon Sep 17 00:00:00 2001
From: Dominique Brunet <dominique.brunet@ec.gc.ca>
Date: Tue, 24 Mar 2026 18:29:16 +0000
Subject: [PATCH 23/32] Capitalize CV2_IMPORTED everywhere in
 motion/farneback.py

---
 pysteps/motion/farneback.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/pysteps/motion/farneback.py b/pysteps/motion/farneback.py
index 4131dc8f3..84ac734f5 100644
--- a/pysteps/motion/farneback.py
+++ b/pysteps/motion/farneback.py
@@ -161,7 +161,7 @@ def farneback(
         print("Computing the motion field with the Farneback method.")
         t0 = time.time()
 
-    if not cv2_imported:
+    if not CV2_IMPORTED:
         raise MissingOptionalDependency(
             "OpenCV (cv2) is required for the Farneback optical flow method, but it is not installed"
         )

From 1cfc15bdd3a7674d35249cc78725f1710266cbdf Mon Sep 17 00:00:00 2001
From: Dominique Brunet <dominique.brunet@ec.gc.ca>
Date: Tue, 24 Mar 2026 18:32:38 +0000
Subject: [PATCH 24/32] Fix Black formatting

---
 pysteps/tests/test_motion.py | 14 ++++++++++----
 1 file changed, 10 insertions(+), 4 deletions(-)

diff --git a/pysteps/tests/test_motion.py b/pysteps/tests/test_motion.py
index 147a8e2eb..d46cbef47 100644
--- a/pysteps/tests/test_motion.py
+++ b/pysteps/tests/test_motion.py
@@ -393,8 +393,14 @@ def test_vet_cost_function():
     assert (returned_values[0] - 1548250.87627097) < 0.001
 
 
-@pytest.mark.parametrize("method", ["LK", "farneback"])
-def test_motion_masked_array(method):
+@pytest.mark.parametrize(
+    "method,kwargs",
+    [
+        ("LK", {"fd_kwargs": {"buffer_mask": 20}, "verbose": False}),
+        ("farneback", {"verbose": False}),
+    ],
+)
+def test_motion_masked_array(method, kwargs):
     """
     Passing a ndarray with NaNs or a masked array should produce the same results.
     Tests for both LK and Farneback motion estimation methods.
@@ -410,11 +416,11 @@ def test_motion_masked_array(method):
     np.ma.set_fill_value(precip_obs, -15)
     ndarray = precip_obs.filled()
     ndarray[ndarray == -15] = np.nan
-    uv_ndarray = motion_method(ndarray, fd_kwargs={"buffer_mask": 20}, verbose=False)
+    uv_ndarray = motion_method(ndarray, **kwargs)
 
     # masked array
     mdarray = np.ma.masked_invalid(ndarray)
     mdarray.data[mdarray.mask] = -15
-    uv_mdarray = motion_method(mdarray, fd_kwargs={"buffer_mask": 20}, verbose=False)
+    uv_mdarray = motion_method(mdarray, **kwargs)
 
     assert np.abs(uv_mdarray - uv_ndarray).max() < 0.01

From 37bea59774764d730ffd60687e82b730cc01cc40 Mon Sep 17 00:00:00 2001
From: dominique-brunet-eccc
 <136009011+dominique-brunet-eccc@users.noreply.github.com>
Date: Tue, 24 Mar 2026 14:43:13 -0400
Subject: [PATCH 25/32] Update test_motion.py to test farneback method

---
 pysteps/tests/test_motion.py | 6 +++++-
 1 file changed, 5 insertions(+), 1 deletion(-)

diff --git a/pysteps/tests/test_motion.py b/pysteps/tests/test_motion.py
index d46cbef47..d7bf20231 100644
--- a/pysteps/tests/test_motion.py
+++ b/pysteps/tests/test_motion.py
@@ -164,6 +164,8 @@ def _create_observations(input_precip, motion_type, num_times=9):
     (reference_field, "proesmans", "linear_y", 2, 0.45),
     (reference_field, "darts", "linear_x", 9, 20),
     (reference_field, "darts", "linear_y", 9, 20),
+    (reference_field, "farneback", "linear_x", 2, 50),
+    (reference_field, "farneback", "linear_y", 2, 50),
 ]
 
 
@@ -256,6 +258,7 @@ def test_optflow_method_convergence(
     ("vet", 3),
     ("darts", 9),
     ("proesmans", 2),
+    ("farneback", 2),
 ]
 
 
@@ -296,6 +299,7 @@ def test_no_precipitation(optflow_method_name, num_times):
     ("vet", 2, 3),
     ("darts", 9, 9),
     ("proesmans", 2, 2),
+    ("farneback", 2, np.inf),
 ]
 
 
@@ -303,7 +307,7 @@ def test_no_precipitation(optflow_method_name, num_times):
 def test_input_shape_checks(
     optflow_method_name, minimum_input_frames, maximum_input_frames
 ):
-    if optflow_method_name == "lk":
+    if optflow_method_name in ("lk", "farneback"):
         pytest.importorskip("cv2")
     image_size = 100
     motion_method = motion.get_method(optflow_method_name)

From 67db9a7b39d01f921f7505f1b2dab3d1f854dd70 Mon Sep 17 00:00:00 2001
From: dominique-brunet-eccc
 <136009011+dominique-brunet-eccc@users.noreply.github.com>
Date: Tue, 24 Mar 2026 14:55:44 -0400
Subject: [PATCH 26/32] Update test_motion.py tuning the relative RMSE for
 Farneback

---
 pysteps/tests/test_motion.py | 36 ++++++++++++++++++++++++------------
 1 file changed, 24 insertions(+), 12 deletions(-)

diff --git a/pysteps/tests/test_motion.py b/pysteps/tests/test_motion.py
index d7bf20231..63d4b4179 100644
--- a/pysteps/tests/test_motion.py
+++ b/pysteps/tests/test_motion.py
@@ -152,20 +152,32 @@ def _create_observations(input_precip, motion_type, num_times=9):
 )
 
 convergence_arg_values = [
-    (reference_field, "lk", "linear_x", 2, 0.1),
-    (reference_field, "lk", "linear_y", 2, 0.1),
-    (reference_field, "lk", "linear_x", 3, 0.1),
-    (reference_field, "lk", "linear_y", 3, 0.1),
-    (reference_field, "vet", "linear_x", 2, 0.1),
-    # (reference_field, 'vet', 'linear_x', 3, 9),
-    # (reference_field, 'vet', 'linear_y', 2, 9),
-    (reference_field, "vet", "linear_y", 3, 0.1),
-    (reference_field, "proesmans", "linear_x", 2, 0.45),
-    (reference_field, "proesmans", "linear_y", 2, 0.45),
-    (reference_field, "darts", "linear_x", 9, 20),
-    (reference_field, "darts", "linear_y", 9, 20),
+    #(reference_field, "lk", "linear_x", 2, 0.1),
+    #(reference_field, "lk", "linear_y", 2, 0.1),
+    #(reference_field, "lk", "linear_x", 3, 0.1),
+    #(reference_field, "lk", "linear_y", 3, 0.1),
+    #(reference_field, "vet", "linear_x", 2, 0.1),
+    ## (reference_field, 'vet', 'linear_x', 3, 9),
+    ## (reference_field, 'vet', 'linear_y', 2, 9),
+    #(reference_field, "vet", "linear_y", 3, 0.1),
+    #(reference_field, "proesmans", "linear_x", 2, 0.45),
+    #(reference_field, "proesmans", "linear_y", 2, 0.45),
+    #(reference_field, "darts", "linear_x", 9, 20),
+    #(reference_field, "darts", "linear_y", 9, 20),
     (reference_field, "farneback", "linear_x", 2, 50),
     (reference_field, "farneback", "linear_y", 2, 50),
+    (reference_field, "farneback", "linear_x", 2, 30),
+    (reference_field, "farneback", "linear_y", 2, 30),
+    (reference_field, "farneback", "linear_x", 2, 20),
+    (reference_field, "farneback", "linear_y", 2, 20),
+    (reference_field, "farneback", "linear_x", 2, 15),
+    (reference_field, "farneback", "linear_y", 2, 15),
+    (reference_field, "farneback", "linear_x", 2, 9),
+    (reference_field, "farneback", "linear_y", 2, 9),
+    (reference_field, "farneback", "linear_x", 2, 6),
+    (reference_field, "farneback", "linear_y", 2, 6),
+    (reference_field, "farneback", "linear_x", 2, 4),
+    (reference_field, "farneback", "linear_y", 2, 4),
 ]
 
 

From ace5a927ef623328c921f1ce05804e4f92b2e1f2 Mon Sep 17 00:00:00 2001
From: dominique-brunet-eccc
 <136009011+dominique-brunet-eccc@users.noreply.github.com>
Date: Tue, 24 Mar 2026 14:58:16 -0400
Subject: [PATCH 27/32] Update test_motion.py adjusting rel RMSE for Farneback

---
 pysteps/tests/test_motion.py | 40 +++++++++++++-----------------------
 1 file changed, 14 insertions(+), 26 deletions(-)

diff --git a/pysteps/tests/test_motion.py b/pysteps/tests/test_motion.py
index 63d4b4179..ede3d92d9 100644
--- a/pysteps/tests/test_motion.py
+++ b/pysteps/tests/test_motion.py
@@ -152,32 +152,20 @@ def _create_observations(input_precip, motion_type, num_times=9):
 )
 
 convergence_arg_values = [
-    #(reference_field, "lk", "linear_x", 2, 0.1),
-    #(reference_field, "lk", "linear_y", 2, 0.1),
-    #(reference_field, "lk", "linear_x", 3, 0.1),
-    #(reference_field, "lk", "linear_y", 3, 0.1),
-    #(reference_field, "vet", "linear_x", 2, 0.1),
-    ## (reference_field, 'vet', 'linear_x', 3, 9),
-    ## (reference_field, 'vet', 'linear_y', 2, 9),
-    #(reference_field, "vet", "linear_y", 3, 0.1),
-    #(reference_field, "proesmans", "linear_x", 2, 0.45),
-    #(reference_field, "proesmans", "linear_y", 2, 0.45),
-    #(reference_field, "darts", "linear_x", 9, 20),
-    #(reference_field, "darts", "linear_y", 9, 20),
-    (reference_field, "farneback", "linear_x", 2, 50),
-    (reference_field, "farneback", "linear_y", 2, 50),
-    (reference_field, "farneback", "linear_x", 2, 30),
-    (reference_field, "farneback", "linear_y", 2, 30),
-    (reference_field, "farneback", "linear_x", 2, 20),
-    (reference_field, "farneback", "linear_y", 2, 20),
-    (reference_field, "farneback", "linear_x", 2, 15),
-    (reference_field, "farneback", "linear_y", 2, 15),
-    (reference_field, "farneback", "linear_x", 2, 9),
-    (reference_field, "farneback", "linear_y", 2, 9),
-    (reference_field, "farneback", "linear_x", 2, 6),
-    (reference_field, "farneback", "linear_y", 2, 6),
-    (reference_field, "farneback", "linear_x", 2, 4),
-    (reference_field, "farneback", "linear_y", 2, 4),
+    (reference_field, "lk", "linear_x", 2, 0.1),
+    (reference_field, "lk", "linear_y", 2, 0.1),
+    (reference_field, "lk", "linear_x", 3, 0.1),
+    (reference_field, "lk", "linear_y", 3, 0.1),
+    (reference_field, "vet", "linear_x", 2, 0.1),
+    # (reference_field, 'vet', 'linear_x', 3, 9),
+    # (reference_field, 'vet', 'linear_y', 2, 9),
+    (reference_field, "vet", "linear_y", 3, 0.1),
+    (reference_field, "proesmans", "linear_x", 2, 0.45),
+    (reference_field, "proesmans", "linear_y", 2, 0.45),
+    (reference_field, "darts", "linear_x", 9, 20),
+    (reference_field, "darts", "linear_y", 9, 20),
+    (reference_field, "farneback", "linear_x", 2, 28),
+    (reference_field, "farneback", "linear_y", 2, 28),
 ]
 
 

From 52113d3f1f4e52057b0cb4d402c8a836f7bdaf68 Mon Sep 17 00:00:00 2001
From: Dominique Brunet <dominique.brunet@ec.gc.ca>
Date: Tue, 24 Mar 2026 19:43:59 +0000
Subject: [PATCH 28/32] Black formatting

---
 pysteps/tests/test_motion_farneback.py | 135 +++++++++++++++++++++++++
 1 file changed, 135 insertions(+)
 create mode 100644 pysteps/tests/test_motion_farneback.py

diff --git a/pysteps/tests/test_motion_farneback.py b/pysteps/tests/test_motion_farneback.py
new file mode 100644
index 000000000..674b80c79
--- /dev/null
+++ b/pysteps/tests/test_motion_farneback.py
@@ -0,0 +1,135 @@
+import pytest
+import numpy as np
+
+from pysteps.motion import farneback
+from pysteps.exceptions import MissingOptionalDependency
+from pysteps.tests.helpers import get_precipitation_fields
+
+fb_arg_names = (
+    "pyr_scale",
+    "levels",
+    "winsize",
+    "iterations",
+    "poly_n",
+    "poly_sigma",
+    "flags",
+    "size_opening",
+    "sigma",
+    "verbose",
+)
+
+fb_arg_values = [
+    (0.5, 3, 15, 3, 5, 1.1, 0, 3, 60.0, False),  # default
+    (0.5, 1, 5, 1, 7, 1.5, 0, 0, 0.0, True),  # minimal settings, sigma=0, verbose
+    (
+        0.3,
+        5,
+        30,
+        10,
+        7,
+        1.5,
+        1,
+        5,
+        10.0,
+        False,
+    ),  # maximal settings, flags=1, big opening
+    (0.5, 3, 15, 3, 5, 1.1, 0, 0, 60.0, True),  # no opening, verbose
+]
+
+
+@pytest.mark.parametrize(fb_arg_names, fb_arg_values)
+def test_farneback_params(
+    pyr_scale,
+    levels,
+    winsize,
+    iterations,
+    poly_n,
+    poly_sigma,
+    flags,
+    size_opening,
+    sigma,
+    verbose,
+):
+    """Test Farneback with various parameters and input types."""
+    pytest.importorskip("cv2")
+    # Input: realistic precipitation fields
+    precip, metadata = get_precipitation_fields(
+        num_prev_files=2,
+        num_next_files=0,
+        return_raw=False,
+        metadata=True,
+        upscale=2000,
+    )
+    precip = precip.filled()
+
+    output = farneback.farneback(
+        precip,
+        pyr_scale=pyr_scale,
+        levels=levels,
+        winsize=winsize,
+        iterations=iterations,
+        poly_n=poly_n,
+        poly_sigma=poly_sigma,
+        flags=flags,
+        size_opening=size_opening,
+        sigma=sigma,
+        verbose=verbose,
+    )
+
+    assert isinstance(output, np.ndarray)
+    assert output.shape[0] == 2
+    assert output.shape[1:] == precip[0].shape
+    assert np.isfinite(output).all() or np.isnan(output).any()
+
+
+def test_farneback_invalid_shape():
+    """Test error when input is wrong shape."""
+    pytest.importorskip("cv2")
+    arr = np.random.rand(64, 64)
+    with pytest.raises(ValueError):
+        farneback.farneback(arr)
+
+
+def test_farneback_nan_input():
+    """Test NaN handling in input."""
+    pytest.importorskip("cv2")
+    arr = np.random.rand(2, 64, 64)
+    arr[0, 0, 0] = np.nan
+    arr[1, 10, 10] = np.inf
+    result = farneback.farneback(arr)
+    assert result.shape == (2, 64, 64)
+
+
+def test_farneback_cv2_missing(monkeypatch):
+    """Test MissingOptionalDependency when cv2 is not injected."""
+    monkeypatch.setattr(farneback, "CV2_IMPORTED", False)
+    arr = np.random.rand(2, 64, 64)
+    with pytest.raises(MissingOptionalDependency):
+        farneback.farneback(arr)
+    monkeypatch.setattr(farneback, "CV2_IMPORTED", True)  # restore
+
+
+def test_farneback_sigma_zero():
+    """Test sigma=0 disables smoothing logic."""
+    pytest.importorskip("cv2")
+    arr = np.random.rand(2, 32, 32)
+    result = farneback.farneback(arr, sigma=0.0)
+    assert isinstance(result, np.ndarray)
+    assert result.shape == (2, 32, 32)
+
+
+def test_farneback_small_window():
+    """Test winsize edge case behavior."""
+    pytest.importorskip("cv2")
+    arr = np.random.rand(2, 16, 16)
+    result = farneback.farneback(arr, winsize=3)
+    assert result.shape == (2, 16, 16)
+
+
+def test_farneback_verbose(capsys):
+    """Test that verbose produces output (side-effect only)."""
+    pytest.importorskip("cv2")
+    arr = np.random.rand(2, 16, 16)
+    farneback.farneback(arr, verbose=True)
+    out = capsys.readouterr().out
+    assert "Farneback method" in out or "mult factor" in out or "---" in out

From 78090fb37ac01f245f1329da857900a507822da9 Mon Sep 17 00:00:00 2001
From: dominique-brunet-eccc
 <136009011+dominique-brunet-eccc@users.noreply.github.com>
Date: Tue, 24 Mar 2026 16:56:22 -0400
Subject: [PATCH 29/32] Update farneback.py for a more verbose description.

---
 pysteps/motion/farneback.py | 12 +++++++++++-
 1 file changed, 11 insertions(+), 1 deletion(-)

diff --git a/pysteps/motion/farneback.py b/pysteps/motion/farneback.py
index 84ac734f5..2a1f7984e 100644
--- a/pysteps/motion/farneback.py
+++ b/pysteps/motion/farneback.py
@@ -50,7 +50,17 @@ def farneback(
     sigma=60.0,
     verbose=False,
 ):
-    """Run the Farneback optical flow routine.
+    """Estimate a dense motion field from a sequence of 2D images using
+    the `Farneback`_ optical flow algorithm.
+
+    This function computes dense optical flow between each pair of consecutive
+    input frames using OpenCV's Farneback method. If more than two frames are
+    provided, the motion fields estimated from all consecutive pairs are
+    averaged to obtain a single representative advection field.
+
+    After the pairwise motion fields are averaged, the resulting motion field
+    can optionally be smoothed with a Gaussian filter. In that case, its
+    amplitude is rescaled so that the mean squared flow magnitude is preserved.
 
     .. _OpenCV: https://opencv.org/
 

From 0e9dd453273a7dac370c8ad6d45ea843e3fa2ae0 Mon Sep 17 00:00:00 2001
From: dominique-brunet-eccc
 <136009011+dominique-brunet-eccc@users.noreply.github.com>
Date: Tue, 24 Mar 2026 17:03:00 -0400
Subject: [PATCH 30/32] Update farneback.py add reference for smoothing and
 renormalization of optical flow,

---
 pysteps/motion/farneback.py | 4 +++-
 1 file changed, 3 insertions(+), 1 deletion(-)

diff --git a/pysteps/motion/farneback.py b/pysteps/motion/farneback.py
index 2a1f7984e..770eba17e 100644
--- a/pysteps/motion/farneback.py
+++ b/pysteps/motion/farneback.py
@@ -156,7 +156,9 @@ def farneback(
     ----------
     Farnebäck, G.: Two-frame motion estimation based on polynomial expansion, 
     In Image Analysis, pages 363–370. Springer, 2003.
-
+    Driedger, N., Mahidjiba, A. and Hortal, A.P. (2022, June 1-8): Evaluation of optical flow
+    methods for radar precipitation extrapolation.
+    Canadian Meteorological and Oceanographic Society Congress, contributed abstract 11801.
     """
 
     if len(input_images.shape) != 3:

From 5bb3be8c01f935cc2a1502a7ade85b9411df09b3 Mon Sep 17 00:00:00 2001
From: dominique-brunet-eccc
 <136009011+dominique-brunet-eccc@users.noreply.github.com>
Date: Tue, 24 Mar 2026 17:05:49 -0400
Subject: [PATCH 31/32] Update description of farneback smoothed in
 plot_optical_flow.py

---
 examples/plot_optical_flow.py | 3 ++-
 1 file changed, 2 insertions(+), 1 deletion(-)

diff --git a/examples/plot_optical_flow.py b/examples/plot_optical_flow.py
index 2e5d92da5..cf196db08 100644
--- a/examples/plot_optical_flow.py
+++ b/examples/plot_optical_flow.py
@@ -147,7 +147,8 @@
 #
 # This module implements the pyramidal decomposition method for motion estimation
 # of Farnebäck as implemented in OpenCV, with an option for smoothing and
-# renormalization of the motion fields proposed by N. Driedger.
+# renormalization of the motion fields proposed by Driedger et al.:
+# https://cmosarchives.ca/Congress_P_A/program_abstracts2022.pdf (p. 392).
 
 oflow_method = motion.get_method("farneback")
 R[~np.isfinite(R)] = metadata["zerovalue"]

From f7fc86c0a4edfc8a31dcdf990c072ee73ee8b1ab Mon Sep 17 00:00:00 2001
From: dominique-brunet-eccc
 <136009011+dominique-brunet-eccc@users.noreply.github.com>
Date: Tue, 24 Mar 2026 17:31:29 -0400
Subject: [PATCH 32/32] Update wording in farneback.py method description

---
 pysteps/motion/farneback.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/pysteps/motion/farneback.py b/pysteps/motion/farneback.py
index 770eba17e..79eebeed8 100644
--- a/pysteps/motion/farneback.py
+++ b/pysteps/motion/farneback.py
@@ -60,7 +60,7 @@ def farneback(
 
     After the pairwise motion fields are averaged, the resulting motion field
     can optionally be smoothed with a Gaussian filter. In that case, its
-    amplitude is rescaled so that the mean squared flow magnitude is preserved.
+    amplitude is rescaled so that the mean motion magnitude is preserved.
 
     .. _OpenCV: https://opencv.org/