数字信号处理实验.docx
- 文档编号:7049782
- 上传时间:2023-01-16
- 格式:DOCX
- 页数:24
- 大小:548.12KB
数字信号处理实验.docx
《数字信号处理实验.docx》由会员分享,可在线阅读,更多相关《数字信号处理实验.docx(24页珍藏版)》请在冰豆网上搜索。
数字信号处理实验
实验报告
课程名称:
数字信号处理实验
实验项目:
FFT作谱分析
FIR数字滤波器设计与软件实现有
限IIR数字滤波器设计与软件实现
实验时间:
实验班级:
姓名:
学号:
指导教师:
学院实验室
二〇一五年月日
实验一:
FFT作谱分析
一、实验目的:
1.进一步加深对DFT算法原理和基本性质的理解(因为FFT只是DFT的一种快速算法,所以FFT的运算结果必然满足DFT的基本性质)。
2.熟悉FFT算法原理和FFT子程序的应用。
3.学习用FFT对连续信号和时域离散信号进行谱分析的方法,了解可能出现的分析误差及其原因,以便在实际中正确应用FFT。
二、实验原理:
1,快速傅立叶变换(fft)算法
长度为n的序列x(n)的离散傅立叶变换x(k)为:
N点的DFT可以分解为两个N/2点的DFT,每个N/2点的DFT有可以分解为两个N/4点的DFT,依次类推,当N为2的整数次幂N=2^m,由于分解依次降低一阶幂次,所以通过M次分解,最后全部成为一系列2点的DFT运算,以上就是按时间抽取的快速傅立叶变化(FFT)算法,当需要进行变化的序列长度不是2的整数次方的时候,为了使用以2为基的FFT,可以用末尾补零的方法,使其长度延至2的整数次方。
序列x(k)的离散傅立叶反变换为:
离散傅立叶反变换与正变换区别在于Wn变为
并多了一个1/N,对于按时间抽取的快速傅立叶变换算法并无实质区别,因此将反变换算法合并在同一个程序
2、利用FFT进行频谱分析
若信号本身是有限长的序列,计算序列就是直接对序列进行FFT运算求得x(k),
x(k)就代表序列在[0,2π]之间的频谱值
三、实验器材
微型计算机,实验箱
四、实验程序:
#include"DSP2833x_Device.h"
#include"DSP2833x_Examples.h"
#include
#definepi3.141593//float小数点后6位
#defineNL256//NL与N的值必须是一致的,即NL=N
intN=256;
floathn[64]={0};//滤波系数补32个零后存放
floatSampleTable[NL];//实数信号序列
//窗函数类型:
汉宁窗,32阶低通
//截止频率:
500Hz
//采样频率:
12800Hz
floath[32]=
{
0,0.0006485134778915,0.002574789225444,0.005709506486408,
0.009931872420821,0.01507448912514,0.02093024026384,0.02726090869226,
0.03380715986684,0.04029946500149,0.04646949558203,0.05206149892494,
0.05684316398774,0.06061550768029,0.06322135360436,0.06455203566049,
0.06455203566049,0.06322135360436,0.06061550768029,0.05684316398774,
0.05206149892494,0.04646949558203,0.04029946500149,0.03380715986684,
0.02726090869226,0.02093024026384,0.01507448912514,0.009931872420821,
0.005709506486408,0.002574789225444,0.0006485134778915,0
};
structComplex//定义复数结构体
{
floatreal,imag;
};
structComplexWn;//定义旋转因子
structComplexVn;//每一级第一个旋转因子虚部为0,实部为1
structComplexT;//存放旋转因子与X(k+B)的乘积
//floatRealin[NL]={0.0};//AD采样输入的实数
floatoutput[NL];//输出的FFT幅值(复数的模)
structComplexSample[NL];//AD采样输入的实数转化为复数
structComplexH[64];//补零后h(n)存放的复数数组
structComplexY[64];
structComplexX[64];
floatResult[288];//
structComplexMUL(structComplexa,structComplexb)//定义复乘
{
structComplexc;
c.real=a.real*b.real-a.imag*b.imag;
c.imag=a.real*b.imag+a.imag*b.real;
return(c);
}
voidMYFFT(structComplex*xin,intN)
{
intL=0;//级间运算层
intJ=0;//级内运算层
intK=0,KB=0;//蝶形运算层
intM=1,Nn=0;//N=2^M
floatB=0;//蝶形运算两输入数据间隔
/*以下是为倒序新建的局部变量*/
intLH=0,J2=0,N1=0,I,K2=0;
structComplexT;
/*以下是倒序*/
LH=N/2;//LH=N/2
J2=LH;
N1=N-2;
for(I=1;I<=N1;I++)
{
if(I { T=xin[I]; xin[I]=xin[J2]; xin[J2]=T; } K2=LH; while(J2>=K2) { J2-=K2; K2=K2/2;//K2=K2/2 } J2+=K2; } /*以下为计算出M*/ Nn=N; while(Nn! =2)//计算出N的以2为底数的幂M { M++; Nn=Nn/2; } /*蝶形运算*/ for(L=1;L<=M;L++)//级间 { B=pow(2,(L-1)); Vn.real=1; Vn.imag=0; Wn.real=cos(pi/B); Wn.imag=-sin(pi/B); for(J=0;J { for(K=J;K { KB=K+B; T=MUL(xin[KB],Vn); xin[KB].real=xin[K].real-T.real; xin[KB].imag=xin[K].imag-T.imag; xin[K].real=xin[K].real+T.real; xin[K].imag=xin[K].imag+T.imag; } Vn=MUL(Wn,Vn);//旋转因子做复乘相当于指数相加,得到的结果 //和J*2^(M-L)是一样的,因为在蝶形因子运算 //层中M与L都是不变的,唯一变x化的是级内的J //而且J是以1为步长的,如J*W等效于W+W+W...J个W相加 } } } /* voidFilterDC(structComplex*ADC,intN)//去除数据中的直流成分,否则直流分量将很大 { inti; floatsum=0; for(i=0;i {sum+=ADC[i].real;} sum=sum/N; for(i=0;i {ADC[i].real-=sum;} }*/ voidModelComplex(structComplex*xin,intN,float*out_put) { inti; for(i=0;i { out_put[i]=sqrt(xin[i].real*xin[i].real+xin[i].imag*xin[i].imag)*2/N; } } voidmain(void) { intcounter,i=0,j=0,k=0,full=0,mid0=0,mid1=0,qt=0; InitSysCtrl(); DINT; InitPieCtrl(); IER=0x0000; IFR=0x0000; for(counter=0;counter<32;counter++)//对h(n)从实数变成复数,并且补零 { X[counter].real=0;//X[64]初始化 X[counter].imag=0; X[counter+32].real=0;//X[64]初始化 X[counter+32].imag=0; H[counter].imag=0; H[counter].real=h[counter];//前32点与h[32]一样 H[counter+32].real=0;//后32点为0 H[counter+32].imag=0; } MYFFT(H,64);//完成补零后h(n)的FFT ModelComplex(H,64,hn);//h(n)的频谱hn //可以观察到滤波器的频谱 for(counter=0;counter<288;counter++)//输出信号初始化 { Result[counter]=0; if(counter {output[counter]=0;} } for(i=0;i { SampleTable[i]=sin(2*pi*5*i/(NL-1))+sin(2*pi*i*5*3/(NL-1))/3+sin(2*pi*i*5*5/(NL-1))/5; Sample[i].real=SampleTable[i];//实数转复数 Sample[i].imag=0; //DELAY_US(10); } for(i=0;i<8;i++) { for(j=0;j<32;j++)//取得512点采样的分段L=32,后32点为0 { X[j].real=Sample[full].real;//h[j];//在初始化时X[64]已经全赋值0 full++; } MYFFT(X,64);//分段后的64点采样FFT for(k=0;k<64;k++)//64点的X与H相乘 { Y[k]=MUL(X[k],H[k]); X[k].real=0; X[k].imag=0; } /*************************** 以下做逆FFT(IFFT) ****************************/ for(counter=0;counter<64;counter++)//得到X(k)的共轭 { Y[counter].real=Y[counter].real; Y[counter].imag=-Y[counter].imag; } MYFFT(Y,64);//共轭的FFT for(counter=0;counter<64;counter++)//得到共轭X(k)的共轭 { Y[counter].imag=-Y[counter].imag; } for(counter=0;counter<64;counter++)//除以N { Y[counter].real=Y[counter].real/64; Y[counter].imag=Y[counter].imag/64; if(fabs(Y[counter].imag)<=0.00001) {Y[counter].imag=0;} } Y[0].real=0; Y[0].imag=0; Y[1].real=0; Y[1].imag=0; ModelComplex(Y,64,output);//取复数的模 for(qt=0;qt<64;qt++)//取模后全为正数,需要恢复原来的负数 { if(Y[qt].real<0) {output[qt]=-output[qt];}//output中存放着IFFT的实数结果 }//IFFT结束 qt=0; mid0=i*32; mid1=mid0+64; for(counter=mid0;counter { if(mid0==0) { Result[counter]=output[qt]; }//counter达到543,但output只在前64中 else {Result[counter]+=output[qt];}//所以需要设qt来表示output数组里的 qt++; } qt=0; counter=0; } while (1) { 五、实验的波形 观察到当代码为: MYFFT(H,64);//完成补零后h(n)的FFT ModelComplex(H,64,hn);//h(n)的频谱hn //可以观察到滤波器的频谱 波形如下: 分析: 因为W=2Πf/F,有因为F>=2f,所以w=[0,π],所以判断出h(n)为低通滤波器 因为代码: SampleTable[i]=sin(2*pi*5*i/(NL-1))+sin(2*pi*i*5*3/(NL-1))/3+sin(2*pi*i*5*5/(NL-1))/5; 三个波形的幅度谱如图: 三个波形的频谱如图: 有因为代码: #definepi3.141593//float小数点后6位 #defineNL256//NL与N的值必须是一致的,即NL=N intN=256; //截止频率: 500Hz //采样频率: 12800Hz SampleTable[i]=sin(2*pi*5*i/(NL-1))+sin(2*pi*i*5*3/(NL-1))/3+sin(2*pi*i*5*5/(NL-1))/5; 可以算出Wc截至频率对应的点数是2.5 因为W3=5W1,W2=3W1 有因为64个点对应一个圆周2π,所以第一个波形w1是2*pi*5/(nl-1)带入pi和nl,算出对于点数为1.25,即w1=1.25 有因为W2=3W1所以带入算的结果为3.76,即W2=3.76 W3=5W1,所以带入算的结果为6.27,即W3=6.27 所以各点对应的波形图的位置如下图: 做逆FFT(IFFT)后的输出结果如下波形: 分析: 因为通过低通滤波器后,在截至频率之间的信号被还原出来,如上可知频率为W1的信号被还原出来,其他信号被滤去。 实验二: FIR数字滤波器设计与软件实现 1、实验目的 (1)掌握用窗函数设计FIR数字滤波器的原理和方法; (2)掌握用等波纹最佳逼近法设计FIR数字滤波器的原理与方法; (3)掌握FIR滤波器的快速卷积实现原理; (4)学会调用MATLAB函数设计与实现FIR滤波器。 2、实验原理 (1)窗函数法设计FIR数字滤波器的原理; 将模拟频率转化为数字频率,设取样时间为T(要满足抽样定理) Ωp=2π*fp*T Ωs=2π*fs*T 过渡带宽度△Ω=Ωp-Ωs (2)等波纹最佳逼近法设计FIR数字滤波器的原理. 1.MATLAB函数fir1的功能及其调用格式请查阅教材; 2.采样频率Fs=1000Hz,采样周期T=1/Fs; 3.根据图1(b)和实验要求,可选择滤波器指标参数: 通带截止频率fp=120Hz,阻带截至频率fs=150Hz,换算成数字频率,通带截止频率,通带最大衰为0.1dB,阻带截至频率,阻带最小衰为60dB。 3、实验器材 计算机与实验箱 4、实验程序 #include"DSP2833x_Device.h"//DSP2833xHeaderfileIncludeFile #include"DSP2833x_Examples.h"//DSP2833xExamplesIncludeFile #include"math.h" #definepi3.14159 #definea256//合成信号的采样点数 #defineb32//h(n)的序列长度 #defineca+b-1//y(n)的序列长度=a+b-1 //Uint16Input[256]; floatSample[a];//采样后的信号序列 floaty1[c];//输出序列 //窗函数类型: 汉宁窗,32阶低通 //截止频率: 500Hz //采样频率: 128000Hz floath[b]= { 0,0.0006485134778915,0.002574789225444,0.005709506486408, 0.009931872420821,0.01507448912514,0.02093024026384,0.02726090869226, 0.03380715986684,0.04029946500149,0.04646949558203,0.05206149892494, 0.05684316398774,0.06061550768029,0.06322135360436,0.06455203566049, 0.06455203566049,0.06322135360436,0.06061550768029,0.05684316398774, 0.05206149892494,0.04646949558203,0.04029946500149,0.03380715986684, 0.02726090869226,0.02093024026384,0.01507448912514,0.009931872420821, 0.005709506486408,0.002574789225444,0.0006485134778915,0 }; /*============================================================= 功能: 实现离散线性卷积 算法原理: 对位相乘求和法实现卷积 形参: xn(*x所指向的数组长度)、hn(*h所指向的数组长度)为参与 运算的两个卷积序列数组长度 *x、*h为指向两个数组的指针 *y指向输出序列的数组 注意: 输出序列长度为xn+hn-1;存放内容为输出数组的0~(xn+hn-2)。 另外,该函数实现的无符号整型的卷积运算,如果要实现浮点型,需要 将形参的Uint16hn,Uint16*x,Uint16*h,Uint16*y的Uint16数据格式 改为float或者double型 ==============================================================*/ voidLinearConvolution(Uint16xn,Uint16hn,float*x,float*h,float*y) { Uint16i,j,k,l; Uint16yn;//输出序列y的长度 yn=xn+hn-1; for(i=0;i k=yn-1; for(i=hn-1;i>0;i--)//将*h作为被乘数 { l=k; for(j=xn-1;j>0;j--)//数组x[n]的1~(xn-1)与h[i]逐一相乘 { y[l]+=h[i]*x[j]; l--; } y[l]+=x[0]*h[i]; k--; } l=k; for(j=xn-1;j>0;j--) { y[l]+=h[0]*x[j]; l--; } y[l]+=x[0]*h[0]; } voidmain(void) { Uint16i; InitSysCtrl(); InitPieCtrl(); IER=0x0000; IFR=0x0000;
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 数字信号 处理 实验