VIEWS: 4 PAGES: 6 POSTED ON: 3/3/2011
AN ONLINE MEASURING SYSTEM TO SUPPORT HEART VALVE SURGERY N.Rajanipriya, Deepa.M B.V.Raju Institute of Technology, B.V.Raju Institute of Technology, Narsapur,Medak(dist). Narsapur,Medak(dist). Email : priya_04_n@yahoo.co.in Email: deepa_madathil@yahoo.co.in ABSTRACT of mitral valve for a suitable valve replacement. It In this paper we present an online system was also intended to provide a friendly solution to be designed to assist heart surgeries, in particular to used by surgeons in real time. provide support to make the decisions during implant Segmentation was done using active heart valve operations. The 2-D echo image is contours and dimensions are measured from Global obtained from the output console of the ultrasound Image Binarization. scanner. The heart valves in the image are We have also developed a Graphical User segregated by processing the image using active Interface using Visual Basic to record the suggestions contours. The properties of that valve: the given by the doctor for artificial heart valve dimensional properties such as number of cusps, replacement. area of valve, radius of each cusp etc., and material properties such as the elasticity, viscosity and BACKGROUND chemical interactions with in the body are obtained. The heart is the center of the cardiovascular The computer assisted system consists of the system. Whereas the term „cardio‟ refers to the heart, database which has the dimensional and the material the term vascular refers to blood vessels. The heart properties of the artificial heart valves. When the propels blood through thousands of miles of blood properties of natural heart valve are given as an vessels, and it is magnificently designed for this task. input to the computer assisted system, it compares The interior of the heart is divided into four them with the properties of artificial heart valves in compartments called chambers that receive the the database using the pattern recognition principles. circulating blood. The two superior chambers are The matching of the properties in the system called the right atrium and left atrium. The two determines the selection of artificial heart valve with inferior chambers are the right ventricle and left which the valve of the patient has to be replaced. This ventricle. As each chamber of the heart contracts, it therefore forms the decision making system assisting pushes a portion of blood into a ventricle or out of the the doctor to identify the appropriate artificial heart heart through an artery. To prevent back flow of valve for replacement. blood, the heart has valves. These structures are composed of dense connective tissue covered by INTRODUCTION endocardium. Valves open and close in response to pressure changes as the heart contracts and relaxes. In recent years mostly since the end of the Atrioventricular valves (AV) lie between the atria last decade growing investment in the development and ventricles. The right AV valve between the right and improvement of computer vision techniques has atrium and right ventricle is also called the tricuspid been a common policy of companies and countries. valve because it consists of three cusps (flaps). The Such techniques rely on the unstoppable advance of left AV between the left atrium and left ventricle has technology which provides cheaper and faster two cusps and is called the bicuspid (mitral) valve. hardware support on a regular basis. Both arteries that emerge from the heart valve that Computer vision aims to increase the speed prevents blood from flowing backward into the heart. and quality of some problem solutions, presenting as These are the Semi lunar (SL) valves. The pulmonary its ultimate goal, the creation of reliable systems semi lunar valve lies in the opening where the capable of taking the right automatic decisions pulmonary trunk leaves the right ventricle. The aortic supported by image analysis thus achieving results valve is situated at the opening between the left that until then were only made possible by human ventricle and the aorta. intervention. Medicine one of the most ancient and Heart valve prostheses are intended for use important fields of science is surely no exception and as the replacement of diseased natural valves. The the main problem that led to our work was related in four natural valves (tricuspid, mitral, pulmonic and this area in particular with heart surgery. aortic) become diseased, due, in part, to rheumatic The problem we tackled in this work was to fever, calcification, endocarditis, and congenital automatically segment and measure the dimensions defects. This results in the restriction of the forward flow of blood known as stenosis or regurgitation flow snaxels, connected by line segments. The position of of blood known as insufficiency. The advent in the snake is described by an M 2 matrix X, where the science and technology has led to new and improved ith row contains coordinates of the ith snaxel designs of these heart valves. Heart valve prostheses x0 y0 may be classified into 2 types, mechanical valves and x1 y1 tissue valves. Mechanical valves are classified into 3 X= ……………………… (1) categories, caged ball valves, tilting disk valves and .... .... bileaflet valves. Tissue valves use materials of xM 1 yM 1 biological origin as the valving elements viz, allograft valves, aortic valve xenograft and The contour motion is a result of interplay pericardial xenigraft. between its mechanical properties (deﬁned as its SEGMENTATION AND TRACKING IN inertia, internal dissipation, and elastic stiffness) and ECHOCARDIOGRAPHIC SEQUENCES: its potential energy (which is inversely proportional ACTIVE CONTOURS GUIDED BY OPTICAL to the image gradient). The equation of motion for FLOW ESTIMATES this model is Tracking algorithm that uses a Kalman ﬁlter MX(t)+CX(t)+KX(t)=F(t) …………….. (2) to track the object and estimate true motion has been proposed. The observations for this Kalman ﬁlter are where F is the matrix of x and y components of optical ﬂow estimates obtained by application of the image forces at each vertex, where M = μ.Im, where μ block-wise Horn and Schunk algorithm. With such is the mass assigned to each vertex, C = γ.Im, where γ observations, the tracker should be able to more is the constant damping density and K is an MxM successfully handle sudden movements of the object stiffness matrix. Dots above the X represent the of interest. However, the use of optical ﬂow as the derivatives in time. only input to the algorithm means that the current Two initial conditions are needed for solving position and shape of the object are estimated by the second-order differential equation. If the initial integrating object motion. This results in position and velocity are known, (2) has a unique accumulation of errors in the optical ﬂow estimates solution. This suggests the way to incorporate the over time. This model does not provide a way of optical ﬂow information into the active contour incorporating the current image information into the model. For the ﬁrst frame of the image sequence, a current estimate of object shape and position. user deﬁnes the initial contour, and the initial velocity In the algorithm, the tracking problem is vector, is set to zero. Once the contour settles in the solved by using optical ﬂow estimates as initial ﬁrst frame, the subsequent initializations are done conditions when integrating the equation of motion. automatically. The initial position is given by the In other words, optical ﬂow is used to push the ﬁnal position from the preceding frame, and the contour toward the new position of the object. initial velocities of vertices are given by the optical In the simpler algorithm proposed by Ayache, an ﬂow estimates. energy measure that tends to preserve the matching Optical ﬂow estimates give the average of high curvature points and to enforce a smooth ﬁeld velocity of the object between two frames, not the of displacements vectors between contours in two initial velocity. They can, however, be used as initial consecutive frames is minimized (such an approach velocities for two reasons. First, we are using a small also requires previous edge detection). In the method value for the damping factor, and second, since the we propose, instead of matching pieces of contours ﬁnal contour position is adjusted through the action around high curvature points, through the use of of image forces, we only need rough estimates of optical ﬂow, the best matches for image regions initial velocity. around each vertex of the contour are found in the To estimate optical ﬂow, we used a next frame. The beneﬁts of the latter approach are multiscale version of Singh‟s algorithm, which is notable in the high deformation frames (our capable of detecting large displacements. algorithm required user intervention in only one frame in three sequences showing the mitral valve. METHOD A. Integration of the Equation of Motion The contour is deﬁned as a set of M vertices, 1) Single-Step Time Marching Scheme: The equation of the dynamic equilibrium of the contour (2) can be integrated using a single-step time- M i 0 2 dFI dXi ……………………… (6) marching scheme (SSTMS). In the SSTMS, the time t is divided into M i 0 1 dFI dXI 0 several intervals, or time steps t i . In . where dF M C X kX F is the force residue each time step, the position and the and dX is the position increment. If is less than velocity of the object are updated. 5%, it is assumed that the contour has reached the Calculations performed for each time equilibrium. interval are i) X.. 2 X.n t 2 n 3) Calculation of the Contour Position: The time interval between two successive frames is assumed to X.. 2 X.n n be one. This time interval is then divided into Nt intervals for which STMS integration is executed. To avoid attraction of the contour by noise and undesired t2 objects that can be in its way while moving toward ii) M t C1 n 2 K 2 the new position of the object of interest, we chose a nonuniform time division. The ﬁrst interval t0 is set . F C X n 2 k X n 1 .. at 0.9 and the remaining 0.1 is divided into Nt-1 equal intervals t n . During t0 , image forces are t2 set to zero to prevent attraction of the contour by iii) X n 2 X n t X.n n 2 undesired valleys of the potential surface. In other X.n 2 X.n n t ………… (3) words, we are making the potential surface ﬂat for the ﬁrst 90% of the time. During this time, elastic where n is the time step number, X 0 X(0) and forces act on the contour, ﬁltering inaccuracies in X. X . 0 and 1 and 2 are constants. The 0 velocity estimates. During the remaining 10% of the unconditional stability of the integration is time, image forces are included, making the contour guaranteed if 1 0.5 and 2 1 . Steps i)–iii) are settle on the edges of the object of interest. The contour is resampled after each iteration to repeated for each of the time intervals. The problem compensate for possible clustering of the snaxels in is nonlinear if the stiffness matrix K is not constant the local energy minima along the valley of the and depends on the position of the contour. In such a potential surface. After resampling, distances case, the equation of motion takes a more general between snaxels are equal. Since both positions and form velocities of points from one iteration are used in the M X C X PX F ……………… (4) .. . following one, velocities of new snaxels have to be recalculated as well. The position of a new snaxel is a In this case, the problem is linearized inside a time linear combination of positions of two old snaxels interval and, hence, its velocity is also a linear combination of _ velocities of the same two snaxels. PX P X n K X ……………………. (5) where K is the value of K at some average value _ B. Contour Stiffness and Internal Forces In the case of the discrete contour, internal forces of X inside the current interval. (4) has to be solved depend on positions of vertices. Internal forces can iteratively using steps ii) and iii) of the procedure usually be decomposed into the product of the given in (3). It will be shown that the stiffness matrix stiffness matrix and the position matrix Fint=-K(X)X in our model depends on the position of the contour. Where such decomposition is not possible , the However, the internal forces can be decomposed as equation of motion takes on the form of (4). P(X)=K(X)X,with stiffness matrix K recalculated at the beginning of each iteration. Depending on the shape of the object of interest, different internal forces can be designed. The role of 2) Equilibrium Criterion: The iteration process is such forces is to enforce on the contour a shape terminated when the contour attains the equilibrium feature that is characteristic of the object. For at each time step. In the proposed model, the ratio of example, in the case of heart valve leaﬂets, assigning contour deformation energy in the current and the the contour elastic properties that try to preserve ﬁrst iteration was chosen as an equilibrium measure length will result in failure of the active contour to the developed approach is quite ﬂexible. A priori detect sudden length changes. However, the leaﬂet knowledge about the shape characteristics of the shape feature that seems to be preserved throughout object of interest can be incorporated into the the image sequence is leaﬂet thickness. This fact deﬁnition of internal forces, and the rest of the suggests that in the case of valve leaﬂets, additional algorithm remains unchanged. forces that act to preserve contour thickness should be added. To deﬁne such forces, for each point i on the contour, the distance di from the matching point on the opposite side is calculated (Fig. 1). Also, the i distance d i0 calculated at the optimal position in the ﬁrst frame is stored and considered to be optimal. The internal force that acts on the vertex i is then deﬁned as a force proportional to the relative change di in length d0 di di di 0 1 1 j f int k k d i 0 d 0 di di di i 1 1 = k 0 X J X i d i di Fig.1. Internal force acts to preserve the thickness of 1 1 the leaﬂet.The force at vertex i will be in the direction d = k 0 x j x i , y j y i ……. (7) of di and its magnitude (positive or negative) will i di depend on the degree of change in the length of d i. where k is the constant elastic coefficient. The form C. Image Forces of (7) indicates that in the ith row of the stiffness To calculate image forces, we first need to matrix, all elements are zero, except Kii and Kij define the potential surface. It is computed as Hi, j G * I (i, j) 2 If we deﬁne ……. (8) 1 1 ki k 0 , then Kii=ki and Kij=-ki d where I is the image intensity, G is a two- i di dimensional (2-D) Gaussian mask with standard Introduction of such a force prevents the opposite deviation and * is a 2-D convolution. The sides of the contour from collapsing. After t 0 the standard deviation of the Gaussian mask has to be contour is close to the leaﬂet, but has not reached it chosen carefully. When the object of interest is very completely. At this time, image forces are introduced, small in at least one dimension (such as in heart valve and very often both sides of the contour can be leaﬂets, which are very thin), too large value for attracted by the same side of the leaﬂet, i.e., be on the may completely blur the edges of the object. We used slopes of the same valley of the potential surface. the value 1.2 for the mitral valve, and 1.8–1.9 for Once opposite sides of the contour begin getting other image sequences of both a healthy and diseased closer, the elastic force pushes them apart and brings heart. These numbers were chosen experimentally. one side of the contour under the inﬂuence of the Components of the image force at each vertex are opposite edge of the leaﬂet. calculated as the negative of the partial derivatives of H in corresponding directions For the left ventricle and aortic root image sequences used to test our algorithm, the elastic coefficient k H xi , yi fi x g ; was set to zero. The velocity estimates seemed to be x norm H xi , yi quite accurate and the problem of snaxel clustering was solved by contour resampling in each iteration. fi y g …….(9) y norm The example of tracking the leaﬂet motion shows that where index norm denotes that the variable is third from 150 to 255. Through histogram analysis normalized to values from zero to one. The Sobel the system will then find out which of these classes operator was used to calculate both H and F contain the largest set of points belonging to the .Parameter g is used to control the amplitude of image. Here, three distinct situations can arise. In the image forces. first situation, when the largest set of pixels belongs to the first defined class [0 - 49], the image consists predominantly of dark tones, thus making it harder to find a good threshold value. For these cases, we developed a binarization algorithm that proved to be somehow very efficient. In this simple method, the threshold value is determined by averaging the histogram maximum gray level (level containing the most image points) with the medium histogram gray level (sum of every pixel‟s level divided by the total number of image pixels). HIST(1) HIST(2)......... .. HIST( N) MAX(HIST) Threshold N 2 Fig.2. 2-D echo images of Mitral valve From the formula above, it can be easily understood MEASURING SYSTEM that for very dark pictures, where other algorithms return a low threshold value that results in poor Automatic Global Binarization quality binarizations and subsequent loss of some calibration target points (due to noise), this method As mentioned before, a reliable image settles a higher threshold value located between the binarization method was developed in order to maximum and the average histogram levels. provide good quality binarized image (in real - time), On the other hand, if the largest set of pixels which represent the most important requirement of belongs to the second defined class [50 - 149], then the template matching algorithm used in the system the image is mostly composed by intermediate gray to detect the calibration target set of points. levels and therefore the best method to apply is the To achieve this goal, some known well known Otsu method (which works especially binarization methods were tested. However, as well if there are no major image brightness expected, no method was found to be general enough variations). – the pictures analyzed, mostly heart tissues, had Finally we reach the last case, where the indeed very specific characteristics (in terms of largest set of pixels belongs to the third defined class, brightness, contrast, and color balance) that changed ranging from 150 to 255. Here, the images are in fact greatly from picture to picture, leading to a too light and both methods described above reveal miscellany of good and bad results for each method their lack of efficiency regarding to finding a applied. threshold value that maximizes binarization quality. To cope with this problem, the solution In this case it is used an algorithm introduced in found was to introduce a decision procedure based on O‟Gorman, which is also a global approach, but uses the image histogram, capable of choosing the most a measure of local information, namely connectivity. appropriate binarization method according to image This method maximizes connectivity of the resulting properties. It should be noticed that histograms, thresholded regions and works exceptionally well if although very simple concepts, encapsulate more than 50% of image gray levels belong to the information that turns out to be very helpful in the third defined class. study of certain cases. At this point, after choosing the binarization As the basis for this binarization decision method according to the rules presented above and procedure was precisely the histogram, image gray performing that operation, some cleaning was found levels are grouped into three different classes, necessary in order to remove isolated pixels. This accordingly to its brightness, the first one ranging was achieved through a morphological filter. from 0 to 49, the second one from 50 to 149 and the RESULTS algorithm requires the user‟s interaction only in the first frame of the sequence. Optical flow information is used to estimate the initial position of contour in the subsequent frames. We have also presented a software package to perform online measuring of mitral valve structures to assist heart surgery. The developed application was designed to capture and process mitral valve pictures. Fig. a Fig. b FUTURE WORK We intend to implement the same concept of active contours and generation of database on the Aortic valve. REFERENCES [1]I.Miki´c, S.Krucinski, and J.D.Thomas, Fig. c Fig. d “Segmentation and tracking of mitral valve leaﬂets in echocardiographic sequences: Active contours guided by optical ﬂow estimates,” SPIE Med. Imag., vol. 2710, pp. 311–320, 1996. [2]T.McInerney and D. Terzopoulos, “Deformable models in medical image analysis: A survey,” Med. Image Anal., vol. 1, no. 2, pp. 91–108, 1996. [3]M.Kass, A.Witkin, and D.Terzopoulos, “Snakes: Active contour models,” in Proc. 1st Int. Conf. Computer Vision, 1987, pp. 259–268. [4]L.D.Cohen, “On active contours models and Fig. e balloons,” CVGIP: Image Understanding, vol. 53, no.2,pp.211–218,1991. [5]L.D.Cohen and I.Cohen,“A ﬁnite element method applied to new active contour models and 3-D reconstruction from cross sections,” in Proc. 3rd Int. Conf. Computer Vision,1990, pp. 587–591. [6]L. D. Cohen and I. Cohen, “Finite-element methods for active contour models and balloons for 2-D and 3-D images,” IEEE Trans. Pattern Anal. Machine Intell. vol. 15, Nov. 1993. [7]A.Amini,T.E.Weymouth, and R.C.Jain, “Using dynamic programming for solving variational problems in vision,” IEEE Trans. Pattern Anal. Fig. f Machine Intell. vol. 12, pp. 855–867, Sept,1990. [8]D. Geiger, A. Gupta, L. A. Costa, and J. Vlontzos, “Dynamic programming for detecting, Fig. a. Initial Contour. tracking, and matching deformable contours,” b. Contour moved to one step. IEEE Trans. Pattern Anal. Machine Intell. vol. 17, c. Contour moved to two steps. pp. 294–302, Mar. 1995. d. Contour moved to three steps. e. Tracking of mitral valve by active contour. f. Patient Case History Form. CONCLUSION We have presented an active contour model for tracking objects in image sequences. The model successfully handles large frame to frame displacements of the object of interest. The presented