【OpenCV】直方图计算

之前学过OpenCV，知道一些基本数据结构和图像处理函数，但似乎掌握的不是很好，这次重温OpenCV，理清自己的思路，属于用到哪里就重温哪里，顺便记录下学习笔记，就大胆的贴出来了。

看到一篇自适应阈值二值化算法，里面有讲到将图像分成较小的块，然后分别计算每块的直方图，根据每个直方图的峰值，然后为每个块计算其阈值，每个像素点的阈值根据相邻的块的阈值进行插值获得。突然想起对于直方图，只记得在OpenCV中可以调用哪些函数来实现，更深的就不清楚了，所以再次学习一下直方图的计算，知其然还要知其所以然。

图像直方图是用以表示数字图像中亮度分布的直方图，标绘了图像中每个亮度值的像素数，可以借助观察该直方图了解需要如何调整亮度分布。这里介绍的是灰度直方图计算，灰度直方图是描绘图像中每一个灰度级像素的个数。横坐标表示灰度级，纵坐标表示图像中该灰度级像素出现的个数，它统计一幅图像中各个灰度级出现的次数或概率。

灰度级范围为[0，L-1]的数字图像直方图是离散函数 h(rk) = nk，其中rk表示第k级灰度值，nk是图像中灰度为rk的像素个数。一维直方图可用下式表示

                                           

其中S(xi)表示某像素的个数，分母部分表示图像像素总数。需要注意的是：直方图中有个特征空间子区段的概念，就是我们把灰度值范围分割成子区域(称作bins)，如下图

                                            

所以前面公式中的 h(xi) 应该认为是一个bins区间内的像素情况，当然bins的区间你可以取值256，这样就是统计每一个灰度级像素的情况。然后再统计位于每一个bins中的像素数目。采用这种方法来统计图像中的像素情况。

让我们再来搞清楚直方图的一些具体细节;

a. dims: 需要统计的特征数目。上例中dims= 1，表示只统计灰度值，若统计RGB图的各个通道的情况，则dims= 3；

b. bins：每个特征空间子区段的数目。就是把整个范围划分为多少个区间段；

c. range：每个特征空间的取值范围。在上例中range= [0，255]。

下面就讲讲OpenCV是怎么来实现直方图的绘制的。

//计算直方图
void calcHist(const Mat* images, int nimages, const int* channels, InputArray mask, 
              OutputArray hist, int dims, const int* histSize, const float** ranges,
              bool uniform=true, bool accumulate=false )

//images：输入图像
//nimages：输入图像的数目
//channels：需要统计的通道(dim)索引。输入单幅图像的情况下，从0到images.channels()-1
//mask：掩码
//hist：储存直方图的矩阵
//dims：直方图维数
//histSize：每个维度的bin数目，就是前面所说的bins
//ranges：直方图每一维的取值范围
//uniform：是否对齐。若设置为true，bin大小相同，则ranges[]是一个只包含两元素的数组，起始点和终止点
//         若设置为false，则bin大小是不同的，你需要在ranges[]中指定各个bin区间的大小
//accumulate：如果设置为true，则直方图在使用之前不清除，用于保存多幅图像集中的单一直方图，
//            或及时更新直方图

上面这个函数就可以实现图像的直方图了，但是我们还要将其显示出来，我们要显示直方图的话，就需要另外创建一幅图片，横坐标表示灰度级，纵坐标表示该灰度级像素的个数，为方便显示，我们通常会对直方图进行归一化。不然，万一某个灰度级的像素数目很多，那要是正常显示的话，这创建的图片的行数还是比较吓人的。将直方图数据显示在图片中可以用矩形的形式也可用线条形式，下面会演示这两种形式，先继续分析需要用到的函数。

//直方图归一化，就是将各个灰度级的像素数之间的比例关系，映射到归一化范围[alpha，beta]
//最大值对应到beta值，而不是简单的每级像素数除以范围值，这就是为什么描绘出来的直方图会置顶的原因
void normalize(InputArray src, OutputArray dst, double alpha=1, double beta=0, 
               int norm_type=NORM_L2, int dtype=-1, InputArray mask=noArray() )

//src：输入图像数组
//dst：归一化后的输出图像数组
//alpha和beta：是归一化后的取值极限，将直方图归一化到范围[alpha，beta]
//norm_type：归一化方法
//dtype = -1：表示归一化后的输出图像数组与输入同类型

接下来就是将直方图画出来，就是将那些值在图中标绘出来。这里介绍两个函数

//画直线，在图像img中画一条颜色为color，粗细为thickness，类型为lineType的直线
void line(Mat& img, Point pt1, Point pt2, const Scalar& color, int thickness=1,
          int lineType=8, int shift=0)
//两点确认一条直线。
//lineType：直线类型
//shift：坐标小数点维数

//画一个单一的实矩形
void rectangle(Mat& img, Point pt1, Point pt2, const Scalar& color, int thickness=1,
                 int lineType=8, int shift=0)
//一条对角线的两个顶点可确定一个矩形
//pt1和pt2互为对顶点
//thickness为负值表示矩形为实矩形

好，万事具备，贴代码

int main()
{
	Mat src, gray, hist;                //hist为存储直方图的矩阵
	src = imread("lena.jpg");
	cvtColor(src, gray, CV_BGR2GRAY);   //转换为灰度图

	int histSize = 256;
	float range[] = { 0, 256 };
	const float* histRange = { range };
	int channels[] = {0};
	bool uniform = true; bool accumulate = false;

	/*计算直方图*/
	calcHist(&gray, 1, channels, Mat(), hist, 1, &histSize, 
		     &histRange, uniform, accumulate);

	/*创建描绘直方图的“图像”，和原图大小一样*/
	int hist_w = src.cols; int hist_h = src.rows;
	int bin_w = cvRound((double)hist_w / histSize);

	Mat histImage(hist_w, hist_h, CV_8UC3, Scalar(0, 0, 0));

	/*直方图归一化范围[0，histImage.rows]*/
	normalize(hist, hist, 0, histImage.rows, NORM_MINMAX, -1, Mat());

	/*画直线*/
	for (int i = 1; i < histSize; ++i)
	{
		//cvRound：类型转换。 这里hist为256*1的一维矩阵，存储的是图像中各个灰度级的归一化值
		line(histImage, Point(bin_w*(i - 1), hist_h - cvRound(hist.at<float>(i - 1))),
			 Point(bin_w*(i), hist_h - cvRound(hist.at<float>(i))),
			 Scalar(0, 0, 255), 2, 8, 0);
	}

	imshow("figure_src", src);
	imshow("figure_hist", histImage);
	
	waitKey(0);
	return 0;
}

结果如下图所示

上面实现的是灰度图的直方图，那么若是实现RGB图多通道的直方图勒，注意到计算直方图的函数calcHist里面有个参数channles，我们可进行多通道索引，来得到多通道的直方图数据，这样我们程序需要进行小小的修改

int main()
{
	Mat src;
	Mat b_hist, g_hist, r_hist;                //*_hist为存储直方图的矩阵
	src = imread("lena1.jpg");

	int histSize = 256;
	float range[] = { 0, 256 };
	const float* histRange = { range };
	int channelsB[] = { 0 };
	int channelsG[] = { 1 };
	int channelsR[] = { 2 };
	bool uniform = true; bool accumulate = false;

	/*计算直方图*/
	calcHist(&src, 1, channelsB, Mat(), b_hist, 1, &histSize,
		     &histRange, uniform, accumulate);
	calcHist(&src, 1, channelsG, Mat(), g_hist, 1, &histSize,
		     &histRange, uniform, accumulate);
	calcHist(&src, 1, channelsR, Mat(), r_hist, 1, &histSize,
		     &histRange, uniform, accumulate);

	/*创建描绘直方图的“图像”，和原图大小一样*/
	int hist_w = src.cols; int hist_h = src.rows;
	int bin_w = cvRound((double)hist_w / histSize);

	Mat histImage(hist_w, hist_h, CV_8UC3, Scalar(0, 0, 0));

	/*直方图归一化*/
	normalize(b_hist, b_hist, 0, histImage.rows, NORM_MINMAX, -1, Mat());
	normalize(g_hist, g_hist, 0, histImage.rows, NORM_MINMAX, -1, Mat());
	normalize(r_hist, r_hist, 0, histImage.rows, NORM_MINMAX, -1, Mat());

	/*画直线*/
	for (int i = 1; i < histSize; ++i)
	{
		//cvRound：类型转换。 这里hist为256*1的一维矩阵，存储的是图像中各个灰度级的归一化值
		line(histImage, Point(bin_w*(i - 1), hist_h - cvRound(b_hist.at<float>(i - 1))),
			 Point(bin_w*(i), hist_h - cvRound(b_hist.at<float>(i))),
			 Scalar(255, 0, 0), 2, 8, 0);
		line(histImage, Point(bin_w*(i - 1), hist_h - cvRound(g_hist.at<float>(i - 1))),
			Point(bin_w*(i), hist_h - cvRound(g_hist.at<float>(i))),
			Scalar(0, 255, 0), 2, 8, 0);
		line(histImage, Point(bin_w*(i - 1), hist_h - cvRound(r_hist.at<float>(i - 1))),
			Point(bin_w*(i), hist_h - cvRound(r_hist.at<float>(i))),
			Scalar(0, 0, 255), 2, 8, 0);
	}

	imshow("figure_src", src);
	imshow("figure_hist", histImage);
	
	waitKey(0);
	return 0;
}

我们前面还介绍了一个画矩形的函数，这里我们也把各个通道的直方图显示在同一图像中，这里不同的是，需要扩展图片的宽度。这里贴一下部分修改代码

/*和原图等高，3倍宽，用于将三幅直方图显示在同一图像中*/
	Mat histImage(hist_h, hist_w * 3, CV_8UC3, Scalar(0, 0, 0));
	
	/*矩形的宽度单位*/
	int scale = 1;

	/*画矩形*/
	for (int i = 0; i < histSize; ++i)
	{
		//“坐标原点”在左上角
		rectangle(histImage, Point(i*scale, hist_h - 1), 
			Point((i+1)*scale, hist_h - cvRound(b_hist.at<float>(i))),
			Scalar(255, 0, 0), CV_FILLED);
		rectangle(histImage, Point((i + histSize)*scale, hist_h - 1), 
			Point((i + histSize + 1)*scale, hist_h - cvRound(g_hist.at<float>(i))),
			Scalar(0, 255, 0), CV_FILLED);
		rectangle(histImage, Point((i + histSize * 2)*scale, hist_h - 1), 
			Point((i + histSize * 2 + 1)*scale, hist_h - cvRound(r_hist.at<float>(i))),
			Scalar(0, 0, 255), CV_FILLED);
	}

结果如下图