An elevation calculation method and system based on delaunay triangulation
By triangulating and projecting the elevation point set, combined with KD tree search and parametric equations, the time-consuming problem of elevation calculation in large-scale Delaunay triangulation is solved, and fast and accurate elevation calculation is achieved.
Patent Information
- Application Number
- CN202411343943.4
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2024-09-25
- Publication Date
- 2025-10-10
- Estimated Expiration
- 2044-09-25
AI Technical Summary
Existing technologies take a long time to calculate the elevation of a specified point in a large-scale Delaunay triangulation. The linear search method and the jump and walk method are inefficient when dealing with uneven point sets or poor starting points.
The Delaunay triangulation is generated by obtaining the elevation point set and performing triangulation. The triangulation is projected onto a two-dimensional plane and split into an undirected graph. The nearest point is searched using a KD tree, and the elevation value is calculated using the triangle parametric equation.
It achieves the rapid and accurate calculation of elevation values in large-scale scenes, reduces the computational time complexity and improves efficiency.
Smart Images

Figure CN119359944B_ABST
Abstract
Description
Technical Field
[0001] The present application relates to the field of elevation estimation, and in particular to an elevation calculation method and system based on Delaunay triangulation. Background Art
[0002] In CAD drawing operations, it is necessary to generate a Delaunay triangulation based on specified points, then determine the triangle where the specified point is located, and then calculate the elevation of the specified point in the triangle. In this process, determining the triangle where the specified point is located is very time-consuming.
[0003] Existing technologies generally use linear search and jump and walk methods to estimate the elevation of a specified point. The linear search method needs to traverse every triangle in the triangulation network and check whether the specified point is inside the triangle. This approach is time-consuming and not suitable for scenes with large triangulation networks. The jump and walk method needs to randomly select a starting triangle and determine whether the specified point is inside the triangle in turn. This approach is related to the distribution of Delaunay triangles. If the point set is unevenly distributed or the starting point is not well selected, then a large number of triangles need to be visited. Therefore, when processing scenes with large triangulation networks, it takes relatively long time.
[0004] Application Content
[0005] This application provides an elevation calculation method based on Delaunay triangulation to quickly and accurately calculate elevation values in large-scale scenes.
[0006] In a first aspect, the present application provides an elevation calculation method based on a Delaunay triangulation, comprising:
[0007] Acquire a plurality of elevation point sets to be calculated, and triangulate the elevation point sets to obtain a plurality of Delaunay triangulations;
[0008] Projecting the Delaunay triangulation onto a two-dimensional plane to obtain a triangle set, and splitting each triangle in the triangle set to obtain a triangle undirected graph, wherein the undirected graph includes a vertex set and an edge set;
[0009] Obtain a point to be queried, and use a KD tree to search for the first point in the undirected graph that is closest to the point to be queried;
[0010] Based on the point to be queried and the first point position, a triangle containing the point to be queried is determined, and based on a parametric equation of a plane where the triangle is located, an elevation value corresponding to the point to be queried is obtained.
[0011] The embodiment of the application can quickly generate a Delaunay triangular network, facilitate subsequent determination of a triangle in which a query point is located, by obtaining a plurality of sets of elevation points to be calculated, and performing triangulation on the sets of elevation points to obtain a plurality of Delaunay triangular networks; can project the Delaunay triangular network into a two-dimensional plane to obtain a set of triangles, and split each triangle in the set of triangles to obtain a triangle-free graph, which can project a triangular network in a three-dimensional space into a two-dimensional plane, facilitate splitting of each vertex set and edge set in each triangle, and further facilitate subsequent determination of a triangle containing the query point and elevation calculation; can quickly obtain a first point closest to the query point by obtaining the query point, and searching for the first point closest to the query point in the graph using a KD tree, which facilitates subsequent determination of a triangle containing the query point and elevation calculation; and can determine a triangle containing the query point based on the query point and the first point, and obtain an elevation value corresponding to the query point based on a parametric equation of a plane in which the triangle is located, which can calculate the elevation value in a large-scale scene quickly and accurately.
[0012] Further, the obtaining of the plurality of sets of elevation points to be calculated specifically comprises:
[0013] Collecting a plurality of elevation points by using an RTK or a total station;
[0014] Importing the plurality of elevation points into CAD to obtain the plurality of sets of elevation points to be calculated.
[0015] In this way, a plurality of elevation points are collected in advance by using an RTK or a total station, and the plurality of sets of elevation points to be calculated are obtained by importing the plurality of elevation points into CAD, which facilitates subsequent construction of a Delaunay triangular network.
[0016] Further, the searching for the first point closest to the query point in the graph using the KD tree specifically comprises:
[0017] Comparing the query point with a root node, if a value of the query point in a current dimension is less than a value of the root node and the root node has no sub-tree, determining that the root node is a second point, and calculating a first distance between the second point and the query point;
[0018] Searching for and calculating a second distance between a third point in a parent node corresponding to the root node and the query point, and constantly updating the second point based on the second distance until there is no smaller distance, and determining that the second point is the first point closest to the query point.
[0019] In this way, by using the KD tree to search the first point in the undirected graph that is closest to the point to be queried, the first point that is closest to the point to be queried can be quickly obtained, which facilitates the subsequent determination of the triangle containing the point to be queried and the elevation calculation.
[0020] Furthermore, the determining of a triangle containing the to-be-queried point based on the to-be-queried point and the first point position is specifically as follows:
[0021] Taking the first point as the starting point, draw a first ray through the point to be queried;
[0022] Extracting a plurality of candidate points adjacent to the first point from the undirected graph, using the first candidate point and the second candidate point as the second ray and the third ray respectively to pass through the first point, and calculating the angles between the first ray and the second ray and the first ray, respectively, to obtain a first angle and a second angle;
[0023] Comparing the first angle with the second angle, and if the second angle is smaller than the first angle, determining the second candidate point corresponding to the second angle as the optimal candidate point;
[0024] The optimal candidate point is used as a starting point to continue drawing rays through the point to be queried until a repeated point appears, and a triangle containing the point to be queried is obtained based on the rays between the repeated points.
[0025] In this way, by determining a triangle containing a point to be queried and performing calculations based on the parametric equation of the determined triangle containing the point to be queried, the elevation value in a large-scale scene can be calculated quickly and accurately.
[0026] Furthermore, the parametric equation of the plane where the triangle is located is specifically:
[0027] A(x-x1)+B(y-y1)+C(z-z1)=0;
[0028] In the formula, A is ((y2-y1)(z3-z1)-(z2-z1)(y3-y1)), B is ((x2-x1)(z3-z1)-(z2-z1)(x3-x1)), C is ((x2-x1)(y3-y1)-(y2-y1)(x3-x1)), and x1, y1 and z1 are the positions of the first point.
[0029] In this way, by performing calculation based on the determined parametric equation of the triangle containing the point to be queried, the elevation value in a large-scale scene can be calculated quickly and accurately.
[0030] In a second aspect, the present application provides an elevation calculation system based on a Delaunay triangulation, comprising: a first acquisition module, a first processing module, a second acquisition module, and a second processing module;
[0031] The first acquisition module is used to acquire a plurality of elevation point sets to be calculated, and triangulate the elevation point sets to obtain a plurality of Delaunay triangulations;
[0032] The first processing module is configured to project the Delaunay triangulation onto a two-dimensional plane to obtain a triangle set, and split each triangle in the triangle set to obtain a triangle undirected graph, wherein the undirected graph includes a vertex set and an edge set;
[0033] The second acquisition module is used to obtain a point to be queried, and use a KD tree to search the undirected graph for a first point that is closest to the point to be queried;
[0034] The second processing module is used to determine a triangle containing the point to be queried based on the point to be queried and the first point position, and obtain an elevation value corresponding to the point to be queried based on a parametric equation of the plane where the triangle is located.
[0035] The embodiment of the present application obtains several elevation point sets to be calculated and triangulates the elevation point sets to obtain several Delaunay triangulations, which can quickly generate a Delaunay triangulation, which is convenient for subsequently determining the triangle where the query point is located; by projecting the Delaunay triangulation onto a two-dimensional plane, a triangle set is obtained, and each triangle in the triangle set is split to obtain a triangle undirected graph, which can project the three-dimensional triangulation onto a two-dimensional plane, and conveniently split each vertex set and edge set in each triangle, which is helpful for subsequently determining the triangle containing the query point. Triangle and elevation calculation; by obtaining the point to be queried, and using the KD tree to search the first point position in the undirected graph that is closest to the point to be queried, the first point position closest to the point to be queried can be quickly obtained, which is convenient for the subsequent determination of the triangle containing the point to be queried and the elevation calculation; based on the point to be queried and the first point position, the triangle containing the point to be queried is determined, and based on the parametric equation of the plane where the triangle is located, the elevation value corresponding to the point to be queried is obtained, and calculation can be performed based on the parametric equation of the triangle containing the point to be queried, and the elevation value in large-scale scenes can be quickly and accurately calculated.
[0036] Furthermore, the first acquisition module includes: a collection unit and an import unit;
[0037] The acquisition unit is used to acquire a plurality of elevation points using RTK or a total station;
[0038] The importing unit is used to import the plurality of elevation points into CAD to obtain a plurality of elevation point sets to be calculated.
[0039] Furthermore, the second acquisition module includes: a comparison unit and an update unit;
[0040] The comparison unit is configured to compare the query point with a root node, and if the value of the query point at the current dimension is less than the value of the root node and the root node has no subtree, determine that the root node is a second point position, and calculate a first distance between the second point position and the query point;
[0041] The updating unit is used to search for and calculate the second distance between the third point position in the parent node corresponding to the root node and the point to be queried, and continuously update the second point position based on the second distance until there is no smaller distance, and determine that the second point position is the first point position closest to the point to be queried.
[0042] Furthermore, the second processing module includes: a third processing unit, a calculation unit, a second comparison unit and a construction unit;
[0043] The third processing unit is configured to generate a first ray starting from the first point and passing through the point to be queried;
[0044] The calculation unit is configured to extract a plurality of candidate points adjacent to the first point from the undirected graph, use the first candidate point and the second candidate point as a second ray and a third ray respectively to pass through the first point, and calculate the angles between the first ray and the second ray and the first ray, respectively, to obtain a first angle and a second angle;
[0045] The second comparison unit is configured to compare the first angle with the second angle, and if the second angle is smaller than the first angle, determine that a second candidate point corresponding to the second angle is an optimal candidate point;
[0046] The construction unit is configured to use the optimal candidate point as a starting point to continue drawing rays through the point to be queried until a repeated point appears, and obtain a triangle containing the point to be queried based on the rays between the repeated points.
[0047] Furthermore, the parametric equation of the plane where the triangle is located is specifically:
[0048] A(x-x1)+B(y-y1)+C(z-z1)=0;
[0049] In the formula, A is ((y2-y1)(z3-z1)-(z2-z1)(y3-y1)), B is ((x2-x1)(z3-z1)-(z2-z1)(x3-x1)), C is ((x2-x1)(y3-y1)-(y2-y1)(x3-x1)), and x1, y1 and z1 are the positions of the first point. BRIEF DESCRIPTION OF THE DRAWINGS
[0050] Figure 1 This is a flow chart of an embodiment of a method for calculating elevation based on a Delaunay triangulation provided by the present application;
[0051] Figure 2 It is a schematic diagram of determining the optimal candidate point provided by this application;
[0052] Figure 3 This is a structural diagram of an embodiment of an elevation calculation system based on Delaunay triangulation provided by the present application. DETAILED DESCRIPTION
[0053] The following will be combined with the drawings in the embodiments of this application to clearly and completely describe the technical solutions in the embodiments of this application. Obviously, the embodiments described are only part of the embodiments of this application, not all of the embodiments. Based on the embodiments in this application, all other embodiments obtained by ordinary technicians in this field without making creative efforts are within the scope of protection of this application.
[0054] It should be understood that the step numbers used herein are only for convenience of description and are not intended to limit the order in which the steps are executed.
[0055] It should be understood that the terms used in the present specification are only for the purpose of describing specific embodiments and are not intended to limit the present invention. As used in the present specification and the appended claims, the singular forms "a", "an" and "the" are intended to include the plural forms unless the context clearly indicates otherwise.
[0056] The terms “include” and “comprising” indicate the presence of described features, integers, steps, operations, elements and / or components, but do not preclude the presence or addition of one or more other features, integers, steps, operations, elements, components and / or groups thereof.
[0057] The term "and / or" refers to and includes any and all possible combinations of one or more of the associated listed items.
[0058] In CAD drawing, after constructing a Delaunay triangulation, it is necessary to locate the triangle to which a specified point belongs in order to calculate its elevation, a time-consuming process. Existing elevation estimation techniques include linear search and jump-and-walk methods. The linear search method requires traversing all triangles, making it inefficient and unsuitable for large-scale triangulations. The jump-and-walk method uses a random starting point and relies on a Delaunay distribution. This can be time-consuming when the point set is uneven or the starting point is poorly chosen.
[0059] The Delaunay triangulation (also known as the Delaunay triangulation) is a special triangulation that is constructed according to a series of strict geometric criteria and consists of a set of adjacent and non-overlapping triangles.
[0060] Example 1
[0061] Please refer to Figure 1 , Figure 1 This is a flow chart of an embodiment of a method for calculating elevation based on a Delaunay triangulation provided by the present application, including steps S1 to S4;
[0062] Step S1, obtaining a plurality of elevation point sets to be calculated, and triangulating the elevation point sets to obtain a plurality of Delaunay triangulations;
[0063] It can be understood that several elevation points can be collected using RTK or a total station; and by importing several of the elevation points into CAD, several sets of elevation points to be calculated can be obtained, where the elevation points are point information in a three-dimensional space, including x-axis coordinates, y-axis coordinates, and z-axis coordinates.
[0064] In this way, several elevation points are collected in advance through RTK or total station, and imported into CAD to obtain several elevation point sets to be calculated, which facilitates the subsequent construction of Delaunay triangulation.
[0065] After obtaining several elevation point sets, it is necessary to preprocess the elevation point sets and select the Delaunay triangulation algorithm to triangulate the preprocessed elevation point sets to obtain several Delaunay triangulation networks; wherein the triangulation algorithm includes but is not limited to the Bowyer-Watson algorithm, the Lawson local optimization algorithm (LOP), etc., for example and not limitation, taking the Bowyer-Watson algorithm for triangulation as an example, first, it is necessary to construct a sufficiently large super triangle as the initial triangular network, the vertices of which contain all the points in the elevation point set; secondly, the points in the elevation point set are inserted into the triangular mesh one by one, for each newly inserted point, the influencing triangle containing the point is found, and the common edges of the influencing triangles are deleted, and the new point is connected to the vertices of the influencing triangle to form three new triangles. Again, each newly formed triangle is tested for an empty circumcircle to determine whether it meets the Delaunay condition. If the circumcircle of a triangle contains other elevation points, it means that the triangle does not meet the Delaunay condition. At this time, the triangle needs to be optimized by exchanging diagonals until all triangles meet the empty circumcircle condition; then, repeat the above steps of inserting points, constructing triangles, and detecting and optimizing empty circumcircles until all elevation points are inserted into the triangular mesh; finally, delete the triangles related to the super triangle to obtain the final Delaunay triangulation.
[0066] It should be noted that data preprocessing methods include but are not limited to removing duplicate points, processing noise points and filling missing data, and these data preprocessing methods can all be achieved through existing technologies, but data preprocessing methods are not the focus of this application, so they will not be expanded here.
[0067] It should be noted that the Bowyer-Watson algorithm is an incremental algorithm for generating Delaunay triangulations. Delaunay triangulation is a special triangulation method that generates a triangulated mesh that meets the Delaunay conditions by gradually inserting points and rebuilding the triangulated mesh until all points are added. The Lawson Local Optimization Procedure (LOP) is an algorithm used to optimize the triangulated mesh structure during the triangulation process to ensure that it meets the Delaunay triangulation properties. It optimizes triangles that do not meet the Delaunay conditions by exchanging diagonals.
[0068] Step S2: Projecting the Delaunay triangulation onto a two-dimensional plane to obtain a triangle set, and splitting each triangle in the triangle set to obtain a triangle undirected graph, wherein the undirected graph includes a vertex set and an edge set;
[0069] It is understandable that when a Delaunay triangulation is obtained, it is necessary to simplify the Delaunay triangulation into a two-dimensional plane for further analysis. Specifically, first, it is necessary to predetermine a suitable projection method to better preserve the structural characteristics of the three-dimensional Delaunay triangulation; second, each point in the Delaunay triangulation is projected into a two-dimensional plane according to the selected projection method. During the projection process, the coordinates of each point on the two-dimensional plane need to be calculated, and the relative position relationship between the points needs to be preserved; finally, the two-dimensional plane is re-triangulated based on the projected point set to obtain a triangle set of the two-dimensional plane. The triangulation method has been detailed in step S1, so it will not be repeated here.
[0070] When a triangle set on a two-dimensional plane is obtained, it is necessary to extract all unique vertices from the triangle set and collect all vertices to form a vertex set V, where V = {v1, v2, ..., v n}, where v i Next, we traverse each triangle in the triangle set and add the three edges connected to each triangle to the edge set E. Each edge consists of two vertices, indicating that the two vertices are adjacent in the triangle. For each triangle Δ ijk , whose vertices are v i , v j , v k , its three sides (v i ,v j ),(v i ,v k ),(v j ,v k ) is added to the set of non-phase edges E, that is, (v i ,v j )∈E, if and only if (v j ,v i )∈E. Finally, use the vertex set V and edge set E to construct an undirected graph G=(V,E).
[0071] It should be noted that the projection method includes but is not limited to orthogonal projection, perspective projection and isometric projection, etc. Projecting the Delaunay triangulation into a two-dimensional plane is not the focus of this application and will not be expanded here.
[0072] It's important to note that an undirected graph is a concept in graph theory, where edges have no directionality. That is, an edge connecting two vertices can be considered bidirectional. In an undirected graph, vertices and edges together constitute the basic elements of the graph, and the edge set generally defines the connection between vertices. A triangular undirected graph can be understood as a special type of undirected graph, consisting of multiple triangles formed by three vertices connected end-to-end by three edges.
[0073] Step S3: obtaining a point to be queried, and using a KD tree to search for the first point in the undirected graph that is closest to the point to be queried;
[0074] It can be understood that after obtaining a triangular undirected graph, a KD tree can be established based on the vertex set V to speed up the search efficiency of the nearest point. The process of constructing the KD tree is as follows: (1) arbitrarily select a vertex from the vertex set V as the root node; (2) In each layer of the KD tree, it is necessary to cyclically select different dimensions as the splitting dimensions; (3) For the selected splitting dimension, calculate the median of all vertices on the dimension, and use the median as the root node to divide the remaining points into two subtrees on the left and right; (4) Repeat steps (2) and (3) for the vertices of the left and right subtrees respectively, and recursively construct the left subtree and the right subtree until the subtree is empty.
[0075] After the KD tree is constructed, the query point p can be obtained i , and use KD tree to search the undirected graph G with p i Euclidean distance to the nearest point v s The specific query process is as follows: First, select the point p to be queried i The first dimension is the current dimension d, and the query point is compared with the same dimension d t The root node p c Compare, if the query point p i The value at the current dimension is less than the root node p c When the value of , then enter the root node p c Continue searching the left subtree of the root node p; otherwise, go to the root node p c Continue searching the right subtree of the node until the root node p c If there is no subtree, then determine the root node p c The first distance d between the second point and the query point is calculated as c ; Secondly, search and calculate the root node p c The corresponding parent node p t The second distance d with the smallest distance between the third point and the query point a , and based on the second distance d aThe second point is continuously updated until there is no smaller distance, and the second point is determined to be the first point v closest to the query point. s For example and not limitation, if the second distance d a <First distance d c , then calculate the parent node p t and the query point p i The distance d t , if d t <d c , then its parent node p t For the node with the closest distance, update the second point, i.e. p c =p t , d c =d t , recursively search for p t The left and right subtrees of are repeatedly recursively searched to determine the point closest to the query point.
[0076] In this way, by using the KD tree to search the first point in the undirected graph that is closest to the point to be queried, the first point that is closest to the point to be queried can be quickly obtained, which facilitates the subsequent determination of the triangle containing the point to be queried and the elevation calculation.
[0077] Step S4: Based on the point to be queried and the first point position, determine a triangle containing the point to be queried, and obtain an elevation value corresponding to the point to be queried based on a parametric equation of the plane where the triangle is located.
[0078] It is understandable that after determining the first point closest to the query point, the next path node v can be searched. p (also referred to as the optimal candidate point below) to determine the triangle containing the query point. Specifically, first, take the first point v s Draw the first ray l as the starting point s Passing the query point p i ; Secondly, arbitrarily extract the first point v from the undirected graph G s The adjacent vertex set V includes several candidate points. The first candidate point v is arbitrarily selected from the vertex set V. p , and take the first candidate point v p Draw the second ray l as the starting point p After the first point v s , where the initial optimal angle d is set p =2π, traverse the vertex set V and select the second candidate point v k , and take the second candidate point v k Draw the third ray l as the starting point k After the first point vs , and calculate the second ray l p and the third ray l k With the first ray l s The angle between the two angles is d. p and the second angle d k ; Again, the first angle d p and the second angle d k For comparison, if the second angle d k smaller than the first angle d p , then determine the second angle d k The corresponding second candidate point v k is the optimal candidate point, even if v p =v k , and make d p =d k , otherwise keep v p is the optimal candidate point; finally, the optimal candidate point v p Push it into the path queue, and continue to make rays passing through the query point as the starting point, repeat the above steps until a duplicate point appears in the path queue, and take the last 3 points of the path queue, that is, form a triangle loop based on the rays between the duplicate points, and finally get the triangle containing the query point p i The triangle ΔS, where the schematic diagram of determining the optimal candidate point is as follows Figure 2 shown.
[0079] In this way, by determining a triangle containing a point to be queried and performing calculations based on the parametric equation of the determined triangle containing the point to be queried, the elevation value in a large-scale scene can be calculated quickly and accurately.
[0080] When the query point p is determined to be included i The triangle ΔS can get the query point p i Specifically, find the parametric equation f(x, y, z) of the plane where ΔS is located. The process of finding the parametric equation is as follows: (1) Let the three points of triangle ΔS be a(x1, y1, z1), b(x2, y2, z2), and c(x3, y3, z3), and find the two vectors: The cross product of the two vectors can be used to obtain the normal vector of the plane where ΔS is located, and the plane normal vector is brought into the query point p i The x coordinate i ,y i Get the elevation value z i .in, right Expand the determinant and get: according to The plane parametric equation f(x,y,z) can be written as:
[0081] A(x-x1)+B(y-y1)+C(z-z1)=0;
[0082] In the formula, A is ((y2-y1)(z3-z1)-(z2-z1)(y3-y1)), B is ((x2-x1)(z3-z1)-(z2-z1)(x3-x1)), C is ((x2-x1)(y3-y1)-(y2-y1)(x3-x1)), and x1, y1 and z1 are the positions of the first point.
[0083] Therefore, the query point p i (x i ,y i ) into f(x,y,z) to get the corresponding elevation value z i By repeating this method, the elevation values corresponding to several candidate points in the vertex set V can be calculated.
[0084] In this way, by performing calculation based on the determined parametric equation of the triangle containing the point to be queried, the elevation value in a large-scale scene can be calculated quickly and accurately.
[0085] The following is a comparison of the time complexity of this application with the linear search method and the jump and walk method, to prove that the technical effect achieved by this application is better than that of the prior art, that is, the elevation values in large-scale scenes can be calculated quickly and accurately. For this application, ignoring the time consumption of the Delaunay triangulation process, assuming that the scale of the Delaunay triangulation is m and the scale of the points to be queried is n, in the preprocessing stage, this algorithm takes O(m) to traverse the triangle, and the process of building an undirected graph uses a red-black tree data structure to store edges (that is, to record edge sets). Each triangle has 3 edges stored 2 times, that is, inserted 6 times into the red-black tree, and the overall complexity is O(mlog m). Therefore, the complexity in the KD tree construction stage is O(mlog m); and in the query stage, the KD tree single query complexity is O(log m). Taking the nearest point as the starting point, the path length of the path loop does not exceed 4, so the overall time complexity of querying n points is O(nlog m). Therefore, the overall complexity of this application is: Under the same problem scale, the time complexity of the linear search method is: O(nm); under the same problem scale, the time complexity of the jump and walk method is:
[0086] The embodiment of the present application obtains several elevation point sets to be calculated and triangulates the elevation point sets to obtain several Delaunay triangulations, which can quickly generate a Delaunay triangulation, which is convenient for subsequently determining the triangle where the query point is located; by projecting the Delaunay triangulation onto a two-dimensional plane, a triangle set is obtained, and each triangle in the triangle set is split to obtain a triangle undirected graph, which can project the three-dimensional triangulation onto a two-dimensional plane, and conveniently split each vertex set and edge set in each triangle, which is helpful for subsequently determining the triangle containing the query point. Triangle and elevation calculation; by obtaining the point to be queried, and using the KD tree to search the first point position in the undirected graph that is closest to the point to be queried, the first point position closest to the point to be queried can be quickly obtained, which is convenient for the subsequent determination of the triangle containing the point to be queried and the elevation calculation; based on the point to be queried and the first point position, the triangle containing the point to be queried is determined, and based on the parametric equation of the plane where the triangle is located, the elevation value corresponding to the point to be queried is obtained, and calculation can be performed based on the parametric equation of the triangle containing the point to be queried, and the elevation value in large-scale scenes can be quickly and accurately calculated.
[0087] Example 2
[0088] Please refer to Figure 3 , Figure 3 This is a structural diagram of an embodiment of an elevation calculation system based on a Delaunay triangulation provided by the present application, comprising a first acquisition module 01, a first processing module 02, a second acquisition module 03, and a second processing module 04;
[0089] The first acquisition module 01 is used to acquire a plurality of elevation point sets to be calculated, and triangulate the elevation point sets to obtain a plurality of Delaunay triangulations; the first acquisition module 01 includes: an acquisition unit and an import unit; the acquisition unit is used to acquire a plurality of elevation points using RTK or a total station; the import unit is used to import the plurality of elevation points into CAD to obtain a plurality of elevation point sets to be calculated.
[0090] It can be understood that several elevation points can be collected using RTK or a total station; and by importing several of the elevation points into CAD, several sets of elevation points to be calculated can be obtained, where the elevation points are point information in a three-dimensional space, including x-axis coordinates, y-axis coordinates, and z-axis coordinates.
[0091] In this way, several elevation points are collected in advance through RTK or total station, and imported into CAD to obtain several elevation point sets to be calculated, which facilitates the subsequent construction of Delaunay triangulation.
[0092] After obtaining several elevation point sets, it is necessary to preprocess the elevation point sets and select the Delaunay triangulation algorithm to triangulate the preprocessed elevation point sets to obtain several Delaunay triangulation networks; wherein the triangulation algorithm includes but is not limited to the Bowyer-Watson algorithm, the Lawson local optimization algorithm (LOP), etc., for example and not limitation, taking the Bowyer-Watson algorithm for triangulation as an example, first, it is necessary to construct a sufficiently large super triangle as the initial triangular network, the vertices of which contain all the points in the elevation point set; secondly, the points in the elevation point set are inserted into the triangular mesh one by one, for each newly inserted point, the influencing triangle containing the point is found, and the common edges of the influencing triangles are deleted, and the new point is connected to the vertices of the influencing triangle to form three new triangles. Again, each newly formed triangle is tested for an empty circumcircle to determine whether it meets the Delaunay condition. If the circumcircle of a triangle contains other elevation points, it means that the triangle does not meet the Delaunay condition. At this time, the triangle needs to be optimized by exchanging diagonals until all triangles meet the empty circumcircle condition; then, repeat the above steps of inserting points, constructing triangles, and detecting and optimizing empty circumcircles until all elevation points are inserted into the triangular mesh; finally, delete the triangles related to the super triangle to obtain the final Delaunay triangulation.
[0093] It should be noted that data preprocessing methods include but are not limited to removing duplicate points, processing noise points and filling missing data, and these data preprocessing methods can all be achieved through existing technologies, but data preprocessing methods are not the focus of this application, so they will not be expanded here.
[0094] It should be noted that the Bowyer-Watson algorithm is an incremental algorithm for generating Delaunay triangulations. Delaunay triangulation is a special triangulation method that generates a triangulated mesh that meets the Delaunay conditions by gradually inserting points and rebuilding the triangulated mesh until all points are added. The Lawson Local Optimization Procedure (LOP) is an algorithm used to optimize the triangulated mesh structure during the triangulation process to ensure that it meets the Delaunay triangulation properties. It optimizes triangles that do not meet the Delaunay conditions by exchanging diagonals.
[0095] The first processing module 02 is configured to project the Delaunay triangulation onto a two-dimensional plane to obtain a triangle set, and split each triangle in the triangle set to obtain a triangle undirected graph, wherein the undirected graph includes a vertex set and an edge set;
[0096] It is understandable that when a Delaunay triangulation is obtained, it is necessary to simplify the Delaunay triangulation into a two-dimensional plane for further analysis. Specifically, first, it is necessary to predetermine a suitable projection method to better preserve the structural characteristics of the three-dimensional Delaunay triangulation; second, each point in the Delaunay triangulation is projected into a two-dimensional plane according to the selected projection method. During the projection process, the coordinates of each point on the two-dimensional plane need to be calculated, and the relative position relationship between the points needs to be preserved; finally, the two-dimensional plane is re-triangulated based on the projected point set to obtain a set of triangles in the two-dimensional plane. The triangulation method has been detailed in the first acquisition module, so it will not be repeated here.
[0097] When a triangle set on a two-dimensional plane is obtained, it is necessary to extract all unique vertices from the triangle set and collect all vertices to form a vertex set V, where V = {v1, v2, ..., v n}, where v i Next, we traverse each triangle in the triangle set and add the three edges connected to each triangle to the edge set E. Each edge consists of two vertices, indicating that the two vertices are adjacent in the triangle. For each triangle Δ ijk , whose vertices are v i , v j , v k , its three sides (v i ,v j ),(v i ,v k ),(v j ,v k ) is added to the set of non-phase edges E, that is, (v i ,v j )∈E, if and only if (v j ,v i )∈E. Finally, use the vertex set V and edge set E to construct an undirected graph G=(V,E).
[0098] It should be noted that the projection method includes but is not limited to orthogonal projection, perspective projection and isometric projection, etc. Projecting the Delaunay triangulation into a two-dimensional plane is not the focus of this application and will not be expanded here.
[0099] It's important to note that an undirected graph is a concept in graph theory, where edges have no directionality. That is, an edge connecting two vertices can be considered bidirectional. In an undirected graph, vertices and edges together constitute the basic elements of the graph, and the edge set generally defines the connection between vertices. A triangular undirected graph can be understood as a special type of undirected graph, consisting of multiple triangles formed by three vertices connected end-to-end by three edges.
[0100] The second acquisition module 03 is used to obtain a point to be queried, and use a KD tree to search the undirected graph for a first point that is closest to the point to be queried;
[0101] It can be understood that after obtaining a triangular undirected graph, a KD tree can be established based on the vertex set V to speed up the search efficiency of the nearest point. The process of constructing the KD tree is as follows: (1) arbitrarily select a vertex from the vertex set V as the root node; (2) In each layer of the KD tree, it is necessary to cyclically select different dimensions as the splitting dimensions; (3) For the selected splitting dimension, calculate the median of all vertices on the dimension, and use the median as the root node to divide the remaining points into two subtrees on the left and right; (4) Repeat steps (2) and (3) for the vertices of the left and right subtrees respectively, and recursively construct the left subtree and the right subtree until the subtree is empty.
[0102] Specifically, the second acquisition module 03 includes: a comparison unit and an update unit; the comparison unit is used to compare the point to be queried with the root node. If the value of the point to be queried in the current dimension is less than the value of the root node and the root node has no subtree, the root node is determined to be the second point position, and the first distance between the second point position and the point to be queried is calculated; the update unit is used to search and calculate the second distance between the third point position in the parent node corresponding to the root node and the point to be queried with the smallest distance, and continuously update the second point position based on the second distance until there is no smaller distance, and determine that the second point position is the first point position closest to the point to be queried.
[0103] After the KD tree is constructed, the query point p can be obtained i , and use KD tree to search the undirected graph G with p i Euclidean distance to the nearest point v s The specific query process is as follows: First, select the point p to be queried i The first dimension is the current dimension d, and the query point is compared with the same dimension d t The root node p c Compare, if the query point p i The value at the current dimension is less than the root node p c When the value of , then enter the root node pc Continue searching the left subtree of the root node p; otherwise, go to the root node p c Continue searching the right subtree of the node until the root node p c If there is no subtree, then determine the root node p c The first distance d between the second point and the query point is calculated as c ; Secondly, search and calculate the root node p c The corresponding parent node p t The second distance d with the smallest distance between the third point and the query point a , and based on the second distance d a The second point is continuously updated until there is no smaller distance, and the second point is determined to be the first point v closest to the query point. s For example and not limitation, if the second distance d a <First distance d c , then calculate the parent node p t and the query point p i The distance d t , if d t <d c , then its parent node p t For the node with the closest distance, update the second point, i.e. p c =p t , d c =d t , recursively search for p t The left and right subtrees of are repeatedly recursively searched to determine the point closest to the query point.
[0104] In this way, by using the KD tree to search the first point in the undirected graph that is closest to the point to be queried, the first point that is closest to the point to be queried can be quickly obtained, which facilitates the subsequent determination of the triangle containing the point to be queried and the elevation calculation.
[0105] The second processing module 04 is configured to determine a triangle containing the query point based on the query point and the first point position, and obtain an elevation value corresponding to the query point based on a parametric equation of the plane where the triangle is located.
[0106] The second processing module 04 includes: a third processing unit, a calculation unit, a second comparison unit and a construction unit; the third processing unit is used to make a first ray passing through the point to be queried with the first point as the starting point; the calculation unit is used to extract several candidate points adjacent to the first point from the undirected graph, and make a second ray and a third ray passing through the first point with the first candidate point and the second candidate point respectively, and calculate the angles between the first ray and the second ray and the first ray respectively to obtain the first angle and the second angle; the second comparison unit is used to compare the first angle with the second angle, and if the second angle is smaller than the first angle, determine that the second candidate point corresponding to the second angle is the optimal candidate point; the construction unit is used to continue making rays passing through the point to be queried with the optimal candidate point as the starting point until a repeated point appears, and obtain a triangle containing the point to be queried based on the rays between the repeated points.
[0107] It is understandable that after determining the first point closest to the query point, the next path node v can be searched. p (also referred to as the optimal candidate point below) to determine the triangle containing the query point. Specifically, first, take the first point v s Draw the first ray l as the starting point s Passing the query point p i ; Secondly, arbitrarily extract the first point v from the undirected graph G s The adjacent vertex set V includes several candidate points. The first candidate point v is arbitrarily selected from the vertex set V. p , and take the first candidate point v p Draw the second ray l as the starting point p After the first point v s , where the initial optimal angle d is set p =2π, traverse the vertex set V and select the second candidate point v k , and take the second candidate point v k Draw the third ray l as the starting point k After the first point v s , and calculate the second ray l p and the third ray l k With the first ray l s The angle between the two angles is d. p and the second angle d k ; Again, the first angle d p and the second angle d k For comparison, if the second angle d k smaller than the first angle d p , then determine the second angle d kThe corresponding second candidate point v k is the optimal candidate point, even if v p =v k , and make d p =d k , otherwise keep v p is the optimal candidate point; finally, the optimal candidate point v p Push it into the path queue, and continue to make rays passing through the query point as the starting point, repeat the above steps until a duplicate point appears in the path queue, and take the last 3 points of the path queue, that is, form a triangle loop based on the rays between the duplicate points, and finally get the triangle containing the query point p i The triangle ΔS, where the schematic diagram for determining the optimal candidate point is as follows Figure 2 shown.
[0108] In this way, by determining a triangle containing a point to be queried and performing calculations based on the parametric equation of the determined triangle containing the point to be queried, the elevation value in a large-scale scene can be calculated quickly and accurately.
[0109] When the query point p is determined to be included i The triangle ΔS can get the query point p i Specifically, find the parametric equation f(x, y, z) of the plane where ΔS is located. The process of finding the parametric equation is as follows: (1) Let the three points of triangle ΔS be a(x1, y1, z1), b(x2, y2, z2), and c(x3, y3, z3), and find the two vectors: The cross product of the two vectors can be used to obtain the normal vector of the plane where ΔS is located, and the plane normal vector is brought into the query point p i The x coordinate i ,y i Get the elevation value z i .in, right Expand the determinant and get: according to The plane parametric equation f(x,y,z) can be written as:
[0110] A(x-x1)+B(y-y1)+C(z-z1)=0;
[0111] In the formula, A is ((y2-y1)(z3-z1)-(z2-z1)(y3-y1)), B is ((x2-x1)(z3-z1)-(z2-z1)(x3-x1)), C is ((x2-x1)(y3-y1)-(y2-y1)(x3-x1)), and x1, y1 and z1 are the positions of the first point.
[0112] Therefore, the query point p i (x i ,y i ) into f(x,y,z) to get the corresponding elevation value z i By repeating this method, the elevation values corresponding to several candidate points in the vertex set V can be calculated.
[0113] In this way, by performing calculation based on the determined parametric equation of the triangle containing the point to be queried, the elevation value in a large-scale scene can be calculated quickly and accurately.
[0114] The following is a comparison of the time complexity of this application with the linear search method and the jump and walk method, to prove that the technical effect achieved by this application is better than that of the prior art, that is, the elevation values in large-scale scenes can be calculated quickly and accurately. For this application, ignoring the time consumption of the Delaunay triangulation process, assuming that the scale of the Delaunay triangulation is m and the scale of the points to be queried is n, in the preprocessing stage, this algorithm takes O(m) to traverse the triangle, and the process of building an undirected graph uses a red-black tree data structure to store edges (that is, to record edge sets). Each triangle has 3 edges stored 2 times, that is, inserted 6 times into the red-black tree, and the overall complexity is O(mlog m). Therefore, the complexity in the KD tree construction stage is O(mlog m); and in the query stage, the KD tree single query complexity is O(log m). Taking the nearest point as the starting point, the path length of the path loop does not exceed 4, so the overall time complexity of querying n points is O(nlog m). Therefore, the overall complexity of this application is: Under the same problem scale, the time complexity of the linear search method is: O(nm); under the same problem scale, the time complexity of the jump and walk method is:
[0115] The embodiment of the present application obtains several elevation point sets to be calculated and triangulates the elevation point sets to obtain several Delaunay triangulations, which can quickly generate a Delaunay triangulation, which is convenient for subsequently determining the triangle where the query point is located; by projecting the Delaunay triangulation onto a two-dimensional plane, a triangle set is obtained, and each triangle in the triangle set is split to obtain a triangle undirected graph, which can project the three-dimensional triangulation onto a two-dimensional plane, and conveniently split each vertex set and edge set in each triangle, which is helpful for subsequently determining the triangle containing the query point. Triangle and elevation calculation; by obtaining the point to be queried, and using the KD tree to search the first point position in the undirected graph that is closest to the point to be queried, the first point position closest to the point to be queried can be quickly obtained, which is convenient for the subsequent determination of the triangle containing the point to be queried and the elevation calculation; based on the point to be queried and the first point position, the triangle containing the point to be queried is determined, and based on the parametric equation of the plane where the triangle is located, the elevation value corresponding to the point to be queried is obtained, and calculation can be performed based on the parametric equation of the triangle containing the point to be queried, and the elevation value in large-scale scenes can be quickly and accurately calculated.
[0116] The specific embodiments described above further illustrate the objectives, technical solutions, and beneficial effects of this application. It should be understood that the above descriptions are merely specific embodiments of this application and are not intended to limit the scope of protection of this application. In particular, it should be noted that any modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of this application by those skilled in the art should be included within the scope of protection of this application.
Claims
1. A method for calculating elevation based on Delaunay triangulation, characterized in that: include: Acquire a plurality of elevation point sets to be calculated, and triangulate the elevation point sets to obtain a plurality of Delaunay triangulations, wherein the acquiring the plurality of elevation point sets to be calculated specifically comprises: collecting a plurality of elevation points using an RTK or a total station; importing the plurality of elevation points into a CAD system to obtain a plurality of elevation point sets to be calculated; Projecting the Delaunay triangulation onto a two-dimensional plane to obtain a triangle set, and splitting each triangle in the triangle set to obtain a triangle undirected graph, wherein the undirected graph includes a vertex set and an edge set; Obtain a point to be queried, and use a KD tree to search for the first point in the undirected graph that is closest to the point to be queried. The process of constructing the KD tree is as follows: S1. Select a vertex from the vertex set as the root node; S2. When constructing each layer of the KD tree, cyclically select different dimensions as the splitting dimensions; S3. For the selected splitting dimension, calculate the median of all vertices on the splitting dimension, use the median as the root node, and divide the remaining points into two left and right subtrees; S4. Repeat steps S2 and S3 for the vertices of the left and right subtrees respectively, and recursively construct the left subtree and the right subtree until the subtree is empty; Based on the point to be queried and the first point position, a triangle containing the point to be queried is determined, and based on a parametric equation of a plane where the triangle is located, an elevation value corresponding to the point to be queried is obtained.
2. The elevation calculation method based on Delaunay triangulation according to claim 1, characterized in that: The method of using the KD tree to search the undirected graph for the first point closest to the query point is as follows: Compare the query point with the root node; if the value of the query point at the current dimension is less than the value of the root node and the root node has no subtree, determine that the root node is a second point; and calculate a first distance between the second point and the query point; Search and calculate the second distance between the third point position in the parent node corresponding to the root node and the point to be queried, and continuously update the second point position based on the second distance until there is no smaller distance, and determine that the second point position is the first point position closest to the point to be queried.
3. The elevation calculation method based on Delaunay triangulation according to claim 1 is characterized in that: The determining of a triangle containing the to-be-queried point based on the to-be-queried point and the first point position is specifically as follows: Taking the first point as the starting point, draw a first ray through the point to be queried; Extracting a plurality of candidate points adjacent to the first point from the undirected graph, drawing a second ray and a third ray through the first point respectively using the first candidate point and the second candidate point, and calculating the angles between the second ray and the third ray and the first ray respectively to obtain a first angle and a second angle; Comparing the first angle with the second angle, and if the second angle is smaller than the first angle, determining the second candidate point corresponding to the second angle as the optimal candidate point; The optimal candidate point is used as a starting point to continue drawing rays through the point to be queried until a repeated point appears, and a triangle containing the point to be queried is obtained based on the rays between the repeated points.
4. The elevation calculation method based on Delaunay triangulation according to claim 1 is characterized in that: The parametric equation of the plane where the triangle is located is specifically: ; Where, for , B is , C is , , and is the position of the first point.
5. An elevation calculation system based on Delaunay triangulation, characterized in that: include: a first acquisition module, a first processing module, a second acquisition module, and a second processing module; The first acquisition module is used to acquire a plurality of elevation point sets to be calculated, and triangulate the elevation point sets to obtain a plurality of Delaunay triangulations, wherein the acquisition of the plurality of elevation point sets to be calculated is specifically: using RTK or a total station to collect a plurality of elevation points; importing the plurality of elevation points into CAD to obtain a plurality of elevation point sets to be calculated; The first processing module is configured to project the Delaunay triangulation onto a two-dimensional plane to obtain a triangle set, and split each triangle in the triangle set to obtain a triangle undirected graph, wherein the undirected graph includes a vertex set and an edge set; The second acquisition module is used to obtain a point to be queried, and use a KD tree to search for the first point in the undirected graph that is closest to the point to be queried, wherein the process of constructing the KD tree is as follows: S1, arbitrarily selecting a vertex from a vertex set as a root node; S2, when constructing each layer of the KD tree, cyclically selecting different dimensions as splitting dimensions; S3, for the selected splitting dimension, calculating the median of all vertices on the splitting dimension, and using the median as the root node to divide the remaining points into two left and right subtrees; S4, repeating steps S2 and S3 for the vertices of the left and right subtrees respectively, and recursively constructing the left subtree and the right subtree until the subtree is empty; The second processing module is used to determine a triangle containing the point to be queried based on the point to be queried and the first point position, and obtain an elevation value corresponding to the point to be queried based on a parametric equation of the plane where the triangle is located.
6. The elevation calculation system based on Delaunay triangulation according to claim 5, characterized in that: The second acquisition module includes: a comparison unit and an update unit; The comparison unit is configured to compare the query point with a root node, and if the value of the query point at the current dimension is less than the value of the root node and the root node has no subtree, determine that the root node is a second point position, and calculate a first distance between the second point position and the query point; The updating unit is used to search for and calculate the second distance between the third point position in the parent node corresponding to the root node and the point to be queried, and continuously update the second point position based on the second distance until there is no smaller distance, and determine that the second point position is the first point position closest to the point to be queried.
7. The elevation calculation system based on Delaunay triangulation according to claim 5, characterized in that: The second processing module includes: a third processing unit, a calculation unit, a second comparison unit and a construction unit; The third processing unit is configured to generate a first ray starting from the first point and passing through the point to be queried; The calculation unit is configured to extract a plurality of candidate points adjacent to the first point from the undirected graph, use the first candidate point and the second candidate point as the second ray and the third ray to pass through the first point, and calculate the angles between the second ray and the third ray and the first ray to obtain the first angle and the second angle; The second comparison unit is configured to compare the first angle with the second angle, and if the second angle is smaller than the first angle, determine that a second candidate point corresponding to the second angle is an optimal candidate point; The construction unit is configured to use the optimal candidate point as a starting point to continue drawing rays through the point to be queried until a repeated point appears, and obtain a triangle containing the point to be queried based on the rays between the repeated points.
8. The elevation calculation system based on Delaunay triangulation according to claim 5, characterized in that: The parametric equation of the plane where the triangle is located is specifically: ; Where, for , B is , C is , , and is the position of the first point.
Citation Information
Patent Citations
GPU (Graphics Processing Unit) acceleration-based Voronoi diagram three-dimensional grid reconstruction method
CN108053483A
Method and system for searching positions on topographic map of K-D tree based on curve
CN114547483A