function [filt, Y1, P1] = calc_filter(band_num, Nfilt, Nfilt1, plot_filter)
% calculate band_num FIR filters "filt".
% ================================================
% band_num = number of bands or filters.
% Nfilt = number of points in each filter.
% ================================================
if nargin < 4
    plot_filter = 0;
end

if nargin < 3
    Nfilt1 = Nfilt;
end

Nfilt_vec = (0:Nfilt - 1) - floor(Nfilt/2);
d1 = calc_d(1e-2, 1e-2);	% Normalized transition width.
% d1 = calc_d(5e-3, 5e-3);	% Normalized transition width.
% Fp = 1/band_num - d1/Nfilt;
% Fs = 1/band_num + d1/Nfilt;
Fp = 1/band_num - d1/Nfilt;
Fs = 1/band_num + d1/Nfilt;
F = [0 Fp Fs 1];
m = [1 1 0 0];

% FIR filter "filt"
% ===================
filt_1 = firls(Nfilt - 1, F, m).';     % center band.
fvec = ceil(-band_num/2):ceil(band_num/2) - 1;
fvec = fvec*2*pi/band_num;
filt = filt_1(:).*exp(-1j*Nfilt_vec(:)*fvec);

% frequency response of "filt"
% ======================================
Y1 = zeros(Nfilt1, band_num);
for kk = 1:band_num
    [Y1(:, kk), P1] = freqres(filt(:, kk), Nfilt1);
end

% plot filt
% ==============
if plot_filter > 0
    if mod(Nfilt, 2) == 1
        figure; plot(P1, real(Y1), P1, imag(Y1), 'm'); grid;
        legend(' real ', ' imag ', 'Location', 'northeast');
    else
        figure; plot(P1, abs(Y1)); grid;
    end
    ax = gca;
    ax.XTick = pi*linspace(-1, 1, 5);
    ax.XMinorTick = 'on';
    ax.XTickLabel = {'-\pi', '-0.5\pi', '0', '0.5\pi', '\pi'};
    title([' desired ', num2str(band_num), '-band FIR filters: Nfilt = ', num2str(Nfilt), '. ']);
    xlabel(' \phi (rad) ');
end
