本文内容原是写毕业论文相关工作时整理的内容，但发现后面并不需要在这里进行详细公式说明，又不忍直接删掉，就让它们换个形式保存在这里吧。
⚠️ 注意：本文写于 2022 年，内容较为基础。
光流法基本假设
光流法具有两个基本假设：
-
相邻帧的亮度恒定不变，即物体在视频中运动时，相邻帧的亮度不发生明显变化；
-
时间连续或微小运动，即相邻帧之间物体的位置变化并不剧烈。
基于这两个假设，我们可以得到图像的约束方程为：
\[I(x,y,t)=I(x+\delta x,y+\delta y,t+\delta t) \tag{1}\]其中 $I(x,y,t)$ 是 $t$ 时刻图像在 $(x,y)$ 像素位置的亮度。在经过 $\delta t$ 时间后，像素在新一帧中的位置相对原位置移动了 $(\delta x,\delta y)$。
将公式 (1) 右端在 $(x,y,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 + \epsilon \tag{2}\]其中 $\epsilon$ 表示泰勒公式的高阶无穷小项，可以近似为 0。
将公式 (1) 和 (2) 联立，得到：
\[\frac{\partial I}{\partial x}\delta x + \frac{\partial I}{\partial y}\delta y + \frac{\partial I}{\partial t}\delta t = 0\]两边同除以 $\delta t$ 得到：
\[\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 \tag{3}\]其中 $\frac{\delta x}{\delta t}$ 和 $\frac{\delta y}{\delta t}$ 分别表示像素点速度矢量沿着 $x$ 方向和 $y$ 方向的导数，分别记为 $u,v$。再令 $I_x=\frac{\partial I}{\partial x}$，$I_y=\frac{\partial I}{\partial y}$，$I_t=\frac{\partial I}{\partial t}$ 分别表示像素点的灰度沿 $x,y,t$ 方向的偏导数。
这时，公式 (3) 可以写成：
\[I_x u + I_y v + I_t = 0\]其中 $I_x, I_y, I_t$ 都可以从图像数据中获得，$(u, v)$ 为未知的光流矢量。
在上述结果中，有一个约束方程，但却有两个未知数，这是不可解的，需要引入额外的约束条件。根据约束条件的不同，光流法产生了不同的分支，其中又可以分为两种实现方案：稀疏光流法和稠密光流法。
稀疏光流法
稀疏光流法指的是先从图像中选取一些关键点（一般选择角点），在之后的识别和追踪中，只关注这些关键点的像素运动，从而大大减小了计算量。
稀疏光流法的代表是 Lucas Kanade 算法 (LK 算法)，它在光流法基本假设的基础上提出了另一个假设：
- 空间一致性，场景中同一物体的相邻像素点具有相似的运动，且在投影到二维图像平面上的距离也比较近。
据此，我们假设在一个大小为 $m\times m$ $(n = m^2)$ 的窗口内，图像的光流是一个恒定值，那么可以得到如下方程组：
\[\begin{aligned} I_{x_1}u + I_{y_1}v & = -I_{t_1} \\ I_{x_2}u + I_{y_2}v & = -I_{t_2} \\ & \cdots \\ & \cdots \\ I_{x_n}u + I_{y_n}v & = -I_{t_n} \end{aligned} \tag{4}\]上述方程组的矩阵形式为：
\[\left[ \begin{array}{cc} I_{x_1} & I_{y_1} \\ I_{x_2} & I_{y_2} \\ \cdot & \cdot \\ \cdot & \cdot \\ I_{x_n} & I_{y_n} \end{array} \right] \left[ \begin{array}{c} u \\ v \end{array} \right] = \left[ \begin{array}{c} -I_{t_1} \\ -I_{t_2} \\ \cdot \\ \cdot \\ -I_{t_n} \end{array} \right]\]记为 $\mathbf{A}\mathbf{v} = \mathbf{b}$
使用最小二乘法得到：
\[\mathbf{A}^{\mathsf{T}}\mathbf{A}\mathbf{v} = \mathbf{A}^{\mathsf{T}}\mathbf{b}\]解得速度矢量为：
\[\mathbf{v} = (\mathbf{A}^{\mathsf{T}}\mathbf{A})^{-1}\mathbf{A}^{\mathsf{T}}\mathbf{b} \tag{5}\]该速度矢量就被视为光流。
稀疏光流法的缺点也十分明显，因为只追踪根据之前帧得到的关键点，所以往往会忽视新出现在图像中的物体。由于这一缺点，本文并不使用 LK 算法作为运动检测的算法，而是对其进行简化，输出两帧图像中像素运动的最大距离，以此衡量物体运动速度，并影响视频参数的配置。
稠密光流法
稠密光流法并没有像稀疏光流那样取巧，而是直接计算一帧中全部像素的光流，因此稠密光流更加精确，可以表现物体整体的运动情况，但同时也意味着对计算资源的巨大需求，速度更慢。
稠密光流法的代表是由 Gunner Farneback 提出的 Farnback 算法，该算法的假设与 LK 算法的相同。
Farnback 算法将输入图像灰度化，并对图像进行二次多项式建模，对于每一个像素位置 $\mathbf{x} = (x\ y)^{\mathsf{T}}$，其灰度值是一个关于 $\mathbf{x}$ 的函数 $f(\mathbf{x})$。在以待求像素点为中心的局部坐标系中，对函数进行二项展开，可以近似为：
\[\begin{aligned} f(\mathbf{x}) & = f(x, y) \\ & \approx r_1 + r_2 x + r_3 y + r_4 x^2 + r_5 y^2 + r_6 xy \\ & = (x\ y)^{\mathsf{T}} \begin{pmatrix} r_4 & r_6/2 \\ r_6/2 & r_5 \end{pmatrix} \begin{pmatrix} x \\ y \end{pmatrix} + (r_2\ r_3) \begin{pmatrix} x \\ y \end{pmatrix} + r_1 \\ & = \mathbf{x}^{\mathsf{T}}\mathbf{A}\mathbf{x} + \mathbf{b}^{\mathsf{T}}\mathbf{x} + c \end{aligned} \tag{6}\]其中 $r_1 \sim r_6$ 是只针对于这一个像素点的参数，可以在该像素点的邻域内使用加权最小二乘法进行计算，距离该像素点越近的点赋予更大的权重。
设像素初始位置的多项式为 $f_1(\mathbf{x})=\mathbf{x}^{\mathsf{T}}\mathbf{A}_1\mathbf{x} + \mathbf{b}_1^{\mathsf{T}}\mathbf{x}+c_1$，根据上面提到的假设 2，像素发生微小移动 $\mathbf{d}$ 后，得到新的多项式：
\[\begin{aligned} f_2(\mathbf{x}) & = f_1(\mathbf{x} - \mathbf{d}) \\ & = (\mathbf{x} - \mathbf{d})^{\mathsf{T}}\mathbf{A}_1(\mathbf{x}-\mathbf{d}) + \mathbf{b}_1^{\mathsf{T}}(\mathbf{x}-\mathbf{d}) + c_1 \\ & = \mathbf{x}^{\mathsf{T}}\mathbf{A}_1\mathbf{x}+(\mathbf{b}_1-2\mathbf{A}_1\mathbf{d})^{\mathsf{T}}\mathbf{x}+\mathbf{d}^{\mathsf{T}}\mathbf{A}_1\mathbf{d} - \mathbf{b}_1^{\mathsf{T}}\mathbf{d} + c_1 \\ & = \mathbf{x}^{\mathsf{T}}\mathbf{A}_2\mathbf{x}+\mathbf{b}_2^{\mathsf{T}}\mathbf{x}+c_2 \end{aligned} \tag{7}\]其中满足
\[\begin{aligned} \mathbf{A}_2 & = \mathbf{A}_1 \\ \mathbf{b}_2 & = \mathbf{b}_1 - 2\mathbf{A}_1\mathbf{d} \\ c_2 & = \mathbf{d}^{\mathsf{T}}\mathbf{A}_1\mathbf{d}-\mathbf{b}_1^{\mathsf{T}}\mathbf{d} + c_1 \end{aligned} \tag{8}\]如果 $\mathbf{A}_1$ 非奇异，则可以得到 $\mathbf{d}=-\frac{1}{2}\mathbf{A}_1^{-1}(\mathbf{b}_2-\mathbf{b}_1)$
虽然我们可以从数学推导中得到 $\mathbf{A}_2=\mathbf{A}_1$，但在实际中该式未必满足。使用局部多项式对理论值进行逼近，然后通过平均值来近似，设第一帧图像对应的系数为 $\mathbf{A}_1(\mathbf{x}),\mathbf{b}_1(\mathbf{x}),c_1(\mathbf{x})$，第二帧图像同理，则
\[\begin{aligned} \mathbf{A}(\mathbf{x}) & = \frac{\mathbf{A}_1(\mathbf{x})+\mathbf{A}_2(\mathbf{x})}{2} \\ \Delta\mathbf{b}(\mathbf{x}) & = -\frac{1}{2}(\mathbf{b}_2(\mathbf{x}) - \mathbf{b}_1(\mathbf{x})) \end{aligned} \tag{9}\]因此
\[\mathbf{A}(\mathbf{x})\mathbf{d}(\mathbf{x}) = \Delta\mathbf{b}(\mathbf{x}) \tag{10}\]使用公式 (10) 可以对一帧图像的每个像素点进行运算，但这样是不现实的，可以尽可能地缩小 $\mathbf{x}$ 的邻域 $I$ 的范围，并找到符合的 $\mathbf{d}(\mathbf{x})$：
\[\sum_{\Delta \mathbf{x}\in I} \omega(\Delta \mathbf{x}) \left\| \mathbf{A}(\mathbf{x} + \Delta \mathbf{x})\mathbf{d}(\mathbf{x}) - \Delta \mathbf{b}(\mathbf{x} + \Delta \mathbf{x}) \right\|^{2} \tag{11}\]其中 $\omega(\Delta \mathbf{x})$ 是像素点对应的权重函数，使用最小二乘法可以求得
\[\mathbf{d}(\mathbf{x}) = \left(\sum \omega\mathbf{A}^{\mathsf{T}}\mathbf{A}\right)^{-1} \sum \omega \mathbf{A}^{\mathsf{T}}\Delta \mathbf{b} \tag{12}\]