Page MenuHomec4science
Files F83474416

erkal2swift
No OneTemporary
Actions

File Metadata

Created
Tue, Sep 17, 08:23

erkal2swift
View Options

#!/usr/bin/env python3
import argparse
import numpy as np
from pNbody import *
from astropy import units as u
####################################################################
# option parser
####################################################################
description="convert an ascii file to an hdf5 swift one"
epilog ="""
Examples:
--------
erkal2swift filein.txt --mass 2.5e8 -o fileout.hdf5
"""
parser = argparse.ArgumentParser(description=description,epilog=epilog,formatter_class=argparse.RawDescriptionHelpFormatter)
parser.add_argument(action="store",
dest="file",
metavar='FILE',
type=str,
default=None,
nargs=1,
help='input file')
parser.add_argument("--mass",action="store",
type=float,
default=1e10,
nargs=1,
help='total model mass [Msol]')
parser.add_argument("-o",
action="store",
type=str,
dest="outputfilename",
default=None,
help="Name of the output file")
####################################################################
# main
####################################################################
if __name__ == '__main__':
opt = parser.parse_args()
# read ascii file
data = np.loadtxt(opt.file[0])
# pos [kpc]
x = data[:,0]
y = data[:,1]
z = data[:,2]
pos = np.transpose((x,y,z))
# vel [km/s]
vx = data[:,3]
vy = data[:,4]
vz = data[:,5]
vel = np.transpose((vx,vy,vz))
# ids
num = data[:,6].astype(int)
# mass
n = len(x)
mass = np.ones(n) * opt.mass/n
mass = mass/1e10 # to 1e10 Msol
# create the pNbody object
nb = Nbody(pos=pos,vel=vel,mass=mass,num=num,status="new",ftype="swift")
# set particles to be stars
md = nb.getParticleMatchingDict()
nb.set_tpe(md["stars"])
# 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()

Event Timeline

c4science · Help