Rotolamento cilindro

# Libreria Santanastasio

import matplotlib.pyplot as plt
import numpy as np
from math import log10, floor, sqrt
import scipy as sp
from scipy import stats
from mpl_toolkits.axes_grid1.inset_locator import zoomed_inset_axes
from mpl_toolkits.axes_grid1.inset_locator import mark_inset
#########################################################################
# funzione per arrotondare con un certo numero di cifre significative
#########################################################################
def PrintResult(name,mean,sigma,digits,unit):
    mean = round(mean,digits)
    sigma = round(sigma,digits)
    nu = sigma / mean
    result = (name+" = ({0} +/- {1} ) ".format(mean,sigma)+unit+" [{0:.2f}%]".format(nu*100))
    print (result)
    #return ""

#########################################################################
# funzioni per fare il fit lineare
#########################################################################
#
def my_mean(x, w):
    return np.sum( x*w ) / np.sum( w )

def my_cov(x, y, w):
    return my_mean(x*y, w) - my_mean(x, w)*my_mean(y, w)

def my_var(x, w):
    return my_cov(x, x, w)

def my_line(x, m=1, c=0):
    return m*x + c

def y_estrapolato(x, m, c, sigma_m, sigma_c, cov_mc):
    y = m*x + c
    uy = np.sqrt(np.power(x, 2)*np.power(sigma_m, 2) +
                   np.power(sigma_c, 2) + 2*x*cov_mc )
    return y, uy

def lin_fit(x, y, sd_y, xlabel="x [ux]", ylabel="y [uy]", xm=0., xM=1., ym=0., yM=1.,
            verbose=True, plot=False, setrange=False):

    #pesi
    w_y = np.power(sd_y.astype(float), -2)

    #m
    m = my_cov(x, y, w_y) / my_var(x, w_y)
    var_m = 1 / ( my_var(x, w_y) * np.sum(w_y) )

    #c
    c = my_mean(y, w_y) - my_mean(x, w_y) * m
    var_c = my_mean(x*x, w_y)  / ( my_var(x, w_y) * np.sum(w_y) )

    #cov
    cov_mc = - my_mean(x, w_y) / ( my_var(x, w_y) * np.sum(w_y) )

    #rho
    rho_mc = cov_mc / ( sqrt(var_m) * sqrt(var_c) )

    if (verbose):

        print ('m         = ', m.round(4))
        print ('sigma(m)  = ', np.sqrt(var_m).round(4))
        print ('c         = ', c.round(4))
        print ('sigma(c)  = ', np.sqrt(var_c).round(4))
        print ('cov(m, c) = ', cov_mc.round(4))
        print ('rho(m, c) = ', rho_mc.round(4))

    if (plot):

        # rappresento i dati
        plt.errorbar(x, y, yerr=sd_y, xerr=0, ls='', marker='.',
                     color="black",label='dati')

        # costruisco dei punti x su cui valutare la retta del fit
        xmin = float(np.min(x))
        xmax = float(np.max(x))
        xmin_plot = xmin-.2*(xmax-xmin)
        xmax_plot = xmax+.2*(xmax-xmin)
        if (setrange):
            xmin_plot = xm
            xmax_plot = xM
        x1 = np.linspace(xmin_plot, xmax_plot, 100)
        y1 = my_line(x1, m, c)

        # rappresento la retta del fit
        plt.plot(x1, y1, linestyle='-.', color="green", label='fit')
        plt.xlabel(xlabel)
        plt.ylabel(ylabel)
        plt.title('Fit lineare')
        plt.xlim(xmin_plot,xmax_plot)
        if (setrange):
            plt.ylim(ym,yM)

        # rappresento le incertezze sulla retta
        y1_plus_1sigma = y1+y_estrapolato(x1, m, c, np.sqrt(var_m), np.sqrt(var_c), cov_mc)[1]
        y1_minus_1sigma = y1-y_estrapolato(x1, m, c, np.sqrt(var_m), np.sqrt(var_c), cov_mc)[1]
        plt.plot(x1,y1_plus_1sigma, linestyle='-', color="orange", label=r'fit $\pm 1\sigma$')
        plt.plot(x1,y1_minus_1sigma, linestyle='-', color="orange")

        plt.grid()

        plt.legend()

    return m, np.sqrt(var_m), c, np.sqrt(var_c), cov_mc, rho_mc

**Fit Lineare**

# plots will be shown inline
%matplotlib inline
import matplotlib.pyplot as plt

import numpy
from numpy import sqrt,floor

import numpy as np
import scipy as sp
from scipy import stats

import pandas as pd
import random
import math
# libreria locale
#import my_lib_santanastasio as my

##NOTA IMPORTANTE: se cambi qualcosa in my_lib_santanastasio
# devi fare Kernel --> Restart prima di rigirare il codice
# altrimenti i cambiamenti non saranno applicati.

# set global random seed
random_Gseed = 1112
np.random.seed(random_Gseed)

# modello e valori di x
m_true = 3
c_true = 2
sigma_y_true = 2

x = np.arange(1,6,0.5) #default
#x = np.arange(1,6,0.25) #N=20
#x = np.arange(3,4,0.1) #var ridotta
#x = np.arange(-1,9,1) #var aumentata
#x = np.arange(-2.25,2.75,0.5) #baricentro x = 0
#x = np.array([1.,1.1,1.2,1.3,1.4,5.1,5.2,5.3,5.4,5.5]) # var aumentata 1

ux = np.repeat(0.1,len(x))


#---------------------

y_true = m_true*x + c_true
uy_true = np.repeat(sigma_y_true,len(x))
y = np.random.normal(y_true,uy_true)
uy = uy_true*1

print (x)
print (y_true)
print (uy_true)
print (y.round(2))
print (uy)

m0, sm0, c0, sc0, cov0, rho0 = lin_fit(x, y, uy, "x [ux]", "y [uy]", 0, 7, 0, 25, plot=True, setrange=True)

# nuove y
uy_new = np.sqrt(uy**2+(m0*ux)**2)
m, sm, c, sc, cov, rho = lin_fit(x, y, uy_new, "x [ux]", "y [uy]", 0, 7, 0, 25, plot=True, setrange=True)
#m, sm, c, sc, cov, rho = my.lin_fit(x, y, uy_new, "x [ux]", "y [uy]", -2, 9, -10, 40, plot=True, setrange=True)

# Studio dei residui
y_atteso = m*x + c
d = y - y_atteso
d_norm = d / uy_new

plt.errorbar(x,d,uy_new,marker='.',linestyle="")
plt.ylabel("Residui $d=y-y_{atteso}$ [uy]")
plt.xlabel("$x$ [ux]")
plt.grid()

plt.errorbar(x,d_norm,uy_new/uy_new,marker='.',linestyle="")
plt.ylabel("Residui normalizzati $d/\sigma_y=(y-y_{atteso})/\sigma_y$")
plt.xlabel("$x$ [ux]")
plt.grid()

# Incertezze a posteriori
sigmy_post = math.sqrt( np.sum(d**2)/(d.size-2) )
uy_post = np.repeat(sigmy_post,y.size)
print (sigmy_post)

# Nuovo fit con incertezze a posteriori sulle y
m1, sm1, c1, sc1, cov1, rho1 = lin_fit(x, y, uy_post, "x [ux]", "y [uy]", 0, 7, 0, 25, plot=True, setrange=True)

# Parte 1: Studio della riproducibilità delle misure di velocità angolare ω ed accelerazione angolare α


# Importo degli array

#15°
for i in range(1,6):
  #15 gradi
  exec("df15_{0} = pd.read_csv('Dati 15°/Raw Data{0}.csv')".format(i))
  exec("a15_{0} = []".format(i))
  exec("a15_{0}.append(np.asarray(df15_{0}['Time (s)']))".format(i))
  exec("a15_{0}.append(np.asarray(df15_{0}['Gyroscope y (rad/s)']))".format(i))


  df15_p = pd.read_csv('Dati 15°/Plane.csv')
  a15_p = []
  a15_p.append(np.asarray(df15_p['t (s)']))
  a15_p.append(np.asarray(df15_p['Inclination (deg)']))
  #20gradi
  exec("df20_{0} = pd.read_csv('Dati 20°/Raw Data{0}.csv')".format(i))
  exec("a20_{0} = []".format(i))
  exec("a20_{0}.append(np.asarray(df20_{0}['Time (s)']))".format(i))
  exec("a20_{0}.append(np.asarray(df20_{0}['Gyroscope y (rad/s)']))".format(i))
  df20_p = pd.read_csv('Dati 20°/Plane.csv')
  a20_p = []
  a20_p.append(np.asarray(df20_p['t (s)']))
  a20_p.append(np.asarray(df20_p['Inclination (deg)']))
  #30gradi
  exec("df30_{0} = pd.read_csv('Dati 30°/Raw Data{0}.csv')".format(i))
  exec("a30_{0} = []".format(i))
  exec("a30_{0}.append(np.asarray(df30_{0}['Time (s)']))".format(i))
  exec("a30_{0}.append(np.asarray(df30_{0}['Gyroscope y (rad/s)']))".format(i))
  df30_p = pd.read_csv('Dati 30°/Plane.csv')
  a30_p = []
  a30_p.append(np.asarray(df30_p['t (s)']))
  a30_p.append(np.asarray(df30_p['Inclination (deg)']))

for i in range(1,19):
  #15 gradi
  exec("df25_{0} = pd.read_csv('Dati 25°/Raw Data{0}.csv')".format(i))
  exec("a25_{0} = []".format(i))
  exec("a25_{0}.append(np.asarray(df25_{0}['Time (s)']))".format(i))
  exec("a25_{0}.append(np.asarray(df25_{0}['Gyroscope y (rad/s)']))".format(i))

for i in range(0,5):
  exec("df25_p{0} = pd.read_csv('Dati 25°/Plane{0}.csv')".format(i))
  exec("a25_p{0} = []".format(i))
  exec("a25_p{0}.append(np.asarray(df25_p{0}['t (s)']))".format(i))
  exec("a25_p{0}.append(np.asarray(df25_p{0}['Inclination (deg)']))".format(i))
a30_1

# Singolo $\alpha$

Grafico omega vs tempo

plt.plot(a25_1[0],a25_1[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.85
xb = 1.35
plt.vlines(xa,-10,40,colors="orange",linestyle="--")
plt.vlines(xb,-10,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('25deg_wy_vs_t')

xnew = a25_1[0][( a25_1[0]>xa)&( a25_1[0]<xb)]
ynew = a25_1[1][(a25_1[0]>xa)&( a25_1[0]<xb)]
#questo sigma serve solo per essere inserito nel codice, è stato scientificamente dimostrato che non ha impatto
#dato che i valori di m e c non cambiano nonstante i sigma messi
sa25_1 = np.repeat(np.std(ynew),len(xnew))

$\sigma[y]$ è ignota e supponiamola costante per tutte le y. Lo calcoleremo a posteriori con la formula 14.22 (p.305) delle dispense ricavandoci m e c con verbose=True (questa opzione calcola i coefficienti del fit senza creare il grafico), e inserendo i sigma trovati nel fit

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa25_1, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)

# Zoommiamo perché non si vedono le barre d'errore
axins = zoomed_inset_axes(ax, 5, loc=4) # zoom-factor: 2.5, location: upper-left
axins.set_xlim(0.79, 0.81) # apply the x-limits
axins.set_ylim(0.65, 0.69) # apply the y-limits
axins.errorbar(m_x, T2_y, yerr=T2_errorbar, fmt='o', ecolor='red')
plt.plot(xp, (alpha1[0]*xp+c1[0]), label="fit lineare", color="green")
mark_inset(ax, axins, loc1=2, loc2=4, fc="none", ec="0.5")
plt.xticks(visible=False)
plt.grid()


plt.savefig('Fit lin 25deg_1')

#this part creates a list which will be used to make a table on latex
alphai = [[],[], []]
alphai[0].append(m)
alphai[1].append(sm)
#rho is also included because, quoting the manual:'Si noti inoltre come il coefficiente di correlazione non sia un optional, ma faccia parte integrante del
#risultato (vedi anche 14.4 e discussione nel paragrafo 14.3, in particolare punti 8 e 10)'.
alphai[2].append(rho)
print('alpha=',m,'+/-',sm)

For a single data set, $\alpha$ is the coefficient of the line, with its associated standard deviation. In other words, $\alpha= m \pm \sigma(m)$.
The expected value of m is: $\frac{MgR}{I_{CM}+MR^2}\cdot sin(\theta)$ where M is the total mass of the system, R is the radius of the cylinder, $I_{CM}$ is the moment of inertia of the cylinder, and is thus of the form: $I_{CM}= m_{cylinder}R^2$.

# $\alpha$ per N set di dati

Le celle che seguono creano, a partire dei dati delle N misure, una tabella con gli $\alpha$ e una stima di un 'best' alpha a partire dai precedenti.

N=2

plt.plot(a25_2[0],a25_2[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.5
xb = 1
plt.vlines(xa,-1,40,colors="orange",linestyle="--")
plt.vlines(xb,-1,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('25deg_wy_vs_t-2')


xnew = a25_2[0][( a25_2[0]>xa)&( a25_2[0]<xb)]
ynew = a25_2[1][( a25_2[0]>xa)&( a25_2[0]<xb)]

sa25_2 = np.repeat(a25_2[0][1]/(12**0.5),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa25_2, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 25deg_2')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alphai[0].append(m)
alphai[1].append(sm)
alphai[2].append(rho)
print('alpha=',m,'+/-',sm)

N=3

plt.plot(a25_3[0],a25_3[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.35
xb = 0.82
plt.vlines(xa,-1,40,colors="orange",linestyle="--")
plt.vlines(xb,-1,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('25deg_wy_vs_t-3')


xnew = a25_3[0][( a25_3[0]>xa)&( a25_3[0]<xb)]
ynew = a25_3[1][( a25_3[0]>xa)&( a25_3[0]<xb)]

sa25_3 = np.repeat(a25_3[0][1]/(12**0.5),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa25_3, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 25deg_3')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alphai[0].append(m)
alphai[1].append(sm)
alphai[2].append(rho)
print('alpha=',m,'+/-',sm)

N=4

plt.plot(a25_4[0],a25_4[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.75
xb = 1.25
plt.vlines(xa,-1,40,colors="orange",linestyle="--")
plt.vlines(xb,-1,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('25deg_wy_vs_t-4')


xnew = a25_4[0][( a25_4[0]>xa)&( a25_4[0]<xb)]
ynew = a25_4[1][( a25_4[0]>xa)&( a25_4[0]<xb)]

sa25_4 = np.repeat(a25_4[0][1]/(12**0.5),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa25_4, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 25deg_4')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alphai[0].append(m)
alphai[1].append(sm)
alphai[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 5

plt.plot(a25_5[0],a25_5[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.35
xb = 0.85
plt.vlines(xa,-1,40,colors="orange",linestyle="--")
plt.vlines(xb,-1,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('25deg_wy_vs_t-5')


xnew = a25_5[0][( a25_5[0]>xa)&( a25_5[0]<xb)]
ynew = a25_5[1][( a25_5[0]>xa)&( a25_5[0]<xb)]

sa25_5 = np.repeat(a25_5[0][1]/(12**0.5),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa25_5, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 25deg_5')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alphai[0].append(m)
alphai[1].append(sm)
alphai[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 6

plt.plot(a25_6[0],a25_6[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.65
xb = 1.12
plt.vlines(xa,-1,40,colors="orange",linestyle="--")
plt.vlines(xb,-1,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('25deg_wy_vs_t-6')


xnew = a25_6[0][( a25_6[0]>xa)&( a25_6[0]<xb)]
ynew = a25_6[1][( a25_6[0]>xa)&( a25_6[0]<xb)]

sa25_6 = np.repeat(a25_6[0][1]/(12**0.5),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa25_6, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 25deg_6')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alphai[0].append(m)
alphai[1].append(sm)
alphai[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 7

plt.plot(a25_7[0],a25_7[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.32
xb = 0.82
plt.vlines(xa,-1,40,colors="orange",linestyle="--")
plt.vlines(xb,-1,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('25deg_wy_vs_t-7')


xnew = a25_7[0][( a25_7[0]>xa)&( a25_7[0]<xb)]
ynew = a25_7[1][( a25_7[0]>xa)&( a25_7[0]<xb)]

sa25_7 = np.repeat(a25_7[0][1]/(12**0.5),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa25_7, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 25deg_7')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alphai[0].append(m)
alphai[1].append(sm)
alphai[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 8

plt.plot(a25_8[0],a25_8[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.32
xb = 0.82
plt.vlines(xa,-1,40,colors="orange",linestyle="--")
plt.vlines(xb,-1,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('25deg_wy_vs_t-8')


xnew = a25_8[0][( a25_8[0]>xa)&( a25_8[0]<xb)]
ynew = a25_8[1][( a25_8[0]>xa)&( a25_8[0]<xb)]

sa25_8 = np.repeat(a25_8[0][1]/(12**0.5),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa25_8, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 25deg_8')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alphai[0].append(m)
alphai[1].append(sm)
alphai[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 9

plt.plot(a25_9[0],a25_9[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.9
xb = 1.4
plt.vlines(xa,-1,40,colors="orange",linestyle="--")
plt.vlines(xb,-1,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('25deg_wy_vs_t-9')


xnew = a25_9[0][( a25_9[0]>xa)&( a25_9[0]<xb)]
ynew = a25_9[1][( a25_9[0]>xa)&( a25_9[0]<xb)]
sa25_9 = np.repeat(a25_9[0][1]/(12**0.5),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa25_9, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 25deg_9')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alphai[0].append(m)
alphai[1].append(sm)
alphai[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 10

plt.plot(a25_10[0],a25_10[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 1.55
xb = 2
plt.vlines(xa,-1,40,colors="orange",linestyle="--")
plt.vlines(xb,-1,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('25deg_wy_vs_t-10')


xnew = a25_10[0][( a25_10[0]>xa)&( a25_10[0]<xb)]
ynew = a25_10[1][( a25_10[0]>xa)&( a25_10[0]<xb)]
sa25_10 = np.repeat(a25_10[0][1]/(12**0.5),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa25_10, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 25deg_10')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alphai[0].append(m)
alphai[1].append(sm)
alphai[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 11

plt.plot(a25_11[0],a25_11[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.62
xb = 1.12
plt.vlines(xa,-1,40,colors="orange",linestyle="--")
plt.vlines(xb,-1,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('25deg_wy_vs_t-11')


xnew = a25_11[0][( a25_11[0]>xa)&( a25_11[0]<xb)]
ynew = a25_11[1][( a25_11[0]>xa)&( a25_11[0]<xb)]
sa25_11 = np.repeat(a25_11[0][1]/(12**0.5),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa25_11, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 25deg_11')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alphai[0].append(m)
alphai[1].append(sm)
alphai[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 12

plt.plot(a25_12[0],a25_12[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 1.15
xb = 1.67
plt.vlines(xa,-1,40,colors="orange",linestyle="--")
plt.vlines(xb,-1,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('25deg_wy_vs_t-12')


xnew = a25_12[0][( a25_12[0]>xa)&( a25_12[0]<xb)]
ynew = a25_12[1][( a25_12[0]>xa)&( a25_12[0]<xb)]
sa25_12 = np.repeat(a25_12[0][1]/(12**0.5),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa25_12, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 25deg_12')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alphai[0].append(m)
alphai[1].append(sm)
alphai[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 13

plt.plot(a25_13[0],a25_13[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.78
xb = 1.25
plt.vlines(xa,-1,40,colors="orange",linestyle="--")
plt.vlines(xb,-1,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('25deg_wy_vs_t-13')


xnew = a25_13[0][( a25_13[0]>xa)&( a25_13[0]<xb)]
ynew = a25_13[1][( a25_13[0]>xa)&( a25_13[0]<xb)]
sa25_13 = np.repeat(a25_13[0][1]/(12**0.5),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa25_13, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 25deg_13')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alphai[0].append(m)
alphai[1].append(sm)
alphai[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 14

plt.plot(a25_14[0],a25_14[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.35
xb = 0.82
plt.vlines(xa,-1,40,colors="orange",linestyle="--")
plt.vlines(xb,-1,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('25deg_wy_vs_t-14')


xnew = a25_14[0][( a25_14[0]>xa)&( a25_14[0]<xb)]
ynew = a25_14[1][( a25_14[0]>xa)&( a25_14[0]<xb)]
sa25_14 = np.repeat(a25_14[0][1]/(12**0.5),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa25_14, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 25deg_14')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alphai[0].append(m)
alphai[1].append(sm)
alphai[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 15

plt.plot(a25_15[0],a25_15[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.425
xb = 0.925
plt.vlines(xa,-1,40,colors="orange",linestyle="--")
plt.vlines(xb,-1,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('25deg_wy_vs_t-15')


xnew = a25_15[0][( a25_15[0]>xa)&( a25_15[0]<xb)]
ynew = a25_15[1][( a25_15[0]>xa)&( a25_15[0]<xb)]
sa25_15 = np.repeat(a25_15[0][1]/(12**0.5),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa25_15, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 25deg_15')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alphai[0].append(m)
alphai[1].append(sm)
alphai[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 16

plt.plot(a25_16[0],a25_16[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.25
xb = 0.75
plt.vlines(xa,-1,40,colors="orange",linestyle="--")
plt.vlines(xb,-1,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('25deg_wy_vs_t-16')


xnew = a25_16[0][( a25_16[0]>xa)&( a25_16[0]<xb)]
ynew = a25_16[1][( a25_16[0]>xa)&( a25_16[0]<xb)]
sa25_16 = np.repeat(a25_16[0][1]/(12**0.5),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa25_16, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 25deg_16')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alphai[0].append(m)
alphai[1].append(sm)
alphai[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 17

plt.plot(a25_17[0],a25_17[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 1
xb = 1.5
plt.vlines(xa,-1,40,colors="orange",linestyle="--")
plt.vlines(xb,-1,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('25deg_wy_vs_t-17')


xnew = a25_17[0][( a25_17[0]>xa)&( a25_17[0]<xb)]
ynew = a25_17[1][( a25_17[0]>xa)&( a25_17[0]<xb)]
sa25_17 = np.repeat(a25_17[0][1]/(12**0.5),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa25_17, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 25deg_17')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alphai[0].append(m)
alphai[1].append(sm)
alphai[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 18

plt.plot(a25_18[0],a25_18[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.36
xb = 0.86
plt.vlines(xa,-1,40,colors="orange",linestyle="--")
plt.vlines(xb,-1,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('25deg_wy_vs_t-18')


xnew = a25_18[0][( a25_18[0]>xa)&( a25_18[0]<xb)]
ynew = a25_18[1][( a25_18[0]>xa)&( a25_18[0]<xb)]
sa25_18 = np.repeat(a25_18[0][1]/(12**0.5),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa25_18, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 25deg_18')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alphai[0].append(m)
alphai[1].append(sm)
alphai[2].append(rho)
print('alpha=',m,'+/-',sm)

# Tabella degli $\alpha$

#Prints the data in a latex-format tabular to copy and paste
df = pd.DataFrame( {'alpha': alphai[0]
                   ,'sigma(alpha)': alphai[1]
                   ,'rho': alphai[2]}
                  )

print(df.to_latex(index=False,formatters={'alphai[0]': '{:.2f}'.format,
                                           'alphai[1]': '{:.2f}'.format,
                                              'alphai[2]': '{:.2f}'.format}
                  ))

# $\alpha_{madre di tutti i fottuti alfa}$

How the fuck do I calculate this. There you go: formule 13.18 e 13.19 (capitolo 13, p.281):
$\alpha_{mother} \equiv E[\alpha] = \frac{\sum_i \alpha_i/\sigma_i^2}{\sum_i \sigma_i^2}$
e
$1/\sigma_{mother}^2 \equiv 1/\sigma_{\alpha}^2 = \frac{1}{\sum_i \sigma_i^2}$

sum_a = 0
sum_s = 0
for i in range(0,19):
    sum_a += alphai[0][i]*(alphai[1][i]**2)
    sum_s += alphai[1][i]**2
alpha_mother = sum_a/sum_s
sigma_mother = np.sqrt(sum_s)
print('alpha_m=', alpha_mother,'+/-', sigma_mother)

# Parte 2: Misura dell’accelerazione di gravità


# Import of data 2: Electric boogaloo

# Misura di $\alpha$ per fit lineare

**Important to note that while it absolutely isn't required, I did the graphs. In fact it would be better if I only showed the results**

## 15°

N = 1

plt.plot(a15_1[0],a15_1[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.28
xb = 1.05
plt.vlines(xa,-10,40,colors="orange",linestyle="--")
plt.vlines(xb,-10,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('15deg_wy_vs_t')

xnew = a15_1[0][( a15_1[0]>xa)&( a15_1[0]<xb)]
ynew = a15_1[1][(a15_1[0]>xa)&( a15_1[0]<xb)]
#questo sigma serve solo per essere inserito nel codice, è stato scientificamente dimostrato che non ha impatto
#dato che i valori di m e c non cambiano nonstante i sigma messi
sa15_1 = np.repeat(np.std(ynew),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa15_1, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 15deg_1')

alpha15 = [[],[],[]]
#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alpha15[0].append(m)
alpha15[1].append(sm)
alpha15[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 2

plt.plot(a15_2[0],a15_2[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.28
xb = 1.05
plt.vlines(xa,-10,40,colors="orange",linestyle="--")
plt.vlines(xb,-10,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('15deg_wy_vs_t2')

xnew = a15_2[0][( a15_2[0]>xa)&( a15_2[0]<xb)]
ynew = a15_2[1][(a15_2[0]>xa)&( a15_2[0]<xb)]
#questo sigma serve solo per essere inserito nel codice, è stato scientificamente dimostrato che non ha impatto
#dato che i valori di m e c non cambiano nonstante i sigma messi
sa15_2 = np.repeat(np.std(ynew),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa15_2, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 15deg_2')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alpha15[0].append(m)
alpha15[1].append(sm)
alpha15[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 3

plt.plot(a15_3[0],a15_3[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.34
xb = 1.05
plt.vlines(xa,-10,40,colors="orange",linestyle="--")
plt.vlines(xb,-10,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('15deg_wy_vs_t3')

xnew = a15_3[0][( a15_3[0]>xa)&( a15_3[0]<xb)]
ynew = a15_3[1][(a15_3[0]>xa)&( a15_3[0]<xb)]
#questo sigma serve solo per essere inserito nel codice, è stato scientificamente dimostrato che non ha impatto
#dato che i valori di m e c non cambiano nonstante i sigma messi
sa15_3 = np.repeat(np.std(ynew),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa15_3, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 15deg_3')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alpha15[0].append(m)
alpha15[1].append(sm)
alpha15[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 4

plt.plot(a15_4[0],a15_4[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 1.97
xb = 2.7
plt.vlines(xa,-10,40,colors="orange",linestyle="--")
plt.vlines(xb,-10,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('15deg_wy_vs_t4')

xnew = a15_4[0][( a15_4[0]>xa)&( a15_4[0]<xb)]
ynew = a15_4[1][(a15_4[0]>xa)&( a15_4[0]<xb)]
#questo sigma serve solo per essere inserito nel codice, è stato scientificamente dimostrato che non ha impatto
#dato che i valori di m e c non cambiano nonstante i sigma messi
sa15_3 = np.repeat(np.std(ynew),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa15_3, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 15deg_4')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alpha15[0].append(m)
alpha15[1].append(sm)
alpha15[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 5

plt.plot(a15_5[0],a15_5[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.55
xb = 1.25
plt.vlines(xa,-10,40,colors="orange",linestyle="--")
plt.vlines(xb,-10,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('15deg_wy_vs_t5')

xnew = a15_5[0][( a15_5[0]>xa)&( a15_5[0]<xb)]
ynew = a15_5[1][(a15_5[0]>xa)&( a15_5[0]<xb)]
#questo sigma serve solo per essere inserito nel codice, è stato scientificamente dimostrato che non ha impatto
#dato che i valori di m e c non cambiano nonstante i sigma messi
sa15_3 = np.repeat(np.std(ynew),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa15_3, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 15deg_5')

#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alpha15[0].append(m)
alpha15[1].append(sm)
alpha15[2].append(rho)
print('alpha=',m,'+/-',sm)

#Prints the data in a latex-format tabular to copy and paste
df = pd.DataFrame( {'alpha': alpha15[0]
                   ,'sigma(alpha)': alpha15[1]
                   ,'rho': alpha15[2]}
                  )

print(df.to_latex(index=False,formatters={'alpha15[0]': '{:.2f}'.format,
                                           'alpha15[1]': '{:.2f}'.format,
                                              'alpha15[2]': '{:.2f}'.format}
                  ))

sum_a = 0
sum_s = 0
for i in range(0,6):
    sum_a += alphai[0][i]*(alphai[1][i]**2)
    sum_s += alphai[1][i]**2
alpha_mother = sum_a/sum_s
sigma_mother = np.sqrt(sum_s)
print('alpha15_m=', alpha_mother,'+/-', sigma_mother)

## 20°

plt.plot(a20_1[0],a20_1[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.38
xb = 1.02
plt.vlines(xa,-10,40,colors="orange",linestyle="--")
plt.vlines(xb,-10,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('20deg_wy_vs_t')

xnew = a20_1[0][( a20_1[0]>xa)&( a20_1[0]<xb)]
ynew = a20_1[1][(a20_1[0]>xa)&( a20_1[0]<xb)]
#questo sigma serve solo per essere inserito nel codice, è stato scientificamente dimostrato che non ha impatto
#dato che i valori di m e c non cambiano nonstante i sigma messi
sa15_1 = np.repeat(np.std(ynew),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa15_1, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 20deg_1')

alpha20 = [[],[],[]]
#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alpha20[0].append(m)
alpha20[1].append(sm)
alpha20[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 2

plt.plot(a20_2[0],a20_2[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.5
xb = 1.15
plt.vlines(xa,-10,40,colors="orange",linestyle="--")
plt.vlines(xb,-10,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('20deg_wy_vs_t2')

xnew = a20_2[0][( a20_2[0]>xa)&( a20_2[0]<xb)]
ynew = a20_2[1][(a20_2[0]>xa)&( a20_2[0]<xb)]
#questo sigma serve solo per essere inserito nel codice, è stato scientificamente dimostrato che non ha impatto
#dato che i valori di m e c non cambiano nonstante i sigma messi
sa15_1 = np.repeat(np.std(ynew),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa15_1, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 20deg_2')


#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alpha20[0].append(m)
alpha20[1].append(sm)
alpha20[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 3

plt.plot(a20_3[0],a20_3[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.55
xb = 1.18
plt.vlines(xa,-10,40,colors="orange",linestyle="--")
plt.vlines(xb,-10,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('20deg_wy_vs_t3')

xnew = a20_3[0][( a20_3[0]>xa)&( a20_3[0]<xb)]
ynew = a20_3[1][(a20_3[0]>xa)&( a20_3[0]<xb)]
#questo sigma serve solo per essere inserito nel codice, è stato scientificamente dimostrato che non ha impatto
#dato che i valori di m e c non cambiano nonstante i sigma messi
sa15_1 = np.repeat(np.std(ynew),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa15_1, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 20deg_3')


#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alpha20[0].append(m)
alpha20[1].append(sm)
alpha20[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 4

plt.plot(a20_4[0],a20_4[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.75
xb = 1.35
plt.vlines(xa,-10,40,colors="orange",linestyle="--")
plt.vlines(xb,-10,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('20deg_wy_vs_t4')

xnew = a20_4[0][( a20_4[0]>xa)&( a20_4[0]<xb)]
ynew = a20_4[1][(a20_4[0]>xa)&( a20_4[0]<xb)]
#questo sigma serve solo per essere inserito nel codice, è stato scientificamente dimostrato che non ha impatto
#dato che i valori di m e c non cambiano nonstante i sigma messi
sa15_1 = np.repeat(np.std(ynew),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa15_1, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 20deg_4')


#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alpha20[0].append(m)
alpha20[1].append(sm)
alpha20[2].append(rho)
print('alpha=',m,'+/-',sm)

N = 5

plt.plot(a20_5[0],a20_5[1],label="wy", marker='.', color='green')
plt.xlabel('Tempo[s]')
plt.ylabel('Velocità angolare in y[rad/s]')
plt.grid()
xa = 0.52
xb = 1.15
plt.vlines(xa,-10,40,colors="orange",linestyle="--")
plt.vlines(xb,-10,40,colors="orange",linestyle="--")
plt.legend()
plt.savefig('20deg_wy_vs_t5')

xnew = a20_5[0][( a20_5[0]>xa)&( a20_5[0]<xb)]
ynew = a20_5[1][(a20_5[0]>xa)&( a20_5[0]<xb)]
#questo sigma serve solo per essere inserito nel codice, è stato scientificamente dimostrato che non ha impatto
#dato che i valori di m e c non cambiano nonstante i sigma messi
sa15_1 = np.repeat(np.std(ynew),len(xnew))

#Qui calcoleremo m e c
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sa15_1, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose = True, plot=False, setrange=False)

#Questa parte del codice applica la formula (applicata solo nel intervallo dove c'è moto)
#Il numeratore:
num = 0
for i in range(0, len(ynew)):
    num += (ynew[i] - m*xnew[i] - c)**2
#E il valore cercato
sigma1 = np.sqrt(num/(len(ynew)-2))
sigma = []
sigma = np.repeat(sigma1, len(ynew))

#Qui graficheremo finalmente il fit lineare
m, sm, c, sc, cov, rho = lin_fit(xnew, ynew, sigma, "tempo [s]", "Vel. angolare $\omega$ [rad/s]", verbose=False, plot=True, setrange=False)
plt.savefig('Fit lin 20deg_5')


#E qui inseriamo i nuovi valori di m e c in alphai per la tabella
alpha20[0].append(m)
alpha20[1].append(sm)
alpha20[2].append(rho)
print('alpha=',m,'+/-',sm)

#Prints the data in a latex-format tabular to copy and paste
df = pd.DataFrame( {'alpha': alpha20[0]
                   ,'sigma(alpha)': alpha20[1]
                   ,'rho': alpha20[2]}
                  )

print(df.to_latex(index=False,formatters={'alpha20[0]': '{:.2f}'.format,
                                           'alpha20[1]': '{:.2f}'.format,
                                              'alpha20[2]': '{:.2f}'.format}
                  ))

sum_a = 0
sum_s = 0
for i in range(0,6):
    sum_a += alphai[0][i]*(alphai[1][i]**2)
    sum_s += alphai[1][i]**2
alpha_mother = sum_a/sum_s
sigma_mother = np.sqrt(sum_s)
print('alpha20_m=', alpha_mother,'+/-', sigma_mother)



# $\alpha$ in funzione di $sin(\theta)$

Per questo quesito, $\alpha$ verrà espressa (con la sua incertezza, quindi in plt.errorbar) in funzione dei diversi angoli $\theta$. (Perciò è importante avere un numero relativamente grande di angoli $\theta$). Da questo si potrà estrarre un coefficiente angolare $m$.

#al = np.asarray(alphai[0])
#sal = np.asarray(alphai[1])
al = [0.3, 0.45, 0.66, 0.73, 0.99, 1.3]
al = np.asarray(al)
sal = [0.03, 0.045, 0.066, 0.073, 0.099, 0.13]
sal = np.asarray(sal)
#The angles for theta are also obviously placeholders, whatsmore it's also not unlikely it will be in degrees and not rads
theta = [0, np.pi/6, np.pi/4, np.pi/3, np.pi/2, np.pi]
theta = np.asarray(theta)
m, sm, c, sc, cov, rho = lin_fit(theta, al, sal, "sin(theta)", "alpha", plot=True, setrange=False)


**Importante: non ho considerato la parte in italico nei quesiti, perciò va ancora fatta quella parte**

# Misura di $I_{CM}$ e di $g$

Questa parte presuppone che le misure dei vari corpi del sistema sono già state fatte. Useremo la formula per il momento d'inerzia del cilindro.


La formula usata per calcolare il momento d'inerzia può variare in funzione di come consideriamo il cilindro: se è pieno, allora $I_{CM}=1/2*MR^2$, mentre se lo consideriamo vuoto $I_{CM}=MR^2$.

#Assumendo che essa sarà considerata piena e omogenea:
M = #kg
sM = #deviazione standard di M
R = #m
sR = #deviazione standard di R
I = M*R^2
#La deviazione standard su I sarà calcolata con la formula di propagazione delle incertezze per i monomi:
sI = I*((sM/M)**2+4*(sR/R)**2)**0.5
print('I_cm=', I,'+/-',sI)

Il calcolo di g, leggermente più complicato di quello di I, si basa sulla formula (5) del PP4: $m = \frac{MgR}{I_{CM}+M\cdot R^2}$ ossia $g = m\frac{I_{CM}+M\cdot R^2}{MR}$ dove $m$ è il coefficiente angolare del grafico di $\alpha$ vs $sin(\theta)$.


Il calcolo di g verrà fatto  rimpiazzando $I_{CM}$ con la sua formula per cilindri pieni.

#Rimpiazzando con la formula. Damn coffee is effective. Staying awake is child's play
g1 = m*3/2*R
#La deviazione standard è calcolata con la formula di propagazione delle incertezze per monomi
sg1 = g1*((sm/m)**2+(sR/R)**2)**0.5
print('g1=',g1,'+/-',sg1)