yt.analysis_modules.absorption_spectrum.absorption_spectrum module

AbsorptionSpectrum class and member functions.

class yt.analysis_modules.absorption_spectrum.absorption_spectrum.AbsorptionSpectrum(lambda_min, lambda_max, n_lambda)[source]

Bases: object

Create an absorption spectrum object.

Parameters:
  • lambda_min (float) – lower wavelength bound in angstroms.
  • lambda_max (float) – upper wavelength bound in angstroms.
  • n_lambda (int) – number of wavelength bins.
add_continuum(label, field_name, wavelength, normalization, index)[source]

Add a continuum feature that follows a power-law.

Parameters:
  • label (string) – label for the feature.
  • field_name (string) – field name from ray data for column densities.
  • wavelength (float) – line rest wavelength in angstroms.
  • normalization (float) – the column density normalization.
  • index (float) – the power-law index for the wavelength dependence.
add_line(label, field_name, wavelength, f_value, gamma, atomic_mass, label_threshold=None)[source]

Add an absorption line to the list of lines included in the spectrum.

Parameters:
  • label (string) – label for the line.
  • field_name (string) – field name from ray data for column densities.
  • wavelength (float) – line rest wavelength in angstroms.
  • f_value (float) – line f-value.
  • gamma (float) – line gamma value.
  • atomic_mass (float) – mass of atom in amu.
make_spectrum(input_file, output_file=None, line_list_file=None, output_absorbers_file=None, use_peculiar_velocity=True, subgrid_resolution=10, observing_redshift=0.0, njobs='auto')[source]

Make spectrum from ray data using the line list.

Parameters:
  • input_file (string or dataset) – path to input ray data or a loaded ray dataset
  • output_file (optional, string) – Option to save a file containing the wavelength, flux, and optical depth fields. File formats are chosen based on the filename extension. .h5 for hdf5, .fits for fits, and everything else is ASCII. Default: None
  • output_absorbers_file (optional, string) – Option to save a text file containing all of the absorbers and corresponding wavelength and redshift information. For parallel jobs, combining the lines lists can be slow so it is recommended to set to None in such circumstances. Default: None
  • use_peculiar_velocity (optional, bool) – if True, include peculiar velocity for calculating doppler redshift to shift lines. Requires similar flag to be set in LightRay generation. Default: True
  • subgrid_resolution (optional, int) – When a line is being added that is unresolved (ie its thermal width is less than the spectral bin width), the voigt profile of the line is deposited into an array of virtual wavelength bins at higher resolution. The optical depth from these virtual bins is integrated and then added to the coarser spectral wavelength bin. The subgrid_resolution value determines the ratio between the thermal width and the bin width of the virtual bins. Increasing this value yields smaller virtual bins, which increases accuracy, but is more expensive. A value of 10 yields accuracy to the 4th significant digit in tau. Default: 10
  • observing_redshift (optional, float) – This is the redshift at which the observer is observing the absorption spectrum. Default: 0
  • njobs (optional, int or "auto") – the number of process groups into which the loop over absorption lines will be divided. If set to -1, each absorption line will be deposited by exactly one processor. If njobs is set to a value less than the total number of available processors (N), then the deposition of an individual line will be parallelized over (N / njobs) processors. If set to “auto”, it will first try to parallelize over the list of lines and only parallelize the line deposition if there are more processors than lines. This is the optimal strategy for parallelizing spectrum generation. Default: “auto”