《实验三DFT和FFT频谱分析》由会员分享，可在线阅读，更多相关《实验三DFT和FFT频谱分析（12页珍藏版）》请在金锄头文库上搜索。
1、实验三 DFT 和 FFT 频谱分析一、实验目的1.掌握 DFT 频谱分析的原理与编程方法。2理解FFT算法的编程思想。2熟练掌握利用FFT对信号作频谱分析，包括正确地进行参数选择、画频谱及读频谱图。3利用FFT频谱分析进行快速卷积和太阳黑子周期性检测。二、实验环境I.Windows xp以上操作系统2安装 MATLAB2007a 软件三、实验原理1.离散傅里叶变换(DFT) 设序列为x(n),长度为N,则N -1X(ejrok)=DFTx(n)=为x(n) e-jwk n,n=02 n其中叭=k (k=0,1,2,M-1),通常MN,以便观察频谱的细节。IX(ejsk)|-x(n)的幅频谱。
2、kM2谱分析参数选择1) 设信号x(t)最高频率为fc,对其进行取样得x(n),根据取样定理，取样频率fs必须满足：fs=2fc。2) 设谱分辨率为F则最小记录时间；取样点数N2fC/F；为使用快速傅里叶变换(FFT)进行谱分析，N还须满足：N=2E (E为整数)。3用FFT计算信号x(n)的频谱。设x(n)为实信号快速傅里叶变换(FFT)是DFT的一种快速算法，其使得DFT的运算速度大为加快。1) 对信 号x(n)作N 点 FFT,得 频谱 X(k)(k=0N-1)X(k)=XR(k)+jXI(k) (k=0N/2-1)， XR(k) X(k)的实部；XI(k) X(k)的虚部。 Matla
3、b语句：Y=fft(xW 其中：x-x(n);Y-X(k)2) 幅频谱:IX(k)I=x2(K) + x2(K)，由于 x(n)为实信号，因此 IX(k)I对称，Matlab语句：abs(Y)iii)功率谱:PSD(k)=IX(k)I2/N=X(k)X*(k)/NMatlab语句：PSD=Y.*conj(Y)/N其中：conj(Y)- X*(k)X(k)的共轭4读频谱图频谱图中任意频率点k对应实际频率为:fk=fs/N*k。5用FFT实现线性卷积运算用FFT实现y (n)=x(n)*h(n)的步骤为：1)设x(n)及h(n)的长度分别为N和N2。为使循环卷积等于线性卷积，用补0的方 法使x(n
4、),h(n)长度均为N,则N须满足NN1+N2-1；为用FFT计算DFT,则N还须满足7=。2) 用 FFT计算X(k),H(k) (N点)。3) Y(k)=_ffl； y(n)= ifftY(K)四、实验内容1根据公式设计DFT原理程序, %DFT/IDFT程序 DFTm clc clear xn=input(x(n)= );M=length(xn); N=input(变换区间 N=); xn=xn zeros(1,N-M); n=0:N-1;k=0:N-1; nk=n*k; wn=exp(-j*2*pi/N); wnK=wn.Ank; xk=xn*wnK并计算：x(n)=1,1,1,1的4
5、,16,64点DFT并绘图。%输入序列x(n)=1 1 1 1%x(n)的长度M%变换区间N%补0,使xn长度为N%旋转因子wn%作区(血)的DFT=xksubplot(211);stem(k,abs(xk), );grid on; %显示xk的幅频谱(离散曲线)subplot(212);plot(k,abs(xk);grid on;%显示xk的幅频谱(连续曲线) 运行结果：问：由此得出怎样的结论？答：n越大越接近原来的dft2理解DIT-FFT算法原理程序，并用它计算X(k)=FFTR/n),分别取N=4,8,16和64,绘出幅频谱IX(k)l。%程序 DIT.m clear clc x=i
6、nput(x= ); N=input(N= ); x(length(x)+1:N)=zeros(1,N-length(x); l=log2(N); x1=zeros(1,N); for j1=1:Nx1(j1)=x(bin2dec(fliplr(dec2bin(j1-1,l)+1); end %FFT(DIT)%M=2; while(M=N)W=exp(-2*j*pi/M);V=1;for k=0:1:M/2-1for i=0:M:N-1 p=k+i; q=p+M/2; A=x1(p+1); B=x1(q+1)*V; x1(p+1)=A+B; x1(q+1)=A-B;end V=V*W;end
7、M=2*M; end % subplot(211);stem(x,.);grid on; title(x(n);%输入序列%做ft的占数% 补0 x(1:N)% 倒序%旋转因子W%k每级蝶形运算旋转因子的个数%伪各群的首序号%本级蝶形运算，x1最终存放X(k)%旋转因子W的变化%第“级%画图%题% 画图%题subplot(212);stem(abs(x1),.);grid on; title(lX(k)l);x(n)|X(k)|x(n)|X(k)|x(n)|X(k)|x(n)|X(k)|3.FFT谱分析设信号为x(t)=sin(2nf1t)+sin(2nf2t)+随机噪声，f1=50Hz, f
8、2=120Hz，以取样频率fs=1kHz对x(t)进行取样，样本长度tp=025s，得x(n),对x(n)作256点FFT,得频谱X(k),画原信号 x(n),幅频谱IX(k)l以及功率谱PSD(k)对信号进行谱分析。%程序 pufenxi.mclearclc fs=1000;t=0:1/fs:0.25;N=256; f1=50;f2=120; s=sin(2*pi*f1*t)+sin(2*pi*f2*t); x=s+randn(size(t);Y=fft(x,N);PSD=Y.*conj(Y)/N;f=fs/N*(0:N/2-1);subplot(311);plot(x); subplot(
9、312);plot(f,abs(Y(1:N/2); subplot(313);plot(f,PSD(1:N/2);%时间范围%做ft的占数%信号频率%产生 x(n)%信号+噪声x(n)%对x做N占ft做功率谱将频率占转化为实际频率画原信号画幅度谱(N/2点)画功率谱(N/2点) 1501005050-5050100150 200 2503000 50 100 150 200 250 300 350 400 450 50002)回答下列问题：i) 观察幅频谱图，可以发现，信号x(n)含有的两个频率分量分别是50.8 Hz和 121.1 Hz。ii) 在该程序中的“f=fs/N*(0:N/2-l)
10、” 下添加“k=0:N/2-1”， “plot(bs(Y(1:N/2)” 改为“plot(k,abs(Y(1:N/2)”重新运行该程序并观察幅频谱图，图中两峰值对应的下标分别是13和，31。它们的含义为主频点。再将该程序中的N改为512,重新运行该程序并观察幅频谱图，这时图中1) 画出图形窗口显示的图形，并两峰值对应的下标分别是_26_和61。结果是否和上面的相同？不同为什么？N不同。iii) 本例的频谱分辨率F是一3.9 z,改变f2=60Hz,问：在幅频谱中，能否分辨和f2对应的频率分量？不能。为什么？间隔小干频谱分辨率。再改变f2=52Hz,问：在幅频谱中，能否分辨勺和孑?对应的频率分量
11、？不能 。为什么？间隔小干频谱分辨率。再改变f2=600Hz,在幅频谱中，f2对应的频率分量出现 398.45Hz；问：在fs=1000Hz的情况下，能否正确检测f2对应的频率分量？ 。为什么？ 不符合采样定理。为了正确检测f2对应的频率分量，则fs至少取多少Hz?1200Hz。在该程序中改变fs,验证你的结论。iv) 比较幅频谱和功率谱，可以发现功率谱具有的特性。4. FFT实现任意两个序列的快速卷积。%程序 fftjuanji.mclearclc x1=input(x1=);x2=input(x2=); N1=length(x1);N2=length(x2); E=ceil(log2(N1
12、+N2-1);N=2AE;x1=x1,zeros(1,N-N1);x2=x2,zeros(1,N-N2);X1=fft(x1,N);X2=fft(x2,N);Y=X1.*X2;y=ifft(Y,N)结果分析：1)回到MATLAB窗口，键入： x1=1 1 1, x2=1 2,回车。结果:y=%输入序列%序列x1(n),x2(n)的长度%ceil-向+8方向取整%做FFT的点数%补零使x1,x2长度为N%对x1做N点的fft% 数列X1和X2的乘积%对y做N点的ifft13322)问：可用Matlab中的什么函数验算上述卷积结果？ Y=conv(x1,x2)5. 利用谱分析观察太阳黑子周期性。以
13、100年中记录到的太阳黑子出现次数为信号x(n),对x(n)作功率谱，从中观察太阳黑子周期 性。% 程序 taiyangheizi.mclear clcx=101 82 66 35 31 7 20 92 154 125 85 68 38 23 10 24 83 . 132 131 118 90 67 60 47 41 21 16 6 4 7 14 34 45 43 48 .42 28 10 8 2 0 1 5 12 14 35 46 41 30 24 16 7 4 2 8 .17 36 50 62 67 71 48 28 8 13 57 122 138 103 86 63 37 24 .11
14、15 40 62 98 124 96 66 64 54 39 21 7 4 23 55 94 96 .77 59 44 47 30 16 7 37 74;%100年中太阳黑子出现的次数subplot(211);plot(x)%画区(血)N=128; fs=1;%fs=1Hz,N=128 点s=x-mean(x);%对区作零均值化处理(去除直流分量)Y= fft(s,n);%对$做N 点 ftPSD=_Y*conj(Y)/N;% 做功率谱 PSDf= fs/N*(0:N/2-1);%将频率定标为实际频率fsubplot(212); plot(f,PSD(1:N/2); %画功率谱(N/2 点)1)填写空格中的画图语句并绘出结果图形。10 20 30 40
