数字图像的直方图处理
实验四 数字图像的直方图处理
一、实验目的
- 掌握灰度直方图的概念及其计算方法;
- 熟练使用matlab实现直方图均衡化的计算过程;
二、实验环境
- PC计算机
- MatLab软件/语言包括图像处理工具箱(Image Processing Toolbox)
- 实验所需要的图片
- 实验原理
自己写直方图均衡化的原理!
四、实验图像
五、实验步骤和结果
- 用:imhist(I,64)和histeq(I) 进行直方图均衡化并显示图像I的直方图!
clear
close all
clc
I=imread('E:\ins风\hai.tif');
subplot(221)
I=rgb2gray(I);
imshow(I)
title('灰度图像')
subplot(222)
imhist(I,64)
title('原始图像直方图')
I1=histeq(I);
subplot(223)
imshow(I1)
title('图像均衡化')
subplot(224)
imhist(I1,64)
title('均衡化h后直方图')
2. 用:J = histeq(I,hgram); %函数完成直方图规定化,其中I是原始图像矩阵,J是变换后所得的图像矩阵;hgram是由用户制定的矢量。
例如:hgram=0:255得到的规定直方图为:
注意:将每一步的函数执行语句和实验结果写入实验报告。
clear
close all
clc
I=imread('E:\ins风\yun.tif');
I=rgb2gray(I);
subplot(221)
imshow(I);
title('灰度图像')
hgram=0:255;%规定化函数
J=histeq(I,hgram);
subplot(222)
imshow(J)
title('图像规定化')
subplot(223)
imhist(I)
title('原始图像直方图')
subplot(224)
imhist(J)
title('规定化后直方图')
- 进入执行直方图均衡化函数histeq,将理论课程学习的直方图均衡化算法步骤和代码对应起来,即,在该代码核心部分进行注释,写明对应理论算法中的哪个步骤。
打开的函数如下图所示:
分析如下:
if nargin == 1%当输入的参数只有一个图像
%HISTEQ(I)
validateattributes(a,{'uint8','uint16','double','int16','single'}, ...
{'nonsparse'}, mfilename,'I',1);%判读参数的输入
n = 64; % Default n
hgram = ones(1,n)*(numel(a)/n);%初始化hgram向量
isIntensityImage = true;
仅向 histeq 函数传入一个参数时,histeq 函数初始化包括 n,hgram 在内的若干个相关变量。
elseif nargin == 2%当输入的参数有两个时
if numel(cm) == 1%两个参数为图像和灰度级
%HISTEQ(I,N)
validateattributes(a,{'uint8','uint16','double','int16','single'}, ...
{'nonsparse'}, mfilename,'I',1);
validateattributes(cm, {'single','double'},...
{'nonsparse','integer','real','positive','scalar'},...
mfilename,'N',2);%判断参数输入
m = cm;
向 histeq 函数传入两个参数时,初始化包括 n,hgram 在内的若干个相关变量。
else
%HISTEQ(I,HGRAM)
validateattributes(a,{'uint8','uint16','double','int16','single'}, ...
{'nonsparse'}, mfilename,'I',1);
validateattributes(cm, {'single','double'},...
{'real','nonsparse','vector','nonempty'},...
mfilename,'HGRAM',2);
hgram = cm;
isIntensityImage = true;
end
向 histeq 函数传入两个参数:一个是图像,一个是用户指定的 hgram 向量
% Intensity image or indexed image
if isIntensityImage%判断图像是否增强
classChanged = false;
if isa(a,'int16')
classChanged = true;
a = im2uint16(a);
end
[nn,cum] = computeCumulativeHistogram(a,NPTS);%计算均匀化后的矩阵
T = createTransformationToIntensityImage(a,hgram,m,NPTS,nn,cum);
% Mex call is equivalent to:
% b = uint8((255.0*T(a+1));
% or uint16((65535.0*T(a+1));
b = images.internal.builtins.grayxform(a, T);%颜色匹配
判断是否是图像增强,计算累计直方图,进行图像增强,最后进行颜色匹配。
function T = createTransformationToIntensityImage(a,hgram,m,n,nn,cum)
% Create transformation to an intensity image by minimizing the error
% between desired and actual cumulative histogram.
% Generate cumulative hgram
cumd = cumsum(hgram);
% Calc error
% tol = nn w/ 1st and last element set to 0, then divide by 2 and tile to MxN
tol = ones(m,1)*min([nn(1:n-1),0;0,nn(2:n)])/2;
% Calculate errors btw cumulative histograms
err = (cumd(:)*ones(1,n)-ones(m,1)*cum(:)')+tol;%计算直方图误差
% Find which combo yielded errors above tolerance
d = find(err < -numel(a)*sqrt(eps));
if ~isempty(d)
% Set to max err
err(d) = numel(a)*ones(size(d));%超出允许的误差组的组合
end
% Get min error
% T will be the bin mapping of a to hgram
% T(oldbinval) = newbinval
[dum,T] = min(err); %#ok 得到最小差值
% Normalize T
T = (T-1)/(m-1);%灰度级变化
- 实验思考
- 为什么直方图均衡化后的图像的直方图不是完全均衡的?
答:因为像素的取值跨度大,像素灰度值的动态范围大。
因为直方图是PDF(概率密度函数)的近似,而且在处理中,不允许造成新的灰度级,所以在实际的直方图均衡应用中,很少见到完美平坦的直方图。因此,直方图均衡技术不能保证直方图的均匀分布,但是却可以扩展直方图的分布范围,也就意味着在直方图上,偏向左的暗区和偏向右的亮区都有像素分布,只是不能保证每个灰度级上都有像素分布。
直方图均衡化处理的“中心思想”是把原始图像的灰度直方图从比较集中的某个灰度区间变成在全部灰度范围内的均匀分布。直方图均衡化就是对图像进行非线性拉伸,重新分配图像像素值,使一定灰度范围内的像素数量大致相同。直方图均衡化就是把给定图像的直方图分布改变成“均匀”分布直方图分布。由于离散图象的直方图也是离散的,其灰度的累积分布函数是一个不减的阶梯函数。如果映射后的图象仍能取到所有256级灰度,那一定是原图象没有任何改变,这种情况只可能发生在原图象的直方图已经是一条水平线的情况下。一般情况下映射后所得到的图象只能取到少于256级灰度,这样在变换后的直方图中会有某些灰度级空缺,当然这些空缺应该均匀分布在0到255之间。于是问题就变成了将原有的256个值,即各灰度的概率,按顺序分成n(n<256)份,每份的概率总和应该相等。显然这个问题是不一定有解的,因此我们只能找到一个近似解。其结果就是最后得到一幅有空缺且不太平坦的直方图。
- 直方图均衡化后的图像的灰度级可能会比原始图像减少,这称之为什么现象?
答:伪轮廓现象;直方图均衡化会造成灰度级的合并,均衡化处理后的图象只能是近
似均匀分布。均衡化图象的动态范围扩大了,但其本质是扩大了量化间隔,而量化级别反而减少了,因此,原来灰度不同的象素经处理后可能变的相同,形成了一片的相同灰度的区域,各区域之间有明显的边界,从而出现了伪轮廓。
直方图均衡化的机理,是将对应像素较少的几个连续灰度级合并成一个灰度级,通过减少灰度级来实现的均衡化。所以直方图均衡化后的图像灰度级可能比原始图像减少。