當(dāng)前位置：首頁 > EDA > 電子設(shè)計(jì)自動化

用matlab設(shè)計(jì)的IIR濾波器源程序-IIR一階低通/高通

時(shí)間：2012-03-22 13:32:00
關(guān)鍵字： 低通    源程序    MATLAB   
手機(jī)看文章

掃描二維碼
隨時(shí)隨地手機(jī)看文章

[導(dǎo)讀]用matlab設(shè)計(jì)的IIR濾波器源程序 （1）IIR一階低通濾波器 P576 clear; fi=1;fs=10;Gc2=0.9; wc=2*pi*fi/fs; omegac=tan(wc/2); alpha=(sqrt(Gc2)/sqrt(1-Gc2))*omegac; a=(1-alpha)/(1+alpha); b=(1-a)/

用matlab設(shè)計(jì)的IIR濾波器源程序
（1）IIR一階低通濾波器 P576
clear;
fi=1;fs=10;Gc2=0.9;
wc=2*pi*fi/fs;    
omegac=tan(wc/2);
alpha=(sqrt(Gc2)/sqrt(1-Gc2))*omegac;
a=(1-alpha)/(1+alpha);
b=(1-a)/2;
w=0:pi/300:pi;    
Hw2=alpha^2./(alpha^2+(tan(w/2)).^2);
plot(w/pi,Hw2);
grid;
hold on;
 
（2）一階高通濾波器 P581
clear;
fi=1;fs=10;Gc2=0.5;
wc=2*pi*fi/fs;
omegac=tan(wc/2);
alpha=(sqrt(1-Gc2)/(sqrt(Gc2)))*omegac;
a=(1-alpha)/(1+alpha);
b=(1+a)/2;
w=0:pi/300:pi;
Hw2=(tan(w/2).^2)./(alpha^2+(tan(w/2)).^2);
plot(w/pi,Hw2);
grid;
hold on;

（3）Notch 嵌波濾波器
clear;
Gb2=0.5;
w0=0.35*pi;
deltaw=0.1*pi;
b=1/(1+tan(deltaw/2)*(sqrt(1-Gb2)/sqrt(Gb2)));
B=[1 -2*cos(w0) 1].*b;
A=[1 -2*b*cos(w0) (2*b-1)];
w=0:pi/500:pi;
H=freqz(B,A,w);
plot(w/pi,abs(H));
grid;

（4）Peak 濾波器
clear;
Ac=3;
Gb2=10^(-Ac/10);
w0=0.35*pi;
deltaw=0.1*pi;
b=1/(1+tan(deltaw/2)*(sqrt(Gb2)/sqrt(1-Gb2)));
B=[1 0 -1].*(1-b);
A=[1 -2*b*cos(w0) (2*b-1)];
w=0:pi/500:pi;
H=freqz(B,A,w);
plot(w/pi,abs(H));
grid;

（5）IIR低通濾波（Butterworth）
% IIR Lowpass Use Butterworth
% copyright by Etual
clear;
fs=20;fpass=4;fstop=5;
Ap=0.5;As=10;
wp=2*pi*fpass/fs;ws=2*pi*fstop/fs;
omegap=tan(wp/2);omegas=tan(ws/2);
ep=sqrt(10^(Ap/10)-1);
es=sqrt(10^(As/10)-1);
N=ceil(log(es/ep)/log(omegas/omegap));
omega0=omegap/ep^(1/N);
K=floor(N/2);
for i=1:K
    theta(i)=pi*(N-1+2*i)/(2*N);
end
for i=1:K
    G(i)=omega0^2/(1-2*omega0*cos(theta(i))+omega0^2);
end
for i=1:K
    a1(i)=2*(omega0^2-1)/(1-2*omega0*cos(theta(i))+omega0^2);
end
for i=1:K
    a2(i)=(1+2*omega0*cos(theta(i))+omega0^2)/(1-2*omega0*cos(theta(i))+omega0^2);
end
if K<(N/2)
    G0=omega0/(omega0+1);a0=(omega0-1)/(omega0+1);
end
w=0:pi/300:pi;
Hw2=1./(1+(tan(w/2)/omega0).^(2*N));
plot(w/pi,Hw2);
grid;

（6）IIR高通濾波（Butterworth）
% IIR Hightpass Use Butterworth
% copyright by Etual
clear;
fs=20;fpass=5;fstop=4;
Ap=0.5;As=10;
wp=2*pi*fpass/fs;ws=2*pi*fstop/fs;
omegap=cot(wp/2);omegas=cot(ws/2);
ep=sqrt(10^(Ap/10)-1);
es=sqrt(10^(As/10)-1);
N=ceil(log(es/ep)/log(omegas/omegap));
omega0=omegap/ep^(1/N);
K=floor(N/2);
for i=1:K
    theta(i)=pi*(N-1+2*i)/(2*N);
end
for i=1:K
    G(i)=omega0^2/(1-2*omega0*cos(theta(i))+omega0^2);
end
for i=1:K
    a1(i)=-2*(omega0^2-1)/(1-2*omega0*cos(theta(i))+omega0^2);
end
for i=1:K
    a2(i)=(1+2*omega0*cos(theta(i))+omega0^2)/(1-2*omega0*cos(theta(i))+omega0^2);
end
if K<(N/2)
    G0=omega0/(omega0+1);a0=-(omega0-1)/(omega0+1);
end
w=(0+eps):pi/300:pi;
Hw2=1./(1+(cot(w/2)/omega0).^(2*N));
plot(w/pi,Hw2);
grid;

（7）IIR帶通濾波（Butterworth）
% IIR Bandpass Use Butterworth
% copyright by Etual
clear;
fs=20;fpa=2;fpb=4;fsa=1.5;fsb=4.5;
Ap=0.0877;As=16.9897;
wpa=2*pi*fpa/fs;wpb=2*pi*fpb/fs;wsa=2*pi*fsa/fs;wsb=2*pi*fsb/fs;
c=sin(wpa+wpb)/(sin(wpa)+sin(wpb));
omegap=abs((c-cos(wpb))/sin(wpb));
omegasa=(c-cos(wsa))/sin(wsa);omegasb=(c-cos(wsb))/sin(wsb);
omegas=min(abs(omegasa),abs(omegasb));
ep=sqrt(10^(Ap/10)-1);es=sqrt(10^(As/10)-1);
N=ceil(log(es/ep)/log(omegas/omegap));
omega0=omegap/ep^(1/N);
K=floor(N/2);
for i=1:K
    theta(i)=pi*(N-1+2*i)/(2*N);
end
for i=1:K
    G(i)=omega0^2/(1-2*omega0*cos(theta(i))+omega0^2);
end
for i=1:K
    a1(i)=4*c*(omega0*cos(theta(i))-1)/(1-2*omega0*cos(theta(i))+omega0^2);
end
for i=1:K
    a2(i)=2*(2*c^2+1-omega0^2)/(1-2*omega0*cos(theta(i))+omega0^2);
end
for i=1:K
    a3(i)=-(4*c*(omega0*cos(theta(i))+1))/(1-2*omega0*cos(theta(i))+omega0^2);
end
for i=1:K
    a4(i)=(1+2*omega0*cos(theta(i))+omega0^2)/(1-2*omega0*cos(theta(i))+omega0^2);
end
if K<(N/2)
    G0=omega0/(1+omega0);a0(1)=-2*c/(1+omega0);a0(2)=(1-omega0)/(1+omega0);
end
w=(0+eps):pi/300:pi;
Hw2=1./(1+((c-cos(w))./(omega0*sin(w))).^(2*N));
plot(w/pi,Hw2);
grid;

（8）IIR帶阻濾波（Butterworth）
% IIR Bandstop Use Butterworth
% copyright by Etual
clear;
fs=20;fpa=1.5;fpb=4.5;fsa=2;fsb=4;
Ap=0.5;As=10;
wpa=2*pi*fpa/fs;wpb=2*pi*fpb/fs;wsa=2*pi*fsa/fs;wsb=2*pi*fsb/fs;
c=sin(wpa+wpb)/(sin(wpa)+sin(wpb));
omegap=abs(sin(wpb)/(c-cos(wpb)));
omegasa=sin(wsa)/(cos(wsa)-c);omegasb=sin(wsb)/(cos(wsb)-c);
omegas=min(abs(omegasa),abs(omegasb));
ep=sqrt(10^(Ap/10)-1);es=sqrt(10^(As/10)-1);
N=ceil(log(es/ep)/log(omegas/omegap));
omega0=omegap/ep^(1/N);
K=floor(N/2);
theta=zeros(1,K);
for i=1:K
    theta(i)=pi*(N-1+2*i)/(2*N);
end
G=zeros(1,K);a1=zeros(1,K);a2=zeros(1,K);
for i=1:K
    G(i)=omega0^2/(1-2*omega0*cos(theta(i))+omega0^2);
end
for i=1:K
    a1(i)=2*(omega0^2-1)/(1-2*omega0*cos(theta(i))+omega0^2);
end
for i=1:K
    a2(i)=(1+2*omega0*cos(theta(i))+omega0^2)/(1-2*omega0*cos(theta(i))+omega0^2);
end
if K<(N/2)
    G0=omega0/(omega0+1);a0=(omega0-1)/(omega0+1);
end
w=(0+eps):pi/300:pi;
Hw2=1./(1+(sin(w)./(omega0*(c-cos(w)))).^(2*N));
plot(w/pi,Hw2);
grid;

（9）IIR低通濾波（chebyshev 1）
% IIR Lowpass Use Chebyshev Type 1
% copyright by Etual
clear;
fs=20;fpass=4;fstop=5;
Ap=0.5;As=10;
wp=2*pi*fpass/fs;ws=2*pi*fstop/fs;
omegap=tan(wp/2);omegas=tan(ws/2);
ep=sqrt(10^(Ap/10)-1);
es=sqrt(10^(As/10)-1);
e=es/ep;w=omegas/omegap;
N=ceil(log(e+sqrt(e^2-1))/log(w+sqrt(w^2-1)));
a=log(1/ep+sqrt(1/ep^2+1))/N;
omega0=omegap*sinh(a);
K=floor(N/2);
theta=zeros(1,K);omega=zeros(1,K);
for i=1:K
    theta(i)=pi*(N-1+2*i)/(2*N);
end
for i=1:K
    omega(i)=omegap*sin(theta(i));
end
G=zeros(1,K);a1=zeros(1,K);a2=zeros(1,K);
for i=1:K
    G(i)=(omega0^2+omega(i)^2)/(1-2*omega0*cos(theta(i))+omega0^2+omega(i)^2);
end
for i=1:K
    a1(i)=2*(omega0^2+omega(i)^2-1)/(1-2*omega0*cos(theta(i))+omega0^2+omega(i)^2);
end
for i=1:K
    a2(i)=(1+2*omega0*cos(theta(i))+omega0^2+omega(i)^2)/(1-2*omega0*cos(theta(i))+omega0^2+omega(i)^2);
end
if K<(N/2)
    G0=omega0/(omega0+1);a0=(omega0-1)/(omega0+1); 
else
    H0=sqrt(1/(1+ep^2));
end
f=0:1/300:10;
Hf2=1./(1+ep^2*(cheby(N,tan(pi*f/fs)/omegap)).^2);
plot(f,abs(Hf2));
grid;
 （9）IIR低通濾波（chebyshev 1）
% IIR Lowpass Use Chebyshev Type 2
% copyright by Etual
clear;
fs=20;fpass=4;fstop=5;
Ap=0.5;As=10;
wp=2*pi*fpass/fs;ws=2*pi*fstop/fs;
omegap=tan(wp/2);omegas=tan(ws/2);
ep=sqrt(10^(Ap/10)-1);
es=sqrt(10^(As/10)-1);
e=es/ep;w=omegas/omegap;
N=ceil(log(e+sqrt(e^2-1))/log(w+sqrt(w^2-1)));
a=log(es+sqrt(es^2+1))/N;
omega0=omegas/sinh(a);
K=floor(N/2);
for i=1:K
    theta(i)=pi*(N-1+2*i)/(2*N);
end
for i=1:K
    omega(i)=omegas/sin(theta(i));
end
for i=1:K
    G(i)=(1+omega(i)^-2)/(1-2*omega0^-1*cos(theta(i))+omega0^-2+omega(i)^-2);
end
for i=1:K
    a1(i)=2*(1-omega0^-2+omega(i)^-2)/(1-2*omega0^-1*cos(theta(i))+omega0^-2+omega(i)^-2);
end
for i=1:K
    a2(i)=(1+2*omega0^-1*cos(theta(i))+omega0^-2+omega(i)^-2)/(1-2*omega0^-1*cos(theta(i))+omega0^-2+omega(i)^-2);
end
for i=1:K
    b1(i)=2*(1-omega(i))/(1+omega(i));
end
if K<(N/2)
    G0=omega0/(omega0+1);a0=(omega0-1)/(omega0+1); 
else
    H0=sqrt(1/(1+ep^2));
end
f=(0+eps):1/100:10;
Hf2=(cheby(N,omegas./tan(pi*f/fs))).^2./((cheby(N,omegas./tan(pi*f/fs))).^2+es^2);
plot(f,abs(Hf2));
grid;

（10）chebyshev 中用到的函數(shù) cheby.m
function CN=cheby(N,x)
if x<=1
    CN=cos(N*acos(x));
else
    CN=cosh(N*log(x+sqrt(x.^2-1)));
end
 

作者：xielinghui

本站聲明： 本文章由作者或相關(guān)機(jī)構(gòu)授權(quán)發(fā)布，目的在于傳遞更多信息，并不代表本站贊同其觀點(diǎn)，本站亦不保證或承諾內(nèi)容真實(shí)性等。需要轉(zhuǎn)載請聯(lián)系該專欄作者，如若文章內(nèi)容侵犯您的權(quán)益，請及時(shí)聯(lián)系本站刪除。
換一批
延伸閱讀

9月2日消息，不造車的華為或?qū)⒋呱龈蟮莫?dú)角獸公司，隨著阿維塔和賽力斯的入局，華為引望愈發(fā)顯得引人矚目。

關(guān)鍵字： 阿維塔 塞力斯 華為

加利福尼亞州圣克拉拉縣2024年8月30日 /美通社/ -- 數(shù)字化轉(zhuǎn)型技術(shù)解決方案公司Trianz今天宣布，該公司與Amazon Web Services （AWS）簽訂了...

關(guān)鍵字： AWS AN BSP 數(shù)字化

倫敦2024年8月29日 /美通社/ -- 英國汽車技術(shù)公司SODA.Auto推出其旗艦產(chǎn)品SODA V，這是全球首款涵蓋汽車工程師從創(chuàng)意到認(rèn)證的所有需求的工具，可用于創(chuàng)建軟件定義汽車。 SODA V工具的開發(fā)耗時(shí)1.5...

關(guān)鍵字： 汽車 人工智能 智能驅(qū)動 BSP

北京2024年8月28日 /美通社/ -- 越來越多用戶希望企業(yè)業(yè)務(wù)能7×24不間斷運(yùn)行，同時(shí)企業(yè)卻面臨越來越多業(yè)務(wù)中斷的風(fēng)險(xiǎn)，如企業(yè)系統(tǒng)復(fù)雜性的增加，頻繁的功能更新和發(fā)布等。如何確保業(yè)務(wù)連續(xù)性，提升韌性，成...

關(guān)鍵字： 亞馬遜 解密 控制平面 BSP

8月30日消息，據(jù)媒體報(bào)道，騰訊和網(wǎng)易近期正在縮減他們對日本游戲市場的投資。

關(guān)鍵字： 騰訊 編碼器 CPU

8月28日消息，今天上午，2024中國國際大數(shù)據(jù)產(chǎn)業(yè)博覽會開幕式在貴陽舉行，華為董事、質(zhì)量流程IT總裁陶景文發(fā)表了演講。

關(guān)鍵字： 華為 12nm EDA 半導(dǎo)體

8月28日消息，在2024中國國際大數(shù)據(jù)產(chǎn)業(yè)博覽會上，華為常務(wù)董事、華為云CEO張平安發(fā)表演講稱，數(shù)字世界的話語權(quán)最終是由生態(tài)的繁榮決定的。

關(guān)鍵字： 華為 12nm 手機(jī) 衛(wèi)星通信

要點(diǎn)： 有效應(yīng)對環(huán)境變化，經(jīng)營業(yè)績穩(wěn)中有升 落實(shí)提質(zhì)增效舉措，毛利潤率延續(xù)升勢 戰(zhàn)略布局成效顯著，戰(zhàn)新業(yè)務(wù)引領(lǐng)增長 以科技創(chuàng)新為引領(lǐng)，提升企業(yè)核心競爭力 堅(jiān)持高質(zhì)量發(fā)展策略，塑強(qiáng)核心競爭優(yōu)勢...

關(guān)鍵字： 通信 BSP 電信運(yùn)營商 數(shù)字經(jīng)濟(jì)

北京2024年8月27日 /美通社/ -- 8月21日，由中央廣播電視總臺與中國電影電視技術(shù)學(xué)會聯(lián)合牽頭組建的NVI技術(shù)創(chuàng)新聯(lián)盟在BIRTV2024超高清全產(chǎn)業(yè)鏈發(fā)展研討會上宣布正式成立。 活動現(xiàn)場 NVI技術(shù)創(chuàng)新聯(lián)...

關(guān)鍵字： VI 傳輸協(xié)議 音頻 BSP

北京2024年8月27日 /美通社/ -- 在8月23日舉辦的2024年長三角生態(tài)綠色一體化發(fā)展示范區(qū)聯(lián)合招商會上，軟通動力信息技術(shù)（集團(tuán)）股份有限公司（以下簡稱"軟通動力"）與長三角投資（上海）有限...

關(guān)鍵字： BSP 信息技術(shù)
關(guān)閉
關(guān)閉