%% FIR陷波器设计
%% 等波纹最佳逼近法Equiripple Best Approximation Method===%%
% 目的:设计一个陷波器阻带在50±1.5Hz以内,采样频率为400Hz的滤波器,
% 并要求通带最大衰减为0.1dB,阻带最小衰减为80dB。
clc;
clear;close all;
rp=0.1; % 通带波纹最大衰减为0.1dB
rs=80; % 阻带衰减为80dB
fs=400; % 采样频率为400Hz
f=[48.5,49.9,50.1,51.5]; % 根据阻带:49.9Hz50.1Hz,通带:048.5Hz和51.5Hz~pi,w=mo_wT=mo_w/fs=2pi*f/fs
m=[1 0 1]; % 期望幅度
dat1=(10(rp/20)-1)/(10(rp/20)+1); % 根据(7.4.4)和(7.4.5)式计算通带和阻带波纹幅度dat1,dat2和rip
dat2=10^(-rs/20);
rip=[dat1 dat2 dat1];
[M,fo,mo,wc]=remezord(f,m,rip,fs);
hn=remez(M,fo,mo,wc); % 得到滤波器单位脉冲响应
[hw,w]=freqz(hn);
%% 产生信号== %%
Nlen=2038400;
n=0:Nlen-1;
Ts=1/fs;
xn=sin(2pi.48nTs)+sin(2pi.48.5nTs)+sin(2pi.49nTs)+sin(2pi50.nTs)+sin(2pi.51nTs)+sin(2pi.51.5nTs)+sin(2pi.52nTs); % 信号叠加
k=2*(0:Nlen-1)/Nlen;
yn=fftfilt(hn,xn,Nlen);
xk_T=fft(xn,Nlen);
yk=fft(yn,Nlen);
%% 波形绘制== %%
figure;
subplot(3,2,1);
plot(wfs/(2pi),abs(hw));
title(‘滤波器幅频特性’);xlabel(‘f/Hz’);ylabel(‘H(ejf)’);
axis([45,55,0,1.2]);
grid on;
subplot(3,2,2);
plot(wfs/(2pi),20log10(abs(hw))/max(abs(hw)));
title(‘损耗函数’);xlabel(‘f/Hz’);ylabel(‘A(s)/dB’);
axis([45,55,min(20log10(abs(hw))/max(abs(hw))),10]);
grid on;
subplot(3,2,3);
plot(n,xn);
title(‘原信号波形图’);xlabel(‘f/Hz’);ylabel(‘幅度’);
axis([0,2000,min(xn),max(xn)]);
grid on;
subplot(3,2,4);
plot(kfs/2,abs(xk_T));
title(‘滤波前输入信号x(n)的频谱’);xlabel(‘f/Hz’);ylabel(‘H(ejf)’);
axis([47,53,0,max(abs(xk_T))]);
grid on;
subplot(3,2,5);
plot(n,yn);
title(‘滤波后波形图’);xlabel(‘f/Hz’);ylabel(‘幅度’);
axis([0,2000,min(yn),max(yn)]);
grid on;
subplot(3,2,6);
plot(kfs/2,abs(yk));
title(‘滤波后输出信号y(n)的频谱’);xlabel(‘f/Hz’);ylabel(‘H(ejf)’);
axis([47,53,0,max(abs(yk))]);
grid on;