Abstract

I. INTRODUCTION

POINT tracking over an image sequence is an essential problem which is common to several domains such as computer vision and biomedical imaging and is involved in many high level tasks such as motion estimation, 3-D reconstruction, surveillance or robot vision. Point trajectory along a given image sequence when non-linear motion is involved, is a pretty-quite difficult task. Unlike structured shape tracking, no shape priors can be imposed and the only possibility is to rely on local characteristics of the point. In addition, tracking a given point over time assumes the invariance of some local features along the trajectory. Another difficulty is often related to the setting of a prior dynamic model of the point motion when no prior knowledge on the evolution of the surrounding of the objects is available. Hence methods reported in the literature often rely on local techniques based on geometric and luminance invariants. One of the most popular is based on correlation techniques [1] but they turn out to be quickly inadequate in the case of large motion, intensity variation or occlusions. Optical flow methods [2], based on assumption of the conservation of intensity properties throughout their temporal evolution are considered well performing for small motion velocities. The brightness constancy constraint equation often requires some averaging in the neighbourhood to reduce the variance of the motion vector estimate. The feature-based optical flow estimation matches small patches over images but its performances may be also limited when the image background is locally highly structured or when dealing with strong non-transactional motions and deformations of the 3-D objects. Techniques using block matching (based on a similarity measure) operate better for larger displacements [3].

The application, we are dealing with, concerns the 3-D temporal tracking of coronaries in multi-slices computed tomography angiography (MSCTA) and is related to the real time observation of the catheter trajectory in interventional surgery. The subject copes with dynamic volume sequence analysis and makes use of a voxel tracking strategy to locally obtain the displacement of a vessel on a cardiac cycle. The proposed method is based on a region matching with a distance criterion that combines geometric moment-based descriptors with intensity information: It is described section II. Preliminary results on simulated and real data are then presented section III. A discussion is finally provided section IV with some perspectives on the work in progress.

II. VOXEL TRACKING PROCESS

The movement of coronaries is directly related to the cardiac beating which induces non-linear motion and deformations of the structures. Because geometric moments are integral-based features and are therefore robust against noise, we used them in association with local intensity information to track a voxel in the sequence. We first locally modeled the vessel by a cylinder of center of gravity P_i, diameter d_i and orientation α_i,- and β_i,- in the 3-D space. Analytical expressions of these parameters were obtained from the 3-D geometric moments of up to order 2 [4]. Central moments were computed on an isotropic grid using a spherical window. Each pixel inside this window was weighted according to its memberships to the sphere, by considering a sub-pixel decomposition. This way, the computation was reduced to mask convolutions.

Thus the tracking process was performed into three steps:

1. Estimation of the local characteristics of the vessel around the considered voxel at time t

2. Region matching based on a distance combining geometric moment descriptors with vessel and background intensities to recover the vessel at time t+1.

3. Adjustment of the position of the voxel inside the vessel using a spatial tracking.

B. Temporal tracking with region matching

Let us consider a given point P_i,t in the volume t. The problem of tracking this feature in a sequence can be defined as estimating the position of the point P_i,t in the volume t+1. We considered that no knowledge on the dynamic of the surrounding object was available. We tried thus to solve the correspondence through a region matching based on a criterion of minimal distance associating both geometric and density information: the central moments up to order 2 and the intensities of the vessel and background, computed for each candidate voxel P_r,_t+1. The maximum displacement magnitude that was observed on a coronary branch, has been evaluated to about 1 cm between two temporal instants, in two recent studies [5–6]. We designed thus a spherical search space of radius Rs = 1 cm, centered on the position of the projected voxel P_i,t, in the volume t+1. This space was then split into N overlapping regions whose center positions were chosen to optimize the covered volume by a set of spheres of radius R equal to the estimated diameter of the vessel at time t (Fig. 1).

Open in a separate window

Fig. 1

Search space for a vessel diameter estimated to R = 8. This space includes 81 spheres

For each spherical region r of radius R, centered at location P_r,t+1, r = 1 …N, we applied the multiresolution moment algorithm to make the center of the window converge towards the brightest structures (cardiac cavities, myocardium, arteries) located in the surrounding of P_r,_t+1. We assumed then that at least one of the spheres was attracted by the vessel we were looking for. The central moments as the vessel and background intensities were then computed for each sphere. We evaluated then the distance between the point P_i,t and the set of points $P_{r,t+1}^G$, characterizing the center of gravity of each sphere (located at the center of the sphere) to select the point likely to belong to the same vessel, in the volume t+1:

$$D_{t,t+1}^{i,r}={\left[\begin{array}{l}{\left(\sum _{p,q,l=0}^2{\alpha }_{p+q+l}{\mu }_{\mathit{pql}}(P_{i,t})-{\alpha }_{p+q+l}{\mu }_{\mathit{pql}}(P_{r,t+1}^G)\right)}^2\\ +{\beta }^2\times {\left(I_v(P_{i,t})-I_v(P_{r,t+1}^G)\right)}^2\\ +{\gamma }^2\times {\left(I_b(P_{i,t})-I_b(P_{r,t+1}^G)\right)}^2\end{array}\right]}^{\frac{1}{2}}$$

(1)

where and ${\mu }_{\mathit{pql}}(P_{r,t+1}^G)$ (p,q,l = 1,2; t = 1,2, …,10; r = 1,2, …, N) represent the 3-D central moments computed at the voxel P_i,t in the volume t and $P_{r,t+1}^G$ the center of gravity of each candidate sphere r in the volume t+1.;respectively. The region that provided the minimal value $D_{t,t+1}^{i,r}$, was retained and its center of gravity $P_{r,t+1}^G$ was assimilated to the position of the voxel P_i,t in the volume t+1:

$${\widehat{P}}_{i,t+1}=\left\{P_{r,t+1}^G\mid r=\underset{r=1\dots N}{argmin}{D}_{t,t+1}^{i,r}\right\}$$

(2)

III. RESULTS

The algorithm was tested on two dynamic cardiac sequences. One was acquired on a Siemens Somatom volume zoom 4 detectors rings and the other on a GE LightSpeed system with 16 detectors rings. These sequences included 10 volumes reconstructed from slices acquired every 10% of the cardiac cycle and on several cycles. The dimension of each volume was 512×512×219 and 512×512×300 slices with a voxel resolution of 0.35×0.35× 0.5 mm³ and 0.47×0.47×0.6 mm³ respectively. A preliminary interpolation was performed to make the datasets isotropic in each direction.

The arteries were divided into proximal, middle and distal segments. Three points were selected on the right coronary: RCA origine, RCA_PM (the point of intersection between the proximal and middle thirds) and RCMMD (the point of intersection between the middle and distal thirds). In the same way, we identified seven points on the left coronary: one on the Left main coronary (LMC), three on the proximal, middle and distal left anterior descending artery (LAD) respectively, one on the left circumflex artery (LCx), one on the second Diagonal (SDA) and one on the Obtuse Marginal artery (PMA).

A first stage aimed at evaluating the behavior of the algorithm with a simulated displacement of the structures. We applied thus a spatial translation of growing magnitude up to 12 mm, on several of the volumes of each sequence. We used then a greedy algorithm to evaluate the tracking results for different values of the weighting coefficients including in the computation of the distance $D_{t,t+1}^{i,r}$. We obtained the best results for α₀ = 1, α₁ = 0, α₂ = 7 and β = 1, γ = 7. We computed each time the mean position error on the set of the 10 tracked points. Fig. 2 illustrates the steps of the tracking process for a spatial translation of magnitude 10 mm, i.e. 16 voxels in each direction x,y and z, on one of the volumes of one sequence. Results are shown on MIP images (Maximum Intensity Projections) computed in the axial plane.

Open in a separate window

Fig. 2

Tracking process applied on the second diagonal artery, a) A lozenge shows the position of a voxel P_u in the volume t. b) Projection of the voxel P_i,t in the volume t+1 and search space centered on the projected voxel P_i,t+1. The spherical space includes 683 spheres of radius R = 4. c) Positions of the centers of the spheres after the temporal region matching. The selected point P_u+i appears in red. d) The lozenge depicts the final position of the estimated point P_i,t+1. Spatial tracking was not useful there since the estimated point was found very close in the sense of the criterion.

The mean position error, computed for a given translation of a volume for each of the sequence, is given Fig. 3. The first stage allows the right vessel to be localized. Nevertheless, if the local characteristics of the vessel, at the estimated point _i,t+1, are relatively far from those measured at the point P_i,t, the spatial tracking, inside the vessel, aimed to refine this position. We have then applied the algorithm to search for the trajectory of the 10 points on the two dynamic volume sequences. We computed also this trajectory on two neighbors located up and down each point to derive a mean trajectory (in order to deal with the noise and motion artifacts). We traced then the mean displacement of each vascular segment (Fig. 4). Fig. 5 illustrates a result of the tracking on the origin of the RCA.

Open in a separate window

Fig. 3

Mean displacement error for a simulated displacement on a volume of each sequence (Siemens and GE): This error has been computed for the 10 points located on the different branches

Open in a separate window

Fig. 4

Mean displacement of the RCA origin during a cardiac cycle. The plot shows results on two points situated at the origin of the RCA (thin lines) and the mean displacement (thick line), a–c) decomposition of the displacement into three 1-D displacements: left/right, post/anterior and inferior/superior, d) 3-D displacement

Open in a separate window

Fig. 5

Temporal tracking of a voxel located at the origin of the RCA. Results are shown on MIP images. The lozenge indicates the estimated position at time t+1 and the star, the position of the point at time t.

IV. DISCUSSION

The performances of the algorithm depend on both the characteristics of the images and the structural complexity of the scene. The images are spatially anisotropic (a linear interpolation between slices need to be applied), not corrected during the reconstruction in term of motion artefacts and the vessel contrast evolves in time. The scene includes different structures such as the myocardium, the cavities and the venous and arterial trees, which are very close together and have a similar contrast. The size of a coronary branch evolves from the proximal to the distal segment and the shape can vary (due to the pathology).

The errors in the estimation of the displacement include two components: the motion of the coronaries but also a spatial translation along the vessel because the matching may not be exact between the two estimated models. Background information and the local orientation have been exploited to differentiate the vessels having similar characteristics. The first temporal process always allows finding the right vessel. Nevertheless, the estimated point is not always well positioned. The second process which corresponds to a tracking process inside the vessel aims at refining the position of the point. Tests on simulated motion have provided a mean error less than 2 mm. This simulation doesn’t include change of contrast with time, the deformation of the structures and motion artefact that were often present in real data. Thus the next stage will be to take into account these parameters to better analyse the limit of the algorithm. Preliminary results on real data appear consistent with two studies reported in the literature [5] [6] on the coronary motion. The computation time to track a point between two volumes takes less than 2 seconds.

V. CONCLUSION

We have proposed a method to track coronaries in dynamic volume sequences in MSCTA. The results on data with linear simulated motion showed that the method works despite the complexity of the structures inside the volume. The first results on real data confirmed these promising results.

Acknowledgments

References