Skip to content
RDS_piezo.py 18 KiB
Newer Older
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
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 barriendo angulo del TISA y viendo kicking de resonancias oscuras

#C:\Users\Usuario\Documents\artiq\artiq_experiments\analisis\plots\20220106_CPT_DosLaseres_v08_TISA_DR\Data

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


"""
en este codigo ploteo espectros CPT de resonancias D-D para configuracion +2/+2 y +2/-2 (usando pentaprisma)
"""

def find_nearest(array, value):
    array = np.asarray(array)
    idx = (np.abs(array - value)).argmin()
    return idx


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


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()))


Piezo1Counts = []
Piezo1Frequencies = []

PIEZO1_FILES = np.arange(791, 826,1)

for i in PIEZO1_FILES:
    #print(str(i) + ' - ' + fname)
    data = h5py.File(f'VaryingBeamlocation/Piezo/000014{i}-IR_Scan_withcal_optimized'+'.h5', 'r')
    Piezo1Counts.append(np.array(data['datasets']['counts_spectrum']))
    Piezo1Frequencies.append(np.array(data['datasets']['IR1_Frequencies']))


Piezo2Counts = []
Piezo2Frequencies = []

PIEZO2_FILES = list(np.arange(834, 841,1))+list(np.arange(842, 872,1))

for i in PIEZO2_FILES:
    #print(str(i) + ' - ' + fname)
    data = h5py.File(f'VaryingBeamlocation/Piezo2/000014{i}-IR_Scan_withcal_optimized'+'.h5', 'r')
    Piezo2Counts.append(np.array(data['datasets']['counts_spectrum']))
    Piezo2Frequencies.append(np.array(data['datasets']['IR1_Frequencies']))


def ErrorDRdepth(p, f, b):
    ep = np.sqrt(p)
    ef = np.sqrt(f)
    eb = np.sqrt(b)
    derivadap = 1/((f-b)**2)
    derivadaf = ((p-b)/((f-b)**2))**2
    derivadab = ((p-f)/((f-b)**2))**2
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
    return 1*np.sqrt(derivadap*ep*ep + derivadaf*ef*ef + derivadab*eb*eb)
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed

Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
def Lorentzian( x, A, B, x0, gam,C):
    #C=0
    return A * gam**2 / ( gam**2 + ( x - x0 )**2) + B - C*(x - x0)


Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
#%%

import seaborn as sns
"""
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
Resonancias DD configuracion -2/-2 colineal variando la ubicacion del ion en los haces
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed

TODO EL OAM
"""

palette = sns.color_palette("tab10")

pmlocmedvec = list(np.arange(0,12,1))+[13,12]+list(np.arange(15,len(PIEZO1_FILES),1))

"""
Hay que invertir la 12 con la 13, y la 14 es la misma que la 12, por las dudas
"""


#pmlocmedvec = [21]

idxvecdr = [85,185, 185,182,182,232, 162,162,175,177,217, 217,184,184,188,186,186,195,195,198,195,198,196,192,192,192,178,175,175,170,140,140,190,185]


plt.figure()

bkg = np.min(Piezo1Counts[1])

pmdepthsdr=[]
errorpmdepthsdr=[]

Intensity = []
errorIntensity = []

idxtest = 185
print(idxtest)

jj=0
for med in pmlocmedvec:
    print(med)
    
    if med == 21 or med == 22 or med == 23:
        pmdepthsdr.append(1-(Piezo1Counts[med][1:][idxvecdr[jj]]-bkg)/(np.mean(Piezo1Counts[med][1:][0:20])-bkg))
        errorpmdepthsdr.append(ErrorDRdepth(Piezo1Counts[med][1:][idxvecdr[jj]],np.mean(Piezo1Counts[med][1:][0:20]), bkg))
            
        Intens = np.mean(Piezo1Counts[med][1:][0:10])-bkg
        
        Intensity.append(Intens)
        errorIntensity.append(2*np.sqrt(np.mean(Piezo1Counts[med][1:][0:10]))+np.sqrt(bkg))

    else:
        pmdepthsdr.append(1-(Piezo1Counts[med][1:][idxvecdr[jj]]-bkg)/(np.mean(Piezo1Counts[med][1:][0:20])-bkg))
        errorpmdepthsdr.append(ErrorDRdepth(Piezo1Counts[med][1:][idxvecdr[jj]],np.mean(Piezo1Counts[med][1:][0:20]), bkg))
            
        Intens = np.mean(Piezo1Counts[med][1:][0:20])-bkg
        
        Intensity.append(Intens)
        errorIntensity.append(2*np.sqrt(np.mean(Piezo1Counts[med][1:][0:20]))+np.sqrt(bkg))
        
    
        
    plt.plot([2*f*1e-6 for f in Piezo1Frequencies[med][1:]], [c for c in Piezo1Counts[med][1:]], '-o', markersize=2, alpha=0.7)
    plt.plot([2*f*1e-6 for f in Piezo1Frequencies[med][1:]][idxtest], [c for c in Piezo1Counts[med][1:]][idxtest], 'o',markersize=14)
    jj=jj+1
# plt.xlabel('Frecuencia (MHz)')
# plt.ylabel('Counts')
# plt.xlim(432, 446.5)
# plt.grid()
# plt.legend()
# #plt.title('Espectros para distintas geometrías')


Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
plt.figure()
plt.errorbar(np.arange(0,len(Intensity),1), [i/np.max(Intensity) for i in Intensity], yerr=[i/np.max(Intensity) for i in errorIntensity], fmt='o',capsize=2, markersize=8)
plt.errorbar(np.arange(0,len(Intensity),1), [p for p in pmdepthsdr], yerr=errorpmdepthsdr, fmt='o',capsize=2, markersize=8)
plt.xlabel('Ion position')
plt.ylabel('Intensity / DR Relative depth')
#plt.xticks([1,2,3,4,5])
#plt.xlim(200,3200)
plt.ylim(-0.1,1.1)
plt.grid()
#plt.axvline(3, color='salmon')
plt.legend()


#%%

import seaborn as sns
"""
Resonancias DD configuracion +2/-2 colineal variando la ubicacion del ion en los haces

Mismo analisis pero con fiteos
"""

palette = sns.color_palette("tab10")

pmlocmedvec = list(np.arange(0,12,1))+[13,12]+list(np.arange(15,len(PIEZO1_FILES),1))

"""
Hay que invertir la 12 con la 13, y la 14 es la misma que la 12, por las dudas
"""


#pmlocmedvec = [30]


plt.figure()

bkg = np.min(Piezo1Counts[1])

pmdepthsdr=[]
errorpmdepthsdr=[]

Intensity = []
errorIntensity = []

idxtest = 185
print(idxtest)

Gamas = []
ErrorGamas = []

jj=0
for med in pmlocmedvec:
    print(med)
    
    Freqs = [2*f*1e-6 for f in Piezo1Frequencies[med][1:]]
    Counts = [c for c in Piezo1Counts[med][1:]]

    if med==30:
        Freqs=Freqs[90:]
        Counts=Counts[90:]
        popt, pcov = curve_fit(Lorentzian, np.array(Freqs), np.array(Counts), p0=(-200,2100,435.8,0.05,0.1), bounds=((-10000,0,435.5,0,0),(0,1e4, 436.1, 1,1e4)))
    else:
        popt, pcov = curve_fit(Lorentzian, np.array(Freqs), np.array(Counts), p0=(-200,2100,435.8,0.05,0.1), bounds=((-10000,0,435.5,0,0),(0,1e4, 436.1, 1,1e4)))

    pmdepthsdr.append(1-(np.min(Lorentzian(np.array(Freqs),*popt))-bkg)/(popt[1]-bkg))
    #pmdepthsdr.append(1-(Piezo1Counts[med][1:][idxvecdr[jj]]-bkg)/(np.mean(Piezo1Counts[med][1:][0:20])-bkg))
    #errorpmdepthsdr.append(ErrorDRdepth(Piezo1Counts[med][1:][idxvecdr[jj]],np.mean(Piezo1Counts[med][1:][0:20]), bkg))
    errorpmdepthsdr.append(ErrorDRdepth(np.min(Lorentzian(Freqs,*popt)),popt[1], bkg))
    

    
    Intens = popt[1]
    Gamas.append(popt[3])
    ErrorGamas.append(np.sqrt(pcov[3,3]))
    Intensity.append(Intens)
    errorIntensity.append(2*np.sqrt(np.mean(Piezo1Counts[med][1:][0:20]))+np.sqrt(bkg))
        
    
        
    plt.plot([2*f*1e-6 for f in Piezo1Frequencies[med][1:]], [c for c in Piezo1Counts[med][1:]], '-o', markersize=2, alpha=0.7)
    plt.plot([f*1 for f in Freqs],Lorentzian(Freqs,*popt))

    #plt.plot([2*f*1e-6 for f in Piezo1Frequencies[med][1:]][idxtest], [c for c in Piezo1Counts[med][1:]][idxtest], 'o',markersize=14)
    jj=jj+1
# plt.xlabel('Frecuencia (MHz)')
# plt.ylabel('Counts')
# plt.xlim(432, 446.5)
# plt.grid()
# plt.legend()
# #plt.title('Espectros para distintas geometrías')


Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
plt.figure()
plt.errorbar(np.arange(0,len(Intensity),1), [i/np.max(Intensity) for i in Intensity], yerr=[i/np.max(Intensity) for i in errorIntensity], fmt='o',capsize=2, markersize=8)
plt.errorbar(np.arange(0,len(Intensity),1), [p for p in pmdepthsdr], yerr=errorpmdepthsdr, fmt='o',capsize=2, markersize=8)
plt.xlabel('Ion position')
plt.ylabel('Intensity / DR Relative depth')
#plt.xticks([1,2,3,4,5])
#plt.xlim(200,3200)
plt.ylim(-0.1,1.1)
plt.grid()
#plt.axvline(3, color='salmon')
plt.legend()

#%%

import seaborn as sns
"""
Resonancias DD configuracion +2/+2 colineal variando la ubicacion del ion en los haces

TODO EL OAM
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed

esto es hecho a mano... se puede ignorar
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
"""

palette = sns.color_palette("tab10")

mmlocmedvec = list(np.arange(0,10,1))+[11,10]+list(np.arange(12,len(PIEZO2_FILES),1))

"""
s
"""

#mmlocmedvec = [18]
#idxvecdr = [126,129,129,129,128,128,216,215,215,215,242,215,162,162,162,162,162,132,138,129,242,182,182,208,205,205,205,181,140,170,140,140,140,138,128,128,126,]
idxvecdr = [126,129,129,129,128,128,216,215,215,215,242,215,162,162,162,162,162,132,132,129,242,182,182,208,205,205,205,181,140,170,140,140,140,138,128,128,126,]
#idxtest = idxvecdr[18]

plt.figure()

bkg = np.min(Piezo2Counts[1])

mmdepthsdr=[]
errormmdepthsdr=[]

Intensity2 = []
errorIntensity2 = []


print(idxtest)

jj=0
for med in mmlocmedvec:
    print(med)
    
    mmdepthsdr.append(1-(Piezo2Counts[med][1:][idxvecdr[jj]]-bkg)/(np.mean(Piezo2Counts[med][1:][0:20])-bkg))
    errormmdepthsdr.append(ErrorDRdepth(Piezo2Counts[med][1:][idxvecdr[jj]],np.mean(Piezo2Counts[med][1:][0:20]), bkg))
        
    Intens = np.mean(Piezo2Counts[med][1:][0:20])-bkg
    
    Intensity2.append(Intens)
    errorIntensity2.append(2*np.sqrt(np.mean(Piezo2Counts[med][1:][0:20]))+np.sqrt(bkg))
      
    
        
    plt.plot([2*f*1e-6 for f in Piezo2Frequencies[med][1:]], [c for c in Piezo2Counts[med][1:]], '-o', markersize=2, alpha=0.7)
    plt.plot([2*f*1e-6 for f in Piezo2Frequencies[med][1:]][idxtest], [c for c in Piezo2Counts[med][1:]][idxtest], 'o',markersize=14)
    jj=jj+1
# plt.xlabel('Frecuencia (MHz)')
# plt.ylabel('Counts')
# plt.xlim(432, 446.5)
# plt.grid()
# plt.legend()
# #plt.title('Espectros para distintas geometrías')


plt.figure()
#plt.errorbar(np.arange(0,len(Intensity2),1), [i/np.max(Intensity2) for i in Intensity2], yerr=[i/np.max(Intensity2) for i in errorIntensity2], fmt='o',capsize=2, markersize=8)
plt.errorbar(np.arange(0,len(Intensity2),1), [p for p in mmdepthsdr], yerr=errormmdepthsdr, fmt='o',capsize=2, markersize=8)
plt.xlabel('Ion position')
plt.ylabel('Intensity / DR Relative depth')
#plt.xticks([1,2,3,4,5])
#plt.xlim(200,3200)
plt.ylim(-0.1,1.1)
plt.grid()
#plt.axvline(3, color='salmon')
plt.legend()


#%%

Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
import seaborn as sns
"""
Resonancias DD configuracion +2/+2 colineal variando la ubicacion del ion en los haces

TODO EL OAM

AHORA CON FITEOS LPM
"""

palette = sns.color_palette("tab10")

mmlocmedvec = list(np.arange(0,10,1))+[11,10]+list(np.arange(12,len(PIEZO2_FILES),1))

"""
s
"""
plotcurvita=1
#mmlocmedvec = [37]

plt.figure()

bkg = np.min(Piezo2Counts[1])

mmdepthsdr=[]
errormmdepthsdr=[]

Intensity2 = []
errorIntensity2 = []


print(idxtest)

jj=0
for med in mmlocmedvec:
    print(med)
    

    Freqs = [2*f*1e-6 for f in Piezo2Frequencies[med][1:]]
    Counts = [c for c in Piezo2Counts[med][1:]]


    if med==18:
        popt, pcov = curve_fit(Lorentzian, np.array(Freqs), np.array(Counts), p0=(-200,210,435.8,0.05,0.1), bounds=((-10000,0,435.5,0,0),(0,1e4, 436.1, 1,1e4)))
    
    else:
        popt, pcov = curve_fit(Lorentzian, np.array(Freqs), np.array(Counts), p0=(-200,2100,435.8,0.05,0.1), bounds=((-10000,0,435.5,0,0),(0,1e4, 436.1, 1,1e4)))

    mmdepthsdr.append(1-(np.min(Lorentzian(np.array(Freqs),*popt))-bkg)/(popt[1]-bkg))
    errormmdepthsdr.append(ErrorDRdepth(np.min(Lorentzian(Freqs,*popt)),popt[1], bkg))
    

    
    Intens = popt[1]
    Gamas.append(popt[3])
    ErrorGamas.append(np.sqrt(pcov[3,3]))
    Intensity2.append(Intens)
    errorIntensity.append(2*np.sqrt(np.mean(Piezo2Counts[med][1:][0:20]))+np.sqrt(bkg))
        
    
        
    plt.plot([2*f*1e-6 for f in Piezo2Frequencies[med][1:]], [c for c in Piezo2Counts[med][1:]], '-o', markersize=2, alpha=0.7)
    plt.plot([f*1 for f in Freqs],Lorentzian(Freqs,*popt))


    
        
    # plt.plot([2*f*1e-6 for f in Piezo2Frequencies[med][1:]], [c for c in Piezo2Counts[med][1:]], '-o', markersize=2, alpha=0.7)
    # plt.plot([2*f*1e-6 for f in Piezo2Frequencies[med][1:]][idxtest], [c for c in Piezo2Counts[med][1:]][idxtest], 'o',markersize=14)
    jj=jj+1
# plt.xlabel('Frecuencia (MHz)')
# plt.ylabel('Counts')
# plt.xlim(432, 446.5)
# plt.grid()
# plt.legend()
# #plt.title('Espectros para distintas geometrías')

if plotcurvita:
    plt.figure()
    #plt.errorbar(np.arange(0,len(Intensity2),1), [i/np.max(Intensity2) for i in Intensity2], yerr=[i/np.max(Intensity2) for i in errorIntensity2], fmt='o',capsize=2, markersize=8)
    plt.errorbar(np.arange(0,len(Intensity2),1), [p for p in mmdepthsdr], yerr=errormmdepthsdr, fmt='o',capsize=2, markersize=8)
    plt.xlabel('Ion position')
    plt.ylabel('Intensity / DR Relative depth')
    #plt.xticks([1,2,3,4,5])
    #plt.xlim(200,3200)
    plt.ylim(-0.1,1.1)
    plt.grid()
    #plt.axvline(3, color='salmon')
    plt.legend()

#%%
"""
graficos de:
    1) intensidad del haz
    2) profundidad de las resonancias en caso +2/+2 y en caso +2/-2
"""

import seaborn as sns
paletadenico=sns.color_palette('rocket')
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed

plt.figure()
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
plt.errorbar(np.arange(-16,18,1), [i/np.max(Intensity) for i in Intensity], yerr=[i/np.max(Intensity) for i in errorIntensity], fmt='o',capsize=2, markersize=8,color=paletadenico[1])
#plt.errorbar(np.arange(0,len(Intensity),1), [p for p in pmdepthsdr], yerr=errorpmdepthsdr, fmt='o',capsize=2, markersize=8)
#plt.errorbar(np.arange(0,len(Intensity),1), [p for p in mmdepthsdr[3:]], yerr=[0.5*m for m in errormmdepthsdr[3:]], fmt='o',capsize=2, markersize=8)
plt.xlabel('Distance to center (motor steps)', fontname='STIXgeneral',fontsize=15)
plt.ylabel('Beam intensity', fontname='STIXgeneral',fontsize=15)
plt.xticks([-15,-10,-5,0,5,10,15], fontname='STIXgeneral',fontsize=15)
plt.yticks([0,0.2,0.4,0.6,0.8,1], fontname='STIXgeneral',fontsize=15)
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
#plt.xlim(200,3200)
plt.ylim(-0.1,1.1)
plt.grid()
#plt.axvline(3, color='salmon')
plt.legend()
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
plt.tight_layout()
plt.savefig('/home/nico/Nextcloud/Nico/Doctorado/Charlas/2023 Europe/beamintensityradial.pdf')

Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed

Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
#%%
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed

Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
plt.figure()
#plt.errorbar(np.arange(0,len(Intensity),1), [i/np.max(Intensity) for i in Intensity], yerr=[i/np.max(Intensity) for i in errorIntensity], fmt='o',capsize=2, markersize=8)
plt.errorbar(np.arange(-16,18,1), [p for p in pmdepthsdr], yerr=errorpmdepthsdr, fmt='o',capsize=2, markersize=8,color=paletadenico[3])
plt.errorbar(np.arange(-16,18,1), [p for p in mmdepthsdr[2:-1]], yerr=[1*m for m in errormmdepthsdr[2:-1]], fmt='o',capsize=2, markersize=8,color=paletadenico[0])
plt.xlabel('Distance to center (motor steps)', fontname='STIXgeneral',fontsize=15)
plt.ylabel('DR Relative depth', fontname='STIXgeneral',fontsize=15)
plt.xticks([-15,-10,-5,0,5,10,15], fontname='STIXgeneral',fontsize=15)
plt.yticks([0,0.2,0.4,0.6,0.8,1], fontname='STIXgeneral',fontsize=15)
#plt.xlim(200,3200)
plt.ylim(-0.1,1.1)
plt.grid()
#plt.axvline(3, color='salmon')
plt.legend()
plt.tight_layout()
plt.savefig('/home/nico/Nextcloud/Nico/Doctorado/Charlas/2023 Europe/depthsradial.pdf')
Nicolas Nunez Barreto's avatar
Nicolas Nunez Barreto committed
#%%
"""
Intento ajustar un modelo para la profundidad
"""

def FunctionTarget(x,x0,A):
    return (A/(1+(x0/x)**2))**1

lim=10

xvec =  np.arange(0,len(Intensity),1)-15.5
xveclong = np.arange(np.min(xvec)-lim, np.max(xvec)+lim,0.01)
popt,pcov=curve_fit(FunctionTarget,xvec,pmdepthsdr)

plt.figure()
plt.errorbar(xvec, pmdepthsdr, yerr=errorpmdepthsdr, fmt='o',capsize=2, markersize=8)
#plt.plot(xveclong,FunctionTarget(xveclong,*popt))
plt.plot(xveclong,FunctionTarget(xveclong,5.93,0.83))
plt.xlabel('Ion position')
plt.ylabel('Intensity / DR Relative depth')
#plt.ylim(-0.1,1.1)
plt.grid()
#plt.axvline(3, color='salmon')
plt.legend()

print(popt)


#%%
"""
Ahora voy a intentar ajustarlas con una lorentziana que es mejor
"""
import seaborn as sns
"""
Resonancias DD configuracion +2/-2 colineal variando la ubicacion del ion en los haces

TODO EL OAM
"""
def Lorentzian( x, A, B, x0, gam ):
    return A * gam**2 / ( gam**2 + ( x - x0 )**2) + B


palette = sns.color_palette("tab10")

pmlocmedvec = list(np.arange(0,12,1))+[13,12]+list(np.arange(15,len(PIEZO1_FILES),1))

"""
Hay que invertir la 12 con la 13, y la 14 es la misma que la 12, por las dudas
"""


#pmlocmedvec = [35]


plt.figure()

bkg = np.min(Piezo1Counts[1])

pmdepthsdr=[]
errorpmdepthsdr=[]

Intensity = []
errorIntensity = []

idxtest = 185
print(idxtest)

jj=0
for med in pmlocmedvec:

    Freqs = [2*f*1e-6 for f in Piezo1Frequencies[med][1:]]
    Counts = [c for c in Piezo1Counts[med][1:]]

    if med==30:
        Freqs = Freqs[100:]
        Counts = Counts[100:]
        popt, pcov = curve_fit(Lorentzian, Freqs, Counts, p0=(-200,2100,435.8,0.05), bounds=((-10000,0,435.7,0),(0,1e4, 436.1, 1)))
        
    else:
        popt, pcov = curve_fit(Lorentzian, Freqs, Counts, p0=(-200,2100,435.8,0.05), bounds=((-10000,0,435.5,0),(0,1e4, 436.1, 1)))

    pmdepthsdr.append(1-(np.min(Lorentzian(Freqs,*popt))-bkg)/(popt[1]-bkg))
    errorpmdepthsdr.append(ErrorDRdepth(np.min(Lorentzian(Freqs,*popt)),popt[1], bkg))
        
    Intens = popt[1]
    
    Intensity.append(Intens)
    # errorIntensity.append(2*np.sqrt(np.mean(Piezo1Counts[med][1:][0:20]))+np.sqrt(bkg))
        
    
    if med in [8,21,25]:
        plt.plot([2*f*1e-6 for f in Piezo1Frequencies[med][1:]], [c for c in Piezo1Counts[med][1:]], '-o', markersize=2, alpha=0.7)
        plt.plot(Freqs,Lorentzian(Freqs,*popt))
    
    jj=jj+1
# plt.xlabel('Frecuencia (MHz)')
# plt.ylabel('Counts')
plt.xlim(435.2, 436.5)
plt.grid()
# plt.legend()
# #plt.title('Espectros para distintas geometrías')


plt.figure()
plt.plot(np.arange(0,len(Intensity),1), [i/np.max(Intensity) for i in Intensity], '-o',markersize=8)
plt.plot(np.arange(0,len(Intensity),1), [p for p in pmdepthsdr], 'o',markersize=8)
plt.xlabel('Ion position')
plt.ylabel('Intensity / DR Relative depth')
#plt.xticks([1,2,3,4,5])
#plt.xlim(200,3200)
plt.ylim(-0.1,1.1)
plt.grid()
#plt.axvline(3, color='salmon')
plt.legend()

#%%
"""
Intento ajustar un modelo para la profundidad
"""

def FunctionTarget(x,x0,A):
    return (A/(1+(x0/x)**2))**1

lim=10

xvec =  np.arange(0,len(Intensity),1)-15.5
xveclong = np.arange(np.min(xvec)-lim, np.max(xvec)+lim,0.01)
popt,pcov=curve_fit(FunctionTarget,xvec,pmdepthsdr)

plt.figure()
plt.errorbar(xvec, pmdepthsdr, yerr=errorpmdepthsdr, fmt='o',capsize=2, markersize=8)
#plt.plot(xveclong,FunctionTarget(xveclong,*popt))
plt.plot(xveclong,FunctionTarget(xveclong,5.93,0.83))
plt.xlabel('Ion position')
plt.ylabel('Intensity / DR Relative depth')
#plt.ylim(-0.1,1.1)
plt.grid()
#plt.axvline(3, color='salmon')
plt.legend()

print(popt)