Retrieval of sulfur dioxide from a ground-based thermal infrared imaging camera

Recent advances in uncooled detector technology now offer the possibility of using relatively inexpensive thermal (7 to 14 μm) imaging devices as tools for studying and quantifying the behaviour of hazardous gases and particulates in atmospheric plumes. An experimental fast-sampling (60 Hz) ground-based uncooled thermal imager (Cyclops), operating with four spectral channels at central wavelengths of 8.6, 10, 11 and 12 μm and one broadband channel (7– 14 μm) has been tested at several volcanoes and at an industrial site, where SO2 was a major constituent of the plumes. This paper presents new algorithms, which include atmospheric corrections to the data and better calibrations to show that SO2 slant column density can be reliably detected and quantified. Our results indicate that it is relatively easy to identify and discriminate SO 2 in plumes, but more challenging to quantify the column densities. A full description of the retrieval algorithms, illustrative results and a detailed error analysis are provided. The noise-equivalent temperature difference (NE1T ) of the spectral channels, a fundamental measure of the quality of the measurements, lies between 0.4 and 0.8 K, resulting in slant column density errors of 20 %. Frame averaging and improved NE 1T ’s can reduce this error to less than 10 %, making a stand-off, day or night operation of an instrument of this type very practical for both monitoring industrial SO2 emissions and for SO 2 column densities and emission measurements at active volcanoes. The imaging camera system may also be used to study thermal radiation from meteorological clouds and the atmosphere.


Introduction
The thermal infrared (3 to 15 µm) region of the electromagnetic spectrum contains several sub-regions which can be exploited for studying atmospheric gases (e.g.Esler et al., 2000).Notable among these are the window regions between 3 and 4 µm, which is often referred to as the mid-infrared (MIR), and 7 and 14 µm, which is referred to as the thermal infrared (TIR).The MIR is used for identifying "hot spots", localised regions of anomalously hot pixels in satellite measurements (Wright et al., 2004).The MIR can also be used from the ground or on airborne platforms to image the heat from forest fires (Lentile et al., 2006) or hot gases rising from volcanic vents (Francis et al., 1995) and to map temperatures in plumes (Sawyer and Burton, 2006) and on lava fields (Realmuto et al., 1992).The TIR has been used less frequently to study volcanic processes.This is largely due to the fact that sensitivity in this region peaks at terrestrial temperatures of 300 K, much lower than the temperature of a typical "hot spot" or volcanic heat source, and because, until recently, thermal imagers operating in the TIR required expensive active detector-cooling systems (nitrogen Dewars or Stirling cycle coolers) to achieve good signal-to-noise performance (Derniak and Boremann, 1996).TIR instruments on satellites do use active cooling systems and in these cases the image data are used to monitor volcanic eruption clouds and discriminate them from meteorological clouds for aviation hazard warnings and for gas measurements (Prata, 2009).Pugnaghi et al. (2002) used the Multi-spectral Infrared and Visible Imaging Spectrometer (MIVIS) on board an aircraft to map the SO 2 emissions from Etna.Their algorithm was based on a split-window formulation using channels centred at 8.74 and 9.56 µm to eliminate the effects of water vapour and determine SO 2 abundance.Realmuto et al. (1994Realmuto et al. ( , 1997) ) showed that SO 2 could be determined from the multi-channel TIR imager Advanced Spaceborne Thermal Emission And Reflection Radiometer (ASTER), on board the Earth Observing System (EOS) Terra satellite, by using detailed radiative transfer calculations to account for water vapour and surface emissivity variations.
All of the studies described above have used passive thermal sensing, relying on emission or absorption by the gas to provide a signal to measure.Measurements can also be made in absorption mode by using the Sun as a source or by providing an artificial source of radiation (typically a globar and retroreflector).In these applications single field-ofview (FoV), medium-spectral-resolution (6-0.5 cm −1 ) interferometers are used to gather quantitative information on multiple gas species simultaneously.Fourier-transform interferometers (FT-IRs) have become a very valuable device for volcanic gas studies (Love et al., 1998;Oppenheimer et al., 1998;Burton et al., 2000;Horrocks, 2001), including measurements of gas ratios reported by Oppenheimer et al. (2002).Systems using ultraviolet light as a source have recently been developed for volcanic SO 2 measurements (Mc-Gonigle, 2005;Horton et al., 2006), for volcanic BrO measurements (Bobrowski et al., 2003), and also for CO 2 slantpath columns (Goff et al., 2001).More recently Stremme et al. (2013) and Krueger et al. (2013) presented measurements of volcanic emissions using a scanning FT-IR, showing two-dimensional visualisations of SO 2 based on thermal emission spectroscopy.Kinoshita et al. (2003) used a groundbased CCD imager together with a near-infrared filter to study volcanic plumes, but they did not attempt a quantitative retrieval of the gases or particulates.Notsu et al. (2003) demonstrated the feasibility of using the 8.6 µm waveband for the measurement of volcanic SO 2 slant column density using a portable spectral infrared radiometer.
This paper presents the first detailed study of the use of a ground-based, uncooled thermal imaging microbolometer radiometer to detect and quantify SO 2 gas from volcanic and industrial sources.The intention for this work was to develop a multi-filter TIR imaging camera capable of sensing gases and particles, principally for applications in volcanology.The details concerning the methods for detecting volcanic ash particles have been provided in a separate paper (Prata and Bernardo, 2009); here we concentrate on the SO 2 gas retrieval methodology.The capability to acquire frequent, real-time images from a fixed platform (e.g.located at a volcanological observatory near to an active volcano, or during a field deployment) day or night offers a practical and safe tool for understanding some aspects of volcanic activity.
The organisation of the paper is as follows: we briefly describe the principal characteristics of uncooled microbolometer thermal imaging devices and then show how such cameras can be adapted for use in detecting and quantifying SO 2 gas emissions.The design of the camera system is described, and the basic theory presented, for SO 2 slant column density (hereafter referred to as SCD) retrieval and then illustrated by showing how SO 2 emissions from an industrial stack can be derived.This is followed by a detailed error analysis of the retrieval scheme.Measurements made at two volcanoes, Etna (Sicily, Italy) and Stromboli (Aeolian Islands, Italy), are provided to show how estimates of volcanic SO 2 emission rates can be estimated.We conclude with comments on how this technology might be improved by integrating it with other remote sensing instruments, for example ultraviolet (UV) spectrometers, and used for quantitative studies of volcanic emissions, for detecting hazards from an airborne platform, and for alerting authorities of volcanic activity during the day and night for hazard warnings.

Thermal imagers
In the last 10-15 years great advances have been made in manufacturing bolometers of high sensitivity (Kruse, 2001).The detectivity of these devices is background-limited, and they are often referred to as background-limited infrared photodetectors (BLIP) devices.The use of silicon semiconductors (silicon nitride substrate with vanadium oxide detecting material) for manufacturing arrays of bolometric detectors has greatly reduced the cost of the production of thermal imaging cameras.These microbolometers, typically consisting of 10 4 -10 6 elements, are sensitive to radiation in the wavelength range of 7-14 µm and operate at 30-60 Hz (Kruse, 2001).Thermal cameras are commercially available with temperature sensitivities of ∼ 50 mK (7-14 µm), array sizes of 320 × 240 pixels (or larger), F1.0 optics and 60 Hz operation.Thus, in principle, a camera of this kind can acquire images showing temperature changes of less than 0.1 K at a rate of tens of frames per second.In practice this is difficult to achieve because of the presence of noise (1/f , background and internal temperature fluctuations and Johnson noise), non-uniformity of the array, the need for calibration, and frame integration.Other factors may also limit achieving the ideal image capture rate: for example, extracting the image frame data rapidly requires fast electronics and a good microprocessor and communications hardware and software.Shaw et al. (2005) describe an uncooled thermal imaging camera for use in atmospheric studies.This camera has a single passband (∼ 8-14 µm) and is used to view the sky overhead for studies of clouds.They report the calibration error of this instrument to be 0.5 W m −2 sr −1 , or about 2 % of the ambient radiance, and also show that the microbolometer is sensitive at low temperatures (< 240 K) by observing changes in signals for very high thin clouds (cirrus).The camera developed here incorporates wavelength selection (filters), and this decreases the sensitivity and adds time delays to the image capture.The camera needs to be sensitive at temperatures of 230-300 K, which cover the typical range of atmospheric plume temperatures.In this paper we describe, in general terms, the thermal imaging camera and its modifications for use as an apparatus for quantifying atmospheric SO 2 gas emissions, and provide an error analysis with some illustrative examples of data captured at an industrial stack and at two erupting volcanoes.

Cyclops -a multifilter thermal infrared camera system
A commercial off-the-shelf (COTS) thermal infrared (IR) camera with 50 mK noise-equivalent temperature difference (NE T ) and a single broadband filter covering the IR wavelength region from 7 to about 14 µm was adapted for use in this work.A schematic of the principal components of the camera, "Cyclops", is provided in Fig. 1.
Radiation enters through the fore-optics (1), which defines the FoV and is focussed onto the microbolometer array (2).Readout electronics (3) convert the sensed radiation signals into voltages, and a digitiser (4), microprocessor (5) and communications port (6) deliver the IR microbolometer array output to a computer for image display and further analysis.For use in detecting and quantifying atmospheric gases at typical atmospheric temperatures, several modifications to the COTS camera are needed.Figure 2 shows photographs of the camera and its main components.Figure 2a illustrates the filter set-up; Fig. 2b shows the camera body -the larger diameter housing holds the filter wheel, filters and the blackbody shutter; Fig. 2c shows the camera mounted on the deck of a ship with a calibration rig attached and an external shutter used to verify the internal shutter calibration.It is important to note that the design concept requires that the filter be placed behind the lens and the shutter (and any external blackbodies) be placed in front of the lens.This arrangement ensures that radiation from the lens is properly accounted for in the calibration.
The two most important modifications of the COTS camera are described below and have been incorporated into an operational camera dubbed "Cyclops".

Filtering
Spectral selection of radiation into narrow bands (0.5-1.0 µm) is achieved by placing a filter wheel between the fore-optics and detector.The filters are carefully selected to match pre-determined specifications for optimal sensing of SO 2 gas and particles.Figure 3 shows the line intensities from the HITRAN-2000 database (Rothman, 2003) illustrating the main absorption features of SO 2 in the region 6.8-10 µm.
The strongest feature at 7.3 µm is not suitable for groundbased sensing of SO 2 because water vapour absorption dominates in this region.The atmospheric transmittance for slant paths with a zenith angle of 75 • and ranges of ∼ 38 and ∼ 6 km, calculated using MODTRAN at a resolution of 5 cm −1 , are also shown on the figure.The 7.3 µm channel (C1) is opaque and hence unsuitable for ground-based use.FIGURE 3. HITRAN line intensities of the two main SO2 absorption bands shown together with the relative response functions of the five "Cyclops" channels (filters) used.The slant-path transmittance between the camera and target at ranges of ∼38 km (green line) and ∼6 km (red line) at 5 cm −1 resolution calculated from MODTRAN over the region 700-1400 cm −1 (7-14 µm) for a standard US atmosphere are also shown.Note that the 7.3 µm channel (C1) is opaque and the influence of the 9.6 µm O3 band on the 10.1 µm channel (C4) decreases for shorter ranges.

35
Figure 3. HITRAN line intensities of the two main SO 2 absorption bands shown together with the relative response functions of the five "Cyclops" channels (filters) used.The slant-path transmittances between the camera and target at ranges of ∼ 38 km (green line) and ∼ 6 km (red line) at 5 cm −1 resolution calculated from MODTRAN over the region 700-1400 cm −1 (7-14 µm) for a standard US atmosphere are also shown.Note that the 7.3 µm channel (C1) is opaque and the influence of the 9.6 µm O 3 band on the 10.1 µm channel (C4) decreases for shorter ranges.
The 10.1 µm channel (C4) is affected by ozone (absorption centre at 9.6 µm), but this effect diminishes with distance to the target (the plume).The feature at 8.6 µm, although less strong, is better suited for SO 2 sensing because water vapour absorption is much reduced compared to at 7.3 µm 1 .Cyclops is restricted to measuring gases that have broad (∼ 1 µm or larger) absorption features within the region 7-14 µm, because of signal-to-noise considerations.Another volcanic gas that meets this criterion is CO 2 , but because of the relatively high abundance of CO 2 in the ambient atmosphere, it is problematic to measure this gas using thermal IR ground-based radiometry.
The design of Cyclops was heavily influenced by knowledge of atmospheric gas and particle absorption characteristics (see for example Gangale et al., 2010), and constrained by current technology.Table 1 shows the Cyclops channels (or filters) chosen for detecting SO 2 and volcanic ash from the ground, and Fig. 3 shows the filter response functions for these channels.

Calibration
Gas and particle discrimination and quantification requires high-fidelity thermal images from Cyclops.To achieve reliability and accuracy, the camera must be calibrated.The procedure is a linear calibration requiring an estimate of the gain and intercept that converts the digital numbers (DNs) to radiances and then to brightness temperatures.A twostep process is implemented: Cyclops is first calibrated in the laboratory under controlled conditions using a blackbody source.Estimates of the gains and intercepts for all channels are determined for a variety of environmental and target (source) conditions.The temperature of the focal plane array (FPA) is also recorded and stored with the data.The FPA temperature is used as a surrogate to correct for radiation from the camera itself and a radiance correction is added to the calibration equation.In the field, environmental conditions cannot be measured with sufficient accuracy to allow sole use of these calibration coefficients.Thus a second step is employed that compensates for changes in the environmental conditions -specifically, the temperatures of the instrument, fore-optics and outer housing.This second step requires the addition of a blackbody shutter, placed in front of the fore-optics, filter wheel and detector.The temperaturecontrolled shutter moves in front of the camera on computer command so as to allow for a single calibration point on the DN-radiance calibration line.The calibration can be repeated as frequently as required and is performed for each of the five filters separately.This two-step procedure gives temperature precisions of 0.2 to 0.7 K at 280 K, depending on channel.Water vapour is typically the largest absorber and emitter of radiation within the Cyclops waveband.Viewing from the ground exacerbates the problem of water vapour absorption and emission because the concentration is largest near the surface and decreases rapidly (exponentially) with increasing height above the surface.At low-elevation viewing angles (high zenith angles), the water vapour pathlength, the product of the water vapour amount and geometrical pathlength, can be large and hence have a significant effect on the measured IR radiation.Furthermore, water vapour absorbs differentially across the waveband, with greater absorption (and emission) occurring at 12 µm than at 11 µm.Since Cyclops views the water vapour against a sky background that is usually colder than the foreground, in the absence of other absorbers (e.g.clouds), Cyclops measures more radiation at 12 µm than at 11 µm.As an example, Fig. 4 shows a series of Cyclops images obtained at a location where no ash or SO 2 was present.The images consist of raw, uncalibrated measurements and their respective histograms (Fig. 4a, b; left-most panels), calibrated temperature images and their respective histograms, and the right-most panels (Fig. 4c, d,  e) show temperature difference histograms for various combinations of Cyclops channels.These images confirm the general comments above: measured radiation increases with wavelength within the 11-12 µm waveband and decreases with increasing camera elevation.The histograms show two distinct peaks; the broad peak covering 230-260 K is due to sky radiation (water vapour and CO 2 ) and the smaller peak centred near 280-290 K is due to radiation from trees captured in the lower left-hand corner of the images.The difference histograms also show that radiation at 10 µm is larger than at 11 µm and at 12 µm, and larger still at 8.6 µm.This is due to the general shape of the water vapour absorption curve between 8 and 12 µm with absorption highest at 8 and 12 µm.The difference histograms for natural objects (e.g.trees and vegetation) are centred near to 0 K difference, the main effect being due to emissivity effects of trees and vegetation.
The images in Fig. 4 also indicate the general trend of decreasing radiation (at all wavelengths) with increasing viewing elevation angle.The rate of decrease with elevation angle is not the same at all wavelengths, and the atmosphere induces a differential absorption effect that depends on viewing angle.The importance of calibrating the images is also apparent.The significant warm patch appearing in the centre of the filtered raw images is caused by unwanted radiation from the lens and housing of Cyclops.There is also a "blooming" effect apparent at the left and right edges of the uncalibrated data, which has been largely removed in the calibrated data in the right edge, but is still partially apparent at the left edge.Finally, it can be seen that image noise is higher at 8.6 µm and lowest in the broadband (lowest panel) image.
These general observations lead to two very significant conclusions regarding the subsequent processing of the Cyclops data.Firstly, raw, uncalibrated data are virtually of no value for identifying gases or particulates in these filtered thermal IR images.Since much of the useful information is contained in difference images, reducing noise and applying a consistent and accurate calibration appear to be fundamental to transforming the data into information.Secondly, we note the strong affect water vapour has on the measurements.Applying an atmospheric correction is crucial to correctly identifying gases and particulates in the images.Furthermore the correction must be applied with a dependence on viewing angle and preferably on a pixel-by-pixel basis.The methodology and results presented here are new and are an improvement to the methodology previously reported by Prata et al. (2004).The atmospheric correction and retrieval procedures are described next.

Quantifying SO 2
The Cyclops camera system was designed to use up to five spectral filters, chosen to optimise the detection of specific atmospheric gases.To quantify SO 2 SCDs from the ground, a filter with a narrow waveband centred near 8.6 µm was selected.The filter response function is plotted in Fig. 5 together with the SO 2 absorption coefficient measured by NIST (National Institute of Standards and Technology) (Chu et al., 1999).
The ground-based thermal imager can view a plume from a volcanic source or from an industrial stack at elevation angles of 10 • or less (zenith angles > 80 • ).The preferred arrangement for Cyclops is with a high elevation angle in order to reduce the effects of water vapour absorption along the path.The camera has a FoV of ∼ 32 • , and the total azimuthal angular variation is similar to the total zenithal variation.In the following analysis each pixel is treated independently of all others and there is a simple mapping between image column and line numbers and azimuth angles and image elevation.
The radiation measured at the imager can be described by three terms, where θ is elevation angle, φ is azimuth angle, i is channel number, and the superscripts refer to foreground radiance (f), background (b), and plume radiance (p).The plume radiance consists of emitted radiation, as well as radiation from the atmosphere that has been attenuated as it traverses through the plume.Scattering is ignored.The plume is considered to be sufficiently opaque that most of the background radiation is blocked by the plume, but in the retrieval scheme it is necessary to consider regions outside the plume where the sum of the background and foreground radiation is denoted as I o i (see Eq. 26 later in the text).The channel radiances are integrations over the channel filter response functions for each pixel within the two-dimensional (2-D) image space.Background radiance refers to radiance from the sky, behind the plume; foreground radiance refers to radiance emanating from the atmosphere between the plume and the imager.In general it is a difficult task to estimate the atmospheric terms I f i and I b i from observations.The goal of this analysis is to isolate the plume radiance term and then estimate the product of the gas concentration and plume thickness.The model used assumes no scattering and that variations in the absorption coefficient of the medium are invariant along the absorption path.Furthermore, the plume is assumed to be plane-parallel and governed by Schwarzschild's radiative transfer equation.The next section provides the mathematical details of the analysis.The resulting equation that is used to retrieve the pathlength concentration amount m * , the product of the absorber density with the pathlength, is stated here and some general remarks are made.
where i,j is an effective emissivity of the plume and is given by and k is the absorption coefficient averaged over the response function of the measurement channel; all other terms are brightness temperature differences ( ) and are defined in the Appendix.The retrieval procedure uses three of the imager's five channels: the 8.6, 10 and 12 µm channels.The information regarding SO 2 in the plume is contained in the 8.6 µm channel, while the 12 µm channel is used to correct for atmospheric effects and the 10 µm channel, which is the most transparent to water vapour absorption, is used to estimate the plume temperature.The retrieval scheme uses temperature differences.Most important of these are the thermal contrast, the temperature difference between the plume and the background atmosphere, and terms involving differences between the spectral brightness temperature, with and without the plume present, and brightness temperature differences between the 8.6 and 12 µm channels.For highly opaque plumes, these differences may be small and the retrieval scheme becomes unstable.For very thin plumes the thermal contrast is low and the retrieval becomes noiselimited.
The sub-section on error analysis (Sect.4.2) provides details on the accuracy of the retrieval scheme and the section following that illustrates the results of using the scheme at several different sites.

Retrieval algorithm
We consider a plane-parallel plume (slab) with thickness d consisting of a homogeneous mixture of two gases with densities ρ 1 and ρ 2 .The absorption coefficients of the gases, k 1 and k 2 , are assumed not to vary within the slab, and radiation is assume to be attenuated by absorption and emitted at a constant plume temperature, T p , but not scattered.In the infrared between 7 and 13 µm wavelengths, scattering is typically much less important than absorption and emission.The camera views the plume in up to five narrowband channels, denoted by i, i = 1,5, and we assume that all quantities (e.g. the radiances and the absorption coefficients) are averages over the channel filter responses.The measurements are also regarded as simultaneous: although this is not strictly true, the actual time difference between images varies depending on which channels are being acquired but is typically a few minutes.The coordinate system adopted is Cartesian, with the leading side of the plume placed at y = 0, the camera placed at x = 0, y = L, z = 0 and the coordinates x and y represent the horizontal axes and z the vertical axis as shown in Fig. 6.
The camera views the plume from a distance R, measured from the centre of the detector to the side of the plume closest to the camera, and at an elevation angle θ n and azimuth angle φ n , which vary with camera pixel number n.In this coordinate system the camera line C l and column C c numbers are related to the camera elevation and azimuth angles through where L is the distance to the plume measured in the x-y plane (z = 0), ζ is the elevation of the camera measured from ground level (height above mean sea level) to the first line of the image, s n is the size of image pixel n, and the image has N c columns by N l lines (320 × 240 in the current setup).The camera is oriented such that an azimuth angle of φ n = 0 corresponds to the centre of the image, or column number N c /2. Pixel numbers are counted from the bottom left of the image, with line 1, column 1 corresponding to pixel number 1, and the last column of the top line corresponding to pixel number N c N l .The pixel size varies with distance from camera to target and can be determined from where F is the focal length of the camera, χ is the pitch of the pixel on the detector chip (∼ 45 µm), and l,c is the FoV of the microbolometer detector array in the vertical ( l ) or horizontal ( c ), and we use subscripts for the pixel size s l,c to denote the size in the line or column directions.The radiation, I i (θ n ), measured by the camera for pixel n and channel i in the direction θ n for this situation is governed by Eq. ( 1).The measured radiation, assumed to arise from radiation along the path R, from the background and from plume radiation, is due to variations in absorbers ρ 1 (θ n , φ n ), ρ 2 (θ n , φ n ) and temperature T (θ n , φ n ), as well as absorption and emission by other well-mixed gases (e.g.CO 2 , CH 4 , N 2 O, O 3 ) that are assumed invariant.The foreground and background radiation can be calculated from the MODTRAN-4 radiative transfer model (Berk et al., 1999) using a nearby radiosonde profile for water vapour (ρ 1 ) and temperature and assuming climatological values for the well-mixed gases.However, an alternate procedure which makes better use of the camera measurements has been adopted.The retrieval uses the difference between the radiation measured by the camera in a channel centred at 12 µm, where there is no SO 2 absorption and some H 2 O absorption, and a channel centred at 8.6 µm, where there is considerable SO 2 absorption and some H 2 O absorption.The 12 µm channel is chosen in preference to a channel at 11 µm or 10 µm because of the concave shape of the water vapour absorption curve from 8 to 12 µm, with absorption greatest at 8 and 12 µm, and lowest at 10 µm.In the ideal case, when the absorption is the same at 8.6 and 12 µm, taking a difference leaves only the contributions from absorbers ρ 2 (SO 2 ) and a smaller contribution from ρ 1 (H 2 O) within the plume.The radiative transfer is divided into two parts: first we analyse the radiation through the plume and treat this as an absorption-emission process.
Next we treat the radiation from the foreground as equivalent to a blackbody radiating at a representative temperature, and attenuated by equivalent transmission functions due to the absorbers.In the case of an opaque plume, the background radiation can be ignored, but we treat this later when radiation near the plume, but not obstructed by it, is considered.This simplified treatment is justified on the basis that we are not interested in the details of the structure of the foreground and background radiation fields but only in their effects as a perturbation on the plume radiance, which is of much greater interest.Schwarzschild's equation for the azimuthally independent plume radiance for one pixel and one channel may be written as and where r is distance along the plume in the direction of θ n , B i is the Planck function, i is the channel number, and T p is the plume temperature (assumed not to vary along the path).This equation can be integrated along the path to yield and where I o i is the radiation from the atmosphere in the direction r, τ i (r 1 , r) is the optical thickness of the plume between r and r 1 , and |r 1 − r| is the pathlength traversed by the radiation within the plume in the direction r.We now assume that the path is homogeneous, k does not vary with position in the plume, and the plume is in thermodynamic equilibrium.Equation (10) shows that the plume radiation measured by channel i consists of terms representing absorption attenuation by the plume and emission from the plume along the path.For two absorbers, Let m j = ρ j d sec θ n sec φ n .
The plume thickness in the r direction is d sec θ n sec φ n , and hence We may write a similar equation for a channel which is unaffected by absorber ρ 2 , The radiances (the measurements) are made at different wavelengths and converted to brightness temperatures so that channel differences can be taken.We use a Taylor series approximation to linearise these equations and then combine them to solve for m 2 .Linearisation of the radiances around a mean temperature has been used by other authors (McMillin and Crosby, 1984) and is a reliable approach provided the radiances i and I o j are all similar.For a plume in thermodynamic equilibrium with the atmospheric environment and for viewing at low elevation angles (θ n < 60 • ), the radiances will be similar.Linearising around an atmospheric radiance (I o i ) unaffected by the plume, and denoting brightness temperatures by the symbol (to avoid confusion with the Planck function), Similarly, Using Eqs. ( 15)-( 17) and substituting for the radiances gives, after some manipulation, The same is the case for the channel with two absorbers: Let e −k i,1 m 1 = e −k j,1 m 1 .This assumption requires that the transmission by water vapour is equal at the two wavelengths chosen, viz.8.6 and 12.0 µm.The section on error analysis examines the efficacy of this approximation.Using this approximation we have Subtracting Eq. ( 18) from Eq. ( 20), and after some tedious algebra, we have where Equation (22) shows that the retrieval of the plume emissivity depends mainly on the plume temperature difference between the two channels but also on the thermal contrast between the plume and the atmosphere outside the plume ( .We now show how the plume and atmosphere brightness temperatures are related to the plume, foreground and background radiances, and how the brightness temperatures are determined for use in Eq. ( 22).
Consider two measurements, one made through the plume and the other without the plume in the FoV.Assuming that the atmosphere does not change appreciably between these two measurements, for the first measurement (dropping reference to angles) we may write and for the second measurement The superscript o refers to atmospheric radiation "outside" the plume.Each of these quantities may be determined by solving integrals of the form In general we do not have information on the path variation of the absorption coefficient, the absorber or the temperature.Let the transmittance of each path be designated τ f i,q , τ p i,q and τ b i,q for the foreground, plume and background, respectively.As before i represents channel and q absorber type (q=1,2).Let the temperatures of the layers be T f , T p , and T b , respectively, and we replace the path integrals with mean radiances, denoted by an overbar.Then, Note that we have assumed that the foreground and background atmospheres have not changed between the measurements and that they contain no SO 2 (absorber q = 2).Subtracting Eq. ( 28) from Eq. ( 29), A similar equation can be obtained for a second channel j , which has no absorption due to absorber q = 2, Subtracting Eq. ( 31) from Eq. ( 30), where The quantities in Eq. ( 33) are all measurable, and hence Eq. ( 32) can be solved after the correction δI o i,j has been applied and the brightness temperature analogues calculated.In this analysis, the reference to the elevation angle θ was dropped for notational convenience, but this is an important variation and must be accounted for.Since the required quantities are temperature differences (viz.o i,j ), the vertical variation is removed by processing the differences.We also need to estimate the quantities p i,j , p i and p j .These quantities are obtained by processing each image to remove the vertical variation of brightness temperature along each image column.A linear least-squares fit is obtained for each image column using data several lines above the plume up to several lines below the top of the image.The plume is discernible in the image data because it has a different temperature to the background sky and the camera viewing orientation can be arranged to completely view the plume, while allowing some clear sky to be imaged.Since each image is 240 lines high, the fit typically uses between 100 and 150 lines.Variations in the number of lines used in the fit occur because the plume is sometimes elevated and because some images contain noisy data towards the top of the image.In general the fit is very good (see Fig. 7).
The linear fit2 removes zenithal variations and provides estimates of o i and o j .Since each column of the image is treated differently, account is taken of any azimuthal variations in the atmosphere.Once this procedure has been applied, the plume temperature is estimated from the 10 µm image (the most transparent channel) after applying a correction for water vapour based on MODTRAN-4 radiative transfer calculations.Figure 8 illustrates the size of the atmospheric correction for the 10 µm channel as a function of the slant range, for three different plume temperatures: a cold plume with T p = 270 K, a plume close to the background atmospheric temperature with T p = 280 K, and a warm plume with T p = 290 K.

Error analysis
The SO 2 retrieval scheme makes several simplifying assumptions that can lead to error in the final results.The scheme depends mainly on the temperature measurements and measurement differences, but also on a few parameters (e.g.absorption coefficients, viewing angles).The sources of error are considered to fall into three distinct groups: -type 1 errors due to measurement noise; -type 2 errors, arising from assumptions and approximations used in the retrieval scheme; -type 3 errors due to inaccurate or incomplete specification of parameters required in the scheme.

Type 1 errors
The theoretical formula for the noise-equivalent temperature difference (NE T ) that produces a signal-to-noise ratio (SNR) of unity for a single microbolometer pixel may be written as (Derniak and Boremann, 1996) where F # is the F number of the camera, f is the sampling frequency, I is the radiance, A d is the area of the detector, and D * is the normalised detectivity or figure of merit of the detector.There are several sources of noise for thermal imager detectors including, Johnson noise, 1/f noise, and noise due to temperature fluctuations.The last of these noise sources is usually the limiting factor.For the Cyclops camera, D * ∼ 2.5 × 10 6 cm Hz 1/2 mW −1 , A d = 45 µm, f = 60 Hz, and F # = 1.Inserting these values into Eq.( 36 . The NE T 's (in mK) for the camera were calculated for a given set of scene brightness temperatures using the derivative of the Planck function at the central wavelengths of the channels, and these are shown in Table 2.
It can be seen by comparing the values in Table 1 with the theoretical noise temperatures of Table 2 that the camera meets the requirements for scene temperatures down to 250 K but not down to 220 K.In practice we have found that the theoretical limits are not met unless some averaging is done.Frame averaging can reduce the noise by √ N f , where N f is the number of frames.However, there is a limit to this as the fixed-pattern noise (FPN) is not reduced by adding more frames.The FPN is reduced by the use of the blackened shutter.Laboratory and field experiments were conducted to establish performance metrics for the thermal imaging camera.These trials suggested that 24-frame images were considerably more noisy than the theoretical results suggest 3 .The measured NE T 's ranged from 0.1 K at 290 K for the broadband channel to 1.8 K at 220 K for the 8.6 µm channel.A least-squares polynomial (third degree) fit to the laboratory data was performed for each channel so that the NE T at any arbitrary scene temperature (T s ) could be obtained.The fit is given by The coefficients for all channels, including the broadband channel, are given in Table 3.At 260 K the NE T = 0.80 K for the 8.6 µm channel, and 0.41 K for the 12 µm channel.The trials also showed that ∼ 0.5 % of the pixels were "dead pixels" -that is, these pixels were constantly off and registering no signal.Once these pixels had been identified they were flagged and not included in any further analyses.
Temperature differences are used in the retrieval scheme.Thus errors due to noisy measurements are increased by NE T 2 i + NE T 2 j , where i and j are the channel numbers.The noise in the measurements represents a large source of uncertainty in the retrieval scheme.We evaluate this by performing a large number of simulations where we specify the temperatures in Eq. ( 22) and include a Gaussian distribution of noise with the mean given by the NE T 's for each channel with a spread of 2σ .A perfect measurement is determined as the result when the NE T 's are zero.The result of these simulations gives an impact of 9-10 % on the retrieved SCD.Reducing the NE T 's by a factor 2 reduces the error to 6-7 %.
3 An improved camera made by FLIR Inc. has a lower NE T .
Calibration data suggest that the absolute errors are between 0.5 and 2 K, depending on the scene temperature, the environmental temperature and the channel used.Since the retrieval scheme uses temperature differences, as long as the channels behave in a similar manner, the actual impact of absolute temperature error is not great.The main impact arises through the estimate of the plume temperature made using the 10 µm channel.The random error associated with the estimate of the plume temperature is given as type 2 error, and here we assume only the component of the calibration error that contributes to bias.The bias error is close to zero when the environmental, scene and camera housing temperatures are the same.Thus the bias error is likely to be variable and may change sign, depending on whether the scene is warmer or colder than the instrument.Temperature offset calibrations are carried out every 5-6 min using a blackened shutter, attempting to minimise the effects of environmental temperature changes.The source of error for these calibrations arises from the non-blackness of the calibration shutter.The performance of the shutter was measured by comparing it to a laboratory blackbody of emissivity ∼ 0.99.It was established that the shutter emissivity was ∼0.98 ± 0.005, with a slight wavelength dependence.An error of ±0.005 in emissivity results in a temperature error of < 0.1 K, which is much smaller than the NE T of the filtered camera channels.These considerations suggest that an absolute calibration accuracy of ± 0.5 K is reasonable.While this is a bias error, the sign of the bias is likely to be variable and difficult to establish unless measurements of the environmental, camera housing and scene temperatures are available.This calibration error translates to an error in the SCD of ∼±5 %.

Type 2 errors
These errors are due to assumptions made in the derivation of the retrieval scheme.These assumptions include 1. plane-parallel radiative transfer, no scattering, radiative transfer (RT) model; 2. linearisation of the radiances to brightness temperatures; 3. constant plume temperature; 4. no spatial variation of the SO 2 absorption coefficient; 5. equivalence of the water vapour absorption coefficients at 8.6 and 12 µm; 6. invariance of the atmospheric structure with or without the plume present.
Assumption 1 includes commonly made assumptions for solving radiative transfer problems in the infrared region.For geometries where the plume is small in comparison to the curvature of the Earth, the radiation paths are almost identical www.atmos-meas-tech.net/7/2807/2014/ to real paths.The effect of assuming that the plume is planeparallel is inconsequential since the retrieval determines the SCD, and hence the actual geometry of the plume is irrelevant.Other aspects of the radiative transfer include use of the MODTRAN-4 code, which has undergone detailed scrutiny and intercomparison (Berk et al., 1998).It is difficult to make a precise estimate of the likely impact of errors in the radiative transfer modelling on the retrieval, but RT models suggest errors of 0.2-0.5 K are possible (Strow et al., 2003).We take ±2 % as an estimate for modelling errors (±0.2K in a 10 K temperature difference).Assumption 2 has been discussed in by McMillin and Crosby (1984), who show that a necessary condition for this approximation to be valid is that the radiances should be similar.This is easily examined by comparing the radiance calculated directly through the Planck function with the radiance calculated using a firstorder Taylor series approximation, e.g.Eq. ( 17). Figure 9 shows the radiance error (in %) dependence on the departure of the temperature from a mean value, and demonstrates that the error is less than 2.5 % in radiance for departures from the mean temperature of up to ±10 K.This radiance error results in an SCD error of less than half of that due to the measurement NE T or ∼ 5 %.
The impact of assuming that the plume temperature is constant (assumption 3) could be significant because the thermal contrast of the plume contributes significantly to the SO 2 signal through Eq. ( 22).In the early stages of generation, the plume is likely to be very inhomogeneous and in thermal disequilibrium.When the plume has been generated from a large explosive eruption, it may remain inhomogeneous for tens of minutes 4 .An idea of the plume temperature variation can be obtained from an analysis of the broadband (7-14 µm) channel data.These data are the least noisy and the variation can be used as a proxy for the variation in the thermodynamic temperature structure.The coefficient of variation for the stable plumes studied here is ∼ 0.01, and the typical temperature variability along the axis of the plume is ±3 K.If it is assumed that these metrics also apply to the thermodynamic temperature and that the magnitude of the variability does not change with position within the plume, then use of 4 We only consider eruptions where the Volcanic Explosivity Index, VEI,< 3 Eq.( 22) with the plume temperature perturbed by ±3 K gives SCD retrieval errors of 12-14 %.
Information about the spatial variation of the SO 2 absorption coefficient (assumption 4) is not available.There is a small pressure and temperature dependence of the absorption coefficient, but given that the range of variability of pressure and temperature is small for the observing conditions, this dependence may be neglected.
Assumption 5 has been examined by use of the water vapour transmission model of Davis and Viezee (1964).The model asserts that the water vapour transmission (τ λ ) within the window region 8-12 µm is governed by where λ is wavelength, w is the precipitable water amount (in cm), P * is the effective pressure, P * = P /P s , P = pressure (mb), P s is the surface pressure, k λ are the absorption coefficients, and a λ are coefficients determined by comparing the model with experimental measurements.The coefficients k λ and a λ are tabulated at 25 cm −1 intervals from 800 to 1200 cm −1 .The model was used to compute the transmission over the 8.6 µm and 12 µm filter response functions5 as a function of water vapour amount, up to 5.5 cm of precipitable water.A measure of the difference between absorption at 8.6 and 12 µm is computed as Err = (τ 8.6 − τ 12 )/τ 8.6 × 100 %.Largest error (Err) is found for greatest precipitable water amounts and reaches about 10 % at 5 cm.We put an upper bound on the error due to assumption 5 as 10 %, and the impact of this error on the retrieved SCD is at most a 3 % positive bias; that is, higher SCDs are recovered under this assumption.
The assumption that the atmosphere is the same whether or not the plume is present seems intuitively reasonable as the atmospheric path under consideration is much larger than the path within the plume.Also, the atmospheric radiance is calculated by a linear or quadratic interpolation of the atmospheric radiance above and below the plume.The very linear nature of the fit obtained demonstrates that this is a good approximation.Nevertheless, there is error involved.This is estimated from the 1-σ uncertainty estimate obtained from the least-squares fit.The uncertainty is evaluated for the 8.6 and 12 µm channels and for the difference.The 1-σ uncertainty for the difference was ±0.3 K, which translates to an SCD error of ±3 %.

Type 3 errors
Several of the parameters used in the retrieval scheme need to be specified.These include k SO 2 , geometry (elevation of the camera and FoV size of the camera), channel filter response functions, and the use of radiosonde data in the RT model.The absorption coefficient was obtained by integration over the filter response function using NIST values of the absorption coefficient measured at 0.125 cm −1 resolution.The likely error incurred is small compared to other errors.An error in the absorption coefficient translates directly into an error in the retrieved SCD.We take this error as 1 %.
Errors in the geometry arise from incorrect specification of the FoV of the instrument, and inaccuracies in measuring the camera elevation.These errors are all small and affect the retrieval only through cos θ and via the RT calculation, which uses radiosonde data and requires specification of the geometry of the calculation.The geometry error is less than ±0.5 %, which corresponds to an error in measuring the angles of ±1 • .
Errors in the radiosonde data (temperature and water vapour profile errors) affect the retrieval through inaccurate calculation of the 10 µm plume temperature.This error has already been incorporated as a type 2 error for the estimation of the plume temperature.The errors arising from all sources of error considered are summarised in Table 4.The final error is the root-meansquared sum of all of the individual random errors -that is, excluding the absolute calibration and transmission approximation errors.Thus the error on the retrieval is estimated to be ∼20 % with a bias of −5 to +6 %.

Field trials, detection and quantification
The retrieval scheme described above is quite complex, and so here we analyse some of the thermal imagery to illustrate the main parts of the scheme.Experiments were conducted at Mt Etna, Sicily (37.755 • N, 14.995 • E; 3330 m a.s.l.), and at Stromboli (38.789 • N, 15.213 • E; 920 m a.s.l.), Aeolian Islands, north of Sicily.Figure 10a shows the temperature difference ( 12,11 ) image between the 12 and 11 µm channels for data acquired at Mt Etna, with a 12,11 height profile shown for a single image column, indicated by the continuous vertical line drawn over the image (profile in Fig. 10a).Above the terrain and vegetation, there is a noticeable decrease in 12,11 which coincides with the plume from Etna.This decrease in 12,11 is likely caused by water vapour in the plume.By contrast, Fig. 10b shows the temperature difference 12,8.6 , which is negative everywhere, and there is also a noticeable anomaly in the vicinity of the Etna plume.This anomaly is due to the presence of both water vapour and SO 2 .Since the absorption by water vapour is slightly greater at 12 µm than at 8.6 µm, if no water vapour were present in the plume, then 12,8.6 would be less negative.There is also water vapour present along the path from the camera lens to the leading edge of the plume, and hence in regions of the atmosphere away from the plume, 12,8.6 is still negative.If the atmosphere were completely absent of water vapour, then 12,8.6 would depend on the temperature profile and the absorption by the uniformly mixed gases, of which CO 2 is the most important in this waveband.The 12,8.6 profile in Fig. 10b also exhibits a marked decrease with height in the atmosphere.Since the SO 2 signal that we wish to recover is masked by these other features due to water vapour and its height variation, it is necessary to try to remove them so that a background value (SO 2 -free atmosphere) for 12,8.6 can be found.This is the purpose of the fitting procedure described earlier.
In the discussion so far we have not looked at the influence of meteorological clouds.Correcting for their effects and attempting to retrieve SO 2 in the presence of clouds is extremely difficult using thermal data as it is necessary to know the microphysics (particle size, shapes and size distributions) as well as the thermodynamic phase of the clouds.The approach taken here is to try to detect clouds and other interfering substances (e.g.volcanic ash) and flag these image pixels as erroneous.Figure 11a and b show an example of cloud detection in Cyclops imagery.Figure 11a shows the 12,8.6 as before, with a single 12,8.6 -height profile taken through what appears to be a small meteorological cloud.The profile shows that the anomaly due to this feature is less negative than the rest of the profile and the difference approaches 0 K.When the corrections for the vertical variation of water vapour are taken into account, this feature appears as a positive anomaly and would be retrieved as a  negative SCD and hence is flagged as erroneous.In Fig. 11b we illustrate how an SO 2 anomaly in the same image appears to cause an opposite effect to that of meteorological water clouds.Ash can also interfere with the retrieval scheme, and ash clouds are often encountered with SO 2 gas emissions.Figure 12 illustrates the effect of an ash plume eruption on the 12-8.6 µm temperature difference.The ash plume eruption was identified in consecutive image frames (different spectral channels) separated by ∼ 1 s that captured the rapid evolution of the cloud when compared to an SO 2 gas emission.The ash cloud is also clearly discerned against the background atmosphere and the SO 2 gas, through its positive temperature difference anomaly.As with meteorological clouds, an ash cloud anomaly is easily identified and removed from the analyses (see also Fig. 17a).Having established that SO 2 can be identified and discriminated from other features, we now turn to the quantification of SO 2 retrieval and begin with a simple case where SO 2 is the only emission.

Port Pirie, South Australia
In order to test the ability of the camera to measure SO 2 , it was taken to a smelter and pointed towards a tall stack known to be emitting an SO 2 plume.The Port Pirie lead smelter, in South Australia (33.18 • S, 138.02 • E), is the single largest lead smelter in Australia, with mean SO 2 emissions of 1 kg s −1 (80-130 t d −1 ; see http://www.epa.sa.gov.au/).The plume is invisible to the eye (low water content) and emanates from a ∼ 200 m tall stack.The camera was placed ∼ 570 m from the stack and viewed it from the ground, looking upwards at an elevation angle of 15 • with a clear blue sky background.Measurements were made continuously, which provided SO 2 estimates at intervals of 4-6 min.The length of time between samples is determined principally by the speed of data transfer and to a lesser degree by the need for capturing images at several wavelengths (different filters) and for acquiring calibration data.A typical sequence consisted of five measurements of the blackbody shutter (one measurement for each filter), followed by five measurements of the of the plume -it is not possible to tell whether the bifurcating plumes are coplanar or whether one plume is moving away from or towards the camera; (e) reduced activity with little or no gas plume visible; and (f) vigorous gas pulse with an indication of a column and some horizontal dispersion.Note that the apparent high SO 2 SCDs just below 5.52 km height and near ∼−4.9 km are probably artefacts due to difficulties in calibrating the images on one part of the focal plane array.scene (the SO 2 plume), followed by a further five measurements of the blackbody shutter.Radiosonde profiles from Adelaide International Airport (about 30 km distant) were acquired for use in calculating the water vapour corrections; however the corrections were small and below the noise limit of the camera and were not applied in the retrieval.The SO 2 signal was very large and clear in the data; however it was necessary to use a second-degree polynomial fit to the brightness-temperature-height profiles (Fig. 13).The final fits and retrieval were robust.Figure 14a-c show a sequence of SO 2 retrievals illustrating the behaviour of the gas plume.At the start of the sequence (Fig. 14a), the plume rose ∼ 50 m above the stack and then became bent over in the light winds.Later, the plume fumigated (Fig. 14b), and eventually, with a change in wind speed and direction, the plume became stronger and was carried away from the viewing site (Fig. 14c).It should be noted that, with one camera, it is not possible to discern the direction of travel of these gas plumes in the plane aligned with the camera viewing direction.For quantitative studies of gas plumes it would be preferable to use three cameras spaced at 120 • to each other.The mean SCD for the Port Pirie plume on this day was ∼ 3 × 10 19 molecules cm −2 , with instantaneous maximum SCD near the stack exit exceeding 10 20 molecules cm −2 .It is possible to estimate the average SO 2 emission rate from these data using estimates of the wind speed at stack height and the effective plume dimensions.An estimate of the SO 2 emission rate can be found from where F is the emission rate (in kg s −1 ), ρ is the concentration (in kg m −3 ), A is the cross-sectional area of the plume (m 2 ), and u is the wind speed (in m s −1 ) at plume height.
Wind speeds at 200 m were ∼ 3-5 m s −1 and the plume width (measured at half maximum) was taken as ∼ 20 m (see Fig. 14).These values give emission rates of ∼ 1.5-2.5 kg s −1 , slightly higher than the mean emissions reported.In principle it is also possible to estimate the plume speed by tracking features in the plume (e.g.Bluth et al., 2007); however in the current configuration of the camera, the data capture and calibration cycles require ∼ 5 min to complete and thus feature tracking is difficult.The success of this field trial at a site where SO 2 could be independently identified and measured gave us confidence to test Cyclops at active volcanoes.

Etna volcano, Italy
In September 2003 the camera was taken to Mt Etna on the island of Sicily to conduct SO 2 measurements under field conditions.Measurements were made from several locations, in most cases more than 10 km from the active vent.The retrieval of SO 2 from Etna is illustrated in Fig. 15.At one site, the camera was mounted on a rooftop in the village of Nicolosi, approximately 17 km from Etna, and viewed the plume almost due north (350 • azimuth) at an elevation angle of about 20 • .At this low angle and distance, the water vapour path was significant and we regard this viewing configuration as being at the limit of the camera's capability.The data were acquired at 4-6 min intervals throughout the evening and into the following morning with no operator intervention and utilising automatic shutter calibration.
The raw images were converted to brightness temperatures using pre-computed laboratory calibrations and adjusted using the off-set shutter calibration procedure.During the sequence of measurements, the plume was blown in a NW direction and was confined to the boundary layer, remaining below ∼ 5 km (a.s.l.) most of the time.In the morning, with the break-up of the nocturnal inversion layer, the plume was observed to rise (Fig. 15b).Some variability in the SO 2 gas emission rate was observed over the period with quiescent periods (Fig. 15e), strong puffing activity (Fig. 15f) and plume bifurcation (Fig. 15c, d).
Emission rates can be determined, as before, from Eq. (39).Values for A and u are not known accurately, but assuming the plume to be symmetric, the data suggest an average plume depth of ∼ 500 m.The mean plume speed was estimated by running a trajectory model -HYSPLIT (Draxler and Rolph, 2003) -starting from the summit elevation at 23:00 LT on 22 September 2003 and run forwards for 8 h.The trajectory of the plume found this way was towards the NW with a mean wind speed (over 8 h) of ∼ 2 m s −1 .Using these values, we find F =∼ 10-20 kg s −1 , and the variation with time over 7 h of continuous measurements is shown in Fig. 16.There are many (unsystematic) measurements of Etna SO 2 emission rates reported in the literature based on different measurement techniques (e.g.Jaeschke et al., 1982;Teggi et al., 1999;Barrancos et al., 2008;Oppenheimer et al., 2006;Bobrowski et al., 2006).These report in situ, remotely sensed UV and IR, ground, aircraft and satellite platform-based retrievals from different years and different months.The variability is high, depending on the degassing phase of activity with emission rates varying from 11 kg s −1 (Oppenheimer et al., 2006), to 82.2 kg s −1 (Teggi et al., 1999).A proper, statistical evaluation and intercomparison of the IR camera retrievals is beyond the scope of this paper, but new work resulting from a volcanic plume workshop, where several UV cameras and the IR camera are compared, has been submitted for publication (Kern et al., 2014;Prata et al., 2014;Lopez et al., 2014).

Stromboli volcano Italy
Measurements at Stromboli were made on two separate occasions in late September 2003.Stromboli is an active stratovolcano which has been erupting and degassing SO 2 throughout historical time.The effusive activity is observed from four vents near the summit and usually consists of small explosions followed by a period of quiescence which lasts from 10 min to a few hours (Andronico et al., 2008).Very little ash was observed during the activity in September 2003.The Cyclops camera was used from two locations: near sea level from the rooftop of a hotel (site A) and ∼ 2.3 km ENE from the active crater, and nearer the volcano at Semaforo Labronzo (site B), 120 m a.s.l., and ∼ 1.7 km north of the crater.At both locations the camera elevation was high (>25 • ).
Long sequences of images were captured at both sites.The SO 2 plume was often mixed with water vapour (judged by its white appearance) and tended to erupt in puffs and disperse in the light winds (< 5 m s −1 ).The retrievals indicate that total SO 2 SCDs varied from 1.1 ± 0.2 × 10 18 to 3.1 ± 0.6 × 10 18 molecules cm −2 at site B and 1.4 ± 0.3 × 10 18 to 2.3 ± 0.5 × 10 18 molecules cm −2 at site A. Individual plumes had variations from ∼ 2 × 10 17 to ∼ 3 × 10 18 molecules cm −2 , but these values are difficult to interpret in terms of concentrations because the plume depth variable and unknown.examples of the retrievals at both sites are shown in Fig. 17.Measurements from site A were made in the late afternoon, and the gas emissions appeared to be continuous.At site B, measurements were made in the evening after the Sun had set, and the gas emissions occurred frequently as discrete puffs.Explosions were also heard and several were imaged by the camera.Since the presence of ash can confound the retrieval scheme, it is important that the algorithm be insensitive to ash or be able to flag regions of the sky contaminated by ash. Figure 17a shows an occasion when an explosion occurred during the imaging.In this case the algorithm has rejected pixels that are ash-contaminated, and this is indicated on the image by the grey-coloured region.Within this region no SO 2 can be retrieved, and consequently the total SCD for the whole plume will be underestimated.Kern et al. (2014) report measurements from a several UV camera systems during an inter-comparison experiment on Stromboli in June 2013.The UV retrievals (SCDs) ranged from 2 to 4 × 10 18 molecules cm −2 .

Conclusions
An instrument for measuring atmospheric gases using passive thermal infrared imaging radiometry has been successfully tested, and a scheme for retrieving the SCD utilising multispectral imagery has been derived.The camera system was tested at an industrial site and at two volcanic sites where plumes of SO 2 were present.The instrument proved reliable and was able to detect SO 2 in the presence of water vapour.Cloud features and ash particles that interfere with the SO 2 measurements could be detected using the multispectral nature of the imagery and removed from the analyses.The retrieval scheme proposed relies on spectral temperature difference and is sensitive to channel NE T 's estimated to be 0.1-0.8K after frame averaging.The error on the retrieval was estimated to be ±20 %, due mostly to the NE T error, but with a significant error to due to inaccurate measurement of the plume temperature.The bias error is estimated to be variable within the range −5 to +6 %, which is slightly poorer than the calibration errors reported by Shaw et al. (2005) for their broadband camera.
The infrared camera is capable of providing 320 × 240 pixel images at frame rates as high as 60 Hz.However, it was found that noise considerations and data capture rates reduced this sampling frequency to several minutes.The limitation in sampling frequency is dominated by the slow data transfer rates, which can be easily overcome6 .Faster sampling would allow for measurements of the dynamic evolution of plumes, and feature tracking could then be used as a means to determine gas emission rates.A more fundamental limitation is the NE T of the spectrally filtered channels.There are cameras available commercially with NE T 's of < 20 mK and 60 Hz frame rates that can provide retrieval errors in SCDs below 10 %.Many improvements to the system can be envisaged.By viewing a target using three cameras arranged with an angular spacing of 120 • , a three-dimensional image could be acquired and quantitative measures of plume dimensions and plume morphology derived.Addition of filters centred at different wavelengths would also permit a range of other gases to be measured.The camera could also be used in atmospheric research for studies of the radiative effects of clouds on the Earth's radiation balance (Smith and Toumi, 2008) and to image toxic gases from industrial accidents or from deliberate gas releases, where personal safety is a major issue.
The system described here has been operated from the ground, but it is quite feasible to use the system from an airborne platform.In this case, operation from higher altitude would permit the use of spectral filters at wavelengths where water vapour is a problem in ground-based use.A filter situated near the 7.3 µm band would have 3 to 5 times the sensitivity to SO 2 as the 8.6 µm filter used here.One application for this technology in airborne use would be to mount the instrument to view forwards from a high-altitude passenger jet aircraft.In this case it would be necessary to remove the filter wheel and use multiple cameras in order to achieve faster sampling rates.The cameras would offer the potential as an on-board early warning device for hazards ahead of the aircraft (Prata and Barton, 1993).Hazards include volcanic ash and potentially small (∼ 1-20 µm particle radii) ice crystals and clear-air turbulence that may be detected through imaging water vapour anomalies.Enhanced night-time viewing capability is another feature of this technology that might be useful for jet aircraft.
Integration of the camera with other instruments is feasible.For example, infrasound arrays, ground-based lidars, ultraviolet cameras and spectrometers, and FT-IRs all offer complementary information which would enhance the ability of a system to detect a suite of gases and measure their concentrations and emission rates (e.g.Lopez et al., 2013).Further improvements to the system have been made, including integration of a webcam, low-light imager, Wi-Fi, and weather proofing.These are described in Prata et al. (2014).Stand-off, 24 h, autonomous operation of the Cyclops camera has been demonstrated at two active volcanoes, and plans are in place to deploy the system for long periods to test the durability of the instrument and the reliability of the detector calibration methodology employed.Absorption coefficient for channel i and absorber q L Path distance from camera to leading side of plume m q Slant column density (SCD) (= ρ q d) for absorber q n Pixel number

Figure 1 .
Figure 1.Schematic showing the main components of the "Cyclops" thermal imaging infrared camera.Note that the filter wheel, containing up to 5 filters, is placed behind the lens.

Figure 2 .
Figure 2. (a) Filters mounted on filter wheel in the arrangement when used for measuring SO 2 gas emissions (central wavelengths in microns are given).(b) "Cyclops" camera mounted on a tripod for field operation.(c) Ship-mounted camera undergoing calibration tests with two moveable blackbodies and an external blackened shutter.

Figure 4 .
Figure 4. Cyclops spectral images obtained at an SO 2 -free, particulate-free site in Australia.(a) Panels showing uncalibrated data (DNs or counts), (b) their respective histograms, (c) panels showing calibrated images, (d) their histograms, and (e) histograms of selected temperature differences.The order of the images starting from the top is 8.6, 10, 11, 12 µm and the broadband (7-14 µm) channel.

FIGURE 5 .Figure 5 .
FIGURE 5. Filter response function (smooth line) for the 8.6 µm Cyclops channel and the variation of the SO2 absorption coefficient with wavenumber as measured by NIST.The integrated absorption coefficient over the waveband is 4.3235 x 10 −5 µmol mol −1 m −1 .
i , T p and θ n , and specification of the absorption coefficient k j,2 .The measurements consist of the plume radiances (I p i , I p j ), the foreground radiances (I f i , I f j ) and the background radiances (I b i , I b j )

Figure 7 .
Figure 7. (a) Brightness temperature versus height variation for the 12 µm filter ( 12 ).(b) Brightness temperature versus height variation for the 8.6 µm filter ( 8.6 ).(c) Brightness temperature difference versus height variation for the 8.6-12 µm fiters ( 8.6,12 ).The straight lines are least-squares linear fits based on profile data above the plume, and extrapolated through and below the plume.

FIGURE 8 . 40 Figure 8 .
FIGURE 8. Atmospheric correction, Tp-T B10 (in K) as a function of the slant range to a plume at three different (uniform) temperatures (Tp).Calculations were performed using MODTRAN-4 for the 10 µm channel.

FIGURE 9 .Figure 9 .
FIGURE 9. Radiance error (in %) versus departure from the mean temperature (K) caused by approximating the radiances using a 1st order Taylor series expansion about a mean temperature.

Figure 10 .
Figure 10.(a) 12-11 µm brightness temperature difference ( 12,11 ) image of the Mt Etna plume.The panel to the right shows a temperature-difference-height profile for one image column, indicated by the vertical line on the image.(b) As for (a) but for the temperature difference between the 12 and 8.6 µm channels ( 12,8.6 ).The height profile for the same column as (a) is shown to the right of this image.

Figure 11 .
Figure 11.(a) 12-8.6 µm brightness temperature difference ( 12,8.6 ) image of the Mt Etna plume.The panel to the right shows a temperature-difference-height profile for an image column that intersects a small meteorological cloud.(b) As for (a) but the height profile now intersects a portion of the Etna SO 2 plume.

Figure 12 .
Figure 12. 12-8.6 µm brightness temperature difference image of the Stromboli plume acquired during a small explosive eruption.The panel to the right shows a temperature-difference-height profile for an image column that intersects the ash cloud eruption.

Figure 13 .
Figure 13.(a) Brightness temperature versus height variation for the 12 µm filter ( 12 ) (b) As for (a) but for the 8.6 µm filter ( 8.6 ).(c) The 8.6-12 µm brightness temperature difference ( 8.6,12).The curved lines are least-squares second-degree polynomial fits based on profile data above the plume, and extrapolated through and below the plume.

Figure 14 .Figure 15 .
Figure 14.Cyclops measurements of the SO 2 plume from the industrial stack at the Port Pirie lead smelter, showing different behaviours of the plume.(a) Lofted plume, (b) fumigation and (c) grounding.

FIGURE 16 .Figure 16 .
FIGURE 16.Variation of the mean SCD (molecules cm −2 ) as a function of time for the Etna plume over a period of 7 hours, starting from midnight until 07:00 LT the following day.

FIGURE 17 .Figure 17 .
FIGURE 17. Stromboli SO 2 observed from two different sites, ∼1.7 km and ∼2.3 km from the active crater.(a) shows an image where there was a small ash eruption obscuring the SO 2 .The 'ashy' parts of the cloud can be masked out by using different channels.Local Time (LT)=UTC+2.

Table 1 .
Channel number, central wavelength, bandwidth, purpose and required noise-equivalent temperature difference (NE T ) for Cyclops.

Table 2 .
Theoretical NE T 's (mK) for four narrow-band channels of the thermal infrared imaging camera and for four different scene temperatures.

Table 3 .
Polynomial fit coefficients for computing NE T as a function of scene temperature and channel.

Table 4 .
Summary of error types and estimated error magnitudes.
l Number of lines in the image N c Number of columns in the image r Radiation path in the direction θ, φ r 1Pathlength of plume radiation in the direction θ n , φ n s n