The detection of carbon dioxide leaks using quasi-tomographic laser absorption spectroscopy measurements in variable wind

Laser absorption spectroscopy (LAS) has been used over the last several decades for the measurement of trace gasses in the atmosphere. For over a decade, LAS measurements from multiple sources and tens of retroreflectors have been combined with sparse-sample tomography methods to estimate the 2-D distribution of trace gas concentrations and underlying fluxes from point-like sources. In this work, we consider the ability of such a system to detect and estimate the position and rate of a single point leak which may arise as a failure mode for carbon dioxide storage. The leak is assumed to be at a constant rate giving rise to a plume with a concentration and distribution that depend on the wind velocity. We demonstrate the ability of our approach to detect a leak using numerical simulation and also present a preliminary measurement.


Introduction
Carbon capture and geological storage (Bachu, 2008;Leung et al., 2014) is one candidate to reduce the concentration of carbon dioxide in the atmosphere.It is well known that the total annual increase in CO 2 is in excess of 1 part per million in volume (Keeling et al., 1995) which is 7.8 Pg (7.8 gigatonnes) of carbon dioxide.Geological sequestration at 30 kg m −3 of CO 2 is possible (Silva et al., 2015).To have a significant effect on the course of climate change, on the order of 1 part per million in volume of CO 2 must be removed annually; hence, an enormous volume scaled to 260 km 3 of storage would need to be added every year for at least the next few decades.
One obvious objection to geological repositories is that the CO 2 could leak.How much leakage is acceptable?Because of the various geological processes for removing CO 2 from the atmosphere that take place over a timescale of hundreds to thousands of years, it is not possible to specify a single time constant (Archer et al., 2009).Still, leakage below 0.1 % per year is required so that even slow geological processes can act to avoid most of the environmental effects of CO 2 in the atmosphere (Hepple and Benson, 2005).An inexpensive method of monitoring storage sites is required.
Remote sensing of atmospheric gasses, frequently for pollution control goes back several decades (Sabins Jr., 1997).One of the principal methods is differential optical absorption spectroscopy (DOAS).As implemented by Pundt et al. (2005), two optical sources send light across a field in various directions to allow for spatially resolved measurements, making what we call quasi-tomographic measurements.A similar technique to understand a vertically oriented slice of the atmosphere has been studied in the context of the plumes of volcanoes (Johansson et al., 2009).
Long-baseline DOAS has been used to measure trace gasses in the atmosphere tomographically.Whereas in medical tomography, the number of individual projections can be in the millions, the practice in tomographic DOAS has been to make tens of measurements in two dimensions (Pundt et al., 2005).For order of magnitude estimation, the spa-tial resolution is proportional to both the square root of the number of measurements and the area of the sample region (Natterer, 1986).In particular, using the filtered backprojection algorithm (Kak and Slaney, 2001) with a first generation tomographic set up in two dimensions (i.e., parallel rays entering a circular region), the spatial resolution is a factor of √ π/2 worse than one would obtain by taking the square root of the number of measurements and applying the Nyquist sampling criterion.
In the case of the verification of the integrity of structures dedicated to the sequestration of carbon dioxide, we may suspect a single point leak in the presence of a reasonably steady wind.The sequestration region might be on the order of 1 km square, so using roughly 36 measurements as considered by Hartl et al. (2006) leads to an estimated spatial resolution just below 200 m.Measurements are typically made a few meters above the ground.Hartl et al. (2006) concluded that measurements with more than two light sources were highly advantageous.Olaguer (2011) argued that the use of two light sources could be sufficient if plume models were used in a reconstruction.Modeling studies include those of Verkruysse and Todd (2005) and Chang and Wu (2012).
Typically, it is assumed that the measurement is taken quickly enough so that the plume does not vary due to changing winds.Here, we consider a different constraint: if we assume that the leak rate Q is constant during the course of the measurement, and that the wind velocity is measured during the measurement, we may seek a plume which shifts direction with the wind and whose density is inversely proportional to the wind speed.Recently, other authors have conducted similar studies (Humphries et al., 2012;Luhar et al., 2014).Other proposals for finding small leaks include consideration of the N 2 , O 2 , and CO 2 concentrations including isotopic ratios, as carbon dioxide sources are not associated with oxygen depletion.

Ermak's plume model
If we assume that the concentration of gas is due to a single steady leak, the strong constraints of a plume model can be imposed.If the speed and direction of the wind are known throughout the measurement then a given leak strength will give rise to a characteristic plume.We adopt coordinates in which x p is downwind, z p is up, and y p is chosen to form a right-handed coordinate system.We also make use of the frame of the Earth, with x being east, y being north, and z = z p being up.The coordinates x, y, and z are fixed, but x p and y p depend on the wind direction.
If the emission is continuous at a constant rate, making the standard approximation of the advection-diffusion equation, i.e., ignoring diffusion parallel to the wind direction (Stockie, 2011), we may assume that the plume amplitude is proportional to the plume basis function specialized to ground-level emission (Ermak, 1977): where is the Heaviside step function.The normalization in z p is over the half-space z p ≥ 0; z p = 0 is the level of the ground which is taken to be flat.We represent the plume as a basis function times a constant c which depends on the leak rate Q and the wind speed.The coordinate origin is at the point of the leak.The width parameter is given by where γ is a positive constant and σ 0 , v 0 , and x 0 are redundant constants which we also refer to, more compactly, with the variables In the plume model we adopt, a single width governs the dispersion in both y p and z p .However, some models suggest using two different functions for y p and z p (Turner, 1994).Three constants are introduced to emphasize scaling from a reference condition: the width is σ 2 0 at the downwind position x p = x 0 with a wind speed of v = v 0 .The constants k 0 and γ are positive and depend on the weather conditions.Although ideal molecular diffusion implies γ = 1, in practice γ takes on a value near 1, depending on the atmospheric stability (Ermak, 1977), with 0.9 being a typical value.The scaling with v occurs in Ermak's model because the transverse diffusion depends on time but is independent of the wind speed v.
The plume basis function is normalized such that reflecting the conservation of mass as the plume undergoes lateral and upward diffusion.The dimensions of ψ are thus inverse area, which may also be seen from Eq. (1).The actual concentration is given by c ψ.The physical meaning of c is the following: there is a certain leak rate Q which may be expressed in the SI units of kg s −1 .Because longitudinal diffusion is neglected in Ermak's plume model, the gas crossing each plane with constant x p must have an area-integrated concentration of Q/v so that the gas is "swept out" at the equilibrium rate.Given that the physical plume function is c ψ(x p , y p , z p ), Eq. (3) leads to the identification with units of kg m −1 .(The dimensions of c ψ are kg m −3 ; it is a concentration.)In our actual case, we will have several values of v, each representing a plume formed from the same leak.We will estimate a single Q and derive values for c for each wind condition based on the known v using Eq. ( 4).In practice, CO 2 concentrations are usually quoted as µmol mol −1 of the dry atmosphere, with 400 µmol mol −1 of CO 2 corresponding to 0.786 g m −3 of CO 2 .(The unit µmol mol −1 is the same as the unit ppmv; here µmol refers to the CO 2 concentration and mol refers to the concentration of dry air.)At a given distance downwind, the total concentration in a line of observation orthogonal to both the wind direction and the direction of gravity is given by Examining Eq. ( 5), there is a singularity as x p → 0 if z p = 0 (physically, the point of emission is at ground level), but lim x p→0 F (x p , z p ) = 0 for z p > 0. Downwind, for small values of x p there is a very rapid increase for x p (z 2 p /k 0 ) 1/γ followed by a gentle decrease for x p (z 2 p /k 0 ) 1/γ .There is a peak at x p = (2 z 2 p /k 0 ) 1/γ .The physical implication is that the largest concentration of observable gas will be located slightly downwind of the gas source.Just above the leak, the gas has not yet entered the line of observation; far downwind from the leak, the gas has dissipated upward too much.
Another property of interest is the peak observed value of the concentration for a given value of z p .By inspection of Eq. ( 1), this maximum occurs for y p = 0.The problem reduces to a 1-D maximization of Eq. (1) specialized by Eq. ( 2).The position of the maximum is while its value is given by The standard deviation at the maximum concentration obeys the simple relation In Ermak's formulation, the transverse diffusion is independent of the wind speed when considered as a function of time.The proportionality constant between time and downwind position is simply the wind speed.If downwind distances which are comparable to the length scale for pressure and density differences in the atmosphere are considered, this approximation will break down.However, for distances up to 500 m, it is a reasonable approximation (Turner, 1994).

Maximum likelihood formulation
The famous filtered back-projection method of tomography requires regular sampling.In Bayesian tomography (Sauer and Bouman, 1993) there is no such restriction.There, we calculate the posterior probability distribution of the model parameters, which indicate likely parameter values for the given data set by regions of high probability.Because we are using a formulation from tomography with far fewer measurements, we characterize the method as "quasitomographic".We take the prior distribution on the reconstructed parameters to be flat; thus, the mode of the posterior distribution corresponds to the maximum likelihood estimate.
Assuming the measurements, denoted n J , are independent, the likelihood is where n = (n 1 , . . ., n N ) are the observed counts whose expected values are (µ 1 , . . ., µ N ).Since we assume a flat prior distribution, Eq. ( 9) is proportional to the posterior distribution.Thus, normalizing Eq. ( 9) so that it integrates to unity (over the reconstructed parameters) yields the posterior distribution.In what follows, we refer to the unnormalized posterior and the likelihood interchangeably.Assuming each source follows Poisson statistics, the measurements, i.e., the counts of light intensity after undergoing some mean attenuation, also follow Poisson statistics (Haus, 2000;Hu et al., 2007).In this case, the likelihood is given by a product of independent Poisson distributions: Below, each µ J is determined by the set of parameters f .We next introduce the log-likelihood L = ln L given by Later in the paper, we will take f = (Q, x, y), three variable plume parameters described below.The function µ J (f ) relates the plume parameters to the observations of the column density of CO 2 .We take the distance-transmission relation to be Beer's Law: where µ (0) J is the mean of a Poisson distribution giving the source strength for observation J , and P J is the projection through the concentration of gas, i.e., the integral of the concentration times the absorption cross section (Robinson et al., 2014) The constant terms n J ln µ (0) J − ln n J ! do not contribute to the derivatives ∂L/∂f i and will be normalized away when a probability density is formed.
The projections are related to the concentration c ψ(r; f ) by The constant c is the attenuation length for a given concentration of gas and r is a position.The units of c are thus m 2 kg −1 , i.e., the inverse of a concentration times a length.The path r(s) is given by the domain of integration J , a line integral, also known as a projection in tomography.The posterior distribution p(f ) for the parameters f given the data n is given by i.e., the likelihood normalized by its integral over all possible values of the parameters.Let the parameters which maximize p be denoted by f 0 , which is assumed to be unique.We take f 0 to be the estimator of f .In performing the maximization, the most time-consuming step is the calculation of the projections in Eq. ( 14).We adopt a strategy of maximizing over the leak rate Q with the location of the leak fixed and maximizing again over Q whenever the leak position is varied.Here, we present the singlesource case, although for several sources, several concentrations may be optimized over without recalculating the projections.Using Eq. (4), we may write We take Q ≡ f 1 (and also x ≡ f 2 and y ≡ f 3 ).This allows us to write the derivative as We use Eq. ( 17) in maximizing L by picking a given point (x, y) for the gas source and setting ∂L/∂Q = 0. We maximize over the leak rate Q for each candidate position (x, y).Maximizing over Q for each (x, y) allows us to make the most efficient use of the projections which are relatively expensive to compute.The log likelihood maximized over Q is itself then maximized over (x, y) by Mathematica's NMaximize with default parameters.Up to 25 draws from a random search drawn uniformly from a square slowly expanding about the candidate maximum point are made until a higher value of the log likelihood was found.Then, Mathematica's local optimizer FindMaximum was used and the maximum found by it was subject to additional sets of up to 25 draws seeking further improvement.After one invocation of NMaximize, 10 iterations of the further process were permitted.While we have formulated Eq. ( 11) in terms of counts, we may obtain the information in the form of a transmission factor T J and an estimate of the signal-to-noise ratio S J .These representations are connected by the simple relations We consider the log likelihood to be a continuous function of the real numbers n J , although strictly speaking the derivation requires each n J to be an integer.In practice, to obtain simulated measurements, we use the model parameters to set the forward model.Next, given the signal-to-noise ratio, we derive expected mean counts for each measurement using Eq. ( 19).A particular measurement is simulated by taking that mean value to be a parameter in a Poisson distribution and taking a single sample.

Experimental setup
Field measurements from the Greenhouse gas Laser Imaging Tomography Experiment or GreenLITE (Dobler et al., 2015) motivated our tomographic reconstruction approach.GreenLITE is a prototype system of hardware and software intended to measure, monitor, and estimate 2-D distributions of CO 2 concentration at carbon sequestration sites (or elsewhere).The hardware component consists of two laser-based differential absorption transceivers and tens of corner-cube retroreflectors used to collect CO 2 transmission measurements in a crisscross pattern such that lines of sight from transceivers to reflectors intersect as shown in Fig. 1.The software component includes algorithms to convert transmission measurements into CO 2 concentrations in µmol mol −1 along each line of sight.These concentration values can be used to validate our approach.
The GreenLITE system was deployed at a farm outside of Fort Wayne, IN, in February of 2015 in the configuration shown in Fig. 1.A full set of transmission measurements -for all transceiver-reflector pairs -were collected approximately every 5 min on 5 February from 13:00 to 9 m plume parameter (Turner, 1994) x 0 100 m plume parameter (Turner, 1994) v 0 5 m s −1 plume parameter (Turner, 1994) z 0 1 or 10 m plume parameter of the two light sources starts 27 line segments which end on each of the 27 reflectors.The pale blue ellipse represents the region in which carbon dioxide sources were simulated.The red dot represents the position of the leak found in the experiment as calculated by maximizing the likelihood using the Ermak plume model.21:00 EST (local).Transmission measurements were then converted to CO 2 concentration in µmol mol −1 and averaged, per transceiver-reflector pair, over 20 min time intervals to minimize the impact of short-term variations in wind velocity and variability in local source and sink contributions, such as automobile traffic, that are not representative of a persistent leak at a sequestration site.Wind was also measured and reported during this time at 2 min intervals.Accordingly, wind was averaged over 20 min periods to coincide with µmol mol −1 averaging and is shown in Fig. 2. The 8 h time interval and 20 min averaging results in 24 sets of concentration values and corresponding wind measurements.Table 1 shows that 1162 transmission measurements were used in our study, and Table 2 shows the break down of transmission measurements per 20 min time window.Some measurements are missing due to routine optimization of the northern transceiver scan pattern during the early portion of the 8 h interval, as well as low return signal on multiple measurements due to snow that had accumulated on some of the reflectors.
Because no source of CO 2 was placed in the field, we were not expecting to detect a source.As reported in Sect.4.3, if all the data were taken into account, in particular the shifting winds, there was an apparent detection.This motivated us to see whether such a detection was plausible in simulation.

Simulations
First, we simulate measurements assuming that the true leak position is exactly at the estimated position from the experimental data.Simulated measurements suggest that the signal can be recovered well at both the measurement height of z 0 = 1 m as well as z 0 = 10 m.Due to the complex structure of the log-likelihood function at z 0 = 1 m, we shift to the value of z 0 = 10 m for most of the following discussion.This allowed us to find the maxima in hundreds of cases without searching interactively or developing a maximization routine which takes account of the cell structure which appears with the 1 m observation-height reconstruction.From Eq. ( 8), σ (x A glance at Fig. 1 shows that 7 m is closer to the scale length of the irregular grid formed by the light paths than 0.7 m.However, the value of z 0 = 10 m is somewhat arbitrary and does not represent an attempt to optimize the height at which measurements are to be taken.
We simulated 389 cases of a gas leak with strength Q = 10 5 µmol mol −1 s −1 and 109 cases with Q = 5 × 10 4 µmol mol −1 s −1 , with the leak position given by sampling from a uniform distribution of the pale blue ellipse shown in Fig. 1.Four cases with Q = 0 were also simulated.
We were able to find the maximum of the log likelihood in 380 of 389 cases with the higher leak strength and 103 of 109 cases with the lower leak strength using the algorithm described in Sect.3. Cases for which the reported maximum was less than the log likelihood of the known true value were discarded.However, the known true value was not otherwise used in the analysis.
A typical example of the log likelihood and its maximum is given in Fig. 3.The function is smooth enough for au- tomatic optimization yet the alignment of the peak of the log likelihood with the light paths shows that the detailed pattern of measurement exerts a considerable influence over the function.The peak lies thousands of log-likelihood units above the log likelihood evaluated at 0 leak rate.A typical case of the log-likelihood function computed without a gas leak is shown in Fig. 4.Although there is a superficial resemblance to Fig. 3, the key point is that the full range of the log-likelihood function is reduced by 3 orders of magnitude -enough that the hypothesis of no leak is consistent with the data.
Our simulated measurements are created from the projection of the CO 2 density derived from the plume model and the attenuation constant from Table 1 and the signal-tonoise ratio also from Table 1 to form Poisson means, which are sampled.For each set of simulated measurements, we may form a credible region for the true parameter values.Since posterior distributions may be approximated by Gaussian distributions centered at the mode of the posterior distribution (Gelman et al., 2014), credible regions may be approximated by the set of parameters with a log-likelihood value within half of the appropriate χ 2 quantile of the maximum value of the log likelihood.In the simulation, we use the convention that the value of the log likelihood at the true parameters is 0, so a plot of the ordered maximum values of the log likelihood vs. the quantiles of a chi-square distribution with 3 degrees of freedom (corresponding to the three parameters) will ideally form a straight line.Figure 5 depicts such a plot.The 0 of the scale is the log likelihood at the true values of (Q, x, y).The ordered log-likelihood values are shown on the y axis.The x axis is obtained as follows: Since the 95th percentile of a 1 2 χ -square distribution with 3 degrees of freedom is about 3.9, Q = 0 is contained in the 95 % credible region for the parameters.Thus, the data are consistent with the hypothesis of no leak.1. CDF refers to the cumulative distribution function.The red line represents ordered random deviates that align perfectly with their corresponding χ 2 quantiles.If the points wander outside of the bounds given by the green lines, the assumption that the random deviates follow a χ 2 distribution is less tenable.the ith point is assigned the probability value from the quantile rule P i = (2i − 1)/(2n).The value of the inverse function [χ 2 3 ] −1 (P i ) is given.If the procedure for creating credible regions based on the χ 2 distribution leads to credible regions that envelope the true parameters with the frequency we expect, the points will form a straight line.Since the plot is Positions of 380 simulations with Q = 10 5 µmol mol −1 m 3 s −1 measured at a height of z 0 = 10 m, with the true positions sampled uniformly from the pale blue ellipse in Fig. 1 with the x, y positions given there.The reconstructed position is shown with a dot.Each line segment runs from the reconstructed position to the true position.(Many of these line segments are very short and appear as small black dots.)The blue dots represent the lowest one-sixth predicted leak rates, and the red dots represent the highest one-sixth predicted leak rates.The green dots represent the middle two-thirds.
subject to noise (because of the Poisson sampling), the 95 % simultaneous bounds help to judge deviations from the ideal.Since the points of that plot all fall within 95 % simultaneous bounds, we are able to make two conclusions: first, the approximate credible regions work well in the sense that they cover the true parameters the frequency of times we expect; second, the maximization performed by our algorithm is substantially complete.
The ability to predict the position of the leak is shown in Figs. 6 and 7.The localization, typically to within 1 or a few meters, is quite good.Moreover, although the distance of the predicted leak position to the true leak position has an observable tail, there are no outliers.No obvious pattern emerges from considering separately points with a low or high predicted leak rate.The means and standard deviations of these discrepancies are summarized in Table 3.There is Table 3. Results of calculations.Units of Q are µmol mol −1 m 3 s −1 .N calc is the number of cases calculated for each leak rate, and N solv is the number for which likelihood maximization was successfully performed.Q is the mean predicted leak rate and σ Q is the standard deviation.is the mean distance from the predicted to the true values in the simulation, and σ is the standard deviation.no detectable bias in the estimates of the leak rates.Not surprisingly, the spatial localization is better at the higher leak rate.

Experimental results
The data were reconstructed with the model given above.Before reconstruction, we subtracted a background of µmol mol −1 (766 mg m −3 of CO 2 ) at the beginning of the measurement, rising linearly to 397.5 µmol mol −1 (781 mg m −3 ) at the end.Because the measurements were taken from 13:00 to 21:00 EST (local time) on 5 February 2015, given the sun set at 18:02 EST it is reasonable that the diminution of photosynthesis would lead to an increase in the CO 2 background.For the experimental data, we found the maximum manually from the likelihood shown in Fig. 8.
We find, as shown in Figs. 8 and 9, there is a very large peak in the log likelihood.The peak is 9537 log-likelihood units above the value for a leak rate of 0. This allows us to estimate that the presence of a leak at that location is roughly e 9537 ≈ 10 4142 times more likely than the absence of a leak considering only statistical uncertainties.We emphasize that we did not verify the presence of this source in this preliminary study.
We find that with a measurement height of z 0 = 1 m, the log-likelihood distribution is strongly affected by the exact positions of the light paths again as illustrated in Figs. 8  and 9. Also, in Fig. 10 the likelihood is shown considering the results of the first complete set of 54 measurements taken with a single value for the wind.While a source is definitely detected, the localization region is on the order of 100 times larger than if all measurements are taken into account.Hence, we conclude that following the wind is an effective way to increase the sensitivity of the experiment.

Conclusions
The challenge of carbon sequestration required the development of a technology to ensure that sequestration sites are not leaking substantial amounts of gas into the atmosphere.The required measurements have certain undesirable characteristics, namely the detection of a small signal against a much larger background and the variability of that background.Laser absorption spectroscopy, combined with the Ermak plume model, can be used to observe relatively small isolated sources of CO 2 in the presence of variable wind.Whereas much of the discussion of quasi-tomographic LAS observations has assumed that the measurements could be made faster than the wind could shift, here we measure the wind velocity and use that as part of the model.It would be difficult if not impossible for an external gas source to emulate the shifts required by the changing wind which greatly suppresses the effect of these external sources.
We have made a preliminary measurement of a gas source, giving both its position and strength.While we did not independently verify the nature and strength of the source, we did show that the strength of the detection was consistent with what would be expected given the correctness of the model and the experimental signal-to-noise value.Using actual wind data and a physically realized experimental protocol, we have found through simulation that the ability to localize both the strength and the spatial position of the leak is quite good.

Data availability
Experimental observations and simulated observations used in this paper are available in the Supplement as text files which are described in the file README.txt.The estimated leak location from Fig. 8 is shown with a dot.The result is in marginal disagreement with these measurements, lying between the 95th and 99th percentile contours of given by the χ 2 2 distribution divided by 2.
The Supplement related to this article is available online at doi:10.5194/amt-9-1627-2016-supplement.

Figure 1 .
Figure 1.The lines of observation in both the experiment and the simulation are shown here, with up being north and right being east.of the two light sources starts 27 line segments which end on each of the 27 reflectors.The pale blue ellipse represents the region in which carbon dioxide sources were simulated.The red dot represents the position of the leak found in the experiment as calculated by maximizing the likelihood using the Ermak plume model.

Table 2 .Figure 2 .
Figure 2. The measured wind velocity during the experiment at 20 min intervals.The first and last measurements are labeled "1" and "24", respectively, which represent 5 February 2015 from 13:00 to 21:00 EST (local time).The average velocity is also given.

Figure 3 .
Figure 3. Contour plot of log likelihood for a source at (−10.30, 99.25) m with Q = 10 5 µmol mol −1 m 3 s −1 at a height z 0 = 10 m with x, y positions as in Fig. 1.In this figure, the log likelihood is referenced to that of the known true position.The maximum is shown with a black dot.The light paths are also shown.

Figure 4 .
Figure 4. Typical contour plot of the log likelihood when Q = 0 with x, y positions as in Fig. 1.The values range from 0 to 3.55.Since the 95th percentile of a 1 2 χ -square distribution with 3 degrees of freedom is about 3.9, Q = 0 is contained in the 95 % credible region for the parameters.Thus, the data are consistent with the hypothesis of no leak.

Figure 5 .
Figure5.Maximum log-likelihood values were acquired from n = 380 of 389 simulations with a leak rate of Q = 10 5 µmol mol −1 m 3 s −1 and a position chosen with uniform probability from the pale blue ellipse in Fig.1.CDF refers to the cumulative distribution function.The red line represents ordered random deviates that align perfectly with their corresponding χ 2 quantiles.If the points wander outside of the bounds given by the green lines, the assumption that the random deviates follow a χ 2 distribution is less tenable.
Figure6.Positions of 380 simulations with Q = 10 5 µmol mol −1 m 3 s −1 measured at a height of z 0 = 10 m, with the true positions sampled uniformly from the pale blue ellipse in Fig.1with the x, y positions given there.The reconstructed position is shown with a dot.Each line segment runs from the reconstructed position to the true position.(Many of these line segments are very short and appear as small black dots.)The blue dots represent the lowest one-sixth predicted leak rates, and the red dots represent the highest one-sixth predicted leak rates.The green dots represent the middle two-thirds.

Figure 8 .
Figure 8. Contour plot of the log likelihood of a source at a given position (in meters, with +x representing east and +y representing north) with the detectors at the height z 0 = 1 m, maximized over the leak rate Q, and referenced to Q = 0. (This convention implies that all values will be non-negative.)The maximum is in the red region at (−124, −85).The secondary maximum near (−190, 100) is hundreds of log-likelihood units lower, even if the region to the left of the displayed area is considered.The structure in the figure is aligned with the light paths shown in Fig. 1.

Figure 9 .
Figure 9. close-up of the plot shown in Fig. 8.The maximum is shown with a black dot.The colors are scaled the minimum and maximum in the graphed region, namely 5960 to 9535 loglikelihood units.The light paths which cross the region are also shown.

Figure 10 .
Figure10.Contour plot of the log likelihood with detectors at z 0 = 1 m based on considering only the data acquired in the eighth 20 min period with x, y positions as in Fig.8.These results are representative of how well the plume parameters can be determined with a traditional measurement.Although localization is good in one dimension, there is an uncertainty of roughly 100 m in the other.The estimated leak location from Fig.8is shown with a dot.The result is in marginal disagreement with these measurements, lying between the 95th and 99th percentile contours of given by the χ 2 2 distribution divided by 2.
along the straight-line path of the light.

Table 1 .
Parameters used in the calculation.Units are omitted for dimensionless quantities.