免費(fèi)AI編程助手
Verilog FFT設(shè)計(jì)
FFT（Fast Fourier Transform），快速傅立葉變換，是一種 DFT（離散傅里葉變換）的高效算法。在以時(shí)頻變換分析為基礎(chǔ)的數(shù)字處理方法中，有著不可替代的作用。
FFT 原理
公式推導(dǎo)
DFT 的運(yùn)算公式為：
其中
將離散傅里葉變換公式拆分成奇偶項(xiàng)，則前 N/2 個(gè)點(diǎn)可以表示為：
同理，后 N/2 個(gè)點(diǎn)可以表示為：
由此可知，后 N/2 個(gè)點(diǎn)的值完全可以通過(guò)計(jì)算前 N/2 個(gè)點(diǎn)時(shí)的中間過(guò)程值確定。對(duì) A[k] 與 B[k] 繼續(xù)進(jìn)行奇偶分解，直至變成 2 點(diǎn)的 DFT，這樣就可以避免很多的重復(fù)計(jì)算，實(shí)現(xiàn)了快速離散傅里葉變換（FFT）的過(guò)程。
算法結(jié)構(gòu)
8 點(diǎn) FFT 計(jì)算的結(jié)構(gòu)示意圖如下。
由圖可知，只需要簡(jiǎn)單的計(jì)算幾次乘法和加法，便可完成離散傅里葉變換過(guò)程，而不是對(duì)每個(gè)數(shù)據(jù)進(jìn)行繁瑣的相乘和累加。
重要特性
- 級(jí)的概念
- 蝶形單元
- 組的概念
-
因子分布
因子存在于 m 級(jí)，其中
。 - 碼位倒置
每分割一次，稱為一級(jí)運(yùn)算。
設(shè) FFT 運(yùn)算點(diǎn)數(shù)為 N，共有 M 級(jí)運(yùn)算，則它們滿足：
每一級(jí)運(yùn)算的標(biāo)識(shí)為 m = 0, 1, 2, ..., M-1。
為了便于分割計(jì)算，F(xiàn)FT 點(diǎn)數(shù) N 的取值經(jīng)常為 2 的整數(shù)次冪。
FFT 計(jì)算結(jié)構(gòu)由若干個(gè)蝶形運(yùn)算單元組成，每個(gè)運(yùn)算單元示意圖如下：
蝶形單元的輸入輸出滿足：
其中
每一個(gè)蝶形單元運(yùn)算時(shí)，進(jìn)行了一次乘法和兩次加法。
每一級(jí)中，均有 N/2 個(gè)蝶形單元。
故完成一次 FFT 所需要的乘法次數(shù)和加法次數(shù)分別為：
每一級(jí) N/2 個(gè)蝶形單元可分為若干組，每一組有著相同的結(jié)構(gòu)與
因子分布。
例如 m=0 時(shí)，可以分為 N/2=4 組。
m=1 時(shí)，可以分為 N/4=2 組。
m=M-1 時(shí)，此時(shí)只能分為 1 組。
在 8 點(diǎn) FFT 第二級(jí)運(yùn)算中，即 m=1 ，蝶形運(yùn)算因子可以化簡(jiǎn)為：
對(duì)于 N=8 點(diǎn)的 FFT 計(jì)算，X(0) ~ X(7) 位置對(duì)應(yīng)的 2 進(jìn)制碼為：
X(000), X(001), X(010), X(011), X(100), X(101), X(110), X(111)
將其位置的 2 進(jìn)制碼進(jìn)行翻轉(zhuǎn)：
X(000), X(100), X(010), X(110), X(001), X(101), X(011), X(111)
此時(shí)位置對(duì)應(yīng)的 10 進(jìn)制為：
X(0), X(4), X(2), X(6), X(1), X(5), X(3), X(7)
恰好對(duì)應(yīng) FFT 第一級(jí)輸入數(shù)據(jù)的順序。
該特性有利于 FFT 的編程實(shí)現(xiàn)。
FFT 設(shè)計(jì)
設(shè)計(jì)說(shuō)明
為了利用仿真簡(jiǎn)單的說(shuō)明 ?FFT ?的變換過(guò)程，數(shù)據(jù)點(diǎn)數(shù)取較小的值 8。
如果數(shù)據(jù)是串行輸入，需要先進(jìn)行緩存，所以設(shè)計(jì)時(shí)數(shù)據(jù)輸入方式為并行。
數(shù)據(jù)輸入分為實(shí)部和虛部共 2 部分，所以計(jì)算結(jié)果也分為實(shí)部和虛部。
設(shè)計(jì)采用流水結(jié)構(gòu)，暫不考慮資源消耗的問(wèn)題。
為了使設(shè)計(jì)結(jié)構(gòu)更加簡(jiǎn)單，這里做一步妥協(xié)，乘法計(jì)算直接使用乘號(hào)。如果 ?FFT ?設(shè)計(jì)應(yīng)用于實(shí)際，一定要將乘法結(jié)構(gòu)換成可以流水的乘法器，或使用官方提供的效率較高的乘法器 IP。
蝶形單元設(shè)計(jì)
蝶形單元為定點(diǎn)運(yùn)算，需要對(duì)旋轉(zhuǎn)因子進(jìn)行定點(diǎn)量化。
借助 matlab 將旋轉(zhuǎn)因子擴(kuò)大 8192 倍（左移 13 位），可參考附錄。
為了防止蝶形運(yùn)算中的乘法和加法導(dǎo)致位寬逐級(jí)增大，每一級(jí)運(yùn)算完成后，要對(duì)輸出數(shù)據(jù)進(jìn)行固定位寬的截位，也可去掉旋轉(zhuǎn)因子倍數(shù)增大而帶來(lái)的影響。 代碼如下：
`timescale 1ns/100ps
/**************** butter unit *************************
Xm(p) ------------------------> Xm+1(p)
- ->
- -
-
- -
- ->
Xm(q) ------------------------> Xm+1(q)
Wn -1
*//////////////////////////////////////////////////////
module butterfly
(
input clk,
input rstn,
input en,
input signed [23:0] xp_real, // Xm(p)
input signed [23:0] xp_imag,
input signed [23:0] xq_real, // Xm(q)
input signed [23:0] xq_imag,
input signed [15:0] factor_real, // Wnr
input signed [15:0] factor_imag,
output valid,
output signed [23:0] yp_real, //Xm+1(p)
output signed [23:0] yp_imag,
output signed [23:0] yq_real, //Xm+1(q)
output signed [23:0] yq_imag);
reg [4:0] en_r ;
always @(posedge clk or negedge rstn) begin
if (!rstn) begin
en_r <= 'b0 ;
end
else begin
en_r <= {en_r[3:0], en} ;
end
end
//=====================================================//
//(1.0) Xm(q) mutiply and Xm(p) delay
reg signed [39:0] xq_wnr_real0;
reg signed [39:0] xq_wnr_real1;
reg signed [39:0] xq_wnr_imag0;
reg signed [39:0] xq_wnr_imag1;
reg signed [39:0] xp_real_d;
reg signed [39:0] xp_imag_d;
always @(posedge clk or negedge rstn) begin
if (!rstn) begin
xp_real_d <= 'b0;
xp_imag_d <= 'b0;
xq_wnr_real0 <= 'b0;
xq_wnr_real1 <= 'b0;
xq_wnr_imag0 <= 'b0;
xq_wnr_imag1 <= 'b0;
end
else if (en) begin
xq_wnr_real0 <= xq_real * factor_real;
xq_wnr_real1 <= xq_imag * factor_imag;
xq_wnr_imag0 <= xq_real * factor_imag;
xq_wnr_imag1 <= xq_imag * factor_real;
//expanding 8192 times as Wnr
xp_real_d <= {{4{xp_real[23]}}, xp_real[22:0], 13'b0};
xp_imag_d <= {{4{xp_imag[23]}}, xp_imag[22:0], 13'b0};
end
end
//(1.1) get Xm(q) mutiplied-results and Xm(p) delay again
reg signed [39:0] xp_real_d1;
reg signed [39:0] xp_imag_d1;
reg signed [39:0] xq_wnr_real;
reg signed [39:0] xq_wnr_imag;
always @(posedge clk or negedge rstn) begin
if (!rstn) begin
xp_real_d1 <= 'b0;
xp_imag_d1 <= 'b0;
xq_wnr_real <= 'b0 ;
xq_wnr_imag <= 'b0 ;
end
else if (en_r[0]) begin
xp_real_d1 <= xp_real_d;
xp_imag_d1 <= xp_imag_d;
//提前設(shè)置好位寬余量，防止數(shù)據(jù)溢出
xq_wnr_real <= xq_wnr_real0 - xq_wnr_real1 ;
xq_wnr_imag <= xq_wnr_imag0 + xq_wnr_imag1 ;
end
end
//======================================================//
//(2.0) butter results
reg signed [39:0] yp_real_r;
reg signed [39:0] yp_imag_r;
reg signed [39:0] yq_real_r;
reg signed [39:0] yq_imag_r;
always @(posedge clk or negedge rstn) begin
if (!rstn) begin
yp_real_r <= 'b0;
yp_imag_r <= 'b0;
yq_real_r <= 'b0;
yq_imag_r <= 'b0;
end
else if (en_r[1]) begin
yp_real_r <= xp_real_d1 + xq_wnr_real;
yp_imag_r <= xp_imag_d1 + xq_wnr_imag;
yq_real_r <= xp_real_d1 - xq_wnr_real;
yq_imag_r <= xp_imag_d1 - xq_wnr_imag;
end
end
//(3) discard the low 13bits because of Wnr
assign yp_real = {yp_real_r[39], yp_real_r[13+23:13]};
assign yp_imag = {yp_imag_r[39], yp_imag_r[13+23:13]};
assign yq_real = {yq_real_r[39], yq_real_r[13+23:13]};
assign yq_imag = {yq_imag_r[39], yq_imag_r[13+23:13]};
assign valid = en_r[2];
endmodule頂層例化
根據(jù) ?FFT ?算法結(jié)構(gòu)示意圖，將蝶形單元例化，完成最后的 ?FFT ?功能。
可根據(jù)每一級(jí)蝶形單元的輸入輸出對(duì)應(yīng)關(guān)系，依次手動(dòng)例化 12 次，也可利用 ?generate ?進(jìn)行例化，此時(shí)就需要非常熟悉 ?FFT ?中"組"和"級(jí)"的特點(diǎn)：
- 8 點(diǎn) FFT 設(shè)計(jì)，需要 3 級(jí)運(yùn)算，每一級(jí)有 4 個(gè)蝶形單元，每一級(jí)的組數(shù)目分別是 4、2、1。
- 每一級(jí)的組內(nèi)一個(gè)蝶形單元中兩個(gè)輸入端口的距離恒為
（m 為級(jí)標(biāo)號(hào)，對(duì)應(yīng)左移運(yùn)算 1<<m），組內(nèi)兩個(gè)蝶形單元的第一個(gè)輸入端口間的距離為 1。 - 每一級(jí)相鄰組間的第一個(gè)蝶形單元的第一個(gè)輸入端口的距離為
（對(duì)應(yīng)左移運(yùn)算 2<<m）。
例化代碼如下。
其中，矩陣信號(hào) xm_real（xm_imag）的一維、二維地址是代表級(jí)和組的標(biāo)識(shí)。
在判斷信號(hào)端口之間的連接關(guān)系時(shí)，使用了看似復(fù)雜的判斷邏輯，而且還帶有乘號(hào)，其實(shí)最終生成的電路和手動(dòng)編寫(xiě)代碼例化 12 個(gè)蝶形單元的方式是完全相同的。因?yàn)?nbsp;generate 中的變量只是輔助生成實(shí)際的電路，相關(guān)值的計(jì)算判斷都已經(jīng)在編譯時(shí)完成。這些變量更不會(huì)生成實(shí)際的電路，只是為更快速的模塊例化提供了一種方法。
timescale 1ns/100ps
module fft8 (
input clk,
input rstn,
input en,
input signed [23:0] x0_real,
input signed [23:0] x0_imag,
input signed [23:0] x1_real,
input signed [23:0] x1_imag,
input signed [23:0] x2_real,
input signed [23:0] x2_imag,
input signed [23:0] x3_real,
input signed [23:0] x3_imag,
input signed [23:0] x4_real,
input signed [23:0] x4_imag,
input signed [23:0] x5_real,
input signed [23:0] x5_imag,
input signed [23:0] x6_real,
input signed [23:0] x6_imag,
input signed [23:0] x7_real,
input signed [23:0] x7_imag,
output valid,
output signed [23:0] y0_real,
output signed [23:0] y0_imag,
output signed [23:0] y1_real,
output signed [23:0] y1_imag,
output signed [23:0] y2_real,
output signed [23:0] y2_imag,
output signed [23:0] y3_real,
output signed [23:0] y3_imag,
output signed [23:0] y4_real,
output signed [23:0] y4_imag,
output signed [23:0] y5_real,
output signed [23:0] y5_imag,
output signed [23:0] y6_real,
output signed [23:0] y6_imag,
output signed [23:0] y7_real,
output signed [23:0] y7_imag
);
//operating data
wire signed [23:0] xm_real [3:0] [7:0];
wire signed [23:0] xm_imag [3:0] [7:0];
wire en_connect [15:0] ;
assign en_connect[0] = en;
assign en_connect[1] = en;
assign en_connect[2] = en;
assign en_connect[3] = en;
//factor, multiplied by 0x2000
wire signed [15:0] factor_real [3:0] ;
wire signed [15:0] factor_imag [3:0];
assign factor_real[0] = 16'h2000; //1
assign factor_imag[0] = 16'h0000; //0
assign factor_real[1] = 16'h16a0; //sqrt(2)/2
assign factor_imag[1] = 16'he95f; //-sqrt(2)/2
assign factor_real[2] = 16'h0000; //0
assign factor_imag[2] = 16'he000; //-1
assign factor_real[3] = 16'he95f; //-sqrt(2)/2
assign factor_imag[3] = 16'he95f; //-sqrt(2)/2
//輸入初始化，和碼位有關(guān)倒置
assign xm_real[0][0] = x0_real;
assign xm_real[0][1] = x4_real;
assign xm_real[0][2] = x2_real;
assign xm_real[0][3] = x6_real;
assign xm_real[0][4] = x1_real;
assign xm_real[0][5] = x5_real;
assign xm_real[0][6] = x3_real;
assign xm_real[0][7] = x7_real;
assign xm_imag[0][0] = x0_imag;
assign xm_imag[0][1] = x4_imag;
assign xm_imag[0][2] = x2_imag;
assign xm_imag[0][3] = x6_imag;
assign xm_imag[0][4] = x1_imag;
assign xm_imag[0][5] = x5_imag;
assign xm_imag[0][6] = x3_imag;
assign xm_imag[0][7] = x7_imag;
//butter instantiaiton
//integer index[11:0] ;
genvar m, k;
generate
//3 stage
for(m=0; m<=2; m=m+1) begin: stage
for (k=0; k<=3; k=k+1) begin: unit
butterfly u_butter(
.clk (clk ) ,
.rstn (rstn ) ,
.en (en_connect[m*4 + k] ) ,
//是否再組內(nèi)？組編號(hào)+組內(nèi)編號(hào)：下組編號(hào)+新組內(nèi)編號(hào)
.xp_real (xm_real[ m ] [k[m:0] < (1<<m) ?
(k[3:m] << (m+1)) + k[m:0] :
(k[3:m] << (m+1)) + (k[m:0]-(1<<m))] ),
.xp_imag (xm_imag[ m ] [k[m:0] < (1<<m) ?
(k[3:m] << (m+1)) + k[m:0] :
(k[3:m] << (m+1)) + (k[m:0]-(1<<m))] ),
.xq_real (xm_real[ m ] [(k[m:0] < (1<<m) ?
(k[3:m] << (m+1)) + k[m:0] :
(k[3:m] << (m+1)) + (k[m:0]-(1<<m))) + (1<<m) ]), //增加蝶形單元兩個(gè)輸入端口間距離
.xq_imag (xm_imag[ m ] [(k[m:0] < (1<<m) ?
(k[3:m] << (m+1)) + k[m:0] :
(k[3:m] << (m+1)) + (k[m:0]-(1<<m))) + (1<<m) ]),
.factor_real(factor_real[k[m:0]<(1<<m)?
k[m:0] : k[m:0]-(1<<m) ]),
.factor_imag(factor_imag[k[m:0]<(1<<m)?
k[m:0] : k[m:0]-(1<<m) ]),
//output data
.valid (en_connect[ (m+1)*4 + k ] ),
.yp_real (xm_real[ m+1 ][k[m:0] < (1<<m) ?
(k[3:m] << (m+1)) + k[m:0] :
(k[3:m] << (m+1)) + (k[m:0]-(1<<m))] ),
.yp_imag (xm_imag[ m+1 ][(k[m:0]) < (1<<m) ?
(k[3:m] << (m+1)) + k[m:0] :
(k[3:m] << (m+1)) + (k[m:0]-(1<<m))] ),
.yq_real (xm_real[ m+1 ][(k[m:0] < (1<<m) ?
(k[3:m] << (m+1)) + k[m:0] :
(k[3:m] << (m+1)) + (k[m:0]-(1<<m))) + (1<<m) ]),
.yq_imag (xm_imag[ m+1 ][((k[m:0]) < (1<<m) ?
(k[3:m] << (m+1)) + k[m:0] :
(k[3:m] << (m+1)) + (k[m:0]-(1<<m))) + (1<<m) ])
);
end
end
endgenerate
assign valid = en_connect[12];
assign y0_real = xm_real[3][0] ;
assign y0_imag = xm_imag[3][0] ;
assign y1_real = xm_real[3][1] ;
assign y1_imag = xm_imag[3][1] ;
assign y2_real = xm_real[3][2] ;
assign y2_imag = xm_imag[3][2] ;
assign y3_real = xm_real[3][3] ;
assign y3_imag = xm_imag[3][3] ;
assign y4_real = xm_real[3][4] ;
assign y4_imag = xm_imag[3][4] ;
assign y5_real = xm_real[3][5] ;
assign y5_imag = xm_imag[3][5] ;
assign y6_real = xm_real[3][6] ;
assign y6_imag = xm_imag[3][6] ;
assign y7_real = xm_real[3][7] ;
assign y7_imag = xm_imag[3][7] ;
endmoduletestbench
testbench 編寫(xiě)如下，主要用于 16 路實(shí)、復(fù)數(shù)據(jù)的連續(xù)輸入。因?yàn)槊看?nbsp;FFT 只有 8 點(diǎn)數(shù)據(jù)，所以送入的數(shù)據(jù)比較隨意，并不是正弦波等規(guī)則的數(shù)據(jù)。
`timescale 1ns/100ps
module test ;
reg clk;
reg rstn;
reg en ;
reg signed [23:0] x0_real;
reg signed [23:0] x0_imag;
reg signed [23:0] x1_real;
reg signed [23:0] x1_imag;
reg signed [23:0] x2_real;
reg signed [23:0] x2_imag;
reg signed [23:0] x3_real;
reg signed [23:0] x3_imag;
reg signed [23:0] x4_real;
reg signed [23:0] x4_imag;
reg signed [23:0] x5_real;
reg signed [23:0] x5_imag;
reg signed [23:0] x6_real;
reg signed [23:0] x6_imag;
reg signed [23:0] x7_real;
reg signed [23:0] x7_imag;
wire valid;
wire signed [23:0] y0_real;
wire signed [23:0] y0_imag;
wire signed [23:0] y1_real;
wire signed [23:0] y1_imag;
wire signed [23:0] y2_real;
wire signed [23:0] y2_imag;
wire signed [23:0] y3_real;
wire signed [23:0] y3_imag;
wire signed [23:0] y4_real;
wire signed [23:0] y4_imag;
wire signed [23:0] y5_real;
wire signed [23:0] y5_imag;
wire signed [23:0] y6_real;
wire signed [23:0] y6_imag;
wire signed [23:0] y7_real;
wire signed [23:0] y7_imag;
initial begin
clk = 0; //50MHz
rstn = 0 ;
#10 rstn = 1;
forever begin
#10 clk = ~clk; //50MHz
end
end
fft8 u_fft (
.clk (clk ),
.rstn (rstn ),
.en (en ),
.x0_real (x0_real),
.x0_imag (x0_imag),
.x1_real (x1_real),
.x1_imag (x1_imag),
.x2_real (x2_real),
.x2_imag (x2_imag),
.x3_real (x3_real),
.x3_imag (x3_imag),
.x4_real (x4_real),
.x4_imag (x4_imag),
.x5_real (x5_real),
.x5_imag (x5_imag),
.x6_real (x6_real),
.x6_imag (x6_imag),
.x7_real (x7_real),
.x7_imag (x7_imag),
.valid (valid),
.y0_real (y0_real),
.y0_imag (y0_imag),
.y1_real (y1_real),
.y1_imag (y1_imag),
.y2_real (y2_real),
.y2_imag (y2_imag),
.y3_real (y3_real),
.y3_imag (y3_imag),
.y4_real (y4_real),
.y4_imag (y4_imag),
.y5_real (y5_real),
.y5_imag (y5_imag),
.y6_real (y6_real),
.y6_imag (y6_imag),
.y7_real (y7_real),
.y7_imag (y7_imag));
//data input
initial begin
en = 0 ;
x0_real = 24'd10;
x1_real = 24'd20;
x2_real = 24'd30;
x3_real = 24'd40;
x4_real = 24'd10;
x5_real = 24'd20;
x6_real = 24'd30;
x7_real = 24'd40;
x0_imag = 24'd0;
x1_imag = 24'd0;
x2_imag = 24'd0;
x3_imag = 24'd0;
x4_imag = 24'd0;
x5_imag = 24'd0;
x6_imag = 24'd0;
x7_imag = 24'd0;
@(negedge clk) ;
en = 1 ;
forever begin
@(negedge clk) ;
x0_real = (x0_real > 22'h3F_ffff) ? 'b0 : x0_real + 1 ;
x1_real = (x1_real > 22'h3F_ffff) ? 'b0 : x1_real + 1 ;
x2_real = (x2_real > 22'h3F_ffff) ? 'b0 : x2_real + 31 ;
x3_real = (x3_real > 22'h3F_ffff) ? 'b0 : x3_real + 1 ;
x4_real = (x4_real > 22'h3F_ffff) ? 'b0 : x4_real + 23 ;
x5_real = (x5_real > 22'h3F_ffff) ? 'b0 : x5_real + 1 ;
x6_real = (x6_real > 22'h3F_ffff) ? 'b0 : x6_real + 6 ;
x7_real = (x7_real > 22'h3F_ffff) ? 'b0 : x7_real + 1 ;
x0_imag = (x0_imag > 22'h3F_ffff) ? 'b0 : x0_imag + 2 ;
x1_imag = (x1_imag > 22'h3F_ffff) ? 'b0 : x1_imag + 5 ;
x2_imag = (x2_imag > 22'h3F_ffff) ? 'b0 : x2_imag + 3 ;
x3_imag = (x3_imag > 22'h3F_ffff) ? 'b0 : x3_imag + 6 ;
x4_imag = (x4_imag > 22'h3F_ffff) ? 'b0 : x4_imag + 4 ;
x5_imag = (x5_imag > 22'h3F_ffff) ? 'b0 : x5_imag + 8 ;
x6_imag = (x6_imag > 22'h3F_ffff) ? 'b0 : x6_imag + 11 ;
x7_imag = (x7_imag > 22'h3F_ffff) ? 'b0 : x7_imag + 7 ;
end
end
//finish simulation
initial #1000 $finish ;
endmodule仿真結(jié)果
大致可以看出，F(xiàn)FT 結(jié)果可以流水輸出。
用 matlab 自帶的 FFT 函數(shù)對(duì)相同數(shù)據(jù)進(jìn)行運(yùn)算，前 2 組數(shù)據(jù) FFT 結(jié)果如下。
可以看出，第一次輸入的數(shù)據(jù)信號(hào)只有實(shí)部有效時(shí)，F(xiàn)FT 結(jié)果是完全一樣的。
但是第二次輸入的數(shù)據(jù)復(fù)部也有信號(hào)，此時(shí)兩者之間的結(jié)果開(kāi)始有誤差，有時(shí)誤差還很大。
用 matlab 對(duì) Verilog 實(shí)現(xiàn)的 FFT 過(guò)程進(jìn)行模擬，發(fā)現(xiàn)此過(guò)程的 FFT 結(jié)果和 Verilog 實(shí)現(xiàn)的 FFT 結(jié)果基本一致。
將有誤差的兩種 FFT 結(jié)果取絕對(duì)值進(jìn)行比較，圖示如下。
可以看出，F(xiàn)FT 結(jié)果的趨勢(shì)大致相同，但在個(gè)別點(diǎn)有肉眼可見(jiàn)的誤差。
設(shè)計(jì)總結(jié)：
就如設(shè)計(jì)蝶形單元時(shí)所說(shuō)，旋轉(zhuǎn)因子量化時(shí)，位寬的選擇就會(huì)引入誤差。
而且每個(gè)蝶形單元的運(yùn)算結(jié)果都會(huì)進(jìn)行截取，也會(huì)引入誤差。
matlab 計(jì)算 FFT 時(shí)不用考慮精度問(wèn)題，以其最高精度對(duì)數(shù)據(jù)進(jìn)行 FFT 計(jì)算。
以上所述，都會(huì)導(dǎo)致最后兩種 FFT 計(jì)算方式結(jié)果的差異。
感興趣的學(xué)者，可以將旋轉(zhuǎn)因子和輸入數(shù)據(jù)位寬再進(jìn)行一定的增加，F(xiàn)FT 點(diǎn)數(shù)也可以增加，然后再進(jìn)行仿真對(duì)比，相對(duì)誤差應(yīng)該會(huì)減小。
附錄：matlab 使用
生成旋轉(zhuǎn)因子
8 點(diǎn) FFT 只需要用到 4 個(gè)旋轉(zhuǎn)因子。旋轉(zhuǎn)因子擴(kuò)大倍數(shù)為 8192。
clear all;close all;clc;
%=======================================================
% Wnr calcuting
%=======================================================
for r = 0:3
Wnr_factor = cos(pi/4*r) - j*sin(pi/4*r) ;
Wnr_integer = floor(Wnr_factor * 2^13) ;
if (real(Wnr_integer)<0)
Wnr_real = real(Wnr_integer) + 2^16 ; %負(fù)數(shù)的補(bǔ)碼
else
Wnr_real = real(Wnr_integer) ;
end
if (imag(Wnr_integer)<0)
Wnr_imag = imag(Wnr_integer) + 2^16 ;
else
Wnr_imag = imag(Wnr_integer);
end
Wnr(2*r+1,:) = dec2hex(Wnr_real) %實(shí)部
Wnr(2*r+2,:) = dec2hex(Wnr_imag) %虛部
endFFT 結(jié)果對(duì)比
matlab 模擬 Verilog 實(shí)現(xiàn) FFT 的過(guò)程如下，也包括 2 種 FFT 結(jié)果的對(duì)比。
clear all;close all;clc;
%=======================================================
% 8dots fft
%=======================================================
for r=0:3
Wnr(r+1) = cos(pi/4*r) - j*sin(pi/4*r) ;
end
x = [10, 20, 30, 40, 10, 20 ,30 ,40];
step = [1+2j, 1+5j, 31+3j, 1+6j, 23+4j, 1+8j, 6+11j, 1+7j];
x2 = x + step;
xm0 = [x2(0+1), x2(4+1), x2(2+1), x2(6+1), x2(1+1), x2(5+1), x2(3+1), x2(7+1)] ;
%% stage1
xm1(1) = xm0(1) + xm0(2)*Wnr(1) ;
xm1(2) = xm0(1) - xm0(2)*Wnr(1) ;
xm1(3) = xm0(3) + xm0(4)*Wnr(1) ;
xm1(4) = xm0(3) - xm0(4)*Wnr(1) ;
xm1(5) = xm0(5) + xm0(6)*Wnr(1) ;
xm1(6) = xm0(5) - xm0(6)*Wnr(1) ;
xm1(7) = xm0(7) + xm0(8)*Wnr(1) ;
xm1(8) = xm0(7) - xm0(8)*Wnr(1) ;
floor(xm1(:))
%% stage2
xm2(1) = xm1(1) + xm1(3)*Wnr(1) ;
xm2(3) = xm1(1) - xm1(3)*Wnr(1) ;
xm2(2) = xm1(2) + xm1(4)*Wnr(2) ;
xm2(4) = xm1(2) - xm1(4)*Wnr(2) ;
xm2(5) = xm1(5) + xm1(7)*Wnr(1) ;
xm2(7) = xm1(5) - xm1(7)*Wnr(1) ;
xm2(6) = xm1(6) + xm1(8)*Wnr(2) ;
xm2(8) = xm1(6) - xm1(8)*Wnr(2) ;
floor(xm2(:))
%% stage3
xm3(1) = xm2(1) + xm2(5)*Wnr(1) ;
xm3(5) = xm2(1) - xm2(5)*Wnr(1) ;
xm3(2) = xm2(2) + xm2(6)*Wnr(2) ;
xm3(6) = xm2(2) - xm2(6)*Wnr(2) ;
xm3(3) = xm2(3) + xm2(7)*Wnr(3) ;
xm3(7) = xm2(3) - xm2(7)*Wnr(3) ;
xm3(4) = xm2(4) + xm2(8)*Wnr(4) ;
xm3(8) = xm2(4) - xm2(8)*Wnr(4) ;
floor(xm3(:))
%% fft
fft1 = fft(x)
fft2 = fft(x2)
plot(1:8, abs(fft2))
hold on
plot(1:8, abs(xm3), 'r')點(diǎn)擊這里下載源碼
更多建議：