Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F100631224
gplot_r-sdens_cw
No One
Temporary
Actions
Download File
Edit File
Delete File
View Transforms
Subscribe
Mute Notifications
Award Token
Subscribers
None
File Metadata
Details
File Info
Storage
Attached
Created
Sat, Feb 1, 07:58
Size
6 KB
Mime Type
text/x-python
Expires
Mon, Feb 3, 07:58 (2 d)
Engine
blob
Format
Raw Data
Handle
24002735
Attached To
rGTOOLS Gtools
gplot_r-sdens_cw
View Options
#!/usr/bin/env python
'''
Plot velocity dispertion in z as a function of the radius
Yves Revaz
Mon Jun 4 16:05:18 CEST 2007
'''
from numpy import *
from pNbody import *
from pNbody import cosmo
import Ptools as pt
import string
import sys
import os
from pNbody.libutil import histogram
from optparse import OptionParser
from Gtools import *
from Gtools import vanderwaals as vw
unit_mass = 0.774
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_units_options(parser)
parser = add_log_options(parser)
parser = add_reduc_options(parser)
parser = add_gas_options(parser)
parser = add_center_options(parser)
parser = add_select_options(parser)
parser = add_cmd_options(parser)
parser = add_display_options(parser)
parser.add_option("-t",
action="store",
dest="ftype",
type="string",
default = None,
help="type of the file",
metavar=" TYPE")
parser.add_option("--nbins",
action="store",
dest="nbins",
type="int",
default = 50,
help="number of bins",
metavar=" INT")
parser.add_option("--rmax",
action="store",
dest="rmax",
type="float",
default = 50.,
help="max radius of bins",
metavar=" INT")
parser.add_option("--mode",
action="store",
dest="mode",
type="string",
default = 'sr',
help="mode",
metavar=" STRING")
parser.add_option("--add_legend",
action="store_true",
dest="add_legend",
default = False,
help="add legend")
(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
#######################################
# MakePlot
#######################################
def MakePlot(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
reduc = options.reduc
ftype = options.ftype
gamma = options.gamma
mu = options.mu
av = options.av
bv = options.bv
center = options.center
select = options.select
cmd = options.cmd
display = options.display
nbins = options.nbins
rmax = options.rmax
mode = options.mode
localsystem = Set_SystemUnits_From_Options(options)
#######################################
# define output system of unit
#######################################
outputunits = units.UnitSystem('mks',[units.Unit_km, units.Unit_Ms, units.Unit_s , units.Unit_K, units.Unit_mol, units.Unit_C])
#######################################
# init graph
#######################################
labelfont = 16
ax = pt.gca()
pt.setp(ax.get_xticklabels(),fontsize=labelfont)
pt.setp(ax.get_yticklabels(),fontsize=labelfont)
palette = pt.GetPalette()
colors = pt.SetColorsForFiles(files,palette)
fig = pt.gcf()
fig.subplots_adjust(left=0.11)
fig.subplots_adjust(top=1)
fig.subplots_adjust(bottom=0.1)
fig.subplots_adjust(right=1)
fig.subplots_adjust(wspace=0.15)
fig.subplots_adjust(hspace=0.1)
fac_sdens = 2.23e9/(1000**2) * unit_mass
#######################################
# LOOP
#######################################
i = 0
lss = ['-','--',':']
labels = ['N=1','N=3','N=5']
# read files
for file in files:
nb = Nbody(file,ftype=ftype)
# set model units
nb.localsystem_of_units = localsystem
# center the model
if center=='hdcenter':
print "centering %s"%(center)
nb.hdcenter()
elif center=='histocenter':
print "centering %s"%(center)
nb.histocenter()
# select
#if select!=None:
# print "select %s"%(select)
# nb = nb.select(select)
# reduc
if reduc!= None:
print "reducing %s"%(reduc)
nb = nb.reduc(reduc)
# cmd
if (cmd!= None) and (cmd!= 'None'):
print nb.nbody
print "exec : %s"%cmd
exec(cmd)
# display
if (display!= None) and (display!= 'None'):
nb.display(obs=None,view=display,marker='cross')
################
# get values
################
nb_wg = nb.select('wg')
nb_cg = nb.select('cg')
nb_ba = nb.select('gas','stars','disk','bulge')
r = getr(nr=31,nt=32,rm=rmax)
r = None
r,sdens = nb_ba.msdens(r=r,nb=nbins,rm=rmax)
x0 = r.astype(float)
y0 = sdens.astype(float)*fac_sdens
r = getr(nr=31,nt=32,rm=rmax)
r = None
r,sdens = nb_wg.msdens(r=r,nb=nbins,rm=rmax)
x1 = r.astype(float)
y1 = sdens.astype(float)*fac_sdens
r = getr(nr=31,nt=32,rm=rmax)
r = None
r,sdens = nb_cg.msdens(r=r,nb=nbins,rm=rmax)
x2 = r.astype(float)
y2 = sdens.astype(float)*fac_sdens
x = concatenate((x0,x1,y1))
y = concatenate((y0,x2,y2))
if file == files[0]:
xmin,xmax,ymin,ymax = pt.SetLimits(options.xmin,options.xmax,options.ymin,options.ymax,x,y,options.log)
pt.plot(x0,y0,c='k',ls=lss[i],label=r"$\rm{%s : total\,\,baryons}$"%labels[i])
pt.plot(x1,y1,c='r',ls=lss[i],label=r"$\rm{%s : visible\,\,gas}$"%labels[i])
pt.plot(x2,y2,c='b',ls=lss[i],label=r"$\rm{%s : dark\,\, gas}$"%labels[i])
i = i + 1
pt.SetAxis(xmin,xmax,ymin,ymax,options.log)
pt.xlabel(r'$\rm{Radius}\,\left[ kpc \right]$',fontsize=labelfont)
pt.ylabel(r'$\rm{Surface\,\,Density}\,\left[ M_\odot/pc^2 \right]$',fontsize=labelfont)
# add legend
if options.add_legend:
#labels = [r'$N=1$ : visible gas',r'$N=3$ : visible gas',r'$N=3$ : dark gas']
#pt.legend(crv_lst,'upper right',shadow=True)
pt.legend(loc='upper right',shadow=True,labelsep=0.0001)
####################################
if __name__ == '__main__':
files,opt = parse_options()
pt.InitPlot(files,opt)
#pt.figure(figsize=(8*2,6*2))
#pt.figure(dpi=10)
pt.pcolors
MakePlot(files,opt)
pt.EndPlot(files,opt)
Event Timeline
Log In to Comment