Unable to fully filter out a tone with a notch filter.

6 views (last 30 days)
I am working on coding a matlab notch filter that is supposed to fitler out a tone in the audio file listed in the code. The Notch Filter works as it is supposed to except at the very beginning and at the very end the tone is still there and not filtered out and I need to have the tone completely gone. I am unable to figure out what to change in my code I have tried changing the frequencies I have tried changing where my filter starts and my notch width but nothing has solved my issue. Any help or suggestions is appreciated Thanks!
Code here below
% Load input sound
[Sound, fs] = audioread('SunshineSquare.wav');
figure(1);
% Plot specgram of input sound
specgram(Sound, 512, fs);
title('Spectrogram of Input Sound (SunshineSquare.wav)');
% Parameters of Notch Filter
frequencies = [1574.99, 3149.97, 4724.96];
notchwidth = 10;
% Tone Removal
a_total = 1;
b_total = 1;
for i = 1:length(frequencies)
Wn = [(frequencies(i) - notchwidth/2) / (fs/2), (frequencies(i) + notchwidth/2) / (fs/2)];
Wn = max(min(Wn, 0.99), 0.01);
[b, a] = butter(2, Wn, 'stop');
b_total = conv(b_total, b);
a_total = conv(a_total, a);
end
filteredAudio = filter(b_total, a_total, Sound);
% Play filtered sound
soundsc(filteredAudio, fs);
figure(2);
% Plot specgram of output sound after using notch filter
specgram(filteredAudio, 512, fs);
title('Specgram after Notch Filter');
% Plot SunshineSquare sound
figure(3);
t = (0:length(Sound)-1)/fs;
plot(t, Sound);
xlabel('Time (sec)');
ylabel('Amplitude');
title('Notch-Filtered SunshineSquare Sound Plot');
figure(4);
plot (ww,abs(HH));
zplane(bb,aa)
xlabel('Real');
ylabel('Imaginary');
title('Magnitude IIR Frequency Response (Z-Plane)');
grid on
figure(5);
% Plot Notch-filtered sound
t = (0:length(Sound)-1)/fs;
plot(t, filteredAudio);
xlabel('Time (sec)');
ylabel('Amplitude');
title('Notch-Filtered Z-Plot');
figure(6);
sgtitle('Redo Tone Removal Step 2: Magnitude Response for Notch Filter');
f1 = 1574.99/fs * 2 * pi;
f2 = 3149.97/fs * 2 * pi;
f3 = 4724.96/fs * 2 * pi;
% Plot nulling filter frequency response and Z-plane representation
ww = -pi:pi/100:pi;
HH = freqz(bb,aa,ww);
bb = poly([exp(1j*f1) exp(-1j*f1) exp(1j*f2) exp(-1j*f2) exp(1j*f3) exp(-1j*f3)]);
aa = poly(0.80*[exp(1j*f1) exp(-1j*f1) exp(1j*f2) exp(-1j*f2) exp(1j*f3) exp(-1j*f3)]);
y = filteredAudio;
% Display the nulling filter frequency response
subplot(1,1,1);
plot(ww, abs(HH));
title('Notch Filter Frequency Response');
xlabel('Frequency (Hz)');
ylabel('Magnitude');
fout = filter(bb,aa,y);
% Calculate the DC gain of the IIR Notch filter
ww_low = -pi:pi/10000:pi; % Frequencies close to zero
HH_low = freqz(b_total, a_total, ww_low, fs);
% Display the DC gain
dc_gain_iir = abs(HH_low(1));
disp(['DC Gain of IIR Notch Filter: ', num2str(dc_gain_iir)]);
Specgrams of the input and output of the file from matlab for reference are attached.

Answers (1)

Star Strider
Star Strider on 1 May 2024
It is not possible to run your code.
The filter may be working, however note that it is only order 2. That may be appropriate for some designs (for example an elliptic filter), however it may not be enough for a Butterworth design. Each filter design has a function to estimate the appropriate order and return other values, for Butterworth filters it is buttord. Using it might help with your filter design, and returning the zero-pole-gain output from butter and then using zp2sos and then filtfilt might also produce better results. So long as you are designing the filters in a loop, an equivalent approach could be to initially filter the signal with the first filter and then the output of that filter with the subsequent filters. It might be easier to use the bandstop function, actually (include 'ImpulseResponse','iir'), unless you have to design your own filters.
Onee problem however is how you are analysing the results. The spectrogram arguments are first the signal, second the window length, and third the overlap. In your call to it, the overlap is equivalent to the sampling frequency. (The sampling frequency is actually the fifth argument to spectrogram.) With a smaller overlap (for example using the default value), and using square brackets [] for the arguments you do not want to specify (in that event they become the function default values), the filtered results might be mor to your liking.

Tags

Products


Release

R2022b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!