Analysis of Phantom output
Visualisation of Phantom output
Dump files
That’s what splash is for! Use splash to interactively inspect the dump files produced by phantom, e.g.:
splash dump_0*
or for a column density rendering:
splash -r 6 dump_0*
and press Enter for “Hollywood mode”. You can also use splash to convert to ascii:
splash to ascii dump_0*
or to interpolate to a 3D grid:
splash to grid dump_0*
To make a movie, just give “/png” as the output:
splash -r 6 dump_0* -dev /png
Then use the ffmpeg script in the splash/scripts directory to convert the png files into an mp4 movie:
~/splash/scripts/movie.sh mymovie.mp4 splash
More details in the splash documentation
For more detailed analysis of Phantom dump files, write yourself an analysis module for the phantomanalysis utility as described below. Analysis modules exist for many common tasks, including interpolating to a 3D grid (both fixed and AMR), computing PDFs, structure functions and power spectra, getting disc surface density profiles, and converting to other formats, and it is simple to write your own.
Global quantities
The .ev files, which are just ascii files containing global quantities as a function of time, can be visualised using asplash with the -ev (or -e) option:
splash -ev dump*.ev
Global quantities not in the .ev file can also be obtained using the splash calc utility, e.g.:
splash calc max dump_0*
which produces a file containing the maximum of each quantity in the dump files as a function of time.
Customised analysis
Reading phantom dump files
The short answer is, do not under any circumstances attempt to do this yourself! (If you need convincing, just have a quick look at how long the read_data_sphNG.f90 file in splash is). The best way to read/analyse phantom dumps, aside from using splash to visualise them, is to use the built-in phantomanalysis utility (described below), or the sarracen python package. A full description of the data format and how to read it can be found here.
Sarracen
Phantomanalysis
Compile the phantomanalysis utility using:
make analysis
which compiles the phantomanalysis binary using the analysis module you specified in build/Makefile_setups:
ifeq ($(SETUP), isodisc)
...
SETUPFILE= setup_disc.f90
ANALYSIS= analysis_disc.f90
...
giving
$ ls
phantomanalysis*
which you can then run on a series of snapshots
$ ./phantomanalysis dump_0*
Phantomanalysis is a simple wrapper that reads all of the dump files on the command line in sequence and calls the analysis routine specified in the ANALYSIS variable, in this case analysis_disc.f90. For a list of pre-built analysis tools, see the list of Phantom utilities.
Compiling your own phantomanalysis module
You can also select the module on the command line using, for example
make analysis ANALYSIS=analysis_blah.f90
You can then write an analysis_blah.f90 to do whatever it is you want, even if you what you want is:
something completely trivial (see for example analysis_dtheader.f90 which just compares the time from each dump file with the time in the previous dump file); or
conversion to another format; or
actually performing some analysis (e.g. analysis_disc.f90 which bins particles into rings for comparison with 1D alpha-disc evolution calculations).
The call to analysis passes the most useful information on the particles (positions, velocities, thermal energy, particle masses and numbers of particles). Any remaining information can also be accessed via the usual phantom modules. For example, you can access sink particle arrays using:
use part, only:xyzmh_ptmass,vxyz_ptmass
Converting to another format
Apart from writing a short analysis module, you can also use the convert utility in splash. For example, to convert all files to ascii format (not recommended, they’ll be huge):
splash to ascii blast_0*
You can also convert to other code formats, e.g.:
splash to gadget blast_0*
Analysis with pyanalysis
An alternative method for analysis in python is to compile the phantom pyanalysis utility using:
make pyanalysis
which compiles the libphantom library giving
$ ls
libphantom.so* libanalysis.py
Now you can import the PhantomAnalysis class from libanalysis.py. This can be done interactively in iPython or in a Python script
In [1]: from libanalysis import PhantomAnalysis as pa
and create an instance of this class with a phantom dumpfile
In [2]: dumpfile = 'blast_00000'
In [3]: dump = pa(dumpfile)
This loads the dumpfile and places particle quantities into numpy arrays. These quantities are accessible as attributes of the PhantomAnalysis class. For example
In [4]: print dump.npart
125000
In [5]: print dump.xyzh
[[-0.49 -0.47 -0.45 ..., 0.45 0.47 0.49 ]
[-0.49 -0.49 -0.49 ..., 0.49 0.49 0.49 ]
[-0.49 -0.49 -0.49 ..., 0.49 0.49 0.49 ]
[ 0.024 0.024 0.024 ..., 0.024 0.024 0.024]]
In [6]: print dump.vxyz
[[ 0. 0. 0. ..., 0. 0. 0.]
[ 0. 0. 0. ..., 0. 0. 0.]
[ 0. 0. 0. ..., 0. 0. 0.]]
In [7]: print dump.utherm
[ 0. 0. 0. ..., 0. 0. 0.]
List of variables
time
hfact
massofgas
units (dictionary) {‘udist’, ‘umass’, ‘utime’, ‘udens’, ‘umagfd’}
npart
xyzh
vxyz
utherm
nptmass
ptmass_xyzmh
ptmass_vxyz
ptmass_spinxyz
Loading phantom HDF5 dumps into python
To get yourself HDF5 dumpfiles, have a look at Running phantom with HDF5 output.
Import h5py and load the dumpfile
In [1]: import h5py
In [2]: f = h5py.File('disc_00000.h5','r')
List the main containers in the file
In [3]: list(f.keys())
Out[3]: ['header', 'particles', 'sinks']
List the particle arrays that are available
In [4]: list(f['particles'].keys())
Out[4]: ['divv', 'dt', 'h', 'itype', 'pressure', 'vxyz', 'xyz']
Extract the xyz array from the file
In [5]: f['particles']['xyz'].value
Out[5]:
array([[ -6.05266606, -6.66164664, -0.34922808],
[ 2.55540523, 17.91264485, 0.52264339],
[ 15.26729989, -6.75512839, -0.70489168],
...,
[ -9.45331138, 1.34188609, 0.69513828],
[ 12.67824199, 3.35761305, -0.39397658],
[-11.34601204, 0.75837632, 0.6858956 ]])
See h5py docs for more information