Automatic extraction of vessel centerlines is viewed as an essential process in majority of image-guided diagnostic and therapeutic applications. Among a considerable number of methods, direct exploratory tracing method is known as an efficient solution for reliable extraction of vessel features from two-dimensional fluoroscopic images. The first step of most automatic exploratory tracing algorithms is to collect some candidate initial seed points as well as their initial tracing directions. To detect reliable initial points, a validation procedure is required to filter out the false candidates and avoid unnecessary tracing. Starting from reliable initial points, the algorithm efficiently extracts the centerline points along the initial direction until certain pre-defined criteria are met. However, most of these algorithms suffer from incomplete results due to inappropriate selection of the initial seed points. The conventional seed point selection algorithms either rely merely on signal-to-noise ratio analysis, which results in a large number of false traces, or impose a set of strict geometrical validation rules which leads to more false negatives and as a consequence more time shall be spent on computation. This paper presents a new method for efficient selection of initial points for exploratory tracing algorithms. The proposed method improves the performance upon existing methods by employing a combination of geometrical and intensity-based approaches. Moreover, it provides a tunable trade-off between the strictness of the validation procedure and computational efficiency. The results of comparative performance with other proposed techniques are included.

Keywords:    Centerline extraction,   coronaries,   exploratory tracing,   seed point detection,   tracking

Imperfect nature of visual assessment for estimating the geometric measures of arterial lesions has led to the emersion of computerized techniques for automated quantitative analysis of vessels in fluoroscopic images. An important step in these techniques is the extraction of vessel centerlines out of an individual image. Through centerline extraction process, other vessel features such as border line positions, vessel diameter profiles, video-densitometric profiles and the length of the stenotic region can automatically be calculated^1-3.

The existing centerline extraction methods can be categorized into interactive, semi-automatic and fully- automatic techniques. In the former techniques^4-6, a set of centerline positions is visually identified by the user clicks. This type of identification takes advantage of highly intelligent human vision in detecting true vessel segments and avoiding background traces. In the second technique, an individual initial point and an initial direction are manually specified for each vessel segment and the extraction algorithm automatically tracks the centerline over the entire artery^7,8. However, in addition to high levels associated with manual estimation of vessel orientation⁹, both techniques are far from ideal for automated image analysis tasks. In fully-automatic techniques, the need for user-interaction is eliminated and reproducibility is enhanced by detecting initial seed points located inside the boundaries of the vessels which define the vessel to be tracked as well as the starting points for centerline extraction algorithm.

Among many automatic centerline extraction approaches, the vessel tracing approach^10-18 holds several attractive properties that facilitates real-time and intra-operative image processing tasks. They have found applications in medical and biomedical areas such as in image guided coronary intervention, computer-assisted laser retinal surgery, neurite outgrowth detection and many others. These algorithms are capable of robustly extracting the linear features in poor quality images with varying contrast and artifacts as well as handling the partial views, i.e., the vessels are not required to be necessarily connected. In addition, they are computationally efficient because they avoid low level image-wide processing and only process the pixels that are close to the vasculature.

The tracing algorithms start with collecting a number of starting points by performing a search for potential starting points (candidate points) along the lines of a sparse grid over the image. Since some of the collected points may be pertinent due to noise, the algorithm must verify the collected points against a set of validation criteria. In the next step, the tracing algorithm starts from an initial point and recursively extracts the center line points of the blood vessels until the complete artery is traced. This process is repeated for every seed point and the traces are combined to constitute the final result.

However, most of the vessel tracing algorithms produce incomplete tracing results when inappropriate parameter values are selected for seed point detection algorithm¹⁸. In fact, the quality of tracing results depends on the performance of the seed point collection and seed point validation in the first stage of the tracing algorithm. Even though increasing the number of grid lines contributes to the number of detected seed points and more complete traces can be achieved, it decreases the computational efficiency of the algorithm. On the other hand, there is a tradeoff between the strictness of the validation algorithm and the tracing performance. Imposing strict validation rules decreases the number of background traces, but results in more false negatives relative to total validated seed points. In order to explain the motivation of the current work, a review of some related prominent studies is presented in the next section.

Collorec and Coatrieux¹³ proposed a recursive algorithm for automatic tracing of vascular network in digital subtraction angiography. In this method the seed points are detected by seeking local maxima intensity values over a coarse grid. After collecting the seed points, each seed point is validated by comparing average gray values of the pixels along different directions with a global threshold. The results achieved using this merely intensity based validation algorithm suffer from relatively large number of false detections.

To overcome this problem, Can¹², et al. proposed a fast and effective method based on the utilization of 2-D differentiator kernels in seed point validation and tracing stages to trace the vessel center lines in ophthalmic images. In the first step, the seed points are collected by exploring the image along a grid one-pixel-wide lines, seeking out local gray-level minima. In the second step, false seed points are filtered out by a set of strict validation rules. This method performs well in rejecting the seed points that are not located on or near any center lines but due to the strictness of the validation rules, the rate of false negatives is relatively high. Moreover, this verification requires much computation time, since for each candidate seed point the existence of a pair of strong parallel edges around the seed point must be tested. This is performed by applying a set of directional kernels to the seed points’ neighboring points along many directions (16 line searches for left and 16 line searches for right kernels).

This idea was then used for seed point detection in many subsequent related papers^10,11,16,18. In the study conducted by Al-Kofahi¹⁰, et al. adaptive 2-D kernels were employed for automatic tracing of neurons from noisy digital confocal microscope images. This enabled their algorithm to trace over apparent discontinuities. They showed that using several 1-D differentiators as a variable length kernel instead of fixed 2-D correlation operators can improve the robustness of the tracing algorithm. In addition, more robustness was achieved by simultaneously detecting both edges of the linear structures.

Al-Kofahi¹¹, et al. used adaptive intensity thresholding for further filtering the candidate points that are due to noise or imaging artifacts. Similarly, in the work of Zhang¹⁸, et al. the reliable candidates are chosen based on the signal-to-noise ratio analysis. They improved earlier methods by utilizing global thresholding and local signal-to-noise ratio analysis to detect reliable seed points in microscopy neuron images. They showed that their method outperforms the methods proposed¹² and balancing between the number of false negatives and false positives¹³. Nevertheless, their method does not take into account the geometrical features of the linear objects and would have a large number of false detections when applied to other image modalities with complex vessel curvatures and imaging artifacts (e.g. coronary angiograms). In addition, their method suffers from the high computational cost of calculating global thresholds over the whole image.

These problems have motivated us to develop a flexible and parametric method that relies on direct detection of valid boundary points and exploiting their location to estimate the location of centerline seed points. The proposed method efficiently validate the seed points using directional features of gradient vectors calculated at the seed point and its neighbors that are located on the same vessel boundary. We present results of applying our method to extract coronary vessels from 2-D X-ray angiogram images (512 × 512, 8-bit). The method exhibits improvements upon a well known prior method¹² in reducing the computational burden and enhancing the performance of the seed point detection algorithm.

In this section we present those aspects of the tracing approach necessary to compare our seed point detection and validation algorithms with those proposed in other methods. We do not discuss other steps of the tracing algorithm which are previously described in the literature, such as recursive tracing of the vessels and branch/crossover point detection.

It should be observed that, although the main objective of the proposed method is to detect reliable center line seed points, the seed points are not directly collected and validated. Instead, the algorithm first detects the valid boundary points which are only meant for selecting the right region of interest, and not for using as input for centerline extraction algorithm. In the next step, the validated left and right boundary points are used for calculating the location of desired center line seed points. The details of the algorithm are explained in the following section.

3.1 Automatic Collection of Boundary Points

In this step a number of candidate boundary points with high probability of being on the vessel boundaries are collected. For computational efficiency, the image is sampled by a coarse grid of 2N lines (N horizontal and N vertical lines) as shown in Fig 1(a). Since the candidate points are supposed to be located on the edges of the linear feature, they can be identified by convolving the gray-level values on each line of the grid with first derivate of one-dimensional Gaussian low-pass filter g′_σ:

${g^{\prime }}_{\text{σ}}=\frac{d}{dx}\left(\frac{1}{\sqrt{2\mathrm{\pi }\mathrm{\sigma }}}{e}^{-\frac{x^2}{2{\mathrm{\sigma }}^2}}\right)$ (1)

For computational efficiency, instead of using the continuous kernel in Eqn (1), a discrete approximation of the form [1, 2, 0, -2, -1]^T is used. By convolving the 1-D kernel, local peak values of the filter response are identified on each line within a small neighborhood of size Ns. The value of Ns determines the number of collected boundary points and affects the performance of the seed point detection process. Figure 1(b) shows a small region of an angiogram with white points indicating the location of collected boundary points. Selecting a small value for Ns improves the detection of both boundary points in small vessels but results in collecting more undesirable edge points and increases the computational burden [shown in Fig 2(d)]. Therefore, the optimal value for neighborhood size Ns is chosen to be equal or greater than half of the maximum expected vessel width M. The number of horizontal and vertical grid lines, denoted by N, is also determined such that the number of detected boundary points as a function of the grid size N is maximized. As shown in Figs 2(a) and 2(b), larger numbers of grid lines increase the number of seed points at the expense of computational time required for boundary point selection and validation processes. In order to tradeoff between the computational time and the number of grid lines, an overall performance measure, denoted by P, is defined as the ratio of the number of true boundary points and computational time. Figure 2(c) plots the average of this ratio for 10 angiograms, exhibiting a broad peak between 8 and 13. This graph indicates that for N values in this range, the cost-performance tradeoff is reasonable. For the computational resources available to us, N was chosen as 13, which yields an average computational time of about 420 ms for seed point collection and validation.

3.2 Boundary Point Validation

Since our method seeks out the boundary points in the previous step, the validation procedure is only required to check the symmetry of linear features in a particular neighborhood of the candidate points and their corresponding points located on the opposite edge. After all the boundary points are validated, the locations of the actual center line seed points are estimated. To examine the symmetry of the linear features, the algorithm uses the directional properties of the gradient vectors around the vessel borders. In the mathematical point of view, the gradient vector is perpendicular to the surface of the level curves. This property can be used for discarding the points that have neighbors with different directions. Figure 3 illustrates a small part of a gradient vector field calculated for a sample linear object.

This figure depicts the following features of the gradient vectors located on or near the borders of the linear feature:

1. The gradient vectors of the points that are located on the same boundary within a small distance have nearly equal directions.
2. The gradient vectors calculated at each boundary point and its diametrically opposite point located at the opposite boundary must have nearly opposite directions.

Before discussing the details of the validation process, we shall mention that due to low signal to noise ratio and since vessel boundaries might cut across the pixel array at any angle, two 5 × 5-gradient masks are used for calculating the x and y gradients centered at point p(x, y). The gradient components ∇_x and ∇y are calculated by correlation of the pixel values with the following masks:

${\nabla }_x=\left(\begin{array}{lllll}-1\hfill & -2\hfill & 0\hfill & 2\hfill & 1\hfill \\ -1\hfill & -2\hfill & 0\hfill & 2\hfill & 1\hfill \\ -1\hfill & -2\hfill & 0\hfill & 2\hfill & 1\hfill \\ -1\hfill & -2\hfill & 0\hfill & 2\hfill & 1\hfill \\ -1\hfill & -2\hfill & 0\hfill & 2\hfill & 1\hfill \end{array}\right){\nabla }_y=\left(\begin{array}{lllll}1\hfill & 1\hfill & 1\hfill & 1\hfill & 1\hfill \\ 2\hfill & 2\hfill & 2\hfill & 2\hfill & 2\hfill \\ 0\hfill & 0\hfill & 0\hfill & 0\hfill & 0\hfill \\ -2\hfill & -2\hfill & -2\hfill & -2\hfill & -2\hfill \\ -1\hfill & -1\hfill & -1\hfill & -1\hfill & -1\hfill \end{array}\right)$ (2)

This estimation is neither unique nor optimal. In fact, more optimal estimations were avoided in favor of suboptimal but computationally beneficial estimate. The above estimates can also be used for calculating the angle between the gradient vector ${\left[{\nabla }_{x}p{\nabla }_{y}p\right]}^T$ calculated at point p and the x-axis:

${\mathrm{\theta }}_p\cong {\tan }^{-1}\frac{{\nabla }_{y}p}{{\nabla }_{x}p}$ (3)

To avoid expensive calculations for computing trigonometric function and divisions, a two-dimensional lookup table data structure is used. By observing the numerical coefficients in Eqn (2) and considering integer values for pixel intensities in the image, the gradient components ∇_x and ∇_y would have pure integer values which make them suitable for indexing the look-up table. The entries of the look-up table are pre-computed and rounded off to integer degree values. This eliminates the need for floating point calculations and entails small rounding errors. The process of validating candidate boundary points consists of 4 steps. The steps are checked for each candidate boundary point in the order and a given candidate point is considered as rejected if it fails in any one of the steps.

Step 1: Detecting unidirectional neighbors

Figure 4 illustrates a schematic view of the validation steps. At the current candidate boundary point pⁱ, the gradient vector $\overrightarrow{{\nabla }_{p^i}}={\left[{\nabla }_x{p}^i{\nabla }_y{p}^i\right]}^T$ is calculated using the gradient masks in Eqn (2). Then, from the point pⁱ, two scan lines along two direction vectors $\overrightarrow{S_L}$ and $\overrightarrow{S_R}$ , both parallel with vessel border direction and perpendicular to vector $\nabla p^i$ are defined. Starting from point pⁱ, the gradient vector is calculated for each point along the scan lines within a certain distance. In this step, the current candidate point is considered valid if for all points m along the both scan lines S_L and S_R:

$\max \left\{\left|{\mathrm{\theta }}_{n^i}-{\mathrm{\theta }}_m\right|\right\}<\mathrm{\phi }m\text{=1,2,}\mathrm{...}\text{,}d$ (4)

where φ is needed to account for quantization effects and possible roughness along the vessel borders. In fact, the parameter φ provides an angular tolerance for the gradient directions between the neighboring points. The parameter is the distance which defines the number of neighboring points that must be verified.

Step 2: Detecting the opposite point

Another scan line is required to explore the corresponding opposite point at the other boundary of the vessel. This scan line is defined as a set of points along the vector $\overrightarrow{-{\nabla }_{p^i}}$ and within the maximum expected vessel width M. Let S_O denote the opposite scan line and m be a point on S_O, the point qⁱ considered as the corresponding opposite point for pⁱ if both of the following conditions are satisfied:

$q^i=\arg \max \left\{\left|{\nabla }_{x}m\right|+\left|{\nabla }_{y}m\right|\right\}{}_{m=1,2,\mathrm{...},M}$ (5)

and

$\left|\mathrm{\pi }-\left|{\mathrm{\theta }}_{p_i}-{\mathrm{\theta }}_{q_i}\right|\right|<\mathrm{\phi }$ (6)

where notation $\left|{\nabla }_{x}m\right|+\left|{\nabla }_{y}m\right|$ represents the estimation of local contrast at point m and is calculated as sum of absolute response values of the gradient masks in (2). The condition (6) accounts for verifying the similar gradient directions of the opposite points pⁱ and qⁱ. As a result, the proposed method largely addresses the problem of detecting true boundary points on the vessels with complex curvatures where there is no guarantee that the vessel borders are necessarily parallel. With this formulation the algorithm returns the first matched point as the corresponding opposite point qⁱ. Otherwise, this step would fail and another candidate point is verified.

Step 3: Detecting unidirectional neighbors of the opposite point

To ensure the proper selection of the opposite point, the unidirectional feature of the gradient vectors calculated at point qⁱ and its neighbors along the both scan lines S_L and S_R must be verified. In fact, the verification process in Step 1 is repeated for the opposite point detected in Step 2.

Step 4: Intensity and contrast thresholding

Although the previous steps can effectively verify the linearity of the features on which the candidate points are located, some detected points could still be due to noise or changes of illumination in the background. Therefore the reliable boundary points are selected based on a local contrast threshold.

The threshold is calculated based on the intensity values of the local minima or local maxima points on the grid lines^11,12,18. However, since the vessel structures only cover a small part of the image, most of the local maxima would be located on the background and their intensity values are contributed to calculate the threshold valuesand can influence the optimal value of the threshold¹². In this step, in order to reduce the effect of false points on the threshold value, the calculation of the threshold value is limited to those points validated in step 3. Therefore, the threshold value is calculated adaptively based on contrast values of the points which are almost validated.

Based on the above explanation, a seed point is considered as a valid seed point if the following condition is satisfied:

$C_p>=T_c$ (7)

Where C_p denotes the contrast value at point p and T_c is the contrast threshold calculated as follows:

$T_c={\mathrm{\mu }}_c-\left(\mathrm{\tau }{\mathrm{.\sigma }}_{{\mathrm{\mu }}_c}\right)$ (8)

in which μ_c is the median contrast value of the candidate points validated in step 3 and ${\mathrm{\sigma }}_{{\mathrm{\mu }}_c}$ is the standard deviation calculated base on μ_c. Parameter τ is considered to control the sensitivity of the validation algorithm to noise and contrast variations. Assuming normal distribution for the contrast values at candidate boundary points, the value of this parameter lies between 0 and 3. However, the optimal value for parameter τ is determined based on an experiment presented in section 4.2.

3.3 Estimating the Center Line Seed Points

After validating the boundary points, for each validated boundary point the location of its corresponding center line seed point is estimated by:

$\overrightarrow{c^i}=\left[\begin{array}{c}{c}^i{}_x\\ c^i{}_y\end{array}\right]=\frac{1}{2}\left[\begin{array}{c}{p}^i{}_x+q^i{}_x\\ p^i{}_y+q^i{}_y\end{array}\right]$ (9)

where cⁱ is the center line seed point and qⁱ is the corresponding opposite point for boundary point pⁱ. Since the boundary point pⁱ is already validated, there is no need to validate the estimated center line seed point.

The initial direction for the tracing algorithm can be calculated as follows:

$\begin{array}{l}\overrightarrow{u^i}=\arg \max \left\{\left|{\nabla }_{x}X\right|+\left|{\nabla }_{y}X\right|\right\}\\ X\in \{p_i,q_i\}\end{array}$ (10)

The tracing algorithm is initiated from center line seed points, once in direction $\overrightarrow{u_i}$ and once along $\overrightarrow{-u_i}$ . Other steps of the tracing algorithm are described in exploratory tracing methods^11,12,18.

A total of 10 angiograms images were selected from 80 angiogram sequences which were acquired by a GE-Innova 2100^IQ C-arm system. The selected images have spatial resolution 512 × 512 and 8-bit quantisation.

In this section, we define the performance measures necessary to discuss how the parameter values are determined and how they affect the detection results. In order to show the effectiveness of our method, we have compared the performance of our algorithm with the method proposed by Can¹², et al.

4.1 Performance Measures

Let TP denote the true positives, i.e., the total number of detected seed points that are located on or near the true vessel center lines. Let FP be the number of false positives, i.e., the number of detected seed points that are not located on or near any true vessel center lines. Let TN be the number of true negatives, i.e., the number of false seed points that are correctly discarded by the validation algorithm. Let FN be the number of false negatives, i.e., the total number of true seed points that are wrongly discarded by the validation algorithm. To measure the performance the standard terms of precision and recall are defined as follows:

• Precision: The percentage of the detected seed points that are on or near center lines in reality.

$P=\frac{TP}{TP+FP}$ (11)

• Recall: The percentage of actual center line seed points that are detected by the algorithm.

$R=\frac{TP}{TP+FN}$ (12)

Both of the above measures are important for performance evaluation. Nevertheless, to obtain a measure of balance between precision and recall, we use F-value measure to combine them into a single measure:

$F=\frac{2\times R\times P}{R+P}$ (13)

4.2 Setting the Parameter Values

The parameters for validation of the boundary points are angular tolerance φ, neighborhood distance d and contrast threshold scaling factor τ. These parameters are used to enhance the capabilities of the validation algorithm to control the strictness and the preciseness of the validation steps. As shown in Fig. 5, by observing on several images, satisfactory performance results are obtained when φ is chosen between 30 and 50. In this practice, we used φ = 45 for all images and obtained nearly perfect seed point detection performance.

The parameter d controls the strictness of the validation algorithm and can be either 1 or 2 or more. It can be seen from Figs. 6(a)-6(c) that for different values of N, the values of parameter d which are higher than 3, significantly increase the number of false negatives [Fig. 6(a)], but slightly affect the precision of the detection algorithm.

Another parameter is the contrast sensitivity control parameter τ. Low value of this parameter decreases the contrast threshold T_c and makes the seed point validation algorithm less sensitive to the local contrast of the collected candidate seed points. The effect of different values of this parameter is illustrated in Fig. 6(d) in which the average performance of the validation process for two categories of low contrast and high contrast images is measured. The results clearly show that the complete performance for both categories is obtained when τ >= 1.

4.3 Comparison Results

For calculating the performance measures, a set of 10 ground truth center line images was established based on the method described¹⁰. Each image was manually traced five times. The tracing was done by the same person but at different times. Then, for a given image, the correspondences between all of its manual traces were established by calculating the average of each correspondence set. Finally, the discontinuities in the estimated ground truth were eliminated by performing a morphological closing and true location of centerlines are obtained by applying a modified version of sequential thinning algorithm proposed¹⁰.

The closeness of each validated seed point to the nearest center line point is measured with respect to the corresponding ground truth center line within a particular distance δ (i.e. disk radius). In this study, equal values were chosen for common parameters of ‘grid searching’ part of both algorithms. The parameters values are as follows: N = 10, M = 26 and Ns = 26. This parameter selection results in equal number of candidate points for both validation algorithm. The following parameters are set for the proposed algorithm: angular tolerance φ = 23, neighborhood verification distance d = 3, contrast sensitivity control parameter τ = 0.5.

Table 1 lists the average values of precision, recall and F-value for 10 test images. The results show that the method proposed by Can¹², et al. has slightly higher precision. This is due to its strictness in seed point validation and small number of false positives. However, it misses a lot of reliable seed points which yields smaller recall values as the disk radius increases. Using the similarity of gradient vector angles as a validation criterion instead of the magnitude of the kernel responses (utilised by kernel-based methods) provides more complete detection where the contrast variations exist around the vessel boundaries. This advantage is due to the fact that variations of the contrast around a linear feature have less impact on the gradient angle than the gradient magnitude.

In this experiment, although the optimal parameter values are not set for our algorithm, higher amounts of F-values are achieved. This indicate that our method outperforms the other method in balancing between completeness (i.e., recall values) and preciseness (i.e., precision values) of the seed point detection process.

The conformance between the graphs in Figs. 2(a) and 2(b) implies that that the complexity of the proposed algorithm is on the order of K.N where K is a constant and N is the number of collected data points. Figure 7 illustrates the plot of running time of the seed point detection method proposed by Can¹², et al. and our method on different number of candidate points. The results show that our method is about 3 times faster than its predecessor. The figure also reveals that by increasing the number of candidate seed points, the running time of our method grows more gradually than kernel-based approach.

In this paper, a new seed point detection algorithm for tracing the centerlines of the vessels in digital angiograms is presented. The performance results show that the proposed method is capable of balancing between the running time, preciseness and completeness of the seed point detection process.

The comparison results on computation time show that the method is very fast compared with kernel-based methods. Its efficiency comes from 1) directly detecting valid boundary points without requiring expensive directional searches using 2-D differentiator kernels 2) avoiding large angular quantization errors associated with kernel-based methods. This method can be a suitable substitute for seed point detection in most of the exploratory tracing methods. For future improvements, we plan to enhance the precision of our algorithm by exploiting variable length gradient masks and perform comparative experiments with more known methods available in the literature.

The authors wish to thank the medical staff of Cardiology Department of Hospital HUKM Malaysia for their support in the data collection effort.

1. Fleagle, S.R. ; Johnson, M.R.; Wilbricht, C.J.; Skorton, D.J.; Wilson, R.F.; White, C.W.; Marcus, M.L. & Collins, S.M. Automated analysis of coronary arterial morphology in cineangiograms: Geometric and physiologic validation in humans. IEEE Trans. Med. Imaging, 2002, 8(4), 387-400.

2. Klein, A.K.; Lee, F. & Amini, A.A. Quantitative coronary angiography with deformable spline models. IEEE Trans. Med. Imaging, 2002, 16(5), 468-82.

3. Reiber, J.H.C.; Kooijman, C.J.; Slager, C.J.; Gerbrands, J.J.; Schuurbiers, J.C.H.; Den Boer, A.; Wijns, W.; Serruys, P.W. & Hugenholtz, P.G. Coronary artery dimensions from ineangiograms-methodology and validation of a computer-assisted analysis procedure. IEEE Trans. Med. Imaging, 2007, 3(3), 131-41.

4. Friman, O.; Kühnel, C. & Peitgen, H. Coronary artery centerline extraction using multiple hypothesis tracking and minimal paths. Insight Journal, 2008.

5. Hoyos, M.H.; Orkisz, M.; Douek, P.C. & Magnin, I.E. Assessment of carotid artery stenoses in 3D contrast-enhanced magnetic resonance angiography, based on improved generation of the centerline. Machine Graphics Vision Inter. J., 2005, 14(4), 349-78.

6. Zhang, Y.; Chen, K. & Wong, S. 3D Interactive Centerline Extraction. The Insight Journal, 2008.

7. Liu, I. & Sun, Y. Recursive tracking of vascular networks in angiograms based on the detection-deletion scheme. IEEE Trans. Med. Imaging, 2002, 12(2), 334-41.

8. Lu, S. & Eiho, S. Automatic detection of the coronary arterial contours with sub-branches from an X-ray angiogram. Paper presented at Computers in Cardiology, 1993. pp. 575-78.

9. Sandor, T.; D’Adamo, A.; Hanlon, W.B. & Spears, J.R. High precision quantitative angiography. IEEE Trans. Med. Imaging, 2007, 6(3), 258-65.

10. Al-Kofahi, K.A.; Can, A.; Lasek, S.; Szarowski, D.H.; Dowell-Mesfin, N.; Shain, W.; Turner, J.N. & Roysam, B. Median-based robust algorithms for tracing neurons from noisy confocal microscope images. IEEE Trans. Inf. Technol. B., 2003, 7(4), 302-17.

11. Al-Kofahi, K.A.; Lasek, S.; Szarowski, D.H.; Pace, C.J.; Nagy, G.; Turner, J.N. & Roysam, B. Rapid automated three-dimensional tracing of neurons from confocal image stacks. IEEE Trans. Inf. Technol. B., 2002, 6(2), 171-87.

12. Can, A.; Shen, H.; Turner, J.N.; Tanenbaum, H.L. & Roysam, B. Rapid automated tracing and feature extraction from retinal fundus images using direct exploratory algorithms. IEEE Trans. Inf. Technol. B., 1999, 3(2), 125-38.

13. Collorec, R. & Coatrieux, J.L. Vectorial tracking and directed contour finder for vascular network in digital subtraction angiography. Pattern Recogn. Lett., 1988, 8(5), 353-58.

14. Janssen, R.D.T. & Vossepoel, A.M. Adaptive Vectorization of Line Drawing Images. Comput. Vis. Image. Und., 1997, 65(1), 38-56.

15. Poli, R. & Valli, G. An algorithm for real-time vessel enhancement and detection. Comput. Meth. Prog. Bio., 1997, 52(1), 1-22.

16. Shen, H. ; Roysam, B. ; Stewart, C. V. ; Turner, J. N. & Tanenbaum, H. L. Optimal scheduling of tracing computations for real-time vascular landmark extraction from retinal fundus images. IEEE Trans. Inf. Technol. B., 2001, 5(1), 77.

17. Sun, Y. Automated identification of vessel contours in coronaryarteriograms by an adaptive tracking algorithm. IEEE Trans. Med. Imaging, 1989, 8(1), 78-88.

18. Zhang, Y.; Zhou, X.; Degterev, A.; Lipinski, M.; Adjeroh, D.; Yuan, J. & Wong, S.T.C. A novel tracing algorithm for high throughput imaging Screening of neuron-based assays. J. Neurosci. Meth., 2007, 160(1), 149-62.