A method for modifying DEM data by considering surface catchment route features

Abstract

This invention discloses a DEM data correction method that takes into account the characteristics of surface runoff routes, belonging to the field of digital elevation model construction technology. It solves the problem that in some scenarios, region-based DEMs can no longer accurately simulate the surface runoff process, thus preventing DEMs from playing their due role in urban surface runoff simulation. The key technical points are that the method includes the following steps: S1, data preparation: acquiring DEM data of the analysis area and marking the outlet and boundary positions in the DEM raster; S2, marking the flow direction of the raster cell at the outlet location as 0, indicating that the water flow cell of that raster does not flow to other raster cells; S3, finding the raster cells adjacent to the outlet cell and calculating the flow direction of each raster cell according to the principle of maximizing elevation difference. The effect is that this invention covers the key technology of raster elevation value correction based on surface flow direction characteristics, which can improve the application value of DEM data in the field of hydrological analysis.

Description

Technical Field

[0001] This invention relates to the field of digital elevation model (DEM) construction technology, and more specifically, to a DEM data correction method that takes into account the characteristics of surface water runoff routes. Background Technology

[0002] Digital Elevation Models (DEMs) are a core component of national basic geographic data. Currently, my country has preliminarily completed the construction of DEM data at scales of 1:1,000,000, 1:250,000, 1:50,000, and 1:10,000 for some regions. These data have played an indispensable role in national economic development, national defense, and scientific research, with applications spanning surveying and mapping, transportation, military, water conservancy, agriculture, environment, resource management, planning, and tourism, among many other fields.

[0003] DEM-based runoff simulation can extract a region's water system and watersheds, providing fundamental data for hydrological analysis at different scales, making it a significant application area for DEMs. The basic process of DEM-based surface runoff simulation involves first filling depressions in the DEM data to eliminate local depressions, and then determining the direction of water flow on each grid cell based on the principle of maximizing elevation difference, thus simulating the runoff path on a region's surface. It should be noted that this process works well for simulating natural terrain with certain undulations, effectively extracting surface valleys and watersheds and classifying watersheds at different levels. However, in urban areas, various engineering projects transform the surface into a complex of continuous and abrupt terrain features, creating many typical artificial landforms such as roads and plazas. In this context, DEMs based on these areas can no longer accurately simulate surface runoff processes, thus limiting their effectiveness in urban surface runoff simulation.

[0004] The main reason for the aforementioned problems is that urban terrain contains a large amount of flat land and slightly undulating terrain. Additionally, some surface features have water-impermeable functions, such as man-made structures, while others have drainage functions, such as underground pipe outlets along roads. This unique terrain and environment dictate that the direction and path of surface water flow cannot be determined solely by elevation differences. In other words, traditional DEMs, which rely on elevation accuracy as the basic evaluation criterion, cannot accurately simulate surface water runoff processes. To achieve accurate simulation of surface water runoff processes, it is necessary to consider the specific environment of urban areas, identify which surface features are water-impermeable and which are drainage features, and thus determine where surface water in a given area will flow.

[0005] Based on the above analysis, to ensure that special terrains such as cities can accurately simulate surface runoff processes, this invention proposes a DEM data correction method that takes into account the characteristics of surface runoff routes, while fully understanding the influence of different land features on the runoff process. This method first identifies the impermeable features and drainage areas in the analysis area, and then uses these impermeable features and drainage areas as constraints to correct the DEM elevation values ​​within the area, enabling the corrected DEM to accurately simulate surface runoff processes. The data sources required for this method are readily available, and the process is easy to implement, demonstrating significant application potential for enhancing the use of DEM in urban surface runoff simulation and urban stormwater simulation research. Summary of the Invention

[0006] To address the aforementioned problems, this invention proposes a method for correcting DEM elevation values ​​in special terrains. The aim is to ensure that the corrected DEM data can accurately simulate surface water runoff processes. This method primarily involves correcting the elevation values ​​of the original DEM data under constraints such as impermeable features and drainage areas, effectively solving the technical challenge of DEMs failing to effectively simulate surface water runoff processes in special areas such as urban areas.

[0007] The technical solution adopted in this invention is as follows:

[0008] A DEM data correction method that takes into account the characteristics of surface water runoff routes includes the following steps:

[0009] S1, Data preparation: Obtain DEM data of the analysis area and mark the outlet location and boundary location in the DEM grid;

[0010] S2, mark the flow direction of the grid cell at the outlet position as 0, indicating that the water flow cell of this grid does not flow to other grids;

[0011] S3, find the grid cells adjacent to the grid cells where the outlet is located, calculate the flow direction of each grid cell according to the principle of maximum elevation difference, and correct the elevation value of the relevant grid cell cells based on the final flow direction to the outlet, and mark its flow direction;

[0012] S4, perform elevation value correction and flow direction marking on the remaining grid cells;

[0013] S5. After processing, the calculation results will be output in raster data format.

[0014] In step S1, DEM data of the analysis area is obtained. Specifically, high-resolution surface image data can be collected by drones, and then an oblique photogrammetry model is constructed based on Context Capture software to generate DEM data.

[0015] Furthermore, the specific process of step S3 is as follows:

[0016] S31, determine the adjacent grid cells of the outlet;

[0017] S32. Based on the principle of maximizing the elevation difference over horizontal distance, calculate the flow direction of all adjacent grids at the outlet. If the flow direction of a grid points to the outlet, mark the flow direction of that grid; otherwise, do not mark the flow direction of that grid.

[0018] S33, Based on the results of S32, if there are still unmarked grids adjacent to the outlet, calculate the minimum value that the grid needs to be adjusted to point to the outlet or to a grid that has been marked with the flow direction. Correct the grid elevation value according to this minimum value, modify the DEM grid elevation value, and mark the flow direction.

[0019] S34, continuously execute S33 until all grid cells have been marked with the flow direction, at which point the elevation values ​​of adjacent grid cells at the outlet are adjusted.

[0020] Furthermore, the specific process of step S4 is as follows:

[0021] S41, Traverse the DEM raster cells of the processing area and record all raster cells that have not yet been marked with flow direction and are adjacent to raster cells that have been marked with flow direction;

[0022] S42, for all the grates recorded in S41, calculate the flow direction of each grates according to the principle of the maximum elevation difference in horizontal distance. If the flow direction of a certain grates points to a grates whose flow direction has been marked, then mark the flow direction of that grates as well.

[0023] S43, for all the rasters recorded in S41, determine whether each raster has been marked with the flow direction. If not, calculate the minimum value that the raster needs to be adjusted to point to the raster with the marked flow direction, correct the raster elevation value according to this minimum value, and mark the flow direction.

[0024] S44, Repeat step S43 until all rasters recorded in S41 have been marked with flow direction;

[0025] S45, Repeat steps S41-S44 until all grids in the entire processing area have completed flow direction marking.

[0026] A DEM data correction device that takes into account the characteristics of surface water runoff routes includes:

[0027] The acquisition module is used to acquire DEM data of the analysis area and mark the outlet location and boundary location in the DEM raster;

[0028] The first marking module is used to mark the flow direction of the grid cell at the outlet position as 0, indicating that the water flow cell of that grid does not flow to other grids;

[0029] The second marking module is used to find the grid cells adjacent to the grid cell where the outlet is located, calculate the flow direction of each grid cell according to the principle of maximum elevation difference, and correct the elevation value of the relevant grid cell cells based on the final flow direction to the outlet, and mark its flow direction.

[0030] The third marking module is used to correct the elevation value and mark the flow direction of the remaining grid cells;

[0031] The output module is used to output the processed calculation results in raster data format.

[0032] A computer-readable storage medium having a computer program stored thereon, the computer program executing the above-described method when run by a computer.

[0033] Technical features and beneficial effects of the present invention:

[0034] (1) The DEM data correction method proposed in this invention takes into account the characteristics of surface water runoff routes. Its main purpose is to ensure that the DEM of special surface areas can correctly simulate the surface water runoff process and enhance the application potential of DEM in the fields of surface runoff simulation and urban stormwater analysis.

[0035] (2) The method of the present invention is ingeniously designed and easy to implement. The corrected DEM does not have a significant change in elevation accuracy compared with the original DEM data, and the corrected DEM data can also guarantee its application value in other fields. Attached Figure Description

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

[0037] Figure 2 This is a schematic diagram of the flow direction markings in the initial stage processing area of ​​the method of the present invention.

[0038] Figure 3 This is a schematic diagram of the drainage outlet and its adjacent grid in the method of the present invention.

[0039] Figure 4 This is a schematic diagram showing the marked flow direction grid and its adjacent grid markings in the method of the present invention.

[0040] Figure 5 This is a schematic diagram of the original DEM data generated from images collected by a drone in the experimental area according to an embodiment of the present invention.

[0041] Figure 6 This is an image of the experimental area collected by a drone in an embodiment of the present invention.

[0042] Figure 7 This is a schematic diagram of the corrected DEM of the experimental area in an embodiment of the present invention.

[0043] Figure 8 This is a schematic diagram of the water catchment route extracted from the original road DEM in the experimental area in an embodiment of the present invention.

[0044] Figure 9 This is a schematic diagram of the water catchment routes extracted from the modified road DEM in the experimental area in an embodiment of the present invention.

[0045] Figure 10 This is a block diagram of the modules in an embodiment of the present invention. Detailed Implementation

[0046] The principles and features of the present invention are described below with reference to the accompanying drawings. The examples given are only for explaining the present invention and are not intended to limit the scope of the present invention.

[0047] like Figure 1 As shown, the specific steps of the DEM data correction method that takes into account the characteristics of surface water runoff routes according to the present invention are as follows:

[0048] Step 1: Read in the DEM data of the study area, as well as the boundary data and drainage outlet location data of the area to be processed in the study area;

[0049] Step 2: Generate a flow direction marker matrix and overlay it with the DEM raster cells of the study area. Initialize the flow direction matrix to null. Then, assign -1 to the matrix elements corresponding to the raster on the boundary to be processed, indicating that the flow direction is not calculated. Assign 0 to the matrix elements corresponding to the raster where the drain outlet is located, indicating that the raster at the drain outlet does not flow to other raster cells. Figure 2 );

[0050] Step 3: In the DEM data, find all the grates adjacent to the drain outlet and label them sequentially as G1, G2, ..., G n The grid corresponding to the drain outlet is sequentially labeled as P1, P2, ..., P n If n = 1, it means that the drain outlet corresponds to only one grid cell. Figure 3 );

[0051] Step 4: For grid G i The flow direction is calculated using the D8 algorithm. If the flow direction points to the grid corresponding to the drain outlet, then grid G ​​is set to [value missing]. i In the flow direction marking matrix, the corresponding element is assigned a value of 1, indicating that the flow direction of the raster has been marked; otherwise, no processing is performed for the time being.

[0052] Step 5: For grids G that have not yet been marked with flow direction i The flow direction is calculated using the D8 algorithm. If the flow direction points to a grid G ​​that has already been marked with a flow direction... j Then the grid G i The corresponding element in the flow direction marking matrix is ​​assigned a value of 1, indicating that the flow direction of the raster has been marked; otherwise, no processing is performed for the time being.

[0053] Step 6: Repeat step 5 until there are no new raster marker flows, then proceed to step 7;

[0054] Step 7: For grids G that have not yet been marked with flow direction i Determine the elevation relationship between the grid and its surrounding grids. If the elevation of the grid is less than or equal to the elevation of the surrounding grids, proceed to step 8; otherwise, proceed to step 9.

[0055] Step 8: At this point, the grid G ​​needs to be corrected. i To determine the elevation value, first calculate the elevation difference between the grid and the drainage outlet grid, or between its adjacent grids that have been marked with the flow direction. The calculation formula is as follows:

[0056]

[0057] In the formula, H(G) i ) represents the raster G i Elevation value, H(O) j ) indicates the grid at the drain outlet or the grid G i Elevation values ​​d of adjacent graticles whose flow direction has been marked j The horizontal distance between two grid cells is calculated using the row and column numbers of the two cells and the grid size, as shown in the following formula:

[0058]

[0059] In the formula, r is the DEM grid size, and G im G in They are grids G i row and column numbers, O jm O jn respectively grid O j The row and column numbers.

[0060] In the calculated Δ j If the value with the smallest value is selected, and its corresponding raster is marked as Q, then raster G is selected. i The elevation difference at the horizontal distance from grid Q is marked as Δ. Q The horizontal distance between two grid cells is marked as d. Q Then grid G i The formula for calculating the elevation value that needs to be corrected is:

[0061] Δ=Δ Q d Q +0.001

[0062] Grid G i After the elevation value is corrected, the raster G will be simultaneously... iThe flow direction marker is changed to marked. If there are still grid cells G at this point... i If no flow direction value is marked, return to step 5; otherwise, continue to step 10.

[0063] Step 9: At this point, it is explained that the grid G... i The flow direction points to grid Q, which has not yet been marked with a flow direction. Therefore, the elevation value of grid Q needs to be modified. First, calculate the elevation value of grid G. i The elevation difference at the horizontal distance from grid Q is calculated using the following formula:

[0064]

[0065] In the formula, H(G) i ) represents the raster G i The elevation value, H(Q) represents the elevation value of raster Q, d Q The horizontal distance between two grid cells can be calculated using their row and column numbers and their grid size. The formulas for calculating these values ​​were already disclosed in step 8, and the calculation method here is consistent with those disclosed above.

[0066] Then, calculate the grid G ​​separately. i The elevation difference in horizontal distance between the grid at the drainage outlet and the grid with the flow direction already marked is calculated using the following formula:

[0067]

[0068] In the formula, H(G) i ) represents the raster G i Elevation value, H(O) j ) indicates the grid at the drain outlet or the grid G i Elevation values ​​d of adjacent graticles whose flow direction has been marked j The horizontal distance between two grid cells can be calculated using their row and column numbers and their grid size. The formulas for calculating these values ​​were already disclosed in step 8, and the calculation method here is consistent with those disclosed above.

[0069] In the calculated Δ j If the value with the largest value is selected, and its corresponding raster is marked as R, then raster G is selected. i The elevation difference at the horizontal distance from the grid R is marked as Δ. R The horizontal distance between two grid cells is marked as d. R The formula for calculating the elevation value that needs to be corrected for raster Q is:

[0070] Δ=Δ Q d Q -Δ R d R +0.001

[0071] After the elevation values ​​of raster Q are corrected, raster G will be... i The flow direction marker is changed to marked. If there are still grid cells G at this point... i If no flow direction value is marked, return to step 5; otherwise, continue to step 10.

[0072] Step 10: In the DEM to be processed, find the cells adjacent to the cells with marked flow directions in the unmarked cells, and label them G1, G2, ..., G... n ( Figure 6 ).

[0073] Step 11: Randomly select an unmarked raster G. i Determine the elevation relationship between the grid and its surrounding grids. If the elevation of the grid is less than or equal to the elevation of the surrounding grids, proceed to step 12; otherwise, proceed to step 13.

[0074] Step 12: At this point, the grid G ​​needs to be corrected. i To determine the elevation value, first calculate the elevation difference between the grid cell and its adjacent grid cells that have already been marked with flow direction. The calculation formula is as follows:

[0075]

[0076] In the formula, H(G) i ) represents the raster G i Elevation value, H(O) j ) indicates the relationship with the grid G i Elevation values ​​d of adjacent graticles whose flow direction has been marked j The horizontal distance between two grid cells can be calculated using their row and column numbers and their grid size. The formulas for calculating these values ​​were already disclosed in step 8, and the calculation method here is consistent with those disclosed above.

[0077] In the calculated Δ j If the value with the smallest value is selected, and its corresponding raster is marked as Q, then raster G is selected. i The elevation difference at the horizontal distance from grid Q is marked as Δ. Q The horizontal distance between two grid cells is marked as d. Q The formula for calculating the elevation value that the raster Gi needs to correct is:

[0078] Δ=Δ Q d Q +0.001

[0079] Grid G i After the elevation value is corrected, the raster G will be simultaneously... iThe flow direction marker is changed to marked. If there are still grid cells G at this point... i If no flow direction value is marked, return to step 11; otherwise, continue to step 14.

[0080] Step 13: At this point, it is explained that the grid G... i The flow direction points to grid Q, which has not yet been marked with a flow direction. Therefore, the elevation value of grid Q needs to be modified. First, calculate the elevation value of grid G. i The elevation difference at the horizontal distance from grid Q is calculated using the following formula:

[0081]

[0082] In the formula, H(G) i ) represents the raster G i The elevation value, H(Q) represents the elevation value of raster Q, d Q The horizontal distance between two grid cells can be calculated using their row and column numbers and their grid size. The formulas for calculating these values ​​were already disclosed in step 8, and the calculation method here is consistent with those disclosed above.

[0083] Then, calculate the grid G ​​separately. i The elevation difference in horizontal distance between adjacent and marked flow direction grids is calculated using the following formula:

[0084]

[0085] In the formula, H(G) i ) represents the raster G i Elevation value, H(O) j ) indicates the relationship with the grid G i Elevation values ​​d of adjacent graticles whose flow direction has been marked j The horizontal distance between two grid cells can be calculated using their row and column numbers and their grid size. The formulas for calculating these values ​​were already disclosed in step 8, and the calculation method here is consistent with those disclosed above.

[0086] In the calculated Δ j If the value with the largest value is selected, and its corresponding raster is marked as R, then raster G is selected. i The elevation difference at the horizontal distance from the grid R is marked as Δ. R The horizontal distance between two grid cells is marked as d. R The formula for calculating the elevation value that needs to be corrected for raster Q is:

[0087] Δ=Δ Q d Q -Δ R d R +0.001

[0088] After the elevation values ​​of raster Q are corrected, raster G will be... i The flow direction marker is changed to marked. If there are still grid cells G at this point... i If no flow direction value is marked, return to step 11; otherwise, continue to step 14.

[0089] Step 14: Repeat steps 10-13 until all rasters in the DEM to be processed have been marked with flow direction, then the DEM correction is complete.

[0090] Example 1

[0091] To verify the application effect of the method of the present invention, an experimental area was selected on a section of road within the campus of Chuzhou University in Anhui Province. The elevation of this section of road is slightly higher at both ends and slightly lower in the middle. There are two underground drainage outlets on one side of the road. During rainfall, the surface water of this section of road is discharged through the two drainage outlets and underground pipes, making it suitable as an experimental area to verify the usability of the method of the present invention.

[0092] Step 1: Obtain DEM data for the example area ( Figure 7 The data includes boundary location information and drainage outlet location information. The original DEM data was generated based on UAV oblique photogrammetry data, while the location information of the example area boundary and drainage outlet was obtained from UAV imagery. Figure 8 );

[0093] Step 2: Using the method of this invention, the DEM of the example area is corrected, and the result is as follows. Figure 9 As shown;

[0094] Step 3: Extract the catchment routes from both the original and corrected DEM data. The results are as follows: Figure 8 and Figure 9 As shown in the figure, the white lines represent the water flow routes.

[0095] Comparing the water runoff routes extracted from the original DEM data and the corrected DEM data, it can be found that the distribution of water runoff routes extracted from the original DEM data has no obvious pattern, and most of the water runoff routes do not reach the drainage outlet. In contrast, the water runoff routes extracted from the corrected DEM data converge towards the drainage outlet more effectively, indicating that the corrected DEM can accurately simulate the surface water runoff process in this area.

Claims

1. A method for correcting DEM data that takes into account the characteristics of surface water runoff routes, characterized in that, Includes the following steps: S1, Data preparation: Obtain DEM data of the analysis area and mark the outlet location and boundary location in the DEM grid; S2, mark the flow direction of the grid cell at the outlet position as 0, indicating that the water flow cell of this grid does not flow to other grids; S3, find the grid cells adjacent to the grid cells where the outlet is located, calculate the flow direction of each grid cell according to the principle of maximum elevation difference, and correct the elevation value of the relevant grid cell cells based on the final flow direction to the outlet, and mark its flow direction; S4, perform elevation value correction and flow direction marking on the remaining grid cells; S5. After processing, output the calculation results in raster data format; The specific process of step S4 is as follows: S41, Traverse the DEM raster cells of the processing area and record all raster cells that have not yet been marked with flow direction and are adjacent to raster cells that have been marked with flow direction; S42, for all the grates recorded in S41, calculate the flow direction of each grates according to the principle of the maximum elevation difference in horizontal distance. If the flow direction of a certain grates points to a grates whose flow direction has been marked, then mark the flow direction of that grates as well. S43, for all the rasters recorded in S41, determine whether each raster has been marked with the flow direction. If not, calculate the minimum value that the raster needs to be adjusted to point to the raster with the marked flow direction, correct the raster elevation value according to this minimum value, and mark the flow direction. S44, Repeat step S43 until all rasters recorded in S41 have been marked with flow direction; S45, Repeat steps S41-S44 until all grids in the entire processing area have completed flow direction marking.

2. The DEM data correction method considering surface water runoff characteristics according to claim 1, characterized in that, The specific process of step S3 is as follows: S31, determine the adjacent grid cells of the outlet; S32. Based on the principle of maximizing the elevation difference over horizontal distance, calculate the flow direction of all adjacent grids at the outlet. If the flow direction of a grid points to the outlet, mark the flow direction of that grid; otherwise, do not mark the flow direction of that grid. S33, Based on the results of S32, if there are still unmarked grids adjacent to the outlet, calculate the minimum value that the grid needs to be adjusted to point to the outlet or to a grid that has been marked with the flow direction. Correct the grid elevation value according to this minimum value, modify the DEM grid elevation value, and mark the flow direction. S34, continuously execute S33 until all grid cells have been marked with the flow direction, at which point the elevation values ​​of adjacent grid cells at the outlet are adjusted.

3. A DEM data correction device that takes into account the characteristics of surface water runoff routes, characterized in that, include: The acquisition module is used to acquire DEM data of the analysis area and mark the outlet location and boundary location in the DEM raster; The first marking module is used to mark the flow direction of the grid cell at the outlet position as 0, indicating that the water flow cell of that grid does not flow to other grids; The second marking module is used to find the grid cells adjacent to the grid cell where the outlet is located, calculate the flow direction of each grid cell according to the principle of maximum elevation difference, and correct the elevation value of the relevant grid cell cells based on the final flow direction to the outlet, and mark its flow direction. The third marking module is used to traverse the DEM raster cells of the processing area, record all raster cells that have not yet been marked with flow direction and are adjacent to raster cells with marked flow direction; then, for all recorded raster cells, the flow direction of each raster is calculated based on the principle of maximizing the elevation difference in horizontal distance. If a raster's flow direction points to a raster cell with marked flow direction, then that raster is also marked with flow direction; if no raster is marked, the minimum adjustment required for that raster to point to a raster cell with marked flow direction is calculated, the raster's elevation value is corrected according to this minimum value, and the flow direction is marked; the above process is repeated until the flow direction of all recorded raster cells has been marked; this module is repeated until the flow direction of all raster cells in the entire processing area has been marked. The output module is used to output the processed calculation results in raster data format.

4. A computer-readable storage medium, characterized in that, It stores a computer program that, when executed by a computer, performs the method as described in any one of claims 1-2.

Images