CN102999937A - Curved planar reconstruction method for cardiac scattered-point cloud data - Google Patents
Curved planar reconstruction method for cardiac scattered-point cloud data Download PDFInfo
- Publication number
- CN102999937A CN102999937A CN 201110264828 CN201110264828A CN102999937A CN 102999937 A CN102999937 A CN 102999937A CN 201110264828 CN201110264828 CN 201110264828 CN 201110264828 A CN201110264828 A CN 201110264828A CN 102999937 A CN102999937 A CN 102999937A
- Authority
- CN
- China
- Prior art keywords
- point
- cloud data
- point cloud
- points
- function
- 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.)
- Pending
Links
- 238000000034 method Methods 0.000 title claims abstract description 50
- 230000000747 cardiac effect Effects 0.000 title claims abstract description 16
- 239000013598 vector Substances 0.000 claims abstract description 79
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 20
- 239000011159 matrix material Substances 0.000 claims description 10
- 238000001914 filtration Methods 0.000 claims description 7
- 239000000284 extract Substances 0.000 claims description 2
- 238000010606 normalization Methods 0.000 claims description 2
- 201000010099 disease Diseases 0.000 abstract description 5
- 208000037265 diseases, disorders, signs and symptoms Diseases 0.000 abstract description 5
- 238000004458 analytical method Methods 0.000 abstract description 4
- 238000005520 cutting process Methods 0.000 abstract description 4
- 238000003745 diagnosis Methods 0.000 abstract description 4
- 238000004088 simulation Methods 0.000 abstract description 3
- 206010003119 arrhythmia Diseases 0.000 abstract description 2
- 230000006793 arrhythmia Effects 0.000 abstract description 2
- 208000019622 heart disease Diseases 0.000 abstract description 2
- 230000006870 function Effects 0.000 description 39
- 238000005070 sampling Methods 0.000 description 14
- 238000002474 experimental method Methods 0.000 description 7
- 238000010586 diagram Methods 0.000 description 6
- 230000008569 process Effects 0.000 description 6
- 238000005516 engineering process Methods 0.000 description 5
- 230000000694 effects Effects 0.000 description 4
- 230000000007 visual effect Effects 0.000 description 4
- 238000012545 processing Methods 0.000 description 3
- 238000004364 calculation method Methods 0.000 description 2
- 239000002131 composite material Substances 0.000 description 2
- 238000004519 manufacturing process Methods 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 230000008520 organization Effects 0.000 description 2
- 238000007781 pre-processing Methods 0.000 description 2
- 230000009467 reduction Effects 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- IYLGZMTXKJYONK-ACLXAEORSA-N (12s,15r)-15-hydroxy-11,16-dioxo-15,20-dihydrosenecionan-12-yl acetate Chemical compound O1C(=O)[C@](CC)(O)C[C@@H](C)[C@](C)(OC(C)=O)C(=O)OCC2=CCN3[C@H]2[C@H]1CC3 IYLGZMTXKJYONK-ACLXAEORSA-N 0.000 description 1
- 101100391182 Dictyostelium discoideum forI gene Proteins 0.000 description 1
- 230000001154 acute effect Effects 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000005094 computer simulation Methods 0.000 description 1
- 238000012885 constant function Methods 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 238000001514 detection method Methods 0.000 description 1
- 238000009499 grossing Methods 0.000 description 1
- 239000000203 mixture Substances 0.000 description 1
- 238000005457 optimization Methods 0.000 description 1
- IYLGZMTXKJYONK-UHFFFAOYSA-N ruwenine Natural products O1C(=O)C(CC)(O)CC(C)C(C)(OC(C)=O)C(=O)OCC2=CCN3C2C1CC3 IYLGZMTXKJYONK-UHFFFAOYSA-N 0.000 description 1
- 238000001356 surgical procedure Methods 0.000 description 1
- HOWHQWFXSLOJEF-MGZLOUMQSA-N systemin Chemical compound NCCCC[C@H](N)C(=O)N[C@@H](CCSC)C(=O)N[C@@H](CCC(N)=O)C(=O)N[C@@H]([C@@H](C)O)C(=O)N[C@@H](CC(O)=O)C(=O)OC(=O)[C@@H]1CCCN1C(=O)[C@H]1N(C(=O)[C@H](CC(O)=O)NC(=O)[C@H](CCCN=C(N)N)NC(=O)[C@H](CCCCN)NC(=O)[C@H](CO)NC(=O)[C@H]2N(CCC2)C(=O)[C@H]2N(CCC2)C(=O)[C@H](CCCCN)NC(=O)[C@H](CO)NC(=O)[C@H](CCC(N)=O)NC(=O)[C@@H](NC(=O)[C@H](C)N)C(C)C)CCC1 HOWHQWFXSLOJEF-MGZLOUMQSA-N 0.000 description 1
- 108010050014 systemin Proteins 0.000 description 1
Images
Landscapes
- Processing Or Creating Images (AREA)
Abstract
The invention relates to the field of curved planar reconstruction, in particular to a curved planar reconstruction method for cardiac scattered-point cloud data. The curved planar reconstruction method for cardiac scattered-point cloud data includes the steps of firstly, defining the data representation form of acquired point cloud data, establishing octree topological relation, and defining a space function; secondly, establishing a vector field, selecting local adjacent planes, and determining vertex normal vectors according to an average of planar normal vectors of the local adjacent planes; thirdly, uniformizing the vertex normal vectors; and fourthly, solving a Poisson equation, extracting a contour surface by MC (marching cubes) algorithm according to an indicator function and gradient thereof, and reconstructing a cardiac three-dimensional surface model. By the curved planar reconstruction method for cardiac scattered-point cloud data, cardiac spatial position relation can be displayed visually and accurately, optional display, zooming, rotation and cutting in three dimensions can be achieved, cardiac diseases can be detected, computer diagnostic analysis for diseases such as arrhythmia and surgical accurate simulation are achieved, and safety and reliability in disease diagnosis are improved.
Description
Technical Field
The invention relates to the field of curved surface reconstruction, in particular to a heart curved surface reconstruction method.
Background
Curved surface reconstruction is a process of rapidly, accurately and efficiently acquiring a complex three-dimensional surface model from three-dimensional data acquired from a real object or a sample, and is currently applied to reverse engineering. The reverse engineering is a technology that an electronic instrument is used for collecting original data from a real object or a sample, then computer equipment is used for converting the data into a conceptual model, and relevant information such as analysis, modification and optimization are carried out on a product on the basis, and the technology belongs to a relatively independent category in CAD/CAM. At present, curved surface reconstruction technology has been widely used in many fields, such as production and manufacturing, animation games, mold design, medical simulation, archaeology, and the like.
Since three-dimensional data collected by a measuring device is usually dense, it is called Point Cloud (Point Cloud) visually, and the Point Cloud can be regarded as a set of points in a three-dimensional space, and each Point Cloud has coordinate values in x, y, and z directions. The point cloud data can be divided into ordered point clouds and scattered point clouds according to different organization forms of the data. The reconstruction of the ordered point cloud refers to that a curved surface passing through a series of sampling points sequentially is constructed for a series of sampling points collected on a known curved surface and is used as an approximation of an original curved surface; in practice, due to the limitation of the acquisition equipment and the difference of the acquisition mode, the acquired data point set has no specific organization form, and is called as scattered point cloud. Because the sampling density is large enough, the acquired data has redundancy and noise, so that an accurate curved surface passing through a sampling point cannot be found, and for scattered point clouds, a surface capable of reflecting the shape of the original point cloud is usually reconstructed according to a sampling point set to serve as a reconstruction curved surface.
At present, research on point cloud data surface reconstruction has been widely developed at home and abroad, and two major methods, namely a method based on a combined structure and a method based on an implicit function, are available. Methods based on composite structures, e.g. Delaunay triangulation,
shape, Crust, etc. this kind of method is to build a triangular mesh to interpolate all or most points to reconstruct the curved surface. The method based on implicit functions, such as fast fourier transform, radial basis function and the like, directly reconstructs an approximate surface by defining a piecewise function, setting internal and external thresholds of a model and extracting an isosurface.
Disclosure of Invention
The invention aims to solve the technical problem of providing a curved surface reconstruction method for scattered point cloud data of a heart, which utilizes the technologies of computer graphics, image processing and the like to carry out curved surface reconstruction on a large amount of scattered point cloud data of the heart to obtain a three-dimensional visual model, can intuitively and accurately display the spatial position relation of the heart, and can realize arbitrary display, scaling, rotation and cutting in a three-dimensional space.
The invention is realized by the following steps: a method for reconstructing a curved surface of scattered point cloud data of a heart comprises the following steps of reconstructing a three-dimensional surface model of the heart from scattered point cloud data of the heart, which is acquired by a catheter, contains noise and a large number of external points and is sampled unevenly, and comprises the following steps:
step one, defining a data representation form of collected point cloud data, establishing an octree topological relation, and defining a space function;
secondly, creating a vector field, selecting a local adjacent plane, and solving a vertex normal vector according to the average value of the surface normal vectors of the local adjacent plane;
step three, unifying vertex normal vectors;
and step four, solving the Poisson equation to obtain an indication function, and extracting an isosurface by adopting an MC algorithm according to the indication function and the gradient thereof to complete the reconstruction of the heart three-dimensional surface model.
In the first step, an octree topological relation is established, and a spatial function is defined, and the specific steps are as follows:
1) establishing an octree topological relation for the scattered point cloud data of the heart collected by the catheter;
2) setting the maximum recursion depth of the octree, and adding all scattered point cloud data of the heart into the octree;
3) defining a node function for each node of the constructed octree, and using bounding box filteringnDimensional convolution to select spatial function
：
Wherein,is the coordinate of any point in the point cloud data point set,
is a function of the bounding box filtering,
，nis the order of the filter, here taken to be 3.
The maximum recursion depth of the octree is 8.
Creating a vector field in the second step, selecting a local adjacent plane, and solving a vertex normal vector, wherein the method specifically comprises the following steps:
1) acquiring the nearest K adjacent points by using a KNN algorithm, and calculating a local adjacent plane fitting the point and an adjacent point set thereof by using a least square approximation method;
2) determining the normal vector of the surface of the adjacent surface：
And calculating the vertex normal vector from the mean value
：
Wherein,
、
is an integer which is a function of the number,
for marking points in a point cloud data point set,
is any point, referred to herein as a vertex,
for marking the firstPoints in the K-neighborhood of the point cloud data,
is any point in the K-neighborhood region,is that
Point and point
The local contiguous plane of points is formed,
is composed of
The distance of the point from the origin of coordinates,
is a Gaussian weight function, which is calculated by
Projection point of point on three-dimensional coordinate plane
Is a parameter, symbol
To represent
Dot sum
Distance between points, function
Is expressed in a constraint
The normal vector is obtained byThe minimized coordinates of (a);
3) defining an approximation indicative of a function vector field
：
Wherein,
is a point-cloud data point set,
is the K neighborhood of any point in the point set,
is any node
K in the neighborhood of K is
The number of the 8 nodes of (a),
is a linear coefficient of the linear coefficient,is thatThe node function of the point is then determined,is that
Vertex normal vectors for points.
In the third step, vertex normal vectors are uniformized, and the method specifically comprises the following steps:
1) selecting two adjacent points
、，、
Is composed of
、The dot product of the normal vectors of two adjacent points is calculated
Converting the problem of the consistency of the normal vector into the best solution of the graphA small spanning tree problem;
2) calculating the cost of an edge between two adjacent pointsAnd construct an undirected graph, cost
Calculated from the following formula:
wherein,
is a line segmentA midpoint of
And
is composed of
Distance on straight line
Two points of a unit distance of a point,
and
is composed of
Distance on straight lineTwo points of a unit distance of a point,
is composed of
On line segment
Projection of (2);
3) the normalization of the normal vector is biased towards the direction of propagation in the tangential direction.
Solving the Poisson equation in the fourth step, and extracting the isosurface, wherein the method comprises the following specific steps:
1) solving the solution of the indicating function equation by adopting a GS matrix iteration mode;
2) selecting a proper threshold, extracting an isosurface by using an MC algorithm, and judging whether each triangular surface is an ambiguous surface point by point: if the points greater than the threshold value and less than the threshold value are respectively positioned at the two ends of the diagonal line, ambiguity exists;
3) and eliminating the ambiguous surface, splicing the triangular surface patches, and completing the reconstruction of the heart three-dimensional surface model.
And when solving the indicating function equation by adopting a GS matrix iteration mode, the iteration number is set to be 6-8.
The method for reconstructing the curved surface of the scattered point cloud data of the heart utilizes the technologies of computer graphics, image processing and the like to carry out curved surface reconstruction on a large amount of scattered point cloud data of the heart to obtain a three-dimensional visual model, so that not only can the spatial position relation of the heart be visually and accurately displayed, but also the random display, the scaling, the rotation and the cutting in a three-dimensional space can be realized; the method can greatly reduce the influence of noise and external points in the process of reconstructing the curved surface of the heart, improves the reconstruction speed and precision, and has good robustness. After the heart is subjected to three-dimensional reconstruction, dynamic simulation and analysis can be performed, detection of heart diseases is facilitated, computer diagnosis and analysis of diseases such as arrhythmia and accurate simulation of surgery are realized, and safety and reliability of disease diagnosis are improved.
Drawings
FIG. 1 is a block diagram of a flow chart of a method for reconstructing a curved surface of scattered point cloud data of a heart according to the present invention;
FIG. 2 is a schematic diagram of establishing an octree topology relationship;
FIG. 3 is a schematic diagram of K neighbors solved using KNN algorithm;
FIG. 4 is a schematic diagram of the vertex normal vector solution principle;
FIG. 5 is a schematic diagram of the principle of normal vector uniformization;
FIG. 6 is a block diagram of a process for disambiguating an ambiguous surface when extracting an iso-surface using the MC algorithm.
In the figure:
is a certain point in the point cloud data point set,
Points in the K neighborhood of the current point,
Other points in the point cloud data point set.
Detailed Description
The invention will be further illustrated with reference to the following specific examples. It should be understood that these examples are for illustrative purposes only and are not intended to limit the scope of the present invention. Furthermore, it should be understood that various changes and modifications can be made by one skilled in the art after reading the description of the present invention, and equivalents fall within the scope of the invention defined by the appended claims.
Example 1
As shown in fig. 1, a method for reconstructing a curved surface of cardiac scattered point cloud data, which reconstructs a cardiac three-dimensional surface model from cardiac scattered point cloud data which is acquired by a catheter, contains noise and a large number of external points and is sampled unevenly, includes the following steps:
step one, defining a data representation form of collected point cloud data, establishing an octree topological relation, and defining a space function, so that the searching speed of points can be improved, and the solving process of the space function can be accelerated;
in this step, an octree topological relation is established, and a spatial function is defined, which specifically comprises the following steps:
1) establishing octree topological relation for scattered point cloud data of heart collected by catheter, reading all point cloud data, recording maximum and minimum coordinates of points in x, y and z directions,
、
、
、
、
、
and numbering for each point. Establishing a maximum bounding box, a so-called bag, of point cloud dataThe bounding box points the three-dimensional space occupied by the minimum external cube of the cloud data, and the maximum bounding box contains all the point cloud data which is used as the root node of the tree;
2) as shown in fig. 2, the maximum recursion depth of the octree is set, the maximum recursion depth selected in this embodiment is 8, and all the scattered point cloud data of the heart are added to the octree;
equally dividing the points in three directions to generate 8 sub-bounding boxes respectively, and adding the point sets into each sub-bounding box in sequence: if the number of the points is smaller than the size of the bounding box, recording the number of the point in the node, and stopping subdivision; if the number is larger than the size of the bounding box, the bounding box is recursively subdivided. The generated 8 sub-bounding boxes are represented by 0-7 until all data is added to the octree.
Each child node contains the coordinates and dimensions of the center of the bounding box, as well as the number of each point in the octree.
Let the number of point cloud data in a certain bounding box be
Then the dimensions of the sub-bounding box are as follows:
the spatial coordinates of the center of each bounding box are:
set a certain point
Has the coordinates of
Then its number in the bounding box is
Wherein:
and when the cardiac scattered point cloud data collected by the catheter is added into the octree according to the set maximum recursion depth, the octree is established.
3) For each node of the constructed octree
（
Representing octree) defines node functions：
Wherein,is a node
Is located in the center of the (c),
，
is any point cloud data point in a certain node,is the size of the node.
To make a vector field
Can be expressed as a linear sum of node functions, using bounding box filteringnDimensional convolution to select spatial function
：
Wherein,
is the coordinate of any point in the point cloud data point set,
is a function of the bounding box filtering,
，nis the order of the filter, here taken to be 3.
Secondly, creating a vector field, selecting a local adjacent plane, and solving a vertex normal vector according to the average value of the surface normal vectors of the local adjacent plane;
in the step, a vector field is created, local adjacent planes are selected, and vertex normal vectors are solved, and the method specifically comprises the following steps:
1) the sampling surface is assumed to be smooth everywhere so that the local neighborhood of any point can be well fitted with a plane. As shown in FIGS. 3 and 4, any point in the point cloud data point set is used herein
Expressing, called vertex, using KNN (K-Nearest Neighbor) algorithm to obtain its Nearest K neighbors, using least square approximation method to calculate local adjacent plane fitting said Neighbor and its neighbors set
Surface normal vector of the adjacent surfaceCalculated from the following formula:
wherein,
is that
Point K is close to any of the neighbors,is composed of
Distance of a point to the origin of coordinates;
2) in order to make the neighboring points of the current point have relatively large influence on the estimation result of the normal vector of the point, thereby obtaining a smoother estimation result, a Gaussian (Gaussian) weight function can be given when calculating the normal vector
，
At each point
To its projected point on a three-dimensional coordinate plane
Is a parameter, a more accurate surface normal vector of the abutment surface can be obtained by the following formula
Estimating a formula:
and by the surface normal vector of the abutment surface
Computing vertex normal vector from the mean of
：
Wherein,
、
is an integer which is a function of the number,
for marking points in a point cloud data point set,
is any point, referred to herein as a vertex,for marking the first
Points in the K-neighborhood of the point cloud data,
is any point in the K-neighborhood region,
is that
Point and pointThe local contiguous plane of points is formed,
is composed of
The distance of the point from the origin of coordinates,
is a Gaussian weight function, which is calculated byProjection point of point on three-dimensional coordinate plane
Is a parameter, symbol
To representDot sum
Distance between points, functionIs expressed in a constraint
The normal vector is obtained by
The minimized coordinates of (a);
3) defining an approximation indicative of a function vector field
：
Wherein,is a point-cloud data point set,
is the K neighborhood of any point in the point set,
is any node
K in the neighborhood of K is
The number of the 8 nodes of (a),
is a linear coefficient of the linear coefficient,
is that
The node function of the point is then determined,
is that
Vertex normal vectors for points.
Step three, the vertex normal vectors are uniformized to greatly shorten the time of curved surface reconstruction, the effects of noise reduction and smoothness can be achieved on the non-preprocessed scattered point cloud data, the vivid surface of the heart three-dimensional surface model can be quickly obtained, and the robustness is high; the vertex normal vector is consistent by firstly determining the direction of the normal vector of the current vertex, and then determining the normal vectors of the adjacent vertexes by taking the direction as a standard. This requires that the sampling point set must be dense enough, and assuming that the sampling surface is smooth everywhere, the unit normal vector variation of two adjacent points is small and can be considered to be close to parallel.
In the step, vertex normal vectors are uniformized, and the method specifically comprises the following steps:
1) selecting two adjacent points
、，
、
Is composed of
、The dot product of the normal vectors of two adjacent points is calculated
The dot product of the two-point normal vector is approximately 1, i.e.
And the sign of the dot product is positive, namely the included angle of the adjacent normal vectors is an acute angle; if the dot product is negative, it indicates that the normal vector of the point needs to be inverted, and the normal vector of a certain point can be set as
If, ifThen, thenOtherwise
Converting the vertex normal vector consistency problem into a minimum spanning tree problem of a solution graph;
2) as shown in fig. 5, the cost of an edge between two adjacent points is calculated
Taking a point in contact with the bounding box to which the point belongs as a starting point, constructing an undirected graph, and calculating the cost of an edge in the undirected graph by considering whether the adjacent points are shielded or not and passing through
And
the gradient value of the algebraic sphere at two points and the included angle of the normal vectors of the two points are determined, at the moment,
and
cost of edge between
Calculated from the following formula:
wherein,
is a line segmentA midpoint of
And
is composed ofDistance on straight line
Two points of a unit distance of a point,and
is composed of
Distance on straight line
Two points of a unit distance of a point,
is composed of
On line segment
Projection of (2);
3) the direction constraint of the normal vector is considered, namely the normal vector is biased to the direction of propagation along the tangential direction as much as possible, and the reconstruction time is shortened by only considering the angle difference of the distance and the normal vector.
And step four, solving a Poisson equation, extracting an isosurface by adopting an MC (Marching cubes) algorithm according to the indication function and the gradient thereof, and completing the reconstruction of the heart three-dimensional surface model.
Solving the Poisson equation in the step, and extracting the isosurface, wherein the method comprises the following specific steps:
1) solving the solution of the indicating function equation by adopting a GS (Gauss-Seidel) matrix iteration mode, wherein the iteration times are generally set to be 6-8 times;
obtaining a vector field
Then, solving the indicator function
Make it
And obtaining by applying a divergence operator:
wherein,
is Laplacian (Laplacian), and thus the Poisson equation is obtained
And converting the problem of surface reconstruction into a Poisson equation to solve.
Known vector
By solving a function
To simplify the problem:
at this time, the solution
Is equivalent to make
Projection of the Laplace equation into a spatial functionEach value of (A) is as close as possible to
The matrix form is expressed as:
thereby, it is possible to obtain:。
solving for vectors
Defining a matrix
：
Make it
Return to
Laplace dot product of each value, thus converting the solution of Poisson's equation to a sparsely symmetric system
In (1).
For all
Solving for
Is equal to solving its minimum bungalow approximationAnd a matrix
The solution of the equation is sparse, so the solution of the equation is solved by adopting a Gauss-Seidel matrix iteration mode. Gauss-Seidel iteration is a common method for solving the Poisson equation, because the Poisson equation is a real symmetric positive definite diagonal matrix, the Gauss-Seidel iteration is convergent, and the limit value is the accurate solution of the calculated Poisson equation.
2) As shown in fig. 6, a proper threshold is selected, an iso-surface is extracted by using the MC algorithm, and it is necessary to determine point-by-point whether each triangular surface is an ambiguous surface: if the points greater than the threshold value and less than the threshold value are respectively positioned at the two ends of the diagonal line, ambiguity exists;
the hypothesis model is
To obtain a reconstructed surface
First, a threshold is selected, which is selected such that the extracted iso-surface most closely approximates the true shape of the sampled point cloud data. For any point in the point cloud data
Is provided with
Since the indicator function is a piecewise constant function, direct calculation of the vector field results in the vector field having infinite values at the surface edges, and therefore smoothing filtering is used firstConvolution indicates the function and then considers the vector field of the smoothed function.
There is an integral relationship between the model surface sampling points and the indicator function, as shown in the following equation:
wherein,
，
is that
The normal vector of (2).
Sampling point set of known point cloud dataEach sampling point is defined by a position
Sum normal vector
Composition, a sampling point set is divided into small area blocksThen for each
The integrated sum is approximated as shown in the following equation:
at this point, the mean of the indicator function can be used to estimate the iso-surface as shown below:
and finally, extracting the isosurface by adopting a Marching Cubes algorithm.
3) And eliminating the ambiguous surface, splicing the triangular surface patches, establishing a smooth and vivid model surface, and completing the reconstruction of the heart three-dimensional surface model.
The preferred embodiment of the method for reconstructing the curved surface of the scattered point cloud data of the heart is explained in detail as follows: application of the inventionA large amount of scattered point cloud data of the inner wall of a heart cavity of a human body, which are acquired by a catheter, are subjected to a curved surface reconstruction experiment on a computer with Intel 2.2GHz and memory 2GB and a Visual Studio 2005 platform under a Windows operating system. The number of points in the point cloud data point set is 78376, the point cloud data point set contains a large number of external points and redundant information, the true surface point set only occupies less than one, and the maximum and minimum values in the x, y and z directions are respectively:
= 74.0390， 
= -33.6701， 
= 24.5505， = -49.4262， 
= 31.5073， = -50.9179。
in the experiment, the maximum recursion depth of the octree is 6-8;
it is the focus of the present invention to unify the normal vectors. Although it is not very difficult to select the normal vectors of the local adjacent planes, keeping the normal vectors of all the points consistent is an important control factor of the time complexity of the reconstruction algorithm. The method for the consistency of the normal vectors can greatly shorten the time of surface reconstruction, can achieve the effects of noise reduction and smoothness on the non-preprocessed scattered point cloud data, has stronger robustness, can quickly acquire the vivid surface of the model without preprocessing the original point cloud data in advance, and saves the operation steps and the time possibly spent. The method has quite good processing effect on the scattered point set containing a large amount of noise and external points and the scattered point set with sharp edge characteristics and non-uniform sampling density in the experiment, and the robustness of the experiment result comparison algorithm is shown in table 1.
TABLE 1 comparison of robustness of existing reconstruction algorithms with the reconstruction algorithm of the present invention
| Method of producing a composite material | Noise(s) | Exterior point | Sharp feature | Non-uniform sampling |
| Existing reconstruction algorithms | √ | √ | ||
| Reconstruction algorithm of the invention | √ | √ | √ | √ |
The difficulty of the invention is to solve the solution of the indicating function, select a proper threshold value and extract the isosurface by the MC algorithm, and at the moment, whether each triangular surface is an ambiguous surface needs to be judged point by point. Comparing each vertex indication function value with the isosurface threshold value, marking the point which is larger than the threshold value as 1, and marking the point which is smaller than the threshold value as 0, in the same triangular mesh, if the points which are 1 and 0 are respectively positioned at the two ends of the diagonal line, two possible connection modes exist, and therefore ambiguity exists.
The flow of eliminating ambiguous surface is shown in FIG. 6:
to be parallel to
Voxel plane of axis
For example, set points
、
、
、
Connecting PQ and MN respectively as the intersection point of the isosurface and the ambiguous surface, the equation of the straight line PQ is as follows:
the equation for MN is:
simultaneously obtaining the following two equations:
order to
，
All of which are constants, the coordinates of the intersection O of PQ and MN can be obtained:
calculating the function value at the intersection point O
Comparison of
And determining the state of the intersection point O according to the relation with the isosurface threshold value, thereby judging the connection mode of the isoline and eliminating the ambiguity problem of curved surface reconstruction:
if it is
If the triangular mesh is judged to be 1 when the triangular mesh is larger than the isosurface threshold, the triangular surface and two vertexes marked as 1 are in the same area; on the contrary, if
If the triangular mesh is judged to be 0 when the isosurface threshold is less than the isosurface threshold, the triangular surface and two vertexes marked as 0 are in the same area. In this way, the connection mode of the isosurface can be uniquely determined, and the ambiguity problem is eliminated.
And splicing the triangular patches to obtain a three-dimensional surface model after the curved surface is reconstructed. The triangle patches represent the surface of the heart model, the small triangle patches in each mesh are represented by three vertexes of the triangle and the connection sequence of the vertexes, and the triangle patches are fitted and curved, so that the purpose of converting abstract heart point cloud data into a visual three-dimensional model is achieved, and a smooth three-dimensional curved surface of the heart surface is formed.
The experimental results are as follows:
in an experiment, when the recursion depth of the octree is 6-8, all data nodes can be contained in the octree, and when the depth is too large, although the precision of a point set is increased, the calculation process of a solution in an iteration process is greatly increased; if the depth is too small, some nodes cannot be reflected, and the reconstructed model is rough and not smooth. Table 2 shows the parameter comparison in the point cloud data surface reconstruction results when the octree recursion depths are 6, 7, and 8.
TABLE 2 comparison of parameters of point cloud data surface reconstruction at different octree depths
| Octree recursion depth | 6 | 7 | 8 |
| Number of vertices | 1585 | 1693 | 2346 |
| Number of triangular patches | 3166 | 4382 | 4688 |
| Reconstruction time (seconds) | 7.53 | 12.31 | 18.24 |
As can be seen from the experimental results, as the octree depth increases, although the reconstruction time increases, the number of vertexes and triangular facets in the reconstruction result is more and more, which means that the reconstruction surface is more smooth and vivid.
In order to more clearly and accurately see the mesh structure of the reconstructed curved surface, when the octree recursion depth is 8, comparison experiments are respectively carried out under 6 and 8 different iterations for solving the poisson equation Gauss-Seidel, and table 3 shows the comparison of experiment parameters.
TABLE 3 comparison of parameters at different Gauss-Seidel iterations
| Number of Gauss-Seidel iterations | 6 | 8 |
| Number of vertices | 2346 | 2437 |
| Number of triangular patches | 4688 | 4870 |
| Reconstruction time (seconds) | 18.24 | 18.24 |
According to the experimental results, the difference of Gauss-Seidel iteration times does not increase the reconstruction time, but can obtain a more accurate three-dimensional model of the heart curved surface.
According to the method for reconstructing the curved surface of the scattered point cloud data of the heart, the reconstructed three-dimensional model of the heart is more vivid, the three-dimensional effect of a large amount of scattered point cloud data acquired by a doctor on the inner wall of the heart by using a catheter is well reflected, the spatial position relation of the heart can be intuitively and accurately displayed, and random display, scaling, rotation and cutting in a three-dimensional space can be realized. The improved Poisson reconstruction theory is utilized to carry out the curved surface reconstruction of the heart, the influence of noise and external points can be greatly reduced, the reconstruction speed and precision are improved, and the robustness is good. The method can quickly acquire the realistic surface of the model without preprocessing the original point cloud data in advance, thereby saving operation steps and time which may be spent. Meanwhile, the method can also be applied to three-dimensional reconstruction of other point cloud data. The research solves the problem of reconstructing the three-dimensional curved surface of the heart surface by using point cloud data, is beneficial to more effectively detecting and diagnosing diseases in the aspect of heart by doctors, and improves the accuracy and safety of medical diagnosis.
Claims (7)
1. A method for reconstructing a curved surface of scattered point cloud data of a heart is characterized in that a three-dimensional surface model of the heart is reconstructed from scattered point cloud data of the heart, which contains noise and a large number of external points and is sampled unevenly and is acquired by a catheter, and the method comprises the following steps:
step one, defining a data representation form of collected point cloud data, establishing an octree topological relation, and defining a space function;
secondly, creating a vector field, selecting a local adjacent plane, and solving a vertex normal vector according to the average value of the surface normal vectors of the local adjacent plane;
step three, unifying vertex normal vectors;
and step four, solving the Poisson equation to obtain an indication function, and extracting an isosurface by adopting an MC algorithm according to the indication function and the gradient thereof to complete the reconstruction of the heart three-dimensional surface model.
2. The method for curved surface reconstruction of cardiac scattered point cloud data according to claim 1, wherein the first step of establishing an octree topological relation and defining a spatial function comprises the following specific steps:
establishing an octree topological relation for the scattered point cloud data of the heart collected by the catheter;
setting the maximum recursion depth of the octree, and adding all scattered point cloud data of the heart into the octree;
defining a node function for each node of the constructed octree, and using bounding box filteringnDimensional convolution to select spatial function
：
Wherein,is the coordinate of any point in the point cloud data point set,
is a function of the bounding box filtering,
，nis the order of the filter, here taken to be 3.
3. The method of curved reconstruction of cardiac scattered point cloud data as recited in claim 2, wherein: the maximum recursion depth of the octree is 8.
4. The method for curved surface reconstruction of cardiac scattered point cloud data according to claim 1, wherein the second step is to create a vector field, select local adjacent planes, and calculate a vertex normal vector, and the specific steps are as follows:
acquiring the nearest K adjacent points by using a KNN algorithm, and calculating a local adjacent plane fitting the point and an adjacent point set thereof by using a least square approximation method;
determining the normal vector of the surface of the adjacent surface
：
And calculating the vertex normal vector from the mean value：
Wherein,
、
is an integer which is a function of the number,
for marking points in a point cloud data point set,
is any point, referred to herein as a vertex,
for marking the firstPoints in the K-neighborhood of the point cloud data,
is any point in the K-neighborhood region,
is that
Point and point
The local contiguous plane of points is formed,
is composed of
The distance of the point from the origin of coordinates,
is a Gaussian weight function, which is calculated by
Projection point of point on three-dimensional coordinate plane
Is a parameter, symbol
To represent
Dot sum
Distance between points, function
Is expressed in a constraint
The normal vector is obtained by
The minimized coordinates of (a);
defining an approximation indicative of a function vector field
：
Wherein,
is a point-cloud data point set,
is the K neighborhood of any point in the point set,
is any node
K in the neighborhood of K is
The number of the 8 nodes of (a),
is a linear coefficient of the linear coefficient,
is that
The node function of the point is then determined,
is thatVertex normal vectors for points.
5. The method for curved surface reconstruction of cardiac scattered point cloud data according to claim 1, wherein the step three is to unify vertex normal vectors, and the specific steps are as follows:
selecting two adjacent points、
，
、
Is composed of
、
The dot product of the normal vectors of two adjacent points is calculated
Converting the problem of the consistency of the normal vector into the problem of the minimum spanning tree of the solution graph;
calculating the cost of an edge between two adjacent points
And construct an undirected graph, cost
Calculated from the following formula:
wherein,is a line segment
A midpoint of
And
is composed of
Distance on straight line
Two points of a unit distance of a point,
and
is composed of
Distance on straight line
Two points of a unit distance of a point,
is composed of
On line segmentProjection of (2);
the normalization of the normal vector is biased towards the direction of propagation in the tangential direction.
6. The method for curved surface reconstruction of the cardiac scattered point cloud data according to claim 1, wherein the solving of the poisson equation in the fourth step extracts an isosurface, and the method comprises the following specific steps:
solving the solution of the indicating function equation by adopting a GS matrix iteration mode;
selecting a proper threshold, extracting an isosurface by using an MC algorithm, and judging whether each triangular surface is an ambiguous surface point by point: if the points greater than the threshold value and less than the threshold value are respectively positioned at the two ends of the diagonal line, ambiguity exists;
and eliminating the ambiguous surface, splicing the triangular surface patches, and completing the reconstruction of the heart three-dimensional surface model.
7. The method of curved reconstruction of cardiac scattered point cloud data according to claim 6, wherein: and when solving the indicating function equation by adopting a GS matrix iteration mode, the iteration number is set to be 6-8.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN 201110264828 CN102999937A (en) | 2011-09-08 | 2011-09-08 | Curved planar reconstruction method for cardiac scattered-point cloud data |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN 201110264828 CN102999937A (en) | 2011-09-08 | 2011-09-08 | Curved planar reconstruction method for cardiac scattered-point cloud data |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN102999937A true CN102999937A (en) | 2013-03-27 |
Family
ID=47928466
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN 201110264828 Pending CN102999937A (en) | 2011-09-08 | 2011-09-08 | Curved planar reconstruction method for cardiac scattered-point cloud data |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN102999937A (en) |
Cited By (20)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103793939A (en) * | 2013-07-29 | 2014-05-14 | 北京正安融翰技术有限公司 | Local increasing type curved-surface reconstruction method of large-scale point cloud data |
| CN103839292A (en) * | 2014-03-06 | 2014-06-04 | 中国科学院自动化研究所 | Method for sampling contour surface and generating high-quality triangular mesh |
| CN104282039A (en) * | 2014-09-29 | 2015-01-14 | 樊晓莉 | Skeleton orthosis brace shaping method based on 3D scanning |
| CN104978440A (en) * | 2014-04-03 | 2015-10-14 | 上海联影医疗科技有限公司 | Heart model building method, heart model registration and heart multi-plane reconstruction method |
| CN105404898A (en) * | 2015-11-26 | 2016-03-16 | 福州华鹰重工机械有限公司 | Loose-type point cloud data segmentation method and device |
| CN105844713A (en) * | 2016-04-22 | 2016-08-10 | 中北大学 | Parametric surface reconstruction method based on manifold |
| CN106504332A (en) * | 2016-10-19 | 2017-03-15 | 未来科技（襄阳）有限公司 | The curve reestablishing method and device of three-dimensional point cloud |
| CN106558102A (en) * | 2016-10-10 | 2017-04-05 | 合肥阿巴赛信息科技有限公司 | A kind of three-dimension modeling method rebuild based on Screened Poisson |
| CN106780751A (en) * | 2017-01-19 | 2017-05-31 | 桂林电子科技大学 | Three-dimensional point cloud method for reconstructing based on improved shielding Poisson algorithm |
| CN107221028A (en) * | 2017-05-27 | 2017-09-29 | 电子科技大学 | A kind of geologic body closed surface three-dimensional rebuilding method based on seismic interpretation data |
| CN107545599A (en) * | 2017-08-21 | 2018-01-05 | 上海妙影医疗科技有限公司 | Method of surface reconstruction and computer equipment in kind, storage medium |
| CN108090404A (en) * | 2016-11-23 | 2018-05-29 | 法乐第(北京)网络科技有限公司 | Three dimensional image processing method |
| CN108510443A (en) * | 2018-03-30 | 2018-09-07 | 河北北方学院 | A kind of medical image rebuilds localization method offline |
| CN108605119A (en) * | 2015-08-03 | 2018-09-28 | M·M·赫菲达 | 2D to 3D video frame is converted |
| CN108805261A (en) * | 2017-04-28 | 2018-11-13 | 微软技术许可有限责任公司 | Convolutional neural networks based on Octree |
| CN108961413A (en) * | 2018-07-16 | 2018-12-07 | 清华大学深圳研究生院 | A kind of internal membrane of heart Incremental insertion dynamic surface method for reconstructing in cardiac mapping system |
| CN110375668A (en) * | 2019-07-08 | 2019-10-25 | 西北农林科技大学 | Loess Surface mima type microrelief Surface Reconstruction based on point cloud data |
| CN111402393A (en) * | 2019-12-06 | 2020-07-10 | 温州大学 | Method for generating parameter curved surface simulation point cloud |
| CN113160186A (en) * | 2021-04-27 | 2021-07-23 | 青岛海信医疗设备股份有限公司 | Lung lobe segmentation method and related device |
| CN113838069A (en) * | 2021-09-27 | 2021-12-24 | 上海景吾智能科技有限公司 | Point cloud segmentation method and system based on flatness constraint |
-
2011
- 2011-09-08 CN CN 201110264828 patent/CN102999937A/en active Pending
Cited By (29)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103793939A (en) * | 2013-07-29 | 2014-05-14 | 北京正安融翰技术有限公司 | Local increasing type curved-surface reconstruction method of large-scale point cloud data |
| CN103839292A (en) * | 2014-03-06 | 2014-06-04 | 中国科学院自动化研究所 | Method for sampling contour surface and generating high-quality triangular mesh |
| CN103839292B (en) * | 2014-03-06 | 2017-05-10 | 中国科学院自动化研究所 | Method for sampling contour surface and generating high-quality triangular mesh |
| CN104978440A (en) * | 2014-04-03 | 2015-10-14 | 上海联影医疗科技有限公司 | Heart model building method, heart model registration and heart multi-plane reconstruction method |
| CN104978440B (en) * | 2014-04-03 | 2020-02-07 | 上海联影医疗科技有限公司 | Method for establishing and registering heart model and multi-plane reconstruction |
| CN104282039A (en) * | 2014-09-29 | 2015-01-14 | 樊晓莉 | Skeleton orthosis brace shaping method based on 3D scanning |
| CN108605119A (en) * | 2015-08-03 | 2018-09-28 | M·M·赫菲达 | 2D to 3D video frame is converted |
| CN108605119B (en) * | 2015-08-03 | 2020-12-01 | M·M·赫菲达 | 2D to 3D video frame conversion |
| CN105404898A (en) * | 2015-11-26 | 2016-03-16 | 福州华鹰重工机械有限公司 | Loose-type point cloud data segmentation method and device |
| CN105404898B (en) * | 2015-11-26 | 2018-11-06 | 福州华鹰重工机械有限公司 | A kind of loose type point cloud data segmentation method and equipment |
| CN105844713A (en) * | 2016-04-22 | 2016-08-10 | 中北大学 | Parametric surface reconstruction method based on manifold |
| CN106558102A (en) * | 2016-10-10 | 2017-04-05 | 合肥阿巴赛信息科技有限公司 | A kind of three-dimension modeling method rebuild based on Screened Poisson |
| CN106504332A (en) * | 2016-10-19 | 2017-03-15 | 未来科技（襄阳）有限公司 | The curve reestablishing method and device of three-dimensional point cloud |
| CN108090404A (en) * | 2016-11-23 | 2018-05-29 | 法乐第(北京)网络科技有限公司 | Three dimensional image processing method |
| CN106780751A (en) * | 2017-01-19 | 2017-05-31 | 桂林电子科技大学 | Three-dimensional point cloud method for reconstructing based on improved shielding Poisson algorithm |
| CN108805261A (en) * | 2017-04-28 | 2018-11-13 | 微软技术许可有限责任公司 | Convolutional neural networks based on Octree |
| CN108805261B (en) * | 2017-04-28 | 2021-11-12 | 微软技术许可有限责任公司 | Convolutional neural network based on octree |
| US11704537B2 (en) | 2017-04-28 | 2023-07-18 | Microsoft Technology Licensing, Llc | Octree-based convolutional neural network |
| CN107221028B (en) * | 2017-05-27 | 2019-06-14 | 电子科技大学 | A kind of geologic body closed surface three-dimensional rebuilding method based on seismic interpretation data |
| CN107221028A (en) * | 2017-05-27 | 2017-09-29 | 电子科技大学 | A kind of geologic body closed surface three-dimensional rebuilding method based on seismic interpretation data |
| CN107545599A (en) * | 2017-08-21 | 2018-01-05 | 上海妙影医疗科技有限公司 | Method of surface reconstruction and computer equipment in kind, storage medium |
| CN108510443A (en) * | 2018-03-30 | 2018-09-07 | 河北北方学院 | A kind of medical image rebuilds localization method offline |
| CN108961413A (en) * | 2018-07-16 | 2018-12-07 | 清华大学深圳研究生院 | A kind of internal membrane of heart Incremental insertion dynamic surface method for reconstructing in cardiac mapping system |
| CN110375668A (en) * | 2019-07-08 | 2019-10-25 | 西北农林科技大学 | Loess Surface mima type microrelief Surface Reconstruction based on point cloud data |
| CN111402393A (en) * | 2019-12-06 | 2020-07-10 | 温州大学 | Method for generating parameter curved surface simulation point cloud |
| CN113160186A (en) * | 2021-04-27 | 2021-07-23 | 青岛海信医疗设备股份有限公司 | Lung lobe segmentation method and related device |
| CN113160186B (en) * | 2021-04-27 | 2022-10-25 | 青岛海信医疗设备股份有限公司 | Lung lobe segmentation method and related device |
| CN113838069A (en) * | 2021-09-27 | 2021-12-24 | 上海景吾智能科技有限公司 | Point cloud segmentation method and system based on flatness constraint |
| CN113838069B (en) * | 2021-09-27 | 2024-02-06 | 杭州景吾智能科技有限公司 | Point cloud segmentation method and system based on flatness constraint |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN102999937A (en) | Curved planar reconstruction method for cardiac scattered-point cloud data | |
| Saha et al. | A survey on skeletonization algorithms and their applications | |
| Jones et al. | 3D distance fields: A survey of techniques and applications | |
| Huang et al. | Manifoldplus: A robust and scalable watertight manifold surface generation method for triangle soups | |
| Gagvani et al. | Parameter-controlled volume thinning | |
| Bouix et al. | Divergence-based medial surfaces | |
| Shen et al. | Skeleton pruning as trade-off between skeleton simplicity and reconstruction error | |
| CN111696210A (en) | Point cloud reconstruction method and system based on three-dimensional point cloud data characteristic lightweight | |
| CN111382777B (en) | Extracting feature trees from a grid | |
| Zhou et al. | Topology repair of solid models using skeletons | |
| Li et al. | On computing mapping of 3d objects: A survey | |
| CN105654483B (en) | The full-automatic method for registering of three-dimensional point cloud | |
| Kalvin et al. | Constructing topologically connected surfaces for the comprehensive analysis of 3-D medical structures | |
| Campen et al. | A Practical Guide to Polygon Mesh Repairing. | |
| Arakeri et al. | An effective and efficient approach to 3D reconstruction and quantification of brain tumor on magnetic resonance images | |
| CN111105476B (en) | Three-dimensional reconstruction method for CT image based on Marching Cubes | |
| Li et al. | On surface reconstruction: A priority driven approach | |
| Dvořák et al. | As-rigid-as-possible volume tracking for time-varying surfaces | |
| Sharma et al. | Signed distance based 3D surface reconstruction from unorganized planar cross-sections | |
| Patel et al. | 3D modeling and rendering of 2D medical image | |
| Khaleel et al. | 3D surface reconstruction of coronary artery trees for vessel locations’ detection | |
| Zou et al. | Area-preserving surface flattening using lie advection | |
| Azernikov et al. | Surface reconstruction of freeform objects based on hierarchical space decomposition | |
| Shi et al. | An integration method for scanned multi-view range images (MRIs) based on local weighted least squares (LWLS) surface fitting | |
| Sharma et al. | 3D Surface Reconstruction from Unorganized Sparse Cross Sections. |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| C06 | Publication | ||
| PB01 | Publication | ||
| C02 | Deemed withdrawal of patent application after publication (patent law 2001) | ||
| WD01 | Invention patent application deemed withdrawn after publication |
Application publication date: 20130327 |