实验方法信号处理Matlab源代码.docVIP

实验方法及数字信号分析处理 课程大作业  第一次作业： 题目：x=sin(2*pi*5*t)+sin(2*pi*50*t)+sin(2*pi*200*t)+0.1*rand(1,m); ①画原信号的曲线； ②从上述信号中滤出5Hz信号，再画曲线；（分别用卷积和递推低通滤波器） ③从上述信号中滤出50Hz信号，再画曲线；（分别用卷积和递推带通滤波器） Matlab源代码： dt=0.002; %dt1/2*fmax=0.0025,满足采样定理 t=0:dt:1; %采样时间序列 m=length(t); %求取时间采样点数 x=sin(2*pi*5*t)+sin(2*pi*50*t)+sin(2*pi*200*t)+0.1*rand(1,m); %原始信号 %figure Subplot(321);plot(t,x); %画出原始信号图形 title(原始信号) %卷积低通滤波器设计 F=10; %截止频率 N=50; %根据经验取N for i=-N:N f(1,N+1+i)=sin(2*pi*5*i*dt)/(pi*i);%根据设计公式求fi end f(1,N+1)=F*2*dt; %求取f与x的卷积 M=size(f,2);NN=size(x,2); %求取矩阵列长 X=zeros(1,M+NN-1);A=rot90(f); B=[zeros(1,M-1),x,zeros(1,M-1)]; for i=1:M+NN-1 for j=1:M temp=A(j)*B(j+i-1); X(i)=X(i)+temp; end end %figure Subplot(322) ;plot(X) %画出低滤波后频率为5Hz信号图 title(卷积低通滤波器) %卷积带通滤波器设计 F2=100;F1=10; %取上下截止频率 for i=-N:N ff(1,N+1+i)=2*sin(pi*(F2-F1)*i*dt)*cos(pi*(F2+F1)*i*dt)/(pi*i); %根据带通滤波器系数公式求ff end ff(1,N+1)=2*(F2-F1)*dt; Y=conv(ff,x);%直接利用卷积公式求ff与信号x的卷积 %M=size(ff,2);NN=size(x,2); %利用上述源程序求卷积 %X=zeros(1,M+NN-1);A=rot90(ff); %B=[zeros(1,M-1),x,zeros(1,M-1)]; %for i=1:M+NN-1 %for j=1:M %temp=A(j)*B(j+i-1); %X(i)=X(i)+temp; %end %end %figure Subplot(323) ;plot(Y) %画出经带通滤波器滤波后的频率为50Hz信号图形 title(卷积带通滤波器) %递推低通滤波器设计 w=tan(dt*pi*F); %固有频率计算 f0=w*w/(1+2^0.5*w+w^2); %根据递推低通滤波器系数设计公式计算各系数 f1=2*w*w/(1+2^0.5*w+w^2); f2=w*w/(1+2^0.5*w+w^2); g1=-2*(1-w*w)/(1+2^0.5*w+w^2); g2=(1-2^0.5*w+w^2)/(1+2^0.5*w+w^2); %对y1，y2初始化 y(1)=f0*x(1); y(2)=f0*x(2)+f1*x(1)-g1*y(1); for k=3:length(x) %循环滤波 y(k)=f0*x(k)+f1*x(k-1)+f2*x(k-2)-g1*y(k-1)-g2*y(k-2); end %figure Subplot(324); plot(t,y) %画出递推低通滤波器滤波后频率为5Hz信号图 title(递推低通滤波器) [nn,wn]=buttord(0.001,0.2 , 0.3,50); %直接利用MatLAB中ButterWorth公式滤波 [B,A]=butter(nn,wn); yy=filter(B,A,x); %figure Subplot(325); plot(t,yy) title(ButterWort