A CT-CBCT Image Deformation Registration Method

Abstract

The present invention provides a CT-CBCT image deformation registration method. Belonging to the technical field of medical image processing, the specific steps are as follows: preprocess the CBCT image and the CT image; establish local control points for images A and B and determine the neighborhood and size of the control points; compare the average gray value E Ai of the neighborhood of the i-th control point Q in image A A with the average gray value E Bi of the neighborhood of the i-th control point Q in image B B to complete the deformation registration of a set of CT images and a set of CBCT images. The present invention uses two-level evaluation for the gray correction registration of the control point neighborhood to judge the dark area and the artifact area, so as to perform local gray correction registration on the two groups of images specifically according to whether it is a dark area and whether it is an artifact area; in addition, the present invention establishes local control points based on the gradient change characteristics of the two groups of images participating in the registration and adopts a grid division mechanism. The local gray correction registration based on the control point neighborhood avoids global blind registration.

Description

Technical Field

[0001] The present invention relates to the technical field of medical image processing, and in particular to a CT-CBCT image deformation registration method, and in particular to a CT-CBCT image deformation registration method based on control point neighborhood grayscale correction. Background Art

[0002] Medical image registration is one of the key technologies in precision radiotherapy and provides a basis for disease diagnosis, target delineation, precise positioning and precise treatment. Multimodal image deformation registration can provide doctors and physicists with richer medical image information.

[0003] Traditional deformable registration methods often cause CT-CBCT image mismatching or mismatching due to problems such as insufficient image contrast and grayscale discontinuity when dealing with CT-CBCT deformable registration problems. To solve these problems, domestic and foreign researchers have proposed a global intensity correction method to preprocess CBCT images, with the aim of converting the CBCT-CT multimodal deformable registration into a homomodal deformable registration problem. However, in actual CT-CBCT deformable registration research, it is found that artifacts always exist locally, which will destroy the original intensity or grayscale distribution correspondence. Therefore, it is necessary to perform local grayscale correction on CBCT images with artifacts, and match CT-CBCT images while correcting the CBCT images to complete the final CT-CBCT registration.

[0004] The present invention overcomes the drawback that the intensity distribution or grayscale distribution caused by CBCT artifacts does not correspond, resulting in mismatch or mismatch of CT-CBCT images. Instead of blindly performing global registration on CT-CBCT images, the present invention performs local grayscale correction registration on CT-CBCT images. First, local control points are constructed based on the grayscale features of CBCT images and by using a grid division mechanism. Then, a two-level evaluation is used to determine the dark area and the artifact area for grayscale correction registration of the control point neighborhood. Thus, local grayscale correction registration is performed on the two groups of images in a targeted manner according to whether there is a dark area and whether there is an artifact area. Finally, local grayscale correction registration is performed on the CBCT images of all layers to be registered to obtain the deformation matrix of all layers of CBCT images, thereby completing the deformation registration of the entire group of CT and CBCT images. Summary of the invention

[0005] The technical problem to be solved by the present invention is to overcome the deficiency that the intensity distribution or gray-scale distribution does not correspond due to CBCT artifacts, resulting in misregistration or incorrect registration of CT-CBCT images. A CT-CBCT image deformation registration method based on gray-scale correction in the neighborhood of control points is proposed. The gray scale in the neighborhood of the control points is evaluated at two levels to identify the dark area and the artifact area. Then, corresponding local gray-scale correction is performed according to whether it is a dark area and whether it is an artifact area, and the deformation matrix of the CBCT image is obtained to complete the deformation registration of the CT-CBCT images.

[0006] The technical solution of the present invention is: a CT-CBCT image deformation registration method, and the specific operation steps are as follows:

[0007] Step (1.1), preprocess the CBCT image and the CT image;

[0008] Step (1.2), establish local control points for images A and B and determine the neighborhood and size of the control points;

[0009] Step (1.3), perform the first-level evaluation on the neighborhood of the local control points of image A, that is, the dark area evaluation, and compare the gray-scale mean value E Ai of the neighborhood of the i-th control point Q A of image A with the gray-scale mean value E Bi of the neighborhood of the i-th control point Q B of image B;

[0010] If E A ≥m E E B , where m E is a constant, and the value is: 0.5 ≤ m E ≤ 1.0, then it is determined that the neighborhood of the i-th control point Q Ai of image A is a non-dark area, and histogram correction registration is performed on the neighborhood of the i-th control point of image A with reference to the neighborhood of the i-th control point of image B, and the deformation matrix M 1i of the neighborhood of the local control points of image A is output, and the matrix size is the same as D Ai ;

[0011] Finally, the deformation matrices of all the control point neighborhoods of image A are aggregated and mapped to the corresponding positions in image A to obtain the deformation matrix of image A and at the same time, the deformation registration of image A and image B is completed;

[0012] If E A <m E E B , where m E is a constant, and the value is: 0.5 ≤ m E ≤ 1.0, then it is determined that the neighborhood of the i-th control point Q Ai of image A is a dark area, and at this time, the control point Q of image A needs to beAi Perform a second-level evaluation on the neighborhood, that is, artifact evaluation, and then compare the i-th control point Q in image A Ai The mean square deviation of gray levels in the neighborhood with the i-th control point Q in image B Bi The square root of the mean square deviation of gray levels in the neighborhood in terms of magnitude,

[0013] If where m σ is a constant, and the value range is: 0.1 ≤ m σ ≤ 0.5, then it is determined that the neighborhood of the i-th control point Q in image A Ai is a non-artifact area. Perform linear gray-level mapping on the neighborhood of the i-th control point in image A with reference to the neighborhood of the i-th control point in image B, and output the deformation matrix M of the local control point neighborhood in image A 2i , and the matrix size is the same as that of D Ai ;

[0014] If where m σ is a constant, 0.1 ≤ m σ ≤ 0.5, then it is determined that the neighborhood of the i-th control point Q in image A Ai is an artifact area. First, locally magnify the gray levels in the neighborhood of Q Ai with the magnification factor being the ratio of the mean gray level in the neighborhood of image B to the mean gray level in the neighborhood of image A, and then perform histogram correction registration on the neighborhood of Q Ai with reference to the neighborhood of Q Bi , and output the deformation matrix M of the local control point neighborhood in image A 3i , and the matrix size is the same as that of D Ai ;

[0015] Finally, aggregate the deformation matrices of all control point neighborhoods in image A and map them to the corresponding positions in image A to obtain the deformation matrix of image A Meanwhile, complete the deformation registration between image A and image B;

[0016] Step (1.4), after obtaining the deformation matrix of image A by mapping the deformation matrices of all control point neighborhoods in image A to the corresponding positions in image A , then perform local correction registration on all CBCT images to be registered simultaneously, and finally complete the deformation registration of a set of CT and a set of CBCT images.

[0017] Further, in step (1.1), the specific operation method for preprocessing the CBCT image and the CT image is as follows: Take the intersection region of the CBCT image and the CT image, and respectively crop the CBCT image and the CT image based on this intersection region. At the same time, sample the CBCT image with reference to the CT image to obtain CBCT images and CT images with the same resolution and size. Take all the gray-scale information of the CBCT image A and the CT image B of the layer to be registered as the input. Denote the width and height of image A as W and H respectively, then the width and height of image B are both W and H, where 1 ≤ i ≤ H, 1 ≤ j ≤ W.

[0018] Further, in step (1.2), the specific operation steps for establishing local control points of images A and B and determining the neighborhood and size of the control points are as follows:

[0019] (1.2.1), Construct local control points based on image A to obtain all control point sets Q of image A A , Denote the control point set Q A There are N points, and the i-th control point is denoted as Q Ai , where 1 ≤ i ≤ N. At the same time, based on the gray-scale features of each local control point of image A, obtain the local control point set Q of image B B , Denote the control point set Q B There are N points, and the i-th control point is denoted as Q Bi , where 1 ≤ i ≤ N;

[0020] (1.2.2), Select the long-direction region centered on the local control point Q of image A as the neighborhood gray-scale range of the control point Q Ai . Denote the neighborhood size of the control point Q of image A Ai as D Ai , where 1 ≤ i ≤ N. Calculate the horizontal displacement difference and vertical displacement difference between Q Ai and other control points in image A, and take the absolute values of the minimum horizontal displacement difference and vertical displacement difference as the length and width of the size D Ai respectively. Similarly, select the long-direction region centered on the local control point Q of image B as the neighborhood gray-scale range of the control point Q Ai . Denote the neighborhood size of the control point Q of image B Bi as D Bi , where 1 ≤ i ≤ N. D Bi and D Bi take the same size. Bi and D Ai take the same size.

[0021] Further,

[0022] Establish a uniform grid of R×R for image A, where 3 ≤ R ≤ 10 and R is a positive integer. Find the points with obvious gradient changes in each grid as the control points of image A, and compare the deformation matrices obtained when R takes the eight values of 3, 4, 5, 6, 7, 8, 9, and 10 respectively. Compare the gray-scale mutual correlation coefficient NMI of the mapped image A and image B in the eight cases. The larger the NMI coefficient, the better the registration result, and thus the value of R is selected accordingly.

[0023] Further, in step (1.3), the m E and m σ are both constants. m E takes the values of 0.5, 0.6, 0.7, 0.8, 0.9, and 1.0 respectively, and m σ takes the values of 0.1, 0.2, 0.3, 0.4, and 0.5 respectively to obtain the deformation matrix. Similarly, compare the gray-scale mutual correlation coefficient NMI of the mapped image A and image B. The larger the NMI coefficient, the better the registration result, and thus the values of m E and m σ are selected accordingly.

[0024] Further, in step (1.4), obtain the deformation matrix P of the current layer CBCT image, that is, image A, and then perform local gray-scale correction registration on the CBCT images of all layers to be registered to obtain the deformation matrices of all layers of CBCT images. Among them, 0 < w < m, m is the CBCT images of all layers to be registered, and the deformation registration of a group of CT and a group of CBCT images is completed.

[0025] The beneficial effects of the present invention are: (1) Overcome the deficiency that the intensity distribution or gray-scale distribution is not corresponding caused by CBCT artifacts, resulting in misregistration or incorrect registration of CT-CBCT images. For the local gray-scale correction registration of the control point neighborhood, two-level evaluation is used to judge the dark area and the artifact area, so as to perform local gray-scale correction registration on the two groups of images specifically according to whether it is a dark area and whether it is an artifact area, avoiding global blind correction of CBCT images; (2) Based on the two groups of images participating in the registration, local control points are established by using the gradient change characteristics and the grid division mechanism. The local gray-scale correction registration based on the control point neighborhood avoids global blind registration. Brief Description of the Drawings

[0026] Figure 1 is the structural flowchart of the present invention;

[0027] Figure 2 Among them, (a) is the CT and CBCT images; (b) is the output image of the deformation registration of the present invention. Detailed Embodiments

[0028] To more clearly illustrate the technical solution of the present invention, the technical solution of the present invention will be further described in detail below with reference to the accompanying drawings:

[0029] As Figure 1 described, a CT-CBCT image deformation registration method, the specific operation steps are as follows:

[0030] Step (1.1), preprocess the CBCT image and the CT image;

[0031] Step (1.2), establish local control points for images A and B and determine the neighborhood and size of the control points;

[0032] Step (1.3), perform a first-level evaluation on the neighborhood of the local control points of image A, that is, dark area evaluation, and compare the average gray value E Ai of the neighborhood of the i-th control point Q A of image A with the average gray value E Bi of the neighborhood of the i-th control point Q B of image B;

[0033] If E A ≥m E E B , where m E is a constant, and the value is: 0.5 ≤ m E ≤1.0, then it is determined that the neighborhood of the i-th control point Q Ai of image A is a non-dark area, and histogram correction registration is performed on the neighborhood of the i-th control point of image A with reference to the neighborhood of the i-th control point of image B, and the deformation matrix M 1i of the neighborhood of the local control points of image A is output, and the matrix size is consistent with D Ai ;

[0034] Finally, summarize the deformation matrices of all the control point neighborhoods of image A and map them to the corresponding positions in image A to obtain the deformation matrix of image A and at the same time complete the deformation registration of image A and image B;

[0035] If E A <m E E B , where m E is a constant, and the value is: 0.5 ≤ m E ≤1.0, then it is determined that the neighborhood of the i-th control point Q Ai of image A is a dark area. At this time, a second-level evaluation of the neighborhood of the control point Q Ai of image A needs to be performed, that is, artifact evaluation, and then compare the mean square deviation of the gray value Ai of the neighborhood of the i-th control point Q of image A with the square root of the mean square deviation of the gray value Bi of the neighborhood of the i-th control point Q The size of

[0036] If where m σ is a constant, and its value range is: 0.1 ≤ m σ ≤ 0.5, then it is determined that the neighborhood of the i-th control point Q of image A Ai is a non-artifact area. Perform linear gray value mapping on the neighborhood of the i-th control point of image A with reference to the neighborhood of the i-th control point of image B, and output the deformation matrix M of the local control point neighborhood of image A 2i , and the matrix size is the same as that of D Ai ;

[0037] If where m σ is a constant, 0.1 ≤ m σ ≤ 0.5, then it is determined that the neighborhood of the i-th control point Q of image A Ai is an artifact area. First, locally magnify the gray value of the neighborhood of Q Ai , and the magnification factor is the ratio of the average gray value of the neighborhood of image B to the average gray value of the neighborhood of image A. Then, perform histogram correction registration on the neighborhood of Q Ai with reference to the neighborhood of Q Bi , and output the deformation matrix M of the local control point neighborhood of image A 3i , and the matrix size is the same as that of D Ai ;

[0038] Finally, summarize the deformation matrices of all control point neighborhoods of image A and map them to the corresponding positions in image A to obtain the deformation matrix of image A At the same time, complete the deformation registration between image A and image B;

[0039] Step (1.4), after obtaining the deformation matrix of image A by mapping the deformation matrices of all control point neighborhoods of image A to the corresponding positions in image A , then perform local correction registration on all CBCT images of the layers to be registered at the same time, and finally complete the deformation registration of a group of CT and a group of CBCT images.

[0040] Furthermore, in step (1.1), the specific operation method for preprocessing the CBCT image and the CT image is as follows: Take the intersection area of the CBCT image and the CT image and crop the CBCT image and the CT image respectively based on this intersection area. At the same time, sample the CBCT image with reference to the CT image to obtain CBCT images and CT images with the same resolution and size. Use all the gray information of the CBCT image A and the CT image B of the layer to be registered as the input. Denote the width and height of image A as W and H respectively, then the width and height of image B are both W and H, where 1 ≤ i ≤ H, 1 ≤ j ≤ W.

[0041] Further, in step (1.2), the specific operation steps for establishing local control points of images A and B and determining the neighborhood and size of the control points are as follows:

[0042] (1.2.1) Based on image A, construct local control points to obtain all control point sets Q of image A A , denote the control point set Q A has N points, and the i-th control point is denoted as Q Ai , where 1 ≤ i ≤ N. At the same time, based on the gray feature of each local control point of image A, obtain the local control point set Q of image B B , denote the control point set Q B has N points, and the i-th control point is denoted as Q Bi , where 1 ≤ i ≤ N;

[0043] (1.2.2) Select the long-direction area centered on the local control point Q of image A Ai as the neighborhood gray range of the control point Q Ai . Denote the neighborhood size of the image A control point Q Ai as D Ai , where 1 ≤ i ≤ N. Calculate the horizontal displacement difference and vertical displacement difference between Q Ai and other control points in image A, and take the absolute values of the minimum horizontal displacement difference and vertical displacement difference as the length and width of the size D Ai . Similarly, select the long-direction area centered on the local control point Q of image B Bi as the neighborhood gray range of the control point Q Bi . Denote the neighborhood size of the image B control point Q Bi as D Bi , where 1 ≤ i ≤ N. D Bi and D Ai take the same size.

[0044] Further,

[0045] Establish an R×R uniform grid for image A, where 3 ≤ R ≤ 10 and R is a positive integer. Find the points with obvious gradient changes in each grid as the control points of image A. Compare the deformation matrices obtained when R takes 3, 4, 5, 6, 7, 8, 9, and 10 respectively Compare the size of the gray mutual correlation coefficient NMI between the mapped image A and image B in the eight cases. The larger the NMI coefficient, the better the registration result, and thus the value of R is selected accordingly.

[0046] Further, in step (1.3), the m E and m σ are both constants. m E are taken as 0.5, 0.6, 0.7, 0.8, 0.9, 1.0 respectively, and mσ The deformation matrix is ​​obtained when taking 0.1, 0.2, 0.3, 0.4, and 0.5 respectively Similarly, by comparing the grayscale mutual correlation coefficient NMI of the mapped image A and image B, the larger the NMI coefficient, the better the registration result, and this is used to determine m. E With m σ The value of .

[0047] Furthermore, in step (1.4), the deformation matrix P of the current layer CBCT image, i.e., image A, is obtained, and then the CBCT images of all layers to be registered are locally grayscale corrected and registered to obtain the deformation matrix of all layers of CBCT images: Wherein, 0<w<m, m is the CBCT images of all layers to be registered, and the deformation registration of a group of CT images and a group of CBCT images is completed.

[0048] like Figure 2 The above-mentioned (a) is a CT and CBCT image (the left image is CT, and the right image is CBCT); (b) is the output image of the deformation registration of the present invention; Figure 2 As described above, the present invention adopts a two-level evaluation for the grayscale correction and registration of the control point neighborhood to determine the dark area and the artifact area, thereby performing local grayscale correction and registration on the two groups of images in a targeted manner according to whether it is a dark area and whether it is an artifact area; in addition, the present invention establishes local control points based on the gradient change characteristics of the two groups of images participating in the registration and adopts a grid division mechanism, and performs local grayscale correction and registration based on the control point neighborhood to avoid global blind matching.

[0049] Finally, it should be understood that the embodiments described in the present invention are only used to illustrate the principles of the embodiments of the present invention; other variations may also fall within the scope of the present invention; therefore, as examples rather than limitations, alternative configurations of the embodiments of the present invention may be considered consistent with the teachings of the present invention; accordingly, the embodiments of the present invention are not limited to the embodiments explicitly introduced and described in the present invention.

Claims

1. A CT-CBCT image deformation registration method, characterized in that, The specific operation steps are as follows: Step (1.1): Preprocess the CBCT image and the CT image, take the intersection region of the CBCT image and the CT image, and respectively crop the CBCT image and the CT image based on this intersection region. At the same time, sample the CBCT image with reference to the CT image to obtain CBCT images and CT images with the same resolution and the same size. Use all the gray information of the CBCT image A and the CT image B of the layer to be registered as the input. Denote the CBCT image as image A, denote the CT image as image B, and denote the width and height of image A as W and H respectively. Then the width and height of image B are both W and H; Step (1.2): Establish local control points for images A and B and determine the neighborhood and size of the control points. The specific operation steps are as follows: (1.2.1) Construct local control points based on Image A to obtain all control point sets Q of Image A A , denote the control point set Q A has N points, and the i-th control point is denoted as Q Ai , where 1 ≤ i ≤ N. Meanwhile, based on the gray feature of each local control point of Image A, obtain the local control point set Q of Image B B , denote the control point set Q B has N points, and the i-th control point is denoted as Q Bi , where 1 ≤ i ≤ N; (1.2.2) Select the local control point Q of Image A Ai The long-direction area centered on it as the control point Q Ai The neighborhood gray range of Image A's control point Q Ai The neighborhood size of Image A's control point Q is denoted as D Ai , where 1 ≤ i ≤ N, calculate the horizontal displacement difference and vertical displacement difference between Q Ai And other control points in Image A. Take the absolute values of the minimum horizontal displacement difference and vertical displacement difference as the length and width of D Ai Similarly, select the local control point Q of Image B Bi The long-direction area centered on it as the control point Q Bi The neighborhood gray range of Image B's control point Q Bi The neighborhood size of Image B's control point Q is denoted as D Bi , where 1 ≤ i ≤ N, D Bi And D Ai Take the same size; Step (1.3): Perform the first-level evaluation on the local control point neighborhood of Image A, i.e., the dark area evaluation, and compare the gray-scale mean value E Ai of the neighborhood of the i-th control point Q A in Image A with the gray-scale mean value E Bi of the neighborhood of the i-th control point Q B in Image B; If E A ≥m E E B , where m E is a constant with a value range of: 0.5 ≤ m E ≤ 1.0, then it is determined that the neighborhood of the i-th control point Q Ai of image A is a non-dark area. The histogram correction registration is performed on the neighborhood of the i-th control point of image A with reference to the neighborhood of the i-th control point of image B, and the deformation matrix M 1i of the local control point neighborhood of image A is output. The size of the matrix is consistent with D Ai ; If E A <m E E B , where m E is a constant, and its value is: 0.5 ≤ m E ≤ 1.0, then it is determined that the neighborhood of the i-th control point Q Ai of the image A is a dark area. At this time, a second-level evaluation of the neighborhood of the control point Q Ai of the image A is required, that is, artifact evaluation. Then, compare the square root of the mean square error of the gray levels in the neighborhood of the i-th control point Q Ai of the image A with the square root of the mean square error of the gray levels in the neighborhood of the i-th control point Q Bi of the image B in terms of magnitude. If where m σ is a constant, and the value range is: 0.1 ≤ m σ ≤ 0.5, then it is determined that the neighborhood of the i-th control point Q of image A Ai is a non-artifact area. Perform linear gray value mapping on the neighborhood of the i-th control point of image A with reference to the neighborhood of the i-th control point of image B, and output the local control point neighborhood deformation matrix M 2i of image A. The size of the matrix is the same as that of D Ai ; If where m σ is a constant, 0.1 ≤ m σ ≤ 0.5, then it is determined that the neighborhood of the i-th control point Q Ai of image A is an artifact area. First, locally magnify the gray level of the neighborhood of Q Ai , and the magnification factor is the ratio of the average gray level of the neighborhood of image B to the average gray level of the neighborhood of image A. Then, perform histogram correction registration on the neighborhood of Q Ai with reference to the neighborhood of Q Bi to output the local control point neighborhood deformation matrix M 3i of image A. The size of the matrix is the same as that of D Ai ; Finally, the deformation matrices of all the control point neighborhoods of image A are mapped to the corresponding positions in image A to obtain the deformation matrix of image A. At the same time, the deformation registration between image A and image B is completed. Step (1.4): Perform deformation registration on images A and B of all layers to be registered, and finally complete the deformation registration of a group of CT images and a group of CBCT images.

2. The CT-CBCT image deformation registration method according to claim 1, wherein In step (1.2), an R×R uniform grid is established for image A, where 3 ≤ R ≤ 10 and R is a positive integer. Points with obvious gradient changes are found as control points of image A in each grid, and the deformation matrices obtained in eight cases where R takes 3, 4, 5, 6, 7, 8, 9, and 10 respectively are compared. The gray-scale mutual correlation coefficient NMI of the mapped image A and image B in the eight cases is compared. The larger the NMI coefficient, the better the registration result, and thus the value of R is selected accordingly.

3. A CT-CBCT image deformation registration method according to claim 1, wherein In step (1.3), the m E and m σ are both constants. m E takes values of 0.5, 0.6, 0.7, 0.8, 0.9, 1.0 respectively, and m σ takes values of 0.1, 0.2, 0.3, 0.4, 0.5 respectively, and the deformation matrix is obtained. Similarly, by comparing the magnitude of the gray-scale cross-correlation coefficient NMI between the mapped image A and the image B, the larger the NMI coefficient, the better the registration result, and thus the values of m E and m σ are selected.

4. A CT-CBCT image deformation registration method according to claim 1, characterized in that, In step (1.4), the deformation matrix P of the current layer CBCT image, i.e., image A, is obtained, and then local gray correction registration is performed on the CBCT images of all layers to be registered to obtain the deformation matrices of the CBCT images of all layers. Among them, 1 ≤ w ≤ m, and m is the number of all layers to be registered of the CBCT images.

Images