1.FFT的必要索引变换
基2算法需要位顺序的反转位逆序,而基4算法需要首先构成一个2位的数字,再反转这些数字,称为数字逆序。
1.1 位逆序和数字逆序
2.FFT的复数乘法转实数乘法
\[X(K)=∑_{n=0}^{N-1}x(n)W_N^{kn}, k=0,1,2…..N-1 \]\[=∑_{n=0}^{N-1}(Re[x(n)]+jIm[x(n)] )∙(Re[W_N^{kn} ]+jIm[W_N^{kn}]) \]\[=∑_{n=0}^{N-1}(Re[x(n)]Re[W_N^{kn} )-Im[x(n)]Im[W_N^{kn}])+j(Re[x(n)]Im[W_N^{kn} ]+Im[x(n)]Re[W_N^{kn}]) \]由式知,计算一个X(K)需要N次复数乘法和N-1次复数加法,计算N点X(K)需要\(N^2\)次复数乘法和N(N-1)次复数加法,计算一次复数乘法需要4次实数乘法和2次实数加法,一次复数的加法需要2次实数的加法。所以每计算一个X(K)需要4N次实数乘法和2N+2(N-1)=2(2N-1)次实数加法, 计算N点X(K)需要\({4N}^2\)次实数乘法和2N(2N-1)次复实加法。
2.基2FFT算法
=∑_{r=0}{r=N/2-1}▒〖x(2r)W_(N/2)kr 〗 +W_N^k ∑_{r=0}{r=N/2-1}▒〖x(2r+1)W_(N/2)kr 〗
令X_0 (k)=∑_{r=0}{r=N/2-1}▒〖x(2r)W_(N/2)kr 〗,X_1 (k)=∑_{r=0}{r=N/2-1}▒〖x(2r+1)W_(N/2)kr 〗所以
X(K)= X_0 (k)+W_N^k X_1 (k) k=0,1,2,... N/2-1
X(K+ N/2)= X_0 (k)-W_N^k X_1 (k) k=0,1,2,... N/2-1,
此时一个N点的DFT可转换为2个N/2DFT,再经过第二次分解,又将N/2点DFT分解成2个N/4DFT,以此类推,经过M次分解,最后将N点FFT分解成N个1点DFT和M级蝶形运算,而1点DFT就是时域序列本身。
1.4点的基2FFT算法
X(K)=∑_{n=0}{N-1}▒〖x(n)W_N 〗=∑_{n=0}3▒〖x(n)W_4 〗
=∑_{r=0}(r=1)▒〖x(2r)W_Nk2r 〗 +∑_{r=0}(r=1)▒〖x(2r+1)W_N(k(2r+1)) 〗
=∑_{r=0}(r=1)▒〖x(2r)W_2kr 〗 +W_4^k ∑_{r=0}(r=1)▒〖x(2r+1)W_2kr 〗
X_0 (k)=∑_{r=0}(r=1)▒〖x(2r)W_2kr 〗, X_1 (k)=∑_{r=0}(r=1)▒〖x(2r+1)W_2kr 〗
X(K)= X_0 (k)+W_N^k X_1 (k) k=0,1 X(K+ 2)= X_0 (k)-W_N^k X_1 (k) k=0,1
X_0 (0)= x(0)+x(2),X_0 (1)= x(0)-x(2)
X_1 (0)= x(1)+x(3),X_1 (1)= x(1)-x(3)
X(0)= X_0 (0)+W_4^0 X_1 (0) = X_0 (0)+X_1 (0)
X(1)= X_0 (1)+W_4^1 X_1 (1)=X_0 (0)-jX_1 (1)
X(2)= X_0 (0)-W_4^0 X_1 (0) = X_0 (0)-X_1 (0)
X(3)= X_0 (1)-W_4^1 X_1 (1)=X_0 (0)+jX_1 (1)
2.8点的基2FFT算法
X(K)=∑_{n=0}{N-1}▒〖x(n)W_N 〗=∑_{n=0}7▒〖x(n)W_4 〗=∑_{r=0}(r=3)▒〖x(2r)W_Nk2r 〗+∑_{r=0}(r=3)▒〖x(2r+1)W_N(k(2r+1)) 〗
=∑_{r=0}(r=1)▒〖x(2r)W_2kr 〗+W_8^k ∑_{r=0}(r=3)▒〖x(2r+1)W_2kr 〗=X_0 (k)+W_8^k X_1 (k) k=0,1,2,3,4
X(K+4)= X_0 (k)-W_8^k X_1 (k)
X_0 (k)=∑_{r=0}^(r=3)▒〖x(2r) W_4kr=∑_(l=0)1▒〖x(4l) W_2kl+W_4k ∑_(l=0)^1▒〖x(4l+2) W_2^kl 〗〗〗.
=X_2 (k)+W_4^k X_3 (k) k=0,1
X_0 (k+2)=X_2 (k)-W_4^k X_3 (k) k=0,1
X_2 (k)=∑_(l=0)^1▒〖x(4l) W_2^kl 〗 X_2 (0)=x(0)+x(4),同理计算X_3 (0)= x(2)+x(6)
X_1 (k)=∑_{r=0}(r=3)▒〖x(2r+1)W_4kr 〗=∑_(l=0)^1▒〖x(4l+1) W_2kl+W_4k ∑_(l=0)^1▒〖x(4l+3) W_2^kl 〗〗
= X_4 (k)+W_4^k X_5 (k) k=0,1
X_1 (k+2)=X_4 (k)-W_4^k X_5 (k) k=0,1
X_4 (0)= x(1)+x(5),X_5 (0)= x(3)+x(7)
3.DIT-FFT的运算规律及编程思想
为了最终写出DIT-FFT运算程序或设计出硬件实现电路,下面介绍DIT-FFT的运算规律和编程思想。
1.原位计算
由图4.2.4可以看出,DIT-FFT的运算过程很有规律。N=2^M点的FFT共进行M级运算,每级由N/2个叠形运算组成。同一级中,每个螺形的两个输入数据只对计算本蝶形有用,而且每个蝶形的输人、输出数据结点又同在一条水平线上,这就意味着计算完1个蝶形后,所得输出数据可立即存人原输入数据所占用的存储单元(数组元素)这样,经过M级运算后,原来存放输入序列数据的N个存储单元(数组A)中便依次存放X(k)的N个值。这种利用同一存储单元存储蝶形计算输入、输出数据的方法称为原位(址)计算,原位计算可节省大量内存,从而使设备成本降低。
2.旋转因子的变化规律
如上所述,N点 DIT-FFT运算流图中,每级都有N/2个蝶形。每个蝶形都要乘以因子W_NP,称其为旋转因子,p为旋转因子的指数。但各级的旋转因子和循环方式都有所不同,为了编写计算程序,应先找出旋转因子W_NP与运算级数的关系,用L表示从左到右的运算级数(L=1,2・…,M)。观察基2的8FFT运算图不难发现,第L级共有2(L-1)个不同的旋转因子,N=23=8时的各级旋转因子表示如下,
L=1时,W_NP=W_(N/4)J=W_(2L)J,J=0
L=2时,W_NP=W_(N/2)J=W_(2L)J,J=0,1
L=3时,W_NP=W_NJ=W_(2L)J,J=0,1,2,3
对于N=2M的一般情况,第L级的旋转因子为W_NP=W_(2L)J,J=0,1,2,.... 2(L-1)-1.因为2L=2M×2(L-M)=N∙2^(L-M),所以:
W_NP=W_(N∙2(L-M))J=W_N(J∙2^(M-L) ), J=0,1,2,.... 2^(L-1)-1, (1)
P= J∙2^(M-L), (2)
这样就可以按照(1)和(2)式确定第L级运算的的旋转因子(实际编程时L为最外层循环变量)。
3.蝶形运算规律
设序列x(n)经时域抽选(倒序)后,按图4.2.4所示的次序(倒序)存人数组A中,如果蝶形运算的两个输入数据相距B个点,应用原位计算,则蝶形运算可表示成如下形式:
〖A 〗L(J)= 〖A 〗(L-1)(J)+ 〖A 〗(L-1)(J+B) W_N^P
〖A 〗L(J+B)= 〖A 〗(L-1)(J)- 〖A 〗(L-1)(J+B) W_N^P,P= J∙2^(M-L), J=0,1,2,.... 2^(L-1)-1;L=1,2,3,…M.
下标L表示第L级运算、〖A 〗L(J)则表示第L级运算后的数组元素A(J)的值(即第L级蝶形的输出数据)。而〖A 〗(L-1)(J)表示第L级运算前A(J)的值(即第L级蝶形的输入数据).
4.编程思想及程序框图
仔细观察基2的8点DITFFT运算图,还可以归纳出一些对编程有用的运算规律:第L级中,每个蝶形的两个输入数据相距B=2(L-1)的力个点,每级有B个不同的旋转因子;同一旋转因子对应着间隔为2L点的2^(M-L)个蝶形。
总结上述运算规律,便可采用下述运算方法:先从输入端(第1级)开始,逐级进行,共进行M级运算。在进行第L级运算时,依次求出B个不同的旋转因子,每求出一个旋转因子,计算完它对应的所有2^(M-L)个蝶形,这样,我们可用三重循环程序实现 DIT FFT运算,程序框图如上图所示。
另列FFT算法运算流图的输出X(K)为自然顺序但为了适应原位计算、其输入序列不是按x(n)的自然顺序排列,这种经过M次偶奇抽选后的排序称为序列x(n)的倒序(倒位),因此,在运算M级蝶形之前先对序列x(n)进行倒序,下面介绍倒序算法。
显而易见,只要将顺序数(n2n1n)的二进制位倒置,则得对应的二进制倒序值。按这一规律,用硬件电路和汇编语言程序产生倒序数很容易。但用有些高级语言程序实现时,直接倒置二进制数位是不行的,因此必须找出产生倒序数的十进制运算规律。自然顺序数增加1,是在顺序数的二进制数最低位加1,逢2向高位进位,而倒序数则是在M位二进数最高位加1,逢2向低位进位。
3.基4FFT算法
X(K)=∑_{n=0}{N-1}▒〖x(n)W_N 〗
=∑_{r=0}(r=N/4-1)▒〖x(4r)W_Nk4r 〗 +∑_{r=0}(r=N/4-1)▒〖x(4r+1)W_N(k(4r+1)) 〗
+∑_{r=0}(r=N/4-1)▒〖x(4r+2)W_N(k(4r+2)) 〗 +∑_{r=0}(r=N/4-1)▒〖x(4r+3)W_N(k(4r+3)) 〗
=∑_{r=0}(r=N/4-1)▒〖x(4r)W_Nk4r 〗 +W_N^k ∑_{r=0}(r=N/4-1)▒〖x(4r+1)W_Nk4r 〗
+W_N^2k ∑_{r=0}(r=N/4-1)▒〖x(4r+2)W_Nk4r 〗 +W_N^3k ∑_{r=0}(r=N/4-1)▒〖x(4r+3)W_Nk4r 〗
=∑_{r=0}(r=N/4-1)▒〖x(4r)W_(N/4)kr 〗 +W_N^k ∑_{r=0}(r=N/4-1)▒〖x(4r+1)W_(N/4)kr 〗
+W_N^2k ∑_{r=0}(r=N/4-1)▒〖x(4r+2)W_(N/4)kr 〗 +W_N^3k ∑_{r=0}(r=N/4-1)▒〖x(4r+3)W_(N/4)kr 〗
令X_0 (k)=∑_{r=0}(r=N/4-1)▒〖x(4r)W_(N/4)kr 〗, X_1 (k)=∑_{r=0}(r=N/4-1)▒〖x(4r+1)W_(N/4)kr 〗
X_2 (k)=∑_{r=0}(r=N/4-1)▒〖x(4r+2)W_(N/4)kr 〗, X_3 (k)=∑_{r=0}(r=N/4-1)▒〖x(4r+3)W_(N/4)kr 〗
X(k)= X_0 (k)+W_N^k X_1 (k)+W_N^2k X_2 (k)+W_N^3k X_3 (k) k=0,1,2,... N/4-1
X(k+N/4)= X_0 (k)+W_N^(k+N/4) X_1 (k) +W_N^(2(k+N/4)) X_2 (k)+W_N^(3(k+N/4)) X_3 (k)
=X_0 (k)-jW_N^k X_1 (k)-W_N^2k X_2 (k)+jW_N^3k X_3 (k)
X(k+N/2)= X_0 (k)+W_N^(k+N/2) X_1 (k) +W_N^(2(k+N/2)) X_2 (k)+W_N^(3(k+N/2)) X_3 (k)
=X_0 (k)-W_N^k X_1 (k) +W_N^2k X_2 (k)-W_N^3k X_3 (k)
X(k+3N/4)= X_0 (k)+W_N^(k+3N/4) X_1 (k) +W_N^(2(k+3N/4)) X_2 (k)+W_N^(3(k+3N/4)) X_3 (k)
=X_0 (k)+〖jW〗_N^k X_1 (k)-W_N^2k X_2 (k)-jW_N^3k X_3 (k)
所以总结如下:
X(k)= X_0 (k)+W_N^k X_1 (k)+W_N^2k X_2 (k)+W_N^3k X_3 (k) k=0,1,2,... N/4-1
X(k+N/4) =X_0 (k)-jW_N^k X_1 (k)-W_N^2k X_2 (k)+jW_N^3k X_3 (k)
X(k+N/2)= X_0 (k)-W_N^k X_1 (k) +W_N^2k X_2 (k)-W_N^3k X_3 (k)
X(k+3N/4) =X_0 (k)+〖jW〗_N^k X_1 (k)-W_N^2k X_2 (k)-jW_N^3k X_3 (k)
U_0(K)= X_0 (k)+W_N^2k X_2 (k), U_1(K)= W_N^k X_1 (k)+W_N^3k X_3 (k)
U_2(K)= X_0 (k)-W_N^2k X_2 (k), U_3(K)= W_N^k X_1 (k)-W_N^3k X_3 (k)
X(k)= U_0(K)+ U_1(K), X(k+N/4)= U_2(K)-jU_3(K),
X(k+N/2)= U_0(K)- U_1(K),X(k+3N/4) =U_2(K)+jU_3(K)
此时一个N点的DFT可转换为N/4DFT。对于一个4点FFT的基4蝶形运算单元为:
对于一个16点FFT的基4运算:
X_0 (k)=∑_{r=0}(r=3)▒〖x(4r)W_(N/4)kr 〗, X_1 (k)=∑_{r=0}(r=3)▒〖x(4r+1)W_(N/4)kr 〗
X_2 (k)=∑_{r=0}(r=3)▒〖x(4r+2)W_(N/4)kr 〗, X_3 (k)=∑_{r=0}(r=3)▒〖x(4r+3)W_(N/4)kr 〗
X_0 (0)=x(0)+x(4)+x(8)+x(12) X_0 (1)=x(0)-jx(4)-x(8)+jx(12)
X_0 (2)=x(0)-x(4)+x(8)-x(12) 〖 X〗_0 (3)=x(0)+jx(4)-x(8)-jx(12)
X(k)= X_0 (k)+W_16^k X_1 (k)+W_16^2k X_2 (k)+W_16^3k X_3 (k) k=0,1,2,3
X(k+4)= X_0 (k)-jW_16^k X_1 (k)-W_16^2k X_2 (k)+jW_16^3k X_3 (k)
X(k+8)=X_0 (k)-W_16^k X_1 (k) +W_16^2k X_2 (k)-W_16^3k X_3 (k)
X(k+12)=X_0 (k)+〖jW〗_16^k X_1 (k)-W_16^2k X_2 (k)-jW_16^3k X_3 (k)
1.N为16的基4FFT分解过程
一个4点序列的DFT运算流图:
对于基4运算,每个蝶形运算有4点输入和4点输出,因此总共有N/4个蝶形运算单元。由于每一个蝶形运算仅需3次复数乘法和8次复数加法,因此,对每一次分解,蝶式运算所需的复数乘法次数和复数加法次数分别为3N/4次和2N次。
在N=4^M(M为正整数)的一般情形下,N点序列分解成一点序列需要作M=log_4N次分解,因此,基4时分FFT算法所需要的复数乘法次数为3N/4 log_4N次,复数加法次数为2Nlog_4N.
2.DIT-FFT的运算规律及编程思想
为了最终写出DIT-FFT运算程序或设计出硬件实现电路,下面介绍DIT-FFT的运算规律和编程思想。
1.原位计算
由基4的16点FFT蝶形图可知,DIT-FFT的运算过程很有规律。N=4^M点的FFT共进行M级运算,每级由N/4个叠形运算组成。同一级中,每个螺形的四个输入数据只对计算本蝶形有用,而且每个蝶形的输人、输出数据结点又同在一条水平线上,这就意味着计算完1个蝶形后,所得输出数据可立即存人原输入数据所占用的存储单元(数组元素)这样,经过M级运算后,原来存放输入序列数据的N个存储单元(数组A)中便依次存放X(k)的N个值。这种利用同一存储单元存储蝶形计算输入、输出数据的方法称为原位(址)计算,原位计算可节省大量内存,从而使设备成本降低。
2.旋转因子的变化规律
如上所述,N点 DIT-FFT运算流图中,每级都有N/4个蝶形。每个蝶形都要乘以因子W_NP,称其为旋转因子,p为旋转因子的指数。但各级的旋转因子和循环方式都有所不同,为了编写计算程序,应先找出旋转因子W_NP与运算级数的关系,用L表示从左到右的运算级数(L=1,2・…,M)。观察蝶形运算图不难发现,第L级共有4(L-1)个不同的旋转因子,N=43=64时的各级旋转因子表示如下,
L=1时,W_NP=W_(N/16)J=W_(4L)J, W_N2P=W_(N/16)2J=W_(4L)2J, W_N3P=W_(N/16)3J=W_(4L)3J,J=0
L=2时,W_NP=W_(N/2)J=W_(4L)J, W_N2P=W_(N/2)2J=W_(4L)2J,W_N3P=W_(N/2)3J=W_(4L)3J,J=0,1,2,3
L=3时,W_NP=W_NJ=W_(4L)J, W_N2P=W_N2J=W_(4L)2J,W_N3P=W_N3J=W_(4L)3J,J=0,1,2,...15
对于N=4M的一般情况,第L级的旋转因子为W_NP=W_(4L)J,W_N2P=W_(4L)2J,W_N3P=W_(4L)3J,J=0,1,2,.... 4(L-1)-1.因为4L=4M×4(L-M)=N∙4^(L-M),所以:
W_NP=W_(N∙4(L-M))J=W_N(J∙4^(M-L) ), J=0,1,2,.... 4^(L-1)-1, (1)
P= J∙4^(M-L), (2)
这样就可以按照(1)和(2)式确定第L级运算的的旋转因子(实际编程时L为最外层循环变量)。
3.蝶形运算规律???
设序列x(n)经时域抽选(倒序)后,按图4.2.4所示的次序(倒序)存人数组A中,如果蝶形运算的4个输入数据分别相距B个点,应用原位计算,则蝶形运算可表示成如下形式:
〖A 〗L(J)= 〖A 〗(L-1)(J)+ 〖A 〗(L-1)(J+B) W_N^P
〖A 〗L(J+B)= 〖A 〗(L-1)(J)- 〖A 〗(L-1)(J+B) W_N^P,P= J∙2^(M-L), J=0,1,2,.... 4^(L-1)-1;L=1,2,3,…M.
X(k)= X_0 (k)+W_16^k X_1 (k)+W_16^2k X_2 (k)+W_16^3k X_3 (k) k=0,1,2,3
〖A 〗_L(J)=
X(k+4)= X_0 (k)-jW_16^k X_1 (k)-W_16^2k X_2 (k)+jW_16^3k X_3 (k)
X(k+8)=X_0 (k)-W_16^k X_1 (k) +W_16^2k X_2 (k)-W_16^3k X_3 (k)
X(k+12)=X_0 (k)+〖jW〗_16^k X_1 (k)-W_16^2k X_2 (k)-jW_16^3k X_3 (k)
下标L表示第L级运算、〖A 〗L(J)则表示第L级运算后的数组元素A(J)的值(即第L级蝶形的输出数据)。而〖A 〗(L-1)(J)表示第L级运算前A(J)的值(即第L级蝶形的输入数据).
4.基4算法的模块实现
FFT模块的构成:
旋转因子乘法器;2.基4蝶形运算器;3.地址产生器;4.控制信号产生器;5中间变量存放组。
1.旋转因子乘法器的实现
FFT设计一个旋转因子存储单元可以提高FFT的计算速度,因此需要合理的设计旋转因子的存储单元,在FFT便件设计中占有重要地位。可以利用高效的算法来实现存储旋转因子。
3.地址产生单元生成 RAM读写地址, 写使能信号以及相关模块的启动、 控制信号, 是系统的控制核心;4 点蝶形运算单元的最后一级输出不是顺序的 ;旋转因子产生单元生成复乘运算中的旋转因子的角度数据;旋转因子ROM中预置了每一级运算中所需的旋转因子
1.运算模块为基4运算模块,控制模块产生所有的控制信号,存储器1和2的读写地址,写使能,运算模块的启动信号因子表的读写地址等信号。存储器1作为当前输入标志对应输入N点数据缓冲器,存储器2作为中间变量结果存储器,用于存储运算模块计算出各个通道的结果。
举例:外部输入为N点数据段流和启动信号,一方面,外部数据存入储器1中,同时通过控制模块的控制,读出存储器1中的前段N点数据和ROM表中的因子及相关控制信号送入运算核心模块进行各个Pass运算,每个Pass运算的输出都存入存储器2中,最后一个Pass的计算结果存入存储器2中,启动并在下一个启动头到来后,输出计算结果。
4.旋转因子与数据相乘实现
//通过旋转因子的地址查找表查找相位因子对应的旋转因子与输入
//数据相乘。
`timescale 1ns/1ps
//name: phase_factor
//function: multiply with input data by phase factor
module phase_factor(line1,line2,line3,line4,
in_pathway1,in_pathway2,in_pathway3,in_pathway4,
lut_phase,CLK20M);
input CLK20M;
input [5:0] lut_phase; // address for lut of phase factor,
//最多2^6=64 phase factor
input [31:0] line1,line2,line3,line4; // input data
output [31:0] in_pathway1,in_pathway2,in_pathway3,in_pathway4;
// result/output data
reg [31:0] in_pathway1,in_pathway2,in_pathway3,in_pathway4;
reg [31:0] delay_line1;
wire [31:0] dout1,dout2,dout3,dout4,
dout5,dout6,dout7,dout8,
dout9,dout10,dout11,dout12;
always @(posedge CLK20M) // delay for wait other 3 line delay multiply
delay_line1 <= line1;
//first16 real,last16image
wire [15:0] temp1_r1 = line1 [31:16];
wire [15:0] temp1_i1 = line1 [15: 0];
wire [15:0] temp1_r2 = line2 [31:16];
wire [15:0] temp1_i2 = line2 [15: 0];
wire [15:0] temp1_r3 = line3 [31:16];
wire [15:0] temp1_i3 = line3 [15: 0];
wire [15:0] temp1_r4 = line4 [31:16];
wire [15:0] temp1_i4 = line4 [15: 0];
wire [15:0] temp2_r2; //= fac2 [31:16];
wire [15:0] temp2_i2; //= fac2 [15: 0];
wire [15:0] temp2_r3; //= fac3 [31:16];
wire [15:0] temp2_i3; //= fac3 [15: 0];
wire [15:0] temp2_r4; //= fac4 [31:16];
wire [15:0] temp2_i4; //= fac4 [15: 0];
//d_rom_fft_phase2 m2(.A(lut_phase),.SPO(fac2));
//d_rom_fft_phase3 m3(.A(lut_phase),.SPO(fac3));
//d_rom_fft_phase4 m4(.A(lut_phase),.SPO(fac4));
d_rom_fft_phase2r m2r(.A(lut_phase),.SPO(temp2_r2));
d_rom_fft_phase2i m2i(.A(lut_phase),.SPO(temp2_i2));
d_rom_fft_phase3r m3r(.A(lut_phase),.SPO(temp2_r3));
d_rom_fft_phase3i m3i(.A(lut_phase),.SPO(temp2_i3));
d_rom_fft_phase4r m4r(.A(lut_phase),.SPO(temp2_r4));
d_rom_fft_phase4i m4i(.A(lut_phase),.SPO(temp2_i4));
// m16 is 16 multiply 16 ,has one delay for pipline
// for multiply (a+bj) * (c+dj) = ac-bd+(ad+bc)j
m16 u1(.clk(CLK20M),.a(temp1_r2),.b(temp2_r2),.o(dout1)); //ac
m16 u2(.clk(CLK20M),.a(temp1_i2),.b(temp2_i2),.o(dout2)); //bd
m16 u3(.clk(CLK20M),.a(temp1_r2),.b(temp2_i2),.o(dout3)); //ad
m16 u4(.clk(CLK20M),.a(temp1_i2),.b(temp2_r2),.o(dout4)); //bc
wire [15:0] dot1 =dout1[29:14]; //= dout1[30] ? dout1[30:15] : dout1[30:15] +1'b1;
wire [15:0] dot2 =dout2[29:14];//= dout2[30] ? dout2[30:15] : dout2[30:15] +1'b1;
wire [15:0] dot3 =dout3[29:14];//= dout3[30] ? dout3[30:15] : dout3[30:15] +1'b1;
wire [15:0] dot4 =dout4[29:14];//= dout4[30] ? dout4[30:15] : dout4[30:15] +1'b1;
wire [15:0] temp3_r2; //= dout1 [30:15] + (~dout2[30:15] +1'b1);
wire [15:0] temp3_i2; //= dout3 [30:15] + dout4 [30:15];
subtract16 sub1(.A(dot1),.B(dot2),.S(temp3_r2));
adder16 add1(.A(dot3),.B(dot4),.S(temp3_i2));
m16 u5(.clk(CLK20M),.a(temp1_r3),.b(temp2_r3),.o(dout5));
m16 u6(.clk(CLK20M),.a(temp1_i3),.b(temp2_i3),.o(dout6));
m16 u7(.clk(CLK20M),.a(temp1_r3),.b(temp2_i3),.o(dout7));
m16 u8(.clk(CLK20M),.a(temp1_i3),.b(temp2_r3),.o(dout8));
wire [15:0] dot5 =dout5[29:14];//= dout5[30] ? dout5[30:15] : dout5[30:15] +1'b1;
wire [15:0] dot6 =dout6[29:14];//= dout6[30] ? dout6[30:15] : dout6[30:15] +1'b1;
wire [15:0] dot7 =dout7[29:14];//= dout7[30] ? dout7[30:15] : dout7[30:15] +1'b1;
wire [15:0] dot8 =dout8[29:14];//= dout8[30] ? dout8[30:15] : dout8[30:15] +1'b1;
wire [15:0] temp3_r3; //= dout5 [30:15] + (~dout6[30:15] +1'b1);
wire [15:0] temp3_i3; //= dout7 [30:15] + dout8 [30:15];
subtract16 sub2(.A(dot5),.B(dot6),.S(temp3_r3));
adder16 add2(.A(dot7),.B(dot8),.S(temp3_i3));
m16 u9(.clk(CLK20M),.a(temp1_r4),.b(temp2_r4),.o(dout9));
m16 u10(.clk(CLK20M),.a(temp1_i4),.b(temp2_i4),.o(dout10));
m16 u11(.clk(CLK20M),.a(temp1_r4),.b(temp2_i4),.o(dout11));
m16 u12(.clk(CLK20M),.a(temp1_i4),.b(temp2_r4),.o(dout12));
wire [15:0] dot9 = dout9[29:14];//= dout9[30] ? dout9[30:15] : dout9[30:15] +1'b1;
wire [15:0] dot10 =dout10[29:14];//= dout10[30] ? dout10[30:15] : dout10[30:15] +1'b1;
wire [15:0] dot11 =dout11[29:14];//= dout11[30] ? dout11[30:15] : dout11[30:15] +1'b1;
wire [15:0] dot12 =dout12[29:14];//= dout12[30] ? dout12[30:15] : dout12[30:15] +1'b1;
wire [15:0] temp3_r4; //= dout9 [30:15] + (~dout10[30:15] +1'b1);
wire [15:0] temp3_i4; //= dout11 [30:15] + dout12 [30:15];
subtract16 sub3(.A(dot9),.B(dot10),.S(temp3_r4));
adder16 add3(.A(dot11),.B(dot12),.S(temp3_i4));
always @(posedge CLK20M) begin
in_pathway1 <= delay_line1;
in_pathway2 <= {temp3_r2,temp3_i2};
in_pathway3 <= {temp3_r3,temp3_i3};
in_pathway4 <= {temp3_r4,temp3_i4};
end
endmodule
//function finish
5.蝶形运算单元
`timescale 1ns/1ps
//name: butterfly
//function: the butterfly unit for radix 4 structure,ever step is same
//add reg_rx reg_ix for high speed pipline; 添加了rx和ix寄存器用于高速流水线实现。
module ifft_butterfly (in_pathway1,in_pathway2,in_pathway3,in_pathway4,
out_pathway1,out_pathway2,out_pathway3,out_pathway4,
CLK20M);
input CLK20M;
input [31:0] in_pathway1,in_pathway2,in_pathway3,in_pathway4;//4 pathway input butterfly
output[31:0]out_pathway1,out_pathway2,out_pathway3,out_pathway4;//4pathwayout butterfly
wire [15:0] in_real1,in_real2,in_real3,in_real4,
in_imag1,in_imag2,in_imag3,in_imag4;
reg [15:0] reg_r1,reg_r2,reg_r3,reg_r4; //中间实部寄存器
reg [15:0] reg_i1,reg_i2,reg_i3,reg_i4; //中间虚部寄存器
//高16位存储的是实部,低16位存储的是虚部
assign in_real1 = in_pathway1 [31:16];// first 16 for real and last for image
assign in_imag1 = in_pathway1 [15: 0];
assign in_real2 = in_pathway2 [31:16];
assign in_imag2 = in_pathway2 [15: 0];
assign in_real3 = in_pathway3 [31:16];
assign in_imag3 = in_pathway3 [15: 0];
assign in_real4 = in_pathway4 [31:16];
assign in_imag4 = in_pathway4 [15: 0];
// begin caculate
// first step for butterfly
wire [15:0] temp1_r1; //= in_real1 + in_real3;
wire [15:0] temp1_i1; //= in_imag1 + in_imag3;
adder16 add4(.A(in_real1),.B(in_real3),.S(temp1_r1));
adder16 add5(.A(in_imag1),.B(in_imag3),.S(temp1_i1));
wire [15:0] temp1_r2; //= in_real2 + in_real4;
wire [15:0] temp1_i2; //= in_imag2 + in_imag4;
adder16 add6(.A(in_real2),.B(in_real4),.S(temp1_r2));
adder16 add7(.A(in_imag2),.B(in_imag4),.S(temp1_i2));
wire [15:0] temp1_r3; //= in_real1 + (~in_real3 +1'b1);// 补码计算,取反加1,in_real1-in_real3
wire [15:0] temp1_i3; //= in_imag1 + (~in_imag3 +1'b1);//
subtract16 sub4(.A(in_real1),.B(in_real3),.S(temp1_r3));
subtract16 sub5(.A(in_imag1),.B(in_imag3),.S(temp1_i3));
wire [15:0] temp1_r4; //= in_real2 + (~in_real4 +1'b1);// in_real2-in_real4
wire [15:0] temp1_i4; //= in_imag2 + (~in_imag4 +1'b1);//
subtract16 sub6(.A(in_real2),.B(in_real4),.S(temp1_r4));
subtract16 sub7(.A(in_imag2),.B(in_imag4),.S(temp1_i4));
always @(posedge CLK20M) begin // reg for pipline
reg_r1 <= temp1_r1;
reg_r2 <= temp1_r2;
reg_r3 <= temp1_r3;
reg_r4 <= temp1_r4;
reg_i1 <= temp1_i1;
reg_i2 <= temp1_i2;
reg_i3 <= temp1_i3;
reg_i4 <= temp1_i4;
end
// second step for butterfly
// (a+bj)*j = -b+aj
wire [15:0] temp2_r1; //= temp1_r1 + temp1_r2; // wire1 = 1 +2;
wire [15:0] temp2_i1; //= temp1_i1 + temp1_i2;
adder16 add8(.A(reg_r1),.B(reg_r2),.S(temp2_r1));
adder16 add9(.A(reg_i1),.B(reg_i2),.S(temp2_i1));
wire [15:0] temp2_r2; //= temp1_r3 + temp1_i4; // wire2 = 3 -4*j;
wire [15:0] temp2_i2; //= temp1_i3 + (~temp1_r4 +1'b1); //
adder16 add10(.A(reg_r3),.B(reg_i4),.S(temp2_r2));
subtract16 sub8 (.A(reg_i3),.B(reg_r4),.S(temp2_i2));
wire [15:0] temp2_r3; //= temp1_r1 + (~temp1_r2 +1'b1); // wire3 = 1 -2;
wire [15:0] temp2_i3; //= temp1_i1 + (~temp1_i2 +1'b1); //
subtract16 sub9 (.A(reg_r1),.B(reg_r2),.S(temp2_r3));
subtract16 sub10(.A(reg_i1),.B(reg_i2),.S(temp2_i3));
wire [15:0] temp2_r4; //= temp1_r3 + (~temp1_i4 +1'b1); // wire4 = 3 +4*j;
wire [15:0] temp2_i4; //= temp1_i3 + temp1_r4;
subtract16 sub11(.A(reg_r3),.B(reg_i4),.S(temp2_r4));
adder16 add11(.A(reg_i3),.B(reg_r4),.S(temp2_i4));
// end caculate
assign out_pathway1 = {temp2_r1,temp2_i1};
assign out_pathway2 = {temp2_r2,temp2_i2};
assign out_pathway3 = {temp2_r3,temp2_i3};
assign out_pathway4 = {temp2_r4,temp2_i4};
// function finish
endmodule
//analysis: first 16 for real and last for image
in_real1 = in_pathway1 [31:16]; in_imag1 = in_pathway1 [15: 0];
in_real2 = in_pathway2 [31:16]; in_imag2 = in_pathway2 [15: 0];
in_real3 = in_pathway3 [31:16]; in_imag3 = in_pathway3 [15: 0];
in_real4 = in_pathway4 [31:16]; in_imag4 = in_pathway4 [15: 0];
// first step for butterfly
temp1_r1= in_real1 + in_real3; temp1_i1= in_imag1 + in_imag3;
temp1_r2= in_real2 + in_real4; temp1_i2= in_imag2 + in_imag4;
temp1_r3= in_real1 - in_real3 ; temp1_i3= in_imag1 - in_imag3;
temp1_r4= in_real2 - in_real4; temp1_i4=in_imag2 - in_imag4;
// second step for butterfly, // (a+bj)*j = -b+aj
temp2_r1= temp1_r1+temp1_r2= in_real1 + in_real3+ in_real2 + in_real4;
temp2_i1=temp1_i1 + temp1_i2= in_imag1 + in_imag3+in_imag2 + in_imag4;
temp2_r2= temp1_r3 + temp1_i4= in_real1 - in_real3+ in_imag2 - in_imag4;
temp2_i2= temp1_i3-temp1_r4= in_imag1 - in_imag3-(in_real2 - in_real4);
temp2_r3= temp1_r1 -temp1_r2= in_real1 + in_real3-( in_real2 + in_real4);
temp2_i3= temp1_i1 -temp1_i2= in_imag1 + in_imag3-( in_imag2 + in_imag4);
temp2_r4= temp1_r3-temp1_i4= in_real1 - in_real3-( in_imag2 - in_imag4);
temp2_i4= temp1_i3 + temp1_r4= in_imag1 - in_imag3+ in_real2 - in_real4;