数字图像处理第五章----图像复原与重建

1 图像退化

图像在传输过程中，由于成像系统、传输介质和设备的不完善，使图像的质量变坏，我们称这种现象为图像退化。图像发生退化后，一般来说我们就无法准确的获得图像的信息，为了得到与原图像尽可能一致的图像，我们就要用到图像复原的技术。
图像发生退化的原因主要有两种，一是图像经过滤波器的等电子装置，二是传输过程中噪声的干扰。
图像退化的一般模型为：g(x,y)=H[f(x,y)]+n(x,y)

1.1 噪声模型

图像发生退化现象的一个重要原因是噪声的干扰，一般来说，图像在传输过程中混入的主要是加性噪声，但也有乘性噪声，这种时候，我们一般对信号取log就可以将乘性噪声转为加性噪声。
1、高斯噪声
所谓高斯噪声是指它的概率密度函数服从高斯分布（即正态分布）的一类噪声。如果一个噪声，它的幅度分布服从高斯分布，而它的功率谱密度又是均匀分布的，则称它为高斯白噪声。由于高斯噪声在空间和频域中数学上的易处理性，所以经常被用于实践中。
它的函数为$p(z)=\frac{1}{\sqrt{2\pi }\sigma }{e}^{-\frac{(z-\mu {)}^2}{2{\sigma }^2}}$
其波形如下图所示:

2、瑞利噪声
瑞利分布是最常见的用于描述平坦衰落信号接收包络或独立多径分量接受包络统计时变特性的一种分布类型。通常，瑞利噪声的函数为：$p(z)=\frac{z}{{\sigma }^2}{e}^{-\frac{z^2}{2{\sigma }^2}}$
其分布函数为：

3、伽马噪声
其分布函数为：

噪声波形为：

4、指数噪声
其分布函数为：

波形为：

5、均匀噪声 $p(z)=\frac{1}{b-a}$
其波形为：

6、椒盐噪声

它的曲线为：

从椒盐噪声的概率密度函数可以看出，其函数是两个独立的脉冲，因此也被成为脉冲噪声。
噪声脉冲可以是正值也可以是负值，负脉冲在图像中相当于一个黑点，正脉冲在图像中相当于一个白点。

1.2 噪声直方图

加入噪声后，图像的直方图与噪声的函数分布是很相似的，有些时候，我们观察图像的直方图后就可以判断图像中混入的是什么噪声。

clear
I=imread('coins.png');
J1=imnoise(I,'salt & pepper',0.3);
J2=imnoise(I,'gaussian',0.2);
[m,n]=size(I); 
a=50;
b=180;
K=double(I);
J3=K+a+(b-a)*rand(m,n);%均匀分布噪声
a=0.04;
k=-1/a;
J4=k*log(1-rand(m, n));%指数分布噪声
figure;
subplot(2,3,1),imhist(I);
title('原图直方图');
subplot(2,3,2),imhist(J1);
title('加椒盐噪声');
subplot(2,3,3),imhist(J2);
title('加高斯噪声');
subplot(2,3,4),imhist(uint8(J3));
title('加均匀噪声');
subplot(2,3,5),imhist(uint8(J4));
title('加指数噪声');

上面是给图像加噪后图像的直方图，通过它们的直方图，我们可以很容易的判断出噪声的类型。
当图像混入的噪声较多，通过图像直方图无法判断噪声类型时，我们可以缩小范围，只对图像的一小块进行观察，一般来说，通过这种方法我们就能够判断出噪声的类型。

2 图像复原

2.1 空间滤波

一般来说，对带有噪声的图像，我们对其进行空间域滤波后就可以达到图像复原的目的。
方法：在已知噪声模型的基础上，根据噪声的类型对图像进行空间域滤波，如：均值滤波、中值滤波等。详见

数字图像处理第三章-----灰度变换与空间滤波

2.2 动态模糊

由于人眼存在视觉暂留效应，所以在看一个运动的物体时，看到的画面相当于当前画面不断发生轻微移动后叠加而成。其退化函数为:
$H(u,v)=\frac{Tsin[\pi (ua+vb)]}{\pi (ua+vb)}{e}^{-j\pi (ua+vb)}$
其中,T表示曝光时间，a和b分别表示水平和垂直方向上的移动量.
下面将一幅正常图像变为运动模糊图像：

I=imread('coins.png'); 
figure;
subplot(1,2,1),imshow(I);
LEN=25; 
THETA=11;
PSF=fspecial('motion',LEN,THETA);
Blurred=imfilter(I,PSF,'circular','conv');
subplot(1,2,2),imshow(Blurred);

2.3 逆滤波

假设图像传输过程中没有噪声的干扰，那么根据公式$F(u,v)=\frac{G(u,v)}{H(u,v)}$我们就可以完美的得到原图像的傅里叶变换，然后对其进行傅里叶反变换即可得到f(x,y)。

clear;
x=imread('coins.png'); 
x1=x(:,:,1);
x1=double(x1); 
y1=fftshift(fft2(x1));
[M,N]=size(y1); 
figure;
subplot(1,3,1);
imshow(x1,[]);
title('原始的图像') ;
m=1:M; 
n=1:N; 
[m,n]=meshgrid(m,n);
a=double(16/100);
b=double(16/100);
t=double(88/100);
f=ones(M,N);   
g=(m-M/2-1).*a+(n-N/2-1).*b+eps;   
f=t.*sin(pi.*g).*exp(-1j.*pi.*g)./(pi.*g);
h=f'.*y1; 
tu=ifft2(ifftshift(h));  
subplot(1,3,2);
imshow(tu,[]);
title('退化图像')
y1=h./f'; 
subplot(1,3,3);
imshow(abs(ifft2(y1)),[]);
title('逆滤波的结果');

可以看出，对不含噪声的图像，我们逆滤波的结果是非常完美的，逆滤波后的图像与原图像几乎一模一样。
但如果H(u,v)有零点的话，在该点进行复原时，结果就会变得非常大，得到的复原效果也会很差。而且我们现实生活中很难存在不含噪声的图像，因此逆滤波的应用范围非常小。

2.4 维纳滤波

采用拉格朗日乘数法，在有噪声的条件下，从退化图像g（x,y）复原出f(x,y)的估计值。
$F(f)=\frac{1}{H(f)}[\frac{|H(f){|}^2}{|H(f){|}^2+\frac{S_i(f)}{S_f(f)}}]$
其中，F代表重构信号的傅里叶变换，$S_i(f)$ 为噪声的能量谱 ，$S_f(f)$是信号的能量谱。
采用维纳滤波要求：未退化图像和噪声的功率必须是已知的。

clear;
x=imread('coins.png'); 
x1=x(:,:,1);
x1=double(x1); 
y1=fftshift(fft2(x1));
[M,N]=size(y1); 
figure;
subplot(1,3,1);
imshow(x1,[]);
title('原始的图像') ;
m=1:M; 
n=1:N; 
[m,n]=meshgrid(m,n);
noise=20.*imnoise(zeros(M,N),'gaussian',0,0.008);  
a=double(21/100);%x方向的最大移动量为ra的0.21倍,可调 
b=double(21/100);%y方向的最大移动量为ca的0.21倍,可调 
t=double(88/100);%移动到最大所需的时间默认为0.88
f=ones(M,N);   
g=(m-M/2-1).*a+(n-N/2-1).*b+eps;   
f=t.*sin(pi.*g).*exp(-j.*pi.*g)./(pi.*g);
h=f'.*y1; 
tu=ifft2(h);  
tu=abs(tu)+noise; 
subplot(1,3,2);
imshow(tu,[]);
title('退化的图像')%原图傅立叶变换估计值 
y1=h./f'; 
h=fftshift(fft2(tu));
x=fftshift(fft2(noise));
K=x.*conj(x)./(y1.*conj(y1));%计算K值  
w=(f.*conj(f))'.*h./(f.*(f.*conj(f)+K'))'; 
weina=abs(ifft2(w)); 
subplot(1,3,3);
imshow(weina,[]);
title('维纳滤波的结果');

可以看出，经过维纳滤波复原的图像效果已经是很不错了，基本没有丢失图像的信息。

3 总结

本章主要学习了图像退化的原因及图像复原的方法，图像退化的原因包括噪声干扰、滤波装置等，介绍了一些噪声的类型，并了解了加噪后图像直方图的特点，并了解了一个典型的图像退化模型—动态模糊。然后学习了一些图像复原方式，如空间滤波、逆滤波、维纳滤波等。图像复原是图像处理中非常重要的内容，如果我们能够把它学的很好的话，毫无疑问它可以帮助我们处理很多问题。