【OpenCV】直方图应用：直方图均衡化，直方图匹配，对比直方图

看着写的不错就给转载了。

转载：http://blog.csdn.net/xiaowei_cqu/article/details/7606607

直方图均衡化

直方图均衡化（Histogram Equalization）是直方图最典型的应用，是图像点运算的一种。对于一幅输入图像，通过运算产生一幅输出图像，点运算是指输出图像的每个像素点的灰度值由输入像素点决定，即：

直方图均衡化是通过灰度变换将一幅图像转换为另一幅具有均衡直方图，即在每个灰度级上都具有相同的象素点数过程。从分布图上的理解就是希望原始图像中y轴的值在新的分布中尽可能的展开。变换过程是利用累积分布函数对原始分布进行映射，生成新的均匀拉伸的分布。因此对应每个点的操作是寻找原始分布中y值在均匀分布中的位置，如下图是理想的单纯高斯分布映射的示意图：

（图片来源：《Learnning OpenCV》 p189）

OpenCV中的cvEqualizeHist

OpenCV中有灰度直方图均衡化的函数cvEqualizeHist，接口很明朗：

[cpp] view plain copy print ?

1. void cvEqualizeHist( const CvArr* src, CvArr* dst ); 

void cvEqualizeHist( const CvArr* src, CvArr* dst );

注意此函数只能处理单通道的灰色图像，对于彩色图像，我们可以把每个信道分别均衡化，再Merge为彩色图像。

实践：图像直方图均衡化

[cpp] view plain copy print ?

1. int main() 
2. { 
3.     IplImage * image= cvLoadImage("baboon.jpg"); 
4.     //显示原图及直方图 
5.     myShowHist("Source",image); 
6.      
7.     IplImage* eqlimage=cvCreateImage(cvGetSize(image),image->depth,3); 
8.     //分别均衡化每个信道 
9.     IplImage* redImage=cvCreateImage(cvGetSize(image),image->depth,1); 
10.     IplImage* greenImage=cvCreateImage(cvGetSize(image),image->depth,1); 
11.     IplImage* blueImage=cvCreateImage(cvGetSize(image),image->depth,1); 
12.     cvSplit(image,blueImage,greenImage,redImage,NULL); 
13.  
14.     cvEqualizeHist(redImage,redImage); 
15.     cvEqualizeHist(greenImage,greenImage);  
16.     cvEqualizeHist(blueImage,blueImage);  
17.     //均衡化后的图像 
18.     cvMerge(blueImage,greenImage,redImage,NULL,eqlimage); 
19.     myShowHist("Equalized",eqlimage); 
20. } 

int main()
{
	IplImage * image= cvLoadImage("baboon.jpg");
	//显示原图及直方图
	myShowHist("Source",image);
	
	IplImage* eqlimage=cvCreateImage(cvGetSize(image),image->depth,3);
	//分别均衡化每个信道
	IplImage* redImage=cvCreateImage(cvGetSize(image),image->depth,1);
	IplImage* greenImage=cvCreateImage(cvGetSize(image),image->depth,1);
	IplImage* blueImage=cvCreateImage(cvGetSize(image),image->depth,1);
	cvSplit(image,blueImage,greenImage,redImage,NULL);

	cvEqualizeHist(redImage,redImage);
	cvEqualizeHist(greenImage,greenImage); 
	cvEqualizeHist(blueImage,blueImage); 
	//均衡化后的图像
	cvMerge(blueImage,greenImage,redImage,NULL,eqlimage);
	myShowHist("Equalized",eqlimage);
}

原始图像及灰度直方图如下：

均衡化后的直方图：

直方图匹配

直方图匹配又叫直方图规定化（Histogram Normalization/Matching）是指对一副图像进行变换，使其直方图与另一幅图像的直方图或特定函数形式的直方图进行匹配。应用场景如不同光照条件下的两幅图像，我们可以在比较两幅图像前先进行匹配变化。

参考shlkl99上传的直方图匹配代码，将图像规定化为高斯分布函数。

[cpp] view plain copy print ?

1. //将图像与特定函数分布histv[]匹配 
2. void myHistMatch(IplImage *img,double histv[]) 
3. { 
4.     int bins = 256; 
5.     int sizes[] = {bins}; 
6.     CvHistogram *hist = cvCreateHist(1,sizes,CV_HIST_ARRAY); 
7.     cvCalcHist(&img,hist); 
8.     cvNormalizeHist(hist,1); 
9.     double val_1 = 0.0; 
10.     double val_2 = 0.0; 
11.     uchar T[256] = {0}; 
12.     double S[256] = {0}; 
13.     double G[256] = {0}; 
14.     for (int index = 0; index<256; ++index) 
15.     { 
16.         val_1 += cvQueryHistValue_1D(hist,index); 
17.         val_2 += histv[index]; 
18.         G[index] = val_2; 
19.         S[index] = val_1; 
20.     } 
21.  
22.     double min_val = 0.0; 
23.     int PG = 0; 
24.     for ( int i = 0; i<256; ++i) 
25.     { 
26.         min_val = 1.0; 
27.         for(int j = 0;j<256; ++j) 
28.         { 
29.             if( (G[j] - S[i]) < min_val && (G[j] - S[i]) >= 0) 
30.             { 
31.                 min_val = (G[j] - S[i]); 
32.                 PG = j; 
33.             } 
34.  
35.         } 
36.         T[i] = (uchar)PG; 
37.     } 
38.  
39.     uchar *p = NULL; 
40.     for (int x = 0; x<img->height;++x) 
41.     {  
42.         p = (uchar*)(img->imageData + img->widthStep*x); 
43.         for (int y = 0; y<img->width;++y) 
44.         { 
45.             p[y] = T[p[y]]; 
46.         } 
47.     } 
48. } 
49.  
50. // 生成高斯分布 
51. void GenerateGaussModel(double model[]) 
52. { 
53.     double m1,m2,sigma1,sigma2,A1,A2,K; 
54.     m1 = 0.15; 
55.     m2 = 0.75; 
56.     sigma1 = 0.05; 
57.     sigma2 = 0.05; 
58.     A1 = 1; 
59.     A2 = 0.07; 
60.     K = 0.002; 
61.  
62.     double c1 = A1*(1.0/(sqrt(2*CV_PI))*sigma1); 
63.     double k1 = 2*sigma1*sigma1; 
64.     double c2 = A2*(1.0/(sqrt(2*CV_PI))*sigma2); 
65.     double k2 = 2*sigma2*sigma2; 
66.     double p = 0.0,val= 0.0,z = 0.0; 
67.     for (int zt = 0;zt < 256;++zt) 
68.     { 
69.         val = K + c1*exp(-(z-m1)*(z-m1)/k1) + c2*exp(-(z-m2)*(z-m2)/k2); 
70.         model[zt] = val; 
71.         p = p +val; 
72.         z = z + 1.0/256; 
73.     } 
74.     for (int i = 0;i<256; ++i) 
75.     { 
76.         model[i] = model[i]/p; 
77.     } 
78. }  

//将图像与特定函数分布histv[]匹配
void myHistMatch(IplImage *img,double histv[])
{
	int bins = 256;
	int sizes[] = {bins};
	CvHistogram *hist = cvCreateHist(1,sizes,CV_HIST_ARRAY);
	cvCalcHist(&img,hist);
	cvNormalizeHist(hist,1);
	double val_1 = 0.0;
	double val_2 = 0.0;
	uchar T[256] = {0};
	double S[256] = {0};
	double G[256] = {0};
	for (int index = 0; index<256; ++index)
	{
		val_1 += cvQueryHistValue_1D(hist,index);
		val_2 += histv[index];
		G[index] = val_2;
		S[index] = val_1;
	}

	double min_val = 0.0;
	int PG = 0;
	for ( int i = 0; i<256; ++i)
	{
		min_val = 1.0;
		for(int j = 0;j<256; ++j)
		{
			if( (G[j] - S[i]) < min_val && (G[j] - S[i]) >= 0)
			{
				min_val = (G[j] - S[i]);
				PG = j;
			}

		}
		T[i] = (uchar)PG;
	}

	uchar *p = NULL;
	for (int x = 0; x<img->height;++x)
	{ 
		p = (uchar*)(img->imageData + img->widthStep*x);
		for (int y = 0; y<img->width;++y)
		{
			p[y] = T[p[y]];
		}
	}
}

// 生成高斯分布
void GenerateGaussModel(double model[])
{
	double m1,m2,sigma1,sigma2,A1,A2,K;
	m1 = 0.15;
	m2 = 0.75;
	sigma1 = 0.05;
	sigma2 = 0.05;
	A1 = 1;
	A2 = 0.07;
	K = 0.002;

	double c1 = A1*(1.0/(sqrt(2*CV_PI))*sigma1);
	double k1 = 2*sigma1*sigma1;
	double c2 = A2*(1.0/(sqrt(2*CV_PI))*sigma2);
	double k2 = 2*sigma2*sigma2;
	double p = 0.0,val= 0.0,z = 0.0;
	for (int zt = 0;zt < 256;++zt)
	{
		val = K + c1*exp(-(z-m1)*(z-m1)/k1) + c2*exp(-(z-m2)*(z-m2)/k2);
		model[zt] = val;
		p = p +val;
		z = z + 1.0/256;
	}
	for (int i = 0;i<256; ++i)
	{
		model[i] = model[i]/p;
	}
} 

实践：直方图匹配

对示例图片每个信道分别进行匹配处理

对比直方图

OpenCV中提供了cvCompareHist函数用以对比两个直方图的相似度：

[cpp] view plain copy print ?

1. double cvCompareHist(  
2.              const CvHistogram* hist1, //直方图1 
3.              const CvHistogram* hist2, //直方图2 
4.              int method//对比方法 
5. ); 

double cvCompareHist( 
             const CvHistogram* hist1, //直方图1
             const CvHistogram* hist2, //直方图2
             int method//对比方法
);

method有CV_COMP_CORREL， CV_COMP_CHISQR，CV_COMP_INTERSECT，CV_COMP_BHATTACHARYYA四种方法，对应公式如下：

实践：对比不同光照条件的两幅图像直方图

直方图的对比主要用以判断两幅图像的匹配度，我们试验以下两幅图像直方图对比的结果：

     

[cpp] view plain copy print ?

1. int main() 
2. { 
3.     IplImage * image= cvLoadImage("myhand1.jpg"); 
4.     IplImage * image2= cvLoadImage("myhand2.jpg"); 
5.     int hist_size=256; 
6.     float range[] = {0,255}; 
7.     float* ranges[]={range}; 
8.      
9.     IplImage* gray_plane = cvCreateImage(cvGetSize(image),8,1); 
10.     cvCvtColor(image,gray_plane,CV_BGR2GRAY); 
11.     CvHistogram* gray_hist = cvCreateHist(1,&hist_size,CV_HIST_ARRAY,ranges,1); 
12.     cvCalcHist(&gray_plane,gray_hist,0,0); 
13.  
14.  
15.     IplImage* gray_plane2 = cvCreateImage(cvGetSize(image2),8,1); 
16.     cvCvtColor(image2,gray_plane2,CV_BGR2GRAY); 
17.     CvHistogram* gray_hist2 = cvCreateHist(1,&hist_size,CV_HIST_ARRAY,ranges,1); 
18.     cvCalcHist(&gray_plane2,gray_hist2,0,0); 
19.  
20.     //相关：CV_COMP_CORREL     
21.     //卡方：CV_COMP_CHISQR 
22.     //直方图相交：CV_COMP_INTERSECT 
23.     //Bhattacharyya距离：CV_COMP_BHATTACHARYYA 
24.     double  com=cvCompareHist(gray_hist,gray_hist2,CV_COMP_BHATTACHARYYA); 
25.  
26.     cout<<com<<endl; 
27. } 

int main()
{
	IplImage * image= cvLoadImage("myhand1.jpg");
	IplImage * image2= cvLoadImage("myhand2.jpg");
	int hist_size=256;
	float range[] = {0,255};
	float* ranges[]={range};
	
	IplImage* gray_plane = cvCreateImage(cvGetSize(image),8,1);
	cvCvtColor(image,gray_plane,CV_BGR2GRAY);
	CvHistogram* gray_hist = cvCreateHist(1,&hist_size,CV_HIST_ARRAY,ranges,1);
	cvCalcHist(&gray_plane,gray_hist,0,0);


	IplImage* gray_plane2 = cvCreateImage(cvGetSize(image2),8,1);
	cvCvtColor(image2,gray_plane2,CV_BGR2GRAY);
	CvHistogram* gray_hist2 = cvCreateHist(1,&hist_size,CV_HIST_ARRAY,ranges,1);
	cvCalcHist(&gray_plane2,gray_hist2,0,0);

	//相关：CV_COMP_CORREL    
	//卡方：CV_COMP_CHISQR
	//直方图相交：CV_COMP_INTERSECT
	//Bhattacharyya距离：CV_COMP_BHATTACHARYYA
	double  com=cvCompareHist(gray_hist,gray_hist2,CV_COMP_BHATTACHARYYA);

	cout<<com<<endl;
}

输出结果为：0.396814

cvCompareHist的结果为【0,1】的浮点数，越小表示两幅图匹配度越高，0.0时两幅图精确匹配。（可以试验两幅完全一样的图即为0.0）。

针对上面两幅图，我们分别先进行直方图匹配变化：

之后使用cvCompareHist()对比两幅图像的直方图，输出结果为0.267421

表明两幅图的匹配度变高了。

注意method用不同的方法对比结果是不同的。

应用

通过对比我们可以设置EMD的阈值来判定皮肤或手的ROI。《Learnning OpenCV》后面有相应的练习题：收集手的肤色直方图，对比室内，室外手的直方图的EMD距离，利用这些测量值设置一个距离阈值。

a.利用该阈值检测第三幅图（如室外阴影），看能能否很好的检测肤色直方图。

b.随机选择不是肤色的背景块直方图，观测EMD变化，试验与真实肤色对比时能否很好的拒绝背景。

如上也即是直方图对比可以应用的场景

实验代码下载：http://download.csdn.net/detail/xiaowei_cqu/4332914