浙江大学学报(工学版), 2021, 55(6): 1199-1207 doi: 10.3785/j.issn.1008-973X.2021.06.021

信息与电子工程

基于改进旋转因子的高性能FFT硬件设计

骆阳,, 张为,

天津大学 微电子学院,天津 300072

Hardware efficient FFT design based on improved rotation factor

LUO Yang,, ZHANG Wei,

School of Microelectronic, Tianjin University, Tianjin 300072, China

通讯作者: 张为,男,教授. orcid.org/0000-0002-2601-3198. E-mail: tjuzhangwei@tju.edu.cn

收稿日期: 2020-06-7  

基金资助: 光电信息控制和安全技术重点实验室资助项目(JCKY2019210C053)

Received: 2020-06-7  

Fund supported: 光电信息控制和安全技术重点实验室资助项目(JCKY2019210C053)

作者简介 About authors

骆阳(1995—),女,硕士生,从事数字集成电路的研究.orcid.org/0000-0002-0691-9478.E-mail:2018232032@tju.edu.cn , E-mail:2018232032@tju.edu.cn

摘要

针对FFT硬件实现中旋转因子模块占用资源较多的问题,设计高性能单路延时反馈结构的基22快速傅里叶变换. 采用CORDIC与MCM混合的方法设计旋转因子模块,实现了无需常规乘法器的FFT架构,不必占用DSP48E资源. 对于旋转角度数量较少的W16旋转因子模块,采用基于三输入加法器的MCM方法设计,将加法器数量降到最低. 对于旋转角度数量较多的W64W256W1 024模块,采用CORDIC方法设计. 依据旋转角度的数学规律,设计旋转角度实时生成模块,与传统的CORDIC方法相比,不需要占用ROM资源,避免了复杂的寻址逻辑和时序控制. 与其他构架相比,设计的16 bit 64点快速傅里叶变换在Xilinx Virtex-7上将单位slice吞吐率提高了35.20%,256点FFT在Virtex-5上提高了30.37%,1 024点FFT在Virtex-7上提高了25.38%.

关键词: 快速傅里叶变换 ; 单路延迟反馈架构 ; 常数乘法器 ; 坐标旋转数字计算方法

Abstract

A hardware efficient Radix-22 fast Fourier transform based on a single-path delay feedback architecture was designed aiming at the problem that the rotation factor module takes up more resources in FFT hardware implementation. The method of mixing CORDIC and MCM was adopted to design the rotation factor module to realize FFT architecture without conventional multiplier and DSP48E resource. MCM method based on ternary adders was used to minimize the number of adders for the W16 rotation factor module with less rotation angles. CORDIC method was adopted for the rotation factor modules of W64, W256 and W1 024 with more rotation angles. The real-time generation module of rotation angles was designed according to the mathematical law of rotation angles. The method does not need to occupy ROM resources and avoids complex addressing logic and timing control compared with the traditional CORDIC method. The designed 16 bit 64-point FFT improves the throughput per slice by 35.20% on Xilinx Virtex-7, the 256-point FFT improves by 30.37% on Virtex-5, and the 1 024-point FFT improves by 25.38% on Virtex-7 compared with other architectures.

Keywords: fast Fourier transform ; single-path delay feedback architecture ; constant multiplier ; coordinate rotation digital computer method

PDF (1304KB) 元数据 多维度评价 相关文章 导出 EndNote| Ris| Bibtex  收藏本文

本文引用格式

骆阳, 张为. 基于改进旋转因子的高性能FFT硬件设计. 浙江大学学报(工学版)[J], 2021, 55(6): 1199-1207 doi:10.3785/j.issn.1008-973X.2021.06.021

LUO Yang, ZHANG Wei. Hardware efficient FFT design based on improved rotation factor. Journal of Zhejiang University(Engineering Science)[J], 2021, 55(6): 1199-1207 doi:10.3785/j.issn.1008-973X.2021.06.021

快速傅里叶变换(fast Fourier transform,FFT)作为离散傅里叶变换(discrete Fourier transformation, DFT)的快速算法,成为信号处理中最重要的算法之一,在通信、滤波及数字频谱分析等领域都有着广泛的应用. 为了满足数字信号处理的实时性需求,必须采用硬件电路,以提升计算速度. 设计高效的快速傅里叶变换的硬件结构具有重要意义.

基2单路延迟反馈架构(single-path delay feedback,SDF)是Groginsky等[1]提出的单路流水线FFT架构. He等[2]提出基22算法,结合基2与基4算法的优势,以资源占用少、控制简单、运算量小等优势,得到了广泛的应用. 近年来,Garrido等[3-4]对FFT的硬件实现进行了大量研究. 在蝶形单元架构的改进方面,Garrido等[5]提出串行架构,将蝶形单元数量减少一半,利用率提高到了100%,但计算速度降低了一半. Ingemarsson等[6]将传统SDF架构中的加法器和选择器换位,LUT资源占用减少了近1/3. Ingemarsson等[7]指出,若将小点数FFT架构中的一个加法器置换为短移位寄存器,则LUT资源占用减少了近43%. 目前对蝶形单元的改进已趋于完善,可以从占用资源很多的旋转因子乘法器模块继续进行改进,进一步改善FFT硬件架构的性能,研究方向主要有基于内存的架构和采用坐标旋转数字计算方法(coordinate rotation digital computer,CORDIC)的架构. Ma等[8]提出基于内存的FFT架构,将旋转因子存储于ROM中,用DSP完成复数乘法器运算. 这种基于内存的架构,对于大点数的FFT,需要大量的内存单元和DSP资源,设计更加复杂,硬件成本更高,因此采用CORDIC的架构逐渐得到更多的研究. Tang等[9]将旋转角度存储于ROM中,利用CORDIC计算旋转因子. Shi等[10]提前计算出每种旋转角度在各级旋转时的方向系数并存储在ROM中,计算时只需寻址读出相应的方向系数. Mankar等[11]提出基于CORDIC-新分布式算法(new distributed arithmetic,NEDA)的折叠FFT架构. Zhang等[12]提出自适应编码CORDIC(adaptive recoding CORDIC, ARC)FFT架构. Mahdavi等[13]提出基于二进制符号数编码(binary-signed-digit, BSD)CORDIC的FFT架构. 目前的研究主要是用CORDIC实现全部的旋转因子并着重于硬件实现方面的改进,鲜少将其他常数乘法器应用于FFT架构中.

针对FFT硬件实现中旋转因子乘法器模块占用资源较多的问题,在本文的设计中,结合常数乘法器和CORDIC算法,对于旋转角度较多的旋转因子模块采用CORDIC方法设计,且旋转角度依据数学规律实时生成,不再占用ROM资源,控制时序更加简单. 对于旋转角度较少的旋转因子模块,采用多元常数乘法(multiple constant multiplication,MCM),将加法器数量减少到最低. 对比其他现有架构,在保证计算精度的前提下,本文设计的基22SDF FFT架构使得硬件效率大大提升.

1. 基本原理

1.1. 基22 FFT算法的基本原理

xn)的离散傅里叶变换为

$X\left( k \right) = {\rm{DFT}}\left[ {x\left( n \right)} \right] = \sum\limits_{n = 0}^{N - 1} {x\left( n \right)W_N^{kn}} ;\;0 \leqslant k \leqslant N - 1.$

xn)和Xk)分解为4个子序列,令

$\left. { \begin{split} &n = \frac{N}{2}{n_1} + \frac{N}{4}{n_2} + {n_3};\;{n_1} = 0,1,\;{n_2} = 0,1,\;\\ &{n_3} = 0,1, \cdots ,{N}/{4} - 1. \end{split}} \right\} $

$\left. {\begin{split} &k = {k_1} + 2{k_2} + 4{k_3};\\ &{k_1} = 0,1,\;{k_2} = 0,1,\;{k_3} = 0,1, \cdots ,{N}/{4} - 1. \end{split}} \right\}$

将式(2)、(3)代入式(1),整理可得

$\begin{split} & X\left( {{k_1} + 2{k_2} + 4{k_3}} \right) =\\ & \sum\limits_{{n_3} \!=\! 0}^{\frac{N}{4}\! -\! 1} {\sum\limits_{{n_2} \!=\! 0}^1 {\sum\limits_{{n_1} \!= 0}^1 {x\left( {\frac{N}{2}{n_1}\! +\! \frac{N}{4}{n_2} \!+ \!{n_3}} \right)} } } W_N^{\left( {\frac{N}{2}{n_1} \!+ \!\frac{N}{4}{n_2}\! +\! {n_3}} \right)\left( {{k_1} + 2{k_2} + 4{k_3}} \right)}= \\ & \sum\limits_{{n_3} = 0}^{\frac{N}{4} - 1} {\left[ {H\left( {{k_1},{k_2},{n_3}} \right)W_N^{{n_3}\left( {{k_1} + 2{k_2}} \right)}} \right]W_{{N}/{4}}^{{n_3}{k_3}}} .\\[-23pt] \end{split} $

式中:

$H\left( {{k_1},{k_2},{n_3}} \right) = \left[ {x\left( {{n_3}} \right) + {{\left( { - 1} \right)}^{{k_1}}}x\left( {{n_3} + \frac{N}{2}} \right)} \right] + {\left( { - {\rm{j}}} \right)^{\left( {{k_1} + 2{k_2}} \right)}}\left[ {x\left( {{n_3} + \frac{N}{4}} \right) + {{\left( { - 1} \right)}^{{k_1}}}x\left( {{n_3} + \frac{3}{4}N} \right)} \right].$

1.2. CORDIC算法的基本原理

直角坐标系中的向量旋转可以表示为

$\begin{split} &\left[ {\begin{array}{*{20}{c}} {{x_{\rm{b}}}}\\ {{y_{\rm{b}}}} \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} {\cos \theta }&{ - \sin \theta }\\ {\sin \theta }&{\cos \theta } \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {{x_{\rm{a}}}}\\ {{y_{\rm{a}}}} \end{array}} \right] = \\ &\cos \theta \left[ {\begin{array}{*{20}{c}} 1&{ - \tan \theta }\\ {\tan \theta }&{1} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {{x_{\rm{a}}}}\\ {{y_{\rm{a}}}} \end{array}} \right]. \end{split}$

CORDIC算法通过一系列固定的、与运算基数有关的角度不断偏转,以逼近所需的旋转角度[14]. 只需要移位和加减运算,就可以完成矢量的旋转,即三角函数的乘法. 设总共需要n次旋转,其中第i次旋转角度为θi,且tanθi=si2isi取值为±1,+1表示逆时针旋转,−1表示顺时针旋转. 第i次旋转的表达式为

$\left[ {\begin{array}{*{20}{c}} {{x_{i + 1}}} \\ {{y_{i + 1}}} \end{array}} \right] = {m_i}\left[ {\begin{array}{*{20}{c}} 1&{ - {s_i}{2^{ - i}}} \\ {{s_i}{2^{ - i}}}&1 \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {{x_i}} \\ {{y_i}} \end{array}} \right].$

式中:mi为每次旋转的增益因子. 随着旋转次数n的增加,式(7)收敛为一个常数M,在n次旋转之后须乘以M校正.

$M = \prod\limits_{i = 0}^{n - 1} {{{{m}}_{{i}}}} = \prod\limits_{i = 0}^{n - 1} {{\rm{cos}}\;\theta } = \prod\limits_{i = 0}^{n - 1} {\frac{1}{{\sqrt {{\rm{1}} + {\rm{ta}}{{\rm{n}}^{\rm{2}}}\theta } }}} \approx 0.607\;252\;3{\rm{.}}$

1.3. MCM基本原理

在信号处理中,MCM可以将常数乘法分解为加减法和移位操作,用来实现一个变量与多个常数相乘. 移位操作可以通过硬件连接实现,资源消耗忽略不计,因此MCM的硬件资源消耗取决于加减法器的数量[15],RPAGT(reduced pipelined adder graph)算法是目前资源消耗较少的求解方法之一. 现代FPGA(Virtex 5-7、Spartan 6、Kintex 7和Artix 7)可以实现三输入加法器,使用三输入加法器的流水线MCM架构相比于传统二输入加法器实现的架构平均可以减少27%的运算量和15%的slice[16],因此将常数分解为三输入加法,可以进一步减小硬件资源占用和流水线级数. 算法流程[17]如下.

1  $Y = \{ {\rm{MSD}}(t)|t \in T\} $

2  $S = \mathop {\max }\limits_{t \in T} {\rm{AD}}_{\min }^3(t) = \mathop {\max }\limits_{y \in Y}\; \left[ {{{\log }_3}(nz(y))} \right]$

3  ${X_S} = \{ {\rm{odd}}(t)|t \in T\} \backslash \{ 0\} $

4  ${\rm{for }}\;s = S,S - 1, \cdots ,{\rm{2}}$

5   ${{W}} = {X_S}$

6   ${{P}} =\phi$

7  do

8    ${{p}}\! \leftarrow \!{\rm{single\_predecessor\; by \;formula \;11}}$

9   If p ≠ 0

10     ${{P}} \leftarrow {{P}} \cup \{ p\}$

11   else

12     ${{P'}} \leftarrow {\rm{multiple}}\_{\rm{predecessor}}$

13     ${{P}} \leftarrow {{P}} \cup p'$

14    ${{W}} \leftarrow {{W}}\backslash \{ w\}$

15   ${\rm{While}}\;|{{W}}| \ne \phi $

16   ${{{X}}_{{{s}} - {\rm{1}}}} \leftarrow {{P}}$

将目标系数集T用最小有符号数(minimum signed digit, MSD)方法表示,计算整体架构的流水线级数S. 将T化简为正奇数集并存入集合XS中. 将XS复制到当前工作集W中,并将前一级元素集P初始化为空集. 1个常数可以拆分为2个或3个常数的加减操作,如下所示:

$A_q^2({{u}},v) = |{2^{{l_u}}}u + {( - 1)^{{s_v}}}{2^{{l_v}}}v|{2^{ - r}}.$

$A_q^3({{u}},v,w) = |{2^{{l_u}}}u + {( - 1)^{{s_v}}}{2^{{l_v}}}v + {( - 1)^{{s_w}}}{2^{{l_w}}}w|{2^{ - r}}.$

式中:lu/v/wN0为左移位数,rN0为右移位数,sv/w∈{−1,1}为vw的符号. 利用下式计算W中所有元素的单一前级元素:

$\left. \begin{array}{l} {P_{\rm{a}}} = \{ w \in W|{\rm{AD}}_{\min }^3(w) < s\}; \\ {P_{\rm{b}}} = \{ w/({2^k} \pm 1)|w \in W,k \in {N_0}\} \cap N\backslash \{ w\}; \\ {P_{\rm{c}}} = \mathop \cup \limits_{\mathop {w \in W}\limits_{p' \in P} } \{ p = A_q^2(w,p')|q\;{\rm{valid}},{\rm{AD}}_{\min }^3(p) < s\} ;\\ {P_{\rm{d}}} = \{ w/({2^k} \pm {2^l} \pm 1)|w \in W,k,l \in {N_0}\} \cap N\backslash {P_{\rm{a}}};\\ {P_{\rm{e}}} = \mathop \cup \limits_{\mathop {w \in W}\limits_{p' \in P} } \{ p = A_q^2(w,p')/({2^k} \pm 1)|q\;{\rm{valid}},\\ \qquad {\rm{AD}}_{\min }^3(p) <s,k \in {N_0}\} \cap N\backslash {P_{\rm{a}}};\\ {P_{\rm{f}}} = \mathop \cup \limits_{\mathop {w \in W}\limits_{p',p'' \in P} } \{ p = A_q^3(w,p',p'')|q\;{\rm{valid}},\\ \qquad {\rm{AD}}_{\min }^3(p) < s,k \in {N_0}\} \backslash {P_{\rm{a}}}. \end{array} \right\}$

若单一前级元素不存在,则计算多元前级元素,将w用MSD方法表示,依据非零位将w分解为2个或3个元素之和,选择加法器深度最小且已求出前级元素p含量最多的分解方式. 将计算出的前级元素添加到集合P中,并将本级中对应的元素wW中移除,当W为空时本级计算结束;将级数s递减,重复上述过程,直到s = 2时终止,因为X2的前级元素p为{1}.

2. 基22 SDF流水线FFT结构

图1所示为1 024点基22 SDF FFT结构图. 该结构共有10级蝶型单元(BF)、10级移位寄存器(SR)、5个简单旋转因子乘法器模块(TR)和4个复杂旋转因子乘法器模块(WN).

图 1

图 1   1 024点基22 SDF FFT结构图

Fig.1   Radix 22 SDF FFT architecture for 1 024-Point


2.1. 控制模块

控制模块主要由5个10 bit计数器组成. 因为数据经过复杂旋转因子乘法器会产生延迟,W1 024W256W64需要24个时钟周期,W16需要3个时钟周期,所以计数器2~4比前一级计数推迟24个时钟周期,计数器5推迟3个时钟周期. 各模块的控制信号为相应计数器的某一位,如表1所示.

表 1   控制信号表

Tab.1  Table of control signal

计数器 计数延迟 计数器位 控制信号 控制模块
Counter_1 0clk Counter_1[9] Control_10 BF10、TR10
Counter_1 0clk Counter_1[8] Control_9 BF9、W1 024
Counter_2 24clk Counter_2[7] Control_8 BF8、TR8
Counter_2 24clk Counter_2[6] Control_7 BF7、W256
Counter_3 24clk Counter_3[5] Control_6 BF6、TR6
Counter_3 24clk Counter_3[4] Control_5 BF5、W64
Counter_4 24clk Counter_4[3] Control_4 BF4、TR4
Counter_4 24clk Counter_4[2] Control_3 BF3、W16
Counter_5 3clk Counter_5[1] Control_2 BF2、TR2
Counter_5 3clk Counter_5[0] Control_1 BF1

新窗口打开| 下载CSV


图2所示为1 024点FFT的时序图. 经过1 024个时钟周期完成数据的输入,在第513个时钟周期开始BF10的计算并输出计算结果,TR10的时序与BF10相同且没有数据延迟,在256个时钟周期后开始BF9模块的计算并输出计算结果,再经过24个时钟周期后数据完成W1 024模块的计算并输入下一级. 后续模块的流水线时序以此类推.

图 2

图 2   1 024点基22 SDF FFT时序图

Fig.2   Radix 22 SDF FFT sequence diagram for 1 024-point


2.2. 蝶形单元模块

传统的SDF蝶形结构如图3(a)所示. 当使用N输入的LUT进行2个数的加法运算时,只会占用2个输入端口,其他几个端口会闲置. 若将加法器输入端一侧的其他逻辑输入合并到该LUT的闲置端口,则可以更加充分地利用LUT资源[18]. 将传统蝶形结构中的加法器和选择器交换位置,使得选择器位于加法器输入端一侧,就可以将选择器和加法器综合进同一个 LUT单元中,从而节省硬件资源. 改进后的碟形单元结构如图3(b)所示[6]. 现代Xilinx FPGA具有6输入LUT,允许将2个长度相同且长度≤16的1 bit位宽短移位寄存器映射到同一个LUT6资源中,输出为1 bit位宽的加法器需要映射到一个完整的LUT6资源中. 这意味着短移位寄存器的LUT成本是加法器的一半. 对于所需移位寄存器长度≤16的情况,将传统蝶形结构中的一个加法器换成一个短移位寄存器,可以减少硬件资源占用. 改进后的碟形单元结构如图3(c)所示[7]. 基于上述2种思想,第6~10级蝶形单元采用图3(b)的改进结构,第1~5级蝶形单元采用图3(c)的改进结构,可以减少硬件资源占用,且这2种改进结构具有相同的数据顺序和控制结构,因此组合起来控制时序很简单. 分别以第6级和第2级蝶形单元为例说明工作时序,如图4所示.

图 3

图 3   蝶形单元结构图

Fig.3   Structure of butterfly elements


图 4

图 4   蝶形单元时序图

Fig.4   Butterfly element sequence diagram


2.3. 旋转因子乘法器模块

2.3.1. 简单旋转因子乘法器模块和复杂旋转因子乘法器模块W16

奇数级蝶形单元之前的简单旋转因子乘法器需要对数据进行乘−j操作:

$\left( {a + b{\rm{j}}} \right){\rm{\cdot}}\left( { - {\rm{j}}} \right) = b - a{\rm{j}}.$

由式(12)可知,只需将实部取反后再将实部、虚部换位即可完成乘法操作,结构如图5所示.

图 5

图 5   简单旋转因子(-j)乘法器

Fig.5   Multiplier of trivial rotation (-j)


W16的系数如表2所示. 当θ=0时,输入数据不需要经过乘法器,可以直接送入下一级;当θ=π/2时,处理方法与2.3.1节相同. 利用旋转因子的对称性,只需要设计出0.923 9、0.382 7、0.707 1这3种数值的常数乘法器. 整数的常数乘法器设计更加简单,因此可以先将其扩大28倍后取整,采用 1.3节的MCM方法,用移位操作和三输入加减法器来代替常规的乘法器,设计出236、97、181的常数乘法器,且最大限度地复用组成这3个常数的中间结果,从而减少硬件资源占用,提高计算速度. 这3个常数的MCM结构如图6所示. 如图7所示为W16旋转因子乘法器结构. 可知,输入数据的实部和虚部分别进入MCM模块后,得到扩大236、97、181倍后的3个值,再将其缩小28倍后得到乘0.923 9、0.382 7、0.707 1的结果. 通过乘法结果选择器中的Counter_4信号,从中选出相应的2个结果相加减,可得乘W16旋转因子的最终值. 该设计方法共使用8个加法器,且在加法器处插入寄存器,关键路径约为1个加法器,有效提高了FFT的系统性能.

表 2   W16的旋转因子系数

Tab.2  Coefficients of rotation W16

Counter_4 ${ {\rm{exp} }\;({ {\rm{ - j} }2{\text{π}} nk/N}) }$ n θ 实部 虚部
0、4、8、12~15 ${{\rm{e}}^0}$ 0 0 1 0
5 ${ {\rm{exp} }\;({ {\rm{ - j} }{\text{π}} k/8}) }$ 1 π/8 0.923 9 −0.382 7
1、6 ${ {\rm{exp} }\;({ {\rm{ - j} }{\text{π}} k/4} )}$ 2 π/4 0.707 1 −0.707 1
7、9 ${ {\rm{exp} }\;({ {\rm{ - j3} }{\text{π}} k/8}) }$ 3 3π/8 0.382 7 −0.923 9
2 ${ {\rm{exp} }\;({ {\rm{ - j} }{\text{π}} k/2}) }$ 4 π/2 0 −1
3、10 ${ {\rm{exp} }\;({ {\rm{ - j3} }{\text{π}} k/4}) }$ 6 3π/4 −0.707 1 −0.707 1
11 ${ {\rm{exp} }\;({ {\rm{ - j9} }{\text{π}} k/8}) }$ 9 9π/8 −0.923 9 0.382 7

新窗口打开| 下载CSV


图 6

图 6   236、97、181的MCM 流水线加法器图

Fig.6   Pipelined adder graph of MCM with coefficients {236、97、181}


图 7

图 7   W16旋转因子乘法器结构图

Fig.7   Multiplier of general rotation W16


2.3.2. 复杂旋转因子乘法器模块W64W256W1 024

 FFT中的旋转因子复数乘法可以表示为

$\left( {{x_{{\rm{re}}}} + {\rm{j}}{x_{{\rm{im}}}}} \right)W_N^n = \left( {{x_{{\rm{re}}}} + {\rm{j}}{x_{{\rm{im}}}}} \right){{\rm{exp}}\;\left( {{ - {\rm{j}}\frac{{2\pi n}}{N}}} \right)}.$

$\begin{array}{*{20}{l}} {\left[ {\begin{array}{*{20}{c}} {{x_{{\rm{re}} }^{'}}}\\ {{x_{{\rm{im}}}^{'}}} \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} {\cos \;\left( { - \dfrac{{2{\text{π}} n}}{N}} \right)}&{ - \sin\; \left( { - \dfrac{{2{\text{π}} n}}{N}} \right)}\\ {\sin \;\left( { - \dfrac{{2{\text{π}} n}}{N}} \right)}&{\cos\; \left( { - \dfrac{{2{\text{π}} n}}{N}} \right)} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {x_{{\rm{re}}}^{}}\\ {x_{{\rm{im}}}^{}} \end{array}} \right]{\rm{ = }}}\\ {\cos \left( { - \dfrac{{2{\text{π}} n}}{N}} \right)\left[ {\begin{array}{*{20}{c}} 1&{ - \tan \;\left( { - \dfrac{{2{\text{π}} n}}{N}} \right)}\\ {\tan \;\left( { - \dfrac{{2{\text{π}} n}}{N}} \right)}&1 \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {x_{{\rm{re}}}^{}}\\ {x_{{\rm{im}}}^{}} \end{array}} \right].} \end{array} $

将式(14)、(6)进行比较可知,两者具有相同的形式,因此将输入序列的实部和虚部代入X0Y0,将旋转角度代入Z0,可以通过CORDIC算法的移位和加减运算,完成输入序列与旋转因子的复数乘法功能. 20级流水线CORDIC算法的硬件结构如图8所示.

图 8

图 8   流水线型CORDIC算法结构图

Fig.8   Structure of pipelined CORDIC unit


CORDIC算法所能获得的最大旋转角度为 $ \mathop {{\rm{lim}}}\limits_{n \to \infty } \displaystyle\sum\limits_{{{i}} = 0}^n {{\rm{arctan}}\;{{\rm{2}}^{{{ - i}}}}} \approx {\rm{99}}^\circ $,故设定旋转区间为[0,π/2]. 对于其他区间的旋转角度,可以利用三角函数的性质,将旋转角度转换到第一象限内进行运算,如表3所示. 表中,Data_re、Data_im分别为信号的实部和虚部.

表 3   旋转角度预处理表格

Tab.3  Pretreatment of rotation angle

旋转角度Z0 角度预处理 X0 Y0
[0, π/2] Z Data_im Data_re
[π/2, π] Z − π/2 −Data_re Data_im
[π, 3π/2] Z − π −Data_im −Data_re

新窗口打开| 下载CSV


根据旋转因子角度的规律,用简单的累加器实时生成旋转角度,设计方法如图9所示.

图 9

图 9   旋转角度生成模块

Fig.9   Angle generator module


每一级旋转因子模块的角度如表4所示. 旋转角度都是2π/N的倍数,因此将2π/N作为累加器的初始角度,累加结果依次为0, 2π/N,···,2π/N×(N/4−1),最终的旋转角度在这一组值的基础上根据k1k2的不同分别扩大0倍、2倍、1倍和3倍.

表 4   旋转角度规律表

Tab.4  Mathematical law of rotation angles

k1 k2 n3 $W_N^{{n_3}\left( {{k_1} + 2{k_2}} \right)}$ 旋转角度
0 0 0, 1,···, N/4−1 $W_N^0$ 0
0 1 0, 1,···, N/4−1 $W_N^{2{n_3}}$ $0,\dfrac{ { {\rm{4} }{\text{π} } } }{ {{N} } },\dfrac{ { {\rm{8} }{\text{π} } } }{ {{N} } }, \ldots ,\dfrac{ { {\rm{4} }{\text{π} } } }{ {{N} } }\left( {\dfrac{ {{N} } }{ {\rm{4} } } - 1} \right)$
1 0 0, 1,···, N/4−1 $W_N^{{n_3}}$ $0,\dfrac{ { {\rm{2} }{\text{π} } } }{ {{N} } },\dfrac{ {4{\text{π} } } }{ {{N} } }, \ldots ,\dfrac{ { {\rm{2} }{\text{π} } } }{ {{N} } }\left( {\dfrac{ {{N} } }{ {\rm{4} } } - 1} \right)$
1 1 0, 1,···, N/4−1 $W_N^{3{n_3}}$ $0,\dfrac{ { {\rm{6} }{\text{π} } } }{ {{N} } },\dfrac{ { {\rm{12} }{\text{π} } } }{ {{N} } }, \ldots ,\dfrac{ { {\rm{6} }{\text{π} } } }{ {{N} } }\left( {\dfrac{ {{N} } }{ {\rm{4} } } - 1} \right)$

新窗口打开| 下载CSV


对CORDIC计算的结果乘增益因子M进行校正,M采用单常数乘法器(single constant multiplication,SCM)的方式实现,如图10所示.

图 10

图 10   增益因子M的SCM流水线加法器图

Fig.10   Pipelined adder graph of SCM with coefficient M


3. 实验结果对比分析

3.1. 功能仿真

测试数据形式如式(15)所示,采用常见的直流分量与余弦分量叠加的信号形式作为输入数据[19]. 如图11所示为Matlab和本文的FFT处理器的结果,测试数据如表5所示. 将FFT处理器的结果与信号的频域特性对比可知,相对误差均小于0.2%. 设计的FFT处理器可以正确地提取信号的频域特性,实现FFT的功能.

图 11

图 11   直流分量与余弦分量叠加信号的FFT结果图

Fig.11   FFT result of signal superimposed by DC component and cosine component


表 5   信号频域特性测试数据表

Tab.5  Table of signal frequency domain characteristic test data

FFT结果 结果幅值 信号幅值 相对误差/% 相位 相对误差/%
第0点 (64, 0) 64 0.062 5 0
第100点 (55.418, −32.031) 64.008 9 0.125 02 0.013 9 −30.027 4 0.091 5
第300点 (0.031, 127.949) 127.949 0 0.249 9 0.039 8 89.986 1 0.015 4
第724点 (0.059, −127.910) 127.910 0 0.249 8 0.070 3 −89.973 6 0.02936
第924点 (55.523, 32.082) 64.125 3 0.125 24 0.195 8 30.019 9 0.066 5

新窗口打开| 下载CSV


$\begin{split} {{S}} =& {\rm{0}}{\rm{.062\; 5}} + 0.125 {\rm{cos}}\left( {{\rm{2}} {\text{π}} \times 100 {{t - 30{\text{π}}}} /{\rm{18}}0} \right){\rm{ + }} \\ & {\rm{ 0}}{\rm{.25}} {\rm{cos}}\left( {{\rm{2}} {\text{π}} \times 300 {{t}} + 90{\text{π}}/{\rm{18}}0} \right). \end{split}\!\!\!\!\!\!\!\! $

3.2. 电路性能分析

本文实现了数据位宽为16 bit的64点、256点和1 024点基22 SDF FFT架构,得到在Virtex-5、Virtex-7 FPGA芯片上的电路性能. 该设计没有使用DSP48E资源,因此在比较硬件性能时,采用文献[20]的方法,将slice regs、LUTs、DSP48E都转换为slice数量进行统一度量. 在Xilinx Virtex-5和Virtex-7 FPGA中,每个DSP48E相当于占用大约500个slice. 电路性能和硬件资源的对比情况如表6所示. 表中,N 为计算点数,f为频率,Tp为吞吐率,Tp,u为单位slice吞吐率. Valencia等[21]将旋转因子储存到ROM中,通过状态机控制读取旋转因子的地址,实现了紧凑的基2/4 FFT架构. 本文设计的64点、1 024点FFT架构吞吐率虽然仅提升10%,但因为硬件资源占用更少,单位面积吞吐率分别提高了35.20%和66.29%. Wang等[22]提出组合SDC-SDF FFT架构,通过共享算术资源,减少了近50%的复数乘法器,但需要使用DSP48E资源,提高了控制时序的复杂度;本文设计的256点FFT单位面积吞吐率提高了30.37%. Nguyen等[23]使用优化的CORDIC算法,通过设置误差阈值来减少CORDIC模块的迭代次数,但吞吐率较低;本文设计的单位面积吞吐率提高了115%. Nguyen等[20]提出基于CORDIC算法的并行双路延迟换向结构FFT处理器,提升了吞吐率且不占用BRAM资源,但整序换向模块会占用较多的LUT资源;本文设计的FFT架构单位面积吞吐率提高了25.38%.

表 6   FPGA实现的电路性能和硬件资源比较表

Tab.6  FPGA implementation results of circuit performance and hardware resources for different architectures

方法 FPGA型号 N Slices LUTs Regs DSP48Es slice总量 f /MHz Tp /(MS·s−1 Tp,u /(kS·s−1
实际数 等效slice数
文献[21]方法 V7 64 848 626 8 4 000 5 474 335 335 61.20
本文方法 V7 64 667 2 423 1 924 0 0 4 347 359.69 359.69 82.74
文献[22]方法 V5 256 1 733 1 073 12 6 000 8 806 297 297 33.72
本文方法 V5 256 1 093 3 925 3 258 0 0 7 183 315.77 315.77 43.96
文献[23]方法 V7 1 024 11 865 1 393 0 0 13 258 200 200 15.08
文献[21]方法 V7 1 024 1 671 2 065 25 12 500 16 236 317 317 19.52
文献[20]方法 V7 1 024 12 737 2 715 0 0 15 452 200 400 25.89
本文方法 V7 1 024 1 645 6 098 4 785 0 0 10 883 353.21 353.21 32.46

新窗口打开| 下载CSV


3.3. 误差分析

FFT处理器的精度通过测量1 024个点的平均相对误差(mean relative error, MRE)来评估,如下所示:

${\rm{MRE}} = {{{N}^{-1}}}{{\displaystyle\sum\limits_{i = 1}^N {{{\left| {A(i) - B(i)} \right|}}/{{\left| {A(i)} \right|}}} }} \times 100{\text{%}} .$

式中:A为Matlab的FFT函数计算结果,B为FFT处理器的计算结果.

本文的FFT架构MRE为0.94%,误差来源主要是旋转因子乘法器模块和10级流水线运算的误差累积. 其中W1 024模块的MRE约为4.3×10−4W256模块约为9.5×10−5W64模块约为2.7×10−5,这3部分采用20级流水线CORDIC算法来实现,由于旋转角度由初始角度2π/N通过累加的方式实时生成,输入角度本身会产生误差. 20级运算的分辨率为10−4度,造成CORDIC算法的误差. W16模块采用MCM方法实现,先将旋转因子左移扩大后量化取整,完成乘法计算后再将结果进行截断取整,左移位数越多,误差越小,但资源消耗越多. 折中考虑这2个因素,选定扩大倍数为28,误差如表7所示.

表 7   W16旋转角度误差表

Tab.7  Rotation angle error of W16

理想系数 扩大取整后系数 实际旋转因子系数 相对误差/%
0.923 9 236 0.921 8 0.227
0.382 7 98 0.382 8 0.026
0.707 1 181 0.707 0 0.014

新窗口打开| 下载CSV


传统的基于内存的FFT处理器MRE为0.52%,文献[20]的MRE为0.81%,文献[23]的MRE为0.72%,文献[24] 的MRE为0.314%. 虽然本文的MRE略高于上述文献的架构,但相对误差小于1%,且硬件性能更具优势,以较少的硬件资源消耗得到了较好的计算结果.

4. 结 语

本文针对FFT硬件实现中旋转因子乘法器消耗硬件资源较多的问题,采用CORDIC与MCM混合的方法,实现了无需常规乘法器的FFT架构,不必占用DSP48E资源. 利用现代FPGA可以实现三输入加法器的特点,设计旋转角度数量较少的W16旋转因子模块,减少了硬件资源占用. 对于旋转角度数量较多的模块,采用流水线CORDIC架构实现,且旋转角度实时生成,避免了存储单元的使用. 蝶形单元中的存储模块均使用移位寄存器实现,整体架构只占用分布式逻辑资源,不需要BRAM资源. 本文提出的FFT架构可以应用于处理数据较多、计算较复杂的雷达成像处理领域,因为它可以把节省下来的硬件资源和ROM、DSP48E用于实现其他功能模块,以减少整体系统的资源消耗. 通过与Nguyen等[20-23]提出的FFT架构对比可知,利用本文的设计方法大大提高了硬件效率.

参考文献

GROGINSKY H L, WORKS G A

A pipeline fast Fourier transform

[J]. IEEE Transactions on Computers, 1970, 19 (11): 1015- 1019

[本文引用: 1]

HE S, TORKELSON M. A new approach to pipeline FFT processor[C]// International Parallel Processing Symposium. Hawaii: IEEE, 1996: 766-770.

[本文引用: 1]

GARRIDO M, QURESHI F, TAKALA J, et al. Hardware architectures for the fast Fourier transform[M]//SHUVRA S, LEUPERS R, TAKALA J, et al. Handbook of signal processing systems. Switzerland: Springer, 2018: 613-648.

[本文引用: 1]

QURESHI F. Optimization of rotations in FFTs[D]. Linköping: Linköping University, 2012.

[本文引用: 1]

GARRIDO M, HUANG S, CHEN S, et al

The serial commutator FFT

[J]. IEEE Transactions on Circuits and Systems II: Express Briefs, 2016, 63 (10): 974- 978

DOI:10.1109/TCSII.2016.2538119      [本文引用: 1]

INGEMARSSON C, KALLSTROM P, GUSTAFSSON O. Using DSP block pre-adders in pipeline SDF FFT implementations in contemporary FPGAs[C]// 22nd International Conference on Field Programmable Logic and Applications. Oslo: IEEE, 2012: 71-74.

[本文引用: 2]

INGEMARSSON C, GUSTAFSSON O

SFF: the single-stream FPGA-optimized feedforward FFT hardware architecture

[J]. Journal of Signal Processing Systems, 2018, 90 (11): 1583- 1592

DOI:10.1007/s11265-018-1370-y      [本文引用: 2]

MA Z G, YIN X B, YU F

A novel memory-based FFT architecture for real-valued signals based on a radix-2 decimation-in-frequency algorithm

[J]. IEEE Transactions on Circuits and Systems II: Express Briefs, 2015, 62 (9): 876- 880

DOI:10.1109/TCSII.2015.2435522      [本文引用: 1]

TANG A M, YU L, HAN F J, et al. CORDIC-based FFT real-time processing design and FPGA implementation[C]// 12th International Colloquium on Signal Processing and its Applications. Malacca: IEEE, 2016: 233-236.

[本文引用: 1]

SHI J Y, TIAN Y H, WANG M X, et al. A novel design of 1024-point pipelined FFT processor based on Cordic algorithm[C]// 2nd International Conference on Intelligent System Design and Engineering Application. Sanya: IEEE, 2012: 80-83.

[本文引用: 1]

MANKAR A, PRASAD N, DAS A D, et al

Multiplier: less VLSI architectures for radix‐22 folded pipelined complex FFT core

[J]. International Journal of Circuit Theory and Applications, 2015, 43 (11): 1743- 1758

DOI:10.1002/cta.2038      [本文引用: 1]

ZHANG J F, LIU H Z, CHEN T, et al

Enhanced hardware efficient FFT processor based on adaptive recoding CORDIC

[J]. Electronics and Electrical Engineering, 2013, 19 (4): 97- 103

[本文引用: 1]

MAHDAVI H, TIMARCHI S

Area-time-power efficient FFT architectures based on binary-signed-digit CORDIC

[J]. IEEE Transactions on Circuits and Systems I: Regular Papers, 2019, 66 (10): 3874- 3881

DOI:10.1109/TCSI.2019.2922988      [本文引用: 1]

MEYER-BASE U, MEYER-BASE A, HILBERG W. Coordinate rotation digital computer (CORDIC) synthesis for FPGA[C]// International Workshop on Field-programmable Logic. Berlin: Springer, 1994: 397–408.

[本文引用: 1]

VORONENKO Y, PUSCHEL M

Multiplierless multiple constant multiplication

[J]. ACM Transactions on Algorithms, 2007, 3 (2): 11- 49

DOI:10.1145/1240233.1240234      [本文引用: 1]

KUMM M, HARDIECK M, WILLKOMM J, et al. Multiple constant multiplication with ternary adders[C]// 23rd International Conference on Field Programmable Logic and Applications. Porto: IEEE, 2013: 1-8.

[本文引用: 1]

KUMM M. Multiple constant multiplication optimizations for programmable gate arrays [M]. Wiesbaden: Springer, 2016.

[本文引用: 1]

EHLIAR A. Optimizing Xilinx designs through primitive instantiation[C]// 7th FPGA World Conference. Copenhagen: ACM, 2010: 20-27.

[本文引用: 1]

MA Y K, LIANG H H. Implementation of a pipeline large-FFT processor based on the FPGA[C]// International Conference in Communications, Signal Processing and Systems. Harbin: Springer, 2017: 638-644.

[本文引用: 1]

NGUYEN H N, KHAN S A, KIM C H, et al

A high-performance, resource-efficient, reconfigurable parallel-pipelined FFT processor for FPGA platforms

[J]. Microprocessors and Microsystems, 2018, 60: 96- 106

DOI:10.1016/j.micpro.2018.04.003      [本文引用: 5]

VALENCIA D, ALIMOHAMMAD A

Compact and high-throughput parameterizable architectures for memory-based FFT algorithms

[J]. IET Circuits, Devices and Systems, 2019, 13 (5): 696- 703

DOI:10.1049/iet-cds.2018.5556      [本文引用: 3]

WANG Z, LIU X, HE B, et al

A combined SDC-SDF architecture for normal I/O pipelined radix-2 FFT

[J]. IEEE Transactions on Very Large Scale Integration Systems, 2015, 23 (5): 973- 977

DOI:10.1109/TVLSI.2014.2319335      [本文引用: 2]

NGUYEN H N, KHAN S A, KIM C H, et al

A pipelined FFT processor using an optimal hybrid rotation scheme for complex multiplication: design, FPGA implementation and analysis

[J]. Electronics, 2018, 7 (8): 137

DOI:10.3390/electronics7080137      [本文引用: 4]

VINAY K M, DAVID S A, SOBHA P M. Area and frequency optimized 1024 point radix-2 FFT processor on FPGA[C]// 2015 International Conference on VLSI Systems, Architecture, Technology and Applications. Bengaluru: IEEE, 2015: 1-6.

[本文引用: 1]

/