Rigid registration CT motion parameter estimation and motion compensation method and system

Abstract

The invention discloses a rigid registration CT (Computed Tomography) motion parameter estimation and motion compensation method and system. The 3D-2D rigid registration CT motion parameter estimation and motion compensation method based on exhaustion search strategy optimization is carried out through six steps. According to the method, aiming at the problems of huge computing power demand and storage occupation caused by an exhaustive search mechanism, the motion parameter interference decoupling method is used for optimizing the exhaustive search strategy, so that the demand of the exhaustive projection search strategy on computing resources is effectively reduced while the motion parameter estimation precision is remarkably improved. Besides, in consideration of motion artifact characterization differences of different stages in the iteration process, the invention provides a time-varying data consistency measurement method which is used for more accurately measuring data consistency between re-projection and measurement projection, so that the accuracy of parameter estimation is further improved.

Description

Technical Field

[0001] The present invention relates to the technical field of 3D-2D rigid registration and computer tomography, and in particular to a method and system for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization. Background Art

[0002] With the continuous advancement of dental technology, cone-beam computed tomography (CBCT) has become an indispensable diagnostic tool in modern dentistry. CBCT accurately captures detailed three-dimensional images of tooth structure, soft tissue, and bone while minimizing radiation exposure. These high-fidelity images provide valuable imaging support for dental implant planning, periodontal diagnosis, maxillofacial surgery, and other procedures. However, the relatively slow scanning speed of CBCT increases the likelihood of subject motion during scanning, posing a significant challenge. High-quality CBCT reconstruction depends on precise acquisition geometry. Subject motion and gantry vibration often occur, leading to misalignment of measured signals during backprojection and resulting in motion artifacts in the reconstructed images. These artifacts manifest as streaks, double-contoured tissue, or blurred structures, severely hindering clinical diagnosis and treatment. Therefore, estimating the motion trajectory and updating the acquisition geometry for image reconstruction are key steps in compensating for subject motion.

[0003] Existing motion estimation and compensation methods can be divided into two main categories: marker-free methods and marker-assisted methods. Marker-free methods mainly include techniques based on projection consistency metrics, 3D-2D registration, and image quality metrics. Methods based on projection consistency metrics fine-tune the acquisition geometry by maximizing the consistency of projection data. 3D-2D registration methods simulate CT forward models, generate virtual projections, and carefully compare them with actual 2D projections to determine the optimal acquisition geometry. On the other hand, methods based on image quality metrics quantitatively evaluate image quality and optimize acquisition parameters to enhance reconstruction quality. Marker-assisted methods rely on external motion tracking systems or fiducial markers to record the motion trajectory of the object, enabling the acquisition geometry parameters to be adjusted during image reconstruction for motion compensation. However, these methods incur additional hardware costs, and fiducial markers may cause local tissue deformation during scanning.

[0004] Therefore, in view of the shortcomings of the existing technology, it is necessary to provide a method and system for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization to solve the shortcomings of the existing technology. Summary of the Invention

[0005] The first object of the present invention is to overcome the shortcomings of the prior art and provide a method for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on an exhaustive search strategy optimization. This method for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on an exhaustive search strategy optimization can improve the accuracy of motion parameter estimation.

[0006] The above-mentioned purpose of the present invention is achieved through the following technical measures:

[0007] A method for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization is provided, which is performed in the following steps:

[0008] S1. Obtain CT image x of rigid motion artifact and corresponding measurement projection y measure and CT Geometry A geo At the same time, the translation distances along the three coordinate axes (w, u, v) in the rotation coordinate system of the CT image x are defined as t u , t w and t v , the rotation angles along the three coordinate axes (w, u, v) in the rotating coordinate system are defined as r u 、r w and r v , initialize the CT geometry A geo At each CT scan angle t u , t w , t v 、r u 、r w and r v , corresponding to 6 initial motion curves, namely t u Motion curve, t w Motion curve, t v Motion curve, r u Motion curve, r w Motion curve and r v Motion curve;

[0009] S2. Update t through exhaustive search strategy v Motion curve;

[0010] S3. Update r through exhaustive search strategy u Motion curve and r v Motion curve;

[0011] S4. Update t through exhaustive search strategy u Motion curve, t w Motion curve and r v Motion curve;

[0012] S5. According to the CT geometry Ageo and updated t u Motion curve, updated t w Motion curve, updated t v Motion curve, updated r u Motion curve, updated r w Motion curve and updated r v The six motion curves update the reconstruction geometry, reconstruct the CT image and assign it as the updated CT image x; then determine whether the iteration stop criterion is met. If not, the updated CT image x is replaced by the updated CT image x, and the updated t u The motion curve replaces the t before the update u Motion curve, updated t w The motion curve replaces the t before the update w Motion curve, updated t v The motion curve replaces the t before the update v Motion curve, updated r u The motion curve replaces the r before the update u Motion curve, updated r w The motion curve replaces the r before the update w Motion curve, updated r v The motion curve replaces the r before the update v Motion curve, then return to S2; if yes, enter S6;

[0013] S6. Use the current CT image x as the final CT image x, and the current t u Motion curve, t w Motion curve, t v Motion curve, r u Motion curve, r w Motion curve and r v The motion curves are used as the final six motion curves.

[0014] In the CT geometry A geo Each CT scanning angle is set to t u , t w , t v 、r u 、r w and r v .

[0015] Preferably, the above-mentioned translation distance and the above-mentioned rotation angle are both motion parameters.

[0016] In the S1, the CT geometry A geo t at each CT scanning angle u , t w , tv 、r u 、r w and r v are all set to 0 to obtain the initial t u Motion curve, initial t w Motion curve, initial t v Motion curve, initial r u Motion curve, initial r w Motion curve and initial r v Motion curve.

[0017] Preferably, the above S2 is specifically: t for each scanning angle v The parameter estimation value update method is: sampling t v Motion parameters, and exhaustively enumerate possible t v Motion parameter example, according to the CT geometry A geo Generate the forward projection geometry corresponding to each instance; then use the generated forward projection geometry to perform a forward projection operation on the CT image x to obtain the reprojection y corresponding to each instance reproj , search for the measurement projection y measure The best similarity corresponds to t v Example, use the tgv instance parameter value to update the t corresponding to the CT scan angle v Parameter estimation value; when t for all scanning angles v After all parameter estimates are updated, the updated t v Motion curve.

[0018] Preferably, the above S3 is specifically: r for each scanning angle example u Parameter estimates and r v The parameter estimation value update method is: sampling by r u and r v The two-dimensional (r u , r w ) motion parameter group, exhaustively enumerate possible (r u , r w ) motion parameter group instance, according to CT geometry A geo Generate the forward projection geometry corresponding to each instance; then use the generated forward projection geometry to perform a forward projection operation on the CT image x to obtain the reprojection y corresponding to each instance reproj , search for the measurement projection y measure The best similarity corresponds to (r u , r w ) instance, use (r u , r w ) instance parameter value updates the corresponding CT scan angle r u Parameter estimates and r vParameter estimation value; when r of all scanning angles u Parameter estimates and r v After all parameter estimates are updated, the updated r u Motion curve and r v Motion curve.

[0019] Preferably, the above S4 is specifically: t for each scanning angle example u Parameter estimates, t w Parameter estimates and r v The parameter estimation value update method is: sampling by tr, t w and r v The three-dimensional (t u , t w , r v ) motion parameter group, exhaustive possible (t u , t w , r v ) Motion parameter group instance, according to CT geometry A geo Generate the forward projection geometry corresponding to each instance; then use the generated forward projection geometry to perform a forward projection operation on the CT image x to obtain the reprojection y corresponding to each instance reproj , search for the reprojection y reproj With the measured projection y measure The best similarity corresponds to (t u , t w , r v ) instance, use (t u , t w , r v ) instance parameter value updates corresponding to the CT scan angle t u Parameter estimates, t w Parameter estimates and r v Parameter estimation value; when t for all scanning angles u Parameter estimates, t w Parameter estimates and r v After all parameter estimates are updated, the updated t u Motion curve, t w Motion curve and r v Motion curve.

[0020] Reprojection y reproj With the measured projection y measure The similarity measurement formulas are shown in formulas (1) to (3):

[0021] Loss = λ j L1+(1-λ j )L2……Formula (1);

[0022] L1=||ymeasure -y reproj ||1……Formula (2);

[0023] L2=||LOG(y measure )-LOG(y reproj )||1……Formula (3);

[0024] Among them, λ is a hyperparameter that changes with the iteration index j, LOG is the Gaussian Laplace operator, L1 and L2 are norms, Loss is the measure of similarity, and Loss is the cost function.

[0025] Preferably, the above CT geometry A geo It includes all parameters that describe the spatial position relationship between the X-ray source, the scanned object and the detector.

[0026] Preferably, the above CT geometry A geo are the source-detector distance (DSD), source-rotation center distance, detector unit size, detector tilt angle, and number of detector arrays.

[0027] In S1, the u-axis of the rotating coordinate system at each scanning angle is aligned with the CT geometry A. geo The vector from the detector pixel (0,0) to (0,1) is aligned; the v axis of the rotating coordinate system at each scan angle follows the direction of the vector from the pixel (0,0) to (1,0) in the CT geometry v; the w axis of the rotating coordinate system at each scan angle is the same as the direction of the normal vector of the detector.

[0028] Preferably, the above-mentioned measurement projection y measure It is the real projection data of the CT scan part of the object.

[0029] In the S2, S3 and S3, according to the CT geometry A geo The forward projection geometry corresponding to each instance is generated by converting the rigid motion of the scanned part of the object into the spatial pose transformation of the source-detector pair.

[0030] Preferably, the above-mentioned front projection geometry is performed by a front projection algorithm, and the front projection algorithm is a ray tracing algorithm or a 3D Gaussian sputtering 3DGS algorithm.

[0031] In S5 , a CT image is reconstructed using a reconstruction algorithm.

[0032] Preferably, the above reconstruction algorithm is an FDK algorithm, an algebraic reconstruction algorithm or a statistical reconstruction algorithm.

[0033] Preferably, the center of the rotating coordinate system coincides with the center of the world coordinate system.

[0034] The second object of the present invention is to overcome the shortcomings of the prior art and provide a system for rigid registration CT motion parameter estimation and motion compensation, which can improve the accuracy of motion parameter estimation.

[0035] The above-mentioned purpose of the present invention is achieved through the following technical measures:

[0036] A system for rigid registration CT motion parameter estimation and motion compensation is provided, characterized by executing the method for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization as described in any one of claims 1 to 9.

[0037] The present invention provides a method and system for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization, wherein the rigid registration CT motion parameter estimation and motion compensation are performed by the following steps: S1, obtaining a CT image x of rigid motion artifacts, and a corresponding measurement projection y measure and CT Geometry A geo At the same time, the translation distances along the three coordinate axes (w, u, v) in the rotation coordinate system of the CT image x are defined as t u , t w and t x , the rotation angles along the three coordinate axes (w, u, v) in the rotating coordinate system are defined as r u 、r w and r v , initialize the CT geometry A geo At each CT scan angle t u , t w , t v 、r u 、r w and r v , corresponding to 6 initial motion curves, namely t u Motion curve, t w Motion curve, t v Motion curve, r u Motion curve, r w Motion curve and r v Motion curve; S2, update t through exhaustive search strategy v Motion curve; S3, update r through exhaustive search strategy u Motion curve and r v Motion curve; S4, update t through exhaustive search strategy u Motion curve, t w Motion curve and r v Motion curve; S5, according to the CT geometry A geo and updated tu Motion curve, updated t w Motion curve, updated t v Motion curve, updated r u Motion curve, updated r w Motion curve and updated r v The six motion curves update the reconstruction geometry, reconstruct the CT image and assign it as the updated CT image x; then determine whether the iteration stop criterion is met. If not, the updated CT image x is replaced by the updated CT image x, and the updated t u The motion curve replaces the t before the update u Motion curve, updated t w The motion curve replaces the t before the update w Motion curve, updated t v The motion curve replaces the t before the update v Motion curve, updated r u The motion curve replaces the r before the update u Motion curve, updated r w The motion curve replaces the r before the update w Motion curve, updated r v The motion curve replaces the r before the update v motion curve, and then return to S2; when yes, enter S6; S6, take the current CT image x as the final CT image x, and take the current t u Motion curve, t w Motion curve, t v Motion curve, r u Motion curve, r w Motion curve and r v The motion curves are used as the final six motion curves. The beneficial effects of the present invention are: 1. The present invention increases the parameter estimation sequence, thereby decoupling the mutual interference between the motion parameters, thereby increasing the accuracy of motion parameter estimation. 2. The present invention customizes a time-varying cost function to accurately measure the similarity between the reprojection and the measured projection under different motion artifact characterizations. 3. The present invention uses an exhaustive search strategy, which can effectively avoid falling into local minima and increase the accuracy of motion parameter estimation. In summary, the present invention addresses the huge computing power requirements and storage occupancy problems brought about by the exhaustive search mechanism, and optimizes the exhaustive search strategy by decoupling the motion parameter interference method. While significantly improving the accuracy of motion parameter estimation, it effectively reduces the demand for computing resources by the exhaustive reprojection search strategy. In addition, the present invention takes into account the differences in motion artifact characterization at different stages of the iterative process, and proposes a time-varying data consistency measurement method for more accurately measuring the data consistency between the reprojection and the measured projection, thereby further improving the accuracy of parameter estimation. BRIEF DESCRIPTION OF THE DRAWINGS

[0038] The present invention is further described with reference to the accompanying drawings, but the contents in the accompanying drawings do not constitute any limitation to the present invention.

[0039] Figure 1 Flowchart of a method for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization.

[0040] Figure 2 Schematic diagram of the world coordinate system and the rotating coordinate system used in the present invention.

[0041] Figure 3 This is a comparison chart of the results before and after the algorithm of the present invention is corrected. DETAILED DESCRIPTION

[0042] The technical solution of the present invention is further described with reference to the following examples.

[0043] Example 1

[0044] A method for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization, such as Figure 1 , proceed as follows:

[0045] S1. Obtain CT image x of rigid motion artifact and corresponding measurement projection y measure and CT Geometry A geo At the same time, the translation distances along the three coordinate axes (w, u, v) in the rotation coordinate system of the CT image x are defined as t u , t w and t v , the rotation angles along the three coordinate axes (w, u, v) in the rotating coordinate system are defined as r u 、r w and r v , initialize CT geometry A geo At each CT scan angle t u , t w , t v 、r u 、r w and r v , corresponding to 6 initial motion curves, namely t u Motion curve, t w Motion curve, t v Motion curve, r u Motion curve, r w Motion curve and r v Motion curve; where CT geometry A is geo t at each CT scanning angle u , tw , t v 、r u 、r w and r v are all set to 0 for initialization; in CT geometry A geo Each CT scanning angle is set to t u , t w , t v 、r u 、r w and r v ; The translation distance and rotation angle are both motion parameters;

[0046] S2. Update t through exhaustive search strategy v Motion curve;

[0047] S3. Update r through exhaustive search strategy u Motion curve and r v Motion curve;

[0048] S4. Update t through exhaustive search strategy u Motion curve, t w Motion curve and r v Motion curve;

[0049] S5, according to CT geometry A geo and updated t u Motion curve, updated t w Motion curve, updated t v Motion curve, updated r u Motion curve, updated r w Motion curve and updated r v The six motion curves update the reconstruction geometry, reconstruct the CT image and assign it as the updated CT image x; then determine whether the iteration stop criterion is met. If not, the updated CT image x is replaced by the updated CT image x, and the updated t u The motion curve replaces the t before the update u Motion curve, updated t w The motion curve replaces the t before the update w Motion curve, updated t v The motion curve replaces the t before the update v Motion curve, updated r u The motion curve replaces the r before the update u Motion curve, updated r w The motion curve replaces the r before the update w Motion curve, updated r v The motion curve replaces the r before the update vMotion curve, then return to S2; if yes, enter S6;

[0050] S6. Use the current CT image x as the final CT image x, and the current t u Motion curve, t w Motion curve, t v Motion curve, r u Motion curve, r w Motion curve and r v The motion curves are used as the final six motion curves.

[0051] Among them, the center of the rotating coordinate system coincides with the center of the world coordinate system, such as Figure 2 shown.

[0052] It should be noted that the rigid motion artifacts of the present invention refer to the image artifacts caused by the rigid motion of the scanned part of the object during the CT scanning process. For the first update, in S5, the updated CT image x of S1 is replaced with the updated one, and the updated corresponding t u Motion curve, t w Motion curve, t v Motion curve, r u Motion curve, r w Motion curve and r v Motion curve, corresponding to the replacement of the initial t u Motion curve, t w Motion curve, t v Motion curve, r u Motion curve, r w Motion curve and r v Motion curve; when it is not the first iteration, the updated CT image x will be used to replace the previous one, and the updated corresponding t u Motion curve, t w Motion curve, t v Motion curve, r u Motion curve, r w Motion curve and r v Motion curve, corresponding to the replacement of the last t u Motion curve, t w Motion curve, t v Motion curve, r u Motion curve, r w Motion curve and r v Motion curve.

[0053] The iteration stopping criterion of the present invention can be that the sum of the differences between the currently updated CT image x and the last updated CT image x is less than a threshold, or can be a maximum number of iterations. The threshold or the maximum number of iterations can be selected according to actual conditions. For example, the threshold is 1, 2, 3, 5, 10, etc., and the maximum number of iterations is 5, 10, 20, etc. In this embodiment, the threshold is 5, and the maximum number of iterations is 10.

[0054] Where S2 is specifically: t for each scanning angle v The parameter estimation value update method is: sampling t v Motion parameters, and exhaustively enumerate possible t v Motion parameter example, according to CT geometry A geo Generate the forward projection geometry corresponding to each instance; then use the generated forward projection geometry to perform a forward projection operation on the CT image x to obtain the reprojection y corresponding to each instance reproj , search and measure projection y measure The best correspondence t of similarity v Example, using t v The instance parameter value updates the corresponding CT scan angle t v Parameter estimation value; when t for all scanning angles v After all parameter estimates are updated, the updated t v Motion curve.

[0055] Among them, S3 is specifically: r for each scanning angle example u Parameter estimates and r v The parameter estimation value update method is: sampling by r u and r v The two-dimensional (r u , r w ) motion parameter group, exhaustively enumerate possible (r u , r w ) motion parameter group instance, according to CT geometry A geo Generate the forward projection geometry corresponding to each instance; then use the generated forward projection geometry to perform a forward projection operation on the CT image x to obtain the reprojection y corresponding to each instance reproj , search and measure projection y measure The best similarity corresponds to (r u , r w ) instance, use (r u , r w ) instance parameter value updates the corresponding CT scan angle r u Parameter estimates and r v Parameter estimation value; when r of all scanning angles u Parameter estimates and r v After all parameter estimates are updated, the updated ru Motion curve and r v Motion curve.

[0056] Where S4 is specifically: t for each scanning angle u Parameter estimates, t w Parameter estimates and r v The parameter estimation value update method is: sampling by t u , t w and r v The three-dimensional (t u , t w , r v ) motion parameter group, exhaustive possible (t u , t w , r v ) Motion parameter group instance, according to CT geometry A geo Generate the forward projection geometry corresponding to each instance; then use the generated forward projection geometry to perform a forward projection operation on the CT image x to obtain the reprojection y corresponding to each instance reproj , search for the reprojection y reproj With the measured projection y measure The best similarity correspondence (t u , t w , r v ) instance, use (t u , t w , r v ) instance parameter value updates corresponding to the CT scan angle t u Parameter estimates, t w Parameter estimates and r v Parameter estimation value; when t for all scanning angles u Parameter estimates, t w Parameter estimates and r v After all parameter estimates are updated, the updated t u Motion curve, t w Motion curve and r v Motion curve.

[0057] The present invention's reprojection y in S2, S3 and S2 reproj With the measured projection y measure The similarity measurement formula is shown in formula (1) to formula (3):

[0058] Loss = λ j L1+(1-λ j )L2……Formula (1);

[0059] L1=||y measure -y reproj ||1……Formula (2);

[0060] L2=||LOG(y measure )-LOG(y reproj )||1……Formula (3);

[0061] Among them, λ is a hyperparameter that changes with the iteration index j, LOG is the Gaussian Laplace operator, L1 and L2 are norms, Loss is the measure of similarity, and Loss is the cost function.

[0062] CT geometry A of the present invention geo It includes all parameters that describe the spatial relationship between the X-ray source, the scanned object, and the detector. geo are the source-detector distance (DSD), source-rotation center distance, detector unit size, detector tilt angle, and number of detector arrays.

[0063] In S1, the u-axis of the rotation coordinate system at each scanning angle is aligned with the coordinate system of the CT geometry A. geo The vector alignment of the detector pixel (0,0) to (0,1) in the CT image is as follows: the v-axis of the rotation coordinate system at each scan angle follows the vector alignment from the CT geometry A geo The vector direction from pixel (0,0) to (1,0) in the image is the same as the normal vector of the detector.

[0064] where the measured projection y measure It is the real projection data of the CT scan part of the object.

[0065] The present invention is based on CT geometry A in S2, S3 and S2 geo The front projection geometry corresponding to each instance is generated by converting the rigid motion of the scanned part of the object into the spatial pose transformation of the source-detector pair; the front projection geometry is performed by a front projection algorithm, which is a ray tracing algorithm or a 3D Gaussian sputtering 3DGS algorithm.

[0066] In S5 , a CT image is reconstructed by a reconstruction algorithm; and the reconstruction algorithm is an FDK algorithm, an algebraic reconstruction algorithm, or a statistical reconstruction algorithm.

[0067] The beneficial effects of the method for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization are as follows: 1. The present invention increases the parameter estimation sequence, thereby decoupling the mutual interference between motion parameters, thereby increasing the accuracy of motion parameter estimation. 2. The present invention customizes a time-varying cost function to accurately measure the similarity between reprojection and measured projection under different motion artifact characterization conditions. 3. The present invention uses an exhaustive search strategy, which can effectively avoid falling into local minima and increase the accuracy of motion parameter estimation. In summary, the present invention addresses the huge computing power requirements and storage occupancy problems brought about by the exhaustive search mechanism, and optimizes the exhaustive search strategy by decoupling the motion parameter interference method. While significantly improving the accuracy of motion parameter estimation, it effectively reduces the demand for computing resources for the exhaustive reprojection search strategy. In addition, the present invention takes into account the differences in motion artifact characterization at different stages of the iterative process, and proposes a time-varying data consistency measurement method for more accurately measuring the data consistency between reprojection and measured projection, thereby further improving the accuracy of parameter estimation.

[0068] Example 2

[0069] A method for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization as in Example 1 is applied. In this example, a CT image x0 of a subject's head scanned by a Fussen commercial dental CBCT (Fussen Matrix5000) is used to simulate the rigid motion of the subject's head during the CT scanning process. Specifically, t is randomly generated at each CT scanning angle. u , t w , t v 、r u 、r w and r v Motion amplitude value, combined with CT geometry A geo Generate the corresponding front projection geometry and perform a front projection operation on the CT image x0 to simulate the measured projection y obtained when the subject's head undergoes rigid motion during the CT scan. measure Without considering the motion situation, directly use CT geometry A geo Reconstruct the CT image x with motion artifacts. measure and CT Geometry A geo The data is input into a CT rigid motion artifact correction algorithm optimized based on an exhaustive search strategy, and the algorithm is run to perform rigid motion parameter estimation and motion compensation, ultimately obtaining a corrected CT image and an estimated rigid motion curve.

[0070] Specifically, CT geometry A geoThe settings are: the image size of the CT image is 800×800×400, the source-to-detector distance (DSD) is set to 635 mm, the source-to-rotation center distance (DSO) is set to 425 mm, the detector unit size is 0.2 mm×0.2 mm, the detector array is 500×1500, the reconstructed image pixel size is 0.2 mm×0.2 mm×0.2 mm, and the number of projections is 300.

[0071] Specific u , t w , t v 、r u 、r w and r v The motion parameter search range is set to [-3 mm, 3 mm] and [-3°, 3°], and the motion parameter search step is 0.2.

[0072] pass Figure 3 This paper demonstrates the effectiveness of a CT rigid motion artifact correction algorithm optimized using an exhaustive search strategy. CT images with rigid motion artifacts exhibit significant double contouring and anatomical distortion. The present invention effectively recovers details such as teeth and bone tissue from severely distorted CT images with rigid motion artifacts, eliminating double contouring artifacts. This result validates the feasibility and effectiveness of the proposed CT rigid motion artifact correction method based on an exhaustive search strategy.

[0073] Example 3

[0074] A system for rigid registration CT motion parameter estimation and motion compensation, wherein the method for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization of embodiment 1 is implemented.

[0075] The beneficial effects of the rigid registration CT motion parameter estimation and motion compensation system are as follows: 1. The present invention increases the parameter estimation sequence, thereby decoupling the mutual interference between motion parameters, thereby increasing the accuracy of motion parameter estimation. 2. The present invention customizes a time-varying cost function to accurately measure the similarity between the reprojection and the measured projection under different motion artifact characterization conditions. 3. The present invention uses an exhaustive search strategy, which can effectively avoid falling into local minima and increase the accuracy of motion parameter estimation. In summary, the present invention addresses the huge computing power requirements and storage occupancy problems brought about by the exhaustive search mechanism, and optimizes the exhaustive search strategy by decoupling the motion parameter interference method. While significantly improving the accuracy of motion parameter estimation, it effectively reduces the demand for computing resources for the exhaustive reprojection search strategy. In addition, the present invention takes into account the differences in motion artifact characterization at different stages of the iterative process, and proposes a time-varying data consistency measurement method for more accurately measuring the data consistency between the reprojection and the measured projection, thereby further improving the accuracy of parameter estimation.

[0076] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and are not intended to limit the scope of protection of the present invention. Although the present invention has been described in detail with reference to 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 essence and scope of the technical solutions of the present invention.

Claims

1. A method for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization, characterized in that: Proceed as follows: S1. Obtain CT image x of rigid motion artifact and corresponding measurement projection y measure and CT Geometry A geo At the same time, the translation distances along the three coordinate axes (w, u, v) in the rotation coordinate system of the CT image x are defined as t u , t w and t v , the rotation angles along the three coordinate axes (w, u, v) in the rotating coordinate system are defined as r u 、r w and r v , initialize the CT geometry A geo At each CT scan angle t u , t w , t v 、r u 、r w and r v , corresponding to 6 initial motion curves, namely t u Motion curve, t w Motion curve, t v Motion curve, r u Motion curve, r w Motion curve and r v Motion curve; S2. Update t through exhaustive search strategy v Motion curve; S3. Update r through exhaustive search strategy u Motion curve and r v Motion curve; S4. Update t through exhaustive search strategy u Motion curve, t w Motion curve and r v Motion curve; S5. According to the CT geometry A geo and updated t u Motion curve, updated t w Motion curve, updated t v Motion curve, updated r u Motion curve, updated r w Motion curve and updated r v The six motion curves update the reconstruction geometry, reconstruct the CT image and assign it as the updated CT image x; then determine whether the iteration stop criterion is met. If not, the updated CT image x is replaced by the updated CT image x, and the updated t u The motion curve replaces the t before the update u Motion curve, updated t w The motion curve replaces the t before the update w Motion curve, updated t v The motion curve replaces the t before the update v Motion curve, updated r u The motion curve replaces the r before the update u Motion curve, updated r w The motion curve replaces the r before the update w Motion curve, updated r v The motion curve replaces the r before the update v Motion curve, then return to S2; if yes, enter S6; S6. Use the current CT image x as the final CT image x, and the current t u Motion curve, t w Motion curve, t v Motion curve, r u Motion curve, r w Motion curve and r v The motion curves are used as the final six motion curves.

2. The method for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization according to claim 1, characterized in that: In the CT geometry A geo Each CT scanning angle is set to t u , t w , t v 、r u 、r w and r v ; The translation distance and the rotation angle are both motion parameters; In the S1, the CT geometry A geo t at each CT scanning angle u , t w , t v 、r u 、r w and r v are all set to 0 to obtain the initial t u Motion curve, initial t w Motion curve, initial t v Motion curve, initial r u Motion curve, initial r w Motion curve and initial r v Motion curve.

3. The method for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization according to claim 1, characterized in that: The S2 is specifically: t for each scanning angle v The parameter estimation value update method is: sampling t v Motion parameters, and exhaustively enumerate possible t v Motion parameter example, according to the CT geometry A geo Generate the forward projection geometry corresponding to each instance; then use the generated forward projection geometry to perform a forward projection operation on the CT image x to obtain the reprojection y corresponding to each instance reproj , search for the measurement projection y measure The best correspondence t of similarity v Example, using t v The instance parameter value updates the corresponding CT scan angle t v Parameter estimation value; when t for all scanning angles v After all parameter estimates are updated, the updated t v Motion curve.

4. The method for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization according to claim 1, characterized in that: The S3 is specifically: r for each scanning angle example u Parameter estimates and r v The parameter estimation value update method is: sampling by r u and r v The two-dimensional (r u , r w ) motion parameter group, exhaustively enumerate possible (r u , r w ) motion parameter group instance, according to CT geometry A geo Generate the forward projection geometry corresponding to each instance; then use the generated forward projection geometry to perform a forward projection operation on the CT image x to obtain the reprojection y corresponding to each instance reproj , search for the measurement projection y measure The best similarity corresponds to (r u , r w ) instance, use (r u , r w ) instance parameter value updates the corresponding CT scan angle r u Parameter estimates and r v Parameter estimation value; when r of all scanning angles u Parameter estimates and r v After all parameter estimates are updated, the updated r u Motion curve and r v Motion curve.

5. The method for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization according to claim 1, characterized in that: The S4 is specifically: t for each scanning angle example u Parameter estimates, t w Parameter estimates and r v The parameter estimation value update method is: sampling by t u , t w and r v The three-dimensional (t u , t w , r v ) motion parameter group, exhaustive possible (t u , t w , r v ) Motion parameter group instance, according to CT geometry A geo Generate the forward projection geometry corresponding to each instance; then use the generated forward projection geometry to perform a forward projection operation on the CT image x to obtain the reprojection y corresponding to each instance reproj , search for the reprojection y reproj With the measured projection y measure The best similarity correspondence (t u , t w , r v ) instance, use (t u , t w , r v ) instance parameter value updates corresponding to the CT scan angle t u Parameter estimates, t w Parameter estimates and r v Parameter estimation value; when t for all scanning angles u Parameter estimates, t w Parameter estimates and r v After all parameter estimates are updated, the updated t u Motion curve, t w Motion curve and r v Motion curve.

6. The method for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization according to any one of claims 3 to 5, characterized in that: Reprojection y reproj With the measured projection y measure The similarity measurement formula is shown in formula (1) to formula (3): Loss=λ j L1+(1 - λ j )L2……Equation (1); L1 = ||y measure -y reproj ||1……Equation (2); L2 = ||LOG(y measure ) - LOG(y reproj )||1……Equation (3); Among them, λ is a hyperparameter that changes with the iteration index j, LOG is the Gaussian Laplace operator, L1 and L2 are norms, Loss is the measure of similarity, and Loss is the cost function.

7. The method for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization according to any one of claims 1 to 5, characterized in that: The CT geometry A geo It includes all parameters that describe the spatial position relationship between the X-ray source, the scanned object and the detector.

8. The method for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization according to claim 7, characterized in that: The CT geometry A geo are the source-detector distance (DSD), source-rotation center distance, detector unit size, detector tilt angle, and number of detector arrays; In S1, the u-axis of the rotating coordinate system at each scanning angle is aligned with the CT geometry A. geo The vector alignment of the detector pixel (0,0) to (0,1) in the CT image is as follows: the v axis of the rotating coordinate system at each scan angle follows the vector alignment of the detector pixel (0,0) to (0,1) in the CT image; ... geo The vector direction of the pixel (0,0) to (1,0) in the image; the w axis of the rotating coordinate system at each scanning angle is the same as the normal vector direction of the detector; The measured projection y measure It is the real projection data of the CT scan part of the object.

9. The method for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization according to any one of claims 3 to 5, characterized in that: According to CT geometry A geo Generate the front projection geometry corresponding to each instance to convert the rigid motion of the scanned part of the object into the spatial pose transformation of the source-detector pair; The front projection geometry is performed by a front projection algorithm, which is a ray tracing algorithm or a 3D Gaussian sputtering 3DGS algorithm; In S5, a CT image is reconstructed by a reconstruction algorithm; and the reconstruction algorithm is an FDK algorithm, an algebraic reconstruction algorithm, or a statistical reconstruction algorithm; The center of the rotating coordinate system coincides with the center of the world coordinate system.

10. A system for rigid registration CT motion parameter estimation and motion compensation, characterized by: Execute the method for 3D-2D rigid registration CT motion parameter estimation and motion compensation based on exhaustive search strategy optimization as described in any one of claims 1 to 9.

Images