In coronary angiography the coronary blood flow velocity can be assessed by tracking a contrast agent in an image sequence. The contrast agent flow velocity is used to estimate the functional behavior of the coronary arteries using TIMI frame count (TFC). In this paper we descibe the 3D reconstruction algorithm used in our method towards the automation of TFC. The method creates a two dimensional map of the contrast agent in which the opacification of the vessel centerline is plotted against time. This map is used to find the velocity of the contrast agent and subsequently the TFC. The vessel centerline is obtained using the Fast Marching Method to find the minimum cost path between\ud
the catheter point and the end of the vessel. The determination of the start and endpoint is estimated using a 3D model reconstructed from bi-plane 2D image\ud
data. The final reconstructed coronary model only includes the segments matching our error criterion

Brinke, G.A. ten

Slump, C.H.

Storm, C.J.

English

In coronary angiography the coronary blood flow velocity can be assessed by tracking a contrast agent in an image sequence. The contrast agent flow velocity is used to estimate the functional behavior of the coronary arteries using TIMI frame count (TFC). In this paper we descibe the 3D reconstruction algorithm used in our method towards the automation of TFC. The method creates a two dimensional map of the contrast agent in which the opacification of the vessel centerline is plotted against time. This map is used to find the velocity of the contrast agent and subsequently the TFC. The vessel centerline is obtained using the Fast Marching Method to find the minimum cost path between the catheter point and the end of the vessel. The determination of the start and endpoint is estimated using a 3D model reconstructed from bi-plane 2D image data. The final reconstructed coronary model only includes the segments matching our error criterion

ten Brinke, G.A.

Slump, Cornelis H.

NARCIS 

Bi-Plane X-ray coronary 3D reconstruction for flow velocity assessment.

University of Twente Research Information

BI-PLANE X-RAY CORONARY 3D RECONSTRUCTION FORFLOW VELOCITY ASSESSMENTGA ten Brinke1, CH Slump1, CJ Storm21Signals and Systems, University of Twente, Netherlands2Cardiology, Hospital Walcheren, NetherlandsAbstractIn coronary angiography the coronary blood flow ve-locity can be assessed by tracking a contrast agent inan image sequence. The contrast agent flow velocityis used to estimate the functional behavior of the coro-nary arteries using TIMI frame count (TFC). In this pa-per we descibe the 3D reconstruction algorithm usedin our method towards the automation of TFC. Themethod creates a two dimensional map of the contrastagent in which the opacification of the vessel center-line is plotted against time. This map is used to findthe velocity of the contrast agent and subsequently theTFC. The vessel centerline is obtained using the FastMarching Method to find the minimum cost path be-tween the catheter point and the end of the vessel. Thedetermination of the start and endpoint is estimatedusing a 3D model reconstructed from bi-plane 2D im-age data. The final reconstructed coronary model onlyincludes the segments matching our error criterion.1 IntroductionTo measure coronary flow velocity, we have proposedthe creation of a flowmap [1]. The flowmap is a twodimensional map of the contrast agent in which theopacification of the vessel centerline is plotted againsttime. This flowmap is used to find the velocity of thecontrast agent and subsequently the TFC. An exam-ple flowmap is shown in figure 1. Currently, the vesselcenterline is obtained using the Fast Marching Methodto find the lowest cost path between the catheter pointand the manually selected endpoint of the vessel. Fit-ting a 3D model of the coronary arteries can automatethe manual point selection. In this paper, we will focuson the creation of a 3D model of the coronary arteriesusing 2D images.2 MethodCanero et al. [2] proposed a method for 3D recon-struction using bi-plane snakes. In our algorithm weFlow map. Dataset: 7Contrast intensity along vessel centerline (s)TFC= 30.0Framenumber (n)0.2 0.4 0.6 0.8 120406080100120Figure 1: A typical flowmap.are also using a B-spline parameterization of vesselsegments, but we use a different feature map for costminimization.The geometry of the acquisition device can be de-scibed by a pinhole camera model. The transformationfrom world Xw = [xw, yw, zw] to camera coordinatesX = [x, y, z] can be described byX = [R|t]xwywzw1 (1)where [R|t] are the extrinsic camera parameters inwhich R is a rotation matrix and t a translation vec-tor.The perspective projection is described by the intrinsicparameters as matrix K:K = fκ 0 u00 fκ v00 0 1 (2)Now we can obtain a projection using the camera po-sition:Xp = K[R|t]xwywzw1 (3)in which K[R|t] is the camera matrix P .IEEE-EMBS Benelux Chapter Symposium November 9-10, 2009X = φ−1(X1,X2,P1,P2) (4)φ−1 is the retroprojection operator. This operator re-constructs a three dimensional point from two givenpoints Xn in two projection planes n = 1, 2 describedby camera matrix Pn. The method solves the triangu-lation problem using the direct linear transform. Thismethod is described in detail in Hartley and Zisserman[3].We will use this function to describe the force requiredfor the deformation of the bi-plane snake. the snakeis described by a 3D-curve ν(s). The deformation iscontrolled by an energy minimization function:E(ν) =∫ 10Eint(ν) + Eext(ν)ds (5)where Eint(ν) is the internal energy and Eext(ν) isthe external energy. The internal energy preservessmoothness and the external energy attracts the snaketo image features. Eint can be described by:Eint(ν) = α∣∣∣∣∂ν∂s∣∣∣∣2 + β ∣∣∣∣∂2ν∂s2∣∣∣∣2 (6)where α and β are constants controlling the tensionand rigidity respectively.The external energy Eext can be described by:Eext(ν) = φ−1 (q1 −∇F1(q1),q2 −∇F2(q2))−X (7)where qn are the projections of ν on view n = 1, 2 andFn is the used feature map. The numerical implemen-tation of the snake algorithm requires approximation ofthe derivatives with finite differences and converting tovector notation with νi = (xi, yi, zi) [4] :Eint(ν) = α |νi − νi−1|2 + β |νi−1 − 2νi + νi+1|2 (8)α∣∣∂νi∂s∣∣2 ≈ |νi − νi−1|2= (xi − xi−1)2+ (yi − yi−1)2+ (zi − zi−1)2(9)β∣∣∣∂2νi∂s2 ∣∣∣2 ≈ |νi−1 − 2νi + νi+1|2= (xi−1 − 2xi + xi+1)2+ (yi−1 − 2yi + yi+1)2+ (zi−1 − 2zi + zi+1)2(10)Notice that the distance between points is not presentin these formulas, because the points are redistributedso that they are equispaced in 3D (euclidian distance).The redistribution of the points is done by cubic B-spline interpolation. Discretization of the integral inequation 5 gives:E =n−1∑i=0Eint(νi) + Eext(νi) (11)Minimization of E allows us to rewrite 11 to be solvedusing dynamic programming:E(ν1, ν2, . . . , νn) = E1(ν1, ν2)+ E2(ν2, ν3) + . . .+ En−1(νn−1, νn)(12)In case of n = 5 nodes we can calculate the minimalenergy using subfunctions sk(νk+1):s1(ν2) = minν1E1(ν1, ν2)s2(ν3) = minν2(s1(ν2) + E2(ν2, ν3))s3(ν4) = minν3(s2(ν3) + E3(ν3, ν4))minν1,...,ν5E(ν1, . . . , ν5) = minν4(s3(ν4) + E4(ν4, ν5))(13)The recurrence relation is now stated as (omitting thesecond order term):sk(νk+1) =minνk{sk−1(νk) + Eext(νk) + |νk+1 − νk|2 . . .}(14)Each stage consists of 27 possibilities, the neighbors[3x3x3], to calculate. The internal energy is stored ateach stage. The indices of node position with the mini-mum energy cost are stored in the position matrix. Theminimum energy can now be found by backtracing inthe position matrix.2.1 FeaturemapThe featuremap is specifically constructed for eachvessel segment. This is done to prevent interferencefrom other vessel segments during the minization pro-cess of the 3D-snake. Since we use annotation forthe begin and endpoint of each vessel segment, wecan generate a featuremap from each segment. Thefeaturemap consists of a vesselcenterline trace usinga minimal cost path algorithm in a landscape gener-ated by the Fast Marching Method. From the map withonly the vesselcenterline, a (euclidean) distance mapis calculated. This distancemap is used as featuremapin the snake energy minimization algorithm.The vessel centerline is found using the Fast MarchingMethod (FMM) [5]. The FMM is typically used to builda map containing the travel time between the startpointIEEE-EMBS Benelux Chapter Symposium November 9-10, 2009Figure 2: Annotation of two images at the same car-diac phase. A vessel segment is highlighted: on topthe view of the featuremap calculated using FMM andminimal cost path, at the bottom, the original image.Figure 3: Reconstruction of phantom. See table 1 forthe reconstruction errors.and all other points using a speed function. The FMMsolves the Eikonal equation:‖∇T (x, y)‖F (x, y) = 1 (15)in which T (x, y) is the arrival time of the front andF (x, y) is the speed function. The speed function usedas input for the Fast Marching method is based on thecontrast enhanced image Ic(x, y):F (x, y) = (1− (Gσ ⊗ Ic))γ (16)segment µ1 σ1 µ2 σ21 8.21 ± 8.81 78.06 ± 37.272 10.98 ± 8.35 5.25 ± 3.333 39.10 ± 32.98 4.41 ± 2.134 64.21 ± 48.23 6.29 ± 6.185 8.75 ± 8.73 2.26 ± 1.106 35.94 ± 31.19 14.88 ± 7.757 3.65 ± 2.56 87.46 ± 58.228 6.89 ± 5.30 57.99 ± 33.459 2.37 ± 1.40 1.85 ± 1.0510 2.06 ± 0.96 2.57 ± 1.3611 2.94 ± 1.60 3.68 ± 2.0612 3.24 ± 1.95 3.83 ± 2.0813 7.62 ± 6.59 3.72 ± 1.84Table 1: Errors in reconstruction of the phantom. Seefigure 3 for the segment locations.Now the vessel centerline is the minimum cost pathfrom the starting point to the end point. The centerlineis defined by ν(s) = (x(s), y(s)) where x and y are thecoordinate functions and s ∈ [0, 1] is the parametric do-main which discribes the vessel from start to endpoint.2.2 CalibrationCalibration is required since the two projections are notacquired simultaneously and table motion may occurbetween successive acquisitions. In this case, we arerelying on the position of the catheter endpoint. Theposition of the catheter is reconstructed in 3D spaceusing the catheter endpoints in both projections. Theaccuracy of this reconstruction is determined by thedistance between the two lines originating from the 2Dplane to the focal spot. By minimizing the distance be-tween these lines, we obtain a better reconstruction. Inthis case, we can alter this distance by using the x andy parameter of the translation t vector in the projectionmatrix. This is equal to moving the image in-plane toboth x and y directions. The z component of t is lockedto guarantee convergence.3 ResultsCreation of a 3D model is required to measure bloodflow. Using the method described above, we are ableto create a 3D model using a sequence of images ob-tained at two different angles. The complete cardiaccycle is covered by the model. Table 1 shows theerrors of the reconstruction results from the phantomdisplayed in figure 3, and table 2 shows the errors ofthe reconstruction results from a clinical example dis-played in figure 4.The errors are expressed as the mean distance µn withstandard deviation σn between the projection (pi) ofthe reconstruction of the segment and the centerline ofIEEE-EMBS Benelux Chapter Symposium November 9-10, 2009Figure 4: Reconstruction of a clinical case. See table2 for the reconstruction errors.segment µ1 σ1 µ2 σ21 54.34 ± 29.79 3.25 ± 1.622 7.27 ± 5.13 3.89 ± 2.273 2.18 ± 0.96 4.06 ± 3.584 2.25 ± 1.13 3.14 ± 3.265 2.30 ± 1.36 3.90 ± 2.926 3.79 ± 1.58 7.31 ± 5.367 4.25 ± 3.25 16.99 ± 11.568 9.25 ± 6.28 12.56 ± 13.609 7.24 ± 7.43 3.07 ± 1.9110 3.60 ± 1.88 6.54 ± 6.6311 2.74 ± 1.46 4.99 ± 6.2112 12.00 ± 10.90 7.17 ± 8.80Table 2: Reconstruction errors of an clinical examplecase. See figure 4 for the segment locations.the vessel (qj) calculated by the FMM on plane n = 1, 2according to equation 17.E =1NN−1∑i=0M−1minj=0||pi − qj || (17)with N the number of points in the projection of thesegment and M the number of points from the calcu-lated centerline. The equation calculates the minimaldistance from each point in the p to all points in set qand calculates the mean. Using this error table we canassume that segments with a small error (µ < 5 pixels)in both planes can be used for our model.4 ConclusionThe reconstruction method described in this paper de-livers 3D models from coronary arteries. Based on theerror measure we discriminate between good and badreconstructed segments. The good segments are in-corperated in our model used for future flow measure-ments.5 References[1] ten Brinke GA, Slump CH, Storm CJ. Digital densitomet-ric determination of clinical relative coronary flow distri-butions. In Proceedings of the SPIE, volume 6143. Inter-national Society for Optical Engineering, 2006; .[2] Canero C, no FV, Mauri J, Radeva P. Predictive (un) dis-tortion model and 3-D reconstruction by biplane snakes.IEEE transactions on medical imaging 2002;21(9):1188–1201.[3] Hartley R, Zisserman A. Multiple View Geometry inComputer Vision. Cambridge University Press, 2003.[4] Amini A, Weymouth TE, Jain RC. Using dynamic pro-gramming for solving variational problems in vision. IEEEtransactions on pattern analysis and machine intelli-gence 1990;12(9):855–867.[5] Sethian JA. Evolution, implementation, and applica-tion of level set and fast marching methods for ad-vancing fronts. Journal of Computational Physics 2001;169(2):503–555.IEEE-EMBS Benelux Chapter Symposium November 9-10, 2009

Digital densitometric determination of clinical relative coronary ﬂow distributions.

Evolution, implementation, and application of level set and fast marching methods for advancing fronts.

Multiple View Geometry in Computer Vision.

Predictive (un) distortion model and 3-D reconstruction by biplane snakes.

Using dynamic programming for solving variational problems in vision. IEEE transactions on pattern analysis and machine intelligence 1990;12(9):855–867.

https://ris.utwente.nl/ws/files/5411176/g.a.tenbrinke_biomedical_imaging_paper_20091110.pdf

Bi-Plane X-ray coronary 3D reconstruction for flow velocity assessment

Abstract

Similar works

Full text

Available Versions

NARCIS

University of Twente Research Information