Commit 7cb31688 authored by Matheus Müller's avatar Matheus Müller

Openflatus versao alpha

parent c71439ae
...@@ -2,7 +2,9 @@ ...@@ -2,7 +2,9 @@
import Tkinter as tk import Tkinter as tk
from PTE import * from PTE import *
from numpy import ones from pylab import *
from numpy import zeros
import os
Edict = dict() Edict = dict()
Ebuttons = list() Ebuttons = list()
...@@ -28,7 +30,7 @@ def elem_callback(x, edict , ebutton, Labelmass, LabelZ, Labelname, Labelsymbol, ...@@ -28,7 +30,7 @@ def elem_callback(x, edict , ebutton, Labelmass, LabelZ, Labelname, Labelsymbol,
############################################################################## ##############################################################################
# Element selection # Element selection
def elem_select(i, button, Labelmass, LabelZ, Labelname, Labelsymbol, Entrymass, EntryZ, Entryname, Entrysymbol): def elem_select(i, button, Labelmass, LabelZ, Labelname, Labelsymbol, LabelLS, Entrymass, EntryZ, Entryname, Entrysymbol,EntryLS):
global Econtrol, Edict global Econtrol, Edict
button[Econtrol]['relief']='raised' button[Econtrol]['relief']='raised'
...@@ -39,34 +41,33 @@ def elem_select(i, button, Labelmass, LabelZ, Labelname, Labelsymbol, Entrymass, ...@@ -39,34 +41,33 @@ def elem_select(i, button, Labelmass, LabelZ, Labelname, Labelsymbol, Entrymass,
LabelZ.configure(text = 'Z - %f' %Edict[i]['Z']) LabelZ.configure(text = 'Z - %f' %Edict[i]['Z'])
Labelname.configure(text = '%s' %Edict[i]['name']) Labelname.configure(text = '%s' %Edict[i]['name'])
Labelsymbol.configure(text = 'Symbol - %s' %Edict[i]['symbol']) Labelsymbol.configure(text = 'Symbol - %s' %Edict[i]['symbol'])
LabelLS.configure(text = 'Line Shape α - %f' %Edict[i]['LineShape'])
Entrymass.delete (0, len(Entrymass.get())) Entrymass.delete (0, len(Entrymass.get()))
Entryname.delete (0, len(Entryname.get())) Entryname.delete (0, len(Entryname.get()))
EntryZ.delete (0, len(EntryZ.get())) EntryZ.delete (0, len(EntryZ.get()))
Entrysymbol.delete (0, len(Entrysymbol.get())) Entrysymbol.delete (0, len(Entrysymbol.get()))
EntryLS.delete (0, len(EntryLS.get()))
Entrymass.insert(0, Edict[int(Econtrol)]['mass']) Entrymass.insert(0, Edict[int(Econtrol)]['mass'])
Entryname.insert(0, Edict[int(Econtrol)]['name']) Entryname.insert(0, Edict[int(Econtrol)]['name'])
EntryZ.insert(0, Edict[int(Econtrol)]['Z']) EntryZ.insert(0, int(Edict[int(Econtrol)]['Z']))
Entrysymbol.insert(0, Edict[int(Econtrol)]['symbol']) Entrysymbol.insert(0, Edict[int(Econtrol)]['symbol'])
EntryLS.insert(0, Edict[int(Econtrol)]['LineShape'])
#if Edict[i]['name'] == 'New':
# butPTE.invoke()
############################################################################## ##############################################################################
# Init # Init
def ewin_build(window, xmax, xstep): def ewin_build(window, xmax, xstep):
############################################################################## ##############################################################################
# Element addition # Element addition
def create(): def create():
global i, Edict, Ebuttons global i, Edict, Ebuttons
Edict[i] = dict(name = "New", symbol = "Hf", mass = 178.00, Z = 72, dist = ones(int(xmax/xstep)), control = 0, LineShape = 200.) Edict[i] = dict(name = "Element X", symbol = "Hf", mass = 178.00, Z = 72, dist = zeros(int(xmax/xstep)), control = 0, LineShape = 200.)
Ebuttons.insert(i, tk.Button(Eframe, text=i, width=1, height=1, command = lambda i=i : elem_select(i,Ebuttons,Labelmass,LabelZ,Labelname,Labelsymbol,Entrymass,EntryZ,Entryname,Entrysymbol))) Ebuttons.insert(i, tk.Button(Eframe, text=i, width=1, height=1, command = lambda i=i : elem_select(i,Ebuttons,Labelmass,LabelZ,Labelname,Labelsymbol,LabelLS,Entrymass,EntryZ,Entryname,Entrysymbol,EntryLS)))
Ebuttons[i].pack(side='left') Ebuttons[i].grid(column = ((len(Ebuttons)-1)%15), row = int(math.floor((len(Ebuttons)-1)/15)))
i=i+1 i=i+1
############################################################################## ##############################################################################
...@@ -86,16 +87,19 @@ def ewin_build(window, xmax, xstep): ...@@ -86,16 +87,19 @@ def ewin_build(window, xmax, xstep):
############################################################################## ##############################################################################
# Element Frames # Element Frames
Labelprop = tk.LabelFrame(window, text = 'Properties', relief='raised', bd=2)
Labelprop.pack(side = 'top')
Label1 = tk.LabelFrame(window, text = 'Elements', relief='raised', bd=2) Label1 = tk.LabelFrame(window, text = 'Elements', relief='raised', bd=2)
Label1.pack() Label1.pack(side = 'top')
Eframe = tk.Canvas(Label1) Eframe = tk.Canvas(Label1)
Eframe.pack(fill='both',expand=0) Eframe.pack(fill='both',expand=0)
############################################################################## ##############################################################################
# Init # Init
create() create()
############################################################################## ##############################################################################
# Properties # Properties
...@@ -109,14 +113,11 @@ def ewin_build(window, xmax, xstep): ...@@ -109,14 +113,11 @@ def ewin_build(window, xmax, xstep):
Labelname.configure(text = '%s' %Edict[i]['name']) Labelname.configure(text = '%s' %Edict[i]['name'])
Edict[i]['symbol'] = str(Entrysymbol.get()) Edict[i]['symbol'] = str(Entrysymbol.get())
Labelsymbol.configure(text = 'Symbol - %s' %Edict[i]['symbol']) Labelsymbol.configure(text = 'Symbol - %s' %Edict[i]['symbol'])
#Edict[i]['LineShape'] = float(EntryLS.get()) Edict[i]['LineShape'] = float(EntryLS.get())
#LabelLS.configure(text = 'Line Shape α - %f' %Edict[i]['LineShape']) LabelLS.configure(text = 'Line Shape α - %f' %Edict[i]['LineShape'])
Edict[i]['dist'] = np.loadtxt(Edict[i]['symbol']+".prof")
print Edict[i] print Edict[i]
Labelprop = tk.LabelFrame(window, text = 'Properties', relief='raised', bd=2)
#Labelprop.pack()
Labelprop.place(anchor='w', y=100, x=5)
LpFrame1 = tk.Frame(Labelprop) LpFrame1 = tk.Frame(Labelprop)
LpFrame1.pack(side='left') LpFrame1.pack(side='left')
LpFrame2 = tk.Frame(Labelprop) LpFrame2 = tk.Frame(Labelprop)
...@@ -134,32 +135,27 @@ def ewin_build(window, xmax, xstep): ...@@ -134,32 +135,27 @@ def ewin_build(window, xmax, xstep):
EPframee2 = tk.Frame(LpFrame2) EPframee2 = tk.Frame(LpFrame2)
EPframee2.pack(side='right') EPframee2.pack(side='right')
Labelname = tk.Label(EPframel, width=20, pady=2, text='%s' %Edict[int(Econtrol)]['name']) Labelname = tk.Label(EPframel, width=26, pady=2, text='%s' %Edict[int(Econtrol)]['name'])
Labelmass = tk.Label(EPframel, width=20, pady=2, text='Mass - %f' %Edict[int(Econtrol)]['mass']) Labelmass = tk.Label(EPframel, width=26, pady=2, text='Mass - %f' %Edict[int(Econtrol)]['mass'])
LabelZ = tk.Label(EPframel, width=20, pady=2, text='Z - %f' %Edict[int(Econtrol)]['Z']) LabelZ = tk.Label(EPframel, width=26, pady=2, text='Z - %f' %Edict[int(Econtrol)]['Z'])
Labelsymbol = tk.Label(EPframel, width=20, pady=2, text= 'Symbol - %s' %Edict[int(Econtrol)]['symbol']) Labelsymbol = tk.Label(EPframel, width=26, pady=2, text= 'Symbol - %s' %Edict[int(Econtrol)]['symbol'])
#LabelLS = tk.Label(EPframel, width=30, pady=2, text='Line Shape α - %f' %Edict[int(Econtrol)]['LineShape']) LabelLS = tk.Label(EPframel, width=26, pady=2, text='Line Shape α - %f' %Edict[int(Econtrol)]['LineShape'])
Labelname.pack() Labelname.pack()
Labelmass.pack() Labelmass.pack()
LabelZ.pack() LabelZ.pack()
Labelsymbol.pack() Labelsymbol.pack()
#LabelLS.pack() LabelLS.pack()
Entryname = tk.Entry(EPframee, width=15) Entryname = tk.Entry(EPframee, width=18)
Entrymass = tk.Entry(EPframee, width=15) Entrymass = tk.Entry(EPframee, width=18)
EntryZ = tk.Entry(EPframee, width=15) EntryZ = tk.Entry(EPframee, width=18)
Entrysymbol = tk.Entry(EPframee, width=15) Entrysymbol = tk.Entry(EPframee, width=18)
#EntryLS = tk.Entry(EPframee, width=22) EntryLS = tk.Entry(EPframee, width=18)
Entryname.pack() Entryname.pack()
Entrymass.pack() Entrymass.pack()
EntryZ.pack() EntryZ.pack()
Entrysymbol.pack() Entrysymbol.pack()
#EntryLS.pack() EntryLS.pack()
BUpdate = tk.Button(EPframee2, text ='Update properties', command=lambda i=i :update(Econtrol), width=15)
BUpdate.pack()#side='right')
BUpdateD = tk.Button(EPframee2, text ='Update Distribution', command=lambda i=i :update(Econtrol), width=15)
BUpdateD.pack()#side='left')
but = tk.Button(EPframel2, command=create, text='Add element', bd=1, height=1,width=15) but = tk.Button(EPframel2, command=create, text='Add element', bd=1, height=1,width=15)
but.pack()#side='left') but.pack()#side='left')
...@@ -170,55 +166,15 @@ def ewin_build(window, xmax, xstep): ...@@ -170,55 +166,15 @@ def ewin_build(window, xmax, xstep):
butPTE = tk.Button(EPframel2, text='Show PTE', bd=1, width=15, command=lambda i=i :elem_callback(Econtrol,Edict,Ebuttons,Labelmass,LabelZ,Labelname,Labelsymbol,Entrymass,EntryZ,Entryname,Entrysymbol)) butPTE = tk.Button(EPframel2, text='Show PTE', bd=1, width=15, command=lambda i=i :elem_callback(Econtrol,Edict,Ebuttons,Labelmass,LabelZ,Labelname,Labelsymbol,Entrymass,EntryZ,Entryname,Entrysymbol))
butPTE.pack()#side='left') butPTE.pack()#side='left')
############################################################################## BUpdate = tk.Button(EPframel2, text ='Update properties', command=lambda i=int(Econtrol) :update(Econtrol), width=15)
# Init BUpdate.pack()#side='right')
Entrymass.insert(0, Edict[int(Econtrol)]['mass'])
Entryname.insert(0, Edict[int(Econtrol)]['name'])
EntryZ.insert(0, Edict[int(Econtrol)]['Z'])
Entrysymbol.insert(0, Edict[int(Econtrol)]['symbol'])
#EntryLS.insert(0, Edict[int(Econtrol)]['LineShape'])
############################################################################## BUpdateD = tk.Button(EPframel2, text ='Distribution', command =lambda i=int(Econtrol) : os.system("python profiler.py %s &" %Edict[int(Econtrol)]['symbol']), width=15)
# Canvas callback BUpdateD.pack()#side='left')
def point(event):
Distcanvas.create_oval(event.x, event.y, event.x+1, event.y+1, fill="black")
points.append(event.x)
points.append(event.y)
return points
def canxy(event):
print event.x, event.y
def graph(event):
global theline
Distcanvas.create_line(points, tags="theline")
def toggle(event):
global spline
if spline == 0:
Distcanvas.itemconfigure(tag1, smooth=1)
spline = 1
elif spline == 1:
Distcanvas.itemconfigure(tag1, smooth=0)
spline = 0
print event.x, event.y
return spline
############################################################################## ##############################################################################
# Drawing Canvas # Init
Ebuttons[0].invoke()
Label2 = tk.LabelFrame(window, text = 'Distribution', bd=1)
#Label2.pack()
Label2.place(anchor='w', y=320, x=5)
Distcanvas = tk.Canvas(Label2, bg="white", width=800, height= 300)
#Distcanvas.configure(cursor="crosshair")
Distcanvas.grid_columnconfigure(999, weight = 1)
Distcanvas.grid_rowconfigure(999, weight = 1)
Distcanvas.pack()
Distcanvas.bind("<Button-1>", point)
Distcanvas.bind("<Button-3>", graph)
Distcanvas.bind("<Button-2>", toggle)
############################################################################## ##############################################################################
...@@ -7,7 +7,7 @@ import sys ...@@ -7,7 +7,7 @@ import sys
############################################################################## ##############################################################################
# Initial parameters # Initial parameters
param = dict(dedx=192. ,Theta_out=70., dW2dx=20740., E0= 100000., Theta_in=0.) #, FWHM0=180.) param = dict(dedx=192. ,Theta_out=70., dW2dx=20740., E0= 100000., Theta_in=0., FWHM0=180.)
ionb = dict(Z=1, mass=1.0079) ionb = dict(Z=1, mass=1.0079)
############################################################################## ##############################################################################
...@@ -21,9 +21,9 @@ def init_calc(): ...@@ -21,9 +21,9 @@ def init_calc():
param['Theta_out'] = float(EntryT_out.get()) param['Theta_out'] = float(EntryT_out.get())
param['Theta_in'] = float(EntryT_in.get()) param['Theta_in'] = float(EntryT_in.get())
param['E0'] = float(EntryE0.get()) param['E0'] = float(EntryE0.get())
#param['FWHM0'] = float(EntryFWHM0.get()) param['FWHM0'] = float(EntryFWHM0.get())
ionb['mass'] = float(Entryionmass.get()) ionb['mass'] = float(Entryionmass.get())
ionb['Z'] = float(EntryionZ.get()) ionb['Z'] = float(EntryionZ.get())
spectro(param, str(modelvar.get()), ionb, float(EntryEmin.get()),float(EntryEmax.get()), float(EntryEstep.get()), float(EntryDstep.get()), int(LSvar.get()) ) spectro(param, str(modelvar.get()), ionb, float(EntryEmin.get()),float(EntryEmax.get()), float(EntryEstep.get()), float(EntryDstep.get()), int(LSvar.get()) )
...@@ -31,29 +31,27 @@ def init_calc(): ...@@ -31,29 +31,27 @@ def init_calc():
# Main window # Main window
root = tk.Tk() root = tk.Tk()
root.minsize(375,270) root.minsize(730,460)
root.title('Open Flatus') root.title('Open Flatus')
root.geometry('375x270+200+400') root.geometry('375x270+200+400')
##############################################################################
# Elements window
ewin = tk.Tk()
ewin.minsize(810,550)
ewin.title('Elements')
ewin.geometry('600x400+600+400')
############################################################################## ##############################################################################
# Frames # Frames
MasterFrame1 = tk.Frame(root) MasterFrame1 = tk.Frame(root)
MasterFrame1.pack(side='left') MasterFrame1.pack(side='left')
MasterFrame2 = tk.Frame(root)
MasterFrame2.pack(side='top')
Label2 = tk.LabelFrame(MasterFrame1, text = 'Parameters', relief='raised', bd=2) Label2 = tk.LabelFrame(MasterFrame1, text = 'Parameters', relief='raised', bd=2)
Label2.pack(side='top') Label2.pack(side='top')
Label3 = tk.LabelFrame(MasterFrame1, text = 'Ion parameters', relief='raised', bd=2) Label3 = tk.LabelFrame(MasterFrame1, text = 'Ion parameters', relief='raised', bd=2)
Label3.pack(side='left') Label3.pack(side='top')
Label4 = tk.LabelFrame(root, text = 'Energy loss model', relief='raised', bd=2) Label4 = tk.LabelFrame(MasterFrame1, text = 'Energy loss model', relief='raised', bd=2)
Label4.pack(side='right') Label4.pack(side='top')
# Elements frame
Label1 = tk.Frame(MasterFrame2)
Label1.pack(side='top')
Pframel = tk.Frame(Label2) Pframel = tk.Frame(Label2)
Pframel.pack(side='left') Pframel.pack(side='left')
...@@ -72,7 +70,7 @@ Labeldedx = tk.Label(Pframel, width=9, text='dƐ/dx', pady=2) ...@@ -72,7 +70,7 @@ Labeldedx = tk.Label(Pframel, width=9, text='dƐ/dx', pady=2)
LabeldW2dx = tk.Label(Pframel, width=9, text='dω²/dx', pady=2) LabeldW2dx = tk.Label(Pframel, width=9, text='dω²/dx', pady=2)
LabelT_out = tk.Label(Pframel, width=9, text='θ out', pady=2) LabelT_out = tk.Label(Pframel, width=9, text='θ out', pady=2)
LabelT_in = tk.Label(Pframel, width=9, text='θ in', pady=2) LabelT_in = tk.Label(Pframel, width=9, text='θ in', pady=2)
#LabelFWHM0 = tk.Label(Pframel, width=9, text='FWHM0') LabelFWHM0 = tk.Label(Pframel, width=9, text='FWHM0')
LabelEmin = tk.Label(Pframel, width=9, text='Emin', pady=2) LabelEmin = tk.Label(Pframel, width=9, text='Emin', pady=2)
LabelEmax = tk.Label(Pframel, width=9, text='Emax', pady=2) LabelEmax = tk.Label(Pframel, width=9, text='Emax', pady=2)
LabelEstep = tk.Label(Pframel, width=9, text='Estep', pady=2) LabelEstep = tk.Label(Pframel, width=9, text='Estep', pady=2)
...@@ -81,7 +79,7 @@ Labeldedx.pack() ...@@ -81,7 +79,7 @@ Labeldedx.pack()
LabeldW2dx.pack() LabeldW2dx.pack()
LabelT_out.pack() LabelT_out.pack()
LabelT_in.pack() LabelT_in.pack()
#LabelFWHM0.pack() LabelFWHM0.pack()
LabelEmin.pack() LabelEmin.pack()
LabelEmax.pack() LabelEmax.pack()
LabelEstep.pack() LabelEstep.pack()
...@@ -92,7 +90,7 @@ Entrydedx = tk.Entry(Pframee, width=11) ...@@ -92,7 +90,7 @@ Entrydedx = tk.Entry(Pframee, width=11)
EntrydW2dx = tk.Entry(Pframee, width=11) EntrydW2dx = tk.Entry(Pframee, width=11)
EntryT_out = tk.Entry(Pframee, width=11) EntryT_out = tk.Entry(Pframee, width=11)
EntryT_in = tk.Entry(Pframee, width=11) EntryT_in = tk.Entry(Pframee, width=11)
#EntryFWHM0 = tk.Entry(Pframee, width=11) EntryFWHM0 = tk.Entry(Pframee, width=11)
EntryEmin = tk.Entry(Pframee, width=11) EntryEmin = tk.Entry(Pframee, width=11)
EntryEmax = tk.Entry(Pframee, width=11) EntryEmax = tk.Entry(Pframee, width=11)
EntryEstep = tk.Entry(Pframee, width=11) EntryEstep = tk.Entry(Pframee, width=11)
...@@ -101,7 +99,7 @@ Entrydedx.pack() ...@@ -101,7 +99,7 @@ Entrydedx.pack()
EntrydW2dx.pack() EntrydW2dx.pack()
EntryT_out.pack() EntryT_out.pack()
EntryT_in.pack() EntryT_in.pack()
#EntryFWHM0.pack() EntryFWHM0.pack()
EntryEmin.pack() EntryEmin.pack()
EntryEmax.pack() EntryEmax.pack()
EntryEstep.pack() EntryEstep.pack()
...@@ -127,7 +125,7 @@ EntrydW2dx.insert(0,param['dW2dx']) ...@@ -127,7 +125,7 @@ EntrydW2dx.insert(0,param['dW2dx'])
EntryT_out.insert(0,param['Theta_out']) EntryT_out.insert(0,param['Theta_out'])
EntryT_in.insert(0,param['Theta_in']) EntryT_in.insert(0,param['Theta_in'])
EntryE0.insert(0,param['E0']) EntryE0.insert(0,param['E0'])
#EntryFWHM0.insert(0,param['FWHM0']) EntryFWHM0.insert(0,param['FWHM0'])
EntryEmin.insert(0,70000) EntryEmin.insert(0,70000)
EntryEmax.insert(0,100000) EntryEmax.insert(0,100000)
EntryEstep.insert(0,20) EntryEstep.insert(0,20)
...@@ -155,22 +153,22 @@ C1.select() ...@@ -155,22 +153,22 @@ C1.select()
############################################################################## ##############################################################################
CalcButton = tk.Button(Label4, text = 'Calculate', command = init_calc, bd=4, width=20, height=4) CalcButton = tk.Button(Label4, text = 'Calculate', command = init_calc, bd=4, width=17, height=2)
CalcButton.pack() CalcButton.pack()
ExitButton = tk.Button(Label4, text = 'Exit', command = sys.exit, bd=4, width=20, height=2) ExitButton = tk.Button(Label4, text = 'Exit', command = sys.exit, bd=4, width=17, height=2)
ExitButton.pack(side='bottom') ExitButton.pack(side='bottom')
############################################################################## ##############################################################################
# Elements callback # Elements callback
ewin_build(ewin, 15, float(EntryDstep.get())) ewin_build(Label1, 15, float(EntryDstep.get()))
############################################################################## ##############################################################################
spectro(param, 'gaussiana', ionb, float(EntryEmin.get()), float(EntryEmax.get()), float(EntryEstep.get()), float(EntryDstep.get()), int(LSvar.get()) ) spectro(param, 'gaussiana', ionb, float(EntryEmin.get()), float(EntryEmax.get()), float(EntryEstep.get()), float(EntryDstep.get()), int(LSvar.get()) )
root.mainloop() root.mainloop()
ewin.mainloop() Label1.mainloop()
############################################################################## ##############################################################################
...@@ -39,51 +39,29 @@ def convoluiF0(espectro, alpha,passoE): ...@@ -39,51 +39,29 @@ def convoluiF0(espectro, alpha,passoE):
######################################################### #########################################################
######## Função de convolução com Gaussiana #######
def convoluiGauss(espectro, Sigma,passoE):
areaantes=espectro.sum()
temp=zeros(len(espectro),float)
wg4_step=4.0*Sigma/passoE;
wg2=Sigma**2
eixo=arange(-6*Sigma,6*Sigma,passoE)
eixo=arange(0,len(temp)*passoE,passoE)
print len(eixo)
for i in range(0,len(eixo)):
for j in range(0,i):
temp[i]+=espectro[j]*exp(-1.0*pow(passoE*(i-j-wg4_step),2.0)/(wg2*2.0));
return temp*areaantes/temp.sum(), wg4_step*passoE
# return temp/temp.max(), wg4_step*passoE # Incluido para teste com Gaussiana -> Altura = 1
#########################################################
#### Parâmetros iniciais #### #### Parâmetros iniciais ####
# Feixe: # Feixe:
E0=99000 #eV E0=100000 #eV
mi=1. mi=1.
Zi=1. Zi=1.
# Alvo: # Alvo:
dedx=215. # eV/nm dedx=192. # eV/nm
dw2dx=20000. # eV^2/nm dw2dx=20740. # eV^2/nm
Theta_in=7. Theta_in=0.
Theta_out=58.67 Theta_out=70.
Theta_s = PI - (Theta_in + Theta_out)*PI/180 Theta_s = PI - (Theta_in + Theta_out)*PI/180
espessura=19 espessura=2
# Cálculo: # Cálculo:
deltaE=20000 deltaE=20000
passoE=20 passoE=20
passox=0.01 passox=0.02
e = arange(E0-deltaE, E0, passoE) e = arange(E0-deltaE, E0, passoE)
totalbessel=e*0 totalbessel=e*0
totalbessel_ls=e*0 totalbessel_ls=e*0
totalgaussiana=e*0 totalgaussiana=e*0
resolucao = 250 # Largura a meia algura da Gaussiana que representa a resolução experimental (FWHM) em eV
## TODO: - Otimizar os cálculos de cada componente: evitar cálculos (de zeros) em energias desnecessários ## TODO: - Otimizar os cálculos de cada componente: evitar cálculos (de zeros) em energias desnecessários
## - Definir o intervalo de plotagem do espectro independente do deltaE ## - Definir o intervalo de plotagem do espectro independente do deltaE
...@@ -91,9 +69,7 @@ resolucao = 250 # Largura a meia algura da Gaussiana que representa a resoluçã ...@@ -91,9 +69,7 @@ resolucao = 250 # Largura a meia algura da Gaussiana que representa a resoluçã
## Define a composição da amostra ## ## Define a composição da amostra ##
elementos=[[57.,138.9,192.,"la.prof"],[38.,87.62,165.,"sr.prof"],[22.,47.867,142.,"sr.prof"]] # La, Sr,ti elementos=[[72.,178.49,217.],[14.,28.,121.],[8.,16.,102]] # Hf, Si e O
#elementos=[[72.,178.49,217.,"hf.txt"],[14.,28.,121.,"si.txt"],[8.,16.,102,"o.txt"]] # Hf, Si e O
#elementos=[[72.,178.49,17.],[14.,28.,21.],[8.,16.,12]] # Hf, Si e O
#elementos=[[72.,178.49,217.]]#,[14.,28.,121.],[8.,16.,102]] # Hf, Si e O #elementos=[[72.,178.49,217.]]#,[14.,28.,121.],[8.,16.,102]] # Hf, Si e O
# Esta lista deve incluir a distribuição # Esta lista deve incluir a distribuição
...@@ -106,8 +82,6 @@ for componente in elementos: ...@@ -106,8 +82,6 @@ for componente in elementos:
mt=componente[1]; mt=componente[1];
Zt=componente[0]; Zt=componente[0];
alpha_lineshape=componente[2] alpha_lineshape=componente[2]
perfil = loadtxt(componente[3], unpack=True)
# plot(perfil)
k = ((sqrt (mt**2+mi**2*(sin(Theta_s)**2)) + mi*cos(Theta_s) ) / (mi+mt) )**2 k = ((sqrt (mt**2+mi**2*(sin(Theta_s)**2)) + mi*cos(Theta_s) ) / (mi+mt) )**2
dedxef=dedx*(k/cos(Theta_in*PI/180.) + 1/cos(Theta_out*PI/180.)) dedxef=dedx*(k/cos(Theta_in*PI/180.) + 1/cos(Theta_out*PI/180.))
...@@ -134,7 +108,7 @@ for componente in elementos: ...@@ -134,7 +108,7 @@ for componente in elementos:
# #
# CSRf(Z1, Z2, Ein, M1, M2, Theta): # CSRf(Z1, Z2, Ein, M1, M2, Theta):
secchoque=CSRf(Zi, Zt,k*E0-dedxef*x,mi,mt, Theta_s)*passox*perfil[int(x/passox)] secchoque=CSRf(Zi, Zt,k*E0-dedxef*x,mi,mt, Theta_s)*passox
# É multiplicada pelo passo para manter a escala do gráfico # É multiplicada pelo passo para manter a escala do gráfico
# independente do mesmo. # independente do mesmo.
...@@ -158,7 +132,7 @@ for componente in elementos: ...@@ -158,7 +132,7 @@ for componente in elementos:
# Gaussiana: # Gaussiana:
Em=k*E0-x*dedxef Em=k*E0-x*dedxef
sigma=(sqrt(x*dw2dxef)+resolucao/2.35482)**2 sigma=x*dw2dxef
gaussianacamada = exp((-(e-Em)**2)/(2.*sigma))/sqrt(2.*PI*sigma) gaussianacamada = exp((-(e-Em)**2)/(2.*sigma))/sqrt(2.*PI*sigma)
sinalgaussiana=gaussianacamada*secchoque+sinalgaussiana sinalgaussiana=gaussianacamada*secchoque+sinalgaussiana
#plt.plot(e,gaussianacamada) #plt.plot(e,gaussianacamada)
...@@ -172,56 +146,19 @@ for componente in elementos: ...@@ -172,56 +146,19 @@ for componente in elementos:
totalbessel=sinalbessel+totalbessel totalbessel=sinalbessel+totalbessel
totalbessel_ls=sinalbessel_ls+totalbessel_ls totalbessel_ls=sinalbessel_ls+totalbessel_ls
totalbessel_ls_gauss,gshift=convoluiGauss(totalbessel_ls,resolucao/2.35482,passoE)
totalgaussiana=sinalgaussiana+totalgaussiana totalgaussiana=sinalgaussiana+totalgaussiana
##################################### #####################################
# Mostra o Gráfico # Mostra o Gráfico
#bes=plt.plot(e,totalbessel,lw=2,color="red") bes=plt.plot(e,totalbessel,lw=2,color="red")
#besls=plt.plot(e,totalbessel_ls,lw=2,color="blue") besls=plt.plot(e,totalbessel_ls,lw=2,color="blue")
beslsgaus=plt.plot(e-gshift,totalbessel_ls_gauss,lw=2,color="orange")
gaus=plt.plot(e,totalgaussiana,lw=2,color="green") gaus=plt.plot(e,totalgaussiana,lw=2,color="green")
#plt.plot(e,totalbessel*10,lw=2, color="blue") #plt.plot(e,totalbessel*10,lw=2, color="blue")
#plt.plot(e,totalgaussiana*10,lw=2, color="green") #plt.plot(e,totalgaussiana*10,lw=2, color="green")
#plt.legend([bes,besls,beslsgaus,gaus],["Bessel","Bessel+lineshape","Bessel+ls+expres","Gaussiana"],loc= "upper left", shadow=True) plt.legend([bes,besls,gaus],["Bessel","Bessel+lineshape","Gaussiana"],loc= "upper left", shadow=True)
plt.legend([beslsgaus,gaus],["Bessel+ls+expres","Gaussiana"],loc= "upper left", shadow=True)
#plt.legend([bes,besls,gaus],["Bessel","Bessel+lineshape","Gaussiana"],loc= "upper left", shadow=True)
plt.xlabel("Energia dos Ions (eV)") plt.xlabel("Energia dos Ions (eV)")
plt.ylabel("Rendimento (contagens)") plt.ylabel("Rendimento (contagens)")
#print totalbessel
################################################
# Testando convolução com deltas de Dirac
#figure()
#xxy=e*0
#xxy[00]=1. # Correspondente a uma delta em 80000
#xxy[100]=1. # Correspondente a uma delta em 82000
#xxy[300]=1.
#xxy[400]=1.
#xxy[500]=1.
# Convoluindo as deltas com uma gaussiana de FWHM = 550 eV
#xxyz,g_shift=convoluiGauss(xxy,resolucao/2.35482,passoE)
# Convoluindo a convolução novamente -> FWHM Res = 777.81 eV = sqrt(550^2+550^2)
#xxyz2,g_shift2=convoluiGauss(xxyz,resolucao/2.35482,passoE)
#deltas = plot(e,xxy,lw=2)
#convo1 = plot(e-g_shift,xxyz+0.001,lw=2)
#convo2 = plot(e-g_shift-g_shift2,xxyz2+0.002,lw=2)
#convo2_analitica=plot(e,exp(-(e-82000)**2.0/(2.*(777.81/2.35482)**2))+0.002,lw=2)
#plt.legend([deltas,convo1,convo2,convo2_analitica],["Deltas","Conv Gauss 1x 550 eV","Conv Gauss 2x 550 eV","Conv Analitica (777.81)"],loc= "upper right", shadow=True)
#print e-g_shift, xxyz
#print str(g_shift) + " <- g_shift"
#print len(e)
#print len(xxy)
#print len(xxyz)
show() show()
################## ##################
...@@ -44,7 +44,7 @@ def spectro(p, modelo, ion, Emin, Emax, EPasso, passox, LScontrol): ...@@ -44,7 +44,7 @@ def spectro(p, modelo, ion, Emin, Emax, EPasso, passox, LScontrol):
e = arange(Emin, Emax, EPasso) e = arange(Emin, Emax, EPasso)
Y = e*0. Y = e*0.
dose = 10e32 dose = 2e34
Theta_s = PI - (p['Theta_in'] + p['Theta_out'])*PI/180 Theta_s = PI - (p['Theta_in'] + p['Theta_out'])*PI/180
mi = ion['mass'] mi = ion['mass']
...@@ -108,8 +108,8 @@ def spectro(p, modelo, ion, Emin, Emax, EPasso, passox, LScontrol): ...@@ -108,8 +108,8 @@ def spectro(p, modelo, ion, Emin, Emax, EPasso, passox, LScontrol):
xlabel("Enegy (eV)") xlabel("Enegy (eV)")
ylabel("Yield") ylabel("Yield")
amostra = np.loadtxt('gustavson.xy') amostra = np.loadtxt('H2.py')
xxx= amostra[:,0] xxx= amostra[:,0]*1000
ccc=amostra[:,1] ccc=amostra[:,1]
plot (xxx, ccc, 'r.') plot (xxx, ccc, 'r.')
......
1 1
2 0.5
3 1
4.3 0.7
4.5 1
5 0.9
6 .1
# draggable rectangle with the animation blit techniques; see
# http://www.scipy.org/Cookbook/Matplotlib/Animations
import numpy as np import numpy as np
import matplotlib.pyplot as plt import matplotlib.pyplot as plt
from scipy.interpolate import interp1d
from scipy import interpolate
from pylab import * from pylab import *
from concentracao import * from concentracao import *
import sys import sys
import os
listay = list() listay = list()
listax = list() listax = list()
...@@ -21,8 +18,9 @@ dxstep = 0.01 ...@@ -21,8 +18,9 @@ dxstep = 0.01
profundidademax = 20.