Voronoi diagram and winding number-based isotropic triangular mesh generation method

Abstract

The invention belongs to the technical field of computer graphics geometric processing, and particularly relates to an isotropic triangular mesh generation method based on a Voronoi graph and a winding number, which comprises the following steps: acquiring a triangular mesh model, and generating three sampling source points for each triangular surface through a gravity center and vertex direction interpolation method; obtaining a semi-plane set of each triangular surface; intersecting all half planes of each triangular surface in pairs to obtain a plurality of common edges; projecting each common edge to the surface of the triangular mesh to obtain a Voronoi diagram line segment; obtaining a region point set and a candidate line segment set of each vertex according to line segments of the Voronoi graph; determining a Voronoi unit area of each vertex by using a fast winding number method; performing CVT re-gridding according to the Voronoi unit area of each vertex to obtain an isotropic triangular mesh; the method can guarantee the quality of the generated grid and improve the operation efficiency at the same time.

Description

Technical Field

[0001] The invention belongs to the technical field of computer graphics geometry processing, and in particular relates to an isotropic triangle mesh generation method based on a Voronoi diagram and a winding number. Background Art

[0002] Triangular mesh optimization is an important task in computer graphics. The triangular faces in an isotropic triangular mesh have similar shapes and sizes to equilateral triangles, and are evenly distributed in anisotropy. Isotropic triangular meshes are stable in numerical calculations and are mainly used in finite element analysis, physical simulation, rendering, etc. As the industry's requirements for rendering and simulation accuracy and efficiency continue to increase, it is necessary to face increasingly complex geometric shapes with a surge in the number of triangular faces. Therefore, it is particularly important to improve the efficiency and quality of triangular mesh optimization operations.

[0003] At present, the basic step of triangular mesh optimization is to generate a Voronoi diagram on the triangular mesh surface, and then perform CVT re-meshing based on the Voronoi diagram. There are two ways to achieve this. One is to parameterize the triangular mesh, map the three-dimensional space to the two-dimensional plane, and use the Euclidean distance on the two-dimensional plane to calculate the Voronoi diagram. Calculating the Voronoi diagram on the two-dimensional plane can ensure time efficiency. However, the introduction of the mesh parameterization step will cause the mesh shape to deform, and the computational efficiency of the mesh parameterization cannot be guaranteed. Another idea is to directly use the Euclidean distance or geodesic distance in the three-dimensional space to calculate the Voronoi diagram of the triangular mesh surface. The use of the Euclidean distance method can ensure time efficiency, but the generated Voronoi diagram is not accurate enough. The use of the geodesic distance method can ensure the quality of the Voronoi diagram, but the calculation time is very long.

[0004] Various methods of generating Voronoi diagrams directly on three-dimensional surfaces have brought new ideas to triangular mesh optimization. The Surface Voronoi method can generate half-planes through the Euclidean distance of the source point to the triangular face vertices. The intersection lines formed between the half-planes are projected onto the mesh surface to generate the Voronoi diagram boundary. This can omit the previous Delaunay triangulation process for the triangular mesh and speed up the calculation efficiency. However, when using Surface Voronoi to optimize the triangular mesh into an isotropic mesh, the source point needs to be placed at the vertex position, which will cause the generated Voronoi diagram to have redundant edges passing through the vertex, which will interfere with the accuracy of the CVT process. Summary of the invention

[0005] To solve the above problems, the present invention provides an isotropic triangular mesh generation method based on Voronoi diagram and winding number, which ensures that the triangular mesh is optimized to an isotropic mesh quality while improving the computational efficiency.

[0006] The specific plan includes the following steps:

[0007] S1. Obtain a triangular mesh model and record the topological relationship data of the triangular mesh model;

[0008] S2. For each triangle face in the triangular mesh model, three sampling source points are generated by the centroid and vertex direction interpolation method;

[0009] S3. Obtain a set of contributing source points for each triangular face, and obtain a set of half-planes for each triangular face based on the set of contributing source points;

[0010] S4. intersect all the half planes of each triangular face in pairs to obtain a plurality of common edges; project each common edge onto the triangular mesh surface to obtain a Voronoi diagram line segment;

[0011] S5. Obtain the region point set and candidate segment set of each vertex according to the Voronoi diagram segment;

[0012] S6. Determine the Voronoi cell area of ​​each vertex using a fast winding number method based on the regional point set, the candidate line segment set, and the 1-neighborhood triangle area of ​​the vertex;

[0013] S7. Perform CVT remeshing according to the Voronoi cell area of ​​each vertex to obtain an isotropic triangular mesh.

[0014] Beneficial effects of the present invention:

[0015] The invention provides an isotropic triangular mesh generation method based on a Voronoi diagram and a winding number, which can ensure the quality of the generated mesh and improve the operation efficiency at the same time. BRIEF DESCRIPTION OF THE DRAWINGS

[0016] Figure 1 It is a flow chart of the isotropic triangular mesh generation method based on the Voronoi diagram and the winding number according to the present invention;

[0017] Figure 2 A schematic diagram of a source point close to a vertex for generating a triangular face according to the present invention;

[0018] Figure 3 A schematic diagram of the Voronoi region determined using the fast winding number according to the present invention;

[0019] Figure 4Schematic diagram of a triangular mesh formed after re-meshing using CVT according to the present invention. DETAILED DESCRIPTION

[0020] The following will be combined with the drawings in the embodiments of the present invention to clearly and completely describe the technical solutions in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0021] As the application of winding number in the field of geometric processing gradually matures, the winding number from a point to a closed curve can be used to determine whether a point is inside or outside the closed range. For this reason, the winding number can be used to determine whether a vertex is inside or outside the region in surface simplification and mesh reconstruction, thereby removing noise and improving mesh quality. The fast winding number extends the winding number to sparse data such as point clouds, and ensures the efficiency of computing time.

[0022] Based on this, the present invention provides an isotropic triangular mesh generation method based on Voronoi diagram and winding number, such as Figure 1 As shown, the following steps are included:

[0023] S1. Obtain a triangular mesh model and record the topological relationship data of the triangular mesh model.

[0024] S2. For each triangle face in the triangular mesh model, three sampling source points are generated by the centroid and vertex direction interpolation method.

[0025] Specifically, step S2 uses the centroid and vertex direction interpolation method to generate sampling source points for any triangular surface, including:

[0026] S21. Obtain the coordinates of the vertices A, B, and C of the triangle surface, and calculate the coordinates of the centroid D of the triangle surface;

[0027] S22. Perform interpolation operations on the direction vectors from the centroid D to vertex A, vertex B, and vertex C respectively to obtain three sampling source points; the sampling source points generated on the direction vector from the centroid D to vertex A are expressed as

[0028]

[0029] Among them, t represents the interpolation parameter, (x A ,y A ,z A ) represents the coordinates of vertex A, (x D ,y D ,zD ) represents the coordinates of the center of gravity D, and p represents the coordinates of the sampling source point obtained by interpolating the direction vector of the center of gravity D and the vertex A.

[0030] In the embodiment of the present invention, t=0.01 is set so that the generated sampling source point is as close to the vertex as possible, such as Figure 2 If the sampling source point is generated at the vertex position, the sampling source point will affect other 1-neighboring triangles of the vertex, thereby introducing redundant noise points; therefore, the sampling source point needs to be limited to the current triangle range, and the sampling source point should be generated as close to the vertex position as possible, so that the vertex can affect the sampling source point in the smallest range, and the Voronoi boundary is confirmed based on the affected sampling source point.

[0031] S3. Obtain a set of contributing source points of each triangular face, and obtain a set of half-planes of each triangular face based on the set of contributing source points.

[0032] Specifically, the present invention uses the Surface Voronoi method based on Euclidean distance to traverse each triangular face, and uses the mark scanning method to progressively select the contributing source point set of each triangular face. Specifically, in step S3, taking any triangular face as the target triangular face, the process of obtaining the contributing source point set of the target triangular face includes:

[0033] S31. For each vertex of the target triangle, calculate the Euclidean distance between it and the three sampling source points generated by the target triangle, record the shortest Euclidean distance as the shortest distance, and take the sampling source point corresponding to the shortest Euclidean distance as the contributing source point; obtain all the sampling source points in the 1-neighborhood triangle of the target triangle, and execute step S32;

[0034] S32. Calculate the Euclidean distance between each sampling source point and each vertex of the target triangle. If the current Euclidean distance is less than the shortest distance recorded for the corresponding vertex, update the shortest distance of the corresponding vertex and use the current Euclidean distance corresponding to the sampling source point as the contributing source point.

[0035] S33. Determine whether the sampling source points of all triangles in the triangular mesh model have been traversed. If so, obtain all contributing source points. If not, obtain all sampling source points within the 1-neighborhood triangle of the current 1-neighborhood triangle, and return to execute step S32.

[0036] Specifically, obtaining a half-plane set according to all contributing source points of the target triangle includes:

[0037] Define the coordinates of the three vertices of the target triangle as (x1, y1, z1), (x2, y2, z2), and (x3, y3, z3);

[0038] For each contributing source point, calculate the square distance d1, d2, d3 between it and the three vertices of the target triangle, expressed as

[0039]

[0040] (x p ,y p ,z p ) is the contributing source point coordinate; thereby, three new coordinates (x1, y1, d1), (x2, y2, d2), (x3, y3, d3) are obtained; according to the three new coordinates, a plane is determined as the half-plane.

[0041] S4. Intersect all the half-planes of each triangular face with each other to obtain multiple common edges; project each common edge onto the triangular mesh surface to obtain a Voronoi diagram segment.

[0042] In the embodiment of the present invention, Figure 3 As shown in the result figure on the left, according to the interpolated sampling source points, the Surface Voronoi method based on Euclidean distance is used to traverse each triangular face, and the contributing source points of the triangular face are progressively selected using the marked scanning method. The half-plane corresponding to each contributing source point is obtained according to the distance from the contributing source point to the vertex. The half-planes of multiple contributing source points intersect to form multiple common edges, which are projected onto the triangular face to obtain the Voronoi diagram line segments of the triangular face, and then a corresponding Voronoi unit area will be formed around each vertex. At this time, there are Voronoi diagram line segments passing through the vertex. This is because each triangular face has three sampling source points. In addition to the sampling source points corresponding to the vertex itself, there are two other sampling source points corresponding to the vertices that will affect the generation of the Voronoi unit area.

[0043] S5. Obtain a region point set and a candidate segment set for each vertex according to the Voronoi diagram segment.

[0044] Specifically, the specific process of obtaining a regional point set of a vertex in step S5 includes:

[0045] S51. Get all Voronoi diagram segments in the 1-neighborhood triangle of the vertex;

[0046] S52. Given a distance threshold threshold, for each Voronoi diagram segment, calculate the Euclidean distances from its two end points to the vertex respectively. As long as there is a Euclidean distance less than the distance threshold threshold, the endpoint corresponding to the maximum value of the two Euclidean distances is selected as the region point; if both Euclidean distances are not less than the distance threshold threshold, the Voronoi diagram segment is recorded as a candidate segment.

[0047] S6. Determine the Voronoi cell area of ​​each vertex using the fast winding number method based on the regional point set, the candidate line segment set and the 1-neighborhood triangle area.

[0048] Specifically, step S6 uses the fast winding number method to determine the Voronoi unit area of ​​a vertex, including:

[0049] S61. Taking the two endpoints of each candidate line segment as two query points;

[0050] S62. For each candidate line segment, the winding numbers of the two query points are calculated respectively. If there is a query point with a winding number greater than 0, it is considered that all query points of the candidate line segment are within the Voronoi unit area of ​​the vertex, and the two query points of the candidate line segment are used as the boundary points of the Voronoi unit area of ​​the vertex; the calculation formula is:

[0051]

[0052] Among them, w(q) represents the winding number of the query point q, m represents the number of regional points corresponding to the vertex, and p i represents the i-th region point, a i Represents the area of ​​the triangle where the i-th region point is located, n i Represents the normal vector of the i-th region point.

[0053] In the embodiment of the present invention, the result after denoising by the fast winding number method is as follows: Figure 3 Shown on the right.

[0054] S7. Perform CVT remeshing according to the Voronoi cell area of ​​each vertex to obtain an isotropic triangular mesh.

[0055] Specifically, Figure 4 As shown in the figure, the original mesh on the left is transformed into the optimized triangular mesh on the right after CVT. CVT can optimize the triangular mesh by Lloyd iteration or minimizing the energy function. However, in order to optimize the triangular mesh in three-dimensional space, CCVT (Constrained Centroidal Voronoi Tessellations) is used to minimize the energy function to optimize the triangular mesh, find the target centroid of the Voronoi unit area corresponding to the vertex, and transform the vertex coordinates to the centroid coordinates to optimize the triangular mesh.

[0056] Step S7 includes:

[0057] Calculate the energy function E(X), expressed as

[0058]

[0059] Among them, Ω i represents the ith Voronoi unit area, ρ(·) represents the density function, x i represents the sampling point in the i-th Voronoi unit area (the sampling point can be obtained in a variety of ways, as long as the sampling point is guaranteed to be in the Voronoi unit area, and the sampling source point close to the vertex is selected to speed up the fitting); k represents the number of Voronoi unit areas, q i Represents the boundary point of the i-th Voronoi unit area; the energy function here uses an integral form to calculate the sum of distances, is to calculate the sampling point x i The sum of the distances to all boundary points of the i-th Voronoi cell region.

[0060] The density function is based on the curvature:

[0061] ρ(x)=1+β·|σ(x)|

[0062] Where β is the tuning parameter and σ is the local curvature of the point.

[0063] Minimize the energy function E(X) to obtain the gradient Expressed as

[0064]

[0065] Among them, c i represents the target centroid of the ith Voronoi cell area, q i,j represents the jth boundary point of the i-th Voronoi cell region, and h represents the number of boundary points in the Voronoi cell region.

[0066] Update the sampling point position according to the gradient to make it close to the target centroid position, which is expressed as

[0067]

[0068] Here, α is the gradient descent rate.

[0069] The sampling point coordinates are continuously updated iteratively, and the vertex coordinates are moved to the sampling point coordinates, and finally an isotropic triangular mesh is obtained.

[0070] In the present invention, unless otherwise clearly stipulated and limited, the terms such as "installation", "setting", "connection", "fixation", "rotation" and the like should be understood in a broad sense. For example, it can be a fixed connection, a detachable connection, or an integral one; it can be a mechanical connection or an electrical connection; it can be a direct connection or an indirect connection through an intermediate medium; it can be the internal connection of two elements or the interaction relationship between two elements. Unless otherwise clearly defined, for ordinary technicians in this field, the specific meanings of the above terms in the present invention can be understood according to the specific circumstances.

[0071] Although embodiments of the present invention have been shown and described, it will be appreciated by those skilled in the art that various changes, modifications, substitutions and variations may be made to the embodiments without departing from the principles and spirit of the present invention, and that the scope of the present invention is defined by the appended claims and their equivalents.

Claims

1. A method for generating isotropic triangular meshes based on Voronoi diagram and winding number, characterized in that: The following steps are involved: S1. Obtain a triangular mesh model and record the topological relationship data of the triangular mesh model; S2. For each triangle face in the triangular mesh model, three sampling source points are generated by the centroid and vertex direction interpolation method; S3. Obtain a set of contributing source points for each triangular face, and obtain a set of half-planes for each triangular face based on the set of contributing source points; S4. intersect all the half planes of each triangular face in pairs to obtain a plurality of common edges; project each common edge onto the triangular mesh surface to obtain a Voronoi diagram line segment; S5. Obtain the region point set and candidate segment set of each vertex according to the Voronoi diagram segment; S6. Determine the Voronoi cell area of ​​each vertex using a fast winding number method based on the regional point set, the candidate line segment set, and the 1-neighborhood triangle area of ​​the vertex; S7. Perform CVT remeshing according to the Voronoi cell area of ​​each vertex to obtain an isotropic triangular mesh.

2. The isotropic triangular mesh generation method based on Voronoi diagram and winding number according to claim 1, characterized in that: Step S2 uses the centroid and vertex direction interpolation method to generate sampling source points for any triangular face, including: S21. Obtain the coordinates of the vertices A, B, and C of the triangle surface, and calculate the coordinates of the centroid D of the triangle surface; S22. Perform interpolation operations on the direction vectors from the centroid D to vertex A, vertex B, and vertex C respectively to obtain three sampling source points; the sampling source points generated on the direction vector from the centroid D to vertex A are expressed as Among them, t represents the interpolation parameter, (x A ,y A ,z A ) represents the coordinates of vertex A, (x D ,y D ,z D ) represents the coordinates of the center of gravity D, and p represents the coordinates of the sampling source point obtained by interpolating the direction vector of the center of gravity D and the vertex A.

3. The isotropic triangular mesh generation method based on Voronoi diagram and winding number according to claim 1, characterized in that: In step S3, the process of taking any triangular face as the target triangular face and obtaining the contributing source point set of the target triangular face includes: S31. For each vertex of the target triangle, calculate the Euclidean distance between it and the three sampling source points generated by the target triangle, record the shortest Euclidean distance as the shortest distance, and take the sampling source point corresponding to the shortest Euclidean distance as the contributing source point; obtain all the sampling source points in the 1-neighborhood triangle of the target triangle, and execute step S32; S32. Calculate the Euclidean distance between each sampling source point and each vertex of the target triangle. If the current Euclidean distance is less than the shortest distance recorded for the corresponding vertex, update the shortest distance of the corresponding vertex and use the current Euclidean distance corresponding to the sampling source point as the contributing source point. S33. Determine whether the sampling source points of all triangles in the triangular mesh model have been traversed. If so, obtain all contributing source points. If not, obtain all sampling source points within the 1-neighborhood triangle of the current 1-neighborhood triangle, and return to execute step S32.

4. The isotropic triangular mesh generation method based on Voronoi diagram and winding number according to claim 3, characterized in that: Obtaining a half-plane set based on all contributing source points of the target triangle includes: Define the coordinates of the three vertices of the target triangle as (x1, y1, z1), (x2, y2, z2), and (x3, y3, z3); For each contributing source point, calculate the square distances d1, d2, and d3 between it and the three vertices of the target triangle to obtain three new coordinates (x1, y1, d1), (x2, y2, d2), and (x3, y3, d3); determine a plane as the half-plane based on the three new coordinates.

5. The isotropic triangular mesh generation method based on Voronoi diagram and winding number according to claim 1, characterized in that: The specific process of obtaining a regional point set of a vertex in step S5 includes: S51. Get all Voronoi diagram segments in the 1-neighborhood triangle of the vertex; S52. For each Voronoi diagram segment, calculate the Euclidean distances from its two end points to the vertex respectively. As long as there is a Euclidean distance less than the distance threshold, select the endpoint corresponding to the maximum value of the two Euclidean distances as the regional point; if both Euclidean distances are not less than the distance threshold, record the Voronoi diagram segment as a candidate segment.

6. The isotropic triangular mesh generation method based on Voronoi diagram and winding number according to claim 1, characterized in that: Step S6 uses a fast winding number method to determine the Voronoi unit area of ​​a vertex, including: S61. Taking the two endpoints of each candidate line segment as two query points; S62. For each candidate line segment, the winding numbers of the two query points are calculated respectively. If the winding number of a query point is greater than 0, the two query points of the candidate line segment are used as the boundary points of the Voronoi unit area of ​​the vertex; the calculation formula is: Among them, w(q) represents the winding number of the query point, m represents the number of regional points corresponding to the vertex, and p i represents the i-th region point, a i Represents the area of ​​the triangle where the i-th region point is located, n i Represents the normal vector of the i-th region point.

7. The isotropic triangular mesh generation method based on Voronoi diagram and winding number according to claim 1, characterized in that: Step S7 includes: S71. Calculate the energy function E(X), expressed as Among them, Ω i represents the ith Voronoi unit area, ρ(·) represents the density function, x i represents the sampling point in the i-th Voronoi unit area; k represents the number of Voronoi unit areas, q i Represents the boundary point of the i-th Voronoi cell area; S72. Minimize the energy function E(X) to obtain the gradient Expressed as Among them, c i represents the target centroid of the ith Voronoi cell area, q i,j represents the jth boundary point of the i-th Voronoi cell area, and h represents the number of boundary points in the Voronoi cell area; S73. Update the sampling point position according to the gradient to make it close to the target centroid position, expressed as Where α is the gradient descent rate; S74. Continuously iterate and update the sampling point coordinates to eventually obtain an isotropic triangular mesh.

Images