1
$\begingroup$

I'm stuck trying to filter a two-channel signal using sosfilt() with the zi initial state. I need to use the zi parameter so I can carry the filter state over multiple buffers.

Omitting zi works as expected

import numpy as np
from scipy.signal import sosfilt

>>> x = np.arange(20, dtype=np.float32).reshape(2, -1)
>>> x
array([[ 0.,  1.,  2.,  3.,  4.,  5.,  6.,  7.,  8.,  9.],
       [10., 11., 12., 13., 14., 15., 16., 17., 18., 19.]], dtype=float32)
>>> sos = np.array([2, 0, 0, 1, 0, 0], dtype=np.float32)  # a very simple filter :)
>>> sos
array([2., 0., 0., 1., 0., 0.], dtype=float32)
>>> y = sosfilt(sos, x)
>>> y
array([[ 0.,  2.,  4.,  6.,  8., 10., 12., 14., 16., 18.],
       [20., 22., 24., 26., 28., 30., 32., 34., 36., 38.]], dtype=float32)

... but what does sosfilt() need for zi?

The sosfilt() documentation for zi says:

Initial conditions for the cascaded filter delays. It is a (at least 2D) vector of shape (n_sections, ..., 2, ...), where ..., 2, ... denotes the shape of x, but with x.shape[axis] replaced by 2.

The filter has two sections, so n_sections = 2. I assume that x.shape[axis] = 2 (because I can't imagine that zi needs to know that the source is 10 frames long). So I tried a few things, but none of them worked:

>>> zi = np.zeros((2, 2, 2), dtype=np.float32)
>>> y, zf = sosfilt(sos, x, zi)
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "...\_signaltools.py", line 4277, in sosfilt
    x_zi_shape[axis] = 2
TypeError: only integer scalar arrays can be converted to a scalar index
>>> zi = np.zeros((1, 2, 2), dtype=np.float32) # maybe it means one second-order section?
>>> y, zf = sosfilt(sos, x, zi)
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "...\_signaltools.py", line 4277, in sosfilt
    x_zi_shape[axis] = 2
TypeError: only integer scalar arrays can be converted to a scalar index
>>> zi = np.zeros((2, 2), dtype=np.float32) # at least 2D?
>>> y, zf = sosfilt(sos, x, zi)
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "...\_signaltools.py", line 4277, in sosfilt
    x_zi_shape[axis] = 2
TypeError: only integer scalar arrays can be converted to a scalar index

What am I missing?

$\endgroup$
2
  • $\begingroup$ Does sosfilt_zi tell you anything? $\endgroup$
    – Peter K.
    Feb 28, 2023 at 18:06
  • $\begingroup$ I should have looked at sosfilt_zi output, but there was another issue hiding the answer. See below... $\endgroup$ Feb 28, 2023 at 18:08

1 Answer 1

3
$\begingroup$

After more fiddling, I figured it out:

  1. need to specify axis=1
  2. n_sections is the number of second-order sections, not filter order.

So if N is the number of input channels, the correct initialization for zi and the correct call to sosfilt() is:

zi = np.zeros((1, N, 2), dtype=np.float32)
y, zf = sosfilt(sos, x, zi=zi, axis=1)
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.