lambda_gg.py
1.42 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
#!/usr/bin/python
from sys import argv
import matplotlib.pyplot as plt
from matplotlib import gridspec
from numpy import *
from scipy.integrate import quad
from Modules.Analytic import Ethreshold_gg
c=2.99792458e10 # cm.s-1
Mpc=(3.0856776e+16)*1e8 # Mpc to cm
H0=67.8*1e5/(Mpc) # s-1
omegaM = 0.3
omegaK = 0
omegaL = 0.7
zlim=-0.
def properIntegrand(z):
return -c/(H0*(1+z)*sqrt(omegaM*(1+z)**3+omegaK*(1+z)**2+omegaL))
def comobileIntegrand(z):
return -c/(H0*sqrt(omegaM*(1+z)**3+omegaK*(1+z)**2+omegaL))
color=['b','r','g','m','c','y']
# Proper distance figure
#========================
fig1 = plt.figure()
ax11 = fig1.add_subplot(111)
labels=["0","0.5","1","2","3"]
theory=[1,2,3,4,5]
ind=0
for lab in labels:
# theoritical curve ===================================================
# f: without cosmo
# g: with cosmo
e,f,g=loadtxt('lambda_e.dat',unpack=True,usecols=[0,int(theory[ind]),int(theory[ind])+5])
e = e*511.e3/1.e9 #GeV
cond= (e>100) & (e<1e5)
e=e[cond]
g=g[cond]
f=f[cond]
p = ax11.plot(e,g,"-"+color[ind],label="z = "+lab)
ax11.plot(e,f,color=p[0].get_color(),linestyle='--')
ind=ind+1
ax11.axvline(x=Ethreshold_gg(), ymin=0., ymax = 1., color='r', linewidth=2)
ax11.legend(loc="best")
ax11.grid(b=True,which='major')
ax11.set_ylim([1e-4,1e4])
ax11.set_xscale('log')
ax11.set_yscale('log')
ax11.set_xlabel("energy [GeV]")
ax11.set_ylabel("$\lambda_{\gamma\gamma}$ [Mpc]")
plt.show()