I'm building a bandstop filter but the bandstop function takes a long long time to run, so I decided to try filtering with a Butterworth filter, but I don't get the same response as with bandstop. Can you help me please?
%% Frequency range
lfreq = 170; % Lower freq
hfreq = 180; % Highest freq
%% Extract audio signal
[f Fs] = audioread('Schumann.wav'); % Extract signal
info = audioinfo('Schumann.wav'); % Extract audio info
n = length(f); % Length of signal
t = 0:seconds(1/Fs):seconds(info.Duration); % Create a x-axis of time in s
t = t(1:end-1); % Set time
freq = 1/time2num(t(2)-t(1))/n*(0:n); % Create a x-axis of frequencies in Hz
freq = freq(1:end-1); % Set frequency
L = 1:floor(n/2); % Only take the first half of freq
%% Compute the Fast Fourier Transform
fhat = fft(f,n); % Compute the FFT
PSD = fhat.*conj(fhat)/n; % Power spectum (energy per freq)
%% Use "bandstop" to filter frequency range
% ffilt = bandstop(f,[lfreq hfreq],Fs); % Remove all freq in range
% fhatfilt = fft(ffilt,n); % Compute the FFT for filtered signal
% PSDfilt = fhatfilt.*conj(fhatfilt)/n; % Filtered power spectum
%% Butterworth filter <----
Wn = [lfreq hfreq]/n;[b a] = butter(5,Wn,'stop');ffilt = filtfilt(b,a,f);%% Plots (original and filtered signal)
figure(1)plot(freq(L),PSD(L),'r','LineWidth',1.5); hold onplot(freq(L),PSDfilt(L),'b','LineWidth',1.5);legend('Original','Filtered');xlabel('Frequency [Hz]'); ylabel('PSD [J/Hz]');figure(2)plot(time2num(t),f,'r','LineWidth',1.5); hold onplot(time2num(t),ffilt,'b','LineWidth',1.5); legend('Original','Filtered');xlabel('Time [s]'); ylabel('Amplitude [Unknown]');%sound(ffilt,Fs)
Best Answer