Discussion:
[SciPy-User] how to treat an invalid value, in signal/filter_design.py
R Schumacher
2016-10-21 19:34:45 UTC
Permalink
In an attempt to computationally invert the effect of an analog RC
filter on a data set and reconstruct the signal prior to the analog
front end, a co-worker suggested: "Mathematically, you just reverse
the a and b parameters. Then the zeros become the poles, but if the
new poles are not inside the unit circle, the filter is not stable."

So then to "stabilize" the poles' issue seen, I test for the DIV/0
error and set it to 2./N+0.j in scipy/signal/filter_design.py ~ line 244
d = polyval(a[::-1], zm1)
if d[0]==0.0+0.j:
d[0] = 2./N+0.j
h = polyval(b[::-1], zm1) / d

- Question is, is this a mathematically valid treatment?
- Is there a better way to invert a Butterworth filter, or work with
the DIV/0 that occurs without modifying the signal library?

I noted d[0] > 2./N+0.j makes the zero bin result spike low; 2/N
gives a reasonable "extension" of the response curve.
The process in general causes a near-zero offset however, which I
remove with a high pass now; In an full FFT of a ~megasample one can
see that the first 5 bins have run away.

An example attached...


Ray Schumacher
Programmer/Consultant

Loading...