~c-tunnell1/+junk/vlenf_plots

« back to all changes in this revision

Viewing changes to create_simple_flux.py

  • Committer: Christopher Tunnell
  • Date: 2011-09-30 10:14:01 UTC
  • Revision ID: c.tunnell1@physics.ox.ac.uk-20110930101401-hdqun0c3ux6y9gvk
flux check

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
 
import numpy as np
2
 
 
3
 
import time
4
 
import math
5
 
import sys
6
 
 
7
 
def frange(startValue, endValue, steps):
8
 
    incrementValue = (float(endValue) - startValue) / steps
9
 
    for i in range(steps):
10
 
        yield float(startValue) + i*incrementValue
11
 
 
12
 
polarization = 0.0
13
 
pi = 3.1415
14
 
parent_energy = 3.0 # GeV 
15
 
 
16
 
E_min = 0.10  # GeV
17
 
E_max = parent_energy # GeV
18
 
E_steps = 500
19
 
 
20
 
mass_muon = 0.1056 # GeV
21
 
 
22
 
stored_muons = 1e17
23
 
 
24
 
theta = 0
25
 
 
26
 
distance = 1 # km, GLoBES normalization
27
 
 
28
 
flux_file_plus = open('flux_file_plus.dat', 'w')
29
 
flux_file_minus = open('flux_file_minus.dat', 'w')
30
 
 
31
 
for E in frange (E_min, E_max, E_steps):
32
 
    cos_theta = 1
33
 
    y = float(E) / parent_energy
34
 
    beta = math.sqrt(1 - mass_muon**2 / parent_energy**2)
35
 
 
36
 
    N = stored_muons
37
 
    N /= pi * distance**2 * mass_muon**6
38
 
    N *= parent_energy**4 * y**2 * (1 - beta * math.cos(theta))
39
 
 
40
 
    flux_mu = 4 * N * (3 * mass_muon**2 - 4 * parent_energy**2 * y * (1 - beta * math.cos(theta)))
41
 
    flux_e = 24 * N * (1 * mass_muon**2 - 2 * parent_energy**2 * y * (1 - beta * math.cos(theta)))
42
 
    
43
 
    
44
 
    flux_e  *= (2*pi) / parent_energy  # solid angle -> cos(theta), y -> E
45
 
    flux_mu *= (2*pi) / parent_energy  # solid angle -> cos(theta), y -> E  
46
 
 
47
 
    flux_file_plus.write('%0.9g %0.9g %0.9g 0 %0.9g %0.9g 0\n' % (E, flux_e, 0, 0, flux_mu))
48
 
    flux_file_minus.write('%0.9g %0.9g %0.9g 0 %0.9g %0.9g 0\n' % (E, 0, flux_mu, flux_e, 0))
49
 
 
50
 
 
51
 
flux_file_plus.close()
52
 
flux_file_minus.close()
53