forked from parejkoj/dust
-
Notifications
You must be signed in to change notification settings - Fork 0
/
constants.py
129 lines (100 loc) · 2.66 KB
/
constants.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
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
import math
import numpy as np
import scipy as sp
def h0():
return 75. #km/s/Mpc
def rho_crit():
return np.float64(1.1e-29)
def omega_d():
return 1e-5
def omega_m():
return 0.3
def omega_l():
return 0.7
def intz(x, y):
from scipy import integrate
return sp.integrate.trapz(y,x)
# Note that scipy calls integration in reverse order as I do
def int(x, y):
dx = x[1:] - x[:-1]
dy = y[1:] - y[:-1]
return np.sum( y[:-1]*dx + 0.5*dx*dy )
def c():
return 3.e10 # cm/s
def h():
return np.float64( 4.136e-18 ) # keV s
def re():
return 2.83e-13 # cm (electron radius)
def mp():
return np.float64( 1.673e-24 ) # g (mass of proton)
def micron2cm():
return 1e-6 * 100
def pc2cm():
return 3.09e18
def cperh0():
return c()*1.e-5 / h0() * (1.e6 * pc2cm() ) #cm
def kev2lam():
return 1.240e-7 # cm keV
def arcs2rad():
return 2.0*np.pi / 360. / 60. / 60. # rad/arcsec
def arcm2rad():
return 2.0*np.pi / 360. / 60. # rad/arcmin
#------- Save and restore functions, similar to IDL -------#
# http://idl2python.blogspot.com/2010/10/save-and-restore-2.html
# Updated April 20, 2012 to store objects
# http://wiki.python.org/moin/UsingPickle
def save( file, varnames, values):
"""
Usage: save('mydata.pysav', ['a','b','c'], [a,b,c] )
"""
import cPickle
f =open(file,"wb")
super_var =dict( zip(varnames,values) )
cPickle.dump( super_var, f )
f.close
def restore(file):
"""
Read data saved with save function.
Usage: data = restore('mydata.pysav')
a = data['a']
b = data['b']
c = data['c']
"""
import cPickle
f=open(file,"rb")
result = cPickle.load(f)
f.close
return result
#------- Read ascii tables --------#
# June 11, 2013
# needed for computers that don't have access to asciidata (hotfoot)
def read_table( filename, ncols, ignore='#' ):
"""
Read data saved in an ascii table
Assumes data is separated by white space
Assumes all the data are floats
Ignores lines that start with the ignore character / sequence
---------------
Usage : read_table( filename, ncols, ignore='#' )
Returns : A dictionary with the column index as keys and the column data as lists
"""
# Initialize result dictionary
result = {}
for i in range(ncols):
result[i] = []
try : f = open( filename, 'r' )
except:
print 'ERROR: file not found'
return
end_of_file = False
while not end_of_file:
try:
temp = f.readline()
if temp[0] == ignore : pass # Ignore the ignore character
else:
temp = temp.split()
for i in range(ncols) : result[i].append( np.float(temp[i]) )
except:
end_of_file = True
f.close()
return result