fft方法的matlab實(shí)現(xiàn)

1、FFTUniversity of Science and Technology of Beijing沈政偉v一，F(xiàn)ourier 級(jí)數(shù)v二，連續(xù)FourierTransformv三，一維離散傅立葉變換（Discrete Fourier Transform）.v一，F(xiàn)ourier 級(jí)數(shù)v法國著名科學(xué)家傅立葉在1807年向法國國家科學(xué)院提交的一篇報(bào)告中提出：“任何周期函數(shù)都可以用一系列正弦波（諧波）來線性表示”Fourier 級(jí)數(shù)考慮正弦波，顯然該函數(shù)周期為 ，對(duì)應(yīng)的頻率為 而一般的樂器發(fā)出的聲音以及電壓等信號(hào)都可以通過具有不同頻率的正弦波函數(shù)疊加來表示。比如：信號(hào)在持續(xù)為 的時(shí)間內(nèi)分別振動(dòng)次數(shù)為

2、1,20,100次，而其中頻率為1的分量振幅最大，達(dá)到100，從而具有決定性作用。 )sin(ktk2k100sin( )3sin(20 )0.5sin(100 )ttt2v v在一般情形下，信號(hào)可以用下面正弦波的無窮和形式來進(jìn)行分解。v在這個(gè)公式中，通過由Riemann-Lebesgue引理 知道：（為什么？）v這表明一般實(shí)際信號(hào)（能量有限）中的高頻成分會(huì)隨著頻率的增大，相應(yīng)的變小，信號(hào)中的主要成分為少數(shù)系數(shù)（頻率低的部分）所控制。而信號(hào)的高頻成分對(duì)應(yīng)著信號(hào)中的細(xì)節(jié)部分，而低頻成分對(duì)應(yīng)著信號(hào)的主要的信息)(tfkkkktbktaatf)sin()cos()(00limlimkkkkbav二，

3、連續(xù)FourierTransformv前面我們討論的函數(shù)傅立葉級(jí)數(shù)都是在具有有限周期（尤其是以 為周期的函數(shù)）情形下進(jìn)行展開的，但是無論是在理論上還是在實(shí)踐中，對(duì)于非周期函數(shù)性質(zhì)的討論都是非常有必要的。v如何從周期函數(shù)的傅立葉級(jí)數(shù)展開擴(kuò)展到非周期函數(shù)的傅立葉“級(jí)數(shù)”性質(zhì)？或者是周期函數(shù)可以展開成為傅立葉級(jí)數(shù)的形式，那么是否非周期函數(shù)也可以展開成為傅立葉級(jí)數(shù)的形式呢？2v一維傅立葉變換的定義v定理：設(shè)函數(shù)設(shè)函數(shù) ，即上的分即上的分段光滑函數(shù)，那么這個(gè)函數(shù)的傅立葉變換定義為：段光滑函數(shù)，那么這個(gè)函數(shù)的傅立葉變換定義為：v或者或者v條件部分也可以這樣描述：條件部分也可以這樣描述：設(shè)為實(shí)變量設(shè)為實(shí)變量

4、x的的連續(xù)且可積函數(shù)連續(xù)且可積函數(shù)1( )( )f xL Rdxxf）（1( )( )2ixff x edx2( )( )ixff x edx ( )f xv這種定義是有意義的，并且函數(shù)可以通過其這種定義是有意義的，并且函數(shù)可以通過其傅立葉反演而得到：傅立葉反演而得到：v或者或者v在上面的式子中知道：，并且會(huì)注意到實(shí)在上面的式子中知道：，并且會(huì)注意到實(shí)函數(shù)的傅立葉變換通常是復(fù)數(shù)形式函數(shù)的傅立葉變換通常是復(fù)數(shù)形式v)(xf1( )()2ixf xfed2( )( )ixf xfed1i )(xf( )( )( )fRiIv指數(shù)形式為：v其中：為幅值函數(shù)，稱為函數(shù)（信號(hào)）的傅立葉譜而 稱為相角傅立

5、葉譜的平方，稱為是能量譜或者是功率譜()()()iffe22()()()fRI()fx1( )( )tan( )IR 22( )( )( )ERI題目：求門函數(shù)的傅立葉變換)(xf0( )0AxXf x其他v結(jié)果是：v其中用到了歐拉公式所以門函數(shù)所對(duì)應(yīng)的傅立葉譜是：( )sin()iXAfX e)(xfsin()( )XfAXXcossini tetit-TT0Atx(t)X()02ATv二維傅立葉變換的定義v設(shè)函數(shù)設(shè)函數(shù) 是連續(xù)可積的，即即那么這個(gè)函數(shù)的傅立葉變換定義為：那么這個(gè)函數(shù)的傅立葉變換定義為：( , )f x y,Rf x y dxdy （）2()(,)(,)juxvyF u vf

6、x y edxdy 2()1222122(,)(,)(,)(,)(,)(,)(,)ta n(,)(,)(,)(,)ju xv yfxyFu ved u d vFu vRu vIu vIu vu vRu vEu vRu vIu v 幅 值相 角能 量v傅立葉變換的性質(zhì)v（），線性性v（），可分離性v（），平移性v（），共軛性v（），尺度變換特性v（），卷積定理v（），arseval定理221( )()2f xdtfdReal Part, Imaginary Part, Magnitude, Phase, SpectrumReal part:Imaginary part:Magnitude-pha

7、se representation:Magnitude(spectrum):Phase(spectrum):PowerSpectrum:Mean of image/ DC component:Highest frequency component:“Half-shifted”Image:Conjugate Symmetry:Magnitude Symmetry:2D DFT Properties2D DFT PropertiesSpatial domain differentiation:Frequency domain differentiation:Laplacian:v三，一維離散傅立葉

8、變換（Discrete Fourier Transform）.v離散傅立葉變換的優(yōu)點(diǎn)：（）比在時(shí)域直接對(duì)數(shù)字信號(hào)進(jìn)行處理所需要的運(yùn)算量要?。唬ǎ┚哂锌焖偎惴╲如果將一維連續(xù)函數(shù)用取個(gè)間隔取樣增量的方法進(jìn)行離散化，變?yōu)殡x散函數(shù)v或者是：)( xfx)(xf0000 ( ), (), (2 ),(1) )f xf xx f xxf xNx 0( )()0,1 ,21f nf x n xnN v那么取樣之后一維離散函數(shù)的離散傅立葉變換定義為：v反變換定義為：()fn1012( )( ) expNniunF ufnNN102( )( ) expNuiunfxF uN0,1, 2,10,1, 2,1nN

9、uN其 中 ，v注意（）上面兩個(gè)公式是對(duì)離散后函數(shù)的準(zhǔn)確的離散傅立葉變換；（）也是一個(gè)取個(gè)等量間隔取樣之后的離散函數(shù)，它可以表示成為并且可以證明函數(shù)在空間域和頻率域取樣間隔和之間的關(guān)系是：（）離散傅立葉變換總是存在的，它并不需要考慮連續(xù)傅立葉變換所需要的可積的條件( )f n0()()FuFunuu( )F u)(xfx1uN x uv一維的v如果對(duì)：v令v那么，上述的公式變?yōu)椋簐1012( )( )expNniunF uf nNN2expiWN102( )( ) expNuiunfnF uN101()()Nu nnFufn WN101( )( )NunufnF u WNv寫成矩陣的形式為：2

10、00000121012(1)(1)(0)(0)(1)(1)1(2)(2)(1)(1)NNNNN NWWWWFfWWWWFfFfNF Nf NWWWW2000001210(1)2(1)(1)(0)(0)(1)(1)(2)(2)(1)(1)NNNNN NWWWWfFWWWWfFfFf NF NWWWW A 4x4 imagejjjjjjjj111111111111336632452889863111111111111144XFFX Compute its 2D-DFT:3366324528898631Xjjjjjjjjjjjj111111111111554213463795542134161921

11、21jjjjjjjjjjjj81174459413611136134574811945235277MATLAB function: fft2lowest frequency componenthighest frequency componentjjjjjjjjjjjj81174459413611136134574811945235277XReal part:114546116135411423277realX874913013047895050imagX60.1306. 840. 685. 932.141132.14134 . 606. 860.1385. 939. 5339. 577mag

12、nitudeX51. 209. 247. 215. 100. 214. 300. 214. 347. 209. 251. 215. 119. 1019. 10phaseXImaginary part:Magnitude:Phase:Computation of 2D-DFT: Examplejjjjjjjjjjjjjjjjjjjj1111111111118117445941361113613457481194523527711111111111141244*FXF Compute the inverse 2D-DFT:X3366324528898631jjjjjjjjjjjj554213463

13、7955421341619212111111111111141MATLAB function: ifft2Centered RepresentationFrom Gonzalez & Woodsvuhighhighhighhigh(-N/2, N/2)(N/2, N/2)low(N/2, -N/2)(-N/2, -N/2)From Prof. Al BovikMATLAB function: fftshiftExample:Log-Magnitude VisualizationFrom Gonzalez & Woods)1log(rcs2D-DFTcenteredApply t

14、o Images2D-DFT centered log intensity transformationv并不是一種新的變換方式，它是離散傅立葉變換的一種算法，這種方法是建立在分析離散傅立葉變換中的多余運(yùn)算的基礎(chǔ)上，進(jìn)而消除這些重復(fù)多余運(yùn)算的思想下得到的，從而減少了運(yùn)算，節(jié)省了時(shí)間，達(dá)到快速運(yùn)算的目的v是以的組成狀況而分為許多種算法，我們?cè)谶@里只介紹為的整數(shù)次冪的算法v令v那么，一維離散傅立葉變換公式為：2e x pNiWN101()()Nu nNnFufn WN2 (0,1, 2,),0 1nNnn u其 中當(dāng) 然分 別 為 ， ， 2，N-1vFFT是利用 的兩個(gè)特點(diǎn)來提高計(jì)算效率的。v（

15、1） 的對(duì)稱性 （2） 的周期性nNWnNW()*()NnnnNNNWWWnNW()nn NNNWWv在此基礎(chǔ)上，將分解成為和對(duì)應(yīng)的偶數(shù)和奇數(shù)兩部分，n的取值范圍由原來的到變?yōu)榈絭令v因此，離散傅立葉變換可以寫成為下面的形式：(2 )fn(21)fn (0,1,2,3,1)2( )(21)gNnh nfn（n）=f(2n)( )f n1112200011222(21)0011220022( )( )( )( )(2 )(21)(2 )(21).NNNunununNNNnnnNNununNNnnNNununuNNNnnF uf n Wg n Wh n Wfn WfnWfn WfnWW111220

16、0011222(21)0011220022( )( )( )( )(2 )(21)(2 )(21).( ).( )NNNunununNNNnnnNNu nunNNnnNNununuNNNnnuNF uf n Wg n Wh n Wfn WfnWfn WfnWWG uWH uv其中的v因?yàn)樗詖所以，我們得到2NuuNNWW 21NNW0 ,1,12Nu ()( ).( )2uNNF uG uWH u()()2NGuGu()( )2NH uH uv因此，一個(gè)求點(diǎn)的離散傅立葉變換可以被分為兩個(gè)求點(diǎn)的離散傅立葉變換如果取，則上式變?yōu)椋?818283808182838(0 )(0 ).(0 )(1)(

17、1).(1)( 2 )( 2 ).( 2 )(3)(3).(3)( 4 )(0 ).(0 )(5 )(1).(1)(6 )( 2 ).( 2 )(7 )(3).(3)FGWHFGWHFGWHFGWHFGWHFGWHFGWHFGWH這是上面公式的蝶型運(yùn)算圖v 由于 和 的過程仍然是常規(guī)的離散傅立葉變換,因此,也可以按照上述規(guī)則進(jìn)行分解,記 為計(jì)算 的DFT運(yùn)算,記 為計(jì)算 的DFT計(jì)算,用 表示 ， 為計(jì)算 的DFT計(jì)算 (4)H(4)G( )A u(2 )Gn( )B u(21)Gn ( )C u(2 )Hn( )D u(21)Hn08280828(0)(0). (0)(1)(1). (1)(

18、2)(0). (0)(3)(1). (1)GAW BGAW BGAW BGAW B08280828(0)(0). (0)(1)(1). (1)(2)(0). (0)(3)(1). (1)HCW DHCW DHCW DHCW Dv由于 現(xiàn)在都是兩點(diǎn)的DFT了，因此不需要再分解，可以直接從原始數(shù)據(jù)計(jì)算得到( )A u( )Bu( )C u( )Du0808(0)(0). (4)(1)(0). (4)AfWfAfWf0808(0)(2). (6)(1)(2). (6)BfWfBfWf0808(0)(1). (5)(1)(1). (5)CfWfCfWf0808(0)(3). (7)(1)(3). (7

19、)DfWfDfWfvFFT的實(shí)現(xiàn)v二維離散Fourier變換v Mean of image/ DC component:Highest frequency component:“Half-shifted”Image:Conjugate Symmetry:Magnitude Symmetry:2D DFT PropertiesvBasics of filtering in the frequency domainv f (x,y)h(x,y)g(x,y)( , )( , )( , )g x yf x yh x yG(u,v)H(u,v)F(u,v) =DFTIDFTDFTIDFTDFTIDFTi

20、nput imageimpulse response(filter)output imagev Ideal lowpass, bandpass and highpasslow-frequency maskmid-frequency maskhigh-frequency maskvSome basic filters and their propertiesFrom Gonzalez & WoodsIdeal lowpass filtering with cutoff frequencies set at radii values of 5, 15, 30, 80, and 230, r

21、espectively2D-DFT Domain Filter DesignGaussian Lowpass FiltervD(u,v): distance from the origin of FTvParameter: =D0 (cutoff frequency)vThe inverse FT of the Gaussian filter is also a GaussianH(u,v) eD2(u,v)/22Image Enhancement in theFrequency Domain using GLPFImage Enhancement in theFrequency Domain

22、Butterworth Filter (lowpass)vHigh-frequency emphasis: Adding a constant to a highpass filter to preserve the low-frequency components.nvuDDvuH20),(/11),(（A）幅頻特性 （B）相頻特性由圖可見，隨n的增大，幅頻特性在截止頻率處下降得越快，則越接近于理想低通濾波器。Sharpening (Highpass) FilteringvImage sharpening can be achieved by a highpass filtering pro

23、cess, which attenuates the low-frequency components without disturbing high-frequency information.v Zero-phase-shift filters: radially symmetric and completely specified by a cross section.Hhp(u,v) 1 Hlp(u,v)Ideal Filter (Highpass)vThis filter is the opposite of the ideal lowpass filter.00),( if 1),( if 0),(DvuDDvuDvuHButterworth Filter (Highpas