MATLAB窗函数实现带通滤波器

%本函数利用窗函数法设计带通滤波器，主要用来滤出单一频率，即中心频率
%data是输入的数据, centerFre是带通的中心频率, offsetFre是频偏,最终带通为centerFre +- offsetFre/2
%,sampFre是采样率
function y = BPassFilter(data, centerFre, offsetFre, sampFre)
    %设计I型带通滤波器
    M = 0 ;    %滤波器阶数（必须是偶数）
    Ap = 0.82; %通带衰减
    As = 45;   %阻带衰减
    Wp1 = 2*pi*(centerFre - offsetFre)/sampFre;  %算出下边频
    Wp2 = 2*pi*(centerFre + offsetFre)/sampFre;  %算出上边频
   
  
  
    % (1)矩形窗
    N = ceil(3.6*sampFre/offsetFre);             %计算滤波器阶数,采用矩形窗，3dB截频在中心频率到上下边频的中点
    M = N - 1;
    M = mod(M,2) + M ; %使滤波器为I型(偶数)
   
    %单位脉冲响应的下脚标
    h = zeros(1,M+1);  %单位冲击响应变量赋初值
    for k = 1:(M+1);
        if (( k -1 - 0.5*M)==0)
            h(k) = Wp2/pi - Wp1/pi;
        else
            h(k) = Wp2*sin(Wp2.*(k - 1 - 0.5*M))/(pi*(Wp2*(k -1 - 0.5*M))) - Wp1*sin(Wp1*(k - 1 - 0.5*M))/(pi*(Wp1*(k -1 - 0.5*M)));
        end
    end
   
   
  % (2) Hann Window
%   N = ceil(12.4*sampFre/offsetFre);             %计算滤波器阶数,采用矩形窗，3dB截频在中心频率到上下边频的中点
%   M = N - 1;
%   M = mod(M,2) + M ; %使滤波器为I型(偶数)
%       h = zeros(1,M+1);  %单位冲击响应变量赋初值
%     for k = 1:(M+1);
%         if (( k -1 - 0.5*M)==0)
%             h(k) = Wp2/pi - Wp1/pi;
%         else
%             h(k) = Wp2*sin(Wp2.*(k - 1 - 0.5*M))/(pi*(Wp2*(k -1 - 0.5*M))) - Wp1*sin(Wp1*(k - 1 - 0.5*M))/(pi*(Wp1*(k -1 - 0.5*M)));
%         end
%     end
%   K = 0:M;
%   w = 0.5 - 0.5*cos(2*pi*K/M);
%   h = h.*w;
 

  % (3)Hamming Window
%   N = ceil(14*sampFre/offsetFre);             %计算滤波器阶数,采用矩形窗，3dB截频在中心频率到上下边频的中点
%   M = N - 1;
%   M = mod(M,2) + M ; %使滤波器为I型(偶数)
%     h = zeros(1,M+1);  %单位冲击响应变量赋初值
%     for k = 1:(M+1);
%         if (( k -1 - 0.5*M)==0)
%             h(k) = Wp2/pi - Wp1/pi;
%         else
%             h(k) = Wp2*sin(Wp2.*(k - 1 - 0.5*M))/(pi*(Wp2*(k -1 - 0.5*M))) - Wp1*sin(Wp1*(k - 1 - 0.5*M))/(pi*(Wp1*(k -1 - 0.5*M)));
%         end
%     end
%   K = 0:M;
%   w = 0.54 - 0.46*cos(2*pi*k/M);
%   h = h.*w;

% (4)Blackman window
%   N = ceil(22.8*sampFre/offsetFre);             %计算滤波器阶数,采用矩形窗，3dB截频在中心频率到上下边频的中点
%   M = N - 1;
%   M = mod(M,2) + M ; %使滤波器为I型(偶数)
%     h = zeros(1,M+1);  %单位冲击响应变量赋初值
%     for k = 1:(M+1);
%         if (( k -1 - 0.5*M)==0)
%             h(k) = Wp2/pi - Wp1/pi;
%         else
%             h(k) = Wp2*sin(Wp2.*(k - 1 - 0.5*M))/(pi*(Wp2*(k -1 - 0.5*M))) - Wp1*sin(Wp1*(k - 1 - 0.5*M))/(pi*(Wp1*(k -1 - 0.5*M)));
%         end
%     end
%   K = 0:M;
%   w = 0.42 - 0.5*cos(2*pi*K/M) + 0.08*cos(4*pi*K/M);
%   h = h.*w;
 
  y = filter(h,[1],data);