From 438bad314769e282f3f3c70a408058cebfe2b428 Mon Sep 17 00:00:00 2001
From: Surabhi Saketh <saketh.surabhi@bitsilica.com>
Date: Tue, 20 Feb 2024 10:59:31 +0000
Subject: [PATCH]
---
Fir_filter-matlab-code.md | 146 ++++++++++++++++++++++++++++++++++++++
1 file changed, 146 insertions(+)
create mode 100644 Fir_filter-matlab-code.md
diff --git a/Fir_filter-matlab-code.md b/Fir_filter-matlab-code.md
new file mode 100644
index 0000000..f59caf4
--- /dev/null
+++ b/Fir_filter-matlab-code.md
@@ -0,0 +1,146 @@
+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