spectrum module

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

The modules is separated by the type of spectrum you wish to model. For now, the CIE, NEI and PShock models are implemented.

Roughly speaking, you load a [Modelname]Session, and you can then obtain:

linelist :
list of lines in a certain wavelength interval)
line_emissivity :
the emissivity of a specific line as a function of temperature, line intensity and more.
spectrum :
the emissivity in photons cm^3 bin^-1 s^-1. This requires a response be set first
class pyatomdb.spectrum.CIESession(linefile='$ATOMDB/apec_line.fits', cocofile='$ATOMDB/apec_coco.fits', elements=False, 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 : iterable of int

Elements to include, listed by atomic number. if not set, include 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 : iterable of int

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

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 (default True)

docont : bool

Calculate Continuum emission (default True)

dopseudo : bool

Calculate PseudoContinuum (weak line) emission (default True)

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.

adjust_line(Z=0, z1=0, z1_drv=0, upper=0, lower=0, quantity='Epsilon', method='Replace')

Change the emissivity or wavelength of a line. Integer parameters set to 0 mean “all”. Note this all happens in memory and does not edit the underlying files.

Parameters:
change : float or str

If float, set the new value to this. If string

Z : int

Element

z1 : int

Ion

z1_drv : int

Driving ion

upper : int

Upper level

lower : int

Lower level

quantity : string

Change “Epsilon” or “Lambda” - emissivity or wavelength - by change

method : string

“Replace”: replace existing value with change “Multiply” : multiply existing value with change “Divide” : divide existing value by change “Add” : add change to existing value “Subtract” : subtract change from existing

Returns:
None
return_line_emissivity(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

Perform linear interpolation on a logT/logEpsilon grid (default), or linear. 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(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.

nearest : bool

Return spectrum at nearest tabulated temperature, without interpolation

apply_binwidth : bool

Divide the line emissivity by the width of the bin they occupy to give emissivity per angstrom or per keV.

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(te, teunit='keV', nearest=False, get_nearest_t=False, log_interp=True, dolines=True, docont=True, dopseudo=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)

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.

log_interp : bool

Perform linear interpolation on a logT/logEpsilon grid, instead of linear.

dolines : bool

Calculate line emission (default True)

docont : bool

Calculate Continuum emission (default True)

dopseudo : bool

Calculate PseudoContinuum (weak line) emission (default True)

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 get_nearest_t is set, return the actual temperature this corresponds to. Units are same as teunit

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)
set_abundset(abundstring=None)

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_broadening(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_eebrems(do_eebrems)

Set whether to do the electron-electron bremsstrahlung in the spectrum

Parameters:
do_eebrems : bool

If True, include the electron-electron bremsstrahlung calculation If False, don’t (this is default, and matches XSPEC behaviour)

set_response(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.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 Non-equilibrium ionization (NEI) 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 : iterable of int

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

abundset : string

The abundance set to use. Default AG89.

Examples

Create a session instance:

>>> s = NEISession()

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, Ne *t = 1e11 cm^-3 s

>>> spec = s.return_spectrum(1.0, 1e11)

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 : NEISpectra

Object storing the actual spectral data

elements : iterable of int

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

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 (default True)

docont : bool

Calculate Continuum emission (default True)

dopseudo : bool

Calculate PseudoContinuum (weak line) emission (default True)

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.

adjust_line(Z=0, z1=0, z1_drv=0, upper=0, lower=0, quantity='Epsilon', method='Replace')

Change the emissivity or wavelength of a line. Integer parameters set to 0 mean “all”. Note this all happens in memory and does not edit the underlying files.

Parameters:
change : float or str

If float, set the new value to this. If string

Z : int

Element

z1 : int

Ion

z1_drv : int

Driving ion

upper : int

Upper level

lower : int

Lower level

quantity : string

Change “Epsilon” or “Lambda” - emissivity or wavelength - by change

method : string

“Replace”: replace existing value with change “Multiply” : multiply existing value with change “Divide” : divide existing value by change “Add” : add change to existing value “Subtract” : subtract change from existing

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

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 or array(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

Perform linear interpolation on a logT/logEpsilon grid (default), or linear.

init_pop : string or float

if ‘ionizing’ : all ionizing from neutral (so [1,0,0,0…]) if ‘recombining’: all recombining from ionized (so[…0,0,1]) if dict of arrays : the acutal fractional populations (so init_pop[6] is the array for carbon) if single float : the temperature (same units as Te)

freeze_ion_pop : bool

If True, skip the ion population calculation, use init_pop as the final pop instead.

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. If Te or Tau was supplied as a scalar, then that index is removed

return_linelist(Te, tau, specrange, init_pop='ionizing', specunit='A', teunit='keV', apply_aeff=False, by_ion_drv=False, nearest=False, apply_binwidth=False, log_interp=True, freeze_ion_pop=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.

by_ion_drv : bool

If true, keep lines which are the same but have different ion_drv separate. Otherwise, merge them.

nearest : bool

calculate spectrum at nearest temperature in linelist, no interpolation. Ionization fraction calculation still exact.

apply_binwidth :
log_interp : bool

Perform linear interpolation on a logT/logEpsilon grid (default), or linear.

freeze_ion_pop : bool

If True, skip the ion population calculation, use init_pop as the final pop instead.

Returns:
linelist : array(linelist)

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(Te, tau, init_pop='ionizing', teunit='keV', nearest=False, log_interp=True, freeze_ion_pop=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

tau : float

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

init_pop : string or float

if ‘ionizing’ : all ionizing from neutral (so [1,0,0,0…]) if ‘recombining’: all recombining from ionized (so[…0,0,1]) if dict of arrays : the acutal fractional populations (so init_pop[6] is the array for carbon) if single float : the temperature (same units as Te)

teunit : {‘keV’ , ‘K’}

Units of te (kev or K, default keV)

nearest : bool

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

log_interp : bool

Perform linear interpolation on a logT/logEpsilon grid (default), or linear.

freeze_ion_pop : bool

If True, skip the ion population calculation, use init_pop as the final pop instead.

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(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(abundstring=None)

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_broadening(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_eebrems(do_eebrems)

Set whether to do the electron-electron bremsstrahlung in the spectrum

Parameters:
do_eebrems : bool

If True, include the electron-electron bremsstrahlung calculation If False, don’t (this is default, and matches XSPEC behaviour)

set_response(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.PShockSession(linefile='$ATOMDB/apec_nei_line.fits', cocofile='$ATOMDB/apec_nei_comp.fits', elements=False, abundset='AG89')

Bases: pyatomdb.spectrum.NEISession

Load and generate a Parallel Shock (PShock) 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 : iterable of int

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

abundset : string

The abundance set to use. Default AG89.

Examples

Create a session instance:

>>> s=PShockSession()

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, Tau_u = 1e11, Tau_l = 0.0 (default)

>>> spec = s.return_spectrum(1.0, 1e11)

Same with Tau_l = 1e9

>>> spec = s.return_spectrum(1.0, 1e11, tau_l = 1e9)

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 : PShockSpectra

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 (default True)

docont : bool

Calculate Continuum emission (default True)

dopseudo : bool

Calculate PseudoContinuum (weak line) emission (default True)

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.

adjust_line(Z=0, z1=0, z1_drv=0, upper=0, lower=0, quantity='Epsilon', method='Replace')

Change the emissivity or wavelength of a line. Integer parameters set to 0 mean “all”. Note this all happens in memory and does not edit the underlying files.

Parameters:
change : float or str

If float, set the new value to this. If string

Z : int

Element

z1 : int

Ion

z1_drv : int

Driving ion

upper : int

Upper level

lower : int

Lower level

quantity : string

Change “Epsilon” or “Lambda” - emissivity or wavelength - by change

method : string

“Replace”: replace existing value with change “Multiply” : multiply existing value with change “Divide” : divide existing value by change “Add” : add change to existing value “Subtract” : subtract change from existing

Returns:
None
return_line_emissivity(Telist, tau_ulist, Z, z1, up, lo, tau_llist=0.0, 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

tau_ulist : 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

tau_llist : float

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

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

Perform linear interpolation on a logT/logEpsilon grid (default), or 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(Te, tau_u, specrange, tau_l=0.0, init_pop='ionizing', specunit='A', teunit='keV', apply_aeff=False, by_ion_drv=False, nearest=False, apply_binwidth=False)

Get the list of line emissivities vs wavelengths

Parameters:
Te : float

Temperature in keV or K

tau_u : float

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

specrange : [float, float]

Minimum and maximum values for interval in which to search

tau_l :

lower limit of ionization timescale, ne * t (cm^-3 s) (defaults to 0)

init_pop : string or float

if ‘ionizing’ : all ionizing from neutral (so [1,0,0,0…]) if ‘recombining’: all recombining from ionized (so[…0,0,1]) if dict of arrays : the acutal fractional populations (so init_pop[6] is the array for carbon) if single float : the temperature (same units as Te)

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.

by_ion_drv : bool

If true, keep lines which are the same but have different ion_drv separate. Otherwise, merge them.

nearest :
apply_binwidth :
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(Te, tau_u, tau_l=0.0, init_pop='ionizing', teunit='keV', nearest=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

tau_u : float

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

tau_l :

lower limit of ionization timescale, ne * t (cm^-3 s) (defaults to 0)

init_pop : string or float

if ‘ionizing’ : all ionizing from neutral (so [1,0,0,0…]) if ‘recombining’: all recombining from ionized (so[…0,0,1]) if dict of arrays : the acutal fractional populations (so init_pop[6] is the array for carbon) if single float : the temperature (same units as Te)

teunit : {‘keV’ , ‘K’}

Units of te (kev or K, default keV)

nearest : bool

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

log_interp : bool

Perform linear interpolation on a logT/logEpsilon grid (default), or linear.

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.

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)
set_abundset(abundstring=None)

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_broadening(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_eebrems(do_eebrems)

Set whether to do the electron-electron bremsstrahlung in the spectrum

Parameters:
do_eebrems : bool

If True, include the electron-electron bremsstrahlung calculation If False, don’t (this is default, and matches XSPEC behaviour)

set_response(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
pyatomdb.spectrum.calc_ee_brems_spec(ebins, Te, dens, teunit='keV')

Calculate the electron-electron bremstrahlung per-bin emissivity

Parameters:
ebins : array(float)

The spectral bin edges in energy order (keV)

Te : float

Electron temperature (default, keV)

dens : float

The electron density (cm^-3)

teunit : string

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

Returns:
eebrems : array(float)

electron-electron bremsstrahlung in ph cm^3 s^-1 bin^-1