from numpy import *

mpr = 938.27208816  #MeV
mp = 1.6726219e-27  #kg
B = 5.e-9           #B in T
E = 1.0             #in keV
qe = 1.60217662e-19 #C
c = 299792458.      #m/s
pa = 30.            #Pitchwinkel in Grad

def deg_to_rad(angle):
    return angle/180.*pi

def eV_to_J(E_eV):
    return E_eV * qe

def E_to_v(E):
    beta_sq = 1. - 1./(1.+E/mp/c/c)**2
    return sqrt(beta_sq)*c

def v_par(v,pa):
    return cos(deg_to_rad(pa))*v

def v_perp(v,pa):
    return sin(deg_to_rad(pa))*v

def gamma(v):
    return 1./sqrt(1.-(v/c*v/c))

#Geschwindigkeit des Protons:
vp = E_to_v(eV_to_J(E*1.e3))
print 'Die Geschwindigkeit des Protons betraegt {} km/s'.format(vp/1.e3)

#Wie weit kommt es in 5 Minuten entlang von B
t = 5*60.       #5 Minuten
s_perp = vp*sin(deg_to_rad(pa))*t
print 'Das Proton fliegt in {} Sekunden {} km senkrecht zu B'.format(t,s_perp/1.e3)

#Wie weit kommt es in 5 Minuten senkrecht zu B
s_par = vp*cos(deg_to_rad(pa))*t
print 'Das Proton fliegt in {} Sekunden {} km entlang von B'.format(t,s_par/1.e3)

#Wie weit ist es effektiv geflogen?
s_tot = vp*t
print 'Das Proton ist in dieser Zeit aber {} km weit geflogen'.format(s_tot/1.e3)

#Gyroradius
rg = mp*v_perp(vp,pa)/qe/B*gamma(vp)
print 'Der Gyroradius betraegt {} m'.format(rg)

#Wie viele Umdrehungen hat es um das Magnetfeld vollfuehrt?
#n = qe*B/2/pi/mp/gamma(v_perp(vp,pa))*t
n = v_perp(vp,pa)*t/2/pi/rg
print 'Das Proton hat in dieser Zeit {} Umdrehungen vollfuehrt'.format(n)

