MATLAB实验二--傅里叶分析及应用

1实验二傅里叶分析及应用一、实验目的（一）掌握使用Matlab进行周期信号傅里叶级数展开和频谱分析1、学会使用Matlab分析傅里叶级数展开，深入理解傅里叶级数的物理含义2、学会使用Matlab分析周期信号的频谱特性（二）掌握使用Matlab求解信号的傅里叶变换并分析傅里叶变换的性质1、学会运用Matlab求连续时间信号的傅里叶变换2、学会运用Matlab求连续时间信号的频谱图3、学会运用Matlab分析连续时间信号的傅里叶变换的性质（三）掌握使用Matlab完成信号抽样并验证抽样定理1、学会运用MATLAB完成信号抽样以及对抽样信号的频谱进行分析2、学会运用MATLAB改变抽样时间间隔，观察抽样后信号的频谱变化3、学会运用MATLAB对抽样后的信号进行重建二、实验条件Win7系统，MATLABR2015a三、实验内容1、分别利用Matlab符号运算求解法和数值计算法求下图所示信号的FT，并画出其频谱图（包括幅度谱和相位谱）[注：图中时间单位为：毫秒(ms)]。2符号运算法数值运算法Code:ft=sym('(t+2)*(heaviside(t+2)-heaviside(t+1))+(heaviside(t+1)-heaviside(t-1))+(2-t)*(heaviside(t-1)-heaviside(t-2))');fw=simplify(fourier(ft));subplot(2,1,1);ezplot(abs(fw));gridon;title('ampspectrum');phi=atan(imag(fw)/real(fw));subplot(2,1,2);ezplot(phi);gridon;title('phasespectrum');Code:dt=0.01;t=-2:dt:2;ft=(t+2).*(uCT(t+2)-uCT(t+1))+(uCT(t+1)-uCT(t-1))+(2-t).*(uCT(t-1)-uCT(t-2));N=2000;k=-N:N;w=pi*k/(N*dt);fw=dt*ft*exp(-i*t'*w);fw=abs(fw);plot(w,fw),gridon;axis([-2*pi2*pi-13.5]);300.511.522.5012345t(20exp(-3t)heaviside(t)-8exp(-5t)heaviside(t))/(2)2、试用Matlab命令求j54-j310)F(j的傅里叶反变换，并绘出其时域信号图。两个单边指数脉冲的叠加3、已知门函数自身卷积为三角波信号，试用Matlab命令验证FT的时域卷积定理。Code：symst;fw=sym('10/(3+i*w)-4/(5+i*w)');ft=ifourier(fw,t);ezplot(ft),gridon;Code：f=sym('heaviside(t+1)-heaviside(t-1)');fw=simplify(fourier(f));F=fw.*fw;subplot(211);ezplot(abs(F),[-9,9]),gridontitle('FW^2')tri=sym('(t+2)*heaviside(t+2)-2*t*heaviside(t)+(t-2)*heaviside(t-2)');Ftri=fourier(tri);F=simplify(Ftri);subplot(212);ezplot(abs(F),[-9,9]),gridon;title('triFT')44、设有两个不同频率的余弦信号，频率分别为Hzf1001，Hzf38002；现在使用抽样频率Hzfs4000对这三个信号进行抽样，使用MATLAB命令画出各抽样信号的波形和频谱，并分析其频率混叠现象Code：f1=100;%f1=100hzts=1/4000;%sample=4000hzdt=0.0001;t1=-0.007:dt:0.007;ft=cos(2*f1*pi*t1);subplot(221);plot(t1,ft),gridon;axis([-0.0060.006-1.51.5])xlabel('Time/s'),ylabel('f(t)')title('Cosinecurve');N=5000;k=-N:N;w=2*pi*k/((2*N+1)*dt);fw=ft*dt*exp(-1i*t1'*w);subplot(222);plot(w,abs(fw));gridon;axis([-200002000000.005]);xlabel('\omega'),ylabel('f(w)')title('Cosfreqspectrum');t2=-0.007:ts:0.007;fst=cos(2*f1*pi*t2);subplot(223);plot(t1,ft,':'),holdonstem(t2,fst),gridon;axis([-0.0060.006-1.51.5])xlabel('Time/s'),ylabel('fs(t)')title('Samplesignal');holdofffsw=ts*fst*exp(-1i*t2'*w);subplot(224);plot(w,abs(fsw)),gridonaxis([-200002000000.006])xlabel('\omega'),ylabel('fsw')title('Samplefreqspectrum');5-505x10-3-101Time/sf(t)Cosinecurve-2-1012x104012345x10-3f(w)Cosfreqspectrum-505x10-3-101Time/sfs(t)Samplesignal-2-1012x1040246x10-3fswSamplefreqspectrumf1=100Hz将代码中f1设为3800即可↓f2=3800Hz-505x10-3-101Time/sf(t)Cosinecurve-2-1012x104012345x10-3f(w)Cosfreqspectrum-505x10-3-101Time/sfs(t)Samplesignal-2-1012x1040246x10-3fswSamplefreqspectrum6-505-0.500.51Sa(t)-2002000.511.5Sa(t)freqspectrum-505-0.500.51Samplingsignal-5005000.511.5spectrumofSamplingsignal5、结合抽样定理，利用MATLAB编程实现)(tSa信号经过冲激脉冲抽样后得到的抽样信号tfs及其频谱[建议：冲激脉冲的周期分别取4*pi/3s、pis、2*pi/3s三种情况对比]，并利用tfs构建)(tSa信号。（**改动第一行代码即可）冲激脉冲的周期=4*pi/3sTs=4/3;%impulseperiod=4*pi/3t1=-5:0.01:5;ft=sinc(t1);subplot(2,2,1)plot(t1,ft),gridonaxis([-66-0.51.2])title('Sa(t)')N=500;k=-N:N;W=pi*k/(N*0.01);Fw=0.01*ft*exp(-1i*t1'*W);subplot(2,2,2)plot(W,abs(Fw)),gridonaxis([-3030-0.051.5])title('Sa(t)freqspectrum')t2=-5:Ts:5;fst=sinc(t2);subplot(2,2,3)plot(t1,ft,':'),holdonstem(t2,fst),gridonaxis([-66-0.51.2])title('Samplingsignal')Fsw=Ts*fst*exp(-1i*t2'*W);subplot(2,2,4)plot(W,abs(Fsw)),gridonaxis([-5050-0.051.5])title('spectrumofSamplingsignal')7冲激脉冲的周期=pis冲激脉冲的周期=2*pi/3s-505-0.500.51Sa(t)-2002000.511.5Sa(t)freqspectrum-505-0.500.51Samplingsignal-5005000.511.5spectrumofSamplingsignal-505-0.500.51Sa(t)-2002000.511.5Sa(t)freqspectrum-505-0.500.51Samplingsignal-5005000.511.5spectrumofSamplingsignal8-4-3-2-101234-1-0.500.511.52Originalwave-4-202400.20.40.60.811classH-wave-4-202400.20.40.60.813classH-wave-4-202400.20.40.60.8113classH-wave-4-202400.20.40.60.8149classH-wave6、已知周期三角信号如下图所示[注：图中时间单位为：毫秒(ms)]：（1）试求出该信号的傅里叶级数[自己求或参见课本P112或P394]，利用Matlab编程实现其各次谐波[如1、3、5、13、49]的叠加，并验证其收敛性；a0=12;an=4(nπ)2sin2(nπ2);bn=0谐波幅度收敛速度：1n2原始波形：第k阶谐波波形9-4-202400.20.40.60.81The1timessuperpose-4-202400.20.40.60.81The3timessuperpose-4-202400.20.40.60.81The13timessuperpose-4-202400.20.40.60.81The49timessuperpose前K次谐波的叠加Code：figure(1);t=-2*pi:0.001:2*pi;f=abs(sawtooth(0.5*pi*t,0.5));plot(t,f),gridon;axis([-4,4,-1,2])title('Originalwave');nclass=[1,3,13,49];figure(2);N=4;a0=1/2;fork=1:Nn=nclass(k);an=4./((n*pi).^2);ft=an*cos(pi*n'*t);ft=ft+a0;subplot(2,2,k);plot(t,ft);axis([-4,4,0,1])title([num2str(nclass(k)),'classH-wave']);endfigure(3);N=4;a0=1/2;fork=1:Nn=1:2:nclass(k);an=4./((n*pi).^2);ft=an*cos(pi*n'*t);ft=ft+a0;subplot(2,2,k);plot(t,ft);axis([-4,4,0,1])title(['The',num2str(nclass(k)),'timessuperpose']);end10-20-15-10-50510152000.20.40.60.81-20-15-10-50510152000.20.40.60.81（2）用Matlab分析该周期三角信号的频谱[三角形式或指数形式均可]。当周期三角信号的周期（如由2ms1ms或由2ms4ms）和宽度（如2ms1ms）分别变化时，试观察分析其频谱的变化。周期为2ms周期为1msdt=0.01;t=-4:dt:4;ft=(t=-1&t0).*(t+1)+(t0&t=1).*(1-t);%subplot(2,1,1)%plot(t,ft);gridonn=2000;k=-n:n;w=pi*k/(n*dt);f=dt*ft*exp(-i*t'*w);f=abs(f);%subplot(2,1,2)plot