Electronic device for processing image-data, for simulating the behaviour of a deformable object

[0002] The invention relates to the field of the processing of digital image data from a set representative of a three-dimensional (3D) image, for simulating the deformable behavior of an object.

[0003] The invention applies more particularly, but nonexclusively, to the processing of a set of image data from a so-called medical image.

[0004] In numerous fields, it is very beneficial to be able to simulate interventions by an operator, with the aid of one or more known tools, on one or more deformable objects. Here, the term intervention is understood to mean either a manipulation, with a view for example to a displacement, or a local transformation, such as for example, in the case of a surgical intervention, incision or extraction of a part of an organ.

[0005] Simulation consists in displaying the image of an object and possibly of the region in which it customarily lies, and the representation of one tool at least whose <<virtual>> (in this document, words which appear within double-arrowhead brackets reflect the fact that a concept, in the context of this document, is designated with such words) displacement, relative to the object, is defined by a user interface of which a harness is maneuvered by an operator, with a view to simulating the handling of the said tool. In order to be able to simulate the reaction of the object on the tool, the user interface is capable of generating a force feedback, in accordance with the reactions of the tool. The term reaction force of an object is understood to mean force feedback.

[0006] In known devices, a reaction module makes it possible to determine this reaction force of the object on the basis of an estimated deformation of this object. This deformation is obtained with the aid of an internal forces module and of an image refresh module. The internal forces module is capable, on designation of a 3D. object appearing in a set of image data, of establishing a field of internal forces, which is representative of the deformation of the object, between nodes of a volume meshing dependent on a surface meshing of this object, on the basis of a deformation law and of an action defined by the user interface and representative of a maneuver of the tool.

[0007] The refresh module then makes it possible to calculate new image data of the object, in the presence of the estimated deformations supplemented with the representation of the tool. These new image data which form the new image of the object and possibly that of the region which surrounds it, are then displayed on a display device so that the operator can see in real time the result of the manipulation of the harness which simulates the action on the tool.

[0008] Such a device must allow the training of an operator or else the tailoring of new techniques of intervention on the object. In a field such as surgery, and more particularly still in the field of laparoscopic surgery, this type of device may make it possible to save human lives. To do this, it is imperative that the simulation makes it possible to reproduce the operator's gesture (or in other words his action on a tool, here virtual) as faithfully as possible. This requires real-time processing of the image data, coupled with reconstruction of the forces induced by the object in response to the deformation generated by the <<tool>>.

[0009] Now, on account of the calculation techniques used by known devices, estimation of the internal forces requires considerable calculation times which are incompatible with continuous dynamic simulation. In other words, contemporary devices do not make it possible to display, in a manner which is continuous in respect of a human eye, the entire action of a tool on a deformable object.

[0010] Moreover, no contemporary device makes it possible to simulate in real time an action such as incision, or tearing, or the removal of material from a deformable object.

[0011] The aim of the present invention is therefore to solve all or some of the aforesaid drawbacks in the field of the processing of digital image data of a 3D object.

[0012] It therefore proposes an electronic device for processing image data of the type described in the introduction, in which, on the one hand, there is provision for a <<collision>> module capable of estimating a point of intersection between a straight line embodying a displacement derived from the defined action and the surface meshing, and on the other hand, the internal forces module is devised so as to estimate the internal force exerted on each node of a first part at least of the volume meshing of the object on the basis of the displacement derived from the action, and applied to the nodes belonging to the surface mesh cell containing the point of intersection, of boundary conditions, and of node tensors and link tensors emanating respectively for each node and each link of at least the first part at least, from stiffness matrices specific to each volume mesh cell of at least the first part and dependent on the deformation law.

[0013] Of course, the first part of the volume mesh cell to which the above technique is applied, which will subsequently be referred to as <<masses/tensors>>, can be equal to the complete volume mesh cell. In the contrary case (when dealing in fact with a part of this volume mesh cell), the internal forces applied to the nodes of the part complementary to this first part (referred to for example as the second part) are determined on the basis of another technique, such as for example that of finite elements relying on precalculations which are stored so as to allow real-time calculations. Such a technique is taught in particular in the article by S. Cotin, H. Delingette, M. Bro-Nielsen and N. Ayache, <<Geometric and physical representations for a simulator of hepatic surgery>>, published in the proceedings of the conference Medicine meets with virtual reality, January 1996. In what follows, the dual technique of calculating internal forces and the deformation of the object will be referred to as a hybrid model.

[0014] This so-called masses/tensors technique used for calculating the internal forces and the deformation of the volume mesh cell of the object permits continuous simulation of an action exerted on a virtual tool at least. It is clear that the smaller the dimension of the first part of the volume mesh cell, the less will be the calculation time.

[0015] According to another characteristic of the invention, the device can comprise a meshing module allowing it to designate the 3D object(s) on which the simulation is to be performed by determination of an external envelope, then to decompose this or these envelopes into surface mesh cells, preferably of triangular form, and lastly to decompose the internal volume of each envelope into volume mesh cells on the basis of the corresponding surface meshing so as to provide the volume meshing of the associated object. It is clear that in the preferred case of a triangular surface meshing the volume mesh cells will be tetrahedral in shape. These shapes are currently preferred since they allow accurate modeling of an object of complex shape. However, of course, other types of meshing may be used.

[0016] Such external envelopes and volume mesh cells may be obtained by methods of segmentation (for example by extracting iso-surfaces) and of the Delaunay-Voronoï type respectively. All these methods are well known to the person skilled in the art.

[0017] In one embodiment of the device, its internal forces module is capable itself of calculating the stiffness matrices of each volume mesh cell, as well as the node tensors and link tensors.

[0018] These calculations are, as was stated earlier, performed on the basis of a deformation law which is preferably of volume linear elastic type. In other words, the force exerted on a node depends on the displacements respectively of this node and of the nodes to which it is connected, relative to their respective positions of equilibrium. Of course, other more complex deformation laws could be used, in particular non-linear laws.

[0019] The internal forces module could also be devised so as to determine the internal forces exerted on some at least of the nodes of the first part of the volume meshing on the basis of the deformation law and of auxiliary surface forces dependent on stored, chosen parameters of the object, such as for example the texture of the object, the presence of structures and underlying substructures, etc.

[0020] Here, the term auxiliary surface forces should be understood to mean for example surface tensions which, in certain situations such as an incision, will make it possible to amplify a visual effect at display level.

[0021] Likewise, the internal forces module may be devised so as to estimate the displacements of the nodes of the volume meshing (at least its first part) on the basis of the displacement derived from the defined action and from external forces, in particular of gravitational force type and/or forces of interaction between objects of one and the same region. This makes it possible to take into account, on the one hand, the partial sagging of an object under its own weight, and on the other hand the presence of neighboring objects and of the ties which exist with these neighboring objects.

[0022] Preferably, the estimated displacements of the nodes, other than those of the said surface mesh cell comprising the said point of intersection, are calculated on the basis of the internal forces by successive integrations with the aid of a method chosen from among at least the Euler method and the Runge-Kutta method, and more preferably still by the so-called <<order 4>> Runge-Kutta method. Of course, other methods of integration may be envisaged.

[0023] According to yet another characteristic of the invention, the internal forces module may be capable of simulating deformations not only of geometrical type, but also of incision and/or removal of material and/or tearing type.

[0024] To allow the simulation of cutting (or incision) and/or of tearing (or fracture), the internal forces module is able, after determining the estimated displacements of the nodes, to delete at least one link between neighboring nodes as a function of a first criterion, then to update the node tensors and the link tensors as a function of the deleted link(s), and lastly to recalculate the internal forces of the nodes of at least the first part of the volume meshing.

[0025] Preferably, the first criterion pertains to at least one parameter chosen from among at least one cue transmitted by the user interface and relating to the type of tool maneuvered, a volume variation of the volume mesh cell comprising the link to be deleted, and a length variation of a link of the volume mesh cell comprising the link to be deleted.

[0026] Here, the term cue is understood to mean for example an item of data specifying that the tool is maneuvered with a view to an incision or a destruction of material.

[0027] Likewise, to allow the simulation of the removal of material, the internal forces module is able, after determining the estimated displacements of the nodes, to delete a node in the event of detecting the deletion of all the links which join the said node to the neighboring nodes or as a function of the first criterion, then to update the node tensors and the link tensors as a function of the node and of the deleted links, and lastly to recalculate the internal forces of the nodes of at least the first part of the volume meshing.

[0028] The tools (here virtual) capable of cutting (or of incising) and/or of removing material are for example scalpels, cutting forceps, or else mechanical or electrical bistoury, or alternatively lasers.

[0029] Moreover, the internal forces module is preferably devised so as, in the event of the deletion of a link and/or of a node and before updating the link tensors and node tensors, to add new mutually independent nodes and new links in such a way as to locally remesh the volume meshing subsequent to the deletion.

[0030] When the device is devised so as to work according to the aforesaid hybrid model, its internal forces module is capable of determining the internal forces exerted on the nodes of at least a second part of the volume meshing on the basis of boundary conditions defined by so-called connection nodes placed at the interface between the first and second parts, and of a table of deformation tensors, each tensor of which is representative of the influence of an elementary displacement of each node of at least the second part on each other node of at least this second part.

[0031] The boundary conditions serving in the calculation of the internal forces of the second part are preferably defined by the internal forces calculated for the connection nodes when calculating the internal forces of the nodes of the first part.

[0032] In the hybrid model, the internal forces module preferably proceeds by successive iterations until a position of so-called <<equilibrium>> of the internal forces of the connection nodes is obtained. To do this, this module is devised so as to deduce from the values of the internal forces exerted on the nodes of the second part of the volume meshing, values of displacement of the connection nodes in such a way as to provide boundary conditions which in turn make it possible to calculate the internal forces of the nodes of the first part.

[0033] The subdivision of the volume meshing into parts (at least two) is determined on the basis of a predetermined criterion pertaining at least to a parameter of the image data of the object chosen from among physical parameters and anatomical parameters, in particular to an intensity. This subdivision can be performed during the formation of the volume meshing of the 3D object. Accordingly, when the device does not comprise a meshing module, the subdivision can be performed by an external processing. The device can comprise a partitioning module, for example forming part of the meshing module if the latter exists, intended to provide the subdivision of the volume meshing.

[0034] Preferably, the <<collision>> determination module is devised so as to determine a collision between at least two tools managed by the user interface. This makes it possible to manage conflicts when one or more operators maneuver at least two tools at the same time.

[0035] For this purpose, each tool is represented by at least one point embodying its end interacting with the object, and a multiplicity of points joined to one another as well as to the end, by segments, embodying its <<shank>>.

[0036] Preferably, the collision determination module estimates the coordinates of the point of intersection between the <<tool>> and a surface mesh cell as follows:

[0037] Firstly, it creates a three-dimensional space encompassing the external envelope of the object, then it decomposes this space into volume blocks, the number of which is chosen so that each block comprises a number of node of the volume meshing of the object substantially equal to the number of nodes contained in the other volume blocks, each block intersecting the external surface comprising at least one node, and lastly it stores in multiplets the coordinates of each node with reference to the volume block which encompasses it. The presence of a point of the tool in the space is then effected, preferably, through a comparison between the multiplets and the coordinates of the point. One of these multiplets makes it possible to designate the volume block of the space in which the point lies.

[0038] Then, preferably, the collision detection module determines the (Euclidean) distance which separates the point of the tool from the node(s) encompassed in the designated volume block so as to determine the smallest of these distances, termed the minimum distance. Thereafter it determines the distance (for example Euclidean) which separates the point from the node(s) encompassed in a predetermined number of volume blocks neighboring the volume block in which it lies so as to compare its distances with the minimum distance. Then, it determines the collection of surface mesh cells adjacent to the node associated with the minimum distance so as to determine whether a segment defined by the position of the point of the tool and by its previous position intersects one at least of these adjacent surface mesh cells. Lastly, from this it deduces the (barycentric) coordinates of the point of intersection (or of collision) between the object and the tool, with a view to their transmission to the internal forces module.

[0039] Preferably, the collision detection module is capable of modifying the contents of the multiplets between two determinations of presence of points of the tool inside the volume blocks, in the event of detection of a deformation of the mesh by the internal forces module. This makes it possible to improve the accuracy of collision detection and hence of the calculation of the deformation.

[0040] In order to allow the most realistic possible simulation of action on one or more tools, the user interface comprises a harness maneuverable by at least one operator hand so as to simulate the maneuvering of each tool. The operator thus maneuvers the harness, which may be a <<joystick>> possibly fitted with actuator(s) or else an articulated <<glove>> fixed on his hand, then the user interface defines the displacement of the associated virtual tool on the basis of this maneuver.

[0041] According to yet another characteristic of the invention, the device comprises display means making it possible to display in real time (and continuously) the image (formed from the image data) of the object and of a representation of the tool.

[0042] The invention applies most particularly to the sets of image data representing a three-dimensional digital image of a region comprising at least one 3D object, including the designated object, and more particularly still to the sets of image data representing a three-dimensional digital image of a region of a living being (animal or human) comprising deformable anatomical structures such as the liver, the kidney, the gall bladder, or alternatively the heart.

[0043] The invention also proposes a process for processing the digital image data for implementing the device described above, comprising the following known steps:

[0044] provide a user interface capable of generating a force feedback, in accordance with the reactions of a tool,

[0045] establish, on the basis of a deformation law and of an action defined by the user interface and representative of a maneuver of the tool, a field of internal forces between nodes of a volume meshing dependent on a surface meshing of a 3D object appearing in a set of image data,

[0046] determine the reaction force of the object which corresponds to its deformation estimated on the basis of the internal forces, so that the force generated by the user interface is substantially balanced by this reaction force,

[0047] calculate new image data of the object, in the presence of the estimated deformations supplemented with the representation of the said tool,

[0048] and characterized in that there is provision for a step in which a point of intersection between a straight line embodying a displacement derived from the defined action and the surface meshing is estimated, and in that the internal force exerted on the nodes of a first part at least of the volume meshing of the object is estimated, on the basis of the displacement applied to the nodes belonging to the surface mesh cell containing the point of intersection, of boundary conditions, and of node tensors and link tensors emanating respectively for each node and each link of this part at least, from stiffness matrices specific to each volume mesh cell of at least the first part and dependent on the deformation law.

[0049] Other characteristics and advantages of the invention will become apparent on examining the detailed description which follows, as well as the appended drawings in which:

[0050] FIG. 1 is a perspective view of part of a simulation device according to the invention;

[0051] FIG. 2 is a functional diagram illustrating the architecture of a simulation device according to the invention in an embodiment with two processing modules;

[0052] FIG. 3 diagrammatically illustrates the mode of displacing a tool controlled by a type of user interface of the device according to the invention;

[0053] FIG. 4 illustrates the external envelope (or external surface) of a human liver, obtained on the basis of a segmentation technique;

[0054] FIG. 5 illustrates the external envelope of the liver of FIG. 4, furnished with a surface mesh of triangular type;

[0055] FIG. 6 is a diagram illustrating a tetrahedral volume mesh obtained from a triangular surface mesh of the type of that of FIG. 5;

[0056] FIG. 7 is a diagram illustrating the displacements at an instant t of two nodes with respect to their respective positions of equilibrium and of the positions of equilibrium of their neighboring nodes;

[0057] FIG. 8 is a diagram illustrating a tetrahedral volume mesh cell and the designation of the node tensors and of the link tensors respectively associated with its nodes and links;

[0058] FIG. 9 is an example of a stiffness matrix of a 3D object represented by 15 nodes, the crosses (X) signifying non-zero values and the dots (.) zero values;

[0059] FIG. 10 is a mechanism illustrating the steps for determining a final mesh of a 3D object at equilibrium;

[0060] FIG. 11 is a mechanism illustrating the steps for determining a deformation and the image of the deformed object;

[0061] FIG. 12 is a diagram illustrating the modification of the topology of a tetrahedral volume mesh following the rupturing of a link;

[0062] FIGS. 13A to 13D illustrate a simulation over four successive images of the action of the end of a virtual tool on a liver lobe;

[0063] FIG. 14 is a diagram illustrating a decomposition of a triangular surface mesh into two domains (or parts), in a hybrid type processing;

[0064] FIG. 15 is a mechanism illustrating a loop for determining positions of equilibrium in a hybrid mesh;

[0065] FIG. 16 is a diagram illustrating a space encompassing the liver of FIG. 4 and decomposed into parallelepipedal volume blocks;

[0066] FIG. 17 is a diagram illustrating part of a tetrahedral volume mesh contained in a parallelepipedal volume block, as well as a point of collision between a tool and one of the triangular surface mesh cells of the part of the said volume mesh; and

[0067] FIG. 18 is a diagram illustrating the embodying of a tool in the form of points and segments, with reference to a part of a triangular surface mesh.

[0068] The drawings are, in essence, of definite character. Accordingly, they form an integral part of the present description. They may therefore serve not only to provide a better understanding of the invention, but also contribute to the definition of the latter.

[0069] In what follows, reference will be made to a processing of a set of digital image data forming a medical three-dimensional (3D) image, and more particularly, but solely by way of example, to images of regions of the liver of the type of that illustrated partially in FIGS. 1 and 4, which have been obtained in a human subject.

[0070] It is however clear that the invention is not limited to the processing of medical images, and still less to that of images of the liver. It applies generally to the processing of digital images of three-dimensional (3D) objects, in particular deformable objects, with a view to simulation of their deformation by at least one virtual tool. The word deformation is understood here within its most general acceptance, namely both of the geometrical (surface or volume) type and also of the incision, or tearing, or alternatively removal of material type.

[0071] In the medical field, a set of image data forming a three-dimensional (3D) image is obtained with the aid of apparatus such as X-ray scanners or nuclear magnetic resonance apparatus (MRI), or more generally any type of apparatus capable of acquiring variable-intensity images. Each elementary part of a region represented by a 3D image is defined by three spatial coordinates and at least one measured physical quantity, such as for example an intensity.

[0072] In the case of MRI, the 3D image of an observed region consists of a multiplicity of stacked 2D sections, within which the intensity variations represent the proton density of the tissues.

[0073] Illustrated partially in FIG. 1 is a device according to the invention, intended for simulating the dynamic behavior of a three-dimensional (3D) object. Here, as indicated earlier, the 3D object is a human liver.

[0074] In the example illustrated, the device comprises a control and calculating unit 1 coupled to a graphics station 2 comprising display means, such as for example a monitor 3, as well as to a user interface 4 having force feedback.

[0075] In this example, the control unit 1 and the graphics station 2 are two separate computers, but it is clear that they could form a single calculating unit.

[0076] In the embodiment illustrated, the control unit 1 is preferably a PC type microcomputer, whereas the graphics station 2 is preferably a DEC ALPHA 233 MHz type workstation equipped with a 3D graphics accelerator card. Here it is preferable to use the most powerful calculating unit to manage the displaying of images, owing to the need to calculate the deformation and to display in real time both the image of the liver and of the representation of a tool, or several tools, this requiring the refreshing of the images at a video type frequency of greater than 20 Hz.

[0077] Moreover, the link between the control unit 1 and the force feedback user interface 4 will preferably be an ISA bus, which permits data transfer rates of 10 kHz which are well above the frequency required for the force feedback, which generally lies between 300 and 500 Hz.

[0078] The effort feedback user interface 4 comprises a manual control 5 (or harness) of the <<joystick>> type capable of reproducing the behavior of a minimally invasive surgical tool. It could in particular be a tool for laparoscopic or coelioscopic surgery. In what follows, the joystick will be regarded as simulating the movement of microscopic cutting forceps intended to be introduced into the body of a patient via a fixed point termed a trocar with a view to an intervention on the liver. Of course, the invention is not limited to this type of surgery or tool alone.

[0079] The effort feedback user interface 4 could for example be the <<laparoscopic impulse engine>> (LIE) from the company Immersion Corporation, or else the PHANToM device developed by the Massachussets Institute of Technology. The <<laparoscopic impulse engine>> (LIE) possesses five degrees of freedom, three to represent the displacement of the tool shank in space, a fourth for the rotation of this shank about an axis passing through the trocar, and a fifth associated with the opening or closing of the forceps located at the end of the shank. The effort feedback can be applied to three of these degrees of freedom with the aid of motors. The forces are transmitted via a system of cables and capstans thus delivering a high and progressive torque.

[0080] Illustrated diagrammatically in FIG. 3 is the mode of displacing the shank 6 of the cutting forceps controlled by the LIE. The end 7 of the shank 6 is embodied here by a sphere which diagrammatically represents the active end (for example the cutting end) of the virtual tool. The trocar is embodied by a fixed point 8 through which the shank 6 of the tool passes, which shank is connected to two motors 9 and 10 which allow it displacements over a spherical cap whose radius varies according to the length of the shank 6, which is managed by a third motor (not represented). The position of the end 7 of the tool on the spherical cap is determined by the values of the angles α and β.

[0081] When, with the aid of his hand 11 (see FIG. 1), an operator maneuvers the grippable part of the joystick 5, he can move the axis 6 in such a way as to vary the length which separates the fixed point 8 (or trocar) from the cutting end 7 thus making it possible to simulate the approaching of the tool towards the 3D object, but he can also move this end 7 over a spherical cap whose radius depends on the trocar/end distance.

[0082] Of course, the user interface 4 can comprise two LIEs so as to allow the simulation of the simultaneous handling of two tools.

[0083] A purpose of the device according to the invention is, as has already been stated, to process, in real time, digital image data forming successive three-dimensional images of a 3D object (here a human liver), which are obtained, for example, by X-ray scanning, or by MRI, depending on the object or organ concerned.

[0084] It is firstly necessary for the device to possess the set of image data forming the image of the 3D object, the external envelope of this object, the surface mesh of this envelope, and the primary volume mesh of the object formed on the basis of the surface mesh. All these data may be obtained outside of the device according to the invention. They could in particular originate from an image data bank coupled to a meshing module. However, the device can comprise a meshing module capable of determining the external envelope and the surface and primary volume meshes on the basis of an image data set provided by an external means.

[0085] The external envelope 12 (see FIG. 4) could be obtained by any technique known to the person skilled in the art, and for example with the aid of a segmentation technique of the iso-surface extraction type, such as the so-called <<3D snakes>> or the <<marching cubes>> techniques which are well known to the person skilled in the art.

[0086] In the case of a meshing module built into the device, the latter will preferably be housed in the graphics station 2.

[0087] The general shape of the surface mesh cells 13 into which the external envelope is decomposed is preferably triangular (see FIG. 5). A criterion for optimizing the surface area of the surface mesh cells could be implemented so as to obtain a regular mesh. Of course, other types of surface meshes may be used, such as in particular the <<simplex>> meshes developed by H. Delingette and disclosed in INRIA technical report No. 2214 entitled <<Simplex meshes: a general representation for 3D shape reconstruction>>, 1994.

[0088] A volume mesh of the 3D object, here the liver, is determined from the surface mesh cells 13. It is clear that the general shape of the mesh cells depends on the shape of the surface mesh cells 13. In this example where the surface mesh cells 13 are triangular, the volume mesh cells (see FIG. 6) are advantageously tetrahedra 14 at the four vertices of which are located nodes 15 over which the object's mass, which is otherwise known, is distributed.

[0089] This decomposition of the volume encompassed by the external envelope 12 into mesh cells, here tetrahedra, can be obtained with the aid of a Delaunay-Voronoï type algorithm using, for example, the Simail module marketed by the French company SIMULOG S.A. Of course, other algorithms could be used with a view to obtaining the volume mesh.

[0090] As was indicated earlier, the device according to the invention does not necessarily comprise a meshing module intended for determining the external envelope of the 3D object and/or the decomposition into surface mesh cells and/or the decomposition into volume mesh cells. All or some of this processing of the image data set of the 3D object can be performed in a manner totally separate of the device according to the invention.

[0091] Knowing the set of image data of the 3D object, and its decomposition into volume mesh cells 14, here tetrahedral, the device will determine the so-called <<internal>> forces which are exerted on each node 15 of the volume mesh. It is clear that the sum of the internal forces exerted on a given node is different from zero only when this node is shifted from its equilibrium position. To do this, an external force must be applied to at least one of the nodes, for example by a tool.

[0092] The internal forces may be determined, either over the entire volume mesh with the aid of a first method using a single so-called <<masses/tensors>> technique which will be described below, or with the aid at least of a second method using the first masses/tensors technique as well as a second, different technique, each technique being applied to at least one domain arising from a decomposition of the volume mesh. This second so-called <<hybrid>> method will be described later, with reference to FIGS. 14 and 15.

[0093] Nevertheless, regardless of the method, the determination of the internal forces is performed by virtue of an internal forces module 16 (see FIG. 2) which is preferably built into the graphics station 2 (but it could be built into the control unit 1), and which requires the assistance of a collision detection module 18 built, preferably, into the control unit 1.

[0094] The purpose of such a module is to accurately determine the locus of collision (or of intersection) between a point embodying a part of the tool and the external envelope of the 3D object, or more precisely a surface mesh cell of this envelope (see FIGS. 16 to 18). This makes it possible to provide the internal forces module 16 with the designation of at least one surface mesh cell node on which the tool acts, and which will make it possible to determine the internal forces exerted at each other node of the volume mesh cell.

[0095] The level of accuracy of the point of intersection (or of collision) can be managed on the basis of the number of points used to embody, on the one hand, the end 7 of the tool, and on the other hand, the latter's shank 6. In other words, the collision detection performed by the collision detection module 18 can be applied at one or more points embodying the tool. An example of discretizing the shape of the tool (shank and end) is given in FIG. 18, in which each point other than the points embodying the nodes of a surface mesh cell 13 are joined together by segments so as to reconstruct the shape of the tool and more particularly that of its shank and of its end, here cutting forceps. It is possible to apply the collision detection to the shank 6 of the tool so as to determine a collision between two shanks of simultaneously actuated tools.

[0096] The collision detection module 18 operates as described below. Initially, the tool is discretized into a certain number of points p_k, which depends both on the instrument used and on the accuracy desired. Next, a space is defined which encompasses the 3D object concerned. This space is then subdivided into volume blocks, for example, of parallelepipedal shapes. The number of volume blocks in each direction in space may be parametrizable, and is preferably chosen as a function of the geometry of the 3D object. Preferably, the parametrization is chosen so that the number of nodes of the volume mesh of the 3D object is substantially identical in all the volume blocks of the space defined by the collision detection module and containing the object. Likewise preferably, a constraint is imposed which consists in insisting that each volume block intersecting the external surface of the object encompass at least one node of the mesh.

[0097] On the basis of this decomposition of the space which encompasses the 3D object into volume blocks, a table is generated in which are stored multiplets comprising the coordinates of each node of the volume mesh with reference to the definition of the volume block which encompasses it. The determination of a point of collision at the level of a surface mesh cell 13 of the external envelope of the 3D object is performed with the aid of a comparison between the position of a moving point p_kof the tool, which position is deduced from the information provided by the force feedback user interface 4, as well as on the basis of the multiplets stored in the table.

[0098] To do this, a code function is defined which makes it possible to calculate an entry inside the table and to find the list of nodes of the volume mesh which correspond to this entry, the said nodes then belonging to the same volume block as the moving point p_k. Such a code function is given by the formulae referenced (1) in the appendix, by way of example.

[0099] Then, for each node N_isituated in the same volume block as the moving point p_k, the Euclidean distance between p_kand each of the N_iis determined. The smallest of the distances thus calculated is then denoted d_min. Next, the distance d, which separates the moving point p_kfrom each node N_jsituated in a predetermined number of neighboring volume blocks, for example 26, is calculated. If one of these distances d, associated with a node N, is less than d_min, then d_min=d. It is clear that the larger the number of nodes inside a volume block, the longer the calculation time will be.

[0100] Next, the number of surface mesh cells 13 (here triangles) which are adjacent to the node N associated with the minimum distance d_min=d is determined. Then, for each surface mesh cell adjacent to N, a check is made as to whether it is intersected by the segment [P_t, P_t−1] defined by the positions of the moving point p_kat the instants t and t−1, as illustrated in FIG. 17. Thus, regardless of the rate of displacement of the tool, it will always be possible to determine the point of intersection (or of collision), given that only the length of the segment [P_t, P_t−1] varies as a function of the rate.

[0101] The position of the point of collision is sub-triangular. Consequently, the coordinates of this point of collision are determined by a barycentric method. These barycentric coordinates are then transmitted to the internal forces module 16 so that the force exerted by the tool on the volume mesh is estimated accurately.

[0102] Knowing the surface mesh cell intersected, the internal forces module 16 deduces therefrom the nodes on which the action of the tool is exerted. Then, from the barycentric coordinates, it deduces the distribution of the displacement of the virtual tool over each of the nodes of the designated surface mesh cell. At the same time, the user interface 4 provides the internal forces module 16 with the information relating to the vector displacement of the tool, either directly via a bus, or via the control unit 1.

[0103] The device according to the invention uses a technique referred to as <<masses/tensors>> in which the nodes belonging to the surface mesh cell, the locus of the collision with the tool, are subjected to the vector displacement provided (distributed according to the barycentric coordinates of the said nodes), then the internal forces applied at each node of the volume mesh are estimated. This estimation is preferably iterative, each iteration making it possible to calculate the internal forces applied at a first <<level>> of nodes, and those of a next <<level>> made up of the nodes neighboring the nodes of this first <<level>>, and so on and so forth until all the nodes have been processed. This constitutes a kind of node-to-node <<constraint propagation>> mechanism. In fact, during the first iteration, only the internal forces of a first node level are non-zero, those of the other nodes being zero. It is in the course of the succeeding iterations that progressively (<<level>> after <<level>>) the (non-zero) internal forces of the other nodes of the mesh are obtained.

[0104] This technique makes it possible to estimate the internal forces independently of the direction of the links. Thus, the deformation induced by the virtual tool on the 3D object is independent of its volume mesh. Such a characteristic of independence relative to the initial mesh is particularly beneficial in the field of simulation, especially of surgical intervention, where the deformations induced by a tool are generally volumic as in the case of an incision, or of a tear, or alternatively of a removal of material.

[0105] In an approach where the deformation law is of volume linear elastic type, the internal force exerted on a node N_isubjected to a displacement u is given by the formula referenced (2) in the appendix.

[0106] This formula (2) can be decomposed into two parts, or two forces, one associated with the actual displacement of the node (or mass) N_i, with respect to its equilibrium position, the other associated with the displacement of a neighboring node N₃. Let ν be the collection of nodes N_jneighboring node N_i, formula 2 for the internal forces applied to node N_ican then be rewritten as indicated in the formula (3) given in the appendix.

[0107] In this formula (3), the expressions between square brackets, [T_ii] and [T_ij], represent matrices of dimension 3×3 which will subsequently be referred to respectively as node tensors and link tensors, on account of their being defined similarly to force tensors which are well known to the person skilled in the art.

[0108] To determine these node and link tensors, it is necessary, firstly, to evaluate the elementary stiffness matrix [K^e] associated with each volume mesh cell (here a tetrahedron) adjacent to each link arising from node N_i.

[0109] As is represented in FIG. 7, the vector displacement of a node N_ibetween an instant indexed 0 and the instant t may be denoted Δx_i. This displacement comprises three components along the three directions in space.

[0110] A stiffness matrix represents, in respect of a 3D object, in particular a deformable one, the various relations that exist between its constituent elements. In other words, it conveys the influence of a node on other nodes, this being given by the deformation law. The numerical values of the components of a stiffness matrix consequently depend on the elasticity parameters of the object, which are known. The deformation law can be stored in the internal forces module 16 of the graphics station 2.

[0111] If T denotes a tetrahedral volume mesh cell, then the elementary stiffness matrix of this tetrahedron is given by the formula (4) given in the appendix, by way of example.

[0112] Once the stiffness matrix [K^e] of each tetrahedron T has been calculated, a node tensor T_ijor a link tensor T_ijis associated with each tetrahedron vertex (or each node), as well as with each of its links, as illustrated in FIG. 8. This association is performed on the basis of implicit formulae given in the appendix under the reference (5).

[0113] The overall stiffness matrix [K] of the 3D object can be determined from the set of elementary stiffness matrices [K^e]. Such an overall matrix is illustrated, by way of example, in FIG. 9, in the case of a 3D object represented by 15 nodes referenced from 1 to F. Here, the crosses represent non-zero values, whilst the dots represent zero values.

[0114] In order to simplify the calculations of the internal forces, it is possible to limit the number of links in the volume mesh cell by filtering the values of the components of the elementary stiffness matrices [K^e], for example by comparing them with a threshold. Thus, any value below this threshold is regarded as zero, this amounting to deleting links between two nodes. This technique for simplifying the matrices makes it possible to determine only the links which are vital to the transmission of the constraints (displacement constraints) through the volume mesh.

[0115] Moreover, it is also possible to simplify the calculation of the link tensors by taking into account the symmetry property of the overall stiffness matrix [K]. Thus, a single, unique link tensor is sufficient to characterize each link. In other words, let l_ijbe the link defined between nodes N_iand N_j, the link tensor associated with this link will be [T_ij] if i>j and [T_ij]^tif i<j. The exponent t denotes the transpose of the matrix concerned.

[0116] Having reached this juncture, the following are available: firstly, the coordinates of each node N_i, secondly, the external envelope, the surface mesh of this external envelope and the final volume mesh of the entire 3D object, thirdly the links l_ijbetween two nodes Ni and Nj of the set of tetrahedra T making up the volume mesh of the 3D object, and fourthly, the set of node tensors [T_ii] and link tensors [T_ij] associated with each node and each link of each tetrahedron T.

[0117] A data structure (or base) will then preferably be defined, making it possible to ascertain, firstly, the index of the two nodes situated at the ends of each link, secondly, the reference of a tetrahedron (volume mesh cell) adjacent to a triangle (surface mesh cell), thirdly, the state of a tetrahedron (active or inactive), and fourthly, the state of a link (active or inactive). This information could be useful for taking into account deformations of the incision, tear or removal of material type, both at graphical and physical level. The various steps for obtaining all these parameters are summarized in FIG. 10.

[0118] In order to take into account auxiliary surface forces, making it possible for example to amplify an effect felt during an incision or when cutting, the internal forces module 16 of the device according to the invention can combine the surface and volume link forces described earlier with tensile forces F_i^tensionof the type of those given in formula (6) in the appendix. These forces may be chosen by varying l₀the unloaded length of the associated spring, as well as its stiffness, depending on the magnitude of the amplification desired. This magnitude may be defined, for example, by the operator, or else by the user interface 4 depending on the type of tool and/or the type of object on which the intervention takes place.

[0119] When such a consideration is taken into account, the formula (3) for the internal force exerted on a node N_iturns into the formula (7), given in the appendix.

[0120] Moreover, it is also possible to take account in the calculation of the internal forces exerted on each node of the volume mesh, of one or more types of exterior forces due to the presence of other objects in the neighborhood of the 3D object concerned, as well as to the gravitational forces which cause a 3D object to tend to sag under its own weight owing to gravity. Such exterior forces may for example be the weight exerted at a part of the 3D object by a neighboring object, or else the force of attachment exerted by a tie, such as for example a ligament, between the 3D object and a neighboring object.

[0121] It is clear that a force of the attachment type may be regarded as a boundary condition akin, for example, to a restoring force in a damped spring possessing a fixed end, the other end being joined to a node of a surface mesh cell of the 3D object.

[0122] The determination of the deformation of a 3D object, or more precisely of its volume mesh, requires the knowledge at each instant of the displacements of each node of this volume mesh. This real-time determination relies on a principle of dynamic animation which relates the sum of the forces applied to the deformable 3D object to its acceleration. Here, the deformable 3D object is formed by a collection of point masses (or nodes N_i) joined together by link forces. It is consequently possible to apply the general equation of mechanics to a volume mesh in such a way as to simulate the dynamic behavior of the 3D object which it represents when the latter is deformed by exterior forces generated by a tool whose displacement is defined by the force feedback user interface 4, in reaction to a maneuvering of the joystick 5 by an operator.

[0123] At a given instant t, the force applied to a node N_iof mass m_iis therefore given by the general equation of dynamics given in formula (8) in the appendix.

[0124] In what follows, it will be assumed that all the nodes of the volume mesh possess the same mass (m_i=μ). Of course, it would be possible to proceed differently, each node possessing its own mass, which can be determined by knowing the makeup of the region which it represents. This makeup information may form part of the set of image data forming the 3D digital image of the object. However, it may also be provided separately.

[0125] It is therefore possible to determine the acceleration a_iof each mass or node N_i, and consequently to retrieve firstly its velocity v_i, then secondly its position x_i+1at an instant t+Δt by performing integrations, for example numerical integrations. It is this calculation of all the positions x_i+1at the instant t+Δt which will provide the displacements of each node and consequently the deformation of the volume mesh cell. To do this, the following types of method could be used: Euler, or Runge-Kutta of order 2, or more preferably of order 4.

[0126] The Euler method is the one most commonly employed. It consists, for example, in starting from a differential equation dx/dt=f(x,t), with initial condition x=x₀at the instant t=0. A time step At is then chosen, such that T_i=iΔt, with i=0,1,2 . . . . Knowing the position of a node x_iat the instant ti, it is then possible to deduce the position of this same node at the instant t_i+1=t+Δt (see formula (9) in the appendix).

[0127] Two additions, two multiplications and an evaluation of the forces are then necessary to obtain the position of a node N together with its velocity at an instant t+Δt, from its acceleration at the instant t. The velocity and the position of the node are given by way of example in the formulae (10) in the appendix.

[0128] The Runge-Kutta formula of order 4 is currently preferred on account in particular of its numerical accuracy. It relies, for example, on the formulae (11) given in the appendix.

[0129] In this method, the external and internal forces at each node must be evaluated four times.

[0130] Of course, other integration techniques, in particular numerical ones, may be used.

[0131] In the left part of FIG. 11 is represented the sequencing of the steps which make it possible to determine the displacement of each node of the volume mesh cell, or of the selected part thereof, from the vector displacement of the end 7 of the shank 6 of the virtual tool relative to the 3D object, due to the action of the operator on the joystick 5. This vector displacement is provided by the force feedback user interface 4, via the ISA bus, preferably at the same time as the position of the end 7.

[0132] As soon as the collision detection module 18 has determined a point of intersection between the <<tool>> and the surface mesh, the internal forces module 16 can begin determining the internal forces exerted on the nodes of the first part of the volume mesh.

[0133] To do this, it applies a vector displacement, equal to that provided by the user interface 4 which defines the displacement of the tool, to each node N_cof the collision surface mesh cell detected (here there are three owing to the triangular shape of the mesh cell). Next, it determines the internal force exerted on each node of the first part of the mesh by using formulae (3) or (7) depending on the external forces chosen by the operator and/or the user interface 4. This determination calls upon the node tensors and link tensors, [T_ii] and [T_ij], previously determined and stored. In fact, during the first iteration, although the internal forces applied at each node are calculated, only those of the nodes neighboring the <<collision>> nodes N_care non-zero, given that the displacements of the other nodes are zero at this juncture in the calculation. From this is deduced the acceleration, the velocity and the displacement of the first-neighbor nodes, and consequently their new positions.

[0134] These new positions of the displaced nodes are then frozen for the remainder of the calculation, and the internal forces module 16 repeats its calculation so as to determine at each node of the first part (here the whole mesh cell), the internal force (non-zero) which is applied to it, and consequently its displacement. The internal forces are therefore calculated according to a model of the <<propagation of constraints>> type. The deformation induced by the vector displacement of the virtual <<tool>> on the volume mesh cell representing the 3D object is then known.

[0135] Between two loops for calculating the internal forces and the displacements due to a vector displacement of the <<tool>>, that is to say while propagating constraints, it is particularly advantageous to take into account information transmitted by the user interface 4, or defined by the operator, pertaining to the type of tool maneuvered and consequently the type of intervention performed. Indeed, such information may indicate to the internal forces module 16 that the action of the tool is intended to cut, or incise, or eliminate material. Thus, after an iteration during which internal forces (non-zero) exerted on certain nodes have been determined, the internal forces module 16 can delete the links which are in the cutting plane of a tool for incising or removing material (embodied by volume mesh cells, here tetrahedra) which lies in the zone of interaction of the virtual <<tool>>, on the basis of the action information with which it is provided by the collision detection module 18 and the user interface 4.

[0136] In this case, a certain number of links and/or of nodes are deleted. The internal forces module must then update the node tensors and the link tensors as a function of the deleted link(s) and/or node(s). This taking into account of the aforesaid information serves as a deletion criterion.

[0137] Virtual tools capable of cutting (or of incising) and/or of removing material are for example scalpels, cutting forceps, or else mechanical or electrical bistoury, or alternatively lasers.

[0138] Moreover, it is also particularly advantageous for the deletion criterion to be capable, between two calculation loops, of performing a test pertaining to the integrity of the links between nodes whose displacements have just been determined and are henceforth frozen.

[0139] The purpose of this test is to determine whether certain displacements are such that in reality they would have caused the rupture of one or more links and/or the deletion of one or more nodes and hence of one or more tetrahedra (or volume mesh cells). The test consists in comparing the displacements with a predetermined link rupture threshold. This test could pertain, for example, to the volume variation of the volume mesh cell comprising the mesh cell to be deleted, and/or the length variation of the links of the volume mesh cell comprising the link to be deleted.

[0140] If no link is to be deleted, the internal forces module 16 goes to the next iteration (next loop) so as to continue the calculation of the internal forces (continue the propagation of constraints).

[0141] On the other hand, if the test indicates that at least one link is to be destroyed, that is to say if the intensity of the displacement exceeds the predetermined threshold, then the internal forces module 16 updates the node tensors [T_ii] and the link tensors [T_ij] of the volume mesh cell or of the first part of the latter, by taking account of the calculated displacements. The value of certain node tensors, and link tensors, may then become zero. This is the case in particular when the four links which join a node to the other four nodes of a given tetrahedron break. In this case, a node of the volume mesh cell, or several if necessary, must be deleted together with all the associated links.

[0142] The operation does not stop at this juncture, given that, for reasons concerned with the equilibrium of the deformed volume mesh cell, the latter must retain its general starting structure at the level of its surface and volume meshes. Now, the removal of a node, or the rupturing of a link give rise to modifications to one or more mesh cells, these no longer being, in the example illustrated, tetrahedral or triangular depending on whether a volume mesh cell or a surface mesh cell is being spoken of. Such a case is illustrated in FIG. 12 in which the action of the end 7 of the tool has given rise to the rupturing of four links in a particular sectional plane. Each of the four links being divided into two, it is then necessary to add nodes, and consequently links, so as to reconstruct tetrahedral mesh cells in the incised zone. In other words, the volume mesh must be locally remeshed.

[0143] Given that the result of adding nodes is an increase in the total mass of the 3D object, it is then necessary to modify the mass of all the other nodes of the volume mesh, and this may be done by a homogeneous distribution.

[0144] Moreover, adding nodes entails adding links so that the mesh cells retain their general shape, here triangular or tetrahedral. Consequently, the internal forces module 16 must create and/or destroy elementary stiffness matrices [K^e] and node and link tensors, and recalculate the elementary stiffness matrices and the node tensors and the link tensors previously calculated, so that the calculations of the internal forces exerted on the nodes can be continued while taking account of these modifications.

[0145] It is clear that it is particularly advantageous for the collision module 18 to be able to update its multiplets alongside the alterations (deformation) in the part of the volume mesh cell.

[0146] As soon as the calculation of the internal forces and displacements is regarded by the internal forces module as having terminated for the value provided of the vector displacement of the tool, the refresh module 17 proceeds to the determination of the new image data of the object and of the embodying of the tool with a view to their display, whilst at the same time, the internal forces module recommences a new calculation with a new vector displacement of the tool, possibly applied to nodes other than those to which they were initially applied.

[0147] Substantially simultaneously with the determination of the new image data of the object and of the <<tool>>, a reaction module 20, preferably built into the control unit 1, determines the reaction force of the object, which corresponds to its deformation (displacements of the nodes) estimated on the basis of the internal forces. For each node displaced by the tool, its initial position (at t=0) and its current position (at the instant t) are known. Moreover, the tensors associated with this node and the links which are associated with it are known, thus making it possible to apply formula (3) given in the appendix.

[0148] The reaction force is then transmitted to the force feedback user interface (LIE) so that it generates on the joystick [lacuna] force (or force feedback) which is substantially balanced by the said reaction force. This makes it possible to transmit physical reactions of the object to the operator, this being essential for good control of the intervention maneuver.

[0149] As emerges from the above description, the masses/tensors calculation technique makes it possible to simulate all types of deformations, be they of geometrical type, of incision or cutting type, of tearing or fracturing type, or else of removal of material type.

[0150] Represented in FIG. 13 are four images of a human liver subjected to the action of a virtual tool of the cutting forceps type. Here, the action is a simple compression.

[0151] In order to improve the quality of the image, it is possible to envisage, before each update of the scene displayed, repeating the complete calculation of the deformation of the volume mesh cell a number of times, and consequently the tests pertaining to the node and link deletions, as well as those pertaining to the values of the components of the tensors. Thus, as long as the number of iterations is below a chosen threshold, and as long as the calculation time remains below the time required for the refreshing of the images on the display means, then the loop can be recommenced. On the other hand, if this is not the case, the updating of the scene displayed is carried out.

[0152] All the main steps of the masses/tensors technique are illustrated in the chart of FIG. 11.

[0153] The invention also proposes an electronic process for processing image data for implementing the aforesaid device, in which the following main steps are provided:

[0154] firstly, provide a user interface (4) capable of giving a force feedback, in accordance with the reactions of a tool,

[0155] then, estimate a point of intersection (or of collision) between a straight line embodying a displacement derived from the action defined by the user interface and the surface mesh,

[0156] then, establish a field of internal forces between nodes of a first part at least of a volume mesh dependent on a surface mesh of a 3D object appearing in a set of image data, on the basis of a deformation law, for example of volume linear elastic type, of a displacement, induced by an action defined by the user interface and representative of a maneuvering of the tool, applied to the nodes belonging to the surface mesh cell containing the point of intersection, of boundary conditions, and of node tensors and link tensors emanating respectively for each node and each link of this first part at least, from stiffness matrices specific to each volume mesh cell of at least the first part and dependent on the deformation law,

[0157] next, determine the reaction force of the object which corresponds to its deformation estimated on the basis of the internal forces, so that the force generated by the user interface is substantially balanced by this reaction force,

[0158] and lastly, calculate new image data of the object, in the presence of the estimated deformations supplemented with the representation of the tool.

[0159] In certain complex cases, where the number of image data forming the three-dimensional digital image of a 3D object is very large, the time required for calculating the internal forces and the deformation of the entire volume mesh may become greater than the time for refreshing the images on the display means (here a monitor). It follows that the application of the technique described above cannot be effected on the entire 3D object. This is why the device according to the invention may apply a second so-called hybrid method to the image data, in which the masses/tensors technique will be applied only to a first part of the volume mesh, whilst the complementary part, termed the second part, of this volume mesh will be processed with the aid of a second real-time calculation technique relying for example on a finite element type approach calling upon stored precalculations. This second technique is described in particular in the already cited article by S. Cotin, H. Delingette, M. Bro-Nielsen and N. Ayache, <<Geometric and physical representations for a simulator of hepatic surgery>>, published in the proceedings of the conference Medecine meets with virtual reality of January 1996.

[0160] In what follows, that part of the volume mesh to which the masses/tensors technique is applied will be called the first part, and that part which is complementary to the first part, to which the finite element technique is applied via stored precalculations, will be called the second part.

[0161] In the hybrid model, the volume mesh (initial or final) is subdivided into at least two domains. It is clear that it would be possible to envisage a subdivision into three domains, or even more, by virtue of a processing, for example in parallel, of the various domains. This subdivision into domains may be performed by the device according to the invention, when it comprises a special-purpose meshing module of the type of that described earlier, but it may also be performed outside the device in the converse case.

[0162] Whether it is done by an external module or by a module internal to the device, the subdivision into domains results from a selection by an operator or from an automatic detection resulting from a predetermined criterion applied to the image data of the 3D object. The criterion could in particular be a comparison of intensity between the various image data or else between these image data and a threshold, for example in the case of a tumor, or a detection of a physical or anatomical parameter contained in the said image data of the 3D object, and designating complex zones in which there are substructures such as for example tumors.

[0163] It is possible to envisage the mesh being subdivided into a multiplicity of regions, for example six or seven, and for that in which the collision with the <<tool>> is detected to be automatically designated as the first part, the complementary part then forming the second part of the volume mesh.

[0164] The device could comprise a specific partitioning module, different from the meshing module, for determining the subdivision into domains. This module will preferably be built into the internal forces module 16.

[0165] In the hybrid model, once the subdivision into parts has been performed, the internal forces module 16 begins the calculation of the internal forces of the entire volume mesh of the 3D object. This calculation can begin either with the first part of the mesh, with the aid of the masses/tensors technique, or with the second part of the mesh, with the aid of the finite element technique via precalculations (the so-called <<precalculations>> technique) which relies on two principles termed superposition and linearity.

[0166] The vector displacement of the tool induces a vector displacement of certain nodes of the volume mesh, which is regarded as a multiple of an elementary vector displacement. Then, the total displacement of each other node i (those on which the action of the tool is not exerted directly) is determined by summing (superposition principle) the contributions from the displacements of all the nodes j having undergone a displacement.

[0167] In order to implement this principle of superposition, it is necessary to determine a pre-calculation table for n×n deformation tensors, where n represents the total number of nodes of the mesh, or of a part thereof (for example only the nodes of the surface). To do this, an elementary displacement is applied at a given node, then the deformation which this displacement induces on all the other nodes of the mesh is calculated. From this are deduced n first elementary deformation tensors, one for each node. Then the two aforesaid operations are recommenced, applying the elementary displacement to another node of the mesh, this providing n second deformation tensors. These two operations are repeated for the n nodes of the mesh, until the n×n deformation tensors are obtained and which are stored in the form of a so-called pre calculation table.

[0168] A source of information is then available which makes it possible instantaneously to ascertain the deformation induced on the mesh cell by an elementary displacement of a given node. It is then merely necessary to determine the ratio of the vector displacement of the node to the elementary displacement, then to multiply the displacement of the nodes of the mesh stored in the table by this ratio (linearity principle). Stated otherwise, the total displacement of a node N_iis obtained by summing the contributions from all the deformation tensors T′_ijdenoted by the index i in the table (T′_i0, T′_i1, T′_i2, . . . , T′_in) without neglecting to multiply these contributions by the ratio.

[0169] More detailed descriptive elements of the second so-called precalculations technique are indicated, in particular, in the publication by S. Cotin, H. Delingette and N. Ayache, <<Real-time non linear elastic deformations of soft tissues for surgery simulation>>, INRIA research report, which will be made public after the present Patent Application has been filed.

[0170] Of course, as in the masses/tensors technique, the second technique takes account of the boundary conditions, and other auxiliary and/or external forces.

[0171] Regardless of which part of the mesh the calculation of the internal forces begins with, it is necessary to take account of the <<connection>> nodes 19 (see FIG. 14) which are located at the interface between the two parts, and which consequently belong to both these parts at the same time. According to the invention, the connection nodes define boundary conditions and consequently they are parameters of the deformation of the complete mesh.

[0172] If we argue on the basis of beginning the determination of the internal forces on the first part (with the aid of the masses/tensors technique), but it would be possible to do the opposite, we determine the internal forces, and consequently the displacements, of all the nodes of the first part of the mesh, including those of the connection nodes. The displacements of these nodes will then make it possible to determine, with the aid of the table of n×n precalculated and stored deformation tensors T′_ij, according to the second technique described above, the displacements of all the nodes of the second part of the mesh. This provides new boundary conditions (external forces) for the first part of the mesh, and these will be taken into account for recalculating the internal forces applied at each node of the first part, and consequently the displacements of its nodes, including those of the connection nodes, this providing new boundary conditions for the second part of the mesh which will in turn make it possible to recalculate the displacements of the nodes of the second part.

[0173] This procedure is repeated, preferably, until equilibrium is obtained between the two parts of the mesh, as regards the forces exerted on the connection nodes. This equilibrium can be fixed by a convergence threshold, or else by a <<floating>> number of iterations depending on the calculation time permitted, which must necessarily be less than the time required for the refreshing of an image (frequency greater than 20 Hz). The animation loop of the hybrid model is illustrated in FIG. 15 where the meshes M1 and M2 respectively represent the first and second parts of the volume mesh.

[0174] It is clear that the internal forces module 16 takes account of the deletions and additions of links and/or of nodes when it uses the hybrid model. Consequently, the node and link tensors used in the first masses/tensors technique are updated in real time.

[0175] As soon as the calculation of the internal forces and displacements is regarded by the internal forces module 16 as having terminated for a value of the vector displacement of the tool, the said refresh module 17 proceeds to the determination of the new image data of the object and of the embodying of the tool with a view to their display, whilst at the same time, the internal forces module 16 recommences a new calculation with a new vector displacement of the tool, possibly applied to nodes other than those to which they were initially applied.

[0176] It is clear that in the case of a hybrid model, it is particularly advantageous for the collision module to define subtables of multiplets for each part of the volume mesh cell. Moreover, it is particularly advantageous, in a desire for accuracy when detecting the points of collision, for the multiplets of the subtables to be updated alongside the alterations (deformation) of the various parts of the volume mesh cell.

[0177] The device according to the invention can be implanted in a memory, for example a mass memory, of one or more means of calculation of the work station and/or computer type, in the form of software modules.

[0178] For all useful purposes, it is specified that more detailed descriptive elements are indicated in the doctoral thesis by S. Cotin, <<Modèles anatomiques déformables en temps-réel—Application à la simulation de chirurgie avec retour d'effort>> [Real-time deformable anatomical models—Application to the simulation of surgery with force feedback], submitted on Nov. 19, 1997, as well as in the already cited publication <<Real-time non linear elastic deformations of soft tissues for surgery simulation>>, by S. Cotin, H. Delingette and N. Ayache. These documents will be made public after the present Patent Application has been filed.

[0179] The invention is not limited to the embodiments described above, merely by way of example, but it encompasses all the variants which could be envisaged by the person skilled in the art within the scope of the Claims below.

[0180] Thus, a 3D image processing device, and the associated process, have been described, in which the action on a tool was induced by the maneuvering of a <<joystick>> by an operator, then defined by the user interface. However, it would be possible to envisage the action on the tool being defined directly by the user interface, this action being stored and decomposable into elementary subactions, the operator now intervening only so as to decide on the choice of certain elementary subactions, for example with a view to training (or learning).

[0181] Moreover, devices and the associated process have been described, which are intended for processing images, especially medical, and in particular of the liver. However, it is clear that the invention is not limited to this field alone. It applies also to sectors other than the medical field, as well as to non-medical fields in which the real-time processing of images of 3D objects on which a simulation of intervention, in the widest sense of the term, is of particular benefit equally in regard to the teaching of intervention techniques as in regard to the improving or tailoring of new intervention techniques.

[0182] Additionally, processing models have been described in which the matrices and tensors were calculated on the basis of a volume linear elastic deformation law. However, other types of deformation laws may be used, and in particular laws of the non-linear type.

[0183] Lastly, a <<hybrid>> model has been described which calls upon two different techniques, a first technique of real-time calculation referred to as <<masses/tensors>>, the second technique being, in the example chosen, a so-called <<precalculations>> real-time calculation technique relying on a finite element approach. However, it is clear that any other second technique could be envisaged, whether or not it calls upon precalculations, provided that it allows image data processing compatible with a dynamic simulation of a real-time intervention.

[0184] code=x+y*N_{blocks x}+z*N_{blocks y}*N_{blocks x} $x=\mathrm{int}\left(\frac{\mathrm{Px}-{b\min }_x}{{\mathrm{block}}_{\mathrm{width}x}}\right)$ $y=\mathrm{int}\left(\frac{\mathrm{Py}-{b\min }_y}{{\mathrm{block}}_{\mathrm{width}y}}\right)$ $z=\mathrm{int}\left(\frac{\mathrm{Pz}-{b\min }_z}{{\mathrm{block}}_{\mathrm{width}z}}\right)$

[0185] Formulae (1):

[0186] <<int>>: function giving the integer part of its argument,

[0187] bmin_i(i=x,y,z): lower bounds on the space encompassing the object,

[0188] block_{width i}(i=x,y,z): dimensions of a volume block,

[0189] N_{blocks i}(i=x,y): number of volume blocks in the x and y directions in space.

[0190] Formula (2):

[0191] F_i^int(t)=[M]u

[0192] u: vector of dimension 3*(N_v+1) representing the displacement of node N_iand of its neighbors,

[0193] N_v: number of neighbors of node N_i, and

[0194] [M]: deformation matrix of dimension 3*3*(NV+1) representing the deformation law of the object (known otherwise).

[0195] Formula (3): $F_i^{\mathrm{int}}=\left[T_{\mathrm{ii}}\right]\Delta x_i+\sum _{j\in v}\left[T_{\mathrm{ij}}\right]\Delta x_j$

[0196] Formula (4):

[0197] [K^e]=Volume(T)[Γ_T]^t[E_T][Γ_T]

[0198] [E_T] is, for example, given by the matrix: $\left[E_T\right]=\left|\begin{array}{ccccccccc}\lambda +2\mu & 0& 0& 0& \lambda & 0& 0& 0& \lambda \\ 0& \mu & 0& \mu & 0& 0& 0& 0& 0\\ 0& 0& \mu & 0& 0& 0& \mu & 0& 0\\ 0& \mu & 0& \mu & 0& 0& 0& 0& 0\\ \lambda & 0& 0& 0& \lambda +2\mu & 0& 0& 0& \lambda \\ 0& 0& 0& 0& 0& \mu & 0& \mu & 0\\ 0& 0& \mu & 0& 0& 0& \mu & 0& 0\\ 0& 0& 0& 0& 0& \mu & 0& \mu & 0\\ \lambda & 0& 0& 0& \lambda & 0& 0& 0& \lambda +2\mu \end{array}\right|$

[0199] λ and u: LAMÉcoefficients,

[0200] [Γ_T]^tis given, for example, by the matrix: ${\left[T_T\right]}^t=\left|\begin{array}{ccccccccc}\mathrm{a11}& \mathrm{a12}& \mathrm{a22}& 0& 0& 0& 0& 0& 0\\ \mathrm{a21}& \mathrm{a22}& \mathrm{a23}& 0& 0& 0& 0& 0& 0\\ \mathrm{a31}& \mathrm{a32}& \mathrm{a33}& 0& 0& 0& 0& 0& 0\\ \mathrm{a41}& \mathrm{a42}& \mathrm{a43}& 0& 0& 0& 0& 0& 0\\ 0& 0& 0& \mathrm{a11}& \mathrm{a12}& \mathrm{a22}& 0& 0& 0\\ 0& 0& 0& \mathrm{a21}& \mathrm{a22}& \mathrm{a23}& 0& 0& 0\\ 0& 0& 0& \mathrm{a31}& \mathrm{a32}& \mathrm{a33}& 0& 0& 0\\ 0& 0& 0& \mathrm{a41}& \mathrm{a42}& \mathrm{a43}& 0& 0& 0\\ 0& 0& 0& 0& 0& 0& \mathrm{a11}& \mathrm{a12}& \mathrm{a22}\\ 0& 0& 0& 0& 0& 0& \mathrm{a21}& \mathrm{a22}& \mathrm{a23}\\ 0& 0& 0& 0& 0& 0& \mathrm{a31}& \mathrm{a32}& \mathrm{a33}\\ 0& 0& 0& 0& 0& 0& \mathrm{a41}& \mathrm{a42}& \mathrm{a43}\end{array}\right|$

[0201] α_jifirst partial derivatives of the barycentric coordinates λ_j(x) of a point x situated inside a tetrahedron T: $\frac{\delta {\lambda }_j\left(x\right)}{\delta x_i}={\alpha }_{\mathrm{ji}},i=1,2,3;j=1,2,3,4$

[0202] Formulae (5):

[0203] a) for the force tensor [T_ij]: $T_{k,l}\leftarrow T_{k,l}+\sum _{i=1}^4{K}_{s_i+4k,s_i+4l}^e;k=1,2,3;l=1,2,3$

[0204] s_i=local (m_i): local index of node N_i,

[0205] <<local>>: transfer function assigning a local index (lying between 1 and 4 in the case of a tetrahedron and corresponding to one of the four vertices of T) to node N_irepresented here by its mass m_i.

[0206] b) for the link tensor [T_ij]: $T_{k,l}\leftarrow T_{k,l}+\sum _{i,j=1}^4{K}_{s_i+4k,s_j+4l}^e;k=1,2,3;l=1,2,3$

[0207] s_i=local (m_i) and s_j=local (m_j)

[0208] Formula (6): $F_i^{\mathrm{tension}}=\sum _{\left(i,j\right)\in v}\{-K_{\left(i,j\right)}\left(l\left(t\right)-{\alpha }_0\right)\frac{l_{\mathrm{ij}}}{l\left(t\right)}+k_{\left(i,j\right)}^d\left(\frac{\partial l}{\partial t}\right)\frac{l_{\mathrm{ij}}}{l\left(t\right)})$

[0209] if i and j s

[0210] F_i^tension=0 otherwise

[0211] α₀<l₀, where l₀is the unloaded length of a spring whose stiffness is chosen according to the magnitude of the amplification desired.

[0212] Formula (7): $F_i^{\mathrm{int}}=\left[T_{\mathrm{ij}}\right]\Delta x_i+\sum _{j\in v}\left[T_{\mathrm{ij}}\right]\Delta x_j+F_i^{\mathrm{tension}}$

[0213] Formula (8):

[0214] F_i(t)=m_ia_i(t)=F_i^int(t)+F_i^ext(t)

[0215] Formula (9):

[0216] x_i+1x_i+Δt*f(x_i,t_i)

[0217] Formulae (10):

[0218] v_i+1=v_i+Δt*a_i

[0219] P_i+1=P_i+Δt*v_i

[0220] P_i: position of node N at the instant t_i.

[0221] Formulae (11):

[0222] δx₁=Δt·f(x_i,t_i)

[0223] δx₂=Δt·f(x_i½δx₁,t_i+½Δt)

[0224] δx₃=Δt·f(x_i½δx₂,t_i+½Δt)

[0225] δx₄=Δt·f(x_i+δx₃,t_i+Δt)

[0226] x_i+1=x_i+⅙δx₁+⅓δx₂+⅓δx₃+⅙δx₄

[0227] the velocity v(t+Δt)=v_i+1and the position P(t+Δt)=P_i+1are deduced from the acceleration a_ithrough the following equations:

[0228] δv₁=Δt·a_i

[0229] δP₁=Δt·v_i

[0230] δv₂=Δt·a′_i

[0231] δP₂=Δt·v′_i

[0232] δv₃=Δt·a″_i

[0233] δP₃=Δt·v″_i

[0234] δv₄=Δt·a′″_i

[0235] δP₄=Δt·v′″_i

[0236] v_i+1=v_i+⅙δv₁+⅓δv₂+⅓δv₃+⅙δv₄

[0237] P_i+1=P_i+⅙δP₁+⅓δP₂+⅓δP₃+⅙δP₄

[0238] a′_i, a″_iand a′″_i: accelerations calculated respectively from the <<intermediate>> positions δP₁, δP₂and δP₃.