ChiantiPy.core package¶
Submodules¶
ChiantiPy.core.Continuum module¶
Continuum module
-
class
ChiantiPy.core.Continuum.
continuum
(ion_string, temperature, abundance=None, em=None)[source] [edit on github]¶ Bases:
object
The top level class for continuum calculations. Includes methods for the calculation of the free-free and free-bound continua.
Parameters: ionStr :
str
CHIANTI notation for the given ion, e.g. ‘fe_12’ that corresponds to the Fe XII ion.
temperature : array-like
In units of Kelvin
abundance :
float
orstr
, optionalElemental abundance relative to Hydrogen or name of CHIANTI abundance file, without the ‘.abund’ suffix, e.g. ‘sun_photospheric_1998_grevesse’.
em : array-like, optional
Line-of-sight emission measure (\(\int\mathrm{d}l\,n_en_H\)), in units of \(\mathrm{cm}^{-5}\), or the volumetric emission measure (\(\int\mathrm{d}V\,n_en_H\)) in units of \(\mathrm{cm}^{-3}\).
Notes
The methods for calculating the free-free and free-bound emission and losses return their result to an attribute. See the respective docstrings for more information.
Examples
>>> import ChiantiPy.core as ch >>> import numpy as np >>> temperature = np.logspace(4,9,20) >>> cont = ch.continuum('fe_15',temperature) >>> wavelength = np.arange(1,1000,10) >>> cont.freeFree(wavelength) >>> cont.freeBound(wavelength, include_abundance=True, include_ioneq=False) >>> cont.calculate_free_free_loss() >>> cont.calculate_free_bound_loss()
-
calculate_free_bound_loss
(**kwargs)[source] [edit on github]¶ Calculate the free-bound energy loss rate of an ion. The result is returned to the
free_bound_loss
attribute.The free-bound loss rate can be calculated by integrating the free-bound emission over the wavelength. This is difficult using the expression in
calculate_free_bound_emission
so we instead use the approach of [R1] and [R2]. Eq. 1a of [R2] can be integrated over wavelength to get the free-bound loss rate,\[\frac{dW}{dtdV} = C_{ff}\frac{k}{hc}T^{1/2}G_{fb},\]in units of erg \(\mathrm{cm}^3\,\mathrm{s}^{-1}\) where \(G_{fb}\) is the free-bound Gaunt factor as given by Eq. 15 of [R2] (see
mewe_gaunt_factor
for more details) and \(C_{ff}\) is the numerical constant as given in Eq. 4 of [R1] and can be written in terms of the fine structure constant \(\alpha\),\[C_{ff}\frac{k}{hc} = \frac{8}{3}\left(\frac{\pi}{6}\right)^{1/2}\frac{h^2\alpha^3}{\pi^2}\frac{k_B}{m_e^{3/2}} \approx 1.43\times10^{-27}\]References
[R1] (1, 2, 3) Gronenschild, E.H.B.M. and Mewe, R., 1978, A&AS, 32, 283 [R2] (1, 2, 3, 4) Mewe, R. et al., 1986, A&AS, 65, 511
-
freeBound
(wavelength, includeAbundance=True, includeIoneq=True, useVerner=True, **kwargs)[source] [edit on github]¶ Calculate the free-bound emission of an ion. The result is returned as a 2D array to the
free_bound_emission
attribute.The total free-bound continuum emissivity is given by,
\[\frac{dW}{dtdVd\lambda} = \frac{1}{4\pi}\frac{2}{hk_Bc^3m_e\sqrt{2\pi k_Bm_e}}\frac{E^5}{T^{3/2}}\sum_i\frac{\omega_i}{\omega_0}\sigma_i^{bf}\exp\left(-\frac{E - I_i}{k_BT}\right)\]where \(E=hc/\lambda\) is the photon energy, \(\omega_i\) and \(\omega_0\) are the statistical weights of the \(i^{\mathrm{th}}\) level of the recombined ion and the ground level of the recombing ion, respectively, \(\sigma_i^{bf}\) is the photoionization cross-section, and \(I_i\) is the ionization potential of level \(i\). This expression comes from Eq. 12 of [R5]. For more information about the free-bound continuum calculation, see Peter Young’s notes on free-bound continuum.
The photoionization cross-sections are calculated using the methods of [R4] for the transitions to the ground state and [R3] for all other transitions. See
verner_cross_section
andkarzas_cross_section
for more details.The free-bound emission is in units of erg \(\mathrm{cm}^3\mathrm{s}^{-1}\mathrm{\mathring{A}}^{-1}\mathrm{str}^{-1}\). If the emission measure has been set, the units will be multiplied by \(\mathrm{cm}^{-5}\) or \(\mathrm{cm}^{-3}\), depending on whether it is the line-of-sight or volumetric emission measure, respectively.
Parameters: wavelength : array-like
In units of angstroms
include_abundance :
bool
, optionalIf True, include the ion abundance in the final output.
include_ioneq :
bool
, optionalIf True, include the ionization equilibrium in the final output
use_verner :
bool
, optionalIf True, cross-sections of ground-state transitions using [R4], i.e.
verner_cross_section
Raises: ValueError
If no .fblvl file is available for this ion
References
[R3] (1, 2) Karzas and Latter, 1961, ApJSS, 6, 167 [R4] (1, 2, 3) Verner & Yakovlev, 1995, A&AS, 109, 125 [R5] (1, 2) Young et al., 2003, ApJSS, 144, 135
-
freeBoundLoss
(**kwargs)[source] [edit on github]¶ Calculate the free-bound energy loss rate of an ion. The result is returned to the
free_bound_loss
attribute.The free-bound loss rate can be calculated by integrating the free-bound emission over the wavelength. This is difficult using the expression in
calculate_free_bound_emission
so we instead use the approach of [R6] and [R7]. Eq. 1a of [R7] can be integrated over wavelength to get the free-bound loss rate,\[\frac{dW}{dtdV} = C_{ff}\frac{k}{hc}T^{1/2}G_{fb},\]in units of erg \(\mathrm{cm}^3\,\mathrm{s}^{-1}\) where \(G_{fb}\) is the free-bound Gaunt factor as given by Eq. 15 of [R7] (see
mewe_gaunt_factor
for more details) and \(C_{ff}\) is the numerical constant as given in Eq. 4 of [R6] and can be written in terms of the fine structure constant \(\alpha\),\[C_{ff}\frac{k}{hc} = \frac{8}{3}\left(\frac{\pi}{6}\right)^{1/2}\frac{h^2\alpha^3}{\pi^2}\frac{k_B}{m_e^{3/2}} \approx 1.43\times10^{-27}\]References
[R6] (1, 2, 3) Gronenschild, E.H.B.M. and Mewe, R., 1978, A&AS, 32, 283 [R7] (1, 2, 3, 4) Mewe, R. et al., 1986, A&AS, 65, 511
-
freeFree
(wavelength, include_abundance=True, include_ioneq=True, **kwargs)[source] [edit on github]¶ Calculates the free-free emission for a single ion. The result is returned as a dict to the
FreeFree
attribute. The dict has the keywordsintensity
,wvl
,temperature
,em
.The free-free emission for the given ion is calculated according Eq. 5.14a of [R8], substituting \(\nu=c/\lambda\), dividing by the solid angle, and writing the numerical constant in terms of the fine structure constant \(\alpha\),
\[\frac{dW}{dtdVd\lambda} = \frac{c}{3m_e}\left(\frac{\alpha h}{\pi}\right)^3\left(\frac{2\pi}{3m_ek_B}\right)^{1/2}\frac{Z^2}{\lambda^2T^{1/2}}\exp{\left(-\frac{hc}{\lambda k_BT}\right)}\bar{g}_{ff},\]where \(Z\) is the nuclear charge, \(T\) is the electron temperature in K, and \(\bar{g}_{ff}\) is the velocity-averaged Gaunt factor. The Gaunt factor is estimated using the methods of [R9] and [R10], depending on the temperature and energy regime. See
itoh_gaunt_factor
andsutherland_gaunt_factor
for more details.The free-free emission is in units of erg \(\mathrm{cm}^3\mathrm{s}^{-1}\mathrm{\mathring{A}}^{-1}\mathrm{str}^{-1}\). If the emission measure has been set, the units will be multiplied by \(\mathrm{cm}^{-5}\) or \(\mathrm{cm}^{-3}\), depending on whether it is the line-of-sight or volumetric emission measure, respectively.
Parameters: wavelength : array-like
In units of angstroms
include_abundance :
bool
, optionalIf True, include the ion abundance in the final output.
include_ioneq :
bool
, optionalIf True, include the ionization equilibrium in the final output
References
[R8] (1, 2) Rybicki and Lightman, 1979, Radiative Processes in Astrophysics, (Wiley-VCH) [R9] (1, 2) Itoh, N. et al., 2000, ApJS, 128, 125 [R10] (1, 2) Sutherland, R. S., 1998, MNRAS, 300, 321
-
freeFreeLoss
(**kwargs)[source] [edit on github]¶ Calculate the free-free energy loss rate of an ion. The result is returned to the
FreeFreeLoss
attribute.The free-free radiative loss rate is given by Eq. 5.15a of [R11]. Writing the numerical constant in terms of the fine structure constant \(\alpha\),
\[\frac{dW}{dtdV} = \frac{4\alpha^3h^2}{3\pi^2m_e}\left(\frac{2\pi k_B}{3m_e}\right)^{1/2}Z^2T^{1/2}\bar{g}_B\]where where \(Z\) is the nuclear charge, \(T\) is the electron temperature, and \(\bar{g}_{B}\) is the wavelength-averaged and velocity-averaged Gaunt factor. The Gaunt factor is calculated using the methods of [R12]. Note that this expression for the loss rate is just the integral over wavelength of Eq. 5.14a of [R11], the free-free emission, and is expressed in units of erg \(\mathrm{cm}^3\,\mathrm{s}^{-1}\).
References
[R11] (1, 2, 3) Rybicki and Lightman, 1979, Radiative Processes in Astrophysics, (Wiley-VCH) [R12] (1, 2) Karzas and Latter, 1961, ApJSS, 6, 167
-
free_free_loss
(**kwargs)[source] [edit on github]¶ Calculate the free-free energy loss rate of an ion. The result is returned to the
free_free_loss
attribute.The free-free radiative loss rate is given by Eq. 5.15a of [R13]. Writing the numerical constant in terms of the fine structure constant \(\alpha\),
\[\frac{dW}{dtdV} = \frac{4\alpha^3h^2}{3\pi^2m_e}\left(\frac{2\pi k_B}{3m_e}\right)^{1/2}Z^2T^{1/2}\bar{g}_B\]where where \(Z\) is the nuclear charge, \(T\) is the electron temperature, and \(\bar{g}_{B}\) is the wavelength-averaged and velocity-averaged Gaunt factor. The Gaunt factor is calculated using the methods of [R14]. Note that this expression for the loss rate is just the integral over wavelength of Eq. 5.14a of [R13], the free-free emission, and is expressed in units of erg \(\mathrm{cm}^3\,\mathrm{s}^{-1}\).
References
[R13] (1, 2, 3) Rybicki and Lightman, 1979, Radiative Processes in Astrophysics, (Wiley-VCH) [R14] (1, 2) Karzas and Latter, 1961, ApJSS, 6, 167
-
ioneqOne
()[source] [edit on github]¶ Provide the ionization equilibrium for the selected ion as a function of temperature. Similar to but not identical to ion.ioneqOne() returned in self.IoneqOne
-
ioneq_one
(stage, **kwargs)[source] [edit on github]¶ Calculate the equilibrium fractional ionization of the ion as a function of temperature.
Uses the
ChiantiPy.core.ioneq
module and does a first-order spline interpolation to the data. An ionization equilibrium file can be passed as a keyword argument,ioneqfile
. This can be passed through as a keyword argument to any of the functions that uses the ionization equilibrium.Parameters: stage : int
Ionization stage, e.g. 25 for Fe XXV
-
itoh_gaunt_factor
(wavelength)[source] [edit on github]¶ Calculates the free-free gaunt factors of [R15].
An analytic fitting formulae for the relativistic Gaunt factor is given by Eq. 4 of [R15],
\[g_{Z} = \sum^{10}_{i,j=0}a_{ij}t^iU^j\]where,
\[t = \frac{1}{1.25}(\log_{10}{T} - 7.25),\ U = \frac{1}{2.5}(\log_{10}{u} + 1.5),\]\(u=hc/\lambda k_BT\), and \(a_{ij}\) are the fitting coefficients and are read in using
ChiantiPy.tools.io.itohRead
and are given in Table 4 of [R15]. These values are valid for \(6<\log_{10}(T)< 8.5\) and \(-4<\log_{10}(u)<1\).See also
ChiantiPy.tools.io.itohRead
- Read in Gaunt factor coefficients from [R15]
References
[R15] (1, 2, 3, 4, 5) Itoh, N. et al., 2000, ApJS, 128, 125
-
karzas_cross_section
(photon_energy, ionization_potential, n, l)[source] [edit on github]¶ Calculate the photoionization cross-sections using the Gaunt factors of [R16].
The free-bound photoionization cross-section is given by,
\[\sigma_i^{bf} = 1.077294\times8065.54\times10^{16}\left(\frac{I_i}{hc}\right)^2\left(\frac{hc}{E}\right)^3\frac{g_{bf}}{n_i},\]where \(I_i\) is the ionization potential of the \(i^{\mathrm{th}}\) level, \(E\) is the photon energy, \(g_{bf}\) is the Gaunt factor calculated according to [R16], and \(n_i\) is the principal quantum number of the \(i^{\mathrm{th}}\) level. \(\sigma_i^{bf}\) is units of \(\mathrm{cm}^{2}\). This expression is given by Eq. 13 of [R17]. For more information on the photoionization cross-sections, see Peter Young’s notes on free-bound continuum.
Parameters: photon_energy : array-like
ionization_potential :
float
n :
int
l :
int
References
[R16] (1, 2, 3) Karzas and Latter, 1961, ApJSS, 6, 167 [R17] (1, 2) Young et al., 2003, ApJSS, 144, 135
-
klgfbInterp
(wvl, n, l)[source] [edit on github]¶ A Python version of the CHIANTI IDL procedure karzas_xs.
Interpolates free-bound gaunt factor of Karzas and Latter, (1961, Astrophysical Journal Supplement Series, 6, 167) as a function of wavelength (wvl).
-
mewe_gaunt_factor
(**kwargs)[source] [edit on github]¶ Calculate the Gaunt factor according to [R18] for a single ion \(Z_z\).
Using Eq. 9 of [R18], the free-bound Gaunt factor for a single ion can be written as,
\[G_{fb}^{Z,z} = \frac{E_H}{k_BT}\mathrm{Ab}(Z)\frac{N(Z,z)}{N(Z)}f(Z,z,n)\]where \(E_H\) is the ground-state potential of H, \(\mathrm{Ab}(Z)\) is the elemental abundance, \(\frac{N(Z,z)}{N(Z)}\) is the fractional ionization, and \(f(Z,z,n)\) is given by Eq. 10 and is approximated by Eq 16 as,
\[f(Z,z,n) \approx f_2(Z,z,n_0) = 0.9\frac{\zeta_0z_0^4}{n_0^5}\exp{\left(\frac{E_Hz_0^2}{n_0^2k_BT}\right)} + 0.42\frac{z^4}{n_0^{3/2}}\exp{\left(\frac{E_Hz^2}{(n_0 + 1)^2k_BT}\right)}\]where \(n_0\) is the principal quantum number, \(z_0\) is the effective charge (see Eq. 7 of [R18]), and \(\zeta_0\) is the number of vacancies in the 0th shell and is given in Table 1 of [R18]. Here it is calculated in the same manner as in fb_rad_loss.pro of the CHIANTI IDL library. Note that in the expression for \(G_{fb}\), we have not included the \(N_H/n_e\) factor.
Raises: ValueError
If no .fblvl file is available for this ion
References
[R18] (1, 2, 3, 4, 5) Mewe, R. et al., 1986, A&AS, 65, 511
-
sutherland_gaunt_factor
(wavelength)[source] [edit on github]¶ Calculates the free-free gaunt factor calculations of [R19].
The Gaunt factors of [R19] are read in using
ChiantiPy.tools.io.gffRead
as a function of \(u\) and \(\gamma^2\). The data are interpolated to the appropriate wavelength and temperature values usingmap_coordinates
.References
[R19] (1, 2, 3) Sutherland, R. S., 1998, MNRAS, 300, 321
-
verner_cross_section
(photon_energy)[source] [edit on github]¶ Calculates the photoionization cross-section using data from [R20] for transitions to the ground state.
The photoionization cross-section can be expressed as \(\sigma_i^{fb}=F(E/E_0)\) where \(F\) is an analytic fitting formula given by Eq. 1 of [R20],
\[F(y) = ((y-1)^2 + y_w^2)y^{-Q}(1 + \sqrt{y/y_a})^{-P},\]where \(E\) is the photon energy, \(n\) is the principal quantum number, \(l\) is the orbital quantum number, \(Q = 5.5 + l - 0.5P\), and \(\sigma_0,E_0,y_w,y_a,P\) are fitting paramters. These can be read in using
ChiantiPy.tools.io.vernerRead
.References
[R20] (1, 2, 3) Verner & Yakovlev, 1995, A&AS, 109, 125
-
ChiantiPy.core.Ion-2016-8-16 module¶
ChiantiPy.core.Ion-2017-5-7 module¶
ChiantiPy.core.Ion-diffs-2016-8-16 module¶
ChiantiPy.core.Ion module¶
Ion class
-
class
ChiantiPy.core.Ion.
ion
(ionStr, temperature=None, eDensity=None, pDensity='default', radTemperature=None, rStar=None, abundance=None, setup=True, em=None, verbose=0)[source] [edit on github]¶ Bases:
ChiantiPy.base.ionTrails
,ChiantiPy.base.specTrails
The top level class for performing spectral calculations for an ion in the CHIANTI database.
Parameters: ionStr :
str
CHIANTI notation for the given ion, e.g. ‘fe_12’ that corresponds to the
Fe XII
ion.temperature :
float
orndarray
, optionalTemperature array (Kelvin)
eDensity :
float
orndarray
, optionalElectron density array (\(\mathrm{cm^{-3}}\) )
pDensity :
float
orndarray
, optionalProton density (\(\mathrm{cm}^{-3}\) )
radTemperature :
float
orndarray
, optionalRadiation black-body temperature (in Kelvin)
rStar :
float
orndarray
, optionalDistance from the center of the star (in stellar radii)
abundance :
float
orstr
, optionalElemental abundance relative to Hydrogen or name of CHIANTI abundance file to use, without the ‘.abund’ suffix, e.g. ‘sun_photospheric_1998_grevesse’.
If True, run ion setup function. Otherwise, provide a limited number of attributes of the selected ion
em :
float
orndarray
, optionalEmission Measure, for the line-of-sight emission measure (\(\mathrm{\int \, n_e \, n_H \, dl}\)) (\(\mathrm{cm}^{-5}\).), for the volumetric emission measure \(\mathrm{\int \, n_e \, n_H \, dV}\) (\(\mathrm{cm^{-3}}\)).
Notes
The keyword arguments temperature, eDensity, radTemperature, rStar, em must all be either a float or have the same dimension as the rest if specified as lists, tuples or arrays.
The
Defaults
dict should have the following keys:- abundfile, the elemental abundance file, unless specified in chiantirc this defaults to sun_photospheric_1998_grevesse.
- ioneqfile, the ionization equilibrium file name. Unless specified in ‘chiantirc’ this is defaults to chianti. Other choices are availble in $XUVTOP/ioneq
- wavelength, the units of wavelength (Angstroms, nm, or kev), unless specified in the ‘chiantirc’ this is defaults to ‘angstrom’.
- flux, specified whether the line intensities are give in energy or photon fluxes, unless specified in the ‘chiantirc’ this is defaults to energy.
- gui, specifies whether to use gui selection widgets (True) or to make selections on the command line (False). Unless specified in the ‘chiantirc’ this is defaults to False.
Attributes
IonStr ( str
) Name of element plus ion, e.g.fe_12
for Fe XIIZ ( int
) the nuclear charge, 26 forfe_12
.Ion ( int
) the ionization stage, 12 forfe_12
.Dielectronic ( bool
) true if the ion is a ‘dielectronic’ ion where the levels are populated by dielectronic recombination.Spectroscopic ( str
) the spectroscopic notation for the ion, such asFe XII
forfe_12
.Filename ( str
) the complete name of the filegeneric
filename in the CHIANTI database, such as$XUVTOP/fe/fe_12/fe_12
.Ip ( float
) the ionization potential of the ionFIP ( float
) the first ionization potential of the elementDefaults ( dict
) these are specified by the software unless a chiantirc file is found in ‘$HOME/.chianti’:-
boundBoundLoss
(allLines=1)[source] [edit on github]¶ Calculate the summed radiative loss rate for all spectral lines of the specified ion.
Parameters: allLines :
bool
If True, include losses from both observed and unobserved lines. If False, only include losses from observed lines.
includes elemental abundance and ionization fraction.
Returns: creates the attribute:
BoundBoundLoss :
dict
with the keys below.rate : the radiative loss rate (\(\mathrm{erg \, cm^{-3}} \, \mathrm{s}^{-1}\)) per unit emission measure.
temperature : (K).
eDensity : electron density (\(\mathrm{cm^{-3}}\))
-
cireclvlDescale
(lvlType)[source] [edit on github]¶ Interpolate and extrapolate cilvl and reclvl rates. lvltype must be either ‘reclvl’, ‘cilvl’ or ‘rrlvl’ Used in level population calculations.
-
diCross
(energy=None, verbose=False)[source] [edit on github]¶ Calculate the direct ionization cross section as a function of the incident electron energy in eV, puts values into DiCross
-
diRate
()[source] [edit on github]¶ Calculate the direct ionization rate coefficient as a function of temperature (K)
-
drRate
()[source] [edit on github]¶ Provide the dielectronic recombination rate coefficient as a function of temperature (K).
-
drRateLvl
(verbose=0)[source] [edit on github]¶ to calculate the level resolved dielectronic rate from the higher ionization stage to the ion of interest rates are determined from autoionizing A-values the dictionary self.DrRateLvl contains rate = the dielectronic rate into an autoionizing level effRate = the dielectronic rate into an autoionizing level mutilplied by the branching ratio for a stabilizing transition totalRate = the sum of all the effRates
-
eaCross
(energy=None, verbose=False)[source] [edit on github]¶ Provide the excitation-autoionization cross section.
Energy is given in eV.
-
eaDescale
()[source] [edit on github]¶ Calculates the effective collision strengths (upsilon) for excitation-autoionization as a function of temperature.
-
eaRate
()[source] [edit on github]¶ Calculate the excitation-autoionization rate coefficient.
-
emiss
(allLines=True)[source] [edit on github]¶ Calculate the emissivities for lines of the specified ion.
units: ergs s^-1 str^-1
Does not include elemental abundance or ionization fraction
Wavelengths are sorted
set allLines = True to include unidentified lines
-
emissList
(index=-1, wvlRange=None, wvlRanges=None, top=10, relative=0, outFile=0)[source] [edit on github]¶ List the emissivities.
wvlRange, a 2 element tuple, list or array determines the wavelength range
Top specifies to plot only the top strongest lines, default = 10
normalize = 1 specifies whether to normalize to strongest line, default = 0
-
emissPlot
(index=-1, wvlRange=None, top=10, linLog='lin', relative=0, verbose=0, plotFile=0, saveFile=0)[source] [edit on github]¶ Plot the emissivities.
wvlRange, a 2 element tuple, list or array determines the wavelength range
Top specifies to plot only the top strongest lines, default = 10
linLog specifies a linear or log plot, want either lin or log, default = lin
normalize = 1 specifies whether to normalize to strongest line, default = 0
-
emissRatio
(wvlRange=None, wvlRanges=None, top=10)[source] [edit on github]¶ Plot the ratio of 2 lines or sums of lines. Shown as a function of density and/or temperature. For a single wavelength range, set wvlRange = [wMin, wMax] For multiple wavelength ranges, set wvlRanges = [[wMin1,wMax1],[wMin2,wMax2], …] A plot of relative emissivities is shown and then a dialog appears for the user to choose a set of lines.
-
gofnt
(wvlRange=0, top=10, verbose=0, plot=True)[source] [edit on github]¶ Calculate the ‘so-called’ G(T) function.
Given as a function of both temperature and eDensity.
Only the top( set by ‘top’) brightest lines are plotted. the G(T) function is returned in a dictionary self.Gofnt
-
intensity
(allLines=1)[source] [edit on github]¶ Calculate the intensities for lines of the specified ion.
units: ergs cm^-3 s^-1 str^-1
includes elemental abundance and ionization fraction.
the emission measure ‘em’ is included if specified
-
intensityRatioInterpolate
(data, scale='lin', plot=0, verbose=0)[source] [edit on github]¶ to take a set of date and interpolate against the IntensityRatio the scale can be one of ‘lin’/’linear’ [default], ‘loglog’, ‘logx’, ‘logy’,
-
ioneqOne
()[source] [edit on github]¶ Provide the ionization equilibrium for the selected ion as a function of temperature. returned in self.IoneqOne
-
ionizCross
(energy=None)[source] [edit on github]¶ Provides the total ionization cross section.
Notes
-
ionizRate
()[source] [edit on github]¶ Provides the total ionization rate.
Calls diRate and eaRate.
-
p2eRatio
()[source] [edit on github]¶ Calculates the proton density to electron density ratio using Eq. 7 of [R21].
Notes
Uses the abundance and ionization equilibrium.
References
[R21] (1, 2) Young, P. R. et al., 2003, ApJS, 144, 135
-
popPlot
(top=10, plotFile=0, outFile=0, pub=0)[source] [edit on github]¶ Plots populations vs temperature or eDensity.
top specifies the number of the most highly populated levels to plot if pub is set, the want publication plots (bw, lw=2).
-
populate
(popCorrect=1, verbose=0, **kwargs)[source] [edit on github]¶ Calculate level populations for specified ion. possible keyword arguments include temperature, eDensity, pDensity, radTemperature and rStar
-
recombRate
()[source] [edit on github]¶ Provides the total recombination rate coefficient.
-
rrRate
()[source] [edit on github]¶ Provide the radiative recombination rate coefficient as a function of temperature (K).
-
setup
(alternate_dir=None, verbose=False)[source] [edit on github]¶ Setup various CHIANTI files for the ion including .wgfa, .elvlc, .scups, .psplups, .reclvl, .cilvl, and others.
Parameters: alternate_dir :
str
directory cotaining the necessary files for a ChiantiPy ion; use to setup an ion with files not in the current CHIANTI directory
verbose :
bool
Notes
If ion is initiated with
setup=False
, call this method to do the setup at a later point.
-
setupIonrec
(alternate_dir=None, verbose=False)[source] [edit on github]¶ Setup method for ion recombination and ionization rates.
Notes
Allows a bare-bones ion object to be setup up with just the ionization and recombination rates. For ions without a complete set of files - one that is not in the MasterList.
-
spectrum
(wavelength, filter=(<function gaussianR>, 1000.0), label=0, allLines=1)[source] [edit on github]¶ Calculates the line emission spectrum for the specified ion.
Convolves the results of intensity to make them look like an observed spectrum the default filter is the gaussianR filter with a resolving power of 1000. Other choices include chianti.filters.box and chianti.filters.gaussian. When using the box filter, the width should equal the wavelength interval to keep the units of the continuum and line spectrum the same.
includes ionization equilibrium and elemental abundances
can be called multiple times to use different filters and widths uses label to keep the separate applications of spectrum sorted by the label for example, do .spectrum( …. labe=’test1’) and do .spectrum( …. label = ‘test2’) then will get self.Spectrum.keys() = test1, test2 and self.Spectrum[‘test1’] = {‘intensity’:aspectrum, ‘wvl’:wavelength, ‘filter’:useFilter.__name__, ‘filterWidth’:useFactor}
Notes
scipy.ndimage.filters also includes a range of filters.
-
twoPhoton
(wvl, verbose=False)[source] [edit on github]¶ to calculate the two-photon continuum - only for hydrogen- and helium-like ions includes the elemental abundance and the ionization equilibrium includes the emission measure if specified
-
twoPhotonEmiss
(wvl)[source] [edit on github]¶ To calculate the two-photon continuum rate coefficient - only for hydrogen- and helium-like ions
-
twoPhotonLoss
()[source] [edit on github]¶ to calculate the two-photon energy loss rate - only for hydrogen- and helium-like ions includes the elemental abundance and the ionization equilibrium does not include the emission measure
-
upsilonDescale
(prot=0)[source] [edit on github]¶ Provides the temperatures and effective collision strengths (upsilons) set prot for proton rates otherwise, ce will be set for electron collision rates uses the new format “scups” files
-
upsilonDescaleSplups
(prot=0, diel=0)[source] [edit on github]¶ Provides the temperatures and effective collision strengths (upsilons) set prot for proton rates otherwise, ce will be set for electron collision rates
ChiantiPy.core.Ioneq module¶
Ionization equilibrium class
-
class
ChiantiPy.core.Ioneq.
ioneq
(el_or_z)[source] [edit on github]¶ Bases:
object
Calculation ion fractions as a function of temperature assuming ionization equilibrium.
Parameters: Atomic number or symbol
-
calculate
(temperature)[source] [edit on github]¶ Calculate ion fractions for given temperature array using the total ionization and recombination rates.
-
load
(ioneqName=None)[source] [edit on github]¶ Read temperature and ion fractions from a CHIANTI “.ioneq” file.
-
plot
(stages=0, xr=0, yr=0, oplot=0, label=1, title=1, bw=0, semilogx=0)[source] [edit on github]¶ Plots the ionization equilibria.
self.plot(xr=None, yr=None, oplot=False) stages = sequence of ions to be plotted, neutral == 1, fully stripped == Z+1 xr = temperature range, yr = ion fraction range
for overplotting: oplot=”ioneqfilename” such as ‘mazzotta’ or if oplot=True or oplot=1 and a widget will come up so that a file can be selected.
-
ChiantiPy.core.IpyMspectrum module¶
-
ChiantiPy.core.IpyMspectrum.
doAll
(inpt)[source] [edit on github]¶ to process ff, fb and line inputs
-
class
ChiantiPy.core.IpyMspectrum.
ipymspectrum
(temperature, eDensity, wavelength, filter=(<function gaussianR>, 1000.0), label=None, elementList=None, ionList=None, minAbund=None, keepIons=0, doLines=1, doContinuum=1, allLines=1, em=None, abundanceName=0, verbose=0, timeout=0.1)[source] [edit on github]¶ Bases:
ChiantiPy.base.ionTrails
,ChiantiPy.base.specTrails
this is the multiprocessing version of spectrum for using inside an IPython Qtconsole or notebook.
be for creating an instance, it is necessary to type something like the following into a console
> ipcluster3 start –n=3 or > ipcluster start –n=3
this is the way to invoke things under the IPython 6 notation
Calculate the emission spectrum as a function of temperature and density.
temperature and density can be arrays but, unless the size of either is one (1), the two must have the same size
the returned spectrum will be convolved with a filter of the specified width on the specified wavelength array
the default filter is gaussianR with a resolving power of 100. Other filters, such as gaussian, box and lorentz, are available in ChiantiPy.filters. When using the box filter, the width should equal the wavelength interval to keep the units of the continuum and line spectrum the same.
A selection of elements can be make with elementList a list containing the names of elements that are desired to be included, e.g., [‘fe’,’ni’]
A selection of ions can be make with ionList containing the names of the desired lines in Chianti notation, i.e. C VI = c_6
Both elementList and ionList can not be specified at the same time
a minimum abundance can be specified so that the calculation can be speeded up by excluding elements with a low abundance. With solar photospheric abundances -
setting minAbund = 1.e-4 will include H, He, C, O, Ne setting minAbund = 2.e-5 adds N, Mg, Si, S, Fe setting minAbund = 1.e-6 adds Na, Al, Ar, Ca, Ni
Setting em will multiply the spectrum at each temperature by the value of em.
em [for emission measure], can be a float or an array of the same length as the temperature/density. allLines = 1 will include lines with either theoretical or observed wavelengths. allLines=0 will include only those lines with observed wavelengths
proc = the number of processors to use timeout - a small but non-zero value seems to be necessary
ChiantiPy.core.Mspectrum module¶
-
class
ChiantiPy.core.Mspectrum.
mspectrum
(temperature, eDensity, wavelength, filter=(<function gaussianR>, 1000.0), label=0, elementList=None, ionList=None, minAbund=None, keepIons=0, abundance=None, doLines=1, doContinuum=1, allLines=1, em=None, proc=3, verbose=0, timeout=0.1)[source] [edit on github]¶ Bases:
ChiantiPy.base.ionTrails
,ChiantiPy.base.specTrails
this is the multiprocessing version of spectrum set proc to the desired number of processors, default=3
Calculate the emission spectrum as a function of temperature and density.
temperature and density can be arrays but, unless the size of either is one (1), the two must have the same size
the returned spectrum will be convolved with a filter of the specified width on the specified wavelength array
the default filter is gaussianR with a resolving power of 100. Other filters, such as gaussian, box and lorentz, are available in chianti.filters. When using the box filter, the width should equal the wavelength interval to keep the units of the continuum and line spectrum the same.
A selection of elements can be make with elementList a list containing the names of elements that are desired to be included, e.g., [‘fe’,’ni’]
A selection of ions can be make with ionList containing the names of the desired lines in Chianti notation, i.e. C VI = c_6
Both elementList and ionList can not be specified at the same time
a minimum abundance can be specified so that the calculation can be speeded up by excluding elements with a low abundance. With solar photospheric abundances -
setting minAbund = 1.e-4 will include H, He, C, O, Ne setting minAbund = 2.e-5 adds N, Mg, Si, S, Fe setting minAbund = 1.e-6 adds Na, Al, Ar, Ca, Ni
Setting em will multiply the spectrum at each temperature by the value of em.
em [for emission measure], can be a float or an array of the same length as the temperature/density. allLines = 1 will include lines with either theoretical or observed wavelengths. allLines=0 will include only those lines with observed wavelengths
proc = the number of processors to use timeout - a small but non-zero value seems to be necessary
ChiantiPy.core.RadLoss module¶
-
class
ChiantiPy.core.RadLoss.
radLoss
(temperature, eDensity, elementList=0, ionList=0, minAbund=0, doContinuum=1, abundance=None, verbose=0, allLines=1, keepIons=0)[source] [edit on github]¶ Bases:
ChiantiPy.base.specTrails
Calculate the emission spectrum as a function of temperature and density.
includes elemental abundances or ionization equilibria
temperature and density can be arrays but, unless the size of either is one (1), the two must have the same size
the returned spectrum will be convolved with a filter of the specified width on the specified wavelength array
the default filter is gaussianR with a resolving power of 1000. Other filters, such as gaussian, box and lorentz, are available in ChiantiPy.filters. When using the box filter, the width should equal the wavelength interval to keep the units of the continuum and line spectrum the same.
A selection of ions can be make with ionList containing the names of the desired lines in Chianti notation, i.e. C VI = c_6
a minimum abundance can be specified so that the calculation can be speeded up by excluding elements with a low abundance. With solar photospheric abundances -
setting minAbund = 1.e-4 will include H, He, C, O, Ne setting minAbund = 2.e-5 adds N, Mg, Si, S, Fe setting minAbund = 1.e-6 adds Na, Al, Ar, Ca, Ni
Setting em will multiply the spectrum at each temperature by the value of em.
em [for emission measure], can be a float or an array of the same length as the temperature/density.
- abundance: to select a particular set of abundances, set abundance to the name of a CHIANTI abundance file,
- without the ‘.abund’ suffix, e.g. ‘sun_photospheric_1998_grevesse’ If set to a blank (‘’), a gui selection menu will popup and allow the selection of an set of abundances
-
radLossPlot
(title=0)[source] [edit on github]¶ to plot the radiative losses vs temperature
ChiantiPy.core.Spectrum module¶
-
class
ChiantiPy.core.Spectrum.
bunch
(temperature, eDensity, wvlRange, elementList=None, ionList=None, minAbund=None, keepIons=0, em=None, abundanceName=None, verbose=0, allLines=1)[source] [edit on github]¶ Bases:
ChiantiPy.base.ionTrails
,ChiantiPy.base.specTrails
Calculate the emission line spectrum as a function of temperature and density.
‘bunch’ is very similar to ‘spectrum’ except that continuum is not calculated and the spectrum is not convolved over a filter. However, this can be done with the inherited convolve method
one of the convenient things is that all of the instantiated ion classes, determined through such keywords as ‘elementList’, ‘ionList’, and ‘minAbund’ are kept in a dictionary self.IonInstances where self.IonInstances[‘mg_7’] is the class instance of ChiantiPy.core.ion for ‘mg_7’. All its methods and attributes are available.
includes elemental abundances and ionization equilibria
the set of abundances, a file in $XUVTOP/abundance, can be set with the keyword argument ‘abundanceName’
temperature and density can be arrays but, unless the size of either is one (1), the two must have the same size
Inherited methods include ‘intensityList’, ‘intensityRatio’ (between lines of different ions), and ‘intensityRatioSave’ and ‘convolve’.
A selection of elements can be make with elementList a list containing the names of elements that are desired to be included, e.g., [‘fe’,’ni’]
A selection of ions can be make with ionList containing the names of the desired lines in Chianti notation, i.e. C VI = c_6
Both elementList and ionList can not be specified at the same time
a minimum abundance can be specified so that the calculation can be speeded up by excluding elements with a low abundance. With solar photospheric abundances -
setting minAbund = 1.e-4 will include H, He, C, O, Ne setting minAbund = 2.e-5 adds N, Mg, Si, S, Fe setting minAbund = 1.e-6 adds Na, Al, Ar, Ca, Ni
At least one of elementList, ionList, or minAbund must be set in order for ‘bunch’ to include any ions.
Setting em will multiply the spectrum at each temperature by the value of em.
em [for emission measure], can be a float or an array of the same length as the temperature/density
-
class
ChiantiPy.core.Spectrum.
spectrum
(temperature, eDensity, wavelength, filter=(<function gaussianR>, 1000.0), label=None, elementList=None, ionList=None, minAbund=None, doLines=1, doContinuum=1, em=None, keepIons=0, abundance=None, verbose=0, allLines=1)[source] [edit on github]¶ Bases:
ChiantiPy.base.ionTrails
,ChiantiPy.base.specTrails
Calculate the emission spectrum as a function of temperature and density.
one of the convenient things is that all of the instantiated ion classes, determined through such keywords as ‘elementList’, ‘ionList’, and ‘minAbund’ are kept in a dictionary self.IonInstances where self.IonInstances[‘mg_7’] is the class instance of ChiantiPy.core.ion for ‘mg_7’. All its methods and attributes are available.
includes elemental abundances and ionization equilibria
the set of abundances, a file in $XUVTOP/abundance, can be set with the keyword argument ‘abundanceName’
temperature and density can be arrays but, unless the size of either is unity (1), the two must have the same size
the returned spectrum will be convolved with a filter of the specified width on the specified wavelength array
the default filter is gaussianR with a resolving power of 1000. Other filters, such as gaussian, box and lorentz, are available in ChiantiPy.tools.filters. When using the box filter, the width should equal the wavelength interval to keep the units of the continuum and line spectrum the same.
Inherited methods include ‘intensityList’, ‘intensityRatio’ (between lines of different ions), ‘intensityRatioSave’ and ‘convolve’
A selection of elements can be make with elementList a list containing the names of elements that are desired to be included, e.g., [‘fe’,’ni’]
A selection of ions can be make with ionList containing the names of the desired lines in CHIANTI notation, i.e. C VI = c_6
Both elementList and ionList can not be specified at the same time
a minimum abundance can be specified so that the calculation can be speeded up by excluding elements with a low abundance. The default of minAbund is 1.e-6
It is necessary to specify at least an elementList, an ionList, or a minAbund to select any ions for a spectrum calculation
With solar photospheric abundances -
setting minAbund = 1.e-4 will include H, He, C, O, Ne setting minAbund = 2.e-5 adds N, Mg, Si, S, Fe setting minAbund = 1.e-6 adds Na, Al, Ar, Ca, Ni
Setting doLines = 0 will skip the calculation of spectral lines. Setting doContinuum =0 will skip the continuum calculation.
Setting em [for emission measure] will multiply the spectrum at each temperature by the value of em.
em [for emission measure] can be a float or an array of the same length as the temperature/density
keepIons: set this to keep the ion instances that have been calculated in a dictionary self.IonInstances with the keywords being the CHIANTI-style ion names
abundance: to select a particular set of abundances, set abundance to the name of a CHIANTI abundance file, without the ‘.abund’ suffix, e.g. ‘sun_photospheric_1998_grevesse’
If set to a blank (‘’), a gui selection menu will popup and allow the selection of an set of abundances
Module contents¶
chianti.core - contains the main classes for ChiantiPy users.
This software is distributed under the terms of the ISC Software License that is found in the LICENSE file