DSP課程設(shè)計(jì)基于MATLAB的FFT算法實(shí)現(xiàn)

1、1 引言 12基于MATLAB勺FFT算法實(shí)現(xiàn) 21.1 系統(tǒng)總體流程圖 21.2 FFT 運(yùn)算規(guī)律及編程思想 31.2.1 語(yǔ)音信號(hào)的采集 31.2.2 DIT-FFT 算法的基本原理 31.2.3 DIT-FFT 算法的運(yùn)算規(guī)律及編程思想 53 Matlab 程序?qū)崿F(xiàn) 104 系統(tǒng)人機(jī)對(duì)話界面 134.1 GUI 簡(jiǎn)介 134.2 界面設(shè)計(jì) 134.3 運(yùn)行調(diào)試 145 心得體會(huì) 16參考文獻(xiàn) 17附錄I 18附錄H 211 引言MATLA是矩陣實(shí)驗(yàn)室(Matrix Laboratory )的簡(jiǎn)稱，是美國(guó) MathWorks公 司出品的商數(shù)學(xué)軟件，用于算法開發(fā)、 數(shù)據(jù)可視化、 數(shù)據(jù)分析以及

2、數(shù)值計(jì)算的高 級(jí)技術(shù)計(jì)算語(yǔ)言和交互式環(huán)境，主要包括MATLAB口 Simulink兩大部分。MATLAB勺應(yīng)用范圍非常廣，包括信號(hào)和圖像處理、通訊、控制系統(tǒng)設(shè)計(jì)、 測(cè)試和測(cè)量、 財(cái)務(wù)建模和分析以及計(jì)算生物學(xué)等眾多應(yīng)用領(lǐng)域。 附加的工具箱 (單 獨(dú)提供的專用MATLAB函數(shù)集)擴(kuò)展了 MATLAB環(huán)境，以解決這些應(yīng)用領(lǐng)域內(nèi)特 定類型的問題。它以矩陣運(yùn)算為基礎(chǔ)，把計(jì)算、可視化、程序設(shè)計(jì)融合在一個(gè)簡(jiǎn) 單易用的交互式工作環(huán)境中， 是一款數(shù)據(jù)分析和處理功能都非常強(qiáng)大的工程適用 軟件。 它可以將聲音文件變換為離散的數(shù)據(jù)文件， 然后利用其強(qiáng)大的矩陣運(yùn)算能 力處理數(shù)據(jù)，如數(shù)據(jù)濾波、傅立葉變換、時(shí)域和頻域分析、

3、聲音回放以及各種圖 的呈現(xiàn)等， 它的信號(hào)處理與分析工具箱位語(yǔ)音信號(hào)分析提供了十分豐富的功能函 數(shù)， 利用這些功能函數(shù)可以快捷而又方便的完成語(yǔ)音信號(hào)的處理和分析以及信號(hào) 的可視化。數(shù)字信號(hào)處理是MATLABt要應(yīng)用的領(lǐng)域之一。對(duì)于有限長(zhǎng)序列x(n),若要求其N點(diǎn)的傅里葉變換(DFT需要經(jīng)過M次復(fù) 數(shù)乘法運(yùn)算和N*(N-1)次復(fù)數(shù)加法運(yùn)算。隨著N的增加，運(yùn)算量將急劇增加，而 在實(shí)際問題中，N往往是較大的，如當(dāng)N=1024時(shí)，完成復(fù)數(shù)乘法和復(fù)數(shù)加法的 次數(shù)分別為百萬(wàn)以上，無論是用通用計(jì)算機(jī)還是用 DSPE片，都需要消耗大量的 時(shí)間和機(jī)器內(nèi)存，不能滿足實(shí)時(shí)的要求。因此， DFT的這種運(yùn)算只能進(jìn)行理論上

4、 的計(jì)算，不適合對(duì)實(shí)時(shí)處理要求高的場(chǎng)合。因此，研究作為 DSP 的快速算法的 FFT是相當(dāng)必要的，快速傅里葉變換(FFT)是為提高DFT運(yùn)算速度而采用的一 種算法， 快速算法的種類很多，而且目前仍在改進(jìn)和提高， 它是根據(jù)離散傅里葉 變換的奇、偶、虛、實(shí)等特性，對(duì)離散傅立葉變換的算法進(jìn)行改進(jìn)獲得的?；?本學(xué)期所學(xué)的 DIT-FFT 的運(yùn)算規(guī)律和編程思想以及Matlab 的學(xué)習(xí)和使用，本課設(shè)要求在Matlab 環(huán)境下編寫基2 DIT-FFT 算法實(shí)現(xiàn)對(duì)離散信號(hào)的快速傅里葉變換，再與Matlab 軟件自帶的FFT 函數(shù)實(shí)現(xiàn)對(duì)離散信號(hào)的傅里葉變換進(jìn)行比較，如果得到的頻譜相同，那么我們編寫的程序就是正

5、確的。其中離散信號(hào)是通過PC自帶的錄音機(jī)錄制一段 wav語(yǔ)音信號(hào)，用Matlab采樣得到離散序列x1。如果 有能力可以選做系統(tǒng)人機(jī)對(duì)話界面。 用GUI界面完成人機(jī)交互方便使用的。 本課 程設(shè)計(jì)主要是對(duì)數(shù)字信號(hào)的分析。2基于MATLAB勺FFT算法實(shí)現(xiàn)2.1 系統(tǒng)總體流程圖本設(shè)計(jì)要求錄制一段個(gè)人自己的語(yǔ)音信號(hào)，并對(duì)錄制的信號(hào)進(jìn)行采樣；畫出采樣后語(yǔ)音信號(hào)的時(shí)域波形和頻譜圖；在Matlab環(huán)境下編寫基2 DIT-FFT算法;利用自己編寫的算法對(duì)已采集的語(yǔ)音信號(hào)進(jìn)行頻譜分析，并畫出語(yǔ)音信號(hào)的時(shí)域與頻譜圖，并與Matlab數(shù)字信號(hào)處理工具箱中的fft函數(shù)進(jìn)行對(duì)比研究，驗(yàn)證 自編算法的正確性。所以得到系統(tǒng)

6、總體流程圖如圖1所示。語(yǔ)音信號(hào)采集編寫fft程序，畫出信號(hào)頻譜圖完成信號(hào)時(shí)域圖兀 成 信 號(hào) 頻 率 響 應(yīng)實(shí) 現(xiàn) 輸 入 信 號(hào) 的 倒 序?qū)崿F(xiàn)一級(jí)中不同種蝶形算運(yùn)實(shí) 現(xiàn)級(jí) 中 相 同 種 蝶 形 運(yùn)Matlab自帶的FFT比較圖1系統(tǒng)總體流程圖2.2 FFT運(yùn)算規(guī)律及編程思想2.2.1 語(yǔ)音信號(hào)的采集利用PM自帶的錄音機(jī)，錄制一段語(yǔ)音信號(hào)，保存格式為 wave的文件，并 將其保存在電腦中。在 MATLAB 中， fn=input(' Enter WAV filename:','s'); x,fs,nb=wavread(fn,n1 n2);用于讀取語(yǔ)音，采樣值

7、放在向量 x中，fs表示采樣頻率(Hz) , nb表示采樣位數(shù)。n1 n2表示讀取從n1點(diǎn)到n2點(diǎn)的值(若 只有一個(gè)n的點(diǎn)則表示讀取前n點(diǎn)的采樣值)。sound(x,fs,nb); 用于對(duì)聲音的 回放。向量x則就代表了一個(gè)信號(hào)(也即一個(gè)復(fù)雜的“函數(shù)表達(dá)式”)也就是說可以像處理一個(gè)信號(hào)表達(dá)式一樣處理這個(gè)聲音信號(hào)。采集到語(yǔ)音信號(hào)之后，需要對(duì)語(yǔ)音信號(hào)進(jìn)行分析，如語(yǔ)音信號(hào)的時(shí)域分析、 頻譜分析、譜圖分析。2.2.2 DIT-FFT算法的基本原理快速傅里葉變換(FFT)是為提高DFT運(yùn)算速度而采用的一種算法。 對(duì)一個(gè)有限長(zhǎng)度序列x(n)的N點(diǎn)的DFT為：N 1X k xn 01 Nx n 一WNnk ,

8、 k 0, 1, 2,., N 1k Wn nk,n 0, 1, 2,., N 1所以，要求N點(diǎn)的DFT,需要2次的復(fù)數(shù)乘法運(yùn)算，N*(N-1)次復(fù)數(shù)乘法運(yùn)算算。隨著N的增加，運(yùn)算量將急劇增加，而在實(shí)際問題中，N往往是較大的，如當(dāng)N=1024時(shí)，完成復(fù)數(shù)乘法和復(fù)數(shù)加法的次數(shù)分別為百萬(wàn)以上，無論是用通用 計(jì)算機(jī)還是用DSP芯片，都需要消耗大量的時(shí)間，不能滿足實(shí)時(shí)的要求不適合于對(duì)實(shí)時(shí)處理要求高的場(chǎng)合。為了能實(shí)時(shí)處理DFT,要想減少DFT的運(yùn)算量可以有兩個(gè)途徑：第一是降N, N的值減小了，運(yùn)算量就減少了；第二是利用旋轉(zhuǎn) 因子的周期性，對(duì)稱性和可約性。利用這兩個(gè)途徑實(shí)現(xiàn)DFT的快速傅里葉變換(FFT)

9、, FFT算法基本上可分為按時(shí)間抽取的 FFT算法(DIT-FFT)和按頻率抽 設(shè)序列x(n)的長(zhǎng)度為N,且N滿足N=2,M為正整數(shù)。若N不能滿足上述關(guān)系， 可以將序列x(n)補(bǔ)零實(shí)現(xiàn)。按時(shí)間抽取基 2-FFT算法的基本思路是將 N點(diǎn)序列 按時(shí)間下標(biāo)的奇偶分為兩個(gè) N/2點(diǎn)序列，計(jì)算這兩個(gè)N/2點(diǎn)序列的N/2點(diǎn)DFI 計(jì)算量可減小約一半；每一個(gè)N/2點(diǎn)序列按照同樣的劃分原則，可以劃分為兩個(gè) N/4點(diǎn)序列，最后，將原序列劃分為多個(gè) 2點(diǎn)序列，將計(jì)算量大大降低。取的FFT算法(DIF-FFT)旋轉(zhuǎn)因子的性質(zhì)：(1)周期性(2)共腕對(duì)稱性(3)可約性本次課設(shè)要求用用基WW;n WWN啊*WNWN)n

10、)* W；n wN k)n* W；(n)*Wkn WCWkn WNk肝2的按時(shí)間抽取的FFT算法(DIT-FFT)實(shí)現(xiàn)FFT功能,按時(shí)間下標(biāo)的奇偶將N點(diǎn)x(n)分別抽取組成兩個(gè)N/2點(diǎn)序列，分別記為x1(n) 和x2(n),將x(n)的DFTW化為x1(n)和x2(n)的DFT的計(jì)算。r 0,1,x(2r) xi(r)x(2r 1)x2(r)n 0,2,4kx nx n WNnkn 1,3,5N12x 2rW2rkr 0,1N12x 2r 1 Wn2r 1 kr 0,12rkxr Wnr 0,1N 122r 1 kx2r Wnr 0,1利用旋轉(zhuǎn)因子的可約性，即:W；*X kj2_2rk e N

11、N 12 .rkx r Wn r 02.2jNrkTN2 k WnrWn：x20r WN"X 1(k) W；X2(k), 0 k N 1用蝶形運(yùn)算可表示為如圖2所示:(k) 人 一_ Y X (#)=國(guó) + 際X" (k)M+幻=X") -際與圖2 DIT-FFT蝶形運(yùn)算流圖符號(hào)以此類推，還可以把x1(n)和x2(n)按n值得奇偶分為兩個(gè)序列，這樣就達(dá)到了降N得目的，從而減少了運(yùn)算量FFT對(duì)DFT的數(shù)學(xué)運(yùn)算量改進(jìn):直接采用DFT進(jìn)行計(jì)算，運(yùn)算量為 M次復(fù)數(shù)乘法和N*(N-1)次復(fù)數(shù)乘法。當(dāng)采用M次FFT時(shí)，由N=皆求得M=logN,運(yùn)算流圖有M級(jí)蝶形，每一級(jí)都 由

12、N/2個(gè)蝶形運(yùn)算構(gòu)成，這樣每一級(jí)蝶形運(yùn)算都需要 N/2次復(fù)數(shù)乘法和N次復(fù)數(shù) 加法。M級(jí)運(yùn)算共需要復(fù)數(shù)乘法次數(shù)為 C=N/2*M,復(fù)數(shù)加法次數(shù)為C=N*M當(dāng)N值較大時(shí)，F(xiàn)FT減少運(yùn)算量的特點(diǎn)表現(xiàn)的越明顯。2.2.3 DIT-FFT 算法的運(yùn)算規(guī)律及編程思想為了編寫 DIT-FFT 算法的運(yùn)算程序， 首先要分析其運(yùn)算規(guī)律， 總結(jié)編程思想 并繪出程序框圖。1. 原位計(jì)算對(duì)N 2M點(diǎn)的FFT共進(jìn)行M級(jí)運(yùn)算，每級(jí)由N/2個(gè)蝶形運(yùn)算組成。在同一級(jí) 中，每個(gè)蝶的輸入數(shù)據(jù)只對(duì)本蝶有用，且輸出節(jié)點(diǎn)與輸入節(jié)點(diǎn)在同一水平線上，這就意味著每算完一個(gè)蝶后， 所得數(shù)據(jù)可立即存入原輸入數(shù)據(jù)所占用的數(shù)組元素 ( 存儲(chǔ)單元 )

13、 ，這種原位( 址 ) 計(jì)算的方法可節(jié)省大量?jī)?nèi)存。2. 蝶形運(yùn)算實(shí)現(xiàn)FFT運(yùn)算的核心是蝶形運(yùn)算，找出蝶形運(yùn)算的規(guī)律是編程的基礎(chǔ)。蝶形 運(yùn)算是分級(jí)進(jìn)行的； 每級(jí)的蝶形運(yùn)算可以按旋轉(zhuǎn)因子的指數(shù)大小排序進(jìn)行； 如果 指數(shù)大小一樣則可從上往下依次蝶算。對(duì) N 2M點(diǎn)的FFT共有M級(jí)運(yùn)算，用L表 示從左到右的運(yùn)算級(jí)數(shù)(L=1, 2,M)。第L級(jí)共有B 2L1個(gè)不同指數(shù)的旋轉(zhuǎn) 因子，用R表示這些不同指數(shù)旋轉(zhuǎn)因子從上到下的順序(R=0, 1,，B-1) 0第R個(gè)旋轉(zhuǎn)因子的指數(shù)P 2MLR,旋轉(zhuǎn)因子指數(shù)為P的第一個(gè)蝶的第一節(jié)點(diǎn)標(biāo)號(hào)k 從R開始，由于本級(jí)中旋轉(zhuǎn)因子指數(shù)相同的蝶共有2ML個(gè)，且這些蝶的相鄰間距為

14、2L ， 故旋 轉(zhuǎn) 因子指數(shù)為 P 的最后 一個(gè)蝶 的 第 一節(jié) 點(diǎn)標(biāo)號(hào) k 為 ：MLLL(21) 2 R N 2 R ，本級(jí)中各蝶的第二個(gè)節(jié)點(diǎn)與第一個(gè)節(jié)點(diǎn)都相距B 點(diǎn)。應(yīng)用原位計(jì)算，蝶形運(yùn)算可表示成如下形式：AL(J)= AL-1(J)+ AL-1 (J+B)* WNPAL(J+B)= AL-1(J)- AL-1(J+B)* WNP總結(jié)上述運(yùn)算規(guī)律， 可采用如下運(yùn)算方法進(jìn)行DIT-FFT 運(yùn)算。 首先讀入數(shù)據(jù)，根據(jù)數(shù)據(jù)長(zhǎng)度確定運(yùn)算級(jí)數(shù) M,運(yùn)算總點(diǎn)數(shù)N 2M ,不足補(bǔ)0處理。然后對(duì)讀入 數(shù)據(jù)進(jìn)行數(shù)據(jù)倒序操作。數(shù)據(jù)倒序后從第1級(jí)開始逐級(jí)進(jìn)行，共進(jìn)行 M級(jí)運(yùn)算。在進(jìn)行第 L 級(jí)運(yùn)算時(shí)，先算出該

15、級(jí)不同旋轉(zhuǎn)因子的個(gè)數(shù)B 2L1 (也是該級(jí)中各個(gè)蝶形運(yùn)算兩輸入數(shù)據(jù)的間距)，再?gòu)腞=0開始按序計(jì)算，直到R=B-1結(jié)束。每個(gè) R對(duì)應(yīng)的旋轉(zhuǎn)因子指數(shù)P 2MLR，旋轉(zhuǎn)因子指數(shù)相同的蝶從上往下依次逐個(gè)運(yùn)算， 各個(gè)蝶的第一節(jié)點(diǎn)標(biāo)號(hào)k都是從R開始，以2L為步長(zhǎng)，到N 2L R(可簡(jiǎn)取極值N-2）結(jié)束?？紤]到蝶形運(yùn)算有兩個(gè)輸出，且都要用到本級(jí)的兩個(gè)輸入數(shù)據(jù)，故第 一個(gè)輸出計(jì)算完畢后，輸出數(shù)據(jù)不能立即存入輸入地址，要等到第二個(gè)輸出計(jì)算 調(diào)用輸入數(shù)據(jù)完畢后才能覆蓋。這樣數(shù)據(jù)倒序后的運(yùn)算可用三重循環(huán)程序?qū)崿F(xiàn) 整個(gè)蝶形運(yùn)算流程圖如圖3所示。3. 序列倒序?yàn)榱吮ＷC運(yùn)算輸出的X(k)按順序排列，要求序列x(n)倒序

16、輸入，即在運(yùn)算 前要先對(duì)輸入的序列進(jìn)行位序顛倒。如果總點(diǎn)數(shù)為n 2M的x(n)的順序數(shù)是用M位二進(jìn)制數(shù)表示，則倒序數(shù)只需將順序數(shù)的二進(jìn)制位倒置即可，按照這一規(guī)律用硬件電路和匯編語(yǔ)言很容易產(chǎn)生倒序數(shù)。但用MATLA舒高級(jí)語(yǔ)言實(shí)現(xiàn)倒序時(shí)，直接倒置二進(jìn)制數(shù)位的方法不可取，還須找出產(chǎn)生倒序的十進(jìn)制規(guī)律。 將十進(jìn)制順序數(shù)用I表示，與之對(duì)應(yīng)的二進(jìn)制數(shù)用 舊表示。十進(jìn)制倒序數(shù)用J表示，與 之對(duì)應(yīng)的二進(jìn)制數(shù)用JB表示。JB是舊的位倒置結(jié)果，十進(jìn)制順序數(shù)I增加1, 相當(dāng)于 舊最低位加1且逢2向高位進(jìn)1,即相當(dāng)于JB最高位加1且逢2向低位 進(jìn)10JB的變化規(guī)律反映到J的變化分二種情況：如果JB的最高位是0(J n

17、/2), 則直接由加1(J J n/2)得到下一個(gè)倒序值；如果JB的最高位是1(J n/2),則 要先將最高位變0(J J n/2),再在次高位加1(J J n/4)。但次高位加1時(shí)， 同樣要判斷0、1值，如果是0 (J N/4),則直接加1(J J N/4),否則要先將 次高位變0(J J n/4),再判斷下一位。依此類推，直到完成最高位加1,逢2向右進(jìn)位的運(yùn)算。利用這一算法可按順序數(shù)I的遞增順序，依次求得與之對(duì)應(yīng)的 倒序數(shù)Jo為了節(jié)省內(nèi)存，數(shù)據(jù)倒序可原址進(jìn)行，當(dāng)I = J時(shí)不需要交換，當(dāng)I* J時(shí)需要交換數(shù)據(jù)。另外，為了避免再次調(diào)換前面已經(jīng)調(diào)換過的一對(duì)數(shù)據(jù)，只對(duì) I <J的情況進(jìn)行數(shù)

18、據(jù)交換即可實(shí)現(xiàn)數(shù)據(jù)倒序操作。圖 3中數(shù)據(jù)倒序的程序流程圖 如圖4所示。例如，N=8時(shí)，序列倒序結(jié)果如表1所示。表1碼位倒序(N=8)00 0 00 0 0010 0 11 0 0420 1 00 1 0230 1 11 1 0641 0 00 0 1151 0 11 0 1561 1 00 1 1371 1 11 1 17圖4倒序流程圖3 Matlab程序?qū)崿F(xiàn)用ringin.wav作為例子，運(yùn)行調(diào)試程序。1.程序運(yùn)行開始時(shí)，要求輸入采樣頻率fs (fs=1024),所需要變換的起點(diǎn)N1 (N1=100O和終點(diǎn)N2 (N2=509§以及要采樣白語(yǔ)音文件(ringin.wav );其 中

19、N1和N2的差值必須在語(yǔ)音信號(hào)長(zhǎng)度范圍內(nèi)， 但不能太小，否則聽不到較為清 晰的語(yǔ)音。具輸入窗口顯示如圖 5所示。圖5程序運(yùn)行開始的輸入窗口將fs=1024, N1=100Q N2=5095和語(yǔ)音信號(hào)ringin.wav 輸入后得到采樣后 的語(yǔ)音信號(hào)x1的時(shí)域波形，如圖6所示，其頻率響應(yīng)特性如圖7所示，Matlab 自帶的FFT函數(shù)實(shí)現(xiàn)的x1的頻譜與編寫的FFT程序?qū)崿F(xiàn)的x1的頻譜的比較如圖 8所示。語(yǔ)音信號(hào)時(shí)域波形0.80.60.40.20-0.2-0.4-0.6-0.8 匚1LJ050010001500200025003000350040004500圖6采樣后的語(yǔ)音信號(hào)的時(shí)域波形60x1的頻

20、率響應(yīng)特性50403020100.511.522.5頻率3.5圖7 x1的頻率響應(yīng)特性Matlab自帶的FFT函數(shù)實(shí)現(xiàn)的x1的頻譜200150 -k 100 一-Y|50 -I ,-JhLI11 JH 10仁；. N'一 2、:-中, 一050010001500200025003000350040004500KK圖8 FFT函數(shù)實(shí)現(xiàn)與編寫的 FFT程序?qū)崿F(xiàn)的比較取不同的點(diǎn)數(shù)進(jìn)行FFT變換，經(jīng)觀察，編寫FFT程序得到的語(yǔ)音信號(hào)的頻譜 圖與Matlab中自帶的FFT函數(shù)得到的語(yǔ)音信號(hào)頻譜圖總是基本一致，但是如果 輸入的N1和N2差值加1不是2的整數(shù)次幕就會(huì)有細(xì)微差別，這是因?yàn)榫帉懙目?速傅

21、里葉變換計(jì)算信號(hào)的 N點(diǎn)傅里葉變換要求N為2的整數(shù)次幕，不夠的話信 號(hào)將會(huì)被補(bǔ)零后運(yùn)算，即參與運(yùn)算的信號(hào)已經(jīng)不同（差別很?。?，所以FFT運(yùn)算后的結(jié)果也不盡相同，所示頻譜圖自然有細(xì)微的差別。由用MATLAB自帶FFT函數(shù)實(shí)現(xiàn)的頻譜圖與用 MATLA踹寫的FFT程序?qū)崿F(xiàn) 的頻譜圖相比較，可知，兩個(gè)算法計(jì)算后的結(jié)果幾乎相同，驗(yàn)證了自編算法的正 確性。如果改變輸入的值和語(yǔ)音信號(hào)，那么得到不一樣的時(shí)域波形，頻率響應(yīng)，通 過FFT得到頻譜都會(huì)不同。需要注意的是自由輸入的語(yǔ)音信號(hào)不同，其長(zhǎng)度也不同，所以根據(jù)要輸入的語(yǔ)音信號(hào)輸入 N1, N2的值。4系統(tǒng)人機(jī)對(duì)話界面4.1 GUI簡(jiǎn)介圖形用戶界面(GUI),是

22、一種提供人機(jī)交互的工具和方法。GUI是包含圖形對(duì)象，如窗口、圖標(biāo)、菜單和文本等圖文并茂的用戶界面。4.2 界面設(shè)計(jì)用MATLA圖形用戶界面開發(fā)環(huán)境設(shè)計(jì) GUI點(diǎn)的一般步驟是：1 .進(jìn)行界面設(shè)計(jì)。2 .設(shè)計(jì)控件屬性。3 .進(jìn)行M語(yǔ)言編程。以本設(shè)計(jì)要求為例介紹。第一步，該選擇本圖形用戶界面需要的控件：三個(gè)文本框(Edit Text) 用來輸入N1, N2和fs。三個(gè)靜態(tài)文本(static text)用來提示我們需要人為輸入的內(nèi)容。兩個(gè)推按鈕(Push button),用來運(yùn)行和退出。四個(gè)軸對(duì)象(axes)用來顯示信號(hào)時(shí)域波形，頻譜圖，和兩種不同F(xiàn)FT實(shí)現(xiàn)的頻譜 圖。完成人機(jī)設(shè)計(jì)界面如圖9所示：圖9

23、完成人機(jī)設(shè)計(jì)界面第二步，設(shè)置控件屬性：雙擊組件可以設(shè)置文本框，推按鈕的屬性，如顯示大小，名稱和默認(rèn)值等。第三步，編寫回調(diào)函數(shù)。組件事件的發(fā)生是通過回調(diào)函數(shù)進(jìn)行工作的。控件設(shè)置完成后保存，然后運(yùn)行GUI （操彳乍為ctrl+T ）,就會(huì)進(jìn)入editor窗口，加入各個(gè)控件功能的函數(shù)代碼。 完成后保存即可。第四步，運(yùn)行GUI。運(yùn)行editor窗口的程序后，會(huì)彈出已經(jīng)激活的人機(jī)對(duì)話界面。系統(tǒng)人機(jī)對(duì)話界面如圖10所示。圖10系統(tǒng)人機(jī)對(duì)話界面4.3運(yùn)行調(diào)試運(yùn)行GUI,已經(jīng)彈出圖10的系統(tǒng)人機(jī)對(duì)話界面，根據(jù)提示“Input N1:"，“Input N2:"和“Input fs :&quo

24、t;分別輸入N1, N2和fs ,然后點(diǎn)擊“運(yùn)行”推按鈕就會(huì)得 到本設(shè)計(jì)的要求。運(yùn)行結(jié)果如圖11所示。如果運(yùn)行結(jié)束，可以通過點(diǎn)擊“退出” 推按鈕退出該人機(jī)對(duì)話界面，返回 Matlab。匚慶圖11運(yùn)行結(jié)果5 心得體會(huì)本次實(shí)習(xí)的主要內(nèi)容是通過用 Matlab 實(shí)現(xiàn)FFT 的設(shè)計(jì)，可以實(shí)現(xiàn)對(duì)一段自己錄制的語(yǔ)音信號(hào)進(jìn)行分析，并畫出采樣信號(hào)的時(shí)域與頻域圖。把自己編寫的FFT 算法與 Matlab 自帶 FFT 算法進(jìn)行比較。程序運(yùn)行調(diào)試時(shí)，自己選擇輸入要采樣的語(yǔ)音信號(hào)， 采樣頻率以及要變換的范圍， 可以實(shí)現(xiàn)對(duì)不同信號(hào)的信號(hào)采樣和進(jìn)行不同點(diǎn)的FFT運(yùn)算。在之前數(shù)字信號(hào)處理的學(xué)習(xí)以及完成實(shí)驗(yàn)的過程中，已經(jīng)使

25、用過Matlab ，對(duì)其有了一些基礎(chǔ)的了解和認(rèn)識(shí)， 通過這次的課程設(shè)計(jì)使我進(jìn)一步了解了信號(hào)的產(chǎn)生， 采樣及頻譜分析的方法， 以及其中產(chǎn)生信號(hào)和繪制信號(hào)的基本命令和一些基礎(chǔ)編程語(yǔ)言。 讓我感受到只有在了解課本知識(shí)的前提下， 才能更好的應(yīng)用這個(gè)工具，并且熟練的應(yīng)用 Matlab 也可以很好的加深我對(duì)課程的理解，方便我的思維。這次課程設(shè)計(jì)使我了解了 Matlab 的使用方法，提高了自己的分析和動(dòng)手實(shí)踐能力。同時(shí)我相信，進(jìn)一步加強(qiáng)對(duì)MATLAB勺學(xué)習(xí)與研究對(duì)我今后的學(xué)習(xí)將會(huì)起到很大的幫助。這次的課程設(shè)計(jì)是對(duì)本學(xué)期所學(xué)知識(shí)的一次重要鞏固， 使得在課堂上掌握的知識(shí)得到了真正的運(yùn)用。 在學(xué)習(xí)的過程中和同學(xué)討

26、論， 更明白了理論知識(shí)與實(shí)踐的聯(lián)系。 書到用時(shí)方恨少， 有些知識(shí)學(xué)會(huì)是一回事， 掌握是一回事， 但應(yīng)用起來，確實(shí)不是那么簡(jiǎn)單的，需要很多知識(shí)的融會(huì)貫通。程序運(yùn)行調(diào)試初期， 曾經(jīng)多次出現(xiàn)錯(cuò)誤、不能產(chǎn)生圖形等問題，但在我翻閱資料認(rèn)真改正及老師同學(xué)的幫助下基本功能還是完成了， 經(jīng)過 1 個(gè)星期的上機(jī)實(shí)習(xí)，程序已得到一些完善，能完成基本的要求的功能。最后經(jīng)過努力，又深入學(xué)習(xí)了圖形用戶界面（GUI）,完成了選做要求的人機(jī)對(duì)話界面。學(xué)習(xí)就是一個(gè)了解， 疑惑，進(jìn)而解惑的過程，這次實(shí)習(xí)就是提供了這樣一個(gè)發(fā)現(xiàn)自己知識(shí)漏洞，與同學(xué)老師探討進(jìn)行解惑的的機(jī)會(huì)。通過這次課程設(shè)計(jì)實(shí)習(xí)， 我更深刻的了解了 Matlab 的

27、運(yùn)用， 重新復(fù)習(xí)了 FFT中的重要的序列倒序和蝶形變換的程序， 對(duì)課本上的知識(shí)有了更深的理解， 使我對(duì)數(shù)字信號(hào)處理有了系統(tǒng)的認(rèn)知。在這里特別感謝董老師和李老師， 他們給了我們很大的發(fā)揮空間， 讓我們真正自己動(dòng)手真正掌握了知識(shí)，感謝他們細(xì)心指導(dǎo)。也非常感謝我的同學(xué)， 他們解開了我在實(shí)習(xí)中出現(xiàn)的諸多知識(shí)死角，謝謝大家！參考文獻(xiàn)1范壽康DSP技術(shù)與DSP芯片.北京：電子工業(yè)出版社2 程佩青 . 數(shù)字信號(hào)處理教程. 北京：清華大學(xué)出版社出版， 20013 高西全 丁玉美等 . 數(shù)字信號(hào)處理. 北京：電子工業(yè)出版社， 20094李勇徐震.MATLAB輔助現(xiàn)代工程數(shù)字信號(hào)處理.西安電子科技大學(xué)出版社5 陳

28、杰 .Matlab 寶典 . 電子工業(yè)出版社6蘇金明張蓮花劉波.MATLAB工具箱應(yīng)用，電子工業(yè)出版社附錄Ifs=input(' 輸入采樣頻率fs='); % 語(yǔ)音信號(hào)采樣頻率為 fsN1=input(' 輸入所需變換的起點(diǎn)N1=');N2=input(' 輸入所需變換的終點(diǎn)N2=');fn=input(' Enter WAV filename:','s'); %獲取一個(gè) *.wav 的文件x,fs,nb=wavread(fn,N1 N2); % 讀取語(yǔ)音信號(hào)的數(shù)據(jù)sound(x,fs,nb); % 播放語(yǔ)音信號(hào)%

29、n=N2-N1+1;臉語(yǔ)音信號(hào)文件較大時(shí)用這兩條%x1=reshape(x,1,2*n);% 語(yǔ)句替換 x1=x'x1=x'y1=fft(x1);figure(1)plot(x1) % 做原始語(yǔ)音信號(hào)的時(shí)域圖形title(' 語(yǔ)音信號(hào)時(shí)域波形 ')xlabel('n');ylabel(' 幅值 ');M=nextpow2(x1); % 求x的長(zhǎng)度對(duì)應(yīng)的2的最低幕次mN=2AM;if length(x1)<Nx1=x1,zeros(1,N-length(x1);%若x的長(zhǎng)度不是2的幕，補(bǔ)零到2的整數(shù)冪end%數(shù)據(jù)倒序操作J=0

30、;%給倒序數(shù)賦初值for I=0:N-1;% 按序交換數(shù)據(jù)和算倒序數(shù)if I<J;% 條件判斷及數(shù)據(jù)交換T=x1(I+1);x1(I+1)=x1(J+1);x1(J+1)=T;end% 算下一個(gè)倒序數(shù)K=N/2;while J>=K;J=J-K;K=K/2;endJ=J+K;end %x1;y=x1; %將x倒序排列作為y的初始值WN=exp(-i*2*pi/N);for L=1:MB=2AL/2;% 第L級(jí)中，每個(gè)蝶形的兩個(gè)輸入數(shù)據(jù)相距B個(gè)點(diǎn)，每級(jí)有Bt不同的旋轉(zhuǎn)因子p=J*2A(M-L);WNp=WNAp;for k=J+1:2AL:N %kp=k+B; %t=y(kp)*WN

31、p; % y(kp)=y(k)-t; %for J=0:B-1 % J 代表了不同的旋轉(zhuǎn)因子本次蝶形運(yùn)算的跨越間隔為 2AL蝶形運(yùn)算的兩個(gè)因子對(duì)應(yīng)單元下標(biāo)的關(guān)系蝶形運(yùn)算的乘積項(xiàng)蝶形運(yùn)算， 注意必須先進(jìn)行減法運(yùn)算，然后進(jìn)行加法運(yùn)算，否則要使用中間變量來傳遞y(k)y(k)=y(k)+t; %蝶形運(yùn)算endendend%yfigure(2)x1,w1=freqz(x1,1); % 繪制原始語(yǔ)音信號(hào)的頻率圖plot(w1,abs(x1);title('x1 的頻率響應(yīng)特性')xlabel(' 頻率');ylabel(' 幅度');figure(3)su

32、bplot(2,1,1);plot(abs(y1) %Matlab 自帶的FFT®數(shù)實(shí)現(xiàn)的語(yǔ)音信號(hào)的FF砌譜圖title('Matlab 自帶的FFT®數(shù)實(shí)現(xiàn)的x1的頻譜')xlabel('K');ylabel('Y1(k)');subplot(2,1,2);plot(abs(y); % 編寫的FFTS序?qū)崿F(xiàn)的語(yǔ)音信號(hào)的FFT1譜圖 title(' 編寫的FFTS序?qū)崿F(xiàn)的的x1的頻譜') xlabel('K');附錄nfunction varargout = gaoli208(varargin)

33、% GAOLI208 M-file for gaoli208.fig% GAOLI208, by itself, creates a new GAOLI208or raises the existing% singleton*.% H = GAOLI208 returns the handle to a new GAOLI208 or the handleto% the existing singleton*.% GAOLI208('CALLBACK',hObject,eventData,handles,.) calls the localfunction named CALL

34、BACKin GAOLI208.M with the given inputarguments.%GAOLI208('Property','Value',.)creates a newGAOLI208or raisesthe% existing singleton*. Starting from the left, property value pairs are% applied to the GUI before gaoli208_OpeningFunction gets called.An% unrecognized property name or in

35、valid value makes property application% stop. All inputs are passed to gaoli208_OpeningFcn via varargin.% *See GUI Options on GUIDE's Tools menu. Choose "GUI allows only one% instance to run (singleton)".% See also: GUIDE, GUIDATA, GUIHANDLES % Edit the above text to modify the respons

36、e to help gaoli208% Last Modified by GUIDE v2.5 11-Jan-2010 14:14:19 % Begin initialization code - DO NOT EDITgui_Singleton = 1;gui_State = struct('gui_Name', mfilename, .'gui_Singleton', gui_Singleton, .'gui_OpeningFcn', gaoli208_OpeningFcn, .'gui_OutputFcn', gaoli20

37、8_OutputFcn, .'gui_LayoutFcn', , .'gui_Callback', );if nargin && ischar(varargin1)gui_State.gui_Callback = str2func(varargin1);endif nargoutvarargout1:nargout = gui_mainfcn(gui_State, varargin:);elsegui_mainfcn(gui_State, varargin:);end% End initialization code - DO NOT EDIT%

38、 - Executes just before gaoli208 is made visible.function gaoli208_OpeningFcn(hObject, eventdata, handles, varargin)% This function has no output args, see OutputFcn.% hObject handle to figure% eventdata reserved - to be defined in a future version of MATLAB% handles structure with handles and user

39、data (see GUIDATA)% varargin command line arguments to gaoli208 (see VARARGIN)% Choose default command line output for gaoli208handles.output = hObject;% Update handles structureguidata(hObject, handles);% UIWAIT makes gaoli208 wait for user response (see UIRESUME) % uiwait(handles.figure1);% - Outp

40、uts from this function are returned to the command line. function varargout = gaoli208_OutputFcn(hObject, eventdata, handles) % varargout cell array for returning output args (see VARARGOUT);% hObject handle to figure% eventdata reserved - to be defined in a future version of MATLAB% handles structu

41、re with handles and user data (see GUIDATA)% Get default command line output from handles structure varargout1 = handles.output;function N1_Callback(hObject, eventdata, handles)% hObject handle to N1 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB% handles structure with

42、 handles and user data (see GUIDATA)% Hints: get(hObject,'String') returns contents of N1 as text% str2double(get(hObject,'String') returns contents of N1 as a double% - Executes during object creation, after setting all properties.function N1_CreateFcn(hObject, eventdata, handles)%

43、hObject handle to N1 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB%handles empty - handles not created until after all CreateFcns called% Hint: edit controls usually have a white background on Windows.% See ISPC and COMPUTER.ifispc&& isequal(get(hObject,'Ba

44、ckgroundColor'),get(0,'defaultUicontrolBackgroundColor')set(hObject,'BackgroundColor','white');end function N2_Callback(hObject, eventdata, handles)% hObject handle to N2 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB% handles structure with

45、handles and user data (see GUIDATA)% Hints: get(hObject,'String') returns contents of N2 as text% str2double(get(hObject,'String') returns contents of N2 as a double% - Executes during object creation, after setting all properties.function N2_CreateFcn(hObject, eventdata, handles)% h

46、Object handle to N2 (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB%handles empty - handles not created until after all CreateFcns called% Hint: edit controls usually have a white background on Windows.% See ISPC and COMPUTER.ifispc&& isequal(get(hObject,'Bac

47、kgroundColor'),get(0,'defaultUicontrolBackgroundColor')set(hObject,'BackgroundColor','white');end function fs_Callback(hObject, eventdata, handles)% hObject handle to fs (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB% handles structure with h

48、andles and user data (see GUIDATA)% Hints: get(hObject,'String') returns contents of fs as text% str2double(get(hObject,'String') returns contents of fs as a double% - Executes during object creation, after setting all properties. function fs_CreateFcn(hObject, eventdata, handles)% h

49、Object handle to fs (see GCBO)% eventdata reserved - to be defined in a future version of MATLAB%handles empty - handles not created until after all CreateFcns called% Hint: edit controls usually have a white background on Windows.% See ISPC and COMPUTER.ifispc&& isequal(get(hObject,'Bac

50、kgroundColor'),get(0,'defaultUicontrolBackgroundColor') set(hObject,'BackgroundColor','white');end% - Executes on button press in begin.function begin_Callback(hObject, eventdata, handles)% hObject handle to begin (see GCBO)% eventdata reserved - to be defined in a future

51、 version of MATLAB% handles structure with handles and user data (see GUIDATA)N1=str2double(get(handles.N1,'String');N2=str2double(get(handles.N2,'String');fs=str2double(get(handles.fs,'String');x,fs,bits=wavread('D:Matlabworkringin.wav',N1 N2);sound(x,fs,bits); % 播放語(yǔ)

52、音信號(hào)%n=N2-N1+1;臉語(yǔ)音信號(hào)文件較大時(shí)用這兩條%x1=reshape(x,1,2*n);% 語(yǔ)句替換 x1=x'x1=x'y1=fft(x1);axes(handles.tu1);plot(x1) % 做原始語(yǔ)音信號(hào)的時(shí)域圖形title(' 語(yǔ)音信號(hào)時(shí)域波形 ')xlabel('n');ylabel(' 幅值 ');M=nextpow2(x1); % 求x的長(zhǎng)度對(duì)應(yīng)的2的最低幕次mN=2AM;if length(x1)<Nx1=x1,zeros(1,N-length(x1);%若x的長(zhǎng)度不是2的幕，補(bǔ)零到2的整數(shù)冪

53、end%數(shù)據(jù)倒序操作J=0;%給倒序數(shù)賦初值for I=0:N-1;% 按序交換數(shù)據(jù)和算倒序數(shù)if I<J;% 條件判斷及數(shù)據(jù)交換T=x1(I+1);x1(I+1)=x1(J+1);x1(J+1)=T;end% 算下一個(gè)倒序數(shù)K=N/2;while J>=K;J=J-K;K=K/2;endJ=J+K;end%x1;y=x1; %將x倒序排列作為y的初始值WN=exp(-i*2*pi/N);for L=1:MB=2AL/2;% 第L級(jí)中，每個(gè)蝶形的兩個(gè)輸入數(shù)據(jù)相距B個(gè)點(diǎn)，每級(jí)有Bt不同的旋轉(zhuǎn)因子for J=0:B-1 % J 代表了不同的旋轉(zhuǎn)因子WNp=WNAp;for k=J+1:

54、2AL:N %kp=k+B; %t=y(kp)*WNp; % y(kp)=y(k)-t; %p=J*2A(M-L);本次蝶形運(yùn)算的跨越間隔為 2AL蝶形運(yùn)算的兩個(gè)因子對(duì)應(yīng)單元下標(biāo)的關(guān)系蝶形運(yùn)算的乘積項(xiàng)蝶形運(yùn)算， 注意必須先進(jìn)行減法運(yùn)算，然后進(jìn)行加法運(yùn)算，否則要使用中間變量來傳遞y(k)y(k)=y(k)+t; %蝶形運(yùn)算endendend%yaxes(handles.tu2);x1,w1=freqz(x1,1); % 繪制原始語(yǔ)音信號(hào)的頻率圖plot(w1,abs(x1);title('x1 的頻率響應(yīng)特性')xlabel(' 頻率');ylabel(' 幅度');axes(handles.tu3);plot(abs(y1) %Matlab自帶的FFT®數(shù)實(shí)現(xiàn)的語(yǔ)音信號(hào)的FF砌譜圖title('Matlab自帶的FFT®數(shù)實(shí)現(xiàn)的x1的頻譜')xlabel('K');ylabel('Y1(k)');axes(handles.tu4);plot(abs(y); % 編寫的FFTS序?qū)崿F(xiàn)的語(yǔ)音信號(hào)的FF砌譜