A Framework to Survey a Region for Gas Leaks Using an Unmanned Aerial Vehicle

To assist gas distribution companies to effectively monitor their infrastructure and locate gas leaks, this paper considers the use of an Unmanned Aerial Vehicle (UAV) carrying a methane sensor to survey a region for gas leaks. As the UAV collects in situ measurements, it gathers evidence regarding the presence or absence of leaks in the region. To relate the UAV measurement to a region on the ground, we propose to use Upwind Survey Regions (USRs). If the UAV collects a measurement of high gas concentration, the USR represents a region that is likely to contain the leak. Likewise, if the UAV collects a measurement of low gas concentration, the USR represents a region that is likely to be clear of gas leaks. We propose a framework to process the measurements and produce a survey map indicating areas that can be reliably cleared of gas leaks and areas that may contain gas leaks. Our framework is composed of two steps: 1) mapping UAV measurements into USRs on the ground; and 2) fusing the various mapped USRs to produce the survey map. We discuss how USRs can be estimated and we test our framework using both simulated and real UAV flights.


I. INTRODUCTION
Natural gas production, transmission, and distribution sectors continuously deal with fugitive gas emissions.Finding and repairing the sources of these unintended emissions would mitigate economic loss [1] and help reduce the negative environmental impact of fugitive emissions on climate and air quality [2].
Traditionally, gas companies search and repair gas leaks using foot patrols [3].Using a portable methane sensor, technicians visit several areas of the region looking for gas leaks and repairing them.Because foot patrols are expensive and time-consuming [4], [5], gas operators and researchers have been considering alternative methods to detect gas leaks.As discussed in Section II, leaks can be detected The associate editor coordinating the review of this manuscript and approving it for publication was Manuel Rosa-Zurera.
by deploying an array of stationary sensors, by carrying a sensor while driving a vehicle or flying an Unmanned Aerial Vehicle (UAV) or aircraft, or by analyzing satellite images.These approaches are intended to assist gas companies in their repair efforts, guiding the foot patrols to the most important areas.
In this paper, we focus on one of the approaches being considered in the literature: surveying a region using an UAV equipped with a methane gas sensor flying inside the gas plume.As illustrated in Fig. 1, an UAV equipped with a GPS receiver, an anemometer, and a methane sensor collects gas concentration measurements at various points in the atmosphere (in situ measurements) while flying downwind of the region being surveyed.For safety and logistical regions, the UAV flies above the treeline, at heights ranging from 5 to 30 meters.For this reason, this framework is most suitable for surveying rural and suburban areas.Each concentration measurement is tagged and stored with its collection time, UAV location, wind direction and wind speed.From such information, the goal is to infer the presence of a gas leak in the region.Many previous works [6], [7], [8], [9], [10], [11] have proposed and studied this approach because methane disperses in the atmosphere and is carried by wind over long distances.This allows a UAV to sense a gas leak hundreds of meters away, enabling the survey of large areas.
Two important questions motivate our paper: when the UAV methane detector measures enhancement above a background concentration level, what can be inferred regarding the location of the leak source?Similarly, when the UAV methane detector does not measure enhancement above the background, which upwind regions can be declared free of gas leaks and to what confidence?Answering these questions is not easy because detection of high gas concentration can be highly intermittent due to the coupling between the filamentous nature of the gas plume and atmospheric turbulence [12], [13].A high methane concentration measured by the UAV only indicates that the UAV has encountered a gas filament and inferring the upwind gas leak location from such a measurement is not trivial [9], [10], [14].
To address this challenge, we use the concept of Upwind Survey Regions (USRs).A USR represents a region on the ground that relates to a UAV measurement.Its key purpose is to map a high or low concentration measurement onto a region: if the UAV collects a measurement of high gas concentration, the USR represents a region that is likely to contain a leak.Likewise, if the UAV collects a measurement of low gas concentration, the USR represents a region that is likely to be clear of gas leaks.Based on the results from existing literature [15], [16], we propose to use ellipsoidal USRs and we discuss in this paper how to estimate the USR dimensions using simulation tools.The concept of using a USR rather than an indication point to describe the presence or absence of leaks is an important step forward in doing auditable surveys where surveyed regions can be unioned together and the total region surveyed quantified.In particular, USRs are important for UAV surveys above the treeline where the upwind fetch for highly sensitive instruments can be many hundreds of meters.The concept of USRs has been used in the industry (e.g.Picarro Leak Indication Search Area (LISA) [17]); however, industry algorithms are often proprietary. 1he main goal of this paper is to propose a framework to enable the survey of a region using the USR concept.As illustrated in Fig. 2, our framework has two components: USR mapping and USR fusion.
• Regarding USR mapping, we use the dimensions of the USR to map each UAV concentration measurement (C mn ) of high or low gas concentration (detections or non-detections) into a region on the ground, as illustrated in Fig. 3.More precisely, we determine the set of locations in the ground that are within the USR relative to the UAV location.At the end of this process, each location in the ground will be covered by multiple USRs: some USRs will be associated with detections while others will be associated with non-detections (see Fig. 4 for an illustration).
• Regarding USR fusion, based on the information provided by the various USRs covering each location on the ground, we propose a statistical fusion procedure to classify each location as containing a gas leak or not.The final product of our USR fusion procedure is a survey map, where the region is partitioned in three areas: the areas that most likely contain a leak (R gas ), the areas that are likely free of gas leaks (R clear ), and the areas from which not enough information is available for a reliable decision (R unknown ).Ultimately, the survey map can be used by gas companies for guidance in infrastructure maintenance efforts and serve as an auditable instrument to document efforts to reduce fugitive gas emissions.It is important to highlight that the mapping of UAV measurements into USRs on the ground is not trivial because USR dimensions change with environmental and flight conditions.Thus, to apply our framework, the surveyor needs a database of USR models, one model for each set of environmental and flight conditions.A USR model is given by the dimensions of the ellipsoidal USR.Before starting the survey, the surveyor assesses the current environmental and flight conditions and selects the associated USR model from the database.
Our vision is to formulate a framework founded on analytical and simulation models to facilitate and guide aerial gas leak surveys.We envision a future in which companies  use our analytical and simulation models to deploy lowcost, small, agile UAVs that can fly just above the ground and reliably complete the survey of the region.As a first step toward this vision, we focus on a single UAV and use simulated data to test and guide our fusion framework.

II. LITERATURE REVIEW
Methane is a greenhouse gas used in industry and many households for various energy-related purposes such as cooking and heating.Methane is transported from production facilities using pipes and other infrastructure that is subject to aging and faults, causing unintended gas leaks.
Such leaks are not a rare occurrence.The current estimated amount of methane leaked in the atmosphere is significant: the Environment Protection Agency calculates that 140 Mmt (millions of metric tons) escape from natural gas systems each year [19] and some authors argue that the actual amount can be 60% higher than this estimate.
Verifying the presence of gas leaks is essential for gas distribution companies not only because of its safety and environmental impact, but also because of the financial loss incurred by gas leaks.It is estimated that the oil and gas industry loses US$2B every year only due to gas leaks [1].To reduce the environmental impact and losses, gas distribution companies periodically survey entire urban areas [1], [6], [20], [21].
To increase the leak surveying efficiency, additional solutions employing moving platforms have been tested and implemented, as we describe in the next subsections.They are complementary approaches intended to facilitate the survey process.Traditional methods using handheld gas sensors are still required to pinpoint leaks [4], [5].Likely, an ideal solution may be a tiered approach, combining multiple technologies to detect and repair leaks [3].
Stationary sensors offer continuous monitoring and low operational costs [27].However, monitoring the entire production, transmission, and/or distribution system with stationary sensors requires a large quantity of sensors.Furthermore, sensors close to the ground may fail to detect 1388 VOLUME 12, 2024 Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.leaks when heat and wind disperse the gas to higher altitudes, being out of reach of sensors [14], [35].For such reasons, stationary sensor solutions are typically deployed close to gas production facilities or high risk areas [5], [28].

B. VEHICLE APPROACHES
Another approach is to place highly sensitive (parts per billion) sensors on vehicles [17], [36], [37], [38], [39] to more rapidly cover large areas, such as urban and suburban distribution systems.Several operators have been exploring this approach to complement their surveys [40], [41].Tests of vehicle approaches demonstrate detection of leak rates below 100 g/hr [4]; however, these technologies often observe the same leak multiple times during a transit, leading to false positives, and may miss leaks that are far from a road [5], [42].For instance, the authors of [37] reported difficulties in detecting leaks more than 20m away from the vehicle.Because of this, the coverage of vehicular-based approaches is typically limited to areas close to roads [4], [17].The logistics of wind transport of the natural gas coupled with roadways from which to sample also imposes difficulty in consistent surveying methodology.

C. SATELLITE-BASED APPROACHES
Remote sensing from satellites is becoming more common to detect larger gas leaks from orbit [43], [44], [45], [46].The hyperspectral imagers implemented on satellites receive sunlight reflected by the Earth.Because methane absorbs certain wavelengths measured by these imagers, the satellite instrument can infer the presence of a leak within the region.
Satellites allow wide coverage and are effective in detecting large gas leaks (e.g.around 117 kg/hr [47]), making them ideal to detect super-emitters at regional or global scales [4], [45].On the other hand, satellite solutions suffer from cloud and snow cover and can detect only large leaks [5].Recent experiments have reported the ability to detect leaks as low as 100 kg/hr [44], [48]; however, gas operators are interested in detecting leaks as low as 36 g/hr [11].

D. AIRBONE-BASED APPROACHES
To improve detection of leaks, researchers have explored the option of performing sensing using aircraft or UAVs.The immediate benefit is that aircraft and UAVs fly at lower altitudes than a satellite, being able to sense gas leaks from a closer distance at much lower emission rates.There are two types of sensing in this category: detecting gas leaks while flying above the gas plume (remote sensing) and detecting gas leaks while flying inside the gas plume (in situ sensing).

1) REMOTE AIRBORNE SENSING
In this approach, the aircraft or UAV carries either a shortwave infrared (SWIR) passive imaging spectrometer measuring reflected sunlight [48], [49], [50], [51], a mid-IR infrared camera [52], or an active laser-based sensor [31], [53], [54], [55], [56], [57], [58], [59].Based on the amount of energy reflected, the sensor estimates the integrated gas concentration in the vertical column of air between the ground and the sensor.To enable detection, the aircraft or UAV needs to fly above the target so that the gas plume is within the field of view of the sensor [57].
Several authors have reported experiments in which this technology was able to detect leaks on the order of 1 kg/hr when flying aircraft at 200m altitude [56], [57] and leaks of 155 kg/hr when flying at 5km altitude [49].In the Stanford/EDF Methane Challenge, an aircraft flying at 850m altitude using a LIDAR sensor was able to reliably detect leaks above 9 kg/hr [4].
Drawbacks of using aircraft are high flying cost, lack of maneuverability [12], need for skilled pilots to fly the aircraft above targets [57], and reduced detection performance when flying at required high altitudes [5].
UAVs allow lower flying altitudes (e.g.Class G airspace not available to manned aircraft [14]) and authors have reported the ability to detect with LIDAR sensors leak rates as low as 5 kg/hr when flying the UAV at altitudes ranging from 3m to 15m above the leak [53] and 13 kg/hr when flying at an altitude of 50m [31].
While reduced altitudes improve detection of small leaks, it also reduces the area within the field of view of the sensor, requiring the UAV to fly right above the gas plume for reliable detection [51].

2) IN SITU AIRBORNE SENSING
An alternative approach, which we follow in this paper, is to fly the aircraft or UAV inside the gas plume [6], [7], [8], [9], [10], [11], [12], [42].In this approach, the aircraft or UAV carries a sensor that measure in situ concentration at the aircraft/UAV location as it flies downwind from the leak.Open or closed path laser absorption spectroscopy are often used to measure the concentration level, enabling detections of concentration at the part-per-billion (ppb) level [8], [11], [60].
Many authors have performed experiments highlighting the potential of this approach.The authors in [11] reported the ability to detecting leak rates as low as 1 kg/hr from a distance of 280m downwind of a leak source.The authors in [9] reported the ability to sense enhancement as low as 20 ppb above the background methane concentration up to 4km downwind of known landfills and flare stacks.Other authors [6] reported being able to measure the plume as far as 10km downwind of known landfills.In the Stanford/EDF Methane Challenge, UAVs using in situ measurements were able to reliably detect leaks above 40 g/hr while flying at an altitude of 1-3m near a test facility [4].
The ability to detect gas at low leak rates, from a distance, and at low altitudes are some of the advantages of flying inside the gas plume, making it a viable solution to quickly screen a region for gas leaks.
As in other approaches, there are also some challenges in implementing this solution.The authors in [8] and [42] reported that in situ measurements are highly intermittent, requiring multiple passes to collect additional data to confirm the existence of a plume.The authors in [7] reported difficulties in detecting the plume: in a controlled experiment with a single gas leak source and considering the UAV flying at 40-50m, the authors evaluated the probability of detecting the gas plume when flying 750m to 2000m downwind from the source.Even with a leak rate of 1 ton/hr, the authors reported detection probabilities as low as 0.21, highlighting the complexity of plume detection.
When deploying this solution in UAVs, the propellers may dilute the gas before it reaches the sensor [10].To mitigate this signal reduction problem, one can place the sensor in front of a UAV with a wide enough frame [12], [35].The authors in [35] further claimed that the UAV needs to fly into the plume with a speed of at least 2m/s to make the UAV propellers' influence in the concentration negligible.
Given the turbulent and intermittent nature of concentration measurements, elaborate techniques were proposed to detect and localize the leak [14].The authors in [6], [9], and [10] used the in situ measurements together with wind information and environmental characteristics in a sequential Bayesian Markov Chain Monte Carlo method to estimate the number, strength, and location of leaks.In [6] and [9], the measurements are used in the Metropolis-Hasting algorithm to update the posterior distribution for these variables.In [10], the posterior distribution is updated with particle filters.
There has also been research work on search algorithms to guide the UAV towards the gas source [61], [62].The authors in [62] propose that the UAV choose the direction that provides the highest expected gain in information.The authors in [61] propose that the UAV follow a three-phase approach: plume acquisition, plume tracking, and source declaration, after which the search task is complete.Each phase involves a different flight pattern and a different way of processing the UAV measurements to determine the next UAV position.The authors in [63] proposed a UAV path-planning algorithm to cover a region and build a map of the gas concentration in the region.Their random-tree based algorithm aims at directing the UAV towards the points of highest concentration, focusing on the scenario without wind.

E. OUR CONTRIBUTIONS
In this paper, we focus on how to perform surveys using UAVs flying inside the gas plume.Although this concept is not new, our paper adds to the existing body of knowledge by providing a formal framework to process the UAV measurements.More precisely, three important aspects differentiate our paper from the existing literature: 1) Our framework is not restricted to locating a gas leak among a small set of candidate leak locations; instead, our method can determine regions with possible gas leaks wherever they occur in the region.
2) Instead of attempting to estimate the leak location as quickly as possible, our framework aims at guiding gas distribution companies by producing a survey map.Such a map would not only indicate the areas that are likely to contain a gas leak, but also indicate areas that can be reliably cleared of gas leaks.With the survey map, the gas company can concentrate its foot patrols or other surveying methods in a smaller region to pinpoint the exact location of the leak.Furthermore, approaches that attempt to estimate the leak location typically assume that the number of leaks is known.By producing a survey map, our approach can identify multiple regions that contain possible leaks.3) We consider a new statistical procedure to fuse the USRs covering each location in the ground in order to determine whether it can be reliably cleared of gas leaks or not.
We conclude this section by highlighting that our approach is not intended to replace foot patrols or other technologies.Instead, our approach should be seen as complementary to other technologies.The main goal of our approach is to screen areas both for the presence of leaks as well as the absence (to a certain confidence level) of leaks and facilitate the work of foot patrols by guiding them to the area most likely to contain gas leaks.

III. MODEL DESCRIPTION AND METHODOLOGY OVERVIEW
Consider a region R to be surveyed for gas leaks.For convenience, we assume R is the plane (x, y, 0) within a threedimensional space.
If K ≥ 0 gas leaks are present in R, we assume that they release unpressured gas that propagates in the atmosphere over time, producing an instantaneous gas concentration field C(p, t) at each location p := (p x , p y , p z ) and time t.
To survey R for gas leaks, we propose to use an UAV.Along with a GPS, anemometer, and communication equipment, the UAV carries an open path laser spectrometer (OPLS) sensor.As the UAV flies over R, the OPLS collects multiple in situ measurements of the instantaneous gas concentration field tagged with 3D location and time; i.e., when the UAV is at position p at a time t, the OPLS obtains the value of C(p, t).As discussed in [12] and [35], it is possible to mount the OPLS sensor in the front of a UAV with a wide enough frame such that the UAV propellers' influence in C(p, t) is negligible when then UAV flies into the plume with a speed of at least 2m/s.We assume such conditions in this paper.
At each measurement, the OPLS sensor compares each C(p, t) against a threshold C min .If C(p, t) > C min , we say that a detection occurred; otherwise, we say that a non-detection occurred.We note that C min should be set high enough to detect enhancement above the background concentration level.
The end goal of our method is to use the sequence of detections and non-detections at various UAV positions to partition R in three subsets: • R clear : set of points that are likely to be clear of leaks; • R gas : set of points that may contain a leak; and • R unknown : set of points for which no decision is made.The principle is to use each measurement of instantaneous gas concentration to gather evidence about the absence or presence of a leak in each part of the region R; and the first issue to address is: how to map each detection and non-detection into a subset of R? In other words, a detection at UAV position p represents an evidence that there is a leak somewhere upwind of p, but where?Likewise, which parts of R can be cleared of leaks when a non-detection occurs?
Mapping detections and non-detections into subsets of R that are likely to contain a leak and subsets that can be cleared is not trivial because the UAV measurements are collected downwind of the region being surveyed.For this, we need a practical way to map UAV measurements into subsets of R.This is the motivation for the concept of Upwind Survey Regions.

A. THE UPWIND SURVEY REGION (USR)
To define the USR, we use concepts from the theory of footprint modeling [16].Meteorologists are often interested in studying gas exchanges close to the surface.For this, various sensors are placed in the top of towers looking over a forest or a region.Since the sensor measurements represent the concentrations at the sensor location, meteorologists are interested in estimating the relative measurement contribution coming from the sources distributed in the underlying region.For this, they define the concept of footprint function: 2 The footprint function f p (x, y) represents the transfer function between a source at location (x, y, 0) and a sensor located at p [65].The footprint function f p can vary not only with the sensor height (p z ), but also with various meteorological parameters, such as wind speed and the inverse Monin-Obukhov (MO) parameter [66,Chp. 4].
Letting Q ss (x, y) be the distribution of source strength, the long-term average concentration measured at p is obtained from It is important to highlight that C avg (p) does not represent the instantaneous concentration C(p, t) collected by the UAV.Rather, it is the long-term average concentration over a period of about 10 minutes [66,Chp. 4].The instantaneous concentration C(p, t) collected by the UAV varies because the gas disperses in filaments with random paths [12], [13].
Considering that gas leaks are point sources, Q ss (x, y) can be modeled with impulses.For instance, if a leak of strength s (2) We use the concept of footprint function to define the USR for a location p.Assuming that the surveyor needs to detect leaks with strength s min or higher, the Upwind Survey Region (USR) for a location p is defined as Intuitively, the USR(p) is used to make inferences regarding the presence or absence of leaks within the region covered by the USR(p).If the footprint function f p were known and the UAV's OPLS sensor collected a long-term average measurement C avg (p) ≤ C min , then the surveyor could infer that the region covered by the USR(p) does not contain a gas leak of strength s min or higher.On the other hand, if the UAV collected a long-term average measurement C avg (p) > C min , then the surveyor could infer that a gas leak with strength s min or higher may be present inside the region covered by USR(p).
Defining the USR in terms of the footprint function f p is useful because results from the footprint modeling literature provide guidance regarding its shape.More precisely, previous results suggest an ellipsoidal shape for USR(p) when p is above the roughness sublayer [15], [16].Also, when R is homogeneous enough, f p does not depend on p x or p y [15], [16].
Unfortunately, analytical forms for the footprint function f p are difficult to obtain, making it difficult to analytically estimate the dimensions of the ellipsoidal USR; however, it is possible use simulation tools to derive models for the footprint function [15].In our work, we use the Quick Urban and Industrial Complex (QUIC) simulation software [67] from Los Alamos National Laboratory to estimate the dimensions of the USR under various environmental and flight conditions.For each set of conditions, we estimate and store a USR model consisting of the fetch start (d start ), fetch end (d end ), and width (w) for the ellipsoidal USR, as illustrated in Fig. 3.We provide details about the QUIC simulator in Appendix A and we describe a procedure to estimate the USR dimensions and build a database of USR models in Appendix B. We note that the procedure to build the database of USR models is done only once: after the database is built, USR models can be applied indefinitely when performing surveys in similar environmental and flight conditions.
Even with the database of USR models, using the USR to perform surveys is not trivial because the UAV does not collect samples of the long-term average C avg (p); instead, it collects samples of the instantaneous concentration C(p, t); however, our methodology addresses this issue to enable the production of survey maps from UAV measurements.

B. METHODOLOGY OVERVIEW
Our framework is illustrated in Fig. 2 and consists of two components: USR mapping, and USR fusion.
• USR mapping: This is the component that maps high concentration measurements (detections) or low concentration measurements (non-detections) collected by the UAV into the USR region.Each measurement has an associated GPS location for the UAV and the USR region is defined with respect to the GPS location.
Multiple USRs are mapped into the ground; and, for each subregion, the procedure tracks how many USR associated with detections and non-detections cover the subregion.Details are given in Section IV.
• USR fusion: Because gas disperses in filaments with random paths [13], UAV measurements are random and USR mapping may provide conflicting information; i.e., a ground subregion may be covered by both USRs associated with detections and USRs associated with non-detections.The USR fusion procedure statistically combines the various information to decide whether a ground location should be included into R gas , R clear , or R unknown .Details are given in Section V.

IV. USR MAPPING
The first step in translating UAV measurements into a survey map is USR mapping; i.e., each UAV measurement collected at a location p is mapped into the region given by the USR(p).
As USRs are mapped on the ground, we track which USRs covered each subregion being surveyed.
Start by dividing R with a uniform grid with any desired granularity, producing the set of subregions {{r ij } I i=1 } J j=1 .Assume the UAV flying at a constant altitude p z in a sweeping flight pattern: as illustrated in Fig. 1, the UAV flies in a direction perpendicular to the wind direction while periodically changing its direction by 180 • and periodically moving upwind.Assume that the UAV changes direction when it reaches the boundary of the space being surveyed.After performing a number of flights in the same perpendicular line, the UAV moves d upwind meters parallel to and against the wind direction and performs a new set of perpendicular flights.This process is repeated as many times as desired.We shall refer to each of these perpendicular flights as a crossflight.
Let M be the total number of crossflights.For each crossflight m, assume that N instantaneous gas concentration measurements are collected.Let p mn , C mn , and t mn respectively represent the UAV position, the instantaneous gas concentration, and time the measurement was collected at the n th measurement in the m th crossflight; i.e., C mn := C(p mn , t mn ).
For every position p mn of the UAV, the USR mapping process determines the USR(p mn ) by considering the wind direction and meteorological conditions.More precisely, based on the UAV height, the inverse MO length, and the wind speed, we first interpolate the USR's fetch start d start , fetch end d end , and minor axis width w from the existing USR models in the database.Let (u x , u y , u z ) be the unitary wind direction vector, and use this vector to determine USR(p mn ) as the ellipsoidal region in the ground plane with major axis starting at (p mn,x − u x • d start , p mn,y − u y • d start ) and ending at (p mn,x − u x • d end , p mn,y − u y • d end ); and minor axis perpendicular to (u x , u y ) and with width w.When the inverse MO length, the wind speed, the wind direction, and the UAV height are constant, the mapped USR(p mn ) have the same ellipsoidal dimensions and same orientation but shifted in space by the UAV location, as illustrated in Fig. 4.
It should be noted that the USR mapping process can also be used when meteorological conditions change during the survey.For instance, if the wind speed changes in between UAV measurements, then USRs with different dimensions would be mapped into the ground.Likewise, if the wind direction changes in between UAV measurements, then the mapped USRs' orientation would also change accordingly.
The USR mapping process associates each USR(p mn ) with its instantaneous concentration measurement C mn ; and compares them against the threshold C min to determine detections (C mn > C min ) and non-detections (C mn ≤ C min ).Each detection represents evidence for a gas leak at subregions r ij within the USR(p mn ).Likewise, each non-detection represents evidence that subregions r ij within the USR(p mn ) are clear of gas leaks.
We note that, given that each USR provides information about subregions located in between d start and d end from the crossflight line, each crossflight must be spaced by no more than d end − d start to ensure that all desired subregions r ij are surveyed; i.e., whenever the UAV moves d upwind parallel and against the wind direction to start a new crossflight, it must move d upwind < d end − d start .
We further note that, if the surveyor knew the strength s of the emitter and if s were greater than the smallest leak strength s min used to define USR(p) in (3), then a larger USR could be used to map measurements on the ground, enabling the survey of larger areas; however, the emitter strength is generally unknown and the conservative option is to use USR(p) corresponding to s min .Using such a USR(p) to map measurements on the ground limits the number of subregions for which the USR provides evidence for or against the presence of gas leaks; however, it results in more conservative survey maps.
After all measurements are collected and all USRs are mapped, we determine if each subregion r ij should be placed into R clear , R gas , or R unknown by fusing the evidence provided by detections and non-detections from all USRs covering r ij .The fusion process that we propose is discussed precisely in the next section.

V. USR FUSION
If all USRs that cover a subregion r ij correspond to detections, then it would be reasonable to place r ij in R gas .Likewise, it would be reasonable to place r ij in R clear if all USR that cover r ij correspond to non-detections.
The challenge is that USRs may provide conflicting information; i.e., it is often the case that one or more USRs corresponding to detections cover r ij while other USRs corresponding to non-detections cover the same r ij , as illustrated in Fig. 4. To understand this, recall from Section III that the concept of USR is defined based on long-term average gas concentrations; however, UAVs collect measurements of instantaneous gas concentration.Instantaneous values of concentration are random and can be significantly different from long-term averages because of the random filament nature of gas dispersion [13].As a result, it may happen that a UAV at p mn is inside the gas plume produced by a leak in (x * , y * ); i.e., the long-term average concentration gas plume is such that C avg (p mn ) > C min , but the instantaneous concentration collected by the UAV may be weak; i.e., C mn < C min .In this case, C mn would correspond to a non-detection, providing incorrect evidence that USR(p mn ) does not contain a leak.Likewise, it may also happen that a UAV at p mn is outside the gas plume produced by a leak in (x * , y * ) but still collects an instantaneous measurement C mn > C min due to a gas filament moving out of the average gas plume.In this case, C mn would correspond to a detection, providing incorrect evidence that USR(p mn ) contains a leak.These phenomena are illustrated in Fig. 5.
Due to the randomness in C mn , we propose to fuse the evidence provided by the various USRs containing r ij to decide whether to place r ij in R gas or R clear using a statistical decision framework.
A first motivation for the statistical decision framework is that it has mechanisms to control the different types of error.In our context, there are two types of decision errors when surveying each point r ij in R: ) is free of leaks.
• False Negative: placing r ij in R clear when r ij should be placed in R gas ; i.e., declaring r ij to be clear of leaks when a gas leak exists in r ij ; or • False Positive: placing r ij in R gas when r ij should be placed in R clear ; i.e., declaring a gas leak is present in r ij when r ij is clear of leaks; and we argue that, among these two types of errors, any survey should aim to keep the probability of any false negative below a small value α (e.g., α = 0.05); i.e., we propose that the survey method clears a subregion r ij of gas leaks only under strong evidence.
A second motivation for the statistical decision framework is that it combines random evidence.In our context, if the USR includes r ij and is associated with a detection, then the surveyor has evidence for placing r ij in R gas ; but if the USR includes r ij and is associated with non-detection, then the surveyor has evidence for placing r ij in R clear .When multiple USRs include r ij , some with high and some with low concentration measurements, a statistical framework provides guidance in how to conciliate the contradictory evidence.
To decide if each r ij should be placed in R gas or R clear , let M ij be the number of crossflights with at least one USR covering r ij .More precisely, letting 1{condition} = 1 when the condition is true and 1{condition} = 0 when the condition is false, Let e (4) Because of the randomness in C mn , we have that M ij , e (m) ij and, therefore, s ij are all random outcomes.Let S ij represent the random variable whose outcomes is s ij and define which is the 'p-value' for the outcome s ij .P-values are typically used by statisticians to indicate the level of evidence against a hypothesis [68].Lower values for a p-value indicates that an outcome is unlikely to occur when the hypothesis is true.In here, P ij indicates the level of evidence that r ij contains a leak and a very low value for P ij provides strong evidence to place r ij in R clear . 4We discuss how to compute P ij in the next section.
From M ij and P ij , we test each r ij individually by comparing its p-value P ij against a threshold α ind : where M min defines the minimum number of crossflights below which we refrain to make a decision given the low amount of information.We highlight that, because of the multiplicity of tests [69], α ind needs to be adjusted to keep the probability of placing any r ij with a leak in R clear below α.
We discuss how to adjust α ind in Section V-C.After testing each subregion, the surveyor builds the survey map, indicating R gas , R leak , and R unknown areas.

A. COMPUTING P-VALUES
To compute P ij , we consider that each e i.e., q is the probability that the UAV collects at least one measurement C mn > C min while crossing the gas plume produced by a leak at r ij in any of the crossflights.We assume that the region is homogeneous such that q does not depend on r ij .This is a simplified model because q does not take into account the distance between the crossflight and r ij ; however, we adopt this model for this initial study.In this model, = 1 among all crossflights that cover any r ij .
We assume that {E m=1 are statistically independent random variables conditioned on the M ij crossflights both when a leak is present at r ij and when a leak is absent at r ij .This assumption is reasonable when a UAV performs the sweeping flight pattern because the configuration of gas particles in the flight path changes from one crossflight to the next.
It then follows that, when conditioned on M ij crossflights, ij has a Binomial distribution with parameters M ij and q [68]; and P ij used in the decision rule ( 6) is obtained from the cumulative distribution function of the Binomial(M ij , q) distribution computed at s ij .

B. ESTIMATING q
Given the difficulty in estimating q analytically, we propose two possible ways to estimate q.
One way to estimate q is to use UAV flight measurements collected when the leak location is known.For instance, the surveyor would place a gas canister in a region r ij , fly a large number M ij of UAV crossflights covering r ij .From the e (m) ij obtained, estimate q with q = (1/M ij ) m e (m) ij .Such a q would be included in a georeferenced database and used in future UAV measurements to survey the same region for actual leaks in unknown locations.
Another way to estimate q is to consider it to be a nuisance parameter in the detection problem [70].In this approach, one estimates q using the same UAV flight measurements used to perform the survey.Since the leak location is unknown, we propose that the surveyor first estimate q as if the leak were present in each subregion r ij .Letting qij represent such estimate, we propose to use qij = (1/M ij ) m e (m) ij , where M ij is the number of USRs covering r ij .Since the subregion containing a leak will have qij among the highest values, we propose to adopt q = max ij qij for the fusion procedure.Such a procedure should provide low-variance estimates when M ij is high.

C. ADJUSTING α IND
To understand how to adjust α ind in the decision rule (6) to keep the probability of placing any r ij with a leak in R clear below α, it is important to discuss how our survey procedure relates to the multiple hypotheses testing problem [69].
Consider testing G hypotheses {H g } G g=1 using G independent tests.The multiple hypotheses testing problem refers to the probability of rejecting any H g when H g is true.Assume that each of the G tests is performed such that P[reject H g |H g is true] = α.It then follows that P[reject any true H g |G 0 hypotheses H g are true] = 1− (1−α) G 0 , which grows to 1 as the number G 0 of true hypotheses H g grows [69].To avoid this problem, one can apply the Bonferroni correction [71]: assuming the number of true hypotheses H g can be as large as G, perform each individual test such that P[reject H g |H g is true] = α/G.Although the 1394 VOLUME 12, 2024 Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.
Bonferroni correction ensures that the probability of rejecting any true H g is below α, it reduces the probability of rejecting H g when H g is false.
The multiple hypothesis testing problem is also present when performing surveys because the surveyor is testing whether each one of I • J subregions r ij contains a leak or not.More precisely, letting H ij refer to the hypothesis that r ij contains a leak, to keep P[reject any true H ij |I • J hypotheses H ij are true] ≤ α, one would need to set α ind = α/(I • J ).However, it is not necessary to use such a small α ind because it is known that the number of leaks (K ) is supposed to be much lower than I •J .This makes the multiple hypothesis testing problem much less severe in our survey problem.
If the surveyor knows that no more than K max leaks exist in the region; i.e., K ≤ K max , then adopting α ind = α/K max ensures that To see this, follow the same proof for the Bonferroni correction [71, pg.350]: let K be the set of subregion indexes (i, j) that contain a gas leak and P[any r ij with a leak placed in If the surveyor cannot assume a value for K max , then K max = I • J could be used to define α ind ; however, α ind would be too small and the survey would likely generate a very large R gas region.
Given that the final purpose of the survey map is to provide guidance, we suggest that the surveyor adopt any reasonable value for K max (even K max = 1) and accept that the probability of placing any r ij with a leak in R clear may be greater than α when the actual number of leaks is greater than K max .

VI. EVALUATING SURVEY PERFORMANCE
Since the survey map will guide maintenance decisions, it is important that gas leaks be present in R gas so that maintenance crews do not miss any leak.This motivates our first performance metric: the probability of placing any r ij with a leak in R clear .As explained in the previous section, this probability depends on the maximum number of leaks assumed in the region.Therefore, we evaluate our survey with P[any r ij with a leak placed in R clear |K leaks in R] (9) for any desired value for K .We propose to estimate (9) by producing survey maps in various test flights with leaks placed in various random locations.In each test, we fly the UAV, collect measurements C mn , and execute the procedure to produce the survey map.Based on the various survey maps, we estimate (9) with the ratio of the number of survey maps containing any leak in R clear to the total number of survey maps.
However, it is important to highlight that ( 9) is not sufficient to evaluate the quality of the survey.To see this, consider an extremely conservative framework that places all subregions in R gas .In this case, the probability that any r ij with a leak is placed in R clear is 0; however, such a survey map does not provide valuable information to a surveyor.
Therefore, to increase the efficiency of maintenance crews, it is important that R gas be as small as possible.For this reason, we also evaluate the quality of the survey map through the size of the R gas region.More precisely, letting | • | represent the number of subregions within a set, we also compute |R gas |/|R|.
We note that |R gas |/|R| is closely related to the accuracy (Acc) of the survey map.Recall that Acc = (TP + TN )/(TP + TN + FP + FN ), where TP, TN , FP and FN respectively represent 'true positives', 'true negatives', 'false positives', or 'false negatives'.In our context, a TP is a subregion that contains a leak and was placed into R gas ; a TN is a subregion that does not contain a leak and was placed into R clear ; a FP is a subregion that does not contain a leak and was placed into R gas ; and a FN is a subregion that contains a leak and was placed into Ideally, it is desirable that survey methods have both the probability of placing any r ij with a leak in R clear and |R gas |/|R| close to 0.

VII. RESULTS
In this section, we provide examples showing how our framework can be used to produce gas leak survey maps.We evaluate our framework using both simulated data and real data.

A. TESTING SURVEY FRAMEWORK USING QUIC SIMULATED DATA
For these tests, we configured QUIC to simulate gas fields spreading in an open field of 200m of width (x-dimension), 1000m of depth (y-dimension), and 20m of height (z-dimension).One, two, or three gas leaks releasing unpressured gas at a constant rate of 100 g/hr were placed at ground locations with a southbound wind spreading the gas particles in the 3D region.The wind speed during simulations was 4m/s and the inverse-MO was −0.11/m.In each simulation, we configured QUIC to produce gas fields with 6 million particles, producing gas fields that approximate the long-term average concentration at each point in the space.Before collecting UAV measurements, we ran the simulation for 800 seconds to ensure that it reached the steady-state condition.
We simulated the UAV flying at a constant height of z = 5m in various lines perpendicular to the wind, with each line at a different distance to the leak.From the average gas concentration produced by QUIC, we simulated the instantaneous concentration measurements C(p, t) collected by the UAV using a Poisson model, in a manner similar to [62], [72], [73], and [74].More details can be found in Appendix C. Each instantaneous measurement C(p, t) was compared against C min = 6.6 • 10 −6 g/m 3 to determine detections and non-detections.
Using a database of USR models estimated as described in Appendix B-A and shown in Tables 1 and 2 of this appendix, the ellipsoidal USR model interpolated for the conditions above had 1.75m of fetch start, 106.78m of fetch end and 35.00m of minor axis.This USR model was used to map the USR(p mn ) for the detections and non-detections obtained at the various UAV locations in all the test simulations that follow.

1) EXAMPLES WITH 1 LEAK
For this first example, we executed our procedure (USR mapping and USR fusion) to produce a survey map while considering a single gas leak placed at (x * = 70, y * = 706, z * = 0).We highlight that our procedure produced the survey map without knowledge of such a leak location.
We simulated the UAV flying in the 3D gas field at a constant height of z = 5m in lines perpendicular to the wind.Starting at the location p 0 = (20, 20, 5), the simulated UAV collected 160 measurements, one per meter, as it moved to the location (180, 20, 5), completing a crossflight.After 10 crossflights at a fixed y, the simulated UAV moved upwind by 100m before doing a new set of 10 crossflights, until it completed the last crossflight at y = 920m.
To produce the final survey map, the 200m-by-1000m ground surface was partitioned using a regular grid, producing 7m-by-7m subregions r ij .
For each measurement, we used our USR mapping procedure of Section IV to map each detection and non-detection USR into the region, tracking the number of detection and non-detection USRs covering each subregion r ij .After completing all crossflights, we applied the fusion procedure of Section V to place each r ij in either R gas , R clear , or R unknown considering α = 0.05 for the maximum probability of any false negative.For this, we first computed for each subregion r ij .Using K max = 1 to define the decision threshold α ind in (6) and without using the knowledge of the simulated leak location, we applied the decision rule (6) to place each r ij in either R gas , R clear , or R unknown .To compute each P ij , we estimated the fusion parameter q using the first approach described in Section V-B; i.e., by using separate simulations under the same conditions as above with the leak at a known location.Such a procedure resulted in the estimate q = 0.715.Fig. 6 illustrates the resulting and the corresponding p-value P ij for each subregion; and, using the decision rule (6) with M min = 6, the resulting survey map is shown in Fig. 7.
1396 VOLUME 12, 2024 Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.Fig. 7 indicates that the procedure was able to include the actual gas leak location into R gas and was able to clear many subregions r ij that do not contain gas leaks.As discussed in Section VI, we also computed |R gas | = 7, 007 and |R| = 183, 211, resulting in |R gas |/|R| = 0.038 for the survey map.
We performed the procedure above 60 times, placing the leak at 60 random locations (x * , y * , z * = 0) within the rectangle given by the vertices (30, 400) and (170, 750).As before, the leak location during the simulation is not used by the procedure to produce the final survey map.For each of the 60 experiments, we inspected the survey map computing its |R gas |/|R| and checking if R gas contained the leak location.
Among the 60 experiments, the final survey map had R clear containing the leak location in 4 of them, which means P[any r ij with a leak placed in R clear |1 leak in R] ≈ 0.067.We also computed the mean and standard deviation of |R gas |/|R| using R gas and R obtained from each of the 60 experiments, resulting in |R gas |/|R| = 0.05 ± 0.02.The experiments resulting in the lowest and highest |R gas |/|R| ratios (0.015 and 0.08) are illustrated in Fig. 8.

2) EXAMPLES WITH 2 AND 3 LEAKS
We also tested our procedure when the region contains more than 1 leak.We again used the QUIC simulation tool to produce simulated gas leaks in the same 3D region of (200mx1000mx20m) but now with 2 or 3 leaks.
Considering the same flight pattern and the same fusion procedure as in the example with 1 leak but using K max = 2 and K max = 3 to define the α ind of (6) for the case with 2 and 3 leaks respectively, we obtained the survey maps shown in the top two displays of Fig. 9.In these figures, the UAV performed 100 crossflights.As can be seen, the procedure was able to produce a survey map with R gas including the leaks; however, the procedure was not able to clear many areas that do not contain a leak.This result can also be seen in the resulting |R gas |/|R|: 0.176 for the 2-leak scenario and 0.292 for the 3-leak scenario.
One way to reduce |R gas |/|R| is to gather more evidence by performing additional crossflights.The middle and bottom displays of Fig. 9 respectively show the survey maps obtained when the UAV performs 190 and 380 crossflights.As can be seen, the additional crossflights are able to reduce |R gas |/|R| significantly.For the 2-leak scenario, the |R gas |/|R| ratio for 190 and 380 crossflights reduces to 0.072 and 0.026 respectively; and for the 3-leak scenario, the |R gas |/|R| ratio for 190 and 380 crossflights reduces to 0.201 and 0.085 respectively.

B. TESTING SURVEY FRAMEWORK USING REAL DATA 1) EXPERIMENTAL METHODOLOGY
We also evaluated our framework with real UAV measurements.In this experiment, 60.8 g/hr (3 standard cubic feet per hour) of methane gas was released unpressured from a known location in Oxnard, California, close to the ground, in an open field with low grass vegetation.An actual UAV manufactured by DJI, model M600 Pro, was equipped with an RKI Open-Path Laser Spectrometer (OPLS) methane-only sensor with a sensitivity of 10 ppb/s.The OPLS sensor was mounted away from the propellers, as indicated in [12], and [35] and as illustrated in Fig. 10.The UAV was also equipped with a GPS and an anemometer to measure the wind direction and wind speed.A ground anemometer was also used to confirm the wind direction and speed measurements.The UAV was flown at an average height of 4.8m, performing 50-meter crossflights perpendicular to the wind at various distances to the leak.The wind direction pointed towards northeast and its average wind speed was measured to be 4.7m/s.The OPLS collected measurements at a rate of 2 measurements/second and transmitted the measurements using wireless telemetry to a ground base station, which stored the measurements for posterior analysis.Fig. 11 shows an aerial photograph of the field with the leak location and the UAV flight path.The figure also shows the instances of detections by the OPLS sensor.Similarly to our simulation experiments, the UAV in these real experiments flew in a line perpendicular to the wind direction; however, differently from our simulation experiments, the UAV did not fly past the leak.Also, the UAV performed a total of 10 crossflights before repeating the flight pattern 3 times.Two crossflights were discarded for being too short, resulting in a total of 28 crossflights.
The UAV instantaneous concentration measurements were filtered to remove the background methane concentration as follows: the smallest value among the 10 past measurements was subtracted from every measurement.We note that this filtering process allows for the removal of background concentration levels that change during the survey.Using the filtered UAV measurements, we considered C min = 2 • 10 −7 ppb to determine detections and non-detections.
Using the USR models estimated from QUIC simulations and considering the average wind speed of 4.7m/s, we obtained the USR model through linear interpolation, 1398 VOLUME 12, 2024 Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.resulting in the following USR dimensions: fetch start of 1.4m, fetch end of 84.6m, and minor axis of 32.6m.We estimated the inverse-MO to be −0.11/m.Given the approximately constant wind conditions and UAV height, all USRs mapped into the ground had the same orientation and dimensions.
Using such a USR model, USRs were mapped from each UAV location and we computed ij for each subregion r ij .Without using the knowledge of the leak location, we applied the fusion procedure of Section V considering K max = 1 to place each r ij in either R gas , R clear , or R unknown .For this, we estimated the fusion parameter q using the second approach described in Section V-B; i.e., using the 28 crossflights, we estimated q as if the leak were present in each r ij and adopted q to be the highest among all r ij .Such an estimation resulted in q = 0.667.We highlight however that a better estimate could have been obtained if the number of crossflights were higher than 28. 2) DISCUSSION OF RESULTS From Figs. 12 and 13, we note that, although R gas contains the actual leak location, the ratio |R gas |/|R| of the survey map was very high (0.27), meaning that close to one-third of the region was placed in R gas .The high |R gas |/|R| ratio is due to the small number of crossflights and their short length relative to the USR width.Lower |R gas |/|R| ratios could be obtained if the UAV were to perform additional or longer crossflights.Furthermore, if the UAV were to fly past the leak, as in the flight patterns used in the QUIC experiments, then the procedure would be able to restrict R gas and reduce the |R gas |/|R| ratio.Note that, when the UAV flies past the leak, the UAV is able to collect low concentration measurements past the leak and, therefore, add areas past the leak into R clear .
From Fig. 12, we can observe that the highest p-value is in a subregion r ij that does not contain the leak.This could be justified by the low number of crossflights or by errors in the estimation of the inverse-MO, wind speed, or wind direction.

VIII. SUMMARY, CONCLUSION, AND FUTURE STUDIES
Focusing on the approach that uses UAVs flying inside the gas plume to assist in the survey of a region for gas leaks, we proposed a framework using Upwind Survey Regions (USRs) to translate UAV measurements into a survey map.Our framework has two components: USR mapping and USR fusion.After formally defining USRs based on long-term concentration averages, we discussed how USRs can be used to map UAV measurements collected at each location p into a region USR(p) on the ground.Given that UAVs collect instantaneous concentration measurements, we proposed a statistical procedure to fuse USRs and produce a survey map that partitions the surveyed region in 3 areas: R gas (area likely to contain a gas leak), R clear (area clear of gas leaks), and R unknown (area from which not enough information is available).Since the USR dimensions vary with environmental and flight conditions, we discussed how a database of USR models can be built using the QUIC simulation tool.We further discussed ways to evaluate the performance of the survey and we showed through both simulation and real UAV measurements that the framework is able to produce satisfactory survey maps.
We note that surveying a region for gas leaks is more than just finding the location of the leak.Gas companies need tools to also determine that an area is clear of leaks.Furthermore, methods to estimate the location of a gas leak assume that there is just one leak in the region and their performance is unclear when multiple or no gas leaks are present.Our framework works when the region has zero, one, or multiple gas leaks.Our framework produces a survey map that gas companies can use to document periodic surveys and show to regulatory agencies that it is performing steps  to reduce fugitive gas leaks, as suggested in recent regulatory efforts [75].
A first conclusion of this paper is that surveying a region with UAVs and USRs needs statistical tools.Gas dispersion is a highly random process: instead of producing a homogeneous plume, a gas leak produces gas filaments that are subject to turbulence, producing random paths.This means that UAVs in situ instantaneous concentration measurements are random.Since USRs are defined based on long-term concentration averages, a surveyor should not infer that the USR associated with a single high concentration measurement contains a gas leak and the surveyor should not infer that the USR associated with a single low concentration measurement is clear of leaks.In fact, a single location in the ground may be covered by multiple USRs, some associated with high concentration measurements and some associated with low concentration measurements.The surveyor should see each UAV concentration measurement as providing evidence for the presence or absence of a gas leak, and should fuse this evidence using a statistical procedure.In this paper, we proposed and justified a possible statistical procedure for this fusion.
A second conclusion of this paper is that QUIC is a valuable simulation tool to design procedures for gas leak surveys.Simulations are needed because collecting measurements from real gas leaks is difficult, expensive, and time-consuming.Simulations can be used to estimate the performance of a survey procedure by testing it under various random conditions.In addition to providing guidance to the design of survey procedures, QUIC can also be used to build the database of USR models.
We highlight the following as future avenues for research: • Given that real UAV experiments are expensive and time-consuming, we evaluated the performance of our procedure with QUIC simulations and a single set of real UAV measurements.We plan to perform further real UAV measurements, with the UAV flying past the leak, to better evaluate the performance of the procedure.
• We evaluated our framework considering static wind conditions; i.e., considering that the UAV height, the wind speed, and the wind direction were constant during the survey.As mentioned in Section IV, our approach allows for USRs of different sizes and orientation to be mapped onto the ground.Future studies should evaluate the framework performance in dynamic wind conditions.
• Our fusion procedure relies on the probability that the UAV collects a high concentration measurement during the crossflight (parameter q defined in Section V-A); and the procedure that we studied in this paper assumes an average value for q, which enables the use of the Binomial model for reaching decisions.We plan to study how the performance of the survey maps would improve if we allowed for values of q that vary with the distance between the crossflight and the subregion being evaluated.
• In this initial study, we considered a fixed sweeping flight pattern that does not attempt to seek the gas leak.This is important because collecting measurements from outside the gas plume is necessary to reliably determine R clear ; however, it is also possible to consider the fusion of USRs generated from multiple sweeping flights: after the first survey map is generated from a first sweeping flight, a second sweeping flight can be designed and flown to refine the first survey map.USRs mapped from the second flight would then be fused with the first set of USRs, generating a second survey map.This process is repeated until the final survey map has a small enough R gas area for further investigation.
• We plan to continue using QUIC to evaluate our framework in areas containing buildings and vegetation in various configurations.
Lastly, we highlight that our framework can also be used to detect other gases.For this, the surveyor needs to equip the UAV with an OPLS sensor tuned to detect the target gas and needs to generate a database of USR models specific to the target gas.
x and y = p y − p ′ y ; and using (10), it follows that (x + x , y + y ) ∈ USR(p x , p y , p z ).This process is illustrated in Fig. 15: if the surveyed region is homogeneous and a leak at (x, y) produces C avg (p 1 ) > C min in the right display of Fig. 15, then we can conclude that a leak at (x 1 , y 1 ) in the left display of Fig. 15 also produces C avg (p) > C min ; i.e., (x 1 , y 1 ) ∈ USR(p) because the relative displacement between p and (x 1 , y 1 ) is the same as the relative displacement between p 1 and (x, y).Thus, by observing {C avg (p 2 ), C avg (p 3 ), . ..} for various points from a single QUIC simulation with a leak at (x, y), we can determine USR(p).
The homogeneous condition (10) also implies that, maintaining all meteorological parameters constant, the dimensions of the USR(p) remain constant for all p with a common height p z , meaning that it is enough to estimate the dimensions of a common USR for each height p z .
In this paper, we consider the survey of open fields; i.e., regions without buildings or large vegetation.In this case, the homogeneous condition applies; and the following process is used to estimate the USR.Fig. 14 illustrates the process.
1) Fix the meteorological conditions: wind direction, wind speed, wind profile, and inverse MO parameter.2) Place a leak with strength s min in a known position (x * , y * , z * = 0).3) Execute the QUIC simulation for a large enough time period until it reaches steady-state.4) From the last 3D concentration matrix A, smooth the average concentration at each height.More precisely, for each discrete height z in the 3D matrix, extract from A the slice at height z, resulting in a 2D matrix A z ; and convolve A z with a 2D Gaussian kernel.The purpose of this step is to smooth the average concentration matrix.Let Ĉavg (p) represent the 3D concentration field after smoothing it at the various heights.Fig. 16 illustrates the 2D concentration matrix A z and the 2D average concentration matrix after convolving it with a 2D Gaussian kernel.5) Estimate the common USR for each of the discrete heights p z as follows: a) Find all points {p l } L l=1 with height p z and Ĉavg (p) > C min .In the example in Fig. 16, such points are located inside the dashed level curve in the right panel when considering C min = 6.6 • 10 −6 g/m 3 .b) For each p l , compute the distance d l := ∥(p l,x , p l,y ) − (x * , y * )∥ and the angle θ l between the wind direction and the line containing (p l,x , p l,y ) and (x * , y * ).Note that d l and θ l define a point that would be inside USR(p l ) if the UAV were to collect a long-term measurement at p l .Note that, because of the homogeneous condition (10), the same {d l , θ l } L l=1 would be obtained if we were to observe the concentration at a single point p and perform multiple QUIC simulations with leaks at various locations.c) Considering p = (0, 0, p z ), fit an ellipsoid containing most of the points located {d l , θ l } L l=1 from the origin.More precisely, an enclosing rectangle with major sides aligned with the wind direction is defined around {d l , θ l } L l=1 by considering a high percentile of points to enclose.The ellipsoid is then defined from the start, end, and width of the rectangle.We refer to the distances to the start and end of the ellipsoid as the fetch start (d start ) and fetch end (d end ) and we refer to the minor axis of the ellipsoid as the width (w) of the ellipsoidal USR, as illustrated in Fig. 14.Such an ellipsoid represents the USR for any p with the same height p z .6) Repeat the above procedure for various sets of meteorological conditions.At the end of the procedure, the surveyor has a database of USR models, one for each set of meteorological conditions and UAV heights.

A. EXAMPLE
To build the database of USR models, we configured QUIC to simulate gas fields in an open field with 200m of width (x-dimension), 1000m of depth (y-dimension), and 20m of height (z-dimension).A gas leak releasing unpressured gas at a constant rate of 100 g/hr was placed at a known location (x * = 100, y * = 900, z * = 0) with a southbound wind spreading the gas particles in the 3D region.
We ran simulations at 3 different wind speeds: 2m/s, 4m/s, and 6m/s; and 4 different inverse MO conditions: −0.02/m, −0.05/m, −0.08/m, and −0.11/m, resulting in 12 different simulations.In each simulation, the wind speed and direction were constant throughout the survey.In each simulation, we configured QUIC to produce gas fields with 6 million particles, producing gas fields that approximate the long-term average concentration at each point in the space.Before estimating the USR model, we ran the simulation for 800 seconds to ensure that it reached the steady-state condition.
For each simulation, we used the procedure described in Appendix B to obtain the ellipsoidal USR dimensions at various UAV heights (p z ∈ {3, 4, 5, 6, 7, 8}).We note that all QUIC simulations used the logarithmic wind profile, which means that the wind speed varies with the height.QUIC was configured so that the nominal wind speed is set at height p z = 5.5m.Recall that the ellipsoidal USR is defined by its fetch start (d start ), fetch end (d end ), and width (w).Table 1 illustrates the USR dimensions obtained for a few of the heights and inverse MO -0.11/m; and Table 2 illustrates the USR dimensions obtained for the various inverse MO values for p z = 5m.
From Table 1, it is possible to note that the USR fetch start (d start ) tends to increase and the USR fetch end (d end ) and the width (w) tend to decrease with the wind speed.Since the topology simulated in QUIC was the open field, this behavior was expected from the Gaussian plume model [66].Since the wind pushes the particles before they can gain height, a UAV at a fixed altitude will start sensing the plume at higher fetch starts when the wind speed increases.Also, since the plume concentration is inversely proportional to the wind speed, higher wind speeds for a fixed minimum concentration threshold (C min ) mean shorter fetch ends (d end ).
From Table 2, it is possible to note that the USR fetch start (d start ), the USR fetch end (d end ), and the width (w) all tend to decrease with the inverse MO.This behavior was also expected because lower values for inverse MO mean more unstable conditions [66], causing particles to disperse more easily, and therefore, reducing the concentration.

APPENDIX C SIMULATING INSTANTANEOUS CONCENTRATION MEASUREMENTS
The QUIC simulation tool is able to simulate average concentrations (C avg (p)) for each position p in the 3D space; however, it does not simulate the coupling between the filamentous nature of the gas plume and atmospheric turbulence that produces the intermittent instantaneous concentration C(p, t) measured by the UAV.
To simulate C(p, t), we follow an approach similar to [62], [72], [73], and [74] in that we model the number of filaments crossing the UAV path using a Poisson process with average proportional to C avg (p).More precisely, after using QUIC to produce C avg (p) for each position p in the 3D space, we simulate the UAV flying in the 3D average concentration field.For each time t mn , the UAV is at position p mn and we collect the long-term average concentration value C avg (p mn ).To obtain the instantaneous concentration C(p mn , t mn ), we model the number of filaments crossing the position p mn at time t mn as a Poisson-distributed random variable F mn with parameter C avg (p mn )/C filament , where C filament is the expected concentration in a filament of the turbulent flow.For each position p mn , we draw a random realization f mn of the random variable F mn ; and use f mn to produce a sample of the instantaneous concentration collected by the UAV: C(p mn , t mn ) = f mn • C filament .In our simulations, we considered C filament = 2 • 10 −4 g/m 3 .Considering such a C filament and the constant emission rate of 100 g/hr, the resulting Poisson parameters are always less than 1, producing a sparse sequence of UAV instantaneous concentration measurements, similar to those observed in real UAV flights, as illustrated in Fig. 17

FIGURE 1 .
FIGURE 1. Illustration of a UAV carrying a methane sensor to survey a region for gas leaks.The UAV flies downwind of a possible emission source following a sweeping flight pattern.The dashed lines respectively illustrate possible trajectories of methane (red) and other molecules (green) as they travel in the atmosphere.

FIGURE 2 .
FIGURE 2. Overview of the proposed USR-based framework to perform survey for gas leaks.The framework has two components: USR mapping and USR fusion.After selecting a USR from a database of USR models, we map the UAV concentration obtained at each location into a region on the ground.We fuse USRs covering each subregion to produce the final survey map.The database of USR models may be built as described in Appendix B.

FIGURE 3 .
FIGURE 3. Illustration of USR mapping.Each instantaneous UAV concentration measurement is mapped into a region on the ground.Each USR is associated with a high concentration (detection) or low concentration (non-detection) measurement.

FIGURE 4 .
FIGURE 4. Illustration of five USR(p mn ) mapped from five UAV measurements collected at positions {p m1 , . . ., p m5 } during crossflight m.USRs encircled with red and green lines respectively represent USRs associated with detections and non-detections.Note that USRs overlap and a subregion r ij may be covered by multiple USRs.If during the crossflight m any of the USRs covering a subregion r ij is associated with a detection, then e (m) ij = 1.

FIGURE 5 .
FIGURE 5. Illustration of instantaneous UAV measurements providing incorrect evidence regarding the presence or absence of leaks.Due to the random path taken by gas filaments, an instantaneous UAV detection at p m1 may provide incorrect evidence that the region covered by USR(p m1) contains a leak.Likewise, although p m2 is inside the average gas plume of a leak, it might happen that the instantaneous UAV measurement is below C min , causing a non-detection and providing incorrect evidence that the region covered by USR(p m2 ij = 1 if at least one of the USR(p mn ) of the crossflight m that covers r ij is associated with a detection (see Fig. 4 for an illustration); otherwise, e (m) ij = 0: e (m) ij = 1 ⇔ N n=1 1 r ij ∈ USR(p mn ) and C mn > C min > 0; i.e., e (m) ij = 1 represents evidence that r ij should be placed in R gas ; and fuse all e R clear .Since TP ≪ TN , Acc ≈ TN /(TP + TN + FP + FN ) = 1 − (TP + FN + FP)/(TP + TN + FP + FN ).And since FN ≪ TN , Acc ≈ 1 − (TP + FP)/(TP + TN + FP + FN ).Lastly, since |R gas | represents the number of regions placed in R gas , which equals TP + FP, Acc ≈ 1 − |R gas |/|R|.

FIGURE 6 .
FIGURE 6. Number of crossflights with high concentration USRs ( M ij m=1 e (m) ij ) and the resulting p-value P ij for each subregion r ij .The actual leak location (not known by the procedure) is indicated with a circle.The wind direction is from top to bottom.The UAV flight path is indicated with cyan lines.

FIGURE 7 .
FIGURE 7. Final survey map.R clear is shown in green, representing areas for which the evidence of a gas leak is small.R gas is shown in red, representing areas that could not be cleared.R unknown is shown in white.The actual leak location (not known by the procedure) is indicated with a circle.The wind direction is from top to bottom.The UAV flight path is indicated with cyan lines.

FIGURE 8 .
FIGURE 8. Final survey maps for the experiments that resulted in the highest (left) and lowest (right) |R gas |/|R|.R clear is shown in green, R gas is shown in red,, and R unknown is shown in white.The actual leak location (not known by the procedure) is indicated with a circle.The wind direction is from top to bottom.The UAV flight path is indicated with cyan lines.

FIGURE 9 .
FIGURE 9. Final survey maps when UAV performs 100, 190, and 380 crossflights when the region has 2 or 3 leaks.R clear is shown in green, R gas is shown in red, and R unknown is shown in white.The actual leak location (not known by the procedure) is indicated with a circle.The wind direction is from top to bottom.The UAV flight path is indicated with cyan lines.

FIGURE 10 .
FIGURE 10.Photograph of UAV and OPLS sensor used in experiments.

FIGURE 11 .
FIGURE 11.Aerial photograph of region where experiment was performed with a real UAV.White lines indicate the flight path of the UAV.The X mark on the bottom left of the photograph indicates the gas leak location.The wind direction points towards northeast (top right of the photograph).Red lines indicate instances of detections by the OPLS sensor.
ij , the corresponding p-values P ij for each subregion and the final survey map are shown in Figs. 12 and 13.

FIGURE 12 .
FIGURE 12. Number of crossflights with high concentration USRs ( M ij m=1 e (m) ij ) and the resulting p-value P ij for each subregion r ij .The actual leak location (not known by the procedure) is indicated with a circle.The wind direction points to the top right corner of the figure.The UAV flight path is indicated with cyan lines.

FIGURE 13 .
FIGURE 13.Final survey map when considering real data collected by a UAV.R clear is shown in green, R gas is shown in red, and R unknown is shown in white.The actual leak location (not known by the procedure) is indicated with a circle.The wind direction points to the top right corner of the figure.The UAV flight path is indicated with cyan lines.

FIGURE 14 .
FIGURE 14. Overview of method to build the database of USR models from QUIC simulations.

FIGURE 15 .
FIGURE 15.Two approaches to determine USR(p): (left) using multiple simulations, each with a leak at different positions; (right) using a single simulation with the leak at a single position, applicable when the surveyed region is homogeneous.

FIGURE 16 .
FIGURE 16. (left) Illustration of 2D concentration matrix A z obtained from a QUIC simulation with leak at fixed location, wind direction from north to south, wind speed of 4m/s, inverse MO −0.05/m, and height z = 5m.(right) Illustration of 2D average concentration matrix, obtained after convolving A z with a Gaussian kernel of standard deviation 12m.The dashed line represents the level curve at C min = 6.6 • 10 −6 g/m 3 .Each of the points inside such a curve produces a distance d k and angle θ k from the gas leak.After referring the various d k and θ k to the origin, we define the ellipsoidal USR.

TABLE 1 .
Database of USR models: Fetch start (d start ), fetch end (d end ), and width (w ) of the USR ellipsoid under various open field conditions and inverse MO −0.11/m.

FIGURE 17 .
FIGURE 17. (left) Illustration of instantaneous concentration measurements obtained in a real UAV flight.(right) Illustration of instantaneous concentration measurements obtained in simulated QUIC flight with Poisson-distributed filament encounters.In both cases, the UAV was flying in a field with a gas leak present.It is possible to observe that the sequence of instantaneous concentration measurements in the simulated UAV flight is sparse, similar to the real UAV flight.

TABLE 2 .
Database of USR models: Fetch start (d start ), fetch end (d end ), and width (w ) of the USR ellipsoid under various inverse MOs for UAV height of 5m. .