劉文斐
(齊魯師范學(xué)院 物理與電子工程學(xué)院,山東 濟(jì)南 250200)
濾波器,顧名思義,是用于“濾波”的,即將某些頻率范圍的信號(hào)進(jìn)行保留,抑制其他頻率的信號(hào).按照處理信號(hào)種類的不同,濾波器分為數(shù)字濾波器和模擬濾波器兩大類,其中,數(shù)字濾波器是通過(guò)程序運(yùn)算將數(shù)據(jù)進(jìn)行處理[1],對(duì)數(shù)字信號(hào)進(jìn)行濾波處理得到期望響應(yīng)特性的離散時(shí)間系統(tǒng).實(shí)際上,數(shù)字濾波器不僅可以處理數(shù)字信號(hào),如果在濾波器的輸入端接入模-數(shù)轉(zhuǎn)換器,它也可對(duì)模擬信號(hào)進(jìn)行處理.與模擬濾波器相比,數(shù)字濾波器具有高精度、高可靠性、可程控改變特性或復(fù)用、便于集成等優(yōu)點(diǎn),廣泛應(yīng)用于語(yǔ)音、圖像、醫(yī)學(xué)生物信號(hào)處理等領(lǐng)域.?dāng)?shù)字濾波器按照其系統(tǒng)中沖激響應(yīng)的特點(diǎn)可分為無(wú)限長(zhǎng)單位沖激響應(yīng)(infinite impulse response,IIR)濾波器和有限長(zhǎng)單位沖激響應(yīng)(finite impulse response,F(xiàn)IR)濾波器.IIR濾波器具有結(jié)構(gòu)簡(jiǎn)單、效率高、與模擬濾波器有對(duì)應(yīng)關(guān)系、易于解析控制及計(jì)算機(jī)輔助設(shè)計(jì)等優(yōu)點(diǎn)[2],與FIR濾波器相比,它的相同階次濾波效果更好.
IIR數(shù)字濾波器的設(shè)計(jì)方法[3]有間接法和直接法,一般采用間接法進(jìn)行設(shè)計(jì).間接法是通過(guò)模擬濾波器設(shè)計(jì)數(shù)字濾波器,這是因?yàn)槟M濾波器的設(shè)計(jì)方法較為成熟,有完整的設(shè)計(jì)公式以及數(shù)據(jù)表格可以直接使用,可使數(shù)字濾波器的設(shè)計(jì)更為簡(jiǎn)便快捷.直接法設(shè)計(jì)IIR數(shù)字濾波器,如果給定的設(shè)計(jì)要求比較復(fù)雜,需通過(guò)計(jì)算機(jī)輔助設(shè)計(jì)求解線性或非線性的聯(lián)立方程組,計(jì)算量大.
間接法設(shè)計(jì)IIR數(shù)字濾波器有3種設(shè)計(jì)方案[3].第1種方案分3步:首先通過(guò)給定的技術(shù)指標(biāo)設(shè)計(jì)歸一化樣本的模擬低通濾波器;然后通過(guò)模擬頻帶變換成其他模擬各種(低通、高通、帶通、帶阻)濾波器;最后將模擬濾波器數(shù)字化成相應(yīng)的數(shù)字濾波器.該方法的實(shí)現(xiàn)需要經(jīng)過(guò)“模擬—模擬頻帶變換”.第2種方案是直接由設(shè)計(jì)出的歸一化樣本模擬低通濾波器變換成各種通帶的數(shù)字濾波器,再將頻帶變換與數(shù)字化結(jié)合在一起,一步變換完成.該方法雖然簡(jiǎn)化了運(yùn)算步驟但實(shí)質(zhì)也是經(jīng)過(guò)“模擬—模擬頻帶變換”.第3種方案也分3步:首先通過(guò)給定的技術(shù)指標(biāo)設(shè)計(jì)歸一化樣本的模擬低通濾波器;然后直接數(shù)字化成數(shù)字低通濾波器;最后用數(shù)字頻帶變換成其他各種通帶的數(shù)字濾波器.該方法的實(shí)現(xiàn)需經(jīng)過(guò)“數(shù)字—數(shù)字頻帶變換”.這3種設(shè)計(jì)方案的根本區(qū)別在于頻帶變換的信號(hào)領(lǐng)域不同,前2種是在模擬信號(hào)域中進(jìn)行,且對(duì)模擬濾波器的運(yùn)算設(shè)計(jì)有可直接使用的公式和表格,所以是較常用的方法;第3種是在數(shù)字信號(hào)域中進(jìn)行頻帶變換,需找到不同數(shù)字頻帶間的函數(shù)關(guān)系.
在間接法設(shè)計(jì)數(shù)字濾波器的3種設(shè)計(jì)方案中,無(wú)論用哪一種對(duì)數(shù)字濾波器進(jìn)行設(shè)計(jì),都有一個(gè)共同的關(guān)鍵步驟——模擬濾波器的數(shù)字化過(guò)程[3],也就是將模擬濾波器的系統(tǒng)函數(shù)H(s)變換成所需的數(shù)字濾波器的系統(tǒng)函數(shù)H(z),即把s平面映射到z平面上.?dāng)?shù)字化方法可采用沖激響應(yīng)不變法、階躍響應(yīng)不變法和雙線性變換法,其中較常用的是沖激響應(yīng)不變法和雙線性變換法.
沖激響應(yīng)不變法是從濾波器的單位抽樣響應(yīng)出發(fā),使數(shù)字濾波器的單位抽樣響應(yīng)h(n)逼近模擬濾波器的單位抽樣響應(yīng)h(t),頻率間的變化是線性變換關(guān)系,克服了雙線性變換法中非線性頻率變換帶來(lái)的臨界頻率點(diǎn)的頻率畸變,所以是最簡(jiǎn)單、最直接的把s平面映射到z平面的映射方法.但是該方法要求模擬濾波器是嚴(yán)格帶限于抽樣頻率的1/2,如果不滿足該要求,數(shù)字濾波器的頻率響應(yīng)將產(chǎn)生混疊失真.
混疊失真是指在對(duì)連續(xù)信號(hào)進(jìn)行等間隔采樣時(shí)[4],如果不能滿足奈奎斯特采樣定理,即信號(hào)的最高頻率成分大于采樣頻率的1/2,那么采樣后信號(hào)的頻譜會(huì)發(fā)生頻譜分量的交疊現(xiàn)象,高于和低于采樣頻率的信號(hào)混雜在一起.如果發(fā)生頻譜混疊失真,原信號(hào)將無(wú)法從采樣信號(hào)中恢復(fù)還原.
如果通過(guò)“模擬—模擬頻帶變換”方法設(shè)計(jì)IIR數(shù)字高通或者帶阻濾波器,沖激響應(yīng)法確實(shí)會(huì)產(chǎn)生混疊失真現(xiàn)象,但如果通過(guò)“數(shù)字—數(shù)字頻帶變換”方法則不存在該問(wèn)題,因?yàn)樵摲椒ǖ臄?shù)字化過(guò)程是將模擬低通濾波器的系統(tǒng)函數(shù)映射為數(shù)字低通濾波器的系統(tǒng)函數(shù),而模擬低通濾波器是嚴(yán)格帶限于抽樣頻率的1/2、是抗混疊的,所以不會(huì)出現(xiàn)頻率混疊失真現(xiàn)象.
設(shè)計(jì)一個(gè)切貝雪夫Ⅰ型數(shù)字高通濾波器,要求通帶的截止頻率為4.5 kHz,阻帶上限截止頻率為3 kHz,通帶內(nèi)衰減不大于2 dB,阻帶內(nèi)衰減不小于15 dB,抽樣頻率為20 kHz.
如果用“模擬—模擬頻帶變換”進(jìn)行設(shè)計(jì),其設(shè)計(jì)流程如圖1所示.用MATLAB輔助設(shè)計(jì)、仿真設(shè)計(jì)結(jié)果,主要的設(shè)計(jì)步驟如圖2所示.
圖1 用“模擬—模擬頻帶變換”進(jìn)行設(shè)計(jì)的流程
圖2 用MATLAB輔助設(shè)計(jì)的仿真結(jié)果
用“模擬—模擬頻帶變換”方法、通過(guò)沖激響應(yīng)不變法完成模擬濾波器的數(shù)字化過(guò)程設(shè)計(jì)的高通濾波器,確實(shí)存在頻譜混疊失真,不符合設(shè)計(jì)參數(shù),達(dá)不到設(shè)計(jì)要求(圖3).為消除頻譜失真,用“數(shù)字—數(shù)字頻帶變換”方法進(jìn)行設(shè)計(jì)(圖4).
圖3 發(fā)生頻譜混疊失真的數(shù)字高通濾波器頻譜
由上述可知,“數(shù)字—數(shù)字頻帶變換”的實(shí)質(zhì)就是從數(shù)字低通濾波器的Z平面映射到另一個(gè)待求所需類型數(shù)字濾波器的z平面的變化關(guān)系,關(guān)鍵點(diǎn)是找到Z到z的映射函數(shù)[3]Z-1=G(z-1),則有
設(shè)Z=ejθ,z=ejω,則數(shù)字低通映射到數(shù)字高通的映射函數(shù)為
若低通濾波器的系統(tǒng)函數(shù)Hl(z)的分子分母皆為Z的多項(xiàng)式,即
Hl(Z)=
由此可得
用該表達(dá)式的形式在MATLAB輔助設(shè)計(jì)中自定義映射函數(shù),將數(shù)字低通濾波器系統(tǒng)函數(shù)Hl(z)映射為數(shù)字高通系統(tǒng)函數(shù)H(z)的形式(圖4),用MATLAB輔助設(shè)計(jì)的設(shè)計(jì)步驟如圖5所示.
圖4 用“數(shù)字—數(shù)字頻帶變換”進(jìn)行設(shè)計(jì)
圖5 用MATLAB輔助設(shè)計(jì)的仿真結(jié)果
關(guān)鍵程序:
%自定義mapping函數(shù):
function[bz,az]=mapping(bZ,aZ,Nz,Dz)
bzord=(length(bZ)-1); %原系統(tǒng)函數(shù)分子階數(shù)
bzord1=(length(bZ)-1)*(length(Nz)-1); %變換后系統(tǒng)函數(shù)分子階數(shù)
azord=(length(aZ)-1); %原系統(tǒng)函數(shù)分母階數(shù)
azord1=(length(aZ)-1)*(length(Dz)-1); %變換后系統(tǒng)函數(shù)分母階數(shù)
bz=zeros(1,bzord1+1); %分子系數(shù)向量初始化為零
for k=0:bzord %求各多項(xiàng)式乘積結(jié)果
pln=[1];
for l=0:k-1
pln=conv(pln,Nz); %求N(z)的k次乘積
end
pld=[1];
for l=0:bzord-k-1
pld=conv(pld,Dz); %求N(z)的bzord-k次乘積
end
bZZ=bZ(k+1)*conv(pln,pld);
bz=bz+bZ(k+1)*conv(pln,pld); %分子多項(xiàng)式系數(shù)向量求和
end
az=zeros(1,azord1+1)
for k=0:azord
pln=[1];
for l=0:k-1
pln=conv(pln,Nz);
end
pld=[1];
for l=0:azord-k-1
pld=conv(pld,Dz);
end
aZZ=aZ(k+1)*conv(pln,pld);
az=az+aZ(k+1)*conv(pln,pld); %分母多項(xiàng)式系數(shù)向量求和
end
az1=az(1);
az=az/az1;
bz=bz/az1;
%主程序
Fs=20000; %設(shè)置采樣頻率
rp=2; %設(shè)置通帶衰減
rs=15; %設(shè)置阻帶衰減
wp=2*pi*3000/Fs; %設(shè)置阻帶截止頻率
ws=2*pi*4500/Fs; %設(shè)置通帶截止頻率
OmegaP=2*Fs*tan(wp/2); %計(jì)算模擬低通通帶截止頻率
OmegaS=2*Fs*tan(ws/2); %計(jì)算模擬低通阻帶截止頻率
[N,OmegaC]=cheb1ord(wp,ws,rp,rs,'s'); %得到切貝雪夫Ⅰ型模擬低通濾波器階數(shù)和通帶截止頻率
[b,a]=cheby1(N,rp,OmegaC,'s'); %得到模擬低通濾波器系統(tǒng)函數(shù)
[bZ,aZ]=impinvar(b,a,Fs); %沖激響應(yīng)不變法得到數(shù)字低通濾波器系統(tǒng)函數(shù)
[bhp,ahp]=mapping(bZ,aZ,Nz,Dz); %調(diào)用mapping函數(shù)得到數(shù)字高通濾波器系統(tǒng)函數(shù)
設(shè)計(jì)結(jié)果如圖6所示.由圖6可知,該數(shù)字高通濾波器通帶的截止頻率約為4.5 kHz,通帶內(nèi)衰減小于5 dB,阻帶的截止頻率為3 kHz,阻帶衰減大于15 dB,設(shè)計(jì)指標(biāo)基本達(dá)到設(shè)計(jì)要求.所以用“數(shù)字—數(shù)字頻帶變換”能夠解決沖激響應(yīng)不變法帶來(lái)的頻譜混疊失真問(wèn)題,實(shí)現(xiàn)高通濾波器的設(shè)計(jì).
圖6 消除頻譜混疊失真的數(shù)字高通濾波器頻譜
設(shè)計(jì)一個(gè)切貝雪夫Ⅰ型數(shù)字帶阻濾波器,其技術(shù)指標(biāo):通帶的截止頻率分別為6,14 kHz,阻帶截止頻率分別為8,12 kHz,通帶內(nèi)衰減不大于2 dB,阻帶內(nèi)衰減不小于20 dB,抽樣頻率為400 kHz.
用“模擬—模擬頻帶變換”方法、通過(guò)沖激響應(yīng)不變法,完成模擬濾波器的數(shù)字化過(guò)程設(shè)計(jì)帶阻濾波器,頻譜混疊現(xiàn)象嚴(yán)重,無(wú)法達(dá)到設(shè)計(jì)要求.如果用“數(shù)字—數(shù)字頻帶變換”方法,則數(shù)字低通映射為數(shù)字帶阻的映射函數(shù)為[3]
θp為數(shù)字低通濾波器通帶的截止頻率;ωp1,ωp2為數(shù)字帶阻濾波器通帶的截止頻率.
用MATLAB輔助設(shè)計(jì),其設(shè)計(jì)步驟與數(shù)字高通濾波器的設(shè)計(jì)相似:
1)確定模擬低通濾波器的階次及截止頻率.
%設(shè)置帶阻濾波器參數(shù)
Fs=40000; %設(shè)置采樣頻率
S1=2*pi*8000/Fs; %設(shè)置阻帶截止頻率
S2=2*pi*12000/Fs;
P1=2*pi*6000/Fs; %設(shè)置阻帶截止頻率
P2=2*pi*14000/Fs;
rp=2; %設(shè)置通帶衰減
rs=20; %設(shè)置阻帶衰減
[N,Wc]=cheb1ord(P1,S1,rp,rs,'s'); %計(jì)算模擬低通濾波器的最小階數(shù)及通帶截止頻率
運(yùn)行得到模擬低通濾波器階數(shù)(N=6),通帶截止頻率為7 kHz.
2)創(chuàng)建模擬低通濾波器原型.
[b,a]=cheby1(N,rp,Wc,'s'); %由階數(shù)和通帶截止頻率得到模擬低通濾波器系統(tǒng)函數(shù)
3)沖激響應(yīng)不變法將模擬低通轉(zhuǎn)化為數(shù)字低通濾波器.
[bZ,aZ]=impinvar(b,a,Fs); %得到數(shù)字低通濾波器系統(tǒng)函數(shù)
4)調(diào)用函數(shù)mapping()實(shí)現(xiàn)數(shù)字帶阻濾波器系統(tǒng)函數(shù).
alpha=cos(p0); %設(shè)置參數(shù)α
klpha=(tan((P2-P1)/2))*(tan(Wc/2)); %設(shè)置參數(shù)k
d1=-2*alpha/(1+klpha);
d2=(1-klpha)/(1+klpha);
Nz=[d2,d1,1];
Dz=[1,d1,d2];
[bhp,ahp]=mapping(bZ,aZ,Nz,Dz);
設(shè)計(jì)結(jié)果如圖7所示.由圖7可知,頻率小于6 kHz和大于14 kHz的范圍屬于通帶范圍,通帶衰減小于5 dB,接近指標(biāo)值;頻率為8~12 kHz,屬于阻帶,衰減大于20 dB,基本滿足設(shè)計(jì)要求.
圖7 消除頻譜混疊失真的數(shù)字帶阻濾波器頻譜
將設(shè)計(jì)實(shí)現(xiàn)的帶阻濾波器應(yīng)用于雜波濾波,驗(yàn)證其濾波效果.
實(shí)現(xiàn)程序:
%引入雜波信號(hào)
t=(1:Fs)/Fs;
y=sin(2*pi*4000*t)+sin(2*pi*7000*t)+
sin(2*pi*10000*t)+
sin(2*pi*13000*t)+
sin(2*pi*15000*t);
%對(duì)雜波信號(hào)濾波,其中bZ,aZ為設(shè)計(jì)實(shí)現(xiàn)帶阻濾波器系統(tǒng)函數(shù)分式中分子、分母的系數(shù)
x=filter(bZ,aZ,y).
雜波波形如圖8所示.由圖8可知,雜波信號(hào)中包含5個(gè)頻譜分量,頻率中心分別為4,9,10,13,15 kHz.進(jìn)入帶阻濾波器進(jìn)行濾波,由于帶阻濾波器的通帶范圍為小于6 kHz和高于14 kHz,阻帶為8~12 kHz.所以經(jīng)過(guò)濾波后,雜波信號(hào)中10 kHz頻率處的頻譜分量被濾掉,9,13 kHz頻率的信號(hào)由于處于帶阻濾波器的過(guò)渡帶范圍內(nèi),頻譜被衰減,處于通帶范圍內(nèi)的信號(hào)沒(méi)有受到影響.濾波后信號(hào)如圖9所示.由圖8~9可知,通過(guò)“數(shù)字—數(shù)字頻帶變換”進(jìn)行設(shè)計(jì)實(shí)現(xiàn)的濾波器能夠達(dá)到設(shè)計(jì)要求,可實(shí)現(xiàn)濾波的功能.
圖8 雜波信號(hào)
圖9 濾波后信號(hào)
提出用“數(shù)字—數(shù)字頻帶變換”設(shè)計(jì)數(shù)字濾波器的方法.該方法在保證實(shí)現(xiàn)設(shè)計(jì)指標(biāo)的同時(shí),消除由沖激響應(yīng)不變法在將模擬濾波器系統(tǒng)函數(shù)映射為數(shù)字濾波器系統(tǒng)函數(shù)時(shí)帶來(lái)的頻譜混疊失真問(wèn)題.通過(guò)數(shù)字高通和數(shù)字帶阻濾波器的設(shè)計(jì)實(shí)例,驗(yàn)證該方法的有效性.用“數(shù)字—數(shù)字頻帶變換”方法設(shè)計(jì)數(shù)字濾波器,可滿足設(shè)計(jì)要求,是一種有效的抗混疊的設(shè)計(jì)方法.