FIR数字滤波器分布式算法的原理及FPGA实现

1、FIR数字滤波器分布式算法的原理及FPGA实现摘要：在利用FPGA实现数字信号处理方面，分布式算法发挥着关键作用，与传统的乘积-积结构相比，具有并行处理的高效性特点。详细研究了基于FPGA、采用分布式算法实现FIR数字滤波器的原理和方法，并通过Xilinx ISE在Modelsim下进行了仿真。     关键词：分布式算法 DALUT FPGA FIR数字滤波器正在迅速地代替传统的由R、L、C元件和运算放大器组成的模块滤波器并且日益成为DSP的一种主要处理环节。FPGA也在逐渐取代ASIC和PDSP，用作前端数字信号处理的运算（如：FIR滤波、CORDIC算法或F

2、FT）。乘累加运算是实现大多数DSP算法的重要途径，而分布式算法则能够大大提高乘累加运算的效能。1 传统的乘累加结构FIR数字滤波器基本理论FIR滤波器被称为有限长脉冲响应滤波器，与IIR数字滤波器相对应，它的单位脉冲响应h(n)只有有限个数据点。输入信号经过线性时不变系系统输出的过程是一个输入信号与单位脉冲响应进行线性卷积的过程，即：式中，x(n)是输入信号，y(n)是卷积输出，h(n)是系统的单位脉冲响应。可以看出，每次采样y(n)需要进行L次乘法和L-1次加法操作实现乘累加之和，其中L是滤波器单位脉冲响应h(n)的长度。可以发现，当L很大时，每计算一个点，则需要很长的延迟时间。2 乘累加

3、运算的位宽分配DSP算法最主要的就是进行乘累加运算。假设采样信号的位宽用N来表示，则N位与N位的乘累结果需要2N位的寄存器来保存；如果两个操作数都是有符号数，则乘积只有2N-1个有效位，因为产生了两个符号位。为了使累加器的结果不产生溢出，需要对累加器进行冗余设计，也就是说要在累加器2N的位宽上多设计出K位，累加器的长度M计算方式如下（L为滤波器的长度）：对于无符号数：M=2N+K=2N+log2 L对于有符号数：M=2N=K=2N+log2 L-13 乘累加运算的分布式算法原理分析得益于Xilinx FPGA查找表结构的潜能，分布式算法在滤波器设计方面显示出了很高的效率，自20世纪90年代初以

4、来越来越受到人们的重要。分布式算法是基于查找表的一种计算方法，在利用FPGA实现数字信号处理方面发挥着重要的作用，可以大大提高信号的处理效率。它主要应用于数字滤波、频率转换等数字信号处理的乘累加运算。分布式算法推导如下：设Ak是已知常数（如滤波器系数、FFT中的正弦/余弦基本函数等），xk(n)是变量，可以看作是n时刻的第k个采样输入数据，y(n)代表n时刻的系统响应。那么它们的内积为：其中，xk(n)变量可以写成下面的格式：式中，B为数据格式的字长,xkb是变量的二进制位，只有“0”和“1”两种状态。将（2）式代入（1）式得：4 FPGA实现过程中查找表的构造方法根据以上论述，括号中的每一乘

5、积项代表着输入变量的某一位与常量的二进制“与”操作，加号代表着算术和操作，指数因子对括号中的值加权。如果事先构造一个查找表，该表存储着括号中所有可能的组合值，就可以通过所有输入变量相对应位的组合向量（XNb,X(N-1)b,.x1b）对该表进行寻址，该查找表称为DALUT。DALUT的构造规则如表1所示。5 采用分布式算法实现FIR数字滤波器为了说明问题，以一个三个系数的FIR数字滤波器为例设计分布式算法，字宽也设置为三位。设FIR数字滤波器系数为：h(0)=5,h(1)=2,h(2)=3。在进行FPGA设计时，该表以组件Component形式构建，设置为ROM结构，提供输入寻址端口table

6、_in2.0,输出端口table_out3.0。FPGA算法的结构图如图2所示。    算法实现中的几个关键问题为：（1）采用状态机实现分布式算法的状态转移状态机的实现如图3所示，设置三个状态s0、s1、s2 。状态s0完成数据的装入，数据寄存器需要成对出现，一个完成数据的延迟，另一个完成数据的移位，并将状态转移到s1；状态s1完成查找表功能、数据移位和分布式算法的乘累加运算，数据移位一个数据宽带后将状态转移到s2;状态s2完成数据的输出，并将状态转移到s0。利用状态机可以条理清楚地简化计算过程，在算法实现时发挥着关键的作用。   

7、（2）系统时钟与数据输入时钟的关系根据上述的状态转移关系，可以得出：每输入一个数据，在下一次数据输入之前，需要在状态s1停留一个数据宽带（三位）的时钟时间，在s2停留一个时钟的数据输出时间。也就是说，系统时钟频率应是数据输入频率的5倍，即fclkock=5fxin。（3）分布式算法中的乘累加式公推导及核心代表实现设B是数据的字宽，Pn是分布式算法第n位的结果，则有：有了该关系式，就可以通过for.loop循环，使用一条语句完成分布式乘累加算法。具体如下：for n in 0 to B-1 loopP:=p/2+tableout(n)*2B-1;End loop;6 算法仿真验证与结论本文实现的

8、FIR滤波器在Xilinx的集成开发环境ISE下利用ModelSim进行了仿真。当输入数据为7，3，1.时，仿真输出依次为35,29,32,16.,与乘累加方式FIR滤波算法得出的结果完全一致。假设查找表和PDSP的通用乘法器延时时间相同，分布式算法的等待时间是Br,通用乘法器的等待时间是N1。可见，对于位宽较小的数据来说，分布式算法的执行速度远高于乘累加运算。可见，利用FPGA实现分布式计算大大提高了计算的速度，在高速信号处理中发挥着重要作用。适合单片机实时处理的简单FIR滤波器设计l 引 言        近年来，心脏病的发病率不

9、断上升，为了使病人能够随时随地得到诊治，科研人员研制了多种便携式移动心电监护设备。在研制便携式心电监护设备时，既要考虑使 用者携带及使用方便，又要保证采集到高质量的心电信号。解决这一矛盾的惟一途径就是用尽量优化的软、硬件资源来完成相应功能。心电信号的频率在0 05100 Hz之间，一个正常的心电图由P波、QRS波群、T波等组成。由频谱分析可知，QRS波群的中心频率在17 Hz带宽约为10 Hz，P波、T波等的频带都在此频带的低端以外1。心电信号采集过程中的噪声源主要有：50 Hz工频干扰及其各次谐波，呼吸运动和电极移动所产生的频率小于5 Hz的干扰，肌肉收缩产生的O100 Hz肌电(EMG)噪

10、声，用于信号处理的电子设备所产生的噪声等2，3。其中最明显的是50 Hz工频干扰及其各次谐波和由于呼吸运动所产生的频率小于5 Hz的干扰。要得到清晰、准确的心电信号，必须在不影响正常心电波形的基础上将这些干扰滤除。这就要求滤波器具有良好的幅频特性和严格的线性相位特性。本 文借助Matlab 704，采用窗函数法结合零、极点调整法1和试探法设计了一种滤除工频50 Hz干扰的FIR滤波器。 2基本原理 21 窗函数设计滤波器6，7原理简述        设滤波器的理想频率响应函数为Hd(ejw)，则其对应的单位脉冲响应为：  &

11、#160;     窗函数法设计滤波器的基本原理就是用有限长度单位脉冲响应序列h(n)去逼近hd(n)。通常hd(n)是无限长序列，且是非因果的，所以用窗函数w (n)将hd(n)进行截断并作加权处理得到： h(n)=hd(n)w(n) h(n)就是实际设计的FIR数字滤波器的单位脉冲响应序列，其频率响应函数为：        其中，N为所选窗函数w(n)的长度。同时，如果要求滤波器具有线性相位特性，则h(n)必须满足： h(n)=±h(N1n) 22方法改进及实现  

12、60;     设计中选用FIR I型滤波器8,9，即满足h(n)=+h(N1n)且N为奇数。选用Hamming窗对理想滤波器的单位脉冲响应进行截断，窗函数为： w(n)=054046 cos(2nn/N)，N=0，1，2，N1        按常规设计方法直接设计50 Hz FIR陷波器，要达到较好的滤波效果，滤波器的阶数较高且滤波器系数多为小数，不利于通用单片机进行实时处理。设定采样频率Fs=200 Hz，先用窗函数法设计一个低频特性较好的30 Hz低通滤波器，然后利用Matlab FDAT

13、OOL工具箱里的量化工具对滤波器的系数进行量化和调整。量化后滤波器的传递函数为： Hl(Z)=0125+025 z1+0375 z2+ 025 z3+0125 z4 3滤波器波形特性比较        由以上滤波器的幅频特性，结合零、极点调整法和试探法反复对滤波器系数进行调整，使其具有更好的低频特性和50 Hz衰减特性。系数调整后滤波器的传递函数为： H2(z)=025+z105 z2+z3025 z4        由滤波器幅频特性曲线可以看出，经系数调整后滤波器具有更好的

14、低频特性，3 dB截止频率为35 Hz，同时对50 Hz工频干扰的衰减达到38 dB。但是，系数调整也使滤波器的高频特性变差，在688 Hz处的增益超过+3 dB。为了滤除心电信号中的高频噪声，必须对50 Hz以上频段范围内的信号进行衰减。        同样，可以用窗函数法结合试探法设计了一个低频特性较好且具有良好高频衰减特性的40 Hz低通滤波器。利用Matlab FDATOOL工具箱里的量化工具对滤波器的系数进行量化和调整。量化后滤波器的传递函数为： H3(z)=025+05 z1+025 z-2   

15、60;     由幅频特性曲线可以看出，量化后滤波器的一3 dB截止频率为363 Hz且具有良好的低频和高频特性。综合以上两种滤波器的幅频特性，将这两种滤波器进行级联。级联后滤波器的传递函数为：           滤波器的幅、相频率特性由图4可以看出，级联后滤波器的一3 dB截止频率为281 Hz，在50 Hz频率处的衰减为-442 dB，高频最小衰减为-10 dB，且滤波器在通带内具有严格的线性相位特性。该滤波器的群延时为： t=(N-1)/2*TS=3*5=1

16、5 ms        由心电信号的功率谱可知，28 dB的截止频率能够满足心电监护的需要。 4滤波效果        将上述设计的滤波器用于实际心电信号滤波取得了良好的滤波效果。5 结 语        由滤波器的幅、相频率特性和实际滤波效果可以看出，该滤波器不但能够很好地滤除工频干扰而且具有良好的低频特性。同时，滤波器的系统函数中分子、分母的各 相系数均为整数，便于乘法运算。并且分母正好是2的4次幂，在进行除法运

17、算时可以直接由一条移位指令完成。这无疑给编程和运算都带来了极大的方便。        本滤波器结构简单、占用资源少、滤波效果好，非常适用于嵌入式便携设备中的信号滤波。同时，该滤波器的计方法也可以推广到其他相似性能指标的滤波器设计中在科研实践中一定会有较大的发展前途基于DA算法的FIR滤波器设计与实现1引言    在数字信号处理系统中，FIR数字滤波器多采用专用DSP芯片(如TMS320CXX系列)，这种基于DSP的处理系统存在很多优点，比如方案灵活、可操作性强、程序易于移植。但这种结构的滤波器多是根据FI

18、R 滤波器的数据移位相乘累加的算法编写相应软件，利用软、硬件相互结合完成滤波器的设计。由于软件运行时，指令都是串行执行的，这严重制约了系统的运行速率，不能满足高传输速率，大数据吞吐量的数字信号的实时性处理要求。而基于DA算法的FPGA滤波器则是一种采用纯硬件的方式实现FIR 滤波器的方式，这种方法突出的优点是运算速度快，特别适合在高速、实时、快变的数字信号处理要求。2DA算法的基本思想    第一个讨论分布式算法(distributed arithmetic,DA)的可以追溯到1973年Cr0i Croisier的一篇论文，并由Peled和Liu来推广。但是直到F

19、PGA出现以后，才在FPG计算乘积和中被广泛地应用。在FPGA 芯片设计中，分布式算法(distributed arithmetic,DA)是一种很重要的技术。它在乘积和的计算中已被广泛应用，为了理解DA算法设计原理，我们考虑有一个线性时不变网络的输出可以用下式表示： 假设系数cn是已知常数，xn是变量，在有符号DA系统中假设变量xn的表达式如下： 其中xbn表示 xn的第b位，即xn是x第n个样点。所以，内积y 可以表示为： 我们重新分别求和次序(这就是分布式算法的由来)，结果是：y=(c0xB-10+c1+cN-1xB-1N-1)2B-1+(c0xB-20+c

20、1+cN-1xB-2N-1)2B-2     (2-4).+(c0x00+c1x01+cN-1x0N-1)20    用更紧凑的形式表示为：     从(24)式可以发现，分布式算法是一种以实现乘加运算为目的的运算方法。它与传统算法实现乘加运算的不同在于执行部分积运算的先后顺序不同。分布式算法在实现乘加功能时，是通过将各输入数据的每一对应位产生的部分积预先进行相加形成相应的部分积，然后再对各个部分积累加而得到最终的结果，而传统算法是等到所有乘积已经产生之后再来相加完成乘加运算的。与传统

21、串行算法相比，分布式算法可极大地减少硬件电路的规模，提高电路的执行速度。3基于DA算法的FIR数字滤波器的硬件实现    由上面分析可以知道，对于任何一个线性时不变网络总可以将它转换成式(2-5)那样的形式，FIR滤波器是一个很典型的线性时不变网络，它的表达式是： 我们可以将它转化成下面这样的形式： 在本系统中，采用的是8 位输入，所以B = 8 则：     线性相位FIR滤波器满足系数对称条件，本系统设计的系统为16阶，那么它的系数关于h7偶对称，即有：hn= h16-1-n= h15-n 

22、  (3-4)    由式(3-3)和(3-4)可得到实现系统硬件框图如图3-1： 4系统的VHDL描述4.1顶层模块设计Library ieee;Use ieee.std_logic_1164.a11;Use ieee.std_logic_unsigned.all;Entity FIR isPort (x: in std_logic_vector(7 downto 0);clk : in std_logic;y : out std_logic_vector(7 downto 0);end FIR;architecture behave of

23、FIR iscomponent lpfirport (in : in std_logic_vector(7 downto 0);clk : in std_logic;out : out std_logic_vector(7 downto 0);end component ;beginprocess(clk)beginu1：lpfir port map (x,clk,y);end process;end behave;4.2 LUT查询表ROM的建立用MIF文件格式编辑的ROM初始化值WIDTH = 16DEPTH = 16ADDRESS_RADIX = HEX ;DATA_RADIX = HEX ;CONTENT BEGIN0 : 0000;1 : 0045;2 : 00E6;F : 0000;END;4.3 滤波器系数的设计    我们利用matlab 自带的滤波器设计工具FDAtool和信号处理工具Pstool设计出一个16阶窗函数FIR低通滤波器