PyAtomDB Spectrum module

This module contains codes for creating spectra from the AtomDB emissivity files

This module contains methods for creating spectra from the AtomDB files. Some are more primitive than others…

class pyatomdb.spectrum.CIESession(linefile='$ATOMDB/apec_line.fits', cocofile='$ATOMDB/apec_coco.fits', elements=[1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30], abundset='AG89')

Bases: object

Load and generate a collisional ionization equilibrium spectrum

Parameters:
linefile : string or HDUList, optional

The line emissivity data file (either name or already open)

cocofile : string or HDUList, optional

The continuum emissivity data file (either name or already open)

elements : arraylike(int), optional

The atomic number of elements to include (default all)

abundset : string

The abundance set to use. Default AG89.

Examples

Create a session instance:

>>> s=CIESession()

Set up the responses, in this case a dummy response from 0.1 to 10 keV, with area 1cm^2 in each bin

>>> ebins = numpy.linspace(0.1,10,1000)
>>> s.set_response(ebins, raw=True)

(Alternatively, for a real response file, s.set_response(rmffile, arf=arffile)

Turn on thermal broadening

>>> s.set_broadening(True)
Will thermally broaden lines with emissivity > 1.000000e-18 ph cm3 s-1

Return spectrum at 1.0keV

>>> spec = s.return_spectrum(1.0)

spec is in photons cm^5 s^-1 bin^-1; ebins are the bin edges (so spec is 1 element shorter than ebins)

Attributes:
datacache : dict

Any Atomdb FITS files which have to be opened are stored here

spectra : CIESpectra

Object storing the actual spectral data

elements : list(int)

Nuclear charge of elements to include.

default_abundset : string

The abundance set used for the original emissivity file calculation

abundset : string

The abundance set to be used for the returned spectrum

abundsetvector : array_like(float)

The relative abundance between default_abundset and abundset for each element

response_set : bool

Have we loaded a response (or set a dummy response)

dolines : bool

Calculate line emission

docont : bool

Calculate continuum emission

dopseudo : bool

Calculate pseudocontinuum emission

broaden_limit : float

Apply broadening to lines with epsilon > this value (ph cm3 s-1)

thermal_broadening : bool

Apply thermal broadening to lines (default = False)

velocity_broadening : float

Apply velocity broadening with this velocity (km/s). If <=0, do not apply.

apply_response(self, spectrum)

Apply a response to a spectrum

Parameters:
spectrum : array(float)

The spectrum, in counts/bin/second, to have the response applied to. Must be binned on the same grid as the rmf.

Returns
——-
array(float)

spectrum folded through the response

return_line_emissivity(self, Te, Z, z1, up, lo, specunit='A', teunit='keV', apply_aeff=False, apply_abund=True, log_interp=True)

Get line emissivity as function of Te.

Parameters:
Te : float or array(float)

Temperature(s) in keV or K

Z : int

nuclear charge of element

z1 : int

ion charge +1 of ion

up : int

upper level for transition

lo : int

lower level for transition

specunit : {‘Angstrom’,’keV’}

Units for wavelength or energy (a returned value)

teunit : {‘keV’ , ‘K’}

Units of Te (kev or K, default keV)

apply_aeff : bool

If true, apply the effective area to the line emissivity in the linelist to modify their intensities.

apply_abund : bool

If true, apply the abundance set in the session to the result.

log_interp : bool

Interpolate between temperature on a log-log scale (default). Otherwise linear

Returns:
ret : dict

Dictionary containing: Te, tau, teunit: as input wavelength : line wavelength (A) energy : line energy (keV) epsilon : emissivity in ph cm^3 s-1 (or ph cm^5 s^-1 if apply_aeff=True)

first index is temperature, second is tau.

return_linelist(self, Te, specrange, specunit='A', teunit='keV', apply_aeff=False, nearest=False, apply_binwidth=False)

Get the list of line emissivities vs wavelengths

Parameters:
Te : float

Temperature in keV or K

specrange : [float, float]

Minimum and maximum values for interval in which to search

specunit : {‘Angstrom’,’keV’}

Units for specrange

teunit : {‘keV’ , ‘K’}

Units of te (kev or K, default keV)

apply_aeff : bool

If true, apply the effective area to the lines in the linelist to modify their intensities.

Returns:
linelist : array(dtype)

The list of lines with lambda (A), energy (keV), epsilon (ph cm3 s-1), epsilon_aeff (ph cm5 s-1) ion (string) and upper & lower levels.

return_spectrum(self, te, teunit='keV', nearest=False, get_nearest_t=False)

Get the spectrum at an exact temperature. Interpolates between 2 neighbouring spectra

Finds HDU with kT closest to desired kT in given line or coco file.

Opens the line or coco file, and looks for the header unit with temperature closest to te. Use result as index input to make_spectrum

Parameters:
te : float

Temperature in keV or K

teunit : {‘keV’ , ‘K’}

Units of te (kev or K, default keV)

raw : bool

If set, return the spectrum without response applied. Default False.

nearest : bool

If set, return the spectrum from the nearest tabulated temperature in the file, without interpolation

get_nearest_t : bool

If set, and nearest set, return the nearest tabulated temperature as well as the spectrum.

Returns:
spectrum : array(float)

The spectrum in photons cm^5 s^-1 bin^-1, with the response, or photons cm^3 s^-1 bin^-1 if raw is set.

nearest_T : float, optional

If nearest is set, return the actual temperature this corresponds to. Units are same as teunit

set_abund(self, elements, abund)

Set the elemental abundance, relative to the abundset. Defaults to 1.0 for everything

Parameters:
elements : int or array_like(int)

The elements to change the abundance of

abund : float or array_like(float)

The new abundances. If only 1 value, set all elements to this abundance Otherwise, should be of same length as elements.

Returns:
None

Examples

Set the abundance of iron to 0.5

>>> myspec.set_abund(26, 0.5)

Set the abundance of iron and nickel to 0.1 and 0.2 respectively

>>> myspec.set_abund([26, 28], [0.1,0.2])

Set the abundance of oxygen, neon, magnesium and iron to 0.1

>>> myspec.set_abund([8,10,12,26],0.1)
set_abundset(self, abundstring)

Set the abundance set.

Parameters:
abundstring : string

The abundance string (e.g. “AG89”, “uniform”. Case insensitive. See atomdb.get_abundance for list of possible abundances

Returns:
none

updates self.abundset and self.abundsetvector.

set_apec_files(self, linefile='$ATOMDB/apec_line.fits', cocofile='$ATOMDB/apec_coco.fits')

Set the apec line and coco files, and load up their data

Parameters:
linefile : str or HDUList

The filename of the line emissivity data, or the opened file.

cocofile : str or HDUList

The filename of the continuum emissivity data, or the opened file.

Returns:
None

Notes

Updates self.linefile, self.linedata, self.cocofile and self.cocodata

set_broadening(self, thermal_broadening, broaden_limit=False, velocity_broadening=0.0, velocity_broadening_units='km/s')

Turn on or off thermal broadening, and the emissivity limit for thost lines

Parameters:
thermal_broadening : bool

If true, turn on broadening. If False, turn it off.

broaden_limit : float

The emissivity limit for lines to be broadened. If False, this value will not be updated.

velocity_broadening : float

velocity broadening to apply. If <=0, not applied

velocity_broadening_units : string

Units of velocity_broadening. ‘km/s’ is default and only value so far.

Notes

Updates attributes thermal_broadening, broaden_limit, velocity_broadening, velocity_broadening_units

set_response(self, rmf, arf=False, raw=False)

Set the response. rmf, arf can either be the filenames or the opened files (latter is faster if called repeatedly)

Amends the following items:

self.rmffile : string
The rmf file name
self.rmf : string
The response matrix
self.arffile : string
The arf file name
self.arf : string
The arf data
self.specbins : array(float)
The spectral bins on which to calculate the spectrum (keV or A)
self.specbin_units : string [‘A’,’keV’]
Units of specbins
self.ebins : array(float)
The spectral bins on which to calculate the spectrum (keV).
self.ebins_out : array(float)
The spectral bins on which to return the spectrum (keV). Can be different from specbins depending on the spectrum
self.response_set : bool
A resonse has been loaded
self.specbins_set : bool
The spectral bins are set
self.ebins_checksum : string
The md5checksum of the specbins
Parameters:
rmf: string or HDUlist or array

The response matrix file or energy bins (see raw k/w)

arf: string or HDUlist

The ancillary response file

raw : bool

If true, the rmf variable contains the energy bin edges (keV) for all the bins, and each bin has a perfect response. This is effectively a dummy response

Returns:
none
class pyatomdb.spectrum.CIESpectrum(linedata, cocodata)

Bases: object

A class holding the emissivity data for CIE emission, and returning spectra

Parameters:
linefile : string or HDUList, optional

The line emissivity data file (either name or already open)

cocofile : string or HDUList, optional

The continuum emissivity data file (either name or already open)

elements : arraylike(int), optional

The atomic number of elements to include (default all)

abundset : string

The abundance set to use. Default AG89.

Attributes:
session : CIESession

The parent CIESession

SessionType : string

“CIE”

spectra : dict of ElementSpectra

a dictionary containing the emissivity data for each HDU, subdivided by element (spectra[12][18] is an ElementSpectrum object containing the argon data for the 12th HDU)

kTlist : array

The temperatures for each emissivity HDU, in keV

logkTlist : array

log of kTlist

get_nearest_Tindex(self, Te, teunit='keV', nearest=False, log_interp=True)

Return the nearest temperature index in the emissivity file, or, alternatively, the array of fractions to sum

Parameters:
Te : float

Temperature (keV by default)

teunit : string

Units of kT

nearest : bool

If true, return only nearest. Otherwise, return nearest 2 and fractions

Returns:
ikT : list[int]

Index of temperature in HDU file (from 0, not 2)

f : list[float]

fractional weight to apply to each ikT. Should sum to 1.

return_line_emissivity(self, Te, Z, z1, up, lo, specunit='A', teunit='keV', abundance=1.0, log_interp=True)

Return the emissivity of a line at kT, tau. Assumes ionization from neutral for now

Parameters:
Te : float

Temperature in keV or K

Z : int

nuclear charge of element

z1 : int

ion charge +1 of ion

up : int

upper level for transition

lo : int

lower level for transition

specunit : {‘Angstrom’,’keV’}

Units for wavelength or energy (a returned value)

teunit : {‘keV’ , ‘K’}

Units of Telist (kev or K, default keV)

abundance : float

Abundance to multiply the emissivity by

Returns:
Emissivity : float

Emissivity in photons cm^3 s^-1

spec : float

Wavelength or Energy of line, depending on specunit

return_linelist(self, Te, teunit='keV', nearest=False, specrange=False, specunit='A', elements=False, abundances=False)

Return the linelist of the element

Parameters:
Te : float

Electron temperature (default, keV)

teunit : string

Units of kT (keV by default, K also allowed)

nearest : bool

If True, return spectrum for the nearest temperature index. If False, use the weighted average of the (log of) the 2 nearest indexes. default is False.

specrange : [float, float]

Minimum and maximum values for interval in which to search

specunit : {‘Ansgtrom’,’keV’}

Units for specrange (default A)

return_spectrum(self, Te, teunit='keV', nearest=False, elements=False, abundances=False, log_interp=True, broaden_object=False)

Return the spectrum of the element on the energy bins in self.session.specbins

Parameters:
Te : float

Electron temperature (default, keV)

teunit : string

Units of kT (keV by default, K also allowed)

nearest : bool

If True, return spectrum for the nearest temperature index. If False, use the weighted average of the (log of) the 2 nearest indexes. default is False.

Returns:
spec : array(float)

The element’s emissivity spectrum, in photons cm^3 s^-1 bin^-1

class pyatomdb.spectrum.ContinuumData(cocoentry, docont=True, dopseudo=True)

Bases: object

A class holding the continuum data for an element in one HDU

Parameters:
cocoentry : array(cocodatatype)

A single row from the continuum data in an AtomDB file.

parentElementSpectrum : ElementSpectrum

Parent ElementSpectrum object

Attributes:
ECont : array(float)

The continuum energies (keV)

EPseudo : array(float)

The pseudocontinuum energies (keV)

Cont : array(float)

The continuum emissivities (ph cm^3 s^-1 keV^-1)

Pseudo : array(float)

The pseudocontinuum energies (ph cm^3 s^-1 keV^-1)

parentElementSpectrum : ElementSpectrum

Parent ElementSpectrum object

spectrum_calculated : bool

True if spectrum has already been calculated, otherwise false

ebins_checksum : string

md5sum of the ebins the spectrum was last calculated on. Used to identify if new calculations are required or can just return the previous value.

docont : bool

Calculate the continuum emission

dopseudo : bool

Calculate the pseudocontinuum emission

return_spec(self, eedges, ebins_checksum=False)
class pyatomdb.spectrum.ElementSpectrum(linedata, cocodata, Z, z1_drv=0)

Bases: object

A class holding the emissivity data for an element in one HDU

Parameters:
linedata : array(linedatatype)

array of line wavelengths and emissivities, from AtomDB files. Should already be filtered to only be from one element.

cocodata : array(cocodatatype)

array of continuum wavelengths and emissivities, from AtomDB files Should already be filtered to only be from one element.

Z : int

The atomic number of the element

z1_drv : int

The charge + 1 for the ion. 0 = whole element.

parent : CIESpectrum

Parent CIESpectrum object

Attributes:
lines : LineData

A LineData object containing all the line information

continuum : ContinuumData

A ContinuumData object containing all the contrinuum information

parent : CIESpectrum

Parent CIESpectrum object

session : CIESession

Parent Session of parent CIESpectrum object

return_linelist(self, specrange, specunit='A', teunit='keV')

Return the list of lines in specrange

Parameters:
specrange : array

spectral range to look for lines in

specunits : string

units of spectral range (A or keV)

teunit : string

units of Te (keV, eV, K)

Returns:
linelist : array

list of lines and epsilons

return_linematch(self, Z, z1, up, lo, z1_drv=0)

Return the line(s) which match the transition

Parameters:
Z : int

nuclear charge

z1 : int

ion charge + 1

up : int

upper level of transition

lo : int

lower level of transition

z1_drv : int

if provided, also filter on driving ion charge (NEI only)

Returns
——-
linelist : array

list of lines and epsilons

return_spectrum(self, eedges, Te, ebins_checksum=False, thermal_broadening=False, broaden_limit=False, velocity_broadening=0.0, teunit='keV', broaden_object=False)

Calculate the spectrum

Parameters:
eedges : array

bin edges (keV)

Te : float

temperature (keV by defult)

ebins_checksum : string

the md5 checksum of eedges

thermal_broadening : bool

true to apply thermal broadening

broaden_limit : float

only broaden lines stronger than this.

velocity_broadening : float

velocity broadening to apply, km/s. Set <=0 for none (default)

teunit : string

Temperature unit (K, keV, eV)

Returns:
spectrum : array(float)

The spectrum in ph cm^3 s^-1 bin^-1

class pyatomdb.spectrum.Gaussian_CDF

Bases: object

For fast interpolation, pre-calculate the CDF and interpolate it when broadening lines

Parameters:
None

Examples

Create a CDF instance:

>>> s=Gaussian_CDF()

Broaden a line on ebins grid, with centroid and width.

>>> cdf = a.broaden(centroid, width, ebins)

Convert to flux in each bin

>>> flux= cdf[1:]-cdf[:-1]
broaden(self, centroid, width, ebins)

Broaden a line, return CDF

Parameters:
centroid : float

The line energy (keV)

width : float

The sigma of the normal distribution, in keV

ebins : array(float)

Energy grid to return CDF on

Returns:
CDF : array(float)

cumulative flux distribution of linen at each bin edge.

class pyatomdb.spectrum.LineData(linelist)

Bases: object

A class holding the line data for an element in one HDU

Parameters:
linelist : array(linedatatype)

array of line wavelengths and emissivities, from AtomDB files.

Attributes:
lines : array(linedatatype)

List of lines, wavelength and emissivities

lineenergies : array(float)

list of line energies

spectrum_calculated : bool

True if spectrum has already been calculated, otherwise false

T: float

Temperature (K) last spectrum was calculated at (for broadening)

v : float

Velocity (km/s) last spectrum was calculated at (for broadening)

ebins_checksum : string

md5sum of the ebins the spectrum was last calculated on. Used to identify if new calculations are required or can just return the previous value.

return_spec(self, eedges, T, ebins_checksum=False, thermal_broadening=False, velocity_broadening=0.0, broaden_limit=1e-18, broaden_object=False)

return the line emission spectrum at tempterature T

Parameters:
eedges : array

energy bin edges, keV

T : float

temperature in Kelvin

ebins_checksum : string

the md5 checksum of eedges

thermal_broadening : bool

true to apply thermal broadening

velocity_broadening : float

velocity broadening to apply, km/s. Set <=0 for none (default)

broaden_limit : float

only broaden lines stronger than this.

Returns:
spectrum : array(float)

Emissivity on eedges spectral bins of the lines, in ph cm^3 s^-1 bin^-1

class pyatomdb.spectrum.NEISession(linefile='$ATOMDB/apec_nei_line.fits', cocofile='$ATOMDB/apec_nei_comp.fits', elements=False, abundset='AG89')

Bases: pyatomdb.spectrum.CIESession

Load and generate a collisional ionization equilibrium spectrum

Parameters:
linefile : string or HDUList, optional

The line emissivity data file (either name or already open)

cocofile : string or HDUList, optional

The continuum emissivity data file (either name or already open)

elements : arraylike(int), optional

The atomic number of elements to include (default all)

abundset : string

The abundance set to use. Default AG89.

Examples

Create a session instance:

>>> s=CIESession()

Set up the responses, in this case a dummy response from 0.1 to 10 keV

>>> ebins = numpy.linspace(0.1,10,1000)
>>> s.set_response(ebins, raw=True)

Turn on thermal broadening

>>> s.set_broadening(True)
Will thermally broaden lines with emissivity > 1.000000e-18 ph cm3 s-1

Return spectrum at 1.0keV

>>> spec = s.return_spectrum(1.0)

spec is in photons cm^3 s^-1 bin^-1; ebins are the bin edges (so spec is 1 element shorter than ebins)

Attributes:
datacache : dict

Any Atomdb FITS files which have to be opened are stored here

spectra : CIESpectra

Object storing the actual spectral data

elements : list(int)

Nuclear charge of elements to include.

default_abundset : string

The abundance set used for the original emissivity file calculation

abundset : string

The abundance set to be used for the returned spectrum

abundsetvector : array_like(float)

The relative abundance between default_abundset and abundset for each element

response_set : bool

Have we loaded a response (or set a dummy response)

dolines : bool

Calculate line emission

docont : bool

Calculate continuum emission

dopseudo : bool

Calculate pseudocontinuum emission

broaden_limit : float

Apply broadening to lines with epsilon > this value (ph cm3 s-1)

thermal_broadening : bool

Apply thermal broadening to lines (default = False)

velocity_broadening : float

Apply velocity broadening with this velocity (km/s). If <=0, do not apply.

apply_response(self, spectrum)

Apply a response to a spectrum

Parameters:
spectrum : array(float)

The spectrum, in counts/bin/second, to have the response applied to. Must be binned on the same grid as the rmf.

Returns
——-
array(float)

spectrum folded through the response

return_line_emissivity(self, Telist, taulist, Z, z1, up, lo, specunit='A', teunit='keV', apply_aeff=False, apply_abund=True, log_interp=True, init_pop='ionizing')

Get line emissivity as function of Te, tau. Assumes ionization from neutral.

Parameters:
Telist : float or array(float)

Temperature(s) in keV or K

taulist : float

ionization timescale(s), ne * t (cm^-3 s).

Z : int

nuclear charge of element

z1 : int

ion charge +1 of ion

up : int

upper level for transition

lo : int

lower level for transition

specunit : {‘Angstrom’,’keV’}

Units for wavelength or energy (a returned value)

teunit : {‘keV’ , ‘K’}

Units of Telist (kev or K, default keV)

apply_aeff : bool

If true, apply the effective area to the line emissivity in the linelist to modify their intensities.

apply_abund : bool

If true, apply the abundance set in the session to the result.

log_interp : bool

Interpolate between temperature on a log-log scale (default). Otherwise linear

Returns:
ret : dict

Dictionary containing: Te, tau, teunit: as input wavelength : line wavelength (A) energy : line energy (keV) epsilon : emissivity in ph cm^3 s-1 (or ph cm^5 s^-1 if apply_aeff=True)

first index is temperature, second is tau.

return_linelist(self, Te, tau, specrange, specunit='A', teunit='keV', apply_aeff=False)

Get the list of line emissivities vs wavelengths

Parameters:
Te : float

Temperature in keV or K

tau : float

ionization timescale, ne * t (cm^-3 s).

specrange : [float, float]

Minimum and maximum values for interval in which to search

specunit : {‘Angstrom’,’keV’}

Units for specrange

teunit : {‘keV’ , ‘K’}

Units of te (kev or K, default keV)

apply_aeff : bool

If true, apply the effective area to the lines in the linelist to modify their intensities.

Returns:
linelist : array(dtype)

The list of lines with lambda (A), energy (keV), epsilon (ph cm3 s-1), epsilon_aeff (ph cm5 s-1) ion (string) and upper & lower levels.

return_spectrum(self, Te, tau, init_pop=False, Te_init=False, teunit='keV', nearest=False, get_nearest_t=False, log_interp=True)

Get the spectrum at an exact temperature. Interpolates between 2 neighbouring spectra

Finds HDU with kT closest to desired kT in given line or coco file.

Opens the line or coco file, and looks for the header unit with temperature closest to te. Use result as index input to make_spectrum

Parameters:
Te : float

Temperature in keV or K

teunit : {‘keV’ , ‘K’}

Units of te (kev or K, default keV)

raw : bool

If set, return the spectrum without response applied. Default False.

nearest : bool

If set, return the spectrum from the nearest tabulated temperature in the file, without interpolation

get_nearest_t : bool

If set, and nearest set, return the nearest tabulated temperature as well as the spectrum.

Returns:
spectrum : array(float)

The spectrum in photons cm^5 s^-1 bin^-1, with the response, or photons cm^3 s^-1 bin^-1 if raw is set.

nearest_T : float, optional

If nearest is set, return the actual temperature this corresponds to. Units are same as teunit

set_abund(self, elements, abund)

Set the elemental abundance, relative to the abundset. Defaults to 1.0 for everything

Parameters:
elements : int or array_like(int)

The elements to change the abundance of

abund : float or array_like(float)

The new abundances. If only 1 value, set all elements to this abundance Otherwise, should be of same length as elements.

Returns:
None

Examples

Set the abundance of iron to 0.5

>>> myspec.set_abund(26, 0.5)

Set the abundance of iron and nickel to 0.1 and 0.2 respectively

>>> myspec.set_abund([26, 28], [0.1,0.2])

Set the abundance of oxygen, neon, magnesium and iron to 0.1

>>> myspec.set_abund([8,10,12,26],0.1)
set_abundset(self, abundstring)

Set the abundance set.

Parameters:
abundstring : string

The abundance string (e.g. “AG89”, “uniform”. Case insensitive. See atomdb.get_abundance for list of possible abundances

Returns:
none

updates self.abundset and self.abundsetvector.

set_apec_files(self, linefile='$ATOMDB/apec_nei_line.fits', cocofile='$ATOMDB/apec_nei_comp.fits')

Set the apec line and coco files, and load up their data

Parameters:
linefile : str or HDUList

The filename of the line emissivity data, or the opened file.

cocofile : str or HDUList

The filename of the continuum emissivity data, or the opened file.

Returns:
None

Notes

Updates self.linefile, self.linedata, self.cocofile and self.cocodata

set_broadening(self, thermal_broadening, broaden_limit=False, velocity_broadening=0.0, velocity_broadening_units='km/s')

Turn on or off thermal broadening, and the emissivity limit for thost lines

Parameters:
thermal_broadening : bool

If true, turn on broadening. If False, turn it off.

broaden_limit : float

The emissivity limit for lines to be broadened. If False, this value will not be updated.

velocity_broadening : float

velocity broadening to apply. If <=0, not applied

velocity_broadening_units : string

Units of velocity_broadening. ‘km/s’ is default and only value so far.

Notes

Updates attributes thermal_broadening, broaden_limit, velocity_broadening, velocity_broadening_units

set_response(self, rmf, arf=False, raw=False)

Set the response. rmf, arf can either be the filenames or the opened files (latter is faster if called repeatedly)

Amends the following items:

self.rmffile : string
The rmf file name
self.rmf : string
The response matrix
self.arffile : string
The arf file name
self.arf : string
The arf data
self.specbins : array(float)
The spectral bins on which to calculate the spectrum (keV or A)
self.specbin_units : string [‘A’,’keV’]
Units of specbins
self.ebins : array(float)
The spectral bins on which to calculate the spectrum (keV).
self.ebins_out : array(float)
The spectral bins on which to return the spectrum (keV). Can be different from specbins depending on the spectrum
self.response_set : bool
A resonse has been loaded
self.specbins_set : bool
The spectral bins are set
self.ebins_checksum : string
The md5checksum of the specbins
Parameters:
rmf: string or HDUlist or array

The response matrix file or energy bins (see raw k/w)

arf: string or HDUlist

The ancillary response file

raw : bool

If true, the rmf variable contains the energy bin edges (keV) for all the bins, and each bin has a perfect response. This is effectively a dummy response

Returns:
none
class pyatomdb.spectrum.NEISpectrum(linedata, cocodata)

Bases: pyatomdb.spectrum.CIESpectrum

A class holding the emissivity data for NEI emission, and returning spectra

Parameters:
linefile : string or HDUList, optional

The line emissivity data file (either name or already open)

cocofile : string or HDUList, optional

The continuum emissivity data file (either name or already open)

elements : arraylike(int), optional

The atomic number of elements to include (default all)

abundset : string

The abundance set to use. Default AG89.

Attributes:
session : CIESession

The parent CIESession

SessionType : string

“CIE”

spectra : dict of ElementSpectra

a dictionary containing the emissivity data for each HDU, subdivided by element (spectra[12][18] is an ElementSpectrum object containing the argon data for the 12th HDU)

kTlist : array

The temperatures for each emissivity HDU, in keV

logkTlist : array

log of kTlist

get_nearest_Tindex(self, Te, teunit='keV', nearest=False, log_interp=True)

Return the nearest temperature index in the emissivity file, or, alternatively, the array of fractions to sum

Parameters:
Te : float

Temperature (keV by default)

teunit : string

Units of kT

nearest : bool

If true, return only nearest. Otherwise, return nearest 2 and fractions

Returns:
ikT : list[int]

Index of temperature in HDU file (from 0, not 2)

f : list[float]

fractional weight to apply to each ikT. Should sum to 1.

return_line_emissivity(self, Te, tau, Z, z1, up, lo, specunit='A', teunit='keV', abundance=1.0, log_interp=True, init_pop='ionizing')

Return the emissivity of a line at kT, tau. Assumes ionization from neutral for now

Parameters:
Te : float

Temperature in keV or K

tau : float

ionization timescale, ne * t (cm^-3 s).

Z : int

nuclear charge of element

z1 : int

ion charge +1 of ion

up : int

upper level for transition

lo : int

lower level for transition

specunit : {‘Angstrom’,’keV’}

Units for wavelength or energy (a returned value)

teunit : {‘keV’ , ‘K’}

Units of Telist (kev or K, default keV)

abundance : float

Abundance to multiply the emissivity by

init_pop : string or float
If string:

if ‘ionizing’ : all ionizing from neutral (so [1,0,0,0…]) if ‘recombining’: all recombining from ionized (so[…0,0,1]) if array of length (Z+1) : the acutal fractional populations if single float : the temperature (same units as Te)

Returns:
Emissivity : float

Emissivity in photons cm^3 s^-1

spec : float

Wavelength or Energy of line, depending on specunit

return_linelist(self, Te, tau, Te_init=False, teunit='keV', nearest=False, specrange=False, specunit='A', elements=False, abundances=False, init_pop='ionizing')

Return the linelist of the element

Parameters:
Te : float

Electron temperature (default, keV)

teunit : string

Units of kT (keV by default, K also allowed)

nearest : bool

If True, return spectrum for the nearest temperature index. If False, use the weighted average of the (log of) the 2 nearest indexes. default is False.

specrange : [float, float]

Minimum and maximum values for interval in which to search

specunit : {‘Ansgtrom’,’keV’}

Units for specrange (default A)

return_spectrum(self, Te, tau, init_pop=False, Te_init=False, teunit='keV', nearest=False, elements=False, abundances=False, log_interp=True)

Return the spectrum of the element on the energy bins in self.session.specbins

Parameters:
Te : float

Electron temperature (default, keV)

tau : float

ionization timescale, ne * t (cm^-3 s).

teunit : string

Units of kT (keV by default, K also allowed)

nearest : bool

If True, return spectrum for the nearest temperature index. If False, use the weighted average of the (log of) the 2 nearest indexes. default is False.

Returns:
spec : array(float)

The element’s emissivity spectrum, in photons cm^3 s^-1 bin^-1

class pyatomdb.spectrum.TestContinuumData(cocoentry, docont=True, dopseudo=True)

Bases: object

A class holding the continuum data for an element in one HDU

Parameters:
cocoentry : array(cocodatatype)

A single row from the continuum data in an AtomDB file.

Attributes:
ECont : array(float)

The continuum energies (keV)

EPseudo : array(float)

The pseudocontinuum energies (keV)

Cont : array(float)

The continuum emissivities (ph cm^3 s^-1 keV^-1)

Pseudo : array(float)

The pseudocontinuum energies (ph cm^3 s^-1 keV^-1)

spectrum_calculated : bool

True if spectrum has already been calculated, otherwise false

ebins_checksum : string

md5sum of the ebins the spectrum was last calculated on. Used to identify if new calculations are required or can just return the previous value.

docont : bool

Calculate the continuum emission

dopseudo : bool

Calculate the pseudocontinuum emission

return_spec(self, eedges, ebins_checksum=False)
class pyatomdb.spectrum.TestSpectrum(linefile='/export1/atomdb_latest/apec_v3.0.9_line.fits', cocofile='/export1/atomdb_latest/apec_v3.0.9_coco.fits')

Bases: object

makeSpectrum(self, ihdu, ebins)
pyatomdb.spectrum.add_lines(Z, abund, lldat, ebins, z1=False, z1_drv=False, broadening=False, broadenunits='A')

Add lines to spectrum, applying gaussian broadening.

Add the lines in list lldat, with atomic number Z, to a spectrum delineated by ebins (these are the edges, in keV). Apply broadening to the spectrum if broadening != False, with units of broadenunits (so can do constant wavelength or energy broadening)

Parameters:
Z : int

Element of interest (e.g. 6 for carbon)

abund : float

Abundance of element, relative to AG89 data.

lldat : dtype linelist

The linelist to add. Usually the hdu from the apec_line.fits file, often with some filters pre-applied.

ebins : array of floats

Energy bins. Will return spectrum with nbins-1 data points.

z1 : int

Ion charge +1 of ion to return

z1_drv : int

Driving Ion charge +1 of ion to return

broadening : float

Apply spectral broadening if > 0. Units of A of keV

broadenunits : {‘A’ , ‘keV’}

The units of broadening, Angstroms or keV

Returns:
array of float

broadened emissivity spectrum, in photons cm^3 s^-1 bin^-1. Array has len(ebins)-1 values.

pyatomdb.spectrum.apply_response(spectrum, rmf, arf=False)

Apply a response to a spectrum

Parameters:
spectrum : array(float)

The spectrum, in counts/bin/second, to have the response applied to. Must be binned on the same grid as the rmf.

rmf : string or pyfits.hdu.hdulist.HDUList

The filename of the rmf or the opened rmf file

arf : string or pyfits.hdu.hdulist.HDUList

The filename of the arf or the opened arf file

Returns
——-
array(float)

energy grid (keV) for returned spectrum

array(float)

spectrum folded through the response

pyatomdb.spectrum.broaden_continuum(bins, spectrum, binunits='keV', broadening=False, broadenunits='keV')

Apply a broadening to the continuum

Parameters:
bins : array(float)

The bin edges for the spectrum to be calculated on, in units of keV or Angstroms. Must be monotonically increasing. Spectrum will return len(bins)-1 values.

spectrum : array(float)

The emissivities in each bin in the unbroadened spectrum

binunits : {‘keV’ , ‘A’}

The energy units for bins. “keV” or “A”. Default keV.

broadening : float

Broaden the continuum by gaussians of this width (if False, no broadening is applied)

broadenunits : {‘keV’ , ‘A’}

Units for broadening (kev or A)

Returns:
array(float)

spectrum broadened by gaussians of width broadening

pyatomdb.spectrum.convert_spec(spec, specunit, specunitout)

Convert spectral ranges from specunit to specunitout

Parameters:
spec : array

The units to return

specunit : string

The input spectral unit (‘keV’, ‘A’)

specunitout : string

The output spectral unit (‘keV’, ‘A’)

Returns:
specout : array

spec, converted to specunitout

pyatomdb.spectrum.convert_temp(Te, teunit, teunitout)

Convert temperature (Te) from units teunit to teunitout

Parameters:
Te : float

The temperature

teunit : string

units of Te

teunitout : string

output temperature units

pyatomdb.spectrum.expand_E_grid(eedges, n, Econt_in_full, cont_in_full)

Code to expand the compressed continuum onto a series of bins.

Parameters:
eedges : float(array)

The bin edges for the spectrum to be calculated on, in units of keV

n : int

The number of good data points in the continuum array

Econt_in_full: float(array)

The compressed continuum energies

cont_in_full: float(array)

The compressed continuum emissivities

Returns:
float(array)

len(bins)-1 array of continuum emission, in units of photons cm^3 s^-1 bin^-1

pyatomdb.spectrum.get_effective_area(rmf, arf=False)

Get the effective area of a response file

Parameters:
rmf : string or pyfits.hdu.hdulist.HDUList

The filename of the rmf or the opened rmf file

arf : string or pyfits.hdu.hdulist.HDUList

The filename of the arf or the opened arf file

Returns
——-
array(float)

energy grid (keV) for returned response

array(float)

effective area for the returned response

pyatomdb.spectrum.get_index(te, filename='$ATOMDB/apec_line.fits', teunits='keV', logscale=False)

Finds HDU with kT closest ro desired kT in given line or coco file.

Opens the line or coco file, and looks for the header unit with temperature closest to te. Use result as index input to make_spectrum

Parameters:
te : float

Temperature in keV or K

teunits : {‘keV’ , ‘K’}

Units of te (kev or K, default keV)

logscale : bool

Search on a log scale for nearest temperature if set.

filename : str or hdulist

line or continuum file, already opened or filename.

Returns:
int

Index in HDU file with nearest temperature to te.

pyatomdb.spectrum.get_nei_line_emissivity(Z, z1, up, lo)

Return the line emissivity for a single line, separated out by the ion driving it

Parameters:
Z : int

nuclear charge

z1 : int

ion charge +1

up : int

the upper level

lo : int

the lower level

Returns:
emiss : dict

a dictionary containing an array, one for each ion_drv, with the emissivity in it. E.g. emiss[6] is a 51 element array, with the emissivity due to z1=6 as a fn of temperature Also emiss[‘Te’] is the tempearture in keV

pyatomdb.spectrum.get_response_ebins(rmf)

Get the energy bins from the rmf file

Parameters:
rmf : string or pyfits.hdu.hdulist.HDUList

The filename of the rmf or the opened rmf file

Returns:
specbins_in : array(float)

input energy bins used. nbins+1 length, with the last item being the final bin This is the array on which the input spectrum should be calculated

specbins_out : array(float)

output energy bins used. nbins+1 length, with the last item being the final bin This is the array on which the output spectrum will be returned. Often (but not always) the same as specbins_in

pyatomdb.spectrum.list_lines(specrange, lldat=False, index=False, linefile=False, units='angstroms', Te=False, teunit='K', minepsilon=1e-20)

Gets list of the lines in a given spectral range

Note that the output from this can be passed directly to print_lines

Parameters:
specrange : [float,float]

spectral range [min,max] to return lines on

lldat : see notes

line data

index : int

index in lldat, see notes

linefile : see notes

line data file, see notes

units : {‘A’ , ‘keV’}

units of specrange (default A)

Te : float

electron temperature (used if index not set to select appropriate data HDU from line file)

teunit : {‘K’,’keV’,’eV’}

units of Te

minepsilon : float

minimum epsilon for lines to be returned, in ph cm^3 s^-1

Returns:
linelist : dtype=([(‘Lambda’, ‘>f4’), (‘Lambda_Err’, ‘>f4’), (‘Epsilon’, ‘>f4’), (‘Epsilon_Err’, ‘>f4’), (‘Element’, ‘>i4’), (‘Ion’, ‘>i4’), (‘UpperLev’, ‘>i4’), (‘LowerLev’, ‘>i4’)])

A line list filtered by the various elements.

Notes

The actual line list can be defined in one of several ways:

specrange = [10,100]

  1. lldat as an actual list of lines:

    a = pyfits.open('apec_line.fits')
    llist = a[30].data
    l = list_lines(specrange, lldat=llist)
    
  2. lldat as a numpy array of lines:

    a = pyfits.open('apec_line.fits')
    llist = numpy.array(a[30].data)
    l = list_lines(specrange, lldat=llist)
    
  3. lldat is a BinTableHDU from pyfits:

    a = pyfits.open('apec_line.fits')
    llist = numpy.array(a[30])
    l = list_lines(specrange, lldat=llist)
    
  4. lldat is a HDUList from pyfits. In this case index must also be set:

    a = pyfits.open('apec_line.fits')
    index = 30
    l = list_lines(specrange, lldat=a, index=index)
    
  5. lldat NOT set, linefile contains apec_line.fits file location, index identifies the HDU:

    linefile = 'mydir/apec_v2.0.2_line.fits'
    index = 30
    l = list_lines(specrange, linefile=linefile, index=index)
    
  6. lldat NOT set & linefile NOT set, linefile is set to $ATOMDB/apec_line.fits. index identifies the HDU:

    index = 30
    l = list_lines(specrange, index=index)
    
pyatomdb.spectrum.list_nei_lines(specrange, Te, tau, Te_init=False, lldat=False, linefile=False, units='angstroms', teunit='K', minepsilon=1e-20, datacache=False)

Gets list of the lines in a given spectral range for a given NEI plasma

For speed purposes, this takes the nearest temperature tabulated in the linefile, and applies the exact ionization balance as calculated to this. This is not perfect, but should be good enough.

Note that the output from this can be passed directly to print_lines

Parameters:
specrange : [float,float]

spectral range [min,max] to return lines on

Te : float

electron temperature

tau : float

electron density * time (cm^-3 s)

Te_init : float

initial ionization balance temperature

lldat : see notes

line data

linefile : see notes

line data file, see notes

units : {‘A’ , ‘keV’}

units of specrange (default A)

teunit : {‘K’ , ‘keV’}

units of temperatures (default K)

minepsilon : float

minimum emissivity (ph cm^3 s^{-1}) for inclusion in linelist

Returns:
linelist : dtype=([(‘Lambda’, ‘>f4’), (‘Lambda_Err’, ‘>f4’), (‘Epsilon’, ‘>f4’), (‘Epsilon_Err’, ‘>f4’), (‘Element’, ‘>i4’), (‘Elem_drv’, ‘>i4’), (‘Ion’, ‘>i4’), (‘Ion_drv’, ‘>i4’), (‘UpperLev’, ‘>i4’), (‘LowerLev’, ‘>i4’)])

A line list filtered by the various elements.

Notes

The actual line list can be defined in one of several ways:

specrange = [10,100]

  1. lldat as an actual list of lines:

    a = pyfits.open('apec_nei_line.fits')
    llist = a[30].data
    l = list_nei_lines(specrange, lldat=llist)
    
  2. lldat as a numpy array of lines:

    a = pyfits.open('apec_nei_line.fits')
    llist = numpy.array(a[30].data)
    l = list_nei_lines(specrange, lldat=llist)
    
  3. lldat is a BinTableHDU from pyfits:

    a = pyfits.open('apec_nei_line.fits')
    llist = numpy.array(a[30])
    l = list_nei_lines(specrange, lldat=llist)
    
  4. lldat is a HDUList from pyfits. In this case index must also be set:

    a = pyfits.open('apec_nei_line.fits')
    index = 30
    l = list_nei_lines(specrange, lldat=a, index=index)
    
  5. lldat NOT set, linefile contains apec_line.fits file location, index identifies the HDU:

    linefile = 'mydir/apec_v3.0.2_nei_line.fits'
    index = 30
    l = list_nei_lines(specrange, linefile=linefile, index=index)
    
  6. lldat NOT set & linefile NOT set, linefile is set to $ATOMDB/apec_line.fits. index identifies the HDU:

    index = 30
    l = list_nei_lines(specrange, Te, tau)
    
pyatomdb.spectrum.list_nei_lines_accurate(specrange, Te, tau, Te_init=False, lldat=False, linefile=False, units='angstroms', teunit='K', minepsilon=1e-20, datacache=False)

Gets list of the lines in a given spectral range for a given NEI plasma

For speed purposes, this takes the nearest temperature tabulated in the linefile, and applies the exact ionization balance as calculated to this. This is not perfect, but should be good enough.

Note that the output from this can be passed directly to print_lines

Parameters:
specrange : [float,float]

spectral range [min,max] to return lines on

Te : float

electron temperature

tau : float

electron density * time (cm^-3 s)

Te_init : float

initial ionization balance temperature

lldat : see notes

line data

linefile : see notes

line data file, see notes

units : {‘A’ , ‘keV’}

units of specrange (default A)

teunit : {‘K’ , ‘keV’}

units of temperatures (default K)

minepsilon : float

minimum emissivity (ph cm^3 s^{-1}) for inclusion in linelist

Returns:
linelist : dtype=([(‘Lambda’, ‘>f4’), (‘Lambda_Err’, ‘>f4’), (‘Epsilon’, ‘>f4’), (‘Epsilon_Err’, ‘>f4’), (‘Element’, ‘>i4’), (‘Elem_drv’, ‘>i4’), (‘Ion’, ‘>i4’), (‘Ion_drv’, ‘>i4’), (‘UpperLev’, ‘>i4’), (‘LowerLev’, ‘>i4’)])

A line list filtered by the various elements.

Notes

The actual line list can be defined in one of several ways:

specrange = [10,100]

  1. lldat as an actual list of lines:

    a = pyfits.open('apec_nei_line.fits')
    llist = a[30].data
    l = list_nei_lines(specrange, lldat=llist)
    
  2. lldat as a numpy array of lines:

    a = pyfits.open('apec_nei_line.fits')
    llist = numpy.array(a[30].data)
    l = list_nei_lines(specrange, lldat=llist)
    
  3. lldat is a BinTableHDU from pyfits:

    a = pyfits.open('apec_nei_line.fits')
    llist = numpy.array(a[30])
    l = list_nei_lines(specrange, lldat=llist)
    
  4. lldat is a HDUList from pyfits. In this case index must also be set:

    a = pyfits.open('apec_nei_line.fits')
    index = 30
    l = list_nei_lines(specrange, lldat=a, index=index)
    
  5. lldat NOT set, linefile contains apec_line.fits file location, index identifies the HDU:

    linefile = 'mydir/apec_v3.0.2_nei_line.fits'
    index = 30
    l = list_nei_lines(specrange, linefile=linefile, index=index)
    
  6. lldat NOT set & linefile NOT set, linefile is set to $ATOMDB/apec_line.fits. index identifies the HDU:

    index = 30
    l = list_nei_lines(specrange, Te, tau)
    
pyatomdb.spectrum.make_ion_index_continuum(bins, element, index=False, cocofile='$ATOMDB/apec_coco.fits', binunits='keV', fluxunits='ph', no_coco=False, no_pseudo=False, ion=0, broadening=False, broadenunits='keV')

Creates the continuum for a given ion.

Parameters:
bins : array(float)

The bin edges for the spectrum to be calculated on, in units of keV or Angstroms. Must be monotonically increasing. Spectrum will return len(bins)-1 values.

element : int

Atomic number of element to make spectrum of (e.g. 6 for carbon)

binunits : {‘keV’ , ‘A’}

The energy units for bins. “keV” or “A”. Default keV.

fluxunits : {‘ph’ , ‘erg’}

Whether to return the emissivity in photons (‘ph’) or ergs (‘erg’). Defaults to photons

no_coco : bool

If true, do not include the compressed continuum

no_pseudo : bool

If true, do not include the pseudo continuum (weak lines)

ion : int

Ion to calculate, e.g. 4 for C IV. By default, 0 (whole element).

index : int

The index to generate the spectrum from. Note that the AtomDB files the emission starts in hdu number 2. So for the first block, you set index=2. Only required if cocofile is a filename or an HDULIST

cocofile : HDUList, HDU or str

The continuum file, either already open (HDULIST) or filename. alternatively, provide the HDU itself, and then do not need to define the index

broadening: float

Broaden the continuum by gaussians of this width (if False, no broadening is applied)

broadenunits: {‘keV’ , ‘A’}

Units for broadening (kev or A)

Returns:
array(float)

len(bins)-1 array of continuum emission, in units of photons cm^3 s^-1 bin^-1 (fluxunits = ‘ph’) or ergs cm^3 s^-1 bin^-1 (fluxunits = ‘erg’)

pyatomdb.spectrum.make_ion_spectrum(bins, index, Z, z1, linefile='$ATOMDB/apec_nei_line.fits', cocofile='$ATOMDB/apec_nei_comp.fits', binunits='keV', broadening=False, broadenunits='keV', abund=False, dummyfirst=False, nei=True, dolines=True, docont=True, dopseudo=True)

make_spectrum is the most generic “make me a spectrum” routine.

It returns the emissivity in counts cm^3 s^-1 bin^-1.

Parameters:
bins : array(float)

The bin edges for the spectrum to be calculated on, in units of keV or Angstroms. Must be monotonically increasing. Spectrum will return len(bins)-1 values.

index : int

The index to plot the spectrum from. note that the AtomDB files the emission starts in hdu number 2. So for the first block, you set index=2

Z : int

Element of spectrum (e.g. 6 for carbon)

z1 : int

Ion charge +1 for the spectrum (e.g. 3 for C III)

linefile : str or HDUList

The file containing all the line emission. Defaults to “$ATOMDB/apec_line.fits”. Can also pass in the opened file, i.e. “linefile = pyatomdb.pyfits.open(‘apec_nei_line.fits’)”

cocofile : str or HDUList

The file containing all the continuum emission. Defaults to “$ATOMDB/apec_coco.fits”. Can also pass in the opened file, i.e. “cocofile = pyatomdb.pyfits.open(‘apec_nei_comp.fits’)”

binunits : {‘keV’,’A’}

The energy units for bins. “keV” or “A”. Default keV.

broadening : float

Line broadening to be applied

broadenunits : {‘keV’,’A’}

Units of line broadening “keV” or “A”. Default keV.

elements : iterable of int

Elements to include, listed by atomic number. if not set, include all.

abund : iterable of float, length same as elements.

If set, and array of length (elements) with the abundances of each element relative to the Andres and Grevesse values. Otherwise, assumed to be 1.0 for all elements

dummyfirst : bool

If true, add a “0” to the beginning of the return array so it is of the same length as bins (can be useful for plotting results)

nei : bool

If set, return the spectrum from the driving ion being Z, rmJ. If not set, return the spectrum for the collisional ionization equilibrium BUT note that the continuum will be wrong, as it is provided for each element as a whole.

dolines : bool

Include lines in the spectrum

docont : bool

Include the continuum in the spectrum

dopseudo : bool

Include the pseudocontinuum in the spectrum.

Returns:
array of floats

Emissivity in counts cm^3 s^-1 bin^-1.

pyatomdb.spectrum.make_spectrum(bins, index, linefile='$ATOMDB/apec_line.fits', cocofile='$ATOMDB/apec_coco.fits', binunits='keV', broadening=False, broadenunits='keV', elements=False, abund=False, dummyfirst=False, dolines=True, docont=True, dopseudo=True)

make_spectrum is the most generic “make me a spectrum” routine.

It returns the emissivity in counts cm^3 s^-1 bin^-1.

Parameters:
bins : array(float)

The bin edges for the spectrum to be calculated on, in units of keV or Angstroms. Must be monotonically increasing. Spectrum will return len(bins)-1 values.

index : int

The index to plot the spectrum from. note that the AtomDB files the emission starts in hdu number 2. So for the first block, you set index=2

linefile : str

The file containing all the line emission. Defaults to “$ATOMDB/apec_line.fits”

cocofile : str

The file containing all the continuum emission. Defaults to “$ATOMDB/apec_coco.fits”

binunits : {‘keV’,’A’}

The energy units for bins. “keV” or “A”. Default keV.

broadening : float

Line broadening to be applied

broadenunits : {‘keV’,’A’}

Units of line broadening “keV” or “A”. Default keV.

elements : iterable of int

Elements to include, listed by atomic number. if not set, include all.

abund : iterable of float, length same as elements.

If set, and array of length (elements) with the abundances of each element relative to the Andres and Grevesse values. Otherwise, assumed to be 1.0 for all elements

dummyfirst : bool

If true, add a “0” to the beginning of the return array so it is of the same length as bins (can be useful for plotting results)

dolines : bool

Include lines in the spectrum

docont : bool

Include the continuum in the spectrum

dopseudo : bool

Include the pseudocontinuum in the spectrum.

Returns:
array of floats

Emissivity in counts cm^3 s^-1 bin^-1.

pyatomdb.spectrum.print_lines(llist, specunits='A', do_cfg=False)

Prints lines in a linelist to screen

This routine is very primitive as things stand. Plenty of room for refinement.

Parameters:
llist: dtype(linelist)

list of lines to print. Typically returned by list_lines.

specunits: {‘A’ , ‘keV’}

units to list the line positions by (A or keV, default A)

do_cfg: bool

Show full configuration information for each level

Returns:
Nothing, though prints data to standard out.