function custom_filter_design_tool_4() f = input(‘Enter the frequency (in Hz): ‘);

if f >= 1000
    [Fpass, Fstop, Wpass, Wstop, N, dens] = set_values_forxxx(f);
    filter_type = 'firpm';
else
    [Fpass, Fstop, Wpass, Wstop, N, dens] = set_values_yy(f);
    filter_type = 'firls';
end

Fs = 40000;

% Filter design based on parameters
if f < 1000
    b = firls(N, [0 Fpass Fstop Fs/2]/(Fs/2), [1 1 0 0], [Wpass Wstop]);
else
    b = firpm(N, [0 Fpass Fstop Fs/2]/(Fs/2), [1 1 0 0], [Wpass Wstop], dens);
end

% Display the filter type
disp(['Filter Type: ' filter_type]);

% Initialize the frequency step
step_size = input('Enter the frequency step size (50 or 100 Hz): ');
if step_size ~= 50 && step_size ~= 100
    error('Invalid step size. Please enter either 50 or 100 Hz.');
end

% Compute magnitude at the given frequency and iterate until the nearest multiple of the step size
while mod(f, step_size) ~= 0
    if f > (f - mod(f, step_size))
        f = f - step_size;
    else
        f = f + step_size;
    end
end

while true
    [magnitude_dB,last_Fpass,Fstop_last, Wpass_Last,Wstop_last, c] = compute_magnitude_at_frequency(b, f, Fs, N, Fpass, Fstop, Wpass, Wstop, dens);
    disp(['Magnitude at ' num2str(f) ' Hz: ' num2str(magnitude_dB) ' dB']);
    disp(['Last value of Wpass: ', num2str(Wpass_Last)]);
    disp(['Last value of Wstop: ', num2str(Wstop_last)]);
    disp(['Last value of Fstop(Hz): ', num2str(Fstop_last)]);
    disp(['Last value of Fpass(Hz): ', num2str(last_Fpass)]);

    % Plot magnitude response
    figure;
    freqz(c, 1, 1024, Fs);
    filename = sprintf('filter_coefficients_6_%dHz_%s.txt', f, filter_type);
    save_filter_coefficients(filename, c);

    % Update frequency for the next iteration
    f = f - step_size;
    if f <= 0
        break;
    end
end

end

function save_filter_coefficients(filename, c) fid = fopen(filename, ‘wt’); fprintf(fid, ‘%f\n’, c); fclose(fid); disp([‘Filter coefficients saved to ' filename]); end

function [Fpass, Fstop, Wpass, Wstop, N, dens] = set_values_forxxx(f) % For frequencies above 1 kHz Fpass = f - 100; Fstop = f + 100; Wstop = 755; Wpass = 2; N = 8; dens = 20; end

function [Fpass, Fstop, Wpass, Wstop, N, dens] = set_values_yy(f) % For frequencies below 1 kHz Fstop = f + 10; Fpass = f - 10; Wstop = 1; Wpass = 1; N = 22; dens = 20; end

function [magnitude_dB, last_Fpass, Fstop_last, Wpass_Last,Wstop_last, c] = compute_magnitude_at_frequency(b, f, Fs, N, Fpass, Fstop, Wpass, Wstop, dens) while true if f < 1000 c = firls(N, [0 Fpass Fstop Fs/2]/(Fs/2), [1 1 0 0], [Wpass Wstop]); else c = firpm(N, [0 Fpass Fstop Fs/2]/(Fs/2), [1 1 0 0], [Wpass Wstop], dens); end

    [H, freq] = freqz(c, 1, 1024, Fs);
    [~, idx] = min(abs(freq - f));
    magnitude_dB = 20*log10(abs(H(idx)));

    % Adjust Fpass and Fstop based on the magnitude
    if f < 1000
        if magnitude_dB < -48
            Wpass = Wpass + 1;
        elseif magnitude_dB > -45
            Fstop = Fstop - 10;
            Fpass = Fpass - 10;
            Wstop = Wstop + 10;
        else
            break; % Exit loop if magnitude is within desired range
        end
    else
        if magnitude_dB < -46 && f >= 3000
            Fpass = Fpass + 100;
            Fstop = Fstop + 100;
            Wstop = Wstop - 10;
        elseif magnitude_dB > -40 && f >= 3000
            Fpass = Fpass - 100;
            Fstop = Fstop - 100;
            Wstop = Wstop + 10;
        elseif magnitude_dB < -47 && f < 3000
            Fpass = Fpass + 10;
            Fstop = Fstop + 10;
            Wstop = Wstop - 1;
        elseif magnitude_dB > -40 && f < 3000
            Fpass = Fpass - 10;
            Fstop = Fstop - 10;
            Wstop = Wstop + 10;
        else
            break; % Exit loop if magnitude is within desired range
        end
    end
end

% Return the last value of Wpass
% last_Wpass = Wpass;
% Fstop_last = Fstop;
% Fpass_Last = Fpass;
% Wstop_last = Wstop;
    last_Fpass = Fpass;
    Fstop_last = Fstop;
    Wpass_Last = Wpass;
    Wstop_last = Wstop;
Hd = dsp.FIRFilter('Numerator', c);

freqz(Hd);

end