Astronomy & Astrophysics manuscript no. aa6820 c© ESO 2013
February 12, 2013
A computer program for fast non-LTE analysis of
interstellar line spectra
With diagnostic plots to interpret observed line intensity ratios
F. F. S. van der Tak1,2, J. H. Black3, F. L. Schöier4, D. J. Jansen5, and E. F. van Dishoeck5
1 Netherlands Institute for Space Research (SRON), Landleven 12, 9747 AD Groningen, The Netherlands
e-mail: vdtak@sron.rug.nl
2 Max-Planck-Institut für Radioastronomie, Auf dem Hügel 69, 53121 Bonn, Germany
3 Onsala Space Observatory, Chalmers University of Technology, 43992 Onsala, Sweden
4 Stockholm Observatory, AlbaNova University Center, 10691 Stockholm, Sweden
5 Leiden University Observatory, P.O. Box 9513, 2300 RA Leiden, The Netherlands
Received 27 November 2006; accepted 27 March 2007
ABSTRACT
Aims. The large quantity and high quality of modern radio and infrared line observations require efficient modeling techniques
to infer physical and chemical parameters such as temperature, density, and molecular abundances.
Methods. We present a computer program to calculate the intensities of atomic and molecular lines produced in a uniform
medium, based on statistical equilibrium calculations involving collisional and radiative processes and including radiation from
background sources. Optical depth effects are treated with an escape probability method. The program is available on the World
Wide Web at http://www.sron.rug.nl/∼vdtak/radex/index.shtml . The program makes use of molecular data files
maintained in the Leiden Atomic and Molecular Database (LAMDA), which will continue to be improved and expanded.
Results. The performance of the program is compared with more approximate and with more sophisticated methods. An
Appendix provides diagnostic plots to estimate physical parameters from line intensity ratios of commonly observed molecules.
Conclusions. This program should form an important tool in analyzing observations from current and future radio and infrared
telescopes.
Key words. Radiative transfer – Methods: numerical – Radio lines – Infrared lines – Submillimeter
1. Introduction
Observations of spectral lines at radio, (sub)millimeter
and infrared wavelengths are a powerful tool to in-
vestigate the physical and chemical conditions in the
dilute gas of astronomical sources where thermody-
namic equilibrium is a poor approximation (e.g., Genzel
1991; Black 2000). To extract astrophysical parame-
ters from the data, the excitation and optical depth of
the lines need to be estimated, for which various meth-
ods may be used, depending on the available observa-
tions (Van Dishoeck & Hogerheijde 1999; Van der Tak
2005).
If only one or two lines of a molecule1 have been ob-
served, the excitation must be deduced from observa-
tions of other species or from theoretical considerations.
An example is the assumption that the excitation tem-
perature equals the kinetic temperature, a case known as
Local Thermodynamic Equilibrium (LTE) which holds
at high densities.
If many lines have been observed, a popular method
is the ‘rotation diagram’, also called ‘Boltzmann
plot’ or ‘population diagram’ (e.g., Blake et al. 1987;
Helmich et al. 1994; Goldsmith & Langer 1999). This
method describes the excitation by a single tempera-
ture, obtained by a fit to the line intensities as a func-
tion of upper level energy. Provided that beam sizes are
similar and optical depths are low, or that appropriate
corrections are made, this method yields estimates of
the excitation temperature and column density of the
molecule. The excitation temperature approaches the ki-
netic temperature in the high-density limit, but generally
depends on both kinetic temperature and volume den-
sity. Spectral line surveys are often analyzed with ro-
tation diagrams, although more advanced methods are
also used (Helmich & van Dishoeck 1997; Comito et al.
2005).
More sophisticated methods retain the assumption of
a local excitation, but solve for the balance of ex-
citation and de-excitation rates from and to a given
http://arxiv.org/abs/0704.0155v1
http://www.sron.rug.nl/~vdtak/radex/index.shtml
2 Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra
probability method and the Large Velocity Gradient
(LVG) method (Sobolev 1960; De Jong et al. 1975;
Goldreich & Scoville 1976). These ‘intermediate-level’
methods require knowledge of molecular collisional
data, whereas the previous ‘basic-level’ methods only
required spectroscopic and dipole moment information.
This extra requirement limits the use of these methods
to some extent, because collisional data do not exist
for all astrophysically relevant species. The advantage
is that column density, kinetic temperature and volume
density can be constrained, if accurate collision rates
are known. As with rotation diagrams, this method can
be used to compute synthetic spectra to be compared
with data with a χ2 statistic (Jansen 1995; Leurini et al.
2004).
The most advanced methods drop the local approxi-
mation and solve for the intensities (or the radiative
rates) as functions of depth into the cloud, as well as of
velocity. Such methods are usually of the Accelerated
Lambda Iteration (ALI) or Monte Carlo (MC) type, al-
though hybrids also exist. The performance and con-
vergence of such programs have recently been tested
by Van Zadelhoff et al. (2002). Using such programs
one can constrain temperature, density, and velocity
gradients within sources (e.g., Van der Tak et al. 1999;
Tafalla et al. 2002; Jakob et al. 2007), and, if enough
observations are available, even molecular abundance
profiles (e.g., Van der Tak et al. 2000a; Schöier et al.
2002; Maret et al. 2005), especially when coupled to
chemical networks (e.g., Doty et al. 2004; Evans et al.
2005; Goicoechea et al. 2006).
This paper presents the public version of a radiative
transfer code at the ‘intermediate’ level. The assump-
tion of a homogeneous medium limits the number of
free parameters and makes the program a useful tool
in rapidly analyzing a large set of observational data,
in order to provide constraints on physical conditions,
such as density and kinetic temperature (Jansen 1995).
The program can be used for any molecule for which
collisional rate coefficients are available. The input for-
mat for spectroscopic and collisional data is that of the
LAMDA database (Schöier et al. 2005)2 where an on-
line calculator for molecular line intensities3, based on
our program, can also be found4.
The paper is set up as follows. Section 2 describes the
radiative transfer formalism and introduces our notation
of the key quantities. Section 3 describes the formalism
which the program actually uses, and discusses its im-
plementation. Section 4 compares the results of the pro-
gram to those of other programs. The paper concludes
in § 5 with suggested future directions of astrophysical
radiative transfer modeling.
2 http://www.strw.leidenuniv.nl/∼moldata
3 In this paper, the ‘strength’ of a line is an intrinsic quantity
2. Radiative transfer and molecular
excitation
This section summarizes the formalism to analyze
molecular line observations which our program adopts.
For more detailed discussions of radiative transfer see,
e.g., Cannon (1985) or Rybicki & Lightman (1979).
2.1. Basic formalism
Describing the transfer of radiation requires a quantity
which is conserved along its path as long as no local
absorption or emission processes take place, and which
includes the direction of travel. The quantity that sat-
isfies this requirement is the specific intensity Iν, de-
fined as the amount of energy passing through a surface
normal to the path, per unit time, surface, bandwidth
(measured here in frequency units), and solid angle. The
transfer equation for radiation propagating a distance ds
can then be written as
= jν − αν Iν, (1)
where jν and αν are the local emission and extinction
coefficients, respectively. The two terms on the right-
hand side may be combined into the source function,
defined by
S ν ≡
. (2)
Writing the transport equation in its integral form and
defining the optical depth, dτν ≡ αν ds, measured along
the ray5 one arrives at
Iν = Iν (0)e
−τν +
S ν(τ
−(τν−τ′ν) dτ′ν, (3)
where Iν is the radiation emerging from the medium and
Iν(0) is the ‘background’ radiation entering the medium.
The above equations hold both for continuum radiation,
which is emitted over a large bandwidth, and for spectral
lines, which arise when the local emission and absorp-
tion properties change drastically over a very small fre-
quency interval, due to the presence of molecules. From
this point the discussion will focus on bound-bound
transitions within a multi-level molecule consisting of
N levels with spontaneous downward rates Aul, Einstein
coefficients for stimulated transitions Bul and Blu, and
collisional rates Cul and Clu, between upper levels u and
lower levels ℓ.
The rate of collision is equal to
Cul = ncolγul, (4)
where ncol is the number density of the collision partner
(in cm−3) and γul is the downward collisional rate coeffi-
cient (in cm3 s−1). The rate coefficient is the Maxwellian
average of the collision cross section, σ,
γul =
8kTkin
)−1/2 (
kTkin
σEe−E/kTkin dE, (5)
http://www.strw.leidenuniv.nl/~moldata
http://www.sron.rug.nl/~vdtak/radex.php
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra 3
where E is the collision energy, k is the Boltzmann con-
stant, Tkin is the kinetic temperature, and µ is the re-
duced mass of the system. The upward rates are ob-
tained through detailed balance
γlu = γul
e−hν/kTkin , (6)
where gi is the statistical weight of level i.
The local emission in transition u→l with laboratory
frequency νul, can be expressed as
nu Aul φν, (7)
where nu is the number density of molecules in level u
and φν is the frequency-dependent line emission profile.
The absorption coefficient reads
(nl Blu ϕν − nu Bul χν) , (8)
where φν and χν are the line profiles for absorption and
stimulated emission (counted as negative extinction),
respectively.
From here on we assume complete angular and fre-
quency redistribution of the emitted photons, so that
φν=ϕν=χν, which is strictly only valid when collisional
excitation dominates. This assumption allows the source
function to be written as
S νul =
nuAul
nlBlu − nuBul
, (9)
where we have used the Einstein relations. It is com-
mon to introduce an excitation temperature Tex defined
through the Boltzmann equation
exp [−(Eu − El)/kTex] , (10)
where Ei is the energy of level i, such that S νul =
Bν(Tex), the specific intensity of a blackbody radiating
at Tex.
In the interstellar medium, the dominant line broadening
mechanism is Doppler broadening. Except in very cold
and dark cloud cores, observed line widths are much
larger than expected from the kinetic temperature: this
effect is commonly ascribed to random macroscopic gas
motions or ‘turbulence’. The result is a Gaussian line
profile
ν − νul − v · n
, (11)
where νD is the Doppler width, v is the velocity vector
of the moving gas at the position of the scattering, n is
a unit vector in the direction of the propagating beam of
radiation, and c is the speed of light. The Doppler width
is the 1/e half-width of the profile, equal to ∆V/2
where ∆V is its full width at half-maximum.
If the level populations ni are known, the radiative trans-
fer equation can be solved exactly. In particular, under
circumstellar media, the density is too low to attain LTE,
but statistical equilibrium (SE) can often be assumed:
= 0 =
n jP ji − ni
Pi j = Fi − niDi , (12)
where Pi j, the destruction rate coefficient of level i, and
its formation rate coefficient P ji are given by
Pi j =
Ai j + Bi j J̄ν + Ci j (i > j)
Bi j J̄ν +Ci j (i < j).
In Eq. 13,
Bi j J̄ν = Bi j
Jν φ(ν) dν (14)
is the number of induced radiative (de-)excitations from
state i to state j per second per particle in state i, and
Iν dΩ (15)
is the specific intensity Iν integrated over solid angle dΩ
and averaged over all directions. The SE equations thus
include the effects of non-local radiation.
This discussion assumes that the state-specific rates of
formationFi [cm3 s−1] and destructionDi [s−1] are zero
to ensure that the radiative transfer is solved indepen-
dently of assumptions about chemical processes. In gen-
eral, formation and destruction processes should be in-
cluded explicitly to be able to deal with situations in
which the chemical time scales are very short or the
radiative lifetimes very long. For example, the forma-
tion temperature (in Fi) affects the rotational excitation
of C3 (Roueff et al. 2002) and the vibrational excitation
of H2 (Black & van Dishoeck 1987; Burton et al. 1990;
Takahashi & Uehara 2001), systems for which line ra-
diation only occurs as slow electric quadrupole transi-
tions. The rotational excitation of reactive ions like CO+
(Fuente et al. 2000; Black 1998) is also sensitive to Fi
and Di because the rates of reactions with H and H2
rival the inelastic collision excitation rates. Similar con-
siderations apply to the excitation of H+3 in the Sgr A
region close to the Galactic Center (Van der Tak 2006),
where electron impact excitation competes with disso-
ciative recombination.
2.2. Molecular line cooling
Once the radiative transfer problem has been solved and
the level populations are known, the cooling (or heating)
from molecular line emission can be estimated. Since
the level populations contain all the information of the
radiative transfer, a general expression for the cooling is
obtained from considering all possible collisional tran-
sitions
(nlγlu − nuγul)hνul, (16)
4 Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra
are in detailed balance at the kinetic temperature; there-
fore it is possible for net heating to occur (Λ < 0) in
cases where the crucial level populations have Tex>Tkin,
owing to strong radiative excitation in a hot external ra-
diation field.
2.3. Escape probability
The difficulty in solving radiative transfer problems is
the interdependence of the molecular level populations
and the local radiation field, requiring iterative solution
methods. In particular, for inhomogeneous or geometri-
cally complex objects, extensive calculations with many
grid points are required. However, if only the global
properties of an interstellar cloud are of interest, the cal-
culation can be greatly simplified through the introduc-
tion of a geometrically averaged escape probability β,
the probability that a photon will escape the medium
from where it was created. This probability depends
only on the optical depth τ and is related to the inten-
sity within the medium, ignoring background radiation
and any local continuum, through
Jνul = S νul (1 − β). (17)
Several authors have developed detailed relations be-
tween β and τ for specific geometrical assumptions.
Our program offers the user a choice of three such ex-
pressions. The first is the expression derived for an ex-
panding spherical shell, the so-called Sobolev or large
velocity gradient (LVG) approximation (Sobolev 1960;
Castor 1970; Elitzur 1992, p. 42-44). This method is
also widely applied for moderate velocity gradients, to
mimic turbulent motions. Our program uses the formula
by Mihalas (1978) and De Jong et al. (1980) for this ge-
ometry:
βLVG =
dτ′ =
1 − e−τ
. (18)
Second, in the case of a static, spherically symmet-
ric and homogeneous medium the escape probability is
(Osterbrock & Ferland 2006, Appendix 2)6
βsphere =
. (19)
Third, for a plane-parallel ‘slab’ geometry, applicable
for instance to shocks,
βslab =
1 − e−3τ
is derived (De Jong et al. 1975). Figure 1 plots the be-
haviour of β as a function of τ for these three cases; for
more detailed comparisons see Stutzki & Winnewisser
(1985) and Ossenkopf et al. (2001). Users of our pro-
gram can select either expression for their calculations.
The on-line version of the program uses the formula for
the uniform sphere, Eq. (19).
 0  1  2  3  4  5  6  7  8  9  10
Optical depth τ
Sphere
Fig. 1. Escape probability β as a function of optical
depth τ for three different geometries: uniform sphere
(solid line), expanding sphere (dotted line) and plane-
parallel slab (dashed line).
3. The program
RADEX is a non-LTE radiative transfer code, written
originally by J. H. Black, that uses the escape proba-
bility formulation assuming an isothermal and homo-
geneous medium without large-scale velocity fields.
With the current increase of observational possibilities
in mind, we have developed a version of this program
which is suitable for public use. A guide for using the
code in practice is provided in Appendix A and on-line7;
Appendix B describes the adopted coding style. This
section focuses on the implementation of the formalism
of § 2 in the program.
3.1. Basic capabilities
For a homogeneous medium with no global velocity
field, the optical depth at line centre can be expressed
using Eqs. (2, 7, 9, 11), as
AulNmol
1.064∆V
, (21)
where Nmol is the total column density,∆V the full width
at half-maximum of the line profile in velocity units,
and xi the fractional population of level i. The formal-
ism is analogous to the LVG method, with the global
n/(dV/dR) replaced by the local N/∆V , as in microtur-
bulent codes (Leung & Liszt 1976). The program itera-
tively solves the statistical equilibrium equations start-
ing from optically thin statistical equilibrium (§ 3.4) for
the initial level populations.
The program can handle up to seven collision partners
simultaneously. In dense molecular clouds, H2 is the
main collision partner for most species, but in some
cases, separate cross sections may exist for collisions
with the ortho and para forms of H2, and electron col-
lisions may be important for ionic species. In diffuse
molecular clouds and PDRs, excitation by atomic H be-
comes important, particularly for fine structure lines,
while for comets, H2O is the main collision partner. We
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra 5
refer to Flower (1989) for the basic theory of molecular
collisions, and to Dubernet (2005) for an update of the
latest results.
The output of the program is the background-subtracted
line intensity in units of the equivalent radiation temper-
ature in the Rayleigh-Jeans limit. The background sub-
traction follows traditional cm- and mm-wave spectro-
scopic observations where the differences between on-
source and off-source measurements are recorded, such
Iemν − I
. (22)
The radiation peak temperature TR can be directly com-
pared to the observed antenna temperature corrected
for the optical efficiency of the telescope. However, it
should be emphasized that RADEX contains no informa-
tion about the geometry or length scale and that it is
assumed that the source fills the antenna beam. If the
source is expected to be smaller than the observational
beam, computed line fluxes must be corrected before
comparing to observed fluxes.
In other types of observations, the continuum may
not be subtracted from the data. In (sub-)millimeter
and THz observations, for example with ESA’s fu-
ture Herschel space observatory, the dust continuum
of many sources will be much stronger than any in-
strumental error, and baseline subtraction may not be
needed. The same is true for interferometer data, where
the instrumental passband is well characterized.
3.2. Background radiation field
The average Galactic background (interstellar radia-
tion field, ISRF) adopted in RADEX consists of several
components. The main contribution is the cosmic mi-
crowave background (CMB) whose absolute tempera-
ture is taken to be TCBR = 2.725±0.001K based on the
full COBE data set as analyzed by Fixsen & Mather
(2002). This model of the microwave background rep-
resents the broadband continuum only and does not in-
clude the strong emission lines, several of which contain
significant power in the far-infrared and (sub-)millime-
ter part of the spectrum (see, e.g., Fixsen et al. 1999).
The ultraviolet/visible/near-infrared part of the spec-
trum is based on the model of average Galactic starlight
in the solar neighborhood of Mathis et al. (1983). The
far-infrared and (sub-)millimeter part of the spectrum is
based on the single-temperature fit to the Galactic ther-
mal dust emission of Wright et al. (1991). At frequen-
cies below 10 cm−1 (30 GHz), there is a background
contribution from non-thermal radiation in the Galaxy.
A tabulation of this spectrum in ASCII format is avail-
able on-line8, and a graphical representation is shown in
Black (1994).
One subtle aspect of the calculation is the distinction
between the background seen by the observer and the
background seen by the molecules. The continuum con-
tribution to the rate equations may be composed of (1)
an external component which arises outside the emitting
region and (2) an internal continuum that arises within
the emitting region. The CMB and ISRF are examples
of external continuum components; dust emission from
the line-emitting region is an example of an internal
continuum. While an external continuum always fills
the entire sky, an internal continuum may only fill a frac-
tion of it, for example in the case of a circumstellar disk.
With this distinction in mind, the internal intensity be-
comes
Jintν = β[Bν(TCBR)+ηI
ν ]+(1−β)[Bν(Tex)+θ(1−η)I
where Iuserν is the continuous spectrum defined by the
user. The factor η, is the fraction of local continuum
which arises outside the line emitting region, and the
factor θ is the fraction of local sky filled by the internal
continuum.
3.3. Chemical formation and destruction rates
The equations of statistical equilibrium (12) include
source and sink terms. By default, RADEX sets the
destruction rates equal to the same small value,
Di ≡ D = 10−15 s−1, appropriate for cosmic-ray
ionization plus cosmic-ray induced photodissociation
(Prasad & Tarafdar 1983; Gredel et al. 1989). The cor-
responding formation rates are
Fi = 10−24ntotalgi exp(−Ei/kTform)/Q(Tform) (24)
where ntotal is the sum of the densities of all collision
partners, Tform is a formation temperature (default value
300 K), and
Q(T ) =
gi exp(−Ei/kTform) (25)
is the partition function. These assumptions imply a
nominal fractional abundance of every molecule
ntotal
ntotalD
= 10−9 . (26)
The value of the nominal abundance is inconsequential
because the results in RADEX depend on Nmol/∆V , but
not on the fractional abundance. For most molecules
currently in the associated database (LAMDA) and for
the most commonly encountered interstellar conditions,
these choices will not affect the observable excitation.
The formation and destruction rates are computed in a
subroutine that can be modified by the user to provide
a more realistic description of chemical processes. For
example, users may treat the combined ortho/para forms
of molecules by introducing a realistic Tform, especially
in cases where no o/p interchange processes is likely
to be effective. Other cases of potential interest include
the photodissociation of large molecules into smaller
molecules, or the evaporation of icy grain mantles into
the gas phase. Our formulation in terms of a volume
rate of formation is chosen to be independent of the de-
tails of the formation process. In general, formation and
destruction processes are important for molecules that
http://www.oso.chalmers.se/~jblack/RESEARCH/
6 Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra
3.4. Calculation
The input parameters of RADEX and its output are de-
scribed in Appendix A. Calculations with RADEX pro-
ceed as follows. A first guess of the populations of the
molecular energy levels is produced by solving statisti-
cal equilibrium in the optically thin case. The only ra-
diation taken into account is the unshielded background
radiation field; internally produced radiation is not yet
available. The solution for the level populations allows
calculation of the optical depths of all the lines, which
are then used to re-calculate the molecular excitation.
The new calculation treats the background radiation in
the same manner as the internally produced radiation.
The program iteratively finds a consistent solution for
the level populations and the radiation field. When the
optical depths of the lines with τ > 10−2 are stable from
one iteration to the next to a given tolerance (default
10−6), the program writes output and stops.
3.5. Results
There are several ways in which RADEX can be used
to analyze molecular line observations. In most of
these applications, the modeled quantity is the velocity-
integrated line intensity, as the excitation is assumed
to be independent of velocity. As a consequence, self-
absorbed lines cannot be modeled satisfactorily with
RADEX. In the simplest case, the temperature and den-
sity are known from other observations and only the col-
umn density of the molecule under consideration needs
to be varied to get the best agreement with the observed
line intensity. If the H2 column density is known from
other observations, for example from an optically thin
CO isotopic line, the ratio of the two column densi-
ties gives the molecular abundance, averaged over the
source. The RADEX distribution contains a Python script
to automate this procedure which is further described in
Appendix D.
Another often-used application of RADEX is to deter-
mine temperatures and densities from the observed in-
tensity ratios of lines of the same molecule. If the abun-
dance of the molecule is constant throughout the source,
the ratios should give source-averaged physical con-
ditions independent of the specific chemistry of the
molecule. Appendix C presents illustrative plots of line
ratios for commonly observed molecules and lines in
the optically thin case. For higher optical depths, the
qualitative trends remain the same but there are quanti-
tative differences. RADEX can readily be used to generate
similar plots for moderately thick cases. Again, Python
scripts are made available to automate this procedure
(Appendix D).
To illustrate the use of RADEX on actual observations, we
take the observations of the HCO+ 1–0 and 3–2 lines to-
ward a relatively simple source, the photon-dominated
region IC 63 (Jansen et al. 1994). The observed 1–0/3–2
ratio corrected for beam dilution is 5.5±1.5. The kinetic
temperature of the source is constrained from CO ob-
servations to be ∼50 K. Fig. C.3 shows that the inferred
Fig. 2. Comparison of the predicted line strengths for
the 10 lowest rotational transitions of CO for a homo-
geneous isothermal sphere, with nH2 = 10
5 cm−3 and
Tkin = 50 K, using different methods. Upper panel: The
total optical depth through the sphere at line centre τ
and excitation temperature Tex as a function of the up-
per rotational level J involved in the transition. Middle
panel: The radiation temperature TR obtained for each
transition using RADEX with different prescriptions of
the escape probability β and compared with the result
from the Monte-Carlo code (MCC) of Schöier (2000).
Also shown are the results for optically thin emission in
LTE. Lower panel: TR obtained from RADEX compared
with the results from MCC, δTR = T
R − TR.
inferred column density from the absolute intensities is
8 × 1012 cm−2, which, together with the overall H2 col-
umn density of 5×1021 cm−2, gives an HCO+ abundance
with respect to H2 of 1.6 × 10−9.
A slightly more complicated situation arises for the
Orion Bar PDR (Hogerheijde et al. 1995). For this
source, both HCO+ 1–0, 3–2 and 4–3 lines have been
observed. The 1–0/3–2 ratio gives an order of magni-
tude lower density than the 3–2/4–3 ratio. This differ-
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra 7
density tracers such as CO are usually much less sensi-
tive to density. One possible solution is a clumpy PDR
model in which the 1–0 line is mostly produced in the
low-density interclump gas containing 90% of the mate-
rial and the 4–3 line in the high-density clumps. Within
this clumpy model, a single column density fits all three
lines and an accurate abundance can be derived. Note
that this technique of adding the results of two models
is only applicable at low optical depths.
If only two lines of a molecule have been observed,
the line ratio can be used as indicator of tempera-
ture or density, depending on molecule and transition
(Appendix C). A single line ratio is never enough
to constrain both temperature and density, though.
For multi-line observations, a comparison of data and
models in terms of χ2 is preferred. See for example
Van der Tak et al. (2000b), Schöier et al. (2002), and
Leurini et al. (2004) for details of such calculations.
3.6. Limitations of the program
The current version of the program does not in-
clude a contribution from continuous (dust or free-free)
opacity to the escape probability, as for example in
Takahashi et al. (1983). Continuum radiation from dust
is generally negligible at long wavelengths ( >∼ 1 mm)
but becomes important for regions with very high col-
umn densities (such as protoplanetary disks) and at
far-infrared and shorter wavelengths ( <∼ 100µm). Free-
free radiation may become important for the calcula-
tion of atomic fine structure lines from H II regions;
other programs such as CLOUDY (Ferland 2003)9 may
be more suitable for this purpose. The absence of con-
tinuous opacity limits the applicability of the program
particularly in situations where infrared pumping is im-
portant, either directly through rotational transitions or
via vibrational transitions (Carroll & Goldsmith 1981;
Hauschildt et al. 1993).
Another limitation of the program is that only one
molecule is treated at a time, so that the effects of line
overlap are not taken into account. Such overlaps may
occur both at radio and at infrared wavelengths (e.g.,
Expósito et al. 2006). In special cases, overlap between
lines of the same molecule may influence their excita-
tion, for example the hyperfine components of HCN or
+ (Daniel et al. 2006).
For certain molecules under certain physical conditions
(especially low density and/or strong radiation field),
population inversions occur, which cause negative op-
tical depth and hence nonlinear amplification of the
incoming radiation (Elitzur 1992). This phenomenon,
known as ‘maser’ action, requires non-local treatment
of the radiative transfer, in particular a fine sampling of
directions, for which RADEX is not set up. Generally, the
escape probability approximation is justified until the
masers saturate, which occurs at τ ≈ −1. In practice,
the computed intensities of lines with τ <∼ − 0.1 are not
as accurate as those of other lines, and the intensities of
lines with τ <∼ − 1 should be disregarded altogether. If
Fig. 3. Excitation temperature of the HCO+ 1–0 tran-
sition as a function of radius for the model cloud of
§ 4.1.2, calculated with RADEX assuming static spher-
ical, expanding spherical, and slab geometry (dashed
/ dash-dotted / dotted lines), with a multi-zone es-
cape probability program (long/short dashes) and with
a Monte Carlo code (solid lines). The panels are for dif-
ferent column densities, hence optical depths. Note the
different vertical scales.
of non-maser lines may also be affected. While special-
ized programs should be used to calculate the intensities
of maser lines (e.g., Spaans & van Langevelde 1992;
Gray & Field 1995; Yates et al. 1997), RADEX may well
be used to predict which lines of a molecule may display
http://www.nublado.org/
8 Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra
which are pumped by infrared radiation (Leurini et al.
2004).
4. Comparison with other methods
This section shows a comparison of RADEX with other
programs, first for the case of constant physical con-
ditions (§ 4.1) and second for variable conditions
(§ 4.2). Comparison is with the analytical rotation di-
agram method and with Monte Carlo methods, which
have been benchmarked to high accuracy, both for
the case of HCO+ (Van Zadelhoff et al. 2002)10 and of
H2O (Van der Tak et al. 2005)
11. Throughout this sec-
tion, molecular data have been taken from the LAMDA
database (Schöier et al. 2005).
4.1. Homogeneous models
4.1.1. The case of CO
To test the RADEX code, we have compared its output
both to an optically thin LTE analysis (rotation diagram
method) and a full radiative transfer analysis using a
Monte-Carlo method (Schöier 2000). The test problem
consists of a spherically symmetric cloud with a con-
stant density, n(H2), of 1 × 105 cm−3 within a radius of
100 AU. In this example only the CO emission is treated
using a fractional abundance of 1 × 10−4 relative to H2
yielding a central CO column density of NCO = 3 ×
1016 cm−2 and an average value of N=2×1016 cm−2. The
kinetic temperature is set to 50 K, the background tem-
perature to 2.73 K, and the line width to ∆V=1.0 km s−1.
Fig. 2 presents the results of the calculations for the
ten lowest rotational transitions. The excitation temper-
atures of the lines vary from being close to thermalized
for transitions involving low J-levels, to sub-thermally
excited for the higher-lying lines. The optical depth in
the lines is moderate (∼ 1 − 2) to low. It is seen that the
expressions of the escape probability for the uniform
sphere and the expanding sphere give almost identical
solutions which are close to that obtained from the full
radiative transfer (MCC in Fig. 2). The slab geometry
gives slightly higher intensities, in particular for high-
lying lines. The optically thin approximation, where the
gas is assumed to be in LTE at 50 K, produces much
larger discrepancies, up to a factor of∼ 2, and only gives
the correct intensity for the J = 1 → 0 line, where the
LTE conditions are met.
4.1.2. The case of HCO+
To further verify the performance of the RADEX
program, we have compared its results to that
of another program that does not use the local
approximation: the Monte Carlo program RATRAN
(Hogerheijde & van der Tak 2000). We also compare
the results to those from the multi-zone escape proba-
bility program by Poelman & Spaans (2005). The test
case is a cloud with n(H2) = 1×104 cm−3, Tkin=10 K,
Tbg=2.73 K, and a line width of ∆V=1.0 km s
−1, equiv-
alent to bD=0.6 km s
−1. The pure rotational emission
spectrum of HCO+ was calculated for column densi-
ties of 1012, 1013, 1014 and 1015 cm−2, which for RADEX
were given directly as input parameters. For the multi-
zone programs, a cloud radius of 1018 cm was specified
along with abundances of 10−10–10−7, distributed over
50 cells.
Figure 3 shows the calculated excitation temperature
of the HCO+ 1–0 transition as a function of radius
for these physical conditions. For N(HCO+) <∼ 10
12 cm−2,
the excitation is independent of radius and the calcu-
lations for the various geometries agree to ≈10%. The
dependence of the excitation on radius and on geom-
etry increases with increasing column density, and for
N(HCO+) >∼ 10
15 cm−2, the curvature of the Tex distri-
bution becomes too large to ignore. The corresponding
line optical depth is ≈100, with ≈20% spread between
the various estimates (Fig. 4). The curvature arises be-
cause at the cloud center, photon trapping thermalizes
the excitation, while at the edge, the emission can es-
cape the cloud (Bernes 1979). We do not recommend to
use RADEX at line optical depths >∼ 100, because the cal-
culated excitation temperature may not be representa-
tive of the emitting region. However, even if some lines
are highly optically thick, RADEX may well be used to
analyze other lines which are optically thin. For exam-
ple for H2O, the ground state lines often have τ ∼ 1000,
but RADEX is well capable of computing intensities for
higher-lying transitions which are not as optically thick.
At low optical depth, variations in Tex translate directly
into changes in emergent line intensity. Thus, differ-
ences as large as 20% in calculated line flux can arise
depending on the choice of escape probability descrip-
tion, even for moderately thick cases. At high optical
depth, the direct connection between Tex and line flux
is lost because of the dependence on the adopted veloc-
ity field. The assumption in the program that the opti-
cal depth is independent of velocity breaks down in this
case. In this limit, the peak line temperature TR gives
the value of Tex at the τ = 1 surface of the cloud in this
specific transition.
The results shown in this section do not translate easily
to other HCO+ lines such as J=3→2, because the ex-
citation is governed by several competing effects. The
optical depth of the J=3→2 line may be higher or lower
than that of the J=1→0 line, depending on temperature
and density. Observers are encouraged to use RADEX to
study the excitation of their lines as a function of these
parameters, and also consider geometric variations.
4.2. Observations of a Young Stellar Object
To compare a typical RADEX analysis with other meth-
ods for a situation which varying physical conditions,
we choose a molecule for which many lines can be
observed: para-formaldehyde, p-H2CO. Figure 5 shows
observations of p-H2CO (sub-)millimeter emission lines
http://www.strw.leidenuniv.nl/astrochem/
http://www.sron.rug.nl/~vdtak/H2O/
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra 9
Fig. 4. Optical depth of the HCO+ 1–0 transition for
the model of § 4.1.2, calculated with RADEX assuming
static spherical, expanding spherical, and slab geometry
(filled triangles/filled squares/open triangles), with the
multi-zone escape probability code (open circles) and
with the Monte Carlo code (filled circles).
Fig. 5. Line strengths of p-H2CO observed towards
the embedded low-mass protostar IRAS 16293–2422
(squares with error bars), modeled assuming LTE (solid
line), using RADEX (triangles), and using a Monte Carlo
program assuming a constant abundance (solid circles)
and an abundance varying with radius (open circles).
Fig. 6. Distributions of the χ2 parameter correspond-
ing to the models in Figure 5. The RADEX results are
for n(H2) = 10
6 cm−3 as found by Van Dishoeck et al.
(1995).
Van Dishoeck et al. (1995). The data are analyzed us-
ing three methods: assuming LTE (with a rotation dia-
gram), assuming SE (using RADEX), and using a Monte
Carlo program. The free parameters for the LTE fit are
the excitation temperature Tex and the column density
Figure 6 shows the distributions of the χ2 parameter for
the LTE and SE fits, calculated in the standard way (see,
e.g., Van der Tak et al. 2000b; Schöier et al. 2002) as-
suming a 20% uncertainty for all observed points except
the line with the highest upper level energy where a 30%
uncertainty was used. As seen from the figure, the non-
LTE method gives a better fit to the data, as quantified
by the lower minimum χ2 value. This result is not nec-
essarily surprising given that more free parameters are
available. A more important difference is that the esti-
mates of temperature and column density between the
two methods are substantially different, in particular the
temperature (50 vs 150 K). Since the non-LTE method
involves fewer assumptions about the physical state of
the cloud, its results are to be preferred.
These results illustrate that rotation diagrams may give
misleading results when determining physical proper-
ties of interstellar gas clouds (cf. Johnstone et al. 2003
for the case of CH3OH). Figure 5 also demonstrates that
temperatures and column densities derived from rota-
tion diagrams tend to depend on which lines happen to
have been observed (cf. the HCO+ case in § 4.1.2). From
other data, IRAS 16293–2422 is actually known to have
a gradient in temperature and density throughout its en-
velope, which cannot be modelled properly with either
technique. For such situation, a full Monte Carlo radia-
tive transfer method is needed in which both the physi-
cal conditions and the abundances can vary with radius
(Fig. 5, circles). Nevertheless, the column densities and
abundances inferred with RADEX using the physical con-
ditions inferred from the line ratios differ by only a fac-
tor of a few from those found with the more sophis-
ticated analysis, at least for the particular zone of the
source to which those conditions apply (Schöier et al.
2002).
5. Conclusions
We have presented a computer program to analyze spec-
tral line observations at radio and infrared wavelengths,
based on the escape probability approximation. The pro-
gram can be used for any molecule for which collisional
data exist; such input data are available in the required
format from the LAMDA database. The program can be
used for optical depths from ≈–0.1 to ≈100.
The limited number of free parameters makes RADEX
very useful to rapidly analyze large datasets. As an
example, observed line intensity ratios may be com-
pared with the plots in Appendix C to estimate density
and kinetic temperature. Ratios of other lines and other
molecules may be easily computed using the Python
scripts included in the RADEX distribution. The program
may also be used to create synthetic spectra. This capa-
bility will be important to model the THz line surveys
from the HIFI instrument onboard the Herschel space
observatory.
In the future, we plan to incorporate a multi-zone es-
cape probability formalism (Poelman & Spaans 2005;
Elitzur & Asensio Ramos 2006) which will enable
10 Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra
depths, the calculation may also start from LTE condi-
tions rather than from optically thin statistical equilib-
rium. Robust convergence may be achieved by starting
from either initial condition and requiring the two an-
swers to be equal. For the modeling of crowded spec-
tra, the effects of line overlap will also need to be con-
sidered, for instance in the ‘all or nothing’ approach
(Cesaroni & Walmsley 1991). Such spectra will be rou-
tinely observed with the superb resolution and sensitiv-
ity of ALMA.
Our program is free for anybody to use for science, pro-
vided that appropriate reference is made to this paper.
For any other purpose such as to incorporate the pro-
gram into other packages which may be distributed to
the public, prior agreement with the authors is needed.
Acknowledgements. The authors wish to thank Huib Jan van
Langevelde for his efforts in documenting RADEX, and Erik
Deul for computing support at Leiden Observatory. JHB and
FLS acknowledge the Swedish Research Council for financial
support. FvdT and EvD thank the Netherlands Organization
for Scientific Research (NWO) and the Netherlands Research
School for Astronomy (NOVA). Finally we thank Volker
Ossenkopf, Marco Spaans, and an anonymous referee for
helpful comments on the manuscript.
References
Bernes, C. 1979, A&A, 73, 67
Black, J. H. 1994, in ASP Conf. Ser. 58: The First
Symposium on the Infrared Cirrus and Diffuse
Interstellar Clouds, ed. R. M. Cutri & W. B. Latter,
Black, J. H. 1998, in Chemistry and Physics of
Molecules and Grains in Space. Faraday Discussions
No. 109, 257
Black, J. H. 2000, in IAU Symposium 197 –
Astrochemistry: From Molecular Clouds to Planetary
Systems, ed. Y. C. Minh & E. F. van Dishoeck, 81
Black, J. H. & van Dishoeck, E. F. 1987, ApJ, 322, 412
Blake, G. A., Sutton, E. C., Masson, C. R., & Phillips,
T. G. 1987, ApJ, 315, 621
Blake, G. A., van Dishoek, E. F., Jansen, D. J.,
Groesbeck, T. D., & Mundy, L. G. 1994, ApJ, 428,
Burton, M. G., Hollenbach, D. J., & Tielens, A. G. G. M.
1990, ApJ, 365, 620
Cannon, C. J. 1985, The transfer of spectral line radia-
tion (Cambridge: University Press)
Carroll, T. J. & Goldsmith, P. F. 1981, ApJ, 245, 891
Castor, J. I. 1970, MNRAS, 149, 111
Cesaroni, R. & Walmsley, C. M. 1991, A&A, 241, 537
Comito, C., Schilke, P., Phillips, T. G., et al. 2005, ApJS,
156, 127
Daniel, F., Cernicharo, J., & Dubernet, M.-L. 2006, ApJ,
648, 461
De Jong, T., Boland, W., & Dalgarno, A. 1980, A&A,
91, 68
De Jong, T., Dalgarno, A., & Chu, S.-I. 1975, ApJ, 199,
Dubernet, M. L. 2005, in IAU Symposium, ed. D. C.
Lis, G. A. Blake, & E. Herbst, 235
Elitzur, M. 1992, Astronomical masers (Kluwer
Academic Publishers)
Elitzur, M. & Asensio Ramos, A. 2006, MNRAS, 365,
Evans, II, N. J., Lee, J.-E., Rawlings, J. M. C., & Choi,
M. 2005, ApJ, 626, 919
Expósito, J. P. F., Agúndez, M., Tercero, B., Pardo, J. R.,
& Cernicharo, J. 2006, ApJ, 646, L127
Ferland, G. J. 2003, ARA&A, 41, 517
Fixsen, D. J., Bennett, C. L., & Mather, J. C. 1999, ApJ,
526, 207
Fixsen, D. J. & Mather, J. C. 2002, ApJ, 581, 817
Flower, D. R. 1989, Physics Reports, 174, 1
Fuente, A., Black, J. H., Martı́n-Pintado, J., et al. 2000,
ApJ, 545, L113
Genzel, R. 1991, in NATO ASIC Proc. 342: The Physics
of Star Formation and Early Stellar Evolution, ed.
C. J. Lada & N. D. Kylafis, 155
Goicoechea, J. R., Pety, J., Gerin, M., et al. 2006, A&A,
456, 565
Goldreich, P. & Scoville, N. 1976, ApJ, 205, 144
Goldsmith, P. F. & Langer, W. D. 1999, ApJ, 517, 209
Gray, M. D. & Field, D. 1995, A&A, 298, 243
Gredel, R., Lepp, S., Dalgarno, A., & Herbst, E. 1989,
ApJ, 347, 289
Hauschildt, H., Güsten, R., Phillips, T. G., et al. 1993,
A&A, 273, L23
Helmich, F. P., Jansen, D. J., de Graauw, T., Groesbeck,
T. D., & van Dishoeck, E. F. 1994, A&A, 283, 626
Helmich, F. P. & van Dishoeck, E. F. 1997, A&AS, 124,
Hogerheijde, M. R., Jansen, D. J., & van Dishoeck, E. F.
1995, A&A, 294, 792
Hogerheijde, M. R. & van der Tak, F. F. S. 2000, A&A,
362, 697
Jakob, H., Kramer, C., Simon, R., et al. 2007, A&A,
461, 999
Jansen, D. J. 1995, Ph.D. Thesis, Leiden University
Jansen, D. J., van Dishoeck, E. F., & Black, J. H. 1994,
A&A, 282, 605
Johnstone, D., Boonman, A. M. S., & van Dishoeck,
E. F. 2003, A&A, 412, 157
Leung, C.-M. & Liszt, H. S. 1976, ApJ, 208, 732
Leurini, S., Schilke, P., Menten, K. M., et al. 2004,
A&A, 422, 573
Mangum, J. G. & Wootten, A. 1993, ApJS, 89, 123
Maret, S., Ceccarelli, C., Tielens, A. G. G. M., et al.
2005, A&A, 442, 527
Mathis, J. S., Mezger, P. G., & Panagia, N. 1983, A&A,
128, 212
Mihalas, D. 1978, Stellar atmospheres (2nd edition)
(San Francisco, W. H. Freeman and Co.)
Müller, H. S. P., Thorwirth, S., Roth, D. A., &
Winnewisser, G. 2001, A&A, 370, L49
Ossenkopf, V., Trojan, C., & Stutzki, J. 2001, A&A,
378, 608
Osterbrock, D. E. & Ferland, G. J. 2006, Astrophysics
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra 11
Poelman, D. R. & Spaans, M. 2005, A&A, 440, 559
Prasad, S. S. & Tarafdar, S. P. 1983, ApJ, 267, 603
Roueff, E., Felenbok, P., Black, J. H., & Gry, C. 2002,
A&A, 384, 629
Rybicki, G. B. & Lightman, A. P. 1979, Radiative
processes in astrophysics (New York, Wiley-
Interscience)
Schöier, F. L. 2000, Ph.D. Thesis, Stockholm University
Schöier, F. L., Jørgensen, J. K., van Dishoeck, E. F., &
Blake, G. A. 2002, A&A, 390, 1001
Schöier, F. L., van der Tak, F. F. S., van Dishoeck, E. F.,
& Black, J. H. 2005, A&A, 432, 369
Sobolev, V. 1960, Moving envelopes of stars (Harvard
University Press)
Spaans, M. & van Langevelde, H. J. 1992, MNRAS,
258, 159
Stutzki, J. & Winnewisser, G. 1985, A&A, 144, 13
Tafalla, M., Myers, P. C., Caselli, P., Walmsley, C. M.,
& Comito, C. 2002, ApJ, 569, 815
Takahashi, J. & Uehara, H. 2001, ApJ, 561, 843
Takahashi, T., Silk, J., & Hollenbach, D. J. 1983, ApJ,
275, 145
Van der Tak, F., Neufeld, D., Yates, J., et al. 2005, in
The Dusty and Molecular Universe: A Prelude to
Herschel and ALMA, ed. A. Wilson, 431–432
Van der Tak, F. F. S. 2005, in IAU Symposium
227: Massive Star Birth, ed. R. Cesaroni, M. Felli,
E. Churchwell, & M. Walmsley (Cambridge:
University Press), 70–79
Van der Tak, F. F. S. 2006, Phil. Trans. R. Soc. Lond.,
364, 3101
Van der Tak, F. F. S., van Dishoeck, E. F., & Caselli, P.
2000a, A&A, 361, 327
Van der Tak, F. F. S., van Dishoeck, E. F., Evans, II,
N. J., Bakker, E. J., & Blake, G. A. 1999, ApJ, 522,
Van der Tak, F. F. S., van Dishoeck, E. F., Evans, II,
N. J., & Blake, G. A. 2000b, ApJ, 537, 283
Van Dishoeck, E. F., Blake, G. A., Jansen, D. J., &
Groesbeck, T. D. 1995, ApJ, 447, 760
Van Dishoeck, E. F. & Hogerheijde, M. R. 1999, in
NATO ASIC Proc. 540: The Origin of Stars and
Planetary Systems, ed. C. J. Lada & N. D. Kylafis,
Van Zadelhoff, G.-J., Dullemond, C. P., van der Tak,
F. F. S., et al. 2002, A&A, 395, 373
Wright, E. L., Mather, J. C., Bennett, C. L., et al. 1991,
ApJ, 381, 200
Yates, J. A., Field, D., & Gray, M. D. 1997, MNRAS,
285, 303
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra, Online Material p 1
Online Material
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra, Online Material p 2
Appendix A: Program input and output
A.1. Program input
The input parameters to RADEX are the following:
1. The name of the molecular data file to be used.
2. The name of the file to write the output to.
3. The frequency range for the output file [GHz]. All
transitions from the molecular data file are always
taken into account in the calculation, but often it is
practical to write only a limited set of lines to the
output.
4. The kinetic temperature of the cloud [K].
5. The number of collision partners to be used. Most
users will want H2 as only collision partner, but in
more specialized cases, additional collisions with
H or electrons may for instance play a role. See
the molecular datafiles for details. For some species
(CO, atoms) separate collision data for ortho and
para H2 exist; the program then uses the thermal or-
tho/para ratio unless the user specifies otherwise.
6. The name (case-insensitive) and the density [cm−3]
of each collision partner. Possibilities are H2, p-H2,
o-H2, electrons, atomic H, He, and H
7. The temperature of the background radiation field
– If >0, a black body at this temperature is used.
Most users will adopt the cosmic microwave
background at TCMB=2.725(1+z) K for a galaxy
at redshift z.
– If =0, the average interstellar radiation field
(ISRF) is used, taken from Black (1994) with
modifications described in § 3.2. This spectrum
is not adjustable by a scale factor because it con-
sists of several components that are not expected
to scale linearly with respect to each other.
– If <0, a user-defined radiation field is used, spec-
ified by values of frequency [cm−1], intensity
[Jy nsr−1], and dilution factor [dimensionless].
Spline interpolation and extrapolation are ap-
plied to this table. The intensity need not be
specified at all frequencies of the line list, but
a warning message will appear if extrapolation
(rather than interpolation) is required.
8. The column density of the molecule [cm−2].
9. The FWHM line width [km s−1].
A.2. Program output
The output file written by RADEX first replicates the in-
put parameters, and then lists the following quantities
for each spectral line within the specified frequency
range.
1. Quantum numbers, upper state energy [K], fre-
quency [GHz], and wavelength [µm]. These num-
bers are just copied from the molecular data file,
which usually comes from the LAMDA database.
Frequencies from this database are generally of
line catalogs such as CDMS (Müller et al. 2001)12
should be consulted.
2. The excitation temperature [K] as defined in
Eq. (10). In general, different lines have different
excitation temperatures. Lines are thermalized if
Tex=Tkin; in LTE, all lines are thermalized.
3. The line optical depth, defined as the optical depth
of the equivalent rectangular line shape (φν = 1/∆ν).
4. The line intensity, defined as the Rayleigh-Jeans
equivalent temperature TR [K].
5. The line flux, defined as the velocity-integrated
intensity, both in units of K km s−1 (common in
radio astronomy) and of erg cm−2 s−1 (common
in infrared astronomy). The line flux is calcu-
lated as 1.0645TR∆V, where the factor 1.0645 =√
ln 2) converts the adopted rectangular line
profile into a Gaussian profile with an FWHM of
∆V. The integrated profile is useful to estimate the
total emission in the line, but it has limited mean-
ing at high optical depths, because the change of
optical depth over the line profile is not taken into
account. Proper modeling of optically thick lines re-
quires programs that resolve the source both spec-
trally and spatially (see § 4.1.2 for further discus-
sion).
Auxiliary output files can be generated, for example to
display the adopted continuum spectrum.
Appendix B: Coding standards
The original version of RADEX was written in such a
way as to minimize the use of machine memory which
was expensive until a decade ago. Nowadays, clarity
and easy maintenance are more important requirements,
which is why the source code has been re-written fol-
lowing the rules below. We hope that these rules will be
useful for the development of other ‘open source’ as-
tronomical software. For further guidelines on scientific
programming we recommend the Software Carpentry13
on-line course.
1. All the action is in subroutines; the sole purpose of
the main program is to show the structure of the pro-
gram. The subroutines are grouped into several files
for a better overview; compilation instructions for
automated builds on a variety of platforms are in a
Makefile.
2. The program text is interspersed with comments at
a ratio of ≈1:1. In particular, each subroutine starts
with a description of its contents, its input and out-
put, where incoming calls come from, and which
calls go out. Then the properties of each variable are
described: contents, units and type.
3. Variables and subroutines have descriptive names
with a length of 5–10 letters. Names of integer vari-
ables start with the letters i..n; names of floating-
point variables (always of double precision) with
a..h or o..z. There are no specific namings for vari-
ables of character or logical type, as for example
http://cdms.de
http://www.swc.scipy.org/
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra, Online Material p 3
in CLOUDY. Names are always based on the English
language. We do not use upper-, lower-, or camel-
case to distinguish types of identifiers, as some pro-
grams do.
4. Loops are marked by indenting the program text.
The loop variables are always called ilev, iline ..,
never just i.
5. Subroutines start with a check whether the input pa-
rameters have reasonable values. Such checks force
soft landings if necessary, and avoid runtime errors.
6. Statements with calculations use spaces around
the =, + and – symbols, but not around others.
Calculations that consist of multiple steps are split
over as many program lines. Multiple assignments
in a row are aligned at the = sign.
7. The program text avoids “magic numbers” both in
calculations and in definitions. Numbers that are of-
ten used such as Planck’s constant are defined at one
central place in the program. Similarly, often-used
variables such as physical parameters are stored in
shared memory rather than passed on via subroutine
calls.
Appendix C: Diagnostic plots of molecular
line ratios
We have used the RADEX program and the LAMDA
database to calculate line ratios of several commonly
observed molecules for a range of kinetic tempera-
tures and H2 densities. The plots in this Appendix may
be used by observers to estimate physical conditions
from their data. Line ratios have the advantage of be-
ing less sensitive to calibration errors than absolute line
strengths, especially when the two lines have been mea-
sured with the same telescope, receiver and spectrome-
The calculations assume a column density of 1012 cm−2,
a line width of 1.0 km s−1, and use a 2.73 K blackbody
as background radiation field. Under these conditions,
the lines are optically thin, so that the line ratios do not
depend on column density. The calculations also assume
that the emission in both lines fills the telescope beams
equally, which may be the case if the lines are close
in frequency. However, lines are generally measured in
beams of different sizes, and the observations need to
be corrected to account for this effect, if the source is
known to be compact.
Linear molecules such as CO (Fig. C.1) are tracers
of density at low densities, when collisions compete
with radiative decay. At higher densities, the excitation
becomes thermalized and the line ratios are sensitive
to temperature. For a given molecule, moving up the
J−ladder means probing higher temperatures and densi-
ties. Note that for the column densities of typical dense
interstellar clouds, the CO lines are optically thick, and
observations of 13CO or even rarer isotopologues must
be used to probe physical conditions.
The critical densities of molecular lines scale as µ2ν3,
where µ is the permanent dipole moment of the
molecule and ν is the frequency of the line. Indeed, the
CS molecule (Fig. C.2) has a larger dipole moment than
CO, and its line ratios are mainly probes of the den-
sity. The small frequency spacing between the lines of
CS makes this molecule very useful to probe density
structure (e.g., Van der Tak et al. 2000b). The HCO+
and HCN molecules (Fig. C.3) display similar trends to
CS, although their line spacing is not as small.
Non-linear molecules such as H2CO (Figs. C.4 and C.5)
have the advantage that both temperature and density
may be probed within the same frequency range. Ratios
of lines from different J−states tend to be density tracers
(left panels), while ratios of lines from the same J−state
but different K−states are mostly temperature probes
(right panels). The lines of H2CO are often quite strong,
making this molecule a favourite tracer of temperature
and density (Mangum & Wootten 1993). Other asym-
metric molecules have also been used, such as H2CS
(Blake et al. 1994) and CH3OH (Leurini et al. 2004) al-
though abundance variations from source to source or
even within sources often complicate the interpretation
(Johnstone et al. 2003).
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra, Online Material p 4
Appendix D: The Python scripts
The RADEX distribution comes with two scripts,
radex line.py and radex grid.py, to automate
standard modeling procedures. The scripts are written
in Python and are run from the Unix shell command
line after manual editing of parameters.
The first script, radex line.py, calculates the column
density of a molecule from an observed line intensity,
given estimates of kinetic temperature and H2 volume
density. The input parameters are:
1. The kinetic temperature [K]
2. The number density of H2 molecules [cm
3. The temperature of the background radiation field
[K], usually 2.73 (CMB).
4. The name of the molecule (or molecular data file)
5. The frequency of the line [GHz]
6. The observed line intensity [K]
7. The observed line width [km s−1]
Furthermore, two numerical parameters have good de-
fault values but will need to be changed occasionally:
1. The free spectral range around the line (default
10%): this number must be smaller for molecules
with many line close in frequency, such as CH3OH.
The program uses this parameter to find the ob-
served line from the list of lines in the molecular
model.
2. The required accuracy (default 10%): The default
corresponds to the calibration uncertainty of most
telescopes.
The script iterates on column density until the observed
and modeled line fluxes agree to within the desired ac-
curacy. The best-fit column density is directly written
to the screen. The file radex.out gives details of the
best-fit model.
The second script, radex grid.py, runs a series of
RADEX models to estimate the kinetic temperature
and/or the volume density from an observed line ratio.
The user needs to set the following input parameters:
1. The grid boundaries: minimum and maximum ki-
netic temperature [K] and minimum and maximum
H2 volume density [cm
2. The temperature of the background radiation field
[K], usually 2.73 (CMB).
3. The molecular column density [cm−2]. For the
illustrative plots in Appendix C, a low value
(N=1012 cm−2) was used, so that the line ratios are
independent of column density (optically thin limit).
However, in modeling specific observations, it is
worth varying this parameter to assess the sensitivity
of the line ratio to column density.
4. The observed line width [km s−1], usually an aver-
age of the widths of the two lines.
1. The number of grid points along the temperature and
density axes.
2. The free spectral range around the line (see above)
The name of the molecule and a list of observed line
ratios and names of associated output files are given at
the start of the main program. The script produces a file
radex.out which is a tabular listing of temperature,
log density, and line ratio. This results may be plotted
with the user’s favourite plotting program.
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra, Online Material p 5
Fig. C.1. Line ratios of CO in the optically thin limit
as a function of kinetic temperature and H2 density.
Contours are spaced linearly and some contours are la-
beled for easy identification.
Fig. C.2. Line ratios of CS in the optically thin limit
as a function of kinetic temperature and H2 density.
Contours are spaced linearly and some contours are la-
beled for easy identification.
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra, Online Material p 6
Fig. C.3. Line ratios of HCO+ and HCN in the optically
thin limit as a function of kinetic temperature and H2
density. Contours are spaced linearly and some contours
are labeled for easy identification.
Fig. C.4. Line ratios of o-H2CO in the optically thin
limit as a function of kinetic temperature and H2 den-
sity. Contours are spaced linearly and some contours are
labeled for easy identification.
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra, Online Material p 7
Fig. C.5. Line ratios of p-H2CO in the optically thin
limit as a function of kinetic temperature and H2 den-
sity. Contours are spaced linearly and some contours are
labeled for easy identification.
	Introduction
	Radiative transfer and molecular excitation
	Basic formalism
	Molecular line cooling
	Escape probability
	The program
	Basic capabilities
	Background radiation field
	Chemical formation and destruction rates
	Calculation
	Results
	Limitations of the program
	Comparison with other methods
	Homogeneous models
	The case of CO
	The case of HCO+
	Observations of a Young Stellar Object
	Conclusions
	Program input and output
	Program input
	Program output
	Coding standards
	Diagnostic plots of molecular line ratios
	The Python scripts

  The large quantity and high quality of modern radio and infrared line
observations require efficient modeling techniques to infer physical and
chemical parameters such as temperature, density, and molecular abundances. We
present a computer program to calculate the intensities of atomic and molecular
lines produced in a uniform medium, based on statistical equilibrium
calculations involving collisional and radiative processes and including
radiation from background sources. Optical depth effects are treated with an
escape probability method. The program is available on the World Wide Web at
http://www.sron.rug.nl/~vdtak/radex/index.shtml . The program makes use of
molecular data files maintained in the Leiden Atomic and Molecular Database
(LAMDA), which will continue to be improved and expanded. The performance of
the program is compared with more approximate and with more sophisticated
methods. An Appendix provides diagnostic plots to estimate physical parameters
from line intensity ratios of commonly observed molecules. This program should
form an important tool in analyzing observations from current and future radio
and infrared telescopes.

Introduction
Observations of spectral lines at radio, (sub)millimeter
and infrared wavelengths are a powerful tool to in-
vestigate the physical and chemical conditions in the
dilute gas of astronomical sources where thermody-
namic equilibrium is a poor approximation (e.g., Genzel
1991; Black 2000). To extract astrophysical parame-
ters from the data, the excitation and optical depth of
the lines need to be estimated, for which various meth-
ods may be used, depending on the available observa-
tions (Van Dishoeck & Hogerheijde 1999; Van der Tak
2005).
If only one or two lines of a molecule1 have been ob-
served, the excitation must be deduced from observa-
tions of other species or from theoretical considerations.
An example is the assumption that the excitation tem-
perature equals the kinetic temperature, a case known as
Local Thermodynamic Equilibrium (LTE) which holds
at high densities.
If many lines have been observed, a popular method
is the ‘rotation diagram’, also called ‘Boltzmann
plot’ or ‘population diagram’ (e.g., Blake et al. 1987;
Helmich et al. 1994; Goldsmith & Langer 1999). This
method describes the excitation by a single tempera-
ture, obtained by a fit to the line intensities as a func-
tion of upper level energy. Provided that beam sizes are
similar and optical depths are low, or that appropriate
corrections are made, this method yields estimates of
the excitation temperature and column density of the
molecule. The excitation temperature approaches the ki-
netic temperature in the high-density limit, but generally
depends on both kinetic temperature and volume den-
sity. Spectral line surveys are often analyzed with ro-
tation diagrams, although more advanced methods are
also used (Helmich & van Dishoeck 1997; Comito et al.
2005).
More sophisticated methods retain the assumption of
a local excitation, but solve for the balance of ex-
citation and de-excitation rates from and to a given
http://arxiv.org/abs/0704.0155v1
http://www.sron.rug.nl/~vdtak/radex/index.shtml
2 Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra
probability method and the Large Velocity Gradient
(LVG) method (Sobolev 1960; De Jong et al. 1975;
Goldreich & Scoville 1976). These ‘intermediate-level’
methods require knowledge of molecular collisional
data, whereas the previous ‘basic-level’ methods only
required spectroscopic and dipole moment information.
This extra requirement limits the use of these methods
to some extent, because collisional data do not exist
for all astrophysically relevant species. The advantage
is that column density, kinetic temperature and volume
density can be constrained, if accurate collision rates
are known. As with rotation diagrams, this method can
be used to compute synthetic spectra to be compared
with data with a χ2 statistic (Jansen 1995; Leurini et al.
2004).
The most advanced methods drop the local approxi-
mation and solve for the intensities (or the radiative
rates) as functions of depth into the cloud, as well as of
velocity. Such methods are usually of the Accelerated
Lambda Iteration (ALI) or Monte Carlo (MC) type, al-
though hybrids also exist. The performance and con-
vergence of such programs have recently been tested
by Van Zadelhoff et al. (2002). Using such programs
one can constrain temperature, density, and velocity
gradients within sources (e.g., Van der Tak et al. 1999;
Tafalla et al. 2002; Jakob et al. 2007), and, if enough
observations are available, even molecular abundance
profiles (e.g., Van der Tak et al. 2000a; Schöier et al.
2002; Maret et al. 2005), especially when coupled to
chemical networks (e.g., Doty et al. 2004; Evans et al.
2005; Goicoechea et al. 2006).
This paper presents the public version of a radiative
transfer code at the ‘intermediate’ level. The assump-
tion of a homogeneous medium limits the number of
free parameters and makes the program a useful tool
in rapidly analyzing a large set of observational data,
in order to provide constraints on physical conditions,
such as density and kinetic temperature (Jansen 1995).
The program can be used for any molecule for which
collisional rate coefficients are available. The input for-
mat for spectroscopic and collisional data is that of the
LAMDA database (Schöier et al. 2005)2 where an on-
line calculator for molecular line intensities3, based on
our program, can also be found4.
The paper is set up as follows. Section 2 describes the
radiative transfer formalism and introduces our notation
of the key quantities. Section 3 describes the formalism
which the program actually uses, and discusses its im-
plementation. Section 4 compares the results of the pro-
gram to those of other programs. The paper concludes
in § 5 with suggested future directions of astrophysical
radiative transfer modeling.
2 http://www.strw.leidenuniv.nl/∼moldata
3 In this paper, the ‘strength’ of a line is an intrinsic quantity
2. Radiative transfer and molecular
excitation
This section summarizes the formalism to analyze
molecular line observations which our program adopts.
For more detailed discussions of radiative transfer see,
e.g., Cannon (1985) or Rybicki & Lightman (1979).
2.1. Basic formalism
Describing the transfer of radiation requires a quantity
which is conserved along its path as long as no local
absorption or emission processes take place, and which
includes the direction of travel. The quantity that sat-
isfies this requirement is the specific intensity Iν, de-
fined as the amount of energy passing through a surface
normal to the path, per unit time, surface, bandwidth
(measured here in frequency units), and solid angle. The
transfer equation for radiation propagating a distance ds
can then be written as
= jν − αν Iν, (1)
where jν and αν are the local emission and extinction
coefficients, respectively. The two terms on the right-
hand side may be combined into the source function,
defined by
S ν ≡
. (2)
Writing the transport equation in its integral form and
defining the optical depth, dτν ≡ αν ds, measured along
the ray5 one arrives at
Iν = Iν (0)e
−τν +
S ν(τ
−(τν−τ′ν) dτ′ν, (3)
where Iν is the radiation emerging from the medium and
Iν(0) is the ‘background’ radiation entering the medium.
The above equations hold both for continuum radiation,
which is emitted over a large bandwidth, and for spectral
lines, which arise when the local emission and absorp-
tion properties change drastically over a very small fre-
quency interval, due to the presence of molecules. From
this point the discussion will focus on bound-bound
transitions within a multi-level molecule consisting of
N levels with spontaneous downward rates Aul, Einstein
coefficients for stimulated transitions Bul and Blu, and
collisional rates Cul and Clu, between upper levels u and
lower levels ℓ.
The rate of collision is equal to
Cul = ncolγul, (4)
where ncol is the number density of the collision partner
(in cm−3) and γul is the downward collisional rate coeffi-
cient (in cm3 s−1). The rate coefficient is the Maxwellian
average of the collision cross section, σ,
γul =
8kTkin
)−1/2 (
kTkin
σEe−E/kTkin dE, (5)
http://www.strw.leidenuniv.nl/~moldata
http://www.sron.rug.nl/~vdtak/radex.php
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra 3
where E is the collision energy, k is the Boltzmann con-
stant, Tkin is the kinetic temperature, and µ is the re-
duced mass of the system. The upward rates are ob-
tained through detailed balance
γlu = γul
e−hν/kTkin , (6)
where gi is the statistical weight of level i.
The local emission in transition u→l with laboratory
frequency νul, can be expressed as
nu Aul φν, (7)
where nu is the number density of molecules in level u
and φν is the frequency-dependent line emission profile.
The absorption coefficient reads
(nl Blu ϕν − nu Bul χν) , (8)
where φν and χν are the line profiles for absorption and
stimulated emission (counted as negative extinction),
respectively.
From here on we assume complete angular and fre-
quency redistribution of the emitted photons, so that
φν=ϕν=χν, which is strictly only valid when collisional
excitation dominates. This assumption allows the source
function to be written as
S νul =
nuAul
nlBlu − nuBul
, (9)
where we have used the Einstein relations. It is com-
mon to introduce an excitation temperature Tex defined
through the Boltzmann equation
exp [−(Eu − El)/kTex] , (10)
where Ei is the energy of level i, such that S νul =
Bν(Tex), the specific intensity of a blackbody radiating
at Tex.
In the interstellar medium, the dominant line broadening
mechanism is Doppler broadening. Except in very cold
and dark cloud cores, observed line widths are much
larger than expected from the kinetic temperature: this
effect is commonly ascribed to random macroscopic gas
motions or ‘turbulence’. The result is a Gaussian line
profile
ν − νul − v · n
, (11)
where νD is the Doppler width, v is the velocity vector
of the moving gas at the position of the scattering, n is
a unit vector in the direction of the propagating beam of
radiation, and c is the speed of light. The Doppler width
is the 1/e half-width of the profile, equal to ∆V/2
where ∆V is its full width at half-maximum.
If the level populations ni are known, the radiative trans-
fer equation can be solved exactly. In particular, under
circumstellar media, the density is too low to attain LTE,
but statistical equilibrium (SE) can often be assumed:
= 0 =
n jP ji − ni
Pi j = Fi − niDi , (12)
where Pi j, the destruction rate coefficient of level i, and
its formation rate coefficient P ji are given by
Pi j =
Ai j + Bi j J̄ν + Ci j (i > j)
Bi j J̄ν +Ci j (i < j).
In Eq. 13,
Bi j J̄ν = Bi j
Jν φ(ν) dν (14)
is the number of induced radiative (de-)excitations from
state i to state j per second per particle in state i, and
Iν dΩ (15)
is the specific intensity Iν integrated over solid angle dΩ
and averaged over all directions. The SE equations thus
include the effects of non-local radiation.
This discussion assumes that the state-specific rates of
formationFi [cm3 s−1] and destructionDi [s−1] are zero
to ensure that the radiative transfer is solved indepen-
dently of assumptions about chemical processes. In gen-
eral, formation and destruction processes should be in-
cluded explicitly to be able to deal with situations in
which the chemical time scales are very short or the
radiative lifetimes very long. For example, the forma-
tion temperature (in Fi) affects the rotational excitation
of C3 (Roueff et al. 2002) and the vibrational excitation
of H2 (Black & van Dishoeck 1987; Burton et al. 1990;
Takahashi & Uehara 2001), systems for which line ra-
diation only occurs as slow electric quadrupole transi-
tions. The rotational excitation of reactive ions like CO+
(Fuente et al. 2000; Black 1998) is also sensitive to Fi
and Di because the rates of reactions with H and H2
rival the inelastic collision excitation rates. Similar con-
siderations apply to the excitation of H+3 in the Sgr A
region close to the Galactic Center (Van der Tak 2006),
where electron impact excitation competes with disso-
ciative recombination.
2.2. Molecular line cooling
Once the radiative transfer problem has been solved and
the level populations are known, the cooling (or heating)
from molecular line emission can be estimated. Since
the level populations contain all the information of the
radiative transfer, a general expression for the cooling is
obtained from considering all possible collisional tran-
sitions
(nlγlu − nuγul)hνul, (16)
4 Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra
are in detailed balance at the kinetic temperature; there-
fore it is possible for net heating to occur (Λ < 0) in
cases where the crucial level populations have Tex>Tkin,
owing to strong radiative excitation in a hot external ra-
diation field.
2.3. Escape probability
The difficulty in solving radiative transfer problems is
the interdependence of the molecular level populations
and the local radiation field, requiring iterative solution
methods. In particular, for inhomogeneous or geometri-
cally complex objects, extensive calculations with many
grid points are required. However, if only the global
properties of an interstellar cloud are of interest, the cal-
culation can be greatly simplified through the introduc-
tion of a geometrically averaged escape probability β,
the probability that a photon will escape the medium
from where it was created. This probability depends
only on the optical depth τ and is related to the inten-
sity within the medium, ignoring background radiation
and any local continuum, through
Jνul = S νul (1 − β). (17)
Several authors have developed detailed relations be-
tween β and τ for specific geometrical assumptions.
Our program offers the user a choice of three such ex-
pressions. The first is the expression derived for an ex-
panding spherical shell, the so-called Sobolev or large
velocity gradient (LVG) approximation (Sobolev 1960;
Castor 1970; Elitzur 1992, p. 42-44). This method is
also widely applied for moderate velocity gradients, to
mimic turbulent motions. Our program uses the formula
by Mihalas (1978) and De Jong et al. (1980) for this ge-
ometry:
βLVG =
dτ′ =
1 − e−τ
. (18)
Second, in the case of a static, spherically symmet-
ric and homogeneous medium the escape probability is
(Osterbrock & Ferland 2006, Appendix 2)6
βsphere =
. (19)
Third, for a plane-parallel ‘slab’ geometry, applicable
for instance to shocks,
βslab =
1 − e−3τ
is derived (De Jong et al. 1975). Figure 1 plots the be-
haviour of β as a function of τ for these three cases; for
more detailed comparisons see Stutzki & Winnewisser
(1985) and Ossenkopf et al. (2001). Users of our pro-
gram can select either expression for their calculations.
The on-line version of the program uses the formula for
the uniform sphere, Eq. (19).
 0  1  2  3  4  5  6  7  8  9  10
Optical depth τ
Sphere
Fig. 1. Escape probability β as a function of optical
depth τ for three different geometries: uniform sphere
(solid line), expanding sphere (dotted line) and plane-
parallel slab (dashed line).
3. The program
RADEX is a non-LTE radiative transfer code, written
originally by J. H. Black, that uses the escape proba-
bility formulation assuming an isothermal and homo-
geneous medium without large-scale velocity fields.
With the current increase of observational possibilities
in mind, we have developed a version of this program
which is suitable for public use. A guide for using the
code in practice is provided in Appendix A and on-line7;
Appendix B describes the adopted coding style. This
section focuses on the implementation of the formalism
of § 2 in the program.
3.1. Basic capabilities
For a homogeneous medium with no global velocity
field, the optical depth at line centre can be expressed
using Eqs. (2, 7, 9, 11), as
AulNmol
1.064∆V
, (21)
where Nmol is the total column density,∆V the full width
at half-maximum of the line profile in velocity units,
and xi the fractional population of level i. The formal-
ism is analogous to the LVG method, with the global
n/(dV/dR) replaced by the local N/∆V , as in microtur-
bulent codes (Leung & Liszt 1976). The program itera-
tively solves the statistical equilibrium equations start-
ing from optically thin statistical equilibrium (§ 3.4) for
the initial level populations.
The program can handle up to seven collision partners
simultaneously. In dense molecular clouds, H2 is the
main collision partner for most species, but in some
cases, separate cross sections may exist for collisions
with the ortho and para forms of H2, and electron col-
lisions may be important for ionic species. In diffuse
molecular clouds and PDRs, excitation by atomic H be-
comes important, particularly for fine structure lines,
while for comets, H2O is the main collision partner. We
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra 5
refer to Flower (1989) for the basic theory of molecular
collisions, and to Dubernet (2005) for an update of the
latest results.
The output of the program is the background-subtracted
line intensity in units of the equivalent radiation temper-
ature in the Rayleigh-Jeans limit. The background sub-
traction follows traditional cm- and mm-wave spectro-
scopic observations where the differences between on-
source and off-source measurements are recorded, such
Iemν − I
. (22)
The radiation peak temperature TR can be directly com-
pared to the observed antenna temperature corrected
for the optical efficiency of the telescope. However, it
should be emphasized that RADEX contains no informa-
tion about the geometry or length scale and that it is
assumed that the source fills the antenna beam. If the
source is expected to be smaller than the observational
beam, computed line fluxes must be corrected before
comparing to observed fluxes.
In other types of observations, the continuum may
not be subtracted from the data. In (sub-)millimeter
and THz observations, for example with ESA’s fu-
ture Herschel space observatory, the dust continuum
of many sources will be much stronger than any in-
strumental error, and baseline subtraction may not be
needed. The same is true for interferometer data, where
the instrumental passband is well characterized.
3.2. Background radiation field
The average Galactic background (interstellar radia-
tion field, ISRF) adopted in RADEX consists of several
components. The main contribution is the cosmic mi-
crowave background (CMB) whose absolute tempera-
ture is taken to be TCBR = 2.725±0.001K based on the
full COBE data set as analyzed by Fixsen & Mather
(2002). This model of the microwave background rep-
resents the broadband continuum only and does not in-
clude the strong emission lines, several of which contain
significant power in the far-infrared and (sub-)millime-
ter part of the spectrum (see, e.g., Fixsen et al. 1999).
The ultraviolet/visible/near-infrared part of the spec-
trum is based on the model of average Galactic starlight
in the solar neighborhood of Mathis et al. (1983). The
far-infrared and (sub-)millimeter part of the spectrum is
based on the single-temperature fit to the Galactic ther-
mal dust emission of Wright et al. (1991). At frequen-
cies below 10 cm−1 (30 GHz), there is a background
contribution from non-thermal radiation in the Galaxy.
A tabulation of this spectrum in ASCII format is avail-
able on-line8, and a graphical representation is shown in
Black (1994).
One subtle aspect of the calculation is the distinction
between the background seen by the observer and the
background seen by the molecules. The continuum con-
tribution to the rate equations may be composed of (1)
an external component which arises outside the emitting
region and (2) an internal continuum that arises within
the emitting region. The CMB and ISRF are examples
of external continuum components; dust emission from
the line-emitting region is an example of an internal
continuum. While an external continuum always fills
the entire sky, an internal continuum may only fill a frac-
tion of it, for example in the case of a circumstellar disk.
With this distinction in mind, the internal intensity be-
comes
Jintν = β[Bν(TCBR)+ηI
ν ]+(1−β)[Bν(Tex)+θ(1−η)I
where Iuserν is the continuous spectrum defined by the
user. The factor η, is the fraction of local continuum
which arises outside the line emitting region, and the
factor θ is the fraction of local sky filled by the internal
continuum.
3.3. Chemical formation and destruction rates
The equations of statistical equilibrium (12) include
source and sink terms. By default, RADEX sets the
destruction rates equal to the same small value,
Di ≡ D = 10−15 s−1, appropriate for cosmic-ray
ionization plus cosmic-ray induced photodissociation
(Prasad & Tarafdar 1983; Gredel et al. 1989). The cor-
responding formation rates are
Fi = 10−24ntotalgi exp(−Ei/kTform)/Q(Tform) (24)
where ntotal is the sum of the densities of all collision
partners, Tform is a formation temperature (default value
300 K), and
Q(T ) =
gi exp(−Ei/kTform) (25)
is the partition function. These assumptions imply a
nominal fractional abundance of every molecule
ntotal
ntotalD
= 10−9 . (26)
The value of the nominal abundance is inconsequential
because the results in RADEX depend on Nmol/∆V , but
not on the fractional abundance. For most molecules
currently in the associated database (LAMDA) and for
the most commonly encountered interstellar conditions,
these choices will not affect the observable excitation.
The formation and destruction rates are computed in a
subroutine that can be modified by the user to provide
a more realistic description of chemical processes. For
example, users may treat the combined ortho/para forms
of molecules by introducing a realistic Tform, especially
in cases where no o/p interchange processes is likely
to be effective. Other cases of potential interest include
the photodissociation of large molecules into smaller
molecules, or the evaporation of icy grain mantles into
the gas phase. Our formulation in terms of a volume
rate of formation is chosen to be independent of the de-
tails of the formation process. In general, formation and
destruction processes are important for molecules that
http://www.oso.chalmers.se/~jblack/RESEARCH/
6 Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra
3.4. Calculation
The input parameters of RADEX and its output are de-
scribed in Appendix A. Calculations with RADEX pro-
ceed as follows. A first guess of the populations of the
molecular energy levels is produced by solving statisti-
cal equilibrium in the optically thin case. The only ra-
diation taken into account is the unshielded background
radiation field; internally produced radiation is not yet
available. The solution for the level populations allows
calculation of the optical depths of all the lines, which
are then used to re-calculate the molecular excitation.
The new calculation treats the background radiation in
the same manner as the internally produced radiation.
The program iteratively finds a consistent solution for
the level populations and the radiation field. When the
optical depths of the lines with τ > 10−2 are stable from
one iteration to the next to a given tolerance (default
10−6), the program writes output and stops.
3.5. Results
There are several ways in which RADEX can be used
to analyze molecular line observations. In most of
these applications, the modeled quantity is the velocity-
integrated line intensity, as the excitation is assumed
to be independent of velocity. As a consequence, self-
absorbed lines cannot be modeled satisfactorily with
RADEX. In the simplest case, the temperature and den-
sity are known from other observations and only the col-
umn density of the molecule under consideration needs
to be varied to get the best agreement with the observed
line intensity. If the H2 column density is known from
other observations, for example from an optically thin
CO isotopic line, the ratio of the two column densi-
ties gives the molecular abundance, averaged over the
source. The RADEX distribution contains a Python script
to automate this procedure which is further described in
Appendix D.
Another often-used application of RADEX is to deter-
mine temperatures and densities from the observed in-
tensity ratios of lines of the same molecule. If the abun-
dance of the molecule is constant throughout the source,
the ratios should give source-averaged physical con-
ditions independent of the specific chemistry of the
molecule. Appendix C presents illustrative plots of line
ratios for commonly observed molecules and lines in
the optically thin case. For higher optical depths, the
qualitative trends remain the same but there are quanti-
tative differences. RADEX can readily be used to generate
similar plots for moderately thick cases. Again, Python
scripts are made available to automate this procedure
(Appendix D).
To illustrate the use of RADEX on actual observations, we
take the observations of the HCO+ 1–0 and 3–2 lines to-
ward a relatively simple source, the photon-dominated
region IC 63 (Jansen et al. 1994). The observed 1–0/3–2
ratio corrected for beam dilution is 5.5±1.5. The kinetic
temperature of the source is constrained from CO ob-
servations to be ∼50 K. Fig. C.3 shows that the inferred
Fig. 2. Comparison of the predicted line strengths for
the 10 lowest rotational transitions of CO for a homo-
geneous isothermal sphere, with nH2 = 10
5 cm−3 and
Tkin = 50 K, using different methods. Upper panel: The
total optical depth through the sphere at line centre τ
and excitation temperature Tex as a function of the up-
per rotational level J involved in the transition. Middle
panel: The radiation temperature TR obtained for each
transition using RADEX with different prescriptions of
the escape probability β and compared with the result
from the Monte-Carlo code (MCC) of Schöier (2000).
Also shown are the results for optically thin emission in
LTE. Lower panel: TR obtained from RADEX compared
with the results from MCC, δTR = T
R − TR.
inferred column density from the absolute intensities is
8 × 1012 cm−2, which, together with the overall H2 col-
umn density of 5×1021 cm−2, gives an HCO+ abundance
with respect to H2 of 1.6 × 10−9.
A slightly more complicated situation arises for the
Orion Bar PDR (Hogerheijde et al. 1995). For this
source, both HCO+ 1–0, 3–2 and 4–3 lines have been
observed. The 1–0/3–2 ratio gives an order of magni-
tude lower density than the 3–2/4–3 ratio. This differ-
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra 7
density tracers such as CO are usually much less sensi-
tive to density. One possible solution is a clumpy PDR
model in which the 1–0 line is mostly produced in the
low-density interclump gas containing 90% of the mate-
rial and the 4–3 line in the high-density clumps. Within
this clumpy model, a single column density fits all three
lines and an accurate abundance can be derived. Note
that this technique of adding the results of two models
is only applicable at low optical depths.
If only two lines of a molecule have been observed,
the line ratio can be used as indicator of tempera-
ture or density, depending on molecule and transition
(Appendix C). A single line ratio is never enough
to constrain both temperature and density, though.
For multi-line observations, a comparison of data and
models in terms of χ2 is preferred. See for example
Van der Tak et al. (2000b), Schöier et al. (2002), and
Leurini et al. (2004) for details of such calculations.
3.6. Limitations of the program
The current version of the program does not in-
clude a contribution from continuous (dust or free-free)
opacity to the escape probability, as for example in
Takahashi et al. (1983). Continuum radiation from dust
is generally negligible at long wavelengths ( >∼ 1 mm)
but becomes important for regions with very high col-
umn densities (such as protoplanetary disks) and at
far-infrared and shorter wavelengths ( <∼ 100µm). Free-
free radiation may become important for the calcula-
tion of atomic fine structure lines from H II regions;
other programs such as CLOUDY (Ferland 2003)9 may
be more suitable for this purpose. The absence of con-
tinuous opacity limits the applicability of the program
particularly in situations where infrared pumping is im-
portant, either directly through rotational transitions or
via vibrational transitions (Carroll & Goldsmith 1981;
Hauschildt et al. 1993).
Another limitation of the program is that only one
molecule is treated at a time, so that the effects of line
overlap are not taken into account. Such overlaps may
occur both at radio and at infrared wavelengths (e.g.,
Expósito et al. 2006). In special cases, overlap between
lines of the same molecule may influence their excita-
tion, for example the hyperfine components of HCN or
+ (Daniel et al. 2006).
For certain molecules under certain physical conditions
(especially low density and/or strong radiation field),
population inversions occur, which cause negative op-
tical depth and hence nonlinear amplification of the
incoming radiation (Elitzur 1992). This phenomenon,
known as ‘maser’ action, requires non-local treatment
of the radiative transfer, in particular a fine sampling of
directions, for which RADEX is not set up. Generally, the
escape probability approximation is justified until the
masers saturate, which occurs at τ ≈ −1. In practice,
the computed intensities of lines with τ <∼ − 0.1 are not
as accurate as those of other lines, and the intensities of
lines with τ <∼ − 1 should be disregarded altogether. If
Fig. 3. Excitation temperature of the HCO+ 1–0 tran-
sition as a function of radius for the model cloud of
§ 4.1.2, calculated with RADEX assuming static spher-
ical, expanding spherical, and slab geometry (dashed
/ dash-dotted / dotted lines), with a multi-zone es-
cape probability program (long/short dashes) and with
a Monte Carlo code (solid lines). The panels are for dif-
ferent column densities, hence optical depths. Note the
different vertical scales.
of non-maser lines may also be affected. While special-
ized programs should be used to calculate the intensities
of maser lines (e.g., Spaans & van Langevelde 1992;
Gray & Field 1995; Yates et al. 1997), RADEX may well
be used to predict which lines of a molecule may display
http://www.nublado.org/
8 Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra
which are pumped by infrared radiation (Leurini et al.
2004).
4. Comparison with other methods
This section shows a comparison of RADEX with other
programs, first for the case of constant physical con-
ditions (§ 4.1) and second for variable conditions
(§ 4.2). Comparison is with the analytical rotation di-
agram method and with Monte Carlo methods, which
have been benchmarked to high accuracy, both for
the case of HCO+ (Van Zadelhoff et al. 2002)10 and of
H2O (Van der Tak et al. 2005)
11. Throughout this sec-
tion, molecular data have been taken from the LAMDA
database (Schöier et al. 2005).
4.1. Homogeneous models
4.1.1. The case of CO
To test the RADEX code, we have compared its output
both to an optically thin LTE analysis (rotation diagram
method) and a full radiative transfer analysis using a
Monte-Carlo method (Schöier 2000). The test problem
consists of a spherically symmetric cloud with a con-
stant density, n(H2), of 1 × 105 cm−3 within a radius of
100 AU. In this example only the CO emission is treated
using a fractional abundance of 1 × 10−4 relative to H2
yielding a central CO column density of NCO = 3 ×
1016 cm−2 and an average value of N=2×1016 cm−2. The
kinetic temperature is set to 50 K, the background tem-
perature to 2.73 K, and the line width to ∆V=1.0 km s−1.
Fig. 2 presents the results of the calculations for the
ten lowest rotational transitions. The excitation temper-
atures of the lines vary from being close to thermalized
for transitions involving low J-levels, to sub-thermally
excited for the higher-lying lines. The optical depth in
the lines is moderate (∼ 1 − 2) to low. It is seen that the
expressions of the escape probability for the uniform
sphere and the expanding sphere give almost identical
solutions which are close to that obtained from the full
radiative transfer (MCC in Fig. 2). The slab geometry
gives slightly higher intensities, in particular for high-
lying lines. The optically thin approximation, where the
gas is assumed to be in LTE at 50 K, produces much
larger discrepancies, up to a factor of∼ 2, and only gives
the correct intensity for the J = 1 → 0 line, where the
LTE conditions are met.
4.1.2. The case of HCO+
To further verify the performance of the RADEX
program, we have compared its results to that
of another program that does not use the local
approximation: the Monte Carlo program RATRAN
(Hogerheijde & van der Tak 2000). We also compare
the results to those from the multi-zone escape proba-
bility program by Poelman & Spaans (2005). The test
case is a cloud with n(H2) = 1×104 cm−3, Tkin=10 K,
Tbg=2.73 K, and a line width of ∆V=1.0 km s
−1, equiv-
alent to bD=0.6 km s
−1. The pure rotational emission
spectrum of HCO+ was calculated for column densi-
ties of 1012, 1013, 1014 and 1015 cm−2, which for RADEX
were given directly as input parameters. For the multi-
zone programs, a cloud radius of 1018 cm was specified
along with abundances of 10−10–10−7, distributed over
50 cells.
Figure 3 shows the calculated excitation temperature
of the HCO+ 1–0 transition as a function of radius
for these physical conditions. For N(HCO+) <∼ 10
12 cm−2,
the excitation is independent of radius and the calcu-
lations for the various geometries agree to ≈10%. The
dependence of the excitation on radius and on geom-
etry increases with increasing column density, and for
N(HCO+) >∼ 10
15 cm−2, the curvature of the Tex distri-
bution becomes too large to ignore. The corresponding
line optical depth is ≈100, with ≈20% spread between
the various estimates (Fig. 4). The curvature arises be-
cause at the cloud center, photon trapping thermalizes
the excitation, while at the edge, the emission can es-
cape the cloud (Bernes 1979). We do not recommend to
use RADEX at line optical depths >∼ 100, because the cal-
culated excitation temperature may not be representa-
tive of the emitting region. However, even if some lines
are highly optically thick, RADEX may well be used to
analyze other lines which are optically thin. For exam-
ple for H2O, the ground state lines often have τ ∼ 1000,
but RADEX is well capable of computing intensities for
higher-lying transitions which are not as optically thick.
At low optical depth, variations in Tex translate directly
into changes in emergent line intensity. Thus, differ-
ences as large as 20% in calculated line flux can arise
depending on the choice of escape probability descrip-
tion, even for moderately thick cases. At high optical
depth, the direct connection between Tex and line flux
is lost because of the dependence on the adopted veloc-
ity field. The assumption in the program that the opti-
cal depth is independent of velocity breaks down in this
case. In this limit, the peak line temperature TR gives
the value of Tex at the τ = 1 surface of the cloud in this
specific transition.
The results shown in this section do not translate easily
to other HCO+ lines such as J=3→2, because the ex-
citation is governed by several competing effects. The
optical depth of the J=3→2 line may be higher or lower
than that of the J=1→0 line, depending on temperature
and density. Observers are encouraged to use RADEX to
study the excitation of their lines as a function of these
parameters, and also consider geometric variations.
4.2. Observations of a Young Stellar Object
To compare a typical RADEX analysis with other meth-
ods for a situation which varying physical conditions,
we choose a molecule for which many lines can be
observed: para-formaldehyde, p-H2CO. Figure 5 shows
observations of p-H2CO (sub-)millimeter emission lines
http://www.strw.leidenuniv.nl/astrochem/
http://www.sron.rug.nl/~vdtak/H2O/
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra 9
Fig. 4. Optical depth of the HCO+ 1–0 transition for
the model of § 4.1.2, calculated with RADEX assuming
static spherical, expanding spherical, and slab geometry
(filled triangles/filled squares/open triangles), with the
multi-zone escape probability code (open circles) and
with the Monte Carlo code (filled circles).
Fig. 5. Line strengths of p-H2CO observed towards
the embedded low-mass protostar IRAS 16293–2422
(squares with error bars), modeled assuming LTE (solid
line), using RADEX (triangles), and using a Monte Carlo
program assuming a constant abundance (solid circles)
and an abundance varying with radius (open circles).
Fig. 6. Distributions of the χ2 parameter correspond-
ing to the models in Figure 5. The RADEX results are
for n(H2) = 10
6 cm−3 as found by Van Dishoeck et al.
(1995).
Van Dishoeck et al. (1995). The data are analyzed us-
ing three methods: assuming LTE (with a rotation dia-
gram), assuming SE (using RADEX), and using a Monte
Carlo program. The free parameters for the LTE fit are
the excitation temperature Tex and the column density
Figure 6 shows the distributions of the χ2 parameter for
the LTE and SE fits, calculated in the standard way (see,
e.g., Van der Tak et al. 2000b; Schöier et al. 2002) as-
suming a 20% uncertainty for all observed points except
the line with the highest upper level energy where a 30%
uncertainty was used. As seen from the figure, the non-
LTE method gives a better fit to the data, as quantified
by the lower minimum χ2 value. This result is not nec-
essarily surprising given that more free parameters are
available. A more important difference is that the esti-
mates of temperature and column density between the
two methods are substantially different, in particular the
temperature (50 vs 150 K). Since the non-LTE method
involves fewer assumptions about the physical state of
the cloud, its results are to be preferred.
These results illustrate that rotation diagrams may give
misleading results when determining physical proper-
ties of interstellar gas clouds (cf. Johnstone et al. 2003
for the case of CH3OH). Figure 5 also demonstrates that
temperatures and column densities derived from rota-
tion diagrams tend to depend on which lines happen to
have been observed (cf. the HCO+ case in § 4.1.2). From
other data, IRAS 16293–2422 is actually known to have
a gradient in temperature and density throughout its en-
velope, which cannot be modelled properly with either
technique. For such situation, a full Monte Carlo radia-
tive transfer method is needed in which both the physi-
cal conditions and the abundances can vary with radius
(Fig. 5, circles). Nevertheless, the column densities and
abundances inferred with RADEX using the physical con-
ditions inferred from the line ratios differ by only a fac-
tor of a few from those found with the more sophis-
ticated analysis, at least for the particular zone of the
source to which those conditions apply (Schöier et al.
2002).
5. Conclusions
We have presented a computer program to analyze spec-
tral line observations at radio and infrared wavelengths,
based on the escape probability approximation. The pro-
gram can be used for any molecule for which collisional
data exist; such input data are available in the required
format from the LAMDA database. The program can be
used for optical depths from ≈–0.1 to ≈100.
The limited number of free parameters makes RADEX
very useful to rapidly analyze large datasets. As an
example, observed line intensity ratios may be com-
pared with the plots in Appendix C to estimate density
and kinetic temperature. Ratios of other lines and other
molecules may be easily computed using the Python
scripts included in the RADEX distribution. The program
may also be used to create synthetic spectra. This capa-
bility will be important to model the THz line surveys
from the HIFI instrument onboard the Herschel space
observatory.
In the future, we plan to incorporate a multi-zone es-
cape probability formalism (Poelman & Spaans 2005;
Elitzur & Asensio Ramos 2006) which will enable
10 Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra
depths, the calculation may also start from LTE condi-
tions rather than from optically thin statistical equilib-
rium. Robust convergence may be achieved by starting
from either initial condition and requiring the two an-
swers to be equal. For the modeling of crowded spec-
tra, the effects of line overlap will also need to be con-
sidered, for instance in the ‘all or nothing’ approach
(Cesaroni & Walmsley 1991). Such spectra will be rou-
tinely observed with the superb resolution and sensitiv-
ity of ALMA.
Our program is free for anybody to use for science, pro-
vided that appropriate reference is made to this paper.
For any other purpose such as to incorporate the pro-
gram into other packages which may be distributed to
the public, prior agreement with the authors is needed.
Acknowledgements. The authors wish to thank Huib Jan van
Langevelde for his efforts in documenting RADEX, and Erik
Deul for computing support at Leiden Observatory. JHB and
FLS acknowledge the Swedish Research Council for financial
support. FvdT and EvD thank the Netherlands Organization
for Scientific Research (NWO) and the Netherlands Research
School for Astronomy (NOVA). Finally we thank Volker
Ossenkopf, Marco Spaans, and an anonymous referee for
helpful comments on the manuscript.
References
Bernes, C. 1979, A&A, 73, 67
Black, J. H. 1994, in ASP Conf. Ser. 58: The First
Symposium on the Infrared Cirrus and Diffuse
Interstellar Clouds, ed. R. M. Cutri & W. B. Latter,
Black, J. H. 1998, in Chemistry and Physics of
Molecules and Grains in Space. Faraday Discussions
No. 109, 257
Black, J. H. 2000, in IAU Symposium 197 –
Astrochemistry: From Molecular Clouds to Planetary
Systems, ed. Y. C. Minh & E. F. van Dishoeck, 81
Black, J. H. & van Dishoeck, E. F. 1987, ApJ, 322, 412
Blake, G. A., Sutton, E. C., Masson, C. R., & Phillips,
T. G. 1987, ApJ, 315, 621
Blake, G. A., van Dishoek, E. F., Jansen, D. J.,
Groesbeck, T. D., & Mundy, L. G. 1994, ApJ, 428,
Burton, M. G., Hollenbach, D. J., & Tielens, A. G. G. M.
1990, ApJ, 365, 620
Cannon, C. J. 1985, The transfer of spectral line radia-
tion (Cambridge: University Press)
Carroll, T. J. & Goldsmith, P. F. 1981, ApJ, 245, 891
Castor, J. I. 1970, MNRAS, 149, 111
Cesaroni, R. & Walmsley, C. M. 1991, A&A, 241, 537
Comito, C., Schilke, P., Phillips, T. G., et al. 2005, ApJS,
156, 127
Daniel, F., Cernicharo, J., & Dubernet, M.-L. 2006, ApJ,
648, 461
De Jong, T., Boland, W., & Dalgarno, A. 1980, A&A,
91, 68
De Jong, T., Dalgarno, A., & Chu, S.-I. 1975, ApJ, 199,
Dubernet, M. L. 2005, in IAU Symposium, ed. D. C.
Lis, G. A. Blake, & E. Herbst, 235
Elitzur, M. 1992, Astronomical masers (Kluwer
Academic Publishers)
Elitzur, M. & Asensio Ramos, A. 2006, MNRAS, 365,
Evans, II, N. J., Lee, J.-E., Rawlings, J. M. C., & Choi,
M. 2005, ApJ, 626, 919
Expósito, J. P. F., Agúndez, M., Tercero, B., Pardo, J. R.,
& Cernicharo, J. 2006, ApJ, 646, L127
Ferland, G. J. 2003, ARA&A, 41, 517
Fixsen, D. J., Bennett, C. L., & Mather, J. C. 1999, ApJ,
526, 207
Fixsen, D. J. & Mather, J. C. 2002, ApJ, 581, 817
Flower, D. R. 1989, Physics Reports, 174, 1
Fuente, A., Black, J. H., Martı́n-Pintado, J., et al. 2000,
ApJ, 545, L113
Genzel, R. 1991, in NATO ASIC Proc. 342: The Physics
of Star Formation and Early Stellar Evolution, ed.
C. J. Lada & N. D. Kylafis, 155
Goicoechea, J. R., Pety, J., Gerin, M., et al. 2006, A&A,
456, 565
Goldreich, P. & Scoville, N. 1976, ApJ, 205, 144
Goldsmith, P. F. & Langer, W. D. 1999, ApJ, 517, 209
Gray, M. D. & Field, D. 1995, A&A, 298, 243
Gredel, R., Lepp, S., Dalgarno, A., & Herbst, E. 1989,
ApJ, 347, 289
Hauschildt, H., Güsten, R., Phillips, T. G., et al. 1993,
A&A, 273, L23
Helmich, F. P., Jansen, D. J., de Graauw, T., Groesbeck,
T. D., & van Dishoeck, E. F. 1994, A&A, 283, 626
Helmich, F. P. & van Dishoeck, E. F. 1997, A&AS, 124,
Hogerheijde, M. R., Jansen, D. J., & van Dishoeck, E. F.
1995, A&A, 294, 792
Hogerheijde, M. R. & van der Tak, F. F. S. 2000, A&A,
362, 697
Jakob, H., Kramer, C., Simon, R., et al. 2007, A&A,
461, 999
Jansen, D. J. 1995, Ph.D. Thesis, Leiden University
Jansen, D. J., van Dishoeck, E. F., & Black, J. H. 1994,
A&A, 282, 605
Johnstone, D., Boonman, A. M. S., & van Dishoeck,
E. F. 2003, A&A, 412, 157
Leung, C.-M. & Liszt, H. S. 1976, ApJ, 208, 732
Leurini, S., Schilke, P., Menten, K. M., et al. 2004,
A&A, 422, 573
Mangum, J. G. & Wootten, A. 1993, ApJS, 89, 123
Maret, S., Ceccarelli, C., Tielens, A. G. G. M., et al.
2005, A&A, 442, 527
Mathis, J. S., Mezger, P. G., & Panagia, N. 1983, A&A,
128, 212
Mihalas, D. 1978, Stellar atmospheres (2nd edition)
(San Francisco, W. H. Freeman and Co.)
Müller, H. S. P., Thorwirth, S., Roth, D. A., &
Winnewisser, G. 2001, A&A, 370, L49
Ossenkopf, V., Trojan, C., & Stutzki, J. 2001, A&A,
378, 608
Osterbrock, D. E. & Ferland, G. J. 2006, Astrophysics
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra 11
Poelman, D. R. & Spaans, M. 2005, A&A, 440, 559
Prasad, S. S. & Tarafdar, S. P. 1983, ApJ, 267, 603
Roueff, E., Felenbok, P., Black, J. H., & Gry, C. 2002,
A&A, 384, 629
Rybicki, G. B. & Lightman, A. P. 1979, Radiative
processes in astrophysics (New York, Wiley-
Interscience)
Schöier, F. L. 2000, Ph.D. Thesis, Stockholm University
Schöier, F. L., Jørgensen, J. K., van Dishoeck, E. F., &
Blake, G. A. 2002, A&A, 390, 1001
Schöier, F. L., van der Tak, F. F. S., van Dishoeck, E. F.,
& Black, J. H. 2005, A&A, 432, 369
Sobolev, V. 1960, Moving envelopes of stars (Harvard
University Press)
Spaans, M. & van Langevelde, H. J. 1992, MNRAS,
258, 159
Stutzki, J. & Winnewisser, G. 1985, A&A, 144, 13
Tafalla, M., Myers, P. C., Caselli, P., Walmsley, C. M.,
& Comito, C. 2002, ApJ, 569, 815
Takahashi, J. & Uehara, H. 2001, ApJ, 561, 843
Takahashi, T., Silk, J., & Hollenbach, D. J. 1983, ApJ,
275, 145
Van der Tak, F., Neufeld, D., Yates, J., et al. 2005, in
The Dusty and Molecular Universe: A Prelude to
Herschel and ALMA, ed. A. Wilson, 431–432
Van der Tak, F. F. S. 2005, in IAU Symposium
227: Massive Star Birth, ed. R. Cesaroni, M. Felli,
E. Churchwell, & M. Walmsley (Cambridge:
University Press), 70–79
Van der Tak, F. F. S. 2006, Phil. Trans. R. Soc. Lond.,
364, 3101
Van der Tak, F. F. S., van Dishoeck, E. F., & Caselli, P.
2000a, A&A, 361, 327
Van der Tak, F. F. S., van Dishoeck, E. F., Evans, II,
N. J., Bakker, E. J., & Blake, G. A. 1999, ApJ, 522,
Van der Tak, F. F. S., van Dishoeck, E. F., Evans, II,
N. J., & Blake, G. A. 2000b, ApJ, 537, 283
Van Dishoeck, E. F., Blake, G. A., Jansen, D. J., &
Groesbeck, T. D. 1995, ApJ, 447, 760
Van Dishoeck, E. F. & Hogerheijde, M. R. 1999, in
NATO ASIC Proc. 540: The Origin of Stars and
Planetary Systems, ed. C. J. Lada & N. D. Kylafis,
Van Zadelhoff, G.-J., Dullemond, C. P., van der Tak,
F. F. S., et al. 2002, A&A, 395, 373
Wright, E. L., Mather, J. C., Bennett, C. L., et al. 1991,
ApJ, 381, 200
Yates, J. A., Field, D., & Gray, M. D. 1997, MNRAS,
285, 303
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra, Online Material p 1
Online Material
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra, Online Material p 2
Appendix A: Program input and output
A.1. Program input
The input parameters to RADEX are the following:
1. The name of the molecular data file to be used.
2. The name of the file to write the output to.
3. The frequency range for the output file [GHz]. All
transitions from the molecular data file are always
taken into account in the calculation, but often it is
practical to write only a limited set of lines to the
output.
4. The kinetic temperature of the cloud [K].
5. The number of collision partners to be used. Most
users will want H2 as only collision partner, but in
more specialized cases, additional collisions with
H or electrons may for instance play a role. See
the molecular datafiles for details. For some species
(CO, atoms) separate collision data for ortho and
para H2 exist; the program then uses the thermal or-
tho/para ratio unless the user specifies otherwise.
6. The name (case-insensitive) and the density [cm−3]
of each collision partner. Possibilities are H2, p-H2,
o-H2, electrons, atomic H, He, and H
7. The temperature of the background radiation field
– If >0, a black body at this temperature is used.
Most users will adopt the cosmic microwave
background at TCMB=2.725(1+z) K for a galaxy
at redshift z.
– If =0, the average interstellar radiation field
(ISRF) is used, taken from Black (1994) with
modifications described in § 3.2. This spectrum
is not adjustable by a scale factor because it con-
sists of several components that are not expected
to scale linearly with respect to each other.
– If <0, a user-defined radiation field is used, spec-
ified by values of frequency [cm−1], intensity
[Jy nsr−1], and dilution factor [dimensionless].
Spline interpolation and extrapolation are ap-
plied to this table. The intensity need not be
specified at all frequencies of the line list, but
a warning message will appear if extrapolation
(rather than interpolation) is required.
8. The column density of the molecule [cm−2].
9. The FWHM line width [km s−1].
A.2. Program output
The output file written by RADEX first replicates the in-
put parameters, and then lists the following quantities
for each spectral line within the specified frequency
range.
1. Quantum numbers, upper state energy [K], fre-
quency [GHz], and wavelength [µm]. These num-
bers are just copied from the molecular data file,
which usually comes from the LAMDA database.
Frequencies from this database are generally of
line catalogs such as CDMS (Müller et al. 2001)12
should be consulted.
2. The excitation temperature [K] as defined in
Eq. (10). In general, different lines have different
excitation temperatures. Lines are thermalized if
Tex=Tkin; in LTE, all lines are thermalized.
3. The line optical depth, defined as the optical depth
of the equivalent rectangular line shape (φν = 1/∆ν).
4. The line intensity, defined as the Rayleigh-Jeans
equivalent temperature TR [K].
5. The line flux, defined as the velocity-integrated
intensity, both in units of K km s−1 (common in
radio astronomy) and of erg cm−2 s−1 (common
in infrared astronomy). The line flux is calcu-
lated as 1.0645TR∆V, where the factor 1.0645 =√
ln 2) converts the adopted rectangular line
profile into a Gaussian profile with an FWHM of
∆V. The integrated profile is useful to estimate the
total emission in the line, but it has limited mean-
ing at high optical depths, because the change of
optical depth over the line profile is not taken into
account. Proper modeling of optically thick lines re-
quires programs that resolve the source both spec-
trally and spatially (see § 4.1.2 for further discus-
sion).
Auxiliary output files can be generated, for example to
display the adopted continuum spectrum.
Appendix B: Coding standards
The original version of RADEX was written in such a
way as to minimize the use of machine memory which
was expensive until a decade ago. Nowadays, clarity
and easy maintenance are more important requirements,
which is why the source code has been re-written fol-
lowing the rules below. We hope that these rules will be
useful for the development of other ‘open source’ as-
tronomical software. For further guidelines on scientific
programming we recommend the Software Carpentry13
on-line course.
1. All the action is in subroutines; the sole purpose of
the main program is to show the structure of the pro-
gram. The subroutines are grouped into several files
for a better overview; compilation instructions for
automated builds on a variety of platforms are in a
Makefile.
2. The program text is interspersed with comments at
a ratio of ≈1:1. In particular, each subroutine starts
with a description of its contents, its input and out-
put, where incoming calls come from, and which
calls go out. Then the properties of each variable are
described: contents, units and type.
3. Variables and subroutines have descriptive names
with a length of 5–10 letters. Names of integer vari-
ables start with the letters i..n; names of floating-
point variables (always of double precision) with
a..h or o..z. There are no specific namings for vari-
ables of character or logical type, as for example
http://cdms.de
http://www.swc.scipy.org/
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra, Online Material p 3
in CLOUDY. Names are always based on the English
language. We do not use upper-, lower-, or camel-
case to distinguish types of identifiers, as some pro-
grams do.
4. Loops are marked by indenting the program text.
The loop variables are always called ilev, iline ..,
never just i.
5. Subroutines start with a check whether the input pa-
rameters have reasonable values. Such checks force
soft landings if necessary, and avoid runtime errors.
6. Statements with calculations use spaces around
the =, + and – symbols, but not around others.
Calculations that consist of multiple steps are split
over as many program lines. Multiple assignments
in a row are aligned at the = sign.
7. The program text avoids “magic numbers” both in
calculations and in definitions. Numbers that are of-
ten used such as Planck’s constant are defined at one
central place in the program. Similarly, often-used
variables such as physical parameters are stored in
shared memory rather than passed on via subroutine
calls.
Appendix C: Diagnostic plots of molecular
line ratios
We have used the RADEX program and the LAMDA
database to calculate line ratios of several commonly
observed molecules for a range of kinetic tempera-
tures and H2 densities. The plots in this Appendix may
be used by observers to estimate physical conditions
from their data. Line ratios have the advantage of be-
ing less sensitive to calibration errors than absolute line
strengths, especially when the two lines have been mea-
sured with the same telescope, receiver and spectrome-
The calculations assume a column density of 1012 cm−2,
a line width of 1.0 km s−1, and use a 2.73 K blackbody
as background radiation field. Under these conditions,
the lines are optically thin, so that the line ratios do not
depend on column density. The calculations also assume
that the emission in both lines fills the telescope beams
equally, which may be the case if the lines are close
in frequency. However, lines are generally measured in
beams of different sizes, and the observations need to
be corrected to account for this effect, if the source is
known to be compact.
Linear molecules such as CO (Fig. C.1) are tracers
of density at low densities, when collisions compete
with radiative decay. At higher densities, the excitation
becomes thermalized and the line ratios are sensitive
to temperature. For a given molecule, moving up the
J−ladder means probing higher temperatures and densi-
ties. Note that for the column densities of typical dense
interstellar clouds, the CO lines are optically thick, and
observations of 13CO or even rarer isotopologues must
be used to probe physical conditions.
The critical densities of molecular lines scale as µ2ν3,
where µ is the permanent dipole moment of the
molecule and ν is the frequency of the line. Indeed, the
CS molecule (Fig. C.2) has a larger dipole moment than
CO, and its line ratios are mainly probes of the den-
sity. The small frequency spacing between the lines of
CS makes this molecule very useful to probe density
structure (e.g., Van der Tak et al. 2000b). The HCO+
and HCN molecules (Fig. C.3) display similar trends to
CS, although their line spacing is not as small.
Non-linear molecules such as H2CO (Figs. C.4 and C.5)
have the advantage that both temperature and density
may be probed within the same frequency range. Ratios
of lines from different J−states tend to be density tracers
(left panels), while ratios of lines from the same J−state
but different K−states are mostly temperature probes
(right panels). The lines of H2CO are often quite strong,
making this molecule a favourite tracer of temperature
and density (Mangum & Wootten 1993). Other asym-
metric molecules have also been used, such as H2CS
(Blake et al. 1994) and CH3OH (Leurini et al. 2004) al-
though abundance variations from source to source or
even within sources often complicate the interpretation
(Johnstone et al. 2003).
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra, Online Material p 4
Appendix D: The Python scripts
The RADEX distribution comes with two scripts,
radex line.py and radex grid.py, to automate
standard modeling procedures. The scripts are written
in Python and are run from the Unix shell command
line after manual editing of parameters.
The first script, radex line.py, calculates the column
density of a molecule from an observed line intensity,
given estimates of kinetic temperature and H2 volume
density. The input parameters are:
1. The kinetic temperature [K]
2. The number density of H2 molecules [cm
3. The temperature of the background radiation field
[K], usually 2.73 (CMB).
4. The name of the molecule (or molecular data file)
5. The frequency of the line [GHz]
6. The observed line intensity [K]
7. The observed line width [km s−1]
Furthermore, two numerical parameters have good de-
fault values but will need to be changed occasionally:
1. The free spectral range around the line (default
10%): this number must be smaller for molecules
with many line close in frequency, such as CH3OH.
The program uses this parameter to find the ob-
served line from the list of lines in the molecular
model.
2. The required accuracy (default 10%): The default
corresponds to the calibration uncertainty of most
telescopes.
The script iterates on column density until the observed
and modeled line fluxes agree to within the desired ac-
curacy. The best-fit column density is directly written
to the screen. The file radex.out gives details of the
best-fit model.
The second script, radex grid.py, runs a series of
RADEX models to estimate the kinetic temperature
and/or the volume density from an observed line ratio.
The user needs to set the following input parameters:
1. The grid boundaries: minimum and maximum ki-
netic temperature [K] and minimum and maximum
H2 volume density [cm
2. The temperature of the background radiation field
[K], usually 2.73 (CMB).
3. The molecular column density [cm−2]. For the
illustrative plots in Appendix C, a low value
(N=1012 cm−2) was used, so that the line ratios are
independent of column density (optically thin limit).
However, in modeling specific observations, it is
worth varying this parameter to assess the sensitivity
of the line ratio to column density.
4. The observed line width [km s−1], usually an aver-
age of the widths of the two lines.
1. The number of grid points along the temperature and
density axes.
2. The free spectral range around the line (see above)
The name of the molecule and a list of observed line
ratios and names of associated output files are given at
the start of the main program. The script produces a file
radex.out which is a tabular listing of temperature,
log density, and line ratio. This results may be plotted
with the user’s favourite plotting program.
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra, Online Material p 5
Fig. C.1. Line ratios of CO in the optically thin limit
as a function of kinetic temperature and H2 density.
Contours are spaced linearly and some contours are la-
beled for easy identification.
Fig. C.2. Line ratios of CS in the optically thin limit
as a function of kinetic temperature and H2 density.
Contours are spaced linearly and some contours are la-
beled for easy identification.
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra, Online Material p 6
Fig. C.3. Line ratios of HCO+ and HCN in the optically
thin limit as a function of kinetic temperature and H2
density. Contours are spaced linearly and some contours
are labeled for easy identification.
Fig. C.4. Line ratios of o-H2CO in the optically thin
limit as a function of kinetic temperature and H2 den-
sity. Contours are spaced linearly and some contours are
labeled for easy identification.
Van der Tak et al.: Fast non-LTE analysis of interstellar line spectra, Online Material p 7
Fig. C.5. Line ratios of p-H2CO in the optically thin
limit as a function of kinetic temperature and H2 den-
sity. Contours are spaced linearly and some contours are
labeled for easy identification.
	Introduction
	Radiative transfer and molecular excitation
	Basic formalism
	Molecular line cooling
	Escape probability
	The program
	Basic capabilities
	Background radiation field
	Chemical formation and destruction rates
	Calculation
	Results
	Limitations of the program
	Comparison with other methods
	Homogeneous models
	The case of CO
	The case of HCO+
	Observations of a Young Stellar Object
	Conclusions
	Program input and output
	Program input
	Program output
	Coding standards
	Diagnostic plots of molecular line ratios
	The Python scripts
