数字信号处理实验一、二、三 、实验四2

1、用MATLAB计算产生x(n)=sin(pi*n/5),n=-5:5,并绘出波形图。

解：n=-5:5;

y=sin(pi*n/5); plot(n,y ,'b*-'); xlabel('-5-5'); ylabel('-1-1');

10.80.60.40.2-1-10-0.2-0.4-0.6-0.8-1-5-4-3-2-10-5-512345

实验二、离散时间系统的时域分析

2、编写教材P35(1.12t)的y(n)运算的MATLAB程序。

（1）、x1=0.5^n*u(n),x2=R5(n)的线性卷积； 解：n=0:Nh-1;

h1=(0.5).^n; %产生h1（n） h2=h1;

nx=0: Nx-1;x1=ones(1,Nx); %产生x1(n) x2= zeros(1,Nx+m);

for k=m+1:m+Nx; %产生x2(n)= x1(n-m) x2(k)= x1(k-m);

end%产生x2(n)

y1=conv(x1,h1); %计算y1(n)=x1(n)*h1(n) y2=conv(x2,h2); %计算y2(n)=x2(n)*h2(n)

1 / 10

subplot(3,2,1) stem(nx,x1,'.')

axis([0 30 0 1.2]),title('x1(n)')

subplot(3,2,3) %绘图 stem(n, h1)

axis([0 30 0 1.2]),title('h1(n)') subplot(3,2,5)

stem(y1),title('y1(n)') subplot(3,2,2) stem(x2,'.')

axis([0 40 0 1.2]),title('x2(n)') subplot(3,2,4) stem(n,h2,'.')

axis([0 40 0 1.2]),title('h2(n)') subplot(3,2,6)

stem(y2),title('y2(n)')

x1(n)10.50010h1(n)10.50210010y1(n)210203010.5001020y2(n)3040203010.5001020h2(n)3040x2(n)01020300102030

（2）、x1=2^R4(n),x2(n)=δ(n)-δ(n-2)的卷积。

解：Nh=20;Nx=3; m=1;%设定Nh=20;Nx=3和位移值m=1 n=0:Nh-1;

h1=2.^n; %产生h1（n） h2=h1;

nx=0: Nx-1;x1=[1 0 -1]; %产生x1(n) x2= zeros(1,Nx+m);

for k=m+1:m+Nx; %产生x2(n)= x1(n-m)

2 / 10

x2(k)= x1(k-m);

end%产生x2(n)

y1=conv(x1,h1); %计算y1(n)=x1(n)*h1(n) y2=conv(x2,h2); %计算y2(n)=x2(n)*h2(n) subplot(3,2,1) stem(nx,x1,'.')

axis([0 30 -1 16]),title('x1(n)')

subplot(3,2,3) %绘图 stem(n, h1)

axis([0 30 -1 16]),title('h1(n)') subplot(3,2,5)

stem(y1),title('y1(n)') subplot(3,2,2) stem(x2,'.')

axis([0 40 -1 16]),title('x2(n)') subplot(3,2,4) stem(n,h2,'.')

axis([0 40 -1 16]),title('h2(n)') subplot(3,2,6)

stem(y2),title('y2(n)')

x1(n)1510500151050010-1x 106x2(n)15105010h1(n)203015105001020h2(n)304010y1(n)203010-10x 1061020y2(n)304001020300102030

3 / 10

实验三、 离散系统的Z域分析

*练习

①例

3：

k=256;

num=[0.8 -0.44 0.36 0.02]; den=[1 0.7 -0.45 -0.6]; w=0:pi/k:pi;

h=freqz(num,den,w); subplot(2,2,1);

plot(w/pi,real(h));grid

title('实部')

xlabel('\\omega/\\pi');ylabel('幅度') subplot(2,2,2);

plot(w/pi,imag(h));grid

title('虚部')

xlabel('\\omega/\\pi');ylabel('Amplitude') subplot(2,2,3);

plot(w/pi,abs(h));grid

title('幅度谱')

xlabel('\\omega/\\pi');ylabel('幅值') subplot(2,2,4);

plot(w/pi,angle(h));grid

title('相位谱')

xlabel('\\omega/\\pi');ylabel('弧度')

4 / 10

实部66虚部20-200.5/幅度谱1Amplitude4420-200.5/相位谱1幅度8632幅值42000.5/1弧度10-100.5/1

例5：

b=[1 2 3 ];

a=[1 -0.5 -0.005 0.3]; figure(1); zplane(a,b); num=[0 1 2 1];

den=[1 -0.5 -0.005 0.3]; h=impz(num,den); figure(2); stem(h);

xlabel('h');

ylabel('幅度');

title('Impulse Respone')

5 / 10

1.510.5Imaginary Part0-0.5-1-1.5-1.5-1-0.50Real Part0.511.5

Impulse Respone2.521.51幅度0.50-0.5-1051015h20253035

②、用MATLAB

1z1H(z),110.2z的频率特性图。 编程，画出

6 / 10

k=256;

num=[1 -1]; den=[1 -0.2]; w=0:pi/k:pi;

h=freqz(num,den,w); subplot(2,2,1);

plot(w/pi,real(h));grid

title('实部')

xlabel('\\omega/\\pi');ylabel('幅度') subplot(2,2,2);

plot(w/pi,imag(h));grid

title('虚部')

xlabel('\\omega/\\pi');ylabel('Amplitude') subplot(2,2,3);

plot(w/pi,abs(h));grid

title('幅度谱')

xlabel('\\omega/\\pi');ylabel('幅值') subplot(2,2,4);

plot(w/pi,angle(h));grid

title('相位谱')

xlabel('\\omega/\\pi');ylabel('弧度')

实部21虚部10.5000.5/幅度谱1Amplitude1.5幅度0.5000.5/相位谱121.521.5幅值10.5000.5/1弧度10.5000.5/1

7 / 10

3、根据MATLAB求解

X(z)1(10.2z1)(10.5z1)2(10.6z1)2z0.6的结果，写出序列x(n)的表达式。

解： b=1;

a=poly([-0.2 0.5 0.5 -0.6 -0.6]); [r,p,k]=residuez(b,a) r =

0.1826 0.4463 0.2032 0.1476 0.0204 p =

-0.6000 -0.6000 0.5000 0.5000 -0.2000 k = []

X(Z)所以

0.18260.44630.20320.14760.0204 10.6000Z-110.6000Z-110.5000Z-110.5000Z-110.2000Z-1x(n)0.1826(0.6)n0.4463(0.6)n0.2032(0.5)n0.1476(0.5)n0.0204(0.2)n实验四、离散傅里叶变换及其快速算法

设x(n)R4(n),X(ej)FT[x(n)]。分别计算X(ej)在频率区间 [0,2]上的16点和32点等间隔采样，并绘制X(ej)采样的幅频特性图和相频特性图。 解：xn=[1 1 1 1]; %输入时域序列向量 Xk16=fft(xn,16); %计算xn的16点DFT Xk32=fft(xn,32); %计算xn的32点DFT %以下为绘图部分

k=0:15;wk=2*k/16; %产生16点DFT对应的采样点频率(关于π归一化值) subplot(2,2,1);

8 / 10

stem(wk,abs(Xk16),'.'); %绘制16点DFT的幅频特性图 title('(a)16点DFT的幅频特性图'); xlabel('ω/π'); ylabel('幅度') subplot(2,2,3);

stem(wk,angle(Xk16),'.'); %绘制16点DFT的相频特性图 line([0,2],[0,0]);

title('(b)16点DFT的相频特性图') xlabel('ω/π'); ylabel('相位');

axis([0,2,-3.5,3.5]);

k=0:31;wk=2*k/32; %产生32点DFT对应的采样点频率(关于π归一化值) subplot(2,2,2);stem(wk,abs(Xk32),'.'); %绘制32点DFT的幅频特性图 title('(c)32点DFT的幅频特性图'); xlabel('ω/π'); ylabel('幅度'); subplot(2,2,4);

stem(wk,angle(Xk32),'.'); %绘制32点DFT的相频特性图 line([0,2],[0,0]);

title('(d)32点DFT的相频特性图'); xlabel('ω/π'); ylabel('相位');

axis([0,2,-3.5,3.5])

(a)16点DFT的幅频特性图4343(c)32点DFT的幅频特性图幅度210011.52ω/π(b)16点DFT的相频特性图0.5幅度210011.5ω/π(d)32点DFT的相频特性图0.5222相位0-200.51ω/π1.52相位0-200.51ω/π1.52

9 / 10

2、编写序列x(n)con(n)con(n)的DFT运算程序。

48解：

N=16;

n=0:1:N-1;%时域采样

xn=cos(n*pi/4)+ cos(n*pi/8); k=0:1:N-1; %频域采样 WN=exp(-j*2*pi/N); nk=n'*k;

WNnk=WN.^nk; %点乘方 Xk=xn*WNnk; Subplot(2,1,1) Stem(n,xn)

title('xn的时域序列'); xlabel('n');ylabel('x (n)'); Subplot(2,1,2) Stem(k,abs(Xk));

title('xn的幅频特性'); xlabel('(N=16 wk=2pik/N)k');ylabel('[X (k)]');

xn的时域序列21x (n)0-1-205nxn的幅频特性101015[X (k)]5005(N=16 wk=2pik/N)k1015

10 / 10