A silt dam extraction method based on high-resolution DEM

Abstract

The application provides a silt dam extraction method based on high-resolution DEM, comprising the following steps: step one, data preprocessing; step two, extracting a river network; step three, performing improved LSD detection based on a mountain shadow map; step four, performing angle filtering on the contour straight line obtained in step three based on the river network obtained in step two; and step five, constructing a cross model according to the contour straight line filtered in step four to identify the silt dam. The application extracts the river network based on flow, narrows the research area of silt dam detection, avoids the influence of terraces and roads with the same linear characteristics in the mountainous area, and effectively reduces the misjudgment in the detection process. The application comprehensively considers the morphological characteristics of the silt dam itself and the morphological characteristics of the position where the silt dam is located, and applies the LSD algorithm to silt dam detection, and improves the LSD algorithm to separate contour straight lines with different slopes.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of topography and relates to a silt dam integrated calculation method, in particular to a silt dam extraction method based on high-resolution DEM. BACKGROUND

[0002] The line detection algorithm (LSD) is widely used in linear object detection due to its high efficiency and robustness, avoiding the influence of most nonlinear noise. Calculating the image gradient is an important step in LSD detection of local straight line contours, where the contour is an area where the gray level changes rapidly from dark to light. The gradient is pseudo-sorted, and the pixels with a gradient amplitude greater than the threshold and error threshold are divided into the same rectangular region, and the main axis of the rectangular region is selected as the detected straight line.

[0003] However, the current silt dam detection method is often based on multiple data to calculate the maximum water storage area or sediment area, and then locate the silt dam. These methods pay less attention to the morphological characteristics of silt dams, resulting in low accuracy in complex areas. However, using a series of features including silt dam morphological characteristics may also have a negative impact. For example, terraces in mountainous areas have similar agricultural functions and linear features to silt dams, which interfere with the detection of silt dams. The roads on both sides of the river network also have similar linear features. Previous studies have focused on calculating storage capacity using different methods to indirectly infer the location of silt dams, and the accurate extraction of silt dams has not been effectively solved. SUMMARY

[0004] In view of the deficiencies in the prior art, the purpose of the present application is to provide a silt dam extraction method based on high-resolution DEM to solve the technical problems that the efficiency and accuracy of the existing silt dam extraction method need to be further improved.

[0005] To solve the above technical problems, the technical scheme is adopted as follows:

[0006] A silt dam extraction method based on high-resolution DEM, the method comprising the following steps:

[0007] Step 1: Data preprocessing

[0008] Step 101: Read DEM data.

[0009] Step 102: Smooth preprocessing.

[0010] Step 2: Extract the river network.

[0011] Step 3: Improved LSD detection based on mountain shadow map.

[0012] Based on the shadow map of mountain, the improved LSD algorithm is used to separate the straight lines with different slopes, and the detection results of contour straight lines are obtained.

[0013] The improved LSD algorithm is as follows: assuming that a straight line L has been output before a straight line R is output; starting from the end point F of the straight line R, it is determined whether there are other straight line pixels except the straight line in the eight neighborhoods; if not, the pixel point F is kept; on the contrary, if yes, the least square method is used to calculate the slope and the difference value of the two straight lines; if the difference value is within the range of 0-5°, it is considered that the two straight lines belong to the same straight line, and the pixel point F is kept; then the next pixel point is checked to meet the condition; if the pixel point is encountered, it is discarded, and two straight lines with different slopes are separated; after all the straight lines are determined, only the contour straight lines intersecting the river network are kept.

[0014] Step four: based on the river network obtained in step two, the contour straight lines obtained in step three are angle filtered.

[0015] Step five: according to the contour straight lines filtered by angle in step four, a cross model is constructed to identify the silt dam:

[0016] Step 501: merging the contour straight lines belonging to the same silt dam:

[0017] The contour straight lines belonging to the same silt dam are merged, and the judgment basis for merging is that if the difference in the inclination angle of two straight lines is less than , the two straight lines are merged; the center point of the merged straight line is determined by formula I ; the gentle value SA is set to determine the extension length of the straight line perpendicular to the cross section of the silt dam.

[0018] ;

[0019] In the formula:

[0020] represents the th straight line, ;

[0021] represents the th contour straight line of the silt dam;

[0022] represents the number of contour straight lines;

[0023] C represents the total number of groups of silt dam contour straight lines;

[0024] represents the th straight line, x directional coordinate of the midpoint of the straight line;

[0025] represents the midpoint of the th line; y directional coordinate;

[0026] represents the slope of the th line;

[0027] represents the midpoint of the x directional coordinate;

[0028] represents the midpoint of the y directional coordinate;

[0029] represents the slope of the

[0030] Step 502, construct a cross model to identify the silt dam:

[0031] Step 50201, the center point as the starting point of the vertical silt dam cross section line, stored in the first node of the linked list; and the slope of the vertical silt dam cross section line is calculated by formula II .

[0032] ;

[0033] In the formula:

[0034] represents the slope of the vertical silt dam cross section line;

[0035] represents the slope of the silt dam cross section.

[0036] Step 50202, according to the Bresenham algorithm, determine the coordinate u1 of the next point along the slope k; judge whether the elevation value of u1 and the standard elevation value ch are equal; if they are equal and the flatness value SA=0, stop extending in this direction and execute step 50203; otherwise, store u1 in the linked list, and update the flatness value ch to the elevation value of u1; if the two elevations are equal and the flatness value SA>0, then SA-1; then repeat step 50202 to calculate the next point u2.

[0037] Step 50203, according to the inverse slope -k of the vertical silt dam cross section, repeat steps 50201 to 50202.

[0038] Step 50204, according to the elevation data corresponding to the line perpendicular to the silt dam, determine the highest point coordinate M(x, y).

[0039] Step 50205, according to the dam crest point M extends to both ends along The extension stop condition is formula III, and the number of grids extended to both sides of the cross section of the silt dam is recorded respectively And The cross model is completed.

[0040] ;

[0041] In the formula,

[0042] H represents the threshold of the height difference of the grid. o

[0043] H represents the threshold of the height difference of the grid.

[0044] H represents the threshold of the height difference of the grid.

[0045] Step 50206, according to each grid on the vertical silt dam cross section straight line, respectively, extend And To construct the corresponding rectangle of the silt dam.

[0046] Step 50207, using morphological closing operation to fill the missed grid due to grid-based extension.

[0047] The application also has the following technical features:

[0048] The specific process of step 101 is:

[0049] Step 10101, read 0.15m resolution DEM data.

[0050] Step 10102, resample the 0.15m resolution DEM data into 2m resolution and 10m resolution DEM data.

[0051] The specific process of step 102 is: using a 3 × 3 window median filter to filter the 10m resolution DEM data.

[0052] The specific steps of step two are:

[0053] Step 201, according to the relationship between the inclination and the elevation, the river network is extracted by region growing, wherein the inclination S is defined as:

[0054] ;

[0055] In the formula,

[0056] h x H represents the threshold of the height difference of the grid. x ​the elevation value of the grid

[0057] h y the distance between the grid and the grid

[0058] the distance between the grid and the grid when the grid and the grid are diagonal , otherwise .

[0059] Step 202, set the growth state for each grid, and set the flow grid which is not marked as river network but greater than the flow accumulation threshold r flow as seed point; formula V and formula VI as the growth criterion of region growth, if the growth criterion of formula V or formula VI is met, it is marked as river network pixel, until the position of the demarcation point is determined.

[0060] ;

[0061] In the formula:

[0062] a the center grid

[0063] b the adjacent grid a of the grid

[0064] the elevation value of the grid a ;

[0065] the elevation value of the grid b ;

[0066] S the slope

[0067] E the slope threshold value, E is 0.9

[0068] Step 203, obtain binary river network data with width.

[0069] The specific process of step 203 is:

[0070] Step 20301, establish a queue and traverse the flow data; the first flow data (x0, y0) which is not marked as river network grid element and greater than the flow accumulation threshold r flow is set as growth point, marked as river network grid, and entered into the queue.

[0071] Step 20302, take the head queue grid element, judge whether the 8-neighbor grid of the head queue grid element satisfies formula V or formula VI, if yes, enter the queue.

[0072] Step 20303, check whether the 8 adjacent grids in the head queue grid element are located at the cutoff point and stop growing.

[0073] Step 20304, sequentially execute step 20301 to step 20303 until the queue is empty.

[0074] Step 20305, when the growable grid is not completely marked as a river network, return to step 20301.

[0075] The specific process of step four is as follows:

[0076] The slope angle of the contour straight line obtained in step three is calculated by using the least square method and the slope angle of the river network straight line K obtained in step two is respectively and If the two straight lines satisfy formula VII, it is considered that the contour straight line belongs to the contour straight line of the silt dam; otherwise, the contour straight line is noise and is discarded.

[0077]

[0078] In the formula:

[0079] indicates the slope angle of the contour straight line ;

[0080] indicates the slope angle of the river network straight line K ;

[0081] α indicates the angle threshold.

[0082] Compared with the prior art, the present application has the following technical effects:

[0083] (I) The present application extracts the river network based on flow for region growing, reduces the research area of silt dam detection, avoids the influence of terraces and roads with the same linear characteristics in the mountainous area, and effectively reduces the misjudgment in the detection process.

[0084] (II) The present application comprehensively considers the morphological characteristics of the silt dam itself and the morphological characteristics of the position where the silt dam is located, and applies the LSD algorithm to the silt dam detection, and improves the LSD algorithm to separate contour straight lines with different slopes.

[0085] (III) The application realizes silt dam detection based on high-resolution DEM, and the cross model detection can detect silt dams of different morphological scales. From the final result, the new silt dam detection method is more simple, the overall recognition rate is higher, and the detection result is acceptable. BRIEF DESCRIPTION OF DRAWINGS

[0086] Figure 1 is a process schematic diagram of the silt dam extraction method based on high-resolution DEM.

[0087] Figure 2 is a comparison diagram of median filter noise removal under different resolutions.

[0088] Figure 3 is a three-dimensional schematic diagram of the silt dam.

[0089] Figure 4 is a comparison diagram of LSD line separation effect.

[0090] Figure 5 is a curve diagram of the number of silt dam contour lines, noise lines and lost silt dams under different alpha values.

[0091] Figure 6 is a cross model schematic diagram for identifying silt dams.

[0092] Figure 7 is a comparison diagram of different river network extraction methods.

[0093] Figure 8 is a spatial distribution diagram of silt dam detection results.

[0094] The specific content of the application will be further explained in detail in combination with the embodiments. DETAILED DESCRIPTION

[0095] It should be noted that all the algorithms in the application, if not specially stated, all use the algorithms known in the prior art.

[0096] How to avoid the influence of terraces and roads with similar morphological characteristics, and efficiently and accurately extract silt dams in complex terrain areas by using the characteristics of silt dams is the innovation direction of the application. The application is based on only one kind of data of high-resolution digital elevation model (DEM), and the river network is obtained by region growing calculation, and the noise influence is reduced according to the nearly vertical position relationship between silt dam and river network. The cross model established based on the characteristics of silt dam effectively avoids the interference of terraces and roads in mountains on silt dam detection, and quickly and accurately detects silt dams of different morphological scales, which is helpful to build a complete silt dam database, and plays an important promoting role for silt dam related scientific research and soil and water loss management.

[0097] This invention improves the LSD algorithm to ensure that straight lines are separated from each other. It comprehensively considers the positional relationship between silt-retaining dams and the river network, as well as the topographic features of the silt-retaining dams themselves and the area where they are located, and proposes a cross model for detecting silt-retaining dams of different sizes from DEM data.

[0098] Since each contour line is modeled using a cross-shaped model, and multiple lines represent the same silt-retaining dam, merging contour lines belonging to the same silt-retaining dam avoids redundant construction of the cross-shaped model. Assume there are C groups of silt-retaining dam contour lines, the... Group contains A straight outline. The first The coordinates of the midpoints and the corresponding slopes of the first two lines are given. Equation I is used to obtain the coordinates of the midpoints and the slopes of the merged lines. The coordinates of the center point and the slope of the silt-retaining dam are as follows: .

[0099] .

[0100] The basic assumption of the cross model is that each silt-retaining dam is represented by only one straight line, without considering complex situations such as dam damage. The center point U after merging the straight lines may be located in the upstream or downstream region of the silt-retaining dam. Therefore, the key to constructing the cross model is to determine the top position of the silt-retaining dam, that is, the center point of the cross model.

[0101] Figure 6 This is a schematic diagram of the cross model, based on satisfying... The slope k of the vertical line perpendicular to the cross section of the silt-retaining dam is calculated as follows:

[0102] .

[0103] in, This represents the slope of the silt-retention dam's cross-section. The extension length of the straight line perpendicular to the silt-retention dam's cross-section is determined using the flatness value (SA). Based on the flat terrain of the dam crest and upstream and downstream areas, SA is set only to determine the dam crest by checking if the elevations of two adjacent grid cells are equal. Finally, the Bresenham algorithm is used to extend the straight line, based on the center point U, stopping when SA reaches 0. The highest elevation corresponding to the extended straight line is the dam crest M.

[0104] The final step in constructing the cross model is to start from grid M on the dam crest and move along... It extends in both directions. For example... Figure 3 As shown, because the crest of the silt-retention dam is a gently sloping road, while the two sides of the cross-section are steep mountains, the elevation changes abruptly and is much higher than the crest elevation. This elevation relationship can be used as the cutoff condition for extending the cross-section of the silt-retention dam:

[0105] .

[0106] where grid o is the point of judgment whether to extend or not, H represents the threshold of grid height difference, and the experimental height difference threshold H is set to 1.5 m. The number of grids extending to both sides of the cross section of the silt dam is recorded respectively and Subsequently, for each grid on the vertical silt dam cross section straight line, the corresponding direction of the cross section is extended respectively and to construct the corresponding rectangle of the silt dam. And the morphological closing operation is used to fill the missed grid due to the grid extension.

[0107] Combined with the above cross model, the silt dam can be detected in the grid DEM data. The DEM is to discretize the continuous ground surface into regular square grids, the grid size is the DEM resolution, and the number of rows and columns of the grid is the DEM data volume.

[0108] According to the above technical solution, the specific embodiments of the present application are given below. It should be noted that the present application is not limited to the following specific embodiments, and any equivalent transformation based on the technical solution of the present application falls within the protection scope of the present application.

[0109] Example 1:

[0110] This embodiment gives a kind of silt dam extraction method based on high-resolution DEM, as shown in Figure 1 The method comprises the following steps:

[0111] Step 1, data preprocessing:

[0112] Step 101, read DEM data:

[0113] Step 10101, read 0.15m resolution DEM data.

[0114] Step 10102, resample 0.15m resolution DEM data into 2m resolution and 10m resolution DEM data.

[0115] In this step, too much detailed terrain data greatly increases the calculation complexity and the existence of a large number of noises.

[0116] Step 102, smoothing preprocessing:

[0117] In this step, the noise existing in the DEM data may have a negative impact on the identification. In order to remove many relatively small and isolated targets, and retain the river network and silt dam terrain features.

[0118] A 3x3 window median filter is used to filter the 10m resolution DEM data.

[0119] In this step, the filtering effect under different resolutions is shown in Figure 2 10m resolution is sufficient to meet the requirements of river network extraction.

[0120] Step two, region growing to extract river network:

[0121] In this step, first fill the 10m resolution DEM, and then calculate the flow direction and flow rate using the steepest slope, and set the threshold as the seed point of region growing. The region growing is performed according to the characteristics of the river network with high sides and low middle, and after median filtering and filling, the overall terrain of the river network is relatively flat.

[0122] Step 201, region growing to extract river network according to the relationship between inclination and elevation, wherein the inclination S is defined as:

[0123] ;

[0124] In the formula:

[0125] h x represents the elevation value of the grid x ;

[0126] h y represents the elevation value of the grid in the 8-neighborhood of grid x;

[0127] represents the distance between grid and grid , when grid and grid are in a diagonal relationship , otherwise .

[0128] Step 202, set the growth state for each grid, and take the flow grid which is not marked as river network but greater than the flow accumulation threshold r flow as the seed point; formula V and formula VI are used as the growth criteria for region growing, if the growth criteria of formula V or formula VI are met, it is marked as river network pixel, until the position of the demarcation point is determined.

[0129] ;

[0130] In the formula:

[0131] a represents the center grid;

[0132] b represents the adjacent grid of grid a ;

[0133] Represents grid a Elevation value;

[0134] Represents grid b Elevation value;

[0135] S Indicates the degree of tilt;

[0136] E This represents the tilt threshold, with E set to 0.9.

[0137] In this step, setting the region growth state avoids repeated growth and improves algorithm efficiency.

[0138] Step 203: Obtain binary river network data with width:

[0139] Step 20301: Create a queue and iterate through the flow data. Elements not marked as river network raster elements and exceeding the cumulative flow threshold r are considered. flow The first flow data (x0, y0) is used as a growth point, marked as a river network grid, and entered into the queue.

[0140] Step 20302: Take the head grid element and determine whether its 8 neighboring grids satisfy equation V or equation VI. If they do, add it to the queue.

[0141] Step 20303: Check if the 8 adjacent cells in the head grid element are at the cutoff point and stop growth.

[0142] Step 20304: Execute steps 20301 to 20303 sequentially until the queue is empty.

[0143] Step 20305: If the growable grid is not fully marked as a river network, return to step 20301.

[0144] In this step, morphological closing operations are used to fill in and extract small, undeveloped areas within the river channel. At a 10m resolution, the topographic features of silt-retaining dams cannot be effectively detected, so the river channel detection results are upsampled to high-resolution river channel data. This data is then nested with a high-resolution DEM to locate the silt-retaining dams, narrowing the study area to the river channel and avoiding impact on terraced fields located on the hillside. In this study, binary river channel data is upsampled to a 2-meter resolution, with each grid cell at the 10-meter resolution divided into 5×5 smaller grids. A thinning algorithm is then used to generate a single-grid-sized river channel centerline, corresponding to the river's flow direction.

[0145] Step 3: Improve LSD detection based on mountain shadow map:

[0146] In this embodiment, asFigure 3 Two features of check dam are summarized: check dam location feature and check dam itself feature. The check dam location feature is that the check dam is built on the river channel with steep mountains on both sides of the cross section and is perpendicular to the river channel. The check dam itself feature is that the check dam is a linear engineering structure, the check dam top is a flat road, and the check dam upstream is usually used to build a flat farmland. The check dam is detected according to the five terrain features.

[0147] The LSD algorithm is applied to the mountain shadow map with a resolution of 2 meters. As shown in Figure 4 , since the check dam is built across the river, the check dam contour straight line and the river bank straight line can be adjacent. The improved LSD algorithm is used in this patent to separate straight lines with different slopes to ensure accurate calculation of the slope of the check dam contour straight line.

[0148] For the mountain shadow map, the improved LSD algorithm is used to separate straight lines with different slopes to obtain the contour straight line detection result.

[0149] The improved LSD algorithm is as follows: Figure 4 Since the detected straight lines are output in order, it is assumed that straight line L has been output before straight line R is ready to be output. Starting from the end point F of straight line R, it is determined whether there are other straight line pixel points in the eight neighborhoods except the straight line. If not, the pixel point F is kept. On the contrary, if there are, the slope of the two straight lines and their difference are calculated using the least squares method. If the difference is within the range of 0-5°, the two straight lines are considered to belong to the same straight line, and the pixel point F is kept. Then the next pixel point is checked to meet the condition. If the pixel point is encountered, it is discarded, and two straight lines with different slopes are separated. After all the straight lines are determined, only the contour straight lines intersecting the river network are kept.

[0150] Step four, based on the river network obtained in step two, the contour straight lines obtained in step three are angle filtered:

[0151] In order to reduce erosion as much as possible, the check dam is built approximately perpendicular to the river. Therefore, the river network and the contour straight line detection result are used for angle filtering to keep the straight lines intersecting the river and almost perpendicular.

[0152] The slope angle of the contour straight line obtained in step three and the straight line K of the river network obtained in step two is calculated using the least squares method, which is and respectively. If the two straight lines satisfy formula VII, the contour straight line is considered to belong to the check dam contour straight line; otherwise, the contour straight line is considered to be noise and discarded.

[0153]

[0154] In the formula:

[0155] represents the slope angle of the contour line ;

[0156] represents the slope angle of the river network line K ;

[0157] α represents the angle threshold value.

[0158] In this step, since the silt dam and the river are not exactly vertical, there is an error in the slope angle calculation based on the grid. In order to minimize the influence of the noise line and avoid the loss of the silt dam contour line, as shown in Figure 5 , the best result is obtained when the value of α is 30°.

[0159] Step five, according to the contour line filtered by step four angle, construct cross model to identify silt dam:

[0160] Step 501, merge the contour lines belonging to the same silt dam:

[0161] In this step, as shown in Figure 6 , since multiple lines represent the same silt dam, the contour lines belonging to the same silt dam need to be merged to avoid repeated construction of the cross model.

[0162] The contour lines belonging to the same silt dam are merged, and the judgment basis for merging is that if the difference between the two lines is less than , then they are merged; the center point of the merged line is determined by formula I ; the gentle value SA is set to determine the extension length of the vertical silt dam cross section line.

[0163] ;

[0164] In the formula:

[0165] represents the first line, ;

[0166] represents the first group of silt dam contour lines;

[0167] represents the number of contour lines;

[0168] C represents the total number of groups of silt dam contour lines;

[0169] the midpoint of the first line; x direction coordinate;

[0170] the midpoint of the first line; y direction coordinate;

[0171] the slope of the first line;

[0172] the midpoint of the merged x line;

[0173] the midpoint of the merged y line;

[0174] the slope of the merged line.

[0175] In this step, according to the flat and gentle value SA of the dam top and the terrain of the upstream and downstream areas, whether the elevations of the adjacent two grids are equal is determined to determine the top of the silt dam.

[0176] Step 502, construct a cross model to identify the silt dam:

[0177] Step 50201, the center point is the starting point of the vertical silt dam cross-section line, stored in the first node of the linked list. And the slope of the vertical silt dam cross-section line is calculated by formula II .

[0178] ;

[0179] In the formula:

[0180] represents the slope of the vertical silt dam cross-section line;

[0181] represents the slope of the silt dam cross-section.

[0182] In this step, the standard elevation value ch is initialized as the elevation of the center point U, and the gentle value SA is initialized as 3 under the resolution of 2m, which can be adjusted according to the resolution.

[0183] Step 50202, according to the Bresenham algorithm along the slope k to determine the coordinates of the next point u1. Determine whether the elevation value of u1 and the standard elevation value ch are equal. If they are equal and the flatness value SA = 0, stop extending in this direction and perform step 50203. Otherwise, store u1 in the linked list, and update the flatness value ch to the elevation value of u1. If the two elevations are equal and the flatness value SA > 0, then SA - 1. Then repeat step 50202 to calculate the next point u2.

[0184] Step 50203, according to the inverse slope -k of the vertical dike cross section, repeat steps 50201 to 50202.

[0185] Step 50204, according to the elevation data corresponding to the straight line perpendicular to the dike, determine the highest point coordinates M (x, y) of the elevation.

[0186] Step 50205, according to the dike top point M, extend to both ends along the direction, and the extension stop condition is formula III, and record the number of grids extended to both sides of the dike cross section and respectively, and the cross model is completed.

[0187] ;

[0188] In the formula:

[0189] represents the elevation of the point o determined whether it can be extended;

[0190] represents the elevation of the dike top point M;

[0191] H represents the threshold value of the grid height difference, which is set to 1.5m through experiments.

[0192] If the cross-sectional elevation difference of the constructed cross model is too large, it is considered as noise and discarded. Otherwise, perform step 50206.

[0193] Step 50206, according to each grid on the straight line perpendicular to the dike cross section, extend and to the corresponding direction of the cross section respectively to construct the rectangle corresponding to the dike.

[0194] Step 50207, use morphological closing operation to fill in the missed grids due to grid-based extension.

[0195] Comparative Example 1:

[0196] This comparative example gives a method for detecting a dike, which is performed according to the following steps:

[0197] Step S1: Preprocessing DEM data.

[0198] Step S2: Extracting river network by using the method of filling the depression.

[0199] Step S3: Improving LSD detection according to the mountain shadow map.

[0200] Step S4: Angle filtering according to the results of contour line detection and river network.

[0201] Step S5: Constructing cross model according to the contour line to detect silt dam.

[0202] Comparative Example 2:

[0203] The comparative example gives a silt dam detection method, which is different from the embodiment in that the embodiment extracts the river network by region growing, while the comparative example extracts the river network by the method of splitting the depression, and the remaining steps are the same as the embodiment.

[0204] Comparative Example 3:

[0205] The comparative example gives a silt dam detection method, which is different from the embodiment in step 3, and the remaining steps are the same as the embodiment. In step 3, the comparative example directly uses the results of the LSD detection algorithm without improving the processing.

[0206] Comparative Example 4:

[0207] The comparative example gives a silt dam detection method, which is different from the embodiment in that the comparative example discards the angle filtering of step 4 in the embodiment, and the remaining steps are the same as the embodiment.

[0208] Application Example:

[0209] The application example is based on the comparison of the above-mentioned embodiments and comparative examples 1 to 4. The continuous river network and silt dam positioning method of the application example is tested in Wangmaogou area and part of Zhoutungou area.

[0210] (1) Experimental data:

[0211] Wangmaogou in Suide County, Shaanxi Province, DEM data resolution 0.15m, elevation 940~1200m.

[0212] Part of Zhoutungou area in Ansai District, Shaanxi Province, DEM data resolution 0.15m, elevation 1026~1439m.

[0213] (2) Experimental scheme:

[0214] According to the specific process of the present application, DEM-based silt dam detection is performed. In combination with actual terrain data and the data resolution used, the inclination threshold E is 0.9; the flow accumulation threshold r of Wangmaogou and Zhoutungou is 250 and 440, respectively. The flat value SA at 2m resolution is 3, and the grid height difference threshold H is 1.5m. flow flow

[0215] (3) Experimental results and analysis:

[0216] The spatial distribution of silt dams detected using the method of the present application is shown in Figure 8 Fig. 6. By observing the terrain details, it can be found that in the complex terrain area caused by water and soil erosion and human factors, the river network obtained by region growing can effectively avoid the influence of mountain terraces, can effectively locate the silt dam, and fully reflects the effectiveness of the present application in detecting silt dams. Existing extraction methods mainly focus on calculating the maximum water storage area or sediment area based on multiple data, and then indirectly locating the silt dam. These methods pay less attention to the morphological characteristics of silt dams, resulting in low accuracy in complex terrain areas. The present application does not continue this method, but determines the location of the silt dam based on only one type of high-resolution DEM data and the characteristics of the silt dam. Moreover, both silt dams with "silted" and "not silted" different morphologies can be detected.

[0217] Table 1. Silt dam detection results in Wangmaogou area by different methods

[0218]

[0219] As shown in Table 1, the effect of detecting silt dams by the traditional method of filling depressions in Comparative Example 1 is the worst. As shown in Figure 7 , because the adjacent grids in the flattened depression have the same elevation, it leads to the problem of unrealistic parallel river network. Such river channels that do not conform to the actual river flow will discard a large number of silt dam contour lines during angle filtering, resulting in the worst detection result. The detection effect of the depression splitting method in Comparative Example 2 is improved compared to the filling method, but it is still not ideal. This is because in the silt dam area, the river channel often has a 45°-90° turn, and the silt dam is parallel until a suitable grid is found to pass through the silt dam. This leads to the situation that many contour lines are discarded during angle filtering. In addition, the river channels extracted in Comparative Examples 1 and 2 are mostly close to the two sides of the river channel. This leads to the situation that at the junction of the main river channel and the tributary river channel, the contour lines on both sides of the river channel may be perpendicular to the river channel. The river channel obtained by region growing in the embodiment is located at the center line of the river channel, reducing the occurrence of such problems. The precision and recall of Comparative Example 3 decrease compared to the embodiment. This is because when only the detection result of the LSD is angle filtered, the situation as shown in Figure 4The silt dam contour line intersects with the river contour line. This causes errors in calculating the slope of the silt dam contour line, resulting in the silt dam contour line being determined as a noise line, and the noise line may also be determined as a valid line due to the connection with the silt dam contour line. In Comparative Example 4, only the step of angle filtering is removed compared with the embodiment. Although the silt dam contour line is retained in its entirety, resulting in an increase in the recall rate, the complex terrain contour line caused by erosion inside the river is also retained, causing some noise lines to be misjudged as silt dams. According to the F1 synthesis, the detection effect is poorer than that of the embodiment. Compared with Comparative Examples 1 and 2, it is found that although Comparative Example 4 does not perform angle filtering on the noise lines inside the river, the number of misjudged silt dams is still less than that of Comparative Examples in which different methods are used to extract the river. This shows that the river extracted by the embodiment is more prominent in detecting silt dams.

[0220] The method of the present application is compared with Comparative Examples 1 and 2. The river extracted by the region growing method used in the present application is more in line with the actual situation, and the nearly vertical position relationship between the river and the silt dam can better detect the silt dam. Compared with Comparative Example 3, the improved LSD algorithm of the present application separates straight lines with different slopes to avoid the influence of noise lines as much as possible, thereby improving the recall rate and precision rate of detection. In addition, the comparison of the present application with Comparative Example 4 can find that angle filtering can reduce the influence of noise lines as much as possible while ensuring the precision rate, thereby greatly improving the recall rate of detecting silt dams.

Claims

1. A method for extracting silt dams based on high-resolution DEMs, characterized by, The method comprises the following steps: Step one, data preprocessing: Step 101, reading DEM data; Step 102, smoothing preprocessing; Step two, extracting river network; Step three, improved LSD detection based on mountain shadow map: Based on the mountain shadow map, the improved LSD algorithm is used to separate straight lines with different slopes, and the contour line detection result is obtained; The improved LSD algorithm is: set in preparation output straight line R has been output before the straight line L; from the end point F of straight line R, determine whether there are other straight line pixel points in eight neighborhoods except the straight line; if not, keep the pixel point F; on the contrary, if there are, use least square method to calculate the slope of the two straight lines and their difference; if the difference is within 0-5°, it is considered that the two straight lines belong to the same straight line, and the pixel point F is kept; then check the next pixel point to meet the condition; if the pixel point is encountered, it is discarded, and two straight lines with different slopes are separated; after all the straight lines are determined, only the contour straight lines intersecting the river network are kept; Step four, angle filtering of the contour line obtained in step three based on the river network obtained in step two; Step five, constructing a cross model to identify silt dam according to the contour line filtered in step four: Step 501, merging contour lines belonging to the same silt dam: The contour straight lines belonging to the same check dam are combined, and the combination basis is that if the inclination angle difference of two straight lines is less than , then the two straight lines are combined; the combined straight line center point is determined by formula I ; the extension length of the straight line of the vertical check dam cross section is determined by setting the gentle value SA. ； In the formula: representing a first straight line, ; represents the first group of contour lines of the dam represents the number of contour lines; C total number of sets of straight lines representing the contour of the check dam; representing the midpoint of the line segment x directional coordinate; representing the midpoint of a straight line representing the midpoint of a straight line y directional coordinates; represents the slope of the straight line of the first line; represents the center of the straight line after merging x directional coordinates; represents the center of the straight line after merging y directional coordinates; represents the slope of the straight line after merging; Step 502, constructing a cross model to identify silt dam: Step 50201, center point As the starting point of the straight line of the vertical check dam cross section, the first node of the linked list is stored; and the slope of the straight line of the vertical check dam cross section is calculated by formula II ; ； In the formula: represents the slope of the straight line representing the cross section of the vertical dike; S represents the slope of the cross section of the check dam; Step 50202, determining the coordinates u1 of the next point according to the Bresenham algorithm along the slope k; judging whether the elevation value of u1 is equal to the standard elevation value ch; if they are equal and the flatness value SA is 0, stop extending in this direction and execute step 50203; otherwise, store u1 in the linked list, and update the flatness value ch to the elevation value of u1; if the two elevations are equal and the flatness value SA>0, then SA-1; then repeat step 50202 to calculate the next point u2; Step 50203, repeating steps 50201 to 50202 according to the inverse slope -k of the cross section perpendicular to the silt dam; Step 50204, determining the highest point coordinate M(x, y) according to the elevation data corresponding to the straight line perpendicular to the silt dam; Step 50205, according to the dam crest point M along the direction to both ends, the extension of the cutoff condition is formula III, and the number of grids extending to both sides of the cross section of the check dam is recorded respectively and and , the cross model is completed ； In the formula: point indicating whether or not extension is possible o elevation of Indicates the elevation of the dam's apex M; H represents the threshold value of grid height difference; Step 50206, according to each grid on the straight line of the vertical check cross section, respectively, to the corresponding direction of the cross section to extend and The corresponding rectangle of the check dam is constructed; Step 50207, using morphological closing operation to fill the missed grid due to grid extension based on the grid.

2. The high resolution DEM based silt dam extraction method of claim 1, wherein, The specific process of step 101 is: Step 10101, reading 0.15m resolution DEM data; Step 10102, resampling the 0.15m resolution DEM data into 2m resolution and 10m resolution DEM data.

3. The high resolution DEM based silt dam extraction method of claim 1, wherein, The specific process of step 102 is to use a 3 × 3 window median filter to filter the 10m resolution DEM data.

4. The high resolution DEM based silt dam extraction method of claim 1, wherein, The specific steps of step two are: Step 201, extracting river network by region growing according to the inclination and elevation relationship, wherein the inclination S is defined as: ； In the formula: h x representing a grid x elevation values; h y representing the elevation values of the 8-neighborhood of the grid x of the grid x represents the distance between the grid and the grid , when the grid and the grid are diagonal , otherwise ; Step 202: Set the growth state for each grid cell, and include those not marked as river networks but exceeding the flow accumulation threshold r. flow The flow grid is used as the seed point; Equations V and VI are used as growth criteria for region growth. If the growth criteria of Equation V or VI are met, it is marked as a river network pixel until the location of the boundary point is determined. ； In the formula: a represents the center grid; b representing a grid a of adjacent grids; representing a grid a of elevation values; representing a grid b of elevation values; S represents the inclination; E represents a threshold for the inclination, E takes 0.9; Step 203, obtaining binary river network data with width.

5. The high resolution DEM based silt dam extraction method of claim 4, wherein, The specific process of step 203 is: Step 20301, establish a queue, traverse the flow data; not marked as river network grid elements and greater than the flow accumulation threshold r flow The first flow data (x0, y0) is taken as a growth point, marked as a river network grid, and entered into the queue; Step 20302, taking the head grid element, judging whether its 8-neighborhood grid satisfies formula V or formula VI, and entering the queue if it satisfies; Step 20303, checking whether the 8 adjacent grids in the head grid element are located at the cutoff point and stop growing; Step 20304, sequentially executing steps 20301 to 20303 until the queue is empty; Step 20305, when the growable grid is not completely marked as river network, return to step 20301.

6. The high resolution DEM based silt dam extraction method of claim 1, wherein, The specific process of step four is: The contour straight line obtained in step three is calculated using the least square method and the slope angle of the river network straight line K obtained in step two, respectively and If the two straight lines satisfy formula VII, it is considered that the contour straight line belongs to the contour straight line of the silt dam; otherwise, the contour straight line is noise and discarded; In the formula: representing the slope angle of the contour straight line ; representing the slope angle of the straight line of the river network K of the river course α represents the angle threshold value.

Images