1、23栗国齐 数字信号处理实验报告指导老师：杨永双,何燕,张彦,李娜院 系：计算机与通信工程学院班 级：通信工程09-1班姓 名：栗国齐学 号：540907040123数字信号处理实验报告实验一 信号、系统及系统响应一、 实验目的：1. 熟悉理想采样的性质，了解信号采样前后的频谱变化，加深对时域采样定理的理解。2. 熟悉离散信号和系统的时域特性。3. 熟悉线性卷积的计算编程方法：利用卷积的方法，观察，分析系统响应的时域特性。4. 掌握序列傅立叶变换的计算机实验方法，利用序列的傅里叶变换对离散信号、系统及系统响应进行频域分析。二、 实验代码： 信号产生子程序：代码：clc;close all;cl
2、ear all;A=444.128;a=50*pi*sqrt(2.0);W0=50*pi*sqrt(2.0);n=-50:1:50;T=1/1000;Xa=A*(exp(a*n*T).*(sin(W0*n*T); %产生理想采样信号序列subplot(3,1,1);stem(n,Xa,r);title(理想采样信号序列);xlabel(n);ylabel(Xa);Xb=(n=0); %产生单位脉冲序列subplot(3,1,2);stem(n,Xb,r);title(单位脉冲序列);xlabel(n);ylabel(Xb);Xc=(n=0)-(n=9); %产生矩形序列subplot(3,1,
3、3);stem(n,Xc,r); title(矩形序列);xlabel(n);ylabel(Xc);运行结果： 分析理想采样信号序列的特性：代码：clc;close all;clear all;A=444.128; a=50*pi*sqrt(2.0);W0=50*pi*sqrt(2.0);n=-50:50;T1=1/1000;Xa=A*(exp(a*n*T1).*(sin(W0*n*T1); % T1=1/1000时得Xa序列stem(n,Xa,r);title(Xa序列);xlabel(n);ylabel(Xa);k=-25:25;X1=Xa*(exp(-j*pi/12.5).(n*k);
4、%X1=FTXapause %暂停，观察波形subplot(3,2,1); stem(k,abs(X1),r); % Xa的幅度谱及相位谱title(Xa的幅度谱);ylabel(幅度);subplot(3,2,2);stem(k,angle(X1),r); title(Xa的相位谱);ylabel(相位);T2=1/300; %T2=1/300时得Xb序列Xb=A*(exp(a*n*T2).*(sin(W0*n*T2);X2=Xb*(exp(-j*pi/12.5).(n*k); %X2=FTXbsubplot(3,2,3);stem(k,abs(X2),r); % Xb的幅度谱及相位谱tit
5、le(Xb的幅度谱);ylabel(幅度);subplot(3,2,4);stem(k,angle(X2),r); title(Xb的相位谱);ylabel(相位);T3=1/200; % T3=1/200时得Xc序列Xc=A*(exp(a*n*T3).*(sin(W0*n*T3);X3=Xc*(exp(-j*pi/12.5).(n*k); %X3=FTXc subplot(3,2,5); stem(k,abs(X3),r); % Xc的幅度谱及相位谱title(Xc的幅度谱);ylabel(幅度);subplot(3,2,6);stem(k,angle(X3),r); title(Xc的相位
6、谱);ylabel(相位);运行结果：图表 1 图表 2 离散信号、系统和系统响应的分析：代码：1) clc;close all;clear all;n=-50:1:50;Xb=(n=0); hb=(n=0)+2.5*(n=1)+2.5*(n=2)+(n=3); Yb=conv(Xb,hb);subplot(3,1,1); stem(n,Xb,r); %产生单位脉冲序列title(单位脉冲序列);xlabel(n);ylabel(Xb);k=-25:25;X=Xb*(exp(-j*pi/12.5).(n*k); %X=FTXbsubplot(3,1,2); stem(k,abs(X),r);
7、%产生Xb的幅度谱及相位谱title(Xb的幅度谱); ylabel(幅度);subplot(3,1,3); stem(k,angle(X),r); title(Xb的相位谱); ylabel(相位);pause %暂停观察波形subplot(3,1,1);stem(n,hb,r); title(hb序列);xlabel(n);ylabel(hb);H=hb*(exp(-j*pi/12.5).(n*k); %H=FThbsubplot(3,1,2);stem(k,abs(H),r); %hb的幅度谱及相位谱title(hb的幅度谱); ylabel(幅度);subplot(3,1,3);ste
8、m(k,angle(H),r); title(hb的相位谱); ylabel(相位);pause %暂停观察波形m=-100:100; %输出信号的定义域的定义subplot(3,1,1);stem(m,Yb,r); %输出信号时域波形title(Yb序列);xlabel(n);ylabel(Yb);Y=Yb*(exp(-j*pi/12.5).(m*k); %Y=FTYb subplot(3,1,2);stem(k,abs(Y),r); %产生Yb的幅度谱及相位谱 title(Yb的幅度谱); ylabel(幅度);subplot(3,1,3);stem(k,angle(Y),r);title
9、(hb的相位谱); ylabel(相位);2) clc;close all;clear all;n=-50:1:50;m=-100:100; Xc=(n=0)-(n=10); ha=Xc;Yb=conv(Xc,ha);k=-25:25;X=Xc*(exp(-j*pi/12.5).(n*k); %X=FTXcH=ha*(exp(-j*pi/12.5).(n*k); %H=FThaY=Yb*(exp(-j*pi/12.5).(m*k); %Y=FTYbfigure(1) %信号时域波形subplot(3,1,1); stem(n,Xc,r); title(单位脉冲序列);xlabel(n);yla
10、bel(Xc);subplot(3,1,2);stem(n,ha,r);title(ha序列);xlabel(n);ylabel(ha);subplot(3,1,3);stem(m,Yb,r); title(Yb序列);xlabel(n);ylabel(Yb);figure(2)subplot(3,1,1); stem(k,abs(X),r); %Xc的幅度谱title(Xc的幅度谱); ylabel(幅度);subplot(3,1,2);stem(k,abs(H),r); %ha的幅度谱title(ha的幅度谱); ylabel(幅度);subplot(3,1,3);stem(k,abs(Y
11、),r); % Yb的幅度谱 title(Yb的幅度谱); ylabel(幅度);3)clc;close all;clear all;n=0:50;m=0:100;T=1;A=1;a=0.4;W0=2.0734;k=-25:25;K=-50:50;ha=(n=0)-(n=9);Xa1=A*(exp(a*n*T).*(sin(W0*n*T);a=0.1;Xa2=A*(exp(a*n*T).*(sin(W0*n*T);WO=1.2516;Xa3=A*(exp(a*n*T).*(sin(W0*n*T);Ya1=conv(Xa1,ha);Ya2=conv(Xa2,ha);Ya3=conv(Xa3,ha
12、);figure(1) %信号时域波形subplot(3,1,1);stem(n,Xa1,r); title(a=0.4，W0=2.0734时Xa序列);xlabel(n);ylabel(Xa);subplot(3,1,2);stem(n,Xa2,r);title(a=0.1，W0=2.0734时Xa序列);xlabel(n);ylabel(Xa);subplot(3,1,3);stem(n,Xa3,r);title(a=0.1，WO=1.2516时Xa序列);xlabel(n);ylabel(Xa);figure(2) %响应时域波形subplot(3,1,1);stem(m,Ya1,r);
13、 title(a=0.4，W0=2.0734时Ya序列);xlabel(n);ylabel(Ya);subplot(3,1,2);stem(m,Ya1,r); title(a=0.4，W0=2.0734时Ya序列);xlabel(n);ylabel(Ya);subplot(3,1,3);stem(m,Ya1,r); title(a=0.4，W0=2.0734时Ya序列);xlabel(n);ylabel(Ya);figure(3) subplot(3,1,1); stem(k,abs(fft(Xa1),r); %Xa的幅度谱title(a=0.4，W0=2.0734时Xa的幅度谱); ylab
14、el(幅度);subplot(3,1,2);stem(k,abs(fft(Xa2),r);title(a=0.1，W0=2.0734时Xa的幅度谱); ylabel(幅度);subplot(3,1,3);stem(k,abs(fft(Xa3),r);title(a=0.1，WO=1.2516时Xa的幅度谱); ylabel(幅度);figure(4)subplot(3,1,1); stem(K,abs(fft(Ya1),r); %Ya的幅度谱title(a=0.4，W0=2.0734时Ya的幅度谱); ylabel(幅度);subplot(3,1,2);stem(K,abs(fft(Ya2),
15、r);title(a=0.1，W0=2.0734时Ya的幅度谱); ylabel(幅度);subplot(3,1,3);stem(K,abs(fft(Ya3),r);title(a=0.1，WO=1.2516时Ya的幅度谱); ylabel(幅度);运行结果：1) 图表 1图表 22) 图表 1图表 23) 图表 1 图表 2 图表 3 图表 4三、 实验思考：1. 在分析理想采样信号序列的特性实验中，利用不痛采样频率所得到的采样信号序列的傅氏变化频谱，数字频谱度量是否相同？它们所对应的模拟频率是否都相同？答：由模拟角频率与数字角频率之间的关系可知，若采样频率不同，则其周期T不同，相应的数字角
16、频率也不相同。而对于同一个信号，其模拟角频率保持不变。实验二 应用FFT对信号进行频谱分析一、 实验目的：1. 在理论学习的基础上，通过本次实验，加深对快速傅立叶变换的理解，熟悉FFT算法及其程序的编写。2. 熟悉应用FFT对典型信号进行频谱分析的方法。3. 了解应用FFT进行信号频谱分析过程中可能出现的问题，以便在实际中正确应用FFT。二、 实验代码： 信号产生子程序：代码：clc;close all;clear all;n=0:15;a=0.1;f=0.0625;p=8;q=2;Xa=exp(-(n-p).2/q); Xb=exp(-a*n).*sin(2*pi*f*n); Xc=(n+1
17、).*(n=0)-(n3)+(8-n).*(n=4)-(n7); Xd=(4-n).*(n=0)-(n3)+(n-3).*(n=4)-(n7); subplot(2,2,1);stem(n,Xa,r);title(高斯序列);xlabel(n);ylabel(Xa);subplot(2,2,2);stem(n,Xb,r);title(衰减正弦序列);xlabel(n);ylabel(Xb);subplot(2,2,3);stem(n,Xc,r); title(三角波序列);xlabel(n);ylabel(Xc);subplot(2,2,4);stem(n,Xd,r); title(反三角序列
18、);xlabel(n);ylabel(Xc);运行结果： 观察高斯序列的时域和频域特性：代码：clc;clear all;close all;n=0:15;p=8;q=2;xa1=exp(-(n-p).2/q);q=4;xa2=exp(-(n-p).2/q);q=6;xa3=exp(-(n-p).2/q);figure(1)subplot(3,1,1);stem(n,xa1,r); % q=2时xa信号波形title(q=2);xlabel(n);ylabel(xa); subplot(3,1,2);stem(n,xa2,r); %q=4时xa信号波形title(q=4);xlabel(n);
19、ylabel(xa);subplot(3,1,3);stem(n,xa3,r); %q=6时xa信号波形title(q=6);xlabel(n);ylabel(xa); figure(2)k=-25:25;X=xa1*(exp(-j*pi/12.5).(n*k); %xa傅立叶变换subplot(3,1,1); stem(k,abs(X),r); %q=2时xa的幅度谱title(xa1的幅度谱); ylabel(幅度);X=xa2*(exp(-j*pi/12.5).(n*k); %q=4时xa傅立叶变换subplot(3,1,2);stem(k,abs(X),r); %q=4时xa的幅度谱t
20、itle(xa2的幅度谱); ylabel(幅度); X=xa3*(exp(-j*pi/12.5).(n*k); subplot(3,1,3);stem(k,abs(X),r); % q=6时xa3的幅度谱title(xa3的幅度谱); ylabel(幅度);运行结果：图表 1图表 2 观察衰减正弦序列的时域和幅域特性：代码：clc;clear all;close all;n=0:15;q=8;k=-25:25;p=8;xa1=exp(-(n-p).2/q);p=13; xa2=exp(-(n-p).2/q);p=14;xa3=exp(-(n-p).2/q);X1=xa1*(exp(-j*pi
21、/12.5).(n*k); % p=8时傅立叶变换X2=xa2*(exp(-j*pi/12.5).(n*k); %p=13时傅立叶变换X3=xa3*(exp(-j*pi/12.5).(n*k); %p=14时傅立叶变换figure(1)subplot(3,1,1);stem(n,xa1,r);title(p=8时xa时域曲线);xlabel(n);ylabel(xa);subplot(3,1,2);stem(n,xa2,r);title(p=13时xa时域曲线);xlabel(n);ylabel(xa); subplot(3,1,3);stem(n,xa3,r);title(p=14时xa时域
22、曲线);xlabel(n);ylabel(xa);figure(2)subplot(3,1,1); stem(k,abs(X1),r); %p=8时Xa1的幅度谱title(p=8时xa的幅度谱); ylabel(幅度);subplot(3,1,2); stem(k,abs(X2),r); %p=13时Xa2的幅度谱title(p=13时xa的幅度谱); ylabel(幅度);subplot(3,1,3);stem(k,abs(X3),r);title(p=14时xa的幅度谱); ylabel(幅度); %p=14时Xa3的幅度谱运行结果：图表 1 图表 2 观察三角波序列和反三角序列的时域和
23、幅域特性：代码：clc;clear all;close all;n=0:15;a=0.1;f=0.0625;xb1=exp(-a*n).*sin(2*pi*f*n); X1=fft(xb1);f=0.4375;xb2=exp(-a*n).*sin(2*pi*f*n);X2=fft(xb2);f=0.5625;xb3=exp(-a*n).*sin(2*pi*f*n);X3=fft(xb3); figure(1)subplot(3,1,1);stem(xb1,r);title(f=0.0625时xb的时域波形);xlabel(n);ylabel(xb);subplot(3,1,2);stem(xb
24、2,r);title(f=0.4375时xb的时域波形);xlabel(n);ylabel(xb);subplot(3,1,3);stem(xb3,r);title(f=0.5625时xb的时域波形);xlabel(n);ylabel(xb);figure(2)subplot(3,1,1);stem(abs(fft(X1),r); %f=0.0625时Xb的幅度谱title(f=0.0625时xb的幅度谱);ylabel(幅度);subplot(3,1,2);stem(abs(fft(X2),r); %f=0.4375时Xb的幅度谱title(f=0.4375时xb的幅度谱);ylabel(幅
25、度);subplot(3,1,3);stem(abs(fft(X3),r); %f=0.5625时Xb的幅度谱title(f=0.5625时xb的幅度谱);ylabel(幅度);运行结果：图表 1图表 2三、 实验思考：1. 对一个有限长序列进行离散傅立叶变换（DFT），等价于讲该序列周期延拓后进行傅立叶级数（DFS）展开。因为DFS也只是取其中一个周期来运算，所以FFT在一定条件下也可以用以分析周期信号序列。如果实正弦信号sin(2fn)，f=0.1，用16点的FFT来做DFS运算，得到的频谱是信号本身的真是谱吗？答：不是，因取的点不是一个完整的周期，所以不是真实谱。实验三 用双线性变换法设
26、计IIR滤波器一、 实验目的：1. 了解两种工程上最常用的变换方法：脉冲响应不变法和双线性变换法。2. 掌握双线性变换法设计IIR滤波器的原理及具体设计方法，熟悉用双线性设计法设计低通、带通和高通IIR数字滤波器的计算机程序。3. 观察用双线性变换法设计的滤波器的频域特性，并与脉冲响应不变法相比较，了解双线性变化法的特点。4. 熟悉用双线性变化法设计数字butterworth和chebyshev滤波器的全过程。5. 了解多项式乘积和多项式乘方运算的计算机编程方法。二、 实验代码： 设计一个Chebyshev高通数字滤波器：代码：clc;clearclose allFS=1;t=1/FS;fp=
27、0.3;fs=0.2;wp=fp/FS*2*pi;ws=fs/FS*2*pi;OmegaP=2*FS*tan(wp/2);OmegaS=2*FS*tan(ws/2);n,Wn=cheb1ord(OmegaP,OmegaS,0.8,20,s);b,a=cheby1(n,0.8,Wn,high,s);freqs(b,a);pausebz,az=bilinear(b,a,FS);freqz(bz,az,512,FS);bz,az运行结果：bz = 0.0262 -0.1047 0.1570 -0.1047 0.0262az = 1.0000 1.5289 1.6537 0.9452 0.2796图表
28、1 图表 2 设计一个数字低通滤波器：代码：clc;clear all;close all;FS=1;Fl=0.2;Fh=0.3;wp=(Fl/FS)*2*pi;ws=(Fh/FS)*2*pi;OmegaP=2*FS*tan(wp/2);OmegaS=2*FS*tan(ws/2);n,Wn=buttord(OmegaP,OmegaS,1,25,s);b,a=butter(n,Wn,s);bz,az=bilinear(b,a,FS);freqz(bz,az,512,FS,whole)运行结果： 设计一个Butterworth带通数字滤波器：代码：clc;clear all;close all;T=1;Fs=1/T;Wp1=tan(0.4/2/Fs*pi);Wph=tan(0.6/2/Fs*pi);Ws1=tan(0.3/2/Fs*pi);Wsh=tan(0.7/2/Fs*pi);Rp=1;As=40;OmigaP=Wph-Wp1;OmigaS=Wsh-Ws1;N,Wn=buttord(OmigaP,OmigaS,Rp,As,s);B,A=butter(N,Wn,s);BT,AT=lp2bp(B,A,sqrt(Wph*Wp1),Wph
