直方图均衡化的原理及其实现

    直方图均衡化在图像处理领域中运用非常广泛,而且非常简单易实现。

    首先我们了解一下什么是图像的直方图：

    设图像的灰度范围为[a,b]，r为此灰度范围内的任一灰度级，p(r)为这幅图像中灰度级为r的像素出现的频率，可以看出p(r)是r的函数。该函数的图形称为这幅图像的直方图。

     p(r)=灰度为r的像素数/图像上的总像素数

    
                                                                      图像灰度级统计信息
                                                                      图像灰度直方图
                                 
          
       原灰度图                     灰度直方图
    

    可以很清楚地看出，灰度直方图抛弃了原灰度图像的空间位置信息，反映了某一像素值在灰度图中出现的频率或者概率信息。如上图所示，可以很清楚地看出0-100灰度值所出现的频率非常大，且越趋近于0概率越高；而灰度值大于200后，出现频率大大降低。由此可以判断出，该图像整体较暗，细节部分不够突出。

    于是这里我们引入---直方图均衡化的方法，希望按照一定的变换公式，将原图映射到新图，使得新图在原图的基础上，直方图分布更加均匀，这样图像的明暗分布更加均匀，给人的视觉效果就是对比度好，细节清晰。

    那么问题的关键是，如何确定这个变换公式呢？

    均衡化方法中，使用直方图的累积分布函数作为变换公式：

    其中，

    实际上就是用某灰度级的累积概率来代替其原出现概率，得到映射后新的灰度值（累积概率乘以最大灰度值）。
    举个例子来说：

      rk代表原图的八个灰度级；nk代表每个灰度级出现的频数；Prk代表每个灰度级出现的概率；Sk代表累积概率；Ps代表新图中rk所对应的出现概率。

     于是，很容易得到，例如：

  原图rk=0，累积概率Sk=1/7，于是其对应的新图灰度值为1*1/7=1/7。

  经过意义映射之后，可以发现，新图中灰度级为0的出现概率为0，于是在rk=0时Ps=0;新图中灰度级为1/7的出现概率为0.19，于是在rk=1/7时Ps=0.19···

  以此类推，很容易换算出新图的灰度级。

  了解了直方图均衡化的原理，那么编程实际上就是小菜一碟了。无论是使用Matlab还是OpenCV，直方图均衡化都很好实现。这里要提一点的就是Matlab相对来说更加地方便，对于直方图的获取与绘制，都有现成的函数可供调用，而OpenCV相对来说麻烦一点，需要自己写绘制函数。

  另外，对输入进来的图像，首先要进行通道数的转换，比如说传入进来的是.jpg格式的图像，那么无论从外表上看它是不是一个灰度图像，读进来的都是三通道图像，所以，往往第一步是需要将多通道图像转换成单通道图像。

 

效果：
           
             原图                           均衡化后

核心函数（来源于网络，我是大自然的搬运工）

1.直方图均衡化

 

void MyEqualizeHist_color(IplImage* src, int color)
{
	if (color == 1)
	{//多通道
		IplImage* imgChannel[4] = { 0, 0, 0, 0 };
		IplImage* dst = cvCreateImage(cvGetSize(src), IPL_DEPTH_8U, 3);
		int i;
		for (i = 0; i < src->nChannels; i++)
		{
			imgChannel[i] = cvCreateImage(cvGetSize(src), IPL_DEPTH_8U, 1);  //要求单通道图像才能直方图均衡化  
		}
		//通道分离  
		cvSplit(src, imgChannel[0], imgChannel[1], imgChannel[2], imgChannel[3]);//BGRA  
		for (i = 0; i < dst->nChannels; i++)
		{
			//直方图均衡化，原始图像和目标图像必须是单通道  
			cvEqualizeHist(imgChannel[i], imgChannel[i]);
		}
		//通道组合  
		cvMerge(imgChannel[0], imgChannel[1], imgChannel[2], imgChannel[3], dst);
		cvNamedWindow("src", 1);
		cvShowImage("src", src);
		cvNamedWindow("Equalize", 1);
		cvShowImage("Equalize", dst);
	}
	else
	{//单通道
		IplImage* dst = cvCreateImage(cvGetSize(src), IPL_DEPTH_8U, 1);
		cvEqualizeHist(src, dst);
		cvNamedWindow("src", 1);
		cvShowImage("src", src);
		cvNamedWindow("Equalize", 1);
		cvShowImage("Equalize", dst);
	}
}

 2.绘制灰度直方图（单通道灰度图像）

 

IplImage* src = cvLoadImage("orange.JPG");//读取图像
IplImage* gray = cvCreateImage(cvGetSize(src), IPL_DEPTH_8U, 1);
cvCvtColor(src, gray, CV_BGR2GRAY);//彩色图像灰度化
CvHistogram *pcvHistogram = CreateGrayImageHist(&gray);//创建一个灰度直方图
int nHistImageWidth = 255; //0到255个灰度级
int nHistImageHeight = 200;  //直方图图像高度  
int nScale = 2;
IplImage* pHistImage = 
CreateHisogramImage(nHistImageWidth, nScale,nHistImageHeight, pcvHistogram);//得到灰度直方图
cvNamedWindow("Hist", 1);
cvShowImage("Hist", pHistImage);
/****
FillWhite是一个填充函数，将pImage填充成白色
*****/
void FillWhite(IplImage *pImage)
{
cvRectangle(pImage, cvPoint(0, 0), cvPoint(pImage->width, pImage->height), CV_RGB(255, 255, 255), CV_FILLED);
}
/****
创建灰度图像的直方图  
*****/
CvHistogram* CreateGrayImageHist(IplImage **ppImage)
{
	int nHistSize = 256;//灰度数
	float fRange[] = { 0, 255 };  //灰度级的范围    
	float *pfRanges[] = { fRange };
	CvHistogram *pcvHistogram = cvCreateHist(1, &nHistSize, CV_HIST_ARRAY, pfRanges);//创建直方图
	cvCalcHist(ppImage, pcvHistogram);//计算ppImage的直方图
	return pcvHistogram;
}
/****
根据直方图创建直方图图像、
*****/  
IplImage * CreateHisogramImage(int nImageWidth, int nScale, int nImageHeight, CvHistogram *pcvHistogram)
{
IplImage *pHistImage = cvCreateImage(cvSize(nImageWidth * nScale, nImageHeight), IPL_DEPTH_8U, 1);
FillWhite(pHistImage);
//统计直方图中的最大直方块  
float fMaxHistValue = 0;
cvGetMinMaxHistValue(pcvHistogram, NULL, &fMaxHistValue, NULL, NULL);
//分别将每个直方块的值绘制到图中  
int i;
for (i = 0; i < nImageWidth; i++)
{
float fHistValue = cvQueryHistValue_1D(pcvHistogram, i); //像素为i的直方块大小  
int nRealHeight = cvRound((fHistValue / fMaxHistValue) * nImageHeight);  //要绘制的高度  
cvRectangle(pHistImage,
		  cvPoint(i * nScale, nImageHeight - 1),
		  cvPoint((i + 1) * nScale - 1, nImageHeight - nRealHeight),
		  cvScalar(i, 0, 0, 0),
		  CV_FILLED
			);
	}
return pHistImage;
}