A point cloud registration algorithm based on head features

Abstract

With the blossoming development and application of computer imaging technology in medicine, medical imaging technology has become one of the most concentrated and fastest-growing fields in medical technology. Magnetic resonance imaging (MRI) technology is a new type of medical imaging technology developed by using computer technology and image reconstruction technology based on the discovery of nuclear magnetic resonance phenomenon in the field of physics. Brainstorm software can be used to process head MRI images and generate point clouds on the head surface. OptiTrack is an optical system that realizes high-precision three-dimensional reconstruction of marker points. The main hardware is an intelligent camera equipped with a processor, with a built-in high-precision marker processing algorithm. The three-dimensional information of a single marker point can be used to create a rigid body and obtain six-degree-of-freedom data. The present invention proposes a point cloud registration algorithm based on head features. By using the point cloud features of the head, high-precision head point cloud registration can be completed when only a small number of head points are collected, which reduces the waiting time of patients and improves the accuracy of instrument treatment.

Description

Technical Field

[0001] The present invention relates to the field of computer vision, and in particular to a point cloud registration algorithm based on head features. Background Art

[0002] With the blossoming development and application of computer imaging technology in medicine, medical imaging technology has become one of the most concentrated and fastest-growing fields in medical technology. Magnetic resonance imaging (MRI) technology is a new type of medical imaging technology developed by using computer technology and image reconstruction technology based on the discovery of nuclear magnetic resonance phenomenon in the field of physics. It has the advantages of being safe and non-invasive, high resolution of brain and other soft tissues, and multi-directional and multi-parameter imaging, and is widely used in the field of brain science. Brainstorm software can be used to process head MRI images and generate point clouds of the head surface. OptiTrack is an optical system that realizes high-precision three-dimensional reconstruction of marker points. The main hardware is an intelligent camera equipped with a processor, with a built-in high-precision marker processing algorithm. The three-dimensional information of a single marker point can be used to create a rigid body and obtain six-degree-of-freedom data.

[0003] The iterative closest point algorithm (ICP) proposed by Besl and McKay in 1992 is currently the most widely used point cloud registration algorithm. However, the algorithm has the following disadvantages: (1) It has high requirements for the initial position and requires a good initial posture between the two point clouds. Otherwise, the iteration will not converge or fall into a local optimal solution, which will eventually lead to mismatching or non-convergence. (2) The traditional ICP algorithm does not eliminate data with large noise during the iteration process, so some points with large noise may have a greater impact on the final registration result. (3) When there is a subset relationship between the two sets of point clouds, the traditional ICP algorithm may cause the iteration result to not converge. Therefore, the traditional ICP algorithm cannot meet the needs of head point cloud registration. Summary of the invention

[0004] In order to overcome the shortcomings of the prior art, the present invention proposes a point cloud registration algorithm based on head features. By using the point cloud features of the head, high-precision head point cloud registration can be completed by collecting only a small number of head points, reducing the patient's waiting time and improving the accuracy of instrument treatment. Specifically, it includes:

[0005] Step 1: Get source point cloud and target point cloud

[0006] Use Brainstorm software to process head MRI images and generate source point clouds. Use OptiTrack cameras and a pen with a marker ball to collect 3D point cloud data in the world coordinate system and generate target point clouds.

[0007] Step 2: Roughly align the source point cloud and the target point cloud

[0008] Take the nasal root point NAS, the left ear preauricular point LPA and the right ear preauricular point RPA in the source point cloud and the target point cloud respectively, and perform rough registration to obtain the rough registration result of the source point cloud and the target point cloud: rotation matrix R c and the translation vector t c ;

[0009] Step 3: Adaptive threshold screening;

[0010] Step 4: Use the rough registration result and the screened threshold to perform fine registration to obtain a fine registration result;

[0011] Step 5: Optimize the precise registration result to minimize the cost;

[0012] Step 6: Calculate the final registration result.

[0013] 2. The point cloud registration algorithm based on head features as claimed in claim 1, characterized in that: the step 3 specifically comprises:

[0014] 3.1) Use the coarse registration result to project the target point cloud to obtain the coarse registration point cloud X c ;

[0015] 3.2) Calculate X c The mutual distances of NAS, LPA and RPA within the point cloud are denoted as d1, d2 and d3;

[0016] 3.3) Traverse the source point cloud and find X in the source point cloud c The set X of the closest points of NAS, LPA and RPA in the point cloud cNASs , X cLPAs and X cRPAs ;

[0017] 3.4) Set the threshold screening range and screening interval, start screening threshold: traverse X cNASs , X cLPAs and X cRPAs , let X cNASs Any point in x cNAS , X cLPAs Any point in x cLPA , X cRPAs Any point in x cRPA , calculate x respectively cNAS and x cLPA The distance is d c1 , calculate x cNAS and x cRPA The distance is d c2 , calculate x cLPA and xcRPA The distance is d c3 ; if |d1-d c1 |、|d2-d c2 | and |d3-d c3 | values ​​are all less than the threshold, then it is considered that x cNAS 、x cLPA and x cRPA A set of candidate point pairs can be formed, and all candidate point pairs under the threshold are grouped together:

[0018]

[0019] If the set X cset If the number of elements contained in is greater than a certain number, the threshold is recorded as an alternative threshold. After recording all the alternative thresholds, it is determined whether the number of alternative thresholds is greater than 0. If so, the smallest alternative threshold is selected as the screening result, and the threshold is recorded as thrd, and step 4 is executed; if not, the program ends and the registration fails.

[0020] 3. The head feature-based point cloud registration algorithm according to claim 2, wherein step 4 specifically comprises:

[0021] 4.1) Select candidate point pairs under the threshold thrd and calculate X cset Candidate point pairs and coarse registration point cloud in the set X c The rotation matrix R of NAS, LPA and RPA cf and the translation vector t cf , forming a set Rt:

[0022]

[0023] 4.2) Take any rotation matrix in the set Rt and translation vector pass Get candidate precise registration point cloud

[0024] 4.3) Traverse the candidate precision registration point cloud All points within the target point cloud X are calculated t If the closest distance is less than 3 mm, the point is considered to be an inner point after registration.

[0025] 4.4) Count all candidate precision registration point clouds The number of inliers is the number of candidate points, and the candidate fine registration point cloud with the largest number of inliers is selected. The corresponding rotation matrix R f and the translation vector t f As the result of fine registration.

[0026] 4. The head feature-based point cloud registration algorithm according to claim 3, wherein step 5 specifically comprises:

[0027] 5.1) Transform the point cloud X obtained by coarse registration c The rotation matrix R obtained by fine alignment f and the translation vector t f Carry on X f =R f X c +t f Transform to obtain the precisely aligned point cloud X f ;

[0028] 5.2) Traverse the precise registration point cloud X f and target point cloud X t , if the distance between a point in the precise registration point cloud and a point in the target point cloud is less than 3 mm, the two points are considered to be a matching point pair, and the set of matching point pairs between the two sets of point clouds is calculated;

[0029] 5.3) Use ICP point cloud registration algorithm to calculate the precise registration point cloud X f and target point cloud X t The rotation matrix R between s and the translation vector t s , making the cost minimum.

[0030] 5. The head feature-based point cloud registration algorithm according to claim 1, wherein step 6 specifically comprises: using the rotation matrix R c , R f and R s and the translation vector t c ,t f and t s Calculate the rotation matrix R of the final registration result n and the translation vector t n :

[0031] R n =(R s R f R c ) T

[0032]

[0033] Make the source point cloud X s and target point cloud X t Between X t =R n X s +t n , and obtain the final registration result.

[0034] The beneficial effects of the present invention are:

[0035] (1) To address the problem that the traditional point cloud registration algorithm fails when the overlap between two groups of point clouds is low and the characteristics of the overlapping parts are not obvious, the present invention combines head features to effectively realize the registration of two groups of head point clouds with low overlap and unclear overlapping features.

[0036] (2) Using the NAS, LPA and RPA points of the human head as constraints, accurate registration results can still be obtained when the number of points in the source point cloud is significantly different from the number of points in the target point cloud and the target point cloud is a subset of the source point cloud.

[0037] (3) Using the NAS, LPA and RPA points of the human head as constraints, the erroneous points in the registration point set can be effectively screened and eliminated, thereby improving the accuracy and robustness of the registration. BRIEF DESCRIPTION OF THE DRAWINGS

[0038] Figure 1 A flow chart of a point cloud registration algorithm based on head features provided in an embodiment of the present invention.

[0039] Figure 2 This is the source point cloud image obtained by processing the head MRI image using Brainstorm software.

[0040] Figure 3 A point cloud of a target obtained by collecting 3D points in world coordinates using an OptiTrack camera and a pen with a marker ball.

[0041] Figure 4 Schematic diagram of the result of point cloud registration using the traditional ICP algorithm.

[0042] Figure 5 This is a schematic diagram of the result of point cloud registration performed by the present invention. DETAILED DESCRIPTION

[0043] The specific implementation of the present invention is described below in conjunction with the accompanying drawings so that those skilled in the art can better understand the present invention. It should be noted that the following description is only used to explain the present invention and is not used to limit the present invention.

[0044] A point cloud registration algorithm based on head features includes the following steps:

[0045] Step 1: Get the source point cloud and the target point cloud.

[0046] 1.1) Use Brainstorm software to process the head MRI image and generate a point cloud of the head surface, which is called the source point cloud. s To express.

[0047]

[0048] in, For any point in the source point cloud, Brainstorm software can simultaneously calculate the rough positions of NAS, LPA, and RPA.

[0049] 1.2) Use the OptiTrack camera and the positioning pen with a marker ball to collect the head point cloud data in the world coordinate system. The first three points collected are the nasal root point (NAS), the left pre-auricular points (LPA) and the right pre-auricular points (RPA). The positions of the subsequent points are not required. This part of the point cloud is called the target point cloud. t To express.

[0050] Step 2: Roughly align the source point cloud and the target point cloud. Use Brainstorm software to calculate NAS, LPA, and RPA as the source point cloud. Take the three points NAS, LPA, and RPA of the source point cloud. Set the coordinates of these three points to be x sNAS 、x sLPA and x sRPA Take the first three points of the target point cloud, NAS, LPA and RPA, and set the coordinates of these three points to be x tNAS 、x tLPA and x tRPA Use the three points of the source point cloud and the target point cloud for coarse registration to obtain the rotation matrix R between the source point cloud and the target point cloud. c and the translation vector t c .

[0051] The steps of three-point registration are as follows: Use the source point cloud three points x sNAS 、x sLPA and x sRPA Calculate the direction vectors of the x-axis, y-axis, and z-axis of the source point cloud. The expressions are as follows:

[0052]

[0053] In the above formula, · represents vector dot product, × represents vector cross product, and the above obtained and are all 3×1 unit column vectors, so the rotation matrix R of the source point cloud relative to the world coordinate system is s It can be expressed as:

[0054]

[0055] Use the first three points x of the target point cloud tNAS 、x tLPAand x tRPA Calculate the direction vectors of the x-axis, y-axis, and z-axis of the target point cloud and Therefore, the rotation matrix R of the target point cloud relative to the world coordinate system t It can be expressed as:

[0056]

[0057] The final three-point registration result is calculated using the following expression:

[0058]

[0059] t c =x sNAS -R c x tNAS .

[0060] Step 3: Adaptive screening threshold.

[0061] 3.1) The target point cloud X t The rotation matrix R obtained by coarse registration c and the translation vector t c To transform:

[0062] X c =R c X t +t c ,

[0063] Among them, X c Represents the target point cloud X t The point cloud is obtained by coarse registration transformation of the points.

[0064] 3.2) Calculate X c The distance between NAS and LPA in the point cloud is recorded as d1; calculate X c The distance between NAS and RPA in the point cloud is recorded as d2; calculate X c The distance between LPA and RPA in the point cloud is recorded as d3.

[0065] 3.3) Traverse the source point cloud X s , calculate the source point cloud X s And coarse registration transform point cloud X c The 30 points with the closest NAS distance are used as the candidate point set for the actual NAS, and X cNASs To represent; calculate the source point cloud X s And coarse registration transform point cloud X c The 30 points closest to LPA in the LPA are used as the candidate point set for the actual LPA, and X cLPAs To represent; calculate the source point cloud X s And coarse registration transform point cloud Xc The 30 points closest to the RPA distance in the RPA are used as the candidate point set for the actual RPA, and X cRPAs To express.

[0066] 3.4) Set the threshold screening range and screening interval, and screen the minimum threshold that meets the conditions. cNASs , X cLPAs and X cRPAs , let X cNASs Any point in x cNAS , X cLPAs Any point in x cLPA , X cRPAs Any point in x cRPA , calculate x respectively cNAS and x cLPA The distance is d c1 , calculate x cNAS and x cRPA The distance is d c2 , calculate x cLPA and x cRPA The distance is d c3 If |d1-d c1 |、|d2-d c2 | and |d3-d c3 | values ​​are all less than the threshold, then xcNAS, xcLPA, and x cRPA A set of candidate point pairs can be formed. All candidate point pairs under this threshold can form a set:

[0067]

[0068] If the set X cset If the number of elements contained in is greater than 10, then the threshold is recorded as an alternative threshold. After recording all the alternative thresholds, determine whether the number of alternative thresholds is greater than 0. If so, select the smallest alternative threshold as the screening result, which is recorded as thrd, and execute step 4; if not, the program ends and the registration fails.

[0069] Step 4: Use the rough registration result and the screened threshold to perform fine registration.

[0070] 4.1) Select the candidate point pairs under the threshold thrd and use the three-point registration algorithm to calculate X cset Candidate point pairs in the set and coarse registration transformation point cloud X c The rotation matrix R of NAS, LPA and RPA cf and the translation vector t cf , these rotation matrices R cf and the translation vector t cf We can form a set Rt:

[0071]

[0072] 4.2) Take any set of rotation matrices in the set Rt and translation vector Using the rotation matrix and translation vector Calculate candidate fine registration results:

[0073]

[0074] in, is the coarse registration point cloud Xtc after the candidate rotation matrix and translation vector The candidate fine registration point cloud obtained by transformation.

[0075] 4.3) Traverse the candidate precision registration point cloud All points within the target point cloud X are calculated. t The closest distance to the point. If the closest distance is less than 3mm, the point is considered to be an inner point after alignment.

[0076] 4.4) Statistics The number of internal points included can form a set N:

[0077] N={n1,n2,…,n i}

[0078] Let Max(N)=n in the set N. j , then the rotation matrix corresponding to this value is and translation vector It is the result of precise alignment.

[0079]

[0080] Among them, R f is the rotation matrix obtained by fine registration, t f is the translation vector obtained by fine registration.

[0081] Step 5: Optimize the precise registration results.

[0082] 5.1) Transform the point cloud X obtained by coarse registration c The rotation matrix R obtained by fine alignment f and the translation vector t f To transform:

[0083] X f =R f X c +t f

[0084] Among them, Xf Represents point cloud X c The point cloud is obtained by precise registration and transformation of the points.

[0085] 5.2) Get the precise registration point cloud X f Any point inside Traverse the target point cloud X t , calculate the target point cloud X t and The closest point is If the point cloud is precisely aligned and the target point cloud If the distance is less than 3 mm, the two points are considered to be a matching point pair, and the precise registration point cloud X is traversed. f , we can get a set of matching point pairs.

[0086]

[0087] The set of matching points can be represented by P and P′:

[0088] P = {p1, p2, p3, ..., p i},

[0089] P′={p′1, p′2, p′3,…,p′ i},

[0090] Among them, p i ∈X f , p′ i ∈X t , and p i and p′ i are matching point pairs.

[0091] 5.3) For the pose estimation problem of two sets of point clouds with known matching points, the traditional ICP algorithm can obtain good results. Use the SVD decomposition method to solve the problem of f Transform to X t The rotation matrix R s and the translation vector t s , so that the cost J is minimized, the expression of the cost function J is as follows:

[0092]

[0093] Step 6: Calculate the final registration result. Use the rotation matrix R calculated in the previous step c , R f and R s and the translation vector t c ,t f and t s The rotation matrix R of the final registration result can be calculated nand the translation vector t n :

[0094] R n =(R s R f R c ) T ,

[0095]

[0096] Make the source point cloud X s and target point cloud X t The following relationship is satisfied:

[0097] X t =R n X s +t n .

[0098] The present invention is introduced in detail above, but the description of the specific implementation mode is only used to explain the method of the present invention and its core idea, so that technicians in this technical field can understand the present invention. However, it should be clear that the present invention is not limited to the scope of the specific implementation mode. For ordinary technicians in this technical field, as long as various changes are within the spirit and scope of the present invention defined and determined by the attached claims, these changes are obvious, and all inventions and creations using the concept of the present invention are protected.

Claims

1. A point cloud registration algorithm based on head features, characterized by: The following steps are involved: Step 1: Get source point cloud and target point cloud Use Brainstorm software to process head MRI images and generate source point clouds. Use OptiTrack cameras and a pen with a marker ball to collect 3D point cloud data in the world coordinate system and generate target point clouds. Step 2: Roughly align the source point cloud and the target point cloud Take the nasal root point NAS, the left ear preauricular point LPA and the right ear preauricular point RPA in the source point cloud and the target point cloud respectively, and perform rough registration to obtain the rough registration result of the source point cloud and the target point cloud: rotation matrix R c and the translation vector t c ; Step 3: Adaptive threshold screening; Step 4: Use the rough registration result and the screened threshold to perform fine registration to obtain a fine registration result; Step 5: Optimize the precise registration result to minimize the cost; Step 6: Calculate the final registration result; The step 3 specifically includes: 3.1) Use the coarse registration result to project the target point cloud to obtain the coarse registration point cloud X c ; 3.2) Calculate X c The mutual distances of NAS, LPA and RPA within the point cloud are denoted as d1, d2 and d3; 3.3) Traverse the source point cloud and find X in the source point cloud c The set X of the closest points of NAS, LPA and RPA in the point cloud cNASs , X cLPAs and X cRPAs ; 3.4) Set the threshold screening range and screening interval, start screening threshold: traverse X cNASs , X cLPA s and X cRPAs , let X cNASs Any point in x cNAS , X cLPAs Any point in x cLPA , X cRPAs Any point in x cRPA , calculate x respectively cNAS and x cLPA The distance is d c1 , calculate x cNAS and x cRPA The distance is d c2 , calculate x cLPA and x cRPA The distance is d c3 ; if |d1-d c1 |、|d2-d c2 | and |d3-d c3 | values ​​are all less than the threshold, then it is considered that x cNAS 、x cLPA and x cRPA A set of candidate point pairs can be formed, and all candidate point pairs under the threshold are grouped together: If the set X cset If the number of elements contained in is greater than a certain number, the threshold is recorded as an alternative threshold. After recording all the alternative thresholds, it is determined whether the number of alternative thresholds is greater than 0. If so, the smallest alternative threshold is selected as the screening result, and the threshold is recorded as thrd, and step 4 is executed; if not, the program ends and the registration fails; The step 5 specifically includes: 5.1) Transform the point cloud X obtained by coarse registration c The rotation matrix R obtained by fine alignment f and the translation vector t f Carry on X f =R f X c +t f Transform to obtain the precisely aligned point cloud X f ; 5.2) Traverse the precise registration point cloud X f and target point cloud X t , if the distance between a point in the precise registration point cloud and a point in the target point cloud is less than 3 mm, the two points are considered to be a matching point pair, and the set of matching point pairs between the two sets of point clouds is calculated; 5.3) Use ICP point cloud registration algorithm to calculate the precise registration point cloud X f and target point cloud X t The rotation matrix R between s and the translation vector t s , making the cost minimum.

2. The head feature-based point cloud registration algorithm according to claim 1, characterized in that: The step 4 specifically includes: 4.1) Select candidate point pairs under the threshold thrd and calculate X cset Candidate point pairs and coarse registration point cloud in the set X c The rotation matrix R of NAS, LPA and RPA cf and the translation vector t cf , forming a set Rt: 4.2) Take any rotation matrix in the set Rt and translation vectors pass Get candidate precise registration point cloud 4.3) Traverse the candidate precision registration point cloud All points within the target point cloud X are calculated t If the closest distance is less than 3 mm, the point is considered to be an inner point after registration. 4.4) Count all candidate precision registration point clouds The number of inliers is the number of candidate points, and the candidate fine registration point cloud with the largest number of inliers is selected. The corresponding rotation matrix R f and the translation vector t f As the result of fine registration.

3. The point cloud registration algorithm based on head features as claimed in claim 1, characterized in that: The step 6 specifically includes: using the rotation matrix R c , R f and R s and the translation vector t c ,t f and t s Calculate the rotation matrix R of the final registration result n and the translation vector t n : R n =(R s R f R c ) T Make the source point cloud X s and target point cloud X t Between X t =R n X s +t n , and obtain the final registration result.

Images