Page MenuHomec4science
Files F100640131

gorbit
No OneTemporary
Actions

File Metadata

Created
Sat, Feb 1, 10:05

gorbit
View Options

#!/home/revaz/local/bin/python
'''
Extract orbit from
a set of nbody files.
default file extension = .orb
Yves Revaz
mar jui 11 15:54:41 CEST 2006
'''
import sys,os
from numarray import *
import SM
from Nbody import *
from Nbody import io
try:
from optparse import OptionParser
except ImportError:
from optik import OptionParser
def parse_options():
usage = "usage: %prog [options] file"
parser = OptionParser(usage=usage)
parser.add_option("-t",
action="store",
dest="ftype",
type="string",
default = 'gadget',
help="file type",
metavar=" FTYPE")
parser.add_option("-n",
action="store",
dest="number",
type="int",
default = 0,
help="particle number",
metavar=" INT")
parser.add_option("--numfile",
action="store",
dest="numfile",
type="string",
default = None,
help="file of particle number",
metavar=" FILE")
(options, args) = parser.parse_args()
return args,options
##########################################################
#
# MAIN
#
#########################################################
# get options
files,options = parse_options()
ftype = options.ftype
number = options.number
numfile = options.numfile
#############################
# read numfile if it exists
#############################
if numfile != None:
numbers = io.read_ascii(numfile,[0])
numbers = numbers[0].astype(Int)
else:
numbers = array([number],Int)
#############################
# open output files
#############################
fs = []
ofiles = []
for number in numbers:
ofile = "%08d.orb"%(number)
#f = open(ofile,'w')
#fs.append(f)
if os.path.exists(ofile):
os.remove(ofile)
ofiles.append(ofile)
#############################
# loops
#############################
for file in files:
nb = Nbody(file,ftype=ftype)
nb = nb.selectp(numbers)
i=0
for number in numbers:
index = argmax(nb.num==number)
x = nb.pos[index,0]
y = nb.pos[index,1]
z = nb.pos[index,2]
vx = nb.vel[index,0]
vy = nb.vel[index,1]
vz = nb.vel[index,2]
print "%s, %10.5f %10.5f %10.5f %10.5f %10.5f %10.5f %10.5f"%(os.path.basename(file),x,y,z,vx,vy,vz,nb.tnow)
#fs[i].write('%10.5f %10.5f %10.5f %10.5f %10.5f %10.5f %10.5f\n'%(x,y,z,vx,vy,vz,nb.tnow))
#fs[i].flush()
f = open(ofiles[i],'a')
f.write('%10.5f %10.5f %10.5f %10.5f %10.5f %10.5f %10.5f\n'%(x,y,z,vx,vy,vz,nb.tnow))
f.close()
i = i + 1
#############################
# close output files
#############################
#for f in fs:
# f.close()

Event Timeline

c4science · Help