Landfill methane source level emission automatic identification and adjacent source correction quantification method and system

Abstract

The present application provides landfill methane source level emission automatic identification and adjacent source correction quantification method and system, concentration and wind field data are collected through underway observation, emission sources are automatically identified through gridding, morphological processing and connected domain clustering, baseline emission rate is obtained by integration with multiple cross sections, and adjacent source interference contribution is deducted based on inverse distance weighting, combined with Monte Carlo uncertainty propagation of five parameters of wind speed, wind direction, concentration, background and cross section distance and four internal diagnostic self-checking, automatic identification of source level emission, adjacent source correction, uncertainty quantification and confidence discrimination in multi-source scene of landfill are realized, subjective bias of manual source identification is eliminated, systematic overestimation caused by plume overlap is inhibited, reliable confidence interval is provided, and implicit quantitative failure results are automatically filtered.

Description

Technical Field

[0001] This invention belongs to the field of environmental monitoring technology and relates to a method and system for quantifying methane source-level emissions from landfills. Background Technology

[0002] For methane emissions from landfills, commonly used top-down emission quantification methods include: mass balance methods (including columnar / box mass balance, upwind / downwind cross-sectional pair integration, and other implementations), tracer gas correlation methods, eddy covariance methods, and satellite inversion methods. Among these, the upwind / downwind cross-sectional pair mass balance implementation (hereinafter referred to as "cross-sectional method," also known as transect-based mass balance in English literature) is widely used in monitoring large-scale industrial emission sources due to its simple equipment, high spatial resolution, and ability to quantify different sub-regions within the site. Its basic principle is to arrange concentration sampling cross-sections perpendicular to the wind direction on both the upwind and downwind sides of the emission source. The emission rate is obtained by multiplying the downwind concentration increase relative to the upwind concentration by the normal wind component and integrating along the vertical direction.

[0003] However, existing cross-sectional methods have the following technical problems when applied to multi-source, dense emission scenarios such as municipal solid waste landfills: First, source identification relies on manual experience, typically requiring operators to manually define emission source boundaries based on concentration contour lines. This approach suffers from strong subjectivity, lack of reproducibility across different runs, and a tendency to overlook small sources, making it difficult to meet the requirements of operational and automated systems. Second, when multiple sources are adjacent, downwind plumes overlap. Multiple sources, such as enclosed historical landfill units, active work areas, and leachate treatment zones, are often distributed adjacently on a scale of tens to hundreds of meters. The downwind cross-section of the target source will simultaneously capture contributions from neighboring sources, resulting in a systematic overestimation of single-source quantification results. Existing methods often treat the entire site as a single source. First, the quantitative analysis lacks a mechanism to deduct inter-source interference. Second, the uncertainty assessment is simplistic; existing methods often only perform sensitivity analysis on a single parameter such as wind speed or concentration, lacking an uncertainty propagation framework for joint disturbances of multiple parameters such as wind speed, wind direction, concentration, background, and cross-sectional location, resulting in unreliable confidence intervals for single-source emission rate values. Third, the confidence level of the quantitative results lacks internal diagnosis; when the cross-sectional spacing in the upwind and downwind directions is improperly selected or local meteorological conditions deviate from steady state, the cross-sectional method may give seemingly reasonable but actually distorted positive results. Existing methods lack a self-checking mechanism for results based on internal indicators such as cross-sectional spacing sensitivity curves, flux decomposition, and connected domain morphology, thus failing to identify such hidden failures. Summary of the Invention

[0004] To address the problems of existing cross-sectional methods described in the background art, which rely on human experience for source identification, overlap of downwind plumes in nearby multi-source, dense emission scenarios such as municipal solid waste landfills, single uncertainty assessment, and lack of internal diagnostics for the confidence level of quantification results, this invention provides an automatic identification and neighboring source correction quantification method and system for methane source-level emissions from landfills. This system can achieve automatic identification of methane source-level emissions from landfills, correction of neighboring source interference, multi-parameter uncertainty propagation, and internal diagnostic self-checking.

[0005] The method of the present invention includes the following steps: S1. Data Acquisition and Matching: The methane concentration, effective path length, coordinates, meteorological parameters, and wind speed and direction of the measuring points are collected synchronously by mobile surveying. The horizontal wind component calculated from the wind speed and direction and the background concentration are spatially interpolated to the concentration measuring points to form a point-by-point matching dataset. S2. Automatic source identification: The dataset is projected onto the UTM coordinate system and rasterized. The threshold is determined according to the statistical characteristics of the positive concentration raster to generate a high-value mask. After morphological processing, connected component clustering and minimum raster number filtering, candidate source connected components and attributes are obtained. S3. Cross-section layout and baseline quantification: Using candidate sources as target sources, the average wind direction and normal direction are determined based on the surrounding coordinated wind field. Multiple sets of upwind and downwind cross-sections are laid out within a preset spacing range. The cross-section concentration is reconstructed by interpolation. The median effective flux is taken as the baseline emission rate by integrating the normal wind component. S4. Neighboring source inverse distance weighted correction: Screen neighboring sources within the downwind half-space of the target source and whose distance does not exceed the set radius, estimate their concentration contribution according to the square weight of inverse distance and deduct it from the downwind enhancement amount of the target source, and re-integrate to obtain the corrected emission rate. S5. Multi-parameter Monte Carlo uncertainty propagation: Perturbations are applied to wind speed, wind direction, cross-sectional reconstructed concentration, background concentration and cross-sectional spacing. Repeat S3 and S4 to obtain quantiles and confidence intervals. S6. Internal diagnostic self-check: The confidence level is determined by the sensitivity of cross-sectional spacing, the sign of upwind and downwind flux, the anomaly of the area of ​​the connected domain, and the consistency of the magnitude of cross-flights. Low confidence results are not included in the site-level summation.

[0006] Furthermore, in step S1, an unmanned aerial vehicle (UAV) platform equipped with a laser methane telemetry instrument and an ultrasonic anemometer is used to conduct near-ground navigation above the target landfill along a preset route, simultaneously recording the following second-by-second data: longitude. lon ,latitude lat ,altitude alt methane volume fraction conc ppm and path integral concentration conc ppm ·m, path length h Wind speed ws ,wind direction wd ,temperature T air pressure P ; Corresponding to the wind observation section, its point-by-point horizontal wind component u , v Calculate according to formula (1): (1), in, u , v These represent the easterly wind component and the northerly wind component, respectively, in m / s; ws Wind speed, in m / s; wd Wind direction, in degrees, using a meteorological coordinate system, with clockwise from true north as positive; π / 180 is the angle-to-radian coefficient; Using a two-dimensional linear interpolation method, the concentration measurement points ( lon , lat Using ) as the query point, match and generate point-by-point. u point , v point For concentration points located outside the convex hull of the wind tunnel, nearest neighbor interpolation is used for extrapolation and filling, followed by... u point , v point The point-by-point wind speed is calculated by reverse calculation. point ws With wind direction point wd Simultaneously, background concentration is also statistically analyzed; the output is a concentration dataset with point-by-point collaborative wind field and background concentration.

[0007] Furthermore, the specific steps of S2 include: S21. Coordinate Projection: Project the latitude and longitude of the WGS84 data set into the local UTM coordinate system to obtain the planar coordinates (x, y). S22. Rasterization: The navigation area is divided into regular square grids according to the resolution dx. The dx value is set according to actual needs. The concentration of all measuring points in each grid is calculated by arithmetic mean to obtain the grid concentration field. C ( i , j Empty rasters are marked with missing values; S23. Background and Threshold Determination: The 5th percentile of all positive concentration grid cells is taken as the background value. C bg The positive concentration grid refers to a grid with a positive concentration value, and the standard deviation of the concentration of all positive concentration grids is calculated. σ bg ; Calculate the detection threshold according to formula (2)T : (2), in, T The concentration detection threshold is expressed in ppm·m. C bg The background concentration is the 5th percentile of all positive concentration grids, expressed in ppm·m. σ bg The standard deviation of the background concentration grid is given in ppm·m. Will C ( i , j ) ≥ T The grid is denoted as a high-value binarization mask. M ; S24. Morphological operations: High-value binarization mask M Three morphological closure operations are performed sequentially to eliminate voids caused by the density of the travel paths within a single source; then one morphological opening operation is performed to remove isolated high-value burrs; the structuring element uses a 3×3 fully connected kernel; S25. Connected Component Clustering and Filtering: Clustering is performed on the morphologically processed mask using the 8-neighbor connected component algorithm, and cells with fewer than 8 connected components are removed. N min The tiny connected regions, N min Based on actual needs, the set of candidate connected components is obtained. D k}; for each candidate source connected component D k Calculate its centroid coordinates ( x k, y k ),area A k (Number of grids × dx) 2 Peak concentration C peak,k With mean concentration C mean,k The above results are then used as the source identification output data.

[0008] Furthermore, S3 identifies each candidate source connected component in the source identification output data. D k For the target source, the emission rate is quantified according to the following steps: S31. Determine the average wind direction: using candidate source connected components D k Measurement points within the surrounding preset buffer zone upoint , v point By performing a vector average, the average wind direction can be obtained. wd k and let the normal unit vector n along wd k Pointing downwind; S32. Layout of upwind and downwind cross-sections: with D k centroid coordinates ( x k, y k For reference, along - n with + n The directions are translated by distance respectively d Determine the center of the cross-section in the upwind and downwind directions; the cross-section direction is... n Vertically, the length covers the projection of the connected domain in the vertical wind direction and extends to both sides by P, where P is set according to actual needs, and the cross-sectional spacing is... d Within a preset range, values ​​are selected one by one according to the step size to form... d Scan sequence; S33. Concentration Field Reconstruction: Sampling points are generated at step sizes on each cross-section, and the concentration is reconstructed using two-dimensional linear interpolation. C (s); Samples outside the convex hull are compensated by nearest neighbor interpolation and multiplication by the attenuation factor. S34. Emission rate calculation: The normal wind component is taken at the center point of the cross-section in the upwind and downwind directions. u point , v point and n dot product average U ⊥, d The methane emission rate at the cross-sectional spacing d Q d Calculate according to formula (3): (3), in h (s) represents the effective path length at the sampling point, which is the local mean of the measured path length during navigation. ρ corr The temperature and pressure standard condition correction factor is calculated according to formula (4), converting the volume fraction ppm to the mass concentration under standard conditions: (4), in, P This refers to the on-site air pressure. P 0 Standard atmospheric pressure; T0 = 273.15 K; T The ambient temperature; S35, Baseline Q base,k Selection: Select the effective range within the preset spacing. Q d The median emission rate was used as the baseline emission rate. Q base,k and record Q d Statistics and follow-up d The change curve is for internal diagnostic use.

[0009] Furthermore, in S4, when multiple sources are distributed adjacently within the site, the downwind cross section of the target source will simultaneously capture the contributions of neighboring sources. The following method is used to subtract these neighboring source contributions: S41. Neighbor source filtering: For target sources D k traverse the remaining recognition sources D m Calculate the Euclidean distance between sources r km Only retain r km ≤ R max The source is used as a candidate neighbor; S42, Directional Gating: Further requirements for neighboring sources D m Located at the target source D k Within the downwind half-space, that is, satisfying ( x m - x k, y m - y k ) · n k > 0; if this condition is not met, the contribution from neighboring sources will not be counted. S43, Inverse Distance Weighted Weight: Neighborhood D m The concentration contribution weight of the downwind cross section of the target source is calculated according to formula (5): (5), in, w m Neighboring source D m For the target source D k The inverse distance-weighted weight of the downwind cross-sectional concentration contribution is dimensionless. rkm For target source D k and neighboring sources D m Euclidean distance in the UTM coordinate system, in meters; w min The minimum weight lower bound is dimensionless and used to avoid long-distance neighbor factors. r km The square decays too quickly and loses its physical meaning; S44. Contribution Estimation: Neighboring Sources D m The contribution to the concentration in the downwind cross section of the target source is estimated by multiplying the downwind concentration increase by the upwind concentration increase in its own cross section by a weight: (6), Where, Δ C m ( s ) are neighboring sources D m The cross section downwind of the target source in coordinates s Concentration contribution at a given location, in ppm·m; Δ C m_self ( s ) are neighboring sources D m The concentration increase from downwind to upwind obtained in its own quantification, in ppm·m; in ΔC m_self (s) are neighboring sources D m The downwind increase minus the upwind increase obtained from its own quantification; if ΔC m_self (s) If not reconstructed within the same coordinate range of the target cross-section, then by ΔC m_self Substitute the mean along its normal direction; S45. Corrected emission rate: Sum the contributions of all candidate neighboring sources and subtract the downwind enhancement from the target source, then re-integrate to obtain the corrected emission rate of the neighboring sources, and calculate according to formulas (7) and (8): (7), (8), Where, Δ C corr ( s The target source's downwind cross section after neighbor source correction is shown in coordinates. sThe concentration increase at a given location, expressed in ppm·m; the first term on the right-hand side of the equation represents the uncorrected downwind minus upwind increase from the target source; summation symbol Σ. m Indicates distance from the target source k The central Euclidean distance does not exceed R max And all candidate neighboring sources located in its downwind half-space D m Summation; Q corr,k For target source D k Emission rates corrected for adjacent sources, in kg / h; U ⊥ The normal wind component at the center of the cross-section downwind of the target source, in m / s; h ( s ) represents the effective path length at the sampling point, in meters (m). ρ corr This is a dimensionless correction factor for temperature and pressure under standard conditions. ds The length of the infinitesimal element along the cross-section is expressed in meters (m). S46. Correction Amplitude Report: Define the correction amplitude δ k Calculate according to formula (9): (9), in, δ k For target source D k The neighbor source correction magnitude is dimensionless and can be reported as a percentage. Q corr,k Emission rate after correction, unit: kg / h; Q base,k Baseline emission rates before correction, in kg / h; the more significant the contribution from neighboring sources... δ k The further away from 0; δ k A value close to 0 indicates that interference from neighboring sources is negligible.

[0010] Furthermore, in step S5, to obtain the confidence interval for the emission rate of each source, multi-parameter Monte Carlo sampling is embedded based on steps S3 and S4. The specific steps are as follows: S51. Set the disturbance parameters and default amplitude: Wind speed: Obey N (1, 0.10 2 The relative perturbation of ), that is, each iteration will U Multiply the whole by a random coefficient that follows this distribution; Wind direction: Obey N (0, 7.5 2 An absolute perturbation of 0° with respect to the normal direction n Rotate; Reconstructed concentrations at each cross-sectional sampling point: follow the rules N (1, 0.10 2 The relative perturbation of ) affects the reconstructed concentration of each cross-sectional sampling point obtained in the S3 cross-sectional layout and baseline quantization. C (s); Background concentration: obeys N (1, 0.05 2 The relative perturbation of ) acts on the background concentration determined in S2 automatic source identification. C bg ; Section spacing d :obey N (1, 0.05 2 The relative disturbance of ) acts on the cross-sectional spacing of the upwind and downwind cross-sections in the S3 cross-section layout and baseline quantization; S52. Fixed random seed: A pseudo-random sequence with seed = 42 is used; S53, Number of iterations: Execution N mc = 500 to 1000 iterations; in each iteration, re-execute steps S3 and S4 with the above perturbation to obtain a set of { Q corr,k}_ i , i = 1… N mc ; S54. Output statistics: P5, P50, and P95 quantiles and 90% confidence interval width for the output emission rate of each source.

[0011] Furthermore, in S6, the four indicators for judging the credibility of the quantification results of each source emission rate are as follows: First, the calculation of the coefficient of variation of the cross-sectional spacing sensitivity curve: (The text abruptly ends here, so the translation stops as well.) d Calculate according to formula (10) within the preset value range. Q d coefficient of variation CV d ,like CV d A value > 0.50 is marked as low confidence. (10), in, CV d The target source emission rate varies with the cross-sectional spacing. d The coefficient of variation, dimensionless; σ( Q d ) for d All valid values ​​within the preset range Q d Standard deviation, in kg / h; mean ( Q d Within this range Q d The arithmetic mean of the values, in kg / h; |·| indicates taking the absolute value; CV d Larger means Q Spacing of sections d The stronger the selected dependence, the more it suggests that the plume pattern deviates from the steady-state advection assumption; Second, check the signs of the upwind and downwind flux components: separately count the signed fluxes in the upwind and downwind cross sections. F up , F down And compare their absolute values; if |F up | > |F down | or F up and F down If the numbers are the same and the magnitudes are close, it indicates that the transmission direction assumption has failed, meaning that the upwind direction has been contaminated by a non-target source or the wind direction has reversed, and it is marked as low confidence. Third, diagnosis of abnormal connected area: When the connected area of ​​the same source in different flights at the same site expands by more than 2 times and corresponds to the decrease in the average wind speed of the site, it indicates that the plume accumulates nearby under low wind speed and wind direction change conditions, and the result of this source in this flight is marked as low confidence. Fourth, cross-flight consistency diagnosis: If the same source can be quantified in more than two flights, and the deviation of the results in each flight exceeds 5 times, it indicates that at least one of them failed, and further diversion is carried out according to internal diagnostic indicators.

[0012] This invention also provides an automatic identification and neighboring source correction quantification system for methane source-level emissions from landfills, used to perform the methods described above, characterized in that it includes: Data acquisition and matching module: used to collect methane concentration, effective path length, coordinates, meteorological parameters and wind speed and direction at the measuring points, and interpolate to form a point-by-point matching dataset; Automatic source identification module: connected to the data acquisition and matching module, used to project and rasterize the dataset, generate a high-value mask according to the statistical characteristics of positive concentration raster, and obtain candidate source connected components and attributes through morphological processing and connected component filtering; Cross-section layout and baseline quantization module: connected to the automatic source identification module, used to lay upwind and downwind cross-section pairs based on the cooperative wind field around the candidate source, interpolate and reconstruct the cross-section concentration and integrate to obtain the baseline emission rate; Neighboring source inverse distance weighted correction module: connected to the cross-sectional layout and baseline quantification module, used to screen downwind neighboring sources, estimate their contribution according to inverse distance weight, and subtract their contribution to obtain the corrected emission rate; Multi-parameter Monte Carlo uncertainty propagation module: connected to the cross-section layout and baseline quantization module and the neighboring source inverse distance weighted correction module, used to perturb wind speed, wind direction, cross-section reconstruction concentration, background concentration and cross-section spacing and output quantiles and confidence intervals; Internal diagnostic self-test module: connected to the cross-section layout and baseline quantification module, the neighboring source inverse distance weighted correction module and the multi-parameter Monte Carlo uncertainty propagation module, used to determine the confidence level based on cross-section spacing sensitivity, flux sign, connected domain area anomaly and cross-flight consistency, and mark those that fail as low confidence and not included in the site-level emission summation.

[0013] The present invention also provides an electronic device, comprising: a memory and a processor, wherein the memory and the processor are communicatively connected to each other, the memory stores computer instructions, and the processor executes the computer instructions to realize the automatic identification and neighboring source correction quantification method for landfill methane source-level emissions as described above.

[0014] The present invention also provides a computer-readable storage medium storing a computer program that, when executed by a processor, implements the method for automatic identification and neighboring source correction quantification of methane source-level emissions from landfills as described above.

[0015] Compared with the prior art, the present invention has the following advantages: (1) Fully automated source identification, elimination of manual dependence and guarantee of reproducibility across flights: Traditional cross-sectional methods rely on operators to manually define emission source boundaries based on concentration contour lines, which has problems such as strong subjectivity, non-reproducibility of results across flights, and easy omission of small-scale emission sources. This invention adopts a three-stage algorithm pipeline of rasterization-morphological processing-connected domain clustering, which automatically determines the detection threshold and generates candidate source connected domains based on the statistical characteristics of positive concentration raster, completely replacing the manual selection process, ensuring that the source identification results under different flights at the same site are objective and consistent, and meeting the automation requirements of continuous business operation; (2) Quantitative deduction of neighboring source interference and suppression of single-source systematic overestimation caused by plume overlap: In the scenario of landfills with multiple densely distributed sources, each independent emission source is distributed in close proximity within a scale of tens to hundreds of meters. The downwind cross section of the target source inevitably captures the contribution of neighboring sources, resulting in a systematic overestimation of the single-source emission rate when the traditional method treats the site as a single source for total quantification. Based on the principle of inverse distance weighting, this invention quantitatively estimates and deducts the interference contribution of each neighboring source to the downwind concentration enhancement of the target source by screening neighboring sources in the downwind half-space of the target source and assigning weights according to the inverse square distance. Thus, under the premise of maintaining the conservation of the total emission of the site, the source-level emission intensity of each independent emission source is obtained as truly comparable. (3) Five-parameter joint uncertainty propagation and provision of emission rate confidence interval: Existing methods often only perform sensitivity analysis on a single parameter such as wind speed or concentration, which cannot reflect the comprehensive impact of the coordinated variation of multiple parameters such as wind speed, wind direction, concentration, background and cross-sectional location on the quantification results, resulting in unreliable confidence intervals for single-source emission rates; This invention establishes a Monte Carlo uncertainty propagation framework that simultaneously applies perturbations to five parameters: wind speed, wind direction, reconstructed concentration at each cross-sectional sampling point, background concentration and cross-sectional spacing. By repeatedly executing the cross-sectional flux integration and neighboring source correction process, the P5, P50, P95 quantiles and 90% confidence interval width of each source emission rate are output, providing quantifiable uncertainty basis for emission inventory verification and emission reduction decisions; (4) Internal diagnostic self-checking mechanism, automatic identification and filtering of hidden quantification failures: When the cross-sectional spacing in the upwind and downwind directions is not properly selected or the local meteorological conditions deviate from the steady-state assumption, the cross-sectional method may still output a positive numerical value but the actual emission rate is distorted. Existing methods lack diagnostic means to identify such hidden failures. This invention is based on a combination of four internal indicators: cross-sectional spacing sensitivity variation coefficient, upwind and downwind flux sign consistency, connected domain area rationality, and cross-flight order consistency. It performs credibility classification on the emission quantification results of each source, automatically marks the sources that fail the diagnosis as low confidence and excludes them from the site-level emission summation, effectively preventing distorted data from entering the emission inventory, and significantly improving the robustness and data reliability of the system under complex field conditions.

[0016] In summary, this invention collects concentration and wind field data through mobile observation, automatically identifies emission sources through rasterization, morphological processing, and connected component clustering, obtains baseline emission rates by integrating multiple cross-sections, and deducts neighboring source interference contributions based on inverse distance weighting. It combines Monte Carlo uncertainty propagation of five parameters (wind speed, wind direction, concentration, background, and cross-section spacing) with four internal diagnostic self-checks to achieve automatic source-level emission identification, neighboring source correction, uncertainty quantification, and confidence level discrimination in multi-source scenarios in landfills. This eliminates subjective bias in manual source identification, suppresses systematic overestimation caused by plume overlap, provides reliable confidence intervals, and automatically filters implicit quantification failure results. Attached Figure Description

[0017] Figure 1 This is a flowchart of the method of the present invention.

[0018] Figure 2 This is a schematic diagram of the three-stage processing of rasterization, morphology, and connected components in automatic source identification.

[0019] Figure 3 Cross-sectional layout and spacing d Scan diagram.

[0020] Figure 4 Typical source emission rate Q With cross-sectional spacing d Sensitivity curves and upwind / downwind flux component decomposition diagrams.

[0021] Figure 5 This is a schematic diagram of the geometric relationship for neighbor-source inverse distance weighted correction.

[0022] Figure 6 This is a schematic diagram of the source-by-source error bar output for the five-parameter Monte Carlo uncertainty propagation.

[0023] Figure 7 This is a system architecture diagram of the present invention. Detailed Implementation

[0024] To make the technical problems to be solved, the technical solutions, and the beneficial effects of the present invention clearer, the present invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative of the present invention and are not intended to limit the present invention.

[0025] Example 1 The flowchart of the automatic identification and neighboring source correction quantification method for methane source-level emissions from landfills is as follows: Figure 1 As shown, the specific steps are as follows.

[0026] S1. Data Acquisition and Matching: The methane concentration, effective path length, coordinates, meteorological parameters, and wind speed and direction at each measuring point are collected synchronously by mobile surveying. The calculated horizontal wind component and background concentration are then spatially interpolated to each concentration measuring point to form a point-by-point matching dataset.

[0027] Specifically, methane concentration, effective path length, meteorological parameters, spatial coordinates, wind speed, and wind direction are simultaneously collected at each measuring point through mobile observation. The horizontal wind component is calculated from the wind speed and direction, and the background concentration is statistically analyzed. The wind field data and background concentration are spatially interpolated and matched to each concentration measuring point to obtain a concentration dataset with point-by-point coordinated wind field and background concentration. This concentration dataset is the point-by-point matched dataset. Mobile observation refers to the method of a measurement platform (such as a drone) moving along a preset route above the target area to conduct measurements and simultaneously collect spatial distribution data.

[0028] Specifically, an unmanned aerial vehicle (UAV) platform equipped with a laser methane telemetry instrument and an ultrasonic anemometer was used to conduct near-ground navigation above the target landfill along a preset route, simultaneously recording the following second-by-second data: longitude lon ,latitude lat ,altitude alt methane volume fraction conc ppm and path integral concentration conc ppm ·m, path length h Wind speed ws ,wind direction wd ,temperature T air pressure P .

[0029] Corresponding to the wind observation section, its point-by-point horizontal wind component u , v Calculate according to formula (1): (1), in, u , v These represent the easterly wind component and the northerly wind component, respectively, in m / s; ws Wind speed, in m / s; wd Wind direction, in degrees, using a meteorological coordinate system, with clockwise from true north as positive; π / 180 is the angle-to-radian coefficient; Using a two-dimensional linear interpolation method, the concentration measurement points ( lon , lat Using ) as the query point, match and generate point-by-point. u point , v point For concentration points located outside the convex hull of the wind tunnel, nearest neighbor interpolation is used for extrapolation and filling, followed by... u point , v point The point-by-point wind speed is calculated by reverse calculation. point ws With wind direction point wdSimultaneously, background concentration is also statistically analyzed; the output is a concentration dataset with point-by-point collaborative wind field and background concentration. In this embodiment, the above interpolation operation can be implemented using linear interpolation based on Delaunay triangulation. The above refers to the smallest convex polygon region that encompasses all wind observation points within the convex hull of the wind segment space.

[0030] S2. Automatic source identification: The point-by-point matching dataset is projected onto the UTM coordinate system and rasterized. The detection threshold is determined based on the statistical characteristics of the positive concentration raster and a high-value binary mask is generated. After morphological closure, opening, connected component clustering and minimum raster number filtering, candidate source connected components and their spatial and concentration attributes are obtained.

[0031] Specifically, the concentration dataset, i.e. the point-by-point matching dataset, is projected onto the UTM coordinate system and then rasterized at a preset resolution. The detection threshold is determined based on the statistical characteristics of the positive concentration raster, and a high-value binary mask is generated. Morphological closure and opening operations are sequentially performed on the binary mask, and then connected component clustering and minimum raster number filtering are performed to obtain the candidate source connected components and their spatial and concentration attributes, which are used as source identification output data.

[0032] This step aims to detect independent emission sources from a continuous mobile concentration field in a reproducible, non-human-interventional manner. The image processing terms used in this step are defined as follows: (1) Rasterization: The irregularly distributed concentration measurement points collected by mobile survey are spatially merged into a regular square grid, and the equidistant concentration field is obtained by the arithmetic mean within the cell, which is convenient for subsequent pixel-level algorithms (threshold determination, morphological operation, connected component analysis, etc.) to process.

[0033] (2) Morphological operations: A small structuring element (3×3 fully connected kernel by default in this invention) is used to slide a window on the binary image. The geometric shape of the binary object is normalized by a combination of two basic operations: dilation (taking the union of neighborhoods) and erosion (taking the intersection of neighborhoods). The operation of dilation followed by erosion is called "closing" and is used to fill small holes and depressions inside the object; the operation of erosion followed by dilation is called "opening" and is used to remove isolated bright spots outside the object.

[0034] (3) 8-connected component: In a binary image, if two pixels with a value of 1 are adjacent in the horizontal, vertical or diagonal direction (each pixel has a maximum of 8 neighbors), they are considered connected; all interconnected pixels are grouped into a connected component, and each connected component corresponds to a candidate emission source. D k .

[0035] A schematic diagram of the three-stage processing of rasterization, morphology, and connected components in automatic source identification is shown below. Figure 2 As shown, the specific steps include: S21. Coordinate Projection: Project the latitude and longitude of the WGS84 data set into the local UTM coordinate system to obtain the planar coordinates (x, y). S22. Rasterization: The navigation area is divided into regular square grids according to the resolution dx. The dx value is set according to actual needs. The concentration of all measuring points in each grid is calculated by arithmetic mean to obtain the grid concentration field. C ( i , j Empty grid cells are marked with missing values; in this embodiment, the resolution dx = 4 m; S23. Background and Threshold Determination: The 5th percentile of all positive concentration grid cells is taken as the background value. C bg The positive concentration grid refers to a grid with a positive concentration value, and the standard deviation of the concentration of all positive concentration grids is calculated. σ bg ; Calculate the detection threshold according to formula (2) T : (2), in, T The concentration detection threshold is expressed in ppm·m. C bg The background concentration is the 5th percentile of all positive concentration grids, expressed in ppm·m. σ bg The standard deviation of the background concentration grid is given in ppm·m. Will C ( i , j ) ≥ T The grid is denoted as a high-value binarization mask. M ; S24. Morphological operations: High-value binarization mask M Three morphological closure operations are performed sequentially to eliminate voids caused by the density of the travel paths within a single source; then one morphological opening operation is performed to remove isolated high-value burrs; the structuring element uses a 3×3 fully connected kernel; S25. Connected Component Clustering and Filtering: Clustering is performed on the morphologically processed mask using the 8-neighbor connected component algorithm, and cells with fewer than 8 connected components are removed. N min The tiny connected regions, N min Based on actual needs, the set of candidate connected components is obtained. Dk}; for each candidate source connected component D k Calculate its centroid coordinates ( x k, y k ),area A k (Number of grids × dx) 2 Peak concentration C peak,k With mean concentration C mean,k The above results are then used as source identification output data. In this embodiment, N min Default value 3 corresponds to an area < 48 m² 2 .

[0036] The aforementioned automatic source identification process is robust to changes in flight altitude, flight path density, and single-point peak noise. It can consistently identify dominant sources (such as closed historical landfill areas and active work surfaces) across different flights at the same site, and also has the capability to detect isolated small sources (such as concentrated areas of gas wells). During implementation, dx and [other parameters] can be adjusted according to the overall concentration level of the site. N min To balance spatial resolution and false detection rate.

[0037] S3. Cross-section layout and baseline quantization: Using the candidate source connected domain as the target source, calculate the average wind direction and lay out the normal unit vector based on the wind field of the surrounding buffer zone; within the preset cross-section spacing range, lay out multiple sets of upwind and downwind cross-sections perpendicular to the normal direction at step size; on each cross-section, obtain the reconstructed concentration and the upwind and downwind enhancement distribution of each cross-section sampling point through two-dimensional linear interpolation; integrate the normal wind component to obtain the flux at each spacing, and take the median of the effective flux as the baseline emission rate.

[0038] Specifically, taking the candidate source connectivity as the target source, the average wind direction is calculated based on the cooperative wind field of each measuring point in the surrounding buffer zone, and a normal unit vector is deployed; within the preset cross-sectional spacing range, multiple sets of upwind and downwind cross-sections perpendicular to the normal are deployed in step size; on each cross-section, the reconstructed concentration and its upwind and downwind enhancement distribution of each cross-section sampling point are obtained by two-dimensional linear interpolation; the flux of a single cross-section is obtained by integrating the normal wind component for each spacing, and the median flux of all effective spacings is used as the baseline emission rate of the source.

[0039] Specifically, the cross-sectional layout and cross-sectional spacing d Scan diagram as shown Figure 3 As shown, each candidate source connected component in the source identification output data is represented. D kFor the target source, the emission rate is quantified according to the following steps: S31. Determine the average wind direction: using D k All concentration measurement points within a pre-defined buffer zone (e.g., 100 m) u point , v point The average wind direction of the target source is obtained by vector averaging. wd k and let the normal unit vector n along wd k Pointing downwind; S32. Layout of upwind and downwind cross-sections: with D k centroid coordinates ( x k, y k For reference, along - n with + n The directions are translated by distance respectively d The center point of the upwind cross section is obtained. P up Center point of the downwind cross section P down ; cross-sectional direction and n Vertically, the length is taken to completely cover the projected width of the connected domain in the vertical wind direction and extends to the left and right by a distance P, where P is set according to actual needs, and the cross-sectional spacing is... d Within a preset range, values ​​are selected one by one according to the step size to form... d Scanning sequence; in this embodiment, the extension distance P = 20m; cross-sectional spacing d The values ​​are selected one by one within the range of {10, 15, 20, 25, 30, 35, 40, 45, 50, 55, 60} m, that is, the preset value range is [10, 60], and the step size is 5m; S33. Concentration Field Reconstruction: Sampling points are generated along the vertical direction on each cross-section with a step size of length p, where p is set according to actual needs. The concentration is reconstructed at these sampling points using two-dimensional linear interpolation. C (s); For sampling points located outside the convex hull of the mobile data (the smallest convex polygon region enclosing the corresponding point set), the nearest neighbor interpolation result is multiplied by an attenuation factor of 0.5 to compensate for the difference, so as to avoid unreasonable concentration increments in the extrapolated region; In this embodiment, the length p = 1 m; S34. Emission rate calculation: The normal wind component is taken at the center point of the cross-section in the upwind and downwind directions. u point , vpoint and n dot product average U ⊥, d The emission rate at the cross-sectional spacing d is calculated according to formula (3): (3), in, Q d Spacing between sections d The methane emission rate under the specified conditions, in kg / h; s The coordinates are integrals along the cross-section, in meters (m). C down ( s ), C up ( s The cross-sections on the downwind and upwind sides are respectively located on the coordinate system. s The reconstructed path integral concentration, in ppm·m; U ⊥,d The normal wind component of the cross section about the center is taken at the center point of the cross section in the upwind and downwind directions. u point , v point With the normal unit vector n The dot product average, in m / s; h ( s ) represents the effective path length at the sampling point, in meters (m). ρ corr The temperature and pressure standard condition correction factor is dimensionless and calculated according to formula (4); the integral variable ds The length of the infinitesimal element along the cross-section is expressed in meters (m). in h (s) represents the effective path length at the sampling point, which is the local mean of the measured path length during navigation. ρ corr The temperature and pressure standard condition correction factor is calculated according to formula (4), converting the volume fraction ppm to the mass concentration under standard conditions: (4), in, ρ corr This is a dimensionless temperature and pressure standard condition correction factor used to convert volume fraction (ppm) to mass concentration under standard conditions. P The pressure is the ambient air pressure, in Pa. P 0 = 101325 Pa is the standard atmospheric pressure; T 0 = 273.15K is the standard temperature;T The ambient temperature is expressed in Kelvin (K). S35, Baseline Q Select: In d Take all values ​​within the range [A, B]m Q d The median is used as the baseline emission rate for that source. Q base,k Where [A, B] is the preset value range in S32, and is recorded. Q d The mean, standard deviation, and upwind and downwind flux components vary with d The change curve is used for internal diagnostic self-testing of S6. In this embodiment, d ∈ [10, 60].

[0040] Among them, typical source emission rates Q With cross-sectional spacing d The sensitivity curve and the decomposition diagram of upwind and downwind flux components are shown below. Figure 4 As shown.

[0041] S4. Neighboring Source Inverse Distance Weighted Correction: Other identified sources located in the downwind half-space and within a set radius from the center of the target source are selected as neighboring sources. The weights are calculated based on the inverse square of the distance and the minimum weight limit. The contribution is estimated by multiplying the difference between the downwind and upwind concentration increases of the neighboring sources by the weights. After deducting the downwind increase from the target source, the emission rate is re-integrated to obtain the corrected emission rate.

[0042] Specifically, based on the distribution of upwind and downwind concentration enhancement at each cross-section spacing, other identified sources located within the downwind half-space of the target source center without exceeding a set radius are selected as neighboring sources. They are assigned weights inversely proportional to the square of the distance and a minimum weight limit is set. The contribution of each neighboring source is estimated by multiplying its downwind and upwind concentration enhancement by the corresponding weight. The contributions of each neighboring source are summed and subtracted from the downwind enhancement of the target source before being re-integrated to obtain the emission rate corrected for the neighboring sources.

[0043] Specifically, the geometric relationship diagram of neighbor-source inverse distance weighted correction is as follows: Figure 5 As shown, when multiple sources are distributed close to each other within the site, the downwind cross section of the target source will simultaneously capture the contributions of neighboring sources. The following method is used to subtract these neighboring source contributions: S41. Neighbor source filtering: For target sources D k traverse the remaining recognition sources D m Calculate the Euclidean distance between sources r km Only retain r km ≤ Rmax The source is used as a candidate neighboring source; in this embodiment, R max The default value is 80 m; S42, Directional Gating: Further requirements for neighboring sources D m Located at the target source D k Within the downwind half-space, that is, satisfying ( x m - x k, y m - y k ) · n k > 0; if this condition is not met, the contribution from neighboring sources will not be counted. S43, Inverse Distance Weighted Weight: Neighborhood D m The concentration contribution weight of the downwind cross section of the target source is calculated according to formula (5): (5), in, w m Neighboring source D m For the target source D k The inverse distance-weighted weight of the downwind cross-sectional concentration contribution is dimensionless. r km For target source D k and neighboring sources D m Euclidean distance in the UTM coordinate system, in meters; w min The minimum weight lower bound is dimensionless and used to avoid long-distance neighbor factors. r km The square decays too quickly and loses its physical meaning; in, w min The default value is 0.05, used to avoid the influence of distant neighbor sources. r km The rapid decay renders it meaningless physically, and it also serves as an approximate compensation for the lateral expansion of the actual plume.

[0044] S44. Contribution Estimation: Neighboring Sources D m The contribution to the concentration in the downwind cross section of the target source is estimated by multiplying the downwind concentration increase by the upwind concentration increase in its own cross section by a weight: (6), Where, Δ C m ( s ) are neighboring sources D m The cross section downwind of the target source in coordinates s Concentration contribution at a given location, in ppm·m; Δ C m_self ( s ) are neighboring sources D m The concentration increase from downwind to upwind obtained in its own quantification, in ppm·m; in ΔC m_self (s) are neighboring sources D m The downwind increase minus the upwind increase obtained from its own quantification; if ΔC m_self (s) If not reconstructed within the same coordinate range of the target cross-section, then by ΔC m_self Substitute the mean along its normal direction; S45. Corrected emission rate: Sum the contributions of all candidate neighboring sources and subtract the downwind enhancement from the target source, then re-integrate to obtain the corrected emission rate of the neighboring sources, and calculate according to formulas (7) and (8): (7), (8), Where, Δ C corr ( s The target source's downwind cross section after neighbor source correction is shown in coordinates. s The concentration increase at a given location, expressed in ppm·m; the first term on the right-hand side of the equation represents the uncorrected downwind minus upwind increase from the target source; summation symbol Σ. m Indicates distance from the target source k The central Euclidean distance does not exceed R max And all candidate neighboring sources located in its downwind half-space D m Summation; Q corr,k For target source D k Emission rates corrected for adjacent sources, in kg / h; U ⊥ The normal wind component at the center of the cross-section downwind of the target source, in m / s; h ( s ), ρ corrSame as formula (3); ds The length of the infinitesimal element along the cross-section is expressed in meters (m). In practice, the median of all cross-sectional spacings d is deducted; if the deduction results in Δ C corr ( s If a negative value appears at a sampling point, the contribution of that point is set to zero to avoid single-point inversion caused by excessive deduction from neighboring sources; S46. Correction Amplitude Report: Define the correction amplitude δ k Calculate according to formula (9): (9), in, δ k For target source D k The neighbor source correction magnitude is dimensionless and can be reported as a percentage. Q corr,k Emission rate after correction, unit: kg / h; Q base,k Baseline emission rates before correction, in kg / h. The more significant the contribution from neighboring sources... δ k The further away from 0; δ k A value close to 0 indicates that interference from neighboring sources is negligible. The values ​​are then given in the output. Q base,k , Q corr,k and δ k This makes it easier for users to determine the significance of neighbor source correction on that source.

[0045] S5. Multi-parameter Monte Carlo uncertainty propagation: Under fixed random seed conditions, apply relative or absolute random perturbations of preset amplitude to five parameters: wind speed, wind direction, reconstructed concentration at each cross-section sampling point, background concentration, and cross-section spacing. Repeat steps S3 and S4 multiple times to output the P5, P50, and P95 quantiles and the 90% confidence interval width for each source emission rate.

[0046] Specifically, to obtain the confidence interval for the emission rate of each source, multi-parameter Monte Carlo sampling is embedded based on steps S3 and S4. A schematic diagram of the source-by-source error bar output for the five-parameter Monte Carlo uncertainty propagation is shown below. Figure 6 As shown, the specific steps are as follows: S51. Set the disturbance parameters and default amplitude: Wind speed: Obey N (1, 0.10 2The relative perturbation of ), that is, each iteration will U Multiply the whole by a random coefficient that follows this distribution; Wind direction: Obey N (0, 7.5 2 An absolute perturbation of 0° with respect to the normal direction n Rotate; Reconstructed concentrations at each cross-sectional sampling point: follow the rules N (1, 0.10 2 The relative perturbation of ) affects the reconstructed concentration of each cross-sectional sampling point obtained in the S3 cross-sectional layout and baseline quantization. C (s); Background concentration: obeys N (1, 0.05 2 The relative perturbation of ) acts on the background concentration determined in S2 automatic source identification. C bg ; Section spacing d :obey N (1, 0.05 2 The relative disturbance of ) acts on the cross-sectional spacing of the upwind and downwind cross-sections in the S3 cross-section layout and baseline quantization; S52. Fixed random seed: A pseudo-random sequence with seed = 42 is used; S53, Number of iterations: Execution N mc = 500 to 1000 iterations; in each iteration, re-execute steps S3 and S4 with the above perturbation to obtain a set of { Q corr,k}_ i , i = 1… N mc ; S54. Output statistics: P5, P50, and P95 quantiles and 90% confidence interval width for the output emission rate of each source.

[0047] Therefore, in subsequent operations, the total site emission can be successively summed based on the corresponding sampling of each source emission rate to output the site-level total emission distribution. In actual measurements, the relative uncertainty of a single source is about 25-27%, and the relative uncertainty of the site level is about 19-22%, with the latter being lower due to the partial hedging effect of multi-source aggregation.

[0048] S6. Internal Diagnostic Self-Check: Based on four indicators, namely, cross-sectional spacing sensitivity variation coefficient, upwind and downwind flux sign, connected domain area anomaly, and cross-flight quantity consistency, the reliability of the emission quantification results of each source is judged; those that fail are marked as low confidence and are not included in the site-level emission summation.

[0049] The definition of the "positive but low confidence" failure mode is as follows: In the engineering practice of the cross-section method, when the site average wind speed is low, the wind direction dispersion is large, or there is a non-target source contribution between the upwind and downwind cross sections, or the plume has not yet fully developed in the downwind direction, the value calculated according to formula (3) is as follows. Q d The value may still be positive and fall within a reasonable historical range, making it visually indistinguishable from a true and valid result. However, such results actually violate one or more of the three basic assumptions upon which the cross-sectional method relies: "steady-state advection dominance, a clean upwind background, and complete downwind plume capture." Therefore, their value does not represent the true emission rate of the target source. This invention refers to such results—"positive values, seemingly reasonable magnitudes, but with invalid underlying assumptions"—as "positive but low-confidence" failure modes. Unlike explicit failures such as negative emission rates or non-physically large values ​​(which can be filtered out by simple thresholds), "positive but low-confidence" failures are latent and must be identified using a combination of the following four internal diagnostic indicators.

[0050] Specifically, to prevent "positive but low-confidence" failure modes from entering the emissions inventory, the following four indicators are used to determine the credibility of the quantification results of the emission rate for each source: First, the calculation of the coefficient of variation of the cross-sectional spacing sensitivity curve: In d Calculate according to formula (10) within the range [A, B]m. Q d coefficient of variation CV d ,like CV d A value > 0.50 is marked as low confidence. (10), in, CV d The target source emission rate varies with the cross-sectional spacing. d The coefficient of variation, dimensionless; σ( Q d ) for d All valid ranges from A to B Q d Standard deviation, in kg / h; mean ( Q d Within this range Q d The arithmetic mean of the values, in kg / h; |·| indicates taking the absolute value; CV d Larger means Q Spacing of sections dThe stronger the selected dependency, the more it suggests that the plume pattern deviates from the steady-state advection assumption; in this embodiment, d ∈ [10, 60] m;σ( Q d ) for d All effective ranges from 10 to 60 m Q d Standard deviation; Second, check the signs of the upwind and downwind flux components: separately count the signed fluxes in the upwind and downwind cross sections. F up , F down And compare their absolute values; if |F up | > |F down | or F up and F down If the numbers are the same and the magnitudes are close, it indicates that the transmission direction assumption has failed, meaning that the upwind direction has been contaminated by a non-target source or the wind direction has reversed, and it is marked as low confidence. Third, diagnosis of abnormal connected area: When the connected area of ​​the same source in different flights at the same site expands by more than 2 times and corresponds to the decrease in the average wind speed of the site, it indicates that the plume accumulates nearby under low wind speed and wind direction change conditions, and the result of this source in this flight is marked as low confidence. Fourth, cross-flight consistency diagnosis: If the same source can be quantified in more than two flights, and the deviation of the results in each flight exceeds 5 times, it indicates that at least one of them failed, and further diversion is carried out according to internal diagnostic indicators.

[0051] The aforementioned internal diagnostic self-check can determine the reliability of single-flight results without relying on external references. Sources marked as low confidence are not included in the site-level emission summation, but their emission rate values ​​and diagnostic indicators are retained in the results file for retrospective purposes.

[0052] Example 2 The architecture diagram of the automatic identification and neighboring source correction quantification system for methane source-level emissions from landfills is as follows: Figure 7 As shown, it consists of a data acquisition and matching module, an automatic source identification module, a cross-sectional layout and baseline quantization module, a neighboring source inverse distance weighted correction module, a multi-parameter Monte Carlo uncertainty propagation module, and an internal diagnostic self-test module.

[0053] Data acquisition and matching module: used to synchronously collect methane concentration, effective path length, meteorological parameters, spatial coordinates, wind speed and wind direction at each measuring point through mobile observation; calculate the horizontal wind component from wind speed and wind direction, and statistically analyze the background concentration; and match the wind field data and background concentration to each concentration measuring point through spatial interpolation to obtain a concentration dataset with point-by-point coordinated wind field and background concentration. Automatic source identification module: connected to the data acquisition and matching module, used to project the concentration dataset onto the UTM coordinate system and rasterize it at a preset resolution, determine the detection threshold based on the statistical characteristics of the positive concentration raster and generate a high-value binary mask; perform morphological closure and opening operations on the binary mask in sequence, and then perform connected component clustering and minimum raster number filtering to obtain candidate source connected components and their spatial and concentration attributes, which are used as source identification output data; Cross-section layout and baseline quantization module: connected to the automatic source identification module, used to calculate the average wind direction and lay out the normal unit vector based on the cooperative wind field of each measuring point in the surrounding buffer zone, with the candidate source connectivity as the target source; lay out multiple sets of upwind and downwind cross-section pairs perpendicular to the normal direction in step size within the preset cross-section spacing range; obtain the reconstructed concentration and upwind and downwind enhancement distribution of each cross-section sampling point through two-dimensional linear interpolation on each cross-section; calculate the single-section flux by integrating the normal wind component for each spacing, and use the median flux of all effective spacings as the baseline emission rate of the source; Neighboring source inverse distance weighted correction module: connected to the cross-sectional layout and baseline quantification module, used to screen other identified sources as neighboring sources based on the upwind and downwind concentration enhancement distribution of each cross-sectional spacing, which are no more than a set radius from the center of the target source and located in its downwind half-space. The neighboring sources are assigned weights according to the inverse square distance ratio and a minimum weight lower limit is set. The contribution of each neighboring source is estimated by multiplying its downwind and upwind concentration enhancement by the corresponding weight. The contributions of each neighboring source are summed and subtracted from the downwind enhancement of the target source, and then integrated again to obtain the emission rate after neighboring source correction. Multi-parameter Monte Carlo uncertainty propagation module: connected to the cross-section layout and baseline quantization module and the neighboring source inverse distance weighted correction module, used to simultaneously apply relative or absolute random perturbations of preset amplitude to five parameters—wind speed, wind direction, reconstructed concentration of each cross-section sampling point, background concentration, and cross-section spacing—under fixed random seed conditions, repeatedly call the cross-section layout and baseline quantization module and the neighboring source inverse distance weighted correction module multiple times, and output the P5, P50, and P95 quantiles and the 90% confidence interval width of the emission rate of each source; Internal diagnostic self-test module: connected to the cross-sectional layout and baseline quantification module, the adjacent source inverse distance weighted correction module, and the multi-parameter Monte Carlo uncertainty propagation module, it is used to judge the credibility of the emission quantification results of each source based on four indicators: cross-sectional spacing sensitivity variation coefficient, upwind and downwind flux sign, connected domain area anomaly, and cross-flight order consistency; those that fail are marked as low confidence and are not included in the site-level emission summation.

[0054] The specific implementation methods of each module in this system are the same as those described in Example 1, and will not be repeated here.

[0055] To enable the implementation of the system in this embodiment, the present invention may provide the following structures or components: (1) Mobile measurement platform: A multi-rotor or fixed-wing UAV equipped with a laser methane telemetry instrument and an ultrasonic anemometer is used to conduct near-ground mobile surveys above the target landfill along a preset route.

[0056] (2) Laser methane telemetry instrument: installed at the bottom of the mobile measurement platform, it emits a near-infrared laser beam downwards / to the side to the ground surface / opposite reflector, and integrates the methane concentration based on the absorption spectrum inversion path. conc ppm ·m, and simultaneously output the effective path length of the measuring point. h And by dividing the path integral concentration by h Obtain the equivalent volume fraction conc ppm .

[0057] (3) Ultrasonic anemometer: synchronously samples with laser methane telemetry instrument and outputs wind speed. ws With wind direction wd The ultrasonic anemometer collects data synchronously during a dedicated wind observation flight segment consisting of a horizontal straight flight segment at a relatively high wind altitude on the target site, and then interpolates and matches the data to the concentration measurement point (see step S1 of Example 1 for details).

[0058] (4) Meteorological and positioning unit: temperature data collection T air pressure P Longitude lon, latitude lat, and altitude alt are used for the standard condition correction in step S3 of Example 1 and the coordinate projection in step S2 of Example 1.

[0059] (5) Data processing unit: can be an airborne embedded computer or a ground station workstation, storing a computer program that implements the entire process of Embodiment 1, and simultaneously outputting data for each identification source. Q base , Q corrThe computer program determines the correction amplitude δ, P5 / P50 / P95 quantiles, 90% confidence interval width, and internal diagnostic indicators. It can also be stored on a separate computer-readable storage medium for execution by other processors.

[0060] Example 3 An electronic device includes a memory and a processor, the memory and the processor being communicatively connected to each other, the memory storing computer instructions, and the processor executing the computer instructions to implement the automatic identification and neighboring source correction quantification method for landfill methane source-level emissions as described in Embodiment 1 above, and the automatic identification and neighboring source correction quantification system for landfill methane source-level emissions as described in Embodiment 2.

[0061] Example 4 A computer-readable storage medium storing a computer program that, when executed by a processor, implements the automatic identification and neighboring source correction quantification method for landfill methane source-level emissions as described in Embodiment 1 above, and the automatic identification and neighboring source correction quantification system for landfill methane source-level emissions as described in Embodiment 2.

[0062] Example 5 Taking an operational municipal solid waste landfill as an example, the low-altitude mobile observation data of the same site are processed according to the method of the present invention. This includes at least one flight with good weather conditions (high average wind speed and stable wind direction) and one flight with poor weather conditions (low average wind speed and large wind direction dispersion) to cover the applicability and failure of the method.

[0063] After processing through steps S1 to S2, multiple candidate emission sources can be automatically detected in both the closed historical landfill unit (divided into several sub-areas) and the active operating area. The main source of the operating area and the main source of the closed area are consistently identified across different flights, and the source identification is repeatable across flights without any manual intervention.

[0064] The baseline emission rates of each source are obtained through cross-sectional quantization in step S3. Q base Then, step S4, neighbor source correction, is performed to obtain... Q corr In multi-source scenarios, the neighbor source correction magnitude δ for isolated large sources (such as the main source on the work surface far from other sources) is close to zero; while the neighbor source correction magnitude for small sources close to strong neighboring sources can reach tens of percentage points, fully verifying the necessity and effectiveness of the neighbor source correction method of this invention in multi-source dense scenarios. The change in site-level total emissions before and after neighbor source correction is usually small (on the order of several percentage points), because neighbor source deduction mainly redistributes emissions among sources, and has a limited impact on the total site emissions.

[0065] Through uncertainty propagation in step S5, the emission rate is obtained for each source. Q The P5, P50, and P95 quantiles and the 90% confidence interval widths were determined. Among comparable flights with favorable weather conditions, the 90% confidence interval widths at the site level overlapped, indicating cross-flight consistency.

[0066] For data from flights with poor weather conditions, the internal diagnosis of step S6 of this invention is applied: During the source identification stage, the connected domain area shows abnormal expansion (the area of ​​the same region can be several times that of flights with good weather conditions), triggering the indicator: abnormal connected domain area diagnosis; during the cross-sectional quantization stage, some sources show upwind and downwind flux sign reversal. Q A non-physical negative value triggers the following indicators: upwind and downwind flux component sign check; although some sources... Q A positive value, but differing from other flights in the same area by more than an order of magnitude, triggers the indicator: cross-flight order-of-magnitude consistency diagnosis, i.e., the "positive value but low confidence" latent failure defined in this invention. All of the above sources are automatically marked as low confidence by the internal diagnostic module of this invention and excluded from the site-level aggregate. This embodiment demonstrates that the internal diagnostic module of this invention can reliably identify failure modes of existing cross-sectional methods under low wind speed and strong wind direction changes, including the most difficult-to-detect "positive value but low confidence" situation, without relying on any external references.

[0067] It must be noted that the raster resolution dx = 4 m and the minimum connected component given in this specific embodiment are... N min = 3, threshold multiple 3σ, morphological closure 3 times and opening 1 time, cross-sectional spacing d Range 10-60 m, step size 5 m, neighboring source radius R max = 80 m, neighbor-source weighted power = 2, minimum weight w min The specific values, such as 0.05, the number of Monte Carlo iterations (500-1000), and the amplitude of the five-parameter perturbation, are preferred values ​​that have yielded good results when this invention was implemented in a municipal solid waste landfill, and do not constitute a limitation on the scope of this invention. Those skilled in the art can adaptively adjust the above parameters based on the specific site size, source density, navigation route density, measurement platform performance, and other conditions, without departing from the technical concept of this invention.

[0068] Those skilled in the art will understand that embodiments of the present invention can be provided as methods, systems, or computer program products. Therefore, the present invention can take the form of a completely hardware embodiment, a completely software embodiment, or an embodiment combining software and hardware aspects. Furthermore, the present invention can take the form of a computer program product implemented on one or more computer-usable storage media (including but not limited to disk storage, CD-ROM, optical storage, etc.) containing computer-usable program code. The solutions in the embodiments of the present invention can be implemented in various computer languages, such as object-oriented programming languages ​​like Java, C++, Python, and interpreted scripting languages ​​like JavaScript.

[0069] This invention is described with reference to flowchart illustrations and / or block diagrams of methods, electronic devices (systems), and computer program products according to embodiments of the invention. It will be understood that each block of the flowchart illustrations and / or block diagrams, and combinations of blocks in the flowchart illustrations and / or block diagrams, can be implemented by computer program instructions. These computer program instructions can be provided to a processor of a general-purpose computer, special-purpose computer, embedded processor, or other programmable data processing electronic device to produce a machine, such that the instructions, which execute via the processor of the computer or other programmable data processing electronic device, generate instructions for implementing the flowchart illustrations and / or block diagrams. Figure 1 One or more processes and / or boxes Figure 1 A device that provides the functions specified in one or more boxes.

[0070] These computer program instructions may also be stored in a computer-readable storage medium that can direct a computer or other programmable data processing electronic device to operate in a particular manner, such that the instructions stored in the computer-readable storage medium produce an article of manufacture including instruction means, which are implemented in a process Figure 1 One or more processes and / or boxes Figure 1 The function specified in one or more boxes.

[0071] These computer program instructions can also be loaded onto a computer or other programmable data processing electronic device to cause a series of operational steps to be performed on the computer or other programmable electronic device to produce a computer-implemented process, thereby providing instructions that execute on the computer or other programmable electronic device for implementing the process. Figure 1 One or more processes and / or boxes Figure 1 The steps of the function specified in one or more boxes.

[0072] Although preferred embodiments of the invention have been described, those skilled in the art, upon learning the basic inventive concept, can make other changes and modifications to these embodiments. Therefore, the appended claims are intended to be interpreted as including both the preferred embodiments and all changes and modifications falling within the scope of the invention.

[0073] Obviously, those skilled in the art can make various modifications and variations to this invention without departing from its spirit and scope. Therefore, if these modifications and variations fall within the scope of the claims of this invention and their equivalents, this invention also intends to include these modifications and variations.

Claims

1. A method for automatic identification and neighboring source correction quantification of methane source-level emissions from landfills, characterized in that, Includes the following steps: S1. Data Acquisition and Matching: The methane concentration, effective path length, coordinates, meteorological parameters, and wind speed and direction of the measuring points are collected synchronously by mobile surveying. The horizontal wind component calculated from the wind speed and direction and the background concentration are spatially interpolated to the concentration measuring points to form a point-by-point matching dataset. S2. Automatic source identification: The dataset is projected onto the UTM coordinate system and rasterized. The threshold is determined according to the statistical characteristics of the positive concentration raster to generate a high-value mask. After morphological processing, connected component clustering and minimum raster number filtering, candidate source connected components and attributes are obtained. S3. Cross-section layout and baseline quantification: Using candidate sources as target sources, the average wind direction and normal direction are determined based on the surrounding coordinated wind field. Multiple sets of upwind and downwind cross-sections are laid out within a preset spacing range. The cross-section concentration is reconstructed by interpolation. The median effective flux is taken as the baseline emission rate by integrating the normal wind component. S4. Neighboring source inverse distance weighted correction: Screen neighboring sources within the downwind half-space of the target source and whose distance does not exceed the set radius, estimate their concentration contribution according to the square weight of inverse distance and deduct it from the downwind enhancement amount of the target source, and re-integrate to obtain the corrected emission rate. S5. Multi-parameter Monte Carlo uncertainty propagation: Perturbations are applied to wind speed, wind direction, cross-sectional reconstructed concentration, background concentration and cross-sectional spacing. Repeat S3 and S4 to obtain quantiles and confidence intervals. S6. Internal diagnostic self-check: The confidence level is determined by the sensitivity of cross-sectional spacing, the sign of upwind and downwind flux, the anomaly of the area of ​​the connected domain, and the consistency of the magnitude of cross-flights. Low confidence results are not included in the site-level summation.

2. The method according to claim 1, characterized in that: In step S1, an unmanned aerial vehicle (UAV) platform equipped with a laser methane telemetry instrument and an ultrasonic anemometer is used to conduct near-ground navigation above the target landfill along a preset route, simultaneously recording the following second-by-second data: longitude. lon ,latitude lat ,altitude alt ; methane volume fraction conc ppm and path integral concentration conc ppm·m Path length h ; wind speed ws ,wind direction wd ,temperature T air pressure P ; Corresponding to the wind observation section, its point-by-point horizontal wind component u , v Calculate using the following formulas respectively: ， in, u , v These are the easterly wind component and the northerly wind component, respectively; ws Wind speed; wd To indicate wind direction, a meteorological coordinate system is used, with clockwise rotation from due north as positive. Using a two-dimensional linear interpolation method, the concentration measurement points ( lon , lat Using ) as the query point, match and generate point-by-point. u point , v point For concentration points located outside the convex hull of the wind tunnel, nearest neighbor interpolation is used for extrapolation and filling, followed by... u point , v point The point-by-point wind speed is calculated by reverse calculation. point ws With wind direction point wd Simultaneously, background concentration is also statistically analyzed; the output is a concentration dataset with point-by-point collaborative wind field and background concentration.

3. The method according to claim 2, characterized in that: The specific steps of S2 include: S21. Coordinate Projection: Project the latitude and longitude of the WGS84 data set into the local UTM coordinate system to obtain the planar coordinates (x, y). S22. Rasterization: The navigation area is divided into regular square grids according to the resolution dx. The dx value is set according to actual needs. The concentration of all measuring points in each grid is calculated by arithmetic mean to obtain the grid concentration field. C ( i , j Empty rasters are marked with missing values; S23. Background and Threshold Determination: The 5th percentile of all positive concentration grid cells is taken as the background value. C bg The positive concentration grid refers to a grid with a positive concentration value, and the standard deviation of the concentration of all positive concentration grids is calculated. σ bg Calculate the detection threshold using the following formula. T : ， Will C ( i , j ) ≥ T The grid is denoted as a high-value binarization mask. M ; S24. Morphological operations: High-value binarization mask M Three morphological closure operations are performed sequentially to eliminate voids caused by the density of the travel paths within a single source; then one morphological opening operation is performed to remove isolated high-value burrs; the structuring element uses a 3×3 fully connected kernel; S25. Connected Component Clustering and Filtering: Clustering is performed on the morphologically processed mask using the 8-neighbor connected component algorithm, and cells with fewer than 8 connected components are removed. N min The tiny connected regions, N min Based on actual needs, the set of candidate connected components is obtained. D k }; for each candidate source connected component D k Calculate its centroid coordinates ( x k, y k ),area A k (Number of grids × dx) 2 Peak concentration C peak,k With mean concentration C mean,k The above results are then used as the source identification output data.

4. The method according to claim 3, characterized in that: S3 includes: S31. Determine the average wind direction: using candidate source connected components D k Measurement points within the surrounding preset buffer zone u point , v point By performing a vector average, the average wind direction can be obtained. wd k and let the normal unit vector n along wd k Pointing downwind; S32. Layout of upwind and downwind cross-sections: with D k centroid coordinates ( x k, y k For reference, along - n with + n The directions are translated by distance respectively d Determine the center of the cross-section in the upwind and downwind directions; the cross-section direction is... n Vertical, the length covers the projection of the connected domain in the vertical wind direction and extends to both sides P, the cross-sectional spacing d Within a preset range, values ​​are selected one by one according to the step size to form... d Scan sequence; S33. Concentration Field Reconstruction: Sampling points are generated at step sizes on each cross-section, and the concentration is reconstructed using two-dimensional linear interpolation. C (s); Samples outside the convex hull are compensated by nearest neighbor interpolation and multiplication by the attenuation factor. S34. Emission rate calculation: The normal wind component is taken at the center point of the cross-section in the upwind and downwind directions. u point , v point and n dot product average U ⊥, d The methane emission rate at the cross-sectional spacing d Q d Calculate using the following formula: ， in h (s) represents the effective path length at the sampling point, which is the local mean of the measured path length during navigation. ρ corr The temperature and pressure standard condition correction factor is calculated using the following formula to convert the volume fraction (ppm) to the mass concentration under standard conditions: ， in, P This refers to the on-site air pressure. P 0 Standard atmospheric pressure; T 0 = 273.15 K; T The ambient temperature; S35, Baseline Q base,k Selection: Select the effective range within the preset spacing. Q d The median emission rate was used as the baseline emission rate. Q base,k and record Q d Statistics and follow-up d The change curve is for internal diagnostic use.

5. The method according to claim 4, characterized in that: In step S4, when multiple sources are distributed close to each other within the site, the downwind cross section of the target source will simultaneously capture the contribution of neighboring sources. The following method is used to subtract this contribution: S41. Neighbor source filtering: For target sources D k traverse the remaining recognition sources D m Calculate the Euclidean distance between sources r km Only retain r km ≤ R max The source is used as a candidate neighbor; S42, Directional Gating: Further requirements for neighboring sources D m Located at the target source D k Within the downwind half-space, that is, satisfying ( x m - x k, y m - y k ) · n k > 0; if this condition is not met, the contribution from neighboring sources will not be counted. S43, Inverse Distance Weighted Weight: Neighborhood D m Concentration contribution weight of the downwind cross section of the target source Calculate using the following formula: ， S44. Contribution Estimation: Neighboring Sources D m The contribution to the concentration in the downwind cross section of the target source is estimated by multiplying the downwind concentration increase by the upwind concentration increase in its own cross section by a weight: ， S45. Corrected emission rate: Summing all candidate neighboring source contributions and subtracting the downwind enhancement from the target source, then re-integrating, yields the corrected methane emission rate from the neighboring source. Calculate using the following formula: ， ， S46. Correction Amplitude Report: Define the correction amplitude δ k Calculate using the following formula: 。 6. The method according to claim 5, characterized in that: In step S5, to obtain the confidence interval for the emission rate of each source, multi-parameter Monte Carlo sampling is embedded based on steps S3 and S4. The specific steps are as follows: S51. Set the disturbance parameters and default amplitude: Wind speed: Obey N (1, 0.10 2 The relative perturbation of ), that is, each iteration will U Multiply the whole by a random coefficient that follows this distribution; Wind direction: Obey N (0, 7.5 2 An absolute perturbation of 0° with respect to the normal direction n Rotate; Reconstructed concentrations at each cross-sectional sampling point: follow the rules N (1, 0.10 2 The relative perturbation of ) affects the reconstructed concentration of each cross-sectional sampling point obtained in the S3 cross-sectional layout and baseline quantization. C (s); Background concentration: obeys N (1, 0.05 2 The relative perturbation of ) acts on the background concentration determined in S2 automatic source identification. C bg ; Section spacing d :obey N (1, 0.05 2 The relative disturbance of ) acts on the cross-sectional spacing of the upwind and downwind cross-sections in the S3 cross-section layout and baseline quantization; S52. Fixed random seed: A pseudo-random sequence with seed = 42 is used; S53, Number of iterations: Execution N mc = 500 to 1000 iterations; In each iteration, steps S3 and S4 are re-executed with the above perturbation to obtain a set of { Q corr,k }_ i , i = 1… N mc ; S54. Output statistics: P5, P50, and P95 quantiles and 90% confidence interval width for the output emission rate of each source.

7. The method according to claim 6, characterized in that: In step S6, the four indicators for judging the credibility of the quantification results of each source emission rate are as follows: First, the calculation of the coefficient of variation of the cross-sectional spacing sensitivity curve: (The text abruptly ends here, so the translation stops as well.) d Calculate within the preset value range using the following formula Q d coefficient of variation CV d ,like CV d A value > 0.50 is marked as low confidence. ， in, CV d The target source emission rate varies with the cross-sectional spacing. d The coefficient of variation, dimensionless; σ( Q d ) for d All valid values ​​within the preset range Q d Standard deviation, in kg / h; mean ( Q d Within this range Q d The arithmetic mean of the values, in kg / h; |·| indicates taking the absolute value; CV d Larger means Q For the cross-sectional spacing d The stronger the selected dependence, the more it suggests that the plume pattern deviates from the steady-state advection assumption; Second, check the signs of the upwind and downwind flux components: separately count the signed fluxes in the upwind and downwind cross sections. F up , F down And compare their absolute values; if |F up | > |F down | or F up and F down If the numbers are the same and the magnitudes are close, it indicates that the transmission direction assumption has failed, meaning that the upwind direction has been contaminated by a non-target source or the wind direction has reversed, and it is marked as low confidence. Third, diagnosis of abnormal connected area: When the connected area of ​​the same source in different flights at the same site expands by more than 2 times and corresponds to the decrease in the average wind speed at the site, it indicates that the plume accumulates nearby under low wind speed and wind direction change conditions, and the result of this source in this flight is marked as low confidence. Fourth, cross-flight consistency diagnosis: If the same source can be quantified in more than two flights, and the deviation of the results in each flight exceeds 5 times, it indicates that at least one of them failed, and further diversion is carried out according to internal diagnostic indicators.

8. An automatic identification and neighboring source correction quantification system for methane source-level emissions from landfills, used to execute the method described in any one of claims 1-7, characterized in that, include: Data acquisition and matching module: used to collect methane concentration, effective path length, coordinates, meteorological parameters and wind speed and direction at the measuring points, and interpolate to form a point-by-point matching dataset; Automatic source identification module: used to project and rasterize datasets, generate high-value masks based on positive concentration raster statistical features, and obtain candidate source connected components and attributes through morphological processing and connected component filtering; Cross-section layout and baseline quantization module: used to lay upwind and downwind cross-section pairs based on the coordinated wind field around the candidate source, interpolate and reconstruct the cross-section concentration and integrate to obtain the baseline emission rate; Neighboring source inverse distance weighted correction module: used to screen downwind neighboring sources, estimate their contribution according to inverse distance weight, and subtract their contribution to obtain the corrected emission rate; Multi-parameter Monte Carlo uncertainty propagation module: used to perturb wind speed, wind direction, cross-sectional reconstructed concentration, background concentration and cross-sectional spacing and output quantiles and confidence intervals; Internal diagnostic self-test module: used to determine confidence level based on cross-sectional spacing sensitivity, flux sign, connected domain area anomaly and cross-flight consistency, and mark those that fail as low confidence and not included in the site-level emission summation.

9. An electronic device, characterized in that, include: The system includes a memory and a processor, which are interconnected. The memory stores computer instructions, and the processor executes these computer instructions to implement the automatic identification and neighboring source correction quantification method for landfill methane source-level emissions as described in any one of claims 1-7.

10. A computer-readable storage medium storing a computer program, characterized in that: When the computer program is executed by the processor, it implements the automatic identification and neighboring source correction quantification method for methane source-level emissions from landfills as described in any one of claims 1-7.

Images