lolo2

analisis/plots/20231123_CPTconmicromocion3/CPT_plotter_20231123.py

@@ -1970,10 +1970,3 @@ for selectedcurve in SelectedCurveVec:
 
     print(f'listo med {selectedcurve}')
     print(popt_3_SA)
-
-
-
-
-
-
-

analisis/plots/20231123_CPTconmicromocion3/analisis_superajuste.py

+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+"""
+Ploteo de datos y ajustes
+@author: lolo
+"""
+
+
+import h5py
+import matplotlib.pyplot as plt
+import numpy as np
+# import sys
+# import re
+# import ast
+from scipy.optimize import curve_fit
+# import os
+# from scipy import interpolate
+
+
+
+
+#%% Importaciones extra
+
+# /home/lolo/Dropbox/marce/LIAF/Trampa_anular/artiq_experiments/analisis/plots/20231123_CPTconmicromocion3/Data/EITfit/MM_eightLevel_2repumps_AnalysisFunctions.py
+
+from Data.EITfit.lolo_modelo_full_8niveles import PerformExperiment_8levels_MM
+
+
+PARAMETROS = np.load('PARAMETROS.npz',allow_pickle=True)
+for var_name in PARAMETROS.keys():
+    globals()[var_name] = PARAMETROS[var_name]
+    print(f'loaded: {var_name}')
+
+
+#%%
+
+"""
+Primero tengo mediciones de espectros cpt de un ion variando la tension dc_A
+"""
+
+#C:\Users\Usuario\Documents\artiq\artiq_experiments\analisis\plots\20220106_CPT_DosLaseres_v08_TISA_DR\Data
+
+folder = '../20231123_CPTconmicromocion3/Data/'
+CPT_FILES = f"""
+{folder}/000016262-IR_Scan_withcal_optimized
+{folder}/000016239-IR_Scan_withcal_optimized
+{folder}/000016240-IR_Scan_withcal_optimized
+{folder}/000016241-IR_Scan_withcal_optimized
+{folder}/000016244-IR_Scan_withcal_optimized
+{folder}/000016255-IR_Scan_withcal_optimized
+{folder}/000016256-IR_Scan_withcal_optimized
+{folder}/000016257-IR_Scan_withcal_optimized
+"""
+
+
+def SeeKeys(files):
+    for i, fname in enumerate(files.split()):
+        data = h5py.File(fname+'.h5', 'r') # Leo el h5: Recordar que nuestros datos estan en 'datasets'
+        print(fname)
+        print(list(data['datasets'].keys()))
+
+print(SeeKeys(CPT_FILES))
+
+#carpeta pc nico labo escritorio:
+#C:\Users\Usuario\Documents\artiq\artiq_experiments\analisis\plots\20211101_CPT_DosLaseres_v03\Data
+
+Counts = []
+Freqs = []
+
+AmpTisa = []
+UVCPTAmp = []
+No_measures = []
+Voltages = []
+
+for i, fname in enumerate(CPT_FILES.split()):
+    print(str(i) + ' - ' + fname)
+    #print(fname)
+    data = h5py.File(fname+'.h5', 'r') # Leo el h5: Recordar que nuestros datos estan en 'datasets'
+
+    # Aca hago algo repugnante para poder levantar los strings que dejamos
+    # que además tenian un error de tipeo al final. Esto no deberá ser necesario
+    # cuando se solucione el error este del guardado.
+    Freqs.append(np.array(data['datasets']['IR1_Frequencies']))
+    Counts.append(np.array(data['datasets']['data_array']))
+    #AmpTisa.append(np.array(data['datasets']['TISA_CPT_amp']))
+    UVCPTAmp.append(np.array(data['datasets']['UV_CPT_amp']))
+    No_measures.append(np.array(data['datasets']['no_measures']))
+    Voltages.append(np.array(data['datasets']['scanning_voltages']))
+
+def Split(array,n):
+    length=len(array)/n
+    splitlist = []
+    jj = 0
+    while jj<length:
+        partial = []
+        ii = 0
+        while ii < n:
+            partial.append(array[jj*n+ii])
+            ii = ii + 1
+        splitlist.append(partial)
+        jj = jj + 1
+    return splitlist
+
+
+CountsSplit = []
+CountsSplit.append(Split(Counts[0],len(Freqs[0])))
+
+
+CountsSplit_2ions = []
+CountsSplit_2ions.append(Split(Counts[4],len(Freqs[4])))
+
+
+
+#%%
+"""
+AHORA INTENTO SUPER AJUSTES O SEA CON OFFSETXPI Y DETDOPPLER INCLUIDOS
+
+La 0 no ajusta bien incluso con todos los parametros libres
+De la 1 a la 11 ajustan bien
+"""
+
+
+#from EITfit.MM_eightLevel_2repumps_AnalysisFunctions import PerformExperiment_8levels
+from scipy.optimize import curve_fit
+import time
+
+"""
+SUPER AJUSTE (SA)
+"""
+
+phidoppler, titadoppler = 0, 90
+phirepump, titarepump = 0,  0
+phiprobe = 0
+titaprobe = 90
+
+Temp = 0.5e-3
+
+sg = 0.544
+sp = 4.5
+sr = 0
+DetRepump = 0
+
+
+lw = 0.1
+DopplerLaserLinewidth, RepumpLaserLinewidth, ProbeLaserLinewidth = lw, lw, lw #ancho de linea de los laseres
+
+
+u = 32.5e6
+
+#B = (u/(2*np.pi))/c
+
+correccion = 13
+
+#DetDoppler = -11.5-correccion
+
+gPS, gPD, = 2*np.pi*21.58e6, 2*np.pi*1.35e6
+alpha = 0
+
+
+drivefreq = 2*np.pi*22.135*1e6
+
+
+SelectedCurveVec = [1,2,3,4,5,6,7,8,9,10,11]
+#SelectedCurveVec = [10]
+
+# if not 'popt_SA_vec' in globals().keys() or len(popt_SA_vec)==0:
+
+popt_SA_vec = []
+pcov_SA_vec = []
+Detuningsshort_vec = []
+Counts_vec = []
+Detuningslong_vec = []
+FittedCounts_vec = []
+
+Betas_vec = []
+ErrorBetas_vec = []
+Temp_vec = []
+ErrorTemp_vec = []
+
+DetuningsUV_vec = []
+ErrorDetuningsUV_vec = []
+
+for selectedcurve in SelectedCurveVec: 
+
+#selectedcurve = 2 #IMPORTANTE: SELECCIONA LA MEDICION
+
+    FreqsDR = Freqs[0]
+    CountsDR = CountsSplit[0][selectedcurve]
+
+
+    if selectedcurve==1:
+        CountsDR[100]=0.5*(CountsDR[99]+CountsDR[101])
+        CountsDR[105]=0.5*(CountsDR[104]+CountsDR[106])
+    if selectedcurve==2:
+        CountsDR[67]=0.5*(CountsDR[66]+CountsDR[68])
+        CountsDR[71]=0.5*(CountsDR[70]+CountsDR[72])
+    if selectedcurve==6:
+        CountsDR[1]=0.5*(CountsDR[0]+CountsDR[2])
+        CountsDR[76]=0.5*(CountsDR[75]+CountsDR[77])
+    if selectedcurve==7:
+        CountsDR[117]=0.5*(CountsDR[116]+CountsDR[118])
+        
+    freqslong = np.arange(min(FreqsDR), max(FreqsDR)+FreqsDR[1]-FreqsDR[0], 0.1*(FreqsDR[1]-FreqsDR[0]))
+    CircPr = 1
+    alpha = 0
+    
+
+    def FitEIT_MM_single(Freqs, offset, DetDoppler, SG, SP, SCALE1, OFFSET, BETA1, TEMP, U, plot=False):
+    #def FitEIT_MM(freqs, SG, SP, SCALE1, OFFSET, BETA1):
+        #BETA = 1.8
+        # SG = 0.6
+        # SP = 8.1
+        # TEMP = 0.2e-3
+    
+        freqs = [2*f*1e-6-offset for f in Freqs]
+      
+        Detunings, Fluorescence1 = PerformExperiment_8levels_MM(SG, SP, gPS, gPD, DetDoppler, U, DopplerLaserLinewidth, ProbeLaserLinewidth, TEMP, alpha, phidoppler, titadoppler, phiprobe, titaprobe,  BETA1, drivefreq, min(freqs), max(freqs)+(freqs[1]-freqs[0]), freqs[1]-freqs[0], circularityprobe=CircPr, plot=False, solvemode=1, detpvec=None)
+    
+        ScaledFluo1 = np.array([f*SCALE1 + OFFSET for f in Fluorescence1])
+        if plot:
+            return ScaledFluo1, Detunings
+        else:
+            return ScaledFluo1
+        #return ScaledFluo1
+    
+    do_fit = True
+    if do_fit:
+        popt_3_SA, pcov_3_SA = curve_fit(FitEIT_MM_single, FreqsDR, CountsDR, p0=[430, -25, 0.9, 6.2, 3e4, 1.34e3, 2, (np.pi**2)*1e-3, 32e6], bounds=((0, -50, 0, 0, 0, 0, 0, 0, 25e6), (1000, 0, 2, 20, 5e4, 5e4, 10, (np.pi**2)*10e-3, 40e6)))
+        
+        popt_SA_vec.append(popt_3_SA)
+        pcov_SA_vec.append(pcov_3_SA)
+    
+        FittedEITpi_3_SA_short, Detunings_3_SA_short = FitEIT_MM_single(FreqsDR, *popt_3_SA, plot=True)
+        freqslong = np.arange(min(FreqsDR), max(FreqsDR)+FreqsDR[1]-FreqsDR[0], 0.1*(FreqsDR[1]-FreqsDR[0]))
+        FittedEITpi_3_SA_long, Detunings_3_SA_long = FitEIT_MM_single(freqslong, *popt_3_SA, plot=True)
+    
+    DetuningsUV_vec.append(popt_3_SA[1])
+    ErrorDetuningsUV_vec.append(np.sqrt(pcov_3_SA[1,1]))
+    
+    Betas_vec.append(popt_3_SA[6])
+    ErrorBetas_vec.append(np.sqrt(pcov_3_SA[6,6]))
+    Temp_vec.append(popt_3_SA[7])
+    ErrorTemp_vec.append(np.sqrt(pcov_3_SA[7,7]))
+    
+    Detuningsshort_vec.append(Detunings_3_SA_short)
+    Counts_vec.append(CountsDR)
+    Detuningslong_vec.append(Detunings_3_SA_long)
+    FittedCounts_vec.append(FittedEITpi_3_SA_long)
+    
+    
+    plt.figure()
+    plt.errorbar(Detunings_3_SA_short, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
+    plt.plot(Detunings_3_SA_long, FittedEITpi_3_SA_long, color='darkolivegreen', linewidth=3, label=f'med {selectedcurve}')
+    #plt.title(f'Sdop: {round(popt[0], 2)}, Spr: {round(popt[1], 2)}, T: {round(popt[2]*1e3, 2)} mK, detDop: {DetDoppler} MHz')
+    plt.xlabel('Detuning (MHz)')
+    plt.ylabel('Counts')
+    plt.legend(loc='upper left', fontsize=20)
+    plt.grid()
+    
+    print(f'listo med {selectedcurve}')
+    print(popt_3_SA)
+
+    
+#%%
+"""
+Grafico distintas variables que salieron del SUper ajuste
+"""
+
+import seaborn as sns
+paleta = sns.color_palette("rocket")
+
+medfin = 12
+
+
+voltages_dcA = Voltages[0][1:medfin]
+
+def lineal(x,a,b):
+    return a*x+b
+
+def hiperbola(x,a,b,c,x0):
+    return a*np.sqrt(((x-x0)**2+c**2))+b
+
+hiperbola_or_linear = True
+
+if hiperbola_or_linear:
+    popthip,pcovhip = curve_fit(hiperbola,voltages_dcA,Betas_vec[:medfin-1],p0=(100,0.1,1,-0.15))
+    
+    xhip = np.linspace(-0.23,0.005,200)
+    
+    plt.figure()
+    plt.errorbar(voltages_dcA,Betas_vec[0:medfin-1],yerr=ErrorBetas_vec[:medfin-1],fmt='o',capsize=5,markersize=5,color=paleta[1])
+    plt.plot(xhip,hiperbola(xhip,*popthip))
+    plt.xlabel('Endcap voltage (V)')
+    plt.ylabel('Modulation factor')
+    plt.grid()
+
+else:
+    poptini,pcovini = curve_fit(lineal,voltages_dcA[0:3],Betas_vec[0:3])
+    poptfin,pcovfin = curve_fit(lineal,voltages_dcA[4:],Betas_vec[4:])
+    
+    minimum_voltage = -(poptini[1]-poptfin[1])/(poptini[0]-poptfin[0]) #voltaje donde se intersectan las rectas, es decir, donde deberia estar el minimo de micromocion
+    minimum_modulationfactor = lineal(minimum_voltage,*poptini) #es lo mismo si pongo *poptfin
+    
+    xini = np.linspace(-0.23,-0.13,100)
+    xfin = np.linspace(-0.15,0.005,100)
+    
+    plt.figure()
+    plt.errorbar(voltages_dcA,Betas_vec,yerr=ErrorBetas_vec,fmt='o',capsize=5,markersize=5,color=paleta[1])
+    plt.plot(xini,lineal(xini,*poptini))
+    plt.plot(xfin,lineal(xfin,*poptfin))
+    plt.axvline(minimum_voltage,linestyle='dashed',color='grey')
+    plt.xlabel('Endcap voltage (V)')
+    plt.ylabel('Modulation factor')
+    plt.grid()
+
+
+print([t*1e3 for t in Temp_vec])
+
+plt.figure()
+plt.errorbar(voltages_dcA,[t*1e3 for t in Temp_vec[:medfin-1]],yerr=[t*1e3 for t in ErrorTemp_vec[:medfin-1]],fmt='o',capsize=5,markersize=5,color=paleta[3])
+#plt.axvline(minimum_voltage,linestyle='dashed',color='grey')
+plt.axhline(0.538)
+plt.xlabel('Endcap voltage (V)')
+plt.ylabel('Temperature (mK)')
+plt.grid()
+#plt.ylim(0,2)
+
+#%%
+"""
+Ahora hago un ajuste con una hiperbola porque tiene mas sentido, por el hecho
+de que en el punto optimo el ion no esta en el centro de la trampa
+sino que esta a una distancia d
+"""
+def hiperbola(x,a,b,c,x0):
+    return a*np.sqrt(((x-x0)**2+c**2))+b
+
+popthip,pcovhip = curve_fit(hiperbola,voltages_dcA,Betas_vec[:10],p0=(100,0.1,1,-0.15))
+
+xhip = np.linspace(-0.23,0.005,200)
+
+plt.figure()
+plt.errorbar(voltages_dcA,Betas_vec[:10],yerr=ErrorBetas_vec[:10],fmt='o',capsize=5,markersize=5,color=paleta[1])
+plt.plot(xhip,hiperbola(xhip,*popthip))
+plt.xlabel('Endcap voltage (V)')
+plt.ylabel('Modulation factor')
+plt.grid()
+
+
+
+
+#%%
+from scipy.special import jv
+
+
+def expo(x,tau,A,B):
+    return A*np.exp(x/tau)+B
+
+def cuadratica(x,a,c):
+    return a*(x**2)+c
+
+def InverseMicromotionSpectra(beta, A, det, x0, gamma, B):
+    ftrap=22.1
+    #gamma=30
+    P = ((jv(0, beta)**2)/((((det-x0)**2)+(0.5*gamma)**2)**2))*(-2*(det-x0))
+    i = 1
+    #print(P)
+    while i <= 5:
+        P = P + (-2*(det-x0))*((jv(i, beta))**2)/(((((det-x0)+i*ftrap)**2)+(0.5*gamma)**2)**2) + (-2*(det-x0))*(((jv(-i, beta))**2)/((((det-x0)-i*ftrap)**2)+(0.5*gamma)**2)**2)
+        i = i + 1
+        #print(P)
+    #return 1/(A*P+B)
+    return 1/(A*P+B)
+
+
+def InverseMicromotionSpectra_raw(beta, A, det, B):
+    ftrap=22.1
+    gamma=21
+    P = ((jv(0, beta)**2)/((((det)**2)+(0.5*gamma)**2)**2))*(-2*(det))
+    i = 1
+    #print(P)
+    while i <= 3:
+        P = P + (-2*(det))*((jv(i, beta))**2)/(((((det)+i*ftrap)**2)+(0.5*gamma)**2)**2) + (-2*(det))*(((jv(-i, beta))**2)/((((det)-i*ftrap)**2)+(0.5*gamma)**2)**2)
+        i = i + 1
+        #print(P)
+    return A/P+B
+
+
+"""
+Temperatura vs beta con un ajuste exponencial 
+"""
+
+popt_exp, pcov_exp = curve_fit(expo,Betas_vec[:10],[t*1e3 for t in Temp_vec[:10]])
+popt_quad, pcov_quad = curve_fit(cuadratica,Betas_vec[:10],[t*1e3 for t in Temp_vec[:10]],p0=(1,10))
+#popt_rho22, pcov_rho22 = curve_fit(InverseMicromotionSpectra,Betas_vec,[t*1e3 for t in Temp_vec],p0=(10,10,-10,1,20)) #esto ajusta muy bien
+#popt_rho22, pcov_rho22 = curve_fit(InverseMicromotionSpectra,Betas_vec, [t*1e3 for t in Temp_vec],p0=(-10,-10,10,1,20)) #esto ajusta muy bien
+popt_rho22_raw, pcov_rho22_raw = curve_fit(InverseMicromotionSpectra_raw,Betas_vec[:10], [t*1e3 for t in Temp_vec[:10]],p0=(-10, -10, 1)) #esto ajusta muy bien
+
+
+print(popt_rho22_raw)
+
+betaslong = np.arange(0,2*2.7,0.01)
+
+print(f'Min temp predicted: {InverseMicromotionSpectra_raw(betaslong,*popt_rho22_raw)[100]}')
+
+plt.figure()
+plt.errorbar(Betas_vec[:10],[t*1e3 for t in Temp_vec[:10]],xerr=ErrorBetas_vec[:10], yerr=[t*1e3 for t in ErrorTemp_vec[:10]],fmt='o',capsize=5,markersize=5,color=paleta[3])
+#plt.plot(betaslong,expo(betaslong,*popt_exp),label='Ajuste exponencial')
+#plt.plot(betaslong,cuadratica(betaslong,*popt_quad),label='Ajuste cuadratico')
+#plt.plot(betaslong,InverseMicromotionSpectra(betaslong,*popt_rho22),label='Ajuste cuadratico')
+plt.plot(betaslong,InverseMicromotionSpectra_raw(betaslong,*popt_rho22_raw),label='Ajuste cuadratico')
+
+#plt.axvline(minimum_voltage,linestyle='dashed',color='grey')
+#plt.axhline(0.538)
+plt.xlabel('Modulation factor')
+plt.ylabel('Temperature (mK)')
+plt.grid()
+
+
+#%%
+"""
+Esto no es del super ajuste sino de los ajustes anteriores en donde DetDoppler y offset son puestos a mano
+Aca grafico los betas con su error en funcion de la tension variada.
+Ademas, hago ajuste lineal para primeros y ultimos puntos, ya que espero que
+si la tension hace que la posicion del ion varie linealmente, el beta varia proporcional a dicha posicion.
+"""
+
+import seaborn as sns
+
+def lineal(x,a,b):
+    return a*x+b
+
+paleta = sns.color_palette("rocket")
+
+betavector = [beta1,beta2,beta3,beta4,beta5,beta6,beta7,beta8,beta9]
+errorbetavector = [errorbeta1,errorbeta2,errorbeta3,errorbeta4,errorbeta5,errorbeta6,errorbeta7,errorbeta8,errorbeta9]
+voltages_dcA = Voltages[0][1:10]
+
+poptini,pcovini = curve_fit(lineal,voltages_dcA[0:3],betavector[0:3])
+poptfin,pcovfin = curve_fit(lineal,voltages_dcA[4:],betavector[4:])
+
+minimum_voltage = -(poptini[1]-poptfin[1])/(poptini[0]-poptfin[0]) #voltaje donde se intersectan las rectas, es decir, donde deberia estar el minimo de micromocion
+minimum_modulationfactor = lineal(minimum_voltage,*poptini) #es lo mismo si pongo *poptfin
+
+xini = np.linspace(-0.23,-0.13,100)
+xfin = np.linspace(-0.15,0.005,100)
+
+plt.figure()
+plt.errorbar(voltages_dcA,betavector,yerr=errorbetavector,fmt='o',capsize=5,markersize=5,color=paleta[1])
+plt.plot(xini,lineal(xini,*poptini))
+plt.plot(xfin,lineal(xfin,*poptfin))
+plt.axvline(minimum_voltage,linestyle='dashed',color='grey')
+plt.xlabel('Endcap voltage (V)')
+plt.ylabel('Modulation factor')
+plt.grid()
+
+
+#%%
+"""
+Aca veo la temperatura del ion en funcion del voltaje del endcap, ya que
+al cambiar la cantidad de micromocion, cambia la calidad del enfriado
+"""
+tempvector = np.array([temp1,temp2,temp3,temp4,temp5,temp6,temp7,temp8,temp9])*1e3
+errortempvector = np.array([errortemp1,errortemp2,errortemp3,errortemp4,errortemp5,errortemp6,errortemp7,errortemp8,errortemp9])*1e3
+
+voltages_dcA = Voltages[0][1:10]
+
+plt.figure()
+plt.errorbar(voltages_dcA,tempvector,yerr=errortempvector,fmt='o',capsize=5,markersize=5,color=paleta[3])
+plt.axvline(minimum_voltage,linestyle='dashed',color='grey')
+plt.xlabel('Endcap voltage (V)')
+plt.ylabel('Temperature (mK)')
+plt.grid()
+plt.ylim(0,2)
+
+
+#%%
+"""
+Por las dudas, temperatura en funcion de beta
+"""
+plt.figure()
+plt.errorbar(betavector,tempvector,yerr=errortempvector,xerr=errorbetavector,fmt='o',capsize=5,markersize=5)
+plt.xlabel('Modulation factor')
+plt.ylabel('Temperature (mK)')
+plt.grid()
+
+
+
+#%%
+"""
+Si quiero ver algun parametro del ajuste puntual. el orden es: 0:SG, 1:SP, 2:SCALE1, 3:OFFSET
+"""
+ki=2
+plt.errorbar(np.arange(0,9,1),[popt_1[ki],popt_2[ki],popt_3[ki],popt_4[ki],popt_5[ki],popt_6[ki],popt_7[ki],popt_8[ki],popt_9[ki]],yerr=[np.sqrt(pcov_1[ki,ki]),np.sqrt(pcov_2[ki,ki]),np.sqrt(pcov_3[ki,ki]),np.sqrt(pcov_4[ki,ki]),np.sqrt(pcov_5[ki,ki]),np.sqrt(pcov_6[ki,ki]),np.sqrt(pcov_7[ki,ki]),np.sqrt(pcov_8[ki,ki]),np.sqrt(pcov_9[ki,ki])], fmt='o',capsize=3,markersize=3)
+
+
+#%%
+
+if False:
+    GUARDAR = {}
+    for var in [ kk for kk in globals().keys() if kk.startswith('pop') ]:
+        print(var)
+        GUARDAR[var] = globals()[var]
+    print('')
+    for var in [ kk for kk in globals().keys() if kk.startswith('pcov') ]:
+        print(var)
+        GUARDAR[var] = globals()[var]
+
+    print('')
+    for var in [ kk for kk in globals().keys() if kk.startswith('Fitted') ]:
+        print(var)
+        GUARDAR[var] = globals()[var]
+    print('')
+    for var in [ kk for kk in globals().keys() if kk.endswith('_vec') ]:
+        print(var)
+        GUARDAR[var] = globals()[var]
+
+
+    np.savez('PARAMETROS.npz', **GUARDAR )
+
+

analisis/plots/20231123_CPTconmicromocion3/lolo_CPT_plotter_20231123.py

@@ -920,6 +920,9 @@ plt.grid()
 #%%
 """
 AHORA INTENTO SUPER AJUSTES O SEA CON OFFSETXPI Y DETDOPPLER INCLUIDOS
+
+La 0 no ajusta bien incluso con todos los parametros libres
+De la 1 a la 11 ajustan bien
 """
 
 
@@ -963,8 +966,8 @@ alpha = 0
 drivefreq = 2*np.pi*22.135*1e6
 
 
-SelectedCurveVec = [1,2,3,4,5,6,7,8,9]
-#SelectedCurveVec = [9]
+SelectedCurveVec = [1,2,3,4,5,6,7,8,9,10,11]
+#SelectedCurveVec = [10]
 
 if not 'popt_SA_vec' in globals().keys() or len(popt_SA_vec)==0:
 
@@ -1152,9 +1155,42 @@ plt.grid()
 
 
 #%%
+from scipy.special import jv
+
+
 def expo(x,tau,A,B):
     return A*np.exp(x/tau)+B
 
+def cuadratica(x,a,c):
+    return a*(x**2)+c
+
+def InverseMicromotionSpectra(beta, A, det, x0, gamma, B):
+    ftrap=22.1
+    #gamma=30
+    P = ((jv(0, beta)**2)/((((det-x0)**2)+(0.5*gamma)**2)**2))*(-2*(det-x0))
+    i = 1
+    #print(P)
+    while i <= 5:
+        P = P + (-2*(det-x0))*((jv(i, beta))**2)/(((((det-x0)+i*ftrap)**2)+(0.5*gamma)**2)**2) + (-2*(det-x0))*(((jv(-i, beta))**2)/((((det-x0)-i*ftrap)**2)+(0.5*gamma)**2)**2)
+        i = i + 1
+        #print(P)
+    #return 1/(A*P+B)
+    return 1/(A*P+B)
+
+
+def InverseMicromotionSpectra_raw(beta, A, det, B):
+    ftrap=22.1
+    gamma=21
+    P = ((jv(0, beta)**2)/((((det)**2)+(0.5*gamma)**2)**2))*(-2*(det))
+    i = 1
+    #print(P)
+    while i <= 3:
+        P = P + (-2*(det))*((jv(i, beta))**2)/(((((det)+i*ftrap)**2)+(0.5*gamma)**2)**2) + (-2*(det))*(((jv(-i, beta))**2)/((((det)-i*ftrap)**2)+(0.5*gamma)**2)**2)
+        i = i + 1
+        #print(P)
+    return A/P+B
+
+
 """
 Temperatura vs
 """

analisis/plots/20231123_CPTconmicromocion3/lolo_analisis_superajuste.py

+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+"""
+Ploteo de datos y ajustes
+@author: lolo
+"""
+
+
+import h5py
+import matplotlib.pyplot as plt
+import numpy as np
+# import sys
+# import re
+# import ast
+from scipy.optimize import curve_fit
+# import os
+# from scipy import interpolate
+
+from numba import jit,njit
+
+from time import time
+
+#%% Importaciones extra
+
+# /home/lolo/Dropbox/marce/LIAF/Trampa_anular/artiq_experiments/analisis/plots/20231123_CPTconmicromocion3/Data/EITfit/MM_eightLevel_2repumps_AnalysisFunctions.py
+
+from Data.EITfit.lolo_modelo_full_8niveles import PerformExperiment_8levels_MM
+
+
+# PARAMETROS = np.load('PARAMETROS.npz',allow_pickle=True)
+# for var_name in PARAMETROS.keys():
+#     globals()[var_name] = PARAMETROS[var_name]
+#     print(f'loaded: {var_name}')
+
+
+#%%
+
+"""
+Primero tengo mediciones de espectros cpt de un ion variando la tension dc_A
+"""
+
+#C:\Users\Usuario\Documents\artiq\artiq_experiments\analisis\plots\20220106_CPT_DosLaseres_v08_TISA_DR\Data
+# os.chdir('../20231123_CPTconmicromocion3/Data/')
+
+folder = '../20231123_CPTconmicromocion3/Data/'
+CPT_FILES = f"""
+{folder}/000016262-IR_Scan_withcal_optimized
+{folder}/000016239-IR_Scan_withcal_optimized
+{folder}/000016240-IR_Scan_withcal_optimized
+{folder}/000016241-IR_Scan_withcal_optimized
+{folder}/000016244-IR_Scan_withcal_optimized
+{folder}/000016255-IR_Scan_withcal_optimized
+{folder}/000016256-IR_Scan_withcal_optimized
+{folder}/000016257-IR_Scan_withcal_optimized
+"""
+
+
+def SeeKeys(files):
+    for i, fname in enumerate(files.split()):
+        data = h5py.File(fname+'.h5', 'r') # Leo el h5: Recordar que nuestros datos estan en 'datasets'
+        print(fname)
+        print(list(data['datasets'].keys()))
+
+print(SeeKeys(CPT_FILES))
+
+#carpeta pc nico labo escritorio:
+#C:\Users\Usuario\Documents\artiq\artiq_experiments\analisis\plots\20211101_CPT_DosLaseres_v03\Data
+
+Counts = []
+Freqs = []
+
+AmpTisa = []
+UVCPTAmp = []
+No_measures = []
+Voltages = []
+
+for i, fname in enumerate(CPT_FILES.split()):
+    print(str(i) + ' - ' + fname)
+    #print(fname)
+    data = h5py.File(fname+'.h5', 'r') # Leo el h5: Recordar que nuestros datos estan en 'datasets'
+
+    # Aca hago algo repugnante para poder levantar los strings que dejamos
+    # que además tenian un error de tipeo al final. Esto no deberá ser necesario
+    # cuando se solucione el error este del guardado.
+    Freqs.append(np.array(data['datasets']['IR1_Frequencies']))
+    Counts.append(np.array(data['datasets']['data_array']))
+    #AmpTisa.append(np.array(data['datasets']['TISA_CPT_amp']))
+    UVCPTAmp.append(np.array(data['datasets']['UV_CPT_amp']))
+    No_measures.append(np.array(data['datasets']['no_measures']))
+    Voltages.append(np.array(data['datasets']['scanning_voltages']))
+
+def Split(array,n):
+    length=len(array)/n
+    splitlist = []
+    jj = 0
+    while jj<length:
+        partial = []
+        ii = 0
+        while ii < n:
+            partial.append(array[jj*n+ii])
+            ii = ii + 1
+        splitlist.append(partial)
+        jj = jj + 1
+    return splitlist
+
+
+CountsSplit = []
+CountsSplit.append(Split(Counts[0],len(Freqs[0])))
+
+
+CountsSplit_2ions = []
+CountsSplit_2ions.append(Split(Counts[4],len(Freqs[4])))
+
+
+
+#%% Definiciones de Numba
+
+@jit
+def FitEIT_MM_single(Freqs, offset, DetDoppler, SG, SP, SCALE1, OFFSET, BETA1, TEMP, plot=False):
+    #BETA = 1.8
+    # SG = 0.6
+    # SP = 8.1
+    # TEMP = 0.2e-3
+
+    freqs = [2*f*1e-6-offset for f in Freqs]
+
+    Detunings, Fluorescence1 = PerformExperiment_8levels_MM(SG, SP, gPS, gPD, DetDoppler, u, 
+                                                            DopplerLaserLinewidth, ProbeLaserLinewidth, 
+                                                            TEMP, alpha, phidoppler, titadoppler, 
+                                                            phiprobe, titaprobe,  BETA1, drivefreq, 
+                                                            min(freqs), max(freqs)+(freqs[1]-freqs[0]), 
+                                                            freqs[1]-freqs[0], circularityprobe=CircPr, 
+                                                            plot=False, solvemode=1, detpvec=None)
+
+    ScaledFluo1 = np.array([f*SCALE1 + OFFSET for f in Fluorescence1])
+
+    return ScaledFluo1
+
+
+
+@jit
+def FitEIT_MM_single_plot(Freqs, offset, DetDoppler, SG, SP, SCALE1, OFFSET, BETA1, TEMP):
+    #BETA = 1.8
+    # SG = 0.6
+    # SP = 8.1
+    # TEMP = 0.2e-3
+
+    freqs = [2*f*1e-6-offset for f in Freqs]
+
+    Detunings, Fluorescence1 = PerformExperiment_8levels_MM(SG, SP, gPS, gPD, DetDoppler, u, 
+                                                            DopplerLaserLinewidth, ProbeLaserLinewidth, 
+                                                            TEMP, alpha, phidoppler, titadoppler, 
+                                                            phiprobe, titaprobe,  BETA1, drivefreq, 
+                                                            min(freqs), max(freqs)+(freqs[1]-freqs[0]), 
+                                                            freqs[1]-freqs[0], circularityprobe=CircPr, 
+                                                            plot=False, solvemode=1, detpvec=None)
+
+    ScaledFluo1 = np.array([f*SCALE1 + OFFSET for f in Fluorescence1])
+    return ScaledFluo1, Detunings
+
+
+
+#%%
+"""
+AHORA INTENTO SUPER AJUSTES O SEA CON OFFSETXPI Y DETDOPPLER INCLUIDOS
+
+La 0 no ajusta bien incluso con todos los parametros libres
+De la 1 a la 11 ajustan bien
+"""
+
+
+#from EITfit.MM_eightLevel_2repumps_AnalysisFunctions import PerformExperiment_8levels
+from scipy.optimize import curve_fit
+
+"""
+SUPER AJUSTE (SA)
+"""
+
+phidoppler, titadoppler = 0, 90
+phirepump, titarepump = 0,  0
+phiprobe = 0
+titaprobe = 90
+
+Temp = 0.5e-3
+
+sg = 0.544
+sp = 4.5
+sr = 0
+DetRepump = 0
+
+
+lw = 0.1
+DopplerLaserLinewidth, RepumpLaserLinewidth, ProbeLaserLinewidth = lw, lw, lw #ancho de linea de los laseres
+
+
+u = 32.5e6
+
+#B = (u/(2*np.pi))/c
+
+correccion = 13
+
+#DetDoppler = -11.5-correccion
+
+gPS, gPD, = 2*np.pi*21.58e6, 2*np.pi*1.35e6
+alpha = 0
+
+
+drivefreq = 2*np.pi*22.135*1e6
+
+
+SelectedCurveVec = [1,2,3,4,5,6,7,8,9,10,11]
+#SelectedCurveVec = [10]
+
+CircPr = 1
+alpha = 0
+
+t0 = time()
+
+if not 'popt_SA_vec' in globals().keys() or len(popt_SA_vec)==0:
+
+    popt_SA_vec = []
+    pcov_SA_vec = []
+    Detuningsshort_vec = []
+    Counts_vec = []
+    Detuningslong_vec = []
+    FittedCounts_vec = []
+
+    Betas_vec = []
+    ErrorBetas_vec = []
+    Temp_vec = []
+    ErrorTemp_vec = []
+
+    DetuningsUV_vec = []
+    ErrorDetuningsUV_vec = []
+
+    for selectedcurve in SelectedCurveVec:
+        print(f"{round(time()-t0,1):6.1f}: Procesando la curva {selectedcurve}")
+        #selectedcurve = 2 #IMPORTANTE: SELECCIONA LA MEDICION
+
+        FreqsDR = Freqs[0]
+        CountsDR = CountsSplit[0][selectedcurve]
+
+
+        if selectedcurve==1:
+            CountsDR[100]=0.5*(CountsDR[99]+CountsDR[101])
+            CountsDR[105]=0.5*(CountsDR[104]+CountsDR[106])
+        if selectedcurve==2:
+            CountsDR[67]=0.5*(CountsDR[66]+CountsDR[68])
+            CountsDR[71]=0.5*(CountsDR[70]+CountsDR[72])
+        if selectedcurve==6:
+            CountsDR[1]=0.5*(CountsDR[0]+CountsDR[2])
+            CountsDR[76]=0.5*(CountsDR[75]+CountsDR[77])
+        if selectedcurve==7:
+            CountsDR[117]=0.5*(CountsDR[116]+CountsDR[118])
+
+        freqslong = np.arange(min(FreqsDR), max(FreqsDR)+FreqsDR[1]-FreqsDR[0], 0.1*(FreqsDR[1]-FreqsDR[0]))
+
+
+        if True:
+            popt_3_SA, pcov_3_SA = curve_fit(FitEIT_MM_single, FreqsDR, CountsDR, p0=[430, -25, 0.9, 6.2, 3e4, 1.34e3, 2, (np.pi**2)*1e-3], bounds=((0, -50, 0, 0, 0, 0, 0, 0), (1000, 0, 2, 20, 5e4, 5e4, 10, (np.pi**2)*10e-3)))
+
+            popt_SA_vec.append(popt_3_SA)
+            pcov_SA_vec.append(pcov_3_SA)
+
+            FittedEITpi_3_SA_short, Detunings_3_SA_short = FitEIT_MM_single_plot(FreqsDR, *popt_3_SA)
+            freqslong = np.arange(min(FreqsDR), max(FreqsDR)+FreqsDR[1]-FreqsDR[0], 0.1*(FreqsDR[1]-FreqsDR[0]))
+            FittedEITpi_3_SA_long, Detunings_3_SA_long = FitEIT_MM_single_plot(freqslong, *popt_3_SA)
+
+        DetuningsUV_vec.append(popt_3_SA[1])
+        ErrorDetuningsUV_vec.append(np.sqrt(pcov_3_SA[1,1]))
+
+        Betas_vec.append(popt_3_SA[6])
+        ErrorBetas_vec.append(np.sqrt(pcov_3_SA[6,6]))
+        Temp_vec.append(popt_3_SA[7])
+        ErrorTemp_vec.append(np.sqrt(pcov_3_SA[7,7]))
+
+        Detuningsshort_vec.append(Detunings_3_SA_short)
+        Counts_vec.append(CountsDR)
+        Detuningslong_vec.append(Detunings_3_SA_long)
+        FittedCounts_vec.append(FittedEITpi_3_SA_long)
+
+
+#%%
+
+# tmp_datos=(Detuningsshort_vec,Counts_vec,Detuningslong_vec,FittedCounts_vec,SelectedCurveVec)
+
+# for Detunings_3_SA_short,CountsDR,Detunings_3_SA_long,FittedEITpi_3_SA_long,selectedcurve in zip(*tmp_datos):
+#     plt.figure()
+#     plt.errorbar(Detunings_3_SA_short, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
+#     plt.plot(Detunings_3_SA_long, FittedEITpi_3_SA_long, color='darkolivegreen', linewidth=3, label=f'med {selectedcurve}')
+#     #plt.title(f'Sdop: {round(popt[0], 2)}, Spr: {round(popt[1], 2)}, T: {round(popt[2]*1e3, 2)} mK, detDop: {DetDoppler} MHz')
+#     plt.xlabel('Detuning (MHz)')
+#     plt.ylabel('Counts')
+#     plt.legend(loc='upper left', fontsize=20)
+#     plt.grid()
+
+#     print(f'listo med {selectedcurve}')
+#     print(popt_3_SA)
+
+
+fig, axx = plt.subplots( 3,4, figsize=(13,8) ,  constrained_layout=True, sharex=True , sharey=True )
+fig.set_constrained_layout_pads(w_pad=2/72, h_pad=2/72, hspace=0, wspace=0)
+
+tmp_datos=(Detuningsshort_vec,Counts_vec,Detuningslong_vec,FittedCounts_vec,SelectedCurveVec,axx.flatten())
+
+for Detunings_3_SA_short,CountsDR,Detunings_3_SA_long,FittedEITpi_3_SA_long,selectedcurve,ax in zip(*tmp_datos):
+    ax.errorbar(Detunings_3_SA_short, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.3, capsize=2, markersize=2)
+    ax.plot(Detunings_3_SA_long, FittedEITpi_3_SA_long, color='black', linewidth=2, label=f'med {selectedcurve}', alpha=0.7)
+    #plt.title(f'Sdop: {round(popt[0], 2)}, Spr: {round(popt[1], 2)}, T: {round(popt[2]*1e3, 2)} mK, detDop: {DetDoppler} MHz')
+    # ax.set_xlabel('Detuning (MHz)')
+    # ax.set_ylabel('Counts')
+    ax.legend(loc='upper left', fontsize=12)
+    ax.grid(True, ls=":")
+
+    print(f'listo med {selectedcurve}')
+    print(popt_3_SA)
+
+
+for ax in axx[:,0]:
+    ax.set_ylabel('Counts')
+
+for ax in axx[-1,:]:
+    ax.set_xlabel('Detuning (MHz)')
+
+#%%
+"""
+Grafico distintas variables que salieron del SUper ajuste
+"""
+
+import seaborn as sns
+paleta = sns.color_palette("rocket")
+
+voltages_dcA = Voltages[0][1:10]
+
+def lineal(x,a,b):
+    return a*x+b
+
+def hiperbola(x,a,b,c,x0):
+    return a*np.sqrt(((x-x0)**2+c**2))+b
+
+hiperbola_or_linear = True
+
+if hiperbola_or_linear:
+    popthip,pcovhip = curve_fit(hiperbola,voltages_dcA,Betas_vec,p0=(100,0.1,1,-0.15))
+
+    xhip = np.linspace(-0.23,0.005,200)
+
+    plt.figure()
+    plt.errorbar(voltages_dcA,Betas_vec,yerr=ErrorBetas_vec,fmt='o',capsize=5,markersize=5,color=paleta[1])
+    plt.plot(xhip,hiperbola(xhip,*popthip))
+    plt.xlabel('Endcap voltage (V)')
+    plt.ylabel('Modulation factor')
+    plt.grid()
+
+else:
+    poptini,pcovini = curve_fit(lineal,voltages_dcA[0:3],Betas_vec[0:3])
+    poptfin,pcovfin = curve_fit(lineal,voltages_dcA[4:],Betas_vec[4:])
+
+    minimum_voltage = -(poptini[1]-poptfin[1])/(poptini[0]-poptfin[0]) #voltaje donde se intersectan las rectas, es decir, donde deberia estar el minimo de micromocion
+    minimum_modulationfactor = lineal(minimum_voltage,*poptini) #es lo mismo si pongo *poptfin
+
+    xini = np.linspace(-0.23,-0.13,100)
+    xfin = np.linspace(-0.15,0.005,100)
+
+    plt.figure()
+    plt.errorbar(voltages_dcA,Betas_vec,yerr=ErrorBetas_vec,fmt='o',capsize=5,markersize=5,color=paleta[1])
+    plt.plot(xini,lineal(xini,*poptini))
+    plt.plot(xfin,lineal(xfin,*poptfin))
+    plt.axvline(minimum_voltage,linestyle='dashed',color='grey')
+    plt.xlabel('Endcap voltage (V)')
+    plt.ylabel('Modulation factor')
+    plt.grid()
+
+
+print([t*1e3 for t in Temp_vec])
+
+plt.figure()
+plt.errorbar(voltages_dcA,[t*1e3 for t in Temp_vec],yerr=[t*1e3 for t in ErrorTemp_vec],fmt='o',capsize=5,markersize=5,color=paleta[3])
+# plt.axvline(minimum_voltage,linestyle='dashed',color='grey')
+print(f'\n\nTE FALTA DEFINIR LA VARIABLE minimum_voltage\n\n')
+plt.axhline(0.538)
+plt.xlabel('Endcap voltage (V)')
+plt.ylabel('Temperature (mK)')
+plt.grid()
+#plt.ylim(0,2)
+
+#%%
+"""
+Ahora hago un ajuste con una hiperbola porque tiene mas sentido, por el hecho
+de que en el punto optimo el ion no esta en el centro de la trampa
+sino que esta a una distancia d
+"""
+def hiperbola(x,a,b,c,x0):
+    return a*np.sqrt(((x-x0)**2+c**2))+b
+
+popthip,pcovhip = curve_fit(hiperbola,voltages_dcA,Betas_vec,p0=(100,0.1,1,-0.15))
+
+xhip = np.linspace(-0.23,0.005,200)
+
+plt.figure()
+plt.errorbar(voltages_dcA,Betas_vec,yerr=ErrorBetas_vec,fmt='o',capsize=5,markersize=5,color=paleta[1])
+plt.plot(xhip,hiperbola(xhip,*popthip))
+plt.xlabel('Endcap voltage (V)')
+plt.ylabel('Modulation factor')
+plt.grid()
+
+
+
+
+#%%
+from scipy.special import jv
+
+
+def expo(x,tau,A,B):
+    return A*np.exp(x/tau)+B
+
+def cuadratica(x,a,c):
+    return a*(x**2)+c
+
+def InverseMicromotionSpectra(beta, A, det, x0, gamma, B):
+    ftrap=22.1
+    #gamma=30
+    P = ((jv(0, beta)**2)/((((det-x0)**2)+(0.5*gamma)**2)**2))*(-2*(det-x0))
+    i = 1
+    #print(P)
+    while i <= 5:
+        P = P + (-2*(det-x0))*((jv(i, beta))**2)/(((((det-x0)+i*ftrap)**2)+(0.5*gamma)**2)**2) + (-2*(det-x0))*(((jv(-i, beta))**2)/((((det-x0)-i*ftrap)**2)+(0.5*gamma)**2)**2)
+        i = i + 1
+        #print(P)
+    #return 1/(A*P+B)
+    return 1/(A*P+B)
+
+
+def InverseMicromotionSpectra_raw(beta, A, det, B):
+    ftrap=22.1
+    gamma=21
+    P = ((jv(0, beta)**2)/((((det)**2)+(0.5*gamma)**2)**2))*(-2*(det))
+    i = 1
+    #print(P)
+    while i <= 3:
+        P = P + (-2*(det))*((jv(i, beta))**2)/(((((det)+i*ftrap)**2)+(0.5*gamma)**2)**2) + (-2*(det))*(((jv(-i, beta))**2)/((((det)-i*ftrap)**2)+(0.5*gamma)**2)**2)
+        i = i + 1
+        #print(P)
+    return A/P+B
+
+
+"""
+Temperatura vs beta con un ajuste exponencial 
+"""
+
+popt_exp, pcov_exp = curve_fit(expo,Betas_vec[:10],[t*1e3 for t in Temp_vec[:10]])
+popt_quad, pcov_quad = curve_fit(cuadratica,Betas_vec[:10],[t*1e3 for t in Temp_vec[:10]],p0=(1,10))
+#popt_rho22, pcov_rho22 = curve_fit(InverseMicromotionSpectra,Betas_vec,[t*1e3 for t in Temp_vec],p0=(10,10,-10,1,20)) #esto ajusta muy bien
+#popt_rho22, pcov_rho22 = curve_fit(InverseMicromotionSpectra,Betas_vec, [t*1e3 for t in Temp_vec],p0=(-10,-10,10,1,20)) #esto ajusta muy bien
+popt_rho22_raw, pcov_rho22_raw = curve_fit(InverseMicromotionSpectra_raw,Betas_vec[:10], [t*1e3 for t in Temp_vec[:10]],p0=(-10, -10, 1)) #esto ajusta muy bien
+
+
+print(popt_rho22_raw)
+
+betaslong = np.arange(0,2*2.7,0.01)
+
+print(f'Min temp predicted: {InverseMicromotionSpectra_raw(betaslong,*popt_rho22_raw)[100]}')
+
+plt.figure()
+plt.errorbar(Betas_vec[:10],[t*1e3 for t in Temp_vec[:10]],xerr=ErrorBetas_vec[:10], yerr=[t*1e3 for t in ErrorTemp_vec[:10]],fmt='o',capsize=5,markersize=5,color=paleta[3])
+#plt.plot(betaslong,expo(betaslong,*popt_exp),label='Ajuste exponencial')
+#plt.plot(betaslong,cuadratica(betaslong,*popt_quad),label='Ajuste cuadratico')
+#plt.plot(betaslong,InverseMicromotionSpectra(betaslong,*popt_rho22),label='Ajuste cuadratico')
+plt.plot(betaslong,InverseMicromotionSpectra_raw(betaslong,*popt_rho22_raw),label='Ajuste cuadratico')
+
+#plt.axvline(minimum_voltage,linestyle='dashed',color='grey')
+#plt.axhline(0.538)
+plt.xlabel('Modulation factor')
+plt.ylabel('Temperature (mK)')
+plt.grid()
+
+
+#%%
+"""
+Esto no es del super ajuste sino de los ajustes anteriores en donde DetDoppler y offset son puestos a mano
+Aca grafico los betas con su error en funcion de la tension variada.
+Ademas, hago ajuste lineal para primeros y ultimos puntos, ya que espero que
+si la tension hace que la posicion del ion varie linealmente, el beta varia proporcional a dicha posicion.
+"""
+
+import seaborn as sns
+
+def lineal(x,a,b):
+    return a*x+b
+
+paleta = sns.color_palette("rocket")
+
+betavector = [beta1,beta2,beta3,beta4,beta5,beta6,beta7,beta8,beta9]
+errorbetavector = [errorbeta1,errorbeta2,errorbeta3,errorbeta4,errorbeta5,errorbeta6,errorbeta7,errorbeta8,errorbeta9]
+voltages_dcA = Voltages[0][1:10]
+
+poptini,pcovini = curve_fit(lineal,voltages_dcA[0:3],betavector[0:3])
+poptfin,pcovfin = curve_fit(lineal,voltages_dcA[4:],betavector[4:])
+
+minimum_voltage = -(poptini[1]-poptfin[1])/(poptini[0]-poptfin[0]) #voltaje donde se intersectan las rectas, es decir, donde deberia estar el minimo de micromocion
+minimum_modulationfactor = lineal(minimum_voltage,*poptini) #es lo mismo si pongo *poptfin
+
+xini = np.linspace(-0.23,-0.13,100)
+xfin = np.linspace(-0.15,0.005,100)
+
+plt.figure()
+plt.errorbar(voltages_dcA,betavector,yerr=errorbetavector,fmt='o',capsize=5,markersize=5,color=paleta[1])
+plt.plot(xini,lineal(xini,*poptini))
+plt.plot(xfin,lineal(xfin,*poptfin))
+plt.axvline(minimum_voltage,linestyle='dashed',color='grey')
+plt.xlabel('Endcap voltage (V)')
+plt.ylabel('Modulation factor')
+plt.grid()
+
+
+#%%
+"""
+Aca veo la temperatura del ion en funcion del voltaje del endcap, ya que
+al cambiar la cantidad de micromocion, cambia la calidad del enfriado
+"""
+tempvector = np.array([temp1,temp2,temp3,temp4,temp5,temp6,temp7,temp8,temp9])*1e3
+errortempvector = np.array([errortemp1,errortemp2,errortemp3,errortemp4,errortemp5,errortemp6,errortemp7,errortemp8,errortemp9])*1e3
+
+voltages_dcA = Voltages[0][1:10]
+
+plt.figure()
+plt.errorbar(voltages_dcA,tempvector,yerr=errortempvector,fmt='o',capsize=5,markersize=5,color=paleta[3])
+plt.axvline(minimum_voltage,linestyle='dashed',color='grey')
+plt.xlabel('Endcap voltage (V)')
+plt.ylabel('Temperature (mK)')
+plt.grid()
+plt.ylim(0,2)
+
+
+#%%
+"""
+Por las dudas, temperatura en funcion de beta
+"""
+plt.figure()
+plt.errorbar(betavector,tempvector,yerr=errortempvector,xerr=errorbetavector,fmt='o',capsize=5,markersize=5)
+plt.xlabel('Modulation factor')
+plt.ylabel('Temperature (mK)')
+plt.grid()
+
+
+
+#%%
+"""
+Si quiero ver algun parametro del ajuste puntual. el orden es: 0:SG, 1:SP, 2:SCALE1, 3:OFFSET
+"""
+ki=2
+plt.errorbar(np.arange(0,9,1),[popt_1[ki],popt_2[ki],popt_3[ki],popt_4[ki],popt_5[ki],popt_6[ki],popt_7[ki],popt_8[ki],popt_9[ki]],yerr=[np.sqrt(pcov_1[ki,ki]),np.sqrt(pcov_2[ki,ki]),np.sqrt(pcov_3[ki,ki]),np.sqrt(pcov_4[ki,ki]),np.sqrt(pcov_5[ki,ki]),np.sqrt(pcov_6[ki,ki]),np.sqrt(pcov_7[ki,ki]),np.sqrt(pcov_8[ki,ki]),np.sqrt(pcov_9[ki,ki])], fmt='o',capsize=3,markersize=3)
+
+
+#%%
+
+if False:
+    GUARDAR = {}
+    for var in [ kk for kk in globals().keys() if kk.startswith('pop') ]:
+        print(var)
+        GUARDAR[var] = globals()[var]
+    print('')
+    for var in [ kk for kk in globals().keys() if kk.startswith('pcov') ]:
+        print(var)
+        GUARDAR[var] = globals()[var]
+
+    print('')
+    for var in [ kk for kk in globals().keys() if kk.startswith('Fitted') ]:
+        print(var)
+        GUARDAR[var] = globals()[var]
+    print('')
+    for var in [ kk for kk in globals().keys() if kk.endswith('_vec') ]:
+        print(var)
+        GUARDAR[var] = globals()[var]
+
+
+    np.savez('PARAMETROS.npz', **GUARDAR )