32阶FIR滤波器的FPGA实现

范晓东 蔡德林 桂 岳 梁本仁

摘 要:阐述了有限冲击响应(FIR)低通滤波器的窗函数设计方法,利用并行分布式算法在现场可编程门阵列上实现了32阶FIR低通滤波器。采用Altera公司中Stratix系列芯片内部的ROM实现了一种基于查找表结构的FIR数字滤波器,从而将卷积运算变换成一种查表后的加法运算,提高了运算速度,节省了逻辑单元。仿真结果表面,基于并行分布式算法的FIR滤波器最大处理速度可以达到200 MHz。

关键词:FIR滤波器;FPGA;并行分布式算法;Matlab;QuartusⅡ

中图分类号:TN713 文献标识码:A

文章编号:1004-373X(2009)21-186-03

Implement of 32 Orders FIR Filter on FPGA

FAN Xiaodong1,CAI Delin2,GUI Yue1,LIANG Benren1

(1.Electronic Science and Technology Institute,Anhui University,Hefei,230039,China;

2.No.38 Institute,China Electronics Technology Group Corporation,Hefei,230031,China)

Abstract:Windows function design method of FIR digital filter is introduced,and the parallel distributed arithmetic is used to implement a 32 orders FIR digital filter.Using the Altera′s Stratix series FPGA to design a new structure of FIR.Using the ROM in FPGA,this design can convert convolution to summation.Thus,a high process speed is improved and the Logic Elements(LE) is saved.The result of simulation shows that the max speed of FIR filter can arrive at 200MHz based on parallel distributed arithmetic.

Keywords:FIR filter;FPGA;parallel distributed arithmetic;Matlab;QuartusⅡ

随着软件无线电的发展,对于滤波器的处理速度要求越来越高。传统的FIR滤波器一般采用通用DSP处理器,但是DSP处理器采用的是串行运算,而FPGA是现场可编程阵列,可以实现专用集成电路,另外还可以采用纯并行结构及考虑流水线结构,因此在处理速度上可以明显高于DSP处理器。本文采用并行分布式算法在FPGA上设计并实现了高速处理的32阶FIR低通滤波器[1],在此过程中利用Matlab的数值计算与分析功能来提高设计效率。

1 FIR低通滤波器的窗函数实现

理想的滤波器频率响应中傅里叶反变换hd(n)一定是无限长的序列,而且是非因果的,而实际要设计的滤波器h(n)是有限长的,因此要用有限长来逼近无限长的,其方法就是用一个有限长度的窗口函数序列w(n)来截取,即:

h(n)=w(n)hd(n)

(1)

常见的窗函数有矩形窗、巴特利特窗、汉宁窗、哈明窗、布莱克曼窗、凯泽窗。其中,凯泽窗提供了可变的过渡带宽。本文采用凯泽窗对FIR滤波器进行设计,其窗函数表达式为:

w(n)=I0β1-1-2nM-12I0[β],0≤n≤M-1

(2)

I0[•]为第一类变形零阶贝赛尔函数,形状参数β为依赖于滤波器阶数M的参数,用来调整主瓣宽度与旁瓣衰减,选择M可产生各种过渡带宽和接近最优的阻带衰减。

给定通带截止频率ωp,阻带起始频率ωs,阻带衰减As,凯泽窗设计中有经典公式[2]可供使用,如下:

过渡带宽:

Δω=ωs-ωp

(3)

滤波器阶数:

M霢s-7.952.286Δω

(4)

形状参数:

β=0.110 2(As-8.7), As≥50

0.584 2(As-21)0.4+0.078 86(As-21),As<50

0,As≤21

(5)

对于低通滤波器:

hd(n)=ωpπsin[ωp(n-N-12)]ωp(n-N-12), 0≤n≤N-1

0,n为其他值

(6)

假设低通数字滤波器设计指标如下:

ωp=0.2π;ωs=0.4π;As=50 dB

采用上面介绍的凯泽窗,利用Matlab编程[4]计算得到32阶FIR低通滤波器参数如下:

h(0)=h(31)=0.001 0

h(1)=h(30)=0.001 9

h(2)=h(29)=0.000 5

h(3)=h(28)=-0.003 8

h(4)=h(27)=-0.007 6

h(5)=h(26)=-0.004 9

h(6)=h(25)=0.006 7

h(7)=h(24)=0.019 3

h(8)=h(23)=0.018 3

h(9)=h(22)=-0.005 3

h(10)=h(21)=-0.039 8

h(11)=h(20)=-0.053 1

h(12)=h(19)=-0.012 8

h(13)=h(18)=0.085 4

h(14)=h(17)=0.205 7

h(15)=h(16)=0.288 4

32阶FIR低通滤波器幅频特性图如图1所示。

图1 低通FIR滤波器的幅频特性

上述求得的系数是浮点型的,而在FPGA设计中使用的数据是定点型的,所以在设计滤波器之前要将系数转化为定点型,即系数的量化。在本文中采用数字信号处理(DSP) 技术中的Q值法[5]对系数进行量化。为了兼顾精度和所占用的资源,本文的系数用12位二进制来量化,得到的整数系数结果如下:

h(0)=h(31)= 2 h(1)=h(30)=4

h(2)=h(29)=1h(3)=h(28)=-8

h(4)=h(27)=-16 h(5)=h(26)=-10

h(6)=h(25)=14 h(7)=h(24)=40

h(8)=h(23)=37 h(9)=h(22)=-11

h(10)=h(21)=-81 h(11)=h(20)=-109

h(12)=h(19)=-26 h(13)=h(18)=175

h(14)=h(17)=421 h(15)=h(16)=591

2 并行分布式算法原理及FPGA设计

32阶FIR滤波器的差分方程表达式为:

y(n)=∑31m=0x(n-m)h(m)

(7)

式中:x(n)为输入;y(n)为输出;h(n)为滤波器系数。

设x(n)用二进制可表示为:

x(n)=x0(n)+21x1(n)+22x2(n)+…+210x10(n)-211x11(n)

(8)

其中,最高位为符号位。则式(7)可写为:

y(31)=h(0)x(31)+h(1)x(30)+…+

h(30)x(1)+h(31)x(0)

=h(0)[x0(31)+21x1(31)+…+

210x2(31)-211x11(31)]+

h(1)[x0(30)+21x1(30)+…+

210x2(30)-211x11(30)]+…+

h(31)[x0(0)+21x1(0)+…+

210x2(0)-211x11(0)]

(9)

转换得到:

y(31)=[h(0)x0(31)+h(1)x0(30)+…+

h(30)x0(1)+h(31)x0(0)]+

[h(0)x1(31)+h(1)x1(30)+…+

h(30)x1(1)+h(31)x1(0)]21+…+

[h(0)x10(31)+h(1)x10(30)+…+

h(30)x10(1)+h(31)x10(0)]210-

[h(0)x11(31)+h(1)x11(30)+…+

h(30)x11(1)+h(31)x11(0)]211

(10)

式(10)为并行分布式算法[3],由上可以看出并行分布式算法是将滤波器表达式重新排列,分别加权求和。与传统算法最大的不同之处是在FPGA设计过程中以查找表[6]代替乘法器,即根据输入数据的不同,将对应的滤波器系数预先求和保存在ROM中,也就是将每一项的乘法求和通过并行结构查表寻值完成,提高运行速度。

具体FPGA实现时,首先将12位的输入数据并行输入到12列32位移位寄存器分别寄存,然后以寄存器中的值为地址,对应于查找表的结果,按照式(10),每列进行相应二次幂加权,最后各列累加,在第32个数据完全输入之后得到正确的滤波器输出。由于输入数据的延迟,在此之前滤波器输出会延迟或者产生不正确的结果,可以在实现过程中加入控制信号进行输出控制。

由于查找表的规模是随着地址的增加呈指数增加的,可以将32位的查找表划分为四个8位的查找表,从而降低对ROM的需求[7]。

在本设计中可采用多级流水线技术[8],也就是将在明显制约系统速度的长路径上插入几级寄存器,虽然流水线会影响器件资源的使用量,但它降低了寄存器间的传播时延,允许维持高的系统时钟速率。

3 FPGA仿真与验证

由于直接将大量数据进行硬件仿真验证很不方便,因此利用Matlab产生一个采样频率为100 MHz,频率分别为1 MHz与30 MHz的两个正弦信号相加后,作为输入信号。同样,浮点变为定点,将此信号进行12位量化,并将负数转化为补码形式,按照一定格式[9]保存为.vec文件,导入到QuartusⅡ中进行仿真,时序功能仿真结果如图2所示。

图2 FIR低通滤波器的时序功能仿真

其中,clk为时钟信号,x_in为滤波器输入信号,y为滤波器输出信号。图2并不能很直观地看出并行分布式算法产生的滤波效果,可以将QuartusⅡ中.vwf文件转化为.tbl文件[10],在Matlab中按照一定形式[9]编程可以得到时域及频域波形图,如图3,图4所示。

图3,图4中,软件仿真是直接在Matlab中用输入信号与滤波系数卷积得到的,在时域波形中软件仿真输出信号与理想信号相比有一定时间延迟,而QuartusⅡ仿真与软件仿真结果中幅度的差别是由于硬件输入量化产生的。

从时域或者频域波形图可以看出,频率为30 MHz的信号被滤除掉,只有频率为1 MHz的信号通过滤波器,达到了滤波的目的。

图3 信号时域波形图

图4 信号频域波形图

4 结 语

本设计选用Stratix系列芯片,最大处理速度可以达到200 MHz以上。本文没有考虑线性相位的滤波器对称性,在考虑线性相位的基础之上结合一些其他算法[4]可以降低器件数量和进一步提高处理速度。由于FPGA器件的可编程特性,在本设计中可以修改滤波器参数,得到高速处理的高通或者带通数字滤波器,具有一定实用价值。另外,本文利用QuartusⅡ与Matlab联合仿真,极大地提高了FPGA的设计效率。

参考文献

[1]Chi-Tsong Chen.Digital Signal Processing Spectral Computation and Filter Design[M].北京:电子工业出版社,2002.

[2]程佩青.数字信号处理教程[M].北京:清华大学出版社,2001.

[3]Baese U M.数字信号处理的FPGA实现[M].刘凌,译.北京:清华大学出版社,2003.

[4]周伟林.FIR 滤波器的软件仿真与硬件实现[J].微计算机信息,2009,25(4):222-224.