数字信号处理实验报告-快速傅立叶变换

数字信号处理实验报告快速傅立叶变换实验一快速傅立叶变换(信息工程专业)一实验目的在理论学习的基础上，通过本实验加深对快速傅立叶变换的理解；熟悉并掌握按时间抽取FFT算法的程序；了解应用FFT进行信号频谱分析过程中可能出现的问题，例如混淆、泄漏、栅栏效应等，以便在实际中正确应用FFT。二实验内容仔细分析教材第六章‘时间抽取法FFT’的算法结构，编制出相应的用FFT进行信号分析的C语言（或MATLAB语言）程序；用FFT程序计算有限长度正弦信号信号频率f＝50Hz，采样点数N=32,采样间隔T=0.000625s信号频率f＝50Hz，采样点数N=32,采样间隔T=0.005s信号频率f＝50Hz，采样点数N=32,采样间隔T=0.0046875s信号频率f＝50Hz，采样点数N=32,采样间隔T=0.004s信号频率f＝50Hz，采样点数N=64,采样间隔T=0.000625s信号频率f＝250Hz，采样点数N=32,采样间隔T=0.005sg)将c）信号后补32个0，做64点FFT三实验步骤记录下实验内容中各种情况下的X(k)值，做出频谱图并深入讨论结果，说明参数的变化对信号频谱产生哪些影响。频谱只做模特性，模的最大值＝1，全部归一化；打印出用C语言（或MATLAB语言）编写的FFT源程序，并且在每一小段处加上详细的注释说明；用C语言（或MATLAB语言）编写FFT程序时，要求采用人机界面形式：N,T,f变量均由键盘输入，补零或不补零要求设置一开关来选择。四实验结果1、人机界面形式的FFT程序运行示意2、各个情况下的X(k)值记录与对应的频谱图a信号频率f＝50Hz，采样点数N=32,采样间隔T=0.000625seq\o\ac(○,1)谱线特征信号频率f＝50Hz，采样间隔T=0.000625s，采样率1/T=1600Hz，故图中k=1谱线为信号频率对应的谱线，k=31为频域上下一个周期负频率部分eq\o\ac(○,2)频谱混叠1600Hz>2*50Hz，无频谱混叠eq\o\ac(○,3)栅栏效应无栅栏效应eq\o\ac(○,4)频率泄露信号频率被频域采样间隔整除，无法观察到频率泄露xk=Columns1through100.00001.00000.00000.00000.00000.00000.00000.00000.00000.0000Columns11through200.00000.00000.00000.00000.00000.00000.00000.00000.00000.0000Columns21through300.00000.00000.00000.00000.00000.00000.00000.00000.00000.0000Columns31through320.00001.0000b信号频率f＝50Hz，采样点数N=32,采样间隔T=0.005seq\o\ac(○,1)谱线特征信号频率f＝50Hz，采样间隔T=0.005s，采样率1/T=200Hz，故图中k=8谱线为信号频率对应的谱线，k=24为频域上下一个周期负频率部分eq\o\ac(○,2)频谱混叠200Hz>2*50Hz，无频谱混叠eq\o\ac(○,3)栅栏效应无栅栏效应eq\o\ac(○,4)频率泄露信号频率被频域采样间隔整除，无法观察到频率泄露xk=Columns1through100.00000.00000.00000.00000.00000.00000.00000.00001.00000.0000Columns11through200.00000.00000.00000.00000.00000.00000.00000.00000.00000.0000Columns21through300.00000.00000.00000.00001.00000.00000.00000.00000.00000.0000Columns31through320.00000.0000c信号频率f＝50Hz，采样点数N=32,采样间隔T=0.0046875seq\o\ac(○,1)谱线特征信号频率f＝50Hz，采样间隔T=0.0046875s，采样率1/T=213.33Hz，图中无法显示出k=7.5与k=24.5处的谱线，因而产生了栅栏效应。eq\o\ac(○,2)频谱混叠213.33Hz>2*50Hz，无频谱混叠eq\o\ac(○,3)栅栏效应图中无法显示出k=7.5与k=24.5处的谱线，因而产生了栅栏效应。eq\o\ac(○,4)频率泄露栅栏效应所呈现出的周期sinc型谱线就是频率泄露产生的谱，可以明显的观察出频率泄露。xk=Columns1through100.10760.11000.11750.13240.15940.21220.34101.00000.99040.3312Columns11through200.20190.14850.12060.10440.09500.09000.08840.09000.09500.1044Columns21through300.12060.14850.20190.33120.99041.00000.34100.21220.15940.1324Columns31through320.11750.1100d信号频率f＝50Hz，采样点数N=32,采样间隔T=0.004seq\o\ac(○,1)谱线特征信号频率f＝50Hz，采样间隔T=0.004s，采样率1/T=250Hz，图中无法显示出k=6.4与k=25.6处的谱线，因而产生了栅栏效应。eq\o\ac(○,2)频谱混叠250Hz>2*50Hz，无频谱混叠eq\o\ac(○,3)栅栏效应图中无法显示出k=6.4与k=25.6处的谱线，因而产生了栅栏效应。eq\o\ac(○,4)频率泄露栅栏效应所呈现出的周期sinc型谱线就是频率泄露产生的谱，可以明显的观察出频率泄露。相比于上一种情况，因为sinc主瓣相差的频率变小，所以频率泄露情况有所减弱。xk=Columns1through100.07950.08280.09360.11600.16200.27961.00000.68890.26920.1739Columns11through200.13260.11030.09680.08840.08320.08040.07950.08040.08320.0884Columns21through300.09680.11030.13260.17390.26920.68891.00000.27960.16200.1160Columns31through320.09360.0828e信号频率f＝50Hz，采样点数N=64,采样间隔T=0.000625seq\o\ac(○,1)谱线特征信号频率f＝50Hz，采样间隔T=0.000625s，采样率1/T=1600Hz，故图中k=2谱线为信号频率对应的谱线，k=62为频域上下一个周期负频率部分eq\o\ac(○,2)频谱混叠1600Hz>2*50Hz，无频谱混叠eq\o\ac(○,3)栅栏效应无栅栏效应eq\o\ac(○,4)频率泄露信号频率被频域采样间隔整除，无法观察到频率泄露xk=Columns1through100.00000.00001.00000.00000.00000.00000.00000.00000.00000.0000Columns11through200.00000.00000.00000.00000.00000.00000.00000.00000.00000.0000Columns21through300.00000.00000.00000.00000.00000.00000.00000.00000.00000.0000Columns31through400.00000.00000.00000.00000.00000.00000.00000.00000.00000.0000Columns41through500.00000.00000.00000.00000.00000.00000.00000.00000.00000.0000Columns51through600.00000.00000.00000.00000.00000.00000.00000.00000.00000.0000Columns61through640.00000.00001.00000.0000f信号频率f＝250Hz，采样点数N=32,采样间隔T=0.005seq\o\ac(○,1)谱线特征信号频率f＝250Hz，采样间隔T=0.005s，采样率1/T=200Hz，故图中k=8谱线为信号频率对应的谱线，k=24为频域上下一个周期负频率部分eq\o\ac(○,2)频谱混叠200Hz<2*250Hz,所以发生了频谱混叠。此时250Hz的信号频率由于采样率200Hz的作用，被搬移到了50Hz，-50Hz的信号频率由于采样率200Hz的作用，被搬移到了150Hz，从而观察到了与第二种情况完全相同的频谱。eq\o\ac(○,3)栅栏效应无栅栏效应eq\o\ac(○,4)频率泄露信号频率被频域采样间隔整除，无法观察到频率泄露xk=Columns1through100.00000.00000.00000.00000.00000.00000.00000.00001.00000.0000Columns11through200.00000.00000.00000.00000.00000.00000.00000.00000.00000.0000Columns21through300.00000.00000.00000.00001.00000.00000.00000.00000.00000.0000Columns31through320.00000.0000g将c）信号后补32个0，做64点FFTeq\o\ac(○,1)谱线特征信号频率f＝50Hz，采样间隔T=0.0046875s，采样率1/T=213.33Hz，故图中k=15谱线为信号频率对应的谱线，k=49为频域上下一个周期负频率部分.虽然k=15处的谱线采到了sinc主瓣的峰值，但由于时域补零引起了频域分辨率提高，使得原本第三种情况的频谱中为零的那些频率中间的旁瓣分量被展示出来，从而观察到了频率泄露。eq\o\ac(○,2)频谱混叠213.33Hz>2*50Hz，无频谱混叠eq\o\ac(○,3)栅栏效应无栅栏效应eq\o\ac(○,4)频率泄露虽然k=15处的谱线采到了sinc主瓣的峰值，但由于时域补零引起了频域分辨率提高，使得原本第三种情况的频谱中为零的那些频率中间的旁瓣分量被展示出来，从而观察到了频率泄露。xk=Columns1through100.06900.00000.07050.00000.07530.00000.08480.00000.10210.0000Columns11through200.13590.00000.21850.00000.64071.00000.63460.00000.21220.0000Columns21through300.12940.00000.09520.00000.07730.00000.06690.00000.06090.0000Columns31through400.05770.00000.05660.00000.05770.00000.06090.00000.06690.0000Columns41through500.07730.00000.09520.00000.12940.00000.21220.00000.63461.0000Columns51through600.64070.00000.21850.00000.13590.00000.10210.00000.08480.0000Columns61through640.07530.00000.07050.0000五实验小结1、实验思考对于时域上的正弦信号，当采样频率，信号频率，采样点数符合一定规律时，可以明显的削弱混淆、泄露现象。具体猜想为频率采样采到单频信号的主频时，即K=N*f/fs,K=N*(1-f/fs)，其中K为整数时，混淆与泄露现象不明显2、参数变化对信号频谱产生大的影响eq\o\ac(○,1)采样率（采样间隔），信号频率：采样率决定了频域的截止频率。采样率越低，信号频率越高，越容易产生频率混叠。eq\o\ac(○,2)时域采样点数（时域窗的大小）：时域采样点数决定了时域加矩形窗后频域卷积的sinc的形式，即影响频率泄露的情况。eq\o\ac(○,3)频域采样点数：决定了频域分辨率，即影响栅栏效应的出现与否。3、关于频谱混叠采样序列的频谱是采样信号频谱的周期延拓。当采样频率不满足奈奎斯特定律，即采样频率小于两倍的信号频率时，经过采样就会发生频谱混叠。这会导致采样后的信号序列频谱不能真实的反映原信号的频谱，而且会导致得到的混叠频谱无法经过内插准确恢复到原信号。4、关于栅栏效应DFT是单位圆上均匀采样，它不可能将频谱视为一个连续函数。这样就产生了栅栏效应，用DFT来观察频谱只能在离散的采样点上看到真实的频谱。这样的话，就会有一些频谱的峰点或谷点被‘栅栏’挡住，不能被我们观察到。减小栅栏效应的一个方法是在原信号序列的末端补一些零值，这种方法的实质是改变了对真实频谱采样的点数和位置，相当于搬动了栅栏的位置，从而使得原来被挡住的一些频谱的峰点或谷点显露出来。5、关于频谱泄露实际中的信号序列往往很长，甚至是无限长序列。为了方便，我们经常用截短的序列（即给原信号序列乘以一个矩形窗函数）来近似它们。而矩形窗函数的频谱不是有限带宽的，从而它和原信号的频谱卷积以后会扩展原信号的频谱。泄露是不能和混淆完全分开的。因为泄露会导致频谱的扩展，从而造成混淆。为了减小泄露的影响，可以选择适当的窗函数使频谱的扩散减小到最小。例如采用升余弦窗，二阶升余弦窗，凯塞窗等。六实验程序clcclearf=input('Inputfrequencyofthesignal:f\n');N=input('Inputnumberofpoints:N\n');T=input('Inputsamplingtime:T\n');flag=input('Addzerotoosamplingsignalornot?yes=1no=0\n');if(flag)ZeroNum=input('Inputnmberofzeros\n');elseZeroNum=0;endn=[0:1:N-1];x=sin(2*pi*f*T*n);%采样ifZeroNum>0x=[xzeros(1,ZeroNum)];endN=length(x);%求出序列长度NV2=N/2;%倒置NM1=N-1;I=0;J=0;whileI<NM1ifI<J%若I<J则交换自然顺序x(I)和倒序x(J)的存储单元T=x(J+1);x(J+1)=x(I+1);x(I+1)=T;endK=NV2;whileK<=JJ=J-K;K=K/2;endJ=J+K;I=I+1;endm=log2(N);