# Third octave bandpass filter with python

I am very new to signal processing and coding which is why my questions might be really basic. I have a signal of the acoustic pressure p'(t) and I would like to use a third octave bandpass filter in Python. I came across these two approaches to filtering with scipy:

At the moment I am able to use one bandpass filter with the following code for a simple example (mostly taken from link 2):

import scipy
import numpy as np
from matplotlib import pyplot as plt
from scipy.signal import butter, lfilter
from scipy.signal import freqz

def butter_bandpass(lowcut, highcut, fs, order=5):
nyq = 0.5 * fs
low = lowcut / nyq
high = highcut / nyq
b, a = butter(order, [low, high], btype='band')
return b, a

def butter_bandpass_filter(data, lowcut, highcut, fs, order=5):
b, a = butter_bandpass(lowcut, highcut, fs, order=order)
y = lfilter(b, a, data)
return y

t = np.linspace(0,1,1000)
s = np.sin(2*np.pi*100*t) + np.sin(2*np.pi*10*t)

T = t - t
fs = 1/T
print(fs)
y = butter_bandpass_filter(s, 90,110,fs,order=5)

plt.figure(1)
plt.plot(t,s,'b',t,y,'r')

plt.figure(2)
N = s.size
f = np.linspace(0, 1/T, N)
s -= np.mean(s)
fft = np.fft.fft(s)
plt.bar(f[:N // 2], np.abs(fft)[:N // 2] * 1 / N, width=5)

plt.figure(3)
N = y.size
f = np.linspace(0, 1/T, N)
y -= np.mean(y)
fft = np.fft.fft(y)
plt.bar(f[:N // 2], np.abs(fft)[:N // 2] * 1 / N, width=5)

plt.show()


The dft shows the results I expected. So my first question would be whether this procedure for a single filter is reasonable?

My second question: How could I use this code to build a third octave bandpass filter with the following center frequencies (for example):

f_one_third = [10,12.5,16,20,25,31.5,40,50,63,80,100,125,160,200,250,315,400,500,630,800,1000,1250,1600,2000,2500,3150,4000,5000,6300,8000,10000,12500,16000,20000]


Simply using one filter after the other seems to be very inefficient and I guess that there is a much more efficient. The first link gives an approach for a third octave filter, but I am honestly not sure how to "use". I am basically missing the equivalent of the def butter_bandpass_filter function from the code above.

I would appreciate any help immensely!

• you should generally use sosfilt, not lfilter. we even added this in the docstring :) Aug 10 '20 at 15:10

I recently developed a function to easily perform octave and fractional octave filtering, it is available on github: PyOctaveBand

It uses the SOS coefficients and performs downsampling to filter correctly at low frequencies.

The dft shows the results I expected. So my first question would be whether this procedure for a single filter is reasonable?

Sort of. For very low frequencies and/or high orders, using transfer function representation (b,a) runs into numerical problems. Use second order section representation and sosfilt() instead.

My second question: How could I use this code to build a third octave bandpass filter with the following center frequencies (for example):

Depends a bit on what you are planning to do with the resulting signals. Just doing lots of parallel filters is not an unreasonable approach. The alternative would be a frequency domain method, but this is unlikely to be cheaper if you need all results in the time domain. However, if you just want to build a 3rd octave spectrum analyzer, you just need the time averaged band energy as the function of time, which can be done in the frequency domain pretty easily.

• Thanks! Alright, I will use sosfilt() then. My desired final ouput is a frequency spectrum in which the amplitudes are assigned to the third octave center frequencies, I think that this is a 3rd octave spectrum analyzer. Does that mean that I don't even need a bandpass filter like the one in the code above? Mar 5 '19 at 19:39