A method for designing terraced fields based on high-precision DEM data

Abstract

This invention belongs to the field of ecological reclamation and terraced field construction technology in mining areas, specifically involving a method for designing slope-to-terraced plots based on high-precision DEM data. Addressing the issues of insufficient accuracy of satellite remote sensing data and the time-consuming and labor-intensive nature of traditional surveying, this invention uses high-precision mountainous DEM data obtained from UAV aerial flights as a basis for automated slope-to-terraced plot design. Computer-automated plot design replaces manual plot design, and the design considers the elevation difference (i.e., ridge height), width, and ridge slope as variables to meet different regulatory requirements.

Description

Technical Field

[0001] This invention belongs to the field of ecological reclamation and terraced field construction technology in mining areas, specifically involving a method for designing slope-to-terraced fields based on high-precision DEM data. Background Technology

[0002] The division of terraced fields should follow the principle of contour lines: fields should be distributed along contour lines to adapt to the terrain, reduce the amount of construction work required, and help conserve water and soil. The current "Technical Specifications for Comprehensive Soil and Water Conservation" also states in its section on sloping farmland management that "in arid and rain-scarce areas of northern my country, the direction of cultivation should generally follow contour lines to facilitate water and soil conservation. In rainy and heavy-soil areas of southern my country, the direction of cultivation should have a 1% to 2% gradient with the contour lines to facilitate drainage and prevent erosion." Modern terraced field design, to accommodate mechanized agricultural operations, requires both a reasonable height for the ridges and specific width requirements for the field surface. The "Technical Specifications for Comprehensive Soil and Water Conservation" provides reference values ​​for both in its section on sloping farmland management.

[0003] Based on these principles and requirements, traditional slope-to-terraced field design often employs a combination of on-site manual surveying and remote sensing image interpretation. This method has two major limitations: first, traditional surveying and manual design of large-scale mountain terraces are time-consuming and labor-intensive; second, satellite remote sensing image data is easily affected by vegetation cover, and publicly available DEM data has low accuracy, which can easily lead to distortion in the representation of the terrain. Summary of the Invention

[0004] To address the issues of insufficient accuracy of satellite remote sensing data and the time-consuming and labor-intensive nature of traditional surveying, this invention provides a method for designing terraced fields based on high-precision DEM data. It utilizes high-precision mountainous DEM data obtained from UAV aerial surveys to automate the design of terraced fields, replacing manual design with computer-automated field design. Furthermore, the design considers the field surface elevation difference (i.e., ridge height), field surface width, and ridge slope as variables to meet different regulatory requirements.

[0005] To achieve the above objectives, the present invention employs the following technical solution:

[0006] This invention provides a method for designing terraced fields based on high-precision DEM data, comprising the following steps:

[0007] S1. Process the UAV laser point cloud data to generate DEM data;

[0008] S2. Use the Gaussian smoothing function to smooth the generated DEM data;

[0009] S3. Generate contour lines based on the smoothed DEM data;

[0010] S4. Check the spacing between the generated contour lines and generate the boundary lines for dividing the terraced fields.

[0011] S5. Generate terrace DEM based on the boundary lines of terraced fields, i.e., the field design results;

[0012] S6. Calculate the earthwork volume variation using the generated DEM data and terrace DEM;

[0013] S7. Calculate the earthwork factor based on the changes in earthwork volume.

[0014] S8. Repeat S4~S5 and S6~S7 to select the final terraced field design scheme based on different field design results and different earthwork factors.

[0015] Furthermore, the specific steps in S1 for processing UAV laser point cloud data to generate DEM data are as follows:

[0016] The drone carries multi-echo laser pulses that penetrate ground vegetation to capture real terrain. Dense ground point cloud data is extracted through filtering algorithms, and an irregular triangular network model is constructed based on the ground points. Finally, high-precision DEM data is generated by interpolation.

[0017] Furthermore, the process of smoothing the generated DEM data using a Gaussian smoothing function in step S2 is as follows:

[0018] Optimize the Gaussian smoothing function:

[0019] ;

[0020] ;

[0021] ;

[0022] ;

[0023] in, , This represents the two-dimensional row and column coordinates (dimensionless index) of a pixel in the DEM raster data. This is the elevation value of the pixel (in meters). , The row and column coordinates of the center cell of the smoothing window. Elevation value (in meters) of the center pixel of the smoothing window. The horizontal standard deviation of the Gaussian kernel represents the smoothing radius (in meters) in geographic space. The vertical Gaussian kernel standard deviation, i.e., the elevation sensitivity parameter (unit: meters), is used to control the impact of elevation differences on the weights. , Spatial resolution of DEM data in the x and y directions (unit: meters / pixel).

[0024] Represents the geographic distance (in meters) of a pixel relative to the center pixel in the x-direction.

[0025] Represents the geographic distance (in meters) of a pixel relative to the center pixel in the y-direction.

[0026] Represents the elevation difference (in meters) between a pixel and the center pixel.

[0027] The DEM data generated from point cloud processing is smoothed using an optimized Gaussian smoothing function.

[0028] Furthermore, the generation of contour lines based on the smoothed DEM data in step S3 specifically involves:

[0029] Based on the elevation value of each pixel in the smoothed DEM data and the set contour interval h, the moving square contour line tracing algorithm is adopted. Taking each grid pixel as the basic unit, the intersection point of the contour value on the grid boundary is calculated by bilinear interpolation. The contour points are traced and connected pixel by pixel to form a discrete polyline. After curve smoothing, topology checking and broken line removal, a continuous, smooth contour line with elevation attributes is finally generated.

[0030] Furthermore, in step S4, the spacing between the generated contour lines is checked, and the boundary lines for dividing the terraced fields are generated as follows:

[0031] S41. Element Input: Input the boundary range of the terraced fields, contour lines, and minimum width requirement for the fields.

[0032] S42. Contour Line Sorting: Take any contour line within the input terrace boundary as the starting line, and start numbering from the starting line. The starting line is the No. 1 contour line. The contour lines on its left and right sides are arranged in a cross pattern and numbered. If there are multiple contour lines with the same elevation, they are numbered with adjacent numbers. The total number of contour line elements is recorded as n.

[0033] S43. Inspection of contour lines with unqualified spacing: Check according to the minimum field width requirement. Check the buffer zone generated by the current contour line A using the minimum field width requirement + basic micro tolerance. If the buffer zone generated by the current contour line A using the minimum field width requirement + basic micro tolerance intersects with a certain contour line B, it means that there is a part between contour line A and contour line B with a width smaller than the minimum field width requirement + basic micro tolerance. Then, contour line B is regarded as an unqualified contour line and needs to be processed in S44.

[0034] During the inspection, all contour lines are sequentially identified by their numbers in ascending order and are used as the current contour line for inspection. Only contour lines with numbers greater than the current contour line number are inspected. After inspecting all contour lines with numbers greater than the current contour line number for each current contour line, it proceeds to step S44 for further processing.

[0035] S44. Handling of contour lines with unqualified spacing: If contour line B is found to be unqualified in S43 by checking the buffer generated by the minimum field width requirement + basic micro tolerance for the current contour line A, then the same check is performed on contour line B here by using the buffer generated by the minimum field width requirement of contour line A + graded micro increment * m, where m is the number of contour line B.

[0036] If contour line B fails both inspections (S43 and S44), proceed as follows; otherwise, skip this step:

[0037] (1) Divide the unqualified contour line B into the first part and the second part:

[0038] The first part is the portion located within the buffer zone generated by the minimum field width requirement plus the graded micro-increment * m at contour line A. The distance between this part and contour line A is not up to standard.

[0039] The second part is the portion outside the buffer zone generated by the minimum field width requirement plus the graded micro-increment * m located at contour line A. The distance between this portion and contour line A is acceptable.

[0040] (2) Divide the buffer zone boundary generated by the minimum field width requirement of contour line A + hierarchical micro-increment * m into the third and fourth parts:

[0041] The third part is the section between contour line A and the unqualified contour line B. The distance from the buffer zone boundary of this part to contour line A is smaller than the distance between contour line A and unqualified contour line B at the corresponding position.

[0042] The fourth part is the portion not located between contour line A and the unqualified contour line B. The distance from the boundary of this buffer zone to contour line A is greater than the distance between contour line A and the unqualified contour line B at the corresponding location.

[0043] (3) Processing Part 1, Part 2, Part 3 and Part 4:

[0044] Delete the first and third parts, keep the second and fourth parts, and merge the remaining second and fourth parts into a new contour line. The new contour line meets the minimum field width requirement between contour line A and contour line B.

[0045] S45. Generate terraced field boundary lines: Use the original qualified contour lines and the processed qualified contour lines to form the terraced field boundary lines, and retain the elevation and number of the original qualified contour lines.

[0046] Furthermore, the specific steps in S5 for generating a terraced field DEM based on the boundary lines of terraced fields are as follows:

[0047] S51. Design of field embankments: Assume a profile perpendicularly cuts the terrain through the profile line. Calculate: the horizontal width of the field embankment. The boundary lines of the terraced fields are shifted to both sides by half the horizontal width d of the ridge, resulting in the upper and lower ridge lines. The upper ridge line retains the original elevation H of the terraced field boundary line, while the lower ridge line is set to elevation Hh.

[0048] The upper and lower ridge lines obtained by offsetting the boundary line of the terraced field, and the slope between them, together form the ridge at the boundary line of the terraced field.

[0049] The upper field surface is composed of the lower field ridge line obtained by offsetting the boundary line of the previous terraced field, the upper field ridge line obtained by offsetting the boundary line of the previous terraced field, and the plane between the two; the lower field surface is composed of the lower field ridge line obtained by offsetting the boundary line of the previous terraced field, the upper field ridge line obtained by offsetting the boundary line of the next terraced field, and the plane between the two.

[0050] S52. Generate Terraced Field DEM: Based on the designed field ridges and terrace boundaries, extract the field surface with elevation information, create an irregular triangular network model of the terraced fields using the field surface with elevation information, and then convert the irregular triangular network model into raster data, which is the terraced field DEM.

[0051] Furthermore, the calculation of earthwork volume changes in S6 using the generated DEM data and terrace DEM specifically involves:

[0052] The difference DEM is the difference between the terrace DEM of S5 and the DEM data generated by S1, that is, terrace DEM - original terrain DEM = difference DEM;

[0053] In the difference DEM, negative value rasters are cut rasters, positive value rasters are fill rasters, and zero value rasters are stationary rasters;

[0054] The absolute value of the excavation volume equals the volume of earthwork.

[0055] =Fill volume;

[0056] MAX (excavation volume, fill volume) = earthwork volume;

[0057] Excavation volume - Fill volume = Required soil volume. A positive value indicates surplus soil, while a negative value indicates shortage soil.

[0058] Furthermore, the calculation of the earthwork factor based on the change in earthwork volume in S7 is specifically as follows:

[0059] Earthwork factor T = volume of earth moved;

[0060] Earthwork factor T = -Required earth volume;

[0061] Earthwork factor T= - Soil required;

[0062] The above three methods for calculating earthwork factors should be selected according to the actual situation.

[0063] Compared with the prior art, the present invention has the following advantages:

[0064] 1. An optimized 3D anisotropic Gaussian smoothing algorithm separates horizontal geographic distance and vertical elevation difference as independent weighting factors, introducing DEM spatial resolution to convert raster indexes into actual geographic distances. The horizontal smoothing radius and elevation sensitivity can be independently adjusted, effectively removing DEM spikes while preserving key micro-topographic features such as steep slopes and gullies to the greatest extent possible, avoiding terrain flattening distortion. By directly using actual geographic distance and elevation difference to calculate weights, rather than raster indexes, the smoothing effect is independent of DEM resolution, adapting to data sources of varying precision.

[0065] 2. A method for generating field boundary lines using contour lines with a buffer zone is employed. This method involves "initial screening - graded verification - local splitting and modification" of contour lines, locally cutting, deleting, and modifying unqualified sections rather than completely eliminating contour lines. While ensuring that field width and ridge height meet standards, the original contour line orientation is preserved as much as possible, adhering to the contour-driven principles of soil and water conservation. A graded micro-increment mechanism and a cross-numbering contour line sorting method jointly ensure the correct topological relationship of contour lines during processing. This guarantees that the terraced fields generated from the processed contour lines are distributed in a stepped manner according to the original terrain, fundamentally reducing changes in earthwork volume before and after design.

[0066] 3. Three earthwork factors are designed, focusing on the amount of earth moved, the amount of earth required, and the weighted value of the amount of earth moved and the amount of earth required, respectively. Then, through an iterative mechanism, the optimal design scheme for the target of interest can be selected.

[0067] 4. By using the height, slope, and width of the field ridges as variable parameters, it is possible to quickly adapt to the requirements of the "Technical Specifications for Comprehensive Soil and Water Conservation" in different regions, as well as the operating width of different agricultural machinery. Attached Figure Description

[0068] Figure 1 This is a flowchart of a slope-to-terraced field design method based on high-precision DEM data according to the present invention;

[0069] Figure 2 This is a schematic diagram for checking contour lines with non-compliant spacing;

[0070] Figure 3 This is a schematic diagram of the first and second parts;

[0071] Figure 4 This is a schematic diagram of Part 3 and Part 4;

[0072] Figure 5 A schematic diagram of the new qualified contour lines;

[0073] Figure 6 A diagram illustrating point 1 in detail;

[0074] Figure 7 A diagram illustrating detail 2;

[0075] Figure 8 A schematic diagram of the design for the paddy field ridges;

[0076] Wherein, 1 is contour line A, 2 is the buffer zone generated by the minimum field width requirement + basic micro-tolerance for the current contour line A, 3 is contour line B, 4 is the first part, 5 is the second part, 6 is the third part, 7 is the fourth part, 8 is the buffer zone boundary generated by the minimum field width requirement + hierarchical micro-increment * m for contour line A, 9 is the first original contour line, 10 is the second original contour line, 11 is the third original contour line, 12 is the field width buffer zone boundary of the first original contour line with the minimum field width requirement + basic micro-tolerance, 13 is the field width buffer zone boundary of the first original contour line with the minimum field width requirement, 14 is the buffer zone of the right side of the first original contour with the minimum field width requirement + basic micro-tolerance, 15 is the modified second contour line, 16 is the modified third contour line, 17 is the boundary line of the terraced field, 18 is the field ridge, 19 is the upper field ridge line, 20 is the lower field ridge line, 21 is the profile line, 22 is the upper field surface, and 23 is the lower field surface. Detailed Implementation

[0077] To further illustrate the technical solution of the present invention, the present invention will be further described below through embodiments.

[0078] like Figure 1 As shown in the figure, a method for designing terraced fields based on high-precision DEM data in this embodiment includes the following steps:

[0079] S1. Process point cloud data into DEM:

[0080] The process of generating DEM data from UAV laser point cloud data is as follows:

[0081] The multi-echo laser pulses carried on the drone can penetrate the ground vegetation to a certain extent, capture the real terrain, extract dense ground point cloud data through filtering algorithms, and build an irregular triangular network model (TIN) based on the ground points. Finally, high-precision DEM data (digital elevation model) is generated by interpolation.

[0082] S2 and DEM optimization and smoothing:

[0083] The generated DEM data is smoothed using a Gaussian smoothing function, specifically as follows:

[0084] Two-dimensional Gaussian smoothing function:

[0085] ;

[0086] , : Row and column coordinates of the corresponding cell in the raster DEM;

[0087] : The horizontal coordinates of the Gaussian kernel center point, i.e., the coordinates of the center cell of the raster DEM smoothing window;

[0088] Gaussian standard deviation is a parameter that controls the smoothness of the surface.

[0089] Isotropic 3D Gaussian smoothing function:

[0090] ;

[0091] , : Row and column coordinates of the corresponding cell in the raster DEM;

[0092] : Elevation values ​​of pixels in a raster DEM;

[0093] : The horizontal coordinates of the Gaussian kernel center point are taken from the center cell coordinates of the raster DEM smoothing window;

[0094] : The mean of Gaussian distribution in the elevation direction, taken as the elevation value of the center cell of the smoothed window of the raster DEM;

[0095] Gaussian standard deviation is a parameter that controls the smoothness of the surface.

[0096] Anisotropic 3D Gaussian smoothing function:

[0097] ;

[0098] , : Row and column coordinates of the corresponding cell in the raster DEM;

[0099] : Elevation values ​​of pixels in a raster DEM;

[0100] : The horizontal coordinates of the Gaussian kernel center point are taken from the center cell coordinates of the raster DEM smoothing window;

[0101] : The mean of Gaussian distribution in the elevation direction, taken as the elevation value of the center cell of the smoothed window of the raster DEM;

[0102] , , Independent Gaussian standard deviations in the x, y, and z directions; parameters that control the smoothness.

[0103] In the smoothing process of DEM data, the purpose of smoothing is to eliminate small surface abrupt changes without altering the actual terrain, so as to facilitate the subsequent generation of contour lines that conform to the actual terrain and are relatively smooth.

[0104] The above three Gaussian functions are used for smoothing DEM data:

[0105] The two-dimensional Gaussian function is used to smooth DEM data. It only considers the planar distance between raster cells and does not consider the difference in the elevation values ​​of raster cells. Therefore, it cannot preserve terrain features such as cliff boundaries where the elevation changes abruptly within a short planar distance.

[0106] The isotropic 3D Gaussian function is used to smooth DEM data. The smoothing effect is the same in the x, y and z directions. However, in actual smoothing, the same horizontal distance and elevation difference have different effects on terrain features. Therefore, two independent sets of parameters should be used for smoothing in the horizontal and elevation directions.

[0107] Anisotropic 3D Gaussian functions are used to smooth DEM data. The smoothing effect is different in the x, y and z directions. This results in the smoothing effect of the terrain being related to the planar direction. The smoothing effect is not the same in different planar directions, which may cause deviation or distortion of the smoothed terrain from the original terrain in the direction.

[0108] Optimize the Gaussian smoothing function:

[0109] ;

[0110] ;

[0111] ;

[0112] ;

[0113] in, , This represents the two-dimensional row and column coordinates (dimensionless index) of a pixel in the DEM raster data. This is the elevation value of the pixel (in meters). , The row and column coordinates of the center cell of the smoothing window. Elevation value (in meters) of the center pixel of the smoothing window. The horizontal standard deviation of the Gaussian kernel represents the smoothing radius (in meters) in geographic space. The vertical Gaussian kernel standard deviation, i.e., the elevation sensitivity parameter (unit: meters), is used to control the impact of elevation differences on the weights. , Spatial resolution of DEM data in the x and y directions (unit: meters / pixel).

[0114] Represents the geographic distance (in meters) of a pixel relative to the center pixel in the x-direction.

[0115] Represents the geographic distance (in meters) of a pixel relative to the center pixel in the y-direction.

[0116] Represents the elevation difference (in meters) between a pixel and the center pixel.

[0117] when As the value approaches 0, the three-dimensional Gaussian smoothing formula approximates the two-dimensional Gaussian smoothing formula. Increase The overall size is reduced, thus decreasing the weight of each pixel in the smoothing calculation. This can be achieved through proper design. The size can be adjusted to more effectively eliminate minor abrupt changes on the surface while preserving the characteristics of the terrain.

[0118] Use the same in both the x and y directions. The parameters ensure that within the smooth window, the weight distribution in the horizontal direction is only related to distance and not to direction.

[0119] The DEM data generated from point cloud processing is smoothed using an optimized Gaussian smoothing function.

[0120] S3 and DEM data generate contour lines:

[0121] Contour lines are generated based on the smoothed DEM data, specifically as follows:

[0122] Based on the elevation value of each pixel in the smoothed DEM data and the set contour interval h (which is the height of the field ridge in the slope-to-terrace design), the moving square contour line tracing algorithm is adopted. Taking each grid pixel as the basic unit, the intersection of contour values ​​on the grid boundary is calculated by bilinear interpolation. The contour points are traced and connected pixel by pixel to form discrete polylines. After curve smoothing, topology checking and broken line removal, continuous, smooth contour lines with elevation attributes are finally generated.

[0123] S4. Contour line processing generates the boundary lines for dividing terraced fields:

[0124] The generated contour lines are checked for spacing, and the boundary lines for dividing the terraced fields are generated, specifically as follows:

[0125] S41. Element Input: Input the boundary range of the terraced fields (surface file), the contour line (line file), and the minimum width requirement of the field (in meters).

[0126] S42. Contour Line Sorting: Take any contour line within the input terrace boundary as the starting line, and start numbering from the starting line. The starting line is the No. 1 contour line. The contour lines on its left and right sides are arranged in a cross pattern and numbered. If there are multiple contour lines with the same elevation, they are numbered with adjacent numbers. The total number of contour line elements is recorded as n.

[0127] S43. Inspection of contour lines with unqualified spacing ( Figure 2 The process begins by checking against the minimum field width requirement. A buffer zone (2) generated from the minimum field width requirement plus the basic micro-tolerance is used to check the current contour line A. If this buffer zone intersects with contour line B3, it indicates that there is a portion between contour line A1 and contour line B3 with a width less than the minimum field width requirement plus the basic micro-tolerance. Contour line B3 is then considered unqualified and requires processing in step S44. Here, the basic micro-tolerance (BMT) is an abbreviation for Base Micro Tolerance, a relatively small value. Similarly, the graded micro-increment (GMI) mentioned below is an abbreviation for Graded Micro Increment, also a relatively small value.

[0128] During the inspection, all contour lines are sequentially identified by their numbers in ascending order and are inspected. Only contour lines with numbers greater than the current contour line number are inspected. After inspecting all other contour lines with numbers greater than the current contour line number, each current contour line proceeds to step S44 for further processing.

[0129] S44. Handling of non-compliant contour lines: If, in S43, the minimum field width requirement + basic micro-tolerance is used to check the buffer 2 generated by the current contour line A and find that contour line B3 is a non-compliant contour line, then the same check is performed on contour line B3 here using the buffer generated by the minimum field width requirement of contour line A1 + hierarchical micro-increment * m (width + GMI * m), where m is the number of contour line B3; BMT is greater than or equal to n times GMI, where n is the total number of contour line elements counted in S42; ;

[0130] If contour line B3 fails both checks by S43 and S44, proceed as follows; otherwise, skip this step:

[0131] (1) Divide the unqualified contour line B3 into part 4 and part 5. Figure 3 ):

[0132] Part 4 is the portion located within the buffer zone generated by the minimum field width requirement plus the graded micro-increment * m at contour line A1. The distance between this portion and contour line A1 is not up to standard.

[0133] Part 2, section 5, is the portion outside the buffer zone generated by the minimum field width requirement plus the graded micro-increment * m located at contour line A1. The distance between this portion and contour line A1 is acceptable.

[0134] (2) Divide the minimum field width requirement of contour line A + the buffer boundary 8 generated by the hierarchical micro-increment *m into the third part 6 and the fourth part 7. Figure 4 ):

[0135] Part 3, section 6, is the section between contour line A1 and the unqualified contour line B3. The distance from the buffer zone boundary of this section to contour line A1 is smaller than the distance between contour line A1 and the unqualified contour line B3 at the corresponding position.

[0136] Part 4, section 7, is the portion not located between contour line A1 and the unqualified contour line B3. The distance from the buffer zone boundary of this portion to contour line A1 is greater than the distance between contour line A1 and the unqualified contour line B3 at the corresponding position.

[0137] (3) Process Part 1, Part 4, Part 2, Part 5, Part 3, Part 6, and Part 4, Part 7. Figure 5 ):

[0138] Delete part 4 and part 6, keep part 5 and part 7, and merge the remaining part 5 and part 7 into a new qualified contour line. The new qualified contour line meets the minimum field width requirement between contour line A1 and contour line B3.

[0139] S45. Generate terraced field boundary lines: Use the original qualified contour lines and the processed qualified contour lines to form the terraced field boundary lines, and retain the elevation and number attributes of the original qualified contour lines.

[0140] S46. Detailed Explanation

[0141] The BMT described in S43 and S44 must be greater than the product of the total number of contour line features n and the GMI, that is... It refers to the maximum number of contour line elements;

[0142] BMT ensures that the width of the inspection buffer is greater than the minimum field width requirement. BMT and GMI together ensure that when there is more than one non-compliant contour line on one side of the buffer 2 generated by the current contour line A, the non-compliant contour lines are arranged in the range from the minimum field width requirement to the minimum field width requirement plus the basic micro-tolerance, with GMI as the interval, according to the original topological relationship.

[0143] For example: Figure 6 In the middle, the second original contour line 10 and the third original contour line 11 are located within the buffer zone 14 of the minimum field width requirement + basic micro tolerance (width + BMT) to the right of the first original contour line, that is, the second original contour line 10 and the third original contour line 11 are unqualified contour lines. Figure 7 In the original contour, the modified second contour line 15 and the modified third contour line 16 are located to the right of the minimum field width requirement buffer boundary 13 of the first original contour line at a distance of two times GMI and three times GMI, respectively, and to the left of the minimum field width requirement + basic micro-tolerance width buffer boundary 12 of the first original contour line; this ensures that the topological relationship between contour lines remains unchanged before and after the modification, and avoids the phenomenon of contour lines crossing or overlapping.

[0144] S5. Generate terraced field DEM based on the boundary lines of terraced field blocks:

[0145] The terraced field DEM, i.e., the field design result, is generated based on the boundary lines of the terraced fields. Specifically:

[0146] Based on the requirements for the height and slope of the field ridges given in the "Technical Specification for Comprehensive Soil and Water Conservation Management" for sloping farmland management, and the contour interval h set in S3, the field ridges are designed.

[0147] The contour interval h is the height of the field ridge. The appropriate field ridge slope α is selected according to the specifications.

[0148] S51, Design of field ridges: such as Figure 8 As shown, suppose a profile cuts the terrain perpendicularly through profile line 21. Calculate: the horizontal width of the field ridge. Offset the terraced field boundary line 17 to both sides by half the horizontal width d of the ridge, resulting in the upper ridge line 19 and the lower ridge line 20. The upper ridge line 19 retains the original elevation H of the terraced field boundary line 17, while the lower ridge line 20 is set to elevation Hh.

[0149] The upper ridge line 19, the lower ridge line 20, and the slope between them, obtained by offsetting the boundary line 17 of the terraced field, together form the ridge 18 at the boundary line 17 of the terraced field.

[0150] The lower ridge line obtained by offsetting the boundary line of one terraced field, the upper ridge line obtained by offsetting the boundary line of the next terraced field, and the plane between them together constitute the field surface. For example... Figure 8 In the middle, the upper field surface 22 is composed of the lower field ridge line obtained by offsetting the boundary line of the previous terraced field and the upper field ridge line 19 obtained by offsetting the boundary line of the terraced field 17, and the plane between the two. The lower field surface 23 is composed of the lower field ridge line 20 obtained by offsetting the boundary line of the terraced field 17 and the upper field ridge line obtained by offsetting the boundary line of the next terraced field, and the plane between the two.

[0151] S52. Generate Terraced Field DEM: Based on the designed field ridges and terraced field boundaries, extract the field surface with elevation information. In software such as ArcGIS, an irregular triangular mesh model of the terraced field can be created from the field surface with elevation information. Then, convert the irregular triangular mesh model into raster data, which is the terraced field DEM.

[0152] S6. Calculate the change in earthwork volume:

[0153] The earthwork volume variation was calculated using DEM data generated from UAV laser point cloud data processing and terraced DEM data. Specifically:

[0154] The difference DEM is the difference between the terrace DEM of S5 and the DEM data generated by S1, that is, terrace DEM - original terrain DEM = difference DEM;

[0155] In the difference DEM, negative value rasters are cut rasters, positive value rasters are fill rasters, and zero value rasters are stationary rasters;

[0156] The absolute value of the excavation volume equals the volume of earthwork.

[0157] =Fill volume;

[0158] MAX (excavation volume, fill volume) = earthwork volume;

[0159] Excavation volume - Fill volume = Required soil volume. A positive value indicates surplus soil, while a negative value indicates shortage soil.

[0160] S7. Calculate the earthwork factor T:

[0161] The earthwork factor is calculated based on the change in earthwork volume, specifically as follows:

[0162] Earthwork factor T = Earthwork volume (the less the earthwork changes, the smaller the T factor).

[0163] Earthwork factor T = -Earthwork requirement (the less earthwork required, the smaller the T factor).

[0164] Earthwork factor T= - Soil demand (the less the change in earthwork, the less earthwork demand, and the smaller the T factor).

[0165] The above three methods for calculating earthwork factors should be selected according to the actual situation.

[0166] S8. Select the final terraced field design scheme:

[0167] Repeat steps S4-S5 and S6-S7, and select the final terraced field design scheme based on different field design results and different earthwork factors, as follows:

[0168] In S4, there are n contour lines. Taking each contour line as the starting line, repeating S4~S5 can yield n different field design results (including field ridges). Repeating S6~S7 can calculate n different earthwork factors T. Finally, the terraced field design scheme is selected manually by combining all field design results and earthwork volume.

[0169] The foregoing has shown and described the main features and advantages of the present invention. It will be apparent to those skilled in the art that the present invention is not limited to the details of the exemplary embodiments described above, and that the invention can be implemented in other specific forms without departing from its spirit or essential characteristics. Therefore, the embodiments should be considered exemplary and non-limiting in all respects, and the scope of the invention is defined by the appended claims rather than the foregoing description. Thus, it is intended that all variations falling within the meaning and scope of equivalents of the claims be included within the present invention.

[0170] Furthermore, it should be understood that although this specification describes embodiments, not every embodiment contains only one independent technical solution. This narrative style is merely for clarity. Those skilled in the art should consider the specification as a whole, and the technical solutions in each embodiment can also be appropriately combined to form other embodiments that can be understood by those skilled in the art.

Claims

1. A method for designing terraced fields on slopes based on high-precision DEM data, characterized in that, Includes the following steps: S1. Process the UAV laser point cloud data to generate DEM data; S2. Use the Gaussian smoothing function to smooth the generated DEM data; S3. Generate contour lines based on the smoothed DEM data; S4. Check the spacing between the generated contour lines and generate the boundary lines for dividing the terraced fields. S5. Generate terrace DEM based on the boundary lines of terraced fields, i.e., the field design results; S6. Calculate the earthwork volume variation using the generated DEM data and terrace DEM; S7. Calculate the earthwork factor based on the changes in earthwork volume. S8. Repeat S4~S5 and S6~S7 to select the final terraced field design scheme based on different field design results and different earthwork factors.

2. The method for designing terraced fields based on high-precision DEM data according to claim 1, characterized in that, The specific steps in S1 to process UAV laser point cloud data into DEM data are as follows: The drone carries multi-echo laser pulses that penetrate ground vegetation to capture real terrain. Dense ground point cloud data is extracted through filtering algorithms, and an irregular triangular network model is constructed based on the ground points. Finally, high-precision DEM data is generated by interpolation.

3. The method for designing terraced fields based on high-precision DEM data according to claim 2, characterized in that, The Gaussian smoothing function is used in S2 to smooth the generated DEM data, specifically as follows: Optimize the Gaussian smoothing function: ； ； ； ； In the formula, , Let the coordinates of a pixel be the two-dimensional row and column coordinates of a cell in the DEM raster data. This is the elevation value of that pixel. v represents the row and column coordinates of the center pixel of the smoothing window. To smooth the elevation value of the center pixel of the window, The horizontal standard deviation of the Gaussian kernel represents the smoothing radius in geographic space. This represents the standard deviation of the Gaussian kernel in the vertical direction, i.e., the elevation sensitivity parameter. , This refers to the spatial resolution of the DEM data in the x and y directions. Represents: the geographical distance of a pixel relative to the center pixel in the x-direction; Represents: the geographical distance of a pixel relative to the center pixel in the y-direction; Representative: The elevation difference between a pixel and the center pixel; The DEM data generated from point cloud processing is smoothed using an optimized Gaussian smoothing function.

4. The method for designing terraced fields based on high-precision DEM data according to claim 3, characterized in that, The specific steps in S3 for generating contour lines based on the smoothed DEM data are as follows: Based on the elevation value of each pixel in the smoothed DEM data and the set contour interval h, the moving square contour line tracing algorithm is adopted. Taking each grid pixel as the basic unit, the intersection point of the contour value on the grid boundary is calculated by bilinear interpolation. The contour points are traced and connected pixel by pixel to form a discrete polyline. After curve smoothing, topology checking and broken line removal, a continuous, smooth contour line with elevation attributes is finally generated.

5. The method for designing terraced fields based on high-precision DEM data according to claim 4, characterized in that, In step S4, the spacing between the generated contour lines is checked, and the boundary lines for dividing the terraced fields are generated as follows: S41. Element Input: Input the boundary range of the terraced fields, contour lines, and minimum requirements for the width of the fields; S42. Contour Line Sorting: Take any contour line within the input terrace boundary as the starting line, and start numbering from the starting line. The starting line is the No. 1 contour line. The contour lines on its left and right sides are arranged in a cross pattern and numbered. If there are multiple contour lines with the same elevation, they are numbered with adjacent numbers. The total number of contour line elements is recorded as n. S43. Inspection of contour lines with unqualified spacing: Check according to the minimum field width requirement. Check the buffer zone generated by the current contour line A using the minimum field width requirement + basic micro tolerance. If the buffer zone generated by the current contour line A using the minimum field width requirement + basic micro tolerance intersects with a certain contour line B, it means that there is a part between contour line A and contour line B with a width smaller than the minimum field width requirement + basic micro tolerance. Then, contour line B is regarded as an unqualified contour line and needs to be processed in S44. During the inspection, all contour lines are used as the current contour lines in ascending order of their numbers, and only the other contour lines with numbers greater than the current contour line number are inspected. After each current contour line has been inspected and all other contour lines with numbers greater than the current contour line number have been inspected, it enters S44 for processing. S44. Handling of unqualified contour lines: If contour line B is found to be unqualified in S43 by checking the buffer generated by the minimum field width requirement + basic micro tolerance for the current contour line A, then the same check is performed on contour line B here by using the buffer generated by the minimum field width requirement of contour line A + graded micro increment * m, where m is the number of contour line B. If contour line B fails both checks by S43 and S44, proceed as follows; otherwise, skip this step: The substandard contour line B is divided into two parts: Part 1 and Part 2. The first part is the portion located within the buffer zone generated by the minimum field width requirement plus the graded micro-increment * m at contour line A. The distance between this part and contour line A is not up to standard. The second part is the portion outside the buffer zone generated by the minimum field width requirement plus the graded micro-increment * m located at contour line A. The distance between this portion and contour line A is acceptable. The minimum field width requirement of contour line A, plus the hierarchical incremental *m, is used to divide the buffer boundary into the third and fourth parts: The third part is the section between contour line A and the unqualified contour line B. The distance from the buffer zone boundary of this part to contour line A is smaller than the distance between contour line A and unqualified contour line B at the corresponding position. The fourth part is the portion not located between contour line A and the unqualified contour line B. The distance from the boundary of this buffer zone to contour line A is greater than the distance between contour line A and the unqualified contour line B at the corresponding location. Process Part 1, Part 2, Part 3, and Part 4: Delete the first and third parts, keep the second and fourth parts, and merge the remaining second and fourth parts into a new contour line. The new contour line meets the minimum field width requirement between contour line A and contour line B. S45. Generate terraced field boundary lines: Use the original qualified contour lines and the processed qualified contour lines to form the terraced field boundary lines, and retain the elevation and number of the original qualified contour lines.

6. The method for designing terraced fields based on high-precision DEM data according to claim 5, characterized in that, The specific steps in S5 for generating a terraced field DEM based on the boundary lines of terraced fields are as follows: S51. Design of field embankments: Assume a profile perpendicularly cuts the terrain through the profile line. Calculate: the horizontal width of the field embankment. The boundary lines of the terraced fields are shifted to both sides by half the horizontal width d of the ridge, resulting in the upper and lower ridge lines. The upper ridge line retains the original elevation H of the terraced field boundary line, while the lower ridge line is set to elevation Hh. The upper and lower ridge lines obtained by offsetting the boundary line of the terraced field, and the slope between them, together form the ridge at the boundary line of the terraced field. The upper field surface is composed of the lower field ridge line obtained by offsetting the boundary line of the previous terraced field, the upper field ridge line obtained by offsetting the boundary line of the previous terraced field, and the plane between the two; the lower field surface is composed of the lower field ridge line obtained by offsetting the boundary line of the previous terraced field, the upper field ridge line obtained by offsetting the boundary line of the next terraced field, and the plane between the two. S52. Generate Terraced Field DEM: Based on the designed field ridges and terrace boundaries, extract the field surface with elevation information, create an irregular triangular network model of the terraced fields using the field surface with elevation information, and then convert the irregular triangular network model into raster data, which is the terraced field DEM.

7. The method for designing terraced fields based on high-precision DEM data according to claim 6, characterized in that, The calculation of earthwork volume changes in S6 using the generated DEM data and terrace DEM specifically involves: The difference DEM is the difference between the terrace DEM of S5 and the DEM data generated by S1, that is, terrace DEM - original terrain DEM = difference DEM; In the difference DEM, negative value rasters are cut rasters, positive value rasters are fill rasters, and zero value rasters are stationary rasters; The absolute value of the excavation volume equals the volume of earthwork. =Fill volume; MAX (excavation volume, fill volume) = earthwork volume; Excavation volume - Fill volume = Required soil volume. A positive value indicates surplus soil, while a negative value indicates shortage soil.

8. The method for designing terraced fields based on high-precision DEM data according to claim 7, characterized in that, The calculation of the earthwork factor based on the change in earthwork volume in S7 is specifically as follows: Earthwork factor T = volume of earth moved; Earthwork factor T = -Required earth volume; Earthwork factor T= - Soil required; The above three methods for calculating earthwork factors should be selected according to the actual situation.

Images