scripts marce

analisis/plots/20231123_CPTconmicromocion3/Data/EITfit/lolo_modelo_full_8niveles.py

+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+"""
+20 de dic 2023
+
+@author: lolo
+
+reingenieria del código que anda
+
+MAPA de FUNCIONES
+
+CPTspectrum8levels_MM
+ |--> FullL_MM => ndarray(64,64,np.complex_)
+       |--> dopplerBroadening => float
+       |--> EffectiveL => ndarray(8,8,np.complex_)
+       |--> H0matrix => ndarray(8,8,np.complex_)
+       |--> HImatrix => ndarray(8,8,np.complex_)
+       |--> CalculateSingleMmatrix => ndarray(64,64,np.complex_)
+       |--> LtempCalculus => ndarray,ndarray
+       |--> GetL1 => ndarray
+
+"""
+
+# pylint: disable=C0301,R0913,R0914,W0621
+
+import time
+import random
+import numpy as np
+import matplotlib.pyplot as plt
+from scipy.signal import savgol_filter as sf
+
+from  numba import jit,njit
+
+@njit
+def PerformExperiment_8levels_MM(sg, sp, gPS, gPD, DetDoppler, u, DopplerLaserLinewidth, ProbeLaserLinewidth, T, alpha, phidoppler, titadoppler, phiprobe, titaprobe, beta, drivefreq, freqMin, freqMax, freqStep, circularityprobe=1, plot=False, solvemode=1, detpvec=None):
+    """
+    solvemode=1: resuelve con np.linalg.solve
+    solvemode=2: resuelve invirtiendo L con la funcion np.linalg.inv
+    """
+    rta = CPTspectrum8levels_MM(sg, sp, gPS, gPD, DetDoppler, u,
+                                DopplerLaserLinewidth, ProbeLaserLinewidth,
+                                T, alpha, phidoppler, titadoppler, phiprobe,
+                                titaprobe, circularityprobe, beta, drivefreq,
+                                freqMin=freqMin, freqMax=freqMax, freqStep=freqStep,
+                                plot=False, solvemode=1)
+    ProbeDetuningVectorL, Fluovector = rta
+    return ProbeDetuningVectorL, Fluovector
+
+
+def GenerateNoisyCPT_MM(sg, sp, gPS, gPD, DetDoppler, u, DopplerLaserLinewidth, ProbeLaserLinewidth, T, alpha, phidoppler, titadoppler, phiprobe, titaprobeVec, kg, kr, v0, drivefreq, freqMin, freqMax, freqStep, circularityprobe=1, plot=False, solvemode=1, detpvec=None, noiseamplitude=0.001):
+    """
+    Genera un resultado de PerformExperiment_8levels_MM con ruido normal agregado
+    """
+    nFrequencyvector, Fluovector = PerformExperiment_8levels_MM(sg, sp, gPS, gPD, DetDoppler, u, DopplerLaserLinewidth, ProbeLaserLinewidth, T, alpha, phidoppler, titadoppler, phiprobe, titaprobeVec, kg, kr, v0, drivefreq, freqMin, freqMax, freqStep, circularityprobe, plot=False, solvemode=1, detpvec=None)
+    NoisyFluovector = [fluo+noiseamplitude*(2*random.random()-1) for fluo in Fluovector]
+    return nFrequencyvector, NoisyFluovector
+
+def GenerateNoisyCPT_MM_fit(sg, sp, gPS, gPD, DetDoppler, u, DopplerLaserLinewidth, ProbeLaserLinewidth, T, alpha, phidoppler, titadoppler, phiprobe, titaprobeVec, beta, drivefreq, freqs, circularityprobe=1, plot=False, solvemode=1, detpvec=None, noiseamplitude=0.001):
+    """
+    Este no se qué hace
+    """
+    Frequencyvector, Fluovector = PerformExperiment_8levels_MM(sg, sp, gPS, gPD, DetDoppler, u, DopplerLaserLinewidth, ProbeLaserLinewidth, T, alpha, phidoppler, titadoppler, phiprobe, titaprobeVec, beta, drivefreq, freqs[0], freqs[-1], freqs[1]-freqs[0], circularityprobe, plot=False, solvemode=1, detpvec=None)
+    #NoisyFluovector = [fluo+noiseamplitude*(2*random.random()-1) for fluo in Fluovector]
+    return Frequencyvector, Fluovector
+
+
+def SmoothNoisyCPT(Fluo, window=11, poly=3):
+    SmoothenFluo = sf(Fluo, window, poly)
+    return SmoothenFluo
+
+
+
+#%% Estos son los auxiliares ###################################################
+
+"""
+Esta parte es la del modelo
+"""
+
+@njit
+def make_diag(vec):
+    "Construye matris diagonal desde una lista o vector"
+    return np.eye(len(vec))*np.array(vec).reshape(-1,1)
+
+# @njit
+# def kron(a,b):
+#     "Hago el producto de Kronecker a mano"
+#     return np.vstack( [  np.hstack( [ a[k,j]*b for j in range(a.shape[1]) ] )  for k in range(a.shape[0])])
+
+# @njit
+# def kron(A, B):
+#     cola = A.shape[1]
+#     rowa = A.shape[0]
+#     colb = B.shape[1]
+#     rowb = B.shape[0]
+#
+#     C = [[0] * (cola * colb)  for _ in range(rowa * rowb) ]
+#
+#     for i in range(rowa):
+#         for k in range(cola):
+#             for j in range(rowb):
+#                 for l in range(colb):
+#                     C[i * rowb + k][j * colb + l]  = A[i][j] * B[k][l]
+#     return np.array(C)
+import numba
+
+@jit
+def kron(A,B):
+    out=np.empty((A.shape[0],B.shape[0],A.shape[1],B.shape[1]),dtype=A.dtype)
+    for i in numba.prange(A.shape[0]):
+        for j in range(B.shape[0]):
+            for k in range(A.shape[1]):
+                for l in range(B.shape[1]):
+                    out[i,j,k,l]=A[i,k]*B[j,l]
+    return out
+
+
+@njit
+def H0matrix(Detg, Detp, u):
+    """
+    Calcula la matriz H0 en donde dr es el detuning del doppler, dp es el retuning
+    del repump y u es el campo magnético en Hz/Gauss.
+    Para esto se toma la energía del nivel P como 0
+    """
+
+    eigenEnergies = (Detg-u, Detg+u, -u/3, u/3, Detp-6*u/5,
+                    Detp-2*u/5, Detp+2*u/5, Detp+6*u/5)
+    #pagina 26 de Oberst. los lande del calcio son iguales a Bario.
+    # H0 = np.diag(eigenEnergies)
+    # H0 = np.eye(len(eigenEnergies))*np.array(eigenEnergies).reshape(-1,1)
+    H0 = make_diag(eigenEnergies)
+    return H0
+
+@njit
+def HImatrix(rabG, rabP, phidoppler, titadoppler, phiprobe, titaprobe, circularityprobe=1):
+    """
+    Calcula la matriz de interacción Hsp + Hpd, en donde rabR es la frecuencia de rabi de
+    la transición Doppler SP, rabP es la frecuencia de rabi de la transición repump DP,
+    y las componentes ei_r y ei_p son las componentes de la polarización
+    del campo eléctrico incidente de doppler y repump respectivamente. Deben
+    estar normalizadas a 1
+    """
+    HI = np.zeros((8, 8), dtype=np.complex_)
+
+
+    i, j = 1, 3
+    HI[i-1, j-1] = (rabG/np.sqrt(3)) * np.cos(titadoppler)
+    HI[j-1, i-1] = np.conjugate(HI[i-1, j-1])
+
+    i, j = 1, 4
+    HI[i-1, j-1] = -(rabG/np.sqrt(3)) * np.sin(titadoppler)*np.exp(1j*phidoppler)
+    HI[j-1, i-1] = np.conjugate(HI[i-1, j-1])
+
+    i, j = 2, 3
+    HI[i-1, j-1] = -(rabG/np.sqrt(3)) * np.sin(titadoppler)*np.exp(-1j*phidoppler)
+    HI[j-1, i-1] = np.conjugate(HI[i-1, j-1])
+
+    i, j = 2, 4
+    HI[i-1, j-1] = -(rabG/np.sqrt(3)) * np.cos(titadoppler)
+    HI[j-1, i-1] = np.conjugate(HI[i-1, j-1])
+
+    i, j = 3, 5
+    HI[i-1, j-1] = -(rabP/2) * np.sin(titaprobe)*(np.cos(phiprobe)-1j*np.sin(phiprobe)*circularityprobe)
+    HI[j-1, i-1] = np.conjugate(HI[i-1, j-1])
+
+    i, j = 3, 6
+    HI[i-1, j-1] = -(rabP/np.sqrt(3)) * np.cos(titaprobe)
+    HI[j-1, i-1] = np.conjugate(HI[i-1, j-1])
+
+    i, j = 3, 7
+    HI[i-1, j-1] = rabP/np.sqrt(12) * np.sin(titaprobe)*(np.cos(phiprobe)+1j*np.sin(phiprobe)*circularityprobe)
+    HI[j-1, i-1] = np.conjugate(HI[i-1, j-1])
+
+
+    i, j = 4, 6
+    HI[i-1, j-1] = -(rabP/np.sqrt(12)) * np.sin(titaprobe)*(np.cos(phiprobe)-1j*np.sin(phiprobe)*circularityprobe)
+    HI[j-1, i-1] = np.conjugate(HI[i-1, j-1])
+
+    i, j = 4, 7
+    HI[i-1, j-1] = -(rabP/np.sqrt(3)) * np.cos(titaprobe)
+    HI[j-1, i-1] = np.conjugate(HI[i-1, j-1])
+
+    i, j = 4, 8
+    HI[i-1, j-1] = (rabP/2) * np.sin(titaprobe)*(np.cos(phiprobe)+1j*np.sin(phiprobe)*circularityprobe)
+    HI[j-1, i-1] = np.conjugate(HI[i-1, j-1])
+
+
+    return HI
+
+# @jit
+# def LtempCalculus(beta:float, drivefreq:float, forma=1):
+#     Hint = np.zeros((8, 8), dtype=np.complex_)
+
+#     ampg=beta*drivefreq
+#     ampr=beta*drivefreq*(397/866)
+#     Hint[0,0] = ampg
+#     Hint[1,1] = ampg
+#     Hint[4,4] = ampr
+#     Hint[5,5] = ampr
+#     Hint[6,6] = ampr
+#     Hint[7,7] = ampr
+
+#     if forma==1:
+#         Ltemp = np.zeros((64, 64), dtype=np.complex_)
+
+#         for r in range(8):
+#             for q in range(8):
+#                 if r!=q:
+#                     Ltemp[r*8+q][r*8+q] = (-1j)*(Hint[r,r] - Hint[q,q])
+
+#     if forma==2:
+#         # deltaKro = np.diag([1, 1, 1, 1, 1, 1, 1, 1])
+#         deltaKro = make_diag([1., 1., 1., 1., 1., 1., 1., 1.]).astype(np.complex_)
+#         # Ltemp = (-1j)*(np.kron(Hint, deltaKro) - np.kron(deltaKro, Hint))
+#         Ltemp = (-1j)*(kron(Hint, deltaKro) - kron(deltaKro, Hint))
+
+#     Omega = np.zeros((64, 64), dtype=np.complex_)
+
+#     for i in range(64):
+#         Omega[i, i] = (1j)*drivefreq
+
+#     return Ltemp, Omega
+
+@njit
+def LtempCalculus(beta:float, drivefreq:float, forma=1):
+    Hint = np.zeros((8, 8), dtype=np.complex_)
+
+    ampg=beta*drivefreq
+    ampr=beta*drivefreq*(397/866)
+    Hint[0,0] = ampg
+    Hint[1,1] = ampg
+    Hint[4,4] = ampr
+    Hint[5,5] = ampr
+    Hint[6,6] = ampr
+    Hint[7,7] = ampr
+
+    Ltemp = np.zeros((64, 64), dtype=np.complex_)
+
+    for r in range(8):
+        for q in range(8):
+            if r!=q:
+                Ltemp[r*8+q][r*8+q] = (-1j)*(Hint[r,r] - Hint[q,q])
+
+    Omega = np.zeros((64, 64), dtype=np.complex_)
+
+    for i in range(64):
+        Omega[i, i] = (1j)*drivefreq
+
+    return Ltemp, Omega
+
+# LtempCalculus(0,1)
+# raise ValueError('aaa')
+
+
+
+
+
+
+@njit
+def GetL1(Ltemp, L0, Omega, nmax):
+    """
+    Devuelve Splus0 y Sminus0
+    """
+
+    # Sp = (-1)*(np.matrix(np.linalg.inv(L0 - (nmax+1)*Omega))*0.5*np.matrix(Ltemp))
+    # Sm = (-1)*(np.matrix(np.linalg.inv(L0 + (nmax+1)*Omega))*0.5*np.matrix(Ltemp))
+
+    Sp = (-1)*np.linalg.inv(L0 - (nmax+1)*Omega).dot(0.5*Ltemp)
+    Sm = (-1)*np.linalg.inv(L0 + (nmax+1)*Omega).dot(0.5*Ltemp)
+
+    for n in list(range(nmax+1))[(nmax+1)::-1][0:len(list(range(nmax+1))[(nmax+1)::-1])-1]: #jaja esto solo es para que vaya de nmax a 1 bajando. debe haber algo mas facil pero kcio
+        # Sp = (-1)*(np.matrix(np.linalg.inv(L0 - n*Omega + (0.5*Ltemp*np.matrix(Sp))))*0.5*np.matrix(Ltemp))
+        # Sm = (-1)*(np.matrix(np.linalg.inv(L0 + n*Omega + (0.5*Ltemp*np.matrix(Sm))))*0.5*np.matrix(Ltemp))
+
+        Sp = (-1)*np.linalg.inv(L0 - n*Omega + (0.5*Ltemp.dot(Sp))).dot(0.5*Ltemp)
+        Sm = (-1)*np.linalg.inv(L0 + n*Omega + (0.5*Ltemp.dot(Sm))).dot(0.5*Ltemp)
+
+    # L1 = 0.5*np.matrix(Ltemp)*(np.matrix(Sp) + np.matrix(Sm))
+    L1 = 0.5*Ltemp.dot(Sp + Sm)
+
+    return L1
+
+@njit
+def EffectiveL(gPS, gPD, lwg, lwp):
+
+    """
+    Siendo Heff = H + EffectiveL, calcula dicho EffectiveL que es (-0.5j)*sumatoria(CmDaga*Cm) que luego sirve para calcular el Liouvilliano
+    """
+    Leff = np.zeros((8, 8), dtype=np.complex_)
+
+    Leff[0, 0] = 2*lwg
+    Leff[1, 1] = 2*lwg
+    Leff[2, 2] = ((2/3)+(1/3))*gPS + ((1/2) + (1/6) + (1/3))*gPD
+    Leff[3, 3] = ((2/3)+(1/3))*gPS + ((1/2) + (1/6) + (1/3))*gPD
+    Leff[4, 4] = 2*lwp
+    Leff[5, 5] = 2*lwp
+    Leff[6, 6] = 2*lwp
+    Leff[7, 7] = 2*lwp
+
+    return (-0.5j)*Leff
+
+@njit
+def CalculateSingleMmatrix(gPS, gPD, lwg, lwp):
+    """
+    Si tomamos el Liuvilliano como L = (-j)*(Heff*deltak - Heffdaga*deltak) + sum(Mm),
+    esta funcion calcula dichos Mm, que tienen dimensión 64x64 ya que esa es la dimensión del L. Estas componentes
+    salen de hacer la cuenta a mano conociendo los Cm y considerando que Mm[8*(r-1)+s, 8*(k-1)+j] = Cm[r,l] + Cmdaga[j,s] = Cm[r,l] + Cm[s,j]
+    ya que los componentes de Cm son reales.
+    Esta M es la suma de las 8 matrices M.
+
+    """
+
+    M = np.zeros((64, 64), dtype=np.complex_)
+
+    M[0,27] = (2/3)*gPS
+
+    M[9,18] = (2/3)*gPS
+
+    M[0,18] = (1/3)*gPS
+    M[1,19] = -(1/3)*gPS
+    M[8,26] = -(1/3)*gPS
+    M[9,27] = (1/3)*gPS
+
+    M[36,18] = (1/2)*gPD
+    M[37,19] = (1/np.sqrt(12))*gPD
+    M[44,26] = (1/np.sqrt(12))*gPD
+    M[45,27] = (1/6)*gPD
+
+    M[54,18] = (1/6)*gPD
+    M[55,19] = (1/np.sqrt(12))*gPD
+    M[62,26] = (1/np.sqrt(12))*gPD
+    M[63,27] = (1/2)*gPD
+
+    M[45,18] = (1/3)*gPD
+    M[46,19] = (1/3)*gPD
+    M[53,26] = (1/3)*gPD
+    M[54,27] = (1/3)*gPD
+
+    M[0,0] = 2*lwg
+    M[1,1] = 2*lwg
+    M[8,8] = 2*lwg
+    M[9,9] = 2*lwg
+
+    for k in [36, 37, 38, 39, 44, 45, 46, 47, 52, 53, 54, 55, 60, 61, 62, 63]:
+        M[k,k]=2*lwp
+
+    return M
+
+@njit
+def dopplerBroadening(wlg, wlp, alpha, T, mcalcio = 6.655e-23*1e-3):
+    """
+    Calcula el broadening extra semiclásico por temperatura considerando que el ion atrapado se mueve.
+    wlg es la longitud de onda doppler, wlp la longitud de onda repump, T la temperatura del ion en kelvin, y alpha (en rads) el ángulo
+    que forman ambos láseres.
+    """
+
+    kboltzmann = 1.38e-23 #J/K
+
+    gammaD = (2*np.pi)*np.sqrt((1/(wlg*wlg)) + (1/(wlp*wlp)) - 2*(1/(wlg*wlp))*np.cos(alpha))*np.sqrt(kboltzmann*T/(2*mcalcio))
+
+    return gammaD
+
+@njit
+def FullL_MM(rabG, rabP, gPS = 0, gPD = 0, Detg = 0, Detp = 0, u = 0, lwg = 0, lwp = 0,
+          phidoppler=0, titadoppler=0, phiprobe=0, titaprobe=0, beta=0,
+          drivefreq=2*np.pi*22.135*1e6, T = 0, alpha = 0, circularityprobe=1):
+
+    """
+    Calcula el Liouvilliano total de manera explícita índice a índice. Suma aparte las componentes de las matrices M.
+    Es la más eficiente hasta ahora.
+    """
+
+
+    db = dopplerBroadening(0.397e-6, 0.866e-6, alpha, T)
+
+    lwg = np.sqrt(lwg**2 + db**2)
+    lwp = np.sqrt(lwp**2 + db**2)
+
+    CC = EffectiveL(gPS, gPD, lwg, lwp)
+
+    Heff = H0matrix(Detg, Detp, u) + HImatrix(rabG, rabP, phidoppler, titadoppler, phiprobe, titaprobe, circularityprobe) + CC
+
+    # Heffdaga = np.matrix(Heff).getH()
+    Heffdaga = np.conj(np.transpose(Heff))
+
+    Lfullpartial = np.zeros((64, 64), dtype=np.complex_)
+
+    for r  in range(8):
+        for q in range(8):
+            for k in range(8):
+                for j in range(8):
+                    if j!=q and r!=k:
+                        pass
+                    elif j==q and r!=k:
+                        if (r < 2 and k > 3) or (k < 2 and r > 3) or (r > 3 and k > 3) or (r==0 and k==1) or (r==1 and k==0) or (r==2 and k==3) or (r==3 and k==2):
+                            #todo esto sale de analizar explicitamente la matriz y tratar de no calcular cosas de más que dan cero
+                            pass
+                        else:
+                            Lfullpartial[r*8+q][k*8+j] = (-1j)*(Heff[r,k])
+                    elif j!=q and r==k:
+                        if (j < 2 and q > 3) or (q < 2 and j > 3) or (j > 3 and q > 3) or (j==0 and q==1) or (j==1 and q==0) or (j==2 and q==3) or (j==3 and q==2):
+                            pass
+                        else:
+                            Lfullpartial[r*8+q][k*8+j] = (-1j)*(-Heffdaga[j,q])
+                    else:
+                        if Heff[r,k] == Heffdaga[j,q]:
+                            pass
+                        else:
+                            Lfullpartial[r*8+q][k*8+j] = (-1j)*(Heff[r,k]-Heffdaga[j,q])
+
+    M = CalculateSingleMmatrix(gPS, gPD, lwg, lwp)
+    # L0 = np.array(np.matrix(Lfullpartial) + M)
+    L0 = Lfullpartial + M
+
+    #ESTA PARTE ES CUANDO AGREGAS MICROMOCION
+    nmax = 3
+
+    #print(nmax)
+    Ltemp, Omega = LtempCalculus(beta, drivefreq)
+    #print(factor)
+    L1 = GetL1(Ltemp, L0, Omega, nmax)
+    Lfull = L0 + L1 #ESA CORRECCION ESTA EN L1
+    #HASTA ACA
+
+    #NORMALIZACION DE RHO
+    i = 0
+    while i < 64:
+        if i%9 == 0:
+            Lfull[0, i] = 1
+        else:
+            Lfull[0, i] = 0
+        i = i + 1
+
+    return Lfull
+
+
+"""
+Scripts para correr un experimento y hacer el análisis de los datos
+"""
+
+@njit
+def CPTspectrum8levels_MM(sg, sp, gPS, gPD, Detg, u, lwg, lwp, Temp, alpha, phidoppler, titadoppler, phiprobe, titaprobe, Circularityprobe, beta, drivefreq, freqMin=-100, freqMax=100, freqStep=1e-1, plot=False, solvemode=1):
+
+    """
+    ESTA ES LA FUNCION QUE ESTAMOS USANDO
+    Hace un experimento barriendo ángulos de repump con el angulo de doppler fijo.
+
+    solvemode=1: resuelve con np.linalg.solve
+    solvemode=2: resuelve invirtiendo L con la funcion np.linalg.inv
+    """
+    phidoppler, titadoppler = phidoppler*(np.pi/180), titadoppler*(np.pi/180)
+    phiprobe, titaprobe = phiprobe*(np.pi/180), titaprobe*(np.pi/180)
+
+    DetProbeVector = 2*np.pi*np.arange(freqMin*1e6, freqMax*1e6+0*freqStep*1e6, freqStep*1e6)
+    Detg = 2*np.pi*Detg*1e6
+    #lwg, lwr, lwp = 2*np.pi*lwg*1e6, 2*np.pi*lwr*1e6, 2*np.pi*lwp*1e6
+    lwg, lwp = lwg*1e6, lwp*1e6
+
+    rabG = sg*gPS
+    rabP = sp*gPD
+    #u = 2*np.pi*u*1e6
+
+    Fluovector = []
+
+    # tinicial = time.time()
+
+    for Detp in DetProbeVector:
+        L = FullL_MM(rabG, rabP, gPS, gPD, Detg, Detp, u, lwg, lwp, phidoppler, titadoppler, phiprobe, titaprobe, beta, drivefreq, Temp, alpha, Circularityprobe)
+        # if solvemode == 1:
+        rhovectorized = np.linalg.solve(L, np.array([int(i==0) for i in range(64)],dtype=np.complex_))
+        Fluo = np.real(rhovectorized[18] + np.real(rhovectorized[27]))
+        Fluovector.append(Fluo)
+
+        # if solvemode == 2:
+        #     Linv = np.linalg.inv(L)
+        #     rhovectorized = [Linv[j][0] for j in range(len(Linv))]
+        #     Fluo = np.real(rhovectorized[18] + np.real(rhovectorized[27])) #estos son los rho33 + rho44
+        #     Fluovector.append(Fluo)
+
+    # tfinal = time.time()
+
+    # print('Done, Total time: ', round((tfinal-tinicial), 2), "s")
+
+    DetProbeVectorMHz = np.arange(freqMin, freqMax, freqStep)
+
+    # if plot:
+    #     plt.xlabel('Probe detuning (MHz)')
+    #     plt.ylabel('Fluorescence (A.U.)')
+    #     plt.plot(DetProbeVectorMHz, [100*f for f in Fluovector], label=str(titaprobe) + 'º, T: ' + str(Temp*1e3) + ' mK')
+    #     plt.legend()
+
+    return DetProbeVectorMHz, Fluovector
+
+
+# @njit
+# def lolo():
+#     L = FullL_MM(100,200,12,123,14)
+#     return np.linalg.solve(L, np.array([int(i==0) for i in range(64)],dtype=np.complex_))
+
+# lolo()
+
+# raise ValueError('áaa')
+
+#%%
+if __name__ == "__main__":
+    ub = 9.27e-24
+    h = 6.63e-34
+    c = (ub/h)*1e-4 #en unidades de MHz/G
+
+    B = 25 #campo magnetico en gauss
+
+    u = c*B
+
+
+    sg, sr, sp = 0.5, 1.5, 4 #parámetros de saturación del doppler y repump
+    gPS, gPD, = 2*np.pi*21.58e6, 2*np.pi*1.35e6 #anchos de linea de las transiciones
+    rabG, rabR, rabP = sg*gPS, sr*gPD, sp*gPD #frecuencias de rabi
+
+    lwg, lwr, lwp = 0.3, 0.3, 0.3 #ancho de linea de los laseres
+
+    Detg = -25
+    Detr = 20 #detuning del doppler y repump
+
+    Temp = 0.0e-3 #temperatura en K
+    alpha = 0*(np.pi/180) #angulo entre los láseres
+
+    phidoppler, titadoppler = 0, 90
+    phirepump, titarepump = 0, 90
+
+    phiprobe, titaprobe = 0, 90
+
+
+    plotCPT = False
+
+    freqMin = -50
+    freqMax = 50
+    freqStep = 5e-2
+
+    # Frequencyvector, Fluovector = CPTspectrum8levels_MM(rabG, rabR, rabP, gPS, gPD, Detg, Detr,   u,  lwg,   lwr,        lwp, Temp, alpha, phidoppler, titadoppler, phiprobe, titaprobe,    phirepump, titarepump,             freqMin=freqMin, freqMax=freqMax, freqStep=freqStep, plot=plotCPT, solvemode=1)
+    Frequencyvector, Fluovector = PerformExperiment_8levels_MM(0.9,6.2,135591138.92893547,8482300.164692441,-24.5,32500000.0,0.1,0.1,0.001,0,0,90,0,90,2.0,139078306.77442014,-54.39999999999998,26.26666666666671,0.6666666666666856,circularityprobe=1,plot=False,solvemode=1,detpvec=None)
+
+
+    plt.plot(Frequencyvector, [100*f for f in Fluovector], label=str(titaprobe) + 'º, T: ' + str(Temp*1e3) + ' mK')
+    plt.xlabel('Probe detuning (MHz)')
+    plt.ylabel('Fluorescence (A.U.)')

analisis/plots/20231123_CPTconmicromocion3/lolo_CPT_plotter_20231123.py

+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+"""
+Ploteo de datos y ajustes
+@author: lolo
+"""
+
+
+import h5py
+import matplotlib.pyplot as plt
+import numpy as np
+# import sys
+# import re
+# import ast
+from scipy.optimize import curve_fit
+# import os
+# from scipy import interpolate
+
+
+
+
+#%% Importaciones extra
+
+# /home/lolo/Dropbox/marce/LIAF/Trampa_anular/artiq_experiments/analisis/plots/20231123_CPTconmicromocion3/Data/EITfit/MM_eightLevel_2repumps_AnalysisFunctions.py
+
+from Data.EITfit.lolo_modelo_full_8niveles import PerformExperiment_8levels_MM
+
+
+PARAMETROS = np.load('PARAMETROS.npz',allow_pickle=True)
+for var_name in PARAMETROS.keys():
+    globals()[var_name] = PARAMETROS[var_name]
+    print(f'loaded: {var_name}')
+
+
+#%%
+
+"""
+Primero tengo mediciones de espectros cpt de un ion variando la tension dc_A
+"""
+
+#C:\Users\Usuario\Documents\artiq\artiq_experiments\analisis\plots\20220106_CPT_DosLaseres_v08_TISA_DR\Data
+# 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
+"""
+
+
+#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]
+#SelectedCurveVec = [9]
+
+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()
+
+
+
+
+#%%
+def expo(x,tau,A,B):
+    return A*np.exp(x/tau)+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 )
+    
+