直方图的概念理解
图像直方图是反映一个图像像素分布的统计表，其实横坐标代表了图像像素的种类，可以是灰度的，也可以是彩色的。纵坐标代表了每一种颜色值在图像中的像素总数或者占所有像素个数的百分比。图像是由像素构成，因为反映像素分布的直方图往往可以作为图像一个很重要的特征。
图像灰度直方图：一副数字图像有[0~255]灰度级，直方图定义如下：
h(gk)=nk
其中，gk是第k个灰度级（如：255），nk是该灰度级的个数。
归一化直方图定义如下:
p(gk)=h(gk)n=nkn
即：第k个灰度级出现的数量，比上所有灰度级数量总和。
直方图均衡化
均衡化的基本思想是：尽量使得每个灰度级的像素数量相等。
图1-1. 直方图均衡化的效果（截自《数字图像处理》）
直方图均衡化的作用：使得图像灰度级跨越更宽的灰度级范围，从而提高图像对比度；
直方图均衡化的另一个优势是：不需要额外参数，整个过程是自动的；
直方图均衡化的缺点是：拉伸后有些灰度级可能不被映射到，造成图像观感上的颗粒感；
算法步骤：
1）计算图像灰度直方图，并进行归一化（除以总像素数）
2）计算归一化直方图的累积直方图 fHist[ L-1 ]
3）对原始图灰度S，变换后的灰度D = fHist[S] * (L-1)
MATLAB中图像灰度直方图操作：
matlab中直方图函数为：imhist()，h=imhist(img,b);
其中:img是输入图像，b是容器，就是分几块显示，默认为256，即256个灰度级。 若要得到归一化后的图像，使用：
h_normal=imhist(img,b)/numel(img);
其中numel(img)是求出所有灰度级数量总和，即图像总像素数。
%% 直方图均衡化
clear all,close all,clc;
H= imread('./image/6.jpg');
if length(size(H))>2
H=rgb2gray(H);
end
subplot(3,2,1);
imshow(H); title('原图');
subplot(3,2,2);
imhist(H); title('原图直方图');
subplot(3,2,3);
H1=adapthisteq(H);
imshow(H1); title('adapthisteq均衡后图');
subplot(3,2,4);
imhist(H1);title('adapthisteq均衡后直方图');
subplot(3,2,5);
H2=histeq(H);
imshow(H2); title('histeq均衡后图');
subplot(3,2,6);
imhist(H1); title('histeq均衡后直方图');
当然按照直方图及均衡的定义，也可以自己完成，方案1：
%% 直方图自己绘制
clc,clear all,close all;
H= imread('./image/6.jpg');
if length(size(H))>2
H=rgb2gray(H);
end
[m,n]=size(H);
p=zeros(1,256);
for i=0:255
p(i+1)=length(find(H==i))/(m*n);
end
subplot(2,2,1);
imshow(H);
title('原图');
subplot(2,2,2);
bar(0:255,p,'b');
title('原图直方图');
s=zeros(1,256);
for i=1:256
for j=1:i
s(i)=p(j)+s(i);
end
end
a=round(s*255);
b=H;
for i=0:255
b(find(H==i))=a(i+1);
end
subplot(2,2,3);
imshow(b)
title('均衡化后图像');
for i=0:255
GPeq(i+1)=sum(p(find(a==i)));
end
subplot(2,2,4);
bar(0:255,GPeq,'b'); title('均衡化后的直方图');
方案2：
clc,clear all,close all;
%读入并将其转化为灰度图，提取图像的长和宽。
Img= imread('./image/6.jpg');
if length(size(Img))>2
Img=rgb2gray(Img);
end
%绘制原始图像的直方图
[height,width]=size(Img);
[counts1, x] = imhist(Img,256);
counts2 = counts1/height/width;
figure(1),
subplot(1,2,1),
imshow(Img);title('原始图像');
subplot(1,2,2),
stem(x, counts2); title('原始图像直方图');
%统计每个灰度的像素值累计数目
NumPixel = zeros(1,256);%统计各灰度数目，共256个灰度级
for i = 1:height
for j = 1: width
%对应灰度值像素点数量+1
%NumPixel的下标是从1开始，而图像像素的取值范围是0~255，所以用NumPixel(Img(i,j) + 1)
NumPixel(Img(i,j) + 1) = NumPixel(Img(i,j) + 1) + 1;
end
end
%将频数值算为频率
ProbPixel = zeros(1,256);
for i = 1:256
ProbPixel(i) = NumPixel(i) / (height * width * 1.0);
end
%函数cumsum来计算cdf，并将频率（取值范围是0.0~1.0）映射到0~255的无符号整数
CumuPixel = cumsum(ProbPixel);
CumuPixel = uint8(255 .* CumuPixel + 0.5);
%直方图均衡。赋值语句右端，Img(i,j)被用来作为CumuPixel的索引
for i = 1:height
for j = 1: width
Img(i,j) = CumuPixel(Img(i,j)+1);
end
end
%显示更新后的直方图
figure(2),
subplot(1,2,1),
imshow(Img); title('直方图均衡化图像');
[counts1, x] = imhist(Img,256);
counts2 = counts1/height/width;
subplot(1,2,2),
stem(x, counts2); title('直方图均衡化后图像的直方图');
除了对灰度图像进行处理，当然也可以对彩色图像进行直方图均衡，只不过需要分通道均衡罢了。
%% 彩色图像直方图均衡
clc,clear all,close all;
Img= imread('./image/6.jpg');
OutImg=Img;
R = Img(:,:,1);
G = Img(:,:,2);
B = Img(:,:,3);
R = histeq(R, 256);
G = histeq(G, 256);
B = histeq(B, 256);
OutImg(:,:,1) = R;
OutImg(:,:,2) = G;
OutImg(:,:,3) = B;
figure,
subplot(1,2,1),
imshow(Img);title('原始图像');
subplot(1,2,2),
imshow(OutImg); title('均衡化后结果');
另外，在还提到了先将RGB空间转换到HSV空间，然后再进行直方图均衡操作，我这里也做了尝试：
%% 转换到HSV空间均衡
clc,clear all,close all;
Img= imread('./image/6.jpg');
hsvImg = rgb2hsv(Img);
V=hsvImg(:,:,3);
[height,width]=size(V);
V = uint8(V*255);
NumPixel = zeros(1,256);
for i = 1:height
for j = 1: width
NumPixel(V(i,j) + 1) = NumPixel(V(i,j) + 1) + 1;
end
end
ProbPixel = zeros(1,256);
for i = 1:256
ProbPixel(i) = NumPixel(i) / (height * width * 1.0);
end
CumuPixel = cumsum(ProbPixel);
CumuPixel = uint8(255 .* CumuPixel + 0.5);
for i = 1:height
for j = 1: width
V(i,j) = CumuPixel(V(i,j)+1); %注意，这里需要+1，要不然会出问题
end
end
V = im2double(V);
hsvImg(:,:,3) = V;
outputImg = hsv2rgb(hsvImg);
figure,
subplot(1,2,1),
imshow(Img);title('原始图像');
subplot(1,2,2),
imshow(outputImg); title('在HSV空间均衡化后结果');
Opencv中图像直方图操作：
opencv中函数为：
void calcHist(const Mat* images, int nimages, const int* channels, InputArray mask, OutputArray hist, int dims, const int* histSize, const float** ranges, booluniform=true, bool accumulate=false )其中：Mat为输入图像；int nimages为计算直方图的图像数，通常情况下都只作用于单一图像，所以通常nimages=1;
const int* channels：图像的通道，它是一个数组，如果是灰度图像则channels[1]={0};如果是彩色图像则channels[3]={0,1,2}；如果是是求彩色图像某一个通道的直方图，则channels[1]={i};
IuputArray mask:用于确定哪些点参与计算，默认情况设置为一个空图像，即：Mat();
OutArray hist：计算得到的直方图;
int dims：直方图的维数，灰度图像为1维，彩色图像为3维;
const int* histSize：直方图横坐标的区间数。如果是10，则它会横坐标分为10份，然后统计每个区间的像素点总和。
const float** ranges：这是一个二维数组，用来指出每个区间的范围。
直方图绘制实例：
#include <opencv2\opencv.hpp>
#include <iostream>
using namespace cv;
using namespace std;
//灰度图像直方图显示
Mat getHistImg(const MatND& hist)
{
double maxVal = 0;
double minVal = 0;
minMaxLoc(hist, &minVal, &maxVal, 0, 0);
int histSize = hist.rows;
Mat histImg(histSize, histSize, CV_8U, Scalar(255));
// 设置最大峰值为图像高度的90%
int hpt = static_cast<int>(0.9*histSize);
for (int h = 0; h < histSize; h++)
{
float binVal = hist.at<float>(h);
int intensity = static_cast<int>(binVal*hpt / maxVal);
line(histImg, Point(h, histSize), Point(h, histSize - intensity), Scalar::all(0));
}
return histImg;
}
int main()
{
Mat Image = imread("D:\\fcq_proMatlab\\vessel_edge_extration\\image\\Penguins.jpg");
Mat imhist;
//灰度图像直方图
cvtColor(Image, Image, CV_BGR2GRAY);
const int channels[1] = { 0 };
const int histSize[1] = { 256 };
float hranges[2] = { 0, 255 };
const float* ranges[1] = { hranges };
MatND hist;
calcHist(&Image, 1, channels, Mat(), hist, 1, histSize, ranges);
imhist = getHistImg(hist);
imshow("直方图", imhist);
彩色图像直方图
//const int channels[3] = { 0, 1, 2 };
//const int histSize[3] = { 256, 256, 256 };
//float hranges[2] = { 0, 255 };
//const float* ranges[3] = { hranges, hranges, hranges };
//MatND hist;
//calcHist(&Image, 1, channels, Mat(), hist, 3, histSize, ranges);
cvWaitKey();
system("pause");
return 0;
}opencv中直方图均衡化的源码程序如下：
CV_IMPL void cvEqualizeHist( const CvArr* srcarr, CvArr* dstarr )
{
CvMat sstub, *src = cvGetMat(srcarr, &sstub);
CvMat dstub, *dst = cvGetMat(dstarr, &dstub);
CV_Assert( CV_ARE_SIZES_EQ(src, dst) && CV_ARE_TYPES_EQ(src, dst) &&
CV_MAT_TYPE(src->type) == CV_8UC1 );
CvSize size = cvGetMatSize(src);
if( CV_IS_MAT_CONT(src->type & dst->type) )
{
size.width *= size.height;
size.height = 1;
}
int x, y;
const int hist_sz = 256;//0到255，一共256个灰度值
int hist[hist_sz];
memset(hist, 0, sizeof(hist));
for( y = 0; y < size.height; y++ )
{
const uchar* sptr = src->data.ptr + src->step*y;
for( x = 0; x < size.width; x++ )
hist[sptr[x]]++; //这里实现了hist中存储各灰度值出现的次数
}
float scale = 255.f/(size.width*size.height);
int sum = 0;
uchar lut[hist_sz+1];
for( int i = 0; i < hist_sz; i++ )
{
sum += hist[i]; //逐级累计
int val = cvRound(sum*scale);
lut[i] = CV_CAST_8U(val);
}
lut[0] = 0;//这是通过程序计算出的灰度值映射关系,不管你怎么映射,0肯定是0
for( y = 0; y < size.height; y++ )
{
const uchar* sptr = src->data.ptr + src->step*y;
uchar* dptr = dst->data.ptr + dst->step*y;
for( x = 0; x < size.width; x++ )
dptr[x] = lut[sptr[x]];
}
}
下面是实例程序及结果：
#include <opencv2\opencv.hpp>
using namespace cv;
using namespace std;
int main()
{
灰度图像均衡
//Mat srcImage = imread("D:\\fcq_proMatlab\\vessel_edge_extration\\image\\6.jpg");
//if (!srcImage.data)
// return 1;
//Mat srcGray;
//cvtColor(srcImage, srcGray, CV_BGR2GRAY);
//imshow("srcGray", srcGray);
直方图均衡化
//Mat heqResult;
//equalizeHist(srcGray, heqResult);
//imshow("heqResult", heqResult);
//彩色图像均衡
Mat srcImage = imread("D:\\fcq_proMatlab\\vessel_edge_extration\\image\\6.jpg");
if (!srcImage.data)
return 1;
// 存储彩色直方图及图像通道向量
Mat colorHeqImage;
vector<cv::Mat> BGR_plane;
// 对BGR通道进行分离,分别对各个通道进行直方图均衡化
split(srcImage, BGR_plane);
for (int i = 0; i < BGR_plane.size(); i++)
cv::equalizeHist(BGR_plane[i], BGR_plane[i]);
merge(BGR_plane, colorHeqImage);
imshow("srcImage", srcImage);
imshow("colorHeqImage", colorHeqImage);
waitKey(0);
system("pause");
}效果图：
自适应直方图均衡化(Adaptive histgram equalization/AHE)
限制对比度自适应直方图均衡（Contrast Limited Adaptive histgram equalization/CLAHE)
关于自适应均衡方法，还是单独拿出来写吧，东西太多了看起开也不方便~
以上算法实现都亲自测试过，如果有问题还请告知啊。
参考：
《冈萨雷斯--图像处理》
http://blog.jobbole.com/84215/
本文章为转载内容，我们尊重原作者对文章享有的著作权。如有内容错误或侵权问题，欢迎原作者联系我们进行内容更正或删除文章。
