update graficos

.gitignore

 **/__pycache__/
 *.pyc
 *.pyon
+**/__*.py
+**/__*.md
+analisis/plots/20231123_CPTconmicromocion3/grafico_*.png
+analisis/plots/20231123_CPTconmicromocion3/grafico_*.pdf

analisis/plots/20231123_CPTconmicromocion3/lolo_graficos_pub.py

@@ -710,3 +710,348 @@ 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/ 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')
+