快速傅立叶变换及其应用.ppt
关于快速傅立叶变换及其应用现在学习的是第1页,共28页一、实验目的一、实验目的 了解计算了解计算DFTDFT算法存在的问题及改进途径。算法存在的问题及改进途径。掌握几种掌握几种DFTDFT算法(时间抽取算法算法(时间抽取算法DITDIT算法,频率抽取算算法,频率抽取算法法DIFDIF算法,线性调频算法,线性调频Z Z变换即变换即CZTCZT法)。法)。学习并掌握学习并掌握FFTFFT的应用。的应用。现在学习的是第2页,共28页二、实验原理二、实验原理有限长序列通过离散傅里叶变换有限长序列通过离散傅里叶变换(DFT)DFT)将其频域离散化将其频域离散化成有限长序列成有限长序列.但其计算量太大但其计算量太大(与与N N的平方成正比)的平方成正比),很很难实时地处理问题难实时地处理问题,因此引出了快速傅里叶变换因此引出了快速傅里叶变换(FFT)FFT)。FFTFFT并不是一种新的变换形式并不是一种新的变换形式,它只是它只是DFTDFT的一种快速算的一种快速算法法.并且根据对序列分解与选取方法的不同而产生了并且根据对序列分解与选取方法的不同而产生了FFTFFT的多种算法的多种算法.现在学习的是第3页,共28页DFTDFT的快速算法的快速算法FFTFFT是数字信号处理的基本方法和基本技是数字信号处理的基本方法和基本技术,是必须牢牢掌握的。术,是必须牢牢掌握的。时间抽选时间抽选FFTFFT算法的理论推导和流图详见数字信号处理教材。算法的理论推导和流图详见数字信号处理教材。该算法遵循两条准则:该算法遵循两条准则:(1 1)对时间奇偶分;()对时间奇偶分;(2 2)对频率前后分。)对频率前后分。这种算法的流图特点是:这种算法的流图特点是:(1 1)基本运算单元都是蝶形)基本运算单元都是蝶形 任何一个长度为任何一个长度为N=2N=2M M的序列,总可通过的序列,总可通过M M次分解最后成为次分解最后成为2 2点的点的DFTDFT计算。如图所示:计算。如图所示:现在学习的是第4页,共28页W WN Nk k称为旋转因子称为旋转因子 计算方程如下:计算方程如下:X Xm+1m+1(p)=X(p)=Xm m(p)+W(p)+WN Nk kX Xm m(q)(q)X Xm+1m+1(q)=X(q)=Xm m(p)-W(p)-WN Nk kX Xm m(q)(q)现在学习的是第5页,共28页 (2 2)同址(原位)计算)同址(原位)计算 这是由蝶形运算带来的好处,每一级蝶形运算的结果这是由蝶形运算带来的好处,每一级蝶形运算的结果X Xm+1m+1(p)(p)无须另外存储,只要再存入无须另外存储,只要再存入X Xm m(p)(p)中即可,中即可,X Xm+1m+1(q)(q)亦然。这样将大大节省存储单元。亦然。这样将大大节省存储单元。(3 3)变址计算)变址计算 输入为输入为“混序混序”(码位倒置)排列,输出按自然序排(码位倒置)排列,输出按自然序排列,因而对输入要进行列,因而对输入要进行“变址变址”计算(即码位倒置计算)。计算(即码位倒置计算)。“变址变址”实际上是一种实际上是一种“整序整序”的行为,目的是保证的行为,目的是保证“同址同址”。现在学习的是第6页,共28页FFTFFT的应用的应用凡是利用付里叶变换来进行分析、综合、变换的地凡是利用付里叶变换来进行分析、综合、变换的地方,都可以利用方,都可以利用FFTFFT算法来减少其计算量。算法来减少其计算量。FFTFFT主要应用在主要应用在 1 1、快速卷积、快速卷积 2 2、快速相关、快速相关 3 3、频谱分析、频谱分析现在学习的是第7页,共28页快速傅立叶变换的快速傅立叶变换的MATLABMATLAB实现实现提供提供fftfft函数计算函数计算DFTDFT格式格式 X=fft(x)X=fft(x)X=fft(x,N)X=fft(x,N)如果如果x x的长度小于的长度小于N N,则在其后填零使其成为,则在其后填零使其成为N N点序列,若省点序列,若省略变量略变量N N,则,则DFTDFT的长度即为的长度即为x x的长度。的长度。如果如果N N为为2 2的幂,则得到高速的基的幂,则得到高速的基-2FFT-2FFT算法;若算法;若N N不是不是2 2的乘的乘方,则为较慢的混合算法。方,则为较慢的混合算法。如果如果x x是矩阵,则是矩阵,则X X是对矩阵的每一列向量作是对矩阵的每一列向量作FFTFFT。现在学习的是第8页,共28页由题目可得由题目可得x=0.5*sin(2*pi*15*t)+2*sin(2*pi*40*t)fs=100N=128/1024例:已知信号由例:已知信号由15Hz15Hz幅值幅值0.50.5的正弦信号和的正弦信号和40Hz40Hz幅值幅值2 2的正弦信号组成,数据采样频率为的正弦信号组成,数据采样频率为100Hz,100Hz,试绘制试绘制N=128N=128点点DFTDFT的幅频图。的幅频图。现在学习的是第9页,共28页fs=100;N=128;n=0:N-1;t=n/fs;x=0.5*sin(2*pi*15*t)+2*sin(2*pi*40*t);y=fft(x,N);f=(0:length(y)-1)*fs/length(y);mag=abs(y);stem(f,mag);title(N=128点)现在学习的是第10页,共28页现在学习的是第11页,共28页利用利用FFTFFT进行功率谱的噪声分析进行功率谱的噪声分析已知带有测量噪声信号已知带有测量噪声信号 其中其中f1=50Hz,f2=120Hz,为均值为零、方差为为均值为零、方差为1的的随机信号,采样频率为随机信号,采样频率为1000Hz,数据点数,数据点数N=512。试绘制信号的频谱图和功率谱图。试绘制信号的频谱图和功率谱图。现在学习的是第12页,共28页t=0:0.001:0.6;t=0:0.001:0.6;x=sin(2*pi*50*t)+sin(2*pi*120*t);x=sin(2*pi*50*t)+sin(2*pi*120*t);y=x+2*randn(1,length(t);y=x+2*randn(1,length(t);Y=fft(y,512);Y=fft(y,512);P=Y.*conj(Y)/512;%P=Y.*conj(Y)/512;%求功率求功率f=1000*(0:255)/512;f=1000*(0:255)/512;subplot(2,1,1);subplot(2,1,1);plot(y);plot(y);subplot(2,1,2);subplot(2,1,2);plot(f,P(1:256);plot(f,P(1:256);现在学习的是第13页,共28页现在学习的是第14页,共28页序列长度和序列长度和FFTFFT的长度对信号频谱的影响。的长度对信号频谱的影响。已知信号已知信号 其中其中f1=15Hz,f2=40Hz,f1=15Hz,f2=40Hz,采样频率为采样频率为100100Hz.Hz.在下列情况下绘制其幅频谱。在下列情况下绘制其幅频谱。Ndata=32,Nfft=32;Ndata=32,Nfft=32;Ndata=32,Nfft=128;Ndata=32,Nfft=128;现在学习的是第15页,共28页fs=100;Ndata=32;Nfft=32;n=0:Ndata-1;t=n/fs;x=0.5*sin(2*pi*15*t)+2*sin(2*pi*40*t);y=fft(x,Nfft);mag=abs(y);f=(0:length(y)-1)*fs/length(y);subplot(2,1,1)plot(f(1:Nfft/2),mag(1:Nfft/2)title(Ndata=32,Nfft=32)现在学习的是第16页,共28页Nfft=128;n=0:Ndata-1;t=n/fs;x=0.5*sin(2*pi*15*t)+2*sin(2*pi*40*t);y=fft(x,Nfft);mag=abs(y);f=(0:length(y)-1)*fs/length(y);subplot(2,1,2)plot(f(1:Nfft/2),mag(1:Nfft/2)title(Ndata=32,Nfft=128)现在学习的是第17页,共28页现在学习的是第18页,共28页快速傅立叶逆变换(快速傅立叶逆变换(IFFTIFFT)函数调用格式函数调用格式 y=ifft(x)y=ifft(x)y=ifft(x,N)y=ifft(x,N)当当N N小于小于x x长度时,对长度时,对x x进行截断,当进行截断,当N N大于大于x x长度时,长度时,对对x x进行补零。进行补零。现在学习的是第19页,共28页对信号对信号 进行进行DFTDFT,对其结果进行,对其结果进行IDFTIDFT,并将,并将IDFTIDFT的结果和原信号进行比较。的结果和原信号进行比较。f1=40Hz f1=40Hz f2=15Hzf2=15HzFs=100Hz Fs=100Hz 现在学习的是第20页,共28页fs=100;N=128;n=0:N-1;t=n/fs;fs=100;N=128;n=0:N-1;t=n/fs;x=sin(2*pi*40*t)+sin(2*pi*15*t);x=sin(2*pi*40*t)+sin(2*pi*15*t);subplot(2,2,1)subplot(2,2,1)plot(t,x)plot(t,x)title(original signal)title(original signal)y=fft(x,N);y=fft(x,N);mag=abs(y);mag=abs(y);f=(0:length(y)-1)*fs/length(y);f=(0:length(y)-1)*fs/length(y);subplot(2,2,2)subplot(2,2,2)plot(f,mag)plot(f,mag)title(FFT to original signal)title(FFT to original signal)现在学习的是第21页,共28页xifft=ifft(y);xifft=ifft(y);magx=real(xifft);magx=real(xifft);ti=0:length(xifft)-1/fs;ti=0:length(xifft)-1/fs;subplot(2,2,3)subplot(2,2,3)plot(ti,magx);plot(ti,magx);title(signal from IFFT)title(signal from IFFT)yif=fft(xifft,N);yif=fft(xifft,N);mag=abs(yif);mag=abs(yif);subplot(2,2,4)subplot(2,2,4)plot(f,mag)plot(f,mag)title(FFT to signal from IFFT)title(FFT to signal from IFFT)现在学习的是第22页,共28页现在学习的是第23页,共28页线性卷积的线性卷积的FFTFFT算法算法在在MATLABMATLAB实现卷积的函数为实现卷积的函数为CONVCONV,对于,对于N N值较小的向量,这是值较小的向量,这是十分有效的。对于十分有效的。对于N N值较大的向量卷积可用值较大的向量卷积可用FFTFFT加快计算速度。加快计算速度。由由DFTDFT性质可知,若性质可知,若DFTxDFTx1 1(n)=X(n)=X1 1(k),DFTx(k),DFTx2 2(n)=X(n)=X2 2(n)(n)则则 若若DFTDFT和和IDFTIDFT均采用均采用FFTFFT和和IFFTIFFT算法,可提高卷积速度。算法,可提高卷积速度。现在学习的是第24页,共28页计算计算x1(n)x1(n)和和x2(n)x2(n)的线性卷积的的线性卷积的FFTFFT算法可由下算法可由下面步骤实现面步骤实现计算计算X1(k)=FFTx1(n);计算计算X2(k)=FFTx2(n);计算计算Y(k)=X1(k)X2(k);计算计算x1(n)*x2(n)=IFFTY(k).现在学习的是第25页,共28页用函数用函数convconv和和FFTFFT计算同一序列的卷积,比较计算同一序列的卷积,比较其计算时间。其计算时间。L=5000;N=L*2-1;n=1:L;L=5000;N=L*2-1;n=1:L;x1=0.5*n;x2=2*n;x1=0.5*n;x2=2*n;t0=clock;yc=conv(x1,x2);t0=clock;yc=conv(x1,x2);conv_time=conv_time=etime(clock,t0)etime(clock,t0)t0=clock;t0=clock;yf=ifft(fft(x1,N).*fft(x2,N);yf=ifft(fft(x1,N).*fft(x2,N);fft_time=etime(clock,t0)fft_time=etime(clock,t0)clockclock函数读取瞬时时钟函数读取瞬时时钟etime(t1,t2)etime(t1,t2)函数计算时刻函数计算时刻t1,t2t1,t2间所经历的时间。间所经历的时间。现在学习的是第26页,共28页四、实验报告要求四、实验报告要求 简述实验目的、原理简述实验目的、原理对于对于8 8点点FFTFFT的显示,讨论其特点。的显示,讨论其特点。与离散卷积结果相比较,讨论快速卷积方法的优越性。与离散卷积结果相比较,讨论快速卷积方法的优越性。现在学习的是第27页,共28页感感谢谢大大家家观观看看现在学习的是第28页,共28页