## TP Filtre de Rauch - Simulation Python

from math import pi, sqrt
import numpy as np
import matplotlib.pyplot  as plt
from cmath import phase # calculer l'argument d'un nombre complexe

# I  filtre de Rauch 2 tracé de e(t) et s(t)

# Script 1
def H1(w):    
    return -1/(1+1j*Q*(w/wo-wo/w)) # wo=8514 rd/s ,fo=1355 Hz
wo=        # A COMPLETER
Q=         # A COMPLETER
e0=    #  amplitude de l'entrée  A COMPLETER
t=np.arange( , , ) # A COMPLETER
f=1300 # fréquence 
w=2*np.pi*f  #  pulsation
s0=e0*abs(H1(w)) # amplitude en sortie
phi=np.angle(H1(w)) # avance de phase de la sortie sur l'entrée
plt.figure(1)
plt.clf()
e=e0*np.sin(w*t)
s=s0*np.sin(      ) # signal en sortie A COMPLETER
plt.plot(t,e,'r' ,label="e(t)")
plt.plot(t,s,'b',label="s(t)")
plt.xlabel("temps (s)")
plt.ylabel("signaux e(t) et s(t)")
plt.legend()
plt.title("Filtre passe-bande de Rauch, fréquence : "+str(f)+" Hz")
plt.show()

## II filtre de Rauch diagramme de Bode

# Script 2  Filtre de Rauch  Bode en amplitude

from math import log10
wmin=                  # A COMPLETER 
wmax=                  # A COMPLETER 
w=np.linspace(   ,   ,   )   # A COMPLETER 
# Bode avec répartition régulière des w en échelle log
logw=np.linspace(log10(wmin),log10(wmax),100); w=10**logw

# calcul du gain
GdB=20*np.log10(         ) # A COMPLETER

# diagramme de gain
plt.figure(2)
plt.clf()
plt.title('Filtre de Rauch: diagramme de Bode d\'amplitude')
plt.plot(w,GdB,'-sb',lw=2)
plt.xlabel('pulsation (rad/s)')
plt.ylabel('Gain en tension  (dB)  ')
plt.xscale('log')
plt.grid()
#plt.show()

# script 3 :  Bode amplitude ajout du tracé asymptotique (GdB)

w_asymptote=[ , , ]# A COMPLETER 
GdB_asymptote=[ , , ] # A COMPLETER
plt.figure(2)
plt.plot(w_asymptote,GdB_asymptote,'r')
plt.show()

# script 4 : Bode tracé du diagramme de phase 
phi=np.angle(    )*180/pi  # A COMPLETER
plt.figure(3)
plt.clf()
plt.title('Filtre de Rauch : diagramme de phase')
plt.plot(w,phi,'-sb',lw=2)
plt.xlabel('pulsation (rad/s)')
plt.ylabel('avance de phase de s(t) sur e(t) en °')
plt.xscale('log')
plt.grid()
plt.show()


## III Analyse spectrale  
# Script 5 Tracé du spectre d'amplitude pour f= 452 Hz et n=5 harmoniques
def H2(w):     
    return -1/(1+1j*Q*(w/wo-wo/w))
wo=8514
f=452 
w= 2840          # 452*2*Pi
Q=20 
n=5 # nombre d'harmoniques considérés, 5 suffisent
A=np.zeros(n)
B=np.zeros(n)
for p in range (n):
    A[p]=        # tableau des amplitudes d'harmoniques A compléter
    B[p]=        # tableau des fréquences des harmoniques A compléter
    e0=A[p]
    s0=e0*abs(H2(B[p]))
    plt.figure(4)
    plt.xlabel('fréquence (Hz)')
    plt.ylabel('amplitude des harmoniques en Volt')
    plt.title("Filtre de Rauch : spectre fréquentiel pour f="+str(f)+"Hz et Q="+str(Q))
    plt.plot([B[p]/2/np.pi,B[p]/2/np.pi],[0,abs(s0)])
plt.show()