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

Ryan May rmay31 at gmail.com
Thu Nov 18 22:53:08 EST 2010


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



More information about the SciPy-User mailing list