信号与系统课程设计——FFT的计算机实现 快速傅立叶变换FFT的计算机实现.doc

信号与系统课程设计——FFT的计算机实现 快速傅里叶变换（FFT）的计算机实现 摘要：本文是信号与系统课程的课程设计，旨在熟悉FFT的计算过程，结合DFT物理意义和实验结果加深对傅立叶变换的理解。文章首先用MATLAB对一个简单信号进行FFT仿真，得出频谱图；其次完成了FFT的C语言实现，结合MATLAB作图及数据处理功能得出了C实现下的FFT结果；最后，讨论分析实验结果。 关键词：DFT、基--2按时间抽取FFT算法、MATLAB、C、频谱、物理意义 算法描述 DFT的运算量 减少运算的方法： 化长序列为短序列。如将长度为N的序列分解为两个长度为N/2的序列 利用的性质（注：本文中的C程序未用到此性质） C程序采用基--2按时间抽取的FFT算法 设输入序列长度为 (M为正整数)，将该序列按时间顺序的奇偶分解为越来越短的子序列，称为基2按时间抽取的FFT算法,也称为Coolkey-Tukey算法。若N不满足条件，则人为地加上若干零值，使。 实验设计与步骤 为简单起见,同时不失一般性,本实验采用三个余弦成分和一个支流偏置成分叠加所得的信号作为信号源。 MATLAB的FFT仿真和C的FFT实现 关于频谱的理论分析： MATLAB的FFT仿真 MATLAB程序代码 function [ output_args ] = FFT2( x0,x1,f1,w1,x2,f2,w2,x3,f3,w3,fs ) %这是一个自定义函数，输入被采样的信号的参数和采样频率，然后输出原信号波形、采样信号序列、采样序列幅度谱和相位谱。函数规定信号由一个直流成分（大小为x0），三个余弦成分（各自频率分别为f1,f2,f3,幅值为x1,x2,x3,初相位为w1,w2,w3）,采样频率为fs，采样时间为1秒。x0,x1,f1,w1,x2,f2,w2,x3,f3,w3,fs作为参数输入 t=0:1/fs:1-1/fs;%定义采样时刻 N=length(t);%采样序列长度 s=x0+x1*cos(2*pi*f1*t+pi/180*w1)+x2*cos(2*pi*f2*t+pi/180*w2)+x3*cos(2*pi*f3*t+pi/180*w3);%采样信号 y=fft(s);%快速傅立叶变换 tt=0:1/10000:1;%原信号描点 ss=x0+x1*cos(2*pi*f1*tt+pi/180*w1)+x2*cos(2*pi*f2*tt+pi/180*w2)+x3*cos(2*pi*f3*tt+pi/180*w3);%原信号 figure; plot(tt,ss); grid; title(原始信号);%原信号波形 figure; subplot(3,1,1); plot(0:N-1,s(1:N),-o); xlim([0 N-1]); grid; title(采样信号); subplot(3,1,2); plot(0:N-1,abs(y(1:N)),-o); xlim([0 N-1]); grid; xlabel(k); ylabel(幅度); title(理想采样信号的幅度谱); subplot(3,1,3); plot(0:N-1,angle(y(1:N)),-o); grid; xlabel(k); ylabel(相位); axis([0 N-1 -pi pi]); title(理想采样信号的相位谱); end; 改变采样频率 依次键入：%其中FFT2()是自定义函数,对信号在1秒内以频率fs进行采样 FFT2(1,1,1,90,2,2,180,3,3,180,32) FFT2(1,1,1,90,2,2,180,3,3,180,16) FFT2(1,1,1,90,2,2,180,3,3,180,8) FFT2(1,1,1,90,2,2,180,3,3,180,4) 原始信号如图： 图1  图2、fs =32Hz 图3、fs =16Hz 图4、fs =8Hz 图5、fs =4Hz  注意到： 1、幅度谱，除却k=0外，图形呈轴对称分布； 2、相位谱，N为偶数时除却k=0和N/2外（N为奇数时，除却k=0外），图形呈中心对称分布。 理论上,由于复指数的周期性，长度为N（假设N为偶数，N为奇数时类似）的的DFT频谱分析，可把k=N/2,N/2+1,…,N-1看作是负频率-N/2,-N/2+1,-N/2+2,…,-2,-1,特别的，对于实序列而言，正负频率成分的幅值必须相等，而初相位必须相反。 下面取fs=8Hz和4Hz部分计算结果分析，其他情况可类似处理。  可看到：fs=8Hz时，由频谱结合物理意义计算的得到的余弦成分的幅值和初相位与原信号一致；fs=4Hz时，由频谱结合物理意义计算的