Coronary artery motion modeling
A method for tracking coronary artery motion includes constructing ( 11 ) a centerline model of a vascular structure in a base phase image in a sequence of 2D images of coronary arteries acquired over a cardiac phase, computing ( 12 ), for each pixel in a region-of-interest in each subsequent image, a velocity vector that represent a change in position between the subsequent image and base phase image, calculating ( 13 ) positions of control points in each phase using the velocity vectors, and applying ( 14 ) PCA to a P×2N data matrix X T constructed from position vectors (x, y) of N centerline control points for P phases to identify d eigenvectors corresponding to the largest eigenvalues of XX T to obtain a d-dimensional linear motion model {circumflex over (α)} p , in which a centerline model for a new image at phase p+1 is estimated by adding {circumflex over (α)} p to each centerline control point of a previous frame at phase p.
1. A method of tracking coronary artery motion in a sequence of 2D images, comprising the steps of:
constructing a centerline model of a segmented vascular structure in a base phase image in a sequence of 2D images of coronary arteries acquired over one or more cardiac phases;
computing, for each pixel in a region-of-interest (ROI) in each subsequent image in the sequence of 2D images, a velocity vector that represent a change in position between the subsequent image and base phase image in the sequence of 2D images;
calculating positions of control points i in each phase p using the velocity vectors v i T =(v i,x , v i,y ) T ; and
applying PCA to a P×2N data matrix X T constructed from position vectors (x, y) of N centerline control points {V 0 (p) , V 1 (p) , . . . , V N−1 (p) } for P phases to identify d eigenvectors corresponding to the largest eigenvalues of XX T to obtain d-dimensional linear motion model
α
^
p
←
m
_
p
+
∑
i
=
1
d
α
p
,
i
v
i
,
where m p is a mean value of 2N control points positions (x, y) at phase p of the cardiac cycle, v i is an eigenvector of U=cov(XX T ), and the α p,i 's are determined for each phase p, wherein a centerline model for a new image at phase p+1 is estimated by adding {circumflex over (α)} p to each centerline control point of a previous frame at phase p.
2. The method of claim 1 , wherein said centerline model is parametrically represented by a set of vessel segments connected by a set of control points, wherein each vessel segment is approximated by a 2D B-spline curve parameterized by chord length.
3. The method of claim 2 , wherein the 2D B-spline curve for a k th centerline segment is defined as
C
k
(
u
)
=
∑
i
=
0
n
k
-
1
B
i
,
q
(
u
)
V
i
k
,
wherein uε[0,1], the {v 1 k }εR 2 are n k control points of the centerline and the {B i,q } are q th degree B-spline basis functions, wherein end values of each curve C k are constrained so that
C k (0)= V 0 k =( s 0 εS k ), C k (1)= V n k k =( s m k εS k )
wherein S k denotes a k th segment, and S m k denotes an m th point in the k th segment.
4. The method of claim 1 , wherein constructing said centerline model comprises processing said base phase image with a Hessian-based vessel enhancement filter, and computing centerlines by numerically integrating a directional vector field obtained from the Hessian-based vessel enhancement filter.
5. The method of claim 4 , further comprising receiving manual adjustments to the centerline model.
6. The method of claim 1 , wherein computing a velocity vector comprises calculating v i =− S i · Q −1 ·q, wherein V i =(v x i , v y i ) is the velocity vector for control point i,
Q
=
∑
i
Q
i
=
∑
i
S
i
T
T
i
S
i
=
(
Q
_
q
q
T
α
)
,
S
i
=
(
x
i
y
i
1
0
0
0
0
0
0
0
x
i
y
i
1
0
0
0
0
0
0
0
1
)
,
S
_
i
=
(
x
i
y
i
1
0
0
0
0
0
0
x
i
y
i
1
)
,
and T i is an isotropy compensated orientation tensor for control point i.
7. The method of claim 6 , wherein the isotropy compensated orientation tensor T i is calculated by stacking all images of the sequence of 2D images onto each other to form a 3D image f, and defining T i as T i =A i A i T +γb i b i T −λ i , min I, wherein I is the identity matrix, λ min is the smallest eigenvalue of {tilde over (T)} i =A i A i T +γb i b i T , A i and b i are found by fitting the image intensity f(x) of the control points x i to a second degree polynomial f(x)˜x T Ax+b T x+c, and γ is a non-negative weight factor between even and odd parts of the image.
8. The method of claim 1 , wherein each image in the sequence of 2D images of coronary arteries is mapped to a cardiac phase by a recorded electrocardiogram signal.
9. The method of claim 1 , further comprising performing reverse-PCA analysis to estimate a vector {circumflex over (m)} p of centerline control points at cardiac phase p from PCA parameter vector {circumflex over (α)} p .
10. The method of claim 1 , wherein the ROI in the subsequent image is determined by projecting a region around the centerline in the base phase image onto the subsequent image, and determining a centerline in the subsequent image by registering the centerline in the base phase image to the ROI in the subsequent image.
11. A non-transitory program storage device readable by a computer, tangibly embodying a program of instructions executed by the computer to perform the method steps for tracking coronary artery motion in a sequence of 2D images, the method comprising the steps of:
constructing a centerline model of a segmented vascular structure in a base phase image in a sequence of 2D images of coronary arteries acquired over one or more cardiac phases;
computing, for each pixel in a region-of-interest (ROI) in each subsequent image in the sequence of 2D images, a velocity vector that represent a change in position between the subsequent image and base phase image in the sequence of 2D images;
calculating positions of control points i in each phase p using the velocity vectors v i T =(v i,x , v i,y ) T ; and
applying PCA to a P×2N data matrix X T constructed from position vectors (x, y) of N centerline control points {V 0 (p) , V 1 (p) , . . . , V N−1 (p) } for P phases to identify d eigenvectors corresponding to the largest eigenvalues of XX T to obtain d-dimensional linear motion model
α
^
p
←
m
_
p
+
∑
i
=
1
d
α
p
,
i
v
i
,
where m p is a mean value of 2N control points positions (x, y) at phase p of the cardiac cycle, v i is an eigenvector of U=cov(XX T ), and the α p,i 's are determined for each phase p, wherein a centerline model for a new image at phase p+1 is estimated by adding {circumflex over (α)} p to each centerline control point of a previous frame at phase p.
12. The computer readable program storage device of claim 11 , wherein said centerline model is parametrically represented by a set of vessel segments connected by a set of control points, wherein each vessel segment is approximated by a 2D B-spline curve parameterized by chord length.
13. The computer readable program storage device of claim 12 , wherein the 2D B-spline curve for a k th centerline segment is defined as
C
k
(
u
)
=
∑
i
=
0
n
k
-
1
B
i
,
q
(
u
)
V
i
k
,
wherein uε[0,1], the {v i k }εR 2 are n k control points of the centerline and the {B i,q } are q th -degree B-spline basis functions, wherein end values of each curve C k are constrained so that
C k (0)= V 0 k =( s 0 εS k ), C k (1)= V n k k =( s m k εS k )
wherein S k denotes a k th segment, and S m k denotes an m th point in the k th segment.
14. The computer readable program storage device of claim 11 , wherein constructing said centerline model comprises processing said base phase image with a Hessian-based vessel enhancement filter, and computing centerlines by numerically integrating a directional vector field obtained from the Hessian-based vessel enhancement filter.
15. The computer readable program storage device of claim 14 , the method further comprising receiving manual adjustments to the centerline model.
16. The computer readable program storage device of claim 11 , wherein computing a velocity vector comprises
calculating v i =− S i · Q −1 ·q, wherein V i =(v x i , v y i ) is the velocity vector for control point i,
Q
=
∑
i
Q
i
=
(
Q
_
q
q
T
α
)
,
S
i
=
(
x
i
y
i
1
0
0
0
0
0
0
0
x
i
y
i
1
0
0
0
0
0
0
0
1
)
,
S
_
i
=
(
x
i
y
i
1
0
0
0
0
0
0
x
i
y
i
1
)
,
and T i is an isotropy compensated orientation tensor for control point i.
17. The computer readable program storage device of claim 16 , wherein the isotropy compensated orientation tensor T i is calculated by stacking all images of the sequence of 2D images onto each other to form a 3D image f, and defining T i as T i =A i A i T +γb i b i T −λ i , min I, wherein I is the identity matrix, λ min is the smallest eigenvalue of {tilde over (T)} i =A i A i T +γb i b i T , A i and b i are found by fitting the image intensity f(x i ) of the control points x i to a second degree polynomial f(x)˜x T Ax+b T x+c, and γ is a non-negative weight factor between even and odd parts of the image.
18. The computer readable program storage device of claim 11 , wherein each image in the sequence of 2D images of coronary arteries is mapped to a cardiac phase by a recorded electrocardiogram signal.
19. The computer readable program storage device of claim 11 , the method further comprising performing reverse-PCA analysis to estimate a vector {circumflex over (m)} p of centerline control points at cardiac phase p from PCA parameter vector {circumflex over (α)} p .
20. The computer readable program storage device of claim 11 , wherein the ROI in the subsequent image is determined by projecting a region around the centerline in the base phase image onto the subsequent image, and determining a centerline in the subsequent image by registering the centerline in the base phase image to the ROI in the subsequent image.
21. The computer readable program storage device of claim 11 , the method further comprising determining model parameters {right arrow over (α)}=(α 1 , . . . , α d ) by minimizing
∑
k
=
1
2
C
k
(
I
k
,
Y
0
,
M
d
(
α
→
)
)
wherein C k (I k , Y 0 , {right arrow over (m)}) is a similarity measure between a 2D image I of camera k and a preoperative 3-D centerline segmentation Y 0 transformed according to a rigid mapping M d and projected into a camera coordinate frame.
22. The computer readable program storage device of claim 11 , the method further comprising fitting polynomial curves for each of the parameters α i and the position/orientation of a sensor that reports its position and orientation using
α
p
,
i
=
∑
j
=
1
6
P
p
,
j
q
(
S
j
)
,
wherein p is the cardiac phase, P is a q th order polynomial functions, and S j are the position and orientation of the sensor.