An isotropic triangular mesh generation method based on voronoi diagram and winding number

Abstract

The application belongs to the technical field of computer graphics geometry processing, and particularly relates to an isotropic triangle mesh generation method based on Voronoi diagram and winding number, which comprises the following steps: obtaining a triangle mesh model; for each triangular face, generating three sampling source points through a barycenter and vertex direction interpolation method; obtaining a half-plane set of each triangular face; intersecting all half-planes of each triangular face to obtain multiple common edges; projecting each common edge to the surface of the triangular mesh to obtain a Voronoi diagram line segment; obtaining a regional point set and a candidate line segment set of each vertex according to the Voronoi diagram line segment; determining the Voronoi unit region of each vertex by using a fast winding number method; and performing CVT re-meshing according to the Voronoi unit region of each vertex to obtain an isotropic triangle mesh.

Description

TECHNICAL FIELD

[0001] The application belongs to the technical field of computer graphics geometry processing, and particularly relates to an isotropic triangular mesh generation method based on a Voronoi diagram and a winding number. BACKGROUND

[0002] Triangular mesh optimization is an important task in computer graphics, in which the triangular faces in the isotropic triangular mesh have similar shapes and sizes to equilateral triangles, and are uniformly distributed in all directions. Isotropic triangular mesh is stable in numerical calculation and is mainly applied to finite element analysis, physical simulation, rendering and the like. With the increasing requirements of the industry on the accuracy and efficiency of rendering and simulation, more and more complex geometries with a sharp increase in the number of triangular faces need to be faced, and therefore it is particularly important to improve the operation efficiency and quality of triangular mesh optimization.

[0003] At present, the basic steps of triangular mesh optimization are to generate a Voronoi diagram on the surface of the triangular mesh, and then to perform CVT re-meshing based on the Voronoi diagram. There are two ways to achieve this. One is to perform mesh parameterization on the triangular mesh, map the three-dimensional space to a two-dimensional plane, and use Euclidean distance to calculate the Voronoi diagram on the two-dimensional plane. 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 operation efficiency of the mesh parameterization cannot be guaranteed. The other way is to directly calculate the Voronoi diagram on the surface of the triangular mesh using Euclidean distance or geodesic distance in three-dimensional space. Using Euclidean distance can ensure time efficiency, but the generated Voronoi diagram is not accurate enough. Using geodesic distance can ensure the quality of the Voronoi diagram, but the operation time is very long.

[0004] Various methods of directly generating a Voronoi diagram on a three-dimensional surface provide new ideas for triangular mesh optimization. The Surface Voronoi method can generate half-planes by the Euclidean distance of the source points to the vertices of the triangular faces, and the intersection lines formed between the half-planes are projected onto the surface of the mesh to generate the boundary of the Voronoi diagram, which can omit the previous Delaunay triangulation process for the triangular mesh, thereby accelerating the operation efficiency. However, using Surface Voronoi to optimize the triangular mesh into an isotropic mesh requires the source points to be placed at the vertex positions, which will cause the generated Voronoi diagram to have extra edges passing through the vertices, which will interfere with the accuracy of the CVT process. SUMMARY

[0005] To solve the above problems, the application provides an isotropic triangular mesh generation method based on a Voronoi diagram and winding number, which ensures that the triangular mesh is optimized to be an isotropic mesh quality and improves the operation efficiency.

[0006] The specific scheme comprises 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 triangular face in the triangular mesh model, generate three sampling source points through a barycenter and vertex direction interpolation method;

[0009] S3. Obtain a contribution source point set of each triangular face, and obtain a half-plane set of each triangular face according to the contribution source point set;

[0010] S4. Intersect all half-planes of each triangular face two by two to obtain a plurality of common edges; and project each common edge to the surface of the triangular mesh to obtain a Voronoi diagram line segment;

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

[0012] S6. According to the region point set, the candidate line segment set and the area of the 1-neighborhood triangular face of the vertex, determine the Voronoi unit region of each vertex by using a fast winding number method;

[0013] S7. Perform CVT re-meshing according to the Voronoi unit region of each vertex to obtain an isotropic triangular mesh.

[0014] The application has the following beneficial effects:

[0015] The isotropic triangular mesh generation method based on the Voronoi diagram and the winding number can ensure the generated mesh quality and improve the operation efficiency. BRIEF DESCRIPTION OF DRAWINGS

[0016] Figure 1 The application provides a flowchart of the isotropic triangular mesh generation method based on the Voronoi diagram and the winding number;

[0017] Figure 2 The application provides a schematic diagram of the source point close to the vertex of the triangular face;

[0018] Figure 3 The application provides a schematic diagram of the Voronoi region determined by using the fast winding number;

[0019] Figure 4A schematic view of a triangular mesh formed after re-meshing using the CVT described in the present application. DETAILED DESCRIPTION

[0020] The technical solutions in the embodiments of the present application will be apparently and completely described below with reference to the drawings in the embodiments of the present application. Obviously, the described embodiments are only part of the embodiments of the present application, rather than all the embodiments. Based on the embodiments in the present application, all the other embodiments obtained by a person of ordinary skill in the art without creative work fall within the protection scope of the present application.

[0021] With the gradual maturity of the application of winding number in the field of geometric processing, the winding number of a point to a closed curve can be used to determine whether the point is inside or outside the closed range. Therefore, the winding number can be applied to determine whether a vertex is inside or outside a region in the surface simplification and mesh reconstruction work, and then remove the noise points and improve the mesh quality. The fast winding number (Fast Winding Number) extends the winding number to sparse data such as point cloud, and ensures the operation time efficiency.

[0022] Based on this, the present application provides an isotropic triangular mesh generation method based on Voronoi diagram and winding number, as shown in Figure 1 , comprising the following steps:

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

[0024] S2. For each triangular face in the triangular mesh model, generate 3 sampling source points by the barycenter and vertex direction interpolation method.

[0025] Specifically, the step S2 generates the sampling source points for any triangular face by the barycenter and vertex direction interpolation method, comprising:

[0026] S21. Obtain the coordinates of the vertices A, B and C of the triangular face, and calculate the coordinates of the barycenter D of the triangular face;

[0027] S22. Perform interpolation operation on the direction vectors from the barycenter D to the vertices A, B and C respectively, to obtain 3 sampling source points; the sampling source point generated on the direction vector from the barycenter D to the vertex A is represented as

[0028]

[0029] Wherein, t represents an interpolation parameter, (x A ,y A ,z A ) represents the coordinates of the 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 on the direction vector from the center of gravity D to the vertex A.

[0030] In the embodiment of the present application, t=0.01 is set, so that the generated sampling source point is as close to the vertex as possible, as shown in the figure. Figure 2 If the sampling source point is generated at the vertex position, the sampling source point will affect other 1-neighbor triangle faces of the vertex, thereby introducing redundant noise points; therefore, the sampling source point needs to be limited within the current triangle face range, and the sampling source point is generated as close to the vertex as possible, so that the vertex can minimally affect the sampling source point, and the Voronoi boundary is confirmed according to the affected sampling source point.

[0031] S3. Obtain a set of contribution source points of each triangle face, and obtain a set of half-planes of each triangle face according to the set of contribution source points.

[0032] Specifically, the Surface Voronoi method based on Euclidean distance is used to traverse each triangle face, and the marked scanning method is used to progressively select the set of contribution source points of each triangle face. Specifically, in step S3, taking any one triangle face as a target triangle face, the process of obtaining the set of contribution source points of the target triangle face includes:

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

[0034] S32. Calculate the Euclidean distance between each sampling source point and each vertex of the target triangle face, if the current Euclidean distance is smaller than the shortest distance recorded for the corresponding vertex, update the shortest distance of the corresponding vertex, and take the sampling source point corresponding to the current Euclidean distance as the contribution source point;

[0035] S33. Determine whether all sampling source points of the triangle mesh model have been traversed, if yes, obtain all contribution source points, if not, obtain all sampling source points in the 1-neighbor triangle face of the current 1-neighbor triangle face, and return to execute step S32.

[0036] Specifically, obtaining the set of half-planes according to all the contribution source points of the target triangle face includes:

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

[0038] For each contribution source point, calculate its squared distance d1, d2, d3 with the 3 vertices of the target triangle, denoted as

[0039]

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

[0041] S4. Intersect every two half planes of each triangle to obtain multiple common edges; project each common edge to the triangle mesh surface to obtain a Voronoi diagram line segment.

[0042] In the embodiment of the present application, as Figure 3 As shown in the left result diagram, according to the interpolated sampling source points, use the Surface Voronoi method based on Euclidean distance to traverse each triangle, and use the marking scanning method to progressively select the contribution source points of the triangle. According to the distance of the contribution source points to the vertices, obtain the half plane corresponding to each contribution source point, intersect the half planes of multiple contribution source points to form multiple common edges, project the common edges to the triangle to obtain the Voronoi diagram line segment of the triangle, and then a corresponding Voronoi unit region will be formed around each vertex. At this time, there is a Voronoi diagram line segment passing through the vertex, because each triangle has three sampling source points, in addition to the sampling source point corresponding to the vertex itself, there are two sampling source points corresponding to the other two vertices, which will affect the generation of the Voronoi unit region.

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

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

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

[0046] S52. Given a distance threshold threshold, for each Voronoi diagram line segment, calculate the Euclidean distance from the two endpoints to the vertex, as long as there is a Euclidean distance less than the distance threshold threshold, select the endpoint corresponding to the maximum of the two Euclidean distances as the region point; if both Euclidean distances are not less than the distance threshold threshold, record the Voronoi diagram line segment as a candidate line segment.

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

[0048] Specifically, step S6 determines the Voronoi cell region of a vertex using the fast winding number method, including:

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

[0050] S62. For each candidate line segment, respectively calculate the winding number of the two query points, if the winding number of one query point is greater than 0, it is considered that all query points of the candidate line segment are in the Voronoi cell region of the vertex, and the two query points of the candidate line segment are taken as the boundary points of the Voronoi cell region of the vertex; the calculation formula is

[0051]

[0052] where w(q) represents the winding number of the query point q, m represents the number of region points corresponding to the vertex, p i represents the ith region point, a i represents the area of the triangle face where the ith region point is located, n i represents the normal vector of the ith region point.

[0053] In the embodiment of the present application, the result after denoising by the fast winding number method is as shown on the right side. Figure 3

[0054] S7. Perform CVT re-meshing according to the Voronoi cell region of each vertex to obtain an isotropic triangular mesh.

[0055] Specifically, as shown in the left side of FIG. 7, the original mesh is changed into the optimized triangular mesh on the right side after CVT. The CVT can be optimized by the Lloyd iteration method or the energy function minimization method, but in order to optimize the triangular mesh in the three-dimensional space, the CCVT (Constrained Centroidal Voronoi Tessellations) is used to minimize the energy function to optimize the triangular mesh, and the vertex coordinates are transformed to the centroid coordinates to optimize the triangular mesh. Figure 4 Step S7 includes:

[0056] Calculate the energy function E(X), which is expressed as

[0057]

[0058]

[0059] ​​where Ω i represents the i-th Voronoi cell region, p(·) represents the density function, x i represents the sampling point in the i-th Voronoi cell region (the sampling point can be obtained in various ways, as long as the sampling point is in the Voronoi cell region, the sampling source point close to the vertex is taken, and the fitting is accelerated) ; k represents the number of Voronoi cell regions, q i represents the boundary point of the i-th Voronoi cell region; the energy function here uses the integral form to calculate the distance sum, is the distance sum of the sampling point x i to all boundary points of the i-th Voronoi cell region.

[0060] The density function is a curvature-based density function:

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

[0062] where β is a tuning parameter, and σ is the local curvature of the point.

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

[0064]

[0065] where c i represents the target centroid of the i-th Voronoi cell region, q i,j represents the j-th boundary point of the i-th Voronoi cell region, and h represents the number of boundary points of the Voronoi cell region.

[0066] Update the sampling point position according to the gradient to approximate the target centroid position, which is represented as

[0067]

[0068] where α is the gradient descent rate.

[0069] Continuously update the sampling point coordinates, and move the vertex coordinates to the sampling point coordinates, to finally obtain the isotropic triangular mesh.

[0070] In this application, unless otherwise clearly specified and limited, the terms "mounting", "setting", "connecting", "fixing", "rotating" and other terms should be understood in a broad sense, for example, can be fixedly connected, can also be detachably connected, or integrated; can be mechanically connected, can also be electrically connected; can be directly connected, can also be indirectly connected through an intermediate medium, can be the internal communication of two elements or the interaction relationship of two elements, unless otherwise clearly limited, the above terms in this application can be understood according to the specific meaning of the above terms by the person skilled in the art.

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

Claims

1. A method for generating an isotropic triangular mesh based on Voronoi diagram and winding number, characterized in that, The method comprises the following steps: S1. Obtain a triangular mesh model and record the topological relationship data of the triangular mesh model; S2. For each triangular face in the triangular mesh model, generate three sampling source points by using the barycenter and vertex direction interpolation method; S3. Obtain a set of contribution source points of each triangular face, and obtain a set of half-planes of each triangular face according to the set of contribution source points; In step S3, the process of obtaining the set of contribution source points of the target triangular face includes: S31. For each vertex of the target triangular face, calculate the Euclidean distance between the vertex and the three sampling source points generated by the target triangular face respectively, record the shortest Euclidean distance as the shortest distance, and take the sampling source point corresponding to the shortest Euclidean distance as the contribution source point; obtain all sampling source points in the 1-neighborhood triangular face of the target triangular face, and execute step S32; S32. Calculate the Euclidean distance between each sampling source point and each vertex of the target triangular face, if the current Euclidean distance is smaller than the shortest distance recorded by the corresponding vertex, update the shortest distance of the corresponding vertex, and take the sampling source point corresponding to the current Euclidean distance as the contribution source point; S33. Determine whether all sampling source points of the triangular mesh model have been traversed, if yes, obtain all contribution source points, if not, obtain all sampling source points in the 1-neighborhood triangular face of the current 1-neighborhood triangular face, and return to execute step S32; According to all the contribution source points of the target triangular face, the set of half-planes includes: Define the coordinates of the three vertices of the target triangular face as (x1, y1, z1), (x2, y2, z2), and (x3, y3, z3); For each contribution source point, calculate the squared distances d1, d2, and d3 between the contribution source point and the three vertices of the target triangular face, and obtain three new coordinates (x1, y1, d1), (x2, y2, d2), and (x3, y3, d3); determine a plane as a half-plane according to the three new coordinates; S4. Intersect all half-planes of each triangular face pairwise to obtain multiple 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 the candidate line segment set of each vertex according to the Voronoi diagram line segment; S6. Determine the Voronoi cell region of each vertex by using the fast winding number method according to the region point set, the candidate line segment set, and the area of the 1-neighborhood triangular face of the vertex; S7. Perform CVT re-meshing according to the Voronoi cell region of each vertex to obtain an isotropic triangular mesh.

2. The method according to claim 1, wherein, In step S2, the barycenter and vertex direction interpolation method is used to generate sampling source points for any triangular face, which includes: S21. Obtain the coordinates of the vertices A, B, and C of the triangular face, and calculate the coordinates of the barycenter D of the triangular face; S22. Perform interpolation operation on the direction vectors from the barycenter D to the vertices A, B, and C respectively to obtain three sampling source points; the sampling source point generated on the direction vector from the barycenter D to the vertex A is represented as where t represents an interpolation parameter, (x A ,y A ,z A ) represents the coordinates of the vertex A, (x D ,y D ,z D ) represents the coordinates of the barycenter D, and p represents the coordinates of the sampling source point obtained by performing an interpolation operation on the direction vectors of the barycenter D and the vertex A.

3. The method of claim 1, wherein, The specific process of obtaining the region point set of a vertex in step S5 includes: S51. Obtain all Voronoi diagram line segments in the 1-neighborhood triangle of the vertex; S52. For each Voronoi diagram line segment, calculate the Euclidean distance from each end point to the vertex, and select the end point corresponding to the maximum Euclidean distance as the region point if there is one Euclidean distance less than the distance threshold; if both Euclidean distances are not less than the distance threshold, record the Voronoi diagram line segment as a candidate line segment.

4. The method of claim 1, wherein, Step S6 uses the fast winding number method to determine the Voronoi cell region of a vertex, including: S61. Take the two end points of each candidate line segment as two query points; S62. For each candidate line segment, calculate the winding number of the two query points, respectively, and if there is one query point with a winding number greater than 0, take the two query points of the candidate line segment as the boundary points of the Voronoi cell region of the vertex; the calculation formula is where w(q) represents the winding number of the query point, m represents the number of region points corresponding to the vertex, p i represents the i-th region point, a i represents the area of the triangular face where the i-th region point is located, n i represents the normal vector of the i-th region point.

5. The method of claim 1, wherein, Step S7 includes: S71. Calculate the energy function E(X), represented as where Ω i represents the i-th Voronoi cell region, p(·) represents the density function, x i represents the sampling point within the i-th Voronoi cell region; k represents the number of Voronoi cell regions, q i represents the boundary point of the i-th Voronoi cell region; S72. Minimize the energy function E(X) to obtain the gradient is represented as where c i represents the target centroid of the i-th Voronoi cell region, q i,j represents the j-th boundary point of the i-th Voronoi cell region, and h represents the number of boundary points within the Voronoi cell region. S73. Update the sampling point position according to the gradient to make it approximate the target centroid position, represented as Wherein, α is the gradient descent rate; S74. Continuously update the sampling point coordinates, and finally obtain the isotropic triangular mesh.

Images