Skip to content

O
OpenFlatus
Commit 369633e2 authored by Matheus Müller's avatar Matheus Müller
Browse files
Options

Novas funçoes e melhorias no metodo de produzir o grafico do sinal

parent 08454278
Expand all Hide whitespace changes
Inline Side-by-side
Showing with 2676 additions and 2531 deletions
Ewindow.py
View file @ 369633e2
......@@ -168,10 +168,10 @@ def ewin_build(window, OFlabentries, mainwindow, calcbutton, distvar, modevar, L
global i, Econtrol
arquivo = tkFileDialog.askopenfile()
for n in arange(9,18):
OFlabentries[n].delete(0,len(OFlabentries[n].get()))
OFlabentries[n].delete(0,'end')
OFlabentries[n].insert(0,float(arquivo.readline()))
for n in arange(21,24):
OFlabentries[n].delete(0,len(OFlabentries[n].get()))
for n in arange(21,25):
OFlabentries[n].delete(0,'end')
OFlabentries[n].insert(0,float(arquivo.readline()))
f = loadtxt(arquivo.name + '.elm', dtype="string")
els = f[:, 0]
......@@ -205,10 +205,10 @@ def ewin_build(window, OFlabentries, mainwindow, calcbutton, distvar, modevar, L
elist = open(config.name + '.elm', 'w+')
for k in arange(9,18):
config.write( (OFlabentries[k].get())+'\n' )
for j in arange(21,24):
for j in arange(21,25):
config.write( (OFlabentries[j].get())+'\n' )
config.write( '##############################################################################\n' )
config.write( '# dƐ/dx\n# dω²/dx\n# θ out\n# θ in\n# FWHM\n# E min\n# E max\n# E step\n# Depth step\n# Ion energy\n# Ion mass\n# Ion Z\n' )
config.write( '# dƐ/dx\n# dω²/dx\n# θ out\n# θ in\n# FWHM\n# E min\n# E max\n# E step\n# Depth step\n# Ion energy\n# Ion mass\n# Ion Z\n # Dose\n' )
for n in range(len(Edict)):
elist.write(Edict[n]['symbol']+' ')
elist.write(Edict[n]['name']+' ')
......
Ewindow.pyc
View file @ 369633e2
No preview for this file type
Hf.dat
View file @ 369633e2
0 0
1.9523573201 0.00520408163265
3.10669975186 0.345255102041
5.20794044664 0.0691326530602
2.84863523573 0.842704081633
7.2129032258 0.932270408162
10.0 0
Hf.prof
View file @ 369633e2
This diff is collapsed.
OpenFlatus.py
View file @ 369633e2
......@@ -3,6 +3,14 @@
import Tkinter as tk
from Ewindow import *
from numerical import *
import tkFileDialog
fig = plt.figure()
ax = fig.add_subplot(111)
fig.canvas.set_window_title("Signal")
xx = arange(0,1, 0.1)
temp, = plot(xx,xx*0)
doseold = 12e34
##############################################################################
# Arbitrary parameters to declare dictionaries
......@@ -10,6 +18,27 @@ param = dict(dedx=192. ,Theta_out=70., dW2dx=20740., E0= 100000., Theta_in=0., F
ionb = dict(Z=1, mass=1.0079)
OFlabentries = list()
##############################################################################
# Load experiment data
def loadexp():
arquivo = tkFileDialog.askopenfilename()
amostra = np.loadtxt(arquivo)
xxx= amostra[:,0]*1000
ccc=amostra[:,1]
plot (xxx, ccc, 'r.')
show()
##############################################################################
# Adjust graphic
def adjust():
temp.set_xdata(arange(OFlabentries[14].get(),OFlabentries[15].get(),OFlabentries[16].get()))
temp.set_ydata(temp.ydata/doseold*float(OFlabentries[24].get()))
doseold = float(OFlabentries[24].get())
show()
##############################################################################
# Calculation callback
# Get parameters
......@@ -29,18 +58,18 @@ def init_calc():
if str(distvar.get()) == 'drawing':
for i in range(len(Edict)):
Edict[i]['dist'] = np.loadtxt(Edict[i]['symbol']+".prof")
spectro(param, str(modelvar.get()), ionb, OFlabentries, int(LSvar.get()), int(CGvar.get()))
spectro(param, str(modelvar.get()), ionb, OFlabentries, int(LSvar.get()), int(CGvar.get()), temp)
elif str(distvar.get()) == 'layer':
spectrolayer(param, str(modelvar.get()), ionb, OFlabentries, int(LSvar.get()), int(CGvar.get()))
spectrolayer(param, str(modelvar.get()), ionb, OFlabentries, int(LSvar.get()), int(CGvar.get()), temp)
elif str(modevar.get()) == 'RRNA':
return 0
spectroRNA(param, str(modelvar.get()), ionb, OFlabentries, int(LSvar.get()), int(CGvar.get()), temp)
##############################################################################
# Main window
root = tk.Tk()
root.minsize(670,550)
root.minsize(670,610)
root.title('Open Flatus')
root.geometry('670x500+200+400')
......@@ -57,6 +86,8 @@ Label2 = tk.LabelFrame(MasterFrame1, text = 'Parameters', relief='raised', bd=2)
Label2.pack(side='top')
Label3 = tk.LabelFrame(MasterFrame1, text = 'Ion parameters', relief='raised', bd=2)
Label3.pack(side='top')
Label5 = tk.LabelFrame(MasterFrame1, text = 'Graphical parameters', relief='raised', bd=2)
Label5.pack(side='top')
Label4 = tk.LabelFrame(MasterFrame1, text = 'Energy loss model', relief='raised', bd=2)
Label4.pack(side='top')
......@@ -69,6 +100,11 @@ Pframel.pack(side='left')
Pframee = tk.Frame(Label2)
Pframee.pack(side='right')
Gframel = tk.Frame(Label5)
Gframel.pack(side='left')
Gframee = tk.Frame(Label5)
Gframee.pack(side='right')
Iframel = tk.Frame(Label3, width=9)
Iframel.pack(side='left')
Iframee = tk.Frame(Label3, width=11)
......@@ -92,17 +128,20 @@ OFlabentries.insert(1, tk.Label(Pframel, width=9, text='dω²/dx', pady=2))
OFlabentries.insert(2, tk.Label(Pframel, width=9, text='θ out', pady=2))
OFlabentries.insert(3, tk.Label(Pframel, width=9, text='θ in', pady=2))
OFlabentries.insert(4, tk.Label(Pframel, width=9, text='FWHM0', pady=2))
OFlabentries.insert(5, tk.Label(Pframel, width=9, text='Emin', pady=2))
OFlabentries.insert(6, tk.Label(Pframel, width=9, text='Emax', pady=2))
OFlabentries.insert(5, tk.Label(Gframel, width=9, text='Emin', pady=2))
OFlabentries.insert(6, tk.Label(Gframel, width=9, text='Emax', pady=2))
OFlabentries.insert(7, tk.Label(Pframel, width=9, text='Estep', pady=2))
OFlabentries.insert(8, tk.Label(Pframel, width=9, text='Depth step', pady=2))
for n in range(9):
OFlabentries[n].pack()
OFlabentries.insert(24, tk.Entry(Gframee, width=11))
OFlabentries.insert(25, tk.Label(Gframel, width=9, text='Dose', pady=2))
# Entries
for n in arange(9,18):
OFlabentries.insert(n, tk.Entry(Pframee, width=11))
OFlabentries[n].pack()
if n == 14 or n == 15:
OFlabentries.insert(n, tk.Entry(Gframee, width=11))
else:
OFlabentries.insert(n, tk.Entry(Pframee, width=11))
#Ion
# Labels // Entries numbers:
......@@ -112,15 +151,15 @@ for n in arange(9,18):
OFlabentries.insert(18, tk.Label(Iframel, width=9, text='Energy'))
OFlabentries.insert(19, tk.Label(Iframel, width=9, text='Mass'))
OFlabentries.insert(20, tk.Label(Iframel, width=9, text='Z'))
for n in arange(18,21):
OFlabentries[n].pack()
for n in arange(21,24):
OFlabentries.insert(n, tk.Entry(Iframee, width=11))
for n in range(26):
OFlabentries[n].pack()
##############################################################################
#Initial values
# Initial values
OFlabentries[9].insert(0,param['dedx'])
OFlabentries[10].insert(0,param['dW2dx'])
......@@ -134,6 +173,16 @@ OFlabentries[17].insert(0,0.01)
OFlabentries[21].insert(0,param['E0'])
OFlabentries[22].insert(0, ionb['mass'])
OFlabentries[23].insert(0, ionb['Z'])
OFlabentries[24].insert(0, 12e34)
doseold = float(OFlabentries[24].get())
##############################################################################
# Grahical parameters buttons
Adjbut = tk.Button(Gframee, text = 'Adjust', command = init_calc, width=8)
Databut = tk.Button(Gframel, text = 'Exp Data', command = loadexp, width=6)
Adjbut.pack()
Databut.pack()
##############################################################################
# Energy loss models
......@@ -141,6 +190,7 @@ OFlabentries[23].insert(0, ionb['Z'])
modevar = tk.StringVar()
modevar.set('IONS')
modelvar = tk.StringVar()
modelvar.set('bessel')
R1 = tk.Radiobutton(Label4, text="Gaussian", variable=modelvar, value='gaussiana')
R1.pack()
R2 = tk.Radiobutton(Label4, text="Bessel", variable=modelvar, value='bessel')
......@@ -151,7 +201,9 @@ R2.select()
# Line shape/Convolve gaussian checkbox
LSvar = tk.IntVar()
LSvar.set(1)
CGvar = tk.IntVar()
CGvar.set(1)
C1 = tk.Checkbutton(Label4, text="Use Line Shape", variable=LSvar)
C1.pack()
C1.select()
......
numerical.py
View file @ 369633e2
......@@ -8,10 +8,7 @@ from numpy import nan_to_num, zeros
from scipy.special import i1
from Ewindow import Edict
from layers import Ldict
PI=math.pi
#xx = arange(0,1, 0.1)
#tempor, = plot(xx,xx*0)
####################################################################################################################################################
# Seção de choque Rutherford
......@@ -82,14 +79,14 @@ def convoluiGaussdev(espectro, Sigma,passoE):
####################################################################################################################################################
####################################################################################################################################################
def spectro(p, modelo, ion, OFlabentries, LScontrol, CGausscontrol):
def spectro(p, modelo, ion, OFlabentries, LScontrol, CGausscontrol, tempor):
EPasso = float(OFlabentries[16].get())
passox = float(OFlabentries[17].get())
e = arange(float(OFlabentries[14].get()), float(OFlabentries[15].get()), EPasso)
resol = float(OFlabentries[13].get())
X = e*0.
dose = 12e34
dose = float(OFlabentries[24].get())
gshift=0
Theta_s = PI - (p['Theta_in'] + p['Theta_out'])*PI/180
......@@ -153,13 +150,15 @@ def spectro(p, modelo, ion, OFlabentries, LScontrol, CGausscontrol):
if CGausscontrol == 1:
X, gshift = convoluiGauss(X,resol/2.35482,EPasso)
plt.plot(e-gshift,X*dose, linecolor)
#tempor.set_xdata(e-gshift)
#tempor.set_ydata(X*dose)
#tempor.set_color(linecolor)
#plt.plot(e-gshift,X*dose, linecolor)
tempor.set_xdata(e-gshift)
tempor.set_ydata(X*dose)
#tempor.set_color(linecolor)
#draw()
ylim([0,max(X*dose)*1.1])
xlim([0,max(e-gshift)*1.05])
suptitle('Signal received', fontsize=12)
#suptitle('Signal received', fontsize=12)
xlabel("Enegy (eV)")
ylabel("Yield")
......@@ -172,20 +171,20 @@ def spectro(p, modelo, ion, OFlabentries, LScontrol, CGausscontrol):
xxx= amostra[:,0]*1000
ccc=amostra[:,1]
#plot (xxx, ccc, 'r.')
plot (xxx, ccc, 'r.')
show()
####################################################################################################################################################
####################################################################################################################################################
def spectrolayer(p, modelo, ion, OFlabentries, LScontrol, CGausscontrol):
def spectrolayer(p, modelo, ion, OFlabentries, LScontrol, CGausscontrol, tempor):
EPasso = float(OFlabentries[16].get())
passox = float(OFlabentries[17].get())
e = arange(float(OFlabentries[14].get()), float(OFlabentries[15].get()), EPasso)
resol = float(OFlabentries[13].get())
X = e*0.
dose = 12e34
dose = float(OFlabentries[24].get())
gshift=0.
curdepht=passox #Avoid zero division
......@@ -266,3 +265,97 @@ def spectrolayer(p, modelo, ion, OFlabentries, LScontrol, CGausscontrol):
####################################################################################################################################################
####################################################################################################################################################
def spectroRNA(p, modelo, ion, OFlabentries, LScontrol, CGausscontrol, tempor):
EPasso = float(OFlabentries[16].get())
passox = float(OFlabentries[17].get())
e = arange(float(OFlabentries[14].get()), float(OFlabentries[15].get()), EPasso)
resol = float(OFlabentries[13].get())
X = e*0.
dose = float(OFlabentries[24].get())
gshift=0
Theta_s = PI - (p['Theta_in'] + p['Theta_out'])*PI/180
mi = ion['mass']
## Calcula o espectro o elemento da amostra
mt = Edict[0]['mass']
c = Edict[0]['dist']
alpha_lineshape = Edict[0]['LineShape']
k = ((sqrt (mt**2+mi**2*(sin(Theta_s)**2)) + mi*cos(Theta_s) ) / (mi+mt) )**2
dedxef=p['dedx']*(k/cos(p['Theta_in']*PI/180.) + 1/cos(p['Theta_out']*PI/180.))
dw2dxef=p['dW2dx']*(k**2./cos(p['Theta_in']*PI/180.)+1./cos(p['Theta_out']*PI/180.))
sinalbessel = e*0
sinalgaussiana = e*0
# Calcula a contribuição de cada camada
for x in arange(passox, Edict[0]['profundidademax'], passox):
# Seção de choque - CSRf(Z1, Z2, Ein, M1, M2, Theta):
CS = c[int(x/passox)]*k*CSRf(ion['Z'], Edict[0]['Z'], k*p['E0']-dedxef*x, mi, mt, Theta_s)*passox
# É multiplicada pelo passo para manter a escala do gráfico independente do mesmo.
# Bessel:
# Ainda falta incluir a delta de Dirac na origem
# Existe um valor máximo em energia no qual a função de Bessel pode ser computada em função
# da precisão numérica. Uma possível solução pode ser trabalhar em outras unidades que não
# resultem em argumentos tão grandes para a função exponencial.
# Também temos um problema com a primeira camada, que tem profundidade zero.
# Perceba que os cálculos de Bessel correspondentes às camadas mais fundas estão truncados
# para uma perda de energia maior que 25 keV. (Observado com Theta_in=70, para x > 30)
if modelo == 'bessel':
alpha=dedxef*(2./dw2dxef)
m=alpha*dedxef
lbd=m*x*alpha
besselcamada = lbd*exp(-m*x-alpha*(k*p['E0']-e))*i1(2.*sqrt(lbd*(k*p['E0']-e)))/(sqrt(lbd*(k*p['E0']-e)))
sinalbessel = besselcamada*CS+sinalbessel
# Gaussiana:
elif modelo == 'gaussiana':
Em=k*p['E0']-x*dedxef
sigma=x*dw2dxef
gaussianacamada = exp((-(e-Em)**2)/(2.*sigma))/sqrt(2.*PI*sigma)
sinalgaussiana = gaussianacamada*CS+sinalgaussiana
if modelo == 'bessel':
besself = nan_to_num(sinalbessel)
if LScontrol == 1:
besself = convoluiF0(besself, alpha_lineshape, EPasso)
X = besself + X
linecolor='r'
elif modelo == 'gaussiana':
if LScontrol == 1:
sinalgaussiana = convoluiF0(sinalgaussiana, alpha_lineshape, EPasso)
X = sinalgaussiana + X
linecolor='b'
if CGausscontrol == 1:
X, gshift = convoluiGauss(X,resol/2.35482,EPasso)
tempor.set_xdata(e-gshift)
tempor.set_ydata(X*dose)
#tempor.set_color(linecolor)
ylim([0,max(X*dose)*1.1])
xlim([0,max(e-gshift)*1.05])
xlabel("Enegy (eV)")
ylabel("Yield")
data2 = open('temp.sim', 'w')
for i in range(len(X)):
print>>data2, e[i]-gshift,X[i]*dose
data2.close()
amostra = np.loadtxt('H2.py')
xxx= amostra[:,0]*1000
ccc=amostra[:,1]
#plot (xxx, ccc, 'r.')
show()
####################################################################################################################################################
####################################################################################################################################################
numerical.pyc
View file @ 369633e2
No preview for this file type
temp.sim
View file @ 369633e2
This diff is collapsed.
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment

Os conteúdos dos repositórios estão sob licenças livres e são responsabilidade dos próprios autores, não representando as opiniões e posicionamento da UFRGS ou do CTA.