A multi-source remote sensing collaborative screening, risk zoning and sampling site method for water-soil-air combined pollution in industrial parks

Abstract

A kind of multi-source remote sensing collaborative screening, risk zoning and sampling point method for water-soil-air composite pollution of industrial park, relate to ecological environment remote sensing, industrial park pollution risk screening and geographic spatial information processing technical field.Firstly, core area and influence area are determined and unified spatial grid is constructed;After mapping multi-source remote sensing and spatial data to grid, respectively generate soil / surface, water body and atmospheric pollution risk indicator, and pollution source intensity index and migration channel index;According to the structure of "medium anomaly-source item support-migration condition", build comprehensive screening index and divide risk level;According to sample type, risk level, key migration node, background control, minimum distance and field constraint, generate sub-type priority sampling point.The present application unifies multi-source heterogeneous data to the same grid calculation, reduces single abnormal misjudgment risk by coupling pollution source intensity and migration channel condition, and directly converts screening result into executable sampling point layout scheme.

Description

Technical Field

[0001] This invention relates to the fields of ecological environment remote sensing, industrial park pollution risk screening and geospatial information processing, and in particular to a method for multi-source remote sensing collaborative screening, risk zoning and sampling point layout for water-soil-air complex pollution in industrial parks. Background Technology

[0002] Industrial parks have a variety of enterprises and diverse pollution emissions. Pollutants can form a complex pollution risk in water, soil and air through exhaust gas deposition, wastewater discharge, rainwater runoff, surface erosion, ditch transport and atmospheric diffusion.

[0003] Existing risk screening technologies have the following main shortcomings:

[0004] 1. Insufficient multi-media coordination. Existing methods mainly focus on monitoring a single pollutant medium or identifying a single remote sensing anomaly, making it difficult to coordinately express the spatial coupling relationship between water, soil, and air pollution. There is a lack of a unified indicator system and collaborative calculation method for water-soil-air complex pollution.

[0005] 2. Lack of a unified spatial grid computing platform. Multiple data sources, such as remote sensing imagery, DEM, wind fields, water systems, enterprise source data, and land use data, differ in spatial resolution, coordinate systems, and data formats. Existing methods mostly remain at the level of layer overlay or thematic mapping, making it difficult to support the gridded construction of comprehensive screening indices.

[0006] 3. Insufficient representation of pollution migration pathways. Existing risk assessments rely heavily on the intensity of pollution sources and the degree of media anomalies, and do not adequately consider whether pollutants can migrate through confluence paths, water system connectivity, ditch drainage, low-lying areas, and sedimentation in the prevailing and leeward directions.

[0007] 4. Insufficient integration between remote sensing screening and sampling point deployment. Existing methods struggle to directly translate remote sensing risk zoning results into actionable sampling point deployment plans. Risk maps often only display high-value areas without further consideration of sample type, risk level, migration nodes, background comparison, minimum spacing, and on-site accessibility, thus failing to generate a categorized candidate sampling point list that can be delivered to on-site personnel.

[0008] Therefore, there is an urgent need for a remote sensing collaborative screening, risk zoning, and sampling point layout method for water-soil-air complex pollution in industrial parks. This method would express multi-media remote sensing anomalies of water, soil, and atmosphere through a unified spatial grid, comprehensively consider the influence of pollution source intensity and migration channels, construct a comprehensive screening index, classify risk levels, and output categorized candidate sampling points that can be used for the preparation of on-site sampling plans. Summary of the Invention

[0009] This invention addresses the problems in existing industrial park remote sensing applications, such as fragmented results from multiple media, difficulty in direct calculation of different data scales, lack of inclusion of migration processes in risk zoning, and difficulty in converting remote sensing anomalies into candidate sampling points. It provides an executable and verifiable spatial screening and sampling point placement method.

[0010] To achieve the above objectives, the present invention adopts the following technical solution:

[0011] A multi-source remote sensing collaborative screening, risk zoning, and sampling point layout method for water-soil-air complex pollution in industrial parks is proposed. First, the core area and affected area are determined based on the industrial park boundary, concentrated industrial activity area, hydrological influence range, and prevailing wind direction, and a unified spatial grid is constructed. Second, multi-temporal surface remote sensing, air pollution remote sensing, DEM, wind field, water system, enterprise source terms, land use, and background reference data are mapped to the spatial grid. Then, soil / surface pollution risk indicators, water pollution risk indicators, air pollution transport risk indicators, pollution source intensity indicators, and pollution migration channel indicators are generated on the spatial grid. Based on this, effective indicators are normalized and weighted by confidence, and a comprehensive screening index is constructed according to the structure of "media anomaly—source term support—migration conditions," and risk levels are classified. Finally, candidate risk patches are extracted, and priority sampling points are generated based on sample type, risk level, key migration nodes, background control, minimum spacing, and field constraints.

[0012] Compared with the prior art, the present invention has the following beneficial effects:

[0013] 1. Reduce the overlay error of multi-source heterogeneous data. By unifying the spatial grid and using confidence weighting, data from different media and scales can be compared within the same computational unit, avoiding the misinterpretation of missing values ​​as low risk.

[0014] 2. Improve the interpretability of risk screening. The comprehensive screening index can be decomposed into the contributions of media anomalies, pollution source intensity, and migration channel conditions, making it easier to determine which factors drive the candidate area.

[0015] 3. Narrow the candidate sampling range and improve the representativeness of the sampling points. By coupling the "medium anomaly - source term support - migration conditions" and using directional migration constraints, the areas that need to be sampled are prioritized within a large park area; the candidate sampling points are classified to cover composite risk areas, single medium anomaly areas, migration channel nodes and background control areas, and the minimum spacing rule is used to reduce the clustering of sampling points.

[0016] 4. Establish a traceable and deliverable results system. The calculation process saves key parameters to make the results verifiable, and the output candidate sampling points include coordinates, priority, and rationale for setting them, which can be directly converted into a point table and task list for on-site sampling plans. Attached Figure Description

[0017] Figure 1 This is a flowchart of the overall process of the present invention;

[0018] Figure 2 A schematic diagram illustrating the unified gridding process for multi-source data;

[0019] Figure 3 A schematic diagram illustrating the remote sensing collaborative identification of water, soil, and air pollution.

[0020] Figure 4 This is a schematic diagram illustrating the coupling relationship between pollution source, media anomaly, and migration channel.

[0021] Figure 5 A schematic diagram illustrating the construction of the comprehensive screening index;

[0022] Figure 6 A schematic diagram illustrating the classification and output results of complex pollution risk levels;

[0023] Figure 7 This is a schematic diagram illustrating the conversion of risk zoning results into categorized candidate sampling points.

[0024] Figure 8 This is a schematic diagram of the spatial organization of the study area in the example embodiment;

[0025] Figure 9 This is a spatial partitioning diagram of the industrial park research area, industrial park source area, CORE, and BUFFER for an example implementation. The diagram shows the spatial relationship and functional division of STUDY_AREA, PARK, CORE, and BUFFER. Among them, STUDY_AREA is the example implementation research area formed by expanding PARK 3000m outward; CORE is the core influence area expanded 500m outward from PARK; and BUFFER is the outer verification area of ​​STUDY_AREA after deducting CORE.

[0026] Figure 10 This is a preliminary remote sensing map of land and water anomalies in an industrial park, illustrating the spatial relationships between Land_Anomaly, Water_Anomaly, and CPRI_wind critical patches. (a), (b), and (c) represent Land_Anomaly (land anomaly), Water_Anomaly (water anomaly), and CPRI_wind critical patches, respectively. Land and water anomalies are remote sensing risk indicators and are not directly equivalent to pollution concentration or exceedance conclusions. CPRI_wind critical patches are used to indicate priority areas for on-site verification.

[0027] Figure 11Example: A proxy map of pollution source intensity, migration pathways, and atmospheric transport risks in an industrial park. The map shows the spatial patterns of three proxy indicators: PSI_proxy, MTI_proxy, and ATR_wind. Among them, (a), (b), and (c) represent PSI_proxy, MTI_proxy, and ATR_wind, respectively. The three proxy indicators express the source term impact, migration pathway conditions, and wind-corrected atmospheric transport risks, respectively. The color scale is used to express the relative high and low levels within each indicator.

[0028] Figure 12 This example uses a risk zoning map of the CPRI_wind comprehensive screening index for an industrial park. The map shows the five comprehensive screening levels, the CORE / BUFFER boundary, and major medium- and high-risk candidate patches within the study area. Specifically, the area within the CORE with a CPRI_wind > 0.60 is 0.2437 km². 2 The value was 1.84% of the core; 0% within the buffer; there were no high-risk areas >0.80 in either the core or buffer; the results in the figure are used for on-site verification and sampling priority ranking.

[0029] Figure 13 This example illustrates the sampling priority index and categorized candidate sampling point map for an industrial park. The map shows the spatial relationship between the Sampling Priority Index and soil / surface dust, water / sediment, atmospheric deposition / road dust, perimeter checkpoints, and background control points. A total of 180 sampling points are included: 40 for soil / surface dust within the core, 40 for water / sediment, and 30 for atmospheric deposition / road dust; 40 for perimeter checkpoints within the buffer, and 30 for background controls.

[0030] Figure 14 This is a spatial distribution map of candidate areas for different types of industrial parks, as shown in the example. The map illustrates the spatial distribution of candidate areas for terrestrial sources, hydrological sources, atmospheric sources, peripheral verification, and background control. A, B, C, D, and E represent candidate areas for terrestrial sources, hydrological sources, atmospheric sources, peripheral verification, and background control, respectively. These candidate areas are used to organize on-site verification and sampling tasks and are not considered as pollution confirmation conclusions. Detailed Implementation

[0031] To facilitate understanding by those skilled in the art, the present invention will be further described below with reference to embodiments. The content mentioned in the embodiments is not intended to limit the present invention.

[0032] This invention aims to solve the following technical problems:

[0033] The first technical challenge is to convert surface remote sensing, atmospheric remote sensing, topography, water systems, wind fields, pollution sources, and background reference data with different spatial resolutions, temporal scales, and physical meanings into comparable indicators within the same grid cell, while preserving the validity and reliability of each indicator.

[0034] The second technical challenge is to avoid directly equating a single spectral anomaly with a confirmed pollution result. This invention defines remote sensing results for soil or surface, water bodies, and atmosphere as risk indicators, and reduces misjudgments caused by non-pollution factors by considering background reference, temporal differences, upstream and downstream relationships, downwind relationships, and proximity to pollution sources.

[0035] The third technical challenge is to incorporate the spatial connectivity between pollution sources, media anomalies, and migration channels into the same computational chain, so that high risk values ​​not only reflect "where the anomaly is," but also whether the anomaly has source support and migration conditions.

[0036] The fourth technical challenge is converting continuous risk indices and risk levels into directly executable sampling point deployment results. Based on risk zoning, this invention generates categorized candidate sampling points according to sample type, risk level, migration nodes, background control, minimum spacing, and on-site accessibility.

[0037] The fifth technical challenge is to establish a complete technical chain from remote sensing data input, indicator construction, risk zoning to candidate sampling point output, so that the deliverables can directly serve the preparation of on-site sampling plans.

[0038] The technical scope of this invention is regional-scale risk screening and candidate sampling point determination. Without legal monitoring or laboratory testing confirmation, the output results will not be presented as specific pollutant concentrations, environmental quality compliance status, or pollution exceedance conclusions.

[0039] To address this, this invention provides a multi-source remote sensing collaborative screening, risk zoning, and sampling point deployment method for complex water-soil-air pollution in industrial parks. First, it acquires surface remote sensing images, air pollution remote sensing products, topographic and hydrological data, wind fields, pollution sources, land use data, and background reference data of the industrial park and its affected areas, mapping them to a unified spatial grid. Second, it generates soil or surface pollution risk indicators, water pollution risk indicators, air pollution transport risk indicators, pollution source intensity indicators, and migration channel indicators, respectively. Then, it performs layered coupling according to the risk screening chain of "pollution source—medium anomaly—migration channel," obtaining the comprehensive screening index and risk level of each grid. Finally, it integrates composite risk, single-medium anomaly, key migration nodes, background comparison, spatial discretization, and on-site accessibility constraints to generate categorized priority sampling points. This method outputs pollution risk candidate areas and sampling priority areas, and does not replace statutory pollution concentration monitoring or exceedance determination; the technical process ends with the output of candidate sampling points.

[0040] It should be noted that this method does not collect or alarm on water, soil and air pollution monitoring values ​​in parallel, nor does it use multi-source remote sensing images to determine whether a pollution accident has occurred. Instead, it expresses remote sensing anomalies of soil or surface pollution, remote sensing anomalies of water pollution and air pollution transport risks in a unified manner, and further couples pollution source intensity and pollution migration channels to form a collaborative screening result that can be used for risk level classification, candidate area identification and priority sampling point placement.

[0041] The overall technical approach of this invention is as follows:

[0042] This invention uses a unified spatial grid as the computational carrier, multi-media remote sensing risk indicators as the observation layer, pollution sources and directional migration channels as the constraint layer, and composite risk zoning and categorized candidate sampling points as the output layer. The steps sequentially correspond to data unification, media anomaly extraction, source term and channel coupling, risk zoning, and sampling point generation.

[0043] The method of the present invention will now be described in detail.

[0044] Step S1: Define the research scope and construct a unified spatial grid

[0045] Step S1 includes: determining the core area and the area of ​​influence based on the boundaries of the industrial park, the concentrated area of ​​industrial activities, the scope of hydrological influence, and the prevailing wind direction, and constructing a unified spatial grid.

[0046] Specifically as follows:

[0047] The research scope of the industrial park was determined, the core area and buffer zone were divided, a unified spatial grid was established, and basic spatial attributes were assigned to each grid unit.

[0048] First, acquire spatial data related to the sampling points, including the boundaries of industrial parks, key enterprises, industrial land areas, park roads, rivers and waterways, ditches, ponds, reservoirs, sewage outlets, storage yards, hazardous waste temporary storage areas, sewage treatment facilities, accident pools, land use, background reference areas, and other spatial data.

[0049] The core area is defined as the boundary of the industrial park or the area with a concentrated distribution of industrial land. Based on the industrial park's potential hydrological migration distance, atmospheric transport distance, surrounding topographic and hydrological conditions, and sampling point requirements, one or more buffer zones are established around the core area to form the research scope. These buffer zones can be determined using fixed-distance buffering, buffering extending along river systems, buffering extending along the prevailing wind direction, or heterogeneous buffering methods combined with topographic confluence directions.

[0050] In one embodiment, the scope of the study includes:

[0051] Core area: Determined by the boundaries of industrial parks, key enterprises, or areas with concentrated industrial land, it is used to characterize the concentrated distribution range of major pollution sources;

[0052] Near-field impact zone: Formed by the outward expansion of the core area, used to characterize the impact of pollutants on surface runoff, ditch transport, dust deposition, or near-field atmospheric diffusion.

[0053] Far-field influence zone: Extends outward from the near-field influence zone, or extends along downstream water systems, prevailing wind directions, and low-lying collection areas. It is used to characterize the peripheral screening range where pollutants may migrate to the outside of the park, and to provide a spatial basis for the deployment of background control points and peripheral verification points.

[0054] After defining the research scope, the research scope is divided into regular spatial grids. These spatial grids can be square, rectangular, hexagonal, or other regular grids, with square grids being preferred as they match the spatial resolution of the remote sensing imagery and the environmental management scale. The grid scale can be determined based on the spatial resolution of the data source, the area of ​​the study area, the density of pollution sources, and the accuracy requirements of risk zoning; for example, it can be set to 10 m × 10 m, 20 m × 20 m, 30 m × 30 m, 50 m × 50 m, 100 m × 100 m, or other scales.

[0055] In a preferred embodiment, when Sentinel-2, GF, or hyperspectral remote sensing images are mainly used to identify surface and water anomalies, the spatial grid scale can be set to 10 m × 10 m or 30 m × 30 m; when the study area is large and mainly used for macro-risk zoning, the spatial grid scale can be set to 50 m × 50 m or 100 m × 100 m; when it is necessary to match the accuracy of management units, enterprise plots, sewage outlets, or sampling points, the grid size can be adjusted according to the actual management scale.

[0056] Each spatial grid, denoted as Gi, serves as the basic unit for subsequent calculations, where i is the grid number. A grid attribute table is established for each grid Gi, including at least the grid number, grid center coordinates, region type, land use type, distance to the nearest enterprise boundary, distance to the nearest sewage outlet, distance to the nearest storage yard, distance to the nearest water body, distance to the nearest surrounding environmental object, slope, aspect, runoff accumulation, prevailing wind direction influence, and whether it is located in the core area, near-field influence area, or far-field influence area.

[0057] Through the above processing, multi-source remote sensing imagery, air pollution products, DEM topographic data, wind field data, water system data, enterprise spatial data, land use data, and background reference data can be uniformly mapped to the same grid cell. Each grid cell can simultaneously obtain soil / surface pollution risk indicators, water pollution risk indicators, air pollution transport risk, pollution source intensity, pollution migration channel impact, and background reference and sampling constraint information, thus providing a unified spatial basis for subsequent comprehensive screening index construction and candidate sampling point determination.

[0058] Step S2: Acquire and preprocess multi-source data

[0059] Step S2 includes: acquiring multi-temporal surface remote sensing, atmospheric pollution remote sensing, DEM, wind field, water system, enterprise source terms, land use and background reference data, and completing quality control, projection unification, time window matching and grid mapping.

[0060] Specifically as follows:

[0061] The process involves acquiring multi-source remote sensing data, acquiring auxiliary geospatial data, performing remote sensing image preprocessing, performing auxiliary data standardization, and mapping all data to a unified spatial grid.

[0062] 2.1 Remote Sensing Data Acquisition

[0063] The remote sensing data includes multi-source remote sensing data used to identify remote sensing anomalies in soil pollution, water pollution, and atmospheric pollution transport risks. The multi-source remote sensing data includes, but is not limited to, the following data types:

[0064] First, multispectral remote sensing imagery. The multispectral remote sensing imagery may include Sentinel-2, Landsat, GF series, or other remote sensing imagery with visible light, near-infrared, and shortwave infrared bands, used to extract information on bare land, vegetation, water bodies, built-up land, storage yards, roads, riverbanks, water color anomalies, turbidity anomalies, and surface disturbance.

[0065] Second, hyperspectral remote sensing imagery. The hyperspectral remote sensing imagery may include GF-5, ZY-1, EnMAP, PRISMA, EMIT, airborne hyperspectral, UAV hyperspectral, or other remote sensing imagery with continuous narrow-band spectral information, used to extract soil spectral anomalies, mineral or heavy metal-related spectral responses, water body spectral anomalies, vegetation stress characteristics, and differences in surface materials in industrial parks.

[0066] Third, remote sensing products for air pollution. These remote sensing products may include NO2, SO2, AOD, aerosols, CO, or other remote sensing products that can reflect industrial emissions and the transport of air pollution, used to identify air pollution anomalies and transport risks in and around industrial parks.

[0067] Fourth, thermal infrared or nighttime light remote sensing data. This thermal infrared or nighttime light remote sensing data can be used to assist in identifying the intensity of industrial activity, abnormal heat sources, changes in production activities, nighttime emissions, or the level of activity on industrial land.

[0068] Fifth, UAV remote sensing data. The UAV remote sensing data may include visible light, multispectral, hyperspectral, thermal infrared, or oblique photography data, which are used for high-resolution verification of high-risk areas, pollution hotspots, and sampling sites identified by satellite remote sensing.

[0069] Sixth, near-surface auxiliary data. This data may include UAV imagery, historical patrol photos, existing ground feature records, GNSS positioning records, road traffic conditions, and safety constraint information, etc., used to assist in identifying ground feature types, interpreting remote sensing anomalies, and determining the accessibility of candidate sampling points. The aforementioned near-surface auxiliary data can serve as background reference and sampling constraint data, but does not constitute necessary input for generating candidate sampling points in this method.

[0070] The aforementioned remote sensing data can be used individually or in combination. Differences in spatial resolution, temporal resolution, spectral resolution, and imaging time from different data sources can be uniformly expressed through resampling, time-matching, image fusion, normalization, or hierarchical calculation.

[0071] 2.2 Auxiliary Geospatial Data Acquisition

[0072] The auxiliary geospatial data includes spatial data used to represent pollution sources, pollution migration pathways, and surrounding environmental objects. Auxiliary geospatial data includes, but is not limited to:

[0073] First, DEM topographic data is used to calculate slope, aspect, surface runoff direction, runoff accumulation, runoff path, low-lying catchment areas, and hydrological connectivity.

[0074] Second, wind field and meteorological data, including wind direction, wind speed, prevailing wind direction, rainfall, temperature, humidity, atmospheric stability or other meteorological factors, are used to analyze the direction of air pollutant transport, diffusion intensity and deposition impact.

[0075] Third, water system data such as rivers, ditches, sewage ditches, stormwater pipe networks, ponds, reservoirs, lakes and wetlands are used to analyze water pollution anomalies, upstream and downstream relationships, hydrological migration paths and water body sensitivity;

[0076] Fourth, enterprise spatial distribution data, including enterprise boundaries, enterprise type, industry category, production scale, sewage outlet location, storage yard location, hazardous waste temporary storage point, sewage treatment facilities, accident pool, loading and unloading area, storage area and historical pollution records, are used to construct pollution source intensity indicators.

[0077] Fifth, land use data, including industrial land, construction land, bare land, cultivated land, forest land, water bodies, roads, residential land and unused land, are used to assist in identifying exposed or disturbed land surfaces, background control areas and risk zone boundaries;

[0078] Sixth, background reference and sampling constraint data, including residential areas, farmland, water bodies, ecological protection objects, roads and accessibility constraints, are used to determine the research scope, background control points, peripheral verification points and sampling accessibility constraints;

[0079] Seventh, existing regulatory and historical data, including historical monitoring data, enterprise emission regulatory data, historical law enforcement records, drone patrol records, historical image interpretation records, and on-site photos, can be used as background references, auxiliary data for setting sample point types, and sampling constraints.

[0080] 2.3 Remote Sensing Image Preprocessing

[0081] The acquired remote sensing images are preprocessed to make them suitable for participation in unified spatial grid computing. The preprocessing includes, but is not limited to, radiometric calibration, atmospheric correction, geometric correction, orthorectification, image cropping, cloud and shadow removal, stripe noise removal, bad band removal, band registration, resampling, unified projection, and image mosaicking.

[0082] For multispectral remote sensing images, surface reflectance data is obtained after atmospheric correction, and pixels that are not suitable for soil or water anomaly identification are removed based on cloud masking, shadow masking, water body masking, vegetation masking, and building shadow masking.

[0083] For hyperspectral remote sensing images, in addition to radiometric calibration, atmospheric correction and geometric correction, bad band removal, water vapor absorption band removal, noise band removal, spectral smoothing, continuous spectrum removal, spectral derivative transformation, standardization or other spectral preprocessing can be performed to enhance the ability to identify soil and water pollution anomalies.

[0084] For remote sensing products of air pollution, the data is cropped according to the study area and screened based on product resolution, time scale, and quality control indicators. For air pollution products with low spatial resolution, spatial interpolation, area weighting, nearest neighbor assignment, or statistical summarization can be used to map them to a unified spatial grid, and the air pollution transport risk can be further calculated by combining wind field data and the spatial distribution of pollution sources.

[0085] 2.4 Standardization of Auxiliary Data

[0086] The auxiliary geospatial data undergoes coordinate system unification, topology checking, spatial clipping, attribute standardization, outlier checking, and spatial gridding.

[0087] For point pollution source data, calculate its distance, density, and impact range with each grid cell; for linear water systems, ditches, sewage channels, and roads, calculate their distance, connectivity, and upstream / downstream relationships with each grid cell; for area data such as enterprise boundaries, industrial land, storage yards, residential areas, farmland, and water sources, calculate the intersection, proximity, and coverage ratio between each grid cell and the area object.

[0088] For DEM data, calculate slope, aspect, surface runoff direction, runoff accumulation, runoff path, low-lying catchment area, and hydrological connectivity path, and assign the above hydro-topographic factors to the corresponding grid cells.

[0089] For wind field data, the prevailing wind direction, wind speed frequency, and wind direction stability during the statistical study period are analyzed, and the downwind influence relationship is established based on the azimuth relationship between pollution sources and grid cells. For multi-temporal wind field data, wind direction weights can be calculated by day, month, season, or pollution event time window.

[0090] For background reference and sampling constraint data, the distance from the grid to residential areas, farmland, water bodies, ecological protection objects, roads, or accessible areas is calculated, and their suitability as background control points, peripheral verification points, or accessibility-restricted points is marked. This information is not a necessary factor in the CPRI calculation but is used for defining the study scope and selecting candidate sampling points.

[0091] 2.5 Data is uniformly mapped to spatial grids

[0092] After completing the preprocessing of remote sensing images and auxiliary data, all data are uniformly mapped to the spatial grid constructed in step S1. For raster data, the data are summarized according to the mean, maximum, median, standard deviation, proportion of outlier pixels, or other statistical measures within the grid area; for vector data, values ​​are assigned according to the intersection area between the grid and vector objects, the number of objects, distance, density, or proximity relationship; for time series data, the data are summarized according to the mean, maximum, cumulative value, outlier frequency, rate of change, or time-weighted value within the study period.

[0093] Through this step, each spatial grid Gi obtains the following basic data fields:

[0094] Remote sensing image fields are used to calculate soil / surface pollution risk indicators and water pollution risk indicators;

[0095] The air pollution field is used to calculate air pollution transport risk indicators;

[0096] The pollution source field is used to calculate the pollution source intensity index;

[0097] Topographic and hydrological fields and wind field fields are used to calculate pollution migration pathway indicators;

[0098] The land use field, background reference field, and sampling constraint field are used to determine the background comparison area, the peripheral verification area, and the sampling accessibility conditions.

[0099] Historical monitoring or existing field record fields are used to assist in setting background references and sampling constraints;

[0100] Near-ground auxiliary fields are used to record existing UAV images, historical photos, positioning information, road conditions and safety constraints within or near each grid cell, and are used for remote sensing anomaly interpretation, background reference setting and accessibility determination of candidate sampling points.

[0101] The aforementioned unified mapping process enables data from different sources, scales, and types to be comprehensively calculated within the same spatial grid. This provides a data foundation for the subsequent construction of soil / surface pollution risk indicators, water pollution risk indicators, air pollution transport risk indicators, pollution source intensity indicators, pollution migration channel indicators, and comprehensive screening indices. It also provides background references and on-site constraints for the selection of candidate sampling points.

[0102] Step S3: Construct a soil / surface pollution risk indicator

[0103] Step S3 includes generating the Soil / Surface Pollution Risk Indicator (TSAI) for bare soil, stockpiles, road dust deposition areas, and riverbanks.

[0104] Specifically as follows:

[0105] Step S3 is used to construct a soil / surface pollution risk indicator, denoted as TSAI, within a unified spatial grid. TSAI stands for Soil Anomaly Index. The TSAI is used to characterize the degree of soil / surface pollution risk indication that may exist in bare soil areas, around storage yards, around factory boundaries, road dust deposition areas, around ditches, riverbanks, low-lying catchment areas, around accident ponds, around hazardous waste temporary storage areas, and other exposed surface areas within industrial parks and their surrounding impact zones.

[0106] Step S3 specifically includes soil candidate region identification, removal of interfering features, extraction of soil spectral features, determination of background reference areas, calculation of soil outliers, and gridding assignment.

[0107] 3.1 Soil Candidate Region Identification

[0108] First, based on multispectral or hyperspectral remote sensing images, land use data, industrial land boundaries, enterprise boundaries, roads, storage yards, riverbanks, ditches, and historically disturbed areas, candidate areas that may be exposed to soil or surface sediments are identified.

[0109] The soil candidate regions include, but are not limited to:

[0110] 1. Exposed ground, construction-disturbed ground, and unpaved surfaces within the industrial park;

[0111] 2. The area surrounding the enterprise's factory boundary, storage yard, temporary hazardous waste storage site, and accident pool;

[0112] 3. Dust accumulation areas on both sides of the road, loading and unloading areas, and the area surrounding vehicle transport routes;

[0113] 4. The area surrounding rivers, ditches, sewage ditches, rainwater runoff channels, and low-lying areas;

[0114] 5. Farmland edges, riverbanks, abandoned land, fallow land, and other surface areas that may receive pollutant deposition or runoff input.

[0115] In one implementation, the study area can be preliminarily classified using Normalized Difference Vegetation Index (NDVI), Water Index (MNDWI), Building Index (NDBI), brightness characteristics, shortwave infrared characteristics, and land use data. Water bodies, dense vegetation, building shadows, cloud shadows, and areas that are clearly unsuitable for soil spectral analysis can be removed to obtain a soil candidate area mask.

[0116] The soil candidate region mask can be represented as:

[0117] This indicates that the pixel or grid belongs to the soil candidate region;

[0118] This indicates that the pixel or grid does not belong to the soil candidate area.

[0119] In a preferred embodiment, a pixel or grid is identified as a candidate soil region when it meets the following conditions: NDVI is less than the vegetation threshold, MNDWI is less than the water body threshold, and it does not belong to building shadows, cloud shadows, or deep water areas. The thresholds can be determined based on image features of the study area, land use data, historical data, or sample training results.

[0120] 3.2 Removal of Interfering Features

[0121] Because industrial parks contain buildings, roads, roofs, paved surfaces, shadows, water bodies, dense vegetation, plastic mulch, metal roofs, and other artificial features, these features may interfere with the identification of soil pollution anomalies by remote sensing. Therefore, further screening of candidate soil areas is necessary.

[0122] In one implementation, interfering features can be removed by combining the following rules:

[0123] 1. Use NDVI to remove areas with dense vegetation;

[0124] 2. Use MNDWI or NDWI to remove water bodies;

[0125] 3. Use shadow index, brightness threshold, or solar altitude angle information to remove building shadows and cloud shadow areas;

[0126] 4. Use NDBI, land use data, or building vector data to remove buildings and large areas of hardened ground;

[0127] 5. Use road vector data, image texture features, or linear feature recognition results to eliminate large areas of asphalt and cement pavement;

[0128] 6. Correct misjudged areas using existing field data, UAV imagery, or high-resolution remote sensing imagery.

[0129] After removing interfering features, the effective soil area for calculating soil / surface pollution risk indicators is obtained.

[0130] 3.3 Extraction of Soil Spectral Features

[0131] Spectral features are extracted from multispectral or hyperspectral remote sensing images within the effective soil area. The soil spectral features include, but are not limited to, the following types:

[0132] First, the original band reflectance characteristics, including visible light, near infrared, red edge, shortwave infrared, or other band reflectances that can be used to characterize soil color, iron oxides, organic matter, water content, mineral composition, and differences in surface sediments.

[0133] Second, band ratio characteristics, including ratios, differences, normalized differences, or combined operation characteristics between different bands, are used to enhance the spectral differences of pollution-related land surfaces.

[0134] Third, normalized index characteristics, including remote sensing indices related to bare soil, iron oxides, clay minerals, organic matter, salinity, moisture, or vegetation stress.

[0135] Fourth, red-edge and near-infrared features are used to help identify sparsely vegetated areas, vegetation-stressed areas, and vegetation-soil mixed areas affected by pollution deposition.

[0136] Fifth, shortwave infrared characteristics are used to characterize the spectral response related to soil mineral composition, moisture, clay minerals, carbonates, sulfates, or industrial stockpile materials.

[0137] Sixth, hyperspectral continuous spectrum removal characteristics, including absorption peak position, absorption depth, absorption area, absorption width, symmetry, and band characteristics after continuous spectrum removal.

[0138] Seventh, spectral derivative features, including first derivative, second derivative, smoothing derivative, or other spectral transformation features used to enhance weak absorption features and reduce background effects.

[0139] Eighth, texture and spatial neighborhood features, including local mean, standard deviation, coefficient of variation, texture contrast, edge features, and anomalous patch neighborhood features, are used to express the spatial heterogeneity of storage yards, bare land, road dust deposition areas, and disturbed surfaces.

[0140] For each grid The soil spectral features within the effective soil area are extracted, and a soil spectral feature vector is formed:

[0141]

[0142] in, Represents a grid Soil spectral feature vectors, to It represents different bands, exponents, ratios, continuous spectrum features, derivative features, or spatial texture features.

[0143] 3.4 Determination of Background Reference Area

[0144] To reduce the impact of non-polluting factors such as soil parent material, surface moisture, vegetation cover, imaging conditions, and terrain shadows on remote sensing anomaly identification, this invention sets up a background reference area for calculating the relative anomaly degree.

[0145] The background reference area can be selected from within the research area or the area surrounding the research area, preferably meeting the following conditions:

[0146] 1. Located outside the industrial park or far from major pollution sources;

[0147] 2. The land use type, soil type, topographic conditions, or surface cover characteristics are comparable to the area to be evaluated;

[0148] 3. Not located in a clear downstream confluence path, downwind settling path, or known pollution impact area;

[0149] 4. The remote sensing image quality is good and is not severely disturbed by clouds, shadows, water bodies, buildings or dense vegetation;

[0150] 5. The pollution risk is considered low based on historical monitoring, existing on-site data, land use data, or expert judgment.

[0151] In one implementation, the background reference area can be determined through manual selection, rule-based screening, cluster analysis, low-risk area screening, or stepwise removal of outliers. For industrial parks with significant topographic, water system, or wind direction influences, the background reference area should preferably be selected from upstream, upwind, or similar surface areas far from pollution sources.

[0152] 3.5 Calculation of Soil Remote Sensing Anomalies

[0153] Based on the spectral characteristics of soil candidate areas and the statistics of background reference areas, a grid is calculated. The Soil / Surface Pollution Risk Indicator (TSAI). For the k-th feature, first determine the direction in which the feature increases with risk. ,in Take 1 or -1.

[0154] When a feature has a clear risk direction, a one-way standardized anomaly is used: When only the degree of deviation from the background is judged without pre-setting the direction, a two-way anomaly is used: .

[0155] in, For the k-th feature of grid Gi, and These are the mean and standard deviation of the features corresponding to the background reference region, respectively; when When the value is 0, this feature is not included in the current batch calculation.

[0156] Integrate effective features into : ,in The weights are non-negative and apply to the effective features involved in the calculation. .

[0157] As a subordinate implementation method, Mahalanobis distance can be used. ,in The background feature covariance matrix; when When irreversibility is not possible, regularize the covariance matrix or remove collinear features.

[0158] Used to express the relative screening priority of soil / surface anomalies, spectral anomalies are not presented as pollutant concentration inversion results; if existing historical data exists, it can be used as auxiliary information to interpret remote sensing anomalies and set background references.

[0159] 3.6 Normalization of Soil Anomaly Indicators

[0160] In order to To enable comprehensive calculations in conjunction with subsequent water pollution risk indicators, air pollution transport risk indicators, pollution source intensity indicators, and pollution migration channel indicators, it is necessary to... Normalize to the 0-1 interval.

[0161] In one implementation, a minimum-maximum normalization method can be used:

[0162]

[0163] When the upper limit of normalization equals the lower limit, the indicator is not included in the fusion in the current batch. When comparing multiple periods, the same reference period and fixed normalization parameters are used for each period; the results of independent normalization for each period are only used for internal ranking within each period.

[0164] in, This represents the soil / surface pollution risk indicator after grid Gi normalization. and Representing the research scope The minimum and maximum values.

[0165] As a subordinate implementation method, quantile normalization, robust normalization, standard deviation normalization, or logistic function normalization methods can be used to reduce the impact of abnormal terrain values, extreme confluence values, abnormal wind field periods, or artificial ditch data errors on the results.

[0166] Normalized The closer a value is to 1, the higher the risk level of soil / surface pollution in that grid; the closer a value is to 0, the lower the risk level of soil / surface pollution in that grid.

[0167] 3.7 Spatial Expression of Soil Anomalies

[0168] The normalized TSAI values ​​are assigned to a unified spatial grid, and remote sensing anomaly maps of soil or surface pollution are generated. For grids without effective soil pixels, whether to assign values ​​can be determined based on the proportion of effective soil pixels in their neighborhood, land use type, proximity of pollution sources, or existing field data; alternatively, they can be marked as having no effective soil observations and not directly involved in TSAI calculation.

[0169] In one implementation, soil pollution remote sensing anomalies can be classified into low anomalies, low-to-medium anomalies, medium anomalies, medium-to-high anomalies, and high anomalies based on the magnitude of the TSAI value. The classification can be performed using the natural breakpoint method, quantile method, standard deviation method, or fixed threshold method.

[0170] Through the above steps, the Soil / Surface Pollution Risk Indicator (TSAI) for each spatial grid can be obtained, and it will be used as one of the components of the subsequent Comprehensive Screening Index (CPRI). TSAI is not used alone as the final determination of pollution risk in industrial parks, but rather participates in the comprehensive screening and evaluation together with the Water Pollution Risk Indicator, Air Pollution Transport Risk Indicator, Pollution Source Intensity Indicator, and Pollution Migration Channel Indicator.

[0171] Step S4: Construct a water pollution risk indicator

[0172] Step S4 includes generating a Water Pollution Risk Indicator (WAI) for rivers, ditches, sewage channels, ponds, and downstream receiving water bodies.

[0173] Specifically as follows:

[0174] Step S4 is used to construct a water pollution risk indicator, denoted as WAI, within a unified spatial grid. WAI stands for Water Anomaly Index. The WAI is used to characterize the degree of anomalies in water color, turbidity, suspended solids, organic pollution, eutrophication, or industrial emissions that may exist in rivers, ditches, sewage channels, stormwater runoff channels, ponds, reservoirs, emergency ponds, downstream water bodies, and other surface water bodies within the industrial park and its surrounding impact area.

[0175] Step S4 specifically includes water body identification, water body object classification, water body spectral feature extraction, background water body or reference water segment determination, upstream and downstream difference analysis, water body anomaly value calculation, and gridded value assignment.

[0176] 4.1 Water Body Identification

[0177] First, based on multispectral or hyperspectral remote sensing images, water system vector data, land use data, DEM confluence paths, sewage channels, ditches, ponds, reservoirs, and existing field data, surface water bodies within the study area are identified.

[0178] The surface water bodies include, but are not limited to:

[0179] 1. Rivers and their tributaries within or around the park;

[0180] 2. Drainage ditches, rainwater ditches, sewage ditches, and open channels in industrial parks;

[0181] 3. Ponds, emergency ponds, sedimentation tanks, water storage tanks, and temporary wastewater storage bodies surrounding the enterprise's factory boundary;

[0182] 4. Downstream river sections, confluences, river bends, low-velocity water bodies, and low-lying catchment areas within the park;

[0183] 5. Ponds, reservoirs, wetlands or other surface water bodies that may receive inputs of industrial wastewater, rainwater runoff, surface erosion and leached material from storage sites.

[0184] In one implementation, a water mask can be constructed using the Normalized Water Index (NDWI), the Improved Normalized Water Index (MNDWI), the Automatic Water Extraction Index (AWEI), near-infrared and short-wave infrared thresholds, topographic depressions, river system vector data, and artificial correction results.

[0185] A water mask can be represented as:

[0186] This indicates that the pixel or grid belongs to a water area;

[0187] This indicates that the pixel or grid does not belong to the water area.

[0188] In a preferred embodiment, for narrow ditches, sewage ditches, and rainwater ditches whose width is less than the spatial resolution of remote sensing images or which are severely mixed with surrounding buildings, roads, and vegetation, corrections can be made by combining high-resolution images, UAV images, water system vector data, DEM confluence paths, and existing field data to reduce the omission of small water bodies.

[0189] 4.2 Classification of Water Bodies

[0190] Because the types of water bodies in industrial parks are complex and the pollution risks of different water bodies are different, the water bodies are classified after identification.

[0191] In one implementation, water bodies can be classified into the following types:

[0192] The first category consists of the main channels and major tributaries of rivers, which are used to characterize the main hydrological channels through which pollutants migrate to the outside of the park.

[0193] The second category includes internal ditches, rainwater ditches, sewage ditches, and open channels within the park, used to characterize short-distance transport channels for pollutants within the park.

[0194] The third category includes ponds, accident ponds, sedimentation ponds, water storage ponds, and low-lying water accumulation areas around enterprises, which are used to characterize the risk of local pollution accumulation or temporary storage.

[0195] The fourth category includes downstream receiving water bodies, water bodies near the inlet, and water bodies in the surrounding impact zone, used to characterize the potential impact of pollutants migrating along the water system to the outside of the park;

[0196] The fifth category is background water bodies or reference water sections, used to calculate the relative changes in remote sensing anomalies of water pollution.

[0197] Different importance weights can be assigned to different water body types. Water bodies that are close to sewage outlets, ditches, enterprise boundaries, storage yards, accident ponds, downstream river sections, or water bodies in the surrounding impact zone can be assigned a higher risk concern weight.

[0198] 4.3 Extraction of Spectral Features of Water Bodies

[0199] Water spectral features are extracted from multispectral or hyperspectral remote sensing images of water bodies. These water spectral features include, but are not limited to, the following types:

[0200] First, water reflectance characteristics, including blue light, green light, red light, red edge, near infrared, shortwave infrared, or other bands of reflectance that can be used to characterize water color, suspended matter, turbidity, chlorophyll, colored dissolved organic matter, and water pollution anomalies.

[0201] Second, water body index characteristics, including NDWI, MNDWI, AWEI, normalized turbidity index, suspended matter related index, chlorophyll related index, water color index, red-green ratio, red-blue ratio, green-blue ratio, or other water quality remote sensing indices.

[0202] Third, the characteristics of water brightness and turbidity, including total visible light brightness, enhanced red light reflectance, abnormal near-infrared reflectance, abnormal short-wave infrared response, and enhanced water turbidity.

[0203] Fourth, abnormal water color characteristics, including the spectral expression of water color changes from normal dark, blue-green or naturally turbid states to abnormal yellow, reddish-brown, grayish-white, black and odorous water hues or other abnormal water color changes.

[0204] Fifth, hyperspectral water body characteristics, including specific absorption peaks, reflection peaks, continuous spectrum removal characteristics, spectral derivative characteristics, chlorophyll absorption characteristics, suspended matter scattering characteristics, and spectral anomalies under the influence of industrial pollution emissions.

[0205] Sixth, spatial neighborhood and morphological characteristics, including the area of ​​abnormal patches in the water body, the length of abnormal patches, the proportion of abnormal pixels, the continuous abnormal length along the river, the location of local abnormal abrupt changes, the abnormal gradient near the sewage outlet, and the abnormal differences between upstream and downstream.

[0206] For each grid Gi, extract the water spectral features within its water region and form a water spectral feature vector:

[0207]

[0208] in, This represents the water spectral eigenvector of grid Gi. to It represents different spectral bands, indices, ratios, water color, turbidity, suspended matter, chlorophyll, hyperspectral features, or spatial neighborhood features.

[0209] 4.4 Determination of background water body or reference water section

[0210] To reduce the impact of non-polluting factors such as seasonal changes, rainfall processes, water level changes, natural sediment, algae growth, solar altitude angle, shadows, and sensor differences on the identification of water anomalies, this invention sets up a background water body or a reference water segment to calculate relative water anomalies.

[0211] The background water body or reference water section can be determined in the following ways:

[0212] 1. Select the water body upstream of the industrial park as the background water section;

[0213] 2. Select water bodies that are far away from sewage outlets, enterprise boundaries, storage yards, ditch inlets, and industrial activity areas as background water bodies;

[0214] 3. Select a stable upstream section of the same river as a reference for downstream anomaly analysis;

[0215] 4. Select water bodies with relatively stable water quality from historical multi-temporal images as the temporal background;

[0216] 5. Determine the background water body by combining upstream water bodies, surrounding reference water bodies, historical data, or expert judgment;

[0217] 6. After clustering or removing outliers from the water body characteristics of the entire study area, low-anomaly water bodies are used as background references.

[0218] For river-type water bodies, it is preferable to use the upstream section of the same water system as a reference to reduce natural differences between different water body types. For ponds, sedimentation tanks, emergency ponds, or isolated water bodies, images of similar water bodies or historical normal conditions can be used as a reference.

[0219] 4.5 Analysis of upstream and downstream differences and anomalies of adjacent sewage sources

[0220] Water pollution in industrial parks exhibits distinct spatial directionality and localized abrupt changes. Pollutants typically migrate downstream along ditches, sewage channels, stormwater runoff channels, river systems, and low-lying terrain. Therefore, this invention incorporates upstream-downstream differences and anomaly analysis of nearby pollution sources when calculating the water pollution index (WAI).

[0221] In one implementation, for rivers, ditches, or sewage channels, the upstream reference grid and downstream influence grid of each water body grid or water body object are determined based on the DEM confluence direction, water system flow direction, river topology, or artificially corrected flow direction.

[0222] For grid Gi, the differences in water anomalies between it and the upstream reference area can be calculated:

[0223]

[0224] in, This represents the water body spectral characteristics, water quality remote sensing index, or comprehensive water body anomaly characteristics of the grid Gi. This represents the corresponding feature value of its upstream reference region. This indicates differences between upstream and downstream sectors.

[0225] when If the value exceeds a preset threshold and the grid is located downstream of a sewage outlet, enterprise boundary, ditch inlet, stockpile runoff path, or rainwater runoff channel, the grid can be determined to have a high probability of remote sensing anomalies related to water pollution.

[0226] As a subordinate implementation method, anomaly features of nearby sewage sources can be constructed for water bodies near sewage outlets, enterprise effluent outlets, ditch inlets, emergency pool outlets, or suspected discharge points. These features comprehensively consider the intensity of water body anomalies, distance from the sewage source, and spatial relationship along the water flow direction. Grids that are closer to the sewage source, located downstream, and exhibit significant spectral anomalies are assigned higher water pollution anomaly values.

[0227] 4.6 Calculation of Remote Sensing Anomalies in Water Bodies

[0228] Based on water spectral characteristics, background water bodies, upstream and downstream differences, and proximity of pollution sources, a computational grid is calculated. Water pollution risk indicator (WAI) is used. Turbidity, water color, suspended solids, chlorophyll, upstream and downstream differences, and pollution source proximity are first converted to the [0, 1] interval to avoid direct addition of different dimensions.

[0229] For the k-th water body characteristic with a clear risk direction, calculate When only judging the degree of deviation from the background, calculate .

[0230] in, As the characteristic direction, For the water features of the grid Gi, and These are the mean and standard deviation of the background water body or the upstream reference water section, respectively; when When the value is 0, this feature is not included in the current batch calculation.

[0231] Preferably, All components have been normalized. to The weights are non-negative and their sum is 1; Represents a grid Normalized values ​​of water turbidity anomalies; This represents the normalized value of water color anomaly; This represents the normalized value of suspended solids anomalies; This represents the normalized value of chlorophyll a or algal abnormalities. The difference between upstream and downstream operations after normalization. The normalized impact of the pollution source's proximity.

[0232] WAI is used to express the priority of water body anomaly screening, and is not directly expressed as water quality level or specific pollutant concentration; if existing historical data exists, it can be used as auxiliary information to explain water body anomalies and set background references.

[0233] 4.7 Normalization of water body anomaly indicators

[0234] In order to enable WAI to be comprehensively calculated with soil / surface pollution risk indicators, air pollution transport risk indicators, pollution source intensity indicators, and pollution migration channel indicators, WAI needs to be normalized to the 0-1 range.

[0235] In one implementation, a minimum-maximum normalization method can be used:

[0236]

[0237] When the upper limit of normalization equals the lower limit, the indicator is not included in the fusion in the current batch. When comparing multiple periods, the same reference period and fixed normalization parameters are used for each period; the results of independent normalization for each period are only used for internal ranking within each period.

[0238] in, Represents a grid Normalized water pollution risk indicator and Representing the research scope The minimum and maximum values.

[0239] As a subordinate implementation method, quantile normalization, robust normalization, standard deviation normalization, or logistic function normalization methods can be used to reduce the impact of abnormal terrain values, extreme confluence values, abnormal wind field periods, or artificial ditch data errors on the results.

[0240] Normalized The closer a value is to 1, the higher the risk level of water pollution in that grid; the closer a value is to 0, the lower the risk level of water pollution in that grid.

[0241] 4.8 Spatial Expression of Water Body Anomalies

[0242] Normalized The values ​​are assigned to a unified spatial grid, and a remote sensing anomaly map of water pollution is generated. For grids containing water body pixels or intersecting with water bodies, calculations can be performed directly. For grids that are adjacent to water bodies but do not contain water body pixels, the determination of whether to assign a water body influence value can be made based on their distance from the water body, the water body anomaly level, upstream and downstream relationships, proximity to sewage sources, and confluence path relationships.

[0243] For narrow ditches, sewage ditches, or storm drains whose width is less than the spatial resolution of remote sensing images, supplementary values ​​can be assigned using water system vector line buffers, UAV imagery, or existing water body records. For water bodies obscured by clouds, shadows, bridges, buildings, or dense vegetation, corrections can be made or the water bodies can be marked as having no valid observations based on adjacent water sections, historical imagery, or existing data.

[0244] In one implementation, it can be based on The numerical value of the anomalies is used to classify remote sensing anomalies of water pollution into low anomalies, low-medium anomalies, medium anomalies, medium-high anomalies, and high anomalies. The classification can be performed using the natural breakpoint method, quantile method, standard deviation method, fixed threshold method, or upstream-downstream abrupt change threshold method.

[0245] By following the steps described above, the water pollution risk indicator for each spatial grid can be obtained. and use it as a subsequent comprehensive screening index. One of the components of the index. It is not used as the sole final determination of water pollution in industrial parks, but rather participates in a comprehensive screening and evaluation together with soil / surface pollution risk indicators, air pollution transport risk indicators, pollution source intensity indicators, and pollution migration channel indicators.

[0246] Step S5: Constructing air pollution transport risk indicators

[0247] Step S5 includes: coupling atmospheric remote sensing anomalies with pollution source location, wind direction consistency, and distance attenuation to generate an atmospheric pollution transport risk index (ATR).

[0248] Specifically as follows:

[0249] Step S5 is used to construct an atmospheric pollution transport risk index, denoted as ATR, within a unified spatial grid. ATR stands for Atmospheric Transport Risk. The ATR is used to characterize the likelihood that each spatial grid within the industrial park and its surrounding impact zone will be affected by the diffusion, transport, deposition, or accumulation of atmospheric pollutants.

[0250] This step is not merely about monitoring or alarming the concentration of air pollutants in industrial parks at a single point, nor is it simply about identifying air pollution anomalies at a particular moment. Instead, it couples and calculates air pollution remote sensing products, the spatial location of potential emission sources, prevailing wind direction, wind speed, distance attenuation, downwind relationship, and grid spatial location to form an air pollution transport risk indicator that can be used in a comprehensive screening and evaluation process together with soil / surface pollution risk indicators, water pollution risk indicators, pollution source intensity indicators, and pollution migration channel indicators.

[0251] Step S5 specifically includes the extraction of remote sensing anomalies of air pollution, identification of potential air pollution sources, calculation of wind field influence relationships, calculation of distance attenuation, calculation of downwind grid influence, fusion of air pollution transport risks, and grid assignment.

[0252] 5.1 Extraction of Remote Sensing Anomalies in Air Pollution

[0253] First, remote sensing products of atmospheric pollution in the study area and its surrounding areas are acquired. These remote sensing products include, but are not limited to, products related to NO2, SO2, AOD, CO, aerosols, particulate matter, smoke plume remote sensing identification results, or other remote sensing products that can reflect industrial emissions, atmospheric pollution diffusion, and regional pollution accumulation.

[0254] For remote sensing products of atmospheric pollution, screening, cropping, and quality control are carried out according to the study area and study time window. The quality control includes, but is not limited to, removing low-quality pixels, removing pixels with severe cloud pollution, removing invalid observations, performing time averaging or time synthesis, checking for outliers, and unifying spatial projection.

[0255] Since the spatial resolution of atmospheric pollution remote sensing products is usually lower than that of surface multispectral or hyperspectral remote sensing images, this invention does not directly use them as a single basis for pollution determination, but rather as an important input for identifying atmospheric pollution background anomalies, regional transport trends, and downwind risks.

[0256] In one implementation, standardized outliers can be calculated for a given atmospheric pollutant concentration or optical thickness product:

[0257]

[0258] in, This represents the remote sensing anomaly value of air pollution in the j-th remote sensing pixel or grid. This represents the atmospheric pollutant concentration, column concentration, optical thickness, or related remote sensing inversion value for the j-th pixel or grid. This represents the mean of the background region or historical reference period. This represents the standard deviation of the background region or historical reference period.

[0259] As a subordinate implementation method, multiple remote sensing products for air pollution can be weighted and fused to obtain a comprehensive remote sensing anomaly value for air pollution:

[0260]

[0261] in, Represents a grid The comprehensive atmospheric pollution remote sensing anomaly values, This indicates an outlier for NO2. This indicates an outlier SO2 value. Indicates an abnormal AOD value. This indicates an outlier in the CO value. This indicates anomalies related to particulate matter or aerosols. to These are the weighting coefficients, and .

[0262] When some remote sensing products for air pollution are unavailable, one or more products can be selected for calculation based on the industry type of the industrial park, the characteristics of pollutants, and the availability of data.

[0263] 5.2 Identification of Potential Air Pollution Sources

[0264] Based on the spatial distribution of enterprises, industry type, location of sewage outlets, location of chimneys, waste gas treatment facilities, combustion facilities, boiler rooms, storage yards, loading and unloading areas, road transport channels, historical emission records, existing on-site data, and remote sensing image identification results, potential sources of air pollution are identified.

[0265] The potential sources of air pollution include, but are not limited to:

[0266] 1. Production areas for enterprises in chemical, metallurgical, building materials, new materials, new energy, equipment manufacturing, warehousing and logistics, etc.

[0267] 2. Chimneys, exhaust outlets, fugitive emission areas, combustion facilities, and boiler rooms;

[0268] 3. Dust storage yards, raw material storage yards, solid waste storage areas, and hazardous waste temporary storage areas;

[0269] 4. Areas prone to road dust, loading and unloading areas, vehicle transport routes, and areas affected by construction disturbances;

[0270] 5. Smoke plumes, thermal anomalies, dust anomalies, exposed stockpiles, or areas of abnormal production activity identified by remote sensing or UAV imagery.

[0271] For each potential air pollution source, record its spatial location, source type, industry category, emission characteristics, possible pollutant types, source area, production activity intensity, and historical regulatory information. The potential air pollution source can be represented as Ps, where s is the pollution source number.

[0272] 5.3 Calculation of Wind Field Influence Relationship

[0273] Obtain wind speed components, wind direction frequency, or prevailing wind direction during the study period, and standardize the wind direction to a pollutant transport azimuth angle increasing clockwise from true north (0°). Meteorological data provides the direction of wind origin. At that time, the azimuth angle of the transmission .

[0274] When the wind field is eastward component and northward component When indicating, the azimuth angle is as follows: calculate; and Calm periods with a wind speed of 0 are not included in the direction weight calculation.

[0275] For pollution sources and grid Calculate the azimuth angle from the source point to the center of the grid in the projected coordinate system. The included angle between their circumferences is... .

[0276] Wind direction consistency weight is When the grid is downwind, the weight is close to 1; when it is crosswind or upwind, the weight decreases to 0.

[0277] When there are wind fields with multiple time periods, the frequency is based on the effective time period. Weighted: ,in are non-negative and .

[0278] 5.4 Distance Attenuation Calculation

[0279] As air pollutants are transported from their source outwards, the intensity of their impact typically decreases with increasing distance. Therefore, this invention introduces a distance attenuation function to characterize the impact of the pollution source on different grid cells.

[0280] For each pollution source and each spatial grid Calculate the spatial distance ds,i between the two. The distance can be a planar Euclidean distance, a distance projected along the prevailing wind direction, a distance along the terrain channel, or a distance corrected for terrain obstruction.

[0281] In one implementation, an exponential decay function can be used:

[0282]

[0283] Distance is calculated in the projected coordinate system, and the distance value and the corresponding attenuation parameter use the same length unit.

[0284] in, Indicates pollution source For the grid Distance decay weight, Indicates pollution source To grid The distance is λ, where λ represents the atmospheric transport distance attenuation parameter.

[0285] As a subordinate implementation method, an inverse distance attenuation function can be used:

[0286]

[0287] As a subordinate implementation method, segmented distance weights can be set according to the type of pollution source, the type of pollutant, wind speed conditions, and the management needs of the industrial park. For example, a certain distance range around the pollution source can be designated as a high-impact zone, a more distant downwind area as a medium-impact zone, and grids far from the pollution source and not in the downwind area as low-impact zones.

[0288] 5.5 Calculation of the Intensity of Air Pollution Source Impact

[0289] For each potential source of air pollution The Atmospheric Pollution Source Impact Intensities (ASIs) are constructed. These ASIs can be determined based on factors such as enterprise industry type, emission characteristics, source area area, production activity intensity, number of discharge outlets, chimney height, storage yard area, road dust intensity, historical emission records, remote sensing thermal anomalies, nighttime light intensity, and existing on-site data or regulatory records.

[0290] In one implementation method, the intensity of the impact of air pollution sources can be expressed as:

[0291]

[0292] in, Indicates the industry's emission potential. Indicates the impact factor of emission outlets or chimneys. Indicates historical emissions or regulatory record factor, Indicates thermal anomalies, phosphorescence, or intensity factors of production activities. This indicates the influencing factors of unorganized emissions, dust from storage yards or roads. to These are the weighting coefficients, and .

[0293] It should be noted that, The pollution source intensity index used to characterize the source term impact in the calculation of atmospheric pollution transport risk is used in subsequent step S6. This is used to comprehensively express the overall impact of pollution sources on the combined pollution risks of water, soil, and air. Both can use some of the same basic pollution source data, but they differ in their calculation purpose and their position in the calculation process.

[0294] 5.6 Calculation of the influence of the downwind grid

[0295] By combining remote sensing anomalies of atmospheric pollution, the intensity of pollution source impact, wind direction consistency weight, and distance attenuation weight, the pollution source is calculated. For the grid Atmospheric transport impact value:

[0296]

[0297] in, Indicates pollution source For the grid The atmospheric transport impact value, Indicates pollution source The intensity of the impact of air pollution sources Indicates the weight of wind direction consistency. This represents the distance decay weight.

[0298] When multiple potential sources of air pollution exist, all pollution sources can be mapped to a grid. The influence values ​​are superimposed, the maximum value is taken, or a weighted fusion is performed to obtain the mesh. Atmospheric transport potential:

[0299]

[0300] or:

[0301]

[0302] in, Represents a grid Atmospheric transport potential.

[0303] In one implementation, for areas with significant remote sensing anomalies in air pollution, the remote sensing anomaly values ​​of air pollution can be further... Atmospheric transport potential Coupling yields atmospheric pollution transport risk indicators. :

[0304]

[0305] in, Represents a grid air pollution transport risk indicators Represents a grid Remote sensing anomalies of air pollution Represents a grid Atmospheric transport potential, and These are the weighting coefficients, and .

[0306] As a subordinate implementation method, a product coupling approach can be adopted:

[0307]

[0308] This approach emphasizes that when a grid simultaneously exhibits high remote sensing anomalies of atmospheric pollution and high downwind transport potential, its atmospheric pollution transport risk is relatively high.

[0309] As a subordinate implementation method, a combination of weighted superposition and product coupling can be adopted:

[0310]

[0311] in, , , These are the weighting coefficients, and .

[0312] 5.7 Normalization of air pollution transport risks

[0313] In order to To achieve comprehensive calculations with soil / surface pollution risk indicators, water pollution risk indicators, pollution source intensity indicators, and pollution migration pathway indicators, it is necessary to... Normalize to the 0-1 interval.

[0314] In one implementation, a minimum-maximum normalization method can be used:

[0315]

[0316] When the upper limit of normalization equals the lower limit, the indicator is not included in the fusion in the current batch. When comparing multiple periods, the same reference period and fixed normalization parameters are used for each period; the results of independent normalization for each period are only used for internal ranking within each period.

[0317] in, Represents a grid Normalized air pollution transport risk indicators and Representing the research scope The minimum and maximum values.

[0318] As a subordinate implementation method, quantile normalization, robust normalization, standard deviation normalization, or logistic function normalization methods can be used to reduce the impact of abnormal terrain values, extreme confluence values, abnormal wind field periods, or artificial ditch data errors on the results.

[0319] Normalized The closer a value is to 1, the higher the risk that the grid will be affected by the transport, diffusion, or deposition of air pollution; the closer a value is to 0, the lower the risk that the grid will be affected by the transport of air pollution.

[0320] 5.8 Spatial Expression of Air Pollutant Transport Risks

[0321] Normalized The data is assigned to a unified spatial grid, and an atmospheric pollution transport risk map is generated. This map is used to express the potential impact of emission sources from industrial parks on the transport of atmospheric pollution to different areas under the combined effects of prevailing wind direction, wind speed, distance attenuation, and remote sensing anomalies of atmospheric pollution.

[0322] For grids located downwind of the pollution source, close to the pollution source, with high superimposed atmospheric pollution remote sensing anomalies, or in areas with unfavorable terrain and ventilation. The values ​​are relatively high; for grids located upwind, far from pollution sources, with low atmospheric remote sensing anomalies, or significantly obstructed by terrain, The value is low.

[0323] In one implementation, it can be based on The magnitude of the numerical value is used to classify the risk of air pollution transport into low risk, low-to-medium risk, medium risk, medium-to-high risk, and high risk levels. This classification can be achieved using the natural breakpoint method, quantile method, fixed threshold method, or standard deviation method.

[0324] Through the above steps, the atmospheric pollution transport risk index for each spatial grid can be obtained. and use it as a subsequent comprehensive screening index. One of the components of the index. It is not used as the sole final determination of air pollution in industrial parks, but rather participates in a comprehensive screening and evaluation together with soil / surface pollution risk indicators, water pollution risk indicators, pollution source intensity indicators, and pollution migration channel indicators.

[0325] Step S6: Construct pollution source intensity indicators

[0326] Step S6 includes: generating pollution source intensity indicators based on information on enterprises, pollution discharge points, storage yards, and industrial activities. .

[0327] Specifically as follows:

[0328] Step S6 is used to construct a pollution source intensity index within a unified spatial grid, denoted as... , This represents Pollution Source Intensity, an index indicating the intensity of pollution sources. It is used to characterize the intensity of the influence of industrial enterprises, sewage outlets, storage yards, hazardous waste temporary storage sites, sewage treatment facilities, accident pools, road transport channels, historically polluted sites or other potential pollution sources on each spatial grid within the industrial park and its surrounding influence area.

[0329] This step is not simply an aggregation based on the number or type of enterprises, nor is it a separate statistical analysis of pollution sources. Instead, it comprehensively calculates pollution source intensity indicators by considering factors such as pollution source type, spatial location, impact range, industry pollution potential, sewage discharge facilities, storage yards and hazardous waste temporary storage areas, historical pollution records, proximity of remote sensing anomalies, and distance attenuation relationships with grid cells. This results in a pollution source intensity index that can be used in conjunction with soil / surface pollution risk indicators, water pollution risk indicators, air pollution transport risk indicators, and pollution migration channel indicators for comprehensive screening and evaluation.

[0330] Step S6 specifically includes pollution source object identification, pollution source classification and weighting, determination of the spatial influence range of pollution sources, distance attenuation calculation, remote sensing anomaly proximity correction, pollution source intensity fusion, and gridded value assignment.

[0331] 6.1 Identification of Pollution Sources

[0332] First, based on enterprise spatial distribution data, industrial land data, pollutant discharge permit data, environmental supervision data, existing on-site data, remote sensing image recognition results, drone patrol results, and historical monitoring data, potential pollution sources within the research scope are identified.

[0333] The sources of pollution include, but are not limited to:

[0334] 1. Industrial enterprises that may generate pollution emissions, including chemical, metallurgical, new materials, new energy, equipment manufacturing, building materials, warehousing and logistics, sewage treatment, solid waste disposal, or other industries.

[0335] 2. Wastewater discharge outlets, rainwater discharge outlets, exhaust gas discharge outlets, chimneys, fugitive emission areas, and factory boundary emission channels;

[0336] 3. Raw material storage yards, coal piles, ore powder piles, slag yards, waste residue storage areas, solid waste temporary storage areas, and hazardous waste temporary storage areas;

[0337] 4. Wastewater treatment plant, emergency pool, sedimentation tank, evaporation tank, equalization tank, wastewater storage tank and leachate collection facilities;

[0338] 5. Loading and unloading areas, storage areas, vehicle transport lanes, areas with high incidence of road dust, and material transfer areas;

[0339] 6. Areas involving historically contaminated sites, suspected contaminated sites, historical accident sites, historical exceedance sites, and regulatory penalty records;

[0340] 7. Abnormal bare land, abnormal storage yards, abnormal water bodies, thermal anomaly areas, plume anomaly areas, or abnormal surface disturbance areas identified by remote sensing imagery, UAV imagery, or existing field data.

[0341] For each pollution source object, record its spatial location, object type, area or length, affiliated enterprise, industry category, possible pollutant types, emission medium, historical records, and spatial relationship with water systems or surrounding environmental objects. A pollution source object can be represented as... , where s is the pollution source number.

[0342] 6.2 Classification and empowerment of pollution sources

[0343] Different types of pollution sources contribute differently to the risk of combined pollution of water, soil, and air, therefore it is necessary to classify and assign weights to pollution sources.

[0344] In one implementation, a basic weight can be set according to the type of pollution source. The basic weights can be determined according to the potential impact of pollution sources on environmental media. For example:

[0345] For entities with high pollution potential, such as chemical industry, metallurgy, hazardous waste temporary storage, sewage treatment, and solid waste storage, a higher basic weight can be set.

[0346] For objects with moderate pollution potential or related to specific processes, such as new materials, new energy, equipment manufacturing, warehousing and logistics, and road transportation, a medium basic weight can be set.

[0347] For general construction land, ordinary warehousing, low-pollution industrial land, or areas that have completed remediation, a lower basic weight can be set.

[0348] As a subordinate implementation method, media weights can be set according to the main influencing media of the pollution source, including soil influence weights, water body influence weights, and atmospheric influence weights. For example:

[0349] Stockpiles, slag yards, exposed material areas, and road dust deposition areas are primarily assigned a higher soil impact weight;

[0350] Sewage outlets, ditch inlets, sewage treatment facilities, emergency pools, and sedimentation tanks are primarily assigned a higher weighting for their impact on water bodies.

[0351] Chimneys, exhaust outlets, fugitive emission areas, dust dumps, and combustion facilities are primarily assigned a higher weight for atmospheric impact.

[0352] pollution source The comprehensive pollution potential weight can be expressed as:

[0353]

[0354] in, Indicates pollution source The overall basic weight, Indicates the industry's pollution potential weight. Indicates the influence weight of the medium. Indicates the weight of historical pollution or regulatory records. This indicates the weight of facilities such as sewage treatment facilities, storage yards, hazardous waste temporary storage areas, or emergency pools. to These are the weighting coefficients, and .

[0355] 6.3 Determination of the Spatial Impact Range of Pollution Sources

[0356] The impact of pollution sources on surrounding grids exhibits spatial decay characteristics. Different types of pollution sources have different impact ranges; therefore, it is necessary to determine the impact range of pollution sources based on the type of pollution source, the pollutant medium, and the spatial conditions of the industrial park.

[0357] In one implementation, for point pollution sources, including sewage outlets, chimneys, accident pools, hazardous waste storage sites, historical accident sites, etc., a buffer zone or distance attenuation zone centered on the location of the pollution source can be established.

[0358] As a subordinate implementation method, for linear pollution sources, including sewage ditches, rainwater ditches, road transport channels, loading and unloading transport routes, buffer zones or distance attenuation zones can be established along the linear objects.

[0359] As a subordinate implementation method, for non-area pollution sources, including enterprise plant areas, storage yards, sewage treatment facilities, solid waste storage areas, industrial land patches, historically polluted sites, etc., the scope of influence can be established outward from the boundary of the non-area object, and modified in combination with the downstream confluence direction, downwind direction and the location of surrounding environmental objects.

[0360] The influence range of the pollution source can be determined by fixed distance buffer, graded buffer, distance attenuation function, extension along water system, extension along prevailing wind direction, extension along DEM confluence path, or a combination of multiple methods.

[0361] 6.4 Distance Attenuation Calculation

[0362] For each pollution source and each spatial grid Calculate the distance between the two. The distance can be the distance from the grid center to the pollution source location, the shortest distance from the grid center to the pollution source boundary, the hydrological distance along the water system or ditch, the projected distance along the prevailing wind direction, the distance along the road transport corridor, or the distance corrected based on the actual pollution migration path.

[0363] In one implementation, an exponential decay function can be used to calculate the impact of pollution sources on the grid:

[0364]

[0365] Distance is calculated in the projected coordinate system, and the distance value and the corresponding attenuation parameter use the same length unit.

[0366] in, Indicates pollution source For the grid Distance decay weight, Indicates pollution source With grid The distance between them Indicates pollution source The influence of distance parameters.

[0367] As a subordinate implementation method, an inverse distance attenuation function can be used:

[0368]

[0369] As a subordinate implementation method, a piecewise attenuation function can be used. For example, a high weight can be assigned to the core influence range of the pollution source, a medium weight to the secondary influence range, a low weight to the far-field influence range, and a value of 0 can be assigned after the maximum influence distance.

[0370] For different types of pollution sources, Ls can be set with different values. For example, for pollution sources such as sewage outlets and ditches, a longer influence distance can be set along the direction of water flow; for pollution sources such as stockpiles, bare land, and road dust, a higher influence weight can be set within the short distance range; and for exhaust gas emission sources, the downwind influence distance can be set in combination with the prevailing wind direction and wind speed.

[0371] 6.5 Migration Direction Correction

[0372] The impact of pollution source intensity on spatial grids depends not only on distance but also on the potential migration direction of pollutants. Therefore, hydrological and atmospheric orientation corrections can be introduced when calculating PSI.

[0373] For pollution sources that may migrate via surface runoff, ditches, sewage channels, or river systems, the grid is determined based on the DEM's confluence direction, water flow direction, and ditch connectivity. Is it located near a pollution source? Downstream, confluence path, or low-lying collection area. If the grid If the pollution source is located downstream or on a confluence path, its impact weight is increased; if it is located upstream and there is no obvious hydrological connectivity, its impact weight is decreased.

[0374] For pollution sources that may affect the surrounding area through exhaust emissions, dust diffusion, or dust deposition, the grid is determined based on the prevailing wind direction, wind speed, and the azimuth angle from the pollution source to the grid. Is it located near a pollution source? The downwind influence area. If the grid If the pollution source is located in the downwind influence area, the weight of the pollution source influence is increased; if it is located upwind or has a weak relationship with the prevailing wind direction, the weight of the pollution source influence is decreased.

[0375] The pollution source orientation correction factor can be expressed as:

[0376]

[0377] in, Indicates pollution source For the grid Direction correction factor, Indicates the correction factor for hydrological migration direction. This represents the atmospheric transport direction correction factor. Indicates the nearest neighbor correction factor. , , These are the weighting coefficients, and .

[0378] For pollution sources that primarily affect a single medium, only one or two directional correction factors may be used; for pollution sources that may simultaneously affect multiple media such as water, soil, and air, a comprehensive directional correction factor may be used.

[0379] 6.6 Remote Sensing Anomaly Proximity Correction

[0380] To integrate pollution source intensity indicators with remote sensing identification results, this invention introduces a proximity correction for remote sensing anomalies. For pollution sources located in or near remote sensing anomaly areas for soil pollution, water pollution, or high-risk areas for air pollution transport, their contribution to the pollution source intensity of surrounding grids can be increased.

[0381] In one implementation, if the pollution source High A grid indicates that there may be abnormal soil or surface sedimentation around the pollution source, which can increase the soil impact weight of the pollution source.

[0382] As a subordinate implementation method, if the pollution source High levels of pollutants exist downstream, near sewage outlets, or near ditch confluences. A grid indicates that the pollution source may have a spatial correlation with water anomalies, which can increase the water body impact weight of the pollution source.

[0383] As a subordinate implementation method, if the pollution source High in the downwind area If the results of grid analysis, remote sensing anomalies of atmospheric pollution, or plume identification indicate that the pollution source may be spatially associated with atmospheric transport risks, the atmospheric impact weight of the pollution source can be increased.

[0384] The remote sensing anomaly proximity correction factor can be expressed as:

[0385]

[0386] in, Indicates pollution source Remote sensing anomaly proximity correction factor, Indicates pollution source Remote sensing anomaly proximity values ​​of surrounding soil pollution Indicates pollution source Remote sensing anomaly values ​​of pollution in downstream or adjacent water bodies Indicates pollution source Downwind atmospheric pollution transport anomaly proximity value, , , These are the weighting coefficients, and .

[0387] 6.7 Calculation of Pollution Source Intensity Index

[0388] The pollution source is calculated by combining the basic weight of the pollution source, distance attenuation, migration direction correction, and remote sensing anomaly proximity correction. For the grid Pollution source impact value:

[0389]

[0390] in, Indicates pollution source For the grid The pollution source impact value, Indicates pollution source The overall basic weight, Indicates the distance decay weight. Indicates the migration direction correction factor. This represents the correction factor for proximity to remote sensing anomalies.

[0391] When multiple pollution sources exist, the grid can be calculated using methods such as superposition, maximum value, weighted superposition, or saturation function. The comprehensive pollution source intensity index.

[0392] In one implementation, an overlay method is used:

[0393]

[0394] in, Represents a grid Pollution source intensity index.

[0395] As a subordinate implementation method, the maximum value method is adopted:

[0396]

[0397] This approach is suitable for situations where the impact of the dominant pollution source is highlighted.

[0398] As a subordinate implementation method, a saturated function approach is adopted:

[0399]

[0400] This method is suitable for situations where multiple pollution sources have overlapping effects, and can avoid the problem of too many pollution sources causing... It can grow indefinitely.

[0401] 6.8 Normalization of pollution source intensity indicators

[0402] In order to To achieve comprehensive calculations with soil / surface pollution risk indicators, water pollution risk indicators, air pollution transport risk indicators, and pollution migration pathway indicators, it is necessary to... Normalize to the 0-1 interval.

[0403] In one implementation, a minimum-maximum normalization method can be used:

[0404]

[0405] When the upper limit of normalization equals the lower limit, the indicator is not included in the fusion in the current batch. When comparing multiple periods, the same reference period and fixed normalization parameters are used for each period; the results of independent normalization for each period are only used for internal ranking within each period.

[0406] in, Represents a grid Normalized pollution source intensity index and Representing the research scope The minimum and maximum values.

[0407] As a subordinate implementation method, quantile normalization, robust normalization, standard deviation normalization, or logistic function normalization methods can be used to reduce the impact of abnormal terrain values, extreme confluence values, abnormal wind field periods, or artificial ditch data errors on the results.

[0408] Normalized The closer a value is to 1, the stronger the influence of pollution sources on the grid; the closer a value is to 0, the weaker the influence of pollution sources on the grid.

[0409] 6.9 Spatial Expression of Pollution Source Intensity

[0410] Normalized The values ​​are assigned to a unified spatial grid, and a pollution source intensity distribution map is generated. This map is used to express the degree of influence of different pollution sources on the source terms of each spatial grid within the industrial park and its surrounding impact zone.

[0411] For grids located inside or around enterprises with high pollution potential, near sewage outlets, around storage yards, around hazardous waste temporary storage areas, around accident pools, around sewage treatment facilities, in areas with high incidence of road dust, or in downstream confluence paths or downwind areas. The value is relatively high; for grids that are far from pollution sources, have no obvious migration connectivity, or have low basic weights for pollution sources, The value is low.

[0412] In one implementation, it can be based on The intensity of pollution sources is classified into low, low-medium, medium, medium-high, and high levels based on the numerical value of the pollution source. This classification can be achieved using the natural breakpoint method, quantile method, fixed threshold method, or standard deviation method.

[0413] Through the above steps, the pollution source intensity index for each spatial grid can be obtained. and use it as a subsequent comprehensive screening index. One of the components of the index. It is not used as the sole final determination of pollution risk in industrial parks, but rather participates in a comprehensive screening and evaluation together with soil / surface pollution risk indicators, water pollution risk indicators, air pollution transport risk indicators, and pollution migration channel indicators.

[0414] Step S7: Constructing indicators for pollution migration pathways

[0415] Step S7 includes: generating pollution migration channel indicators based on DEM confluence, water system topology, ditch drainage, and downwind relationships. .

[0416] Specifically as follows:

[0417] Step S7 is used to construct pollution migration channel indicators within a unified spatial grid, denoted as... , This refers to the Migration Transport Index, which is an indicator of pollution migration pathways. This is used to characterize the potential of each spatial grid within an industrial park and its surrounding impact zone to serve as a pollutant migration path, accumulation area, transport corridor, or diffusion impact zone.

[0418] This step does not rely solely on the distance from the pollution source for risk assessment, nor does it simply involve conventional topographic or water system analysis. Instead, it couples the DEM confluence path, river system connectivity, ditches and sewage channels, low-lying catchment areas, surface runoff direction, prevailing wind transport path, and the spatial location of the pollution source to construct a spatial channel indicator that can express the possibility of pollutants migrating from the source area to downstream, downwind, and surrounding impact areas.

[0419] Step S7 includes hydrological migration channel identification, atmospheric migration channel identification, migration channel importance weighting, pollution source connectivity calculation, comprehensive migration channel index calculation, and gridded value assignment.

[0420] 7.1 Hydrological Migration Channel Indicators

[0421] Hydrological migration channels are marked as , This refers to the Hydrological Transport Index, which is an index of hydrological migration channels. It is used to characterize the likelihood of pollutants migrating through surface runoff, rivers, ditches, sewage channels, stormwater runoff channels, low-lying catchment areas, and water system connectivity pathways.

[0422] First, slope, aspect, surface runoff direction, runoff accumulation, runoff path, and low-lying catchment areas are calculated based on DEM data. For industrial parks with obvious artificial ditches, sewage ditches, stormwater pipe networks, river channel modifications, or factory drainage systems, the DEM runoff results can be corrected by combining water system vector data, sewage ditches data, stormwater pipe network data, UAV imagery, existing field data, and high-resolution remote sensing imagery.

[0423] In one implementation, hydrological migration channels include the following types:

[0424] 1. The main channel of the river, its tributaries, the riverbank zone, and the downstream receiving water bodies;

[0425] 2. Ditches, sewage ditches, rainwater ditches, and open channels within the industrial park;

[0426] 3. Surface runoff confluence paths around storage yards, plant boundaries, roads, bare ground, and accident pools;

[0427] 4. Low-lying areas, waterlogged areas, ponds, sedimentation tanks, and areas suspected of accumulating pollutants;

[0428] 5. Potential hydrological connectivity pathways from key enterprises, sewage outlets, storage yards, hazardous waste temporary storage sites, accident ponds, or exposed ground surfaces to downstream water bodies, ditch inlets, and surrounding impact areas.

[0429] For each spatial grid Hydrological migration channel indicators can be calculated based on factors such as the cumulative amount of runoff, whether it is located in the runoff path, whether it intersects with rivers or ditches, whether it is located in a low-lying catchment area, whether it is downstream of a pollution source, and whether it is connected to a downstream receiving water body.

[0430] In one implementation, It can be represented as:

[0431]

[0432] in, Represents a grid Hydrological migration channel indicators This represents the normalized value of the cumulative flow. The channel factor indicates whether it is located on a confluence path or a ditch path. This indicates the hydrological connectivity factor with rivers, ditches, drainage channels, or stormwater channels. Indicates the low-lying accumulation zone factor. This indicates the influence factor located downstream of the pollution source. Indicates the influencing factors related to connectivity with downstream receiving water bodies or the water systems of surrounding influence areas. to These are the weighting coefficients, and .

[0433] As a subordinate implementation method, for areas within industrial parks where artificial drainage systems are more prominent, the weights of ditches, sewage ditches, rainwater pipe networks, and inlets can be increased; for areas where natural slope runoff has a significant impact, the weights of DEM runoff accumulation, slope aspect, and low-lying collection areas can be increased; for industrial parks where rivers flow through or are adjacent to rivers, the weights of river connectivity and upstream-downstream relationships can be increased.

[0434] 7.2 Calculation of Hydrological Connectivity

[0435] In order to identify the spatial connections between pollution sources and water bodies and surrounding environmental objects, this invention further calculates hydrological connectivity relationships.

[0436] For each pollution source and grid Based on the DEM confluence direction, water system topology, ditch connectivity, and sewage ditch flow direction, the grid is determined. Is it located near a pollution source? On the downstream confluence path. If the grid Located at the pollution source Potential hydrological migration paths to water bodies or surrounding environmental objects are assigned higher hydrological connectivity weights; if the grid... With pollution sources If there is no obvious hydrological connection between them, then a lower weight is assigned.

[0437] Hydrological connectivity weights can be expressed as:

[0438]

[0439] in, Indicates pollution source With grid Hydrological connectivity weights between them Represents the connectivity discriminant factor, when the grid... Located at the pollution source When the downstream confluence path, ditch connection path, or river system connection path is on the river, take the higher value; otherwise take the lower value or 0. Indicates the distance along the hydrological channel; This represents the hydrological migration distance attenuation parameter.

[0440] When multiple pollution sources exist, each pollution source and grid can be monitored. The hydrological connectivity weights are superimposed, the maximum value is taken, or a saturation function is used to fuse them to correct the HTI.

[0441] 7.3 Atmospheric migration channel indicators

[0442] Atmospheric migration channels are marked as , This refers to the Atmospheric Transport Index, which is an index of atmospheric transport routes. It is used to characterize the likelihood that pollutants will affect the spatial grid through prevailing winds, downwind diffusion, dust transport, particulate matter deposition, or plume paths.

[0443] It should be noted that in step S5 above... The key focus is on characterizing the atmospheric pollution transport risk resulting from the combined effects of remote sensing anomalies and downwind transport from pollution sources; this section... The focus is on characterizing the channel properties of spatial grids as atmospheric pollution migration channels or downwind diffusion paths. Both can use some of the same wind field and pollution source baseline data, but... More biased towards pollution risk intensity It is more inclined to consider migration pathway conditions.

[0444] In one implementation, for each potential air pollution source and spatial grid Calculate the angle between the azimuth of the pollution source pointing to the grid center and the prevailing wind direction, and combine this with wind speed, wind direction frequency, terrain obstruction, and distance attenuation to determine the grid. Is it located near a pollution source? The downwind direction affects the passageway.

[0445] It can be represented as:

[0446]

[0447] in, Represents a grid Atmospheric migration channel indicators, This represents the downwind consistency factor. This indicates the frequency of wind direction or the stability factor of the prevailing wind direction. This represents the atmospheric transport distance attenuation factor. Indicates the terrain's ventilation or barrier factor. Indicates the influencing factors that are spatially connected to air pollution sources. to These are the weighting coefficients, and .

[0448] As a secondary implementation method, the fan-shaped area downwind of the pollution source can be used as an atmospheric migration channel. (When the grid...) When located within a fan-shaped area downwind of the pollution source, and close to the pollution source or located in the overlapping downwind area of ​​multiple pollution sources, Higher values; when the grid When located upwind, far from pollution sources, or significantly obstructed by terrain, The value is low.

[0449] 7.4 Atmospheric Channel Topographic Correction

[0450] The transport of air pollutants from industrial parks is affected not only by wind direction and speed, but also by topography, building clusters, valleys, river valleys, open areas, and low-lying areas. Therefore, this invention can introduce a topographic correction factor to modify the atmospheric migration channels.

[0451] In one implementation, the impact of topography on atmospheric pollution transport can be determined by calculating topographic openness, slope aspect and wind direction consistency, valley passages, topographic barrier height, or relative elevation difference based on the DEM. When the grid is located in an open passage, valley passage, or low-barrier area consistent with the prevailing wind direction, its atmospheric migration channel value is increased; when the grid is significantly blocked by mountains, building complexes, or elevation differences, its atmospheric migration channel value is decreased.

[0452] The topographically corrected atmospheric migration channel index can be expressed as:

[0453]

[0454] in, This indicates the atmospheric migration channel index after terrain correction. This represents the terrain correction factor. The determination can be made based on the terrain openness, slope consistency, relative elevation difference, building density, or other ventilation conditions.

[0455] 7.5 Calculation of Comprehensive Migration Channel Indicators

[0456] Comprehensive hydrological migration channel indicators and atmospheric migration channel indicators Construct comprehensive pollution migration channel indicators .

[0457] In one implementation, A weighted fusion method can be used:

[0458]

[0459] in, Represents a grid Comprehensive pollution migration channel indicators, Represents a grid Hydrological migration channel indicators Represents a grid Atmospheric migration channel indicators, and These are the weighting coefficients, and .

[0460] As a subordinate implementation method, when the influence of water systems, ditches, sewage channels, and surface runoff in the study area on pollution migration is more significant, it can improve... Weighting; when exhaust emissions, dust deposition, downwind influence, and atmospheric diffusion risks are more pronounced in the study area, the weighting can be increased. Weighting; when both hydrological migration and atmospheric migration are significant, a relatively balanced weighting can be used.

[0461] As a subordinate implementation method, calculations can be performed separately according to different pollution source types. For pollution sources that are primarily likely to migrate through water bodies or surface runoff, Mainly composed of Identify; for pollution sources that are primarily likely to migrate through exhaust gas, plumes, dust, or particulate matter deposition, Mainly composed of It is determined that for comprehensive pollution sources that simultaneously pose risks of wastewater, waste gas, solid waste storage, and dust generation, Depend on and To be determined jointly.

[0462] 7.6 Pollution Source-Migration Corridor Connectivity

[0463] In order to To better serve comprehensive screening and sampling site deployment, this invention further calculates the connectivity between pollution sources and migration channels.

[0464] For each spatial grid Determine whether it simultaneously meets one or more of the following conditions:

[0465] 1. Located downstream of the pollution source's confluence path, ditch connection path, or downstream path of a river system;

[0466] 2. Located downwind of the pollution source or within the fan-shaped area influenced by the prevailing wind direction;

[0467] 3. Located on potential migration paths between pollution sources and water systems, ditches, road dust deposition zones, low-lying catchment areas, or downwind subsidence zones;

[0468] 4. Located in an area where the migration paths of multiple pollution sources overlap;

[0469] 5. Located in the transitional area where areas with high risk of soil or surface remote sensing anomalies, water body remote sensing anomalies, or air pollution transport intersect with migration channels;

[0470] 6. Located in low-lying catchment areas, river inlets, ditch confluences, road dust deposition zones, downwind subsidence zones, or other key migration nodes.

[0471] If grid Having both pollution source connectivity and migration channel attributes improves its Value. This process can identify sampling candidate nodes that are not particularly strong in terms of simple pollution anomalies, but are located on critical migration paths.

[0472] In one implementation, a connectivity enhancement factor may be introduced:

[0473]

[0474] in, This represents the comprehensive migration channel index after enhanced connectivity. This represents the connectivity enhancement factor. According to the grid The number, intensity, and superposition relationships between the pollution source, ditches, water system, low-lying catchment area, downwind subsidence zone, and other migrating elements are determined.

[0475] 7.7 Normalization of Pollution Migration Channel Indicators

[0476] In order to To achieve comprehensive calculations with soil / surface pollution risk indicators, water pollution risk indicators, air pollution transport risk indicators, and pollution source intensity indicators, it is necessary to... Normalize to the 0-1 interval.

[0477] In one implementation, a minimum-maximum normalization method can be used:

[0478]

[0479] When the upper limit of normalization equals the lower limit, the indicator is not included in the fusion in the current batch. When comparing multiple periods, the same reference period and fixed normalization parameters are used for each period; the results of independent normalization for each period are only used for internal ranking within each period.

[0480] in, Represents a grid Normalized pollution migration pathway indicators and Representing the research scope The minimum and maximum values.

[0481] As a subordinate implementation method, quantile normalization, robust normalization, standard deviation normalization, or logistic function normalization methods can be used to reduce the impact of abnormal terrain values, extreme confluence values, abnormal wind field periods, or artificial ditch data errors on the results.

[0482] Normalized The closer a value is to 1, the higher the likelihood that the grid will serve as a pollution migration channel, pollution accumulation area, or pollution transport impact area; the closer a value is to 0, the lower the likelihood that the grid will serve as a pollution migration channel.

[0483] 7.8 Spatial Representation of Pollution Migration Channels

[0484] Normalized The data is assigned to a unified spatial grid, and a pollution migration channel map is generated. This map is used to express the spatial pattern of pollutants from industrial parks that may migrate along surface runoff, ditches, river systems, low-lying catchment areas, prevailing wind directions, downwind diffusion paths, and particulate matter deposition paths.

[0485] For grids located on rivers, ditches, sewage channels, stormwater runoff paths, low-lying catchment areas, downstream inlets, downstream paths of pollution sources, or downwind paths of pollution sources, The values ​​are relatively high; for grids located upstream of the terrain, far from water systems, without obvious confluence and connectivity, or upwind of pollution sources or with weak migration channel attributes, The value is low.

[0486] In one implementation, it can be based on The numerical value of the pollution migration channel is used to divide the pollution migration channel into low-channel influence zone, medium-low-channel influence zone, medium-channel influence zone, medium-high-channel influence zone, and high-channel influence zone. This classification can be achieved using the natural breakpoint method, quantile method, fixed threshold method, or standard deviation method.

[0487] By following the steps above, the pollution migration pathway indicators for each spatial grid can be obtained. and use it as a subsequent comprehensive screening index. One of the components of the index. It is not used as the sole final determination of pollution risk in industrial parks, but rather participates in risk zoning and sampling point selection together with soil / surface pollution risk indicators, water pollution risk indicators, air pollution transport risk indicators, and pollution source intensity indicators.

[0488] Step S8: Construct the Comprehensive Screening Index (CPRI)

[0489] Step S8 includes: normalizing and weighting the effective indicators based on their credibility, and constructing a comprehensive screening index according to "media anomaly - source term support - migration conditions". .

[0490] Specifically as follows:

[0491] 8.1 Indicator Validity, Reliability, and Normalization

[0492] For the grid , obtain , , , and Each indicator is located in the interval [0, 1], and its validity is recorded separately. and credibility .

[0493] Missing data due to inapplicable media type is not treated as zero-risk; missing data due to clouds, shadows, or low-quality observations can be supplemented with historical observations from the same season or nearby valid observations, and the risk level will be reduced. .

[0494] When using minimum-maximum normalization, if the upper limit equals the lower limit, the indicator will not be included in the fusion in the current batch. When used for multi-period comparisons, each period uses the same background period, the same upper and lower limits, and the same weights; the results of normalization according to the minimum and maximum values ​​of each period are only used for relative ranking within the period.

[0495] 8.2 Multi-media anomaly comprehensive index

[0496] make For grid Set of effective media indicators, In order to represent , and The comprehensive index of multi-media anomalies is: .

[0497] in, Non-negative medium weights. Denominator recalibration ensures that non-aqueous meshes are not affected by... Not applicable and mechanically assigned zero; when the denominator is 0, This is recorded as no valid observation.

[0498] 8.3 Source Term Support and Migration Sub-Indices

[0499] Source term support sub-index It is used to express the combined effect of abnormal media and pollution sources.

[0500] Migration index It is used to enhance anomalous grids located downstream of pollution sources, downwind, or in low-lying convergence paths. , All coefficients in each group are non-negative and the sum of the coefficients in each group is 1.

[0501] In this methodology, no independent sub-index for exposure to surrounding objects is set; information on surrounding residential areas, farmland, water bodies, and ecological protection objects is only used as a reference for the research scope, background control points, and sampling accessibility constraints, and is not included in the methodology. Necessary calculations.

[0502] 8.4 Comprehensive Screening Index

[0503] Preferably, ,in to The weights are non-negative and their sum is 1. , , , and All are located in the interval [0, 1], therefore It is also located in the interval [0, 1].

[0504] The aforementioned items respectively retain high-value information on media anomalies, pollution source intensity, and migration channel conditions, and through... and Enhance the mesh that simultaneously possesses "medium anomaly - source term support - migration conditions".

[0505] It is used for relative risk zoning within the study area, candidate patch extraction, and sampling priority ranking, but not for directly determining the concentration of specific pollutants or the environmental quality compliance status.

[0506] 8.5 Weight determination, sensitivity analysis, and parameter tracking

[0507] The weights can be determined using preset initial weights, expert prior weights, or the results of stability comparisons of multiple sets of weights, and the calculation caliber should be kept consistent within the same batch before the candidate sampling points are output.

[0508] Each calculation saves the data time window, spatial reference, normalization upper and lower limits, background reference, weights, thresholds, distance units, and result version. All charts are regenerated after weight or threshold changes, and different parameter versions are not mixed within the same result chain.

[0509] 8.6 Spatial Neighborhood Correction and Credibility Output

[0510] Isolated noisy meshes can be corrected using connected components, minimum patch area, and neighborhood consistency, but local hotspots supported by both high-confidence single anomalies and explicit source terms are not removed. Except... In addition, it simultaneously outputs the overall credibility and missing data identifiers.

[0511] Step S9: Determine the risk level of the comprehensive screening

[0512] Step S9 includes: based on Risk levels are classified using fixed thresholds, quantile thresholds, or natural breakpoint thresholds.

[0513] Specifically as follows:

[0514] Step S9 is used to determine the comprehensive screening index. Risk levels are classified for the industrial park and its surrounding impact zone. This risk level classification is used to categorize continuous... Converted into spatial hierarchy results that can be used for candidate region extraction and sampling point placement.

[0515] This step does not directly determine whether environmental quality meets or exceeds standards based on the concentration of a single pollutant, nor does it replace the results of statutory environmental monitoring. Instead, it uses a comprehensive screening index formed by multi-source remote sensing anomalies, pollution source intensity, and pollution migration channels to spatially classify and express candidate risks of water-soil-air composite pollution in industrial parks.

[0516] 9.1 Risk Level System

[0517] In one implementation, according to Based on the numerical values, the study area is divided into five risk levels: 1. Low-risk area; 2. Low-to-medium risk area; 3. Medium-risk area; 4. Medium-to-high risk area; 5. High-risk area.

[0518] Among them, low-risk areas indicate that the overall pollution remote sensing anomalies, pollution source impacts, and migration channel conditions in the area are relatively weak; medium-low risk areas indicate that the area has certain risk factors, but the comprehensive screening value is low; medium-risk areas indicate that the area has obvious pollution anomalies, pollution source impacts, or migration channel conditions, and can be used as general candidate sampling areas; medium-high risk areas indicate that the area has a high comprehensive screening value and should be used as priority sampling candidate areas; high-risk areas indicate that the area has strong media anomalies, source term support, and migration conditions, and should be used as the highest priority candidate sampling areas.

[0519] As a subordinate implementation method, the risk level can also be divided into three, four, or more levels according to the sampling site layout requirements. For example, it can be divided into low-risk, medium-risk, and high-risk areas, or into general candidate areas, key candidate areas, and priority sampling areas. The number of risk levels can be determined based on the study area, data accuracy, and the availability of sampling resources.

[0520] 9.2 Fixed Threshold Division Method

[0521] In one implementation, a fixed threshold method can be used for... Risk levels are classified. When When the data has been normalized to the 0-1 interval, the following threshold can be used:

[0522] When 0 ≤ When the value is < 0.2, it is classified as a low-risk area;

[0523] When 0.2 ≤ When the value is less than 0.4, it is classified as a medium-to-low risk area;

[0524] When 0.4 ≤ When the value is < 0.6, it is classified as a medium-risk area;

[0525] When 0.6 ≤ If the value is less than 0.8, the area is classified as a medium- to high-risk area.

[0526] When 0.8 ≤ When the value is ≤ 1.0, it is classified as a high-risk area.

[0527] The fixed threshold method is only used for annual or multi-park comparisons when the same background reference, normalization parameters, weights and thresholds are used in each batch; when the parameters are inconsistent, the five-level threshold only represents the screening level of the current batch.

[0528] 9.3 Natural Breakpoint Method for Partitioning

[0529] As a subordinate implementation method, the natural breakpoint method can be used for... Risk levels are classified. The natural breakpoint method is based on... The natural clustering characteristics of the numerical distribution determine the grading threshold, ensuring that data within the same grade are clustered together. The differences are small, between different levels The differences are significant.

[0530] The natural breakpoint method is applicable within the study area. In cases of uneven numerical distribution, high-risk areas show obvious clustering or skewed distribution. For industrial parks with concentrated pollution sources, obvious pollution hotspots, or prominent migration channels, the natural discontinuity method can effectively highlight relatively high-risk patches.

[0531] 9.4 Quantile Method

[0532] As a subordinate implementation method, the quantile method can be used to... Risk levels are classified. The quantile method is used according to... The numerical ranking results divide the study area grid into several risk levels with similar proportions.

[0533] For example, it can be After sorting, the areas were divided into low-risk, low-to-medium-risk, medium-risk, medium-to-high-risk, and high-risk zones based on 20%, 40%, 60%, and 80% quantiles. The quantile method is suitable for situations where it is necessary to identify relatively high-risk areas within the study area, ensure a relatively balanced area ratio across different risk levels, or for initial screening of sampling sites.

[0534] It should be noted that the quantile method yields a relative risk level, which is suitable for expressing the ranking of risk levels within a study area, but should not be used as the sole basis for comparing absolute risk levels in different parks or at different times.

[0535] 9.5 Standard Deviation Method for Division

[0536] As a subordinate implementation method, the standard deviation method can be used to... Risk levels are classified. The standard deviation method is used to classify risks. Based on the mean and standard deviation, The degree to which something is above or below the average level is used as the basis for grading.

[0537] For example, let The mean is The standard deviation is Then it can be based on and The degree of deviation is used to classify risk levels. Significantly higher than... The area is divided into medium- and high-risk areas or high-risk areas, close to The area was classified as a medium-risk area, significantly lower than... The area is divided into medium- and low-risk areas or low-risk areas.

[0538] Standard deviation method is applicable The distribution is approximately continuous, and the goal is to identify areas that significantly deviate from the average risk level.

[0539] 9.6 Threshold Classification Method Based on Management Objectives

[0540] As a subordinate implementation method, when there is historical regulatory data, existing monitoring sections, park management requirements, or upper limits of sampling resources, a threshold division method based on management objectives can be adopted.

[0541] This method is based on the proportion of the area to be focused on, the number of candidate sampling points, the quota for different sample types, and the internal structure of the study area. Distribution, and determination of low risk, medium risk, medium-high risk, and high risk thresholds.

[0542] In one implementation, it is possible to select an area where the proportion of high-risk zones is controllable and can cover the main high-value patches. The value is used as a high-risk threshold; a value is selected that can cover the main candidate sampling area without excessively expanding the range. The value is used as a medium-to-high risk threshold; based on the background reference area or low value area. The distribution determines the low-risk threshold.

[0543] The threshold division method based on sampling point requirements can be expressed as:

[0544]

[0545] in, This indicates the risk level threshold obtained through screening. This represents an evaluation function constructed with the objectives of sampling resource constraints, the proportion of risk area, and the coverage of high-value patches.

[0546] 9.7 Comprehensive Grading Method

[0547] As a subordinate implementation method, a comprehensive classification method can be used to classify risk levels. The comprehensive classification method can first use a fixed threshold method to form an initial level, and then combine the natural breakpoint method, quantile method and sampling resource constraints to adjust the local thresholds.

[0548] For example, when the high-risk area delineated by the fixed threshold method is too small and may miss candidate hotspots, the quantile method can be combined to expand the candidate range of medium- and high-risk areas; when the high-risk area delineated by the natural discontinuity method is too fragmented, it can be corrected by combining the minimum patch area and neighborhood consistency; when some low-risk areas are located in key migration channels or downstream or downwind paths of pollution sources, it can be adjusted according to... High values ​​are considered as sampling candidate areas.

[0549] The comprehensive classification method is suitable for situations where industrial parks have complex pollution sources, diverse remote sensing anomaly types, and require the selection of priority sampling areas with limited sampling resources.

[0550] 9.8 Risk Level Spatial Adjustment

[0551] To make the risk zoning results more in line with the sampling point layout requirements of industrial parks, this invention can spatially correct the initially classified risk levels.

[0552] In one implementation, isolated high-risk grids that are too small, have a fragmented shape, and lack support from pollution sources or migration channels can have their risk level reduced or be marked as low-priority candidate areas.

[0553] As a subordinate implementation method, for medium-risk grids located inside high-risk patches, surrounded by high-risk areas, or continuously distributed with high-risk areas, their risk level can be increased according to neighborhood relationships to maintain the spatial continuity of risk patches.

[0554] As a subordinate implementation method, for grids located downstream of pollution sources, along sewage channels, at ditch confluences, river inlets, low-lying catchment areas, or downwind channels of pollution sources, even if Even if the high-risk threshold is not reached, it can still be based on High values ​​are considered as sampling candidate areas.

[0555] As a subordinate implementation method, grids whose remote sensing anomalies may be caused by shadows, buildings, construction sites, natural sediment, seasonal water level changes, or non-pollution factors can be marked as uncertain candidate areas, and their priority can be reduced or they can be used as verification candidate points when screening sampling points.

[0556] 9.9 Meaning of Risk Level

[0557] In one implementation, the five risk levels can be interpreted as follows:

[0558] Low-risk areas: The values ​​are low, indicating weak pollution remote sensing anomalies, pollution source influence, and migration pathway conditions. This area can be used as a background control point or a low-priority candidate area.

[0559] Medium and low risk areas: The value is slightly higher than the background level, which may indicate a slight remote sensing anomaly or the influence of a weak pollution source. This area can be used as a general candidate area or a reference for the layout of background gradient sampling points.

[0560] Medium-risk areas: The value is at a moderate level, indicating that there is some pollution anomaly, pollution source influence, or migration channel conditions in the area. This area is suitable as a candidate area for supplementary sampling.

[0561] Medium- and high-risk areas: High values ​​typically indicate abnormal pollution remote sensing activity, high pollution source intensity, or high levels of multiple indicators in the migration pathway. This area should be considered a key candidate area for sampling.

[0562] High-risk areas: The region with the highest value typically exhibits a high degree of media anomaly, strong source term support, and high migration potential. This region should be considered the highest priority candidate sampling area.

[0563] 9.10 Risk Level Results Output

[0564] After completing the risk level classification, the risk level of each spatial grid is assigned to a unified spatial grid, and a comprehensive screening risk level map of the industrial park is generated. This map should include the study area boundary, core area and buffer zone, risk level zones, major pollution sources, river system, migration channels, and necessary legend information.

[0565] In one implementation, the risk level map can serve as the basis for subsequent identification of pollution candidate hotspots, priority sampling site selection, and sampling route design.

[0566] As a subordinate implementation method, risk level results can be overlaid with administrative boundaries, enterprise plots, control units, land use patches, or sampling management units to form risk control zones for environmental management.

[0567] Through the above steps, the comprehensive screening index can be obtained. The results are transformed into risk level results such as low-risk area, medium-low-risk area, medium-risk area, medium-high-risk area, and high-risk area, which provide a foundation for subsequent output of risk zoning results, pollution hotspot results, migration channel results, and priority sampling site results.

[0568] Step S10: Output risk zoning results and candidate sampling point results

[0569] Step S10 includes: extracting risk candidate patches and generating categorized priority sampling points based on sample type, risk level, key migration nodes, background control, minimum spacing, and field constraints.

[0570] Specifically as follows:

[0571] Step S10 is used to output the risk zoning results of industrial parks, candidate risk patches, migration channel results, typed candidate areas and priority sampling point results based on the soil / surface pollution risk indicators, water pollution risk indicators, air pollution transport risk indicators, pollution source intensity indicators, pollution migration channel indicators, comprehensive screening index and risk level formed in the aforementioned steps.

[0572] This step does not simply output a single remote sensing anomaly map or a single pollution index map, but rather expresses the results of multi-media remote sensing anomaly identification, spatial relationships of pollution sources, migration channel analysis results, and comprehensive risk levels, forming a result system that can serve the risk screening and sampling point layout of industrial parks.

[0573] 10.1 Comprehensive Screening Risk Level Chart

[0574] Based on the risk level classification results obtained in step S9, a composite pollution risk level map of the industrial park is output. This composite pollution risk level map uses a unified spatial grid as the basic unit of expression, dividing the study area into low-risk, low-to-medium-risk, medium-risk, medium-to-high-risk, and high-risk zones.

[0575] In one implementation, the composite pollution risk level map includes the following:

[0576] 1. Boundaries of the core area, buffer zone, and research area of ​​the industrial park;

[0577] 2. Low-risk areas, low-to-medium-risk areas, medium-risk areas, medium-to-high-risk areas, and high-risk areas;

[0578] 3. Locations of pollution sources such as key enterprises, sewage outlets, storage yards, hazardous waste temporary storage sites, sewage treatment facilities, and emergency pools;

[0579] 4. Rivers, ditches, sewage channels, stormwater runoff channels, and other water systems;

[0580] 5. Necessary reference information such as background reference area, outer verification area, major residential areas, farmland, water bodies and ecological protection objects;

[0581] 6. Legend, scale, north arrow, coordinate system, data source, and mapping date.

[0582] The comprehensive screening risk level map is used to visually represent the spatial distribution of candidate risks of water-soil-air composite pollution in industrial parks. It can serve as a basic map for conducting sampling site selection, candidate area screening, and key area identification.

[0583] 10.2 Remote sensing anomaly map of soil or surface pollution

[0584] Based on the soil / surface pollution risk indicator obtained in step S3 Output a remote sensing anomaly map of soil pollution. This map is used to indicate the degree of soil / surface pollution risk in bare land in industrial parks, around storage yards, around factory boundaries, road dust deposition areas, around ditches, riverbanks, low-lying catchment areas, and other exposed surface areas.

[0585] In one implementation, the remote sensing anomaly map of soil pollution includes:

[0586] 1. Continuous value distribution;

[0587] 2. Low anomaly, low-to-medium anomaly, medium anomaly, medium-to-high anomaly, and high anomaly levels;

[0588] 3. Areas with effective soil observation and areas without effective soil observation;

[0589] 4. Key pollution sources, storage yards, roads, ditches, riverbanks, and historically polluted sites;

[0590] 5. Candidate anomalous patches for subsequent on-site soil sampling.

[0591] This map can be used to help identify potential soil pollution hotspots, dust deposition areas, landfill impact areas, and areas where soil sampling points should be prioritized.

[0592] 10.3 Remote sensing anomaly map of water pollution

[0593] Based on the water pollution risk indicator obtained in step S4 Output a remote sensing anomaly map of water pollution. This map is used to express the degree of anomalies in water color, turbidity, suspended solids, organic pollution, or industrial emissions affecting rivers, ditches, sewage channels, stormwater runoff channels, ponds, emergency ponds, sedimentation ponds, downstream water bodies, and other surface water bodies.

[0594] In one implementation, the remote sensing anomaly map of water pollution includes:

[0595] 1. Continuous value distribution;

[0596] 2. Water body anomaly levels: low, low-medium, medium, medium-high, and high.

[0597] 3. Boundaries of rivers, ditches, sewage channels, rainwater runoff channels, and water bodies;

[0598] 4. Upstream reference water section, downstream abnormal water section, inlet and outlet adjacent abnormal areas;

[0599] 5. Locations of water bodies that require priority for water quality sampling or drone verification.

[0600] This map can be used to help identify sewage channels, ditches, downstream river sections, and suspected polluted water bodies in industrial parks, providing a basis for water quality sampling and sewage outlet verification.

[0601] 10.4 Air Pollutant Transport Risk Map

[0602] Based on the air pollution transport risk indicators obtained in step S5 Output an atmospheric pollution transport risk map. This map is used to express the potential impacts of atmospheric pollution sources in industrial parks on the transport, diffusion, or deposition of atmospheric pollution in each spatial grid under the combined effects of remote sensing anomalies, wind direction and speed, distance attenuation, and downwind relationship.

[0603] In one implementation, an air pollution transport risk map includes:

[0604] 1. Continuous value distribution;

[0605] 2. Air pollution transport zones are categorized into low-risk, low-medium-risk, medium-risk, medium-high-risk, and high-risk areas;

[0606] 3. Potential sources of exhaust gas emissions, chimneys, areas of unorganized emissions, dust storage areas, road dust areas, and areas with thermal anomalies;

[0607] 4. Dominant wind direction, downwind influence sector, wind frequency, or wind rose diagram;

[0608] 5. Downwind subsidence impact area, road dust impact area, outer verification area, and background reference gradient location.

[0609] This map can be used to help identify downwind areas affected by air pollution, areas affected by dust deposition, and areas where candidate sampling points for atmospheric deposition / road dust / particulate matter need to be established.

[0610] 10.5 Distribution map of pollution source intensity

[0611] Based on the pollution source intensity index obtained in step S6 Output a pollution source intensity distribution map. This map is used to express the degree of influence of industrial enterprises, sewage outlets, storage yards, hazardous waste temporary storage areas, sewage treatment facilities, emergency ponds, road transport channels, and historically contaminated sites on the spatial grid.

[0612] In one implementation, the pollution source intensity distribution map includes:

[0613] 1. Continuous value distribution;

[0614] 2. Pollution source intensity level;

[0615] 3. Enterprise type, boundaries of key enterprises, sewage outlets, storage yards, hazardous waste temporary storage sites, emergency pools, and sewage treatment facilities;

[0616] 4. The impact range of the pollution source, the distance attenuation zone, and the area of ​​proximity correction for remote sensing anomalies;

[0617] 5. Key pollution sources and high-value areas of source term impact that require priority attention.

[0618] This map can be used to help identify the areas affected by major source terms in industrial parks, providing spatial clues for the selection of candidate sampling points, the focus on key source terms, and the formulation of sampling tasks.

[0619] 10.6 Pollution Migration Channel Map

[0620] Based on the pollution migration channel indicators obtained in step S7 Output a pollution migration channel map. This map is used to express the spatial channels through which pollutants may migrate along surface runoff, river systems, ditches, sewage channels, stormwater runoff channels, low-lying catchment areas, prevailing wind directions, downwind diffusion paths, and particulate matter deposition paths.

[0621] In one implementation, the pollution migration pathway map includes:

[0622] 1. Continuous value distribution;

[0623] 2. Hydrological migration channels, low-lying catchment areas, confluence paths, river inlets, and ditch junctions;

[0624] 3. Atmospheric migration channels, prevailing wind paths, downwind influence areas, and subsidence influence areas;

[0625] 4. Potential migration paths from pollution sources to downstream water bodies, ditch inlets, downwind subsidence zones, and surrounding affected areas;

[0626] 5. Key migration nodes, channel intersection nodes, and high-value areas of migration channels where candidate sample points need to be deployed.

[0627] This map can be used to help identify key pathways for pollutants to migrate from the source area to downstream, downwind, and surrounding impact areas, avoiding the omission of migration channel nodes due to field sampling only being placed around the pollution source.

[0628] 10.7 Candidate Hotspot Map of Complex Pollution

[0629] based on High value area High value area High value area High value area High value area and High-value areas are identified, and candidate hotspots of compound pollution are generated, and a candidate hotspot map of compound pollution is output.

[0630] In one implementation, hotspots of complex pollution can be determined by the following rules:

[0631] 1. Located in a high-risk or medium-to-high-risk area;

[0632] 2. The simultaneous existence of two or more types of high-value indicators, for example... High value and High value superposition, High value and High value superposition, High value and or High values ​​superimposed;

[0633] 3. Located in the vicinity of the pollution source, at a key node of the migration channel, in the peripheral verification area, or at the location of the background reference gradient;

[0634] 4. It exhibits spatial consistency with historical pollution records, existing on-site anomaly records, or remote sensing multi-temporal anomalies;

[0635] 5. High-value patches reach the set threshold area, or although the area is small, they are located in key locations such as sewage outlets, inlets, ditch intersections, and downwind settlement zones.

[0636] Composite pollution candidate hotspot maps are used to identify key areas where sampling points should be prioritized. Unlike single pollution anomaly maps, this map focuses on representing comprehensive candidate hotspots formed by the overlay of multiple indicators.

[0637] 10.8 Results of Priority Sampling Point Layout

[0638] Based on a comprehensive screening of risk levels, candidate hotspots, migration pathways, and data gaps, priority sampling site layout results are generated. These priority sampling site layout results are used to guide the determination of candidate sites for soil, water bodies, sediments, road dust or settling, and background control samples.

[0639] Priority sampling points are not randomly selected, nor are they selected solely around known pollution sources. Instead, they are determined based on the comprehensive screening results of risk zoning, following the principles of prioritizing high-risk areas, migration pathways, anomaly overlap, background comparison, and spatial discrete constraints.

[0640] In one implementation, priority sampling points include the following types:

[0641] 1. High-risk area sampling points: deployed in High-value areas and high-risk areas are used to prioritize the acquisition of samples from key candidate areas;

[0642] 2. Pollution candidate hotspot sampling points: deployed in , , , or High-value areas of multiple indicators are superimposed to cover candidate hotspots supported by multiple factors;

[0643] 3. Sampling points near pollution sources: These are set up around key enterprises, sewage outlets, storage yards, hazardous waste temporary storage areas, accident pools, and sewage treatment facilities to verify the impact of pollution sources;

[0644] 4. Migration channel sampling points: deployed at key migration nodes such as ditches, sewage channels, river inlets, low-lying catchment areas, downstream paths, and downwind subsidence zones;

[0645] 5. Peripheral verification sampling points: deployed around the CORE, in the extension area of ​​the migration channel, or in the downwind influence area to cover the peripheral candidate influence range;

[0646] 6. Background control sampling points: located far from pollution sources, far from migration routes, and in low-risk areas. Areas with comparable surface conditions to the study area are used to establish background references;

[0647] 7. Uncertainty sampling points: These are placed in areas with poor remote sensing image quality, missing indicators, or significant changes in risk level to improve the representativeness of the sampling scheme.

[0648] 10.9 Priority Sampling Point Selection Rules

[0649] Construct a sampling priority index for each spatial grid. : .

[0650] in, The intensity of multiple indicator hotspots superimposed, Due to missing data or spatial uncertainty, For on-site accessibility, all three are normalized to the interval [0, 1] and The larger the value, the easier it is to reach; to The weights are non-negative and their sum is 1.

[0651] Preferably, candidate grid sets are established according to sample types such as soil, water, sediment, road dust or settling, and background control, and the following screening is performed on each set:

[0652] (1) According to Sort the grids from highest to lowest and select the grids that have the highest current scores and meet the on-site safety and land access requirements;

[0653] (2) After recording the selected grid, remove candidate grids whose distance to it is less than the minimum spacing of the corresponding sample type;

[0654] (3) Repeat the selection and elimination process until the preset number of samples of that type is reached or the candidate set is empty;

[0655] (4) Water samples should cover at least the upstream reference, the vicinity of potential emission nodes, and the downstream impact location; soil samples should cover the vicinity of the source area and the migration path; atmospheric or deposition samples should cover the upwind background, the vicinity of the source area, and the downwind location.

[0656] (5) Other low-risk background control sites that are far from pollution sources and migration routes and have comparable surface conditions are reserved.

[0657] The candidate set for each category, the number of sample points, the minimum spacing, and the gradient constraints are all stored as parameters, making the sampling point selection process recalculated.

[0658] 10.10 Presentation of Sampling Point Layout Results

[0659] Priority sampling results may include a sampling point location map, a sampling point location table, a sampling priority list, and a sampling task list.

[0660] In one implementation, the sampling point table includes at least the following fields:

[0661] 1. Sample point numbering;

[0662] 2. Sampling point types, including soil sampling points, water sampling points, sediment sampling points, road dust or sediment sampling points, peripheral verification points or background control points;

[0663] 3. Latitude and longitude or plane coordinates of the sample point;

[0664] 4. Risk level;

[0665] 5. Grid number;

[0666] 6. Value and , , , , Indicator value;

[0667] 7. Reasons for setting up sample points;

[0668] 8. Information on nearby pollution sources and migration routes;

[0669] 9. Suggested sample types or categories of indicators to focus on;

[0670] 10. On-site accessibility and safety precautions.

[0671] As a subordinate implementation method, different thematic sampling maps can be generated according to the sampling point type, including soil priority sampling map, water priority sampling map, sediment or road dust map, peripheral verification point map and background control point map.

[0672] 10.11 Output Format

[0673] The above results can be output in the form of maps, tables, vector files, raster files, databases, or geographic information services.

[0674] In one implementation, map results can be output as TIFF, GeoTIFF, PNG, JPEG, PDF or other formats; vector results can be output as Shapefile, GeoJSON, GPKG, KML or other formats; tabular results can be output as Excel, CSV or database tables; spatial database results can be stored in PostGIS, File Geodatabase, GeoPackage or other geospatial databases.

[0675] As a subordinate implementation method, the comprehensive screening risk level map, pollution candidate hotspot map, migration channel map and priority sampling point map can be connected to the industrial park environmental supervision platform, remote sensing monitoring platform or geographic information system for results display, query statistics and sampling task compilation.

[0676] Through the above steps, the present invention can transform the results of remote sensing collaborative screening and comprehensive risk index calculation into operable spatial results, realizing the transformation from "remote sensing anomaly identification" to "risk zoning" and "priority sampling point determination".

[0677] It should be noted that this method ends with the classification of candidate sampling points and the corresponding chart output; contamination confirmation, laboratory testing evaluation, or parameter inversion are not necessary steps in this main line.

[0678] The foregoing , and These are all risk indicators and do not constitute independent confirmation of pollution. It is used for spatial priority ranking and does not replace statutory monitoring and evaluation for specific pollutants.

[0679] The output of this invention is a risk candidate area, a risk level map, a categorized candidate area, and a priority sampling point. It does not directly output pollutant concentration, exceedance judgment, or detection calibration results.

[0680] The key technical features of this invention are described below.

[0681] The core of this invention lies in transforming multi-media remote sensing anomalies from independent thematic maps into comprehensive screening results constrained by pollution sources and directional migration channels, and further generating candidate sampling points according to categorized candidate sets, spatial spacing, and gradient constraints. Specific technical points are as follows.

[0682] 1. A unified grid representation with validity and credibility

[0683] Data with different spatial resolutions, time windows, and applicable media are mapped to a unified grid, and the validity and reliability of indicators are recorded simultaneously. For inapplicable media, effective weights are used for recalibration to avoid mistaking "no observations" for "no risk."

[0684] 2. Layered coupling of source term—medium anomaly—migration path

[0685] First, risk indicators from soil or surface, water, and atmosphere are used to form multi-media anomalies. These anomalies are then coupled with the relationships between pollution source intensity and migration channels to form source term support, migration, and comprehensive screening indices. This structure differs from a simple weighted overlay of multiple layers.

[0686] 3. Expression of migration risk constrained by wind direction and hydrological direction

[0687] By utilizing the consistency between the pollution source-grid azimuth and the time-period wind direction, distance attenuation, DEM confluence direction, water system topology, and upstream-downstream relationships, directional constraints are imposed on atmospheric and hydrological migration, reducing spatial misjudgments caused by diffusion based solely on Euclidean distance.

[0688] 4. Conversion from risk zoning to categorized sampling task

[0689] By combining composite risk, single-medium anomaly, critical migration nodes, background control, sample point type quota, minimum spacing, and on-site accessibility for sampling point screening, remote sensing risk results can be directly used to generate soil, water, sediment, road dust or settling and background sample point tasks.

[0690] The beneficial effects of this invention are as follows:

[0691] 1. Reduce the direct superposition error of multi-source heterogeneous data. By unifying the grid, validity identification, credibility, and effective weight recalibration, data from different media and scales can be compared within the same computational unit, and missing values ​​are avoided from being mistaken for low-risk data.

[0692] 2. Enhance the process interpretability of risk candidate regions. High-risk values ​​can be decomposed into contributions from media anomalies, pollution sources, and migration pathways, facilitating the determination of which factors drive the candidate region and providing a basis for selecting sampling point types.

[0693] 3. Narrowing the candidate sampling range. By coupling source terms, medium anomalies, and migration channels, directional migration constraints, and key patch extraction, areas requiring sampling are prioritized within a larger campus area, reducing omissions and duplications caused by relying solely on experience-based sampling.

[0694] 4. Improve the representativeness of sampling points. Classified sampling points simultaneously cover composite risk areas, single-medium anomaly areas, migration channel nodes, peripheral verification areas, and background areas, and reduce sampling point clustering through minimum spacing and spatial discretization rules.

[0695] 5. Supports result tracking and sampling task verification. The calculation process saves data time windows, spatial references, normalization parameters, weights, thresholds, sample point type quotas, and minimum intervals, enabling the risk zoning and candidate sampling point generation process to be verified, recalculated, and traced.

[0696] 6. Define the deliverables and interfaces for subsequent work. This invention not only outputs risk maps but also simultaneously outputs candidate patches, candidate areas categorized by type, candidate sampling points, sample point types, coordinates, priorities, and rationales for their placement. This allows remote sensing screening results to be directly converted into base maps, point locations, and task lists for on-site sampling plans. Subsequent work can then be based on this information to plan sampling routes, sample types, personnel assignments, access coordination, and background control point setup, thereby reducing the uncertainty of repeated pre-sampling reconnaissance and empirical point selection.

[0697] The above effects all point to the optimization of risk screening and sampling tasks; without confirmation by statutory monitoring or laboratory testing, no claims shall be made regarding the accuracy of pollution concentration inversion, environmental quality compliance status, or conclusions of pollution exceeding standards.

[0698] The following is a brief description of the accompanying drawings used in the embodiments of the present invention. It should be understood that the following drawings are only used to illustrate the technical approach, data relationships, indicator construction, and result output methods of the present invention, and do not constitute a limitation on the scope of protection of the present invention.

[0699] like Figure 1 As shown, this invention focuses on industrial parks and their surrounding impact areas. First, it determines the spatial organization relationship of STUDY_AREA, PARK, CORE, and BUFFER, and constructs a unified spatial grid. Second, it acquires and preprocesses multi-source remote sensing data and auxiliary geospatial data. Then, it constructs soil / surface pollution risk indicators, water pollution risk indicators, air pollution transport risk indicators, pollution source intensity indicators, and pollution migration channel indicators, respectively. Furthermore, it constructs a comprehensive screening index and classifies risk levels. Finally, it outputs risk zoning results and priority sampling points by type. Figure 1 The technical endpoint of the process shown is to output risk zones, candidate patches, categorized candidate sampling points, and charts; this process does not include post-sampling detection, pollution confirmation, or model updates.

[0700] like Figure 2 As shown, data from various sources, scales, and types, including multi-source remote sensing images, atmospheric pollution remote sensing products, DEM topographic data, wind field data, river system data, enterprise spatial data, land use data, and background reference data, are uniformly mapped onto a spatial grid within the study area. Each grid cell can obtain remote sensing image features, atmospheric pollution features, pollution source attributes, topographic and hydrological attributes, wind field attributes, and land use attributes, providing a unified spatial computational foundation for subsequent comprehensive screening index construction. Furthermore, Depend on , , , , The results are generated through normalization, confidence weighting, and coupled calculations; background reference and accessibility information are used for candidate sampling point selection.

[0701] like Figure 3 As shown, this invention constructs remote sensing anomaly identification processes for three types of environmental media: soil, water, and atmosphere. Soil pollution remote sensing anomaly identification primarily targets areas such as bare land, storage yards, factory boundaries, road dust deposition areas, ditches, and riverbanks; water pollution remote sensing anomaly identification primarily targets areas such as rivers, ditches, sewage channels, stormwater runoff channels, ponds, and downstream water bodies; and atmosphere pollution transport risk identification mainly combines remote sensing products for NO2, SO2, and AOD, the location of potential emission sources, and wind field data. The three types of remote sensing anomaly indicators are respectively formulated... , and They will also participate in subsequent comprehensive screening and evaluation.

[0702] like Figure 4 As shown, industrial enterprises, sewage outlets, storage yards, hazardous waste temporary storage sites, sewage treatment facilities, emergency ponds, road transport routes, and historically contaminated sites are considered potential pollution sources; soil or surface anomalies, water body anomalies, and atmospheric pollution transport anomalies are considered medium anomalies; surface runoff paths, river systems, ditches, sewage channels, low-lying catchment areas, and prevailing and leeward transport paths are considered pollution migration channels. This invention expresses the spatial relationship between pollution sources, medium anomalies, and migration channels through pollution source intensity indicators and pollution migration channel indicators. Furthermore, background references, information on residential areas / farmland / ecological protection objects, etc., are used for scope delineation, control points, and site layout constraints.

[0703] like Figure 5 As shown, this invention uses soil / surface pollution risk indicators. Water pollution risk indicator air pollution transport risk indicators Pollution source intensity indicators and pollution migration channel indicators After unification and normalization, a multi-media anomaly comprehensive index, source term support sub-index, and migration sub-index are constructed, and a comprehensive screening index is further formed. The aforementioned Used to comprehensively characterize the candidate risk levels of water-soil-air composite pollution for each spatial grid.

[0704] like Figure 6 As shown, this invention is based on a comprehensive screening index. The study area was divided into low-risk, low-to-medium-risk, medium-risk, medium-to-high-risk, and high-risk zones. Furthermore, risk level maps, soil or surface remote sensing anomaly maps, water pollution remote sensing anomaly maps, air pollution transport risk maps, pollution source intensity distribution maps, pollution migration channel maps, candidate hotspot maps for complex pollution, and priority sampling site maps were generated. These results are used to support pollution risk screening and sampling point determination in industrial parks.

[0705] like Figure 7 As shown, after completing remote sensing collaborative identification, comprehensive screening index calculation, and risk zoning output, this invention converts high-risk areas, medium-to-high-risk areas, pollution candidate hotspots, key migration channel areas, peripheral verification areas, and background reference areas into different types of candidate sampling points. Candidate points are screened according to sample type, risk level, minimum spacing, background comparison, and on-site accessibility constraints to form a deliverable sampling point map and sampling point table.

[0706] like Figure 8 As shown, in a specific embodiment, the core area, near-field influence area, and far-field influence area of ​​the industrial park can be delineated based on the industrial park boundary, distribution of key enterprises, water system, roads, surrounding land use, and background reference requirements. This figure is used to illustrate the spatial scope setting method of the present invention in a specific industrial park, as well as the basis for constructing a unified spatial grid.

[0707] The following specific examples will provide further explanation.

[0708] The following example, using an industrial park, illustrates the process of multi-source remote sensing collaborative screening, risk zoning, and priority sampling point placement at the industrial park scale. This embodiment utilizes publicly available remote sensing images, topographic data, wind field data, and park spatial boundaries to generate comprehensive screening indices for land anomalies, water anomalies, pollution source intensity proxies, migration channel proxies, atmospheric transport risk proxies, wind direction correction, and sampling priority indices on a unified 30m spatial grid.

[0709] The purpose of this embodiment is to demonstrate that the present invention can organize multi-source remote sensing anomalies, pollution source impacts, migration channel conditions, atmospheric transport directions, and sampling point deployment requirements into an executable spatial screening process. In this embodiment, HighLand, High Water, PSI_proxy, MTI_proxy, ATR_wind, CPRI_wind, and Sampling PriorityIndex are all remote sensing collaborative screening or sampling priority indicators, and are not equivalent to specific pollutant concentrations, pollution exceedance conclusions, or statutory environmental quality assessment results.

[0710] I. Spatial Partitioning and Functional Positioning in the Implementation Example

[0711] This embodiment adopts a spatial organization method of "STUDY_AREA—PARK—CORE—BUFFER". Here, PARK represents the source area of ​​the industrial park or the concentrated area of ​​industrial activity; CORE represents the core influence area formed by expanding 500m outward from PARK, used to accommodate the main risk zoning and sampling point analysis; STUDY_AREA represents the embodiment study area formed by expanding 3000m outward from PARK; and BUFFER represents the outer verification area of ​​STUDY_AREA after deducting CORE, used to identify potential diffusion, migration, and background control relationships in the periphery.

[0712] This embodiment preserves the spatial partitioning of CORE and BUFFER, enabling the method to perform detailed screening within concentrated industrial activity areas as well as identify verification points and background control points within the peripheral impact area, avoiding the simple exclusion or simple determination of high risk in the peripheral area.

[0713] Figure 9 The diagram illustrates the spatial zoning of the industrial park research area, industrial park source area, CORE, and BUFFER areas in this embodiment. Table 1 shows the spatial zoning and functional descriptions of this embodiment.

[0714] Table 1 Spatial Partitioning and Functional Description of This Embodiment

[0715]

[0716] II. Data Sources and Processing Environment

[0717] This embodiment uses multi-source remote sensing and auxiliary geospatial data as input, and employs a remote sensing processing platform and geospatial computing environment to complete data filtering, projection transformation, cropping, normalization, raster overlay, and vectorization output. To avoid limiting this invention to a specific computing platform, the processing can be implemented in local remote sensing processing software, remote sensing cloud platforms, unit servers, dedicated environmental monitoring platforms, or other equivalent computing environments.

[0718] The core deliverables of this embodiment include index results and sampling point distribution results. The index results output CPRI_wind, PSI_proxy, MTI_proxy, ATR_wind, and key threshold statistics; the sampling point distribution results output categorized candidate regions, Sampling Priority Index, and candidate sampling points. All major maps are represented using the WGS 84 / UTM 49N (EPSG:32649) projected coordinate system and a 30m spatial grid.

[0719] Table 2 shows the main deliverables and usage descriptions of this embodiment.

[0720] Table 2. Main Deliverables and Usage Descriptions of This Embodiment

[0721]

[0722] III. Preliminary Screening of Remote Sensing Anomalies in Land and Water Bodies

[0723] Within a unified spatial grid, two types of remote sensing preliminary screening results are first generated: Land_Anomaly and Water_Anomaly. Land_Anomaly mainly reflects land anomaly features such as exposed surface, construction disturbance, road dust deposition, vegetation stress, and proximity to industrial activities; Water_Anomaly mainly reflects water-related features such as turbidity, water color, water body edges, and anomalies in local drainage channels.

[0724] It should be noted that both land anomalies and water anomalies are remote sensing risk indicators. Specifically, high anomaly areas on land are not equivalent to soil heavy metal content inversion results, and high anomaly areas in water are not equivalent to legally mandated water quality exceedance results. They are used to narrow down the candidate sampling range and serve as inputs for subsequent CPRI_wind and sampling priority indices.

[0725] Figure 10 The example shows a preliminary remote sensing map of land and water anomalies in an industrial park. Table 3 shows statistical data on high anomaly areas in land and water.

[0726] Table 3 Statistical Table of High Anomaly Zones in Land and Water Bodies

[0727]

[0728] IV. Construction of Risk Proxy Layer for Pollution Source Intensity, Migration Channels, and Atmospheric Transport

[0729] In this embodiment, PSI_proxy characterizes the intensity of industrial activities and potential pollution source impacts, MTI_proxy characterizes terrain depressions, low slopes, water channels, and potential migration paths, and ATR_wind characterizes the atmospheric transport risk after the atmospheric pollution remote sensing anomalies are corrected for the mean wind field direction. These three indices are not interchangeable but rather respectively characterize source terms, channels, and atmospheric transport conditions.

[0730] Figure 11 The data shows that high PSI_proxy values ​​are mainly concentrated in PARK and its adjacent areas, high MTI_proxy values ​​are distributed along potential confluence and channel structures, and high ATR_wind values ​​exhibit a directional spatial gradient after wind direction correction. These results indicate that a single industrial activity intensity, a single topographic channel, or a single atmospheric anomaly is insufficient to independently explain the risk of complex pollution; a comprehensive screening index needs to be formed through multi-factor coupling.

[0731] Figure 11The example diagram illustrates the pollution source intensity, migration pathways, and atmospheric transport risk proxy map for the industrial park. Table 4 shows the statistical data for key thresholds of PSI, MTI, and ATR_wind.

[0732] Table 4. Statistics of Key Thresholds for PSI, MTI, and ATR_wind

[0733]

[0734] V. Construction and Zoning of the Comprehensive Screening Index CPRI_wind for Wind Direction Correction

[0735] Land anomalies, water anomalies, pollution source intensity proxies, migration channel proxies, and wind-direction-corrected atmospheric transport risks are coupled in a unified grid to generate a wind-direction-corrected comprehensive screening index, CPRI_wind. In this embodiment, CPRI_wind is used as a preliminary comprehensive screening index for risk candidate area ranking and sampling priority design, and is not presented as a pollution exceedance assessment or pollutant concentration confirmation result.

[0736] In this embodiment, CPRI_wind is partitioned using five threshold levels: 0.00-0.20, 0.20-0.40, 0.40-0.60, 0.60-0.80, and 0.80-1.00. CPRI_wind > 0.60 is used as a candidate threshold for medium- to high-risk conditions. Figure 12 The red border indicates the primary candidate patches with CPRI_wind > 0.60.

[0737] Figure 12 The example shows a risk zoning map of the industrial park using the CPRI_wind comprehensive screening index (wind direction correction). Table 5 shows the area statistics for CPRI_wind risk classification.

[0738] Table 5. CPRI_wind Risk Classification Area Statistics Table

[0739]

[0740] Statistical results show that the average CPRI_wind value within the CORE is 0.3358, and within the BUFFER it is 0.1323. The area of ​​medium-to-high risk candidate zones with CPRI_wind > 0.60 within the CORE is 0.2437 km², accounting for 1.84% of the CORE. No areas with CPRI_wind > 0.60 were identified within the BUFFER, and no high-risk areas with CPRI_wind > 0.80 were identified in either the CORE or the BUFFER. These results indicate that this embodiment can narrow down the key candidate zones to localized patches within concentrated industrial activity areas, while avoiding classifying the entire surrounding area as high-risk.

[0741] VI. On-site sampling site layout guided by comprehensive screening results

[0742] Based on the CPRI_wind risk zoning and categorized candidate regions, this embodiment further constructs a Sampling Priority Index and generates candidate sampling points by combining sample type, spatial zoning, minimum spacing, and background control requirements. Sampling points are divided into core points for soil / surface dust within the CORE, core points for water / sediment, and core points for atmospheric deposition / road dust, as well as verification points and background control points within and around the BUFFER. Categorized candidate regions are used to constrain the spatial sources of different sample types, and the Sampling Priority Index is used to further rank and generate candidate sampling points within the corresponding candidate regions. Therefore, candidate sampling points are the result of the combined effect of categorized candidate regions, comprehensive screening index, and field constraints.

[0743] Figure 13 The example shows the sampling priority index and categorized candidate sampling point map for the industrial park. Table 6 shows the statistical data on the types and quantities of candidate sampling points.

[0744] Table 6. Statistics on the Types and Quantities of Candidate Sampling Points

[0745]

[0746] This embodiment outputs a total of 180 candidate sampling points. The CORE includes 40 core points for soil / surface dust, 40 core points for water / sediment, and 30 core points for atmospheric deposition / road dust. The BUFFER includes 40 peripheral verification points and 30 background control points. The average Sampling Priority Index for the CORE is 0.6436, higher than the 0.3641 for the BUFFER. This result indicates that sampling priority is mainly concentrated in the core impact area of ​​the industrial park, while retaining peripheral verification and background calibration capabilities.

[0747] The above sampling points are candidate locations. Before actual sampling, adjustments should be made based on factors such as site accessibility, safety conditions, land ownership, enterprise access requirements, weather and hydrological conditions, and regulatory needs. Figure 13 This indicates the sampling site layout plan, but does not mean that on-site sampling has been completed or that pollutant detection results have been obtained.

[0748] VII. Analysis of Key Risk Candidate Regions and Classified Candidate Regions

[0749] To illustrate the controlling factors for medium- and high-risk candidate patches, this embodiment overlays regions with CPRI_wind > 0.60 with high PSI, high source term influence, high MTI, High Land, and High Water. The results show that regions within the CORE with CPRI_wind > 0.60 spatially overlap with high PSI, high source term influence, and High Land, while direct overlap with high MTI and High Water is weak. This indicates that the main medium- and high-risk candidate patches identified in this version are more inclined towards regions supported by industrial source terms and coupled with land disturbances.

[0750] Furthermore, the candidate regions are divided into five categories: terrestrial source region, hydrological source region, atmospheric source region, peripheral verification, and background control. These categorized candidate regions are not mutually exclusive; the same spatial grid can be influenced by multiple types of factors simultaneously. Any_Typed_Key represents the union of any type of candidate region.

[0751] Figure 14 A spatial distribution map of candidate zones for different types of industrial parks is shown in the example. Table 7 shows the area statistics of the candidate zones for different types.

[0752] Table 7. Statistics on the Area of ​​Candidate Areas by Category

[0753]

[0754] The categorized statistics show that the Any_Typed_Key area within the CORE is 11.8183 km², accounting for 89.01% of the CORE; among which, the candidate areas for land source terms and atmospheric source terms account for 58.78% and 38.54% of the CORE, respectively. Within the BUFFER, the Any_Typed_Key area is 9.2300 km², accounting for 16.59% of the BUFFER; the background control candidate area is 17.0151 km², accounting for 30.59% of the BUFFER. This indicates that the main function of the BUFFER is peripheral verification and background calibration, rather than being interpreted as a high-risk area as a whole.

[0755] VIII. Interim Conclusions of this Embodiment

[0756] Based on the figures and statistical results of this embodiment, the following preliminary conclusions are reached:

[0757] First, this invention can simultaneously express land anomalies, water anomalies, pollution source intensity, migration channels, and atmospheric transport risks on a unified 30m spatial grid, forming verifiable multi-source remote sensing collaborative screening results.

[0758] Second, the high-risk candidate areas in CPRI_wind were narrowed down to local patches within the CORE, and no areas with CPRI_wind > 0.60 were identified within the BUFFER. This indicates that the method can distinguish between the core influence area and the peripheral review area, thus avoiding the overall high-risk classification of the peripheral area.

[0759] Third, the classification of candidate regions and sampling priority index can convert risk partitioning results into executable candidate sampling points.

[0760] Fourth, neither CORE nor BUFFER showed a high-risk area with CPRI_wind > 0.80, indicating that the risk expression in this embodiment is relatively restrained, which is in line with the technical positioning of remote sensing pre-screening and sampling point determination.

[0761] The purpose of this embodiment is to prove that the method chain is executable, the results can be mapped, the statistics can be verified, and the locations can be converted, which can support the technical main line of this invention: "multi-source remote sensing collaborative screening - risk zoning - candidate area extraction - priority sampling point layout".

[0762] IX. Applicable Conditions and Boundaries of this Embodiment

[0763] This embodiment is a preliminary screening version based on publicly available remote sensing, topographic, wind field, and spatial boundary data. High-value areas of CPRI_wind, candidate areas for terrestrial sources, candidate areas for hydrological sources, candidate areas for atmospheric sources, peripheral verification candidate areas, and background control candidate areas are all used as the basis for determining candidate sampling points, but not as the basis for confirming pollution exceedances or pollutant concentrations.

[0764] The technical endpoint of this embodiment is the classification of candidate sampling points and their graphical output. Formal sampling implementation, sample testing, contamination confirmation, evaluation of test results, and parameter inversion are subsequent tasks and are not included in the main technical line of this embodiment.

[0765] The study area, PARK extension distance, 30m grid scale, 0.60 and 0.80 thresholds, categorized candidate area rules, and 180 candidate sampling points in this embodiment are specific implementation conditions of the industrial park in this embodiment and do not constitute a limitation on the scope of protection of this invention. In other industrial parks or under different data conditions, the study area, grid scale, index weights, thresholds, and number of sampling points can be adjusted according to the park area, pollution source type, prevailing wind direction, hydrological conditions, background reference requirements, on-site accessibility, and the quantity of sampling resources.

[0766] As can be seen from this embodiment, the present invention can organize multi-source remote sensing anomalies, pollution source impacts, migration channel conditions, and sampling constraints into verifiable spatial results before on-site sampling. These results include not only the CPRI_wind risk zoning map, but also candidate risk patches, categorized candidate areas, sampling priority indices, candidate sampling point maps, and candidate sampling point tables.

[0767] In summary, this invention completes the conversion from multi-source remote sensing data and auxiliary spatial data into executable candidate sampling task drafts. Specifically, it includes: uniformly mapping remote sensing anomalies related to water bodies, soil / surface, and atmospheric transport to the same spatial grid; incorporating pollution source intensity and migration channel conditions into a comprehensive screening index to avoid judging risk based solely on a single remote sensing anomaly; converting continuous indices into risk levels and candidate patches to clarify the spatial locations requiring priority attention; further converting the risk zoning results into categorized candidate sampling points such as soil / surface dust, water / sediment, atmospheric deposition / road dust, peripheral verification, and background control; and generating exportable maps, tables, and spatial data, enabling the screening process and the basis for point generation to be verified.

[0768] This invention provides the following interfaces for subsequent work: On-site sampling plan preparation can directly determine sampling locations, sample types, priorities, and rationales based on the candidate sampling point list; on-site reconnaissance can focus on candidate patches and migration pathways for verification, reducing the scope of untargeted reconnaissance; sampling organization can arrange personnel, routes, access coordination, and sample container preparation according to the categorized sampling points; laboratory testing can pre-classify sample types such as soil, water, sediment, road dust, or settling based on the sample point type; and regulatory or scientific research verification can retain the corresponding CPRI, TSAI, WAI, ATR, PSI, and MTI indicators for each point.

[0769] The convenience of this invention lies in the fact that, before formal sampling is implemented, a wide range of dispersed and heterogeneous risk clues are converged into a limited number of candidate sites with clear types and traceable reasons. This provides a spatial prior and task organization basis for subsequent sampling, sample testing, and contamination confirmation. Subsequent sampling and testing evaluation can be carried out based on the output results of this invention, but this is not the technical endpoint of this invention in this embodiment.

[0770] The embodiments described above are merely preferred embodiments of the present invention and are not exhaustive examples of all possible implementations of the present invention. Any obvious modifications made by those skilled in the art without departing from the principles and spirit of the present invention should be considered to be included within the scope of protection of the claims of the present invention.

Claims

1. A method for multi-source remote sensing collaborative screening, risk zoning, and sampling point layout for water-soil-air complex pollution in industrial parks, characterized in that, Includes the following steps: S1. Determine the core area and the area of ​​influence based on the boundaries of the industrial park, the concentrated area of ​​industrial activities, the scope of hydrological influence, and the prevailing wind direction, and construct a unified spatial grid; S2. Map the multi-source remote sensing and spatial data to the spatial grid; S3. Generate the Soil / Surface Pollution Risk Indicator (TSAI), the Water Pollution Risk Indicator (WAI), and the Atmospheric Pollution Transport Risk Indicator (ATR) on the spatial grid, respectively. S4. Generate the pollution source intensity index PSI based on pollution source information, and generate the pollution migration channel index MTI based on topography, water system and wind field; S5. Normalize and weight the TSAI, WAI, ATR, PSI and MTI, construct the comprehensive screening index CPRI according to "media anomaly - source term support - migration conditions" and classify the risk level; S6. Extract risk candidate patches and generate categorized priority sampling points based on sample type, risk level, key migration nodes, background control, minimum spacing, and on-site constraints.

2. The method according to claim 1, characterized in that, The core area and influence area mentioned in step S1 include a core area, a near-field influence area, and a far-field influence area. The core area is determined by the boundary of the industrial park or the area where industrial activities are concentrated. The near-field influence area is formed by the outward expansion of the core area and is used to characterize the area affected by pollutants through surface runoff, ditch transport, dust deposition, or near-distance atmospheric diffusion. The far-field influence area is formed by the continued outward expansion of the near-field influence area or by extending along downstream water systems, prevailing wind directions, and low-lying collection areas, and provides a spatial basis for the layout of background control points and peripheral verification points.

3. The method according to claim 1, characterized in that, Step S2 also includes providing the following basic data fields for each spatial grid: Remote sensing image fields are used to calculate soil / surface pollution risk indicators and water pollution risk indicators; The air pollution field is used to calculate air pollution transport risk indicators; The pollution source field is used to calculate the pollution source intensity index; Topographic and hydrological fields and wind field fields are used to calculate pollution migration pathway indicators; The land use field, background reference field, and sampling constraint field are used to determine the background comparison area, the peripheral verification area, and the sampling accessibility conditions. Historical monitoring or existing field record fields are used to assist in setting background references and sampling constraints; Near-ground auxiliary fields are used to record existing UAV images, historical photos, positioning information, road conditions and safety constraints within or near each grid cell, and are used for remote sensing anomaly interpretation, background reference setting and accessibility determination of candidate sampling points.

4. The method according to claim 1, characterized in that, The generation of the Soil / Surface Pollution Risk Indicator (TSAI) in step S3 includes: Based on multispectral or hyperspectral remote sensing images, land use data, industrial land boundaries, enterprise boundaries, roads, storage yards, riverbanks, ditches and historical disturbance area data, candidate areas that may be exposed to soil or surface sediments are identified. After removing interfering features, the effective soil area for calculating the soil / surface pollution risk indicator is obtained; Spectral features are extracted from multispectral or hyperspectral remote sensing images within the effective soil area; Set the background reference area; Based on the spectral characteristics of the soil candidate area and the statistics of the background reference area, the soil / surface pollution risk indicator (TSAI) for each spatial grid is calculated. The Soil / Surface Pollution Risk Indicator (TSAI) was normalized and assigned to a unified spatial grid.

5. The method according to claim 1, characterized in that, The generation of the Water Pollution Risk Indicator (WAI) in step S3 includes: Based on multispectral or hyperspectral remote sensing images, water system vector data, land use data, DEM confluence paths, sewage channels, ditches, ponds, reservoirs, and existing field data, identify surface water bodies within the core area and the affected area; After water bodies are identified, they are classified into different categories. Extracting water spectral features from multispectral or hyperspectral remote sensing images of water bodies; Set a background water body or reference water section; Incorporate analysis of upstream and downstream differences and anomalies from nearby sewage sources; Based on the spectral characteristics of water bodies, background water bodies, upstream and downstream differences, and proximity of pollution sources, the water pollution risk indicator (WAI) for each spatial grid is calculated. The water pollution risk indicator (WAI) is normalized and then assigned to a unified spatial grid.

6. The method according to claim 1, characterized in that, The generation of the Atmospheric Pollutant Transport Risk (ATR) index in step S3 includes: Calculate remote sensing anomalies of air pollution; Potential sources of air pollution were identified based on enterprise source data, remote sensing imagery, and existing field data. To calculate the influence of wind field, a distance attenuation function is introduced. Construct the impact intensity of air pollution sources; By combining remote sensing anomalies of air pollution, intensity of pollution source impact, wind direction consistency weight, and distance attenuation weight, the air pollution transport risk index (ATR) for each spatial grid is calculated. The ATR values ​​are normalized and then assigned to a unified spatial grid.

7. The method according to claim 1, characterized in that, The generation of the pollution source intensity index (PSI) in step S4 includes: Potential pollution sources are identified and classified based on enterprise source data, remote sensing images, and existing field data; Determine the impact range of the pollution source and calculate distance attenuation; Hydrological orientation correction, atmospheric orientation correction, and remote sensing anomaly proximity correction are introduced; The pollution source intensity index (PSI) for each spatial grid is calculated by integrating the basic pollution source weight, distance attenuation, migration direction correction, and remote sensing anomaly proximity correction. The PSI values ​​are normalized and then assigned to a unified spatial grid.

8. The method according to claim 1, characterized in that, The generation of the pollution migration pathway index (MTI) in step S4 includes: Based on DEM data, slope, aspect, confluence direction, confluence accumulation, confluence path, and low-lying catchment area are calculated, and the hydrological migration channel index (HTI) is calculated in combination with water system vector data. Calculate hydrological connectivity; The Atmospheric Migration Indicator (ATI) is calculated based on wind field data and pollution source locations. By combining HTI and ATI, a comprehensive pollution migration pathway index MTI is constructed. After normalizing the MTI, the values ​​are assigned to a unified spatial grid.

9. The method according to claim 1, characterized in that, The comprehensive screening index CPRI mentioned in step S5 is: in: For spatial grid The comprehensive screening index; For spatial grid Multi-media anomaly comprehensive index; For spatial grid The normalized pollution source intensity index; For spatial grid Normalized pollution migration pathway indicators; For spatial grid The source term supporting sub-index; For spatial grid The migration index; to The weights are non-negative and their sum is 1; , , , and All are located in the interval [0, 1].

10. The method according to claim 1, characterized in that, The generation of categorized priority sampling points in step S6 includes: identifying candidate hotspots of compound pollution based on high-value areas of CPRI, TSAI, WAI, ATR, PSI, and MTI, and generating categorized priority sampling points by combining risk level, migration channel, and data missing information. Preferably, the categorized priority sampling points include source term near-field land anomaly sampling points, source term associated hydrological migration sampling points, source term associated atmospheric transport sampling points, peripheral verification sampling points, and background control sampling points. Specifically, source term near-field land anomaly sampling points are determined based on PSI, TSAI, CPRI, and high-value areas of surface disturbance; source term associated hydrological migration sampling points are determined based on WAI, MTI, water anomalies, ditch connectivity, and hydrological convergence conditions; source term associated atmospheric transport sampling points are determined based on ATR, prevailing wind direction, leeward relationship, and distance attenuation relationship; peripheral verification sampling points are located in the migration channel extension area, leeward influence area, or peripheral anomaly area within the peripheral verification area; and background control sampling points are located in areas far from pollution sources, far from migration channels, low-risk areas, and with surface conditions comparable to the study area.

Images