数字信号作业1.doc

1、作业 1：x(t)=sin(2*pi*480*t)+2*sin(2*pi*490*t)，(I) 若要在 x(t)的频谱上分辨出这两个频率，且保证频谱不混叠，采样频率和信号时长至少要取多大？(II) 以 1KHz 的采样率采 128 个点，做 128 点 FFT，画出频谱图； 以 1KHz 的采样率采 256 个点，做 256 点 FFT，画出频谱图； 以 1KHz 的采样率采 128 个点，后面添加 128 个 0，做 256 点 FFT，画频谱图； 以 2KHz 的采样率采 128 个点，做 128 点 FFT，画出频谱图； 以 2KHz 的采样率采 256 个点，做 256 点 FFT，画

2、出频谱图； 以 2KHz 的采样率采 512 个点，做 512 点 FFT，画出频谱图；要求：频谱图横轴坐标以 Hz 为单位，比较这 5 种情况下信号的频谱有何异同， 并说说你从这些现象中明白了什么道理。(I) 解：若要在 x(t)的频谱上分辨出这两个频率，且保证频谱不混叠，采样频率 fs=980Hz；信号时长 tp=1/F=1/10Hz=0.1s;(II) 以 1KHz 的采样率采 128 个点，做 128 点 FFT，画出频谱图；用 Matlab 编写的程序为clf;clear;clc;fs=1000;%采样频率和数据点数N=128; n=0:N-1;t=n/fs; %时间序列x=sin(

3、2*pi*480*t)+2*sin(2*pi*490*t); %信号y=fft(x,N); %对信号进行快速 Fourier 变换mag=abs(y); %求得 Fourier 变换后的振幅f=n*fs/N; %频率序列subplot(1,1,1)plot(f,mag); %绘出随频率变化的振幅xlabel(频率/Hz);ylabel(振幅);title(N=128);grid on;绘出的图像如图：0 100 200 300 400 500 600 700 800 900 1000020406080100120140位 位 /Hz位位N=128 以 1KHz 的采样率采 256 个点，做 2

4、56 点 FFT，画出频谱图用 Matlab 编写的程序为clf;clear;clc;fs=1000;%采样频率和数据点数N=256; n=0:N-1;t=n/fs; %时间序列x=sin(2*pi*480*t)+2*sin(2*pi*490*t); %信号y=fft(x,N); %对信号进行快速 Fourier 变换mag=abs(y); %求得 Fourier 变换后的振幅f=n*fs/N; %频率序列subplot(1,1,1)plot(f,mag); %绘出随频率变化的振幅xlabel(频率/Hz);ylabel(振幅);title(N=256);grid on;绘出的图像如图：0 1

5、00 200 300 400 500 600 700 800 900 1000020406080100120140160180位 位 /Hz位位N=256 以 1KHz 的采样率采 128 个点，后面添加 128 个 0，做 256 点 FFT，画频谱图；用 Matlab 编写的程序为clf;clear;clc;fs=1000;%采样频率和数据点数N=128; n=0:N-1;t=n/fs; %时间序列x=sin(2*pi*480*t)+2*sin(2*pi*490*t); %信号x1=zeros(1,256);x1(1:N)=x;y=fft(x1,256); %对信号进行快速 Fourier

6、 变换mag=abs(y); %求得 Fourier 变换后的振幅n1=0:255;f=n1*fs/256; %频率序列subplot(1,1,1)plot(f,mag); %绘出随频率变化的振幅xlabel(频率/Hz);ylabel(振幅);title(N=128 N1=256);grid on; 绘出的图像为：0 100 200 300 400 500 600 700 800 900 1000020406080100120140位 位 /Hz位位N=128 N1=256 以 2KHz 的采样率采 128 个点，做 128 点 FFT，画出频谱图；用 Matlab 编写的程序为：clf;c

7、lear;clc;fs=2000;%采样频率和数据点数N=128; n=0:N-1;t=n/fs; %时间序列x=sin(2*pi*480*t)+2*sin(2*pi*490*t); %信号y=fft(x,N); %对信号进行快速 Fourier 变换mag=abs(y); %求得 Fourier 变换后的振幅f=n*fs/N; %频率序列subplot(1,1,1)plot(f,mag); %绘出随频率变化的振幅xlabel(频率/Hz);ylabel(振幅);title(N=128);grid on;绘出的图像为：0 200 400 600 800 1000 1200 1400 1600

8、1800 20000102030405060708090100位 位 /Hz位位N=128 以 2KHz 的采样率采 256 个点，做 256 点 FFT，画出频谱图；用 Matlab 编写的程序为：clf;clear;clc;fs=2000;%采样频率和数据点数N=256; n=0:N-1;t=n/fs; %时间序列x=sin(2*pi*480*t)+2*sin(2*pi*490*t); %信号y=fft(x,N); %对信号进行快速 Fourier 变换mag=abs(y); %求得 Fourier 变换后的振幅f=n*fs/N; %频率序列subplot(1,1,1)plot(f,mag

9、); %绘出随频率变化的振幅xlabel(频率/Hz);ylabel(振幅);title(N=256);grid on;绘出的图像为：0 200 400 600 800 1000 1200 1400 1600 1800 2000050100150200250位 位 /Hz位位N=256 以 2KHz 的采样率采 512 个点，做 512 点 FFT，画出频谱图；用 Matlab 编写的程序为：clf;clear;clc;fs=2000;%采样频率和数据点数N=512; n=0:N-1;t=n/fs; %时间序列x=sin(2*pi*480*t)+2*sin(2*pi*490*t); %信号y=

10、fft(x,N); %对信号进行快速 Fourier 变换mag=abs(y); %求得 Fourier 变换后的振幅f=n*fs/N; %频率序列subplot(1,1,1)plot(f,mag); %绘出随频率变化的振幅xlabel(频率/Hz);ylabel(振幅);title(N=512);grid on;绘出的图像为：0 200 400 600 800 1000 1200 1400 1600 1800 2000050100150200250300350400位 位 /Hz位位N=512从对这 6 幅图的比较中可以看出：在采样频率 fs 一定的情况下，采样有效点数 N 越大，信号在频域

11、上所占的脉宽越窄，频率混叠程度越小，频谱分辨率越大。若对采样点后面补零，可以使采样间隔更小，频谱形状更平滑，但不能提高分辨率。作业 2：x(t)= sin(2*pi*100*t)+ 0.5*cos(2*pi*1500*t) +0.5*sin(2*pi*2900*t)，(I) 分别设计巴特沃思模拟低通、带通、高通滤波器，将信号中的三个分量分离出来； (II) 分别设计切比雪夫 I、II 型模拟低通、带通、高通滤波器，将信号中的三个分量分离出来；画出信号分离前后的时域波形图和频谱图，画出滤波器的幅频、相频特性曲线（幅频曲线横轴用普通坐标，纵轴用对数坐标） ；自行确定滤波器的带宽和波纹，不允许用直接

12、设计的方法，用低通滤波器转换函数一，巴特沃兹模拟低通滤波器clc;clear;clf;Wp=800*2*pi;Ws=1200*2*pi;%通带阻带截止频率Rp=3;Rs=20;%通带波纹和阻带衰减N,Wc=buttord(Wp,Ws,Rp,Rs,s); %求得滤波器的最小阶数和截止频率 z,p,k=buttap(N); %设计Butterworth滤波器b,a=zp2tf(z,p,k); %将零点极点增益形式转换为传递 函数形式bt,at=lp2lp(b,a,Wc);w=linspace(1,5000,1000)*2*pi; %设置绘制频率响应的频率点H=freqs(bt,at,w); %计算

13、给定频率点的复数频率响应magH=abs(H);%-figure(1)subplot(2,1,1),plot(w/2/pi,20*log10(magH); %绘制幅频响应xlabel(频率 /Hz);ylabel(振幅/dB );grid on;title(幅频图);subplot(2,1,2),plot(w/2/pi,unwrap(angle(H);%绘制相频响应xlabel(频率 /Hz);ylabel(相位 /o);grid on;title(相频图)%-figure(2)dt=1/10000; %模拟信号采样间隔f1=100;f2=1500;f3=2900;%输入信号的三个频率成分t=

14、0:dt:0.04; %给定模拟时间段x=sin(2*pi*f1*t)+0.5*cos(2*pi*f2*t)+0.5*sin(2*pi*f3*t); %输入信号h=tf(bt,at); %滤波器在MATLAB系统中的表示y,t1=lsim(h,x,t); %模拟输出subplot(2,1,1),plot(t,x),title(输入信号 ) %绘出输入信号subplot(2,1,2),plot(t1,y) %绘制输出信号title(输出信号 ),xlabel(时间/s)%-figure(3)Fy=abs(fft(y);w=(-(length(Fy)-1)/2:(length(Fy)-1)/2)/

15、(length(Fy)*dt);plot(w(length(Fy)-1)/2+1:(length(Fy)-1),Fy(1:(length(Fy)-1)/2)*2/length(Fy) title(输出信号 频谱 )xlabel(频率 /Hz) 0 500 1000 1500 2000 2500 3000 3500 4000 4500 5000-100-500位 位 /Hz位位/dB位 位 位0 500 1000 1500 2000 2500 3000 3500 4000 4500 5000-10-50位 位 /Hz位位/o位 位 位0 0.005 0.01 0.015 0.02 0.025 0

16、.03 0.035 0.04-2-1012 位 位 位 位0 0.005 0.01 0.015 0.02 0.025 0.03 0.035 0.04-2-1012 位 位 位 位位 位 /s0 500 1000 1500 2000 2500 3000 3500 4000 4500 500000.20.40.60.81 位 位 位 位 位 位位 位 /Hz二、巴特沃兹模拟带通滤波器clc;clear;clf;Wp=1200,1800*2*pi;Ws=1000,2000*2*pi;%通带阻带截止频率Rp=3;Rs=20;%通带波纹和阻带衰减N,Wc=buttord(Wp,Ws,Rp,Rs,s);

17、%求得滤波器的最小阶数和截止频率 z,p,k=buttap(N); %设计Butterworth滤波器b,a=zp2tf(z,p,k); %将零点极点增益形式转换为传递 函数形式Wo=1500*2*pi; %中心 频率Bw=600*2*pi; %频带宽度bt,at=lp2bp(b,a,Wo,Bw); %频率转换w=linspace(1,5000,1000)*2*pi; %设置绘制频率响应的频率点H=freqs(bt,at,w); %计算给定频率点的复数频率响应magH=abs(H);%-figure(1)subplot(2,1,1),plot(w/2/pi,20*log10(magH); %绘

18、制幅频响应xlabel(频率 /Hz);ylabel(振幅/dB );grid on;title(幅频图);subplot(2,1,2),plot(w/2/pi,unwrap(angle(H);%绘制相频响应xlabel(频率 /Hz);ylabel(相位 /o);grid on;title(相频图)%-figure(2)dt=1/10000; %模拟信号采样间隔f1=100;f2=1500;f3=2900;%输入信号的三个频率成分t=0:dt:0.04; %给定模拟时间段x=sin(2*pi*f1*t)+0.5*cos(2*pi*f2*t)+0.5*sin(2*pi*f3*t); %输入信号

19、h=tf(bt,at); %滤波器在MATLAB系统中的表示y,t1=lsim(h,x,t); %模拟输出subplot(2,1,1),plot(t,x),title(输入信号 ) %绘出输入信号subplot(2,1,2),plot(t1,y) %绘制输出信号title(输出信号 ),xlabel(时间/s)%-figure(3)Fy=abs(fft(y);w=(-(length(Fy)-1)/2:(length(Fy)-1)/2)/(length(Fy)*dt);plot(w(length(Fy)-1)/2+1:(length(Fy)-1),Fy(1:(length(Fy)-1)/2)*2

20、/length(Fy) title(输出信号 频谱 )xlabel(频率 /Hz) 0 500 1000 1500 2000 2500 3000 3500 4000 4500 5000-400-2000200位 位 /Hz位位/dB位 位 位0 500 1000 1500 2000 2500 3000 3500 4000 4500 5000-20-10010位 位 /Hz位位/o位 位 位0 0.005 0.01 0.015 0.02 0.025 0.03 0.035 0.04-2-1012 位 位 位 位0 0.005 0.01 0.015 0.02 0.025 0.03 0.035 0.0

21、4-0.500.51 位 位 位 位位 位 /s0 500 1000 1500 2000 2500 3000 3500 4000 4500 500000.10.20.30.40.5 位 位 位 位 位 位位 位 /Hz三、巴特沃兹模拟高通滤波器clc;clear;clf;Wp=2600*2*pi;Ws=2000*2*pi;%通带阻带截止频率Rp=3;Rs=25;%通带波纹和阻带衰减N,Wc=buttord(Wp,Ws,Rp,Rs,s); %求得滤波器的最小阶数和截止频率 z,p,k=buttap(N); %设计Butterworth滤波器b,a=zp2tf(z,p,k); %将零点极点增益形式

22、转换为传递 函数形式Wo=sqrt(Wp*Ws); %中心 频率bt,at=lp2hp(b,a,Wo); %频率转换w=linspace(1,5000,1000)*2*pi; %设置绘制频率响应的频率点H=freqs(bt,at,w); %计算给定频率点的复数频率响应magH=abs(H);%-figure(1)subplot(2,1,1),plot(w/2/pi,20*log10(magH); %绘制幅频响应xlabel(频率 /Hz);ylabel(振幅/dB );grid on;title(幅频图);subplot(2,1,2),plot(w/2/pi,unwrap(angle(H);%

23、绘制相频响应xlabel(频率 /Hz);ylabel(相位 /o);grid on;title(相频图)%-figure(2)dt=1/10000; %模拟信号采样间隔f1=100;f2=1500;f3=2900;%输入信号的三个频率成分t=0:dt:0.04; %给定模拟时间段x=sin(2*pi*f1*t)+0.5*cos(2*pi*f2*t)+0.5*sin(2*pi*f3*t); %输入信号h=tf(bt,at); %滤波器在MATLAB系统中的表示y,t1=lsim(h,x,t); %模拟输出subplot(2,1,1),plot(t,x),title(输入信号 ) %绘出输入信号

24、subplot(2,1,2),plot(t1,y) %绘制输出信号title(输出信号 ),xlabel(时间/s)%-figure(3)Fy=abs(fft(y);w=(-(length(Fy)-1)/2:(length(Fy)-1)/2)/(length(Fy)*dt);plot(w(length(Fy)-1)/2+1:(length(Fy)-1),Fy(1:(length(Fy)-1)/2)*2/length(Fy) title(输出信号 频谱 )xlabel(频率 /Hz) 0 500 1000 1500 2000 2500 3000 3500 4000 4500 5000-800-6

25、00-400-2000位 位 /Hz位位/dB位 位 位0 500 1000 1500 2000 2500 3000 3500 4000 4500 5000-10-505位 位 /Hz位位/o位 位 位0 0.005 0.01 0.015 0.02 0.025 0.03 0.035 0.04-2-1012 位 位 位 位0 0.005 0.01 0.015 0.02 0.025 0.03 0.035 0.04-0.500.5 位 位 位 位位 位 /s0 500 1000 1500 2000 2500 3000 3500 4000 4500 500000.10.20.30.4 位 位 位 位

26、位 位位 位 /Hz四、切比雪夫I型模拟低通滤波器clc;clear;clf;Wp=800*2*pi;Ws=1200*2*pi;%通带阻带截止频率Rp=3;Rs=20;%通带波纹和阻带衰减N,Wc=cheb1ord(Wp,Ws,Rp,Rs,s) %计算ChebyshevI滤波器的最小阶数和截止频率z,p,k=cheb1ap(N,Rp); %设计Chebyshev I型滤波器b,a=zp2tf(z,p,k); %将零点极点增益形式转换为传递 函数形式bt,at=lp2lp(b,a,Wc);w=linspace(1,5000,1000)*2*pi; %设置绘制频率响应的频率点H=freqs(bt,

27、at,w); %计算给定频率点的复数频率响应magH=abs(H);%-figure(1)subplot(2,1,1),plot(w/2/pi,20*log10(magH); %绘制幅频响应xlabel(频率 /Hz);ylabel(振幅/dB );grid on;title(幅频图);subplot(2,1,2),plot(w/2/pi,unwrap(angle(H);%绘制相频响应xlabel(频率 /Hz);ylabel(相位 /o);grid on;title(相频图)%-figure(2)dt=1/10000; %模拟信号采样间隔f1=100;f2=1500;f3=2900;%输入信

28、号的三个频率成分t=0:dt:0.04; %给定模拟时间段x=sin(2*pi*f1*t)+0.5*cos(2*pi*f2*t)+0.5*sin(2*pi*f3*t); %输入信号h=tf(bt,at); %滤波器在MATLAB系统中的表示y,t1=lsim(h,x,t); %模拟输出subplot(2,1,1),plot(t,x),title(输入信号 ) %绘出输入信号subplot(2,1,2),plot(t1,y) %绘制输出信号title(输出信号 ),xlabel(时间/s)%-figure(3)Fy=abs(fft(y);w=(-(length(Fy)-1)/2:(length(

29、Fy)-1)/2)/(length(Fy)*dt);plot(w(length(Fy)-1)/2+1:(length(Fy)-1),Fy(1:(length(Fy)-1)/2)*2/length(Fy) title(输出信号 频谱 )xlabel(频率 /Hz) 0 500 1000 1500 2000 2500 3000 3500 4000 4500 5000-100-500位 位 /Hz位位/dB位 位 位0 500 1000 1500 2000 2500 3000 3500 4000 4500 5000-8-6-4-20位 位 /Hz位位/o位 位 位0 0.005 0.01 0.015

30、 0.02 0.025 0.03 0.035 0.04-2-1012 位 位 位 位0 0.005 0.01 0.015 0.02 0.025 0.03 0.035 0.04-1-0.500.51 位 位 位 位位 位 /s0 500 1000 1500 2000 2500 3000 3500 4000 4500 500000.20.40.60.8 位 位 位 位 位 位位 位 /Hz五、切比雪夫I型模拟带通滤波器clc;clear;clf;Wp=1200,1800*2*pi;Ws=1000,2000*2*pi;%通带阻带截止频率Rp=3;Rs=20;%通带波纹和阻带衰减N,Wc=cheb1o

31、rd(Wp,Ws,Rp,Rs,s) %计算ChebyshevI滤波器的最小阶数和截止频率z,p,k=cheb1ap(N,Rp); %设计Chebyshev I型滤波器b,a=zp2tf(z,p,k); %将零点极点增益形式转换为传递 函数形式Wo=1500*2*pi; %中心 频率Bw=600*2*pi; %频带宽度bt,at=lp2bp(b,a,Wo,Bw); %频率转换w=linspace(1,5000,1000)*2*pi; %设置绘制频率响应的频率点H=freqs(bt,at,w); %计算给定频率点的复数频率响应magH=abs(H);%-figure(1)subplot(2,1,1

32、),plot(w/2/pi,20*log10(magH); %绘制幅频响应xlabel(频率 /Hz);ylabel(振幅/dB );grid on;title(幅频图);subplot(2,1,2),plot(w/2/pi,unwrap(angle(H);%绘制相频响应xlabel(频率 /Hz);ylabel(相位 /o);grid on;title(相频图)%-figure(2)dt=1/10000; %模拟信号采样间隔f1=100;f2=1500;f3=2900;%输入信号的三个频率成分t=0:dt:0.04; %给定模拟时间段x=sin(2*pi*f1*t)+0.5*cos(2*pi

33、*f2*t)+0.5*sin(2*pi*f3*t); %输入信号h=tf(bt,at); %滤波器在MATLAB系统中的表示y,t1=lsim(h,x,t); %模拟输出subplot(2,1,1),plot(t,x),title(输入信号 ) %绘出输入信号subplot(2,1,2),plot(t1,y) %绘制输出信号title(输出信号 ),xlabel(时间/s)%-figure(3)Fy=abs(fft(y);w=(-(length(Fy)-1)/2:(length(Fy)-1)/2)/(length(Fy)*dt);plot(w(length(Fy)-1)/2+1:(length

34、(Fy)-1),Fy(1:(length(Fy)-1)/2)*2/length(Fy) title(输出信号 频谱 )xlabel(频率 /Hz) 0 500 1000 1500 2000 2500 3000 3500 4000 4500 5000-400-300-200-1000位 位 /Hz位位/dB位 位 位0 500 1000 1500 2000 2500 3000 3500 4000 4500 5000-15-10-50位 位 /Hz位位/o位 位 位0 0.005 0.01 0.015 0.02 0.025 0.03 0.035 0.04-2-1012 位 位 位 位0 0.005

35、 0.01 0.015 0.02 0.025 0.03 0.035 0.04-0.500.5 位 位 位 位位 位 /s0 500 1000 1500 2000 2500 3000 3500 4000 4500 500000.10.20.30.4 位 位 位 位 位 位位 位 /Hz六、切比雪夫I型模拟高通滤波器clc;clear;clf;Wp=2600*2*pi;Ws=2300*2*pi;%通带阻带截止频率Rp=3;Rs=25;%通带波纹和阻带衰减N,Wc=cheb1ord(Wp,Ws,Rp,Rs,s) %计算ChebyshevI滤波器的最小阶数和截止频率z,p,k=cheb1ap(N,Rp

36、); %设计Chebyshev I型滤波器b,a=zp2tf(z,p,k); %将零点极点增益形式转换为传递 函数形式Wo=sqrt(Wp*Ws); %中心 频率bt,at=lp2hp(b,a,Wo); %频率转换w=linspace(1,5000,1000)*2*pi; %设置绘制频率响应的频率点H=freqs(bt,at,w); %计算给定频率点的复数频率响应magH=abs(H);%-figure(1)subplot(2,1,1),plot(w/2/pi,20*log10(magH); %绘制幅频响应xlabel(频率 /Hz);ylabel(振幅/dB );grid on;title(

37、幅频图);subplot(2,1,2),plot(w/2/pi,unwrap(angle(H);%绘制相频响应xlabel(频率 /Hz);ylabel(相位 /o);grid on;title(相频图)%-figure(2)dt=1/10000; %模拟信号采样间隔f1=100;f2=1500;f3=2900;%输入信号的三个频率成分t=0:dt:0.04; %给定模拟时间段x=sin(2*pi*f1*t)+0.5*cos(2*pi*f2*t)+0.5*sin(2*pi*f3*t); %输入信号h=tf(bt,at); %滤波器在MATLAB系统中的表示y,t1=lsim(h,x,t); %

38、模拟输出subplot(2,1,1),plot(t,x),title(输入信号 ) %绘出输入信号subplot(2,1,2),plot(t1,y) %绘制输出信号title(输出信号 ),xlabel(时间/s)%-figure(3)Fy=abs(fft(y);w=(-(length(Fy)-1)/2:(length(Fy)-1)/2)/(length(Fy)*dt);plot(w(length(Fy)-1)/2+1:(length(Fy)-1),Fy(1:(length(Fy)-1)/2)*2/length(Fy) title(输出信号 频谱 )xlabel(频率 /Hz) 0 500 1

39、000 1500 2000 2500 3000 3500 4000 4500 5000-400-300-200-1000位 位 /Hz位位/dB位 位 位0 500 1000 1500 2000 2500 3000 3500 4000 4500 5000-10-50位 位 /Hz位位/o位 位 位0 0.005 0.01 0.015 0.02 0.025 0.03 0.035 0.04-2-1012 位 位 位 位0 0.005 0.01 0.015 0.02 0.025 0.03 0.035 0.04-0.4-0.200.20.4 位 位 位 位位 位 /s0 500 1000 1500 2

40、000 2500 3000 3500 4000 4500 500000.10.20.30.4 位 位 位 位 位 位位 位 /Hz七、切比雪夫II型模拟低通滤波器clc;clear;clf;Wp=800*2*pi;Ws=1200*2*pi;%通带阻带截止频率Rp=1;Rs=150;%通带波纹和阻带衰减N,Wc=cheb2ord(Wp,Ws,Rp,Rs,s) %计算ChebyshevII滤波器的最小阶数和截止频率z,p,k=cheb2ap(N,Rs); %设计Chebyshev II型滤波器b,a=zp2tf(z,p,k); %将零点极点增益形式转换为传递 函数形式bt,at=lp2lp(b,a

41、,Wc);w=linspace(1,5000,1000)*2*pi; %设置绘制频率响应的频率点H=freqs(bt,at,w); %计算给定频率点的复数频率响应magH=abs(H);%-figure(1)subplot(2,1,1),plot(w/2/pi,20*log10(magH); %绘制幅频响应xlabel(频率 /Hz);ylabel(振幅/dB );grid on;title(幅频图);subplot(2,1,2),plot(w/2/pi,unwrap(angle(H);%绘制相频响应xlabel(频率 /Hz);ylabel(相位 /o);grid on;title(相频图)

42、%-figure(2)dt=1/10000; %模拟信号采样间隔f1=100;f2=1500;f3=2900;%输入信号的三个频率成分t=0:dt:0.04; %给定模拟时间段x=sin(2*pi*f1*t)+0.5*cos(2*pi*f2*t)+0.5*sin(2*pi*f3*t); %输入信号h=tf(bt,at); %滤波器在MATLAB系统中的表示y,t1=lsim(h,x,t); %模拟输出subplot(2,1,1),plot(t,x),title(输入信号 ) %绘出输入信号subplot(2,1,2),plot(t1,y) %绘制输出信号title(输出信号 ),xlabel(

43、时间/s)%-figure(3)Fy=abs(fft(y);w=(-(length(Fy)-1)/2:(length(Fy)-1)/2)/(length(Fy)*dt);plot(w(length(Fy)-1)/2+1:(length(Fy)-1),Fy(1:(length(Fy)-1)/2)*2/length(Fy) title(输出信号 频谱 )xlabel(频率 /Hz) 0 500 1000 1500 2000 2500 3000 3500 4000 4500 5000-200-1000100位 位 /Hz位位/dB位 位 位0 500 1000 1500 2000 2500 3000 3500 4000 4500 5000-30-20-100位 位 /Hz位位/o位 位 位0 0.005 0.01 0.015 0.02 0.025 0.03 0