Aprimorado o programa de comparação das perdas de energias de Bessel e Gaussian

	Ainda apresenta alguns problemas como trucamento da perda de Energia de Bessel.
	Falta adicionar a delta de Dirac à perda de energia de Bessel.

gaussbesseltest.py

-# Teste para a perda de energia baseada na funcao modificada de Bessel do primeiro tipo de ordem 1.
+#-*- coding: utf-8 -*-
+# Teste para a perda de energia baseada na função modificada de Bessel do primeiro tipo de ordem 1.
 # Apresenta os graficos para a perda de energia para algumas espessuras na geometria definida,
 # comparadas com as correspondentes perdas de energias Gaussianas nas mesmas espessuras.
 # Referencia: R. P. Pezzi, et al. Applied Physics Letters, 92 164102 (2008)
@@ -7,35 +8,51 @@ from pylab import *
 from scipy.special import i1
 PI=math.pi
 
+E0=100000 #eV
+deltaE=25000
+passo=20
+
 dedx=192. # eV/nm
 dw2dx=20740. # eV^2/nm
 Theta_in=0.
 Theta_out=70.
 Theta_s = PI - (Theta_in + Theta_out)*PI/180
 mt = 178.49
-mi = 1
-k = ( (sqrt (mt**2+mi**2*(sin(Theta_s)**2)) + mi*cos(Theta_s) ) / (mi+mt) )**2	
+mi = 1.
+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.))
 dw2dxef=dw2dx*(k**2./cos(Theta_in*PI/180.)+1./cos(Theta_out*PI/180.))
-e = arange(0,23000,20)
-print dedx
-print dedxef/dedx
-print (k/cos(Theta_in*PI/180.) + 1/cos(Theta_out*PI/180.))
+e = arange(0,deltaE, passo)
+
+for x in [.1,1,5,10,15,20,25]:
 
-for x in [1,5,10,15,20,25]:
     # Bessel:  
+    # Ainda falta incluir a delta de Dirac na origem e a convolução com a resolução experimental
+    # 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.
+
+    # Perceba que o cálculo de Bessel correspondente à camada mais funda está trucado.
+
     alpha=dedxef*(2./dw2dxef)
     m=alpha*dedxef
     lbd=m*x*alpha
     bessel = lbd*exp(-m*x-alpha*e)*i1(2.*sqrt(lbd*e))/(sqrt(lbd*e))
-    plt.plot(e,bessel)
+    plt.plot(e+k*E0-deltaE, bessel[::-1])  
+     # [::-1] inverte a ordem dos elementos da matriz. Este é um problema que não
+     # existe na perda de energia Gaussiana pela sua simetria.
+     # +k*E0-deltaE ajusta o vetor energia para a posição correta de kE0 na 
+     # superficie, levando em conta a inversão dos elementos da matriz bessel.
+     # Este artifício vai dificultar a soma dos espectros de cada elemento.
+
+
 
     # Gaussiana:
     Em=x*dedxef
     sigma=x*dw2dxef
     gaussian = exp((-(e-Em)**2)/(2.*sigma))/sqrt(2.*PI*sigma) 
-    plt.plot(e,gaussian)
+    plt.plot(e+k*E0-deltaE,gaussian[::-1])
 
 show()