实验三 用FFT进行谱分析
一、实验目的
1.进一步加深对DFT算法原理和基本性质的理解（因为FFT只是DFT的一种快速算法，所以FFT的运算结果必然满足DFT的基本性质）。
2.熟悉FFT算法原理和FFT子程序的应用。
3.学习用FFT对连续信号和时域离散信号进行谱分析的方法，了解可能出现的分析误差及其原因，以便在实际中正确应用FFT。
二、实验原理
1．快速傅立叶变换(FFT)算法
长度为N的序列x(n)的离散傅立叶变换X(k)为：
nkX(k)??x(n)WN，k?0，??，N?1
n?0N?1
N点的DFT可以分解为两个N/2点的DFT，每个N/2点的DFT又可以分解为两个N/4点的DFT。依此类推，当N为2的整数次幂时(N?2)，由于每分解一次降低一阶幂次，所以通过M次的分解，最后全部成为一系列2点DFT运算。以上就是按时间抽取的快速傅立
叶变换(FFT)算法。当需要进行变换的序列的长度不是2的整数次方的时候，为了使用以2为基的FFT，可以用末尾补零的方法，使其长度延长至2的整数次方。
序列X(k)的离散傅立叶反变换为
M
1N?1?nkx(n)??X(k)WN，n?0，??，N?1
Nk?0?1离散傅立叶反变换与正变换的区别在于WN变为WN，并多了一个1/N的运算。因为?1对于推导按时间抽取的快速傅立叶变换算法并无实质性区别，因此可将FFT和快WN和WN速傅立叶反变换（IFFT）算法合并在同一个程序中。 2．利用FFT进行频谱分析
若信号本身是有限长的序列，计算序列的频谱就是直接对序列进行FFT运算求得X(k)，
X(k)就代表了序列在[0,2?]之间的频谱值。
幅度谱 相位谱
X(k)?2XR(k)?XI2(k)
?(k)?arctanXI(k) XR(k)若信号是模拟信号，用FFT进行谱分析时，首先必须对信号进行采样，使之变成离散信号，然后就可按照前面的方法用FFT来对连续信号进行谱分析。按采样定理，采样频率fs应
6
大于2倍信号的最高频率，为了满足采样定理，一般在采样之前要设置一个抗混叠低通滤波器。用FFT对模拟信号进行谱分析的方框图如下所示。
三、主要实验仪器及材料
微型计算机、Matlab6.5教学版、TC编程环境。
四、实验内容
1．知识准备
实验前学生应认真复习DFT和FFT有关的知识，掌握快速傅里叶变换的基本原理以及如何用FFT等计算信号频谱。
2．离散时间信号（序列）的产生
利用MATLAB或C语言编程产生和绘制下列两有限长序列：
cos(?sin(nT)nT)、、sin(0.25?nT)?2?cos(0.5?nT) 8nTcos(0.125?nT)?2cos(0.25?nT)、sin(0.25?nT)
3．用一种语言编写FFT的通用程序块
4．画主程序实现框图并编写主程序，实现信号的谱分析。 5．记录下实验内容中各信号x(n)的X(k)值，作出频谱图。
五、思考题
1. 根据实验中各x(n)的X(k)值以及频谱图，说明参数的变化对信号频谱产生哪些影响？
2.如果周期信号的周期预先不知道，如何用FFT进行分析？
六、实验报告要求
1.简述实验原理及目的。
2.结合实验中所得给定典型序列幅频特性曲线，与理论结合比较，并分析说明误差产生的原因以及用FFT作谱分析时有关参数的选择方法。
3.总结实验所得主要结论。 4.简要回答思考题。 FFT参考子程序
void fft(double *dr,double *di,int N) {
int lh,m; lh=N/2;
m=int(log(N)/log(2)+0.9999); int j=lh; int k;
for(int i=1;i<(N-1);i++)
7
{
if(i double tr,ti; tr=dr[i];ti=di[i]; dr[i]=dr[j];di[i]=di[j]; dr[j]=tr;di[j]=ti; } k=lh; while(j>=k) {j=j-k;k=k/2;} j=j+k; } for(i=1;i<=m;i++) { int b=int(pow(2.0,(i-1))); for(j=0;j double p; p=j*pow(2.0,(m-i))*2.0*PI/N; for(k=j;k double tr,ti; tr=dr[k+b]*cos(p)+di[k+b]*sin(p); ti=di[k+b]*cos(p)-dr[k+b]*sin(p); dr[k+b]=dr[k]-tr; di[k+b]=di[k]-ti; dr[k]=dr[k]+tr; di[k]=di[k]+ti; k=int(k+pow(2.0,i)); } } } 8 } 实验四 利用FFT实现快速卷积 一、实验目的 1．加深理解FFT在实现数字滤波（或快速卷积）中的重要作用，更好的利用FFT进行数字信号处理。 2．掌握循环卷积和线性卷积两者之间的关系。 二、实验原理 数字滤波器根据系统的单位脉冲响应h(n)是有限长还是无限长可分为有限长单位脉冲响应（Finite Impulse Response）系统（简记为FIR系统）和无限长单位脉冲响应（Infinite Impulse Response）系统（简记为IIR系统）。 对于FIR滤波器来说，除了可以通过数字网络来实现外，也可以通过FFT的变换来实现。 首先我们知道，一个信号序列x(n)通过FIR滤波器时，其输出应该是x(n)与h(n)的卷积： y(n)?x(n)?h(n)?m????x(m)h(n?m) ?当h(n)是一个有限长序列，即h(n)是FIR滤波器，且0?n?N?1时 y(n)??h(m)x(n?m) m?0N?1在数字网络类的FIR滤波器中，普遍使用的横截型结构就是按这个卷积公式构成的。 应用FFT实现数字滤波器实际上就是用FFT来快速计算有限长度序列的线性卷积。 这种方法就是先将输入信号x(n)通过FFT变换为它的频谱采样值X(k)，然后再和FIR滤波器的频响采样值H(k)相乘，H(k)可事先存放在存储器中，最后再将乘积H(k)X(k)通过快速傅里叶变换（简称IFFT）还原为时域序列，即得到输出y(n)。 现以FFT求有限长序列间的卷积及求有限长度序列与较长序列间的卷积为例来讨论FFT的快速卷积方法。 1．序列x(n)和h(n)的长差不多。设x(n)的长为N1，h(n)的长为N2，要求 y(n)?x(n)?y(n)??h(m)x(n?m) m?0N?1用FFT完成这一卷积的具体步骤如下： ①为使两有限长序列的线性卷积可用其循环卷积代替而不发生混叠，必须选择循环卷积长度N?N1?N2?1，若采用基2-FFT完成卷积运算，要求N?2（m为整数）。 ②用补零方法使x(n)和h(n)变成列长为N的序列。 m 9 ?x(n)0?n?N1?1x(n)??N1?n?N?1?0?h(n)0?n?N2?1h(n)??N2?n?N?1?0 ③用FFT计算x(n)和h(n)的N点离散傅里叶变换 ④完成X(k)和H(k)乘积，Y(k)?x(k)H(k) ⑤用FFT计算Y(k)的离散傅里叶反变换得 ?1??nk?N?1?1??nk?y(n)???Y(k)?WN????Y(k)?WN???k?0?N?k?0?N?N?1? 2．当x(n)长度很长时，即，通常不允许等x(n)全部采集齐后再进行卷积，否则使输出相 对于输入有较长的延时，另外，若21NN21NN太大，h(n)要补上太多的零点，很不经济，且FFT的计算时间也要很长。为此，采用分段卷积的方法，即把x(n)分成长度与h(n)相仿的一段段，分别求出每段卷积的结果，然后用相应的方式把它们结合起来，便是总的输出。分段卷积方法主要有两种，即重叠相加法和重叠保留法。具体内容请参考数字信号处理教材中“快速离散傅里叶变换”一章中的线性卷积的FFT算法部分，本实验这部分不作重点要求。 三、主要实验仪器及材料 微型计算机、Matlab6.5教学版、TC编程环境。 四、实验内容 nh(n)?(1/2)RN2(n)，N取8。 1．数字滤波器的脉冲响应为 2 输入序列x(n)可选下列几种情况 ①x(n)?RN1(n) ，N1取8。 ②) x(n)?cos2?RN(n)，N1取8。 N11n?1?③) x(n)???RN1(n)，N1取8。 ?3?2．实验前，预先计算好x(n)?h(n)的值。 3．实验前，预先编制一个应用FFT实现数字滤波器的通用程序。 4．上机独立调试，并打印或记录实验结果。 5．将实验结果与预先笔算的结果比较，验证其正确性。 五、思考题 1．N同时取8、16、32时，线性卷积和循环卷积的结果有何不同，为什么？ 10