Page MenuHomec4science

gplotqa
No OneTemporary

File Metadata

Created
Sat, Feb 1, 09:52
#!/usr/bin/env python
'''
Extract and plot data from output from jeans fortran program
Yves Revaz
Tue Aug 16 14:05:26 CEST 2005
'''
from numarray import *
import SM
import string
import sys
from libjeans import *
from scipy.interpolate import splrep,splev
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("--minx",
action="store",
dest="mnx",
default=0,
help="min value in x")
parser.add_option("--maxx",
action="store",
dest="mxx",
default=30,
help="max value in x")
parser.add_option("--miny",
action="store",
dest="mny",
default=0,
help="min value in y")
parser.add_option("--maxy",
action="store",
dest="mxy",
default=220,
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",
dest="colors",
type="string",
default = None,
help="colors",
metavar=" 0,64,192")
(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
rmin = options.mnx
rmax = options.mxx
vmin = options.mny
vmax = options.mxy
col = options.colors
# set colors
colors = {}
i = 0
for file in files:
if col!=None:
colors[file]=col[i]
else:
colors[file] = i*255/len(files)
i = i + 1
#######################################
# 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, 17750, 31000)
g.limits(rmin,rmax,-1,3)
g.box(0,2)
g.relocate(-200,1)
g.draw(200,1)
g.ylabel('Q')
#######################################
# second
g.location(3500, 31000, 3500, 16750)
g.limits(rmin,rmax,0.1,1.)
g.box()
g.relocate(0,0.293)
g.draw(100,0.293)
g.xlabel('R [kpc]')
g.ylabel('A')
#######################################
# LOOP
#######################################
for file in files:
R,Sdens,Omega,Kappa,OmK2,Nu,Vr,Vt,Vz,Vc,Sr,St,Sz,H,Q=read_jeans(file)
r = R
q = Q
sr = Sr
sz = Sz
nr1 = 5
nr2 = 25
r = r[nr1:nr2]
q = q[nr1:nr2]
sr = sr[nr1:nr2]
sz = sz[nr1:nr2]
# color
g.ctype(colors[file])
# first
# interpolate
a = splrep(r,q,s=0)
ri = arange(min(r),max(r),0.1).astype(Float)
q = splev(ri,a)
g.location(3500, 31000, 17750, 31000)
g.limits(rmin,rmax,-1,3)
g.connect(ri,q)
# second
# interpolate
a = splrep(r,sr,s=0)
b = splrep(r,sz,s=0)
r = arange(min(r),max(r),0.1).astype(Float)
sr = splev(ri,a)
sz = splev(ri,b)
szr = sz/sr
g.location(3500, 31000, 3500, 16750)
g.limits(rmin,rmax,0.1,1.)
g.connect(ri,szr)
# colors
if len(files)>1:
x0 = 0.05
y0 = 0.95
x = x0
y = y0
dx = 0.5/len(files)
g.location(3500,31000,3500,31000)
g.limits(0,1,0,1)
g.ptype(10,3)
g.expand(1.5)
for i in range (len(colors)):
xp=array([x],Float)
yp=array([y],Float)
g.ctype(colors[files[i]])
g.points(xp,yp)
x = x + dx
# -- end ---
if ps==None:
g.show()
else:
g.write()
g.clean()

Event Timeline