Page MenuHomec4science

gvel
No OneTemporary

File Metadata

Created
Tue, May 21, 21:26
#!/usr/bin/env python
'''
Extract and plot velocities dispersions and
mean azimuthal velocities from an Nbody file
default file extension = .vel
Yves Revaz
Tue May 31 10:18:56 CEST 2005
'''
import sys
from numarray import *
import SM
from Nbody import *
from Nbody import io
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("-l", "--label",
# action="store_true",
# dest="label",
# default=0,
# help="display the movie")
parser.add_option("--minx",
action="store",
dest="mnx",
default=None,
help="min value in x")
parser.add_option("--maxx",
action="store",
dest="mxx",
default=None,
help="max value in x")
parser.add_option("--miny",
action="store",
dest="mny",
default=None,
help="min value in y")
parser.add_option("--maxy",
action="store",
dest="mxy",
default=None,
help="max value in y")
parser.add_option("-p",
action="store",
dest="ps",
type="string",
default = None,
help="postscript filename",
metavar=" FILE")
parser.add_option("-c",
action="store_true",
dest="center",
default = 1,
help="center using histocenter")
parser.add_option("-t",
action="store",
dest="ftype",
type="string",
default = 'gadget',
help="file type",
metavar=" FTYPE")
parser.add_option("-i",
action="store",
dest="ifile",
type="string",
default = None,
help="input file",
metavar=" FILE")
parser.add_option("-o",
action="store",
dest="ofile",
type="string",
default = None,
help="output file",
metavar=" FILE")
parser.add_option("-s",
action="store",
dest="select",
type="string",
default = None,
help="select particles of type NAME",
metavar=" NAME")
(options, args) = parser.parse_args()
# if len(args) == 0:
# print "you must specify a filename"
# sys.exit(0)
if len(args)==0:
file = None
else:
file = args[0]
return file,options
##########################################################
#
# MAIN
#
#########################################################
# get options
file,options = parse_options()
ps = options.ps
mnx = options.mnx
mxx = options.mxx
mny = options.mny
mxy = options.mxy
center = options.center
ftype = options.ftype
ifile = options.ifile
ofile = options.ofile
select = options.select
f = 978.
#############################
# open file
#############################
if ifile==None:
nb = Nbody(file,ftype=ftype)
# histocenter
if center:
print "histocenter"
nb.histocenter()
# select
if select:
print "select"
nb = nb.select(select)
# compute velocities
r,sr,st,sz,mt = nb.sigma()
sr = sr*f
st = st*f
sz = sz*f
mt = mt*f
else:
r,sr,st,sz,mt = io.read_ascii(ifile,[0,1,2,3,4])
#############################
# output file
#############################
if ofile!=None:
fd = open(ofile,'w')
for i in range(len(r)):
fd.write('%8.3f %8.3f %8.3f %8.3f %8.3f\n'%(r[i],sr[i],st[i],sz[i],mt[i]))
fd.close()
#############################
# graph
#############################
# 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')
# set limits
x = concatenate((sr,st,sz,mt))
if mnx==None:
mnx = min(r)
if mxx==None:
mxx = max(r)
if mny==None:
mny = min(x)
if mxy==None:
mxy = max(x)
g.location(3500, 30000, 3500, 19000)
g.limits(mnx,mxx,mny,mxy)
g.box()
# draw lines
g.ctype('cyan')
g.connect(r,mt)
g.ctype('green')
g.connect(r,sr)
g.ctype('red')
g.connect(r,st)
g.ctype('blue')
g.connect(r,sz)
#g.ctype('black')
#g.connect(r,vc)
g.ctype('black')
# labels
g.xlabel('Radius [kpc]')
g.ylabel('velocities [km/s]')
#------- labels -------
g.location(3500, 30000, 19000, 23000)
g.limits(0,13,1,3)
#
#g.ltype(0)
#g.ctype('black')
#g.relocate(5,1.5)
#g.draw(7,1.5)
#g.relocate(8,1.5)
#g.putlabel(5,'\\ v_c')
#
g.ltype(1)
g.ctype('cyan')
g.relocate(1,1.5)
g.draw(3,1.5)
g.relocate(4,1.5)
g.putlabel(5,'\\ \\bar{v_\\theta}')
#
g.ltype(2)
g.ctype('green')
g.relocate(1,2)
g.draw(3,2)
g.relocate(4,2)
g.putlabel(5,'\\ \\sigma_r')
#
g.ltype(3)
g.ctype('red')
g.relocate(5,2)
g.draw(7,2)
g.relocate(8,2)
g.putlabel(5,'\\ \\sigma_\\theta')
#
g.ltype(4)
g.ctype('blue')
g.relocate(9,2)
g.draw(11,2)
g.relocate(12,2)
g.putlabel(5,'\\ \\sigma_z')
g.ltype(0)
g.ctype('black')
# old limits
g.location(3500, 30000, 3500, 19000)
g.limits(mnx,mxx,mny,mxy)
# -- end ---
if ps==None:
g.show()
else:
g.write()
g.clean()

Event Timeline