带阻滤波器matlab代码

合集下载

滤波器代码

滤波器代码
hd[i]=hd[i]*(0.42-0.5*cos(2*PI*i/(N-1))+0.08*cos(4*PI*i/(N-1)));
printf("%lf\t",hd[i]);
}
printf("\n\n\n");
}
}
void HPF() //高通滤波器
{
int N,M,i,k;
double wp,ws,as,wd,Bt;
}
}
void BPF() //带通滤波器
{
int N,M,i,k;
double wpl,wpu,wsl,wsu,as,wd1,wd2,Bt1,Bt2,Bt;
double hd[10000]={0.0};
printf("\n请输入设计指标:\n\n");
printf("*********************************\n\n");
wd=0.5*(wp+ws);
Bt=wp-ws;
printf("以下窗ห้องสมุดไป่ตู้类型可供你选择\n\n");
printf("***********************************************\n");
printf("标号窗口类型阻带最小衰减*\n");
printf(" 1:矩形窗21dB *\n");
}
printf("\n\n\n");
}
else if(k=4) //选用汉明窗
{
N=(int)(8/Bt+0.5);
if(fmod(N,2.)==0)

使用MATLAB设计FIR滤波器

使用MATLAB设计FIR滤波器

使⽤MATLAB设计FIR滤波器1. 采⽤fir1函数设计,fir1函数可以设计低通、带通、⾼通、带阻等多种类型的具有严格线性相位特性的FIR滤波器。

语法形式:b = fir1(n, wn)b = fir1(n, wn, ‘ftype’)b = fir1(n, wn, ‘ftype’, window)b = fir1(n, wn, ‘ftype’, window, ‘noscale’)参数的意义及作⽤:b:返回的FIR滤波器单位脉冲响应,脉冲响应为偶对称,长度为n+1;n:滤波器的介数;wn:滤波器的截⽌频率,取值范围为0<wn<1,1对应信号采样频率⼀半。

如果wn是单个数值,且ftype参数为low,则表⽰设计截⽌频率为wn的低通滤波器,如果ftype参数为high,则表⽰设计截⽌频率为wn的⾼通滤波器;如果wn是有两个数组成的向量[wn1wn2],ftype为stop,则表⽰设计带阻滤波器,ftype为bandpass,则表⽰设计带通滤波器;如果wn是由多个数组成的向量,则根据ftype的值设计多个通带或阻带范围的滤波器,ftype为DC-1,表⽰设计的第⼀个频带为通带,ftype为DC-0,表⽰设计的第⼀个频带为阻带;window:指定使⽤的窗函数,默认为海明窗;noscale:指定是否归⼀化滤波器的幅度。

⽰例:N=41; %滤波器长度fs=2000; %采样频率%各种滤波器的特征频率fc_lpf=200;fc_hpf=200;fp_bandpass=[200 400];fc_stop=[200 400];%以采样频率的⼀半,对频率进⾏归⼀化处理wn_lpf=fc_lpf*2/fs;wn_hpf=fc_hpf*2/fs;wn_bandpass=fp_bandpass*2/fs;wn_stop=fc_stop*2/fs;%采⽤fir1函数设计FIR滤波器b_lpf=fir1(N-1,wn_lpf);b_hpf=fir1(N-1,wn_hpf,'high');b_bandpass=fir1(N-1,wn_bandpass,'bandpass');b_stop=fir1(N-1,wn_stop,'stop');%求滤波器的幅频响应m_lpf=20*log(abs(fft(b_lpf)))/log(10);m_hpf=20*log(abs(fft(b_hpf)))/log(10);m_bandpass=20*log(abs(fft(b_bandpass)))/log(10);m_stop=20*log(abs(fft(b_stop)))/log(10);%设置幅频响应的横坐标单位为Hzx_f=0:(fs/length(m_lpf)):fs/2;%绘制单位脉冲响应%绘制单位脉冲响应subplot(421);stem(b_lpf);xlabel('n');ylabel('h(n)');subplot(423);stem(b_hpf);xlabel('n');ylabel('h(n)');subplot(425);stem(b_bandpass);xlabel('n');ylabel('h(n)');subplot(427);stem(b_stop);xlabel('n');ylabel('h(n)');%绘制幅频响应曲线subplot(422);plot(x_f,m_lpf(1:length(x_f)));xlabel('频率(Hz)','fontsize',8);ylabel('幅度(dB)','fontsize',8);subplot(424);plot(x_f,m_hpf(1:length(x_f)));xlabel('频率(Hz)','fontsize',8);ylabel('幅度(dB)','fontsize',8);subplot(426);plot(x_f,m_bandpass(1:length(x_f)));xlabel('频率(Hz)','fontsize',8);ylabel('幅度(dB)','fontsize',8);subplot(428);plot(x_f,m_stop(1:length(x_f)));xlabel('频率(Hz)','fontsize',8);ylabel('幅度(dB)','fontsize',8);2. 采⽤fir2函数设计,函数算法是:⾸先根据要求的幅频响应向量形式进⾏插值,然后进⾏傅⾥叶变换得到理想滤波器的单位脉冲响应,最后利⽤窗函数对理想滤波器的单位脉冲响应激进型截断处理,由此得到FIR滤波器系数。

matlab滤波器代码

matlab滤波器代码

数字信号处理:已知通带截止频率fp=5kHz,通带最大衰减ap=2dB,阻带截止频率fs=2kHz,阻带最小衰减as=30dB,按照以上技术指标设计巴特沃斯低通滤波器:wp=2*pi*5000;ws=2*pi*12000;Rp=2;As=30;[N,wc]=buttord(wp,ws,Rp,As,'s');[B,A]=butter(N,wc,'s');k=0:511;fk=0:14000/512:14000;wk=2*pi*fk;Hk=freqs(B,A,wk);subplot(2,2,1);plot(fk/1000,20*log10(abs(Hk)));grid onxlabel('频率(kHz)');ylabel('幅度(dB)')axis([0,14,-40,5])切比雪夫1型低通滤波器:wp=2*pi*3000;ws=2*pi*12000;Rp=0.1;As=60;[N1,wpl]=cheb1ord(wp,ws,Rp,As,'s');%cheb1ord,里面的是1,不是L[B1,A1]=cheby1(N1,Rp,wpl,'s');subplot(2,2,1);fk=0:12000/512:12000;wk=2*pi*fk;Hk=freqs(B1,A1,wk);plot(fk/1000,20*log10(abs(Hk)));grid onxlabel('频率(kHz)');ylabel('幅度(dB)')axis([0,12,-70,5])椭圆模拟低通滤波器:wp=2*pi*3000;ws=2*pi*12000;Rp=0.1;As=60;[N,wpo]=ellipord(wp,ws,Rp,As,'s');[B,A]=ellip(N,Rp,As,wpo,'s');subplot(2,2,1);fk=0:12000/512:12000;wk=2*pi*fk;Hk=freqs(B1,A1,wk);plot(fk/1000,20*log10(abs(Hk)));grid onxlabel('频率(kHz)');ylabel('幅度(dB)')axis([0,12,-70,5])p195-14wp=2*4/80;ws=2*20/80;rp=0.5;rs=45;[N,wc]=buttord(wp,ws,rp,rs);[B,A]=butter(N,wc);[hk,wk]=freqz(B,A);fk=wk/pi*40;plot(fk,20*log10(abs(hk)));axis([0,30,-100,0])xlabel('频率(kHZ)');ylabel('幅度(db)');grid on P195-16wp=2*325/2500;ws=2*225/2500;rp=1;rs=40;[N,wc]=ellipord(wp,ws,rp,rs);[B,A]=ellip(N,rp,rs,wc);[hk,wk]=freqz(B,A);fk=wk/pi*40;plot(fk,20*log10(abs(hk)));axis([0,30,-100,0])xlabel('频率(kHZ)');ylabel('幅度(db)');grid onP195-15wp=2*4/80;ws=2*20/80;rp=0.5;rs=45;[N,wc]=cheb1ord(wp,ws,rp,rs);[B,A]=cheby1(N,rp,wc);[hk,wk]=freqz(B,A);fk=wk/pi*40;plot(fk,20*log10(abs(hk)));axis([0,30,-100,0])xlabel('频率(kHZ)');ylabel('幅度(db)');grid on 切比雪夫低通滤波器wp=2*pi*3000;ws=2*pi*12000;rp=0.1;as=60;[N1,wp1]=cheb1ord(wp,ws,rp,as,'s');[B1,A1]=cheby1(N1,rp,wp1,'s');subplot(2,2,1);fk=0:12000/512:12000;wk=2*pi*fk;hk=freqs(B1,A1,wk);plot(fk/1000,20*log10(abs(hk)));grid onxlabel('频率(kHZ)');ylabel('幅度(db)');axis([0,12,-70,5])双音频检测audiofile='test.wav'[in_audio,fs,bits]=wavread(audiofile); [b,a]=cheby1(5,0.1,0.3);out_audio=filter(b,a,in_audio);sound(out_audio,fs,bits);wavwrite(out_audio,fs,bits,'test_out'); xk1=fft(in_audio,512);xk2=fft(out_audio,512);subplot(2,1,1);stem(abs(xk1));subplot(2,1,2);stem(abs(xk2));巴特沃斯模拟高通滤波器。

(完整word版)巴特沃斯带阻数字滤波器设计matlab程序及仿真图 - 副本

(完整word版)巴特沃斯带阻数字滤波器设计matlab程序及仿真图 - 副本

fs=15000;T= 1/fs;rp=1;rs=40;wp1=0.11*pi;wp2=0.81*pi;ws1=0.31*pi;ws2=0.61*pi;%数字带阻滤波器技术指标wc1=(2/T)*tan(wp1/2);%频率预畸变wc2=(2/T)*tan(wp2/2);wr1=(2/T)*tan(ws1/2);wr2=(2/T)*tan(ws2/2);w0=sqrt(wc1*wc2);B=wc2-wc1;wp=1;%归一化通带截止频率ws=wp*(wr1*B) / (w0^2-wr1^2) ; %归一化阻带截止频率[N,wc]=buttord(wp,ws,rp,rs,'s')%求滤波器阶数和3dB截止频率[Z,P,K]=buttap(N)%设计模拟低通滤波器[Md,Nd]=zp2tf(Z,P,K)%将零极点形式转换为传输函数形式[M,N]=lp2bs(Md,Nd,w0,B)%对低通滤波器进行频率变换,转换为带阻滤波器[h,w]=freqs(M,N);%模拟带阻滤波器的幅频响应plot(w/(2*pi),abs(h));grid;xlabel('频率Hz');ylabel('幅度');title('模拟带阻滤波器');[b,a]=bilinear(M,N,15000)%对模拟滤波器双线性变换figure(1);freqz(b,a);[H,W]=freqz(b,a); %绘出频率响应;axis([0,1,-100,20]);figure(2);plot(W*fs/(2*pi),abs(H));grid on;xlabel('频率/Hz');ylabel('幅值');n=0:199;t=n/fs;x=sin(2*pi*400*t)+3*sin(2*pi*3000*t)+2*sin(2*pi*5000*t);figure(3);subplot(311);plot(t,x);axis([0,0.01,-5,5]);title('输入信号');grid on;y=filter(b,a,x);subplot(312);stem(y,'.');title('输出序列');grid on;ya=y*sinc(fs*(ones(length(n),1)*t-(n/fs)'*ones(1,length(t))));subplot(313);plot(t,ya);axis([0,0.01,-3,3]);title('输出波形');grid on;t=(0:100)/fs;figure(4)fs=1.5*10000;n=(0:100)/fs;f=sin(2*pi*400*t)+3*sin(2*pi*3000*t)+2*sin(2*pi*5000*t);y=fftfilt(b,x);[H1,f1]=freqz(f,[1]);[H2,f2]=freqz(y,[1]);f1=f1/pi*fs/2;f2=f2/pi*fs/2;subplot(2,1,1);plot(f1,abs(H1));title('输入信号的频谱');subplot(2,1,2);plot(f2,abs(H2));title('输出信号的频谱');基于Matlab 的带阻滤波器设计.10.20.30.40.50.60.70.80.91-800-600-400-2000N o r m a l i z e d Fre q u⨯π r a d /s a m p l e Ph a se(d e g r e e s )00.10.20.30.40.50.60.70.80.91-100-50N o r m a l i z e d Fr e q u⨯π r a d /s a m p l e M a g n i tu d e1000200030004000500060007000800000.20.40.60.811.21.4频率/Hz幅值00.0010.0020.0030.0040.0050.0060.0070.0080.0090.01-505输入信号020406080100120140160180200-22输出序列0.0010.0020.0030.0040.0050.0060.0070.0080.0090.01-202输出波形01000200030004000500060007000800050100150200输入信号的频谱010002000300040005000600070008000102030输出信号的频谱N =4wc =1.7947b =0.0186 -0.0410 0.1082 -0.1355 0.1810 -0.1355 0.1082 -0.0410 0.0186a =1.0000 -0.6707 -1.3750 0.5678 1.1964 -0.2996 -0.4631 0.0496 0.0762>。

用MATLAB设计低通带通高通和带阻FIR数字滤波器

用MATLAB设计低通带通高通和带阻FIR数字滤波器

抽样频率为f 为1KHZ的数字带通滤波器,性能要求为:通带范围从200HZ到250HZ,在此两频率处衰减不大于3dB,在100HZ和400HZ频率处频率衰减不得小于20dB,采用巴特沃斯滤波器4阶Butterworth滤波器源代码n = 2;Wn = [200 250]/500;[b,a] = butter(n,Wn);freqz(b,a,512,1000);这个滤波器100Hz、400Hz处衰减在-30db以上。

可以通过增大n增加衰减。

用MATLAB设计低通带通高通和带阻FIR数字滤波器(1)低通滤波器的技术指标:H(ejw)=1,0<w<0.3pi;H(ejw)=0,0.3pi<w<pi.(2)带通滤波器的技术指标:H(ejw=0,0<w<0.4pi;H(ejw)=1;0.4pi<w<0.6pi;H(ejw)=0,0.6<w<pi;(3)高通滤波器的技术指标:H(ejw)=0,0<w<0.7pi;H(ejw)=1,0.7pi<w<pi.(4)带阻滤波器的技术指标:H(ejw)=1,0<w<0.3pi;H(ejw)=0;0.3pi<w<0.7pi;H(ejw)=1,0.7<w<pi;低通采样定理实验1.1 实验目的1.了解数字信号处理系统的一般构成;2.掌握奈奎斯特抽样定理。

1.2 实验仪器1.YBLD智能综合信号源测试仪1台2.双踪示波器1台3.MCOM-TG305数字信号处理与现代通信技术实验箱1台4.PC机(装有MATLAB、MCOM-TG305配套实验软件)1台1.3 实验原理一个典型的DSP系统除了数字信号处理部分外,还包括A/D和D/A两部分。

这是因为自然界的信号,如声音、图像等大多是模拟信号,因此需要将其数字化后进行数字信号处理,模拟信号的数字化即称为A/D转换。

matlab滤波器代码要点

matlab滤波器代码要点

数字信号处理:已知通带截止频率fp=5kHz,通带最大衰减ap=2dB,阻带截止频率fs=2kHz,阻带最小衰减as=30dB,按照以上技术指标设计巴特沃斯低通滤波器:wp=2*pi*5000;ws=2*pi*12000;Rp=2;As=30;[N,wc]=buttord(wp,ws,Rp,As,'s');[B,A]=butter(N,wc,'s');k=0:511;fk=0:14000/512:14000;wk=2*pi*fk;Hk=freqs(B,A,wk);subplot(2,2,1);plot(fk/1000,20*log10(abs(Hk)));grid onxlabel('频率(kHz)');ylabel('幅度(dB)')axis([0,14,-40,5])切比雪夫1型低通滤波器:wp=2*pi*3000;ws=2*pi*12000;Rp=0.1;As=60;[N1,wpl]=cheb1ord(wp,ws,Rp,As,'s'); %cheb1ord,里面的是1,不是L[B1,A1]=cheby1(N1,Rp,wpl,'s');subplot(2,2,1);fk=0:12000/512:12000;wk=2*pi*fk;Hk=freqs(B1,A1,wk);plot(fk/1000,20*log10(abs(Hk)));grid onxlabel('频率(kHz)');ylabel('幅度(dB)')axis([0,12,-70,5])椭圆模拟低通滤波器:wp=2*pi*3000;ws=2*pi*12000;Rp=0.1;As=60;[N,wpo]=ellipord(wp,ws,Rp,As,'s');[B,A]=ellip(N,Rp,As,wpo,'s');subplot(2,2,1);fk=0:12000/512:12000;wk=2*pi*fk;Hk=freqs(B1,A1,wk);plot(fk/1000,20*log10(abs(Hk)));grid onxlabel('频率(kHz)');ylabel('幅度(dB)')axis([0,12,-70,5])p195-14wp=2*4/80;ws=2*20/80;rp=0.5;rs=45;[N,wc]=buttord(wp,ws,rp,rs);[B,A]=butter(N,wc);[hk,wk]=freqz(B,A);fk=wk/pi*40;plot(fk,20*log10(abs(hk)));axis([0,30,-100,0])xlabel('频率(kHZ)');ylabel('幅度(db)');grid on P195-16wp=2*325/2500;ws=2*225/2500;rp=1;rs=40;[N,wc]=ellipord(wp,ws,rp,rs);[B,A]=ellip(N,rp,rs,wc);[hk,wk]=freqz(B,A);fk=wk/pi*40;plot(fk,20*log10(abs(hk)));axis([0,30,-100,0])xlabel('频率(kHZ)');ylabel('幅度(db)');grid onP195-15wp=2*4/80;ws=2*20/80;rp=0.5;rs=45;[N,wc]=cheb1ord(wp,ws,rp,rs);[B,A]=cheby1(N,rp,wc);[hk,wk]=freqz(B,A);fk=wk/pi*40;plot(fk,20*log10(abs(hk)));axis([0,30,-100,0])xlabel('频率(kHZ)');ylabel('幅度(db)');grid on 切比雪夫低通滤波器wp=2*pi*3000;ws=2*pi*12000;rp=0.1;as=60;[N1,wp1]=cheb1ord(wp,ws,rp,as,'s');[B1,A1]=cheby1(N1,rp,wp1,'s');subplot(2,2,1);fk=0:12000/512:12000;wk=2*pi*fk;hk=freqs(B1,A1,wk);plot(fk/1000,20*log10(abs(hk)));grid onxlabel('频率(kHZ)');ylabel('幅度(db)');axis([0,12,-70,5])双音频检测audiofile='test.wav'[in_audio,fs,bits]=wavread(audiofile); [b,a]=cheby1(5,0.1,0.3);out_audio=filter(b,a,in_audio);sound(out_audio,fs,bits);wavwrite(out_audio,fs,bits,'test_out'); xk1=fft(in_audio,512);xk2=fft(out_audio,512);subplot(2,1,1);stem(abs(xk1));subplot(2,1,2);stem(abs(xk2));巴特沃斯模拟高通滤波器。

数字滤波器的设计及其MATLAB实现

数字滤波器的设计及其MATLAB实现

设计低通数字滤波器,要求在通带内频率低于0.2pirad时,允许幅度误差在1dB以内,在频率0.3pi rad~pi rad之间的阻带衰减大于15dB,用脉冲响应不变法设计数字滤波器,T=1: 切比雪夫I型模拟滤波器的设计子程序:function [b,a]=afd_chb1(Omegap,Omegar,Ar)if Omegap<=0error('通带边缘必须大于0')endif(Dt<=0)|(Ar<0)error('通带波动或阻带衰减必须大于0');endep=sqrt(10^(Dt/10)-1);A=10^(Ar/20);OmegaC=Omegap;OmegaR=Omegar/Omegap;g=sqrt(A*A-1)/ep;N=ceil(log10(g+sqrt(g*g-1))/log10(OmegaR+sqrt(OmegaR*OmegaR-1)));fprintf('\n***切比雪夫I型模拟低通滤波器阶数=%2.0f\n',N);[b,a]=u_chblap(N,Dt,OmegaC);设计非归一化切比雪夫I型模拟低通滤波器原型程序:function [b,a]=u_chblap(N,Dt,OmegaC)[z,p,k]=cheb1ap(N,Dt);a=real(poly(p));aNn=a(N+1);p=p*OmegaC;a=real(poly(p));aNu=a(N+1);k=k*aNu/aNn;b0=k;B=real(poly(z));b=k*B;直接形式转换成级联形式子程序:function [C,B,A]=sdir2cas(b,a)Na=length(a)-1;Nb=length(b)-1;b0=b(1);b=b/b0;a0=a(1);a=a/a0;C=b0/a0;p=cplxpair(roots(a));K=floor(Na/2);if K*2==NaA=zeros(K,3);for n=1:2:NaArow=p(n:1:n+1,:);Arow=poly(Arow);A((fix(n+1)/2),:)=real(Arow);elseif Na==1A=[0 real(poly(p))];elseA=zeros(K+1,3);for n=1:2:2*KArow=p(n:1:n+1,:);Arow=poly(Arow);A((fix(n+1)/2),:)=real(Arow);endA(K+1,:)=[0 real(poly(p(Na)))];endz=cplxpair(roots(b));K=floor(Nb/2);if Nb==0B=[0 0 poly(z)];elseif K*2==NbB=zeros(K,3);for n=1:2:NbBrow=z(n:1:n+1,:);Brow=poly(Brow);B((fix(n+1)/2),:)=real(Brow);endelseif Nb==1B=[0 real(poly(z))];elseB=zeros(K+1,3);for n=1:2:2*KBrow=z(n:1:n+1,:);Brow=poly(Brow);B((fix(n+1)/2),:)=real(Brow);endB(K+1,:)=[0 real(poly(z(Nb)))];End计算系统函数的幅度响应和相位响应子程序:function [db,mag,pha,w]=freqs_m(b,a,wmax)w1=0:500;w=w1*wmax/500;h=freqs(b,a,w);mag=abs(h);db=20*log10((mag+eps)/max(mag));pha=angle(h);脉冲响应不变法程序:function [b,a]=imp_invr(c,d,T)[R,p,k]=residue(c,d);p=exp(p*T);[b,a]=residuez(R,p,k);b=real(b).*T;数字滤波器响应子程序:function [db,mag,pha,grd,w]=freqz_m(b,a);[H,w]=freqz(b,a,1000,'whole');H=(H(1:501))';w=(w(1:501))';mag=abs(H);db=20*log10((mag+eps)/max(mag));pha=angle(H);grd=grpdelay(b,a,w);直接转换成并联型子程序:function [C,B,A]=dir2par(b,a)M=length(b);N=length(a);[r1,p1,C]=residuez(b,a);p=cplxpair(p1,10000000*eps);x=cplxcomp(p1,p);r=r1(x);K=floor(N/2);B=zeros(K,2);A=zeros(K,3);if K*2==Nfor i=1:2:N-2br=r(i:1:i+1,:);ar=p(i:1:i+1,:);[br,ar]=residuez(br,ar,[]);B((fix(i+1)/2),:)real(br');A((fix(i+1)/2),:)real(ar');end[br,ar]=residuez(r(N-1),p(N-1),[]);B(K,:)=[real(br') 0];A(K,:)=[real(ar') 0];elsefor i=1:2:N-1br=r(i:1:i+1,:);ar=p(i:1:i+1,:);[br,ar]=residuez(br,ar,[]);B((fix(i+1)/2),:)real(br);A((fix(i+1)/2),:)real(ar);endEnd比较两个含同样标量元素但(可能)有不同下标的复数对及其相位留数向量子程序:function I=cplxcomp(p1,p2)I=[];for i=1:length(p2)for j=1:length(p1)if(abs(p1(j)-p2(i))<0.0001)I=[I,j];endendendI=I';双线性变换巴特沃斯低通滤波器设计:巴特沃思模拟滤波器的设计子程序:function [b,a]=afd_butt(wp,ws,Rp,rs)if wp<=0error('通带边缘必须大于0');endif ws<=wperror('阻带边缘必须大于通带边缘');endif(Rp<=0)|(Rs<0)error('通带波动或阻带衰减必须大于0');endN=ceil((log10((10^(Rp/10)-1)/(10^(Rs/10)-1)))/(2*log10(wp/ws))); fprintf('\n***Butterworth Filter Order=%2.0f\n',N);OmegaC=wp/((10^(Rp/10)-1)^(1/(2*N)));[b,a]=u_buttap(N,OmegaC)设计非归一化巴特沃思模拟低通滤波器原型子程序:function [b,a]=u_buttap(N,OmegaC)[z,p,k]=buttap(N);p=p*OmegaC;k=k*OmegaC^N;B=real(poly(z));b0=k;b=k*B;a=real(poly(p));直接型到级联型形式的转换:function [b0,B,A]=dir2cas(b,a)b0=b(1);b=b/b0;a0=a(1);a=a/a0;b0=b0/a0;M=length(b);N=length(a);if N>Mb=[b,zeros(1,N-M)];a=[a,zeros(1,M-N)];elseNM=0;endk=floor(N/2);B=zeros(k,3);A=zeros(k,3);if k*2==Nb=[b,0];a=[a,0];endbroots=cplxpair(roots(b));aroots=cplxpair(roots(a));for i=1:2:2*kbr=broots(i:1:i+1,:);br=real(polt(br));B((fix(i+1)/2),:)=br;ar=aroots(i:1:i+1,:);ar=real(polt(ar));A((fix(i+1)/2),:)=ar;Endfunction [db,mag,pha,grd,w]=freqz_m(b,a)[h,w]=freqz(b,a,1000,'whole');h=(h(1:501))';w=(w(1:501))';mag=abs(h);db=20*log10((mag+eps)/max(mag));pha=angle(h);grd=grdelay(b,a,w);设计一个巴特沃思高通滤波器,要求通带截止频率为0.6pi,通带内衰减不大于1dB,阻带·起始频率为0.4pi,阻带内衰减不小于15dB,T=1:>> wp=0.6*pi;ws=0.4*pi;>> Rp=1;Rs=15;T=1;>> [N,wn]=buttord(wp/pi,ws/pi,Rp,Rs) 计算巴特沃思滤波器阶数和截止频率N =4wn =>> [b,a]=butter(N,wn,'high'); 频率变换法计算巴特沃思高通滤波器>> [C,B,A]=dir2cas(b,a)C =0.0751B =1.0000 -2.0000 1.00001.0000 -2.0000 1.0000A =1.0000 0.1562 0.44881.0000 0.1124 0.0425>> [db,mag,pha,grd,w]=freqz_m(b,a);>> subplot(2,1,1);plot(w/pi,mag);>> subplot(2,1,2);plot(w/pi,db);椭圆带通滤波器的设计--ellip函数的应用:>> ws=[0.3*pi 0.75*pi]; 数字阻带边缘频率>> wp=[0.4*pi 0.6*pi]; 数字通带边缘频率>> Rp=1;Rs=40;>> Ripple=10^(-Rp/20); 通带波动>> Attn=10^(-Rs/20); 阻带衰减>> [N,wn]=ellipord(wp/pi,ws/pi,Rp,Rs) 计算椭圆滤波器参数N =4wn =0.4000 0.6000>> [b,a]=ellip(N,Rp,Rs,wn); 数字椭圆滤波器的设计>> [b0,B,A]=dir2cas(b,a) 级联形式实现b0 =0.0197B =1.0000 1.5066 1.00001.0000 0.9268 1.00001.0000 -0.9268 1.00001.0000 -1.5066 1.0000A =1.0000 0.5963 0.93991.0000 0.2774 0.79291.0000 -0.2774 0.79291.0000 -0.5963 0.9399>> figure(1);>> [db,mag,pha,grd,w]=freqz_m(b,a);>> subplot(2,2,1);plot(w/pi,mag);>> grid on;>> subplot(2,2,3);plot(w/pi,db);grid on;>> subplot(2,2,2);plot(w/pi,pha/pi);grid on;>> subplot(2,2,4);plot(w/pi,grd);设计一个巴特沃思带阻滤波器,要求通带上下截止频率为0.8pi、0.2pi,通带内衰减不大于1dB,阻带上起始频率为0.7pi、0.4pi,阻带内衰减不小于30dB:>> wp=[0.2*pi 0.8*pi];>> ws=[0.4*pi 0.7*pi];>> Rp=1;Rs=30;>> [N,wn]=buttord(wp/pi,ws/pi,Rp,Rs);>> [b,a]=butter(N,wn,'stop');>> [C,B,A]=dir2cas(b,a)C =0.0394B =1.0000 0.3559 0.99941.0000 0.3547 1.00401.0000 0.3522 0.99541.0000 0.3499 1.00461.0000 0.3475 0.99601.0000 0.3463 1.0006A =1.0000 1.3568 0.79281.0000 1.0330 0.46331.0000 0.6180 0.17751.0000 -0.2493 0.11131.0000 -0.6617 0.37551.0000 -0.9782 0.7446>> [db,mag,pha,grd,w]=freqz_m(b,a); >> subplot(2,1,1);plot(w/pi,mag);>> subplot(2,1,2);plot(w/pi);数字低通---数字带阻:function [bz,az]=zmapping(bZ,aZ,Nz,Dz) bzord=(length(bZ)-1)*(length(Nz)-1); azord=(length(aZ)-1)*(length(Dz)-1);bz=zeros(1,bzord+1);for k=0:bzordpln=[1];for i=0:k-1pln=conv(pln,Nz);endpld=[1];for i=0:bzord-k-1pld=conv(pld,Dz);endbz=bz+bZ(k+1)*conv(pln,pld); endfor k=0:azordpln=[1];for i=0:k-1pln=conv(pln,Nz);endpld=[1];for i=0:azord-k-1pld=conv(pld,Dz);endaz=az+aZ(k+1)*conv(pln,pld); endall=az(1);az=az/az1;bz=bz/az1;线性相位FIR滤波器的幅度特性:function pzkplot(num,den)hold on;axis('square');x=-1:0.01:1;y=(1-x.^2).^0.5;y1=-(1-x.^2).^0.5;plot(x,y,'b',x,y1,'b');num1=length(num);den1=length(den);if(num1>1)z=roots(num);elsez=0;endif(den1>1)p=roots(den);elsep=0;endif(num>1&den1>1)r_max_z=max(abs(real(z)));i_max_z=max(abs(imag(z)));a_max_z=max(r_max_z,i_max_z);r_max_p=max(abs(real(p)));i_max_p=max(abs(imag(p)));a_max_p=max(r_max_p,i_max_p);a_max=max(a_max_z,a_max_p);elseif (num1>1)r_max_z=max(abs(real(z)));i_max_z=max(abs(imag(z)));a_max=max(r_max_z,i_max_z);elser_max_p=max(abs(real(p)));i_max_p=max(abs(imag(p)));a_max=max(r_max_p,i_max_p);endaxis([-a_max a_max -a_max a_max]);plot([-a_max a_max],[0 0],'b');plot([0 0],[-a_max a_max],'b');plot([-a_max a_max],[a_max a_max],'b');plot([a_max a_max],[-a_max a_max],'b');Lz=length(z);for i=1:Lz;plot(real(z(i)),imag(z(i)),'bo');endLp=length(p);for j=1:Lpplot(real(p(j)),imag(p(j)),'bx');endtitle('The zeros-pole plot');xlabel('虚部');ylabel('实部');function [Hr,w,a,L]=Hr_Type1(h)M=length(h);L=(M-1)/2;a=[h(L+1) 2*h(L:-1:1)];n=[0:1:L];w=[0:1:500]'*pi/500;Hr=cos(w*n)*a';设计I型线性相位FIR滤波器:>> h=[-4 1 -1 -2 5 6 5 -2 -1 1 -4];>> M=length(h);n=0:M-1;>> [Hr,w,a,L]=Hr_Type1(h);>> amax=max(a)+1;>> amin=min(a)-1;>> subplot(2,2,1);stem(n,h);>> axis([-1 2*L+1 amin amax]);text(2*L+1.5,amin,'n'); >> xlabel('n');ylabel('h(n)');title('脉冲响应');>> subplot(2,2,3);stem(0:L,a);>> axis([-1 2*L+1 amin amax]);>> xlabel('n');ylabel('a(n)');title('a(n) 系数');>> subplot(2,2,2);plot(w/pi,Hr);>> grid on;text(1.05,-20,'频率pi');>> xlabel('频率');ylabel('Hr');title('I 型振幅响应');>> subplot(2,2,4);pzkplot(h,1);>> title('零极点分布');function [hr,w,b,L]=Hr_Type2(h)M=length(h);L=M/2;b=2*h(L:-1:1);n=[1:1:L];n=n-0.5;w=[0:1:500]'*pi/500;hr=cos(w*n)*b';II型线性相位FIR滤波器:>> h=[-4 1 -1 -2 5 6 5 -2 -1 1 -4];>> M=length(h);n=0:M-1;>> [Hr,w,b,L]=Hr_Type2(h);Warning: Integer operands are required for colon operator when used as index. > In Hr_Type2 at 2>> bmax=max(b)+1;bmin=min(b)-1;>> subplot(2,2,1);stem(n,h);axis([-1 2*L+1 bmin bmax]);text(2*L+1.5,bmin,'n');xlabel('n');ylabel('h(n)');title('脉冲响应');>> subplot(2,2,3);stem(1:L,b);axis([-1 2*L+1 bmin bmax]);xlabel('n');ylabel('b(n)');title('b(n) 系数');>> subplot(2,2,2);plot(w/pi,Hr);grid on;text(1.05,-20,'频率pi');xlabel('频率');ylabel('Hr');title('II 型振幅响应');>> subplot(2,2,4);pzkplot(h,1);title('零极点分布');function [hr,w,c,L]=Hr_Type3(h)M=length(h);L=(M-1)/2;b=2*h(L+1:-1:1);n=[1:1:L];w=[0:1:500]'*pi/500;hr=cos(w*n)*c';用MA TLAB编程绘制各种窗函数的形状。

滤波matlab代码

滤波matlab代码

滤波matlab代码滤波是信号处理中常用的技术,用于去除信号中的噪声或者滤波信号以得到感兴趣的频率成分。

在MATLAB中,有多种滤波函数可以使用,例如`filter`、`designfilt`和`fir1`等。

本文将介绍这些函数的用法和原理,并通过实例说明如何使用MATLAB进行滤波。

我们来介绍一下`filter`函数。

该函数可以用于实现各种滤波器,如低通滤波器、高通滤波器和带通滤波器等。

其基本语法为:```Matlaby = filter(b,a,x)```其中,`b`和`a`是滤波器的系数,`x`是输入信号的向量。

这个函数将输出滤波后的信号`y`。

接下来,我们来看一个实例。

假设我们有一个包含噪声的信号`x`,我们希望通过低通滤波器来去除噪声。

我们可以使用`filter`函数来实现这个功能。

首先,我们需要设计一个低通滤波器的系数。

可以使用`fir1`函数来设计一个FIR滤波器的系数。

例如,我们可以使用以下代码来设计一个阶数为10的低通滤波器:```Matlaborder = 10; % 滤波器阶数cutoff = 0.2; % 截止频率b = fir1(order, cutoff);```然后,我们可以使用这个滤波器对信号进行滤波:```Matlaby = filter(b, 1, x);```这样,我们就得到了滤波后的信号`y`。

除了`filter`函数,MATLAB还提供了`designfilt`函数用于设计各种类型的滤波器。

该函数可以设计IIR滤波器、带通滤波器、带阻滤波器等。

使用`designfilt`函数需要指定滤波器的类型、阶数以及其他参数。

例如,我们可以使用以下代码来设计一个IIR低通滤波器:```Matlaborder = 6; % 滤波器阶数cutoff = 0.2; % 截止频率d = designfilt('lowpassiir', 'FilterOrder', order, 'PassbandFrequency', cutoff);```然后,我们可以使用这个滤波器对信号进行滤波:```Matlaby = filter(d, x);```同样地,我们得到了滤波后的信号`y`。

  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

带阻滤波器(Bandstop Filter)是一种滤波器,它的作用是在频域上阻止特定频率的信号通过,而允许其他频率的信号通过。

在MATLAB中,你可以使用`iirnotch`函数设计带阻滤波器。

下面是一个简单的例子,说明如何使用`iirnotch`函数设计一个带阻滤波器:
1. 设定滤波器的参数。

在这个例子中,我们将设计一个中心频率为1000Hz,半带宽为500Hz的带阻滤波器:
```matlab
% 设定参数
fs = 1000; % 采样频率,单位Hz
fpass = 1000; % 通带频率,单位Hz
fstop = 1500; % 阻带频率1,单位Hz
Bandwidth = 500; % 阻带宽度,单位Hz
```
2. 设计滤波器:
```matlab
% 设计带阻滤波器
[N, Wn] = iirnotch(fstop, fpass, Bandwidth, fs);
```
3. 绘制滤波器的频率响应:
```matlab
% 绘制频率响应
[H, f] = freqz(N, Wn, fs);
plot(f, 20 * log10(abs(H)));
title('Bandstop Filter Frequency Response');
xlabel('Frequency (Hz)');
ylabel('Amplitude (dB)');
```
这个例子中的代码将设计一个带阻滤波器,并绘制其频率响应。

你可以根据你的具体需求调整参数。

相关文档
最新文档