blob: c1236818a0bef40a733fcf50c60575070f528ebd (
plain)
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
|
#!/usr/bin/env python
import sys
import numpy
import scipy.interpolate
center_entropy_file=open('entropy_center.qdp')
entropy_file=open('summary_entropy.qdp')
confidence_level=.68
rout=float(sys.argv[1])
center_s=0
for i in center_entropy_file:
r,s=i.split()
r=float(r)
s=float(s)
if r>rout:
center_s=s
break
new_data=True
s_list=[]
for i in entropy_file:
if i[0]=='n':
new_data=True
continue
if new_data==False:
continue
r,s=i.split()
r=float(r)
s=float(s)
if r>rout:
new_data=False
s_list.append(s)
s_idx=-1
s_list.sort()
for i in range(len(s_list)-1):
if (center_s-s_list[i])*(center_s-s_list[i+1])<=0:
m_idx=i
break
slidx=int(s_idx*(1-confidence_level))
suidx=s_idx-1+int((len(s_list)-s_idx)*confidence_level)
serr1=s_list[slidx]-center_s
serr2=s_list[suidx]-center_s
print("S=\t%e\t %e/+%e keV cm^2 (1 sigma)"%(center_s,serr1,serr2))
|