A cartilage repair method based on point cloud completion

Abstract

The present invention discloses a cartilage repair method based on point cloud completion, which constructs a cartilage point cloud repair dataset by preprocessing a cartilage model; adopts a point cloud completion Transformer, first downsamples part of the input point cloud to obtain the center point; then uses a lightweight DGCNN to extract local features around the center point; after adding position embedding to the local features, uses a transformer encoder-decoder architecture to predict the point proxy of the missing part; finally, uses a simple MLP and FoldingNet to complete the point cloud based on the predicted point proxy in a coarse-to-fine manner; the present invention breaks through the key technology of repairing and completing the missing part of the shape while effectively maintaining the fine structural information of the input point cloud shape; solves the small sample problem of joint bones by adopting the transfer learning method; and achieves accurate restoration by removing redundant completion points using precise extraction technology, thereby realizing the application of 3D shape completion in cartilage repair.

Description

Technical Field

[0001] The present invention relates to a cartilage repair method based on point cloud completion. Background Art

[0002] Currently, injuries caused by mechanical trauma, aging, or congenital diseases are increasing year by year. Because cartilage tissue has a very limited self-repair capacity, even a minor injury can cause joints to slowly wear out and degenerate during exercise. Existing clinical treatments can only relieve pain and slow further cartilage degeneration but are ineffective in promoting the regeneration of damaged cartilage. The biochemical and mechanical properties of newly formed cartilage differ significantly from those of normal cartilage, making it far from sufficient to meet the body's demands for exercise.

[0003] Therefore, cartilage repair surgery and cartilage transplantation are necessary. Cartilage repair surgery is highly operator-dependent, especially in restoring the healthy cartilage morphology, which requires a surgeon with extensive clinical experience. Artificial intelligence (AI) has the potential to revolutionize the field of radiology by improving the efficiency and accuracy of both interpretive and non-interpretive tasks, thereby resolving the question of "what kind of repair is needed" and assisting surgeons in better performing cartilage repair surgeries.

[0004] Predicting and restoring healthy cartilage shape aims to repair damaged cartilage and estimate a complete point cloud from the missing point cloud to obtain a point cloud of healthy cartilage shape, also known as point cloud completion. Currently, the vast majority of point cloud completion research focuses on objects, prioritizing learning general features of the genus / class rather than the local details of a specific object. Consequently, they may alter the positions of known points and suffer from categorical distortion, resulting in noise and loss of detailed geometry. 3D shape completion for biological skeletons should leverage relational structural information (such as geometric symmetry, regular arrangement, and surface smoothness) that existing methods cannot capture to reconstruct a reasonably fine-grained complete point cloud. Furthermore, because point clouds of 3D shapes are highly unstructured, models trained on 3D shapes typically produce point clouds with significant inhomogeneity. This inhomogeneity can lead to unwanted holes in the shape, compromising the integrity of our predictions, which is unacceptable in the task of bone completion.

[0005] In addition, due to the complex point distribution and the difficulty of training deep learning networks, the current network model can only produce smaller-scale point sets and cannot produce large-scale point sets, making it impossible to finely complete the precise morphology of biological bones. Summary of the Invention

[0006] The purpose of the present invention is to provide a cartilage repair method based on point cloud completion to solve the problem that existing cartilage defect repair methods are unable to finely complete biological bones with precise morphology.

[0007] To solve the above technical problems, the present invention provides a cartilage repair method based on point cloud completion, comprising:

[0008] S1: Preprocess the knee cartilage model to obtain residual defects of the knee cartilage;

[0009] S2: First, downsample the defect cloud to be repaired to obtain the center point; then extract the local features around the center point, and then add the local feature point proxy obtained by position embedding to the local features;

[0010] S3: Use the point cloud completion model to predict missing proxy points based on local feature point proxies, and then restore the missing point cloud corresponding to the missing proxy points from coarse to fine based on the predicted missing proxy points;

[0011] S4: splicing the missing point cloud obtained in step S3 with the femoral residual point cloud;

[0012] S5: Extract boundary points of missing point cloud;

[0013] S6: Use RANSAC to perform 3D circle fitting in the 3D point cloud space to obtain the spherical equation;

[0014] S7: Determine whether all points in the missing point cloud are within the sphere, delete the points that are not within the sphere, and leave the points within the sphere to form the point cloud of the area to be repaired. Combine the point cloud of the area to be repaired with the defective point cloud to obtain the healthy cartilage point cloud.

[0015] Furthermore, step S1 specifically includes:

[0016] S11: obtaining an STL model of the knee joint and femur, and converting the STL model of the knee joint and femur into a point cloud format to obtain a point cloud model of the knee joint and femur;

[0017] S12: using a straight-through filter to filter out all points in the Z-axis dimension that are not within a given value range, and segmenting the knee joint femur point cloud model;

[0018] S13: Randomly sample the segmented knee joint femoral point cloud model to obtain the original femoral point cloud set;

[0019] S14: The original point cloud set is cut into the upper surface point cloud of the femur according to the position and size of the cartilage defect in the actual surgery, and each original femoral point cloud is cut into several incomplete femoral point clouds to obtain a incomplete point cloud set.

[0020] Furthermore, step S2 specifically includes:

[0021] S21: Sampling the farthest point of the femoral residual defect cloud to obtain N fixed center points {q1,…,q i ,…,q N};

[0022] S22: Use lightweight DGCNN to extract features from the local area around the center point and obtain N local area features {F1′,…,F i ,…,F N ′}, where F i ′ corresponds to q i The local area features of the center point;

[0023] S23: Through MLP network Extract each local feature F i ′ is embedded in the position, and the point proxy is obtained after addition, that is,

[0024] Furthermore, the point cloud completion model is a Transformer encoder-decoder architecture, including:

[0025] An encoder, configured to encode the local feature point proxy obtained in step S2;

[0026] Geometric perceptron, used to capture geometric relationships in point clouds using the kNN model; given a query coordinate p Q , according to the bond coordinate p k The nearest key feature is queried and adjacent point proxies in space are spliced ​​together. Then, the local geometric structure is learned by performing a maximum pooling operation through feature aggregation using a linear layer. Finally, the geometric features and the semantic features of the self-attention mechanism are concatenated and mapped to the original dimension to form the output.

[0027] Query generator, which is used to first obtain the global feature f through the output feature of the encoder g , and predict the rough missing point cloud center coordinates {c1,…,c i ,…c M}, after concatenating the coordinates of the missing point cloud center point with the global features, the cartilage repair network model is used to generate the query features, i.e.

[0028]

[0029] Multi-scale point cloud generation module, translating query features into missing proxy points Use the Query Builder to replace the missing point cloud center coordinates {c1,…, c i ,…c M} as the local center of the missing point cloud, and use FoldingNet to recover the missing proxy points Detailed local shape centered on:

[0030]

[0031] Among them, C i The coordinates of the center point c of the missing point cloud are i The set of adjacent points centered at is the missing point cloud corresponding to the missing proxy point.

[0032] Furthermore, the cartilage repair network model construction method includes:

[0033] A cartilage repair network model is pre-trained using the existing point cloud completion dataset as the source model;

[0034] Create a new cartilage repair network model as the target model and migrate all the model designs and parameters of the source model to the target model;

[0035] The target model was trained using the knee cartilage dataset and fine-tuned based on the parameters of the source model to obtain a cartilage repair network model based on the fine-grained structure of cartilage.

[0036] Furthermore, the cartilage repair network model adopts the chamfer distance loss function as its loss function.

[0037] Furthermore, step S5 specifically includes:

[0038] S51: For any point P in the point cloud P i Set the search radius r, and record the neighborhood points within the search radius r as set N(P i ),Right now

[0039] N(P i )={P j |P j ∈P,‖P i -P j ‖ <r}

[0040] S52: Set the surface equation ax+by+cz=d(d≥0),a 2 +b 2 +c 2 =1, take P i The corresponding set N(P i ), calculate N(P i ) to the surface

[0041] d i =|ax+by+cz-d|

[0042] Solve N(P i ) The eigenvector corresponding to the minimum sum of the distances from all points to the surface is the normal vector of the point

[0043]

[0044] According to Pi Make the tangent plane Ψ of the point with the normal vector n;

[0045] S53: Set N(P i ) is projected onto the tangent plane Ψ, denoted as N(P i ′). In N(P i ′) select a point P j ',by is the u axis, n is the w axis, u×w is the v axis, and P i Construct a local coordinate system for the coordinate center, denoted as (P i ,u,v,w);

[0046] S54: Calculate the set N(P i Other points P in n ' to point P i The vector P i P n ′ with the coordinate axis u clockwise angle Ω=(θ1,θ2,θ3,…θ n ), make the difference between the adjacent angles and get the angle set Ω′=(θ′1,θ′2,θ′3···θ′ n-1 ),in

[0047] θ′ n-1 =θ n -θ n-1

[0048] Arrange the elements in the set Ω′ in descending order and find the largest angle θ′ max .

[0049] θ′ max =max(θ′1,θ′2,θ′3…θ′ n-1 )

[0050] When θ′ max When it is greater than the threshold, the point is changed to a boundary point.

[0051] Furthermore, step S6 specifically includes:

[0052] S61: Initialize the number of loops k = 0 and the local point set Inputs = NULL;

[0053] S62: Randomly select 3 points from the boundary point set D, and use these 3 points to solve the parameters of the circle, that is, the center (x0, y0, z0) and the radius R;

[0054] S63: Calculate the distance d from each boundary point to the center of the obtained circle. If dR≤ε, the point is counted as an internal point set; otherwise, it is considered an external point.

[0055] S64: Calculate the number of points inside the circle, record it as M, if M is greater than the threshold M min , then the estimation is considered successful and the process goes to step S65; otherwise, the process goes to step S26;

[0056] S65: recalculate the parameter model of the circle using the minimum multiplication method for all the parameters in the set Inputs to obtain the final result;

[0057] S66: k=k+1, if k>k max , then end; otherwise go to S62.

[0058] S67: Using the parametric model of a circle to obtain the equation of a sphere

[0059] (x-x0) 2 +(y-y0) 2 +(z-z0) 2 =R 2 .

[0060] The beneficial effects of the present invention are: the present invention breaks through the key technology of repairing and completing the missing parts of the shape while effectively maintaining the fine structural information of the input point cloud shape; at the same time, it also achieves accurate restoration by removing redundant completion points through the use of precise extraction technology, thereby realizing the application of 3D shape completion in knee cartilage repair. BRIEF DESCRIPTION OF THE DRAWINGS

[0061] The drawings described herein are used to provide a further understanding of the present application and constitute a part of the present application. The same reference numerals are used in these drawings to represent the same or similar parts. The exemplary embodiments of the present application and their descriptions are used to explain the present application and do not constitute an improper limitation on the present application. In the drawings:

[0062] Figure 1 The flowchart of one embodiment of the present invention is shown in FIG. DETAILED DESCRIPTION

[0063] like Figure 1 The cartilage repair method based on point cloud completion shown in the figure comprises:

[0064] S1: Preprocess the knee cartilage model to obtain residual defects of the knee cartilage;

[0065] S2: First, downsample the defect cloud to be repaired to obtain the center point; then extract the local features around the center point, and then add the local feature point proxy obtained by position embedding to the local features;

[0066] S3: Use the point cloud completion model to predict missing proxy points based on local feature point proxies, and then restore the missing point cloud corresponding to the missing proxy points from coarse to fine based on the predicted missing proxy points;

[0067] S4: splicing the missing point cloud obtained in step S3 with the femoral residual point cloud;

[0068] S5: Extract boundary points of missing point cloud;

[0069] S6: Use RANSAC to perform 3D circle fitting in the 3D point cloud space to obtain the spherical equation;

[0070] S7: Determine whether all points in the missing point cloud are within the sphere, delete the points that are not within the sphere, and leave the points within the sphere to form the point cloud of the area to be repaired. Combine the point cloud of the area to be repaired with the defective point cloud to obtain the healthy cartilage point cloud.

[0071] The above step S1 specifically includes:

[0072] S11: obtaining an STL model of the knee joint and femur, and converting the STL model of the knee joint and femur into a point cloud format to obtain a point cloud model of the knee joint and femur;

[0073] S12: using a straight-through filter to filter out all points in the Z-axis dimension that are not within a given value range, and segmenting the knee joint femur point cloud model;

[0074] S13: Randomly sample the segmented knee joint femoral point cloud model, limit the number of points of all femoral point clouds (e.g., limit it to 16384), and obtain the original femoral point cloud set;

[0075] S14: Cutting the original point cloud set into the upper surface point cloud of the femur according to the position and size of the cartilage defect in the actual surgery, cutting each original femoral point cloud into several (e.g., 8) incomplete femoral point clouds to obtain a incomplete point cloud set.

[0076] The above step S2 specifically includes:

[0077] S21: Sampling the farthest point of the femoral residual defect cloud to obtain N fixed center points {q1,…,q i ,…,q N};

[0078] S22: Use lightweight DGCNN to extract features from the local area around the center point and obtain N local area features {F1′,…,F i ,…,F N ′}, where F i ′ corresponds to q i The local area features of the center point;

[0079] S23: Through MLP network Extract each local feature F i ′ is embedded in the position, and the point proxy F is obtained after addition. i,Right now

[0080] The point cloud completion model in step S3 above is a Transformer encoder-decoder architecture, including:

[0081] An encoder, configured to encode the local feature point proxy obtained in step S2;

[0082] Geometric perceptron, used to capture geometric relationships in point clouds using the kNN model; given a query coordinate p Q , according to the bond coordinate p k The nearest key feature is queried and adjacent point proxies in space are spliced ​​together. Then, the local geometric structure is learned by performing a maximum pooling operation through feature aggregation using a linear layer. Finally, the geometric features and the semantic features of the self-attention mechanism are concatenated and mapped to the original dimension to form the output.

[0083] Query generator, which is used to first obtain the global feature f through the output feature of the encoder g , and predict the rough missing point cloud center coordinates {c1,…,c i ,…c M}, after concatenating the coordinates of the missing point cloud center point with the global features, the cartilage repair network model is used to generate the query features, i.e.

[0084]

[0085] Multi-scale point cloud generation module, translating query features into missing proxy points In order to reduce redundant calculations, the Query generator is used to replace the missing point cloud center coordinates {c1,…, c i ,…c M} as the local center of the missing point cloud, and use FoldingNet to recover the missing proxy points Detailed local shape centered on:

[0086]

[0087] Among them, P i The coordinates of the center point c of the missing point cloud are i The set of adjacent points centered at is the missing point cloud corresponding to the missing proxy point.

[0088] The above-mentioned method for constructing the cartilage repair network model includes:

[0089] A cartilage repair network model is pre-trained using the existing point cloud completion dataset as the source model;

[0090] Create a new cartilage repair network model as the target model and migrate all the model designs and parameters of the source model to the target model;

[0091] The target model was trained using the knee cartilage dataset and fine-tuned based on the parameters of the source model to obtain a cartilage repair network model based on the fine-grained structure of cartilage.

[0092] The above cartilage repair network model uses the chamfer distance loss function as its loss function, and the loss function is expressed as:

[0093]

[0094] Where C represents n C local center points, P represents the final generated complete point cloud n P points, G represents the n of the original point cloud G points.

[0095] The predicted local center and the center of the input point cloud are connected to form the local center of the entire object C. The high-resolution point cloud G is used to supervise the sparse point cloud C. The final objective function is J = αJ0 + βJ1. Among them, α and β are parameters for adjustment and weight.

[0096] The loss function of point cloud completion provides a quantitative measurement of output quality, enabling the entire cartilage point cloud restoration network to be gradually optimized during training.

[0097] The above step S5 specifically includes:

[0098] S51: For any point P in the point cloud P i Set the search radius r, and record the neighborhood points within the search radius r as set N(P i ),Right now

[0099] N(P i )={P j |P j ∈P,‖P i -P j ‖ <r}

[0100] S52: Set the surface equation ax+by+cz=d(d≥0),a 2 +b 2 +c 2 =1, take P i The corresponding set N(P i ), calculate N(P i ) to the surface

[0101] d i =|ax+by+cz-d|

[0102] Solve N(P i ) The eigenvector corresponding to the minimum sum of the distances from all points to the surface is the normal vector of the point

[0103]

[0104] According to P i Make the tangent plane Ψ of the point with the normal vector n;

[0105] S53: Set N(P i ) is projected onto the tangent plane Ψ, denoted as N(P i ′). In N(P i ′) select a point P j ',by is the u axis, n is the w axis, u×w is the v axis, and P i Construct a local coordinate system for the coordinate center, denoted as (P i ,u,v,w);

[0106] S54: Calculate the set N(P i Other points P in n ' to point P i The vector P i P n ′ with the coordinate axis u clockwise angle Ω=(θ1,θ2,θ3,…θ n ), make the difference between the adjacent angles and get the angle set Ω′=(θ′1,θ′2,θ′3···θ′ n-1 ),in

[0107] θ′n-1=θ n -θ n-1

[0108] Arrange the elements in the set Ω′ in descending order and find the largest angle θ′ max .

[0109] θ′ max =max(θ′1,θ′2,θ′3…θ′ n-1 )

[0110] When θ′ max When it is greater than the threshold, the point is changed to a boundary point.

[0111] Furthermore, step S6 specifically includes:

[0112] S61: Initialize the number of loops k = 0 and the local point set Inputs = NULL;

[0113] S62: Randomly select 3 points from the boundary point set D, and use these 3 points to solve the parameters of the circle, that is, the center (x0, y0, z0) and the radius R;

[0114] S63: Calculate the distance d from each boundary point to the center of the obtained circle. If dR≤ε, the point is counted as an internal point set; otherwise, it is considered an external point.

[0115] S64: Calculate the number of points inside the circle, record it as M, if M is greater than the threshold M min , then the estimation is considered successful and the process goes to step S65; otherwise, the process goes to step S26;

[0116] S65: recalculate the parameter model of the circle using the minimum multiplication method for all the parameters in the set Inputs to obtain the final result;

[0117] S66: k=k+1, if k>k max , then end; otherwise go to S62.

[0118] S67: Using the parametric model of a circle to obtain the equation of a sphere

[0119] (x-x0) 2 +(y-y0) 2 +(z-z0) 2 =R 2 .

[0120] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and are not limiting. Although the present invention has been described in detail with reference to the preferred embodiments, those skilled in the art should understand that the technical solutions of the present invention may be modified or replaced by equivalents without departing from the purpose and scope of the technical solutions of the present invention, which should all be included in the scope of the claims of the present invention.

Claims

1. A cartilage repair method based on point cloud completion, characterized in that: include: S1: Preprocessing the knee cartilage model in the system to construct a knee cartilage dataset; the knee cartilage dataset includes a defective femoral point cloud set consisting of a plurality of defective femoral point clouds; S2: First, the defect cloud of the repaired knee cartilage is downsampled to obtain the center point; then, local features around the center point are extracted, and the local feature point proxy obtained by position embedding is added to the local features; S3: Use the point cloud completion model to predict missing proxy points based on local feature point proxies, and then restore the missing point cloud corresponding to the missing proxy points from coarse to fine based on the predicted missing proxy points; The point cloud completion model is a Transformer encoder-decoder architecture, including: An encoder, configured to encode the local feature point proxy obtained in step S2; Geometric perceptron, used to capture geometric relationships in point clouds using the kNN model; given a query coordinate p Q , according to the bond coordinate p k The nearest key feature is queried and adjacent point proxies in space are spliced ​​together. Then, the local geometric structure is learned by performing a maximum pooling operation through feature aggregation using a linear layer. Finally, the geometric features and the semantic features of the self-attention mechanism are concatenated and mapped to the original dimension to form the output. Query generator, which is used to first obtain the global feature f through the output feature of the encoder g , and predict the rough missing point cloud center coordinates {c1,…,c i ,…c M }, after concatenating the coordinates of the missing point cloud center point with the global features, the cartilage repair network model is used to generate the query features, i.e. Multi-scale point cloud generation module, translating query features into missing proxy points Use the Query Builder to replace the missing point cloud center coordinates {c1,…, c i ,…c M } as the local center of the missing point cloud, and use FoldingNet to recover the missing proxy points Detailed local shape centered on: Among them, C i The coordinates of the center point c of the missing point cloud are i The set of adjacent points centered at is the missing point cloud corresponding to the missing proxy point; S4: splicing the missing point cloud obtained in step S3 with the incomplete femoral point cloud; S5: Extract boundary points of missing point cloud; S6: Use RANSAC to perform 3D circle fitting in the 3D point cloud space to obtain the spherical equation; S7: Determine whether all points in the missing point cloud are within the sphere, delete the points that are not within the sphere, and leave the points within the sphere to form the point cloud of the area to be repaired. Combine the point cloud of the area to be repaired with the defective point cloud to obtain the healthy cartilage point cloud.

2. The cartilage repair method based on point cloud completion according to claim 1, characterized in that: Step S1 specifically includes: S11: obtaining an STL model of the knee joint and femur, and converting the STL model of the knee joint and femur into a point cloud format to obtain a point cloud model of the knee joint and femur; S12: using a straight-through filter to filter out all points in the Z-axis dimension that are not within a given value range, and segmenting the knee joint femur point cloud model; S13: Randomly sample the segmented knee joint femoral point cloud model to obtain the original femoral point cloud set; S14: Cutting the upper surface point cloud of the femur according to the position and size of the cartilage defect in the actual surgery, cutting each original femoral point cloud into several incomplete femoral point clouds, and obtaining an incomplete femoral point cloud set.

3. The cartilage repair method based on point cloud completion according to claim 2, characterized in that: Step S2 specifically includes: S21: Sampling the farthest point of the incomplete femoral point cloud to obtain N fixed center points {q1,…,q i ,…,q N }; S22: Use lightweight DGCNN to extract features from the local area around the center point and obtain N local area features {F1′,…,F i ,…,F N ′}, where F i ′ corresponds to q i The local area features of the center point; S23: Through MLP network Extract each local feature F i ′ is embedded in the position, and the point proxy is obtained after addition, that is, 4. The cartilage repair method based on point cloud completion according to claim 1, characterized in that: The method for constructing a cartilage repair network model comprises: A cartilage repair network model is pre-trained using the existing point cloud completion dataset as the source model; Create a new cartilage repair network model as the target model and migrate all the model designs and parameters of the source model to the target model; The target model was trained using the knee cartilage dataset and fine-tuned based on the parameters of the source model to obtain a cartilage repair network model based on the fine-grained structure of cartilage.

5. The cartilage repair method based on point cloud completion according to claim 4, characterized in that: The cartilage repair network model adopts the chamfer distance loss function as its loss function.

6. The cartilage repair method based on point cloud completion according to claim 1, characterized in that: The step S5 specifically includes: S51: For any point P in the point cloud P i Set the search radius r, and record the neighborhood points within the search radius r as set N(P i ),Right now N(P i )={P j |P j ∈P,‖P i -P j ‖<r} S52: Set the surface equation ax+by+cz=d(d≥0),a 2 +b 2 +c 2 =1, take P i The corresponding set N(P i ), calculate N(P i ) to the surface d i =|ax+by+cz-d| Solve N(P i ) The eigenvector corresponding to the minimum sum of the distances from all points to the surface is the normal vector of the point According to P i Make the tangent plane Ψ of the point with the normal vector n; S53: Set N(P i ) is projected onto the tangent plane Ψ, denoted as N(P i ′); in N(P i ′) select a point P j ',by is the u axis, n is the w axis, u×w is the v axis, and P i Construct a local coordinate system for the coordinate center, denoted as (P i ,u,v,w); S54: Calculate the set N(P i Other points P in n ' to point P i The vector P i P n ′ with the coordinate axis u clockwise angle Ω=(θ1,θ2,θ3,…θ n ), make the difference between the adjacent angles and get the angle set Ω′=(θ′1,θ′2,θ′3…θ′ n-1 ),in θ′ n-1 =θ n -θ n-1 Arrange the elements in the set Ω′ in descending order and find the largest angle θ′ max ; θ′ max =max(θ′1,θ′2,θ′3…θ′ n-1 ) When θ′ max When it is greater than the threshold, the point is changed to a boundary point.

7. The cartilage repair method based on point cloud completion according to claim 6, characterized in that: The step S6 specifically includes: S61: Initialize the number of loops k = 0 and the local point set Inputs = NULL; S62: Randomly select 3 points from the boundary point set D, and use these 3 points to solve the parameters of the circle, that is, the center (x0, y0, z0) and the radius R; S63: Calculate the distance d from each boundary point to the center of the obtained circle. If dR≤ε, the point is counted as an internal point set; otherwise, it is considered an external point. S64: Calculate the number of points inside the circle, record it as M, if M is greater than the threshold M min , then the estimation is considered successful and the process goes to step S65; otherwise, the process goes to step S26; S65: recalculate the parameter model of the circle using the minimum multiplication method for all the parameters in the set Inputs to obtain the final result; S66: k=k+1, if k>k max , then end; otherwise go to S62; S67: Using the parametric model of a circle to obtain the equation of a sphere (x-x0) 2 +(y-y0) 2 +(z-z0) 2 =R 2 。

Images