# Design a 15th-order bandpass filter in SOS format.

from scipy import signal
sos = signal.ellip(15, 0.5, 60, (0.2, 0.4), btype='bandpass',
                   output='sos')

# Compute the frequency response at 1500 points from DC to Nyquist.

w, h = signal.sosfreqz(sos, worN=1500)

# Plot the response.

import matplotlib.pyplot as plt
plt.subplot(2, 1, 1)
db = 20*np.log10(np.abs(h))
plt.plot(w/np.pi, db)
plt.ylim(-75, 5)
plt.grid(True)
plt.yticks([0, -20, -40, -60])
plt.ylabel('Gain [dB]')
plt.title('Frequency Response')
plt.subplot(2, 1, 2)
plt.plot(w/np.pi, np.angle(h))
plt.grid(True)
plt.yticks([-np.pi, -0.5*np.pi, 0, 0.5*np.pi, np.pi],
           [r'$-\pi$', r'$-\pi/2$', '0', r'$\pi/2$', r'$\pi$'])
plt.ylabel('Phase [rad]')
plt.xlabel('Normalized frequency (1.0 = Nyquist)')
plt.show()

# If the same filter is implemented as a single transfer function,
# numerical error corrupts the frequency response:

b, a = signal.ellip(15, 0.5, 60, (0.2, 0.4), btype='bandpass',
                   output='ba')
w, h = signal.freqz(b, a, worN=1500)
plt.subplot(2, 1, 1)
db = 20*np.log10(np.abs(h))
plt.plot(w/np.pi, db)
plt.subplot(2, 1, 2)
plt.plot(w/np.pi, np.angle(h))
plt.show()
