phantom icon indicating copy to clipboard operation
phantom copied to clipboard

python reader for native file format

Open danieljprice opened this issue 2 years ago • 2 comments

this is posted here as a reference, currently the supported solution for this is to use libphantom to call the read_dump routine directly, or to convert to another format (e.g. using splash or phantom2hdf5)

A better solution would be to have a pure python reader for the native format. A preliminary version of this was kindly contributed by Lionel Siess and Ward Homan, based on work by Stéven Toupin. Just needs to be cleaned up and put into a separate repository (i.e. as a standalone python package)

This issue is just a reminder that this needs to be done, and to flag that this exists for anyone trying to do this.

danieljprice avatar Oct 15 '21 01:10 danieljprice

The file scripts/readPhantomDump.py is now available in the master branch. It works most of the time. Note that handling hdf5 files has never been tested. Here is a simple example on how to use it :

from readPhantomDump import *
from pylab import *

fdump='/my/path/to/run_00041'
data3D = read_dump(fdump)
fpar='/my/path/to/run.in'
data_in=read_infile(fpar)

tmax=data_in['tmax']
npart = data3D['part_numbers']['npartoftype'][0]
udist, umass, utime = [data3D['units'][u] for u in 'udist', 'umass', 'utime']
hfact = data3D['quantities']['hfact']
mpart = data3D['quantities']['massoftype'][0]*umass
x, y, z, vx, vy, vz, h, u = [array(data3D['blocks'][0]['data'][c],dtype=float)*conv for c, conv in [('x', udist), ('y', udist), ('z', udist), ('vx', udist/utime), ('vy', udist/utime), ('vz', udist/utime), ('h', udist), ('u', udist**2/utime**2)]]
r = sqrt(x**2+y**2+z**2)
v = sqrt(vx**2+vy**2+vz**2)
rho = mpart * (hfact/h)**3

# Plot velocity
figure(1)
clf()
plot(r,v/1.e5, '.', color='b', markersize=1.592, rasterized=True, label='3D')
tight_layout()
show()

Cheers Lionel

lsiess avatar Dec 02 '21 08:12 lsiess

Hi @lsiess, just to let you know: for reading Phantom HDF5 files there's also Plonk. If you use it and have any issues, you're welcome to raise over at the Plonk issue tracker.

dmentipl avatar Dec 02 '21 10:12 dmentipl