利用汉宁窗设计Ⅰ型数字高通滤波器
clear all; Wp=0.6*pi; Ws=0.4*pi;
tr_width=Wp-Ws; %过渡带宽度 N=ceil(6.2*pi/tr_width) %滤波器长度 n=0:1:N-1;
Wc=(Ws+Wp)/2; %理想低通滤波器的截止频率 hd=ideal_hp1(Wc,N); %理想低通滤波器的单位冲激响应 w_han=(hanning(N))'; %汉宁窗
h=hd.*w_han; %截取得到实际的单位脉冲响应 [db,mag,pha,w]=freqz_m2(h,[1]); %计算实际滤波器的幅度响应 delta_w=2*pi/1000;
Ap=-(min(db(Wp/delta_w+1:1:501))) %实际通带纹波 As=-round(max(db(1:1:Ws/delta_w+1))) %实际阻带纹波 subplot(221) stem(n,hd)
title('理想单位脉冲响应hd(n)') subplot(222) stem(n,w_han) title('汉宁窗w(n)') subplot(223)
精彩文档
实用标准文案
stem(n,h)
title('实际单位脉冲响应hd(n)') subplot(224) plot(w/pi,db) title('幅度响应(dB)') axis([0,1,-100,10]) clear all; Wp=0.6*pi; Ws=0.4*pi;
tr_width=Wp-Ws; %过渡带宽度 N=ceil(6.2*pi/tr_width) %滤波器长度 n=0:1:N-1;
Wc=(Ws+Wp)/2; %理想低通滤波器的截止频率 hd=ideal_hp1(Wc,N); %理想低通滤波器的单位冲激响应 w_han=(hanning(N))'; %汉宁窗
h=hd.*w_han; %截取得到实际的单位脉冲响应 [db,mag,pha,w]=freqz_m2(h,[1]); %计算实际滤波器的幅度响应 delta_w=2*pi/1000;
Ap=-(min(db(Wp/delta_w+1:1:501))) %实际通带纹波 As=-round(max(db(1:1:Ws/delta_w+1))) %实际阻带纹波 subplot(221) stem(n,hd)
精彩文档
实用标准文案
title('理想单位脉冲响应hd(n)') subplot(222) stem(n,w_han) title('汉宁窗w(n)') subplot(223) stem(n,h)
title('实际单位脉冲响应hd(n)') subplot(224) plot(w/pi,db) title('幅度响应(dB)') axis([0,1,-100,10])
基于切比雪夫一致逼近法设计FIR数字低通滤波器
clear all;
f=[0 0.6 0.7 1]; %给定频率轴分点
A=[1 1 0 0]; %给定在这些频率分点上理想的幅频响应 weigh=[1 10]; %给定在这些频率分点上的加权
b=remez(32,f,A,weigh); %设计出切比雪夫最佳一致逼近滤波器 [h,w]=freqz(b,1,256,1); h=abs(h); h=20*log10(h); subplot(211) stem(b,'.');
精彩文档
实用标准文案
grid;
title('切比雪夫逼近滤波器的抽样值') subplot(212) plot(w,h); grid;
title('滤波器幅频特性(dB)')
利用汉宁窗设计Ⅰ型数字带阻滤波器
clear all; Wpl=0.2*pi; Wph=0.8*pi; Wsl=0.4*pi; Wsh=0.6*pi;
tr_width=min((Wsl-Wpl),(Wph-Wsh)); %过渡带宽度 N=ceil(6.2*pi/tr_width) %滤波器长度 n=0:1:N-1;
Wcl=(Wsl+Wpl)/2; %理想低通滤波器的截止频率 Wch=(Wsh+Wph)/2;
hd=ideal_bs(Wcl,Wch,N); %理想低通滤波器的单位冲激响应 w_hann=(hanning(N))'; %汉宁窗
h=hd.*w_hann; %截取得到实际的单位脉冲响应 [db,mag,pha,w]=freqz_m2(h,[1]); %计算实际滤波器的幅度响应 delta_w=2*pi/1000;
精彩文档
实用标准文案
Ap=-(min(db(1:1:Wpl/delta_w+1))) %实际通带纹波 As=-round(max(db(Wsl/delta_w+1:1:Wsh/delta_w+1))) %实际阻带纹波 subplot(221) stem(n,hd)
title('理想单位脉冲响应hd(n)') subplot(222) stem(n,w_hann) title('汉宁窗w(n)') subplot(223) stem(n,h)
title('实际单位脉冲响应hd(n)') subplot(224) plot(w/pi,db) title('幅度响应(dB)') axis([0,1,-100,10])
利用三角窗设计Ⅲ型数字带通滤波器
clear all; Wpl=0.4*pi; Wph=0.6*pi; Wsl=0.2*pi; Wsh=0.8*pi;
tr_width=min((Wpl-Wsl),(Wsh-Wph)); 精彩文档
%过渡带宽度
实用标准文案
N=ceil(6.1*pi/tr_width) %滤波器长度 n=0:1:N-1;
Wcl=(Wsl+Wpl)/2; %理想低通滤波器的截止频率 Wch=(Wsh+Wph)/2;
hd=ideal_bp2(Wcl,Wch,N); %理想低通滤波器的单位冲激响应
w_tri=(triang(N))'; %三角窗
h=hd.*w_tri; %截取得到实际的单位脉冲响应 [db,mag,pha,w]=freqz_m2(h,[1]); %计算实际滤波器的幅度响应
delta_w=2*pi/1000;
Ap=-(min(db(Wpl/delta_w+1:1:Wph/delta_w+1))) %实际通带纹波 As=-round(max(db(Wsh/delta_w+1:1:501))) %实际阻带纹波 subplot(221) stem(n,hd)
title('理想单位脉冲响应hd(n)') subplot(222) stem(n,w_tri) title('三角窗w(n)') subplot(223) stem(n,h)
title('实际单位脉冲响应hd(n)')
精彩文档
实用标准文案
subplot(224) plot(w/pi,db) title('幅度响应(dB)') axis([0,1,-100,10])
利用布拉克曼窗设计Ⅱ型数字带通滤波器
clear all; Wpl=0.4*pi; Wph=0.6*pi; Wsl=0.2*pi; Wsh=0.8*pi;
tr_width=min((Wpl-Wsl),(Wsh-Wph)); %过渡带宽度 N=ceil(11*pi/tr_width)+1 %滤波器长度 n=0:1:N-1;
Wcl=(Wsl+Wpl)/2; %理想低通滤波器的截止频率 Wch=(Wsh+Wph)/2;
hd=ideal_bp1(Wcl,Wch,N); %理想低通滤波器的单位冲激响应
w_bman=(blackman(N))'; %布拉克曼窗
h=hd.*w_bman; %截取得到实际的单位脉冲响应 [db,mag,pha,w]=freqz_m2(h,[1]); %计算实际滤波器的幅度响应
delta_w=2*pi/1000;
精彩文档
实用标准文案
Ap=-(min(db(Wpl/delta_w+1:1:Wph/delta_w+1))) %实际通带纹波 As=-round(max(db(Wsh/delta_w+1:1:501))) %实际阻带纹波 subplot(221) stem(n,hd)
title('理想单位脉冲响应hd(n)') subplot(222) stem(n,w_bman) title('布拉克曼窗w(n)') subplot(223) stem(n,h)
title('实际单位脉冲响应hd(n)') subplot(224) plot(w/pi,db) title('幅度响应(dB)') axis([0,1,-100,10])
利用海明窗设计Ⅱ型数字低通滤波器
clear all; Wp=0.2*pi; Ws=0.4*pi;
tr_width=Ws-Wp; %过渡带宽度 N=ceil(6.6*pi/tr_width)+1 %滤波器长度 n=0:1:N-1;
精彩文档
实用标准文案
Wc=(Ws+Wp)/2; %理想低通滤波器的截止频率 hd=ideal_lp1(Wc,N); %理想低通滤波器的单位冲激响应 w_ham=(hamming(N))'; %海明窗
h=hd.*w_ham; %截取得到实际的单位脉冲响应 [db,mag,pha,w]=freqz_m2(h,[1]); %计算实际滤波器的幅度响应 delta_w=2*pi/1000;
Ap=-(min(db(1:1:Wp/delta_w+1))) %实际通带纹波 As=-round(max(db(Ws/delta_w+1:1:501))) %实际阻带纹波 subplot(221) stem(n,hd)
title('理想单位脉冲响应hd(n)') subplot(222) stem(n,w_ham) title('海明窗w(n)') subplot(223) stem(n,h)
title('实际单位脉冲响应hd(n)') subplot(224) plot(w/pi,db) title('幅度响应(dB)') axis([0,1,-100,10])
%--------------------------------------------------------
精彩文档
实用标准文案
function[db,mag,pha,w]=freqz_m2(b,a) %滤波器的幅值响应(相对、绝对)、相位响应 %db:相对幅值响应 %mag：绝对幅值响应 %pha: 相位响应 %w 采样频率；
%b 系统函数H(z)的分子项(对FIR，b=h) %a 系统函数H(z)的分母项(对FIR，a=1) [H,w]=freqz(b,a,1000,'whole'); H=(H(1:1:501))'; w=(w(1:1:501))';
mag=abs(H); %绝对幅值响应
db=20*log10((mag+eps)/max(mag)); %相对幅值响应 pha=angle(H); %相位响应
利用模拟Butterworth滤波器设计数字低通滤波器
% exa4-8_pulseDF for example4-8
% using Butterworth analog lowpass filter to design digital lowpass filter %利用模拟Butterworth滤波器设计数字低通滤波器 %脉冲响应不变法 wp=0.2*pi; ws=0.3*pi;
精彩文档
实用标准文案
Rp=1; As=15; T=1; %性能指标 Rip=10^(-Rp/20); Atn=10^(-As/20); OmgP=wp*T; OmgS=ws*T;
[N,OmgC]=buttord(OmgP,OmgS,Rp,As,'s');%选取模拟滤波器的阶数 [cs,ds]=butter(N,OmgC,'s'); %设计出所需的模拟低通滤波器 [b,a]=impinvar(cs,ds,T); %应用脉冲响应不变法进行转换 %求得相对、绝对频响及相位、群迟延响应 [db,mag,pha,grd,w]=freqz_m(b,a); %下面绘出各条曲线 subplot(2,2,1); plot(w/pi,mag); title('幅频特性'); xlabel('w(/pi)'); ylabel('|H(jw)|'); axis([0,1,0,1.1]);
set(gca,'XTickMode','manual','XTick',[0 0.2 0.3 0.5 1]); set(gca,'YTickMode','manual','YTick',[0 Atn Rip 1]);
精彩文档
实用标准文案
grid
subplot(2,2,2); plot(w/pi,db); title('幅频特性(dB)'); xlabel('w(/pi)'); ylabel('dB'); axis([0,1,-40,5]);
set(gca,'XTickMode','manual','XTick',[0 0.2 0.3 0.5 1]); set(gca,'YTickMode','manual','YTick',[-40 -As -Rp 0]); grid
subplot(2,2,3); plot(w/pi,pha/pi); title('相频特性'); xlabel('w(/pi)'); ylabel('pha(/pi)'); axis([0,1,-1,1]);
set(gca,'XTickMode','manual','XTick',[0 0.2 0.3 0.5 1]); grid
subplot(2,2,4); plot(w/pi,grd); title('群延迟'); xlabel('w(/pi)');
精彩文档
实用标准文案
ylabel('Sample'); axis([0,1,0,12]);
set(gca,'XTickMode','manual','XTick',[0 0.2 0.3 0.5 1]); grid
function[db,mag,pha,grd,w]=freqz_m(b,a) %滤波器幅值响应(绝对、相对)、相位响应及群延迟
%Usage: [db,mag,pha,grd,w]=freqz_m(b,a) %500点对应[0,pi]
%db 相对幅值响应；mag 绝对幅值响应；pha 相位响应；grd 群延迟响应 %w 采样频率；b 系统函数H(z)的分子项(对FIR，b＝h) %a 系统函数H(z)的分母项(对FIR，a＝1) [H,w]=freqz(b,a,500);%500点的复频响应 mag=abs(H);
db=20*log10((mag+eps)/max(mag)); pha=angle(H); grd=grpdelay(b,a,w);
基于频域抽样法的FIR数字带阻滤波器设计
clear all; N=41; T1=0.598; alpha= (N-1)/2; l=0:N-1; wl= (2*pi/N)*l;
精彩文档
实用标准文案
Hrs=[ones(1,6),T1,zeros(1,7),T1,ones(1,11),T1,zeros(1,7),T1,ones(1,6)]; %理想振幅采样响应 Hdr=[1,1,0,0,1,1]; wdl=[0,0.3,0.3,0.7,0.7,1]; k1=0:floor((N-1)/2); k2=floor((N-1)/2)+1:N-1;
angH=[pi/2-alpha*(2*pi)/N*(k1+0.5),-pi/2+alpha*(2*pi)/N*(N-k2-0.5)]; %相位约束条件
Hdk=Hrs.*exp(j*angH); %构成Hd(k) h1=ifft(Hdk,N); n=0:1:N-1;
h=real(h1.*exp(j*pi*n/N)); %实际单位冲激响应 [db,mag,pha,w]=freqz_m2(h,[1]);
[Hr,ww,a,L]=hr_type3(h); %实际振幅响应 subplot(221)
plot(wl/pi+1/N,Hrs,'.',wdl,Hdr) title('频率样本Hd(k) :N=41') axis([0 1 -0.1 1.2]) subplot(222) stem(l,h)
title('实际单位脉冲响应h(n)') subplot(223)
精彩文档
实用标准文案
plot(ww/pi,Hr,wl/pi+1/N,Hrs,'.') title('实际振幅响应H(w)') axis([0 1 -0.1 1.2]) subplot(224) plot(w/pi,db) title('幅度响应(dB)') axis([0 1 -80 10])
function [db,mag,pha,w] = freqz_m(b,a); %滤波器的幅值响应(相对、绝对)、相位响应 %db:相对幅值响应 %mag：绝对幅值响应 %pha: 相位响应 %w 采样频率；
%b 系统函数H(z)的分子项(对FIR，b=h) %a 系统函数H(z)的分母项(对FIR，a=1) [H,w] = freqz(b,a,1000,'whole'); H = (H(1:1:501))'; w = (w(1:1:501))'; mag = abs(H);
db = 20*log10((mag+eps)/max(mag)); pha = angle(H);
% pha = unwrap(angle(H));
function [Hr,w,c,L]=hr_type3(h);
精彩文档
实用标准文案
%计算所设计的3型滤波器的振幅响应 %Hr=振幅响应 %b=3型滤波器的系数 %L=Hr的阶次
%h=3型滤波器的单位冲击响应 M=length(h); L=(M-1)/2; c= [2*h(L+1:-1:1)]; n=[0:1:L];
w=[0:1:500]'*2*pi/500; Hr=sin(w*n)*c';
基于频域抽样法的FIR数字带通滤波器设计
wsl=0.12*pi;%低阻带边缘 wsh=0.82*pi;%高阻带边缘 wpl=0.32*pi;%低通带边缘 wph=0.62*pi;%高通带边缘 delta=(wpl-wsl);%过度带 M=ceil(2*pi*3/delta);%抽样点数 al=(M-1)/2;
wl=(2*pi/M); %抽样间隔 k=0:M-1;
精彩文档
实用标准文案
T1=0.12; T2=0.6;%过渡带样本点
Hrs=[zeros(1,ceil(0.12*pi/wl)+1),T2,T1,ones(1,ceil(0.3*pi/wl)),T1,T2,zeros(1,ceil(0.3734*pi/wl)),T2,T1,ones(1,ceil(0.3*pi/wl)),T1,T2,zeros(1,ceil(0.12*pi/wl)+1)]; wdl=[0 0.12 0.32 0.62 0.82 1];
k1=0:floor((M-1)/2);k2=floor((M-1)/2)+1:M-1; angH=[-al*(2*pi)/M*k1,al*(2*pi)/M*(M-k2)]; H=Hrs.*exp(j*angH); h=real(ifft(H));%傅立叶反变换 figure(1);%冲击响应图 stem(k,h);
title('impulse response'); xlabel('n');ylabel('h(n)'); grid;
figure(2);%幅频曲线图 Hf=abs(H);w=k*wl/pi; plot(w,Hf,'*b-') axis([0 1 -0.1 1.1]); title('amplitude response');
xlabel('frequency in pi units');ylabel('Hr(w)'); set(gca,'xtickmode','manual','xtick',wdl);
set(gca,'ytickmode','manual','ytick',[0 0.12 0.6 1]); grid;
精彩文档
实用标准文案
figure(3); fs=15000; [c,f3]=freqz(h,1); f3=f3/pi*fs/2;
plot(f3,20*log10(abs(c))); title('频谱特性');
xlabel('频率/HZ');ylabel('衰减/dB'); grid; t=(0:100)/fs;
x=sin(2*pi*t*700)+sin(2*pi*t*3200)+sin(2*pi*t*6200); q=filter(h,1,x); [a,f1]=freqz(x); f1=f1/pi*fs/2; [b,f2]=freqz(q); f2=f2/pi*fs/2; figure(4); subplot(2,1,1); plot(f1,abs(a)); title('输入波形频谱图'); xlabel('频率');ylabel('幅度') subplot(2,1,2); plot(f2,abs(b));
精彩文档
实用标准文案
title('输出波形频谱图'); xlabel('频率');ylabel('幅度')
基于汉宁窗的FIR数字高通滤波器设计
function s2 Fs=15000; t=(0:100)/Fs;
x=sin(2*pi*500*t)+sin(2*pi*3000*t) subplot(245); stem(x); title('原始信号'); axis([0,100,-2,2]); Ws=7*pi/30; Wp=13*pi/30;
tr_wid=Wp-Ws; %过渡带宽度 N=ceil(11*pi/tr_wid) %滤波器长度 n=0:1:N-1;
Wc=(Ws+Wp)/2; %理想高通滤波器的截止频率 hd=ideal_hp1(Wc,N); %理想高通滤波器的单位冲激响应 w_bla=(blackman(N))'; %布拉克曼
h=hd.*w_bla; %截取得到实际的单位脉冲响应
精彩文档
实用标准文案
[db,mag,pha,grd,w]=freqz_m(h,[1]); %计算实际滤波器的幅度响应 delta_w=2*pi/1000;
As=-round(max(db(1:1:Ws/delta_w+1))) %实际阻带纹波,round是取整函数
y=filter(h,1,x) subplot(246) plot(y)
title('滤波后的信号'); axis([0,100,-1,1]) subplot(241) stem(n,hd)
title('理想单位脉冲响应hd(n)') subplot(242) stem(n,w_bla) title('布拉克满窗w(n)') subplot(243) stem(n,h)
title('实际单位脉冲响应hd(n)') subplot(244) plot(w/pi,db) title('幅度响应(dB)') axis([0,1,-100,10])
精彩文档
实用标准文案
function [db,mag,pha,grd,w] = freqz_m(b,a); [H,w] = freqz(b,a,1000,'whole'); H = (H(1:1:501))'; w = (w(1:1:501))'; mag = abs(H);
db = 20*log10((mag+eps)/max(mag)); pha = angle(H); grd = grpdelay(b,a,w); subplot(247); plot(pha) title('相频响应')
function hd=ideal_hp1(Wc,N) alp=(N-1)/2; n=0:1:N-1;
m=n-alp+eps; %eps是一个很小很小的数 hd=[sin(pi*m)-sin(Wc*m)]./(pi*m);
用双线性法设计巴特沃斯高通数字滤波器
clear all; clc; close all fs=120; T=1/fs; rp=1; rs=30;
Wp=0.35*pi; Ws=0.65*pi; %数字滤波器指标 wp=2*tan(Wp/2)/T; ws=2*tan(Ws/2)/T; %预畸变，将数字滤波器的指标变为模拟滤波器的指标
精彩文档
实用标准文案
[N,w]=buttord(wp,ws,rp,rs,'s'); %求滤波器阶数和3dB截止频率 [Z,P,K]=buttap(N); %设计模拟低通滤波器
[Md,Nd]=zp2tf(Z,P,K); %将零极点形式转换为传输函数形式 [M,N]=lp2hp(Md,Nd,w); %对低通滤波器进行频率变换 [h,w]=freqs(M,N,512); %模拟滤波器的幅频响应 subplot(2,1,1);plot(w,abs(h)); grid; xlabel('Hz');ylabel('幅度'); title('模拟高通滤波器');
[Mh,Nh]=bilinear(M,N,1/T); %对模拟滤波器双线性变换 [h1,w1]=freqz(Mh,Nh); %数字滤波器的幅频响应 subplot(2,1,2);
plot(w1/pi,20*log10(abs(h1))); grid;
xlabel('ω/π');ylabel('幅度(dB)'); title('数字高通滤波器');
%图-5 模拟滤波器与设计的滤波器的单位冲击响应 k=0:2000; k2=1:1001;
x=10*sin(pi/10*k/fs)+5*sin(10*pi*k/fs)+3*sin(30*pi*k/fs); figure subplot(2,1,1) X=fft(x)*2/2001; y=filter(Mh,Nh,x);
plot(k,y); ylim([-5 5]); title('高通数字滤波器输出'); Y=fft(y)*2/2001;
df=fs/2001; ff=(k2-1)*df;
精彩文档
实用标准文案
subplot(2,1,2); plot(ff,abs(X(k2)),'r','linewidth',2); hold on plot(ff,abs(Y(k2)),'b'); title('输入输出频谱比较'); grid;
精彩文档
因篇幅问题不能全部显示，请点此查看更多更全内容