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.CXSession(linefile='$ATOMDB/apec_line.fits', cocofile='$ATOMDB/apec_coco.fits', elements=False, abundset='AG89', collisionunits='kev/amu')

Bases: pyatomdb.spectrum.Session

A Charge Exchange session using the same line and coco files, and/or responses

Attributes

linefile (string) The line emissivity data file
cocofile (string) The continuum emissivity data file
linedata: HDUList The line emissivity data
cocodata: HDUList The line emissivity data
elements (array_like, int) The atomic number of the elements to include. Defaults to all.
abundset (string) The elemental abundances to be used. Defaults to Anders and Grevesse 1989.
collisionunits (string) Whether the units are given in energy (kev/amu) or velocity (cm/s)
ready (bool) Set when line, continuum and spectral bin data has been read in, and a spectrum can be calculated.
default_abundset (string) The abundance set used in line and continuum files
abundset (string) The abundance set to be used in calculating the spectra.
response_set (bool) If a response (rmf & arf) have been loaded, set to true
spectra (dict of array_like) Holds the spectra at each temperature.
rmffile (string) Filename of RMF file
arffile (string) Filename of ARF file
rmf (HDUList) RMF data
arf (HDUList) ARF data
ionbal (dict of array like) ionization balance of each ion, normalized to 1 e.g. ionbal[6]=numpy.array([0.5,0.4,0.1,0,0,0,0]) for Carbon
class IonSpec(session, index, Z, z1)

An individual ion spectrum within a session, from a specifically tabulated temperature in a line/coco file.

Attributes

energy (float) The energy of this spectrum, in keV/amu
index (int) The index in the line file for this spectrum
calc_spectrum(session, dolines=True, docont=True, dopseudo=True)

Calculates the spectrum for each ion on a single energy

Parameters:

session : Session

The parent Session

dolines : bool

Include lines in the spectrum

docont : bool

Include continuum in the spectrum

dopseudo : bool

Include pseudocontinuum in the spectrum

Outputs

——-

none

Notes

Modifies:

dict : self.spectrum the spectrum of the ion

dict : self.spectrum_withresp the spectrum of the ion, folded through response

Then calls recalc() to update the spectra

recalc(session)

Recalculate the spectrum - just for changing abundances etc. Does not recalculate spectrum fully, just changes the multipliers. Does nothing if self.ready is False, should be run after calc_spectrum.

Parameters:

session : Session

The parent session

Returns:

none

Notes

Modifies:

self.spectrum : array_like (float)

self.spectrum_withresp : array_like (float)

set_index(E, Eunit='kev/amu', logscale=False)

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:

E : float

Energy in keV/amu

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

Units of E (kev/amu default)

logscale : bool

Search on a log scale for nearest temperature if set.

Returns:

none

Notes

modifies self.index : int Index in HDU file with nearest temperature to te.

class CXSession.Spec(session, index)

An individual spectrum within a session, from a specifically tabulated temperature in a line/coco file.

Attributes

temperature (float) The temperature of this spectrum, in keV
index (int) The index in the line file for this spectrum
calc_spectrum(session, dolines=True, docont=True, dopseudo=True)

Calculates the spectrum for each element on a single temperature

Parameters:

session : Session

The parent Session

dolines : bool

Include lines in the spectrum

docont : bool

Include continuum in the spectrum

dopseudo : bool

Include pseudocontinuum in the spectrum

Outputs

——-

none

Notes

Modifies:

dict : self.spectrum_by_Z the spectrum of each element

dict : self.spectrum_by_Z_withresp the spectrum of each element, folded through response

Then calls recalc() to update the spectra

recalc(session)

Recalculate the spectrum - just for changing abundances etc. Does not recalculate spectrum fully, just changes the multipliers. Does nothing if self.ready is False, should be run after calc_spectrum.

Parameters:

session : Session

The parent session

Returns:

none

Notes

Modifies:

self.spectrum : array_like (float)

self.spectrum_withresp : array_like (float)

set_index(T, teunit='K', logscale=False)

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

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

Units of te (kev or K, default keV)

logscale : bool

Search on a log scale for nearest temperature if set.

Returns:

none

Notes

modifies self.index : int Index in HDU file with nearest temperature to te.

CXSession.recalc()

Recalculate the spectrum - just for changing abundances etc. Does not recalculate spectrum fully, just changes the multipliers. Does nothing if self.ready is False, should be run after calc_spectrum.

Parameters:none
Returns:none

Notes

modifies self.spectrum

CXSession.return_spectra(collision, raw=False, nearest=False)

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

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:

collision : float

The energy (kev/amu) or velocity (cm/s) of the collision

raw : bool

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

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.

CXSession.set_abund(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)
CXSession.set_abundset(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.

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

Set the apec line and coco files

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.

elements : array_like(int)

The atomic numbers of the elements to include. Defaults to all (1-30)

abundset : string

The abundance set to use. Defaults to AG89. See atomdb.set_abundance

Returns:

None

Notes

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

CXSession.set_ionbal_temperature(te, teunit='keV')

Set the ionization balance to that of a given electron temperature

Parameters:

te : float

Electron Temperature

teunit : string

Units for the temperature. keV or A.

Notes

Modifies self.ionbal

CXSession.set_response(rmf, arf=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
Parameters:

rmf: string or HDUlist

The response matrix file

arf: string or HDUlist

The ancillary response file

Returns:

none

CXSession.set_specbins(specbins, specunits='A')

Set the energy or wavelength bin for the raw spectrum

Note that this is overridden if a response is loaded

Parameters:

ebins : array(float)

The edges of the spectral bins (for n bins, have n+1 edges)

specunits : {‘a’,’kev’}

The spectral bin units to use. Default is angstroms

Returns:

None

Notes

updates self.specbins, self.binunits, self.specbins_set

class pyatomdb.spectrum.Session(linefile='$ATOMDB/apec_line.fits', cocofile='$ATOMDB/apec_coco.fits', elements=False, abundset='AG89')

A session using the same line and coco files, and/or responses

Attributes

linefile (string) The line emissivity data file
cocofile (string) The continuum emissivity data file
linedata: HDUList The line emissivity data
cocodata: HDUList The line emissivity data
elements (array_like, int) The atomic number of the elements to include. Defaults to all.
abundset (string) The elemental abundances to be used. Defaults to Anders and Grevesse 1989.
ready (bool) Set when line, continuum and spectral bin data has been read in, and a spectrum can be calculated.
default_abundset (string) The abundance set used in line and continuum files
abundset (string) The abundance set to be used in calculating the spectra.
response_set (bool) If a response (rmf & arf) have been loaded, set to true
spectra (dict of array_like) Holds the spectra at each temperature.
rmffile (string) Filename of RMF file
arffile (string) Filename of ARF file
rmf (HDUList) RMF data
arf (HDUList) ARF data
class Spec(session, index)

An individual spectrum within a session, from a specifically tabulated temperature in a line/coco file.

Attributes

temperature (float) The temperature of this spectrum, in keV
index (int) The index in the line file for this spectrum
calc_spectrum(session, dolines=True, docont=True, dopseudo=True)

Calculates the spectrum for each element on a single temperature

Parameters:

session : Session

The parent Session

dolines : bool

Include lines in the spectrum

docont : bool

Include continuum in the spectrum

dopseudo : bool

Include pseudocontinuum in the spectrum

Outputs

——-

none

Notes

Modifies:

dict : self.spectrum_by_Z the spectrum of each element

dict : self.spectrum_by_Z_withresp the spectrum of each element, folded through response

Then calls recalc() to update the spectra

recalc(session)

Recalculate the spectrum - just for changing abundances etc. Does not recalculate spectrum fully, just changes the multipliers. Does nothing if self.ready is False, should be run after calc_spectrum.

Parameters:

session : Session

The parent session

Returns:

none

Notes

Modifies:

self.spectrum : array_like (float)

self.spectrum_withresp : array_like (float)

set_index(T, teunit='K', logscale=False)

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

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

Units of te (kev or K, default keV)

logscale : bool

Search on a log scale for nearest temperature if set.

Returns:

none

Notes

modifies self.index : int Index in HDU file with nearest temperature to te.

Session.recalc()

Recalculate the spectrum - just for changing abundances etc. Does not recalculate spectrum fully, just changes the multipliers. Does nothing if self.ready is False, should be run after calc_spectrum.

Parameters:none
Returns:none

Notes

modifies self.spectrum

Session.return_spectra(te, teunit='keV', raw=False, nearest=False, get_nearest_t=False)

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

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

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

Session.set_abund(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)
Session.set_abundset(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.

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

Set the apec line and coco files

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.

elements : array_like(int)

The atomic numbers of the elements to include. Defaults to all (1-30)

abundset : string

The abundance set to use. Defaults to AG89. See atomdb.set_abundance

Returns:

None

Notes

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

Session.set_response(rmf, arf=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
Parameters:

rmf: string or HDUlist

The response matrix file

arf: string or HDUlist

The ancillary response file

Returns:

none

Session.set_specbins(specbins, specunits='A')

Set the energy or wavelength bin for the raw spectrum

Note that this is overridden if a response is loaded

Parameters:

ebins : array(float)

The edges of the spectral bins (for n bins, have n+1 edges)

specunits : {‘a’,’kev’}

The spectral bin units to use. Default is angstroms

Returns:

None

Notes

updates self.specbins, self.binunits, self.specbins_set

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.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_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:

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

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.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 monotonicallyincreasing. Spectrum will return len(bins)-1 values.

index : int

The index to plot the spectrum from. note that the AtomDB filesthe emission starts in hdu number 2. So for the first block, youset 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 eachelement relative to the Andres and Grevesse values. Otherwise, assumed tobe 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 monotonicallyincreasing. Spectrum will return len(bins)-1 values.

index : int

The index to plot the spectrum from. note that the AtomDB filesthe emission starts in hdu number 2. So for the first block, youset 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 eachelement relative to the Andres and Grevesse values. Otherwise, assumed tobe 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')

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)

Returns:

Nothing, though prints data to standard out.