信号与系统实验2

信号与系统实验--系统的时域分析实验目的1、熟悉和掌握常用的用于系统时域分析的MATLAB函数。2、掌握连续时间系统零状态响应的求解3、掌握离散时间系统零状态响应的求解4、掌握时间系统单位脉冲响应的求解。实验器材软件、计算机实验原理系统仿真实质上就是对系统模型的求解，对控制系统来说，一般模型可以转化成某个微分方程或养分方程表示，因此在仿真过程中，一般以某种数值算法从初态出发，逐步计算系统的响应，最后绘制出系统的响应曲线，进而可分析系统的性能。控制系统最常用的时域分析方法是当输入信号为单位阶跃和单位冲激函数时，求出系统的输出响应，分别称为单位阶跃响应和单位冲激响应。在MATLAB中提供了求取连续系统的单位阶跃响应函数step，单位冲激响应函数impulse，零输入响应函数initial等等。实验内容及步骤例3-34图3-21所示力学系统中物体位移用y（t）与外力f(t)的关系md2y(t)/dt2+fddy(t)/dt+ksy(t)=f(t)质量m=1kg，ks=100N/m，fd=2Ns/m。ts=0;te=5;dt=0.01;sys=tf([1],[12100]);t=ts:dt:te;f=10*sin(2*pi*t);y=lsim(sys,f,t);plot(t,y);xlabel('Time(sec)')ylabel('y(t)')例3-35在例3-34所述力学系统中，若外力f(t)是强度为10的冲激信号，求物体的位移y(t)%program3_2ImpulsereponseofLTIsystemts=0;te=5;dt=0.01;sys=tf([10],[12100]);t=ts:dt:te;y=impulse(sys,t);plot(t,y);xlabel('Time(sec)')ylabel('y(t)')例3-36受噪声干扰的信号为f[k]=s[k]+d[k]，其中s[k]=(2k)0.9^k是原始信号，d[k]是噪声。已知M点滑动平均系统的输入输出关系为y[k]=≡1/M∑f[k-n]试编程实现用M点滑动平均系统对受干扰的信号去噪%program3_3SignalSmoothingbyMovingAverageFilterR=51;d=rand(1,R)-0.5;k=0:R-1;s=2*k.*(0.9.^k);f=s+d;figure(1);plot(k,d,'r-.',k,s,'b--',k,f,'g-');xlabel('Timeindexk');legend('d[k]','s[k]','f[k]');M=5;b=ones(M,1)/M;a=1;y=filter(b,a,f);figure(2);plot(k,s,'b--',k,y,'r-')xlabel('Timeindexk');legend('s[k]','y[k]');例3-37用impz函数求离散时间系统y[k]+3y[k-1]+2y[k-2]=f[k]的单位脉冲响应h[k]，并与理论值h[k]=-(-1)^k+2(-2)^k,k=0比较。%program3_4Impulseresponseofdiscretesystemk=0:10;a=[132];b=[1];h=impz(b,a,k);subplot(2,1,1)stem(k,h)title('jishizhi');hk=-(-1).^k+2*(-2).^k;subplot(2,1,2)stem(k,hk)title('lilunzhi')例3-38已知序列x[k]={1,2,3,4;k=0,1,2,3},y[k]={1,1,1,1,1;k=0,1,2,3,4},计算x[k]*y[k]并画出卷积结果。%program3_5sequenceconvolutionx=[1,2,3,4];y=[1,1,1,1,1];z=conv(x,y);N=length(z);stem(0:N-1,z);M3-1一系统满足的微分方程为y''(t)+5y'(t)+6y(t)=u(t)-u(t-1)(1)求出该系统的零状态响应yf(t).(2)用lsim求出该系统的零状态响应的数值解。利用（1）所求的结果，比较不同的抽样间隔对数值解精度的影响。（1）t=-6:0.001:6;sys=tf([1],[156]);ft2=heaviside(t)-heaviside(t-1);y1=lsim(sys,ft2,t);plot(t,y1)(2)t=-6:0.1:6;sys=tf([1],[156]);ft2=heaviside(t)-heaviside(t-1);y1=lsim(sys,ft2,t);plot(t,y1)精度越小越精确。M3-2在题M3-2图所示电路中，L=1H，C=1F，R1=1Ω，R2=2Ω，f(t)是输入信号，y（t）是输出响应。（1）建立描述该系统的微分方程（2）用impulse函数求系统的冲激响应；（3）用step函数求系统的阶跃响应。（1）1.5y''(t)+1.5y'(t)+0.5y(t)=f'(t)(2)、（3）t=0:0.1:10;sys=tf([10],[1.5,1.5,0.5])y=impulse(sys,t);subplot(1,2,1)plot(t,y)subplot(1,2,2)stem(t,y)M3-3下列系统分别为一阶、二阶和三阶BW模拟低通滤波器，用impulse函数分别求出各系统的冲激响应，并比较它们的特征。（1）y'(t)+y(t)=f(t)t=0:0.01:10;sys=tf([1],[11]);y=impulse(sys,t);plot(t,y（2）y''(t)+√2y'(t)+y(t)=f(t)t=0:0.01:10;sys=tf([1],[1sqrt(2)1]);y=impulse(sys,t);plot(t,y)(3)y'''(t)+2y''(t)+2y'(t)+y(t)=f(t)t=0:0.01:10;sys=tf([1],[1221]);y=impulse(sys,t);plot(t,y)M3-4下列系统分别为BW模拟低通、高通、带通、带阻滤波器，用impulse函数分别求出各系统的冲激响应，并比较特征。（1）y''(t)+√2y'(t)+y(t)=f(t)(2)y''(t)+√2y'(t)+y(t)=f''(t)(3)y''(t)+y'(t)+y(t)=f'(t)(4)y''(t)+y'(t)+y(t)=f''(t)+f(t)t=0:0.01:10;sys=tf([1],[1sqrt(2)1]);y=impulse(sys,t);subplot(2,2,1)plot(t,y)title('(1)');t=0:0.01:10;sys=tf([100],[1sqrt(2)1]);y=impulse(sys,t);subplot(2,2,2)plot(t,y)title('(2)');t=0:0.01:10;sys=tf([10],[111]);y=impulse(sys,t);subplot(2,2,3)plot(t,y)title('(3)');t=0:0.01:10;sys=tf([101],[111]);y=impulse(sys,t);subplot(2,2,4)plot(t,y)title('(4)');M3-5利用conv函数验证卷积和的交换律、分配律、结合律。交换律：x=[1,2,3,4,5];y=[1,1,1,1,1,1];z=conv(x,y);z1=conv(y,x);N=length(z);M=length(z1);subplot(1,2,1)stem(0:N-1,z)subplot(1,2,2)stem(0:M-1,z1)分配律：t=0:10;x=exp(2*t);y=t;n=x+y;z=conv(n,h);z1=conv(y,h)+conv(x,h);N=length(z);M=length(z1);subplot(1,2,1)stem(0:N-1,z)subplot(1,2,2)stem(0:M-1,z1)结合律：t=0:10;x=exp(2*t);y=t;n=t-1;m=conv(x,y);z=conv(m,n);m1=conv(y,n);z1=conv(x,m1);N=length(z);M=length(z1);subplot(1,2,1)stem(0:N-1,z)subplot(1,2,2)stem(0:M-1,z1)M3-6已知f[k]=+_f[N-1-k],h[k]=+_h[N-1-k],用conv函数计算f[k]*h[k],并总结f[k]*h[k]的对称关系。M3-7两个连续信号的卷积定义为用y(t)=∫f(τ)h(t-τ)dτ为了进行数值计算，需对连续信号进行抽样。记f[k]=f(kΔ)，h(k)=h(kΔ),Δ为进行数值计算的抽样间隔，则……..(1)(2)fork=1:3m=[0.01,0.1,1];t=0:m(k):3;ft=heaviside(t)-heaviside(t-1);ht=conv(ft,ft);y=conv(ft,ht);n=0:length(y)-1;subplot(1,3,k)plot(n,y)endM3-8利用impz函数，计算系统y[k]+0.7[k-1]-0.45y[k-2]-0.6y[k-3]=0.8f[k]-0.44f[k-1]+0.36f[k-2]+0.02f[k-3]的单位脉冲响应，并画出前31点的图。k=0:30;a=[10.7-0.45-0.6];b=[0.8-0.440.360.02];h=impz(b,a,k);stem(k,h)M3-9利用MATLAB的filter函数，求出下列单位脉冲响应（1）y[k]-1.845y[k-1]+0.850586y[k-2]=f[k]（2）y[k]-1.85y[k-1]+0.85y[k-2]=f[k]k=-10:250;b=[1];a=[1,-1.845,0.850586];f=[k==0];h1=filter(b,a,f);subplot(2,1,1),stem(k,h1),axis([-20,300,0,5.5])xlabel('k'),ylabel('h[k]'),title('(1)')c=[1,-1.85,0.85];h2=filter(b,c,f);subplot(2,1,2),stem(k,h2),axis([-20,300,0,10])xlabel('k'),ylabel('h[k]'),title('(2)')M3-10题M3-10图所示力学系统中物体质量为m,弹簧的弹性系数分别为ks1和ks2，物体与地面的摩擦系数为fd,物体m在外力f(t)作用下的位移y(t)……(1)md^2y(t)/dt^2+fd*dy(t)/dt+(ks1+ks2)y(t)=f(t)(2)t=0:0.01:15;m=2;ks1=4;ks2=6;fd=20;b=[1];a=[m,fd,(ks1+ks2)];sys=tf(b,a);y=impulse(sys,t);plot(t,y),xlabel('t'),ylabel('y(t)')实验思考与心得本次实验用到了连续信号的微积分与卷积的求解，通过MATLAB中强大的运算符号，可以很方便算出数值和画出图形来，实现了现实中算术的麻烦。