I am trying to design a bandpass butterworth filter and am getting an error import scipy.signal as sig Fs = 400. Nyq = Fs/2. # low pass stop freq, high pass corner freq, etc... lpsf = 5. lpcf = 7. hpcf = 12. hpsf = 18. Nyq = Fs/2. wp = [lpcf/Nyq, hpcf/Nyq] ws = [lpsf/Nyq, hpsf/Nyq] gpass = 3. gstop = 15. ord, Wn = sig.buttord(wp, ws, gpass, gstop) mybutt = sig.butter(ord, Wn, btype='bandpass') # pun intended With scipy version 0.3.2 I get the following the traceback Traceback (most recent call last): File "<console>", line 1, in ? File "/usr/lib/python2.4/site-packages/scipy/signal/filter_design.py", line 482, in butter return iirfilter(N, Wn, btype=btype, analog=analog, output=output, ftype='butter') File "/usr/lib/python2.4/site-packages/scipy/signal/filter_design.py", line 456, in iirfilter b, a = lp2bp(b,a,wo=wo,bw=bw) File "/usr/lib/python2.4/site-packages/scipy/signal/filter_design.py", line 253, in lp2bp aprime[Dp-j] = val TypeError: can't convert complex to float; use abs(z) Any ideas? Thanks, JDH
John Hunter wrote: >I am trying to design a bandpass butterworth filter and am getting an >error > > import scipy.signal as sig > Fs = 400. > Nyq = Fs/2. > # low pass stop freq, high pass corner freq, etc... > lpsf = 5. > lpcf = 7. > hpcf = 12. > hpsf = 18. > Nyq = Fs/2. > wp = [lpcf/Nyq, hpcf/Nyq] > ws = [lpsf/Nyq, hpsf/Nyq] > gpass = 3. > gstop = 15. > ord, Wn = sig.buttord(wp, ws, gpass, gstop) > mybutt = sig.butter(ord, Wn, btype='bandpass') # pun intended > > > > >With scipy version 0.3.2 I get the following the traceback > >Traceback (most recent call last): > File "<console>", line 1, in ? > File "/usr/lib/python2.4/site-packages/scipy/signal/filter_design.py", line 482, in butter > return iirfilter(N, Wn, btype=btype, analog=analog, output=output, ftype='butter') > File "/usr/lib/python2.4/site-packages/scipy/signal/filter_design.py", line 456, in iirfilter > b, a = lp2bp(b,a,wo=wo,bw=bw) > File "/usr/lib/python2.4/site-packages/scipy/signal/filter_design.py", line 253, in lp2bp > aprime[Dp-j] = val >TypeError: can't convert complex to float; use abs(z) > > val must be complex for some reason where it is not expected to be. Not sure what is going on. Which version of Numeric do you have installed? This is what I get. >>> import scipy >>> scipy.__version__ '0.3.3_304.4613' import scipy.signal as sig ord,Wn = sig.buttord([7/200.,12/200.],[5/200.0,18/200.0],3,15) >>> ord 2 >>> Wn array([ 0.034979417860466, 0.060035166350418]) >>> mybutt = sig.butter(ord,Wn,btype='bandpass') >>> mybutt (array([ 0.001466661479123, 0. , -0.002933322958246, 0. , 0.001466661479123]), array([ 1. , -3.848535133256346, 5.594309123067282, -3.640021384921099, 0.894652727053865])) -Travis
"Travis" == Travis Oliphant <oliphant@ee.byu.edu> writes:
Travis> val must be complex for some reason where it is not Travis> expected to be. Travis> Not sure what is going on. Which version of Numeric do Travis> you have installed? 23.7 Here is what I get for the same commands in a plain-ol-python-shell.
uname -a Linux peds-pc311 2.6.10-5-386 #1 Fri Jun 24 16:53:01 UTC 2005 i686 GNU/Linux
python Python 2.4.1 (#2, Mar 30 2005, 21:51:10) [GCC 3.3.5 (Debian 1:3.3.5-8ubuntu2)] on linux2 Type "help", "copyright", "credits" or "license" for more information.
import Numeric as N N.__version__ '23.7' import scipy as s s.__version__ '0.3.2' import scipy.signal as sig ord,Wn = sig.buttord([7/200.,12/200.],[5/200.0,18/200.0],3,15) ord 2 Wn array([ 0.03497942, 0.06003517]) mybutt = sig.butter(ord,Wn,btype='bandpass') Traceback (most recent call last): File "<stdin>", line 1, in ? File "/usr/lib/python2.4/site-packages/scipy/signal/filter_design.py", line 482, in butter return iirfilter(N, Wn, btype=btype, analog=analog, output=output, ftype='butter') File "/usr/lib/python2.4/site-packages/scipy/signal/filter_design.py", line 456, in iirfilter b, a = lp2bp(b,a,wo=wo,bw=bw) File "/usr/lib/python2.4/site-packages/scipy/signal/filter_design.py", line 253, in lp2bp aprime[Dp-j] = val TypeError: can't convert complex to float; use abs(z)
Let me know if there is anything else you want me to try. I can always upgrade to scipy CVS and try that. Thanks! JDH
John Hunter wrote:
import scipy as s s.__version__
'0.3.2'
[...]
Let me know if there is anything else you want me to try. I can always upgrade to scipy CVS and try that.
By all means do. 0.3.2 is quite out of date, so there's a good chance CVS will take care of things. Cheers, f
"Fernando" == Fernando Perez <Fernando.Perez@colorado.edu> writes:
Fernando> By all means do. 0.3.2 is quite out of date, so there's Fernando> a good chance CVS will take care of things. That did the trick. Because I need this fix for software that I distribute to windoze users, I still have some interest in seeing if I can backport the fix the 0.3.2, since this is the current binary release and thus I suspect it will also be the release in the new enthought python. I'm posting the diff betwween scipy 0.3.3_309.4626 and 0.3.2 in the signal module, in case the person(s) who made those changes will have some insight into what is contributing to the fix. If not, it's no big deal because I can always try and build scipy CVS on windows (arggg) or just tell folks they can't use this feature. It would be much easier just to say, "drop this file in site-packages/scipy/signal" if it's just a one or two line fix. Thanks! JDH Only in scipy/Lib/signal/: CVS Only in scipy/Lib/signal/: docs diff -cr SciPy_complete-0.3.2/Lib/signal/filter_design.py scipy/Lib/signal/filter_design.py *** SciPy_complete-0.3.2/Lib/signal/filter_design.py 2004-02-26 02:55:12.000000000 -0600 --- scipy/Lib/signal/filter_design.py 2004-11-22 14:25:00.000000000 -0600 *************** *** 5,10 **** --- 5,11 ---- from scipy_base.fastumath import * from scipy_base import atleast_1d, poly, polyval, roots, imag, real, asarray,\ allclose, Float, resize, pi, concatenate, absolute, logspace, c_ + from scipy_base import mintypecode from scipy import comb, special, optimize, linalg import string, types *************** *** 228,236 **** a,b = map(atleast_1d,(a,b)) D = len(a) - 1 N = len(b) - 1 ! artype = b.typecode() ! if artype not in ['F','D','f','d']: ! artype = Num.Float ma = max([N,D]) Np = N + ma Dp = D + ma --- 229,235 ---- a,b = map(atleast_1d,(a,b)) D = len(a) - 1 N = len(b) - 1 ! artype = mintypecode((a,b)) ma = max([N,D]) Np = N + ma Dp = D + ma *************** *** 261,269 **** a,b = map(atleast_1d,(a,b)) D = len(a) - 1 N = len(b) - 1 ! artype = b.typecode() ! if artype not in ['F','D','f','d']: ! artype = Num.Float M = max([N,D]) Np = M + M Dp = M + M --- 260,266 ---- a,b = map(atleast_1d,(a,b)) D = len(a) - 1 N = len(b) - 1 ! artype = mintypecode((a,b)) M = max([N,D]) Np = M + M Dp = M + M diff -cr SciPy_complete-0.3.2/Lib/signal/firfilter.c scipy/Lib/signal/firfilter.c *** SciPy_complete-0.3.2/Lib/signal/firfilter.c 2003-07-08 20:52:20.000000000 -0500 --- scipy/Lib/signal/firfilter.c 2005-05-31 01:51:32.000000000 -0500 *************** *** 112,118 **** if ((sum = calloc(type_size,2))==NULL) return -3; /* No memory */ value = sum + type_size; ! if (outsize == FULL) {Os[0] = Ns[0]+Nwin[0]-1; Os[1] = Ns[0]+Nwin[1]-1;} else if (outsize == SAME) {Os[0] = Ns[0]; Os[1] = Ns[1];} else if (outsize == VALID) {Os[0] = Ns[0]-Nwin[0]+1; Os[1] = Ns[1]-Nwin[1]+1;} else return -1; /* Invalid output flag */ --- 112,118 ---- if ((sum = calloc(type_size,2))==NULL) return -3; /* No memory */ value = sum + type_size; ! if (outsize == FULL) {Os[0] = Ns[0]+Nwin[0]-1; Os[1] = Ns[1]+Nwin[1]-1;} else if (outsize == SAME) {Os[0] = Ns[0]; Os[1] = Ns[1];} else if (outsize == VALID) {Os[0] = Ns[0]-Nwin[0]+1; Os[1] = Ns[1]-Nwin[1]+1;} else return -1; /* Invalid output flag */ *************** *** 120,125 **** --- 120,128 ---- if ((boundary != PAD) && (boundary != REFLECT) && (boundary != CIRCULAR)) return -2; /* Invalid boundary flag */ + /* Speed this up by not doing any if statements in the for loop. Need 3*3*2=18 different + loops executed for different conditions */ + for (m=0; m < Os[0]; m++) { /* Reposition index into input image based on requested output size */ if (outsize == FULL) new_m = convolve ? m : (m-Nwin[0]+1); diff -cr SciPy_complete-0.3.2/Lib/signal/info_signal.py scipy/Lib/signal/info_signal.py *** SciPy_complete-0.3.2/Lib/signal/info_signal.py 2004-02-26 01:15:10.000000000 -0600 --- scipy/Lib/signal/info_signal.py 2005-05-31 01:51:32.000000000 -0500 *************** *** 6,15 **** --- 6,17 ---- convolve -- N-dimensional convolution. correlate -- N-dimensional correlation. + fftconvolve -- N-dimensional convolution using the FFT. convolve2d -- 2-dimensional convolution (more options). correlate2d -- 2-dimensional correlation (more options). sepfir2d -- Convolve with a 2-D separable FIR filter. + B-splines: bspline -- B-spline basis function of order n. *************** *** 45,52 **** freqs --- Analog filter frequency response freqz --- Digital filter frequency response ! ! matlab-style IIR filter design: butter (buttord) -- Butterworth cheby1 (cheb1ord) -- Chebyshev Type I --- 47,53 ---- freqs --- Analog filter frequency response freqz --- Digital filter frequency response ! Matlab-style IIR filter design: butter (buttord) -- Butterworth cheby1 (cheb1ord) -- Chebyshev Type I *************** *** 68,74 **** tf2ss -- transfer function to state-space. ss2tf -- state-pace to transfer function. zpk2ss -- zero-pole-gain to state-space. ! ss2zpk -- state-space to pole-zero-gain. """ postpone_import = 1 --- 69,88 ---- tf2ss -- transfer function to state-space. ss2tf -- state-pace to transfer function. zpk2ss -- zero-pole-gain to state-space. ! ss2zpk -- state-space to pole-zero-gain. ! ! Waveforms: ! ! sawtooth -- Periodic sawtooth ! square -- Square wave ! gausspulse -- Gaussian modulated sinusoid ! chirp -- Frequency swept cosine signal ! ! Wavelets: ! ! daub -- return low-pass filter for daubechies wavelets ! qmf -- return quadrature mirror filter from low-pass ! cascade -- compute scaling function and wavelet from coefficients """ postpone_import = 1 diff -cr SciPy_complete-0.3.2/Lib/signal/__init__.py scipy/Lib/signal/__init__.py *** SciPy_complete-0.3.2/Lib/signal/__init__.py 2004-04-13 16:23:51.000000000 -0500 --- scipy/Lib/signal/__init__.py 2005-05-31 01:51:32.000000000 -0500 *************** *** 10,13 **** --- 10,14 ---- from filter_design import * from ltisys import * from signaltools import * + from wavelets import * diff -cr SciPy_complete-0.3.2/Lib/signal/ltisys.py scipy/Lib/signal/ltisys.py *** SciPy_complete-0.3.2/Lib/signal/ltisys.py 2003-12-13 07:44:51.000000000 -0600 --- scipy/Lib/signal/ltisys.py 2005-02-24 04:42:24.000000000 -0600 *************** *** 223,230 **** self.__dict__['zeros'], self.__dict__['poles'], \ self.__dict__['gain'] = ss2zpk(*args) self.__dict__['num'], self.__dict__['den'] = ss2tf(*args) ! self.inputs = B.shape[-1] ! self.outputs = C.shape[0] else: raise ValueError, "Needs 2, 3, or 4 arguments." --- 223,230 ---- self.__dict__['zeros'], self.__dict__['poles'], \ self.__dict__['gain'] = ss2zpk(*args) self.__dict__['num'], self.__dict__['den'] = ss2tf(*args) ! self.inputs = self.B.shape[-1] ! self.outputs = self.C.shape[0] else: raise ValueError, "Needs 2, 3, or 4 arguments." diff -cr SciPy_complete-0.3.2/Lib/signal/setup_signal.py scipy/Lib/signal/setup_signal.py *** SciPy_complete-0.3.2/Lib/signal/setup_signal.py 2004-02-06 02:17:52.000000000 -0600 --- scipy/Lib/signal/setup_signal.py 2005-03-24 04:22:16.000000000 -0600 *************** *** 5,10 **** --- 5,11 ---- from scipy_distutils.misc_util import get_path, default_config_dict, dot_join def configuration(parent_package='',parent_path=None): + from scipy_distutils.system_info import get_info package = 'signal' local_path = get_path(__name__,parent_path) config = default_config_dict(package, parent_package) diff -cr SciPy_complete-0.3.2/Lib/signal/signaltools.py scipy/Lib/signal/signaltools.py *** SciPy_complete-0.3.2/Lib/signal/signaltools.py 2004-06-02 17:57:55.000000000 -0500 --- scipy/Lib/signal/signaltools.py 2005-01-17 18:53:50.000000000 -0600 *************** *** 10,17 **** import types import scipy from scipy.stats import mean ! import Numeric ! from Numeric import array, arange, where, sqrt, rank, zeros, NewAxis, argmax from scipy_base.fastumath import * _modedict = {'valid':0, 'same':1, 'full':2} --- 10,17 ---- import types import scipy from scipy.stats import mean ! import scipy_base as Numeric ! from scipy_base import array, arange, where, sqrt, rank, zeros, NewAxis, argmax, product from scipy_base.fastumath import * _modedict = {'valid':0, 'same':1, 'full':2} *************** *** 66,72 **** kernel = asarray(in2) if rank(volume) == rank(kernel) == 0: return volume*kernel ! if (Numeric.product(kernel.shape) > Numeric.product(volume.shape)): temp = kernel kernel = volume volume = temp --- 66,72 ---- kernel = asarray(in2) if rank(volume) == rank(kernel) == 0: return volume*kernel ! if (product(kernel.shape) > product(volume.shape)): temp = kernel kernel = volume volume = temp *************** *** 76,81 **** --- 76,117 ---- return sigtools._correlateND(volume, kernel, val) + def _centered(arr, newsize): + # Return the center newsize portion of the array. + newsize = asarray(newsize) + currsize = array(arr.shape) + startind = (currsize - newsize) / 2 + endind = startind + newsize + myslice = [slice(startind[k], endind[k]) for k in range(len(endind))] + return arr[tuple(myslice)] + + def fftconvolve(in1, in2, mode="full"): + """Convolve two N-dimensional arrays using FFT. SEE convolve + """ + s1 = array(in1.shape) + s2 = array(in2.shape) + if (s1.typecode() in ['D','F']) or (s2.typecode() in ['D', 'F']): + cmplx=1 + else: cmplx=0 + size = s1+s2-1 + IN1 = fftn(in1,size) + IN1 *= fftn(in2,size) + ret = ifftn(IN1) + del IN1 + if not cmplx: + ret = real(ret) + if mode == "full": + return ret + elif mode == "same": + if product(s1) > product(s2): + osize = s1 + else: + osize = s2 + return _centered(ret,osize) + elif mode == "valid": + return _centered(ret,abs(s2-s1)+1) + + def convolve(in1, in2, mode='full'): """Convolve two N-dimensional arrays. Only in scipy/Lib/signal/tests: CVS Only in scipy/Lib/signal/: wavelets.py
participants (3)
-
Fernando Perez
-
John Hunter
-
Travis Oliphant