matlab-空域和频域进行滤波处理

1、图像平滑处理的空域算法和频域分析1 技术要求对已知图像添加高斯白噪声,并分别用低通滤波器(频域法)和邻域平均法(空域法)对图像进行平滑处理(去噪处理),并分析比较两种方法处理的效果。2 基本原理2.1 图像噪声噪声在理论上可以定义为“不可预测，只能用概率统计方法来认识的随机误差”。实际获得的图像一般都因受到某种干扰而含有噪声。引起噪声的原因有敏感元器件的内部噪声、相片底片上感光材料的颗粒、传输通道的干扰及量化噪声等。噪声产生的原因决定了噪声的分布特性及它和图像信号的关系。根据噪声和信号的关系可以将其分为两种形式：（1）加性噪声。有的噪声与图像信号g(x,y)无关，在这种情况下，含噪图像f(x,

2、y)可表示为f(x,y)=g(x,y)+n(x,y)（2）乘性噪声。有的噪声与图像信号有关。这又可以分为两种情况：一种是某像素处的噪声只与该像素的图像信号有关，另一种是某像点处的噪声与该像点及其邻域的图像信号有关，如果噪声与信号成正比，则含噪图像f(x,y)可表示为f(x,y)=g(x,y)+n(x,y)g(x,y)另外，还可以根据噪声服从的分布对其进行分类，这时可以分为高斯噪声、泊松噪声和颗粒噪声等。如果一个噪声，它的幅度分布服从高斯分布，而它的功率谱密度又是均匀分布的，则称它为高斯白噪声，一般为加性噪声。2.2 图像平滑处理技术平滑技术主要用于平滑图像中的噪声。平滑噪声在空间域中进行，其基

3、本方法是求像素灰度的平均值或中值。为了既平滑噪声又保护图像信号，也有一些改进的技术，比如在频域中运用低通滤波技术。（1）空域法在空域中对图像进行平滑处理主要是邻域平均法。这种方法的基本思想是用几个像素灰度的平均值来代替每个像素的灰度。假定有一幅N*N个像素的图像f(x,y)，平滑处理后得到一幅图像g(x,y)。g(x,y)由下式决定 式中，x,y=0,1,2,，N-1；S是(x,y)点邻域中点的坐标的集合，但其中不包括(x,y)点；M是集合内坐标点的总数。上式说明，平滑化的图像g(x,y)中每个像素的灰度值均由包含在(x,y)的预定邻域中的f(x,y)的几个像素的灰度值的平均值来决定。（2）频

4、域法低通滤波法是一种频域处理方法。在分析图像信号的频率特性时，一幅图像的边缘、跳跃部分以及颗粒噪声代表图像信号的高频分量，而大面积的背景区则代表图像信号的低频分量。用滤波的方法滤除其高频部分就能去掉噪声，使图像得到平滑。由卷积定理可知 其中F(u,v)是含有噪声的图像的傅立叶变换，G(u,v)是平滑处理后的图像的傅立叶变换，H(u,v)是传递函数。选择传递函数H(u,v)，利用H(u,v)使F(u,v)的高频分量得到衰减，得到G(u,v)后再经傅立叶反变换后就可以得到所希望的平滑图像g(x,y)了。根据前面的分析，显然H(u,v)应该具有低通滤波特性，所以这种方法叫低通滤波法平滑化处理。常用的

5、低通滤波器有如下几种：a.理想低通滤波器一个理想的二维低通滤波器有一个参数 。它是一个规定的非负的量，叫做理想低通滤波器的截止频率。所谓理想低通滤波器是指以截频 为半径的圆内的所有频率都能无损地通过，而在截频之外的频率分量完全被衰减。理想低通滤波器可以用计算机模拟实现，但是却不能用电子元器件实现。b.布特沃斯(Butterworth)低通滤波器一个n阶布特沃斯低通滤波器的传递函数由下式表示式中， 为截止频率。布特沃斯低通滤波器又称最大平坦低通滤波器。它与理想低通滤波器不同，它的通带与阻带之间没有明显的不连续性。也就是说，在通带和阻带之间有一个平滑的过度带，通常把H(u,v)下降到某一值的那一点

6、定为截止频率 。一般情况下常常采用下降到H(u,v)最大值的 那一点为截止频率点，该点也常称为半功率点，这样上式可修改成 c.指数低通滤波器 在图像处理中常用的另一种平滑滤波器是指数低通滤波器，由于指数低通滤波器有更快的衰减率，所以，经指数低通滤波器处理的图像比布特沃斯低通滤波器处理的图像稍模糊一些。3 建立模型描述图像平滑处理用到两种方法，其中空域法使用邻域平均法处理，其流程图如图1所示Imread()函数读取图像rgb2gray()函数将原图转换为灰度图imnoise()函数为原图加高斯白噪声double()将加了噪声的图转换为双精度型创建3*3模板创建5*5模板创建7*7模板用模板对噪声

7、图邻域平均用模板对噪声图邻域平均用模板对噪声图邻域平均Image()显示处理后的图像Image()显示处理后的图像Image()显示处理后的图像图1 邻域平均法模型图频域法用二阶布特沃斯低通滤波器对图像滤除噪声，其程序流程图如图2所示：imread()函数读取图像 Size（）函数获取图像矩阵大小 Image()显示处理后的图像Image()显示处理后的图像与噪声图的傅立叶变换矩阵相乘后用ifft2()求逆傅立叶变换与噪声图的傅立叶变换矩阵相乘后用ifft2()求逆傅立叶变换设定巴特我司低通滤波器的参数，并求出传递函数，截频20设定巴特我司低通滤波器的参数，并求出传递函数，截频100fft2(

8、)函数对噪声图进行傅立叶变换 将原图转换为灰度图wgn()函数产生高斯白噪声并加入原来图像中fft2()对图像进行傅立叶变换巴特沃斯低通滤波器传递函数与图像傅立叶变换相乘ifft2()函数进行逆傅立叶变换image()函数显示图像图2 低通滤波法模型图4 源程序代码（1）邻域平均法源程序clear;clc;close all;I=imread('tu.jpg'); %读取图像f=rgb2gray(I);%转化成灰度图subplot(2,3,1);% 分割2*3个窗口。取第一个窗口，下面在第一个窗口处显示图像imshow (f);%显示灰度图colormap(gray);titl

9、e('原始图')%给显示的图像命名为“原始图”J=imnoise(f,'gaussian',0,0.09);%给原始图加入参数为0.09的高斯白噪声subplot(2,3,2);% 分割2*3个窗口。取第二个窗口，如果下面有绘图语句，就表示要 %在第二个窗口中绘图imshow(J);%显示加了高斯白噪声的图像Jtitle('噪声图')%命名为“噪声图”m n=size(f); %获取灰度图的大小 f=double(f);%转换f为双精度型c=1/9*1 1 1;1 1 1;1 1 1; %3*3模板for i=1:m for j=1:n L=f(

10、i:i,j:j).*c; %求点积 G(i,j)=sum(sum(L);%求和 endendsubplot(2,3,3);image(G);%取第三个窗口title('3*3模板')%命名为“3*3模板”c=1/25*1 1 1 1 1;1 1 1 1 1;1 1 1 1 1; %5*5模板for i=1:m for j=1:n L=f(i:i,j:j).*c; %求点积 G(i,j)=sum(sum(L);%求和 endendsubplot(2,3,4);image(G);%在第四个窗口中显示图像title('5*5模板')% 命名为“5*5模板”c=1/49

11、*1 1 1 1 1 1 1;1 1 1 1 1 1 1;1 1 1 1 1 1 1; %7*7模板for i=1:m for j=1:n L=f(i:i,j:j).*c; %求点积 G(i,j)=sum(sum(L);%求和 endendsubplot(2,3,5);image(G);%在第五个窗口中显示图像title('7*7模板')% 命名为“7*7模板”（2）低通滤波法源程序clear;clc;close all;I=imread('tu.jpg'); %读取图像m n p=size(I);%获取图像矩阵大小I=double(I);将原来的图像矩阵转换为

12、双精度型的I=I(1:m,1:n,1).*0.3+I(1:m,1:n,2).*0.51+I(1:m,1:n,1).*0.11;%将原图转换为灰度图subplot(2,3,1);image(I);colormap(gray); %分割2*3个窗口。取第一个窗口，显示灰度图 title('原始图')%命名为“原始图”Noise=wgn(m,n,25);%产生25dBm的高斯白噪声New = Noise + I;%将高斯白噪声与原始信号叠加subplot(2,3,2);image(New);%在第二个窗口中显示加了噪声的图像title('噪声图')%命名为“噪声图”g

13、=fft2(New); % 傅立叶变换g=fftshift(g); % 转换数据矩阵M,N=size(g); %获取矩阵大小nn=2; % 定义二阶巴特沃斯(Butterworth)低通滤波器d0=20; %截止频率为20m=fix(M/2); n=fix(N/2);圆整函数for i=1:M for j=1:N d=sqrt(i-m)2+(j-n)2); h=1/(1+0.414*(d/d0)(2*nn); % 计算低通滤波器传递函数 result(i,j)=h*g(i,j);%传递函数与信号傅立叶变换函数相乘 endendresult=ifftshift(result);J2=ifft2(

14、result); %逆傅立叶变换J3=uint8(real(J2);%将输出图像转换为8位数据存储subplot(2,3,3);image(J3); %在第三个窗口中显示滤波后的图像title('截止频率20 ') %命名为“截止频率20”nn=2; % 定义二阶巴特沃斯(Butterworth)低通滤波器d0=100; %截止频率为100m=fix(M/2); n=fix(N/2);圆整函数for i=1:M for j=1:N d=sqrt(i-m)2+(j-n)2); h=1/(1+0.414*(d/d0)(2*nn); % 计算低通滤波器传递函数 result(i,j)

15、=h*g(i,j); %传递函数与信号傅立叶变换函数相乘 endendresult=ifftshift(result);J2=ifft2(result);%二维逆傅立叶变换J3=uint8(real(J2); %转换为8位存储的图像subplot(2,3,4);image(J3); %在第四个窗口中显示绿波后的图像title('截止频率100 ')%命名为“截止频率100”5 调试过程及结论在MatLab中输入邻域平均法程序代码后运行程序，得到输出结果如图3所示： 图3 邻域平均法输出图像由噪声图与滤波后的图像对比可看出，邻域平均法对抑制噪声有明显的效果，但随着邻域的加大，就是

16、随着模板的加大，图像的模糊程度也愈加严重。在MatLab中输入低通滤波法程序代码后运行程序，得到输出结果如图4所示：图4 低通滤波法输出图像巴特沃斯低通滤波器的去噪效果与所选的截止频率有关，由于一幅图像的边缘、跳跃部分以及颗粒噪声代表图像信号的高频分量，而大面积的背景区则代表图像信号的低频分量所以，截频设的太低会使图像变得越模糊，因为图像的许多细节信息也被滤掉了。6 心得体会这次能力强化训练是图像平滑处理的空预算法和频域分析，需要有一定的对数字图像处理方面的知识，由于没有学习过数字图像处理这门课，已开始题目布置下来时感到无从下手,查阅了大量资料，对图像平滑处理有了一定的了解，对于MatLab刚

17、开始接触，许多语法、函数都不会使用，于是和同学互相讨论，通过参考大量资源，站在巨人的肩膀上稍加整合创新，最终得出了自己的程序，程序调试运行后，发现最开始根本无法显示图像，检查原因，原来是对MatLab的符号规则不清楚，MatLab中的标点符号必须是在英文状态下输入，这在编写程序时必须小心！修改后得到了正确的显示，但是图像效果对比不明显，通过不断摸索试探，使用不同的图像显示函数，得到了对比鲜明的图像，对分析不同方法的优缺点提供了很大方便。 这次数字图像处理课程设计虽然时间短暂，但我在短时间内学到了许多以前没有接触过的数字图像处理方面的知识，并最终完成了程序设计，得到正确的结果，学到了许多MatLab的知识。为了完成这次课程设计，我花费了很多精力学习MatLab，也查阅了很多关于MatLab