A method for locating a tracer ball under three-dimensional imagery

Abstract

The application discloses a tracer ball positioning method under three-dimensional image, comprising the following steps: scanning to obtain three-dimensional image of a registration piece on which a marker and a tracer ball are installed; identifying the position of the marker in the three-dimensional image, registering the position with the actual marker position, and calculating the position of the actual tracer ball in the image according to the design parameters of the registration piece, thereby identifying the accurate position of the tracer ball. The application firstly coarsely positions the marker, calculates the area of the tracer ball in the image by using the prior relationship between the marker and the tracer ball, then locally searches the area of the tracer ball in the image obtained previously in a small range, thereby reducing the positioning time of the tracer ball, improving the registration accuracy, and reducing the misregistration.

Description

Technical Field

[0001] This invention relates to the field of image processing, and more particularly to a method for locating a tracer sphere in a three-dimensional image. Background Technology

[0002] Lung cancer is one of the most life-threatening cancers, and lung biopsy is the gold standard for diagnosing lung tumors. Lung biopsy surgery requires a surgeon with extensive clinical experience to remove the lesion for pathological analysis. The main principle of a surgical navigation system is to register the preoperative CT image coordinate system with the intraoperative human body coordinate system, thereby enabling the surgical instruments to be displayed in real-time in the preoperative CT image coordinate space.

[0003] The patent application with application number CN202111181190.8, entitled "Marker Ball Positioning Method and Its Application", provides a marker ball positioning method. However, the marker ball needs to be searched globally in the 3D image, which is very time-consuming. At the same time, the marker ball positioning only uses a threshold, which may lead to misregistration. Summary of the Invention

[0004] Purpose of the invention: To address the above-mentioned shortcomings, this invention proposes a method for locating a tracer sphere under three-dimensional images, which can reduce marker positioning time, improve registration accuracy, and reduce misregistration.

[0005] Technical solution:

[0006] A method for locating a tracer sphere in a 3D image, comprising the following steps:

[0007] The scan yielded a three-dimensional image of the registration document, on which markers and tracer balls were mounted.

[0008] The location of the marker in the 3D image is identified, registered with the actual location of the marker, and the actual location of the tracer ball in the image is calculated based on the design parameters of the registration document. The precise location of the tracer ball is then identified.

[0009] The location of markers in the 3D image is identified by threshold filtering.

[0010] The threshold filtering specifically involves:

[0011] Histogram analysis was performed on the CT values ​​of each voxel in the 3D image.

[0012] The optimal group index on the histogram is determined based on the number of voxel points occupied by the marker, thereby calculating the screening threshold of the marker and selecting the corresponding voxel points accordingly.

[0013] Candidate markers are obtained by clustering, and the number of candidate markers corresponding to the markers on the registration document are selected and matched with the markers on the registration document in terms of relative position. If they match, they are combined as candidate markers.

[0014] The candidate marker combination is registered with the marker combination on the registration document, and the candidate marker combination with the smallest registration error is selected as the final marker combination, thereby identifying the marker position in the 3D image.

[0015] The registration error is specifically the root mean square error between the two.

[0016] The precise location of the tracer ball in the aforementioned area is determined by sequentially performing circle detection, circle screening, circle clustering, centroid determination, tracer ball combination matching, and registration to locate the tracer ball.

[0017] The circle detection specifically involves:

[0018] The region of interest for each tracer sphere in the 3D image is calculated based on the actual position of the tracer sphere in the image and the radius of the tracer sphere.

[0019] Traverse the planar regions of interest of each cross-sectional image corresponding to the region of interest of the tracer sphere on the 3D image;

[0020] Circle detection is performed using the Hough transform on a region of interest in a plane.

[0021] The circular filtering specifically refers to:

[0022] Four rectangular regions with set widths and lengths are set on the edge of a detected planar circle, located at the left, right, top, and bottom boundaries respectively. The proportion of pixels with pixel values ​​less than a fixed threshold in each rectangular region is determined. If the proportion is greater than the set background threshold, the rectangular region is the background; otherwise, it is the foreground.

[0023] If at least two rectangular regions are the background, then the planar circle is considered to be the planar circle of the tracer sphere on the slice.

[0024] The width and length of the rectangular region are specifically as follows:

[0025] The width of the rectangular region The length is 2 * (radius of the plane circle + pad);

[0026] Where sx and sy represent the lengths of the pixels on the plane containing the circle along the x and y axes of the corresponding planar coordinate system.

[0027] The specific matching of the tracer ball combination is as follows:

[0028] Iterate through and select the number of tracer ball combinations that correspond to the actual number of tracers on the registered document, and perform relative positional matching. If a match is found, then use it as a candidate tracer ball combination.

[0029] The specific registration of the tracer sphere combination is as follows:

[0030] The candidate tracer ball combination obtained by matching the tracer ball combination is registered with the actual tracer ball combination on the registration document, and the candidate tracer ball combination with the smallest registration error is selected as the final tracer ball combination.

[0031] It also includes further positioning steps:

[0032] Based on the coordinates of the centroid obtained from the localization of the tracer sphere and the radius of the tracer sphere, the regions of interest are obtained on the cross-sections of the centroids of each tracer sphere in the coronal and sagittal planes, as well as the regions of interest on the two cross-sections before and after it. The outer contour of each tracer sphere on each cross-section is obtained, and then the outer contour points are fitted to calculate the new centroid of the tracer sphere. The registration error is recalculated. If the updated registration error is less than the previous registration error, the position of the tracer sphere is updated.

[0033] The identification of the precise location of the tracer ball also includes:

[0034] Set the initial and final window levels, set the adjustment step size from the initial to the final window level, set a fixed window width, adjust the window level and window width to obtain images with different window levels and window widths, and identify the precise position of the tracer ball;

[0035] If the corresponding tracer ball is identified, the registration error is updated. If this registration error is less than the previous registration error, the tracer ball combination and the corresponding tracer ball centroid are updated.

[0036] After updating the tracer ball combination, determine whether the proportion of the number of times the corresponding tracer ball was detected after the window level and window width adjustment is greater than the proportion of all adjustment times. If it is greater, update the obtained tracer ball combination as the final tracer ball combination; otherwise, readjust the window level and window width.

[0037] The set ratio is 0.6.

[0038] The markers and the tracer spheres are all in at least three units, and are all mounted on the registration piece in a coplanar but not collinear manner.

[0039] Beneficial effects: This invention first performs coarse localization of the markers, calculates the region of the tracer ball in the image using the prior relationship between the markers and the tracer ball, and then performs a small-scale local search on the previously obtained region of the tracer ball in the image, thus reducing the localization time of the tracer ball. Simultaneously, this invention performs fine localization through iterative adjustment of window level and window width, improving registration accuracy and reducing misregistration. Attached Figure Description

[0040] Figure 1 This is a flowchart of the tracer ball positioning method of the present invention;

[0041] Figure 2 A schematic diagram of the structure of a registration document used for navigation and positioning;

[0042] Figure 3 A flowchart illustrating the coarse positioning of markers;

[0043] Figure 4 The CT image after window level and window width adjustment;

[0044] Figure 5 This is a schematic diagram of the human body coordinate system.

[0045] Figure 6 This is a schematic diagram of the positioning process for the tracer ball;

[0046] Figure 7 This is a schematic diagram of the circular clustering algorithm.

[0047] Among them, 1 is the registration document, 11 is the bracket, 12 is the marker, and 13 is the tracer ball. Detailed Implementation

[0048] The present invention will be further explained below with reference to the accompanying drawings and specific embodiments.

[0049] The method for locating a tracer sphere under three-dimensional images in this invention is as follows: Figure 1 As shown, it includes the following steps:

[0050] (1) Scan the three-dimensional image of the registration document used for navigation and positioning;

[0051] Figure 2 This is a schematic diagram of the structure of the registration document for navigation and positioning according to the present invention, as shown below. Figure 2 As shown, the registered component 1 of this invention is installed on the human body surface near the lesion, including a support 11 and markers 12 and tracer balls 13 mounted on the support 11. Specifically, the support 11 is made of polyetheretherketone (PEEK); at least three markers 12 are mounted on the support 11 in a coplanar but not collinear manner; more specifically, the markers 12 are steel balls with a diameter of 2.5 mm, which minimizes artifacts on CT images due to their small diameter. In this invention, four steel balls are mounted on the support 11, with the centers of each ball coplanar but not collinear. The tracer balls 13 are made of frosted material, with CT values ​​similar to those of human soft tissue, and do not produce artifacts on CT images like metal markers, thus not affecting the doctor's diagnosis. They have a diameter of 11-13 mm and are at least three mounted on the support 11 in a coplanar but not collinear manner; in this invention, four tracer balls 13 are mounted on the support 11 in a coplanar but not collinear manner, forming a tracer assembly.

[0052] This invention uses a CT scanner to scan the registration document to obtain its CT image.

[0053] (2) Perform coarse localization of markers in the 3D image;

[0054] This invention employs threshold filtering to locate the position of a marker, achieving coarse localization. Specifically, since the steel ball has a high CT value in CT images, threshold filtering can be used to select voxel points with high CT values, significantly reducing the number of voxel points and allowing for rapid location of the steel ball. However, in areas where the human body experiences significant breathing amplitude, motion artifacts will exist, affecting the positioning accuracy of the steel ball. Therefore, the steel ball can only be used for coarse localization. Figure 3 As shown, the steps are as follows:

[0055] (21) Histogram statistics were performed on the CT values ​​of each voxel in the three-dimensional image;

[0056] If the number of groups in the histogram is set to bin, then the CT value interval (i.e., step size) between adjacent groups is step = (Z) / bin. max -Z min ) / bin, where Z max Z min These represent the maximum and minimum CT values ​​of voxel points in the 3D image, respectively.

[0057] The CT value of each voxel in the 3D image is denoted as Z. i Then, the group index D of the histogram corresponding to each voxel point is calculated using the following formula. i At the same time, the number of voxels in the group corresponding to that voxel is incremented by 1, where the following formula is used to calculate and round down to the nearest integer to obtain D. i ;

[0058] D i =(Z i -Z min ) / bin;

[0059] In this invention, the number of groups in the histogram is set to bin = 1000, and the index of the first group is 0. Therefore, the group index of the histogram is [0, 999].

[0060] (22) Calculate the screening threshold for the markers;

[0061] Traverse the histogram from the largest group index down to index 0, recording the total number of voxels N within each traversed group. When N is greater than a set threshold A0, the optimal group index of the histogram is determined to be D. best ;

[0062] Set the threshold A0 = k * 4 / 3πr1 3 / (sx*sy*sz);

[0063] Where k is the threshold coefficient, with a value range of (0,1), and is set according to actual needs; r1 is the radius of the marker; sx, sy, and sz are the lengths of the space occupied by the voxel points on the three coordinate axes of the image coordinate system, respectively.

[0064] Therefore, the screening threshold for the marker can be determined as L = Z. min +bin*(D best +1);

[0065] (23) Mark object pixel filtering;

[0066] Traverse each voxel in the 3D image. If the CT value of a voxel is greater than or equal to the screening threshold L, then put the voxel's 3D coordinates in the image coordinate system, CT value, and access status (initial access status is false) into set O.

[0067] (24) Identify voxel combinations of suspected markers;

[0068] 1) Randomly select a voxel in set O. If the access status of the voxel is false, set the voxel as a seed point, set the access status of the seed point to true, and put it into the seed point set A.

[0069] 2) Set the current seed point as the center point, and search for surrounding neighboring points with a distance of 1 along the x, y, and z axes of the image coordinate system. Determine whether the corresponding neighboring point is in set O. If the point is in set O, set the neighboring point as the seed point and set its access status to true, and put it into set A. Iterate through the seed points in set A and repeat step 2) to find surrounding neighboring points with a distance of 1 from the seed point in set O, perform clustering to obtain a voxel point combination of a suspected marker, put it into set B, and set the access status of the obtained neighborhood point to true.

[0070] 3) Repeat steps 1) to 2) until all voxel points in set O with a visit status of false are traversed. The resulting set B is the combination of voxel points of each suspected marker.

[0071] (25) Identify candidate markers;

[0072] Since the size of the marker is fixed, the number of voxel points within the marker is also fixed, that is, n = 4 / 3πr1 3 / (sx*sy*sz), then the lower limit of the number of voxel points of the marker is set to n. lower =n / 10, setting the maximum number of voxel points for markers to n. upper=n*10, meaning the number of voxel points in the combination of voxel points of each suspected marker in set B is less than n. lower or greater than n upper The voxel point combinations are deleted and set B is updated to obtain the voxel point combinations of the candidate markers, that is, the candidate markers.

[0073] (26) Calculate the centroid of each candidate marker;

[0074] Iterate through each combination of voxels in set B, calculate the centroid of each combination, which is the centroid of each candidate marker. The specific centroid p of a combination of voxels is obtained. t The calculation is as follows:

[0075]

[0076]

[0077]

[0078] Where, p t .x represents the centroid p of this voxel-point combination. t The x-coordinate value, p t .y represents the centroid p t The y-coordinate value, p t z is the centroid point p t The z-coordinate value, p j .x represents a voxel p in this voxel combination. j The x-coordinate value, p j .y represents a voxel point p j The y-coordinate value, p j .z represents a voxel point p j The z-coordinate value, where n is the number of voxels in each voxel combination.

[0079] (27) Obtain candidate marker combinations;

[0080] The centroid coordinates of each candidate marker are obtained through step (26). Four candidate markers in set B are randomly selected as the set of points to be matched. Their relative positional relationships are matched one-to-one with the relative positional relationships of the actual point set composed of the actual markers in the registration document. The set of points to be matched that meets the requirements of the relative positional relationship with the actual point set in the registration document is traversed and selected as the candidate point set, that is, the candidate marker combination, specifically:

[0081] Randomly select 4 voxel point combinations (i.e. candidate markers) from set B as the set of points to be matched, and determine whether the relative positional relationship of the 4 candidate markers in the set of points to be matched matches the relative positional relationship of the actual point set (using coordinate measuring machine). If they match, put the set of points to be matched into the set of candidate marker combinations. Otherwise, repeat the step until all voxel point combinations in set B have been traversed.

[0082] In this invention, a voxel point combination can be randomly selected first, and then other voxel point combinations with the same positional relationship between the voxel point combination and the relative positional relationship between the four markers in the actual point set can be selected to form a set of points to be matched, and then matched with the actual point set.

[0083] Specifically, whether the relative positional relationship of the point set to be matched matches the relative positional relationship of the actual point set is determined by the following three points:

[0084] a. Whether the deviation between the Euclidean distance between any two points in the set of points to be matched and the Euclidean distance between the corresponding two points in the actual set of points is less than the first threshold;

[0085] b. Whether the deviation between the included angle between any two lines connecting any three points in the set of points to be matched and the actual included angle between the corresponding two lines in the actual set of points is less than the second threshold.

[0086] c. Whether the deviation between the flatness of the plane containing the four points of the point set to be matched and the actual flatness of the plane containing the four points of the actual point set is less than the third threshold;

[0087] When the deviations of a, b, and c are all within the corresponding threshold range, the set of points to be matched is determined to match the actual set of points; otherwise, it is determined to be a mismatch. The first, second, and third thresholds are set as needed.

[0088] (28) Registration yields the final marker combination;

[0089] The candidate marker combination calculated in step (27) is registered with the actual point set formed by the actual marker combination in the registered document under the design coordinate system (i.e., marker coordinate system) to obtain the registration matrix. Specifically, the centroid of each marker is used for registration, and the candidate marker combination with the smallest projection error is selected as the marker combination under the image coordinate system.

[0090] The registration adopts the rigid body registration method (common method SVD decomposition). The error evaluation standard index is to transform the candidate marker combination in the image coordinate system to the design coordinate system, and calculate the root mean square error (RMSE) between it and the actual marker combination in the registered part, which is used as the projection error of the registration.

[0091] Transformation formula: P2 = R * P1 + T, where R is the rotation transformation between the image coordinate system and the design coordinate system, and T is the corresponding translation transformation;

[0092] The root mean square error is:

[0093] P1 represents the centroid coordinates of the candidate marker in the image coordinate system, P2 represents the centroid coordinates of the candidate marker transformed from the image coordinate system to the design coordinate system, P3 represents the centroid coordinates of the actual marker in the design coordinate system, and k represents the centroid of the kth candidate marker.

[0094] (3) Determine the region of the tracer sphere in the 3D image;

[0095] Based on step (2), the transformation relationship T between the marker coordinate system and the three-dimensional image coordinate system is obtained. img1 The actual centroids of the marker material and the tracer sphere of the registered component were measured using a coordinate measuring machine, thereby obtaining the transformation relationship T between the marker coordinate system and the world coordinate system. w1 Transformation relationship T between the tracer sphere coordinate system and the world coordinate system w2 Thus, the transformation relationship between the marker coordinate system and the tracer sphere coordinate system can be calculated: T 12 =T w1 *T w2 -1 Finally, the transformation relationship T between the tracer sphere coordinate system and the image coordinate system is obtained. img2 =T img1 *T 12 :

[0096] This gives us the coordinates of the actual centroid of the tracer sphere in the image coordinate system. By combining this with the actual radius of the tracer sphere, we can determine the area of ​​the tracer sphere in the 3D image. However, due to motion artifacts caused by human breathing, the position of the tracer sphere is not very accurate, so precise positioning is required.

[0097] (4) Coarse positioning of the tracer sphere;

[0098] Tracer sphere localization is achieved through circle detection, circle screening, circle clustering, centroid determination, tracer sphere combination matching, and registration.

[0099] Because the marker is a sphere, along a cross-section parallel to the sphere ( Figure 5 (xy plane) In each slice, the sphere is a planar circle, and the center coordinates of each circle are the same. As long as the center of more than half of the planar circles in all the slice circles are detected, and multiple center circles are connected by a straight line, the midpoint of the straight line is the centroid, and the centroid of the marker can be located, ensuring the robustness of the coarse positioning algorithm. Repeat the above steps to obtain the centroid of all suspected spheres.

[0100] The algorithm flow is as follows Figure 6 As shown, the details are as follows:

[0101] (41) Adjustment of the window level and width of the 3D image;

[0102] The window level and width are determined. Since the window level and width of the tracer sphere are relatively fixed, the window level is set to -100 and the window width to 400. For the original 3D image containing markers, it is first converted to Hu values, and then the window level and width are adjusted to normalize the data to the [0, 255] interval. Figure 4 As shown below:

[0103] p v =(h v -h min )*255 / (h max -h min )

[0104] h min =(2*W L -W W ) / 2+0.5

[0105] h max =(2*W L +W W ) / 2+0.5

[0106] Where, p v h represents the normalized grayscale value of a pixel in the original CT image. v Hu represents the value, h min h represents the minimum Hu value of the original CT image. max W represents the maximum value of the Hu value in the original CT image. L Indicates window level, W W Indicates window width; h min and h max The final value is the actual calculated value rounded up to the nearest integer.

[0107] (42) Planar circle detection;

[0108] Based on step (3), the coordinates of the actual tracer sphere centroid in the image coordinate system are obtained, which is q(x,y,z). The region of interest (ROI) of each tracer sphere in the 3D image is determined by the following calculation method:

[0109] roi.x = qx - 2 * r2;

[0110] roi.y = qy - 2 * r2;

[0111] roi.z = qz - 2 * r2;

[0112] roi.width = 2 * r2;

[0113] roi.height = 2 * r2;

[0114] roi.depth = 2 * r2;

[0115] Where roi.x, roi.y, and roi.z are the coordinates of the starting point of the region of interest in the x, y, and z directions, respectively; qx, qy, and qz are the coordinates of the centroid of the tracer sphere on the x, y, and z axes of the image coordinate system, respectively; r2 is the radius of the tracer sphere; and roi.width, roi.height, and roi.depth are the width, height, and depth of the region of interest roi in the 3D image in the image coordinate system, respectively.

[0116] Traverse the cross-sections corresponding to the ROI regions of the tracer sphere in the 3D image (e.g.) Figure 4 The planar region of interest P of each cross-sectional image of ) roi The details are as follows:

[0117] P roi .x = roi.x / sx;

[0118] P roi .y = roi.y / sy;

[0119] P roi .width = roi.width / sx;

[0120] P roi .height = roi.height / sy;

[0121] Among them, P roi .x、P roi .y represent the coordinates of the tracer sphere's centroid on the x-axis and y-axis of the image coordinate system, respectively. P roi .width represents the planar region of interest P corresponding to the tracer sphere. roi The length P on the x-axis of the corresponding two-dimensional coordinate system on the cross-sectional slice. roi .height represents the planar region of interest P corresponding to the tracer sphere. roi The length on the y-axis of the corresponding two-dimensional coordinate system on the cross-sectional slice;

[0122] In P roi By using the Hough transform to detect circles in a two-dimensional image, this invention can significantly reduce the time required for circle detection compared to the global search for circles across all cross-sections described in the patent mentioned in the background art.

[0123] (43) Planar circle screening;

[0124] If there are too many background planar circles, it will lead to too many subsequent tracer ball combinations, resulting in low tracer ball combination efficiency. Therefore, it is necessary to filter out some planar circles that do not meet the requirements. The filtering rule is to set four rectangular regions on the edge of a detected planar circle, located at the left, right, top, and bottom boundaries respectively. The proportion of pixels with pixel values ​​less than a fixed threshold in each of the four rectangular regions is determined. If the proportion is greater than the set background threshold, the rectangular region is considered the background; otherwise, it is considered the foreground. In this invention, the fixed threshold is 30, and the background threshold is 0.6. If at least two of the four rectangular regions are background, the planar circle is considered to be a planar circle of the tracer ball on the slice.

[0125] In this invention, the rectangular region is specifically defined as a rectangular region with a predetermined width and a predetermined length set on the top, bottom, left, and right boundaries of a planar circle, specifically as follows:

[0126] On the plane containing the circle, the positive x-axis is defined as the horizontal direction to the right, and the positive y-axis is defined as the vertical direction downwards; set the width. The length is set to 2*(c.r+pad), where sx and sy are the lengths of the space occupied by the voxel on the x and y axes of the image coordinate system, respectively. In the plane coordinate system, this is the length of the pixel on the x and y axes of the plane coordinate system.

[0127] The rectangular region on the left boundary of the planar circle is then expanded downwards and to the right from its starting point, as follows:

[0128] l roi .x = cx - cr - pad

[0129] l roi .y = cy-cr-pad

[0130] l roi .width=pad

[0131] l roi .height = 2 * (c.r + pad)

[0132] Where cx and cy are the coordinates of the center of the planar circle on the x-axis and y-axis of the planar coordinate system, cr is the radius of the planar circle, and l roi .x、l roi .y represent the x-coordinates and y-coordinates of the starting point of the rectangular region on the left boundary of the planar circle, respectively, on the planar coordinate system. roi .width、l roi .height represents the width and length of the rectangular region on the left boundary of the planar circle, respectively;

[0133] The rectangular region on the right boundary of the planar circle is obtained by extending it downwards and to the right from its starting point, as follows:

[0134] r roi x = c.x + cr

[0135] r roi .y = cy-cr-pad

[0136] r roi .width=pad

[0137] r roi .height = 2 * (c.r + pad)

[0138] Where, r roi .x、r roi .y represent the x-coordinates and y-coordinates of the starting point of the rectangular region on the right boundary of the planar circle, respectively, on the planar coordinate system. r roi .width、r roi .height represents the width and length of the rectangular region on the right boundary of the plane circle, respectively;

[0139] The rectangular region on the upper boundary of the planar circle is obtained by extending it downwards and to the right from its starting point, as follows:

[0140] u roi .x = cx - cr - pad

[0141] u roi .y = cy-cr-pad

[0142] u roi .width = 2 * (c.r + pad)

[0143] u roi .height=pad

[0144] Among them, u roi .x、u roi .y represent the x-coordinates and y-coordinates of the starting point of the rectangular region on the upper boundary of the plane circle, respectively, on the plane coordinate system. u roi .height、u roi .width represents the width and length of the rectangular region on the upper boundary of the plane circle, respectively;

[0145] The rectangular region at the lower boundary of the planar circle is obtained by extending downwards and to the right from its starting point, as follows:

[0146] d roi .x = cx - cr - pad

[0147] d roi .y=c.y+cr

[0148] d roi .width = 2 * (c.r + pad)

[0149] d roi .height=pad

[0150] Where, d roi .x、d roi .y represent the x-coordinates and y-coordinates of the starting point of the rectangular region at the lower boundary of the plane circle on the plane coordinate system, respectively. d roi .height、d roi .width represents the width and length of the rectangular region at the lower boundary of the plane circle, respectively.

[0151] (44) Circular clustering;

[0152] Based on the property that the centers of the planar circles on each cross-section of the same sphere are theoretically identical in a two-dimensional coordinate system parallel to that cross-section, the planar circles obtained in step (43) are clustered, such as... Figure 7 As shown, the specific circular clustering process is as follows:

[0153] Set the access status of all planar circles to false, and designate one of them as the seed circle. Represent the center and radius of the seed circle with vector b, and set its access status to true. Iterate through and detect other circles, calculating the Euclidean distance between the center of the detected circle and the center of the seed circle. If the Euclidean distance is less than a set threshold of 6mm, represent the center and radius of the detected circle with vector a, set its access status to true, and classify it with the seed circle. Continue this process until the seed circle has traversed all circles. Finally, group all circles with vector b and seed circles with vector a together to obtain the combination of circles in all planes parallel to the cross section for a suspected tracer sphere. Then, iterate through the planar circles with access status false, select one circle as the seed circle, and repeat the above clustering process to obtain the combination of circles for all suspected tracer spheres.

[0154] (45) Find the centroid and radius of the suspected tracer sphere;

[0155] According to step (44), the coordinates of the centers of all circles in the circle combination of each suspected tracer ball obtained by clustering in the same two-dimensional coordinate system parallel to the cross section satisfy a Gaussian distribution. μ is the mean of all samples, σ is the standard deviation of all samples, and the confidence interval of the normal distribution n~N(μ,σ) is [μ-λσ,μ+λσ], where λ depends on the confidence level. In this invention, λ is equal to 3.

[0156] The x and y coordinates of the centroid of a tracer sphere are the average values ​​of the samples within the confidence interval, the z coordinate is the median z coordinate of each circle combination obtained by circle clustering, and the radius of the tracer sphere is the largest radius among all circles in its corresponding circle combination.

[0157] (46) The candidate tracer ball combinations are calculated using the method in step (27);

[0158] (47) The final combination of tracer spheres is calculated using the method in step (28), thereby calculating the three-dimensional coordinates of the centroid of each tracer sphere in the image coordinate system and the transformation relationship of the coordinate system corresponding to each tracer sphere relative to the image coordinate system.

[0159] (48) Further positioning of the tracer ball;

[0160] Because the CT scan slices are too thick, the tracer sphere may not be fully scanned. Furthermore, the tracer sphere's circular planes at the initial and final sections of the transverse section may not be detected, thus affecting the precise positioning of the tracer sphere's centroid in the transverse plane. The tracer sphere's slice images in the coronal and sagittal planes are also circular (the coronal plane is...). Figure 4 The xz plane, the sagittal plane is Figure 4 The yz plane can be used to locate the precise coordinates of the tracer sphere's centroid on the cross section. Based on the three-dimensional coordinates and radius of the centroid of each tracer sphere obtained in step (4), the regions of interest (ROI) of the cross section where the centroid of each tracer sphere is located and the two cross sections before and after it are obtained in the coronal and sagittal planes. The outer contour of each tracer sphere on each cross section is obtained, and its outer contour points are fitted to calculate the new centroid of the tracer sphere. Then, the registration is performed again in step (47) to update the projection error. If the updated projection error is less than the previous projection error, the three-dimensional coordinates of the tracer sphere are updated.

[0161] The fitting method is as follows:

[0162] a. Based on the centroid coordinates of the tracer sphere obtained above, obtain the coordinates of the two-dimensional center O of the tracer sphere on the current sectional plane, and define an ideal circle with the radius of the tracer sphere; remove points on the sectional plane that are more than a set threshold (in this invention, the set threshold is set to 6 mm);

[0163] b. Perform least-squares circle fitting on the remaining outer contour points to obtain the center of the remaining outer contour points. Find the center of the circle with the smallest Euclidean distance from the center O among the three cross planes and use it as the new center of the circle. Use it as the new center of mass of the tracer sphere.

[0164] (5) Precise positioning of the tracer ball

[0165] Since step (41) uses window level and window width adjustment to perform coarse positioning of the area where the tracer ball is located, the positioning of the tracer ball may not be optimal, resulting in the positioning accuracy of the tracer ball not reaching the best. Therefore, this invention iteratively adjusts the window level and window width, and after each adjustment, iterative processing is performed through step (4) to obtain the final tracer balls, as follows:

[0166] (51) Therefore, the present invention sets the initial window level to -200, the final window level to 30, the adjustment step from the initial window level to the final window level to 4, and the window width to be fixed at 400 to adjust the window level and window width, thereby obtaining images with different window levels and window widths;

[0167] (52) Set u1 to the number of iterations, with an initial value of 0; u2 to the number of times no tracer balls were detected, with an initial value of 0; each time the window level and width are adjusted, u1 is incremented by 1, and steps (41) to (48) are repeated; if no tracer balls are obtained in this window level and width adjustment (i.e., no four tracer balls are detected in the image after the window level and width adjustment), then u2 is incremented by 1; if tracer balls are obtained, the centroids of the new four tracer balls are calculated according to the method in step (4), and the corresponding projection error between the registration and the actual marker combination in the registration document is calculated. If the projection error of this registration is less than the projection error obtained in step (2), then the tracer ball combination and the centroids and registration matrix of each tracer ball are updated, thereby improving the positioning accuracy of the marker balls. After the traversal is completed, only when the proportion of detected tracer balls rate = (u1-u2) / u1 > 0.6 is the tracer ball combination considered to be the final tracer ball combination, thereby reducing the phenomenon of misregistration.

[0168] This invention first performs coarse localization on four markers, then calculates the region of the tracer ball in the image using the prior relationship between the markers and the tracer ball. Next, it performs a small-scale local search on the previously obtained region of the tracer ball in the image, thus reducing the localization time. Simultaneously, this invention performs fine localization through iterative adjustment of window level and width, which improves registration accuracy and reduces misregistration.

[0169] The preferred embodiments of the present invention have been described in detail above. However, the present invention is not limited to the specific details of the above embodiments. Within the scope of the technical concept of the present invention, various equivalent transformations (such as quantity, shape, position, etc.) can be made to the technical solution of the present invention, and these equivalent transformations all fall within the protection scope of the present invention.

Claims

1. A method for locating a tracer sphere under three-dimensional imagery, characterized in that: Including the following steps: The scan yielded a three-dimensional image of the registration document, on which markers and tracer balls were mounted. The positions of markers in the 3D image are identified and registered with the actual positions of the markers to obtain the final marker combination. The position of the actual tracer ball in the 3D image is calculated by combining the positional relationship between the markers and the tracer ball calculated based on the design parameters of the registered part. Set the initial window level and window width, and then perform circle detection, circle filtering, circle clustering, centroid finding, tracer sphere combination matching and registration to locate the tracer spheres and obtain the tracer sphere combination. Set the final window position, set the adjustment step size from the initial window position to the final window position, set a fixed window width, adjust the window position to obtain images at different window positions, and reposition the tracer ball to obtain the corresponding tracer ball combination; Update the projection error between the corresponding tracer ball combination and the actual tracer ball combination. If the projection error is less than the previous registration projection error, update the tracer ball combination and the centroid of the corresponding tracer ball.

2. The tracer sphere positioning method according to claim 1, characterized in that: The location of markers in the 3D image is identified by threshold filtering.

3. The tracer sphere positioning method according to claim 2, characterized in that: The threshold filtering specifically involves: Histogram analysis was performed on the CT values ​​of each voxel in the 3D image. The optimal group index on the histogram is determined based on the number of voxel points occupied by the marker, thereby calculating the screening threshold of the marker and selecting the corresponding voxel points accordingly. Candidate markers are obtained by clustering, and the number of candidate markers corresponding to the markers on the registration document are selected and matched with the markers on the registration document in terms of relative position. If they match, they are combined as candidate markers. The candidate marker combination is registered with the marker combination on the registration document, and the candidate marker combination with the smallest registration error is selected as the final marker combination, thereby identifying the marker position in the 3D image.

4. The tracer sphere positioning method according to claim 3, characterized in that: The registration error is specifically the root mean square error between the two.

5. The tracer sphere positioning method according to claim 1, characterized in that: The circle detection specifically involves: The region of interest for each tracer sphere in the 3D image is calculated based on the actual position of the tracer sphere in the image and the radius of the tracer sphere. Traverse the planar regions of interest of each cross-sectional image corresponding to the region of interest of the tracer sphere on the 3D image; Circle detection is performed using the Hough transform on a region of interest in a plane.

6. The tracer sphere positioning method according to claim 1, characterized in that: The circular filtering specifically refers to: Four rectangular regions with set widths and lengths are set on the edge of a detected planar circle, located at the left, right, top, and bottom boundaries respectively. The proportion of pixels with pixel values ​​less than a fixed threshold in each rectangular region is determined. If the proportion is greater than the set background threshold, the rectangular region is the background; otherwise, it is the foreground. If at least two rectangular regions are the background, then the planar circle is considered to be the planar circle of the tracer sphere on the slice.

7. The tracer sphere positioning method according to claim 6, characterized in that: The width and length of the rectangular region are specifically as follows: The width of the rectangular region The length is 2 * (radius of the plane circle + pad); Where sx and sy represent the lengths of the pixels on the plane containing the circle along the x and y axes of the corresponding planar coordinate system.

8. The tracer sphere positioning method according to claim 1, characterized in that: The specific matching of the tracer ball combination is as follows: Iterate through and select the number of tracer ball combinations that correspond to the actual number of tracers on the registered document, and perform relative positional matching. If a match is found, then use it as a candidate tracer ball combination.

9. The tracer sphere positioning method according to claim 5, characterized in that: The specific registration of the tracer ball combination is as follows: The candidate tracer ball combination obtained by matching the tracer ball combination is registered with the actual tracer ball combination on the registration document, and the candidate tracer ball combination with the smallest registration error is selected as the final tracer ball combination.

10. The tracer sphere positioning method according to claim 5, characterized in that: It also includes further positioning steps: Based on the coordinates of the centroid obtained from the localization of the tracer sphere and the radius of the tracer sphere, the regions of interest are obtained on the cross-sections of the centroids of each tracer sphere in the coronal and sagittal planes, as well as the regions of interest on the two cross-sections before and after it. The outer contour of each tracer sphere on each cross-section is obtained, and then the outer contour points are fitted to calculate the new centroid of the tracer sphere. The registration error is recalculated. If the updated registration error is less than the previous registration error, the position of the tracer sphere is updated.

11. The tracer sphere positioning method according to claim 1, characterized in that: After updating the tracer ball combination, determine whether the proportion of the number of times the corresponding tracer ball was detected after the window level and window width adjustment is greater than the proportion of all adjustment times. If it is greater, update the obtained tracer ball combination as the final tracer ball combination; otherwise, readjust the window level and window width.

12. The tracer sphere positioning method according to claim 11, characterized in that: The set ratio is 0.

6.

13. The tracer sphere positioning method according to claim 1, characterized in that: The markers and the tracer spheres are all in at least three units, and are all mounted on the registration piece in a coplanar but not collinear manner.

Images