Agricultural non-point source pollutant into river process tracing early warning method and system

Abstract

The application provides a kind of agricultural non-point source pollutant into river process tracing early warning method and system, it is related to pollution control technical field, the method includes: obtaining farmland runoff catchment area, river section and soil infiltration layer multi-source hydrological monitoring data, constructs and includes rainfall runoff process and pollutant concentration field three-dimensional space-time distribution characteristic data matrix;Tensor decomposition is carried out to three-dimensional space-time distribution characteristic data matrix, and rainfall interference term and runoff coupling term are stripped, to obtain the pure pollutant migration basic tensor;Based on pollutant migration basic tensor, the instantaneous change rate of pollutant concentration of each spatial node in runoff field is calculated point by point, to obtain two-dimensional concentration gradient evolution vector;According to two-dimensional concentration gradient evolution vector, the instantaneous diffusion trajectory of pollutants in complex surface runoff network is tracked, to obtain the tracing path node set.The application realizes the whole chain accurate tracing and early warning from river into river point to farmland plot level initial release source.

Description

Technical Field

[0001] This invention relates to the field of pollution control technology, and in particular to a method and system for tracing and early warning of agricultural non-point source pollutants entering rivers. Background Technology

[0002] The process of agricultural non-point source pollutants entering rivers is driven by the coupling of multiple hydrological processes such as rainfall runoff, surface runoff, and soil infiltration. As pollutants migrate from farmland to rivers along a complex surface runoff network, their concentration distribution is simultaneously affected by both the instantaneous disturbance caused by direct dilution or scouring from rainfall and the spatial migration caused by runoff transport.

[0003] Existing methods for tracing non-point source pollution typically rely on time-series monitoring data of pollutant concentrations at runoff catchment outlet sections or river sections, employing hydrochemical mass balance models or isotope tracing methods to infer the source region of pollutants. However, when rainfall events and agricultural runoff processes highly overlap in time, the disturbance signal of pollutant concentration caused by rainfall itself and the coupling signal of pollutant migration carried by runoff become intertwined in the monitoring data. Existing methods lack decoupling mechanisms for the contributions of different driving factors, leading to misinterpretations of concentration decreases caused by rainfall dilution as a decrease in pollutant release intensity, or misinterpretation of concentration peaks caused by rainfall runoff as concentrated release events in agricultural areas. Consequently, the tracing path deviates from the actual pollutant migration trajectory, ultimately failing to accurately locate the initial release source of pollutants within the agricultural runoff catchment area. Summary of the Invention

[0004] This invention provides a method and system for tracing and early warning of agricultural non-point source pollutants entering rivers, realizing precise tracing and early warning of the entire chain from the point of entry into the river to the initial release source at the farmland plot level.

[0005] To solve the above-mentioned technical problems, the technical solution of the present invention is as follows: Firstly, a method for tracing and early warning of agricultural non-point source pollutants entering rivers, the method comprising: Step 1: Obtain multi-source hydrological monitoring data from farmland runoff collection areas, river cross-sections, and soil infiltration layers; construct a three-dimensional spatiotemporal distribution feature data matrix that includes rainfall generation processes and pollutant concentration fields; perform tensor decomposition on the three-dimensional spatiotemporal distribution feature data matrix to remove rainfall interference terms and runoff coupling terms, and obtain a pure pollutant migration fundamental tensor. Step 2: Based on the pollutant migration fundamental tensor, calculate the instantaneous change rate of pollutant concentration at each spatial node in the runoff field point by point to obtain the two-dimensional concentration gradient evolution vector. Step 3: Based on the two-dimensional concentration gradient evolution vector, track the instantaneous diffusion trajectory of pollutants in the complex surface runoff network to obtain the set of source tracing path nodes; Step 4: Combine the set of source tracing path nodes to calculate the total spatiotemporal cumulative flux of pollutants into the river within the confluence interval, and obtain the continuous curve of river load. Step 5: Extract the peak section of the continuous curve of the inflow load into the river, and combine it with the water flow velocity and riverbed surface roughness of the river cross section to estimate the advection transport distance of pollutants in the river channel. Step 6: Based on the advection transport distance and the set of source tracing path nodes, accurately locate the initial source of pollutant release within the farmland runoff catchment area and obtain the spatial coordinates of the source; Step 7: Deeply correlate the source spatial coordinates with the peak time of the continuous curve of the inflow load to obtain the source tracing and early warning control instructions for non-point source pollution. Step 8: Based on the non-point source pollution tracing and early warning control command, trigger the physical interception facilities in the farmland runoff collection area and output the tracing and early warning results of the process of agricultural non-point source pollutants entering the river.

[0006] Secondly, a source tracing and early warning system for agricultural non-point source pollutants entering rivers includes: The acquisition module is used to acquire multi-source hydrological monitoring data from farmland runoff catchment areas, river cross-sections, and soil permeability layers, and construct a three-dimensional spatiotemporal distribution feature data matrix that includes rainfall generation processes and pollutant concentration fields; tensor decomposition is performed on the three-dimensional spatiotemporal distribution feature data matrix to remove rainfall interference terms and runoff coupling terms, and obtain a pure pollutant migration fundamental tensor. The calculation module is used to calculate the instantaneous rate of change of pollutant concentration at each spatial node in the runoff field based on the pollutant migration tensor, and obtain the two-dimensional concentration gradient evolution vector. The tracking module is used to track the instantaneous diffusion trajectory of pollutants in a complex surface runoff network based on the two-dimensional concentration gradient evolution vector, and obtain the set of source tracing path nodes; The accumulation module is used to combine the set of source tracing path nodes to calculate the total spatiotemporal cumulative flux of pollutants into the river within the confluence interval, and obtain a continuous curve of the inflow load. The extrapolation module is used to extract the peak section of the continuous curve of the inflow load into the river, and combine the water flow velocity and riverbed surface roughness of the river cross section to extrapolate the advection transport distance of pollutants in the river channel. The positioning module is used to accurately locate the initial release source of pollutants in the farmland runoff catchment area based on the advection transport distance and the set of source tracing path nodes, and obtain the spatial coordinates of the source. The correlation module is used to deeply correlate the spatial coordinates of the source with the peak time of the continuous curve of the inflow load to obtain the source tracing, early warning and control instructions for non-point source pollution. The processing module is used to trigger physical interception facilities in the farmland runoff collection area according to the non-point source pollution tracing and early warning control instructions, and output the tracing and early warning results of the process of agricultural non-point source pollutants entering the river.

[0007] The above-described solution of the present invention has at least the following beneficial effects: By acquiring multi-source hydrological monitoring data from farmland runoff catchment areas, river cross-sections, and soil infiltration layers, a three-dimensional spatiotemporal distribution feature data matrix was constructed, encompassing rainfall-runoff generation processes and pollutant concentration fields. Tensor decomposition was performed on this matrix to remove rainfall interference terms and runoff coupling terms, yielding a purified pollutant migration fundamental tensor. Based on this migration fundamental tensor, the instantaneous rate of change of pollutant concentration at each spatial node in the runoff field was calculated point-by-point, and a two-dimensional concentration gradient evolution vector was constructed. This allowed for tracing the instantaneous diffusion trajectory of pollutants within the complex surface runoff network, obtaining a set of source-tracing path nodes. Subsequently, the total spatiotemporal cumulative flux of pollutants into the river within the catchment area was calculated across the entire region, generating a continuous inflow load curve. Combining river flow velocity and riverbed roughness to extrapolate the advection transport distance, the initial pollutant release source within the farmland runoff catchment area was accurately located. Finally, the spatial coordinates of the source were deeply correlated with the peak inflow load time to generate an early warning control command, triggering physical interception facilities to implement targeted interception. This method effectively decouples rainfall disturbance from runoff migration from multi-source monitoring data, eliminates the misleading influence of rainfall disturbance on the source tracing path, and realizes accurate source tracing and early warning of the entire chain from the river inflow point to the initial release source at the farmland plot level. Attached Figure Description

[0008] Figure 1 This is a schematic flowchart of a method for tracing and early warning of agricultural non-point source pollutants entering rivers, provided by an embodiment of the present invention.

[0009] Figure 2 This is a schematic diagram of an agricultural non-point source pollutant river entry tracing and early warning system provided by an embodiment of the present invention. Detailed Implementation

[0010] Exemplary embodiments of the present disclosure will now be described in more detail with reference to the accompanying drawings. While exemplary embodiments of the present disclosure are shown in the drawings, it should be understood that the present disclosure may be implemented in various forms and should not be limited to the embodiments set forth herein. Rather, these embodiments are provided so that this disclosure will be thorough and complete, and will fully convey the scope of the disclosure to those skilled in the art.

[0011] like Figure 1 As shown, an embodiment of the present invention proposes a method for tracing and early warning of agricultural non-point source pollutants entering rivers, the method comprising the following steps: Step 1: Obtain multi-source hydrological monitoring data from farmland runoff collection areas, river cross-sections, and soil infiltration layers; construct a three-dimensional spatiotemporal distribution feature data matrix that includes rainfall generation processes and pollutant concentration fields; perform tensor decomposition on the three-dimensional spatiotemporal distribution feature data matrix to remove rainfall interference terms and runoff coupling terms, and obtain a pure pollutant migration fundamental tensor. Step 2: Based on the pollutant migration fundamental tensor, calculate the instantaneous change rate of pollutant concentration at each spatial node in the runoff field point by point to obtain the two-dimensional concentration gradient evolution vector. Step 3: Based on the two-dimensional concentration gradient evolution vector, track the instantaneous diffusion trajectory of pollutants in the complex surface runoff network to obtain the set of source tracing path nodes; Step 4: Combine the set of source tracing path nodes to calculate the total spatiotemporal cumulative flux of pollutants into the river within the confluence interval, and obtain the continuous curve of river load. Step 5: Extract the peak section of the continuous curve of the inflow load into the river, and combine it with the water flow velocity and riverbed surface roughness of the river cross section to estimate the advection transport distance of pollutants in the river channel. Step 6: Based on the advection transport distance and the set of source tracing path nodes, accurately locate the initial source of pollutant release within the farmland runoff catchment area and obtain the spatial coordinates of the source; Step 7: Deeply correlate the source spatial coordinates with the peak time of the continuous curve of the inflow load to obtain the source tracing and early warning control instructions for non-point source pollution. Step 8: Based on the non-point source pollution tracing and early warning control command, trigger the physical interception facilities in the farmland runoff collection area and output the tracing and early warning results of the process of agricultural non-point source pollutants entering the river.

[0012] In this embodiment of the invention, by stripping away rainfall interference to obtain pure background data on pollutant migration, the interference of complex meteorological environmental noise is effectively eliminated, significantly improving the purity and calculation accuracy of the source tracing data; by calculating the instantaneous rate of concentration change at each spatial node of the runoff field, the dynamic diffusion trend of pollutants in the complex surface runoff network is accurately captured, and the continuous migration physical path and diffusion range of pollutants are realistically restored; by combining hydrodynamic conditions to calculate the spatiotemporal cumulative total flux into the river in the confluence interval, the continuous and accurate calculation of multi-source confluence load is achieved, effectively avoiding the problem of distorted warning thresholds and ensuring the reliability of warning results; finally, by deriving the deep correlation between the advection transport distance of the river channel and the source tracing path, the specific emission source of agricultural non-point source pollution is accurately located, and the physical interception facilities in the agricultural runoff confluence area are triggered in conjunction, realizing the automated execution from pollution source tracing warning to physical control, and improving the timeliness and interception effectiveness of agricultural non-point source pollution into the river.

[0013] In a preferred embodiment of the present invention, step 1 above may include: Step 11: Based on multi-source hydrological monitoring data, the rainfall event coverage period is divided into time step sequences according to equal monitoring intervals, the farmland surface space is divided into spatial node arrays, and the soil permeable layer is divided into layered unit groups along the depth direction. Then, the rainfall intensity value, runoff depth value, and pollutant concentration value of each spatial node at each time step are filled into a three-dimensional data array with time step sequence number, spatial horizontal axis sequence number, and spatial vertical axis sequence number as three dimensions to obtain a three-dimensional spatiotemporal distribution feature data matrix. Specifically, based on multi-parameter water quality monitoring sensors deployed at key nodes of surface runoff in farmland runoff collection areas, integrated velocity and water level monitoring stations deployed at river cross-sections, and soil moisture leaching monitoring probes buried at different depths in the soil permeable layer, the rainfall intensity time sequence, surface runoff depth spatial distribution, soil vertical infiltration rate profile, and total nitrogen and total phosphorus pollutant concentration sampling values ​​at different spatial locations within the same rainfall event coverage period are simultaneously collected at equal monitoring intervals in minutes to obtain the original dataset of multi-source hydrological monitoring.

[0014] Based on the original dataset from multi-source hydrological monitoring, the time periods covered by rainfall events are divided according to the aforementioned equal monitoring intervals. Each time step yields a time step number. The sequence, where The value ranges from 1 to Integers; the farmland surface space is divided into horizontal planes with sides of length . The regular square grid is divided into OK The spatial node array of columns yields the row numbers. sequence and column number The sequence, where The value ranges from 1 to integers, The value ranges from 1 to Integers; the soil permeable layer is divided along the vertical depth direction into layers with a thickness of... The equal-thickness layering method is divided into Each layered unit group yields a layer number. The sequence, where The value ranges from 1 to Integers.

[0015] Based on this, construct a three-dimensional data array. This array is indexed by time step number. The first dimension, based on spatial row numbering As the second dimension, using spatial column numbers For a three-dimensional data array, the third dimension is used. The Middle The time step, the first Line 1 For the spatial grid cell positions corresponding to each column, fill in the data as follows: Construct a rainfall intensity array ,in The Line 1 Liede The element value at column position indicates the position of the first element. Each time step spatial grid unit Rainfall intensity values ​​at the location, in millimeters per hour; construct runoff depth arrays. ,in The Line 1 Liede The element value at column position indicates the position of the first element. Each time step spatial grid unit Surface runoff depth values ​​at the location, in millimeters; construct pollutant concentration arrays. ,in The Line 1 Liede The element value at column position indicates the position of the first element. Each time step spatial grid unit The pollutant concentration at the location is measured in milligrams per liter.

[0016] Therefore, the rainfall intensity array Runoff depth array and pollutant concentration array Stacked along the fourth dimension, they form a three-dimensional spatiotemporal distribution feature data matrix containing three hydrological variables.

[0017] Step 12: Based on the three-dimensional spatiotemporal distribution feature data matrix, perform Tucker decomposition on the three-dimensional data array along the time dimension. This decomposes the array into a rainfall modal factor matrix, a runoff spatial modal factor matrix, a pollutant concentration spatial modal factor matrix, and a core tensor reflecting the interaction between the three modes, to obtain the decomposition result. Specifically, this includes: based on the obtained three-dimensional spatiotemporal distribution feature data matrix, performing Tucker decomposition on the three-dimensional data array constituting the matrix along the first dimension of the matrix, i.e., the time dimension. Specifically, the three-dimensional spatiotemporal distribution feature data matrix is ​​denoted as a third-order tensor. Along the first modal direction, this tensor is expanded into a first-modal expansion matrix, where the rows are time step indices and the columns are spatial grid indices; along the second modal direction, this tensor is expanded into a second-modal expansion matrix, where the rows are spatial row indices and the columns are combined indices of time step and spatial column indices; along the third modal direction, this tensor is expanded into a third-modal expansion matrix, where the rows are spatial column indices and the columns are combined indices of time step and spatial row indices.

[0018] Singular value decomposition was performed on the first, second, and third mode expansion matrices, respectively, and the left singular vectors of the first few principal components of each were extracted to form the rainfall mode factor matrix. Runoff spatial modal factor matrix and pollutant concentration spatial mode factor matrix Among them, the rainfall mode factor matrix The dimension is Multiply List, The cutoff rank of the rainfall modes represents the number of rainfall modes retained; the runoff spatial mode factor matrix. The dimension is Multiply List, The cutoff rank of the runoff spatial mode; the pollutant concentration spatial mode factor matrix. The dimension is Multiply List, It represents the truncated rank of the spatial mode of pollutant concentration.

[0019] The core tensor is calculated using the following formula based on the three modality factor matrices. : ; in: The obtained three-dimensional spatiotemporal distribution feature data matrix is ​​a third-order tensor. Represents the rainfall modal factor matrix; Represents the runoff spatial modal factor matrix; Represents the spatial modal factor matrix of pollutant concentration; superscript Indicates matrix transpose; symbol Represents tensor-matrix modular multiplication along the first mode direction; symbol Represents tensor-matrix modular multiplication along the second mode direction; symbol This represents the tensor-matrix modular multiplication operation along the direction of the third mode; Represents the core tensor, with dimension . Multiply column multiplied by The layer reflects the intensity of the interaction among the three spatial modes of rainfall, runoff, and pollutant concentration.

[0020] Thus, the rainfall mode factor matrix is ​​obtained. Runoff spatial modal factor matrix Spatial mode factor matrix of pollutant concentration and core tensor The tensor decomposition results.

[0021] Step 13: Based on the decomposition results, calculate the Pearson correlation coefficient between each column of the rainfall modal factor matrix and the rainfall intensity time series. Load vectors with correlation coefficients exceeding a preset interference threshold are marked as rainfall interference terms and removed. Accordingly, delete the corresponding sub-tensor slices from the core tensor. Then, resynthesize the pruned factor matrices and the core tensor along their respective modal directions to obtain a pure pollutant migration fundamental tensor stripped of direct rainfall interference. Specifically, this includes: based on the obtained tensor decomposition results, processing the rainfall modal factor matrix... For each column of load vectors, calculate the Pearson correlation coefficient between it and the obtained rainfall intensity time series.

[0022] Specifically, for the rainfall mode factor matrix The Column load vector, put the column load vector in The values ​​at each time step are denoted as follows: The spatial mean of rainfall intensity at the corresponding time step in the time series of rainfall intensity is denoted as follows: Calculate the number using the following formula. Pearson correlation coefficient between column load vector and rainfall intensity time series : ; in: Represents the rainfall mode factor matrix Column number, The value ranges from 1 to Integers; Indicates the first The time step to the first The element values ​​of the column load vector are dimensionless; Indicates the first Column load vectors in all The mean at each time step, dimensionless; Indicates the first The average rainfall intensity across the entire spatial grid at each time step, expressed in millimeters per hour; Indicating the time series of rainfall intensity in the entire The average value at each time step, in millimeters per hour; Indicates the first The Pearson correlation coefficient between the column load vector and the time series of rainfall intensity is dimensionless and ranges from -1 to 1. Indicates the total number of time steps; Symbols indicate from equal to 1 to The summation operation.

[0023] The calculated Pearson correlation coefficient Each with a preset interference threshold Compare the results. Preset interference threshold. The value is 0.75, which means that when the absolute value of the Pearson correlation coefficient between a certain load vector and the time series of rainfall intensity exceeds 0.75, it is considered that the mode represented by the load vector mainly reflects the direct interference effect of rainfall rather than the independent migration signal of pollutants carried by runoff.

[0024] The absolute value will be greater than The column number corresponding to the load vector of the condition Marked as rainfall disturbance terms. From the rainfall modal factor matrix. Remove the marked columns to obtain the trimmed rainfall mode factor matrix. Its dimensions are Multiply Column, among which This indicates the number of non-rainfall modes retained after removing the rainfall interference column. equal Subtract the number of marked columns.

[0025] Based on the column number being marked Along the core tensor Delete the sub-tensor slice corresponding to the column index in the first modal direction to obtain the pruned core tensor. Its dimensions are Multiply column multiplied by layer.

[0026] The pruned core tensor With the trimmed rainfall modal factor matrix The original runoff spatial mode factor matrix and the original pollutant concentration spatial mode factor matrix Reconstruct the tensor using the following formula: ; in: This represents the reconstructed pure three-dimensional tensor, that is, the basic tensor of pollutant migration stripped of direct interference from rainfall; This represents the pruned core tensor; This represents the trimmed rainfall modal factor matrix; This represents the original runoff spatial modal factor matrix; This represents the original spatial mode factor matrix of pollutant concentrations; This represents a modular multiplication operation along the direction of the first mode; This represents the modular multiplication operation along the direction of the second mode; This represents a modular multiplication operation along the direction of the third mode.

[0027] Thus, the pure pollutant migration fundamental tensor is obtained. In this tensor, the instantaneous dilution and scouring interference of rainfall on pollutant concentrations have been removed, and only the spatiotemporal distribution information of pollutants carried independently by runoff is retained.

[0028] In this embodiment of the invention, a three-dimensional spatiotemporal distribution feature data matrix containing the rainfall-runoff process and pollutant concentration field is constructed. Tensor decomposition is then performed along the time dimension to accurately separate rainfall interference terms and runoff coupling terms, effectively eliminating the interference of complex meteorological environmental noise on the source tracing data. By calculating the correlation coefficient between rainfall modal factors and the time series of rainfall intensity, direct interference components from rainfall are accurately identified and eliminated, and a pure pollutant migration fundamental tensor is resynthesized. This process significantly improves the purity and computational accuracy of the source tracing fundamental data, realistically restoring the background state of pollutant migration on farmland surfaces and eliminating the impact of data abrupt changes caused by short-term heavy rainfall.

[0029] In a preferred embodiment of the present invention, step 2 above may include: Step 21: Based on the purified pollutant migration fundamental tensor, extract spatial distribution slices of pollutant concentration from two adjacent time points along the time dimension. For each spatial node, calculate the concentration difference between the previous and subsequent time points and divide it by the time interval to obtain the pollutant concentration temporal change rate. Specifically, this includes: using the purified pollutant migration fundamental tensor obtained in Step 13... Along the time dimension, spatial distribution slices of pollutant concentrations that are adjacent in time are extracted sequentially.

[0030] Specifically, for the first The time step and the first A time step, from Extract the first from each of the following: Spatial distribution slices of pollutant concentrations at each time step and the Spatial distribution slices of pollutant concentrations at each time step Each slice is a OK A matrix of columns, the first column in the matrix Line 1 The element values ​​of a column represent the corresponding spatial grid cells. The pollutant concentration value at the corresponding time, in milligrams per liter.

[0031] For each spatial grid cell in the spatial node array The node at the 1st epoch is calculated using the following formula. The and the first Rate of change of pollutant concentration over time within a time window consisting of several time steps: ; in: Representing spatial nodes In the time window The rate of change of pollutant concentration over time, expressed in milligrams per liter per hour; Slice The Middle Line 1 The pollutant concentration values ​​at each location are listed, in milligrams per liter. Slice The Middle Line 1 The pollutant concentration values ​​at each location are listed, in milligrams per liter. This indicates the monitoring time interval between two adjacent time steps, i.e., the equal monitoring interval, in hours.

[0032] Traverse all Each time window and all The pollutant concentration time change rate field is obtained by calculating the pollutant concentration time change rate field for each node in each time window using a spatial grid cell.

[0033] Step 22: Based on the pollutant concentration change rate over time, extract adjacent spatial node pairs along the main direction of surface runoff and its orthogonal directions within the same time window. Calculate the ratio of the concentration difference to the spatial distance for each pair of nodes to obtain the spatial gradient of concentration in the horizontal and vertical axes. Specifically, this includes: based on the obtained farmland surface spatial digital elevation model data, for each spatial grid cell... Calculate the direction of the surface slope, and determine the direction of the steepest slope as the main direction of surface runoff at that spatial grid cell. Let the spatial grid cell be... The angle between the main direction of surface runoff at a given location and the positive direction of the horizontal axis is... Then the unit direction vector in the horizontal direction is horizontal to the right, and the unit direction vector in the vertical direction is vertical upward.

[0034] For the horizontal axis direction, i.e., the main direction of runoff, within the time window Inner space node neighboring nodes The spatial gradient of concentration along the horizontal axis is calculated using the following formula: ; in: Representing spatial nodes Concentration spatial gradient along the horizontal axis, in milligrams per liter per meter; Representing spatial nodes The average pollutant concentration at two time steps within this time window, expressed in milligrams per liter, is equal to... ; Representing spatial nodes The average pollutant concentration at two time steps within this time window, expressed in milligrams per liter, is equal to... ; This indicates the center-to-center distance between adjacent spatial grid cells, i.e., the grid side length, in meters.

[0035] For the vertical axis, which is orthogonal to the main runoff direction, within the time window Inner space node neighboring nodes The spatial gradient of concentration along the vertical axis is calculated using the following formula: ; in: Representing spatial nodes Concentration spatial gradient along the vertical axis, in milligrams per liter per meter; Representing spatial nodes The average pollutant concentration at two time steps within this time window, expressed in milligrams per liter; the meanings of other symbols are the same as above.

[0036] Traverse all Each time window and all For nodes that can form effective adjacent node pairs in a spatial grid cell, calculate the spatial concentration gradient of each node in the horizontal and vertical directions within each time window.

[0037] Step 23: Based on the concentration-time change rate and the concentration spatial gradients in the horizontal and vertical directions, construct a three-dimensional vector, where the first component is the concentration-time change rate, the second component is the spatial gradient in the horizontal direction, and the third component is the spatial gradient in the vertical direction. Assign values ​​to each time window and spatial node to obtain a two-dimensional concentration gradient evolution vector. Specifically, this includes: based on the pollutant concentration-time change rate of each node within each time window... And the obtained concentration spatial gradient along the horizontal axis. and the spatial gradient of concentration along the vertical axis The three elements are combined into a three-dimensional vector.

[0038] Specifically, for spatial nodes In the time window The two-dimensional concentration gradient evolution vector within is assigned values ​​as follows: ; in: Representing spatial nodes In the time window The two-dimensional concentration gradient evolution vector within the vector; the first component is equal to the obtained pollutant concentration time change rate, in milligrams per liter per hour, reflecting the instantaneous trend of pollutant concentration increasing or decreasing over time; the second component is equal to the obtained spatial gradient of concentration in the horizontal axis, in milligrams per liter per meter, reflecting the spatial distribution trend of pollutant concentration with spatial location along the main direction of runoff; the third component is equal to the obtained spatial gradient of concentration in the vertical axis, in milligrams per liter per meter, reflecting the spatial distribution trend of pollutant concentration with spatial location along the orthogonal direction of runoff.

[0039] Traverse all Each time window and all Each spatial grid cell is assigned a value, resulting in a vector field composed of all two-dimensional concentration gradient evolution vectors, which is the two-dimensional concentration gradient evolution vector field.

[0040] In this embodiment of the invention, based on the pure pollutant migration fundamental tensor, the rate of change of time is obtained by calculating the ratio of the concentration difference between adjacent time points to the time interval, and the spatial gradient is calculated along the main direction of surface runoff and its orthogonal directions to construct a two-dimensional concentration gradient evolution vector.

[0041] In a preferred embodiment of the present invention, step 3 above may include: Step 31: Based on the two-dimensional concentration gradient evolution vector, take the spatial node where the river cross-section enters the river as the tracking starting point, take the opposite direction of the spatial gradient component at the tracking starting point as the instantaneous upward direction, and advance along this direction with a preset tracking step size to obtain the next tracking candidate spatial position. Specifically, this includes: based on the obtained two-dimensional concentration gradient evolution vector field, take the spatial node where the river cross-section enters the river as the tracking starting point, and start the pollutant upward tracking.

[0042] Specifically, let the row number of the river inlet point in the spatial node array be . Column number is Record this node as the starting point of the tracking. Retrieve the starting point of the tracking The two-dimensional concentration gradient evolution vector within the time window corresponding to the first tracking time step. Extract the second component of the vector. and the third component That is, the spatial gradient components in the horizontal and vertical directions.

[0043] Reversing the directions of these two spatial gradient components creates the instantaneous upward direction vector. The x-coordinate component of this vector is... The vertical axis components are The reason for choosing this direction is that after pollutants are released from the upstream source, they migrate downstream along the runoff and their concentration gradually decreases. Therefore, the concentration spatial gradient vector points to the direction of increasing concentration, that is, the direction of the upstream source. Using this gradient direction as the instantaneous upstream direction, and advancing along this direction with a preset tracking step size, the source of pollution can be traced back.

[0044] The instantaneous upward direction vector is normalized to a unit direction vector, and the tracking step is performed along this direction with a preset step size. Proceeding forward, the spatial coordinates of the next candidate tracking location are obtained, where the preset tracking step size is... The value is the side length of the spatial grid cell. In this embodiment, the positive integer multiples are used. Values .

[0045] Step 32: Based on the next tracking candidate spatial location, retrieve the two-dimensional concentration gradient evolution vector on the spatial node to which the location belongs. Perform a vector weighted average of the upward direction of the current tracking step and the opposite direction of the spatial gradient component of the node to obtain the corrected upward direction. Based on this, continue to advance to the next tracking point along the corrected direction with a preset tracking step size to obtain the advancement result. Specifically, this includes: based on the obtained coordinates of the next tracking candidate spatial location, retrieve the corresponding two-dimensional concentration gradient evolution vector on the spatial grid node to which the candidate location belongs.

[0046] Specifically, let the continuous coordinate values ​​of the next tracking candidate spatial position obtained after the advancement be... In the spatial node array, determine the spatial grid cell containing the continuous coordinate value, and let the row number corresponding to the grid cell be . Column number is Then the spatial node corresponding to that grid cell is the next tracking node reached in the current tracking step. .

[0047] Retrieve Node The two-dimensional concentration gradient evolution vector within the time window corresponding to the current tracking step. Extract its second component and the third component This forms the spatial gradient vector at the node, and its opposite direction is used to obtain the local up-tracking direction vector at the node.

[0048] The tracking direction vector used in the current tracking step is compared with that node. The local uplink direction vectors at each point are averaged with the same weight to obtain the corrected uplink direction vector. The weighting coefficients of both vectors are both set to 0.5, which means that the tracking inertia of the previous tracking step is retained, while the new direction information provided by the local spatial gradient change of pollutant concentration at the current node is also incorporated.

[0049] Calculate the corrected upward direction vector using the following formula. : ; in: This represents the tracking direction vector used in the current tracking step; Indicates the next tracking node The local up-tracing direction vector at that location; and These represent the magnitudes of the two vectors, used for normalization; 0.5 represents the weighting coefficients of the two direction vectors. This represents the corrected upward direction vector, which has been normalized again to a unit vector. For example, a preset positive minimum value. This is used to avoid division by zero errors. and If the tracking stops, the corrected upward direction vector will be used. After normalization, continue with the preset tracking step size. Proceed along the corrected direction to reach the next tracking point.

[0050] Step 33: Based on the tracking results, repeatedly perform direction retrieval, vector weighted averaging, and step-size advancement operations until the upstream boundary of the confluence region is reached, thus obtaining the upstream path formed by arranging the spatial nodes traversed in all tracking steps in chronological order. Specifically, this includes: based on the tracking results, repeatedly perform direction retrieval, vector weighted averaging, and step-size advancement operations until the upstream boundary of the confluence region is reached, thus obtaining the upstream path formed by arranging the spatial nodes traversed in all tracking steps in chronological order. After each tracking step is completed, the spatial node reached in the current tracking step is denoted as... And record the row and column numbers of that node. and the time window number corresponding to the tracking step At the beginning of the first During each tracking step, with nodes To begin anew, repeat step 32.

[0051] After each tracking step is completed, it is determined whether the tracking termination condition is met. The tracking termination condition consists of two parts: the first termination condition is that the current tracking node has reached the upstream boundary of the confluence region, i.e., the row and column index of the current tracking node. Located at the upstream boundary row of the pre-defined farmland runoff collection area in the spatial node array; the second termination condition is the currently tracked node. The pollutant concentration at the location is lower than the preset background concentration threshold. Among them, the background concentration threshold The value is set at 1.2 times the average pollutant concentration at the upstream reference node of the farmland runoff catchment area during the monitoring period. This indicates that the pollutant concentration has dropped to near the natural background level, and further upstream tracing is no longer meaningful. If any of the above termination conditions are met, the tracking cycle terminates; otherwise, the next tracking step continues. After tracking is terminated, all [resources / resources] will be [accounted for / etc.]. The spatial nodes traversed by each tracking step Arrange them in the order of tracing to obtain the upstream path sequence.

[0052] Step 34: Based on the upstream path, extract all surface spatial nodes located within the farmland runoff catchment area along the upstream path, arrange them in the tracing order, and remove duplicate nodes due to path intersections to obtain the source path node set. Specifically, this includes: based on the obtained upstream path sequence, extracting all surface spatial nodes located within the farmland runoff catchment area from the sequence. Specifically, the upstream path sequence... The row and column numbers of each node in the array The spatial boundaries of the farmland runoff catchment area are compared one by one. If the node row number Falling into the upstream boundary row number of the confluence area Row number to downstream boundary of the confluence area Within the closed interval, and the column index The sequence number falls from the left boundary of the confluence zone. Column number to the right boundary of the confluence area If a node is within a closed interval, it is determined that the node is located within the farmland runoff collection area and is retained; otherwise, it is discarded.

[0053] Arrange all retained nodes into an ordered set according to the tracing order. Remove duplicate nodes from this ordered set due to path intersections: traverse backwards from the first node, and if the current node's row and column index... If a node has the same row and column number as a node that has already appeared in the set, only the position of its first appearance is retained, and the nodes that appear repeatedly afterward are removed from the set.

[0054] After deduplication, the set of source path nodes is obtained, denoted as set. The nodes in this set are arranged in spatial order from upstream to downstream.

[0055] In this embodiment of the invention, the river inlet point is taken as the starting point for tracing, and the process is progressively deduced in the opposite direction of the spatial gradient. The upstream direction is continuously corrected through vector weighted averaging until the upstream boundary of the runoff area is reached. This accurately extracts the surface spatial nodes within the farmland runoff collection area, resulting in a set of source tracing path nodes. This process enables accurate reverse tracing of pollutants from the river inlet to the upstream of farmland, effectively overcoming the problems of trajectory deviation and divergence caused by the variable flow direction and topographic undulations in complex surface runoff networks. By continuously tracking and eliminating overlapping nodes, the continuous migration physical path of pollutants is realistically restored, and the actual flow area of ​​pollutants is locked.

[0056] In a preferred embodiment of the present invention, step 4 above may include: Step 41: Based on the source tracing path node set, extract the pollutant concentration value of each node at each monitoring time, multiply it by the runoff depth and grid cell area at that node to obtain the instantaneous pollutant mass stock at each monitoring time. Specifically, this includes: based on the obtained source tracing path node set... Extract the pollutant concentration value of each node in the set at each monitoring time, and calculate the instantaneous pollutant mass inventory of each node. Specifically, for the set of source tracing path nodes... The first in Let there be nodes, and their row and column numbers be... Extract the node at the 1st Pollutant concentration values ​​at each time step The unit is milligrams per liter, from the obtained runoff depth array Extract the runoff depth value of the same node at the same time step. The unit is millimeters. The node at the [node number] is calculated using the following formula. Instantaneous pollutant mass inventory at each time step: ; in: Represents the set of trace path nodes The node number in the middle, The value ranges from 1 to integers, For set The total number of nodes in the middle; Indicates the first The source tracing path node is at the... The instantaneous pollutant mass stock at each time step, in grams; This indicates the pollutant concentration value at that node at that time step, in milligrams per liter; This represents the runoff depth value at that node at that time step, in millimeters; This represents the surface area of ​​a single spatial grid cell, expressed in square meters, and its value is equal to... The square of, is the grid side length; dividing by 1000 is the dimension conversion factor for converting runoff depth from millimeters to meters and pollutant mass from milligrams to grams.

[0057] Traverse the set of source path nodes All of them Each node and all Each time step calculates the instantaneous pollutant mass inventory at each node at each monitoring moment.

[0058] Step 42: Based on the instantaneous pollutant mass stock at each node, calculate the mass difference between upstream and downstream adjacent node pairs along the source tracing path from upstream to downstream. Record this difference as the comprehensive loss mass between adjacent node pairs caused by soil infiltration, microbial degradation, and surface volatilization. Specifically, this includes: based on the obtained instantaneous pollutant mass stock at each node at each monitoring time... The quality difference between upstream and downstream nodes is calculated for each adjacent node along the tracing path from upstream to downstream. Specifically, for the set of nodes along the tracing path... The adjacent number The node and the first Each node, at the same monitoring time The difference in pollutant quality between two adjacent nodes is calculated using the following formula: ; in, Indicates the first upstream node at each time step With downstream nodes The instantaneous difference in pollutant mass between the two, expressed in grams; Indicates upstream node The instantaneous mass of pollutants, expressed in grams; Indicates downstream node The instantaneous mass of pollutants, expressed in grams.

[0059] like A positive value indicates that the mass of pollutants lost along the runoff section between adjacent node pairs due to the combined effects of three natural attenuation mechanisms: downward seepage through soil, degradation and transformation by soil microorganisms, and volatilization from the surface to the atmosphere. This positive value represents the total mass loss between the adjacent node pairs. If... A negative value indicates that there may be other inflowing tributaries carrying pollutants between the adjacent node pairs. In this case, the absolute value of the negative value is recorded as the inflow mass, and the total loss mass between the adjacent node pairs is set to zero to conservatively estimate the loss value.

[0060] The calculated combined loss mass between each pair of adjacent nodes and at each time step is recorded as a loss array. , the elements Indicates the first in the tracing path The node and the first Between the nodes, at the node The total loss mass at each time step, expressed in grams.

[0061] Step 43: Based on the comprehensive loss mass, the instantaneous pollutant mass stock at the downstream boundary node adjacent to the river section in the source tracing path is deducted segment by segment from the comprehensive loss mass between each adjacent node pair along the upstream direction to obtain the actual net inflow flux into the river at that boundary node after deducting the loss along the way. Specifically, this includes: based on the obtained loss array... Set of trace path nodes The instantaneous pollutant mass stock at the most downstream node is deducted segment by segment upstream from the comprehensive loss mass between adjacent node pairs to obtain the actual net inflow flux into the river channel at that boundary node after deducting the loss along the way. Specifically, let the set of source-tracing path nodes be defined. The serial number of the downstream border node is In its first The instantaneous pollutant mass stock at each time step is Starting from this border node, the adjacent node pairs passed sequentially upstream along the tracing path are: , ... Each pair of adjacent nodes in the th The comprehensive loss mass at each time step is as follows: , ... .

[0062] The border node is calculated using the following formula at the [number]th [year]. Net inflow flux at each time step: ; in: Indicates the first The actual net inflow flux into the river channel at each time step, expressed in grams per time step. Indicates the border node at the th The instantaneous pollutant mass stock at each time step, in grams; Indicates the first The nth adjacent node pair in the th... The total loss mass at each time step, in grams; Symbols indicate from equal to 1 to Summation of .

[0063] Traverse all For each time step, calculate the net inflow flux at each time step. The time series sequence of net inflow flux was obtained.

[0064] Step 44: Based on the net inflow flux at each monitoring time, with the monitoring time as the horizontal axis and the net inflow flux as the vertical axis, connect all data points in chronological order to form a continuous time-varying curve, thus obtaining a continuous inflow load curve. Specifically, based on the obtained net inflow flux time series, construct a continuous inflow load curve with each monitoring time as the horizontal axis variable and the net inflow flux as the vertical axis variable. Specifically, use the actual time point corresponding to each time step as a discrete coordinate point on the horizontal axis, and use the net inflow flux value of the corresponding time step as a discrete coordinate point on the vertical axis. Use a piecewise linear interpolation method to connect all discrete data points in chronological order to form a continuous time-varying curve.

[0065] For two adjacent time steps and at any time in between Calculate using the following linear interpolation formula Net inflow flux at time: ; in: Indicates in Interpolated net inflow flux at time step, in grams per time step; and They represent the first time. The and the first Net inflow flux at each time step, expressed in grams per time step; and They represent the first The and the first Each time step corresponds to the actual monitoring time; This represents any time interval at which interpolation is to be performed, with values ​​ranging from 1 to 2. and between.

[0066] By merging all interpolated points with the original discrete points, a continuous curve of inflow load in time is obtained, denoted as curve . Its horizontal axis represents the time variable. The vertical axis represents the net flux into the river.

[0067] In this embodiment of the invention, the instantaneous pollutant mass stock at each node is calculated, and the comprehensive loss mass caused by soil infiltration, microbial degradation, and surface volatilization is calculated for each adjacent node. The net inflow flux is obtained by subtracting the loss segment by segment along the upstream direction. Finally, a continuous curve of inflow load is plotted. This process realizes continuous and accurate measurement of multi-source confluence load, fully considers the physical and biochemical loss factors of pollutants during migration, effectively avoids the problem of distortion of the warning threshold caused by direct accumulation, and intuitively reflects the dynamic evolution law of inflow pollution load by constructing a continuous time variation curve, thus ensuring the objectivity of the inflow pollution load assessment results.

[0068] In a preferred embodiment of the present invention, step 5 above may include: Step 51: Based on the continuous curve of the inflow load, scan the load value time by time from the starting point to the ending point, extract the continuous time periods where the load value exceeds the preset peak value judgment threshold as peak segments, and obtain the set of peak segments after traversing the entire curve. Specifically, this includes: based on the obtained continuous curve of the inflow load... The time corresponding to the starting point of the curve The time corresponding to the end of the curve The load values ​​are scanned periodically, specifically by moving gradually along the entire time axis in a sliding window manner, with the time width of each analysis window equal to the monitoring interval. Three times The sliding step size is equal to Within each analysis window, determine whether the load value of the continuous curve of the inflow load at each time point in that window meets the preset peak value determination condition.

[0069] The preset peak value determination criterion is: the net inflow flux value at all time points within this analysis window. The average value exceeds the preset peak value threshold. Furthermore, the net inflow flux at the center time point of the analysis window is the maximum value within that window, where a preset peak value determination threshold is defined. Determine using the following formula: ; in: This represents the average net inflow flux over the entire time period of the continuous curve of inflow load, expressed in grams per time step. The standard deviation of net inflow flux over the entire time period is expressed as grams per time step; the coefficient 1.5 is an empirical coefficient used to set the threshold at 1.5 standard deviations above the mean to achieve a reasonable distinction between peak and normal loads.

[0070] The time period covered by the analysis window that meets the above peak determination criteria is determined as the candidate interval for the peak segment. If two adjacent candidate intervals overlap on the time axis or the interval is less than [a certain value], [further details are needed]. Then, these two candidate intervals are merged into a single continuous peak segment. All the merged peak segments are then arranged in chronological order to form a set of peak segments, denoted as set . .

[0071] Step 52: Based on the set of peak segments, retrieve the monitoring time with the largest load value for each peak segment and mark that time as the peak inflow time corresponding to that segment. Specifically, this includes: based on the obtained set of peak segments... For each peak segment in the set, the monitoring time with the largest load value is retrieved within its covered time range. Specifically, for the peak segment set... The first in Each peak segment has a time range from time 1 to 12. At the time Within this time frame, the continuous curve of inflow load is traversed. For all discrete time points within this interval, compare the net inflow flux values ​​at each time point, and mark the time point with the largest net inflow flux value as the peak inflow time corresponding to that peak segment. Traverse the set of peak segments We obtain a one-to-one correspondence between each peak segment and the time of entry into the river for each peak, denoted as the mapping set. ,in This indicates the total number of peak segments.

[0072] Step 53: Based on the entry time of each peak flow, retrieve the measured flow velocity values ​​of the river cross-section at that time, and find the corresponding riverbed surface roughness parameter values ​​from the preset Manning roughness coefficient reference table according to the riverbed material type. Specifically, this includes: based on the obtained entry times of each peak flow... The measured flow velocity at the cross-section of the river channel at that moment is retrieved, and the corresponding riverbed surface roughness parameter value is looked up from the preset Manning roughness coefficient reference table according to the riverbed material type. Specifically, for the p-th peak inflow moment... From the historical monitoring database of integrated flow velocity and water level monitoring stations deployed along the river cross-section, the closest measured flow velocity record to that moment is retrieved, and the cross-sectional average flow velocity value is extracted from that record. The unit is meters per second.

[0073] Based on the location of the river section where the cross-section is situated, the riverbed material type for that section is determined. The riverbed material type is selected from a pre-defined riverbed material classification system, which includes five types: natural fine sand riverbed, coarse sand and gravel riverbed, pebble riverbed, clay riverbed, and artificially lined riverbed. The Manning roughness coefficient corresponding to each type is determined according to the following mapping relationship: Manning roughness coefficient corresponding to natural fine sand riverbed. The coefficient is 0.020, indicating that this type of riverbed is composed of uniform fine sand, resulting in low flow resistance; the Manning roughness coefficient corresponds to a coarse sand and gravel riverbed. The coefficient is 0.028, indicating that this type of riverbed consists of a mixture of coarse sand and gravel, with moderate flow resistance; the Manning roughness coefficient corresponds to the pebble riverbed. The coefficient is 0.035. This type of riverbed surface is composed of large-diameter pebbles, resulting in high flow resistance; the Manning roughness coefficient corresponding to clay riverbeds is... The coefficient of roughness n is 0.018, indicating that the surface of this type of riverbed is composed of smooth clay, resulting in minimal water flow resistance. The coefficient of roughness n corresponding to the artificially lined riverbed is 0.015, indicating that the surface of this type of riverbed is lined with concrete or masonry, resulting in extremely low water flow resistance.

[0074] Based on the actual riverbed material type of the river section, the corresponding Manning roughness coefficient value can be found from the above mapping relationship. This serves as the surface roughness parameter value of the riverbed at that location.

[0075] Step 54: Based on the measured values ​​of water flow velocity and riverbed surface roughness parameters, the water flow velocity is reduced according to the riverbed surface roughness to obtain the effective transport rate. This effective transport rate is then multiplied by the time span from the peak inflow time to the warning target time to obtain the advection transport distance of the pollutant plume in the river channel corresponding to each peak segment. Specifically, this includes: retrieving the measured water flow velocity values ​​at the river cross-section at each peak inflow time. Based on the riverbed material type, the corresponding riverbed surface roughness parameter value is found from the preset Manning roughness coefficient reference table. Based on the measured water flow velocity Riverbed surface roughness parameter values Calculate the effective transport rate using the following formula. : ; in, Based on Manning roughness coefficient The velocity reduction factor has a value ranging from 0 to 1, and the specific value is determined through a pre-calibration experiment.

[0076] Based on this effective transport rate Calculate the time of entry into the river from the peak using the following formula. From the start of the warning to the target of attention The advection transport distance of the pollutant plume within the river channel over the specified time span: ; in: Indicates the first The advection transport distance of pollutant plumes within the river channel corresponding to each peak segment, in meters; Indicates the effective transport rate, measured in meters per second; The time of the warning focus is the time after the preset warning lead time is calculated from the current time. Indicates the first The peak moment of river entry; This indicates the time span from the peak of the pollutant plume entering the river to the moment of the warning target, expressed in seconds.

[0077] Traverse all For each peak segment, the corresponding advection transport distance is calculated. .

[0078] In this embodiment of the invention, by extracting the peak segment of the continuous curve of the inflow load and the corresponding peak inflow time, the measured value of the water flow velocity at the river cross-section is retrieved, and the surface roughness parameter value of the riverbed is found according to the riverbed material type. The water flow velocity is then reduced according to the roughness to obtain the effective transport rate, and the advection transport distance of the pollutant plume in the river channel is deduced. This process deeply integrates the river hydrodynamic conditions and riverbed physical characteristics, accurately quantifies the actual transport distance of the pollutant plume in the river channel, and effectively eliminates the transport error caused by water flow turbulence and riverbed resistance. By introducing a roughness coefficient to scientifically reduce the flow velocity, the physical mechanism of advection transport of pollutants in complex river cross-sections is truly reflected.

[0079] In a preferred embodiment of the present invention, step 6 above may include: Step 61: Based on the advection transport distance, measure the length of the river segment equidistant from the river cross-section inlet point upstream along the river centerline to determine the current river segment location coordinates of the pollutant plume after advection transport. Specifically, this includes: based on the obtained advection transport distance... Taking the river inlet point of the cross-section as the starting point for measurement, the length of the river segment equidistant from the centerline of the river upstream, i.e., in the opposite direction to the water flow direction, is measured. Specifically, the measurement is performed based on the coordinates of the determined river inlet point. The corresponding latitude and longitude coordinates or projected coordinates in real geographic space Starting from the centerline of the river, the lengths of each segment along the digital vector path are accumulated upstream until the accumulated path length first exceeds the advection transport distance. Stop when the accumulation path is located at a point where the accumulation length is exactly equal to the sum of the sums. Geographic coordinates of the location That is, the first The coordinates of the current river section where the pollutant cluster corresponding to each peak section is located after advection transport.

[0080] Traverse all For each peak segment, determine the coordinates of the current river segment's location. .

[0081] Step 62: Based on the current river segment location coordinates, calculate the spatial Euclidean distance between each node in the obtained source tracing path node set and the river segment location coordinates, and mark the node with the smallest distance as the spatial intersection node of the river channel and farmland runoff. Specifically, this includes: based on the determined current river segment location coordinates... In the obtained set of traceability path nodes The spatial Euclidean distance between each source tracing path node and the current river segment's location coordinates is calculated one by one.

[0082] Specifically, for the first Current river section location coordinates and the set of traceability path nodes The first in There are nodes, and let the true geographic coordinates of each node be... The spatial Euclidean distance between the two is calculated as follows: ; in: Indicates the first The current river segment location coordinates and the first The spatial Euclidean distance between each traceability path node is expressed in meters. and They represent the first The x and y coordinates of the current river section location, in meters; and They represent the first The x and y coordinates of each traceability path node, in meters.

[0083] For each Value, in the set of trace path nodes All Searching within each node The node number that takes the minimum value , this node The spatial intersection of river channels and farmland runoff is marked as a node. This spatial intersection node is located at the end of the river's advection transport path and on the source-tracing path of farmland surface runoff, making it a key spatial connection point for pollutants to enter the river from the farmland runoff network.

[0084] Step 63: Based on the spatial convergence nodes, traverse all upstream nodes in reverse order along the source tracing path upstream of the farmland, extracting the pollutant concentration value of each upstream node at the peak inflow time, and calculating the ratio of the concentration value of each upstream node to the total concentration value of all nodes along that upstream path to obtain the concentration contribution ratio of each node. Specifically, this includes: based on the obtained spatial convergence nodes... Node set along the tracing path The path from this node to the upstream direction of the farmland runoff confluence area is traversed in reverse order, visiting all source-tracing path nodes upstream of this spatial intersection node. Specifically, this involves the set of source-tracing path nodes. In the middle, from the spatial intersection node Serial number Begin by extracting nodes sequentially in descending order of their serial numbers. , forming a subset of upstream nodes .

[0085] for Each upstream node in (in The value ranges from 1 to (integer), extract the node at the th position. Peak times when the river enters the river Pollutant concentration values ​​at the corresponding monitoring time steps The unit is milligrams per liter. In all The peak inflow time at each time step The time step number closest to the current moment.

[0086] Calculate the upstream node subset using the following formula. Concentration contribution percentage of each node: ; in: Indicates the first The proportion of the concentration contribution of each upstream node to the pollutant load entering the river is dimensionless and ranges from 0 to 1. The sum of the ratio values ​​of all upstream nodes is equal to 1. Indicates the first The pollutant concentration values ​​at the time step corresponding to the peak inflow time of each upstream node, in milligrams per liter; Represents all nodes in the upstream subset. The sum of pollutant concentrations at each node, expressed in milligrams per liter.

[0087] Step 64: Based on the concentration contribution ratio of each node, the node with the largest ratio is identified as the dominant release node with the most concentrated contribution to the pollutant load entering the river. The horizontal and vertical coordinates of this node in the farmland surface spatial grid are extracted to obtain the source spatial coordinates. Specifically, this includes: based on the obtained concentration contribution ratio of each upstream node... ,exist All Searching within each node The node number that takes the maximum value .

[0088] The node with the highest ratio Determined to be the first The dominant release node, where the pollutant load contribution to the river is most concentrated in each peak section, is the primary source of pollutant release from the farmland surface into runoff. This dominant release node is then extracted. Row number in the farmland surface spatial grid and column number And map it to the horizontal coordinate in real geographic space. and ordinate , obtained the Source spatial coordinates corresponding to each peak segment Traverse all After identifying the peak segments, the source spatial coordinate set is obtained. .

[0089] In this embodiment of the invention, the current river section coordinates of the pollutant plume are determined based on the advection transport distance. Spatial Euclidean distance is calculated from the set of source-tracing path nodes to pinpoint spatial intersection nodes. Then, upstream nodes are traversed in reverse order, and their concentration contribution percentages are calculated. The dominant release node with the largest contribution percentage is identified as the initial release source, thus obtaining the spatial coordinates of the source. This process achieves precise source location from the macroscopic river cross-section to the microscopic farmland plot. By quantifying the concentration contribution percentage of each node, interference from non-dominant pollution sources is eliminated, accurately pinpointing the specific emission source of farmland non-point source pollution.

[0090] In a preferred embodiment of the present invention, step 7 above may include: Step 71: Based on the source spatial coordinates, perform a spatial point query in the pre-constructed geographic information database of farmland runoff catchment areas to determine the farmland plot number and the runoff sub-catchment identifier to which the coordinates fall. Specifically, this includes: based on the source spatial coordinates corresponding to each peak segment... Spatial point queries are performed in a pre-built geographic information database of farmland runoff catchment areas.

[0091] This geographic information database stores the spatial boundary data of all plots within the farmland runoff catchment area in the form of vector areal features. Each plot isal feature includes a plot number attribute field and a runoff sub-catchment identifier attribute field. The spatial point query process is as follows: For the first... Source space coordinates The database is searched sequentially among all land parcel areal features based on spatial inclusion relationships to determine which land parcel isal feature the coordinate point falls within. The land parcel number attribute value of the land parcel isal feature containing the coordinate point is recorded as follows: Record the identifier attribute value of its respective runoff sub-catchment area as .

[0092] Therefore, the first The farmland plot number where the source spatial coordinates fall. and its associated runoff sub-catchment area markers .

[0093] Step 72: Based on the runoff sub-catchment area identifier, retrieve the equipment identification code of the interception facility bound to the sub-catchment area from the physical interception facility ledger to obtain the correspondence between the land parcel and the interception facility. Specifically, this includes: based on the determined runoff sub-catchment area identifier... The equipment identification code of the physical interception facility, which is bound to the runoff sub-catchment area identifier, is retrieved from the pre-deployed and registered physical interception facility ledger.

[0094] This physical interception facility ledger stores basic information on all deployed physical interception facilities within the farmland runoff catchment area in the form of data records. Each record includes fields such as: facility identification code, facility type, facility geographical coordinates, runoff sub-catchment area identifier to which the facility belongs, and the current operational status of the facility. The retrieval process is as follows: To find the query criteria, search the physical interception facility ledger for the value of the runoff sub-catchment area identifier field. For records that are exact matches, extract the Facility / Equipment Identification Code field value from that record, and denote it as... If multiple matching records exist, the facility and equipment identification codes corresponding to all matching records will be extracted to form a list of interception facility and equipment identification codes corresponding to that plot of land. Thus, the first... The plot number corresponding to the source spatial coordinates With the identification code of the interception facility equipment A one-to-one correspondence between them.

[0095] Step 73: Based on the continuous inflow load curve, extract the peak time of the source tracing path corresponding to the spatial coordinates of the source on the curve. Bind the peak time, plot number, and equipment identification code into a triplet to obtain the source tracing-time association record. Specifically, this includes: based on the obtained continuous inflow load curve... Extract the curve that is related to the first... Peak inflow time of the tracing path corresponding to each source spatial coordinate. The peak time when it enters the river , Determined land parcel number and the identified device identification code Perform triple binding to generate a source-time association record in the form of a data record. The structure of this record includes: the peak time field is assigned a value. The source plot number field is assigned the value The interception facility equipment identification code field is assigned the value. traverse all After each peak segment, we obtain the result from... A set of traceability-time-related records, consisting of traceability-time-related records.

[0096] Step 74: Based on the source-time correlation records and combined with the preset early warning lead time parameter, the peak time is shifted forward by this lead time to obtain the interception trigger time. The interception trigger time, plot number, and equipment identification code are combined into control commands. After sorting the interception trigger times, a sequence of control commands for non-point source pollution source tracing and early warning is obtained. Specifically, this includes: based on each record in the obtained source-time correlation record set, combined with the preset early warning lead time parameter... Generate corresponding early warning control instructions for each record.

[0097] Specifically, for the first The source-time-related records have a peak time of [time value]. The interception trigger time is calculated using the following formula: ; in: Indicates the first The interception trigger time corresponding to each record; Indicates the first The peak inflow time recorded in the records; The preset advance warning time parameter is expressed in seconds. It is a preset proportion of the average runoff confluence time from the farthest source point to the river inlet point in the catchment area. In this embodiment, the preset proportion is 0.5, that is, the advance is half of the average confluence time, so as to ensure that the physical interception facilities have sufficient deployment time before the pollutant cloud reaches the interception site.

[0098] Therefore, , and The three parameters are combined into a control instruction structure, denoted as the first one. A warning control instruction, the semantics of which is: in At any time, the device identification code is The physical interception facility sends an activation signal to the plot numbered Pollutants are intercepted from farmland runoff.

[0099] Traversing the source - all records in the time-related record set Each record is generated one by one. Issue warning control commands and intercept trigger times according to each command. The instructions are sorted in order to obtain the sequence of instructions for tracing and controlling non-point source pollution.

[0100] In this embodiment of the invention, the farmland plot number and runoff sub-catchment area identifier are determined by querying the geographic information database based on the source spatial coordinates. The bound physical interception facility identification code is retrieved, and the peak time, plot number, and equipment identification code are bound as a triple. The interception trigger time is obtained by shifting forward with the preset early warning lead time parameter, and a sequence of control instructions for non-point source pollution tracing early warning is generated. This process realizes the automated and precise linkage between pollution tracing early warning and physical control facilities. By introducing an early warning lead time, it is ensured that physical interception facilities are deployed and prepared before the pollution plume arrives. This effectively eliminates the time delay between early warning signal transmission and facility response, improves the timeliness of early warning response and the execution of control instructions, and achieves seamless connection from data tracing to physical interception.

[0101] In a preferred embodiment of the present invention, step 8 above may include: Step 81: According to the non-point source pollution tracing and early warning control command sequence, when the time reaches the interception trigger time specified in the command, a start signal is sent to the drive controller of the physical interception facility indicated by the command, driving the trash rack or interception net to be lowered from the standby position to the runoff channel interception station to complete the interception action. Specifically, this includes: according to the obtained non-point source pollution tracing and early warning control command sequence, when the current system time reaches the interception trigger time specified by the first command in the command sequence... At that time, the early warning control system automatically executes the interception triggering process. Specifically, the early warning control system transmits the device identification code in the instruction through the wireless communication module. The drive controller of the physical interception facility sends a start signal. This drive controller is an embedded control unit installed next to the interception facility. Upon receiving the start signal, the drive motor or hydraulic device lowers the debris barrier or net from its standby position (the height at which the barrier / net is retracted to one side of the runoff channel without obstructing normal runoff) to the interception station in the runoff channel. This interception station is positioned where the bottom edge of the barrier / net is lowered to the bottom of the runoff channel, and the top edge is above a preset safety margin above the historical maximum runoff level, ensuring that the interception facility completely blocks the cross-section of the runoff channel, thus intercepting pollutants carried in the runoff. After the interception action is completed, the drive controller returns an interception completion confirmation signal to the early warning control system. At the moment of interception triggering... The system then continues to monitor the timeline, and when the time reaches the interception trigger moment for the next instruction in the instruction sequence... When the above triggering process m is repeated until all instructions are executed or the rainfall event ends.

[0102] Step 82: Based on the completed interception action, activate the water quality verification sensor deployed downstream of the interception facility. Within a preset verification time window after the interception trigger time, continuously sample at a preset frequency to obtain the measured value sequence of pollutant concentration within the verification window. Compare the measured value sequence with the background concentration value sequence of the same duration before the interception to obtain interception effect comparison data. Specifically, based on the completed interception action, after each interception action is completed, the early warning control system automatically activates the water quality verification sensor deployed at a location adjacent to the downstream of the physical interception facility.

[0103] Specifically, in the first Upon receiving the confirmation signal for the completion of the interception action corresponding to the instruction, the online water quality monitoring sensor located downstream of the interception facility and at a preset distance along the runoff channel is immediately activated. This sensor automatically samples and analyzes the runoff water at a preset verification sampling frequency, outputting the concentration detection value of the target pollutant for each analysis. The preset verification sampling frequency is one-third of the equal monitoring interval described in step 11 to ensure a sufficient density of sampling points within the interception effect verification window. The preset verification time window length is from the interception trigger time... Peak entry point into the river The time extension is then doubled to ensure coverage of the entire process of the pollutant plume before and after passing through the interception facility. Thus, within this validation time window, a sequence of measured pollutant concentrations from multiple sampling points arranged chronologically is obtained, denoted as... The sequence is further analyzed, and simultaneously, a sequence of background pollutant concentrations within a time period equal to the length of the verification time window, prior to the execution of the interception action, is extracted and denoted as [missing information]. sequence, will Sequence and Each sampling point in the sequence that corresponds one-to-one in time is compared point by point, and the concentration difference at each point is calculated to obtain the interception effect comparison data.

[0104] Step 83: Based on the interception effect comparison data, calculate the percentage decrease in the average concentration after interception relative to the average background concentration before interception. If this percentage exceeds the preset effective interception judgment threshold, the interception is deemed effective; otherwise, generate a facility anomaly alarm and push it to the maintenance management terminal. Specifically, this includes: quantifying the effectiveness of the interception operation based on the obtained interception effect comparison data; specifically, calculating the arithmetic mean of the measured pollutant concentration values ​​at all sampling points within the verification window after interception. And the arithmetic mean of the background pollutant concentrations at all sampling points within the background window before interception. Both are in milligrams per liter. The percentage decrease in the average concentration after interception relative to the average background concentration before interception is calculated using the following formula: ; in: Represents the percentage decrease in concentration, dimensionless, and ranges from negative infinity to 1. A positive value indicates that the concentration decreases after interception, while a negative value indicates that the concentration increases after interception. This represents the average background concentration before interception, expressed in milligrams per liter. This represents the average concentration within the verification window after interception, in milligrams per liter; the calculated value... Compared with the preset effective interception threshold Comparison. Preset effective interception threshold. The value is 0.30, which means that when the average concentration of pollutants after interception decreases by more than 30% relative to the average background concentration before interception, the interception operation is considered to have achieved the expected effect of intercepting pollutants.

[0105] like If the interception operation is deemed valid, an interception validity determination conclusion is generated as valid; if If the interception operation fails to achieve the expected results, an invalid interception validity conclusion is generated, and a facility anomaly alarm message is automatically generated. This alarm message includes: the equipment identification code of the physical interception facility corresponding to this interception. Interception trigger time Peak inflow time Concentration decrease ratio The alarm information includes the actual calculated value and the average concentration before and after interception. This alarm message is pushed to the maintenance management terminal via the wireless communication network, prompting maintenance personnel to conduct an on-site inspection of the interception facility to troubleshoot the cause of the malfunction.

[0106] Step 84: Based on the interception effectiveness determination conclusion, summarize and package the source spatial coordinates, peak time, interception trigger time, interception effect comparison data, and determination conclusion to obtain a comprehensive report containing source tracing path information, early warning trigger information, and interception effect information. The output is the source tracing and early warning result, specifically including: summarizing and packaging all key data and results generated and referenced throughout the entire source tracing and early warning process based on the interception effectiveness determination conclusion. Specifically, for the first... The source tracing and early warning process corresponding to each instruction is summarized as follows to form a source tracing and early warning result record: Source spatial coordinates Peak inflow time Interception trigger time The interception effect comparison data includes the background concentration sequence, the verification concentration sequence, and the conclusion on the interception effectiveness.

[0107] Organize the above content into a comprehensive report record by field format, and apply it to all... After each instruction's source tracing and early warning process generates a corresponding comprehensive report record, all of them will be... The report records are summarized in chronological order to generate a complete comprehensive report containing the following three types of information: Source tracing path information: including the source spatial coordinates of each farmland plot and the corresponding source tracing path node name; warning trigger information: including the peak entry time into the river, the corresponding interception trigger time, and the identification code of the triggered interception facility equipment; interception effect information: including the concentration reduction ratio of each interception and the conclusion of the interception effectiveness judgment.

[0108] The comprehensive report will be output as a structured data file to the human-computer interface of the early warning management system for environmental management personnel to review and make decisions. At the same time, the comprehensive report will be archived in the historical source tracing and early warning database as the source tracing and early warning result of agricultural non-point source pollutants entering the river during this rainfall event.

[0109] In this embodiment of the invention, according to a sequence of control commands, the debris barrier or interception net is lowered to the interception station at a specified time, and downstream water quality verification sensors are activated for continuous sampling. The measured concentration sequence is compared point by point with the background concentration sequence to calculate the rate of decrease and determine the effectiveness of the interception. Finally, the data is summarized and packaged to generate a comprehensive report containing information on source tracing, early warning, and interception effectiveness. This process achieves closed-loop management of the entire process from pollution source tracing and early warning to physical interception and effectiveness verification. By verifying the effectiveness of interception through measured data and promptly reporting anomalies, it ensures the actual effectiveness of controlling agricultural non-point source pollution entering rivers, ensures the reliability of physical interception facilities, and provides complete data support and effectiveness evaluation basis for the refined management of the watershed water environment.

[0110] like Figure 2As shown, embodiments of the present invention also provide a source tracing and early warning system for agricultural non-point source pollutants entering rivers, comprising: The acquisition module is used to acquire multi-source hydrological monitoring data from farmland runoff catchment areas, river cross-sections, and soil permeability layers, and construct a three-dimensional spatiotemporal distribution feature data matrix that includes rainfall generation processes and pollutant concentration fields; tensor decomposition is performed on the three-dimensional spatiotemporal distribution feature data matrix to remove rainfall interference terms and runoff coupling terms, and obtain a pure pollutant migration fundamental tensor. The calculation module is used to calculate the instantaneous rate of change of pollutant concentration at each spatial node in the runoff field based on the pollutant migration tensor, and obtain the two-dimensional concentration gradient evolution vector. The tracking module is used to track the instantaneous diffusion trajectory of pollutants in a complex surface runoff network based on the two-dimensional concentration gradient evolution vector, and obtain the set of source tracing path nodes; The accumulation module is used to combine the set of source tracing path nodes to calculate the total spatiotemporal cumulative flux of pollutants into the river within the confluence interval, and obtain a continuous curve of the inflow load. The extrapolation module is used to extract the peak section of the continuous curve of the inflow load into the river, and combine the water flow velocity and riverbed surface roughness of the river cross section to extrapolate the advection transport distance of pollutants in the river channel. The positioning module is used to accurately locate the initial release source of pollutants in the farmland runoff catchment area based on the advection transport distance and the set of source tracing path nodes, and obtain the spatial coordinates of the source. The correlation module is used to deeply correlate the spatial coordinates of the source with the peak time of the continuous curve of the inflow load to obtain the source tracing, early warning and control instructions for non-point source pollution. The processing module is used to trigger physical interception facilities in the farmland runoff collection area according to the non-point source pollution tracing and early warning control instructions, and output the tracing and early warning results of the process of agricultural non-point source pollutants entering the river.

[0111] It should be noted that this system is a system corresponding to the above method. All implementation methods in the above method embodiments are applicable to this embodiment and can achieve the same technical effect.

[0112] The above description represents the preferred embodiments of the present invention. It should be noted that those skilled in the art can make various improvements and modifications without departing from the principles of the present invention, and these improvements and modifications should also be considered within the scope of protection of the present invention.

Claims

1. A method for tracing and early warning of agricultural non-point source pollutants entering rivers, characterized in that, The method includes: Step 1: Obtain multi-source hydrological monitoring data from farmland runoff collection areas, river cross-sections, and soil infiltration layers; construct a three-dimensional spatiotemporal distribution feature data matrix that includes rainfall generation processes and pollutant concentration fields; perform tensor decomposition on the three-dimensional spatiotemporal distribution feature data matrix to remove rainfall interference terms and runoff coupling terms, and obtain a pure pollutant migration fundamental tensor. Step 2: Based on the pollutant migration fundamental tensor, calculate the instantaneous change rate of pollutant concentration at each spatial node in the runoff field point by point to obtain the two-dimensional concentration gradient evolution vector. Step 3: Based on the two-dimensional concentration gradient evolution vector, track the instantaneous diffusion trajectory of pollutants in the complex surface runoff network to obtain the set of source tracing path nodes; Step 4: Combine the set of source tracing path nodes to calculate the total spatiotemporal cumulative flux of pollutants into the river within the confluence interval, and obtain the continuous curve of river load. Step 5: Extract the peak section of the continuous curve of the inflow load into the river, and combine it with the water flow velocity and riverbed surface roughness of the river cross section to estimate the advection transport distance of pollutants in the river channel. Step 6: Based on the advection transport distance and the set of source tracing path nodes, accurately locate the initial source of pollutant release within the farmland runoff catchment area and obtain the spatial coordinates of the source; Step 7: Deeply correlate the source spatial coordinates with the peak time of the continuous curve of the inflow load to obtain the source tracing and early warning control instructions for non-point source pollution. Step 8: Based on the non-point source pollution tracing and early warning control command, trigger the physical interception facilities in the farmland runoff collection area and output the tracing and early warning results of the process of agricultural non-point source pollutants entering the river.

2. The method for tracing and early warning of agricultural non-point source pollutants entering rivers according to claim 1, characterized in that, Step 1: Obtain multi-source hydrological monitoring data from farmland runoff collection areas, river cross-sections, and soil infiltration layers, and construct a three-dimensional spatiotemporal distribution feature data matrix that includes rainfall runoff generation processes and pollutant concentration fields; Tensor decomposition is performed on the three-dimensional spatiotemporal distribution feature data matrix to remove the rainfall disturbance term and runoff coupling term, resulting in a pure pollutant migration fundamental tensor, including: Based on multi-source hydrological monitoring data, the rainfall event coverage period is divided into time step sequences according to equal monitoring intervals, the farmland surface space is divided into spatial node arrays, and the soil permeable layer is divided into layered unit groups along the depth direction. Then, the rainfall intensity value, runoff depth value and pollutant concentration value of each spatial node at each time step are filled into a three-dimensional data array with time step number, spatial horizontal axis number and spatial vertical axis number as three dimensions, to obtain a three-dimensional spatiotemporal distribution feature data matrix. Based on the three-dimensional spatiotemporal distribution feature data matrix, Tucker decomposition is performed on the three-dimensional data array along the time dimension, decomposing the array into a precipitation modal factor matrix, a runoff spatial modal factor matrix, a pollutant concentration spatial modal factor matrix, and a core tensor reflecting the interaction relationship of the three modes, so as to obtain the decomposition result; Based on the decomposition results, the Pearson correlation coefficient between each column of the rainfall modal factor matrix and the rainfall intensity time series is calculated one by one. Load vectors with correlation coefficients exceeding the preset interference threshold are marked as rainfall interference terms and removed. Accordingly, the corresponding sub-tensor slices are deleted from the core tensor. Then, the pruned factor matrices and the core tensor are resynthesized along their respective modal directions to obtain a pure pollutant migration basis tensor stripped of direct rainfall interference.

3. The method for tracing and early warning of agricultural non-point source pollutants entering rivers according to claim 2, characterized in that, Step 2: Based on the fundamental tensor of pollutant migration, calculate the instantaneous rate of change of pollutant concentration at each spatial node in the runoff field point by point to obtain the two-dimensional concentration gradient evolution vector, including: Based on the pure pollutant migration fundamental tensor, the spatial distribution slices of pollutant concentration at two adjacent time points are extracted sequentially along the time dimension. For each spatial node, the concentration difference between the next time point and the previous time point is calculated and divided by the time interval to obtain the pollutant concentration temporal change rate. Based on the pollutant concentration change rate over time, adjacent spatial node pairs are extracted along the main direction of surface runoff and its orthogonal direction within the same time window. The ratio of the concentration difference to the spatial distance of each pair of nodes is calculated to obtain the spatial gradient of concentration in the horizontal and vertical directions. Based on the concentration-time change rate and the concentration spatial gradient in the horizontal and vertical directions, a three-dimensional vector is constructed, where the first component is the concentration-time change rate, the second component is the spatial gradient in the horizontal direction, and the third component is the spatial gradient in the vertical direction. By traversing all time windows and spatial nodes and assigning values ​​one by one, a two-dimensional concentration gradient evolution vector is obtained.

4. The method for tracing and early warning of agricultural non-point source pollutants entering rivers according to claim 3, characterized in that, Step 3: Based on the two-dimensional concentration gradient evolution vector, trace the instantaneous diffusion trajectory of pollutants in the complex surface runoff network to obtain the set of source tracing path nodes, including: Based on the two-dimensional concentration gradient evolution vector, the spatial node where the river cross-section enters the river is located is taken as the tracking starting point. The opposite direction of the spatial gradient component at the tracking starting point is taken as the instantaneous upward direction. The next tracking candidate spatial position is obtained by advancing along this direction with a preset tracking step size. Based on the next candidate spatial location, retrieve the two-dimensional concentration gradient evolution vector on the spatial node to which the location belongs. Perform a vector weighted average of the upward direction of the current tracking step and the opposite direction of the spatial gradient component of the node to obtain the corrected upward direction. Based on this, continue to advance to the next tracking point along the corrected direction with a preset tracking step size to obtain the advancement result. Based on the progress results, the direction retrieval, vector weighted average and step size advancement operations are performed repeatedly until the upstream boundary of the confluence zone is reached, and the tracing path is obtained by arranging the spatial nodes traversed by all tracking steps in chronological order. Based on the upstream path, all surface spatial nodes located in the farmland runoff catchment area along the upstream path are extracted, arranged in the tracing order, and duplicate nodes due to path intersections are removed to obtain the source path node set.

5. The method for tracing and early warning of agricultural non-point source pollutants entering rivers according to claim 4, characterized in that, Step 4: Combining the set of source tracing path nodes, calculate the total spatiotemporal cumulative flux of pollutants into the river within the confluence interval across the entire region, obtaining a continuous curve of river load, including: Based on the set of source tracing path nodes, the pollutant concentration value of each node at each monitoring time is extracted, and multiplied by the runoff depth and grid cell area at that node to obtain the instantaneous pollutant mass stock at each monitoring time. Based on the instantaneous pollutant mass stock at each node, the mass difference between upstream and downstream adjacent node pairs is calculated along the source tracing path from upstream to downstream. This difference is recorded as the comprehensive loss mass between adjacent node pairs caused by soil infiltration, microbial degradation and surface volatilization. Based on the comprehensive loss mass, the instantaneous pollutant mass stock of the downstream boundary node adjacent to the river section in the source tracing path is deducted segment by segment from the comprehensive loss mass between each adjacent node pair along the upstream direction to obtain the actual net inflow flux into the river at the boundary node after deducting the loss along the way. Based on the net inflow flux at each monitoring time, with the monitoring time as the horizontal axis and the net inflow flux as the vertical axis, all data points are connected in chronological order to form a continuous time variation curve, thus obtaining the continuous inflow load curve.

6. The method for tracing and early warning of agricultural non-point source pollutants entering rivers according to claim 5, characterized in that, Step 5: Extract the peak segment of the continuous curve of the inflow load into the river. Combine this with the water flow velocity and riverbed surface roughness at the river cross-section to estimate the advection transport distance of pollutants within the river channel, including: Based on the continuous curve of the load entering the river, the load value is scanned from the starting point to the end point in time intervals. The continuous time intervals in which the load value exceeds the preset peak value judgment threshold are extracted as peak segments. After traversing the entire curve, a set of peak segments is obtained. Based on the set of peak segments, the monitoring time with the largest load value is retrieved for each peak segment, and this time is marked as the peak inflow time corresponding to that segment. Based on the peak entry time into the river, retrieve the measured value of the water flow velocity at the river cross section at that time, and find the corresponding riverbed surface roughness parameter value from the preset Manning roughness coefficient reference table according to the riverbed material type. Based on the measured values ​​of water flow velocity and riverbed surface roughness parameters, the effective transport rate is obtained by subtracting the water flow velocity from the riverbed surface roughness. Then, the effective transport rate is multiplied by the time span from the peak entry time to the warning target time to obtain the advection transport distance of the pollutant cluster in the river channel corresponding to each peak segment.

7. The method for tracing and early warning of agricultural non-point source pollutants entering rivers according to claim 6, characterized in that, Step 6: Based on the advection transport distance and the set of source tracing path nodes, accurately locate the initial source of pollutant release within the farmland runoff catchment area, and obtain the spatial coordinates of the source, including: Based on the advection transport distance, the length of the river segment equidistant from the river cross-section entry point upstream along the river centerline is measured to determine the current river segment location coordinates of the pollutant plume after advection transport. Based on the current river section location coordinates, calculate the spatial Euclidean distance between each node and the river section location coordinates in the obtained source tracing path node set, and mark the node with the smallest distance as the spatial intersection node of the river channel and farmland runoff. Based on the spatial intersection nodes, all upstream nodes are traversed in reverse from the source tracing path upstream of the farmland. The pollutant concentration values ​​of each upstream node at the peak entry time into the river are extracted one by one. The ratio of the concentration value of each upstream node to the total concentration values ​​of all nodes in that upstream path is calculated to obtain the concentration contribution ratio of each node. Based on the proportion of concentration contribution of each node, the node with the largest proportion is identified as the dominant release node that contributes the most to the load of pollutants entering the river. The horizontal and vertical coordinates of this node in the farmland surface spatial grid are extracted to obtain the source spatial coordinates.

8. The method for tracing and early warning of agricultural non-point source pollutants entering rivers according to claim 7, characterized in that, Step 7: Perform a deep correlation between the source spatial coordinates and the peak time of the continuous curve of the inflow load to obtain non-point source pollution tracing and early warning control instructions, including: Based on the source spatial coordinates, a spatial point query is performed in the pre-constructed geographic information database of farmland runoff catchment areas to determine the farmland plot number and the runoff sub-catchment area identifier to which the coordinates fall. Based on the runoff sub-catchment area identification, the equipment identification code of the interception facility bound to the sub-catchment area is retrieved from the physical interception facility ledger to obtain the correspondence between the land parcel and the interception facility; Based on the continuous curve of the inflow load, the peak time of the source tracing path corresponding to the spatial coordinates of the source is extracted on the curve. The peak time, plot number and equipment identification code are bound together as a triple to obtain the source tracing-time association record. Based on the source-time correlation record and combined with the preset early warning lead time parameter, the peak time is shifted forward by this lead time to obtain the interception trigger time. The interception trigger time, plot number and equipment identification code are combined into control instructions, and the sequence of control instructions for non-point source pollution source tracing and early warning is obtained by sorting them according to the order of the interception trigger times.

9. The method for tracing and early warning of agricultural non-point source pollutants entering rivers according to claim 8, characterized in that, Step 8: Based on the non-point source pollution tracing and early warning control command, trigger the physical interception facilities in the farmland runoff collection area, and output the tracing and early warning results of the agricultural non-point source pollutant entry process into the river, including: According to the sequence of control instructions for tracing and early warning of non-point source pollution, when the time reaches the interception trigger time specified in the instruction, a start signal is sent to the drive controller of the physical interception facility pointed to by the instruction, driving the trash rack or interception net to be lowered from the standby position to the interception station in the runoff channel to complete the interception action; Based on the completed interception action, the water quality verification sensor deployed downstream of the interception facility is activated. It continuously samples at a preset frequency within a preset verification time window after the interception trigger time to obtain the measured value sequence of pollutant concentration within the verification window. The measured value sequence is then compared point by point with the background concentration value sequence of the same duration before the interception is executed to obtain the interception effect comparison data. Based on the interception effect comparison data, calculate the percentage decrease of the average concentration after interception relative to the average background concentration before interception. If the percentage exceeds the preset effective interception judgment threshold, the interception is deemed effective; otherwise, generate facility abnormality alarm information and push it to the maintenance management terminal. Based on the interception effectiveness determination, the source spatial coordinates, peak time, interception trigger time, interception effect comparison data and determination conclusion are summarized and packaged to obtain a comprehensive report containing source tracing path information, early warning trigger information and interception effect information, which is output as source tracing early warning results.

10. A source tracing and early warning system for agricultural non-point source pollutants entering rivers, the system implementing the method as described in any one of claims 1 to 9, characterized in that, include: The acquisition module is used to acquire multi-source hydrological monitoring data from farmland runoff catchment areas, river cross-sections, and soil permeable layers, and to construct a three-dimensional spatiotemporal distribution feature data matrix that includes rainfall runoff generation processes and pollutant concentration fields. Tensor decomposition is performed on the three-dimensional spatiotemporal distribution feature data matrix to remove the rainfall interference term and runoff coupling term, resulting in a pure pollutant migration fundamental tensor. The calculation module is used to calculate the instantaneous rate of change of pollutant concentration at each spatial node in the runoff field based on the pollutant migration tensor, and obtain the two-dimensional concentration gradient evolution vector. The tracking module is used to track the instantaneous diffusion trajectory of pollutants in a complex surface runoff network based on the two-dimensional concentration gradient evolution vector, and obtain the set of source tracing path nodes; The accumulation module is used to combine the set of source tracing path nodes to calculate the total spatiotemporal cumulative flux of pollutants into the river within the confluence interval, and obtain a continuous curve of the inflow load. The extrapolation module is used to extract the peak section of the continuous curve of the inflow load into the river, and combine the water flow velocity and riverbed surface roughness of the river cross section to extrapolate the advection transport distance of pollutants in the river channel. The positioning module is used to accurately locate the initial release source of pollutants in the farmland runoff catchment area based on the advection transport distance and the set of source tracing path nodes, and obtain the spatial coordinates of the source. The correlation module is used to deeply correlate the spatial coordinates of the source with the peak time of the continuous curve of the inflow load to obtain the source tracing, early warning and control instructions for non-point source pollution. The processing module is used to trigger physical interception facilities in the farmland runoff collection area according to the non-point source pollution tracing and early warning control instructions, and output the tracing and early warning results of the process of agricultural non-point source pollutants entering the river.

Images