Cone beam CT geometry correction method suitable for irregular trajectories

Abstract

The application discloses a cone beam CT geometry correction method suitable for irregular tracks, adopts a new optimization constraint iteration method based on a back projection model to solve geometry parameters under each projection view, and most of existing correction methods all assume that the internal geometry relation of a reference point in a phantom is accurate and known, and rarely deeply study whether the phantom precision is suitable for the correction model. Through the disclosure of the application scheme, the correction method is a high-noise robust correction method, and high-precision correction can be realized even when the resolution of a system to be calibrated is twice that of a phantom measurement system.

Description

Technical Field

[0001] This invention relates to the field of image processing technology in computed tomography, and more specifically to a geometric correction method for high-resolution cone-beam computed tomography with irregular scanning trajectories. Background Technology

[0002] Cone-beam CT is currently widely used in many fields, such as medical diagnosis and treatment, preclinical imaging and research, and structural analysis in materials science. Typically, different applications have different resolution requirements for cone-beam CT. The resolutions for clinical, preclinical, and materials science cone-beam CT are 100 micrometers, 10 micrometers, and 1 micrometer, respectively. With increasing resolution, cone-beam CT systems are expected to become more stable and accurate. Therefore, the corresponding geometric correction steps, a prerequisite for obtaining ideal CT images, become more challenging.

[0003] To address calibration problems encountered in various scenarios, numerous calibration methods have been proposed in the literature over the past few decades. One type is online calibration, which does not rely on a specific phantom but rather on loss functions based on different image quality metrics such as image entropy, image gradient, or a combination of multiple metrics. This type of method is suitable for acquisitions of regular but non-repeatable geometries. These methods are typically computationally intensive and prone to degrading image quality. Most commercial CT systems are mechanically repeatable and employ another method called offline calibration. Offline calibration methods are applicable to both regular and irregular geometries. For regular geometric acquisition systems, a simple phantom, with very few parameters, can be used by embedding one or a few small spheres to obtain all geometric parameters. In these calibration cases, the model is not complex, manual operation is easy, and the corresponding methods are also very simple. For irregular systems where the geometric parameters differ across views, a simple phantom with low accuracy requirements cannot work. This is because, unless additional conditions are given, such as the spatial relationships between point or line markers, the complete projection matrix cannot be solved. A search revealed that most literature assumes that the geometric relationships between markers are known. Some of these methods rely on specially fabricated phantoms in which reference points are precisely combined in a specific pattern, such as a spiral, circle, or orthogonal line. The geometry of the specific pattern is considered known, so manufacturing accuracy is a key prerequisite for obtaining accurate geometric parameters. Because manufacturing errors are not easily quantifiable, previous methods have not paid much attention to this issue. In fact, if the system does not require high spatial resolution, fabricating a dedicated phantom with low accuracy requirements may not be too difficult. However, as resolution increases, this becomes a challenging task. Many researchers may encounter situations where they cannot reproduce known calibration methods or achieve the same calibration accuracy based on self-made phantoms. One reason is that the accuracy of the phantom does not meet the requirements of the system being calibrated, but unfortunately, due to the lack of suitable detection methods, it will ultimately be impossible to know whether this is truly the case. Measuring a roughly fabricated phantom seems more practical than fabricating a high-precision phantom; therefore, Li et al. proposed a calibration method based on projection matrices. However, they did not provide detailed measurement methods or further investigate the degree of matching between the measured model and the calibrated model. In theory, high-precision phantoms and robust calibration methods are mutually restrictive and supportive; both are crucial for calibration accuracy. Therefore, it is necessary to investigate the conditions under which the calibration model and phantom measurements can be matched to achieve optimal calibration accuracy for irregular trajectory CBCT. Furthermore, since the reference points are not arranged regularly, a calibration model based on the projection matrix is ​​currently the best choice. In practical applications, however, this model has limited noise robustness, especially when the number of reference points is low. This necessitates a more reliable calibration method to compensate for this deficiency.Therefore, with good phantom measurement methods and more reliable calibration models, we can use relatively rough, handcrafted phantoms to calibrate high-resolution CT systems with irregular trajectories. This is extremely helpful for researchers who want to build high-resolution cone-beam CT systems but lack advanced machining equipment to manufacture high-precision phantoms. Summary of the Invention

[0004] Purpose of the Invention: The purpose of this invention is to address the shortcomings of existing technologies by providing a geometric correction method for cone-beam CT with irregular trajectories. First, a measurement model is introduced to obtain the intrinsic geometric relationship of the reference points in the phantom. The measurement model is implemented by extending the previous high-precision geometric correction model, which is suitable for cone-beam CT with circular trajectories. Then, a new optimization constraint iteration method based on a back-projection model is proposed to solve the geometric parameters under each projection view.

[0005] Technical solution: This invention provides a cone-beam CT geometric correction method suitable for irregular trajectories, comprising the following steps:

[0006] S1: After obtaining the spatial coordinates of the randomly arranged balls through calculation, an iterative model is established to achieve geometric correction of a single projected view;

[0007] S2: There are two coordinate systems in the model, namely and ;in This is called the detector coordinate system, where U and V are along the detector's horizontal and vertical axes, respectively, and W is unique to the detector plane. And is called the object coordinate system;

[0008] S3: Assume point B is in Coordinates in a coordinate system are Point B is located at... The corresponding coordinates in the coordinate system can be represented as

[0009]

[0010] In the formula: T is the translation vector.

[0011] R is the rotation matrix.

[0012] R is written as R = R1(η)R2(θ)R3(φ), where η, θ and φ are Euler angles; R1, R2 and R3 represent yaw, pitch and roll angles, respectively;

[0013] S4: In the model, the S-coordinate of the x-ray source is S(w s , u s , v s ) T Point Bw The coordinates of the projection point Q are Q = (0, u, v). T There are k points in the target coordinate system, called B. k K∈[1,k], their coordinates are given by b k = (x k , y k , z k ) T The corresponding projection point in the detector coordinate system is Q. k Its coordinates are Q k = (0, u k , v k ) T The intermediate solution of the geometric parameters in the iteration is (w) s , u s , v s , t1, t2, t3, η, φ, θ), B is obtained through equation (1) k Convert to the detector coordinate system to obtain its corresponding point B. w-k Its coordinates are written as B w-k = [R|T] B k ;

[0014] S5: Let point B be... k to line Q k The foot of the perpendicular line from S is H. k By establishing a system consisting of two straight lines Q k S and B k H k The system of equations formed yields H k coordinates h k Then deduce B k to line Q k The distance to S is as follows:

[0015]

[0016] Combining all K points, the total distance can be expressed as:

[0017]

[0018] When all geometric parameters are accurate, the corresponding =0; therefore, the basic objective function can be written as:

[0019] (4).

[0020] Furthermore, the aforementioned cone-beam CT geometric correction method applicable to irregular trajectories introduces a variance penalty term in the aforementioned equation (4):

[0021]

[0022] in With this term added, the objective function becomes:

[0023]

[0024] In the formula, λ is a weighting factor with an empirical value of 10.

[0025] Furthermore, the aforementioned geometric correction method for cone-beam CT applicable to irregular trajectories is applicable in most cone-beam CT systems where the X-ray source and detector remain relatively stationary during the scanning process; by solving equation (6), the parameter vector for each projection view can be obtained. Here, N is the total number of projections; if there is no noise, all They should be the same; therefore, all average This will be very close to the true coordinates of the X-ray source; typically, a system has hundreds of projected views that need to be corrected, thus ensuring... The accuracy; after determining the source location, the objective function with 9 variables is simplified to a function with 6 variables, as shown below:

[0026] (7).

[0027] Beneficial Effects: Compared with existing technologies, the advantages of this invention are as follows: This invention proposes an easily implemented geometric correction scheme that is independent of phantoms with complexly arranged reference points and exhibits good robustness to noise and measurement errors. The phantoms used in this embodiment are handmade, thus requiring almost no processing equipment or precision. The precise coordinates of the spheres are obtained using the NOC method proposed in previously published schemes, serving as the benchmark for subsequent irregular trajectory correction. Simulations in the embodiment show that the overall coordinate error of each sphere is on the order of 0.01 pixels. When Gaussian noise with a standard deviation of 0.5 is added, the overall error remains at the same magnitude. This result demonstrates that the NOC method can lay a good foundation for the next step of irregular trajectory correction. Simulation verification shows that the simulated annealing algorithm can effectively solve for the global optimum even when the initial values ​​are far from the true values. This means that solving this non-convex model is not a problem in practical applications. Compared with analytical methods, iterative methods are generally more flexible and stable. Therefore, the flexibility of iterative methods can be fully utilized in research to improve robustness to noise and measurement errors. Attached Figure Description

[0028] Figure 1 It is the geometry of the projected view of this invention;

[0029] Figure 2This is the back projection model of the present invention;

[0030] Figure 3 This is the convergence curve of the iterative process in Example 2;

[0031] Figure 4 (a) is the projection metric. (b) is the distribution of the proposed method and the back projection metric based on the GGC method and the proposed method. Distribution;

[0032] Figure 5 These are corrected phantom reconstruction images based on geometric parameters obtained using the GGC method (a1-a4) and the method presented in this paper (b1-b4);

[0033] Figure 6 The images are real images and their MTF curves obtained based on the GGC correction method and the correction method presented in this paper.

[0034] Figure 7 The images are reconstructed based on geometric parameters obtained from the competing GGC method (a1-a4) and the method in this paper (b1-b4);

[0035] Figure 8 This is the first low-resolution mouse femur reconstruction image in Example 4;

[0036] Figure 9 This is the second high-resolution mouse femur reconstruction image in Example 4. Detailed Implementation

[0037] The technical solution of the present invention will be described in detail below with reference to the accompanying drawings, but the scope of protection of the present invention is not limited to the embodiments described.

[0038] First, after obtaining the spatial coordinates of the randomly arranged spheres, an iterative model is established to achieve geometric correction of a single projected view. The geometric structure of the projected view is as follows: Figure 1 As shown, there are two coordinate systems, namely... and . This is called the detector coordinate system. U and V are along the detector's horizontal and vertical axes, respectively. W is unique to the detector plane. This is called the object's coordinate system. Assume point B is in... Coordinates in a coordinate system are This point is at The corresponding coordinates in the coordinate system can be represented as

[0039] .

[0040] In the formula, T is the translation vector, and R is the rotation matrix. R can be written as R = R1(η)R2(θ)R3(φ), where η, θ, and φ are Euler angles, representing yaw, pitch, and roll angles, respectively. Let the coordinates of the X-ray source S be S(w s , u s , v s ) T B w The coordinates of the projection point Q are Q = (0, u, v). T If we define all 9 geometric parameters (w) of this geometry s , u s , v s Given that t1, t2, t3, η, φ, θ are all known, then point B... w S and Q are collinear. Based on this observation, the following back-projection model is established: Figure 2 As shown, these 9 geometric parameters are obtained through iteration. Assume there are k points in the target coordinate system, denoted as B. k K∈[1,K], their coordinates are given by b k = (x k , y k , z k ) T The corresponding projection point in the detector coordinate system is Q. k Its coordinates are Q k =(0, u k , v k ) T If the intermediate solution of the geometric parameters in the iteration is (w) s , u s , v s , t1, t2, t3, η, φ, θ), B is obtained through equation (1) k Convert to the detector coordinate system to obtain its corresponding point B. w-k Its coordinates can be written as B w-k = [R|T] B k Let point B be... k to line Q k The foot of the perpendicular line from S is H. k By establishing a system consisting of two straight lines Q k S and B k H k The system of equations formed by these equations can easily yield H. k coordinates h k Derive B k to line Q k The distance to S is as follows:

[0041]

[0042] Combining all K points, the total distance can be expressed as...

[0043]

[0044] Theoretically, when all nine geometric parameters are accurate, the corresponding = 0. Therefore, the basic objective function can be written as:

[0045]

[0046] Iterative methods are more flexible than analytical methods because various penalty functions or constraints can be added to improve convergence and reduce the impact of metric errors on the final solution. In the model, a variance penalty term is first introduced.

[0047]

[0048] in After adding this term, the objective function becomes...

[0049]

[0050] Here, λ is a weighting factor with an empirical value of 10. The variance penalty is a constraint that makes all... Uniform convergence. When measurement error exists, the objective function value of the optimal parameters is no longer minimized. Therefore, even if the global minimum of the objective function can be obtained, the optimal geometric parameters may not be obtained. In order to reduce the influence of measurement error, traditional analysis methods introduce more marker points. According to the experience of papers in the field of computer vision, about 28 marker points are needed. A famous correction method proposed by Li et al. is to complete the correction of the tomographic synthesis system with 44 points. However, in high-resolution cone-beam CT, due to the limitation of the field of view, it is not easy to fit so many marker points in a small model. Therefore, a new strategy based on the flexibility of the iterative method is introduced instead of increasing the number of marker points to reduce the influence of measurement error. The flexibility of the iterative method is that if any variable in the function is known, then these variables can remain unchanged during the iteration process. When correcting this system, the gantry of the X-ray source and detector can be rotated synchronously by moving / rotating the object or rotating, and the relative position of the X-ray source and detector can be kept unchanged. For example, in most cone-beam CT systems, the X-ray source and detector remain relatively stationary during the scanning process. Therefore, by solving equation (6), the parameter vector of each projection view can be obtained. Here, N is the total number of projected views. Theoretically, if there is no noise, all They should be the same. Therefore, all average This will be very close to the true coordinates of the X-ray source. Typically, a system has hundreds of projected views that need to be corrected, thus ensuring... The accuracy. After determining the source location, the objective function with 9 variables is simplified to a function with 6 variables, as shown below:

[0051]

[0052] Therefore, an approximate solution for s is used to improve the accuracy of the remaining variables. (Function) and Since these are all non-convex functions, the simulated annealing algorithm can be used to solve them. The simulated annealing algorithm combines the simplex simulated annealing method with the traditional simulated annealing method, which can significantly avoid getting trapped in local minima and improve iterative convergence.

[0053] The specific experiment in this embodiment is as follows:

[0054] Example 1. Numerical simulation of phantom measurement

[0055] To verify the feasibility and robustness of the measurement method, two simulation experiments were designed in this embodiment. The scanning trajectory of the cone-beam CT in the simulation was circular. Geometric parameters are shown in Table 1. The simulated detector had 2016×1344 pixels, with a pixel size of 10µm. Since the ideal coordinates of the spheres obtained through the measurement model are the original coordinates, multiplied by the magnification of the measurement system, to intuitively verify the accuracy of the measurement model, both SRD and SDD were set to 8 cm, thus the magnification was 1.0. The simulated model was a cylinder with eight spheres arranged in a spiral around its circumference. In the first simulation, noise-free analytical projection was performed on the spheres to obtain their projection data. A total of 360 projection data points within 360° were obtained. The following metrics were used to evaluate the accuracy of the estimated spatial coordinates of the spheres:

[0056]

[0057] Where K is the number of balls. This is the distance between the actual coordinates and the estimated coordinates of the i-th ball. In the second simulation, the geometric parameters are the same as in the first simulation. The diameter of the ball is set to 400 μm to simulate a realistic projection environment. In addition, 360 projections are collected within a 360° range.

[0058] Table 1: Geometric parameters of the simulation experiment

[0059]

[0060] Table 2: True and measured 3D coordinates of the ball in the noise-free experiment. X, Y, and Z represent the true values ​​in the three directions, respectively. estY est Z est These are the measured values ​​in three directions. The unit is pixels.

[0061]

[0062] The centroid of the projected sphere can be extracted using a threshold segmentation method. To verify the robustness of the measurement method, Gaussian noise with a standard deviation of 0.5 pixels was added to the projected centroid. The entire simulation process was repeated eight times. By randomly moving the model, the final three-dimensional spatial coordinates were obtained using a publicly available measurement method. The results of the first simulation are shown in Table 2, where the first three rows are the true values ​​of the sphere's spatial coordinates in the three directions. The units in the table are pixels. Rows four through six are the estimated coordinates. It can be seen that the estimated spatial coordinates of the sphere are very accurate, with errors in all three directions less than 0.01 pixels, and an overall error precision (DE) of 0.013 pixels. The results of the second simulation experiment containing 50% Gaussian noise are shown in Table 3. The first three rows are the estimated coordinates of the sphere in the three directions. The last three rows are the errors of the sphere. It can be seen that most of the errors in the three directions are much smaller than 0.1 pixels, with an overall error precision (DE) of 0.039 pixels. Based on these experimental results, it can be seen that the phantom measurement method has good robustness to noise and has great potential for practical application.

[0063] Table 3: Three-dimensional coordinates of the ball measured in the noise experiment and the corresponding errors. X est Y est Z est These are the measured values ​​in three directions. X error Y error and Z error These represent the errors in three directions, expressed in pixels.

[0064]

[0065] Example 2. Convergence of the Iterative Method

[0066] To verify the convergence of the proposed geometric correction method, a simulation experiment was conducted using the aforementioned phantom and arbitrarily setting a geometric parameter and its value range. The simulated detector was 2016×1344 pixels in size, with a pixel pitch of 10µm. The geometric relationship between the phantom coordinate system XYZ and the detector coordinate system WUV under a specific projection view is shown in the first row of Table 4. Using the standard three-dimensional spatial coordinates of the eight spheres as known conditions, a non-convex objective function was constructed, i.e., equation (6), and then solved iteratively using the simulated annealing algorithm. The simulation experiment was only to verify the convergence performance of the iterative algorithm, so measurement error was not considered in this simulation. The initial value set for the experiment is located in the second row of Table 4. To prove the convergence performance of the algorithm, this value was set to be sufficiently different from the true value. The errors of the nine geometric parameters obtained are shown in the last row of Table 4. It can be seen that the error is almost equal to 0. The convergence curve of the iterative process is shown in the figure. Figure 3 As shown, the simulated annealing algorithm can effectively converge the objective function to the global minimum even when the initial value is far from the true value.

[0067] Table 4: Actual values, initial values, and errors of the geometric parameters of the scanned view

[0068]

[0069] Example 3. Geometric Correction of a Cone-Beam CT System Simulating Non-Circular Trajectories

[0070] To verify the robustness of the proposed iterative correction method, this embodiment simulates a CBCT system with an elliptical trajectory. The system has an SDD of 8 cm, a detector with 2016 × 1344 pixels, and a pixel pitch of 20 µm. The equation of the elliptical trajectory is as follows: η and φ are 1° and 3° respectively, and the rotation angle of each projection view is 2.0 + k degrees, k∈[0,359) (k is an integer). The phantom used is the one from the simulation experiment above, containing 8 small balls with a diameter of 400 μm. The true coordinates of the small balls are shown in the first three rows of Table 2. In this simulation experiment, the coordinates measured in Table 2 (from row 4 to row 6) are used for geometric correction. To further verify the noise robustness of the proposed correction method, Gaussian noise with a standard deviation of 0.5 is added to the obtained measured coordinates. To quantitatively evaluate the correction accuracy of the algorithm, we simulated a test phantom containing 360 uniformly distributed marker points and performed a simulated scan of the phantom. In addition, we introduced two indices: average forward projection error and average back projection error. The average forward projection error is expressed as follows:

[0071]

[0072] in, Let xi be the ideal projection point of spatial point xi on the m-th projection view. Estimate the projection points, where K is the number of test points and m is the number of projected views. The average backprojection error is:

[0073]

[0074] here, Test point The back projection estimation. The following loss function can be obtained:

[0075]

[0076] in It is the back-projection ray from spatial point x to the m-th projected view, passing through the i-th test point. The Euclidean distance. (Measurement) It reflects the degree of divergence of the back-projected light rays, and its value directly affects the quality of the reconstructed image.

[0077] This embodiment uses the well-known analytical correction algorithm GGC as a benchmark for comparison, illustrating that significant noise reduces correction accuracy and consequently degrades the quality of the reconstructed image, while the proposed method mitigates the impact of noise. Two sets of geometric parameters, obtained using both the benchmark method and the proposed method, are used to scan and reconstruct a calibrated phantom and a volumetric dataset. The reconstruction method employed is the general-purpose FDK algorithm provided by ASTRA Toolbox. The resolution of the reconstructed phantom image is quantitatively analyzed using the MTF (modulation transfer function). The two methods are compared... and The distribution is as follows Figure 4 As shown in the figure, the forward projection error and back projection error of the proposed method can be observed. and The values ​​are highly concentrated, and the average forward projection error and average back projection error obtained by the method in this application are much smaller than those obtained by the comparative method. These quantitative analysis results demonstrate that the proposed iterative correction model can better suppress noise and measurement errors, and obtain more accurate geometric parameters. Since the back projection error is directly related to the quality of the reconstructed image, it can be inferred that the reconstructed image obtained based on the geometric parameters obtained by this method has better quality and higher resolution. The reconstructed images of the corrected phantom obtained by the two methods are as follows: Figure 5As shown, the first to third columns of images are axial images, and the fourth column is a sagittal image. The second and fourth rows are magnified images of the block-marked regions in the first and third rows, respectively. The image in row (a) is obtained based on the geometric parameters obtained using the contrast method, and the image in row (b) is the reconstruction result based on the geometric parameters obtained using the proposed method. Clearly, the edges of the image in row (b) are sharper. Since the tomographic images of the spheres in the phantom resemble planar disks, the MTF curve can be easily obtained by extracting the boundaries of these disks and fitting the edge spread function (ESF) using the Gaussian error function. Figure 6 Three MTF curves are shown, derived from a real image, a reconstructed image based on the method of this application, and a comparative GGC method, respectively. It can be seen that the MTF curve of the reconstructed image obtained by the correction method of this application is closer to the real image, a result consistent with intuitive observation. These results all indicate that, under conditions of high noise, the correction method of this application can obtain more accurate geometric parameters. The comparison results of volume data reconstructed based on two sets of geometric parameters are shown below. Figure 7 As shown, (a) row images are reconstructed based on the geometric parameters obtained using the contrastive GGC method, and (b) row images are the result of obtaining the geometric parameters based on the proposed method. These images contain many fine structural features. From the magnified areas in the images, it can be seen that the images reconstructed using the geometric parameters obtained by the correction method based on this application have clearer boundaries (indicated by arrows).

[0078] Example 4. Validation in high-resolution micro-CT based on gantry structure

[0079] Currently, there are generally two types of scanning devices for micro-CT. The first is sample rotation-based, where the sample rotates while the source and detector remain stationary. This device typically has better spatial resolution and is often used for scanning ex vivo samples. Therefore, it is usually called ex vivo micro-CT. The second is gantry rotation-based, where the sample is stationary while the source and detector mounted on the gantry rotate. The advantage of this device is that it can be used to scan live small animals. This device is usually called in vivo micro-CT. The disadvantage of this setup is that it has lower spatial resolution compared to ex vivo micro-CT. One major reason is that in ex vivo micro-CT, the gantry carrying the X-ray source and detector is too heavy, especially in high-resolution scanning, making it more difficult to guarantee the roundness and eccentricity of the gantry. Generally, for micrometer-level resolution ex vivo micro-CT, the roundness of the trajectory is relatively easy to guarantee, so ex vivo micro-CT can be used to measure handcrafted phantoms. In this experiment, the phantom body was made of polymethyl methacrylate (PMMA) and contained eight metal spheres, each approximately 400 μm in diameter, spirally embedded within the body. We used an in vitro Micro-CT scanner with a resolution of 9 µm to perform eight scans of the phantom. During each scan, the phantom's position was randomly changed, uniformly acquiring 360 projections within a 360° range. The measured three-dimensional coordinates of the eight spheres were then registered to obtain the final three-dimensional coordinates of the spheres. The mean and variance of the sphere's three-dimensional coordinates are shown in Table 5. The table shows that the variances of the three components of the sphere coordinates are all at the sub-pixel level, indicating that the measurement results are stable and the measurement accuracy is very high in the actual system.

[0080] Table 5: Mean and variance of the three-dimensional coordinates of the eight spheres in the real experiment

[0081]

[0082] After obtaining the coordinate estimates of the spheres in the model, we used the calibration model of this application to calibrate a high-resolution Micro-CT based on gantry rotation (living body) to verify the accuracy of the coordinate estimation and the performance of the proposed calibration method in practical applications. This Micro-CT is equipped with a 5µm focal length X-ray tube and a 4008×2672 pixel CCD camera with a pixel pitch of 9µm. The X-ray source and CCD camera are mounted at opposite ends of the gantry, 20cm apart, with a magnification of approximately 2. The operating voltage is 60kV and the operating current is 130uA. The CCD camera has two binning modes: 1 × 1 and 2 × 2. In the 2 × 2 binning mode, the system resolution is approximately 9µm, which is essentially the same as the resolution of the ex vivo Micro-CT used to measure the calibration phantom. The first experiment below is based on this binning mode. The calibration phantom is placed on the scanning bed, and 400 projections are obtained within a 360° range. Due to the large load and high resolution, the trajectory of this Micro-CT is no longer a standard circle. To demonstrate this, we used the NOC method. The NOC method is a high-precision calibration method suitable for circular scan orbit geometry calibration. This embodiment also uses the classic GGC method from the simulation experiment to visually demonstrate the measurement accuracy of the sphere. This method can also be used as a comparison method to compare performance in real-world scenarios. After obtaining the geometric parameters using these three calibration methods, we collected mouse femoral projection data at 720 angles to reconstruct tomographic images. Assuming the trajectory is circular, we used the traditional FDK reconstruction method; for the other two cases, we used the general FDK algorithm provided by ASTRAToolbox for image reconstruction. In the second experiment, the binning mode was set to 1 × 1. At this point, the resolution of the Micro-CT is 4.5 µm, which is twice the resolution of the ex vivo Micro-CT used to measure the sphere's coordinates. This significantly amplifies the measurement error of the sphere, and these increased errors may have a greater negative impact on calibration accuracy. Therefore, the 1 × 1 mode can be used to further test the performance of the geometry calibration method under extreme conditions of more severe noise or measurement errors. In the second experiment, we again used the classic GGC method as a comparison method. The test sample was a mouse femur sample, with a total of 720 projections. The reconstruction method was consistent with the previous experiment. The reconstructed image of the mouse femur from the first experiment is shown below. Figure 8As shown, the resolution is 9µm. The images in columns (a), (b), and (c) are obtained based on the geometric parameters obtained using the methods of NOC, GGC, and this application, respectively. Each column contains one axial image and one sagittal image. The first column shows a very blurry image with many ghosting and artifacts. Since the NOC method is a high-precision geometric correction method suitable for circular scan trajectories, the only reason for these artifacts and blurring is that the trajectory is not a standard circle. The second column shows better image quality with no obvious artifacts. Combined with the conclusions drawn from the first column, this further proves that the system trajectory is not circular and that the measured ball coordinates are accurate, ensuring that the GGC method achieves relatively good correction results based on 8 balls. Comparing the images in the second and third columns reveals that the image quality is similar. This result indicates that when the phantom measurement system and the system to be corrected have similar resolutions, the correction method of this application has performance comparable to the compared GGC method. In the magnified areas selected from the reconstructed images in columns 2 and 3, subtle differences in details can be seen between the two images; for example, the magnified area in the third column shows sharper edges. This suggests that the calibration method in this application has higher calibration accuracy, and the reason for the insignificant performance improvement may be that the measurement of the microsphere is accurate enough that the calibration accuracy of the two methods is comparable. In the second experiment, the measured microsphere was used to calibrate a Micro-CT with a resolution of 4.5µm, and the measurement error was at least twice that of the previous method. The mouse femur reconstruction images from the second experiment are shown below. Figure 9 As shown, the first and second rows of images are the reconstruction results of geometric parameters obtained by the comparative GGC method and the correction method of this application, respectively. It can be intuitively seen that these images have a higher resolution than the images obtained in the first experiment, and the difference between the first and second rows is very obvious; that is, the first row of images has many artifacts, but the second row shows almost no artifacts. Furthermore, the details in the second row of images are very clear. We know that the resolution of the system to be corrected is much higher than that of the phantom measurement system. Therefore, noise will increase when extracting the projected centroid of the sphere obtained from the system to be corrected, and the measurement error of the sphere will be at least twice the original. However, we can still use the proposed correction method to obtain sufficiently accurate geometric parameters to generate high-quality reconstructed images. These results verify that the proposed iterative correction model has a better ability to suppress noise and measurement errors than classical methods.

[0083] As described above, although the invention has been shown and described with reference to specific preferred embodiments, it should not be construed as limiting the invention itself. Various changes in form and detail may be made without departing from the spirit and scope of the invention as defined in the appended claims.

Claims

1. A cone-beam CT geometric correction method applicable to irregular trajectories, characterized in that... Includes the following steps: S1: After obtaining the spatial coordinates of the randomly arranged balls through calculation, an iterative model is established to achieve geometric correction of a single projected view; S2: There are two coordinate systems in the model, namely and ;in This is called the detector coordinate system. and These are along the horizontal and vertical axes of the detector, respectively. It is the line-of-sight direction unique to the detector plane; This is called the object coordinate system; S3: Assume point B is in Coordinates in a coordinate system are Point B is located at The corresponding coordinates in the coordinate system can be represented as ； In the formula: T is the translation vector. , , This represents the components of the translation vector T along the X, Y, and Z axes; R is the rotation matrix. R is written as R = R1(η)R2(θ)R3(φ), where η, θ and φ are Euler angles, and R1, R2 and R3 represent yaw, pitch and roll angles respectively; S4: In the model, the S-coordinate of the x-ray source is S(w s , u s , v s ) T Point B w The coordinates of the projection point Q are Q = (0, u, v). T There are k points in the target coordinate system, called B. k K∈[1,k], their coordinates are given by b k = (x k , y k , z k ) T The corresponding projection point in the detector coordinate system is Q. k Its coordinates are Q k = (0, u k , v k ) T The intermediate solution of the geometric parameters in the iteration is (w) s , u s , v s , t1, t2, t3, η, φ, θ), B is obtained through equation (1) k Convert to the detector coordinate system to obtain its corresponding point B. w-k Its coordinates are written as = ; S5: Let point B be... k to line Q k The foot of the perpendicular line from S is H. k By establishing a system consisting of two straight lines Q k S and B k H k The system of equations formed yields H k coordinates h k Then deduce B k to line Q k Distance of S as follows: ； Combining all K points, the total distance It can be represented as: ； When all geometric parameters are accurate, the corresponding total distance =0; therefore, the basic objective function can be written as: ； Introduce a variance penalty term into equation (4). : ； In equation (5) ; After adding this term, the objective function becomes: ； In the formula, λ is a weighting factor with an empirical value of 10; In most cone-beam CT systems, the X-ray source and detector remain relatively stationary during the scan; by solving equation (6), the parameter vector for each projection view can be obtained. N is the total number of projected views; If there is no noise, all They should be the same; therefore all average The coordinates of the X-ray source are very close to the actual coordinates. After determining the source location, the objective function with 9 variables is simplified to a function with 6 variables, as shown below: （7）。

Images