David Stavens的光流法检测代码理解

1、为了进一步了解光流法的性质和使用方法，学习了Stanford Artificial Intelligence Lab的David Stavens的课件，从一个规范的角度去了解Optical Flow在OPENCV中的使用，并调试及运行了其提供的这一版代码。环境是VS2013+OPENCV2.4.11。此代码的作用是稀疏光流的检测，并用箭头指明方向。

 

2、接下来对David Stavens这一版OpenCV_OpticalFlow代码进行逐句分析，并查阅相关材料添加注释。

 

#include<stdio.h>#include<cv.h>#include<highgui.h>#include<math.h>static const doublepi = 3.14159265358979323846; inline static doublesquare(int a){         return a * a;} /* 该函数为图像img分配空间大小size，并设定格式位深depth和通道数channels */inline static voidallocateOnDemand( IplImage **img, CvSize size, int depth, int channels){         if ( *img != NULL ) return;                   *img = cvCreateImage( size,depth, channels ); // 创建首地址并分配存储空间         if ( *img == NULL )         {                   fprintf(stderr, "Error:Couldn't allocate image. Out of memory?\n");                   exit(-1);         }} int main(void){         CvCapture *input_video =cvCaptureFromFile("video1.avi"); /* CvCapture结构体，用来保存图像捕获的信息，cvCaptureFromFile(argv[1])从文件中打开一个视频，另一个是cvCaptureFromCAM( int index )则是从摄像头捕获视频，这里打开video1.avi */          if (input_video == NULL)         {         fprintf(stderr, "Error: Can't openvideo.\n");         return -1;         }         cvQueryFrame( input_video ); // 从摄像头或者文件中抓取并返回一帧         CvSize frame_size; /* CvSize，OpenCV的基本数据类型之一，表示矩阵框大小，以像素为精度,数据成员是integer类型的width和height */          frame_size.height =  (int) cvGetCaptureProperty( input_video,CV_CAP_PROP_FRAME_HEIGHT ); /* doublecvGetCaptureProperty( CvCapture* capture, int property_id )获取视频流的各种属性,capture为视频流，property_id 为属性名，CV_CAP_PROP_FRAME_HEIGHT视频流中的帧高度, CV_CAP_PROP_FRAME_HEIGHT视频流中的帧宽度 */         frame_size.width =   (int) cvGetCaptureProperty( input_video,CV_CAP_PROP_FRAME_WIDTH );          long number_of_frames; // 第几帧         cvSetCaptureProperty( input_video,CV_CAP_PROP_POS_AVI_RATIO, 1. );/* intcvSetCaptureProperty( CvCapture* capture, int property_id, double value )获取视频流的各种属性,当property_id 为CV_CAP_PROP_POS_MSEC从文件开始的位置（单位为毫秒）、CV_CAP_PROP_POS_FRAMES单位为帧数的位置（只对视频文件有效）、CV_CAP_PROP_POS_AVI_RATIO视频文件的相对位置（0 - 影片的开始，1 - 影片的结尾)时可设置value值，此时value为1表示到视频文件的结尾 */          number_of_frames = (int) cvGetCaptureProperty(input_video, CV_CAP_PROP_POS_FRAMES );// 获取结尾帧位置         cvSetCaptureProperty( input_video,CV_CAP_PROP_POS_FRAMES, 0. ); // 视频到开头         cvNamedWindow("Optical Flow",CV_WINDOW_AUTOSIZE);// 打开一个自适应大小的窗口Optical Flow         long current_frame = 0; // 当前帧设置为0          while(true) //第一层循环         {                   static IplImage *frame =NULL, *frame1 = NULL, *frame1_1C = NULL, *frame2_1C =                   NULL, *eig_image = NULL,*temp_image = NULL, *pyramid1 = NULL, *pyramid2 = NULL;/* IplImage结构体用来声明图像，eig_image、temp_image、pyramid1、pyramid2为临时空间 */                    cvSetCaptureProperty(input_video, CV_CAP_PROP_POS_FRAMES, current_frame );// 将视频设置到指定的当前帧                   frame = cvQueryFrame(input_video );                   if (frame == NULL)                   {                            fprintf(stderr,"Error: Hmm. The end came sooner than we thought.\n");                            return -1;                   }                   allocateOnDemand(&frame1_1C, frame_size, IPL_DEPTH_8U, 1 );/* 为frame1_1C分配空间，IPL_DEPTH_8U意思是8位无符号整型，即2的8次方256色图像深度，通道数为1。通道数：单通道为灰度，3通道为RGB，4通道为RGB+透明度 */                   cvConvertImage(frame,frame1_1C, CV_CVTIMG_FLIP); // 翻转freme到frame1_1C？                   allocateOnDemand(&frame1, frame_size, IPL_DEPTH_8U, 3 ); // 3通道                   cvConvertImage(frame, frame1,CV_CVTIMG_FLIP); // 翻转freme到frame1？                   frame = cvQueryFrame(input_video );                   if (frame == NULL)                   {                            fprintf(stderr,"Error: Hmm. The end came sooner than we thought.\n");                            return -1;                   }                   allocateOnDemand(&frame2_1C, frame_size, IPL_DEPTH_8U, 1 );                   cvConvertImage(frame,frame2_1C, CV_CVTIMG_FLIP);                   allocateOnDemand(&eig_image, frame_size, IPL_DEPTH_32F, 1 ); //32bits float型灰度图像                   allocateOnDemand( &temp_image,frame_size, IPL_DEPTH_32F, 1 );                   CvPoint2D32fframe1_features[400]; /* frame1_features为frame_1的特征点，CvPoint2D32f 表示32bits float型二维坐标下的点*/                   int number_of_features; // 特征点数量                   number_of_features = 400;                    cvGoodFeaturesToTrack(frame1_1C,eig_image, temp_image, frame1_features, &number_of_features, .01, .01,NULL);/* void cvGoodFeaturesToTrack(const CvArr* image, CvArr* eig_image, CvArr* temp_image,CvPoint2D32f* corners, int* corner_count, doublequality_level, double min_distance,const CvArr* mask=NULL,int block_size =NULL, int use_harris = 0, double k = 0.4);检测角点作为特征点，输出到frame1_features和number_of_featuresimage输入图像，8-位或浮点32-比特，单通道eig_image临时浮点32-位图像，尺寸与输入图像一致temp_image另外一个临时图像，格式与尺寸与 eig_image 一致corners输出参数，检测到的角点corner_count输出参数，检测到的角点数目quality_level最大最小特征值的乘法因子。定义可接受图像角点的最小质量因子。min_distance限制因子。得到的角点的最小距离。使用 Euclidian 距离maskROI:感兴趣区域。函数在ROI中计算角点，如果 mask 为 NULL，则选择整个图像。 */                    CvPoint2D32fframe2_features[400]; // frame2中和frame1对应的特征点                   charoptical_flow_found_feature[400]; //frame2中发现的frame1中的特征点                   floatoptical_flow_feature_error[400];                   CvSize optical_flow_window =cvSize(3,3); //避免孔径问题                   CvTermCriteria optical_flow_termination_criteria= cvTermCriteria( CV_TERMCRIT_ITER | CV_TERMCRIT_EPS, 20, .3 ); /* 停止条件，当迭代次数达到20次或结果精确性ε优于0.3时停止【1】CV_TERMCRIT_ITER---在完成最大的迭代次数之后,停止算法  【2】CV_TERMCRIT_EPS----当算法的精确度小于参数double epsilon指定的精确度时，停止算法  【3】CV_TERMCRIT_ITER+CV_TERMCRIT_EPS--无论是哪一个条件先达到,都停止算法 */                    allocateOnDemand(&pyramid1, frame_size, IPL_DEPTH_8U, 1 ); // 临时空间                   allocateOnDemand(&pyramid2, frame_size, IPL_DEPTH_8U, 1 ); /*金字塔Lucas-kanade检测光流frame1_1C：第一帧frame2_1C：第二帧pyramid1：第一帧临时缓存pyramid2：第二帧临时缓存frame1_features：第一帧中指定的特征点集frame2_features：第二帧中计算出来的特征点集number_of_features：特征点的数量optical_flow_window：每个金字塔层的搜索窗口尺寸5：最大的金字塔层数，这里为6层，0的话为1层optical_flow_found_feature：如果第二帧中对应位置的特征点被发现，则该元素被设置为1，否则设置为0optical_flow_feature_error：第二帧中计算出的特征点和第一帧中对应特征点的差值optical_flow_termination_criteria：在每个金字塔层，为某点寻找光流的迭代过程的终止条件0 */                  cvCalcOpticalFlowPyrLK(frame1_1C,frame2_1C, pyramid1, pyramid2, frame1_features,                   frame2_features,number_of_features, optical_flow_window, 5,                   optical_flow_found_feature,optical_flow_feature_error,                   optical_flow_termination_criteria,0 );                    for(int i = 0; i <number_of_features; i++) // 第二层循环，对于每个特征点                   {                            if (optical_flow_found_feature[i] == 0 ) continue; // 没找到则下一个点                            int line_thickness;line_thickness = 1; // 画线的粗度为1                            CvScalar line_color;line_color = CV_RGB(255,0,0); // 线的颜色为RGB（255,0,0）红色                            CvPoint p,q;                            p.x = (int)frame1_features[i].x; // 第一帧中的特征点                            p.y = (int)frame1_features[i].y;                            q.x = (int)frame2_features[i].x; // 第二帧中找到的对应特征点                            q.y = (int)frame2_features[i].y;                            double angle; angle= atan2( (double) p.y - q.y, (double) p.x - q.x );// 计算两点连线的角度                           double hypotenuse; hypotenuse = sqrt( square(p.y- q.y) + square(p.x - q.x) );// 计算斜边长度，即两点连线的长度                            q.x = (int) (p.x - 3* hypotenuse * cos(angle)); // 将连线延长3倍，便于看清                            q.y = (int) (p.y - 3* hypotenuse * sin(angle));                            cvLine( frame1, p,q, line_color, line_thickness, CV_AA, 0 ); // 画线                            p.x = (int) (q.x + 9* cos(angle + pi / 4));                            p.y = (int) (q.y + 9* sin(angle + pi / 4));                            cvLine( frame1, p,q, line_color, line_thickness, CV_AA, 0 ); // 画箭头的一翼                            p.x = (int) (q.x + 9* cos(angle - pi / 4));                            p.y = (int) (q.y + 9* sin(angle - pi / 4));                            cvLine( frame1, p,q, line_color, line_thickness, CV_AA, 0 ); // 画箭头的另一翼                   }                   cvShowImage("OpticalFlow", frame1); // 将画好的图像显示在窗口Optical Flow中                  // 等待按键控制，b或B为上一帧，其他为下一帧int key_pressed;                   key_pressed = cvWaitKey(0);                   if (key_pressed == 'b' ||key_pressed == 'B') current_frame--;                   else current_frame++;// 防止视频帧溢出                  if (current_frame< 0) current_frame = 0;                   if (current_frame >=number_of_frames - 1) current_frame = number_of_frames - 2;         }}

 

3、运行结果

如下图所示，各帧都可以较好地计算出视频中的光流变化，表现在具体事物上来看就是道路上的车辆往右行驶的时候，计算出其上的特征点的对应关系，箭头反映出，车辆是如箭头所指方向向右行驶的。