FIR滤波器的DSP设计 FFT的DSP实验

1. 代码被修改的部分： #define S2_Freq 160 // 信号 2 频率 void init_fft_tab(void) //输入波形的初始化
{ float wt1; float wt2; inti; for (i=0;i<Sample_Numb;i++) { wt1=2*pi*i*S1_Freq; wt1=wt1/SampleFreq; wt2=2*pi*i*S2_Freq; wt2=wt2/SampleFreq; //此处实现 60Hz 和 160Hz 两个信号相乘。 SignalInput[i]=cos(wt1)*cos(wt2)*100; //SignalInput[i]=(cos(wt1)+cos(wt2))/2*100; } } 由时域相乘频域卷积可以知道 60Hz 和 160Hz 的余弦信号相乘将得到有 100Hz 和 220Hz 频率分量的信号。 2. CCS 仿真图形如下： A、时域信号波形为（图形入口地址：SignalInput） ：
b) 用 SPtool 验证所设计的滤波器的功能 在 FDATool 里设计了一个 FIR 低通滤波器后，可以用 SPTool 来分析所设计的滤波器是否 实现了功能与 FDATool 一样。SPTool 可以对信号进行时域与频域的分析。测试滤波器是否实 现功能的主要步骤： 在 MATLAB 环境下产生一个混叠波信号，在 MATLAB 命令窗口输入： Fs = 20000;
} void main() { int p=0; n=0; wavein(); for(n=0;n<N;n++ ) { p=0; s=0; while((p<order)&&(p<=n)) { s=fir[p]*input[n-p]+s; p++; } output[n]=s; } }
的频谱：
c)
输出波形的频谱：
一、 通过代码修改、仿真，理解定点运算和浮点运算相比的优缺点。 定点运算的优点是：效率高，硬件单元比较简单，适合在低成本或者一些高速的 DSP 上应用；缺点主要就是丢失了精度，并且计算过程容易溢出。 浮点运算的优点很明显，浮点数有更大的数字表示范围，数据运算时不容易产生溢出， 运算精度也比定点运算高；但是其缺点是运算复杂，硬件较复杂，带有浮点运算单元的 DSP 明显比定点 DSP 价格高，增加了成本。 二、 通过仿真波形的分析、理解滤波器阶数与滤波器延时以及滤波器的阻带衰减之间 的关系。 在利用窗函数设计的 FIR 滤波器中，用同一种窗函数设计的 FIR 数字滤波器的阶数 N 不 同，所得到的滤波器的过渡带会不一样，N 越大，过渡带越窄。滤波器的阻带衰减只有窗形 状决定，不受滤波器阶数的影响；但是滤波器的阶数太高，会导致过渡带很窄，暂态增加， 延时很长。 三、 实验程序： #include <math.h> #include <stdlib.h> #define order 14 #define N 128 #define pi 3.14159265358979323846 float fir[order] = {0,0.00254885266601421,0.0147553639227699, 0.0448623472097668,0.0949715275858371,0.152016225973152,0.19084568264246, 0.19084568264246,0.152016225973152,0.0949715275858371,0.0448623472097668, 0.0147553639227699,0.00254885266601421,0}; float s; inti,n; float input[N],output[N]; int fs=20000; int f1=1000; int f2=5000; voidwavein() { for(n=0;n<N;n++ ) input[n]= sin(2*pi*n*f2/fs)+sin(2*pi*n*f1/fs);
B、对应谱为（图形入口地址：SignalInput） ：
C、FFT 输出结果（图形入口地址：FFT_W） ：
分析采样率、FFT 点数与频谱分辨率之间的关系 由采样定理知：采样率 Fs 应当满足������ 。频率 ������ ≥ 2������ℎ （������ℎ 为信号的最高频率） 分辨率为：������0 = ������ ������ ������ = 1 ������������������ = 1 ������ （N 为 FFT 的点数，Ts 为采样间隔，NTs 就是采样前模拟信号的时间长度 T） ， 可以得到信号长度越长， 频率分辨率越好。
实验内容：
一、 完成所给例程对应实验，验收结果如下： A、时域信号波形为（图形入口地址：SignalInput） ：
B、对应谱为（图形入口地址：SignalInput） ：
C、FFT 输出结果（图形入口地址：FFT_W） ：
二、
修改代码为时域信号 60Hz 和 160Hz 两个信号相乘，分析其频谱，并说明理由。
FIR 滤波器的 DSP 设计
作业要求：
一、 修改代码：滤波器设计要求：线性相位低通滤波器，模拟信号通带为 0~1kHz，阻 带衰减在 5kHz 处为-60dB。采样率为 20kHz。 1. FIR 滤波器设计流程为，根据给出的滤波器设计（通带带宽、阻带衰减、相位要求）要 求， 在 matlab 环境下仿真得到滤波器系数， 得到了滤波器系数后在 CCS 下编程设计 DSP 的 FIR 滤波器。 a) 使用 matlab 中的 FDAtool 设计 FIR 滤波器 MATLAB 工具箱中的 FDATool 是通过指定滤波器的性能指标来快速设计数字滤波器, 它 是一种交互式工具，只要在 FDATool 工具里选择相应的滤波器设计方法，并填入需要的技术 指标，如通带阻带截止频率、阶数等，即可生成所需的滤波器。并可以在界面上可以看到幅 频相频、零极点分布等图形。 在本设计中：选择 fir 低通滤波器，阶数为 13，采样率为 20KHz，截止频率 fc=1KHz，窗 口设置为 Blackman 窗。
t = (0:500)/Fs; s = sin(2*pi*t*1000)+sin(2*pi*t*5000); subplot(1,1,1); plot(t, s); axis([0 0.005 -2 2]); grid
产生了 1000Hz 和 5000Hz 的混叠波，其波形显示图形为：
c)
在 FDATool 的界面里输出到 SPTool，然后在将导入 SPTool 的滤波器命名为 filter，点击
“Apply” ，这样，设计的滤波器就导入了 SPTool 中。在 SPtool 主界面里面，将输入信号 S，采样 频率 Fs 导入进去，并生成名为 output 的输出信号。通过点击该界面下的“view”可以得到输入 输出波形，点击“create”可以得到输出的频谱，这三个图形如下所示：
通过观察输入输出的波形图及频谱图验证所设计的滤波器的功能。 可以看出 5000HZ 的 频率分量被滤，所设计的滤波器实现了低通滤波功能。 2. 使用 CCS 编程实现 DSP 的 FIR 滤波 启动 CCS 软件，编写及修改程序，我们得到的仿真波形如下： a) 输入与输出波形：
四、
B、对应谱为（图形入口地址：SignalInput） ：
由频谱图可以看出，60Hz 和 160Hz 相乘产生的信号的频率分量有 100Hz 和 220Hz。
C、FFT 输出结果（图形入口地址：FFT_W） ：
三、 修改实验代码为 1024 点 FFT，进行程序的编译、运行和对应结果的查看。 1、代码被修改的部分（时域信号为 60Hz 和 180Hz 的两个波形相加） ： #define Sample_Numb 1024 //设置 FFT 点数为 1024 点 voidfft(float datar[Sample_Numb],float datai[Sample_Numb]) { int x0,x1,x2,x3,x4,x5,x6,x7,x8,x9,xx; inti,j,k,b,p,L; floatTR,TI,temp; //倒位序；由 1024 点可以得到顺序数可以由 10 位数表示，这里取 0-9； for(i=0;i<Sample_Numb;i++) { x0=x1=x2=x3=x4=x5=x6=x7=x8=0; x0=i&0x01;x1=(i/2)&0x01;x2=(i/4)&0x01;x3=(i/8)&0x01; x4=(i/16)&0x01;x5=(i/32)&0x01;x6=(i/64)&0x01;x7=(i/128)&0x01; x8=(i/256)&0x01;x9=(i/512)&0x01; xx=x0*512+x1*256+x2*128+x3*64+x4*32+x5*16+x6*8+x7*4+x8*2+x9; datai[xx]=datar[i]; } for(i=0;i<Sample_Numb;i++) { datar[i]=datai[i];datai[i]=0; } //10 级蝶形运算 for(L=1;L<=10;L++) { b=1;i=L-1; while(i>0) { b=b*2;i--; } for(j=0;j<=b-1;j++) {//
p=1;i=10-L; while(i>0) { p=p*2;i--;} p=p*j; for(k=j;k<1024;k=k+2*b) { TR=datar[k];TI=datai[k];temp=datar[k+b]; datar[k]=datar[k]+datar[k+b]*cos_tab[p]+datai[k+b]*sin_tab[p]; datai[k]=datai[k]-datar[k+b]*sin_tab[p]+datai[k+b]*cos_tab[p]; datar[k+b]=TR-datar[k+b]*cos_tab[p]-datai[k+b]*sin_tab[p]; datai[k+b]=TI+temp*sin_tab[p]-datai[k+b]*cos_tab[p]; } } } //计算 FFT 输出的功率 for(i=0;i<Sample_Numb/2;i++) { FFT_W[i]=sqrt(datar[i]*datar[i]+datai[i]*datai[i]); } } A、时域信号波形为（图形入口地址：SignalInput） ：