MedianFlow代码 & 程序示例

 TLD开源代码MedianFlow部分，改造成一个运行程序示例，包含测试主程序及测试视频序列。运行环境：vs2008 + opencv2.4.2

代码及测试视频序列下载地址：https://download.csdn.net/download/qq_28584889/12012477

//头文件包含：
    MFSystemStruct.h
    MedianFlow.h
    OpticalFlow.h

//cpp文件包含
    MedianFlow.cpp
    OpticalFlow.cpp
    main.cpp

//头文件：MFSystemStruct.h

#ifndef MedianFlow_systemStruct_h
#define MedianFlow_systemStruct_h

#include <opencv2/opencv.hpp>
#include <algorithm>
#include <vector>

using namespace cv;
using namespace std;

typedef float TYPE_OF_COORD;
typedef TYPE_OF_COORD TYPE_MF_COORD;

typedef Point_<TYPE_OF_COORD> TYPE_OF_PT;
typedef TYPE_OF_PT TYPE_MF_PT;
typedef Rect_<TYPE_MF_COORD> TYPE_MF_BB;

typedef pair<Mat, char> TYPE_TRAIN_DATA;
typedef vector<TYPE_TRAIN_DATA> TYPE_TRAIN_DATA_SET;

typedef pair<Point2f, Point2f> TYPE_FERN_LEAF; // save pixel comparision
typedef vector<vector<TYPE_FERN_LEAF> > TYPE_FERN_FERNS; // save all ferns

static const TYPE_OF_PT PT_ERROR = TYPE_OF_PT(-1, -1);
static const TYPE_MF_BB BB_ERROR = TYPE_MF_BB(PT_ERROR, PT_ERROR);

typedef Rect TYPE_BBOX;

static const bool OF_USE_OPENCV = 1;

static const int MF_HALF_PATCH_SIZE = 4; // NNC patch size
static const int MF_NPTS = 12; // number of points in the patch(both vertical and horizontal)
static const int MF_FB_ERROR_DIST = 10; // threshold of detecting confusing condition

//跟踪状态反馈
static const int MF_TRACK_SUCCESS = 0;
static const int MF_TRACK_F_PTS = -1; // number of points after filtering is too little
static const int MF_TRACK_F_BOX = -2; // result box is out of bounds
static const int MF_TRACK_F_CONFUSION = -3; // tracking result is disordered
static const int MF_TRACK_F_BOX_SMALL = -4; // input box is too small

static const int MF_REJECT_OFERROR = 1 << 0; // filtered by OF error
static const int MF_REJECT_NCC = 1 << 1; // filtered by NCC
static const int MF_REJECT_FB = 1 << 2; // filtered by Forward-Backward

static const int TLD_TRACK_SUCCESS = 1;
static const int TLD_TRACK_FAILED = 0;

static const Scalar COLOR_GREEN = Scalar(156, 188, 26);
static const Scalar COLOR_BLUE = Scalar(219, 152, 52);
static const Scalar COLOR_BLACK = Scalar(94, 73, 52);
static const Scalar COLOR_WHITE = Scalar(241, 240, 236);
static const Scalar COLOR_YELLOW = Scalar(15, 196, 241);
static const Scalar COLOR_RED = Scalar(60, 76, 231);
static const Scalar COLOR_PURPLE = Scalar(182, 89, 155);


static const bool NCC_USE_OPENCV = 0; // 1(lower speed): use matchTemplate(), 0(faster)
static const bool NCC_FAST = 1; // 1 : my own implementation
static const bool RND_SHUFFLE_STD = 1;
static const bool QUIET_MODE = 0;
static const bool SHOW_NEW_NN_SAMPLES = 1;

static const float RF_FEA_SHIFT = 1.f / 5;
static const float RF_FEA_OFF = RF_FEA_SHIFT;

static void outputInfo(const string module ,const string info)
{
    if(QUIET_MODE) return;
    
    cerr << "[" << module << "] " << info << endl;
}

#endif

//头文件：MedianFlow.h

#ifndef __MedianFlow__MedianFlow__
#define __MedianFlow__MedianFlow__

#include "MFSystemStruct.h"

#include <cmath>
#include <iostream>

#include "OpticalFlow.h"

using namespace std;
using namespace cv;

class MedianFlow
{
private:
    Mat prevImg, nextImg;
    
    OpticalFlow *opticalFlow, *opticalFlowSwap;
    
    bool isPointInside(const TYPE_MF_PT &pt, const TYPE_MF_COORD border = 0);
    bool isBoxUsable(const TYPE_MF_BB &rect);
    
    void generatePts(const TYPE_MF_BB &box, vector<TYPE_MF_PT> &ret);
    
    float calcNCC(const Mat &img0, const Mat &img1);
    
    void filterOFError(const vector<TYPE_MF_PT> &pts, const vector<uchar> &retF, vector<int> &rejected);
    void filterFB(const vector<TYPE_MF_PT> &initialPts, const vector<TYPE_MF_PT> &FBPts, vector<int> &rejected);
    void filterNCC(const vector<TYPE_MF_PT> &initialPts, const vector<TYPE_MF_PT> &FPts, vector<int> &rejected);
    
    TYPE_MF_BB calcRect(const TYPE_MF_BB &rect, const vector<TYPE_MF_PT> &pts, const vector<TYPE_MF_PT> &FPts,  const vector<TYPE_MF_PT> &FBPts, const vector<int> &rejected, int &status);
    
public:
    
    MedianFlow();
 
    // prevImg & nextImg should be CV_8U
    MedianFlow(const Mat &prevImg, const Mat &nextImg);
    
    ~MedianFlow();
    
    static bool compare(const pair<float, int> &a, const pair<float, int> &b);
    
    TYPE_MF_BB trackBox(const TYPE_MF_BB &inputBox, int &status);
};

#endif /* defined(__MedianFlow__MedianFlow__) */

//头文件：OpticalFlow.h

#ifndef __MedianFlow__OpticalFlow__
#define __MedianFlow__OpticalFlow__

#include <vector>
#include <iostream>
#include <opencv2/opencv.hpp>

#include "MFSystemStruct.h"

using namespace std;
using namespace cv;

class OpticalFlow
{
private:
    Mat prevImg, nextImg;
    
public:
    OpticalFlow();
    
    // prevImg & nextImg should be CV_8U
    OpticalFlow(const Mat &prevImg, const Mat &nextImg);
    
    ~OpticalFlow();
    
    void trackPts(vector<TYPE_OF_PT> &pts, vector<TYPE_OF_PT> &retPts, vector<uchar> &status);

};

#endif /* defined(__MedianFlow__OpticalFlow__) */

//MedianFlow.cpp

#include "MedianFlow.h"

MedianFlow::MedianFlow() {}

MedianFlow::MedianFlow(const Mat &prevImg, const Mat &nextImg)
{
    this->prevImg = prevImg;

    this->nextImg = nextImg;
    
    opticalFlow = new OpticalFlow(this->prevImg, this->nextImg);
    opticalFlowSwap = new OpticalFlow(this->nextImg, this->prevImg);
}

MedianFlow::~MedianFlow()
{
    delete opticalFlow;
    opticalFlow = NULL;
    
    delete opticalFlowSwap;
    opticalFlowSwap = NULL;
}

void MedianFlow::generatePts(const TYPE_MF_BB &_box, vector<TYPE_MF_PT> &ret)
{
    TYPE_MF_PT tl(max(0.f, _box.tl().x), max(0.f, _box.tl().y));
    TYPE_MF_PT br(min((float)prevImg.cols, _box.br().x), min((float)prevImg.rows, _box.br().y));
    TYPE_MF_BB box(tl, br);
    
    float stepX = (float)(box.width - 2 * MF_HALF_PATCH_SIZE) / (MF_NPTS - 1);
    float stepY = (float)(box.height - 2 * MF_HALF_PATCH_SIZE) / (MF_NPTS - 1);
    int x0 = box.x + MF_HALF_PATCH_SIZE;
    int y0 = box.y + MF_HALF_PATCH_SIZE;
    
    if(!ret.empty()) ret.clear();
    
    for(int fx = 0; fx < MF_NPTS; fx++)
    {
        for(int fy = 0; fy < MF_NPTS; fy++)
        {
            ret.push_back(TYPE_MF_PT(x0 + fx * stepX, y0 + fy * stepY));
        }
    }
}

bool MedianFlow::compare(const pair<float, int> &a, const pair<float, int> &b)
// caution : prefix static can only be specified inside the class definition
{
    return a.first < b.first;
}

bool MedianFlow::isPointInside(const TYPE_MF_PT &pt, const TYPE_MF_COORD alpha)
{
    int width = prevImg.cols, height = prevImg.rows;
    return (pt.x >= 0 + alpha) && (pt.y >= 0 + alpha) && (pt.x <= width - alpha) && (pt.y <= height - alpha);
}

bool MedianFlow::isBoxUsable(const TYPE_MF_BB &rect)
{
    int width = prevImg.cols, height = prevImg.rows;
 
    // bounding box is too large
    if(rect.width > width || rect.height > height) return false;
    
    // intersection between rect and img is too small
    TYPE_MF_COORD tlx = max((TYPE_MF_COORD)rect.tl().x, (TYPE_MF_COORD)0);
    TYPE_MF_COORD tly = max((TYPE_MF_COORD)rect.tl().y, (TYPE_MF_COORD)0);
    TYPE_MF_COORD brx = min((TYPE_MF_COORD)rect.br().x, (TYPE_MF_COORD)width);
    TYPE_MF_COORD bry = min((TYPE_MF_COORD)rect.br().y, (TYPE_MF_COORD)height);
    
    TYPE_MF_BB bb(tlx, tly, brx - tlx, bry - tly);
    if(bb.width <= 2 * MF_HALF_PATCH_SIZE || bb.height <= 2 * MF_HALF_PATCH_SIZE) return false;
 
    // otherwise
    return true;
}

void MedianFlow::filterOFError(const vector<TYPE_MF_PT> &pts, const vector<uchar> &status, vector<int> &rejected)
{
    for(int i = 0; i < pts.size(); i++)
    {
        if(status[i] == 0) rejected[i] |= MF_REJECT_OFERROR;
    }
}

void MedianFlow::filterFB(const vector<TYPE_MF_PT> &initialPts, const vector<TYPE_MF_PT> &FBPts, vector<int> &rejected)
{
    int size = int(initialPts.size());
    vector<pair<float, int> > V;
    
    for(int i = 0; i < size; i++)
    {
        if(rejected[i] & MF_REJECT_OFERROR) continue;
        
        float dist = norm(Mat(initialPts[i]), Mat(FBPts[i]));
        V.push_back(make_pair(dist, i));
    }
    
    sort(V.begin(), V.end(), compare);
    
    for(int i = (int)V.size() / 2; i < V.size(); i++)
    {
        rejected[V[i].second] |= MF_REJECT_FB;
    }
}

float MedianFlow::calcNCC(const cv::Mat &img0, const cv::Mat &img1)
{
    if(NCC_USE_OPENCV)
    {
        Mat nccMat;
        matchTemplate(img0, img1, nccMat, CV_TM_CCORR_NORMED);
        
        return nccMat.at<float>(0);
    }
    else
    {
        Mat v0, v1; // convert image to 1 dimension vector
        
        img0.convertTo(v0, CV_32F);
        img1.convertTo(v1, CV_32F);
        
        v0 = v0.reshape(0, v0.cols * v0.rows);
        v1 = v1.reshape(0, v1.cols * v1.rows);
        
        Scalar mean, stddev;
        meanStdDev(v0, mean, stddev);
        v0 -= mean.val[0];
        meanStdDev(v1, mean, stddev);
        v1 -= mean.val[0];
        
        Mat v01 = v0.t() * v1;
        
        float norm0, norm1;
        
        norm0 = norm(v0);
        norm1 = norm(v1);
        
        return v01.at<float>(0) / norm0 / norm1;
    }
}

void MedianFlow::filterNCC(const vector<TYPE_MF_PT> &initialPts, const vector<TYPE_MF_PT> &FPts, vector<int> &rejected)
{
    int size = int(initialPts.size());
    vector<pair<float, int> > V;
    
    for(int i = 0; i < size; i++)
    {
        if(rejected[i] & MF_REJECT_OFERROR) continue;
        
        if(!isPointInside(initialPts[i], MF_HALF_PATCH_SIZE)) continue;
        if(!isPointInside(FPts[i], MF_HALF_PATCH_SIZE)) continue;
        
        Point2d win(MF_HALF_PATCH_SIZE, MF_HALF_PATCH_SIZE);
        Point2d pt1(initialPts[i].x, initialPts[i].y);
        Point2d pt2(FPts[i].x, FPts[i].y);
        
        // must be int
        Rect_<int> rect0(pt1 - win, pt1 + win);
        Rect_<int> rect1(pt2 - win, pt2 + win);
        
        float ncc = calcNCC(this->prevImg(rect0), this->nextImg(rect1));
        
        V.push_back(make_pair(ncc, i));
    }
    
    sort(V.begin(), V.end(), compare);
    
    //for(int i = int(V.size()) / 2; i < V.size(); i++)   //原代码：感觉这样不对，跟TLD原文的NCC的过滤方法不一致,故修改如下
	for(int i = 0; i < int(V.size()) / 2; i++)            //应该剔除相似性小于中值的
    {
        rejected[V[i].second] |= MF_REJECT_NCC;
    }
}

TYPE_MF_BB MedianFlow::calcRect(const TYPE_MF_BB &rect, const vector<TYPE_MF_PT> &pts, const vector<TYPE_MF_PT> &FPts, const vector<TYPE_MF_PT> &FBPts, const vector<int> &rejected, int &status)
{
    const int size = int(pts.size());
    
    vector<TYPE_MF_COORD> dxs, dys;
    
    for(int i = 0; i < size; i++)
    {
        if(rejected[i]) continue;
        
        dxs.push_back(FPts[i].x - pts[i].x);
        dys.push_back(FPts[i].y - pts[i].y);
    }
    
    if(dxs.size() <= 1)
    {
        status = MF_TRACK_F_PTS;
        outputInfo("Tracker", "Error : Too little points after filter.");
        return BB_ERROR;
    }
    
    sort(dxs.begin(), dxs.end());
    sort(dys.begin(), dys.end());
    
    TYPE_MF_COORD dx = dxs[dxs.size() / 2];
    TYPE_MF_COORD dy = dys[dys.size() / 2];
    TYPE_MF_PT delta(dx, dy);
    
    vector<float> ratios;
    vector<float> absDist;
    
    for(int i = 0; i < size; i++)
    {
        if(rejected[i]) continue;
        
        for(int j = i + 1; j < size; j++)
        {
            if(rejected[j]) continue;
            
            float dist0 = norm(Mat(pts[i]), Mat(pts[j]));
            float dist1 = norm(Mat(FPts[i]), Mat(FPts[j]));
            float ratio = dist1 / dist0;
            
            ratios.push_back(ratio);
        }
    }
    
    sort(ratios.begin(), ratios.end());
    float ratio = ratios[ratios.size() / 2];
    
    TYPE_MF_BB ret(delta + rect.tl(), delta + rect.br());
    
    TYPE_MF_PT center((float)(ret.tl().x + ret.br().x) / 2, (float)(ret.tl().y + ret.br().y) / 2);
    TYPE_MF_PT tl(center.x - ret.width / 2 * ratio, center.y - ret.height / 2 * ratio);
    TYPE_MF_PT br(center.x + ret.width / 2 * ratio, center.y + ret.height / 2 * ratio);
    
    ret = TYPE_MF_BB(tl, br);
    
    for(int i = 0; i < size; i++)
    {
        if(rejected[i] == MF_REJECT_OFERROR) continue;
        
        float dist = norm(Mat(pts[i]), Mat(FBPts[i]));
        
        absDist.push_back(dist);
    }

    sort(absDist.begin(), absDist.end());
    
    float medianAbsDist = absDist[(int)absDist.size() / 2];
    
    //for(auto &i : absDist) //caution : must add '&'
    //    i = abs(i - medianAbsDist);

	for(int i = 0; i < absDist.size(); i++)
		absDist[i] = abs(absDist[i] - medianAbsDist);
    
    sort(absDist.begin(), absDist.end());

    if(medianAbsDist > MF_FB_ERROR_DIST)
    {
        status = MF_TRACK_F_CONFUSION;
        outputInfo("Tracker", "Error : Large foward-backward distance.");
        return BB_ERROR;
    }
    
    if(!isBoxUsable(ret))
    {
        status = MF_TRACK_F_BOX;
        outputInfo("Tracker", "Error : Result bounding box is unusable.");
        return BB_ERROR;
    }

	ret.x = max(0.0f, ret.x), ret.y = max(0.0f, ret.y);//修改：保证最终结果跟踪框的左上角坐标不越界
    status = MF_TRACK_SUCCESS;
    outputInfo("Tracker", "Tracked successfully.");
    return ret;
}

TYPE_MF_BB MedianFlow::trackBox(const TYPE_MF_BB &inputBox, int &status)
{
    if(!isBoxUsable(inputBox))
    {
        status = MF_TRACK_F_BOX;
        outputInfo("Tracker", "Error : Input bounding box is unusable.");
        return BB_ERROR;
    }
    
    vector<TYPE_MF_PT> pts;
    generatePts(inputBox, pts);
    
    vector<TYPE_MF_PT> retF, retFB;
    vector<uchar> statusF, statusFB;
    
    retF = retFB = pts;
    
    opticalFlow->trackPts(pts, retF, statusF);
    opticalFlowSwap->trackPts(retF, retFB, statusFB);

    vector<int> rejected(MF_NPTS * MF_NPTS);
    
    filterOFError(retF, statusF, rejected);
    filterOFError(retFB, statusFB, rejected);
    
    filterFB(pts, retFB, rejected);
    filterNCC(pts, retF, rejected);
    
    TYPE_MF_BB ret;
    
    ret = calcRect(inputBox, pts, retF, retFB, rejected, status);
    
    if(status != MF_TRACK_SUCCESS)
    {
        return BB_ERROR;
    }
    
    return ret;
}

//OpticalFlow.cpp

#include "OpticalFlow.h"

OpticalFlow::OpticalFlow() {}

OpticalFlow::OpticalFlow(const Mat &prevImg, const Mat &nextImg)
{
    this->prevImg = prevImg;
    this->nextImg = nextImg;
}

OpticalFlow::~OpticalFlow() {}

void OpticalFlow::trackPts(vector<TYPE_OF_PT> &pts, vector<TYPE_OF_PT> &retPts, vector<uchar> &status)
{
    vector<float> err;
        
    calcOpticalFlowPyrLK(prevImg, nextImg, pts, retPts, status, err, Size(21, 21), 5, cvTermCriteria(CV_TERMCRIT_ITER | CV_TERMCRIT_EPS, 30, 0.01), OPTFLOW_USE_INITIAL_FLOW);
}

//main.cpp
/* 主程序功能：
       输入一组视频图片序列，使用MedianFlow算法跟踪目标框 
	   初始目标框选择：在第一帧手动选择目标框，按空格键继续跟踪
 */
#undef UNICODE  //使用多字节字符集，也可以在项目配置属性-常规中更改字符集
#include <opencv2/core/core.hpp>
#include <opencv2/highgui/highgui.hpp>
#include "opencv2/video/tracking.hpp"
#include <iostream>
#include <fstream>
#include <sstream>
#include <stdio.h>
#include <string.h>
#include "MedianFlow.h"

using namespace cv;
using namespace std;

void drawObjectRectangle(int event, int x, int y, int flags, void*);
void readImageSequenceFiles(char* ImgFilePath,vector <string> &imgNames);

Mat firstFrame;
Point previousPoint, currentPoint;
Rect initObjectRect;
bool lButtonDownFlag = false;

int main(int argc, char * argv[])
{
	char imgFilePath[100];
	memset(imgFilePath, 0, sizeof(imgFilePath));
	strcpy(imgFilePath, "./kitesurf");

	char tmpDirPath[MAX_PATH+1];
	
	Rect rect; // [x y width height] tracking position

	vector <string> imgNames;
    
	readImageSequenceFiles(imgFilePath, imgNames);

	Mat frame;
	Mat grayImg, lastImg;

	sprintf(tmpDirPath, "%s/", imgFilePath);
	imgNames[0].insert(0,tmpDirPath);
	frame = imread(imgNames[0], 1);

	namedWindow("MedianFlow", WINDOW_AUTOSIZE);
    imshow("MedianFlow", frame);
	firstFrame = frame;
    setMouseCallback("MedianFlow", drawObjectRectangle, 0);
    waitKey(0); //等待在第一帧：左键单击画一个矩形框作为追踪目标区域，按空格键继续
	rect = initObjectRect;
	cout << "#" << 0 << " " << rect.x << "  " << rect.y << "  " << rect.width << "  " << rect.height << endl; 

    cvtColor(frame, grayImg, CV_RGB2GRAY);    

	char strFrame[20];

    FILE* resultStream;
	resultStream = fopen("TrackingResults.txt", "w");
	fprintf (resultStream,"#0  %i %i %i %i\n",(int)rect.x,(int)rect.y,(int)rect.width,(int)rect.height);

	for(int i = 1; i < imgNames.size() - 1; i++)
	{
        imgNames[i].insert(0,tmpDirPath);
        
		lastImg = grayImg.clone();//记录上一帧
		frame = imread(imgNames[i], 1);// get frame
		cvtColor(frame, grayImg, CV_RGB2GRAY);

		//MedianFlow framework
		MedianFlow mf(lastImg, grayImg);
		int trackResStatus = MF_TRACK_SUCCESS;
		Rect rectTmp = mf.trackBox(rect, trackResStatus); 
		if(trackResStatus == MF_TRACK_SUCCESS) rect = rectTmp;

		rectangle(frame, rect, Scalar(200,0,0),2);// Draw rectangle
		cout << "#" << i << " " << rect.x << "  " << rect.y << "  " << rect.width << "  " << rect.height << endl; 

		fprintf (resultStream,"#%d  %i %i %i %i\n", i, (int)rect.x,(int)rect.y,(int)rect.width,(int)rect.height);

		sprintf(strFrame, "#%d ",i);

		putText(frame,strFrame,cvPoint(0,20),2,1,CV_RGB(25,200,25));
		
		imshow("MedianFlow", frame);// Display
		waitKey(100);		
	}
	fclose(resultStream);

	system("pause");
	return 0;
}

void drawObjectRectangle(int event, int x, int y, int flags, void*)
{
	if(event == EVENT_LBUTTONDOWN)
	{
		previousPoint = Point(x, y);
		lButtonDownFlag = true;
	}else if(event == EVENT_MOUSEMOVE && lButtonDownFlag == true){
		Mat tmp;
        firstFrame.copyTo(tmp);
        currentPoint = Point(x, y);
        rectangle(tmp, previousPoint, currentPoint, Scalar(200, 0, 0), 2);
        imshow("MedianFlow", tmp);
	}else if(event == EVENT_LBUTTONUP){
		initObjectRect.x = previousPoint.x;
		initObjectRect.y = previousPoint.y;
		initObjectRect.width = abs(currentPoint.x - previousPoint.x);
		initObjectRect.height = abs(currentPoint.y - previousPoint.y);
		lButtonDownFlag = false;
	}
}

//将./data目录下的图片名字全部存入imageNames
void readImageSequenceFiles(char* imgFilePath,vector <string> &imgNames)
{	
	imgNames.clear();

	char tmpDirSpec[MAX_PATH+1];
	sprintf (tmpDirSpec, "%s/*", imgFilePath);
	 
	WIN32_FIND_DATA f;
	HANDLE h = FindFirstFile(tmpDirSpec , &f); //read . --modified by Lumengru
	if(h != INVALID_HANDLE_VALUE)
	{
		FindNextFile(h, &f);	//read ..  --modified by Lumengru
		FindNextFile(h, &f);	//read the first picture  --modified by Lumengru
		do
		{
			imgNames.push_back(f.cFileName);
		} while(FindNextFile(h, &f)); 
	}
	FindClose(h);	
}