Photon Lists API

Classes for generating lists of photons

class pyxsim.photon_list.PhotonList(filespec)[source]
write_spectrum(specfile, emin, emax, nchan, overwrite=False)[source]

Bin photon energies into a spectrum and write it to a FITS binary table. This is for an unconvolved spectrum.

Parameters:
  • specfile (string) – The name of the FITS file to be written.

  • emin (float) – The minimum energy of the spectral bins in keV.

  • emax (float) – The maximum energy of the spectral bins in keV.

  • nchan (integer) – The number of channels.

  • overwrite (boolean, optional) – Set to True to overwrite a previous file.

pyxsim.photon_list.make_photons(photon_prefix, data_source, redshift, area, exp_time, source_model, point_sources=False, parameters=None, center=None, dist=None, cosmology=None, velocity_fields=None, bulk_velocity=None, observer='external', fields_to_keep=None)[source]

Write a photon list dataset to disk from a yt data source and assuming a source model for the photons. The redshift, collecting area, exposure time, and cosmology are stored in the parameters dictionary which is passed to the source_model function.

Parameters:
  • photon_prefix (string) – The prefix of the filename(s) to contain the photon list. If run in serial, the filename will be “{photon_prefix}.h5”, if run in parallel, the filenames will be “{photon_prefix}.{mpi_rank}.h5”.

  • data_source (YTSelectionContainer) – The data source from which the photons will be generated.

  • redshift (float) – The cosmological redshift for the photons.

  • area (float, (value, unit) tuple, YTQuantity, or Quantity) – The collecting area to determine the number of photons. If units are not specified, it is assumed to be in cm^2.

  • exp_time (float, (value, unit) tuple, YTQuantity, or Quantity) – The exposure time to determine the number of photons. If units are not specified, it is assumed to be in seconds.

  • source_model (SourceModel) – A source model used to generate the photons.

  • point_sources (boolean, optional) – If True, the photons will be assumed to be generated from the exact positions of the cells or particles and not smeared around within a volume. Default: False

  • parameters (dict, optional) – A dictionary of parameters to be passed for the source model to use, if necessary.

  • center (string or array_like, optional) – The origin of the photon spatial coordinates. Accepts “c”, “max”, or a coordinate. If array-like and without units, it is assumed to be in units of kpc. If not specified, pyxsim attempts to use the “center” field parameter of the data_source.

  • dist (float, (value, unit) tuple, YTQuantity, or Quantity, optional) – The angular diameter distance, used for nearby sources. This may be optionally supplied instead of it being determined from the redshift and given cosmology. If units are not specified, it is assumed to be in kpc. To use this, the redshift must be set to zero.

  • cosmology (Cosmology, optional) – Cosmological information. If not supplied, we try to get the cosmology from the dataset. Otherwise, LCDM with the default yt parameters is assumed.

  • velocity_fields (list of fields) – The yt fields to use for the velocity. If not specified, the following will be assumed: [‘velocity_x’, ‘velocity_y’, ‘velocity_z’] for grid datasets [‘particle_velocity_x’, ‘particle_velocity_y’, ‘particle_velocity_z’] for particle datasets

  • bulk_velocity (array-like, optional) – A 3-element array or list specifying the local velocity frame of reference. If not a YTArray, it is assumed to have units of km/s. Default: [0.0, 0.0, 0.0] km/s.

  • observer (string, optional) – Controls whether the source is at a large distance, in which case the observer sits far outside the source (“external”), or if the observer is inside the X-ray-emitting source (“internal”), such as a galaxy. Default: “external”

  • fields_to_keep (list of tuples) – A list of fields to add to the HDF5 photon list in addition to the cell or particle positions, velocities, and sizes. Default: None

Returns:

  • A tuple of two integers, the number of photons and the number of cells with

  • photons

Examples

>>> thermal_model = pyxsim.CIESourceModel("apec", 0.1, 10.0, 1000, 0.3)
>>> redshift = 0.05
>>> area = 6000.0 # assumed here in cm**2
>>> exp_time = 2.0e5 # assumed here in seconds
>>> sp = ds.sphere("c", (500., "kpc"))
>>> fields_to_keep = [("gas", "density"), ("gas", "temperature")]
>>> n_photons, n_cells = pyxsim.make_photons(
...     sp, redshift, area, exp_time, thermal_model,
...     fields_to_keep=fields_to_keep
... )
pyxsim.photon_list.project_photons(photon_prefix, event_prefix, normal, sky_center, absorb_model=None, nH=None, abund_table='angr', no_shifting=False, north_vector=None, sigma_pos=None, flat_sky=False, kernel='top_hat', save_los=False, prng=None)[source]

Projects photons onto an image plane given a line of sight, and stores them in an HDF5 dataset which contains an event list.

Parameters:
  • photon_prefix (string) – The prefix of the filename(s) containing the photon list. If run in serial, the filename will be “{photon_prefix}.h5”, if run in parallel, the filenames will be “{photon_prefix}.{mpi_rank}.h5”.

  • event_prefix (string) – The prefix of the filename(s) which will be written to contain the event list. If run in serial, the filename will be “{event_prefix}.h5”, if run in parallel, the filename will be “{event_prefix}.{mpi_rank}.h5”.

  • 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]

  • sky_center (array-like) – Center RA, Dec of the events in degrees.

  • absorb_model (string) – A model for foreground galactic absorption, to simulate the absorption of events before being detected. Known options for are “wabs” and “tbabs”.

  • nH (float, optional) – The foreground column density in units of 10^22 cm^{-2}. Only used if absorption is applied.

  • abund_table (string) – The abundance table to be used for abundances in the absorption model (only used for TBabs). Default is set in the SOXS configuration file, the default for which is “angr”. Built-in options are: “angr” : from Anders E. & Grevesse N. (1989, Geochimica et Cosmochimica Acta 53, 197) “aspl” : from Asplund M., Grevesse N., Sauval A.J. & Scott P. (2009, ARAA, 47, 481) “feld” : from Feldman U. (1992, Physica Scripta, 46, 202) “wilm” : from Wilms, Allen & McCray (2000, ApJ 542, 914 except for elements not listed which are given zero abundance) “lodd” : from Lodders, K (2003, ApJ 591, 1220) “cl17.03” : the default abundance table in Cloudy 17.03

  • no_shifting (boolean, optional) – If set, the photon energies will not be Doppler shifted. Default: False

  • 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 perpendicular to the normal is chosen. Ignored in the case where a particular axis (e.g., “x”, “y”, or “z”) is explicitly specified.

  • sigma_pos (float, optional) – Apply a gaussian smoothing operation to the sky positions of the events. This may be useful when the binned events appear blocky due to their uniform distribution within simulation cells. However, this will move the events away from their originating position on the sky, and so may distort surface brightness profiles and/or spectra. Should probably only be used for visualization purposes. Supply a float here to smooth with a standard deviation with this fraction of the cell size. Default: None

  • flat_sky (boolean, optional) – If set, we assume that the sky is “flat” and RA, Dec positions are computed using simple linear offsets

  • kernel (string, optional) – The kernel used when smoothing positions of X-rays originating from SPH particles, “gaussian” or “top_hat”. Default: “top_hat”.

  • save_los (boolean, optional) – If True, save the line-of-sight positions along the projection axis in units of kpc to the events list. Default: False

  • prng (integer or RandomState object) – 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 to use the numpy.random module.

Return type:

A integer for the number of events created

Examples

>>> L = np.array([0.1,-0.2,0.3])
>>> n_events = pyxsim.project_photons("my_photons.h5", "my_events.h5", L,
...                                   [30., 45.], absorb_model='tbabs',
...                                   nH=0.04)
pyxsim.photon_list.project_photons_allsky(photon_prefix, event_prefix, normal, absorb_model=None, nH=None, abund_table='angr', no_shifting=False, center_vector=None, kernel='top_hat', save_los=False, prng=None)[source]

Projects photons onto the sky sphere given a normal vector (“z” or “up” in spherical coordinates), and stores them in an HDF5 dataset which contains an event list.

Parameters:
  • photon_prefix (string) – The prefix of the filename(s) containing the photon list. If run in serial, the filename will be “{photon_prefix}.h5”, if run in parallel, the filenames will be “{photon_prefix}.{mpi_rank}.h5”.

  • event_prefix (string) – The prefix of the filename(s) which will be written to contain the event list. If run in serial, the filename will be “{event_prefix}.h5”, if run in parallel, the filename will be “{event_prefix}.{mpi_rank}.h5”.

  • normal (array-like) – The vector determining the “z” or “up” vector for the spherical coordinate system for the all-sky projection, something like [1.0, 2.0, -3.0]. It will be normalized before use.

  • absorb_model (string) – A model for foreground galactic absorption, to simulate the absorption of events before being detected. Known options are “wabs” and “tbabs”.

  • nH (float, optional) – The foreground column density in units of 10^22 cm^{-2}. Only used if absorption is applied.

  • abund_table (string) – The abundance table to be used for abundances in the absorption model (only used for TBabs). Default is set in the SOXS configuration file, the default for which is “angr”. Built-in options are: “angr” : from Anders E. & Grevesse N. (1989, Geochimica et Cosmochimica Acta 53, 197) “aspl” : from Asplund M., Grevesse N., Sauval A.J. & Scott P. (2009, ARAA, 47, 481) “feld” : from Feldman U. (1992, Physica Scripta, 46, 202) “wilm” : from Wilms, Allen & McCray (2000, ApJ 542, 914 except for elements not listed which are given zero abundance) “lodd” : from Lodders, K (2003, ApJ 591, 1220) “cl17.03” : the default abundance table in Cloudy 17.03

  • no_shifting (boolean, optional) – If set, the photon energies will not be Doppler shifted. Default: False

  • center_vector (a sequence of floats) – A vector defining what direction will be placed at the center of the lat/lon coordinate system. If not set, an arbitrary grid-aligned center_vector perpendicular to the normal is chosen.

  • kernel (string, optional) – The kernel used when smoothing positions of X-rays originating from SPH particles, “gaussian” or “top_hat”. Default: “top_hat”.

  • save_los (boolean, optional) – If True, save the line-of-sight radii in units of kpc to the events file. Default: False

  • prng (integer or RandomState object) – 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 to use the numpy.random module.

Return type:

A integer for the number of events created

Examples

>>> L = np.array([0.1,-0.2,0.3])
>>> n_events = pyxsim.project_photons_allsky("my_photons.h5", "my_events.h5", L)
pyxsim.utils.merge_files(input_files, output_file, overwrite=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.

  • overwrite (boolean, default False) – If 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 pyxsim import merge_files
>>> merge_files(["events_0.h5","events_1.h5","events_3.h5"], "events.h5",
...             overwrite=True, add_exposure_times=True)

Notes

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