PROJECT 03-02

PROJECT 03-02 [Multiple Uses]
Histogram Equalization  
(a) Write a computer program for computing the histogram of an image.
(b) Implement the histogram equalization technique discussed in Section 3.3.1.
(c) DownloadFig. 3.8(a) and perform histogram equalization on it.

As a minimum, your report should include the original image, a plot of its histogram, a plot of the histogram-equalization transformation function, the enhanced image, and a plot of its histogram. Use this information to explain why the resulting image was enhanced as it was.

利用MATLAB自带的函数处理非常简单

%显示图像直方图
f=imread('F:\大学学习\大三上\数字图像处理\实验\PROJECT\DIP3E_Original_Images_CH03\Fig0308(a)(fractured_spine).tif');

figure
imshow(f);%显示原图像

%h=imhist(f);
%h1=h(1:10:256);
%horz=1:10:256;
%figure(3)
%bar(horz,h1);%绘制条形图，将水平轴等分为若干单元，降低水平轴分辨率
%figure(4)
%stem(horz,h1,'fill');

figure
imhist(f);%缺省的灰度级b=256
title('imhist函数绘制的原图像直方图') 

%Matlab函数进行直方图均衡化
J=histeq(f);
figure,imhist(J)
title('histeq函数均衡化后的直方图') 

figure,imshow(J)
title('histeq函数均衡化后的图像') 

也可以自己编写。算法在冈萨雷斯的数字图像处理里有。（参考网上程序改编，原地址找不到了）

%自己编写的直方图均衡化程序

Clear all 
%一，图像的预处理，读入彩色图像将其灰度化 
f=imread('F:\大学学习\大三上\数字图像处理\实验\PROJECT\DIP3E_Original_Images_CH03\Fig0308(a)(fractured_spine).tif');
%二，绘制直方图 
[m,n]=size(f); %测量图像尺寸参数 
GP=zeros(1,256); %预创建存放灰度出现概率的向量 
for k=0:255 
GP(k+1)=length(find(f==k))/(m*n); %计算每级灰度出现的概率，将其存入GP中相应位置 
end 

figure,bar(0:255,GP) %绘制直方图 
title('原图像直方图') 
xlabel('灰度值') 
ylabel('出现概率') 

%三，直方图均衡化 
S1=zeros(1,256); 
for i=1:256 
for j=1:i 
S1(i)=GP(j)+S1(i); %计算Sk 
end 
end 
S2=round((S1*256)+0.5); %将Sk归到相近级的灰度 
for i=1:256 
GPeq(i)=sum(GP(find(S2==i))); %计算现有每个灰度级出现的概率 
end 
figure,bar(0:255,GPeq,'b') %显示均衡化后的直方图 
title('均衡化后的直方图') 
xlabel('灰度值') 
ylabel('出现概率') 

%四，图像均衡化 
PA=f; 
for i=0:255 
PA(find(f==i))=S2(i+1); %将各个像素归一化后的灰度值赋给这个像素 
end 
figure,imshow(PA) %显示均衡化后的图像 
title('均衡化后图像') 
imwrite(PA,'PicEqual.bmp'); 

%a=tabulate(f(:));%统计数组中元素频数a(:,2)、频率a(:,3) /100

改写为函数形式

%http://blog.sina.com.cn/s/blog_550f161701000820.html

function myhist = myimhist(ima)
%绘制直方图
%输入参数为ima，输出参数为myhist
myhist = zeros(256,1);
imau = double(ima);
[xs,ys] = size(imau);
for ix = 1:xs
    for iy = 1:ys
        myhist(imau(ix,iy)+1) = myhist(imau(ix,iy)+1)+1;
    end
end
return;

function myheqf = myhisteqf(hist)
%生成直方图均衡化变换函数
%输入参数为hist，输出参数为myheqf
myheqd = zeros(256,1);
myhist = double(hist)/sum(hist);
myheqd(1) = myhist(1);
for idx = 2:256
    myheqd(idx) = myheqd(idx-1) + myhist(idx);
end
myheqf = round(myheqd*255);
return;

function proj03_02(ima)
%直方图均衡化增强
WH = 2; WL = 3;
subplot(WH,WL,1),imshow(ima),title('原图像');
myhist = myimhist(ima);
subplot(WH,WL,2),bar(myhist),title('原图的直方图');
myheqf = myhisteqf(myhist);
subplot(WH,WL,3),plot(myheqf),title('直方图均衡化变换函数图');
imeq = double(ima);
[xs,ys] = size(imeq);
for ix = 1:xs;
    for iy = 1:ys;
        imeq(ix,iy) = myheqf(imeq(ix,iy)+1);
    end
end
imeq = uint8(imeq);
subplot(WH,WL,4),imshow(imeq),title('增强后的图像');
myhist = myimhist(imeq);
subplot(WH,WL,5),bar(myhist),title('增强后图像的直方图');