x = np.linspace(0,np.pi/2,100)
plt.plot(x,np.sin(x),label='sin(x)')
plt.plot(x,x,label='x')
plt.plot(x,x-x**3/6,label='x-x**3/6')
plt.xlabel('Angle [rad]')
plt.legend()
plt.show()

round(conv_rad_deg(0.5))

x = np.linspace(0,np.pi/2,100)
plt.plot(x,np.cos(x),label='cos(x)')
plt.plot(x,1-x**2/2,label='1-x²/2')
plt.plot(x,1-x**2/2+x**4/24,label='1-x²/2+x**4/24')
plt.xlabel('Angle [rad]')
plt.legend()
plt.show()

T_nonlin = np.array([])
T_nonlin = T_nonlin[:len(T_nonlin)-10]
#Sans masse
T_nonlin = np.array([])
#Avec masse
T_nonlin = np.array([])
T_nonlin = T_nonlin[::2]

n = np.linspace(0,len(T_nonlin)-1,len(T_nonlin))

print(n)

popc, pcov = curve_fit(affine, n[:10], T_nonlin[:10])

a,b = popc

plt.plot(T_nonlin,'x') # Permet de montrer la baisse de la période avec le temps dû aux frottements
plt.plot(n, a*n + b)
plt.xlabel("Numéro d'interruption")
plt.ylabel("Temps d'interruption (en s)")

print(n[:3])

def affine(x,a,b):
    return a*x + b

theta = np.array([10, 20, 30, 40, 50, 60, 70, 80, 90, 100, 110])*(np.pi/180)
utheta = 1*(np.pi/180)
theta2 = theta**2

T10 = np.array([])
T10 = T10[::2]

T20 = np.array([])
T20 = T20[::2]

T30 = np.array([])
T30 = T30[::2]

T40 = np.array([])
T40 = T40[::2]

T50 = np.array([])
T50 = T50[::2]

T60 = np.array([])
T60 = T60[::2]

T70 = np.array([])
T70 = T70[::2]

T80 = np.array([])
T80 = T80[::2]

T90 = np.array([])
T90 = T90[::2]

T100 = np.array([])
T100 = T100[::2]

T110 = np.array([])
T110 = T110[::2]

lin = 5
T10 = T10[:lin]
T20 = T20[:lin]
T30 = T30[:lin]
T40 = T40[:lin]
T50 = T50[:lin]
T60 = T60[:lin]
#T65 = T65[:3]
T70 = T70[:lin]
T80 = T80[:lin]
T90 = T90[:lin]
T100 = T100[:lin]
T110 = T110[:lin]
#T120 = T120[:3]
#T130 = T130[:3]
#T140 = T140[:3]

T = np.zeros(len(theta))

n = np.linspace(0,len(T10) - 1,len(T10))

T_list = [T10, T20, T30, T40, T50, T60, T70, T80, T90, T100, T110]


# On récupère le coefficient directeur de la droite T=f(n) qui va correspondre à la période
for i in range(len(theta)):
    popc, pcov = curve_fit(affine,n,T_list[i])
    T[i] = popc[0]*2
print(T)

def func(x,a,b,c):
    return a*(1 + b*x + c*x**2)

popc,pcov = curve_fit(func, theta2, T, p0=[0.3, 6.25e-2, 3.6e-3])

a,b,c = popc
ub = np.sqrt(pcov[1,1])
uc = np.sqrt(pcov[2,2])

plt.plot(theta2, T, 'x', label='data')
plt.plot(theta2, func(theta2, *popc), label='fit')
plt.xlabel('theta**2')
plt.ylabel('Période T')
plt.legend()
plt.show()

print("La valeur théorique de b est 6.25e-2")
print("La valeur expérimentale de b est (", np.round(b*1e2,1), '+-', np.round(ub*1e2,1), ') e-2')
#print("La valeur théorique de c est 3.6e-3")
#print("La valeur expérimentale de c est (", np.round(c*1e3,1), '+-', np.round(uc*1e3,1), ') e-3')

data = pd.read_csv('gyro_.csv', delimiter = ',', decimal='.')

t = data['Time (s)'].values
vwz = data['Gyroscope z (rad/s)'].values # Données du gyroscope, vitesse angulaire en z

d1 = 200
f1 = 1400
t_fft_r = t[d1:f1]
vwz_fft_r = vwz[d1:f1]

d2 = 3200
f2 = 4500
t_fft_g = t[d2:f2]
vwz_fft_g = vwz[d2:f2]

d3 = 8000
f3 = 9200
t_fft_k = t[d3:f3]
vwz_fft_k = vwz[d3:f3]

plt.subplot(2,1,1)
plt.plot(t,vwz)
plt.plot(t_fft_r, vwz_fft_r, 'r-')
plt.plot(t_fft_g, vwz_fft_g, 'g-')
plt.plot(t_fft_k, vwz_fft_k, 'k-')
plt.title("Evolution des oscillations au cours du temps (theta = 90°)")
plt.xlabel('Temps (en s)')
plt.ylabel('Vitesse angulaire (en rad/s)')

# Number of sample points
N = len(t_fft_r)
# sample spacing
T = t[1] - t[0]
#yf_r = sp.fft.fft(vwz_fft_r)
#yf_g = sp.fft.fft(vwz_fft_g)
#yf_k = sp.fft.fft(vwz_fft_k)
#xf = sp.fft.fftfreq(N, T)[:N//2]

# si bug avec fft de scipy
from scipy.fftpack import fft
from scipy.fftpack import fftfreq

yf_r = fft(vwz_fft_r)
yf_g = fft(vwz_fft_g)
yf_k = fft(vwz_fft_k)
xf = fftfreq(N, T)[:N//2]

plt.subplot(2,1,2)
plt.plot(xf, 2.0/N * np.abs(yf_r[0:N//2]), 'r')
plt.plot(xf, 2.0/N * np.abs(yf_g[0:N//2]), 'g')
plt.plot(xf, 2.0/N * np.abs(yf_k[0:N//2]), 'k')
plt.xlim(0,3)
plt.xlabel('Fréquence (en Hz)')
plt.ylabel('Amplitude')