Problems using a Butterworth filter instead of a bandstop to filter frequencies.

Question

Manuel Cruz 2020-11-29

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/668038-problems-using-a-butterworth-filter-instead-of-a-bandstop-to-filter-frequencies

评论： Star Strider 2020-11-29

Schumann.zip

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 on
plot(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 on
plot(time2num(t),ffilt,'b','LineWidth',1.5); 
legend('Original','Filtered');
xlabel('Time [s]'); ylabel('Amplitude [Unknown]');
%sound(ffilt,Fs)

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Star Strider 2020-11-29

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/668038-problems-using-a-butterworth-filter-instead-of-a-bandstop-to-filter-frequencies#answer_559538

在 MATLAB Online 中打开

You are not designing the filter correctly.

Try this instead:

lfreq = 170;                                                        % Lower freq
hfreq = 180;                                                        % Highest freq
Fs = 44100;                                                         % Sampling Frequency
Fn = Fs/2;                                                          % Nyquist Frequency
Ws = [lfreq hfreq]/Fn;                                              % Stopband Frequency (Normalised)
Wp = [0.99 1.01].*Ws;                                               % Passband Frequency (Normalised)
Rp =  1;                                                            % Passband Ripple
Rs = 90;                                                            % Passband Ripple (Attenuation)
[n,Wp] = ellipord(Wp,Ws,Rp,Rs);                                     % Elliptic Order Calculation
[z,p,k] = ellip(n,Rp,Rs,Wp,'stop');                                 % Elliptic Filter Design: Zero-Pole-Gain 
[sos,g] = zp2sos(z,p,k);                                            % Second-Order Section For Stability
figure
freqz(sos, 2^18, Fs)                                                % Filter Bode Plot
set(subplot(2,1,1), 'XLim',Wp*Fn)                                   % Optional
set(subplot(2,1,2), 'XLim',Wp*Fn)                                   % Optional