实验2_fft的实现

《实验2_fft的实现》由会员分享，可在线阅读，更多相关《实验2_fft的实现（13页珍藏版）》请在金锄头文库上搜索。

1、电 子 科 技 大 学实 验 报 告学生姓名： 学 号： 指导教师： 一、实验室名称：数字信号处理实验室二、实验项目名称：FFT 的实现三、实验原理：一FFT 算法思想：1DFT 的定义：对于有限长离散数字信号xn，0 n N-1，其离散谱xk可以由离散付氏变换（DFT）求得。DFT 的定义为：，k=0,1,N-1210NjkNnXkxe通常令 ，称为旋转因子。2jNeW2直接计算 DFT 的问题及 FFT 的基本思想：由 DFT 的定义可以看出，在 xn为复数序列的情况下，完全直接运算 N 点DFT 需要（N-1） 2次复数乘法和 N（N-1）次加法。因此，对于一些相当大的 N值（如 102

2、4）来说，直接计算它的 DFT 所作的计算量是很大的。FFT 的基本思想在于，将原有的 N 点序列分成两个较短的序列，这些序列的 DFT 可以很简单的组合起来得到原序列的 DFT。例如，若 N 为偶数，将原有的 N 点序列分成两个（N/2）点序列，那么计算 N 点 DFT 将只需要约(N/2) 2 2=N2/2 次复数乘法。即比直接计算少作一半乘法。因子（N/2） 2表示直接计算（N/2）点 DFT 所需要的乘法次数，而乘数 2 代表必须完成两个 DFT。上述处理方法可以反复使用，即（N/2）点的 DFT 计算也可以化成两个（N/4）点的DFT（假定 N/2 为偶数） ，从而又少作一半的乘法。

3、这样一级一级的划分下去一直到最后就划分成两点的 FFT 运算的情况。3基 2 按时间抽取（DIT）的 FFT 算法思想：设序列长度为 ，L 为整数（如果序列长度不满足此条件，通过在后2N面补零让其满足） 。将长度为 的序列 ，先按 n 的奇偶分成两组：2LN(0,1.)xnN，r=0,1,N/2-12rDFT 化为： 1/21/21(21)000/21/212200/21/21/2/200NNNnkrk rkNnrNrk kNNr rrkkrkr rXkDFTxWxxWxx 上式中利用了旋转因子的可约性，即： 。又令/2rkrkNW，则上式可以写成：/21/211/2/200,NNrk rkN

4、r rXkxWXx(k=0,1,N/2-1)12kNX可以看出， 分别为从 中取出的 N/2 点偶数点和奇数点序列12,k的 N/2 点 DFT 值，所以，一个 N 点序列的 DFT 可以用两个 N/2 点序列的 DFT 组合而成。但是，从上式可以看出，这样的组合仅表示出了 前 N/2 点的 DFTXk值，还需要继续利用 表示 的后半段本算法推导才完整。利用旋12,XkXk转因子的周期性，有： ，则后半段的 DFT 值表达式：(/2)/rrNNW，同样， /21/21()21/ /2100Nrkrkr rXkxx22NXk（k=0,1,N/2-1） ，所以后半段（k=N/2,N-1）的 DFT

5、 值可以用前半段 k 值表达式获得，中间还利用到 ，得到后半段的 值表达()22NkkkW式为： （k=0,1,N/2-1） 。1kXkWX这样，通过计算两个 N/2 点序列 的 N/2 点 DFT ，可以12,xn12,Xk组合得到 N 点序列的 DFT 值 ，其组合过程如下图所示：k 1Xk 12kNXW-1 2XknkNW12kNX比如，一个 N = 8 点的 FFT 运算按照这种方法来计算 FFT 可以用下面的流程图来表示：W0W0W2W0W2W0W1W2W3x(0)x(1)x(2)x(3)x(4)x(5)x(6)x(7) X(7)X(6)X(5)X(4)X(3)X(2)X(1)X(0

6、)W0W0W04基 2 按频率抽取（DIF）的 FFT 算法思想：设序列长度为 ，L 为整数（如果序列长度不满足此条件，通过在后2N面补零让其满足） 。在把 按 k 的奇偶分组之前，把输入按 n 的顺序分成前后两半：X1/21100/2/21/21()200/2120,0,1.NNNnknknkNnN nknkNNn knnDFTxxWxxWx g因为 ，则有 ，所以：2NW2(1)Nkk/210(),01,.2NknkNnXkxxWg按 k 的奇偶来讨论，k 为偶数时： /2120,.NrnNnrxkk 为奇数时：/21(21)0,0,.1rnnXWN g前面已经推导过 ，所以上面的两个等式

7、可以写为：2/2rkrkNW/21 /20,01,./2rnNnrxg/21 /20,./1NnrnXWN 通过上面的推导， 的偶数点值 和奇数点值 分别可以XkX2Xr由组合而成的 N/2 点的序列来求得，其中偶数点值 为输入 xn的前半段和后半段之和序列的 N/2 点 DFT 值，奇数点值 为输入 xn的前半段和21r后半段之差再与 相乘序列的 N/2 点 DFT 值。nNW令 ， ，则有：12xx22nNxnxWg/2 /11/22/20 0,0,1.2N Nrn rnNn nXrX g这样，也可以用两个 N/2 点 DFT 来组合成一个 N 点 DFT，组合过程如下图所示：x 2xn-

8、1 2NxnnNW2nNxW二在 FFT 计算中使用到的 MATLAB 命令：函数 fft(x)可以计算 R 点序列的 R 点 DFT 值；而 fft(x,N)则计算 R 点序列的 N 点 DFT，若 RN，则直接截取 R 点 DFT 的前 N 点，若 RN，则 x 先进行补零扩展为 N 点序列再求 N 点 DFT。函数 ifft(X)可以计算 R 点的谱序列的 R 点IDFT 值；而 ifft(X,N)同 fft(x,N)的情况。四、实验目的：离散傅氏变换（DFT）的目的是把信号由时域变换到频域，从而可以在频域分析处理信息，得到的结果再由逆 DFT 变换到时域。FFT 是 DFT 的一种快速

9、算法。在数字信号处理系统中，FFT 作为一个非常重要的工具经常使用，甚至成为 DSP 运算能力的一个考核因素。本实验通过直接计算 DFT，利用 FFT 算法思想计算 DFT，以及使用 MATLAB函数中的 FFT 命令计算离散时间信号的频谱，以加深对离散信号的 DFT 变换及FFT 算法的理解。五、实验内容：a) 计算实数序列 的 256 点 DFT。5()cos,0261xnnb) 计算周期为 1kHz 的方波序列（占空比为 50，幅度取为/-512，采样频率为 25kHz，取 256 点长度） 256 点 DFT。六、实验器材（设备、元器件）：安装 MATLAB 软件的 PC 机一台，DS

10、P 实验演示系统一套。七、实验步骤：（1） 先利用 DFT 定义式，编程直接计算 2 个要求序列的 DFT 值。（2） 利用 MATLAB 中提供的 FFT 函数，计算 2 个要求序列的 DFT 值。（3） (拓展要求)不改变序列的点数，仅改变 DFT 计算点数（如变为计算 1024 点 DFT 值） ，观察画出来的频谱与前面频谱的差别，并解释这种差别。通过这一步骤的分析，理解频谱分辨力的概念，解释如何提高频谱分辨力。（4） 利用 FFT 的基本思想（基 2DIT 或基 2DIF） ，自己编写 FFT计算函数，并用该函数计算要求序列的 DFT 值。并对前面 3 个结果进行对比。（5） (拓展要

11、求)尝试对其他快速傅立叶变换算法（如 Goertzel 算法）进行 MATLAB 编程实现，并用它来计算要求的序列的 DFT 值。并与前面的结果进行对比。（6） （拓展要求）在提供的 DSP 实验板上演示要求的 2 种序列的 FFT算法（基 2DIT） ，用示波器观察实际计算出来的频谱结果，并与理论结果对比。八、实验数据及结果分析：注：本次实验在寝室电脑上完成，所用 MATLAB 版本为 MATLAB R2010b程序：（1） 对要求的 2 种序列直接进行 DFT 计算的程序clc,clearclose all%正弦序列Fs = 1; %设采样频率为 1HzTs = 1/Fs; %采样周期N

12、= 256; %采样点数n = 0:N-1;x = cos(5*pi/16*n);A = myDFT(x,256); %直接法计算DFT%周期方波序列N2 = 256;Fs2 = 25*1e3;Ts2 = 1/Fs2;f2 = linspace(0,Fs2,N2);n2 = 1:N2;x2 = 512*square(2*pi*1000/Fs2*n2, 50);A2 = myDFT(x2,256);%-作图-figure,subplot(211), plot(n,x)xlabel(Time index n), ylabel(Amplitude)title(Sinusoidal, time-dom

13、ain sequence)subplot(212), stem(n,abs(A)xlabel(Frequency index k), ylabel(Magnitude)title(Magnitude of DFT samples(direct method)figure, subplot(211), plot(n2,x2)xlabel(Time index n), ylabel(Amplitude)title(Periodic square wave, time-domain sequence)subplot(212), stem(n2,abs(A2)xlabel(Frequency inde

14、x k), ylabel(Magnitude)title(Magnitude of DFT samples(direct method)函数文件：myDFT.mfunction X = myDFT(x,N)%利用定义式计算序列的DFT%输入参数：% x-离散时间信号% N-计算的DFT点数%输出参数：% X-N点DFT序列WN = exp( -1i*(2*pi/N) ); %旋转因子for t1 = 1:Nk = t1-1;X(t1) = 0;for t2 = 1:Nn = t2-1;X(t1) = X(t1) + x(t2)*WN( n*k );endend（2） 对要求的 2 种序列进行基

15、 2DIT 和基 2DIF FFT 算法程序序列生成的代码同（1） 。这里给出 FFT 算法的程序：基 2-DIT FFT 算法程序：函数文件：myfitfft.mfunction X = myditfft(x)%按时间抽选的基2-FFT算法%输入参数：% x-离散时间信号% %输出参数：% X-序列x的N点DFT （N是序列长度，必须是2 的整数次幂）M = nextpow2(length(x); % x的长度对应的2的最低次幂N = 2M; % 序列长度if length(x) Nx = x,zeros(1,N-length(x); % 若x的长度不是2的整数次幂，则补零直到长度为Nendnxd = b