实验二FFT谱分析

熟悉傅里叶变换的性质熟悉常见信号的傅里叶变换了解傅里叶变换的MATLAB实现方法连续时间信号直接调用专用函数法傅里叶变换的数值计算实现法离散时间信号正/反Z变换离散系统的频率特性直接调用专用函数法傅里叶变换：F=fourier(f)对f(t)进行傅里叶变换，其结果为F(w)F=fourier(f,v)对f(t)进行傅里叶变换，其结果为F(v)F=fourier(f,u,v)对f(u)进行傅里叶变换，其结果为F(v)傅里叶反变换f=ifourier(F)对F(w)进行傅里叶反变换，其结果为f(x)f=ifourier(F,U)对F(w)进行傅里叶反变换，其结果为f(u)f=ifourier(F,v,u)对F(v)进行傅里叶反变换，其结果为f(u)注意：在调用函数fourier()及ifourier()之前，要用syms命令对所有需要用到的变量进行说明，即要将这些变量说明成符号变量。对fourier()中的f及ifourier()中的F也要用符号定义符sym将其说明为符号表达式。注意：采用fourier()及ifourier()得到的返回函数，仍然为符号表达式。在对其作图时要用ezplot()函数，而不能用plot()函数。注意：fourier()及ifourier()函数的应用有很多局限性，如果在返回函数中含有δ(ω)等函数，则ezplot()函数也无法作出图来。另外，在用fourier()函数对某些信号进行变换时，其返回函数如果包含一些不能直接表达的式子，则此时当然也就无法作图了。另一个局限是在很多场合，尽管原时间信号f(t)是连续的，但却不能表示成符号表达式，此时只能应用下面介绍的数值计算法来进行傅氏变换了，当然，大多数情况下，用数值计算法所求的频谱函数只是一种近似值。例1：求函数的傅里叶反变换f(t)symstw%定义两个符号变量t,wFw=sym('1/(1+w^2)');%定义频谱函数F(jw)ft=ifourier(Fw,w,t);%对频谱函数F(jw)进行傅氏反变换运行结果：ft=1/2*exp(-t)*Heaviside(t)+1/2*exp(t)*Heaviside(-t)21()1Fj傅里叶变换的数值计算实现法例2：用数值计算法实现门函数的傅里叶变换，并画出幅度频谱图.()(1)(1)ftttMATLAB程序如下：R=0.02;%取样间隔τ=0.02t=-2:R:2;%t为从-2到2，间隔为0.02的行向量%有201个样本点ft=[zeros(1,50),ones(1,101),zeros(1,50)];%产生f(t)的样值矩阵（即f(t)%的样本值组成的行向量）W1=10*pi;%取要计算的频率范围M=500;k=0:M;w=k*W1/M;%频域采样数为M,w为频率正半轴的采样点Fw=ft*exp(-j*t'*w)*R;%求傅氏变换FRw=abs(Fw);%取振幅W=[-fliplr(w),w(2:501)];%形成负半轴和正半轴的2M+1个频率点WFW=[fliplr(FRw),FRw(2:501)];%形成对应于2M+1个频率点的值subplot(2,1,1);plot(t,ft);grid;%画出原时间函数f(t)的波形，并加网格xlabel('t');ylabel('f(t)');%坐标轴标注title('f(t)=u(t+1)-u(t-1)');%文本标注subplot(2,1,2);plot(W,FW);grid;%画出振幅频谱的波形,并加网格xlabel('W');ylabel('F(W)');%坐标轴标注title('f(t)的振幅频谱图');%文本标注在MATLAB语言中有专门对信号进行正反Z变换的函数ztrans()和itrans()。其调用格式分别如下：F=ztrans(f)对f(n)进行Z变换，其结果为F(z)F=ztrans(f,v)对f(n)进行Z变换，其结果为F(v)F=ztrans(f,u,v)对f(u)进行Z变换，其结果为F(v)f=itrans(F)对F(z)进行Z反变换，其结果为f(n)f=itrans(F,u)对F(z)进行Z反变换，其结果为f(u)f=itrans(F,v,u)对F(v)进行Z反变换，其结果为f(u)注意：在调用函数ztran()及iztran()之前，要用syms命令对所有需要用到的变量（如t,u,v,w）等进行说明，即要将这些变量说明成符号变量。例3．用MATLAB求出离散序列的Z变换MATLAB程序如下：symskzf=0.5^k;%定义离散信号Fz=ztrans(f)%对离散信号进行Z变换运行结果如下：Fz=2*z/(2*z-1)()(0.5)()kfkk例4．已知一离散信号的Z变换为，求出它所对应的离散信号f(k)MATLAB程序如下：symskzFz=2*z/(2*z-1);%定义Z变换表达式fk=iztrans(Fz,k)%求反Z变换运行结果如下：fk=(1/2)^k2()21zFzz离散系统的频率特性MATLAB为我们提供了专门用于求解离散系统频率响应的函数fregz()，其调用格式如下：[H,w]=fregz(B,A,N)其中，B和A分别是表示待分析的离散系统的系统函数的分子，分母多项式的向量，N为正整数，返回向量H则包含了离散系统频率响应函数在范围内的N个频率等分点的值。向量w则包含范围内的N个频率等分点。在默认情况下N=512。[H,w]=fregz(B,A,N,'whole')其中，B，A和N的意义同上，而返回向量H包含了频率响应函数在范围内N个频率等分点的值。()jHe()jHe0~2离散系统的频率特性由于调用fregz()函数只能求出离散系统频率响应的数值，不能直接绘制曲线图，因此，我们可以先用fregz()函数求出系统频率响应的值，然后再利用MATLAB的abs()和angle()函数以及plot()命令，即可绘制出系统在或范围内的幅频特性和相频特性曲线。例5．用MATLAB计算前面离散系统在频率范围内200个频率等分点的频率响应值，并绘出相应的幅频特性和相频特性曲线。MATLAB程序如下：A=[10];B=[1-0.5];[H,w]=freqz(B,A,200);[H,w]=freqz(B,A,200,'whole');%求出对应范围内200个频率点的频率响%应样值HF=abs(H);%求出幅频特性值HX=angle(H);%求出相频特性值subplot(2,1,1);plot(w,HF)%画出幅频特性曲线subplot(2,1,2);plot(w,HX)%画出相频特性曲线运行结果如下：运行结果分析：从该系统的幅频特性曲线可以看出，该系统呈高通特性，是一阶高通滤波器。如果用FFT对模拟信号进行谱分析，首先要把模拟信号转换成数字信号，转换时要求知道模拟信号的最高截止频率，以便选择满足采样定理的采样频率。一般选择采样频率是模拟信号中最高频率的3～4倍。另外要选择对模拟信号的观测时间，如果采样频率和观测时间确定，则采样点数也确定了。这里观测时间和对模拟信号进行谱分析的分辨率有关，最小的观测时间和分辨率成倒数关系。要求选择的采样点数和观测时间大于它的最小值。用FFT作谱分析时，要求做FFT的点数服从２的整数幂，这一点在上面选择采样点数时可以考虑满足，即使满足不了，可以通过在序列尾部加０完成。如果要进行谱分析的模拟信号是周期信号，最好选择观测时间是信号周期的整数倍。如果不知道信号的周期，要尽量选择观测时间长一些，以减少截断效应的影响。用DFT对连续信号作谱分析。已知xa(t)=cos(200*pi*t)+sin(100*pi*t)+cos(50*pi*t);选取不同的截取长度Tp，观察用DFT进行频谱分析时存在的截取效应（频谱泄漏和谱间干扰）。在计算机上用DFT对模拟信号进行谱分析时，只能以有限大的采样频率fs对模拟信号采样。对有限点样本序列（等价于截取模拟信号一段进行采样）作DFT变换得到模拟信号的近似频谱clear;closeall;fs=400;T=1/fs;Tp=0.04;N=Tp*fs;N1=[N,4*N,8*N];%三种长度0.04s4*0.04s8*0.04s%矩形窗截断form=1:3n=1:N1(m);xn=cos(200*pi*n*T)+sin(100*pi*n*T)+cos(50*pi*n*T);Xk=fft(xn,4096);fk=fs*[0:4095]/4096;subplot(3,2,2*m-1);plot(fk,abs(Xk)/max(abs(Xk)));ifm==1title('矩形窗截断');endend%加海明窗截断form=1:3n=1:N1(m);wn=hamming(N1(m));xn=(cos(200*pi*n*T)+sin(100*pi*n*T)+cos(50*pi*n*T)).*wn';Xk=fft(xn,4096);fk=fs*[0:4095]/4096;subplot(3,2,2*m)plot(fk,abs(Xk)/max(abs(Xk)));ifm==1title('hamming窗截断');end1、求出下列离散序列的Z变换1122()()cos()()kkfkk1122()()cos()()kkfkk1122()()cos()()kkfkk11()()cos()()22kkfkk2()()(5)fkkk2、分别以变换区间N=8，N=16进行FFT，画出相应的幅频特性曲线。14,()3,0,nxnnnnn其他74302()cos4xnn