Skip to content
CPT_plotter_20220520.py 29.4 KiB
Newer Older
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

#Mediciones de CPT para campo magnetico bajo y terrestre

#/home/nico/Documents/artiq_experiments/analisis/plots/20220526_CPTvariandoB_v2/Data

Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
# os.chdir('/home/nico/Documents/artiq_experiments/analisis/plots/20220526_CPTvariandoB_v2/Data')

os.chdir("C://Users//nicon//Doctorado//artiq_experiments//analisis//plots//20220526_CPTvariandoB_v2//Data")


ALL_FILES = """000007697-IR_Scan_withcal_optimized
000007696-IR_Scan_withcal_optimized
000007695-IR_Scan_withcal_optimized
000007699-IR_Scan_withcal_optimized
000007700-IR_Scan_withcal_optimized
000007701-IR_Scan_withcal_optimized
000007702-IR_Scan_withcal_optimized
000007703-IR_Scan_withcal_optimized
000007712-IR_Scan_withcal_optimized
000007713-IR_Scan_withcal_optimized
000007714-IR_Scan_withcal_optimized
000007715-IR_Scan_withcal_optimized
000007716-IR_Scan_withcal_optimized
000007717-IR_Scan_withcal_optimized
000007718-IR_Scan_withcal_optimized
000007719-IR_Scan_withcal_optimized
000007720-IR_Scan_withcal_optimized
000007727-IR_Scan_withcal_optimized
000007730-IR_Scan_withcal_optimized
000007731-IR_Scan_withcal_optimized
000007732-IR_Scan_withcal_optimized
000007733-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(ALL_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 = []

for i, fname in enumerate(ALL_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']['IR_Frequencies']))
    Counts.append(np.array(data['datasets']['counts_spectrum']))
    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']))

Counts_B = []
Freqs_B = []


for i, fname in enumerate(ALL_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_B.append(np.array(data['datasets']['IR_Frequencies']))
    Counts_B.append(np.array(data['datasets']['counts_spectrum']))



#%% Estos bloques fitean curvas cpt pero se ven feas, no vale la pena mostrarlas
"""
#Barriendo angulo del IR con tisa apagado

jvec = [0, 1, 2, 3, 4, 5, 6, 7]
Bvec = [3, 2, 1, 0, -1, -1.3, -1.5, -1.7]

jsel = [4, 5, 6, 7]
plt.figure()
i = 0
for j in jvec:
    if j in jsel:
        plt.errorbar([2*f*1e-6 for f in Freqs[j]], Counts[j], yerr=np.sqrt(Counts[j]), fmt='o', label='B: {Bvec[i]}', capsize=2, markersize=2)
        #plt.plot([2*f*1e-6 for f in Freqs[j]], Counts[j], 'o-', label=f'Amp Tisa: {AmpTisa[i]}', mb  arkersize=3)
    i = i + 1
#plt.grid()
#plt.xlabel('Frecuencia (MHz)')
#plt.ylabel('counts')
#plt.legend()

#Barriendo angulo del IR con tisa apagado

#60 uW de IR

jvec = [17, 18, 19, 20, 21]
Bvec = [-2.2, -2, -1.75, -1.5, -1.2, -1, 0, 2]
jsel = jvec
plt.figure()
i = 0
for j in jvec:
    if j in jsel:
        plt.errorbar([2*f*1e-6 for f in Freqs[j]], Counts[j], yerr=np.sqrt(Counts[j]), fmt='o', label=f'B: {Bvec[i]}', capsize=2, markersize=2)
        #plt.plot([2*f*1e-6 for f in Freqs[j]], Counts[j], 'o-', label=f'Amp Tisa: {AmpTisa[i]}', mb  arkersize=3)
    i = i + 1
#plt.grid()
plt.xlabel('Frecuencia (MHz)')
plt.ylabel('counts')
plt.legend()
"""

#%%

#Barriendo angulo del IR con tisa apagado

#60 uW de IR

jvec = [10, 9, 11, 12, 13, 14, 15, 16]
Ibobinavec = [2, 0, -1, -1.2, -1.5, -1.75, -2, -2.2]

jsel = [9,11,12,13,14,15,16]
plt.figure()
i = 0
for j in jvec:
    if j in jsel:
        plt.errorbar([2*f*1e-6 for f in Freqs[j]], Counts[j], yerr=np.sqrt(Counts[j]), fmt='o', label=f'I={Ibobinavec[i]} A', capsize=2, markersize=2)
        #plt.plot([2*f*1e-6 for f in Freqs[j]], Counts[j], 'o-', label=f'Amp Tisa: {AmpTisa[i]}', mb  arkersize=3)
    i = i + 1
#plt.grid()
plt.xlabel('Frecuencia (MHz)')
plt.ylabel('counts')
plt.legend()


#%%
"""
Ahora empezamos a intentar fitear de a una mejor
"""
#Ajustamos una curva para obtener escala y offset

phidoppler, titadoppler = 0, 90
phirepump, titarepump = 0,  0
phiprobe = 0
titaprobe = 90

T = 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 = -5 #con 8 fitea bien

offsetxpi = 440+1+correccion
DetDoppler = -20.0-correccion

#fitsvec = [9, 11, 12, 13, 14, 15, 16] son las posibilidades

curvafiteada = 13

FreqsDRpi_3 = [2*f*1e-6-offsetxpi+14 for f in Freqs_B[curvafiteada]]
CountsDRpi_3 = Counts_B[curvafiteada]

freqslongpi_3 = np.arange(min(FreqsDRpi_3), max(FreqsDRpi_3)+FreqsDRpi_3[1]-FreqsDRpi_3[0], 0.1*(FreqsDRpi_3[1]-FreqsDRpi_3[0]))


def FitEITpiplotter(freqs, scale, offset, SG, SP, temp):
    U=285732
    #SG = 1.3
    #SP = 3.39e0
    #temp = 3.7e-4
    #temp=4e-3
    MeasuredFreq, MeasuredFluo = GenerateNoisyCPT_fit(SG, sr, SP, gPS, gPD, DetDoppler, DetRepump, U, DopplerLaserLinewidth, RepumpLaserLinewidth, ProbeLaserLinewidth, temp, alpha, phidoppler, titadoppler, phiprobe, [titaprobe], phirepump, titarepump, freqs, plot=False, solvemode=1, detpvec=None, noiseamplitude=noiseamplitude)
    FinalFluo = [f*scale + offset for f in MeasuredFluo]
    return FinalFluo

popt_tisaoffprueba, pcov_tisaoffprueba = curve_fit(FitEITpiplotter, FreqsDRpi_3, CountsDRpi_3, p0=[1e5, 1e2, 1, 4, 1e-3], bounds=((0,0,0,0,0), (1e6, 1e4, 10, 20, 10e-3)))


FittedEITpi_3 = FitEITpiplotter(freqslongpi_3, *popt_tisaoffprueba)

plt.figure()
plt.errorbar(FreqsDRpi_3, CountsDRpi_3, yerr=2*np.sqrt(CountsDRpi_3), fmt='o', capsize=2, markersize=2)
plt.plot(freqslongpi_3, FittedEITpi_3)
#plt.title(f'Sdop: {round(popt[0], 2)}, Spr: {round(popt[1], 2)}, T: {round(popt[2]*1e3, 2)} mK, detDop: {DetDoppler} MHz')

#%%
#ESTE CODIGO AJUSTA UNA DE LAS CURVAS, LA 12!!!

phidoppler, titadoppler = 0, 90
phirepump, titarepump = 0,  0
phiprobe = 0
titaprobe = 90

T = 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 = -5 #con 8 fitea bien

offsetxpi = 440+1+correccion
DetDoppler = -20.0-correccion

#jvec = [10, 9, 11, 12, 13, 14, 15, 16]
#Bvec = [2, 0, -1, -1.2, -1.5, -1.75, -2, -2.2]

curvafiteada = 12

FreqsDRpi_12 = [2*f*1e-6-offsetxpi+14 for f in Freqs_B[curvafiteada]]
CountsDRpi_12 = Counts_B[curvafiteada]

freqslongpi_12 = np.arange(min(FreqsDRpi_12), max(FreqsDRpi_12)+FreqsDRpi_12[1]-FreqsDRpi_12[0], 0.1*(FreqsDRpi_12[1]-FreqsDRpi_12[0]))


def FitEITpiplotter(freqs, U, SG, SP, DETDOPPLER):
    #U=285732
    #SG = 1.3
    #SP = 3.39e0
    #temp = 3.7e-4
    temp=0
    MeasuredFreq, MeasuredFluo = GenerateNoisyCPT_fit(SG, sr, SP, gPS, gPD, DETDOPPLER, DetRepump, U, DopplerLaserLinewidth, RepumpLaserLinewidth, ProbeLaserLinewidth, temp, alpha, phidoppler, titadoppler, phiprobe, [titaprobe], phirepump, titarepump, freqs, plot=False, solvemode=1, detpvec=None, noiseamplitude=noiseamplitude)
    FinalFluo = [f*2.0917e5 + 7.53e2 for f in MeasuredFluo]
    return FinalFluo

popt_12, pcov_12 = curve_fit(FitEITpiplotter, FreqsDRpi_12, CountsDRpi_12, p0=[285000, 1, 4, -15], bounds=((0,0,0,-20), (500000, 10, 20, -10)))


FittedEITpi_12 = FitEITpiplotter(freqslongpi_12, *popt_12)

plt.figure()
plt.errorbar(FreqsDRpi_12, CountsDRpi_12, yerr=2*np.sqrt(CountsDRpi_12), fmt='o', capsize=2, markersize=2)
plt.plot(freqslongpi_12, FittedEITpi_12)
#plt.title(f'Sdop: {round(popt[0], 2)}, Spr: {round(popt[1], 2)}, T: {round(popt[2]*1e3, 2)} mK, detDop: {DetDoppler} MHz')



#ESTE CODIGO AJUSTA UNA DE LAS CURVAS, LA 13!!! pero ahora ajustando demas parametros

phidoppler, titadoppler = 0, 90
phirepump, titarepump = 0,  0
phiprobe = 0
titaprobe = 90

T = 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 = -5 #con 8 fitea bien

offsetxpi = 440+1+correccion
DetDoppler = -20.0-correccion

#fitsvec = [9, 11, 12, 13, 14, 15, 16] son las posibilidades
#las corrientes son [2, 0, -1, -1.5, -1.75, -2, -2.2]

curvafiteada = 13

FreqsDRpi_13 = [2*f*1e-6-offsetxpi+14 for f in Freqs_B[curvafiteada]]
CountsDRpi_13 = Counts_B[curvafiteada]

freqslongpi_13 = np.arange(min(FreqsDRpi_13), max(FreqsDRpi_13)+FreqsDRpi_13[1]-FreqsDRpi_13[0], 0.1*(FreqsDRpi_13[1]-FreqsDRpi_13[0]))


def FitEITpiplotter(freqs, U, SG, SP, DETDOPPLER):
    #U=285732
    #SG = 1.3
    #SP = 3.39e0
    temp = 3.7e-10
    #temp=4e-3
    MeasuredFreq, MeasuredFluo = GenerateNoisyCPT_fit(SG, sr, SP, gPS, gPD, DETDOPPLER, DetRepump, U, DopplerLaserLinewidth, RepumpLaserLinewidth, ProbeLaserLinewidth, temp, alpha, phidoppler, titadoppler, phiprobe, [titaprobe], phirepump, titarepump, freqs, plot=False, solvemode=1, detpvec=None, noiseamplitude=noiseamplitude)
    FinalFluo = [f*2.0917e5 + 7.53e2 for f in MeasuredFluo]
    return FinalFluo

popt_13, pcov_13 = curve_fit(FitEITpiplotter, FreqsDRpi_13, CountsDRpi_13, p0=[285000, 1, 4, -15], bounds=((0,0,0,-20), (500000, 10, 20, 0)))


FittedEITpi_13 = FitEITpiplotter(freqslongpi_13, *popt_13)

plt.figure()
plt.errorbar(FreqsDRpi_13, CountsDRpi_13, yerr=2*np.sqrt(CountsDRpi_13), fmt='o', capsize=2, markersize=2)
plt.plot(freqslongpi_13, FittedEITpi_13)

#ESTE CODIGO AJUSTA UNA DE LAS CURVAS, LA 14!!!

phidoppler, titadoppler = 0, 90
phirepump, titarepump = 0,  0
phiprobe = 0
titaprobe = 90

T = 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 = -5 #con 8 fitea bien

offsetxpi = 440+1+correccion
DetDoppler = -20.0-correccion

#fitsvec = [9, 11, 12, 13, 14, 15, 16] son las posibilidades
#las corrientes son [2, 0, -1, -1.5, -1.75, -2, -2.2]

curvafiteada = 14

FreqsDRpi_14 = [2*f*1e-6-offsetxpi+14 for f in Freqs_B[curvafiteada]]
CountsDRpi_14 = Counts_B[curvafiteada]

freqslongpi_14 = np.arange(min(FreqsDRpi_14), max(FreqsDRpi_14)+FreqsDRpi_14[1]-FreqsDRpi_14[0], 0.1*(FreqsDRpi_14[1]-FreqsDRpi_14[0]))


def FitEITpiplotter(freqs, U, SG, SP, DETDOPPLER):
    #U=285732
    #SG = 1.3
    #SP = 3.39e0
    temp = 3.7e-10
    #temp=4e-3
    MeasuredFreq, MeasuredFluo = GenerateNoisyCPT_fit(SG, sr, SP, gPS, gPD, DETDOPPLER, DetRepump, U, DopplerLaserLinewidth, RepumpLaserLinewidth, ProbeLaserLinewidth, temp, alpha, phidoppler, titadoppler, phiprobe, [titaprobe], phirepump, titarepump, freqs, plot=False, solvemode=1, detpvec=None, noiseamplitude=noiseamplitude)
    FinalFluo = [f*2.0917e5 + 7.53e2 for f in MeasuredFluo]
    #FinalFluo = [f*scale + offset for f in MeasuredFluo]
    return FinalFluo

popt_14, pcov_14 = curve_fit(FitEITpiplotter, FreqsDRpi_14, CountsDRpi_14, p0=[285000, 1, 2.2, -15], bounds=((0,0,1,-20), (500000, 2.5, 20, 0)))


FittedEITpi_14 = FitEITpiplotter(freqslongpi_14, *popt_14)
#FittedEITpi_14 = FitEITpiplotter(freqslongpi_14, popt_14[0],popt_14[1],2.1,popt_14[3])

plt.figure()
plt.errorbar(FreqsDRpi_14, CountsDRpi_14, yerr=2*np.sqrt(CountsDRpi_14), fmt='o', capsize=2, markersize=2)
plt.plot(freqslongpi_14, FittedEITpi_14)
#plt.title(f'Sdop: {round(popt[0], 2)}, Spr: {round(popt[1], 2)}, T: {round(popt[2]*1e3, 2)} mK, detDop: {DetDoppler} MHz')

Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed

#ESTE CODIGO AJUSTA UNA DE LAS CURVAS, LA 15!!!

phidoppler, titadoppler = 0, 90
phirepump, titarepump = 0,  0
phiprobe = 0
titaprobe = 90

T = 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 = -5 #con 8 fitea bien

offsetxpi = 440+1+correccion
DetDoppler = -20.0-correccion

#fitsvec = [9, 11, 12, 13, 14, 15, 16] son las posibilidades
#las corrientes son [2, 0, -1, -1.5, -1.75, -2, -2.2]

curvafiteada = 15

FreqsDRpi_15 = [2*f*1e-6-offsetxpi+14 for f in Freqs_B[curvafiteada]]
CountsDRpi_15 = Counts_B[curvafiteada]

freqslongpi_15 = np.arange(min(FreqsDRpi_15), max(FreqsDRpi_15)+FreqsDRpi_15[1]-FreqsDRpi_15[0], 0.1*(FreqsDRpi_15[1]-FreqsDRpi_15[0]))


def FitEITpiplotter(freqs, U, SG, SP, DETDOPPLER):
    #U=285732
    #SG = 1.3
    #SP = 3.39e0
    temp = 3.7e-10
    #temp=4e-3
    MeasuredFreq, MeasuredFluo = GenerateNoisyCPT_fit(SG, sr, SP, gPS, gPD, DETDOPPLER, DetRepump, U, DopplerLaserLinewidth, RepumpLaserLinewidth, ProbeLaserLinewidth, temp, alpha, phidoppler, titadoppler, phiprobe, [titaprobe], phirepump, titarepump, freqs, plot=False, solvemode=1, detpvec=None, noiseamplitude=noiseamplitude)
    FinalFluo = [f*2.0917e5 + 7.53e2 for f in MeasuredFluo]
    #FinalFluo = [f*2.6917e5 + 3.53e2 for f in MeasuredFluo]
    return FinalFluo

popt_15, pcov_15 = curve_fit(FitEITpiplotter, FreqsDRpi_15, CountsDRpi_15, p0=[285000, 1, 4, -15], bounds=((0,0,0,-20), (500000, 10, 20, 0)))

print('B:')
print(popt_15[0])

FittedEITpi_15 = FitEITpiplotter(freqslongpi_15, *popt_15)

plt.figure()
plt.errorbar(FreqsDRpi_15, CountsDRpi_15, yerr=2*np.sqrt(CountsDRpi_15), fmt='o', capsize=2, markersize=2)
plt.plot(freqslongpi_15, FittedEITpi_15)
#plt.title(f'Sdop: {round(popt[0], 2)}, Spr: {round(popt[1], 2)}, T: {round(popt[2]*1e3, 2)} mK, detDop: {DetDoppler} MHz')


#%%

Detuningsfits = [popt_12[-1], popt_13[-1], popt_14[-1], popt_15[-1]]
PotIRfits = [popt_12[-2], popt_13[-2], popt_14[-2], popt_15[-2]]
PotUVfits = [popt_12[-3], popt_13[-3], popt_14[-3], popt_15[-3]]

ErrorDetuningfits = [np.sqrt(pcov_12)[3,3], np.sqrt(pcov_13)[3,3],np.sqrt(pcov_14)[3,3], np.sqrt(pcov_15)[3,3]]
ErrorPotIRfits = [np.sqrt(pcov_12)[2,2], np.sqrt(pcov_13)[2,2],np.sqrt(pcov_14)[2,2], np.sqrt(pcov_15)[2,2]]
ErrorPotUVfits = [np.sqrt(pcov_12)[1,1], np.sqrt(pcov_13)[1,1],np.sqrt(pcov_14)[1,1], np.sqrt(pcov_15)[1,1]]

Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
ErrorLarmorsvecs = [np.sqrt(pcov_12)[0,0], np.sqrt(pcov_13)[0,0],np.sqrt(pcov_14)[0,0], np.sqrt(pcov_15)[0,0]]

print(Detuningsfits)
print(PotIRfits)
print(PotUVfits)

plt.figure()
plt.errorbar([1,2,3,4],Detuningsfits, yerr=ErrorDetuningfits, fmt='o', capsize=2, markersize=15)

plt.figure()
plt.errorbar([1,2,3,4],PotIRfits, yerr=ErrorPotIRfits, fmt='o', capsize=2, markersize=15)

plt.figure()
plt.errorbar([1,2,3,4],PotUVfits, yerr=ErrorPotUVfits, fmt='o', capsize=2, markersize=15)

print('Media de detuning: ', np.mean(Detuningsfits), 'MHz', 'con media de error:', np.mean(ErrorDetuningfits),'MHz')
print('Media de potencia UV: ', np.mean(PotUVfits), 'con media de error:', np.mean(ErrorPotUVfits))
print('Media de potencia IR ', np.mean(PotIRfits), 'con media de error:', np.mean(ErrorPotIRfits))



#%%
#corrijo glitches xdxd

#12
CountsDRpi_12[97]=0.5*(CountsDRpi_12[96]+CountsDRpi_12[98])

#13
CountsDRpi_13[126]=(CountsDRpi_13[125]-CountsDRpi_13[128])*0.33+CountsDRpi_13[125]
CountsDRpi_13[127]=(CountsDRpi_13[125]-CountsDRpi_13[128])*0.66+CountsDRpi_13[125]

#14 ok

#15
CountsDRpi_15[65]=(CountsDRpi_15[64]-CountsDRpi_15[67])*0.33+CountsDRpi_15[64]
CountsDRpi_15[66]=(CountsDRpi_15[64]-CountsDRpi_15[67])*0.66+CountsDRpi_15[64]
CountsDRpi_15[83]=(CountsDRpi_15[82]+CountsDRpi_15[84])*0.5

#%%

#plot raw sin ajustar escalas

plt.figure()
plt.errorbar(FreqsDRpi_12, CountsDRpi_12, yerr=2*np.sqrt(CountsDRpi_12), fmt='o', capsize=2, markersize=2)
plt.plot(freqslongpi_12, FittedEITpi_12)
plt.errorbar(FreqsDRpi_13, CountsDRpi_13, yerr=2*np.sqrt(CountsDRpi_13), fmt='o', capsize=2, markersize=2)
plt.plot(freqslongpi_13, FittedEITpi_13)
plt.errorbar(FreqsDRpi_14, CountsDRpi_14, yerr=2*np.sqrt(CountsDRpi_14), fmt='o', capsize=2, markersize=2)
plt.plot(freqslongpi_14, FittedEITpi_14)
plt.errorbar(FreqsDRpi_15, CountsDRpi_15, yerr=2*np.sqrt(CountsDRpi_15), fmt='o', capsize=2, markersize=2)
plt.plot(freqslongpi_15, FittedEITpi_15)


#%%
#filtro y mejoro el plot
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
"""
Figura 4 del paper (apendice)
"""

Ufields = [popt_12[0], popt_13[0], popt_14[0], popt_15[0]]
ErrorUfields = [np.sqrt(pcov_12)[0,0], np.sqrt(pcov_13)[0,0],np.sqrt(pcov_14)[0,0], np.sqrt(pcov_15)[0,0]]

Bfields = [u/(2*np.pi)/1398190 for u in Ufields]
ErrorBfields = [eu/(2*np.pi)/1398190 for eu in ErrorUfields]


from scipy.signal import savgol_filter as sf
import seaborn as sns

offs = 7.53e2
scal = 2.0917e5

winl = 11
po = 3

FiltCounts12 = np.array([(c-offs)*100/scal for c in sf(CountsDRpi_12, winl, po)])
#ErrorCounts12 = np.array([0.5*np.sqrt(c/scal) for c in CountsDRpi_12])
ErrorCounts12 = np.array([0.5*np.sqrt(c/scal) for c in sf(CountsDRpi_12, winl, po)])

FiltCounts13 = np.array([(c-offs)*100/scal for c in sf(CountsDRpi_13, winl, po)])
#ErrorCounts13 = np.array([0.5*np.sqrt(c/scal) for c in CountsDRpi_13])
ErrorCounts13 = np.array([0.5*np.sqrt(c/scal) for c in sf(CountsDRpi_13, winl, po)])

FiltCounts14 = np.array([(c-offs)*100/scal for c in sf(CountsDRpi_14, winl, po)])
#ErrorCounts14 = np.array([0.5*np.sqrt(c/scal) for c in CountsDRpi_14])
ErrorCounts14 = np.array([0.5*np.sqrt(c/scal) for c in sf(CountsDRpi_14, winl, po)])

FiltCounts15 = np.array([(c-offs)*100/scal for c in sf(CountsDRpi_15, winl, po)])
#ErrorCounts15 = np.array([0.5*np.sqrt(c/scal) for c in CountsDRpi_15])
ErrorCounts15 = np.array([0.5*np.sqrt(c/scal) for c in sf(CountsDRpi_15, winl, po)])


#plot con escalas atomicas y lindo

Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
colors=sns.color_palette("mako", 10)

color1=colors[1]
color2=colors[2]
color3=colors[4]
color4=colors[7]

ms = 4

plt.figure(figsize=(3.5, 3))

Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
plt.plot(FreqsDRpi_12, FiltCounts12, 'o', markersize=ms, color=color1, label=fr'$B=${round(1e3*Bfields[0])}({round(1e3*ErrorBfields[0])}) mG', alpha=0.3)
plt.plot(FreqsDRpi_13, FiltCounts13, 'o', markersize=ms, color=color2, label=fr'$B=${round(1e3*Bfields[1])}({round(1e3*ErrorBfields[1])}) mG', alpha=0.3)
plt.plot(FreqsDRpi_14, FiltCounts14, 'o', markersize=ms, color=color3, label=fr'$B=${round(1e3*Bfields[2])}({round(1e3*ErrorBfields[2])}) mG', alpha=0.3)
plt.plot(FreqsDRpi_15, FiltCounts15, 'o', markersize=ms, color=color4, label=fr'$B=${round(1e3*Bfields[3])}({round(1e3*ErrorBfields[3])}) mG', alpha=0.3)

plt.plot(freqslongpi_12, [(c-offs)*100/scal for c in FittedEITpi_12], color=color1, linewidth=3)
plt.plot(freqslongpi_13, [(c-offs)*100/scal for c in FittedEITpi_13], color=color2, linewidth=3)
plt.plot(freqslongpi_14, [(c-offs)*100/scal for c in FittedEITpi_14], color=color3, linewidth=3)
plt.plot(freqslongpi_15, [(c-offs)*100/scal for c in FittedEITpi_15], color=color4, linewidth=3)

plt.fill_between(FreqsDRpi_12, FiltCounts12+ErrorCounts12, FiltCounts12-ErrorCounts12, color=color1, alpha=0.2)
plt.fill_between(FreqsDRpi_13, FiltCounts13+ErrorCounts13, FiltCounts13-ErrorCounts13, color=color2, alpha=0.2)
plt.fill_between(FreqsDRpi_14, FiltCounts14+ErrorCounts14, FiltCounts14-ErrorCounts14, color=color3, alpha=0.2)
plt.fill_between(FreqsDRpi_15, FiltCounts15+ErrorCounts15, FiltCounts15-ErrorCounts15, color=color4, alpha=0.2)

plt.xlim(-40,30)
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
plt.ylim(-0.1,1.4)
plt.xlabel('IR laser detuning (MHz)', fontsize=11, fontname='STIXgeneral')
plt.ylabel('Normalized fluorescence', fontsize=11, fontname='STIXgeneral')
plt.xticks([-40, -20, 0, 20], fontsize=11, fontname='STIXgeneral')
plt.yticks([0, 0.3, 0.6, 0.9, 1.2], fontsize=11, fontname='STIXgeneral')
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
visible_ticks = {
   "top": False,
   "right": False
}
plt.tick_params(axis="x", which="both", **visible_ticks)
plt.tick_params(axis="y", which="both", **visible_ticks)
#plt.legend(loc='upper left', frameon=True, fontsize=7.6, handletextpad=0.1)
plt.grid()
plt.tight_layout()
#plt.savefig('/home/nico/Nextcloud/G_liaf/Publicaciones/Work/2022 B vs k race/Figuras/Figuras jpg trabajadas/CPT_exp.png',dpi=500)
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
# plt.savefig('/home/nico/Nextcloud/G_liaf/Publicaciones/Papers/2022 B vs K eigenbasis/Figuras_finales/Finalesfinales/Fig4_final.pdf')

#%%
from scipy.optimize import curve_fit

def LinearLarmortoCurrent(I, a, b):
    Larmor = a*I+b
    return Larmor

def ConvertLarmortoBfield(u):
    c = 1398190.0452488689
    return u/(2*np.pi)/c

"""
mediciones = [10, 9, 11,  12,   13,    14, 15,   16]
corrientes = [2, 0, -1, -1.2, -1.5, -1.75, -2, -2.2]
"""

"""
#Esto anda:
IFit = [-1, -1.5, -1.75]
LarmorFit = [popt_13[0], popt_14[0], popt_15[0]]
"""

IFit = [-1.2, -1.5, -1.75, -2]
LarmorFit = [popt_12[0], popt_13[0], popt_14[0], popt_15[0]]
#LarmorFit = [287282.374931893, 203998.09641511342, 158507.0255951109] por si fallan los ajustes

Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
ErrorLarmorFit = [np.sqrt(pcov_12[0,0]),np.sqrt(pcov_13[0,0]),np.sqrt(pcov_14[0,0]),np.sqrt(pcov_15[0,0])]
MeanError = 11139.353180216529 #por si fallan los ajustes

Ilong = np.arange(2, -3, -0.01)
popt_larmor, pcov_larmor = curve_fit(LinearLarmortoCurrent, IFit, LarmorFit)
LarmorLong = LinearLarmortoCurrent(Ilong, *popt_larmor)

print(popt_larmor)

plt.figure()
plt.plot(IFit, LarmorFit, 'o', markersize=5)
plt.plot(Ilong, LarmorLong)

Bfitted = [ConvertLarmortoBfield(u) for u in LarmorFit]
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
ErrorBfitted = [ConvertLarmortoBfield(u) for u in ErrorLarmorFit]
BLong = [ConvertLarmortoBfield(u) for u in LarmorLong]

Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed

#ESTE GRAFICO TIENE QUE IR EN LA TESIS CUANDO PONGA
#COMO CALIBRE EL CAMPO MAGNETICO CON LA CORRIENTE
plt.figure()
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
plt.plot(Ilong, [b*1e3 for b in BLong])
plt.plot(IFit, [b*1e3 for b in Bfitted], 'o', markersize=8)
plt.xlim(-2.2,-1)
plt.ylim(0.012*1e3,0.045*1e3)
plt.xlabel('Corriente (A)')
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
plt.ylabel('Campo magnetico (mG)')
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
#%%
plt.figure()
plt.plot(Ilong, [b*1e-6 for b in LarmorLong])
plt.errorbar(IFit, [b*1e-6 for b in LarmorFit], yerr=1e-6*MeanError, fmt='o', capsize=4, markersize=4)
plt.xlim(-2.2,-1)
plt.ylim(0.1,0.43)
plt.xlabel('Coil current (A)')
plt.ylabel('Larmor frequency (MHz)')
plt.grid()
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed



#%%
#filtro y mejoro el plot
"""
PLOTS TESIS (LOS DEL PAPER SON LOS DE ANTES)
"""

Ufields = [popt_12[0], popt_13[0], popt_14[0], popt_15[0]]
ErrorUfields = [np.sqrt(pcov_12)[0,0], np.sqrt(pcov_13)[0,0],np.sqrt(pcov_14)[0,0], np.sqrt(pcov_15)[0,0]]

Bfields = [u/(2*np.pi)/1398190 for u in Ufields]
ErrorBfields = [eu/(2*np.pi)/1398190 for eu in ErrorUfields]


from scipy.signal import savgol_filter as sf
import seaborn as sns

offs = 7.53e2
scal = 2.0917e5

winl = 11
po = 3

FiltCounts12 = np.array([(c-offs)*100/scal for c in sf(CountsDRpi_12, winl, po)])
#ErrorCounts12 = np.array([0.5*np.sqrt(c/scal) for c in CountsDRpi_12])
ErrorCounts12 = np.array([0.5*np.sqrt(c/scal) for c in sf(CountsDRpi_12, winl, po)])

FiltCounts13 = np.array([(c-offs)*100/scal for c in sf(CountsDRpi_13, winl, po)])
#ErrorCounts13 = np.array([0.5*np.sqrt(c/scal) for c in CountsDRpi_13])
ErrorCounts13 = np.array([0.5*np.sqrt(c/scal) for c in sf(CountsDRpi_13, winl, po)])

FiltCounts14 = np.array([(c-offs)*100/scal for c in sf(CountsDRpi_14, winl, po)])
#ErrorCounts14 = np.array([0.5*np.sqrt(c/scal) for c in CountsDRpi_14])
ErrorCounts14 = np.array([0.5*np.sqrt(c/scal) for c in sf(CountsDRpi_14, winl, po)])

FiltCounts15 = np.array([(c-offs)*100/scal for c in sf(CountsDRpi_15, winl, po)])
#ErrorCounts15 = np.array([0.5*np.sqrt(c/scal) for c in CountsDRpi_15])
ErrorCounts15 = np.array([0.5*np.sqrt(c/scal) for c in sf(CountsDRpi_15, winl, po)])


#plot con escalas atomicas y lindo

colors=sns.color_palette("mako", 10)

color1=colors[1]
color2=colors[2]
color3=colors[4]
color4=colors[7]

ms = 4

plt.figure(figsize=(3.6, 3))

# plt.plot(FreqsDRpi_12, FiltCounts12, 'o', markersize=ms, color=color1, label=fr'$B=${round(1e3*Bfields[0])}({round(1e3*ErrorBfields[0])}) mG', alpha=0.3)
# plt.plot(FreqsDRpi_13, FiltCounts13, 'o', markersize=ms, color=color2, label=fr'$B=${round(1e3*Bfields[1])}({round(1e3*ErrorBfields[1])}) mG', alpha=0.3)
# plt.plot(FreqsDRpi_14, FiltCounts14, 'o', markersize=ms, color=color3, label=fr'$B=${round(1e3*Bfields[2])}({round(1e3*ErrorBfields[2])}) mG', alpha=0.3)
# plt.plot(FreqsDRpi_15, FiltCounts15, 'o', markersize=ms, color=color4, label=fr'$B=${round(1e3*Bfields[3])}({round(1e3*ErrorBfields[3])}) mG', alpha=0.3)

# plt.plot(freqslongpi_12, [(c-offs)*100/scal for c in FittedEITpi_12], color=color1, linewidth=3)
# plt.plot(freqslongpi_13, [(c-offs)*100/scal for c in FittedEITpi_13], color=color2, linewidth=3)
# plt.plot(freqslongpi_14, [(c-offs)*100/scal for c in FittedEITpi_14], color=color3, linewidth=3)
# plt.plot(freqslongpi_15, [(c-offs)*100/scal for c in FittedEITpi_15], color=color4, linewidth=3)

plt.plot(freqslongpi_12, [c for c in FittedEITpi_12], color=color1, linewidth=3)
plt.plot(freqslongpi_13, [c for c in FittedEITpi_13], color=color2, linewidth=3)
plt.plot(freqslongpi_14, [c for c in FittedEITpi_14], color=color3, linewidth=3)
plt.plot(freqslongpi_15, [c for c in FittedEITpi_15], color=color4, linewidth=3)

plt.errorbar(FreqsDRpi_12, [c*scal/100+offs for c in FiltCounts12], yerr=np.sqrt(np.array([c*scal/100+offs for c in FiltCounts12])),  capsize=2,fmt='o', markersize=2, color=color1, label=fr'$B=${round(1e3*Bfields[0])}({round(1e3*ErrorBfields[0])}) mG', alpha=0.7)
plt.errorbar(FreqsDRpi_13, [c*scal/100+offs for c in FiltCounts13], yerr=np.sqrt(np.array([c*scal/100+offs for c in FiltCounts13])),  capsize=2,fmt='o', markersize=2, color=color2, label=fr'$B=${round(1e3*Bfields[1])}({round(1e3*ErrorBfields[1])}) mG', alpha=0.7)
plt.errorbar(FreqsDRpi_14, [c*scal/100+offs for c in FiltCounts14], yerr=np.sqrt(np.array([c*scal/100+offs for c in FiltCounts14])),  capsize=2,fmt='o', markersize=2, color=color3, label=fr'$B=${round(1e3*Bfields[2])}({round(1e3*ErrorBfields[2])}) mG', alpha=0.7)
plt.errorbar(FreqsDRpi_15, [c*scal/100+offs for c in FiltCounts15], yerr=np.sqrt(np.array([c*scal/100+offs for c in FiltCounts15])),  capsize=2,fmt='o', markersize=2, color=color4, label=fr'$B=${round(1e3*Bfields[3])}({round(1e3*ErrorBfields[3])}) mG', alpha=0.7)


# plt.plot(FreqsDRpi_13, [c*scal/100+offs for c in FiltCounts13], 'o', markersize=ms, color=color2, label=fr'$B=${round(1e3*Bfields[1])}({round(1e3*ErrorBfields[1])}) mG', alpha=0.3)
# plt.plot(FreqsDRpi_14, [c*scal/100+offs for c in FiltCounts14], 'o', markersize=ms, color=color3, label=fr'$B=${round(1e3*Bfields[2])}({round(1e3*ErrorBfields[2])}) mG', alpha=0.3)
# plt.plot(FreqsDRpi_15, [c*scal/100+offs for c in FiltCounts15], 'o', markersize=ms, color=color4, label=fr'$B=${round(1e3*Bfields[3])}({round(1e3*ErrorBfields[3])}) mG', alpha=0.3)



# plt.fill_between(FreqsDRpi_12, FiltCounts12+ErrorCounts12, FiltCounts12-ErrorCounts12, color=color1, alpha=0.2)
# plt.fill_between(FreqsDRpi_13, FiltCounts13+ErrorCounts13, FiltCounts13-ErrorCounts13, color=color2, alpha=0.2)
# plt.fill_between(FreqsDRpi_14, FiltCounts14+ErrorCounts14, FiltCounts14-ErrorCounts14, color=color3, alpha=0.2)
# plt.fill_between(FreqsDRpi_15, FiltCounts15+ErrorCounts15, FiltCounts15-ErrorCounts15, color=color4, alpha=0.2)

plt.xlim(-40,30)
# plt.ylim(-0.1,1.4)
plt.xlabel('Desintonía IR (MHz)', fontsize=10)
plt.ylabel('Cuentas', fontsize=10)
plt.xticks([-40, -20, 0, 20], fontsize=10)
# plt.yticks([0, 0.3, 0.6, 0.9, 1.2], fontsize=10)
# visible_ticks = {
#    "top": False,
#    "right": False
# }
# plt.tick_params(axis="x", which="both", **visible_ticks)
# plt.tick_params(axis="y", which="both", **visible_ticks)
# plt.legend(loc='upper left', frameon=True, fontsize=7.6, handletextpad=0.1)
plt.grid(alpha=0.5)
plt.tight_layout()
#plt.savefig('/home/nico/Nextcloud/G_liaf/Publicaciones/Work/2022 B vs k race/Figuras/Figuras jpg trabajadas/CPT_exp.png',dpi=500)
# plt.savefig('/home/nico/Nextcloud/G_liaf/Publicaciones/Papers/2022 B vs K eigenbasis/Figuras_finales/Finalesfinales/Fig4_final.pdf')

plt.savefig('C://Users//nicon//Nextcloud//Nico//Doctorado//Tesis//Tesis_doctorado//Chapters//figures//Cap6//fig2test.pdf')


#%%
from scipy.optimize import curve_fit

def LinearLarmortoCurrent(I, a, b):
    Larmor = a*I+b
    return Larmor

def ConvertLarmortoBfield(u):
    c = 1398190.0452488689
    return u/(2*np.pi)/c

"""
mediciones = [10, 9, 11,  12,   13,    14, 15,   16]
corrientes = [2, 0, -1, -1.2, -1.5, -1.75, -2, -2.2]
"""

"""
#Esto anda:
IFit = [-1, -1.5, -1.75]
LarmorFit = [popt_13[0], popt_14[0], popt_15[0]]
"""

IFit = [-1.2, -1.5, -1.75, -2]
LarmorFit = [popt_12[0], popt_13[0], popt_14[0], popt_15[0]]
#LarmorFit = [287282.374931893, 203998.09641511342, 158507.0255951109] por si fallan los ajustes

ErrorLarmorFit = [np.sqrt(pcov_12[0,0]),np.sqrt(pcov_13[0,0]),np.sqrt(pcov_14[0,0]),np.sqrt(pcov_15[0,0])]
MeanError = 11139.353180216529 #por si fallan los ajustes

Ilong = np.arange(2, -3, -0.01)
popt_larmor, pcov_larmor = curve_fit(LinearLarmortoCurrent, IFit, LarmorFit)
LarmorLong = LinearLarmortoCurrent(Ilong, *popt_larmor)

print(popt_larmor)

plt.figure()
plt.plot(IFit, LarmorFit, 'o', markersize=5)
plt.plot(Ilong, LarmorLong)

Bfitted = [ConvertLarmortoBfield(u) for u in LarmorFit]
ErrorBfitted = [ConvertLarmortoBfield(u) for u in ErrorLarmorFit]
BLong = [ConvertLarmortoBfield(u) for u in LarmorLong]


#ESTE GRAFICO TIENE QUE IR EN LA TESIS CUANDO PONGA
#COMO CALIBRE EL CAMPO MAGNETICO CON LA CORRIENTE
plt.figure(figsize=(3.6,3))
plt.plot(Ilong, [b*1e3 for b in BLong])
plt.plot(IFit, [b*1e3 for b in Bfitted], 'o', markersize=8)
plt.xlim(-2.2,-1)
plt.ylim(0.012*1e3,0.045*1e3)
plt.xlabel('Corriente bobina superior (A)')
plt.ylabel('Campo magnetico (mG)')
plt.grid()