Notch filter parts 5 and 6

6th-Semester-Spring-2024/DSP/Labs/NotchFilter/notch_filter_5.py

@@ -0,0 +1,84 @@
+import numpy as np
+import matplotlib.pyplot as plt
+import scipy.signal
+from  matplotlib import patches
+from matplotlib.figure import Figure
+from matplotlib import rcParams
+    
+def zplane(b,a,filename=None):
+    """Plot the complex z-plane given a transfer function.
+    """
+
+    # get a figure/plot
+    ax = plt.subplot(111)
+
+    # create the unit circle
+    uc = patches.Circle((0,0), radius=1, fill=False,
+                        color='black', ls='dashed')
+    ax.add_patch(uc)
+
+    # The coefficients are less than 1, normalize the coeficients
+    if np.max(b) > 1:
+        kn = np.max(b)
+        b = b/float(kn)
+    else:
+        kn = 1
+
+    if np.max(a) > 1:
+        kd = np.max(a)
+        a = a/float(kd)
+    else:
+        kd = 1
+        
+    # Get the poles and zeros
+    p = np.roots(a)
+    z = np.roots(b)
+    k = kn/float(kd)
+    
+    # Plot the zeros and set marker properties    
+    t1 = plt.plot(z.real, z.imag, 'o', ms=10, label="zero")
+    plt.setp(t1, markersize=10.0, color='tab:blue', fillstyle='none')
+
+    # Plot the poles and set marker properties
+    t2 = plt.plot(p.real, p.imag, 'x', ms=10, label="pole")
+    plt.setp(t2, markersize=10.0, color='tab:blue') 
+
+    ax.spines['left'].set_position('center')
+    ax.spines['bottom'].set_position('center')
+    ax.spines['right'].set_visible(False)
+    ax.spines['top'].set_visible(False)
+
+    # set the ticks
+    r = 1.5; plt.axis('scaled'); plt.axis([-r, r, -r, r])
+    ticks = [-1, -.5, .5, 1]; plt.xticks(ticks); plt.yticks(ticks)
+    return z, p, k
+
+
+f_s = 2         # Sampling frequency of 2[Hz]
+f_0 = 0.25      # Notch frequency at 0.25[Hz]
+r = 0.99        # Poles at radius 0.99
+
+b_0 = 1
+b_1 = -2*np.cos(2*np.pi*f_0/f_s)
+b_2 = 1
+b = [b_2, b_1, b_0]
+
+a_0 = 1
+a_1 = -2*r*np.cos(2*np.pi*f_0/f_s)
+a_2 = r**2
+a = [a_2, a_1, a_0]
+
+omega, h scipy.signal.freqz(b,a, fs=f_s)
+
+h_f = np.abs(h)
+
+plt.figure(figsize=(16,9))
+plt.plot(omega*f_s / (2*np.pi), h_f)
+plt.title("Magnitude Response")
+plt.xlabel("Frequency in Hz")
+plt.savefig("IIR_099r_2fs_025f0_freq.png")
+
+plt.figure(figsize=(16,9))
+zplane(b,a)
+plt.legend()
+plt.savefig("IIR_099r_2fs_025f0_pz.png")

6th-Semester-Spring-2024/DSP/Labs/NotchFilter/notch_filter_6.py

@@ -0,0 +1,77 @@
+import numpy as np
+import matplotlib.pyplot as plt
+import scipy.signal
+from  matplotlib import patches
+from matplotlib.figure import Figure
+from matplotlib import rcParams
+    
+def zplane(b,a,filename=None):
+    """Plot the complex z-plane given a transfer function.
+    """
+
+    # get a figure/plot
+    ax = plt.subplot(111)
+
+    # create the unit circle
+    uc = patches.Circle((0,0), radius=1, fill=False,
+                        color='black', ls='dashed')
+    ax.add_patch(uc)
+
+    # The coefficients are less than 1, normalize the coeficients
+    if np.max(b) > 1:
+        kn = np.max(b)
+        b = b/float(kn)
+    else:
+        kn = 1
+
+    if np.max(a) > 1:
+        kd = np.max(a)
+        a = a/float(kd)
+    else:
+        kd = 1
+        
+    # Get the poles and zeros
+    p = np.roots(a)
+    z = np.roots(b)
+    k = kn/float(kd)
+    
+    # Plot the zeros and set marker properties    
+    t1 = plt.plot(z.real, z.imag, 'o', ms=10, label="zero")
+    plt.setp(t1, markersize=10.0, color='tab:blue', fillstyle='none')
+
+    # Plot the poles and set marker properties
+    t2 = plt.plot(p.real, p.imag, 'x', ms=10, label="pole")
+    plt.setp(t2, markersize=10.0, color='tab:blue') 
+
+    ax.spines['left'].set_position('center')
+    ax.spines['bottom'].set_position('center')
+    ax.spines['right'].set_visible(False)
+    ax.spines['top'].set_visible(False)
+
+    # set the ticks
+    r = 1.5; plt.axis('scaled'); plt.axis([-r, r, -r, r])
+    ticks = [-1, -.5, .5, 1]; plt.xticks(ticks); plt.yticks(ticks)
+    return z, p, k
+
+
+f_s = 2         # Sampling frequency of 2[Hz]
+f_0 = 0.25      # Notch frequency at 0.25[Hz]
+r = 0.99        # Poles at radius 0.99
+
+b_0 = 1
+b_1 = -2*np.cos(2*np.pi*f_0/f_s)
+b_2 = 1
+b = [b_2, b_1, b_0]
+
+a_0 = 2 - 2*np.cos(2*np.pi*f_0/f_s)
+a = [a_0]
+
+omega, h = scipy.signal.freqz(b,a, fs=f_s)
+
+h_f = np.abs(h)
+
+plt.figure(figsize=(16,9))
+plt.plot(omega, h_f)
+plt.title("Magnitude Response")
+plt.xlabel("Frequency in Hz")
+plt.savefig("FIR_099r_2fs_025f0_freq.png")