《TMS320C55x算法编程.ppt》由会员分享，可在线阅读，更多相关《TMS320C55x算法编程.ppt（53页珍藏版）》请在三一办公上搜索。
1、定标与溢出处理FIR滤波器IIR滤波器快速傅里叶变换(FFT),第7讲 DSP算法编程,7.1定标与溢出处理,数的定标溢出的处理方法常用信号处理算法中的定标方法,7.1.1 数的定标,小数定标的概念设定一个16位数的小数点处于该数中的哪一位可以表示不同大小和不同精度的小数Q表示法表7-1列出了一个16位数的16种Q表示及它们所能表示的十进制数值范围,表7-1 Q表示及数值范围,不同的Q所表示的数不仅范围不同，而且精度也不相同Q越大，数值范围越小，但精度越高Q越小，数值范围越大，但精度就越低例如，Q0的数值范围是-32768到+32767，其精度为1；Q15的数值范围为-1到0.9999695，
2、精度为 1/32768=0.00003051,对定点数而言，数值范围与精度是一对矛盾。一个变量要想能够表示比较大的数值范围，必须以牺牲精度为代价；而想提高精度，则数的表示范围就相应地减小。在实际的定点算法中，应该根据具体问题进行折衷处理,以达到最佳效果。,在C55x中，16位整数采用补码形式表示。每个采用Qi定标的16位数用1个符号位、i个小数位和15-i个整数位来表示。,表7-2 同样的数在不同定标方式下所表示的具体数值,同样一个16位数，若小数点设定的位置不同，它所表示的数也就不同。,7.1.2 溢出的处理方法,如果算术运算结果超出寄存器所能表示的最大数就会出现溢出因为16位定点DSP的动
3、态范围有限，所以在使用时必须注意动态范围以防溢出溢出还与输入信号的特性和运算法则有关,1.溢出,C55x有以下几种硬件特性可以处理溢出：保护位 C55x的每个累加器都有8个保护位（3932位），允许连续256次乘加操作而累加器不溢出溢出标志位 C55x的每个累加器都有相关的溢出标志位，当累加器操作结果出现溢出时，这个标志位就会置位,2.C55x的溢出处理机制,饱和方式位SATD和SATASATD控制D单元的操作，SATA控制A单元的操作。如果SATD=1，当D单元发生溢出时，对D单元的结果进行饱和处理。不管饱和方式位的值是什么，当累加器发生溢出时，相应的溢出标志位都会被置位A单元没有溢出标志位
4、，但如果SATA=1，发生溢出时，结果也会进行饱和处理,饱和处理是用最近的边界值代替溢出结果。例如，16位寄存器的范围是8000h（最小负数）7FFFh（最大正数），饱和处理就是用7FFFh代替比7FFFh大的结果；用8000h代替比8000h小的结果。,饱和。饱和是一种处理溢出的方法，但是饱和会剪掉部分输出信号，可能会引起信号失真和引起系统非线性。输入定标。分析所要使用的系统，假定最坏的情况，然后对输入信号定标，以防止溢出。但是这种方法会极大地降低输出信号的精确度。固定定标。假定最坏的情况，对中间结果定标。这种方法可以防止溢出，同时增加了系统的信噪比。动态定标。可以监测中间结果的范围，只在需
5、要的时候对中间结果定标。这种方法可以防止溢出但会增加计算量。,3.溢出的处理方法,7.1.3 常用信号处理 算法中的定标方法,FIR滤波器的定标方法在FIR滤波器中处理溢出的最好方法是设计时使滤波器的增益小于1，这样就不需要对输入信号定标。这种方法和累加器的保护位结合起来，可以有效地防止溢出。由于对信号处理的负面影响，在FIR滤波器中不使用固定定标和输入定标。如果不考虑计算量，在FIR滤波器中可以使用动态定标。对一些类型的音频信号，饱和处理也是一种常用的方法。,IIR滤波器的定标方法IIR滤波器的定点实现推荐使用多个二阶基本节级联组成，这样可以减小高阶滤波器频率响应灵敏度。由于滤波器系数的量化
6、引入误差，因此避免溢出对IIR滤波器非常重要。可以通过把中间结果保存在处理器累加器来避免节间数据溢出。为防止在第k阶内部发生数据溢出，需要用增益系数对滤波器的单位脉冲响应（前馈通道）定标。动态标定方法。在每个阶段滤波器内部状态都被减半，以提高指令周期换取为代价提高了结果的精度。,FFT的定标方法,在FFT操作里，每次蝶形运算后数据平均增加一位。输入定标需要移位（FFT长度为N），这会导致在计算FFT之前就衰减 6 dB。在固定定标中，每级蝶形运算输出除以2，这是最常用的FFT定标方法，因为它简单而且有比较好的信噪比。但是，对于大的FFT，这种定标可能会使信息丢失。另一种方法是动态定标，即在输出
7、溢出时再除以2。在这种情况下，会在这个过程中指定一个变量，每定标一次变量的值加1，计算结束后根据变量的值把结果乘以一个系数。动态定标的信噪比最好，但会增加FFT循环次数。,7.3 FIR滤波器,数字滤波器是DSP的基本应用，有2种基本类型：有限冲激响应滤波器FIR无限冲激滤波器IIR一般来说，如果需要线性相位则选择用FIR滤波器，对于相位要求不敏感的场合可以选用IIR滤波器。本节主要讨论FIR滤波器的DSP实现方法，有关IIR滤波器的实现将在下一节中介绍。,7.3.1 FIR滤波器的基本结构,一个FIR滤波器的输出序列和输入序列之间的关系，满足差分方程：传递函数为,FIR滤波器的结构：FIR滤
8、波器的单位冲激响应是一个有限长序列。若为实数，且满足偶对称或奇对称的条件，则FIR滤波器具有线性相位特性。偶对称线性相位FIR滤波器的差分方程为：,7.3.2 FIR滤波器的C语言编程实现,例7-9，直接型FIR滤波器的C语言编程实现。,/*fir.c 该程序用于实现FIR滤波器*L滤波器的阶数*bi滤波器的系数,i=0,1,L-1*xi输入信号向量,i=0,1,L-1;x0对应于当前值,x1对应于上一采样值*x_in 输入信号的当前值*y_out 输出信号的当前值*/,float fir(float x_in,float*x,float*b,int L)float y_out;int i;/
9、*-*/*把上一个采样时间的输入信号向量延迟一个单元,得到当前采样时间的输入信号向量*/for(i=L-1;i0;i-)xi=xi-1;x0=x_in;/*-*/*完成FIR滤波*/y_out=0.0;for(i=0;iL;i+)y_out=y_out+bi*xi;return y_out;,直接型FIR滤波器的实现涉及到两个基本操作，一个是输入信号向量与滤波器系数向量的内积计算，另一个是输入信号向量的更新处理。在每个采样周期信号缓冲器都要更新一次，最老的采样被抛弃，而其他的信号则向缓冲器的右方移动一个单元，一个新的采样被插入存储单元，并被标记。如果这个操作过程不用DSP硬件完成，那么它需要很
10、多的时间。,7.3.3 FIR滤波器的汇编语言编程实现,处理信号缓冲器的最有效方法，是把信号采样加载到循环缓冲器中。在循环缓冲器中，采取数据保持固定、反时针方向移动地址的方式，代替保持缓冲器地址固定且正方向移动数据。信号采样的起点由指针x(n)指定，其它诸采样则沿着顺时针方向，从起点开始依次顺序加载。,图7-4 FIR滤波器的循环缓冲区,(a)信号循环缓冲区(b)系数循环缓冲区,例7-10，FIR滤波器的C55x汇编语言实现。,(1)主程序 fir_test.c/*fir_test.c*/#includemath.h#define L 64/*Number of FIR filter coef
11、ficients*/#define Fs 8000/*8000 Hz sampling frequency*/#define T 1/Fs#define f1 800/*800 Hz frequency*/#define f2 1800/*1800 Hz frequency*/#define f3 3300/*3300 Hz frequency*/#define PI 3.1415926#define w1(2*PI*f1*T)/*2*pi*f1/Fs*/#define w2(2*PI*f2*T)/*2*pi*f2/Fs*/#define w3(2*PI*f3*T)/*2*pi*f3/Fs*/
12、#define a1 0.333/*Magnitude for wave 1*/#define a2 0.333/*Magnitude for wave 2*/#define a3 0.333/*Magnitude for wave 3*/,extern int fir(int*,int*,unsigned int,int);/*Low-pass FIR filter coefficients*/int coeffL=-26,-13,14,36,31,-8,-58,-71,-15,83,139,76,-90,-231,-194,50,331,383,78,-405,-654,-347,403,
13、1024,863,-228,-1577,-1972,-453,2910,6836,9470,9470,6836,2910,-453,-1972,-1577,-228,863,1024,403,-347,-654,-405,78,383,331,50,-194,-231,-90,76,139,83,-15,-71,-58,-8,31,36,14,-13,-26;int inL;/*input buffer*/int outL;/*Output buffer*/,main()unsigned int i;float signal;unsigned int n=0;int index=0;for(i
14、=0;iL;i+)ini=0;outi=0;while(1)signal=a1*cos(float)w1*n);signal+=a2*cos(float)w2*n);signal+=a3*cos(float)w3*n);n+;inindex=(int)(0 x7fff*signal)+0.5);outindex=fir(in,coeff,L,index);index-;if(index=-1)index=L-1;,(2)汇编语言整数fir滤波器函数：fir.asm,;fir.asm 该程序用于实现FIR滤波器,可被C语言程序调用;int fir(int*,int*,unsigned int，i
15、nt);参数0:AR0 输入信号缓冲区指针;参数1:AR1-FIR滤波器系数向量指针;参数2:T0-FIR 滤波器的阶数L；参数3:T1-输入信号当前值在循环缓冲区的序数;返回值:T0-输出信号当前值.def _fir _fir pshm ST1_55;现场ST1,ST2和ST3入栈 pshm ST2_55 pshm ST3_55or#0 x340,mmap(ST1_55);设置FRCT,SXMD,SATDbset SMUL;置位SMUL,mov mmap(AR0),BSA01;AR0=输入信号循环缓冲区的起始地址mov mmap(AR1),BSA23;AR1=滤波器系数循环缓冲区的起始地址m
16、ov mmap(T0),BK03;设置循环缓冲区大小or#0 x5,mmap(ST2_55);AR0和AR2为循环缓冲区指针mov T1,AR0;AR0从index偏移量开始mov#0,AR2;AR2从0偏移量开始sub#2,T0;T0=L-2 mov T0,CSR;设置外部循环次数为L-1 mpym*AR0+,*AR2+,AC0;执行第一次运算|rpt CSR;启动循环macm*AR0+,*AR2+,AC0mov hi(AC0),T0;用Q15格式存放结果popm ST3_55;恢复ST1,ST2和 ST3popm ST2_55 popm ST1_55 ret.end,7.4 IIR滤波器,
17、IIR滤波器的优点：结构简单，运算量小，可以用较少的阶数获得很高的选择性。IIR滤波器的缺点：具有相位特性差，存在稳定性问题。高阶IIR滤波器经常以串联或并联二阶环节的形式予以实现。,7.4.1 二阶IIR滤波器的结构,二阶IIR滤波器，又称为二阶基本节，分为直接型、标准型和变换型。二阶IIR滤波器传递函数为：,图7-5 二阶IIR滤波器的直接型的实现,图7-6 的信号流图,图7-7 二阶IIR滤波器的直接型的实现,由于直接型对于给定的传递函数具有最小可能的延迟数、加法器数和乘法器数，所以被称为标准型。,7.4.2 高阶IIR滤波器的结构,高阶IIR滤波器的差分方程和系统函数分别为：,图7-8
18、 高阶IIR滤波器的直接型实现（L=M+1）,图7-9 高阶IIR滤波器的串联型结构,图7-10 第k个二阶节,对于，有,其中，,注意：,7.4.3 IIR滤波器的C语言实现,例7-11,采用二维数组编写的IIR滤波器C语言程序。,temp=xin;/*xin为IIR滤波器的输入*/for(k=0;kN_IIR;k+)/*N_IIR为IIR滤波器二阶节的个数w0k=temp-a1k*w1k-a2k*w2k;/*这里,temp为本二阶节的输入,也是上一级二阶节的输出*/temp=b0k*w0k+b1k*w1k+b2k*w2k;/*这里,temp为本二阶节的输出,也是下一级二阶节的输入*/w2k=
19、w1k;w1k=w0k;xoutput=temp;/*xoutput为IIR滤波器的输出*/,aik（i=1,2）和bik（i=0,1,2）分别为第k（k=0,1,K）个二阶IIR滤波器环节的系数。,wik(i=0,1,2)是第k个二阶IIR滤波器环节中对应于图7-10中(m=0,1,2)的中间信号。,7.5 快速傅里叶变换FFT,FFT算法原理库利一图基算法FFT算法的实现,7.5.1 FFT算法原理,快速傅里叶变换（FFT）是离散傅里叶变换（DFT）的一种快速算法。通过FFT算法，DFT的计算量大大减少，运算时间缩短12个数量级。DFT的变换公式为 其中 为旋转因子。,FFT之所以减少运算
20、量，主要是利用了旋转因子的以下3点特性：对称性周期性可约性,利用这些特性可以使DFT运算中有些项进行合并，将长序列的DFT分解为短序列的DFT。DFT从算法上分为按时间抽选（DIT）和按频率抽选(DIF）。基2的DIT又被称为库利一图基算法。基2的DIF又称为桑德图基算法。,7.5.2 库利一图基算法,信号流图比特反转蝶形运算,1.信号流图,2比特反转,图7-11的输入信号的顺序是按照比特反转排列的，输出序列是按照自然顺序的。比特反转就是将序列下标用二进制表示，然后将二进制数按照相反的方向排列，即得到这个序列的实际位置。按照自然排序的时域信号数据是x(0)、x(1)、x(2)、x(3)、x(4
21、)、x(5)、x(6)、x(7)，其序号写成二进制数分别为000b、001b、010b、011b、100b、101b、110b、111b，将这些二进制数前后倒转，即得到进行FFT前数据所对应的实际二进制数地址：000b、100b、010b、110b、001b、101b、011b、111b，对应的十进制数是：0、4、2、6、1、5、3、7。序号为3的存储单元，按照自然排序应该存放x(3)，但由于FFT计算规则的要求，现在应该存放x(6)。,3.蝶形运算,基2DIT FFT算法，共由M级构成，每级计算由N/2个蝶形运算构成。基本运算单元为以下所谓蝶形运算：蝶形运算中上下两个节点p、q的间距为：,7
22、.5.3 FFT算法的实现,为了叙述简单，本书给出采用C语言编写的FFT程序，相应的汇编程序请读者自行完成。,（1）主程序：fft_test.c#include#include fcomplex.h/*包含浮点复数结构体定义头文件fcomplex.h*/extern void bit_rev(complex*,unsigned int);/*位反转函数声明*/extern void fft(complex*,unsigned int,complex*,unsigned int);/*fft函数声明*/,例7-13，基2DIT FFT算法的C语言实现。,extern void generator
23、(int*,unsigned int)#define N 128/*FFT的数据个数*/#define M 7/*M=log2(N)*/#define PI 3.1415926 complex XN;/*说明输入信号数组,为复数*/complex WM;/*说明旋转因子数组e(-j2PI/N),为复数*/complex temp;/*说明临时复数变量*/float xinN;float spectrumN;/*说明功率谱信号数组,为实数*/float re1N,im1N;/*说明临时变量数组,为实数*/,void main()unsigned int i,j,L,LE,LE1;/*-*/*产生
24、旋转因子表*/for(L=1;L1;/*子FFT中的蝶形运算数目*/WL-1.re=cos(PI/LE1);WL-1.im=-sin(PI/LE1);/*-*/generator(xin,N);,for(;)for(i=0;iN;i+)/*构造输入信号样本*/Xi.re=xini;Xi.im=0;/*复制到参考缓冲器*/re1i=Xi.re;im1i=Xi.im;/*启动 FFT*/bit_rev(X,M);/*以倒位次序排列X*/fft(X,M,W,1);/*执行 FFT*/*计算功率谱，验证FFT结果*/for(i=0;iN;i+)temp.re=Xi.re*Xi.re;temp.im=X
25、i.im*Xi.im;spectrumi=(temp.re+temp.im)*4;,（2）浮点复数基2 DIT FFT函数：fft_float.c#include fcomplex.hvoid fft(complex*X,unsigned int M,complex*W,unsigned int SCALE)complex temp;/*复变量临时存储器*/complex U;/*旋转因子Wk*/unsigned int i,j;unsigned int id;/*蝶形运算中下位节点的序号*/unsigned int N=11;U.re=1.0;U.im=0.;,for(j=0;jLE1;j+
26、)for(i=j;iN;i+=LE)/*进行蝶形计算*/id=i+LE1;temp.re=(Xid.re*U.re-Xid.im*U.im)*scale;temp.im=(Xid.im*U.re+Xid.re*U.im)*scale;Xid.re=Xi.re*scale-temp.re;Xid.im=Xi.im*scale-temp.im;Xi.re=Xi.re*scale+temp.re;Xi.im=Xi.im*scale+temp.im;/*递推计算Wk*/temp.re=U.re*WL-1.re-U.im*WL-1.im;U.im=U.re*WL-1.im+U.im*WL-1.re;U.
27、re=temp.re;,（3）位反转函数：bit_rev.c#include fcomplex.hvoid bit_rev(complex*X,unsigned int M)complex temp;unsigned int i,j,k;unsigned int N=11;for(j=0,i=1;i=1;j+=k;if(ij)temp=Xj;Xj=Xi;Xi=temp;,（4）信号发生器函数：generator.c/*generator.c-该程序用于产生一组信号样本*/#include math.h#define Fs=8000;采样频率设为8000Hz#define T=1/Fs;采样时间为0.25ms#define f1=500;信号源频率1取为500Hz#define a1=0.5;信号源幅度1取为0.5#define w1=2*PI*f1*T void generator(int*x,unsigned int N)unsigned int i;for(i=0;iN;i+)Xi=a1*cos(float)w1*i);,（5）复数结构定义头文件：fcomplex.hstruct cmpxfloat re;float im;typedef struct cmpx complex;,
