Page MenuHomec4science

gplot_rho
No OneTemporary

File Metadata

Created
Mon, Aug 19, 18:57

gplot_rho

#!/home/revaz/local/bin/python
'''
Plot density of the model as a function of radius
(isotropique density)
Yves Revaz
Mon Nov 21 14:54:44 CET 2005
'''
from numarray import *
from Nbody import *
import SM
import string
import sys
import os
from libjeans import *
from Nbody.libutil import histogram
try:
from optparse import OptionParser
except ImportError:
from optik import OptionParser
def parse_options():
usage = "usage: %prog [options] file"
parser = OptionParser(usage=usage)
parser.add_option("-t",
action="store",
dest="ftype",
type="string",
default = None,
help="type of the file",
metavar=" TYPE")
parser.add_option("--xmin",
action="store",
dest="xmin",
type="float",
default=0,
help="min value in x")
parser.add_option("--xmax",
action="store",
dest="xmax",
type="float",
default=50,
help="max value in x")
parser.add_option("--nb",
action="store",
dest="nb",
type="int",
default=25,
help="number of bins")
parser.add_option("--ymin",
action="store",
dest="mny",
type="float",
default=0,
help="min value in y")
parser.add_option("--ymax",
action="store",
dest="mxy",
type="float",
default=None,
help="max value in y")
parser.add_option("-p",
action="store",
dest="ps",
type="string",
default = None,
help="postscript filename",
metavar=" FILE")
(options, args) = parser.parse_args()
if len(args) == 0:
print "you must specify a filename"
sys.exit(0)
file = args[0]
return file,options
#############################
# graph
#############################
# get options
file,options = parse_options()
ps = options.ps
xmin = options.xmin
xmax = options.xmax
nb = options.nb
ymin = options.mny
ymax = options.mxy
ftype = options.ftype
#######################################
# open sm
#######################################
if ps==None:
g = SM.plot("x11 -bg white -fg black ")
else:
g = SM.plot("postencap %s"%ps)
# some init
g.palette('bgyrw')
g.expand(0.999)
g.setvariable('TeX_strings', '1')
#######################################
# first
g.location(3500, 31000, 3500, 31000)
#######################################
# LOOP
#######################################
vs = array([],Float)
ts= array([],Float)
nbody = Nbody(file,ftype=ftype)
r,dens = nbody.dens(nb=nb,rm=xmax)
dens = log10(dens)
print r,dens
if ymin == None:
ymin = min(dens)
if ymax == None:
ymax = max(dens)
g.location(3500, 31000, 3500, 31000)
g.limits(xmin,xmax,ymin,ymax)
g.box()
g.connect(r,dens)
g.histogram(r,dens)
#rho = log10( 1000./(1.+r)**2 ) + 0.1
#g.ctype(192)
#g.connect(r,rho)
#g.ctype(0)
g.ctype(0)
g.xlabel('R')
g.ylabel('N')
# -- end ---
if ps==None:
g.show()
else:
g.write()
g.clean()

Event Timeline