# coding: utf-8

# Rutimas básicas para trabajar en el laboratorio. 
# Carga de datos y Ajuste lineal

# In[3]:


# Paso 1. Cargo las librerias que usaré.

import numpy as np
from matplotlib import pyplot as plt
from scipy.optimize import curve_fit
#%matplotlib inline   
# si agrego %matplotlib inline muestra los gráficos en la misma ventana sino los abre en otra aparte
get_ipython().run_line_magic('matplotlib', 'inline')


# In[4]:


# CASO 1: MIS DATOS ESTÀN APARTE EN UN TXT
# Experimento de laboratorio. Mido la variable x, la variable y, estimo el error dy y armo una planilla con esas 3 columnas
# Guardo esa planilla en un archivo txt 

# Leo los datos 
Misdatos = np.loadtxt('Misdatos.txt')
np.shape(Misdatos)   #con esto veo las dimensiones y verifico que se cargaron todos los datos


# In[5]:


# Para mayor comodidad asigno Misdatos a las variables x, y, dy
x = Misdatos[:,0]
y = Misdatos[:,1]
dy = Misdatos[:,2]

# significado de los ':'.  Misdatos es una matriz de 35 x 3. Si le pido [:,0] significa que quiero todas las filas 
# de la primer col, si le pido [0,:] significa que quiero de la primer fila todas las columnas.


# In[9]:


# ploteo los datos para chequearlos
plt.plot(x,y,'bo', label = 'datos')
plt.errorbar(x,y,dy, linestyle = 'None')
plt.title('X vs Y')
plt.xlabel('Mediciones X')
plt.ylabel('Valores y')
plt.legend(loc = 'best')  #pùede ser upper center, upper right, upper left, (idem pero con:) lower 
plt.show()


# In[6]:


# si quiero reescalar el plot agrego los comandos xlim o ylim
plt.plot(x,y,'bo', label = 'datos')
plt.errorbar(x,y,dy, linestyle = 'None')
plt.title('X vs Y')
plt.xlabel('Mediciones X')
plt.ylabel('Valores y')
plt.legend(loc = 'best')

plt.xlim(600,1000)
plt.ylim(2,4)

plt.show()


# In[7]:


# Hago un ajuste lineal, para ello defino la función lineal que usará mi variable x
# los parámetros a y b se inicializan solos en 1, sino se puede cambiar
def func(x, a, b):
        return a*x + b

param, cov = curve_fit(func, x, y, sigma = dy)

# param y cov son las salidas del ajuste
# en par están los coeficientes de ajuste, en este caso a y b
# en cov está la matriz de covarianza, donde los elementos de la diagonal me permiten caluclar el error de a y b como: 
da = np.sqrt(cov[0,0])
db = np.sqrt(cov[1,1])


# In[8]:


plt.plot(x,y,'bo', label = 'datos')
plt.plot(x,func(x, param[0], param[1]), 'r-', label = 'Ajuste') 
plt.errorbar(x, y, dy, linestyle = 'None')
plt.title('X vs Y')
plt.xlabel('Mediciones X')
plt.ylabel('Valores y')
plt.legend()
plt.show()


# In[10]:


# le pido los parámetros de ajuste cofs y el valor estadístico p
print('a=', param[0], 'b=', param[1])
print('da=', da, 'db=', db)

# Función para calcular el valor de R2. Estas lineas no se cambian
def rsquared(popt1,popt2, pcov, f, xdata, ydata):
  residuals = ydata- f(xdata, popt1,popt2)
  ss_res = np.sum(residuals**2)
  ss_tot = np.sum((ydata-np.mean(ydata))**2)
  r_squared = 1 - (ss_res / ss_tot)
  return r_squared
# fin

# Acá le paso mis parámetros a la función rsqared
rsquared(param[0],param[1],cov,func,x,y)


# In[20]:


# EJEMPLO: veo en el plot que hay datos que quiero sacar del ajuste. 
# por ejemplo quiero quedarme solo con los x desde 600 

x_min = 600

# uso este comando de lógica para redefinir un nuevo vectos x_sub que contenga los valores de > 600
# luego busco los y que se corresponden con esos mismos x
x_sub = x[x > x_min]
y_sub = y[x > x_min]
dy_sub =dy[x > x_min]


# In[19]:


def func_sub(x_sub, a, b):
        return a*x_sub + b

param_sub, cov_sub = curve_fit(func_sub, x_sub, y_sub, sigma = dy_sub)

plt.plot(x_sub,y_sub,'bo', label = 'datos')
plt.plot(x_sub,func(x_sub, param_sub[0], param_sub[1]), 'r-', label = 'Ajuste') 
plt.errorbar(x_sub, y_sub, dy_sub, linestyle = 'None')
plt.title('X vs Y subplot')
plt.xlabel('Mediciones X')
plt.ylabel('Valores y')
plt.legend()
plt.show()

print('a=', param_sub[0], 'b=', param_sub[1])
print('pa=', np.sqrt(cov_sub[0,0]), 'pb=', np.sqrt(cov_sub[1,1]))


# In[4]:


# CASO 2: voy cargando los datos acá mismo

# los voy cargando en estos vectores a medida que los mido
x = np.array([1, 2, 3])   
y = np.array([10, 15, 17])
dy = np.array([0.1, 0.2, 0.4])


# In[5]:


plt.plot(x,y, 'bo', label = 'datos')
plt.errorbar(x,y,dy, linestyle = 'None')
plt.show()