稠密光流----Farneback

作者lowkeyway​

我们在洋洋洒洒介绍过KLT时，有提到：

在光流方程中，有两个未知数，一个方程，两个未知数，初中生都知道是没办法求解的。但是这可难不倒Bruce D. Lucas 和Takeo Kanade，他们提出来稀疏光流跟踪算法。

其实有个概念隐含在原理中，为什么称KLT为稀疏光流跟踪算法？

这一切都可以从稀疏矩阵说起。

什么是稀疏矩阵？

在矩阵中，若数值为0的元素数目远远多于非0元素的数目，并且非0元素分布没有规律时，则称该矩阵为稀疏矩阵；与之相反，若非0元素数目占大多数时，则称该矩阵为稠密矩阵。

看不懂定义没关系，看看图再反过来看就好理解了：

看我的0多不多

稀疏矩阵有什么好的呢？当然是存储可压缩、计算方便简单啦！

那为什么称KLT为稀疏矩阵呢？

还记得我们运用KLT的时候必须要先计算角点吗？如果我把角点位置当做1，非角点位置当做0，那么一张图像是不是就像极了稀疏矩阵？

KLT的光流图

黑点就是对应的角点位置

你就说：稀疏不稀疏？

补充：

LK稀疏算法只需要每个感兴趣点周围小窗口的局部信息，这样就有一个问题：

• 较大的运动会将点移除这个小窗口，从而造成算法无法再找到这些点。（所以才有了KLT的三个前提假设）

1. 相邻帧之间的亮度恒定；
2. 相邻帧之间物体的运动比较“微小”；
3. 保持空间一致性；即，相邻像素点具有相同的运动

不过，金字塔的LK算法可以改善这个问题，即从金字塔的最高层（细节最少）开始向金字塔的最低层（丰富的细节）进行跟踪。跟踪图像金字塔允许小窗口部或较大的运动。

还记得API接口calcOpticalFlowPyrLK中的maxLevel参数吗？

难道还有稠密光流算法？

没有金刚钻别揽瓷器活，既然有人提出稀疏光流，就有人来解决稠密光流。

OpenCV中支持的稠密光流算法是由Gunner Farneback在2003年提出来的，它是基于前后两帧所有像素点的移动估算算法，其效果要比稀疏光流算法更好。

从逻辑上讲，Gunnar Farneback算法计算全局性的稠密光流算法（即图像上所有像素点的光流都计算出来），由于要计算图像上所有点的光流，故计算耗时，速度慢。需要使用某种插值方法在比较容易跟踪的像素之间进行插值以解决那些运动不明确的像素。

宏观上理解到时很容易，稀疏的关键点是计算角点，稠密的关键点是找图像上所有移动的点。可是，从原理上能够自洽吗？准确的问，有可计算性吗？

当然，这是算法专家的灵魂拷问，关我们工程狗什么事？可以去看论文“Two-Frame Motion Estimation Based on Polynomial Expansion”、可以去参考：

​​http://vision.middlebury.edu/flow/​​

利用多项式对每个像素的邻域信息进行近似表示。

可是。。。

Farneback算法原理

非常简单，只有两步。

1. 将图像视为二维信号的函数（输出图像是灰度图像），因变量是二维坐标位置  ，并且利用二次多项式对于图像进行近似建模的话，会得到:

其中，

• A 是一个2×2的对称矩阵（是通过像素的邻域信息的最小二乘加权拟合得到的，权重系数与邻域的像素大小和位置有关）
• b是一个2×1的矩阵向量；
• c为标量；

因此系数化之后，以上公式等号右侧可以写为（看不懂系列）：

2. 如果将原有（笛卡尔坐标系）图像的二维信号空间，转换到以 (  )作为基函数的空间，则表示图像需要一个六维向量作为系数，代入不同像素点的位置 x,y求出不同像素点的灰度值。

至于A中的最小二乘法的、权重的选择

OpenCV API

老规矩，OpenCV已经给我们实现了这一切，下面我们只要明白原理（而不是推导奥），熟知参数含义，就可以用起来了

CV_EXPORTS_W void calcOpticalFlowFarneback( InputArray prev, InputArray next, InputOutputArray flow,
                                            double pyr_scale, int levels, int winsize,
                                            int iterations, int poly_n, double poly_sigma,
                                            int flags );

• prev：前一帧图像
• next：后一帧图像
• flow：输出的光流矩阵。矩阵大小同输入的图像一样大，但是矩阵中的每一个元素可不是一个值，而是两个值，分别表示这个点在x方向与y方向的运动量（偏移量）。
• pyr_scale：金字塔上下两层之间的尺度关系
• levels：金字塔层数
• winsize：均值窗口大小，越大越能denoise并且能够检测快速移动目标，但会引起模糊运动区域
• iterations：迭代次数
• poly_n：像素领域大小，一般为5，7等
• poly_sigma：高斯标注差，一般为1-1.5
• flags：计算方法

C++

以下C++代码我修改自官方demo。

注意点：1. 最重要的一点，这里引入了UMat的概念，简单的说，UMat是给OpenCL定制的GPU处理定制的，OpenCV3.0之后可以根据底层是否支持OpenCL来决定图像处理是用CPU还是GPU。

2. calcOpticalFlowFarneback运行的结果flow是包含x/y两个方向的偏移量，需要分别提取出来。
3. 本例子中虽然FB处理了全图，但是绘制的时候，是根据step的间隔决定绘制那些关键点的。

#include <iostream>#include <opencv2/opencv.hpp>
using namespace cv;
using namespace std;

static void drawOptFlowMap(const Mat& flow, Mat& cflowmap, int step,
double, const Scalar& color)
{
for (int y = 0; y < cflowmap.rows; y += step)
for (int x = 0; x < cflowmap.cols; x += step)
    {
const Point2f& fxy = flow.at<Point2f>(y, x);
      line(cflowmap, Point(x, y), Point(cvRound(x + fxy.x), cvRound(y + fxy.y)),
        color);
      circle(cflowmap, Point(x, y), 2, color, -1);
    }
}

int main(int argc, char** argv)
{
//VideoCapture cap(0);
  VideoCapture cap;
  cap.open("vtest.avi");
if (!cap.isOpened())
return -1;

  Mat flow, cflow, frame, preframe;
//Mat gray, prevgray, uflow;
  UMat gray, prevgray, uflow;
  namedWindow("flow", 1);

for (;;)
  {
//cap >> frame;
    bool ret = cap.read(frame);
    cvtColor(frame, gray, COLOR_BGR2GRAY);

if (!prevgray.empty())
    {
      calcOpticalFlowFarneback(prevgray, gray, uflow, 0.5, 3, 15, 3, 5, 1.2, 0);
//cvtColor(prevgray, cflow, COLOR_GRAY2BGR);
      uflow.copyTo(flow);
      drawOptFlowMap(flow, preframe, 16, 1.5, Scalar(0, 255, 0));
      imshow("flow", preframe);
    }
if (waitKey(30) >= 0)
break;
    std::swap(prevgray, gray);
    std::swap(preframe, frame);
  }
return 0;
}

运行结果

Python:

改自互联网。

注意事项：
这里的思路比较奇特，他把所有的位置转换成了极坐标，然后通过极坐标中的半径和角度对应到HSV上，通过HSV再转换成RGB显示出来。这么做的好处是非常直观，而且炫彩！

#!/usr/bin/env python# -*- coding: utf-8 -*-'''# author：lowkeyway time:11/2/2019# This file creat for test calcOpticalFlowFarneback API'''import sysimport cv2 as cvimport numpy as npfb_params = dict(pyr_scale = 0.5,
                 levels = 3,
                 winsize = 15,
                 iterations = 3,
                 poly_n = 5,
                 poly_sigma = 1.2,
                 flags = 0)class app:
    def __init__(self, src):
        self.cap = cv.VideoCapture(src)

    def run(self):
        ret, preFrame = self.cap.read()
        if ret is not True:
            return
        preFrameGary = cv.cvtColor(preFrame, cv.COLOR_BGR2GRAY)
        hsv = np.zeros_like(preFrame)
        hsv[..., 1] = 255

        while True:
            ret, curFrame = self.cap.read()
            if ret is not True:
                break
            curFrameGary = cv.cvtColor(curFrame, cv.COLOR_BGR2GRAY)
            flow = cv.calcOpticalFlowFarneback(preFrameGary,curFrameGary, None, **fb_params)
            mag, ang = cv.cartToPolar(flow[..., 0], flow[..., 1])
            hsv[..., 0] = ang * 180 / np.pi / 2
            hsv[..., 2] = cv.normalize(mag, None, 0, 255, cv.NORM_MINMAX)

            bgr = cv.cvtColor(hsv, cv.COLOR_HSV2BGR)
            cv.imshow("Frame", bgr)

            ch = cv.waitKey(30) & 0xff
            if ch == 27:
                break

            preFrameGary = curFrameGary

        self.cap.release()def main_func(argv):
    try:
        videoSrc = sys.argv[1]
    except:
        videoSrc = "vtest.avi"

    app(videoSrc).run()if __name__ == '__main__':
    print(__doc__)
    main_func(sys.argv)

运行结果

上述内容，如有侵犯版权，请联系作者，会自行删文。