Skip to content
UV_CPT_spectrums.py 5.82 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

# Solo levanto algunos experimentos
Muriel Bonetto's avatar
Muriel Bonetto committed
Calib_Files = """000007324-UV_Scan_withcalib_Haeffner
000007325-IR_Scan_withcal_optimized
000007326-UV_Scan_withcalib_Haeffner
000007327-IR_Scan_withcal_optimized
000007328-UV_Scan_withcalib_Haeffner
000007327-IR_Scan_withcal_optimized
000007197-IR_Scan_withcal_optimized
000007198-UV_Scan_withcalib_Haeffner""" 
#directory = '/home/liaf-murib/Documents/Artiq/artiq_experiments/analisis/plots/20220520_EspectrosUVyCPT_descompensacion/Data/'
directory = '/home/nico/Documents/artiq_experiments/analisis/plots/20220520_EspectrosUVyCPT_descompensacion/Data/'

#carpeta pc nico labo escritorio:
#C:\Users\Usuario\Documents\artiq\artiq_experiments\analisis\plots\20220503_EspectrosUVnuevos\Data

def SeeKeys(files,directory = ''):
    for i, fname in enumerate(files.split()):
        data = h5py.File(directory + fname+'.h5', 'r') # Leo el h5: Recordar que nuestros datos estan en 'datasets'
        print(fname)
        print(list(data['datasets'].keys()))

#%%    

Amps = []
Freqs = []
Counts = []
T_readouts = []

for i, fname in enumerate(Calib_Files.split()):
    print(SeeKeys(Calib_Files,directory = directory))
    print(i)
    print(fname)
    data = h5py.File(directory + fname+'.h5', 'r') # Leo el h5: Recordar que nuestros datos estan en 'datasets'
    print(list(data['datasets'].keys()))
Muriel Bonetto's avatar
Muriel Bonetto committed
    try:
        Amps.append(np.array(data['datasets']['UV_Amplitudes']))
    except KeyError:
        Amps.append(np.array(data['datasets']['IR_Amplitudes']))
    try:
        Freqs.append(np.array(data['datasets']['UV_Frequencies']))
    except KeyError:
        Freqs.append(np.array(data['datasets']['IR_Frequencies']))
    Counts.append(np.array(data['datasets']['counts_spectrum']))
    T_readouts.append(np.array(data['datasets']['t_readout']))    


#def GetBackground(countsper100ms, )



#%%
from scipy.special import jv

def Lorentzian(f, A, x0, gamma, offset):
    return (A/np.pi)*0.5*gamma/(((f-x0)**2)+((0.5*gamma)**2)) + offset #40 es el piso de ruido estimado


Muriel Bonetto's avatar
Muriel Bonetto committed
jvec = [4] #UV_cooling en 90 MHz

plt.figure()

for j in jvec:

    FreqsChosen = [2*f*1e-6 for f in Freqs[j]]
    CountsChosen = Counts[j]
    
    portion = 0.
    
Muriel Bonetto's avatar
Muriel Bonetto committed
    #popt, pcov = curve_fit(Lorentzian, FreqsChosen[int(portion*len(FreqsChosen)):], CountsChosen[int(portion*len(FreqsChosen)):], p0=[12000, 208, 30, 30])
    
    freqslong = np.arange(min(FreqsChosen)-10, max(FreqsChosen)+10, (FreqsChosen[1]-FreqsChosen[0])*0.01)
    
    plt.errorbar(FreqsChosen, CountsChosen, yerr=np.sqrt(np.array(CountsChosen)), fmt='o', capsize=5, markersize=5)
    #plt.plot(freqslong, Lorentzian(freqslong, popt[0], popt[1], popt[2]), label=f'FWHM {round(popt[1])} MHz')
Muriel Bonetto's avatar
Muriel Bonetto committed
    #plt.plot(freqslong, Lorentzian(freqslong, popt[0], popt[1], popt[2], popt[3]), label=f'FWHM 30 MHz')
    #plt.axvline(popt[2]-22.1, linestyle='--', linewidth=1)
    #plt.axvline(popt[2]+22.1, linestyle='--', linewidth=1)
    plt.xlabel('Frecuencia (MHz)')
    plt.ylabel('Cuentas')
    plt.legend()
    
Muriel Bonetto's avatar
Muriel Bonetto committed
    #print(f'Ancho medido: {round(popt[2])} MHz') 

#%%




#%%
from scipy.special import jv
from scipy.optimize import curve_fit

def Lorentzian(f, A, gamma, x0):
    return (A/np.pi)*0.5*gamma/(((f-x0)**2)+((0.5*gamma)**2))

def MicromotionSpectra(det, A, beta, x0, offset):
    ftrap=22.1
Muriel Bonetto's avatar
Muriel Bonetto committed
    gamma= 29
    P = A*(jv(0, beta)**2)/(((det-x0)**2)+(0.5*gamma)**2)+ offset
    i = 1
    #print(P)
    while i <= 5:
        P = P + A*((jv(i, beta))**2)/((((det-x0)+i*ftrap)**2)+(0.5*gamma)**2) + A*((jv(-i, beta))**2)/((((det-x0)-i*ftrap)**2)+(0.5*gamma)**2) 
        i = i + 1
        #print(P)
    return P


jvec = [7] #UV_cooling en 90 MHz

"""
plt.figure()

for j in jvec:

    FreqsChosen = [2*f*1e-6 for f in Freqs[j]]
    CountsChosen = Counts[j]
    
    portion = 0.
    
    popt, pcov = curve_fit(Lorentzian, FreqsChosen[int(portion*len(FreqsChosen)):], CountsChosen[int(portion*len(FreqsChosen)):], p0=[12000, 25, 208, 30])
    
    freqslong = np.arange(min(FreqsChosen)-10, max(FreqsChosen)+10, (FreqsChosen[1]-FreqsChosen[0])*0.01)
    
    plt.errorbar(FreqsChosen, CountsChosen, yerr=np.sqrt(np.array(CountsChosen)), fmt='o', capsize=5, markersize=5)
    #plt.plot(freqslong, Lorentzian(freqslong, popt[0], popt[1], popt[2]), label=f'FWHM {round(popt[1])} MHz')
    plt.plot(freqslong, Lorentzian(freqslong, popt[0], popt[1], popt[2], popt[3]), label=f'FWHM 30 MHz')
    plt.axvline(popt[2]+2*22.1, linestyle='--', linewidth=1)
    plt.axvline(popt[2]+22.1, linestyle='--', linewidth=1)
    plt.xlabel('Frecuencia (MHz)')
    plt.ylabel('Cuentas')
    plt.legend()
    
    print(f'Ancho medido: {round(popt[1])} MHz') 
"""    
   
plt.figure()

for j in jvec:

    FreqsChosen = [2*f*1e-6 for f in Freqs[j]]
    CountsChosen = Counts[j]
    
    portion = 0.
    
    popt, pcov = curve_fit(MicromotionSpectra, FreqsChosen[int(portion*len(FreqsChosen)):], CountsChosen[int(portion*len(FreqsChosen)):],p0=[100000,1.5,220,200], bounds=((7000,0.1,200,-500),(1000000,10,300,500)))
    
    freqslong = np.arange(min(FreqsChosen)-10, max(FreqsChosen)+10, (FreqsChosen[1]-FreqsChosen[0])*0.01)
    
    plt.errorbar(FreqsChosen, CountsChosen, yerr=np.sqrt(np.array(CountsChosen)), fmt='o', capsize=5, markersize=5)
    #plt.plot(freqslong, Lorentzian(freqslong, popt[0], popt[1], popt[2]), label=f'FWHM {round(popt[1])} MHz')
    plt.plot(freqslong, MicromotionSpectra(freqslong, *popt), label='Beta ={:0.2f}'.format(popt[1]))
    #plt.plot(freqslong, MicromotionSpectra(freqslong, 70000,0.2,215,220), label=f'FWHM 30 MHz')
    #plt.axvline(popt[2]+2*22.1, linestyle='--', linewidth=1)
    #plt.axvline(popt[2]+22.1, linestyle='--', linewidth=1)
    plt.xlabel('Frecuencia (MHz)')
    plt.ylabel('Cuentas')
    plt.legend()
    
    print(f'Beta medido: {popt[1]}') 
Muriel Bonetto's avatar
Muriel Bonetto committed
    print(popt)