Distributed model-oriented small and medium-sized basin sub-basin division method

Abstract

The invention discloses a distributed model-oriented medium and small watershed sub-watershed division method, which comprises the following steps of S1, collecting digital elevation model data, performing depression filling processing, and extracting preliminary river network information; s2, extracting an actual water system, and calculating a river network main stream grid elevation difference; s3, adjusting digital elevation model data according to the elevation difference, and optimizing a river network extraction result; s4, basin morphological parameters are calculated, and basin classification is carried out through a manifold learning method; s5, constructing a sub-basin optimal threshold calculation model, and dynamically adjusting a sub-basin division threshold; s6, flow direction calculation, confluence analysis and sub-basin division are carried out, and sub-basin boundaries are optimized; and S7, verification is carried out in combination with the distributed hydrological model, the influence of the division result is analyzed, and sub-basin division parameters are optimized. According to the method, multi-scale features, dynamic threshold optimization and hydrological network topology are fused, sub-basin division precision, boundary continuity and hydrological simulation stability are improved, and intelligent basin division is realized.

Description

Technical Field

[0001] The present invention relates to the technical field of water resource management, and particularly to a method for dividing sub-watersheds of small and medium-sized watersheds for distributed models. Background Art

[0002] In the fields of hydrological simulation and watershed management, sub-watershed division is one of the key technical links in water resource management, flood prediction, ecological environment assessment, etc. The rationality of sub-watershed division directly affects the calculation accuracy of hydrological models and determines the reliability of watershed hydrodynamic simulation. At present, traditional sub-watershed division methods mainly rely on Digital Elevation Model (DEM) data to determine sub-watershed boundaries through flow direction calculation, catchment area analysis, and setting of hydrological characteristic parameters. However, these methods have certain limitations and are difficult to adapt to the complex terrain characteristics and dynamic hydrological processes of watersheds, resulting in the lack of scientificity and adaptability of the division results.

[0003] Existing sub-watershed division methods usually rely on a fixed catchment area threshold for division, that is, an empirical value is set as the minimum catchment area threshold of the sub-watershed. All areas larger than this threshold are identified as sub-watersheds, and areas smaller than this threshold are merged into adjacent sub-watersheds. Although this method can provide a certain division reference under regular terrain, for watersheds with complex geomorphic features, such as areas with uneven distribution of mountains, hills, and lakes, it is difficult to accurately depict the watershed morphology. In addition, the setting of the fixed threshold lacks self-adaptability and cannot be flexibly adjusted in different types of watersheds, resulting in a large deviation between the division result and the actual hydrological process.

[0004] During the watershed division process, traditional methods mainly rely on flow direction calculation based on DEM data and use the D8 (unidirectional flow direction) algorithm for water flow direction calculation. Although these methods can describe the water flow direction to a certain extent, when encountering flat areas, depressions, or areas with large elevation errors, unreasonable flow direction results may occur, leading to discontinuous sub-watershed boundaries and even the phenomenon of "reverse water flow", affecting the accuracy of hydrological simulation. At the same time, existing flow direction calculation methods usually ignore the hydrodynamic connection relationship between sub-watersheds and only consider terrain features, failing to dynamically optimize and adjust sub-watersheds in combination with the watershed hydrological process, further reducing the rationality of the division.

[0005] Another key challenge in sub - watershed division lies in boundary optimization. Traditional sub - watershed division methods usually directly perform threshold segmentation based on catchment area, without considering the smoothness and coherence of sub - watershed boundaries, resulting in a large number of "fragmented" sub - watersheds in the division results, that is, regions with too small area but forming independent domains. These small sub - watersheds may cause numerical instability in hydrological simulations and affect the stability of the overall hydrodynamic calculation. Some studies have tried to optimize the division results by manually adjusting sub - watershed boundaries or merging small sub - watersheds, but these methods rely on artificial settings and are difficult to adapt to different types of watershed characteristics.

[0006] In recent years, with the development of Geographic Information System (GIS) technology, artificial intelligence, and machine learning technology, some scholars have tried to introduce data - driven methods to optimize the sub - watershed division process. For example, some studies combine clustering algorithms for sub - watershed optimization, but most of these methods rely on empirical parameter adjustment and still have certain limitations. In addition, traditional methods have not effectively utilized watershed morphological characteristics, such as curvature, fractal dimension, complexity of confluence paths, etc., resulting in the failure to fully consider the overall morphology of the watershed during the sub - watershed division process and affecting the scientific nature of the division.

[0007] Therefore, how to provide a sub - watershed division method for small and medium - sized watersheds for distributed models is an urgent problem that needs to be solved by those skilled in the art. Summary of the Invention

[0008] An object of the present invention is to propose a sub - watershed division method for small and medium - sized watersheds for distributed models. The present invention makes full use of technologies such as watershed morphological feature analysis, multi - scale feature fusion, dynamic threshold optimization, hydrological network topology construction, and neural network optimization of sub - watershed division, details the calculation method of intelligent sub - watershed division, and has the advantages of strong adaptability, high calculation accuracy, good continuity of boundary optimization, and strong stability of hydrological simulation.

[0009] The sub - watershed division method for small and medium - sized watersheds for distributed models according to the embodiments of the present invention includes the following steps:

[0010] S1. Obtain the digital elevation model data and remote sensing image data of the study area, perform depression filling on the digital elevation model data, extract preliminary river network information, and identify the main stream area;

[0011] S2. Extract the actual water system based on the remote sensing image data, perform overlay analysis with the preliminary river network information, identify the deviation area of the main stream grid, and calculate the elevation difference between the elevation value of the main stream grid of the river network and the elevation of the surrounding depression grids;

[0012] S3. Based on the calculated elevation difference, perform local elevation adjustment on the digital elevation model data, re - execute depression filling, flow direction calculation, and cumulative runoff calculation, and optimize the river network extraction result;

[0013] S4. Calculate the basin morphological parameters based on the optimized river network extraction results, classify the basins using the manifold learning method, and identify different basin types.

[0014] S5. Construct a calculation model for the optimal threshold of sub-basins based on the basin morphological parameters, calculate the threshold of catchment area, and dynamically adjust the threshold of sub-basin division.

[0015] S6. Based on the threshold of catchment area, perform flow direction calculation and confluence analysis on the adjusted digital elevation model data, divide the sub-basins, optimize the boundaries in combination with the basin morphological parameters, and perform merging processing on small-area sub-basins.

[0016] S7. Import the sub-basin division results into the distributed hydrological model, analyze the impact of sub-basin division on basin runoff simulation, and adjust the sub-basin division parameters according to the hydrological simulation results.

[0017] Optionally, the S2 specifically includes:

[0018] S21. Extract the actual water system based on remote sensing image data, use a water body recognition model that fuses graph attention network and basin hydrodynamic characteristics for water body detection, and calculate the water body classification confidence S by combining spectral characteristics, terrain characteristics, and hydrodynamic simulation information w :

[0019]

[0020] where σ is the normalization activation function, is the adjacency weight calculated based on the attention mechanism, is the water body feature of node j at the k-th scale, ω is the weight of hydrodynamic information, and F(i, j) is the hydrodynamic feature;

[0021] S22. Perform spatial interpolation processing on the extracted preliminary river network information, and use the high-order interpolation method with adaptive weights to calculate the elevation value Z r (i, j) of the main river network raster, and store it as the interpolated river network raster data R(i, j):

[0022]

[0023] where Z(m, n) is the known elevation value of the surrounding raster, μ is the adjustment coefficient, is the second-order Laplacian smoothing term, w m,n is the adaptive weight, M and N are the number of rasters, and Z(i, j) is the elevation value of the digital elevation model data at the raster coordinates (i, j);

[0024] S23. Binarize the confidence of water body classification to generate water body distribution raster data W(i, j), and perform raster overlay operation with the interpolated river network raster data R(i, j) obtained, and calculate the consistency error using cross-entropy error:

[0025]

[0026] where E is the consistency error value. If E > T e , it is determined as the main stream raster deviation area, and T e is the error threshold, and ln is the logarithmic function;

[0027] S24. Calculate the elevation difference for the identified main stream raster deviation area to solve the elevation difference ΔZ max :

[0028] ΔZ max = max(|Z r (i, j) - Z w (i, j)|) + γ·σ Z ;

[0029] where σ Z is the elevation standard deviation of the adjacent area, γ is the adaptive correction coefficient, max is the maximum value operation, Z r (i, j) is the elevation value of the river network main stream raster, and Z w (i, j) is the elevation value of the surrounding depression raster;

[0030] S25. Store the calculated main stream raster deviation area and elevation difference data.

[0031] Optionally, the specific content of S3 includes:

[0032] S31. Based on the calculated local elevation difference ΔZ max , adjust the elevation of the deviation area raster to construct a local elevation adjustment model:

[0033] Z′ r (i, j) = Z r (i, j) - α·ΔZ max ;

[0034] where Z′ r (i, j) is the adjusted elevation value of the river network main stream raster, α is the local elevation adjustment coefficient, and Z r (i, j) is the elevation value of the river network main stream raster before adjustment;

[0035] S32. Use the Markov random field model to smooth the elevation adjustment area and construct a local elevation optimization model:

[0036]

[0037] Among them, is the optimized raster elevation value of the main river network, and argmin Z′ is to solve the candidate elevation adjustment value Z' to minimize the objective function value, and Z' r (m, n) represents the adjusted elevation value of the neighborhood raster (m, n), λ is the regularization parameter, H(Z') is the regularization term of elevation change, N(i, j) is the neighborhood set of the raster, Z' is the candidate elevation adjustment value, m and n are the row index and column index of the neighborhood raster, and β is the neighborhood smoothing factor;

[0038] S33. Based on the optimized digital elevation model data Perform a depression filling operation to fill all internal depressions and calculate the digital elevation model data after elevation adjustment:

[0039]

[0040] Among them, Z f (i, j) is the elevation value of the digital elevation model data after depression filling, and Z w (i, j) is the elevation value of the adjacent depression raster;

[0041] S34. Based on the generated digital elevation model data Z f (i, j) after depression filling, use the directional weighted flow direction calculation method to calculate the flow direction matrix D(i, j):

[0042]

[0043] Among them, Z f (m, n) is the elevation value of the neighborhood raster (m, n), d m,n is the Euclidean distance between the raster (i, j) and the neighborhood raster (m, n), ε is the flow direction adjustment index, and argmax (m,n)∈N(i,j) is to determine the adjacent raster (m, n) to which the water flow from the current raster (i, j) flows;

[0044] Among them, the catchment area of the river network raster is calculated using the adaptive confluence weight model:

[0045]

[0046] Among them, A(i, j) is the catchment area of the raster (i, j), exp is the exponential function, ∈ is the confluence weight adjustment parameter, A(m, n) is the catchment area of the neighborhood raster (m, n), and d p,q is the Euclidean distance between the raster (i, j) and the neighborhood raster (p, q);

[0047] S35. Optimize the river network extraction result based on the calculated flow direction matrix and cumulative catchment area:

[0048]

[0049] Among them, is the optimized raster data for river network extraction, and T A is the minimum catchment area threshold for river network extraction.

[0050] Optionally, the S4 specifically includes:

[0051] S41. Calculate the basin morphological parameters based on the optimized raster data for river network extraction, including the curvature C, fractal dimension D, complexity of the confluence path P, and topological centrality T c , and construct a basin morphological feature set:

[0052]

[0053] Among them, L r is the actual channel length of the main stream of the basin, and L s is the shortest straight-line distance of the main stream of the basin, is the minimum number of grids required to cover the basin boundary at scale , is the grid scale, Q is the total number of confluence paths, and L i is the length of the i-th confluence path in the sub-basin, ζ is the hydrodynamic adjustment factor, Q max and Q min are the maximum and minimum flow rates in the sub-basin, N(i) is the set of sub-basins connected to sub-basin i, and W ij is the confluence weight between sub-basin i and sub-basin j, and W ik is the confluence weight between sub-basin i and sub-basin K;

[0054] S42. Construct a basin morphological feature matrix M f :

[0055]

[0056] Among them, C n , D n , P n , T cn are the curvature, fractal dimension, complexity of the confluence path, and topological centrality of the n-th sub-basin, respectively, and N is the total number of sub-basins;

[0057] S43. Based on the basin morphological feature matrix, use a graph neural network based on an adaptive attention mechanism for basin classification to construct a morphological classification model:

[0058]

[0059] Among them, is the feature vector of basin i in the (l + 1)-th layer network, is the feature vector of basin j in the l-th layer network, and σ is a non-linear activation function. is the adjacency weight calculated based on the attention mechanism, exp is the exponential function, τ is the weight vector in the attention mechanism, and W (l( is the transformation matrix of the l-th layer network, and LeakyReLU is the leaky rectified linear unit activation function. is the feature vector of basin i in the l-th layer network, is the feature vector of basin k in the l-th layer network;

[0060] S44. Based on the calculated basin morphology classification results, construct a basin classification label matrix based on soft clustering and neighborhood optimization, and use a joint optimization method to calculate the final classification label.

[0061]

[0062] Among them, d(T n , M f ) is the Euclidean distance between the morphological features of sub-basin n and the centroid of category T n , K is the total number of classification categories, is the temperature adjustment parameter, d(T j , M f ) is the Euclidean distance between sub-basin j and the centroid of category T f in the morphological feature matrix M j , d(T k , M f ) is the Euclidean distance between sub-basin k and the centroid of category T f in the morphological feature matrix M k , ω n,j is the hydrological connectivity weight, T n is the classification label, and μ is the neighborhood influence factor. is to select a category from the classification label T n that makes the total score the highest;

[0063] S45. Based on the obtained optimized basin classification label, generate a visualized basin classification map and store the classification results.

[0064] Optionally, the specific content of S5 includes:

[0065] S51. Based on the stored classification results and the basin morphological feature matrix M f, construct a multi-scale feature fusion model, introduce basin morphology, topographic gradient, and hydrodynamic characteristics, and calculate the threshold of the catchment area of the sub-basins

[0066]

[0067] Among them, F k is the comprehensive morphological feature of sub-basin k, β k is the optimized weight coefficient, N(n) is the set of neighboring sub-basins of sub-basin n, ω n,j is the neighborhood influence weight, F j is the comprehensive morphological feature of sub-basin j, α o is the neighborhood weighted influence factor;

[0068] S52. Based on the calculated threshold of the catchment area of the sub-basins, combine the information entropy of the sub-basins to calculate the dynamic threshold finally used for sub-basin division

[0069]

[0070] Among them, η and α e are the dynamic adjustment coefficients, E A is the information entropy of the sub-basin, is the average catchment area threshold, is the standard deviation of the catchment area threshold;

[0071] S53. Based on the calculated dynamic threshold, optimize the sub-basin division:

[0072]

[0073] Among them, ω n,j is the neighborhood clustering weight, N s is the number of sub-basin splits, T A,min is the minimum area threshold of the sub-basin, T A,max is the maximum area threshold of the sub-basin, is the catchment area threshold of neighboring sub-basin j;

[0074] S54. Based on the calculated optimized sub-basin data, store the final sub-basin division result and output the sub-basin division matrix M sub :

[0075]

[0076] Among them, X N and Y N are the geographical coordinates of sub-basin n, is the final catchment area threshold of sub-basin n;

[0077] Store the calculated sub-basin division data to support the sub-basin division calculation in step S6.

[0078] Optionally, the specific steps of S6 include:

[0079] S61. Calculate the multi-scale hydrological topological flow matrix D(i,j) based on the catchment area threshold, and construct the basin topological network:

[0080] D(i,j) = argmax (m,n)∈N(i,j) (Z f (i,j) - Z f (m,n)) + γ 1 ·Φ(i,j);

[0081] Where, Z f (i,j) is the elevation value of the grid (i,j), Z f (m,n) is the elevation value after depression filling of the neighborhood grid (m,n), γ 1 is the hydrological topological optimization factor, Φ(i,j) is the water flow clustering topological optimization function, argmax (m,n)∈N(i,j) is to select the most likely flow direction of the water flow among all neighborhood candidate flow directions, and N(i,j) is the neighborhood of the grid (i,j);

[0082] S62. Based on the calculated flow matrix, adopt the cumulative catchment area calculation method optimized by information entropy:

[0083] A(i,j) = ∑ (m,n)∈N(i,j) A(m,n) + W(i,j)·(1 + γ 2 ·E A );

[0084] Where, A(i,j) is the catchment area of the grid (i,j), W(i,j) is the basin contribution rate weight, γ 2 is the information entropy adjustment factor, E A is the information entropy of the sub-basin, and A(m,n) is the cumulative catchment area of the upstream neighborhood grid (m,n);

[0085] S63. Based on the calculated catchment area, combined with the calculated dynamic catchment area threshold Adopt a dynamic neural network to optimize the sub-basin division:

[0086]

[0087] Where, M sub (i,j) is the sub-basin division result matrix, σ is the adaptive classification activation function, ω n,j is the neighborhood influence weight, γ 3 is the neighborhood smoothing adjustment parameter;

[0088] S64. Optimize the sub - basin boundaries based on the calculated sub - basin division results and adjust the merger relationships of small - area sub - basins:

[0089]

[0090] Among them, B(i, j) is the optimized sub - basin boundary matrix, and δ(M sub (m, n), M sub (i, j)) is the sub - basin merger function, and T A,min is the minimum area threshold of the sub - basin;

[0091] S65. Store the final optimized sub - basin data based on the calculated sub - basin division results.

[0092] The beneficial effects of the present invention are as follows:

[0093] By introducing basin - form constraints, the present invention proposes an adaptive sub - basin division method. Compared with the traditional fixed - threshold division method, the present invention can adaptively adjust the catchment area threshold of sub - basins, making the division results more in line with the actual hydrological characteristics of different types of basins. The present invention adopts a dynamic - threshold optimization method based on multi - scale feature fusion, comprehensively considering basin curvature, fractal dimension, complexity of confluence paths, topological centrality, slope, hydrodynamic characteristics, etc., so that the sub - basin division not only depends on a single hydrological parameter, but combines the overall morphological characteristics of the basin, thereby improving the scientificity and adaptability of the division.

[0094] By optimizing the flow - direction calculation method and adopting a hydrological network topology optimization technology based on water - flow clustering, the present invention enhances the hydrodynamic connectivity between sub - basins and reduces the calculation errors caused by sudden changes in flow direction or unreasonable flow directions. Traditional flow - direction calculation methods are prone to problems such as flow - direction breaks or incorrect divisions when dealing with depressions, flat areas or complex terrains. However, the present invention introduces a flow - direction optimization matrix and combines information entropy to optimize the calculation of cumulative catchment area, enabling the sub - basin boundaries to more reasonably adapt to terrain changes and improving the stability and accuracy of hydrological simulation calculations.

[0095] The present invention adopts an optimized sub - basin division method based on a dynamic neural network, introduces the Wasserstein distance to calculate the neighborhood influence weight, ensures the smoothness of the sub - basin division boundary, and avoids the problem of fragmented sub - basins commonly found in traditional methods. By introducing a dynamic activation function during the sub - basin division process, the continuity and self - adaptability of the sub - basin division are ensured, reducing the numerical instability problems caused by boundary jumps during the calculation process and improving the interpretability of the calculation results.

[0096] In addition, in the sub-basin optimization stage, the present invention combines the Wasserstein distance merging method, enabling small-area sub-basins to be intelligently merged into the optimal adjacent sub-basins, avoiding the influence of overly small sub-basins on the stability of hydrological calculations. Meanwhile, for sub-basins with overly large areas, the present invention uses the minimum cut maximum flow method to subdivide the sub-basins, making the sub-basin division conform to the true distribution of the basin hydrological process and improving the calculation accuracy of hydrological simulation.

[0097] Generally speaking, the present invention provides a more intelligent and adaptive sub-basin division method, overcoming the limitations of the fixed threshold division method in the prior art and improving the rationality of sub-basin boundaries and calculation stability. The present invention can be applied to the basin division of different terrains and hydrological conditions, providing a more accurate division scheme for water resource management, flood prediction, hydrological model optimization, etc., and enhancing the reliability and practicality of basin hydrological calculations. BRIEF DESCRIPTION OF THE DRAWINGS

[0098] The drawings are used to provide a further understanding of the present invention and constitute a part of the specification. Together with the embodiments of the present invention, they are used to explain the present invention and do not constitute a limitation to the present invention. In the drawings:

[0099] Figure 1 is a flowchart of the sub-basin division method for small and medium-sized basins for a distributed model proposed by the present invention;

[0100] Figure 2 is a schematic diagram of the process of extracting the basin morphological characteristics of the sub-basin division method for small and medium-sized basins for a distributed model proposed by the present invention. DETAILED DESCRIPTION OF THE EMBODIMENTS

[0101] Now, the present invention will be further described in detail with reference to the drawings. These drawings are all simplified schematic diagrams, only illustrating the basic structure of the present invention in a schematic manner, and therefore only showing the components related to the present invention.

[0102] Refer to Figure 1-2 , the sub-basin division method for small and medium-sized basins for a distributed model includes the following steps:

[0103] S1. Obtain the digital elevation model data and remote sensing image data of the study area, perform depression filling on the digital elevation model data, extract preliminary river network information, and identify the main stream area;

[0104] S2. Extract the actual water system based on the remote sensing image data, perform overlay analysis with the preliminary river network information, identify the deviation area of the main stream grid, and calculate the elevation difference between the elevation values of the main stream grid of the river network and the surrounding depression grids;

[0105] S3. Based on the calculated elevation difference, locally adjust the digital elevation model data, re - perform the depression filling process, flow direction calculation, and cumulative flow accumulation calculation to optimize the river network extraction result;

[0106] S4. Based on the optimized river network extraction result, calculate the basin morphological parameters, and classify the basin based on the manifold learning method to identify different basin types;

[0107] S5. Construct a calculation model for the optimal threshold of sub - basins based on the basin morphological parameters, calculate the threshold of catchment area, and dynamically adjust the threshold of sub - basin division;

[0108] S6. Based on the threshold of catchment area, perform flow direction calculation and confluence analysis on the adjusted digital elevation model data, divide sub - basins, optimize the boundaries in combination with the basin morphological parameters, and perform merging processing on small - area sub - basins;

[0109] S7. Import the sub - basin division result into the distributed hydrological model, analyze the impact of sub - basin division on basin runoff simulation, and adjust the sub - basin division parameters according to the hydrological simulation results.

[0110] In this embodiment, the S2 specifically includes:

[0111] S21. Extract the actual water system based on remote sensing image data, use a water body recognition model that fuses graph attention network and basin hydrodynamic characteristics for water body detection, and combine spectral characteristics, topographic characteristics, and hydrodynamic simulation information to calculate the water body classification confidence S w :

[0112]

[0113] where σ is the normalization activation function, is the adjacency weight calculated based on the attention mechanism, is the water body feature of node j at the k - th scale, ω is the weight of hydrodynamic information, and F(i, j) is the hydrodynamic feature;

[0114] S22. Perform spatial interpolation processing on the extracted preliminary river network information, and use the high - order interpolation method with adaptive weights to calculate the elevation value Z r (i, j) of the river network main stream raster, and store it as the interpolated river network raster data R(i, j):

[0115]

[0116] where Z(m, n) is the known elevation value of the surrounding raster, μ is the adjustment coefficient, is the second - order Laplacian smoothing term, w m,nis the adaptive weight, M and N are the number of grids, and Z(i, j) is the elevation value of the digital elevation model data at the grid coordinates (i, j);

[0117] S23. Binarize the water body classification confidence to generate the water body distribution grid data W(i, j), and perform a grid overlay operation with the interpolated river network grid data R(i, j) obtained, and calculate the consistency error using the cross-entropy error:

[0118]

[0119] where E is the consistency error value. If E > T e , it is determined as the main stream grid deviation area, and T e is the error threshold, and ln is the logarithmic function;

[0120] S24. Calculate the elevation difference for the identified main stream grid deviation area to solve the elevation difference ΔZ max :

[0121] ΔZ max = max(|Z r (i, j) - Z w (i, j)|) + γ·σ Z ;

[0122] where σ Z is the elevation standard deviation of the adjacent area, γ is the adaptive correction coefficient, max is the maximum value operation, Z r (i, j) is the elevation value of the river network main stream grid, and Z w (i, j) is the elevation value of the surrounding depression grid;

[0123] S25. Store the calculated main stream grid deviation area and elevation difference data.

[0124] In this embodiment, the specific steps of S3 include:

[0125] S31. Based on the calculated local elevation difference ΔZ max , adjust the elevation of the deviation area grid to construct a local elevation adjustment model:

[0126] Z′ r (i, j) = Z r (i, j) - α·ΔZ max ;

[0127] where Z′ r (i, j) is the adjusted elevation value of the river network main stream grid, α is the local elevation adjustment coefficient, and Z r (i, j) is the elevation value of the river network main stream grid before adjustment;

[0128] S32. Use the Markov random field model to smooth the elevation adjustment area and construct a local elevation optimization model:

[0129]

[0130] Among them, is the elevation value of the main river network raster after optimization, argmin Z′ is to solve the candidate elevation adjustment value Z' to minimize the objective function value, Z' r (m,n) is the adjusted elevation value of the neighborhood raster (m,n), λ is the regularization parameter, H(Z') is the regularization term of elevation change, N(i,j) is the neighborhood set of the raster, Z' is the candidate elevation adjustment value, m and n are the row index and column index of the neighborhood raster, and β is the neighborhood smoothing factor;

[0131] S33. Based on the optimized digital elevation model data Perform a depression filling operation to fill all internal depressions and calculate the digital elevation model data after elevation adjustment:

[0132]

[0133] Among them, Z f (i,j) is the elevation value of the digital elevation model data after depression filling, Z w (i,j) is the elevation value of the adjacent depression raster;

[0134] S34. Based on the generated digital elevation model data Z f (i,j) after depression filling, use the direction weighted flow calculation method to calculate the flow direction matrix D(i,j):

[0135]

[0136] Among them, Z f (m,n) is the elevation value of the neighborhood raster (m,n), d m,n is the Euclidean distance between the raster (i,j) and the neighborhood raster (m,n), ε is the flow direction adjustment index, argmax (m,n)∈N(i,j) is to determine the adjacent raster (m,n) to which the water flow from the current raster (i,j) flows;

[0137] Among them, the catchment area of the river network raster is calculated using the adaptive confluence weight model:

[0138]

[0139] Among them, A(i,j) is the catchment area of the raster (i,j), exp is the exponential function, ∈ is the confluence weight adjustment parameter, A(m,n) is the catchment area of the neighborhood raster (m,n), dp,q is the Euclidean distance between the grid (i, j) and the neighboring grid (p, q);

[0140] S35. Optimize the river network extraction result based on the calculated flow direction matrix and cumulative catchment area:

[0141]

[0142] Among them, is the optimized river network extraction grid data, and T A is the minimum catchment area threshold for river network extraction.

[0143] In this embodiment, the S4 specifically includes:

[0144] S41. Calculate the basin morphological parameters based on the optimized river network extraction grid data, including the curvature C, fractal dimension D, confluence path complexity P, and topological centrality T c , and construct a basin morphological feature set:

[0145]

[0146] Among them, L r is the actual channel length of the main stream of the basin, L s is the shortest straight-line distance of the main stream of the basin, is the minimum number of grids required to cover the basin boundary at the scale , is the grid scale, Q is the total number of confluence paths, L i is the length of the i-th confluence path in the sub-basin, ζ is the hydrodynamic adjustment factor, Q max and Q min are the maximum and minimum flow rates in the sub-basin, N(i) is the set of sub-basins connected to the sub-basin i, W ij is the confluence weight between the sub-basin i and the sub-basin j, and W ik is the confluence weight between the sub-basin i and the sub-basin K;

[0147] S42. Construct a basin morphological feature matrix M f :

[0148]

[0149] Among them, C n , D n , P n , T cn are the curvature, fractal dimension, confluence path complexity, and topological centrality of the n-th sub-basin, and N is the total number of sub-basins;

[0150] S43. Based on the basin morphological feature matrix, use a graph neural network with an adaptive attention mechanism to classify basins and construct a morphological classification model:

[0151]

[0152] Among them, is the feature vector of basin i in the (l + 1)-th layer network, is the feature vector of basin j in the l-th layer network, σ is a non-linear activation function, is the adjacency weight calculated based on the attention mechanism, exp is the exponential function, τ is the weight vector in the attention mechanism, W (l) is the transformation matrix of the l-th layer network, LeakyReLU is the leaky rectified linear unit activation function, is the feature vector of basin i in the l-th layer network, is the feature vector of basin k in the l-th layer network;

[0153] S44. Based on the calculated basin morphological classification results, construct a basin classification label matrix based on soft clustering and neighborhood optimization, and use a joint optimization method to calculate the final classification label

[0154]

[0155] Among them, d(T n , M f ) is the Euclidean distance between the morphological features of sub-basin n and the centroid of category T n , K is the total number of classification categories, is the temperature adjustment parameter, d(T j , M f ) is the Euclidean distance between sub-basin j and the centroid of category T f in the morphological feature matrix M j , d(T k , M f ) is the Euclidean distance between sub-basin k and the centroid of category T f in the morphological feature matrix M k , ω n,j is the hydrological connectivity weight, T n is the classification label, μ is the neighborhood influence factor, is to select a category from the classification label T n that makes the total score the highest;

[0156] S45. Based on the obtained optimized basin classification label, generate a visualized basin classification map and store the classification results.

[0157] In this embodiment, the specific content of S5 includes:

[0158] S51. Based on the stored classification results and the basin morphological feature matrix M f , construct a multi-scale feature fusion model, introduce basin morphology, topographic gradient, and hydrodynamic characteristics, and calculate the threshold of the catchment area of sub-basins

[0159]

[0160] where F k is the comprehensive morphological feature of sub-basin k, β k is the optimized weight coefficient, N(n) is the set of neighboring sub-basins of sub-basin n, ω n,j is the neighboring influence weight, F j is the comprehensive morphological feature of sub-basin j, and α o is the neighboring weighted influence factor;

[0161] S52. Based on the calculated threshold of the catchment area of sub-basins, combine the information entropy of sub-basins to calculate the final dynamic threshold for sub-basin division

[0162]

[0163] where η and α e are dynamic adjustment coefficients, E A is the information entropy of sub-basins, is the average threshold of the catchment area, is the standard deviation of the threshold of the catchment area;

[0164] S53. Based on the calculated dynamic threshold, optimize sub-basin division:

[0165]

[0166] where ω n,j is the neighboring clustering weight, N s is the number of sub-basin splits, T A,min is the minimum area threshold of sub-basins, T A,max is the maximum area threshold of sub-basins, is the threshold of the catchment area of neighboring sub-basin j;

[0167] S54. Based on the calculated optimized sub-basin data, store the final sub-basin division results and output the sub-basin division matrix M sub :

[0168]

[0169] where X N and Y N are the geographical coordinates of sub-basin n, is the final catchment area threshold for sub - basin n;

[0170] S55. Store the calculated sub - basin division data to support the sub - basin division calculation in step S6.

[0171] In this embodiment, the said S6 specifically includes:

[0172] S61. Based on the catchment area threshold, calculate the multi - scale hydrological topological flow matrix D(i,j) and construct the basin topological network:

[0173] D(i,j) = argmax (m,n)∈N(i,j) (Z f (i,j) - Z f (m,n)) + γ 1 ·Φ(i,j);

[0174] Among them, Z f (i,j) is the elevation value of grid (i,j), Z f (m,n) is the elevation value of the filled - depression neighborhood grid (m,n), γ 1 is the hydrological topological optimization factor, Φ(i,j) is the water flow clustering topological optimization function, argmax (m,n)∈N(i,j) is to select the most likely flow direction of water flow among all neighborhood candidate flow directions, and N(i,j) is the neighborhood of grid (i,j);

[0175] S62. Based on the calculated flow matrix, adopt the cumulative catchment area calculation method optimized by information entropy:

[0176] A(i,j) = ∑ (m,n)∈N(i,j) A(m,n) + W(i,j)·(1 + γ 2 ·E A );

[0177] Among them, A(i,j) is the catchment area of grid (i,j), W(i,j) is the basin contribution rate weight, γ 2 is the information entropy adjustment factor, E A is the information entropy of the sub - basin, and A(m,n) is the cumulative catchment area of the upstream neighborhood grid (m,n);

[0178] S63. Based on the calculated catchment area and combined with the calculated dynamic catchment area threshold adopt a dynamic neural network to optimize the sub - basin division:

[0179]

[0180] Among them, M sub (i,j) is the sub - basin division result matrix, σ is the adaptive classification activation function, ω n,jis the neighborhood influence weight, γ 3 is the neighborhood smoothing adjustment parameter;

[0181] S64. Based on the calculated sub-basin division results, optimize the sub-basin boundaries and adjust the merging relationships of small-area sub-basins:

[0182]

[0183] Among them, B(i, j) is the optimized sub-basin boundary matrix, and δ(M sub (m, n), M sub (i, j)) is the sub-basin merging function, and T A,min is the minimum area threshold of the sub-basin;

[0184] S65. Based on the calculated sub-basin division results, store the final optimized sub-basin data.

[0185] Example 1:

[0186] To verify the feasibility of the present invention in implementation, the present invention is applied to a hilly terrain basin for experimental analysis, and comparative tests are carried out from four aspects: sub-basin division accuracy, boundary smoothness, calculation efficiency, and hydrological simulation effect.

[0187] In the experiment, a hilly basin was selected for sub-basin division testing. During the experiment, first, the elevation information of the basin was obtained using DEM data, and the method of the present invention was used for sub-basin division. The basin morphological feature analysis method was adopted to extract key parameters such as basin curvature, fractal dimension, and confluence path complexity. Based on these features, the dynamic catchment area threshold was calculated to perform adaptive optimization of the sub-basins. The hydrological network topology optimization flow direction calculation was adopted to ensure hydrodynamic connectivity and improve the rationality of flow direction calculation.

[0188] During the experiment, in order to compare the advantages and disadvantages of the present invention with traditional methods, the fixed threshold method, the GIS-based division method, and the method of the present invention were respectively used for sub-basin division to test their calculation efficiency, sub-basin boundary continuity, and hydrological simulation accuracy. Multiple sub-basins were selected in each basin for detailed analysis, and experimental data were collected for comparison.

[0189] The experimental data are from recently measured hydrological data, and the time range is from 2023 to 2024. The data cover indicators such as rainfall, flow changes in different periods, and hydrological simulation errors after sub-basin division.

[0190] Table 1 Comparison table of sub-basin division effects

[0191]

[0192] The experimental results show that the present invention is superior to the traditional methods in terms of the sub - basin division accuracy, calculation efficiency, and hydrological simulation stability. In terms of the optimization of the number of sub - basins, the average number of sub - basins divided by the method of the present invention is 187, which is 18.5% less than that of the traditional method. This indicates that through the adaptive dynamic threshold optimization, redundant sub - basins can be effectively reduced, making the division results more compact and avoiding the over - division problem caused by the unreasonable setting of the fixed threshold in the traditional method. At the same time, although the GIS semi - automatic division method has been optimized, it still cannot achieve the effect of the present invention, indicating the limitations of the traditional method in complex terrains.

[0193] In terms of the optimization of the sub - basin boundary, the number of boundary mutation points of the present invention is reduced to 735, which is 41.7% lower than that of the traditional method and 32.1% lower than that of the GIS semi - automatic method. This shows that through the Wasserstein distance optimization and the analysis of the basin morphological characteristics, the present invention improves the continuity of the sub - basin boundary, reduces the boundary mutation caused by the calculation error of the flow direction, and optimizes the coherence between sub - basins. The reduction of the boundary mutation points means a smoother division, ensuring the stability of the hydrological simulation calculation and avoiding the error accumulation caused by the unreasonable boundary in the traditional method.

[0194] In terms of the calculation efficiency, the method of the present invention only takes 46 minutes, which is 51.6% shorter than that of the traditional method (95 minutes) and 37% higher than that of the GIS semi - automatic method (73 minutes). This result shows that through the optimized calculation of the flow direction, the adjustment of the information entropy, and the optimization of the dynamic neural network, the present invention significantly reduces the calculation redundancy, improves the processing ability of large - scale basin data, and at the same time increases the degree of calculation automation, avoiding the disadvantages of the need for manual intervention in the GIS semi - automatic method.

[0195] In terms of the hydrological simulation accuracy, the hydrological simulation error of the present invention is only 8.5%, which is 47.5% lower than that of the traditional method (16.2%) and 34.1% lower than that of the GIS semi - automatic method (12.9%). This result shows that through the dynamic optimization of the sub - basin threshold, the optimization of the flow - direction calculation, and the optimization of the boundary division by the neural network, the present invention effectively reduces the basin - division error, ensures the accuracy of the hydrological simulation calculation, and improves the matching degree of the hydrological characteristics of the sub - basins. The lower hydrological simulation error means that the division results are more in line with the actual hydrological process, providing reliable support for water - resource management, flood prediction, and ecological assessment.

[0196] Generally speaking, the present invention has achieved significant improvements in terms of the sub - basin division accuracy, boundary smoothness, calculation efficiency, and hydrological simulation stability. Compared with the traditional method, the present invention can not only adaptively adjust the sub - basin threshold according to the basin morphological characteristics, but also optimize the flow - direction calculation, improve the calculation efficiency, ensure the rationality of the sub - basin division, and provide more accurate and reliable technical support for basin hydrological calculation.

[0197] As described above, it is only the preferred specific implementation manner of the present invention, but the protection scope of the present invention is not limited thereto. Any person skilled in the art within the technical scope disclosed by the present invention, according to the technical solution and inventive concept of the present invention, makes equivalent substitutions or changes, and all should be covered within the protection scope of the present invention.

Claims

1. A method for dividing small and medium-sized watersheds into sub-basins for distributed models, characterized in that: The steps include: S1. Obtain digital elevation model data and remote sensing image data of the study area, perform depression filling processing on the digital elevation model data, extract preliminary river network information, and identify the main stream area; S2. Extract the actual water system based on remote sensing image data, overlay and analyze it with the preliminary river network information, identify the mainstream grid deviation area, and calculate the elevation difference between the mainstream grid of the river network and the surrounding depression grid; S3. Based on the calculated elevation difference, the local elevation adjustment of the digital elevation model data is performed, and the depression filling process, flow direction calculation and cumulative flow calculation are re-executed to optimize the river network extraction results; S4. Based on the optimized river network extraction results, the basin morphological parameters are calculated, and the basins are classified based on the manifold learning method to identify different basin types; S5. Construct a sub-basin optimal threshold calculation model based on basin morphological parameters, calculate the catchment area threshold, and dynamically adjust the sub-basin division threshold; S6. Based on the catchment area threshold, the adjusted digital elevation model data is used for flow direction calculation and confluence analysis, sub-basins are divided, the boundaries are optimized in combination with the basin morphological parameters, and small-area sub-basins are merged; S7. Import the sub-basin division results into the distributed hydrological model, analyze the impact of sub-basin division on basin runoff simulation, and adjust the sub-basin division parameters according to the hydrological simulation results.

2. The method for dividing small and medium-sized watersheds into sub-basins for distributed models according to claim 1 is characterized in that: The S2 specifically includes: S21. Extract the actual water system based on remote sensing image data, use the water body recognition model based on the fusion of graph attention network and watershed hydrodynamic characteristics to detect water bodies, and calculate the water body classification confidence S by combining spectral characteristics, terrain characteristics and hydrodynamic simulation information. w : Among them, σ is the normalized activation function, is the adjacency weight calculated based on the attention mechanism, is the water body characteristic of node j at the kth scale, ω is the hydrodynamic information weight, and F(i,j) is the hydrodynamic characteristic; S22, perform spatial interpolation processing on the extracted preliminary river network information, and use the high-order interpolation method with adaptive weights to calculate the elevation value Z of the river network mainstream grid r (i, j), stored as interpolated river network grid data R(i, j): Among them, Z(m,n) is the known elevation value of the surrounding grid, μ is the adjustment coefficient, is the second-order Laplace smoothing term, w m,n is the adaptive weight, M and N are the grid numbers, and Z(i, j) is the elevation value of the digital elevation model data at the grid coordinate (i, j); S23. Binarize the water body classification confidence to generate water body distribution grid data W(i, j), and perform grid overlay operation with the interpolated river network grid data R(i, j). Use cross entropy error to calculate consistency error: Where E is the consistency error value, if E>T e , it is determined to be the mainstream grid deviation area, T e is the error threshold, ln is the logarithmic function; S24, calculate the elevation difference of the identified mainstream grid deviation area to solve the elevation difference ΔZ max : ΔZ max =max(|Z r (i,j)-Z w (i,j)|)+γ·σ Z ; Among them, σ Z is the elevation standard deviation of the neighboring area, γ is the adaptive correction coefficient, max is the maximum value operation, Z r (i, j) is the grid elevation value of the main stream of the river network, Z w (i, j) is the elevation value of the surrounding depression grid; S25. Store the calculated mainstream grid deviation area and elevation difference data.

3. The method for dividing small and medium-sized watersheds into sub-basins for distributed models according to claim 1 is characterized in that: The S3 specifically includes: S31, based on the calculated local elevation difference ΔZ max , adjust the elevation of the deviation area grid and build a local elevation adjustment model: Z′ r (i,j)=Z r (i,j)-α·ΔZ max ; Among them, Z′ r (i, j) is the adjusted grid elevation value of the main stream of the river network, α is the local elevation adjustment coefficient, Z r (i, j) is the grid elevation value of the main stream of the river network before adjustment; S32. Use the Markov random field model to smooth the elevation adjustment area and build a local elevation optimization model: in, is the optimized river network main stream grid elevation value, argmin Z′ To solve the candidate elevation adjustment value Z′ so that the objective function value is minimized, Z′ r (m,n) represents the adjusted elevation value of the neighborhood grid (m,n), λ is the regularization parameter, H(Z′) is the regularization term of the elevation change, N(i,j) is the neighborhood set of the grid, Z′ is the candidate elevation adjustment value, m and n are the row and column indices of the neighborhood grid, and β is the neighborhood smoothing factor; S33, based on optimized digital elevation model data Perform a fill operation to fill all internal depressions and calculate the digital elevation model data after elevation adjustment: Among them, Z f (i, j) is the elevation value of the digital elevation model data after filling, Z w (i, j) is the elevation value of the adjacent depression grid; S34, based on the generated digital elevation model data Z after filling f (i, j), the flow matrix D(i, j) is calculated using the direction-weighted flow calculation method: Among them, Z f (m,n) is the elevation value of the neighborhood grid (m,n), d m,n is the Euclidean distance between grid (i, j) and its neighboring grid (m, n), ε is the flow adjustment index, argmax (m,n)∈N(i,j) To determine which adjacent grid (m,n) the water flows from the current grid (i,j); Among them, the adaptive confluence weight model is used to calculate the catchment area of ​​the river network grid: Where A(i,j) is the catchment area of ​​grid (i,j), exp is the exponential function, ∈ is the catchment weight adjustment parameter, A(m,n) is the catchment area of ​​the neighboring grid (m,n), and d p,q is the Euclidean distance between grid (i, j) and its neighboring grid (p, q); S35. Optimize river network extraction results based on the calculated flow direction matrix and cumulative flow: in, Extract raster data for the optimized river network, T A The minimum catchment area threshold for river network extraction.

4. The method for dividing small and medium-sized watersheds into sub-basins for distributed models according to claim 1, characterized in that: The S4 specifically includes: S41. Extract raster data based on the optimized river network and calculate the basin morphological parameters, including curvature C, fractal dimension D, confluence path complexity P and topological centrality T. c , construct a watershed morphological feature set: Among them, L r is the actual length of the main stream in the basin, L s is the shortest straight-line distance to the main stream of the basin, N(θ) is the minimum number of grids required to cover the basin boundary at scale θ, θ is the grid scale, Q is the total number of confluence paths, and L i is the length of the i-th confluence path in the sub-basin, ζ is the hydrodynamic adjustment factor, Q max and Q min is the maximum and minimum flow in the sub-basin, N(i) is the set of sub-basins connected to sub-basin i, and W ij is the confluence weight between sub-basin i and sub-basin j, W ik is the confluence weight between sub-basin i and sub-basin K; S42. Based on the calculated basin morphological parameters, construct the basin morphological characteristic matrix M f : Among them, C n ,D n ,P n ,T cn The tortuosity, fractal dimension, confluence path complexity and topological centrality of the nth sub-basin, where N is the total number of sub-basins; S43. Based on the watershed morphological feature matrix, a graph neural network based on an adaptive attention mechanism is used to classify watersheds and construct a morphological classification model: in, is the feature vector of basin i in the l+1th layer network, is the feature vector of basin j in the l-th layer network, σ is the nonlinear activation function, is the neighbor weight calculated based on the attention mechanism, exp is the exponential function, τ is the weight vector in the attention mechanism, W (l) is the transformation matrix of the l-th layer network, LeakyReLU is the leakage rectified linear unit activation function, is the feature vector of basin i in the l-th layer network, is the feature vector of basin k in the l-th layer network; S44. Based on the calculated watershed morphology classification results, a watershed classification label matrix based on soft clustering and neighborhood optimization is constructed, and the final classification label is calculated using a joint optimization method. Among them, d(T n ,M f ) is the morphological characteristics and category T of sub-basin n n The Euclidean distance of the centroid, K is the total number of classification categories, is the temperature adjustment parameter, d(T j ,M f ) is the matrix M of sub-basin j and morphological characteristics f Medium Category T j The Euclidean distance of the centroid, d(T k ,M f ) is the sub-basin k and the morphological feature matrix M f Medium Category T k Euclidean distance to the centroid, ω n,j is the hydrological connectivity weight, T n is the classification label, μ is the neighborhood influence factor, From the classification label T n Choose a category that gives you the highest overall score; S45. Based on the obtained optimized watershed classification labels, a visual watershed classification map is generated, and the classification results are stored.

5. The method for dividing small and medium-sized watersheds into sub-basins for distributed models according to claim 1, characterized in that: The S5 specifically includes: S51, based on the stored classification results and watershed morphological feature matrix M f , construct a multi-scale feature fusion model, introduce watershed morphology, terrain gradient, and hydrodynamic characteristics, and calculate the sub-basin catchment area threshold Among them, F k is the comprehensive morphological characteristics of sub-basin k, β k To optimize the weight coefficient, N(n) is the set of neighboring sub-basins of sub-basin n, ω n,j is the neighborhood influence weight, F j is the comprehensive morphological characteristics of sub-basin j, α o is the neighborhood weighted impact factor; S52. Based on the calculated sub-basin catchment area threshold, the dynamic threshold for sub-basin division is calculated in combination with the sub-basin information entropy. Among them, η and α e is the dynamic adjustment coefficient, E A is the sub-basin information entropy, is the average catchment area threshold, is the standard deviation of the catchment area threshold; S53. Optimize sub-basin division based on calculated dynamic thresholds: Among them, ω n,j is the neighborhood clustering weight, N s is the number of sub-basin splits, T A,min is the minimum area threshold of the sub-basin, T A,max is the maximum area threshold of the sub-basin, is the catchment area threshold of the neighborhood sub-basin j; S54, based on the calculated optimized sub-basin data, store the final sub-basin division result and output the sub-basin division matrix M sub : Among them, X N and Y N is the geographical coordinate of sub-basin n, is the final catchment area threshold of sub-basin n; S55. Store the calculated sub-basin division data to support the sub-basin division calculation in step S6.

6. The method for dividing small and medium-sized watersheds into sub-basins for distributed models according to claim 1, characterized in that: The S6 specifically includes: S61. Based on the catchment area threshold, calculate the multi-scale hydrological topological flow matrix D(i,j) and construct the watershed topological network: D(i,j)=argmax (m,n)∈N(i,j) (Z f (i,j)-Z f (m,n))+γ1·Φ(i,j); Among them, Z f (i,j) is the elevation value of grid (i,j), Z f (m,n) is the elevation value of the neighborhood grid (m,n) after filling, γ1 is the hydrological topology optimization factor, Φ(i,j) is the water flow clustering topology optimization function, argmax (m,n)∈N(i,j) To select the most likely direction of water flow among all candidate flow directions in the neighborhood, N(i,j) is the neighborhood of grid (i,j); S62. Based on the calculated flow direction matrix, the cumulative confluence area calculation method based on information entropy optimization is adopted: A(i,j)=∑ (m,n)∈N(i,j) A(m,n)+W(i,j)·(1+γ2·E A ); Among them, A(i,j) is the catchment area of ​​grid (i,j), W(i,j) is the watershed contribution weight, γ2 is the information entropy adjustment factor, E A is the information entropy of the sub-basin, A(m,n) is the cumulative catchment area of ​​the upstream neighboring grid (m,n); S63, based on the calculated catchment area, combined with the calculated dynamic catchment area threshold Using dynamic neural network to optimize sub-basin division: Among them, M sub (i, j) is the sub-basin partition result matrix, σ is the adaptive classification activation function, ω n,j is the neighborhood influence weight, γ3 is the neighborhood smoothing adjustment parameter; S64. Based on the calculated sub-basin division results, optimize the sub-basin boundaries and adjust the merging relationship of small-area sub-basins: Among them, B(i,j) is the optimized sub-basin boundary matrix, δ(M sub (m,n),M sub (i,j)) is the sub-basin merging function, T A,min is the minimum area threshold of the sub-basin; S65. Based on the calculated sub-watershed division results, store the final optimized sub-watershed data.

Images