[SciPy-Dev] signal.sosfilt filter representation

[Eric Moore]

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
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/scipy-dev/attachments/20150605/d6136cab/attachment.html>