Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F78441764
gplot_r-mr
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
Tue, Aug 20, 22:01
Size
3 KB
Mime Type
text/x-python
Expires
Thu, Aug 22, 22:01 (2 d)
Engine
blob
Format
Raw Data
Handle
20047082
Attached To
rGTOOLS Gtools
gplot_r-mr
View Options
#!/home/revaz/local/bin/python
'''
Plot the mass of the model as a function of radius,
assuming a sperical density distribution.
Yves Revaz
Fri Apr 27 09:30:11 CEST 2007
'''
from
numarray
import
*
from
pNbody
import
*
from
pNbody
import
cosmo
import
SM
import
string
import
sys
import
os
from
libjeans
import
*
from
pNbody.libutil
import
histogram
from
optparse
import
OptionParser
from
Gtools
import
*
from
Gtools
import
vanderwaals
as
vw
def
parse_options
():
usage
=
"usage: %prog [options] file"
parser
=
OptionParser
(
usage
=
usage
)
parser
=
add_postscript_options
(
parser
)
parser
=
add_color_options
(
parser
)
parser
=
add_limits_options
(
parser
)
parser
=
add_units_options
(
parser
)
parser
=
add_log_options
(
parser
)
parser
=
add_reduc_options
(
parser
)
parser
=
add_center_options
(
parser
)
parser
=
add_select_options
(
parser
)
parser
.
add_option
(
"-t"
,
action
=
"store"
,
dest
=
"ftype"
,
type
=
"string"
,
default
=
None
,
help
=
"type of the file"
,
metavar
=
" TYPE"
)
parser
.
add_option
(
"--rmax"
,
action
=
"store"
,
dest
=
"rmax"
,
type
=
"float"
,
default
=
50
,
help
=
"max radius"
,
metavar
=
" FLOAT"
)
parser
.
add_option
(
"--nbins"
,
action
=
"store"
,
dest
=
"nbins"
,
type
=
"int"
,
default
=
25
,
help
=
"number of bins"
)
parser
.
add_option
(
"--cgs"
,
action
=
"store_true"
,
dest
=
"cgs"
,
default
=
0
,
help
=
"invert into cgs units"
)
(
options
,
args
)
=
parser
.
parse_args
()
if
options
.
colors
!=
None
:
exec
(
"options.colors = array([
%s
])"
%
(
options
.
colors
))
if
len
(
args
)
==
0
:
print
"you must specify a filename"
sys
.
exit
(
0
)
files
=
args
return
files
,
options
#############################
# graph
#############################
# get options
files
,
options
=
parse_options
()
ps
=
options
.
ps
col
=
options
.
colors
xmin
=
options
.
xmin
xmax
=
options
.
xmax
ymin
=
options
.
ymin
ymax
=
options
.
ymax
log
=
options
.
log
reduc
=
options
.
reduc
center
=
options
.
center
select
=
options
.
select
cgs
=
options
.
cgs
ftype
=
options
.
ftype
rmax
=
options
.
rmax
nbins
=
options
.
nbins
localsystem
=
Set_SystemUnits_From_Options
(
options
)
localsystem
.
get_UnitDensity_in_cgs
()
#######################################
# define output system of unit
#######################################
outputunits
=
units
.
UnitSystem
(
'mks'
,[
units
.
Unit_kpc
,
units
.
Unit_Ms
,
units
.
Unit_yr
,
units
.
Unit_K
,
units
.
Unit_mol
,
units
.
Unit_C
])
#######################################
# open sm
#######################################
g
=
Graph_Init
(
ps
)
Graph_SetDefaultsGraphSettings
(
g
)
colors
=
Graph_SetColorsForFiles
(
files
,
col
)
#######################################
# LOOP
#######################################
# read files
for
file
in
files
:
nb
=
Nbody
(
file
,
ftype
=
ftype
)
# set unit parameter
nb
.
localsystem_of_units
=
localsystem
# center the model
if
center
==
'hdcenter'
:
print
"centering
%s
"
%
(
center
)
nb
.
hdcenter
()
elif
center
==
'histocenter'
:
print
"centering
%s
"
%
(
center
)
nb
.
histocenter
()
# select
if
select
!=
None
:
print
"select
%s
"
%
(
select
)
nb
=
nb
.
select
(
select
)
if
reduc
!=
None
:
print
"reducing
%s
"
%
(
reduc
)
nb
=
nb
.
reduc
(
reduc
)
# get values
x
,
y
=
nb
.
mr
(
nb
=
nbins
,
rm
=
rmax
)
# use log
if
log
!=
None
:
x
,
y
=
Graph_UseLog
(
x
,
y
,
log
)
if
file
==
files
[
0
]:
Graph_SetLimits
(
g
,
xmin
,
xmax
,
ymin
,
ymax
,
x
,
y
)
g
.
box
()
# plot points
g
.
ctype
(
colors
[
file
])
g
.
connect
(
x
,
y
)
g
.
ctype
(
0
)
if
log
==
'xy'
or
log
==
'yx'
:
g
.
xlabel
(
'log R'
)
g
.
ylabel
(
'log Velocity'
)
elif
log
==
'x'
:
g
.
xlabel
(
'log R'
)
g
.
ylabel
(
'Velocity'
)
elif
log
==
'y'
:
g
.
xlabel
(
'R'
)
g
.
ylabel
(
'log Velocity'
)
else
:
g
.
xlabel
(
'R'
)
g
.
ylabel
(
'Velocity'
)
# -- end ---
Graph_End
(
g
,
ps
)
Event Timeline
Log In to Comment