用雙線性變換法設(shè)計(jì)原型低通為橢圓型的數(shù)字IIR帶通濾波器重點(diǎn)

1、 課 程 設(shè) 計(jì)課程設(shè)計(jì)名稱：用雙線性變換法設(shè)計(jì)原型低通為 橢圓型的數(shù)字IIR帶通濾波器 專 業(yè) 班 級(jí) ： xxxxxx 學(xué) 生 姓 名 ： xxx 學(xué) 號(hào) ： xxxxxxxxxxxx 指 導(dǎo) 教 師 ： xxx 課程設(shè)計(jì)時(shí)間： 2014.6.16至2014.6.20 電子信息工程 專業(yè)課程設(shè)計(jì)任務(wù)書(shū)學(xué)生姓名XXX專業(yè)班級(jí)XXXX 學(xué)號(hào)XXXXXX題 目用雙線性變換法設(shè)計(jì)原型低通為橢圓型的數(shù)字IIR低通濾波器課題性質(zhì)其他課題來(lái)源自擬指導(dǎo)教師XX同組姓名主要內(nèi)容用雙線性變換法設(shè)計(jì)原型低通為橢圓型的數(shù)字IIR低通濾波器，要求通帶邊界頻率為400Hz，阻帶邊界頻率分別為500Hz，通帶最大衰減1

2、dB，阻帶最小衰減40dB，抽樣頻率為2000Hz，用MATLAB畫(huà)出幅頻特性，畫(huà)出并分析濾波器傳輸函數(shù)的零極點(diǎn)；信號(hào)經(jīng)過(guò)該濾波器，其中300Hz，600Hz，濾波器的輸出是什么？用Matlab驗(yàn)證你的結(jié)論并給出的圖形。任務(wù)要求1、掌握用雙線性變換法設(shè)計(jì)原型低通為橢圓型的數(shù)字IIR低通濾波器的原理和設(shè)計(jì)方法。2、求出所設(shè)計(jì)濾波器的Z變換。3、用MATLAB畫(huà)出幅頻特性圖。4、驗(yàn)證所設(shè)計(jì)的濾波器。參考文獻(xiàn)1、程佩青著，數(shù)字信號(hào)處理教程，清華大學(xué)出版社，20012、Sanjit K. Mitra著，孫洪，余翔宇譯，數(shù)字信號(hào)處理實(shí)驗(yàn)指導(dǎo)書(shū)（MATLAB版），電子工業(yè)出版社，2005年1月3、郭仕劍等

3、，MATLAB 7.x數(shù)字信號(hào)處理，人民郵電出版社，2006年4、胡廣書(shū)，數(shù)字信號(hào)處理 理論算法與實(shí)現(xiàn)，清華大學(xué)出版社，2003年審查意見(jiàn)指導(dǎo)教師簽字：教研室主任簽字： 2014 年06 月 12 日 說(shuō)明：本表由指導(dǎo)教師填寫(xiě)，由教研室主任審核后下達(dá)給選題學(xué)生，裝訂在設(shè)計(jì)（論文）首頁(yè)一需求分析和技術(shù)要求 1、需求分析濾波器從廣義上來(lái)說(shuō)對(duì)特定的頻點(diǎn)或頻點(diǎn)以外的頻率進(jìn)行有效濾波的電路，這種電路保留輸入信號(hào)中的有用信息，濾除不需要的信息，從而達(dá)到信號(hào)的檢測(cè)、提取、識(shí)別等不同的目的。如果處理的信號(hào)是時(shí)域離散信號(hào)，那么相應(yīng)的處理系統(tǒng)就稱為數(shù)字濾波器，由于在實(shí)際工作中被處理的信號(hào)都是幅度量化的數(shù)字信號(hào)，因

4、此，數(shù)字濾波器實(shí)際上是用有限精度的算法實(shí)現(xiàn)一個(gè)線性時(shí)不變的時(shí)域離散系統(tǒng)。數(shù)字濾波器的種類(lèi)很多，分類(lèi)方法也不同，可以從功能上分類(lèi)，也可以從實(shí)現(xiàn)方法上及設(shè)計(jì)方法上來(lái)分類(lèi)等等。濾波器在功能上總的可分為四類(lèi)，即低通(LP)、高通(HP)、帶通(BP)、帶阻(BS)濾波器等，從實(shí)現(xiàn)方法上，由有限長(zhǎng)沖激響應(yīng)的數(shù)字濾波器被稱為FIR濾波器，具有無(wú)限長(zhǎng)沖激響應(yīng)的數(shù)字濾波器稱作IIR濾波器。橢圓濾波器（Elliptic filter）又稱考爾濾波器（Causer filter），是在通帶和阻帶等波紋的一種濾波器。橢圓濾波器相比其他類(lèi)型的濾波器，在階數(shù)相同的條件下有著最小的通帶和阻帶波動(dòng)。它在通帶和阻帶的波動(dòng)相同

5、，這一點(diǎn)區(qū)別于在通帶和阻帶都平坦的巴特沃斯濾波器，以及阻帶平坦、通帶等波紋的切比雪夫?yàn)V波器。利用雙線性變換法將模擬傳輸信號(hào)Ha(s)變換為數(shù)字傳輸函數(shù)G(z)，從而是z域的數(shù)字傳輸函數(shù)保留s域的模擬傳輸函數(shù)的基本性質(zhì)。設(shè)計(jì)成的IIR數(shù)字低通濾波器能夠去掉信號(hào)中不必要的高頻成分，降低采樣頻率，避免頻率混淆，去掉高頻干擾。2、技術(shù)要求 1、掌握用雙線性變換法設(shè)計(jì)原型低通為橢圓型的數(shù)字IIR低通濾波器的原理和設(shè)計(jì)方法。2、求出所設(shè)計(jì)濾波器的Z變換。3、用MATLAB畫(huà)出幅頻特性圖。4、驗(yàn)證所設(shè)計(jì)的濾波器。二設(shè)計(jì)原理及設(shè)計(jì)思路1、設(shè)計(jì)原理數(shù)字濾波器是一種用來(lái)過(guò)濾時(shí)間離散信號(hào)的數(shù)字系統(tǒng)，通過(guò)對(duì)抽樣數(shù)據(jù)進(jìn)

6、行數(shù)學(xué)處理來(lái)達(dá)到頻域?yàn)V波的目的。數(shù)字濾波器也是具有一定傳輸選擇特性的數(shù)字信號(hào)處理裝置，其輸入、輸出均為數(shù)字信號(hào)，實(shí)質(zhì)上是一個(gè)由有限精度算法實(shí)現(xiàn)的線性時(shí)不變離散系統(tǒng)。IIR數(shù)字濾波器采用遞歸型結(jié)構(gòu)，即結(jié)構(gòu)上帶有反饋環(huán)路。IIR濾波器運(yùn)算結(jié)構(gòu)通常由延時(shí)、乘以系數(shù)和相加等基本運(yùn)算組成，可以組合成直接型、正準(zhǔn)型、級(jí)聯(lián)型、并聯(lián)型四種結(jié)構(gòu)形式，都具有反饋回路。數(shù)字濾波器根據(jù)其沖激響應(yīng)函數(shù)的時(shí)域特性，可分為兩種，即無(wú)限長(zhǎng)沖激響應(yīng)(IIR)數(shù)字濾波器和有限長(zhǎng)沖激響應(yīng)(FIR)數(shù)字濾波器。IIR 數(shù)字濾波器的特征是，具有無(wú)限持續(xù)時(shí)間沖激響應(yīng)，需要用遞歸模型來(lái)實(shí)現(xiàn)，其差分方程為： （2.1） （2.2）設(shè)計(jì)IIR

7、濾波器的任務(wù)就是尋求一個(gè)物理上可實(shí)現(xiàn)的系統(tǒng)函數(shù)H(z)，使其頻率響應(yīng)H(z)滿足所希望得到的頻域指標(biāo)，即符合給定的通帶截止頻率、阻帶截止頻率、通帶衰減系數(shù)和阻帶衰減系數(shù)。本次課程設(shè)計(jì)采用雙線性變換法。雙線性變換法設(shè)計(jì)IIR數(shù)字濾波器雙線性變換法主要是采用非線性頻率壓縮方法，將整個(gè)頻率軸上的頻率范圍壓縮到-/T/T之間，再用z=esT轉(zhuǎn)換到Z平面上。也就是說(shuō)，第一步先將整個(gè)S平面壓縮映射到S1平面的-/T/T一條橫帶里；第二步再通過(guò)標(biāo)準(zhǔn)變換關(guān)系z(mì)=es1T將此橫帶變換到整個(gè)Z平面上去。這樣就使S平面與Z平面建立了一一對(duì)應(yīng)的單值關(guān)系，消除了多值變換性，也就消除了頻譜混疊現(xiàn)象，映射關(guān)系如圖1所示。

8、圖1雙線性變換的映射關(guān)系為了將S平面的整個(gè)虛軸j壓縮到S1平面j1軸上的-/T到/T段上，可以通過(guò)以下的正切變換實(shí)現(xiàn) (2.3)式中,T仍是采樣間隔。當(dāng)1由-/T經(jīng)過(guò)0變化到/T時(shí)，由-經(jīng)過(guò)0變化到+，也即映射了整個(gè)j軸。將式（2.3）寫(xiě)成(2.4)將此關(guān)系解析延拓到整個(gè)S平面和S1平面，令j=s，j1=s1，則得 （2.5）再將S1平面通過(guò)以下標(biāo)準(zhǔn)變換關(guān)系映射到Z平面z=es1T (2.6) 從而得到S平面和Z平面的單值映射關(guān)系為： （2.7） （2.8） 式（2.7）與式（2.8）是S平面與Z平面之間的單值映射關(guān)系，這種變換都是兩個(gè)線性函數(shù)之比，因此稱為雙線性變換式（2.6）與式（2.7）

9、的雙線性變換符合映射變換應(yīng)滿足的兩點(diǎn)要求。首先,把z=ej，可得 （2.9）即S平面的虛軸映射到Z平面的單位圓。其次，將s=+j代入式（2.8），得 (2.10)因此(2.11)由此看出，當(dāng)<0時(shí)，|z|<1；當(dāng)>0時(shí)，|z|>1。也就是說(shuō)，S平面的左半平面映射到Z平面的單位圓內(nèi)，S平面的右半平面映射到Z平面的單位圓外，S平面的虛軸映射到Z平面的單位圓上。因此，穩(wěn)定的模擬濾波器經(jīng)雙線性變換后所得的數(shù)字濾波器也一定是穩(wěn)定的。2、設(shè)計(jì)思路具體實(shí)現(xiàn)方法如下：（1）確定數(shù)字低通濾波器的技術(shù)指標(biāo)：通帶截止頻率p、通帶衰減Rp、阻帶截止頻率s、阻帶衰減Rs。（2）將數(shù)字低通濾波器的

10、技術(shù)指標(biāo)轉(zhuǎn)換成模擬低通濾波器的技術(shù)指標(biāo)。如果采用脈沖響不變法，邊界頻率的轉(zhuǎn)換關(guān)系為: (2.12)如果采用雙線性變換法，邊界頻率的轉(zhuǎn)換關(guān)系為(2.13)（3）按照模擬低通濾波器的技術(shù)指標(biāo)設(shè)計(jì)模擬低通濾波器。（4）利用雙線性變換法將模擬濾波器Ha(s)，從s平面轉(zhuǎn)換到z平面，得到數(shù)字低通濾波器系統(tǒng)函數(shù)H(z)。（5）數(shù)字低通技術(shù)指標(biāo)為：p=0.4rad， Rp=1dB； s=0.5rad， Rs=40dB（6）模擬低通的技術(shù)指標(biāo)為：通帶截止頻率p=(2/T)*tan(p/2)=tan(0.392)=0.0068阻帶截止頻率s=(2/T)*tan(s/2)=tan(0.628)=0.011通帶最小

11、衰減Rp=1dB和阻帶最大衰減Rs=40dB；歸一化截止角頻率p=2pi*Fs/Ft； s=2pi*Fs/Ft（7） 利用模擬橢圓濾波器設(shè)計(jì)數(shù)字濾波器。通帶截止頻率為：p=0.4*pi；阻帶截止頻率為：s=0.5*pi；通帶最大衰減為：Rp=1；阻帶最大衰減為：Rs=40；設(shè)定周期為1s；模擬低通濾波器的生成：num1,den1 =ellip(N,Rp,Rs,Wn,'s')；最后實(shí)現(xiàn)輸入輸出、幅頻特性、相頻特性的圖形。三程序流程圖確定數(shù)字低通濾波器的各項(xiàng)性能指標(biāo)開(kāi)始數(shù)字低通濾波器的技術(shù)指標(biāo)轉(zhuǎn)換成模擬低通濾波器的技術(shù)指標(biāo)按照模擬低通濾波器的技術(shù)指標(biāo)設(shè)計(jì)模擬低通濾波器利用雙線性變換

12、法將模擬濾波器轉(zhuǎn)變?yōu)閿?shù)字低通濾波器利用模擬橢圓濾波器設(shè)計(jì)數(shù)字濾波器結(jié)束圖2 濾波器設(shè)計(jì)流程方框圖四程序源代碼MATLAB 程序如下：clear all close all Coeff_LF = 1; % 低頻信號(hào)系數(shù) F1(300 Hz)Coeff_HF = 1; % 高頻信號(hào)系數(shù) F2(600 Hz) Coeff_NOISE = 0; % 噪聲系數(shù), 如0(沒(méi)有噪聲), 1, sqrt(2), 3等等. F1 = 300; % HzF2 = 600; % HzFS = 2000; % Hz 抽樣頻率 % 時(shí)間定義EndTime = 1/FS * 1023;Time = 0:1/FS:EndT

13、ime;N = length(Time); % 產(chǎn)生正弦信號(hào) x1(t),x2(t)和 x(t)SinSignal_x1= Coeff_LF*sin(2*pi*F1*Time)% x1(t)SinSignal_x2= Coeff_HF*sin(2*pi*F2*Time)% x2(t)SinSignalSum = SinSignal_x1 + SinSignal_x2+ Coeff_NOISE*randn(1, N);SinSignal = SinSignalSum;% x(t) % 繪制信號(hào) x1(t)figuresubplot(3,1,1)plot(Time, SinSignal_x1)xl

14、abel('時(shí)間')xlim(0 0.5)title('時(shí)域信號(hào) x1(t)') subplot(3,1,2)plot(Time, SinSignal_x1)xlabel('時(shí)間')xlim(0 0.05)title('時(shí)域信號(hào) x1(t)') FFTSinSignal_x1 = 2/N*fft(SinSignal_x1)subplot(3,1,3)FnSin_x1 = abs(FFTSinSignal_x1); FW = 2*pi*(Time)/EndTime; plot(FW(2:length(FW)*FS/pi/2, FnS

15、in_x1(2:length(FW)hold onplot(FW(2:length(FW)-2*pi)*FS/pi/2, FnSin_x1(2:length(FW)xlim(-FS/2 FS/2) xlabel('頻率') % 繪制 x2(t)figuresubplot(3,1,1)plot(Time, SinSignal_x2)xlabel('時(shí)間')xlim(0 0.5)title('時(shí)域信號(hào) x(2)') subplot(3,1,2)plot(Time, SinSignal_x2)xlabel('時(shí)間')xlim(0 0.0

16、5)title('時(shí)域信號(hào) x2(t)') FFTSinSignal_x2 = 2/N*fft(SinSignal_x2)subplot(3,1,3)FnSin_x2 = abs(FFTSinSignal_x2); plot(FW(2:length(FW)*FS/pi/2, FnSin_x2(2:length(FW)hold onplot(FW(2:length(FW)-2*pi)*FS/pi/2, FnSin_x2(2:length(FW)xlim(-FS/2 FS/2) xlabel('頻率') % 繪制信號(hào) x(t)figuresubplot(3,1,1)

17、plot(Time, SinSignal)xlabel('時(shí)間')xlim(0 0.5)title('時(shí)域信號(hào)') subplot(3,1,2)plot(Time, SinSignal)xlabel('時(shí)間')xlim(0 0.05)title('時(shí)域信號(hào) x(t)') FFTSinSignal = 2/N*fft(SinSignal)subplot(3,1,3)FnSin = abs(FFTSinSignal); plot(FW(2:length(FW)*FS/pi/2, FnSin(2:length(FW)hold onplo

18、t(FW(2:length(FW)-2*pi)*FS/pi/2, FnSin(2:length(FW)xlim(-FS/2 FS/2) xlabel('頻率') %-設(shè)計(jì)低通濾波器 - Rp = 1; % 通帶衰減 /dBRs = 40; % 阻帶衰減 /dBOmegaP1_1=400; % 通帶邊界頻率OmegaS1_1=500; % 阻帶邊界頻率FS=2000; % 抽樣頻率 Wp1=2*pi*OmegaP1_1/FS; % 將模擬頻率變?yōu)閿?shù)字角頻率Ws1=2*pi*OmegaS1_1/FS; % 將模擬頻率變?yōu)閿?shù)字角頻率OmegaP1=tan(Wp1/2); % 非線性變

19、換OmegaS1=tan(Ws1/2); % 非線性變換Eta_P=OmegaP1/OmegaP1; % 歸一化Eta_S=OmegaS1/OmegaP1; % 歸一化 Lemta_P_EquivalentLowPass=1/Eta_P; % 變?yōu)榈刃У屯▍?shù)Lemta_S_EquivalentLowPass=1/Eta_S; % 變?yōu)榈刃У屯▍?shù) % 估計(jì)濾波器階數(shù)N, Wn = Ellipord(Lemta_P_EquivalentLowPass, Lemta_S_EquivalentLowPass,Rp, Rs,'s'); % 時(shí)間濾波器num1,den1 = ellip

20、(N,Rp,Rs,Wn,'s');num2,den2=lp2lp(num1,den1,OmegaP1);num,den=bilinear(num2,den2,0.5); % 計(jì)算增益響應(yīng)w = 0:pi/255:pi;h = freqz(num,den,w);g = 20*log10(abs(h); % 繪制增益響應(yīng)figure;plot(w/pi,g);grid axis(0 1 -60 5);xlabel('omega /pi'); ylabel('增益 /dB');title('橢圓低通濾波器的增益響應(yīng)'); % 繪制零點(diǎn)和

21、極點(diǎn)figure;z,p,k=tf2zp(num,den);zplane(z,p);title('零點(diǎn)和極點(diǎn)') figuresubplot(2,1,1)YLP = filter(num,den, SinSignal);plot(Time, YLP)xlim(0 0.05)xlabel('時(shí)間')title('低通濾波后的時(shí)域信號(hào)') FFTSinSignalSum = 2 / N * fft(YLP)subplot(2,1,2)FnSinSum = abs(FFTSinSignalSum); plot(FW(2:length(FW)*FS/pi/2, FnSinSum(2:length(FW)hold onplot(FW(2:length(FW)-2*pi)*FS/pi/2, FnSinSum(2:length(FW)xlim(-FS/2 FS/2) ylim(0 1)xlabel('頻率') %return%-End Program-五仿真結(jié)果圖 六 參考資料1程佩青著.數(shù)字信號(hào)處理教程