Page MenuHomec4science
Files F74274150

pplotdmp
No OneTemporary
Actions

File Metadata

Created
Fri, Jul 26, 20:40

pplotdmp
View Options

#!/usr/bin/env python
from optparse import OptionParser
import Ptools as pt
from pNbody import *
from pNbody import units
from pNbody import ctes
from pNbody import thermodyn
import string
from scipy import optimize
def parse_options():
usage = "usage: %prog [options] file"
parser = OptionParser(usage=usage)
parser = pt.add_postscript_options(parser)
parser = pt.add_limits_options(parser)
parser = pt.add_log_options(parser)
parse = pt.add_units_options(parser)
parser.add_option("--x",
action="store",
dest="x",
type="string",
default = 'r',
help="x value to plot",
metavar=" STRING")
parser.add_option("--y",
action="store",
dest="y",
type="string",
default = 'T',
help="y value to plot",
metavar=" STRING")
parser.add_option("--z",
action="store",
dest="z",
type="string",
default = None,
help="z value to plot",
metavar=" STRING")
parser.add_option("--legend",
action="store_true",
dest="legend",
default = False,
help="add a legend")
parser.add_option("--colorbar",
action="store_true",
dest="colorbar",
default = False,
help="add a colorbar")
parser.add_option("--nopoints",
action="store_true",
dest="nopoints",
default = False,
help="do not plot points")
parser.add_option("--map",
action="store_true",
dest="map",
default = False,
help="plot map instead of points")
parser.add_option("--accumulate",
action="store_true",
dest="accumulate",
default = False,
help="integrate histogramm")
parser.add_option("--data",
action="store",
type="string",
dest="data",
default = None,
help="data")
parser.add_option("--mode",
action="store",
type="string",
dest="mode",
default = 'Velocities',
help="mode")
(options, args) = parser.parse_args()
pt.check_files_number(args)
files = args
return files,options
def get_value_and_label(vals,val):
if val == 'R':
label = r'$\rm{Radius}\,\left[ kpc \right]$'
out_units = units.UnitSystem('local',[units.Unit_kpc,units.Unit_Ms,units.Unit_Myr,units.Unit_K])
val = vals['R']
return val,label
if val == 'vc':
label = r'$\rm{Radius}\,\left[ kpc \right]$'
out_units = units.UnitSystem('local',[units.Unit_kpc,units.Unit_Ms,units.Unit_Myr,units.Unit_K])
val = vals['vc']
return val,label
else:
label = r'%s'%val
print "val = %s"%vals[val]
exec("val = %s"%vals[val])
return val, label
#######################################
# MakePlot
#######################################
def MakePlot(files,opt):
# some inits
colors = pt.Colors(n=len(files))
datas = []
# read files
for file in files:
vals = pt.io.read_dmp(file)
################
# units
################
# define local units
unit_params = pt.do_units_options(opt)
if unit_params!=None:
localsystem_of_units = units.Set_SystemUnits_From_Params(unit_params)
else:
if vals.has_key('localsystem_of_units')==0:
localsystem_of_units=units.cgs
else:
localsystem_of_units = vals['localsystem_of_units']
################
# some info
################
print "---------------------------------------------------------"
localsystem_of_units.info()
print "---------------------------------------------------------"
################
# get values
################
if opt.mode=='Numbers':
xlabel = r'$\rm{Radius}\,\left[ kpc \right]$'
ylabel = r'$\rm{Number}\,\rm{of}\,\rm{points}\,\rm{per}\,\rm{bins}$'
x = vals['R']
y = vals['nm']
data = pt.DataPoints(x,y,color=colors.get(),label=r'nn',tpe='line')
datas.append(data)
if opt.mode=='Dispersions':
xlabel = r'$\rm{Radius}\,\left[ kpc \right]$'
ylabel = r'$\rm{Velocity}\,\left[ km/s \right]$'
x = vals['R']
out_units_y = units.UnitSystem('local',[units.Unit_km,units.Unit_Ms,units.Unit_s,units.Unit_K])
tokms = localsystem_of_units.convertionFactorTo(out_units_y.UnitVelocity)
# vct
y = vals['vct']*tokms
data = pt.DataPoints(x,y,color='k',label=r'$V_{c,tot}$',tpe='line')
datas.append(data)
# vc
#y = vals['vc']*tokms
#data = pt.DataPoints(x,y,color='c',label=r'$V_c$',tpe='line')
#datas.append(data)
# vm
y = vals['vm']*tokms
data = pt.DataPoints(x,y,color='r',label=r'$V_m$',tpe='line')
datas.append(data)
# sr
y = vals['sr']*tokms
data = pt.DataPoints(x,y,color='b',label=r'$\sigma_r$',tpe='line')
datas.append(data)
# sp
y = vals['sp']*tokms
data = pt.DataPoints(x,y,color='y',label=r'$\sigma_\theta$',tpe='line')
datas.append(data)
# sz
y = vals['sz']*tokms
data = pt.DataPoints(x,y,color='g',label=r'$\sigma_z$',tpe='line')
datas.append(data)
if opt.mode=='Velocities':
xlabel = r'$\rm{Radius}\,\left[ kpc \right]$'
ylabel = r'$\rm{Velocity}\,\left[ km/s \right]$'
x = vals['R']
out_units_y = units.UnitSystem('local',[units.Unit_km,units.Unit_Ms,units.Unit_s,units.Unit_K])
tokms = localsystem_of_units.convertionFactorTo(out_units_y.UnitVelocity)
if file==files[0]: # skip black
colors.get()
# vct
y = vals['vct']*tokms
data = pt.DataPoints(x,y,color='k',label=r"%s"%'total',tpe='line')
datas.append(data)
# vc
y = vals['vc']*tokms
data = pt.DataPoints(x,y,color=colors.get(),label=r"%s"%file,tpe='line')
datas.append(data)
if opt.mode=='Q':
xlabel = r'$\rm{Radius}\,\left[ kpc \right]$'
ylabel = r'$Q$'
x = vals['R']
y = vals['Q']
data = pt.DataPoints(x,y,color=colors.get(),label=r"%s"%file,tpe='line')
datas.append(data)
if opt.mode=='Energy':
xlabel = r'$\rm{Radius}\,\left[ kpc \right]$'
ylabel = r'$\rm{Energy}$'
x = vals['R']
out_units_y = units.UnitSystem('local',[units.Unit_km,units.Unit_Ms,units.Unit_s,units.Unit_K])
tokms = localsystem_of_units.convertionFactorTo(out_units_y.UnitVelocity)
# vct
y = vals['vct']**2
data = pt.DataPoints(x,y,color='k',label=r'$V_{c,tot}^2$',tpe='line')
datas.append(data)
# total
y = (vals['vm'])**2 + (vals['sr'])**2 + (vals['sp'])**2 + (vals['sz'])**2 + vals['u']
data = pt.DataPoints(x,y,color='c',label=r'$u+V_m^2+\sigma_r^2+\sigma_\theta^2+\sigma_z^2$',tpe='line')
datas.append(data)
if opt.mode=='srQ1':
xlabel = r'$\rm{Radius}\,\left[ kpc \right]$'
ylabel = r'$\sigma$'
out_units_y = units.UnitSystem('local',[units.Unit_km,units.Unit_Ms,units.Unit_s,units.Unit_K])
tokms = localsystem_of_units.convertionFactorTo(out_units_y.UnitVelocity)
x = vals['R']
y = vals['srQ1'] * tokms
data = pt.DataPoints(x,y,color=colors.get(),label=r"%s"%file,tpe='line')
datas.append(data)
if opt.mode=='lambdaj':
xlabel = r'$\rm{Radius}\,\left[ kpc \right]$'
ylabel = r'$\lambda_{\rm{J}}\,\rm{\left[ kpc \right]}$'
x = vals['R']
y = vals['lambdaj']
data = pt.DataPoints(x,y,color=colors.get(),label=r"%s"%file,tpe='line')
datas.append(data)
norm = None
# now, plot
for d in datas:
if d.tpe=='points' or d.tpe=='both':
if opt.log!=None:
if string.find(opt.log,'z')!=-1:
#norm = pt.colors.LogNorm(opt.zmin,opt.zmax)
norm = pt.colors.LogNorm()
else:
norm = None
cmap = pt.GetColormap('rainbow4',revesed=True)
pt.scatter(d.x,d.y,c=d.color,edgecolor=d.color,s=1,linewidths=1,marker='o',vmin=opt.zmin,vmax=opt.zmax,norm=norm,cmap=cmap)
pt.plot(d.x,d.y,c=d.color)
if d.tpe=='line' or d.tpe=='both':
pt.plot(d.x,d.y,color=d.color)
# set limits and draw axis
xmin,xmax,ymin,ymax = pt.SetLimitsFromDataPoints(opt.xmin,opt.xmax,opt.ymin,opt.ymax,datas,opt.log)
# plot axis
pt.SetAxis(xmin,xmax,ymin,ymax,log=opt.log)
pt.xlabel(xlabel,fontsize=pt.labelfont)
pt.ylabel(ylabel,fontsize=pt.labelfont)
pt.grid(False)
if opt.legend:
pt.LegendFromDataPoints(datas)
if opt.colorbar:
#lev = array([-2,0,15])
#l_f = pt.LogFormatter(10, labelOnlyBase=False)
#cb = pt.colorbar(ticks=lev, format = l_f)
cb = pt.colorbar()
cb.set_label(zlabel)
########################################################################
# MAIN
########################################################################
if __name__ == '__main__':
files,opt = parse_options()
pt.InitPlot(files,opt)
pt.pcolors
MakePlot(files,opt)
pt.EndPlot(files,opt)

Event Timeline

c4science · Help