High-order continuity maintenance and error control method in STL curved surface reconstruction process

Abstract

The invention provides a high-order continuity maintenance and error control method in an STL curved surface reconstruction process, and relates to the technical field of three-dimensional modeling, and the method comprises the steps: obtaining the connection information of discrete data, constructing a continuous normal vector field based on spherical harmonic decomposition, carrying out the frequency domain processing of a vertex coordinate and a normal vector, and recognizing high-frequency details, and calculating a main curvature direction to determine an insertion control point, and iteratively optimizing the normal deviation until convergence. According to the method, the high-order continuity of the STL curved surface can be effectively kept, the reconstruction error is accurately controlled, and the geometric accuracy and quality of the reconstructed curved surface are improved.

Description

Technical Field

[0001] This invention relates to the field of 3D modeling technology, and in particular to a method for maintaining high-order continuity and controlling errors in the STL surface reconstruction process. Background Technology

[0002] STL, a widely used 3D model representation format in additive manufacturing, numerical simulation, and computer-aided design, is typically composed of discretized triangular facets and can approximate complex surfaces. In high-precision manufacturing and analysis applications, it is necessary to reconstruct discrete triangular facets into parametric surfaces with high-order continuity to meet the requirements of product design and precision machining.

[0003] Traditional surface reconstruction methods mainly include point cloud-based surface fitting, mesh-based parametric methods, and subdivision surface techniques. These methods convert discrete data into continuous representations, such as NURBS surfaces and B-spline surfaces, using different mathematical models. The key challenges in surface reconstruction are maintaining geometric continuity, controlling approximation errors, and handling complex topological structures. However, existing STL surface reconstruction techniques still struggle to simultaneously guarantee high-order geometric continuity and accurate local feature representation, resulting in continuity breaks or over-smoothing. They typically employ global optimization strategies, making it difficult to effectively control errors for local high-frequency features, leading to insufficient reconstruction accuracy. Furthermore, they cannot automatically adjust the control point density based on the local geometric characteristics of the surface, resulting in wasted computational resources and low representation efficiency. Summary of the Invention

[0004] This invention provides a method for maintaining high-order continuity and controlling errors in the STL surface reconstruction process, which can at least solve some of the problems existing in the prior art.

[0005] A first aspect of this invention provides a method for preserving high-order continuity and controlling errors during STL surface reconstruction, comprising: Discretized geometric data is obtained and the shared edges and vertices of triangular facets are identified to obtain connection information. Based on the connection information, the vertex coordinates of the common boundary on adjacent triangular facets are determined and projected onto the tangent plane of the corresponding triangular facets to obtain parameter coordinates. Extract the vertex normal vectors on both sides of the common boundary and perform spherical harmonic decomposition to obtain spherical harmonic coefficients. Based on the spherical harmonic coefficients and the parameter coordinates, construct a continuous normal vector field and perform path integration to obtain the surface patch control points and curvature tensor components. The vertex coordinates and the vertex normal vector are transformed in the frequency domain to obtain frequency domain coefficients. High-frequency coefficients are determined by combining the frequency domain coefficients with a preset coefficient threshold. The high-frequency coefficients are then transformed inversely to obtain high-frequency vertex coordinates and determine high-frequency triangular facets. Calculate the principal curvature direction of the parameter coordinate region corresponding to the high-frequency triangular patch based on the curvature tensor components, determine the insertion position and additional control points based on the principal curvature direction, and merge them with the control points of the surface patch to obtain the merged control points; Based on the merged control points and the parameter coordinates, a parameterized surface patch is constructed and the surface sampling points are determined. The normal projection distance from the surface sampling points to the surface of the triangular patch is calculated to obtain the normal deviation and the deviation spectrum is determined. The merged control points are iteratively updated according to the amplitude distribution of the deviation spectrum until the normal deviation converges, and the converged control points are obtained. The reconstructed surface is generated by combining the parameter coordinates.

[0006] In one alternative implementation, Discretized geometric data is acquired, and shared edges and vertices of triangular facets are identified to obtain connection information. Based on this connection information, the vertex coordinates of the common boundary on adjacent triangular facets are determined and projected onto the tangent plane of the corresponding triangular facets to obtain parameter coordinates, including: Obtain discrete geometric data containing multiple triangular facets and traverse the vertex coordinates of each triangular facet in the discrete geometric data. Calculate the vertex overlap between any two triangular facets. Identify triangular facet pairs with vertex overlap greater than zero as adjacent triangular facets and extract the overlapping vertices of the adjacent triangular facets as shared vertices. Based on the shared vertices, the edges of the adjacent triangular faces are traversed, and edges containing the same pair of shared vertices are identified as shared edges. An index relationship between the shared edges and the adjacent triangular faces is established to obtain the connection information. Based on the connection information, the shared edge of each pair of adjacent triangular facets is extracted, the shared vertex on the shared edge is determined as the common boundary vertex, and the vertex coordinates of the common boundary vertex are extracted as the common boundary vertex coordinates. Calculate the centroid coordinates of the three vertices of the triangular facet, construct a tangent plane with the centroid coordinates as the origin based on the vertex normal vectors of the triangular facet, project the coordinates of the common boundary vertices onto the tangent plane to obtain the projected coordinates, construct a two-dimensional coordinate system on the tangent plane, and represent the projected coordinates in the two-dimensional coordinate system to obtain the parameter coordinates.

[0007] In one alternative implementation, Extracting the vertex normal vectors on both sides of the common boundary and performing spherical harmonic decomposition yields spherical harmonic coefficients. Based on the spherical harmonic coefficients and the parameter coordinates, a continuous normal vector field is constructed and path integration is performed to obtain the surface patch control points and curvature tensor components, including: Based on the connection information, the adjacent triangular facets on both sides of the common boundary are identified, the vertex normal vectors of the adjacent triangular facets are extracted as the vertex normal vectors of the common boundary, and the vertex normal vectors of the common boundary are mapped to a unit sphere to obtain spherical coordinates. The spherical coordinates are expanded using a spherical harmonic function, and the function value of the spherical harmonic function at the spherical coordinates is calculated. Based on the function value and the normal vector of the common boundary vertex, the least squares fitting is performed to obtain the spherical harmonic coefficients of each order. A set of grid points is established based on the parameter coordinates. For each grid point, the reconstructed value of the spherical harmonic function is calculated based on the spherical harmonic coefficient. The reconstructed value of the spherical harmonic function is mapped from the unit sphere to the three-dimensional space to obtain the grid point normal vector. The continuous normal vector field is constructed based on the grid point normal vector. Extract the normal vector components of the continuous normal vector field along the parameter coordinate direction, perform line integration on the normal vector components along the corresponding parameter coordinate direction to obtain the integral curve, and determine the control points of the surface patch based on the intersection coordinates of the integral curves in different parameter coordinate directions. Calculate the partial derivatives of the continuous normal vector field with respect to each parameter coordinate direction, and construct the curvature tensor components based on the partial derivatives.

[0008] In one alternative implementation, Performing a frequency domain transformation on the vertex coordinates and the vertex normal vector to obtain frequency domain coefficients, and determining high-frequency coefficients by combining them with a preset coefficient threshold, and then performing an inverse transformation on the high-frequency coefficients to obtain high-frequency vertex coordinates and determine high-frequency triangular patches includes: The vertex coordinates of all triangular facets in the discretized geometric data are extracted and combined with the connection information to construct a vertex adjacency matrix. Based on the vertex adjacency matrix, a normalized graph Laplacian matrix is ​​constructed and spectral decomposition is performed to obtain a set of eigenvectors. The vertex coordinates and the vertex normal vectors are jointly embedded in the spectral space formed by the set of eigenvectors. The coordinate frequency domain coefficients and normal frequency domain coefficients are calculated by minimizing the reconstruction error. The gradient correlation between the coordinate frequency domain coefficients and the normal frequency domain coefficients is calculated, and the magnitudes of the coordinate frequency domain coefficients and the normal frequency domain coefficients are adaptively weighted based on the gradient correlation to obtain the spectral energy distribution. The spectral energy distribution is then subjected to multi-scale wavelet decomposition to extract energy abrupt change points, and the spectral energy value corresponding to the energy abrupt change point is used as the coefficient threshold. The coordinate frequency domain coefficients whose spectral energy exceeds the coefficient threshold are identified as the high-frequency coefficients, and the corresponding feature vectors are extracted to construct a high-frequency basis. The high-frequency coefficients are then back-projected onto the high-frequency basis and optimized and reconstructed in combination with preset local geometric constraints to obtain the high-frequency vertex coordinates. Traverse all triangular facets, identify triangular facets containing the high-frequency vertex coordinates and calculate the corresponding discrete Gaussian curvature and mean curvature, construct a curvature tensor based on the discrete Gaussian curvature and mean curvature and calculate the principal curvature, and mark triangular facets with principal curvature exceeding a preset curvature threshold as the high-frequency triangular facets.

[0009] In one alternative implementation, The principal curvature direction of the parametric coordinate region corresponding to the high-frequency triangular patch is calculated based on the curvature tensor components. The insertion position and additional control points are determined based on the principal curvature direction and merged with the control points of the surface patch to obtain merged control points, including: Extract the parameter coordinate region corresponding to the high-frequency triangular patch, obtain the curvature tensor components of each parameter grid point in the parameter coordinate region and perform feature decomposition to obtain eigenvalues ​​and eigenvectors, take the eigenvector corresponding to the maximum eigenvalue as the maximum principal curvature direction, and take the eigenvector corresponding to the minimum eigenvalue as the minimum principal curvature direction. Curvature change rate fields are constructed along the direction of maximum principal curvature and the direction of minimum principal curvature, respectively. Morphological gradient operations are performed on the curvature change rate fields to identify extreme points of curvature change rate. The parameter coordinates of the extreme points of curvature change rate in the parameter coordinate region are used as the insertion positions. At the insertion position, the covariant derivatives of the curvature tensor are calculated along the direction of maximum principal curvature and the direction of minimum principal curvature, respectively. A local curvature compensation vector is constructed based on the covariant derivatives. The parameter coordinates of the insertion position are offset along the local curvature compensation vector to obtain the parameter coordinates of the additional control point. The spatial coordinates of the additional control points are obtained by performing surface mapping calculation on the parameter coordinates of the additional control points, and the surface patch control points are obtained. The spatial coordinates of the additional control points and the surface patch control points are merged in the order of parameter coordinates, and the merged control points are topologically reorganized to obtain the merged control points.

[0010] In one alternative implementation, Based on the merged control points and the parameter coordinates, a parameterized surface patch is constructed and the surface sampling points are determined. The normal projection distance from the surface sampling points to the surface of the triangular patch is calculated to obtain the normal deviation and the deviation spectrum is determined, including: Obtain the spatial coordinates of the merged control points and the parameter boundary of the parameter coordinate region. Construct a T-spline basis function based on the spatial coordinates of the merged control points. Insert node vectors into the parameter coordinate region and generate a non-uniform parameter grid using a local refinement algorithm. Calculate the spatial coordinates of the surface sampling points based on the nodes in the non-uniform parameter grid and the T-spline basis function. An octree spatial index is constructed based on the spatial coordinates of the surface sampling point, and the nearest triangular facet is queried. The normal vector of the tangent plane of the parameterized surface facet at the surface sampling point is calculated. The surface sampling point is projected onto the surface of the nearest triangular facet along the direction of the tangent plane normal vector to obtain the projection point. The directed distance between the surface sampling point and the projection point is calculated as the normal projection distance. The parameter coordinates corresponding to the nodes in the non-uniform parameter grid are used as the sampling positions, and the normal projection distance is used as the deviation amplitude at the sampling positions. A non-uniform sampling deviation field is constructed based on the sampling positions and the deviation amplitude. The non-uniform sampling deviation field is subjected to a non-uniform fast Fourier transform to obtain a frequency domain representation, and the power spectral density of each frequency component in the frequency domain representation is extracted to obtain the deviation spectrum.

[0011] In one alternative implementation, Iteratively updating the merged control points based on the amplitude distribution of the deviation spectrum until the normal deviation converges, obtaining converged control points, and generating a reconstructed surface by combining the parameter coordinates includes: The amplitude of each frequency component in the deviation spectrum is extracted and frequency-weighted to obtain the frequency domain deviation energy. A control point sensitivity matrix is ​​constructed based on the frequency domain deviation energy. The dominant sensitivity direction is obtained by performing singular value decomposition on the control point sensitivity matrix. The position adjustment amount of the merged control point is calculated along the dominant sensitivity direction. The spatial coordinates of the merged control point are updated along the position adjustment amount to obtain the optimized merged control point. Based on the optimized merged control points and the parameter coordinates, the parameterized surface patch is reconstructed and the updated surface sampling points are determined. The normal projection distance from the updated surface sampling points to the surface of the triangular patch is calculated to obtain the optimized normal deviation. Calculate the rate of change of the deviation between the optimized normal deviation and the normal deviation before the update, and determine whether the rate of change of the deviation is less than a preset convergence threshold. If it is less, the optimized merged control point is taken as the convergence control point. If it is not less, the optimized merged control point is taken as the new merged control point, and the process of constructing the control point sensitivity matrix and calculating the normal deviation is repeated until the rate of change of the deviation is less than the convergence threshold. Based on the convergence control point and the parameter coordinates, the parameterized surface patch is reconstructed to obtain the reconstructed surface.

[0012] A second aspect of the present invention provides a high-order continuity preservation and error control system for STL surface reconstruction, comprising: The information extraction unit is used to acquire discretized geometric data and identify the shared edges and vertices of triangular facets to obtain connection information. Based on the connection information, the vertex coordinates of the common boundary on adjacent triangular facets are determined and projected onto the tangent plane of the corresponding triangular facets to obtain parameter coordinates. The normal field construction unit is used to extract the vertex normal vectors on both sides of the common boundary, perform spherical harmonic decomposition to obtain spherical harmonic coefficients, construct a continuous normal vector field based on the spherical harmonic coefficients and the parameter coordinates, and perform path integration to obtain the surface patch control points and curvature tensor components. The high-frequency identification unit is used to perform frequency domain transformation on the vertex coordinates and the vertex normal vector to obtain frequency domain coefficients and determine high-frequency coefficients by combining them with a preset coefficient threshold. The high-frequency coefficients are then inversely transformed to obtain high-frequency vertex coordinates and determine high-frequency triangular facets. The control point solving unit is used to calculate the principal curvature direction of the parameter coordinate region corresponding to the high-frequency triangular patch based on the curvature tensor components, determine the insertion position and additional control points based on the principal curvature direction, and merge them with the control points of the surface patch to obtain the merged control points. The surface reconstruction unit is used to construct a parameterized surface patch based on the merged control points and the parameter coordinates, determine the surface sampling points, calculate the normal projection distance from the surface sampling points to the surface of the triangular patch to obtain the normal deviation and determine the deviation spectrum, iteratively update the merged control points according to the amplitude distribution of the deviation spectrum until the normal deviation converges, obtain the converged control points, and generate the reconstructed surface by combining the parameter coordinates.

[0013] A third aspect of the present invention provides an electronic device, comprising: A processor and a memory for storing processor-executable instructions, wherein the processor is configured to invoke instructions stored in the memory to perform the aforementioned method.

[0014] A fourth aspect of the present invention provides a computer-readable storage medium having stored thereon computer program instructions that, when executed by a processor, implement the aforementioned method.

[0015] In this invention, by performing spherical harmonic decomposition on vertex normal vectors and constructing a continuous normal vector field, high-order continuity between adjacent surface patches is preserved. By performing frequency domain transformation on vertex coordinates and vertex normal vectors and extracting high-frequency information, complex detail regions of the surface can be accurately identified, improving the ability of the reconstructed surface to express local geometric features. The method of calculating the principal curvature direction and determining the insertion position and additional control points based on curvature tensor components makes the distribution of control points more reasonable and better adaptable to the geometric change characteristics of the surface. By adopting a normal deviation convergence mechanism and deviation spectrum analysis method, the error in the surface reconstruction process is precisely controlled, ensuring the geometric accuracy of the reconstructed surface and the original STL model. Attached Figure Description

[0016] Figure 1 This is a flowchart illustrating the method for preserving high-order continuity and controlling errors during STL surface reconstruction according to an embodiment of the present invention. Figure 2 This is a flowchart of the frequency domain deviation optimization for surface reconstruction in the STL surface reconstruction process according to an embodiment of the present invention. Detailed Implementation

[0017] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0018] The technical solution of the present invention will be described in detail below with reference to specific embodiments. These specific embodiments can be combined with each other, and the same or similar concepts or processes may not be described again in some embodiments.

[0019] Figure 1 This is a flowchart illustrating the method for preserving high-order continuity and controlling errors during STL surface reconstruction according to an embodiment of the present invention. Figure 1 As shown, the method includes: Discretized geometric data is obtained and the shared edges and vertices of triangular facets are identified to obtain connection information. Based on the connection information, the vertex coordinates of the common boundary on adjacent triangular facets are determined and projected onto the tangent plane of the corresponding triangular facets to obtain parameter coordinates. Extract the vertex normal vectors on both sides of the common boundary and perform spherical harmonic decomposition to obtain spherical harmonic coefficients. Based on the spherical harmonic coefficients and the parameter coordinates, construct a continuous normal vector field and perform path integration to obtain the surface patch control points and curvature tensor components. The vertex coordinates and the vertex normal vector are transformed in the frequency domain to obtain frequency domain coefficients. High-frequency coefficients are determined by combining the frequency domain coefficients with a preset coefficient threshold. The high-frequency coefficients are then transformed inversely to obtain high-frequency vertex coordinates and determine high-frequency triangular facets. Calculate the principal curvature direction of the parameter coordinate region corresponding to the high-frequency triangular patch based on the curvature tensor components, determine the insertion position and additional control points based on the principal curvature direction, and merge them with the control points of the surface patch to obtain the merged control points; Based on the merged control points and the parameter coordinates, a parameterized surface patch is constructed and the surface sampling points are determined. The normal projection distance from the surface sampling points to the surface of the triangular patch is calculated to obtain the normal deviation and the deviation spectrum is determined. The merged control points are iteratively updated according to the amplitude distribution of the deviation spectrum until the normal deviation converges, and the converged control points are obtained. The reconstructed surface is generated by combining the parameter coordinates.

[0020] In one alternative implementation, Discretized geometric data is acquired, and shared edges and vertices of triangular facets are identified to obtain connection information. Based on this connection information, the vertex coordinates of the common boundary on adjacent triangular facets are determined and projected onto the tangent plane of the corresponding triangular facets to obtain parameter coordinates, including: Obtain discrete geometric data containing multiple triangular facets and traverse the vertex coordinates of each triangular facet in the discrete geometric data. Calculate the vertex overlap between any two triangular facets. Identify triangular facet pairs with vertex overlap greater than zero as adjacent triangular facets and extract the overlapping vertices of the adjacent triangular facets as shared vertices. Based on the shared vertices, the edges of the adjacent triangular faces are traversed, and edges containing the same pair of shared vertices are identified as shared edges. An index relationship between the shared edges and the adjacent triangular faces is established to obtain the connection information. Based on the connection information, the shared edge of each pair of adjacent triangular facets is extracted, the shared vertex on the shared edge is determined as the common boundary vertex, and the vertex coordinates of the common boundary vertex are extracted as the common boundary vertex coordinates. Calculate the centroid coordinates of the three vertices of the triangular facet, construct a tangent plane with the centroid coordinates as the origin based on the vertex normal vectors of the triangular facet, project the coordinates of the common boundary vertices onto the tangent plane to obtain the projected coordinates, construct a two-dimensional coordinate system on the tangent plane, and represent the projected coordinates in the two-dimensional coordinate system to obtain the parameter coordinates.

[0021] Obtain discretized geometric data containing multiple triangular facets. The discretized geometric data can be 3D model data obtained through 3D scanning, computer-aided design systems, or other means, containing multiple triangular facets, each defined by three vertex coordinates. For example, the acquired data might contain a spherical model composed of thousands of triangular facets, each with its vertex coordinates represented in 3D space. Specifically, for a given triangular facet in the model, its three vertex coordinates might be: the first vertex coordinates are (10.5, 15.2, 8.7), the second vertex coordinates are (11.2, 15.8, 8.9), and the third vertex coordinates are (10.8, 15.5, 9.4).

[0022] The algorithm iterates through the vertex coordinates of each triangular facet in the discretized geometric data, calculating the vertex overlap between any two triangular facets. Vertex overlap refers to the number of vertices shared by two triangular facets. The overlap is determined by comparing the vertex coordinates of two triangular facets. When the overlap is greater than zero, the two triangular facets are identified as adjacent facets, and their overlapping vertices are extracted as shared vertices. For example, the coordinates of the three vertices of triangle facet A are (10.5, 15.2, 8.7), (11.2, 15.8, 8.9), and (10.8, 15.5, 9.4), while the coordinates of the three vertices of triangle facet B are (10.5, 15.2, 8.7), (10.8, 15.5, 9.4), and (10.2, 15.0, 9.1). By comparing the coordinates, it was found that the two triangular faces share two vertices, namely (10.5, 15.2, 8.7) and (10.8, 15.5, 9.4). Therefore, the degree of overlap between the vertices of the two triangular faces is 2, they are adjacent triangular faces, and the two overlapping vertices are shared vertices.

[0023] By traversing the edges of adjacent triangular faces based on shared vertices, edges containing the same pair of shared vertices are identified as shared edges. An edge of a triangular face is defined by two vertices; when two triangular faces have two identical shared vertices, the edge connecting these two vertices is a shared edge. An index relationship is established between shared edges and adjacent triangular faces to obtain connection information. For example, triangular face A and triangular face B share vertices (10.5, 15.2, 8.7) and (10.8, 15.5, 9.4), therefore, the edge defined by these two vertices is a shared edge. The index relationship records that this shared edge connects triangular face A and triangular face B; this index relationship constitutes the connection information.

[0024] Based on the connection information, the shared edges of each pair of adjacent triangular faces are extracted, and the shared vertices on the shared edges are determined as common boundary vertices. The coordinates of the common boundary vertices are then extracted as the common boundary vertex coordinates. For example, for triangular face A and triangular face B, the two shared vertices on the shared edge (10.5, 15.2, 8.7) and (10.8, 15.5, 9.4) are the common boundary vertices, with coordinates of (10.5, 15.2, 8.7) and (10.8, 15.5, 9.4) respectively.

[0025] Calculate the centroid coordinates of the three vertices of a triangular facet. The centroid coordinates are the arithmetic mean of the three vertex coordinates of the facet. For example, if the three vertices of facet A are (10.5, 15.2, 8.7), (11.2, 15.8, 8.9), and (10.8, 15.5, 9.4), its centroid coordinates are (10.83, 15.50, 9.00). Construct a tangent plane with the centroid coordinates as the origin, based on the vertex normals of the facet. The vertex normals can be calculated by the cross product of two edges of the facet. For example, for facet A, its vertex normals might be (0.18, 0.32, 0.93). Construct a tangent plane with the centroid coordinates (10.83, 15.50, 9.00) as the origin and the vertex normals (0.18, 0.32, 0.93) as the normals.

[0026] Projecting the coordinates of the common boundary vertices onto the tangent plane yields the projected coordinates. The projection process involves projecting the common boundary vertices onto the tangent plane along the normal vector direction of the tangent plane. For example, projecting the common boundary vertex (10.5, 15.2, 8.7) onto the tangent plane of triangle A yields the projected coordinates (10.52, 15.24, 8.74). A two-dimensional coordinate system is then constructed on the tangent plane, and the projected coordinates are expressed in this system to obtain the parametric coordinates. When constructing the two-dimensional coordinate system, two orthogonal directions on the tangent plane can be selected as coordinate axes. For example, (0.98, -0.20, 0.00) and (0.07, 0.34, -0.94) can be selected as the basis vectors of the two-dimensional coordinate system. Transforming the projected coordinates (10.52, 15.24, 8.74) into this two-dimensional coordinate system yields the parametric coordinates (-0.32, 0.27).

[0027] In this embodiment, by introducing a vertex overlap determination mechanism, the true adjacency relationship between any two triangular faces is automatically identified, avoiding the incompleteness and error accumulation problems caused by relying solely on preset topological indexes or manual annotations. This improves the adaptability and robustness to complex, irregular, or reconstruction error-prone discrete geometric models. By accurately extracting shared vertices and shared edges and establishing index relationships between shared edges and adjacent triangular faces, the unique determination of common boundary vertices is achieved, eliminating the problems of boundary repetition, inconsistent boundary directions, or boundary ambiguity between adjacent faces. This improves the consistency and accuracy of the face connection relationship expression, which is beneficial to the stable execution of subsequent geometric analysis and numerical calculations. By constructing a local tangent plane with the centroid of the triangular face as a reference and projecting the coordinates of the common boundary vertices onto this tangent plane, a unified representation of the common boundary in the local two-dimensional space is achieved while maintaining the original three-dimensional geometric relative relationships. This significantly reduces computational complexity and reduces numerical instability problems caused by spatial rotation and pose differences.

[0028] In one alternative implementation, Extracting the vertex normal vectors on both sides of the common boundary and performing spherical harmonic decomposition yields spherical harmonic coefficients. Based on the spherical harmonic coefficients and the parameter coordinates, a continuous normal vector field is constructed and path integration is performed to obtain the surface patch control points and curvature tensor components, including: Based on the connection information, the adjacent triangular facets on both sides of the common boundary are identified, the vertex normal vectors of the adjacent triangular facets are extracted as the vertex normal vectors of the common boundary, and the vertex normal vectors of the common boundary are mapped to a unit sphere to obtain spherical coordinates. The spherical coordinates are expanded using a spherical harmonic function, and the function value of the spherical harmonic function at the spherical coordinates is calculated. Based on the function value and the normal vector of the common boundary vertex, the least squares fitting is performed to obtain the spherical harmonic coefficients of each order. A set of grid points is established based on the parameter coordinates. For each grid point, the reconstructed value of the spherical harmonic function is calculated based on the spherical harmonic coefficient. The reconstructed value of the spherical harmonic function is mapped from the unit sphere to the three-dimensional space to obtain the grid point normal vector. The continuous normal vector field is constructed based on the grid point normal vector. Extract the normal vector components of the continuous normal vector field along the parameter coordinate direction, perform line integration on the normal vector components along the corresponding parameter coordinate direction to obtain the integral curve, and determine the control points of the surface patch based on the intersection coordinates of the integral curves in different parameter coordinate directions. Calculate the partial derivatives of the continuous normal vector field with respect to each parameter coordinate direction, and construct the curvature tensor components based on the partial derivatives.

[0029] Adjacent triangular faces on both sides of a common boundary are identified based on connectivity information. In a 3D geometric model, each common boundary is typically shared by two triangular faces. By querying the index relationships in the connectivity information, adjacent triangular faces on both sides of the common boundary can be identified. For example, for a vertex on the common boundary, the connectivity information can identify two triangular faces associated with that vertex, located on either side of the common boundary. The vertex normal vectors of the adjacent triangular faces are extracted as the vertex normal vectors of the common boundary. For example, if the normal vectors of the two triangular faces associated with a vertex of a common boundary are (0.18, 0.32, 0.93) and (0.15, 0.36, 0.92), these two vectors are the normal vectors of the vertex of the common boundary.

[0030] Mapping the normal vectors of the common boundary vertices to a unit sphere yields spherical coordinates. The normal vectors in 3D space are then standardized to unit vectors and converted to their representation in spherical coordinates. For example, standardizing the normal vector (0.18, 0.32, 0.93) results in the unit vector (0.1824, 0.3243, 0.9419), whose corresponding spherical coordinates can be represented as angle values, such as (1.0585, 0.3399). The first value represents the angle with the positive axis, and the second value represents the angle in the vertical plane; all angle values ​​are in radians.

[0031] Spherical coordinates are expanded using spherical harmonic functions. Spherical harmonic functions are special functions defined on a sphere and can be used to represent data distributed on that sphere. In practice, a series of predefined spherical harmonic function basis functions are typically used to calculate the function values ​​of the basis functions at a given spherical coordinate. For example, for the spherical coordinates (1.0585, 0.3399), different orders of spherical harmonic function values ​​can be calculated, such as the zeroth-order function value of 0.2821, first-order function values ​​of 0.2466, 0.1682, and 0.4897, and second-order function values ​​of 0.1534, 0.2638, 0.0952, 0.3721, and 0.1856, etc.

[0032] The spherical harmonic coefficients of various orders are obtained by least-squares fitting of the function values ​​and the normal vectors of the common boundary vertices. Using the least-squares method, a set of spherical harmonic coefficients can be found such that a linear combination of the spherical harmonic functions best fits the normal vectors of the common boundary vertices. For example, through fitting calculations, the zeroth-order coefficient might be 0.9832, the first-order coefficients might be 0.0324, 0.1278, and 0.0591, and the second-order coefficients might be 0.0097, 0.0156, 0.0082, 0.0129, and 0.0073, etc.

[0033] A set of grid points is established based on parametric coordinates. Grid points can be set uniformly or non-uniformly in the parameter space to form a grid point set. For example, within the parametric coordinate range, a 10×10 uniform grid of points with a grid spacing of 0.1 can be set. For each grid point, a reconstructed spherical harmonic function value is calculated based on the spherical harmonic coefficient. The reconstructed spherical harmonic function value at each grid point is obtained by multiplying the spherical harmonic coefficient by the corresponding spherical harmonic function value and summing the results. For example, for the grid point at parametric coordinates (0.3, 0.4), the calculated reconstructed spherical harmonic function value might be (0.9526, 0.2345).

[0034] The reconstructed values ​​of the spherical harmonic function are inversely mapped from the unit sphere to three-dimensional space to obtain the grid point normal vectors. This step transforms the spherical coordinates back into a vector representation in three-dimensional space. For example, the spherical coordinates (0.9526, 0.2345) are inversely mapped to the three-dimensional normal vectors (0.1953, 0.3452, 0.9178). Based on the grid point normal vectors, a continuous normal vector field is constructed, forming a continuous normal vector distribution across the entire parameter space.

[0035] Extract the normal vector components of the continuous normal vector field along the parametric coordinate directions. Extract the components of the normal vector field in the two-dimensional parameter space along the two parametric directions respectively. For example, for the normal vector (0.1953, 0.3452, 0.9178), the components in the first parametric direction (0.1953, 0.3452, 0) and the components in the second parametric direction (0, 0, 0.9178) can be extracted. Perform line integration on the normal vector components along the corresponding parametric coordinate directions to obtain the integral curves. The line integral calculates the cumulative change of the normal vector components along a path in the parameter space. For example, performing line integration along the first parametric direction from coordinates (0, 0) to (1, 0) yields a series of coordinates, forming the integral curve; similarly, another series of integral curves can be obtained along the second parametric direction.

[0036] Control points for a surface patch are determined by the intersection coordinates of integral curves in different parameter coordinate directions. The intersection points of integral curves are the points where integral curves in different parameter directions intersect, and these intersection points can serve as control points for the surface patch. For example, the intersection coordinates (11.25, 16.34, 9.78) of the integral curves in the first parameter direction and the integral curves in the second parameter direction in three-dimensional space can be used as a control point. By calculating the intersection points of multiple integral curves, a set of control points can be obtained, which define the shape and curvature characteristics of the surface patch.

[0037] The partial derivatives of the continuous normal vector field with respect to each parameter coordinate direction are calculated, and curvature tensor components are constructed based on these partial derivatives. The partial derivatives represent the rate of change of the normal vector field in the parameter direction. For example, at the parameter point (0.3, 0.4), the partial derivative of the normal vector field with respect to the first parameter direction might be (0.0253, 0.0421, 0.0127), and the partial derivative with respect to the second parameter direction might be (0.0185, 0.0376, 0.0098). These partial derivatives can be used to construct curvature tensor components, such as the principal curvature values ​​of 0.0462 and 0.0287, and the principal curvature directions (0.7321, 0.6812) and (-0.6812, 0.7321). The curvature tensor components describe the degree of curvature of the surface in each direction and are important indicators for evaluating the continuity and smoothness of the surface.

[0038] In this embodiment, by mapping the normal vectors of the common boundary vertices to a unit sphere and introducing spherical harmonic function expansion, a smooth fitting of the discrete, unevenly distributed, and noise-sensitive normal vector information under global constraints is performed. This effectively suppresses abrupt changes in normal vectors and local oscillations, significantly improving the continuity and stability of the normal vector field within the common boundary and its neighborhood. By reconstructing the normal vector field using spherical harmonic coefficients obtained from least-squares fitting, a continuous and differentiable representation of normal information in the parametric coordinate domain is achieved. This ensures that the normal vectors are consistent not only at discrete vertices but also maintain a smooth transition of direction changes within and across surface regions. This overcomes the problem that the normal vectors can only be effectively expressed at a limited number of sampling points and cannot support high-order geometric analysis. By using the intersection points of line integrals of the continuous normal vector field in different parametric coordinate directions to determine the control points of the surface patch, the method that originally relied on empirical rules or local geometric heuristics is transformed into a deterministic calculation process based on the overall normal consistency, improving the rationality and geometric consistency of the control point distribution.

[0039] In one alternative implementation, Performing a frequency domain transformation on the vertex coordinates and the vertex normal vector to obtain frequency domain coefficients, and determining high-frequency coefficients by combining them with a preset coefficient threshold, and then performing an inverse transformation on the high-frequency coefficients to obtain high-frequency vertex coordinates and determine high-frequency triangular patches includes: The vertex coordinates of all triangular facets in the discretized geometric data are extracted and combined with the connection information to construct a vertex adjacency matrix. Based on the vertex adjacency matrix, a normalized graph Laplacian matrix is ​​constructed and spectral decomposition is performed to obtain a set of eigenvectors. The vertex coordinates and the vertex normal vectors are jointly embedded in the spectral space formed by the set of eigenvectors. The coordinate frequency domain coefficients and normal frequency domain coefficients are calculated by minimizing the reconstruction error. The gradient correlation between the coordinate frequency domain coefficients and the normal frequency domain coefficients is calculated, and the magnitudes of the coordinate frequency domain coefficients and the normal frequency domain coefficients are adaptively weighted based on the gradient correlation to obtain the spectral energy distribution. The spectral energy distribution is then subjected to multi-scale wavelet decomposition to extract energy abrupt change points, and the spectral energy value corresponding to the energy abrupt change point is used as the coefficient threshold. The coordinate frequency domain coefficients whose spectral energy exceeds the coefficient threshold are identified as the high-frequency coefficients, and the corresponding feature vectors are extracted to construct a high-frequency basis. The high-frequency coefficients are then back-projected onto the high-frequency basis and optimized and reconstructed in combination with preset local geometric constraints to obtain the high-frequency vertex coordinates. Traverse all triangular facets, identify triangular facets containing the high-frequency vertex coordinates and calculate the corresponding discrete Gaussian curvature and mean curvature, construct a curvature tensor based on the discrete Gaussian curvature and mean curvature and calculate the principal curvature, and mark triangular facets with principal curvature exceeding a preset curvature threshold as the high-frequency triangular facets.

[0040] The vertex coordinates of all triangular faces in the discretized geometric data are extracted and combined with connectivity information to construct a vertex adjacency matrix. The vertex adjacency matrix describes the topology of a 3D mesh, where matrix elements indicate whether vertices are connected. For example, for a 3D model with 1000 vertices, a 1000×1000 adjacency matrix can be constructed. If vertex 5 and vertex 8 are connected, the element in the 5th row, 8th column and the 8th row, 5th column of the matrix has a value of 1; otherwise, it has a value of 0. A normalized graph Laplacian matrix is ​​then constructed based on the vertex adjacency matrix. The normalized graph Laplacian matrix is ​​an important tool in graph theory used to analyze the topological properties of mesh structures. For a model with 1000 vertices, the corresponding normalized graph Laplacian matrix is ​​also a 1000×1000 matrix, where the diagonal elements are the degree of the vertex (the number of edges connected to that vertex) divided by the total number of edges, and the off-diagonal elements are the negative values ​​of the corresponding elements in the adjacency matrix divided by the total number of edges.

[0041] Spectral decomposition of the normalized graph Laplacian matrix yields a set of eigenvectors. Spectral decomposition is the process of decomposing a matrix into eigenvalues ​​and eigenvectors. For example, for a 1000×1000 normalized graph Laplacian matrix, spectral decomposition yields 1000 eigenvalues ​​and their corresponding eigenvectors. These eigenvectors form the basis of the spectral space and can be sorted according to the magnitude of the eigenvalues, such as from smallest to largest, resulting in eigenvectors vector1, vector2, ..., vector1000. Vertex coordinates and vertex normal vectors are then jointly embedded in the spectral space formed by the eigenvector set. Joint embedding is the process of simultaneously projecting spatial coordinates and normal vector information onto the spectral space. For vertex coordinates, the three-dimensional coordinates (11.25, 16.34, 9.78) can be multiplied by each eigenvector to obtain their representation in the spectral space; similarly, the normal vector (0.18, 0.32, 0.93) can also be represented in the spectral space.

[0042] The coordinate frequency domain coefficients and normal frequency domain coefficients are obtained by minimizing the reconstruction error. Reconstruction error refers to the error incurred when reconstructing the original data using partial eigenvectors. Through optimized calculation, a set of coordinate and normal frequency domain coefficients can be obtained that minimizes the reconstruction error. For example, the coordinate frequency domain coefficients might yield a vector containing 1000 values, where the first 100 coefficients may have larger amplitudes, such as 0.5632, 0.4215, and 0.3987, while the subsequent coefficients have smaller amplitudes. Similarly, the normal frequency domain coefficients are also a vector containing 1000 values, reflecting the distribution characteristics of the normal vector in the spectral domain.

[0043] The gradient correlation between the coordinate frequency domain coefficients and the normal frequency domain coefficients is calculated. Gradient correlation reflects the degree of association between coordinate changes and normal changes. For example, calculating the correlation coefficient between the gradients of the coordinate frequency domain coefficients and the gradients of the normal frequency domain coefficients yields a correlation vector with element values ​​ranging from -1 to 1, such as 0.8756, 0.7632, and 0.6589. Based on the gradient correlation, the amplitudes of the coordinate and normal frequency domain coefficients are adaptively weighted to obtain the spectral energy distribution. Adaptive weighting is a process of dynamically adjusting the weights according to the degree of correlation. For example, when the correlation is high (e.g., 0.8756), a larger weight (e.g., 0.9) can be assigned; when the correlation is low (e.g., 0.3254), a smaller weight (e.g., 0.4) can be assigned. Through weighted calculation, the spectral energy distribution can be obtained, with energy values ​​of 23.56, 18.92, and 15.67, respectively.

[0044] Multi-scale wavelet decomposition is used to extract energy abrupt change points in the spectral energy distribution. Multi-scale wavelet decomposition is a method for analyzing the characteristics of signal changes at different scales. Through wavelet transform, abrupt change points in the spectral energy distribution can be found, that is, locations where the energy value changes significantly. For example, after performing a 4-level wavelet decomposition, it may be found that the energy value drops sharply from 15.67 to 5.23 at the 237th frequency point; this point is the energy abrupt change point. The spectral energy value corresponding to the energy abrupt change point is used as a coefficient threshold. For example, if the energy value at the 237th frequency point is 5.23, then 5.23 can be set as the coefficient threshold.

[0045] Coefficients in the coordinate frequency domain with spectral energy exceeding a coefficient threshold are identified as high-frequency coefficients. For example, if the coefficient threshold is 5.23, then frequency domain coefficients with energy values ​​greater than 5.23 (such as 23.56, 18.92, 15.67, etc.) are identified as high-frequency coefficients. Corresponding feature vectors are extracted to construct a high-spectral basis. For the identified high-frequency coefficients, their corresponding feature vectors are extracted, such as feature vector 5, feature vector 12, feature vector 28, etc., to form the high-spectral basis. The high-frequency coefficients are then spectrally back-projected onto the high-spectral basis. Spectral back-projection is the process of converting frequency domain information back to the spatial domain. For example, a linear combination of the high-frequency coefficients and their corresponding feature vectors yields a spatial representation of the high-frequency components.

[0046] The coordinates of high-frequency vertices are obtained by optimizing and reconstructing them using pre-defined local geometric constraints. These constraints may include conditions such as curvature continuity and normal continuity. Through optimization calculations, high-frequency vertex coordinates that satisfy the constraints can be obtained. For example, the original coordinates of a vertex are (11.25, 16.34, 9.78), and the reconstructed coordinates after high-frequency optimization may be (11.27, 16.38, 9.82), preserving high-frequency details while satisfying the geometric constraints.

[0047] Traverse all triangular patches and identify those containing high-frequency vertex coordinates. For example, if any of the three vertices of a triangular patch is a high-frequency vertex, that patch requires special handling. Calculate the corresponding discrete Gaussian curvature and mean curvature. Discrete Gaussian curvature and mean curvature are important indicators describing the local geometric properties of a surface. For triangular meshes, these curvature values ​​can be calculated using discrete differential geometry methods. For example, a triangular patch might have a discrete Gaussian curvature of 0.0235 and a mean curvature of 0.0568.

[0048] A curvature tensor is constructed based on discrete Gaussian curvature and mean curvature, and the principal curvatures are calculated. The curvature tensor is a second-order tensor describing the degree of curvature of a surface in various directions. The principal curvatures, i.e., the maximum and minimum curvatures of the surface at a given point, can be calculated from the curvature tensor. For example, the principal curvatures of a triangular facet might be 0.0823 and 0.0314, representing the degree of curvature of the surface at that point in two orthogonal directions, respectively. Triangular facets with principal curvatures exceeding a preset curvature threshold are marked as high-frequency triangular facets. For example, if the preset curvature threshold is 0.05, then a principal curvature of 0.0823 exceeds the threshold, and this triangular facet will be marked as a high-frequency triangular facet, requiring special handling during reconstruction to preserve high-frequency details.

[0049] In this embodiment, by constructing a normalized graph Laplacian matrix based on vertex adjacency and performing spectral decomposition, the geometric information that was originally strongly dependent on spatial location is mapped to a spectral domain consistent with the mesh topology. This makes the analysis of geometric features no longer directly affected by uneven mesh resolution or local topological complexity, improving the consistency of the expression of overall geometric structure changes. By jointly embedding vertex coordinates and vertex normals in the spectral domain and using reconstruction error constraints to calculate frequency domain coefficients, geometric shape changes and normal changes are collaboratively characterized under the same analysis framework. This avoids the false detection or false negative problem caused by high-frequency feature extraction based solely on coordinates or solely on normal information, improving the completeness and reliability of high-frequency geometric feature recognition. By introducing the gradient correlation between coordinate frequency domain coefficients and normal frequency domain coefficients and adaptively weighting the spectral energy, the determination of high-frequency energy can simultaneously reflect the intensity of shape changes and the consistency of normal changes, effectively enhancing the ability to distinguish between real geometric abrupt changes and noise disturbances.

[0050] In one alternative implementation, The principal curvature direction of the parametric coordinate region corresponding to the high-frequency triangular patch is calculated based on the curvature tensor components. The insertion position and additional control points are determined based on the principal curvature direction and merged with the control points of the surface patch to obtain merged control points, including: Extract the parameter coordinate region corresponding to the high-frequency triangular patch, obtain the curvature tensor components of each parameter grid point in the parameter coordinate region and perform feature decomposition to obtain eigenvalues ​​and eigenvectors, take the eigenvector corresponding to the maximum eigenvalue as the maximum principal curvature direction, and take the eigenvector corresponding to the minimum eigenvalue as the minimum principal curvature direction. Curvature change rate fields are constructed along the direction of maximum principal curvature and the direction of minimum principal curvature, respectively. Morphological gradient operations are performed on the curvature change rate fields to identify extreme points of curvature change rate. The parameter coordinates of the extreme points of curvature change rate in the parameter coordinate region are used as the insertion positions. At the insertion position, the covariant derivatives of the curvature tensor are calculated along the direction of maximum principal curvature and the direction of minimum principal curvature, respectively. A local curvature compensation vector is constructed based on the covariant derivatives. The parameter coordinates of the insertion position are offset along the local curvature compensation vector to obtain the parameter coordinates of the additional control point. The spatial coordinates of the additional control points are obtained by performing surface mapping calculation on the parameter coordinates of the additional control points, and the surface patch control points are obtained. The spatial coordinates of the additional control points and the surface patch control points are merged in the order of parameter coordinates, and the merged control points are topologically reorganized to obtain the merged control points.

[0051] Extract the parameter coordinate region corresponding to the high-frequency triangular facet. A high-frequency triangular facet refers to a triangular facet marked as containing high-frequency geometric features in the aforementioned embodiments. For a high-frequency triangular facet, it is necessary to extract the corresponding region in the parameter space. For example, if the coordinates of the three vertices of a certain high-frequency triangular facet in three-dimensional space are (11.27, 16.38, 9.82), (11.56, 16.42, 9.95), and (11.35, 16.63, 9.89), the corresponding parameter coordinate region is the triangular region in the parameter space formed by the points (0.32, 0.45), (0.38, 0.46), and (0.34, 0.51). Obtain the curvature tensor components of each parameter grid point within the parameter coordinate region. Set uniform grid points in the parameter space, for example, setting grid points within the parameter coordinate region with a grid spacing of 0.01, and calculate the curvature tensor components at each grid point. The curvature tensor is a second-order tensor. For a point (0.35, 0.47) in the parameter space, its curvature tensor components may have four values: 0.0723, 0.0256, 0.0256 and 0.0512, which correspond to the four elements of the tensor, respectively.

[0052] Eigenvalues ​​and eigenvectors are obtained by performing eigenvalue decomposition on the components of the curvature tensor. Eigenvalue decomposition is the process of decomposing a tensor into eigenvalues ​​and eigenvectors. For example, for the aforementioned curvature tensor, eigenvalue decomposition yields two eigenvalues, 0.0823 and 0.0412, and corresponding eigenvectors (0.9325, 0.3612) and (-0.3612, 0.9325). The eigenvector corresponding to the largest eigenvalue is taken as the direction of maximum principal curvature. In the example above, the eigenvector (0.9325, 0.3612) corresponding to the eigenvalue 0.0823 is the direction of maximum principal curvature, indicating that the curvature of the surface is greatest in this direction. The eigenvector corresponding to the smallest eigenvalue is taken as the direction of minimum principal curvature. The eigenvector (-0.3612, 0.9325) corresponding to the eigenvalue 0.0412 is the direction of minimum principal curvature, indicating that the curvature of the surface is least in this direction.

[0053] Curvature change rate fields are constructed along the directions of maximum and minimum principal curvature, respectively. The curvature change rate refers to the degree of change in curvature along a certain direction. In parameter space, the curvature change rate along the directions of maximum and minimum principal curvature can be calculated at each grid point. For example, at the point (0.35, 0.47), the curvature change rate along the direction of maximum principal curvature might be 0.0156, and the curvature change rate along the direction of minimum principal curvature might be 0.0082. Morphological gradient operations are performed on the curvature change rate fields to identify extreme points of the curvature change rate. Morphological gradient operations are a method used for edge detection in image processing, used to identify extreme points in the curvature change rate field. By applying the morphological gradient operator to the curvature change rate field, the locations where the change rate changes most significantly can be found. For example, it might be found that the curvature change rate reaches a local maximum of 0.0189 at the point (0.36, 0.48).

[0054] The parametric coordinates of the extreme points of the rate of change of curvature within the parametric coordinate region are used as the insertion positions. In the example above, the point (0.36, 0.48) is determined as the location for inserting additional control points. At the insertion positions, the covariant derivatives of the curvature tensor are calculated along the directions of maximum and minimum principal curvature, respectively. The covariant derivatives take into account the intrinsic geometric properties of the surface and can accurately describe the change in curvature on the surface. For example, at the point (0.36, 0.48), the covariant derivative along the direction of maximum principal curvature might be (0.0223, 0.0087), and the covariant derivative along the direction of minimum principal curvature might be (0.0065, 0.0142).

[0055] A local curvature compensation vector is constructed based on covariant derivatives. This vector is used to adjust the position of control points to adapt to local curvature variations. The local curvature compensation vector can be constructed through a linear combination of covariant derivatives. For example, for the point (0.36, 0.48), the calculated local curvature compensation vector might be (0.0187, 0.0105). The parameter coordinates of the insertion position are then offset along the local curvature compensation vector to obtain the parameter coordinates of the additional control point. For example, offsetting the point (0.36, 0.48) along the vector (0.0187, 0.0105) yields the parameter coordinates of the additional control point (0.3787, 0.4905).

[0056] The spatial coordinates of the additional control points are obtained by performing surface mapping calculations on their parametric coordinates. Surface mapping is the process of mapping points in parametric space to three-dimensional space. Using the surface equation, the parametric coordinates (0.3787, 0.4905) can be mapped to three-dimensional spatial coordinates (11.42, 16.51, 9.88). The surface patch control points are then obtained. These control points are the set of control points determined by the aforementioned method and are used to define the shape of the surface. For example, a surface patch may be defined by 16 control points, one of which has spatial coordinates of (11.35, 16.45, 9.85).

[0057] The spatial coordinates of the additional control points are merged with those of the surface patch control points in parametric coordinate order. The merging process must consider the distribution order of the control points in the parameter space to ensure the merged control point mesh maintains topological consistency. For example, inserting the additional control point (11.42, 16.51, 9.88) into the original 16 control points results in a set of 17 control points. Topological reorganization of the merged control points yields the merged control points. Topological reorganization refers to re-establishing the connections between control points based on their spatial distribution, forming an effective control point mesh structure. For example, it may be necessary to reorganize the original 4×4 control point mesh into a non-uniform control point mesh, where the control point density is higher in high-curvature regions.

[0058] In this embodiment, by performing feature decomposition on the curvature tensor and analyzing the rate of curvature change along the directions of maximum and minimum principal curvature respectively, the dominant direction of local geometric changes on the surface can be accurately identified, improving the perception of the anisotropic features of the surface and more accurately reflecting the true surface morphology. By constructing a rate of curvature change field along the principal curvature direction and automatically identifying the extreme points of the rate of curvature change using morphological gradient operations, the control point insertion position is adaptively matched with the region where curvature changes are most significant, improving the pertinence and effectiveness of the control point distribution for geometric features. By introducing the covariant derivative of the curvature tensor at the insertion position and constructing a local curvature compensation vector, the position of the newly added control point is adjusted in terms of direction and amplitude compensation, so that the additional control point is not only located in the high curvature region, but also conforms to the curvature change trend of the surface in terms of position, effectively reducing the geometric deviation caused by abrupt changes in local curvature and improving the accuracy and smoothness of surface reconstruction.

[0059] In one alternative implementation, Based on the merged control points and the parameter coordinates, a parameterized surface patch is constructed and the surface sampling points are determined. The normal projection distance from the surface sampling points to the surface of the triangular patch is calculated to obtain the normal deviation and the deviation spectrum is determined, including: Obtain the spatial coordinates of the merged control points and the parameter boundary of the parameter coordinate region. Construct a T-spline basis function based on the spatial coordinates of the merged control points. Insert node vectors into the parameter coordinate region and generate a non-uniform parameter grid using a local refinement algorithm. Calculate the spatial coordinates of the surface sampling points based on the nodes in the non-uniform parameter grid and the T-spline basis function. An octree spatial index is constructed based on the spatial coordinates of the surface sampling point, and the nearest triangular facet is queried. The normal vector of the tangent plane of the parameterized surface facet at the surface sampling point is calculated. The surface sampling point is projected onto the surface of the nearest triangular facet along the direction of the tangent plane normal vector to obtain the projection point. The directed distance between the surface sampling point and the projection point is calculated as the normal projection distance. The parameter coordinates corresponding to the nodes in the non-uniform parameter grid are used as the sampling positions, and the normal projection distance is used as the deviation amplitude at the sampling positions. A non-uniform sampling deviation field is constructed based on the sampling positions and the deviation amplitude. The non-uniform sampling deviation field is subjected to a non-uniform fast Fourier transform to obtain a frequency domain representation, and the power spectral density of each frequency component in the frequency domain representation is extracted to obtain the deviation spectrum.

[0060] Obtain the spatial coordinates of the merged control points and the parameter boundaries of the parameter coordinate region. The merged control points are the set of control points obtained by inserting additional control points and performing topological reorganization in the aforementioned method. For example, a set of merged control points may contain 25 spatial points, one of which has spatial coordinates of (11.42, 16.51, 9.88) and corresponding parameter coordinates of (0.3787, 0.4905). The parameter boundaries of the parameter coordinate region define the range of the surface patch to be reconstructed in the parameter space; for example, the parameter range is from 0.25 to 0.45 in the horizontal direction and from 0.40 to 0.60 in the vertical direction. Construct T-spline basis functions based on the spatial coordinates of the merged control points. T-splines are a high-order continuous parametric surface representation method. Compared to traditional non-uniform rational B-splines (NURBS), T-splines allow the insertion of control points in local regions without affecting the overall structure. When constructing T-spline basis functions, a local node vector needs to be specified for each control point. For example, for the control point (11.42, 16.51, 9.88), its corresponding local node vector may have horizontal node values ​​of 0.32, 0.36, 0.40, 0.44, and 0.48, and vertical node values ​​of 0.42, 0.46, 0.50, 0.54, and 0.58. These node values ​​define the support region of the T-spline basis function.

[0061] A non-uniform parametric mesh is generated by inserting node vectors within the parametric coordinate region using a local refinement algorithm. The algorithm dynamically adjusts the node density based on the surface complexity, increasing density in areas of high curvature and decreasing density in smooth areas. For example, in the high curvature region near the parametric coordinates (0.36, 0.48), the node density can be increased with a node spacing of 0.005; while in the smooth region, the node spacing can be set to 0.02. The spatial coordinates of the surface sampling points are calculated based on the nodes within the non-uniform parametric mesh and the T-spline basis functions. For each node in the non-uniform parametric mesh, its corresponding three-dimensional spatial coordinates are calculated using the T-spline basis functions. For example, for a mesh node with parametric coordinates (0.35, 0.47), the calculated spatial coordinates of the surface sampling point are (11.38, 16.48, 9.86).

[0062] An octree spatial index is constructed based on the spatial coordinates of surface sampling points. An octree is a spatial partitioning data structure that efficiently organizes and queries points in three-dimensional space. When constructing the octree, a bounding box is defined for the space, for example, a horizontal range of 10.00 to 12.00, a vertical range of 15.00 to 17.00, and a height range of 9.00 to 10.50. This space is recursively divided into eight subspaces until a preset precision is reached or other termination conditions are met. Using the octree structure, the nearest triangle to a surface sampling point can be quickly queried. For example, for the surface sampling point (11.38, 16.48, 9.86), the nearest triangle might be composed of vertices (11.35, 16.45, 9.85), (11.40, 16.50, 9.87), and (11.37, 16.51, 9.84).

[0063] Calculate the normal vector of the tangent plane of the parametric surface patch at the surface sampling point. The tangent plane normal vector is a unit vector describing the direction of the normal to the surface at a point, and can be obtained by calculating the cross product of the tangent vectors of the parametric surface at that point. For example, at the surface sampling point (11.38, 16.48, 9.86), the calculated tangent plane normal vector is (0.21, 0.34, 0.92). Project the surface sampling point along the direction of the tangent plane normal vector onto the nearest triangular patch surface to obtain the projection point. The projection process involves finding the intersection point with the triangular patch along the direction of the normal vector. For example, projecting the surface sampling point (11.38, 16.48, 9.86) along the direction of the normal vector (0.21, 0.34, 0.92) onto the nearest triangular patch yields the projection point (11.39, 16.50, 9.88). Calculate the directed distance between the surface sampling point and the projection point as the normal projection distance. Directed distance takes directionality into account. If the projection point is in the positive direction of the normal vector of the surface sampling point, the distance is positive; otherwise, it is negative. For example, the directed distance between the surface sampling point (11.38, 16.48, 9.86) and the projection point (11.39, 16.50, 9.88) may be 0.0282 mm.

[0064] The parameter coordinates corresponding to nodes within a non-uniform parameter mesh are used as sampling locations. For example, a non-uniform parameter mesh may contain 400 nodes, each with corresponding parameter coordinates, such as node 25 with parameter coordinates of (0.35, 0.47). The normal projection distance is used as the deviation magnitude at the sampling location. For example, at the sampling location with parameter coordinates of (0.35, 0.47), the deviation magnitude is 0.0282 mm. A non-uniform sampling deviation field is constructed based on the sampling location and the deviation magnitude. The non-uniform sampling deviation field is a two-dimensional function, defined in the parameter space, with the range being the normal projection distance. For example, the value at parameter coordinates (0.35, 0.47) in the deviation field is 0.0282 mm, and the value at parameter coordinates (0.38, 0.52) is -0.0173 mm.

[0065] A non-uniform fast Fourier transform (NFT) is performed on the non-uniform sampling bias field to obtain its frequency domain representation. The NFT is an extension of the Fourier transform, applicable to non-uniformly sampled data. Through the NFT, a non-uniform sampling bias field in the spatial domain can be converted into a frequency domain representation. For example, performing a NFT on a non-uniform sampling bias field containing 400 sampling points yields a 20×20 frequency domain coefficient matrix, where the complex coefficients at frequency point (5, 8) are 0.0152 + 0.0083i. The bias spectrum is obtained by extracting the power spectral density of each frequency component in the frequency domain representation. The power spectral density describes the distribution of signal energy in the frequency domain. The power spectral density of each frequency component can be obtained by calculating the square of the modulus of the frequency domain coefficients. For example, the power spectral density at frequency point (5, 8) is 0.000302, and the power spectral density at frequency point (10, 12) is 0.000125.

[0066] In this embodiment, by constructing T-spline basis functions based on merged control points and performing local refinement within the parametric coordinate region, a non-uniform parametric mesh adaptively matches the geometric complexity of the surface is generated. This enables surface sampling to achieve higher resolution in high-frequency detail regions, effectively reducing redundant computation and improving the ability to express local geometric features. By utilizing octree spatial indexing to quickly locate the nearest triangular facet corresponding to the surface sampling point and projecting it along the normal vector direction of the parametric surface tangent plane, a normal projection distance with clear geometric meaning is obtained, improving the consistency and comparability between error evaluation results and actual geometric deviations. By constructing a non-uniform sampling deviation field at the nodes of the non-uniform parametric mesh, where sampling positions and deviation amplitudes correspond one-to-one, discrete surface deviations are transformed into a field representation in the continuous parameter domain. This allows error distributions at different locations and scales to be comprehensively analyzed in the same parameter domain, significantly improving the diagnostic capability of error analysis and providing a quantitative basis for subsequent targeted adjustments to control point density, spline order, or local refinement strategies.

[0067] In one alternative implementation, Iteratively updating the merged control points based on the amplitude distribution of the deviation spectrum until the normal deviation converges, obtaining converged control points, and generating a reconstructed surface by combining the parameter coordinates includes: The amplitude of each frequency component in the deviation spectrum is extracted and frequency-weighted to obtain the frequency domain deviation energy. A control point sensitivity matrix is ​​constructed based on the frequency domain deviation energy. The dominant sensitivity direction is obtained by performing singular value decomposition on the control point sensitivity matrix. The position adjustment amount of the merged control point is calculated along the dominant sensitivity direction. The spatial coordinates of the merged control point are updated along the position adjustment amount to obtain the optimized merged control point. Based on the optimized merged control points and the parameter coordinates, the parameterized surface patch is reconstructed and the updated surface sampling points are determined. The normal projection distance from the updated surface sampling points to the surface of the triangular patch is calculated to obtain the optimized normal deviation. Calculate the rate of change of the deviation between the optimized normal deviation and the normal deviation before the update, and determine whether the rate of change of the deviation is less than a preset convergence threshold. If it is less, the optimized merged control point is taken as the convergence control point. If it is not less, the optimized merged control point is taken as the new merged control point, and the process of constructing the control point sensitivity matrix and calculating the normal deviation is repeated until the rate of change of the deviation is less than the convergence threshold. Based on the convergence control point and the parameter coordinates, the parameterized surface patch is reconstructed to obtain the reconstructed surface.

[0068] The amplitude of each frequency component in the deviation spectrum is extracted and frequency-weighted to obtain the frequency domain deviation energy. The deviation spectrum is the frequency domain representation obtained by performing a non-uniform fast Fourier transform on the non-uniform sampling deviation field in the aforementioned method. For example, for a 20×20 frequency domain coefficient matrix, the power spectral density at frequency point (5, 8) is 0.000302, and the power spectral density at frequency point (10, 12) is 0.000125. Frequency weighting is the process of assigning different weights to each frequency component according to the importance of the frequency. Generally, high-frequency components correspond to detailed features, and low-frequency components correspond to the overall shape, so different weights can be set according to application requirements. For example, the weight of the low-frequency region (frequency range 0 to 5) can be set to 0.3, the weight of the mid-frequency region (frequency range 6 to 12) to 0.5, and the weight of the high-frequency region (frequency range 13 to 20) to 0.2. The frequency domain bias energy can be obtained by multiplying the power spectral density of each frequency component by its corresponding weight and summing the results. For example, the calculated frequency domain bias energy is 0.000185.

[0069] A control point sensitivity matrix is ​​constructed based on the frequency domain deviation energy. The control point sensitivity matrix describes the degree to which changes in the control point position affect the frequency domain deviation energy. For a surface containing 25 control points, a 25×3 sensitivity matrix is ​​constructed, where each row corresponds to one control point and each column corresponds to the three components of the control point's coordinates. The elements of the sensitivity matrix can be calculated using the perturbation method: a small perturbation (e.g., 0.001 mm) is applied to each coordinate component of each control point, and the rate of change of the frequency domain deviation energy before and after the perturbation is calculated. For example, changing the lateral coordinate of control point 5 from 11.42 to 11.421 results in a frequency domain deviation energy of 0.000182, with a rate of change of -0.0162. Similarly, the sensitivity vector corresponding to control point 5 is calculated as (-0.0162, 0.0253, 0.0127), representing the sensitivity of this control point in the three coordinate directions.

[0070] Singular value decomposition (SVD) is performed on the control point sensitivity matrix to obtain the dominant sensitivity direction. SVD is a matrix factorization method that decomposes a matrix into the product of three matrices to extract its main features. SVD of a 25×3 sensitivity matrix yields three singular values: 0.0412, 0.0287, and 0.0153, along with their corresponding right singular vectors. The right singular vector (0.5731, 0.6824, 0.4536) corresponding to the largest singular value of 0.0412 is the dominant sensitivity direction, indicating that adjusting the control point position along this direction is most effective in reducing deviation energy. The position adjustment amount for merging control points is calculated along the dominant sensitivity direction. The calculation of the position adjustment amount needs to consider the sensitivity magnitude and the adjustment step size. The adjustment step size can be determined using a line search method or set to a fixed value. For example, if the basic adjustment step size is set to 0.05 mm, then for the sensitivity vector (-0.0162, 0.0253, 0.0127), the corresponding position adjustment amount is (-0.0008, 0.0013, 0.0006) mm.

[0071] The spatial coordinates of the merged control points are updated along the position adjustment to obtain optimized merged control points. For example, the original coordinates of control point 5 are (11.42, 16.51, 9.88), and the adjusted coordinates are (11.4192, 16.5113, 9.8806). Similar adjustments are made to all control points to obtain a set of optimized control point coordinates. Based on the optimized merged control points and parameter coordinates, a parametric surface patch is reconstructed, and the updated surface sampling points are determined. The T-spline surface is reconstructed using the optimized control points, and surface sampling points are calculated at the nodes of the non-uniform parametric mesh. For example, for the mesh node with parameter coordinates (0.35, 0.47), the spatial coordinates of the surface sampling point calculated using the optimized control points are (11.3802, 16.4803, 9.8604).

[0072] The optimized normal deviation is obtained by calculating the normal projection distance from the updated surface sampling point to the triangular facet surface. This step is similar to the previous method, calculating the projection distance from the surface sampling point along the normal vector direction to the nearest triangular facet. For example, for the updated surface sampling point (11.3802, 16.4803, 9.8604), the calculated normal projection distance is 0.0251 mm, which is reduced compared to 0.0282 mm before the update. Similar calculations are performed for all sampling points to obtain a set of optimized normal deviation values. The rate of change of deviation between the optimized normal deviation and the normal deviation before the update is calculated. The rate of change of deviation can be measured by calculating the relative change of the root mean square value of the normal deviation. For example, the root mean square value of the normal deviation before the update is 0.0215 mm, and after the update it is 0.0194 mm, the rate of change of deviation is (0.0194-0.0215) / 0.0215=-0.0977, that is, a reduction of 9.77%.

[0073] The algorithm determines whether the rate of change of deviation is less than a preset convergence threshold. The convergence threshold is a parameter set according to accuracy requirements; for example, it can be set to 0.01 or 0.005, indicating that the optimization process is considered converged when the absolute value of the rate of change of deviation is less than this threshold. In the previous example, the rate of change of deviation is -0.0977, and the absolute value is 0.0977, which is greater than the assumed convergence threshold of 0.01, therefore further iterations are needed. If the rate of change of deviation is less than the convergence threshold, the optimized merged control point is used as the convergence control point; if it is not less, the optimized merged control point is used as the new merged control point, and the process of constructing the control point sensitivity matrix and calculating the normal deviation is repeated until the rate of change of deviation is less than the convergence threshold. In practical applications, multiple iterations may be required to reach the convergence condition. For example, after 5 iterations, the rate of change of deviation may drop to 0.0092, which is less than the convergence threshold of 0.01. In this case, the control point obtained in the 5th iteration can be used as the convergence control point.

[0074] Based on the convergence control points and parameter coordinates, a parametric surface patch is reconstructed to obtain the reconstructed surface. The optimized converged control points are then used to construct a T-spline surface that meets the high-order continuity requirements. For example, using 25 convergence control points, one of which has coordinates (11.4173, 16.5095, 9.8789), a third-order T-spline surface is constructed, resulting in a reconstructed surface with second-order geometric continuity. The maximum deviation between this surface and the original STL model is controlled within 0.02 mm, meeting the requirements of high-precision engineering applications.

[0075] In this embodiment, by weighting each frequency component in the deviation spectrum and constructing the frequency domain deviation energy, the contribution of surface error at different scales can be quantified and explicitly introduced into the control point optimization process, improving the targeting of the optimization direction to the real geometric deviation. By constructing the control point sensitivity matrix based on the frequency domain deviation energy and extracting the dominant sensitivity direction through singular value decomposition, the complex coupling relationship between multiple control points and multi-scale errors is compressed into a small number of dominant adjustment directions, improving the numerical stability and convergence efficiency of the control point update process. By performing overall coordinated adjustment of the merged control points along the dominant sensitivity direction, the control point position update has global consistency, reducing the damage to the overall continuity of the surface caused by local overcorrection, and effectively reducing the risk of getting trapped in local optima or generating geometric distortion.

[0076] Figure 2 This is a flowchart of the frequency domain deviation optimization for surface reconstruction in the STL surface reconstruction process according to an embodiment of the present invention.

[0077] A second aspect of the present invention provides a high-order continuity preservation and error control system for STL surface reconstruction, comprising: The information extraction unit is used to acquire discretized geometric data and identify the shared edges and vertices of triangular facets to obtain connection information. Based on the connection information, the vertex coordinates of the common boundary on adjacent triangular facets are determined and projected onto the tangent plane of the corresponding triangular facets to obtain parameter coordinates. The normal field construction unit is used to extract the vertex normal vectors on both sides of the common boundary, perform spherical harmonic decomposition to obtain spherical harmonic coefficients, construct a continuous normal vector field based on the spherical harmonic coefficients and the parameter coordinates, and perform path integration to obtain the surface patch control points and curvature tensor components. The high-frequency identification unit is used to perform frequency domain transformation on the vertex coordinates and the vertex normal vector to obtain frequency domain coefficients and determine high-frequency coefficients by combining them with a preset coefficient threshold. The high-frequency coefficients are then inversely transformed to obtain high-frequency vertex coordinates and determine high-frequency triangular facets. The control point solving unit is used to calculate the principal curvature direction of the parameter coordinate region corresponding to the high-frequency triangular patch based on the curvature tensor components, determine the insertion position and additional control points based on the principal curvature direction, and merge them with the control points of the surface patch to obtain the merged control points. The surface reconstruction unit is used to construct a parameterized surface patch based on the merged control points and the parameter coordinates, determine the surface sampling points, calculate the normal projection distance from the surface sampling points to the surface of the triangular patch to obtain the normal deviation and determine the deviation spectrum, iteratively update the merged control points according to the amplitude distribution of the deviation spectrum until the normal deviation converges, obtain the converged control points, and generate the reconstructed surface by combining the parameter coordinates.

[0078] A third aspect of the present invention provides an electronic device, comprising: A processor and a memory for storing processor-executable instructions, wherein the processor is configured to invoke instructions stored in the memory to perform the aforementioned method.

[0079] A fourth aspect of the present invention provides a computer-readable storage medium having stored thereon computer program instructions that, when executed by a processor, implement the aforementioned method.

[0080] This invention can be a method, apparatus, system, and / or computer program product. The computer program product may include a computer-readable storage medium having computer-readable program instructions loaded thereon for performing various aspects of the invention.

[0081] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, and not to limit them; although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some or all of the technical features; and these modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the scope of the technical solutions of the embodiments of the present invention.

Claims

1. A method for preserving high-order continuity and controlling errors in STL surface reconstruction, characterized in that, include: Discretized geometric data is obtained and the shared edges and vertices of triangular facets are identified to obtain connection information. Based on the connection information, the vertex coordinates of the common boundary on adjacent triangular facets are determined and projected onto the tangent plane of the corresponding triangular facets to obtain parameter coordinates. Extract the vertex normal vectors on both sides of the common boundary and perform spherical harmonic decomposition to obtain spherical harmonic coefficients. Based on the spherical harmonic coefficients and the parameter coordinates, construct a continuous normal vector field and perform path integration to obtain the surface patch control points and curvature tensor components. The vertex coordinates and the vertex normal vector are transformed in the frequency domain to obtain frequency domain coefficients. High-frequency coefficients are determined by combining the frequency domain coefficients with a preset coefficient threshold. The high-frequency coefficients are then transformed inversely to obtain high-frequency vertex coordinates and determine high-frequency triangular facets. Calculate the principal curvature direction of the parameter coordinate region corresponding to the high-frequency triangular patch based on the curvature tensor components, determine the insertion position and additional control points based on the principal curvature direction, and merge them with the control points of the surface patch to obtain the merged control points; Based on the merged control points and the parameter coordinates, a parameterized surface patch is constructed and the surface sampling points are determined. The normal projection distance from the surface sampling points to the surface of the triangular patch is calculated to obtain the normal deviation and the deviation spectrum is determined. The merged control points are iteratively updated according to the amplitude distribution of the deviation spectrum until the normal deviation converges, and the converged control points are obtained. The reconstructed surface is generated by combining the parameter coordinates.

2. The method according to claim 1, characterized in that, Discretized geometric data is acquired, and shared edges and vertices of triangular facets are identified to obtain connection information. Based on this connection information, the vertex coordinates of the common boundary on adjacent triangular facets are determined and projected onto the tangent plane of the corresponding triangular facets to obtain parameter coordinates, including: Obtain discrete geometric data containing multiple triangular facets and traverse the vertex coordinates of each triangular facet in the discrete geometric data. Calculate the vertex overlap between any two triangular facets. Identify triangular facet pairs with vertex overlap greater than zero as adjacent triangular facets and extract the overlapping vertices of the adjacent triangular facets as shared vertices. Based on the shared vertices, the edges of the adjacent triangular faces are traversed, and edges containing the same pair of shared vertices are identified as shared edges. An index relationship between the shared edges and the adjacent triangular faces is established to obtain the connection information. Based on the connection information, the shared edge of each pair of adjacent triangular facets is extracted, the shared vertex on the shared edge is determined as the common boundary vertex, and the vertex coordinates of the common boundary vertex are extracted as the common boundary vertex coordinates. Calculate the centroid coordinates of the three vertices of the triangular facet, construct a tangent plane with the centroid coordinates as the origin based on the vertex normal vectors of the triangular facet, project the coordinates of the common boundary vertices onto the tangent plane to obtain the projected coordinates, construct a two-dimensional coordinate system on the tangent plane, and represent the projected coordinates in the two-dimensional coordinate system to obtain the parameter coordinates.

3. The method according to claim 1, characterized in that, Extracting the vertex normal vectors on both sides of the common boundary and performing spherical harmonic decomposition yields spherical harmonic coefficients. Based on the spherical harmonic coefficients and the parameter coordinates, a continuous normal vector field is constructed and path integration is performed to obtain the surface patch control points and curvature tensor components, including: Based on the connection information, the adjacent triangular facets on both sides of the common boundary are identified, the vertex normal vectors of the adjacent triangular facets are extracted as the vertex normal vectors of the common boundary, and the vertex normal vectors of the common boundary are mapped to a unit sphere to obtain spherical coordinates. The spherical coordinates are expanded using a spherical harmonic function, and the function value of the spherical harmonic function at the spherical coordinates is calculated. Based on the function value and the normal vector of the common boundary vertex, the least squares fitting is performed to obtain the spherical harmonic coefficients of each order. A set of grid points is established based on the parameter coordinates. For each grid point, the reconstructed value of the spherical harmonic function is calculated based on the spherical harmonic coefficient. The reconstructed value of the spherical harmonic function is mapped from the unit sphere to the three-dimensional space to obtain the grid point normal vector. The continuous normal vector field is constructed based on the grid point normal vector. Extract the normal vector components of the continuous normal vector field along the parameter coordinate direction, perform line integration on the normal vector components along the corresponding parameter coordinate direction to obtain the integral curve, and determine the control points of the surface patch based on the intersection coordinates of the integral curves in different parameter coordinate directions. Calculate the partial derivatives of the continuous normal vector field with respect to each parameter coordinate direction, and construct the curvature tensor components based on the partial derivatives.

4. The method according to claim 1, characterized in that, Performing a frequency domain transformation on the vertex coordinates and the vertex normal vector to obtain frequency domain coefficients, and determining high-frequency coefficients by combining them with a preset coefficient threshold, and then performing an inverse transformation on the high-frequency coefficients to obtain high-frequency vertex coordinates and determine high-frequency triangular patches includes: The vertex coordinates of all triangular facets in the discretized geometric data are extracted and combined with the connection information to construct a vertex adjacency matrix. Based on the vertex adjacency matrix, a normalized graph Laplacian matrix is ​​constructed and spectral decomposition is performed to obtain a set of eigenvectors. The vertex coordinates and the vertex normal vectors are jointly embedded in the spectral space formed by the set of eigenvectors. The coordinate frequency domain coefficients and normal frequency domain coefficients are calculated by minimizing the reconstruction error. The gradient correlation between the coordinate frequency domain coefficients and the normal frequency domain coefficients is calculated, and the magnitudes of the coordinate frequency domain coefficients and the normal frequency domain coefficients are adaptively weighted based on the gradient correlation to obtain the spectral energy distribution. The spectral energy distribution is then subjected to multi-scale wavelet decomposition to extract energy abrupt change points, and the spectral energy value corresponding to the energy abrupt change point is used as the coefficient threshold. The coordinate frequency domain coefficients whose spectral energy exceeds the coefficient threshold are identified as the high-frequency coefficients, and the corresponding feature vectors are extracted to construct a high-frequency basis. The high-frequency coefficients are then back-projected onto the high-frequency basis and optimized and reconstructed in combination with preset local geometric constraints to obtain the high-frequency vertex coordinates. Traverse all triangular facets, identify triangular facets containing the high-frequency vertex coordinates and calculate the corresponding discrete Gaussian curvature and mean curvature, construct a curvature tensor based on the discrete Gaussian curvature and mean curvature and calculate the principal curvature, and mark triangular facets with principal curvature exceeding a preset curvature threshold as the high-frequency triangular facets.

5. The method according to claim 1, characterized in that, The principal curvature direction of the parametric coordinate region corresponding to the high-frequency triangular patch is calculated based on the curvature tensor components. The insertion position and additional control points are determined based on the principal curvature direction and merged with the control points of the surface patch to obtain merged control points, including: Extract the parameter coordinate region corresponding to the high-frequency triangular patch, obtain the curvature tensor components of each parameter grid point in the parameter coordinate region and perform feature decomposition to obtain eigenvalues ​​and eigenvectors, take the eigenvector corresponding to the maximum eigenvalue as the maximum principal curvature direction, and take the eigenvector corresponding to the minimum eigenvalue as the minimum principal curvature direction. Curvature change rate fields are constructed along the direction of maximum principal curvature and the direction of minimum principal curvature, respectively. Morphological gradient operations are performed on the curvature change rate fields to identify extreme points of curvature change rate. The parameter coordinates of the extreme points of curvature change rate in the parameter coordinate region are used as the insertion positions. At the insertion position, the covariant derivatives of the curvature tensor are calculated along the direction of maximum principal curvature and the direction of minimum principal curvature, respectively. A local curvature compensation vector is constructed based on the covariant derivatives. The parameter coordinates of the insertion position are offset along the local curvature compensation vector to obtain the parameter coordinates of the additional control point. The spatial coordinates of the additional control points are obtained by performing surface mapping calculation on the parameter coordinates of the additional control points, and the surface patch control points are obtained. The spatial coordinates of the additional control points and the surface patch control points are merged in the order of parameter coordinates, and the merged control points are topologically reorganized to obtain the merged control points.

6. The method according to claim 1, characterized in that, Based on the merged control points and the parameter coordinates, a parameterized surface patch is constructed and the surface sampling points are determined. The normal projection distance from the surface sampling points to the surface of the triangular patch is calculated to obtain the normal deviation and the deviation spectrum is determined, including: Obtain the spatial coordinates of the merged control points and the parameter boundary of the parameter coordinate region. Construct a T-spline basis function based on the spatial coordinates of the merged control points. Insert node vectors into the parameter coordinate region and generate a non-uniform parameter grid using a local refinement algorithm. Calculate the spatial coordinates of the surface sampling points based on the nodes in the non-uniform parameter grid and the T-spline basis function. An octree spatial index is constructed based on the spatial coordinates of the surface sampling point, and the nearest triangular facet is queried. The normal vector of the tangent plane of the parameterized surface facet at the surface sampling point is calculated. The surface sampling point is projected onto the surface of the nearest triangular facet along the direction of the tangent plane normal vector to obtain the projection point. The directed distance between the surface sampling point and the projection point is calculated as the normal projection distance. The parameter coordinates corresponding to the nodes in the non-uniform parameter grid are used as the sampling positions, and the normal projection distance is used as the deviation amplitude at the sampling positions. A non-uniform sampling deviation field is constructed based on the sampling positions and the deviation amplitude. The non-uniform sampling deviation field is subjected to a non-uniform fast Fourier transform to obtain a frequency domain representation, and the power spectral density of each frequency component in the frequency domain representation is extracted to obtain the deviation spectrum.

7. The method according to claim 1, characterized in that, Iteratively updating the merged control points based on the amplitude distribution of the deviation spectrum until the normal deviation converges, obtaining converged control points, and generating a reconstructed surface by combining the parameter coordinates includes: The amplitude of each frequency component in the deviation spectrum is extracted and frequency-weighted to obtain the frequency domain deviation energy. A control point sensitivity matrix is ​​constructed based on the frequency domain deviation energy. The dominant sensitivity direction is obtained by performing singular value decomposition on the control point sensitivity matrix. The position adjustment amount of the merged control point is calculated along the dominant sensitivity direction. The spatial coordinates of the merged control point are updated along the position adjustment amount to obtain the optimized merged control point. Based on the optimized merged control points and the parameter coordinates, the parameterized surface patch is reconstructed and the updated surface sampling points are determined. The normal projection distance from the updated surface sampling points to the surface of the triangular patch is calculated to obtain the optimized normal deviation. Calculate the rate of change of the deviation between the optimized normal deviation and the normal deviation before the update, and determine whether the rate of change of the deviation is less than a preset convergence threshold. If it is less, the optimized merged control point is taken as the convergence control point. If it is not less, the optimized merged control point is taken as the new merged control point, and the process of constructing the control point sensitivity matrix and calculating the normal deviation is repeated until the rate of change of the deviation is less than the convergence threshold. Based on the convergence control point and the parameter coordinates, the parameterized surface patch is reconstructed to obtain the reconstructed surface.

8. A high-order continuity preservation and error control system for STL surface reconstruction, used to implement the method of any one of claims 1-7, characterized in that, include: The information extraction unit is used to acquire discretized geometric data and identify the shared edges and vertices of triangular facets to obtain connection information. Based on the connection information, the vertex coordinates of the common boundary on adjacent triangular facets are determined and projected onto the tangent plane of the corresponding triangular facets to obtain parameter coordinates. The normal field construction unit is used to extract the vertex normal vectors on both sides of the common boundary, perform spherical harmonic decomposition to obtain spherical harmonic coefficients, construct a continuous normal vector field based on the spherical harmonic coefficients and the parameter coordinates, and perform path integration to obtain the surface patch control points and curvature tensor components. The high-frequency identification unit is used to perform frequency domain transformation on the vertex coordinates and the vertex normal vector to obtain frequency domain coefficients and determine high-frequency coefficients by combining them with a preset coefficient threshold. The high-frequency coefficients are then inversely transformed to obtain high-frequency vertex coordinates and determine high-frequency triangular facets. The control point solving unit is used to calculate the principal curvature direction of the parameter coordinate region corresponding to the high-frequency triangular patch based on the curvature tensor components, determine the insertion position and additional control points based on the principal curvature direction, and merge them with the control points of the surface patch to obtain the merged control points. The surface reconstruction unit is used to construct a parameterized surface patch based on the merged control points and the parameter coordinates, determine the surface sampling points, calculate the normal projection distance from the surface sampling points to the surface of the triangular patch to obtain the normal deviation and determine the deviation spectrum, iteratively update the merged control points according to the amplitude distribution of the deviation spectrum until the normal deviation converges, obtain the converged control points, and generate the reconstructed surface by combining the parameter coordinates.

9. An electronic device, characterized in that, include: processor; Memory used to store processor-executable instructions; The processor is configured to invoke instructions stored in the memory to execute the method according to any one of claims 1 to 7.

10. A computer-readable storage medium having computer program instructions stored thereon, characterized in that, When the computer program instructions are executed by the processor, they implement the method described in any one of claims 1 to 7.

Images