[SciPy-User] scipy.signal.butter function different from Matlab

[gwenael.guillaume at aliceadsl.fr]

Thanks Ryan ,

Appreciate your reply.

Best regards

Gwenaël

Selon Ryan May <rmay31 at gmail.com>:

> On Thu, Nov 18, 2010 at 4:55 AM,  <gwenael.guillaume at aliceadsl.fr> wrote:
> > % Parameters:
> > N=5;
> > F_inf=88.38834764831843;
> > F_sup=111.3623397675424;
> > Fs=32768.00013421773;
> >
> > % Filter coefficients computing:
> > [z,p,k]=butter(N,[F_inf F_sup]/(Fs/2));
> >
> > % Result:
> > z=[1;1;1;1;1;-1;-1;-1;-1;-1;]
> > p=[0.999020109086358+0.021203989980732i;...
> >   0.999020109086358-0.021203989980732i;...
> >   0.99789313059316+0.020239448648803i;...
> >   0.99789313059316-0.020239448648803i;...
> >   0.99762168426512+0.018853164086747i;...
> >   0.99762168426512-0.018853164086747i;...
> >   0.998184731408311+0.017659802911797i;...
> >   0.998184731408311-0.017659802911797i;...
> >   0.999249126282689+0.017020854458606i;...
> >   0.999249126282689-0.017020854458606i;]
> > k=5.147424357763888e-014
>
> Just because MATLAB is giving you an answer doesn't mean it's not
> having the same problems. You're asking it to make a bandpass filter
> with a normalized width of 0.0014 using only 5 coefficients--that's
> too tight with too few coefficients and will *always* result in
> problems. That k value there is supposed to be the gain, which is
> almost 0 (140dB of loss!).  If you look at a plot of the transfer
> function you're generating in MATLAB, you'll see it's complete
> garbage. Using sig.buttord() and moving the stopband 0.0001 off of the
> passbands edges, 5db attenuation in the passband and 30db in the stop,
> I get an order of 26, which is >>5.
>
> The code for that is:
>
> Wn = np.array([F_inf F_sup])/(Fs/2)
> sig.buttord(Wn, np.array([Wn[0]-0.0001, Wn[1]+0.0001]), 5, 30)
>
> If we make the passband a bit less crazy:
>
> import matplotlib.pyplot as plt
> z,p,k = sig.butter(5, [0.3, 0.4], btype='bandpass', output='zpk')
> b,a = sig.zpk2tf(z, p, k)
> w,h = sig.freqz(b, a)
> plt.plot(w, 10.0 * np.log10(np.abs(h)))
> plt.grid()
> plt.ylim(-40, 3)
>
> You get a nice plot of the filter response, which is as expected for
> only order 5.
>
> As for the other things, it appears that there is no equivalent for SOS.
>
> Ryan
>
> --
> Ryan May
> Graduate Research Assistant
> School of Meteorology
> University of Oklahoma
> _______________________________________________
> SciPy-User mailing list
> SciPy-User at scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-user
>