数字信号处理实验指导书.docx
《数字信号处理实验指导书.docx》由会员分享,可在线阅读,更多相关《数字信号处理实验指导书.docx(26页珍藏版)》请在冰豆网上搜索。
数字信号处理实验指导书
实验一离散傅立叶变换
一、实验目的
1.在理论学习的基础上,通过本实验,加深对离散傅立叶变换(DFT)的理解。
2.以正弦信号为例,学习和掌握利用离散傅立叶变换(DFT)分析信号频谱的方法。
3.了解应用DFT进行信号频谱分析过程中可能出现的问题,以便在实际中正确应用DFT。
4.学习和熟悉在计算机上作图显示和研究信号的离散频谱
二、实验原理与方法
在各种信号序列中,有限长序列的数字信号处理占有很重要地位,对有限长序列,我们可以使用离散傅立叶变换(DFT),这一变换不但可以很好地反映序列的频谱特性,而且易于用快速算法在计算机上实现,当序列x(n)的长度为N时,它的DFT定义为:
X(k)=
反变换为:
x(n)=
有限长序列的DFT是其Z变换在单位圆上的等距采样,或者说是序列傅立叶变换的等距采样,因此可以用于序列的谱分析。
DFT要求被分析序列的长度为有限长,对无限长和很长的序列则需要截取。
截取的结果使得信号的频谱分析出现频率泄漏。
本实验分析的周期性信号,对周期信号,当截取的序列长度是信号周期的整数倍的时候,DFT分析的频率离散点恰好采在周期信号的谐波分量和其它零值频率点上,能正确反映信号的频谱。
当截取的序列长度不为信号周期的整数倍时,DFT分析的结果就会有频率泄漏的反映。
三、实验内容
1.设计说明
编制程序用DFT分析正弦信号的频谱。
具体包括以下几部分
(1)编程产生正弦抽样输入信号
x(t)=sin(Ωt)=sin(2πft)
x(n)=x(nT)=sin(2πfnT)
正弦信号的频率f、取样点数N和取样间隔T通过人机对话的方式由键盘输入,以便调整。
(2)编程实现DFT的运算
DFT的计算公式为
(3)计算DFT的幅度
(4)作图画出DFT的幅度谱线。
2.实验步骤
将编制的程序输入计算机,调试运行正确后,进行以下实验。
(1)分别在下列几种情况下,输入正弦信号的频率f、采样点数N和采样周期T,观察和研究DFT的结果,并记录谱线的大体形状。
(a)f=50Hz,N=32,T=0.000625s
(b)f=50Hz,N=32,T=0.005s
(c)f=50Hz,N=32,T=0.004s
(d)f=150Hz,N=32,T=0.000625s
(2)修改程序,将输入信号改变成50Hz和150Hz两个正弦信号分量的迭加,并令50Hz分量的幅度为2,150Hz分量的幅度为1。
即
x(n)=2sin(2π50nT)+sin(2π150nT)
取N=32,分别在T=0.000625s、T=0.0025s和T=0.002s情况下进行DFT分析,观察和记录结果。
(3)有余力的同学还可以将输入信号改成50Hz的周期性方波,再对其进行DFT分析。
四、C语言参考程序
/*用DFT分析正弦信号频谱的C语言程序*/
#include
#include
#include
#include
#definePI3.1415926
floatx[1024],y[1024],w[1024];
voiddraw(int);
voidaxis(int,int);
main()
{
intN,f,n,k;
floatT,r,i,c;
/*键盘输入f、N、T*/
printf("ThefrequencyoftheSinewavef=");
scanf("%d",&f);
printf("ThenumberofsamplesN=");
scanf("%d",&N);
printf("ThesamplingperiodT=");
scanf("%f",&T);
c=2*PI/N;
/*产生正弦抽样信号*/
for(n=0;nx[n]=sin(2*PI*f*n*T);
}
/*计算DFT*/
for(k=0;kr=i=0.0;
for(n=0;nr=r+x[n]*cos(c*n*k);
i=i+x[n]*sin(c*n*k);
}
y[k]=r;
w[k]=i;
}
/*计算DFT的幅度*/
for(k=0;kw[k]=sqrt(y[k]*y[k]+w[k]*w[k]);
}
/*画频谱图*/
draw(N);
getch();
closegraph();
}
voiddraw(intN)/*作图*/
{
intk,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,"16");
outtextxy((x+step-19),y+8,
"01234567890123456789012345678901");
/*绘制频谱图*/
for(k=0;kline(x,y,x,(y-floor(w[k]*7)));
x+=step;
}
}
voidaxis(intx,inty)/*绘制横轴与纵轴*/
{
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-260),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运算次数的一种快速算法。
它是对DFT变换式进行一次次分解,使其成为若干小点数的组合,从而减少运算量。
常用的FFT是以2为基数的,其长度N=2M。
它的效率高,程序简单,使用非常方便,当要变换的序列长度不等于2的整数次方时,为了使用以2为基数的FFT,可以用末尾补零的方法,使其长度延长至2的整数次方。
本实验运用时间抽取基2FFT,其原理、信号流图和运算过程可参阅课堂笔记、教材和其它教科书。
FFT的实现要比DFT复杂,通常采用三个嵌套循环来实现。
最外面的循环是分级循环,N=2M点的FFT分为M级计算。
中间一层是分组循环,一级内蝶形系数Wk相同的蝶形构成一组。
最内层为蝶形计算的循环,一组内不同输入数据的蝶形逐个计算。
编程时要注意各级蝶形组之间的间隔、组内蝶形之间间隔、以及系数Wk变化。
有不少书中有FFT程序可供参考,但要注意理解和弄懂,不可一味照搬,要尽量自己去编。
三、实验内容
1.设计说明
编制时间抽取基2FFT程序计算前面DFT程序分析过的正弦信号,与DFT
计算的结果进行比较,以验证所编程序的正确性。
FFT为复数运算,若用实数运算来实现,需要设置两个数组来存放输入输出数据。
其中一个用于存放数据的实部,另一个存放数据的虚部。
正弦输入信号为实数,则令其虚部为零。
同样,碟形运算也要化成实数来进行,分别算出实部和虚部。
当然也可以直接用复数数组和语句实现。
程序设计的难点和重点在于要合理安排和正确设置碟形运算的分级循环、级内分组循环和组内分碟形循环。
要注意乘法系数WNk的变化以及各循环变量的变化和调整。
程序中,在进行FFT计算之前,先要将输入数据排列成二进制倒序的形式,这可通过将有关数据相互换位来实现。
正弦抽样信号的产生和输出频谱图的绘制与前面的DFT程序相同。
2.实验步骤
(1)将编制的程序输入计算机,调试、运行正确。
即输入实验一中的1~2组数
据,看FFT输出谱线的结果是否与前面DFT的结果一致。
(2)在不同的FFT长度下(可分别令N=64,128,256,512,1024),运行FFT
程序,观察、比较和记录运行时间。
由于现在计算机的运算速度很快,计算一次FFT的时间很短(在ms数量级),不便观察和记录。
可在程序中添置循环语句,让其重复计算FFT许多次(可重复100次、1000次或更多,视情况而定),使总的运算时间达到数秒或数十秒钟,从而可以通过人工来观察和记录。
此时注意,为了防止FFT多次循环计算,使输入输出数据越来越大而发生溢出,可令输入数据全部为零。
这样,FFT的运算量保持不变,但输入输出始终为零,不会因多次循环而发生溢出。
在观测FFT的运算时间时,作图部分的语句不要运行(一方面,它不应该计入FFT或DFT的运算时间。
另一方面,当N较大时,原有程序中的作图语句会出错),可将其注释掉。
测量FFT运行时间的一个更好、更精确的方法是在程序中调用读取计算机时钟的函数。
TURBOC中提供多种读取计算机时间的方法和函数,一种比较简单的方法是用
biostime()函数,后面给出了相应的参考程序。
(3)用上述类似的方法,在不同的长度下(N=64,128,256,512,1024),运
行实验一中的DFT程序,记录运行时间,比较它与FFT的速度差距。
四、C语言参考程序
(1)用FFT分析正弦信号频谱的程序
#include
voidfft(x,y,n)/*FFT作为子程序*/
intn;
floatx[1024],y[1024];
{
inti,j,k,l,m,n1,n2;
floatc,s,e,tr,ti;
/*计算FFT的级数M*/
for(j=1,i=1;i{
m=i;
j=2*j;
if(j==n)break;
}
/*改变输入数据的顺序*/
n1=n-1;
for(j=0,i=0;i{
if(i{
tr=x[j];
ti=y[j];
x[j]=x[i];
y[j]=y[i];
x[i]=tr;
y[i]=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;k{
c=cos(k*e);/*计算WNk*/
s=-sin(k*e);
for(i=k;i{
j=i+n2;
tr=c*x[j]-s*y[j];/*计算碟形*/
ti=c*y[j]+s*x[j];
x[j]=x[i]-tr;
y[j]=y[i]-ti;
x[i]=x[i]+tr;
y[i]=y[i]+ti;
}
}
}
}
/*主程序部分*/
#include
#include
#include
#include
#definePI3.1415926
floatx[1024],y[1024],w[1024];
voiddraw(int);
voidaxis(int,int);
main()
{
intN,f,n,k;
floatT;
/*键盘输入f、N、T*/
printf("thefrequencyofthesinewavef=");
scanf("%d",&f);
printf("thenumberofsamplesN=");
scanf("%d",&N);
printf("thetimestepofsamplesT=");
scanf("%f",&T);
/*产生正弦输入数据*/
for(n=0;n{
x[n]=sin(2*PI*f*n*T);
y[n]=0;
}
/*调用子程序计算FFT*/
fft(x,y,N);
/*计算FFT的幅度*/
for(k=0;k{
w[k]=sqrt(x[k]*x[k]+y[k]*y[k]);
}
/*画FFT的幅度谱*/
draw(N);
getch();
closegraph();
}
/*作图子程序*/
voiddraw(intN)
{
intk,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,"16");
outtextxy((x+step-19),y+8,
"01234567890123456789012345678901");
for(k=0;k{
line(x,y,x,(y-floor(w[k]*7)));
x+=step;
}
}
voidaxis(intx,inty)
{
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-260),64,120);
line((x+4),(y-260),64,120);
}
(2)测量FFT运行时间的程序
以下给出的是用于测量FFT运算时间的有关程序语句,应将它们插入前面FFT程序中main()函数内适当的位置。
unsignedintstart,end;
doubletime;inti;/*提示:
这两行放在主程序开头数据类型说明部分*/
/*提示:
以下部分替代主程序中的fft(x,y,N)*/
start=biostime(0,0);/*读取开始时间*/
for(n=0;n<2000;n++)/*循环计算FFT2000次*/
{
fft(x,y,N);
for(i=0;i{
x[i]=0;
y[i]=0;
}
}
end=biostime(0,0);/*读取结束时间*/
time=(end-start)/2000.0*55;/*biostime()读取的一个时间单位约为55ms*/
printf("Ittook%fms\n",time);/*屏幕输出FFT的运算时间*/
/*draw(N);*//*不运行作图语句*/
(3)测量DFT运行时间的程序
以下给出的是用于测量DFT运算时间的有关程序语句,应将它们插入实验一DFT程序中main()函数内适当的位置。
unsignedintstart,end;
floattime;
intj;
start=biostime(0,0);/*读取开始时间*/
for(j=0;j<50;j++){/*循环计算DFT50次*/
for(k=0;kr=i=0.0;
for(n=0;nr=r+x[n]*cos(c*n*k);
i=i+x[n]*sin(c*n*k);
}
y[k]=r;
w[k]=i;
}
}
end=biostime(0,0);/*读取结束时间*/
time=(end-start)/50.0*55;
printf("Ittook%fms\n",time);
/*draw(N);*/
五、实验报告要求
(1)简要说明FFT程序的主要结构和组成部分,各部分的功能和作用如何。
(2)给出各种长度下DFT、FFT的运行时间,比较和分析不同长度下的FFT的运
算速度(相差倍数),比较和分析同样长度下的FFT与DFT的运算速度(相差倍数),并从理论上对这些数据进行分析和讨论。
(3)通过本实验,谈谈你对FFT的进一步认识。
实验三IIR数字滤波器
一、实验目的和要求
1.学习和掌握IIR数字滤波器的计算机实现方法和过程
2.运用计算机模拟验证滤波器的性能,获得对IIR数字滤波器的感性认识。
3.巩固和加深对所学的IIR数字滤波器设计方法的理解和掌握。
二、实验内容
1.设计要求
用双线性变换法设计一巴特沃什数字低通滤波器,要求其技术指标为通带内ω<ωp=0.2π,幅度衰减不大于k1=-1dB;阻带内ω>ωs=0.3π,幅度衰减不小于k2=-15dB。
采用级联型结构,在计算机上编程实现所设计的IIR数字滤波器,验证其滤波性能。
2.设计过程
(1)对数字频率指标进行预畸变,得到模拟频率指标(取T=1)
(2)根据Ωp、Ωs和k1、k2设计巴特沃什模拟低通滤波器
确定滤波器的阶数:
取N=6。
确定滤波器的截止频率:
确定模拟滤波器的系统函数:
由N=6查表,并将Ωc=0.72729代入,经整理后得到
(3)根据双线性变换法,将模拟滤波器转换成数字滤波器:
将
代入上式,经整理得到数字滤波器的系统函数
(4)用二阶级联结构实现上述数字滤波器:
x(n)w1(n)0.08338y1(n)w2(n)0.08338y2(n)w3(n)0.08338y(n)
z-1z-1z-1
1.314320.166761.05410.166760.945920.16676
z-1z-1z-1
-0.714890.08338-0.375340.08338-0.234320.08338
由上述结构,写出系统的输入输出方程:
w1(n)=x(n)+1.31432w1(n–1)–0.71489w1(n–2)
y1(n)=0.08338w1(n)+0.16676w1(n–1)+0.08338w1(n–2)
w2(n)=y1(n)+1.0541w2(n–1)–0.37534w2(n–2)
y2(n)=0.08338w2(n)+0.16676w2(n–1)+0.08338w2(n–2)
w3(n)=y2(n)+0.94592w3(n–1)–0.23422w3(n–2)
y(n)=0.08338w3(n)+0.16676w3(n–1)+0.08338w3(n–2)
初始条件为n<0时,x(n)=w1(n)=w2(n)=w3(n)=0。
3.计算机实现
在理解和掌握以上设计过程的基础上,根据系统的输入输出方程,编制程序实现滤波器的计算,并验证其滤波性能。
滤波器的输入仍采用正弦抽样信号,方法同实验一和实验二。
其频率f、取样间隔T、取样点数N仍通过人机对话方式输入,以便调整。
为了直观地看出系统的滤波性能,程序中需作图画出输入和输出信号的波形。
三、实验步骤
该数字滤波器的截止频率ωc=2tan-1(Ωc/2)=0.69756=0.22204π。
给定系统采样周期T,则对应的模拟截止频率为
若取系统的采样周期T=0.001s,则fc=116Hz。
因此,50Hz的正弦信号应无衰减地通过该滤波器,120Hz左右的信号处于过度带,有一定的衰减;150Hz以上的信号处于阻带,几乎不能通过滤波器。
实验可分三步进行:
(1)令T=0.001s,N=100,分别输入50Hz、120Hz、150Hz和200Hz的正弦信号,观察输出波形,并与输入进行比较,验证滤波器的性能。
(2)将输入信号分别改变成50Hz与80Hz、50Hz与150Hz、50Hz与200Hz两正弦信号的迭加,再观察滤波器的输入输出波形,体会和评价滤波结果。
(3)改变取样间隔T=0.0001s(减小10倍),则滤波器的模拟截止频率应增大10倍,即fc=1160Hz。
分别输入500Hz、1200、1500和2000的正弦信号,观测滤波器输出,验证fc。
由此体会数字频率和模拟频率之间与系统取样间隔T的关系。
(4)有余力的同学,可尝试将该低通数字滤波器转换成高通数字滤波器,再进行实验。
四、C语言参考程序
/*IIR数字滤波器实现程序*/
#include
#include
#include
#include
#definePI3.1415926
floatx[100],y[100];
voiddraw();
voidaxis();
main()
{
intN,f,n;
floatT,w1[100],y1[100],w2[100],y2[100],w3[100];
/*键盘输入f、N、T*/
printf("ThefrequencyoftheSinewavef=");
scanf("%