数字信号处理实验指导书.docx

1、数字信号处理实验指导书实验一 离散傅立叶变换一、实验目的1. 在理论学习的基础上，通过本实验，加深对离散傅立叶变换（DFT）的理解。2. 以正弦信号为例，学习和掌握利用离散傅立叶变换（DFT）分析信号频谱的方法。3. 了解应用DFT进行信号频谱分析过程中可能出现的问题，以便在实际中正确应用DFT。4. 学习和熟悉在计算机上作图显示和研究信号的离散频谱二、实验原理与方法在各种信号序列中，有限长序列的数字信号处理占有很重要地位，对有限长序列，我们可以使用离散傅立叶变换（DFT）, 这一变换不但可以很好地反映序列的频谱特性，而且易于用快速算法在计算机上实现，当序列x(n)的长度为N时，它的DFT定义

2、为： X (k)= 反变换为： x(n)= 有限长序列的DFT是其Z变换在单位圆上的等距采样，或者说是序列傅立叶变换的等距采样，因此可以用于序列的谱分析。 DFT要求被分析序列的长度为有限长，对无限长和很长的序列则需要截取。截取的结果使得信号的频谱分析出现频率泄漏。本实验分析的周期性信号，对周期信号，当截取的序列长度是信号周期的整数倍的时候，DFT分析的频率离散点恰好采在周期信号的谐波分量和其它零值频率点上，能正确反映信号的频谱。当截取的序列长度不为信号周期的整数倍时，DFT分析的结果就会有频率泄漏的反映。三、实验内容1 设计说明编制程序用DFT分析正弦信号的频谱。具体包括以下几部分（1）编程

3、产生正弦抽样输入信号 x(t) = sin(t) = sin(2f t)x(n) = x(nT) = sin(2f n T)正弦信号的频率f、取样点数N和取样间隔T通过人机对话的方式由键盘输入，以便调整。（2）编程实现DFT的运算DFT的计算公式为 （3）计算DFT的幅度（4）作图画出DFT的幅度谱线。2实验步骤将编制的程序输入计算机，调试运行正确后，进行以下实验。（1）分别在下列几种情况下，输入正弦信号的频率f、采样点数N和采样周期T，观察和研究DFT的结果，并记录谱线的大体形状。（a） f = 50 Hz，N = 32，T = 0.000625 s（b） f = 50 Hz，N = 32，

4、T = 0.005 s（c） f = 50 Hz，N = 32，T = 0.004 s（d） f = 150 Hz，N = 32，T = 0.000625 s （2）修改程序，将输入信号改变成50 Hz和150 Hz两个正弦信号分量的迭加，并令50 Hz分量的幅度为2，150 Hz分量的幅度为1。即x(n) = 2 sin(250 n T) + sin(2150 n T)取N = 32，分别在T = 0.000625 s、T = 0.0025 s和T = 0.002 s情况下进行DFT分析，观察和记录结果。（3）有余力的同学还可以将输入信号改成50Hz的周期性方波，再对其进行DFT分析。四、C

5、语言参考程序/* 用DFT分析正弦信号频谱的C语言程序 */#include#include#include#include#define PI 3.1415926float x1024,y1024, w1024;void draw(int);void axis(int,int);main() int N,f,n,k; float T,r,i,c;/* 键盘输入f、N、T */ printf(The frequency of the Sine wave f=); scanf(%d,&f); printf(The number of samples N=); scanf(%d,&N); prin

6、tf(The sampling period T=); scanf(%f,&T); c=2*PI/N;/* 产生正弦抽样信号 */ for(n=0;nN;n+) xn=sin(2*PI*f*n*T); /* 计算DFT */ for(k=0;kN;k+) r=i=0.0; for(n=0;nN;n+) r=r+xn*cos(c*n*k); i=i+xn*sin(c*n*k); yk=r; wk=i; /* 计算DFT的幅度 */ for(k=0;kN;k+) wk=sqrt(yk*yk+wk*wk); /* 画频谱图 */ draw(N); getch(); closegraph();void

7、 draw(int N) /* 作图 */ int k,driver,mode,step,x,y; driver=VGA; mode=VGAHI; step=512/N; x=64; y=400; registerbgidriver(driver); initgraph(&driver,&mode,c:tc);/* 绘制坐标及刻度 */ axis(x,y); outtextxy(x+step+546),y+8,k); outtextxy(x+step-8),y-290,X(k); outtextxy(x+step-8),y-234,32); outtextxy(x+step-8),y-121,

8、16); outtextxy(x+step-19),y+8, 0 1 2 3 4 5 6 7 8 9 0 1 2 3 4 5 6 7 8 9 0 1 2 3 4 5 6 7 8 9 0 1);/* 绘制频谱图 */ for(k=0;kN;k+) line(x,y,x,(y-floor(wk*7); x+=step; void axis(int x,int y) /* 绘制横轴与纵轴 */ line(10,y,630,y); line(610,(y-4),630,y); line(610,(y+4),630,y); line(x,(y+20),x,(y-280); line(x-4),(y-26

9、0),64,120); line(x+4),(y-260),64,120);五、实验报告要求（1）给出几种情况下DFT分析的结果，绘制相应的谱线。（2）对每种情况的结果进行分析，从理论上加以说明。（3） 通过本实验，谈谈你对DFT的进一步认识。实验二 快速傅立叶变换一、实验目的1 学习和掌握快速傅立叶变换（FFT）的实现过程和编程技术2 运用FFT分析正弦信号的频谱3 测试FFT的运算时间，比较FFT与DFT的运算速度，获得对 FFT“快速”的感性认识。4 锻炼和提高数字信号处理的程序设计和调试能力。二、实验原理与方法FFT并不是与DFT不同的另一种变换，而是为了减少DFT运算次数的一种快速算

10、法。它是对DFT变换式进行一次次分解，使其成为若干小点数的组合，从而减少运算量。常用的FFT是以2为基数的，其长度N=2M。它的效率高，程序简单，使用非常方便，当要变换的序列长度不等于2的整数次方时，为了使用以2为基数的FFT，可以用末尾补零的方法，使其长度延长至2的整数次方。本实验运用时间抽取基2 FFT，其原理、信号流图和运算过程可参阅课堂笔记、教材和其它教科书。FFT的实现要比DFT复杂，通常采用三个嵌套循环来实现。最外面的循环是分级循环，N=2M点的FFT分为M级计算。中间一层是分组循环，一级内蝶形系数Wk相同的蝶形构成一组。最内层为蝶形计算的循环，一组内不同输入数据的蝶形逐个计算。编

11、程时要注意各级蝶形组之间的间隔、组内蝶形之间间隔、以及系数Wk变化。有不少书中有FFT程序可供参考，但要注意理解和弄懂，不可一味照搬，要尽量自己去编。三、实验内容1设计说明编制时间抽取基2 FFT程序计算前面DFT程序分析过的正弦信号，与DFT计算的结果进行比较，以验证所编程序的正确性。FFT为复数运算，若用实数运算来实现，需要设置两个数组来存放输入输出数据。其中一个用于存放数据的实部，另一个存放数据的虚部。正弦输入信号为实数，则令其虚部为零。同样，碟形运算也要化成实数来进行，分别算出实部和虚部。当然也可以直接用复数数组和语句实现。程序设计的难点和重点在于要合理安排和正确设置碟形运算的分级循环

12、、级内分组循环和组内分碟形循环。要注意乘法系数WNk的变化以及各循环变量的变化和调整。程序中，在进行FFT计算之前，先要将输入数据排列成二进制倒序的形式，这可通过将有关数据相互换位来实现。正弦抽样信号的产生和输出频谱图的绘制与前面的DFT程序相同。2实验步骤（1）将编制的程序输入计算机，调试、运行正确。即输入实验一中的12组数据，看FFT输出谱线的结果是否与前面DFT的结果一致。（2）在不同的FFT长度下（可分别令N=64，128，256，512，1024），运行FFT程序，观察、比较和记录运行时间。由于现在计算机的运算速度很快，计算一次FFT的时间很短（在ms数量级），不便观察和记录。可在程

13、序中添置循环语句，让其重复计算FFT许多次（可重复100次、1000次或更多，视情况而定），使总的运算时间达到数秒或数十秒钟，从而可以通过人工来观察和记录。此时注意，为了防止FFT多次循环计算，使输入输出数据越来越大而发生溢出，可令输入数据全部为零。这样，FFT的运算量保持不变，但输入输出始终为零，不会因多次循环而发生溢出。在观测FFT的运算时间时，作图部分的语句不要运行（一方面，它不应该计入FFT或DFT的运算时间。另一方面，当N较大时，原有程序中的作图语句会出错），可将其注释掉。 测量FFT运行时间的一个更好、更精确的方法是在程序中调用读取计算机时钟的函数。TURBO C中提供多种读取计算

14、机时间的方法和函数，一种比较简单的方法是用biostime()函数，后面给出了相应的参考程序。（3）用上述类似的方法，在不同的长度下（N=64，128，256，512，1024），运行实验一中的DFT程序，记录运行时间，比较它与FFT的速度差距。四、C语言参考程序（1） 用FFT分析正弦信号频谱的程序 #includevoid fft(x,y,n) /* FFT作为子程序 */int n;float x1024,y1024; int i,j,k,l,m,n1,n2; float c,s,e,tr,ti;/* 计算FFT的级数M */ for(j=1,i=1;in;i+) m=i; j=2*j;

15、 if(j=n) break; /* 改变输入数据的顺序 */ n1=n-1; for(j=0,i=0;in1;i+) if(ij) tr=xj; ti=yj; xj=xi; yj=yi; xi=tr; yi=ti; k=n/2; while(k(j+1) j=j-k; k=k/2; j=j+k; /* 三重嵌套循环进行碟形计算 */ n1=1; for(l=1;l=m;l+) /* 分级循环 */ n1=2*n1; n2=n1/2; e=3.14159265359/n2; for(k=0;kn2;k+) /* 分组循环 */ c=cos(k*e); /* 计算WNk */ s=-sin(k*

16、e); for(i=k;in;i+=n1) /* 碟形循环 */ j=i+n2; tr=c*xj-s*yj; /* 计算碟形 */ ti=c*yj+s*xj; xj=xi-tr; yj=yi-ti; xi=xi+tr; yi=yi+ti; /* 主程序部分 */ #include #include #include #include #define PI 3.1415926 float x1024,y1024,w1024; void draw(int); void axis(int,int); main() int N,f,n,k; float T; /* 键盘输入f、N、T */printf

17、(the frequency of the sine wave f=); scanf(%d,&f); printf(the number of samples N=); scanf(%d,&N); printf(the time step of samples T=); scanf(%f,&T); /* 产生正弦输入数据 */ for (n=0;nN;n+) xn=sin(2*PI*f*n*T); yn=0; /* 调用子程序计算FFT */ fft(x,y,N);/* 计算FFT的幅度 */ for(k=0;kN;k+) wk=sqrt(xk*xk+yk*yk); /* 画FFT的幅度谱 *

18、/ draw(N); getch(); closegraph(); /* 作图子程序 */ void draw(int N) int k,driver,mode,step,x,y; driver=VGA; mode=VGAHI; step=512/N; x=64;y=400; registerbgidriver(driver); initgraph(&driver,&mode,c:tc); axis(x,y); outtextxy(x+step+546),y+8,k); outtextxy(x+step-8),y-290,x(k); outtextxy(x+step-8),y-234,32);

19、 outtextxy(x+step-8),y-121,16);outtextxy(x+step-19),y+8,0 1 2 3 4 5 6 7 8 9 0 1 2 3 4 5 6 7 8 9 0 1 2 3 4 5 6 7 8 9 0 1); for(k=0;kN;k+) line(x,y,x,(y-floor(wk*7); x+=step; void axis(int x,int y) line(10,y,630,y); line(610,(y-4),630,y); line(610,(y+4),630,y); line(x,(y+20),x,(y-280); line(x-4),(y-26

20、0),64,120); line(x+4),(y-260),64,120);（2） 测量FFT运行时间的程序 以下给出的是用于测量FFT运算时间的有关程序语句，应将它们插入前面FFT程序中main()函数内适当的位置。unsigned int start,end;double time;int i; /*提示：这两行放在主程序开头数据类型说明部分*/ /*提示：以下部分替代主程序中的fft(x,y,N)*/start=biostime(0,0); /* 读取开始时间 */ for (n=0;n2000;n+) /* 循环计算FFT 2000次 */ fft(x,y,N);for (i=0;iN

21、;i+) /*清除x,y的值，防止下一次的FFT数据溢出 xi=0; yi=0; end=biostime(0,0); /* 读取结束时间 */ time=(end-start)/2000.0*55; /* biostime()读取的一个时间单位约为55ms */ printf(It took %f msn, time); /* 屏幕输出FFT的运算时间 */ /* draw(N); */ /* 不运行作图语句 */（3） 测量DFT运行时间的程序 以下给出的是用于测量DFT运算时间的有关程序语句，应将它们插入实验一DFT程序中main()函数内适当的位置。 unsigned int star

22、t,end; float time; int j; start=biostime(0,0); /* 读取开始时间 */ for(j=0;j50;j+) /* 循环计算DFT 50次 */ for(k=0;kN;k+) r=i=0.0; for(n=0;nN;n+) r=r+xn*cos(c*n*k); i=i+xn*sin(c*n*k); yk=r; wk=i; end=biostime(0,0); /* 读取结束时间 */ time=(end-start)/50.0*55; printf(It took %f msn,time); /* draw(N); */五、实验报告要求（1）简要说明F

23、FT程序的主要结构和组成部分，各部分的功能和作用如何。（2）给出各种长度下DFT、FFT的运行时间，比较和分析不同长度下的FFT的运算速度（相差倍数），比较和分析同样长度下的FFT与DFT的运算速度（相差倍数），并从理论上对这些数据进行分析和讨论。（3）通过本实验，谈谈你对FFT的进一步认识。实验三 IIR数字滤波器一、实验目的和要求1学习和掌握IIR数字滤波器的计算机实现方法和过程2运用计算机模拟验证滤波器的性能，获得对IIR数字滤波器的感性认识。3巩固和加深对所学的IIR数字滤波器设计方法的理解和掌握。二、实验内容1设计要求用双线性变换法设计一巴特沃什数字低通滤波器，要求其技术指标为通带内

24、s= 0.3，幅度衰减不小于k2= -15dB。采用级联型结构，在计算机上编程实现所设计的IIR数字滤波器，验证其滤波性能。2设计过程（1）对数字频率指标进行预畸变，得到模拟频率指标（取T1） （2）根据p、s和k1 、k2设计巴特沃什模拟低通滤波器确定滤波器的阶数： 取N6。确定滤波器的截止频率： 确定模拟滤波器的系统函数：由N=6查表，并将c= 0.72729代入，经整理后得到 （3）根据双线性变换法，将模拟滤波器转换成数字滤波器： 将代入上式，经整理得到数字滤波器的系统函数（4）用二阶级联结构实现上述数字滤波器：x(n) w1(n) 0.08338 y1(n) w2(n) 0.08338

25、 y2(n) w3(n) 0.08338 y(n) z-1 z-1 z-1 1.31432 0.16676 1.0541 0.16676 0.94592 0.16676 z-1 z-1 z-1-0.71489 0.08338 -0.37534 0.08338 -0.23432 0.08338由上述结构，写出系统的输入输出方程： w1(n)=x(n)+1.31432w1(n1)0.71489w1(n2)y1(n)=0.08338 w1(n )+0.16676w1(n1)+0.08338w1(n2) w2(n)= y1(n)+1.0541w2(n1)0.37534w2(n2) y2(n)=0.08

26、338 w2(n )+0.16676w2(n1)+0.08338w2(n2) w3(n)= y2(n)+0.94592w3(n1)0.23422w3(n2) y(n)=0.08338 w3(n )+0.16676w3(n1)+0.08338w3(n2)初始条件为n 0时，x(n)=w1(n) =w2(n) =w3(n)=0。3计算机实现在理解和掌握以上设计过程的基础上，根据系统的输入输出方程，编制程序实现滤波器的计算，并验证其滤波性能。滤波器的输入仍采用正弦抽样信号，方法同实验一和实验二。其频率f、取样间隔T、取样点数N仍通过人机对话方式输入，以便调整。为了直观地看出系统的滤波性能，程序中需作

27、图画出输入和输出信号的波形。三、实验步骤该数字滤波器的截止频率c=2tan-1(c/2)=0.69756=0.22204。给定系统采样周期T，则对应的模拟截止频率为若取系统的采样周期T=0.001 s，则fc=116 Hz。因此，50 Hz的正弦信号应无衰减地通过该滤波器，120 Hz左右的信号处于过度带，有一定的衰减；150 Hz以上的信号处于阻带，几乎不能通过滤波器。 实验可分三步进行：（1）令T=0.001 s，N=100，分别输入50 Hz、120 Hz、150 Hz和200 Hz的正弦信号，观察输出波形，并与输入进行比较，验证滤波器的性能。（2）将输入信号分别改变成50 Hz与80

28、Hz、50 Hz与150 Hz、50 Hz与200 Hz两正弦信号的迭加，再观察滤波器的输入输出波形，体会和评价滤波结果。（3）改变取样间隔T=0.0001 s（减小10 倍），则滤波器的模拟截止频率应增大10倍，即fc=1160 Hz。分别输入500 Hz、1200、1500和2000的正弦信号，观测滤波器输出，验证fc。由此体会数字频率和模拟频率之间与系统取样间隔T的关系。（4）有余力的同学，可尝试将该低通数字滤波器转换成高通数字滤波器，再进行实验。四、C语言参考程序/* IIR数字滤波器实现程序 */ #include #include #include #include #define PI 3.1415926 float x100,y100; void draw(); void axis(); main() int N,f,n; float T,w1100,y1100,w2100,y2100,w3100;/* 键盘输入f、N、T */ printf(The frequency of the Sine wave f=); scanf(%