1 光流 (Optical Flow) 的定义与介绍

1.1 光流 (Optical Flow) 的提出和定义

光流或视流是观察者与场景之间的相对运动引起的视觉场景中物体、表面和边缘的表观运动模式。光流还可以定义为图像中亮度模式运动的表观速度分布。

光流的概念是由美国心理学家 James J. Gibson 在20世纪40年代提出的，用来描述动物在世界中移动时所接收到的视觉刺激。Gibson 强调了光流对于可利用性感知的重要性，即可辨别环境中行动可能性的能力。

图 1.旋转观察者 (在本例中为苍蝇) 所经历的光流。每个位置的光流方向和大小由每个箭头的方向和长度表示。

直观一些理解。我们可以做如下的总结 :

光流的概念是指在连续的两帧图像中由于图像中的物体移动或者摄像头的移动导致的图像中目标像素的移动。
光流是二维矢量场，表示了一个点从第一帧到第二帧的位移。

为了更直观地理解光流，我们用上面的图表示了一个球在连续的 5 帧图像中的运动。箭头表示了它的位移矢量。

1.2 光流 (Optical Flow) 的估计

传统光流法的工作原理基于如下假设 :

场景的像素强度在相邻帧之间基本不变
相邻像素具有相似的运动

亮度不变假设则是期待两帧图像的同一物体的亮度不变，这个假设通常是成立的，因为环境光照通常不会发生太大的变化。假设 $t$ 时刻，位于 $(x,y)$ 像素位置的物体，在 $t+\Delta t$ 时刻位于 $(x+\Delta x,y+\Delta y)$ 位置，基于亮度不变假设，有 :

I(x,y,t)=I(x+\Delta x,y+\Delta y,t+\Delta t)

因为相邻两帧之间的移动会比较小，因此我们对右边的式子使用泰勒展开，在此只展开到一阶 :

I(x+\Delta x,y+\Delta y,t+\Delta t)=I(x,y,t)+{\frac {\partial I}{\partial x}}\,\Delta x+{\frac {\partial I}{\partial y}}\,\Delta y+{\frac {\partial I}{\partial t}}\,\Delta t+O\left(\left(\Delta t\right)^2\right)

我们忽略高阶项带来的影响，只考虑线性化的泰勒展开，即可得到 :

\begin{aligned} {\frac {\partial I}{\partial x}}\Delta x+{\frac {\partial I}{\partial y}}\Delta y+{\frac {\partial I}{\partial t} }\Delta t&=0\\ {\frac {\partial I}{\partial x}}{\frac {\Delta x}{\Delta t}}+{\frac {\partial I}{\partial y}}{\frac {\Delta y}{\Delta t}}+{\frac {\partial I}{\partial t}}&=0 \end{aligned}

因为两帧之间的时间很短，我们可以认为 $V_x=\frac {\Delta x}{\Delta t}$ 和 $V_y=\frac {\Delta y}{\Delta t}$ 是图像的移动速度。将 $\frac {\partial I}{\partial x}, \frac {\partial I}{\partial y}, \frac {\partial I}{\partial t}$ 记为 $I'_x,I'_y,I'_t$ ，其表示光流函数 $I(x,y,t)$ 在 $(x,y,t)$ 上偏导的方向。式子可以写为 :

\begin{aligned} I'_{x}V_{x}+I'_{y}V_{y}&=-I'_{t}\\ \begin{bmatrix} I'_x, \space I'_y \end{bmatrix}\begin{bmatrix}V_{x} \\ V_{y}\end{bmatrix}&=-I'_{t} \end{aligned}

基于这个公式，以下有几种光流的常见估计方式。

1.3 Lucas-Kanade 方法

Lucas-Kanade 方法假设两个邻近时刻（帧）之间图像内容的位移很小，并且在该点的邻域内近似恒定 $p$ 在考虑中。因此，可以假设光流方程对于以 $p$ 为中心的窗口内的所有像素都成立。即局部图像流 (速度) 向量 $(V_x.V_y)$ 必须满足 :

\begin{aligned} I'_{x}(q_{1})V_{x}+I'_{y}(q_{1})V_{y}&=-I'_{t}(q_{1} )\\ I'_{x}(q_{2})V_{x}+I'_{y}(q_{2})V_{y}&=-I'_{t}(q_{2})\\ &\; \ \vdots \\ I'_{x}(q_{n})V_{x}+I'_{y}(q_{n})V_{y}&=-I'_{t}(q_{n}) \end{aligned}

其中 ${q_{1},q_{2},\dots ,q_{n}}$ 是窗口内的像素，这些式子可以转化成矩阵形式 $Av=b$ ，其中 :

\begin{aligned} A={\begin{bmatrix}I'_{x}(q_{1})&I'_{y}(q_{1})\\[10pt] I'_{x}(q_{2})&I'_{y}(q_{2})\\[10pt] \vdots &\vdots \\[10pt] I'_{x}(q_{n})&I'_{y}(q_{n})\end{bmatrix}}\quad \quad v={\begin{bmatrix}V_{x}\\[10pt]V_{y}\end{bmatrix}}\quad \quad b={\begin{bmatrix}-I'_{t}(q_{1})\\[10pt] -I'_{t}(q_{2})\\[10pt] \vdots \\[10pt] -I'_{t}(q_{n})\end{bmatrix}} \end{aligned}

这样就转化成了一个常见的最小二乘矩阵优化 $A^TAv=A^Tb$ 。该算法所需的工作主要是对于图中进行区域的划分。

\begin{aligned} {\begin{bmatrix}V_{x}\\[10pt] V_{y} \end{bmatrix}} ={\begin{bmatrix} \sum _{i}I'_{x}(q_{i} )^{2}&\sum _{i}I'_{x}(q_{i})I'_{y}(q_{i})\\[10pt]\sum _{i}I'_{y}(q_{ i})I'_{x}(q_{i})&\sum _{i}I'_{y}(q_{i})^{2}\end{bmatrix}}^{-1} {\begin{bmatrix}-\sum _{i}I'_{x}(q_{i})I'_{t}(q_{i})\\[10pt] -\sum _{i}I'_{y}(q_{i})I'_{t}(q_{i}) \end{bmatrix}} \end{aligned}

通常，在实验中，我们给予更中心的像素更多权重 $A^{T}WAv=A^{T}Wb$ ，其中 $W$ 是一个 $n\times n$ 的包含权重的对角矩阵 $W=\operatorname{diag}(w_i)$ ， $w_i$ 为像素 $q_i$ 对应的权重，经常设置为 $q_i$ 和 $p$ 之间的高斯距离，即 :

\begin{aligned} {\begin{bmatrix}V_{x}\\[10pt]V_{y}\end{bmatrix}}={\begin{bmatrix}\sum _{i}w_{i}I'_{x}(q_{i})^{2}&\sum _{i}w_{i}I'_{x}(q_{i})I'_{y}(q_{i})\\[10pt]\sum_{i}w_{i}I'_{x}(q_{i})I'_{y}(q_{i})&\sum _{i}w_{i}I'_{y}(q_{i})^{2} \end{bmatrix}}^{-1}{\begin{bmatrix}-\sum_{i}w_{i}I'_{x}(q_{i})I'_{t}(q_{i})\\[10pt] -\sum_{i}w_{i}I'_{y}(q_{i})I'_{t}(q_{i})\end{bmatrix}} \end{aligned}

实际应用中，因为 Lucas-Kanade 通常会指定一个固定大小的区域，因此其是稀疏的光流计算方法。

1.4 Horn-Schunck 方法

Horn-Schunck 算法假设整个图像上的流是平滑的 (有点奥卡姆剃刀那味了)。因此，它试图尽量减少光流的扭曲，并更喜欢显示更平滑的解，为了更好理解，我们在此将之前的 $\Delta x, \Delta y$ 进行这样的重写 $u=\Delta x, v=\Delta y$ 。

该光流被公式化为全局能量函数，然后寻求最小化。对于二维图像流，其对应的全局能量为 (这里使用拉格朗日乘子法对其进行优化) :

E=\iint\left[(I'_{x}u+I'_{y}v+I'_{t})^{2}+\alpha ^{2}(\lVert \nabla u\rVert ^{2 }+\lVert \nabla v\rVert ^{2})\right]{{\rm {d}}x{\rm {d}}y}=\iint L(x,y) {{\rm {d}}x{\rm {d}}y}

其中 $\overset{\rightarrow}{V}=\left[u(x,y),v(x,y)\right]^T$ 是光流向量，参数 $\alpha$ 用来提升整体的平滑性，对 $u,v$ 求导则有

\begin{aligned} {\frac {\partial L}{\partial u}}-{\frac {\partial }{\partial x}}{\frac {\partial L}{\partial u_{x}}}-{ \frac {\partial }{\partial y}}{\frac {\partial L}{\partial u_{y}}}&=0\\ {\frac {\partial L}{\partial v}}-{\frac {\partial }{\partial x}}{\frac {\partial L}{\partial v_{x}}}-{ \frac {\partial }{\partial y}}{\frac {\partial L}{\partial v_{y}}}&=0 \end{aligned}

同时对于 $x,y$ 求导还能得到，其中 $\Delta ={\frac {\partial ^{2}}{\partial x^{2}}}+{\frac {\partial ^{2}}{\partial y^{2}}}$ 表示拉普拉斯算子:

\begin{aligned} I'_{x}(I'_{x}u+I'_{y}v+I'_{t})-\alpha ^{2}\Delta u&=0\\ I'_{y}(I'_{x}u+I'_{y}v+I'_{t})-\alpha ^{2}\Delta v&=0 \end{aligned}

在实践中，拉普拉斯算子是使用有限差分进行数值近似的，并且可以写成 $\Delta u(x,y)=({\overline {u}}(x,y)-u(x,y))$ ，其中 $\overline{u}(x,y)$ 是 $u$ 在像素 $(x,y)$ 处邻域的加权平均值。使用这个加权估计，可以得到 :

\begin{aligned} ({I'_{x}}^{2}+\alpha ^{2})u+{I'_{x}}{I'_{y}}v&=\alpha ^{2}{\overline {u}}-{I'_{x}}{I'_{t}}\\ {I'_{x}}{I'_{y}}u+({I'_{y}}^{2}+\alpha ^{2})v&=\alpha ^{2}{\overline {v}}-{I'_{y}}{I'_{t}} \end{aligned}

所以按照 Cramer 法则给出以下的迭代解 :

\begin{aligned} u^{k+1}&={\overline{u}}^{k}-{\frac {{I'_{x}}({I'_{x}}{\overline {u}}^{k}+{I'_{y}}{\overline {v}}^{k}+{I'_{t}})}{\alpha ^{2}+{I'_{x}}^{2}+{I'_{y}}^{2}}}\\ v^{k+1}&={\overline{v}}^{k}-{\frac {{I'_{y}}({I'_{y}}{\overline {u}}^{k}+{I'_{y}}{\overline {v}}^{k}+{I'_{t}})}{\alpha ^{2}+{I'_{x}}^{2}+{I'_{y}}^{2}}} \end{aligned}

可以看出，相比 Lucas-Kanade 方法直接给出一个解析解，Horn-Schunck 方法的计算量更大，其同时也是一个稀疏的光流计算方法。

1.5 Farneback 方法

该方法是 opencv 库中给出的稠密的光流计算方法，其使用了多项式展开的思想，即用多项式逼近每个像素的某个邻域。给出的方法中为了降低复杂度使用二次函数模型 :

f(\mathbf{x}) \sim \mathbf{x}^T \mathbf{A} \mathbf{x}+\mathbf{b}^T \mathbf{x}+c

其中 $\mathbf{A}$ 是一个对称矩阵， $\mathbf{b}$ 是向量而 $c$ 是标量。系数是从适合邻域信号值的加权最小二乘估计的。

由于多项式展开的结果是每个邻域都由一个多项式近似，首先分析多项式经历理想的转换会发生什么。考虑精确的二次多项式

f_1(\mathbf{x})=\mathbf{x}^T \mathbf{A}_1 \mathbf{x}+\mathbf{b}_1^T \mathbf{x}+c_1

并通过 $d$ 的全局位移构造一个新的信号 $f_2$ ，

\begin{aligned} f_2(\mathbf{x}) & =f_1(\mathbf{x}-\mathbf{d})=(\mathbf{x}-\mathbf{d})^T \mathbf{A}_1(\mathbf{x}-\mathbf{d})+\mathbf{b}_1^T(\mathbf{x}-\mathbf{d})+c_1 \\ & =\mathbf{x}^T \mathbf{A}_1 \mathbf{x}+\left(\mathbf{b}_1-2 \mathbf{A}_1 \mathbf{d}\right)^T \mathbf{x}+\mathbf{d}^T \mathbf{A}_1 \mathbf{d}-\mathbf{b}_1^T \mathbf{d}+c_1 \\ & =\mathbf{x}^T \mathbf{A}_2 \mathbf{x}+\mathbf{b}_2^T \mathbf{x}+c_2 . \end{aligned}

将二次多项式中的系数相加得到 :

\begin{aligned} \mathbf{A}_2 & =\mathbf{A}_1, \\ \mathbf{b}_2 & =\mathbf{b}_1-2 \mathbf{A}_1 \mathbf{d}, \\ c_2 & =\mathbf{d}^T \mathbf{A}_1 \mathbf{d}-\mathbf{b}_1^T \mathbf{d}+c_1 . \end{aligned}

通过观察关键的 $\mathbf{b}_2 =\mathbf{b}_1-2 \mathbf{A}_1 \mathbf{d}$ ，我们可以求解平移 $d$ ，如果 $A_1$ 是非奇异的，

\begin{aligned} 2 \mathbf{A}_1 \mathbf{d} & =-\left(\mathbf{b}_2-\mathbf{b}_1\right), \\ \mathbf{d} & =-\frac{1}{2} \mathbf{A}_1^{-1}\left(\mathbf{b}_2-\mathbf{b}_1\right) . \end{aligned}

这一观察结果适用于任何信号维度。

1.5.1 Farneback 方法的实际考量

首先，将等式 $f_1(\mathbf{x})=\mathbf{x}^T \mathbf{A}_1 \mathbf{x}+\mathbf{b}_1^T \mathbf{x}+c_1$ 中的全局多项式替换为局部多项式近似。因此，我们首先对两幅图像进行多项式展开，得到第一个图像的展开系数 $\mathbf{A}_1(x)$ 、 $\mathbf{b}_1(x)$ 和 $c_1(x)$ 和第二个图像的 $\mathbf{A}_2(x)$ 、 $\mathbf{b}_2(x)$ 和 $c_2(x)$ 。理想情况下，这应该根据之前的等式给出 $A_1 = A_2$ ，但在实践中必须确定近似值

\mathbf{A}(\mathbf{x})=\frac{\mathbf{A}_1(\mathbf{x})+\mathbf{A}_2(\mathbf{x})}{2}

同时引入 :

\Delta \mathbf{b}(\mathbf{x})=-\frac{1}{2}\left(\mathbf{b}_2(\mathbf{x})-\mathbf{b}_1(\mathbf{x})\right)

以满足主要约束

\mathbf{A}(\mathbf{x}) \mathbf{d}(\mathbf{x})=\Delta \mathbf{b}(\mathbf{x})

1.5.2 邻域估计

至此，这个方程是可以逐像素来求解的，但计算复杂且容易受到噪声影响，因此作者假设位移场只是缓慢变化的。作者尝试在 $x$ 的邻域 $I$ 上找到满足 $\mathbf{A}(\mathbf{x}) \mathbf{d}(\mathbf{x})=\Delta \mathbf{b}(\mathbf{x})$ 的 $\mathbf{d}(\mathbf{x})$ ，或者对下式进行最小化 :

\sum_{\Delta \mathbf{x} \in I} w(\Delta \mathbf{x})\|\mathbf{A}(\mathbf{x}+\Delta \mathbf{x}) \mathbf{d}(\mathbf{x})-\Delta \mathrm{b}(\mathbf{x}+\Delta \mathbf{x})\|^2

可以通过最小二乘法得到最小化的表达 :

\mathbf{d}(\mathbf{x})=\left(\sum w \mathbf{A}^T \mathbf{A}\right)^{-1} \sum w \mathbf{A}^T \boldsymbol{\Delta} \mathbf{b}

剩下部分则是对于光流变化进行的建模以及对于前面部分的补充，在此不过多赘述，可以阅读参考资料中两篇文章 (一篇原文一篇讲解文章，写得还是很不错的)。

下面两张光流图分别是 LK 的稀疏光流方法和 Farneback 的稠密光流方法，可以看出精度以以及目标划分上还是有较大差距的。

2 FlowNet 的论文相关信息

因为是古早时期的 CNN 论文，结构也较为简单也没有太多原理可言，感觉都是纯技巧偏多，直接放结构图和效果吧。(实际是我也没看出什么道理)

图 1.FlowNetSimple 和 FlowNetCorr 的神经网络结构图。

图 2.将粗特征图细化为高分辨率预测。

最终实现的效果如下

3 FlowNet 2 及后续改进

一样，给出 flownet 2 的模型结构图 (还是同一个团队做的，就不放论文信息了)。

图 3.完整架构示意图 : 为了计算大位移光流，作者组合了多个 FlowNet。大括号表示输入的连接。亮度误差是第一幅图像和第二幅图像与先前估计的流之间的差异。为了最佳地处理小位移，作者在 FlowNetS 架构中引入了较小的开始步长和上卷积之间的卷积。最后，作者应用一个小的融合网络来提供最终的估计。

4 思考和总结

总得来说，目前的光流方法还是基于平面图形的假设 (不知道深度学习的方法有没有这种局限性)，从基本的位移考量到后面的仿射和射影，也从一个侧面提供了 scene flow (场景流) 的动机。感觉如果单纯研究平面光流，也很难拿出更好的解决方案，如果在深度学习中以此作为组件，则又存在没有对应的 groundtruth 进行训练的难题。东西是个好东西，但是在学术研究中作用有点有限，不知道工业场景中会不会有更实际的用法。

Optical Flow and FlowNet : 从传统算法到深度模型