forked from ctrendafilova/FisherLens
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathplot_fisherlens_result.py
70 lines (58 loc) · 2.3 KB
/
plot_fisherlens_result.py
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
58
59
60
61
62
63
64
65
66
67
68
69
70
import numpy as np
import matplotlib.pyplot as plot
ANN = True
DECAY = False # not impl
SCATTER = False
c = 2.997e8
f_sky = 0.4
TAU_PRIOR = True
if SCATTER:
data = np.load('CLASS_delens/results/step-25/scatter_m0_n0.pkl', allow_pickle=True)
elif ANN:
data = np.load('CLASS_delens/results/ann.pkl', allow_pickle=True)
elif DECAY:
data = np.load('CLASS_delens/results/decay_g26.0.pkl', allow_pickle=True)
else:
raise Exception('bro?')
print(data.keys())
fish = data['fisherGaussian']['delensed']
if TAU_PRIOR:
fish[4, :] = 0
fish[:, 4] = 0
fish[4,4] = 1/0.0074**2
cov = np.linalg.inv(fish)/f_sky
#print(np.sqrt(cov[4,4]))
if ANN:
# fix pann factors of c
units_conversion = 9.e16 # 5.61e20
cov[-1, :] *= units_conversion#c ** 2
cov[:, -1] *= units_conversion#c ** 2
print('pann 95% =', np.sqrt(cov[-1,-1])*2)
elif SCATTER:
#print('m 95% = 1 ±', np.sqrt(cov[-2,-2])*2, 'GeV')
print('σ 95% =', np.sqrt(cov[-1, -1]) * 2, 'cm^2')
elif DECAY:
#print('m 95% = 1 ±', np.sqrt(cov[-2,-2])*2, 'GeV')
print('Γ 95% =', np.sqrt(cov[-1, -1]) * 2, 's^-1')
#print(cov.tolist())
print(np.sqrt(np.diagonal(cov)))
# graph
import plot_triangle
if ANN:
SCALES = [0, 0, 9, 0, 0, 0, 7]
FID = np.array([0.1197, 0.0222, 2.196e-9, 0.9655, 0.054, 0.010409*100, 0])
LABELS = [r'$\Omega_{DM} h^2$', r'$\Omega_b h^2$', r'$10^{9}A_s$', r'$n_s$', r'$\tau_{reio}$', r'$10^{2}\theta_s$',r'$p_{\mathrm{ann}}$ ($10^{-7}\mathrm{m}^3\mathrm{s}^{-1}\mathrm{kg}^{-1}$)']
elif SCATTER:
SCALES = [0, 0, 9, 0, 0, 0, 26]
FID = np.array([0.1197, 0.0222, 2.196e-9, 0.9655, 0.054, 0.010409*100, 0])
LABELS = [r'$\Omega_{DM} h^2$', r'$\Omega_b h^2$', r'$10^{9}A_s$', r'$n_s$', r'$\tau_{reio}$', r'$10^{2}\theta_s$', r'$\sigma_0\ \ (10^{-26}\mathrm{cm}^2)$'] # r'$\log_{10}(m_\chi)$',
elif DECAY:
SCALES = [0, 0, 9, 0, 0, 0, 26]
FID = np.array([0.1197, 0.0222, 2.196e-9, 0.9655, 0.054, 0.010409*100, 0])
LABELS = [r'$\Omega_{DM} h^2$', r'$\Omega_b h^2$', r'$10^{9}A_s$', r'$n_s$', r'$\tau_{reio}$', r'$10^{2}\theta_s$', r'$\Gamma\ \ (10^{-26}\mathrm{s}^{-1})$'] # r'$\log_{10}(m_\chi)$',
plot_triangle.ONE_SIDED.append(LABELS[-1])
plot_triangle.triplot(LABELS, SCALES, FID, cov)
#plot.gcf().set_dpi(100)
plot.gcf().set_size_inches(12, 9)
plot.gcf().savefig('sigma_cmbs4_planck.png')
plot.show()