CN113450396B - Three-dimensional/two-dimensional image registration method and device based on bone characteristics - Google Patents
Three-dimensional/two-dimensional image registration method and device based on bone characteristics Download PDFInfo
- Publication number
- CN113450396B CN113450396B CN202110682717.9A CN202110682717A CN113450396B CN 113450396 B CN113450396 B CN 113450396B CN 202110682717 A CN202110682717 A CN 202110682717A CN 113450396 B CN113450396 B CN 113450396B
- Authority
- CN
- China
- Prior art keywords
- dimensional
- image
- data
- registration
- gradient
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 58
- 210000000988 bone and bone Anatomy 0.000 title claims abstract description 17
- 238000012549 training Methods 0.000 claims abstract description 42
- 238000011524 similarity measure Methods 0.000 claims abstract description 37
- 230000009466 transformation Effects 0.000 claims abstract description 37
- 230000008569 process Effects 0.000 claims abstract description 25
- 238000009877 rendering Methods 0.000 claims abstract description 25
- 238000005070 sampling Methods 0.000 claims abstract description 16
- 238000005457 optimization Methods 0.000 claims abstract description 15
- 238000012360 testing method Methods 0.000 claims abstract description 11
- 238000004088 simulation Methods 0.000 claims abstract description 10
- 238000007781 pre-processing Methods 0.000 claims abstract description 7
- 230000006870 function Effects 0.000 claims description 42
- 238000013519 translation Methods 0.000 claims description 18
- 239000011159 matrix material Substances 0.000 claims description 13
- 239000013598 vector Substances 0.000 claims description 12
- 238000004422 calculation algorithm Methods 0.000 claims description 9
- 238000013527 convolutional neural network Methods 0.000 claims description 8
- 238000012545 processing Methods 0.000 claims description 7
- 238000004364 calculation method Methods 0.000 claims description 5
- 230000010354 integration Effects 0.000 claims description 4
- 238000005259 measurement Methods 0.000 claims description 4
- 238000009826 distribution Methods 0.000 claims description 3
- 238000011156 evaluation Methods 0.000 claims description 3
- 101100379079 Emericella variicolor andA gene Proteins 0.000 claims description 2
- 108010063499 Sigma Factor Proteins 0.000 claims description 2
- 229940060587 alpha e Drugs 0.000 claims description 2
- 238000005266 casting Methods 0.000 claims description 2
- PCHJSUWPFVWCPO-UHFFFAOYSA-N gold Chemical group [Au] PCHJSUWPFVWCPO-UHFFFAOYSA-N 0.000 claims description 2
- 238000003384 imaging method Methods 0.000 claims description 2
- 238000011176 pooling Methods 0.000 claims description 2
- 230000000694 effects Effects 0.000 description 5
- 238000012804 iterative process Methods 0.000 description 3
- 230000003044 adaptive effect Effects 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000003860 storage Methods 0.000 description 2
- 238000012935 Averaging Methods 0.000 description 1
- 230000006978 adaptation Effects 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 239000000463 material Substances 0.000 description 1
- 230000005405 multipole Effects 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- 230000000704 physical effect Effects 0.000 description 1
- 230000011218 segmentation Effects 0.000 description 1
- 238000001356 surgical procedure Methods 0.000 description 1
- 230000000007 visual effect Effects 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/30—Determination of transform parameters for the alignment of images, i.e. image registration
- G06T7/33—Determination of transform parameters for the alignment of images, i.e. image registration using feature-based methods
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
- G06N3/00—Computing arrangements based on biological models
- G06N3/02—Neural networks
- G06N3/04—Architecture, e.g. interconnection topology
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
- G06N3/00—Computing arrangements based on biological models
- G06N3/02—Neural networks
- G06N3/08—Learning methods
- G06N3/084—Backpropagation, e.g. using gradient descent
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T19/00—Manipulating 3D models or images for computer graphics
- G06T19/20—Editing of 3D images, e.g. changing shapes or colours, aligning objects or positioning parts
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/11—Region-based segmentation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10072—Tomographic images
- G06T2207/10081—Computed x-ray tomography [CT]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20081—Training; Learning
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20084—Artificial neural networks [ANN]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/30—Subject of image; Context of image processing
- G06T2207/30004—Biomedical image processing
- G06T2207/30008—Bone
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2219/00—Indexing scheme for manipulating 3D models or images for computer graphics
- G06T2219/20—Indexing scheme for editing of 3D models
- G06T2219/2004—Aligning objects, relative positioning of parts
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2219/00—Indexing scheme for manipulating 3D models or images for computer graphics
- G06T2219/20—Indexing scheme for editing of 3D models
- G06T2219/2016—Rotation, translation, scaling
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
- Y02T10/00—Road transport of goods or passengers
- Y02T10/10—Internal combustion engine [ICE] based vehicles
- Y02T10/40—Engine management systems
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Software Systems (AREA)
- General Engineering & Computer Science (AREA)
- Computing Systems (AREA)
- Mathematical Physics (AREA)
- Data Mining & Analysis (AREA)
- Evolutionary Computation (AREA)
- General Health & Medical Sciences (AREA)
- Molecular Biology (AREA)
- Biophysics (AREA)
- Biomedical Technology (AREA)
- Artificial Intelligence (AREA)
- Computational Linguistics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Health & Medical Sciences (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Architecture (AREA)
- Computer Graphics (AREA)
- Computer Hardware Design (AREA)
- Image Analysis (AREA)
- Image Processing (AREA)
Abstract
The three-dimensional/two-dimensional image registration method and device based on the skeleton features have small occupied resources, can be used for back propagation optimization of registration parameters, solve the problem that an image similarity function is not convex or concave in the three-dimensional/two-dimensional registration process, enlarge the initial registration value range, realize the registration of the three-dimensional/two-dimensional skeleton features in the TIPS operation, and improve the registration precision. The method comprises the following steps: (1) preprocessing the segmented three-dimensional bone data; (2) Generating a large amount of simulation data for subsequent training and testing; (3) The three-dimensional data is processed by a three-dimensional convolution network jump connection module to obtain parameterized three-dimensional characteristic data, and the superimposed original image and the three-dimensional characteristic data are projected to two dimensions by a differential rendering module based on grid sampling; (4) Calculating similarity measure between the two-dimensional decoding network and the target image by using the two-dimensional decoding network; in order for the similarity measure to have the property of a concave function in European space, a gradient of geodesic distance of a rigid transformation parameter is used as training target.
Description
Technical Field
The invention relates to the technical field of medical image processing, in particular to a three-dimensional/two-dimensional image registration method based on skeleton characteristics and a three-dimensional/two-dimensional image registration device based on skeleton characteristics.
Background
In the field of three-dimensional/two-dimensional registration, in the iterative process of the traditional registration method, re-projection and similarity measure calculation are required to be repeatedly carried out, so that excessive calculation resources are caused, and iterative optimization of parameters of a network model is inconvenient. Moreover, the similarity measure has a great limitation in solving the cross-modal problem due to the presence of the multipole problem.
Disclosure of Invention
In order to overcome the defects of the prior art, the technical problem to be solved by the invention is to provide a three-dimensional/two-dimensional image registration method based on skeleton features, which occupies small resources and can be used for counter-propagating and optimizing registration parameters, so that the problem that an image similarity function is not convex or concave in the three-dimensional/two-dimensional registration process is solved, the registration initial value range is enlarged, the registration of the three-dimensional/two-dimensional skeleton features in a TIPS (tip-in-place) operation is realized, and the method has an effective effect on improving the registration precision.
The technical scheme of the invention is as follows: the three-dimensional/two-dimensional image registration method based on the bone characteristics comprises the following steps:
(1) Preprocessing the segmented three-dimensional bone data, including region-of-interest clipping and downsampling;
(2) Carrying out random translation and rotation on the processed TIPS preoperative CT skeleton data to obtain different initial poses, and then carrying out random translation and rotation again to obtain target poses; projecting the TIPS preoperative CT data of the two poses to generate a large amount of simulation data for subsequent training and testing;
(3) In the training stage, three-dimensional data is firstly processed by a three-dimensional convolution network jump connection module to obtain parameterized three-dimensional characteristic data, and then the superimposed original image and the three-dimensional characteristic data are projected to two dimensions by a differential rendering module based on grid sampling;
(4) According to the two-dimensional image obtained by projection, calculating the similarity measure between the two-dimensional image and the target image by using a two-dimensional decoding network; in order to make the similarity measure obtained by network training have the property of a concave function in European space, a gradient of geodesic distance of a rigid transformation parameter is used as a training target.
Firstly, a traditional three-dimensional/two-dimensional registration frame is built; then constructing a differential rendering module based on grid sampling, wherein the differential rendering module occupies small resources and can be used for back propagation to optimize registration parameters, and replaces rigid transformation and re-projection of the volume data in each iterative process in a registration algorithm; finally, in the training stage of the concave function similarity measure model, the geodesic distance is used as a target, the problem that the image similarity function is not convex or concave in the three-dimensional/two-dimensional registration process is solved, the registration initial value range is enlarged, the registration of three-dimensional/two-dimensional skeleton features in the TIPS operation is realized, and the method has an effective effect on improving the registration precision.
There is also provided a three-dimensional/two-dimensional image registration apparatus based on bone features, comprising:
a data preprocessing module configured to preprocess the segmented three-dimensional bone data, including region of interest clipping and downsampling;
the simulation data generation module is configured to perform random translation and rotation on the processed TIPS preoperative CT skeleton data to obtain different initial poses, and then perform random translation and rotation again to obtain target poses; projecting the TIPS preoperative CT data of the two poses to generate a large amount of simulation data for subsequent training and testing;
the training module is configured to obtain parameterized three-dimensional characteristic data by processing the three-dimensional data through the three-dimensional convolution network jump connection module at the training stage, and then project the superimposed original image and the three-dimensional characteristic data to two dimensions through the differential rendering module based on grid sampling;
the similarity measure calculating module is used for calculating the similarity measure of the two-dimensional image obtained by projection and the target image according to the two-dimensional image obtained by projection by registration; in order to make the similarity measure obtained by network training have the property of a concave function in European space, a gradient of geodesic distance of a rigid transformation parameter is used as a training target.
Drawings
Fig. 1 is a flow chart of a three-dimensional/two-dimensional image registration method based on bone features according to the present invention.
Fig. 2 shows the steps of the covariance matrix adaptive evolution strategy registration algorithm.
Detailed Description
As shown in fig. 1, the three-dimensional/two-dimensional image registration method based on bone features comprises the following steps:
(1) Preprocessing the segmented three-dimensional bone data, including region-of-interest clipping and downsampling;
(2) Carrying out random translation and rotation on the processed TIPS preoperative CT skeleton data to obtain different initial poses, and then carrying out random translation and rotation again to obtain target poses; projecting the TIPS preoperative CT data of the two poses to generate a large amount of simulation data for subsequent training and testing;
(3) In the training stage, three-dimensional data is firstly processed by a three-dimensional convolution network jump connection module to obtain parameterized three-dimensional characteristic data, and then the superimposed original image and the three-dimensional characteristic data are projected to two dimensions by a differential rendering module based on grid sampling;
(4) According to the two-dimensional image obtained by projection, calculating the similarity measure between the two-dimensional image and the target image by using a two-dimensional decoding network; in order to make the similarity measure obtained by network training have the property of a concave function in European space, a gradient of geodesic distance of a rigid transformation parameter is used as a training target.
Firstly, a traditional three-dimensional/two-dimensional registration frame is built; then constructing a differential rendering module based on grid sampling, wherein the differential rendering module occupies small resources and can be used for back propagation to optimize registration parameters, and replaces rigid transformation and re-projection of the volume data in each iterative process in a registration algorithm; finally, in the training stage of the concave function similarity measure model, the geodesic distance is used as a target, the problem that the image similarity function is not convex or concave in the three-dimensional/two-dimensional registration process is solved, the registration initial value range is enlarged, the registration of three-dimensional/two-dimensional skeleton features in the TIPS operation is realized, and the method has an effective effect on improving the registration precision.
Preferably, in the step (4), two three-dimensional/two-dimensional registration frames are included:
(I) The test data is subjected to iterative optimization through concave function similarity measurement comprising a differentiable rendering module so as to obtain a rigid transformation parameter corresponding to the target image;
(II) optimization using a three-dimensional/two-dimensional registration optimizer, a conventional parameterized registration framework is implemented.
In the process of three-dimensional/two-dimensional matching image registration, the angle of three-dimensional projection needs to be continuously transformed, so that three-dimensional data is projected into a two-dimensional space. Such a process may alsoEquivalent expression is inverse rigid transformation of three-dimensional volume data, i.e. projection with the projection environment fixed, which is adopted by this subsection as T (θ x ,θ y ,θ z ,t x ,t y ,t z )·V m Simplified expression was performed. Preferably, in the step (3), the image projected to two dimensions is expressed by formula (4.1):
I m =P·T(θ x ,θ y ,θ z ,t x ,t y ,t z )·V m (4.1)
wherein Vm Is three-dimensional volume data, I m Is a corresponding image projected to two dimensions, wherein T (alpha), alpha E SE (3) is a three-dimensional rigid transformation matrix, the formulas are (4.2) - (4.5), P is an expression of a 3 x 4 point source projection matrix, and alpha= (theta) x ,θ y ,θ z ,t x ,t y ,t z ) Refer to the rigidity transformation parameters
wherein θx ,θ y ,θ z The rotation angles around the x coordinate axis, the y coordinate axis and the z coordinate axis are respectively t x ,t y ,t z To translate along three axes, x 0 ,y 0 ,z 0 Is V m Coordinates of the pixel points on the pixel points, x, y, z are x 0 ,t 0 ,z 0 Coordinates of the rigidly transformed points。
The three-dimensional/two-dimensional registration is aimed at the I after projection m With reference image (Fixed image) I f More closely, the key to the problem of registration is to measure I f and Im Similarity measure between the two. Preferably, in the step (3),
measurement of I f and Im Similarity measure between them, which represents a loss function of formula (4.6)
wherein ,p″CT Representing the projection of three-dimensional data onto a two-dimensional camera internal reference matrix, the gradient difference (Gradient Difference, GD) and the normalized cross-correlation coefficient (Normalized Cross Correlation, NCC) are similarity measure functions commonly used in conventional registration methods;
gradient differenceIs based on image i M and if The similarity calculated for the gradients of (2) is suitable for processing relatively long and narrow structures, as shown in formulas (4.7), (4.8) and (4.9)
Diff width and Diffhigh Representing the difference in gradients of the two images in the vertical and horizontal directions respectively,representing the gradient of the image I in the vertical direction, +.>Representing the gradient in the horizontal direction, s being a gray scale decay factor, sigma width and σhigh Respectively represent I f Variance in the vertical and horizontal directions. It is clear from the formula (4.9) that if gradient difference is used as the evaluation function of registration, the floating image and the reference image need to have higher similarity in gradient information such as edges, and are not affected by the gray level distribution of the image;
normalized cross-correlation coefficientEvaluation of similarity and correlation of two images in pixel distribution, NC and NCC are formulas (4.10) and (4.11)
wherein ,If,i and Im,i Respectively represent I f and Im Is used for the pixel value of the (i) th pixel, andRespectively represent I f and Im Average pixel value of (a); in the process of calculating the normalized cross-correlation coefficient, the image is required to be subjected to centering treatment to obtain the centered cosine similarity, which is equivalent to the normalized cross-correlation coefficient, and as can be seen from a formula (4.11), the robustness of the similarity measure on gray value difference can be effectively increased by the zero-averaging operation of the data;
the gradient normalized cross-correlation coefficient GradNCC combining the gradient information and the normalized cross-correlation coefficient is the formula (4.12)
wherein , andRespectively representing a gradient image of the ith pixel value of the two-dimensional image I in the vertical direction and a gradient image of the ith pixel value in the horizontal direction, and (2)> andRepresenting average pixel values of the gradient image in the vertical direction and the gradient image in the horizontal direction of the two-dimensional image I, respectively. .
The goal of three-dimensional/two-dimensional registration is to have the post-projection I m With reference image I f More closely, for the final purpose of rigid registration, the result of solving a globally optimal rigid transformation parameter. The process of solving for optima can be generally divided into two types: and resolving and searching. The optimization of analytical solution is generally to solve a linear equation by specific coordinates or parameters to obtain transformation parameters, such as PnP algorithm. The search solution is generally that an initial value is given, and the optimal solution is gradually approximated in a plurality of iterations with a certain step length. Among the local optimization algorithms commonly used in registration algorithms, the evolution strategies (Evolutionary Strategies, ES) are better performing, and the main idea of such methods is to simulate the principle of biological evolution. The optimization process of iteratively solving the optimal result can also be used for three-dimensional/two-dimensional registration problems.
The most common method in evolutionary strategies is covariance matrix adaptive evolutionPolicy (Covariance matrix adaptation evolution strategy, CMA-ES). The three-dimensional/two-dimensional registration problem is a typical nonlinear optimization problem, and because three-dimensional data is required to be mapped to two dimensions in a projection mode, information of a plurality of three-dimensional data is missing in the two-dimensional data, and the change in unit data cannot find rules on the two-dimensional effect of projection, especially in the case of large translation and rotation transformation ranges. There are many challenges in nonlinear optimization problems, such as non-smoothness of the function (non-leadership, multi-modal, noise), excessive number of parameters, variable-to-variable dependencies, etc. These are problems that conventional optimization algorithms cannot solve, but CMA-ES does not rely on gradient optimization compared to other methods, and thus can solve the above-mentioned complex problems. On the other hand, the gradient-based optimization method requires gradient calculation of the objective function, namely While the process of improved DRR is designed to be straightforward, it requires very large computing resources, a condition that existing computing resources are difficult to meet.
Preferably, as shown in fig. 2, in the step (3), the evolution strategy implements the spatial search by iterating the following four steps: (a) sampling a new set of candidate solutions; (b) calculating an optimal value for the new candidate solution; (c) selecting an optimal solution based on f (α); (d) updating the parameters according to the selected solution; in the course of the evolutionary algorithm, a complete iteration comprising (a) - (d) is referred to as a generation, each iteration producing new candidate solutions as children for generating the children as parents; for CMA-ES, the average value m according to the current generation after choosing the optimal point in step (c) (g-1) Calculating covariance matrix C of next generation (g) And at step (4) using the updated average value m (g) 。
Differential rendering has strong rendering capability, and can realize reverse drawing (Inverse Rendering) while guaranteeing rendering effect. The three-dimensional scene information required for generating the image is obtained by reversely drawing the information of the two-dimensional image, wherein the information can be one or more of geometry, lamplight, material and visual angle, and even three-dimensional reconstruction is carried out. The common rendering process is to input a three-dimensional scene into a renderer to obtain a two-dimensional image, and different from the common rendering process, the differential rendering which introduces reverse drawing can link the three-dimensional scene and the two-dimensional image, so that the method has become a popular method in recent years.
In the three-dimensional/two-dimensional registration problem, what needs to be solved is the corresponding three-dimensional shooting view angle when a two-dimensional image is shot, but many existing differentiable rendering methods are realized by performing the grid parameterization of the surface-mount processing on the three-dimensional scene, because the grid parameterization provides a relatively proper compromise between spatial resolution and data volume compared with the volumetric representation. For medical applications, in particular data observed during surgery, are mostly obtained based on a volumetric representation method, and each individual pixel has its specific physical properties. But the mesh data cannot represent the internal features and data of each closed surface, so the differentiable rendering of the study volume datamation process is particularly important.
The spatial transformation network proposed by Google deep allows the traditional convolutional neural network to learn characteristics such as clipping, translation, rotation and scaling effectively. In the spatial transformation network architecture, a feature map or an original image is input to a positioning network, which outputs transformation parameters of the image, and a grid generator and a sampler transform the feature map or the original image using the parameters output from the positioning network. Currently, spatial transformation networks have been used to solve the perspective transformation problem of three-dimensional data reconstruction and the estimation problem of three-dimensional registration deformation fields.
Preferably, in the step (3),
for a two-dimensional image I with the size of MxN, in the DRR process based on Ray-casting, the image is formed by integrating MxN rays through volume data to obtain MxN pixels I (m,n) Is formula (4.13); each ray is emitted by source, andand K control points G are uniformly distributed on each light ray (m,n,k) Representing a continuous integration over the light by discrete integration;
where (m, n) is coordinates on the image, G represents a control point in space;
for the rigid transformation matrix T (alpha), alpha epsilon SE (3), the coordinates of the point source are changed into T (alpha) · [0, 0], the parameter representation of the imaging plane image is changed into T (alpha) · [0, SID,1], and the M multiplied by N multiplied by K control points are subjected to rigid transformation and then subjected to difference calculation, and the transformation result is shown as a formula (4.14):
G′=f Interpolation (T(α)·G) (4.14)
wherein fInterpolation Representing bilinear difference function, G' is control point after rigid transformation, and the transformed image is formula (4.15); the sampling process enables the parameter to carry out back propagation training through the module in the process of carrying out network training; deriving G' yields equation (4.16):
for the three-dimensional/two-dimensional registration problem, the presence of the non-linearity condition results in a large number of extreme points in the objective function constructed based on the conventional similarity measure, which makes the range of the initial rigid transformation parameters of the registration process a large limitation. The ideal objective function is in the shape of a concave function, so that the objective function can be effectively ensured not to be affected by the limitation of the range of the excessively small initial value in the registration process.
In order to obtain a similarity measure applicable to a larger initial parameter range, the present invention solves this problem using a concave function similarity measure model based on the segmentation results. The network training counter-propagation is realized by designing a differential rendering module based on grid sampling, so that in the integral concave function similarity measure, only the three-dimensional convolutional neural network and the two-dimensional image decoder need to be optimized for learning parameters.
Preferably, in the step (4), the three-dimensional volume data is added with a 6-layer three-dimensional convolutional neural network before being input into the differential rendering module based on grid sampling, the three-dimensional convolutional neural network is composed of five convolutional layers, the convolutional kernel of each convolutional layer is 3 multiplied by 3, the step length is (1, 1), after the three-dimensional data are subjected to convolutional neural network to obtain three-channel characteristic images, adding the original images with the obtained characteristic images through one jump connection, outputting three-dimensional data of three channels by a network, wherein the size of each channel is the same as that of the original images, and obtaining four-dimensional data containing image characteristics and the original images after superposition;
for the projected two-dimensional image I m Reference image I f Corresponding feature vectors are obtained through two decoders with the same network structure respectively, and the two decoders do not share parameters in the training process; the decoder is formed by connecting a convolution layer with a convolution kernel of 3 multiplied by 3 and a step length of (1, 1), a residual bottleneck module and a pooling layer, and finally outputting a feature vector through one convolution layer;
obtain I m and If Feature vector E of (2) m and Ef Then, the concave function similarity measure will be E m and Ef As output of the network, this is I m and If Is a similarity measure of (1).
Preferably, in the step (4),
for image similarity measure L in three-dimensional/two-dimensional registration process C =MSE(E m ,E f ) For the concave function, the network uses the Euclidean space geodesic distance with concave function characteristics as the training objective function, as formula (4.17):
L C =MSE(E m ,E f )~L G (α,α f ) (4.17)
wherein LG Represents geodesic distance, alpha, in European space f Represents I f The corresponding rigidity transformation parameter is a gold standard in the registration task; in order to make the two functions more similar in terms of trend, the most straightforward way is to make the gradient of the two functions andAs close as possible, obtaining a target of network training;
during one iteration of the registration, the output is denoted as I C (φ;V m ,I f α), wherein φ represents a learning parameter required by the network; during the back propagation, the networks respectively generate andTwo gradient vectors, the loss function is designed as alpha for L C and LG Two gradient vectors of the two> andThe distance measures of (2) are the formulas (4.18), (4.19) and (4.20), ensuring that the network parameter phi is not updated during back propagation, and the rotation parameter and the translation parameter are calculated respectively;
wherein v represents a vector, trans represents a translation parameter, rot represents a rotation parameter, M D Representing a distance measure;
gradient-based optimization of the rigid transformation parameter α during registration testing by back propagation is achieved as equation (4.21)
Optimizing the parameter alpha by iteration and according to I m Calculate L C To solve for the registration parameters.
It will be understood by those skilled in the art that all or part of the steps in implementing the above embodiment method may be implemented by a program to instruct related hardware, where the program may be stored in a computer readable storage medium, where the program when executed includes the steps of the above embodiment method, and the storage medium may be: ROM/RAM, magnetic disks, optical disks, memory cards, etc. Accordingly, the present invention also includes, corresponding to the method of the present invention, a three-dimensional/two-dimensional image registration device based on skeletal features, which is generally represented in the form of functional modules corresponding to the steps of the method. The device comprises:
a data preprocessing module configured to preprocess the segmented three-dimensional bone data, including region of interest clipping and downsampling;
the simulation data generation module is configured to perform random translation and rotation on the processed TIPS preoperative CT skeleton data to obtain different initial poses, and then perform random translation and rotation again to obtain target poses; projecting the TIPS preoperative CT data of the two poses to generate a large amount of simulation data for subsequent training and testing;
the training module is configured to obtain parameterized three-dimensional characteristic data by processing the three-dimensional data through the three-dimensional convolution network jump connection module at the training stage, and then project the superimposed original image and the three-dimensional characteristic data to two dimensions through the differential rendering module based on grid sampling;
the similarity measure calculating module is used for calculating the similarity measure of the two-dimensional image obtained by projection and the target image according to the two-dimensional image obtained by projection by registration; in order to make the similarity measure obtained by network training have the property of a concave function in European space, a gradient of geodesic distance of a rigid transformation parameter is used as a training target.
The present invention is not limited to the preferred embodiments, but can be modified in any way according to the technical principles of the present invention, and all such modifications, equivalent variations and modifications are included in the scope of the present invention.
Claims (5)
1. The three-dimensional/two-dimensional image registration method based on skeleton features is characterized by comprising the following steps of: which comprises the following steps:
(1) Preprocessing the segmented three-dimensional bone data, including region-of-interest clipping and downsampling;
(2) Carrying out random translation and rotation on the processed TIPS preoperative CT skeleton data to obtain different initial poses, and then carrying out random translation and rotation again to obtain target poses; projecting the TIPS preoperative CT data of the two poses to generate a large amount of simulation data for subsequent training and testing;
(3) In the training stage, three-dimensional data is firstly processed by a three-dimensional convolution network jump connection module to obtain parameterized three-dimensional characteristic data, and then the superimposed original image and the three-dimensional characteristic data are projected to two dimensions by a differential rendering module based on grid sampling;
(4) According to the two-dimensional image obtained by projection, calculating the similarity measure between the two-dimensional image and the target image by using a two-dimensional decoding network; in order to enable the similarity measure obtained by network training to have the property of a concave function in European space, a gradient of geodesic distance of a rigid transformation parameter is used as a training target;
in the step (4), two three-dimensional/two-dimensional registration frames are included:
(I) The test data is subjected to iterative optimization through concave function similarity measurement comprising a differentiable rendering module so as to obtain a rigid transformation parameter corresponding to the target image;
(II) optimizing using a three-dimensional/two-dimensional registration optimizer, implementing a traditional parameterized registration framework;
in the step (3), the step of (c),
the projection to two dimensions is expressed by formula (4.1):
I m =P·T(θ x ,θ y ,θ z ,t x ,t y ,t z )·V m (4.1)
wherein Vm Is three-dimensional volume data, I m Is a corresponding image projected to two dimensions, wherein T (alpha), alpha E SE (3) is a three-dimensional rigid transformation matrix, the formulas are (4.2) - (4.5), P is an expression of a 3 x 4 point source projection matrix, and alpha= (theta) x ,θ y ,θ z ,t x ,t y ,t z ) Refer to the rigidity transformation parameters
wherein θx ,θ y ,θ z The rotation angles around the x coordinate axis, the y coordinate axis and the z coordinate axis are respectively t x ,t y ,t z To translate along three axes, x 0 ,y 0 ,z 0 Is V m Coordinates of the pixel points on the pixel points, x, y, z are x 0 ,y 0 ,z 0 Coordinates of the rigidly transformed points;
in the step (3), the step of (c),
measurement of I f and Im Similarity measure between them, which represents a loss function of formula (4.6)
wherein ,p″CT Representing that three-dimensional data are projected to a two-dimensional camera internal reference matrix, wherein the gradient difference GD and the normalized cross correlation coefficient NCC are similarity measure functions commonly used by a traditional registration method;
gradient differenceIs based on image I m and If The similarity calculated for the gradients of (2) is suitable for processing relatively long and narrow structures, as shown in formulas (4.7), (4.8) and (4.9)
Diff width and Diffhigh Representing two images in the vertical and horizontal directions, respectivelyThe difference in the gradient is such that,representing the gradient of the image I in the vertical direction, +.>Representing the gradient in the horizontal direction, s being a gray scale decay factor, sigma width and σhigh Respectively represent I f Variance in the vertical and horizontal directions;
normalized cross-correlation coefficientEvaluation of similarity and correlation of two images in pixel distribution, NC and NCC are formulas (4.10) and (4.11)
wherein ,If,i and Im,i Respectively represent I f and Im Is used for the pixel value of the (i) th pixel, andRespectively represent I f and Im Average pixel value of (a);
the gradient normalized cross-correlation coefficient GradNCC combining the gradient information and the normalized cross-correlation coefficient is the formula (4.12)
wherein , andRespectively representing a gradient image of the ith pixel value of the two-dimensional image I in the vertical direction and a gradient image of the ith pixel value in the horizontal direction, and (2)> andAverage pixel values of the gradient image in the vertical direction and the gradient image in the horizontal direction of the two-dimensional image I are represented respectively;
in the step (3), the evolution strategy realizes the space search by iterating the following four steps: (a) sampling a new set of candidate solutions; (b) calculating an optimal value for the new candidate solution; (c) selecting an optimal solution based on f (α); (d) updating the parameters according to the selected solution; in the course of the evolutionary algorithm, a complete iteration comprising (a) - (d) is referred to as a generation, each iteration producing new candidate solutions as children for generating the children as parents; for CMA-ES, the average value m according to the current generation after choosing the optimal point in step (c) (g-1) Calculating covariance matrix C of next generation (g) And at step (4) using the updated average value m (g) 。
2. The bone feature-based three-dimensional/two-dimensional image registration method of claim 1, wherein: in the step (3), the step of (c),
for a two-dimensional image I with the size of MxN, in the DRR process based on Ray-casting, the image is formed by integrating MxN rays through volume data to obtain MxN pixels I (m,n) Is formula (4.13); each ray is emitted by source, and eachThe light is uniformly distributed with K control points G (m,n,k) Representing a continuous integration over the light by discrete integration;
where (m, n) is coordinates on the image, G represents a control point in space;
for the rigid transformation matrix T (alpha), alpha epsilon SE (3), the coordinates of the point source are changed into T (alpha) · [0, 0], the parameter representation of the imaging plane image is changed into T (alpha) · [0, SID,1], and the M multiplied by N multiplied by K control points are subjected to rigid transformation and then subjected to difference calculation, and the transformation result is shown as a formula (4.14):
G′=f Interpolation (T(α)·G) (4.14)
wherein fInterpolation Representing bilinear difference function, G' is control point after rigid transformation, and the transformed image is formula (4.15); the sampling process enables the parameter to carry out back propagation training through the module in the process of carrying out network training; deriving G' yields equation (4.16):
3. the bone feature-based three-dimensional/two-dimensional image registration method of claim 2, wherein: in the step (4), a three-dimensional convolutional neural network of 6 layers is added to the three-dimensional volume data before the three-dimensional volume data is input to the differentiable rendering module based on grid sampling, the three-dimensional convolutional neural network is composed of five convolutional layers, the convolutional kernel of each convolutional layer is 3 multiplied by 3, the step length is (1, 1), after the three-dimensional data are subjected to convolutional neural network to obtain three-channel characteristic images, adding the original images with the obtained characteristic images through one jump connection, outputting three-dimensional data of three channels by a network, wherein the size of each channel is the same as that of the original images, and obtaining four-dimensional data containing image characteristics and the original images after superposition;
for the projected two-dimensional image I m Reference image I f Corresponding feature vectors are obtained through two decoders with the same network structure respectively, and the two decoders do not share parameters in the training process; the decoder is formed by connecting a convolution layer with a convolution kernel of 3 multiplied by 3 and a step length of (1, 1), a residual bottleneck module and a pooling layer, and finally outputting a feature vector through one convolution layer; obtain I m and If Feature vector E of (2) m and Ef Then, the concave function similarity measure will be E m and Ef As output of the network, this is I m and If Is a similarity measure of (1).
4. A three-dimensional/two-dimensional image registration method based on bone features according to claim 3, characterized in that: in the step (4) of the above-mentioned method,
for image similarity measure L in three-dimensional/two-dimensional registration process C =MSE(E m ,E f ) For the concave function, the network uses the Euclidean space geodesic distance with concave function characteristics as the training objective function, as formula (4.17):
L C =MSE(E m ,E f )~L G (α,α f ) (4.17)
wherein LG Represents geodesic distance, alpha, in European space f Represents I f The corresponding rigidity transformation parameter is a gold standard in the registration task; in order to make the two functions more similar in terms of trend, the most straightforward way is to make the gradient of the two functions andAs close as possible, obtaining a target of network training;
during one iteration of the registration, the output is denoted as L C (φ;V m ,I f α), wherein φ represents a learning parameter required by the network; during the back propagation, the networks respectively generate andTwo gradient vectors, the loss function is designed as alpha for L C and LG Two gradient vectors of the two> andThe distance measures of (2) are the formulas (4.18), (4.19) and (4.20), ensuring that the network parameter phi is not updated during back propagation, and the rotation parameter and the translation parameter are calculated respectively;
wherein v represents a vector, trans represents a translation parameter, rot represents a rotation parameter, M D Representing a distance measure;
gradient-based optimization of the rigid transformation parameter α during registration testing by back propagation is achieved as equation (4.21)
Optimizing the parameter alpha by iteration and according to I m Calculate L C To solve for the registration parameters.
5. The apparatus of the bone feature-based three-dimensional/two-dimensional image registration method according to claim 1, wherein: it comprises the following steps:
a data preprocessing module configured to preprocess the segmented three-dimensional bone data, including region of interest clipping and downsampling;
the simulation data generation module is configured to perform random translation and rotation on the processed TIPS preoperative CT skeleton data to obtain different initial poses, and then perform random translation and rotation again to obtain target poses; projecting the TIPS preoperative CT data of the two poses to generate a large amount of simulation data for subsequent training and testing;
the training module is configured to obtain parameterized three-dimensional characteristic data by processing the three-dimensional data through the three-dimensional convolution network jump connection module at the training stage, and then project the superimposed original image and the three-dimensional characteristic data to two dimensions through the differential rendering module based on grid sampling;
the similarity measure calculating module is used for calculating the similarity measure of the two-dimensional image obtained by projection and the target image according to the two-dimensional image obtained by projection by registration; in order to make the similarity measure obtained by network training have the property of a concave function in European space, a gradient of geodesic distance of a rigid transformation parameter is used as a training target.
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202110682717.9A CN113450396B (en) | 2021-06-17 | 2021-06-17 | Three-dimensional/two-dimensional image registration method and device based on bone characteristics |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202110682717.9A CN113450396B (en) | 2021-06-17 | 2021-06-17 | Three-dimensional/two-dimensional image registration method and device based on bone characteristics |
Publications (2)
Publication Number | Publication Date |
---|---|
CN113450396A CN113450396A (en) | 2021-09-28 |
CN113450396B true CN113450396B (en) | 2023-05-30 |
Family
ID=77811980
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202110682717.9A Active CN113450396B (en) | 2021-06-17 | 2021-06-17 | Three-dimensional/two-dimensional image registration method and device based on bone characteristics |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN113450396B (en) |
Families Citing this family (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN113920178B (en) * | 2021-11-09 | 2022-04-12 | 广州柏视医疗科技有限公司 | Mark point-based multi-vision 2D-3D image registration method and system |
CN114240740B (en) * | 2021-12-16 | 2022-08-16 | 数坤(北京)网络科技股份有限公司 | Bone expansion image acquisition method and device, medical equipment and storage medium |
CN115619835B (en) * | 2022-09-13 | 2023-09-01 | 浙江大学 | Heterogeneous three-dimensional observation registration method, medium and equipment based on depth phase correlation |
CN115223023B (en) * | 2022-09-16 | 2022-12-20 | 杭州得闻天下数字文化科技有限公司 | Human body contour estimation method and device based on stereoscopic vision and deep neural network |
CN117671162B (en) * | 2024-02-01 | 2024-04-19 | 江苏一影医疗设备有限公司 | Three-dimensional imaging method and device for 4D joint standing position |
CN118429403B (en) * | 2024-07-04 | 2024-09-27 | 湘江实验室 | Image registration method, terminal equipment and medium for periacetabular osteotomy |
Citations (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN107507234A (en) * | 2017-08-29 | 2017-12-22 | 北京大学 | Cone beam computed tomography image and x-ray image method for registering |
Family Cites Families (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US7327865B2 (en) * | 2004-06-30 | 2008-02-05 | Accuray, Inc. | Fiducial-less tracking with non-rigid image registration |
US10818019B2 (en) * | 2017-08-14 | 2020-10-27 | Siemens Healthcare Gmbh | Dilated fully convolutional network for multi-agent 2D/3D medical image registration |
CN110009669B (en) * | 2019-03-22 | 2021-12-10 | 电子科技大学 | 3D/2D medical image registration method based on deep reinforcement learning |
CN112509022A (en) * | 2020-12-17 | 2021-03-16 | 安徽埃克索医疗机器人有限公司 | Non-calibration object registration method for preoperative three-dimensional image and intraoperative perspective image |
CN112614169B (en) * | 2020-12-24 | 2022-03-25 | 电子科技大学 | 2D/3D spine CT (computed tomography) level registration method based on deep learning network |
-
2021
- 2021-06-17 CN CN202110682717.9A patent/CN113450396B/en active Active
Patent Citations (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN107507234A (en) * | 2017-08-29 | 2017-12-22 | 北京大学 | Cone beam computed tomography image and x-ray image method for registering |
Also Published As
Publication number | Publication date |
---|---|
CN113450396A (en) | 2021-09-28 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN113450396B (en) | Three-dimensional/two-dimensional image registration method and device based on bone characteristics | |
Xie et al. | Neural fields in visual computing and beyond | |
Meuleman et al. | Progressively optimized local radiance fields for robust view synthesis | |
US10740897B2 (en) | Method and device for three-dimensional feature-embedded image object component-level semantic segmentation | |
CN111902825A (en) | Polygonal object labeling system and method for training object labeling system | |
CN112215050A (en) | Nonlinear 3DMM face reconstruction and posture normalization method, device, medium and equipment | |
CN113516659A (en) | Medical image automatic segmentation method based on deep learning | |
Jiang et al. | H $ _ {2} $-Mapping: Real-time Dense Mapping Using Hierarchical Hybrid Representation | |
Cornells et al. | Real-time connectivity constrained depth map computation using programmable graphics hardware | |
WO2022198684A1 (en) | Methods and systems for training quantized neural radiance field | |
CN113177592B (en) | Image segmentation method and device, computer equipment and storage medium | |
Li et al. | Coarse-to-fine PatchMatch for dense correspondence | |
US12136230B2 (en) | Method for training neural network, system for training neural network, and neural network | |
Rabby et al. | BeyondPixels: A comprehensive review of the evolution of neural radiance fields | |
Zhou et al. | NeRFLix: High-quality neural view synthesis by learning a degradation-driven inter-viewpoint mixer | |
Wang et al. | Super-resolution of multi-observed RGB-D images based on nonlocal regression and total variation | |
Nakashima et al. | Lidar data synthesis with denoising diffusion probabilistic models | |
CN112967373B (en) | Facial image feature coding method based on nonlinear 3DMM | |
Kläser et al. | Deep boosted regression for MR to CT synthesis | |
Liu et al. | Video frame interpolation via optical flow estimation with image inpainting | |
Peng et al. | PDRF: progressively deblurring radiance field for fast scene reconstruction from blurry images | |
Tan et al. | High dynamic range imaging for dynamic scenes with large-scale motions and severe saturation | |
CN115965765A (en) | Human motion capture method in deformable scene based on neural deformation | |
CN115439669A (en) | Feature point detection network based on deep learning and cross-resolution image matching method | |
Tran et al. | Self-supervised learning with multi-view rendering for 3d point cloud analysis |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
CB03 | Change of inventor or designer information | ||
CB03 | Change of inventor or designer information |
Inventor after: Fu Tianyu Inventor after: Yang Jian Inventor after: Fan Jingfan Inventor after: Song Hong Inventor after: Hu Guoyu Inventor before: Li Wenjie Inventor before: Ai Danni Inventor before: Yang Jian Inventor before: Fan Jingfan Inventor before: Song Hong |
|
GR01 | Patent grant | ||
GR01 | Patent grant |