数字信号处理实验50题 清华出版社.docx
《数字信号处理实验50题 清华出版社.docx》由会员分享，可在线阅读，更多相关《数字信号处理实验50题 清华出版社.docx（36页珍藏版）》请在冰点文库上搜索。
数字信号处理实验50题清华出版社
1
clc;
b=[1,1];a=[1,-0.5];
subplot（3,1,1）;zplane（b,a）;title（'因果系统零极图'）;
n=0:
50;
x=3*cos（pi*n/3）;
y=filter（b,a,x）;
subplot（3,1,2）;stem（n,x,'.'）;title（'输入x的波形'）;
subplot（3,1,3）;stem（n,y,'.'）;title（'输出y的波形'）;
1-2
b=[1,1,1];a=[1,0.5,-0.25];
subplot（3,1,1）;zplane（b,a）;title（'因果系统零极图'）;
n=0:
50;
x=3*cos（pi*n/3）;
y=filter（b,a,x）;
subplot（3,1,2）;stem（n,x,'.'）;title（'输入x的波形'）;
subplot（3,1,3）;stem（n,y,'.'）;title（'输出y的波形'）;
2clear;
clc;
b=[0,1];a=[1,-1,-1];
x=impseq（0,-5,50）;n=-5:
50;
h=filter（b,a,x）;
stem（n,h,'.'）;
title（'单位脉冲响应'）
sum（abs（h））
3
b=[2];
a=[1-0.8-0.5];
subplot（4,1,1）;
zplane（b,a）;
title（'系统的零极图'）;
[H,w]=freqz（b,a,100,'whole'）;
magH=abs（H）;
phaH=angle（H）;
subplot（4,1,2）;
plot（w/pi,magH）;
title（'系统的幅频响应'）;
subplot（4,1,3）;
plot（w/pi,phaH/pi）;
title（'系统的相频响应'）;
n=0:
100;
x=impseq（0,0,100）;
h=filter（b,a,x）;
subplot（4,1,4）;
stem（n,h,'.'）;
title（'系统的冲激响应'）;
4
b=[11];
a=[1-0.90.81];
[H,w]=freqz（b,a,400,'whole'）;
magH=abs（H）;
phaH=angle（H）;
subplot（4,1,1）;
plot（w/pi,magH）;
title（'系统的幅频响应'）;
subplot（4,1,2）;
plot（w/pi,phaH/pi）;
title（'系统的相频响应'）;
n=0:
200;
x=sin（pi*n/3）+5*cos（pi*n）;
y=filter（b,a,x）;
subplot（4,1,3）;
plot（n,x）;
title（'输入信号X'）;
subplot（4,1,4）;
plot（n,y）;
title（'输出信号Y'）;grid;
5
x11=[1111];
n=0:
5;
x12=cos（pi*n/4）;
y11=circonvt（x11,x12,8）
y12=conv（x11,x12）
y13=[y11（1:
1:
8）,zeros（1,1）]
e1=y13-y12
x21=[1-11-1];
x22=[10-10];
y21=circonvt（x21,x22,5）
y22=conv（x21,x22）
y23=[y21（1:
1:
5）,zeros（1,2）]
e2=y23-y22
n=0:
15;
x31=cos（2*pi*n/32）;
x32=sin（2*pi*n/32）;
y31=circonvt（x31,x32,32）
y32=conv（x31,x32）
y33=[y32（1:
1:
31）,zeros（1,1）]
e3=y31-y33
n=0:
9;
x41=（0.8）.^n;
x42=（-0.8）.^n;
y41=circonvt（x41,x42,15）
y42=conv（x41,x42）
y43=[y41（1:
1:
15）,zeros（1,4）]
e4=y43-y42
6
clear;
clc;
x1=[2112];
x2=[1-1-11];
n=[0:
8-1];
y11=circonvt（x1,x2,4）
y12=circonvt（x1,x2,7）
y13=circonvt（x1,x2,8）
y2=conv（x1,x2）
%N最小值7
7--1
x=[2,2,2,2,2,2,2,2];w=[0:
1:
500]*2*pi/500;
[H]=freqz（x,1,w）;
magH=abs（H）;phaH=angle（H）;
subplot（2,2,1）;plot（w/pi,magH）;grid
xlabel（''）;ylabel（'|x|'）;title（'DTFT的幅度'）
subplot（2,2,2）;plot（w/pi,phaH/pi*180）;grid
xlabel（'以pi为单位的频率'）;ylabel（'度'）;title（'DTFT的相角'）
N=8;w1=2*pi/N;k=0:
N-1;
X=dft（x,N）;
magX=abs（X）,phaX=angle（X）*180/pi
subplot（2,2,3）;plot（w*N/（2*pi）,magH,'--'）;
axis（[-0.1,8.1,0,20]）;holdon
stem（k,magX）;
ylabel（'|x（k）|'）;title（'DFT的幅度:
N=8'）;text（4.3,-1,'k'）
holdoff
subplot（2,2,4）;plot（w*N/（2*pi）,phaH*180/pi,'--'）;
axis（[-0.1,8.1,-200,200]）;holdon
stem（k,phaX）;
ylabel（'度'）;title（'DFT的相角:
N=8'）;text（4.3,-200,'k'）
7--2
clear;
clc;
x=[2,2,2,2,2,2,2,2];w=[0:
1:
500]*2*pi/500;
[H]=freqz（x,1,w）;
magH=abs（H）;phaH=angle（H）;
subplot（2,2,1）;plot（w/pi,magH）;grid
xlabel（''）;ylabel（'|x|'）;title（'DTFT的幅度'）
subplot（2,2,2）;plot（w/pi,phaH/pi*180）;grid
xlabel（'以pi为单位的频率'）;ylabel（'度'）;title（'DTFT的相角'）
N=16;w1=2*pi/N;k=0:
N-1;
X=fft（x,N）;
magX=abs（X）,phaX=angle（X）*180/pi
subplot（2,2,3）;plot（w*N/（2*pi）,magH,'--'）;
axis（[-0.1,16.1,0,20]）;
holdon
stem（k,magX）;
ylabel（'|x（k）|'）;title（'DFT的幅度:
N=16'）;text（4.3,-1,'k'）
holdoff
subplot（2,2,4）;plot（w*N/（2*pi）,phaH*180/pi,'--'）;
axis（[-0.1,16.1,-200,200]）;
holdon
stem（k,phaX,'.'）;
ylabel（'度'）;title（'DFT的相角:
N=16'）;text（4.3,-250,'k'）
8--1
N=12;w1=2*pi/N;k=0:
N-1;
x=[1,2,3,4,5,6,6,5,4,3,2,1];
X=dft（x,N）;
magX=abs（X）,phaX=angle（X）*180/pi
subplot（2,1,1）;
axis（[-0.1,12.1,0,50]）;
holdon
stem（k,magX）;
ylabel（'|x（k）|'）;title（'DFT的幅度:
N=12'）;holdoff
subplot（2,1,2）;
axis（[-0.1,12.1,-400,400]）;holdon
stem（k,phaX）;
ylabel（'度'）;title（'DFT的相角:
N=12'）;
8--2
x=[1,2,3,4,5,6,6,5,4,3,2,1];w=[0:
1:
500]*2*pi/500;
[H]=freqz（x,1,w）;
magH=abs（H）;phaH=angle（H）;
subplot（2,2,1）;plot（w/pi,magH）;grid
axis（[0,2,0,50]）;
xlabel（''）;ylabel（'|x|'）;title（'DTFT的幅度'）
subplot（2,2,2）;plot（w/pi,phaH/pi*180）;grid
axis（[0,2,-400,400]）;
xlabel（'以pi为单位的频率'）;ylabel（'度'）;title（'DTFT的相角'）
N=12;w1=2*pi/N;k=0:
N-1;
x=[1,2,3,4,5,6,6,5,4,3,2,1];
X=dft（x,N）;
magX=abs（X）,phaX=angle（X）*180/pi
subplot（2,2,3）;plot（w*N/（2*pi）,magH,'--'）;
axis（[-0.1,12.1,0,50]）;
holdon
stem（k,magX）;
ylabel（'|x（k）|'）;title（'DFT的幅度:
N=12'）;
holdoff
subplot（2,2,4）;plot（w*N/（2*pi）,phaH*180/pi,'--'）;
axis（[-0.1,12.1,-400,400]）;
holdon
stem（k,phaX）;
ylabel（'度'）;title（'DFT的相角:
N=12'）;
9
clear;clc;
N1=40;
n=0:
1:
N1-1;
t=0.01*n;
x=2*sin（4*pi*t）+5*cos（16*pi*t）;
x1=fft（x）;
magx1=abs（x1）;
w=2*pi/N1*n;
subplot（3,1,1）;
plot（（w*100）/（2*pi）,magx1）;title（'DFT幅度'）;
axis（[0,25,0,200]）;
N2=60;
n=0:
1:
N2-1;
t=0.01*n;
x=2*sin（4*pi*t）+5*cos（16*pi*t）;
x2=fft（x）;
magx2=abs（x2）;
w=2*pi/N2*n;
subplot（3,1,2）;
plot（（w*100）/（2*pi）,magx2）;title（'DFT幅度'）;
axis（[0,25,0,200]）;
N3=128;
n=0:
1:
N3-1;
t=0.01*n;
x=2*sin（4*pi*t）+5*cos（16*pi*t）;
x3=fft（x）;
magx3=abs（x3）;
w=2*pi/N3*n;
subplot（3,1,3）;
plot（（w*100）/（2*pi）,magx3）;title（'DFT幅度'）;
axis（[0,25,0,400]）;
10
clear;
clc;
N=128;
n=0:
1:
N-1;
t=0.01*n;
x=2*sin（4*pi*t）+5*cos（16*pi*t）;
y=x+0.8*randn（1,length（t））;
x1=fft（x）;
magx1=abs（x1）;
w=2*pi/N*n;
subplot（2,1,1）;plot（（w*100）/（2*pi）,magx1）;title（'DFT幅度'）;
axis（[0,40,0,400]）;
y1=fft（y）;
magy1=abs（y1）;
w=2*pi/N*n;
subplot（2,1,2）;plot（（w*100）/（2*pi）,magy1）;title（'被噪声污染后DFT幅度'）
axis（[0,100,0,400]）;
11
clear;clc;
N=512;n=0:
N-1;t=0.01*n;
x=sin（2*pi*5*t）+sin（2*pi*15*t）+sin（2*pi*30*t）;
X=fft（x,N）;
magx=abs（X）;
k=[0:
1:
N-1];w=2*pi/N*k;
plot（k/N*100,magx）;
title（'FFTN=512'）
xlabel（'频率（单位：
Hz）'）;
ylabel（'|X|'）;grid
axis（[0,100,0,300]）
12
clearall
clc
N1=128;
n1=0:
N1-1;
t1=0.01*n1;
x1=0.5*sin（2*pi*15*t1）+2*sin（2*pi*40*t1）;
k1=0:
1:
127;
w1=2*pi/N1*k1;
X1=fft（x1）;
magX1=abs（X1）;
subplot（2,1,1）;
plot（（w1*100）/（2*pi）,magX1）;
axis（[0,50,0,150]）;
title（'DFTN=128'）;
xlabel（'频率（单位：
pi）'）;
ylabel（'X（k）'）;grid;
N2=1024;
n2=0:
N2-1;
t2=0.01*n2;
x2=0.5*sin（2*pi*15*t2）+2*sin（2*pi*40*t2）;
k2=0:
1:
1023;
w2=2*pi/N2*k2;
X2=fft（x2）;
magX2=abs（X2）;
subplot（2,1,2）;
plot（（w2*100）/（2*pi）,magX2）;
axis（[0,50,0,900]）;
title（'DFTN=1024'）;
xlabel（'频率（单位：
pi）'）;
ylabel（'X（k）'）;grid;
13
t=0:
0.001:
1;
x=sin（2*pi*60*t）+sin（2*pi*200*t）;
subplot（2,1,1）;stem（t,x,'.'）;
title（'signialx（n）'）;grid;
y=x+1.5*randn（1,length（t））;
Y=fft（y,1024）;
p=Y.*conj（Y）/1024;
N=1:
1024;n=N/1000*1024;
subplot（2,1,2）;plot（n,p）;
axis（[0,600,0,280]）;
title（'signialy（n）'）;grid;
xlabel（'频率（单位：
Hz）'）;
ylabel（'p'）;grid
14
n=[0:
1:
9];
x=cos（0.48*pi*n）+cos（0.52*pi*n）;
X=fft（x）;
magx=abs（X（1:
1:
10））;
k=0:
1:
9;
w=2*pi/10*k;
subplot（3,1,1）;
stem（w/pi,magx）;
title（'N=10点DFT幅度'）;
xlabel（'频率（单位:
pi）'）;
axis（[0,1,0,10]）;
n=[0:
1:
9];
y=cos（0.48*pi*n）+cos（0.52*pi*n）;
n1=[0:
1:
99];
x=[y（1:
1:
10）zeros（1,90）];
x1=fft（x）;
magx1=abs（x1（1:
1:
50））;
k1=0:
1:
49;
w1=2*pi/100*k1;
subplot（3,1,2）;
stem（w1/pi,magx1）;
title（'补零到一百点DFT幅度'）;
xlabel（'频率（单位:
pi）'）;
axis（[0,1,0,10]）;
n=[0:
1:
99];
x=cos（0.48*pi*n）+cos（0.52*pi*n）;
X=fft（x）;
magx=abs（X（1:
1:
50））;
k=0:
1:
49;
w=2*pi/100*k;
subplot（3,1,3）;
stem（w/pi,magx）;
title（'N=100点DFT幅度'）;
xlabel（'频率（单位:
pi）'）;axis（[0,1,0,60]）;
15
n=0:
10;
x=10*（0.8.^n）;
x1=fft（x）;
k=0:
10;
y1=x1.*（exp（8*j*pi*k/11））;
y=ifft（y1）;
subplot（2,2,1）;
stem（n,x）;
title（'原序列x（n）'）;
xlabel（'n'）;
axis（[0,10,0,12]）;
subplot（2,2,2）;
stem（n,y）;
title（'移位序列y（n）'）;
axis（[0,10,0,12]）;
n=0:
10;
y=10*（0.8.^n）;
x=[y（1:
1:
11）zeros（1,4）];
n1=0:
14;
subplot（2,2,3）;
stem（n1,x）;
title（'15点序列x（n）'）;
xlabel（'n'）;
axis（[0,14,0,12]）;
x1=fft（x）;
k=0:
14;
y1=x1.*（exp（8*j*pi*k/15））;
y=ifft（y1）;
subplot（2,2,4）;
stem（n1,y）;
title（'15点移位序列y（n）'）;axis（[0,14,0,12]）;
16
clc;
N=31;
n=[0:
N];
x=n.*（stepseq（0,0,N）-stepseq（16,0,N））;
y=stepseq（0,0,N）-stepseq（8,0,N）;
X=fft（x）;
Y=fft（y）;
Z=X.*Y;
z=ifft（Z）;
subplot（3,1,1）;
stem（n,z）;
title（'线性卷积'）;
axis（[0,25,0,100]）;
N1=15;
n1=[0:
31];
x1=n1.*（stepseq（0,0,31）-stepseq（16,0,31））;
y1=stepseq（0,0,N1）-stepseq（8,0,N1）;
X1=fft（x1,16）;
Y1=fft（y1）;
Z1=X1.*Y1;
z1=ifft（Z1）;
subplot（3,1,2）;
n1=[0:
15];
stem（n1,z1）;
title（'16点圆周卷积'）;
axis（[0,20,0,100]）;
N=31;
n=[0:
N];
x=n.*（stepseq（0,0,N）-stepseq（16,0,N））;
y=stepseq（0,0,31）-stepseq（8,0,31）;
X=fft（x）;
Y=fft（y）;
Z=X.*Y;
z=ifft（Z）;
subplot（3,1,3）;
stem（n,z）;
title（'32点圆周卷积'）;axis（[0,25,0,100]）;
17
Rp=0.5;
T=0.001;
ws=200*2*pi*T;
ws1=（2/T）*tan（ws/2）;
[b,a]=cheby1（9,Rp,ws1,'high','s'）;
[bz,az]=bilinear（b,a,1/T）;
[db,mag,pha,grd,w]=freqz_m（bz,az）;
subplot（2,1,1）;
plot（w/pi,db）;
grid;
axis（[0,1,-400,100]）;
title（'系统的幅频响应'）;
subplot（2,1,2）;
plot（w/pi,pha）;
title（'系统的相频响应'）;
18
Wn=2*pi*100;
fs=1000;
[b,a]=butter（6,Wn,'s'）;
[bz,az]=impinvar（b,a,fs）;
[db,mag,pha,grd,w]=freqz_m（bz,az）;
subplot（2,2,1）;
plot（w/pi,db）;
title（'系统的幅频响应'）;
axis（[0,1,-50,5]）;
subplot（2,2,2）;
plot（w/pi,pha）;
title（'系统的相频响应'）;
%Filter
x=[-4,-2,0,-4,-6,-4,-2,-4,-6,-6,-4,-4,-6,-6,...
-2,6,12,8,0,-16,-38,-60,-84,-90,-66,-32,...
-4,2,-4,8,12,12,10,6,6,6,4,0,0,0,0,0,-2,...
-4,0,0,0,-2,-2,0,0,-2,-2,-2,-2,0];
y=filter（bz,az,x）;
N=56;
n=0:
N-1;
subplot（2,2,3）;
plot（n,x）;
title（'输入波形'）;
subplot（2,2,4）;
plot（n,y）;
title（'输出波形'）;
19
%最高f=30Hz，可取fs=100Hz，即t=0.01n
%s（n）=sin（0.1*pi*n）+sin（0.3*pi*n）+sin（0.6*pi*n）;
%s（n）的样本取301点
%注意：
这不是双线性变换法，是完全设计法，不过，效果一样。
参考P
Rp=0.1;
Rs=40;
wp1=0.2*pi;
wp2=0.4*pi;
wn=[wp1,wp2]/pi;
n=4;
[b,a]=ellip（n,Rp,Rs,wn）;%默认时表带通
[db,mag,pha,grd,w]=freqz_m（b,a）;
subplot（3,1,1）;
plot（w/pi,db）;
axis（[0,1,-100,5]）;
n=0:
300;
s=sin（0.1*pi*n）+sin（0.3*pi*n）+sin（0.6*pi*n）;
subplot（312）;
plot（n,s）
holdon
y=filter（b,a,s）;
subplot（313）;
plot（n,y）
20
%可抽象成一低通或带阻滤波器。
抽象成低通来设计
%抽样频率取fs=1000Hz；
%验证看指标
fp=100;
fs=130;
Rp=2;
Rs=50;
T=0.001;
wp=2*pi*fp*T;
ws=2*pi*fs*T;
wp1=（2/T）*tan（wp/2）;
ws1=（2/T）*tan（ws/2）;
[n,wn]=cheb1ord（wp1,ws1,Rp,Rs,'s'）;
[b,a]=cheby1（n,Rp,wn,'low','s'）;
[bz,az]=bilinear（b,a,1/T）;
[db,mag,pha,grd,w]=freqz_m（bz,az）;
gridon
plot（w/pi,db）;
axis（[0,1,-80,5]）;
21
b=[1,1];a=[1,5,6];T=1;
[bz,az]=impinvar（b,a,1/T）
[bz1,az1]=bilinear（b,a,1/T）
22
Rp=2;Rs=30;T=0.001;
wp1=2*pi*
升级会员 