Farneback Optical Flow

光流算法在计算机视觉领域扮演着重要角色，用于估计图像序列中物体的运动。在众多光流算法中，Farneback算法因其在密集光流估计方面的出色表现而备受关注。本文将深入剖析Farneback算法的数学原理、实现细节及其广泛应用。

1. 光流的基本概念

#

在深入Farneback算法之前，我们需要理解光流的基本概念。光流是描述图像亮度模式随时间变化的二维向量场，每个向量代表对应像素点的位移。
详细可看：

1.1 光流的数学表达

#

假设在时间 \(t\) 时，空间位置 \((x,y)\) 处的图像亮度为 \(I(x,y,t)\)，那么光流的基本约束方程可以表示为：

对上式进行泰勒展开并保留一阶项，可得：

整理可得光流约束方程：

其中，\(u = \frac{dx}{dt}\) 和 \(v = \frac{dy}{dt}\) 分别表示水平和垂直方向的光流速度。

1.2 光流的基本假设

#

光流估计通常基于以下假设：

1. 亮度恒定假设：同一物体在不同帧中的亮度保持不变
2. 小位移假设：连续帧之间的物体运动较小
3. 空间一致性假设：相邻像素具有相似的运动

然而，这些假设在实际场景中并不总是成立，这也是光流估计具有挑战性的原因之一。

1.3 光流类别：

#

• 当描述部分像素时，称为：稀疏光流(跟踪角点、固定位置等)
• 当描述全部像素时，称为：稠密光流(全局运动估计、背景建模等)
  以下是稀疏光流（Sparse Optical Flow）和稠密光流（Dense Optical Flow）的对比表格：

2. Farneback

#

Farneback算法由Gunnar Farnebäck于2003年提出，它是一种基于多项式展开的密集光流估计方法。与传统的Lucas-Kanade或Horn-Schunck等基于梯度的方法不同，Farneback算法采用了多项式信号模型来表示图像局部区域，能够有效处理较大范围的运动。

2.1 基本概念

#

更详细的请自行搜索或问AI。

图像金字塔

#

想象你有一张高清照片，现在你不断缩小它，得到一系列越来越模糊、越来越小的版本，就像金字塔一样从底到顶越来越小。这就是图像金字塔。

• 作用：让计算机能在不同尺度（大小）下分析图像，比如检测远处的小目标和近处的大目标。
• 常见类型：
  • 高斯金字塔（不断模糊+缩小）
  • 拉普拉斯金字塔（记录不同尺度下的细节差异）

高斯模糊

#

高斯模糊就像让照片“变柔焦”，让图像看起来更平滑，减少噪点（比如照片上的小颗粒）。

• 怎么做的？
  每个像素的新值 = 周围像素的加权平均，越靠近中心的像素权重越大。
• 为什么叫“高斯”？
  因为权重的分布符合高斯函数（类似钟形曲线）。
• 用途：
  • 去噪（让图像更干净）
  • 降低细节（比如人脸美化）

多项式拟合图像窗口

#

一个图像的5x5的窗口（小方块）一般需要25个像素点来描述。
有时候我们想减少参数来描述这个窗口，比如用一个简单的数学公式（多项式）来表示。

多项式展开

#

举个直观例子🌰： 假设你有一个5×5像素的小方块，它的灰度值分布如下：

100 105 110 115 120  
105 110 115 120 125  
110 115 120 125 130  
115 120 125 130 135  
120 125 130 135 140

你会发现：

• 水平方向：每向右一列，灰度值+5（线性变化）

• 垂直方向：每向下一行，灰度值+5（线性变化）

但如果这个小方块的灰度变化更复杂（比如有弯曲的边缘或纹理），线性函数就不够用了，这时就需要二次多项式：

表示为矩阵乘法形式：

2. 紧凑矩阵表示
定义： • 坐标向量：\(\mathbf{x} = \begin{bmatrix} x \\ y \end{bmatrix}\)

• 二次项系数矩阵：\(\mathbf{A} = \begin{bmatrix} a & c/2 \\ c/2 & b \end{bmatrix}\)

• 一次项系数向量：\(\mathbf{b} = \begin{bmatrix} d \\ e \end{bmatrix}\)

• 常数项：\(c = f\)

二次多项式用6个参数能更好地描述图像的局部结构（边缘、角点、纹理）！

多项式系数估计

#

怎么找到这些系数（a, b, c…）？
用最小二乘法（Least Squares）：让数学公式计算的亮度与实际图像的亮度误差最小。

• 步骤：
  1. 取图像窗口的所有像素点 \((x_i, y_i)\) 和它们的亮度 \(I(x_i, y_i)\)。
  2. 解方程，找到最合适的 \(a, b, c...\) 使误差最小。

2.2 Farneback算法的基本流程

#

2.2.1. 图像金字塔构建

#

Farneback算法首先构建高斯金字塔，以处理不同尺度的运动：

• 对输入图像进行高斯模糊，然后降采样（缩小尺寸），生成多层金字塔。
• 从最顶层（最小分辨率）开始计算光流，逐步向下细化，提高精度。

为什么需要金字塔？

• 大运动在低分辨率下更容易捕捉（小位移≈大运动）。
• 逐层优化，避免直接在大分辨率下计算导致误差过大。

下面的文章由粗糙到精细部分有介绍：

2.2.2. 局部多项式拟合

#

对每一层金字塔图像，用二次多项式拟合每个像素的邻域（如 5×5 窗口）：

其中：

• \( x = \begin{bmatrix} x \\ y \end{bmatrix} \) 是像素坐标（以窗口中心为原点）
• \( A_1 \) 是2×2对称矩阵，表示二阶项
• \( b_1 \) 是2×1向量，表示一阶项
• \( c_1 \) 是常数项

拟合方法：通过最小二乘法计算多项式系数。

2.2.3. 光流方程推导

#

位移假设

#

假设在下一帧图像中，该窗口发生了位移 \( d = \begin{bmatrix} u \\ v \end{bmatrix} \)，则新图像可以表示为：

新帧的多项式表示

#

同样对新帧用二次多项式表示：

等式联立

#

比较两个多项式表达式，可以得到：

光流方程

#

从第二个等式可以得到光流方程：

位移求解

#

对于每个像素点，解这个线性方程组：

解得：

实际计算步骤

#

1. 对两帧图像分别计算每个像素的多项式系数（A,b,c）
2. 计算系数差 Δb = b₂ - b₁
3. 对每个像素解上述方程得到(u,v)

2.2.4. 迭代优化

#

为提高精度，Farneback算法采用迭代式优化：

1. 从粗到精：先在金字塔顶层计算粗略光流，再将结果传递到下一层作为初始值。
2. 局部加权：对每个像素的光流进行加权平滑，减少噪声影响（类似高斯模糊）。

2.2.5. 输出稠密光流场

#

最终，算法为每个像素计算一个位移矢量 \((d_x, d_y)\)，形成稠密光流场。

• 可用箭头图或颜色编码可视化（如OpenCV的cv2.calcOpticalFlowFarneback）。

3. Farneback算法的实现

#

下面我们将介绍如何使用OpenCV库实现Farneback算法。OpenCV提供了calcOpticalFlowFarneback函数，可以直接用于计算密集光流。

import cv2
import numpy as np
import matplotlib.pyplot as plt

def compute_farneback_flow(prev_frame, curr_frame):
    # 转换为灰度图
    prev_gray = cv2.cvtColor(prev_frame, cv2.COLOR_BGR2GRAY)
    curr_gray = cv2.cvtColor(curr_frame, cv2.COLOR_BGR2GRAY)
    
    # 计算光流
    # 参数说明:
    # prev: 前一帧图像
    # next: 当前帧图像
    # flow: 输出的光流
    # pyr_scale: 金字塔尺度，通常设为0.5，表示每一层的尺寸是上一层的一半
    # levels: 金字塔的层数
    # winsize: 窗口大小，用于寻找多项式展开的邻域
    # iterations: 每一层的迭代次数
    # poly_n: 用于多项式展开的邻域大小
    # poly_sigma: 高斯标准差，用于平滑导数
    # flags: 计算标志
    flow = cv2.calcOpticalFlowFarneback(
        prev_gray, 
        curr_gray, 
        None, 
        pyr_scale=0.5, 
        levels=3, 
        winsize=15, 
        iterations=3, 
        poly_n=5, 
        poly_sigma=1.2, 
        flags=0
    )
    
    return flow

3.1 光流可视化

#

光流场的可视化对于理解和分析结果非常重要。常用的可视化方法有两种：颜色编码和箭头表示。

def visualize_flow(flow, type="color"):
    """
    可视化光流场
    flow: 光流场
    type: 可视化类型，"color"为HSV颜色编码，"arrows"为箭头表示
    """
    if type == "color":
        # 计算光流的大小和角度
        mag, ang = cv2.cartToPolar(flow[..., 0], flow[..., 1])
        
        # 创建HSV图像
        hsv = np.zeros((flow.shape[0], flow.shape[1], 3), dtype=np.uint8)
        
        # 将角度映射到色调（H）通道
        hsv[..., 0] = ang * 180 / np.pi / 2
        
        # 将大小映射到饱和度（S）通道，并归一化
        hsv[..., 1] = cv2.normalize(mag, None, 0, 255, cv2.NORM_MINMAX)
        
        # 将明度（V）通道设为最大
        hsv[..., 2] = 255
        
        # 转换HSV到RGB
        rgb = cv2.cvtColor(hsv, cv2.COLOR_HSV2BGR)
        
        return rgb
    
    elif type == "arrows":
        # 创建网格
        step = 16
        h, w = flow.shape[:2]
        y, x = np.mgrid[step/2:h:step, step/2:w:step].reshape(2, -1)
        fx, fy = flow[y.astype(int), x.astype(int)].T
        
        # 创建画布
        plt.figure(figsize=(12, 8))
        plt.imshow(np.zeros((h, w)), cmap='gray')
        plt.quiver(x, y, fx, fy, color='r', angles='xy', scale_units='xy', scale=0.25)
        plt.axis('off')
        
        return plt

3.2 完整实现

#

下面是一个完整的示例，展示如何应用Farneback算法处理视频：

import cv2
import numpy as np
import matplotlib.pyplot as plt

def compute_farneback_flow(prev_frame, curr_frame):
    # 转换为灰度图
    prev_gray = cv2.cvtColor(prev_frame, cv2.COLOR_BGR2GRAY)
    curr_gray = cv2.cvtColor(curr_frame, cv2.COLOR_BGR2GRAY)
    
    # 计算光流
    flow = cv2.calcOpticalFlowFarneback(
        prev_gray, 
        curr_gray, 
        None, 
        pyr_scale=0.5,  # 金字塔尺度，通常设为0.5，表示每一层的尺寸是上一层的一半
        levels=3,       # 金字塔的层数
        winsize=15,     # 窗口大小，用于寻找多项式展开的邻域
        iterations=3,   # 每一层的迭代次数
        poly_n=5,       # 用于多项式展开的邻域大小
        poly_sigma=1.2, # 高斯标准差，用于平滑导数
        flags=0         # 计算标志
    )
    
    return flow

def draw_flow(img, flow, step=16):
    """
    绘制光流场的箭头表示
    img: 输入图像
    flow: 光流场
    step: 网格步长，控制绘制密度
    """
    h, w = img.shape[:2]
    y, x = np.mgrid[step/2:h:step, step/2:w:step].reshape(2, -1).astype(int)
    fx, fy = flow[y, x].T
    
    # 创建线条终点
    lines = np.vstack([x, y, x + fx, y + fy]).T.reshape(-1, 2, 2)
    lines = np.int32(lines)
    
    # 创建图像副本
    vis = img.copy()
    
    # 绘制线条
    for (x1, y1), (x2, y2) in lines:
        cv2.line(vis, (x1, y1), (x2, y2), (0, 255, 0), 1)
        cv2.circle(vis, (x1, y1), 1, (0, 255, 0), -1)
    
    return vis

def flow_to_color(flow):
    """
    将光流场转换为HSV颜色编码图像
    flow: 光流场
    """
    # 计算光流的大小和角度
    mag, ang = cv2.cartToPolar(flow[..., 0], flow[..., 1])
    
    # 创建HSV图像
    hsv = np.zeros((flow.shape[0], flow.shape[1], 3), dtype=np.uint8)
    
    # 将角度映射到色调（H）通道
    hsv[..., 0] = ang * 180 / np.pi / 2
    
    # 将大小映射到饱和度（S）通道，并归一化
    hsv[..., 1] = cv2.normalize(mag, None, 0, 255, cv2.NORM_MINMAX)
    
    # 将明度（V）通道设为最大
    hsv[..., 2] = 255
    
    # 转换HSV到BGR
    bgr = cv2.cvtColor(hsv, cv2.COLOR_HSV2BGR)
    
    return bgr

def main():
    # 打开视频文件或摄像头
    cap = cv2.VideoCapture("video.mp4")  # 或者 cap = cv2.VideoCapture(0) 使用摄像头
    
    # 读取第一帧
    ret, prev_frame = cap.read()
    if not ret:
        print("无法读取视频")
        return
    
    # 创建输出视频
    fourcc = cv2.VideoWriter_fourcc(*'XVID')
    out_color = cv2.VideoWriter('farneback_color.avi', fourcc, 20.0, 
                               (prev_frame.shape[1], prev_frame.shape[0]))
    out_arrows = cv2.VideoWriter('farneback_arrows.avi', fourcc, 20.0, 
                                (prev_frame.shape[1], prev_frame.shape[0]))
    
    frame_count = 0
    
    while True:
        # 读取当前帧
        ret, curr_frame = cap.read()
        if not ret:
            break
        
        # 每隔5帧计算一次光流(可根据需要调整)
        if frame_count % 5 == 0:
            # 计算光流
            flow = compute_farneback_flow(prev_frame, curr_frame)
            
            # 可视化光流
            flow_color = flow_to_color(flow)
            flow_arrows = draw_flow(curr_frame, flow)
            
            # 保存结果
            out_color.write(flow_color)
            out_arrows.write(flow_arrows)
            
            # 显示结果
            cv2.imshow("原始视频", curr_frame)
            cv2.imshow("光流颜色编码", flow_color)
            cv2.imshow("光流箭头", flow_arrows)
            
            # 更新前一帧
            prev_frame = curr_frame.copy()
        
        frame_count += 1
        
        # 按ESC退出
        if cv2.waitKey(5) & 0xFF == 27:
            break
    
    # 释放资源
    cap.release()
    out_color.release()
    out_arrows.release()
    cv2.destroyAllWindows()

if __name__ == "__main__":
    main()

4. Farneback算法的数学推导

#

为了更深入理解Farneback算法，让我们对其核心步骤进行详细的数学推导。这将帮助我们把握算法的本质，并为参数调优提供理论基础。

4.1 多项式展开的详细推导

#

在Farneback算法中，我们使用二次多项式来近似图像的局部区域。对于一个像素点 \(x = (x, y)^T\)，其周围的图像强度可以表示为：

为了求解多项式系数 \(A\)、\(b\) 和 \(c\)，我们需要最小化加权平方误差：

其中，\(\Omega\) 是以 \(x\) 为中心的局部区域，\(w\) 是权重函数，通常选择高斯函数：

设 \(z = y - x\)，那么 \(y = z + x\)，将其代入二次多项式：

因此，多项式系数的计算可以通过解一系列加权最小二乘问题来实现。在实际实现中，这可以通过卷积操作高效完成。

4.2 位移估计的详细推导

#

假设我们有两个连续帧 \(f_1\) 和 \(f_2\)，它们的多项式展开分别为：

若存在一个位移 \(d\)，使得 \(f_2(x) = f_1(x - d)\)，则：

展开上式：

对比两个多项式展开，我们有：

从第二个方程，我们可以解出位移 \(d\)：

若 \(A_1\) 是非奇异矩阵，则：

在实际应用中，为了提高算法的鲁棒性，我们通常采用迭代优化的方法来求解位移场。具体来说，对于每个迭代步骤，我们使用上一步得到的位移场作为初始估计，然后基于残差进行更新。

4.3多尺度策略数学描述

#

真实计算的时候需要金字塔多尺度最优近似解优化结果。
针对一个像素点的计算过程：
以它为中心周围\(\mathbf{poly_n} \times \mathbf{poly_n}\)的小窗口，计算它的多项式系数，即A、b、c。
然后再以它为中心，周围\(\mathbf{winsize} \times \mathbf{winsize}\)的大窗口，即小窗口的移动范围，在这个大窗口中使用最小二乘法找到最优的位移矢量\(d\)或者\(\delta d\)。

• 这里的\(d\)是当前像素点的位移矢量，\(\delta d\)是当前像素点的残差位移矢量。

完整数学流程

1. 金字塔构建
   对两帧图像 \(I_1, I_2\) 分别构建高斯金字塔：

   \[ I_l^k = \text{GaussianBlur}(I_l^{k-1}) \downarrow_2, \quad l \in \{1,2\}, \ k \in \{1,...,L\} \]

2. 顶层粗估计
   在最低分辨率层（尺度 \(L\)）直接计算初始光流：

   \[ d_L = \arg\min_d \sum_{x} w(x) \|A_1(x)d - \frac{1}{2}(b_1(x)-b_2(x))\|^2 \]

3. 逐层优化
   对于每一层 \(k\)（从 \(L-1\) 到 \(0\)）：

   • 上采样光流：\(d_{k+1} \rightarrow 2d_{k+1}↑\)

   • Warp Frame1（反向变形）：

     \[ I_{1\_warp}^k(x) = I_1^k(x - 2d_{k+1}↑) \]

   • 计算残差（基于warp后的Frame1）：

     \[ \delta d_k = \arg\min_{\delta d} \sum_{x} w(x) \|A_{1\_warp}(x)\delta d - (b_{1\_warp}(x)-b_2(x))\|^2 \]

     其中 \(A_{1_warp}\) 是warp后图像的梯度矩阵

   • 合成光流：

     \[ d_k = 2d_{k+1}↑ + \delta d_k \]

4. 输出结果
   最终光流场：

   \[ d = d_0 \]

参数说明：

• 金字塔层数通常3-5层，缩放因子0.5（pyr_scale=0.5）

• 每层迭代次数3-5次（iterations=3）

• 多项式邻域大小5-7像素（poly_n=5）

5. 参数调优与优化

#

Farneback算法的性能在很大程度上取决于参数设置。下面我们将分析关键参数的影响，并提供调优策略。

5.1 金字塔参数

#

• 金字塔层数 (levels)： 控制算法处理的运动范围。层数越多，可以处理的运动范围越大，但计算复杂度也越高。

  对于运动幅度较大的场景，建议设置为3-5层；对于运动幅度较小的场景，2-3层通常已经足够。金字塔层数 \(L\) 与可处理的最大位移 \(d_{max}\) 的关系可以近似为：

  $$ d_{max} \approx \frac{\text{winsize}}{2} \cdot \frac{1}{s^{L-1}} $$

• 金字塔尺度 (pyr_scale)： 控制相邻层之间的缩放比例。较小的缩放比例（如0.5）可以捕获更大范围的运动，但可能会丢失细节；较大的缩放比例（如0.8）可以保留更多细节，但处理大运动的能力受限。

5.2 多项式拟合参数

#

• 窗口大小 (winsize)： 控制局部多项式拟合的邻域大小。较大的窗口提供更平滑的光流场，但可能会模糊运动边界；较小的窗口可以更好地保留边界细节，但对噪声更敏感。

  窗口大小与图像分辨率和运动特性相关。对于高分辨率图像或复杂场景，建议使用较大的窗口（如15-21）；对于低分辨率图像或简单场景，较小的窗口（如7-15）可能更合适。

• 多项式邻域大小 (poly_n)： 控制用于估计多项式展开的像素邻域大小。通常设置为5或7。较大的值可以提供更鲁棒的估计，但计算复杂度更高。

• 高斯标准差 (poly_sigma)： 控制用于平滑导数的高斯核标准差。较大的值会产生更平滑的光流场，但可能会丢失细节；较小的值可以保留更多细节，但对噪声更敏感。

  一般建议设置 \(\text{poly\_sigma}\) = \(\frac{\text{poly\_n} - 1}{2}\)，即 \(\text{poly\_n} = 5\) 时，\(\text{poly\_sigma} \approx 2.0\)；\(\text{poly\_n} = 7\) 时，\(\text{poly\_sigma} \approx 3.0\)。

5.3 迭代参数

#

• 迭代次数 (iterations)： 控制每个金字塔层级上的迭代次数。较多的迭代可以提高精度，但会增加计算时间。通常设置为3-5次。

  迭代次数与误差收敛速度相关。可以通过监控残差变化来自适应地调整迭代次数：

  $$ \text{residual} = \frac{1}{N} \sum_{i=1}^{N} ||f_2(x_i) - f_1(x_i - d_i)||^2 $$

  当残差变化小于预设阈值时，可以提前终止迭代。

5.4 性能优化策略

#

为了在实际应用中提高Farneback算法的性能，可以采用以下优化策略：

1. 区域选择性处理： 只在感兴趣区域(ROI)计算光流，而不是整个图像。这可以显著减少计算量。

2. GPU加速： 利用GPU的并行计算能力可以大幅提升算法性能。OpenCV提供了CUDA版本的Farneback实现。

3. 帧间隔处理： 在帧率较高的视频中，可以考虑每隔几帧计算一次光流，而不是每一帧都计算。

4. 自适应参数： 根据场景特性动态调整算法参数。例如，在运动较大的区域使用更多的金字塔层数和迭代次数。

6. 应用场景与实例分析

#

Farneback算法在多个领域有着广泛的应用。下面我们将介绍几个典型的应用场景，并分析其实现细节。

6.1 运动检测与目标跟踪

#

在视频监控系统中，光流是检测和跟踪移动目标的有效工具。通过分析光流场的模式，可以识别场景中的运动区域，并对移动目标进行分割和跟踪。

实现思路：

1. 使用Farneback算法计算连续帧之间的光流场
2. 计算光流场的幅值，作为运动强度图
3. 对运动强度图进行阈值处理，得到运动掩码
4. 应用形态学操作（如膨胀和腐蚀）来去除噪声和填充空洞
5. 提取连通区域作为移动目标
6. 使用卡尔曼滤波