Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

SciPy "lfilter" returns only NaNs

All -

I am trying to use SciPy's signal.lfilter function to filter a vector of samples - unfortunately, all that is returned is a vector of NaN.

I have plotted the frequency response of the filter, and the filter coefficients look correct; I'm fairly certain the issue is with the actual call to lfilter.

It's a high-pass Chebychev I filter, which I'm creating with:

b,a = signal.iirdesign(wp = 0.11, ws= 0.1, gstop= 60, gpass=1, ftype='cheby1')

I am then filtering data with:

filtered_data = signal.lfilter(b, a, data)

Below, I am printing a selection of 20 samples from the pre-filtered data, and then the filtered data. You can clearly see the issue:

### Printing a small selection of the data before it is filtered:

((-0.003070347011089325+0.0073614344000816345j), (-0.003162827342748642+0.007342938333749771j), (-0.003310795873403549+0.0073614344000816345j), (-0.0031813234090805054+0.007342938333749771j), (-0.003255307674407959+0.007398426532745361j), (-0.003162827342748642+0.007287450134754181j), (-0.003125835210084915+0.007509402930736542j), (-0.003162827342748642+0.007342938333749771j), (-0.0031073391437530518+0.007287450134754181j), (-0.0032368116080760956+0.007398426532745361j), (-0.0030888430774211884+0.007342938333749771j))


### Printing a small selection of the filtered data:

[ nan nanj  nan nanj  nan nanj  nan nanj  nan nanj  nan nanj  nan nanj
  nan nanj  nan nanj  nan nanj  nan nanj  nan nanj  nan nanj  nan nanj
  nan nanj  nan nanj  nan nanj  nan nanj  nan nanj  nan nanj]

Like I said before, the coefficients of the filter look good. They are:

b = [  4.06886235e-02  -7.73083846e-01   6.95775461e+00  -3.94272761e+01
   1.57709105e+02  -4.73127314e+02   1.10396373e+03  -2.05021836e+03
   3.07532754e+03  -3.75873366e+03   3.75873366e+03  -3.07532754e+03
   2.05021836e+03  -1.10396373e+03   4.73127314e+02  -1.57709105e+02
   3.94272761e+01  -6.95775461e+00   7.73083846e-01  -4.06886235e-02]
a = [  1.00000000e+00  -1.27730099e+01   7.81201390e+01  -3.03738394e+02
   8.40827723e+02  -1.75902089e+03   2.88045462e+03  -3.77173152e+03
   3.99609428e+03  -3.43732844e+03   2.38415171e+03  -1.30118368e+03
   5.21654119e+02  -1.18026566e+02  -1.85597824e+01   3.24205235e+01
  -1.65545917e+01   5.02665439e+00  -9.09697811e-01   7.68172820e-02]

So why would lfilter return only NaN? How am I using this function incorrectly?

Thanks in advance for your help!

Edit:

Okay, I solved it.

For anyone that encounters this in the future:

For whatever reason, even though the returned coefficients for the filter looked good, when I then used those coefficients in SciPy's lfilter function, the filtered values were unbounded. Simply changing the passband edge to ANY number other than 0.11 fixed the problem. Even this works:

b,a = signal.iirdesign(wp = 0.119, ws= 0.1, gstop= 60, gpass=1, ftype='cheby1')

Other than manually grepping through the poles and zeros of the filter, I'm not sure how you would detect instability of the filter. Bizarre.

like image 208
bhilburn Avatar asked Jan 10 '12 22:01

bhilburn


1 Answers

An IIR filter is stable if the absolute values of the roots of the denominator of the discrete transfer function a(z) are all less than one. So, you can detect the instability by following code:

from scipy import signal
import numpy as np
b1, a1 = signal.iirdesign(wp = 0.11, ws= 0.1, gstop= 60, gpass=1, ftype='cheby1')
b2, a2 = signal.iirdesign(wp = 0.119, ws= 0.1, gstop= 60, gpass=1, ftype='cheby1')

print "filter1", np.all(np.abs(np.roots(a1))<1)
print "filter2", np.all(np.abs(np.roots(a2))<1)
like image 138
HYRY Avatar answered Sep 25 '22 09:09

HYRY