數(shù)字信號處理-實(shí)驗(yàn)-Matlab實(shí)驗(yàn)FFT算法實(shí)現(xiàn)

實(shí)驗(yàn)二應(yīng)用FFT對信號進(jìn)行頻譜分析實(shí)驗(yàn)?zāi)康?、加深對離散信號的DTFT和DFT的及其相互關(guān)系的理解。2、在理論學(xué)習(xí)的基礎(chǔ)上，通過本次實(shí)驗(yàn)，加深對快速傅里葉變換的理解，熟悉FFT算法極其程序的編寫。3、熟悉應(yīng)用FFT對典型信號進(jìn)行頻譜分析的方法。4、了解應(yīng)用FFT進(jìn)行信號頻譜分析過程中可能出現(xiàn)的問題，以便在實(shí)際中正確應(yīng)用FFT。實(shí)驗(yàn)原理和方法一個(gè)連續(xù)信號的頻譜可以用它的傅里葉變換表示（2—1）如果對該信號進(jìn)行理想采樣，可以得到采樣序列：（2—2）同樣可以對該序列進(jìn)行Z變換，其中T為采樣周期（2—3）當(dāng)?shù)脮r(shí)候，我們就得到了序列的傅里葉變換（2—4）其中稱為數(shù)字頻率，它和模擬頻域的關(guān)系為（2—5）式中的是采樣頻率，上式說明數(shù)字頻率是模擬頻率對采樣頻率的歸一化。同模擬域的情況相似，數(shù)字頻率代表了序列值變化的速率，而序列的傅里葉變換稱為序列的頻譜。序列的傅里葉變換和對應(yīng)的采樣信號頻譜具有下式的對應(yīng)關(guān)系：（2—6）即序列的頻譜是采樣信號頻譜的周期延拓。從式（2—6）可以看出，只要分析采樣序列的頻譜，就可以得到相應(yīng)的連續(xù)信號的頻譜。注意：這里的信號必須是帶限信號，采樣也必須滿足Nyquist定理。在各種信號序列中，有限長序列在數(shù)字信號處理中占有很重要的地位。無限長的序列也往往可以用有限長序列來逼近。對于有限長的序列我們可以使用離散傅里葉變換（DFT），這一變換可以很好的放映序列的頻域特性，并且容易利用快速算法在計(jì)算機(jī)上實(shí)現(xiàn)當(dāng)序列的長度是N時(shí)，我們定義離散傅里葉變化為：（2—7）其中，，它的反變換定義為：（2—8）根據(jù)式（2—3）和（2—7）令，則有（2—9）可以得到，是Z平面單位圓上幅角為的點(diǎn)，就是見單位圓進(jìn)行N等分以后第K個(gè)點(diǎn)。所以，是變換在單位圓上的等距采樣，或者說是序列傅里葉變換的等距采樣。時(shí)域采樣在滿足Nyquist定理時(shí)，就不會(huì)發(fā)生頻譜混淆；同樣地，在頻率域進(jìn)行采樣的時(shí)候，只要采樣間隔足夠小，也不會(huì)發(fā)生時(shí)域序列混淆。DFT時(shí)對序列傅里葉變換的等距采樣，因此可以用于序列的頻譜分析。在運(yùn)用DFT進(jìn)行頻譜分析的時(shí)候可能有三種誤差，分析如下：（1）混淆現(xiàn)象從式（2—6）中可以看出，序列的頻譜是采樣信號頻譜的周期延拓，周期是，因此當(dāng)采樣速率不滿足Nyquist定理，即采樣頻率小于兩倍的信號（這里指的是實(shí)信號）頻率時(shí)，經(jīng)過采樣就會(huì)發(fā)生頻譜混淆。這導(dǎo)致采樣后的信號序列頻譜不能真實(shí)的反映原信號的頻譜。所以，在利用DFT分析連續(xù)信號頻譜的時(shí)候，必須注意這一問題。這就告訴我們，在確定信號的采樣頻率之前，需要對頻譜的性質(zhì)有所了解。在一般的情況下，為了保證高于折疊頻率的分量不會(huì)出現(xiàn)，在采樣之前，先用低通模擬濾波器對信號進(jìn)行濾波。（2）泄露現(xiàn)象實(shí)際中的信號序列往往很長，甚至是無線長序列。為了方便，我們往往用截短的序列來近似它們，這樣可以使用較短的DFT來對信號進(jìn)行頻譜分析。這種截短等價(jià)于給原始信號序列乘以一個(gè)矩形窗函數(shù)，而矩形窗函數(shù)的頻譜不是有限帶寬的，從而它和原信號的頻譜進(jìn)行卷積以后會(huì)擴(kuò)展原信號的頻譜。值得一提的是，泄漏時(shí)不能和混淆完全分離的，因?yàn)樾孤?dǎo)致頻譜的擴(kuò)展，從而造成混淆。為了減小泄漏的影響，可以選擇是黨的窗函數(shù)使頻譜的擴(kuò)散減小到最小。柵欄效應(yīng)因?yàn)镈FT是對單位圓上Z變換的均勻采樣，所以它不可能將頻譜視為一個(gè)連續(xù)函數(shù)。這樣就產(chǎn)生了柵欄效應(yīng)，從某種角度來看，用DFT來觀看頻譜就好像通過一個(gè)柵欄來觀看一幅景象，只能在離散點(diǎn)上看到真實(shí)的頻譜。這樣的話就會(huì)有一些頻譜的峰點(diǎn)或谷點(diǎn)被“柵欄”擋住，不能被我們觀察到。減小柵欄效應(yīng)的一個(gè)方法是在原序列的末端補(bǔ)一些零值，從而變動(dòng)DFT的點(diǎn)數(shù)。這種方法的實(shí)質(zhì)是人為地改變了對真是頻譜采樣的點(diǎn)數(shù)和位置，相當(dāng)于搬動(dòng)了“柵欄”的位置，從而使得原來被擋住的一些峰點(diǎn)或谷點(diǎn)顯露出來。注意，這時(shí)候每根譜線對應(yīng)的頻率和原來的已經(jīng)不相同了。從上面的分析過程可以看出，DFT可以用于信號的頻譜分析，但必須注意可能產(chǎn)生的誤差，在應(yīng)用過程中要盡可能減小和消除這些誤差的影響?？焖俑道锶~變換FFT并不是與DFT不相同的另一種變換，而是為了減少DFT運(yùn)算次數(shù)的一種快速算法。它是對變換式（2—7）進(jìn)行一次次的分解，使其成為若干小數(shù)點(diǎn)DFT的組合，從而減小運(yùn)算量。常用的FFT是以2為基數(shù)，其長度。它的運(yùn)算效率高，程序比較簡單，使用也十分菲娜改變。當(dāng)需要進(jìn)行變換的序列的長度不是2的整數(shù)次方的時(shí)候，為了使用以2為基的FFT，可以用末尾補(bǔ)零的方法，使其長度延長至2的整數(shù)次方。IFFT一般可以通過FFT程序來完成，比較式（2—7）和（2—8），只要對取共軛，進(jìn)行FFT運(yùn)算，然后再取其共軛，并乘以因子，就可以完成IFFT。實(shí)驗(yàn)內(nèi)容及步驟（一）編制實(shí)驗(yàn)用主程序及相應(yīng)子程序1、在實(shí)驗(yàn)之前，認(rèn)真復(fù)習(xí)DFT和FFT有關(guān)的知識，閱讀本實(shí)驗(yàn)原理與方法和實(shí)驗(yàn)附錄部分中和本實(shí)驗(yàn)有關(guān)的子程序，掌握子程序的原理并學(xué)習(xí)調(diào)用方法。2、編制信號產(chǎn)生子程序及本實(shí)驗(yàn)的頻譜分析主程序。試驗(yàn)中需要用到的基本信號包括：（1）高斯序列：（2）衰減正弦序列（3）三角波序列，（4）反三角序列，（二）上機(jī)實(shí)驗(yàn)內(nèi)容1、觀察高斯序列的時(shí)域和頻域特性（1）固定信號的參數(shù)p=8,改變q的值，使q分別等于2，4，8。觀察它們的時(shí)域和幅頻特性，了解q取不同值的時(shí)候，對信號時(shí)域特性和幅頻特性的影響。（2）固定q=8，改變p，使p分別等于8，13，14，觀察參數(shù)p變化對信號序列時(shí)域及幅頻特性的影響。注意p等于多少時(shí)，會(huì)發(fā)生明顯的泄漏現(xiàn)象，繪制相應(yīng)的食欲序列和幅頻特性曲線。1.高斯序列>>n=0:15;p=8;q=8;x=exp(-1*(n-p).^2/q);>>subplot(3,2,2);stem(abs(fft(x)));subplot(3,2,1);stem(n,x);>>p=8;q=4;x=exp(-1*(n-p).^2/q);>>subplot(3,2,3);stem(n,x);subplot(3,2,4);stem(abs(fft(x)));>>p=8;q=2;x=exp(-1*(n-p).^2/q);>>subplot(3,2,5);stem(n,x);subplot(3,2,6);stem(abs(fft(x)));>>n=0:15;p=8;q=8;x=exp(-1*(n-p).^2/q);>>subplot(3,2,2);stem(abs(fft(x)));subplot(3,2,1);stem(n,x);>>p=13;q=8;x=exp(-1*(n-p).^2/q);>>subplot(3,2,3);stem(n,x);subplot(3,2,4);stem(abs(fft(x)));>>p=14;q=8;x=exp(-1*(n-p).^2/q);>>subplot(3,2,5);stem(n,x);subplot(3,2,6);stem(abs(fft(x)));2.衰減正弦序列>>n=0:15;a=0.1;f=0.0625;x=exp(-a*n).*sin(2*pi*f*n);>>subplot(3,2,1);stem(n,x);subplot(3,2,2);stem(n,abs(fft(x)));>>n=0:15;a=0.1;f=0.4375;x=exp(-a*n).*sin(2*pi*f*n);>>subplot(3,2,3);stem(n,x);subplot(3,2,4);stem(n,abs(fft(x)));>>n=0:15;a=0.1;f=0.5625;x=exp(-a*n).*sin(2*pi*f*n);>>subplot(3,2,5);stem(n,x);subplot(3,2,6);stem(n,abs(fft(x)));3．三角波序列>>fori=1:4x(i)=i;end>>fori=5:8x(i)=9-i;end>>subplot(2,1,1);stem(x);subplot(2,1,2);stem(abs(fft(x,16)))4.反三角序列>>fori=1:4x(i)=5-i;end>>fori=5:8x(i)=i-4;end>>subplot(2,1,1);stem(x);subplot(2,1,2);stem(abs(fft(x,16)))四、思考題1.實(shí)驗(yàn)中的信號序列x（n）和x(n),在單位圓上的Z變換頻譜和會(huì)相同嗎？如果不同，你能說出哪一個(gè)低頻分量更多一些嗎？答：不同在單位圓上的Z變換頻譜中，xc(n)的低頻分量比xd(n)的多一些。2.對一個(gè)有限長序列進(jìn)行離散傅里葉變換（DFT）,等價(jià)于將該序列周期延拓后進(jìn)行傅里葉級數(shù)展開。因?yàn)镈FS也只是取其中一個(gè)周期來運(yùn)算，所以FFT在一定條件下也可以用分析周期信號序列。如果正弦信號sin(２ｆｎ)，ｆ＝０.1，用16點(diǎn)的FFT來做DFS運(yùn)算，得到的頻譜是信號本身的真實(shí)譜嗎？答：當(dāng)N與進(jìn)行FFT變換的點(diǎn)數(shù)K一樣的時(shí)候，可以認(rèn)為DFS與FFT的變換時(shí)相等的，這時(shí)我們可以用DFS來分析FFT。但是在N與K不相等的時(shí)候，DFS與F