WO2001043047A2

WO2001043047A2 - Dynamic computing imagery for visceral osteopathy and for articular kinetics

Info

Publication number: WO2001043047A2
Application number: PCT/BE2000/000145
Authority: WO
Inventors: Christian Williame; Georges Finet
Original assignee: Christian Williame; Georges Finet
Priority date: 1999-12-10
Filing date: 2000-12-08
Publication date: 2001-06-14
Also published as: AU1845501A; WO2001043047A3

Abstract

The present invention is related to a computerized method for studying alterations of movements of objects, physiological movements, i.e. visceral dynamics, of abdominal and uro-genital viscera or parts of viscera, said alterations being organized and repeated regularly and having a rhythm determined by respiration, said method comprising the steps of: capturing two images of said viscera in maximum inhalation condition and maximum exhalation condition respectively, improving the quality of said images by using image enhancing techniques, selecting said viscera on a screen, defining the contours of said viscera on both said images, calculating the amplitude of the viscera movement between said inhalation and said exhalation position.

Description

DYNAMIC COMPUTING IMAGERY FOR VISCERAL OSTEOPATHY AD FOR

ARTICULAR KINETICS

Field of the invention

[0001] The present invention is related to a method for visually controlling visceral movement or articular kinetics based on medical imaging techniques. The invention is equally related to a computer program for performing said method.

Aim of the invention

[0002] The aim is to propose a method and computer program to allow an automatic an repeatable definition of the contour of an object, essentially of an organ, and an automatic calculation of the alteration m the condition of said object, between two or more different images of said obj ect .

Summary of the invention

[0003] The starting point of this method and computer software is one or several images m levels of grey (BMP format) containing the viscera m inhalation and exhalation positions. These images have been scanned which may alter them. These images are taken on a film or are directly captured from X-rays or echograms systems. [0004] In order to extract the viscera, the quality of the image obtained has to be improved by increasing the contrast, enhancing the image dynamics and if necessary by smoothing the noise. Once these pre-treatment algorithms have been executed, two other algorithms are run m order to extract the viscera. The first one allows to extract the area of the organ by replenishment of the useful surface (replenishment) . The second one (DDCM) determines the edge by using a technique of active contour. The extraction of the viscera is followed by the determination of the features which provide information about the evolution of the organ between the exhalation and the inhalation.

Vertical and horizontal shifts, variation of the inclination, surface, perimeter, gravity centres, etc. are studied m this respect .

[0005] Both methods are started the same way by the user selecting the viscera on the image containing the viscera to be extracted, i.e. by clicking with a mouse inside the viscera or by defining a number of points on the contour of said viscera.

Replenishment

[0006] Starting form the click point, points that have the same ordmates are browsed and the colour of every pixel is changed until the edge of the organ is reached. After having determined and coloured this segment, the following segments are treated one by one until the whole organ is replenished.

DDCM (Dynamic Discrete Contour Model)

[0007] Starting from an initial shape, which is created with a minimum of user interaction, the dynamic contour model modifies its shape to evolve towards the desired contour. The driving force leading to the shape deformation is calculated from a so called internal force that is derived from the contour model itself and from an external force that is derived from the main characteristics of the image

[0008] In order to determine the features, aside from the organ as a whole, parts of viscera are compared m order to have a more precise idea of the movement. A stomach for example is first compared on the first third and then on the two remaining thirds. It should be noted that the comparisons of the gravity centres which give horizontal and vertical displacement take into account the patient movement correction. Results are indeed, if possible, based on a stationary reference of the patient, determined m advance such as a bone. The shape of the bone is determined by using the same algorithms as for a viscera with a contrast adaptation because bones generate more noise and are less contrasted.

[0009] This reference allows to exactly superpose the organs regardless of the patient's moving or of the image acquisition technique. [0010] When the complete contour is defined by either of the two methods, the mam features of the viscera are determined, such as the perimeter, the area, and the centre of gravity as well as the orientation of the viscera. The positions (inhalation and exhalation) are then compared m terms of coordinates of center of gravity in a given coordinate system and inclination of a given axis, preferably the principal axis of the contour. The movement thus derived is subsequently compared to prescribed physiological movements. [0011] The method of the invention is not restricted to physiological movements. Any object having a contour represented on an image can be subjected to the method of the invention.

[0012] The invention is further related to a computer program for performing the method of the invention. This computer program is a user-friendly, interactive, visual and modular image treatment software which runs on Windows on a PC-type computer. Such basic functions as contrast manipulation, filtering, inter- images operations as well as more complex functions like contour extraction, images superposition, reference search (and specifically bones used as reference considering their stability compared to viscera) are available. All these functions are available from unrolling menus and dialog boxes.

[0013] This software was developed with Visual C++ m a Microsoft environment. This software is multi- documents and multi-views: several images can simultaneously be loaded and displayed. The file format

[0014] Similar search carried out on bones articulations lead to a second computer program which allows the analysis of articular kinematics using medical imaging and will allow measurements of articular mobility and restrictions of this mobility.

Detailed description of the invention

rrect analysis of the viscera ynamics Introduction

[0015] The main difficulty m the treatment of a vision application is the lack of general theory or universal method offering the answer to several problems at the same time. For example, m order to examine the surface quality of continuous plane, products on a high-speed production line [1,2] cannot be used for recognition of characters or for medical diagnosis. Thus, every particular application requires a proper survey.

[0016] We easily admit that the relative volume of information for an application can considerably vary depending on its nature. The treatment will be also quite different according to the requirements of the application (speed, precision...) . In our application, the form of the viscera m normal exhalation and inhalation need to be recognised. It will allow to reduce the visceral dynamic m comparison with the forms.

[0017] Many applications m the field of medical imagery use the concept of contour. Intuitively, a contour can be defined as the ensemble of points that is the place of a sharp transition of a specific feature (luminance, texture,...) m a given direction. [0018] Object contour definition can be accomplished by several ways:

Completely manual : a user draws the contour using some pointing device like a mouse or a graphic tablet. This requires a certain knowledge from the operator about the clinical problem and is expected to consume much time for a significant result. This method is not desirable m the processing of large data set and its reproducibility is reduced. Completely automatic : the available techniques presented for automatic contouring are not yet enough sophisticated for many typical applications. These methods are reserved for the simplest cases of contour detection.

Automatic (followed by manual ) .- the first step (automatic) is usually based on simple techniques of contour detection (threshold or region growing) and the second step (manual) modifies the generated contour on the basis of anatomical and pathological knowledge. This approach is very time consuming and the results are not reproducible.

Manual (Automatic ) -. the manual intervention just serves to establish an initial contour which will define the initial conditions for the automatic phase; the latter procedure will refine the first contour. In this case, the influence of the users is not direct since different starting contours, determined by several users, will lead to the same final result. This method leads to a high degree of reproducibility .

[0019] In this article, we have developed two contour algorithms. They are all based on this last approach.

Osteopathic Study

[0020] Osteopaths sustain that the visceral dynamics can be origin of certain lesion if the visceral movement is not normally accomplished.

[0021] Before 1985, neither scientific study nor viscera (internal organs) displacement statistic m the abdomen under the influence of diaphragm pressure is available.

[0022] An estimate of this movement based on anatomy and physiology argument has been proposed, and some observations were made without statistical study. From an epistemological view, Georges FINET and Christian WILLIAME stress that this study has been carried out quite rigorously, thus neglecting the influence of all preconceived ideas. They ust wanted to reflect as soon as possible reality.

[0023] They started their work m 1985 m studying, with X-ray and echograms the movement m abdomen, stomach, small and large intestines, hepar, lien, reni and pancreas that are induced by diaphragma respiratory movement; note that this visceral movement does not correspond to peristaltism. They defined an optimal study of the organs on an image, established a strict protocol and methodology and took all the necessary and sufficient precautions that are required to achieve a trustworthy scientific work. After they analysed around 3000 images by a fully manual procedure and published them results m [3] unfortunately, this approach is too much time consuming for application daily life. Automatic method is required.

Methodology

[0024] All the images refer to shifts occurring from exnalation to inhalation, with the patient standing upright to reproduce the conditions of everyday life. [0025] Radiographic and echographic examinations are recorded on videotapes. The tapes are visualised on the screen of the echograph equipped with a camera that allows us to capture the pictures (the first during the inhalation, the other during the exhalation) . These pictures are computerised; the computer memorises the horizontal and vertical shifts and the variation of the inclination axis m each space plan.

[0026] The goal is to elaborate a general software for image processing that would allow us to extract the viscera on which we want to apply the treatment . This treatment aim at determining the features that would best provide the evolution of the movement of the organ between the two phases .

[0027] The first stage is related to the detection of edges and to the extraction of the viscera after a possible pre-treatment on the image computerised m 256 grey level. To do so, we developed two algorithms. The first one allows us to extract the area of the organ by replenishment of the useful surface. The second method determines the edge m using technique of active contour. The extraction of the viscera is followed by the determination of its features at the surface, along the perimeter, at the centre of gravity as well as along its principal axes.

Algorithms

[0028] The first X-rays image reports the viscera in normal inhalation and normal exhalation positions. To be processed by the computer, the image has to be described as a matrix or some other discrete data structure. Therefore, X-rays are scanned m order to get the numeric image (at the present time, we do not have the possibility to directly capture the image from X-ray or echograms system) . These images, m levels of grey, represent the starting point of the software that we have developed.

[0029] The basis and the outlines of our process are represented m Fig .1. As the original image is scanned, we loose the "real grey level" information, i.e. the contrast between the different components constituting the image. In order to extract the viscera, we would like first to enhance the contrast of the image. For example, overexposed, underexposed and/or blurred images can be improved with contrast enhancement techniques. [0030] The grey level histogram can be used for image enhancement, since it has provided several information on the overall brightness and contrast of an image and is a valuable tool for image processing at both the qualitative and quantitative levels. [0031] Histogram equalisation can also be used to produce a coarser quantification of an image. We can smooth the noise (reduce various spurious effects of a local nature m the image) or remove it. [0032] Edges detection can be performed using convolution or gradient operator, followed by a threshold operation on the gradient m order to decide whether an edge has been found. The pixels that have been identified, as edges must then be linked to form closed curves surrounding regions. Edges m picture are just a preliminary step m the extraction of the features. In fact it does not provide the separated contour of interest for the viscera.

[0033] To extract the viscera of interests, we apply two different methods that can provide complementary information. These methods allow us to get either the full organ or its edge. The following stage is to determine features that would inform us on the evolution of the organ between the normal inhalation and normal exhalation positions. In that purpose, we study the vertical, horizontal shifts and the variation of inclination axis.

1. Replenishment

[0034] To explain this method, we will refer to the

Fig .2. P₀ which represents the starting point of the method. This point corresponds to the "click" of the device done by an user. We will demonstrate that this does not affect the nature of the final edge. We formulate the hypothesis that P₀ is located inside the viscera that we want to extract . [0035] Starting from the point P₀ with coordinates Xo Yo, we browse points m the same way and change the colour of every pixel until we meet pixels representing the edge of the organ (x₀d) ■ In order to minimise the risk of considering noise as pixels of edge (residues of the pre- treatment), the algorithm introduces some criteria to assess the reliability of the viscera edge. After having determined and coloured the segment limited by the right (x_0d) ^and the left pixels (x_0g) of the contour, we move to an another line of the image. The transition to the following line is made with a step named K. This means, that the following line is distant from the first by K pixels m the vertical of the image. In our case, K is a parameter allowing the adjustment of the treatment required and its accuracy. Speed treatment governs the choice of K. Indeed, small K value increases the computing time. To fill the organ with a maximum of precision, we choose K=l . This condition makes us browsing all lines from y₀. The choice of the following point P_l starting point for the second segment, is chosen such that its abscissa is m the middle of the previous segment (x₀g/ Xod) ■ This means Ε> IS calculated m the following way:

(1) Where K represents the displacement step m pixels.

[0036] The choice of this criterion (1) does not depend on the abscissa of P₀ and allows P and the final contour to be completely independent of the choice of P₀. The right and left parts of contour between y₀ and yi made while joining on one hand x_od to x_id and on the other hand

X₀g t O X_lg .

[0037] The replenishment takes place m a very similar way for the following lines while incrementing Y.

[0038] Once the lower side has been detected, we perform the same treatment for decreasing Y (a negative step K) from P₀. The final edge is obtained while joining all points of two vectors. The right vector containing the points of the right edge sorts out according to the increasing order of Y_x . The left vector of the left edge sorts out according to the same order.

[0039] If the step K is equal to 1, we will then get the contour and the area of the organ. [0040] When we extract the edge and/or the surface of the organ, we can apply a post treatment on the picture in order to smooth the contour of the organ. Algorithms that we used in post treatment are for example the dilation and the erosion. After the post treatment, we examine data to determine the main features of the viscera such as the perimeter, the area, and the centre of gravity as well as the orientation of the organ.

2. Dynamic Discrete Contour Model (DDCM)

[0041] Terzopoulos [4] who introduced the model of the 'Snakes', constructs a contour with connected spline segments and optimizes the contour approximate to a desired form by minimizing an energy function containing internal and external energy. The internal energy corresponds to the energy of curvature of the spline while the external energy is calculated by integrating image features. Additional constraints are also introduced to obtain a better final contour . [0042] Another model proposed by Miller [5] the 'Geometrically Deformed Model' or ^λ GDM ' describes a contour as a set of verticles connected by edges. The energy function is quite different from that used m the Snake model. The energy term depends on an estimation of local curvature and the distance between a vertex and its neighbours; we add an additional term originating from the threshold values for pixels of the picture. This energy function is evaluated only for the vertex positions and not for the trajectory of the connecting edge segments. [0043] Our approach follows Miller's method by adopting the basic structure of the model, i.e., vertices connected by edges. We will describe our dynamic process m term of force, more adequate than the energy of contribution [6]

(potential energy) since force, acceleration and speed are directly related.

[0044] Starting from an initial shape, which is created with a minimum of user interaction, the dynamic contour model modifies its shape to evolve towards the desired contour. The driving force leading to the shape deformation is calculated from internal forces, derived from the shape of the contour model itself, and from an external force field derived from the main characteristics of the image. The internal force will try to minimize local contour curvature while the external forces will try to make the model follow a valley or a ridge through the « landscape ».

[0045] An image feature may be simple, like the pixel or voxel grey value, or the magnitude of the grey value gradient but might also be quite complex, like those derived by means of differential geometry. [0046] The deformation process is organized m a number of discrete steps; after each of them, the situation is analyzed with respect to position, velocity, and acceleration for each of the vertices. In this evaluation, internal and external forces on a vertex are calculated from the position of the vertex and its neighbours. These forces result m an acceleration, which changes the velocity of the vertex. This velocity determines the displacement of the vertex during the next deformation step. After a given number of deformation steps, a final stable situation will be reached and will be characterized by an equilibrium such that the velocity and the acceleration are nulls for each vertex. When described in term of energies, this situation represents a local minimum m the energy function.

[0047] During the deformation process, two undesirable effects may occur: shrinking and clustering. The first one is due to the internal forces while the second originates from the external forces. These two perturbation effects will be eliminated by an improvement of the basic technique

(cf below) .

[0048] In comparison to other models, this method guarantees the connectivity of the contour.

[0049] The method that we describe is not strictly reserved to the medical imagery treatment but can also be implemented in other domains for contour detection.

Force and Force Field

[0050] The Fig. . represents the structure of the model : an organized contour of vertices connected by a line or a segment . The position of the vertex V_L is represented by a vector p_x and the edge between V and V₁₊₁ by a vector d_± (Cartesian coordinates) . The deformation of the model is caused by the combination of forces acting on vertices; the acceleration on the vertex V is noted by a vector a_x . We will label v_x the speed to the V vertex (not shown on the diagram) . Length d on an edge between segments represents the local resolution of the model : if it is large, the model will not be able to follow the weak variations of energy m the image feature.

[0051] In order to preserve a limited variation, the length of the edge is evaluated at regular intervals and if necessary, some vertices are withdrawn or inserted while keeping a resolution m a specified range. This is described m the section "Resampling" .

Internal forces

[0052] The internal forces are used to minimize the local curvature of the contour (Fig.5.) . The introduction of this concept allows for an effect of smoothing on the contour .

[0053] Let us introduce the concept of local curvature m our model. Strictly speaking, the local curvature is null on a right segment between two vertices but is not defined precisely (first-class discontinuity) on a vertex. The solution is to calculate the local curvature on the vertex by the difference between the two directions of the segments joining it. The segment leaving from the V₁ vertex is represented by the vector d_l ; its direction oy the unit vector d_:

c_t = d_l — d_t__χ (3)

[0054] According to the described definition, the local curvature c_t on the vertex V_± is given by (Fig.6.) :

The local curvature has a length as well as a direction and gives a measure unique and usable of the angle between the two segments. Besides, the length of the vector local curvature depends solely on this angle, and is not influenced by the length of the two segments joining the vertex.

SUBSTITUTE SHEET(RULE 25) [0055] We can define the radial vector and the tangential vector with respect to the position of the vertex. In order to define the tangential unit vector t , , we use the normalized sum of the unit vectors joining edge segments (Fig.7.)

The radial unit vector is given by a rotation of

[0056] Vectors /, and ι represent a system of local coordinates (r,t) at the V_x vertex that will be used to calculate the internal and external forces. If we use this system of coordinates to define the local curvature c_x, it has to be noted that this vector is oriented m the direction of the radial vector or m the opposite way. In other words, the vector c_x is oriented with respect to the axis of the radial vector and its length is derived from the product (c_t -r . Let us use the fact that the length of curvature vectors c_L can be either positive or negative c, = (c_t r_t ) r_t (4), to obtain (Fig.8.).

[0057] Now that we have defined the local curvature as a one-dimensional variable (m the direction of the radial vector) , we can define internal strengths that will act on vertices of the model to limit the deformation. [0058] To clearly assess the influence of the internal forces on the deformation of the model, let us consider a situation m which the external forces will be completely absent . [0059] Following shapes (5,6 and 8) we can observe that the internal forces and the curvature vector are closely bound. These are two vectors oriented m the same direction. However, it would be dangerous to define the internal forces proportional to vectors of local curvature ; indeed, we can be convinced by simply looking at the left side of Fig.9. Simple shapes, m the absence of external force, can be distorted m shapes having a local curvature minimum everywhere, e.g. a circle (or a symmetrical polygon due to discretization process) . In this case, the process does not stop and continuously displaces the vertices towards the center of the model (Fig .9a) . During this phase, the model shrinks (shrinking) while the local curvature is not affected, contrary to the role of the internal forces acting to reduce the local curvature.

[0060] In order to avoid this problem, a second constraint is added. This acts like an elastic force that preserves the distance between vertices of the neighbourhood m certain limits. This strength will act on vertices that are migrating towards each others by stopping the shrinking process. The point where the shrinking is observed depends on the weight assigned to the forces of the two internal constraints. In other words, we must balance the effects coming from the elastic internal forces and the minimisation of the local curvature. This is still not possible m the presence of the external strength, as shown below. In order to avoid the introduction of the second constraint (more artificial) , we looked for an alternative definition of the minimization forces of the local curvature.

[0061] As demonstrated before, Fig.9a, m the absence of external force, the contour will be reduced to a single point if we take internal forces proportional to the local curvature . [0062] The contour of Fig.9b is facing the same problem. If we apply a force proportional to the local curvature, the contour will become shorter (m term of length) . Besides, as is the case for Fig.9a, the local curvature is not reduced but simply displaced. Fig 9c does not encounter such a problem since this represents a typical situation for which a local curvature reduction is required.

[0063] On the figures located on the right side, one has used the same shapes but m a local set of coordinates r,t. The analysis of these vector curvatures and the need for the internal forces of these vectors to be determined when taking into account some criteπons, leads us to a simple solution. The internal forces f_{m ι} acting on the V_± vertices must have the same direction (radial) than that of the curvature vectors. This means that the internal forces are derived from the curvature vectors by simply modifying their length. Secondly, m order to ( uce the local curvature without affecting surfaces where the curvature is constant, the length of the internal force vector has to vanished where the contour reaches a constant curvature.

The two conditions can be fulfilled if we consider the sequence (c, ^■ r along the contour as a discrete scalar function, depending on the I position and if we use the product of convolution of this function with a discrete filter k_λ as the representation of the lengths sequence of the internal forces vectors f_{m ι} .

f_^ = (c r_t )® k_{t (5)}

[0064] The first condition can be filled when using r, for the direction of f_{m ι}

J f m, i = J f in i ^■ r i ( 6 )

SUBSTITUTE SHEET (RULE 25) [0065] The second condition is verified when choosing some suitable coefficients for the k_λ . A large class of filter can be used; following the experimental results obtained on test images, we have chosen to apply the constant

k_t = k ^*= 05 ( 7 )

[0066] In a near future, additional values could be used according to new experiments.

[0067] Fig.11 shows internal force vectors for the shapes shown m Fig.9. The left side describes the internal forces m coordinates r,t ;

[0068] In summary, the internal forces such defined lead us to the desired deformation of the contour m the absence of external forces. The next section will define the external forces and will present a solution to the gathering .

External forces

[0069] The role of external forces is to deform the model. They are defined m term of potential energy of distribution that represents « energy » or more precisely the force of an image feature. We will use straightforward but significant features of the images such as the pixel, the value of grey, and the value of the gradient. [0070] If we want, for example, that the model follows the path of the maximum gradient, we can express the gradient like a characteristic of the image and define an energy of distribution that would large when the given characteristic reaches high values . One would represent the path of the maximum gradient like one ridge m the energy of distribution. The normal behaviour of the distortion process is that it tends to drive the vertices to a local minimum of distribution energy. For the entire model, this means that the path follows the weak energy regions or a valley through the landscape of energies. This implies that if we want the model to follow a ridge, we must reverse the energy of distribution. Let's call it £„„the potential energy of resulting distribution. Force field driving an object m the direction of weak energy writes :

f = -V£

[0071] If we consider the case where there is no internal forces, the application of this relation ends up with a contour joining the local minima along a valley. [0072] If we apply the force fields described above, we are always confronted with the same problem : the force f_{un l} (force m Vi) will not only have a component perpendicular to the local direction (radial component) but also along of the path of the model (tangential component) . If there are no restrictions on the curvature of the model, the final situation will be a contour passing by the local minima of the external energy of contribution and the field will be located along the path of the contour. The tangential components will thus make the vertices moving along the contour and forming a Cluster.

[0073] As previously noticed, the introduction of elastic force can be the answer to the clustering. In this situation, the size of the elastic force must be tuned locally to the size of the external force/,-,,, at the V_x top. Experimentally, we didn't succeed m modulating the elastic force to both the internal _m,and external f„„ χ forces. [0074] The solution to the clustering is then obvious and corresponds to the one suggested above with the

SUBSTITUTE SHEET (RULE 25) 'shrinking', namely the displacement of vertices along the contour doesn't bring any contribution to the deformation of the model. As a result, we won't use directly f_{lm i} but we will decompose it m its radial component and its tangential component, only the first being used. If we label the radial component of f_ιm by f_{m ι} then its length is given by the product of f_{m l} and r,

J mi r \J nn X i J i

[0075] The local radial force f_{m r} and f_{m ι} provide the resulting force that acts on the vertex only by deforming the model, without any movement of vertices along the contour.

[0076] In practice, to define the characteristic of the picture, we defined a 5*5 neighbourhood centred on pixel ι(x,y) of the model.

We then apply f_{un r} = A •

y)<8> W] wi th JT

[0077] Where A is a multiplication coefficient, I(x,y) is the level of brightness of the pixel . We carried out a study m order to fix the weight of vector W m order to favour different research direction for the contour. The chosen vector is the one offering the best results. [0078] In the case of our images (cavity m white, black contour) , a high external force is associated to a pixel of the image corresponding to a cavity. An external force of low value is representative of a contour.

[0079] Regions that have the largest external forces encourage the expansion of the curve, m contrast to the zones presenting weak external forces where the expansion has to stop .

[0080 ] We have stil l to convert the external forces into displacement of vert ices we cal led s_t the vector displacement of the l pixel in the current iteration :

-K. where

where Kl , K2 , K3 , K4 are the constants that adjust the speed of expansion of the curve and have the dimensions of the pixel ; FI, F2 , F3 are constants that adjust the sensitivity to partitions and are measured in level of gray. These coefficients have been thoroughly tested to provide the following results:

Kl = 0 FI = 100

K2 = 1 F2 = 150

K3 = 2 F3 = 200

K4 = 3 FI is a gate to which the model is confronted.

Deformation, resampling, stop conditions and locking up

Distortion [0081] The total force/, acting on a vertex is a weighted sum of internal and external forces.

J f i — VV im J f im ,r, + w in J f m , ι

[0082] The weights can be fixed to default values for a given application but can also be fixed as free parameters by the user. One is accentuating the role of external forces, namely to follow the features of the image more precisely, the other accentuates the influence of internal wich is the smoothing process improvement .

Resampling [0083] During expansion, the vertices of the curve have a tendency to move away and to induce a decrease m the geometric resolution. In order to avoid this, we introduced a parameter l_des representing the length wanted for one segment between two consecutive vertices. From l_des , two other parameters are derived I_mn and/_n , which are representing the minimum distance and the authorised maximum between two neighbouring vertices.

[0084] Before every expansion, the resampling is evaluated m two steps : the first one checks along of the contour that there is a segment smaller than the minimal length l_mn . If it is the case, the segment is removed from the model and the two vertices are replaced by a single one located m the middle (see fig 9a) . The second step verifies along the contour that there is a segment larger than the maximal length /_m . When this occurs, this segment is divided m two smaller ones of equal size (to see fig 9b) .

[0085] We will use /_mι„ = V 2, l d„es and /_mma_nx = / ³/, l_des to insure that that the model doesn't have an oscillating behaviour, that is that vertices are withdrawn m a step of resampling and are reintegrated m the following step.

Detection of buckles

[0086] Much often, due to noise induced by the pre- treatment, some isolated points of the cavity correspond to maxima of the gradient, thus possessing a weak external strength that could stop the model m its expansion.

[0087] We observe experimentally that after some iterations a such point can lead to the formation of

SUBSTITUTE SHEET (RULE 25) buckl es .

[0088] To avoid this problem, it is sufficient just to verify after the resampling that tnere are no vertices with mterdistance lower than /_mιn . If one is found and fits the criteria, this means that they cross and that the model forms a buckle. As soon as it is detected, this buckle can be suppressed by eliminating the points of the model. Once the suppression of the buckle is done, the model finds a possibility to go beyond the initial barrier. The expansion of the model is therefore not too much influenced by noise points inside the cavity.

Automatic stop

[0089] In our approach, the insertion of a stopping criterion would improve significantly our model. This criterion can be related to the immobility of points at the tops during a fixed number of iterations. At the present stage, this criterion has not been inserted since one lets to the user the possibility of controlling the number of iterations.

Application

[0090] Following images represent an organ (stomach) to study m expiration and expiration. After adequate pre- treatment, we apply the different algorithm, compare results and finally, we compare same work with manual contour extraction.

[0091] In order to study the displacement of the organ, we need two fundamental notions: the gravi ty centre and the principal axis . Centre of gravity is the point of the object where the gravity force applies. It coincides with the mass if the gravity fields are uniform. The gravity centre C (C_r, C_c) is the pixel having the same amount of object m any direction. C_r is the row centre and C_c is the column centre. If the object is constituted of N pixels, C_r and C_c can be computed by:

is the row position of image pixel where

y, _{l €} l N] ^ls the column position of image pixel

[0092] The principal axis of a binary object is the inertia axis passing trough the gravity centre of this object and having a minimal total distance with all the pixels of the object. This axis divides the surface m two equal parts, and around which the surface is balanced. This axis also passes through one of the contour pixels; its orientation is along the direction of the maximum extent of the object. Sometimes for computational reasons, we use the minus axis (± perpendicular) as principal axis. It is not problematic for our application, since our interest focuses on the variation of the axis between the expiration and inspiration phases.

[0093] In addition of a global organ study, we propose to cut organ m different part (1/3, 2/3) m agreement with osteopath doctor. Let's recall that calculation on upper part (1/3) is based on the surface m pixels. Let's recall also that comparisons of gravity centre (give horizontal and vertical displacement) take m account the patient movement correction. Indeed, when it is possible, Results are based on a stationary reference of the patient i.e. bones. This reference permits us to superimpose organs precisely without manual correction due to a patient movement during X-ray process. [0094] The figures 16, 17, 18, 19 and 20 show the imagery treatment succession. We extract the viscera, we divide it in two parts (1/3, 2/3) then we compare features as showing m Table I .

[0095] After (f) , we apply a post -treatment m order to smooth boundary of stomachs, then we work separately on each of them to extract features as showing m the following figures.

[0096] In order to verify the method validity, we have proceed manually for extract the viscera contour , then we applied the two methods and compared results. [0097] Table I gives the difference between automatic and manual extraction results given by Replenishement . Principal axis is used to determine the angle oc : the inclination axis between exhalation and inhalation. Table II gives the difference between DDCM and manual method Here we use the minus axis instead of principal axis to determine the angle α.

These tables give as difference between replenishment and DDCM and also with manual method. Remark, results are described m Pixel .

Lower 2/3 Lower 2/ 3 Lower 2 '3 .

≠ CGx = -3 ≠ CGκ = - I ≠ CGx = -2

≠ CGy = 0 ≠ CGy = 3 ≠ CGy = -3 α = - 2 6 ° α 0 5 ° α = -3 5 °

Table I: Replenishment: Features comparison between exhalation and inhalation. 300T1 and 300T2: automatic detection, mainl and main2 : manual cuts.

= -36° u. = - I ° α = 46°

Lower 2/3 Lower 2/3 Lower 2'3

≠ CGx = _2 ≠ CGx = .2 ≠ CGx = 1

≠ CGy = 1 ≠ CGy = 5 ≠ CGy = -3 α = -03° α = 3 ° α = -4.5 °

Table II: DDCM : Features comparison between exhalation and inhalation. 300T1 and 300T2: automatic detection, mainl and main2 : manual cuts.

Results and difference between DDCM and Replenishment on other image .

Validation

[0098] The Fig.22. illustrating a comparison between results found automatically by our methods and those used manually by G. Finet and Ch. Williame [3] .

Φl = 1-2 : represents the variation m degree of the angle of the exhalation towards the inhalation between the result (image 1) of the original image 192-1 and the result (image 2) of 192-2.

Ωl : represents tne variation m degree of the angle of the exhalation towards the inhalation determined manually by G. Finet and Ch . Williame. Δl : gives the gap m degree between Φl and Ωl . Φ2, Ω2 and Δ2 give results for images 231-1 and 231-2.

[0099] We can note that gaps Δ are minimal between results that we found and those gotten by osteopaths. The gap Δ is included m the interval of error determined by G. Finet and Ch. Williame at the time of the statistical study of error tests and margin of error [3] . We can conclude therefore that we have obtained accurate result with an interesting consuming machine time (sec Vs Hours) . [0100] Results of this study have shown that an organised dynamics does exist on visceral level : viscera move on a specific way under diaphragmal pressure.

Furthermore, G. Finet and Ch. Williame noticed that this visceral course can be modified or disappears m that case, a new statistical study shows that stomach trouble (for example) corresponds to abnormal dynamics of the same organ.

[0101] In order to re -equilibrate those movements , this study permits to found new visceral manual normalisation; technique is simple, based on the precise knowledge of visceral dynamics discovered m the research.

Conclusions and perspectives

[0102] The study done by G. Finet and Ch. Williame have shown that an organised and repetitive dynamics does exist on visceral level. It means that viscera move on a specific way under diaphragmal pressure. For example, stomach trouble corresponds to abnormal dynamics of the same organ, this visceral course can be modified or disappears with manual test. This permit, after manual normalisation, to control the return to normal dynamics .

[0103] Owing to our algorithms of contour detection and extraction of viscera, the treatment can be made quickly and m a brief time.

Results gotten by these algorithms have been tested in relation to those founded manually m this paper and m [3] . These tests prove the validity and the efficiency of us our two methods to detect quickly and of adequate manner viscera .

[0104] At present, beginning pictures are the X-ray, what provokes a loss of information at the time of the transformation m images. Our objective is that the software runs directly with echograms and X-ray machines, what alleviates the phase of the pre- treatment and permits by the same opportunity to treat several viscera at a time. This software must be completely automatic m order to facilitate the task to users and permits to note anomalies directly at the patient. A last point is to allow the software to animate the dynamics of viscera, detect anomalies and propose solutions. References :

[1] : M. BENJELLOUN, Le parallelisme pour

1' inspection de produits plats en defilement rapide , ph.D. en Sciences Appliquees at Faculte Polytechnique de Mons, 1994. [2] : M. BENJELLOUN and G. LIBERT, Parall el

Compu ting for Defects De tection ' ' , Fourth International Workshop on Parallel Image Analysis, IWPIA' 95, LIP-ENS LYON, December 7-8, 1995.

[3] : G. FINET and Ch . WILLIAME, Bιometrιe de la dynamique viscerale et nouvelles normalisations osteopathiques' ' , Ed. R. JOLLOIS 1992. [4] : M. Kass, A. Witkm ' ' Snakes : Active contour models'' Int. J. Comput . Vision, pp. 321-331, 1988.

[5] : J. V. Miller, D.E. Breen, ''Extracting geometric models through constraint minimization, ' ' Rensselaer Polytechnic Institute, Tech. Rep. N° 90024, 1990. [6] : M. Alonso and E . J. Finn, Fundamental University Physics. Reading, MA : Addison-Wesley, 1968.

[7] : S. Lobregt and Max A. Viergever, A Discrete

Dynamic Contour Model ' ' , IEEE Transactions on Medical

Imaging, Vol. 14, No. 1, March 1995. Pp. 12-24.

[8] : F. Mauro, Detection automatique du contour ventriculaire gauche par images echographiques ' ' T.F.E.

Faculte Polytechnique de Mons, 1997.

Claims

1. A computerized method for studying alterations of physiological movements, i.e. visceral dynamics, of abdominal and uro-gemtal viscera or parts of viscera, said alterations being organized and repeated regularly and having a rhythm determined by respiration, said method comprising the steps of capturing two images of said viscera m maximum inhalation condition and maximum exhalation condition respectively, improving the quality of said images by using image enhancing techniques, selecting said viscera on a screen, - defining the contours of said viscera on both said images, calculating the amplitude of the viscera movement between said inhalation and said exhalation position.

2. Method according to claim 1, wherein said contour definition takes place through a technique of replenishment of the surface delimited by said contour, said contour definition comprising the steps of :

Receiving a starting point, selected by a user on said surface by the use of a pointer, e.g. a mouse or a graphic tablet, by clicking inside said surface onto a first point, browsing through the row of points having the same ordmate as said first point, and changing the colour of each pixel, until the two edges of said surface are reached, repeating said browsing on a plurality of ord ates, until a plurality of horizontal lines on said surface have been coloured,

3. Method according to claim 2, further comprising a post-treatment, for example based on dilation or erosion for smoothing said contour.

4. Method according to claim 1 wherein said contour definition takes place through the technique of the dynamic discrete contour model, comprising the steps of : receiving a number of points on said contour, selected by the user, - creating a dynamic contour model, comprising said points (vertices) connected by lines, calculating the internal and external forces working on said dynamic contour model, calculating the deformations of said vertices as a consequence of a linear combination of said forces, thus creating an improved dynamic contour model, repeating said calculating steps until the desired contour is reached.

5. Method according to claim 4, wherein said internal forces are derived from the shape of said dynamic contour model, m such a way that they have the direction of the curvature vectors of said model m said vertices and become zero when the curvature becomes constant , and wherein said external forces are calculated from the main characteristics of the image, in particular from the brightness of a pixel, and wherein only the radial component of said external force is taken into account.

6. A method according to claim 5, wherein before each recalculation step, a check is done whether the distance between each two consecutive vertices is lying between predefined boundaries [l_mιn, lmaxl ^and wherein two vertices are replaced by a single one m case said distance is smaller than l_mιn, and wherein a vertex is added m case said distance is larger than l_max.

7. A method according to claims 1 to 6, wherein said images are captured from X-rays or echographs systems, magnetic resonance IRM, scanner or any other medical imagery equipment.

8. A method according to claims 1 to 7, wherein said image enhancing technique may be one or a combination of the group consisting of contrast increase, enhancement of image dynamics and smoothing of noise.

9. A method according to claims 1 to 8, wherein said amplitude is calculated as the difference between X-coordmates of the center of gravity of said surface, the difference between Y-coordmates of the center of gravity of said surface, the difference between the inclination of the principal axis of said surface, said differences being taken between said inhalation and said exhalation postions.

10. A method according to claims 1 to 9 , wherein a stationary reference of the patient is taken into account, such as a bone of said patient.

11. A method according to claims 1 to 10, wherein said method is performed on images m each of the space planes of interest, e.g. three orthogonal planes.

12. A computerized method of measuring articulatory restrictions of mobility, comprising : automatically determining an articulation m medical imaging, - automatically observing and calculating an amplitude of the articulation along an axis to determine the mobility, and automatically comparing the mobility with a prescribed mobility to determine any restriction m the articulation.

13. A computer program comprising program code means for performing all of the steps of any one of the claims 1 to 11, when said program is run on a computer.

14. A computer program comprising program code means stored on a computer readable medium, for performing all of the steps of any one of the claims 1 to 11, when said program is run on a computer.

15. A computer program comprising program code means for performing all of the steps of claim 12, when said program is run on a computer.

16. A computer program comprising program code means stored on a computer readable medium for performing all of the steps claim 12, when said program is run on a computer.