第5章〓图像增强
本章思维导图







图像增强(Image Enhancement)主要是为了改善图像的质量以及增强感兴趣部分，改善图像的视觉效果或使图像变得更利于计算机处理。相关的图像增强技术有针对单个像素点的点运算，也有针对像素局部邻域的模板运算，根据模板运算的具体功能还可以分为图像平滑、图像锐化等。本章主要讲解图像增强技术中的灰度级变换、直方图修正法、基于照度反射模型的图像增强和色彩变换等。

5.1灰度级变换


灰度级变换就是借助于变换函数将输入的像素灰度值映射成一个新的输出值，通过改变像素的亮度值来增强图像，如式(51)所示。


g(x,y)=T［f(x,y)］(51)
其中，f(x,y)是输入图像，g(x,y)是变换后的输出图像，T是灰度变换函数。

变换函数T的不同将导致不同的输出，其实现的变换效果也不一样。因此，在实际应用中，可以通过灵活地设计变换函数T来实现各种处理。

根据变换函数的不同，灰度级变换可以分为线性灰度级变换和非线性灰度级变换。

5.1.1线性灰度级变换
1. 基本线性灰度级变换

最基本的线性灰度级变换如图51所示。将输入图像fx,y变换为gx,y，变换结果由倾角α所决定，如式(52)所示。

gx,y=fx,y·tanα(52)


可以看出，当α=45°，图像无变化； 当α<45°，变换后灰度取值范围压缩，图像均匀变暗； 当α>45°，变换后灰度取值范围拉伸变长，图像均匀变亮。

当线性灰度级变换如图52所示，将图像高低灰度值发生翻转，实现对图像求反，简单说来就是使黑变白，使白变黑。

2．  分段线性灰度级变换

分段线性灰度级变换如图53所示，它是将输入图像fx,y的灰度级区间分成两段乃至多段，分别作线性灰度级变换，以获得增强图像gx,y。典型的三段线性灰度级变换如式(53)所示。

gx,y=cafx,y，0≤fx,y<a
d－cb－afx,y－a+c，a≤fx,y<b
L－1－dL－1－bfx,y－b+d，b≤fx,y<L－1(53)






图51基本线性灰度级变换



图52图像求反



图53分段线性灰度级变换


其中，参数a、b、c、d为用于确定三段线段斜率的常数，取值可根据具体变换需求灵活设定。但也存在一些情况，用户仅对某个范围内的灰度感兴趣，只需对其进行线性拉伸，而将其余灰度级保持不变或变为某固定值。

MATLAB也提供了函数imadjust用于实现图像的灰度级映射变换。

【例5．1】编写程序，对图像进行线性灰度级变换。

解： 程序如下。

Image=im2double(rgb2gray(imread('img5_1.bmp')));%读取图像，灰度化并转换为double型

［h,w］=size(Image);   %获取图像尺寸

NewImage1=tand(60).*Image;   %基本线性灰度级变换

NewImage2=1-Image;   %图像求反 

NewImage4=Image; 

a=30/256; b=100/256; c=75/256; d=200/256;  %参数设置         

for x=1:w

for y=1:h

if Image(y,x)<a 

NewImage3(y,x)=Image(y,x)*c/a;

elseif Image(y,x)<b

NewImage3(y,x)=(Image(y,x)-a)*(d-c)/(b-a)+c;  %三段线性灰度级变换

else

NewImage3(y,x)=(Image(y,x)-b)*(1-d)/(1-b)+d;

end  

if Image(y,x)>a && Image(y,x)<b

NewImage4(y,x)=(Image(y,x)-a)*(d-c)/(b-a)+c; %%［a,b］内线性变换，其余值不变

end        

end

end

NewImage5=imadjust(Image,［a b］,［c d］);  %［a,b］内线性变换，其余低端值变为c，高端值变为d

NewImage6=imadjust(Image,［a b］,［］);  %将［a,b］线性拉伸到［0,1］，进行对比度拉伸

imshow(Image);title('原始图像');

figure,

subplot(2,3,1),imshow(NewImage1);title('基本线性灰度级变换');

subplot(2,3,2),imshow(NewImage2);title('图像取反');

subplot(2,3,3),imshow(NewImage3);title('分段线性灰度级变换');

subplot(2,3,4),imshow(NewImage4);title('高低端灰度保持不变');

subplot(2,3,5),imshow(NewImage5);title('高低端灰度赋固定值');

subplot(2,3,6),imshow(NewImage6);title('［a,b］线性拉伸到［0,1］');

线性灰度级变换的处理结果如图54所示。



图54线性灰度级变换效果


5.1.2非线性灰度级变换

当用某些非线性变换函数作为灰度变换的变换函数时，可实现图像灰度的非线性变换。对数变换、指数变换和幂变换是常见的非线性变换。

1. 对数变换

基于对数变换的非线性灰度级变换如式(54)所示。


g(x,y)=c·log［f(x,y)+1］(54)
式中，c是尺度比例常数，其取值可以结合输入图像的范围来定。

对数变换函数示意图如图55所示。当希望对图像的低灰度区作较大拉伸、高灰度区压缩时，可采用这种变换，它能使图像的灰度分布与人的视觉特性相匹配。对数变换一般适用于处理过暗图像。 


2. 指数变换

基于指数变换的非线性灰度变换如式(55)所示。


g(x,y)=bc·［f(x,y)－a］－1(55)

式中，a用于决定指数变换函数曲线的初始位置。当取f(x,y)=a时，g(x,y)=0，曲线与x轴交叉。b是底数，c用于决定指数变换曲线的陡度。

指数变换函数示意图如图56所示。当希望对图像的低灰度区压缩，高灰度区作较大拉伸时，可采用这种变换。指数变换一般适用于处理过亮图像。



图55对数变换函数




图56指数变换函数





图57幂次变换函数

3. 幂次变换

基于幂次变换的非线性灰度变换如式(56)所示。


g(x,y)=c·［f(x,y)］γ(56)
式中，c和γ为正常数。当c取1，γ取不同值时，可以得到一簇变换曲线，如图57所示。


与对数变换的情况类似，幂次变换可以将部分灰度区域映射到更宽的区域中，从而增强图像的对比度。当γ=1时，幂次变换转变为线性正比变换； 当0<γ<1时，幂次变换可以扩展原图中低灰度级、压缩高灰度级，从而使图像变亮，增强原图中暗区的细节； 当γ>1时，幂次变换可以扩展原图中高灰度级，压缩低灰度级，从而使图像变暗，增强原图中亮区的细节。

幂次变换常用于图像获取、打印和显示的各种装置设备的伽马校正中，这些装置设备的光电转换特性都是非线性的，是根据幂次规律产生响应的。幂次变换的指数值就是伽马值，因此幂次变换也称为伽马变换。

【例5.2】编写程序对图像进行非线性灰度变换。

解： 程序如下。

Image=rgb2gray(imread('img5_2.bmp'));  

Image=double(Image);

NewImage1=log(Image+1);                          %对数函数非线性灰度级变换

NewImage2=exp(0.325*(Image-225)/30)+1;       %指数函数非线性灰度级变换

a=0.5; c=1.1;

NewImage3=［(Image/255).^a］*255*c;                 %幂次函数非线性灰度级变换

imshow(Image,［］);title('原图');

figure;imshow(NewImage1,［］);title('对数变换');

figure;imshow(NewImage2,［］);title('指数变换');

figure;imshow(NewImage3,［］);title('幂次变换');

程序运行结果如图58所示。可以看出，原图经过对数变换处理后（如图58(b)），整体亮度增强； 经过指数变换处理后（如图58（c）），整体亮度变暗； 图58(d)为幂次变换处理结果，幂为0.5，图像变亮，但变亮程度低于对数变换。



图58非线性灰度变换处理效果


5.2直方图修正法

在数字图像处理中，灰度直方图是最简单和常用的工具。本节介绍灰度直方图的概念及直方图修正技术。

5.2.1灰度直方图
1. 灰度直方图的定义

灰度直方图是灰度级的函数，表示的是数字图像中每一灰度级与其出现频数(呈现该灰度的像素数目)间的统计关系。通常，用横坐标表示灰度级，纵坐标表示频数或相对频数(呈现该灰度级的像素出现的概率)。灰度直方图的定义如式(57)所示。


p(rk)=nkN(57)

式中，N为一幅数字图像的总像素数，nk是第k级灰度的像素数，rk表示第k个灰度级，p(rk)为该灰度级rk出现的相对频数。

给定一幅6×6的图像，如图59(a)所示，共有0~7八个灰度级。灰度分布统计如表51所示，则可以绘制并显示图像的灰度直方图，如图59(b)所示。


表51灰度分布统计



rk01234567
nk69654321
p(rk)6/369/366/365/364/363/362/361/36





图59数字图像及其直方图


【例5.3】编写程序，统计图像的灰度直方图分布。

解： 统计图像的灰度直方图可以通过扫描图像并统计各个灰度出现的次数，进而计算频数并绘制直方图的方法。

程序如下。

Image=rgb2gray(imread('img5_3.bmp'));

histgram=zeros(256);                                    

［h w］=size(Image);

for x=1:w

for y=1:h                                                     %循环扫描

histgram(Image(y,x)+1)=histgram(Image(y,x)+1)+1;   %统计并累加

end

end

imshow(Image);title('couple灰度图像');

figure;stem(histgram(),'.'); 

axis tight;

也可以直接调用MATLAB提供的函数imhist,程序如下。 

figure;imhist(Image);

axis tight;

程序运行结果如图510所示。



图510原图及其直方图


2. 灰度直方图的性质

一幅图像的灰度直方图通常具有如下性质。

(1) 直方图不具有空间特性。直方图描述了每个灰度级具有的像素的个数，但不能反映图像像素空间位置信息，即不能为这些像素在图像中的位置提供任何线索。

(2) 直方图反映图像的大致描述，如图像灰度范围、灰度级分布、整幅图像平均亮度等。
图511为两幅图像的直方图，可以从中判断出图像的相关特性。图511(a)中，大部分像素值集中在低灰度级区域，图像偏暗； 图511(b)中的图像则相反，大部分像素的灰度集中在高灰度级区域，图像偏亮； 两幅图像都存在动态范围不足的现象。



图511灰度动态范围不足的图像灰度直方图


(3) 一幅图像唯一对应相应的直方图，而不同的图像可以具有相同的直方图。
因直方图只是统计图像中灰度出现的次数，与各个灰度出现的位置无关，因此，不同的图像可能具有相同的直方图。图512(a)为四幅大小相同、空间灰度分布不同的二值图像，图512(b)为它们具有的相同的直方图。



图512不同图像具有相同直方图分布特性


(4) 若一幅图像可分为多个子区，则多个子区直方图之和等于对应的全图直方图。

5.2.2直方图修正法理论

直方图修正法的基本原理就是通过构造灰度级变换函数改造原图像的直方图，使变换后的图像的直方图达到一定的要求。

设变量r代表要增强的图像中像素的灰度级，变量s代表增强后新图像中的灰度级。为了研究方便，将r、s归一化，得


0≤r≤1，0≤s≤1(58)
则直方图修正法变换函数T(·)的定义如式(59)所示。


s=T(r)(59)
其中，每一像素灰度值r对应产生一个s值。

对图像的直方图修正变换过程中，要满足： 

(1) T(r)在0≤r≤1区域内单值单调增加，以保证灰度级从黑到白的次序不变； 

(2)  T(r)在0≤r≤1区域内满足0≤s≤1，以保证变换后的像素灰度级仍在允许的灰度级范围内。

直方图修正法的核心就是寻找满足这两个条件的变换函数T(r)。

5.2.3直方图均衡化

直方图均衡化是采用灰度级r的累积分布函数作为变换函数的直方图修正法。

假设用pr(r)表示原图像灰度级r的灰度级概率密度函数，直方图均衡化变换函数为


s=T(r)=∫r0pr(ω)dω(510)

T(r)是r的累积分布函数，随着r增大，s值单调增加，最大为1，满足直方图修正法的两个条件。

根据概率论知识，用pr(r)和ps(s)分别表示r和s的灰度级概率密度函数，得



ps(s)=pr(r)·drds=pr(r)·1pr(r)=1(511)
即利用r的累积分布函数作为变换函数可产生一幅灰度级分布具有均匀概率密度的图像。

给出一幅数字图像，共有L个灰度等级，总像素个数为N，其中，第j级灰度rj对应的像素数为nj。根据式(57)进行灰度直方图统计，则图像进行直方图均衡化处理的变换函数T(r)为


sk=T(rk)=∑kj=0pr(rj)=∑kj=0njN(512)

对一幅数字图像进行直方图均衡化处理的算法步骤如下。

(1) 由式(57)统计原始图像直方图。

(2) 由式(512)计算新的灰度级。

(3) 修正sk为合理的灰度级。

(4) 计算新的直方图。

(5) 用处理后的新灰度代替处理前的灰度，生成新图像。

【例5.4】假定一幅大小为64×64、灰度级为8级的图像，其灰度级分布如表52所示，对其进行直方图均衡化处理。


表52例5.4中图像的灰度级分布



灰度级rk01/72/73/74/75/76/71
像素数nk790102385065632924512281
pr(rk)0.190.250.210.160.080.060.030.02


解： 由原图的灰度分布统计可以看出，该图像中绝大部分像素灰度值集中在低灰度区，图像整体偏暗。

(1) 计算新的灰度级： 



s0=T(r0)=∑0j=0pr(rj)=pr(r0)=0.19
s1=T(r1)=∑1j=0pr(rj)=pr(r0)+pr(r1)=0.19+0.25=0.44
以此类推，可得到


s2=0.19+0.25+0.21=0.65
s3=0.19+0.25+0.21+0.16=0.81，s4=0.89
s5=0.95，s6=0.98，s7=1

(2) 修正sk为合理的灰度级s′k： 



s0=0.19≈17，s1=0.44≈37，s2=0.65≈57，s3=0.81≈67
s4=0.89≈67，s5=0.95≈1，s6=0.98≈1，s7=1

则新图像对应只有5个不同灰度级别，为1/7、3/7、5/7、6/7、1，即


s′0=17，s′1=37，s′2=57，s′3=67，s′4=1

(3) 计算新的直方图： 



ps(s′0)=pr(r0)=0.19
ps(s′1)=pr(r1)=0.25
ps(s′2)=pr(r2)=0.21
ps(s′3)=pr(r3)+pr(r4)=0.16+0.08=0.24
ps(s′4)=pr(r5)+pr(r6)+pr(r7)=0.06+0.03+0.02=0.11

(4) 生成新图像：
 按照表53中变换前后的灰度对应关系改变像素的灰度，即可生成新的图像。


表53直方图均衡化变换前后灰度级对应关系



变换前灰度级01/72/73/74/75/76/71
变换后灰度级1/73/75/76/76/7111



原始图像直方图和直方图均衡化处理后的图像直方图显示结果如图513所示。可以看出，图513(b)中对应的变换后的新直方图比图513(a)中的原图像的直方图要平坦多了。理想情况下，经过直方图均衡化处理后的图像直方图应是十分均匀平坦的，但实际情况并非如此，和理论分析有差异，这是由于图像在直方图均衡化处理过程中，灰度级作“近似简并”引起的结果。



图513直方图均衡化处理前后的直方图分布对比


【例5.5】编写程序，对图像进行直方图均衡化。

解： 程序如下。

Image=rgb2gray(imread('img5_3.bmp'));

histgram =imhist(Image);                         %统计图像直方图

［h w］=size(Image);       

NewImage=zeros(h,w);

s=zeros(256);       s(1)=histgram(1);

for t=2:256

s(t)=s(t-1)+histgram(t);                      %计算新的灰度值    

end

for x=1:w

for y=1:h

NewImage(y,x)=s(Image(y,x)+1)/(w*h);       %生成新图像

end

end

imshow(Image);title('原始灰度图像');

figure;imhist(Image);title('原始灰度图像的直方图'); 

axis tight;

figure;imshow(NewImage);title('直方图均衡化处理后图像');

figure;imhist(NewImage);title('直方图均衡化处理后图像的直方图');

axis tight;

也可以直接调用MATLAB提供的函数histeq,程序如下。 

NewImage=histeq(Image,256)

程序运行结果如图514所示。



图514直方图均衡化处理前后的图像及直方图


前面介绍的灰度级变换和直方图修正都是以灰度图像为处理目标的，假如是一幅彩色图像比较暗，能不能采用这两种方法进行增强呢？请扫二维码，查看讲解。



微课视频


5.3基于照度反射模型的图像增强

一般情况下，自然景物图像f(x,y)可以表示为光源照度场(照明函数)i(x,y)和场景中物体反射光的反射场(反射函数)r(x,y)的乘积，如式(513)所示。


f(x,y)=i(x,y)·r(x,y)(513)
其中，0<i(x,y)<∞，0<r(x,y)<1。一般把式(513)称为图像的照度反射模型。

近似认为，照明函数i(x,y)描述景物的照明，其性质取决于照射源，与景物无关。反射函数r(x,y)描述景物内容，其性质取决于成像物体的特性，而与照明无关。由于照明亮度一般是缓慢变化的，所以认为照明函数的频谱集中在低频段。由于反射函数随图像细节不同在空间快速变化，所以认为反射函数的频谱集中在高频段。这样，就可根据式(513)将图像理解为高频分量与低频分量的乘积的结果。

基于照度反射模型的处理算法，通常会借助于对数变换，将式(513)中两个相乘分量变成为两个相加分量。这样不仅能够简化计算，而且对数变换接近人眼亮度感知能力，能够增强图像的视觉效果。

下面主要介绍两种基于照度反射模型的处理算法： 同态滤波和Retinex增强。

5.3.1同态滤波

若物体受到照度明暗不匀，则图像上对应照度暗的部分的细节就较难辨别。同态滤波的主要目的就是消除不均匀照度的影响，增强图像细节。

同态滤波的基本原理是根据图像的照度反射模型，对原始图像f(x,y)中的反射分量r(x,y)进行扩展，对光照分量i(x,y)进行压缩，以获得所要求的增强图像。

同态滤波的具体算法步骤如下。

(1) 对图像函数f(x,y)取对数，即进行对数变换处理



z(x,y)=ln［f(x,y)］=ln［i(x,y)·r(x,y)］
=ln［i(x,y)］+ln［r(x,y)］(514)

(2) 进行傅里叶变换处理


Z(u,v)=DFT{z(x,y)}=DFT{ln［i(x,y)］}+DFT{ln［r(x,y)］}
=I(u,v)+R(u,v)(515)

(3) 进行同态滤波处理



S(u,v)=Homo(u,v)·Z(u,v)
=Homo(u,v)·I(u,v)+Homo(u,v)·R(u,v)(516)

如前所述，图像的对数傅里叶变换的高频分量主要对应反射分量，低频分量主要对应照射分量，因此，需要设计合适的同态滤波函数Homo(u,v)，对高频、低频成分产生不同的影响，目的是为了压低照明分量、扩大反射分量。

可以看出，同态滤波函数的功能类似于高通滤波函数，因此，可在传统高通滤波函数的基础上加以改动，逼近同态滤波函数。假设高通滤波转移函数为High(u,v)，则由High(u,v)到Homo(u,v)的映射关系如式(517)所示。

图515同态滤波函数特性曲线




Homo(u,v)=(γH－γL)·High(u,v)+γL
(517)

式中，γH和γL分别表示高频分量频率场和低频分量频率场滤波特性，且当取值γH>1，0<γL<1时，照射分量受到抑制，反射分量得到增强，从而突出图像的轮廓细节。同态滤波函数Homo(u,v)的特性曲线如图515所示。

(4) 求傅里叶反变换



s(x,y)=DFT－1{S(u,v)}
=DFT－1{Homo(u,v)·I(u,v)}+DFT－1{Homo(u,v)·R(u,v)}
=i′(x,y)+r′(x,y)(518)

(5) 求指数变换，得到经同态滤波处理的图像


g(x,y)=e{s(x,y)}=e{i′(x,y)+r′(x,y)}
=i0(x,y)·r0(x,y)(519)
式中，i0(x,y)是处理后的照射分量，r0(x,y)是处理后的反射分量。

需要指出的是，在傅里叶平面上用同态滤波器来增强高频分量以突出轮廓细节的同时，也平滑了低频分量，使得图像中灰度变化平缓区域出现模糊。因此，通常会增加一个后滤波处理补偿低频分量，以使得图像得到很大的改善。

【例5．6】编写程序，对图像进行同态滤波增强。

解： 程序如下。

Image=rgb2gray(imread('img5_6.jpg')); 

［M,N］=size(Image);  

Hmin=0.5;  Hmax=2.0;%可根据需要效果调整参数  

c=1.1;  d0=1800;  

Image1=log(double(Image)+1);%取对数  

FImage=fft2(Image1);%傅里叶变换  

n1=floor(M/2);    n2=floor(N/2);  

for i=1:M  

for j=1:N  

D(i,j)=((i-n1)^2+(j-n2)^2); %高斯同态滤波 

H(i,j)=(Hmax-Hmin)*(1-exp(c*(-D(i,j)./(d0^2))))+Hmin;

end  

end  

Image2=ifft2(H.*FImage);%傅里叶反变换  

result=real(exp(Image2));  

imshow(I),title('原图');

figure,imshow(result,［］);title('同态滤波增强后');

程序运行结果如图516所示。可以看出，原始的光照不均匀图像［见图516(a)］经过同态滤波增强处理后，质量有了很大改善。



图516同态滤波增强


5.3.2Retinex增强


Retinex源于Retina(视网膜)和Cortex(大脑皮层)合成词的缩写，故Retinex理论又被称为“视网膜大脑皮层理论”。Retinex理论的基本原理模型是以人类视觉系统为出发点发展而来的一种基于颜色恒常性的色彩理论，该理论认为： ①人眼对物体颜色的感知与物体表面的反射性质有着密切关系，即由物体对红、绿、蓝三色光线的反射能力来决定，反射率低的物体看上去较暗，反射率高的物体看上去较亮； ②人眼对物体色彩的感知具有一致性，不受光照变化的影响。

基于Retinex理论的增强方法的基本原理就是根据图像的照度反射模型，通过从原始图像中估计光照分量，然后设法去除(或降低)光照分量，获得物体的反射性质，从而获得物体的本来面貌。

根据采用的不同的估计光照分量的方法，产生了各种Retinex算法。这里主要介绍中心环绕Retinex方法。在中心环绕Retinex方法中，估计光照分量的计算如式(520)所示，其物理意义是通过计算被处理像素与其周围区域加权平均值的比值消除照度变化的影响。


i′C(x,y)=F(x,y)fC(x,y)(520)
其中，C∈{R,G,B}，fC(x,y)是图像f(x,y)的第C颜色通道的亮度分量。i′C(x,y)是第C颜色通道的光照分量估计值。F(x,y)是中心环绕函数，一般采用高斯函数形式的环绕函数,即


F(x,y)=K·e－x2+y2σ2(521)
其中，σ为标准差，表示高斯环绕函数的尺度常数，决定了卷积核的作用范围。K为归一化因子，使得



∫∫F(x,y)dxdy=1(522)

中心环绕Retinex方法主要分为单尺度Retinex(SingleScale Retinex, SSR)增强方法、多尺度Retinex(MultiScale Retinex，MSR)增强方法和带色彩恢复的多尺度Retinex（MultiScale Retinex with Color Restoration,MSRCR）增强方法。

1. 单尺度Retinex增强方法

单尺度Retinex增强方法的具体步骤如下： 

(1) 根据式(520)、式(521)计算第C颜色通道的光照分量估计值i′C(x,y)； 

(2) 对fC(x,y)取对数，即进行对数变换处理，有



ln［fC(x,y)］=ln［iC(x,y)·rC(x,y)］
=ln［iC(x,y)］+ln［rC(x,y)］(523)

(3) 在对数域中，用fC(x,y)减去光照分量估计i′C(x,y)，计算反射分量，即获得图像的高频分量，有



RC(x,y)=ln［r′C(x,y)］=ln［fC(x,y)］－ln［i′C(x,y)］
=ln［fC(x,y)］－ln［F(x,y)fC(x,y)］(524)
式中，RC(x,y)是第C颜色通道的单尺度Retinex增强输出图像。

2. 多尺度Retinex增强方法

根据式(521)，高斯环绕函数F(x,y)的尺度σ的大小会直接影响图像增强结果。尺度σ 越小，越能够较好地完成动态范围压缩，但全局照度损失，图像呈现“白化”； 尺度σ 越大，越能够较好地保证图像的色感一致性，但局部细节模糊，强边缘处有明显“光晕”。单尺度Retinex增强方法很难在动态范围压缩和色感一致性上寻找到平衡点。因此，对单尺度改进后，就产生了多尺度Retinex增强方法。

多尺度Retinex增强实质上是多个不同单尺度Retinex增强的加权平均。

多尺度Retinex增强方法的具体步骤如下。

(1) 设置不同尺度σn，n=1,2,…,N。其中，N为设置的不同尺度个数。

(2) 根据式(521)，计算不同尺度的中心环绕函数Fn(x, y)。

(3) 根据式(525)，求图像的不同尺度的Retinex增强输出。


RnC(x,y)={ln［fC(x,y)］－ln［Fn(x,y)*fC(x,y)］}(525)
式中，RnC(x,y)是第C颜色通道第n个尺度的Retinex增强输出。

(4) 对多个不同尺度的Retinex增强输出结果进行加权平均。

RMSRCx,y=∑Nn=1wnRnCx,y(526)

式中，RMSRCx,y是第C颜色通道的多尺度Retinex增强输出。wn是给不同尺度σn分配的权重因子，满足： ∑wn=1。

3．  带色彩恢复的多尺度Retinex增强方法

在上述算法的彩色图像增强过程中，由于是对RGB每一个通道直接进行N次不同尺度滤波，导致三通道增幅不一，恢复的RGB比值与原图比值不太一样，引起颜色失真的问题。因此，产生了带色彩恢复的多尺度Retinex增强方法。

带色彩恢复的多尺度Retinex增强方法的表达式如下： 

RMSRCRCx,y=RMSRCx,y·γCx,y(527)

式中，γCx,y是第C颜色通道的色彩恢复系数，其定义如式(528)所示。

γCx,y=η·lnβ·fCx,y∑C∈R,G,BfCx,y(528)

式中，η为增益常数，β为非线性强度的控制因子。

【例5．7】编写程序，对图像进行Retinex增强方法。

解： 程序如下。

Image=double(imread('img5_6.jpg'));  %打开图像并转换为double数据

［height,width,c］=size(Image);  %获得图像尺寸

RI=Image(:,:,1);   GI=Image(:,:,2);   BI=Image(:,:,3);

alpha=0.4;   beta=125;

CR=alpha*(log(beta*(RI+1))-log(RI+GI+BI+1));%R通道色彩恢复系数

CG=alpha*(log(beta*(GI+1))-log(RI+GI+BI+1));%G通道色彩恢复系数

CB=alpha*(log(beta*(BI+1))-log(RI+GI+BI+1));%B通道色彩恢复系数

sigma=［100 15 80 250］;  filtersize=［height,width］; %设置高斯滤波器参数

Rhigh1=zeros(height,width);Rhigh2=zeros(height,width);

Ghigh1=zeros(height,width);Ghigh2=zeros(height,width); %初始化

Bhigh1=zeros(height,width);Bhigh2=zeros(height,width);

for i=1:4

gaussfilter=fspecial('gaussian',filtersize,sigma(i));  %构造高斯低通滤波器

Rlow=imfilter(RI,gaussfilter,'replicate','conv'); %高斯卷积后，获得R通道低频分量

Glow=imfilter(GI,gaussfilter,'replicate','conv'); %获得G通道低频分量

Blow=imfilter(BI,gaussfilter,'replicate','conv'); %获得B通道低频分量

if i==1

Rhigh=log(RI./Rlow+1);%获得R通道高频分量

Ghigh=log(GI./Glow+1);%获得G通道高频分量

Bhigh=log(BI./Blow+1);%获得B通道高频分量

SSRI=cat(3,Rhigh,Ghigh,Bhigh);  %%SSR增强

else

Rhigh1=1/3*(log(RI./Rlow+1)+Rhigh1);%获得R通道高频分量

Ghigh1=1/3*(log(GI./Glow+1)+Ghigh1);%获得G通道高频分量

Bhigh1=1/3*(log(BI./Blow+1)+Bhigh1);%获得B通道高频分量

Rhigh2=1/3*(CR.*log(RI./Rlow+1)+Rhigh2);%获得R通道高频分量

Ghigh2=1/3*(CG.*log(GI./Glow+1)+Ghigh2);%获得G通道高频分量

Bhigh2=1/3*(CB.*log(BI./Blow+1)+Bhigh2);%获得B通道高频分量

end

end 

MSRI=cat(3,Rhigh1,Ghigh1,Bhigh1);  %MSR增强图像

MSRCRI=cat(3,Rhigh2,Ghigh2,Bhigh2);  %MSRCR增强图像

subplot(2,2,1);imshow(uint8(Image));title('原图');

subplot(2,2,2);imshow(SSRI);title('SSR增强后图像');

subplot(2,2,3);imshow(MSRI);title('MSR增强后图像');

subplot(2,2,4);imshow(MSRCRI);title('MSRCR增强后图像');

程序运行结果如图517所示。



图517Retinex增强效果


启发： 

Retinex理论诠释了同样物体在有差异的光源或光线底下颜色恒定的机理，即认为人眼能够以某种方式过滤掉光照等影响，而只保留反映物体本质特征的信息，如反射率，从而确定颜色。正如“透过现象看本质”，需要突破认知上的束缚，识别出信息背后的本质，做出更有效的决策。

5．4色彩变换

色彩是表征图像信息的重要组成部分。色彩信息使得图像能够提供更丰富的内容、更清晰的细节，从而具有更好的视觉表达力。

5．4．1色彩迁移

图像色彩迁移是一种改变图像颜色基调的方法，实现将一幅图像的色彩信息传递到另一幅图像中，使两幅图像具有相类似的色彩特征。

色彩空间的选取对图像的色彩迁移算法是否有效有很大的影响。Smith等于1975年基于人体视网膜锥状细胞对光的波长相当敏感这一现象提出了LMS颜色空间，其中L、M、S三个通道分别表示长、中、短激发光谱。由于 LMS 颜色空间三个通道间有较大的相关性，给图像处理过程带来一定的困难。针对这种情况，Ruderman等于1998年在LMS颜色空间基础上提出了一种不相关的、近似正交的均匀色彩空间——lαβ颜色空间，其中，l表示非彩色的亮度通道，α表示彩色的黄蓝相关通道，β表示彩色的红绿相关通道，0≤l≤100，-128≤α,β≤127。lαβ颜色空间不仅假设人类视觉系统理想地适应自然基色的处理，建立把光线波长转换为亮度和色相的一套描述色彩数据，而且最大限度地降低了各通道间的相关性，使得每个通道的处理和运算都可以独立进行，不会影响到其他通道的信息，避免了通道交叉的问题，它为图像色彩迁移的发展奠定了基础。

Reinhard等于2001年通过对lαβ颜色空间各分量进行统计分析，首次提出了彩色图像的色彩迁移算法，算法思想是在lαβ颜色空间，通过对目标图像和源图像的色彩特征信息分别进行统计分析，进而得到两者间的变换关系，可使得目标图像和源图像最终具有相近的均值和标准差，这样就将源图像的色彩特征传递给了目标图像，使之看上去更加符合人类的视觉审美要求。

1．  lαβ与RGB相互转换

RGB转换为lαβ的过程如下。

由RGB空间先转换到LMS空间，再转换到lαβ空间，有

LMS=0．38110．57830．04020．19670．72440．07820．02410．12880．8444RGB (529)
lαβ=13000160001211111－21－10log10Llog10Mlog10S(530)


lαβ转换为RGB的过程如下。

由lαβ空间先转换到LMS空间，再转换到RGB空间，有

LMS=11111－11－20130001600012lαβ(531)
RGB=4．4679－3．58730．1193－1．21862．3809－0．16240．0497－0．24391．204510L10M10S (532)

式（530）中选择10为底的对数运算，因此式（532）中进行10为底的指数运算，其他合理底数也可以。

2． Reinhard色彩迁移

Reinhard色彩迁移算法的具体步骤如下。

(1) 将源图像和目标图像由RGB变换到lαβ颜色空间，并分别计算l、α、β三分量的均值和方差。

(2) 利用下式匹配色彩统计量，使其从源图像传递到目标图像。

lres=lin－μlin×σlrefσlin+μlref
αres=αin－μαin×σαrefσαin+μαref
βres=βin－μβin×σβrefσβin+μβref(533)

其中，lres、αres、βres为色彩迁移结果图像的l、α、β三分量值； lin、αin、βin为目标图像的l、α、β三分量值； μlref、μαref、μβref、σlref、σαref、σβref分别为源图像的lαβ三分量的均值和标准差； μlin、μαin、μβin、σlin、σαin、σβin分别为目标图像的l、α、β三分量的均值和标准差。

(3) 将调整后的目标图像分量由lαβ变换到RGB颜色空间，完成目标图像的色彩处理。

【例5．8】编写程序，对彩色图像进行色彩迁移。

解： 程序如下。 

clear,clc,close all;

ImageS=imread('img5_8_s.jpg');                 

ImageT=imread('img5_8_t.jpg');

subplot(131),imshow(ImageS);title('源参考图像');

subplot(132),imshow(ImageT);title('要处理的目标图像');

［sh,sw,sc］=size(ImageS);  ［th,tw,tc］=size(ImageT);     

［sR,sG,sB］=imsplit(double(ImageS));

［tR,tG,tB］=imsplit(double(ImageT));

Trans1=［0.3811,0.5783,0.0402;0.1967,0.7244,0.0782;0.0241,0.1288,0.8444］;

Trans2=［1 1 1;1 1 -2;1 -1 0］;

Trans3=［1/sqrt(3) 0 0;0 1/sqrt(6) 0;0 0 1/sqrt(2)］;   %色彩空间变换矩阵

svR=sR(:)'; svG=sG(:)'; svB=sB(:)';

sLMS=Trans1*［svR;svG;svB］;                   

tvR=tR(:)'; tvG=tG(:)'; tvB=tB(:)';

tLMS=Trans1*［tvR;tvG;tvB］;                         %将RGB数据转换为LMS数据

slab=Trans3*Trans2*［log10(sLMS)］;

tlab=Trans3*Trans2*［log10(tLMS)］;                 %再转到lαβ色彩空间

s_mu=mean(slab,2);    s_sigma=std(slab,0,2);     %源图像lαβ各通道均值和标准差

t_mu=mean(tlab,2);    t_sigma=std(tlab,0,2);     %目标图像lαβ各通道均值和标准差

tlab=(tlab-t_mu).*s_sigma./t_sigma+s_mu;         %目标图像色彩迁移

tLMS=(Trans3*Trans2)＼tlab;                  

tRGB=Trans1＼(10.^tLMS);                             %从lαβ色彩空间变换回RGB空间

tRGB(tRGB>255)=255;    tRGB(tRGB<0)=0;            %限幅

result1=reshape(tRGB(1,:),th,tw);

result2=reshape(tRGB(2,:),th,tw);

result3=reshape(tRGB(3,:),th,tw);

result=uint8(cat(3,result1,result2,result3));

subplot(133),imshow(result);title('Reinhard色彩迁移图像');

程序运行结果如图518所示。



图518Reinhard色彩迁移效果


色彩迁移不仅可以在两幅彩色图像之间进行，也可以在彩色图像和灰度图像之间进行，灰度图像色彩迁移主要是利用源参考图像的彩色信息，对灰度图像进行自动彩色化。关于灰度图像色彩迁移算法的设计实现，请扫二维码查看讲解。



微课视频


5．4．2色彩平衡

当一幅彩色图像数字化后，在显示时，颜色有时会看起来不正常，这是因为颜色通道中不同的敏感度、增光因子、偏移量等，导致数字化中的三个图像分量出现不同的变换，使结果图像的三原色“不平衡”，从而使景物中所有物体的颜色都偏离了其原有的真实色彩。因此，为了获得正常颜色的图像，将有色偏的图像进行颜色校正的过程就是色彩平衡。

下面主要介绍两种色彩平衡处理算法： 白平衡和最大颜色值平衡算法。

1．  白平衡算法

白平衡算法基本原理是，对于原始场景中某些理应为白色的像素点(即R*k=G*k=B*k=255)，由于色偏现象的影响，这些点的RGB三分量的值不再保持相同，通过调整这些点的三分量值，使之达到平衡，由此获得一套色彩平衡映射关系式，将这套映射关系式应用于整幅彩色图像，即可达到图像色彩平衡的目的。

白平衡方法有多种，下列是一种基本算法的步骤。

(1) 对存在色偏的图像，按照下式计算亮度分量： 

Y=0．299R+0．587G+0．114B (534)


(2) 由于现实场景中的白色点，在图像中可能不是理想状态的白色，即Y≠255，但白色点亮度为图像中的最大亮度，则需要求出图像中的最大亮度Ymax和平均亮度Y。

(3) 考虑到环境光照的适应性，寻找图像中所有亮度≤0．95Ymax的像素点构成的运算点集。

(4) 计算白色点集中的所有像素的RGB三分量的均值、、。

(5) 按下式计算色彩平衡调整参数kR、kG、kB。

kR=Y,kG=Y,kB=Y (535)


(6) 按下式对整幅图像的R、G、B三个颜色分量进行色彩平衡调整： 

Rk=kR·R,Gk=kG·G,Bk=kB·B (536)


【例5．9】编写程序，对图像进行白平衡处理。

解： 程序如下。 

clear,clc,close all;

Image=im2double(imread('img5_9.jpg'));

［R,G,B］=imsplit(Image);

Y=0.299*R+0.587*G+0.114*B;

maxY=max(Y(:));    muY=mean(Y(:));

R_set=R(Y<=0.95*maxY);  G_set=G(Y<=0.95*maxY); B_set=B(Y<=0.95*maxY);

%运算点集对应的RGB值

kR=muY/mean(R_set(:)); kG=muY/mean(G_set(:)); kB=muY/mean(B_set(:));

%计算色彩平衡调整参数

result=cat(3,kR*R,kG*G,kB*B);

subplot(121),imshow(Image),title('原图');

subplot(122),imshow(result),title('白平衡后图像');

程序运行结果如图519所示。



图519白平衡效果


2．  最大颜色值平衡算法

白平衡算法对存在白色像素点的图像有较好的色彩平衡效果。但若一幅图像中白色点不存在，或者仅占很少比例，则白平衡方法的处理就不是很有效。最大颜色值平衡算法就是针对这种情况提出的色彩平衡算法。

最大颜色值平衡基本原理是，如果一幅图像存在色偏现象，则RGB三通道中存在某个颜色信息比较强的通道，通过对该颜色通道的抑制，或者对另外两个颜色信息较弱的颜色通道的增强，就可以达到色彩平衡的目的。

最大颜色值平衡算法的具体步骤如下。

(1) 对存在色偏的图像，计算其RGB三通道中每个颜色通道的最大强度值Rmax、Gmax、Bmax。并计算最小值Cmax=minRmax，Gmax，Bmax。

(2) 统计RGB三通道中每个颜色通道强度值≥Cmax的像素的个数NR、NG、NB。并计算最大值Nmax=maxNR，NG，NB，其所对应的为颜色信息最强通道。

(3) 将RGB三通道中每个颜色通道像素值从大到小排序，寻找每个颜色通道中第Nmax个最大的强度值Rth、Gth、Bth。

(4) 按下式计算色彩平衡的调整参数kR、kG、kB： 

kR=CmaxRth,kG=CmaxGth,kB=CmaxBth(537)


(5) 按下式对整幅图像的R、G、B三个颜色分量进行色彩平衡调整。

Rk=kR·R,Gk=kG·G,Bk=kB·B (538)


【例5．10】编写程序，对图像进行最大颜色值平衡处理。

解： 程序如下。 

Image=double(imread('3.jpg'));

IR=Image(:,:,1);IG=Image(:,:,2);IB=Image(:,:,3);   

Cmax=min(［max(IR(:)) max(IG(:)) max(IB(:))］);  %从RGB每个通道的最大强度值中取最小值

numR=length(find(IR>=Cmax));   numG=length(find(IG>=Cmax)); 

numB=length(find(IB>=Cmax));  %统计RGB每个通道中像素值>=Cmax的像素数目

Nmax=max(［numR numG numB］); 

Tr_temp=sort(IR(:),'descend');   Tg_temp=sort(IG(:),'descend');

Tb_temp=sort(IB(:),'descend');    %排序

Tr=Tr_temp(Nmax);  Tg=Tg_temp(Nmax); 

Tb=Tb_temp(Nmax); %寻找每个通道中第Nmax个最大的强度

kr=Cmax/Tr;  kg=Cmax/Tg;   kb=Cmax/Tb; %计算色彩平衡调整参数

res=cat(3,kr*IR,kg*IG,kb*IB);

imshow(uint8(Image));

figure,imshow(uint8(res));title('最大颜色值平衡后图像');

程序运行结果如图520所示。



图520最大颜色值平衡效果


5．5其他图像增强方法

除了前述的常见的图像增强方法外，模糊增强、基于对数图像处理模型的图像增强、图像去雾增强等方法在图像增强处理中也得到了成功应用。

5．5．1模糊增强

模糊性就是事物的性质或类属的一种不分明性。1965年，美国加利福尼亚大学控制论专家L．A．Zadeh教授提出了用模糊集来处理这种模糊性。从此，以模糊集合为基础的模糊数学诞生了。目前，模糊技术已被广泛地应用于自然科学与社会科学的许多领域。

令U为元素(对象)集，u表示U的一类元素，即U=u，则该集合称为论域U。论域U到0,1闭区间的任一映射μA为

μA： U→0,1，u→μAu(539)

都确定U的一个模糊集合A，μA称为模糊集合的隶属函数。μAu称为u对于A的隶属度，取值范围为0,1，其大小反映了u对于模糊集合A的从属程度。

因此，模糊集合A完全可由u值和相应的隶属度函数μAu来表示描述，即

A=u,μAu|u∈U(540)

相应地，模糊集合的运算可由其隶属函数的运算来定义。

1． 图像的模糊特征平面

依照模糊集的概念，一幅灰度级为L的M×N的二维图像X，可以看作一个模糊点阵，记为

X=μ11x11μ12x12…μ1Nx1Nμ21x21μ22x22…μ2Nx2N︙︙︙μM1xM1μM2xM2…μMNxMN

或

X=∪Mm=1∪Nn=1μmnxmn(541)

式中，xmn是图像中像素m,n的灰度，μmnxmn表示图像中像素m,n的灰度xmn相对于某些特定灰度级x的隶属度，且0≤μmn≤1。换句话说，一幅图像X的模糊集合是一个从X到0,1的映射μmn。对于任一点xmn∈X，称μmnxmn为xmn在μmn中的隶属度。

可以看出，隶属度函数μmnxmn可以用于表征图像的模糊特征，可将图像从空间灰度域变换到模糊域。

2．  图像的模糊增强

图像的模糊增强的具体算法步骤如下。

(1) 进行模糊特征提取，将图像从空间灰度域变换到模糊域，有

μmn=Txmn=1+xmax－xmnFd－Fe(542)

式中，xmn是图像中像素m,n的灰度，xmax是图像中最大灰度值，Fd和Fe分别为指数和分数模糊因子，这些模糊因子可以在模糊域内改变μmn的值。一般情况下，指数因子Fe取值为2，分数因子Fd取值为

Fd=xmax－xc21Fe－1(543)

式中，xc为渡越点，其取值需要满足： μc=Txc=0．5且xc∈X。

(2) 在模糊域，对模糊特征进行一定的增强变换处理，有

μ′mn=Irμmn=2μ2mn,0≤μmn<0．5
1－21－μmn2,0．5≤μmn<1
Irμmn=I1Ir－1μmn(544)

式中，μ′mn为增强后的模糊域像素灰度值； r为正整数，表示迭代次数。可以根据不同需求选择迭代次数。

（3) 反变换，得到新的模糊增强后的输出图像。反变换如下： 

zmn=I－1μ′mn=xmax－Fdμ′mn1Fe－1(545)


图像的模糊增强结果如图521所示。



图521模糊增强效果


5.5.2基于对数图像处理模型的图像增强

通常，在对图像进行线性运算处理过程中，会出现许多边际效应。首先，可能会使输出值超出允许的灰度取值范围，产生“超区间值”问题，这类输出值会被修剪掉，从而导致图像信息的严重丢失； 其次，线性运算结果与人类的视觉感官往往有所偏差。所以，一般的标量乘法“×”和加法“+”运算并不通用于图像的形成法则。对数图像处理(Logarithmic Image Processing，LIP)模型的提出，不仅解决了“数值溢出”的缺陷问题，而且该模型符合人类视觉系统的非线性(对数)感知处理过程，能更好地分析图像的细节，突出图像中人们感兴趣的区域，一直是人们研究的热点。后来，人们从经典的对数图像处理模型发展出许多改进的对数图像处理模型，例如GLIP、PLIP、HLIP、SLIP、PSLIP等。本节主要介绍经典的对数图像处理模型及其在对比度增强中的应用。

1. 对数图像处理模型

在对数图像处理模型中，一幅图像通常用灰度色调函数表示。设图像支撑集为DR2，灰度色调函数的定义如下所示。


f(x,y)=M－I(x,y)(546)
式中，I(x,y)为输入的原始图像，f(x,y)为定义在［0,M)区间内的灰度色调函数，M 是严格的正值。对于一般的8比特的图像，M一般取值为256。

灰度色调函数的物理意义是光线通过一个光强滤波器形成了进入人眼的透射光从而成像，本质上就是把原图像转换为对应的光强滤波函数表示，该光强滤波函数就被定义为灰度色调函数。

对数图像处理模型是一个完备的数学理论，它规定了一系列特殊的代数运算和函数使运算处理前后像素的灰度值都在允许的灰度范围［0,M)内。

1) 基础向量运算

对数图像处理模型定义了三种基本向量运算： 对数加法、对数乘法和对数减法。

假设f1(x,y)、f2(x,y)是定义在［0, M)区间内的灰度色调函数，则有如下3种运算。

(1) 对数加法运算： 
f1(x,y)和f2(x,y)的对数加的定义为


f1(x,y)f2(x,y)=f1(x,y)+f2(x,y)－f1(x,y)·f2(x,y)M(547)

(2) 对数乘法运算： 
设λ为一个正实标量，则f1(x,y)和λ的对数乘的定义为


λf1(x,y)=M－M1－f1(x,y)Mλ(548)

(3) 对数减法运算： 
f1(x,y)和f2(x,y)的对数减的定义为


f1(x,y)Θf2(x,y)=M·f1(x,y)－f2(x,y)M－f2(x,y)(549)
可以根据以上基本运算来构造图像处理的一些复杂运算来对图像进行相应的处理。

2) 基本同态函数

对数图像处理模型的基本同态函数的正变换定义为


φ(f)=－M·ln1－f(x,y)M(550)

其反变换定义为


φ－1(f)=M·1－e－f(x,y)M(551)

通过对数图像处理模型的基本同态函数的正变换，灰度色调函数f(x,y)的取值范围由区间［0,M)扩大到(－∞,M)，但小于零的部分没有实际意义，且容易扭曲图像的负值信息。由于对数图像处理模型的运算空间和输入图像的实数空间通过正变换、反变换公式同构，所以可以由此推出对数图像处理模型空间的运算公式在实数空间的表示式。

2. 基于对数图像处理模型的增强

假设f(x,y)是原始输入图像的灰度色调函数，f′(x,y)是增强处理后图像的灰度色调函数，A(x,y)是以像素点(x,y)为中心的n×n的邻域窗口的灰度色调函数的平均值，基于对数图像处理模型的增强算法表达式为


f′(x, y)=αA(x, y)β［f(x, y)ΘA(x, y)］(552)
其中，参数α用来调整图像对比度。当α>1时，图像变亮； 当α<1时，图像变暗。参数β用来调整图像的锐化效果。当β>1时，增强图像边缘； 当β<1时则反之，图像边缘模糊。

若将f(x,y)进行基本正态函数的正变换处理，且令f-(x,y)=1-f(x, y)/M，则可进一步得到化简后的基于对数图像处理模型的增强算法表达式为

lnf-′（x,y）=α·ln （x,y）+β·［lnf-（x,y）-ln （x,y）］
(553)
其中


ln(x, y)=1n×n∑x+n2k=x－n2∑y+n2l=y－n2lnf-(x, y)(554)

简化后的增强算法的计算复杂度大幅度降低，同时增大了动态调节范围及边缘信息，实现了快速增强图像的目的。

基于对数图像处理模型的图像增强的具体算法步骤如下。

(1) 对输入图像I(x,y)进行变换，求得灰度色调函数f(x,y)。

(2) 对f(x,y)进行基本正态函数的正变换φ(f)，且有φ(f)=lnf-。

(3) 如式(554)所示，求以像素点(x,y)为中心的n×n的邻域窗口的平均值ln。

(4) 如式(553)所示，进行基于对数图像处理模型的增强处理。

(5) 进行基本正态函数的反变换φ-1(f)。

基于对数图像处理模型的图像增强结果如图522所示。图522(a)为原始图像； 图522(b)为基于对数图像处理模型增强后的图像，可以看出图像的对比度明显增强。



图522基于对数图像处理模型的图像增强效果


5.5.3图像去雾增强

摄像机所拍摄的图像，由于受到雾的影响，往往出现模糊不清、对比度下降、色彩失真等质量退化的现象，严重影响到后续的处理工作。图像去雾就是采用一定的处理方法或手段，去除图像中雾的影响，提高图像的对比度，增强图像的细节信息，获得视觉效果较好的图像。

目前，图像去雾方法主要分为两类。一类方法是基于非物理模型的去雾方法，是不考虑雾导致图像退化的成因，只通过实现对比度增强方法来达到去雾的目的。该类方法虽然能一定程度上达到去雾效果，但是可能会造成部分信息的损失。例如，子块部分重叠的局部直方图均衡化方法、对比度受限自适应直方图均衡化方法、Retinex增强方法等。另一类方法是基于物理模型的去雾方法，考虑雾导致图像退化的成因，进行数学建模，补偿退化过程造成的失真，恢复无雾图像。该类方法获得的无雾图像，视觉效果自然，一般无信息损失。

1. 图像去雾物理模型

在恶劣的大气环境(霾、雾、灰尘等)中进行拍摄时，光线从物体表面到达观测点(如摄像头)之前，被大气中的悬浮粒子散射，使得拍摄到的图像模糊不清，图像质量下降，雾天图像退化模型如图523所示。这种复杂的散射作用分为两类： 光线通过大气介质时引起的入射光衰减效应和大气介质中的粒子散射的环境光所引起的大气光散射效应。



图523雾天图像退化模型



那么，入射光线到达摄像头的总光强度，为入射光线衰减后到达摄像头的光强和周围环境中的各种散射光线进入摄像头后的附加光强之和，即


I(x)=Ed(d,λ)+Es(d,λ)(555)
其中，


Ed(d,λ)=J(x)·e－β(λ)·d(x)(556)
Es(d,λ)=A(1－e－β(λ)·d(x))(557)

式中，x代表二维空间位置。d(x)为场景深度(图像上的目标到观察者的距离)。Ed(d, λ)为直接衰减，描述了场景目标反射光在介质中衰减的结果。Es(d,λ)为空气光，反映了全局大气光的散射导致杂散光成像的情况。I(x)为观察图像，即输入图像。J(x)为场景目标直接反射光成像的亮度。A为大气光，是沿着观测视线场景物体在无穷远处的光强，与x无关。β为空气的散射系数，由于大雾天气对可见光的衰减几乎不受波长影响，因此可以认为β与波长λ无关，将β假设为一个不变的常量值，得


I(x)=J(x)·e－β·d(x)+A(1－e－β·d(x))(558)
这就是图像去雾的物理模型。

2. 基于暗原色先验的去雾方法

He Kaiming等通过对户外大量清晰无雾的自然图像观察统计得到： 在绝大多数非天空的局部区域内，某一些像素在RGB三色通道中至少有一个通道的像素颜色值比较低。换言之，该区域光强度的最小值是个很小的数。

若给出一幅清晰无雾图像J，暗通道Jdark用数学公式定义为


Jdark(i,j)=minC∈{r,g,b}miny∈Ω(i,j)(JC(y))(559)

式中，JC表示图像J的RGB三色通道中的一个颜色通道，Ω(i,j)表示以点(i,j)为中心的一个局部区域块。对于非天空清晰无雾图像J，Jdark值非常低，趋于0，那么Jdark称为J的暗颜色，并且把以上观察得出的经验性规律称为暗原色先验。

基于图像去雾的物理模型，定义t(x)为介质传播函数，表示透射率，0≤t(x)≤1，且为


t(x)=e－β·d(x)(560)
则将式(558)变换为


I(x)=J(x)·t(x)+A(1－t(x))(561)
暗原色先验去雾的目的就是希望能够通过雾天图像I(x)恢复无雾图像J(x)。

附加的大气光导致图像被雾干扰之后往往要比其本身亮度更大，透射率一般较小，所以被浓雾覆盖的图像的暗原色具有较高的强度值。视觉上看来，暗原色强度值是雾浓度的粗略近似。因此，可利用这一性质来估算大气光A和透射率t(x)。

1) 估计大气光

对于大气光A的估计，由于其亮度值通常较高，所以主要通过对图像最不透明区域像素的亮度值估算而得。通常通过计算图像的暗原色最亮区域来估算大气光A。

估计大气光A的具体算法步骤如下。

(1) 划分雾天图像I的局部块Ω，求取其暗原色图Idark_channel； 

(2) 从暗原色图Idark_channel中按照亮度值从大到小取最亮的前0.1%像素，该部分像素的位置通常对应图像中最不透明的区域； 

(3) 在图像I中定位步骤(2)求取的像素位置区域，求该区域像素的亮度平均值Imean_dc； 

(4) 设置参数Amax=240，计算min(Amax,Imean_dc)作为大气光A的估计值。这里，设置参数Amax表示所允许的大气光A的估算最大值，以避免所计算的Imean_dc值接近于255时可能会造成的处理后图像的色调偏差。

大气光A的估计如图524所示。其中，图524(a)为原始图像，所标记的蓝框区域，代表原始图像中亮度最大区域，值为254。图524(b)为暗原色图，所标记的红框区域，代表暗原色图中的最亮点，值为208，作为估计的大气光A的值。可以看出，估计的大气光A并不一定为图像中亮度的最大值。



图524大气光估计


2) 估计透射率

基于暗原色先验理论和雾天图像成像模型，可以直接估计出成像时刻的雾浓度，从而粗略估算透射率t(x)。首先已知大气光A，假定局部区域Ω内透射率保持一致，对式(562)在RGB三个颜色通道进行最小运算，可求得暗原色图。



minC∈{r,g,b}minx∈Ω(i,j)IC(x)A=t(x)·minC∈{r,g,b}minx∈Ω(i,j)JC(x)A+(1－t(x))(562)

这里，C表示图像RGB三通道中的一个颜色通道。又根据暗原色先验规律，无雾自然图像的暗原色值接近于0，且A>0，则得


minC∈{r,g,b}minx∈Ω(i,j)JC(x)A→0(563)
将式(563)代入式(562)进行变换，可估算透射率为



t(x)=1－minC∈{r,g,b}minx∈Ω(i,j)IC(x)A(564)

又由于考虑空间透视现象、人眼的视觉特性及人观看景物时对于图像真实感和深度感的需求，针对式(564)引入一个调节参数ω∈(0,1)，ω值的选取通过经验获得，则得



t(x)=1－ω·minC∈{r,g,b}minx∈Ω(i,j)IC(x)A(565)

3) 获得无雾图像

将式(565)变形，可得


JC(x)=A+IC(x)－At(x)(566)
将有雾图像I(x)、估计大气光A和透射率t(x)，代入式(566)，计算获得无雾图像J(x)。

基于暗原色先验去雾的图像增强结果如图525所示。




图525暗原色先验去雾增强效果


3. 基于暗原色先验的低照度图像增强

低照度图像通常是指在光照较暗或者夜间所拍摄得到的图像。如果将低照度图像反转后，其图像表征及直方图表征与雾天图像具有很高的相似性，因此，可以采用图像去雾技术来实现低照度图像的增强。该算法的基本思想是： 对输入的低照度图像进行反转，对反转图像进行去雾处理，然后对去雾结果反转得到低照度图像增强效果。

基于暗原色先验的低照度图像增强的具体算法步骤如下： 

(1) 对输入低照度图像IC(x)进行取反操作，求得反转图像RC(x)： 


RC(x)=255－IC(x)(567)
式中，C表示图像的RGB三个颜色通道，x表示图像的坐标点。

(2) 对反转图像RC(x)进行暗原色先验去雾处理，获得无“雾”图像JC(x)。

(3) 将无“雾”图像JC(x)反转，获得低照度增强图像PC(x)： 


PC(x)=255－JC(x)(568)

基于暗原色先验的低照度图像增强结果如图526所示。



图526基于暗原色先验的低照度图像增强结果


习题

5.1试给出把灰度范围［0，10］伸长为［0，15］，把范围［10，20］移到［15，25］并把范围［20，30］压缩为［25，30］的变换方程。

5.2在对图像进行直方图均衡化时，为什么会产生简并现象？

5.3一幅图像的直方图如图527所示，可以对其进行什么处理？ 



图527习题5.3图


5．4一幅大小为64×64的图像，8个灰度级对应像素数及概率prr如表54所示，试对其进行直方图均衡化。


表54习题5.4中图像的灰度级分布


灰度级rk
0
1/7
2/7
3/7
4/7
5/7
6/7
1
像素数nk
560
920
1046
705
356
267
170
72
概率prr
0．14
0．22
0．26
0．17
0．09
0．06
0．04
0．02


5．5已知一幅图像的灰度级为8，图像的左边一半为深灰色，其灰度值为1/7，右边一半为黑色，其灰度值为0。试对此图像进行直方图均衡化，并描述处理后的图像视觉效果。

5．6同态滤波的特点是什么？适用于什么情况？

5．7Retinex增强方法的特点是什么？单尺度和多尺度Retinex增强方法的区别是什么？

5．8Reinhard色彩迁移方法的设计思想是什么？

5．9白平衡的基本原理是什么？色彩补偿的基本原理是什么？

5．10编写程序，实现习题5．1中的灰度级变换。

5．11编写程序，实现习题5．3中所设计的处理方法。

5．12编写程序，实现真彩色图像的同态滤波。

5．13编写程序，实现真彩色图像的基于对数图像处理模型的增强。

5．14编写程序，实现有雾图像的基于暗原色先验的去雾处理，并且与其他基于非物理模型的去雾方法进行结果对比。