数字图像处理作业

基于直方图的图像增强方法姓名：

学号：导师：

基于直方图的图像增强方法

一．直方图均衡化算法

直方图均衡化是图像处理领域中利用图像直方图对对比度进行调整的方法，这种方法通常用来增加许多图像的局部对比度，尤其是当图像的有用数据的对比度相当接近的时候。通过这种方法，亮度可以更好地砸直方图上分布。这样就额可以拥有增强局部的对比度而不影响整体的对比度，直方图均衡化通过有效地扩展常用的亮度来实现这种功能。

1.1算法介绍与仿真

图像的直方图是图像的重要统计特征，是用来表征数字图像的每一灰度级与其出现的频率间的统计关系的方法。其数学公式如式（1）所示：

n k P (r k ) =(k =0, 1, 2,.... L -1) (1) N

式中，P （r k ）为图像f （x ，y ）的第k 级灰度出现的概率；r k 为第k 级灰度的灰度值级；n k 为图像中灰度值为r k 的像素的个数；N 为图像f （x ，y ）的总像素数；L 为总的灰度级数，需要注意的是直方图不表示图像的空间信息，且任一特定图像都有唯一直方图，但反之并不成立在图像直方图中，整体较暗的图像其直方图的组成部分集中在灰度低的一侧，而明亮的图像的直方图组成部分集中在灰度较高的一

侧。当一副图像其像素占据全部可能的灰度范围且分布均匀时，该图像有较高的对比度，其图像也比较清晰。因此，可以通过改变直方图的灰度级分布达到增强图像的目的。

直方图均衡化是一种以累积分布函数变换法为基础的直方图修正方法。其基本思路是将一已知灰度概率分布的图像经过变换，使之成为具有均匀灰度概率分布、输出图像的直方图近似服从均匀分布的变换算法。

其计算步骤如下：

1）列出原始图像的灰度级r k ，k ＝0，1，2，…，L －1，L 为总的灰度级数；

2）统计各灰度级像素数目n k ，k ＝0，1，2，…，L －1；

3）计算原始图像直方图各灰度级的频率P （r k ）＝n k ／N ，k ＝0， 1，2，…，L －1，N 为图像f （x ，y ）的总像素数；

4）计算累计分布函数C (k ) =∑P (r k ) , k =0, 1, 2,... L -1;

0L -1

5）计算映射后的输出图像灰度级g （i ）＝INT ［（gmax －gmin ）C （k ）＋gmin ＋0．5］，i ＝0，1， …，P －1，P 为输出图像灰度等级个数，

INT 为取整符号；

6）统计映射后各灰度等级的像素数目n i ，i ＝0，1，…，P －1；

7）计算输出图像直方图P （g i ）＝n i ／N ，i ＝0，1，…，P －1；

8）用r k 和g i 的映射关系修改原始图像的灰度级，获得均衡化后的输出图像。

为了验证直方图均衡方法是否可以达到增强图像的目的，下面我们就利用在宿舍拍摄的照片来进行验证，图像直方图均衡化前后效果如图1所示。由图1中的a 和b 可看出均衡后的图像对比度明显增强，亮度增大，原本图中的黑暗区域也可以看到了。由c 和d 可看出原图的灰度级集中在数值较低的部分而直方图均衡后使得直方图变得均匀分布了，在整体灰度级上都有分布。

图1

1.2算法存在的问题

由于直方图均衡化理论来源于连续函数而数字图像的的灰度是离散值应用于数字图像处理的变换函数进行了从连续到离散的近似。如积

分运算就变成了累加运算因而在直方图均衡存在以下问题。

1）量化误差

会造成原图像信息的丢失，原来的像素的总的灰度值为连续的而变换为离散的这样就存在量化误差如1.1中步骤5）中得取整，如变换后的灰度值为50.15和50.00由于灰度级只有256个所以只能把变换后的两个灰度值都量化到50，那这样就势必会造成原图像的信息丢失。丢失的一定是数量很少的像素。

2）无法增强局部细节

由于直方图均衡变换是针对整幅图像所有像素值进行统一变换，所以这种算法只能整体增强图像，而不能针对局部细节部分也进行增强。并且当灰度集中在地区域或高区域时这是会造成画面过亮，而整幅图像的灰度的范围没有提升，也达不到突出图像细节的目的。所以这样就提出了自适应直方图均衡化（AHE ）。

二．自适应直方图均衡化

上面的介绍的直方图均衡是图像中的每一点都进行运算, 也就是说, 灰度变换函数运算与像素所处的位置无关, 这种全局性处理的算法, 它具有算法简单, 计算速度较快等优点, 但由于它是对所有像素点都作同样的处理, 忽略了图像的局部特征, 这样导致经过直方图均衡法算法的图像将丢失图像的有用信息, 对图像的去噪处理及边缘检测带来损失。

那么该如何提取图像的局部特征呢? 采用局部法对图像进行处

理, 不同局部采用不同的对比度增强方法, 也就是说, 根据图像象素的局部统计特征来决定处理方法. 每个象素的灰度值都通过一个均衡化变换函数得到的, 而该变换函数是由以该象素中心的一个局部子图像的直方图得到的，称其为局部对比度增强法局部对比度增强法的计算公式为

x i , j = m i , j + k ( x i , j - m i , j ) （2 ）

其中x i, j , x i, j 分布为变换前、后的中心象素, m i , j =1x i , j 为窗W ∑m ⨯n (i , j ) ∈w 内象素的平均灰值, 从（2）可看出, 当k> 1 时, 如果x i, j > mi, j , 那么x i, j>x i, j, 否则x i, j

改变局部窗口内的对比度, 同样达到了细节增强的效果, 优化了图像质量. 下面我们详细地分析局部对比度增强法的具体算法。

为了保存图像中的细节部分, 在实施直方图均衡化之前, 先对原图的细节部分保存, 再到直方图均衡时将这些细节添入算法过程。于是, 图像均衡算法修改为:

x i , j 0≤x i , j ≤255⎧⎪T (x i , j +k (x i , j -m i , j )) 　　　=⎨ （3） T (x ) 　　　　　　　　　其它⎪i , j ⎩

其中x i, j和x i, j表示变换前后的图像灰度值, mi, j表示为以x i, j为中心的窗邻域均值, T表示对x i, j的变换函数。从式（3）中我们可以看出, T 起到了调节直方图动态范围的作用, k ( x i, j- mi, j )则相当于一个高通滤波器, 起到了增强细节( 即局部对比度) 的作用, 强化细节的同时, 也在增加高频噪声。自然想到, 我们该如何即可增强细节又能避免增加高频噪声呢? 其实, 要达到这种效果比较简单, 算法关键就在k 这个

自适应参量的选取上。

考虑到参量k 的自适应性, 我们先对k 进行分析。当窗W 的中心象素x i, j位于无图像细节处时，k 趋向于0，而当窗W 的中心象素x i,

k j 位于图像细节处时，取一个较大正值。基于以上考虑我们选取窗W

内邻域灰度方差作为自适应变量, 那么的表达式可以写为:

⎡σ2i , j ⎤k =k () -1 ⎢ 2n ⎥ ⎣⎦

其中k 为比例系数，σ2n 为整幅图像的噪声方差，σ

度方差。若令r =

方差σ2i , j 2i , j （4）为窗W 内的灰σi , j ，则式（2）可表为k=k[r2-1],由于图像窗W 内的灰度σn 不小于图像的噪声方差，可知r ≥1, 当r=1时，k=0，也即σi , j =σn 那么图像无细节存在，因此该领域中点的像素值只能增强，当r 越大，k 值也越大，由此可知r>1，得到 σi , j >σn ，该领域内有细节存在，而

k 值也变大，细节被增强。

综上所述, 我们得出自适应直方图均衡算法具体的实现过程如下。

1) 给出原始图像的所有灰度级i, 统计原始图像各灰度级的像素数ni

2) 计算原始图像的直方图与累积直方图。

3) 用式（4）计算k 值. 用式（3）进行计算局部灰度值。

4) 用式（2）计算局部对比度, 实现均衡化。

5) 用p (ti ) =

n i 计算新的直方图。 n

三．对比度受限的自适应直方图均衡

AHE 有过度放大图像中相同区域的噪音的问题，另外一种自适应的直方图均衡算法即限制对比度直方图均衡（CLAHE ）算法能有限的限制这种不利的放大。

3.1算法简介

对比度受限自适应直方图均衡法（CLAHE ）通过限制局部直方图的高度来限制局部对比度的增强幅度，从而限制噪声的放大及局部对比度的过增强。

对比度增强的幅度可定义为灰阶映射函数的斜率。假定自适应直方图均衡方法的滑动窗口的大小为M ×M ，则局部映射函数为：

255⨯CDF (i ) m (i ) = (5) M ⨯M

CDF(i)为滑动窗口局部直方图的累积分布函数。累积分布函数CDF(i)的导数为直方图Hist(i)，从而局部映射函数m(i)的斜率S 为：

d (m (i ) )255S ==Hist (i ) ⨯ (6) di M ⨯M

因此限制直方图高度就等效于限制映射函数m(i)的斜率S ，进而限制对比增强度。若限定最大斜率为Smax ，则允许的最大直方图高度为：

H max =S max ⨯M ⨯M (7) 255

从而，对于高度大于H max 的直方图应截去多余的那部分，如图1所示。由图可知，实际上是从阈值T （而非H max ）处对直方图进行截断，然后将截去的部分均匀地分布在整个灰阶范围上，以保证总的直

方图面积不变，从而使整个直方图上升一个高度L 。因此H max 、T 、L 三者之间应满足下面关系：

H max =T+L

(8)

最后，改进的直方图值为：

⎧Hist (i ) +L 　　　　Hist (i )

通过改变最大的映射函数斜率S max 及相应的最大直方图高度H max ，便可以获得不同增强效果的图像。

CLAHE 可以有效抑制局部对比度的增强及噪声放大。然而，在输出图像中仍然会产生大量的人为噪声，尤其是在灰阶突变的交界区域，这是由于灰阶突变交界区域的局部直方图剧烈变化而引起的。基本步骤：

1) 图像分块

2) 产生局部直方图

3) 分别对每个图像块进行限制对比度直方图均衡，每个图像块分别产生独立的变换函数

4) 将多个图像块的灰度查找表拼接成输出图像

其基本思想是构造限制函数，限制灰度级的概率密度，并将超过限制函数的像素点在直方图内进行重整。

具体算法：

1）定义限制函数clipLimit ，计算超过限制函数的总像素点数totalExcess

2）计算totalExcess 平均分到每个灰度级的平均值avgBinIncr

3）定义限制函数与平均值的差upperLimit = clipLimit - avgBinIncr;

4）对直方图进行限高处理

● if 第K 灰度级像素数>限制函数

●降低该灰度级像素数=限制函数

●else

●if 第K 灰度级像素数> upperLimit

●增高该灰度级概率密度=限制函数；

●并在超过限制函数的总像素点数totalExcess 中减去增加的像素点数

●Else 第K 灰度级像素数

●该灰度级像素点数+多余像素点平均值avgBinIncr ；

●并在超过限制函数的总像素点数totalExcess 中减去该平均值。

5）对直方图进行均衡

3.2仿真结果

图2

由上图的仿真结果可以看出来，限制对比度的直方图均衡得到的图片细节更加的突出，没有过亮或过暗的情况出现，弥补了图像直方图均衡在处理低灰度层密集分布的图像时, 不能有效增加图像动态范围与对比度的缺陷. 该算法与直方图匹配算法相比, 算法简单, 计算量小, 对图像处理系统要求低, 容易实现. 在图像实时处理领域中可以使用该算法来取代直方图匹配算法, 以实现图像直方图均衡的处理。

四．总结

通过这次的数字图像处理学习，我也更好的理解了什么是数字图像处理，并且也学到了很多图像处理的相关知识，包括图像增强、图像切割、图像变换、图像编码、图像分析等好多的知识，并通过本次的实验也熟悉了MATLAB 的编程和软件的操作，并且了解了上述三种算法的过程与优劣并且自己仿真，所以说还是收益颇丰的。

五．附录

下面这部分是所涉及到的英文文献的粗略翻译与文章的程序

A Novel Approach for Image Enhancement by Using Contrast Limited Adaptive Histogram Equalization Method（一种通过运用对比度受限自适应直方图均衡方法的图像增强方法）

摘要

一个图像增强的新方法通过运用对比度受限的自适应直方图均衡方法将铲射给一个好的对比度图像如医学图像。在本文中，我们提出了一个新的图像增强方式通过运用对比度首先的自适应直方图均衡的方法。我们提出了一个自适应直方图均衡方法总体的框架。我们将通过对比其他的对比度增强方法来证明本方法的有效性。索引词汇：、、、、、

介绍

图像增强的目的是产生一个重要的步骤在图像处理中，这个步骤是通过编辑原始的图像来观测到更多的图像对比度在一个特定的应用领域中。对比度增强技术将会在图像处理应用中扮演重要的角色，

如移动图像，数字照片，分析医疗图像，遥感，和各种图形科学中。所有的图像会由几点因素产生不好的对比度，因为运用不好质量的成像设备或者光线和环境所以我们的提出了一个新的对比度增强技术意在消除这些问题。不同的对比度增强技术被用在提高图像的对比度如直方图均衡，直方图修正，贪婪算法，自适应直方图均衡等等。本文提出了一个新的一个基于对比度受限自适应直方图均衡的直方图均衡方法。

退出对比度增强方法

现在有很多的对比度增强技术可使用在不清晰图像上，所以我们先讨论下面几中对比度增强技术

直方图均衡

直方图均衡的目的是分配一个图像的灰度级，从而使没给灰度级的发生概率近乎相等。直方图均衡将会增加亮度和一幅暗的和低对比度图像的对比度。使原始图像中观察到的特征不可见。它还用来规范图像的亮度和对比度，直方图均衡进程是找到一个映射函数来映射，图像输入直方图函数和使分布式输出直方图函数相一致。直方图均衡化还被用在生物神经网络上来是神经元输出功率达到最大做为输入统计函数。这已经在fly retina 中被证明了。直方图均衡化是在更一般的直方图映射方法中较为特别的。这些方法是通过调整图像来时期更容易被分析和增强其显示质量。

B) 直方图平滑

为了避免尖峰带来的强烈的排斥固定点问题一个平滑度约束可

以完成这个目标。直方图的向后方差被用来确保其平滑度。平滑可以修改直方图是其尖峰更少因为他们在直方图中基本上都是突然变化的。

C) 自适应直方图均衡

普通的直方图均衡运用的是从图像直方图之得到的一样的变换源来去变换所有得像素。这个在像素值分布相似的图像中可以工作的很好。然而，当图像包含区域比大部分图像部分更暗或更亮时，这些区域的对比度将得不到足够的增强。自适应直方图均衡可以通过运用来自于临近区域的变换函数来变换每一个像素值来改善这一点。它的首次开发是运用在飞机驾驶舱显示器中。当图像区域包含的一个像素的临域是相当均匀的，它的直方图将会特别的尖锐，转移函数将只会映射整个结果图像中的一个狭窄域的像素值。这导致了自适应直方图均衡过于放大了很大程度上同质的的图像上的少量信号。这个方法被用来增强图像的对比度。它通过自适应的方式改变了直方图均衡化而这种方法是通过计算几个直方图来完成的，每个直方图对应着一个图像的不同区域，从而重新分配图像的亮度值。因此它便于增减图像的局部对比度和显示更多的细节。

3提出的算法

（a 对比度受限的自适应直方图均衡

一个提出算法为了医学图像特别的开发出来它提供了很好的图像增强相比于原始图像。CLAHE 算法把图像分割成些相关区域然后对没个区域进行直方图均衡化。这些事件产生了使用过的灰度值并重

新分配因此可以使图像隐藏的特征更加的可见。CLAHE 是AHE 的改进算法。我们运用上述提到的算法来去增强图像，直方图均衡化，直方图平滑，自适应直方图均衡化和加强了的对比度首先的自适应直方图均衡化。这些提到的加强技术馋了的如图一的结果，第一次测试结果的可视化表示（乳房肿瘤）。从观察到的可视化结果分析所有的算法都能够加强其对比度但是我们提到的算法中对比度受限的自适应直方图可以更准切的显示出图像。

（3

对于图像质量来说人类的视觉不被认为是基准，所以估计上述提到的算法的质量指标已经被计算出来了通过对比原始图像和输出图像。下图2表示了直方图增强图像等级的映射。

为了评价上述提到的图像技术增强方法产生了更好的增强图像。表1揭示了CLAHE 对比与其他算法的效果。得出的结果再次更好的揭示了这些增强算法的优劣。其它的算法也产生了不错的图像增强，但是还是没有CLAHE 好。

算法步骤

获得输入：制定行和列的数量然后设定限幅界限和用于直方图的统计特征数量（bin ）用于构建对比度增强的转换。更高的数值引起更大的动态范围代价是处理速度变慢。对比度增强技术的限幅界限是从0到1来限制对比度增强。更高的数值可以得到更高的对比度。

处理输入：制定限幅界限为归一化值如果必要的话，在图像分裂成区域前补充图像（ pad the imagebefore splitting it into regions）

对比每一个过程行和列区域（瓦）从而产生灰度映射值做一个此区域的直方图并用制定的数量的统计特征数量（bin ）用限幅界限来限制此直方图，对比度限制程序（由变换公式推出）必须使用与每一个临近区。CLAHE 为了防止噪声过度的被放大自适应直方图均衡可以解决它。

差值允许效率的显著提高而不影响结果的质量

上图4代表的是运用对比度受限自适应直方图均衡时不同限幅界限所产生的图。CLAHE 的限幅界限归一化从0到1. 更高的数值产生更大的对比度。表1.1代表了不同限幅界限时对比度水平。

四结论

在此文中，对比度受限的自适应直方图均衡已经运用于乳房癌症图像的对比度增强。对比与现在已经提出的流行的的方法直方图均衡，直方图平滑，自适应直方图均衡他它已经被证可以得到更好的结果了。

下面是仿真程序和参考程序

直方图均衡：

i1=imread('f:\1\1.jpg');

i2=rgb2gray(i1);

figure;

subplot(131);imshow(i1);title(' 原图' );

subplot(132);imshow(i2);title(' 灰度图' );

subplot(133);imhist(i2);title(' 原图直方图' );

i3=histeq(i2);

figure;

subplot(221);imshow(i2);title(' ａ灰度图' );

subplot(222);imshow(i3);title('b 直方图均衡后的图' );

subplot(223);imhist(i2);title('c 原图直方图' );

subplot(224);imhist(i3);title('d 直方图均衡后直方图' );

%zhifangtuqunhenghua

clear;

x1=imread('1.tif' );

[m,n]=size(x1);

x2=double(x1);

lenna=zeros(m,n);

lenna_equ=zeros(n,n);

histgram=zeros(256);

cdf=zeros(256);

[lenna,map]=imread('1.tif' );

%get histogram

for i=1:m

for j=1:n

k=lenna(i,j);

histgram(k)=histgram(k)+1;

end

%get cdf

cdf(1)=histgram(1);

for i=2:256

cdf(i)=cdf(i-1)+histgram(i);

end

%run point operation

for i=1:m

for j =1:n

k=lenna(i,j);

lenna_equ(i,j)=cdf(k)*256/(m*n);

end

%生成直方图均衡化后的lenna 图

g=mat2gray(lenna_equ);

figure;

subplot(2,2,1);imshow(lenna);

subplot(2,2,2);imhist(lenna);

ylim('auto' );

subplot(2,2,3);imshow(g);

subplot(2,2,4);imhist(g);

ylim('auto' )

f=imread('1.tif' );

figture;

subplot(2,2,1);imshow(f);

title('the orignal image');

subplot(2,2,2);imhist(f);

ylim('auto' );

g=histeq(f,256);

subplot(2,2,3);imshow(g);

title('image after equalization')

subplot(2,2,4);imhist(g);

ylim('auto' );

%zhifangtupipeichengxu

I = imread('pout.tif' );

Ieq=histeq(I);

[c,x]=imhist(Ieq);

I1=imread('tire.tif' );

J=histeq(I1,c);

subplot(2,2,1);

imshow(I);

title('Image A');

subplot(2,2,2);

imshow(Ieq);

title('Image A after Equalization');

subplot(2,2,3);

imshow(I1);

title('Image B');

subplot(2,2,4);

imshow(J);

title('Image B after Matching');

受限的自适应直方图均衡

oimg=imread('f:\1\1.jpg');

oimg1=rgb2gray(oimg);

i3=adapthisteq(oimg1);

i4=histeq(oimg1);

figure;

subplot(131);imshow(oimg1);title(' ａ灰度图' );

subplot(132);imshow(i4);title('b 直方图均衡' );

subplot(133);imshow(i3);title('c 限制对比度的直方图均衡' );

function [CEImage] = runCLAHE(Image,Min,Max,NrX,NrY,NrBins,Cliplimit)

% "Contrast Limited Adaptive Histogram Equalization"

% by Karel Zuiderveld

% The program was reproduced by Alireza Saberi in April 2010

% These functions implement Contrast Limited Adaptive Histogram Equalization.

% The main routine (CLAHE) expects an input image that is stored contiguously in

% memory; the CLAHE output image overwrites the original input image and has the

% same minimum and maximum values (which must be provided by the user). % This implementation assumes that the X- and Y image resolutions are an integer

% multiple of the X- and Y sizes of the contextual regions. A check on various other

% error conditions is performed.

% Image - The input/output image

% XRes - Image resolution in the X direction

% YRes - Image resolution in the Y direction

% Min - Minimum greyvalue of input image (also becomes minimum of output image)

% Max - Maximum greyvalue of input image (also becomes maximum of output image)

% NrX - Number of contextial regions in the X direction (min 2, max uiMAX_REG_X)

% NrY - Number of contextial regions in the Y direction (min 2, max uiMAX_REG_Y)

% NrBins - Number of greybins for histogram ("dynamic range")

% Cliplimit - Normalized cliplimit (higher values give more contrast)

% The number of "effective" greylevels in the output image is set by

uiNrBins; selecting

% a small value (eg. 128) speeds up processing and still produce an output image of

% good quality. The output image will have the same minimum and maximum value as the input

% image. A clip limit smaller than 1 results in standard (non-contrast limited) AHE.

Image=imread('f:\1\1.jpg');

[XRes,YRes]=size(Image);

% CEimage = Image;

CEImage = zeros(XRes,YRes);

if Cliplimit == 1

return

end

NrBins=max(NrBins,128);

XSize = round(XRes/NrX);

YSize = round(YRes/NrY);

NrPixels = XSize*YSize;

XSize2 = round(XSize/2);

YSize2 = round(YSize/2);

if Cliplimit > 0

ClipLimit = max(1,Cliplimit*XSize*YSize/NrBins);

else

ClipLimit = 1E8;

end

LUT=makeLUT(Min,Max,NrBins);

% avgBin = NrPixels/NrBins;

Bin=1+LUT(round(Image));

Hist = makeHistogram(Bin,XSize,YSize,NrX,NrY,NrBins);

if Cliplimit > 0

Hist = clipHistogram(Hist,NrBins,ClipLimit,NrX,NrY);

end

Map=mapHistogram(Hist,Min,Max,NrBins,NrPixels,NrX,NrY);

% Interpolate

xI = 1;

for i = 1:NrX+1

if i == 1

subX = XSize/2;

xU = 1;

xB = 1;

elseif i == NrX+1

subX = XSize/2;

xU = NrX;

xB = NrX;

else

subX = XSize;

xU = i - 1;

xB = i;

end

yI = 1;

for j = 1:NrY+1

if j == 1

subY = YSize/2;

yL = 1;

yR = 1;

elseif j == NrY+1

subY = YSize/2;

yL = NrY;

yR = NrY;

else

subY = YSize;

yL = j - 1;

yR = j;

end

UL = Map(xU,yL,:);

UR = Map(xU,yR,:);

BL = Map(xB,yL,:);

BR = Map(xB,yR,:);

subImage = Bin(xI:xI+subX-1,yI:yI+subY-1);

subImage = interpolate(subImage,UL,UR,BL,BR,subX,subY); CEImage(xI:xI+subX-1,yI:yI+subY-1) = subImage;

yI = yI + subY;

end

xI = xI + subX;

end

基于直方图的图像增强方法姓名：

学号：导师：

基于直方图的图像增强方法

一．直方图均衡化算法

1.1算法介绍与仿真

图像的直方图是图像的重要统计特征，是用来表征数字图像的每一灰度级与其出现的频率间的统计关系的方法。其数学公式如式（1）所示：

n k P (r k ) =(k =0, 1, 2,.... L -1) (1) N

其计算步骤如下：

1）列出原始图像的灰度级r k ，k ＝0，1，2，…，L －1，L 为总的灰度级数；

2）统计各灰度级像素数目n k ，k ＝0，1，2，…，L －1；

3）计算原始图像直方图各灰度级的频率P （r k ）＝n k ／N ，k ＝0， 1，2，…，L －1，N 为图像f （x ，y ）的总像素数；

4）计算累计分布函数C (k ) =∑P (r k ) , k =0, 1, 2,... L -1;

0L -1

5）计算映射后的输出图像灰度级g （i ）＝INT ［（gmax －gmin ）C （k ）＋gmin ＋0．5］，i ＝0，1， …，P －1，P 为输出图像灰度等级个数，

INT 为取整符号；

6）统计映射后各灰度等级的像素数目n i ，i ＝0，1，…，P －1；

7）计算输出图像直方图P （g i ）＝n i ／N ，i ＝0，1，…，P －1；

8）用r k 和g i 的映射关系修改原始图像的灰度级，获得均衡化后的输出图像。

图1

1.2算法存在的问题

由于直方图均衡化理论来源于连续函数而数字图像的的灰度是离散值应用于数字图像处理的变换函数进行了从连续到离散的近似。如积

分运算就变成了累加运算因而在直方图均衡存在以下问题。

1）量化误差

2）无法增强局部细节

二．自适应直方图均衡化

那么该如何提取图像的局部特征呢? 采用局部法对图像进行处

x i , j = m i , j + k ( x i , j - m i , j ) （2 ）

改变局部窗口内的对比度, 同样达到了细节增强的效果, 优化了图像质量. 下面我们详细地分析局部对比度增强法的具体算法。

x i , j 0≤x i , j ≤255⎧⎪T (x i , j +k (x i , j -m i , j )) 　　　=⎨ （3） T (x ) 　　　　　　　　　其它⎪i , j ⎩

自适应参量的选取上。

考虑到参量k 的自适应性, 我们先对k 进行分析。当窗W 的中心象素x i, j位于无图像细节处时，k 趋向于0，而当窗W 的中心象素x i,

k j 位于图像细节处时，取一个较大正值。基于以上考虑我们选取窗W

内邻域灰度方差作为自适应变量, 那么的表达式可以写为:

⎡σ2i , j ⎤k =k () -1 ⎢ 2n ⎥ ⎣⎦

其中k 为比例系数，σ2n 为整幅图像的噪声方差，σ

度方差。若令r =

k 值也变大，细节被增强。

综上所述, 我们得出自适应直方图均衡算法具体的实现过程如下。

1) 给出原始图像的所有灰度级i, 统计原始图像各灰度级的像素数ni

2) 计算原始图像的直方图与累积直方图。

3) 用式（4）计算k 值. 用式（3）进行计算局部灰度值。

4) 用式（2）计算局部对比度, 实现均衡化。

5) 用p (ti ) =

n i 计算新的直方图。 n

三．对比度受限的自适应直方图均衡

AHE 有过度放大图像中相同区域的噪音的问题，另外一种自适应的直方图均衡算法即限制对比度直方图均衡（CLAHE ）算法能有限的限制这种不利的放大。

3.1算法简介

对比度受限自适应直方图均衡法（CLAHE ）通过限制局部直方图的高度来限制局部对比度的增强幅度，从而限制噪声的放大及局部对比度的过增强。

对比度增强的幅度可定义为灰阶映射函数的斜率。假定自适应直方图均衡方法的滑动窗口的大小为M ×M ，则局部映射函数为：

255⨯CDF (i ) m (i ) = (5) M ⨯M

CDF(i)为滑动窗口局部直方图的累积分布函数。累积分布函数CDF(i)的导数为直方图Hist(i)，从而局部映射函数m(i)的斜率S 为：

d (m (i ) )255S ==Hist (i ) ⨯ (6) di M ⨯M

因此限制直方图高度就等效于限制映射函数m(i)的斜率S ，进而限制对比增强度。若限定最大斜率为Smax ，则允许的最大直方图高度为：

H max =S max ⨯M ⨯M (7) 255

方图面积不变，从而使整个直方图上升一个高度L 。因此H max 、T 、L 三者之间应满足下面关系：

H max =T+L

(8)

最后，改进的直方图值为：

⎧Hist (i ) +L 　　　　Hist (i )

通过改变最大的映射函数斜率S max 及相应的最大直方图高度H max ，便可以获得不同增强效果的图像。

1) 图像分块

2) 产生局部直方图

3) 分别对每个图像块进行限制对比度直方图均衡，每个图像块分别产生独立的变换函数

4) 将多个图像块的灰度查找表拼接成输出图像

其基本思想是构造限制函数，限制灰度级的概率密度，并将超过限制函数的像素点在直方图内进行重整。

具体算法：

1）定义限制函数clipLimit ，计算超过限制函数的总像素点数totalExcess

2）计算totalExcess 平均分到每个灰度级的平均值avgBinIncr

3）定义限制函数与平均值的差upperLimit = clipLimit - avgBinIncr;

4）对直方图进行限高处理

● if 第K 灰度级像素数>限制函数

●降低该灰度级像素数=限制函数

●else

●if 第K 灰度级像素数> upperLimit

●增高该灰度级概率密度=限制函数；

●并在超过限制函数的总像素点数totalExcess 中减去增加的像素点数

●Else 第K 灰度级像素数

●该灰度级像素点数+多余像素点平均值avgBinIncr ；

●并在超过限制函数的总像素点数totalExcess 中减去该平均值。

5）对直方图进行均衡

3.2仿真结果

图2

四．总结

五．附录

下面这部分是所涉及到的英文文献的粗略翻译与文章的程序

A Novel Approach for Image Enhancement by Using Contrast Limited Adaptive Histogram Equalization Method（一种通过运用对比度受限自适应直方图均衡方法的图像增强方法）

摘要

介绍

退出对比度增强方法

现在有很多的对比度增强技术可使用在不清晰图像上，所以我们先讨论下面几中对比度增强技术

直方图均衡

B) 直方图平滑

为了避免尖峰带来的强烈的排斥固定点问题一个平滑度约束可

以完成这个目标。直方图的向后方差被用来确保其平滑度。平滑可以修改直方图是其尖峰更少因为他们在直方图中基本上都是突然变化的。

C) 自适应直方图均衡

3提出的算法

（a 对比度受限的自适应直方图均衡

（3

算法步骤

处理输入：制定限幅界限为归一化值如果必要的话，在图像分裂成区域前补充图像（ pad the imagebefore splitting it into regions）

差值允许效率的显著提高而不影响结果的质量

四结论

下面是仿真程序和参考程序

直方图均衡：

i1=imread('f:\1\1.jpg');

i2=rgb2gray(i1);

figure;

subplot(131);imshow(i1);title(' 原图' );

subplot(132);imshow(i2);title(' 灰度图' );

subplot(133);imhist(i2);title(' 原图直方图' );

i3=histeq(i2);

figure;

subplot(221);imshow(i2);title(' ａ灰度图' );

subplot(222);imshow(i3);title('b 直方图均衡后的图' );

subplot(223);imhist(i2);title('c 原图直方图' );

subplot(224);imhist(i3);title('d 直方图均衡后直方图' );

%zhifangtuqunhenghua

clear;

x1=imread('1.tif' );

[m,n]=size(x1);

x2=double(x1);

lenna=zeros(m,n);

lenna_equ=zeros(n,n);

histgram=zeros(256);

cdf=zeros(256);

[lenna,map]=imread('1.tif' );

%get histogram

for i=1:m

for j=1:n

k=lenna(i,j);

histgram(k)=histgram(k)+1;

end

%get cdf

cdf(1)=histgram(1);

for i=2:256

cdf(i)=cdf(i-1)+histgram(i);

end

%run point operation

for i=1:m

for j =1:n

k=lenna(i,j);

lenna_equ(i,j)=cdf(k)*256/(m*n);

end

%生成直方图均衡化后的lenna 图

g=mat2gray(lenna_equ);

figure;

subplot(2,2,1);imshow(lenna);

subplot(2,2,2);imhist(lenna);

ylim('auto' );

subplot(2,2,3);imshow(g);

subplot(2,2,4);imhist(g);

ylim('auto' )

f=imread('1.tif' );

figture;

subplot(2,2,1);imshow(f);

title('the orignal image');

subplot(2,2,2);imhist(f);

ylim('auto' );

g=histeq(f,256);

subplot(2,2,3);imshow(g);

title('image after equalization')

subplot(2,2,4);imhist(g);

ylim('auto' );

%zhifangtupipeichengxu

I = imread('pout.tif' );

Ieq=histeq(I);

[c,x]=imhist(Ieq);

I1=imread('tire.tif' );

J=histeq(I1,c);

subplot(2,2,1);

imshow(I);

title('Image A');

subplot(2,2,2);

imshow(Ieq);

title('Image A after Equalization');

subplot(2,2,3);

imshow(I1);

title('Image B');

subplot(2,2,4);

imshow(J);

title('Image B after Matching');

受限的自适应直方图均衡

oimg=imread('f:\1\1.jpg');

oimg1=rgb2gray(oimg);

i3=adapthisteq(oimg1);

i4=histeq(oimg1);

figure;

subplot(131);imshow(oimg1);title(' ａ灰度图' );

subplot(132);imshow(i4);title('b 直方图均衡' );

subplot(133);imshow(i3);title('c 限制对比度的直方图均衡' );

function [CEImage] = runCLAHE(Image,Min,Max,NrX,NrY,NrBins,Cliplimit)

% "Contrast Limited Adaptive Histogram Equalization"

% by Karel Zuiderveld

% The program was reproduced by Alireza Saberi in April 2010

% These functions implement Contrast Limited Adaptive Histogram Equalization.

% The main routine (CLAHE) expects an input image that is stored contiguously in

% memory; the CLAHE output image overwrites the original input image and has the

% same minimum and maximum values (which must be provided by the user). % This implementation assumes that the X- and Y image resolutions are an integer

% multiple of the X- and Y sizes of the contextual regions. A check on various other

% error conditions is performed.

% Image - The input/output image

% XRes - Image resolution in the X direction

% YRes - Image resolution in the Y direction

% Min - Minimum greyvalue of input image (also becomes minimum of output image)

% Max - Maximum greyvalue of input image (also becomes maximum of output image)

% NrX - Number of contextial regions in the X direction (min 2, max uiMAX_REG_X)

% NrY - Number of contextial regions in the Y direction (min 2, max uiMAX_REG_Y)

% NrBins - Number of greybins for histogram ("dynamic range")

% Cliplimit - Normalized cliplimit (higher values give more contrast)

% The number of "effective" greylevels in the output image is set by

uiNrBins; selecting

% a small value (eg. 128) speeds up processing and still produce an output image of

% good quality. The output image will have the same minimum and maximum value as the input

% image. A clip limit smaller than 1 results in standard (non-contrast limited) AHE.

Image=imread('f:\1\1.jpg');

[XRes,YRes]=size(Image);

% CEimage = Image;

CEImage = zeros(XRes,YRes);

if Cliplimit == 1

return

end

NrBins=max(NrBins,128);

XSize = round(XRes/NrX);

YSize = round(YRes/NrY);

NrPixels = XSize*YSize;

XSize2 = round(XSize/2);

YSize2 = round(YSize/2);

if Cliplimit > 0

ClipLimit = max(1,Cliplimit*XSize*YSize/NrBins);

else

ClipLimit = 1E8;

end

LUT=makeLUT(Min,Max,NrBins);

% avgBin = NrPixels/NrBins;

Bin=1+LUT(round(Image));

Hist = makeHistogram(Bin,XSize,YSize,NrX,NrY,NrBins);

if Cliplimit > 0

Hist = clipHistogram(Hist,NrBins,ClipLimit,NrX,NrY);

end

Map=mapHistogram(Hist,Min,Max,NrBins,NrPixels,NrX,NrY);

% Interpolate

xI = 1;

for i = 1:NrX+1

if i == 1

subX = XSize/2;

xU = 1;

xB = 1;

elseif i == NrX+1

subX = XSize/2;

xU = NrX;

xB = NrX;

else

subX = XSize;

xU = i - 1;

xB = i;

end

yI = 1;

for j = 1:NrY+1

if j == 1

subY = YSize/2;

yL = 1;

yR = 1;

elseif j == NrY+1

subY = YSize/2;

yL = NrY;

yR = NrY;

else

subY = YSize;

yL = j - 1;

yR = j;

end

UL = Map(xU,yL,:);

UR = Map(xU,yR,:);

BL = Map(xB,yL,:);

BR = Map(xB,yR,:);

subImage = Bin(xI:xI+subX-1,yI:yI+subY-1);

subImage = interpolate(subImage,UL,UR,BL,BR,subX,subY); CEImage(xI:xI+subX-1,yI:yI+subY-1) = subImage;

yI = yI + subY;

end

xI = xI + subX;

end

数字图像处理作业

相关文章