graphcut.h

来自「动态场景中运动目标检测提取与跟踪 对新手很有用」· C头文件 代码 · 共 164 行

#ifndef GRAPHCUT_H
#define GRAPHCUT_H

#include "graph.h"
#include "math.h"

double getE2(int xr,int xs,int Ir, int Is,int bgr,int bgs);
#define K1   5
#define GAMAZ 10

void graphCuts(int height, int width, int** PV,double** w,int* bgV,bool* b)
{
	typedef Graph<double,double,double> GraphType;
	GraphType *g = new GraphType( width*height, (width-1)*height+(height-1)*width);
	//int node_num=0;
    for(int x=0;x<width;x++)
		for(int y=0;y<height;y++)
		{	
			g->add_node();
			//给t-link赋值，0代表source,1 代表sink
			//g->add_tweights(y*width+x,-log(w[y*width+x][0]),-log(w[y*width+x][1]));
		}
	for(x=0;x<width;x++)
		for(int y=0;y<height;y++)
		{
			//g->add_node();
			//给t-link赋值，0代表source,1 代表sink
			g->add_tweights(y*width+x,-log(w[y*width+x][0]),-log(w[y*width+x][1]));
		}

		//ofstream output("e2.txt",ios::ate);
		
	for(x=1;x<width;x++)
		for(int y=1;y<height;y++)
		{
			int index1=y*width+x;
			int index2=(y-1)*width+x;
			double E2=getE2(PV[index1][0],PV[index2][0],
				PV[index1][1],PV[index2][1],
				bgV[index1],bgV[index2]);
			if(E2<0)
			{
				E2=0;
			}
			//output<<x<<","<<y<<","<<PV[index1][1]<<","<<PV[index2][1]<<","<<E2<<",";
			g -> add_edge(index1, index2,  E2, E2 );//与上面的点
			index2=y*width+x-1;
			E2=getE2(PV[index1][0],PV[index2][0],
						PV[index1][1],PV[index2][1],
						bgV[index1],bgV[index2]);
			if(E2<0)
			{
				E2=0;
			}
			//output<<PV[index1][1]<<","<<PV[index2][1]<<","<<E2<<endl;
			g -> add_edge( index1,index2,  E2, E2 );//与左面的点
		}
	for(x=1;x<width;x++)
	{
		int index2=x-1;
		double E2=getE2(PV[x][0],PV[index2][0],
						PV[x][1],PV[index2][1],
						bgV[x],bgV[index2]);
		if(E2<0)
		{
			E2=0;
		}
		g -> add_edge( x,index2, E2, E2 );//与左面的点
	}
	for(int y=1;y<height;y++)
	{
		int index1=y*width;
		int index2=(y-1)*width;
		double E2=getE2(PV[index1][0],PV[index2][0],
						PV[index1][1],PV[index2][1],
						bgV[index1],bgV[index2]);
		if(E2<0)
		{
			E2=0;
		}
		g -> add_edge( index1,index2, E2, E2 );//与上面的点
	}
	double flow = g -> maxflow();
	for(int i=0;i<height;i++)
		for(int j=0;j<width;j++)
		{
			if(g->what_segment(i*width+j) == GraphType::SOURCE)
				b[i*width+j]=false;
			else
				b[i*width+j]=true;
		}
	delete g;
}

double getDrs(int Irb, int Isb);
double getDrs(int Ir, int Is, int Irb, int Irs);
double getE2(int xr,int xs,int Ir, int Is,int bgr,int bgs)
{
	double mm;
	double beita;//=1.0/(2*(getDrs(bgr,bgs)));
	if((mm=getDrs(bgr,bgs))==0)
	{
		beita=1;
	}
	else
	{
		beita=1.0/(2*mm);
	}
	
	//cout<<getDrs(bgr,bgs)<<endl;
	//cout<<"beita"<<beita<<endl;
	return fabs((abs(xr-xs)*exp(-beita*getDrs(Ir,Is,bgr,bgs))));
}

double getDrs(int Irb, int Isb)
{
	Irb &=0x00ffffff;
	int rred=(Irb>>16) & 0xff;
	int rgreen=(Irb>>8) & 0xff;
	int rblue=Irb & 0xff;
	Isb &=0x00ffffff;
	int sred=(Isb>>16) & 0xff;
	int sgreen=(Isb>>8) & 0xff;
	int sblue=Isb & 0xff;
	return ((rred-sred)*(rred-sred)+(rgreen-sgreen)*(rgreen-sgreen)
		+(rblue-sblue)*(rblue-sblue));
}

double getDrs(int Ir, int Is, int Irb, int Isb)
{
	Ir &=0x00ffffff;
	int rred=(Ir>>16) & 0xff;
	int rgreen=(Ir>>8) & 0xff;
	int rblue=Ir & 0xff;
	Is &=0x00ffffff;
	int sred=(Is>>16) & 0xff;
	int sgreen=(Is>>8) & 0xff;
	int sblue=Is & 0xff;
	double drs= ((rred-sred)*(rred-sred)+(rgreen-sgreen)*(rgreen-sgreen)
		+(rblue-sblue)*(rblue-sblue));
	
	Irb &=0x00ffffff;
    int rredb=(Irb>>16) & 0xff;
	int rgreenb=(Irb>>8) & 0xff;
	int rblueb=Irb & 0xff;
	Isb &=0x00ffffff;
	int sredb=(Isb>>16) & 0xff;
	int sgreenb=(Isb>>8) & 0xff;
	int sblueb=Isb & 0xff;
	double drsb= ((rredb-sredb)*(rredb-sredb)+(rgreenb-sgreenb)*(rgreenb-sgreenb)
		+(rblueb-sblueb)*(rblueb-sblueb));
	double zrs1=((rredb-rred)*(rredb-rred)+(rgreenb-rgreen)*(rgreenb-rgreen)
		+(rblueb-rblue)*(rblueb-rblue));
	double zrs2=((sredb-sred)*(sredb-sred)+(sgreenb-sgreen)*(sgreenb-sgreen)
		+(sblueb-sblue)*(sblueb-sblue));
	double zrs= (zrs1>zrs2 ? zrs1 : zrs2);
	
	return (drs/((1+(drsb/(K1*K1)))*exp(-zrs/GAMAZ)));
	
	
	
}

#endif