ChiantiPy.tools package¶
Submodules¶
ChiantiPy.tools.archival module¶
Functions for reading pre-v8 CHIANTI files
-
ChiantiPy.tools.archival.
elvlcRead
(ions, filename=0, verbose=0, useTh=0)¶ read a chianti energy level file and returns {“lvl”:lvl,”conf”:conf,”term”:term,”spin”:spin,”l”:l,”spd”:spd,”j”:j ,”mult”:mult,”ecm”:ecm,”eryd”:eryd,”ecmth”:ecmth,”erydth”:erydth,”ref”:ref,”pretty”:pretty, ‘ionS’:ions} if a energy value for ecm or eryd is zero(=unknown), the theoretical values (ecmth and erydth) are inserted
-
ChiantiPy.tools.archival.
elvlcWrite
(info, outfile=None, addLvl=0)¶ Write Chianti data to .elvlc file.
- Parameters
info (dict) – Information about the Chianti data to write. Should contain the following keys: ionS, the Chianti style name of the ion such as c_4 conf, an integer denoting the configuration - not too essential term, a string showing the configuration spin, an integer of the spin of the state in LS coupling l, an integer of the angular momentum quantum number spd, an string for the alphabetic symbol of the angular momemtum, S, P, D, etc j, a floating point number, the total angular momentum ecm, the observed energy in inverse cm, if unknown, the value is 0. eryd, the observed energy in Rydbergs, if unknown, the value is 0. ecmth, the calculated energy from the scattering calculation, in inverse cm erydth, the calculated energy from the scattering calculation in Rydbergs ref, the references in the literature to the data in the input info
outfile (str) – Output filename. ionS+’.elvlc’ (in current directory) if None
addLvl (int) – Add a constant value to the index of all levels
Notes
For use with files created before elvlc format change in November 2012
See also
ChiantiPy.tools.io.elvlcWrite()
Write .elvlc file using the new format.
-
ChiantiPy.tools.archival.
wgfaRead
(ions, filename=None, elvlcname=-1, total=False, verbose=False)¶ Read CHIANTI data from a .wgfa file.
- Parameters
ions (str) – Ion, e.g. ‘c_5’ for C V
filename (str) – Custom filename, will override that specified by ions
elvlcname (str) – If specified, the lsj term labels are returned in the pretty1 and pretty2 keys of Wgfa
total (bool) – Return the level 2 avalue data in Wgfa
verbose (bool)
- Returns
Wgfa – Information read from the .wgfa file. The dictionary structure is {“lvl1”,”lvl2”,”wvl”,”gf”,”avalue”,”ref”,”ionS”,”filename”}
- Return type
dict
Notes
This is text-wise not different than the v8 version except that it uses the archival elvlcRead in ~ChiantiPy.tools.archival though this has now been commented out. Can this routine be removed? Should the elvlcRead routine be uncommented?
See also
ChiantiPy.tools.io.wgfaRead()
Read .wgfa file with the new format.
ChiantiPy.tools.constants module¶
A set of physical constants, in (mostly) cgs unit system. Most are from 11.
References
Notes
Many of these can be replaced by the ~astropy.constants module. Elemental symbols can be removed in favor of the periodictable module. Spectroscopic roman numerals can be removed in favor of roman module. The Gauss-Laguerre weights can be calculated by ~numpy.polynomial.laguerre.laggauss.
ChiantiPy.tools.data module¶
Module for collecting various top-level Chianti data
Descriptions for keywordArgs:
temperature : temperature (in K)
eDensity : electron density (in cm−3)
hDensity : hydrogen density (in cm−3)
pDensity : proton density (in cm−3)
radTemperature : radiation temperature of central source (in K)
rStar : distance of the plasma from the source (in units of the source’s radius)
distance : distance from the central source
- Parameters
XUVTOP (str) – the root of the CHIANTI database
Defaults (dict) – default values used by ChiantiPy, can also be set by an optional HOME/.chianti/chiantirc file
Ip (np.ndarray) – the ionization potentials for all ionization stages up to Zn, in eV
MasterList (list) – the CHIANTI style names of all ions in the CHIANTI database
IoneqAll (dict) – a dict containing the ionization equilibrium values for the default ionization file
ChiantiVersion (str) – the version of the CHIANTI database
AbundanceDefault (dict) – the elemental abundances in the default abundance file
AbundanceList (list) – the names of all abundance files included in the CHIANTI database
GrndLevels (list) – the number of levels that should be considered in an ionization calculation
ChiantiPy.tools.filters module¶
Line profile filters for creating synthetic spectra.
-
ChiantiPy.tools.filters.
boxcar
(wvl, wvl0, factor=None)¶ Box-car filter
- Parameters
wvl (~numpy.ndarray) – Wavelength array
wvl0 (~numpy.float64) – Wavelength filter should be centered on.
factor (~numpy.float64) – Full width of the box-car filter
-
ChiantiPy.tools.filters.
gaussian
(wvl, wvl0, factor=1.0)¶ A gaussian filter
- Parameters
wvl (~numpy.ndarray) – Wavelength array
wvl0 (~numpy.float64) – Wavelength filter should be centered on.
factor (~numpy.float64) – Gaussian width
integrated value is unity
-
ChiantiPy.tools.filters.
gaussianR
(wvl, wvl0, factor=1000.0)¶ A gaussian filter where the gaussian width is given by wvl0/factor.
- Parameters
wvl (~numpy.ndarray) – Wavelength array
wvl0 (~numpy.float64) – Wavelength filter should be centered on.
factor (~numpy.float64) – Resolving power
-
ChiantiPy.tools.filters.
lorentz
(wvl, wvl0, factor=1.0)¶ Lorentz profile filter with the exception that all factors are in wavelength units rather than frequency as the lorentz profile is usually defined.
- Parameters
wvl (~numpy.ndarray) – Wavelength array
wvl0 (~numpy.float64) – Wavelength filter should be centered on.
factor (~numpy.float64) – Value of the so-called constant gamma
integrated value is unity
the FWHM is 2*gamma
.. math:: – L = frac{1}{pi gamma} frac{ gamma^2}{(lambda - lambda_0)^2 + gamma^2}
-
ChiantiPy.tools.filters.
moffat
(wvl, wvl0, factor=2.5)¶ Moffat profile with parameters suited to Chandra Letg spectra
- Parameters
wvl (~numpy.ndarray) – Wavelength array
wvl0 (~numpy.float64) – Wavelength the filter is centered on.
factor (~numpy.float64) – Resolving power (TODO: correct description)
-
ChiantiPy.tools.filters.
voigt
(wvl, wvl0, factor=(0.5, 1.0))¶ pseudo-Voigt filter the sum of a Gaussian and a Lorentzian
- Parameters
wvl (~numpy.ndarray) – Wavelength array
wvl0 (~numpy.float64) – Wavelength the filter is centered on.
factor (array-type) – contains the following 2 parameters
A (~numpy.float64) – relative size of gaussian and lorentz components must be between 0. and 1. but this is not currently checked
sigma (~numpy.float64) – the gaussian width of the gaussian profile (the standard deviation) also creates the lorentz component with the same fwhm
ChiantiPy.tools.io module¶
Reading and writing functions
TODO: see if klgfbRead is needed or not
-
ChiantiPy.tools.io.
abundanceRead
(abundancename=None, verbose=False)¶ Read abundance file abundancename and return the abundance values relative to hydrogen
- Keyword Arguments
abundancename (str) – the name of the abundance file in the $XUVTOP/abundance directory to be read the default is an empty string and then the ‘default’ abundance values are read
verbose (bool) – if true, prints out some info
-
ChiantiPy.tools.io.
autoRead
(ions, filename=None, total=True, verbose=False)¶ Read CHIANTI autoionization rates from a .auto file.
- Parameters
ions (str) – Ion, e.g. ‘c_5’ for C V
filename (str) – Custom filename, will override that specified by ions
elvlcname (str) – If specified, the lsj term labels are returned in the ‘pretty1’ and ‘pretty2’ keys of ‘Wgfa’ dict
total (bool) – Return the summed level 2 autoionization rates in ‘Auto’
verbose (bool)
- Returns
Auto – Information read from the .wgfa file. The dictionary structure is {“lvl1”, “lvl2”, “avalue”, “pretty1”, “pretty2”, “ref”,”ionS”, “filename”}
- Return type
dict
-
ChiantiPy.tools.io.
autoWrite
(info, filename=None, minBranch=None)¶ Write data to a CHIANTI .wgfa style file
- Parameters
info (dict) – Should contain the following: ionS, the Chianti style name of the ion such as c_4 for C IV, lvl1, the lower level, the ground level is 1, lvl2, the upper level, wvl, the wavelength (in Angstroms), avalue, the autoionization rate, pretty1, descriptive text of the lower level (optional), pretty2, descriptive text of the upper level (optiona), ref, reference text, a list of strings
outfile (str)
minBranch (~numpy.float64) – The transition must have a branching ratio greater than the specified to be written to the file
-
ChiantiPy.tools.io.
cireclvlRead
(ions, filename=None, filetype='cilvl')¶ Read Chianti cilvl, reclvl, or rrlvl files and return data
- Parameters
ions (str) – Ion, e.g. ‘c_5’ for C V
filename (str, optional) – Custom filename, will override that specified by ions
filetype (str) – {‘cilvl’, ‘reclvl’, ‘rrlvl’} Type of the file to read
-
ChiantiPy.tools.io.
convertName
(name)¶ Convert ion name string to Z and Ion and other interesting info
- Parameters
name (str) – a generic name of an ion in the CHIANTI database, such as fe_14 for Fe XIV
Todo
Put in separate module to avoid multiple copies
Notes
A duplicate of the routine in ChiantiPy.tools.util but needed by masterList Info
-
ChiantiPy.tools.io.
defaultsRead
(verbose=False)¶ Read in configuration from .chiantirc file or set defaults if one is not found.
-
ChiantiPy.tools.io.
demRead
(demName='')¶ Read emission measure file emName and return the temperatures, densities and emission measures
- Keyword Arguments
demName (str) – the name of the differential emission measure file to read in the $XUVTOP/dem directory
- Returns
DEM – keywords are temperature, density, dem, em, dt, filename
- Return type
dict
-
ChiantiPy.tools.io.
diRead
(ions, filename=None)¶ Read chianti direct ionization .params files and return data.
- Parameters
ions (str) – Ion, e.g. ‘c_5’ for C V
filename (str, optional) – Custom filename, will override that specified by ions
-
ChiantiPy.tools.io.
drRead
(ions, filename=None)¶ Read CHIANTI dielectronic recombination .drparams files if filename is set, then reads that file
- Parameters
ions (str) – Ion, e.g. ‘c_5’ for C V
filename (str, optional) – Custom filename, will override that specified by ions
-
ChiantiPy.tools.io.
eaRead
(ions, filename=None)¶ Read a CHIANTI excitation-autoionization file and calculate the EA ionization rate data derived from splupsRead.
- Parameters
ions (str) – Ion, e.g. ‘c_5’ for C V
filename (str, optional) – Custom filename, will override that specified by ions
-
ChiantiPy.tools.io.
elvlcRead
(ions, filename=None, getExtended=False, getLatex=False, verbose=False, useTh=True)¶ Reads the new format elvlc files.
- Parameters
ions (str) – Ion, e.g. ‘c_5’ for C V
filename (str, optional) – Custom filename, will override that specified by ions
getExtended (bool) – get energy levels in columns beyond the theoretical energy
getLatex (bool) – get energy level descriptions in latex markup
verbose (bool) – print out the elvlc information
useTh (bool) – If True, the theoretical values (ecmth and erydth) are inserted when an energy value for ecm or eryd is zero(=unknown)
-
ChiantiPy.tools.io.
elvlcWrite
(info, outfile=None, round=0, addLvl=0, includeRyd=False, includeEv=False, verbose=False)¶ Write Chianti data to .elvlc file.
- Parameters
info (dict) – Information about the Chianti data to write. Should contain the following keys: ionS, the Chianti style name of the ion such as c_4 term, a string showing the configuration spin, an integer of the spin of the state in LS coupling l, an integer of the angular momentum quantum number spd, an string for the alphabetic symbol of the angular momemtum, S, P, D, etc j, a floating point number, the total angular momentum ecm, the observed energy in inverse cm, if unknown, the value is 0. eryd, the observed energy in Rydbergs, if unknown, the value is 0. ecmth, the calculated energy from the scattering calculation, in inverse cm erydth, the calculated energy from the scattering calculation in Rydbergs ref, the references in the literature to the data in the input info
outfile (str) – Output filename. ionS + ‘.elvlc’ (in current directory) if None
round (int) – input to ‘np.round’ to round input values to maintain the correct number of significant figures
addLvl (int) – Add a constant value to the index of all levels
includeRyd (bool) – If True, write the Rydberg energies in the extended area, delimited by a comma
includeEv (bool) – If True, write the energies in eV in the extended area, delimited by a comma
Notes
For use with files created after elvlc format change in November 2012
See also
ChiantiPy.tools.archival.elvlcWrite()
Write .elvlc file using the old format.
-
ChiantiPy.tools.io.
emRead
(filename=None, directory=None, verbose=False)¶ Read emission measure file emName and return the temperatures, densities and emission measures
- Keyword Arguments
filename (str) – the name of the emission measure file to read in the $XUVTOP/em directory
- Returns
EM – keywords are temperature, density, em, filename
- Return type
dict
-
ChiantiPy.tools.io.
fblfRead
(ions, filename=None, verbose=False)¶ Read the .fbparams file to return the parameters to calculate the fb factor of Mao, Kaastra, Badnell Astron. Astrophys. 599, A10 (2017)
- Parameters
ions (str) – Ion, e.g. ‘c_5’ for C V
filename (str, optional) – Custom filename, will override that specified by ions
-
ChiantiPy.tools.io.
fblvlRead
(ions, filename=None, verbose=False)¶ Read a Chianti energy level file for calculating the free-bound continuum
- Parameters
ions (str) – Ion, e.g. ‘c_5’ for C V
filename (str, optional) – Custom filename, will override that specified by ions
-
ChiantiPy.tools.io.
gffRead
()¶ Read the free-free gaunt factors of 1.
References
Notes
This function reads the file and reverses the values of g2 and u
-
ChiantiPy.tools.io.
grndLevelsRead
()¶ to read the grndLevels.dat file give the number of ground levels to sum over in populate and drPopulate
-
ChiantiPy.tools.io.
ioneqRead
(ioneqName='', minIoneq=1e-20, verbose=False)¶ Reads an ioneq file ionization equilibrium values less then minIoneq are returns as zeros
- Keyword Arguments
ioneqName (str) – reads the file in the $XUVTOP/ioneq directory, if a blank, the default is read
minIoneq (float) – sets values to zero if less the minIoneq
verbose (bool) – if true, prints into to the terminal
- Returns
{‘ioneqname’,’ioneqAll’,’ioneqTemperature’,’ioneqRef’} – Ionization equilibrium values and the reference to the literature
- Return type
dict
-
ChiantiPy.tools.io.
ipRead
(verbose=False)¶ Reads the ionization potential file
- Returns
ip – Ionization potential (in eV)
- Return type
array-like
-
ChiantiPy.tools.io.
klgbfnRead
(filename)¶ Read CHIANTI files containing the free-bound gaunt factors for n=1-6 from 14. This reads the corrected KL files in CHIANTI version 10+
- Parameters
filename (str) – the filename for the KL data in the continuum directory, such as klgfb_1.dat, n=1,6
- Returns
[{‘pe’, ‘klgfb’}] – Photon energy and the bound-free gaunt factors
- Return type
list
References
-
ChiantiPy.tools.io.
maoParsRead
(filename=None)¶ to read the mao et al par2.dat file for calculating the ratio of rrloss to rrrecomb The electron energy-loss rate due to radiative recombination. Mao J., Kaastra J., Badnell N.R. <Astron. Astrophys. 599, A10 (2017)> =2017A&A…599A..10M 1- 2 I2 — s Isoelectronic sequence number 4- 5 I2 — z Atomic number 7- 16 E10.4 — a0 Primary fitting parameter 18- 27 E10.4 — b0 Primary fitting parameter 29- 38 E10.4 — c0 Addiational fitting parameter 40- 49 E10.4 — a1 Primary fitting parameter 51- 60 E10.4 — b1 Primary fitting parameter 62- 71 E10.4 — a2 Addiational fitting parameter 73- 82 E10.4 — b2 Addiational fitting parameter 84- 86 F3.1 — mdp Maximum deviation in percent
- Keyword Arguments
filename (str) –
- Returns
data
- Return type
dict
-
ChiantiPy.tools.io.
masterListInfo
(force=False, verbose=False)¶ Get information about ions in the CHIANTI masterlist.
- Keyword Arguments
force (bool) – if true, recreates the masterListInfo file
verbose (bool) – if true, prints into to the terminal
- Returns
masterListInfo – {‘wmin’, ‘wmax’, ‘tmin’, ‘tmax’} Minimum and maximum wavelengths in the wgfa file. Minimum and maximum temperatures for which the ionization balance is nonzero.
- Return type
dict
Notes
This function speeds up multi-ion spectral calculations. The information is stored in a pickled file ‘masterlist_ions.pkl’ If the file is not found, one will be created.
-
ChiantiPy.tools.io.
masterListRead
()¶ Read a CHIANTI masterlist file.
- Returns
masterlist – All ions in Chianti database
- Return type
list
-
ChiantiPy.tools.io.
rrLossRead
()¶ to read the Mao 2017 rr loss parameters 12
References
- 12
Mao J., Kaastra J., Badnell N.R., 2017 Astron. Astrophys. 599, A10
-
ChiantiPy.tools.io.
rrRead
(ions, filename=None)¶ Read CHIANTI radiative recombination .rrparams files
- Parameters
ions (str) – Ion, e.g. ‘c_5’ for C V
- Keyword Arguments
filename (str, optional) – Custom filename, will override that specified by ions
- Returns
{‘rrtype’,’params’,’ref’}
- Return type
dict
-
ChiantiPy.tools.io.
scupsRead
(ions, filename=None, verbose=False)¶ Read the new format v8 scups file containing the scaled temperature and upsilons from 8.
- Parameters
ions (str) – Ion, e.g. ‘c_5’ for C V
- Keyword Arguments
filename (str, optional) – Custom filename, will override that specified by ions
verbose (bool) – if True, prints info to the terminal
-
ChiantiPy.tools.io.
splomRead
(ions, ea=False, filename=None)¶ Read chianti .splom files
- Parameters
ions (str) – Ion, e.g. ‘c_5’ for C V
- Keyword Arguments
ea (bool) – if true, reads the .easplom file
filename (str, optional) – Custom filename, will override that specified by ions
- Returns
{‘lvl1’, ‘lvl2’, ‘ttype’, ‘gf’, ‘deryd’, ‘c’, ‘splom’, ‘ref’}
- Return type
dict
Notes
Still needed for ionization cross sections
-
ChiantiPy.tools.io.
splupsRead
(ions, filename=None, filetype='splups')¶ Read a CHIANTI .splups file
- Parameters
ions (str) – Ion, e.g. ‘c_5’ for C V
- Keyword Arguments
filename (str, optional) – Custom filename, will override that specified by ions
filetype (str, optional) – {psplups,`cisplups`,`splups`} Type of file to read
- Returns
{‘lvl1’, ‘lvl2’, ‘ttype’, ‘gf’, ‘de’, ‘cups’, ‘bsplups’, ‘ref’}
- Return type
dict
-
ChiantiPy.tools.io.
twophotonHRead
()¶ Read the two-photon Einstein A values and distribution function for the H sequence.
- Returns
{‘y0’, ‘z0’, ‘avalue’, ‘asum’, ‘psi0’}
- Return type
dict
-
ChiantiPy.tools.io.
twophotonHeRead
()¶ Read the two-photon Einstein A values and distribution function for the He sequence.
- Returns
{‘y0’, ‘avalue’, ‘psi0’}
- Return type
dict
-
ChiantiPy.tools.io.
vernerRead
()¶ Reads the photoionization cross-section data from 6.
- Returns
{‘pqn’,’l’,’eth’,’e0’,’sig0’,’ya’,’p’, yw’} – pqn is the principal quantum number, l is the subshell orbital quantum number, eth (in eV) is the subshell ionization threshold energy; sig0, ya, p, and yw are all fit parameters used in calculating the total photoionization cross-section.
- Return type
dict
References
-
ChiantiPy.tools.io.
versionRead
()¶ Read the version number of the CHIANTI database
-
ChiantiPy.tools.io.
wgfaRead
(ions, filename=None, elvlcname=0, total=False, getLatex=False, verbose=False)¶ Read CHIANTI data from a .wgfa file.
- Parameters
ions (str) – Ion, e.g. ‘c_5’ for C V
- Keyword Arguments
filename (str) – Custom filename, will override that specified by ions
elvlcname (str) – If specified, the lsj term labels are returned in the ‘pretty1’ and ‘pretty2’ keys of ‘Wgfa’ dict
total (bool) – Return the summed level 2 avalue data in ‘Wgfa’
getLatex (bool) – retrieve the level descriptions in latex markup returned in the ‘latex1’ and ‘latex2’ keys of the Wgfa dict
verbose (bool) –
- Returns
Wgfa – Information read from the .wgfa file. The dictionary structure is {“lvl1”,”lvl2”,”wvl”,”gf”,”avalue”,”ref”,”ionS”,”filename”}
- Return type
dict
See also
ChiantiPy.tools.archival.wgfaRead()
Read .wgfa file with the old format.
-
ChiantiPy.tools.io.
wgfaWrite
(info, filename=None, minBranch=1e-05, rightDigits=4, maxLvl1=None, comment=None)¶ Write data to a CHIANTI .wgfa file
- Parameters
info (dict) – Should contain the following keys: ionS, the Chianti style name of the ion such as c_4 for C IV, lvl1, the lower level, the ground level is 1, lvl2, the upper level, wvl, the wavelength (in Angstroms), gf,the weighted oscillator strength, avalue, the A value, pretty1, descriptive text of the lower level (optional), pretty2, descriptive text of the upper level (optiona), ref, reference text, a list of strings
filename (str)
minBranch (~numpy.float64) –
- The transition must have a branching ratio greater than the specified minBranch
to be written to the file. default = 1.e-5
rightDigits (int) – the number of digits to the right of the decimal point in the wavelength default = 5
maxLvl1 (int) – the largest level to be written. default is None
comment (str) – add a comment to the reference section. default = None
-
ChiantiPy.tools.io.
zion2name
(z, ion, dielectronic=False)¶ Convert Z and ion to generic name, e.g. 26, 13 -> fe_13
- Parameters
z (int) – the nuclear charge, for example 26 for Fe XIV
ion (int) – the ion stage, for example, 14 for Fe XIV
- Keyword Arguments
dielectronic (bool) – if True, created the name of a dielectronic ion, with a ‘d’ at the end
Todo
See if dielectronic is still appropriate
Put in separate module to avoid multiple copies
Notes
A duplicate of the routine in ChiantiPy.tools.util but needed by masterList Info
ChiantiPy.tools.mputil module¶
Functions needed for standard Python multiprocessing module mspectrum
-
ChiantiPy.tools.mputil.
doFbLossQ
(inQ, outQ)¶ Multiprocessing helper for ChiantiPy.core.continuum.freeBound
- Parameters
inQ (~multiprocessing.Queue) – Ion free-bound emission jobs queued up by multiprocessing module
outQ (~multiprocessing.Queue) – Finished free-bound emission jobs
-
ChiantiPy.tools.mputil.
doFbQ
(inQ, outQ)¶ Multiprocessing helper for ChiantiPy.core.continuum.freeBound
- Parameters
inQ (~multiprocessing.Queue) – Ion free-bound emission jobs queued up by multiprocessing module
outQ (~multiprocessing.Queue) – Finished free-bound emission jobs
-
ChiantiPy.tools.mputil.
doFfLossQ
(inQ, outQ)¶ Multiprocessing helper for ChiantiPy.core.continuum.freeFreeLoss
- Parameters
inQ (~multiprocessing.Queue) – Ion free-free emission jobs queued up by multiprocessing module
outQ (~multiprocessing.Queue) – Finished free-free emission jobs
-
ChiantiPy.tools.mputil.
doFfQ
(inQ, outQ)¶ Multiprocessing helper for ChiantiPy.core.continuum.freeFree
- Parameters
inQ (~multiprocessing.Queue) – Ion free-free emission jobs queued up by multiprocessing module
outQ (~multiprocessing.Queue) – Finished free-free emission jobs
-
ChiantiPy.tools.mputil.
doIonLossQ
(inQueue, outQueue)¶ Multiprocessing helper for ChiantiPy.core.ion and ChiantiPy.core.ion.twoPhoton
- Parameters
inQueue (~multiprocessing.Queue) – Jobs queued up by multiprocessing module
outQueue (~multiprocessing.Queue) – Finished jobs
-
ChiantiPy.tools.mputil.
doIonQ
(inQueue, outQueue)¶ Multiprocessing helper for ChiantiPy.core.ion and ChiantiPy.core.ion.twoPhoton
- Parameters
inQueue (~multiprocessing.Queue) – Jobs queued up by multiprocessing module
outQueue (~multiprocessing.Queue) – Finished jobs
ChiantiPy.tools.sources module¶
Blackbody temperature calculations
-
class
ChiantiPy.tools.sources.
blackStar
(temperature, radius)¶ Bases:
object
Calculate blackbody radiation
- Parameters
temperature (~numpy.ndarray) – Temperature in Kelvin
radius (~numpy.ndarray) – Stellar radius in cm
- Variables
Temperature (~numpy.ndarray) – Temperature in Kelvin
Radius (~numpy.ndarray) – Stellar radius in cm
Incident (~numpy.ndarray) – Blackbody photon distribution
-
incident
(distance, energy)¶ Calculate photon distribution times the visible cross-sectional area.
- Parameters
distance (~numpy.ndarray) – Distance to the stellar object in cm
energy (~numpy.ndarray) – Energy range in erg
Notes
This function returns the photon distribution instead of the distribution times the cross-sectional area. Is this correct? Why is the incident photon distribution calculated at all?
-
ChiantiPy.tools.sources.
blackbody
(temperature, variable, hnu=1)¶ Calculate the blackbody photon distribution as a function of energy (hnu = 1) or as a function of wavelength (hnu = 0) in units of photonscm−2s−1str−1erg−1
- Parameters
temperature (~numpy.float64) – Temperature at which to calculate the blackbody photon distribution
variable (~numpy.ndarray) – Either energy (in erg) or wavelength (in angstrom)
hnu (int) – If 1, calculate distribution as a function of energy. Otherwise, calculate it as a function of wavelength
- Returns
{‘photons’, ‘temperature’, ‘energy’} or {‘photons’, ‘temperature’, ‘wvl’}
- Return type
dict
ChiantiPy.tools.util module¶
Utility functions
Notes
Some of these functions can be replaced by roman numeral and periodic table lookup libraries. some functions using os.walk can be replaced by os.path
-
ChiantiPy.tools.util.
between
(array, limits)¶ Find the indices of array corresponding to values in the range given by limits
- Parameters
array (list, ndarray) – contains a list of elements such a wavelengths
limits (list, ndarray) – a 2 element array specifying the lower and upper limits of the array to be included
-
ChiantiPy.tools.util.
convertName
(name)¶ Convert ion name string (e.g. ‘fe_13’) to atomic number and ion number.
- Parameters
name (str) – the CHIANTI style name for an ion, such as fe_14
- Returns
{‘Z’, ‘Ion’, ‘Dielectronic’, ‘Element’, ‘higher’, ‘lower’} – higher and lower are the Chianti-style names for the higher and lower ionization stages, respectively.
- Return type
dict
-
ChiantiPy.tools.util.
descale_bt
(bte, btomega, f, ev1)¶ Apply excitation descaling of 3 to energy and collision strength
- Parameters
bte (array-like) – Scaled energy
btomega (array-like) – Scaled collision strength
f (array-like)
ev1 (array-like)
- Returns
[energy,omega] – Descaled energy and collision strength
- Return type
list
Notes
Need more details here. Not clear which equations are being used.
See also
scale_bt()
Apply scaling to excitation energy and cross-section
References
-
ChiantiPy.tools.util.
descale_bti
(bte, btx, f, ev1)¶ Apply ionization descaling of 9 to energy and cross-sections of 7.
- Parameters
bte (array-like) – Scaled energy
btx (array-like) – Scaled cross-section
f (float) – Scaling parameter
ev1 (float) – ionization potential - units determine the output energy units
- Returns
[energy,cross] – Descaled energy and cross-section
- Return type
list
Notes
This is the scaling used and discussed in the Dere (2007) calculation 7 of cross sections. It was derived from similar scalings provided by reference [2]
See also
scale_bti()
Descale ionization energy and cross-section
References
-
ChiantiPy.tools.util.
descale_bti_rate
(btTemperature, btRate, ip, f=2.0)¶ Apply ionization descaling of 7, a Burgess-Tully type scaling to bt scaled ionization rates and temperatures. The result is to return a temperature array and a ionization rate array.
- Parameters
btTemperature (array-like) – the bt scaled temperature
btRate (array-like) – the bt scaled ionization rate
ip (float) – the ionization potential in eV.
f (float (optional)) – the scaling parameter, 1.7 generally works well
-
ChiantiPy.tools.util.
dilute
(radius)¶ Calculate the dilution factor.
- Parameters
radius (array-like) – Distance from the center of a star in units of the stellar radius.
Notes
If radius is less than 1.0, returns 0.
-
ChiantiPy.tools.util.
el2z
(els)¶ Convert elemental symbol to atomic number
- Parameters
els (str) – the abreviated element name
- Returns
z – the atomic number or nuclear charge of the element
- Return type
int
-
ChiantiPy.tools.util.
findFileTypes
(wpath, type='*.txt', verbose=False)¶ to find all the files in wpath and below with file names matching fname
-
ChiantiPy.tools.util.
ion2filename
(ions)¶ Convert ion name string to generic directory-file name. convertName has probably made this redundant
- Parameters
ions (str) – the CHIANTI style name for an ion, such as fe_14
- Returns
fname – the full file name of the ion directory in the CHIANTI database assumes a top directory from the environmental variable XUVTOP
- Return type
str
Todo
this duplicates what ‘convertName’ does
-
ChiantiPy.tools.util.
listFiles
(dir)¶ Walks the path and subdirectories to return a list of files.
- Parameters
dir (str) – the top directory to search subdirectories are also searched
- Returns
listname – a list of files in dir and subdirectories
- Return type
list
Notes
This can be replaced by functions in os.path, as if 3.4, pathlib is probably better. It is not clear that this function is used anywhere in ChiantiPy
-
ChiantiPy.tools.util.
listRootNames
(dir)¶ Walks the path and subdirectories to return a list of file root names.
Notes
This can be replaced by functions in os.path, as if 3.4, pathlib is probably better. Only seems to be used by
-
ChiantiPy.tools.util.
qrp
(z, u)¶ Calculate Q′R(Z,u), where u=ϵ/I is the impact electron energy in threshold units, from Eq. 2.12 of 4.
- Parameters
z (int) – Atomic number
u (array-like) – Impact electron energy in threshold units.
- Returns
q – 1s ionization cross-section, Q′R(Z,u)
- Return type
array-like
Notes
Used for calculations of direct ionization cross-sections of the H and He sequences in ChiantiPy.tools.io.twophotonHRead and ChiantiPy.tools.io.twophotonHeRead, respectively.
References
-
ChiantiPy.tools.util.
scale_bt
(evin, omega, f, ev1)¶ Apply excitation scaling of 8 to energy and collision strength.
- Parameters
evin (array-like)
omega (array-like)
f (array-like)
ev1 (array-like)
- Returns
[bte,btomega] – Scaled energy and collision strength
- Return type
list
Notes
Need more details here. Not clear which equations are being used.
See also
descale_bt()
Descale excitation energy and cross-section
-
ChiantiPy.tools.util.
scale_bt_rate
(inDict, ip, f=2.0)¶ Apply ionization scaling of 7, a Burgess-Tully type scaling to ionization rates and temperatures. The result of the scaling is to return a scaled temperature between 0 and 1 and a slowly varying scaled rate as a function of scaled temperature. In addition, the scaled rates vary slowly along an iso-electronic sequence.
- Parameters
inDict (dict) – the input dictionary should have the following key pairs: temperature, array-like and rate, array-like
ip (float) – the ionization potential in eV.
f (float (optional)) – the scaling parameter, 1.7 generally works well
Notes
btTemperature and btRate keys are added to inDict
-
ChiantiPy.tools.util.
scale_bti
(evin, crossin, f, ev1)¶ Apply ionization scaling of 7,[8]_, to energy and cross-section.
- Parameters
evin (float) – Energy - same units as ev1
crossin (array-like) – Cross-section
f (float - the scale factor)
ev1 (float) – the ionization potential units - the same as evin
- Returns
[bte,btx] – Scaled energy and cross-section
- Return type
list
Notes
This is the scaling used and discussed in the Dere (2007) calculation [1] of cross sections. It was derived from similar scalings derived in reference [2]
See also
descale_bti()
Descale ionization energy and cross-section
References
-
ChiantiPy.tools.util.
scale_classical
(inDict, ip)¶ Apply the classical scaling to the input data
- Parameters
inDict (dictionary) –
- the input dictionary should have the following key pairs
energy and cross or temperature and rate
energy (array-like) – energy values of the cross-section
cross (array-like) – a cross-section
temperature (array-like)
rate (array-like)
ip (float) – the ionization potential. Typically in eV.
Returns – the following keys are added to inDict
——-
{‘csEnergy’, ‘csCross’, ‘ip’} or {‘csTemperature’, ‘csRate’, ‘ip’}
-
ChiantiPy.tools.util.
spectroscopic2name
(el, roman)¶ Convert from spectroscopic notation, e.g. Fe XI to ‘fe_11’.
- Parameters
el (str) – Elemental symbol, e.g. Fe for iron
roman (str) – Roman numeral spectroscopic symbol
-
ChiantiPy.tools.util.
splomDescale
(splom, energy)¶ Calculate the collision strength for excitation-autoionization as a function of energy.
- Parameters
energy (array-like) – In eV
splom (dict) – Structure returned by ChiantiPy.tools.io.splomRead
- Returns
omega – Collision strength
- Return type
array-like
-
ChiantiPy.tools.util.
units
(defaults)¶ to create a set of units compatible with ChiantiPy default values
-
ChiantiPy.tools.util.
z2element
(z)¶ Convert atomic number z to its elemental symbol.
- Parameters
z (int) – The atomic number/nuclear charge
- Returns
element – the abbreviated element name
- Return type
str
-
ChiantiPy.tools.util.
zion2dir
(z, ion, dielectronic=False, xuvtop='')¶ Convert atomic number and ion number to CHIANTI database directory.
- Parameters
z (int) – The atomic number/nuclear charge
ion (int) – the ionization stage, 1 for the neutral, 2 for the first ionization stage, …
dielectronic (bool, optional)
xuvtop (str, optional) – Set different CHIANTI database than the default
- Returns
fname – the CHIANTI directory where the file for the ion specified by z and ion are found
- Return type
str
-
ChiantiPy.tools.util.
zion2experimental
(z, ion, dielectronic=False)¶ Convert atomic number and ion number to spectroscopic notation string
- Parameters
z (int) – The atomic number/nuclear charge
ion (int) – the ionization stage, 0 for neutrals, 1 for singly ionized
dielectronic (bool, optional)
- Returns
expt – the experimental/laboratory notation for the ion, such as Fe 13+
- Return type
str
-
ChiantiPy.tools.util.
zion2filename
(z, ion, dielectronic=False, xuvtop='')¶ Convert nuclear charge/atomic number and ion number to CHIANTI database filename.
- Parameters
z (int) – The atomic number/nuclear charge
ion (int) – the ionization stage, 1 for neutrals, 2 for singly ionized
dielectronic (bool, optional) – whether the ion is the simple dielectronic model
xuvtop (str, optional) – the top directory of the CHIANTI database to be used
- Returns
fname – CHIANTI database filename
- Return type
str
-
ChiantiPy.tools.util.
zion2localFilename
(z, ion, dielectronic=False)¶ Convert atomic number and ion number to generic file name with current directory at top.
- Parameters
z (int) – The atomic number/nuclear charge
ion (int) – the ionization stage, 1 for neutrals, 2 for singly ionized
dielectronic (bool, optional)
-
ChiantiPy.tools.util.
zion2name
(z, ion, dielectronic=False)¶ Convert atomic number and ion number to generic name, e.g. (26,13) to ‘fe_13’
- Parameters
z (int) – The atomic number/nuclear charge
ion (int) – the ionization stage, 1 for the neutral, 2 for the first ionization stage, …
dielectronic (bool, optional)
- Returns
thisone – the CHIANTI style ion name, such as ‘fe_13’
- Return type
str
-
ChiantiPy.tools.util.
zion2spectroscopic
(z, ion, dielectronic=False)¶ Convert atomic number and ion number to spectroscopic notation string
- Parameters
z (int) – The atomic number/nuclear charge
ion (int) – the ionization stage, 1 for neutrals, 2 for singly ionized
dielectronic (bool, optional)
- Returns
spect – the spectroscopic notation for the ion, such as Fe XIV
- Return type
str
Module contents¶
Basic tools and utilities used in ChiantiPy