數(shù)字信號(hào)處理實(shí)驗(yàn)報(bào)告

數(shù)字信號(hào)處理實(shí)驗(yàn)報(bào)告目錄實(shí)驗(yàn)1利用DFT分析信號(hào)頻譜一、實(shí)驗(yàn)?zāi)康亩?shí)驗(yàn)原理三、實(shí)驗(yàn)內(nèi)容四、實(shí)驗(yàn)總結(jié)實(shí)驗(yàn)2利用FFT計(jì)算線性卷積一、實(shí)驗(yàn)?zāi)康亩?shí)驗(yàn)原理三、實(shí)驗(yàn)內(nèi)容四、實(shí)驗(yàn)總結(jié)實(shí)驗(yàn)3IIR數(shù)字濾波器設(shè)計(jì)一、實(shí)驗(yàn)?zāi)康亩?shí)驗(yàn)原理三、實(shí)驗(yàn)內(nèi)容四、實(shí)驗(yàn)總結(jié)實(shí)驗(yàn)4FIR數(shù)字濾波器設(shè)計(jì)一、實(shí)驗(yàn)?zāi)康亩?shí)驗(yàn)原理三、實(shí)驗(yàn)內(nèi)容四、實(shí)驗(yàn)總結(jié)實(shí)驗(yàn)1利用DFT分析信號(hào)頻譜一、實(shí)驗(yàn)?zāi)康?.加深對(duì)DFT原理的理解。2.應(yīng)用DFT分析信號(hào)頻譜。3.深刻理解利用DFT分析信號(hào)頻譜的原理，分析現(xiàn)實(shí)過(guò)程現(xiàn)象及解決辦法。二、實(shí)驗(yàn)原理1、DFT和DTFT的關(guān)系有限長(zhǎng)序列的離散時(shí)間傅里葉變換在頻率區(qū)間的N個(gè)等分點(diǎn)上的N個(gè)取樣值可以由下式表示：由上式可知，序列的N點(diǎn)DFT，實(shí)際上就是序列的DTFT在N個(gè)等間隔頻率點(diǎn)上樣本。2、利用DFT求DTFT方法1：由恢復(fù)出的方法如圖2.1所示：圖2.1.由N點(diǎn)DFT恢復(fù)頻譜DTFT的流程由圖2.1所示流程圖可知：由式2-2可以得到其中為內(nèi)插函數(shù)方法2：然而在實(shí)際MATLAB計(jì)算中，上訴插值公式不見(jiàn)得是最好的方法。由于DFT是DTFT的取樣值，其相鄰的兩個(gè)頻率樣本點(diǎn)的間距為,所以如果我們?cè)黾訑?shù)據(jù)的長(zhǎng)度N，使得得到的DFT譜線就更加精細(xì)，其包絡(luò)就越接近DTFT的結(jié)果，這樣可以利用DFT來(lái)近似計(jì)算DTFT。如果沒(méi)有更多的數(shù)據(jù)，可以通過(guò)補(bǔ)零來(lái)增加數(shù)據(jù)長(zhǎng)度。3、利用DFT分析連續(xù)時(shí)間信號(hào)的頻譜采用計(jì)算機(jī)分析連續(xù)時(shí)間信號(hào)的頻譜，第一步就是把連續(xù)時(shí)間信號(hào)離散化，這里需要進(jìn)行連個(gè)操作：一是采樣，二是截?cái)?。?duì)于連續(xù)非周期信號(hào)，按采樣間隔T進(jìn)行采樣，截取長(zhǎng)度為M，那么對(duì)進(jìn)行N點(diǎn)的頻率采樣，得到因此，可以將利用DFT分析連續(xù)非周期信號(hào)頻譜的步驟歸納如下：（1）確定時(shí)域采樣間隔T，得到離散序列；（2）確定截取長(zhǎng)度M，得到M點(diǎn)離散序列，這里的為窗函數(shù)。（3）確定頻域采樣點(diǎn)數(shù)N，要求。（4）利用FFT計(jì)算離散序列的N點(diǎn)DFT，得到。（5）根據(jù)式（2-6）由計(jì)算采樣點(diǎn)的近似值。采用上訴方法計(jì)算的頻譜，需要注意如下三點(diǎn)問(wèn)題：（1）頻譜混疊。如果不滿足采樣定理的條件，頻譜會(huì)很出現(xiàn)混疊誤差。對(duì)于頻譜無(wú)限寬的信號(hào)，應(yīng)考慮覆蓋大部分主要頻率的范圍。（2）柵欄效應(yīng)和頻譜分辨率。使用DFT計(jì)算頻譜，得到的結(jié)果只是N個(gè)頻譜樣本值，樣本值之間的頻譜是未知的，就像通過(guò)一個(gè)柵欄觀察頻譜，稱(chēng)為“柵欄效應(yīng)”。頻譜分辨率與記錄長(zhǎng)度成正比，提高頻譜分辨率，就要增加記錄時(shí)間。（3）頻譜泄露。對(duì)于信號(hào)截?cái)鄷?huì)把窗函數(shù)的頻譜會(huì)引入到信號(hào)頻譜中，造成頻譜泄露。解決這問(wèn)題的主要辦法是采用旁瓣小的窗函數(shù)，頻譜泄露和窗函數(shù)均會(huì)引起誤差。因此，要合理選取采樣間隔和截取長(zhǎng)度，必要時(shí)還需考慮適當(dāng)?shù)拇啊?duì)于連續(xù)周期信號(hào)，我們?cè)诓捎糜?jì)算機(jī)進(jìn)行計(jì)算時(shí)，也總是要進(jìn)行截?cái)啵蛄锌偸怯邢揲L(zhǎng)的，仍然可以采用上訴方法近似計(jì)算。4、可能用到MATLAB函數(shù)與代碼實(shí)驗(yàn)中的DFT運(yùn)算可以采用MATLAB中提供的FFT來(lái)實(shí)現(xiàn)。DTFT可以利用MATLAB矩陣運(yùn)算的方法進(jìn)行計(jì)算。三、實(shí)驗(yàn)內(nèi)容1.，完成如下要求：（1）計(jì)算其DTFT，并畫(huà)出區(qū)間的波形。（2）計(jì)算4點(diǎn)DFT，并把結(jié)果顯示在（1）所畫(huà)的圖形中。（3）對(duì)補(bǔ)零，計(jì)算64點(diǎn)DFT，并顯示結(jié)果。（4）是否可以由DFT計(jì)算DTFT，如果可以，請(qǐng)編程實(shí)現(xiàn)。代碼及圖形：（1）計(jì)算其DTFT，并畫(huà)出區(qū)間的波形。n=0:3;x=[2-111];w=-pi:0.001:pi;X=x*exp(-j*n'*w);subplot(211);plot(w,abs(X));xlabel('w');title('Magnitude');%畫(huà)出幅度譜>>axistightsubplot(212);plot(w,angle(X)/pi);>>xlabel('w');title('Phase');%畫(huà)出相位譜axistight（2）計(jì)算4點(diǎn)DFT，并把結(jié)果顯示在（1）所畫(huà)的圖形中。%只需在(1)的基礎(chǔ)上加上部分代碼，并稍作更改n=0:3;x=[2-111];w=-pi:0.001:2*pi;X=x*exp(-j*n'*w);Y=fft(x);>>subplot(211);plot(w,abs(X));>>holdon>>stem([0:pi/2:3*pi/2],abs(Y));%畫(huà)出4點(diǎn)DFT>>xlabel('w');title('Magnitude');>>axistight>>subplot(212);plot(w,angle(X)/pi);holdon>>stem([0:pi/2:3*pi/2],angle(Y));%畫(huà)出4點(diǎn)DFT>>xlabel('w');title('Phase');axistight（3）對(duì)補(bǔ)零，計(jì)算64點(diǎn)DFT，并顯示結(jié)果。n=0:3;x=[2-111];Y=fft(x,64);%為了計(jì)算64點(diǎn)DFTk=0:63;subplot(211);stem(k,abs(Y));xlabel('k');title('|X(k)|');axistightsubplot(212);stem(k,angle(Y));xlabel('k');title('φ(X(k))');axistight是否可以由DFT計(jì)算DTFT，如果可以，請(qǐng)編程實(shí)現(xiàn)。由（1）和（2）（3）可知，DFT相當(dāng)于在DTFT的頻譜上取樣。當(dāng)DFT的點(diǎn)數(shù)比較少時(shí)，無(wú)法根據(jù)根據(jù)DFT的頻譜得到理想的DTFT頻譜；如果通過(guò)采取在序列后面補(bǔ)零的方法，則可以提高采樣密度，使DFT的點(diǎn)數(shù)更多，從而得到精細(xì)的DTFT頻譜，且補(bǔ)零越多，頻譜越精細(xì)。2.考察序列（1）時(shí)，用DFT估計(jì)的頻譜；將補(bǔ)零加長(zhǎng)到長(zhǎng)度為100點(diǎn)序列用DFT估計(jì)的頻譜。要求畫(huà)出相應(yīng)波形。（2）時(shí)，用DFT估計(jì)的頻譜，并畫(huà)出波形。（3）根據(jù)實(shí)驗(yàn)結(jié)果，分析怎樣提高頻譜分辨率。代碼及圖形：（1）時(shí)，用DFT估計(jì)的頻譜；將補(bǔ)零加長(zhǎng)到長(zhǎng)度為100點(diǎn)序列用DFT估計(jì)的頻譜。要求畫(huà)出相應(yīng)波形。步驟一：時(shí)，用DFT估計(jì)的頻譜n=0:10;x=cos(0.48*pi*n)+cos(0.52*pi*n);>>Y=fft(x);>>subplot(211);plot(abs(Y));>>xlabel('0<=n<=10,DFT估計(jì)x(n)頻譜');title('Magnitude');>>axistight>>subplot(212);plot(angle(Y));>>xlabel('0<=n<=10,DFT估計(jì)x(n)頻譜');title('Phase');axistight步驟二：將補(bǔ)零加長(zhǎng)到長(zhǎng)度為100點(diǎn)序列用DFT估計(jì)的頻譜n=0:10;x=cos(0.48*pi*n)+cos(0.52*pi*n);>>Y=fft(x,100);>>subplot(211);plot(abs(Y));>>xlabel('0<=n<=10,補(bǔ)零至100點(diǎn)，DFT估計(jì)x(n)頻譜');title('Magnitude');>>axistight>>subplot(212);plot(angle(Y));>>xlabel('0<=n<=10,補(bǔ)零至100點(diǎn)，DFT估計(jì)x(n)頻譜');title('Phase');axistight（2）時(shí)，用DFT估計(jì)的頻譜，并畫(huà)出波形。n=0:100;x=cos(0.48*pi*n)+cos(0.52*pi*n);>>Y=fft(x);>>subplot(211);plot(abs(Y));>>xlabel('0<=n<=100，DFT估計(jì)x(n)頻譜');title('Magnitude');>>axistight>>subplot(212);plot(angle(Y));>>xlabel('0<=n<=100，DFT估計(jì)x(n)頻譜');title('Phase');axistight根據(jù)實(shí)驗(yàn)結(jié)果，分析怎樣提高頻譜分辨率。由以上實(shí)驗(yàn)結(jié)果可知可以通過(guò)增大截取長(zhǎng)度和增加補(bǔ)零的個(gè)數(shù)來(lái)使得序列列長(zhǎng)變長(zhǎng)，從而提高頻譜分辨率。已知信號(hào)，其中，，。從的表達(dá)式可以看出，它包含三個(gè)頻率的正弦波，但是，從其時(shí)域波形來(lái)看，似乎是一個(gè)正弦信號(hào)，利用DFT做頻譜分析，確定適合的參數(shù)，使得到的頻譜的頻率分辨率符合需要。代碼及圖形：t=0:0.1:15;x=0.15*sin(2*pi*t)+sin(2*pi*2*t)-0.1*sin(2*pi*3*t);>>Y=fft(x);>>subplot(211);stem(abs(Y));>>xlabel('DFT估計(jì)x(n)頻譜');title('Magnitude');>>axistight>>subplot(212);stem(angle(Y));>>xlabel('DFT估計(jì)x(n)頻譜');title('Phase');axistight利用DFT近似分析連續(xù)時(shí)間信號(hào)的頻譜（幅度譜）。分析采用不同的采樣間隔和截取長(zhǎng)度進(jìn)行計(jì)算的結(jié)果，并最終確定合適的參數(shù)。代碼及圖形：n=0:2:120;x=exp(-0.1*n);X=fft(x);k=0:60;stem(k,abs(X));xlabel('k');title('Magnitude');%采樣區(qū)間可以為0到120，采樣間隔為2。實(shí)驗(yàn)總結(jié)本實(shí)驗(yàn)研究怎樣用DFT分析信號(hào)頻譜的問(wèn)題。在理論學(xué)習(xí)中，已經(jīng)了解到DFT由DFS取主值序列得到，而DFS可以從DTFT推導(dǎo)而來(lái)，這說(shuō)明DFT與DTFT之間有著密切關(guān)系。實(shí)際上，DTFT的頻譜抽樣后能夠得到DFT的頻譜；即DTF的點(diǎn)數(shù)越多，越能通過(guò)擬合得到精細(xì)的DTFT頻譜。在實(shí)驗(yàn)過(guò)程中，我對(duì)提高頻譜分辨率的方法有了更加深刻的認(rèn)識(shí)，通過(guò)序列補(bǔ)零和選取合適的截取長(zhǎng)度都可以達(dá)到目的。當(dāng)然在實(shí)驗(yàn)中也遇到不少困難，由于編程能力有限，我對(duì)一些函數(shù)的用法不甚了解，在查找資料后才恍然大悟。我相信這是一個(gè)不斷提高能力的過(guò)程。實(shí)驗(yàn)2利用FFT計(jì)算線性卷積一、實(shí)驗(yàn)?zāi)康?.掌握利用FFT計(jì)算線性卷積的原理及具體實(shí)現(xiàn)方法。2.加深理解重疊相加法和重疊保留法。3.考察利用FFT計(jì)算線性卷積各種方法的適用范圍。二、實(shí)驗(yàn)原理1.線性卷積與圓周卷積設(shè)x(n)為L(zhǎng)點(diǎn)序列，h(n)為M點(diǎn)序列，x(n)和h(n)的線性卷積為（3-1）的長(zhǎng)度為L(zhǎng)+M-1x(n)和h(n)的圓周卷積為（3-2）圓周卷積與線性卷積相等而不產(chǎn)生交疊的必要條件為(3-3)圓周卷積定理：根據(jù)DFT性質(zhì)，x(n)和h(n)的N點(diǎn)圓周卷積的DFT等于它們的DFT的乘積：(3-4)2.快速卷積快速卷積發(fā)運(yùn)用圓周卷積實(shí)現(xiàn)線性卷積，根據(jù)圓周卷積定理利用FFT算法實(shí)現(xiàn)圓周卷積?？蓪⒖焖倬矸e運(yùn)算的步驟歸納如下：(1)必須選擇；為了能使用基-2算法，要求。采用補(bǔ)零的辦法使得x(n)和h(n)的長(zhǎng)度均為N。(2)計(jì)算x(n)和h(n)的N點(diǎn)FFT。(3)組成乘積(4)利用IFFT計(jì)算Y(k)的IDFT，得到線性卷積y(n)3.分段卷積我們考察單位取樣響應(yīng)為h(n)的線性系統(tǒng)，輸入為x(n)，輸出為y(n)，則當(dāng)輸入序列x(n)極長(zhǎng)時(shí)，如果要等x(n)全部集齊時(shí)再開(kāi)始進(jìn)行卷積，會(huì)使輸出有較大延時(shí)；如果序列太長(zhǎng)，需要大量存儲(chǔ)單元。為此，我們把x(n)分段，為別求出每段的卷積，合在一起得到最后的總輸出。這稱(chēng)為分段卷積。分段卷積可以細(xì)分為重疊保留法和重疊相加法。重疊保留法：設(shè)x(n)的長(zhǎng)度為，h(n)的長(zhǎng)度為M。把序列x(n)分成多段N點(diǎn)序列，每段雨前一段重寫(xiě)M-1個(gè)樣本。并在第一個(gè)輸入段前面補(bǔ)M-1個(gè)零。計(jì)算每一段與h(n)的圓周卷積，其結(jié)果中前M-1個(gè)不等與線性卷積，應(yīng)當(dāng)舍去，只保留后面N-M+1個(gè)正確的輸出樣本，把它們合起來(lái)得到總的輸出。利用FFT實(shí)現(xiàn)重疊保留法的步驟如下：(1)在x(n)前面填充M-1個(gè)零，擴(kuò)大以后的序列為(2)將x(n)分為若干段N點(diǎn)子段，設(shè)L=N-M+1為每一段的有效長(zhǎng)度，則第i段的數(shù)據(jù)為：(3)計(jì)算每一段與h(n)的N點(diǎn)圓周卷積，利用FFT計(jì)算圓周卷積(4)舍去每一段卷積結(jié)果的前M-1個(gè)樣本，連接剩下的樣本得到卷積結(jié)果y(n)。重疊相加法：設(shè)h(n)長(zhǎng)度為M，將信號(hào)x(n)分解成長(zhǎng)為L(zhǎng)的子段。以表示沒(méi)斷信號(hào)，則：每一段卷積的長(zhǎng)度為L(zhǎng)+M-1，所以在做求和時(shí)，相鄰兩段序列由M-1個(gè)樣本重疊，即前一段的最后M-1個(gè)樣本和下一段前M-1個(gè)樣本序列重疊，這個(gè)重疊部分相加，再與不重疊的部分共同組成y(n)。利用FFT實(shí)現(xiàn)重疊保留法的步驟如下：(1)將x(n)分為若干L點(diǎn)子段。(2)計(jì)算每一段與h(n)的卷積，根據(jù)快速卷積法利用FFT計(jì)算卷積。(3)將各段相加，得到輸出y(n)。4、可能得到的MATLAB函數(shù)實(shí)驗(yàn)中FFT運(yùn)算可采用MATLAB中提供的函數(shù)fft來(lái)實(shí)現(xiàn)。三、實(shí)驗(yàn)內(nèi)容假設(shè)要計(jì)算序列x(n)=u(n)-u(n-L),0≤n≤L和h(n)=cos(0.2πn),0≤n≤M的線性卷積完成以下實(shí)驗(yàn)內(nèi)容。1.設(shè)L=M,根據(jù)線性卷積的表達(dá)式和快速卷積的原理分別編程實(shí)現(xiàn)計(jì)算兩個(gè)序列線性卷積的方法,比較當(dāng)序列長(zhǎng)度分別為8,16,32,64,256,512,1024時(shí)兩種方法計(jì)算線性卷積所需時(shí)間。2．當(dāng)L=2048且M=256時(shí)比較直接計(jì)算線性卷積和快速卷積所需的時(shí)間,進(jìn)一步考察當(dāng)L=4096且M=256時(shí)兩種算法所需的時(shí)間。3.編程實(shí)現(xiàn)利用重疊相加法計(jì)算兩個(gè)序列的線性卷積,考察L=2048且M=256時(shí)計(jì)算線性卷積的時(shí)間,與2題的結(jié)果進(jìn)行比較。4.編程實(shí)現(xiàn)利用重疊保留法計(jì)算兩個(gè)序列的線性卷積,考察L=2048且M=256時(shí)計(jì)算線性卷積的時(shí)間,與2題的結(jié)果進(jìn)行比較。題目一：快速卷積法：L=input('L=');M=input('M=');x=1;h=1;fori=1:L-1x=[1x];h=[hcos(0.2*pi*i)];endNum=pow2(nextpow2(M+L-1));ticxs=fft(x,Num);hs=fft(h,Num);ys=xs.*hs;y=ifft(ys,Num);toc結(jié)果：L=8M=8時(shí)間已過(guò)0.004177秒。L=16M=16時(shí)間已過(guò)0.007234秒。L=32M=32時(shí)間已過(guò)0.012298秒。L=64M=64時(shí)間已過(guò)0.016842秒。L=256M=256時(shí)間已過(guò)0.018824秒。L=512M=512時(shí)間已過(guò)0.020002秒。L=1024M=1024時(shí)間已過(guò)0.024016秒。線性卷積法：L=input('L=');M=input('M=');n=0:L;x=heaviside(n)-heaviside(n-L);m=0:M;h=cos(0.2*pi*m);ticy=conv(x,h);toc結(jié)果：L=8;M=8;時(shí)間已過(guò)0.002205秒。L=16;M=16;時(shí)間已過(guò)0.008691秒。L=32;M=32;時(shí)間已過(guò)0.014703秒。L=64;M=64;時(shí)間已過(guò)0.021067秒。L=256;M=256;時(shí)間已過(guò)0.036894秒。L=512;M=512;時(shí)間已過(guò)0.058499秒。L=1024;M=1024;時(shí)間已過(guò)0.080003秒。結(jié)論：在L和M相等的情況下，使用快速卷積比使用線性卷積運(yùn)算更快。題目二：快速卷積法：L=input('L=');M=input('M=');x=1;h=1;fori=1:L-1x=[1x];h=[hcos(0.2*pi*i)];endNum=pow2(nextpow2(M+L-1));ticxs=fft(x,Num);hs=fft(h,Num);ys=xs.*hs;y=ifft(ys,Num);toc結(jié)果：L=2048M=256時(shí)間已過(guò)0.007841秒。L=4096M=256時(shí)間已過(guò)0.009291秒。線性卷積法：L=2048;M=256;n=0:L;x=heaviside(n)-heaviside(n-L);m=0:M;h=cos(0.2*pi*m);ticy=conv(x,h);toc結(jié)果：時(shí)間已過(guò)0.050958秒。L=4096;M=256;n=0:L;x=heaviside(n)-heaviside(n-L);m=0:M;h=cos(0.2*pi*m);ticy=conv(x,h);toc時(shí)間已過(guò)0.080958秒。結(jié)論：在L=2048，M=256；L=4096，M=256;兩種情況下仍然是快速卷積運(yùn)算更快。題目三：L=input('L=');M=input('M=');N=M+1;x=1;h=1;forj=1:L-1x=[1x];h=[hcos(0.2*pi*j)];endx=[x,zeros(1,N-1)];t=zeros(1,M-1);H=zeros(1,L+M-1);a=floor(L/N);ticfork=0:aA=x(k*N+1:k*N+N);H1=fft(A,L+M-1);H2=fft(h,L+M-1);HH=H1.*H2;q=ifft(HH,L+M-1);H(k*N+1:k*N+M-1)=q(1:M-1)+t(1:M-1);H(k*N+M:k*N+N)=q(M:N);t(1:M-1)=q(N+1:N+M-1);endtoc結(jié)果：L=2048M=256時(shí)間已過(guò)0.003830秒。題目四：L=input('L=');M=input('M=');N=M+1;x=1;y=1;fori=1:L-1x=[1x];y=[ycos(0.2*pi*i)];endLL=N+M-1;h=[h,zeros(1,N-M)];x=[zeros(1,M-1),x,zeros(1,N-1)];a=floor((L+M-2)/(LL))+1;Y=zeros(1,N);ticfork=0:a-1xk=x(k*LL+1:k*LL+N);b=fft(xk,N);C=fft(h,N);Z=b.*C;Y(k+1,:)=ifft(Z,N);endtoc結(jié)果：L=2048M=256時(shí)間已過(guò)0.007472秒。結(jié)論：L=2048，M=256時(shí)使用重疊相加法和重疊保留法比使用快速卷積的方法運(yùn)算更快一點(diǎn)。實(shí)驗(yàn)總結(jié)本實(shí)驗(yàn)研究了利用FFT計(jì)算線性卷積的方法。通過(guò)比較線性卷積和快速卷積兩種方法，我深深地體會(huì)到了快速卷積的優(yōu)越性，它可以大大提高計(jì)算效率?？焖倬矸e的核心思想是“分而治之”，通過(guò)利用旋轉(zhuǎn)因子的對(duì)稱(chēng)性和周期性合并一些運(yùn)算項(xiàng)，并且可以使長(zhǎng)序列的DFT分更小點(diǎn)的DFT.盡管實(shí)驗(yàn)的結(jié)果是對(duì)課堂已有結(jié)論的證明，但我仍在實(shí)踐中找到探索發(fā)現(xiàn)的樂(lè)趣。實(shí)驗(yàn)3IIR數(shù)字濾波器設(shè)計(jì)一、實(shí)驗(yàn)?zāi)康?.掌握利用脈沖響應(yīng)不變法和雙線性變換法設(shè)計(jì)IIR數(shù)字濾波器的原理及具體方法。2.加深理解數(shù)字濾波器和模擬濾波器之間的技術(shù)指標(biāo)轉(zhuǎn)化。3.掌握脈沖響應(yīng)不變法和雙線性變換法設(shè)計(jì)IIR數(shù)字濾波器的優(yōu)缺點(diǎn)及適用范圍。二、實(shí)驗(yàn)原理1、基本原理（1）從時(shí)域響應(yīng)出發(fā)，使數(shù)字濾波器的單位脈沖響應(yīng)h(n)模仿模擬濾波器的單位沖激響應(yīng)，h(n)等于的取樣值。（2）從頻率響應(yīng)出發(fā)，直接使數(shù)字濾波器的頻率響應(yīng)逼近模擬濾波器的頻率響應(yīng)，進(jìn)而求得數(shù)字濾波器的系統(tǒng)函數(shù)H（z）。2、變換方法脈沖響應(yīng)不變法思路：(1)將進(jìn)行部分分式展開(kāi)(2)對(duì)進(jìn)行拉式變換(3)對(duì)時(shí)域采樣得到h(n)(4)對(duì)h(n)進(jìn)行z變換雙線性變換法思路：通過(guò)正切變換Ω=tan(Ω1T/2)將s平面的整個(gè)虛軸jΩ壓縮到s1平面jΩ1軸上的-jπ/T到j(luò)π/T段上。通過(guò)z=es1T變換將s1平面映射到z平面上，得到s平面和z平面的映射關(guān)系為：s=(1-z-1)/(1+z-1)一般來(lái)說(shuō)為了調(diào)節(jié)模擬頻帶與數(shù)字頻帶之間的關(guān)系，可引入常數(shù)c,使映射關(guān)系為s=c(1-z-1)/(1+z-1)。設(shè)計(jì)時(shí)可采用模擬頻率特性與數(shù)字頻率特性在低頻出有確切對(duì)應(yīng)關(guān)系的常數(shù)c=2/T。設(shè)計(jì)步驟(1)確定數(shù)字濾波器的性能指標(biāo)。(2)將數(shù)字濾波器頻率指標(biāo)轉(zhuǎn)換成響應(yīng)的模擬濾波器頻率指標(biāo)(3)根據(jù)指標(biāo)，，和設(shè)計(jì)模擬濾波器。(4)將展成部分分式形式。(5)把模擬極點(diǎn)轉(zhuǎn)換成數(shù)字極點(diǎn)，得到數(shù)字濾波器?？梢?jiàn)至H(z)間的變換關(guān)系為方法1：利用residue函數(shù)和residuez函數(shù)實(shí)現(xiàn)脈沖響應(yīng)不變變換法，實(shí)用方法如下：[r,p,k]=residue(b,a)[b,a]=residue(r,p,k)實(shí)現(xiàn)多項(xiàng)式形式和部分分式形式之間的裝換[r,p,k]=residuez(b,a)[b,a]=residuez(r,p,k)實(shí)現(xiàn)多項(xiàng)式形式和部分分式形式之間的轉(zhuǎn)換方法2：（1）matlab中提供了impinvar函數(shù)采用脈沖響應(yīng)不變法實(shí)現(xiàn)模擬濾波器到數(shù)字濾波器的變換，其使用如下：[bz,az]=impinvar(b,a,fs)采用脈沖響應(yīng)不變法將模擬濾波器系統(tǒng)函數(shù)的系數(shù)向量b和a變換成為數(shù)字濾波器系統(tǒng)函數(shù)的系數(shù)向量bz和az，fs為采樣頻率（默認(rèn)為1）。[bz,az]=impinvar(b,a)采樣頻率默認(rèn)為1的情況下，采用脈沖響應(yīng)不變法將模擬濾波器變換為數(shù)字濾波器。（2）matlab中提供了bilinear函數(shù)采用雙線性變換法實(shí)現(xiàn)模擬濾波器到數(shù)字濾波器的變換，其使用如下：[bz,az]=bilinear(b,a,fs)采用雙線性變換法將模擬濾波器系統(tǒng)函數(shù)的系數(shù)向量b和a變換成為數(shù)字濾波器系統(tǒng)函數(shù)的系數(shù)向量bz和az，fs為采樣頻率。實(shí)驗(yàn)內(nèi)容設(shè)采樣設(shè)采樣頻率為fs=10kHz，設(shè)計(jì)數(shù)字低通濾波器，滿足如下指標(biāo)通帶截止頻率：fp=1kHz，通帶波動(dòng)：Rp=1dB阻帶截止頻率：fst=1.5kHz，阻帶衰減：As=15dB要求分別采用巴特沃斯、切比雪夫I型、切比雪夫II型和橢圓模擬原型濾波器，并分別結(jié)合脈沖響應(yīng)不變法和雙線性變換法進(jìn)行設(shè)計(jì)。結(jié)合實(shí)驗(yàn)結(jié)果，分別討論采用上述設(shè)計(jì)的數(shù)字濾波器是否都能滿足給定指標(biāo)要求，分析脈沖響應(yīng)不變法和雙線性變換法設(shè)計(jì)IIR數(shù)字濾波器的優(yōu)缺點(diǎn)及適用范圍。代碼及圖形：脈沖響應(yīng)不變法：巴特沃斯型：Wp=1000*2*pi;Ws=1500*2*pi;Rp=1;As=15;[NOmegaC]=buttord(Wp,Ws,Rp,As,'s');[ba]=butter(N,OmegaC,'s');[bzaz]=impinvar(b,a,10000);w=[0:500]*pi/500;[Hw]=freqz(bz,az);subplot(221);plot(w/pi,abs(H));gridon;xlabel('\omega(\pi)');ylabel('|H(e^j^\omega)|');subplot(222);plot(w/pi,20*log10(abs(H)/max(abs(H))));gridon;xlabel('\omega(\pi)');ylabel('|H(e^j^\omega)|(dB)');subplot(223);plot(w/pi,angle(H)/pi);gridon;xlabel('\omega(\pi)');ylabel('PhaseofH(e^j^\omega)(\pi)');grd=grpdelay(bz,az,w);%計(jì)算群延時(shí)響應(yīng)subplot(224);plot(w/pi,grd);gridon;xlabel('\omega(\pi)');ylabel('Groupdelay');切比雪夫Ⅰ型：Wp=1000*2*pi;Ws=1500*2*pi;Rp=1;As=15;[NOmegaC]=cheb1ord(Wp,Ws,Rp,As,'s');[b,a]=cheby1(N,Rp,OmegaC,'s');[bzaz]=impinvar(b,a,10000);w=[0:500]*pi/500;[Hw]=freqz(bz,az);subplot(221);plot(w/pi,abs(H));gridon;xlabel('\omega(\pi)');ylabel('|H(e^j^\omega)|');subplot(222);plot(w/pi,20*log10(abs(H)/max(abs(H))));gridon;xlabel('\omega(\pi)');ylabel('|H(e^j^\omega)|(dB)');subplot(223);plot(w/pi,angle(H)/pi);gridon;xlabel('\omega(\pi)');ylabel('PhaseofH(e^j^\omega)(\pi)');grd=grpdelay(bz,az,w);subplot(224);plot(w/pi,grd);gridon;xlabel('\omega(\pi)');ylabel('Groupdelay');切比雪夫Ⅱ型：Wp=1000*2*pi;Ws=1500*2*pi;Rp=1;As=15;[NOmegaC]=cheb2ord(Wp,Ws,Rp,As,'s');[b,a]=cheby2(N,As,OmegaC,'s');[bzaz]=impinvar(b,a,10000);w=[0:500]*pi/500;[Hw]=freqz(bz,az);subplot(221);plot(w/pi,abs(H));gridon;xlabel('\omega(\pi)');ylabel('|H(e^j^\omega)|');subplot(222);plot(w/pi,20*log10(abs(H)/max(abs(H))));gridon;xlabel('\omega(\pi)');ylabel('|H(e^j^\omega)|(dB)');subplot(223);plot(w/pi,angle(H)/pi);gridon;xlabel('\omega(\pi)');ylabel('PhaseofH(e^j^\omega)(\pi)');grd=grpdelay(bz,az,w);subplot(224);plot(w/pi,grd);gridon;xlabel('\omega(\pi)');ylabel('Groupdelay');橢圓型：Wp=1000*2*pi;Ws=1500*2*pi;Rp=1;As=15;[NOmegaC]=ellipord(Wp,Ws,Rp,As,'s');[b,a]=ellip(N,Rp,As,OmegaC,'s');[bzaz]=impinvar(b,a,10000);w=[0:500]*pi/500;[Hw]=freqz(bz,az);subplot(221);plot(w/pi,abs(H));gridon;xlabel('\omega(\pi)');ylabel('|H(e^j^\omega)|');subplot(222);plot(w/pi,20*log10(abs(H)/max(abs(H))));gridon;xlabel('\omega(\pi)');ylabel('|H(e^j^\omega)|(dB)');subplot(223);plot(w/pi,angle(H)/pi);gridon;xlabel('\omega(\pi)');ylabel('PhaseofH(e^j^\omega)(\pi)');grd=grpdelay(bz,az,w);subplot(224);plot(w/pi,grd);gridon;xlabel('\omega(\pi)');ylabel('Groupdelay');雙線性變換法：巴特沃斯型：Wp=1000*2*pi;Ws=1500*2*pi;Rp=1;As=15;[NOmegaC]=buttord(Wp,Ws,Rp,As,'s');[b,a]=butter(N,OmegaC,'s');[bzaz]=bilinear(b,a,10000);w=[0:500]*pi/500;[Hw]=freqz(bz,az);subplot(221);plot(w/pi,abs(H));gridon;xlabel('\omega(\pi)');ylabel('|H(e^j^\omega)|');subplot(222);plot(w/pi,20*log10(abs(H)/max(abs(H))));gridon;xlabel('\omega(\pi)');ylabel('|H(e^j^\omega)|(dB)');subplot(223);plot(w/pi,angle(H)/pi);gridon;xlabel('\omega(\pi)');ylabel('PhaseofH(e^j^\omega)(\pi)');grd=grpdelay(bz,az,w);%計(jì)算群延時(shí)響應(yīng)subplot(224);plot(w/pi,grd);gridon;xlabel('\omega(\pi)');ylabel('Groupdelay');切比雪夫Ⅰ型：Wp=1000*2*pi;Ws=1500*2*pi;Rp=1;As=15;[NOmegaC]=cheb1ord(Wp,Ws,Rp,As,'s');[b,a]=cheby1(N,Rp,OmegaC,'s');[bzaz]=bilinear(b,a,10000);w=[0:500]*pi/500;[Hw]=freqz(bz,az);subplot(221);plot(w/pi,abs(H));gridon;xlabel('\omega(\pi)');ylabel('|H(e^j^\omega)|');subplot(222);plot(w/pi,20*log10(abs(H)/max(abs(H))));gridon;xlabel('\omega(\pi)');ylabel('|H(e^j^\omega)|(dB)');subplot(223);plot(w/pi,angle(H)/pi);gridon;xlabel('\omega(\pi)');ylabel('PhaseofH(e^j^\omega)(\pi)');grd=grpdelay(bz,az,w);subplot(224);plot(w/pi,grd);gridon;xlabel('\omega(\pi)');ylabel('Groupdelay');切比雪夫Ⅱ型：Wp=1000*2*pi;Ws=1500*2*pi;Rp=1;As=15;[NOmegaC]=cheb2ord(Wp,Ws,Rp,As,'s');[b,a]=cheby2(N,As,OmegaC,'s');[bzaz]=bilinear(b,a,10000);w=[0:500]*pi/500;[Hw]=freqz(bz,az);subplot(221);plot(w/pi,abs(H));gridon;xlabel('\omega(\pi)');ylabel('|H(e^j^\omega)|');subplot(222);plot(w/pi,20*log10(abs(H)/max(abs(H))));gridon;xlabel('\omega(\pi)');ylabel('|H(e^j^\omega)|(dB)');subplot(223);plot(w/pi,angle(H)/pi);gridon;xlabel('\omega(\pi)');ylabel('PhaseofH(e^j^\omega)(\pi)');grd=grpdelay(bz,az,w);subplot(224);plot(w/pi,grd);gridon;xlabel('\omega(\pi)');ylabel('Groupdelay');橢圓型：Wp=1000*2*pi;Ws=1500*2*pi;Rp=1;As=15;[NOmegaC]=ellipord(Wp,Ws,Rp,As,'s');[b,a]=ellip(N,Rp,As,OmegaC,'s');[bzaz]=bilinear(b,a,10000);w=[0:500]*pi/500;[Hw]=freqz(bz,az);subplot(221);plot(w/pi,abs(H));gridon;xlabel('\omega(\pi)');ylabel('|H(e^j^\omega)|');subplot(222);plot(w/pi,20*log10(abs(H)/max(abs(H))));gridon;xlabel('\omega(\pi)');ylabel('|H(e^j^\omega)|(dB)');subplot(223);plot(w/pi,angle(H)/pi);gridon;xlabel('\omega(\pi)');ylabel('PhaseofH(e^j^\omega)(\pi)');grd=grpdelay(bz,az,w);subplot(224);plot(w/pi,grd);gridon;xlabel('\omega(\pi)');ylabel('Groupdelay');采用脈沖響應(yīng)不變法時(shí)：巴特沃斯型低通濾波器滿足題目要求：在0到1kHz時(shí)滿足通帶波動(dòng)1dB；在1.5kHz滿足阻帶衰減15dB。切比雪夫I型低通濾波器滿足題目要求：在0到1kHz時(shí)滿足通帶波動(dòng)1dB；在1.5kHz滿足阻帶衰減15dB。切比雪夫II型低通濾波器不滿足題目要求，通帶波動(dòng)明顯大于1dB，阻帶衰減明顯小于15dB。橢圓型低通濾波器不滿足題目要求，通帶波動(dòng)明顯大于1dB，阻帶衰減明顯小于15dB。優(yōu)點(diǎn)：頻率坐標(biāo)的變換是線性的：。因此，如果模擬濾波器的頻響是限帶于折疊頻率以內(nèi)的話，則通過(guò)變換后所得的數(shù)字濾波器的頻響可以不失真的反映原響應(yīng)與頻率的關(guān)系。缺點(diǎn)：脈沖響應(yīng)不變法的最大缺點(diǎn)是有頻譜的周期性延拓效應(yīng)，因此只能用于限帶的頻響特性。對(duì)于高通、帶阻濾波器，由于它們?cè)诟哳l部分不衰減，當(dāng)一定要追求頻率線性關(guān)系而采用脈沖響應(yīng)不變法時(shí)，必須先對(duì)模擬高通和帶阻濾波器加一層保護(hù)濾波器，濾除高于折疊頻率以上的頻率，然后再轉(zhuǎn)變?yōu)閿?shù)字濾波器以避免混疊失真。雙線性變換法：巴特沃斯型低通濾波器滿足題目要求：在0到1kHz時(shí)滿足通帶波動(dòng)1dB；在1.5kHz滿足阻帶衰減15dB。切比雪夫I型低通濾波器滿足題目要求：在0到1kHz時(shí)滿足通帶波動(dòng)1dB；在1.5kHz滿足阻帶衰減15dB。切比雪夫II型低通濾波器不滿足題目要求，橢圓型低通濾波器不滿足題目要求。優(yōu)點(diǎn)：消除了脈沖響應(yīng)不變法所固有的混疊誤差。缺點(diǎn)：頻率Ω與w間的非線性。這種非線性關(guān)系要求被變換的連續(xù)時(shí)間系統(tǒng)的幅度響應(yīng)必須是分段常數(shù)型的，不然所映射出的數(shù)字頻率響應(yīng)相對(duì)于原來(lái)的模擬頻率響應(yīng)會(huì)產(chǎn)生變形。解決雙線性變換中的頻率非線性關(guān)系的方法有：在低頻段，模擬與數(shù)字兩濾波器的頻率關(guān)系處于近似的線性范圍內(nèi)，故可忽略非線性影響；或者采用補(bǔ)償?shù)姆椒?，也可稱(chēng)為預(yù)畸的方法。四、實(shí)驗(yàn)總結(jié)本實(shí)驗(yàn)研究了用脈沖響應(yīng)不變法和雙線性變換法設(shè)計(jì)IIR濾波器的問(wèn)題。在這次實(shí)驗(yàn)中我對(duì)兩種方法的優(yōu)缺點(diǎn)有了更加熟練的掌握，并且比較分析了巴特沃斯低通濾波器、切比雪夫I型低通濾波器、切比雪夫II型低通濾波器、橢圓低通濾波器的不同之處。理論知識(shí)的學(xué)習(xí)在實(shí)踐中得到鞏固。在分析實(shí)驗(yàn)中的問(wèn)題時(shí)，（比如分析脈沖響應(yīng)不變法的優(yōu)缺點(diǎn)），我才體會(huì)到老師上課所反復(fù)講的這方面的知識(shí)在實(shí)際設(shè)計(jì)時(shí)產(chǎn)生的結(jié)果。因?yàn)橛行〇|西只靠做題時(shí)無(wú)法深刻理解的，只有具體操作過(guò)才能明白其中含義。這次實(shí)驗(yàn)讓我明白了DSP的學(xué)習(xí)不能生搬硬套，而應(yīng)該自己多動(dòng)手，多時(shí)間，方可達(dá)到把這門(mén)課學(xué)精、學(xué)通的境界。實(shí)驗(yàn)4FIR數(shù)字濾波器設(shè)計(jì)一、實(shí)驗(yàn)?zāi)康?掌握窗函數(shù)法和頻率取樣法設(shè)計(jì)FIR數(shù)字濾波器的原理及具體方法實(shí)驗(yàn)原理窗函數(shù)法： 1.基本原理 窗函數(shù)設(shè)計(jì)法的基本思想為，首先選擇一個(gè)適當(dāng)?shù)睦硐氲臑V波器，然后用窗函數(shù)截取它的單位脈沖響應(yīng)，得到線性相位和因果的FIR濾波器。這種方法的重點(diǎn)是選擇一個(gè)合適的窗函數(shù)和理想濾波器，使設(shè)計(jì)的濾波器的單位脈沖響應(yīng)逼近理想濾波器的單位脈沖響應(yīng)。 2.設(shè)計(jì)步驟 （1）給定理想濾波器的頻率響應(yīng)，在通帶上具有單位增益和線性相位，在阻帶上具有零響應(yīng)。一個(gè)帶寬為的低通濾波器由下式給定：其中為采樣延遲，其作用是為了得到一個(gè)因果系統(tǒng)。（2）確定這個(gè)濾波器的單位脈沖響應(yīng)為了得到一個(gè)長(zhǎng)度為N的因果的線性相位FIR濾波器，我們令（3）用窗函數(shù)截取得到所設(shè)計(jì)FIR數(shù)字濾波器：3.窗函數(shù)的選擇常用的窗函數(shù)有矩形（Rectangular）窗，漢寧（Hanning）窗，海明（Hamming）窗、布萊克曼（Blackman）窗、凱瑟（Kaiser）窗等。MATLAB提供了一些函數(shù)用于產(chǎn)生窗函數(shù)。如表7-1所示。表7-1MATLAB中產(chǎn)生窗函數(shù)的命令MATLAB函數(shù)窗函數(shù)MATLAB函數(shù)窗函數(shù)Boxcar矩形窗函數(shù)Blackman布萊克曼窗Hanning漢寧窗函數(shù)Kaiser凱瑟窗函數(shù)Hamming海明窗在設(shè)計(jì)過(guò)程中我們需要根據(jù)給定的濾波器技術(shù)指標(biāo)，選擇濾波器長(zhǎng)N和窗函數(shù)。表7-2列出了常用窗函數(shù)的一些特性，可供設(shè)計(jì)參考。表7-2常用窗函數(shù)的特性窗函數(shù)窗函數(shù)頻率特性加窗后濾波器指標(biāo)旁瓣峰值dB主瓣寬度過(guò)渡帶寬最小阻帶衰減dB矩形窗-134π/N1.8π/N-21漢寧窗-318π/N6.2π/N-44海明窗-418π/N6.6π/N-53布萊克曼窗-5712π/N11π/N-74凱瑟窗是一種廣泛在實(shí)際中廣泛應(yīng)用的窗函數(shù)，它由下式給定：其中是修正的零階貝塞爾函數(shù)，參數(shù)控制最小阻帶衰減，這種窗函數(shù)對(duì)于相同的N可以提供不同的過(guò)渡帶寬。由于貝塞爾函數(shù)比較復(fù)雜，這種窗函數(shù)的設(shè)計(jì)方程很難推導(dǎo)，然而幸運(yùn)的是，有一些經(jīng)驗(yàn)設(shè)計(jì)方程可以直接使用。已知給定的指標(biāo)，濾波器長(zhǎng)度N和凱瑟窗參數(shù)可以按如下凱瑟窗方程給出過(guò)渡帶帶寬：頻率取樣法：基本原理：頻率取樣法從頻域出發(fā)，把理想的濾波器Hd(ejw)等間隔采樣得到Hd(k）,將Hd(k）作為實(shí)際設(shè)計(jì)濾波器的H（k）。H(k）=Hd(k）=H(ejw)|w=2πk/Nk=0,1,…N-1得到H(k）以后可以由H(k）唯一確定濾波器的單位脈沖響應(yīng)h(n),H(ejw)也可以由H(k）求得：h(n)=IDFT[H（k）],H(ejw)=如果我們?cè)O(shè)計(jì)的是線性相位FTR濾波器，則H（k)的幅度和相位一定滿足線性相位濾波器的約束條件。我們將H(k)表示成如下形式：H(k)=|H(k)|ejθ(k)=Hr(k)ejθ(k)當(dāng)h(n)為實(shí)數(shù)，則H(k)=H*(N-k),由此得到：Hr(k)=Hr(N-k)三、實(shí)驗(yàn)內(nèi)容1、設(shè)計(jì)一個(gè)數(shù)字低通FIR濾波器，其技術(shù)指標(biāo)如下：分別采用矩形窗、漢寧窗、海明窗、布萊克曼窗、凱瑟窗設(shè)計(jì)該濾波器。結(jié)合實(shí)驗(yàn)結(jié)果，分別討論采用上述方法設(shè)計(jì)的數(shù)字濾波器是否都能滿足給定指標(biāo)要求。矩形窗wp=0.2*pi;wst=0.3*pi;rp=0.25;as=50;tr_width=wst-wp;N=ceil(1.8*pi/tr_width)+1;n=0:(N-1);wc=(wp+wst)/2;alpha=(N-1)/2;hd=(wc/pi)*sinc((wc/pi)*(n-alpha));w_boxcar=boxcar(N)';%產(chǎn)生窗函數(shù)，并對(duì)hd(n)加窗得到h(n)h=hd.*w_boxcar;subplot(221);stem(n,hd,'filled');axistight;xlabel('n');ylabel('hd(n)');[Hr,w1]=zerophase(h);subplot(222);plot(w1/pi,Hr);axistight;xlabel('\omega/\pi');ylabel('H(\omega)');subplot(223);stem(n,h,'filled');axistight;xlabel('n');ylabel('h(n)');[H,w]=freqz(h,1);subplot(224);plot(w/pi,20*log10(abs(H)/max(abs(H))));xlabel('\omega/\pi');ylabel('dB');gridon>>漢寧窗：wp=0.2*pi;wst=0.3*pi;rp=0.25;as=50;tr_width=wst-wp;N=ceil(6.2*pi/tr_width)+1;n=0:(N-1);wc=(wp+wst)/2;alpha=(N-1)/2;hd=(wc/pi)*sinc((wc/pi)*(n-alpha));w_hanning=hanning(N)';%產(chǎn)生窗函數(shù)，并對(duì)hd(n)加窗得到h(n)h=hd.*w_hanning;subplot(221);stem(n,hd,'filled');axistight;xlabel('n');ylabel('hd(n)');[Hr,w1]=zerophase(h);subplot(222);plot(w1/pi,Hr);axistight;xlabel('\omega/\pi');ylabel('H(\omega)');subplot(223);stem(n,h,'filled');axistight;xlabel('n');ylabel('h(n)');[H,w]=freqz(h,1);subplot(224);plot(w/pi,20*log10(abs(H)/max(abs(H))));xlabel('\omega/\pi');ylabel('dB');gridon>>海明窗：wp=0.2*pi;wst=0.3*pi;rp=0.25;as=50;tr_width=wst-wp;N=ceil(6.6*pi/tr_width)+1;n=0:(N-1);wc=(wp+wst)/2;alpha=(N-1)/2;hd=(wc/pi)*sinc((wc/pi)*(n-alpha));w_hamming=hamming(N)';%產(chǎn)生窗函數(shù)，并對(duì)hd(n)加窗得到h(n)h=hd.*w_hamming;subplot(221);stem(n,hd,'filled');axistight;xlabel('n');ylabel('hd(n)');[Hr,w1]=zerophase(h);subplot(222);plot(w1/pi,Hr);axistight;xlabel('\omega/\pi');ylabel('H(\omega)');subplot(223);stem(n,h,'filled');axistight;xlabel('n');ylabel('h(n)');[H,w]=freqz(h,1);subplot(224);plot(w/pi,20*log10(abs(H)/max(abs(H))));xlabel('\omega/\pi');ylabel('dB');gridon>>布萊克曼窗：wp=0.2*pi;wst=0.3*pi;rp=0.25;as=50;tr_width=wst-wp;N=ceil(11*pi/tr_width)+1;n=0:(N-1);wc=(wp+wst)/2;alpha=(N-1)/2;hd=(wc/pi)*sinc((wc/pi)*(n-alpha));w_blackman=blackman(N)';%產(chǎn)生窗函數(shù)，并對(duì)hd(n)加窗得到h(n)h=hd.*w_blackman;subplot(221);stem(n,hd,'filled');axistight;xlabel('n');ylabel('hd(n)');[Hr,w1]=zerophase(h);subplot(222);plot(w1/pi,Hr);axistight;xlabel('\omega/\pi');ylabel('H(\omega)');subplot(223);stem(n,h,'filled');axistight;xlabel('n');ylabel('h(n)');[H,w]=freqz(h,1);subplot(224);plot(w/pi,20*log10(abs(H)/max(abs(H))));xlabel('\omega/\pi');ylabel('dB');gridon>>凱瑟窗：wp=0.2*pi;wst=0.3*pi;rp=0.25;as=50;tr_width=wst-wp;N=ceil((as-7.95)/(2.285*tr_width))+1;n=0:(N-1);wc=(wp+wst)/2;alpha=(N-1)/2;hd=(wc/pi)*sinc((wc/pi)*(n-alpha));w_kaiser=kaiser(N)';%產(chǎn)生窗函數(shù)，并對(duì)hd(n)加窗得到h(n)h=hd.*w_kaiser;subplot(221);stem(n,hd,'filled');axistight;xlabel('n');ylabel('hd(n)');[Hr,w1]=zerophase(h);subplot(222);plot(w1/pi,Hr);axistight;xlabel('\omega/\pi');ylabel('H(\omega)');subplot(223);stem(n,h,'filled');axistight;xlabel('n');ylabel('h(n)');[H,w]=freqz(h,1);subplot(224);plot(w/pi,20*log10(abs(H)/max(abs(H))));xlabel('\omega/\pi');ylabel('dB');gridon>>2、設(shè)計(jì)一個(gè)數(shù)字帶通FIR濾波器，其技術(shù)指標(biāo)如下：下阻帶邊緣：下通帶邊緣：上通帶邊緣：上阻帶邊緣：用海明窗設(shè)計(jì)：wp1=0.35*pi;wst1=0.2*pi;wp2=0.65*pi;wst2=0.8*pi;tr_width=wp1-wst1;N=ceil(6.6*pi/tr_width)+1;n=0:(N-1);alpha=(N-1)/2;wc1=(wp1+wst1)/2;wc2=(wp2+wst2)/2;hd=(wc2/pi)*sinc((wc2/pi)*(n-alpha))-(wc1/pi)*sinc((wc1/pi)*(n-alpha));w_hamming=hamming(N)';h=hd.*w_hamming;subplot(221);stem(n,hd,'filled');axistight;xlabel('n');ylabel('hd(n)');[Hr,w1]=zerophase(h);subplot(222);plot(w1/pi,Hr);axistight;xlabel('\omega/\pi');ylabel('H(\omega)');subplot(223);stem(n,h,'filled');axistight;xlabel('n');ylabel('h(n)');[H,w]=freqz(h,1);subplot(224);plot(w/pi,20*log10(abs(H)/max(abs(H))));xlabel('\omega/\pi');ylabel('dB');gridon3.采用頻率取樣法設(shè)計(jì)FIR數(shù)字低通濾波器，滿足以下指標(biāo)：Wp=0.2π,Rp=0.25dBWst=0.3π,As=50dB取N=20，過(guò)渡帶沒(méi)有樣本取N=40，過(guò)渡帶有一個(gè)樣本，T=0.39取N=60，過(guò)渡帶有兩個(gè)樣本，T1=0.5925，T2=0.1099代碼及圖形：（1）取N=20，過(guò)渡帶沒(méi)有樣本N=20;alpha=(N-1)/2;l=0:N-1;wl=(2*pi/N)*l;Hrs=[1,1,1,zeros(1,15),1,1];%由于過(guò)渡帶無(wú)樣本Hdr=[1,1,0,0];wdl=[0,0.25,0.25,1];k1=0:floor((N-1)/2);k2=floor((N-1)/2)+1:N-1;angH=[-alpha*(2*pi)/N*k1,alpha*(2*pi)/N*(N-k2)];H=Hrs.*exp(j*angH);h=ifft(H,N);w=[0:500]*pi/500;H=freqz(h,l,w);[Hr,wr]=zerophase(h);subplot(221);plot(wdl,Hdr,wl(1:11)/pi,Hrs(1:11),'o');axis([0,1,-0.1,1.1]);xlabel('\omega(\pi)');ylabel('Hr(k)');subplot(222);stem(l,h,'filled');axis([0,N-1,-0.1,0.3]);xlabel('n');ylabel('h(n)');subplot(223);plot(wr/pi,Hr,wl(1:11)/pi,Hrs(1:11),'o');axis([0,1,-0.2,1.2]);xlabel('\omega(\pi)');ylabel('Hr(ω)');subplot(224);plot(w/pi,20*log10((abs(H)/max(abs(H)))));axis([0,1,-50,5]);grid;xlabel('\omega(\pi)');ylabel('dB');（2）取N=40，過(guò)渡帶有一個(gè)樣本，T=0.39N=40;alpha=(N-1)/2;l=0:N-1;wl=(2*pi/N)*l;Hrs=[1,1,1,1,1,0.39,zeros(1,29),0.39,1,1,1,1];%由于過(guò)渡帶有一個(gè)樣本Hdr=[1,1,0.39,0,0];wdl=[0,0.2,0.25,0.3,1];k1=0:floor((N-1)/2);k2=floor((N-1)/2)+1:N-1;angH=[-alpha*(2*pi)/N*k1,alpha*(2*pi)/N*(N-k2)];H=Hrs.*exp(j*angH);h=ifft(H,N);w=[0:500]*pi/500;H=freqz(h,l,w);[Hr,wr]=zerophase(h);subplot(221);plot(wdl,Hdr,wl(1:21)/pi,Hrs(1:21),'o');axis([0,1,-0.1,1.1]);xlabel('\omega(\pi)');ylabel('Hr(k)');subplot(222);stem(l,h,'filled');axis([0,N-1,-0.1,0.3]);xlabel('n');ylabel('h(n)');subplot(223);plot(wr/pi,Hr,wl(1:21