专业综合课程设计说明书图像频率域处理程序设计.doc

目录摘要21 MATLAB软件功能简介32 数字图像处理简介42.1 数字图像处理的目的42.2 数字图像处理的方法52.3 图像频域处理的概述53 二维傅里叶变换73.1 二维连续傅里叶变换73.2 二维离散傅里叶变换83.3 二维离散傅里叶变换的性质83.4 周期延拓在卷积中的作用94 图像频率域处理程序设计步骤134.1 找出两幅大小不一的256级的灰度图像134.2 频率域处理程序设计175 运行结果及结果分析206课程设计心得体会22参考文献24摘要图像的频域处理是指根据一定的图像模型，对图像频谱进行不同程度修改的技术。二维正交变换是图像处理中常用的变换，其特点是变换结果的能量分布向低频成份方向集中，图像的边缘、线条在高频成份上得到反映，因此正交变换在图像处理中得到广泛运用。傅里叶作为一种典型的正交变换，在数学上有比较成熟和快速的处理方法。卷积特性是傅里叶变换性质之一，由于它在通信系统和信号处理中的重要地位应用最广。在用频域方法进行卷积过程中尤其要注意傅里叶变换的周期性，注意周期延拓的重要作用，本次课设将对此作详细的介绍。关键字：频域处理，二维傅里叶变换，卷积，周期延拓1 MATLAB软件功能简介MATLAB的名称源自Matrix Laboratory,1984年由美国Mathworks公司推向市场。它是一种科学计算软件，专门以矩阵的形式处理数据。MATLAB将高性能的数值计算和可视化集成在一起，并提供了大量的内置函数，从而被广泛的应用于科学计算、控制系统和信息处理等领域的分析、仿真和设计工作。MATLAB软件包括五大通用功能，数值计算功能（Nemeric）、符号运算功能（Symbolic）、数据可视化功能（Graphic）、数字图形文字统一处理功能（Notebook）和建模仿真可视化功能（Simulink）。其中，符号运算功能的实现是通过请求MAPLE内核计算并将结果返回到MATLAB命令窗口。该软件有三大特点，一是功能强大；二是界面友善、语言自然；三是开放性强。目前，Mathworks公司已推出30多个应用工具箱。 MATLAB在线性代数、矩阵分析、数值及优化、数值统计和随机信号分析、电路与系统、系统动力学、图像处理、控制理论分析和系统设计、过程控制、建模和仿真、通信系统以及财政金融等众多领域的理论研究和工程设计中得到了广泛应用。MATLAB 的应用范围非常广，包括信号和图像处理、通讯、控制系统设计、测试和测量、财务建模和分析以及计算生物学等众多应用领域。附加的工具箱（单独提供的专用 MATLAB 函数集）扩展了 MATLAB 环境，以解决这些应用领域内特定类型的问题。MATLAB应用：MATLAB 产品族可以用来进行以下各种工作：数值分析 数值和符号计算 工程与科学绘图 控制系统的设计与仿真 数字图像处理 数字信号处理 通讯系统设计与仿真 财务与金融工程 2 数字图像处理简介2.1 数字图像处理的目的数字图像处理(Digital Image Processing)是通过计算机对图像进行去除噪声、增强、复原、分割、提取特征等处理的方法和技术。数字图像处理的产生和迅速发展主要受三个因素的影响：一是计算机的发展；二是数学的发展（特别是离散数学理论的创立和完善）;三是广泛的农牧业、林业、环境、军事、工业和医学等方面应用需求的增长。一般来讲，对图像进行处理（或加工、分析）的主要目的有三个方面：（1）提高图像的视感质量，如进行图像的亮度、彩色变换，增强、抑制某些成分，对图像进行几何变换等，以改善图像的质量。（2）提取图像中所包含的某些特征或特殊信息，这些被提取的特征或信息往往为计算机分析图像提供便利。提取特征或信息的过程是模式识别或计算机视觉的预处理。提取的特征可以包括很多方面，如频域特征、灰度或颜色特征、边界特征、区域特征、纹理特征、形状特征、拓扑特征和关系结构等。（3）图像数据的变换、编码和压缩，以便于图像的存储和传输。不管是何种目的的图像处理，都需要由计算机和图像专用设备组成的图像处理系统对图像数据进行输入、加工和输出。图像增强的目的是采用某种技术手段，改善图像的视觉效果，或将图像转换到更适合人眼观察和机器分析、识别的形式，以便从图像中获取更有用的信息。图像增强的基本方法可分为两大类：空间域和频域方法。空间域是指图像平面自身，这类方法是以对图像的像素直接处理为基础的；而频率域处理技术是以修改图像的傅里叶变换为基础的。两者的具体方法包括以下内容：（1）空间域处理：点处理，模板处理即邻域处理（2）频率域处理：高、低通滤波，同态滤波等。2.2 数字图像处理的方法傅里叶变换能将满足一定条件的某个函数表示成三角函数（正弦和/或余弦函数）或者它们的积分的线性组合。在不同的研究领域，傅里叶变换具有多种不同的变体形式，如连续傅里叶变换和离散傅里叶变换。最初傅里叶分析是作为热过程的解析分析的工具被提出的1。傅里叶变换属于谐波分析。傅里叶变换的逆变换容易求出，而且形式与正变换非常类似。正弦基函数是微分运算的本征函数，从而使得线性微分方程的求解可以转化为常系数的代数方程的求解。在线性时不变的物理系统内，频率是个不变的性质，从而系统对于复杂激励的响应可以通过组合其对不同频率正弦信号的响应来获取。卷积定理指出：傅里叶变换可以化复杂的卷积运算为简单的乘积运算，从而提供了计算卷积的一种简单手段。离散形式的傅里叶变换可以利用数字计算机快速的实现（其算法称为快速傅里叶变换算法（FFT）。MATLAB中提供的变换函数（1）fft2:用于计算二维快速傅立叶变换，语句格式：B=fft2(I,m,n) 按指定的点数计算m，返回矩阵B的大小为mn，不写默认为原图像大小（2）fftn:用于计算n维快速傅立叶变换（3）fftshift:用于将变换后的图像频谱中心从矩阵的原点移到矩阵地中心，语法格式B=fftshift(I) （4）ifft2:用于计算图像的二维傅立叶反变换，语法格式：B=ifft2(i)（5）ifftn:用于计算n维傅立叶变换，快速卷积实验：傅立叶变换一个重要特性是可以实现快速卷积设A为MN矩阵，B为PQ的矩阵，快速卷积方法如下： （1）对A和B补0，使其大小都为(M+P-1)(N+Q-1) （2）利用fft2对矩阵A和B进行二维变换 （3）将两个FFT结果相乘，利用ifft2对得到的乘积进行傅立叶反变换2.3 图像频域处理的概述图像的频率是表征图像中灰度变化剧烈程度的指标，是灰度在平面空间上的梯度。如大面积的沙漠在图像中是一片灰度变化缓慢的区域，对应的频率值很低；而对于地表属性变化剧烈的边缘区域在图像中是一片灰度变化剧烈的区域，对应的频率值较高。频域处理是指根据一定的图像模型，对图像频谱进行不同程度修改的技术，通常作如下假设：1)引起图像质量下降的噪声占频谱的高频段；2)图像边缘占高频段；3)图像主体或灰度缓变区域占低频段。基于这些假设，可以在频谱的各个频段进行有选择性的修改。为什么要在频率域研究图像增强 （1）可以利用频率成分和图像外表之间的对应关系。一些在空间域表述困难的增强任务，在频率域中变得非常普通。（2）滤波在频率域更为直观，它可以解释空间域滤波的某些性质。 （3）可以在频率域指定滤波器，做反变换，然后在空间域使用结果滤波器作为空间域滤波器的指导。 （4）一旦通过频率域试验选择了空间滤波，通常实施都在空间域进行。3 二维傅里叶变换由于图像的频率是表征图像中灰度变化剧烈程度的指标，是灰度在平面空间上的梯度。傅立叶变换在实际中的物理意义，设f是一个能量有限的模拟信号，则其傅立叶变换就表示f的谱。从纯粹的数学意义上看，傅立叶变换是将一个函数转换为一系列周期函数来处理的。从物理效果看，傅立叶变换是将图像从空间域转换到频率域，其逆变换是将图像从频率域转换到空间域。换句话说，傅立叶变换的物理意义是将图像的灰度分布函数变换为图像的频率分布函数，傅立叶逆变换是将图像的频率分布函数变换为灰度分布函数。 3.1 二维连续傅里叶变换 如果二维连续函数f(x,y)满足狄里赫莱条件，则将有下面的傅立叶变换对存在：与一维傅立叶变换类似，二维傅立叶变换的傅立叶谱和相位谱为：3.2 二维离散傅里叶变换一个MN大小的二维函数f(x,y)，其离散傅立叶变换对为 ：在数字图像处理中，图像一般取样为方形矩阵，即NN，则其傅立叶变换及其逆变换为 ：3.3 二维离散傅里叶变换的性质 离散傅里叶变换主要有以下性质：1.平移性质2.分配律3.尺度变换（缩放）4.旋转性5.周期性和共轭对称性6.平均值7.可分性8.卷积9.相关性。这里主要简述周期性，卷积相关内容会在下一节中介绍。离散傅里叶变换有如下周期性性质：反变换也是周期性的：频谱也是关于原点对称的：这些等式的有效性是建立在二维离散傅里叶变换公式基础上的。图像的周期性在图像处理中有非常重要的作用，下面会在卷积部分继续阐述周期性的相关内容。3.4 周期延拓在卷积中的作用基于卷积理论，频率域的乘法相当于空间域的卷积，反之亦然。当处理离散变量和傅里叶变换时，要记住不同函数所包含的周期性。虽然可能不太直观，但周期性是定义离散傅里叶变换对时产生的数学副产品。周期性是处理操作的一部分,不应忽视。图1列举了周期性的重要性。图1 傅里叶变换周期性图1左边（ae）：两个离散函数的卷积，右边（fj）：相同函数的卷积，考虑DFT周期性的应用。图的左边一列是用下式的一维形式计算的卷积：在此详细地解释卷积运算的过程。为简化表示，简单的数字将代替那些表示函数长度和高度的通用符号。图1(a)和(b)是两个要进行卷积的函数。每个函数包含400个点。卷积的第一步是将一个函数关于原点进行镜像映射(倒转)，在本例情况下，对第二个函数进行，在图1(c)中以h(-m)示出。下一步是将h(-m)滑过f(m)。这要增加一个常数x到h(-m)，即变成h(x-m)，如图1(d)所示。注意只有一个置换值。在第一次遇到时这个简单步骤通常是引起混乱的根源。而这恰好是卷积计算的全部关键。换言之，为了执行卷积，倒转了一个函数，并将它滑过另一个函数。在每一个置换点(的每一个值)都要计算式的全部总和。这个总和不比在给定位移处f和h乘积的和更太。位移x的范围为h完全滑过f需要的所有值。图1(e)显示了h完全滑过f后的结果，并在x的每个点计算式。在此例中，为使h(x-m)完全滑过f,x值的范围是从0到799。这幅图是两个函数的卷积，要清楚地记住卷积中的变量是x.从上面介绍的卷积理论可知，由F(u)H(u)的傅里叶反变换能得到同样的准确结果。但是，从前面对周期性的讨论又知离散傅里叶变换自动地将输入函数周期化。换言之，采用DFT允许在频率域进行卷积计算，但函数必须看做周期性的，且周期等于函数的长度。可以通过图1右边一列考察这种隐含的周期性。图1(f)同图1(a)一样，但同样的函数在两个方向上周期性地无限扩展(扩展部分用虚线表示)。从图1(g)到图1(i)同样应用该扩展。现在，可以通过将h(x-m)滑过f(m)进行卷积。如前面一样，变化x完成滑动。然而，h(x-m)的周期性扩展产生了图1左边的计算中所没有的值。例如，在图1(i)中，当x=0时，看到h(x-m)右侧第一个扩展周期的一部分进 入图1(f)中所示的f(m)(从原点开始)的一部分。当h(x-m)向右滑动时，在f(m)中的那部分开始向右侧移出，但被h(x-m)左侧相同部分所取代。这引起卷积产生一个常量值，如图1(j)所示的0，100的一段.从100到4OO的一段是正确的，但周期性是周而复始的，这样就引起卷积函数尾部的一部分丢失，由图1(j)和图1(e)实线部分的比较可以看出这一点。在频率域，该过程需要计算图1(a)和(b)中函数的傅里叶变换。根据卷积理论，两个变换要相乘，再计算傅里叶反变换。结果包含40O个点的卷积，如图1(j)的实线部分所示。简单的解释表明当使用傅里叶变换得出卷积函数时，错误地处理周期性将得到错误的结论。结果，在开头有错误数据，结尾将丢失数据。问题的解决办法很简单。假设f和h分别由A和B个点组成。对两个函数同时添加零，以使它们具有相同的周期，表示为P。这个过程产生扩展的或延拓的函数，如下所示：和可以看出,除非选择PA+B-1，否则卷积的独立周期将会混叠。已经在图1中看到了这种现象的结果，这通常归于缠绕误差。若P=A+B-1，周期便会邻接起来。若PA+B-1，周期将会是分隔开的，分隔的程度等于P与A+B-1的差。 扩展后的卷积结果如图2所示。在这里，选择P=A+B-1(799)，即可知卷积周期是相邻的。遵循与前面的解释相同的过程，得到如图2所示的卷积函数。该结果的一个周期与图1(e)相同，是正确的。这样，如果要在频率域计算卷积，应该：(1)得到两个扩展序列的傅里叶变换(每个序列有8OO个点)；(2)将两个变换相乘；(3)计算傅里叶反变换。结果便得到正确的8OO个点的卷积函数。见图2中周期加重的部分。图2 卷积函数这些概念扩展到二维函数时遵循了相同的前提。假设有f(x,y)和h(x,y)两幅图像，大小分别为AB和CD。如同一维情况，这些行列必须假定在x方向上有相同的周期P，在y方向上有相同的周期Q。二维卷积的混叠可由选择如下周期避免：扩展f(x,y)和h(x,y)形成如下周期性序列：为了简化图例，假设f和h是方形的，且大小相同。图3 二维函数周期延拓图3对二维函数周期延拓的说明。图3（a）没有延拓执行二维卷积的结果；图3（b）合格的函数延拓；图3（c）正确的卷积结果。图3(a)显示了图像没有延拓时得到的滤波结果。这通常是由于没有对一幅输入图像进行延拓就进行傅里叶变换，然后又乘上同样大小的函数(也没有延拓)，计算傅里叶反变换。结果就是与输入图像相同的大小为AB的图像，如图3(a)左上象限所示。如同一维情况，图像前面边沿(阻影部分)由于周期性而引入了错误数据，而在尾部边沿将丢失数据。如图3(b)所示，通过对输入图像和函数进行合适的延拓，将得到正确的、大小为PQ的过滤图像。这幅图像在两个坐标方向上是原始图像的两倍大小，有原始图像4倍数量的像素点。4 图像频率域处理程序设计步骤4.1 找出两幅大小不一的256级的灰度图像选择两幅图片如图4图5 图4 源图片1 图5源图片2（1）检查源图片1和源图片2格式info=imfinfo(1.jpg)info=imfinfo(2.jpg)显示如下info = Filename: 1.jpg FileModDate: 15-Aug-2014 13:44:02 FileSize: 51822 Format: jpg FormatVersion: Width: 440 Height: 737 BitDepth: 24 ColorType: truecolor FormatSignature: NumberOfSamples: 3 CodingMethod: Huffman CodingProcess: Sequential Comment: info = Filename: 2.jpg FileModDate: 08-Jan-2015 08:38:58 FileSize: 16439 Format: jpg FormatVersion: Width: 300 Height: 400 BitDepth: 8 ColorType: grayscale FormatSignature: NumberOfSamples: 1 CodingMethod: Huffman CodingProcess: Sequential Comment: 用MATLAB检查发现是图片1是truecolor 格式，ColorType: truecolor。（2）用MATLAB将其转换为灰度图像a=imread(1.jpg);b=rgb2gray(a);imwrite(b,3.jpg) info=imfinfo(3.jpg)显示如下info = Filename: 3.jpg FileModDate: 08-Jan-2015 08:54:41 FileSize: 43637 Format: jpg FormatVersion: Width: 440 Height: 737 BitDepth: 8 ColorType: grayscale FormatSignature: NumberOfSamples: 1 CodingMethod: Huffman CodingProcess: Sequential Comment: 在 ColorType: grayscale 属性行如此显示，则说明以3.jpg命名的文件为灰度图像如图6所示。图6 图片3（3）检查图像灰度级在命令窗口输入f=imread(2.jpg)whos fg=imread(3.jpg)whos g显示如下，表明图像为256级灰度图像Name Size Bytes Class f 400x300 120000 uint8 arrayName Size Bytes Class g 737x440 324280 uint8 arrayGrand total is 120000 elements using 120000 bytes4.2 频率域处理程序设计MATLAB中提供的变换函数（1）fft2:用于计算二维快速傅立叶变换，语句格式：B=fft2(I,m,n)按指定的点数计算m,返回矩阵B的大小为mn,不写默认为原图像大小。（2）ifft2:用于计算图像的二维傅立叶反变换，语法格式：B=ifft2(i)%直接卷积程序 I1=imread(2.jpg);I2=imread(3.jpg);I5=conv2(I1,I2);figure(3);imshow(I5,);title(直接函数卷积得到的图像 (黄深)%正确的频域处理程序I1=imread(2.jpg);I2=imread(3.jpg);m1,n1=size(I1);m2,n2=size(I2);I1(m1+m2-1,n1+n2-1)=0;I2(m1+m2-1,n1+n2-1)=0;I3=ifft2(fft2(I1).*fft2(I2);I3=I3(1:m1+m2-1,1:n1+n2-1);I3=real(I3);figure(1);imshow(I3,);title(正确延拓频域法得到的卷积图像(黄深)%比较频域方法与直接卷积的结果，显示差矩阵并且显示错误数据数F=minus(I3,I5);figure(4)imshow(F);title(正确延拓差矩阵的二值图像(黄深)s=0;for i=1:m1+m2-1 for j=1:n1+n2-1 if (minus(abs(F(i,j),0.000001)0) s=s+1; end; end;end;disp(sprintf(差错 1: %d,s);%补0不够的频域处理程序I1=imread(2.jpg);I2=imread(3.jpg);m1,n1=size(I1);m2,n2=size(I2);I1(m1+m2-100,n1+n2-100)=0;I2(m1+m2-100,n1+n2-100)=0;I3=ifft2(fft2(I1).*fft2(I2);I3=I3(1:m1+m2-100,1:n1+n2-100);I3=real(I3);I3(m1+m2-1,n1+n2-1)=0;figure(2);imshow(I3,);title(补0不够频域法得到的卷积图像(黄深)%比较频域方法与直接卷积的结果，显示差矩阵并且显示错误数据数F=minus(I3,I5);figure(5)imshow(F);title(补0不够的差矩阵的二值图像(黄深)s=0;for i=1:m1+m2-100 for j=1:n1+n2-100 if (minus(abs(F(i,j),0.000001)0) s=s+1; end; end;end;disp(sprintf(差错 2: %d,s);5 运行结果及结果分析在MATLAB中输入程序后，显示的卷积结果如下: 图7 正确延拓频域法得到的卷积图像 图8补0不够频域法得到的卷积图像 图9直接函数卷积得到的图像 图10正确延拓差矩阵的二值图像图11 补0不够的差矩阵的二值图像比较图7和图9，看不出两个图像有任何区别。通过作差，认为舍入误差小于0.000001的均可作为0来处理，这里S=785，差值矩阵的二值图像全为黑，可以认为两图几乎没有任何区别，即频域方法的卷积结果是完全正确的。比较图8和图9，表面上也看不出两个图像有什么区别，图8的靠左和靠上部分有亮度增加，这部分是叠加错误，而靠下和靠右部分是两条黑杠，这是补零的数据，也就是原来丢失的数据。通过检测差值矩阵，S=99099,错误的有很多，即没有补0的频域方法计算的结果不正确。值得注意的是这里差值矩阵应该四周都是白色，因为左边和上边是混叠错误的地方应该为，行数:100，列数100；同理右边和下边是数据丢失人为补0的地方也有与混叠相同的行数和列数。但因为这里gyy.jpg周围为0，因此正确卷积的结果也为0，因此差矩阵得到的相应区域也为0，显示的2值图像就看不到白色地方了。通过以上分析说明，二维图像或矩阵的线性卷积可以通过补零周期延拓后，经二维傅里叶变换相乘，再做反变换来实现。而不补零或补零不足，用此方法求得卷积图像靠左靠上会有叠加错误和靠下靠右会有数据丢失。6课程设计心得体会通过本次课程设计，我主要了解了用MATLAB处理数字图像中的问题，特别是数字图像的频域分析法，快速傅里叶变换，把数字图像技术和MATLAB等通信类科目的内容应用到本课程设计中来，进一步巩固复习数字图像处理技术，MATLAB等课程，以达到融会贯通的目的。本次课设的主要任务是用MATLAB编程来实现数字图像的傅里叶计算，卷积原理。开始我对数字图像技术了解不是很多，通过查阅相关资料，我熟悉了用MATLAB处理数字图像基础，通过查阅资料学会了如何用MATLAB编程进行图像的处理，感觉非常有成就感。并且加深了对数字图像处理的认识，经