This page was created by the IDL library routine
mk_html_help. For more information on
this routine, refer to the IDL Online Help Navigator
or type:
? mk_html_help
at the IDL command line prompt.
Last modified: Tue Jan 4 12:52:37 2011.
PROJECT: CHIANTI
CHIANTI is an Atomic Database Package for Spectroscopic Diagnostics
of Astrophysical Plasmas. It is a collaborative project involving
the Naval Research Laboratory (USA), the University of Florence
(Italy), the University of Cambridge and the Rutherford Appleton
Laboratory (UK).
NAME:
FB_RAD_LOSS
PURPOSE:
Calculate the total radiative losses of a plasma due to the
free-bound (radiative recombination) continuum.
EXPLANATION
This routine does not use the same method of calculating the ion
continuum emissivities as the FREEBOUND routine. This is because a
modified version of FREEBOUND would be very slow for this purpose.
Instead, we use the method of Mewe et al. (A&AS 65, 511, 1986) which
is outlined in their Sect.2.2. The integration over the quantity
P_c(lambda,T) is very simple as the gaunt factor, G_c, has no
intrinsic lambda dependence other than through the limits.
Comparisons between the wavelength-resolved continuum emission
derived from the Mewe et al. method with the more sophisticated
method employed in FREEBOUND show excellent agreement with at most
10% differences at specific wavelengths.
CALLING SEQUENCE:
fb_rad_loss, temp, int, min_abund=min_abund, /no_setup
OUTPUTS
TEMP Temperatures (K). These are the temperatures at which the
ion fractions are defined (typically 10^4 to 10^8 in 0.1
dex intervals).
INT The emissivity in units of erg cm^3 s^-1.
OPTIONAL INPUTS:
MIN_ABUND Exclude elements whose abundances are less than MIN_ABUND.
(Note that Ab(H)=1.)
KEYWORD PARAMETERS:
NO_SETUP If the procedure setup_elements has already been called
then the keyword /nosetup should be set to avoid
repeating this step
PROGRAMMING NOTES
This routine computes the free-bound gaunt factor following the
prescription set out in Sect.2.2 of Mewe et al. (1986). The
expression for f_2 (Eq. 16 of Mewe et al.) contains several
quantities that need to be computed. Zeta_0 is computed internally
by the function ZETA_0 through a prescription evident from browsing
Table I of Mewe et al. Z_0 is computed from the ionization potential
of the recombined ion which is contained in the .ip file within the
database. The quantity n_0 (also used in deriving Z_0) is derived
using the routine CONF2N which extracts the highest n value from
the configuration description of the ground term of each ion.
The quantity Z is just the charge on the recombined ion; E_0 is the
ionization potential of the recombined ion, while E_n_0+1 is derived
from Mewe's Eq.7 with the prescription that z_n=Z and n=n_0+1.
The ions that are considered for the continuum are those for which
.fblvl files exist.
COMMON BLOCKS:
common elements,abund,abund_ref,ioneq,ioneq_t,ioneq_ref
INTERNAL FUNCTIONS
ZETA_0
CALLS
READ_IP, ZION2FILENAME, SETUP_ELEMENTS, FILE_EXIST,
READ_FBLVL, CONCAT_DIR, GET_IEQ
EXAMPLES
IDL> fb_rad_loss,temp,int
IDL> plot,temp,int,/xlog,/ylog
MODIFICATION HISTORY:
Ver.1, 1-Aug-2002, Peter Young
Completely new version of fb_rad_loss. Incorporates code from
a version of freebound not available in CHIANTI.
V 2, 25-May-2005, GDZ
corrected routine header.
VERSION: 2, 25-May-2005
(See continuum/fb_rad_loss.pro)
NAME
ZETA_0
EXPLANATION
Returns the value of zeta_0 (the number of vacancies in the ion given
by IZ and ION). See Sect. 2.2 of Mewe et al. (1986, A&AS 65, 511).
INPUTS
IZ Atomic number of ion (e.g., 26 -> Fe)
ION Spectroscopic number of ion (e.g., 13 -> XIII)
OUTPUTS
Value of zeta_0.
(See continuum/fb_rad_loss.pro)
PROJECT: CHIANTI
CHIANTI is an atomic database package for the calculation of
continuum and emission line spectra from astrophysical plasmas. It is a
collaborative project involving the Naval Research Laboratory
(Washington DC, USA), the Arcetri Observatory (Firenze, Italy), and the
Cambridge University (United Kingdom).
NAME:
FF_RAD_LOSS
PURPOSE:
Calculate the free-free radiative energy losses losses.
Uses the free-free integrated gaunt factor calculations of
Sutherland, 1998, MNRAS, 300, 321
CALLING SEQUENCE:
FF_RAD_LOSS,Temperature,LossRate
INPUTS:
OPTIONAL INPUTS:
None
KEYWORD PARAMETERS:
NO_SETUP: If the procedure setup_elements has already been called then
the keyword /no_setup should be set to avoid repeating this step
MIN_ABUND: If set, calculates the continuum only from those elements which
have an abundance greater than min_abund. Can speed up the
calculations. For example:
abundance (H) = 1.
abundance (He) = 0.085
abundance (C) = 3.3e-4
abundance (Si) = 3.3e-5
abundance (Fe) = 3.9e-5
OUTPUTS:
Temperature: temperature in degrees Kelvin, can be a 1 D array
LossRate: radiative energy loss rate in erg s^-1 cm^3
(radiative loss rate per emission measure (N_e N_H V)
COMMON BLOCKS:
common elements,abund,abund_ref,ioneq,ioneq_t,ioneq_ref
EXAMPLE:
> ff_rad_loss,t,rad
> ff_rad_loss,t,rad,min_abund=3.e-5
> ff_rad_loss,t,rad,/no_setup,min_abund=1.e-6
MODIFICATION HISTORY:
Written by: Ken Dere
April 2000: Version 3.0
(See continuum/ff_rad_loss.pro)
PROJECT: CHIANTI
CHIANTI is an Atomic Database Package for Spectroscopic Diagnostics of
Astrophysical Plasmas. It is a collaborative project involving the Naval
Research Laboratory (USA), the University of Florence (Italy), the
University of Cambridge and the Rutherford Appleton Laboratory (UK).
NAME
FREEBOUND
PURPOSE:
Calculates the free-bound (radiative recombination) continuum.
A report giving detailed information on how the CHIANTI
free-bound is calculated is available at:
http://solar.bnsc.rl.ac.uk/~young/chianti/freebound.pdf
To calculate the continuum spectrum the user should specify
whether the plasma structure is specified through a differential
emission measure (DEM) curve or through emission measure (EM)
values. For the former, the DEM values are specified through the
DEM_INT keyword and for the latter the EM values are specified
through the EM_INT keyword. If neither are specified EM_INT
values of unity are assumed.
If DEM_INT values are specifed they must be specified on a
temperature scale that has constant intervals in log10(T). EM_INT
values do not have this restriction.
Note that EM_INT must be used for isothermal plasmas.
INPUTS
TEMP Temperature in K (can be an array). If the keyword
/DT is used, the
temperatures must be spaced at equal intervals of
log10(T). E.g., log T=6.0, 6.1, 6.2, 6.3, ...
WVL Wavelength in angstroms (can be an array). If /keV is
set, then WVL is assumed to contain energies in keV units.
OUTPUTS
INT Free-bound continuum intensity in units of
10^-40 erg cm^3/s/sr/Angstrom
( integral(N_H N_e dh) in cm^-5) if a DEM is not defined.
If DEM values are defined, it is assumed that they are given
as N_H N_e dh/dT. The units are 10^-40 erg/cm^2/s/srAngstrom
If T is given as a 1-D array, then the output will be a
2-D array, with one element for each temperature and
wavelength (but also see SUMT).
If the keyword /keV is set, then the units of INT will be
10^-40 erg cm^3/s/sr/keV.
OPTIONAL INPUTS
DEM_INT This should be the same size as TEMP and contain
differential emission measure values. Specifying
DEM_INT places a restriction on TEMP (see above). Note
that DEM_INT values are multiplied by the factor dT
within the routine, whereas EM_INT values are not.
EM_INT This should be the same size as TEMP and contain
emission measure values. Note
that DEM_INT values are multiplied by the factor dT
within the routine, whereas EM_INT values are not. If
neither EM_INT or DEM_INT are specified, then a EM_INT
vector of same size as TEMP with values of unity is
assumed.
IZ Only calculate continuum for the element with atomic
number IZ
ION (To be used in conjunction with IZ.) Calculated continuum
for a single ion (IZ, ION).
KEYWORDS
NO_SETUP If the procedure setup_elements has already been called
then the keyword /nosetup should be set to avoid
repeating this step
MIN_ABUND If set, calculates the continuum only from those
elements which have an abundance greater than
min_abund. Can speed up the calculations. For
example:
abundance (H) = 1.
abundance (He) = 0.085
abundance (C) = 3.3e-4
abundance (Si) = 3.3e-5
abundance (Fe) = 3.9e-5
PHOTONS The output spectrum is given in photon units rather
than ergs.
SUMT When a set of temperatures is given to FREEBOUND, the
default is to output INTENSITY as an array of size
(nwvl x nT). With this keyword set, a summation over
the temperatures is performed.
VERBOSE Output information from FREEBOUND.
KEV If set, then WVL is assumed to contain energies in keV units
rather than wavelengths in Angstroms.
COMMON BLOCKS
ELEMENTS
CALLS
FREEBOUND_ION, SETUP_ELEMENTS, READ_KLGFB, GET_IEQ
HISTORY
Ver.1, 24-Jul-2002, Peter Young
Ver.2, 26-Jul-2002, Peter Young
revised call to freebound_ion; corrected ion fraction problem
V 3, 25-May-2005, GDZ
corrected routine header.
Ver.4, 10-Mar-2006, Peter Young
added keV keyword
Ver.5, 8-Sep-2009, Peter Young
routine now checks to make sure temperatures are evenly
spaced when the DEM_INT keyword is used; if an error is
found due to incorrect temperature specificiations, then
the continuum intensity is returned as zero.
Ver.6, 6-Oct-2009, Peter Young
Introduced EM_INT keyword for specifying emission measure
values for isothermal plasmas.
VERSION : 6, 6-Oct-2009
(See continuum/freebound.pro)
PROJECT: CHIANTI
CHIANTI is an Atomic Database Package for Spectroscopic Diagnostics of
Astrophysical Plasmas. It is a collaborative project involving the Naval
Research Laboratory (USA), the University of Florence (Italy), the
University of Cambridge and the Rutherford Appleton Laboratory (UK).
NAME
FREEBOUND
EXPLANATION
Calculates the free-bound (radiative recombination) continuum
from a single ion. Note that the output does not contain the ion
fraction, element abundance or differential emission measure term.
INPUTS
TEMP Temperature in K (can be an array).
WVL Wavelength in angstroms (can be an array). If the keyword
/keV is set, then wvl should instead be given in units of
keV.
IZ Atomic number of ion (e.g., 26 = Fe)
ION Spectroscopic number of ion (e.g., 13 = XIII)
OUTPUTS
INT Free-bound continuum intensity. Needs to be multiplied by
element abundance, ion fraction and DEM to obtain the final
continuum intensity.
OPTIONAL INPUTS
IP The ionization potential of the ion.
VDATA An array containing the Verner & Yakovlev data array.
PE An array containing the PE data from READ_FBLVL
KLGFB An array containing the KLGFB data from READ_FBLVL
[Note: the above 3 inputs are used when calling freebound_ion from
freebound]
KEYWORDS
NOVERNER If set, then the Verner & Yakovlev cross-sections will
not be used.
KEV If set, then WVL is assumed to contain energies in keV units
rather than wavelengths in Angstroms.
COMMON BLOCKS
None.
CALLS
READ_FBLVL, ZION2FILENAME, VERNER_XS, KARZAS_XS, CONCAT_DIR,
READ_IP, READ_KLGFB, FILE_EXIST
PROGRAMMING NOTES
The way I treat the exponential function in the expression for the
emissivity may seem strange, but it saves a bit of time in the
calculation. Basically calculating exp(E-IP+E_i/T) for each level
was time consuming so I split it into exp(E-IP/T)*exp(E_i/T). The
first term comes out of the for loop, while the second term is the
exponential of a vector rather than an array. This made a ~30%
time-saving for freebound.
HISTORY
Ver.1, 24-Jul-2002, Peter Young
Ver.2, 26-Jul-2002, Peter Young
Added /noverner keyword.
Ver.3, 30-Jul-2002, Peter Young
Speeded up routine by modifying treatment of exponential.
Ver.4, 10-Mar-2006, Peter Young
Added /keV keyword
Ver. 5 8-Jun-2010, Ken Dere
added check for n_params()
(See continuum/freebound_ion.pro)
PROJECT : CHIANTI
CHIANTI is an Atomic Database Package for Spectroscopic Diagnostics of
Astrophysical Plasmas. It is a collaborative project involving the Naval
Research Laboratory (USA), the University of Florence (Italy), the
University of Cambridge and the Rutherford Appleton Laboratory (UK).
NAME
FREEFREE
PURPOSE:
This routine computes the free-free continuum (bremsstrahlung)
using the fitting formulae of Itoh et al. (ApJS 128, 125, 2000)
and Sutherland (MNRAS 300, 321, 1998).
The Itoh et al. data are valid for smaller ranges for temperature
and wavelength than Sutherland and so for points outside of their
ranges we use the data of Sutherland.
INPUTS
TEMP Temperature (in K). They do not have to be in monotonically
increasing order.
WVL Wavelengths in angstroms. Can be a scalar or vector. If /keV is
set, then WVL is assumed to contain energies in keV units.
OUTPUTS
INT Free-free continuum intensity in units of
10^-40 erg cm^3/s/sr/Angstrom
[ integral(N_H N_e dh) in cm^-5 ] if a DEM is not defined.
If DEM values are defined, it is assumed that they are given
as N_H N_e dh/dT. The units are 10^-40 erg/cm^2/s/sr/Angstrom.
If T is given as a 1-D array, then the output will be a 2-D array,
with one element for each temperature and wavelength
(but also see SUMT).
If the keyword /keV is set, then the units of INT will be
10^-40 erg cm^3/s/sr/keV
OPTIONAL INPUTS
DEM_INT This should be the same size as TEMP and contain
differential emission measure values. Specifying
DEM_INT places a restriction on TEMP (see above). Note
that DEM_INT values are multiplied by the factor dT
within the routine, whereas EM_INT values are not.
EM_INT This should be the same size as TEMP and contain
emission measure values. Note
that DEM_INT values are multiplied by the factor dT
within the routine, whereas EM_INT values are not. If
neither EM_INT or DEM_INT are specified, then a EM_INT
vector of same size as TEMP with values of unity is
assumed.
MIN_ABUND This keyword allows the specification of a minimum abundance,
such that any elements with an abundance (relative to
hydrogen) less than MIN_ABUND will not be included in the
calculation. E.g., MIN_ABUND=1e-5.
KEYWORDS
NO_SETUP By default the routine asks the user which ion balance
and abundance files to use via pop-up widgets. If
/no_setup is used then this data is taken from the common
block.
SUMT The default is to output the intensity array as an array
of size (nwvl x nT). Setting this keyword performs a sum
over the temperatures to yield a vector of same size as
the input wavelengths, thus producing the complete
free-free spectrum.
PHOTONS Gives output emissivity in photon units rather than ergs.
KEV If set, then WVL is assumed to contain energies in keV units
rather than wavelengths in Angstroms.
CALLS
SUTHERLAND, ITOH
COMMON BLOCKS
ELEMENTS
PROGRAMMING NOTES
The Itoh fitting formula is only valid for (6.0 LE logT LE 8.5).
For temperatures below this, we thus switch to the Sutherland
fitting formula. There is very little (<1%) between the two at
logT=6.
Itoh also has a constraint on the quantity u=hc/kTl (l=wavelength),
such that (-4 LE log u LE 1.0). The upper limit corresponds to the
continuum being cut-off prematurely at low wavelengths. E.g., for
T=10^6 the cutoff is at 14.39 angstroms. For these low wavelengths
we also use the Sutherland data to complete the continuum. Note that
the continuum at these wavelengths is very weak
MODIFICATION HISTORY
Ver.1, 5-Dec-2001, Peter Young
Completely revised to call the separate itoh.pro and
sutherland.pro routines.
V. 2, 21-May-2002, Giulio Del Zanna (GDZ),
Corrected the description of the units.
Added verbose keyword and a printout.
V. 3, 22-May-2002, Peter Young (PRY)
Added MIN_ABUND optional input.
Changed ioneq_t to ioneq_logt (GDZ).
V 4, 25-May-2005, GDZ
corrected routine header.
V. 5, 9-Mar-2006, Peter Young
Added /kev keyword
V. 6, 9-May-2008, Peter Young & Kevin Tritz
Modifed code so that temperatures can be specified in a unordered
list.
V. 7, 7-Oct-2009, Peter Young
Added EM_INT= keyword; simplified how the routine decides
which calculation to use (Itoh or Sutherland)
v. 8, 8-Jun-2010, Ken Dere
added print statement following check for params
VERSION : 8, 8-Jun-2010
(See continuum/freefree.pro)
PROJECT: CHIANTI
CHIANTI is an Atomic Database Package for Spectroscopic Diagnostics of
Astrophysical Plasmas. It is a collaborative project involving the Naval
Research Laboratory (USA), the University of Florence (Italy), the
University of Cambridge and the Rutherford Appleton Laboratory (UK).
NAME
ITOH
EXPLANATION
Calculates the relativistic free-free continuum using the fitting
formula of Itoh et al. (ApJS 128, 125, 2000).
To calculate the continuum spectrum the user should specify
whether the plasma structure is specified through a differential
emission measure (DEM) curve or through emission measure (EM)
values. For the former, the DEM values are specified through the
DEM_INT keyword and for the latter the EM values are specified
through the EM_INT keyword. If neither are specified EM_INT
values of unity are assumed.
If DEM_INT values are specifed they must be specified on a
temperature scale that has constant intervals in log10(T). EM_INT
values do not have this restriction.
Note that EM_INT must be used for isothermal plasmas.
INPUTS
TEMP Temperature (in K). If the keyword DEM_INT= is set, then the
temperatures must be specified with a fixed spacing in
logT. E.g., TEMP=10.^[4.0,4.1,4.2].
WVL Wavelengths in angstroms. Can be a scalar or vector. If /keV is
set, then WVL is assumed to contain energies in keV units.
OUTPUTS
INT Free-free continuum emissivity in units of
10^-40 erg cm^3 / s / sr / Angstrom per unit emission
measure [ integral(N_e N_H dh) in cm^-5 ]. If T is given as
a 1-D array, then RAD will be output as a 2-D array,
with one element for each temperature and wavelength
(but also see SUMT).
If the keyword /keV is set, then the units of INT will be
10^-40 erg cm^3/s/sr/keV
OPTIONAL INPUTS
DEM_INT This should be the same size as TEMP and contain
differential emission measure values. Specifying
DEM_INT places a restriction on TEMP (see above). Note
that DEM_INT values are multiplied by the factor dT
within the routine, whereas EM_INT values are not.
EM_INT This should be the same size as TEMP and contain
emission measure values. Note
that DEM_INT values are multiplied by the factor dT
within the routine, whereas EM_INT values are not. If
neither EM_INT or DEM_INT are specified, then a EM_INT
vector of same size as TEMP with values of unity is
assumed.
MIN_ABUND This keyword allows the specification of a minimum abundance,
such that any elements with an abundance (relative to
hydrogen) less than MIN_ABUND will not be included in the
calculation. E.g., MIN_ABUND=1e-5.
KEYWORDS
NO_SETUP By default the routine asks the user which ion balance
and abundance files to use via pop-up widgets. If
/no_setup is used then this data is taken from the common
block.
SUMT The default is to output the intensity array as an array
of size (nwvl x nT). Setting this keyword performs a sum
over the temperatures to yield a vector of same size as
the input wavelengths, thus producing the complete
free-free spectrum.
PHOTONS Gives output emissivity in photon units rather than ergs.
KEV If set, then WVL is assumed to contain energies in keV units
rather than wavelengths in Angstroms.
CALLS
SETUP_ELEMENTS
PROGRAMMING NOTES
In the Itoh et al. paper they state that their fitting formula is
valid for ion charges (Z_j) of 1-28. The data-file they provide
actually goes up to Z_j=30, and so you will see that the loop over
z (see below) goes up to 30.
There is no restriction on the elements for which the fitting
formula is valid, and so all elements are allowed to be included
(subject to their abundances being non-zero).
HISTORY
Ver.1, 3-Dec-2001, Peter Young
Ver.2, 22-May-2002, Peter Young
Added MIN_ABUND optional input.
Ver.3, 28-May-2002, Peter Young
Corrected way in which DEM is handled.
Ver.4, 9-Mar-2006, Peter Young
Added /keV keyword.
Ver.5, 10-Sep-2009, Peter Young
Routine now checks that temperatures are evenly spaced
when DEM_INT is specified. Also 0.1 dex spacing is no
longer assumed.
Ver.6, 6-Oct-2009, Peter Young
Introduced EM_INT keyword for specifying emission measure
values for isothermal plasmas.
(See continuum/itoh.pro)
PROJECT: CHIANTI
CHIANTI is an Atomic Database Package for Spectroscopic Diagnostics of
Astrophysical Plasmas. It is a collaborative project involving the Naval
Research Laboratory (USA), the University of Florence (Italy), the
University of Cambridge and the Rutherford Appleton Laboratory (UK).
NAME
KARZAS_XS()
EXPLANATION
Outputs the photoionization cross-section at the wavelengths WVL for
ionization of an N,L electron with ionization energy IONIZ_EN,
calculated using the Karzas & Latter (ApJS 6, 167, 1961) formulation.
The bound-free gaunt factor is derived from the tables of Karzas &
Latter.
INPUTS
WVL Wavelengths (in angstroms) for which cross-sections are
required (1-D array).
N Principal quantum number of the electron being removed in the
ionization.
L The orbital angular momentum of the electron being removed in the
ionization.
IONIZ_EN The ionization energy (in cm^-1) of the electron being
removed.
OPTIONAL INPUTS
PE Allows the PE array from READ_KLGFB to be directly input to
KARZAS_XS thus avoiding the need to read the K&L data
repeatedly for many ions. Requires KLGFB to also be input.
KLGFB Allows the KLGFB array from READ_KLGFB to be directly input to
KARZAS_XS thus avoiding the need to read the K&L data
repeatedly for many ions. Requires PE to also be input.
OUTPUT
The photoionization cross-section for the ionization of the outer
electron in units of mega-barns (Mb = 10^-18 cm^2) at the input
wavelengths. E.g., for Fe XIII (ground configuration
1s2.2s2.2p6.3s2.3p2) it is the cross-section for the removal of the
3p electron.
CALLS
READ_KLGFB
HISTORY
Ver.1, 24-Jul-2002, Peter Young
(See continuum/karzas_xs.pro)
PROJECT: CHIANTI
CHIANTI is an Atomic Database Package for Spectroscopic Diagnostics of
Astrophysical Plasmas. It is a collaborative project involving the Naval
Research Laboratory (USA), the University of Florence (Italy), the
University of Cambridge and the Rutherford Appleton Laboratory (UK).
NAME:
READ_FBLVL
PURPOSE:
To read files containing observed free-bound energy levels
CATEGORY:
Science.
CALLING SEQUENCE:
READ_FBLVL, File, L1, Conf, pqn, ll, spd, mult, Ecm, Ecmth, Ref
INPUTS:
File: the name of the file
e.g., !xuvtop+'/si/si_12/si_12.fblvl' for Si XII
OPTIONAL INPUTS:
None:
OUTPUTS: L1 - level index
CONF - configuration description
PQN - principal quantum number of electron being ionized
LL - L-value of ionized electron
SPD - 'S', 'P', etc to denote L value
MULT - multiplicity 2J+1
ECM - energy (cm^-1)
ECMTH - energy (cm^-1)
REF - reference
EXAMPLE:
> file = !xuvtop+'/si/si_12/si_12.elvlc'
> read_fblvl,file,l1,conf,pqn,ll,mult,ecm,ecmth,ref
>
HISTORY
Ver.1, 31-Jul-2002, Peter Young
(See continuum/read_fblvl.pro)
PROJECT: CHIANTI
CHIANTI is an atomic database package for the calculation of
continuum and emission line spectra from astrophysical plasmas. It is a
collaborative project involving the Naval Research Laboratory
(Washington DC, USA), the Arcetri Observatory (Firenze, Italy), and the
Cambridge University (United Kingdom).
NAME:
READ_KLGFB
PURPOSE:
to read CHIANTI files file containing the free-bound gaunt factors for n=1-6
from Karzas and Latter, 1961, ApJSS, 6, 167
CATEGORY:
science.
CALLING SEQUENCE:
READ_KLGFB,PE,GFB
INPUTS:
File: none
OUTPUTS:
PE: an array of natural log of 41 values of Photon Energy values relative to
the recombination edge i.e., at the recombination edge, PE=0.
GFB: an array containing the natural log of free-bound gaunt factors indexed
by energy, n and l
PROCEDURE: reads the file: '!xuvtop/continuum/klgfb.dat'
EXAMPLE:
> read_klgfb,pe,klgfb
the free-bound gaunt factor vs. energy for recombination to a hydrogen 2p state
(n=2 and l=1) is given by klgfb(*,n-1,l-1)
MODIFICATION HISTORY:
Written by: Ken Dere
July 2002: Version 1.0
Version 2, 8-Aug-02, GDZ
corrected a typo.
VERSION : V.2, 8-Aug-02
(See continuum/read_klgfb.pro)
PROJECT: CHIANTI
CHIANTI is an Atomic Database Package for Spectroscopic Diagnostics of
Astrophysical Plasmas. It is a collaborative project involving the Naval
Research Laboratory (USA), the University of Florence (Italy), the
University of Cambridge and the Rutherford Appleton Laboratory (UK).
NAME:
SUTHERLAND
PURPOSE:
Calculate the free-free continuum from an hot, low density plasma
Uses the free-free gaunt factor calculations of Sutherland, 1998,
MNRAS, 300, 321.
Note that Sutherland's Eq.(15) has units of erg/cm^3/s. Comparing
with Rybicki & Lightman's Eq.5.14(a) (in their book 'Radiative
Processes in Astrophysics'), suggests that Sutherland's units
should be erg/cm^3/s/sr/Hz. We are assuming the latter to be
correct in this routine.
When using the DEM_INT optional input, FREEFREE expects the
differential emission measure to have been derived from a product
of N_e*N_H rather than N_e*N_e. This can be important when dealing
with a regime (typically T < 10^4.5) where H and He are not fully
ionised.
To calculate the continuum spectrum the user should specify
whether the plasma structure is specified through a differential
emission measure (DEM) curve or through emission measure (EM)
values. For the former, the DEM values are specified through the
DEM_INT keyword and for the latter the EM values are specified
through the EM_INT keyword. If neither are specified EM_INT
values of unity are assumed.
If DEM_INT values are specifed they must be specified on a
temperature scale that has constant intervals in log10(T). EM_INT
values do not have this restriction.
Note that EM_INT must be used for isothermal plasmas.
CALLING SEQUENCE:
SUTHERLAND,temperature, wavelength, intensity
INPUTS:
T Temperature in degrees Kelvin, can be a 1 D
array. If the keyword /DT is set, then the
temperatures must be specified with a fixed spacing
in logT. E.g., TEMP=10.^[4.0,4.1,4.2].
WVL Wavelength in Angstroms. If /keV is set, then WVL is
assumed to contain energies in keV units.
OPTIONAL INPUTS:
None
KEYWORD PARAMETERS:
NO_SETUP: If the procedure setup_elements has already been called
then the keyword /no_setup should be set to avoid
repeating this step.
MIN_ABUND: If set, calculates the continuum only from those
elements which have an abundance greater than min_abund.
This can speed up the calculations. ;
DEM_INT This should be the same size as TEMP and contain
differential emission measure values. Specifying
DEM_INT places a restriction on TEMP (see above). Note
that DEM_INT values are multiplied by the factor dT
within the routine, whereas EM_INT values are not.
EM_INT This should be the same size as TEMP and contain
emission measure values. Note
that DEM_INT values are multiplied by the factor dT
within the routine, whereas EM_INT values are not. If
neither EM_INT or DEM_INT are specified, then a EM_INT
vector of same size as TEMP with values of unity is
assumed.
PHOTONS The output spectrum is given in photon units rather
than ergs.
SUMT The default is to output the intensity array as an array
of size (nwvl x nT). Setting this keyword performs a sum
over the temperatures to yield a vector of same size as
the input wavelengths, thus producing the complete
free-free spectrum.
KEV If set, then WVL is assumed to contain energies in keV
units rather than wavelengths in Angstroms.
OUTPUTS:
RAD Free-free continuum intensity in units
10^-40 erg cm^3 s^-1 str^-1 angstrom^-1 per unit
emission measure [ integral (N_e N_H dh) in cm^-5 ].
If T is given as a 1-D array, then RAD will be output
as a 2-D array, with one element for each temperature
and wavelength (but also see SUMT).
If the keyword /keV is set, then the units of INT will be
10^-40 erg cm^3/s/sr/keV
PROGRAMMING NOTES
The gaunt factors from Sutherland (MNRAS 300, 321, 1998) are a
function of uu and gg (see his Eq. 14). uu is a function of both
wavelength and T, while gg is a function of T only.
The gaunt factor (gff) is tabulated for values of uu and gg at
fixed intervals in log(uu) and log(gg). The log(uu) values go from
-4 to 4 in 0.2 steps; the log(gg) values go from -4 to 4 in 0.1 steps.
A particular (input) temperature and wavelength give rise to values
uu0 and gg0, the logs of which lie between -4 and 4. To derive the
corresponding gff0 value, I use the IDL routine BILINEAR.
BILINEAR requires, not uu and gg values as input, but indices.
E.g., the indices corresponding to the tabulated values of uu are
0 (=-4.0), 1 (=-3.8), 2 (=-3.6), etc. Thus, if log(uu0)=-3.76, then
i_uu0=1.20 is the index of uu0.
In order to make significant time-savings, I give BILINEAR all of
the wavelengths and temperatures in the same call for a particular
ion. To do this, I make my i_uu and i_gg values 2-D arrays of size
(nwvl x nT), and BILINEAR then produces a (nwvl x nT) array
containing the gff values.
A problem occurred if nT=1, as BILINEAR will turn the input uu and
gg vectors into 2-D arrays of size (nwvl x nwvl). If there are a
large number of wavelengths, this uses a lot of memory. To solve
this I make a 2 element temperature vector whose values are
identical, and then change this back to a 1 element vector after
BILINEAR has been called. See the parts of the code where I use
'tst1'.
COMMON BLOCKS:
common elements,abund,abund_ref,ioneq,ioneq_t,ioneq_ref
CALLS
READ_IP, READ_GFFGU, SETUP_ELEMENTS
EXAMPLES:
IDL> freefree,1.e+6,wvl,int
IDL> freefree,1.e+6,wvl,int,min_abund=3.e-5
IDL> freefree,1.e+6,wvl,int,/no_setup,min_abund=1.e-6
IDL> wvl=findgen(5001)/10. + 50.
IDL> temp=10.^(findgen(41)/10. +4.)
IDL> freefree,temp,wvl,int
MODIFICATION HISTORY:
Written by: Ken Dere
March 1999: Version 2.0
September 1999: Version 3.0
Ver.3.1, 11-Aug-00, Peter Young
Improved call to bilinear, allowing routine to solve for
all temperatures in one go. This makes the routine quicker,
and also lowers the memory usage of the routine when dealing
with many wavelengths.
Ver.3.2, 16-Aug-00, Peter Young
Corrected expression for 'gg', replacing ip with z^2.
Ver.3.3, 16-Oct-00, Peter Young
Now deals with dem_int correctly
Ver.3.4, 10-Oct-01, Ken Dere
Corrected for labelling errors in Sutherland's gffgu.dat file
No longer reads ionization potential file
Ver.3.5, 5-Dec-01, Peter Young
Corrected expression for gamma^2
Renamed routine sutherland.pro
Restructured code to make it run quicker.
Ver.3.6, 22-May-01, Peter Young
Re-instated the MIN_ABUND optional input.
Changed ioneq_t to ioneq_logt (GDZ).
Ver.3.7, 18-Aug-03, Peter Young
Activated /PHOTONS keyword
Ver.3.8, 5-Nov-03, Peter Young
Corrected bug found by Jim McTiernan when multiple temperatures
were input. The quantity 'newfactor' was not being calculated
correctly due to indexing problems.
Ver. 4, 9-Mar-2006, Peter Young
Added /keV keyword.
Ver.5, 10-Sep-2009, Peter Young
Routine now checks that temperatures are evenly spaced
when DEM_INT is specified. Also 0.1 dex spacing is no
longer assumed.
Ver.6, 6-Oct-2009, Peter Young
Introduced EM_INT keyword for specifying emission measure
values for isothermal plasmas.
(See continuum/sutherland.pro)
PROJECT: CHIANTI
CHIANTI is an Atomic Database Package for Spectroscopic Diagnostics of
Astrophysical Plasmas. It is a collaborative project involving the Naval
Research Laboratory (USA), the University of Florence (Italy), the
University of Cambridge and the Rutherford Appleton Laboratory (UK).
NAME:
TWO_PHOTON
PURPOSE:
Calculate the 2 photon continuum from a hot, low density plasma.
For the hydrogen isoelectronic sequence, A values
Parpia, F. A., and Johnson, W. R., 1982, Phys. Rev. A, 26, 1142.
and distribution function as a function of Z is taken from:
Goldman, S.P. and Drake, G.W.F., 1981, Phys Rev A, 24, 183
For the helium isoelectronic sequence, A values are from:
Drake, G.W.F., 1986, Phys. Rev. A, 34, 2871
and the distribution function as a function of Z is taken from:
Drake, G.W.F., Victor, G.A., Dalgarno, A., 1969, Phys.
Rev. A, 180, 25.
in this paper the distribution is only given up to Z=10 but
extrapolation to higher Z appears to be OK.
Note that, unlike the freefree and freebound routines, two_photon
requies the electron density to be specified. This is because there
is a call to pop_solver
To calculate the continuum spectrum the user should specify
whether the plasma structure is specified through a differential
emission measure (DEM) curve or through emission measure (EM)
values. For the former, the DEM values are specified through the
DEM_INT keyword and for the latter the EM values are specified
through the EM_INT keyword. If neither are specified EM_INT
values of unity are assumed.
If DEM_INT values are specifed they must be specified on a
temperature scale that has constant intervals in log10(T). EM_INT
values do not have this restriction.
Note that EM_INT must be used for isothermal plasmas.
CALLING SEQUENCE:
TWO_PHOTON,temperature, density, wavelength, intensity
INPUTS:
TEMP Temperature (in K). If the keyword /DT is set, then the
temperatures must be specified with a fixed spacing in
logT. E.g., TEMP=10.^[4.0,4.1,4.2].
WAVELENGTH Wavelengths in Angstroms. If /keV is set, then WVL is
assumed to contain energies in keV units.
OPTIONAL INPUTS:
DENSITY Electron density in cm^-3, can also be a 1D array
of the same size as Temperature. If there are several
temperatures specified but only one density, then
density is assumed the same at all temperatures.
If not specified, then default densities of 10^10
electrons/cm^3 are assumed.
DEM_INT This should be the same size as TEMP and contain
differential emission measure values. Specifying
DEM_INT places a restriction on TEMP (see above). Note
that DEM_INT values are multiplied by the factor dT
within the routine, whereas EM_INT values are not.
EM_INT This should be the same size as TEMP and contain
emission measure values. Note
that DEM_INT values are multiplied by the factor dT
within the routine, whereas EM_INT values are not. If
neither EM_INT or DEM_INT are specified, then a EM_INT
vector of same size as TEMP with values of unity is
assumed.
KEYWORD PARAMETERS:
NO_SETUP: If the procedure setup_elements has already been called
then the keyword /no_setup should be set to avoid
repeating this step
MIN_ABUND: If set, calculates the continuum only from those
elements which have an abundance greater than min_abund.
Can speed up the
calculations. For example:
abundance (H) = 1.
abundance (He) = 0.085
abundance (C) = 3.3e-4
abundance (Si) = 3.3e-5
abundance (Fe) = 3.9e-5
SUMT If several temperatures have been specified, then /sumt
will sum the emissivities over the different
temperatures, giving an output INTENSITY that has the
same size as WAVELENGTH.
PHOTONS If set the continuum emissivity is output in photon
units rather than erg units.
VERBOSE If set, information is printed to screen as the routine runs.
KEV If set, then WVL is assumed to contain energies in keV units
rather than wavelengths in Angstroms.
OUTPUTS:
RAD 2 photon continuum intensity in units of
10^-40 erg cm^3/s/sr/Angstrom per unit emission measure
( integral(N_H N_e dh) in cm^-5) if a DEM is not defined.
If DEM values are defined, it is assumed that they are given
as N_H N_e dh/dT. The units are
10^-40 erg/cm^2/s/sr/Angstrom
If T is given as a 1-D array, then the output will be a 2-D
array, with one element for each temperature and wavelength
(but also see SUMT).
If the keyword /keV is set, then the units of INT will be
10^-40 erg cm^3/s/sr/keV
CALLS
POP_SOLVER, SETUP_ION, SETUP_ELEMENTS, READ_MASTERLIST,
CONVERTNAME
COMMON BLOCKS:
ELEMENTS, ELVLC, WGFA, UPSILON, PROTON
EXAMPLE:
> two_photon,1.e+6,3.e+10,wvl,int
> two_photon,1.e+6,3.e+10,wvl,int,min_abund=3.e-5
> two_photon,1.e+6,3.e+10,wvl,int,/no_setup
PROGRAMMING NOTES
For He 2-photon decays, the distribution function is from Table II
of Drake et al. (1969), except that the values have been divided by
the A-value from Drake (1986).
MODIFICATION HISTORY:
Written by: Ken Dere
February 2001: Version 1.0
Ver.2, 19-Dec-2001, Peter Young
Now allows an array of temperatures.
Added /SUMT keyword.
Added DEM_INT= optional input.
Switched to using spl_init & spl_interp for the spline fits.
Corrected a small bug.
Ver.3, 20-Dec-2001, Peter Young
Corrected a bug related to density indices.
Ver.4, 23-Apr-2002, Peter Young
Added /photons keyword.
Ver.5, 28-May-2002, Peter Young
Corrected way in which DEM_INT is handled.
V. 6, 28-May-2002, Giulio Del Zanna (GDZ)
generalized directory concatenation to work for
Unix, Windows and VMS.
Corrected the description of the units and various
inaccuracies in the header.
V.7, 14-Aug-02, GDZ
Added keyword VERBOSE
V.8, 4-Oct-2003, GDZ
modified the input to POP_SOLVER, now it is dealt with an
input structure.
V.9, 8-Jun-2004, EL
modified the input to POP_SOLVER, now it includes ion/rec
V.10, 5-Jul-2005
corrected problems with the input structure for pop_solver
v.11 29-Jul-2005, GDZ
fixed bug, only define ioneq_ionrec when files are
there (it was failing with neutrals)
v.12, 31-Aug-2005, GDZ
missed one correct definition of ioneq_ionrec. Fixed it now.
v.13, 10-Mar-2006, Peter Young
added /keV keyword
v.14, 12-Jun-2009, Enrico Landi
Changed the definition of the temperature array for ion fractions
in the IONREC variable, now taken directly from the output of
READ_IONEQ.PRO
v.15, 11-Sep-2009, Peter Young
Routine now checks that temperatures are evenly spaced
when DEM_INT is specified. Also 0.1 dex spacing is no
longer assumed.
v.16, 6-Oct-2009, Peter Young
Introduced EM_INT keyword for specifying emission measure
values for isothermal plasmas.
v.17, 13-Apr-2010, Enrico Landi
Remove a bug in the definition of temperature for ion/rec effects when
dealing with He-like ions
VERSION : 17, 13-Apr-2010
(See continuum/two_photon.pro)
PROJECT: CHIANTI
CHIANTI is an Atomic Database Package for Spectroscopic Diagnostics of
Astrophysical Plasmas. It is a collaborative project involving the Naval
Research Laboratory (USA), the University of Florence (Italy), the
University of Cambridge and the Rutherford Appleton Laboratory (UK).
NAME
VERNER_XS()
EXPLANATION
Reads the Verner & Yakovlev (A&AS 109, 125, 1995) photoionization
cross-section data and generates the values of the cross-section at
the wavelengths WVL
INPUTS
IZ Atomic number of ion (e.g., 26 = Fe)
ION Spectroscopic number of ion (e.g., 13 = XIII)
WVL Wavelengths (in angstroms) for which cross-sections are
required (1-D array).
OPTIONAL INPUTS
DATA By default VERNER_XS reads the Verner & Yakovlev data file when
it is called. Through the DATA keyword the data array can be
sent to VERNER_XS instead.
OUTPUT
The photoionization cross-section for the ionization of the outer
electron in units of mega-barns (Mb = 10^-18 cm^2) at the input
wavelengths. E.g., for Fe XIII (ground configuration
1s2.2s2.2p6.3s2.3p2) it is the cross-section for the removal of the
3p electron.
HISTORY
Ver.1, 24-Jul-2002, Peter Young
(See continuum/verner_xs.pro)