A slope unit division method combining computer graphics and hydrological principles

Abstract

The application discloses a slope unit division method combining computer graphics and hydrological principles, and is based on digital elevation model data, binary images are constructed by setting a curvature threshold and a standardization threshold; morphological skeleton processing and pretreatment are carried out, and a skeleton graph with consistent topological results is determined; the skeleton graph with consistent topological structures is analyzed based on iterative flow accumulation, and topological sorting results are determined; a valley skeleton line is taken as a catchment boundary, a ridge skeleton line is taken as a watershed boundary, and the topological sorting results are connected into coherent initial watershed units; finally, subunits in the initial watershed units are subjected to cluster analysis through topographic features, hydrological parameters and mechanical parameters, and slope units are obtained; the application combines computer graphics and hydrological principles, guarantees hydrological consistency of the division results, avoids the problem that hydrological boundaries and topographic features do not match, and improves slope unit division efficiency.

Description

Technical Field

[0001] This invention relates to the technical field of geological monitoring, and more specifically, to a method for dividing slope units by combining computer graphics and hydrological principles. Background Technology

[0002] Current methods for slope unit delineation based on Digital Elevation Models (DEMs) mainly include hydrological analysis methods and methods based on unit homogeneity. Hydrological analysis methods construct closed catchment basins through watershed extraction techniques, such as using a two-way DEM processing framework to integrate flow direction and discharge parameters through forward and reverse watershed analysis, and delineate slope units based on the Hydrological Process Analysis (HPAM) method. Furthermore, curvature analysis can be used to extract concave and convex topographic features, leading to watershed segmentation methods based on topographic curvature. While these methods are simple to operate, they still have limitations such as reliance on manual technical correction and significant heterogeneity in mechanical parameters within units. In contrast, methods based on unit homogeneity aim to meet the computational needs of landslide stability analysis. Slope units should possess spatially consistent mechanical properties. By introducing area and azimuth control parameters, iterative subdivision methods are used to delineate secondary watersheds obtained from cumulative discharge analysis until user-defined threshold conditions are met. By developing a hybrid method combining morphological segmentation and azimuth-slope homogeneity principles, slope units are redefined as three-dimensional, coherent, homogeneous, and closed regions. This method automates the algorithm while retaining the flexibility of manual threshold adjustment. Although these homogeneity-oriented methods improve parameter uniformity, they neglect the regulatory role of hydrological processes in landslide formation, leading to systematic differences between the assessment model and the actual hydraulic-mechanical coupling mechanism. Summary of the Invention

[0003] To address at least one of the aforementioned technical problems, the present invention aims to provide a slope unit division method that combines computer graphics and hydrological principles, ensuring the hydrological consistency of the division results and avoiding the mismatch between hydrological boundaries and topographic features.

[0004] This invention provides a method for slope unit division that combines computer graphics and hydrological principles, including: Step 1: Obtain digital elevation model data of the target area, extract the average curvature raster image from the digital elevation model data, and determine the binary image based on the preset curvature threshold and standardization threshold; Step 2: Perform morphological skeletonization on the binary image to extract the topographic skeleton lines representing ridgelines and valleys. After preprocessing, a skeleton map with consistent topological structure is obtained. Step 3: Through iterative flow accumulation analysis, determine the topological sorting results of the nodes and edges in the corresponding skeleton graph with consistent topological structure; use the valley skeleton line as the water catchment boundary and the ridge skeleton line as the water divide boundary, and connect them into a coherent initial watershed unit according to the topological sorting results. Step 4: Obtain topographic features, hydrological parameters, and mechanical parameters; and perform cluster analysis on the sub-units in the initial watershed unit based on the topographic features, hydrological parameters, and mechanical parameters to obtain the slope unit.

[0005] In this scheme, the step of obtaining the average curvature raster image from the digital elevation model data specifically includes: The digital elevation model data is divided according to a preset grid size to obtain a raster image; Construct a grid centered on one of the raster images. The local neighborhood window, where m is an odd number greater than or equal to 5; A quadratic surface equation is fitted using the least squares method based on all elevation values ​​within the neighborhood window. Based on the quadratic surface equation, the average curvature of the raster image at the corresponding center is determined.

[0006] In this scheme, the step of determining the binary image based on a preset curvature threshold and a normalization threshold specifically includes: The average curvature raster image is segmented according to a preset curvature threshold, and the DEM is divided into ridge area raster images and valley area raster images. Centered on the raster image of the ridge area or the raster image of the valley area, extract its The average curvature distribution of the raster image within the neighborhood is calculated, and the average curvature deviation of the raster image within the neighborhood is also calculated. If the average curvature deviation exceeds a preset deviation threshold, the classification result of the central grid is corrected according to the neighborhood majority principle to obtain a corrected ridge area grid image or valley area grid image. Based on a standardized threshold, the ridge area raster is assigned a logical value of 1, and the valley area raster is assigned a logical value of 0, resulting in a binary image.

[0007] In this scheme, the step of performing morphological skeletonization processing on the binary image to extract the topographic skeleton lines representing ridgelines and valley lines specifically includes: Based on a preset distance transformation thinning algorithm, Euclidean distance transformation is performed on the foreground grid in the binary image, and the Euclidean distance value from each foreground grid to the nearest background grid is calculated to generate a distance map; the foreground grid includes ridge area grids and valley area grids, and the background grids are non-ridge area grids and non-valley area grids. Based on the distance map, local maxima points are extracted from the distance map as skeleton candidate points; the local maxima point is a grid whose distance value is greater than that of all neighboring points in a 3×3 neighborhood; The candidate skeleton points are connected to form an initial skeleton line using the 8-neighborhood connectivity criterion, which states that a connection is established when two candidate skeleton points are adjacent in the horizontal, vertical, or diagonal directions. After the initial skeleton lines are generated, boundary pixels that do not meet the preset conditions are deleted to obtain terrain skeleton lines with a width of one pixel. The terrain skeleton lines with a width of one pixel include ridge skeleton lines and valley skeleton lines.

[0008] In this scheme, the step of obtaining a skeleton diagram with a consistent topological structure after preprocessing is specifically as follows: Connect the first and last endpoints of the skeleton line to form a straight line segment, and calculate the distance from the vertex on the skeleton line to the straight line segment; If the distance value is less than the preset initial distance threshold, the vertex corresponding to the distance value is deleted until all vertices on the skeleton line are traversed or the remaining number of vertices on the skeleton line is the preset number, thus obtaining the initial simplified skeleton line. Extract the feature values ​​of vertices on the initially simplified skeleton lines. The feature values ​​of the vertices include at least the local curvature change of the corresponding vertex, the length of the branch where the vertex is located, and the distance between the vertex and the intersection point of the adjacent skeleton lines. After normalizing the feature values ​​of the vertices, a weighted summation is performed to obtain the importance score of the corresponding vertex. If the importance score of the vertex is less than the preset importance threshold, the corresponding vertex is deleted. After traversing the feature values ​​of all vertices, the simplified skeleton line is obtained. Then, a topology consistency check is performed on the simplified skeleton lines to obtain a skeleton diagram with a consistent topology.

[0009] In this scheme, the step of determining the topological sorting of nodes and edges in the corresponding skeleton graph by iterative flow accumulation analysis based on the skeleton graph with consistent topological structure is as follows: Extract grid cells from a skeleton diagram with a consistent topology; Using the grid cell as the center, calculate the elevation difference between the corresponding grid cell and its 8 neighboring grid cells; The direction of the area with the largest elevation difference is selected as the flow direction of the grid cell; Based on the flow direction assignment of all grid cells, starting from the highest grid cell in the target area, the number of grid cells flowing into each grid cell is accumulated in the order from upstream to downstream, so as to obtain the cumulative flow value of each grid cell. Based on the cumulative flow value of each grid cell, determine the topological sorting of nodes and edges in the corresponding grid cell; After traversing all grid cells, determine the topological sorting result of the nodes and edges in the corresponding skeleton graph.

[0010] In this scheme, after obtaining the cumulative flow value of each grid cell, the method further includes: The elevation value of the central grid cell is compared with the elevation values ​​of its eight neighboring grid cells. If the elevation values ​​of all eight neighboring grid cells are higher than the elevation value of the corresponding central grid cell, then the corresponding central grid cell is set as a low-lying grid cell. After raising the elevation of the low-lying grid cell to the elevation value of its lowest neighboring grid cell, the flow direction calculation and flow accumulation calculation are re-executed, and the iteration calculation is recorded. When the total number of iterations exceeds the preset number of iterations or all low-lying grid cells are eliminated, the final cumulative flow value of each grid cell is output, and the topological sorting of nodes and edges in the grid cell is determined accordingly.

[0011] In this scheme, the step of performing cluster analysis on the sub-units of the initial watershed unit based on topographic features, hydrological parameters, and mechanical parameters to obtain the slope unit specifically includes: Each sub-unit in the initial watershed unit is treated as an independent initial cluster; By using a bottom-up aggregation sequence, independent initial clusters are compared and analyzed based on terrain features to determine similarity values; If the similarity value is greater than the preset similarity threshold, the two initial clusters will be merged to obtain a merged subunit. Based on mechanical parameter thresholds and hydrological parameter thresholds, the coefficient of variation of mechanical parameters and the continuity parameters of surface runoff in the merged sub-units are verified respectively. If they meet the requirements, the merged sub-units are correct; if they do not meet the requirements, the merged sub-units are split back into two initial clusters. After traversing all initial clusters, multiple merged sub-units are obtained, which are then merged to form the corresponding slope unit.

[0012] One or more technical solutions proposed in this invention have at least the following technical effects: 1. This invention achieves fully automated slope unit division, significantly reducing human intervention. In step 2, the morphological skeletonization and preprocessing techniques automatically extract and optimize the terrain skeleton lines. In step 3, the iterative flow accumulation analysis and topology sorting algorithm automatically construct the initial watershed units. The entire process, from DEM data input to slope unit output, does not require manual setting of key division parameters, overcoming the shortcomings of existing methods that heavily rely on repeated trial and error and manual parameter adjustment, and greatly improving the division efficiency and objectivity. 2. This invention achieves a deep integration of hydrological principles and computer graphics, ensuring the hydrological consistency of the classification results. In step 1, the invention accurately identifies ridge and valley areas through average curvature calculation; in step 2, it extracts topographic skeleton lines representing ridges and valleys through morphological skeletonization; in step 3, it explicitly uses the valley skeleton lines as catchment boundaries and the ridge skeleton lines as watershed boundaries, connecting the skeleton network into coherent initial watershed units using the D8 unidirectional flow algorithm and topological sorting algorithm, forming closed catchment areas. This design ensures that each slope unit has a complete surface runoff path, avoiding the problem of mismatch between hydrological boundaries and topographic features in traditional methods. In summary, this invention, by combining computer graphics and hydrological principles, ensures the hydrological consistency of the division results, avoids the mismatch between hydrological boundaries and topographic features, and improves the efficiency of slope unit division. Attached Figure Description

[0013] The present invention will now be described in further detail with reference to the accompanying drawings and specific implementation methods.

[0014] Figure 1 A flowchart of a slope unit division method combining computer graphics and hydrological principles according to the present invention is shown; Figure 2 A schematic diagram of the slope element division of the present invention is shown. Detailed Implementation

[0015] To further illustrate the technical means and effects of the present invention in achieving its intended purpose, the following detailed description of the specific implementation methods, structures, features, and effects of the present invention, in conjunction with the accompanying drawings and preferred embodiments, is provided below.

[0016] Figure 1 The flowchart illustrates a slope unit division method that combines computer graphics and hydrological principles according to the present invention.

[0017] like Figure 1 As shown, this invention discloses a slope unit division method combining computer graphics and hydrological principles, comprising: Step 1: Obtain digital elevation model data of the target area, extract the average curvature raster image from the digital elevation model data, and determine the binary image based on the preset curvature threshold and standardization threshold; Step 2: Perform morphological skeletonization on the binary image to extract the topographic skeleton lines representing ridgelines and valleys. After preprocessing, a skeleton map with consistent topological structure is obtained. Step 3: Through iterative flow accumulation analysis, determine the topological sorting results of the nodes and edges in the corresponding skeleton graph with consistent topological structure; use the valley skeleton line as the water catchment boundary and the ridge skeleton line as the water divide boundary, and connect them into a coherent initial watershed unit according to the topological sorting results. Step 4: Obtain topographic features, hydrological parameters, and mechanical parameters; and based on these features, perform cluster analysis on the sub-units of the initial watershed unit to obtain slope units, such as... Figure 2 .

[0018] According to an embodiment of the present invention, before taking the average curvature raster image from the digital elevation model data, the DEM data is filtered by a preset spatial filtering and noise reduction method to eliminate noise interference and local outliers in the DEM data; the average curvature calculation method in digital image processing technology is used to calculate the average curvature value of the terrain surface grid by grid to obtain the average curvature raster image.

[0019] It should be noted that the spatial filtering and noise reduction process employs an adaptive median filtering algorithm. This algorithm dynamically adjusts the size of the filtering window based on the grayscale value distribution characteristics of each raster point and its neighborhood window in the DEM data. Larger windows are used in flat areas to enhance noise suppression, while smaller windows are used in areas with abrupt terrain changes (including steep slopes and cliffs) to preserve terrain details. For example, the initial size of the filtering window is a 3×3 raster cell, with a maximum size not exceeding 11×11 raster cells. The output value of the adaptive median filtering algorithm is the median of all valid elevation values ​​within the filtering window. For raster cells containing holes or outliers, these holes or outliers are not included in the median calculation; instead, they are filled by inverse distance weighted interpolation of the filtered neighborhood valid values. The filtering process iterates until the root mean square change in all raster elevation values ​​between two consecutive iterations is less than a preset convergence threshold (set as one-tenth of the DEM's vertical resolution) to ensure the stability and convergence of the filtering results.

[0020] According to an embodiment of the present invention, the step of obtaining the average curvature raster image from the digital elevation model data specifically includes: The digital elevation model data is divided according to a preset grid size to obtain a raster image; Construct a grid centered on one of the raster images. The local neighborhood window, where m is an odd number greater than or equal to 5; A quadratic surface equation is fitted using the least squares method based on all elevation values ​​within the neighborhood window. Based on the quadratic surface equation, the average curvature of the raster image at the corresponding center is determined.

[0021] It should be noted that, for example, setting the equation of the quadratic surface as Where z is the elevation value, x and y are plane coordinates, and a, b, c, d, e, and f are the surface coefficients to be fitted; then, the average curvature H of the raster image at the corresponding center is determined using the coefficients in the quadratic surface equation, and its formula is: .

[0022] According to an embodiment of the present invention, the step of determining the binary image based on a preset curvature threshold and a normalization threshold specifically includes: The average curvature raster image is segmented according to a preset curvature threshold, and the DEM is divided into ridge area raster images and valley area raster images. Centered on the raster image of the ridge area or the raster image of the valley area, extract its The average curvature distribution of the raster image within the neighborhood is calculated, and the average curvature deviation of the raster image within the neighborhood is also calculated. If the average curvature deviation exceeds a preset deviation threshold, the classification result of the central grid is corrected according to the neighborhood majority principle to obtain a corrected ridge area grid image or valley area grid image. Based on a standardized threshold, the ridge area raster is assigned a logical value of 1, and the valley area raster is assigned a logical value of 0, resulting in a binary image.

[0023] It should be noted that, for example, if the preset curvature threshold is between -0.05 and 0.05, when the average curvature of the average curvature raster image is 0.06, which is greater than the maximum value in the preset curvature threshold, the corresponding raster image is set as a ridge area raster image; if the average curvature is -0.08, which is less than the minimum value in the preset curvature threshold, the corresponding raster image is set as a valley area raster image; the preset deviation threshold is 1.5 times the neighborhood standard deviation; for transition raster in flat areas (curvature close to zero), its periodic 8 neighboring raster images are examined, and the minority follows the majority to forcibly classify them, preventing large areas of "unassigned" regions from appearing.

[0024] It should be noted that for grids in flat terrain areas whose absolute curvature is lower than a preset flatness threshold, the flatness threshold is set to 0.3 times the global standard deviation of the average curvature. These grids are classified as transition zone grids. In the binary image, the transition zone grids are assigned a value based on their spatial adjacency with adjacent ridge or valley areas. Specifically, if the number of ridge grids in the 8-neighborhood of the transition zone grid is greater than the number of valley grids, then the value is assigned a logical value of 1; otherwise, the value is assigned a logical value of 0.

[0025] According to an embodiment of the present invention, the step of performing morphological skeletonization processing on the binary image to extract topographic skeleton lines representing ridgelines and valley lines specifically includes: Based on a preset distance transformation thinning algorithm, Euclidean distance transformation is performed on the foreground grid in the binary image, and the Euclidean distance value from each foreground grid to the nearest background grid is calculated to generate a distance map; the foreground grid includes ridge area grids and valley area grids, and the background grids are non-ridge area grids and non-valley area grids. Based on the distance map, local maxima points are extracted from the distance map as skeleton candidate points; the local maxima point is a grid whose distance value is greater than that of all neighboring points in a 3×3 neighborhood; The candidate skeleton points are connected to form an initial skeleton line using the 8-neighborhood connectivity criterion, which states that a connection is established when two candidate skeleton points are adjacent in the horizontal, vertical, or diagonal directions. After the initial skeleton lines are generated, boundary pixels that do not meet the preset conditions are deleted to obtain terrain skeleton lines with a width of one pixel. The terrain skeleton lines with a width of one pixel include ridge skeleton lines and valley skeleton lines.

[0026] It should be noted that the preset conditions include at least the following: the boundary pixel is not an endpoint, the deletion does not destroy the connectivity of the skeleton line, and there are at least two foreground connected components in the neighborhood of the boundary pixel.

[0027] According to an embodiment of the present invention, the step of obtaining a skeleton diagram with a consistent topological structure after preprocessing specifically includes: Connect the first and last endpoints of the skeleton line to form a straight line segment, and calculate the distance from the vertex on the skeleton line to the straight line segment; If the distance value is less than the preset initial distance threshold, the vertex corresponding to the distance value is deleted until all vertices on the skeleton line are traversed or the remaining number of vertices on the skeleton line is the preset number, thus obtaining the initial simplified skeleton line. Extract the feature values ​​of vertices on the initially simplified skeleton lines. The feature values ​​of the vertices include at least the local curvature change of the corresponding vertex, the length of the branch where the vertex is located, and the distance between the vertex and the intersection point of the adjacent skeleton lines. After normalizing the feature values ​​of the vertices, a weighted summation is performed to obtain the importance score of the corresponding vertex. If the importance score of the vertex is less than the preset importance threshold, the corresponding vertex is deleted. After traversing the feature values ​​of all vertices, the simplified skeleton line is obtained. Then, a topology consistency check is performed on the simplified skeleton lines to obtain a skeleton diagram with a consistent topology.

[0028] It should be noted that the skeleton line has many jagged redundant vertices. First, small bends deviating from the straight line are filtered out using a preset initial distance threshold. Then, important branch points and endpoints are selected based on importance scores. Finally, a topology consistency check is performed. The vertices include branch nodes, merging nodes, and skeleton line endpoints. The preset initial distance threshold ranges from 1 to 3 times the DEM resolution; for example, the preset number is 200, and the preset importance threshold is 60%. The topology consistency check verifies whether the simplified skeleton line maintains the same branch structure and connectivity as the original skeleton line. If a topology inconsistency is found, the process reverts to the previous level of simplification and adjusts the simplification parameters before re-performing the simplification operation.

[0029] According to an embodiment of the present invention, the step of determining the topological sorting result of nodes and edges in the corresponding skeleton graph by iterative flow accumulation analysis of the skeleton graph with consistent topological structure specifically includes: Extract grid cells from a skeleton diagram with a consistent topology; Using the grid cell as the center, calculate the elevation difference between the corresponding grid cell and its 8 neighboring grid cells; The direction of the area with the largest elevation difference is selected as the flow direction of the grid cell; Based on the flow direction assignment of all grid cells, starting from the highest grid cell in the target area, the number of grid cells flowing into each grid cell is accumulated in the order from upstream to downstream, so as to obtain the cumulative flow value of each grid cell. Based on the cumulative flow value of each grid cell, determine the topological sorting of nodes and edges in the corresponding grid cell; After traversing all grid cells, determine the topological sorting result of the nodes and edges in the corresponding skeleton graph.

[0030] It should be noted that, looking at the surrounding eight directions from each grid, water flows towards the lowest neighboring grid, and the amount of upstream water received by each grid is counted to determine the cumulative flow value of each grid unit; the larger the cumulative flow value of each grid unit, the later it is in the topology sort.

[0031] According to an embodiment of the present invention, after obtaining the cumulative flow value of each grid cell, the method further includes: The elevation value of the central grid cell is compared with the elevation values ​​of its eight neighboring grid cells. If the elevation values ​​of all eight neighboring grid cells are higher than the elevation value of the corresponding central grid cell, then the corresponding central grid cell is set as a low-lying grid cell. After raising the elevation of the low-lying grid cell to the elevation value of its lowest neighboring grid cell, the flow direction calculation and flow accumulation calculation are re-executed, and the iteration calculation is recorded. When the total number of iterations exceeds the preset number of iterations or all low-lying grid cells are eliminated, the final cumulative flow value of each grid cell is output, and the topological sorting of nodes and edges in the grid cell is determined accordingly.

[0032] It should be noted that the elevation values ​​of these low-lying grids are raised to their lowest neighboring elevation values ​​to eliminate the blocking effect of closed depressions in the terrain on the water flow path; after the depressions are filled, the flow direction calculation and flow accumulation calculation are re-executed, and iterative execution is performed until all depressions are eliminated or the preset number of iterations is reached, for example, the number of iterations is set to half of the sum of the number of rows and columns of the DEM.

[0033] According to an embodiment of the present invention, the step of performing cluster analysis on the sub-units in the initial watershed unit based on topographic features, hydrological parameters, and mechanical parameters to obtain slope units specifically includes: Each sub-unit in the initial watershed unit is treated as an independent initial cluster; By using a bottom-up aggregation sequence, independent initial clusters are compared and analyzed based on terrain features to determine similarity values; If the similarity value is greater than the preset similarity threshold, the two initial clusters will be merged to obtain a merged subunit. Based on mechanical parameter thresholds and hydrological parameter thresholds, the coefficient of variation of mechanical parameters and the continuity parameters of surface runoff in the merged sub-units are verified respectively. If they meet the requirements, the merged sub-units are correct; if they do not meet the requirements, the merged sub-units are split back into two initial clusters. After traversing all initial clusters, multiple merged sub-units are obtained, which are then merged to form the corresponding slope unit.

[0034] It should be noted that if the topographic skeleton lines in the initial watershed unit are not connected, then the two corresponding unconnected topographic skeleton lines are set as different sub-units; the topographic features include the average slope, average aspect, elevation variance, etc. of the corresponding topographic skeleton lines; by normalizing the topographic features, the similarity value of each topographic feature parameter is determined, and then a weighted sum is performed to determine the similarity value of the two corresponding independent initial clusters.

[0035] It should be noted that the mechanical parameters include the internal friction angle, cohesion, and natural density. The coefficient of variation for each mechanical parameter is calculated through normalization and weighting. If the coefficient of variation is less than or equal to a set threshold (mechanical parameter threshold), the mechanical parameters meet the requirements. A lower coefficient of variation indicates higher mechanical consistency. The hydrological connectivity requirement is quantified by calculating the flow direction consistency index at the boundary of the merged sub-units. This flow direction consistency index is calculated based on the cosine of the D8 flow direction angle between the grids on both sides of the boundary, determining the hydrological connectivity angle. If the hydrological connectivity angle is less than a preset angle threshold (hydrological parameter threshold), the corresponding hydrological connectivity meets the requirements. A smaller angle indicates better hydrological connectivity. When both the mechanical parameters and hydrological connectivity meet the requirements, the merged sub-units are correct.

[0036] This invention discloses a slope unit partitioning method combining computer graphics and hydrological principles. Based on digital elevation model data, a binary image is constructed by setting curvature and standardization thresholds. Morphological skeleton processing and preprocessing are then performed to determine a skeleton map with consistent topological results. The skeleton map with consistent topological structure is analyzed based on iterative flow accumulation to determine the topological sorting result. Valley skeleton lines are used as catchment boundaries, and ridge skeleton lines as watershed boundaries. These are connected according to the topological sorting result to form coherent initial watershed units. Finally, the sub-units in the initial watershed units are clustered using topographic features, hydrological parameters, and mechanical parameters to obtain slope units. This invention, by combining computer graphics and hydrological principles, ensures the hydrological consistency of the partitioning results, avoids the problem of mismatch between hydrological boundaries and topographic features, and improves the efficiency of slope unit partitioning.

[0037] The above description is merely a preferred embodiment of the present invention and is not intended to limit the present invention in any way. Although the present invention has been disclosed above with reference to preferred embodiments, it is not intended to limit the present invention. Any person skilled in the art can make some modifications or alterations to the above-disclosed technical content to create equivalent embodiments without departing from the scope of the present invention. Any simple modifications, equivalent changes and alterations made to the above embodiments based on the technical essence of the present invention without departing from the scope of the present invention shall still fall within the scope of the present invention.

Claims

1. A method for dividing slope units by combining computer graphics and hydrological principles, characterized in that, include: Step 1: Obtain digital elevation model data of the target area, extract the average curvature raster image from the digital elevation model data, and determine the binary image based on the preset curvature threshold and standardization threshold; Step 2: Perform morphological skeletonization on the binary image to extract the topographic skeleton lines representing ridgelines and valleys. After preprocessing, a skeleton map with consistent topological structure is obtained. Step 3: Through iterative flow accumulation analysis, determine the topological sorting results of the nodes and edges in the corresponding skeleton graph with consistent topological structure; use the valley skeleton line as the water catchment boundary and the ridge skeleton line as the water divide boundary, and connect them into a coherent initial watershed unit according to the topological sorting results. Step 4: Obtain topographic features, hydrological parameters, and mechanical parameters; and perform cluster analysis on the sub-units in the initial watershed unit based on the topographic features, hydrological parameters, and mechanical parameters to obtain the slope unit; The step of performing cluster analysis on the sub-units of the initial watershed unit based on topographic features, hydrological parameters, and mechanical parameters to obtain the slope unit specifically includes: Each sub-unit in the initial watershed unit is treated as an independent initial cluster; By using a bottom-up aggregation sequence, independent initial clusters are compared and analyzed based on terrain features to determine similarity values; If the similarity value is greater than the preset similarity threshold, the two initial clusters will be merged to obtain a merged subunit. Based on mechanical parameter thresholds and hydrological parameter thresholds, the coefficient of variation of mechanical parameters and the continuity parameters of surface runoff in the merged sub-units are verified respectively. If they meet the requirements, the merged sub-units are correct; if they do not meet the requirements, the merged sub-units are split back into two initial clusters. After traversing all initial clusters, multiple merged sub-units are obtained, which are then merged to form the corresponding slope unit.

2. The slope unit division method combining computer graphics and hydrological principles according to claim 1, characterized in that, The steps for obtaining the average curvature raster image from the digital elevation model data are as follows: The digital elevation model data is divided according to a preset grid size to obtain a raster image; Construct a grid centered on one of the raster images. The local neighborhood window, where m is an odd number greater than or equal to 5; A quadratic surface equation is fitted using the least squares method based on all elevation values ​​within the neighborhood window. Based on the quadratic surface equation, the average curvature of the raster image at the corresponding center is determined.

3. The slope unit division method combining computer graphics and hydrological principles according to claim 2, characterized in that, The step of determining the binary image based on the preset curvature threshold and standardization threshold is as follows: The average curvature raster image is segmented according to a preset curvature threshold, and the DEM is divided into ridge area raster images and valley area raster images. Centered on the raster image of the ridge area or the raster image of the valley area, extract its The average curvature distribution of the raster image within the neighborhood is calculated, and the average curvature deviation of the raster image within the neighborhood is also calculated. If the average curvature deviation exceeds a preset deviation threshold, the classification result of the central grid is corrected according to the neighborhood majority principle to obtain a corrected ridge area grid image or valley area grid image. Based on a standardized threshold, the ridge area raster is assigned a logical value of 1, and the valley area raster is assigned a logical value of 0, resulting in a binary image.

4. The slope unit division method combining computer graphics and hydrological principles according to claim 3, characterized in that, The step of performing morphological skeletonization processing on the binary image to extract the topographic skeleton lines representing ridgelines and valley lines is as follows: Based on a preset distance transformation thinning algorithm, Euclidean distance transformation is performed on the foreground grid in the binary image, and the Euclidean distance value from each foreground grid to the nearest background grid is calculated to generate a distance map; the foreground grid includes ridge area grids and valley area grids, and the background grids are non-ridge area grids and non-valley area grids. Based on the distance map, local maxima points are extracted from the distance map as skeleton candidate points; the local maxima point is a grid whose distance value is greater than that of all neighboring points in a 3×3 neighborhood; The candidate skeleton points are connected to form an initial skeleton line using the 8-neighborhood connectivity criterion, which states that a connection is established when two candidate skeleton points are adjacent in the horizontal, vertical, or diagonal directions. After the initial skeleton lines are generated, boundary pixels that do not meet the preset conditions are deleted to obtain terrain skeleton lines with a width of one pixel. The terrain skeleton lines with a width of one pixel include ridge skeleton lines and valley skeleton lines.

5. The slope unit division method combining computer graphics and hydrological principles according to claim 4, characterized in that, The step of obtaining a skeleton diagram with a consistent topological structure after preprocessing is as follows: Connect the first and last endpoints of the skeleton line to form a straight line segment, and calculate the distance from the vertex on the skeleton line to the straight line segment; If the distance value is less than the preset initial distance threshold, the vertex corresponding to the distance value is deleted until all vertices on the skeleton line are traversed or the remaining number of vertices on the skeleton line is the preset number, thus obtaining the initial simplified skeleton line. Extract the feature values ​​of vertices on the initially simplified skeleton lines. The feature values ​​of the vertices include at least the local curvature change of the corresponding vertex, the length of the branch where the vertex is located, and the distance between the vertex and the intersection point of the adjacent skeleton lines. After normalizing the feature values ​​of the vertices, a weighted summation is performed to obtain the importance score of the corresponding vertex. If the importance score of the vertex is less than the preset importance threshold, the corresponding vertex is deleted. After traversing the feature values ​​of all vertices, the simplified skeleton line is obtained. Then, a topology consistency check is performed on the simplified skeleton lines to obtain a skeleton diagram with a consistent topology.

6. The slope unit division method combining computer graphics and hydrological principles according to claim 5, characterized in that, The step of determining the topological sorting of nodes and edges in the skeleton graph with consistent topological structure through iterative flow accumulation analysis is as follows: Extract grid cells from a skeleton diagram with a consistent topology; Using the grid cell as the center, calculate the elevation difference between the corresponding grid cell and its 8 neighboring grid cells; The direction of the area with the largest elevation difference is selected as the flow direction of the grid cell; Based on the flow direction assignment of all grid cells, starting from the highest grid cell in the target area, the number of grid cells flowing into each grid cell is accumulated in the order from upstream to downstream, so as to obtain the cumulative flow value of each grid cell. Based on the cumulative flow value of each grid cell, determine the topological sorting of nodes and edges in the corresponding grid cell; After traversing all grid cells, determine the topological sorting result of the nodes and edges in the corresponding skeleton graph.

7. The slope unit division method combining computer graphics and hydrological principles according to claim 6, characterized in that, After obtaining the cumulative flow value for each grid cell, the process further includes: The elevation value of the central grid cell is compared with the elevation values ​​of its eight neighboring grid cells. If the elevation values ​​of all eight neighboring grid cells are higher than the elevation value of the corresponding central grid cell, then the corresponding central grid cell is set as a low-lying grid cell. After raising the elevation of the low-lying grid cell to the elevation value of its lowest neighboring grid cell, the flow direction calculation and flow accumulation calculation are re-executed, and the iteration calculation is recorded. When the total number of iterations exceeds the preset number of iterations or all low-lying grid cells are eliminated, the final cumulative flow value of each grid cell is output, and the topological sorting of nodes and edges in the grid cell is determined accordingly.

Images