光流法：Farneback

光流法：Farnback

---

现实世界中，万物都在在运动，且运动的速度和方向可能均不同，这就构成了运动场。物体的运动投影在图像上反应的是像素的移动。这种像素的瞬时移动速度就是光流。光流法是利用图像序列中的像素在时间域上的变化、相邻帧之间的相关性来找到的上一帧跟当前帧间存在的对应关系，计算出相邻帧之间物体的运动信息的一种方法。

光流法按照不同的实现方式可以分为：基于梯度的方法、基于匹配的方法、基于能量的方法、基于相位的方法等。本文介绍的是一种基于梯度的经典光流方法：Farnback法。光流法的前提假设包括：相邻帧之间亮度恒定；相邻帧之间取时间连续或者运动变化微小；同一子图像中像素点具有相同的运动。

基本假设

假定图像序列记作$I(x,y,t)$,其中$X = [x,y]$。视频中的每个前后帧提取出来之后就是一个图像序列。假设图像亮度恒定，即图像亮度没有变化，则导数为0：

$$\frac{dI( X, t)}{d t} = \frac{\partial I}{\partial X}\frac{\partial X}{\partial t}+\frac{\partial I}{\partial t} = 0$$
或者根据泰勒展开来得出上述式子：
$$I( X , t) = I( X + \Delta X, t + \Delta t) \approx I( X , t) + \frac{\partial I}{\partial X}\Delta X+\frac{\partial I}{\partial t}\Delta t \\ \frac{\partial I}{\partial X}\Delta X+\frac{\partial I}{\partial t}\Delta t = 0\\ \frac{\partial I}{\partial X}\frac{\Delta X}{\Delta t}+\frac{\partial I}{\partial t}=0\\ \frac{\partial I}{\partial x}\frac{\partial x}{\partial t}+\frac{\partial I}{\partial y}\frac{\partial y}{\partial t}+\frac{\partial I}{\partial t}=0$$
其中，在微小时间内 $\frac{\partial X}{\partial t}$或者 $\frac{\Delta X}{\Delta t}$表示速度，可以记为：
$$\frac{\partial X}{\partial t} = [\frac{\partial x}{\partial t}, \frac{\partial y}{\partial t}]=[u, v] = d$$
则有：
$$I_xu+I_yv+I_t=0\\ \begin{bmatrix} I_x \ \ I_y \end{bmatrix} \begin{bmatrix} u \\ v \end{bmatrix} = -I_t$$

Farneback光流法

Farneback是一种基于梯度的方法，假设图像梯度恒定且假设局部光流恒定。局部光流恒定，即对于任意的$y\in N(x), d = \frac{\partial X}{\partial t} 不变$。梯度恒定即：

$$\frac{d}{dt}\nabla I( X, t) = \frac{\partial \nabla I}{\partial X}\frac{\partial X}{\partial t}+\frac{\partial \nabla I}{\partial t} = H(I) \cdot d + ( \nabla I)_t =0$$
假设：
$$E( X, d) = ||(H(I)\cdot d +(\nabla I)_t)||^2$$
上式在最优值处有导数为0：
$$\frac{\partial E}{\partial d} = 0 \\ d= -\left(H^T(I)H(I)\right)^{-1}(H^T(I)(\nabla I)_t)$$

若对时间离散化：
（后向差分）

$$( \nabla I)_t ( X, t) \approx \nabla I(X,t) - \nabla I(X,t-1)$$
（时间中心差分）
$$( \nabla I)_t ( X, t-1/2) \approx \nabla I(X,t) - \nabla I(X,t-1)\\ H(I)(X,t-1/2) \approx \frac{1}{2}(H(I)(X,t)+H(I)(X,t-1))$$

图像模型

图像一般是二维的（灰度图像），那么图像像素点的灰度值可以看成是一个二维变量的函数$f(x,y)$。假设以感兴趣的像素点为中心，构建一个局部坐标系（并不是针对整张图像）。对函数进行二项展开，可以近似为：

$$\begin{align} f(x,y) &\approx r_1+r_2x+r_3y+r_4x^2+r5y^2+r_6xy\\&= \left( \begin{matrix} x & y \end{matrix} \right)^T \left( \begin{matrix} r_{4} & r_{6}/2 \\ r_{6}/2 & r_{5} \end{matrix} \right) \left( \begin{matrix} x \\ y \end{matrix} \right)+ \left( \begin{matrix} r_{2} \\ r_{3} \end{matrix} \right)^T \left( \begin{matrix} x \\ y \end{matrix} \right)+r_{1}\\ &=\mathbf{x^TAx+b^Tx+}c \end{align}$$

其中，$\bf x$为二维列向量， $\bf A$为$2\times 2$的对称矩阵，$\bf b$为$2\times 1$的矩阵。注意，此处的系数确定后只针对在确定点$(x,y)$而言，对于其他点可能并不适用，也就是说，每个像素点对应一组系数。

取该像素点的一个邻域（通常以该像素为中心，大小为$2n+1的方形区域$），利用这些像素点的值和坐标来进行系数的估计，估计的算法可以使用加权最小二乘法。加权是因为在邻域内，距中心越近的像素点与中心像素具有更大的相关性，而越远的点提供的信息越少。其实可以将邻域以外的像素点的权重都看成是0。

位移估计

考虑多项式扩展是在一个像素的邻域内，如果像素经过移动$d$后，则整个多项式应该会发生变化。
原始位置：

$$f_1({\bf x}) = {\bf x^TA_1x+b_1^Tx}+c_1$$
像素移动后：

$$\begin{aligned} f_2({\bf x}) &= f_1({\bf x-d})\\ &= {\bf (x-d)^TA_1(x-d)+b_1^T(x-d)}+c_1\\ &= {\bf x^TA_1x+(b_1-2A_1d)^Tx + d^TA_1d-b_1^Td}+c_1\\ &= {\bf x^TA_2x+b_2^Tx}+c_2 \end{aligned}$$
其中， $${\bf A_2 = A_1 \\b_2= b_1-2A_1d} \\c_2={\bf d^TA_1d-b_1^Td}+c_1$$
如果 $\bf A_1$非奇异，则有上述的第二个式子可以得到：
$${\bf d} = -\frac{1}{2}\bf A_1^{-1}(b_2 - b_1)$$
按照理论推导，其中必定有 $\bf A_1 = A_2$,但实际情况中未必能满足这一项要求，因此可以通过来求平均来近似真实值。如果令：
$${\bf A(x)} = \frac{{\bf A_1(x) + A_2(x)}}{2}\\ \Delta {\bf b(x)} =-\frac{1}{2} \bf(b_2 - b_1)$$
那么：
$$\bf A(x)d(x) =\Delta {\bf b(x)} \\ d = (A^TA)^{-1}(A^T\Delta b)$$

可以构建目标函数来进行优化求得位移：

$$e({\bf x} ) = ||{\bf Ad - \Delta b}||^2$$
实际情况中，这种方法求得的结果中噪声太多，因此可以使用兴趣像素点的邻域，然后使用加权的目标函数：
$$e({\bf x} ) =\sum_{\Delta x \in I} w(\Delta x)||{\bf A(x + \Delta x)d - \Delta b(x+\Delta x)}||^2$$

Reference

[1] 图像分析之光流之经典
[2] Farneback 光流算法详解与 calcOpticalFlowFarneback 源码分析
[3] 光流Optical Flow介绍与OpenCV的实现
[4] 光流法简单介绍
[5] Farneback, 2003, Two-Frame Motion Estimation Based on Polynomial Expansion