
import numpy as np
import matplotlib.pyplot as plt
import pycanum.main as pycan
from scipy.signal import blackman
from numpy.fft import fft

def frequence(t,x):
    N=len(x)
    te=t[1]-t[0]
    zeros=np.zeros(6*N)
    X = np.concatenate((x*blackman(N),zeros))
    NN=len(X)
    spectre = np.absolute(fft(X))*2.0/N/0.42
    freq = np.arange(NN)*1.0/(NN*te)
    k=np.argmax(spectre[0:NN//2])
    return freq[k]

def mesure(can,R,techant,N=50000,Np=100,plot=False):
    can.config_echantillon(techant*1e6,N)
    can.acquerir()
    t=can.temps()[0]
    signaux = can.entrees()
    u=signaux[0]
    ur=signaux[1]
    u=u-u.mean()
    ur=ur-ur.mean()
    i=ur/R
    N = len(t)
    Z = u.std()/i.std()
    cosphi = np.mean(i*u)/(u.std()*i.std())
    phi = np.arccos(cosphi)
    freq = frequence(t,u)
    if plot:
        plt.figure()
        plt.plot(t,u,'b')
        plt.plot(t,R*i,'r')
        plt.grid()
        plt.ylim(-10,10)
        plt.xlim(0,20/freq)
        plt.show()
    return (freq,Z,cosphi,phi)
   
liste_f = []
liste_Z = []
liste_phi = []
R=100.7 # résistance en série
can = pycan.Sysam("SP5")
Vmax = 2
techant = 1e-6
can.config_entrees([0,1],[Vmax,Vmax],diff=[0])
r = 'O'
while r=='O':
    (f,Z,cosphi,phi) = mesure(can,R,techant)
    print("f = %f Hz, Z = %f, cosphi = %f, phi = %f"%(f,Z,cosphi,phi))
    liste_f.append(f)
    liste_Z.append(Z)
    liste_phi.append(phi)
    r = input("mesure (O/N) ?")
can.fermer()
liste_Z = np.array(liste_Z)
liste_phi = np.unwrap(liste_phi)
np.savetxt('impedance-E=2.txt',np.array([liste_f,liste_Z,liste_phi]).T,header='f\t Z\t phi')
plt.figure()
plt.plot(liste_f,liste_Z,'b-')
plt.grid()
plt.xlabel('f (Hz)')
plt.ylabel(r'$|Z|\ (\rm\Omega)$')
plt.figure()
plt.plot(liste_f,liste_phi/np.pi,'b-')
plt.xlabel('f (Hz)')
plt.ylabel(r'$\phi/\pi (\rm rad)$')
plt.grid()
plt.figure()
plt.plot(liste_f,liste_Z*np.cos(liste_phi),'b-')
plt.grid()
plt.xlabel('f (Hz)')
plt.ylabel(r'$Re(Z)\ (\rm\Omega)$')
plt.figure()
plt.plot(liste_f,liste_Z*np.sin(liste_phi),'b-')
plt.xlabel('f (Hz)')
plt.ylabel(r'$Im(Z)\ (\rm\Omega)$')
plt.grid()
plt.show()
                     