光流法介绍1

光流的概念是Gibson在1950年首先提出来的，它是空间运动物体在观察成像平面上的像素运动的瞬时速度，是利用图像序列中像素在时间域上的变化以及相邻帧之间的相关性来找到上一帧跟当前帧之间存在的对应关系，从而计算出相邻帧之间物体的运动信息的一种方法。一般而言，光流是由于场景中前景目标本身的移动、相机的运动，或者两者的共同运动所产生的。
     光流场(Optical Flow Field)法的基本思想：在空间中，运动可以用运动场描述，而在一个图像平面上，物体的运动往往是通过图像序列中不同图像灰度分布的不同体现的，从而，空间中的运动场转移到图像上就表示为光流场。光流场反映了图像上每一点灰度的变化趋势，可看成是带有灰度的像素点在图像平面上运动而产生的瞬时速度场，也是一种对真实运动场的近似估计。光流场是运动场在二维图像的投影，研究光流场也是为了研究运动场。
     IJCV2011有一篇文章，《A Database and Evaluation Methodology for Optical Flow》里面对主流的光流算法做了简要的介绍和对不同算法进行了评估。网址是：vision.middlebury.edu/flow。另外，需要提到的一个问题是，光流场是图片中每个像素都有一个x方向和y方向的位移，所以在上面那些光流计算结束后得到的光流flow是个和原来图像大小相等的双通道图像。那怎么可视化呢？这篇文章用的是Munsell颜色系统来显示。关于孟塞尔颜色系统（MunsellColor System），可以看wikibaike。它是美国艺术家阿尔伯特孟塞尔（Albert H. Munsell，1858－1918）在1898年创制的颜色描述系统。

 
    在上面的那个评估的网站有这个从flow到color显示的Matlab和C++代码。但是感觉C++代码分几个文件，有点乱，然后我自己整理成两个函数了，并配合OpenCV的Mat格式。 
 把C++的一个评估代码提出来共享下：

void makecolorwheel(vector<Scalar> &colorwheel)  

{  

int RY = 15; int YG = 6;  int GC = 4;  int CB = 11;  int BM = 13;  int MR = 6;  

int i;  

for (i = 0; i < RY; i++) colorwheel.push_back(Scalar(255,255*i/RY,0));  

for (i = 0; i < YG; i++) colorwheel.push_back(Scalar(255-255*i/YG,255,0));  

for (i = 0; i < GC; i++) colorwheel.push_back(Scalar(0,255,255*i/GC));  

for (i = 0; i < CB; i++) colorwheel.push_back(Scalar(0,255-255*i/CB, 255));  

for (i = 0; i < BM; i++) colorwheel.push_back(Scalar(255*i/BM,0,255));  

for (i = 0; i < MR; i++) colorwheel.push_back(Scalar(255,0,255-255*i/MR));  

}  

void motionToColor(Mat flow, Mat &color)  

{  

if (color.empty())  color.create(flow.rows,flow.cols,CV_8UC3);  

static vector<Scalar> colorwheel;

if (colorwheel.empty())  makecolorwheel(colorwheel);  

float maxrad = -1;  

for (int i= 0; i < flow.rows; ++i)   

{  

for (int j = 0; j < flow.cols; ++j)   

{  

Vec2f flow_at_point = flow.at<Vec2f>(i, j);  

float fx = flow_at_point[0];  

float fy = flow_at_point[1];  

if ((fabs(fx) >  UNKNOWN_FLOW_THRESH) || (fabs(fy) >  UNKNOWN_FLOW_THRESH))  continue;  

float rad = sqrt(fx * fx + fy * fy);  

maxrad = maxrad > rad ? maxrad : rad;  

}  

}  

for (int i= 0; i < flow.rows; ++i)   

{  

for (int j = 0; j < flow.cols; ++j)   

{  

uchar *data = color.data + color.step[0] * i + color.step[1] * j;  

Vec2f flow_at_point = flow.at<Vec2f>(i, j);  

float fx = flow_at_point[0] / maxrad;  

float fy = flow_at_point[1] / maxrad;  

if ((fabs(fx) >  UNKNOWN_FLOW_THRESH) || (fabs(fy) >  UNKNOWN_FLOW_THRESH))  

{  

data[0] = data[1] = data[2] = 0;  

continue;  

}  

float rad = sqrt(fx * fx + fy * fy);  

float angle = atan2(-fy, -fx) / CV_PI;  

float fk = (angle + 1.0) / 2.0 * (colorwheel.size()-1);  

int k0 = (int)fk;  

int k1 = (k0 + 1) % colorwheel.size();  

float f = fk - k0; 

for (int b = 0; b < 3; b++)   

{  

float col0 = colorwheel[k0][b] / 255.0;  

float col1 = colorwheel[k1][b] / 255.0;  

float col = (1 - f) * col0 + f * col1;  

if (rad <= 1)  col = 1 - rad * (1 - col);

else           col *= .75;

data[2 - b] = (int)(255.0 * col);  

}  

}  

}  

} 

1. Horn_Schunck方法

1981年，Horn和Schunck创造性地将二维速度场与灰度相联系，引入光流约束方程，得到光流计算的基本算法。 此方法是首次使用亮度恒定假设和推导出基本的亮度恒定方程的方法之一。用u(x, y,t)表示视频序列，(x(t), y(t))表示一个点在像平面的轨迹，那么亮度不变假设意味着u(x(t), y(t),t)受下式约束：

对于这个欠定的方程，Horn和Schunck求解该等式假设一个平滑性约束，通过对光流速度分量的二阶导数进行规则化获得。
 opencv提供了函数CalcOpticalFlowHS计算Horn_Schunck(稠密光流)
2.  Lucas-Kanade方法   
参考这篇论文：”Pyramidal Implementation of the Lucas Kanade Feature TrackerDescription of the algorithm” 
 通过金字塔Lucas-Kanade 光流方法计算某些点集的光流（稀疏光流）。
Lucas-Kanade方法有三个假设：
  （1）亮度恒定，就是同一点随着时间的变化，其亮度不会发生改变。这是基本光流法的假定（所有光流法变种都必须满足），用于得到光流法基本方程；

  （2）小运动，这个也必须满足，就是时间的变化不会引起位置的剧烈变化，这样灰度才能对位置求偏导（换句话说，小运动情况下我们才能用前后帧之间单位位置变化引起的灰度变化去近似灰度对位置的偏导数），这也是光流法不可或缺的假定；

  （3）空间一致，一个场景上邻近的点投影到图像上也是邻近点，且邻近点速度一致。这是Lucas-Kanade光流法特有的假定，因为光流法基本方程约束只有一个，而要求x，y方向的速度，有两个未知变量。我们假定特征点邻域内做相似运动，就可以连立n多个方程求取x，y方向的速度（n为特征点邻域总点数，包括该特征点）。
设：速度的y分量为v，x分量为u，则光流约束方程为：

 
对于这个欠定方程： Lucas-Kanade提出：
 

 

 

opencv提供了该函数：calcOpticalFlowPyrLK
3.Gunnar Farneback 的算法
 相关论文是："Two-Frame Motion Estimation Based on PolynomialExpansion"
opencv的函数calcOpticalFlowFarneback 

 4.块匹配的方法
块匹配算法对像素集合进行处理而非单个像素。
opencv的函数：CalcOpticalFlowBM 
5.SimpleFlow 

 opencv的函数：calcOpticalFlowSF
    这一个是2012年欧洲视觉会议的一篇文章的实现："SimpleFlow:A Non-iterative,Sublinear Optical FlowAlgorithm"，网站是：SimpleFlow: A Non-iterative, Sublinear Optical Flow Algorithm - U.C. Berkeley Computer Graphics Research  

上面介绍了几种光流算法，这里给出Gunnar Farneback 的算法的测试：

图像：
 

测试代码： 

void testImage()

{

Mat prevgray, gray, flow, cflow, frame1,frame2;Mat motion2color;

frame1=imread("I0.png");

frame2=imread("I1.png");

cvtColor(frame1,prevgray,CV_BGR2GRAY);

cvtColor(frame2,gray,CV_BGR2GRAY);

calcOpticalFlowFarneback(prevgray, gray,flow,0.5,3,15,3,5,1.2,0);

motionToColor(flow, motion2color);  

imshow("flow", motion2color);  

waitKey(0);

}

 
 视频测试：
 void testVideo()

{

Mat frame;Mat gray, prevGray,flow; Mat motion2color;

VideoCapture cap;  

cap.open("video1.avi");

if (cap.isOpened())  

{

while(true)

{

cap>>frame;

if (frame.empty())  break;

imshow("orl",frame);

cvtColor(frame, gray, COLOR_BGR2GRAY);//转成灰度图像

if (prevGray.empty()) gray.copyTo(prevGray);

//用Gunnar Farneback 的算法计算稠密光流

calcOpticalFlowFarneback(prevGray, gray, flow,0.5,3,15,3,5,1.2,0);

swap(prevGray, gray);

motionToColor(flow, motion2color);  

imshow("flow", motion2color); 

int c = waitKey(100);

if( (char)c == 27 ) break; 

}

}

}

 

摄像头测试：

void testVideo2()

{

Mat gray, prevGray,flow; 

Mat motion2color;

CvCapture* capture=cvCreateCameraCapture(0);

//设置分辨率

cvSetCaptureProperty(capture,CV_CAP_PROP_FRAME_WIDTH ,240);

cvSetCaptureProperty(capture,CV_CAP_PROP_FRAME_HEIGHT,160);

Mat frame;

while(capture)

{

frame = cvQueryFrame(capture);//获得一帧图象

if (frame.empty())  break;

imshow("orl",frame);

cvtColor(frame, gray, COLOR_BGR2GRAY);//转成灰度图像

if (prevGray.empty()) gray.copyTo(prevGray);

//用Gunnar Farneback 的算法计算稠密光流

calcOpticalFlowFarneback(prevGray,gray,flow,0.5,3,15,3,5,1.2,0);

swap(prevGray, gray);

motionToColor(flow, motion2color);  

imshow("flow", motion2color); 

if(cvWaitKey(1)==27) break;

}

cvDestroyAllWindows(); 

}