Method and apparatus for three-dimensional reconstruction of coronary vessels from angiographic images and analytical techniques applied thereto

Abstract

A method for in-room computer reconstruction of a three-dimensional (3-D) coronary vascular tree from routine biplane angiograms acquired at arbitrary angles and without using calibration objects. The method includes seven major steps: (1) acquisition of two standard angiogram sequences by use of a single-plane or biplane imaging system, (2) identification of 2-D coronary vascular trees and feature extractions including bifurcation points, directional vectors, vessel centerlines, and construction of vessel hierarchy in the two images, (3) determination of transformation in terms of a rotation matrix R and a translation vector {right arrow over (t)} based on the extracted vessel bifurcation points and directional vectors, (4) establishment of vessel centerline correspondence, (5) calculation of the skeleton of 3-D coronary vascular trees, (6) rendering of 3-D coronary vascular tree with associated gantry angulation, and (7) calculation of optimal view(s) and 3-D QCA as quantitative measures associated with the selected vascular segment(s) of interest. The calculation of optimal views includes determination of 2-D projections of the reconstructed 3-D vascular tree so as to minimize foreshortening of a selected vascular segment, overlap of a selected vascular segment or both overlap and foreshortening of a selected arterial segment.

Description

BACKGROUND OF THE INVENTION

A portion of the disclosure of this patent document contains material which is subject to copyright protection. The copyright owner has no objection to the facsimile reproduction by anyone of the patent document or the patent disclosure, as it appears in the Patent and Trademark Office patent file or records, but otherwise reserves all copyright rights whatsoever.

1. Field of the Invention

The present invention relates generally to a method for reconstructing images of vascular structures and more specifically to an improved method for three-dimensional (3-D) reconstruction of vascular structures from two two-dimensional biplane projection images and methods and structures for quantitative analysis of such a reconstructed structure.

2. Discussion of Related Art

Several investigators have reported computer assisted methods for estimation of the 3-D coronary arteries from biplane projection data. These known methods are based on the known or standard X-ray geometry of the projections, placement of landmarks, known vessel shape, and on iterative identification of matching structures in two or more views. Such methods are described in a publication entitled “3-D digital subtraction angiography”, IEEE Trans. Med. Imag., vol. MI-1, pp. 152-158, 1982 by H. C. Kim, B. G. Min, T. S. Lee, et. al. and in a publication entitled “Methods for evaluating cardiac wall motion in 3-D using bifurcation points of the coronary arterial tree”, Invest. Radiology, January-February pp. 47-56, 1983 by M. J. Potel, J. M. Rubin, and S. A. Mackay, et al. Because the computation was designed for predefined views only, it is not suitable to solve the reconstruction problem on the basis of two projection images acquired at arbitrary and unknown relative orientations.

Another known method is based on motion and multiple views acquired in a single-plane imaging system. Such a method is described in a publication entitled “Reconstructing the 3-d medial axes of coronary arteries in single-view cineangiograms”, IEEE Trans. MI, vol. 13, no. 1, pp. 48-60, 1994 by T. V. Nguyen and J. Sklansky uses motion transformations of the heart model. However, the motion transformations of the heart model consist only of rotation and scaling. By incorporation of the center-referenced method, initial depth coordinates, and center coordinates, a 3-D skeleton of the coronary arteries was obtained. However, the real heart motion during the contraction involves five specific movements: translation, rotation, wringing, accordion-like motion, and movement toward the center of the ventricular chamber. Therefore, the model employed is not general enough to portray the true motion of the heart, especially toward the end-systole.

Knowledge-based or rule-based systems have been proposed for 3-D reconstruction of coronary arteries by use of a vascular network model. One such knowledge-based system is described in a publication entitled “An expert system for the labeling and 3-D reconstruction of the coronary arteries from two projections”, International Journal of Imaging, Vol. 5, No. 2-3, pp. 145-154, 1990 by Smets, Vandewerf, Suctens, and Oosterlinck. Because the rules or knowledge base were organized for certain specific conditions, it does not generalize the 3-D reconstruction process to arbitrary projection data. In other knowledge-based systems, the 3-D coronary arteries were reconstructed from a set of X-ray perspective projections by use of an algorithm from computed tomography. Due to the motion of the heart and only a limited number of projections (four or six), accurate reconstruction and quantitative measurement are not easily achieved.

Closed-form solutions of the 3-D reconstruction problem using a linear approach was a significant development and is described in, for example, a publication entitled “Determining 3-d motion and structure of a rigid body using the spherical projection”, CVGIP, vol. 21, pp. 21-32, 1983 by B. L. Yen and T. S. Huang. Unfortunately, actual data is always corrupted by noise or errors and the linear approach based techniques may not be sufficiently accurate when using noisy data. Hence, optimal estimation has been explicitly investigated. Additionally, U.S. Pat. No. 4,875,165 entitled Method for Determination of 3-D Structures in Biplane Angiography issued in the name of Fencil et al. also has significant drawbacks.

Use of a two-step method is known for producing an optimal estimation for a 3-D structure based on maximum-likelihood and minimum-variance estimation. In these techniques, for example, two publications entitled “Optimal motion and structure estimation”, IEEE Trans. on PAMI, Vol. 15, no. 9, September 1993, pp. 864-884, and “Structure from motion using the reconstruction and projection technique”, Proc. IEEE Workshop Computer Vision, November 1987, pp. 345-348, image error was employed in the objective function for a non-constricted minimization process. Preliminary estimates computed by a linear algorithm were used as initial estimates for the process of optimal estimation. However, if the initial solution from the linear approach is not sufficient, (e.g., with more than 2 pixels=0.6 mm error in the input 2-D image data), the iterative minimization process at the second step may become trapped in a local minimum due to a lack of prior information concerning the variables to be optimized.

Quantitative coronary analysis (“QCA”) of a reconstruction of an arterial tree was known and developed in the 1970's to quantify vessel geometry and the effects of drugs on the regression and progression of coronary artery disease. In the mid-1980's, digital systems were introduced into the catheterization laboratory to support the angiographer during the interventional procedures. With the advent digital angiographic technology, on-line QCA has been widely used predominantly for the selection of the optimal sizes of the interventional devices, and for the assessment of the efficacy of the procedure. However, current QCA techniques are performed on the 2-D projection views in which the vessel overlap and foreshortening are subjectively minimized in a “trial-and-error” manner by the interventionist. FIGS. 56 ( a )-( d ) demonstrate the problems resulting from use of such 2-D projections in QCA techniques. FIG. 56 c is an RCA image as known in the art used for QCA and resulting in a 30% narrowing measurement as shown in FIG. 56 d. With 2-D projection views there is no way to know or estimate how much error occurs in the QCA process due to foreshortening with respect to the stentotic segment. In FIGS. 56 a and 56 b, the same vessel segment between two bifurcation points as marked by two dots at its proximal and distal ends depicts 77% and 52% foreshortening, respectively.

It is evident from the above discussion that a need exists for improved reconstruction of 3-D images from 2-D image data and that a further need exists for improved QCA techniques utilizing such 3-D reconstruction to provide needed analysis in the intervention process.

SUMMARY OF THE INVENTION

The disadvantages of present methods known for visual reconstruction of vascular structures are substantially overcome with the present invention by providing a novel method for three-dimensional reconstruction of vessels using two-dimensional angiograms and further by providing improved QCA techniques utilizing such 3-D reconstructions of vascular structures.

This patent application teaches the methods and structures of the invention with reference to vascular structures (i.e., coronary arteries or coronary arterial trees). Those skilled in the art will readily recognize that the improved reconstruction methods and the quantitative analysis methods and structures are equally applicable to any vascular structure. References herein to “coronary” or “arterial” applications of the invention are intended merely as exemplary of common vascular structures where such reconstruction and analysis techniques are beneficially applied. Nothing in this patent application should be construed as limiting the reconstruction and analysis features to coronary applications. In particular, reference to “quantitative coronary analysis” (or “QCA”) should be understood to refer broadly to quantitative analysis of any vascular structures. Likewise, reference to an “artery” or “arteries” should be broadly understood to refer to any vessel or vessels. “Artery tree”, “arterial tree”, “coronary arterial tree” or “coronary artery tree” should all be understood to broadly refer to any vascular tree structure. Coronary applications of the present invention are therefore merely one common exemplary application of the features, methods and structures of the present invention.

In the present inventive method, a novel optimization technique minimizes the image point errors and directional vector errors in both imaging systems subject to the constraints derived from the individual intrinsic imaging parameters of the employed imaging system. Given five or more corresponding object points in both views (projection image), a constrained nonlinear optimization algorithm is applied to obtain an optimal estimate (transformation) of the biplane imaging geometry in the form of R and {right arrow over (t)} which characterize the position and orientation of one imaging system (imaging portion of the imaging system) relative to the other. The initial solution is estimated on the basis of the individual intrinsic imaging parameters.

Angiograms of eight hundred patients were analyzed in which two cases are selected for discussion hereinafter. The biplane imaging geometry was first determined without a calibration object, and the 3-D coronary arterial trees were reconstructed, including both left and right coronary artery systems. Various two-dimensional (2-D) projection images of the reconstructed 3-D coronary arterial tree were generated and compared to other viewing angles obtained in the actual patient study. Similarity between the real and reconstructed arterial structures was excellent. The accuracy of this method was evaluated by using a computer-simulated coronary arterial tree. Root-mean-square (RMS) errors in the 3-D position and the 3-D configuration of vessel centerlines and in the angles defining the R matrix and {right arrow over (t)} vector were 0.9-5.5 mm, 0.7-1.0 mm, and less than 1.5 and 2.0 degrees, respectively, when using 2-D vessel centerlines with RMS normally distributed errors varying from 0.4-4.2 pixels (0.25-1.26 mm).

More specifically, the method for three-dimensional reconstruction of a target object from two-dimensional images involves a target object having a plurality of branch-like vessels. The method includes the steps of: a) placing the target object in a position to be scanned by an imaging system, the imaging system having a plurality of imaging portions; b) acquiring a plurality of projection images of the target object, each imaging portion providing a projection image of the target object, each imaging portion disposed at an unknown orientation relative to the other imaging portions; c) identifying two-dimensional vessel centerlines for a predetermined number of the vessels in each of the projection images; d) creating a vessel hierarchy data structure for each projection image, the data structure including the identified two-dimensional vessel centerlines defined by a plurality of data points in the vessel hierarchy data structure; e) calculating a predetermined number of bifurcation points for each projection image by traversing the corresponding vessel hierarchy data structure, the bifurcation points defined by intersections of the two-dimensional vessel centerlines; f) determining a transformation in the form of a rotation matrix and a translation vector utilizing the bifurcation points corresponding to each of the projections images, the rotation matrix, and the translation vector representing imaging parameters corresponding to the relative orientations of the imaging portions of the imaging system; g) utilizing the data points and the transformation to establish a correspondence between the two-dimensional vessel centerlines corresponding to each of the projection images such that each data point corresponding to one projection image is linked to a data point corresponding to the other projection images, the linked data points representing an identical location in the vessel of the target object after the projections; h) calculating three-dimensional vessel centerlines utilizing the two-dimensional vessel centerlines and the correspondence between the data points of the two-dimensional vessel centerlines; and i) reconstructing a three-dimensional visual representation of the target object based on the three-dimensional vessel centerlines and diameters of each vessel estimated along the three-dimensional centerline of each vessel; and j) determining the optimal view of the vessel segments with minimal vessel foreshortening.

Further features of the present invention utilize such 3-D reconstructions to provide improved QCA techniques for analysis of coronary artery images.

BRIEF DESCRIPTION OF THE DRAWINGS

The patent or application file contains at least one drawing executed in color. Copies of this patent or patent application publication with color drawing(s) will be provided by the Office upon request and payment of the necessary fee.

FIG. 1 is a flowchart of a high level view of the reconstruction and analysis procedures of the present invention.

FIG. 2 is a flowchart providing additional details of the collection of gantry system parameters and image acquisition.

FIG. 3 is a diagram of a typical bi-plane imaging system.

FIG. 4 is a diagram of a typical single plane imaging system.

FIG. 5 is a diagram of intrinsic imaging parameters of a typical single plane imaging system with 7″ by 7″ field-of-view and 94 cm SID characterized by mathematical model.

FIG. 6 depicts a mathematical model of two projections, showing a 3-D object point and its projection in two views.

FIG. 7 is exemplary of estimated pixel size on 7″ field-of-view.

FIG. 8 is exemplary of the determination of the image centroid by calculation of the intersection of two lines.

FIG. 9 is a flowchart providing additional details of the identification of 2-D coronary arterial tree and feature extraction.

FIG. 10 depicts a manual indication of initial vessel centerline by selection of a sequence of points.

FIG. 11 depicts a manual indication of pseudo-branches by selection of a sequence of points originating from a desired landmark.

FIG. 12 depicts resultant ridge-point map where each dot corresponds to a local minima grey level.

FIG. 13 depicts final results of deformation process after an initial curve is aligned with the actual vessel centerlines.

FIG. 14 depicts a bifurcation point calculated by finding the intersection between a side branch and the precedent vessel.

FIGS. 15 a-c depict a vessel hierarchical digraph constructed from 2-D coronary tree identified in two images.

FIG. 16 is a flowchart providing additional details of the determination of the transformation defining spatial relationship of the employed two views.

FIG. 17 depicts the actual and initially calculated bi-plane geometry.

FIG. 18 depicts the calculation of bi-plane transformation based on bifurcation points and directional vectors.

FIG. 19 depicts cone-shaped bounding regions used to constrain the calculated 3-D points.

FIG. 20 is a flowchart providing additional details of the establishment of vessel centerline correspondence for each pair of arteries in the coronary trees.

FIG. 21 depicts the selection of corresponding pairs of arteries using vessel hierarchical digraphs.

FIG. 22 depicts epipolar plane geometry constraints.

FIG. 23 depicts epipolar line geometry corresponding to FIG. 22 .

FIGS. 24 a-b depict typical examples of finding corresponding points using epipolar lines.

FIG. 25 depicts a refinement process for establishment of vessel centerline correspondences.

FIGS. 26 a-d depict comparison of previous methods for establishing correspondence with the methods of the present invention.

FIG. 27 is a flowchart providing additional details of calculation of 3-D coronary arterial skeleton in terms of 3-D vessel centerlines and cross-sectional diameters.

FIG. 28 is a flowchart providing additional details of rendering the reconstructed 3-D coronary arterial tree with gantry angulation.

FIG. 29 is an exemplary skeleton image of 3-D coronary arteries.

FIG. 30 is an exemplary surface based modeling image of coronary arteries.

FIGS. 31 a-f depict 3-D reconstruction of RCA images viewing at four different angles of sample acquired images.

FIGS. 32 a-f depict 3-D reconstruction of LCA images viewing at four different angles of sample acquired images.

FIG. 33 is a typical example of the vessel skeleton structure in terms of 3-D curves.

FIGS. 34 a-d depict schematic diagrams of curvature and torsion.

FIGS. 35 a-b depict schematic diagrams of 3-D and 2-D take-off angles.

FIG. 36 depicts a schematic diagram of insertion angle of a bypass graft.

FIG. 37 is a flowchart providing additional details of QCA processing.

FIGS. 38 a-b depict a selected region of interest on a bypass graft and its QCA result.

FIGS. 39 a-b depict a selected region of interest on a bypass graft and its regional QCA result.

FIGS. 40 a-c depict 3-D and 2-D take-off angles of LCA and an example of insertion angle of a bypass graft.

FIGS. 41 a-c depict the percent foreshortening.

FIGS. 42 a-b depict typical vessel foreshortening on LCA segments.

FIGS. 43 a-d depict accumulated percent foreshortening on the entire LCA reconstructed from an image pair.

FIGS. 44 a-d depict a selected single segment of a region of interest and its maps.

FIGS. 45 a-d depict a selected bifurcation region of a region of interest and its maps.

FIGS. 46 a-d depict a selected bypass graft segment of a region of interest and its three regional maps.

FIG. 47 is a flowchart providing additional detail of the processing to assess length of a selected arterial segment.

FIG. 48 is a flowchart providing additional detail of the processing to assess tortuosity of a selected arterial segment.

FIG. 49 is a flowchart providing additional detail of the processing to assess take-off angles of a selected arterial segment.

FIG. 50 is a flowchart providing additional detail of the processing to assess insertion angles of a selected arterial segment.

FIG. 51 is a flowchart providing additional detail of the processing to verify connectivity of a selected arterial segment.

FIG. 52 is a flowchart providing additional detail of the processing to optimal view of a selected arterial segment to minimize foreshortening.

FIG. 53 is a flowchart providing additional detail of the processing to determine accumulated percent foreshortening of a selected arterial segment.

FIG. 54 is a flowchart providing additional detail of the processing to optimal view of a selected arterial segment to minimize overlap.

FIG. 55 is a flowchart providing additional detail of the processing to optimal view of a selected arterial segment to minimize both foreshortening and overlap.

FIGS. 56 a-d depict exemplary vessel foreshortening and existing 2-D QCA process.

DETAILED DESCRIPTION OF PREFERRED EMBODIMENTS

While the invention is susceptible to various modifications and alternative forms, a specific embodiment thereof has been shown by way of example in the drawings and will herein be described in detail. It should be understood, however, that it is not intended to limit the invention to the particular form disclosed, but on the contrary, the invention is to cover all modifications, equivalents, and alternatives falling within the spirit and scope of the invention as defined by the appended claims.

Overview

A mono-plane or bi-plane angiographic system can be employed for cardiac angiographic acquisition. The pair of angiograms are acquired based on any two selected gantry orientation in terms of left/right anterior oblique (RAO/LAO) angle, caudal/cranial (CAUD/CRAN) angle. The angulation is defined with respect to the iso-center to which the gantry arm rotates. The gantry system information including focal spot to image intensifier distance (SID), field of view (FOV, e.g., 5, 7, or 9 inch mode), and gantry orientation as in FIG. 1 , and the acquired images are automatically recorded. The employed 3-D reconstruction protocol involves the following steps:

1. Routine cardiac catheterization is initiaeeed and a standard coronary artenogram is completed in two views; one injection in the biplane system, two injections in the single plane system.

2. A pair of image sequences, consisting of 8 images (including from end-diastole to end-systole frames) at each view, are chosen and transferred to a computer workstation.

3. The 2-D coronary arterial trees are identified in the selected pair of angiograms (e.g., images corresponding to end-diastole) and the vessel features are extracted including vessel centerlines, bifurcation points, and hierarchy.

4. The processing of 3-D reconstruction is completed by use of the information of identified vessel features yielding a 3-D patient-specific coronary arterial tree model.

5. The computer-simulated projection from reconstructed 3-D coronary model along with the corresponding gantry angles are displayed on a computer graphics workstation.

6. The process of 3-D quantitative coronary analysis (QCA) process including optimal-view strategy map (OSM) and assessment of coronary morphology relative to any selected arterial segment(s) is completed. The results are used as guidance for minimizing the vessel foreshortening and overlap, and for facilitating selection of other intra-coronary devices, i.e., stents.

The proposed method of reconstruction consists of seven major steps shown in the flowchart figures and descriptive text to follow. The major steps are:

(1) acquisition of two standard angiogram sequences by use of a single-plane or biplane imaging system,

(2) identification of 2-D coronary arterial trees and feature extractions including bifurcation points, vessel centerlines, and construction of vessel hierarchy in the two images,

(3) determination of transformation in terms of a rotation matrix R and a translation vector {right arrow over (t)} based on the extracted vessel bifurcation points,

(4) establishment of vessel centerline correspondence,

(5) calculation of the skeleton of 3-D coronary arterial trees,

(6) rendering of 3-D coronary arterial tree with associated gantry angulation, and

(7) calculation of optimal view(s) and 3-D QCA associated with the selected arterial segment(s) of interest.

A. Image Acquisition

Details of the methods for acquiring images of the arterial tree are shown in the flowchart of FIG. 2 .

A pair of angiographic views are acquired by a typical single-plane or biplane imaging system as shown in FIG. 3 and FIG. 4 with intrinsic imaging parameters as illustrated in FIG. 5 . The spatial relationship of the acquired two views can be characterized by use of a mathematical model as shown in FIG. 6 . The imaging system parameters including pixel size and the imaging center coordinates are calculated as shown in FIG. 7 and FIG. 8 , respectively. Once the information are calculated, they will be utilized in al the subsequent 3-D reconstruction process. After routine cardiac catheterization is initiated, a standard coronary arteriogram is completed in two standard views based on one injections in a bi-plane imaging system. With the acquired two angiographic sequences, a pair of images were chosen corresponding to the same or close acquisition time instance (e.g., end-diastole) in the respective cardiac cycles. They were transferred to the computer workstation through the fiber optic link for analysis and reconstruction. The transfer rate is 2 seconds per image with the size of 512×512 pixels and one byte storage per pixel.

B. Segmentation and Feature Extraction

Details of the methods for defining segments and extracting features from the arterial tree are shown in the flowchart of FIG. 19 .

A semi-automatic system based on the technique of deformation model and segmentation technique was employed for identification of the 2-D coronary arterial tree in the given angiograms. The required user interaction in the identification process involves only the indication of several points inside the lumen of each artery, near the projection of vessel centerline in the angiogram as shown in FIG. 10 . The artery with stenotic lesion or tortuous shape can be selected as landmarks to improve the reconstruction results as shown in FIG. 11 . Afterwards, a spline-curve is formed based on the selected points, which serves as the initial centerline of the vessel. Meanwhile, an operator was automatically applied to the image to identify those pixels which are local minimum on intensity (i.e., the intensity of coronary arterial structures is darker than that of background). These identified pixels indicate the location of vessel centerline as shown in FIG. 12 . By use of the deformation process, the identified pixels serve as the external forces to act on the initial model curve such that it will be gradually deformed and finally resides on the actual centerline of vessel as shown in FIG. 13 .

The identified centerlines and the branching relationships are used for construction of the vessel hierarchy in each angiogram by a labeling process according to the anatomy and connectivity among the primary and secondary coronary arteries as shown in FIG. 14 . The labeling process on the coronary tree is performed automatically by application of the algorithm of the breadth-first search to traverse identified vessel centerlines and associated bifurcations. From each vessel of the coronary tree that is currently visited, this approach tries to search as broadly as possible by next visiting all of the vessel centerlines that are adjacent to it, and finally results in a vessel hierarchically directed graph (digraph) containing a set of nodes with corresponding depth levels and two types of arcs (descendant and sibling arcs) defining the coronary anatomy. In the hierarchical digraph, each node represents a vessel in the coronary arterial tree, and the descendent arc connects the nodes between the vessel and its branch. If a vessel has multiple branches, the nodes associated with these branches are connected as a linked list by sibling arcs according to the branching locations relative to the proceeding vessel.

With the constructed hierarchical digraph, the remaining vessel features including bifurcation points and diameters as shown in FIGS. 15 a-c are then calculated for all of the coronary arterial segments, and the results are saved into the nodes associated with the branching vessels in the hierarchical digraph.

The deformation model is controlled under the influence or energy E def of the internal forces characterized by smoothness constraint and external forces derived from image features. $$ E def = ∑ i = 1 n ⁢ { ∑ j = 1 m 1 ⁢ E int ( l i ⁡ ( v ⁢   ⁢ ( s j ) ) + E ext ( l i ⁡ ( v ⁢   ⁢ ( s j ) ) }

where n denotes the number of coronary arteries in the coronary arterial tree, and m i , i=1, . . . ,n denotes the number of centerline points v(s j )=(x(s j ), y(s j )) on each artery l i (v),The internal curve energy can be written as

E int ( l i ( v ( s j ))=(α( s j )| v s ( s j )| 2 +β( s j )| v ss ( s j )| 2 )/2,

where α(s) and the first order term v s (s) make the model curve like a membrane and β(s) and the second order term v ss (s) make the it act like a thin plate. The external energy is derived from the ridge points of image I denoting the minimal (or maximal) gray level in a local window (e.g., k pixels by k pixels) region as

E ext ( l i ( v ( s j ))=∇ k *( G k ( v ( s j ))⊕ I ))

where ⊕ denotes the convolution operator, G k (v(s)) is a Gaussian based filter centered at v(s) with k by k pixels of window, ∇ k denotes the ridge point filter defined within a k by k pixels of window.

C. Calculation of Transformation Characterizing the Two Views

Details of the methods for transforming the views of the arterial tree are shown in the flowchart of FIG. 16 .

By use of a mono-plane or biplane imaging system for image acquisitions, the locations and dimensions of imaging plane and x-ray source can be defined by use of a mathematical model as shown in FIG. 16 . The spatial relationship between the employed two views can be characterized by a transformation in the forms of a rotation matrix R and a translation vector {right arrow over (t)} on the basis of the respective locations of x-ray source (or focal spot) served as the origin in 3-D coordinate space.

In the first view, let (u i , v i ) denote the image coordinates of the ith object point, located at position (x i , y i , z i ). We have u i =Dx i /z i , v i =Dy i /z i , where D is the perpendicular distance between the x-ray focal spot and the image plane when the view was acquired. Let (ξ i , η i ) denote the scaled image coordinates, defined as ξ i =u i /D=x i /z i , η i =v i /D=y i /z i . The second projection view of the employed bi-plane imaging system can be described in terms of a second pair of image and object coordinate systems u′v′ and x′y′z′ defined in an analogous manner. Scaled image coordinates (ξ i , η′ i ′) in the second view for the ith object point at position (x i ′, y i ′, z i ′) are given by ξ i ′=u i ′/D′=x i ′/z i ′, η i ′=v i ′/D′=y i ′/z i ′. The spatial relationship between the two views can be characterized by $$ [ x i ′ y i ′ z i ′ ] = R · { [ x i y i z i ] - t → } = ⌊ r 11 r 12 r 13 r 21 r 22 r 23 r 31 r 32 r 33 ⌋ · [ x i - t x y i - t y z i - t z ] .

In our previously proposed technique, the transformation was calculated based on the identified bifurcation points of the 2-D coronary arterial trees in the two views and can be only effective to the images acquired from the biplane systems or a single-plane imaging system without table panning or movement during the acquisition.

In this invention, a new 3-D reconstruction technique is developed by incorporation of the corresponding bifurcation points and the vessel directional vectors of bifurcations. The proposed method can also be employed to obtain accurate results of 3-D reconstruction without the need of ECG-triggered acquisition even when there is table panning or moving during the injection. The required prior information (i.e., the intrinsic parameters of single-plane imaging system) for 3-D reconstruction include: (1) the distance between each focal spot and its image plane, SID (focal-spot to imaging-plane distance), (2) the field of view in terms of the pixel size, p size (e.g., 0.3 mm/pixel), (3) the distance {overscore (ff′)} between the two focal spots or the known 3-D distance between two points in the projection images, and (4) iso-center distance (with respect to which the rotary motion of the gantry arm rotates) relative to the focal spot. Items (1) and (2) are prerequisite and (3) and (4) are optional.

In the single-plane imaging system utilized, the gantry angulation associated with different view is always recorded and accompany with the acquired image file. However, the transformation directly calculated from the recorded gantry angles may not accurately define the spatial relationship between the two views due to (i) images non-simultaneously acquired from two different cardiac cycles and (ii) table panning or movement during the acquisitions. Although such a directly calculated transformation can not truly characterize the two views, the resultant estimate is generally good enough to serve as an initial guess to the non-linear minimization process for obtaining an optimal solution. In the biplane imaging system utilized, the gantry angulations associated with different views may or may not be recorded with the acquired image file. If the gantry angulations of the employed two views are unknown, orthogonal angles will be used as the initial solutions of transformation as shown in FIG. 17 by use of the equations herein below.

If the gantry angulations are known in advance, however, the transformation directly calculated from the recorded gantry angles may not accurately define the spatial relationship between the two views due to (i) patient's respiration, (ii) non-synchronous x-ray source firing on the x-ray tubes, and (iii) mis-alignment of the respective iso-centers corresponding to the individual x-ray gantry instruments of the bi-plane system. Although such a directly calculated transformation can not truly characterize the two views, the resultant estimate is generally good enough to serve as an initial guess to the non-linear minimization process for obtaining an optimal solution. With the initial solution calculated from the pre-defined orthogonal angels or recorded gantry angles, the global minimum resulting from the non-linear optimization can be guaranteed since the initial estimate is close to the actual solution.

On the basis of the initial solution, the final accurate transformation is calculated from an objective function by minimization of (1) image point errors in term of the square of Euclidean distance between the 2-D input data and the projection of calculated 3-D data points and (2) directional vector errors between the 2-D vessel directional vectors and the projection of calculated 3-D vessel directional vectors based on equations below and as shown in FIG. 18 subject the constraints shown below and in FIG. 19 .

The initial estimates of the z and z′ axes of the two imaging systems are taken to be orthogonal in this situation, and the unit translation vector are assumed to lie on the x-z plane. Two different initial solutions ({overscore (R 1 )}, {right arrow over (t)} {overscore (u1)} ) and ({overscore (R 2 )}, {right arrow over (t)} {overscore (u2)} ) are employed as follows: $$ R _ 1 = [ 0 0 - 1 0 1 0 1 0 0 ] ,   ⁢ t → u _ 1 = ⌊ - D ′ α ′ · t d 0 D α · t d ⌋ ,   ⁢ and $$ R _ 2 = [ 0 0 1 0 1 0 - 1 0 0 ] ,   ⁢ t → u _ 2 = [ D ′ α ′ · t d 0 D α · t d ] ,

where D and D′ denote the SID distances, α and α′ are the MF factors, and t d represents the magnitude of {right arrow over (t)}. If the magnitude of {right arrow over (t)} is not available, an approximate measure is taken as: $$ t d = ( D α ) 2 + ( D ′ α ′ ) 2 .

Let α 1 ,β 1 ) and (α 2 ,β 2 ) denote the gantry orientations of the 1st and 2nd views, where α x and β x are the LAO/RAO, and CRAN/CAUD angles, respectively. Let (x 1 , y 1 , z 1 ) and (x 2 , Y 2 , z 2 ) denote the respective position vectors defined in the 1st and 2nd views. The spatial relationship between the two views can be characterized based on the recorded gantry information as $$ [ x 2 y 2 z 2 ] = R · { [ x 1 y 1 z 1 ] - t → } = ⌊ r 11 r 12 r 13 r 21 r 22 r 23 r 31 r 32 r 33 ⌋ · [ x 1 - t x y 1 - t y z 1 - t z ] .

where

  R=R x (−α 1 )· R y (β 1 )· R y (−β 2 )·R x (α 1 ),

and $$ ⌊ t x t y t z ⌋ = ( R x ⁡ ( - α 1 ) · R y ⁡ ( β 1 ) · R y ⁡ ( - β 2 ) · R x ⁡ ( α 2 ) · [ 0 0 - Iso 2 ] ⁢   ) - [ 0 0 - Iso 1 ] .

The matrices R x (θ) and R y (φ) denote the rotation with respect to x-axis and y-axis with θ and φ angles, respectively, and Iso 1 and Iso 2 denote the two iso-center distances associated with the 1st and 2nd views.

The scaled 3-D points ({overscore (x)}′ i , {overscore (y)}′ i , {overscore (z)}′ i ) defined in the x′y′z′ coordinate system are initialized as $$ x _ i ′ = u i ′ α ′ · t d ,   ⁢ y _ i ′ = v i ′ α ′ · t d ,   ⁢ z _ i ′ = D ′ α ′ · t d ,   ⁢ i = 1 , … ⁢   , n ,

where (u′ i , v′ i ) denotes the 2-D input coordinates of the input points. In other words, the initial scaled 3D points are taken to lie in a plane that is parallel to the image plane in the x′y′z′ system and passes through the iso-center.

Let {right arrow over (μ)} i ={μ i x , μ i y } and {right arrow over (μ)} i ′={μ′ i x , μ′ i y } denote the respective identified 2-D vessel directional unit vectors of bifurcations in each views. Let v={{right arrow over (v)} 1 , {right arrow over (v)} 2 , . . . , and {right arrow over (v)} n } and v′={{right arrow over (v)} 1 ′, {right arrow over (v)} 2 ′, . . . , {right arrow over (v)} n ′} denote the calculated 3-D vessel directional unit vectors of bifurcations in two views. The initial estimate of the 3-D directional vectors defined in the two views are characterized as

{right arrow over (μ)} i =[μ i x , μ i y , 0] t , {right arrow over (μ)} i ′=[μ′ i x ,μ′ i y ,0] t

With the known estimates of the scaled image coordinates (ξ i , η i ) and (ξ i , η i ′), transformation in terms of R and {right arrow over (t)}, the object coordinates (x i ,y i ,z i ) can be calculated as follows: $$ ⌊ r 11 - r 31 ⁢ ξ i ′ r 12 - r 32 ⁢ ξ i ′ r 13 - r 33 ⁢ ξ i ′ r 21 - r 31 ⁢ η i ′ r 22 - r 32 ⁢ η i ′ r 23 - r 33 ⁢ η i ′ 1 0 - ξ i 0 1 - η i ⌋ ⁢   [ x i y i z i ] = [ a → · t → b → · t → 0 0 ] ,

where $$ a → = ⌊ r 11 - r 31 ⁢ ξ i ′ r 12 - r 32 ⁢ ξ i ′ r 13 - r 33 ⁢ ξ i ′ ⌋ ,   ⁢ b → = ⌊ r 21 - r 31 ⁢ η i ′ r 22 - r 32 ⁢ η i ′ r 23 - r 33 ⁢ η i ′ ⌋ .

The scaled coordinates of the initial 3D points ({right arrow over (x)} i , {right arrow over (y)} i , {right arrow over (z)} i ) defined in the x′y′z! coordinate system are calculated by ®

([ x i ,y i ,z i t] t −{right arrow over (t)})/| {right arrow over (t)}|.

Let {right arrow over (μ)} i ={μ i x , μ i y } and {right arrow over (μ)}′ i ={μ′ i x , μ′ i y } denote the respective identified 2-D vessel directional unit vectors of bifurcations in each view. Let v={{right arrow over (v)} 1 , {right arrow over (v)} 2 , . . . , {right arrow over (v)} n } and v′={{right arrow over (v)}′ 1 , {right arrow over (v)}′ 2 , . . . , {right arrow over (v)}′ n } denote the calculated 3-D vessel directional unit vectors of bifurcations in two views. The initial estimate of the 3-D directional vectors defined in the two views are characterized as

{right arrow over (μ)} i =[μ i x , μ i y ,0] t ,{right arrow over (μ)}′ i =[μ′ i x ,μ′ i y ,0] t

When there occurs a table panning or movement during the image acquisitions, a translation vector {right arrow over (r)}=[r x ,r y ,r z ] t that defines the initial solution of such a movement can be estimated if the transformation in terms of a rotation matrix and a translation vector is known as $$ min r → ⁢   ⁢ ℱ ⁢ ( r → ) = ∑ i = 1 n ⁢ ( ξ i ′ - x i ′ z i ′ ) 2 + ( η i ′ - y i ′ z i ′ ) 2 + ( ξ i - c → 1 · p → i ′ + t x + r x c → 3 · p → i ′ + t z + r z ) 2 + ( η i - c → 2 · p → i ′ + t y + r y c → 3 · p → i ′ + t z + r z ) 2 ,

where {right arrow over (c)} k denotes the respective kth column vectors of matrix R. The initial estimates of the 3-D bifurcation points {right arrow over (p)}′ i , i=1, . . . , n can be calculated based on the formula as shown in the above equations.

Given the set of 2-D bifurcation points and directional vectors extracted from the pair of images, “optimal” estimates of the transformation and 3-D bifurcation points are obtained by finding the minimum of the following equation: $$ min p , p ′ , v , v ′ ⁢ ℱ 1 ⁢ ( P , P ′ , v , v ′ ) = ∑ i = 1 n ⁢ ( ξ i - x i z i ) 2 + ( η i - y i z i ) 2 + ( ξ i ′ - x i ′ z i ′ ) 2 + ( η i ′ - y i ′ z i ′ ) 2 + &LeftDoubleBracketingBar; μ → i - [ v → i ] xy &RightDoubleBracketingBar; 2 + &LeftDoubleBracketingBar; μ → i ′ - [ v → i ′ ] x ′ ⁢ y ′ &RightDoubleBracketingBar; 2 ,

where n denotes the number of pairs of corresponding points extracted from the two images and P and P′ denote the sets of 3-D object position vectors (i.e., bifurcation points), {right arrow over (p)} i =(x i ,y i ,z i ) and {right arrow over (p)}′ i =(x′ i ,y′ i , z′ i ), respectively, {right arrow over (μ)} i and {right arrow over (μ)}′ i , denote the respective 2-D vessel directional vectors of bifurcations in each views, and v={{right arrow over (v)} 1 , {right arrow over (v)} 2 , . . . , {right arrow over (v)} n } and v′={{right arrow over (v)}′ 1 , {right arrow over (v)}′ 2 , . . . , {right arrow over (v)}′ n } denote the calculated 3-D vessel directional vectors of bifurcations in two views, and [{right arrow over (v)} i ] xy and [{right arrow over (v)}′ i ] x′y′ , i=1, . . . ,n, denote the calculated 2-D projection vectors from {right arrow over (v)} i and {right arrow over (v)} i onto the image planes, respectively.

Since the relationship between the two views can be characterized by a rotation matrix R and a translation vector {right arrow over (t)}=[t x , t y , t z ] t , the above equation can be expressed as $$ min R , t → , p ′ , v ′ ⁢ ℱ 2 ⁢ ( R , t → , p ′ , v ′ ) = ∑ i = 1 n ⁢ ( ξ i ′ - x i ′ z i ′ ) 2 + ( η i ′ - y i ′ z i ′ ) 2 + ( ξ i - c → 1 · p → i ′ + t x c → 3 · p → i ′ + t z ) 2 + ( η i - c → 2 · p → i ′ + t y c → 3 · p → i ′ + t z ) 2 + &LeftDoubleBracketingBar; μ → i - [ R - 1 · v → i ′ ] xy &RightDoubleBracketingBar; 2 + &LeftDoubleBracketingBar; μ → i ′ - [ v → i ′ ] x ′ ⁢ y ′ &RightDoubleBracketingBar; 2 ,

where {right arrow over (c)} k denotes the respective kth column vectors of matrix R and R 1 is the inverse matrix of R. When there occurs a table panning or movement during the image acquisitions, however, the above equation can not truly characterize the spatial relationship of the two views. Instead, a translation vector {right arrow over (r)}=[r x , r y , r z ] t needs to be incorporated to characterize such a situation, and the above equation can be rewritten as $$ min R , t → , p ′ , v ′ , r → ⁢ ℱ 3 ⁢ ( R , t → , p ′ , v ′ ⁢ r → ) = ∑ i = 1 n ⁢ ( ξ i ′ - x i ′ z i ′ ) 2 + ( η i ′ - y i ′ z i ′ ) 2 + ( ξ i - c → 1 · p → i ′ + t x + r x c → 3 · p → i ′ + t z + r z ) 2 + ( η i - c → 2 · p → i ′ + t y + r y c → 3 · p → i ′ + t z + r z ) 2 + &LeftDoubleBracketingBar; μ → i - [ R - 1 · v → i ′ ] xy &RightDoubleBracketingBar; 2 + &LeftDoubleBracketingBar; μ → i ′ - [ v → i ′ ] x ′ ⁢ y ′ &RightDoubleBracketingBar; 2 ,

It is well known that any 3D rigid motion can be decomposed uniquely into a translation and a rotation by an angle θ around a unit vector {right arrow over (v)} u passing through the origin of the coordinate system. In the proposed method, a quaternion denoted by {tilde over (q)}=(s,{right arrow over (w)})=(s, w 1 , w 2 , w 3 ) of norm equal to 1 is employed to represent the rotation transformation as

{right arrow over (w)} =sin(θ/2){right arrow over (v)} u , s =cos(θ/2).

For any quaternion {tilde over (q)}=(s,{right arrow over (w)})=(s, w 1 , w 2 , W 3 ) of norm 1, there exists a rotation matrix defined as: $$ R ≡ 2 ⁢ ⌊ s 2 + ( w 1 ) 2 - 1 2 w 1 ⁢ w 2 - sw 3 sw 2 + w 1 ⁢ w 3 w 1 ⁢ w 2 + sw 3 s 2 + ( w 2 ) 2 - 1 2 w 2 ⁢ w 3 - sw 1 w 1 ⁢ w 3 - sw 2 sw 1 + w 2 ⁢ w 3 s 2 + ( w 3 ) 2 - 1 2 ⌋

With this quaternion representation, the above equations can be rewritten as $$ min q ~ , t → , p ^ ′ ⁢ ℱ 4 ⁡ ( q ~ , t → u ⁢ p ^ ′ ) = ∑ i = 1 n ⁢ ( ξ i ′ - x ^ i ′ z ^ i ′ ) 2 + ( η i ′ - y ^ i ′ z ^ i ′ ) 2 + { ξ i - 2 ⁢ ( s 2 + w 1 2 - 1 / 2 ) ⁢ x ^ i ′ + 2 ⁢ ( w 1 ⁢ w 2 + sw 3 ) ⁢ y ^ i ′ + 2 ⁢ ( w 1 ⁢ w 3 - sw 2 ) ⁢ z ^ i ′ + t u x + r _ x 2 ⁢ ( sw 2 + w 1 ⁢ w 3 ) ⁢ x ^ i ′ + 2 ⁢ ( w 2 ⁢ w 3 - sw 1 ) ⁢ y ^ i ′ + 2 ⁢ ( s 2 + w 3 2 - 1 2 ) ⁢ z ^ i ′ + t u z + r _ z } 2 + { η i - 2 ⁢ ( w 1 ⁢ w 2 - sw 3 ) ⁢ x ^ i ′ + 2 ⁢ ( s 2 + w 2 2 - 1 / 2 ) ⁢ y ^ i ′ + 2 ⁢ ( sw 1 + w 2 ⁢ w 3 ) ⁢ z ^ i ′ + t u y + r _ y 2 ⁢ ( sw 2 + w 1 ⁢ w 3 ) ⁢ x ^ i ′ + 2 ⁢ ( w 2 ⁢ w 3 - sw 1 ) ⁢ y ^ i ′ + 2 ⁢ ( s 2 + w 3 2 - 1 2 ) ⁢ z ^ i ′ + t u z + r _ z } 2 + &LeftDoubleBracketingBar; μ → i - [ R - 1 · v → i ′ ] xy &RightDoubleBracketingBar; 2 + &LeftDoubleBracketingBar; μ → i ′ - [ v → i ′ ] x ′ ⁢ y ′ &RightDoubleBracketingBar; 2 ,

where ({circumflex over (x)}′ i , ŷ′ i , {circumflex over (z)}′ i ) denotes the scaled 3D coordinates, within a scale factor of the magnitude of the translation vector |{right arrow over (t)}| and where [tu x ,tu y ,tu z ] t denotes the unit translation vector that corresponds to {right arrow over (t)}, and {right arrow over (r)}/|{right arrow over (t)}|=[{right arrow over (r)} x , {right arrow over (r)} y , {right arrow over (r)} z ]. $$ C 1 ⁢ : s 2 + w 1 2 + w 2 2 + w 3 2 = 1 C 2 ⁢ : t u x 2 + t u y 2 + t u z 2 = 1 C 3 ⁢ : 0 < z ^ i ′ ,   ⁢ i = 1 , … ⁢   , n C 4 ⁢ : 0 < 2 ⁢ ( sw 2 + w 1 ⁢ w 3 ) ⁢ x ^ i ′ + 2 ⁢ ( w 2 ⁢ w 3 - sw 1 ) ⁢ y ^ 1 ′ + 2 ⁢ ( s 2 + w 3 2 - 1 2 ) ⁢ z ^ i ′ + t u z + r _ z ,     ⁢ i = 1 , … ⁢   , n .

The constraint C 1 characterizes the quaternion norm; the constraint C 2 ensures a unit translation vector; and the constraints C 3 and C 4 force the scaled coordinates {circumflex over (z)}′ i and {circumflex over (z)} i to be in front of the respective focal spots.

If the iso-center distances of the biplane imaging systems or the MF factors are available, the constraints C 3 and C 4 above can be modified as: $$ D ′ / α ′ - δ h &LeftBracketingBar; t → &RightBracketingBar; < z ^ i ′ < D ′ / α ′ + δ h &LeftBracketingBar; t → &RightBracketingBar; , … ⁢   , n $$ D / α - δ h &LeftBracketingBar; t → &RightBracketingBar; < z ^ i = 2 ⁡ [ ( sw 2 + w 1 ⁢ w 3 ) ⁢ x ^ i ′ + ( w 2 ⁢ w 3 - sw 1 ) ⁢ y ^ i ′ + ( s 2 + w 3 2 - 1 2 ) ⁢ z ^ i ′ ] + t u z + r _ z ⁢   ⁢ … ⁢   , n ,

where α and α′ denote the MF factors and δ h (≈12.5±2.0 cm at end-diastole) denotes the maximal length of the heart in a long-axis view. For each 3D object point, the ray connecting that point and the focal spot intersects the image plane near the associated 2D image point even when the input data are corrupted by noise. In addition to the constraints imposed on the z and z′ coordinates, two other constraints are incorporated to confine the x, x′, y, and y′ coordinates of each calculated 3D point as follows: $$ C 5 ⁢ : ⁢   ⁢ ( x ^ i z ^ i - ξ i ) 2 + ( y ^ i z ^ i - η i ) 2 ≤ ( δ c p size ) 2 ,   ⁢ i = 1 , … ⁢   , n , ⁢ C 6 ⁢ : ⁢   ⁢ ( x ^ i ′ z ^ i ′ - ξ i ′ ) 2 + ( y ^ i ′ z ^ i ′ - η i ′ ) 2 ≤ ( δ c p size ) 2 ,   ⁢ i = 1 , … ⁢   , n ,

where δ c defines the radius of a circular disk (e.g., 20 pixels) centered at (ξ i , η i ) or (ξ′ i , η′ i ). The constraints C 3 to C 6 form a cone-shape bounding region for each 3D point with the focal spot as its apex and the a circular disk (of radius δ c ) centered at (ξ i , η i ) as its base to confine the searching space of the calculated 3D object point which is shown in FIG. 19 .

D. Establishment of Vessel Centerline Correspondence

Details of the methods for establishing the vessel centerline correspondence of the arterial tree are shown in the flowchart of FIG. 20 .

After the transformation (R, {right arrow over (t)}) that defines the two views is obtained, this information will be used for establishing the point correspondences on 2-D vessel centerlines in the pair of images and calculating 3-D morphologic structures of coronary arterial tree. The calculated transformation in conjunction with the epi-polar constraints were employed as the framework for establishing the point correspondences on the vessel centerlines based on the two identified 2-D coronary arterial trees and the vessel hierarchy as shown in FIG. 21 .

If the information on relative positions of two cameras is known (i.e., the known locations and orientations of two focal spots F a and F b ), the correspondences of image points can be solved by use of “epi-polar constraints.” Let P denote a 3-D point, and let P a and P b , the projections of P, denote the pair of corresponding points in two images. Such constraints state that these points, F a , F b , P a , P b , and P, are all on one plane called the epi-polar plane as shown in FIG. 22 . The epi-polar plane and each image plane intersect along a straight line called the epi-polar line (see FIG. 23 ). Therefore, the location of point P b in the second image must lie on the epi-polar line resulting from the intersection between the second image plane and the plane defined by point P a and the two focal spots. A typical example of local correspondence establishment by use of epi-polar constraints where points A, B, C, and D at the first view are defined by finding the intersections between the 2-D curve associated with left anterior descending artery (LAD) artery and respective epi-polar lines defined by a, b, c, and d at the second view as shown in FIG. 24 .

If the epi-polar line is “tangential” relative to the 2-D vessel centerline curve, the accuracy of calculated intersection point becomes sensitive to any small errors in epi-polar line derived based on the calculated transformation in terms of R and {right arrow over (t)}. To overcome the problem, a refinement process as shown in FIG. 25 will be employed to calculate the global correspondences using optimal parametric arguments t 1 , t 2 , . . . t n−1 , t n based on the equations below. $$ min t ^ , p ′ , v ′ ⁢ ℱ epi ⁢ ( t ^ , P ′ , v ′ ) = ∑ i = 1 n ⁢ ( ξ i ′ - x i ′ z i ′ ) 2 + ( η i ′ - y i ′ z i ′ ) 2 + ( f x ⁢ ( t i ) - c → 1 · p → i ′ + t x c → 3 · p → i ′ + t z ) 2 + ( f y ⁢ ( t i ) - c → 2 · p → i ′ + t y c → 3 · p → i ′ + t z ) 2 + &LeftDoubleBracketingBar; f ′ ⁢ ( t i ) - [ R - 1 · v → i ′ ] xy &RightDoubleBracketingBar; 2 + &LeftDoubleBracketingBar; μ → i ′ - [ v → i ′ ] x ′ ⁢ y ′ &RightDoubleBracketingBar; 2 ,

where {circumflex over (t)}={t 1 , t 2 , . . . , t n }, P′={{right arrow over (p)}′ i =(x′ 1 ,y′ i ,z′ i ) | i=1, . . . , n}, and v′={{right arrow over (v)} 1 , {right arrow over (v)} 2 , . . . , {right arrow over (v)} n } denote the set of parametric variables, the set of 3-D object position vectors, and the set of 3-D vessel directional vectors to be optimized, respectively, n denotes the number of points on the 2-D vessel centerline, (f x (t i ), f y (t i )) and (ξ i , η′ i ) are the respective spline-curve modeled and extracted 2-D vessel centerline points in two views, f′ (t i )=(f′ x (t i ), f′ y (t i )) and {right arrow over (μ)}′ i denote the respective spline-curve modeled and extracted 2-D directional vectors of vessel centerlines in each view, [R −1 {right arrow over (v)} i ] xy and [{right arrow over (v)} i ] x′y′ denote the projection of the calculated 3-D directional vector on the respective image planes, {right arrow over (c)} k , k=1,2,3, denotes the respective kth column vectors of matrix R, and R −1 is the inverse matrix of R. The 3-D position vectors P′≡{{right arrow over (p)} i =(x′ i ,y′ i ,z′ i ) | i=1, . . . , n} on vessel centerline are computed by use of equations described above where the 2-D pair correspondence (ξ,η i ) are substituted by the 2-D spline-curve modeling function f(t i )=(f x (t i ),f y (t i )). Based on the equations, the vessel centerline correspondences are globally established by incorporating the entire vessel shape in terms of directions and locations that will yield more accurate results than those obtained by only use of epi-polar constraints with local vessel segment points; especially when epi-polar line and vessel segment are tangential. FIGS. 26 ( c ) and ( d ) show the results of 3-D LCA coronary reconstructed from a pair of angiograms as shown in FIGS. 26 ( a ) and ( b ) based on the simple epi-polar technique and the refinement process, respectively. The reconstruction of left main artery apparently illustrates inaccurate results based on the simple epi-polar process which were corrected after employing the refinement process.

E. Calculation of 3-D Skeleton of Coronary Arterial Tree

Details of the methods for calculation of the 3-D skeleton of the arterial tree are shown in the flowchart of FIG. 27 .

With the point correspondences on 2-D vessel centerlines (ξ j , η j ) and (ξ′ j , η′ j ) and the transformation (R, {right arrow over (t)}), the 3-D vessel centerline points of coronary arteries (x j , y j , z j ) can then be calculated as shown in equations below. For each artery, the cross-sectional diameters at the two end-points can be calculated as $$ d i = d _ i · ( x i 2 + y i 2 + z i 2 / ξ i 2 + η i 2 + D 2 ) ,

where {right arrow over (d)} is the extracted 2-D vessel diameter centered at (ξ i , η i ), and D is the distance between the focal spot and the image intensifier. The remaining cross-sectional diameters on a artery are linearly interpolated based on the calculated diameters at two end-points as shown in FIG. 29 .

In the first view, let (u i , v i ) denote the image coordinates of the ith object point, located at position {right arrow over (p)} i =(x i , y i , z i ) relative to the focal spot of the first view. Then $$ u i = Dx i z i ,   ⁢ v i = Dy i z i ,

where D is the SID (focal-spot to imaging-plane distance) in the first imaging system. Let (ξ i , η i ) denote the scaled image coordinates, defined as $$ ξ i = u i D = x i z i ,   ⁢ η i = v i D = y i z i .

The second projection view of the biplane imaging system can be described in terms of the SID D′, and a second pair of image and object coordinate systems u′v′ and x′y′z′ can be defined in an analogous manner. Scaled image coordinates (ξ′ i , η′ i ) in the second view for the ith object point at position {right arrow over (p)} i =(x′ i , y′ i , z′ i ) are then given by $$ ξ i ′ = u i ′ D ′ = x i ′ z i ′ ,   ⁢ η i ′ = v i ′ D ′ = y i ′ z i ′ .

With the known estimates of the scaled image coordinates (ξ i , η i ) and (ξ′ i , η′ i ), transformation in terms of R and {right arrow over (t)}, the object coordinates (x i , y i , z i ) can be calculated as follows: $$ ⌊ r 11 - r 31 ⁢ ξ i ′ r 12 - r 32 ⁢ ξ i ′ r 13 - r 33 ⁢ ξ i ′ r 21 - r 31 ⁢ η i ′ r 22 - r 32 ⁢ η i ′ r 23 - r 33 ⁢ η i ′ 1 0 - ξ i 0 1 - η i ⌋ ⁡ [ x i y i z i ] = [ a → · t → b → · t → 0 0 ] ,

where $$ a → = ⌊ r 11 - r 31 ⁢ ξ i ′ r 12 - r 32 ⁢ ξ i ′ r 13 - r 33 ⁢ ξ i ′ ⌋ , b → = ⌊ r 21 - r 31 ⁢ η i ′ r 22 - r 32 ⁢ η i ′ r 23 - r 33 ⁢ η i ′ ⌋ .

F. Rendering of 3-D Coronary Arterial Tree with Associated Gantry Angulation

Details of the methods for rendering the 3-D arterial tree with associated gantry angulation are shown in the flowchart of FIG. 28 .

To obtain the geometrical continuity of the vessel centerline, a spline-based smoothing process is applied to the originally calculated data sets in accord with the equations below. The surface-based modeling technique is employed to recover the complete morphology of the 3-D arterial lumen. For every two consecutive cross sections centered at (x i , y i , z i ) and (x i+1 , y i+1 , z i+1 ), a sequence of polygonal patches are formed with the corresponding cross sectional diameters as shown in FIG. 30 .

When an arbitrary computer-generated image is produced, the gantry information defining the current projection is calculated in the form of left anterior oblique/right anterior oblique (LAO/RAO) and caudal/cranial (CAUD/CRAN). These gantry angles are defined in a spatial coordinate system with the iso-center as the origin. The LAO/RAO angulation is defined on the y-z plane, while the CAUD/CRAN angulation is defined on the x-z plane. With the spatial coordinates, the position of focal spot (x f , y f , z f ) can be formulated by use of two rotations R x (α) and R y (−β) as presented in equations to follow.

The reconstructed 3-D coronary arterial tree can be rotated to any selected viewing angle yielding multiple computer-generated projections to improve the quality and utility of the images subsequently acquired as shown in FIGS. 31 a-f and 32 a-f.

Due to the poor imaging resolution in the depth or z-axis, the 3-D shape of coronary vessel tree become joggling. To overcome the drawbacks, the technique of Bezier curve for data smoothing is adopted. The Bezier curve can be defined as a parametric function as follows: $$ P ⁢ ( u ) = ∑ k = 0 n ⁢   ⁢ p k ⁢ B k , n ⁢ ( u ) , p k = ( x k , y k , z k ) .

where p k denote the n+1 input control points. Each B k,n (u) is a polynomial function of degree one less than the number of control points used and is defined as   B k,n ( u )= C ( n, k ) u k (1− u ) n−k

and the C(n,k) represent the binomial coefficients $$ C ⁢ ( n , k ) = n ! k ! ⁢ ( n - k ) !

By applying the Bezier equation, a new representation of 3-D coronary vessels can be obtained.

When an arbitrary computer-generated image is produced, the gantry information defining the current projection is calculated in the form of left anterior oblique/right anterior oblique (LAO/RAO)and caudal/cranial (CAUD/CRAN). These gantry angles are defined in a spatial coordinate system with the iso-center as the origin. The LAO/RAO angulation is defined on the y-z plane, while the CAUD/CRAN angulation is defined on the x-z plane. With the spatial coordinates, the position of focal spot (x f , y f , z f ) can be formulated by use of two rotations R x (α) and R y (−β) as, $$ ⌊ x f y f z f ⌋ = R y ⁡ ( - β ) ⁢ R x ⁡ ( α ) · ⌊ x n y n z n ⌋ = [ cos ⁡ ( - β ) sin ⁡ ( α ) ⁢ sin ⁡ ( - β ) cos ⁡ ( α ) ⁢ sin ⁡ ( - β ) 0 cos ⁡ ( α ) - sin ⁡ ( α ) - sin ⁡ ( - β ) sin ⁡ ( α ) ⁢ cos ⁡ ( - β ) cos ⁡ ( α ) ⁢ cos ⁡ ( - β ) ] · [ x n y n z n ] ,

where (x n , Y n , z n )=(0, 0, D fc ) denotes the neutral position of focal spot (or anterior-posterior (AP) view with 0° LAO and 0° CAUD, D fc is the distance between iso-center and focal spot which is provided by the manufacture), R x and R y denotes the rigid rotations with respect to x-axis and y-axis, and α(−α) and β(−β) denote the LAO (RAO) and CAUD (CRAN) angles, respectively.

G. 3-D Coronary Quantitative Analysis: Optimal Views and Morphology Assessment

FIG. 37 provides a flowchart describing the high level operations included in the improved QCA of the present invention. Specifically, the method of FIG. 37 presumes the availability of 3-D reconstruction image data of a coronary arterial tree. The methods presented above are useful for generation of such 3-D reconstructed data from 2-D coronary images. Those skilled in the art will recognize a variety of other sources for such 3-D image data as well as alternate techniques for reconstructing such data from 2-D images.

Element 3700 is first operable to determine whether the user of the QCA method desire to assess the geometry of the coronary artery reconstruction or whether the user wishes to determine optimal view of the 3-D image for particular parameters. If the user is assessing the geometry, element 3702 is operable to determine what aspect of the arterial image geometry is to be assessed. Element 3702 is operable to evaluate take-off angles of each arterial branch in the arterial tree image presently under study. Element 3706 is operable to evaluate the insertion angle a bypass graft selected in the arterial tree image presently under study. Element 3708 is operable to evaluate the length of a selected arterial segment from the arterial tree image presently under study. And element 3710 is operable to evaluate the tortuosity of a selected arterial pathway from the arterial tree image presently under study. Following completion of the selected evaluative task processing is completed or may continue by looping back to element 3700 to await a next selection from the arteriographer.

If the user wishes to determine optimal view for particular parameters as determined by element 3700 , element 3712 is next operable to determine which view the user desires for purposes of optimizing. Element 3714 is operable to present the user with a graphical map useful for viewing and optimizing foreshortening of the coronary arterial tree image presently under study. Element 3716 is operable to present the user with a graphical map useful for viewing and optimizing overlap of the arterial tree image presently under study. Element 3718 is operable to present the user with a graphical map useful for optimizing the combination offoreshortening and overlap in the arterial tree image presently under study. Element 3720 is operable to present to the user information indicative of the accumulated percent foreshortening of each arterial segment divided by bifurcation points in the arterial tree image presently under study. Following completion of the selected view optimization, processing is completed or may continue by looping back to element 3700 to await a next selection from the arteriographer.

Each of the above QCA tasks is further described in detail herein below with reference to additional figures as required for a clear understanding.

G.1. Assessment of the Coronary Geometry

The curvilinear nature of vessels is important in the pathophysiology of coronary artery disease and the performance of interventions. Therefore such a geometric property can not be quantified from 2-D projection data. The skeleton of a 3-D vessel can be mathematically defined as a curve function p(s)=(x(s), y(s), z(s)) connecting all the 3-D centerline points as shown in FIG. 33 , where 0≦s≦1 is the parametric variable.

G.1.1. Evaluation of Length of Arterial Segment

As noted above, element 3708 is operable to evaluate the length of a selected artery segment of the arterial tree image presently under study. FIG. 47 is a flowchart describing additional details of the operation of element 3708 to evaluate the length of a selected artery segment in the arterial tree image presently under study. In general, the 3-D vessel centerline is characterized by a parametric curve function in terms of an ordered sequence of 3-D data points as described previously. The length of any selected artery segment on the vessel centerline can be calculated by the summation of the inter-distance of every pair of consecutive vessel centerline points.

More specifically, element 4700 is operable to accept user input defining two points on the displayed arterial tree. The user provides the two points by clicking a pointer device at each of the desired points. The two selected points define a segment of the displayed arterial tree for which the user wishes to evaluate the length. Element 4702 then determines whether a path exists in the arterial tree from the start of the displayed tree through bot selected points. Calculation of the segment length from the 3 -D data requires a path exists from the start of the image tree through both selected points. Details of this element are provided further herein below. Element 4704 then tests whether the required path is found. If not, processing loops back to element 4700 to permit the user to select different points for a desired segment. If the required path through the elected segment exists, processing continues with element 4706 to calculate the length between the two identified points along the path identified by element 4702 .

FIG. 51 is a flowchart providing additional details of the processing of element 4702 (and 4802 discussed below) to determine the connectivity of a selected arterial pathway defined by two selected end points. In particular, element 5100 uses the vessel hierarchy data structure to start at the distal end point of the selected pathway. Element 5102 then traverses the vessel centerline point sequence from the present point in the upstream direction. Element 5104 then determines whether the selected proximal end point of the selected pathway is encountered in the next point traversed by element 5102 . If the proximal end point is encountered, then the selected pathway is connected and element 5106 is operable to return a “CONNECTED” status to the invoking procedure. If the proximal end point is not yet encountered, element 5108 is operable to determine whether the beginning of the coronary tree is encountered in this next point traversed by element 5102 . If so, element 5110 is operable to return a “NOT CONNECTED” status to the invoking procedure because the beginning of the tree is encountered without traversing through the proximal end point. If the beginning of the tree is not yet encountered, processing loops back to element 5102 to traverse a next point along the selected path and repeat the tests of elements 5104 and 5108 .

G. 1.2. Evaluation of the tortuosity of artery

Element 3710 is operable to evaluate the tortuosity of a selected arterial segment pathway. FIG. 48 is a flowchart providing additional details of the processing to evaluate tortuosity of a selected arterial segment.

G. 1.2.1. Evaluation of the tortuosity of an artery by the F-S theory

Generally, by use of such spline-based curve modeling technique, we are able to apply the theory of differential geometry such as Frenet-Serret apparatus (or F-S theory) or its variation to study the geometrical nature of the 3-D coronary artery shape at any time point in the cardiac cycle. A technique based on the Frenet-Serret (F-S) theory of differential geometry has been developed to study the geometrical nature or tortuosity of the 3-D coronary artery shape at any time point in the cardiac cycle. The F-S theory consists of five components: three vector fields along the given curve (the tangent T(s), the normal N(s), and the bi-normal B(s) vectors) and two scalar valued functions (the curvature κ(s) and the torsion τ(s) where s denotes the parametric variable defining the location of point on the curve function α(s). The curvature κ(s 0 ) measures the rate of change of the angle defined by the two neighboring tangents T(s 0 ) and T(s 0 +δ s ) associated with the two respective points α(s 0 ) and α(s 0 +δ s ). In other word, it characterizes the degree of bending pertinent to a local curve segment in the neighborhood between s 0 and s 0 +δ s (i.e., how rapidly the curve pulls away from the plane Ω n perpendicular to the tangent vector at T(s 0 )) as shown in FIG. 34 a . Similarly, the torsion at τ(s 0 ) measures a rate of change (or a twisting) at a point α(s 0 ) how its trajectory twists out of the plane Ω n perpendicular to the normal vector B(s 0 )) as shown in FIG. 34 b.

In FIG. 34 c , a circle defining on the x-y plane of radius r>0 is illustrated, where the vector fields T(s) and N(s) vary along the curve and the binormal field B(s) is constant. Note that the curvature κ(s) (for all s) of a circle is equal to the reciprocal of its radius (i.e., 1/r). Since the binormal field B is a constant, the torsion at every point is equal to 0. The example illustrated in FIG. 34 d shows a circular helix in which has constant curvature and torsion at every point in the curve. Given the curvature and torsion functions κ(s) and τ(s), they completely describe the curve α(s) geometrically.

The calculation of curvature and torsion at every vessel centerline point s 0 is characterized by the following equations:

T ( s 0 )=α (1) ( s 0 )/|α (1) ( s 0 )|

B ( s 0 )=α (1) ( s 0 )×α (2) ( s 0 )/|α (1) ( s 0 )×α (2) ( s 0 )|

N ( s 0 )= B ( s 0 )× T ( s 0 )

κ( s 0 )=|α (1) ( s 0 )×α (2) ( s 0 )|/|α (1) ( s 0 )| 3

τ( s 0 )=[α (1) ( s 0 )×α (2) ( s 0 )]·α (3) ( s 0 )/|α (1) ( s )×α (2) ( s 0 )| 2 ,

where α (1) (s) denotes the ith derivative with respect to s. In FIG. 38 a , the length of the selected pathway is 10.0 cm. The narrowest diameter of the diseased arterial segment is 3.4 mm. The measurement of tortuosity associated with this pathway was illustrated in FIG. 38 b.

1.2.2. Evaluation of the Tortuosity of Artery by Point-to-Point Method

A new method, which is a variation of F-S theory, is employed to characterize the vessel tortuosity in order to facilitate the study of vessel straightening due to various intra-coronary devices. Three types of vessel tortuosity are defined: (1) point-to-point (or local) tortuosity, (2) point-to-line (or relative) tortuosity, and (3) tortuosity spectrum. Given an arterial segment, the point-to-point (P-to-P) tortuosity at each vessel centerline point is estimated by calculation of the differential of tangent vectors (i.e., the angle in degree between the two tangent vectors) associated with two consecutive arterial points on the 3-D centerline curve. The point-to-line (P-to-L) tortuosity at each vessel centerline point is defined by the differential of the tangent vector and the directional vector of the line connecting the two end points of the given arterial segment (i.e., the angle in degree between the tangent vector of a vessel centerline point and the directional vector of line). Afterwards, the mean values of the P-to-P tortuosity and P-to-L tortuosity are calculated. To evaluate the “global variation” of vessel geometry, the tortuosity spectrum is derived by calculation of the normalized integration of P-to-L tortuosity along the arterial segment.

FIG. 39 a shows the reconstructed 3-D bypass graft. In FIG. 39 b the P-to-P (bottom) and P-to-L (the 2nd graph from the bottom) tortuosity measurements associated with the selected 3 arterial segments defined between the proximal and distal locations in FIG. 39 a are illustrated, respectively. The mean values of regional P-to-P and P-to-Line are (7.6,6.8, 1.9) and (37.5,27.5,9.0) degrees, respectively. The total mean values of P-to-P and P-to-Line are 5.4 and 23.3 degrees, respectively. The tortuosity spectrum in terms of the normalized integration of P-to-L measurements and the diameters along the arterial segments are shown at the 2nd graph from the top and the top graphs in FIG. 39 b.

Referring again to FIG. 48 , elements 4800 . . . 4804 are operable in a manner analogous to elements 4700 . . . 4704 of FIG. 47 . Specifically, the user is prompted to enter two points on the arterial tree image presented to define a segment of interest. When a valid segment is so selected, element 4806 is then operable to determine whether the user wishes to evaluate the tortuosity in accordance with the point-to-point method or the F-S method both described above. If the user selects the point-to-point (p-to-p) method, elements 4808 and 4810 are operable to calculate the point-to-point regional curvature and point-to-line regional curvature, respectively, of each point on the vessel centerline as described above. If the F-S method is selected by the user, elements 4812 and 4814 are operable to calculate the curvature and torsion, respectively, of each vessel centerline point along the selected segment pathway as described above.

G. 1.3. Evaluation of the Take-off Angle

Element 3704 is operable to evaluate the take-off angles of a selected arterial segment pathway. FIG. 49 is a flowchart providing additional details of the processing to evaluate take-off angles of a selected arterial segment. The 3-D take-off angle of a branch vessel is estimated by calculating the angle between the two tangent vectors originating from the bifurcation; i.e., one toward the direction of precedent vessel and the other toward the direction of branching vessel as shown in FIG. 35 a . The 2-D take-off angle of a branch vessel is obtained by first computing the branch vessel and preceding vessel centerlines projected onto the imaging plane based on the gantry orientation (α,β) followed by calculating the angle of the 2-D tangent vectors as shown in FIG. 35 b . When an arbitrary computer-generated image is produced, the gantry information defining the current projection is calculated in the form of left anterior oblique/right anterior oblique (LAO/RAO) and caudal/cranial (CAUD/CRAN). Mathematically, these gantry angles are defined in a spatial coordinate system with the iso-center as the origin. The LAO/RAO angulation is defined on the y-z plane, while the CAUD/CRAN angulation is defined on the x-z plane. With the gantry viewing angle (α,β), the 2-D projection coordinates (u (α,β) ,ν (α,β) of a vessel centerline point can be calculated by use of the following equations: $$ ⌊ u ( α , β ) v ( α , β ) ⌋ = ⌊ x ( α , β ) / z ⁡ ( α , β ) x ( α , β ) / z ⁡ ( α , β ) ⌋ ⁢ [ x ( α , β ) y ( α , β ) z ( α , β ) ] = R y ⁡ ( - β ) ⁢ R x ⁡ ( α ) · [ x y z ] = [ cos ⁡ ( - β ) sin ⁡ ( α ) ⁢ sin ⁡ ( - β ) cos ⁡ ( α ) ⁢ sin ⁡ ( - β ) 0 cos ⁡ ( α ) - sin ⁡ ( α ) - sin ⁡ ( - β ) sin ⁡ ( α ) ⁢ cos ⁡ ( - β ) cos ⁡ ( α ) ⁢ cos ⁡ ( - β ) ] · [ x y z ] ,

where R x and R y denotes the rigid rotations with respect to x-axis, (x (α,β) , γ (α,β) , z (α,β) ) denotes the 3-D coordinates of a vessel centerline point defined at the current gantry coordinate system, (x,y,z) denotes the 3-D coordinates of a vessel centerline point in the AP viewing angle (0 degree LAO and 0 degree CAUD) coordinate system, and α (−α) and β (−β) denote a gantry position at LAO (RAO) and CAUD (CRAN) angles, respectively. In FIGS. 40 ( a ) and ( b ), the 3-D and 2-D take-off angles of individual coronary arteries are shown based on two different gantry viewing angles.

The above equation can be employed to model the majority of existing imaging systems having only two degrees of freedom. Another variable R z (γ), where R z denotes the rigid rotation with respect to z-axis by γ degrees, can be readily added to emulate the imaging system if a rotational C-arm base is employed. With the new third degree of freedom in term of the new angle variable, the above equation can be re-written as follows. $$ ⌊ u ( α , β , γ ) v ( α , β , γ ) ⌋ = ⌊ x ( α , β , γ ) / z ⁡ ( α , β , γ ) x ( α , β , γ ) / z ⁡ ( α , β , γ ) ⌋ ⁢ [ x ( α , β , γ ) y ( α , β , γ ) z ( α , β , γ ) ] = R y ⁡ ( - β ) ⁢ R x ⁡ ( α ) ⁢ R z ⁡ ( γ ) · [ x y z ] = [ C ⁡ ( - β ) ⁢ C ⁡ ( γ ) + S ⁡ ( α ) ⁢ S ⁡ ( - β ) ⁢ S ⁡ ( γ ) - S ⁡ ( γ ) ⁢ C ⁡ ( - β ) + S ⁡ ( α ) ⁢ S ⁡ ( - β ) ⁢ C ⁡ ( γ ) C ⁡ ( α ) ⁢ S ⁡ ( - β ) C ⁡ ( α ) ⁢ S ⁡ ( γ ) C ⁡ ( α ) ⁢ S ⁡ ( γ ) - S ⁡ ( α ) - S ⁡ ( - β ) ⁢ C ⁡ ( γ ) + S ⁡ ( α ) ⁢ C ⁡ ( - β ) ⁢ S ⁡ ( γ ) S ⁡ ( - β ) ⁢ S ⁡ ( γ ) + S ⁡ ( α ) ⁢ C ⁡ ( - β ) ⁢ C ⁡ ( γ ) C ⁡ ( α ) ⁢ C ⁡ ( - β ) ] · [ x y z ] ,

where C(x) and S(x) denote the cos and sin functions respectively.

FIG. 49 is a flowchart providing additional details of the processing of element 3704 of FIG. 37 comprising a method for computing take-off angles in accordance with a 2-D or 3-D model. Element 4900 is first operable to determine which model the user has requested. If a 3-D model is used, elements 4902 and 4904 compute the two required tangent vectors from the bifurcation point downstream and from the bifurcation point upstream, respectively. If a 2-D model is used, elements 4906 and 4908 calculate the required two tangent vectors as originated from the bifurcation point downstream in the direction of the branch and from the bifurcation point upstream in the direction of the preceding vessel, respectively. In both models, element 4910 then computes the take-off angel as the angle between the two tangent vectors as determined in accordance with the requested model.

G. 1.4. Evaluation of the Insertion Angle

Element 3706 is operable to evaluate the insertion angle of a selected arterial segment pathway. FIG. 50 is a flowchart providing additional details of the processing to evaluate the insertion angle of a selected arterial segment. The insertion angles (upstream and downstream angles) of bypass graft at the end-to-side anastomosis are estimated by calculation of the angles between the 3-D tangent vectors at the native vessels and the tangent vector at the graft as shown in FIG. 36 . Element 5000 of FIG. 50 is operable to receive operator input selecting two end points on the arterial tree having the anastomosis located between them. Element 5002 then receives operator input to select one point on the native artery connecting to the anastomosis. Element 5004 then calculates the insertion angle as the upstream and downstream angles defined by the three selected points.

G.2. Determination of Optimal View Strategy

When an arbitrary computer-generated image is produced, the gantry information defining the current projection is calculated in the form of left anterior oblique/right anterior oblique (LAO/RAO)and caudal/cranial (CAUD/CRAN). These gantry angles are defined in a spatial coordinate system with the iso-center as the origin. The LAO/RAO angulation is defined on the y-z plane, while the CAUD/CRAN angulation is defined on the x-z plane. With the spatial coordinates, the position of focal spot (x f , y f , z f ) can be formulated by use of two rotations R x (α) and R y (−β) as, $$ ⌊ x f y f z f ⌋ = R y ⁡ ( - β ) ⁢ R x ⁡ ( α ) · ⌊ x n y n z n ⌋ = [ cos ⁡ ( - β ) sin ⁡ ( α ) ⁢ sin ⁡ ( - β ) cos ⁡ ( α ) ⁢ sin ⁡ ( - β ) 0 cos ⁡ ( α ) - sin ⁡ ( α ) - sin ⁡ ( - β ) sin ⁡ ( α ) ⁢ cos ⁡ ( - β ) cos ⁡ ( α ) ⁢ cos ⁡ ( - β ) ] · [ x n y n z n ] ,

where (x n , y n , z n )=(0, 0, D fc ) denotes the neutral position of focal spot or anterior-posterior (AP) view with 0° LAO and 0° CAUD, D fc is the distance between iso-center and focal spot which is provided by the manufacture), R x and R y denotes the rigid rotations with respect to x-axis andy-axis, and α (−α) and β (−β) denote the LAO (RAO) and CAUD (CRAN) angles, respectively.

The above equation can be employed to model the majority of existing imaging systems having only two degrees of freedom. Another variable R z (γ), where R z denotes the rigid rotation with respect to z-axis by γ degrees, can be readily added to emulate the imaging system if a rotational C-arm base is employed. With the new third degree of freedom in term of the new angle variable, the above equation can be re-written as follows. $$ ⌊ x f y f z f ⌋ = R y ⁡ ( - β ) ⁢ R x ⁡ ( α ) ⁢ R z ⁡ ( γ ) · ⌊ x n y n z n ⌋ = [ C ⁡ ( - β ) ⁢ C ⁡ ( γ ) + S ⁡ ( α ) ⁢ S ⁡ ( - β ) ⁢ S ⁡ ( γ ) - S ⁡ ( γ ) ⁢ C ⁡ ( - β ) + S ⁡ ( α ) ⁢ S ⁡ ( - β ) ⁢ C ⁡ ( γ ) C ⁡ ( α ) ⁢ S ⁡ ( - β ) C ⁡ ( α ) ⁢ S ⁡ ( γ ) C ⁡ ( α ) ⁢ S ⁡ ( γ ) - S ⁡ ( α ) - S ⁡ ( - β ) ⁢ C ⁡ ( γ ) + S ⁡ ( α ) ⁢ C ⁡ ( - β ) ⁢ S ⁡ ( γ ) S ⁡ ( - β ) ⁢ S ⁡ ( γ ) + S ⁡ ( α ) ⁢ C ⁡ ( - β ) ⁢ C ⁡ ( γ ) C ⁡ ( α ) ⁢ C ⁡ ( - β ) ] · [ x n y n z n ] ,

where C(x) and S(x) denote the cos and sin functions respectively.

Given the 3-D character of the coronary artery tree, it is expected that any projection will foreshorten a variety of segments. The visual perception provided by such a projection shows a variety of arterial segments with a degree of foreshortening that may or may not be appreciated. In FIG. 41 a , the graphical example shows that the projected object has 0% degree of foreshortening since the tangential direction of such a line segment is perpendicular to the normal vector of image plane. Note that the projections of the two line segments have various magnifications due to the different locations relative to the image intensifier. If the two line segments are rotated by an arbitrary angle, e.g., 45 degrees relative to their original orientation in FIG. 41 a , there will be approximately 30% degree of foreshortening as shown in FIG. 41 b . In FIG. 41 c , a synthetic example of vessel foreshortening for a typical segment of a curvilinear arterial centerline is illustrated.

Due to the problem of vessel overlap and foreshortening, multiple projections are necessary to adequately evaluate the coronary arterial tree with arteriography. The elimination or at least minimization of foreshortening and overlap is a prerequisite for an accurate quantitative coronary analysis (QCA) such as determination of intra coronary lengths in a 2-D display. Hence, the patient might receive substantial radiation and contrast material during diagnostic and interventional procedures. This traditional trial-and-error method provides views in which overlapping and foreshortening are minimized, but only in terms of the subjective experience-based judgment of the angiographer. In FIG. 42 a and 42 b , the same vessel segment between two bifurcation points as marked by two dots at its two distal ends depicts 77% and 52% foreshortening in the respective angiograms.

The reconstructed 3-D coronary arterial tree can be rotated to any selected viewing angle yielding multiple computer-generated projections to determine for each patient which standard views are useful and which are of no clinical value due to excessive overlap. Therefore, the 3-D computer assistance provides means to improve the quality and utility of the images subsequently acquired.

G.2. 1. Calculation of Percent Foreshortening Map

Element 3714 is operable to present foreshortening information for an arterial tree to permit a user to optimize the view for foreshortening. FIG. 52 is a flowchart providing additional details of the processing to optimize the view of an arterial tree 3-D projection for reducing foreshortening. Let p i , i=0, 1, . . . , m denote the points on the 3-D vessel centerline curve {haeck over (L)}. Let {right arrow over (l)} j =[l j x , l j y , l j z ] t , j=1, 2, . . . , m denote the vector formed by every two consecutive vessel centerline points p j−1 and p j . The percent foreshortening of a vessel segment {haeck over (L)} based on any specific gantry orientation α and β angles can be obtained as follows: $$ ℱ fores ⁢ ( α , β , L ⋓ ) =   ⁢ ∑ j = 1 m ⁢   ⁢ &LeftBracketingBar; l → j &RightBracketingBar; ⁢ cos ⁢ ( θ j ) ∑ j = 1 m ⁢   ⁢ &LeftBracketingBar; l → j &RightBracketingBar; , ⁢ 0 ⁢ ° ≤ θ j ≤ 180 ⁢ ° , - 120 ⁢ ° < α < 120 ⁢ ° , - 60 ⁢ ° < β < 60 ⁢ ° .

where θ is the angle between the directional vector {right arrow over (l)} j and projection vector {right arrow over (z)} p =[−cos(α)sin(β), −sin(α),cos(α)cos(β)] t , and the constraints associated with the angle parameters α and β pertinent to the achievable gantry angles of imaging system. If there are k multiple arterial segments to be evaluated (e.g., a bifurcation region), the aforementioned equation is extended as follows: $$ ℱ fores ⁢ ( α , β , L ⋓ ) = ∑ i = 1 k ⁢   ⁢ ∑ j = 1 m i ⁢   ⁢ &LeftBracketingBar; l → j i &RightBracketingBar; ⁢ cos ⁢ ( θ j ) ∑ i = 1 k ⁢   ⁢ ∑ j = 1 m i ⁢   ⁢ &LeftBracketingBar; l → j i &RightBracketingBar; , ⁢ 0 ⁢ ° ≤ θ j i ≤ 180 ⁢ ° , - 120 ⁢ ° < α < 120 ⁢ ° , - 60 ⁢ ° < β < 60 ⁢ ° .

where {haeck over (L)}={L 1 , L 2 , . . . , L k } denote the set of arterial segments and each arterial segment of interest L i ={{right arrow over (l)} 1 i , {right arrow over (l)} 2 i , . . . , {right arrow over (l)} i m i } consists of (m i +1) centerline points, i=1, 2, . . . , k.

The above equations can be employed to model the majority of existing imaging systems having only two degrees of freedom. Another variable γ, denoting the rigid rotation with respect to z-axis by γ degrees, can be readily added to emulate the imaging system if a rotational C-arm base is employed. With the new third degree of freedom in term of the new angle variable, the above equation can be re-written as follows. $$ 𝔉 fores ⁡ ( α , β , γ , L ˘ ) = ∑ j = 1 m ⁢ &LeftBracketingBar; 1 -> j &RightBracketingBar; ⁢ cos ⁢ ( θ j ) ∑ j = 1 m ⁢ &LeftBracketingBar; 1 -> j &RightBracketingBar; , ⁢ 0 ⁢ ° ≤ θ j ≤ 180 ⁢ ° , - 120 ⁢ ° < α < 120 ⁢ ° , - 60 ⁢ ° < β < 60 ⁢ ° , - 90 ⁢ ° < γ < 90 ⁢ ° . ⁢ 𝔉 fores ⁡ ( α , β , γ , L ˘ ) = ∑ i = 1 k ⁢ ∑ j = 1 m 1 ⁢ &LeftBracketingBar; 1 -> j i &RightBracketingBar; ⁢ cos ⁢ ( θ j i ) ∑ i = 1 k ⁢ ∑ j = 1 m 1 ⁢ &LeftBracketingBar; 1 -> j i &RightBracketingBar; , ⁢ 0 ⁢ ° ≤ θ j i ≤ 180 ⁢ ° , - 120 ⁢ ° < α < 120 ⁢ ° , - 60 ⁢ ° < β < 60 ⁢ ° , - 90 ⁢ ° < γ < 90 ⁢ ° .

Based on the percent foreshortening equation (α, β, {haeck over (L)}), the percent foreshortening map relative to the arterial segments of interest {haeck over (L)} can be generated. In the foreshortening map, the native and positive x-axes denote the RAO and LAO angles. Similarly, the native and positive y-axes denote the CAUD and CRAN angles, respectively. Each element of the map Map f (α,β) contains the calculated percent foreshortening measurement. The contents of Map f (α,β) is further coded by 5 colors with the spectrum from 0 to 255 for each color that is inverse proportional to the magnitude of percent foreshortening. Five colors are employed for representing various ranges of percent foreshortening measurement including: white color (0-5)%, yellow color (6-14)%, green color (15-29)%, blue color (30-49)%, and red color (50-100)%, respectively. FIG. 44 b shows the foreshortening map of a selected arterial segment of interest corresponding to a stenotic lesion in FIG. 44 a. In FIG. 45 b, the map corresponding to a bifurcation lesion of FIG. 46 a as region of interest is produced. FIG. 46 b demonstrates the multiple foreshortening maps with respect to the three arterial segments of interest selected from a bypass graft.

Referring again to FIG. 52 , element 5200 is operable to rotate the 3-D projection of the arterial tree so that a selected region of interest (“ROI”) is visible on the user's display screen. Element 5202 then determines whether the selected ROI represents a single segment of the arterial tree or includes a bifurcation. If the ROI includes a bifurcation, element 5204 is operable to receive user input to select two end points on the displayed projection such that the bifurcation is encompassed within the segment defined by the two selected end points. Element 5206 is operable to receive user input to select an additional point on the branching vessel of the bifurcation. Element 5208 is then operable to verify that the selected three points are correct in that they are within the region of interest and correctly define the bifurcated branches. If the points entered are not correct, processing continues by looping back to element 5204 to await correct input from the user. If the entered points are correct, processing continues with element 5214 discussed below.

If element 5202 determines that the ROI comprises a single segment of the arterial tree, element 5210 is operable to select two end points defining the segment of the ROI to be optimized for foreshortening. Element 5212 then verifies that the selected segment is correct in that the two points define a continuous segment within the ROI. If the points are correct, processing continues with element 5214 below. If not, processing loops back to element 5210 to await entry by a user of correct points defining the segment to be displayed.

Element 5214 is then operable to calculate the accumulated percent vessel foreshortening over all points of the selected segment or bifurcation for all viewing angles. Element 5216 then displays the calculated values as a color coded map as discussed above. The color coding of the foreshortening map allows the user to rapidly assess the optimal view of the segment or bifurcation in the ROI. The user then points to the portion of the map for the desired view of the ROI and segment or bifurcation. The user may view all other (sub-optimal) views by simply pointing to other areas of the displayed map. The 3-D projection of the arterial tree (the ROI) is moved accordingly to present the corresponding view angles.

G.2.2. Calculation of Accumulated Percentage Foreshortening of Entire Coronary Tree

Element 3720 is operable to calculate the accumulated foreshortening percentage for an entire segment or bifurcated segment. FIG. 53 is a flowchart providing additional details of the processing to perform such computations. The calculation of the percent foreshortening of the coronary tree is performed by calculating the accumulated foreshortening (a) on the multiple segments of an artery divided by bifurcation points or (b) on the entire artery if it has no branching vessel. If an artery has N branching vessels, there will be N+1 segments associated with the artery. The calculated results of accumulated foreshortening of each segment is displayed with the true 3-D length in millimeters near the locations of every distal-end and bifurcations of the projection according to the current gantry viewing angle. In FIGS. 43 ( c ) and 43 ( d ), the overall percent foreshortening of the entire coronary tree are illustrated. For example, the true length (4.8 mm) of the arterial segment from the distal-end to the bifurcation on the mid-LAD is foreshortened by 4% on the current gantry angle LAO 43.7, CRAN 34.1 in FIG. 43 ( c ). Similarly, the arterial segment of true length 3.6 mm defined by the bifurcations of the first and second diagonal arteries has 16% and 3% foreshortening as shown in FIGS. 43 ( c ) and 43 ( d ), respectively.

Element 5300 is operable to start with the first artery in the arterial tree using the vessel hierarchy data structure. Starting at that first artery, element 5302 is then operable to determine whether there are any bifurcation points on the current artery. If note, processing skips to element 5306 below. If there are bifurcations on the current artery, element 5304 is next operable to divide the current vessel into N+1 arterial segments based on the N bifurcation points on the current artery. Element 5306 then calculates the percent foreshortening of each arterial segment for all gantry angles. Element 5308 then determines whether any additional arterial segments remain to be processed by further calculations of element 5306 . If additional segments remain to be processed, element 5310 is then operable to select the next segment of the current artery and processing loops back to element 5302 to process the next segment. Elements 5302 through 5310 are therefore iteratively operable on all segments of the selected artery of the arterial tree. When all segments have been processed, element 5312 is operable to display the results as a color coded map as described above.

G.2.3. Calculation of Percent Overlap Map

Element 3716 is operable to present overlap information for an arterial tree to permit a user to optimize the view for overlap. FIG. 54 is a flowchart providing additional details of the processing to optimize the view of an arterial tree 3-D projection for reducing foreshortening. Let p i , i=0, 1, . . . , m denote the points on the 3-D vessel centerline curve {haeck over (L)}. Let {right arrow over (l)} j =[l j x , l j y , l j z ] t , j=1, 2, . . . , m denote the vector of every two consecutive vessel centerline points p j−1 and p j . Based on any specific gantry orientation α and β angles, the percent vessel overlap of a vessel segment {haeck over (L)} relative to all the arteries ζ, i=1, . . . , n, in the coronary tree T after projection can be obtained as follows: $$ 𝔉 ovlp ⁢ ( α , β , L ˘ , T ˘ ) = min ⁡ ( ∑ i = 1 , i ≠ k n ⁢ { Π α , β ⁡ ( Ϛ i ) ⋂ Π α , β ⁡ ( Ϛ 1 -> k ) } , Π α , β ⁡ ( Ϛ 1 -> k ) ) Π α , β ⁡ ( Ϛ 1 -> k ) · 100 ⁢ %

subject to the constraints pertinent to the achievable gantry angles of imaging system

−120°<α<120°, −60°<β<60°,

where ζ l , and ζ k l , denote the 3-D arterial lumen surfaces of the ith artery and of the selected segment-of-interest at the kth artery (or {haeck over (L)}), respectively. Π α,β (O) denotes an operator that counts number of pixels after perspectively projecting the object O onto the image plane based on the gantry angles (α,β), and min(V a , V b ) denotes a binary operator that returns the smaller value between the two variables V a and V b .

The above equations can be employed to model the majority of existing imaging systems having only two degrees of freedom. Another variable γ, denoting the rigid rotation with respect to z-axis by γ degrees, can be readily added to emulate the imaging system if a rotational C-arm base is employed. With the new third degree of freedom in term of the new angle variable, the above equation can be re-written as follows. $$ 𝔉 ovlp ⁢ ( α , β , γ , L ˘ , T ˘ ) = min ⁢ ( ∑ i = 1 , i ≠ k n ⁢ { Π α , β , γ ⁡ ( Ϛ i ) ⋂ Π α , β , γ ⁡ ( Ϛ 1 -> k ) } , Π α , β , γ ⁡ ( Ϛ 1 -> k ) ) Π α , β ⁡ ( Ϛ 1 -> k ) · 100 ⁢ %

subject to the constraints pertinent to the achievable gantry angles of imaging system

−120°<α<120°, −60°<β<60°, −90°<γ<90°,

Based on the percent overlap equation (α,β, {haeck over (L)},T), the map relative to the selected arterial segments ofinterest {haeck over (L)} can be generated. In the overlap map Map o (α, β), the native and postive x-axes denote the RAO and LAO angles. Similarly, the native and postive y-axes denote the CAUD and CRAN angles, respectively. Each element of the map Map o (α, β) contains the calculated percent vessel overlap measurement. Similarly, the contents of Map o (α, β) is further coded by 5 colors with the spectrum from 0 to 255 for each color that is inversely proportional to the magnitude of percent overlap. Five colors are utilized for illustration of various ranges of percent overlap measurement including: white color (0-5)%, yellow color (6-14)%, green color (15-29)%, blue color (30-49)%, and red color (50-100)%, respectively. FIG. 44 c shows the percent overlap map of a selected arterial segment of interest corresponding to a stenotic lesion. In FIG. 45 c, the percent overlap map corresponding to a bifurcation lesion as region of interest is produced. FIG. 46 c demonstrates multiple percent overlap maps with respect to the three arterial segments of interest selected from a bypass graft.

FIG. 54 is a flowchart describing the processing to calculate and display a color coded map for optimizing the view for overlap analysis. The flowchart of FIG. 54 is identical to that of FIG. 52 except the calculations are for overlap as described in the above overlap equations. Elements 5400 through 5416 perform analogous processing to that of elements 5200 through 5216 , respectively.

G.2.4. Calculation of Composite Map

The ultimate goal of optimal view strategy is to minimize the vessel overlap as well as foreshortening on the view(s) acquired from the gantry orientation with respect to any arterial segment(s) of interest. To achieve such a goal, the composite map Map c (α,β) is derived by integration of the respective percent foreshortening map Map f (α,β), and overlap map Map o (α,β). The resultant hybrid information is color coded by use of the percent foreshortening color scheme as background super-imposed by the percent overlap color. The percent foreshortening color scheme is characterized by four different colors with 0 to 255 magnitude including: white (0-7)%, yellow (8-15)%, green (16-31)%, red (32-100)%. The magnitude of percent foreshortening is inversely proportional to that of the color magnitude. The percent overlap color scheme is characterized by a translucent white color toward to opacified blue colors that is proportional to the magnitude of percent overlap (0-100)%. In other word, the background color (i.e., percent foreshortening color) at location (α,β) in the composite map is rendered if there occurs 0 percent overlap. The larger the magnitude of percent overlap at Ma c (α,β), the stronger the opacified blue color is illustrated. FIG. 44 d shows the percent overlap map of a selected arterial segment of interest corresponding to a stenotic lesion. In FIG. 45 d, the composite map corresponding to a bifurcation lesion as region of interest is illustrated. FIG. 46 d illustrates multiple composite maps with respect to the three arterial segments of interest selected from a bypass graft.

FIG. 55 is a flowchart describing the processing to calculate and display a color coded map for optimizing the view for both overlap and foreshortening analysis. The flowchart of FIG. 55 a combination of the processing of FIGS. 52 and 54 discussed above. Specifically, elements 5500 through 5512 perform analogous processing to elements 5200 through 5212 and 5400 through 5412 , respectively. Element 5514 is identical to element 5214 and element 5516 is identical to element 5414 . Element 5518 then displays a composite color coded map presenting information for optimal viewing of both foreshortening and overlap.

While the invention has been illustrated and described in detail in the drawings and foregoing description, such illustration and description is to be considered as exemplary and not restrictive in character, it being understood that only the preferred embodiment and minor variants thereof have been shown and described and that all changes and modifications that come within the spirit of the invention are desired to be protected.

Claims

1. A method for three-dimensional reconstruction of a target object from two-dimensional images, said target object having a plurality of branch-like vessels, the method comprising the steps of:a) placing the target object in a position to be scanned by an imaging system, said imaging system having a plurality of imaging portions; b) acquiring a plurality of projection images of the target object, each imaging portion providing a projection image of the target object, each imaging portion disposed at an unknown orientation relative to the other imaging portions; c) identifying two-dimensional vessel centerlines for a predetermined number of the vessels in each of the projection images; d) creating a vessel hierarchy data structure for each projection image, said data structure including the identified two-dimensional vessel centerlines defined by a plurality of data points in the vessel hierarchy data structure; e) calculating a predetermined number of bifurcation points and directional vectors for each projection image by traversing the corresponding vessel hierarchy data structure, said bifurcation points defined by intersections of the two-dimensional vessel centerlines and said directional vectors defined by tangent vectors on said bifurcation points; f) determining a transformation in the form of a rotation matrix and a translation vector utilizing the bifurcation points and the directional vectors corresponding to each of the projections images, said rotation matrix and said translation vector representing imaging parameters corresponding to the orientation of each imaging portion relative to the other imaging portions of the imaging system; g) utilizing the data points and the transformation to establish a correspondence between the two-dimensional vessel centerlines corresponding to each of the projection images such that each data point corresponding to one projection image is linked to a data point corresponding to the other projection images, said linked data points representing an identical location in the vessel of the target object; h) calculating three-dimensional vessel centerlines utilizing the two-dimensional vessel centerlines and the correspondence between the data points of the two-dimensional vessel centerlines; and i) reconstructing a three-dimensional visual representation of the target object based on the three-dimensional vessel centerlines and diameters of each vessel estimated along the three-dimensional centerline of each vessel.

2. The method of claim 1 wherein the step of identifying vessel centerlines includes the step of:identifying two-dimensional vessel centerlines for a predetermined number of landmarks in each of the projection images.

3. The method of claim 1 wherein the step of utilizing the data points and the transformation to establish a correspondence comprises the step of:utilizing the data points and the transformation to establish local and global correspondences between the two-dimensional vessel centerlines corresponding to each of the projection images such that each data point corresponding to one projection image is linked to a data point corresponding to the other projection images, said linked data points representing an identical location in the vessel of the target object.

4. The method of claim 1 further comprising the step of:identifying a predetermined number of landmarks in each of the projection images wherein said landmarks are man-made objects identified in the projection images, wherein the step of determining said transformation comprises the step of: f) determining a transformation in the form of a rotation matrix and a translation vector utilizing the landmarks and utilizing the bifurcation points and utilizing the directional vectors corresponding to each of the projections images, said rotation matrix and said translation vector representing imaging parameters corresponding to the orientation of each imaging portion relative to the other imaging portions of the imaging system.

5. The method of claim 1 further comprising the step of:identifying a predetermined number of landmarks in each of the projection images wherein said landmarks are points and curves identified in the projection images, wherein the step of determining said transformation comprises the step of: f) determining a transformation in the form of a rotation matrix and a translation vector utilizing the landmarks and utilizing the bifurcation points and utilizing the directional vectors corresponding to each of the projections images, said rotation matrix and said translation vector representing imaging parameters corresponding to the orientation of each imaging portion relative to the other imaging portions of the imaging system.

6. A method for quantitative analysis of a vascular tree comprising the steps of:providing a 3-D vessel hierarchy data structure describing the vascular tree; displaying a 2-D projection of said 3-D vessel hierarchy data structure on a user's display screen; receiving input from the user to select a desired segment of said vascular tree wherein said selected segment includes any number of contiguous segments forming a sub-tree of said vascular tree; receiving input from the user to select a desired quantitative measure of said selected segment wherein said desired quantitative measure is a value determinable from data derived from said 3-D vessel hierarchy data structure; calculating said quantitative measure using said data derived from said 3-D vessel hierarchy; and displaying said quantitative measure on the users display screen.

7. The method of claim 6 further comprising the steps of:receiving input from the user to select an alternate viewing angle of said vascular tree; computing an updated 2-D projection of said 3-D vessel hierarchy data structure based on said alternate viewing angle; and displaying said updated 2-D projection on the user's display screen.

8. The method of claim 6wherein the step of receiving input from the user to select said desired quantitative measure includes the step of receiving user input to evaluate a take-off angle of said selected segment, and wherein the step of calculating includes the step of calculating said take-off angle.

9. The method of claim 8 wherein the step of calculating said take-off angle comprises the step of:calculating said take-off angle of said selected segment wherein said 3-D vessel hierarchy is reconstructed using data generated by an imaging system having two degrees of freedom.

10. The method of claim 8 wherein the step of calculating said take-off angle comprises the step of:calculating said take-off angle of said selected segment wherein said 3-D vessel hierarchy is reconstructed using data generated by an imaging system having three degrees of freedom.

11. The method of claim 6wherein the step of receiving input from the user to select said desired quantitative measure includes the step of receiving user input to evaluate an insertion angle of said selected segment, and wherein the step of calculating includes the step of calculating said insertion angle.

12. The method of claim 6wherein the step of receiving input from the user to select said desired quantitative measure includes the step of receiving user input to evaluate a length of said selected segment, and wherein the step of calculating includes the step of calculating said length.

13. The method of claim 6wherein the step of receiving input from the user to select said desired quantitative measure includes the step of receiving user input to evaluate tortuosity of said selected segment, and wherein the step of calculating includes the step of calculating said tortuosity.

14. The method of claim 13 wherein the step of calculating said tortuosity comprises the step of:calculating said tortuosity in accordance with the F-S theory.

15. The method of claim 13 wherein the step of calculating said tortuosity comprises the step of:calculating said tortuosity in accordance with the Point-to-Point method.

16. The method of claim 6 wherein the step of receiving input from the user to select said desired quantitative measure includes the step of receiving user input to determine an optimal viewing angle of said selected segment relative to said desired quantitative measure.

17. The method of claim 16wherein the step of receiving input from the user to select said desired quantitative measure further includes the step of receiving user input to determine said optimal viewing angle relative to foreshortening of said selected segment, and wherein the step of calculating comprises the step of calculating foreshortening of said selected segments for all possible viewing angles, and wherein the method further comprises the steps of: determining said optimal viewing angle from said all possible viewing angles so as to minimize said foreshortening; computing an updated 2-D projection of said 3-D vessel hierarchy data structure based on said optimal viewing angle; and displaying said updated 2-D projection on the user's display screen.

18. The method of claim 17 wherein the step of calculating said foreshortening comprises the step of:calculating said foreshortening of said selected segment wherein said 3-D vessel hierarchy is reconstructed using data generated by an imaging system having two degrees of freedom.

19. The method of claim 17 wherein the step of calculating said foreshortening comprises the step of:calculating said foreshortening of said selected segment wherein said 3-D vessel hierarchy is reconstructed using data generated by an imaging system having three degrees of freedom.

20. The method of claim 17 wherein the step of determining said optimal view comprises the step of:displaying a color coded map having a plurality of points on the user's display screen wherein the color of each said point indicates the magnitude of said foreshortening at a corresponding viewing angle; receiving input from the user to select a desired to select a point of said plurality of points; and determining said optimal viewing angle as the viewing angle corresponding to the selected point.

21. The method of claim 16wherein the step of receiving input from the user to select said desired quantitative measure further includes the step of receiving user input to determine said optimal viewing angle relative to overlap of said selected segment, and wherein the step of calculating comprises the step of calculating overlap of said selected segments for all possible viewing angles, and wherein the method further comprises the steps of: determining said optimal viewing angle from said all possible viewing angles so as to minimize said overlap; computing an updated 2-D projection of said 3-D vessel hierarchy data structure based on said optimal viewing angle; and displaying said updated 2-D projection on the user's display screen.

22. The method of claim 21 wherein the step of calculating said overlap comprises the step of:calculating said overlap of said selected segment wherein said 3-D vessel hierarchy is reconstructed using data generated by an imaging system having two degrees of freedom.

23. The method of claim 21 wherein the step of calculating said overlap comprises the step of:calculating said overlap of said selected segment wherein said 3-D vessel hierarchy is reconstructed using data generated by an imaging system having three degrees of freedom.

24. The method of claim 21 wherein the step of determining said optimal view comprises the step of:displaying a color coded map having a plurality of points on the user's display screen wherein the color of each said point indicates the magnitude of said overlap at a corresponding viewing angle; receiving input from the user to select a desired to select a point of said plurality of points; and determining said optimal viewing angle as the viewing angle corresponding to the selected point.

25. The method of claim 16wherein the step of receiving input from the user to select said desired quantitative measure further includes the step of receiving user input to determine said optimal viewing angle relative to a combination of foreshortening and overlap of said selected segment, and wherein the step of calculating comprises the steps of: calculating overlap of said selected segments for all possible viewing angles; and calculating foreshortening of said selected segments for all possible viewing angles, and wherein the method further comprises the steps of: determining said optimal viewing angle from said all possible viewing angles so as to minimize the combination of said overlap and said foreshortening; computing an updated 2-D projection of said 3-D vessel hierarchy data structure based on said optimal viewing angle; and displaying said updated 2-D projection on the user's display screen.

26. The method of claim 25wherein the step of calculating said overlap comprises the step of: calculating said overlap of said selected segment wherein said 3-D vessel hierarchy is reconstructed using data generated by an imaging system having two degrees of freedom, and wherein the step of calculating said foreshortening comprises the step of: calculating said foreshortening of said selected segment wherein said 3-D vessel hierarchy is reconstructed using data generated by an imaging system having two degrees of freedom.

27. The method of claim 25wherein the step of calculating said overlap comprises the step of: calculating said overlap of said selected segment wherein said 3-D vessel hierarchy is reconstructed using data generated by an imaging system having three degrees of freedom, and wherein the step of calculating said foreshortening comprises the step of: calculating said foreshortening of said selected segment wherein said 3-D vessel hierarchy is reconstructed using data generated by an imaging system having three degrees of freedom.

28. The method of claim 25 wherein the step of determining said optimal view comprises the step of:displaying a color coded map having a plurality of points on the user's display screen wherein the color of each said point indicates the magnitude of said overlap and the magnitude of said foreshortening at a corresponding viewing angle; receiving input from the user to select a desired to select a point of said plurality of points; and determining said optimal viewing angle as the viewing angle corresponding to the selected point.