Homec4science
Diffusion ARRAKIHS (master) Browse

arrakihs/Workspace/InjectedDwarfs_conflimit/ark_idealized_dwarfs_pR12vsL
d192e801c42bmaster

Last Change
Raw File

ark_idealized_dwarfs_pR12vsL
View Options

#!/usr/bin/python3
import os
import numpy as np
import argparse
import pickle
import matplotlib.pyplot as plt
####################################################################
# option parser
####################################################################
description="plot the stellar metallicity as a function of luminosity"
epilog ="""
Examples:
--------
ark_idealized_dwarfs_pFevsL output.pkl
"""
parser = argparse.ArgumentParser(description=description,epilog=epilog,formatter_class=argparse.RawDescriptionHelpFormatter)
parser.add_argument(action="store",
dest="files",
metavar='FILE',
type=str,
default=None,
nargs='*',
help='a list of files')
parser.add_argument("-o",
action="store",
type=str,
dest="outputfilename",
default=None,
help="Name of the output file")
parser.add_argument("--data",
action="store",
type=str,
dest="data",
default=None,
help="data file")
parser.add_argument('--xmin',
action="store",
dest="xmin",
metavar='FLOAT',
type=float,
default=1e5,
help='x min')
parser.add_argument('--xmax',
action="store",
dest="xmax",
metavar='FLOAT',
type=float,
default=1e9,
help='x max')
parser.add_argument('--ymin',
action="store",
dest="ymin",
metavar='FLOAT',
type=float,
default=5e-2,
help='y min')
parser.add_argument('--ymax',
action="store",
dest="ymax",
metavar='FLOAT',
type=float,
default=5,
help='y max')
parser.add_argument('--log',
action="store",
dest="log",
metavar='STR',
type=str,
default='xy',
help='log scale (None,x,y,xy)')
class sdict(dict):
def selectc(self,c):
import copy
d = copy.deepcopy(self)
for key in d.keys():
d[key] = np.compress(c,d[key])
return d
def readpkl(file):
f = open(file,"rb")
data = pickle.load(f)
f.close()
return data
def readData(file):
f = open(file,"rb")
data = pickle.load(f)
f.close()
odata = sdict()
names = np.array([])
for name in data.keys():
names = np.append(names,name)
odata["name"] = names
attrs = data[name].keys()
# loop over
for attr_name in attrs:
# create a variable with the name attr
vect = np.array([],float)
# loop over all galaxies
for name in names:
vect = np.append( vect, data[name][attr_name] )
odata[attr_name] = vect
return odata
def SetAxis(ax,xmin, xmax, ymin, ymax, log=None):
"""
Set ticks for the axis
"""
#####################################
# first : slightly extend the limits
#####################################
if log is not None:
if str.find(log, 'x') != -1:
if xmin is not None:
xmin = np.log10(xmin)
if xmax is not None:
xmax = np.log10(xmax)
if log is not None:
if str.find(log, 'y') != -1:
if ymin is not None:
ymin = np.log10(ymin)
if ymax is not None:
ymax = np.log10(ymax)
if xmin is not None and xmax is not None:
if xmin == xmax:
xmin = xmin - 0.05 * xmin
xmax = xmax + 0.05 * xmax
else:
xmin = xmin - 0.05 * (xmax - xmin)
xmax = xmax + 0.05 * (xmax - xmin)
if ymin is not None and ymax is not None:
if ymin == ymax:
ymin = ymin - 0.05 * ymin
ymax = ymax + 0.05 * ymax
else:
ymin = ymin - 0.05 * (ymax - ymin)
ymax = ymax + 0.05 * (ymax - ymin)
if log is not None:
if str.find(log, 'x') != -1:
xmin = 10**xmin
xmax = 10**xmax
if log is not None:
if str.find(log, 'y') != -1:
ymin = 10**ymin
ymax = 10**ymax
#####################################
# second : set log log or lin log
#####################################
if log is not None:
if str.find(log, 'x') != -1 and str.find(log, 'y') != -1:
ax.loglog()
elif str.find(log, 'x') != -1:
ax.semilogx()
else:
ax.semilogy()
plt.axis([xmin, xmax, ymin, ymax])
if log is None:
log = 'None'
#####################################
# third : adapt ticks
#####################################
if str.find(log, 'x') == -1:
ax.xaxis.set_major_locator(plt.AutoLocator())
x_major = ax.xaxis.get_majorticklocs()
dx_minor = (x_major[-1] - x_major[0]) / (len(x_major) - 1) / 5.
ax.xaxis.set_minor_locator(plt.MultipleLocator(dx_minor))
if str.find(log, 'y') == -1:
ax.yaxis.set_major_locator(plt.AutoLocator())
y_major = ax.yaxis.get_majorticklocs()
dy_minor = (y_major[-1] - y_major[0]) / (len(y_major) - 1) / 5.
ax.yaxis.set_minor_locator(plt.MultipleLocator(dy_minor))
####################################################################
# main
####################################################################
if __name__ == '__main__':
opt = parser.parse_args()
params = {
"axes.labelsize": 14,
"axes.titlesize": 18,
"font.size": 12,
"legend.fontsize": 12,
"xtick.labelsize": 14,
"ytick.labelsize": 14,
"text.usetex": True,
"figure.subplot.left": 0.15,
"figure.subplot.right": 0.95,
"figure.subplot.bottom": 0.15,
"figure.subplot.top": 0.95,
"figure.subplot.wspace": 0.02,
"figure.subplot.hspace": 0.02,
"figure.figsize" : (8, 6),
"lines.markersize": 6,
"lines.linewidth": 2.0,
}
plt.rcParams.update(params)
# create the plot
fig = plt.gcf()
fig.set_size_inches(8,6)
ax = plt.gca()
# loop over files
for filename in opt.files:
with open(filename, 'rb') as file:
database = pickle.load(file)
xs = []
ys = []
for key in database.keys():
x = database[key]["Lv"]
y = database[key]["R12"]
xs.append(x)
ys.append(y)
ax.scatter(xs,ys)
############################################################################
# add data
############################################################################
if opt.data is not None:
some_names = ["WLM", "Sextans (I)"]
data = readData(opt.data)
data = data.selectc(data["L"] != np.array(None))
data = data.selectc(data["Rh"] != np.array(None))
data = data.selectc(data["L"]>1e5)
data = data.selectc(data["Rh"]<5)
host = data["host"]
data_MW = data.selectc(host=="MW")
data_M31 = data.selectc(host=="M31")
############################
# MW
xn = "L"
yn = "Rh"
data = data_MW
data = data.selectc(data[yn]!=np.array([None]))
names = data["name"]
xs = data[xn]
exp = data["e%sp"%xn]
exm = data["e%sm"%xn]
ys = data[yn]
eyp = data["e%sp"%yn]
eym = data["e%sm"%yn]
# clean none
exp = np.where(exp==np.array(None),0,exp)
exm = np.where(exm==np.array(None),0,exm)
eyp = np.where(eyp==np.array(None),0,eyp)
eym = np.where(eym==np.array(None),0,eym)
ax.errorbar(x=xs,y=ys,xerr=[exm,exp],yerr=[eym,eyp],fmt='s',c='k', alpha=0.1,label=r"$\rm{Milky\,Way\,\,dwarfs}$")
############################
# M31
xn = "L"
yn = "Rh"
data = data_M31
data = data.selectc(data[yn]!=np.array([None]))
names = data["name"]
xs = data[xn]
exp = data["e%sp"%xn]
exm = data["e%sm"%xn]
ys = data[yn]
eyp = data["e%sp"%yn]
eym = data["e%sm"%yn]
# clean none
exp = np.where(exp==np.array(None),0,exp)
exm = np.where(exm==np.array(None),0,exm)
eyp = np.where(eyp==np.array(None),0,eyp)
eym = np.where(eym==np.array(None),0,eym)
ax.errorbar(x=xs,y=ys,xerr=[exm,exp],yerr=[eym,eyp],fmt='v',c='k',alpha=0.1,label=r"$\rm{Andromeda\,\,dwarfs}$")
# set the axis
SetAxis(ax,opt.xmin,opt.xmax,opt.ymin,opt.ymax,log=opt.log)
# labels
ax.set_xlabel(r"V-Band Luminosity")
ax.set_ylabel(r"$r_{1/2}$ [kpc]")
# save or display
if opt.outputfilename:
plt.savefig(opt.outputfilename)
else:
plt.show()
c4science · Help