CN109671155B - Surface grid reconstruction method, system and equipment based on point cloud data - Google Patents
Surface grid reconstruction method, system and equipment based on point cloud data Download PDFInfo
- Publication number
- CN109671155B CN109671155B CN201811570403.4A CN201811570403A CN109671155B CN 109671155 B CN109671155 B CN 109671155B CN 201811570403 A CN201811570403 A CN 201811570403A CN 109671155 B CN109671155 B CN 109671155B
- Authority
- CN
- China
- Prior art keywords
- points
- cross
- section
- point cloud
- point
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 51
- 238000012952 Resampling Methods 0.000 claims abstract description 39
- 230000006870 function Effects 0.000 claims description 38
- 230000011218 segmentation Effects 0.000 claims description 30
- 238000004364 calculation method Methods 0.000 claims description 16
- 238000012545 processing Methods 0.000 claims description 3
- 238000010586 diagram Methods 0.000 description 15
- 238000012937 correction Methods 0.000 description 3
- 230000000694 effects Effects 0.000 description 3
- 230000005484 gravity Effects 0.000 description 3
- 239000003814 drug Substances 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000006467 substitution reaction Methods 0.000 description 2
- 238000013459 approach Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000000052 comparative effect Effects 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 238000002474 experimental method Methods 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T17/00—Three dimensional [3D] modelling, e.g. data description of 3D objects
- G06T17/20—Finite element generation, e.g. wire-frame surface description, tesselation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/11—Region-based segmentation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10028—Range image; Depth image; 3D point clouds
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
- Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Computer Graphics (AREA)
- Geometry (AREA)
- Software Systems (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Image Generation (AREA)
Abstract
The invention relates to the technical field of three-dimensional surface reconstruction, in particular to a surface grid reconstruction method, system and equipment based on point cloud data, and aims to improve reconstruction accuracy. The reconstruction method of the invention comprises the following steps: projecting coordinate information of the point cloud data to a z axis, calculating a z coordinate probability density function, and further calculating a curvature local extremum point of the z coordinate probability density function; drawing a plane as a cross section through each local extreme point and perpendicular to the z-axis; dividing between local extreme points according to a preset dividing layer number to obtain a plurality of cross sections perpendicular to a z axis; projecting points in a specific distance range above and below each cross section onto the cross section, and further calculating boundary points of the cross section; fitting a boundary curve of the corresponding cross section according to the calculated boundary points, and resampling; and connecting resampling boundary points on all adjacent cross sections to obtain the reconstructed surface grid. The invention improves the accuracy of the reconstruction of the point cloud data surface.
Description
Technical Field
The invention relates to the technical field of three-dimensional surface reconstruction, in particular to a surface grid reconstruction method, system and equipment based on point cloud data.
Background
The physical model can be quickly and accurately reconstructed by using the surface grid reconstruction method, so that the method is widely applied to various fields of clinical medicine, digital entertainment, industry, agriculture and the like.
The three-dimensional discrete point cloud surface reconstruction algorithm based on the slicing method can rapidly acquire the profile data of the two-dimensional section from the point cloud data, thereby completing three-dimensional surface reconstruction. However, in the prior art, a slice with a fixed thickness is adopted to fit point cloud data or each slice is simply fitted to calculate a contour surface, so that the characteristics of a part of the measured object may be lost, and accurate three-dimensional surface reconstruction cannot be realized.
Disclosure of Invention
In order to solve the problems in the prior art, the invention provides a surface grid reconstruction method, a system and equipment based on point cloud data, which improve the accuracy of the surface reconstruction of the point cloud data.
In one aspect of the present invention, a method for reconstructing a surface mesh based on point cloud data is provided, the method comprising:
taking three-dimensional coordinate data of the point cloud data as an initial value, projecting coordinate information of the point cloud data to a z axis, and calculating a z coordinate probability density function;
calculating curvature local extremum points of the z coordinate probability density function according to the z coordinate probability density function;
drawing a plane as a cross section through each local extreme point and perpendicular to the z-axis; dividing the local extreme points according to the preset dividing layer number to obtain a plurality of cross sections perpendicular to the z axis;
projecting points in a specific distance range above and below each cross section onto the cross section, and further calculating boundary points of the cross section;
fitting a boundary curve of the corresponding cross section according to the calculated boundary points, and resampling according to the fitted boundary curve to obtain resampled boundary points;
and connecting the resampling boundary points on all adjacent cross sections to obtain the reconstructed surface grid.
Preferably, the step of calculating a curvature local extremum point of the z-coordinate probability density function according to the z-coordinate probability density function comprises:
the probability density is modified as follows:
p'＝corr p *p
wherein,,
max (z) represents the maximum value of the z-axis projection, min (z) represents the minimum value of the z-axis projection, and max (p) represents the maximum value of the probability density p;
the curvature of the discrete points on the z-axis projection is calculated as follows:
wherein,,
Δf k (t)＝f(t+kΔt)-f(t+(k-1)Δt)
c (t) is the calculated curvature value, f (t) is a function, k=1, 2,3, Δt is the unit step size of the parameter t of the function f (t), the symbol "×" represents the cross product, the symbol "·" represents the dot product, the symbol "|||||": representing the mode length;
selecting a point with the curvature value larger than a curvature threshold value from discrete points on the z-axis projection as a to-be-determined extreme point;
and selecting a point with the maximum curvature value from the adjacent undetermined extreme points to obtain a local extreme point.
Preferably, the step of dividing between the local extremum points according to the preset dividing layer number to obtain a plurality of cross sections perpendicular to the z axis includes:
calculating the segmentation step length between each pair of adjacent local extreme points according to the following formula:
equally-spaced segmentation is carried out between each pair of adjacent local extremum points according to the corresponding segmentation step length, so that a plurality of cross sections perpendicular to the z axis are obtained;
wherein τ i And τ i+1 The heights of the ith and the (i+1) th local extreme points are respectively, the ith and the (i+1) th local extreme points are two adjacent local extreme points, i=1, 2, …, m and m are the local extreme pointsThe number of extreme points; Δz is the height of the point cloud, and L is the preset number of division layers.
Preferably, the step of "projecting points within a specific distance range above and below each cross-section onto the cross-section, and further calculating the boundary point of the cross-section" includes:
calculating the average distance between the points of the point cloud data according to the following formula:
wherein,,
ρ is the density of the point cloud in the height direction, N is the number of points in the point cloud data, and S is the projection area of the point cloud in the horizontal direction of the geodetic coordinates; Δz is the point cloud height;
projecting points on each cross section within a specific distance from top to bottom, wherein the specific distance is epsilon z /2；
And calculating the boundary point of each cross section by adopting an alpha-shape algorithm of the boundary contour of the two-dimensional plane for the point cloud data of each cross section.
Preferably, the step of "fitting a boundary curve of the corresponding cross section according to the calculated boundary points and resampling according to the fitted boundary curve to obtain resampled boundary points" includes:
performing third-order Bezier curve fitting on each cross section according to the calculated boundary points to obtain a fitted boundary curve equation;
the barycentric coordinates of the cross section are calculated as follows:
wherein,,
is the barycentric coordinates; o (o) j To calculate the boundary point of the cross section, p j ∈P i J is the serial number of the boundary points, j=1, …, r, r is the calculated number of the boundary points of the cross section;
and starting from an x positive half axis, taking the barycentric coordinate of the cross section as a starting point, taking 2 pi/K as an angle interval step length, and solving the intersection point of each ray and the fitted boundary curve to obtain the resampling boundary points, wherein K is the number of the preset resampling boundary points.
Preferably, the step of connecting said resampled border points on all adjacent cross sections, resulting in a reconstructed surface mesh, comprises:
and connecting all the space quadrilaterals by taking four resampling points which correspond to serial numbers in the cross sections of two adjacent layers and are adjacent to each other as a group as four vertexes of the space quadrilaterals, thereby obtaining the reconstructed surface grid.
Preferably, the curvature threshold is 0.03.
In another aspect of the present invention, a surface mesh reconstruction system based on point cloud data is provided, the system comprising:
the probability density function calculation module is used for taking three-dimensional coordinate data of the point cloud data as an initial value, projecting coordinate information of the point cloud data to a z axis, and calculating a z coordinate probability density function;
the extreme point calculation module is used for calculating curvature local extreme points of the z-coordinate probability density function according to the z-coordinate probability density function;
the segmentation module is used for drawing a plane which passes through each local extreme point and is perpendicular to the z axis as a cross section; dividing the local extreme points according to the preset dividing layer number to obtain a plurality of cross sections perpendicular to the z axis;
the boundary point calculation module is used for projecting points in the upper specific distance range and the lower specific distance range of each cross section onto the cross section, and further calculating boundary points of the cross section;
the fitting and resampling module is used for fitting a boundary curve of the corresponding cross section according to the calculated boundary points and resampling according to the fitted boundary curve to obtain resampled boundary points;
and the surface grid reconstruction module is used for connecting the resampling boundary points on all adjacent cross sections so as to obtain a reconstructed surface grid.
In a third aspect the invention proposes a storage device comprising a memory adapted to store a program adapted to be loaded and executed by a processor for implementing the above-described point cloud data based surface mesh reconstruction method.
In a fourth aspect of the present invention, there is provided a processing apparatus comprising:
a processor adapted to load a program;
a memory adapted to store the program;
the program is adapted to be loaded and executed by the processor to implement the above-described method of reconstructing a surface mesh based on point cloud data.
Compared with the closest prior art, the invention has the following beneficial effects:
according to the surface grid reconstruction method based on the point cloud data, the unequal interval segmentation strategy is given by utilizing the density information of the three-dimensional point cloud in the vertical direction, and meanwhile, a more accurate cross section profile is given by utilizing an alpha-shape algorithm for solving the boundary profile by utilizing a two-dimensional plane and third-order Bezier curve fitting, so that the accuracy of the surface reconstruction of the point cloud data is improved. The surface grid reconstruction method can reconstruct a physical model rapidly and accurately, and can be widely applied to various fields of clinical medicine, digital entertainment, industry, agriculture and the like.
Drawings
FIG. 1 is a schematic diagram of the main steps of an embodiment of a point cloud data based surface mesh reconstruction method of the present invention;
FIG. 2 is a schematic diagram of a probability density function of point cloud data of a crown in an embodiment of the invention;
FIG. 3 is a schematic diagram comparing a prior art segmentation method with the segmentation method of the present invention;
FIG. 4 is a schematic diagram of a method for resampling a boundary curve in an embodiment of the invention;
FIG. 5 is a schematic diagram of a method for reconstructing a surface grid according to the center of gravity of a cross section and resampling points in an embodiment of the invention;
fig. 6 is a schematic diagram of the main components of an embodiment of the point cloud data based surface mesh reconstruction system of the present invention.
Detailed Description
Preferred embodiments of the present invention are described below with reference to the accompanying drawings. It should be understood by those skilled in the art that these embodiments are merely for explaining the technical principles of the present invention, and are not intended to limit the scope of the present invention.
In order to solve the defects of an equidistant slicing method in the traditional method, the invention provides an unequal interval segmentation strategy by utilizing the density information of the three-dimensional point cloud in the vertical direction, and simultaneously provides a more accurate cross section profile by utilizing an alpha-shape algorithm for solving the boundary profile by utilizing a two-dimensional plane and a third-order Bezier curve fitting, thereby ensuring that the reconstruction effect is more ideal.
Fig. 1 is a schematic diagram of main steps of an embodiment of a point cloud data-based surface mesh reconstruction method of the present invention. As shown in fig. 1, the surface mesh reconstruction method of the present embodiment includes steps S1 to S6:
and S1, taking three-dimensional coordinate data of the point cloud data as an initial value, projecting coordinate information of the point cloud data to a z axis, and calculating a z coordinate probability density function.
Fig. 2 is a schematic diagram of a probability density function of point cloud data of a crown in an embodiment of the invention. As shown in fig. 2, the left half is a schematic diagram of a crown, and the right half is a probability density function corresponding to point cloud data of the crown. Wherein the vertical axis is the z-coordinate of the point cloud data, and the horizontal axis is the probability density of the point cloud data projected to the z-axis.
And S2, calculating curvature local extremum points of the z coordinate probability density function according to the z coordinate probability density function.
Due toThe dimensions of the z-axis projection coordinate value and the probability density value are different, and a correction coefficient corr is added on the probability density value p p 。
Step S2 may specifically include steps S21-S24:
step S21, correcting the probability density according to the formula (1):
p'＝corr p *p (1)
wherein,,
max (z) represents the maximum value of the z-axis projection, min (z) represents the minimum value of the z-axis projection, and max (p) represents the maximum value of the probability density p.
Step S22, calculating the curvature of the discrete point on the z-axis projection by using the idea of calculating the curvature of the discrete point by using the difference, as shown in formula (3):
wherein,,
Δf k (t)＝f(t+kΔt)-f(t+(k-1)Δt) (4)
c (t) is the calculated curvature value, f (t) is a function, k=1, 2,3, Δt is the unit step size of the parameter t of the function f (t), the symbol "×" represents the cross product, the symbol "·" represents the dot product, the symbol "|||||": representing the mode length.
And S23, selecting a point with the curvature value larger than a curvature threshold value from discrete points on the z-axis projection as a undetermined extreme point.
When the selected curvature threshold value is smaller, the number of local curvature extreme points is increased, so that the time cost is increased; when the curvature threshold value is large, the number of local curvature extreme points becomes small, and the calculation accuracy is insufficient. According to the invention, through a large number of experiments, for the corrected z-coordinate probability density function diagram, when the curvature threshold is set to be 0.03, a better local extremum point selection effect can be obtained.
In step S24, since the curvature change range of the adjacent points is small, the points with curvature larger than the threshold value are secondarily screened, the point represented by the maximum curvature is selected from the adjacent points, and the points adjacent to the maximum curvature are discarded, thereby obtaining the local extremum point LEz (Local Extremum).
Step S3, drawing a plane as a cross section for each local extreme point passing through the point and perpendicular to the z axis; and dividing between local extreme points according to the preset dividing layer number to obtain a plurality of cross sections perpendicular to the z axis. The method specifically comprises the steps S31-S33:
step S31, for each local extremum point, a plane is drawn through the point and perpendicular to the z-axis as a cross-section.
Step S32, calculating a segmentation step length between each pair of adjacent local extreme points according to a formula (5):
wherein τ i And τ i+1 The heights of the ith and the (i+1) th local extremum points (namely coordinate values of the two points on the z axis) are respectively, the ith and the (i+1) th local extremum points are two adjacent local extremum points, i=1, 2, …, m, m is the number of the local extremum points calculated in the step S2, Δz is the point cloud height, and Δz=max (z) -min (z), wherein max (z) and min (z) respectively represent the maximum value and the minimum value of the point cloud data projected on the z axis; l is a preset number of division layers (20 in this embodiment);
representing an upward rounding.
And S33, equally-spaced segmentation is carried out between each pair of adjacent local extreme points according to the corresponding segmentation step length, so that a plurality of cross sections perpendicular to the z axis are obtained.
According to the above algorithm, the number of division layers obtained finally is usually greater than n because of the rounding-up r 。
FIG. 3 is a prior art segmentation method and the segmentation of the present inventionA comparative schematic of the process. As shown in fig. 3, the left graph is a schematic diagram of the equidistant dividing method using the prior art, and the right graph is a schematic diagram of the unequal dividing method using the method. Wherein z is 1 、z 2 、z 3 、z 4 、z 5 Representing different cross sections, respectively. On a longitudinal sectional view, the intersection points of the cross section and the original contour are sequentially connected to obtain the segmented contour, so that the contour segmented by the method is closer to the original contour, and the error is smaller.
And S4, projecting points in a specific distance range above and below each cross section onto the cross section, and further calculating the boundary point of the cross section. The method specifically comprises the steps S41-S43:
step S41, calculating the average distance between the points of the point cloud data according to the formula (6):
wherein,,
ρ is the density of the point cloud in the height direction, N is the number of points in the point cloud data, and S is the projected area of the point cloud in the geodetic coordinate horizontal direction; Δz is the point cloud height, and Δz=max (z) -min (z), max (z) and min (z) represent the maximum and minimum values, respectively, of the projection of the point cloud data on the z-axis.
In step S42, points within a specific distance range above and below each cross section are projected onto the cross section, so that both the points on each cross section are guaranteed, and the points on the cross section can approximately reflect the cross section shape. In this embodiment, the specific distance is ε z /2。
Step S43, for the point cloud data on each cross section, calculating the boundary points of the cross section by adopting an alpha-shape algorithm of the boundary profile of the two-dimensional plane, and orderly storing the boundary points into a boundary point set P i Wherein i is a pair ofSerial number of cross section. The term "orderly deposit" as used herein means that the boundary points are deposited sequentially in the order of the boundary contour.
And S5, fitting a boundary curve of the corresponding cross section according to the calculated boundary points, and resampling according to the fitted boundary curve to obtain resampled boundary points. The method specifically comprises the steps of S51-S53:
and S51, performing third-order Bezier curve fitting on each cross section according to the calculated boundary points to obtain a fitted boundary curve equation.
Step S52, calculating the barycentric coordinates of the cross section according to the formula (8):
wherein,,
a barycentric coordinate; o (o) j To calculate the boundary point of the cross section, p j ∈P i J is the serial number of the boundary points, j=1, …, r, r is the calculated number of the boundary points of the cross section.
Step S53, starting from the x positive half axis, establishing rays with the barycentric coordinates of the cross section as a starting point and with 2 pi/K as an angle interval step length, and obtaining the intersection point of each ray and the fitted boundary curve to obtain a resampling boundary point
And sequentially storing (sequentially storing in order of resampling boundary points on the boundary curve) the resampling boundary point set P i 'in'.
Wherein K is the number of preset resampling boundary points,
k=1, …, K, which is the kth resampling boundary point of the ith cross section.
FIG. 4 is a method of resampling a boundary curve in the present embodimentSchematic diagram. Wherein the method comprises the steps of
Is the center of gravity of the cross section,
resampling points.
And S6, connecting resampling boundary points on all adjacent cross sections, so as to obtain the reconstructed surface grid.
In this embodiment, resampling points corresponding to serial numbers and adjacent to each other in two adjacent layers of cross sections are taken as a group of four points, and are used as four vertexes of a space quadrilateral to connect all the space quadrilaterals, so as to obtain a reconstructed surface grid.
Fig. 5 is a schematic diagram of a method for reconstructing a surface grid according to the center of gravity of a cross section and resampling points in the present embodiment. As shown in FIG. 5, the i-th layer and the i+1-th layer are adjacent two-layer cross sections, and their barycentric coordinates are respectively
And->
When reconstructing the surface grid, four resampling points with corresponding serial numbers and adjacent positions in the two adjacent layers of cross sections are added>
And
as a group, four vertices of a space quadrilateral, similarly, +.>
And->
As a group, as four vertexes of the space quadrangle, all the space quadrangles are connected to obtain the reconstructed space quadrangleSurface mesh.
Although the steps are described in the above-described sequential order in the above-described embodiments, it will be appreciated by those skilled in the art that in order to achieve the effects of the present embodiments, the steps need not be performed in such order, and may be performed simultaneously (in parallel) or in reverse order, and such simple variations are within the scope of the present invention.
Based on the same technical conception as the method embodiment, the invention also provides a surface grid reconstruction system, and the method is specifically described below.
Fig. 6 is a schematic diagram of the main components of an embodiment of the point cloud data based surface mesh reconstruction system of the present invention. As shown in fig. 6, the point cloud data-based surface mesh reconstruction system 1 of the present embodiment includes: probability density function calculation module 10, extreme point calculation module 20, segmentation module 30, boundary point calculation module 40, fitting and resampling module 50, and surface mesh reconstruction module 60.
The probability density function calculation module 10 is configured to project coordinate information of the point cloud data to a z axis with three-dimensional coordinate data of the point cloud data as an initial value, and calculate a z coordinate probability density function; the extremum point calculating module 20 is configured to calculate a curvature local extremum point of the z-coordinate probability density function according to the z-coordinate probability density function; the segmentation module 30 draws a plane as a cross-section through each local extremum point and perpendicular to the z-axis; dividing between local extreme points according to the preset dividing layer number to obtain a plurality of cross sections perpendicular to the z axis; the boundary point calculating module 40 is configured to project points within a specific distance range above and below each cross section onto the cross section, so as to calculate boundary points of the cross section; the fitting and resampling module 50 is configured to fit a boundary curve of the corresponding cross section according to the calculated boundary points, and resample according to the fitted boundary curve to obtain resampled boundary points; the surface mesh reconstruction module 60 is configured to connect resampled boundary points on all adjacent cross sections to obtain a reconstructed surface mesh.
In this embodiment, the extreme point calculating module 20 may specifically include: the device comprises a correction unit, a curvature calculation unit, a to-be-determined extreme point selection unit and a local extreme point selection unit.
Wherein the correction unit is used for correcting the probability density according to the methods shown in formulas (1) - (2); the curvature calculating unit is used for calculating the curvature of the discrete points on the z-axis projection according to the formulas (3) - (4); the undetermined extreme point selecting unit is used for selecting a point with a curvature value larger than a curvature threshold value from discrete points on the z-axis projection as an undetermined extreme point; the local extremum point selecting unit is used for selecting the point with the maximum curvature value from the adjacent undetermined extremum points to obtain the local extremum point.
In this embodiment, the segmentation module 30 may specifically include: the device comprises an extreme point cross section generating unit, a segmentation step length calculating unit and a segmentation unit.
The extreme point cross section generating unit is used for drawing a plane which passes through each local extreme point and is perpendicular to the z axis as a cross section; the segmentation step length calculation unit is used for calculating the segmentation step length between each pair of adjacent local extremum points according to the formula (5); the segmentation unit is used for carrying out equidistant segmentation between each pair of adjacent local extreme points according to the corresponding segmentation step length, so as to obtain a plurality of cross sections perpendicular to the z axis.
In this embodiment, the boundary point calculating module 40 may specifically include: an inter-point average distance calculation unit, a projection unit, and a boundary point calculation unit.
The point average distance calculation unit is used for calculating the point average distance of the point cloud data according to formulas (6) - (7); the projection unit is used for projecting points in a specific distance range on and below each cross section onto the cross section; the boundary point calculating unit is used for calculating the boundary point of each cross section by adopting an alpha-shape algorithm of the two-dimensional plane boundary contour for the point cloud data of each cross section.
In this embodiment, the fitting and resampling module 50 may specifically include: fitting unit, barycentric coordinates calculating unit and fitting unit.
The fitting unit is used for carrying out third-order Bezier curve fitting on each cross section according to the calculated boundary points to obtain a fitted boundary curve equation; the barycentric coordinate calculation unit is used for calculating barycentric coordinates of the cross section according to a formula (8); the resampling unit is used for establishing rays by taking the barycentric coordinates of the cross section as a starting point and taking 2 pi/K as an angle interval step length from the x positive half axis, solving the intersection point of each ray and the fitted boundary curve, and obtaining a resampling boundary point.
In this application, the surface mesh reconstruction system is divided into modules and units, so that the functions related to the technical solution of the present invention are better understood, and in practice, the functions corresponding to the modules may be loaded and executed by single or multiple hardware.
Based on the above surface mesh reconstruction method, the present invention also proposes a storage device, comprising a memory adapted to store a program, said program being adapted to be loaded and executed by a processor to implement the above surface mesh reconstruction method based on point cloud data.
Further, the present invention also proposes a processing apparatus including: a processor and a memory. Wherein the processor is adapted to load a program and the memory is adapted to store said program, said program being adapted to be loaded and executed by said processor for implementing the above-described point cloud data based surface mesh reconstruction method.
Those of skill in the art will appreciate that the various illustrative method steps, modules, and elements described in connection with the embodiments disclosed herein may be implemented as electronic hardware, computer software, or combinations of both, and that the various illustrative components and steps have been described above generally in terms of functionality in order to clearly illustrate the interchangeability of electronic hardware and software. Whether such functionality is implemented as electronic hardware or software depends upon the particular application and design constraints imposed on the solution. Those skilled in the art may implement the described functionality using different approaches for each particular application, but such implementation is not intended to be limiting.
Thus far, the technical solution of the present invention has been described in connection with the preferred embodiments shown in the drawings, but it is easily understood by those skilled in the art that the scope of protection of the present invention is not limited to these specific embodiments. Equivalent modifications and substitutions for related technical features may be made by those skilled in the art without departing from the principles of the present invention, and such modifications and substitutions will be within the scope of the present invention.
Claims (7)
1. A method for reconstructing a surface mesh based on point cloud data, the method comprising:
taking three-dimensional coordinate data of the point cloud data as an initial value, projecting coordinate information of the point cloud data to a z axis, and calculating a z coordinate probability density function;
calculating curvature local extremum points of the z coordinate probability density function according to the z coordinate probability density function;
the probability density is modified as follows:
p′＝corr p *p
wherein,,
max (z) represents the maximum value of the z-axis projection, min (z) represents the minimum value of the z-axis projection, and max (p) represents the maximum value of the probability density p;
the curvature of the discrete points on the z-axis projection is calculated as follows:
wherein,,
Δf k (t)＝f(t+kΔt)-f(t+(k-1)Δt)
c (t) is the calculated curvature value, f (t) is a function, k=1, 2,3, Δt is the unit step size of the parameter t of the function f (t), the symbol "×" represents the cross product, the symbol "·" represents the dot product, the symbol "|||||": representing the mode length;
selecting a point with the curvature value larger than a curvature threshold value from discrete points on the z-axis projection as a to-be-determined extreme point;
selecting a point with the maximum curvature value from the adjacent undetermined extreme points to obtain a local extreme point;
drawing a plane as a cross section through each local extreme point and perpendicular to the z-axis; dividing the local extreme points according to the preset dividing layer number to obtain a plurality of cross sections perpendicular to the z axis;
calculating the segmentation step length between each pair of adjacent local extreme points according to the following formula:
equally-spaced segmentation is carried out between each pair of adjacent local extremum points according to the corresponding segmentation step length, so that a plurality of cross sections perpendicular to the z axis are obtained;
wherein τ i And τ i+1 The heights of the ith and the (i+1) th local extremum points are respectively, the ith and the (i+1) th local extremum points are two adjacent local extremum points, i=1, 2, & m, m is the number of the local extremum points; Δz is the height of the point cloud, and L is the preset number of division layers;
projecting points within a specific distance range above and below each cross section onto the cross section, and further calculating boundary points of the cross section, wherein the method comprises the following steps: calculating the average distance between the points of the point cloud data according to the following formula:
wherein,,
ρ is the density of the point cloud in the height direction, N is the number of points in the point cloud data, and S is the projection area of the point cloud in the horizontal direction of the geodetic coordinates; Δz is the point cloud height;
projecting points on each cross section within a specific distance from top to bottom, wherein the specific distance is epsilon z /2；
For the point cloud data on each cross section, calculating the boundary point of the cross section by adopting an alpha-shape algorithm of the boundary contour of the two-dimensional plane;
fitting a boundary curve of the corresponding cross section according to the calculated boundary points, and resampling according to the fitted boundary curve to obtain resampled boundary points;
and connecting the resampling boundary points on all adjacent cross sections to obtain the reconstructed surface grid.
2. The method for reconstructing a surface mesh based on point cloud data according to claim 1, wherein the step of fitting a boundary curve of a corresponding cross section based on the calculated boundary points and resampling based on the fitted boundary curve to obtain resampled boundary points comprises:
performing third-order Bezier curve fitting on each cross section according to the calculated boundary points to obtain a fitted boundary curve equation;
the barycentric coordinates of the cross section are calculated as follows:
wherein,,
is the barycentric coordinates; o (o) j For the calculated boundary points of the cross section, j is the serial number of the boundary points, j=1, …, r, r is the calculated number of the boundary points of the cross section;
and starting from an x positive half axis, taking the barycentric coordinate of the cross section as a starting point, taking 2 pi/K as an angle interval step length, and solving the intersection point of each ray and the fitted boundary curve to obtain the resampling boundary points, wherein K is the number of the preset resampling boundary points.
3. The method of reconstructing a surface mesh based on point cloud data according to claim 1, wherein the step of connecting the resampled boundary points on all adjacent cross sections to obtain a reconstructed surface mesh comprises:
and connecting all the space quadrilaterals by taking four resampling points which correspond to serial numbers in the cross sections of two adjacent layers and are adjacent to each other as a group as four vertexes of the space quadrilaterals, thereby obtaining the reconstructed surface grid.
4. The method of point cloud data based surface mesh reconstruction of claim 1, wherein the curvature threshold is 0.03.
5. A point cloud data-based surface mesh reconstruction system, the system comprising:
the probability density function calculation module is used for taking three-dimensional coordinate data of the point cloud data as an initial value, projecting coordinate information of the point cloud data to a z axis, and calculating a z coordinate probability density function;
the extreme point calculation module is used for calculating curvature local extreme points of the z-coordinate probability density function according to the z-coordinate probability density function;
the probability density is modified as follows:
p'＝corr p *p
wherein,,
max (z) represents the maximum value of the z-axis projection, min (z) represents the minimum value of the z-axis projection, and max (p) represents the maximum value of the probability density p;
the curvature of the discrete points on the z-axis projection is calculated as follows:
wherein,,
Δf k (t)＝f(t+kΔt)-f(t+(k-1)Δt)
c (t) is the calculated curvature value, f (t) is a function, k=1, 2,3, Δt is the unit step size of the parameter t of the function f (t), the symbol "×" represents the cross product, the symbol "·" represents the dot product, the symbol "|||||": representing the mode length;
selecting a point with the curvature value larger than a curvature threshold value from discrete points on the z-axis projection as a to-be-determined extreme point;
selecting a point with the maximum curvature value from the adjacent undetermined extreme points to obtain a local extreme point;
the segmentation module is used for drawing a plane which passes through each local extreme point and is perpendicular to the z axis to serve as a cross section; dividing the local extreme points according to the preset dividing layer number to obtain a plurality of cross sections perpendicular to the z axis;
calculating the segmentation step length between each pair of adjacent local extreme points according to the following formula:
equally-spaced segmentation is carried out between each pair of adjacent local extremum points according to the corresponding segmentation step length, so that a plurality of cross sections perpendicular to the z axis are obtained;
wherein τ i And τ i+1 The height of the ith local extreme point and the i+1th local extreme point are respectively, the ith local extreme point and the i+1th local extreme point are two adjacent local extreme points, and i=1, 2, …, m and m are the number of the local extreme points; Δz is the height of the point cloud, and L is the preset number of division layers;
the boundary point calculating module is used for projecting points in a specific distance range above and below each cross section onto the cross section, and further calculating boundary points of the cross section, and comprises the following steps: calculating the average distance between the points of the point cloud data according to the following formula:
wherein,,
ρ is the density of the point cloud in the height direction, N is the number of points in the point cloud data, and S is the projection area of the point cloud in the horizontal direction of the geodetic coordinates; Δz is the point cloud height;
projecting points on each cross section within a specific distance from top to bottom, wherein the specific distance is epsilon z /2；
For the point cloud data on each cross section, calculating the boundary point of the cross section by adopting an alpha-shape algorithm of the boundary contour of the two-dimensional plane;
the fitting and resampling module is used for fitting a boundary curve of the corresponding cross section according to the calculated boundary points and resampling according to the fitted boundary curve to obtain resampled boundary points;
and the surface grid reconstruction module is used for connecting the resampling boundary points on all adjacent cross sections so as to obtain a reconstructed surface grid.
6. A storage device comprising a memory adapted to store a program, wherein the program is adapted to be loaded and executed by a processor to implement the point cloud data based surface mesh reconstruction method of any one of claims 1-4.
7. A processing apparatus, comprising:
a processor adapted to load a program;
a memory adapted to store the program;
characterized in that the program is adapted to be loaded and executed by the processor for implementing the point cloud data based surface mesh reconstruction method of any one of claims 1-4.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201811570403.4A CN109671155B (en) | 2018-12-21 | 2018-12-21 | Surface grid reconstruction method, system and equipment based on point cloud data |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201811570403.4A CN109671155B (en) | 2018-12-21 | 2018-12-21 | Surface grid reconstruction method, system and equipment based on point cloud data |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN109671155A CN109671155A (en) | 2019-04-23 |
| CN109671155B true CN109671155B (en) | 2023-06-23 |
Family
ID=66145781
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201811570403.4A Active CN109671155B (en) | 2018-12-21 | 2018-12-21 | Surface grid reconstruction method, system and equipment based on point cloud data |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN109671155B (en) |
Families Citing this family (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110264565B (en) * | 2019-05-27 | 2021-07-30 | 浙江大学 | Three-dimensional reconstruction method based on half-peak probability density distribution |
| CN110211171A (en) * | 2019-06-18 | 2019-09-06 | 华志微创医疗科技（北京）有限公司 | The method and device of focal area is extracted from medical image |
| CN110664438B (en) * | 2019-10-22 | 2021-09-10 | 深圳瀚维智能医疗科技有限公司 | Ultrasonic scanning track planning method and device, storage medium and computer equipment |
| CN113033270B (en) * | 2019-12-27 | 2023-03-17 | 深圳大学 | 3D object local surface description method and device adopting auxiliary axis and storage medium |
| CN112308761A (en) * | 2020-11-13 | 2021-02-02 | 济南浪潮高新科技投资发展有限公司 | Robot depth camera point cloud down-sampling filtering method based on slicing and interpolation |
| CN112884903A (en) * | 2021-03-22 | 2021-06-01 | 浙江浙能兴源节能科技有限公司 | Driving three-dimensional modeling system and method thereof |
| CN113610916B (en) * | 2021-06-17 | 2024-04-12 | 同济大学 | Irregular object volume determination method and system based on point cloud data |
| CN116309788A (en) * | 2023-03-23 | 2023-06-23 | 中铁武汉电气化局集团有限公司 | 3D point cloud railway contact line cross-sectional area detection method |
Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2011079421A1 (en) * | 2009-12-30 | 2011-07-07 | 中国科学院自动化研究所 | Method for global paremeterization and quadrilateral gridding of point cloud data |
| CN106910216A (en) * | 2017-03-03 | 2017-06-30 | 广东工业大学 | A kind of grid surface based on ship hull plate cloud data determines method and device |
-
2018
- 2018-12-21 CN CN201811570403.4A patent/CN109671155B/en active Active
Patent Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2011079421A1 (en) * | 2009-12-30 | 2011-07-07 | 中国科学院自动化研究所 | Method for global paremeterization and quadrilateral gridding of point cloud data |
| CN106910216A (en) * | 2017-03-03 | 2017-06-30 | 广东工业大学 | A kind of grid surface based on ship hull plate cloud data determines method and device |
Non-Patent Citations (3)
| Title |
|---|
| Parameter free torsion estimation of curves in 3D images;Blankenburg C，et al;《IEEE International Conference on Image Processing》;20161231;全文 * |
| 激光测量果树冠层体积方法与试验研究;黄健等;《中国农业工程学会2011年学术年会论文集》;20111022;第1-6页 * |
| 迭代切片算法在点云曲面拟合中的应用;张崇军等;《地理空间信息》;20180131(第1期);第1-4页 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN109671155A (en) | 2019-04-23 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN109671155B (en) | Surface grid reconstruction method, system and equipment based on point cloud data | |
| JP5222946B2 (en) | How to make object model compatible with 3D point cloud | |
| US9373190B2 (en) | High-quality stereo reconstruction featuring depth map alignment and outlier identification | |
| CN107330903B (en) | Skeleton extraction method of human point cloud model | |
| CN106023298B (en) | Point cloud Rigid Registration method based on local Poisson curve reestablishing | |
| CN109697731B (en) | Plant canopy volume calculation method and system based on point cloud data and related equipment | |
| US7825925B2 (en) | Method and system for repairing triangulated surface meshes | |
| CN107610061B (en) | Edge-preserving point cloud hole repairing method based on two-dimensional projection | |
| Kwok et al. | Efficient optimization of common base domains for cross parameterization | |
| CN107016649A (en) | A kind of vision data complementing method estimated based on local low-rank tensor | |
| GB2502680A (en) | Three-Dimensional Shape Measurement using patterned image projection | |
| CN104134220A (en) | Low-altitude remote sensing image high-precision matching method with consistent image space | |
| CN113155054B (en) | Automatic three-dimensional scanning planning method for surface structured light | |
| CN108122280A (en) | The method for reconstructing and device of a kind of three-dimensional point cloud | |
| Holdstein et al. | Three-dimensional surface reconstruction using meshing growing neural gas (MGNG) | |
| CN109671154A (en) | The curved surface non-iterative gridding method again that triangle gridding indicates | |
| CN109551768A (en) | A kind of data processing method of the 3D printing file based on STL | |
| Mao et al. | Robust surface reconstruction of teeth from raw pointsets | |
| JP3950376B2 (en) | Shape model generation method and apparatus from three-dimensional point group, shape model generation program from three-dimensional point group, and recording medium recording the program | |
| CN107194998B (en) | A kind of method of multi-layer three-dimension point cloud single-layered | |
| KR20200142392A (en) | Method for Estimating Rotation and Movement Position Between Corresponding 3D Marker | |
| JP2002008014A (en) | Method and device for extracting three-dimensional shape, and recording medium | |
| CN111667581B (en) | Processing method, device, equipment and storage medium of 3D model | |
| Kim et al. | Mesh modelling of 3D point cloud from UAV images by point classification and geometric constraints | |
| CN114964176A (en) | Landform surveying and mapping method for permanent shadow area of moon |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |