MATLAB: Butterworth Bandpass filter issues

codefilterMATLAB

Hi there! I am trying to make a butterworth bandpass filter in matlab and am having some trouble with it. The specifications for my filter are a center frequency of 10kHz with limits of +/- 30 Hz. In addition to this, I'd like to have 100dB of attenuation at 1kHz away from the center frequency, so if anyone could shed some light on how I could do this that would be awesome!

Right now I am trying to test my filter by simulating a 10kHz sinusoidal signal in LTSpice, exporting the time-voltage data to matlab and then resampling the signal so that it is periodic. I have then made and applied a butterworth bandpass filter with passband 9970Hz – 10030 Hz in the hope that I can completely allow the 10kHz signal to pass through. Right now this is not the case and I am getting some cutoff which progressively increases with time. Running the provided matlab script with the provided text file will illustrate this. I have attached the matlab code as well as the exported time-voltage text file.

Any help would be greatly appreciated!

Best Answer

That is likely asking too much of a Butterworth design. Use an elliptic filter instead:

Fs = 44100;                                                         % Sampling Frequency
Fn = Fs/2;                                                          % Nyquist Frequency
Wp = [1E+3-30 1E+3+30]/Fn;                                          % Passband Frequency (Normalised)
Ws = [1E+3-100 1E+3+100]/Fn;                                        % Stopband Frequency (Normalised)
Rp =   1;                                                           % Passband Ripple
Rs = 100;                                                           % Passband Ripple (Attenuation)
[n,Wp] = ellipord(Wp,Ws,Rp,Rs);                                     % Elliptic Order Calculation
[z,p,k] = ellip(n,Rp,Rs,Wp,'bandpass');                             % Elliptic Filter Design: Zero-Pole-Gain 
[sos,g] = zp2sos(z,p,k);                                            % Second-Order Section For Stability
figure
freqz(sos, 2^16, Fs)                                                % Filter Bode Plot

Use your own sampling frequency.

.

Related Solutions

MATLAB: Design band pass filter

Here is prototype code I used for a different project that you can adapt, if you want an IIR filter:

Fs = 250;                                                       % Sampling Frequency
Fn = Fs/2;                                                      % Nyquist Frequency
Ts = 1/Fs;                                                      % Estimate Sampling Interval
t = linspace(0, numel(EKG), numel(EKG))/Fs;
Wp = [2  90]/Fn;                                                % Normalised Passband (Passband =  2 Hz To  90 Hz)
Ws = [1 100]/Fn;                                                % Normalised Stopband (Passband =  1 Hz To 100 Hz)
Rp =  1;                                                        % Passband Ripple/Attenuation
Rs = 50;                                                        % Stopband Ripple/Attenuation
[n,Wp] = ellipord(Wp, Ws, Rp, Rs);                              % Calculate Elliptic Filter Optimum Order
[z,p,k] = ellip(n, Rp, Rs, Wp,'bandpass');                      % Elliptic Filter
[sos,g] = zp2sos(z,p,k);                                        % Convert To Second-Order-Sections For Stability
figure
freqz(sos,2^16,Fs)
filtered_signal = filtfilt(sos,g,signal);

There are a number of FIR designs on Answers as well (I have posted many of them), although I usually reserve those for filters with multiple passbands and stopbands.

MATLAB: Filter design parameter: radian/sec or 1/sec

The passbands and stopbands for the Signal Processing Toolbox filter functions are actually normalised to the closed interval [0, pi] radians. It is calculated as the passband and stopband frequency in Hz divided by the Nyquist frequency.

Example for frequency calculations for the filters your signal:

Ts = 0.002;                                             % Sampling Interval (seconds)
Fs = 1/Ts;                                              % Sampling Frequency (Hz)

Fn = Fs/2;                                              % Nyquist Frequency

Note that the highest frequency you can design in your filter is the Nyquist frequency, here 250 Hz.

An example of a Chebyshev Type II filter design for filtering an EKG signal is here:

Fs = 250;                                               % Sampling Frequency (Hz)
Fn = Fs/2;                                              % Nyquist Frequency (Hz)
Wp = [1.0   100]/Fn;                                    % Passband Frequencies (Normalised)
Ws = [0.5   101]/Fn;                                    % Stopband Frequencies (Normalised)
Rp = 10;                                                % Passband Ripple (dB)
Rs = 50;                                                % Stopband Ripple (dB)
[n,Ws] = cheb2ord(Wp,Ws,Rp,Rs);                         % Filter Order
[z,p,k] = cheby2(n,Rs,Ws);                              % Filter Design
[sosbp,gbp] = zp2sos(z,p,k);                            % Convert To Second-Order-Section For Stability
figure(3)
freqz(sosbp, 2^17, Fs)
EKG_filt = filtfilt(sosbp, gbp, EKG_original);          % Filter EKG Signal

The filter passband goes from 1 Hz to 100 Hz.

Related Question