Page MenuHomec4science
Files F73278432

gplot_energy
No OneTemporary
Actions

File Metadata

Created
Sat, Jul 20, 05:58

gplot_energy
View Options

#!/export/revaz/local/bin/python
'''
Extract and plot energy and mass values contained in the
output Gadget file called by default "energy.txt".
Yves Revaz
ven jun 9 10:43:59 CEST 2006
'''
from numpy import *
from pNbody import *
import SM
import string
import sys
import os
from pNbody.libutil import histogram
from optparse import OptionParser
from Gtools import *
from Gtools import io
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_log_options(parser)
(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
#######################################
# open sm
#######################################
g = Graph_Init(ps)
Graph_SetDefaultsGraphSettings(g)
colors = Graph_SetColorsForFiles(files,col)
#######################################
# LOOP
#######################################
# set ltype
ltype = [0,1,2,3,4,5,6]
i = 0
# read files
for file in files:
try:
vals = io.read_new_energy(file)
except "NotNewEgyFileError":
vals = io.read_energy(file,iobs=None)
x = vals['Time']
y = 0
xs = array([],float)
ys = array([],float)
EnergyInt = None
EnergyPot = None
EnergyKin = None
EnergyRadSph = None
EnergyRadSticky = None
EnergySfr = None
EnergyFeedback = None
EnergyBubbles = None
EnergyAGNHeat = None
if vals.has_key('EnergyInt'):
EnergyInt = vals['EnergyInt']
y = y + EnergyInt
xs = concatenate((xs,x))
ys = concatenate((ys,EnergyInt))
if vals.has_key('EnergyPot'):
EnergyPot = vals['EnergyPot']
y = y + EnergyPot
xs = concatenate((xs,x))
ys = concatenate((ys,EnergyPot))
if vals.has_key('EnergyKin'):
EnergyKin = vals['EnergyKin']
y = y + EnergyKin
xs = concatenate((xs,x))
ys = concatenate((ys,EnergyKin))
if vals.has_key('EnergyRadSph'):
EnergyRadSph = vals['EnergyRadSph']
y = y + EnergyRadSph
xs = concatenate((xs,x))
ys = concatenate((ys,EnergyRadSph))
if vals.has_key('EnergyRadSticky'):
EnergyRadSticky = vals['EnergyRadSticky']
y = y + EnergyRadSticky
xs = concatenate((xs,x))
ys = concatenate((ys,EnergyRadSticky))
if vals.has_key('EnergySfr'):
EnergySfr = vals['EnergySfr']
y = y + EnergySfr
xs = concatenate((xs,x))
ys = concatenate((ys,EnergySfr))
if vals.has_key('EnergyFeedback'):
EnergyFeedback = vals['EnergyFeedback']
y = y + EnergyFeedback
xs = concatenate((xs,x))
ys = concatenate((ys,EnergyFeedback))
if vals.has_key('EnergyBubbles'):
EnergyBubbles = vals['EnergyBubbles']
y = y + EnergyBubbles
xs = concatenate((xs,x))
ys = concatenate((ys,EnergyBubbles))
if vals.has_key('EnergyAGNHeat'):
EnergyAGNHeat = vals['EnergyAGNHeat']
y = y + EnergyAGNHeat
xs = concatenate((xs,x))
ys = concatenate((ys,EnergyAGNHeat))
# # use log
# if log != None:
# x,y = Graph_UseLog(x,y,log)
if file == files[0]:
xmin,xmax,ymin,ymax = Graph_SetLimits(g,xmin,xmax,ymin,ymax,xs,ys)
g.box()
# plot points
g.ltype(ltype[i])
t = x
g.ctype('black')
g.connect(t,y)
g.points(t,y)
if EnergyInt != None:
g.ctype(160) # yellow
g.connect(t,EnergyInt)
g.points(t,EnergyInt)
print 'EnergyInt'
if EnergyPot != None:
g.ctype(192) # red
g.connect(t,EnergyPot)
g.points(t,EnergyPot)
print 'EnergyPot'
if EnergyKin != None:
g.ctype(128) # green
g.connect(t,EnergyKin)
g.points(t,EnergyKin)
print 'EnergyKin'
if EnergyRadSph != None:
g.ctype(64) # blue
g.connect(t,EnergyRadSph)
g.points(t,EnergyRadSph)
print 'EnergyRadSph'
if EnergyRadSticky != None:
g.ctype(96) # cyan
g.connect(t,EnergyRadSticky)
g.points(t,EnergyRadSticky)
print 'EnergyRadSticky'
if EnergySfr != None:
g.ctype(132) # green+
g.connect(t,EnergySfr)
g.points(t,EnergySfr)
print 'EnergySfr'
if EnergyFeedback != None:
g.ctype(250) # red+
g.connect(t,EnergyFeedback)
g.points(t,EnergyFeedback)
print 'EnergyFeedback'
if EnergyBubbles != None:
g.ctype(175) # pink
g.connect(t,EnergyBubbles)
g.points(t,EnergyBubbles)
print 'EnergyBubbles'
if EnergyAGNHeat != None:
g.ctype(168) # ?
g.connect(t,EnergyAGNHeat)
g.points(t,EnergyAGNHeat)
print 'EnergyAGNHeat'
i=i+1
# labels
obs = "energy"
g.ltype(0)
g.ctype(0)
if log == 'xy' or log == 'yx':
g.xlabel('log T')
g.ylabel('log %s'%obs)
elif log == 'x':
g.xlabel('log T')
g.ylabel('%s'%obs)
elif log == 'y':
g.xlabel('T')
g.ylabel('log %s'%obs)
else:
g.xlabel('T')
g.ylabel('%s'%obs)
# -- end ---
Graph_End(g,ps)

Event Timeline

c4science · Help