Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F77896101
gplot_histovals
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
Sat, Aug 17, 01:50
Size
4 KB
Mime Type
text/x-python
Expires
Mon, Aug 19, 01:50 (2 d)
Engine
blob
Format
Raw Data
Handle
19943216
Attached To
rGTOOLS Gtools
gplot_histovals
View Options
#!/home/revaz/local/bin/python
'''
Plot histogram of a certain value of the model
Yves Revaz
ven oct 6 10:58:38 CEST 2006
'''
from numarray import *
from Nbody import *
import SM
import string
import sys
import os
import numarray.mlab as mlab
from libjeans import *
from Nbody.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_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("--cgs",
action="store_true",
dest="cgs",
default = 0,
help="invert into cgs units")
parser.add_option("--nbins",
action="store",
dest="nbins",
type="int",
default=25,
help="number of bins")
parser.add_option("--cmd",
action="store",
dest="cmd",
type="string",
default = 'v = nb.T()',
help="command to extract value ex. 'v = nb.T()' ",
metavar=" STRING")
(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
cgs = options.cgs
ftype = options.ftype
select = options.select
nbins = options.nbins
cmd = options.cmd
localsystem = Set_SystemUnits_From_Options(options)
#######################################
# define output system of unit
#######################################
outputunits = UnitSystem('mks',[Unit_kpc, Unit_Ms, Unit_yr , Unit_K, Unit_mol, 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)
# 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
nb.localsystem_of_units = localsystem
exec(cmd)
if xmin==None:
xmin = min(v)
if xmax==None:
xmax = max(v)
x = arange(xmin,xmax,(xmax-xmin)/nbins)
y = histogram(v.astype(Float),x)[:-1]
x = x[:-1]
# print some info
print "vmin = %g vmax =%g"%(v.min(),v.max())
print "mean = %g std =%g"%(v.mean(),mlab.std(v))
print "median = %g "%(mlab.median(nb.rho))
print
print "xmin = %g xmax =%g"%(xmin,xmax)
print "ymin = %g ymax =%g"%(min(y),max(y))
print "number of bins = %d"%(len(y))
print
print "total number of particles = %d\n"%(sum(y))
# use log
if log != None:
x,y = Graph_UseLog(x,y,log)
g.ticksize(-1,0,-1,-1)
if file == files[0]:
xmin,xmax,ymin,ymax = Graph_SetLimits(g,xmin,xmax,ymin,ymax,x,y)
Graph_DrawBox(g,xmin,xmax,ymin,ymax,log)
g.ctype(colors[file])
# plot points
g.histogram(x,y)
g.ctype(0)
if log == 'xy' or log == 'yx':
g.xlabel('log %s'%cmd[4:])
g.ylabel('log N')
elif log == 'x':
g.xlabel('log %s'%cmd[4:])
g.ylabel('N')
elif log == 'y':
g.xlabel('%s'%cmd[4:])
g.ylabel('log N')
else:
g.xlabel('%s'%cmd[4:])
g.ylabel('N')
# -- end ---
Graph_End(g,ps)
Event Timeline
Log In to Comment