Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F83460562
ark_idealized_dwarfs_generate
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
Tue, Sep 17, 06:34
Size
6 KB
Mime Type
text/x-python
Expires
Thu, Sep 19, 06:34 (2 d)
Engine
blob
Format
Raw Data
Handle
20841879
Attached To
rARRAKIHS ARRAKIHS
ark_idealized_dwarfs_generate
View Options
#!/usr/bin/env python3
import numpy as np
from pNbody import *
from pNbody.mass_models import plummer
from pNbody import ic
from astropy import constants as cte
from astropy import units as u
import argparse
from scipy.stats import skewnorm
from scipy.stats import rv_continuous
from scipy.stats.sampling import NumericalInversePolynomial
from tqdm import tqdm
from pNbody import pychem
from pNbody import Isochrones
import matplotlib.pyplot as plt
import matplotlib.colors as mcolors
import copy
import pickle
from astropy.io import fits
from pNbody import Mockimgs
from pNbody.Mockimgs.obs_object import Object
from pNbody.Mockimgs import utils
from pNbody.Mockimgs.los import LineOfSights
from pNbody.Mockimgs import filters
from pNbody.Mockimgs import luminosities
####################################################################
# option parser
####################################################################
description="""generate a dwarf galaxy
"""
epilog ="""
Examples:
--------
ark_idealized_dwarfs_generate --Fe -2.5 --age 3 --M0 1e4 -e 0.1 -o 1e4.hdf5
"""
parser = argparse.ArgumentParser(description=description,epilog=epilog,formatter_class=argparse.RawDescriptionHelpFormatter)
parser.add_argument("-o",
action="store",
type=str,
dest="outputfilename",
default=None,
help="Name of the output file")
parser.add_argument("--AgeMean",
action="store",
dest="AgeMean",
metavar='FLOAT',
type=float,
default=1,
help='mean Age')
parser.add_argument("--AgeFlat",
action="store_true",
default=False,
help="flatten age")
parser.add_argument("--Fe",
action="store",
dest="Fe",
metavar='FLOAT',
type=float,
default=-2.5,
help='Fe mode')
parser.add_argument("--M0",
action="store",
dest="M0",
metavar='FLOAT',
type=float,
default=1e5,
help='stellar mass to generate [Msol]')
parser.add_argument("-e",
action="store",
dest="e",
metavar='FLOAT',
type=float,
default=0.5,
help='scale length [kpc]')
parser.add_argument("--rmax_e_ratio",
action="store",
dest="rmax_e_ratio",
metavar='FLOAT',
type=float,
default=10,
help='rmax to scale length ratio')
parser.add_argument("--ftype",
action="store",
type=str,
dest="ftype",
default="swift",
help="type of file")
parser.add_argument("--do_not_compute_rsp",
action="store_true",
default=False,
help="do not compute rsp")
parser.add_argument("--nngb",
action="store",
type=int,
dest="nngb",
default=5,
help="Number of neighbouring particles to consider to compute RSP(==HSML)")
parser.add_argument("-N",
action="store",
dest="N",
metavar='INT',
type=int,
default=10000,
help='number of particles')
class AgesDF():
def __init__(self,AgeMean=2,AgeMax=13.7,flat=False):
self.AgeMax = AgeMax
self.AgeMean = AgeMean
self.s = AgeMean/AgeMax/2
self.flat = flat
self.x = np.linspace(0,1,10000)
self.y = self.DF(self.x)
self.yc = np.add.accumulate(self.y)
self.yc = self.yc/max(self.yc)
def DF(self,x):
'''
max pour x=s
max = 1
'''
n=1
y = x**n * np.exp(-x**n/self.s**n) #/ (self.s**2*np.exp(-1))
if self.flat:
xpeak = self.x[np.argmax(y)]
ypeak = np.max(y)
y = np.where(self.x>=xpeak,ypeak,y)
return y
def Sample(self,n):
X = np.random.random(n)
x = np.interp(X,self.yc,self.x)
return x*self.AgeMax
####################################################################
# main
####################################################################
def Do(opt):
n = opt.N
#############################################
# set particle mass
#############################################
mass = np.ones(n)*opt.M0/n
mass = mass/1e10 # output units are in 1e10 Msol
#############################################
# set ages
#############################################
df = AgesDF(AgeMean=opt.AgeMean,flat=opt.AgeFlat)
Ages = df.AgeMax - df.Sample(n)
Ages.sort()
print('Ages done')
#############################################
# generate metallicities
#############################################
Fes = opt.Fe + np.log10(-np.log(1-np.random.rand(n)))
Fes.sort()
Fes = np.flip(Fes)
#############################################
# generate a plummer model
#############################################
# set maximal radius
rmax = opt.e * opt.rmax_e_ratio
# generate the model
nb = ic.plummer(n,1,1,1,opt.e,rmax,M=1,irand=1,vel='no',name=opt.outputfilename,ftype=opt.ftype)
# compute Hsml
if opt.do_not_compute_rsp is False:
print("Compute Rsp...")
# remove doublets
#u,idx = np.unique(nb.rxyz(),return_index=True)
#nb = nb.selectp(lst=nb.num[idx])
nb.ComputeRsp(opt.nngb)
print("done.")
# set particles to be stars
md = nb.getParticleMatchingDict()
nb.set_tpe(md["stars"])
# set masses (as mass or luminosities)
nb.mass = mass
# store age
if Ages is not None:
nb.age = Ages # store in Gyr
# store metallicity
if Fes is not None:
nb.mh = Fes
# define units
u_Length = 1* u.kpc
u_Mass = 10**10 * u.M_sun
u_Velocity = 1* u.km/u.s
u_Time = u_Length/u_Velocity
toMsol = u_Mass.to(u.M_sun).value
# add units
nb.UnitLength_in_cm = u_Length.to(u.cm).value
nb.UnitMass_in_g = u_Mass.to(u.g).value
nb.UnitVelocity_in_cm_per_s = u_Velocity.to(u.cm/u.s).value
nb.Unit_time_in_cgs = u_Time.to(u.s).value
# non cosmological
nb.setComovingIntegrationOff()
nb.cosmorun=0
nb.time =0
if opt.outputfilename is not None:
nb.rename(opt.outputfilename)
nb.write()
if __name__ == '__main__':
opt = parser.parse_args()
Do(opt)
Event Timeline
Log In to Comment