Drawing method for making three-dimensional model based on three-dimensional laser point cloud

Abstract

The invention relates to the technical field of three-dimensional models, in particular to a drawing method for making a three-dimensional model based on three-dimensional laser point clouds, which comprises the following steps of: sequentially reading point cloud data files of the three-dimensional laser point clouds, inserting coordinate values of each point into an octree structure, and when octree node buffer areas in a memory are full, drawing a three-dimensional laser point cloud into the octree structure; and writing the node data in the buffer area into a hard disk, and traversing the hard disk octree node set from bottom to top. According to the method, the three-dimensional laser point cloud data is inserted step by step, and the hierarchical aggregated octree structure is utilized, so that the point cloud data processing efficiency is effectively improved, and the memory occupation pressure is reduced; representative points in an octree structure are fused with an average normal to construct a continuous function field, a grid model is adaptively generated in a multi-detail-level mode, and the detail retaining capacity and drawing precision of the model are improved. Furthermore, accurate identification and positioning of topological features are realized by adopting calculation of a cell complex sequence and a topological feature noise rank.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of three-dimensional model, in particular to a drawing method for making three-dimensional model based on three-dimensional laser point cloud. BACKGROUND

[0002] The technical field of three-dimensional model refers to using computer graphics, computer vision, geometric modeling, spatial analysis and other related technologies to create, describe, store and display the three-dimensional geometric structure and attribute characteristics of objects or scenes in the real world or virtual world in a digital form through mathematical expression or data representation, thereby providing basic data and visualization support for various industrial design, architectural design, urban planning, film and television animation production, virtual reality, augmented reality, robot navigation, autonomous driving, cultural heritage protection, medical image analysis and other application fields.

[0003] The prior art often uses direct point cloud data processing and mesh generation method when constructing three-dimensional model for complex or large-scale scenes. This method will cause a sharp increase in memory occupation, reduce data loading and processing efficiency, and frequently cause memory overflow or processing efficiency bottleneck in actual application, which is difficult to meet the application requirements of large-scale data scenes. Therefore, improvement is needed. SUMMARY

[0004] The purpose of the present application is to solve the shortcomings in the prior art and provide a drawing method for making three-dimensional model based on three-dimensional laser point cloud.

[0005] In order to achieve the above purpose, the present application adopts the following technical scheme, a drawing method for making three-dimensional model based on three-dimensional laser point cloud, comprising the following steps: Sequentially read the point cloud data file of three-dimensional laser point cloud, insert the coordinate value of each point into the octree structure, when the octree node buffer in the memory is full, write the node data in the buffer to the hard disk, traverse the hard disk octree node set from bottom to top, and establish a hierarchical aggregated attribute octree; Based on the hierarchical aggregated attribute octree, select the coarse level and call the representative point and average normal in the node, construct the function field covering the three-dimensional model, obtain the basic implicit function field, perform isosurface extraction on the basic implicit function field, and generate an adaptive multi-detail level mesh model according to the viewpoint distance and memory limit; Based on the adaptive multi-detail level mesh model, construct a sequence of cell complex, obtain cell complex basic data, traverse each point in the cell complex basic data, calculate the topological noise rank, and obtain a list of topological feature noise ranks; Based on the topology feature noise rank list, each topology noise rank in the topology feature noise rank list is compared with a preset rank value threshold, a classified topology feature set is obtained, and a topology optimized three-dimensional model is established according to the classified topology feature set.

[0006] Preferably, the obtaining step of the hierarchical aggregated attribute octree is: The point cloud data file of the three-dimensional laser point cloud is sequentially read, the three-dimensional coordinate values of each point in the file are parsed one by one, and the parsed coordinate values of each point are sequentially inserted into the pre-created octree structure. When the insertion process causes the storage amount of the octree node buffer in the memory to reach the memory limit threshold, the buffer writing mechanism is triggered, and all octree node data in the current buffer is written to the hard disk storage space as a whole to form a hard disk octree node set. Based on the hard disk octree node set, the node data stored in the hard disk is traversed layer by layer in a bottom-up order according to the octree structure, the spatial geometric data information contained by all child nodes under each parent node is extracted one by one, the coordinate values and normal vector information of the child nodes are merged, and the representative point coordinate values, average normal values and error measurement values of each parent node are calculated according to the merged geometric data information to generate a parent node geometric information aggregation set. Based on the parent node geometric information aggregation set, the representative point coordinate values, average normal values and error measurement values of all parent nodes in the parent node geometric information aggregation set are re-inserted into the corresponding positions of the parent nodes of the octree structure layer by layer, and the calculation is repeated upward layer by layer until the hierarchical aggregation of the entire octree structure is completed to form a hierarchical aggregated attribute octree.

[0007] Preferably, the obtaining step of the basis implicit function field is: Based on the hierarchical aggregated attribute octree, all nodes of the rough level are located, the representative point coordinate vectors and average normal vectors stored in the nodes are extracted one by one, all representative point coordinate vectors and average normal vectors are traversed, and each query position in the three-dimensional space is mapped to establish an initial function field covering the three-dimensional model to obtain initial data of the function field. According to the initial data of the function field, the basis implicit function field is calculated.

[0008] Preferably, the obtaining step of the adaptive multi-detail level grid model is: Based on the basis implicit function field, isosurface extraction is performed on each query position in the three-dimensional space, and the deeper level node data of the hierarchical aggregated attribute octree is called to replace the implicit function values of the corresponding region under the condition of the viewpoint distance and the memory limit, the local function field is updated and re-subdivided block by block to generate an adaptive multi-detail level grid model.

[0009] Preferably, the step of obtaining the cavity complex basis data is: Based on the adaptive multi-detail level grid model, the geodesic distance of all vertices to the reference point is calculated, the vertices and connected edges and faces are introduced in order of geodesic distance value from small to large, and the complete cavity complex sequence containing vertices, edges and faces is formed by gradually expanding under the order, and the cavity complex basis data is obtained.

[0010] Preferably, the step of obtaining the topological feature noise rank list is: According to the cavity complex basis data, the topological feature noise rank is calculated; Based on the topological feature noise rank, all topological feature points are traversed, the topological feature noise ranks are summarized, and the topological feature noise rank list is generated.

[0011] Preferably, the step of obtaining the classified topological feature set is: Based on the topological feature noise rank list, the noise rank corresponding to each topological feature in the topological feature noise rank list is extracted, the noise rank of each topological feature is compared with the pre-set rank value threshold one by one, if the noise rank of the topological feature is greater than the rank value threshold, it is marked as a noise feature, and if the noise rank of the topological feature is less than or equal to the rank value threshold, it is marked as a significant structure, and a classified topological feature set is generated.

[0012] Preferably, the step of obtaining the topological optimization three-dimensional model is: Based on the classified topological feature set, all topological features marked as noise features are screened, the spatial position of each noise feature in the adaptive multi-detail level grid model is traced one by one, all edge elements and face elements corresponding to each noise feature are located, and a noise feature positioning set is formed; Based on the noise feature positioning set, all edge elements and face elements in the noise feature positioning set are deleted one by one in the adaptive multi-detail level grid model, and the grid elements are reconnected and filled and repaired through the local area boundary after topological cutting, and a topological optimization three-dimensional model is established.

[0013] Compared with the prior art, the application has the advantages and positive effects that: The application effectively improves the point cloud data processing efficiency and reduces the memory occupation pressure by inserting three-dimensional laser point cloud data step by step and using the hierarchical aggregation octree structure; the representative points and average normals in the octree structure are fused to construct a continuous function field, and the grid model is adaptively generated in a multi-detail level manner, which improves the detail retention capability and rendering accuracy of the model; further, the accurate identification and positioning of the topological features are realized by using the calculation of the cell complex sequence and the topological feature noise rank, the noise data and redundant topological structure in the model are removed, and the authenticity and structural stability of the three-dimensional model in the local area are improved; by using the fine topological noise positioning and local reconstruction manner, the local topological defect area is accurately cut and refilled, so that the topological structure of the generated three-dimensional model in the complex area is more reasonable, the overall model quality and visual authenticity are improved, and the reliability and practicality of the model in subsequent application in the fields of precision industrial design, cultural heritage protection, building engineering planning and the like are enhanced. BRIEF DESCRIPTION OF DRAWINGS

[0014] Figure 1 The figure is a schematic diagram of the steps of the application. DETAILED DESCRIPTION

[0015] In order to make the purpose, technical scheme and advantages of the application more clear, the application is further described in detail below in combination with the drawings and examples. It should be understood that the specific examples described herein are only used to explain the application, and are not used to limit the application.

[0016] Please refer to Figure 1 The application provides a technical scheme, a drawing method for making a three-dimensional model based on three-dimensional laser point cloud, comprising the following steps: The point cloud data file of the three-dimensional laser point cloud is sequentially read, the coordinate value of each point is inserted into the octree structure, when the octree node buffer in the memory is full, the node data in the buffer is written to the hard disk, the hard disk octree node set is traversed from bottom to top, and the hierarchical aggregation attribute octree is established; Based on the hierarchical aggregation attribute octree, the representative points and average normals in the nodes are selected and called, the function field covering the three-dimensional model is constructed, the basic implicit function field is obtained, the isosurface of the basic implicit function field is extracted, and the adaptive multi-detail level grid model is generated according to the viewpoint distance and memory limit; Based on the adaptive multi-detail level grid model, the cell complex sequence is constructed, the cell complex basic data is obtained, each point in the cell complex basic data is traversed, the topological noise rank is calculated, and the topological feature noise rank list is obtained; Based on the topology feature noise rank list, each topology noise rank in the topology feature noise rank list is compared with a preset rank value threshold, a classified topology feature set is obtained, and a topology optimization three-dimensional model is established according to the classified topology feature set.

[0017] The obtaining step of the hierarchical aggregation attribute octree is: Sequentially read the point cloud data file of the three-dimensional laser point cloud, parse the three-dimensional coordinate values of each point in the file one by one, sequentially insert the coordinate values of each parsed point into the pre-created octree structure, when the insertion process causes the memory of the octree node buffer in the memory to reach the memory limit threshold, trigger the buffer write mechanism, write all the octree node data in the current buffer to the hard disk storage space as a whole, and form a hard disk octree node set; Based on the hard disk octree node set, the node data stored in the hard disk is traversed layer by layer from bottom to top according to the octree structure, the spatial geometric data information of all child nodes under each parent node is extracted one by one, the coordinate values and normal vector information of the child nodes are merged, and the representative point coordinate values, average normal values and error measurement values of each parent node are calculated according to the merged geometric data information, and a parent node geometric information aggregation set is generated; Based on the parent node geometric information aggregation set, the representative point coordinate values, average normal values and error measurement values of all parent nodes in the parent node geometric information aggregation set are re-inserted into the corresponding positions of the parent nodes of the octree structure layer by layer, and the calculation is repeated layer by layer upwards until the hierarchical aggregation of the entire octree structure is completed, forming a hierarchical aggregation attribute octree.

[0018] Specifically, the point cloud data file of the three-dimensional laser point cloud is read in sequence, specifically, the.xyz format file separated by spaces or commas is scanned line by line, each line of text is parsed, three floating point numbers representing spatial points are extracted, which correspond to X, Y and Z three-dimensional coordinate values respectively, before parsing the first point, the maximum and minimum values of the three-dimensional coordinates of all points are determined by pre-defining a space range or performing a quick pre-scan on the entire point cloud file, so as to establish an initial cubic space that can completely enclose all point cloud data, and based on this cubic space, an octree structure containing only a root node is created in the memory, then the coordinate values of each point parsed, such as point (10.5, 20.3, 5.8), are inserted into the root node of the octree according to their coordinate positions, and recursively inserted down to a leaf node, if there is a data point in the leaf node, the leaf node is split into eight new child nodes, and the data points in the original node and the current data point to be inserted are redistributed to the corresponding new child nodes, the insertion process continues, while the system monitors the special buffer area in the memory for storing octree node data in real time, the size of the buffer area is preset as a memory limit threshold, the setting of the threshold is based on the available physical memory of the current running environment, the specific setting method is to obtain the total physical memory of the current system, subtract the reserved memory occupied by the operating system and other key processes, and take 80% of the remaining available memory as the threshold, for example, if the available memory of the system is 10GB, the threshold is set to 8GB, when the node buffer area used for storage exceeds 8GB for the first time due to the insertion operation, the buffer area write mechanism is triggered, which will pause the new point cloud insertion operation, serialize all node data in the octree node buffer area in the current memory, including the level, space range, parent-child relationship index and point coordinate data contained in the leaf node, and write it as a whole data block to the specified storage area on the hard disk at one time, empty the memory buffer area and continue to process the subsequent point cloud data, repeat the process until all point cloud data files are read and processed, forming one or more ordered node data files on the hard disk, which together constitute a complete hard disk octree node set.

[0019] Based on the hard disk octree node set, a bottom-up aggregation process is started, first locating to the deepest level of the octree structure, that is, the level where the leaf nodes containing the original point cloud data are located, and traversing the node data stored in the hard disk layer by layer in the order of the bottom-up of the octree structure. The specific implementation is to load all the node data of the current processing level into the memory at one time. For a parent node, the spatial geometric data information contained in all eight child nodes under it is extracted one by one. If the child node is a leaf node, the original point cloud coordinate values stored in it are directly extracted. If the normal vector information is not contained in the original point cloud data, the initial normal vector of each point set in the leaf node needs to be calculated first. The calculation method is as follows: for each point, all the points in the leaf node where it is located are selected as the neighborhood, and the covariance matrix of the neighborhood point set is calculated through principal component analysis (PCA). The eigenvector corresponding to the smallest eigenvalue is taken as the initial normal vector of the point. Then, the coordinate values and the calculated normal vector information of all the points contained in all the child nodes under the parent node are merged, and the three core aggregation attributes of the parent node are calculated according to the merged geometric data information. First, the representative point coordinate value is calculated, which is the arithmetic average of the coordinate values of all points contained in the child nodes, that is, the X, Y and Z coordinates of all points are summed and then divided by the total number of points. Second, the average normal value is calculated, which is the vector sum of the normal vectors of all points contained in the child nodes, and then the result vector is normalized to obtain a unit vector. Finally, the error measurement value is calculated, which is used to quantify the deviation of all points in the current node from a fitting plane defined by the representative point and the average normal. The calculation method is as follows: for each point in the child node set, the square of the perpendicular distance from the point to the fitting plane is calculated, then all these squared distance values are accumulated and averaged. The representative point coordinate value, average normal value and error measurement value of each parent node calculated are packaged into an aggregation record. When all the parent nodes of a level are calculated, these aggregation records are collected to generate a parent node geometric information aggregation set.

[0020] Based on the parent node geometry information aggregation set, the representative point coordinate value, the average normal value and the error measurement value of all parent nodes calculated in the previous stage recorded in the aggregation set are updated layer by layer into the hard disk octree structure, specifically, these aggregation attribute data are written into the storage location of the corresponding parent node, replacing the original placeholder information or low-level data of the node, after completing the update of a level, the program automatically moves up a level, taking the parent node layer just completed aggregation as the new "child node layer", and starts a new round of aggregation calculation, in the aggregation process of this higher level, the calculation method is adjusted, when calculating the representative point coordinate value of a higher level parent node (for example, L-2 level), the original point cloud is no longer used, but the representative point coordinates of its direct child node (L-1 level) are used, and a weighted average method is used for calculation, wherein the weight of each child node representative point is the number of original point clouds it represents, similarly, when calculating the average normal value, the same weight is used to perform weighted vector summation on the average normal of the child node, and then normalization is performed, for the upward propagation of the error measurement value, the total error measurement value of a parent node is composed of two parts, one part is the fitting plane error calculated according to the representative points of all child nodes of the parent node, the calculation method is the same as the foregoing, the other part is the maximum value of the error measurement of all child nodes, and the final error measurement value of the parent node is set as the sum of the two parts, the cycle process of repeated calculation and update layer by layer upwardly continues, starting from the deepest non-leaf node level, and processing upwardly in turn, until the aggregation attributes of the direct child node of the root node are calculated and updated, and finally the aggregation attribute calculation of the root node itself is completed, thus, each non-leaf node in the entire octree structure stores the condensed representation of the geometric features in the space range it covers, that is, the representative point, the average normal and the error measurement, thereby forming a hierarchical aggregation attribute octree which is complete in structure, rich in attributes and supports multi-level detail query.

[0021] The acquisition step of the basic implicit function field is: Based on the hierarchical aggregation attribute octree, all nodes of the rough level are located, and the representative point coordinate vector and the average normal vector stored in the nodes are extracted one by one, all representative point coordinate vectors and average normal vectors are traversed, and are mapped to each query position in the three-dimensional space, an initial function field covering the three-dimensional model is established, and initial data of the function field is obtained; According to the initial data of the function field, the basic implicit function field is calculated, and the calculation formula is: ; The calculation formula of the middle weight function is: ; Wherein, is the basic implicit function value of the query point , is a coordinate vector of a position to be calculated in three-dimensional space, is a coordinate vector of a representative point, is a coordinate vector of a representative point, is an average normal vector of a node, is an error metric value of a node, is the largest error metric value among all nodes, is an error influence factor, is a smoothing radius parameter, is the number of adjacent representative points involved in the calculation. In particular, based on the hierarchical aggregation attribute octree, a suitable rough level needs to be determined first. The selection criteria of the rough level is to minimize the number of nodes involved in the calculation under the premise of ensuring the integrity of the initial model profile. The specific positioning method is to start from the next layer of the root node of the hierarchical aggregation attribute octree (i.e. level 1) and traverse layer by layer downward, count the total number of nodes at each level, and set a node number threshold, which is dynamically set according to the memory capacity of the processing hardware and the target processing efficiency. For example, in a workstation with 32GB of physical memory, to ensure the smoothness of the initial calculation, the threshold can be set to 30000. The system will select the level whose total number of nodes first exceeds 30000 as the rough level. After selecting the rough level, the system will traverse all nodes in this level and extract the core aggregation attributes from the storage structure of each node, i.e. the calculated representative point coordinate vector and the average normal vector. Then, to establish the initial function field covering the three-dimensional model, a three-dimensional grid needs to be defined as the query space within the bounding box space of the entire model. The size of the grid is determined according to the approximate size of the model and the required initial resolution. For example, for a scanning scene with a size of 10 meters by 10 meters by 5 meters, a 128 by 128 by 64 query grid can be set. Each center point in the grid is a query position. Subsequently, each query position in the three-dimensional grid is traversed. For each query position, a set of representative points that have an impact on it needs to be determined. The implementation is to use the spatial indexing capability of the hierarchical aggregation attribute octree to quickly retrieve all representative points in the rough level that fall within a preset influence radius range with the current query position as the center. The setting of the influence radius is related to the average size of the octree nodes, and is usually set to 2 to 3 times the average diagonal length of the nodes in the rough level. The coordinate vector, average normal vector and error metric value of each representative point retrieved are extracted together to form a data set associated with the current query position. All query positions and their corresponding associated data sets are summarized to finally obtain the initial data of the function field.

[0022] In particular, based on the hierarchical aggregation attribute octree, a suitable rough level needs to be determined first. The selection criteria of the rough level is to minimize the number of nodes involved in the calculation under the premise of ensuring the integrity of the initial model profile. The specific positioning method is to start from the next layer of the root node of the hierarchical aggregation attribute octree (i.e. level 1) and traverse layer by layer downward, count the total number of nodes at each level, and set a node number threshold, which is dynamically set according to the memory capacity of the processing hardware and the target processing efficiency. For example, in a workstation with 32GB of physical memory, to ensure the smoothness of the initial calculation, the threshold can be set to 30000. The system will select the level whose total number of nodes first exceeds 30000 as the rough level. After selecting the rough level, the system will traverse all nodes in this level and extract the core aggregation attributes from the storage structure of each node, i.e. the calculated representative point coordinate vector and the average normal vector. Then, to establish the initial function field covering the three-dimensional model, a three-dimensional grid needs to be defined as the query space within the bounding box space of the entire model. The size of the grid is determined according to the approximate size of the model and the required initial resolution. For example, for a scanning scene with a size of 10 meters by 10 meters by 5 meters, a 128 by 128 by 64 query grid can be set. Each center point in the grid is a query position. Subsequently, each query position in the three-dimensional grid is traversed. For each query position, a set of representative points that have an impact on it needs to be determined. The implementation is to use the spatial indexing capability of the hierarchical aggregation attribute octree to quickly retrieve all representative points in the rough level that fall within a preset influence radius range with the current query position as the center. The setting of the influence radius is related to the average size of the octree nodes, and is usually set to 2 to 3 times the average diagonal length of the nodes in the rough level. The coordinate vector, average normal vector and error metric value of each representative point retrieved are extracted together to form a data set associated with the current query position. All query positions and their corresponding associated data sets are summarized to finally obtain the initial data of the function field.

[0023] ​formula: Medium weight function The calculation formula is: ; The benefit of the formula is that it constructs a smooth scalar field whose zero-value surface can accurately approximate the surface of the real object represented by the original point cloud, and the weight function The error adjustment term is introduced in , which utilizes the error metrics pre-calculated at each node in the hierarchical aggregated attribute octree , which makes it possible to use the point cloud data in areas with poor quality, high noise or complex details (corresponding to The weight of its representative points in the function field calculation is actively reduced, which effectively suppresses the interference of noise on surface reconstruction and improves the robustness and smoothness of the model. The weight of the representative points is relatively higher when the value is low, ensuring the accurate expression of the model features. This adaptive weighting mechanism enables the model reconstruction process to intelligently distinguish between signals and noise, ultimately generating high-quality 3D models with clear details and smooth surfaces. is the coordinate vector of the position to be calculated in three-dimensional space. This parameter represents any point in the three-dimensional space for function value query. In specific implementation, these query points are usually the grid centers of a predefined three-dimensional grid. The three-dimensional grid completely covers the entire model space defined by the root node of the hierarchical aggregation attribute octree. The resolution of the grid determines the basic accuracy of the final generated model. After obtaining the initial data of the function field, That is, the coordinates of each query location recorded in the data. For example, in the previous step, a 128 by 128 by 64 query grid was created, then The center coordinates of these 128 times 128 times 64 grid points will be taken in turn to calculate the implicit function value of each position in the entire space, forming a complete scalar field. For example, the coordinate vector of a specific query point is , in meters; The number of adjacent representative points involved in the calculation. This parameter defines the number of adjacent representative points involved in the calculation of any query point. When calculating the function value at , it is necessary to consider the number of representative points around it. The value of this parameter directly affects the locality, smoothness and efficiency of the calculation. The setting method is to first randomly select 10% of the sample points from all the representative points of the coarse level, and for each sample point, use the octree structure to quickly query a certain radius around it (for example, three times the smoothing radius). ) and record the number of neighbors within it, then calculate the mean and standard deviation of all these values, and finally The value is set to the mean plus one standard deviation and rounded up. This setting makes It can adapt to the local density of representative points. For example, after statistical calculation, the average value of the number of neighbors is 18.7 and the standard deviation is 4.2. The setting value is , rounded up to get ; For the The coordinate vectors of the representative points, For the The average normal vector, For the The error metric value of each node, these three parameters are directly derived from the initial data of the function field generated in the previous step, for each query point , the function field initial data has stored the associated The nearest neighbor represents the information of the point, so when calculating When the program reads this data structure directly Group Triple data, which are calculated and stored in the aggregation stage of constructing the hierarchical aggregation attribute octree, reflect the local geometric features of the point cloud at different spatial scales. For example, for the query point , the first representative point information read from the initial data of the function field is: coordinate vector , the average normal vector , error metric value ; The maximum error metric value among all nodes. This parameter is used as a normalization factor to scale the error metric values ​​of all nodes to a controllable range. The way to obtain it is to traverse all nodes of the coarse level after selecting the level and read the error metric value stored in each node. , and record the maximum value. This process is performed once before the first calculation of the function field, and all subsequent calculations of the coarse level use the same This can ensure the global consistency of the error impact. For example, the system traverses all 32768 nodes in the coarse level and finds that the error metric value range is between 0.001 and 0.45. The value of is set to 0.45; is the error impact factor, which is a dimensionless hyperparameter used to adjust the error metric The numerical value of the weight function is empirically set and needs to be adjusted according to the noise level of the point cloud data and the smoothness requirements of the target model. The setting process is to select a small area with representative features (such as planes, curved surfaces, and sharp corners) and use different Generate local models with different values ​​(e.g. 1.0, 5.0, 10.0, 20.0) and select the one that best preserves geometric details while effectively suppressing noise by visual inspection or comparison with a high-precision reference model. For example, for a building point cloud from a terrestrial laser scanner, the noise level is low, but sharp features such as wall corners need to be retained. After testing, it was found that The effect is best when ; is the smoothing radius parameter, which defines the weight function The spatial influence range of determines the smoothness of the surface reconstruction. Its setting is closely related to the average spacing of the representative points in the coarse level. The specific setting steps are as follows: randomly select 1000 representative points in the coarse level, and for each representative point, use the octree to calculate the distance from it to its fifth nearest neighbor representative point, and then calculate the average of these 1000 distance values, which is recorded as , and finally the smoothing radius Set to 2.5 times the average distance, that is, , so that the influence radius can cover enough neighboring points to ensure the continuity and smoothness of the function field. For example, the calculation results are is 0.8 meters, then the smoothing radius Set to rice.

[0024] Calculation process: With a query point As an example, to simplify the demonstration, only the three nearest representative points (i.e. , used in actual calculations ) is calculated, and the values ​​of each parameter are as follows: .

[0025] The information of the three representative points is: The first step is to calculate the weight of each representative point : for : Distance squared: Weight: For : Distance square: Weight: For : Distance square: Weight: Second step, calculate the weighted directed distance sum in the molecule: Item 1: Item 2: Item 3: Molecule and: ; Third step, calculate the weight sum in the denominator: Denominator sum: ; Fourth step, calculate the final function value : ; The results show that the implicit function value at the query point is 0.0734, since the value is greater than 0, it indicates that the point is located in the "outside" of the implicit surface defined according to the local representative point, and the numerical value 0.0734 approximately represents the shortest distance from the point to the implicit surface. By repeating this calculation for all grid query points, a complete three-dimensional scalar field, i.e. the basic implicit function field, can be obtained, and the subsequent isosurface extraction algorithm will find the positions of all to construct the final mesh model.

[0026] The steps for obtaining an adaptive multi-detail level mesh model are: Based on the basic implicit function field, perform isosurface extraction for each query position in the three-dimensional space, and combine the viewpoint distance condition and memory limit to replace the implicit function value of the corresponding region with the data of the deeper level node of the hierarchical aggregation attribute octree, update and re-subdivide the local function field block by block, and generate an adaptive multi-detail level mesh model.

[0027] Specifically, based on the basic implicit function field, the marching cubes algorithm is first used to extract the initial isosurface for each query position in the three-dimensional space. Specifically, each cube unit (Voxel) in the three-dimensional grid composed of the basic implicit function field is traversed, and the implicit function values ​​of its eight vertices are read. These values ​​are compared with the isosurface threshold (fixed to 0). According to the positive and negative sign combination of the eight vertex function values, the topological configuration of the triangle face to be generated in the cube is determined from a pre-calculated lookup table containing 256 cases. The position of the triangle face vertex is obtained by linear interpolation of the two endpoints of the edge crossing the isosurface, completing the topology of all After traversing the cube units, the basic mesh covering the entire model is obtained. Then, the system enters a dynamic, adaptive refinement loop based on viewpoint distance and memory limit. The loop monitors the position of the virtual camera, i.e., the viewpoint, in real time and adjusts the level of detail according to the distance between the viewpoint and each part of the model. Specifically, two distance thresholds are set, the near threshold and the far threshold. The setting of these two thresholds is based on the diagonal length of the model bounding box. For example, if the diagonal length of the model bounding box is 100 meters, the near threshold can be set to 15 meters and the far threshold can be set to 60 meters. For any point on the model, if its If the distance from the viewpoint is less than 15 meters, it is judged as a high-detail area. If the distance is between 15 meters and 60 meters, it is a medium-detail area. If the distance is greater than 60 meters, it is a low-detail area. At the same time, the system sets a memory limit, which is specifically reflected in a target total number of triangles, for example, 3 million triangles. During the refinement process, the program will count the total number of faces of the current model in real time. When the total number of faces approaches this limit, it will prioritize the refinement of the near area and may stop or simplify the mesh of the distant area. For the local mesh judged as a high-detail area, the system will locate its corresponding hierarchical aggregation attribute. The node of the octree is retrieved, and the child node data of one or two deeper levels below the node is called. These denser and more accurate representative points and average normal data are used to recalculate the basic implicit function field only in the corresponding high-detail area, thereby updating the implicit function value of the local area. Finally, the marching cube algorithm is called again to re-mesh the updated small-scale function field, and replace the original coarse mesh with the generated high-precision mesh. This series of "judgment-recall-update-mesh" operations will continue as the viewpoint moves, thereby dynamically generating an adaptive multi-level detail mesh model.

[0028] The steps to obtain the basic data of the cell complex are: Based on the adaptive multi-level detail mesh model, the geodesic distances from all vertices to the reference point are calculated, and the vertices and connected edges and faces are introduced in ascending order of geodesic distance values. They are gradually expanded in this order to form a cellular complex, and a complete cellular complex sequence including vertices, edges and faces is constructed to obtain the basic data of the cellular complex.

[0029] Specifically, based on the adaptive multi-detail level grid model, a unique reference point needs to be determined first, the selection of the reference point has a benchmark significance for subsequent calculation, the selection method is to traverse all vertices of the model, find out the vertex with the minimum Z-axis coordinate value in the preset coordinate system, and designate it as the reference point, then, taking the reference point as the starting point, the geodesic distance of all other vertices on the model to the reference point is calculated, that is, the shortest path length along the model surface, the calculation is completed by executing Dijkstra algorithm on the adaptive multi-detail level grid model, the specific process is that the vertices of the grid are regarded as nodes of a graph, the edges of the grid are regarded as edges of the graph, the weight of the edge is set as the Euclidean distance between the two endpoints, the algorithm starts from the reference point, sets the distance of the reference point as 0, and initializes the distance of all other vertices as infinity, then a priority queue is used to iteratively access and update the distance of all vertices until the shortest distance of all vertices to the reference point is calculated, after obtaining the geodesic distance value of each vertex, all components of the model, including vertices, edges and faces, are introduced in a specific order to build a filtration, this order is completely determined by the geodesic distance value, first, all vertices are arranged in ascending order according to their geodesic distance values to form an ordered vertex list, then the cell complex is gradually expanded according to this vertex order, the specific rule is that when all component edges and faces of a vertex have been introduced according to the geodesic distance of the other vertices, the edge or face is added to the sequence, a more explicit introduction time definition is that the introduction time of an edge is the larger geodesic distance value of its two endpoints, and the introduction time of a face is the largest geodesic distance value of its three vertices, finally, all vertices, edges and faces are uniformly sorted according to their introduction time (i.e. geodesic distance value), when the introduction time is the same, the sorting is performed according to the principle of dimension priority (vertices first, edges first, and faces first), thereby a complete cell complex sequence containing all vertices, edges and faces of the model and strictly sorted according to the geodesic distance value is constructed, this sequence is the cell complex basic data.

[0030] The obtaining step of the topological feature noise rank list is: According to the cell complex basic data, the topological feature noise rank is calculated, and the calculation formula is: ; Wherein, is the dimensionless noise rank of the i-th topological feature, is the birth time value of the i-th topological feature, is the death time value of the i-th topological feature, is the birth time value of the i-th topological feature, is the death time value of the i-th topological feature, is the birth time value of the i-th topological feature, is the death time value of the i-th topological feature, The life cycle of a topological feature, To prevent extremely small positive numbers with zero denominators, is the sensitivity factor of the dimensionless ratio; Based on the topological feature noise rank, all topological feature points are traversed, the topological feature noise ranks are summarized, and a topological feature noise rank list is generated.

[0031] Specifically, the formula: , the formula is beneficial in that the topological characteristic noise rank The calculation formula can quantify and distinguish the real topological structure in the 3D model from the false topological features introduced by point cloud noise or reconstruction algorithm, and the birth time of the topological features is Its life cycle The ratio of is used as the core metric to capture the "prematurity" of a feature. Usually, small bumps or holes caused by noise will appear in the early stage of the filtering process (i.e., small geodesic distance values), i.e. It is very small and will soon be filled by the units introduced subsequently and disappear, that is, the life cycle is very short, which results in its noise rank Very large, on the contrary, the real macrostructure on the model, such as a real hole or handle, is born at Relatively large, and life cycle Very long, calculated noise rank will be very small, through the sensitivity factor This difference can be further amplified so that noise and significant structures are clearly separated in terms of rank value. This method does not rely on geometric size, but is judged purely from the perspective of topological persistence. Therefore, it is highly robust to various forms of noise. For the The birth time value of the topological feature, For the The death time value of the topological feature. These two parameters are directly derived from the result of persistent homology calculation of the cell complex basic data obtained in the previous step. The persistent homology algorithm will track the generation and disappearance of topological features (such as connected branches, rings, cavities) of the cell complex sequence during the filtering process. Topological features, whose birth time Corresponds to the filtering value that first causes the feature to appear after the introduction of a unit (vertex, edge or face) in the cell complex sequence, that is, the geodesic distance value, and its death time This corresponds to the filter value that causes the feature to be filled or merged with other features and disappear after the introduction of another unit. These pairs Values constitute the persistent diagram, which is the standard output of persistent homology analysis, for example, when processing the cell complex base data, the persistent homology analysis (such as Gudhi, Dionysus) will output a list, each item in the list records a topological feature pair, for example, for a 1-dimensional ring (topological feature) representing a small indentation, it is recorded as , where 5.2 and 5.4 are both geodesic distance values, in meters; To prevent the denominator from being a very small positive number, this is a numerical stability parameter, used to handle a special case where the birth and death of a topological feature occur at exactly the same filtering value , Can be set to ; Sensitivity factor for dimensionless ratio, this is a hyperparameter for adjusting the sensitivity of noise rank to the ratio of birth time and life cycle, its role is to amplify or reduce the discrimination between noise and significant structure, The value depends on the specific application scenario and the statistical properties of the noise, the setting method is usually based on experience or determined by cross-validation, an effective setting process is, first prepare a small test model with known topological noise and significant structure annotation, then calculate the noise rank of all features under a series of candidate values (for example ), and use a preset rank threshold for classification, and finally select the value that makes the F1 score of the classification result (compared with the annotation) highest as the final parameter, for a typical three-dimensional laser point cloud reconstruction model, the noise usually manifests as high-frequency and small-scale features, experience shows that when the value is 2.0, it can well distinguish these noises, therefore, in this example, it is set to .

[0032] Calculation process: Consider two topological feature pairs obtained from persistent homology analysis, one is a suspected noise feature, and the other is a suspected significant structure feature, and bring them into the formula to calculate their noise ranks.

[0033] Parameter setting: .

[0034] Example 1: Suspected noise feature (denoted as the first feature, k = 1); Persistent pair is .

[0035] Birth time value: ; Death time value: ; Lifetime: ; Compute noise rank : ; Example 2: Feature of suspected significant structure (denoted as the 2nd feature, k = 2); Persistent pair .

[0036] Birth time value: ; Death time value: ; Lifetime: ; Compute noise rank : ; The result shows that the noise rank of the first topological feature is 676.0, which is a very large value, which confirms its characteristics of being born early (k = 1) and having a very short lifetime, which is a typical topological noise. The noise rank of the second topological feature is 0.1704, which is a very small value, reflecting its characteristic of having a lifetime much greater than its birth time, which is a typical significant topological structure. Through this formula, each topological feature in the model will get a noise rank, which will be used in the subsequent classification step.

[0037] Based on the noise rank of the topological feature, the system first traverses all the topological features obtained by persistent homology analysis. These features exist in the form of persistent pairs , and the corresponding topological feature noise rank has been calculated for each feature. In the traversal process, the unique identifier of each topological feature (such as its index number k in the persistent homology calculation result) and its calculated topological feature noise rank are stored as a key-value pair. At the same time, in order to facilitate subsequent analysis and tracing, the dimension of the feature (0-dimensional represents connected components, 1-dimensional represents loops or tunnels, and 2-dimensional represents cavities), the birth time value and the death time value are also recorded. All this information (feature ID, noise rank, dimension, birth time, and death time) is integrated together to form a structured record. When all the topological features are traversed, these records are collected into a single data structure, usually a list or array. This data set containing all the topological features and their quantitative noise measurements, i.e., the topological feature noise rank list, is the topological feature noise rank list.

[0038] The steps to obtain the classified topological feature set are: Based on the topological feature noise rank list, the noise rank corresponding to each topological feature in the topological feature noise rank list is extracted item by item, and the noise rank of each topological feature is numerically compared with the pre-set rank value threshold one by one. If the noise rank of the topological feature is greater than the rank value threshold, it is marked as a noise feature. If the noise rank of the topological feature is less than or equal to the rank value threshold, it is marked as a significant structure to generate a classified topological feature set.

[0039] Specifically, based on the topological feature noise rank list, the system begins to classify each topological feature. The core is to compare its noise rank with a pre-set rank value threshold. The setting of the rank value threshold is the key to distinguishing noise from significant structures. The setting method is as follows: First, all the noise ranks in the topological feature noise rank list are compared. Extract it and perform logarithmic transformation on it, that is, calculate Because the distribution of noise ranks is usually highly skewed, taking the logarithm can make its distribution closer to the normal distribution, which is convenient for analysis. Then, the K-means clustering algorithm is applied to all the logarithmic noise rank data, and the number of clusters is set to 2. The K-means algorithm will automatically divide these logarithmic values ​​into two clusters, one corresponding to the larger logarithmic value (high noise rank) and the other corresponding to the smaller logarithmic value (low noise rank). Then, the boundaries of the two clusters are found, that is, the average of the minimum value in the high-value cluster and the maximum value in the low-value cluster is taken, and this average value is used as the threshold in the logarithmic space. Finally, this logarithmic threshold is exponentially transformed (i.e. ), get the rank value threshold in the original noise rank space. For example, after cluster analysis, the threshold in the logarithmic space is 3.5, then the original rank value threshold is set to After setting the rank value threshold, the program extracts the noise rank corresponding to each topological feature in the topological feature noise rank list one by one, and compares the noise rank of each topological feature with the threshold 33.1. If the noise rank of the topological feature is greater than 33.1, the feature is marked as a "noise feature". If its noise rank is less than or equal to 33.1, it is marked as a "significant structure". All features and their corresponding classification labels (noise features or significant structures) are summarized to generate a classified topological feature set.

[0040] The steps to obtain the topology optimization 3D model are as follows: Based on the classified topological feature set, all topological features marked as noise features are screened, and the spatial position of each noise feature in the adaptive multi-level detail mesh model is traced one by one. All edge units and surface units corresponding to each noise feature are located to form a noise feature location set. Based on the noise feature positioning set, all edge elements and surface elements in the noise feature positioning set are deleted one by one in the adaptive multi-detail level grid model, and the grid element reconnection and filling repair are performed through the local region boundary after topological cutting, to establish a topologically optimized three-dimensional model.

[0041] Specifically, based on the classified topological feature set, the system first starts a screening program, which traverses each record in the classified topological feature set and only extracts those topological features that are explicitly marked as "noise features". For each noise feature screened out, the system will perform a backtracking positioning operation, which aims to accurately find the specific grid elements that constitute this noise feature in the adaptive multi-detail level grid model. The basis for backtracking is the associated information recorded during the persistent homology calculation. When calculating the persistent homology, the system not only outputs the birth and death time pairs of each topological feature, but also records the "generator" and "annihilator" of each topological feature, i.e., which edge element or surface element is introduced to cause the birth and death of the feature. For example, for a 1D ring marked as noise (representing a small hole or recess), the system reads the "generating edge" that causes the formation of this ring and the "annihilating surface" that eventually fills the ring from the associated information. These two elements are the core of the noise feature. Then, in order to completely locate the entire spatial region affected by the noise feature, the system takes all the vertices of these two core elements (the generating edge and the annihilating surface) as the starting point and performs a restricted breadth-first search on the connection relationship graph of the adaptive multi-detail level grid model. The restriction condition of the search range is to only include those adjacent elements whose geodesic distance values are less than or equal to the death time value of the noise feature This process will continue to expand until no more adjacent elements that meet the conditions can be found. All edge elements and surface elements visited in this search process are collectively identified as the complete spatial entities that constitute the noise feature. The indices or pointers of these elements are collected to form an independent set. This tracing and positioning process is repeated for all screened noise features, and finally all these independent element sets are aggregated to form the noise feature positioning set.

[0042] Based on the noise feature positioning set, the system starts the grid optimization process, which first marks all the edge elements and face elements in the adaptive multi-detail level grid model data structure for deletion according to all the edge element and face element indexes recorded in the noise feature positioning set, and after marking all the elements to be deleted, the system performs a one-time batch deletion operation to remove all the marked edge elements and face elements from the model. This operation will generate one or more holes at the positions where the topological noise originally existed, and each hole is surrounded by a ring of boundary edges. The next step is to repair the holes generated by the topological cutting. The system traverses each hole, first extracts the ordered vertex ring that constitutes the boundary of the hole by tracking the connection relationship of its boundary edges, and for each unstructured polygon region defined by a vertex ring, the system uses the constrained Delaunay triangulation algorithm to fill it. This algorithm takes the boundary edges of the hole as a hard constraint condition and generates a series of new triangular patches inside the polygon. These newly generated patches completely cover the hole area and their boundaries match the original boundary edges exactly. During the triangulation process, in order to ensure the quality of the newly generated grid, for holes with large area or long and narrow shape, the algorithm is allowed to add new vertices, i.e. Steiner points, inside the hole to avoid generating poor quality long and narrow triangles. After the initial filling, the system performs local smoothing and quality optimization processing on the newly generated grid patches and their adjacent regions. Specifically, it performs 3 to 5 iterations of Laplace smoothing operation on all newly added Steiner points and vertices on the original boundary ring, each iteration moves the vertex to the geometric center of the positions of all its adjacent vertices. Then, the edge flip operation is performed on the newly generated triangular patches and the original triangular patches that share edges with them. Each internal edge is iteratively checked, and if flipping the edge can increase the minimum internal angle of its two adjacent triangles, the flip is performed. The process continues until no further improvement in grid quality can be achieved through edge flipping. After all the holes have been processed through the deletion and repair process, a topologically optimized three-dimensional model is finally obtained, which has removed the small topological noise, has a smoother surface, and has a more reasonable structure.

[0043] The above is only a preferred embodiment of the present application, and does not limit the present application in other forms. Any skilled person in the art can modify or change the above disclosed technical content to equivalent embodiments applied to other fields, but any simple modification, equivalent change and modification made according to the technical essence of the present application to the above embodiments within the technical solution content of the present application still falls within the protection scope of the present application.

Claims

1. A method for producing a three-dimensional model based on a three-dimensional laser point cloud, characterized in that: The following steps are involved: Sequentially read the point cloud data files of the 3D laser point cloud and insert the coordinate value of each point into the octree structure. When the octree node buffer in the memory is full, write the node data in the buffer to the hard disk, traverse the hard disk octree node set from bottom to top, and establish a hierarchical aggregated attribute octree. Based on the hierarchical aggregation attribute octree, a coarse level is selected and representative points and average normals within the nodes are retrieved to construct a function field covering the three-dimensional model, thereby obtaining a basic implicit function field. Isosurface extraction is performed on the basic implicit function field, and an adaptive multi-level-of-detail mesh model is generated according to viewpoint distance and memory limitations. Based on the adaptive multi-level detail grid model, a cell complex sequence is constructed to obtain cell complex basic data, each point in the cell complex basic data is traversed, a topological noise rank is calculated, and a topological feature noise rank list is obtained; Based on the topological feature noise rank list, each topological noise rank in the topological feature noise rank list is compared with a preset rank value threshold to obtain a classified topological feature set, and a topological optimization three-dimensional model is established according to the classified topological feature set.

2. The method for producing a three-dimensional model based on a three-dimensional laser point cloud according to claim 1, characterized in that: The steps for obtaining the hierarchical aggregated attribute octree are: Sequentially read the point cloud data file of the 3D laser point cloud, parse the 3D coordinate value of each point in the file one by one, and insert the coordinate value of each parsed point into the pre-created octree structure in sequence. When the insertion process causes the storage capacity of the octree node buffer in the memory to reach the memory limit threshold, the buffer write mechanism is triggered, and all the octree node data in the current buffer is written to the hard disk storage space as a whole, forming a hard disk octree node set; Based on the hard disk octree node set, the node data stored in the hard disk is traversed layer by layer in a bottom-up order of the octree structure, the spatial geometric data information contained in all child nodes under each parent node is extracted one by one, the coordinate values ​​and normal vector information of the child nodes are merged, and the representative point coordinate value, average normal value and error metric value of each parent node are calculated based on the merged geometric data information to generate an aggregated set of parent node geometric information; Based on the parent node geometric information aggregation set, the representative point coordinate values, average normal values ​​and error metric values ​​of all parent nodes in the parent node geometric information aggregation set are reinserted into the corresponding positions of the parent nodes of the octree structure layer by layer, and the calculation is repeated layer by layer upwards until the step-by-step aggregation of the entire octree structure is completed to form a hierarchical aggregation attribute octree.

3. The method for producing a three-dimensional model based on a three-dimensional laser point cloud according to claim 1, characterized in that: The steps for obtaining the basic implicit function field are: Based on the hierarchical aggregation attribute octree, all nodes of the coarse level are located, representative point coordinate vectors and average normal vectors stored in the nodes are extracted one by one, all representative point coordinate vectors and average normal vectors are traversed, and mapped to each query position in the three-dimensional space, an initial function field covering the three-dimensional model is established, and initial function field data is obtained; The basic implicit function field is calculated according to the initial function field data.

4. The method for producing a three-dimensional model based on a three-dimensional laser point cloud according to claim 1, characterized in that: The steps for obtaining the adaptive multi-level detail grid model are: Based on the basic implicit function field, isosurface extraction is performed for each query position in the three-dimensional space. Combined with the viewpoint distance condition and memory limitation, the deeper level node data of the hierarchical aggregation attribute octree is retrieved to replace the implicit function value of the corresponding area. The local function field is updated and re-divided block by block to generate an adaptive multi-level detail grid model.

5. The method for producing a three-dimensional model based on a three-dimensional laser point cloud according to claim 1, characterized in that: The steps for obtaining the cell complex basic data are: Based on the adaptive multi-level-of-detail mesh model, the geodesic distances from all vertices to the reference point are calculated, and the vertices and connected edges and faces are introduced in ascending order of geodesic distance values, and gradually expanded in this order to form a cell complex, and a complete cell complex sequence including vertices, edges and faces is constructed to obtain the basic data of the cell complex.

6. The method for producing a three-dimensional model based on a three-dimensional laser point cloud according to claim 1, characterized in that: The steps for obtaining the topological feature noise rank list are: Calculate the topological characteristic noise rank based on the cell complex basic data; Based on the topological feature noise rank, all topological feature points are traversed, the topological feature noise ranks are summarized, and a topological feature noise rank list is generated.

7. The method for producing a three-dimensional model based on a three-dimensional laser point cloud according to claim 1, characterized in that: The steps for obtaining the classified topological feature set are: Based on the topological feature noise rank list, the noise rank corresponding to each topological feature in the topological feature noise rank list is extracted item by item, and the noise rank of each topological feature is numerically compared with the pre-set rank value threshold one by one. If the noise rank of the topological feature is greater than the rank value threshold, it is marked as a noise feature; if the noise rank of the topological feature is less than or equal to the rank value threshold, it is marked as a significant structure, and a classified topological feature set is generated.

8. The method for producing a three-dimensional model based on a three-dimensional laser point cloud according to claim 1, wherein: The steps for obtaining the topology optimized three-dimensional model are: Based on the classified topological feature set, all topological features marked as noise features are screened, the spatial position of each noise feature in the adaptive multi-level-of-detail grid model is traced one by one, and all edge units and surface units corresponding to each noise feature are located to form a noise feature location set; Based on the noise feature positioning set, all edge units and surface units in the noise feature positioning set are deleted one by one in the adaptive multi-level detail grid model, and the grid units are reconnected and filled through the local area boundaries after topological cutting to establish a topologically optimized three-dimensional model.

Images