Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F86500824
gasNvsFe.py
No One
Temporary
Actions
Download File
Edit File
Delete File
View Transforms
Subscribe
Mute Notifications
Award Token
Subscribers
None
File Metadata
Details
File Info
Storage
Attached
Created
Sun, Oct 6, 21:13
Size
5 KB
Mime Type
text/x-python
Expires
Tue, Oct 8, 21:13 (1 d, 23 h)
Engine
blob
Format
Raw Data
Handle
21435016
Attached To
rMTOOLS Mtools
gasNvsFe.py
View Options
#!/usr/bin/env python
from
optparse
import
OptionParser
from
numpy
import
histogram
as
histo
import
Ptools
as
pt
from
pNbody
import
*
def
parse_options
():
usage
=
"usage: %prog [options] file"
parser
=
OptionParser
(
usage
=
usage
)
parser
=
pt
.
add_postscript_options
(
parser
)
parser
=
pt
.
add_ftype_options
(
parser
)
parser
=
pt
.
add_reduc_options
(
parser
)
parser
=
pt
.
add_center_options
(
parser
)
parser
=
pt
.
add_select_options
(
parser
)
parser
=
pt
.
add_cmd_options
(
parser
)
parser
=
pt
.
add_display_options
(
parser
)
parser
=
pt
.
add_info_options
(
parser
)
parser
=
pt
.
add_limits_options
(
parser
)
parser
=
pt
.
add_log_options
(
parser
)
parser
.
add_option
(
"--data"
,
action
=
"store"
,
type
=
"string"
,
dest
=
"data"
,
default
=
None
,
help
=
"data"
)
(
options
,
args
)
=
parser
.
parse_args
()
pt
.
check_files_number
(
args
)
files
=
args
return
files
,
options
def
add_data_FeHistogram
(
data
):
if
data
==
'scl'
:
Fe
,
f
=
pt
.
io
.
read_ascii
(
'Sculptor-Fe-histogram.txt'
,[
0
,
2
])
pt
.
plot
(
Fe
,
f
,
'r--'
,
linewidth
=
1
,
c
=
'r'
)
elif
data
==
'fnx'
:
Fe
,
f
=
pt
.
io
.
read_ascii
(
'Fornax-Fe-histogram.txt'
,[
0
,
2
])
pt
.
plot
(
Fe
,
f
,
'r--'
,
linewidth
=
1
,
c
=
'r'
)
elif
data
==
'car'
:
Fe
,
f
=
pt
.
io
.
read_ascii
(
'Carina-Fe-histogram.txt'
,[
0
,
2
])
pt
.
plot
(
Fe
,
f
,
'r--'
,
linewidth
=
1
,
c
=
'r'
)
elif
data
==
'sex'
:
Fe
,
f
=
pt
.
io
.
read_ascii
(
'Sextans-Fe-histogram.txt'
,[
0
,
2
])
pt
.
plot
(
Fe
,
f
,
'r--'
,
linewidth
=
1
,
c
=
'r'
)
elif
data
==
'leoII'
:
Fe
,
f
=
pt
.
io
.
read_ascii
(
'LeoII-Fe-histogram.txt'
,[
0
,
2
])
pt
.
plot
(
Fe
,
f
,
'r--'
,
linewidth
=
1
,
c
=
'r'
)
#######################################
# MakePlot
#######################################
def
MakePlot
(
files
,
ps
):
class
Option
:
def
__init__
(
self
):
pass
opt
=
Option
()
opt
.
size
=
(
5.12
,
5.12
)
opt
.
cmd
=
"nb.atime=3000;nb.histocenter();nb=nb.selectc(nb.rxyz()<10)"
opt
.
xmin
=
-
4.
opt
.
xmax
=
-
0.5
opt
.
ymin
=
0
opt
.
ymax
=
0.15
opt
.
nxh
=
0
opt
.
select
=
"gas"
opt
.
nopoints
=
False
opt
.
accumulate
=
False
opt
.
density
=
False
opt
.
log
=
None
opt
.
colorbar
=
False
# other default options
opt
.
ftype
=
"gadget"
opt
.
reduc
=
None
opt
.
center
=
None
opt
.
info
=
None
opt
.
display
=
None
opt
.
data
=
None
opt
.
ps
=
ps
# init plot
pt
.
InitPlot
([],
opt
)
pt
.
pcolors
pt
.
labelfont
=
16
fig
=
pt
.
gcf
()
fig
.
subplots_adjust
(
left
=
0.14
)
fig
.
subplots_adjust
(
right
=
0.98
)
fig
.
subplots_adjust
(
bottom
=
0.12
)
fig
.
subplots_adjust
(
top
=
0.98
)
fig
.
subplots_adjust
(
wspace
=
0.0
)
fig
.
subplots_adjust
(
hspace
=
0.0
)
# read files
for
file
in
files
:
#nb = Nbody(file,ftype=opt.ftype)
nb
=
file
# apply options
nb
=
pt
.
do_reduc_options
(
nb
,
opt
)
nb
=
pt
.
do_select_options
(
nb
,
opt
)
nb
=
pt
.
do_center_options
(
nb
,
opt
)
nb
=
pt
.
do_cmd_options
(
nb
,
opt
)
nb
=
pt
.
do_info_options
(
nb
,
opt
)
nb
=
pt
.
do_display_options
(
nb
,
opt
)
################
# get values
################
x
=
nb
.
Fe
()
y
=
nb
.
Mg
()
-
nb
.
Fe
()
db
=
0.1
bins
=
arange
(
-
4
,
0
,
db
)
#n,bins,patches = pt.hist(Fe,bins,normed=True,fc=pt.pcolors[file],lw=0.5)
#c = x > -3
#x = compress(c,x)
h
=
libutil
.
myhistogram
(
x
,
bins
)
h
=
h
.
astype
(
float
)
/
sum
(
h
)
pt
.
bar
(
bins
,
h
,
width
=
0.9
*
(
bins
[
1
]
-
bins
[
0
]),
color
=
'k'
,
alpha
=
0.3
)
x
=
bins
y
=
h
#############################################################
# new method to compute histogram
#############################################################
# ok, ne change quasiment rien, car n_stars nearly identical for all star particles
stellar_mass
=
9.8e-08
N_i
=
stellar_mass
*
2.846012674
*
1.0e10
m_tau
=
(
nb
.
mass
/
stellar_mass
*
(
100.0
**
(
-
0.35
)
-
0.1
**
(
-
0.35
))
+
0.1
**
(
-
0.35
))
**
(
-
1.0
/
0.35
)
n_stars
=
N_i
*
(
m_tau
**
(
-
1.35
)
-
0.1
**
(
-
1.35
))
/
(
100.0
**
(
-
1.35
)
-
0.1
**
(
-
1.35
))
'''
# normal histogram
n = len(Fe)
shape = (len(bins),)
mn = min(bins)
mx = max(bins)
mass = ones(n)
val = ones(n)
# scale between 0 and 1
Fes = (Fe-mn)/(mx-mn)
h = mapping.mkmap1d(Fes.astype(float32),mass.astype(float32),val.astype(float32),shape)
h = h/sum(h)
pt.plot(bins,h,'b')
'''
# weighted histogram
n
=
len
(
x
)
shape
=
(
len
(
bins
),)
mn
=
min
(
bins
)
mx
=
max
(
bins
)
mass
=
n_stars
unity
=
ones
(
n
)
# scale between 0 and 1
Fes
=
(
x
-
mn
)
/
(
mx
-
mn
)
h
=
mapping
.
mkmap1d
(
Fes
.
astype
(
float32
),
mass
.
astype
(
float32
),
unity
.
astype
(
float32
),
shape
)
h
=
h
/
sum
(
h
)
#pt.plot(bins,h,'g')
#############################################################
if
file
==
files
[
0
]:
xmin
,
xmax
,
ymin
,
ymax
=
pt
.
SetLimits
(
opt
.
xmin
,
opt
.
xmax
,
opt
.
ymin
,
opt
.
ymax
,
x
,
y
,
opt
.
log
)
if
file
==
files
[
0
]:
add_data_FeHistogram
(
opt
.
data
)
#if opt.data!=None:
# pt.legend([r'$\rm{%s}$'%opt.data],'upper right',shadow=True)
# plot axis
pt
.
SetAxis
(
xmin
,
xmax
,
ymin
,
ymax
,
opt
.
log
)
pt
.
xlabel
(
r'$\left[ \rm{Fe}/\rm{H} \right]$'
,
fontsize
=
pt
.
labelfont
)
pt
.
ylabel
(
r'$\rm{Stellar\,\,Fraction}$'
,
fontsize
=
pt
.
labelfont
)
pt
.
grid
(
False
)
pt
.
EndPlot
(
files
,
opt
)
########################################################################
# MAIN
########################################################################
#if __name__ == '__main__':
# files,opt = parse_options()
# pt.InitPlot(files,opt)
# pt.pcolors
# MakePlot(files,opt)
# pt.EndPlot(files,opt)
Event Timeline
Log In to Comment