实验四FFT算法的应用

2015-2016下期数字信号处理实验指导书1实验四FFT算法的应用一、实验目的：1、加深对离散信号的DFT的理解及其FFT算法的运用。2、学会利用matlab编译基2-FFT。二、实验原理：1、定义式N点序列的DFT和IDFT变换定义式如下：10[][]NknNnXkxnW，101[][]NknNkxnXkWN][kX为输入序列][nx对应的N个离散频率点的相对幅度。其中knNW为旋转因子，n为原序列里的某一点，k是DFT过后的序列里的某一点，N为变换的点数。利用旋转因子2jnkknNNWe具有周期性，可以得到快速算法（FFT）。在MATLAB中，可以用函数X=fft（x，N）和x=ifft（X，N）计算N点序列x的DFT正、反变换。2、DFT算法（1）DFT的运算量10)1,2,1,0()()]([)(NnknNNkWnxnxDFTkX在所有复指数值knNW的值全部已算好的情况下，要计算一个][kX需要N次复数乘法和N－1次复数加法。算出全部N点][kX共需2N次复数乘法和N（N-1）次复数加法。复数乘：2N次，复数加：2)1(NNN次。（2）减少DFT运算量的方法：a)将长度N变短。b)利用knNW的性质：周期性：)(rNnNnNWW共轭对称性：*)()(nNnNWW2015-2016下期数字信号处理实验指导书2可约性：nNrnrNWW3、DFT的快速算法（FFT，FastFourierTransformation）（1）FFT算法首先由Cooly-Tuky提出了“基-2FFT算法”。要求长度N满足MN2（M为整数），若不满足可将序列补零延长，使其满足长度要求。（2）FFT算法分为时间抽取（DIT）和频率抽取（DIF）两大类，其中时间抽取算法是把时间序列x(n)分解为两个长度为N/2点的序列，即偶数序列和奇数序列：)1(),5(),3(),1()12()()2(),4(),2()2()(21NxxxxrxrxNxxxrxrxr=1,2,…,N/2-1所以DFT运算也分成两部分：10)()(NnknNWnxkX12/0)12(12/02)12()2(NrkrNNrrkNWrxWrx12/02212/021)()(NrkNrkNNrrkNWWrxWrx【注意由x到x1】由可约性：nNrnrNWW可知rkNrkNWW2/2，得：12/02/212/02/1)()()(NrkNrkNNrrkNWWrxWrxkX对于某一个k值来言，kNW为常量提出来得：12/02/212/02/1)()()(NrrkNkNNrrkNWrxWWrxkX)()(21kXWkXkN公式（一）这时，)(1kX和)(2kX分别为)(1nx和)(2nx（n=0，1,2，…，N/2-1）的N/2点DFT，即：12/02/11)()(NnknNWnxkX12/02/22)()(NnknNWnxkX（k=0，1,2，…，N/2-1）2015-2016下期数字信号处理实验指导书3但是这时只是得到)(kX的前N/2个点的DFT。那么后N/2个点的DFT)2/(NkX和)(1kX、)(2kX：由)()()(21kXWkXkXkN推得：)2/()2/()2/(2)2/(1NkXWNkXNkXNkN12/0)2/(2/2)2/(12/0)2/(2/1)()(NnnNkNNkNNnnNkNWnxWWnx由周期性)(rNnNnNWW可知knNnNknNnNkN)2/(2/，以及kNjNkjNNkjNkNWeeeW/2/))2/(2(2/，所以：)()()2/(21kXWkXNkXkN公式（二）这样，N点的DFT就可以由两个2N个点的DFT来计算。三、实验内容：(1)128点实数序列nNnnNnNnx其它,012,...,2,1,0),192cos(21)72cos()(用一个128点DFT程序，一次算出)]([)(nxDFTkX，并绘出)(kX。(2)对（1）式，用两个64点的复数FFT程序，一次算出)]([)(nxDFTkX，并2015-2016下期数字信号处理实验指导书4绘出)(kX。设计流程：A、给出128点的原始序列，并绘图；B、将原始序列分成两个64点的序列，分别是偶数序列和奇数序列；C、分别使用fft函数求这两个64点序列的64点DFT；D、利用公式（一）公式（二）求得原始128点序列的DFT，并绘图，比较参考程序1的图是否一样；E、对D项中的DFT结果，使用ifft函数进行DFT反变换，并绘图。（3）已知某序列)(nx在单位圆上的N=64等分样点的Z变换为63,...,2,1,0,8.011)()(/2kekXzXNkjk用N点IFFT程序计算)]([)(_kXIDFTnx，绘出)(_nx。四、实验结果：内容提要实验结果内容（1）1参考程序1；N=64;k=0:2*N-1;x=cos(2*pi/N*7*k)+1/2*cos(2*pi/N*19*k);F_128=fft(x);stem(k/128,abs(F_128));title('序列x(k)图形');2参考程序1绘图结果；2015-2016下期数字信号处理实验指导书500.10.20.30.40.50.60.70.80.91010203040506070序列x(k)图形内容（2）1程序；N=64;x=zeros(1,128);u=zeros(1,128);y=zeros(1,64);y1=zeros(1,64);y2=zeros(1,64);fork=0:2*N-1;x(k+1)=cos(2*pi/N*7*k)+1/2*cos(2*pi/N*19*k);endfork=1:64y1(k)=x(2*k);y2(k)=x(2*k-1);endt1=fft(y1,N);t2=fft(y2,N);form=1:64y(m)=t2(m)+t1(m)*exp(-i*m*2*pi/128);t(m)=t2(m)-t1(m)*exp(-i*m*2*pi/128);endform=1:642015-2016下期数字信号处理实验指导书6u(m)=y(m);u(m+64)=t(m);endk=1:128subplot(3,1,1);stem(k,x);xlabel('k');ylabel('x');title('序列x(k)图形');subplot(3,1,2);stem(k/128,abs(u));title('序列x(k)FFT图形');v=ifft(u,128);subplot(3,1,3);stem(k,v);title('序列x(k)dft的反变换图形');2绘图结果；020406080100120140-2-1012kx序列x(k)图形00.10.20.30.40.50.60.70.80.91020406080序列x(k)FFT图形020406080100120140-2-1012序列x(k)dft的反变换图形内容（3）1参考程序2；m=0:63;N=64;m2=-j*2*pi*m/N;x=1./(1-0.8*exp(m2));subplot(2,1,1);stem(m,x);2015-2016下期数字信号处理实验指导书7xlabel('m');ylabel('Mangitude');title('X[m]的原图型');h=ifft(x,N);subplot(2,1,2);stem(m,h);xlabel('k');ylabel('Mangitude');title('反变换后x[k]的图形');2参考程序2绘图结果；0102030405060700246mMangitudeX[m]的原图型01020304050607000.511.5kMangitude反变换后x[k]的图形