結(jié)合頻譜移位的二維傅里葉變換FPGA實(shí)現(xiàn)

錢 輝,史 瑤,龔 敏,高 博*

(1.四川大學(xué)物理科學(xué)與技術(shù)學(xué)院微電子學(xué)系,成都 610064;2.微電子技術(shù)四川省重點(diǎn)實(shí)驗(yàn)室,成都 610064)

結(jié)合頻譜移位的二維傅里葉變換FPGA實(shí)現(xiàn)

錢 輝1,2,史 瑤1,2,龔 敏1,2,高 博1,2*

(1.四川大學(xué)物理科學(xué)與技術(shù)學(xué)院微電子學(xué)系,成都 610064;2.微電子技術(shù)四川省重點(diǎn)實(shí)驗(yàn)室,成都 610064)

在三維成像技術(shù)中,對(duì)圖像進(jìn)行傅里葉變換之后需要通過頻譜移位將零頻集中到中心位置。為了提高圖像處理的運(yùn)算速度,提出了一種結(jié)合頻譜移位的二維傅里葉變換的FPGA實(shí)現(xiàn)方法,將頻譜移位結(jié)合到二維傅里葉變換硬件系統(tǒng)中,在實(shí)現(xiàn)圖像二維傅里葉變換的同時(shí)也完成了一半的頻譜移位。采用ALTERA 系列EP4CE115F29C7芯片,針對(duì)256×256的圖像實(shí)現(xiàn)設(shè)計(jì),最高工作頻率分別達(dá)到84.2 MHz;資源消耗為3 849個(gè)LE。采用SignalTap Ⅱ Logic Analyzer工具實(shí)時(shí)驗(yàn)證了模塊的正確性。

FPGA;二維傅里葉變換;頻譜移位;圖像處理

二維傅里葉變換和頻譜移位是目前三維成像系統(tǒng)中非常重要的模塊。二維傅里葉變換可以分解成兩次一維傅里葉變換,先行變換再列變換[1-2]。二維傅里葉變換得到的圖像頻譜的零頻是位于矩陣的4個(gè)頂角,為了方便進(jìn)行頻譜分析,需要進(jìn)行頻譜移位將零頻移到矩陣的中心。頻譜移位需要分別沿著行和列兩個(gè)方向進(jìn)行兩次移位。傅里葉變換和頻譜移位都需要進(jìn)行兩次變換,消耗大量的時(shí)間[3-5]。

本文基于Altera公司的FFT IP 核,提出了一種結(jié)合頻譜移位的二維傅里葉變換的FPGA實(shí)現(xiàn)方法。通過改變二維傅里葉變換的列變換起始位置,從矩陣的中間列開始進(jìn)行列變換,使得二維傅里葉變換頻譜的零頻從矩陣的4個(gè)頂角移到中間列,在完成二維傅里葉變換的同時(shí)完成了一半的頻譜移位。

1 二維傅里葉變換和頻譜移位的介紹

1.1 二維傅里葉變換

數(shù)字圖像是一個(gè)時(shí)域上的二維函數(shù),它包含了周期性成分、非周期性成分、噪聲等信息,這些成分往往緊密交織在一起,很難進(jìn)行分離和分析。傅里葉變換可以將圖像時(shí)域的灰度分布函數(shù)變換為圖像的頻率分布函數(shù),得到了不同頻率的周期函數(shù),可以更加方便地對(duì)這些頻域信號(hào)進(jìn)行處理、分析。

一維的離散傅里葉變換(DFT)可以用式(1)表示,逆傅里葉變換用式(2)表示,f(x)表示時(shí)域上的一維數(shù)字信號(hào),x=0,1,…N-1[5-7]。

(1)

(2)

一維離散的傅里葉變換可以擴(kuò)展到二維,用f(x,y)表示數(shù)字圖像或大小為M×N的二維離散信號(hào),u和v是離散的頻域變量,x和y是離散的時(shí)域變量,M和N是二維矩陣的行數(shù)和列數(shù)。本文研究的是一個(gè)256像素×256像素大小的圖片,所以令M=N=256。二維離散傅里葉變換對(duì)可以用式(3)、式(4)表示。

(3)

(4)

從式(3)可以看出二維的離散傅里葉變換可以分解成兩次一維的傅里葉變換,先對(duì)f(x,y)進(jìn)行一維的行傅里葉變換得到f(u,y),然后再對(duì)f(u,y)進(jìn)行一維的列傅里葉變換得到F(u,v)。同時(shí),每一次的變換都是只有一行數(shù)或者一列數(shù)獨(dú)立地進(jìn)行變換,所以無論從哪一列開始進(jìn)行變換,都不會(huì)改變得到的頻譜信息,只是改變了頻譜的位置。

直接使用DFT的定義計(jì)算的計(jì)算復(fù)雜度為O(N2),即對(duì)N點(diǎn)的f(x)進(jìn)行傅里葉變換需要進(jìn)行N2次的復(fù)數(shù)乘法。實(shí)現(xiàn)一次復(fù)數(shù)乘法需要4次實(shí)數(shù)乘2次實(shí)數(shù)加,當(dāng)N很大時(shí),其計(jì)算量是相當(dāng)可觀,對(duì)于2-D圖像處理,所需計(jì)算量更是大得驚人。本文實(shí)現(xiàn)的256×256像素大小的圖像的二維傅里葉變換,則需要進(jìn)行256×256×256×2=33 554 432次復(fù)數(shù)乘法。本文使用的FFTip是基于快速傅里葉變換算法實(shí)現(xiàn)的,快速傅里葉變換(FFT)是離散傅里葉變換的一種高效算法,可以將復(fù)雜度降低為O(N/2log2N),對(duì)于256×256大小的圖像的二維傅里葉變換只需要進(jìn)行262 144次復(fù)數(shù)乘法。隨著處理圖像的大小的增加,快速傅里葉變換的優(yōu)勢(shì)會(huì)更加明顯。本文設(shè)計(jì)中使用的FFT IP核是高性能、高參數(shù)化的快速傅里葉變換(FFT)處理器[8-9]。

1.2 頻譜移位

數(shù)字圖像的二維傅里葉變換的頻譜的零頻分別位于矩陣的4個(gè)頂角,如圖1(a)所示。為了更方便地對(duì)頻譜信息進(jìn)行處理和分析,需要通過頻譜移位把零頻集中到矩陣的中央,如圖1(c)所示。

圖1 二維圖像頻譜移位圖

通常用軟件對(duì)圖像進(jìn)行處理時(shí)是用MALTAB中的FFT2+fftshift()命令對(duì)原始圖像進(jìn)行二維傅里葉變換和頻譜移位的[4-5]。如圖1所示,把圖像矩陣等分成4個(gè)小矩陣,頻譜移位(fftshift)是實(shí)現(xiàn)矩陣1和4進(jìn)行平移交換,矩陣2和3進(jìn)行平移交換。在硬件上,頻譜移位是通過對(duì)頻譜信息的存儲(chǔ)和地址變換來實(shí)現(xiàn)。首先對(duì)圖1(a)的頻譜的每一行進(jìn)行移位,得到圖1(b)的頻譜圖,然后對(duì)圖1(b)的頻譜圖的每一列進(jìn)行移位,最后得到圖1(c)所示的頻譜圖,完成了頻譜移位。頻譜移位的硬件實(shí)現(xiàn)需要消耗大量的時(shí)間。

通過1.1節(jié)對(duì)二維傅里葉變換的分析可知,改變列變換的順序不會(huì)改變頻譜信息,所以考慮對(duì)圖像的二維傅里葉變換的列變換進(jìn)行改進(jìn),從第N/2列開始依次進(jìn)行第(N/2)列至N列變換,然后從第0列開始依次進(jìn)行第0列至(N/2-1)列變換。這樣完成二維傅里葉變換之后得到的頻譜如圖1(b)所示,零頻已經(jīng)移到了矩陣的中間列。所以在二維傅里葉變換的完成的同時(shí)已經(jīng)完成了一半的頻譜移位,縮短了一半的頻譜移位時(shí)間。

2 結(jié)合頻譜移位的二維傅里葉變換的FPGA實(shí)現(xiàn)

二維傅里葉變換系統(tǒng)的FPGA實(shí)現(xiàn)框圖如圖2所示。clk是系統(tǒng)時(shí)鐘,reset_n是系統(tǒng)的復(fù)位信號(hào),start_ac為二維傅里葉變換的開始信號(hào)。主要的模塊為FFT ip核模塊、SRAM控制模塊和RAM控制模塊。

原始圖片存儲(chǔ)在SRAM中,像素信息按行為單元送入到FFT模塊進(jìn)行行變換,行變換的中間結(jié)果存入到RAM中,然后再從RAM中按列為單元讀出數(shù)據(jù)送入FFT模塊進(jìn)行列變換,通過兩次的一維傅里葉變換實(shí)現(xiàn)了二維傅里葉變換。

圖2 二維傅里葉變換的硬件框圖

2.1 FFT ip核模塊

快速傅里葉變換在現(xiàn)代數(shù)字信號(hào)處理系統(tǒng)中是最常見的處理算法之一。Altera FFT IP核是高性能、高參數(shù)化的快速傅里葉變換(FFT)處理器,可以實(shí)現(xiàn)高性能的復(fù)數(shù)FFT或反FFT(IFFT)。FFT處理器的實(shí)現(xiàn)結(jié)構(gòu)有緩存、突發(fā)、流和可變流結(jié)構(gòu)。為了實(shí)現(xiàn)自然數(shù)順序的輸入和自然數(shù)順序的輸出,本文使用的是定點(diǎn)的可變流結(jié)構(gòu),數(shù)據(jù)變換長(zhǎng)度為256,輸入數(shù)據(jù)的精度為16 bit,旋轉(zhuǎn)因子的精度為12 bit,輸出數(shù)據(jù)的精度為25 bit[10]。旋轉(zhuǎn)因子的精度越高,傅里葉變換的計(jì)算結(jié)果越精確。

FFT IP核接口信號(hào)如圖3所示。左邊部分信號(hào)是輸入控制信號(hào),右邊部分是輸出信號(hào)。

圖3 FFT IP核接口信號(hào)

主要的接口信號(hào)的介紹:

Sink_real,Sink_imag:輸入數(shù)據(jù)的實(shí)部和虛部;

Sink_valid:輸入數(shù)據(jù)有效信號(hào),當(dāng)該信號(hào)為1時(shí),FFT處理器認(rèn)為輸入的Sink_real和Sink_imag信號(hào)有效,所以可以通過控制該信號(hào)的值來控制FFT的運(yùn)算的進(jìn)行和暫停。

source_real,source_imag:輸出數(shù)據(jù)的實(shí)部和虛部;

source_valid:輸出數(shù)據(jù)有效信號(hào),當(dāng)該信號(hào)為1的同時(shí)就要控制RAM寫數(shù)據(jù)。

2.1.1 行傅里葉變換

行傅里葉變換,可以看成是虛部為0的復(fù)數(shù)的傅里葉變換。在50 MHz的時(shí)鐘頻率下從SRAM中按順序讀出原始數(shù)據(jù)的實(shí)部(Sink_real),虛部(Sink_imag)直接賦值為0。行變換結(jié)果的實(shí)部和虛部分別存儲(chǔ)在兩個(gè)RAM中。

2.1.2 列傅里葉變換

從RAM按列的順序從第129列開始依次讀出第129列至255列的數(shù)據(jù),再從讀出第0列至128列的數(shù)據(jù)送入傅里葉變換模塊進(jìn)行列變換。

3 系統(tǒng)設(shè)計(jì)結(jié)果

系統(tǒng)針對(duì)256×256大小的圖像完成二維傅里葉變換的過程,設(shè)計(jì)采用Verilog源代碼以及 Altera公司的EP4CE115F29C7芯片進(jìn)行原型驗(yàn)證,采用50 MHz為系統(tǒng)時(shí)鐘,在Altera開發(fā)環(huán)境Quartus Ⅱ 13.0 SPI中綜合,綜合結(jié)果如表1所示。

表1 綜合分析報(bào)告

采用SignalTap Ⅱ logic Analyser 邏輯分析儀工具,以50 MHz時(shí)鐘作為采樣時(shí)鐘,對(duì)該二維傅里葉變換模塊進(jìn)行了實(shí)時(shí)驗(yàn)證。圖4是硬件驗(yàn)證的結(jié)果,圖5是MATLAB軟件仿真得到的結(jié)果,圖6是兩者的二維傅里葉變換的第129列的頻譜的對(duì)比圖,從圖6可以看到兩者的結(jié)果基本完全重合。誤差在1%以內(nèi),誤差存在的主要原因是MALTAB軟件計(jì)算是用雙精度浮點(diǎn)型數(shù)據(jù)類型進(jìn)行計(jì)算的,而硬件實(shí)現(xiàn)是基于二進(jìn)制計(jì)算的,所以精確度也會(huì)降低,這樣的微小的誤差不可消除,始終存在。

圖4 FPGA實(shí)時(shí)驗(yàn)證得到的二維傅里葉變換結(jié)果(source_real_cut和source_imag_cut分別是得到的頻譜的實(shí)部和虛部)

圖5 MATLAB軟件仿真得到的二維傅里葉變換結(jié)果(黑色方框中數(shù)據(jù)是得到的頻譜的第129列的值)

圖6 FPGA和MATLAB計(jì)算結(jié)果對(duì)比圖

硬件實(shí)現(xiàn)一開始得到的結(jié)果和軟件仿真結(jié)果的第129列一致,說明該硬件系統(tǒng)實(shí)現(xiàn)了從圖像的第129列開始進(jìn)行列傅里葉變換,在完成二維傅里葉變換的同時(shí)把頻譜的零頻移到了矩陣的中間列。設(shè)計(jì)是通過改變列變換的起始列來實(shí)現(xiàn)頻譜移位的,這種設(shè)計(jì)方法在完成一半的頻譜移位的同時(shí)也不會(huì)給原來的系統(tǒng)增加額外的邏輯單元,不會(huì)影響原來系統(tǒng)的時(shí)序。

實(shí)現(xiàn)二維傅里葉變換和頻譜移位均大約需要256×256×2=131 072個(gè)時(shí)鐘周期的時(shí)間,本文設(shè)計(jì)的二維傅里葉變換系統(tǒng)可以節(jié)省一半的頻譜移位的時(shí)間,節(jié)省二維傅里葉變換和頻譜移位系統(tǒng)總耗時(shí)的25%。按50 MHz的頻率計(jì)算,可以節(jié)省1.3 ms。實(shí)際的圖像處理系統(tǒng)中,需要處理的圖像的幀數(shù)非常多,圖像的大小也可能會(huì)很大,如果處理的是大小為1 024×1 024的50幀的圖像,就會(huì)節(jié)省大約1 s的時(shí)間,是非常可觀的。

5 總結(jié)

本文結(jié)合頻譜移位完成了256×256像素大小的圖像的二維傅里葉變換的FPGA實(shí)現(xiàn),并實(shí)時(shí)驗(yàn)證了該系統(tǒng)的正確性,最高工作頻率達(dá)84.2 MHz。通過改變傅里葉變換的起始列,從行傅里葉變換的第129列開始進(jìn)行列傅里葉變換,把頻譜移位結(jié)合到二維傅里葉變換中,使得在完成二維傅里葉變換的同時(shí)也完成了一半的頻譜移位,縮短了一半的頻譜移位的時(shí)間,大大縮短圖像處理系統(tǒng)的運(yùn)算時(shí)間,為圖像處理的硬件實(shí)現(xiàn)技術(shù)中頻譜分析和處理提供了一種新的方法的思路。

[1] 楊軍,于艷艷,陳成.基于FPGA的二維FFT處理器的研究與設(shè)計(jì)[J]. 云南大學(xué)學(xué)報(bào)(自然科學(xué)版),2013,35(06):750-755.

[2] Hojin Kee,Shuvra S. Bhattacharyya,Newton Petersen and Jacob KornerupR esource-Efficient Acceleration of 2-Dimensional Fast Fourier Tran-Form Computation on FPGAs[C]//Third ACM/IEEE International Conferenceon Distributed Smart Cameras.Como:IEEE,2009:1-8.

[3] 孫文昌,李忠科,楊威. 圖像的傅里葉變換頻譜特性分析[J]. 計(jì)算機(jī)與信息技術(shù),2005(6):46-48.

[4] Abdellah M,Saleh S,Eldeib A,et al. High Performance Multi-Dimensional(2D/3D)FFT-Shift Implementation on Graphics Processing Units(GPUs)[C]//2012 Cairo International on Biomedical Engineering Conference. Giza:IEEE,2012:171-174.

[5] Zavadil J,Tuma J,Valicek J,et al. Two Dimensional Fourier Transform Using MATLAB[C]//2013 14th on International Carpathian Control Conference,Rytro:IEEE,2013:432-435.

[6] 劉微,朱明,李向榮,等. 運(yùn)動(dòng)模糊圖像實(shí)時(shí)恢復(fù)的多DSP方案[J]. 電子器件,2005,28(4):722-725.

[7] 胡廣書. 數(shù)字信號(hào)處理:理論、算法與實(shí)現(xiàn)[M]. 北京:清華大學(xué)出版社,2012:166-184.

[8] Anitha T G,Ramachandran S. Novel Algorithms for 2-D FFT and Its Inverse for Image Compression[C]//2013 International Conference on Signal Processing Image Processing and Pattern Recognit-Ion. Coimbatore:IEEE,2013:62-65.

[9] Vinay Kumar M,David Selvakumar A,Sobha P M. Area and Frequency Optimized 1024 Point Radix-2 FFT Processor on FPGA[C]//2015 International Conference on Vlsi Systems,Archit-Ecture,Technology and Applications. Bangalore:IEEE. 2015:1-6.

[10] 劉東華. Altera系列FPGA芯片IP核詳解[M]. 北京:電子工業(yè)出版社,2014:172-286.

ImplementationofTwo-DimensionalFourierTransformCombinedwithSpectrum-ShiftingBasedonFPGA

QIANHui1,2,SHIYao1,2,GONGMin1,2,GAOBo1,2*

(1.Micro-Electronics of College of Physical Science and Technology,Chengdu 610064,China;2.Key Laboratory of Micro-Electronics Technology of Sichuan Province,Chengdu 610064,China)

In 3d imaging technology,after the Fourier transform of the image zero frequency will need to be concentrated to the center position by spectrum-shift. In order to improve the speed of image processing,a method of the implementation of two-dimensional Fourier transform combining with spectrum-shift is proposed based on FPGA. Combing the spectrum-shift with the two-dimensional Fourier transform hardware systems,the two-dimensional Fourier transform of the image and the half of the spectrum-shift can be achieved at the same time. The Altera chip of EP4CE115F29C7 was used in this design and the hardware implementation was completed for 256×256 of the image. The highest working frequency reached to 84.2 MHz. The amount of the resource usage reached to 3 849 LEs. The SignalTap Ⅱ Logic Analyzer has verified the correctness of the module realtimely.

FPGA;two-dimensional Fourier transform;spectrum-shift;image processing

10.3969/j.issn.1005-9490.2017.05.009

2016-08-19修改日期2016-10-23

TN402

A

1005-9490(2017)05-1092-05

錢輝(1992-),女,漢族,安微蕪湖人,四川大學(xué)物理學(xué)院微電子專業(yè),碩士研究生,研究方向?yàn)槌笠?guī)模集成電路設(shè)計(jì);

史瑤(1991-),男,漢族,四川成都人,四川大學(xué)物理學(xué)院微電子專業(yè),碩士研究生,研究方向?yàn)槌笠?guī)模集成電路設(shè)計(jì);

龔敏(1961-),男,教授(導(dǎo)師),博士生導(dǎo)師,從事新型半導(dǎo)體材料與器件工藝、集成電路設(shè)計(jì)和工藝及半導(dǎo)體器件的輻照效應(yīng)研究;

高博(1975-),男,副教授,主要從事CMOS集成電路芯片設(shè)計(jì)和生物醫(yī)學(xué)成像領(lǐng)域的研究。