trident.spectrum_generator.SpectrumGenerator

class trident.spectrum_generator.SpectrumGenerator(instrument=None, lambda_min=None, lambda_max=None, n_lambda=None, dlambda=None, lsf_kernel=None, line_database='lines.txt', ionization_table=None, bin_space='wavelength')[source]

Preferred class for generating, storing, and plotting absorption-line spectra. SpectrumGenerator is a subclass of yt’s AbsorptionSpectrum class with additional functionality like line lists, post-processing to include Milky Way foreground, quasar backgrounds, applying line-spread functions, and plotting.

User first specifies the telescope/instrument used for generating spectra (e.g. ‘COS’ for the Cosmic Origins Spectrograph aboard the Hubble Space Telescope). This can be done by naming the Instrument, or manually setting all of the spectral characteristics including lambda_min, lambda_max, lsf_kernel, and n_lambda or dlambda. If none of these arguments are set, defaults to ‘COS’ as the default instrument covering 1150-1450 Angstroms with a binsize (dlambda) of 0.1 Angstroms.

Once a SpectrumGenerator has been initialized, pass it LightRay objects using make_spectrum to actually generate the spectra themselves. Then one can post-process, plot, or save them using add_milky_way_foreground, add_qso_spectrum, apply_lsf, save_spectrum, and plot_spectrum.

Parameters

Instrument

string, optional

The spectrograph to use. Currently, the only named options are different observing modes of the Cosmic Origins Spectrograph ‘COS’: ‘COS-G130M’, ‘COS-G160M’, and ‘COS-G140L’ as well as ‘COS’ which aliases to ‘COS-G130M’. These automatically set the lambda_min, lambda_max, dlambda and lsf_kernel``s appropriately. If you're going to set ``lambda_min, lambda_max, et al manually, leave this set to None. Default: None

Lambda_min

float, YTQuantity, or ‘auto’

lower wavelength bound in angstroms or velocity bound in km/s (if bin_space set to ‘velocity’). If set to ‘auto’, the lower bound will be automatically adjusted to encompass all absorption lines. The window will not be expanded for continuum features, only absorption lines.

Lambda_max

float, YTQuantity, or ‘auto’

upper wavelength bound in angstroms or velocity bound in km/s (if bin_space set to ‘velocity’). If set to ‘auto’, the upper bound will be automatically adjusted to encompass all absorption lines. The window will not be expanded for continuum features, only absorption lines.

N_lambda

int

The number of bins in the spectrum (inclusive), so if extrema = 10 and 20, and dlambda (binsize) = 1, then n_lambda = 11. Default: None

Dlambda

float

size of the wavelength bins in angstroms or velocity bins in km/s. Default: None

Bin_space

‘wavelength’ or ‘velocity’

Sets the dimension in which spectra are created. If set to wavelength, the resulting spectra are flux (or tau) vs. observed wavelength. If set to velocity, the spectra are flux vs. velocity offset from the rest wavelength of the absorption line. Default: wavelength

Lsf_kernel

string, optional

The filename for the LSF kernel. Files are found in trident.__path__/data/lsf_kernels or current working directory. Only necessary if you are applying an LSF to the spectrum in postprocessing. Default: None

Line_database

string or LineDatabase, optional

A text file listing the various lines to insert into the line database, or a LineDatabase object in memory. The line database provides a list of all possible lines that could be added to the spectrum. For a text file, it should have 4 tab-delimited columns of name (e.g. MgII), wavelength in angstroms, gamma of transition, and f-value of transition. See example datasets in trident.path/data/line_lists for examples. Default: lines.txt

Ionization_table

hdf5 file, optional

An HDF5 file used for computing the ionization fraction of the gas based on its density, temperature, metallicity, and redshift. If set to None, will use the ion table defined in your Trident config file. Default: None

Example

Create a one-zone ray, and generate a COS spectrum from that ray.

>>> import trident
>>> ray = trident.make_onezone_ray()
>>> sg = trident.SpectrumGenerator('COS')
>>> sg.make_spectrum(ray)
>>> sg.plot_spectrum('spec_raw.png')

Create a one-zone ray at redshift 0.5, and generate a spectrum with 1 angstrom spectral bins from 2000-4000 angstroms, then post-process by adding a MW foreground a QSO background at z=0.5 and add a boxcar line spread function of 100 angstroms width. Plot it and save the figure to ‘spec_final.png’.

>>> import trident
>>> ray = trident.make_onezone_ray(redshift=0.5)
>>> sg = trident.SpectrumGenerator(lambda_min=2000, lambda_max=4000,
... dlambda=1)
>>> sg.make_spectrum(ray)
>>> sg.add_qso_spectrum(emitting_redshift=.5)
>>> sg.add_milky_way_foreground()
>>> sg.apply_lsf(function='boxcar', width=100)
>>> sg.plot_spectrum('spec_final.png')

Methods

__init__

add_continuum

Add a continuum feature that follows a power-law.

add_gaussian_noise

Postprocess a spectrum to add gaussian random noise of a given SNR.

add_line

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

add_line_to_database

Adds desired line to the current LineDatabase object.

add_milky_way_foreground

Postprocess a spectrum to add a Milky Way foreground.

add_noise_vector

Add an array of noise to the spectrum.

add_qso_spectrum

Postprocess a spectrum to add a QSO spectrum background.

apply_lsf

Postprocess a spectrum to apply a line spread function.

clear_spectrum

Clear the current spectrum in the SpectrumGenerator.

error_func

Approximate the flux error for a spectrum.

load_spectrum

Load data arrays into an existing spectrum object.

make_spectrum

Make a spectrum from ray data depositing the desired lines.

plot_spectrum

Plot the current spectrum and save to disk.

save_spectrum

Save the current spectrum data to an output file.

Attributes

current_tau_field

This is the optical depth array for the current absorption line being deposited.

lambda_field

The lambda field.

tau_field

This is the total optical depth of all lines and continua.