yt.analysis_modules.photon_simulator.photon_simulator module

Classes for generating lists of photons and detected events

The SciPy Proceeding that describes this module in detail may be found at:

http://conference.scipy.org/proceedings/scipy2014/zuhone.html

The algorithms used here are based off of the method used by the PHOX code (http://www.mpa-garching.mpg.de/~kdolag/Phox/), developed by Veronica Biffi and Klaus Dolag. References for PHOX may be found at:

Biffi, V., Dolag, K., Bohringer, H., & Lemson, G. 2012, MNRAS, 420, 3545 http://adsabs.harvard.edu/abs/2012MNRAS.420.3545B

Biffi, V., Dolag, K., Bohringer, H. 2013, MNRAS, 428, 1395 http://adsabs.harvard.edu/abs/2013MNRAS.428.1395B

class yt.analysis_modules.photon_simulator.photon_simulator.EventList(events, parameters)[source]

Bases: object

filter_events(region)[source]

Filter events using a ds9 region. Requires the pyregion package. Returns a new EventList.

classmethod from_fits_file(fitsfile)[source]

Initialize an EventList from a FITS file with filename fitsfile.

classmethod from_h5_file(h5file)[source]

Initialize an EventList from a HDF5 file with filename h5file.

has_key(key)[source]
items()[source]
keys()[source]
values()[source]
write_fits_file(fitsfile, clobber=False)[source]

Write events to a FITS binary table file with filename fitsfile. Set clobber to True if you need to overwrite a previous file.

write_fits_image(imagefile, clobber=False, emin=None, emax=None)[source]

Generate a image by binning X-ray counts and write it to a FITS file.

Parameters:
  • imagefile (string) – The name of the image file to write.
  • clobber (boolean, optional) – Set to True to overwrite a previous file.
  • emin (float, optional) – The minimum energy of the photons to put in the image, in keV.
  • emax (float, optional) – The maximum energy of the photons to put in the image, in keV.
write_h5_file(h5file)[source]

Write an EventList to the HDF5 file given by h5file.

write_simput_file(prefix, clobber=False, emin=None, emax=None)[source]

Write events to a SIMPUT file that may be read by the SIMX instrument simulator.

Parameters:
  • prefix (string) – The filename prefix.
  • clobber (boolean, optional) – Set to True to overwrite previous files.
  • e_min (float, optional) – The minimum energy of the photons to save in keV.
  • e_max (float, optional) – The maximum energy of the photons to save in keV.
write_spectrum(specfile, bin_type='channel', emin=0.1, emax=10.0, nchan=2000, clobber=False, energy_bins=False)[source]

Bin event energies into a spectrum and write it to a FITS binary table. Can bin on energy or channel. In that case, the spectral binning will be determined by the RMF binning.

Parameters:
  • specfile (string) – The name of the FITS file to be written.
  • bin_type (string, optional) – Bin on “energy” or “channel”. If an RMF is detected, channel information will be imported from it.
  • emin (float, optional) – The minimum energy of the spectral bins in keV. Only used if binning without an RMF.
  • emax (float, optional) – The maximum energy of the spectral bins in keV. Only used if binning without an RMF.
  • nchan (integer, optional) – The number of channels. Only used if binning without an RMF.
  • energy_bins (boolean, optional) – Bin on energy or channel. Deprecated in favor of bin_type.
class yt.analysis_modules.photon_simulator.photon_simulator.PhotonList(photons, parameters, cosmo, p_bins)[source]

Bases: object

classmethod from_file(filename)[source]

Initialize a PhotonList from the HDF5 file filename.

classmethod from_scratch(data_source, redshift, area, exp_time, photon_model, parameters=None, center=None, dist=None, cosmology=None)[source]

Initialize a PhotonList from a photon model. The redshift, collecting area, exposure time, and cosmology are stored in the parameters dictionary which is passed to the photon_model function.

Parameters:
  • data_source (yt.data_objects.data_containers.YTSelectionContainer) – The data source from which the photons will be generated.
  • redshift (float) – The cosmological redshift for the photons.
  • area (float) – The collecting area to determine the number of photons in cm^2.
  • exp_time (float) – The exposure time to determine the number of photons in seconds.
  • photon_model (function) – A function that takes the data_source and the parameters dictionary and returns a photons dictionary. Must be of the form: photon_model(data_source, parameters)
  • parameters (dict, optional) – A dictionary of parameters to be passed to the user function.
  • center (string or array_like, optional) – The origin of the photons. Accepts “c”, “max”, or a coordinate.
  • dist (tuple, optional) – The angular diameter distance in the form (value, unit), used mainly for nearby sources. This may be optionally supplied instead of it being determined from the redshift and given cosmology.
  • cosmology (yt.utilities.cosmology.Cosmology, optional) – Cosmological information. If not supplied, we try to get the cosmology from the dataset. Otherwise, LambdaCDM with the default yt parameters is assumed.

Examples

This is the simplest possible example, where we call the built-in thermal model:

>>> thermal_model = ThermalPhotonModel(apec_model, Zmet=0.3)
>>> redshift = 0.05
>>> area = 6000.0
>>> time = 2.0e5
>>> sp = ds.sphere("c", (500., "kpc"))
>>> my_photons = PhotonList.from_user_model(sp, redshift, area,
...                                         time, thermal_model)

If you wish to make your own photon model function, it must take as its arguments the data_source and the parameters dictionary. However you determine them, the photons dict needs to have the following items, corresponding to cells which have photons:

“x” : the x-position of the cell relative to the source center in kpc, YTArray “y” : the y-position of the cell relative to the source center in kpc, YTArray “z” : the z-position of the cell relative to the source center in kpc, YTArray “vx” : the x-velocity of the cell in km/s, YTArray “vy” : the y-velocity of the cell in km/s, YTArray “vz” : the z-velocity of the cell in km/s, YTArray “dx” : the width of the cell in kpc, YTArray “NumberOfPhotons” : the number of photons in the cell, NumPy array of unsigned 64-bit integers “Energy” : the source rest-frame energies of the photons, YTArray

The last array is not the same size as the others because it contains the energies in all of the cells in a single 1-D array. The first photons[“NumberOfPhotons”][0] elements are for the first cell, the next photons[“NumberOfPhotons”][1] are for the second cell, and so on.

The following is a simple example where a point source with a single line emission spectrum of photons is created. More complicated examples which actually create photons based on the fields in the dataset could be created.

>>> import numpy as np
>>> import yt
>>> from yt.analysis_modules.photon_simulator.api import PhotonList
>>> def line_func(source, parameters):
...
...     ds = source.ds
...
...     num_photons = parameters["num_photons"]
...     E0  = parameters["line_energy"] # Energies are in keV
...     sigE = parameters["line_sigma"]
...     src_ctr = parameters["center"]
...
...     energies = norm.rvs(loc=E0, scale=sigE, size=num_photons)
...
...     # Place everything in the center cell
...     for i, ax in enumerate("xyz"):
...         photons[ax] = (ds.domain_center[0]-src_ctr[0]).in_units("kpc")
...     photons["vx"] = ds.arr([0], "km/s")
...     photons["vy"] = ds.arr([0], "km/s")
...     photons["vz"] = ds.arr([100.0], "km/s")
...     photons["dx"] = ds.find_field_values_at_point("dx", ds.domain_center).in_units("kpc")
...     photons["NumberOfPhotons"] = np.array(num_photons*np.ones(1), dtype="uint64")
...     photons["Energy"] = ds.arr(energies, "keV")
>>>
>>> redshift = 0.05
>>> area = 6000.0
>>> time = 2.0e5
>>> parameters = {"num_photons" : 10000, "line_energy" : 5.0,
...               "line_sigma" : 0.1}
>>> ddims = (128,128,128)
>>> random_data = {"density":(np.random.random(ddims),"g/cm**3")}
>>> ds = yt.load_uniform_grid(random_data, ddims)
>>> dd = ds.all_data()
>>> my_photons = PhotonList.from_user_model(dd, redshift, area,
...                                         time, line_func,
...                                         parameters=parameters)
items()[source]
keys()[source]
project_photons(normal, area_new=None, exp_time_new=None, redshift_new=None, dist_new=None, absorb_model=None, psf_sigma=None, sky_center=None, responses=None, convolve_energies=False, no_shifting=False, north_vector=None, prng=<module 'numpy.random' from '/home/fido/.local/lib64/python3.6/site-packages/numpy/random/__init__.py'>)[source]

Projects photons onto an image plane given a line of sight.

Parameters:
  • normal (character or array_like) – Normal vector to the plane of projection. If “x”, “y”, or “z”, will assume to be along that axis (and will probably be faster). Otherwise, should be an off-axis normal vector, e.g [1.0,2.0,-3.0]
  • area_new (float, optional) – New value for the effective area of the detector. If responses are specified the value of this keyword is ignored.
  • exp_time_new (float, optional) – The new value for the exposure time.
  • redshift_new (float, optional) – The new value for the cosmological redshift.
  • dist_new (tuple, optional) – The new value for the angular diameter distance in the form (value, unit), used mainly for nearby sources. This may be optionally supplied instead of it being determined from the cosmology.
  • absorb_model ('yt.analysis_modules.photon_simulator.PhotonModel`, optional) – A model for galactic absorption.
  • psf_sigma (float, optional) – Quick-and-dirty psf simulation using Gaussian smoothing with standard deviation psf_sigma in degrees.
  • sky_center (array_like, optional) – Center RA, Dec of the events in degrees.
  • responses (list of strings, optional) – The names of the ARF and/or RMF files to convolve the photons with.
  • convolve_energies (boolean, optional) – If this is set, the photon energies will be convolved with the RMF.
  • no_shifting (boolean, optional) – If set, the photon energies will not be Doppler shifted.
  • north_vector (a sequence of floats) – A vector defining the ‘up’ direction. This option sets the orientation of the plane of projection. If not set, an arbitrary grid-aligned north_vector is chosen. Ignored in the case where a particular axis (e.g., “x”, “y”, or “z”) is explicitly specified.
  • prng (NumPy RandomState object or numpy.random) – A pseudo-random number generator. Typically will only be specified if you have a reason to generate the same set of random numbers, such as for a test. Default is the numpy.random module.

Examples

>>> L = np.array([0.1,-0.2,0.3])
>>> events = my_photons.project_photons(L, area_new="sim_arf.fits",
...                                     redshift_new=0.05,
...                                     psf_sigma=0.01)
values()[source]
write_h5_file(photonfile)[source]

Write the photons to the HDF5 file photonfile.

yt.analysis_modules.photon_simulator.photon_simulator.convert_old_file(input_file, output_file, clobber=False)[source]

Helper function for converting old PhotonList or EventList HDF5 files (pre yt v3.3) to their new versions.

Parameters:
  • input_file (list of strings) – The filename of the old-versioned file to be converted.
  • output_file (string) – Name of the new file to be outputted.
  • clobber (boolean, default False) – If a the output file already exists, set this to True to overwrite it.

Examples

>>> from yt.analysis_modules.photon_simulator.api import convert_old_file
>>> convert_old_file("photons_old.h5", "photons_new.h5", clobber=True)
yt.analysis_modules.photon_simulator.photon_simulator.force_unicode(value)[source]
yt.analysis_modules.photon_simulator.photon_simulator.merge_files(input_files, output_file, clobber=False, add_exposure_times=False)[source]

Helper function for merging PhotonList or EventList HDF5 files.

Parameters:
  • input_files (list of strings) – List of filenames that will be merged together.
  • output_file (string) – Name of the merged file to be outputted.
  • clobber (boolean, default False) – If a the output file already exists, set this to True to overwrite it.
  • add_exposure_times (boolean, default False) – If set to True, exposure times will be added together. Otherwise, the exposure times of all of the files must be the same.

Examples

>>> from yt.analysis_modules.photon_simulator.api import merge_files
>>> merge_files(["events_0.h5","events_1.h5","events_3.h5"], "events.h5",
...             clobber=True, add_exposure_times=True)

Notes

Currently, to merge files it is mandated that all of the parameters have the same values, with the possible exception of the exposure time parameter “exp_time” if add_exposure_times=False.

yt.analysis_modules.photon_simulator.photon_simulator.parse_value(value, default_units)[source]
yt.analysis_modules.photon_simulator.photon_simulator.validate_parameters(first, second, skip=[])[source]