The libRadtran software package for radiative transfer calculations ( Version 2 . 0 . 1 )

libRadtran is a widely used software package for radiative transfer calculations. It allows to compute (polarized) radiances, irradiances, and actinic fluxes in the solar and thermal spectral regions. libRadtran has been used for various applications, including remote sensing of clouds, 5 aerosols and trace gases in the Earth’s atmosphere, climate studies, e.g., for the calculation of radiative forcing due to different atmospheric components, for UV-forcasting, the calculation of photolysis frequencies, and for remote sensing of other planets in our solar system. The package has 10 been described in Mayer and Kylling (2005). Since then several new features have been included, for example polarization, Raman scattering, a new molecular gas absorption parameterization, and several new cloud and aerosol scattering parameterizations. Furthermore a graphical user in15 terface is now available which greatly simplifies the usage of the model, especially for new users. This paper gives an overview of libRadtran version 2.0.1 with focus on new features. Applications including these new features are provided as examples of use. A complete description of libRadtran and 20 all its input options is given in the user manual included in the libRadtran software package, which is freely available at http://www.libradtran.org.


Introduction
Radiative transfer modelling is essential not only for remote sensing of planetary atmospheres, but also for many other fields in atmospheric physics: atmospheric chemistry, which is largely influenced by photochemical reactions, calculation of radiative forcing in climate models, and radiatively driven dynamics in numerical weather prediction models.
The libRadtran software package is a versatile toolbox, which has been used for various applications related to atmospheric radiation, a list of publications that have used the package can be found on the website http://www.libradtran.org; currently it includes more than 400 entries.Applications include the following topics (the given references are taken as examples out of the list of publications): -analysis of UV-radiation measurements, from which parameters, e.g.ozone concentrations, aerosol optical thickness, UV-index, are derived.Since the libRadtran package originally was a radiative transfer code for the UV spectral range (the main executable is still called uvspec), the model is well established in this research area and frequently used (e.g.Seckmeyer et al., 2008;Kreuter et al., 2014); Published by Copernicus Publications on behalf of the European Geosciences Union.

C. Emde et al.:
The libRadtran software package -cloud and aerosol remote sensing using measurements in solar and thermal spectral regions.The developed retrieval methods are for ground-based, satellite and airborne instruments that measure (polarized) radiances (e.g.Painemal and Zuidema, 2011;Bugliaro et al., 2011;Zinner et al., 2010;Alexandrov et al., 2012); -volcanic ash studies including remote sensing of ash mass concentrations (e.g.Gasteiger et al., 2011;Kylling et al., 2015) and visibility of ash particles from the pilot's perspective (e.g.Weinzierl et al., 2012); -remote sensing of surface properties: a model like li-bRadtran is particularly important to develop atmospheric correction methods (e.g.Drusch et al., 2012;Schulmann et al., 2015); -trace gas remote sensing: libRadtran used as a forward model for retrievals of O 3 , NO 2 , and BrO from DOAS (Differential Optical Absorption Spectroscopy) measurements (e.g.Theys et al., 2007;Emde et al., 2011); -calculation of actinic fluxes in order to quantify photolysis rates for atmospheric chemistry (e.g.Suminska-Ebersoldt et al., 2012); -determination of solar direct irradiance and global irradiance distributions in order to optimize locations of solar energy platforms (e.g.Lohmann et al., 2006) and calculation of circumsolar irradiance (Reinhardt et al., 2014); -simulation of satellite radiances to be used for data assimilation in numerical weather prediction models (Kostka et al., 2014); -validation of radiation schemes included in climate models (Forster et al., 2011), calculation of radiative forcing of clouds and contrail cirrus (Forster et al., 2012), impacts of aviation on climate (e.g. Lee et al., 2010); -simulation of heating rates in three-dimensional (3-D) atmospheres to develop fast radiation parameterizations for large eddy simulation (LES) models (Klinger and Mayer, 2014); -simulation of solar radiation during a total eclipse (Emde and Mayer, 2007); -rotational Raman scattering explaining the filling-in of Fraunhofer lines in the solar spectrum (Kylling et al., 2011); -Estimation of background radiation affecting lidar measurements (e.g.Ehret et al., 2008); -Remote sensing of planetary atmospheres (e.g.Rannou et al., 2010).
Since the publication of the first libRadtran reference paper (Mayer and Kylling, 2005), the model has been further developed.It includes numerous new features that will be the focus of this paper.
One of the major extensions is the implementation of polarization in the radiative transfer solver MYSTIC (Monte Carlo code for the phYSically correct Tracing of photons In Cloudy atmospheres) (Emde et al., 2010), which is important because an increasing number of polarimetric observations have been performed during the last years and are planned for the future, from ground, satellite, and aircraft.These observations include more information about optical and microphysical properties of atmospheric particles than total radiances alone (Kokhanovsky et al., 2010b;Mishchenko et al., 2007).Another important reason for considering polarization is that in the shortwave spectral region (below about 500 nm), the neglect of polarization can lead to large errors: more than 10 % for a molecular atmosphere and up to 5 % for an atmosphere with aerosol (Mishchenko et al., 1994;Kotchenova et al., 2006).
Moreover libRadtran now includes a solver to calculate rotational Raman scattering (Kylling et al., 2011), which improves the accuracy of trace gas retrievals.Further the Raman-scattering signal can be used to estimate cloud top pressure from satellite measurements and aerosol properties from surface and satellite observations.
Numerous state-of-the-art parameterizations for aerosol and ice cloud optical properties have been included (see Sects. 5 and 6).These new parameterizations provide more accurate radiance calculations.In particular for polarized radiative transfer, which requires not only the scattering phase function but the full scattering phase matrix, new data on optical properties were required.In order to improve the accuracy for highly peaked phase functions -which are typical for ice clouds -an improved intensity correction method has been developed and included into the DISORT solver (Buras et al., 2011), and new variance reduction methods have been developed for the Monte Carlo solver MYSTIC (Buras and Mayer, 2011).libRadtran has also been rewritten to allow for simulations with an arbitrary number of cloud and aerosol types -which can, e.g., be used to take into account detailed particle size distributions (number densities for discretized size bins) that can be different in each layer.In earlier versions it was only possible to take into account parameterized size distributions such as gamma or log-normal distributions.
A new gas absorption parameterization for the solar and thermal spectral ranges has been developed (Gasteiger et al., 2014).It is available in different spectral resolutions and can be applied for the simulation of radiances and irradiance.It is particularly useful for efficient simulations of radiances measured by satellite instruments (see Sect. 4.1).
The DISORT radiative transfer solver has been translated from FORTRAN77 to the C programming language.All variables were transferred from single to double precision.These changes improved the numerical stability of the code and reduced computational time significantly (for details see Buras et al., 2011).
The paper is organized as follows: Sect. 2 provides an overview of the uvspec radiative transfer model, which is the core of the libRadtran package.Section 3 gives a short description of the radiative transfer solvers included in uvspec.Section 4 provides a summary of how molecules are handled and outlines various ways to include molecular absorption.Moreover, Rayleigh-scattering parameterizations are described.Section 5 summarizes the available parameterizations for aerosol microphysical and optical properties.Section 6 gives an overview of the parameterizations for water and ice clouds and also outlines how these were generated.In Sect.7 available surface properties are described, including Lambertian reflection, bi-directional distribution functions and fluorescent surfaces.In Sect.8 we describe code and implementation improvements relevant for users.Section 9 introduces the graphical user interface for uvspec.Section 10 provides a short summary of additional tools that come with the libRadtran package.Finally, Sect.11 shows a few applications as examples of the usage of libRadtran.

The uvspec radiative transfer model
The main tool of the libRadtran package is the uvspec radiative transfer model, which consists of the following parts: 1.The atmospheric state (e.g.trace gas profiles, cloud liquid water content, cloud droplet size, aerosol concentration profiles) needs to be provided as input to the model.quantities are normalized to the source function, i.e. the solar irradiance in the solar spectral region.In order to get physical quantities with corresponding units the output may be post-processed.The uvspec output then corresponds to calibrated radiances or brightness temperatures for a given instrumental filter function.It is also possible to obtain integrated solar or thermal irradiance.
The overall structure of the uvspec model is shown in Fig. 1.
The model was originally designed to compute UV radiation; therefore, its name is uvspec.As said before it now covers the complete solar and thermal spectral range.
The usage of the model is described in the user guide, which comes along with the package.The user guide includes descriptions of the RTE solvers, examples of use as well as detailed documentation of all options and respective parameters.Below uvspec input options are put in teletypefont, for example rte_solver.
The uvspec model may be run either from the command line using uvspec < input_file > output_file or from the graphical user interface (see Sect. 9).

Radiative transfer equation solvers
The RTE for a macroscopically isotropic medium, i.e. randomly oriented particles and molecules, may be written as   (Chandrasekhar, 1950;Mishchenko et al., 2002) where the source function J is Here I = (I, Q, U, V ) is the Stokes vector at location (x, y, z), β the volume extinction coefficient, ω 0 the singlescattering albedo, P( , ) the scattering phase matrix, and B e (T ) = (B(T ), 0, 0, 0) the emission vector including the Planck function B(T ).For most applications in the Earth's atmosphere, thermal emission can be neglected for wavelengths below about 3 µm.Polarization is also often neglected, in this case the Stokes vector in Eqs.(1) and ( 2) is replaced by the radiance L, the phase matrix becomes the scalar phase function p( , ) and the emission vector is just the Planck function B(T ).
The uvspec model includes various methods to solve Eq. ( 1).The list of solvers, which may be selected using the option rte_solver, is shown in Table 1.

DISORT
The solver disort is used by default in libRadtran.DIS-ORT (Stamnes et al., 2000) is based on discrete ordinates and allows one to compute radiances, irradiance, and actinic fluxes in plane-parallel geometry.The original FORTRAN77 version of the algorithm exhibited several numerical instabilities for certain combinations of geometries and optical properties.The FORTRAN77 code has been translated to C-code and is entirely in double precision (the FORTRAN77 version is mostly in single precision) and includes dynamic memory allocation (not possible in FORTRAN77).As such, the C version is numerically stable and also faster than the original FORTRAN77 version.We thus use the C version of the DIS-ORT algorithm by default.The original FORTRAN77 version may still be invoked by fdisort2.Both the C-code and the FORTRAN77 version include the new intensity correction method for peaked phase functions by Buras et al. (2011), which is used by default.
For calculations with rotational Raman scattering, the C version has been generalized so that arbitrary source functions (not only a solar or thermal source function) can be handled (Kylling and Stamnes, 1992;Kylling et al., 2011).Rotational (inelastic) Raman scattering from other wavelengths into the wavelength, for which the RTE is solved, is included into the source term.

MYSTIC
The most comprehensive solver in libRadtran is the Monte Carlo model MYSTIC (Mayer, 2009), which may be used to calculate (polarized) radiances, irradiance, and actinic fluxes in the solar and thermal spectral regions.Within MYSTIC photons are traced through the atmosphere from the source towards the sensor or backwards, from the sensor to the source, which is much more efficient especially in the thermal wavelength region.One of the main applications of MYSTIC is to calculate radiances in cloudy atmospheres.The sharp forward scattering of clouds and aerosols causes numerical problems in Monte Carlo models.In order to avoid these, sophisticated variance reduction methods have been developed (Buras and Mayer, 2011).These are enabled using mc_vroom on.Solar radiation is initially unpolarized and becomes polarized by molecular, aerosol, or cloud scattering in the atmosphere.With the option mc_polarisation (Emde et al., 2010), the full Stokes vector is calculated.For 1-D atmospheres, MYSTIC may also be operated in spherical geometry using the option mc_spherical (Emde and Mayer, 2007).
The public version of MYSTIC allows for calculations in 1-D (plane-parallel or spherical) geometry.A full 3-D version is also available for joint projects.The non-public version includes several other features: complex 3-D topography (Mayer et al., 2010) and efficient high-spectralresolution calculations using absorption lines importance sampling (Emde et al., 2011).

Two-stream solvers
For the calculation of irradiance, two fast two-stream solvers are available.
The first solver, twostr, is described in detail in Kylling et al. (1995).twostr is optimized for calculating actinic fluxes, and hence heating rates.It can be run in plane parallel as well as in pseudo-spherical geometry.
The second two-stream method available in libRadtran is rodents, which is based on the delta-Eddington two stream described, e.g., in Zdunkowski et al. (2007), Sects.6.1-6.41 .Based on a different two-stream approach than twostr, it naturally yields different results.In contrast to twostr, neither the pseudo-spherical approximation is implemented nor is rodents capable of calculating actinic fluxes.
For actinic fluxes and atmospheric heating rates, twostr is the better choice.However, for calculating solar irradiance, we recommend using rodents: for cases where the resulting irradiance is not negligible (larger than 2 % of the extraterrestrial irradiance), the difference between rodents and exact disort calculations is on average 5 % (7 %) for down(up)-welling irradiance.For twostr the values are 9 % (11 %).Especially in case the atmosphere is only weakly absorbing, the average differences at top-ofatmosphere (TOA) and at the surface are only 2 % (1 %) for rodents, whereas they are 5 % at TOA and even 13 % (18 %) at surface for twostr.
For the thermal irradiance, rodents also gives better results at TOA (1.6 %) and at the surface (1 %) than twostr (3 %).For irradiance within the atmosphere, no real preference can be given.

Lidar and radar simulations
In order to complement the instruments that can be simulated by libRadtran, a lidar simulator called sslidar has been implemented.It only takes into account single scattering and reflection and is based on the lidar equation, which is integrated over each range.Note that in order to obtain a smooth signal, a fine vertical resolution of the model atmosphere is required.The vertical resolution should correspond to the range width of the simulated lidar instrument.For radar simulations a stand-alone tool is available (see Sect. 10.2).

Other solvers
The solver tzs (see Appendix B) is based on the zeroscattering approximation in the thermal spectral range.It may be used for clear-sky calculations of radiances at TOA.It also calculates "black cloud" radiances for the application of the CO 2 -slicing algorithm (Smith et al., 1970;Chahine, 1974;Smith and Platt, 1978;Menzel et al., 1983;Eyre and Menzel, 1989), which may be used for the determination of cloud top temperatures from passive remote sensing measurements in the thermal spectral range.
For the solar region a fast single-scattering solver sss is available.These solvers may be used for fast but approximate simulations of satellite measurements.
Several other RTE solvers are included in uvspec for compatibility with earlier releases of the package.These include sdisort (pseudospherical disort), spsdisort (single precision, pseudospherical disort), fdisort1 (version 1 of DISORT), and polradtran (Evans and Stephens, 1991).While they may still be used, we do not recommend their use as the other solvers listed in Table 1 perform better.

Accuracy of solvers
The MYSTIC model has been validated in many international model intercomparison studies, for radiance calculations with highly peaked phase functions (Kokhanovsky et al., 2010a), for polarized radiance calculations (Emde et al., 2015), and for radiances and irradiance in 3-D model domains (Cahalan et al., 2005).In all studies MYSTIC belongs to the core of models that produce equal results within their uncertainty range.MYSTIC agrees perfectly with DIS-ORT for radiances and irradiance with only a few exceptions, e.g. for circum-solar radiation, where the second-order intensity correction included in DISORT is not accurate enough for highly peaked scattering phase functions (Buras et al., 2011).In Emde et al. (2011), a comparison between DIS-ORT and MYSTIC for a radiance spectrum in the O 2 -A band is shown.The relative difference between the solvers is less than 0.05 % here.All other solvers are approximations and hence less accurate: as mentioned before the two-stream solvers are only appropriate for irradiance and the tzs solver only provides radiances in thermal atmospheres and neglects scattering completely.
The accuracy of MYSTIC depends only on the number of traced photons.The standard deviation of MYSTIC is calculated when the option mc_std is enabled.The user may run MYSTIC with many photons as reference for some cases in order to check the accuracy of other solvers for specific applications.for the midlatitude-summer atmosphere of Anderson et al. (1986).All calculations were performed with the MYSTIC solver using the "absorption lines importance sampling" method (Emde et al., 2011).(Top) high spectral resolution calculation, based on line-by-line absorption cross sections calculated using ARTS (Eriksson et al., 2011); (bottom) pseudo-spectral calculations using the representative wavelengths band parameterizations (reptran) with different resolutions and lowtran.For comparison see also Fig. 3 in Mayer and Kylling (2005), which shows transmittances for genln2 line-by-line calculations and lowtran for the same spectral regions.

Molecular absorption parameterizations
Spectral ranges affected by molecular absorption comprising a complex line structure require parameterizations to reduce the computational cost.Molecular absorption parameterizations included in libRadtran are listed in Table 2.By default the reptran parameterization is applied.Using the option mol_abs_param, the user may select the most appropriate parameterization for the specific application.As an example Fig. 2 shows radiance calculations for nadir viewing direction at the top of the atmosphere using the parameterizations reptran and lowtran and line-by-line calculations.
The reptran parameterization (Gasteiger et al., 2014) has recently been included in libRadtran.In reptran integrals over spectral intervals, e.g.integrated over a narrow spectral band or an instrument channel response function, are parameterized as weighted means over representative wavelengths similar to the method described by Buehler et al. (2010).The selection of an optimum set of representative wavelengths is based on accurate line-by-line simulations for top-of-atmosphere radiances of a highly variable set of atmospheric states.The ARTS (Atmospheric Radiative Transfer Simulator) model (Eriksson et al., 2011) including state-of-the-art continuum models and spectroscopic data from HITRAN 2004 (Rothman et al., 2005) were used to calculate the gas absorption properties.For wavelengths below 1130 nm measured absorption cross sections of O 3 (Molina and Molina, 1986), O 4 (Greenblatt et al., 1990), and NO 2 (Burrows et al., 1998) are included, as they are not covered by HITRAN or the continua (see also Sect.4.2).Three-band resolutions (fine: 1 cm −1 ; medium: 5 cm −1 ; and coarse: 15 cm −1 ) are available in the solar and thermal spectral range, as well as a number of instruments on the following satellites: ADEOS (Advanced Earth Observing Satellite), ALOS (Advanced Land Observing Satellite), Earth-CARE (Earth Clouds, Aerosols and Radiation Explorer), Envisat (Environmental Satellite), ERS (European Remote-Sensing Satellite), Landsat, MSG (Meteosat Second Generation), PARASOL (Polarization and Anisotropy of Reflectances for Atmospheric Sciences coupled with Observations from Lidar), Proba (Project for On-Board Autonomy), Sentinel, Seosat (Satélite Espanol de Observación de la Tierra), and SPOT (Satellite Pour l'Observation de la Terre).The parameterization has been validated by comparison to high spectral resolution calculations.For solar and thermal radiation at the top of atmosphere, as well as for solar radiation on the ground, the mean parameterization error is in the range of 1 %.The mean error is slightly larger than 1 % for thermal radiation at the surface.
For the simulation of radiances and irradiance, we recommend to use reptran because it is faster and more accurate than lowtran.
Several correlated-k parameterizations with different numbers of bands, i.e. different accuracy, are included in libRadtran.For the calculation of integrated solar and thermal irradiance and heating rates, the correlated-k parameterizations by Kato et al. (1999) and Fu andLiou (1992, 1993) are recommended.Also for the calculation of heating/cooling rates in the higher atmosphere (above ∼ 20 km), we recommend these parameterizations because reptran and lowtran are affected by large errors.

Molecular absorption cross sections
For the spectral region from 160 to 850 nm, libRadtran includes measured absorption cross sections of various molecules in the atmosphere (see Table 3).Using the option mol_abs_param crs, these cross sections are used instead of the default reptran parameterization.For wavelengths below 500 nm, reptran yields approximately the same results as mol_abs_param crs because the cross sections from HITRAN and the continua are very small at these wavelengths and the same measured cross sections are relevant in both cases.
For O 2 for instance the cross section data include the Schumann-Runge bands between 176 and 192.6 nm and the Herzberg continuum between 205 and 240 nm.Ozone ab-sorption bands are for example the Huggins bands between 320 and 360 nm and the Chappuis bands between 375 and 650 nm.Using the option crs_model the user may specify which cross section data should be used in the simulations.Alternatively with crs_file, the users may specify their own absorption cross section data.

Line-by-line calculations
In the shortwave infrared, thermal infrared, and microwave region, we find a huge number of absorption lines that are due to vibrational or rotational transitions in molecules.A line-by-line model is required in order to calculate spectrally resolved radiances.Line-by-line models take the absorption line positions as well as line strength parameters from spectral databases like HITRAN, calculate line broadening, which depends on pressure and temperature in the atmosphere, and finally obtain absorption optical thickness profiles.libRadtran does not include a line-by-line model but it allows one to specify absorption optical thickness profiles using the option mol_tau_file abs.It is convenient to use the ARTS model (Eriksson et al., 2011) to generate spectrally resolved molecular absorption data because it outputs the format required by libRadtran.ARTS includes a comprehensive line-by-line module, it allows one to use different spectroscopic databases like HITRAN as input and it also includes various state-of-the-art absorption continuum models.The toolbox Py4CATS (Schreier and Böttger, 2003;Schreier, 2006;Schreier and Kohlert, 2008) which can be downloaded from www.libradtran.org,also includes conve- nient command line programs to generate spectrally resolved absorption data.The Py4CATS tools, however, do not include continuum models; hence, it should only be used for simulations where the continua are not relevant.

Rayleigh-scattering cross sections
The Rayleigh-scattering cross sections are by default calculated using Eqs.( 22)-( 23) of Bodhaine et al. (1999).Using the option crs_model rayleigh, the user may select Eq. ( 29) of Bodhaine et al. (1999) or the formulas proposed by Nicolet (1984) and Penndorf (1957), respectively.The analytical Rayleigh-scattering phase matrix P R (Hansen and Travis, 1974) is and δ is the depolarization factor that accounts for the anisotropy of the molecules; δ is also calculated according to Bodhaine et al. (1999).The Rayleigh phase matrix for δ = 0 is shown in Fig. 3.For calculations neglecting polarization only the (1, 1) element of the phase matrix, which corresponds to the scattering phase function, is required.

Aerosols
Besides the models by Shettle (1989), which are described in Mayer and Kylling (2005), libRadtran now includes additional aerosol properties based on the OPAC (Optical Properties of Aerosols and Clouds) database (Hess et al., 1998).
For the soluble aerosols the parameters depend on humidity because the aerosol particles swell in humid air.Relative humidities of 0, 50, 70, 80, 90, 95, 98, and 99 % are included in OPAC.The option aerosol_species_file allows one to define arbitrary mixtures of these basic types or to select pre-defined mixtures from OPAC, such as, e.g., continental_average, for which uvspec automatically uses the optical properties closest to the background humidity profile.
Optical properties of all basic aerosol types were calculated using libRadtran's Mie tool (see Sect. 10.1).For mineral aerosols, which are highly aspherical, we additionally provide optical properties calculated with the T-matrix method (Mishchenko and Travis, 1998), assuming an aspect ratio distribution of prolate spheroids as described by Koepke et al. (2015).
As an example Fig. 3 shows the phase matrix elements of the basic OPAC aerosol types, of liquid cloud droplets with an effective radius of 10 µm and the Rayleigh-scattering phase matrix.Note that for spherical particles only four elements of the 4 × 4 scattering phase matrix are independent whereas for aspherical particles six elements are required (see, e.g., Hansen and Travis, 1974).Figure 4 shows the absorption and the scattering optical thicknesses (integrated from the surface to the top of the atmosphere) for the standard aerosol mixtures in the spectral region from 300 to 800 nm.As expected, the optical thickness of the urban aerosol is the largest and that of the antarctic aerosol the smallest.In general the continental aerosol mixtures show a stronger wavelength dependency than the maritime mixtures.
The users may also provide their own optical properties data, which may be generated using libRadtran's Mie tool or  Phase matrix elements for the basic OPAC aerosol types "water soluble" (waso), "sea salt accumulated mode" (ssam), and soot, for a water cloud with a droplet effective radius of 10 µm, and for Rayleigh scattering (with δ = 0) at a wavelength of 350 nm.θ is the scattering angle, i.e. the angle between incoming and scattered directions.
other external programs; more detailed instructions are provided in the libRadtran user guide.

Water clouds
Table 4 summarizes the parameterizations of water cloud optical properties, which may be selected in libRadtran using the option wc_properties.
For the simulation of irradiance and heating rates, it is normally sufficient to use a simple parameterization to convert from cloud liquid-water content and droplet effective radius to the respective optical properties: extinction coefficient, single-scattering albedo, and asymmetry parameter.For this purpose libRadtran includes the parameterization generated by Hu and Stamnes (1993).
For the simulation of radiances more accurate optical properties are needed and the phase function should not be approximated by a Henyey-Greenstein function as it is done in Hu and Stamnes (1993).Therefore, we have pre-calculated cloud optical properties using libRadtran's Mie tool, assuming that the cloud droplets are gamma distributed: Calculations have been performed for effective radii r eff from 1 to 25 µm with a step width of 1 µm.The effective variance was set to a value of v eff = 0.1 and the constant N was determined by normalization.The size distributions were cut off at a minimum radius of 0.02 • r eff and a maximum radius of 8 • r eff .The size distribution bins are sampled on a size parameter ( 2πr λ ) grid with a resolution of 0.003.This fine resolution is necessary to obtain smooth phase matrices.The pre-calculated data include the wavelength ranges from 250 to 2200 nm (solar) with a resolution of 10 nm and the range from 2.2 to 100 µm (thermal) in 100 steps of equal wavenumbers.The refractive index of water has been taken from Warren (1984).In the solar (thermal) region, the phase matrices are computed from 5000 (500) Legendre polynomials.In the optical properties files, 129 of the Legendre polynomials are stored, as well as the phase matrix elements, which are stored on scattering angle grids θ optimized such that the error of the phase matrix -when interpolated linearly in cos θ between the grid points -is smaller than 1 %.As an example Fig. 3 shows the four phase matrix elements of a cloud droplet distribution with r eff = 10 µm at 350 nm.
Here the cloudbow at θ ≈ 140 • is clearly visible in the P 11 and P 12 /P 11 elements of the phase matrix.P 12 /P 11 corresponds to the degree of polarization in the principal plane after single scattering; it can be seen that the maximum in the cloudbow region is about 80 %.The mystic solver uses the phase matrix stored on the θ-grid, whereas all other solvers use the Legendre polynomials, except for the intensity correction in disort, which uses the phase function (see also Buras et al., 2011).For specific applications, e.g.different size distributions, the user can easily generate optical properties using libRadtran's Mie tool.

Ice clouds
For ice clouds libRadtran includes a variety of parameterizations (see Table 5) from which the user may select the most appropriate one for a specific application by specifying the option ic_properties.Ice clouds are more complex than water clouds because they consist of ice crystals of different shapes.Some of the ice cloud parameterizations allow the crystal habit (ic_habit) to be specified.
As described in the previous section the exact phase matrix is not needed when irradiance are calculated.For this purpose the parameterizations by Fu (1996), Fu et al. (1998), andKey et al. (2002) are included in libRadtran.Fu (1996) and Fu et al. (1998) approximate the phase function by a Henyey-Greenstein function.Key et al. (2002) is slightly more accurate because it uses a double-Henyey-Greenstein function, which represents the backscattering of ice crystals much better.The parameterization is based on single-scattering calculations for various ice crystal habits and on measured size distributions.It is available in the wavelength range from 0.2 to 5 µm.Based on single-scattering data provided by P. Yang and on the size distributions from J. R. Key we have extended the original parameterization by Key et al. (2002) to the thermal wavelength region up to 100 µm.
For accurate radiance calculations the parameterizations by Baum et al. (2005a, b) (baum) and the newer one by Heymsfield et al. (2013), Yang et al. (2013), andBaum et al. (2014) (baum_v36) are available: baum includes full phase functions for a mixture of particle shapes, the parameterization is based on single-scattering properties of smooth ice crystals and on a large number of measured size distributions.baum_v36 includes full phase matrices and three different habit models: a general habit mixture similar to baum but for rough ice crystals, and the single habits solid column and aggregate, both of them severely roughened.
We have generated two further parameterizations (hey and yang2013) for individual habits, which also include the full phase matrices (see Appendix A): hey is available for the wavelength region from 0.2 to 5 µm for smooth particles in the effective radius range from 5 to 90 µm.The full wavelength region from 200 nm to 99 µm is available for yang2013, effective radii may be in the range from 5 to 90 µm and a roughness parameter may also be specified, ranging from smooth to severely rough.For the yang2013 parameterization, the single-scattering properties of nine individual ice crystal habits, which are commonly observed in ice clouds, have been taken from the database by Yang et al. (2013).The hey parameterization was generated before this database existed and it is based on single-scattering data provided by Hong Gang, who used the improved geometrical optics method (IGOM), the same method as used by Yang et al. (2013).
Please refer to the libRadtran user guide for a list of available habits for each parameterization.
Figure 5 shows the phase matrix elements of ice crystal distributions with an effective radius of 40 µm at 550 nm wavelength.The red lines correspond to smooth crystals and the blue lines to severely rough crystals.The individual habits are for the yang2013 parameterization.General habit mixtures, which are available for the hey parameterization based on smooth crystals and for the baum_v36 parameterization based on severely rough crystals, are also shown.For most smooth crystals and also for the general habit mixture ghm of the hey parameterization-scattering features of hexagonal ice crystals, the most prominent being the halo at 22 • scattering angle, are visible in all phase matrix elements.The phase matrices for severely rough crystals do not show halo features and they are relatively similar for all habits.In reality ice clouds are highly variable: There are situations when the halo is visible, in this case obviously there must be regular smooth ice crystals in the cirrus clouds.When no halo is visible, the assumption of severely roughened crystals might be more realistic.

Bi-directional reflectance distribution functions
All solvers included in libRadtran may include Lambertian surfaces, while DISORT and MYSTIC can also handle bidirectional reflectance distribution functions.libRadtran provides a variety of BRDFs (bi-directional reflection distribution function), which are listed in Table 6.
Two parameterizations for land surfaces are available.The first is the "RPV (Rahman, Pinty, and Verstraete)" parameter- (polarized) radiances Single-scattering properties generated by Hong Gang using the code by Yang et al. (2013), Appendix A yang2013 bulk optical properties including phase matrices for 9 habits and 3 degrees of roughness, covers wavelength range from 0.2 to 99 µm.
(polarized) radiances Yang et al. (2013), Appendix A ization by Rahman et al. (1993) with the extension by Degünther and Meerkötter (2000) for modelling snow-covered surfaces.The second is the "RossLi" BRDF first presented by Roujean et al. (1992).The original RossLi BRDF is used in the AMBRALS (the Algorithm for Model Bi-directional Reflectance Anisotropy of the Land Surface) BRDF modelling framework (Wanner et al., 1997), and consists of four different kernel combinations, of which the RossThickLiSparse-Reciprocal combination was identified in several studies to be the model best suited for the operational MODIS BRDF/Albedo algorithm (see Schaaf et al., 2002).An additional factor for simulating the hot spot in vegetation canopies was added by Maignan et al. (2004).The version implemented in libRadtran is the RossThickLiSparse-Reciprocal model as used in MODIS data, as presented in Lucht et al. (2000).The hot spot correction factor can be turned on if required.As already stated in Mayer and Kylling (2005), but repeated here for completeness, a parameterization of the BRDF of water surfaces is also included, which depends mainly on wind speed and to a lesser degree on plankton concentration and salinity.For the MYSTIC solver, also the wind direction can be set.In contrast to vegetation where the typical hot spot occurs in the 180 • backscatter direction, the main feature for water is specular reflection.The parameterization in uvspec was adopted from the 6S code (Vermote et al., 1997) and is based on the measurements of Cox and Munk (1954a, b) and the calculations of Nakajima and Tanaka (1983).A vector version of the ocean parameterization, developed by Tsang et al. (1985) and Mishchenko and Travis (1997), is available for polarization calculations with MYSTIC.The vector version uses only wind speed as a parameter and does not take into account plankton concentration, salinity or wind direction.
Finally, the parameterization of the surfaces of extraterrestrial solid bodies such as the Moon, asteroids, or the inner planets by Hapke (1993) is available.Only the ocean BRDF parameterizations depend directly on the wavelength.For all other BRDF models, the parameterization can either be given as being constant with wavelength (by using, e.g., the option brdf_rpv), or as a file containing the parameters for each wavelength (using, e.g., brdf_rpv_file).

Fluorescence
For vegetation covered surfaces, a weak solar-induced chlorophyll fluorescence signal is emitted in the red and far-red spectral regions.The contribution of fluorescence to the radiance leaving the bottom boundary is where F (λ) is the fluorescence source in the same units as the incoming solar flux at the top of the atmosphere (for example mW (m 2 nm sr) −1 ).The fluorescence source of radiation is included in the disort solver.It may either be constant or vary as a function of wavelength.Additional surface bidirectional reflection of radiation may also be included.The fluorescence source depends on the solar radiation impinging the vegetation and the type of vegetation.Output from vegetation fluorescence canopy models, such as that described by Miller et al. (2005), may readily be used by uvspec.
8 Implementation improvements

Multiple atmospheric constituents
The previous versions of libRadtran were restricted to using at most four types of atmospheric constituents: molecules, aerosols, and water and ice clouds.Any user defined constituent could only be included by replacing, e.g., water clouds with them.Also, it was not possible to use several types of ice cloud habits at the same time.
A recent major internal restructuring of the libRadtran code has now made it possible to use any number of atmospheric constituents for a radiative transfer simulation.The number is only limited by computational memory and time.The new input options needed for loading the additional constituents are profile_file and profile_properties.They work very similar to the cloud input options; merely the name of the constituent needs to be defined.
This option increases the flexibility of libRadtran in many ways; e.g., it can be used to load the optical properties for each size bin of an aerosol or water or ice cloud.This way, the size distribution may differ between the atmospheric layers.An example can be found in Kylling et al. (2013).

Change of nomenclature and backward compatibility
As the number of input options had grown to more than 300 over the years, we decided to restructure the language of the input options.The input options now have a largely consistent naming and their usage follows certain rules, making it more easy to find related input options.
We have included a python script in order to provide backward compatibility for long-established libRadtran users.The script can be found in the directory src_py.By invoking the command python translate.pyinput_file \ > new_input_file input files written in the old nomenclature will be translated to the new nomenclature automatically.Alternatively, the old input file can be sent directly to uvspec with the following command: python translate.pyinput_file | uvspec

Graphical user interface
The large number of input options available in the uvspec model may appear overwhelming.To help the user to cre-ate uvspec input files a graphical user interface (GUI) has been developed.The GUI organizes the input options in logical groups such as "molecular atmosphere", "aerosol", "surface", etc.; see also the grey bar at the top in Fig. 6.Input options that are set by the user and that will be written to the given input files are shown in bold face (for example option rte_solver in Fig. 6).Options that may be set are shown as normal characters, while options that are not compatible with other set options are greyed (for example in Fig. 6 mc_ipa is greyed since it is not possible to combine it with rte_solver set to disort).
Online documentation of the options are available and this is identical to the documentation in the libRadtran user manual.In Fig. 6 the documentation for the option number_of_streams is shown in the lower left corner.The online help is activated by pointing the mouse at the requested input variable.
Input options that refer to input data files, such as wavelength-dependent surface albedo, may be plotted from the GUI.In the example in Fig. 6, the extraterrestrial flux (upper left subplot), the surface fluorescence spectrum (lower left subplot), and surface albedo (lower right subplot) inputs are plotted.Note that the wavelength coverage (x axis) differs reflecting the different wavelength regions included in the input data files.Once all wanted input options are set, they are saved to a user specified file, and uvspec is run from within the GUI.The output from the run may readily be plotted using the GUI.For example, in Fig. 6, the calculated nadir radiance at the top of the atmosphere is shown in the upper right subplot.The GUI includes numerous working examples.Users may add more examples to the GUI specific to their interests.

Other tools
Several additional tools are included in the libRadtran package.An overview is given in Mayer and Kylling (2005, Table 4).New tools are ssradar, a single-scattering radar simulator (see below), and pmom, which calculates Legendre polynomials for a given phase function.

Mie calculations
The tool for Mie calculations (mie) has been extended considerably.The user may select between two Mie codes, MIEV0 by Wiscombe (1980) or bhmie by Bohren and Huffman (1983).The tool allows one to generate input optical properties for uvspec calculations for arbitrary size distributions.It generates full phase matrices, which are stored on optimized angular grids for a user-defined accuracy.The radiative transfer solvers MYSTIC and DISORT with the new intensity correction method (Buras et al., 2011) use the phase functions/matrices rather than Legendre polynomials, which are calculated by the Mie codes.

Single-scattering radar simulator
Single-scattering radar (ssradar) is a stand-alone 1-D pure Rayleigh-scattering cloud radar simulator that handles arbitrary cloud layers and droplet size distributions as well as tilted viewing angles and supercooled water droplets.The radar reflectivity factor is calculated directly from the droplet distribution with Z = i n i D 6 i (Rinehart, 2010) where D is the droplet diameter and n i the distribution number density for the discrete interval D i , D i+1 .Internally available distributions are gamma and log-normal, arbitrary distributions can be entered using input files.

Some applications
The libRadtran package has been used for numerous applications.Many of these are listed under the publications link at http://www.libradtran.org.The examples directory also includes a number input files that may be used especially by new users to create input files.Below some applications of libRadtran are described.

uvspec and ARTS
The high number of absorption lines in the shortwave infrared and the thermal infrared requires a line-by-line approach to resolve the spectral structure.Below it is shown how molecular absorption data from ARTS may be combined with uvspec to perform line-by-line calculations in both the solar and thermal parts of the spectrum.For both examples the spectral resolution, the molecules to be included and the line function properties are specified in the input to ARTS.It is noted that the same ambient atmospheric profile should be used in both, ARTS and uvspec.

Solar source
Solar induced chlorophyll fluorescence is emitted in the 660 to 800 nm spectral region with two broad peaks at about 685 and 740 nm.In this spectral region are the O 2 -A and O 2 -B bands which contain a large number of absorption lines.Although the fluorescence signal is weak, especially the O 2 -B region holds promise for retrieval of vegetation fluorescence from spectrally high-resolution space-borne instruments (Guanter et al., 2010).In this spectral region the surface albedo is typically low while there is a fluorescence peak around 685 nm (see red line right plot Fig. 7).The optical depths from ARTS are input to uvspec, which calcu-lates the top of the atmosphere radiance (blue line, left plot of Fig. 7) including the fluorescence signal (red line, right plot of Fig. 7), surface albedo (green line, right plot of Fig. 7), and molecular scattering.Measurements may be made at a lower spectral resolution.The right plot of Fig. 7 shows radiance spectra convolved with a triangular spectral response function with a full width at half maximum (FWHM) of 0.3 nm using the conv tool of libRadtran.The spectral response function was generated with the make_slitfunction tool.Spectra with (blue line) and without (purple line) fluorescence are presented.It is seen that the fluorescence signal is relatively larger when the surface albedo is low, below about 690 nm, compared to larger wavelengths.

Thermal source
The Infrared Atmospheric Sounding Interferometer (IASI) on board the MetOp satellite measures the radiance from 645 to 2760 cm −1 (15.50-3.6 µm) with a spectral resolution of 0.25 cm −1 .Its main purpose is high-resolution atmospheric sounding of temperature and humidity, and trace gas column retrievals (Clerbaux et al., 2009;Hilton et al., 2011).It may also be used to detect volcanic ash (see Clarisse et al., 2013, and references therein).
The top panel of Fig. 8 shows IASI spectra from a granule covering the ash cloud following the eruption of Mt Kelud, Indonesia, in February 2014.The spectra are classified as cloudless (green), ice cloud (blue), and volcanic ash (red).To investigate the realism of this identification the spectra were simulated with ARTS/uvspec.For all simulated spectra, the surface emissivity was set equal to one which is representative for water.The simulated spectra are shown in the bottom plot of Fig. 8.
The cloudless spectrum has brightness temperatures representative for the ocean at these latitudes.The main molecular absorption features in this part of the spectrum are water vapour lines throughout the spectrum, ozone (broad band feature centred around 1050 cm −1 ), and CO 2 (feature below 800 cm −1 ).The data from ARTS include absorption lines from these molecules.In the cloudless spectrum the ozone band around 1050 cm −1 has a lower brightness temperature than the radiation at lower and higher wavenumber, indicating that the radiation in the ozone band was emitted at a higher altitude with lower temperature than the surface.Overall the ARTS/uvspec cloudless spectrum agrees well with the measured spectrum.
For the simulation with an ice cloud, the ice cloud was located between 12 and 13 km.Ice water content was set to 1 g m −3 .The ice particles were assumed to consist of solid columns with r eff = 40.0µm.The ice cloud parameterization ic_properties yang was selected.The spectrum identified as ice cloud (blue curve in top plot of Fig. 8) appears saturated for nearly all wavenumbers except for the ozone band centred around 1050 cm −1 .The rather low brightness temperature and wavenumber-independent behaviour outside the ozone band, indicates that this is an ice cloud and that it is opaque.The simulation with an ice cloud (blue curve in bottom plot of Fig. 8) agrees well with the measured spectrum.The higher temperatures in the ozone band implies that this radiation was emitted at a higher altitude in the stratosphere where the temperature is higher than at the altitude of the cloud.
The ash simulation included an ash cloud between 17 and 18 km.The ash particles were assumed to be made of andesite, spherical and mono-disperse with a radius of 3 µm.The refractive index of andesite was taken from Pollack et al. (1973) and the optical properties were calculated using the mie tool.The ash density was 1 × 10 −3 g m −3 which corresponds to a mass loading of 1 g m −2 for a 1 km thick cloud.
The red curve in the top plot of Fig. 8 is classified as ash using the difference in brightness temperature method described by Clarisse et al. (2010).This spectrum has a lower brightness temperature than the cloudless spectrum indicating a colder effective emitting temperature overall.The general spectral shape is similar to the cloudless spectrum below 1000 cm −1 .Above about 1200 cm −1 the brightness temperature of the cloudless spectrum generally decreases with increasing wavenumber, while the converse is true for the ash spectrum.The simulated ash cloud spectrum (black curve in bottom plot of Fig. 8) differs from the measured spectrum classified as ash.Both the simulated and measured ash spectra increase in magnitude with increasing wavelength above 1100 cm −1 , but the simulated spectrum increases more.Below about 900 cm −1 the spectral behaviour of the measured and simulated spectra differs.This may be due to either wrong assumptions about the ash type and hence refractive index and/or the mixing of ice with ash.Ice clouds have an opposite effect of ash clouds on the brightness temperature between 800 and 1000 cm −1 , whereas above 1075 cm −1 ice clouds have only a very weak dependence on wavenumber (see Fig. 2 of Gangale et al., 2010).To test if the presence of both ash and ice could reproduce the measured spectrum, simulations were made with both an ash cloud and an ice cloud.The altitude and thickness of the clouds were as above, but the ash cloud density was 2 × 10 −4 g m −3 and the ice water content 1.5 × 10 −2 g m −3 .The resulting spectrum is shown in maroon in the bottom plot of Fig. 8.The mixed scene with both ash and ice is seen to well reproduce the measured ash spectrum in the top plot of Fig. 8.

Simulated satellite image
Figure 9 shows a simulated satellite image (top) and the corresponding observation (bottom).Three visible channels of the SEVIRI (Spinning Enhanced Visible and Infrared Imager) instrument on the MSG satellite were simulated based on input data from the operational COSMO-DE forecast (Baldauf et al., 2011) of Deutscher Wetterdienst for the 15 July 2012, 12:00 UTC.The spatial resolution of the simulation is 2.8 km × 2.8 km; the SEVIRI observation is 3 km × 3 km at the sub-satellite point.A false colour composite was generated using the simulated radiance of the 1.6 µm channel for red, the 0.8 µm radiance for green, and 0.6 µm radiance for blue.The simulations were performed using the 1-D disort solver.The MODIS surface albedo data set was used (Schaaf et al., 2002) to set the Lambertian surface albedo.The effective radii of liquid clouds were parameterized according to Martin et al. (1994), and for the optical properties the mie parameterization was applied.Ice cloud effective radii were parameterized according to Wyser (1998) and for the corresponding optical proper- ties the parameterization baum_v36 was used with the general habit mixture.Molecular absorption was included using the reptran parameterization.In the false colour composite water clouds appear white and ice clouds appear blueish, because ice absorbs in the region of about 1.6 µm.The simulated image looks very similar to the observation.A major difference is that the ice clouds in the observation appear more blueish, the reason is that their real optical thickness is larger than in the COSMO-DE forecast.

Polarization
The MYSTIC solver can be applied to simulate multiangle multi-spectral polarized radiances using the option mc_polarisation (Emde et al., 2010).Polarized radiative transfer using MYSTIC has been validated in extensive model intercomparison projects (Kokhanovsky et al., 2010a;Emde et al., 2015).
Figure 10 shows an example for simulations at wavelengths of 443, 670, and 865 nm; these are measured by the POLDER (Polarization and Directionality of the Earth's Reflectances) instrument onboard PARASOL (Deschamps et al., 1994).All simulations are for a solar zenith angle of 30 • and show the reflected radiances (normalized to incom- ing solar irradiance) at the top of the atmosphere in the solar principal plane.The viewing angle of 30 • corresponds to the exact backscattering direction.The angular resolution is 2 • .All simulations are for the US-standard atmosphere.The figure shows the first and second components of the Stokes vector I and Q; the components U and V are exactly 0 in the principal plane for symmetry reasons.
The first row shows the results for a clear atmosphere, i.e.Rayleigh scattering and molecular absorption.Here I is largest for the shortest wavelength because the Rayleighscattering cross section decreases with λ −4 , where λ is the wavelength.The absolute value of Q also increases with increasing Rayleigh-scattering cross section.A negative Q means that Rayleigh scattering polarizes perpendicular to the scattering plane, which, for single scattering, corresponds to the principal plane for this geometry.
The second row of the figure shows the same simulation but with an underlying ocean surface, which is modelled according to Mishchenko and Travis (1997) (bpdf_tsang).The wind speed was set to 2 m s −1 .I and Q clearly show the sun glint, which has a maximum at a viewing angle of about -30 • and is highly polarized.The intensity of the sun-glint increases with increasing wavelength since the incoming radiance at the surface becomes less diffuse when there is less Rayleigh scattering in the atmosphere.
The third row shows the result for desert aerosol as defined in the OPAC database (aerosol_species_file desert), with an underlying Lambertian surface albedo of 0.3.I shows a backscatter peak at 670 and 865 nm.Q looks similar to Rayleigh scattering; however, there are differences mainly around the backscatter region.At wavelengths of 670 and 865 nm, Q has a minimum in the exact backscatter direction and becomes positive for viewing angles around this direction.
The fourth row shows a simulation including a water cloud (wc_properties mie) in 2-3 km altitude with an optical thickness of 10 and an effective droplet radius of 10 µm.I and Q show the glory about the backscatter direction and the rainbow at a viewing angle of about −10 • corresponding to a scattering angle of 140 • .In Q the rainbow is more pronounced than in I because Q is less affected by multiple scattering.The angular resolution shown here is not sufficient to separate the glory from the backscattering peak in I .The sign of Q in the rainbow region is the same as for Rayleigh scattering, whereas it is opposite in the glory region, which means that the rainbow is polarized perpendicular to the scattering plane, whereas the glory is polarized parallel to the scattering plane.
The last two rows show simulations with ice clouds, where we have used the yang2013 parameterization.An ice cloud layer with an optical thickness of 2 was included at an altitude from 9-10 km.The selected habit was solid_column and we performed simulations for smooth crystals and for severely rough crystals.The effective crystal radius in both simulations is 30 µm.The smooth crystals show a backscatter peak in I and a positive Q about the backscatter direction.Also there are some smaller features in I and Q.The radiances (I and Q) for rough crystals are smooth functions of viewing angle.This different behaviour has been used to determine the fraction of smooth crystals in ice clouds from POLDER measurements (Cole et al., 2014).

nm
Figure 11.Twilight radiance at 500 and 700 nm calculated using fully spherical geometry for the US-standard atmosphere.The lines are for different solar depression angles.The x axis corresponds to the viewing zenith angle.

Fully spherical geometry
MYSTIC can be operated in fully spherical geometry (mc_spherical 1D).The implementation of 1-D spherical geometry is described in Emde and Mayer (2007) where it has been used to simulate radiation in the umbral shadow of a solar eclipse.A comparison to measurements during the total eclipse in Greece in March 2006 (Kazantzidis et al., 2007) showed a very good agreement for modelled and measured UV irradiance, which decreased during totality by 2 to 3 orders of magnitude depending on wavelength.
Fully spherical geometry has also been used to simulate actinic fluxes at high solar zenith angles up to 92 •  (Suminska-Ebersoldt et al., 2012).
Another interesting application is the simulation of polarized radiance at the surface at twilight, because polarized radiance measurements at twilight can be used to retrieve aerosol optical properties (e.g.Saito and Iwabuchi, 2015).
As an example we calculated polarized clear-sky radiances for solar depression angles up to 9 • for the US-standard atmosphere and default Rayleigh-scattering and absorption settings.Figure 11 shows the result as a function of viewing zenith angle.The relative azimuth angle between sun and observer is 0 • , which means that the observer looks into the direction of the sun.We see that the intensity decreases by about 4 orders of magnitude for solar depression angles between 0 • (sun at horizon) and 9 • (sun 9 • below horizon).The degree of polarization (not shown) at a viewing angle of 5 • is more than 90 %.All results agree to published results by Blättner et al. (1974), which indicates that fully spherical geometry works correctly in MYSTIC.

Summary
We have presented the libRadtran software package (version 2.0.1), which is a comprehensive and powerful collection of tools for radiative transfer simulations of the Earth's atmosphere.It is user-friendly, well-documented, and is widely used in the scientific community.We have described various new features and parameterizations, which have been included after the first publication of libRadtran in 2005.New features are for example a vector radiative transfer solver and a solver for rotational Raman scattering.The package includes state-of-the-art parameterizations for aerosol and ice cloud optical properties and a newly developed efficient absorption parameterization.angle, respectively.Planck's function at a given frequency ν is represented by B ν (τ ) and its temperature dependence is contained implicitly in τ .
The first term on the right-hand side in Eq. (B1) represents the contribution of the surface and the second one the contribution of the atmosphere.The surface contribution can be written as

Figure 2 .
Figure2.Nadir top-of-the-atmosphere radiance in the oxygen-A band around 760 nm (left) and in the IR (infra-red) window region (right) for the midlatitude-summer atmosphere ofAnderson et al. (1986).All calculations were performed with the MYSTIC solver using the "absorption lines importance sampling" method(Emde et al., 2011).(Top) high spectral resolution calculation, based on line-by-line absorption cross sections calculated using ARTS(Eriksson et al., 2011); (bottom) pseudo-spectral calculations using the representative wavelengths band parameterizations (reptran) with different resolutions and lowtran.For comparison see also Fig.3inMayer and Kylling (2005), which shows transmittances for genln2 line-by-line calculations and lowtran for the same spectral regions.
Figure3.Phase matrix elements for the basic OPAC aerosol types "water soluble" (waso), "sea salt accumulated mode" (ssam), and soot, for a water cloud with a droplet effective radius of 10 µm, and for Rayleigh scattering (with δ = 0) at a wavelength of 350 nm.θ is the scattering angle, i.e. the angle between incoming and scattered directions.

Figure 4 .
Figure 4.Absorption (left)  and scattering (right) optical thickness for various aerosol mixtures specified using the option aerosol_species_file.The aerosol optical properties as well as the mixtures have been generated based on OPAC(Hess et al., 1998) parameters.

Figure 5 .
Figure5.Phase matrix elements of ice crystal distributions with an effective radius of 40 µm at 550 nm wavelength.The red lines correspond to smooth and the blue lines to severely rough crystals.The individual habits (solid-column, column-8elements and plate) are for the parameterization yang2013, and the general habit mixtures (ghm) are for hey including smooth crystals and baum_v36 including severely rough particles.

Figure 6 .
Figure 6.Screenshot of the graphical user interface for a spectral high-resolution simulation of the O 2 -B band including a fluorescence source.Plots of input and output data are included together with the help information for one option.See text for further explanation.

Figure 7 .
Figure 7. (Left) the transmittance from ARTS output and radiance from uvspec.(Right) the top of the atmosphere nadir viewing radiance in the O 2 -B band with (black line) and without (cyan line with circles) a surface fluorescence source (red line with triangles).The radiances have been convolved with a spectral response function with FWHM of 0.3 nm.

Figure 8 .
Figure 8. (Top plot) brightness temperature spectra for different locations as measured by IASI on 15 February 2014, 02:33 UTC, during the Mt Kelud, Indonesia, eruption.Tentative classification of the spectra is given in the legend.See text for details.(Bottom plot) simulated brightness temperature spectra using ARTS/uvspec.The atmospheric processes included in the simulations are given in the legend.

Figure 9 .
Figure 9. (Top) simulation of MSG-SEVIRI image.False colour composite, where red corresponds to the 1.6 µm channel, green to 0.8 µm, and blue to 0.6 µm.The simulation was performed using the disort solver with input data from the operational COSMO-DE forecast for the 15 July 2012, 12:00 UTC.The axes correspond to SEVIRI pixel.(Bottom) corresponding SEVIRI image.

Figure 10 .
Figure10.Stokes vector components I and Q at wavelengths of 443 nm (blue solid lines), 670 nm (green dashed lines), and 865 nm (red dashed-dotted lines) for various atmospheric setups (see text for details).The radiances are calculated at the top of the atmosphere for viewing angles from −50 to 50 • , where 0 • corresponds to the nadir direction.
I ν (τ * , µ, φ) = s B ν (τ * ) term representing the emission of the surface ( s = surface emissivity) and the second one the reflection at the surface of the radiation emitted by the atmosphere toward the surface.The factor 2 comes from the integration over the azimuth angle φ.Under the approximation of Planck's function B ν (τ ) as a piecewise linear function in τ between two consecutive levels, both integrals can be solved as a function of the exponential integral Ei(x)

Table 1 .
The radiative transfer equation solvers currently implemented in libRadtran.

Table 3 .
Absorption cross section data included in libRadtran; the non-default parameterizations are put in parentheses.

Table 5 .
Ice cloud parameterizations in libRadtran

Table 6 .
The surface reflection models currently implemented in libRadtran.