Skip to content

A
artiq_experiments
Commit b700e139 authored by Nicolas Nunez Barreto's avatar Nicolas Nunez Barreto
Browse files
Options

agrego codgos de analisis de cpt de blade trap

parent f9b573e5
Hide whitespace changes
Inline Side-by-side
Showing with 2 additions and 5949 deletions
analisis/plots/20250516_CPT_bladetrap/CPT_plotter_20231123.py → analisis/plots/20250516_CPT_bladetrap/CPT_bladetrap.py
View file @ b700e139
...@@ -331,7 +331,7 @@ for selectedcurve in SelectedCurveVec: ...@@ -331,7 +331,7 @@ for selectedcurve in SelectedCurveVec:
Counts_vec.append(CountsDR) Counts_vec.append(CountsDR)
Detuningslong_vec.append(Detunings_3_SA_long) Detuningslong_vec.append(Detunings_3_SA_long)
FittedCounts_vec.append(FittedEITpi_3_SA_long) FittedCounts_vec.append(FittedEITpi_3_SA_long)
#%%
plt.figure() 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.errorbar(Detunings_3_SA_short[:], CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
plt.plot([c-corri for c in Detunings_3_SA_long], FittedEITpi_3_SA_long, color='darkolivegreen', linewidth=3, label=f'med {selectedcurve}') plt.plot([c-corri for c in Detunings_3_SA_long], FittedEITpi_3_SA_long, color='darkolivegreen', linewidth=3, label=f'med {selectedcurve}')
...@@ -341,6 +341,7 @@ for selectedcurve in SelectedCurveVec: ...@@ -341,6 +341,7 @@ for selectedcurve in SelectedCurveVec:
# plt.legend(loc='upper left', fontsize=20) # plt.legend(loc='upper left', fontsize=20)
plt.grid() plt.grid()
print(f'listo med {selectedcurve}') print(f'listo med {selectedcurve}')
print(popt_3_SA) print(popt_3_SA)
analisis/plots/20250516_CPT_bladetrap/analisis_superajuste.py deleted 100644 → 0
View file @ f9b573e5
#!/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/20250516_CPT_bladetrap/lolo_CPT_plotter_20231123.py deleted 100644 → 0
View file @ f9b573e5
#!/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
# 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])))
#%%
"""
Ploteo la cpt de referencia / plotting the reference CPT
"""
jvec = [9] # de la 1 a la 9 vale la pena, despues no
drs = [390.5, 399.5, 406, 413.5]
drive=22.1
Frequencies = Freqs[0]
plt.figure()
i = 0
for j in jvec:
plt.errorbar([2*f*1e-6 for f in Frequencies], CountsSplit[0][j], yerr=np.sqrt(CountsSplit[0][j]), fmt='o', capsize=2, markersize=2)
i = i + 1
plt.xlabel('Frecuencia (MHz)')
plt.ylabel('counts')
plt.grid()
#for dr in drs:
# plt.axvline(dr)
#plt.axvline(dr+drive)
plt.legend()
#%%
#from EITfit.MM_eightLevel_2repumps_AnalysisFunctions import PerformExperiment_8levels
from scipy.optimize import curve_fit
import time
"""
MEDICION 1
"""
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
offsetxpi = 419+correccion+3*0.8
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
selectedcurve = 1
FreqsDR = [2*f*1e-6-offsetxpi for f in Freqs[0]]
CountsDR = CountsSplit[0][selectedcurve]
CountsDR[100]=0.5*(CountsDR[99]+CountsDR[101])
CountsDR[105]=0.5*(CountsDR[104]+CountsDR[106])
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, SG, SP, SCALE1, OFFSET, BETA1, TEMP):
#def FitEIT_MM(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#BETA = 1.8
# SG = 0.6
# SP = 8.1
# TEMP = 0.2e-3
t0 = time.time()
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)
# print(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], dict(circularityprobe=CircPr, plot=False, solvemode=1, detpvec=None))
print('Done, Total time: ', round((time.time()-t0), 2), "s")
ScaledFluo1 = np.array([f*SCALE1 + OFFSET for f in Fluorescence1])
return ScaledFluo1
#return ScaledFluo1
if not 'popt_1' in globals().keys():
popt_1, pcov_1 = curve_fit(FitEIT_MM_single, FreqsDR, CountsDR, p0=[0.9, 6.2, 3e4, 1.34e3, 2, 1e-3], bounds=((0, 0, 0, 0, 0, 0), (2, 20, 5e4, 5e4, 10, 10e-3)))
FittedEITpi_1 = FitEIT_MM_single(freqslong, *popt_1)
beta1 = popt_1[4]
errorbeta1 = np.sqrt(pcov_1[4,4])
temp1 = popt_1[5]
errortemp1 = np.sqrt(pcov_1[5,5])
plt.figure()
plt.errorbar(FreqsDR, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
plt.plot(freqslong, FittedEITpi_1, color='darkolivegreen', linewidth=3, label='med 1')
#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()
#%%
#from EITfit.MM_eightLevel_2repumps_AnalysisFunctions import PerformExperiment_8levels
from scipy.optimize import curve_fit
import time
"""
MEDICION 2
"""
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
offsetxpi = 419+correccion+1.6
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
selectedcurve = 2
FreqsDR = [2*f*1e-6-offsetxpi for f in Freqs[0]]
CountsDR = CountsSplit[0][selectedcurve]
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, SG, SP, SCALE1, OFFSET, BETA1, TEMP):
#def FitEIT_MM(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#BETA = 1.8
# SG = 0.6
# SP = 8.1
# TEMP = 0.2e-3
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
#return ScaledFluo1
if not 'popt_2' in globals().keys():
popt_2, pcov_2 = curve_fit(FitEIT_MM_single, FreqsDR, CountsDR, p0=[0.9, 6.2, 3e4, 1.34e3, 2, 1e-3], bounds=((0, 0, 0, 0, 0, 0), (2, 20, 5e4, 5e4, 10, 10e-3)))
FittedEITpi_2 = FitEIT_MM_single(freqslong, *popt_2)
beta2 = popt_2[4]
errorbeta2 = np.sqrt(pcov_2[4,4])
temp2 = popt_2[5]
errortemp2 = np.sqrt(pcov_2[5,5])
plt.figure()
plt.errorbar(FreqsDR, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
plt.plot(freqslong, FittedEITpi_2, color='darkolivegreen', linewidth=3, label='med 2')
#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()
#%%
#from EITfit.MM_eightLevel_2repumps_AnalysisFunctions import PerformExperiment_8levels
from scipy.optimize import curve_fit
import time
"""
MEDICION 3
"""
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
offsetxpi = 419+correccion+0.8
DetDoppler = -11.5-correccion
print(offsetxpi,DetDoppler)
gPS, gPD, = 2*np.pi*21.58e6, 2*np.pi*1.35e6
alpha = 0
drivefreq = 2*np.pi*22.135*1e6
selectedcurve = 3
FreqsDR = [2*f*1e-6-offsetxpi for f in Freqs[0]]
CountsDR = CountsSplit[0][selectedcurve]
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, SG, SP, SCALE1, OFFSET, BETA1, TEMP):
#def FitEIT_MM(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#BETA = 1.8
# SG = 0.6
# SP = 8.1
# TEMP = 0.2e-3
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
#return ScaledFluo1
if not 'popt_3' in globals().keys():
popt_3, pcov_3 = curve_fit(FitEIT_MM_single, FreqsDR, CountsDR, p0=[0.9, 6.2, 3e4, 1.34e3, 2, 1e-3], bounds=((0, 0, 0, 0, 0, 0), (2, 20, 5e4, 5e4, 10, 10e-3)))
FittedEITpi_3 = FitEIT_MM_single(freqslong, *popt_3)
beta3 = popt_3[4]
errorbeta3 = np.sqrt(pcov_3[4,4])
temp3 = popt_3[5]
errortemp3 = np.sqrt(pcov_3[5,5])
plt.figure()
plt.errorbar(FreqsDR, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
plt.plot(freqslong, FittedEITpi_3, color='darkolivegreen', linewidth=3, label='med 3')
#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()
#%%
#from EITfit.MM_eightLevel_2repumps_AnalysisFunctions import PerformExperiment_8levels
from scipy.optimize import curve_fit
import time
"""
MEDICION 4
"""
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
offsetxpi = 419+correccion
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
selectedcurve = 4
FreqsDR = [2*f*1e-6-offsetxpi for f in Freqs[0]]
CountsDR = CountsSplit[0][selectedcurve]
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, SG, SP, SCALE1, OFFSET, BETA1, TEMP):
#def FitEIT_MM(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#BETA = 1.8
# SG = 0.6
# SP = 8.1
# TEMP = 0.2e-3
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
#return ScaledFluo1
if not 'popt_4' in globals().keys():
popt_4, pcov_4 = curve_fit(FitEIT_MM_single, FreqsDR, CountsDR, p0=[0.9, 6.2, 3e4, 1.34e3, 2, 1e-3], bounds=((0, 0, 0, 0, 0, 0), (2, 20, 5e4, 5e4, 10, 10e-3)))
FittedEITpi_4 = FitEIT_MM_single(freqslong, *popt_4)
beta4 = popt_4[4]
errorbeta4 = np.sqrt(pcov_4[4,4])
temp4 = popt_4[5]
errortemp4 = np.sqrt(pcov_4[5,5])
plt.figure()
plt.errorbar(FreqsDR, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
plt.plot(freqslong, FittedEITpi_4, color='darkolivegreen', linewidth=3, label='med 4')
#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()
#%%
#from EITfit.MM_eightLevel_2repumps_AnalysisFunctions import PerformExperiment_8levels
from scipy.optimize import curve_fit
import time
"""
MEDICION 5
"""
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
offsetxpi = 419+correccion-1
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
selectedcurve = 5
FreqsDR = [2*f*1e-6-offsetxpi for f in Freqs[0]]
CountsDR = CountsSplit[0][selectedcurve]
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, SG, SP, SCALE1, OFFSET, BETA1, TEMP):
#def FitEIT_MM(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#BETA = 1.8
# SG = 0.6
# SP = 8.1
#TEMP = 0.2e-3
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
#return ScaledFluo1
if not 'popt_5' in globals().keys():
popt_5, pcov_5 = curve_fit(FitEIT_MM_single, FreqsDR, CountsDR, p0=[0.9, 6.2, 3e4, 1.34e3, 2, 1e-3], bounds=((0, 0, 0, 0, 0, 0), (2, 20, 5e4, 5e4, 10, 10e-3)))
FittedEITpi_5 = FitEIT_MM_single(freqslong, *popt_5)
beta5 = popt_5[4]
errorbeta5 = np.sqrt(pcov_5[4,4])
temp5 = popt_5[5]
errortemp5 = np.sqrt(pcov_5[5,5])
plt.figure()
plt.errorbar(FreqsDR, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
plt.plot(freqslong, FittedEITpi_5, color='darkolivegreen', linewidth=3, label='med 5')
#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()
#%%
#from EITfit.MM_eightLevel_2repumps_AnalysisFunctions import PerformExperiment_8levels
from scipy.optimize import curve_fit
import time
"""
MEDICION 6
"""
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
offsetxpi = 419+correccion-2.2
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
selectedcurve = 6
FreqsDR = [2*f*1e-6-offsetxpi for f in Freqs[0]]
CountsDR = CountsSplit[0][selectedcurve]
CountsDR[76]=0.5*(CountsDR[75]+CountsDR[77])
CountsDR[1]=0.5*(CountsDR[0]+CountsDR[2])
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, SG, SP, SCALE1, OFFSET, BETA1, TEMP):
#def FitEIT_MM(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#BETA = 1.8
# SG = 0.6
# SP = 8.1
# TEMP = 0.2e-3
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
#return ScaledFluo1
if not 'popt_6' in globals().keys():
popt_6, pcov_6 = curve_fit(FitEIT_MM_single, FreqsDR, CountsDR, p0=[0.9, 6.2, 5e4, 1e3, 2, 1e-3], bounds=((0, 0, 0, 0, 0, 0), (2, 20, 5e4, 5e4, 10, 10e-3)))
FittedEITpi_6 = FitEIT_MM_single(freqslong, *popt_6)
beta6 = popt_6[4]
errorbeta6 = np.sqrt(pcov_6[4,4])
temp6 = popt_6[5]
errortemp6 = np.sqrt(pcov_6[5,5])
plt.figure()
plt.errorbar(FreqsDR, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
plt.plot(freqslong, FittedEITpi_6, color='darkolivegreen', linewidth=3, label='med 6')
#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()
#%%
#from EITfit.MM_eightLevel_2repumps_AnalysisFunctions import PerformExperiment_8levels
from scipy.optimize import curve_fit
import time
"""
MEDICION 7
"""
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
offsetxpi = 419+correccion-3.7
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
selectedcurve = 7
FreqsDR = [2*f*1e-6-offsetxpi for f in Freqs[0]]
CountsDR = CountsSplit[0][selectedcurve]
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, SG, SP, SCALE1, OFFSET, BETA1, TEMP):
#def FitEIT_MM(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#BETA = 1.8
# SG = 0.6
# SP = 8.1
# TEMP = 0.2e-3
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
#return ScaledFluo1
if not 'popt_7' in globals().keys():
popt_7, pcov_7 = curve_fit(FitEIT_MM_single, FreqsDR, CountsDR, p0=[0.9, 6.2, 3e4, 1.34e3, 2, 1e-3], bounds=((0, 0, 0, 0, 0, 0), (2, 20, 5e4, 5e4, 10, 10e-3)))
FittedEITpi_7 = FitEIT_MM_single(freqslong, *popt_7)
beta7 = popt_7[4]
errorbeta7 = np.sqrt(pcov_7[4,4])
temp7 = popt_7[5]
errortemp7 = np.sqrt(pcov_7[5,5])
plt.figure()
plt.errorbar(FreqsDR, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
plt.plot(freqslong, FittedEITpi_7, color='darkolivegreen', linewidth=3, label='med 7')
#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()
#%%
#from EITfit.MM_eightLevel_2repumps_AnalysisFunctions import PerformExperiment_8levels
from scipy.optimize import curve_fit
import time
"""
MEDICION 8
"""
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
offsetxpi = 419+correccion-4.9
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
selectedcurve = 8
FreqsDR = [2*f*1e-6-offsetxpi for f in Freqs[0]]
CountsDR = CountsSplit[0][selectedcurve]
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, SG, SP, SCALE1, OFFSET, BETA1, TEMP):
#def FitEIT_MM(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#BETA = 1.8
# SG = 0.6
# SP = 8.1
# TEMP = 0.2e-3
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
#return ScaledFluo1
if not 'popt_8' in globals().keys():
popt_8, pcov_8 = curve_fit(FitEIT_MM_single, FreqsDR, CountsDR, p0=[0.9, 6.2, 3e4, 1.34e3, 2, 1e-3], bounds=((0, 0, 0, 0, 0, 0), (2, 20, 5e4, 5e4, 10, 10e-3)))
FittedEITpi_8 = FitEIT_MM_single(freqslong, *popt_8)
beta8 = popt_8[4]
errorbeta8 = np.sqrt(pcov_8[4,4])
temp8 = popt_8[5]
errortemp8 = np.sqrt(pcov_8[5,5])
print()
plt.figure()
plt.errorbar(FreqsDR, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
plt.plot(freqslong, FittedEITpi_8, color='darkolivegreen', linewidth=3, label='med 8')
#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()
#%%
#from EITfit.MM_eightLevel_2repumps_AnalysisFunctions import PerformExperiment_8levels
from scipy.optimize import curve_fit
import time
"""
MEDICION 9
"""
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 = 16
offsetxpi = 419+correccion-6
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
selectedcurve = 9
FreqsDR = [2*f*1e-6-offsetxpi for f in Freqs[0]]
CountsDR = CountsSplit[0][selectedcurve]
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, SG, SP, SCALE1, OFFSET, BETA1, TEMP):
#def FitEIT_MM(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#BETA = 1.8
# SG = 0.6
# SP = 8.1
#TEMP = 0.2e-3
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
#return ScaledFluo1
if not 'popt_9' in globals().keys():
popt_9, pcov_9 = curve_fit(FitEIT_MM_single, FreqsDR, CountsDR, p0=[0.9, 6.2, 3e4, 1.34e3, 2, 1e-3], bounds=((0, 0, 0, 0, 0, 0), (2, 20, 5e4, 5e4, 10,10e-3)))
FittedEITpi_9 = FitEIT_MM_single(freqslong, *popt_9)
beta9 = popt_9[4]
errorbeta9 = np.sqrt(pcov_9[4,4])
temp9 = popt_9[5]
errortemp9 = np.sqrt(pcov_9[5,5])
plt.figure()
plt.errorbar(FreqsDR, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
plt.plot(freqslong, FittedEITpi_9, color='darkolivegreen', linewidth=3, label='med 9')
#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()
#%%
"""
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, 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
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(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)
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)
#%%
"""
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
"""
popt_exp, pcov_exp = curve_fit(expo,Betas_vec,[t*1e3 for t in Temp_vec])
betaslong = np.arange(0,2.7,0.01)
plt.figure()
plt.errorbar(Betas_vec,[t*1e3 for t in Temp_vec],xerr=ErrorBetas_vec, yerr=[t*1e3 for t in ErrorTemp_vec],fmt='o',capsize=5,markersize=5,color=paleta[3])
plt.plot(betaslong,expo(betaslong,*popt_exp))
#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)
#%%
"""
AHORA VAMOS A MEDICIONES CON MAS DE UN ION!!!
"""
"""
Ploteo la cpt de referencia / plotting the reference CPT
1: 2 iones, -100 mV dcA
2: 2 iones, -150 mV dcA
3: 2 iones, -50 mV dcA
4: 2 iones, 5 voltajes (el ion se va en la 4ta medicion y en la 5ta ni esta)
5, 6 y 7: 3 iones en donde el scaneo esta centrado en distintos puntos
"""
jvec = [3] # desde la 1, pero la 4 no porque es un merge de curvitas
plt.figure()
i = 0
for j in jvec:
plt.errorbar([2*f*1e-6 for f in Freqs[j]], Counts[j], yerr=np.sqrt(Counts[j]), fmt='o', capsize=2, markersize=2)
i = i + 1
plt.xlabel('Frecuencia (MHz)')
plt.ylabel('counts')
plt.grid()
#for dr in drs:
# plt.axvline(dr)
#plt.axvline(dr+drive)
plt.legend()
#%%
"""
Mergeo la 5, 6 y 7
"""
Freqs5 = [2*f*1e-6 for f in Freqs[5]]
Freqs6 = [2*f*1e-6 for f in Freqs[6]]
Freqs7 = [2*f*1e-6 for f in Freqs[7]]
Counts5 = Counts[5]
Counts6 = Counts[6]
Counts7 = Counts[7]
i_1_ini = 0
i_1 = 36
i_2_ini = 0
i_2 = 24
f_1 = 18
f_2 = 30
scale_1 = 0.92
scale_2 = 0.98
#Merged_freqs_test = [f-f_2 for f in Freqs6[i_2_ini:i_2]]+[f-f_1 for f in Freqs5[i_1_ini:i_1]]+Freqs7
#plt.plot(Merged_freqs_test,'o')
Merged_freqs = [f-f_2 for f in Freqs6[0:i_2]]+[f-f_1 for f in Freqs5[0:i_1]]+Freqs7
Merged_counts = [scale_2*c for c in Counts6[0:i_2]]+[scale_1*c for c in Counts5[0:i_1]]+list(Counts7)
Merged_freqs_rescaled = np.linspace(np.min(Merged_freqs),np.max(Merged_freqs),len(Merged_freqs))
#drs = [391.5, 399.5, 405.5, 414]
drs = [370,379,385,391.5]
plt.figure()
i = 0
for j in jvec:
plt.plot([f-f_1 for f in Freqs5[0:i_1]], [scale_1*c for c in Counts5[0:i_1]],'o')
plt.plot([f-f_2 for f in Freqs6[0:i_2]], [scale_2*c for c in Counts6[0:i_2]],'o')
plt.plot(Freqs7, Counts7,'o')
plt.errorbar(Merged_freqs, Merged_counts, yerr=np.sqrt(Merged_counts), fmt='o', capsize=2, markersize=2)
i = i + 1
plt.xlabel('Frecuencia (MHz)')
plt.ylabel('counts')
plt.grid()
for dr in drs:
plt.axvline(dr)
plt.axvline(dr+drive, color='red', linestyle='dashed', alpha=0.3)
plt.axvline(dr-drive, color='red', linestyle='dashed', alpha=0.3)
plt.legend()
#%%
"""
ajusto la mergeada de 3 iones
"""
raise ValueError('STOP')
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 = -20
offsetxpi = 438+correccion
DetDoppler = -35-correccion-22
gPS, gPD, = 2*np.pi*21.58e6, 2*np.pi*1.35e6
alpha = 0
drivefreq = 2*np.pi*22.135*1e6
FreqsDR = [f-offsetxpi for f in Merged_freqs]
CountsDR = Merged_counts
freqslong = np.arange(min(FreqsDR), max(FreqsDR)+FreqsDR[1]-FreqsDR[0], 0.1*(FreqsDR[1]-FreqsDR[0]))
CircPr = 1
alpha = 0
import numba
@numba.jit
def FitEIT_MM(freqs, SG, SP, SCALE1, SCALE2, SCALE3, OFFSET, BETA1, BETA2, BETA3):
#def FitEIT_MM(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#BETA = 1.8
# SG = 0.6
# SP = 8.1
TEMP = 0.1e-3
#BETA1, BETA2, BETA3 = 0, 0, 2
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)
Detunings, Fluorescence2 = PerformExperiment_8levels_MM(SG, SP, gPS, gPD, DetDoppler, u, DopplerLaserLinewidth, ProbeLaserLinewidth, TEMP, alpha, phidoppler, titadoppler, phiprobe, titaprobe, BETA2, drivefreq, min(freqs), max(freqs)+(freqs[1]-freqs[0]), freqs[1]-freqs[0], circularityprobe=CircPr, plot=False, solvemode=1, detpvec=None)
Detunings, Fluorescence3 = PerformExperiment_8levels_MM(SG, SP, gPS, gPD, DetDoppler, u, DopplerLaserLinewidth, ProbeLaserLinewidth, TEMP, alpha, phidoppler, titadoppler, phiprobe, titaprobe, BETA3, 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 for f in Fluorescence1])
ScaledFluo2 = np.array([f*SCALE2 for f in Fluorescence2])
ScaledFluo3 = np.array([f*SCALE3 for f in Fluorescence3])
return ScaledFluo1+ScaledFluo2+ScaledFluo3+OFFSET
#return ScaledFluo1
if not 'popt_3ions' in globals().keys():
popt_3ions, pcov_3ions = curve_fit(FitEIT_MM, FreqsDR, CountsDR, p0=[0.6, 6.2, 3.5e5, 3.5e5, 3.5e5, 2e3, 1, 1, 1], bounds=((0, 0, 0, 0, 0, 0, 0, 0, 0), (2, 20, 5e8, 5e8, 5e8, 7e3, 10, 10, 10)))
#popt, pcov = curve_fit(FitEIT_MM, FreqsDR, CountsDR, p0=[0.8, 8, 4e4, 3.5e3, 0], bounds=((0, 0, 0, 0, 0), (2, 15, 1e5, 1e5, 10)))
#array([7.12876797e-01, 7.92474752e+00, 4.29735308e+04, 1.74240582e+04,
#1.53401696e+03, 1.17073206e-06, 2.53804151e+00])
FittedEITpi_3ions = FitEIT_MM(freqslong, *popt_3ions)
#FittedEITpi_3ions = FitEIT_MM(freqslong, popt_3ions[0],popt_3ions[1],popt_3ions[2],popt_3ions[3],popt_3ions[4],popt_3ions[5],4,2,0)
#FittedEITpi_3ions = FitEIT_MM(freqslong, *popt_3ions)
print(popt_3ions)
plt.figure()
plt.errorbar(FreqsDR, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
plt.plot(freqslong, FittedEITpi_3ions, color='darkgreen', linewidth=3)
#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.title(f'Corr:{correccion},DetD:{DetDoppler}')
plt.grid()
#%%
"""
Veo la medicion de varios voltajes uno atras de otro
Se va en medio de la medicion 4, y en la 5 ni esta
"""
jvec = [2] # desde la 1, pero la 4 no porque es un merge de curvitas
Freqs
plt.figure()
i = 0
for j in jvec:
plt.errorbar([2*f*1e-6 for f in Freqs[4]], CountsSplit_2ions[0][j], yerr=np.sqrt(CountsSplit_2ions[0][j]), fmt='o', capsize=2, markersize=2)
i = i + 1
plt.xlabel('Frecuencia (MHz)')
plt.ylabel('counts')
plt.grid()
#for dr in drs:
# plt.axvline(dr)
#plt.axvline(dr+drive)
plt.legend()
#%%
#from EITfit.MM_eightLevel_2repumps_AnalysisFunctions import PerformExperiment_8levels
from scipy.optimize import curve_fit
import time
"""
AJUSTO LA CPT DE 2 IONES CON UN MODELO EN DONDE SUMO DOS ESPECTROS CON BETAS DISTINTOS
"""
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 = 27
offsetxpi = 421+correccion
DetDoppler = -16-correccion+5
gPS, gPD, = 2*np.pi*21.58e6, 2*np.pi*1.35e6
alpha = 0
drivefreq = 2*np.pi*22.135*1e6
FreqsDR = [2*f*1e-6-offsetxpi for f in Freqs[1]]
CountsDR = Counts[1]
freqslong = np.arange(min(FreqsDR), max(FreqsDR)+FreqsDR[1]-FreqsDR[0], 0.1*(FreqsDR[1]-FreqsDR[0]))
CircPr = 1
alpha = 0
def FitEIT_MM(freqs, SG, SP, SCALE1, SCALE2, OFFSET):
#def FitEIT_MM(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#BETA = 1.8
# SG = 0.6
# SP = 8.1
TEMP = 0.1e-3
BETA1, BETA2 = 3, 0
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)
Detunings, Fluorescence2 = PerformExperiment_8levels_MM(SG, SP, gPS, gPD, DetDoppler, u, DopplerLaserLinewidth, ProbeLaserLinewidth, TEMP, alpha, phidoppler, titadoppler, phiprobe, titaprobe, BETA2, 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])
ScaledFluo2 = np.array([f*SCALE2 + OFFSET for f in Fluorescence2])
return ScaledFluo1+ScaledFluo2
#return ScaledFluo1
if not 'popt_2ions_1' in globals().keys():
popt_2ions_1, pcov_2ions_1 = curve_fit(FitEIT_MM, FreqsDR, CountsDR, p0=[0.9, 6.2, 3.5e3, 2.9e3, 3e3], bounds=((0, 0, 0, 0, 0), (2, 20, 5e8, 5e8, 8e3)))
#popt, pcov = curve_fit(FitEIT_MM, FreqsDR, CountsDR, p0=[0.8, 8, 4e4, 3.5e3, 0], bounds=((0, 0, 0, 0, 0), (2, 15, 1e5, 1e5, 10)))
#array([7.12876797e-01, 7.92474752e+00, 4.29735308e+04, 1.74240582e+04,
#1.53401696e+03, 1.17073206e-06, 2.53804151e+00])
FittedEITpi_2sp = FitEIT_MM(freqslong, *popt_2ions_1)
#FittedEITpi = FitEIT_MM(freqslong, 0.8, 8, 4e4, 3.5e3, 0)
# beta1_2ions = popt_2ions_1[5]
# beta2_2ions = popt_2ions_1[6]
# errbeta1_2ions = np.sqrt(pcov_2ions_1[5,5])
# errbeta2_2ions = np.sqrt(pcov_2ions_1[6,6])
"""
Estos params dan bien poniendo beta2=0 y correccion=0 y son SG, SP, SCALE1, SCALE2, OFFSET, BETA1
#array([9.03123248e-01, 6.25865542e+00, 3.47684055e+04, 2.92076804e+04, 1.34556420e+03, 3.55045904e+00])
"""
"""
Ahora considerando ambos betas, con los parametros iniciales dados por los que se obtuvieron con beta2=0
y correccion=0 dan estos parametros que son los de antes pero con BETA2 incluido:
array([8.52685426e-01, 7.42939084e+00, 3.61998310e+04, 3.40160472e+04, 8.62651715e+02, 3.89756335e+00, 7.64867601e-01])
"""
#arreglito = np.array([8.52685426e-01, 7.42939084e+00, 3.61998310e+04, 3.40160472e+04, 8.62651715e+02, 3.89756335e+00, 7.64867601e-01])
FittedEITpi_2ions_1 = FitEIT_MM(freqslong, *popt_2ions_1)
print(popt_2ions_1)
plt.figure()
plt.errorbar(FreqsDR, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
plt.plot(freqslong, FittedEITpi_2ions_1, color='darkgreen', linewidth=3)
#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.title(f'Corr:{correccion},DetD:{DetDoppler}')
plt.grid()
#%%
"""
SUPER AJUSTE PARA MED DE 2 IONES
"""
#from EITfit.MM_eightLevel_2repumps_AnalysisFunctions import PerformExperiment_8levels
from scipy.optimize import curve_fit
import time
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 = [3]
if not 'popt_SA_vec_2ions' in globals().keys():
popt_SA_vec_2ions = []
pcov_SA_vec_2ions = []
for selectedcurve in SelectedCurveVec:
FreqsDR = Freqs[selectedcurve]
CountsDR = Counts[selectedcurve]
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, SCALE2, OFFSET, BETA1, BETA2, TEMP, 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)
Detunings, Fluorescence2 = PerformExperiment_8levels_MM(SG, SP, gPS, gPD, DetDoppler, u, DopplerLaserLinewidth, ProbeLaserLinewidth, TEMP, alpha, phidoppler, titadoppler, phiprobe, titaprobe, BETA2, 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])
ScaledFluo2 = np.array([f*SCALE2 for f in Fluorescence2])
if plot:
return ScaledFluo1+ScaledFluo2, Detunings
else:
return ScaledFluo1+ScaledFluo2
#return ScaledFluo1
if True:
popt_3_SA_2ions, pcov_3_SA_2ions = curve_fit(FitEIT_MM_single, FreqsDR, CountsDR, p0=[448, -42, 0.6, 8.1, 4e4, 4e4, 6e3, 1, 1.2, 0.5e-3], bounds=((0, -100,0, 0, 0,0,0,0,0, 0), (1000, 0, 2, 20,5e6, 5e6,5e4, 10, 10,10e-3)))
#popt_3_SA_2ions = [448, -42, 8e4, 6e3, 2, 0.5e-3]
popt_SA_vec_2ions.append(popt_3_SA_2ions)
pcov_SA_vec_2ions.append(pcov_3_SA_2ions)
FittedEITpi_3_SA_short, Detunings_3_SA_short = FitEIT_MM_single(FreqsDR, *popt_3_SA_2ions, 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_2ions, plot=True)
raise ValueError('Acá tenes que levantar de nuevo los valores que van')
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_2ions)
#print(f'Detdop:{popt_3_SA[1]},popt_3_SA:{popt[0]}')
#%%
#from EITfit.MM_eightLevel_2repumps_AnalysisFunctions import PerformExperiment_8levels
from scipy.optimize import curve_fit
import time
"""
AJUSTO LA CPT DE 2 IONES CON UN MODELO EN DONDE SUMO DOS ESPECTROS CON BETAS DISTINTOS
"""
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 = 27
offsetxpi = 421+correccion
DetDoppler = -16-correccion+5
gPS, gPD, = 2*np.pi*21.58e6, 2*np.pi*1.35e6
alpha = 0
drivefreq = 2*np.pi*22.135*1e6
FreqsDR = [2*f*1e-6-offsetxpi for f in Freqs[1]]
CountsDR = Counts[1]
freqslong = np.arange(min(FreqsDR), max(FreqsDR)+FreqsDR[1]-FreqsDR[0], 0.1*(FreqsDR[1]-FreqsDR[0]))
CircPr = 1
alpha = 0
def FitEIT_MM(freqs, SG, SP, SCALE1, SCALE2, OFFSET):
#def FitEIT_MM(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#BETA = 1.8
# SG = 0.6
# SP = 8.1
TEMP = 0.1e-3
BETA1, BETA2 = 3, 0
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)
Detunings, Fluorescence2 = PerformExperiment_8levels_MM(SG, SP, gPS, gPD, DetDoppler, u, DopplerLaserLinewidth, ProbeLaserLinewidth, TEMP, alpha, phidoppler, titadoppler, phiprobe, titaprobe, BETA2, 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])
ScaledFluo2 = np.array([f*SCALE2 + OFFSET for f in Fluorescence2])
return ScaledFluo1+ScaledFluo2
#return ScaledFluo1
if not 'popt_2ions_1' in globals().keys():
popt_2ions_1, pcov_2ions_1 = curve_fit(FitEIT_MM, FreqsDR, CountsDR, p0=[0.9, 6.2, 3.5e3, 2.9e3, 3e3], bounds=((0, 0, 0, 0, 0), (2, 20, 5e8, 5e8, 8e3)))
#popt, pcov = curve_fit(FitEIT_MM, FreqsDR, CountsDR, p0=[0.8, 8, 4e4, 3.5e3, 0], bounds=((0, 0, 0, 0, 0), (2, 15, 1e5, 1e5, 10)))
#array([7.12876797e-01, 7.92474752e+00, 4.29735308e+04, 1.74240582e+04,
#1.53401696e+03, 1.17073206e-06, 2.53804151e+00])
FittedEITpi_2sp = FitEIT_MM(freqslong, *popt_2ions_1)
#FittedEITpi = FitEIT_MM(freqslong, 0.8, 8, 4e4, 3.5e3, 0)
# beta1_2ions = popt_2ions_1[5]
# beta2_2ions = popt_2ions_1[6]
# errbeta1_2ions = np.sqrt(pcov_2ions_1[5,5])
# errbeta2_2ions = np.sqrt(pcov_2ions_1[6,6])
"""
Estos params dan bien poniendo beta2=0 y correccion=0 y son SG, SP, SCALE1, SCALE2, OFFSET, BETA1
#array([9.03123248e-01, 6.25865542e+00, 3.47684055e+04, 2.92076804e+04, 1.34556420e+03, 3.55045904e+00])
"""
"""
Ahora considerando ambos betas, con los parametros iniciales dados por los que se obtuvieron con beta2=0
y correccion=0 dan estos parametros que son los de antes pero con BETA2 incluido:
array([8.52685426e-01, 7.42939084e+00, 3.61998310e+04, 3.40160472e+04, 8.62651715e+02, 3.89756335e+00, 7.64867601e-01])
"""
#arreglito = np.array([8.52685426e-01, 7.42939084e+00, 3.61998310e+04, 3.40160472e+04, 8.62651715e+02, 3.89756335e+00, 7.64867601e-01])
FittedEITpi_2ions_1 = FitEIT_MM(freqslong, *popt_2ions_1)
print(popt_2ions_1)
plt.figure()
plt.errorbar(FreqsDR, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
plt.plot(freqslong, FittedEITpi_2ions_1, color='darkgreen', linewidth=3)
#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.title(f'Corr:{correccion},DetD:{DetDoppler}')
plt.grid()
#%%
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/20250516_CPT_bladetrap/lolo_CPT_plotter_20231123_multi.py deleted 100644 → 0
View file @ f9b573e5
#!/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
#%% Funciones auxiliares
from scipy.stats.distributions import t,chi2
def estadistica(datos_x,datos_y,modelo,pcov,parametros,nombres=None,alpha=0.05):
if nombres is None:
nombres = [ f'{j}' for j in range(len(parametros)) ]
# Cantidad de parámetros
P = len(parametros)
# Número de datos
N = len(datos_x)
# Grados de libertas (Degrees Of Freedom)
dof = N-P-1
# Cauculamos coordenadas del modelo
# modelo_x = datos_x if modelo_x_arr is None else modelo_x_arr
# modelo_y = modelo( modelo_x, *parametros )
# Predicción del modelo para los datos_x medidos
prediccion_modelo = modelo( datos_x, *parametros )
# Calculos de cantidades estadísticas relevantes
COV = pcov # Matriz de Covarianza
SE = np.sqrt(np.diag( COV )) # Standar Error / Error estandar de los parámetros
residuos = datos_y - prediccion_modelo # diferencia enrte el modelo y los datos
SSE = sum(( residuos )**2 ) # Resitual Sum of Squares
SST = sum(( datos_y - np.mean(datos_y))**2) # Total Sum of Squares
# http://en.wikipedia.org/wiki/Coefficient_of_determination
# Expresa el porcentaje de la varianza que logra explicar el modelos propuesto
Rsq = 1 - SSE/SST # Coeficiente de determinación
Rsq_adj = 1-(1-Rsq) * (N-1)/(N-P-1) # Coeficiente de determinación Ajustado
# https://en.wikipedia.org/wiki/Pearson_correlation_coefficient#In_least_squares_regression_analysis
# Expresa la correlación que hay entre los datos y la predicción del modelo
r_pearson = np.corrcoef( datos_y , prediccion_modelo )[0,1]
# Reduced chi squared
# https://en.wikipedia.org/wiki/Reduced_chi-squared_statistic
chi2_ = sum( residuos**2 )/N
chi2_red = sum( residuos**2 )/(N-P)
# Chi squared test
chi2_test = sum( residuos**2 / abs(prediccion_modelo) )
# p-value del ajuste
p_val = chi2(dof).cdf( chi2_test )
sT = t.ppf(1.0 - alpha/2.0, N - P ) # student T multiplier
CI = sT * SE # Confidence Interval
print('R-squared ',Rsq)
print('R-sq_adjusted',Rsq_adj)
print('chi2 ',chi2_)
print('chi2_reduced ',chi2_red)
print('chi2_test ',chi2_test)
print('r-pearson ',r_pearson)
print('p-value ',p_val)
print('')
print('Error Estandard (SE):')
for i in range(P):
print(f'parametro[{nombres[i]:>5s}]: ' , parametros[i], ' ± ' , SE[i])
print('')
print('Intervalo de confianza al '+str((1-alpha)*100)+'%:')
for i in range(P):
print(f'parametro[{nombres[i]:>5s}]: ' , parametros[i], ' ± ' , CI[i])
return dict(R2=Rsq,R2_adj=Rsq_adj,chi2=chi2_,chi2_red=chi2_red,
chi2_test=chi2_test,r=r_pearson,pvalue=p_val,
SE=SE,CI=CI)
#%% 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])))
#%%
"""
Ploteo la cpt de referencia / plotting the reference CPT
"""
jvec = [2] # de la 1 a la 9 vale la pena, despues no
drs = [390.5, 399.5, 406, 413.5]
drive=22.1
Frequencies = Freqs[0]
plt.figure()
i = 0
for j in jvec:
plt.errorbar([2*f*1e-6 for f in Frequencies], CountsSplit[0][j], yerr=np.sqrt(CountsSplit[0][j]), fmt='o', capsize=2, markersize=2)
i = i + 1
plt.xlabel('Frecuencia (MHz)')
plt.ylabel('counts')
plt.grid()
#for dr in drs:
# plt.axvline(dr)
#plt.axvline(dr+drive)
plt.legend()
#%%
#%%
#from EITfit.MM_eightLevel_2repumps_AnalysisFunctions import PerformExperiment_8levels
from scipy.optimize import curve_fit
import time
#%%
"""
AHORA VAMOS A MEDICIONES CON MAS DE UN ION!!!
Las mediciones estan buenas, habria que ver de ajustarlas bien, yo no lo logre.
"""
"""
Ploteo la cpt de referencia / plotting the reference CPT
1: 2 iones, -100 mV dcA
2: 2 iones, -150 mV dcA
3: 2 iones, -50 mV dcA
4: 2 iones, 5 voltajes (el ion se va en la 4ta medicion y en la 5ta ni esta)
5, 6 y 7: 3 iones en donde el scaneo esta centrado en distintos puntos
"""
jvec = [3] # desde la 1, pero la 4 no porque es un merge de curvitas
plt.figure()
i = 0
for j in jvec:
plt.errorbar([2*f*1e-6 for f in Freqs[j]], Counts[j], yerr=np.sqrt(Counts[j]), fmt='o', capsize=2, markersize=2)
i = i + 1
plt.xlabel('Frecuencia (MHz)')
plt.ylabel('counts')
plt.grid()
#for dr in drs:
# plt.axvline(dr)
#plt.axvline(dr+drive)
plt.legend()
#%%
"""
Mergeo la 5, 6 y 7
"""
Freqs5 = [2*f*1e-6 for f in Freqs[5]]
Freqs6 = [2*f*1e-6 for f in Freqs[6]]
Freqs7 = [2*f*1e-6 for f in Freqs[7]]
Counts5 = Counts[5]
Counts6 = Counts[6]
Counts7 = Counts[7]
i_1_ini = 0
i_1 = 36
i_2_ini = 0
i_2 = 24
f_1 = 18
f_2 = 30
scale_1 = 0.92
scale_2 = 0.98
#Merged_freqs_test = [f-f_2 for f in Freqs6[i_2_ini:i_2]]+[f-f_1 for f in Freqs5[i_1_ini:i_1]]+Freqs7
#plt.plot(Merged_freqs_test,'o')
Merged_freqs = [f-f_2 for f in Freqs6[0:i_2]]+[f-f_1 for f in Freqs5[0:i_1]]+Freqs7
Merged_counts = [scale_2*c for c in Counts6[0:i_2]]+[scale_1*c for c in Counts5[0:i_1]]+list(Counts7)
Merged_freqs_rescaled = np.linspace(np.min(Merged_freqs),np.max(Merged_freqs),len(Merged_freqs))
#drs = [391.5, 399.5, 405.5, 414]
drs = [370,379,385,391.5]
plt.figure()
i = 0
for j in jvec:
plt.plot([f-f_1 for f in Freqs5[0:i_1]], [scale_1*c for c in Counts5[0:i_1]],'o')
plt.plot([f-f_2 for f in Freqs6[0:i_2]], [scale_2*c for c in Counts6[0:i_2]],'o')
plt.plot(Freqs7, Counts7,'o')
plt.errorbar(Merged_freqs, Merged_counts, yerr=np.sqrt(Merged_counts), fmt='o', capsize=2, markersize=2)
i = i + 1
plt.xlabel('Frecuencia (MHz)')
plt.ylabel('counts')
plt.grid()
for dr in drs:
plt.axvline(dr)
plt.axvline(dr+drive, color='red', linestyle='dashed', alpha=0.3)
plt.axvline(dr-drive, color='red', linestyle='dashed', alpha=0.3)
plt.legend()
#%%
"""
ajusto la mergeada de 3 iones
"""
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 = -20
offsetxpi = 438+correccion
DetDoppler = -35-correccion-22
gPS, gPD, = 2*np.pi*21.58e6, 2*np.pi*1.35e6
alpha = 0
drivefreq = 2*np.pi*22.135*1e6
FreqsDR = [f-offsetxpi for f in Merged_freqs]
CountsDR = Merged_counts
freqslong = np.arange(min(FreqsDR), max(FreqsDR)+FreqsDR[1]-FreqsDR[0], 0.1*(FreqsDR[1]-FreqsDR[0]))
CircPr = 1
alpha = 0
import numba
@numba.jit
def FitEIT_MM1(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#def FitEIT_MM(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#BETA = 1.8
# SG = 0.6
# SP = 8.1
TEMP = 0.1e-3
#BETA1, BETA2, BETA3 = 0, 0, 2
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 for f in Fluorescence1])
return ScaledFluo1+OFFSET
#return ScaledFluo1
@numba.jit
def FitEIT_MM(freqs, SG, SP, SCALE1, SCALE2, SCALE3, OFFSET, BETA1, BETA2, BETA3):
#def FitEIT_MM(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#BETA = 1.8
# SG = 0.6
# SP = 8.1
TEMP = 0.1e-3
freqs = np.array(freqs)
#BETA1, BETA2, BETA3 = 0, 0, 2
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)
Detunings, Fluorescence2 = PerformExperiment_8levels_MM(SG, SP, gPS, gPD, DetDoppler, u, DopplerLaserLinewidth, ProbeLaserLinewidth, TEMP, alpha, phidoppler, titadoppler, phiprobe, titaprobe, BETA2, drivefreq, min(freqs), max(freqs)+(freqs[1]-freqs[0]), freqs[1]-freqs[0], circularityprobe=CircPr, plot=False, solvemode=1, detpvec=None)
Detunings, Fluorescence3 = PerformExperiment_8levels_MM(SG, SP, gPS, gPD, DetDoppler, u, DopplerLaserLinewidth, ProbeLaserLinewidth, TEMP, alpha, phidoppler, titadoppler, phiprobe, titaprobe, BETA3, 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 for f in Fluorescence1])
# ScaledFluo2 = np.array([f*SCALE2 for f in Fluorescence2])
# ScaledFluo3 = np.array([f*SCALE3 for f in Fluorescence3])
# return ScaledFluo1+ScaledFluo2+ScaledFluo3+OFFSET
Fluorescence1 = np.array(Fluorescence1)
Fluorescence2 = np.array(Fluorescence2)
Fluorescence3 = np.array(Fluorescence3)
return SCALE1*Fluorescence1+SCALE2*Fluorescence2+SCALE3*Fluorescence3+OFFSET
if not 'popt_1ions' in globals().keys():
t0 = time.time()
print("Arranamos FIT 1")
par_ini = [0.65, 7.06, 86070, 3917, 1.64]
popt_1ions, pcov_1ions = curve_fit(FitEIT_MM1, FreqsDR, CountsDR, p0=par_ini, bounds=((0, 0, 0, 0, 0), (2, 20, 5e8, 7e3, 10)))
pp1 = estadistica(FreqsDR,CountsDR,FitEIT_MM1,pcov_1ions,popt_1ions,nombres=None,alpha=0.05)
print(f"time: {round(time.time()-t0,1)} seg")
if not 'popt_3ions' in globals().keys():
t0 = time.time()
print("Arranamos FIT 1")
par_ini = [0.65, 7.06, 86070, 3917, 1.64]
popt_1ions, pcov_1ions = curve_fit(FitEIT_MM1, FreqsDR, CountsDR, p0=par_ini, bounds=((0, 0, 0, 0, 0), (2, 20, 5e8, 7e3, 10)))
print(f"time: {round(time.time()-t0,1)} seg")
print("Arranamos FIT 3")
def fun(freqs, SG, SP, SCALE1, SCALE2, OFFSET, BETA1, BETA2, BETA3):
SCALE3 = max(1-SCALE2-SCALE1-OFFSET,0)
Fluorescence3 = FitEIT_MM(freqs, SG, SP, SCALE1, SCALE2, SCALE3, OFFSET, BETA1, BETA2, BETA3)
return Fluorescence3/sum(Fluorescence3)
par_ini = popt_1ions.tolist()[:2]+[1/3]*2 + [0.1] +popt_1ions.tolist()[-1:]*3
bounds = ((0, 0, 0, 0, 0, 0, 0, 0),
(2, 20, 1, 1, 1, 10, 10, 10))
popt_3ions, pcov_3ions = curve_fit(fun, FreqsDR, CountsDR, p0=par_ini,bounds=bounds )
print(f"time: {round(time.time()-t0,1)} seg")
freqs = np.array(FreqsDR)
SG, SP, SCALE1, SCALE2, OFFSET, BETA1, BETA2, BETA3 = popt_3ions
SCALE3 = max(1-SCALE2-SCALE1-OFFSET,0)
# Fluorescence= FitEIT_MM(freqs, SG, SP, SCALE1, SCALE2, SCALE3, OFFSET, BETA1, BETA2, BETA3)
Fluorescence= fun(freqs, SG, SP, SCALE1, SCALE2, OFFSET, BETA1, BETA2, BETA3)
plt.figure()
plt.errorbar(freqs, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
plt.plot(freqs, Fluorescence*sum(CountsDR))
plt.plot(freqs, fun(freqs, *par_ini )*sum(CountsDR))
#%% Acá hay un ajuste del de 3 iones que da razonable
FreqsDR, CountsDR = np.array(FreqsDR) , np.array(CountsDR)
t0 = time.time()
print("Arranamos FIT 1")
par_ini = [0.65, 7.06, 86070, 3917, 1.64]
popt_1ions, pcov_1ions = curve_fit(FitEIT_MM1, FreqsDR, CountsDR, p0=par_ini, bounds=((0, 0, 0, 0, 0), (2, 20, 5e8, 7e3, 10)))
print(f"time: {round(time.time()-t0,1)} seg")
par_ini = popt_1ions.tolist()[:2]+[1/3]*2 + [0.1] +popt_1ions.tolist()[-1:]*3
bounds = ((0, 0, 0, 0, 0, 0, 0, 0),
(2, 20, 1, 1, 1, 10, 10, 10))
# bounds = (-np.inf, np.inf)
plt.figure()
plt.errorbar(freqs, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
l1, =plt.plot(freqs, fun(freqs, *par_ini )*sum(CountsDR))
plt.draw()
print("Arranamos FIT 3")
def fun(freqs, SG, SP, SCALE1, SCALE2, OFFSET, BETA1, BETA2, BETA3):
SCALE3 = max(1-SCALE2-SCALE1-OFFSET,0)
Fluorescence3 = FitEIT_MM(freqs, SG, SP, SCALE1, SCALE2, SCALE3, OFFSET, BETA1, BETA2, BETA3)
l1.set_ydata(Fluorescence3/sum(Fluorescence3)*sum(CountsDR))
print(SG, SP, SCALE1, SCALE2, OFFSET, BETA1, BETA2, BETA3)
plt.pause(0.1)
return Fluorescence3/sum(Fluorescence3)
popt_3ions, pcov_3ions = curve_fit(fun, FreqsDR, CountsDR/CountsDR.sum(), p0=par_ini,bounds=bounds )
print(f"time: {round(time.time()-t0,1)} seg")
freqs = np.array(FreqsDR)
SG, SP, SCALE1, SCALE2, OFFSET, BETA1, BETA2, BETA3 = popt_3ions
SCALE3 = max(1-SCALE2-SCALE1-OFFSET,0)
# Fluorescence= FitEIT_MM(freqs, SG, SP, SCALE1, SCALE2, SCALE3, OFFSET, BETA1, BETA2, BETA3)
Fluorescence= fun(freqs, SG, SP, SCALE1, SCALE2, OFFSET, BETA1, BETA2, BETA3)
plt.plot(freqs, Fluorescence*sum(CountsDR))
if False:
popt_3ions = np.array([0.70790618, 7.41165226, 0.1707309 , 0.13974759, 0.02322333,
3.69298275, 3.68847693, 1.36216489])
#%%
#popt, pcov = curve_fit(FitEIT_MM, FreqsDR, CountsDR, p0=[0.8, 8, 4e4, 3.5e3, 0], bounds=((0, 0, 0, 0, 0), (2, 15, 1e5, 1e5, 10)))
#array([7.12876797e-01, 7.92474752e+00, 4.29735308e+04, 1.74240582e+04,
#1.53401696e+03, 1.17073206e-06, 2.53804151e+00])
if False:
# Ejemplo 3 iones
t0 = time.time()
freqs, SG, SP = np.array(FreqsDR), 0.65, 7.06
SCALE1, SCALE2, SCALE3, OFFSET = 1,1,1,1
BETA1, BETA2, BETA3 = 1.64, 2, 3
TEMP = 0.1e-3
Fluorescence3= FitEIT_MM(freqs, SG, SP, SCALE1, SCALE2, SCALE3, OFFSET, BETA1, BETA2, BETA3)
print(f"time: {round(time.time()-t0,1)} seg")
Detunings, Fluorescence = np.array(freqs), np.array(Fluorescence3)
plt.figure()
plt.errorbar(FreqsDR, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
plt.plot(Detunings, Fluorescence*sum(CountsDR)/sum(Fluorescence))
if False:
# Ejemplo 1 ion
t0 = time.time()
freqs, SG, SP, BETA1 = FreqsDR, 0.65, 7.06, 1.64
TEMP = 0.1e-3
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)
print(f"time: {round(time.time()-t0,1)} seg")
Detunings, Fluorescence1 = np.array(Detunings), np.array(Fluorescence1)
plt.figure()
plt.errorbar(FreqsDR, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
plt.plot(Detunings, Fluorescence1*sum(CountsDR)/sum(Fluorescence1))
FittedEITpi_3ions = FitEIT_MM(freqslong, *popt_3ions)
#FittedEITpi_3ions = FitEIT_MM(freqslong, popt_3ions[0],popt_3ions[1],popt_3ions[2],popt_3ions[3],popt_3ions[4],popt_3ions[5],4,2,0)
#FittedEITpi_3ions = FitEIT_MM(freqslong, *popt_3ions)
print(popt_3ions)
plt.figure()
plt.errorbar(FreqsDR, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
plt.plot(freqslong, FittedEITpi_3ions, color='darkgreen', linewidth=3)
#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.title(f'Corr:{correccion},DetD:{DetDoppler}')
plt.grid()
#%%
"""
Veo la medicion de varios voltajes uno atras de otro
Se va en medio de la medicion 4, y en la 5 ni esta
"""
jvec = [2] # desde la 1, pero la 4 no porque es un merge de curvitas
Freqs
plt.figure()
i = 0
for j in jvec:
plt.errorbar([2*f*1e-6 for f in Freqs[4]], CountsSplit_2ions[0][j], yerr=np.sqrt(CountsSplit_2ions[0][j]), fmt='o', capsize=2, markersize=2)
i = i + 1
plt.xlabel('Frecuencia (MHz)')
plt.ylabel('counts')
plt.grid()
#for dr in drs:
# plt.axvline(dr)
#plt.axvline(dr+drive)
plt.legend()
#%%
#from EITfit.MM_eightLevel_2repumps_AnalysisFunctions import PerformExperiment_8levels
from scipy.optimize import curve_fit
import time
"""
AJUSTO LA CPT DE 2 IONES CON UN MODELO EN DONDE SUMO DOS ESPECTROS CON BETAS DISTINTOS
"""
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 = 27
offsetxpi = 421+correccion
DetDoppler = -16-correccion+5
gPS, gPD, = 2*np.pi*21.58e6, 2*np.pi*1.35e6
alpha = 0
drivefreq = 2*np.pi*22.135*1e6
FreqsDR = [2*f*1e-6-offsetxpi for f in Freqs[1]]
CountsDR = Counts[1]
freqslong = np.arange(min(FreqsDR), max(FreqsDR)+FreqsDR[1]-FreqsDR[0], 0.1*(FreqsDR[1]-FreqsDR[0]))
CircPr = 1
alpha = 0
def FitEIT_MM(freqs, SG, SP, SCALE1, SCALE2, OFFSET):
#def FitEIT_MM(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#BETA = 1.8
# SG = 0.6
# SP = 8.1
TEMP = 0.1e-3
BETA1, BETA2 = 3, 0
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)
Detunings, Fluorescence2 = PerformExperiment_8levels_MM(SG, SP, gPS, gPD, DetDoppler, u, DopplerLaserLinewidth, ProbeLaserLinewidth, TEMP, alpha, phidoppler, titadoppler, phiprobe, titaprobe, BETA2, 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])
ScaledFluo2 = np.array([f*SCALE2 + OFFSET for f in Fluorescence2])
return ScaledFluo1+ScaledFluo2
#return ScaledFluo1
if not 'popt_2ions_1' in globals().keys():
popt_2ions_1, pcov_2ions_1 = curve_fit(FitEIT_MM, FreqsDR, CountsDR, p0=[0.9, 6.2, 3.5e3, 2.9e3, 3e3], bounds=((0, 0, 0, 0, 0), (2, 20, 5e8, 5e8, 8e3)))
#popt, pcov = curve_fit(FitEIT_MM, FreqsDR, CountsDR, p0=[0.8, 8, 4e4, 3.5e3, 0], bounds=((0, 0, 0, 0, 0), (2, 15, 1e5, 1e5, 10)))
#array([7.12876797e-01, 7.92474752e+00, 4.29735308e+04, 1.74240582e+04,
#1.53401696e+03, 1.17073206e-06, 2.53804151e+00])
FittedEITpi_2sp = FitEIT_MM(freqslong, *popt_2ions_1)
#FittedEITpi = FitEIT_MM(freqslong, 0.8, 8, 4e4, 3.5e3, 0)
# beta1_2ions = popt_2ions_1[5]
# beta2_2ions = popt_2ions_1[6]
# errbeta1_2ions = np.sqrt(pcov_2ions_1[5,5])
# errbeta2_2ions = np.sqrt(pcov_2ions_1[6,6])
"""
Estos params dan bien poniendo beta2=0 y correccion=0 y son SG, SP, SCALE1, SCALE2, OFFSET, BETA1
#array([9.03123248e-01, 6.25865542e+00, 3.47684055e+04, 2.92076804e+04, 1.34556420e+03, 3.55045904e+00])
"""
"""
Ahora considerando ambos betas, con los parametros iniciales dados por los que se obtuvieron con beta2=0
y correccion=0 dan estos parametros que son los de antes pero con BETA2 incluido:
array([8.52685426e-01, 7.42939084e+00, 3.61998310e+04, 3.40160472e+04, 8.62651715e+02, 3.89756335e+00, 7.64867601e-01])
"""
#arreglito = np.array([8.52685426e-01, 7.42939084e+00, 3.61998310e+04, 3.40160472e+04, 8.62651715e+02, 3.89756335e+00, 7.64867601e-01])
FittedEITpi_2ions_1 = FitEIT_MM(freqslong, *popt_2ions_1)
print(popt_2ions_1)
plt.figure()
plt.errorbar(FreqsDR, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
plt.plot(freqslong, FittedEITpi_2ions_1, color='darkgreen', linewidth=3)
#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.title(f'Corr:{correccion},DetD:{DetDoppler}')
plt.grid()
#%%
"""
SUPER AJUSTE PARA MED DE 2 IONES
"""
#from EITfit.MM_eightLevel_2repumps_AnalysisFunctions import PerformExperiment_8levels
from scipy.optimize import curve_fit
import time
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 = [3]
if not 'popt_SA_vec_2ions' in globals().keys():
popt_SA_vec_2ions = []
pcov_SA_vec_2ions = []
for selectedcurve in SelectedCurveVec:
FreqsDR = Freqs[selectedcurve]
CountsDR = Counts[selectedcurve]
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, SCALE2, OFFSET, BETA1, BETA2, TEMP, 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)
Detunings, Fluorescence2 = PerformExperiment_8levels_MM(SG, SP, gPS, gPD, DetDoppler, u, DopplerLaserLinewidth, ProbeLaserLinewidth, TEMP, alpha, phidoppler, titadoppler, phiprobe, titaprobe, BETA2, 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])
ScaledFluo2 = np.array([f*SCALE2 for f in Fluorescence2])
if plot:
return ScaledFluo1+ScaledFluo2, Detunings
else:
return ScaledFluo1+ScaledFluo2
#return ScaledFluo1
if True:
popt_3_SA_2ions, pcov_3_SA_2ions = curve_fit(FitEIT_MM_single, FreqsDR, CountsDR, p0=[448, -42, 0.6, 8.1, 4e4, 4e4, 6e3, 1, 1.2, 0.5e-3], bounds=((0, -100,0, 0, 0,0,0,0,0, 0), (1000, 0, 2, 20,5e6, 5e6,5e4, 10, 10,10e-3)))
#popt_3_SA_2ions = [448, -42, 8e4, 6e3, 2, 0.5e-3]
popt_SA_vec_2ions.append(popt_3_SA_2ions)
pcov_SA_vec_2ions.append(pcov_3_SA_2ions)
FittedEITpi_3_SA_short, Detunings_3_SA_short = FitEIT_MM_single(FreqsDR, *popt_3_SA_2ions, 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_2ions, plot=True)
raise ValueError('Acá tenes que levantar de nuevo los valores que van')
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_2ions)
#print(f'Detdop:{popt_3_SA[1]},popt_3_SA:{popt[0]}')
#%%
#from EITfit.MM_eightLevel_2repumps_AnalysisFunctions import PerformExperiment_8levels
from scipy.optimize import curve_fit
import time
"""
AJUSTO LA CPT DE 2 IONES CON UN MODELO EN DONDE SUMO DOS ESPECTROS CON BETAS DISTINTOS
"""
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 = 27
offsetxpi = 421+correccion
DetDoppler = -16-correccion+5
gPS, gPD, = 2*np.pi*21.58e6, 2*np.pi*1.35e6
alpha = 0
drivefreq = 2*np.pi*22.135*1e6
FreqsDR = [2*f*1e-6-offsetxpi for f in Freqs[1]]
CountsDR = Counts[1]
freqslong = np.arange(min(FreqsDR), max(FreqsDR)+FreqsDR[1]-FreqsDR[0], 0.1*(FreqsDR[1]-FreqsDR[0]))
CircPr = 1
alpha = 0
def FitEIT_MM(freqs, SG, SP, SCALE1, SCALE2, OFFSET):
#def FitEIT_MM(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#BETA = 1.8
# SG = 0.6
# SP = 8.1
TEMP = 0.1e-3
BETA1, BETA2 = 3, 0
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)
Detunings, Fluorescence2 = PerformExperiment_8levels_MM(SG, SP, gPS, gPD, DetDoppler, u, DopplerLaserLinewidth, ProbeLaserLinewidth, TEMP, alpha, phidoppler, titadoppler, phiprobe, titaprobe, BETA2, 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])
ScaledFluo2 = np.array([f*SCALE2 + OFFSET for f in Fluorescence2])
return ScaledFluo1+ScaledFluo2
#return ScaledFluo1
if not 'popt_2ions_1' in globals().keys():
popt_2ions_1, pcov_2ions_1 = curve_fit(FitEIT_MM, FreqsDR, CountsDR, p0=[0.9, 6.2, 3.5e3, 2.9e3, 3e3], bounds=((0, 0, 0, 0, 0), (2, 20, 5e8, 5e8, 8e3)))
#popt, pcov = curve_fit(FitEIT_MM, FreqsDR, CountsDR, p0=[0.8, 8, 4e4, 3.5e3, 0], bounds=((0, 0, 0, 0, 0), (2, 15, 1e5, 1e5, 10)))
#array([7.12876797e-01, 7.92474752e+00, 4.29735308e+04, 1.74240582e+04,
#1.53401696e+03, 1.17073206e-06, 2.53804151e+00])
FittedEITpi_2sp = FitEIT_MM(freqslong, *popt_2ions_1)
#FittedEITpi = FitEIT_MM(freqslong, 0.8, 8, 4e4, 3.5e3, 0)
# beta1_2ions = popt_2ions_1[5]
# beta2_2ions = popt_2ions_1[6]
# errbeta1_2ions = np.sqrt(pcov_2ions_1[5,5])
# errbeta2_2ions = np.sqrt(pcov_2ions_1[6,6])
"""
Estos params dan bien poniendo beta2=0 y correccion=0 y son SG, SP, SCALE1, SCALE2, OFFSET, BETA1
#array([9.03123248e-01, 6.25865542e+00, 3.47684055e+04, 2.92076804e+04, 1.34556420e+03, 3.55045904e+00])
"""
"""
Ahora considerando ambos betas, con los parametros iniciales dados por los que se obtuvieron con beta2=0
y correccion=0 dan estos parametros que son los de antes pero con BETA2 incluido:
array([8.52685426e-01, 7.42939084e+00, 3.61998310e+04, 3.40160472e+04, 8.62651715e+02, 3.89756335e+00, 7.64867601e-01])
"""
#arreglito = np.array([8.52685426e-01, 7.42939084e+00, 3.61998310e+04, 3.40160472e+04, 8.62651715e+02, 3.89756335e+00, 7.64867601e-01])
FittedEITpi_2ions_1 = FitEIT_MM(freqslong, *popt_2ions_1)
print(popt_2ions_1)
plt.figure()
plt.errorbar(FreqsDR, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
plt.plot(freqslong, FittedEITpi_2ions_1, color='darkgreen', linewidth=3)
#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.title(f'Corr:{correccion},DetD:{DetDoppler}')
plt.grid()
#%%
"""
AHORA INTENTO SUPER AJUSTES O SEA CON OFFSETXPI Y DETDOPPLER INCLUIDOS
"""
#%%
"""
SUPER AJUSTE (SA)
"""
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]
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]
SelectedCurveVec = [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)
analisis/plots/20250516_CPT_bladetrap/lolo_ajuste_3iones_A.py deleted 100644 → 0
View file @ f9b573e5
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Acá hacemos el ajuste de 3 iones del archivo
artiq_experiments/analisis/plots/20231123_CPTconmicromocion3/CPT_plotter_20231123.py
sacado de la parte que dice:
"ajusto la mergeada de 3 iones"
@author: lolo
"""
import h5py
import matplotlib.pyplot as plt
import numpy as np
from scipy.optimize import curve_fit
from Data.EITfit.lolo_modelo_full_8niveles import PerformExperiment_8levels_MM
import time
#%% Funciones auxiliares para información estadística
from scipy.stats.distributions import t,chi2
import inspect
def estadistica(datos_x,datos_y,modelo,pcov,parametros,nombres=None,alpha=0.05):
if nombres is None:
nombres = inspect.getfullargspec(modelo).args[1:]
# Cantidad de parámetros
P = len(parametros)
# Número de datos
N = len(datos_x)
# Grados de libertas (Degrees Of Freedom)
dof = N-P-1
# Cauculamos coordenadas del modelo
# modelo_x = datos_x if modelo_x_arr is None else modelo_x_arr
# modelo_y = modelo( modelo_x, *parametros )
# Predicción del modelo para los datos_x medidos
prediccion_modelo = modelo( datos_x, *parametros )
# Calculos de cantidades estadísticas relevantes
COV = pcov # Matriz de Covarianza
SE = np.sqrt(np.diag( COV )) # Standar Error / Error estandar de los parámetros
residuos = datos_y - prediccion_modelo # diferencia enrte el modelo y los datos
SSE = sum(( residuos )**2 ) # Resitual Sum of Squares
SST = sum(( datos_y - np.mean(datos_y))**2) # Total Sum of Squares
# http://en.wikipedia.org/wiki/Coefficient_of_determination
# Expresa el porcentaje de la varianza que logra explicar el modelos propuesto
Rsq = 1 - SSE/SST # Coeficiente de determinación
Rsq_adj = 1-(1-Rsq) * (N-1)/(N-P-1) # Coeficiente de determinación Ajustado
# https://en.wikipedia.org/wiki/Pearson_correlation_coefficient#In_least_squares_regression_analysis
# Expresa la correlación que hay entre los datos y la predicción del modelo
r_pearson = np.corrcoef( datos_y , prediccion_modelo )[0,1]
# Reduced chi squared
# https://en.wikipedia.org/wiki/Reduced_chi-squared_statistic
chi2_ = sum( residuos**2 )/N
chi2_red = sum( residuos**2 )/(N-P)
# Chi squared test
chi2_test = sum( residuos**2 / abs(prediccion_modelo) )
# p-value del ajuste
p_val = chi2(dof).cdf( chi2_test )
sT = t.ppf(1.0 - alpha/2.0, N - P ) # student T multiplier
CI = sT * SE # Confidence Interval
return dict(R2=Rsq,R2_adj=Rsq_adj,chi2=chi2_,chi2_red=chi2_red,
chi2_test=chi2_test,r=r_pearson,pvalue=p_val,
SE=SE,CI=CI,parametros=parametros,nombres=nombres,alpha=alpha,N=N)
def print_estadistica(dic):
pars = 'N,pvalue,R2,R2_adj,chi2,chi2_red,chi2_test,r,SE,CI,parametros,nombres,alpha'.split(',')
N,pvalue,R2,R2_adj,chi2,chi2_red,chi2_test,r,SE,CI,parametros,nombres,alpha = [ dic[k] for k in pars ]
P = len(parametros)
print('R-squared ',R2)
print('R-sq_adjusted',R2_adj)
print('chi2 ',chi2)
print('chi2_reduced ',chi2_red)
print('chi2_test ',chi2_test)
print('r-pearson ',r)
print('p-value ',pvalue)
print('')
print('Error Estandard (SE):')
for i in range(P):
print(f'parametro[{nombres[i]:>5s}]: ' , parametros[i], ' ± ' , SE[i])
print('')
print('Intervalo de confianza al '+str((1-alpha)*100)+'%:')
for i in range(P):
print(f'parametro[{nombres[i]:>5s}]: ' , parametros[i], ' ± ' , CI[i])
def fit(model,x_data, y_data, fixed=None, plot=False, **kwargs):
"""Usa los mismos argumentos que curve_fit pero permite fixear parámetros y
visualizar el ajuste
fixed: nombres de los parametros fijos, con espacios
plot: si plotea o no el ajuste
"""
# Extraemos los parámetros del modelo
param_names = inspect.getfullargspec(model).args[1:]
if not 'p0' in kwargs.keys():
kwargs['p0'] = np.ones(len(param_names))
params_initial = kwargs['p0']
# Nos fijamos si hay que fijar algo
if isinstance(fixed, str):
fixed = fixed.split()
if fixed is None:
fixed = []
# Armoa lista de True/False de los parametros fijos
params_fixed = [ True if p in fixed else False for p in param_names ]
if True in params_fixed and 'bounds' in kwargs.keys():
try:
b0 = [ bb for fi,bb in zip(params_fixed,kwargs['bounds'][0]) if fi==False ]
except:
b0 = kwargs['bounds'][0]
try:
b1 = [ bb for fi,bb in zip(params_fixed,kwargs['bounds'][1]) if fi==False ]
except:
b1 =kwargs['bounds'][0]
kwargs['bounds'] = (b0,b1)
print("Parámetros:", param_names)
# print(locals())
print('\n')
if plot:
plt.figure()
plt.plot(x_data, y_data,'.')
ylim=plt.gca().get_ylim()
l1, = plt.plot(x_data, model(x_data,*params_initial),'r-')
plt.ylim(ylim)
# Función auxiliar para el ajuste
def aux_func(x_data,*args):
all_params = np.ones(len(param_names))
jj=0
for ii in range(len(all_params)):
if params_fixed[ii]:
all_params[ii] = params_initial[ii]
else:
all_params[ii] = args[jj]
jj += 1
# print(params_fixed,all_params)
rta = model(x_data,*all_params)
if plot:
l1.set_ydata(rta)
plt.draw()
plt.pause(0.01)
return rta
# Hacemos el curve_fit
t0 = time.time()
print("Arancamos el fit")
kwargs['p0'] = [ p for p,ff in zip(params_initial,params_fixed) if ff is False ]
params, pcov = curve_fit(aux_func, x_data, y_data, **kwargs )
print(f"time: {round(time.time()-t0,1)} seg")
# Rearmamos los parámetros para incluir los fijos
params2 = np.ones(len(param_names))
jj=0
for ii in range(len(params2)):
if params_fixed[ii]:
params2[ii] = params_initial[ii]
else:
params2[ii] = params[jj]
jj += 1
# Rearmo la pcov
ii, a = 0, []
for fixed in params_fixed:
if fixed:
a.append( [0]*len(pcov[0]) )
else:
a.append( pcov[ii].tolist() )
ii+=1
a = np.array(a)
ii, b = 0, []
for fixed in params_fixed:
if fixed:
b.append( [0]*len(a.T[0]) )
else:
b.append( a.T[ii].tolist() )
ii+=1
pcov2 = np.array(b)
# Caluclo información estadística
info = estadistica(x_data,y_data,model,pcov2,params2,nombres=None,alpha=0.05)
print_estadistica(info)
return params2, pcov2 , info
# fit(modelo, FreqsDR, CountsDR/CountsDR.sum(), p0=(1,1,1),bounds=(-1e6,1e6) )
# fit(modelo, FreqsDR, CountsDR/CountsDR.sum(), fixed="", plot=True, p0=(1,2,3),bounds=(-1e6,1e6) )
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
#%% Levantamos los datos %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
"""
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])))
###############################################################################
"""
AHORA VAMOS A MEDICIONES CON MAS DE UN ION!!!
Las mediciones estan buenas, habria que ver de ajustarlas bien, yo no lo logre.
"""
"""
Ploteo la cpt de referencia / plotting the reference CPT
1: 2 iones, -100 mV dcA
2: 2 iones, -150 mV dcA
3: 2 iones, -50 mV dcA
4: 2 iones, 5 voltajes (el ion se va en la 4ta medicion y en la 5ta ni esta)
5, 6 y 7: 3 iones en donde el scaneo esta centrado en distintos puntos
"""
"""
Mergeo la 5, 6 y 7
"""
Freqs5 = [2*f*1e-6 for f in Freqs[5]]
Freqs6 = [2*f*1e-6 for f in Freqs[6]]
Freqs7 = [2*f*1e-6 for f in Freqs[7]]
Counts5 = Counts[5]
Counts6 = Counts[6]
Counts7 = Counts[7]
i_1_ini = 0
i_1 = 36
i_2_ini = 0
i_2 = 24
f_1 = 18
f_2 = 30
scale_1 = 0.92
scale_2 = 0.98
#Merged_freqs_test = [f-f_2 for f in Freqs6[i_2_ini:i_2]]+[f-f_1 for f in Freqs5[i_1_ini:i_1]]+Freqs7
#plt.plot(Merged_freqs_test,'o')
Merged_freqs = [f-f_2 for f in Freqs6[0:i_2]]+[f-f_1 for f in Freqs5[0:i_1]]+Freqs7
Merged_counts = [scale_2*c for c in Counts6[0:i_2]]+[scale_1*c for c in Counts5[0:i_1]]+list(Counts7)
Merged_freqs_rescaled = np.linspace(np.min(Merged_freqs),np.max(Merged_freqs),len(Merged_freqs))
#drs = [391.5, 399.5, 405.5, 414]
drs = [370,379,385,391.5]
if False:
plt.figure()
jvec = [2]
i = 0
drive=22.1
for j in jvec:
plt.plot([f-f_1 for f in Freqs5[0:i_1]], [scale_1*c for c in Counts5[0:i_1]],'o')
plt.plot([f-f_2 for f in Freqs6[0:i_2]], [scale_2*c for c in Counts6[0:i_2]],'o')
plt.plot(Freqs7, Counts7,'o')
plt.errorbar(Merged_freqs, Merged_counts, yerr=np.sqrt(Merged_counts), fmt='o', capsize=2, markersize=2)
i = i + 1
plt.xlabel('Frecuencia (MHz)')
plt.ylabel('counts')
plt.grid()
for dr in drs:
plt.axvline(dr)
plt.axvline(dr+drive, color='red', linestyle='dashed', alpha=0.3)
plt.axvline(dr-drive, color='red', linestyle='dashed', alpha=0.3)
plt.legend()
#%% Funciones para hacer ajustes
"""
ajusto la mergeada de 3 iones
"""
import numba
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 = -20
offsetxpi = 438+correccion
DetDoppler = -35-correccion-22
gPS, gPD, = 2*np.pi*21.58e6, 2*np.pi*1.35e6
alpha = 0
drivefreq = 2*np.pi*22.135*1e6
FreqsDR = [f-offsetxpi for f in Merged_freqs]
CountsDR = Merged_counts
freqslong = np.arange(min(FreqsDR), max(FreqsDR)+FreqsDR[1]-FreqsDR[0], 0.1*(FreqsDR[1]-FreqsDR[0]))
CircPr = 1
alpha = 0
# @numba.jit
def FitEIT_MM1(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#def FitEIT_MM(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#BETA = 1.8
# SG = 0.6
# SP = 8.1
TEMP = 0.1e-3
#BETA1, BETA2, BETA3 = 0, 0, 2
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 for f in Fluorescence1])
return ScaledFluo1+OFFSET
#return ScaledFluo1
# @numba.jit
def FitEIT_MM(freqs, SG, SP, SCALE1, SCALE2, SCALE3, OFFSET, BETA1, BETA2, BETA3):
#def FitEIT_MM(freqs, SG, SP, SCALE1, OFFSET, BETA1):
#BETA = 1.8
# SG = 0.6
# SP = 8.1
TEMP = 0.1e-3
freqs = np.array(freqs)
#BETA1, BETA2, BETA3 = 0, 0, 2
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)
Detunings, Fluorescence2 = PerformExperiment_8levels_MM(SG, SP, gPS, gPD, DetDoppler, u, DopplerLaserLinewidth, ProbeLaserLinewidth, TEMP, alpha, phidoppler, titadoppler, phiprobe, titaprobe, BETA2, drivefreq, min(freqs), max(freqs)+(freqs[1]-freqs[0]), freqs[1]-freqs[0], circularityprobe=CircPr, plot=False, solvemode=1, detpvec=None)
Detunings, Fluorescence3 = PerformExperiment_8levels_MM(SG, SP, gPS, gPD, DetDoppler, u, DopplerLaserLinewidth, ProbeLaserLinewidth, TEMP, alpha, phidoppler, titadoppler, phiprobe, titaprobe, BETA3, 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 for f in Fluorescence1])
# ScaledFluo2 = np.array([f*SCALE2 for f in Fluorescence2])
# ScaledFluo3 = np.array([f*SCALE3 for f in Fluorescence3])
# return ScaledFluo1+ScaledFluo2+ScaledFluo3+OFFSET
Fluorescence1 = np.array(Fluorescence1)
Fluorescence2 = np.array(Fluorescence2)
Fluorescence3 = np.array(Fluorescence3)
return SCALE1*Fluorescence1+SCALE2*Fluorescence2+SCALE3*Fluorescence3+OFFSET
#%% Acá hay un ajuste del de 3 iones que da razonable
# FreqsDR, CountsDR = np.array(FreqsDR) , np.array(CountsDR)
# freqs = np.array(FreqsDR)
# t0 = time.time()
# print("Arranamos FIT 1")
# par_ini = [0.65, 7.06, 86070, 3917, 1.64]
# popt_1ions, pcov_1ions = curve_fit(FitEIT_MM1, FreqsDR, CountsDR, p0=par_ini, bounds=((0, 0, 0, 0, 0), (2, 20, 5e8, 7e3, 10)))
# print(f"time: {round(time.time()-t0,1)} seg")
# par_ini = popt_1ions.tolist()[:2]+[1/3]*2 + [0.1] +popt_1ions.tolist()[-1:]*3
# bounds = ((0, 0, 0, 0, 0, 0, 0, 0),
# (2, 20, 1, 1, 1, 10, 10, 10))
# # bounds = (-np.inf, np.inf)
# print("Arranamos FIT 3")
# def fun(freqs, SG, SP, SCALE1, SCALE2, OFFSET, BETA1, BETA2, BETA3):
# SCALE3 = max(1-SCALE2-SCALE1-OFFSET,0)
# Fluorescence3 = FitEIT_MM(freqs, SG, SP, SCALE1, SCALE2, SCALE3, OFFSET, BETA1, BETA2, BETA3)
# try:
# l1.set_ydata(Fluorescence3/sum(Fluorescence3)*sum(CountsDR))
# except:
# print("Primer corrida")
# finally:
# print(SG, SP, SCALE1, SCALE2, OFFSET, BETA1, BETA2, BETA3)
# plt.pause(0.1)
# return Fluorescence3/sum(Fluorescence3)
# plt.figure()
# plt.errorbar(freqs, CountsDR, yerr=2*np.sqrt(CountsDR), fmt='o', color='darkgreen', alpha=0.5, capsize=2, markersize=2)
# l1, =plt.plot(freqs, fun(freqs, *par_ini )*sum(CountsDR))
# plt.draw()
# popt_3ions, pcov_3ions = curve_fit(fun, FreqsDR, CountsDR/CountsDR.sum(), p0=par_ini,bounds=bounds )
# print(f"time: {round(time.time()-t0,1)} seg")
# SG, SP, SCALE1, SCALE2, OFFSET, BETA1, BETA2, BETA3 = popt_3ions
# SCALE3 = max(1-SCALE2-SCALE1-OFFSET,0)
# # Fluorescence= FitEIT_MM(freqs, SG, SP, SCALE1, SCALE2, SCALE3, OFFSET, BETA1, BETA2, BETA3)
# Fluorescence= fun(freqs, SG, SP, SCALE1, SCALE2, OFFSET, BETA1, BETA2, BETA3)
# plt.plot(freqs, Fluorescence*sum(CountsDR))
# if False:
# popt_3ions = np.array([0.70790618, 7.41165226, 0.1707309 , 0.13974759, 0.02322333,
# 3.69298275, 3.68847693, 1.36216489])
#%% Armo algo para meter parametros fijos
FreqsDR, CountsDR = np.array(FreqsDR) , np.array(CountsDR)
t0 = time.time()
print("Arranamos FIT 1")
par_ini = [0.65, 7.06, 86070, 3917, 1.64]
popt_1ions, pcov_1ions, info1 = fit(FitEIT_MM1, FreqsDR, CountsDR, p0=par_ini, bounds=((0, 0, 0, 0, 0), (2, 20, 5e8, 7e3, 10)))
print(f"time: {round(time.time()-t0,1)} seg")
par_ini = popt_1ions.tolist()[:2]+[1/3]*2 + [0.1] +popt_1ions.tolist()[-1:]*3
bounds = ((0, 0, 0, 0, 0, 0, 0, 0),
(2, 20, 1, 1, 1, 10, 10, 10))
def fun(freqs, SG, SP, SCALE1, SCALE2, OFFSET, BETA1, BETA2, BETA3):
SCALE3 = max(1-SCALE2-SCALE1-OFFSET,0)
Fluorescence3 = FitEIT_MM(freqs, SG, SP, SCALE1, SCALE2, SCALE3, OFFSET, BETA1, BETA2, BETA3)
return Fluorescence3/sum(Fluorescence3)
par_ini = np.array([0.64872536, 7.06495239, 0.9985873 , 0.9985873 , 0.09090475,
1.63941009, 1.63941006, 1.80228481])
SG, SP, SCALE1, SCALE2, OFFSET, BETA1, BETA2, BETA3 = par_ini
SCALE1, SCALE2 = 0.3, 0.3
par_ini = SG, SP, SCALE1, SCALE2, OFFSET, BETA1, BETA2, BETA3
t0 = time.time()
print("Arranamos FIT 1")
popt_3ions, pcov_3ions, info2 = fit(fun, FreqsDR, CountsDR/sum(CountsDR), fixed="SG SP", plot=True, p0=par_ini, bounds=bounds)
print(f"time: {round(time.time()-t0,1)} seg")
par_ini = popt_3ions
t0 = time.time()
print("Arranamos FIT 1")
popt_3ions, pcov_3ions, info2 = fit(fun, FreqsDR, CountsDR/sum(CountsDR), fixed="", plot=True, p0=par_ini, bounds=bounds)
print(f"time: {round(time.time()-t0,1)} seg")
analisis/plots/20250516_CPT_bladetrap/lolo_analisis_superajuste.py deleted 100644 → 0
View file @ f9b573e5
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Ploteo de datos y ajustes del barrido en voltaje del Endcap
equivalente a pasear el ion en algun entorno del punto de compensación ideal
y ver los efectos de micromoción (y tal vez temperatura)
@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}')
# Funciones auxiliares
from scipy.stats.distributions import t,chi2
def estadistica(datos_x,datos_y,modelo,pcov,parametros,nombres=None,alpha=0.05):
if nombres is None:
nombres = [ f'{j}' for j in range(len(parametros)) ]
# Cantidad de parámetros
P = len(parametros)
# Número de datos
N = len(datos_x)
# Grados de libertas (Degrees Of Freedom)
dof = N-P-1
# Cauculamos coordenadas del modelo
# modelo_x = datos_x if modelo_x_arr is None else modelo_x_arr
# modelo_y = modelo( modelo_x, *parametros )
# Predicción del modelo para los datos_x medidos
prediccion_modelo = modelo( datos_x, *parametros )
# Calculos de cantidades estadísticas relevantes
COV = pcov # Matriz de Covarianza
SE = np.sqrt(np.diag( COV )) # Standar Error / Error estandar de los parámetros
residuos = datos_y - prediccion_modelo # diferencia enrte el modelo y los datos
SSE = sum(( residuos )**2 ) # Resitual Sum of Squares
SST = sum(( datos_y - np.mean(datos_y))**2) # Total Sum of Squares
# http://en.wikipedia.org/wiki/Coefficient_of_determination
# Expresa el porcentaje de la varianza que logra explicar el modelos propuesto
Rsq = 1 - SSE/SST # Coeficiente de determinación
Rsq_adj = 1-(1-Rsq) * (N-1)/(N-P-1) # Coeficiente de determinación Ajustado
# https://en.wikipedia.org/wiki/Pearson_correlation_coefficient#In_least_squares_regression_analysis
# Expresa la correlación que hay entre los datos y la predicción del modelo
r_pearson = np.corrcoef( datos_y , prediccion_modelo )[0,1]
# Reduced chi squared
# https://en.wikipedia.org/wiki/Reduced_chi-squared_statistic
chi2_ = sum( residuos**2 )/N
chi2_red = sum( residuos**2 )/(N-P)
# Chi squared test
chi2_test = sum( residuos**2 / abs(prediccion_modelo) )
# p-value del ajuste
p_val = chi2(dof).cdf( chi2_test )
sT = t.ppf(1.0 - alpha/2.0, N - P ) # student T multiplier
CI = sT * SE # Confidence Interval
print('R-squared ',Rsq)
print('R-sq_adjusted',Rsq_adj)
print('chi2 ',chi2_)
print('chi2_reduced ',chi2_red)
print('chi2_test ',chi2_test)
print('r-pearson ',r_pearson)
print('p-value ',p_val)
print('')
print('Error Estandard (SE):')
for i in range(P):
print(f'parametro[{nombres[i]:>5s}]: ' , parametros[i], ' ± ' , SE[i])
print('')
print('Intervalo de confianza al '+str((1-alpha)*100)+'%:')
for i in range(P):
print(f'parametro[{nombres[i]:>5s}]: ' , parametros[i], ' ± ' , CI[i])
return dict(R2=Rsq,R2_adj=Rsq_adj,chi2=chi2_,chi2_red=chi2_red,
chi2_test=chi2_test,r=r_pearson,pvalue=p_val,
SE=SE,CI=CI)
#%%
"""
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])))
#%% Cargo parámetros fiteados de antes
PARAMETROS = np.load('analisis_superajuste_PARAMETROS.npz',allow_pickle=True)
for var_name in PARAMETROS.keys():
globals()[var_name] = PARAMETROS[var_name]
print(f'loaded: {var_name}')
if False:
# Esto es para correr en caso de necesidad de limpiar todos los vectores de parametros
print('Limpio los vectores de parámetros')
for var in [ kk for kk in globals().keys() if kk.endswith('_vec') ]:
print(f'del {var}')
del(globals()[var])
#%% Definiciones de Numba
@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
@jit
def FitEIT_MM_single(Freqs, offset, DetDoppler, SG, SP, SCALE1, OFFSET, BETA1, TEMP):
"Esta verison de la función devuelve sólo el eje y, para usar de modelo en un ajuste"
return FitEIT_MM_single_plot(Freqs, offset, DetDoppler, SG, SP, SCALE1, OFFSET, BETA1, TEMP)[0]
param_names = 'offset DetDoppler SG SP SCALE1 OFFSET BETA1 TEMP'.split()
#%%
"""
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 = []
Estadistica_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)
# estadistica(datos_x,datos_y,modelo,pcov,parametros,nombres=None,alpha=0.05)
est_tmp = estadistica(FreqsDR,CountsDR,FitEIT_MM_single,pcov_3_SA,popt_3_SA,
nombres=param_names,alpha=0.05)
Estadistica_vec.append(est_tmp)
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)
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
#%% Graficamos todos los fiteos
# 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)')
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
#%% Inspección de parámetros
param_names = 'offset DetDoppler SG SP SCALE1 OFFSET BETA1 TEMP'.split()
err_vecs = np.array([ np.sqrt(np.diag(el)) for el in pcov_SA_vec ])
num_med = np.arange(len(pcov_SA_vec)) +1
r2_values = np.array([ el['R2_adj'] for el in Estadistica_vec ])
fig, axx = plt.subplots( len(popt_SA_vec[0])+1,1, figsize=(13,8) , constrained_layout=True, sharex=True , sharey=False )
fig.set_constrained_layout_pads(w_pad=2/72, h_pad=2/72, hspace=0, wspace=0)
for ax,param_vec,err_vec,par_name in zip(axx,popt_SA_vec.T,err_vecs.T,param_names) :
ax.plot(num_med, param_vec, '.-')
ax.errorbar( num_med, param_vec, yerr=err_vec,
fmt='s', mfc='none', elinewidth = 1, capsize=3, ms=1)
ax.grid(True, ls=":", color='lightgray')
ax.set_ylabel(par_name)
ax=axx[-1]
ax.plot( num_med , r2_values, '.-')
ax.set_ylabel(r'$R^2$')
ax.grid(True, ls=":", color='lightgray')
fig.align_ylabels()
ax.set_xticks(num_med)
ax.set_xlabel('Num. de medición')
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
#%% Endcap hiperbola (con residuos)
"""
Veamos cómo varía el Beta con el voltaje del endcap
"""
import seaborn as sns
paleta = sns.color_palette("rocket")
I = slice(None,9)
voltages_dcA = Voltages[0][SelectedCurveVec]
def lineal(x,a,b):
return a*x+b
def hiperbola(x,a,y0,b,x0):
"""
Hiperbola de ecuación:
1 =(y-y0)²/a² - (x-x0)²/b²
"""
return a*np.sqrt(((x-x0)**2+b**2))+y0
es_hiperbola = False
# par_inicial = (100,0.1,1,-0.15)
# a y0 b x0
par_inicial = (12,0.1,1,-0.13)
popthip,pcovhip = curve_fit(hiperbola,voltages_dcA[I],Betas_vec[I],p0=par_inicial)
xhip = np.linspace(-0.23,0.005,200)
# plt.figure()
# plt.errorbar(voltages_dcA[I],Betas_vec[I],yerr=ErrorBetas_vec[I],fmt='o',capsize=5,markersize=5,color=paleta[1])
# plt.plot(xhip,hiperbola(xhip,*popthip))
# # plt.plot(xhip,hiperbola(xhip,*par_inicial),'--',color='red')
# plt.xlabel('Endcap voltage (V)')
# plt.ylabel('Modulation factor')
# plt.grid()
fig, axx = plt.subplots( 2, figsize=(10,7) ,
constrained_layout=True, sharex=True,
gridspec_kw=dict(height_ratios=[10,2]))
fig.set_constrained_layout_pads(w_pad=2/72, h_pad=2/72, hspace=0, wspace=0)
ax = axx[0]
ax.errorbar(voltages_dcA[I],Betas_vec[I],yerr=ErrorBetas_vec[I],fmt='o',capsize=5,markersize=5,color=paleta[1])
ax.plot(xhip,hiperbola(xhip,*popthip))
# plt.plot(xhip,hiperbola(xhip,*par_inicial),'--',color='red')
ax.set_ylabel('Modulation factor')
ax = axx[1]
ax.errorbar(voltages_dcA[I],Betas_vec[I]-hiperbola(voltages_dcA[I],*popthip),
yerr=ErrorBetas_vec[I],fmt='o',capsize=5,markersize=5,color=paleta[1])
ax.set_ylabel('Res.')
ax.set_xlabel('Endcap voltage (V)')
for ax in axx:
ax.grid(True, ls=":", color='lightgray')
print([t*1e3 for t in Temp_vec])
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
#%% Este hay que armarlo aún
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)
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
#%% Ajuste de los betas y la temperatura
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('analisis_superajuste_PARAMETROS.npz', **GUARDAR )
analisis/plots/20250516_CPT_bladetrap/lolo_graficos_pub.py deleted 100644 → 0
View file @ f9b573e5
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Ploteo de datos y ajustes del barrido en voltaje del Endcap
equivalente a pasear el ion en algun entorno del punto de compensación ideal
y ver los efectos de micromoción (y tal vez temperatura)
@author: lolo
"""
import h5py
import matplotlib.pyplot as plt
import numpy as np
from scipy.optimize import curve_fit
from time import time
from numba import jit,njit
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
#%% Carga de datos y librerías auxiliares
## Librerias ##################################################################
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}')
# Funciones auxiliares
from scipy.stats.distributions import t,chi2
import inspect
def estadistica(datos_x,datos_y,modelo,pcov,parametros,nombres=None,alpha=0.05):
if nombres is None:
nombres = inspect.getfullargspec(modelo).args[1:]
# Cantidad de parámetros
P = len(parametros)
# Número de datos
N = len(datos_x)
# Grados de libertas (Degrees Of Freedom)
dof = N-P-1
# Cauculamos coordenadas del modelo
# modelo_x = datos_x if modelo_x_arr is None else modelo_x_arr
# modelo_y = modelo( modelo_x, *parametros )
# Predicción del modelo para los datos_x medidos
prediccion_modelo = modelo( datos_x, *parametros )
# Calculos de cantidades estadísticas relevantes
COV = pcov # Matriz de Covarianza
SE = np.sqrt(np.diag( COV )) # Standar Error / Error estandar de los parámetros
residuos = datos_y - prediccion_modelo # diferencia enrte el modelo y los datos
SSE = sum(( residuos )**2 ) # Resitual Sum of Squares
SST = sum(( datos_y - np.mean(datos_y))**2) # Total Sum of Squares
# http://en.wikipedia.org/wiki/Coefficient_of_determination
# Expresa el porcentaje de la varianza que logra explicar el modelos propuesto
Rsq = 1 - SSE/SST # Coeficiente de determinación
Rsq_adj = 1-(1-Rsq) * (N-1)/(N-P-1) # Coeficiente de determinación Ajustado
# https://en.wikipedia.org/wiki/Pearson_correlation_coefficient#In_least_squares_regression_analysis
# Expresa la correlación que hay entre los datos y la predicción del modelo
r_pearson = np.corrcoef( datos_y , prediccion_modelo )[0,1]
# Reduced chi squared
# https://en.wikipedia.org/wiki/Reduced_chi-squared_statistic
chi2_ = sum( residuos**2 )/N
chi2_red = sum( residuos**2 )/(N-P)
# Chi squared test
chi2_test = sum( residuos**2 / abs(prediccion_modelo) )
# p-value del ajuste
p_val = chi2(dof).cdf( chi2_test )
sT = t.ppf(1.0 - alpha/2.0, N - P ) # student T multiplier
CI = sT * SE # Confidence Interval
return dict(R2=Rsq,R2_adj=Rsq_adj,chi2=chi2_,chi2_red=chi2_red,
chi2_test=chi2_test,r=r_pearson,pvalue=p_val,
SE=SE,CI=CI,parametros=parametros,nombres=nombres,alpha=alpha,N=N)
def print_estadistica(dic):
pars = 'N,pvalue,R2,R2_adj,chi2,chi2_red,chi2_test,r,SE,CI,parametros,nombres,alpha'.split(',')
N,pvalue,R2,R2_adj,chi2,chi2_red,chi2_test,r,SE,CI,parametros,nombres,alpha = [ dic[k] for k in pars ]
P = len(parametros)
print('R-squared ',R2)
print('R-sq_adjusted',R2_adj)
print('chi2 ',chi2)
print('chi2_reduced ',chi2_red)
print('chi2_test ',chi2_test)
print('r-pearson ',r)
print('p-value ',pvalue)
print('')
print('Error Estandard (SE):')
for i in range(P):
print(f'parametro[{nombres[i]:>5s}]: ' , parametros[i], ' ± ' , SE[i])
print('')
print('Intervalo de confianza al '+str((1-alpha)*100)+'%:')
for i in range(P):
print(f'parametro[{nombres[i]:>5s}]: ' , parametros[i], ' ± ' , CI[i])
## Datos ######################################################################
"""
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])))
## Cargo parámetros fiteados de antes #########################################
PARAMETROS = np.load('analisis_superajuste_PARAMETROS.npz',allow_pickle=True)
for var_name in PARAMETROS.keys():
globals()[var_name] = PARAMETROS[var_name]
print(f'loaded: {var_name}')
if False:
# Esto es para correr en caso de necesidad de limpiar todos los vectores de parametros
print('Limpio los vectores de parámetros')
for var in [ kk for kk in globals().keys() if kk.endswith('_vec') ]:
print(f'del {var}')
del(globals()[var])
## Definiciones de Numba ######################################################
@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
@jit
def FitEIT_MM_single(Freqs, offset, DetDoppler, SG, SP, SCALE1, OFFSET, BETA1, TEMP):
"Esta verison de la función devuelve sólo el eje y, para usar de modelo en un ajuste"
return FitEIT_MM_single_plot(Freqs, offset, DetDoppler, SG, SP, SCALE1, OFFSET, BETA1, TEMP)[0]
param_names = 'offset DetDoppler SG SP SCALE1 OFFSET BETA1 TEMP'.split()
#%%
"""
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
voltages_dcA = Voltages[0][SelectedCurveVec]
def hiperbola(x,a,y0,b,x0):
"""
Hiperbola de ecuación:
1 =(y-y0)²/a² - (x-x0)²/b²
"""
return a*np.sqrt(((x-x0)**2+b**2))+y0
def hiperbola2(x,a,b,x0):
"""
Hiperbola de ecuación:
1 =(y)²/a² - (x-x0)²/b²
"""
return a*np.sqrt(((x-x0)**2+b**2))
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
#%% Graficamos todos los fiteos
# 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)')
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
#%% Gráfico central OPCION A
import matplotlib.pyplot as plt
plt.rcParams.update({'font.size': 8})
# Para esto hace falta:
# sudo apt install dvipng
plt.rcParams['text.usetex']=True
# plt.rcParams['text.latex.unicode']=True
# params= {'text.latex.preamble' : [r'\usepackage{amsmath}']}
# plt.rcParams.update(params)
# fig, axx = plt.subplots( 3,4, figsize=(8.6/2.54*2,4) , constrained_layout=True, sharex=True , sharey=True )
# fig.set_constrained_layout_pads(w_pad=2/72, h_pad=2/72, hspace=0, wspace=0)
## Layout ###########################
from matplotlib import gridspec
figsize=(8.6/2.54*2,4)
width_ratios=[1,2,1]
fig, axx = plt.subplots(ncols=3, nrows=3, sharex=True, sharey=False, constrained_layout=True,
figsize=figsize, gridspec_kw=dict(width_ratios=width_ratios))
fig.set_constrained_layout_pads(w_pad=1/72, h_pad=1/72, hspace=0, wspace=0)
fig.set_constrained_layout_pads(w_pad=0, h_pad=0, hspace=0, wspace=0)
gs = axx[0, 0].get_gridspec()
# remove the underlying axes
for ax in axx[:,1]:
ax.remove()
gs2 = gridspec.GridSpec(4,3, width_ratios=width_ratios,left=0.16,right=0.9,bottom=0.12)
ax_central = fig.add_subplot(gs2[:-1, 1])
ax_res = fig.add_subplot(gs2[-1, 1], sharex=ax_central,zorder=-5)
plt.setp(ax_central.get_xticklabels(), visible=False)
for ax in axx[:,2]:
ax.yaxis.tick_right()
ax.yaxis.set_label_position("right")
## CPTs ########################################
selection= (0,1,3,4,7,8)
axes_vec = (axx[0,0],axx[1,0],axx[2,0],axx[2,2],axx[1,2],axx[0,2])
tmp_datos=(Detuningsshort_vec[selection,:],
Counts_vec[selection,:],
Detuningslong_vec[selection,:],
FittedCounts_vec[selection,:],
np.array(selection)+1,
axes_vec,'abcfed')
for Detunings_3_SA_short,CountsDR,Detunings_3_SA_long,FittedEITpi_3_SA_long,selectedcurve,ax,le in zip(*tmp_datos):
ax.errorbar(Detunings_3_SA_short, CountsDR, yerr=2*np.sqrt(CountsDR),
fmt='o', color='darkgreen', alpha=0.2, capsize=2,
markersize=2, label='measured data')
ax.plot(Detunings_3_SA_long, FittedEITpi_3_SA_long,
color='black', linewidth=2, label=f'micromotion model', alpha=0.7)
ax.grid(True, ls=":", color='lightgray')
ax.set_yticklabels([ f"{int(y/1000)}k" for y in ax.get_yticks() ])
ax.set_title(f'({le})', x=0.1, y=0.78, color='gray')
# ax.legend()
# ax.legend(ax.get_legend_handles_labels()[0], f'{selectedcurve}')
print(f'listo med {selectedcurve}')
axx[2,0].set_xlabel('Detuning [MHz]')
axx[2,2].set_xlabel('Detuning [MHz]')
axx[1,0].set_ylabel('Counts')
axx[1,2].set_ylabel('Counts')
## Grafico central ###############################
I = slice(None,9)
par_inicial = (12,0.1,1,-0.13)
param,pcov = curve_fit(hiperbola,voltages_dcA[I],Betas_vec[I],p0=par_inicial)
x_hip = np.linspace(-0.23,0.005,200)
EST=estadistica(voltages_dcA[I],Betas_vec[I],hiperbola,pcov,param,nombres=None,alpha=0.05)
print_estadistica(EST)
I = slice(None,9)
par_inicial = (12,1,-0.13)
param2,pcov2 = curve_fit(hiperbola2,voltages_dcA[I],Betas_vec[I],p0=par_inicial)
x_hip = np.linspace(-0.23,0.005,200)
EST2=estadistica(voltages_dcA[I],Betas_vec[I],hiperbola2,pcov2,param2,nombres=None,alpha=0.05)
print_estadistica(EST2)
ax = ax_central
ax.errorbar(voltages_dcA[I],Betas_vec[I],yerr=ErrorBetas_vec[I],fmt='o',
capsize=4,markersize=3,color='C3', label=r'fitted $\beta$')
ax.plot(x_hip,hiperbola(x_hip,*param),color='C0', label=r'hyperbola model')
ax.plot(x_hip,hiperbola2(x_hip,*param2),'--',color='C4', label=r'hyperbola model')
ax.set_ylabel(r'Modulation factor $\beta$', labelpad=-5)
ax.set_ylim(-0.05,3)
ax.set_xlim(-0.22,0)
ax.set_title(f'(g)', x=0.95, y=0.006, color='gray')
ax = ax_res
ax.errorbar(voltages_dcA[I],Betas_vec[I]-hiperbola(voltages_dcA[I],*param),
yerr=ErrorBetas_vec[I],fmt='o',capsize=4,markersize=3,color='C3')
ax.errorbar(voltages_dcA[I]+0.005,Betas_vec[I]-hiperbola2(voltages_dcA[I],*param2),
yerr=ErrorBetas_vec[I],fmt='o',capsize=4,markersize=3,color='C4')
ax.axhline( 0 , color='C0')
ax.set_ylabel('Res.', labelpad=-5)
ax.set_xlabel('Endcap voltage [V]')
ax.set_title(f'(h)', x=0.95, y=0.72, color='gray')
ax_res.get_xticklabels()[-1].set_visible(False)
for ax in [ax_central,ax_res]:
ax.grid(True, ls=":", color='lightgray')
# print([t*1e3 for t in Temp_vec])
# Anotaciones ##########################3
for jj,ax in zip(selection,axes_vec):
Axes_x = 1 if axx.flatten().tolist().index(ax)%3==0 else 0
ax_central.annotate("",
xy=(ax_central.get_lines()[0].get_xdata()[jj], ax_central.get_lines()[0].get_ydata()[jj]),
xycoords=ax_central.transData,
xytext=(Axes_x, 0.5), textcoords=ax.transAxes,
arrowprops=dict(arrowstyle="<-",connectionstyle="arc3,rad=-0.2", color='gray', alpha=0.5),
zorder=-2)
# Leyenda ##############################
# h1, l1 = ax_central.get_legend_handles_labels()
# h2, l2 = axx[0,0].get_legend_handles_labels()
# ax_central.legend(h1+h2, l1+l2, loc='upper left')
fig.align_ylabels([ax_central,ax_res])
fig.tight_layout()
# fig.savefig('grafico_central_opcion_A.png', dpi=300)
# fig.savefig('grafico_central_opcion_A.pdf')
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
#%% Gráfico central OPCION B
import matplotlib.pyplot as plt
plt.rcParams.update({'font.size': 8})
# Para esto hace falta:
# sudo apt install dvipng
plt.rcParams['text.usetex']=True
# plt.rcParams['text.latex.unicode']=True
# params= {'text.latex.preamble' : [r'\usepackage{amsmath}']}
# plt.rcParams.update(params)
# fig, axx = plt.subplots( 3,4, figsize=(8.6/2.54*2,4) , constrained_layout=True, sharex=True , sharey=True )
# fig.set_constrained_layout_pads(w_pad=2/72, h_pad=2/72, hspace=0, wspace=0)
voltages_dcA = Voltages[0][SelectedCurveVec]*4
## Layout ###########################
from matplotlib import gridspec
figsize=(8.6/2.54*2,4)
width_ratios=[1,2,1]
fig, axx = plt.subplots(ncols=3, nrows=3, sharex=True, sharey=False, constrained_layout=True,
figsize=figsize, gridspec_kw=dict(width_ratios=width_ratios))
fig.set_constrained_layout_pads(w_pad=1/72, h_pad=1/72, hspace=0, wspace=0)
fig.set_constrained_layout_pads(w_pad=0, h_pad=0, hspace=0, wspace=0)
gs = axx[0, 0].get_gridspec()
# remove the underlying axes
for ax in axx[:,1]:
ax.remove()
gs2 = gridspec.GridSpec(4,3, width_ratios=width_ratios,left=0.16,right=0.9,bottom=0.12)
ax_central = fig.add_subplot(gs2[1:-1, 1])
ax_dib = fig.add_subplot(gs2[0, 1], sharex=ax_res, zorder=-20)
ax_dib.spines[:].set_visible(False)
ax_dib.xaxis.set_tick_params(labelbottom=False)
ax_dib.yaxis.set_tick_params(labelleft=False)
ax_dib.set_xticks([])
ax_dib.set_yticks([])
ax_res = fig.add_subplot(gs2[-1, 1], sharex=ax_central,zorder=-5)
plt.setp(ax_central.get_xticklabels(), visible=False)
for ax in axx[:,2]:
ax.yaxis.tick_right()
ax.yaxis.set_label_position("right")
## CPTs ########################################
selection= (0,1,3,4,7,8)
axes_vec = (axx[0,0],axx[1,0],axx[2,0],axx[2,2],axx[1,2],axx[0,2])
tmp_datos=(Detuningsshort_vec[selection,:],
Counts_vec[selection,:],
Detuningslong_vec[selection,:],
FittedCounts_vec[selection,:],
np.array(selection)+1,
axes_vec,'abcfed')
from scipy.signal import savgol_filter as savgol
for Detunings_3_SA_short,CountsDR,Detunings_3_SA_long,FittedEITpi_3_SA_long,selectedcurve,ax,le in zip(*tmp_datos):
# ax.errorbar(Detunings_3_SA_short, CountsDR, yerr=2*np.sqrt(CountsDR),
# fmt='o', color='darkgreen', alpha=0.2, capsize=2,
# markersize=2, label='measured data')
yerr_max = CountsDR+2*np.sqrt(CountsDR)
yerr_min = CountsDR-2*np.sqrt(CountsDR)
yerr_max = savgol(yerr_max,3,2)
yerr_min = savgol(yerr_min,3,2)
ax.fill_between(Detunings_3_SA_short, yerr_max, yerr_min,
color='darkgreen', alpha=0.2)
ax.plot(Detunings_3_SA_short, CountsDR, '.', ms=2,
color='darkgreen', alpha=0.9)
ax.plot(Detunings_3_SA_long, FittedEITpi_3_SA_long,
color='black', linewidth=1, label=f'micromotion model', alpha=0.7)
ax.grid(True, ls=":", color='lightgray')
ax.set_yticklabels([ f"{int(y/1000)}k" for y in ax.get_yticks() ])
ax.set_title(f'({le})', x=0.1, y=0.78, color='gray')
# ax.legend()
# ax.legend(ax.get_legend_handles_labels()[0], f'{selectedcurve}')
print(f'listo med {selectedcurve}')
axx[2,0].set_xlabel('Detuning [MHz]')
axx[2,2].set_xlabel('Detuning [MHz]')
axx[1,0].set_ylabel('Counts')
axx[1,2].set_ylabel('Counts')
## Grafico central ###############################
I = slice(None,9)
par_inicial = (12,1,-0.13)
param,pcov = curve_fit(hiperbola2,voltages_dcA[I],Betas_vec[I],p0=par_inicial)
x_hip = np.linspace(-0.23*4,0.005,200)
EST=estadistica(voltages_dcA[I],Betas_vec[I],hiperbola2,pcov,param,nombres=None,alpha=0.05)
print_estadistica(EST)
ax = ax_central
ax.errorbar(voltages_dcA[I],Betas_vec[I],yerr=ErrorBetas_vec[I],fmt='o',
capsize=4,markersize=3,color='C3', label=r'fitted $\beta$')
ax.plot(x_hip,hiperbola2(x_hip,*param),color='C0', label=r'hyperbola model')
ax.set_ylabel(r'Modulation factor $\beta$', labelpad=-5)
ax.set_ylim(-0.05,3)
ax.set_xlim(-0.22*4,0)
ax.set_title(f'(h)', x=0.95, y=0.006, color='gray')
ax = ax_res
ax.errorbar(voltages_dcA[I],Betas_vec[I]-hiperbola2(voltages_dcA[I],*param),
yerr=ErrorBetas_vec[I],fmt='o',capsize=4,markersize=3,color='C3')
ax.axhline( 0 , color='C0')
ax.set_ylabel('Res.', labelpad=-5)
ax.set_xlabel('Endcap voltage [V]')
ax.set_title(f'(i)', x=0.95, y=0.72, color='gray')
ax_res.get_xticklabels()[-1].set_visible(False)
for ax in [ax_central,ax_res]:
ax.grid(True, ls=":", color='lightgray')
# print([t*1e3 for t in Temp_vec])
# Anotaciones ##########################
for jj,ax in zip(selection,axes_vec):
rad = 0.2 if jj==8 or jj==3 else -0.2
Axes_x = 1 if axx.flatten().tolist().index(ax)%3==0 else 0
ax_central.annotate("",
xy=(ax_central.get_lines()[0].get_xdata()[jj], ax_central.get_lines()[0].get_ydata()[jj]),
xycoords=ax_central.transData,
xytext=(Axes_x, 0.5), textcoords=ax.transAxes,
arrowprops=dict(arrowstyle="<-",connectionstyle=f"arc3,rad={rad}", color='gray', alpha=0.5),
zorder=-2)
# dibujos ##########################
# a, y0, b, x0 = param
a, b, x0 = param
y0 = -0.212663964284234
ax_dib.plot( x_hip , x_hip*0+y0 , color='gray')
ax_dib.set_ylim(-y0*0.1,y0*1.2)
ax_dib.plot( [x0], [0], 'x' )
ax_dib.text(x0+0.01, 0, 'trap center', fontsize = 8, ha='left',va='center', color='C0')
ax_dib.plot( voltages_dcA[I], voltages_dcA[I]*0+y0, '.' )
for Vx in voltages_dcA[I]:
ax_dib.plot( [Vx,x0], [y0,0], ':', color='lightgray', lw=1)
import matplotlib.patches as mpatches
arr = mpatches.FancyArrowPatch((x0, y0*1.1), (x0/2, y0*1.1),
arrowstyle='->,head_width=.15', mutation_scale=10, color='C1')
ax_dib.add_patch(arr)
ax_dib.annotate("ion position", (.5, .5), xycoords=arr, ha='center', va='bottom', color='C1')
# r
arr = mpatches.FancyArrowPatch((x0,0), ( voltages_dcA[5], y0),
arrowstyle='->,head_width=.25', mutation_scale=10, color='black')
ax_dib.add_patch(arr)
ax_dib.annotate(r"$\vec r$", (.5, .5), xycoords=arr, ha='center', va='bottom', color='black')
bbox = ax_dib.get_window_extent().transformed(fig.dpi_scale_trans.inverted())
aspect_ratio = bbox.width / bbox.height
ax_dib.plot(0.1, 0, ">k", transform=ax_dib.transAxes, clip_on=False, ms=3, lw=1)
ax_dib.plot([0,0.1], [0,0], "-k", transform=ax_dib.transAxes, clip_on=False, lw=1)
ax_dib.plot(0, 0.1*aspect_ratio, "^k", transform=ax_dib.transAxes, clip_on=False, ms=3, lw=1)
ax_dib.plot([0,0], [0,0.1*aspect_ratio], "-k", transform=ax_dib.transAxes, clip_on=False, lw=1)
ax_dib.text(0.1, 0.13, "x", transform=ax_dib.transAxes,
fontsize=8, va='top')
ax_dib.text(0.03, 0.13*aspect_ratio, "y", transform=ax_dib.transAxes,
fontsize=8, va='top')
ax_dib.set_title(f'(g)', x=0.95, y=0.5, color='gray')
# Leyenda ##############################
# h1, l1 = ax_central.get_legend_handles_labels()
# h2, l2 = axx[0,0].get_legend_handles_labels()
# ax_central.legend(h1+h2, l1+l2, loc='upper left')
fig.align_ylabels([ax_central,ax_res])
fig.tight_layout()
fig.savefig('grafico_central_opcion_B.png', dpi=300)
fig.savefig('grafico_central_opcion_B.pdf')
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
#%% Grafico de RF Heating y temp
from scipy.special import jv
import seaborn as sns
paleta = sns.color_palette("rocket")
"""
Por si no quiero correr todo de nuevo,
los vectores del grafico son estos
Betas_vec = [1.5063238290734708, 0.901629236713104, 0.5280696724084589, 0.32347492275852474, 0.581352716043054, 0.9655636982224044, 1.6186284658310388, 2.1055294502704047, 2.548008896618126, 3.130666254406969, 3.1665902866111795]
ErrorBetas_vec = [0.034184989382424424, 0.04827624922956116, 0.054704762232730145, 0.06708926613955447, 0.047526038227881616, 0.05210296764897591, 0.03630379898730492, 0.034064021790261405, 0.050405155463964166, 0.10838273932613045, 0.2828605486723366]
Temp_vec = [0.0017384071359855713, 0.0006380145223733723, 0.0007457923288975645, 0.0006442238745101592, 0.0006270680749881928, 0.0008995355804703286, 0.0017799223158665226, 0.002941224610139307, 0.008378628768005558, 0.026250695067608725, 0.09869604401089357]
ErrorTemp_vec = [0.0004229476096548473, 0.00014439375508413987, 0.00013204015146487435, 9.307939678673377e-05, 0.000100129717662808, 0.0001841318633900307, 0.0003595040837509976, 0.0005950353892849986, 0.001866844309182069, 0.012656306714647434, 0.13143081065882864]
"""
Betas_vec = [1.5063238290734708, 0.901629236713104, 0.5280696724084589, 0.32347492275852474, 0.581352716043054, 0.9655636982224044, 1.6186284658310388, 2.1055294502704047, 2.548008896618126, 3.130666254406969, 3.1665902866111795]
ErrorBetas_vec = [0.034184989382424424, 0.04827624922956116, 0.054704762232730145, 0.06708926613955447, 0.047526038227881616, 0.05210296764897591, 0.03630379898730492, 0.034064021790261405, 0.050405155463964166, 0.10838273932613045, 0.2828605486723366]
Temp_vec = [0.0017384071359855713, 0.0006380145223733723, 0.0007457923288975645, 0.0006442238745101592, 0.0006270680749881928, 0.0008995355804703286, 0.0017799223158665226, 0.002941224610139307, 0.008378628768005558, 0.026250695067608725, 0.09869604401089357]
ErrorTemp_vec = [0.0004229476096548473, 0.00014439375508413987, 0.00013204015146487435, 9.307939678673377e-05, 0.000100129717662808, 0.0001841318633900307, 0.0003595040837509976, 0.0005950353892849986, 0.001866844309182069, 0.012656306714647434, 0.13143081065882864]
def MicromotionSpectra(beta,det, gamma):
"""
Espectro de transicion doppler considerando micromocion
"""
ftrap=2*np.pi*22.1e6
#gamma=23
P = (jv(0, beta)**2)/(((det)**2)+(0.5*gamma)**2)
i = 1
#print(P)
while i <= 2: #numero de bandas de micromocion a considerar
P = P + ((jv(i, beta))**2)/((((det)+i*ftrap)**2)+(0.5*gamma)**2) + ((jv(-i, beta))**2)/((((det)-i*ftrap)**2)+(0.5*gamma)**2)
i = i + 1
#print(P)
return (gamma/2*np.pi)*P
#return P
def InverseDerivMicromotionSpectra(beta, det, gamma):
"""
La inversa de la derivada del espectro de micromocion
"""
ftrap=2*np.pi*22.1e6
#gamma=23
#det = -gamma/2
P = ((jv(0, beta)**2)/((((det)**2)+(0.5*gamma)**2)**2))*(-2*(det))
i = 1
#print(P)
while i <= 2:
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 1/((gamma/2*np.pi)*P)
#return 1/P
def FinalTemp(beta,det, C,D):
"""
Funcion para ajustar la tmeperatura final en funcion de beta.
Esta funcion tiene un parametro de ajuste C para estudiar posibles desviacioens de la temperatura
doppler teorica para un sistema de dos niveles.
"""
gamma = 23e6
hbar = 1.05e-34
kb = 1.38e-23
return ((hbar/(kb))*C*MicromotionSpectra(beta,2*np.pi*det,2*np.pi*gamma)+D*(beta**2))*InverseDerivMicromotionSpectra(beta, 2*np.pi*det, 2*np.pi*gamma)
def FinalTemp_withoutC(beta,det, D, forgetmicromotion=False):
"""
Funcion para ajustar la temperatura final en funcion de beta.
Esta funcion no tiene ese parametro de ajuste ya que vi que era muy cercano a 1.
"""
gamma = 23e6
hbar = 1.05e-34
kb = 1.38e-23
C = 1 #justamente aca dejo fijo ese parametro
if forgetmicromotion:
D=0
return ((hbar/(kb))*C*MicromotionSpectra(beta,2*np.pi*det,2*np.pi*gamma)+D*(beta**2))*InverseDerivMicromotionSpectra(beta, 2*np.pi*det, 2*np.pi*gamma)
#return (C*MicromotionSpectra(beta,det,gamma))*InverseDerivMicromotionSpectra(beta, det, gamma)
def FinalTemp_withoutmicromotion(beta,det):
"""
Funcion para ajustar la temperatura final en funcion de beta.
Esta funcion no tiene ese parametro de ajuste ya que vi que era muy cercano a 1.
"""
gamma = 23e6
hbar = 1.05e-34
kb = 1.38e-23
C = 1 #justamente aca dejo fijo ese parametro
D=0
return ((hbar/(kb))*C*MicromotionSpectra(beta,2*np.pi*det,2*np.pi*gamma)+D*(beta**2))*InverseDerivMicromotionSpectra(beta, 2*np.pi*det, 2*np.pi*gamma)
#return (C*MicromotionSpectra(beta,det,gamma))*InverseDerivMicromotionSpectra(beta, det, gamma)
def CoolingLimit_noMM(det):
gamma = 2*np.pi*23e6
hbar = 1.05e-34
kb = 1.38e-23
beta = 0
rta = ((2*hbar/(2*kb))*MicromotionSpectra(beta,det,gamma))
rta *= InverseDerivMicromotionSpectra(beta, det, gamma)
return rta
DetVec = np.arange(-2*np.pi*30e6,-2*np.pi*1e6,2*np.pi*0.1e6)
Temps_limit = []
for detuning in DetVec:
Temps_limit.append(CoolingLimit_noMM(detuning))
print(f'best: {1e3*np.min(Temps_limit)}')
if False:
plt.figure()
plt.plot(DetVec*1e-6/(2*np.pi), np.array(Temps_limit)*1e3)
plt.grid()
plt.xlabel('detuning (MHz)')
plt.ylabel('final temperature')
#Dos pruebas... ignorar
#popt_rho22_balance, pcov_rho22_balance = curve_fit(FinalTemp,list(Betas_vec[:9]), [t for t in Temp_vec[:9]],p0=(-11e6,0.8,4e-28)) #esto ajusta muy bien
#popt_rho22_balance, pcov_rho22_balance = curve_fit(FinalTemp,list(Betas_vec[:10]), [t for t in Temp_vec[:10]],sigma=ErrorTemp_vec[:10],absolute_sigma=False,p0=(-10e6,0.8,5e-20)) #esto ajusta muy bien
#Buen ajuste final:
popt_rho22_balance, pcov_rho22_balance = curve_fit(FinalTemp_withoutC,list(Betas_vec[:10]), [t for t in Temp_vec[:10]],sigma=ErrorTemp_vec[:10],absolute_sigma=False,p0=(-10e6,5e-20)) #esto ajusta muy bien
popt_rho22_balance_noMM, pcov_rho22_balance_noMM = curve_fit(FinalTemp_withoutmicromotion,list(Betas_vec[:10]), [t for t in Temp_vec[:10]],sigma=ErrorTemp_vec[:10],absolute_sigma=False,p0=(-10e6)) #esto ajusta muy bien
betaslong = np.arange(0,2.8,0.01)
print(f'Min temp predicted: {FinalTemp_withoutC(betaslong,*popt_rho22_balance)[0]}')
#print(f'Min temp predicted: {1e3*FinalTemp(betaslong,*popt_rho22_balance)[0]} mK')
print(f'Detuning: {1e-6*popt_rho22_balance[0]} MHz')
print(f'params: {popt_rho22_balance}')
print(f'errores: {np.sqrt(np.diag(pcov_rho22_balance))}')
k_plot = 9
# #plt.plot(betaslong,FinalTemp_fixedall(betaslong,*popt_rho22_balance),label='Ajuste con espectro modulado')
# # plt.xlim(-0.1,1.1)
# #plt.ylim(0,1)
# #plt.axvline(minimum_voltage,linestyle='dashed',color='grey')
# #plt.axhline(0.538)
# plt.xlabel('Modulation factor')
# plt.ylabel('Temperature (mK)')
# plt.legend()
# plt.grid()
import matplotlib.pyplot as plt
plt.rcParams.update({'font.size': 8})
# Para esto hace falta:
# sudo apt install dvipng
plt.rcParams['text.usetex']=True
# Acá estimo un detuning de (10 +- 4) MHz
Detuning_range = ( 10 + 4*np.array([-1,1]) )*1e6
Doppler_limit_area = CoolingLimit_noMM(-2*np.pi*Detuning_range)*1e3
figsize=(8.6/2.54,3.5)
fig, ax = plt.subplots(1,1, constrained_layout=True,
figsize=figsize)
# fig.set_constrained_layout_pads(w_pad=1/72, h_pad=1/72, hspace=0, wspace=0)
# fig.set_constrained_layout_pads(w_pad=0, h_pad=0, hspace=0, wspace=0)
ax.errorbar(Betas_vec[:k_plot],[t*1e3 for t in Temp_vec[:k_plot]],xerr=ErrorBetas_vec[:k_plot],
yerr=[t*1e3 for t in ErrorTemp_vec[:k_plot]],fmt='o',
capsize=3,markersize=3,color='C3',zorder=10,alpha=0.9, label='inferred data')
#plt.plot(betaslong,[t*1e3 for t in FinalTemp(betaslong,*popt_rho22_balance)],label='Ajuste con espectro modulado')
ax.plot(betaslong,[t*1e3 for t in FinalTemp_withoutC(betaslong,*popt_rho22_balance)],
'-', color='C0', label='model w/ RF heating')
ax.plot(betaslong,[t*1e3 for t in FinalTemp_withoutmicromotion(betaslong,*popt_rho22_balance_noMM)],
'--', color='C1',label='model w/o RF heating')
ax.fill_between( [0, 2.8] , np.ones(2)*Doppler_limit_area[0], np.ones(2)*Doppler_limit_area[1],
color='gray', alpha=0.5, label='Doppler limit w/o micromotion')
ax.grid(True, ls=":", color='lightgray')
# ax.set_yticklabels([ f"{int(y/1000)}k" for y in ax.get_yticks() ])
# ax.set_title(f'({le})', x=0.1, y=0.78, color='gray')
ax.set_xlabel(r'Modulation factor $\beta$')
ax.set_ylabel('Temperature [mK]')
ax.set_xlim(0,2.8)
# ax.legend(fontsize=8)
# Leyenda ##############################
h1, l1 = ax.get_legend_handles_labels()
ax.legend(h1[-1:]+h1[:-1], l1[-1:]+l1[:-1], loc='upper left', fontsize=8)
fig.savefig('grafico_temp_A.png', dpi=300)
fig.savefig('grafico_temp_A.pdf')
#%% Versión B
#Buen ajuste final:
popt_rho22_balance, pcov_rho22_balance = curve_fit(FinalTemp_withoutC,list(Betas_vec[:10]), [t for t in Temp_vec[:10]],sigma=ErrorTemp_vec[:10],absolute_sigma=False,p0=(-10e6,5e-20)) #esto ajusta muy bien
popt_rho22_balance_noMM, pcov_rho22_balance_noMM = curve_fit(FinalTemp_withoutmicromotion,list(Betas_vec[:10]), [t for t in Temp_vec[:10]],sigma=ErrorTemp_vec[:10],absolute_sigma=False,p0=(-10e6)) #esto ajusta muy bien
betaslong = np.arange(0,2.8,0.01)
print(f'Min temp predicted: {FinalTemp_withoutC(betaslong,*popt_rho22_balance)[0]}')
#print(f'Min temp predicted: {1e3*FinalTemp(betaslong,*popt_rho22_balance)[0]} mK')
print(f'Detuning: {1e-6*popt_rho22_balance[0]} MHz')
print(f'params: {popt_rho22_balance}')
print(f'errores: {np.sqrt(np.diag(pcov_rho22_balance))}')
k_plot = 9
import matplotlib.pyplot as plt
plt.rcParams.update({'font.size': 8})
plt.rcParams['text.usetex']=True
# Acá estimo un detuning de (10 +- 4) MHz
Detuning_range = ( 10 + 4*np.array([-1,1]) )*1e6
Doppler_limit_area = CoolingLimit_noMM(-2*np.pi*Detuning_range)*1e3
figsize=(8.6/2.54,4)
height_ratios=[10,3]
fig, axx = plt.subplots(2,1, sharex=True, constrained_layout=True,
figsize=figsize, gridspec_kw=dict(height_ratios=height_ratios))
fig.set_constrained_layout_pads(w_pad=1/72, h_pad=0, hspace=0, wspace=0)
ax=axx[0]
ax.errorbar(Betas_vec[:k_plot],[t*1e3 for t in Temp_vec[:k_plot]],xerr=ErrorBetas_vec[:k_plot],
yerr=[t*1e3 for t in ErrorTemp_vec[:k_plot]],fmt='o',
capsize=3,markersize=3,color='C3',zorder=10,alpha=0.9, label='inferred data')
#plt.plot(betaslong,[t*1e3 for t in FinalTemp(betaslong,*popt_rho22_balance)],label='Ajuste con espectro modulado')
ax.plot(betaslong,[t*1e3 for t in FinalTemp_withoutC(betaslong,*popt_rho22_balance)],
'-', color='C0', label='model w/ RF heating')
ax.plot(betaslong,[t*1e3 for t in FinalTemp_withoutmicromotion(betaslong,*popt_rho22_balance_noMM)],
'--', color='C1',label='model w/ RF heating')
ax.fill_between( [0, 2.8] , np.ones(2)*Doppler_limit_area[0], np.ones(2)*Doppler_limit_area[1],
color='gray', alpha=0.5, label='Doppler limit w/o micromotion')
ax.grid(True, ls=":", color='lightgray')
# ax.set_yticklabels([ f"{int(y/1000)}k" for y in ax.get_yticks() ])
# ax.set_title(f'({le})', x=0.1, y=0.78, color='gray')
ax.set_ylabel('Temperature [mK]')
# Leyenda ##############################
h1, l1 = ax.get_legend_handles_labels()
ax.legend(h1[-1:]+h1[:-1], l1[-1:]+l1[:-1], loc='upper left', fontsize=8)
ax=axx[1]
residuos = np.array([t*1e3 for t in Temp_vec[:k_plot]]) - np.array([t*1e3 for t in FinalTemp_withoutC(np.array(Betas_vec[:k_plot]),*popt_rho22_balance)])
ax.errorbar(Betas_vec[:k_plot],residuos,xerr=ErrorBetas_vec[:k_plot],
yerr=[t*1e3 for t in ErrorTemp_vec[:k_plot]],fmt='o',
capsize=3,markersize=3,color='C3',zorder=10,alpha=0.9, label='inferred data')
modelA = np.array([t*1e3 for t in FinalTemp_withoutC(betaslong,*popt_rho22_balance)])
modelB = np.array([t*1e3 for t in FinalTemp_withoutmicromotion(betaslong,*popt_rho22_balance_noMM)])
ax.plot(betaslong,modelB-modelA,
'--', color='C1',label='model w/ RF heating')
ax.axhline( 0 , color='C0')
ax.set_xlim(0,2.8)
ax.set_xlabel(r'Modulation factor $\beta$')
ax.set_ylabel('Res. [mK]')
ax.grid(True, ls=":", color='lightgray')
ax.set_ylim(-5,5)
# ax.legend(fontsize=8)
fig.align_ylabels(axx)
fig.savefig('grafico_temp_B.png', dpi=300)
fig.savefig('grafico_temp_B.pdf')
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment