###### August 27, 2008

# Theoretical Evaluation of a Method for Locating Gaseous Emission Hot Spots

By Hashmonay, Ram A

ABSTRACT This paper describes and theoretically evaluates a recently developed method that provides a unique methodology for mapping gaseous emissions from non-point pollutant sources. The horizontal radial plume mapping (HRPM) methodology uses an open- path, path-integrated optical remote sensing (PI-ORS) system in a horizontal plane to directly identify emission hot spots. The radial plume mapping methodology has been well developed, evaluated, and demonstrated. In this paper, the theoretical basis of the HRPM method is explained in the context of the method's reliability and robustness to reconstruct spatially resolved plume maps. Calculation of the condition number of the inversion's kernel matrix showed that this method has minimal error magnification (EM) when the beam geometry is optimized. Minimizing the condition number provides a tool for such optimization of the beam geometry because it indicates minimized EM. Using methane concentration data collected from a landfill with a tunable diode laser absorption spectroscopy (TDLAS) system, it is demonstrated that EM is minimal because the averaged plume map of many reconstructed plume maps is very similar to a plume map generated by the averaged concentration data. It is also shown in the analysis of this dataset that the reconstructions of plume maps are unique for the optimized HRPM beam geometry and independent of the actual algorithm applied.

(ProQuest: ... denotes formulae omitted.)INTRODUCTION

Optical remote sensing (ORS) is a powerful technique for measuring air contaminant emissions from fugitive area sources.1-3 The radial plume mapping (RPM)-based methodologies use a single ORS monitor to collect path-integrated concentration (PIC) data from multiple beam paths in a plane and combine these with optimization algorithms to map the field of concentration across the plume of contaminant.4-14 Horizontal RPM (HRPM) methodology, which is the subject of this paper, was designed to map pollutant concentrations in a horizontal plane. This methodology is used to locate hot spots close to the ground. The HRPM methodology was first proposed4 and computationally evaluated in 1999, and later that year it was experimentally tested in a controlled chamber environment.5 After introduction and testing of the HRPM method, it was realized that the beam configuration had to be optimized to achieve stable results. The first generation of reconstruction approach relied on smooth basis function minimization (SBFM), in which a superposition of basis functions (bivariate Gaussian) are assumed to describe the field of concentration in the plane of measurement.6 The premise was that the HRPM configuration was inferior to conventional computed tomography (CT) configurations in terms of reconstruction accuracy, but the cost benefits of using one ORS monitor and fewer mirrors would justify the application of HRPM method in various scenarios. It was also assumed that the SBFM was the only algorithm that might provide reasonable reconstruction. When this premise was challenged, it was clear that any typical CT reconstruction algorithm might provide reasonable results.

In a computational simulation study,7 it was shown that a 16- beam optimized HRPM configuration was superior in terms of error magnification (EM, defined below in the theoretical evaluation section of this paper) to an idealized 18-beam CT configuration (would require 18 pairs of source and detector). This simulation compared the uniqueness and sensitivity to error for a wide range of plausible test maps for both beam configurations using the same CT reconstruction algorithm and the same spatial resolution (16 pixels). This led to investigation of the underlying kernel matrices in respect to EM as described below in this study.

Under the auspices of the U.S. Department of Defense's (DoD) Environmental Security Technology Certification Program (ESTCP) and the U.S. Environmental Protection Agency (EPA), an RPM methodology to directly characterize gaseous emissions from area sources has been demonstrated and validated, and a protocol has been developed and peer reviewed. This EPA "Other Test Method" (OTM-10) was made available for use on the EPA website in July 2006 (http:// www.epa.gov/ttn/emc/tmethods.html). The Other Test Methods category of the EPA Emission Measurement Center website includes test methods that have not yet been subject to the federal rulemaking process. Each of these methods, as well as the available technical documentation supporting them, have been peer reviewed and revised by the Emission Measurement Center staff and have been found to be potentially useful to the emission measurement community. This test method currently describes three methodologies, each for a specific use. As part of the development of this test method, the HRPM was validated in a series of outdoor experiments in which simulated hot spots were accurately located by the HRPM methodology.8,9

The methodologies in OTM-10 are independent of the particular path-integrated ORS (PI-ORS) system used to generate the PIC data. This test method can be applied to any scanning PI-ORS system that can provide PIC data. These may include the following: open-path Fourier transform infrared (OP-FTIR) spectroscopy, ultraviolet differential optical absorption spectroscopy (UV-DOAS), open-path tunable diode laser absorption spectroscopy (TDLAS), and path- integrated differential absorption LIDAR (PI-DIAL; light detection and ranging [LIDAR]). The choice of instrument must be made based on its performance relative to the data quality objectives of the study. The OP-FTIR and UV-DOAS technologies are widely used because of their capability of simultaneous chemical detection for a large number of gas species of environmental interest. However, when only a few gas species are of interest, it may be more beneficial to use other PI-ORS instrumentation, such as the TDLAS or PIDIAL. 3,15 For landfills, a scanning TDLAS system (hardware and software) was developed to identify hot spots of methane emissions.

In this paper, the current HRPM procedure is described, and the results of the computational study of the kernel matrices of the various beam geometries are discussed. Using methane PIC data collected in a landfill with a TDLAS system, it is also demonstrated in this paper that EM is minimal, because the averaged plume map of many reconstructed plume maps is very similar to a plume map generated using the averaged PIC data from the dataset. Analysis of this dataset also shows that the reconstructions of plume maps are unique for the optimized HRPM beam geometry and independent of the actual algorithm applied.

HRPM METHODOLOGY

This methodology was applied and demonstrated at various sites for locating emission hot spots in landfills.10-15 PIC data are collected over time with a single PI-ORS instrument, completing many cycles through the defined beam configuration. Because data are acquired sequentially, a moving average is required to reduce errors that originate from temporal variability. Typically, a moving average with a grouping of three cycles is sufficient to provide stable results.9-15 If the grouping of three cycles does not provide the desired quality of fit (discussed below), it is recommended to increase the number of cycles averaged in each group until the quality of fit is satisfactory. Once the PIC for all beam paths are averaged with the predetermined grouping of cycles for the gas species of interest, the RPM calculations make use of the information to reconstruct a plume map.

In landfills, the HRPM methodology is used to locate the source of fugitive emissions or hot spots, primarily methane. A rectangular area is defined around the ground location where the suspected gaseous emissions are originating. Ideally, the HRPM configuration will cover the entire suspected source area; however, this may not always be possible because of equipment limitations or site conditions (topography in landfills). Larger areas may need to be divided into smaller sections and studied separately.

Beam Geometry

The PI-ORS instrument is typically placed at the origin (in the first quadrant of the Cartesian convention) of the rectangular area to be measured (see Figure 1). For ease of presentation, determining components (PDC) are used to denote the component on the other end of the optical path from the PI-ORS instrument. Depending on the instrument selected, this could be a source, detector, retroreflecting mirror, or other reflecting object.

Once the HRPM measurement area and the number of path PDC have been determined, the area is divided into smaller rectangular areas called pixels. The total number of pixels required is smaller or equal to the total number of beam paths.

In Figure 1, the survey area is divided into nine pixels of equal size. It should be noted that the survey area might be asymmetric (e.g., 2 x 4 pixels, 3 x 5 pixels, etc.). Each pixel will have at least one optical beam path that terminates within its boundaries at a PDC (a retroreflecting mirror for most ORS instruments). This geometry maximizes the spread of the optical beams inside the area of emissions by passing one optical beam through the center of each pixel (as explained below). Dwelling time (the time that the ORS instrument dwells aimed on each PDC and collects PIC data) per PDC is determined by the specific study goals and the PI-ORS instrument- specific detection limits of the expected target gases. Typically, a range for dwelling time per PDC is 10-60 sec. Mathematical Formulation

Once the PIC for all beam paths are averaged with the predetermined grouping of cycles for the gas species of interest, the HRPM mathematical inversion procedure reconstructs a plume map over the area of interest. The measured PIC, as a function of the field of concentration, is given by:

... (1)

where K is a kernel matrix that incorporates the specific beam geometry with the pixel dimensions, k is the number index for the beam paths, m is the number index for the pixels, and c is the average concentration in the mth pixel.

Each value in the kernel matrix K is the length of the kth beam within the mth pixel; therefore, the matrix is specific to the beam geometry. Equation 1 is a reduction of the Fredholm linear integral equation of the first kind16 to a system of linear algebraic equations to be solved numerically. In mathematical terms, a solution is sought for f(x) from measured g(y), which is defined by the Fredholm linear integral equation:

g(y)= [integral]^sup b^^sub a^K(y,x)f(x)dx (2)

In earlier simulation studies, the multiplicative algebraic reconstruction technique (MART) was applied to solve for the average concentrations for each pixel.7,8,17 Currently, the HRPM procedure is executed by applying a non-negative least squares (NNLS) algorithm to achieve exactly the same results but more quickly (up to 100 times faster).18

The second stage of the plume reconstruction involves interpolation among the reconstructed pixel's average concentration, providing a peak concentration not limited to the center of the pixels. A triangle-based cubic interpolation procedure (in Cartesian coordinates) is currently used in the HRPM procedure.19 The HRPM procedure provides a plume map and calculates the location of the peak concentrations.

Quality of Fit

As described in earlier studies,4,5 the concordance correlation factor (CCF) is used to represent the level of fit for the reconstruction in the path-integrated domain (predicted vs. measured PIC). CCF is defined as the product of two components:

CCF = rA (3)

where r is the Pearson correlation coefficient, and A is a correction factor for the shift in slope and y-intercept of the linear regression.

This shift is a function of the relationship between the averages and standard deviations of the measured and predicted PIC vectors:

... (4)

where sigma^sub PIC^sub P^^ is standard deviation of the predicted PIC vector, sigma^sub PIC^sub M^^ is standard deviation of the measured PIC vector, ... is the mean of the predicted PIC vector, and ... is the mean of the measured PIC vector.

The Pearson correlation coefficient is a good indicator of the quality of reconstructing the location of hot spots. Typically an r close to 1 will be followed by an A very close to 1 but not vice versa. This means that the averages and standard deviations in the two concentration vectors are very similar, and the mass is conserved (reconstructed concentrations are accurate). When a poor CCF is reported (CCF

THEORETICAL EVALUATION

It is possible to evaluate an inversion method theoretically without running a computer simulation or without performing an experiment with real-world data. Most mathematical inversion problems in remote sensing (e.g., atmospheric temperature structure measurements from satellites, particulate size distribution from multispectral light extinction data, and CT) are described by the integral equation, eq 2, and therefore share some common characteristics. The most important characteristic is that the kernel function or matrix must be nonsingular and "stable" to minimize EM through the process of finding f(x).16 It is well established that the stability of the K matrix can be assessed a priori by calculating the determinant (Delta) when K is squared or by calculating the singular values when K is not squared. If the determinant value is very close to zero (the smallest singular value is very small), the K matrix approaches singularity. A singular matrix (Delta=0) cannot be inverted because the solution is not unique.16 If all singular values are nonzero, there is a unique solution for f for any given g. However, a matrix close to singularity is very unstable and will result in enormous EM in finding f(x). This means that a very small change in g will produce a totally different solution for f. The relative EM is defined as the percent error in f divided by the percent error in g. A good measure for the K matrix stability is the matrix condition number, which is defined as the ratio of the largest to the smallest singular values.18 An estimate for the average EM regardless of the solving algorithm16 is given by:

EM = M^sup -1^[the square root of]kappa (5)

where M is the size of g (number of measurements) and kappa is the condition number. The choice of the algorithm for finding f (including constraints such as non-negative values in f) and the error structure in g will determine the actual EM for finding f.

Table 1 demonstrates the stability of various kernel matrices (of HRPM beam geometries and of the CT beam geometry in the simulated study7 presented in 2000) through their condition number and estimated EM. Figure 2 depicts the two beam geometries compared in the simulation study.7 Generally, it was found that when each beam of the HRPM geometry was approximately passing through the center of the pixel where it terminated, the condition number of the kernel matrix was the lowest. As shown in Table 1, the EM dropped from a value of 0.69 for the 19-beam nonoptimized beam geometry to 0.28 for the 16-beam optimized configurations. This allowed for the optimization of the HRPM beam geometry using a trial-and-error approach.

All HRPM kernel matrices have a relatively small condition number, and therefore very low estimated EM. It is also shown in the table that optimized HRPM configurations have lower EM than nonoptimized HRPM configurations of the same size. Again, an optimized configuration (lowest EM) was found to be such where each beam goes through the center of the pixel where it is terminated (for an example of optimized beam configuration see Figures 1 and 2a). The kernel matrix of the ideal CT 18-beam configuration has a much larger condition number and therefore much larger average estimated EM than any of the HRPM beam configurations. This explains the finding of the earlier computational study7 that the HRPM is superior (in terms of the kernel matrix stability) to conventional CT regardless of the algorithm for finding the solution for f. The selected CT configuration represents the best case and artificially overdetermines (18 beams and 16 pixels) the solution to avoid bias toward the HRPM method. Any other practical configuration would perform much worse in terms of EM.

EVALUATION OF HRPM METHODOLOGY IN A LANDFILL USING TDLAS

As mentioned above, the HRPM method was primarily applied in landfills, mapping hot spots of methane emissions. 10-15 Because methane is the only compound of interest, a scanning TDLAS instrument is typically deployed for such measurement campaigns. The representative example shown below was arbitrarily chosen to demonstrate the method's stability to error and independence from the actual inversion algorithm. The beam geometry in this example is 15 TDLAS beams defined by retroreflecting mirrors. The dwelling time on each mirror was 10 sec, and the duration of each cycle was approximately 4.3 min. The entire run consisted of 34 completed cycles of measurement, which were collected over a period of more than 2 hr. Table 2 provides the angle from the horizontal axis, the path length, and averaged PIC for each of the 15 beams. Figure 3a shows the reconstructed plume map of the averaged PIC data in Table 2. A second averaged surface plume was generated using the same field data. To retrieve this map, 32 reconstructions were performed using a moving average of three cycles grouped together (cycles 1- 3, 2-4, 3-5, etc.). The resultant 32 plume maps were then averaged (by pixel) to yield one plume map, representing the average concentration distribution over the surface during the duration of the survey. This plume map is shown in Figure 3b and is essentially identical to the plume map in Figure 3a. Regardless of whether the averaging is performed on the input PIC data, or is performed on the temporally resolved output plume maps (pixels), the final averaged plume map for the duration of the survey is almost identical. It should be noted that such a successful comparison can be achieved only if all CCF values of the 32 individual reconstructions are close to 1 (>0.8). This demonstrates the stability of the HRPM configuration and proves the uniqueness of each solution (reconstruction). If this beam geometry was prone to error, a small change in the input data would provide significantly different plume maps, and the final result of the averaged 32 plume maps would have differed from the plume map generated directly from the averaged PIC data in Table 2.

A sample of seven plume maps representing a time interval of 17 min is given in Figure 4. It is demonstrated here that along this run the plume maps varied considerably. The seven plume maps on the left side of Figure 4 are reconstructions using the NNLS reconstruction algorithm, and the corresponding plume maps on the right side are reconstructions using the MART algorithm.

Table 3 compares the CCF of the 32 plume maps between the two algorithms. Generally, all CCF values are very good for this dataset, which confirms the grouping of three cycles for the moving average is more than sufficient. In runs where CCF values are marginal, it is strongly recommended to increase the number of cycles grouped together for the moving average. This most likely will improve the values of CCF and brings them closer to 1 on the expense of temporal resolution of the reconstructions. There is almost no difference between the CCF values between the two algorithms used. The CCF values generated using the MART algorithm are always equal or very slightly smaller than the CCF values generated by the NNLS algorithm. Furthermore, the NNLS is approximately 100 times faster and, therefore, is currently used. This comparison emphasizes that when beam geometry is robust, the solving algorithm becomes secondary in importance and that any algorithm that can achieve high CCF values can be applied in using this method. CONCLUSIONS

As shown in the examples, the theoretical evaluation results of the kernel matrices successfully demonstrated the stability of the HRPM beam geometry with regards to EM. Also, it was shown that HRPM configuration can be optimized using the kernel matrix condition number. An optimized configuration (lowest condition number) was found to be such where each beam goes through the center of the pixel where it is terminated. Through examination of experimental PIC data, the HRPM methodology was proven to be very stable because the average of 32 plume maps (34 cycles of measurements over more than 2 hr) was almost identical to a plume map reconstructed from the averaged PIC data from all 34 cycles. This can be achieved only with a very robust underlying inversion problem, and the importance of the reconstruction algorithm becomes secondary. This conclusion is reinforced by this study's demonstration that the NNLS and MART algorithms produced identical plume maps when performed on the same dataset. The HRPM methodology is almost independent of the reconstruction algorithm used, as long as the solution is constrained to positive values.

ACKNOWLEDGMENTS

The author would like to thank Gary Hater and Roger Green of Waste Management, Inc. for providing the TDLAS data used in this computational demonstration. Also, the author would like to recognize Michael Chase of ARCADIS for his development of the plume maps graphical presentation and for his help in computer programming.

IMPLICATIONS

HRPM has been applied mostly at landfills and a few wastewater treatment plants for locating hot spots of methane emissions. However, this method has applications for spatial characterization of other area and fugitive emission sources such as contaminated soil, industrial facilities for leak detection, and large urban areas, using a variety of optical remote sensing ORS technologies. This method is not suitable for mapping air contaminants emitted from the surface a body of water such as waste lagoons because it requires the setup of mirrors within the measurement domain.

REFERENCES

1. Walmsley, H.L.; O'Connor, S.J. The Use of Differential Absorption LIDAR to Measure Atmospheric Emission Rates at Industrial Facilities. In Proceedings of the A&WMA's International Conference on Optical Sensing for Environmental and Process Monitoring; A&WMA: Pittsburgh, PA, 1996; p 127.

2. Gronlund, R.; Sjoholm, M.; Weibring, P.; Edner, H.; Svanberg, S. Elemental Mercury Emissions from Chlor-Alkali Plants Measured by LIDAR Techniques; Atmos. Environ. 2005, 39, 7474-7480.

3. Grant, W.B.; Kagann, R.H.; McClenney, W.A. Optical Remote Measurement of Toxic Gases; J. Air & Waste Manage. Assoc. 1992, 42, 18-30.

4. Hashmonay, R.A.; Yost, M.G.; Wu, C.F. Computed Tomography of Air Pollutants Using Radial Scanning Path-Integrated Optical Remote Sensing; Atmos. Environ. 1999, 33, 267-74.

5. Wu, C.F.; Yost, M.G.; Hashmonay, R.A.; Park, D.Y. Experimental Evaluation of a Radial Beam Geometry for Mapping Air Pollutants Using Optical Remote Sensing and Computed Tomography; Atmos. Environ. 1999, 33, 4709-4716.

6. Drescher, A.C.; Park, D.Y.; Yost, M.G.; Gadgil, A.J.; Levine, S.P.; Nazaroff, W.W. Stationary and Time-Dependent Indoor Tracer- Gas Concentration Profiles Measured by OP-FTIR Remote Sensing and SBFMComputed Tomography; Atmos. Environ. 1997, 31, 727-740.

7. Hashmonay, R.A. Real Time Mapping of Indoors Gaseous Air Contaminants Using Radial Scanned Tunable Diode Laser Absorption Spectroscopy System. In Engineering Solutions to Indoor Air Quality Problem Symposium, VIP 98; A&WMA: Pittsburgh, PA, 2000; pp 332-340.

8. Hashmonay, R.A.; Wagoner, K.; Natschke, D.F.; Harris, D.B.; Thompson, E.L. Radial Computed Tomography of Air Contaminant Using Optical Remote Sensing. In Proceedings of the A&WMA's 95th Annual Conference, A&WMA: Pittsburgh, PA, 2002.

9. Varma, R.; Hashmonay, R.A.; Kagann, R.; Bolch, M.A. Optical Remote Sensing to Determine Strength of Nonpoint Sources: Duke Forest Validation Study; ESTCP Report No. A528044, Contract No. CP- 0214, AFRL-MLTY-TR-2005-4584, November 2005; available at http:// stinet.dtic.mil/ (accessed 2008).

10. Modrak, M.T.; Hashmonay, R.A.; Keagan, R. Measurement of Fugitive Emissions at a Region I Landfill; EPA-600/R-04-001; U.S. Environmental Protection Agency; Office of Research and Development: Washington, DC, 2004.

11. Modrak, M.T.; Hashmonay, R.A.; Varma, R.; Kagann, R. Evaluation of a Former Landfill Site in Fort Collins, Colorado Using Ground-Based Optical Remote Sensing Technology; EPA-600/R-05/042; U.S. Environmental Protection Agency; Office of Research and Development: Washington, DC, 2005.

12. Modrak, M.T.; Hashmonay, R.A.; Varma, R.; Kagann, R. Evaluation of Fugitive Emissions at a Former Landfill Site in Colorado Springs, Colorado Using Ground-Based Optical Remote Sensing Technology; EPA-600/R-05/041; U.S. Environmental Protection Agency; Office of Research and Development: Washington, DC, 2005.

13. Modrak, M.T.; Hashmonay, R.A.; Varma, R.; Kagann, R. Measurement of Fugitive Emissions at a Landfill Practicing Leachate Recirculation and Air Injection; EPA-600/R-05/088; U.S. Environmental Protection Agency; Office of Research and Development: Washington, DC, 2005.

14. Modrak, M.T.; Hashmonay, R.A.; Varma, R.; Kagann, R. Measurement of Fugitive Emissions at a Bioreactor Landfill; EPA-600/ R-05/096; U.S. Environmental Protection Agency; Office of Research and Development: Washington, DC, 2005.

15. Evaluation of Fugitive Emissions Using Ground-Based Optical Remote Sensing Technology; EPA-600/R-07/032; U.S. Environmental Protection Agency; Office of Research and Development: Washington, DC, 2007.

16. Twomey, S. Introduction to the Mathematics of Inversion in Remote Sensing and Indirect Measurements; Dover: Mineola, NY, 1996.

17. Chahine, M.T. Inverse Problems in Radiative Transfer: Determination of Atmospheric Parameters; J. Atmos. Sci. 1970, 17, 960-967.

18. Lawson, C.L.; Janson, R.J. Solving Least Squares Problems; Society for Industrial and Applied Mathematics: Philadelphia, PA, 1995; pp 158-165.

19. Barber, C.B.; Dobkin, D.P.; Huhdanpaa, H.T. The Quickhull Algorithm for Convex Hulls; ACM Trans. Math. Software 1996, 22, 469- 483.

Ram A. Hashmonay

ARCADIS, Inc., Durham, NC

About the Author

Dr. Ram A. Hashmonay is Director of the Advanced Air Monitoring Solutions Service at ARCADIS, Inc. He is also Chair of the Optical Sensing Division at A&WMA's Technical Council. Please address correspondence to: Dr. Ram A. Hashmonay, ARCADIS, 4915 Prospectus Drive, Suite F, Durham, NC 27713; phone: +1-919-544-4535; fax: +1- 919-544-5690; e-mail: [email protected]

Copyright Air and Waste Management Association Aug 2008

(c) 2008 Journal of the Air & Waste Management Association. Provided by ProQuest LLC. All rights Reserved.