A method for vegetation removal in high-slope areas

Abstract

The application discloses a high-slope-area vegetation removing method and belongs to the technical field of information processing. The method comprises the following steps: projecting all points in a point set P to be processed onto an xoy plane, performing 2D Delaunay triangulation to obtain a triangular mesh set T; calculating the slope and slope direction of each point in the point set P; classifying the point set P according to the slope and slope direction; finding a seed triangle t in the triangular mesh set T seed ; gradually changing the z coordinate of a vegetation area point by using a region growing method; regarding a point whose z coordinate value has been changed and a point whose z coordinate value of an adjacent point has been changed as an untrusted point; and deleting the untrusted point to obtain a point set P' after vegetation removing. The method can obtain more accurate results in earthwork calculation, solves the problem of cost control, and has higher accuracy in landslide prediction, so that the loss caused by landslides to agriculture can be better prevented.

Description

Technical Field

[0001] This invention relates to the field of information processing technology, and specifically to a method for removing vegetation in high-slope areas. Background Technology

[0002] The construction of Digital Terrain Models (DTMs) has always been an important research area in the surveying and mapping industry. For example, by obtaining ground model data from different periods and calculating elevation differences to obtain earthwork volumes, it has been widely applied in basic engineering projects such as agriculture, construction, highways, and water conservancy. Accurate and feasible earthwork calculations enable the rational planning of engineering schemes. Another example is landslide monitoring, which requires accurate ground model data of landslide-prone areas to predict the probability of landslides and provide early warnings of landslide disasters. The construction of a DTM first requires obtaining the three-dimensional coordinates of ground points in the relevant geographical area. Current main methods include remote sensing, LiDAR laser scanning, and UAV photogrammetry. UAV photogrammetry, in particular, has become an important method for acquiring topographic point clouds due to its ability to quickly acquire large-area point cloud data, high accuracy, and strong initiative. However, since the areas where data is collected are often covered by vegetation, the three-dimensional point clouds converted from UAV images often contain vegetation points in addition to ground points, thus failing to capture the true shape of the ground. This leads to significant errors in earthwork volume calculations and inaccurate landslide prediction results. Therefore, the ability to quickly and accurately distinguish between surface point cloud data and vegetation data from a large amount of scattered point cloud data has become an important condition for accurate earthwork volume calculation and landslide monitoring.

[0003] In recent years, experts and scholars both domestically and internationally have proposed various vegetation removal methods. Based on different strategies, these can be mainly divided into slope methods, mathematical morphology methods, and surface fitting methods. Slope methods typically assume that the slope within a terrain area is gradual, while the slope between trees and the ground varies significantly. Based on this, Vosselman, San, and Meng et al. proposed different methods for calculating the slope of adjacent points for vegetation filtering, but all struggled to obtain the optimal slope threshold applicable to different terrain features. Subsequent studies proposed various automatic threshold definitions, such as adaptive filters and bidirectional filters, but the filtering accuracy of vegetation decreased with increasing slope. Mathematical morphology methods use a fixed-size structuring element window to perform opening operations on the point cloud, determining point cloud attributes by calculating the elevation difference before and after filtering. Zhang and Chen et al. defined different window size variation methods for vegetation removal using mathematical morphology, but the window size remains difficult to determine, leading to poor vegetation removal results. Surface fitting methods first establish an initial point cloud, then iteratively select ground measurements from the original dataset to gradually approximate the surface. Axelsson proposed an adaptive irregular triangular mesh algorithm, which gradually densifies the sparse triangular mesh generated by selected seed points to obtain a surface point cloud. However, it is difficult to correctly acquire ground points around steep terrain. Zhang et al. proposed a cloth simulation algorithm, which simulates the initial point cloud as cloth laid flat on the z-coordinate inverted surface of the original data by setting rigidity values ​​for the points in the initial point cloud. The parameter selection of this algorithm is difficult, often resulting in low vegetation not being effectively filtered out or surface points in high-slope areas being incorrectly filtered out.

[0004] In conclusion, among the many vegetation removal methods proposed so far, there is still no method that can effectively preserve the surface of high-slope areas while effectively filtering out vegetation. Summary of the Invention

[0005] To address the aforementioned shortcomings, this invention provides a method for vegetation removal in high-slope areas, which can effectively filter out vegetation points while preserving the landform point set in high-slope areas.

[0006] The objective of this invention is achieved through the following technical solutions:

[0007] A method for vegetation removal in high-slope areas includes the following steps:

[0008] S1. Project all points in the point set P to be processed onto the xoy plane and perform 2D Delaunay triangulation to obtain the triangulation set T.

[0009] S2. Calculate the slope and aspect of each point in the point set P;

[0010] S3. Classify the point set P according to slope and aspect;

[0011] S4. Find a seed triangle t in the triangulation set T. seed ;

[0012] S5. Use the region growth method to gradually change the z-coordinate of the vegetation zone points;

[0013] S6. If point a and point b are connected by a line segment in the triangulation set T, then a and b are adjacent points of each other. Points whose z-coordinate values ​​have changed, as well as points whose z-coordinate values ​​have changed, are called untrusted points. Deleting untrusted points yields the point set P' after vegetation removal.

[0014] S7. Output the point set P' after removing untrusted points, which is the surface point.

[0015] Furthermore, step S2 specifically includes the following steps:

[0016] S21. Calculate the uniform normal and area of ​​each triangle in the triangulation set T;

[0017] S22. For each point p, find all the triangles connected to it, forming a set T. p Find set T p The weighted average of the product of the area of ​​the middle triangle and the normalized normal is the normal direction of point p;

[0018] S23. Let the normalized normal line at point p be (x, y, z). Then the formula for calculating the slope at point p is: Slope aspect

[0019] Furthermore, step S3 specifically includes the following steps:

[0020] S31. Find all adjacent points of each point in the input point set P, set the state of each point to unprocessed, and let k = 0;

[0021] S32. Randomly select a point p' whose state is unprocessed and add p' to class B. k Then, set the state of p' to "processed" and let n = 1;

[0022] S33. For each of the adjacent points of p' that is in an unprocessed state, determine whether p" can be merged into the same class as p'. If it can, add p" to class B. k And set the "p" state to processed;

[0023] S34. When n is greater than B k When determining the number of points in B, let n = n + 1 and proceed to the next step; otherwise, leave B.k The nth point is taken as p', and step S33 is repeated;

[0024] S35. Determine if there are still points in P that are in an unprocessed state. If not, the point classification ends; if so, let k = k + 1 and repeat step S32.

[0025] Furthermore, in step S33, can the adjacent point p" be merged into the class B to which p' belongs? k The specific judgment criteria include the following steps:

[0026] p” can be merged into B when it meets all three of the following conditions. k :

[0027] a, p” join B k B can then be made k The difference between the maximum and minimum slope aspects at all points remains less than or equal to π / 2;

[0028] b, p” added to B k B can then be made k The difference between the maximum and minimum slopes at all points remains less than or equal to π / 8;

[0029] The slope difference between c, p' and p” is less than π / 16, and the aspect difference is less than π / 8.

[0030] Furthermore, step S4 specifically includes the following steps:

[0031] S41. Find the point p with the smallest z-coordinate in the input point set P. zmin ;

[0032] S42, Let p in T zmin The triangles with vertices form a set T';

[0033] S43. Let R be the average minimum outer circle radius of the triangle in T. ave If there exists a triangle t in T' that satisfies both of the following conditions, then t can be used as a seed triangle t seed :

[0034] The difference between the maximum and minimum z-coordinates of the three vertices of points a and t is less than or equal to R. ave ;

[0035] The average of the differences between the maximum and minimum z-coordinates of the three vertices of the adjacent triangle of b and t (i.e., the triangle that shares a side with t) is less than or equal to R. ave ;

[0036] If none of the triangles in T' can simultaneously satisfy both of the above conditions, then the z-coordinate value of P will be second only to p. zminSmaller points as new p zmin Repeat step S42 until the seed triangle t can be found. seed .

[0037] Furthermore, step S5 specifically includes the following steps:

[0038] S51. Mark all points in the input point set P as unprocessed points, and mark all triangles in T as unprocessed triangles.

[0039] S52, Place the seed triangle t seed Marked as a processed triangle, t seed The three vertices are marked as processed points;

[0040] S53. The common edge between each processed and unprocessed triangle is called the growing edge. From t seed Initially, the coordinates of unprocessed triangles in T connected by each growing edge are iteratively modified until there are no more growing edges.

[0041] Furthermore, step S53 specifically includes the following steps:

[0042] S531. Sort all unprocessed triangles connected by growing edges in descending order of priority. The priority is the minimum value minus the maximum value of the z-coordinates of the three vertices of the triangle. S53 ends when there are no growing edges.

[0043] S532. Let R be the average value of the minimum circumcircle radius of the triangle in triangle T. ave The unprocessed triangle t0 is processed according to priority from high to low as follows: Let p be the vertex in t0 that is not connected to the current growing edge. t0 If t0 meets any one of the following three conditions, t0 is marked as a processed triangle, and S531 is repeated.

[0044] The difference between the maximum and minimum z-coordinates of the three vertices a and t0 is less than R. ave ,

[0045] b、p t0 The number of blocks is greater than 200.

[0046] c, p t0 The z-coordinate value is not the maximum or minimum z-coordinate value of the three vertices in t0.

[0047] Otherwise, find p t0 In triangle T, a vertex q' is called a reliable neighbor in z-coordinate if q' is a vertex of any processed triangle; find p t0 For all reliable neighboring points with z-coordinates, calculate their average z-coordinate and then use p. t0The z-value of the point is modified to this value, and then t0 is marked as the processed triangle, and S531 is repeated.

[0048] Beneficial effects:

[0049] (1) The vegetation removal method for high slope areas of the present invention obtains the triangle connection relationship of the points by projecting the input points onto the xoy plane using 2DDelaunay triangulation, and then obtains the adjacent relationship of the points; at the same time, the normal direction of the points is estimated by calculating the weighted average of the product of the area of ​​the triangle connected to the point and the normal, which reduces the situation of abrupt changes in the normal of adjacent points and effectively solves the vegetation removal error caused by the deviation of normal calculation in the uneven sampling area.

[0050] (2) The vegetation removal method for high slope areas of the present invention clusters the input points by the difference in slope and aspect between the point and its neighboring points and the difference between the maximum and minimum slope and aspect of each class. Then, by taking advantage of the fact that the slope and aspect of the hillside surface point set are relatively continuous while the slope and aspect of the vegetation area point set are discontinuous, the hillside surface point and the vegetation area point are distinguished, which effectively solves the problem that it is difficult to distinguish between hillside surface point and vegetation point.

[0051] (3) The vegetation removal method for high slope areas of the present invention uses a new seed triangle selection method to select a triangle connected to a point with a small z value and whose z value span between its vertex and the vertex of its adjacent triangle is small as a seed triangle. This effectively avoids the problem of the seed triangle being located in a vegetation area or a noise point area, which is an important prerequisite for the correct removal of vegetation points.

[0052] (4) The high-slope area vegetation removal method of the present invention modifies the coordinates of the vertices of triangles with large z-value spans through a regional growth method. However, when encountering a point with a large number of points of its class, the coordinates of the point are not modified, thereby ensuring the effective preservation of the hillside surface morphology points and the effective filtering of vegetation points.

[0053] (5) The vegetation removal method for high slope areas of the present invention, based on obtaining a triangulation network by performing 2DDelaunay triangulation on the input points, improves the point normal prediction method, seed triangle selection method, and clusters the points according to the slope and aspect, and modifies the z-coordinate of the vegetation points using the region growth method. It can effectively filter out vegetation points while preserving the surface point cloud morphology of high slope areas, and provides a better method for solving the vegetation removal problem of mountain point cloud data collected by UAVs, so that the point cloud data can be better used for terrain modeling.

[0054] (6) The earthwork calculation performed by the vegetation removal method for high-slope areas of the present invention can obtain more accurate results, solving the problem of cost control; the accuracy of landslide prediction is also higher, which can better prevent agricultural losses caused by landslides. After vegetation removal is carried out in high-slope areas by the method of the present invention, more accurate results can be obtained in earthwork calculation, which has positive significance for more accurate project cost accounting, more reasonable earthwork allocation, accelerated project progress, and improved project quality in high-slope vegetation-covered areas; in addition, it can provide more accurate ground model data for high-slope landslide-prone areas, improve the accuracy of landslide prediction in high-slope areas, and further reduce the major casualties and economic losses caused by landslides in high-slope areas. Brief Description of the Drawings

[0055] Figure 1 is a flowchart of the present invention;

[0056] Figure 2 is a flowchart for calculating the slope and aspect of a point set according to the present invention;

[0057] Figure 3 is a flowchart for clustering a point set according to slope and aspect according to the present invention;

[0058] Figure 4 is a flowchart for selecting a seed triangle according to the present invention;

[0059] Figure 5 is a flowchart for modifying the z-coordinate of a point according to the present invention;

[0060] Figure 6 is a schematic diagram of the original data of the test point set according to the present invention;

[0061] Figure 7 is a schematic diagram of the triangulation network set T obtained after Delaunay triangulation according to the present invention;

[0062] Figure 8 is a schematic diagram of estimating the normal vector of a certain point according to the present invention;

[0063] Figure 9 is a schematic diagram of classifying a point set according to slope and aspect according to the present invention, wherein a is Figure 7 a schematic diagram of processed points (represented by solid circles) in a local part of the triangulation network T shown and unprocessed adjacent points of point p' 11 (represented by hollow circles); b is a schematic diagram of the points in class B1 (represented by crosses) after three points p"1, p"2, p"3 that meet the condition among the 5 unprocessed adjacent points of p' in a are added into class B1; c is to select Figure 9 the second point in class B1 in b as point p' 11 , and p' Figure 9 12 12 ​​A schematic diagram of unprocessed adjacent points (represented by hollow circles); d represents the selected... Figure 9 The third point in class B1 of b is p' 13 A schematic diagram; e represents the processed points (represented by solid circles) and a randomly selected point p' in the region when n is greater than the number of points in class B1. 21 (The cross indicates) a schematic diagram of the first point in class B2;

[0064] Figure 10 This is a schematic diagram of the final overall classification result of the example point set in this invention;

[0065] Figure 11 This is a schematic diagram of the seed triangle selected for the present invention, where a is the point p with the smallest z-coordinate in set P. zmin (represented by solid circles) and T with p zmin 'a' is a schematic diagram of the triangle with vertices (represented by bold lines); 'b' is a seed triangle t selected from 'T' that meets the conditions. seed (Bold lines indicate)

[0066] Figure 12 This is a schematic diagram illustrating the gradual modification of the z-coordinate of vegetation zone points using the region growth method of this invention. Here, 'a' represents the seed triangle (thick line and 't') at the beginning of the region growth method. seed A diagram shows the following: (a) the processed points (represented by solid circles), the growing edges (represented by bold lines), and the unprocessed triangles connected to the current growing edge (represented by dashed lines); b is a diagram showing the priority order of the unprocessed triangles connected to the current growing edge in a (represented by 1, 2, 3), the highest priority triangle t0, and the vertices p in t0 that are not connected to the current growing edge. t0 The diagram shows the priority order of the unprocessed triangles connected by the newly grown edge (represented by numbers 1, 2, 3, 4). The triangle with the highest priority is t0, and the vertex p in t0 that is not connected to the current grown edge is... t0 p t0 The 9 adjacent points (represented by a cross shape) and p t0 A schematic diagram of two trusted adjacent points q'1 and q'2;

[0067] Figure 13 This is a schematic diagram of the final vegetation removal result of the point set in this invention;

[0068] Figure 14 This is a schematic diagram of the final vegetation removal result of the cloth algorithm for the same embodiment;

[0069] Figure 15 This is the result obtained by three-dimensional modeling of the vegetation removal results in this embodiment of the present invention;

[0070] Figure 16The result is obtained by 3D modeling of the vegetation removal results in this embodiment using the cloth algorithm. Detailed Implementation

[0071] To clearly illustrate the technical features of this solution, the following detailed implementation method, in conjunction with its accompanying drawings, will be used to describe the solution.

[0072] This implementation case focuses on point cloud data converted from images of Changshan Island collected by drones, such as... Figure 6 As shown, a total of 612,820 points were used to construct a digital terrain model using a vegetation removal method for high-slope areas proposed in this invention. The overall process of this invention is as follows: Figure 1 As shown, it includes the following steps:

[0073] S1. Project the overall point cloud data P input points onto the xoy plane and perform 2D Delaunay triangulation to obtain a triangulation set T. The triangulation set T obtained after Delaunay triangulation of the input points in this embodiment is shown below. Figure 7 As shown.

[0074] S2. Calculate the slope and aspect of each point in the point set P, as follows: Figure 2 As shown, the specific steps include:

[0075] S21. Calculate the uniform normal and area of ​​each triangle in the triangulation set T;

[0076] S22. Set a counter variable i and let i = 1;

[0077] S23. Let the i-th point in the point set P be p. i Find p i All triangles connected by points; such as Figure 8 As shown, this area is a selection. Figure 7 The diagram shows a local case of the triangular network set T, where the bold triangle is p. i Connect the points to form triangles, and form a set T. pi ;

[0078] S24, Calculate set T pi The product of the area of ​​the middle triangle and its respective normal is added together and then divided by the set T. pi The sum of the areas of the middle triangles gives us p. i The direction of the normal to the point;

[0079] S25, Let p i The normalized normal of a point is (x i ,y i ,z i ), then p i The formula for calculating the point slope is: Slope aspect

[0080] S26. Let i = i + 1; if i is less than or equal to the total number of points in the point set P, repeat step S23; otherwise, the slope and aspect calculation for each point is complete.

[0081] S3. Classify the point set P according to slope and aspect, as follows: Figure 3 As shown, the specific steps include:

[0082] S31. Set the state of each point in the input point set P to unprocessed, and set a counter variable k, let k = 0;

[0083] S32. Randomly select one point p' whose state is unprocessed. k1 , will p' k1 Add class B k In the middle, p' k1 Set the state to "processed", and set a counter variable n, letting n = 1; for example Figure 9 As shown in a, this region is a selection. Figure 7 The diagram shows a local case of the triangular network T, where k = 1. Processed points are represented by solid circles, p' 11 Five unprocessed points among the adjacent points are represented by hollow circles;

[0084] S33. Set a counter variable j and let j = 1;

[0085] S34, B k The nth point is p' kn If point a and point b are connected by a line segment in a triangulation network T, then a and b are adjacent points of each other. Let p' kn The number of unprocessed points among the adjacent points is num;

[0086] S35, Let p' kn The j-th unprocessed neighbor is p” j , for p” j Perform the following condition judgment:

[0087] a、p” j Join B k B can then be made k The difference between the maximum and minimum slope aspects at all points remains less than or equal to π / 2;

[0088] b、p” j Join B k B can then be made k The difference between the maximum and minimum slopes at all points remains less than or equal to π / 8;

[0089] c, p' kn and p'' j have a gradient difference less than π / 16 and an aspect difference less than π / 8;

[0090] If p'' j satisfies all three conditions above at the same time, p'' j is added to B k , and the state of p'' j is set to processed; if the conditions are not all satisfied, the next step is executed directly; Figure 9 After each of the 5 unprocessed neighboring points of point p' in a 11 is judged through the above process, if 3 points p''1, p''2, p''3 satisfy all three conditions above at the same time, they are added to class B1, as shown in Figure 9 b, at this time there are 4 points in B1, all represented by crosses;

[0091] S36: Let j = j+1. If j ≤ num at this time, re-perform step S35; otherwise, let n = n+1 and proceed to the next step;

[0092] S37: If n is less than or equal to the number of points in B k , re-perform step S34; otherwise proceed to the next step; after the steps shown in b Figure 9 are completed, n = 1+1 = 2, which is less than or equal to the number 4 of points in B1, then the 2nd point in B1 is selected as p' 12 , as shown in Figure 9 c, re-perform step S34; suppose when Figure 9 the three unprocessed neighboring points of point p' represented by a hollow circle in c 12 are judged in step S35, none of them can satisfy all three conditions in step S35 at the same time, at this time n = 2+1 = 3, which is still less than or equal to the number 4 of points in B1, then the 3rd point in B1 is selected as p' 13 , as shown in Figure 9 d, re-perform step S34; suppose that no unprocessed neighboring points satisfying all three conditions in step S35 at the same time are found for all 4 points in B1 in b Figure 9 , at this time n = 4+1 = 5 which is greater than the number of points in B1, then proceed to the next step;

[0093] S38: Judge whether there is a point with an unprocessed state in P, if yes, let k = k+1 and re-perform step S32, otherwise the classification ends; as shown in Figure 9 e, the unlabeled points in the figure are all unprocessed points, let k = 1+1 = 2, randomly select one point p' 21 , represented by a cross, re-perform step S32, and the classification ends when all points are processed.

[0094] The final classification result of the point set in this implementation is as follows: Figure 10 As shown, different classes are represented by different shades of gray.

[0095] S4. Find a seed triangle t in T. seed The process is as follows Figure 4 As shown, the specific steps include:

[0096] S41. Find the point p with the smallest z-coordinate in the input point set P. zmin ;like Figure 11 As shown in figure a, the solid circle marks the point p with the smallest z-coordinate in the entire point set. zmin ;

[0097] S42, Let p in T zmin The triangles with vertices form a set T', such as Figure 11 The triangle represented by the bold line in a is shown;

[0098] S43. Let R be the average minimum outer circle radius of the triangle in T. ave If there exists a triangle t in T' that satisfies both of the following conditions, then t can be used as a seed triangle t seed :

[0099] The difference between the maximum and minimum z-coordinates of the three vertices of points a and t is less than or equal to R. ave ;

[0100] The average of the differences between the maximum and minimum z-coordinates of the three vertices of the adjacent triangle of b and t (i.e., the triangle that shares a side with t) is less than or equal to R. ave ;

[0101] like Figure 11 If the triangle t highlighted in bold in b meets both of the above conditions, then this triangle is used as the seed triangle t. seed ;

[0102] Otherwise, the z-coordinate value in P is second only to p. zmin The point as the new p zmin Repeat step S42 until the seed triangle t is found. seed .

[0103] S5. Gradually change the z-coordinate of the vegetation zone points using the region growth method, as follows: Figure 5 As shown, the specific steps include:

[0104] S51. Mark all points in the input point set P as unprocessed points, and mark all triangles in T as unprocessed triangles.

[0105] S52, Place the seed triangle t seed Marked as processed triangles, such as Figure 12 The bold line in 'a' represents the triangle; t seed The three vertices are marked as processed points, such as... Figure 12 The solid circles in 'a' represent points; the common edge between each processed and unprocessed triangle is called the growth edge, such as... Figure 12 The bolded border in 'a';

[0106] S53. Find all unprocessed triangles connected by growing edges, and form a queue Q based on the minimum and maximum z-coordinates of the three vertices of the triangle; from t seed Initially, the unprocessed triangles connected by the currently growing edges are in Figure 12 In section 'a', the dashed lines indicate the priority order. Figure 12 b is labeled with numbers;

[0107] S54. Determine if Q is empty. If it is empty, end step 5; otherwise, let the triangle with the highest priority in Q be t0, and let p be the vertex in t0 that is not connected to the current growing edge. t0 ,like Figure 12 As shown in b;

[0108] S55. Let the average value of the minimum circumcircle radius of triangle T be R. ave Determine whether t0 meets any one of the following three conditions:

[0109] The difference between the maximum and minimum z-coordinates of the three vertices a and t0 is less than R. ave ;

[0110] b、p t0 The category score is greater than 200;

[0111] c, p t0 The z-coordinate value is not the maximum or minimum z-coordinate value of the three vertices in t0;

[0112] If so, mark t0 as a processed triangle and repeat step S53; otherwise, proceed to step S56; for example... Figure 12 b. If t0 is the triangle with the highest priority and satisfies condition a, then mark t0 as a processed triangle and repeat step S53. Figure 12 In step c, the priority order of the unprocessed triangles connected by the newly grown edge is marked with numbers. The triangle with the highest priority is t0. The three conditions in step S55 are judged. At this time, none of the three conditions of t0 are met, so step S56 is performed.

[0113] S56. Set up counting variables l and m, let l = 1, m = 0; set up an accumulation variable sum, let sum = 0;

[0114] S57, Let q' l For pt0 In the l-th adjacent node of T, if q' l Let q' be any vertex of a processed triangle. l For p t0 Let m = m + 1, sum = sum + q' be the trusted adjacent points. l If the z-coordinate value is p, then l = l + 1; otherwise, simply set l = l + 1; repeat step S57 until l is greater than p. t0 The number of adjacent nodes in T is determined in step S58; such as... Figure 12 As shown in c, p t0 There are 9 adjacent vertices in T, represented by cross shapes. After performing a loop check, it is found that only two of them are p. t0 Trusted neighboring points are labeled with “q’1” and “q’2” respectively, and step S58 is performed;

[0115] S58, p t0 Modify the z-coordinate of the point to sum / m, mark t0 as the processed triangle, and proceed to step S53; for example... Figure 12 In c, find the two p's, q'1 and q'2. t0 After finding the reliable neighbor points, let the z-coordinate of q'1 be 38.1793 and the z-coordinate of q'2 be 37.7783. Then the average z-coordinate is (38.1793 + 37.7783) / 2 = 37.9788. Then p... t0 After modifying the z-value of the point to 37.9788, mark t0 as the processed triangle and repeat step S53.

[0116] S6. Points whose z-coordinate values ​​have changed, as well as points whose z-coordinate values ​​have changed for all their adjacent points, are called untrusted points; deleting untrusted points yields the point set P' after vegetation removal.

[0117] like Figure 13 The image shows the final vegetation removal results for the point set in this embodiment. Figure 14 The results of the cloth removal algorithm proposed by Zhang et al. on the removal of vegetation points in this embodiment show that although the cloth removal algorithm removes vegetation points, it also removes hillside points incorrectly. The vegetation removal algorithm proposed in this invention can effectively filter out vegetation points while preserving the surface point cloud in high-slope areas.

[0118] S7. Output the point cloud set P' after removing untrusted points, which is the surface point.

[0119] Figure 15 The results of the three-dimensional modeling of the vegetation removal results in this embodiment of the present invention show that the terrain of the high slope area was effectively preserved. Figure 16The 3D modeling results of vegetation removal in this embodiment, obtained using the cloth-laying algorithm, show that topographical features in high-slope areas were largely incorrectly removed. After vegetation removal in high-slope areas using the method of this invention, earthwork volume calculations will yield more accurate results. This is significant for more accurate cost accounting in high-slope vegetated areas, more rational earthwork allocation, faster project progress, and improved project quality. Furthermore, it provides more accurate ground model data for high-slope landslide-prone areas, improving the accuracy of landslide prediction in these areas and thus reducing significant casualties and economic losses caused by landslides.

[0120] For those skilled in the art, any changes, modifications, substitutions, and variations made to the embodiments without departing from the principles and spirit of the present invention, based on the teachings of the present invention, still fall within the protection scope of the present invention.

Claims

1. A method for vegetation removal in high-slope areas, characterized in that, Includes the following steps: S1. Project all points in the point set P to be processed onto the xoy plane and perform 2D Delaunay triangulation to obtain the triangulation set T. S2. Calculate the slope and aspect of each point in the point set P; S3. Classify the point set P according to slope and aspect; S4. Find a seed triangle t in the triangulation set T. seed Specifically, the steps include the following: S41. Find the point p with the smallest z-coordinate in the input point set P. zmin ; S42, Let p in T zmin The triangles with vertices form a set T'; S43. Let R be the average minimum outer circle radius of the triangle in T. ave If there exists a triangle t in T' that satisfies both of the following conditions, then t can be used as a seed triangle t seed : The difference between the maximum and minimum z-coordinates of the three vertices of points a and t is less than or equal to R. ave ; The average of the differences between the maximum and minimum z-coordinates of the three vertices of the adjacent triangles of b and t is less than or equal to R. ave ; If none of the triangles in T' can simultaneously satisfy both of the above conditions, then the z-coordinate value of P will be second only to p. zmin Smaller points as new p zmin Repeat step S42 until the seed triangle t can be found. seed ; S5. Gradually change the z-coordinates of the vegetation zone points using the region growth method; specifically including the following steps: S51. Mark all points in the input point set P as unprocessed points, and mark all triangles in T as unprocessed triangles. S52, Place the seed triangle t seed Marked as a processed triangle, t seed The three vertices are marked as processed points; S53. The common edge between each processed and unprocessed triangle is called the growing edge. From t seed Initially, the coordinates of the unprocessed triangles in T connected by each growing edge are iteratively modified until there are no more growing edges. S6. If point a and point b are connected by a line segment in the triangulation set T, then a and b are adjacent points of each other. Points whose z-coordinate values ​​have changed, as well as points whose z-coordinate values ​​have changed, are called untrusted points. Deleting untrusted points yields the point set P' after vegetation removal. S7. Output the point set P' after removing untrusted points, which is the surface point.

2. The method for vegetation removal in high-slope areas according to claim 1, characterized in that, Step S2 specifically includes the following steps: S21: Calculate the uniform normal and area of ​​each triangle in the triangulation set T; S22: For each point p, find all the triangles connected to it, forming a set T. p Find set T p The weighted average of the product of the area of ​​the middle triangle and the normalized normal is the normal direction of point p; S23: Let the normalized normal line at point p be (x, y, z), then the formula for calculating the slope at point p is: The slope direction is .

3. The method for vegetation removal in high-slope areas according to claim 1, characterized in that, Step S3 specifically includes the following steps: S31. Find all adjacent points of each point in the input point set P, set the state of each point to unprocessed, and let k=0; S32: Randomly select a point p' whose state is unprocessed, and add p' to class B. k Then, set the state of p' to "processed" and let n=1; S33. Let p'' be the unprocessed point among the adjacent points of p'. Check each point in the class whether p'' can be merged with p'. If it can, add p'' to class B. k And set the state of p'' to processed; S34. When n is greater than B k When determining the number of points in B, let n = n + 1 and proceed to the next step; otherwise, leave B. k The nth point is taken as p', and step S33 is repeated; S35. Determine if there are still points in P that are in an unprocessed state. If not, the point classification ends; if so, let k = k + 1 and repeat step S32.

4. The method for vegetation removal in high-slope areas according to claim 3, characterized in that, In step S33, can p'' be merged into the class B to which p' belongs? k The specific judgment criteria include the following steps: p'' can be merged into B when it simultaneously meets the following three conditions. k : a、p''Add B k B can then be made k The difference between the maximum and minimum slope aspects at all points remains less than or equal to π / 2; b, p'' Add B k B can then be made k The difference between the maximum and minimum slopes at all points remains less than or equal to π / 8; The slope difference between c, p', and p'' is less than π / 16, and the aspect difference is less than π / 8.

5. The method for vegetation removal in high-slope areas according to claim 1, characterized in that, Step S53 specifically includes the following steps: S531. Sort all unprocessed triangles connected by growing edges in descending order of priority. The priority is the minimum value minus the maximum value of the z-coordinates of the three vertices of the triangle. S53 ends when there are no growing edges. S532. Let R be the average value of the minimum circumcircle radius of the triangle in triangle T. ave The unprocessed triangle t0 is processed according to priority from high to low as follows: Let p be the vertex in t0 that is not connected to the current growing edge. t0 If t0 meets any one of the following three conditions, t0 is marked as a processed triangle, and S531 is repeated. The difference between the maximum and minimum z-coordinates of the three vertices a and t0 is less than R. ave , b、p t0 The number of blocks is greater than 200. c, p t0 The z-coordinate value is not the maximum or minimum z-coordinate value of the three vertices in t0. Otherwise, find p t0 In triangle T, a vertex q' is called a reliable neighbor in z-coordinate if q' is a vertex of any processed triangle; find p t0 For all reliable neighboring points with z-coordinates, calculate their average z-coordinate and then use p. t0 The z-value of the point is modified to this value, and then t0 is marked as the processed triangle, and S531 is repeated.

Images