欢迎大家关注我的微信公众号：
原文链接：64点FFT处理器（上）
前言
截止2022年2月15日，中国科学院大学《高等数字集成电路分析及设计》课程终于完结，所以我计划分享几个自己完成的实践作业，供大家交流学习。
设计收获
- 对FFT/IFF算法有了清晰的理解
- 因为本设计为结课大作业，所以我进行了比较详细的文档介绍，并在源码中增加了自动化测试脚本，方便读者快速复现。
64点FFT处理器设计报告正文
一、设计内容
设计一个FFT处理器时序逻辑电路，计算64点FFT和IFFT（N = 64）。模块整体采用流水线结构实现，能够处理连续多组输入数据。
顶层模块名为fft_64，输入输出功能定义：
| 名称 | 方向 | 位宽 | 描述 |
|---|---|---|---|
| clk | I | 1 | 系统时钟 |
| rst_n | I | 1 | 系统异步复位，低电平有效 |
| inv | I | 1 | 模式控制，0表示FFT运算，1表示IFFT运算 |
| valid_in | I | 1 | 输入数据有效指示，高电平有效 |
| sop_in | I | 1 | 每组输入数据（64个数）第一个有效数据指示，高电平有效 |
| x_re | I | 16 | 输入数据实部，二进制补码定点格式 |
| x_im | I | 16 | 输入数据虚部，二进制补码定点格式 |
| valid_out | O | 1 | 输出数据有效指示，高电平有效 |
| sop_out | O | 1 | 每组输出数据（64个数）第一个有效数据指示，高电平有效 |
| y_re | O | 16 | 输出数据实部，二进制补码定点格式 |
| y_im | O | 16 | 输出数据虚部，二进制补码定点格式 |
设计要求：
- Verilog实现代码可综合，给出详细设计文档、综合以及仿真结果。
- 在每组数据处理过程中，inv信号值保持不变。
- 计算过程进行适当精度控制，保证输出结果精确度，输出定点格式（精度范围）可以根据需要进行调整，需要对计算结果进行误差分析。
二、FFT原理及设计
2.1 FFT算法原理
FFT（Fast Fourier Transform），快速傅立叶变换，是一种DFT（离散傅里叶变换）的高效算法。在以时频变换分析为基础的数字处理方法中，有着不可替代的作用。
设x(n)为N点有限长序列，则其DFT的运算公式为
X ( k ) = ∑ n = 1 N − 1 x ( n ) W N n k , k = 0 , 1 , . . . , N − 1 (2-1) X(k)=\sum_{n=1}^{N-1}x(n)W_N^{nk}, k=0,1,...,N-1 \tag{2-1} X(k)=n=1∑N−1x(n)WNnk,k=0,1,...,N−1(2-1)
其中：
W
N
=
e
−
j
2
π
/
N
W_N=e^{-j2π/N}
WN=e−j2π/N
将离散傅里叶变换公式拆分成奇偶项，则前 N/2 个点可以表示为：
X
(
k
)
=
∑
r
=
0
N
/
2
−
1
X
(
2
r
)
W
N
2
r
k
+
∑
r
=
0
N
/
2
−
1
X
(
2
r
+
1
)
W
N
(
2
r
+
1
)
k
=
∑
r
=
0
N
/
2
−
1
X
(
2
r
)
W
N
/
2
r
k
+
W
N
k
∑
r
=
0
N
/
2
−
1
X
(
2
r
+
1
)
W
N
/
2
r
k
=
A
(
k
)
+
W
N
k
B
(
k
)
,
k
=
0
,
1
,
…
,
N
/
2
−
1
X(k)=\sum_{r=0}^{N/2-1}X(2r)W_N^{2rk} +\sum_{r=0}^{N/2-1}X(2r+1)W_N^{(2r+1)k} \\ = \sum_{r=0}^{N/2-1}X(2r) W_{N/2}^{rk}+W_N^k \sum_{r=0}^{N/2-1}X(2r+1) W_{N/2}^{rk} \\ = A(k)+W_N^k B(k), k=0,1,…,N/2-1
X(k)=r=0∑N/2−1X(2r)WN2rk+r=0∑N/2−1X(2r+1)WN(2r+1)k=r=0∑N/2−1X(2r)WN/2rk+WNkr=0∑N/2−1X(2r+1)WN/2rk=A(k)+WNkB(k),k=0,1,…,N/2−1
同理，后N/2个点可以表示为：
X ( k + N / 2 ) = A ( k ) − W N k B ( k ) , k = 0 , 1 , … , N / 2 − 1 X(k+N/2)=A(k)-W_N^k B(k), k=0,1,…,N/2-1 X(k+N/2)=A(k)−WNkB(k),k=0,1,…,N/2−1
由此可知，后 N/2 个点的值完全可以通过计算前 N/2 个点时的中间过程值确定。对 A(k) 与 B(k) 继续进行奇偶分解，直至变成 2 点的 DFT，这样就可以避免很多的重复计算，实现了快速离散傅里叶变换（FFT）的过程。
考虑到64点数据量很大，不方便绘图解释原理，这里以8点FFT为例介绍基2FFT算法结构。
说明：图中 W ( r ) = W N r , r = 0 , 1 , 2 , 3 , 4 , 5 , 6 , 7 W(r)=W_N^r, r=0,1,2,3,4,5,6,7 W(r)=WNr,r=0,1,2,3,4,5,6,7
要介绍上面这一计算流图，首先要了解结构中的一些概念：
（1） 级的概念
如图2-1中所以，每一列即为一级运算，设FFT运算点数为N，则需要被分割成
M
=
l
o
g
2
N
M=log_2N
M=log2N级运算，每一级的标识为
m
=
0
,
1
,
2
,
…
,
M
−
1
m=0,1,2,…,M-1
m=0,1,2,…,M−1
（2） 蝶形单元
FFT计算结构由若干个蝶形运算单元组成，每个运算单元示意图如下：
蝶形单元的输入输入满足：
x
m
+
1
(
p
)
=
x
m
(
p
)
+
x
m
(
q
)
⋅
W
N
r
x_{m+1} (p)=x_m (p)+x_m(q) \cdot W_N^r
xm+1(p)=xm(p)+xm(q)⋅WNr
x
m
+
1
(
q
)
=
x
m
(
p
)
−
x
m
(
q
)
⋅
W
N
r
x_{m+1} (q)=x_m (p)-x_m (q) \cdot W_N^r
xm+1(q)=xm(p)−xm(q)⋅WNr
其中 q = p + 2 m q=p+2^m q=p+2m,每个蝶形运算单元运算时，进行一次乘法和两次加法，每一级中，均有N/2个蝶形单元。故完成一次FFT所需要的乘法次数为 N / 2 ⋅ l o g 2 N N/2 \cdot log_2N N/2⋅log2N，所需要的加法次数为： N ⋅ l o g 2 N N \cdot log_2N N⋅log2N
（3） 组的概念
每一级 N/2 个蝶形单元可分为若干组，每一组有着相同的结构与
W
N
r
W_N^r
WNr因子分布。例如
- m=0 时，可以分为 N/2=4 组。
- m=1 时，可以分为 N/4=2 组。
- m=M-1 时，此时只能分为 1 组。
（4）
W
N
r
W_N^r
WNr因子分布
W
2
m
+
1
r
W_{2^{m+1}}^r
W2m+1r因子存在于m级，其中
r
=
0
,
1
,
2
,
…
,
2
m
−
1
r=0,1,2,…,2^m-1
r=0,1,2,…,2m−1。举例说明：
在8点FFT第二级运算中，即m=1，蝶形因子运算为：
W
4
0
W_4^0
W40和
W
4
1
W_4^1
W41，对应的8等分的旋转因子即为
W
8
0
W_8^0
W80和
W
8
2
W_8^2
W82，对应图2-1中的W(0)和W(2)；
在硬件实现时，旋转因子事先会计算好，所以会统一成最大等分数的旋转因子。更易理解的旋转因子分布式可以用如下关系表示：
假设N点FFT第m+1级的旋转因子为
W
N
k
t
W_N^{kt}
WNkt，则其满足如下关系：
t
=
N
/
2
m
+
1
t=N/2^{m+1}
t=N/2m+1
k
=
0
,
1
,
2
,
…
,
2
m
−
1
k=0,1,2,…,2^m-1
k=0,1,2,…,2m−1
（5） 码位倒置
对于 N=8 点的 FFT 计算，X(0) ~ X(7) 位置对应的 2 进制码为：
X
(
000
)
,
X
(
001
)
,
X
(
010
)
,
X
(
011
)
,
X
(
100
)
,
X
(
101
)
,
X
(
110
)
,
X
(
111
)
X(000), X(001), X(010), X(011), X(100), X(101), X(110), X(111)
X(000),X(001),X(010),X(011),X(100),X(101),X(110),X(111)
将其位置的 2 进制码进行翻转：
X
(
000
)
,
X
(
100
)
,
X
(
010
)
,
X
(
110
)
,
X
(
001
)
,
X
(
101
)
,
X
(
011
)
,
X
(
111
)
X(000), X(100), X(010), X(110), X(001), X(101), X(011), X(111)
X(000),X(100),X(010),X(110),X(001),X(101),X(011),X(111)
此时位置对应的 10 进制为：
X
(
0
)
,
X
(
4
)
,
X
(
2
)
,
X
(
6
)
,
X
(
1
)
,
X
(
5
)
,
X
(
3
)
,
X
(
7
)
X(0), X(4), X(2), X(6), X(1), X(5), X(3), X(7)
X(0),X(4),X(2),X(6),X(1),X(5),X(3),X(7)
恰好对应FFT第一级输入数据的顺序。
2.2 FFT设计
（1）蝶形单元设计
端口定义如下：
module butterfly2(
input clk ,//系统时钟
input rstn ,//系统异步复位，低电平有效
input en ,//使能信号，表示输入数据有效
input signed [15:0] xp_real ,//Xm(p)
input signed [15:0] xp_imag ,
input signed [15:0] xq_real ,//Xm(q)
input signed [15:0] xq_imag ,
input signed [15:0] factor_real,//扩大8192 倍（左移 13 位）后的旋转因子
input signed [15:0] factor_imag,
output valid ,//输入数据有效执行信号
output signed [15:0] yp_real ,//Xm+1(p)
output signed [15:0] yp_imag ,
output signed [15:0] yq_real ,//Xm+1(q)
output signed [15:0] yq_imag
);
注意：蝶形单元为定点运算，输入factor_real和factor_imag为扩大 8192 倍（左移 13 位）后的旋转因子，且为了防止蝶形运算中的乘法和加法导致位宽逐级增大，每一级运算完成后，会对输出数据进行固定位宽的截位，也可去掉旋转因子倍数增大而带来的影响。
内部实现逻辑介绍:
蝶形单元内部使用流水线结构，将计算划分成了三级流水，每一级功能如下：
- Xm(q) 与 factor进行复数乘法的乘法运算，Xm§ 左移13位（保证扩大相同的倍数）
- Xm(q) 与 factor进行复数乘法的加法运算 a，Xm§进行寄存
- 进行蝶形运算的加法运算，并将得到的结果指定位宽截位，得到 X m + 1 ( p ) X_{m+1}(p) Xm+1(p)和 X m + 1 ( q ) X_{m+1}(q) Xm+1(q)
（2）输入缓存的设计
因为设计内容中要求输入为串行输入，所以需要对输入数据缓存。这里使用一个计数器来生成输入要寄存的地址信息，当valid_in有效时更新为counter_next的值，而counter_next的取值逻辑如下：
wire [5:0] counter_next = {6{sop_in }} & 6'b0
| {6{~sop_in}} & (counter + 1);
注意，在写入操作时，为了保证数据与地址的一致，写输入缓冲的地址使用counter计数器的次态值，即counter_next的值。其输入地址代码如下：
// store input
reg signed [15:0] mem_xm_real[63:0];
reg signed [15:0] mem_xm_imag[63:0];
reg inv_reg ;
always @(posedge clk or negedge rst_n) begin
if (~rst_n) begin
end
else if (valid_in) begin
mem_xm_real[counter_next] <= x_re;
mem_xm_imag[counter_next] <= {16{~inv}}&x_im | {16{inv}}&(-x_im); // inv==0 -> fft inv==1 -> ifft
inv_reg <= inv;
end
end
（3） 例化蝶形单元，实现fft/ifft
依据2.1 FFT算法原理可知，64点fft的实现，需要
l
o
g
2
64
=
6
log_264=6
log264=6级运算（对应for(m=0; m<=5; m=m+1)），且每一级需要64/2=32个蝶形单元（对应for (k=0; k<=31; k=k+1)）。所以例化的主体代码如下：
//butter instantiaiton
genvar m, k;
generate
//6 stage
for(m=0; m<=5; m=m+1) begin: stage
for (k=0; k<=31; k=k+1) begin: unit
butterfly2 u_butter(
.clk (clk ) ,
.rstn (rst_n ) ,
这里做了省略
);
end
end
endgenerate
另外fft算法中，要求输入数据索引二进制码位倒置，所以对xm_real[0]、xm_imag[0]初始化时，做如下处理：
genvar czw;
generate
for (czw = 0; czw<=63; czw=czw+1) begin
assign xm_real[0][czw] = mem_xm_real[index[czw]]; //输入地址二进制码翻转
assign xm_imag[0][czw] = mem_xm_imag[index[czw]];
end
endgenerate
其中index是一个一维向量，存放了0到63的二进制码位倒置值，这些值是事先用matlab计算出来的：
N=64;
index = bin2dec(fliplr(dec2bin(0:N-1)));
再者，第一级蝶形单元的计算要在输入缓冲满64个时才能开始，所以其使能的生成逻辑为counter==63，且需要第一辑32个蝶形单元同时有效，所以其代码如下：
wire en_connect [223:0] ;
wire en = counter==63;
assign {
en_connect[ 0],en_connect[ 1]
这里做了省略
en_connect[30],en_connect[31]
} = {
en ,en ,
这里做了省略
en ,en
};
（4）输出缓存的设计
因为设计内容中要求输出为串行输出，所以需要对输入数据缓存。这里同样使用一个计数器来生成输出要寄存的地址信息，每个时钟将counter_o更新为counter_next_o的值，而counter_next_o的取值逻辑如下：
wire [5:0] counter_next_o = {6{sop_out_next }} & 6'b0
| {6{~sop_out_next}} & (counter_o + 1);
注意，在写入操作时，为了保证数据与地址的一致，输出缓冲送到输出端口的地址使用counter_o计数器的次态值，即counter_next_o的值。其输出地址代码如下：
always @(posedge clk or negedge rst_n) begin
if (~rst_n) begin
counter_o <= 0;
sop_out <= 0;
valid_out <= 0;
y_re <= 0;
y_im <= 0;
end
else begin
counter_o <= counter_next_o;
y_re <= xm_real[6][counter_next_o];
y_im <= {16{~inv_r[17]}}&xm_imag[6][counter_next_o] | {16{inv_r[17]}}&(-xm_imag[6][counter_next_o]); // inv==0 -> fft inv==1 -> ifft
这里做了省略
end
end
说明
参考链接：https://www.runoob.com/w3cnote/verilog-fft.html
标签：蝶形,counter,15,运算,FFT,源码,64 来源： https://blog.csdn.net/qq_42816434/article/details/123027013
本站声明： 1. iCode9 技术分享网（下文简称本站）提供的所有内容，仅供技术学习、探讨和分享； 2. 关于本站的所有留言、评论、转载及引用，纯属内容发起人的个人观点，与本站观点和立场无关； 3. 关于本站的所有言论和文字，纯属内容发起人的个人观点，与本站观点和立场无关； 4. 本站文章均是网友提供，不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属；如您发现该文章侵犯了您的权益，可联系我们第一时间进行删除； 5. 本站为非盈利性的个人网站，所有内容不会用来进行牟利，也不会利用任何形式的广告来间接获益，纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。