數(shù)字信號課程設(shè)計

1、一、利用DFT分析離散序列頻譜應(yīng)用傅里葉變換DFT，分析各種離散序列x(n)的頻譜。離散周期信號的分析 離散周期信號可以展開成離散傅里葉級數(shù)，其中傅里葉系數(shù)如下所示：式中：N是序列的周期，n為時間離散變量， k為數(shù)字頻率離散變量，是k次諧波的數(shù)字頻率。所以離散周期信號的頻譜X (k)是一個以N為周期的周期性離散頻譜，各譜線之間的間隔為0=2/ N，而且存在著諧波的關(guān)系。 離散周期序列在時域與頻域都是離散的、以N為周期的序列。正確地在一個周期內(nèi)選取N，即可以準(zhǔn)確地求取周期序列的頻譜。離散周期序列頻譜的求解步驟： （1）確定離散周期序列的基本周期N； （2）使用fft命令作N點FFT計算。頻率分辨

2、率F0=2/N（3）令X(k0 )=1/NX（k）。 1利用DFT計算序列的頻譜（1）程序N=100;n=0:N-1;x0=(0.5).n;x1=(n>=0);x=x0.*x1;X=1/N*fft(x,N);omega=2*pi/N*(n-N/2);subplot(2,1,1); stem(omega,abs(fftshift(X); ylabel('Magnitude'); xlabel('Frequency (rad)' );subplot(2,1,2); stem(omega,angle(fftshift(X); ylabel('Phase&

3、#39;); xlabel('Frequency (rad) ' );（2）程序N=256;n=0:N-1; x0=(0.5).n;x1=(n>=0);x=x0.*x1;X=1/N*fft(x,N);omega=2*pi/N*(n-N/2);subplot(2,1,1); stem(omega,abs(fftshift(X); ylabel('Magnitude'); xlabel('Frequency (rad)' );subplot(2,1,2); stem(omega,angle(fftshift(X); ylabel('Ph

4、ase'); xlabel('Frequency (rad) ' );（3）程序N=300;n=0:N-1;x0=(0.5).n;x1=(n>=0);x=x0.*x1;X=1/N*fft(x,N);omega=2*pi/N*(n-N/2);subplot(2,1,1); stem(omega,abs(fftshift(X); ylabel('Magnitude'); xlabel('Frequency (rad)' );subplot(2,1,2); stem(omega,angle(fftshift(X); ylabel('

5、;Phase'); xlabel('Frequency (rad) ' );二、連續(xù)系統(tǒng)的分析零輸入、零狀態(tài)法求解系統(tǒng)的響應(yīng)MATLAB的控制工具箱（control toolbox）里包含了許多可用于分析線性時不變（LTI）系統(tǒng)的函數(shù)，使用命令help control可以查看控制工具箱里的這些函數(shù)。 對于線性時不變的連續(xù)系統(tǒng)，在時域中其數(shù)學(xué)模型用常系數(shù)線性微分方程來描述。如下：在拉普拉斯變換中，可用系統(tǒng)函數(shù)H(s)來描述，且根據(jù)其零，極點的分布情況，可以決定系統(tǒng)的結(jié)構(gòu)及系統(tǒng)的穩(wěn)定性。 另外還可以用狀態(tài)變量分析法進(jìn)行系統(tǒng)描述。此設(shè)計完成用atlab實現(xiàn)幾種分析方法的相互轉(zhuǎn)

6、換。在Matlab中，描述系統(tǒng)的傳遞函數(shù)（系統(tǒng)函數(shù)）型tf(transfer function),零極點型zp(zero pole)以及狀態(tài)變量型ss(state space)三種方式可以方便的轉(zhuǎn)換。Matlab相應(yīng)的語句為：tf2zp系統(tǒng)函數(shù)型轉(zhuǎn)換到零極點型tf2ss系統(tǒng)函數(shù)型轉(zhuǎn)換到狀態(tài)變量型zp2tf零極點型轉(zhuǎn)換到系統(tǒng)函數(shù)型zp2ss零極點型轉(zhuǎn)換到狀態(tài)變量型ss2tf狀態(tài)變量型轉(zhuǎn)換到系統(tǒng)函數(shù)型ss2zp狀態(tài)變量型轉(zhuǎn)換到零極點型1:已知一個因果LTI離散系統(tǒng)的系統(tǒng)函數(shù)為： （1）階躍響應(yīng)M=30; N=32;omega=0:(N-1)*2*pi/N;H=(0.0675+0.1349*exp(

7、-j.*omega)+0.0675*exp(-2*j.*omega)./(1-1.143*exp(-j.*omega)+0.4128*exp(-2*j.*omega);x=ones(1,M);X=fft(x,N);Y=X.*H ;y=ifft(Y,N);figure(1);subplot(2,1,1); stem(0:M-1,x); xlabel('Time index n'); ylabel('xn');title ('Input signal and output signal');subplot(2,1,2); stem(0:M-1,rea

8、l(y(1:M); xlabel('Time index n'); ylabel('yn');（2）單位脈沖響應(yīng)h(n)M=30; N=32;omega=0:(N-1)*2*pi/N;H=(0.0675+0.1349*exp(-j.*omega)+0.0675*exp(-2*j.*omega)./(1-1.143*exp(-j.*omega)+0.4128*exp(-2*j.*omega);x=1,zeros(1,M-1);X=fft(x,N);Y=X.*H ;y=ifft(Y,N);figure(1);subplot(2,1,1); stem(0:M-1,x)

9、; xlabel('Time index n'); ylabel('xn');title ('Input signal and output signal');subplot(2,1,2); stem(0:M-1,real(y(1:M); xlabel('Time index n'); ylabel('yn');三、 IIR數(shù)字低通濾波器的設(shè)計（雙線性變換法）從模擬濾波器設(shè)計IIR數(shù)字濾波器就是要由模擬濾波器的系統(tǒng)函數(shù)進(jìn)一步求得。歸根結(jié)底是一個由s平面到z平面的變換。這個變換應(yīng)遵循兩個基本的目標(biāo)。（1）H(z)的頻

10、響必須要模仿的頻響，也即s平面的虛軸j應(yīng)該映射到z平面的單位圓上。（2）的因果穩(wěn)定性，通過映射后仍應(yīng)在所得到的H(z)中保持，也即s平面的左半平面應(yīng)該映射到z平面單位圓以內(nèi)。從模擬濾波器映射(變換)成數(shù)字濾波器常用的是脈沖響應(yīng)不變變換法和雙線性變換法兩種。其中脈沖響應(yīng)不變法設(shè)計數(shù)字濾波器的優(yōu)點是器頻率坐標(biāo)的轉(zhuǎn)換是線性的，在不存在頻率混疊的情況下，能夠完全逼近模擬濾波器的頻率特性。但是當(dāng)模擬濾波器的傳輸函數(shù)超出時，映射后必然存在頻譜混疊的現(xiàn)象，這是脈沖響應(yīng)不變法不足的地方，所以一般只用來設(shè)計低通、帶通濾波器。對于高通、帶通濾波器，由于它們在高頻部分不衰減，當(dāng)一定要追求頻率響應(yīng)線性而采用脈沖響應(yīng)不

11、變法時，必須先對模擬高通、帶阻濾波器加一保護(hù)濾波器，從而加大了設(shè)計成不，一般不采用這種方法。脈沖響應(yīng)不變法設(shè)計濾波器時，從模擬預(yù)到數(shù)字域的轉(zhuǎn)換是“多對一”的關(guān)系，而雙線性變換法是“一對一”的關(guān)系，從而避免了頻率混疊現(xiàn)象，因此本課程設(shè)計中采用雙線性變換法設(shè)計濾波器，下面主要介紹雙線性變換法的原理。1. 變換原理雙線性變換法是使數(shù)字濾波器的頻率響應(yīng)與模擬濾波器的頻率響應(yīng)相似的一種變換方法。為了克服脈沖響應(yīng)不變法的多值映射這一缺點，我們首先把整個s平面壓縮變換到某一中介的平面的一橫帶里（寬度為，即從到），然后再通過上面討論過的標(biāo)準(zhǔn)變換關(guān)系將此橫帶變換到整個z平面上去，這樣就使s平面與z平面是一一對應(yīng)

12、的關(guān)系，消除了多值變換性，也就消除了頻譜混疊現(xiàn)象，基本原理如圖2-3所示。圖2-3 雙線性變換法的映射關(guān)系將s平面整個平面壓縮到平面的到，可采用以下的變換關(guān)系 （2-3）其中C為常數(shù)；這樣經(jīng)（2-3）式變?yōu)椋優(yōu)?，可將上式寫成令，則可得 （2-4）再將平面通過以下標(biāo)準(zhǔn)變化關(guān)系映射到z平面 （2-5）這樣（2-4）式可表示為 （2-6） （2-7）2變換常數(shù)C的選擇為了使模擬濾波器與數(shù)字濾波器在低頻處有較確切的對應(yīng)關(guān)系，即在低頻處有，當(dāng)較小時有由（2-3）式可知因而得 （2-8）則（2-6）和（2-7）式可重新寫成： （2-9）即 （2-10）3.逼近情況雙線性變換具備模擬域到數(shù)字域映射變換的總

13、要求，現(xiàn)分析如下：（1）將代入到(2-9)式則得或 （2-11）由上式可見，當(dāng)時，；當(dāng)時，；當(dāng)時，。這就是說雙線性變換把s左半平面映射在單位圓的內(nèi)部；把s平面的整個軸映射成單位圓，把s右半平面映射在單位圓的外部。(2)令，則由(2-9)式得所以由此得出模擬濾波器的頻率和數(shù)字濾波器頻率的關(guān)系式為 (2-12)這一公式的關(guān)系如圖2-4所示?？梢钥闯?，當(dāng)時，當(dāng)時，,當(dāng)時，。這就是說：s平面的原點映射為z平面(1，0)點，而s平面的正虛軸和負(fù)虛軸分別映射成z平面單位圓的上半圓和下半圓。圖2-4 雙線性變換的頻率間非線性關(guān)系由上所述，可得如下結(jié)論：a模擬濾波器中最大和最小值將保留在數(shù)字濾波器中，因此模擬

14、濾波器的通帶或阻帶變換成數(shù)字濾波器的通帶或阻帶。b如果模擬濾波器是穩(wěn)定的，則通過雙線性變換后所得的數(shù)字濾波器也一定是穩(wěn)定的。c由于s平面的整個虛軸映射為z平面上的單位圓，因此雙線性變換法確實消除了脈沖響應(yīng)不變變換法所存在的混疊誤差，所以逼近是良好的。但由(2-13)式可見，在頻率與間存在嚴(yán)重的非線性。4. 模擬濾波器的數(shù)字化由于雙線性變換法中，s與z之間有簡單的代數(shù)關(guān)系，故可由模擬系統(tǒng)函數(shù)通過代數(shù)置換直接得到數(shù)字濾波器的系統(tǒng)函數(shù)。即 (2-14)可見數(shù)字濾波器的極點數(shù)等于模擬濾波器的極點數(shù)。頻率響應(yīng)也可用直接置換得到 (2-15)這一公式可用于將濾波器的數(shù)字域指標(biāo)，轉(zhuǎn)換為模擬域指標(biāo)。再者，可在

15、未進(jìn)行雙線性變換前把原模擬系統(tǒng)函數(shù)分解成并聯(lián)或級聯(lián)子系統(tǒng)函數(shù)，然后再對每個子系統(tǒng)函數(shù)分別加以雙線性變換。就是說，所有的分解，都可以就模擬濾波器系統(tǒng)函數(shù)來進(jìn)行，因為模擬濾波器已有大量圖表可供利用，且分解模擬系統(tǒng)函數(shù)比較容易。3.2 IIR濾波器的設(shè)計由前面IIR濾波器的設(shè)計原理確定運用雙線性變換法設(shè)計IIR濾波器的流程如圖3-1所示：數(shù)字濾波器指標(biāo)模擬濾波器指標(biāo)設(shè)計模擬濾波器對應(yīng)的數(shù)字濾波器器雙線性變換開始結(jié)束圖3-1 IIR濾波器設(shè)計流程圖3.3 IIR低通濾波器的設(shè)計1.IIR低通濾波器設(shè)計指標(biāo)：，阻帶最小衰減，通帶最大衰減。2.IIR低通濾波器的設(shè)計在Matlab中設(shè)計巴特沃斯型、切比雪夫

16、I型和橢圓型IIR低通濾波器的程序分別如下：%巴特沃斯濾波器設(shè)計fp=1000;fs=1200;Rp=1;Rs=10;Fs=22050;wp=2*pi*fp/Fs;%確定數(shù)字濾波器指標(biāo)ws=2*pi*fs/Fs;Wp=2*tan(wp/2);%轉(zhuǎn)換為模擬濾波器指標(biāo)Ws=2*tan(ws/2);N,Wn=buttord(Wp,Ws,Rp,Rs,'s');%選擇巴特沃斯濾波器的最小階數(shù)N和截止頻率Wnb,a=butter(N,Wn,'s');%設(shè)計巴特沃斯濾波器b1,a1=bilinear(b,a,1);%運用雙線性變換為數(shù)字濾波器%切比雪夫I型濾波器設(shè)計fp=10

17、00;fs=1200;Rp=1;Rs=10;Fs=22050;wp=2*pi*fp/Fs;%確定數(shù)字濾波器指標(biāo)ws=2*pi*fs/Fs;Wp=2*tan(wp/2);%轉(zhuǎn)換為模擬濾波器指標(biāo)Ws=2*tan(ws/2);N,Wn=cheb1ord(Wp,Ws,Rp,Rs,'s');%選擇切比雪夫I型濾波器的最小階數(shù)N和截止頻率Wnb,a=cheby1(N,Rp,Wn,'s');%設(shè)計切比雪夫I型濾波器b1,a1=bilinear(b,a,1);%運用雙線性變換法變?yōu)閿?shù)字濾波器%橢圓型濾波器設(shè)計fp=1000;fs=1200;Rp=1;Rs=10;Fs=22050

18、;wp=2*pi*fp/Fs;%確定數(shù)字濾波器指標(biāo)ws=2*pi*fs/Fs;Wp=2*tan(wp/2);%轉(zhuǎn)換為模擬濾波器指標(biāo)Ws=2*tan(ws/2);N,Wn=ellipord(Wp,Ws,Rp,Rs,'s');%選擇橢圓型濾波器的最小階數(shù)N和截止頻率Wnb,a = ellip(N,Rp,Rs,Wn,'s');%設(shè)計橢圓型濾波器b1,a1=bilinear(b,a,1);%運用雙線性變換法變?yōu)閿?shù)字濾波器得到巴特沃斯濾波器、切比雪夫I型濾波器和橢圓型濾波器的幅頻特性分別如圖3-2、圖3-3和圖3-4所示。圖3-2 IIR低通巴特沃斯型濾波器幅頻特性圖3-

19、3 IIR低通切比雪夫I型濾波器幅頻特性圖3-4 IIR低通橢圓型濾波器幅頻特性總結(jié)經(jīng)過為期一周的課程設(shè)計實踐，我完全自主地完成了既定的任務(wù)，實現(xiàn)了課程設(shè)計的要求，這使我對DSP課程有了更深一步的認(rèn)識和了解，要想學(xué)好它重在實踐，要通過不斷的實際操作才能更好地學(xué)習(xí)和運用知識，并且基本掌握了運用MATLAB軟件進(jìn)行濾波器的設(shè)計和GUI界面的設(shè)計。通過課程設(shè)計這一實踐環(huán)節(jié)，加深了我對本學(xué)期所學(xué)知識的理解，也提高了我運用所學(xué)知識分析問題和解決實際問題的能力。在此過程中我也發(fā)現(xiàn)自己好多不足之處，首先是自己對所學(xué)的基礎(chǔ)知識掌握不深，還有對理論知識的運用還有欠缺，通過這次課程設(shè)計實踐，使我在這幾個方面都有所提高。由于時間原因，這次課程設(shè)計最大的遺憾是沒能實現(xiàn)對語音信號進(jìn)行實時錄入，然后進(jìn)行采樣分析，接下來我會努力嘗試實現(xiàn)這項功能。在這個過程中，我也曾經(jīng)因為實踐經(jīng)驗的缺乏失落過，也曾經(jīng)因設(shè)計成功而熱情高漲。生活就是這樣，汗水預(yù)示著結(jié)果也見證著收