agrego cosas pendientes soy nico

analisis/plots/20250313_fluorescencecorrelationtests/Transitorios_osci.py

@@ -152,7 +152,7 @@ T0 = -0.4e-6
 
 #files = [20577,20578,20579,20580,20581] #el 19613 es el que muestra el pico en 48.5 Hz
 
-files = [20586] #el 19613 es el que muestra el pico en 48.5 Hz
+files = [20585] #el 19613 es el que muestra el pico en 48.5 Hz
 
 def bin_time_arrivals(arrival_times, tau):
     """
@@ -252,16 +252,23 @@ y_extended = np.tile(hist, 3)
 def sinusoidal(x, A, C, D):
     return A * np.sin(2*np.pi*freqrf * x + C) + D
 
-def lorentzian(x,A,b,C,det,fase,cte):
-    det = -10e6*2*np.pi
-    return A*((C**2)/(C**2 + (det-b*2*np.pi*freqrf*np.sin(2*np.pi*freqrf*x+fase))**2))+cte
+def lorentzian(x,A,beta,C,det,fase,cte):
+    #det = -10e6*2*np.pi
+    return A*((C**2)/(C**2 + (det-beta*2*np.pi*freqrf*np.sin(2*np.pi*freqrf*x+fase))**2))+cte
+
+#def lorentzian(x,A,beta,det,fase,cte):
+#    #det = -10e6*2*np.pi
+#    return A/(1 + (det-beta*2*np.pi*freqrf*np.sin(2*np.pi*freqrf*x+fase))**2)+cte
+
 
 x_extended_dense = np.arange(np.min(x_extended),np.max(x_extended),0.1*(x_extended[1]-x_extended[0]))
 
-#params, _ = curve_fit(sinusoidal, x_extended, y_extended, p0=[max(y_extended),  0, np.mean(y_extended)])
-params, _ = curve_fit(lorentzian, x_extended, y_extended, p0=[1,1,1e7,1e7,1,1])
-amps.append(np.abs(params[0]))
+#params, cov = curve_fit(sinusoidal, x_extended, y_extended, p0=[max(y_extended),  0, np.mean(y_extended)])
+params, cov_matrx = curve_fit(lorentzian, x_extended, y_extended, p0=[200,15,-200,-1,80, 100],bounds = ([0,0,-2*np.pi*200e6,-2*np.pi,0,0],[1000,20,0,2*np.pi,200, 1000]))
 
+amps.append(np.abs(params[0]))
+#params = [ 2.0005516e+02,  15.72948375e+00,  2.34874470e+08, -7.26868536e+08,
+#        -1.5,  80]
 
 # Graficar los datos y el ajuste
 plt.figure(figsize=(8, 4))
@@ -274,7 +281,7 @@ for i in range(3):
 
 plt.xlabel("Tiempo bineado (repetido 3 veces)")
 plt.ylabel("Frecuencia")
-plt.title("Distribución de tiempos bineados con ajuste sinusoidal")
+plt.title("Distribución de tiempos bineados con ajuste completo")
 #plt.ylim(0,np.max(hist)*1.2)
 plt.ylim(0,400)
 plt.legend()
@@ -284,3 +291,17 @@ print(f'Amplitud normalizada: {np.abs(round(100*params[0]/params[2],3))}%')
 
 print(f'Beta: {np.abs(params[1])}')
 
+print(f'Amplitud: {params[0]}')
+
+# Step 1: Calculate the standard deviations (square root of diagonal elements)
+std_devs = np.sqrt(np.diagonal(cov_matrix))
+
+# Step 2: Create a matrix of standard deviations (outer product of std_devs)
+std_dev_matrix = np.outer(std_devs, std_devs)
+
+# Step 3: Calculate the correlation matrix by dividing the covariance matrix by the standard deviation matrix
+correlation_matrix = cov_matrix / std_dev_matrix
+
+# Print the resulting correlation matrix
+print("Correlation Matrix:")
+print(correlation_matrix)
\ No newline at end of file

analisis/plots/20250408_FluorescenceCorrelationsTests/Transitorios_osci_20250408.py

+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
+
+#os.chdir('/home/nico/Documents/artiq_experiments/analisis/plots/20221006_transitoriosv2')
+
+os.chdir('/home/nico/Documents/artiq_experiments/analisis/plots/20250408_FluorescenceCorrelationsTests/Data1')
+
+def expo(T, tau, N0, C):
+    
+    global T0
+    return N0*np.exp(-(T-T0)/tau) + C
+
+def pow_from_amp(amp):
+    """Paso de amplitud urukul a potencia medida por Nico"""
+    # Forma altamente ineficiente de hacer esto, pero me salio asi
+    amplitudes_UV = np.flip(np.array([0.08, 0.10, 0.12, 0.14, 0.16, 0.18, 0.20, 0.22, 0.24, 0.26, 0.28, 0.30]))
+    assert amp in amplitudes_UV
+    potencias_UV = np.flip(np.array([4, 10, 19, 32, 49, 71, 96, 125, 155, 183, 208, 229]))
+    return potencias_UV[np.where(amplitudes_UV == amp)][0]
+
+
+def SP_Bkgr_builder(amp_in, amp_fin, derivadainicio, derivadafin, longbins):
+    CalibCurve = []
+    j=0
+    while j<longbins:
+        if j<=derivadainicio:
+            CalibCurve.append(amp_in)
+        elif j>=derivadainicio and j<=derivadafin:
+            pendiente=(amp_fin-amp_in)/(derivadafin-derivadainicio)
+            CalibCurve.append(amp_in+pendiente*(j-derivadainicio))
+        else:
+            CalibCurve.append(amp_fin)
+        j=j+1
+    return CalibCurve
+
+    
+
+"""
+plt.plot(amplitudes_UV, potencias_UV, 'ko-', lw=0.2)
+plt.xlabel("Amplitud Urukul")
+plt.ylabel("Potencia /uW")
+plt.grid()
+"""
+#%%
+import scipy.fftpack
+
+BINW = 1e-2
+T0 = -0.4e-6
+
+files = [20689, 20690, 20691]
+
+SP_Heigths = []
+SP_Bins = []
+
+for i, fname in enumerate(files):
+    #print(i)
+    #print(fname)
+    #data = h5py.File('0000'+str(fname)+'-MicromotionCompensation.h5', 'r')
+    data = h5py.File('0000'+str(fname)+'-FluorescenceCorrelation.h5', 'r')
+    # counts = np.array(data['datasets']['countsrf'])
+    counts = np.array(data['datasets']['time_stamps'])
+    bines = np.arange(counts.min(), counts.max()+BINW, BINW)
+
+    heigs, binsf  = np.histogram(counts, bines[bines>T0])
+    
+    SP_Heigths.append(heigs)
+    SP_Bins.append(binsf)
+
+
+#freq: 7.262 mhZ
+for k in range(len(counts)):
+    print(counts[k]-counts[k-1])
+    #print(1e-6/(counts[k]-counts[k-1]))
+
+
+#%%
+
+"""
+Esto binea fotones asumienod que la cantidad de fotones que llegan en un intervalo tau
+esta modulada por una funcion periodica. Entonces a los tiempos de llegada entre
+0 y tau no les hace nada, a los tiempos entre tau y 2tau les resta tau,
+a los tiempos entre 2tau y 3tau les resta tau, y asi. Y despues
+binea y ajusta con una sinusoidal y devuelve la frecuencia del ajuste.
+
+Asi anda barbaro para ver la oscilacion a casi 50 Hz
+"""
+
+def bin_time_arrivals(arrival_times, tau):
+    """
+    Binea los tiempos de llegada de los fotones según el periodo tau.
+    
+    Parameters:
+        arrival_times (numpy array): Vector con los tiempos de llegada de los fotones.
+        tau (float): Periodo de bineado.
+        
+    Returns:
+        numpy array: Vector con los tiempos bienados.
+    """
+    return arrival_times - tau * (arrival_times // tau)
+
+#taurf = 1/(22.135e6)
+taurf = 1/(7262766.1)
+
+tiemposarreglados = bin_time_arrivals(counts, taurf)
+
+b = taurf/100  # Ancho de bineo
+bins = np.arange(0, taurf, b)
+hist, bin_edges = np.histogram(tiemposarreglados, bins=bins)
+
+colors = ['b', 'g', 'r']  # Colores para cada repetición
+
+x_extended = np.concatenate([bin_edges[:-1] + i * taurf for i in range(3)])
+y_extended = np.tile(hist, 3)
+
+def sinusoidal(x, A, B, C, D):
+    return A * np.sin(B * x + C) + D
+
+x_extended_dense = np.arange(np.min(x_extended),np.max(x_extended),0.1*(x_extended[1]-x_extended[0]))
+
+params, _ = curve_fit(sinusoidal, x_extended, y_extended, p0=[max(y_extended), 2*np.pi/taurf, 0, np.mean(y_extended)])
+
+# Graficar los datos y el ajuste
+plt.figure(figsize=(8, 4))
+plt.plot(x_extended, y_extended, marker='o', linestyle='-', label="Datos", color='gray')
+plt.plot(x_extended_dense, sinusoidal(x_extended_dense, *params), linestyle='--', label="Ajuste sinusoidal", color='red')
+for i in range(3):
+    plt.plot(bin_edges[:-1] + i * taurf, hist, marker='o', linestyle='-', color=colors[i])
+
+plt.xlabel("Tiempo bineado (repetido 3 veces)")
+plt.ylabel("Frecuencia")
+plt.title("Distribución de tiempos bineados con ajuste sinusoidal")
+plt.legend()
+plt.show()
+
+print(f'Frecuencia del ajuste: {round(params[1]/(2*np.pi),2)} Hz')
+
+
+#%%
+
+"""
+Esto deberia hacerse asi para ver la oscilacion a la RF (micromocion)
+"""
+
+BINW = 1e-2
+T0 = -0.4e-6
+
+#files = [20577,20578,20579,20580,20581] #el 19613 es el que muestra el pico en 48.5 Hz
+
+files = [20689] 
+
+def bin_time_arrivals(arrival_times, tau):
+    """
+    Binea los tiempos de llegada de los fotones según el periodo tau.
+    
+    Parameters:
+        arrival_times (numpy array): Vector con los tiempos de llegada de los fotones.
+        tau (float): Periodo de bineado.
+        
+    Returns:
+        numpy array: Vector con los tiempos bienados.
+    """
+    return arrival_times - tau * (arrival_times // tau)
+
+
+SP_Heigths = []
+SP_Bins = []
+
+
+plt.figure()
+
+for i, fname in enumerate(files):
+    #print(i)
+    #print(fname)
+    data = h5py.File('0000'+str(fname)+'-FluorescenceCorrelation.h5', 'r')
+    counts = np.array(data['datasets']['time_stamps'])
+    bines = np.arange(counts.min(), counts.max()+BINW, BINW)
+
+    heigs, binsf  = np.histogram(counts, bines[bines>T0])
+    
+    SP_Heigths.append(heigs)
+    SP_Bins.append(binsf)
+    
+    
+
+    
+    
+    """
+    Amplitud en funcion de la frecuencia de bineo de la rf
+    """
+    
+    freqrf = 7.2618e6
+    
+    freqrfvec = np.arange(7262765,7262767,0.01)
+    
+    amps = []
+    
+    
+    for freqrf in freqrfvec:
+        taurf = 1/freqrf
+        
+        tiemposarreglados = bin_time_arrivals(counts, taurf)
+        
+        b = taurf/50  # Ancho de bineo
+        bins = np.arange(0, taurf, b)
+        hist, bin_edges = np.histogram(tiemposarreglados, bins=bins)
+        
+        colors = ['b', 'g', 'r']  # Colores para cada repetición
+        
+        x_extended = np.concatenate([bin_edges[:-1] + i * taurf for i in range(3)])
+        y_extended = np.tile(hist, 3)
+        
+        def sinusoidal(x, A, C, D):
+            return A * np.sin(2*np.pi*freqrf * x + C) + D
+        
+        x_extended_dense = np.arange(np.min(x_extended),np.max(x_extended),0.1*(x_extended[1]-x_extended[0]))
+        
+        params, _ = curve_fit(sinusoidal, x_extended, y_extended, p0=[max(y_extended),  0, np.mean(y_extended)])
+        amps.append(np.abs(params[0]))
+
+
+    plt.plot([f*1e-6 for f in freqrfvec],amps,'o')
+plt.xlabel('Frecuencia RF (MHz)')
+plt.ylabel('Amplitud oscilacion')
+
+
+"""
+Grafico temporal replicado 3 veces a la frecuencia dada
+"""
+
+
+freqrf = 7262766.1
+
+taurf = 1/freqrf
+
+tiemposarreglados = bin_time_arrivals(counts, taurf)
+
+b = taurf/50  # Ancho de bineo
+bins = np.arange(0, taurf, b)
+hist, bin_edges = np.histogram(tiemposarreglados, bins=bins)
+
+colors = ['b', 'g', 'r']  # Colores para cada repetición
+
+x_extended = np.concatenate([bin_edges[:-1] + i * taurf for i in range(3)])
+y_extended = np.tile(hist, 3)
+
+def sinusoidal(x, A, C, D):
+    return A * np.sin(2*np.pi*freqrf * x + C) + D
+
+def lorentzian(x,A,b,C,det,fase,cte):
+    det = -10e6*2*np.pi
+    return A*((C**2)/(C**2 + (det-b*2*np.pi*freqrf*np.sin(2*np.pi*freqrf*x+fase))**2))+cte
+
+x_extended_dense = np.arange(np.min(x_extended),np.max(x_extended),0.1*(x_extended[1]-x_extended[0]))
+
+#params, _ = curve_fit(sinusoidal, x_extended, y_extended, p0=[max(y_extended),  0, np.mean(y_extended)])
+params, _ = curve_fit(lorentzian, x_extended, y_extended, p0=[1,1e-3,1e7,1e4,1,1])
+amps.append(np.abs(params[0]))
+
+
+# Graficar los datos y el ajuste
+plt.figure(figsize=(8, 4))
+
+#plt.plot(x_extended_dense, sinusoidal(x_extended_dense, *params), linestyle='--', label="Ajuste sinusoidal", color='black',linewidth=3,zorder=2)
+plt.plot(x_extended_dense, lorentzian(x_extended_dense, *params), linestyle='--', label="Ajuste completo", color='black',linewidth=3,zorder=2)
+
+for i in range(3):
+    plt.plot(bin_edges[:-1] + i * taurf, hist, marker='o', linestyle='-', color=colors[i],markersize=3,zorder=1)
+
+plt.xlabel("Tiempo bineado (repetido 3 veces)")
+plt.ylabel("Frecuencia")
+plt.title("Distribución de tiempos bineados con ajuste completo")
+#plt.ylim(0,np.max(hist)*1.2)
+plt.ylim(0,400)
+plt.legend()
+plt.show()
+
+print(f'Amplitud normalizada: {np.abs(round(100*params[0]/params[2],3))}%')
+
+print(f'Beta: {np.abs(params[1])}')
+
+print(f'Amplitud: {params[0]}')

analisis/plots/20250421_fluorescencorrelation/check_rffrequency.py

+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
+
+#os.chdir('/home/nico/Documents/artiq_experiments/analisis/plots/20221006_transitoriosv2')
+
+# os.chdir('/home/nico/Documents/artiq_experiments/analisis/plots/20250408_FluorescenceCorrelationsTests/Data1')
+
+def expo(T, tau, N0, C):
+    
+    global T0
+    return N0*np.exp(-(T-T0)/tau) + C
+
+def pow_from_amp(amp):
+    """Paso de amplitud urukul a potencia medida por Nico"""
+    # Forma altamente ineficiente de hacer esto, pero me salio asi
+    amplitudes_UV = np.flip(np.array([0.08, 0.10, 0.12, 0.14, 0.16, 0.18, 0.20, 0.22, 0.24, 0.26, 0.28, 0.30]))
+    assert amp in amplitudes_UV
+    potencias_UV = np.flip(np.array([4, 10, 19, 32, 49, 71, 96, 125, 155, 183, 208, 229]))
+    return potencias_UV[np.where(amplitudes_UV == amp)][0]
+
+
+def SP_Bkgr_builder(amp_in, amp_fin, derivadainicio, derivadafin, longbins):
+    CalibCurve = []
+    j=0
+    while j<longbins:
+        if j<=derivadainicio:
+            CalibCurve.append(amp_in)
+        elif j>=derivadainicio and j<=derivadafin:
+            pendiente=(amp_fin-amp_in)/(derivadafin-derivadainicio)
+            CalibCurve.append(amp_in+pendiente*(j-derivadainicio))
+        else:
+            CalibCurve.append(amp_fin)
+        j=j+1
+    return CalibCurve
+
+    
+
+"""
+plt.plot(amplitudes_UV, potencias_UV, 'ko-', lw=0.2)
+plt.xlabel("Amplitud Urukul")
+plt.ylabel("Potencia /uW")
+plt.grid()
+"""
+#%%
+import scipy.fftpack
+
+BINW = 1e-2
+T0 = -0.4e-6
+
+files = [21090]
+
+SP_Heigths = []
+SP_Bins = []
+
+for i, fname in enumerate(files):
+    #print(i)
+    #print(fname)
+    #data = h5py.File('0000'+str(fname)+'-MicromotionCompensation.h5', 'r')
+    data = h5py.File('0000'+str(fname)+'-FluorescenceCorrelation.h5', 'r')
+    # counts = np.array(data['datasets']['countsrf'])
+    counts = np.array(data['datasets']['time_stamps'])
+    bines = np.arange(counts.min(), counts.max()+BINW, BINW)
+
+    heigs, binsf  = np.histogram(counts, bines[bines>T0])
+    
+    SP_Heigths.append(heigs)
+    SP_Bins.append(binsf)
+
+
+#freq: 7.262 mhZ
+for k in range(len(counts)):
+    print(counts[k]-counts[k-1])
+    #print(1e-6/(counts[k]-counts[k-1]))
+
+
+
+#%%
+
+"""
+Esto deberia hacerse asi para ver la oscilacion a la RF (micromocion)
+"""
+
+BINW = 1e-3
+T0 = -0.4e-6
+
+#files = [20577,20578,20579,20580,20581] #el 19613 es el que muestra el pico en 48.5 Hz
+
+files = [21090] 
+
+def bin_time_arrivals(arrival_times, tau):
+    """
+    Binea los tiempos de llegada de los fotones según el periodo tau.
+    
+    Parameters:
+        arrival_times (numpy array): Vector con los tiempos de llegada de los fotones.
+        tau (float): Periodo de bineado.
+        
+    Returns:
+        numpy array: Vector con los tiempos bienados.
+    """
+    return arrival_times - tau * (arrival_times // tau)
+
+
+SP_Heigths = []
+SP_Bins = []
+
+
+plt.figure()
+
+for i, fname in enumerate(files):
+    #print(i)
+    #print(fname)
+    data = h5py.File('0000'+str(fname)+'-FluorescenceCorrelation.h5', 'r')
+    counts = np.array(data['datasets']['time_stamps'])
+    bines = np.arange(counts.min(), counts.max()+BINW, BINW)
+
+    heigs, binsf  = np.histogram(counts, bines[bines>T0])
+    
+    SP_Heigths.append(heigs)
+    SP_Bins.append(binsf)
+    
+    
+
+    
+    
+    """
+    Amplitud en funcion de la frecuencia de bineo de la rf
+    """
+    
+    freqrf = 7.2618e6
+    
+    freqrfvec = np.arange(7262765,7262767,0.01)
+    
+    amps = []
+    
+    
+    for freqrf in freqrfvec:
+        taurf = 1/freqrf
+        
+        tiemposarreglados = bin_time_arrivals(counts, taurf)
+        
+        b = taurf/500  # Ancho de bineo
+        bins = np.arange(0, taurf, b)
+        hist, bin_edges = np.histogram(tiemposarreglados, bins=bins)
+        
+        colors = ['b', 'g', 'r']  # Colores para cada repetición
+        
+        x_extended = np.concatenate([bin_edges[:-1] + i * taurf for i in range(3)])
+        y_extended = np.tile(hist, 3)
+        
+        def sinusoidal(x, A, C, D):
+            return A * np.sin(2*np.pi*freqrf * x + C) + D
+        
+        x_extended_dense = np.arange(np.min(x_extended),np.max(x_extended),0.1*(x_extended[1]-x_extended[0]))
+        
+        params, _ = curve_fit(sinusoidal, x_extended, y_extended, p0=[max(y_extended),  0, np.mean(y_extended)])
+        amps.append(np.abs(params[0]))
+
+
+    plt.plot([f*1e-6 for f in freqrfvec],amps,'o')
+plt.xlabel('Frecuencia RF (MHz)')
+plt.ylabel('Amplitud oscilacion')
+plt.axvline(7262766.13*1e-6,linewidth=4,color='red',linestyle='dashed')
+#%%
+
+"""
+Grafico temporal replicado 3 veces a la frecuencia dada
+"""
+
+
+freqrf = 7262766.13
+
+taurf = 1/freqrf
+
+tiemposarreglados = bin_time_arrivals(counts, taurf)
+
+b = taurf/50  # Ancho de bineo
+bins = np.arange(0, taurf, b)
+hist, bin_edges = np.histogram(tiemposarreglados, bins=bins)
+
+colors = ['b', 'g', 'r']  # Colores para cada repetición
+
+x_extended = np.concatenate([bin_edges[:-1] + i * taurf for i in range(3)])
+y_extended = np.tile(hist, 3)
+
+def sinusoidal(x, A, C, D):
+    return A * np.sin(2*np.pi*freqrf * x + C) + D
+
+def lorentzian(x,A,b,C,det,fase,cte):
+    det = -10e6*2*np.pi
+    return A*((C**2)/(C**2 + (det-b*2*np.pi*freqrf*np.sin(2*np.pi*freqrf*x+fase))**2))+cte
+
+x_extended_dense = np.arange(np.min(x_extended),np.max(x_extended),0.1*(x_extended[1]-x_extended[0]))
+
+#params, _ = curve_fit(sinusoidal, x_extended, y_extended, p0=[max(y_extended),  0, np.mean(y_extended)])
+params, _ = curve_fit(lorentzian, x_extended, y_extended, p0=[1,1e-3,1e7,1e4,1,1])
+amps.append(np.abs(params[0]))
+
+
+# Graficar los datos y el ajuste
+plt.figure(figsize=(8, 4))
+
+#plt.plot(x_extended_dense, sinusoidal(x_extended_dense, *params), linestyle='--', label="Ajuste sinusoidal", color='black',linewidth=3,zorder=2)
+plt.plot(x_extended_dense, lorentzian(x_extended_dense, *params), linestyle='--', label="Ajuste completo", color='black',linewidth=3,zorder=2)
+
+for i in range(3):
+    plt.plot(bin_edges[:-1] + i * taurf, hist, marker='o', linestyle='-', color=colors[i],markersize=3,zorder=1)
+
+plt.xlabel("Tiempo bineado (repetido 3 veces)")
+plt.ylabel("Frecuencia")
+plt.title("Distribución de tiempos bineados con ajuste completo")
+#plt.ylim(0,np.max(hist)*1.2)
+plt.ylim(0,400)
+plt.legend()
+plt.show()
+
+print(f'Amplitud normalizada: {np.abs(round(100*params[0]/params[2],3))}%')
+
+print(f'Beta: {np.abs(params[1])}')
+
+print(f'Amplitud: {params[0]}')