title: Python: IIR Filter response
body: <p>This is related to <a href="https://stackoverflow.com/questions/23341566/python-high-pass-fir-filter-by-windowing">this post</a> but is a question in its own right:
  Sampling at 40MHz, I have created a narrow pass band IIR filter which has a center frequency of 1MHz and a BW of 20kHz. This gives the following coefficients -</p>
  
  <pre><code>Fc = 1e6 /40e6 # Fcenter as a fraction of Fsample 
  BW = 20e3/40e6 # BW as a fraction of Fsample 
  
  a0 = 0.00140
  a2 = 0.00018
  
  b1 = 1.97241
  b2 = -0.9970
  </code></pre>
  
  <p>Applying the signal appears to work as indicated below - My question is: 
  How may I plot the magnitude and phase response in python?
  [As per my previous post this is intentionally idiomatic] </p>
  
  <p><img src="https://i.stack.imgur.com/g8dWu.png" alt="enter image description here"></p>
  
  <pre><code>import numpy as np
  import matplotlib.pyplot as plt
  
  # create an array of 1024 points sampled at 40MHz
  # [each sample is 25ns apart and the key signal is 1MHz]
  Fs = 40e6
  T  = 1/Fs
  t  = np.arange(0,(1024*T),T)
  
  f     = 1e6   
  Omega = 2*np.pi*f
  x  = np.sin(Omega*t) * (t**3) * np.exp(-t/2e-6)
  x /= max(x)
  y  = [0]*len(x)
  
  # create a narrow passband IIR filter with fcentre=1MHz
  # and BW=0.0005
  Fc = 1e6
  Ft = Fc/Fs
  BW = 0.0005
  R  = 1 - (3*BW)
  K  = (1 - 2*R*np.cos(2*np.pi*Ft) + (R*R)) / (2 - 2*np.cos(2*np.pi*Ft))
  
  # coefficients
  a0 = 1 - K
  a1 = 2*(K-R)*np.cos(2*np.pi*Ft)
  a2 = (R*R) - K
  b1 = 2*R*np.cos(2*np.pi*Ft)
  b2 = -(R*R)
  
  for n in range(2, len(x)):
      y[n] = a0*x[n] + a1*x[n-1] + a2*x[n-2] + b1*y[n-1] + b2*y[n-2]
  y /= max(y)
  
  plt.subplot(211)
  plt.plot( x,'r-', linewidth=2)
  plt.xlabel( 'sample length' )
  plt.ylabel( 'ip value' )
  plt.grid()
  
  plt.subplot(212)
  plt.plot( y,'k-', linewidth=2)
  plt.xlabel( 'sample length' )
  plt.ylabel( 'op value' )
  plt.grid()
  plt.show()
  </code></pre>