图像的变换域处理及应用.doc

第9章 图像的变换域处理及应用本章要点： 图像的正交变换 频域低通滤波 频域高通滤波9.1 概述数字图像处理的方法很多，根据它们处理数字图像时所用系统，主要可以归纳为两大类：空间域处理法（空域法）及频域法（或称为变换域法）。前面几章所介绍的几何变换，图像的增强、边缘检测等所用算法都是在空间域中进行图像处理的，本章将着重介绍数字图像处理中一些常见的频域处理方法。数字图像处理经常要用到线性系统，在图像处理中使用空间作为参数来描述，通常用二维系统进行表示，输入函数f(x,y)表示原始图像，输出函数g(x,y)表示经处理后的图像，线性系统可看作是输入函数和输出函数之间的一种映射，反映了各种线性的图像处理方法。关系如同公式：g(x,y)=f(x,y) (9-1)一般数字图像处理的计算方法本质上都是线性的，处理后的输出图像阵列就是输入图像阵列中的各个元素经加权线性组合而得到，通常这种线性空间线性处理要比非线性处理容易理解并且算法简单。线性系统可用传递函数来刻画，将其看作黑箱（Black Box），线性系统的输入信号和输出信号之间的关系，在时域可用卷积运算来表达，在频域可直接用乘积来确定：传递函数(冲击响应)h(t)H(s) x(t) y(t)X(s)Y(s)输入信号 输出信号时域关系：y(t)=h(t)*x(t)=x()h(t-) d (9-2)频域关系：Y(s)=H(s)X(s) (9-3)在图像处理中，图像的锐化与平滑处理可采用空间域处理方式（又称空间滤波）和频域处理方式（又称频域变换）两类。从数学角度看，空间滤波是采用微分、积分、多项式运算、坐标变换等方法对图像进行某种形式的处埋，具有方法直观，制作简便等优点；但当要处理较大的数字图像数据时，由于图像阵列很大，如果没有发现比较高效的算法，计算上会变得很繁琐，存在着滤波的广度和构成方式的模糊，计算时间长，预测性差等缺点，这样就会降低其在现实工作中的实用价值。同样情况下如果采用图像变换的方法，如傅里立叶算法、沃尔夫算法等间接处理技术，就可以获得更为有效的处理方法。所谓的图像频谱变换则是将图像从空间域进行付里叶变换于频谱域，检测和研究图像频谱特性，并进行滤波处理，最终将处理的频谱经傅里叶逆变换恢复图像于空间域。如下图所示。其优点是处理速度快，构成方式清晰，滤波广度大，预测性好，但数学过程复杂，不易理解。f（x,y） F（u,v） H（u,v） G（u,v） g（x,y）付立叶变换滤波付立叶反变换图中，F(u,v)是带噪声的原始图像f(x,y)的付立叶变换，H(u,v)为滤波器的传递函数，经过滤波处理后的G(u,v)=H(u,v)*F(u,v), 再进行付立叶反变换得到增强的图像g(x,y)。当H(u,v)为低通滤波器的传递函数时，经过付立叶反变换会得到去除噪声后的平滑图像g(x,y)。当H(u,v)为高通滤波器的传递函数时，经过付立叶反变换会得边缘增强的图像，衰减图像信号的低频部分能相对增强图像高频部分，从而实现图像锐化的目的。目前，图像变换技术被广泛地运用于图像增强、图像复原、图像压缩、图像特征提取、图像识别以及图像特征提取等领域。本章将重点介绍这些与图像变换相关的算法。9.2 图像的正交变换在将数字图像由空间域变换到频域时，所采用的变换方式一般都是线性正交变换，又称为酉变换。正交变换是信号分析学科中的一个重要部分，它是计算机图像处理的前续课程。多年来，变换理论在图像处理（频域法处理）中起着关键作用。下面我们将介绍使用正交变换的傅立叶变换、离散余弦变换和沃尔什变换。 傅立叶变换基本概念 一维离散傅立叶变换 二维离散傅立叶变换 离散余弦变换9.2.1 傅立叶变换基本概念 傅立叶变换是一种经常被使用的正交变换，尤其是在一维信号处理中被广泛使用。在这里我们将介绍它在数字图像处理中的使用方法。1傅立叶的定义傅立叶变换在数学中的定义非常严格，它的定义如下：设为的函数，如果满足下面的狄里赫莱条件：(1) 具有有限个间断点；(2) 具有有限个极值点；(3) 绝对可积。则定义的傅立叶变换公式为： (9-4)它的逆反变换公式为： (9-5)其中为时域变量，为频域变量。由上面的公式可以看出，傅立叶变换结果是一个复数表达式。设的实部为，虚部为，则： (9-6)或者写成指数形式： 其中： (9-7) (9-8)通常把称作的傅立叶幅度谱，为的相位谱。我们可以把傅立叶变换推广到二维情况。如果二维函数满足狄里赫莱条件，那么将可以导出下面的二维傅里叶变换： (9-9) (9-10)同样，二维傅立叶变换的幅度谱和相位谱为： (9-11) (9-12)可以定义： (9-13)通常称为能量谱。2傅立叶变换的性质傅立叶变换具有很多方便运算处理的性质。下面列出二维傅立叶变换的一些重要性质。(1)线性傅立叶变换是一个线性变换，即： (9-14)(2)可分离性一个二维傅立叶变换可以用二次一维傅立叶变换来实现。推导如下： (9-15)(3)平移性傅立叶变换具有平移特性，推导公式如下： (9-16) (9-17)(4)共轭性如果函数的傅立叶变换为,为傅立叶变换的共轭函数，那么： (9-18)(5)尺度变换特性如果函数的傅立叶变换为，和为两个标量，那么： (9-19) (9-20)(6)旋转不变性如果空间域函数旋转角度为，则在变换域中该函数的傅立叶变换函数也将旋转同样的角度。表达式如下： (9-21)公式（919）中，和为极坐标表达式，其中，的傅立叶变换结果为。(7)对称性如果函数的傅立叶变换为，那么： (9-22)(8)能量保持定理能量保持定理也称帕斯维尔(Parseval)定理，该定理的数学描述如下： (9-23) 该公式表明了傅立叶变换前后能量守恒。(9)相关定理如果和为两个二维时域函数，那么可以定义相关运算如下： (9-24)则： (9-25) (9-26)其中为函数的傅立叶变换，为函数的傅立叶变换；为的共轭，为的共轭。(10)卷积定理如果和为两个二维时域函数，那么可以定义卷积运算*如下： (9-27)则： (9-28) (9-29)其中为函数的傅立叶变换，为函数的傅立叶变换。9.2.2 一维离散傅立叶变换1概念连续函数的傅里叶变换是波形分析的有力工具，但是为了使之用于计算机技术，必须将连续变换转变成离散变换，这样就必须引入离散傅里叶变换(DFT,Discrete Fourier Transform)的概念。离散傅立叶变换在数字信号处理和数字图像处理中都得到了十分广泛的应用，它在离散时域和离散频域之间建立了联系。如果直接应用卷积和相关运算在时域中处理，计算量将随着取样点数N的平方而增加，这使计算机的计算量迅速增大，耗时增多，很难达到对数字图像实际处理的要求。因此，一般可采用离散傅立叶变换方法，将输入的数字信号首先进行频域处理，再利用离散时域与离散频域之间的联系，将在离散频域中处理的效果反馈给离散时域，这样就比在时域中直接对数字图像处理变得更加快捷便利，计算量也会大大减少，同时提高数字图像的处理速度，增强算法的实用性。因此，离散傅立叶变换在数字图像处理领域中有很大的使用价值。离散傅立叶变换还有一个明显的优点就是具有快速算法，即快速傅立叶算法(Fast Fourier Transform)，它可以大大减少计算次数，使计算量减少到只是相当于直接使用离散傅立叶变换所用的一小部分。并且，二维离散傅立叶变换很容易从一维的概念推广得到。在数字图像处理中，二维离散傅立叶被广泛的应用于图像增强、复原、编码和分类中。如果为一长度为N的数字序列，则其离散傅里叶正变换定义由下式来表示： (9-30)傅里叶反变换定义由下式来表示： (9-31)其中：=0，1，2N-1如果令，那么上述公式变成： (9-32) (9-33)公式(932)写成矩阵形式为： (9-34)公式(933)写成矩阵形式为： (9-35)2快速傅立叶变换的实现现在，离散傅里叶变换已成为数字信号处理的重要工具，但是它的计算量较大，运算时间长，在某种程度上限制了它的使用。为了解决这一矛盾，引用了快速傅里叶变换的思想。快速傅里叶变换并不是一种新的变换方式，它是离散傅里叶变换的一种算法，这种方法是建立在分析离散傅里叶变换中的多余运算的基础上，进而消除这些重复工作的思想指导下得到的，从而在运算中节省了大量的计算时间，达到快速运算的目的。快速傅立叶算法以的组成状况可以分成为2的整数幂的算法；为高复合数的算法；为素数的算法三种情况。这里介绍第一种算法。令 (9-36)一维离散傅立叶变换公式变为 (9-37) 分别为。再令 在此基础上，将分解成为和对应的偶数和奇数两部分，的取值范围由原来的0到改为0到。下面我们按照奇偶来将序列进行划分，设： 因此，离散傅立叶变换可以改写成下面的形式： (9-38) 因此，一个求点的离散傅立叶变换可以被转换成为两个求点的离散傅立叶变换。以的离散傅立叶变换为例，利用公式(9-38)可得： (9-39)由于和都是4点的离散傅立叶变换，所以它们均以4为周期。因此，。再加上的对称性，因此公式(9-39)可以改进为： (9-40)如下面的蝶形运算单元图所示，定义由、和所构成的结构为蝶形运算单元，它的左方两个节点为输入节点，代表输入数值；右方两个节拍为输出节点，表示输入数值的叠加。运算由左向右进行。线旁的和为加权系数。 图9-1 蝶形运算单元图图9-1所示运算为： 由于和的过程仍然是常规的离散傅立叶变换，因此也可以按照上述规则进行分解，记为计算的DFT运算，记为计算的DFT运算，记为计算的DFT运算，记为计算的DFT运算。 (9-41)由于、和现在已经是2点DFT了，因此不能继续分解，可直接由原始数据计算得出： (9-42)至此，可以进一步写出8点DFT的完整FFT计算的流程框图， 图9-2 FFT计算的流程框图从上可见和仅仅在运算符上不同,同理可以推出与、与、与在运算符上不同。为了快速实现傅立叶变换，要对进行“逆序”重排。例如： ，其逆序形式为。数组由形式变成形式称为“逆序”。如已知，则的“逆序”排列将如下表9-1所示： 表9-1 时数组的“逆序”排列表有了“逆序”的概念，就可以着手讨论快速傅立叶变换的实现。我们仍以为例，计算时，先将原始数组按宗量值“逆序”顺序，重新排列见上表，然后自左向右，每两个相邻元素位一组，共4组，具体是：、；、；、；、。在每一组中，第一个元素作为偶元素，第二个元素作为奇元素（注意对所有组都一样）。在此基础上完成四个“两点变换”。下一步是利用“两点变换”结果形成两个“四点变换”，最后是建立在“四点变换”基础上的一个“八点变换”。具体见图9-3。 两点变换四点变换八点变换图9-3快速傅立叶变换以上快速傅立叶变换流程简称为逐次加倍法，因为“两点变换”由两个“一点变换”算出，“四点变换”由两个“两点变换”算出，“八点变换”由两个“四点变换”算出，依此类推。通过这种计算方法快速傅立叶变换比二维离散傅立叶变换的计算量大量减少，效率大大提高，增强了算法的实用性，使傅立叶变换在数字图像处理中得到更为广泛的应用。3. 快速傅立叶变换实现步骤(1) 开辟存储空间用以保存加权系数Wi及中间变量。(2) 采用分解法进行蝶形运算。(3) 重新排列序列顺序。(4) 释放存储空间。4. 编程代码struct CNumberdouble re;double im;/* *函数名称：QFC(CNumber* t,CNumber* f,int r)*参 数: t、f分别是指向时域和频域的指针，r是2的幂数 *函数类型：void*功能：此函数实现快速傅立叶变换 */void ZhengJiaoBianHuanDib:QFC(CNumber* t,CNumber* f,int r)long count;/傅立叶变换点数int i,j,k,p,bfsize;CNumber *w,*x,*a,*b;/复数结构类型的指针变量，其中w指向加权系数double angle;/计算加权系数所用的角度值count=1r;/计算傅立叶变换点数/分配所需运算空间w=(CNumber*)malloc(sizeof(CNumber)*count/2);a=(CNumber*)malloc(sizeof(CNumber)*count);b=(CNumber*)malloc(sizeof(CNumber)*count);/计算加权系数for(i=0;icount/2;i+)angle=-i*pi*2/count;wi.re=cos(angle);wi.im=sin(angle);memcpy(a,t,sizeof(CNumber)*count);/采用频率分解法进行蝶形运算for(k=0;kr;k+)for(j=0;j1k;j+)bfsize=1(r-k);for(i=0;ibfsize/2;i+)p=j*bfsize;bi+p=Add(ai+p,ai+p+bfsize/2);bi+p+bfsize/2=Mul(Sub(ai+p,ai+p+bfsize/2),wi*(1k);x=a;a=b;b=x;/将乱序的变换序列重新排序for(j=0;jcount;j+)p=0;for(i=0;ir;i+)if(j&(1i)p+=1GetData();/指向原图像数据区width=this-GetWidth();/得到图像宽度height=this-GetHeight();/得到图像高度 long lLineBytes=WIDTHBYTES(width*8);/计算图像每行的字节数while(w*2=width)/计算进行傅立叶变换的宽度（2的整数次方）w*=2;wp+;while(h*2=height)/计算进行傅立叶变换的高度（2的整数次方）h*=2;hp+; t=(CNumber*)malloc(sizeof(CNumber)*w*h);/分配存储器空间 f=(CNumber*)malloc(sizeof(CNumber)*w*h);for(j=0;jh;j+)for(i=0;iw;i+)p=p_data+lLineBytes*(height-j-1)+i;/指向第i行第j列像素 ti+w*j.re=*(p);/给时域赋值ti+w*j.im=0;for(j=0;jh;j+)/在垂直方向上进行快速傅立叶变换QFC(&tw*j,&fw*j,wp);for(j=0;jh;j+)/转换变换结果for(i=0;iw;i+)tj+h*i=fi+w*j;for(j=0;jw;j+)/水平方向进行快速傅立叶变换QFC(&tj*h,&fj*h,hp);for(j=0;jh;j+)for(i=0;i255)temp=255;p=p_data+lLineBytes*(height-(jh/2?j+h/2:j-h/2)-1)+(iw/2?i+w/2:i-w/2);/将变换后的原点移到中心p=p_data+lLineBytes*(height-(jh/2?j:j)-1)+(iw/2?i:i);*(p)=(BYTE)(temp);free(t);free(f);5效果对比图 (a)原图 (b)二维傅立叶变换效果图图9-7二维傅立叶变换效果图6应用傅立叶变换时应当注意的问题尽管傅立叶变换提供了很多有用的属性，在数字图像处理领域中得到广泛的应用，但是它也有自身的不足，主要表现在两个方面。一是需计算复数，而进行复数运算相对比较费时。如采用其它合适的完备的正交函数来代替傅里叶变换所用的正、余弦函数构成完备的正交函数系，就可避免这种复数运算。如后面介绍的沃尔什(Walsh)函数系,每个函数只取+1和-1两个值，组成二值正交函数。因此，以沃尔什函数为基础所构成的变换，是实数加减运算，其运算速度要比傅立叶变换快。傅立叶变换的另一个缺点是收敛慢，这在图像编码应用中尤为突出。9.2.4 离散余弦变换前面介绍的傅立叶变换作为可分离变换中一种经典的变换，在运算方式上与其它几种可分离变换有相当大的差别：傅立叶变换计算的对象是复数；而采用了其它完备正交函数系的一些可分离变换(如离散余弦变换、沃尔什变换等)则是以实数为对象的余弦函数。虽然此类变换没有傅立叶变换的功能强大，但是离散余弦变换（Discrete Cosine Transform）的计算速度要比对象为复数的离散傅立叶变换块得多，并且已经被广泛应用到图像压缩编码、语音信号处理等众多领域。1 一维离散余弦变换一维离散余弦变换的定义可以用下式表示： (9-83) (9-84)式中是第个余弦变换系数，是广义频率变量，；是时域点序列。一维离散余弦变换的定义由下式表示： (9-85)（1）实现步骤计算离散余弦变换点数。对时域空间进行延拓。将时域点写入已开辟存储空间。调用一维快速傅里叶变换（过程见傅里叶变换的实现）。调整系数。转换变换结果，将变换结果转存回时域存储区。（2）编程代码/* LiSan()* 参数:* double * t- 指向时域值的指针* double * f- 指向频域值的指针* r- 2的幂数* 返回值: 无* 说明:* 该函数用来实现快速离散余弦变换。该函数利用2N点的快速付立叶变换* 调用QFC来实现离散余弦变换。*/void ZhengJiaoBianHuanDib:LiSan(double *t, double *f, int r)/ 离散余弦变换点数longcount;/ 循环变量inti;/ 中间变量doubledTemp;CNumber *X;/ 计算离散余弦变换点数count = 1r;/ 分配内存X=(CNumber*)malloc(sizeof(CNumber)*count*2);/ 赋初值为0memset(X, 0, sizeof(CNumber) * count * 2);/ 将时域点写入数组Xfor(i=0;icount;i+)Xi.re=ti;Xi.im=0;/ 调用快速付立叶变换QFC(X,X,r+1);/ 调整系数dTemp = 1/sqrt(count);/ 求F0f0 = X0.re * dTemp;dTemp *= sqrt(2);/ 求Fufor(i = 1; i count; i+)fi=(Xi.re * cos(i*PI/(count*2) + Xi.im * sin(i*PI/(count*2) * dTemp;/ 释放内存delete X;2二维离散余弦变换二维离散余弦变换的定义由下式表示：= (9-86)= 其中为空间域中二维向量，为变换系数矩阵，。由上可知，二维离散余弦的变换核是可分离的可通过两次一维变换实现。传统的方法是行-列法，即先沿行（列）进行一维离散余弦变换计算，再沿列（行）计算一维离散余弦变换。其离散余弦反变换由下式表示：+ + (9-87)傅立叶变换中指数项通过尤拉公式和进行分解，其傅立叶变换实数部分对应于余数项，其虚数部分对应于正弦项，因此，离散余弦变换可以从傅立叶变换的实数部分求得，即离散余弦变换可以改写成以下形式： (9-88)式中：。其中，对于有，代表花括号项取实数部分，其中求和项就是个点上的离散余弦变换。为了提高计算效率，在实际应用中需要找一种快速算法。（1）离散余弦变换的编程方法首先，我们将进行延拓： (9-89)由于为的点离散傅立叶变换。因此，再做离散余弦变换时，可以把长度为的序列的长度延拓为的序列，然后对延拓的结果进行离散傅立叶变换，最后取离散傅立叶变换的实部便是离散余弦变换的结果。同理，对于离散余弦反变换，也可以按照下式延拓： (9-90)按照公式（9-87）可得： (9-91) 从上式可见离散余弦反变换可以由的点的反傅立叶变换来实现。（2）实现步骤获取复制来的原图像的数据区首地址、图像的高度和图像的宽度。计算进行离散余弦变换的宽度和高度，这两个值必须是2的整数次方；计算变换时所用的迭代次数，包括水平方向的和垂直方向的。依行列顺序依次读取数据区的值，存储到开辟的复数存储区。调用一维离散余弦变换函数进行垂直方向的变换。调用一维离散余弦变换函数，在水平方向上进行离散余弦变换。将计算结果转换成可显示图像。（3）编程代码/* DIBLiSanYuXuan()* 参数:* LPBYTE lpDIBBits - 指向源DIB图像指针* LONG lWidth - 源图像宽度（像素数）* LONG lHeight - 源图像高度（像素数）* 返回值:* BOOL - 成功返回TRUE，否则返回FALSE。* 说明:* 该函数调用前面的LiSan（）函数来实现对图像的离散余弦变换。*/bool ZhengJiaoBianHuanDib:DIBLiSanYuXuan(LPBYTE lpDIBBits, LONG lWidth, LONG lHeight)/ 指向源图像的指针unsigned char* lpSrc;/ 循环变量LONGi;LONGj;/ 进行付立叶变换的宽度和高度（2的整数次方）LONGw;LONGh;/ 中间变量doubledTemp;intwp;inthp;/ 图像每行的字节数LONGlLineBytes;/ 计算图像每行的字节数lLineBytes = WIDTHBYTES(lWidth * 8);/ 赋初值w = 1;h = 1;wp = 0;hp = 0;/ 计算进行离散余弦变换的宽度和高度（2的整数次方）while(w * 2 = lWidth)w *= 2;wp+;while(h * 2 = lHeight)h *= 2;hp+;/ 分配内存double *f = new doublew * h;double *F = new doublew * h;/ 行for(i = 0; i h; i+)/ 列for(j = 0; j w; j+)/ 指向DIB第i行，第j个像素的指针lpSrc = (unsigned char*)lpDIBBits + lLineBytes * (lHeight - 1 - i) + j;/ 给时域赋值fj + i * w = *(lpSrc);for(i = 0; i h; i+)/ 对y方向进行离散余弦变换LiSan(&fw * i, &Fw * i, wp);/ 保存计算结果for(i = 0; i h; i+)for(j = 0; j w; j+)fj