数字信号处理实验报告.docx
- 文档编号:23980325
- 上传时间:2023-05-23
- 格式:DOCX
- 页数:24
- 大小:150.21KB
数字信号处理实验报告.docx
《数字信号处理实验报告.docx》由会员分享,可在线阅读,更多相关《数字信号处理实验报告.docx(24页珍藏版)》请在冰豆网上搜索。
数字信号处理实验报告
华中科技大学武昌分校
实验(课程设计)报告
实验课程名称:
数字信号处理实验
系(部):
信息科学与技术系
专业班:
通信工程0702
姓名:
汪斯琪
学号:
20071181088
实验老师:
吴莉
时间2009年12月30日——2010年1月8日
目录
实验一信号、系统及系统响应1
1.实验目的1
2.实验原理1
3.实验内容2
冲击响应和阶跃响应2
实验题目2
仿真程序及解释3
调试结果及分析3
卷积定理4
实验题目4
仿真程序及解释4
调试结果及分析5
实验二时域采样与频域采样6
1.实验目的6
2.实验原理6
3.实验内容6
时域采样理论的验证。
6
仿真程序及解释7
调试结果及分析8
频域采样理论的验证9
仿真程序及解释9
调试结果及分析10
4.思考题11
实验三用双线性变换法设计IIR数字滤波器12
1.实验目的12
2.实验内容12
3.实验步骤12
仿真程序及解释13
运行结果及分析15
4.思考题16
实验心得17
实验一信号、系统及系统响应
1.实验目的
(1)加深对离散线性移不变(LSI)系统基本理论的理解,明确差分方程与系统函数之
间的关系。
(2)初步了解用MATLAB语言进行离散时间系统研究的基本方法。
(3)掌握求解离散时间系统冲激响应和阶跃响应程序的编写方法,了解常用子函数。
(4)通过实验进一步理解卷积定理,了解卷积的过程。
(5)了解MATLAB中有关卷积的子函数及其应用方法。
2.实验原理
(1)离散LSI系统的响应与激励
由离散时间系统的时域和频域分析方法可知,一个线性移不变离散系统可以用线性常系数差分方程表示:
系统函数H(z)反映了系统响应与激励的关系。
一旦上式中的bm和ak的数据确定了,则系统的性质也就确定了。
其中特别注意:
a0必须进行归一化处理,
即a0=1。
对于复杂信号激励下的线性系统,可以将激励信号在时域中分解为单位脉冲序列或单位阶跃序列,把这些单元激励信号分别加于系统求其响应,然后把这些响应叠加,即可得到复杂信号加于系统的零状态响应。
因此,求解系统的冲激响应和阶跃响应尤为重要。
由图1-1可以看出一个离散LSI系统响应与激励的关系。
同时,图1-1显示了系统时域分析方法和z变换域分析法的关系。
如果已知系统的冲激响应h(n),则对它进行z变换即可求得系统函数H(z);反之,知道了系统函数H(z),对其进行z逆变换,即可求得系统的冲激响应h(n)。
y(n)=x(n)*h(n)
x(n)
h(n)
H(z)
Y(z)=X(z)H(z)
X(z)
图1-1离散LSI系统响应与激励的关系
(2)离散LSI系统的线性卷积
由理论学习我们已知,对于线性移不变离散系统,任意的输入信号x(n)可以用
及其位移的线性组合来表示,即
当输入为时,系统的输出y(n)=h(n),由系统的线性移不变性质可以得到系统对x(n)的响应y(n)为
称为离散系统的线性卷积,简记为
也就是说,如果已知系统的冲激响应,将输入信号与系统的冲激响应进行卷积运算,即可求得系统的响应。
3.实验内容
●冲击响应和阶跃响应
实验题目
(1)已知一个因果系统的差分方程为
满足初始条件y(-1)=0,x(-1)=0,求系统的冲激响应和阶跃响应。
编写仿真程序,并调试得到结果,进行分析。
仿真程序及解释
a=[1,0,1/3,0];%输入函数的归一化增益,数组同时除以6,保证a的系数为1
b=[1/6,1/2,1/2,1/6];%输入函数的归一化增益
N=32;%作图取的点数为32
n=0:
N-1;%步长为1的N个等间隔采样点,n取0~N-1
hn=impz(b,a,n);%调用单位冲激响应函数
gn=dstep(b,a,n);%调用单位阶跃响应函数
subplot(1,2,1),stem(n,hn,'k');%视图以一行两列的形式显示,当前选定的窗体区域的序号为1,并定义横纵坐标及图表属性,画序列图
title('系统的单位冲激响应');%定义视图一的标题
ylabel('h(n)');xlabel('n');%以n为横坐标,h(n)为纵坐标显示视图
subplot(1,2,2),stem(n,gn,'k');%视图以一行两列的形式显示,当前选定的窗体区域的序号为2,并定义横纵坐标及图表属性,画序列图
title('系统的单位阶跃响应');%定义视图二的标题
ylabel('g(n)');xlabel('n');%以n为横坐标,h(n)为纵坐标显示视图二
调试结果及分析
✓系统的冲激响应是在输入序列为冲激序列的情况下所求得的输出序列。
此实验用递推法来编程求出系统的单位冲激响应。
✓差分方程是用来描述时域离散系统系统输出与输入之间关系的。
✓如果其序列点数N不是32,而是64的话,它图中相应的横坐标也会变为只有64个点,但是相较于上图会紧凑很多;反之,N为16,会稀疏很多。
●卷积定理
实验题目
(2)已知两个信号序列:
f1=0.8n(0 f2=u(n)(0 求两个序列的卷积和。 编写仿真程序,并调试得到结果,进行分析。 仿真程序及解释 nf1=0: 20;%定义步长为1的21点等间隔采样点,n为0~20 f1=0.8.^nf1;%序列f1的函数表达式f1=0.8n subplot(2,2,1);stem(nf1,f1,'filled');%视图以两行两列的形式显示,当前选定的窗体区域的序号为1,并定义横纵坐标及图表属性,画序列图 title('f1(n)');%视图一标题为f1(n) nf2=0: 10;%定义步长为1的11点等间隔采样点 lf2=length(nf2);%求序列nf2的长度 f2=ones(1,lf2);%序列nf2的函数表达式 subplot(2,2,2);stem(nf2,f2,'filled');%视图以两行两列的形式显示,当前选定的窗体区域的序号为2,并定义横纵坐标及图表属性,画序列图 title('f2(n)');%视图二标题为f2(n) y=conv(f1,f2);%求函数f1与f2的卷积和 subplot(2,1,2);stem(y,'filled');title('y(n)');%以两行一列的形式显示最后结果的图示,并定义横纵坐标及图表属性,视图标题为y(n) 调试结果及分析 ✓线性时不变系统的输出等于输入序列和该系统的单位脉冲响应的卷积。 y(n)=∑0.8m*u(n-m);卷积长度N=21+11-1=31,所以y(n)长度>=31. 实验二时域采样与频域采样 1.实验目的 (1)掌握模拟信号采样前后频谱的变化,以及如何选择采样频率才能使采样后的信号不丢失信息 (2)掌握频率域采样会引起时域周期化的概念,以及频率域采样定理及其对频域采样点数选择的指导作用。 (3)会用MATLAB语言进行时域抽样与信号重建的方法,以及频域抽样与恢复时程序的编写方法。 2.实验原理 了解时域采样定理的要点,理解理想采样信号 和模拟信号 之间的关系, 了解频域采样定理的要点,掌握这两个采样理论的结论: “时域采样频谱周期延拓,频域采样时域信号周期延拓”。 3.实验内容 (1)时域采样理论的验证。 给定模拟信号, 式中A=444.128, =50 π, =50 πrad/s (2)用DFT(FFT)求该模拟信号的幅频特性,选取三种采样频率,以验证时域采样理论。 (3)编写实验程序,计算 、 和 的幅度特性,并绘图显示。 观察分析频谱混叠失真。 仿真程序及解释 %时域采样理论验证 Tp=64/1000; Fs=1000;T=1/Fs; M=Tp*Fs;n=0: M-1;%产生M=64采样序列x(n) A=444.128;alph=pi*50*2^0.5;omega=pi*50*2^0.5;%给模拟信x(t)的系数赋值A=444.12, =50 π, =50 πrad/s xnt=A*exp(-alph*n*T).*sin(omega*n*T);%进行nt采样,xnt=xa(nt) Xk=T*fft(xnt,M);%M点FFT[xnt],即64个点 yn=xnt;%赋值 subplot(3,2,1);stem(n,yn,’k’);boxon;title('(a)Fs=1000Hz');%第一个图是xa(t)的Fs=1000Hz的时域抽样 k=0: M-1;fk=k/Tp; subplot(3,2,2);plot(fk,abs(Xk));title('(a)T*FT[xa(nT)],Fs=1000Hz');%第二个图是xa(t)的幅频特性曲线 xlabel('f(Hz)');ylabel('幅度');axis([0,Fs,0,1.2*max(abs(Xk))]);%设置坐标比例 Fs=300;T=1/Fs; M=ceil(Tp*Fs);n=0: M-1;%取整 xnt=A*exp(-alph*n*T).*sin(omega*n*T);%xnt=xa(nt) Xk=T*fft(xnt,M);%M点FFT[xnt)] yn=xnt; subplot(3,2,3);stem(n,yn,’k’);boxon;title('(b)Fs=300Hz');%第三个图是xa(t)的Fs=300Hz的时域抽样 k=0: M-1;fk=k/Tp; subplot(3,2,4);plot(fk,abs(Xk));title('(b)T*FT[xa(nT),Fs=300Hz');%第四个图是xa(t)的幅频特性曲线 xlabel('f(Hz)');ylabel('幅度');axis([0,Fs,0,1.2*max(abs(Xk))]);%设置坐标比例 Fs=200;T=1/Fs; M=ceil(Tp*Fs);n=0: M-1; xnt=A*exp(-alph*n*T).*sin(omega*n*T);%xnt=xa(nt) Xk=T*fft(xnt,M);%M点FFT[xnt] yn=xnt;subplot(3,2,5);stem(n,yn,’k’);boxon;title('(c)Fs=200Hz');%第五个图是xa(t)的Fs=200Hz的时域抽样 k=0: M-1;fk=k/Tp; subplot(3,2,6);plot(fk,abs(Xk));title('(c)T*FT[xa(nT),Fs=200Hz');%第六个图是xa(t)的幅频特性曲线 xlabel('f(Hz)');ylabel('幅度');axis([0,Fs,0,1.2*max(abs(Xk))]);%设置坐标比例 调试结果及分析 式中A=444.128, =50 π, =50 πrad/s ✓^xa(t)=∑xa(t)δ(t-nT),pδ(t)=∑δ(t-nT),其中n=-∞~∞,δ(t)是单位冲激信号,只有当t=nT时,才有非零值,故^xa(t)=∑xa(nT)δ(t-nT) 在傅里叶变换中,两信号在时域相乘的傅里叶变换等于两个信号分别的傅里叶变换的卷积。 ˆXa(jΩ)=Xa(jΩ)*Pδ(jΩ)/2π==∑xa(nT)e^(-jΩnt) 在数值上xa(nT)=x(n)故ˆXa(jΩ)=∑x(n)e^(-jωn),n=-∞~∞,ω=ΩT ✓理想采样信号的频谱是原模拟信号的频谱以Ωs为周期,进行周期性延拓而产生的;采样频率Ωs必须大于等于模拟信号最高频率两倍以上,才能使采样信号的频谱不产生频谱混叠。 在采样频率不变的情况下,采样点数越多,频谱分辨率越高。 如图当采样频率为1000Hz时,频谱混叠很小;当采样频率为300Hz时,频谱混叠严重;当采样频率为200Hz时,频谱混叠很严重。 如图所示,频率在0~500,0~150,0~100时周期再现。 (4)频域采样理论的验证 给定信号如下: (5)编写程序分别对频谱函数 在区间 上等间隔采样32 和16点,得到 ,再分别对 进行32点和16点IFFT,得到 。 仿真程序及解释 %频域采样理论验证 M=27;N=32;n=0: M; %产生M长三角波序列x(n) xa=0: floor(M/2);xb=ceil(M/2)-1: -1: 0;xn=[xa,xb]; Xk=fft(xn,1024);%1024点FFT[x(n)],用于近似序列x(n)的TF X32k=fft(xn,32);%32点FFT[x(n)] x32n=ifft(X32k);%32点IFFT[X32(k)]得到x32(n) X16k=X32k(1: 2: N);%隔点抽取X32k得到X16(K) x16n=ifft(X16k,N/2);%16点IFFT[X16(k)]得到x16(n) subplot(3,2,2);stem(n,xn,'.');boxon title('(b)三角波序列x(n)'); xlabel('n');ylabel('x(n)');axis([0,32,0,20]);%横纵轴的长度 k=0: 1023;%取样点0~1023 wk=2*k/1024;%赋值 subplot(3,2,1);plot(wk,abs(Xk));title('(a)FT[x(n)]');%视图显示为三行两列,当前选定的窗体区域序号为1 xlabel('\omega/\pi');ylabel('|X(e^j^\omega)|'); %x轴y轴内容 axis([0,1,0,200]);%横纵轴的取值范围 k=0: N/2-1; subplot(3,2,3);stem(k,abs(X16k),'.');boxon title('(c)16点频域采样'); xlabel('k');ylabel('|X_1_6(k)|');axis([0,8,0,200]) n1=0: N/2-1; subplot(3,2,4);stem(n1,x16n,'.');boxon title('(d)16点IDFT[X_1_6(k)]');xlabel('n');ylabel('x_1_6(n)');axis([0,32,0,20]) k=0: N-1; subplot(3,2,5);stem(k,abs(X32k),'.');boxon title('(e)32点频域采样');xlabel('k');ylabel('|X_3_2(k)|');axis([0,16,0,200]) n1=0: N-1; subplot(3,2,6);stem(n1,x32n,'.');boxon title('(f)32点IDFT[X_3_2(k)]');xlabel('n');ylabel('x_3_2(n)');axis([0,32,0,20]) 调试结果及分析 ✓频域采样时域信号周期延拓,由于实序列的DFT满足共轭对称性,因此频域图仅画在[0,π]上的幅频特性波形。 ✓频率采样点数N必须大于等于时域离散信号的长度M(即N≥M),才能使时域不产生混叠,则N点IDFT[XN(k)]得到的序列xN(n)就是原序列x(n),即xN(n)=x(n)。 如果N>M,xN(n)比原序列尾部多N-M个零点;如果N 如图所示,当采样点数N=16 (6)分别画出 、 的幅度谱,并绘图显示x(n)、 的波形,进行对比和分析,验证总结频域采样理论。 4.思考题: 如果序列x(n)的长度为M,希望得到其频谱 在 上的N点等间隔采样,当N 答: 先对原序列x(n)以N为周期进行周期延拓后取主值区序列,x(n)=[∑x(n+iN)]R(n),再计算N点DFT,则得到N点频域采样: XN(k)=DFT[x(n)]= |w=2πk/Nk=0,1…N-1 实验三用双线性变换法设计IIR数字滤波器 1.实验目的 (1)熟悉用双线性变换法设计IIR数字滤波器的原理与方法。 (2)掌握数字滤波器的计算机仿真方法。 (3)通过观察对实际心电图信号的滤波作用,获得数字滤波的感性知识。 2.实验内容 (1)用双线性变换法设计一个巴特沃斯低通IIR数字滤波器。 设计指标参数为: 在通带内频率低于0.2π时,最大衰减小于1dB;在阻带内[0.3π,π]频率区间上,最小衰减大于15dB。 (2)以0.02π为采样间隔,打印出数字滤波器在频率区间[0,π/2]上的幅频响应特性曲线。 (3)用所设计的滤波器对实际心电图信号采样序列(在本实验后面给出)进行仿真滤波处理,并分别打印出滤波前后的心电图信号波形图,观察总结滤波作用与效果。 3.实验步骤 (1)复习有关巴特沃斯模拟滤波器设计和用双线性变换法设计IIR数字滤波器的内容,按照例6.4.2,用双线性变换法设计数字滤波器系统函数H(z)。 例6.4.2中已求出满足本实验要求的数字滤波器系统函数: (2-1) (2-2) A=0.09036 B1=1.2686,C1=-0.7051 B2=1.0106,C2=-0.3583 B3=0.9044,C3=-0.2155 由(2-1)式和(2-2)式可见,滤波器H(z)由三个二阶滤波器H1(z),H2(z)和H3(z)级联组成,如图2-1所示。 图2-1滤波器H(z)的组成 (2)编写滤波器仿真程序,计算H(z)对心电图信号采样序列x(n)的响应序列y(n)。 设为第k级二阶滤波器的输出序列,为输入序列,如图2-1所示。 由(2-2)式可得到差分方程: (2-3) 当k=1时,所以的总响应序列 可以用顺序迭代算法得到。 即依次对k=1,2,3,求解差分方程(2-3), 最后得到仿真程序就是实现上述求解差分方程和顺序迭代算法的通用程序。 也可以直接调用MATLABfilter函数实现仿真。 (3)在通用计算机上运行仿真滤波程序,并调用通用绘图子程序,完成实验内容 (2)和(3)。 仿真程序及解释 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]; k=1;%给k赋初值,使滤波从第一级开始 closeall;%关闭所有参数 figure (1);%图标 subplot(2,2,1);%将x(n)序列用图形显示出来 n=0: 55;%n值从0取到55 stem(n,x,'.');%画出序列函数离散点为n,幅度响应为x axis([055-10050]);%纵坐标从0到50;横坐标从-100到50 holdon;%保存当前图像不被覆盖 n=0: 60;%n值从0取到60 m=zeros(61);%提供零基准线 plot(n,m);%画平面图 xlabel('n');%x坐标显示n ylabel('x(n)');%y坐标显示x(n) title('心电图信号采样序列x(n)'); B=[0.090362*0.090360.09036];%输入序列的增益系数 A=[1.2686-0.7051];%第一级滤波器输出序列的增益系数 A1=[1.0106-0.3583];%第二级滤波器输出序列的增益系数 A2=[0.9044-0.2155];%第三级滤波器输出序列的增益系数 while(k<4) y=filter(B,A,x); x=y; if(k==2) A=A1; elseif(k==3) A=A2; else; end k=k+1; end;%完成三级滤波 subplot(2,2,3);%绘图 n=0: 55;%n值从0取到55 stem(n,y,'.');%画出序列函数离散点为n,幅度响应为y axis([055-155]);%纵坐标从0到55;横坐标从-15到5 holdon;%保存当前图像不被覆盖 n=0: 60;%n值从0取到60 m=zeros(61);%提供零基准线 plot(n,m);%画平面图 xlabel('n');%x坐标显示n ylabel('y(n)');%y坐标显示y(n) title('三级滤波后的心电图信号');%显示“三级滤波后的心电图信号”的图题 A=[0.09036,0.1872,0.09036];%输入序列的增益系数 B1=[1,-1.2686,0.7051];%第一级滤波器输出序列的增益系数 B2=[1,-1.0106,0.3583];%第二级滤波器输出序列的增益系数 B3=[1,-0.9044,0.2155];%第三级滤波器输出序列的增益系数 [H1,w]=freqz(A,B1,100);%求第一级滤波器的频率响应 [H2,w]=freqz(A,B2,100);%求第二级滤波器的频率响应 [H3,w]=freqz(A,B3,100);%求第三级滤波器的频率响应 H4=H1.*(H2);%第三级的频率响应等于前两级频率响应的卷积 H=H4.*(H3);%一级二级滤波器矩阵相乘结果与三级滤波器矩阵相乘 mag=abs(H);%求幅度变化 db=20*log10((mag+eps)/max(mag));%通带最大衰减幅度 subplot(2,2,2)%绘图 plot(w/pi,db);%画平面图 axis([00.5-5010]);%纵坐标从0到0.5;横坐标从-50到10 title('滤波器的幅频响应曲线');%显示“滤波器的幅频响应曲线”的图题 运行结果及分析 ✓如图所示,图中第1和第3张图中红色部分是由程序中n=0: 60和m=zeros(61)两句实现; 图中第2张图中0.2表示通带截止频率; 程序中的while循环部分完成三级滤波的作用; 程序中A=[0.09036,0.1872,0.09036]中的0.09036是式子开三次方所得。 ✓双线性变换法的特点: 采用非线性频率压缩的方法,将整个模拟频率轴压缩到±π/T之间,再变换到Z平面内,使之不产生频谱混叠现象,但是也不能保真的模仿模拟滤波器的频响曲线形状。 比较图中第1和第3张图
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 数字信号 处理 实验 报告
![提示](https://static.bdocx.com/images/bang_tan.gif)