\documentclass[preprint,12pt]{elsarticle}
%
\journal{J. Quant. Spectrosc. Radiat. Transfer}
%
\bibliographystyle{unsrt}


%- Extra packages
%
\usepackage{url}
%\usepackage{listings}
%\usepackage{color}
\usepackage{lineno}  

%- Special commands
%
\newcommand{\todo}[1]{[* #1 *]}


%- Symbol definitions
%
\def\ARTS{ARTS}       % To easily change how we write ARTS
\def\V1{ARTS-1}       % To easily change how we write ARTS-1.0.x
%
\def\AUG{AUG}         % ARTS user guide
\def\ADG{ADG}         % ARTS developer guide       
\def\ATD{ATD}         % ARTS theory document       
%
\newcommand{\unit}[1]{{\it#1}}
%
% Use the symbol definitions from AUG as far as possible
\input{symbol_defs}
\renewcommand{\PDir}  {\ensuremath{{\bf \hat{n}}}}
\newcommand{\PPos}    {\ensuremath{{\bf r}}}
\newcommand{\Iv}      {\ensuremath{I_v}}
\newcommand{\Ih}      {\ensuremath{I_h}}
\newcommand{\Ipff}    {\ensuremath{I_{+45^\circ}}}
\newcommand{\Imff}    {\ensuremath{I_{-45^\circ}}}
\newcommand{\Irhc}    {\ensuremath{I_{rhc}}}
\newcommand{\Ilhc}    {\ensuremath{I_{lhc}}}
\newcommand{\BT}      {\aTmp{B}}




\begin{document}

\begin{frontmatter}

%--- Titel + authors ----------------------------------------------------------
\title{ARTS, the Atmospheric Radiative Transfer Simulator, Version 2}

\author[chalmers]{P. Eriksson\corref{cors}}
\ead{patrick.eriksson@chalmers.se}

\author[ltu]{S. A. Buehler\corref{cors}}
\ead{sbuehler@ltu.se}

\author[metservice]{C. P. Davis}
\author[lmu]{C. Emde}
\author[ltu]{O. Lemke}

\address[chalmers]{%
  Department of Earth and Space Sciences,
  Chalmers University of Technology,\\
  SE-41296 Gothenburg,
  Sweden
}
\address[ltu]{%
  Department of  Computer Science, Electrical and Space Engineering,
  Lule\aa\ University of Technology,
  Space Campus 1,
  SE-98128 Kiruna,
  Sweden
}
\address[metservice]{%
  Meteorological Service of New Zealand Ltd.,
  Wellington 6012,
  New Zealand
}

\address[lmu]{%
  Meteorological Institute,
  Ludwig-Maximilians-Universit\"at,
  Theresienstr. 37,\\
  D-80333 Munich,
  Germany
}
  
\cortext[cors]{Corresponding authors}




%--- Abstract + keywords -----------------------------------------------------
\begin{abstract}
  The second version of the Atmospheric Radiative Transfer Simulator, \ARTS, is
  introduced. This is a general software package for long wavelength radiative
  transfer simulations, with a focus on passive microwave observations. The
  core part provides a workspace environment, in line with script languages.
  New for this version is an agenda mechanism that gives a high degree of
  modularity. The framework is intended to be as general as possible: the
  polarisation state can be fully described, the model atmosphere can be one-
  (1D), two- (2D) or three-dimensional (3D), a full description of geoid and
  surface is possible, observation geometries from the ground, from satellite,
  and from aeroplane or balloon are handled, and surface reflection can be
  treated in simple or complex manners. Remote sensing applications are
  supported by a comprehensive and efficient treatment of sensor
  characteristics. Jacobians can be calculated for the most important
  atmospheric variables in non-scattering conditions. Finally, the most
  prominent feature is the rigorous treatment of scattering that has been
  implemented in two modules: A discrete ordinate iterative approach mainly
  used for 1D atmospheres, and a Monte Carlo approach which is the preferred
  algorithm for 3D atmospheres. \ARTS\ is freely available, and maintained as
  an open source project.
\end{abstract}

\begin{keyword}
%% keywords here, in the form: keyword \sep keyword
atmospheric radiative transfer \sep
emission \sep
polarisation \sep
scattering \sep
sensor modelling \sep
remote sensing

%% PACS codes here, in the form: \PACS code \sep code
\PACS 42.68.Ay \sep 42.68.Mj \sep 42.68.Wt \sep 92.60.Ta
\end{keyword}
\end{frontmatter}


\linenumbers
%\setpagewiselinenumbers

\section{Introduction}
%-----------------------------------------------------------------------------
\label{sec:introduction}
%
\ARTS\ is a free open-source software program that simulates atmospheric
radiative transfer. It focuses on thermal radiation from the microwave to the
infrared spectral range. Version 1.0 of \ARTS\ \citep{buehler:artst:05}, which
handles cases without scattering, was mainly developed between 2000 and 2005.
It is well validated \citep{melsheimer:inter:05,
  buehler06:_radiat_trans_calcul_passiv_microw_satel_sensor_jgr,
  saunders07:_atmos_infrar_sound_airs_jgr} and still used, primarily for the
analysis of ground-based and satellite-based measurements in the
millimetre/submillimetre spectral region \citep[e.g.,][]{Urban20051, john05:_station_comp_acp, john04:_amsu_ozone_grl}. A large part
of its popularity is due to the retrieval software Qpack
\citep{eriksson:qpack:05}, which uses \ARTS\ as the forward model. But \ARTS\
version 1.0 has also been used for the simulation of clear-sky broadband
energy fluxes in the thermal infrared spectral range \citep{buehler:recen:06,
  john:under:06}. This model version is below denoted as \V1.

From 2002, the \ARTS\ developer community became increasingly interested in the
remote sensing of clouds, particularly ice clouds. A main driver was the ESA
mission proposal CIWSIR \citep{buehler:aconc:07}, a submillimetre instrument
for the characterisation of ice clouds, which required a radiative transfer
model that could simulate the scattering by ice particles, including
polarisation effects
\citep{miao03:_potential_of_polarization_measurements,davis07:_3d_pol_sim_acp}.

Another strong driver was the treatment of microwave limb sounders: firstly for
the analysis of cloud-affected data from the MLS and Odin-SMR satellite
instruments
\citep{davis07:_3d_pol_sim_acp,davis05:_cirrus_induc_polar_ghz_aura,
  eriksson:first:07,rydberg:nonga:09}. Secondly, future limb sensors will
sample the atmosphere more densely in order to increase the `tomographic'
capability. This and the scattering by clouds demand going beyond a 1D
representation of the atmosphere.

The interest in such atmospheric sounding techniques led to an internal fork in
the \ARTS\ program development. Active development shifted to version
1.1.x, which included modules to simulate scattering by cloud particles and
other significant improvements, while \V1\
was maintained to provide a stable version for existing users. The new version
with scattering is now complete and stable enough to fully replace the old
version. We mark this by calling the
latest version \ARTS\ 2.0. The purpose of this article is to present \ARTS\
2.0, and give an overview of its features, strengths, and limitations. In the
remaining text ``\ARTS'' refers to the latest 2.0 version.


\section{Overview}
%-----------------------------------------------------------------------------
\label{sec:overview}

\subsection{Scope}
%
The ambition is to accommodate simulations of any type of passive longwave
observation, and \ARTS\ is designed to have no limitations when it comes
to the representation of polarisation state, atmospheric fields and geometrical
aspects: 
\begin{enumerate}
\item The Stokes formalism is used to describe polarisation
  (Sec.~\ref{sec:rtvar}).
\item The model atmosphere can be represented with a one (1D), two (2D) or
  three (3D) dimensional view (Sec.~\ref{sec:modelatm}).
\item No assumption of a ``flat Earth'', the geoid and the surface are either
  spherical (by definition for 1D), or can be given arbitrary shapes
  (Sec.~\ref{sec:modelatm}).
\item Radiative transfer calculations can be made from any position and along
  any direction, as long as the resulting calculations make sense with respect
  to the model atmosphere (Sec.~\ref{sec:obsgeom}).
\end{enumerate}
Individual functions can be limited
to some configurations, for example, the Monte Carlo scattering module
(Sec.~\ref{sec:scattering}) is restricted to 3D.

As mentioned, handling of scattering is a primary aim of \ARTS,
where the goal is to allow arbitrary complex scattering properties. This goal
has been reached for surface reflection (Sec.~\ref{sec:surfscat}), but not
completely for particle scattering (Sec.~\ref{sec:optprop}). The development
has so far focused on exact algorithms and the model's strongest side is that
complicated simulations can be performed in a stringent manner. \ARTS\ is thus
primarily a research tool. Speed has not been a primary objective, and
extremely fast, but approximative, algorithms like RTTOV
\cite{saunders99:_rttov_qjrms} are not in the scope of \ARTS.

The software is mainly developed for remote sensing applications, and an
extensive support for inclusion of sensor characteristics is provided
(Sec.~\ref{sec:sensor}). In addition, weighting functions (also called
Jacobians) \cite{rodgers90,rodgers:inver:00} can be obtained for a number of
variables in non-scattering conditions.

\ARTS\ comes with a small amount of input data. The purpose of these data is to
provide some usage examples, and allow the developers to perform standardised
tests before committing changes of the code. Normally, the user has to provide
the bulk of input data, such as temperatures, volume mixing ratios and
spectroscopic parameters. A noticeable exception is that a number of
``absorption models'' are built into the model (Sec.~\ref{sec:abstable}).
%Handling of line coupling requires dedicated absorption models.


\subsection{Documentation}
\label{sec:doc}
%
The efforts to document \ARTS\ focus on the practical usage of the software.
This is mainly achieved through the built-in documentation, that provides a
definition and a basic description of individual variables and methods. This
documentation can be browsed on-line at \url{www.sat.ltu.se/arts/docserver}. An
introduction to the usage of \ARTS\ on a system level is given by some example
cases that are distributed along with the source code.

This article provides a compact overview of ARTS. A more detailed description
can be found in the three documents of guide type that are distributed with
\ARTS. Model definitions and algorithms are described in the `\ARTS\ user
guide' (\AUG). The `\ARTS\ development guide' (\ADG) gives practical
information for the source code. Background theory for some core subjects is
provided by the `\ARTS\ theory document' (\ATD). Some parts are described
further, or solely, by dedicated research articles
\cite{emde:apoli:04,davis05:_monte_carlo_i3etgrs,eriksson:effic:06,buehler:effic:10}.
See further \url{www.sat.ltu.se/arts/docs/}. Download instructions and
technical specifications are found at \url{www.sat.ltu.se/arts/getarts/}.


\subsection{Supporting tools}
%
Functions for creating input files and for reading output files (for e.g.\
plotting) are provided for two popular higher-level and interactive
environments, Python through PyARTS and Matlab through Atmlab. These packages
provide also additional features. For example, PyARTS allows the calculation of
particle optical properties using the T-matrix code by
\citep{mishchenko02:_scatt_absor_emiss_light_small_partic} and a new version of
Qpack is being implemented inside Atmlab. The packages can be downloaded
from the \ARTS\ home page and are documented separately.




\section{The workspace}
%-----------------------------------------------------------------------------
\label{sec:concept}
%
\subsection{\ARTS\ as a scripting language}
%
One of the main goals in the \ARTS\ development was to make the program as
flexible as possible, so that it can be used for a wide range of applications
and new features can be added in a relatively simple manner. As a result,
\ARTS\ behaves like a scripting language. An \ARTS\ controlfile contains a
sequence of instructions. When \ARTS\ is executed, the controlfile is parsed,
and then the instructions are executed sequentially. 

This feature is build around the ``workspace'' \citep{buehler:artst:05}. The
basic structure is unchanged from \V1, but some improvements
have been introduced. Regarding the ``workspace variables'', the set of
variables is now not fixed. The user is free to specify new variables, as part
of the controlfile operations. User-defined variables can
replace any of the pre-defined variables, as long as they are of the same type.

The syntax around the ``workspace methods'' is also somewhat changed. This change is
not described here, it should be clear from the on-line documentation and the
example cases (Sec.~\ref{sec:doc}). 
%For example, the on-line documentation of
%the method `\texttt{WriteXML}' is obtained through the command line switch
%`\texttt{-d}', as `\texttt{arts -d WriteXML}'. 
% [FIXME: I commented this out, because I think it is redundant. Also,
% better emphasize the HTML version, as in Section \ref{sec:doc}. - Stefan]

\subsection{Agendas}
\label{sec:agendas}
%
It became increasingly clear that the workspace methods alone do not provide
the flexibility sought. In order to avoid increasingly complex methods, the
concept of agendas was introduced. 
An agenda is a user-defined list of workspace methods, which are executed in sequence to calculate a predefined set of workspace variables from a predefined set of input (workspace) variables.  
As an example, the absorption is handled by an agenda. Several radiative
transfer methods use this agenda as an input variable. When they need local
absorption coefficients for a point in the atmosphere, they execute the agenda
with the local pressure, temperature, and trace gas volume mixing ratio values
as inputs. The agenda then provides absorption coefficients as output. If the
absorption is extracted from a pre-calculated lookup table or is calculated
from basic spectroscopic data (Sec.~\ref{sec:abstable}) depends on
which methods the user has elected to include in the agenda.

\section{Model definitions and input}
%-----------------------------------------------------------------------------
\label{sec:input}
%
This section gives some basic model definitions, and comments on mandatory
and other input of general character required for a model run. Units used for
ARTS specific input and output files, as well as internal definitions of
variables, follow the SI system.

\subsection{Radiative transfer, nomenclature and variables}
\label{sec:rtvar}
%
\ARTS\ describes (spectral) radiances using the Stokes vector, \StoVec. The
calculations can be selected to treat one to four elements of the Stokes
vector, all methods adjust automatically to this choice. The phrase ``scalar
radiative transfer'' refers to the case when just the first Stokes vector
element is considered. The other options are all termed as vector radiative
transfer. The four elements of the Stokes vector, $\StoVec=[I,Q,U,V]^T$, are
defined as:
\begin{eqnarray}
  I &=&   \Iv + \Ih = \Ipff + \Imff = \Ilhc + \Irhc,  \label{eq:I}\\
  Q &=&   \Iv - \Ih,                                  \label{eq:Q}\\
  U &=&   \Ipff - \Imff,                              \label{eq:U}\\
  V &=&   \Ilhc - \Irhc,                              \label{eq:V}
\end{eqnarray}
where \Iv, \Ih, \Ipff, and \Imff\ are the intensity of the component linearly
polarised at the vertical, horizontal, +45\degree\ and -45\degree\ direction,
respectively, and \Irhc, and \Ilhc\ are the intensity for the right- and
left-hand circular components. The definition used here follows
\cite{mishchenko02:_scatt_absor_emiss_light_small_partic}, see also \ATD.

Accordingly, $I$ is the total radiance and the other Stokes
elements give the difference between two orthogonal components. Individual
components are extracted as combinations of $I$ and the other elements, e.g.
\begin{equation}
  \label{eq:Iv}
  \Iv = (I+Q)/2.
\end{equation}

The standard vector radiative transfer equation (VRTE) for cases involving
multiple scattering is \cite{mishchenko02:_scatt_absor_emiss_light_small_partic}
 \begin{eqnarray}
    \label{eq:VRTE}
  \frac {\DiffD\StoVec(\Frq,\PPos,\PDir)}{\DiffD s} =
    &{}-\ExtMat(\Frq,\PPos,\PDir) \StoVec(\Frq,\PPos,\PDir) +
    \AbsVec(\Frq,\PPos,\PDir) \rlap{$B(\Frq,\PPos)$} \\ \nonumber
    &{}+ \int_{4\pi} \DiffD\PDir' \PhaMat(\Frq,\PPos,\PDir,\PDir')
    \StoVec(\Frq,\PPos,\PDir'),
\end{eqnarray}
where \Frq\ is frequency, \PPos\ represents the atmospheric position, \PDir\ is
the propagation direction (at \PPos), $s$ is distance along \PDir, $\ExtMat$ is
the extinction matrix, $\AbsVec$ is the absorption vector, $B$ is the Planck
function and $\PhaMat$ is the phase (or scattering) matrix. This equation
assumes local thermodynamic equilibrium and that the scattering events can be
treated as incoherent.

Equation~\ref{eq:VRTE}, or some simplified version of it, is solved, giving
simulated radiances. The inclusion of sensor characteristics requires that
radiative transfer calculations are performed for a set of monochromatic
frequencies, the frequency grid, and a number of pencil beams
(Sec.~\ref{sec:sensor}). The frequency grid is a primary input variable; it
determines the frequencies for which absorption and radiative transfer are
calculated. The propagation through the atmosphere of the unscattered,
but possibly refracted, pencil beam is below denoted as the propagation path.



\subsection{Atmospheric and surface variables}
\label{sec:modelatm}
%
Atmospheric quantities can be defined to vary in one, two and three dimensions.
The atmospheric dimensionality can thus be 1D, 2D or 3D. Pressure is the
vertical coordinate in all cases. The two horizontal dimensions for 3D coincide
with standard latitude and longitude. The second dimension for 2D is for
simplicity denoted as latitude, but is not demanded to have a direct
geophysical interpretation. This “latitude” can, for example, represent the
angular distance inside the plane of a satellite orbit.

Each (active) atmospheric dimension has an associated grid. This gives an
atmospheric grid mesh, for which temperature, geometrical altitude (above the
geoid) and the volume mixing ratio for the species must be specified. The
basic definition of the model atmosphere is completed by the geoid radius and
the surface altitude, as a function of latitude and longitude. 

The minimum value of the pressure grid sets the upper limit of the model
atmosphere (vacuum assumed above). The lower limit for the calculation is set
by the ground, which constitutes a surface (with arbitrary topography) at the
boundary or inside the atmospheric domain. The atmosphere is undefined outside
the latitude and longitude grid ranges.


\subsection{Radiative transfer domains}
\label{sec:cloubox}
%
The default assumptions are that scattering can be neglected, and that
absorption and emission are unpolarised. More complicated calculation
conditions are restricted to a special domain of the atmosphere,
introduced initially to handle cloud scattering and consequently called
the ``cloudbox''.  For simplicity, the calculations outside the cloudbox are
denoted as ``clear sky''.

The vertical limits of the cloudbox are two pressure surfaces. For 1D, the
cloudbox extends around the model planet, as implied by the spherical symmetry
for this case. For higher atmospheric dimensions, the horizontal limits are
found at latitude and longitude grid points. The cloudbox can extend below the
surface, or be restricted to the atmosphere. The surface is allowed to cause
both scattering and polarisation effects outside the cloudbox.

\subsection{Particle optical properties}
\label{sec:optprop}
%
The optical properties of cloud droplets and ice crystals (\ExtMat, \AbsVec\ and
\PhaMat; see Eq.~\ref{eq:VRTE}) are required as input for scattering
calculations. They have to be pre-calculated outside the \ARTS\ program.

For liquid water clouds the droplets are in good approximation of spherical
shape and the optical properties can be computed using the well known Mie
theory \citep{mie1908}. The Atmlab toolbox includes a Mie program
\citep{maetzler2002} to generate optical properties of spherical particles. Ice
crystals have complex hexagonal shapes like solid columns, plates, aggregates
etc. The PyARTS package provides tools for the calculation of optical
properties of aspherical particles (cylinders, plates, and spheroids) which may
be used as approximations for the complex ice crystal shapes. PyARTS uses the
T-matrix codes by \citep{mishchenko02:_scatt_absor_emiss_light_small_partic}. 

\ARTS\ offers the possibility to define an arbitrary number of particle types.
For each particle type the user needs to define the particle number density
field, so that the desired mixture is obtained. Size and shape are not
specified specifically. Instead, each particle type is defined by its single
scattering properties. A common assumption is that aspherical cloud particles
are randomly oriented, this is one of the options in \ARTS. A special feature
of \ARTS\ is that it also allows to include oriented, more specifically
horizontally aligned, particles. Arbitrarily oriented particles can in
principle easily be implemented in \ARTS, but it is not done yet for the
practical reason that the optical properties for arbitrarily oriented particles
require a huge amount of computational memory. See further \AUG\ and
\cite{battaglia-2006-emiss}.


\subsection{Observation geometry}
\label{sec:obsgeom}
%
There are no intrinsic limitations for the observation geometry. Radiative
transfer can be performed for any position inside and above the model
atmosphere, and with arbitrary observation direction, as long as the radiative
transfer does not reach undefined parts of the atmosphere (Sec.\
\ref{sec:modelatm}). As long as this constraint is met, the observation position
can be outside the horizontal region covered by the latitude and
longitude grids. This option is useful for satellite limb sounding where the
distance between the sensor and the practical atmospheric entry point
can exceed 1500\,km. 

The observation geometry is defined by combinations of sensor position and
line-of-sight (LOS). The term sensor is used here, but this can be a
hypothetical instrument observing monochromatic radiances.
Inclusion of sensor characteristics is discussed in Sec.~\ref{sec:sensor}. The
vertical coordinate used for the sensor position is the radius (distance from
the origin). Horizontal position is defined by latitude and longitude.

The LOS is specified by a zenith angle, and for 3D also an azimuth angle. The
zenith angle is the angle between the observation direction and the radial unit
vector. This angle is inside the range [0\degree,180\degree]. For 2D, zenith
angles down to -180\degree\ are also defined, where a positive / negative value
signifies an observation direction towards higher / lower latitudes. The
azimuth angle is defined as the clockwise angle between the observation
direction and meridional plane north of the observation point. Westward
observations have negative azimuth angles and the allowed range is
[-180\degree,180\degree].

ARTS is designed to handle a complete measurement sequence by default, and the
involved variables can hold a series of position and LOS combinations. Each
combination of position and LOS is denoted as a measurement block. This
reflects that the operations for a single position and LOS combination can
involve numerous radiance calculations, and that the output can correspond
to several measurement spectra. A static sensor is assumed inside each
measurement block and any shift of the observation position requires a new
such block. See further Sec.~\ref{sec:sensor}.




\section{Calculation algorithms}
%-----------------------------------------------------------------------------
\label{sec:calc}
\subsection{Gas absorption}
\label{sec:abstable}
%
The actual gas absorption calculation routines in \ARTS\ are identical to those
in \V1\ \cite{buehler:artst:05}. In particular, \ARTS\ can do line-by-line
absorption calculations, but includes also some predefined complete absorption
models and continua. 
%The complete absorption models include the model for O$_2$
%absorption in the microwave spectral region, including line coupling, by
%\citep{rosenkranz93:_absorption_microwaves}. The pure continuum absorption
%models include the MT\_CKD infrared water vapour continuum by
%\citep{clough05:_atmospheric_jqsrt}. 
The absorption can be calculated explicitly for each position along the
propagation path, that gives highest possible accuracy but slow calculations.
As a more rapid alternative, a lookup table approach has been implemented,
which stores pre-calculated absorption cross-sections as a function of
pressure, temperature, and water vapour concentration \cite{buehler-absor}.


\subsection{Ray tracing}
\label{sec:ray}
%
The radiative transfer equation is solved along a pre-calculated propagation
path. Such a path is basically described by a set of positions and the distance
between these points. The ray inside each grid box is calculated separately.
There are two reasons for this. Firstly, the DOIT algorithm
(Sec.~\ref{sec:scattering}) operates only with such local propagation paths.
Secondly, interpolation tends to cause a smoothing of atmospheric structures
and to decrease this effect it is desirable that the propagation points include
all crossings with the atmospheric grids. A step-wise approach is then
required, considering that these points can not be calculated in an analytical
manner with refraction. The same applies for the crossings with pressure
surfaces, even without refraction, as \ARTS\ allows the radius for each surface
to vary. 
%Propagation path points are also inserted for the
%position of the sensor and the tangent point (when applicable), if these are
%found inside the model atmosphere. A maximum allowed distance between the grid
%points can also be specified, resulting in additional points being inserted to
%fulfil the criterion. The end points of complete paths match either the sensor
%position, the top of the atmosphere, the surface or a cloud box boundary.

The details of the path calculations are described in \AUG, and are not
repeated here. Refraction is so far only handled in a very straightforward, but
inefficient, way, and further work is needed on this point.


\subsection{Radiative transfer algorithms}
\label{sec:csrte}
%
\subsubsection{Clear-sky}
%
As described in Sec.~\ref{sec:cloubox}, the term ``clear sky'' refers in \ARTS\
to the domain outside the cloudbox. For this domain it is assumed that
scattering can be neglected and that absorption (and thus also emission) is
unpolarised. However, the radiative transfer must be performed in a vector
manner, to correctly propagate polarisation effects generated inside the
cloudbox and by the surface.

This part is totally reimplemented but the calculations are basically performed
as in \V1, including the calculation of weighting functions
\citep{buehler:artst:05}. As emission is unpolarised for this domain, only
transmission has to be considered for the $Q$, $U$ and $V$ elements of the
Stokes vector (i.e.\ the Beer-Lambert law). An analytical approach can be used
for the weighting functions of some variables, so far implemented in \ARTS\ for
gas concentrations and atmospheric temperatures (neglecting non-local effects
due to refraction and hydrostatic equilibrium). 
%Other variables are hardly
%treated analytically and pure numerical approaches are used, where the code has
%been restructured to apply the ``perturbations'' with a lower calculation cost.



\subsubsection{Cloud scattering}
\label{sec:scattering}
%
The most unique feature of \ARTS\ is the possibility to handle scattering in a
rigorous manner. In fact, two algorithms that solve the VRTE
(Eq.~\ref{eq:VRTE}) have been implemented as part of the development of \ARTS.
One of the algorithms is based on a Discrete Ordinate ITerative (DOIT) scheme
\cite{emde:apoli:04}. The second algorithm applies Monte Carlo (MC) integration
with importance sampling \cite{davis05:_monte_carlo_i3etgrs}. The DOIT scheme calculates the entire radiation field within the `cloudbox' and is the preferred method for 1-D calculations.  The MC scheme calculates the Stokes' Vector for only a single viewing direction, but, due to the efficiency of Monte Carlo methods for highly dimensioned integration, is the preferred method for 3-D calculations.  The MC algorithm is not confined to the cloudbox, and handles surface effects in a similar fashion to cloud scattering and emission.  Also, if desired, scalar antenna response functions can be handled by Monte Carlo integration over viewing directions.

DOIT and MC make use of the general features of \ARTS\ described in this article, and we refer
to \cite{emde:apoli:04,davis05:_monte_carlo_i3etgrs}, \AUG\ and \ATD\ for 
details of the specific algorithms. Both DOIT and MC have been applied for
theoretical studies, as well as practical retrievals, for example
\cite{hoepfneremde05:_ir_comparison,pietranera07:_cmb_polarization_mnras,buehler:aconc:07,rydberg:predi:07,eriksson:first:07,sreerekha2008}
and
\cite{davis05:_cirrus_induc_polar_ghz_aura,battaglia07:_microw_jqsrt,davis07:_3d_pol_sim_acp,adams2008,rydberg:nonga:09},
respectively.



\subsection{Surface scattering}
\label{sec:surfscat}
%
The Stokes vector for upwelling radiation from the surface, $\aStoVec{u}$, in
the direction of \PDir, is calculated as
\begin{eqnarray}
  \label{eq:surface}
  \aStoVec{u}(\Frq,\PDir) &=& \aStoVec{e}(\Frq,\PDir) + 
  \int_{2\pi} \DiffD\PDir' \MtrStl{R}(\Frq,\PDir,\PDir')\StoVec^d(\Frq,\PDir')
                                                                   \nonumber\\
  &\approx& \aStoVec{e}(\Frq,\PDir) + 
  \sum_{i=1}^n \MtrStl{R}_i(\Frq,\PDir,\PDir')\StoVec^d_i(\Frq,\PDir')
\end{eqnarray}
The first term, $\aStoVec{e}$, is surface emission for the direction of
concern. The second term treats the reflection of down-welling radiation, where we use a discrete approximation.
This equation can be compared to the last term of Eq.~\ref{eq:VRTE}, that
describes scattering into the line-of-sight. The main differences are that this
integration is performed only over a half sphere and \MtrStl{R}\ is
denoted as the bidirectional polarised reflectance distribution function (BPDF).

Accordingly, the down-welling radiation, $\StoVec^d$, is calculated for $n$
directions, giving $\StoVec^d_i$. The set of $\StoVec^d_i$ are weighted with
the matrices $\MtrStl{R}_i$ that are a combination of the BPDF and the solid
beam angle that each direction $i$ represents.

The down-welling term of Eq.~\ref{eq:surface} vanishes if the surface is
treated to be a blackbody. For surfaces that can be treated as lacking
roughness, $n$ is one and $\PDir'$ is the specular direction. The required
value for $n$ and the best selection of the $\PDir'$-directions for other
situations is open for experimentation. Methods for blackbody, specular and
Lambertian surface conditions have been implemented (applied equations found in
\AUG).

As noted in Sec.~\ref{sec:scattering}, the MC algorithm has its own way of handling surface scattering: using Monte Carlo integration to evaluate the integral in Eq.~\ref{eq:surface}.


\subsection{Sensor characteristics}
\label{sec:sensor}
%
Several instrumental effects can be expressed as 
\begin{equation}
  \label{eq:sensorint}
  \int r(x) I(x)\, \DiffD x,
\end{equation}
where $r$ is the instrument's response function, $I$ is the radiance
and $x$ is frequency or some other variable, depending on which
response that is treated. The normal case is that simulations are repeated for
the same sensor characteristics, and a direct implementation of
Eq.~\ref{eq:sensorint} is normally not most efficient. In practise, $I$ is a
discrete quantity, and we have a set of values; $I_i$. The approach taken in
\ARTS\ is based on the observation that the practical calculation of
Eq.~\ref{eq:sensorint} can be written as
\begin{equation}
  \label{eq:sensorsum}
  \sum_i h_i I_i.
\end{equation}
This expression assumes that $r$ is independent of $I$, which is generally
valid. The summation weights $h_i$ are pre-calculated and stored in a matrix
\SnsMtr. The \SnsMtr\ of each sensor component can be calculated separately:
\begin{equation}
  \label{eq:Hs}
  \SnsMtr = \aSnsMtr{n}\dots\aSnsMtr{2}\aSnsMtr{1},
\end{equation}
where $n$ is the number of sensor components considered. The inclusion of sensor
characteristics is then simply made as
\begin{equation}
  \label{eq:Hi}
   \MsrVct = \SnsMtr\MpiVct,
\end{equation}
where \MpiVct\ is a vector, where the Stokes vectors from each monochromatic
radiance calculation are appended, and \MsrVct\ is the final ``measurement
vector''. This approach was introduced by \cite{eriksson:studi:02} and
elaborated further in \cite{eriksson:effic:06}.

%Methods exist for including the following sensor components and observation
%strategies: antenna, mixer, spectrometer, beam switching and frequency
%switching. Instruments having mixers at different frequencies (such as AMSU)
%require special attention and have a dedicated method. The description of the
%responses can be detailed. For example, the response of the antenna can be
%described as a function of polarisation, frequency, zenith angle and azimuth
%angle. For flexibility, for dimensions where the corresponding grid has length
%1, the convention is that the response is taken as constant or is ignored, and
%the antenna response data can be as simple as a vector describing the zenith
%angle variation. 

The method presented in \cite{buehler:effic:10} to efficiently handle broadband
infrared channels is also implemented in \ARTS. The approach can be seen as an
extension of Eq.~\ref{eq:Hi}, where the frequencies of \MpiVct\ and the
``weights'' in \SnsMtr\ are selected in parallel. This in order to approximate
Eq.~\ref{eq:sensorint} over a large range of atmospheric conditions with the
lowest possible number of monochromatic frequencies (length of \MpiVct).


\subsection{Transmission and batch calculations}
%
The standard \ARTS\ case is measurements of direct or scattered emission, but
also pure transmission calculations can be treated. For example,
%The modular approach of \ARTS\ made this addition
%simple. The main step was to make a workspace method for the clear-sky
%radiative transfer agenda that applies the Beer-Lambert law for all Stokes
%components, accompanied with the corresponding changes for the Jacobian
%expressions. This method does only need to consider radiative transfer along a
%single pencil beam, all other operations are handled by the existing methods.
%This was complemented with a method for the cloudbox, including also particle
%extinction.
%These two methods make 
it possible to simulate solar occultation and satellite-to-satellite
transmissions. This includes particle effects, as long as the (re-)scattering
into the line-of-sight can be neglected. 
%The radiation at the top of the
%atmosphere is then set to solar radiation or a unit Stokes vector (e.g.\
%[1,0,0,$\pm$1]$^T$ if circular polarisation is transmitted), instead of cosmic
%background radiation as for the emission case. The agenda for surface
%properties should also be configured to not add any emission.

\ARTS\ includes now a very general mechanism for batch calculations. This is
handled by an agenda (Sec.~\ref{sec:agendas}) that contains the methods that
should be executed for each batch case. Batch calculations are particularly
efficient with absorption lookup tables (Sec.~\ref{sec:abstable}), since the
table has to be calculated (or read from file) only once, and can then be used
for all cases. A typical application of this is to simulate satellite
measurements for a large number of atmospheric scenarios.




\subsection{Radiance units}
\label{sec:units}
%
The flexibility of \ARTS\ has the consequence that there is no fixed unit for
the measurement vector $\MsrVct$. The unit depends primarily on the method
selected to set the emission source term, but the sensor response
matrix (\SnsMtr) can also include operations that change the unit.

The standard definition inside \ARTS\ of the Planck function is
\begin{equation}
  \label{eq:planck}
  B(T) = \frac{2h\Frq^3}{c^2(exp(h\Frq/k_BT)-1)},
\end{equation}
where $h$ is the Planck constant, $c$ the speed of light and $k_B$ the
Boltzmann constant. This expression gives the (total) power per unit frequency
per unit area per solid angle and the resulting unit is W/(m$^2$\,Hz\,sr). As
long as Eq.~\ref{eq:planck} is followed, \ARTS\ supports
conversion to the following units: \vspace{1ex}\\
\noindent
{\bf W/(m$^2$\,m\,sr)}, power per unit wavelength per unit area per solid angle\\
\noindent
{\bf W/(m$^2$\,m$^{-1}$\,sr)}, as above but per unit wavenumber\\
\noindent
{\bf RJBT}, brightness temperature (\BT) following the Rayleigh-Jeans
approximation of the Planck function [K]\\
\noindent {\bf PlanckBT}, brightness temperature following the Planck function
[K]
\vspace{1ex}\\
The two first conversions correspond to linear mappings, and a common rescaling
factor can be applied for all Stokes elements, polarisation components and the
Jacobian. The conversion to brightness temperatures is more complex. In the
text below, all primed quantities ($I'$, $Q'$, $\Iv'$, ...) refer to brightness
temperatures (RJ or Planck), whereas all unprimed quantities ($I$, $Q$, $\Iv$,
...) refer to radiances.


\subsubsection{Stokes element $I$}
%
The first Stokes element is converted to PlanckBT by inverting
Eq.~\ref{eq:planck},
\begin{equation}
  \label{eq:invplanck}
  I' = B^{-1}(I)=\frac{h\Frq}{k_B\ln((p_nh\Frq^3/c^2I)+1)},
\end{equation}
while the conversion to RJBT uses the standard approximative
expression
\begin{equation}
  \label{eq:rjtb}
  I' = \frac{c^2}{p_n\Frq^2k_B} I.
\end{equation}
The factor $p_n$, representing the number of polarisation modes
\cite{matzler-2006-radia}, is introduced for reasons of generality (see below).
For $I$, $p_n=2$ in both equations above (to match Eq.~\ref{eq:planck}).

The conversion of the Jacobian to PlanckBT requires further considerations.
The derivative of a radiance in PlanckBT, with respect to a variable $x$, can
be formulated as
\begin{equation}
  \label{eq:dipdx}
  \frac{\PartD I'}{\PartD x} = \frac{\PartD B^{-1}(I)}{\PartD x} =   
   \frac{\PartD B^{-1}(I)}{\PartD I}\frac{\PartD I}{\PartD x}.
\end{equation}
The term $\PartD I/\PartD x$ is the weighting function for the original
unit, that shall be converted to PlanckBT. The conversion term can be derived
to be
\begin{equation}
  \frac{\PartD B^{-1}}{\PartD I} = \frac{k_B\left[B^{-1}(I)\right]^2}
                     {h\Frq I(1+(c^2I/2h\Frq^3))}.
\end{equation}



\subsubsection{Stokes elements $Q$, $U$ and $V$}
%
The conversion of $Q$, $U$ and $V$ to RJBT is made exactly as for $I$. That is,
Eq.~\ref{eq:rjtb} is applied with $p_n=2$. This deviates from e.g.\
\cite{matzler-2006-radia} (setting $p_n=1$ for these Stokes elements), but is
preferred for reasons of generality. A practical consideration for the Stokes
vector is that the ratio between the elements must be the same independent of
the selected unit. Otherwise it would be needed to adapt optical properties,
e.g.\ \ExtMat\ (Eq.~\ref{eq:VRTE}), to the selected unit. Another way to
express this is that, in the Rayleigh-Jeans limit, the same result shall be
obtained if Eq.~\ref{eq:planck} is used and radiances are converted to RJBT, as
if the emission source term ($B$) is replaced by the physical temperature
(\Tmp). \ARTS\ allows the latter, see \cite{liu-2008-conve} for a discussion of
this choice. (It should be noted that these two options do not generally give
the same \BT).

As Eq.~\ref{eq:invplanck} is a non-linear mapping, it can not be applied
directly on $Q$, $U$ and $V$. To maintain the basic properties of the Stokes
vector, $Q$ is converted to PlanckBT as (cf.\ Eqs.~\ref{eq:Q} and \ref{eq:Iv})
\begin{equation}
  \label{eq:invplanck2}
  Q' = B^{-1}\left(\left[I+Q\right]/2\right) - 
       B^{-1}\left(\left[I-Q\right]/2\right).
\end{equation}
The conversion of weighting functions must be done in a similar manner
\begin{equation}
  \frac{\PartD Q'}{\PartD x} = \left[
      \left.\frac{\PartD B^{-1}}{\PartD I}\right|_{(I+Q)/2}  +
      \left.\frac{\PartD B^{-1}}{\PartD I}\right|_{(I-Q)/2}
     \right] \frac{\PartD Q}{\PartD x}.
\end{equation}
The elements $U$ and $V$ are treated likewise.



\subsubsection{Individual polarisation components}
%
The measurement vector \MsrVct\ can contain either Stokes elements ($I$, $Q$,
...) or individual polarisation components (\Iv , \Ih , ..., see
Sec.~\ref{sec:rtvar}). In the later case this is taken as a calibrated
observation and, as the data correspond to a single polarisation mode, the
conversion to \BT\ must be adapted. The reference for the conversion is then
the blackbody radiation for a single polarisation mode, that is a factor 2
smaller than Eq.~\ref{eq:planck}. The conversion from radiance to \BT\ is thus
made through Eqs.~\ref{eq:invplanck} and \ref{eq:rjtb} with $p_n=1$.

If individual polarisation components are extracted outside \ARTS, it is
important to note that the definitions above have the consequence that
Eq.~\ref{eq:Iv} can not be applied if the data have been converted to \BT. As
example, the brightness temperature for the vertical linear component is
obtained as
\begin{equation}
  \Iv' = I'+Q',
\end{equation}
which differs from Eq.~\ref{eq:Iv} with a factor of two.


\section{Conclusions}
%-----------------------------------------------------------------------------
%
The first version of \ARTS\ (\V1) was a traditional microwave to infrared
clear-sky forward model; it was 1D and had no treatment of
scattering. The main novelty of \V1\ was the introduction of the
workspace. However, the ambition of easily extendable software was not
fully met by \V1, and the concept was for this version extended by an
agenda mechanism. Our experience so far is that the desired degree of
modularity has been reached.

The new \ARTS\ version (2.0) is a state-of-the-art radiative transfer model for
the thermal spectral region, as it combines the following features:
\begin{itemize}
\item The model atmosphere can be 1D, 2D or 3D. Tomographic limb sounding
  retrievals require 2D or 3D, and rigorous cloud scattering simulations are
  only possible in 3D.
\item Spherical geoid and surface are throughout default. For 2D and 3D more
  complex topography are also possible. A `flat Earth' is not a viable
  option for limb sounding.
\item Radiative transfer can be made for 1 to 4 Stokes elements. Polarisation
  effects can thus be fully described.
\item Basically no restriction in complexity of surface reflection (but is currently handled only in a simplistic manner).
\item For particle single scattering properties, not only the standard
  assumption of spherical or completely randomly oriented particles, but also
  the case of horizontally aligned particles is handled. 
\item Two modules for solving radiative transfer with particle scattering: MC
  and DOIT. Both modules lack intrinsic
  approximations, and have been verified by practical retrievals
  \cite{eriksson:first:07,rydberg:nonga:09}.
\item Sensor responses can be incorporated in an efficient manner
  \cite{eriksson:effic:06,buehler:effic:10}.
\end{itemize}
Another way to judge the scientific merits of \ARTS-2.0 is the fact that it has
already been used for a number of scientific publications. Direct usage of 
\ARTS-2.0 includes 
\cite{emde:theef:04,davis05:_cirrus_induc_polar_ghz_aura,hoepfneremde05:_ir_comparison,teichmann06:_understanding_jqsrt,battaglia07:_microw_jqsrt,buehler:aconc:07,davis07:_3d_pol_sim_acp,ekstrom:first:07,eriksson:first:07,pietranera07:_cmb_polarization_mnras,rydberg:predi:07,adams2008,sreerekha2008,Haefele2009,harlow-2009-milli,milz09:_comparison_grl,rydberg:nonga:09,holl-2010-collo}, 
and indirect usage is found in yet more journal articles. 

The main limitations of \ARTS-2.0 are:
\begin{itemize}
\item Physical mechanisms not yet implemented include non-local thermodynamic
  equilibrium and polarised gas absorption.
\item Particle single scattering properties must be calculated externally.
\item Extremely fast calculations are not within the present
  scope of \ARTS. The same applies to calculation of radiative fluxes
  and cooling rates.
\item Weighting functions can be obtained, but so far only for a limited number
  of variables under non-scattering conditions.
\end{itemize}
The web address for \ARTS\ is \url{www.sat.ltu.se/arts}, where the software can
be downloaded freely and additional documentation is found. Please, note the
``code of conduct'' found on the web site, asking users to cite this and the
relevant module specific articles
\cite[at the time of writing:][]{emde:apoli:04,davis05:_monte_carlo_i3etgrs,eriksson:effic:06,buehler:effic:10,buehler-absor}.



\section*{Acknowledgements}

First of all, we thank all others that have contributed to the ARTS
development over the years, particularly 
%
Hermann Berg, 
Mattias Ekstr\"om, 
Axel von Engeln, 
Wolfram Haas,  
Carlos Jim\'enez, 
Viju John,  
Nikolay Koulev,  
Thomas Kuhn,  
Mashrab Kuvatov, 
Christian Melsheimer,  
Mathias Milz, 
Randall Skelton, 
Claas Teichmann,  
and Sreerekha T.R.
%
The second acknowledgement is due to all those in the ARTS user community that
have contributed interesting use cases, applications, and discussions on the
ARTS mailing list (\url{arts-users@www.sat.ltu.se}). Seeing the model being
used is the biggest reward for us who have developed it!

Last but not least, we acknowledge our funding agencies. The ESA study
`Development of a Radiative Transfer Model for Frequencies between 200 and 1000
GHz' (contract 17632/03/NL/FF) was especially valuable. Further, direct or
indirect support has been provided by the Swedish National Space Board (e.g.\
59/07 and 92/08) and by the Swedish Research Council (2007-3720 and 2007-5370),
through different contracts at both the Swedish groups. Development of the
ARTS-MC algorithm was largely funded by the UK Natural Environment Research
Council.



%--- References --------------------------------------------------------------
\bibliography{j_abbr,references}




\end{document}

