西电数字信号处理大作业

第二章

2.25 已知线性时不变系统的差分方程为

y(n)=−1

2

y(n−1)+x(n)+2x(n−2)

若系统的输入序列x(x)={1,2,3,4,2,1}编写利用递推法计算系统零状态响应的MATLAB程序，并计算出结果。

代码及运行结果：

>> A=[1,-0.5];

>> B=[1,0,2];

>> n=0:5;

>> xn=[1,2,3,4,2,1];

>> zx=[0,0,0];zy=0;

>> zi=filtic(B,A,zy,zx);

>> yn=filter(B,A,xn,zi);

>> figure(1)

>> stem(n,yn,'.');

>> grid on;

2.28图所示系统是由四个子系统T1、T2、T3和T4组成的，分别用单位脉冲响应或差分方程描述为

T1：h 1(n )={1,12,14,1

8,1

16,1

32,n =0,1,2,3,4,5

0, 其他

T2: h 2(n )={

1,1,1,1,1,1,n =0,1,2,3,4,5

0, 其他

T3：y 3(n )=1

4x (n )+1

2x (n −1)+1

4x(n −2)

T4:y (n )=0.9y (n −1)−0.81y (n −2)+v (n )+v(n −1)

编写计算整个系统的单位脉冲响应h(n)，0≤n ≤99的MATLAB 程序，并计算结果。

代码及结果如下：

>> a=0.25;b=0.5;c=0.25; >> ys=0;

>> xn=[1,zeros(1,99)]; >> B=[a,b,c]; >> A=1;

>> xi=filtic(B,A,ys);

>> yn1=filter(B,A,xn,xi);

>> h1=[1,1/2,1/4,1/8,1/16,1/32]; >> h2=[1,1,1,1,1,1]; >> h3=conv(h1,h2);

>> h31=[h3,zeros(1,89)];

>> yn2=yn1+h31;

>> D=[1,1];C=[1,-0.9,0.81]; >> xi2=filtic(D,C,yn2,xi); >> xi2=filtic(D,C,ys); >> yn=filter(D,C,yn2,xi); >> n=0:99; >> figure(1)

>> stem(n,yn,'.'); >> title('单位脉冲响应'); >> xlabel('n');ylabel('yn');

2.30 利用MATLAB画出受高斯噪声干扰的正弦信号的波形，表示为

x(n)=10sin(0.02πn)+v(n), 0≪n≪100=N

其中v(n)是均值为零、方差为1的高斯噪声。

代码及结果如下：

>> N=100;

>> n=0:N;

>> xn=10*sin(0.02*pi*n);

>> R=randn(1,N+1);

>> x=xn+R;

>> figure(2);

>>plot(n,x,'.'),title('受高斯噪声干扰的正弦信号'),xlabel('n'),ylabel('x');

第三章

3.47 利用Matlab 工具箱函数zplane(b,a),画出下列Z 变换的零极点分布图，并给出所有可能的收敛域及对应序列的特性（左边序列，右边序列，双边序列）。

（1）12

181533325644162)(2

342341-+-+++++=z z z z z z z z z X

>> b=[2,16,44,56,32]; >> a=[3,3,-15,18,-12]; >> zplane(b,a)

（2）65

610204

.874.2698.1768.84)(2342342+++--+--=z z z z z z z z z X

>> b=[4,-8.68,-17.98,26.74,-8.04]; >> a=[1,-2,10,6,65]; >> zplane(b,a)

3.53 利用Matlab 语言，画出下列无限长脉冲响应系统的幅频响应特性曲线和相频响应特性曲线，并指出系统的类型。

（1）)

61.088.01)(4.01()1()(2

112

111-----+---=z z z z z z H >> syms z;

>> ps=z^-1*(1-z^-1)^2; >> ps1=expand(ps) ps1 =

1/z - 2/z^2 + 1/z^3 >> syms z;

>> ps = (1 - 0.4*z^ - 1)*(1 - 0.88*z^ - 1 + 0.61*z^ - 2); >> ps1=expand(ps) ps1 =

481/(500*z^2) - 32/(25*z) - 61/(250*z^3) + 1 >> a=[1,-32/25,481/500,-61/250]; >> b=[0,1,-2,1];

>> [H,w]=freqz(b,a,'whole');

>> subplot(2,1,1),plot(w/pi,abs(H));

>> xlabel('\omega/\pi');ylabel('|H(e^j^\omega)|') >> subplot(2,1,2),plot(w/pi,angle(H)/pi);

>> xlabel('\omega/\pi');ylabel('\phi(\omega)/\pi')

(2) )

7957.04461.1.1)(683.01()0166.11)(1(0534.0)(2112

2112------+-+-+=z z z z z z z H

>> syms z;

>> ps=0.0534*(1+z^-1)*(1-1.0166*z^-1+z^-2)^2; >> ps1=expand(ps) ps1 =

6676839363/(125000000000*z^2) - 689661/(12500000*z) + 6676839363/(125000000000*z^3) - 689661/(12500000*z^4) + 267/(5000*z^5) + 267/5000