您的当前位置:首页正文

对二维FFT的详细介绍及实验

2021-04-02 来源:榕意旅游网
matlab fft fftshift

(2009-04-15 10:32:51)

转载▼ 标签:

it

在图象处理的广泛应用领域中,傅立叶变换起着非常重要的作用,具体表现在包括图象分析、图象增强及图象压缩等方面。

fftshift的作用正是让正半轴部分和负半轴部分的图像分别关于各自的中心对称。因为直接用fft得出的数据与频率不是对应的,fftshift可以纠正过来

假设f(x,y)是一个离散空间中的二维函数,则该函数的二维傅立叶变换的定义如下:

p=0,1…M-1 q=0,1…N-1 (1)

或 p=0,1…M-1 q=0,1…N-1 (2)

离散傅立叶反变换的定义如下:

m=0,1…M-1 n=0,1…N-1(3)

F(p,q)称为f(m,n)的离散傅立叶变换系数。这个式子表明,函数f(m,n)可以用无数个不同频率的复指数信号和表示,而在频率(w1,w2)处的复指数信号的幅度和相位是F(w1,w2)。

例如,函数f(m,n)在一个矩形区域内函数值为1,而在其他区域为0,如图所示。

了简便起见,假设f(m,n)为一个连续函数,则f(m,n)的傅立叶变换的幅度值(即 )显示为网格图,如图所示。

将傅立叶变换的结果进行可视化的另一种方法是用图象的方式显示变换结果的对数幅值 ,如图所示。

几种简单函数的傅立叶变换的频谱可以直观的表示为图所示的样子。

2、MATLAB提供的快速傅立叶变换函数 (1fft2

fft2函数用于计算二维快速傅立叶变换,其语法格式为: B = fft2(I)

B = fft2(I)返回图象I的二维fft变换矩阵,输入图象I和输出图象B大小相同。

例如,计算图象的二维傅立叶变换,并显示其幅值的结果,如图所示,其命令格式如下

load imdemos saturn2

imshow(saturn2)

B = fftshift(fft2(saturn2));

imshow(log(abs(B)),[],'notruesize')

(2)fftshift

MATLAB提供的fftshift函数用于将变换后的图象频谱中心从矩阵的原点移到矩阵的中心,其语法格式为:

B = fftshift(I)

对于矩阵I,B = fftshift(I)将I的一、三象限和二、四象限进行互换。

(2)ifft2

ifft2函数用于计算图象的二维傅立叶反变换,其语法格式为:

B = ifft2(I)

B = ifft2(A)返回图象I的二维傅立叶反变换矩阵,输入图象I和输出图象B大小相同。其语法格式含义与fft2函数的语法格式相同,可以参考fft2函数的说明。

3、简单低通滤波器的设计

一个图象经过傅立叶变换后,就从空域变到了频域,因此我们可以用信号处理中对于频域信号的处理方法对一幅图象进行处理。比如对图象进行低通滤波等。

一个二维的理想低通滤波器(ILPF),它的传递函数由下式确定:

若 (4)

0 若

式中D0是一个规定的非负的量,称为截止频率,虽然在计算机中必定能够模拟一个锐截止频率的理想低通滤波器,但它们不能用电子元件来实现。实际中比较常用的低通滤波器有:巴特沃思(Butterworth)滤波器、指数滤波器(ELPF)、梯形低通滤波器等。

在实验中我们设计一个理想的低通滤波器。

设计理想的低通滤波器由其定义可知只要设计一个与频域图象大小完全相同的矩阵。在某一个域值内该矩阵的值为1,其余为0即可。

例:若图象的大小为128*128,则可以这样设计一个低通滤波器:

H=zeros(128);

H(32:96,32:96)=1; %此处的范围是人为取定的,可以根据需要更改。

若图象矩阵I的傅立叶变换是B(已经用fftshift将频谱中心移至矩阵的中心),则对这幅图象做低通滤波,再做傅立叶逆变换命令为

LOWPASS=B.*H; %此处设变换后的矩阵为LOWPASS,另注意这儿是矩阵的点乘。

C=ifft2(LOWPASS);

Imshow(abs(C))

*******************************************************************************象富里哀级数,富里哀变换以及它们离散时间相应部分构成了信号处理的基础。为了便于这类问题的分析,MATLAB提供了函数fft,ifft,fft2,ifft2和fftshift。这类函数集执行一维和二维离散富里哀变换及其逆变换。这些函数允许人们完成很多信号处理任务。除此之外,还可在可选的信号处理工具箱中得到其他扩展的信号处理工具。

因为信号处理包含如此广泛的领域,甚至要说明用MATLAB中离散富里哀变换函数可解决的这类小问题,就超出了本书的范围。因此,这里将只介绍用函数fft近似连续时间信号的富里哀变换的一个例子。此外,还将讨论《精通MATLAB工具箱》中处理富里哀级数的函数集。

14.1 快速富里哀变换

在MATLAB中,函数fft计算一个信号的离散富里哀变换。在数据的长度是2的幂次或质因数的乘积的情况下,就用快速富里哀变换(FFT)来计算离散富里哀变换。当数据长度是2的幂次时,计算速度显著增加,因此,只要可能,选择数据

长度为2的幂次或者用零来填补数据,使得数据长度等于2的幂次显得非常重要。在《MATLAB参考指南》中可找到有关该问题的讨论。

MATLAB中实现的快速富里哀变换,是按照工科教材中常使用的方法。

F(k)=FFT{f(n)}

因为MATLAB不允许零下标,所以移动了一个下标值。

相应的逆变换为:

为了说明FFT的使用,考虑估计连续信号的富里哀变换的问题。

解析上,该富里哀变换为:

虽然在这种情况下,由于知道了富里哀变换的解析结果,再运用FFT没有多大的实用价值,但这个例子说明了对不常见的信号,特别是那些解析上难以找到富里哀变换的信号,一个估计富里哀变换的方法。下面的MATLAB语句用FFT估计F(w),并且用图形把所得到结果与上面的解析表达式的结果进行比较:

>>N=128; % choose a power of 2 for speed

>>t=linspace(0, 3, N); % time points for function evaluation

>>f=2*exp(-3*t); % evaluate the function and minimize aliasing:f(3)~0 >>Ts=t(2)-t(1); % the sampling period

>>Ws=2*pi/Ts; % the sampling frequency in rad/sec >>F=fft(f); % compute the fft >>Fp=F(1 : N/2+1)*Ts;

图14.1 富里哀变换两种结果的比较

仅从F中取正频率分量,并且乘以采样间隔计算F(w)。

>>W=Ws*(0 : N/2)/N

它建立了连续频率轴,该轴起始于0,终止于奈魁斯特(Nyquist)频率Ws/2,

>>Fa=2./(3+j*w); % evaluate analytical Fourier transform

>>plot(W, abs(Fa), W, abs(Fp), ‘ + ‘ ) % generate plot, ‘ + ‘ mark fft results

>>xlabel(‘ Frequency, Rad/s ‘),ylabel(‘ |F(w)| ‘)

MATLAB提供了大量的完成一般信号处理任务的函数。它们列于表14.1:

表14.1

信号处理函数 conv 卷积

conv2 2维卷积

fft 快速富里哀变换

fft2 2维快速富里哀变换 ifft 逆快速富里哀变换

ifft2 2维逆快速富里哀变换 filter 离散时间滤波器

filter2 2维离散时间滤波器 abs 幅值

angle 四个象限的相角

unwrap 在360°边界清除相角突变 fftshift 把FFT结果平移到负频率上 nextpow2 2的下一个较高幂次

14.2 富里哀级数

MATLAB本身没有特别关于富里哀级数分析和处理的函数。不过,通过创建M文件函数,可容易加上这些函数。在这一节,将介绍《精通MATLAB工具箱》中富里哀级数函数。在介绍之前,首先定义实周期信号f(t)的富里哀级数表示形式。 给出富里哀级数的复指数形式为:

式中的富里哀级数的系数是:

且基频为 。式中T0满足f(t+ T0)=f(t)。 给出富里哀级数的三角数形式为:

式中的富里哀级数的系数是:

且基频为 。式中T0满足f(t+ T0)=f(t)。

表14.2

精通MATLAB的富里哀级数函数

fsderiv(Kn,Wo) 富里哀级数的微分 fseval(Kn,t,Wo) 计算富里哀级数

因篇幅问题不能全部显示,请点此查看更多更全内容