-
Notifications
You must be signed in to change notification settings - Fork 0
/
moment_like.py
60 lines (45 loc) · 1.53 KB
/
moment_like.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
from cosmosis.datablock import option_section, names
from scipy.interpolate import interp1d
import numpy as np
import csv
mypipe = "moment"
cosmosis_dir = "output/"
def setup(options):
#load in data/make up test data
# fake_theta = np.linspace(0, 2*np.pi*(1./360), num=100)
# fake_moment = []
# for ft in fake_theta:
# fake_moment.append(2*ft)
# err = 1
section = "momentlike"
data_dir = options[section, "data_dir"]
data_name = options[section, "data_name"]
filename = data_dir+data_name
theta, moment, err = np.genfromtxt(filename, unpack=True)
#the return defines config parameter for execute
#return fake_theta, fake_moment, err
return theta, moment, err
def execute(block,config):
theta, moment, err = config
like_name = names.likelihoods
#pull in real data with config
#theoretically computed
theta_theo = block[mypipe, "theta"]
moment_theo = block[mypipe, "moment"]
with open('moments/theo.csv', 'w') as f:
writer = csv.writer(f, delimiter = " ")
writer.writerows(zip(theta_theo, moment_theo))
#compute chi-sq
#chi = (moment - moment_theo(theta))/err
chi = (moment - moment_theo)/err
dof = 2
chi_sq = chi**2
red_chi_sq = chi_sq.sum()/(chi.size-dof)
likeli = -0.5*chi_sq.sum()
print "\t\t\t red_chi_sq", "{:10.1f}".format(red_chi_sq),
print " likelihood {:10.1f} ".format(likeli)
#put chi-sq in block
block[like_name, "MOMENTLIKE_LIKE"] = likeli
return 0
def cleanup(config):
return 0