# Ripple in a SciPy Butterworth SOS Filter

I have a script using SciPy for checking the RMS of various Butterworth Bandpass Filters of varying orders.

I would expect the RMS values to increase and decrease linearly and consistently as you get closer and further from the passband, but this is not what happens at all. The RMS fluctuates dependent on the frequency, for instance, in the example below, 162Hz has a much higher RMS than 158, despite 162 being further from the high-pass cut-off of 133Hz. This appears to be a cyclical pattern, and independent of the order, but I am not good enough with matplotlib to create a pretty chart that visually represents this issue. The exact cycle also appears to be very dependent on the buffer size.

I am using https://www.szynalski.com/tone-generator/ to test tones, and here is my code:

import numpy as np
import scipy.signal
import pyaudio

def normalize(block):
array = np.frombuffer(block, dtype=np.float32)
return array

def get_rms(samples: np.array) -> float:
"""Get the RMS of an array of audio samples

Args:
samples: the samples to get the RMS from

Returns:
float: the RMS
"""
samples_array = np.array(samples)
return np.sqrt(np.mean(samples_array ** 2))  # type:ignore

def design_filter(lowcut, highcut, fs, order=3):
nyq = 0.5 * fs
low = lowcut / nyq
high = highcut / nyq
sos = scipy.signal.butter(order, [low, high], btype="band", output="sos")
return sos

def main():
sample_rate = 44100
buffer_size = 2048

filters = {}
for i in range(10):
sos = design_filter(101, 133, sample_rate, i)
zi = scipy.signal.sosfilt_zi(sos)
filters[i] = [sos, zi]

stream = pyaudio.PyAudio().open(
format=pyaudio.paFloat32,
channels=1,
rate=sample_rate,
input=True,
frames_per_buffer=buffer_size,
)

update_every = 3
update = 0
while True:

if update_every == update:
update = 0
samples = normalize(block)

# blank out terminal
print(chr(27) + "[2J")
# move cursor to to left
print(chr(27) + "[1;1f")

print("rms")
print(int(get_rms(samples) * 200) * "-")

for order, tup in filters.items():
print(f"Order: {order}")
bandpass_samples, zi = scipy.signal.sosfilt(tup[0], samples, zi=tup[1])
tup[1] = zi.copy()
print(int(get_rms(bandpass_samples) * 200) * "-")
else:
update += 1

if __name__ == "__main__":
main()


UPDATE After looking at some more examples and messing with my code, it appears that it is more related to the buffer size than I initially thought.

If buffer_size == sample_rate the ripple (?, oscillation? noise?) does not happen at all (this is wrong, see below). Also, there is seemingly a relationship between factors of buffer_size and sample_rate that relates to the division. Not sure why this is, and I am able to work around it, but if anyone knows why this happens it would be good to know, for sure.

UPDATE 2 It appears that the above is incorrect. There is some transient noise that appears every period. Having a larger period decreases the relative impact of that noise, but it does not eliminate it, and it does not allow for a work around of the issue by increasing the buffer size and then chopping up the filter output later to get the granularity back.

• Hint : the RMS calculation only works if you have a whole number of periods. Think : When does this happen?
– Ben
Oct 20, 2020 at 0:48

I can find 2 errors at a quick glance :

1 - You should discard the transient when measuring the steady-state gain at for a given frequency $$f$$. I haven't checked how long the transient last, but you should factor it in your RMS measurement.

2 - The RMS calculation only works when the number of samples corresponds to a whole number of periods. Otherwise, the calculation is biased.

Edit :

Period in this context means the signal period. As you are in the discrete domain, this is the ratio between the sampling frequency and the signal frequency

$$T = \frac{f_s}{f}$$

For transient, see this site

Edit 2 :

Assuming your signal frequency is 1 kHz

$$T = \frac{44.1 kHz}{1 kHz} = 44.1$$

Since you have 44.1 samples per period, the RMS calculation with only 44 samples will have some bias. However, if you use 441 samples, you will have exactly 10 signal periods and then the RMS calculation will have no bias.

As for the transient, you should plot the time response. You should be able to evaluate how long the transient lasts.

• I apologize in advanced for being big dumb. What is "the transient" in this context? What is a "period" in this context? Oct 20, 2020 at 0:59
• I really appreciate your help. Unfortunately I am even less sure of what you are saying after the update. Oct 20, 2020 at 2:45
• Would the signal period be one second in this case? 44,100 samples? Oct 20, 2020 at 2:45
• If so, are you saying that RMS only works if you are calculating it for a factor of 44,100 samples? AKA only once every second? That doesn't seem right but also I have no idea what I am doing. Oct 20, 2020 at 2:47
• The transient duration will be filter-specific
– Ben
Oct 20, 2020 at 15:20