-
Notifications
You must be signed in to change notification settings - Fork 24
reading and processing outputs
The py_progs folder contains some programs for reading outputs. Add the py_progs folder to your $PYTHONPATH. you might want to import the following modules:
import py_read_output as rd
import py_plot_output as plot
import pylab as p
import py_plot_util as util
To read ASCII outputs, Python mostly uses the astropy.io.ascii module. To install astropy type
pip install astropy
or visit the astropy website
With the new (as of November 2014) file formats which are consistent with astropy, one can read a spectrum, or spec_tot file with the following commands
from astropy.io import ascii
spectrum = ascii.read(filename)
one can then a spectrum with the column headings, e.g.
p.plot(spectrum["Lambda"], spectrum["A10P0.50"])
and view the column headings with
print s.colnames
Alternatively, use the py_read_output module to read into an astropy.table.table.Table object
spectrum = rd.read_spectrum(filename)
one can then plot a spectrum as above, or use py_plot_output to make a standard plot from an astropy.table.table.Table object:
plot.make_spec_plot(spectrum, filename)
To read a wind file, you first need to have compiled py_wind. You can then run this from a python shell to generate a complete summary of wind properties by typing
util.get_pywind_summary(filename)
you can then read the data from filename.complete into an astropy Table
data = rd.read_pywind(filename)
and reshape into a masked array for a given colname, e.g. ne or lum_lines
x,z,ne = util.wind_to_masked (data, "ne")
contour plots can then be created with
p.contourf(z,x,np.log10(ne)) # note that contourf takes z as first argument
Alternatively, running the following will make a standard set of contourf plots for the user
plot.make_wind_plot(data, filename, var = ["ne", "te", "tr", "IP", "nphot", "v", "w", "ionc4"], shape=(4,2) )
and if data is set to Nonetype then the program will actually run py_wind and read in filename.complete for you.
pf = rd.read_pf("sv")
View the disk mdot rate
print pf['disk.mdot(msol/yr)']