[SciPy-Dev] signal.sosfilt filter representation

[Eric Moore]

On Fri, Jun 5, 2015 at 8:29 AM, Eric Moore <ewm at redtetrahedron.org> wrote:

> On Wed, Jun 3, 2015 at 11:06 PM, Eric Larson <larson.eric.d at gmail.com>
> wrote:
>
>> Hey Eric,
>>
>> Your proposal sounds good to me. I can see the value in being able to
>> choose between DF2 and biquad methods based on considerations of speed and
>> quantization error based on the input and filter characteristics. I agree
>> that a keyword argument switch is the way to go. I did some work on the
>> existing SOS filtering implementation -- feel free to have a look at the
>> current pole-zero pairing function and suggest possible improvements if you
>> see them there, too.
>>
>> The feature I'd most like to see next related to SOS is actually
>> `filtfilt`-like support, so if you have ideas for that (or the time to do
>> it!) it would be great to make some progress toward that as well.
>>
>> Cheers,
>> Eric
>>
>>
>> On Wed, Jun 3, 2015 at 11:30 AM, Eric Quintero <ericq at caltech.edu> wrote:
>>
>>> Hi All,
>>>
>>> I’m very happy that second order section filtering is being introduced
>>> in scipy! Looking at the implementation, I see that it was chosen to use
>>> the direct form 2 realization of a second order section, which has the
>>> advantage of speed. Looking at the code, this maybe was also chosen to
>>> simply use signal.lfilter for the actual filtering.
>>>
>>> However, in the LIGO project, we initially used this form in our digital
>>> control systems, and found that this form leads to high levels of
>>> quantization noise in floating point signals, especially when the signals
>>> have a high dynamic range. This is due to the fact that the signal is first
>>> propagated through the filter poles, before going through the zeros.
>>>
>>> We’ve now moved on to using the Biquad form, which has one additional
>>> summation step compared to DF2, but avoids large internal values. The noise
>>> introduced by filtering in this way can be hundreds of times less than the
>>> DF2 results, which for our purposes makes the modest increase in
>>> computational time definitely worth it. Also, as a point of reference, the
>>> sosfilt function in MATLAB also uses this form.
>>>
>>> Long story short, I think it would be to scipy’s benefit to include
>>> biquad SOS filtering, at the very least as a keyword argument option, if
>>> not the default. I am, of course, willing to work on this myself.
>>>
>>> I appreciate any feedback you all may have, thanks for your
>>> consideration,
>>>
>>> Eric Q.
>>>
>>
> This seems like a fine plan.  I think that the way to go here is to make
> sosfilt at least internally into a gufunc.  Moving sosfilt to c or cython
> should also be a win for large data vectors or high order filters.  Since
> the current implementation calls lfilter in a loop, n_sections-1 extra
> input sized arrays are allocated and then destroyed during the call while
> the function shouldn't need any extra copies.
>
> Seems this is an all Eric discussion so far.
>
> -Eric
>

I forgot to mention that the slide deck at
https://dcc.ligo.org/LIGO-G0900928/public shows examples of the value of
the biquad form Eric discussed.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/scipy-dev/attachments/20150605/3e7a61f2/attachment.html>