# Patent

## Title
Method for determining a three-dimensional representation of an object using a sequence of cross-section images, computer program product, and corresponding method for analyzing an object and imaging system

## Abstract

The method comprises, for each cross-section image, determining the position of the object (O) in relation to the cross-section plane at the moment the cross-section image is captured, and determining a three-dimensional representation (V) of the object (O) using cross-section images (X0 . . . Xm) and positions of the object (O) in relation to the cross-section plane at each moment of image capture.
The positions of the object (O) in relation to the cross-section plane are determined using the cross-section images (X0 . . . Xm).

## Background
This is a 371 of PCT/FRO9/050513 filed Mar. 25, 2009 which has a priority of French no. 08 51984 filed Mar. 7, 2008, hereby incorporated by reference.
The present invention concerns a method for determining a three-dimensional representation of an object.
The invention specifically concerns the reconstruction of micro-objects in microscopic imaging systems.
Prior art document EP 1 413 911 A1 describes a method for determining a three-dimensional representation object using a sequence of cross-section images of the image in a cross-section plane, each cross-section image having been taken at a respective moment of image capture, of the type comprising:

    
    
        for each cross-section image, determining the position of the object in relation to the cross-section plane at the moment the cross-section image was captured,
        determining a three-dimensional representation of the object using cross-section images and positions of the object in relation to the cross-section plane at each moment of image capture.
    
    


In this document of the prior art, the object is a cell with dimensions in the vicinity of several micrometers. The cell is marked with a fluorescent compound and is placed in a receptacle of a microscopic imaging system. A confinement field is created in the receptacle to control the position of the cell, without said cell being crushed. The microscopic imaging system also comprises an optical microscope having a focal plane forming the cross-section plane. To obtain the sequence of cross-section images, the cell is made to rotate at a constant angular speed around an axis of the focal plane, and cross-section images are taken at a speed between 1 and 1000 acquisitions per second.
The document of the prior art also proposes to determine the three-dimensional representation of the cell using usual reconstruction techniques used in tomography. These usual techniques consist of determining the three-dimensional representation of an object from cross-section images and positions of the object in relation to the cross-section plane, the positions of the object in relation to the cross-section plane being determined by the knowledge of the adjustments of the radiography machine and/or physical sensors making it possible to know the position of the radiography machine.
The inventor author of the instant invention realized that radiography reconstruction techniques were yielding mediocre results, when applied to the preceding microscopic imaging system.
The invention aims to provide a reconstruction method adapted to a microscopic imaging system.
Thus, the invention concerns a method of the aforementioned type, characterized in that the positions of the object in relation to the cross-section plane are determined using the cross-section images.
According to other features of the invention:

    
    
        the object having a movement in relation to the cross-section plane made up, on one hand, of a regular movement defined by regular movement parameters and, on the other hand, a sequence of displacements disrupting the object between the moment the cross-section image is captured and the following moment of image capture, the disrupting displacements being disrupting translations, of variable direction and value, the determination of the positions of the object in relation to the cross-section plane comprises the determination of the sequence of disrupting displacements using the cross-section images,
        the regular movement being a rotational movement defined by the following movement parameters: a fixed axis of rotation and a constant angular speed, determining the positions of the object in relation to the cross-section plane comprises determining at least part of the regular movement parameters using the cross-section images,
        determining the positions of the object in relation to the cross-section plane comprises determining the absolute value of the angular speed from the determination of a period of revolution such that two cross-section images taken at respective moments separated from each other by a time substantially equal to a non-zero multiple of the period, are substantially similar,
        determining the period comprises determining an initial group of candidate periods; selecting, among the candidate periods, for each cross-section image, a first subset grouping together candidate periods having the highest levels of periodic similarity, the periodic similarity level characterizing, for each candidate period, the similarity level of the cross-section image with substantially periodic cross-section images, for the candidate period; determining a law of probability that a candidate period will be selected, the law of probability having the period as a parameter; selecting the period from the candidate periods of the initial group, as being that which maximizes the likelihood, according to the law of probability, knowing the selected candidate periods,
        determining positions of the object in relation to the cross-section plane comprises, the absolute value of the angular speed being known, for each cross-section image, determining cross-section images spatially neighboring on the cross-section image; for each cross-section image, calculating an orthogonal projection of the cross-section image on a support plane of each spatially neighboring cross-section image, the orthogonal projection being expressed as a function of an axis of rotation value and a disrupting translation sequence value; comparing the value of each of a set of pixels of the orthogonal projection with the value of said pixel of the spatially neighboring cross-section image on which the projection was done; determining the axis of rotation and sequence of disrupting translations whereof the values yield the values of close pixels during the comparison;
        the orthogonal projection of the cross-section image is done by calculating an affine transformation of the cross-section image, the affine transformation having a linear transformation component and a translation component depending on, respectively, the direction of the axis of rotation value and the speed of rotation on one hand, and, on the other hand, the axis of rotation value, the rotational speed and the disrupting translation sequence value,
        determining positions of the object in relation to the cross-section plane comprises, the absolute value of the angular speed being known; for each cross-section image, determining spatially neighboring cross-section images of the cross-section image; for each cross-section image, calculating an orthogonal projection of the cross-section image on a support plane of each spatially neighboring cross-section image, the orthogonal projection being expressed as a function of an axis of rotation value and a disrupting translation sequence value, the orthogonal projection being done by calculating an affine transformation of the cross-section image, the affine transformation having a linear transformation component and a translation component, the linear transformation component depending on the direction of the axis of rotation value and the speed of rotation; comparing the value of each of a set of pixels of the orthogonal projection with the value of that pixel of the spatially neighboring cross-section image on which the projection was done; determining the direction of the axis of rotation of the so-called reference translation components, whereof the values yield close pixel values during the comparison,
        determining an estimate of the position of the object in relation to the cross-section plane comprises, the direction of the axis of rotation and the reference translation components having been determined: the expression of affine transformation translation components as a function of a disrupting translation sequence value and a passage point value of the axis of rotation; determining the sequence of disrupting translations and a passage point of the axis of rotation whereof the values yield translation components close to reference translation components,
        determining positions of the object in relation to the cross-section plane comprises determining the direction of the axis of rotation, and determining the direction of the axis of rotation comprises determining a projection of the axis of rotation on the cross-section plane, using cross-section images,
        each cross-section image comprising a grid of pixels, identical for all of the cross-section images, the method comprises: for each pixel of the grid, determining the histogram of the values said pixel assumes in at least part of the sequence of cross-section images and the histogram of the values assumed by the point symmetrical to the pixel in relation to a projection value of the axis of rotation; determining the projection of the axis of rotation whereof the value yields close histograms between each pixel and its point of symmetry,
        the method comprises: determining, for at least one pixel of the grid, a symmetrization time, necessary for a point of the object to move from said pixel to the point of symmetry in relation to the projection of the axis of rotation; determining an angle between the axis of rotation and its projection on the cross-section plane, using the symmetrization time(s),
        determining positions of the object in relation to the cross-section plane comprises, the angular speed and axis of rotation being known: for each cross-section image, determining cross-section images spatially neighboring on the cross-section image; for each cross-section image, calculating an orthogonal projection of the cross-section image on a support plane of each spatially neighboring cross-section image, the orthogonal projection being expressed as a function of an axis of rotation value and a disrupting translation sequence value, the orthogonal projection being done by calculating an affine transformation of the cross-section image, the affine transformation having a linear transformation component and a translation component, the linear transformation component depending on the direction of the axis of rotation value and the angular speed of rotation; comparing the value of each of a set of pixels of the orthogonal projection with the value of said pixel of the spatially neighboring cross-section image on which the projection was done; determining, the axis of rotation being known, so-called reference translation components, the values of which yield close pixel values during the comparison; expressing affine transformation translation components as a function of a disrupting translation sequence value and a passage point value of the axis of rotation; determining the sequence of disrupting translations and a passage point of the axis of rotation whereof the values yield translation components close to the reference translation components.
    
    


The invention also concerns a computer program that, when implemented on a computer, carries out a method according to the invention.
The invention also concerns a method for analyzing an object to determine a three-dimensional representation thereof, comprising: introducing the object into an imaging system defining a focal plane forming a cross-section plane; configuring the imaging system to rotate the object around a fixed axis of rotation and with a fixed angular speed; taking a sequence of cross-section images of the object; determining a three-dimensional representation of the object from the sequence of cross-section images, using a method according to the invention.
According to other features of the invention:

    
    
        the method comprises: taking a sequence of cross-section images of the object, called adjustment cross-section images; determining the axis of rotation of the object, using the adjustment cross-section images,
        the method comprises adjusting the imaging system to move the cross-section plane in relation to the object such that the cross-section plane passes through a middle of the object,
        the adjustment of the imaging system comprises determining a barycenter of the object, and the middle of the object is the projection of the barycenter on the axis of rotation.
    
    


The invention also concerns an imaging system characterized in that it comprises means making it possible to obtain images in a focal plane, a receptacle for receiving an object, means for setting the object in motion, means for receiving cross-section images taken in the focal plane, adapted to implement a method according to the invention.



These features, as well as others, will appear upon reading the following description of preferred embodiments of the invention. This description is done according to the appended drawings, in which:
FIG. 1 shows an imaging system according to the invention,
FIG. 2 shows the steps of a method according to the invention for analyzing an object,
FIG. 3 shows the steps of a method, implemented in the analysis method of FIG. 2, for processing cross-section images of the object,
FIG. 4 shows the steps for determining parameters of the movement of the object in relation to a cross-section plane,
FIG. 5 illustrates the steps for determining, according to a first alternative, an axis of rotation of the object and a sequence of disrupting translations undergone by the object,
FIG. 6 illustrates the orthogonal projection of a cross-section image on a support plane of a spatially neighboring image,
FIG. 7 shows the steps for determining, according to a second alternative, an axis of rotation of the object and a sequence of disrupting translations that the object undergoes,
FIG. 8 shows the steps for determining, according to a third alternative, an axis of rotation of the object and a sequence of disrupting translations undergone by the object,
FIG. 9 shows a circle cutting the cross-section plane, that describes a point of the object around the axis of rotation,
FIG. 10 shows the steps for beginning a determination of a three-dimensional representation of the object,
FIG. 11 shows the steps for ending, according to a first alternative, the determination of a three-dimensional representation of the object,
FIG. 12 shows the steps for ending, according to a second alternative, the determination of a three-dimensional representation of the object,
FIG. 13 illustrates the steps of an alternative for determining a set of points, to determine the three-dimensional representation of the object of FIG. 10,
FIG. 14 shows the steps of an alternative for selecting groups of points, done in the determination of the set of points of FIG. 13,
FIG. 15 shows the steps, according to a third alternative, for determining a three-dimensional representation of the object, and
FIG. 16 shows the steps for adjusting the imaging system of FIG. 1.



The description that follows concerns the determination of a three-dimensional representation of a living cell. However, those skilled in the art will have no difficulty transposing this description to other types of objects, whether living or inert.
Furthermore, the mathematical relations indicated hereinafter are expressed in a fixed reference XYZ, attached to the cross-section plane P (which will be described later) and whereof the plane XY is combined with said cross-section plane P.
Moreover, in the present description, the term “value” does not refer only to a numerical value, but can also be a vector value (i.e. a specific vector) or a line value (i.e. a specific line), etc., depending on the nature of the mathematical object whereof the value is being considered.
DESCRIPTION OF THE IMAGING SYSTEM
In reference to FIG. 1, a microscopic imaging system 10 first comprises an optical microscope 12. The optical microscope 12 comprises a lens 14 defining a focal plane P. The optical microscope 12 also comprises a camera 16, for example a CCD camera making it possible to obtain images in the focal plane P.
The microscopic imaging system 10 also comprises a receptacle 18 for receiving an object O of microscopic size, such as a cell, in which a fluorescent material has been introduced.
Generally, it is possible to replace the fluorescent material with any marker adapted to the studied object, and capable of being detected by the imaging system used.
The receptacle 18 comprises a chamber 20 designed to contain a fluid microsystem comprising the object O. The chamber 20 is situated opposite the lens 14 of the optical microscope 12. The chamber 20 is defined by a support 24, side walls 26 and a window 28 covering the walls 26 so that the lens can observe the content of the chamber 20. The chamber 20 defines a U-shaped volume.
The side walls comprise microelectrodes 30 to create an electric field, the latter making it possible to position the object O.
The microscopic imaging system 10 also comprises a device 32 for illuminating the marker contained in the object O, such that each point o of the object O emits a luminosity O(o). The environment of the object O, on the contrary, emits a very weak or even non-existent luminosity.
The microscopic imaging system 10 also comprises a control unit 34, acting in particular on the microelectrodes 30 to set the object O in motion, and on the camera 16 to take a sequence of cross-section images X0 . . . Xm in the focal plane P (which thus forms a cross-section plane P of the volume U of the chamber 20, and in particular of the object O), at respective moments of image capture t0 . . . tm.
A reference XYZ is attached to the cross-section plane P. Another reference, called reference of the object O, is also attached to the object O. This other reference is preferably chosen such that it is combined with the reference XYZ at the initial moment t0. Moreover, the term “a point of the object” subsequently means, unless specifically otherwise indicated, a point in the reference of the object O.
Each cross-section image Xk extends over a support plane Pk attached to the object O, this support plane Pk being combined with the cross-section plane P at the moment tk the cross-section image Xk is captured.
Each cross-section image Xk comprise a grid G of pixels s. The grid G is the same for all of the cross-section images X0 . . . Xm. Each pixel s records a value Xk(s) of the marker's luminosity at the position of the pixel s in the volume U of the chamber 20. In the present description, this luminosity value is recorded in monochrome, in the form of a gray level.
Thus, when a point o of the object O is in the position of the pixel s, at the moment tk the cross-section image Xk is captured, the value Xk(s) of that pixel s depends in particular on the luminosity O(o) of the point o. When no point of the object O is in the position of the pixel s, the luminosity of the “vacuum” is recorded (in fact, that of the fluid comprising the object O). Thus, the background of the cross-section images X0 . . . Xm has a low gray level.
Moreover, the value Xk(s) of a pixel s also depends on a point spread function (PSF) introducing blurring. In general, the point spread function has an elongated shape perpendicular to the cross-section plane P.
The microscopic imaging system 10 also comprises a computer image processing device 36, connected to the control unit 30 to receive the cross-section images X0 . . . Xm. A reconstruction computer program 38 is installed on the computer device 36. The computer program 38 is designed to implement a reconstruction method intended to determine a three-dimensional representation V of the object O using the sequence of cross-section images X0 . . . Xm. The computer device 36 is capable of exporting the three-dimensional representation V in the form of a digital file and/or of displaying said three-dimensional representation V on a screen 40.
DESCRIPTION OF THE ANALYSIS METHOD
A method for analyzing an object O, carried out by the imaging system 10, is illustrated in FIG. 2. In reference to this FIG. 2, the method for analyzing an object O comprises a step 50 for introducing the object O into the chamber 20, then a step 52 for configuring the control unit 34 in order to cause the object O to rotate around a fixed axis of rotation L, with a fixed angular speed τ. The axis L is defined by a point u0 on the axis of rotation L—hereinafter called passage point u0—and a direction {right arrow over (a)} of the axis of rotation L, with unit standard: ∥{right arrow over (a)}∥=1. The axis of rotation L is not perpendicular to the cross-section plane P.
The analysis method also comprises a step 54 for acquiring at least one sequence of cross-section images X0 . . . Xm at respective moments of image capture t0 . . . tm, and a step 56 for processing the cross-section images X0 . . . Xm by the computer program 38.
In practice, the axis of rotation L is never exactly that adjusted by the control unit 34. The invention therefore proposes to determine the axis of rotation L using cross-section images X0 . . . Xm, then to determine a three-dimensional representation V using the determined axis of rotation L, rather than using the adjusted axis of rotation from the adjustment by the control unit 34.
Moreover, in practice, the movement of the object O is never completely rotational. The movement error in relation to the fixed axis rotation is represented by a sequence of disrupting translations T1 . . . Tm, each disrupting translation being undergone by the object O between two successive respective cross-section images Xk−1, Xk. The disrupting translations T1 . . . Tm have variable directions and values.
Thus, the position of a point o of the object O at a moment of image capture tk, starting from position u of point o at the preceding moment of image capture tk−1, is:

R{right arrow over (a)},τ(tk−tk−1)(u−u0)+u0+Tk,

where R{right arrow over (a)},τ(tk−tk−1) is the angle rotation matrix τ(tk−tk−1) around the axis of direction {right arrow over (a)}passing through the origin of the reference XYZ. It will be noted that the angle α rotation matrix R{right arrow over (a)},α around an axis of direction {right arrow over (a)} passing through the origin is given, according to Rodrigues' formula, by:

R{right arrow over (a)},α=I+sin α[{right arrow over (a)}]x+(1−cos α)[{right arrow over (a)}]x2,

where I is the 3 by 3 identity matrix, and





  
    
      [
      
        a
        ⇀
      
      ]
    
    x
  
  =
  
    [
    
      
        
          0
        
        
          
            -
            
              a
              3
            
          
        
        
          
            a
            2
          
        
      
      
        
          
            a
            3
          
        
        
          0
        
        
          
            -
            
              a
              1
            
          
        
      
      
        
          
            -
            
              a
              2
            
          
        
        
          
            a
            1
          
        
        
          0
        
      
    
    ]
  




with {right arrow over (a)}=(a1, a2, a3).

It will also be noted that the disrupting translation Tk does not depend on the position of u0, on the axis of rotation L.
Processing of the Acquired Cross-Section Images
In reference to FIG. 3, the processing step 56 carried out by the computer program 38 first comprises a step 100 for extension of each cross-section image X0 . . . Xm, during which point values between the pixels s are calculated, as well as point values outside the grid G. For example, the values of the points between the pixels s are calculated by interpolation or smoothing from values of the pixels s of the grid G, whereas the values of the points outside the grid G are put at a low gray level, for example 0. Subsequently, any point whatsoever of a cross-section image, pixel whereof the value is measured or point whereof the value is calculated, will be noted x, whereas a pixel strictly speaking of a cross-section image will be noted s. The value of a point x will be noted Xk (X), while the value of a pixel s will be noted Xk.
Step 56 carried out by the computer program 38 also comprises a step 200 for determining parameters of the movement of the object O in relation to the cross-section plane P. During this step for determining parameters of the movement 200, the regular movement parameters (angular speed τ, axis of rotation L) are determined, as well as the disrupting movement parameters (sequence of disrupting translations T1 . . . Tm).
The movement parameters (angular speed τ, axis of rotation L and sequence of disrupting movements T1 . . . Tm) determine the position of the object O at each moment of capture tk of a cross-section image Xk: the position at the capture moment tk of a point o located at position u at initial moment t0 is given by:

R{right arrow over (a)},τ(tk−t0)(u−u0)+u0+ Tk,

where Tk=R{right arrow over (a)},τ(tk−t1)T1+R{right arrow over (a)},τ(tk−t2)T2+ . . . +R{right arrow over (a)},τ(tk−tk−1)Tk−1+Tk, the accumulated translation since initial moment t0 until the moment tk when the cross-section image Xk is captured, with T0=0.

The position of the object O at each capture moment tk of cross-section images Xk determines the position of the cross-section plane P at capture moment tk in any reference whatsoever attached to the object O (and vice versa): the position in the reference of the chosen object O (which is combined with the fixed reference XYZ at the initial capture moment t0) of a pixel s of the cross-section plane P at the capture moment tk is:

R{right arrow over (a)},τ(tk−t0)(π3s−u0− Tk)+u0 

with R{right arrow over (a)},τ(tk−t0) the transposed matrix of R{right arrow over (a)},τ(tk−t0),





  
    
      π
      3
    
    =
    
      [
      
        
          
            1
          
          
            0
          
        
        
          
            0
          
          
            1
          
        
        
          
            0
          
          
            0
          
        
      
      ]
    
  
  ,




π3x is the position in the three-dimensional reference XYZ of a point x of a cross-section image, combining the support plane of the cross-section image with plane XY. The position of the cross-section plane P in any other reference attached to the object O at each capture moment tk follows from the position of the cross-section plane P in the reference of the object O chosen at capture moment tk and the relation between the reference of the chosen object and said other reference.

In the description that follows, the movement of the object O in relation to the cross-section P is expressed in reference XYZ attached to cross-section plane P. Of course, the movement of the object O in relation to the cross-section plane P could be expressed in another reference point, which would not necessarily be attached to the cross-section plane. In this case, determining the movement of the object O in relation to the cross-section plane P would also comprise determining the movement of the cross-section plane P in relation to said other reference point.
Step 56 carried out by the computer program 38 also comprises a step 400 for determining a three-dimensional representation V of the object O, from cross-section images X0 . . . Xm and movement parameters π,L, T1 . . . Tm.
Determining the Angular Speed in Absolute Value
The sign of the angular speed τ indicates the direction of rotation of the movement. The sign of the angular speed τ is positive if the rotation occurs in the positive direction in relation to direction {right arrow over (a)}, and negative if the rotation occurs in the negative direction in relation to direction {right arrow over (a)}. The sign of the angular speed τ is known, for example, according to the adjustment of the imaging system 10, once the direction {right arrow over (a)} of the axis of rotation is chosen.
If the sign of the angular rotation τ is not known, it can be chosen arbitrarily: in this case, the three-dimensional representation V of the object O will at worst be the three-dimensional representation V of the mirror image of the object O.
In reference to FIG. 4, step 200 for determining the parameters of the movement comprises a step 210 for determining the angular speed in absolute value |τ|.
This step 210 first comprises a step 212 for determining a period p>0 such that each pair of cross-section images Xk, Xk′ captured at respective moments tk, tk′ separated from each other by a time substantially equal to a (non-zero) multiple of the period p, are substantially similar. This period p is thus the period of revolution of the object O.
More precisely, step 212 for determining the period of revolution p comprises a step 214 for determining an initial group of positive candidate periods p1 . . . pn, and a step 216 for selecting the period of revolution p among the candidate periods p1 . . . pn of the initial group.
In the described example, the step 214 for determining the candidate periods p1 . . . pn consists of choosing the candidate periods p1 . . . pn. Preferably, the candidate periods p1 . . . pn are chosen uniformly spaced.
Simple Determination of the Period
Step 216 for selecting the period of revolution p comprises, in a simple alternative, determining the best period among the candidate periods p1 . . . pn, while maximizing the likelihood according to a chosen probability model.
Evolved Determination of the Period
However, to improve the reliability of the determination of the period of revolution p, one or several prior selections are made on the candidate periods p1 . . . pn, the probability model taking this or these selection(s) into account.
In reference to FIG. 4, step 216 for selecting the period of revolution p comprises two selection steps 218, 220 in order to obtain, on one hand, for each cross-section image Xk, a first respective subset pj(k,1), . . . , p(k,e) of candidate periods, and, on the other hand, for all of the cross-section images X0 . . . Xm, a second subset pj(1) . . . pj(l) of candidate periods.
Selection of the First Subsets
The first step 218 for selecting first subsets pj(k,1), . . . , pj(k,e) comprises, for each cross-section image Xk and for each candidate period pj, a step 222 for determining substantially periodic cross-section images Xk′ (according to the candidate period pj) at the cross-section image Xk.
Step 222 for determining substantially periodic cross-section images Xk′ comprises a step 224 for determining cross-section images Xk′ taken at capture moments tk′ separated from a capture moment tk of the cross-section image Xk by a time close to a (non zero) multiple of the candidate period pj. “Close” means that the difference between the close time and the (non zero) multiple of the candidate period pj, is in a time interval J including 0. Preferably, the time interval J is centered on 0. For example, J=[−ζ,ζ] with small ζ in relation to each candidate period pj, with ζ≦pjβ for the set of candidate periods pj. Alternatively, ζ varies as a function of the candidate period pj by choosing ζ=pj/10 for example for each candidate period pj.
The first step 218 for selecting first subsets {pj(k,1), . . . pj(k,e)} also comprises a step 226 for re-centering each substantially periodic cross-section image Xk′ in relation to the cross-section image Xk.
The re-centering step 226 first comprises a step 228 for selecting light pixels, in the cross-section image Xk and in the substantially periodic cross-section image Xk′. Preferably, the selected pixels are those whereof the gray level is above a predetermined threshold α, this threshold for example being the gray level q-quantile of the cross-section images X0 . . . Xm (which means that the proportion of the pixels of the cross-section images X0 . . . Xm that have a gray level less than or equal to α is substantially equal to q, and that of the pixels that have a gray level greater than α is substantially equal to 1−q), with q for example being between 60% and 95%.
The re-centering step 226 also comprises a step 230 for calculating, on one hand, a first center d(Xk) of the light points selected in the cross-section image Xk and, on the other hand, calculating a second center d(Xk′) of the light points selected in the cross-section image Xk′.
The center of an image X (Xk or Xk′) is given by:




  
    d
    ⁡
    
      (
      X
      )
    
  
  =
  
    
      
        ∑
        s
        
          
        
      
      ⁢
      
        
      
      ⁢
      
        
          1
          
            
              X
              ⁡
              
                (
                s
                )
              
            
            >
            α
          
        
        ⁢
        S
      
    
    
      
        ∑
        s
        
          
        
      
      ⁢
      
        
      
      ⁢
      
        1
        
          
            X
            ⁡
            
              (
              s
              )
            
          
          >
          α
        
      
    
  




where 1A>B is the indicator function: 1A>B=1 if A>B, and 1A>B=0 if not. The re-centering step 226 also comprises the determination 232 of a shift dk,k′ between the centers of the light points of the cross-section image Xk and of the substantially periodic cross-section image

Xk′:dk,k′=d(Xk)−d(Xk′).

The re-centering step 226 also comprises a step 234 for translation of the substantially periodic cross-section image Xk′ of the shift dk,k′ between the centers, in order to obtain a centered substantially periodic cross-section image, noted Trans(Xk′,dk,k′). The centered substantially periodic cross-section image Trans(Xk′, dk,k′) is calculated, for each pixel s, by:

Trans(Xk′,dk,k′)(s)=Xk′(s−dk,k′).

The first step 218 for selecting first subsets pj(k,1), . . . , pj(k,e)) also comprises a step 236 for determining a distance T (k,k′) between the cross-section image Xk and each centered substantially periodic cross-section image Trans(Xk′, dk,k′). Preferably, the distance T (k,k′) is given by the following relation:

∀0≦k,k′≦m,T(k,k′)=χ(Xk,Trans(Xk′,dk,k′)),

with χ a distance function, χ(X,Y) measures the difference between two images X and Y. The distance function χ is for example a quadratic distance of the gray levels of the pixels between the two images, given by:





  
    χ
    ⁡
    
      (
      
        X
        ,
        Y
      
      )
    
  
  =
  
    
      ∑
      s
      
        
      
    
    ⁢
    
      
    
    ⁢
    
      
        
          (
          
            
              X
              ⁡
              
                (
                s
                )
              
            
            -
            
              Y
              ⁡
              
                (
                s
                )
              
            
          
          )
        
        2
      
      .
    
  




The first step 218 for selecting first subsets pj(k,1), . . . , pj(k,e) also comprises a step 238 for determining a periodic similarity level sim(Xk, pj) of the cross-section image Xk using the distances T(k, k′).
The periodic similarity level sim(Xk, pj) characterizes the similarity level of the cross-section image Xk with the substantially periodic cross-section images Xk′, for the candidate period pj. Preferably, the step 238 for determining a periodic similarity level sim(Xk, pj) comprises a step for calculating the inverse similarity by the following relation:




  
    
      
        
          
            sim
            
              -
              1
            
          
          ⁡
          
            (
            
              
                X
                k
              
              ,
              
                p
                j
              
            
            )
          
        
        =
        
        ⁢
        
          {
          
            
              
                
                  
                    
                      
                        ∑
                        
                          r
                          ≠
                          0
                        
                        
                          
                        
                      
                      ⁢
                      
                        
                      
                      ⁢
                      
                        
                          v
                          ⁡
                          
                            (
                            
                              k
                              ,
                              j
                              ,
                              r
                            
                            )
                          
                        
                        /
                        
                          
                            ∑
                            
                              r
                              ≠
                              0
                            
                            
                              
                            
                          
                          ⁢
                          
                            
                          
                          ⁢
                          
                            h
                            ⁡
                            
                              (
                              
                                k
                                ,
                                j
                                ,
                                r
                              
                              )
                            
                          
                        
                      
                    
                    ⁢
                    
                      
                    
                  
                
                
                  
                    
                      si
                      ⁢
                      
                        
                          ∑
                          
                            r
                            ≠
                            0
                          
                          
                            
                          
                        
                        ⁢
                        
                          
                        
                        ⁢
                        
                          h
                          ⁡
                          
                            (
                            
                              k
                              ,
                              j
                              ,
                              r
                            
                            )
                          
                        
                      
                    
                    >
                    0
                  
                
              
              
                
                  ∞
                
                
                  
                    if
                    ⁢
                    
                      
                    
                    ⁢
                    not
                  
                
              
            
            .
          
        
      
    
  
  
    
      
        
          with
          ⁢
          
            :
          
        
        ⁢
        
          

        
        ⁢
        
          
            
              v
              ⁡
              
                (
                
                  k
                  ,
                  j
                  ,
                  r
                
                )
              
            
            =
            
            ⁢
            
              
                ∑
                
                  
                    
                      i
                      ⁢
                      
                        :
                      
                    
                    -
                    ζ
                  
                  ⁢
                  
                    
                  
                  ≤
                  
                    
                      t
                      i
                    
                    -
                    
                      rp
                      j
                    
                  
                  ≤
                  
                    
                  
                  ⁢
                  ζ
                
                
                  
                
              
              ⁢
              
                
              
              ⁢
              
                
                  w
                  ⁡
                  
                    (
                    
                      
                        t
                        i
                      
                      -
                      
                        rp
                        j
                      
                    
                    )
                  
                
                ⁢
                
                  T
                  ⁡
                  
                    (
                    
                      k
                      ,
                      i
                    
                    )
                  
                
              
            
          
          ,
          
            

          
          ⁢
          
            
              h
              ⁡
              
                (
                
                  k
                  ,
                  j
                  ,
                  r
                
                )
              
            
            =
            
            ⁢
            
              
                ∑
                
                  
                    
                      i
                      ⁢
                      
                        :
                      
                    
                    -
                    ζ
                  
                  ⁢
                  
                    
                  
                  ≤
                  
                    
                      t
                      i
                    
                    -
                    
                      rp
                      j
                    
                  
                  ≤
                  
                    
                  
                  ⁢
                  ζ
                
                
                  
                
              
              ⁢
              
                
              
              ⁢
              
                w
                ⁢
                
                  (
                  
                    
                      t
                      i
                    
                    -
                    
                      rp
                      j
                    
                  
                  )
                
              
            
          
          ,
        
      
    
  




where r is a non-zero integer, ζ is defined in the same way as for step 224 for determining substantially periodic cross-section images, and w positive weighting function defined over the interval J=[−ζ,ζ]. Preferably, w is symmetrical in relation to 0 (w(t)=w(−t)), and the high values of w are concentrated around 0. For example w(t)=exp(−ct2) with c a positive constant. The function w makes it possible to decrease the influence of the substantially periodic cross-section images Xk, that move away from the multiple of the candidate period pj.

The first step 218 of selecting first subsets pj(k,1), . . . , pj(k,e) also comprises a step 240 for selecting, among the candidate periods p1 . . . pn, for each cross-section image Xk, a first subset pj(k,1), . . . , pj(k,e) grouping together the candidate periods p1 . . . pn, having the highest similarity levels (i.e. the smallest values of sim−1(Xk, pj)). Preferably, a predetermined number e of candidate periods are selected. Also preferably, this number e is chosen between 1 and 15.
Selecting the Second Subset
The second step 220 for selecting the second subset pj(1) . . . pj(l) of candidate periods comprises, for each candidate period pj, a step 242 for calculating a number of appearances S(pj) of the candidate period pj, corresponding to the number of first subsets pj(k,1), . . . , pj(k,e) in which the candidate period pj appears.
The values of the number of appearances S for the candidate periods p1 . . . pn are for example given by the relation:




  
    ∀
    
      1
      ≤
      j
      ≤
      n
    
  
  ,
  
    
      S
      ⁡
      
        (
        
          p
          j
        
        )
      
    
    =
    
      
        ∑
        
          k
          =
          0
        
        m
      
      ⁢
      
        
      
      ⁢
      
        1
        
          
            p
            j
          
          ⁢
          ε
          ⁢
          
            
          
          ⁢
          
            {
            
              
                p
                
                  j
                  ⁡
                  
                    (
                    
                      k
                      ,
                      1
                    
                    )
                  
                
              
              ,
              ...
              ⁢
              
                
              
              ,
              
                p
                
                  j
                  ⁡
                  
                    (
                    
                      k
                      ,
                      e
                    
                    )
                  
                
              
            
            }
          
        
      
    
  
  ,




where 1AεB is the indicator function: 1AεB=1 if AεB, and 1AεB=0 if not.

The second step 220 for selecting the second subset pj(1) . . . pj(l) of candidate periods also comprises, for each candidate period pj, a step 244 for calculating a dispersion I of the values of the number of appearances S, around each multiple (greater than or equal to 1) of the candidate period pj. The dispersion I for a candidate period pj indicates whether high values of the appearance number S are concentrated (low dispersion) or dispersed (high dispersion) around each multiple (greater than or equal to 1) of the candidate period pj.

Preferably, the dispersion I of a candidate period pj is calculated such that, the further another candidate period pj, gets from the closest multiple (greater than or equal to 1) of the candidate period pj, the more the value of the appearance number S of said other candidate period pj′ contributes to the dispersion.
Preferably, the dispersion I of the candidate period p is given by the relation:




  
    ∀
    
      1
      ≤
      j
      ≤
      n
    
  
  ,
  
    
            ⁢
      
        (
        
          p
          j
        
        )
      
    
    =
    
      
        ∑
        
          i
          =
          1
        
        n
      
      ⁢
      
        
      
      ⁢
      
        
          
            
            
              
                
                  [
                  
                    
                      p
                      i
                    
                    
                      p
                      j
                    
                  
                  ]
                
                ⁢
                
                  p
                  j
                
              
              -
              
                p
                i
              
            
            
          
          
            p
            j
          
        
        ⁢
        
          S
          ⁡
          
            (
            
              p
              i
            
            )
          
        
      
    
  
  ,




where [x] is the integer closest to x.

The second step 220 for selecting the second subset pj(1) . . . pj(l) of candidate periods also comprises a step 246 for selecting, among the candidate periods pj, the second subset pj(1) . . . pj (l) of candidate periods having the lowest dispersions I. Preferably, a predetermined number of candidate periods are selected. Also preferably, this number is chosen between 4 and 40.
Selection of the Period from the First Subsets, Second Subset, and a Law of Probability.
The selection 216 of the period of revolution p also comprises a step 248 for determining a law of probability Pp describing the probability that a period pj will be chosen, during step 218, in at least one of the first sets pj(k,1), . . . , pj(k,e). The law of probability Pp is defined with the period of revolution p as parameter.
The selection 216 of the period of revolution p also comprises a step 249 for calculating a histogram h of the candidate periods selected in step 218. It has been noted that large values of the histogram h are generally concentrated around the true value of the period of revolution p and multiples of the true value of the period of revolution p.
Preferably, it is considered that the selections to obtain the candidate periods in step 218 are approximately independent, and that the choice of each of the candidate periods in step 218 is done randomly according to a law of probability Pp.
For this, the law of probability Pp, having as parameter the period of revolution p, is such that the shape of the probability function Pp can “fit” that of the histogram h to within about an expansion factor, making the value of the parameter p vary.
More precisely, it is considered that each period pj obtained in step 218 is selected either randomly, with a low probability β, among all of the candidate periods p1 . . . pn according to uniform distribution, or with probability 1−β by an alternative choice consisting of choosing a period around a multiple (greater than or equal to 1) of the period of revolution p. This alternative choice is made by first choosing a multiple ip (i an integer and i≦1) with a probability bi and then choosing pj around ip (ip already being chosen) with probability vip(pj).
In this case, the law of probability Pp is defined by the relation:




  
    ∀
    
      1
      ≤
      j
      ≤
      n
    
  
  ,
  
    
      
        P
        p
      
      ⁡
      
        (
        
          p
          j
        
        )
      
    
    =
    
      
        β
        n
      
      +
      
        
          (
          
            1
            -
            β
          
          )
        
        ⁢
        
          
            ∑
            
              i
              =
              1
            
            
              l
              ⁡
              
                (
                p
                )
              
            
          
          ⁢
          
            
          
          ⁢
          
            
              b
              i
            
            ⁢
            
              
                
                  v
                  _
                
                ip
              
              ⁡
              
                (
                
                  p
                  j
                
                )
              
            
          
        
      
    
  
  ,




where, l(p) is the number of multiples of p.

Law Pp is a mixture of the uniform distribution over all of the candidate periods p1 . . . pn, and laws vip with 1≦i≦l(p).
The probability vip(pj) is preferably the translation of a support function v of the value ip: vip(x)=v(x−ip). Preferably, the support function v is chosen finite and positive to define bi and vip. Also preferably, v is symmetrical around 0 (i.e. v(−x)=v(x)) and centered on 0. For example, one will take v(x)∝e−dx21|x|≦δ, or v(x)∝(δ−|x|)1|x|≦δ for d and δ constant positive givens.
Preferably l(p)=max{iεN*:{p1, . . . , pn}∩supp(vip)≠Ø}, with supp(vip)={x:vip(x)≠0}, and vip(pj)=vip(pj)/ci,




  
    
      b
      i
    
    =
    
      
        c
        i
      
      
        
          ∑
          
            k
            =
            1
          
          
            l
            ⁡
            
              (
              p
              )
            
          
        
        ⁢
        
          
        
        ⁢
        
          c
          k
        
      
    
  
  ,
  
    
      with
      ⁢
      
        
      
      ⁢
      
        c
        i
      
    
    =
    
      
        ∑
        
          k
          =
          1
        
        n
      
      ⁢
      
        
      
      ⁢
      
        
          
            v
            ip
          
          ⁡
          
            (
            
              p
              k
            
            )
          
        
        .
      
    
  




In practice β is chosen between 0 and 25%.

The selection 216 of the period of revolution p also comprises a step 250 for determining the period of revolution p, as being the period of the second subset of candidate periods pj(1) . . . pj(l) that maximizes the likelihood Vrais(p) (or, equivalently, the log-likelihood log Vrais(p)) associated with the preceding law of probability Pp, given the first subsets of selected periods pj(k,1), . . . , pj(k,e),0≧k≧m. The log-likelihood is given by:




  
    log
    ⁢
    
      
    
    ⁢
    
      Vrais
      ⁡
      
        (
        p
        )
      
    
  
  =
  
    
      ∑
      
        j
        =
        1
      
      n
    
    ⁢
    
      
    
    ⁢
    
      
        S
        ⁡
        
          (
          
            p
            j
          
          )
        
      
      ⁢
      log
      ⁢
      
        
      
      ⁢
      
        
          
            P
            p
          
          ⁡
          
            (
            
              p
              j
            
            )
          
        
        .
      
    
  




Alternatively, the step 212 for determining the period of revolution p does not comprise step 220 for determining the second subset, and, in step 250, the period of revolution p is determined as being the period of the set of candidate periods p1 . . . pn, maximizing the preceding likelihood Vrais(p).
Step 210 for determining the angular speed in absolute value |τ| then comprises a step 252 for calculating the absolute value of the angular speed |τ| from the period of revolution p:

|τ|=2π/p. 

Once the direction {right arrow over (a)} of the axis of rotation is chosen, the angular speed τ is determined from the known (or assumed) direction of rotation in relation to the direction {right arrow over (a)} of the axis of rotation, by the following relation:




  τ
  =
  
    {
    
      
        
          
            2
            ⁢
            
              
            
            ⁢
            
              π
              /
              p
            
          
        
        
          
            for
            ⁢
            
              
            
            ⁢
            positive
            ⁢
            
              
            
            ⁢
            rotations
          
        
      
      
        
          
            
              -
              2
            
            ⁢
            
              
            
            ⁢
            
              π
              /
              p
            
          
        
        
          
            for
            ⁢
            
              
            
            ⁢
            negative
            ⁢
            
              
            
            ⁢
            
              rotations
              .
            
          
        
      
    
  




Step 200 for determining the parameters of the movement also comprises a step 300 for determining the axis of rotation L and the sequence of disrupting translations T1 . . . Tm.
First Alternative to Determine the Axis of Rotation and the Sequence of Disrupting Translations
FIG. 5 illustrates a first alternative of step 300 for determining the axis of rotation L and the sequence of disrupting translations T1 . . . Tm.
Step 300 for determining the axis of rotation L and the sequence of disrupting translations T1 . . . Tm is carried out by determining an axis of rotation L and a sequence of disrupting translations T1 . . . Tm such that the value Xk(s) of each pixel s of each cross-section image Xk, is substantially equal, for each spatially neighboring cross-section image Xk′ (i.e. whereof the support plane Pk′ is spatially neighboring on the support plane Pk of the cross-section image Xk), to the value Xk′(x′) of a point x′ “close” to said pixel s, on the spatially neighboring cross-section image Xk′. The closest point x′ is determined according to a given proximity criterion.
More specifically, step 300 for determining the axis of rotation L and the sequence of disrupting translations T1 . . . Tm comprises a step 304 for determining, for each cross-section image Xk, spatially neighboring cross-section images Xk′ of said cross-section image Xk.
A cross-section image Xk′ is spatially neighboring on the cross-section image Xk, when the angle separating their respective support planes Pk, Pk′, is below a predetermined threshold value ΔI. The angle is preferably determined from the absolute value of the angular speed |τ| and moments tk, tk′ two cross-section images Xk, Xk′ are captured. This predetermined value ΔI is preferably taken less than or equal to 12 degrees. For a cross-section image Xk, the spatially neighboring cross-section images Xk′ are therefore the cross-section images verifying:

an integer q exists such that |τ|(tk−tk′)=2πq+r, with |r|≦Δ1.

The preceding condition on the angle separating the cross-section images Xk and Xk′ is equivalent to a condition on the time interval separating the capture of two cross-section images.
Step 300 for determining the axis of rotation L and the sequence of disrupting translations T1 . . . Tm also comprises a step 306 for choosing a proximity criterion.
In reference to FIG. 6, preferably, this proximity criterion is that the point x′ of a spatially neighboring cross-section image Xk′, the closest to a pixel s of the cross-section image Xk, is the position Projk,k′(s)=x′ of the orthogonal projection o′ of a point o of the object O on the support plane Pk of the cross-section image Xk, situated on the pixel s of the cross-section image Xk, on the support plane Pk′ of the spatially neighboring image Xk′.
This choice of this proximity criterion is particularly adapted in the case of a point spread function having an elongated shape perpendicular to the cross-section plane P.
Thus, returning to FIG. 5, step 300 for determining the axis of rotation L and the sequence of disrupting translations T1 . . . Tm also comprises a step 308 for calculating the orthogonal projection Fk,k′(Xk) of the cross-section image Xk on the support plane Pk′ of each spatially neighboring cross-section image Xk′ of the cross-section image Xk.
Preferably, step 308 for calculating the orthogonal projection Fk,k′(Xk) is done by calculating an affine transformation Aff(Qk,k′,vk,k′, Xk) of the cross-section image Xk, the affine transformation Aff(Qk,k′,vk,k′,Xk) having a linear transformation component Qk,k′({right arrow over (a)},τ) and a translation component vk,k′(u0,{right arrow over (a)},τ,T1 . . . Tm) depending, respectively, on one hand, on the direction {right arrow over (a)} of the axis of rotation value and the rotational speed τ and, on the other hand, the axis of rotation value L, the speed of rotation τ and the disrupting translation sequence value T1 . . . Tm. Thus:

Fk,k′(X)=Aff(Qk,k′({right arrow over (a)},τ),vk,k′(u0,{right arrow over (a)},τ,T1 . . . Tm),X),

with Aft (Q,v,X) the affine transform of X defined by:

∀x,Aff(Q,v,X)(x)=X(Q−1(x−v))

with Q a 2×2 invertible matrix and v a translation vector, and

Qi,j({right arrow over (a)},τ)=π2R{right arrow over (a)},τ(tj−ti)π3, vi,j(u0,{right arrow over (a)},τ,T1 . . . Tm)=(π2−π2R{right arrow over (a)},τ(tj−ti))u0+π2 Tj−π2R{right arrow over (a)},τ(tj−ti) Ti,

where





  
    
      π
      3
    
    =
    
      [
      
        
          
            1
          
          
            0
          
        
        
          
            0
          
          
            1
          
        
        
          
            0
          
          
            0
          
        
      
      ]
    
  
  ,




is the projection on plane XY,





  
    π
    2
  
  =
  
    [
    
      
        
          1
        
        
          0
        
        
          0
        
      
      
        
          0
        
        
          1
        
        
          0
        
      
    
    ]
  




π3x is the position in the three-dimensional reference XYZ of a point x of a cross-section image, combining the support plane of the cross-section image with the plane XY, and Tk=R{right arrow over (a)},τ(tk−t1)T1+R{right arrow over (a)},τ(tk−t2)T2+ . . . +R{right arrow over (a)},τ(tk−tk−1)Tk−1+Tk, the accumulated translation from the initial moment t0 to the moment tk the cross-section image Xk is taken with T0=0, and one notes that for each possible direction {right arrow over (a)} of the axis of rotation, the sign of the angular speed τ is determined from the known (or assumed) direction of rotation in relation to the direction {right arrow over (a)} of the axis of rotation, for example according to the adjustment of the imaging system 10, by the following relation:





  τ
  =
  
    {
    
      
        
          
            
            τ
            
          
        
        
          
            for
            ⁢
            
              
            
            ⁢
            positive
            ⁢
            
              
            
            ⁢
            rotations
          
        
      
      
        
          
            -
            
              
              τ
              
            
          
        
        
          
            
              for
              ⁢
              
                
              
              ⁢
              negative
              ⁢
              
                
              
              ⁢
              rotations
            
            ,
          
        
      
    
  




(if the direction of rotation is chosen arbitrarily: in this case, the three-dimensional representation V of the object O will at worst be the three-dimensional representation V of the mirror image of the object O).

Step 300 for determining the axis of rotation L and the sequence of disrupting translations T1 . . . Tm also comprises a step 310 for comparing the values Xk(s) of the pixels s of each cross-section image Xk with the values Xk′(x′) of the points x′ closest to the spatially neighboring cross-section images Xk′. Preferably, comparison step 310 is done by calculating a distance χ(Aff (Qk,k′, vk,k′, Xk), Xk′) between the affine transform Aff(Qk,k′, vk,k′, Xk) of the cross-section image Xk and the neighboring cross-section image Xk′.
Step 300 for determining the axis of rotation L and the sequence of disrupting translations T1 . . . Tm also comprises a step 312 for determining an axis of rotation L and a sequence of disrupting translations T1 . . . Tm whereof the values minimize a cost function E comprising a first part depending on the calculated distances, such that reducing the calculated distances causes the reduction of the cost function E.
Preferably, the cost function E also comprises a second portion that depends, for a sequence of disrupting translations T1 . . . Tm, on the amplitude of the translations of the sequence of disrupting translations T1 . . . Tm. This second portion thus comprises a regularization function Reg, which yields a low value when the translations of the sequence of disrupting translations T1 . . . Tm, have a low amplitude, and yields a high value in the contrary case. The regularization function Reg is preferably a quadratic norm, for example given by the relation:

Reg(T1 . . . Tm)=∥T∥M2≡TtMT, 

with M a symmetrical and positive definite matrix, and T=T1 . . . Tm, is written in the form of a column vector with the translations Tk one after the other. For example, M=I3m, matrix identity 3 m×3 m.

Preferably, the cost function E is given by the following relation:




  
    
      E
      ⁡
      
        (
        
          
            a
            ⇀
          
          ,
          
            x
            0
          
          ,
          
            
              T
              1
            
            ⁢
            
              
            
            ⁢
            ...
            ⁢
            
              
            
            ⁢
            
              T
              m
            
          
        
        )
      
    
    =
    
      
        
          ∑
          k
          
            
          
        
        ⁢
        
          
        
        ⁢
        
          
            ∑
            
              
                k
                ′
              
              ⁢
              ε
              ⁢
              
                
              
              ⁢
              
                N
                k
              
            
            
              
            
          
          ⁢
          
            
          
          ⁢
          
            χ
            ⁢
            
              
            
            ⁢
            
              (
              
                
                  
                    F
                    
                      k
                      ,
                      
                        k
                        ′
                      
                    
                  
                  ⁡
                  
                    (
                    
                      X
                      k
                    
                    )
                  
                
                ,
                
                  X
                  
                    k
                    ′
                  
                
              
              )
            
          
        
      
      +
      
        λ
        ⁢
        
          
        
        ⁢
        
          Reg
          ⁡
          
            (
            
              
                T
                1
              
              ⁢
              
                
              
              ⁢
              ...
              ⁢
              
                
              
              ⁢
              
                T
                m
              
            
            )
          
        
      
    
  
  ,




with Nk being the set of indices of the spatially neighboring cross-section images of the cross-section image Xk, χ(Fk,k′(Xk), Xk′) a measurement of the distance between, on one hand, the projection Fk,k′(Xk) of the cross-section image Xk on the support plane Pk′ of the spatially neighboring cross-section image Xk′ and, on the other hand, the spatially neighboring cross-section image Xk′, and λ≧0 is a compromise parameter between, on one hand, the distance between the projections Fk,k′(Xk) of the cross-section images Xk and the spatially neighboring cross-section images Xk′ of the cross-section images Xk, and, on the other hand, the regularization of translations T1 . . . Tm.

Second Alternative to Determine the Axis of Rotation and the Sequence of Disrupting Translations
A second alternative of step 300 for determining the axis of rotation L and the sequence of disrupting translations T1 . . . Tm, is illustrated in FIG. 7. This second alternative of step 300 repeats certain steps of the first alternative of step 300.
In reference to this figure, step 300 for determining the axis of rotation L and the sequence of disrupting translations T1 . . . Tm, comprises, the absolute value of the angular speed |τ| being known, a step 314 for determining the direction {right arrow over (a)} of the axis of rotation L, before a step 326 for determining a passage point u0 of the axis of rotation L and the sequence of disrupting translations T1 . . . Tm.
Determining the Direction of the Axis of Rotation
Step 314 for determining the direction {right arrow over (a)} of the axis of rotation L comprises a step 316, identical to the preceding step 304, for determining, for each cross-section image Xk, spatially neighboring cross-section images Xk′ of said cross-section image Xk.
Step 314 for determining the direction {right arrow over (a)} of the axis of rotation L also comprises a step 318, identical to the preceding step 306, for choosing a proximity criterion, according to which the point x′ of a spatially neighboring cross-section image Xk′, closest to a pixel s of the cross-section image Xk, is the position Projk,k′(s) of the orthogonal projection o′ of the point of the object o combined with the pixel s of the cross-section image Xk, on the support plane Pk′ of the spatially neighboring image Xk′.
Thus, step 314 for determining the direction {right arrow over (a)} of the axis of rotation L also comprises a step 320 for calculating the orthogonal projection Fk,k′(Xk) of the cross-section image Xk on the support plane Pk′ of each spatially neighboring cross-section image Xk′ of the cross-section image Xk.
Preferably, step 320 for calculating the orthogonal projection Fk,k′(Xk) is done by calculating an affine transformation Aff(Qk,k′,vk,k′,Xk) of the cross-section image Xk, the affine transformation Aff(Qk,k, vk,k′, Xk) having a linear transformation component Qk,k′({right arrow over (a)},τ) and a translation component vk,k′. The linear transformation component Qk,k′({right arrow over (a)},τ) is a function of an axis of rotation direction {right arrow over (a)} and the rotational speed τ, whereas, unlike the first alternative, the translation component vk,k′ is considered a variable, and is therefore not expressed as a function of the movement parameters of the object O. The family of translation vectors in the affine transformations Fk,k′ will hereinafter be noted v=(vk,k′)0≦k≦m,k′εNk. The orthogonal projection Fk,k′ (Xk) is therefore expressed as follows:

Fk,k′(X)=Aff(Qk,k′({right arrow over (a)},τ),vk,k′,X),

with Aff (Q,v,X) the affine transform of X defined by:

∀x,Aff(Q,v,X)(x)=X(Q−1(x−v))

with Q a 2×2 invertible matrix and v a translation vector, and

Qi,j({right arrow over (a)},τ)=π2R{right arrow over (a)},τ(ti−ti)π3,

where





  
    π
    2
  
  =
  
    [
    
      
        
          1
        
        
          0
        
        
          0
        
      
      
        
          0
        
        
          1
        
        
          0
        
      
    
    ]
  




is the projection on plane XY,





  
    
      π
      3
    
    =
    
      [
      
        
          
            1
          
          
            0
          
        
        
          
            0
          
          
            1
          
        
        
          
            0
          
          
            0
          
        
      
      ]
    
  
  ,




π3x is the position in the three-dimensional reference XYZ of a point x of a cross-section image, combining the support plane of the cross-section image with plane XY, and it is noted that for each possible direction {right arrow over (a)} of the axis of rotation, the sign of the angular speed τ is determined from the known (or assumed) direction of rotation in relation to the direction {right arrow over (a)} of the axis of rotation, for example according to the adjustment of the imaging system 10, by the following relation:





  τ
  =
  
    {
    
      
        
          
            
            τ
            
          
        
        
          
            for
            ⁢
            
              
            
            ⁢
            positive
            ⁢
            
              
            
            ⁢
            rotations
          
        
      
      
        
          
            -
            
              
              τ
              
            
          
        
        
          
            
              for
              ⁢
              
                
              
              ⁢
              negative
              ⁢
              
                
              
              ⁢
              rotations
            
            ,
          
        
      
    
  




(if the direction of rotation is chosen arbitrarily: in this case, the three-dimensional representation V of the object O will at worst be the three-dimensional representation V of the mirror image of the object O).

Step 314 for determining the direction {right arrow over (a)} of the axis of rotation L also comprises a step 322 for comparing the values Xk(s) of the pixels s of each cross-section image Xk with the values Xk′(x′) of the points x′ closest to the spatially neighboring cross-section images Xk′. Preferably, comparison step 322 is done by calculating a distance χ(Aff(Qk,k′,vk,k′,Xk),Xk′) between the affine transform Aff(Qk,k′, vk,k′, Xk) of the cross-section image Xk and the neighboring cross-section image Xk′.
Step 314 for determining the direction {right arrow over (a)} of the axis of rotation L also comprises a step 324 for minimizing a cost function U, in which the translation components vk,k′ are not expressed as a function of the parameters of the movement, but are left as variables of the cost function U.
Step 324 for minimizing the cost function U therefore amounts to finding the axis of rotation direction {right arrow over (a)} and the family of translation components vk,k′ that minimize the cost function U.
Preferably, the cost function U is given by the relation:




  
    
      U
      ⁡
      
        (
        
          
            a
            →
          
          ,
          v
        
        )
      
    
    =
    
      
        ∑
        k
        
          
        
      
      ⁢
      
        
      
      ⁢
      
        
          ∑
          
            
              k
              ′
            
            ⁢
            ε
            ⁢
            
              
            
            ⁢
            
              N
              k
            
          
          
            
          
        
        ⁢
        
          
        
        ⁢
        
          χ
          ⁢
          
            
          
          ⁢
          
            (
            
              
                Aff
                ⁡
                
                  (
                  
                    
                      
                        Q
                        
                          k
                          ,
                          
                            k
                            ′
                          
                        
                      
                      ⁡
                      
                        (
                        
                          
                            a
                            →
                          
                          ,
                          τ
                        
                        )
                      
                    
                    ,
                    
                      v
                      
                        k
                        ,
                        
                          k
                          ′
                        
                      
                    
                    ,
                    
                      X
                      k
                    
                  
                  )
                
              
              ,
              
                X
                
                  k
                  ′
                
              
            
            )
          
        
      
    
  
  ,




with Nk being the set of indices of the cross-section images neighboring on the cross-section image Xk, χ(Aff (Qk,k′({right arrow over (a)},τ),vk,k′,Xk), Xk′) a measurement of the distance between the affine transform and the cross-section image Xk and the spatially neighboring cross-section image Xk′.

In the preceding definition of the cost function U, the translation components vk,k′ are variables independent of u0,{right arrow over (a)},τ,T1 . . . Tm.
Determination of a Passage Point of the Axis of Rotation and the Sequence of Disrupting Translations
The direction {right arrow over (a)} of the axis of rotation L and of the translation components vk,k′, hereinafter called reference translation components vk,k′ having been determined, step 326 for determining a passage point u0 of the axis of rotation L and the sequence of disrupting translations T1 . . . Tm comprises a step 328 for expressing translation components vk,k′(u0,{right arrow over (a)},τ,T1 . . . Tm) as a function of a sequence of disrupting translations T1 . . . Tm and a passage point u0, the direction {right arrow over (a)} of the axis of rotation L and the rotational speed τ being known:

vk,k′(u0,{right arrow over (a)},τ,T1 . . . Tm)=(π2−π2R{right arrow over (a)},τ(tk′−tk))u0+π2 Tk′−π2R{right arrow over (a)},τ(tk′−tk) Tk,

where Tk=R{right arrow over (a)},τ(tk−t1)T1+R{right arrow over (a)},τ(tk−t2)T2+ . . . +R{right arrow over (a)},τ(tk−tk−1)Tk−1+Tk, the accumulated translation from the initial moment t0 until the moment tk when the cross-section image Xk is captured with T0=0.

Step 326 for determining a passage point u0 of the axis of rotation L and the sequence of disrupting translations T1 . . . Tm, also comprises a step 330 for determining the sequence of disrupting translations T1 . . . Tm and the passage point u0, such that the translation components vk,k′(u0,{right arrow over (a)},τ,T1 . . . Tm) approach the reference translation components vk,k′.
Preferably, step 330 for determining the sequence of disrupting translations T1 . . . Tm and the passage point u0 comprises a step 332 for minimizing a cost function K that comprises a first portion representing, for a sequence of disrupting translations T1 . . . Tm and a passage point u0, a distance between the reference translation components vk,k′ and the translation components vk,k′(u0,{right arrow over (a)},τ,T1 . . . Tm) expressed as a function of the sequence of disrupting translations T1 . . . Tm and the passage point u0.
Also preferably, the cost function K comprises a second portion representing the regularization Reg of the value T1 . . . Tm of disrupting translation sequences whereof the minimization reduces the amplitude of the translations T1 . . . Tm. The regularization function Reg is preferably a quadratic norm, for example given by the relation:

Reg(T1 . . . Tm)=∥T∥m2≡TtMT, 

with M a symmetrical and positive definite matrix, and T=T1 . . . Tm, is written in the form of a column vector with the translations Tk one after the other. For example, M=I3m, matrix identity 3 m×3 m.

Preferably, the cost function K is given by the relation:




  
    
      K
      ⁡
      
        (
        
          
            x
            0
          
          ,
          
            
              T
              1
            
            ⁢
            
              
            
            ⁢
            ...
            ⁢
            
              
            
            ⁢
            
              T
              m
            
          
        
        )
      
    
    =
    
      
        
          ∑
          k
          
            
          
        
        ⁢
        
          
        
        ⁢
        
          
            ∑
            
              
                k
                ′
              
              ⁢
              ε
              ⁢
              
                
              
              ⁢
              
                N
                k
              
            
            
              
            
          
          ⁢
          
            
          
          ⁢
          
            
              
              
                
                  v
                  
                    k
                    ,
                    
                      k
                      ′
                    
                  
                
                -
                
                  
                    v
                    
                      k
                      ,
                      
                        k
                        ′
                      
                    
                  
                  ⁡
                  
                    (
                    
                      
                        a
                        →
                      
                      ,
                      τ
                      ,
                      
                        
                          T
                          1
                        
                        ⁢
                        
                          
                        
                        ⁢
                        ...
                        ⁢
                        
                          
                        
                        ⁢
                        
                          T
                          m
                        
                      
                      ,
                      
                        x
                        0
                      
                    
                    )
                  
                
              
              
            
            2
          
        
      
      +
      
        α
        ⁢
        
          
        
        ⁢
        
          Reg
          ⁡
          
            (
            
              
                T
                1
              
              ⁢
              
                
              
              ⁢
              ...
              ⁢
              
                
              
              ⁢
              
                T
                m
              
            
            )
          
        
      
    
  
  ,




with α≧0 a compromise between the quadratic deviation between vk,k′ and vk,k′(u0,{right arrow over (a)},τ,T1 . . . Tm), and the regularization of T1 . . . Tm to control the amplitude of the translations T1 . . . Tm.

Third Alternative to Determine the Axis of Rotation and the Sequence of Disrupting Translations
A third alternative of step 300 for determining the axis of rotation L and the sequence of disrupting translations T1 . . . Tm is illustrated in FIG. 8.
In reference to FIG. 8, step 300 for determining the axis of rotation L and the sequence of disrupting translations T1 . . . Tm comprises a step 334 for determining the direction {right arrow over (a)} of the axis of rotation L, during which the disrupting translations T1 . . . Tm are considered to be negligible (which means that the movement of the object O is considered a stable rotation around the axis of rotation L),
Determining the Direction of the Axis of Rotation
Step 334 for determining the direction {right arrow over (a)} comprises a step 336 for determining a projection LXY of the axis of rotation L on the cross-section plane P, using cross-section images X0 . . . Xm.
Determining a Projection of the Axis of Rotation on the Cross-Section Plane
The method for determining the projection LXY uses the fact that each pixel s and its symmetrical point LXY(s) on the cross-section plane P in relation to the projection LXY, have histograms hs and hLXY(s) that are close. This is due to the fact that, during a rotation of the object O, the pixel s and its symmetrical point L(XY)(s) assume substantially the same values, just with a temporal shift corresponding to the time the point o of the object O takes to move from the pixel s to its symmetrical point LXY(s).
During step 336, the projection LXY of the axis of rotation L on the cross-section plane P is determined by choosing the line lXY of the cross-section plane P having the highest histogram symmetry level, i.e. such that the histograms hs and hlXY(s) of each pixel s and its symmetrical point lXY(s) are close.
Thus, step 336 for determining a projection LXY comprises, for each pixel s of the grid G, a step 338 for calculating a histogram hs of the gray levels assumed by said pixel s in at least part of the cross-section images X0, . . . , Xm. The histogram hs illustrates the count of the values of the pixel s, without taking into account the order in which those values appear in the sequence of cross-section images X0, . . . , Xm.
Preferably, the cross-section images part corresponds to the cross-section images X0, . . . , Xm′ (with m′≦m) acquired while the object O performs a whole number of rotations around the axis of rotation L. This translates to:




  
    
      t
      
        m
        ′
      
    
    -
    
      t
      0
    
  
  ≈
  
    
      
        2
        ⁢
        
          
        
        ⁢
        π
        ⁢
        
          
        
        ⁢
        r
      
      
        
        τ
        
      
    
    ⁢
    
      
    
    ⁢
    for
    ⁢
    
      
    
    ⁢
    an
    ⁢
    
      
    
    ⁢
    integer
    ⁢
    
      
    
    ⁢
    
      r
      .
    
  




Step 336 for determining a projection LXY also comprises a step 340 for determining, for each pixel s of the grid G and for a line lXY of the cross-section plane P, a histogram hlXY(s) of the gray levels assumed by the symmetrical point lXY(s) of said pixel s in the preceding part of the cross-section images.
Step 336 for determining a projection LXY also comprises a step 342 for determining distances Δ(hs,hlXY(s)) between the histogram hs of each pixel s and the histogram hlXY(s) of the symmetrical point lXY(s). The distance Δ is for example the Kolmogorov-Smirnov distance:




  
    
      Δ
      ⁡
      
        (
        
          h
          ,
          
            h
            ′
          
        
        )
      
    
    =
    
      
        sup
        x
      
      ⁢
      
        
        
          
            
              h
              _
            
            ⁡
            
              (
              x
              )
            
          
          -
          
            
              
                h
                _
              
              ′
            
            ⁡
            
              (
              x
              )
            
          
        
        
      
    
  
  ,
  
    
      with
      ⁢
      
        
      
      ⁢
      
        
          h
          _
        
        ⁡
        
          (
          x
          )
        
      
    
    =
    
      
        ∑
        
          y
          ≤
          x
        
        
          
        
      
      ⁢
      
        
      
      ⁢
      
        
          h
          ⁡
          
            (
            y
            )
          
        
        .
      
    
  




Step 336 for determining a projection LXY also comprises a step 344 for minimizing a cost function Ψ, which represents, for a line lXY of the cross-section plane P, the differences between the histogram hs of each pixel s and the histogram hlXY(s) of the symmetrical point lXY(s) in relation to the line lXY. Preferably, the cost function Ψ is given by the relation:




  
    Ψ
    ⁡
    
      (
      
        l
        XY
      
      )
    
  
  =
  
    
      ∑
      
        s
        ⁢
        
          
        
        ⁢
        ε
        ⁢
        
          
        
        ⁢
        G
      
      
        
      
    
    ⁢
    
      
    
    ⁢
    
      Δ
      ⁢
      
        
      
      ⁢
      
        
          (
          
            
              h
              s
            
            ,
            
              h
              
                
                  l
                  XY
                
                ⁡
                
                  (
                  s
                  )
                
              
            
          
          )
        
        .
      
    
  




LXY is determined as being the line lXY, on the plane XY that minimizes the cost function Ψ.
Determining the Angle Between the Axis of Rotation and its Projection
Step 334 for determining the direction {right arrow over (a)} also comprises a step 346 for determining the angle between the axis of rotation L and its projection LXY on the cross-section plane P. This angle between the axis of rotation L and its projection LXY is determined by calculating the angle φ between the axis Z of the reference XYZ (perpendicular to the cross-section plane P) and the axis of rotation L.
In reference to FIG. 9, step 346 for determining the angle φ uses the fact that points o of the object O describe, over time, a respective circle, centered substantially on the axis of rotation L, this circle cutting the cross-section plane P, on one hand, in a pixel s and, on the other, in its symmetrical point LXY (s) in relation to the projection LXY.
Returning to FIG. 8, step 346 for determining the angle φ comprises a step 348 for determining, for each pixel s of the grid G, a time t(s), called symmetrization time. The symmetrization time t(s) of a pixel s is the time necessary for a point o of the object O to move from said pixel s to the symmetrical point LXY (s) in relation to the projection LXY.
This translates to the fact that the value Xk(s) of a pixel s on a cross-section image Xk captured at a moment tk is substantially equal to the value Xk′(LXY (s)) of the symmetrical point LXY (s) on a cross-section image Xk′ captured at a moment tk′ substantially shifted from the symmetrization time t(s) in relation to the capture moment tk′.
Thus, the symmetrization time t(s) is preferably determined by determining the best temporal retiming between the two vectors (X0 (s), X1(s), . . . , Xm(s)) and (X0(LXY(s)), X1(LXY(s)), . . . , Xm(LXY (s))). Thus, the vector retimed by time μ of the vector (X0(LXY(s)), X1(LXY(s)), . . . , Xm(LXY (s))) is defined by v(μ)=(Xq(t0+μ)(LXY(s)),Xq(t1+μ)(LXY(s)), . . . , Xq(tj(μ)+μ)(LXY(s))), with j(μ) such that j(μ)+1 is the number of images captured between moment t0 and moment tm+μ, and q(t) such that capture moment tq(t) is the closest moment to the moment t among all of the moments cross-section images are captured. Preferably, the retiming error is given by

Err(μ)=κ((X0(s),X1(s), . . . ,Xj(μ)(s)),v(μ)),

with κ(u,v) a distance measurement between two vectors u and v, K is chosen for example as being a normalized quadratic distance:





  
    κ
    ⁢
    
      
    
    ⁢
    
      (
      
        u
        ,
        v
      
      )
    
  
  =
  
    
      1
      l
    
    ⁢
    
      
        ∑
        i
        
          
        
      
      ⁢
      
        
      
      ⁢
      
        
          (
          
            
              u
              i
            
            -
            
              v
              i
            
          
          )
        
        2
      
    
  




with l the length of the vectors. Thus t(s) is obtained by minimizing the retiming error Err(μ):





  
    
      t
      ⁡
      
        (
        s
        )
      
    
    =
    
      arg
      ⁢
      
        
      
      ⁢
      
        
          min
          
            0
            ≤
            μ
            <
            p
          
        
        ⁢
        
          Err
          ⁡
          
            (
            μ
            )
          
        
      
    
  
  ,




with p=2π/|τ| the period of revolution.

Step 346 for determining the angle φ also comprises a step 350 for determining the angle φ from symmetrization times t(s).
Step 350 for determining the angle φ using symmetrization times t(s) comprises a step 352 for determining a direction {right arrow over (b)} of the projection LXY with ∥{right arrow over (b)}∥=1, as well as a passage point y0 of the projection LXY. The projection LXY is determined by its direction and the passage point (y0, {right arrow over (b)}). The direction {right arrow over (b)} of the projection LXY is chosen as being, to within a multiplicative positive constant, the projection on the plane XY of the direction of the axis of rotation {right arrow over (a)} that one subsequently determines:

{right arrow over (b)}=πXY({right arrow over (a)})/∥πXY({right arrow over (a)})∥,

with





  
    π
    XY
  
  =
  
    [
    
      
        
          1
        
        
          0
        
        
          0
        
      
      
        
          0
        
        
          1
        
        
          0
        
      
      
        
          0
        
        
          0
        
        
          0
        
      
    
    ]
  




the projection matrix on the plane XY. Thus we assume that the direction of rotation of the object O in relation to the direction of the axis of rotation a that is subsequently determined, is known (or assumed), for example according to the adjustment of the imaging system 10, once the direction {right arrow over (b)} is chosen (such that {right arrow over (b)}=πXY({right arrow over (a)})/∥πXY({right arrow over (a)})∥). The angular speed of rotation τ is thus determined:





  τ
  =
  
    {
    
      
        
          
            
            τ
            
          
        
        
          
            for
            ⁢
            
              
            
            ⁢
            positive
            ⁢
            
              
            
            ⁢
            rotations
          
        
      
      
        
          
            -
            
              
              τ
              
            
          
        
        
          
            for
            ⁢
            
              
            
            ⁢
            negative
            ⁢
            
              
            
            ⁢
            
              rotations
              .
            
          
        
      
    
  




If the direction of rotation of the object is unknown, it is chosen arbitrarily: in this case, the three-dimensional representation V of the object O will at worst be the three-dimensional representation V of the mirror image of the object O.

Step 350 for determining the angle φ using the symmetrization times t(s) also comprises a step 354 for determining, for each pixel s, a distance z(s) between, on one hand, the middle




  
    s
    +
    
      
        L
        XY
      
      ⁡
      
        (
        s
        )
      
    
  
  2




of the segment between the pixel s and its symmetrical LXY(s) in relation to the projection LXY, and, on the other hand, the center c(s) of the circle around which a point o of the object O passing through the pixel s turns. This distance is positive if the center c(s) is above the plane XY, it is zero if c(s) is in the plane XY, and negative otherwise.

Preferably, the distance z(s) is given by the relation:




  
    
      z
      ⁡
      
        (
        s
        )
      
    
    =
    
      
        d
        ⁡
        
          (
          s
          )
        
      
      ⁢
      
        
          
          
            s
            -
            
              
                L
                XY
              
              ⁡
              
                (
                s
                )
              
            
          
          
        
        
          2
          ⁢
          
            
          
          ⁢
          
            tan
            ⁡
            
              (
              
                
                  
                    t
                    ⁡
                    
                      (
                      s
                      )
                    
                  
                  ⁢
                  τ
                
                2
              
              )
            
          
        
      
    
  
  ,




with t(s)τ the angle of the circle arc from s to the symmetrical point LXY(s) following the direction of rotation,





  
    d
    ⁡
    
      (
      s
      )
    
  
  =
  
    {
    
      
        
          1
        
        
          
            
              if
              ⁢
              
                
              
              ⁢
              s
            
            ∈
            
              P
              R
            
          
        
      
      
        
          
            -
            1
          
        
        
          
            
              
                if
                ⁢
                
                  
                
                ⁢
                s
              
              ∈
              
                P
                L
              
            
            ,
          
        
      
    
  




PR and PL being the two half-planes of the cross-section plane P separated by the projection LXY, PR and PL are defined by:

PL={sεP:<{right arrow over (b)}(s−y0),{right arrow over (e)}3>≧0}

PR={sεP:<{right arrow over (b)}(s−y0),{right arrow over (e)}3><0},

with  the vectorial product, <. , . > the scalar product, and





  
    
      e
      →
    
    3
  
  =
  
    [
    
      
        
          0
        
      
      
        
          0
        
      
      
        
          1
        
      
    
    ]
  




the direction of the axis Z of the reference. PR and PL do not depend on the position of y0 on the projection LXY.

In theory, there is an affine relation between the distance z(s) and the position on projection LXY of the axis of rotation, of the point




  
    s
    +
    
      
        L
        XY
      
      ⁡
      
        (
        s
        )
      
    
  
  2




(which is also the projection of son the projection LXY of the axis of rotation) in relation to the point y0:

z(s)=<s−y0,{right arrow over (b)}>cos φ+z0,

with z0≧0 if the projection of the point y0 on the axis of rotation L is above the plane XY or in the plane XY and z0<0 otherwise, and |z0| is the distance between y0 and the axis of rotation L.

Thus, preferably, step 350 for determining the angle φ using the symmetrization times t(s) also comprises a step 356 for determining the angle φ using the regression of the preceding affine relation:




  
    
      (
      
        φ
        ,
        
          z
          0
        
      
      )
    
    =
    
      arg
      ⁢
      
        
          min
          
            
              φ
              ^
            
            ,
            
              
                z
                ^
              
              0
            
          
        
        ⁢
        
          
            ∑
            
              s
              ∈
              
                G
                
                  δ
                  ,
                  
                    σ
                    min
                  
                
              
            
          
          ⁢
          
            
              (
              
                
                  
                    
                      z
                      ⁡
                      
                        (
                        s
                        )
                      
                    
                    -
                  
                  <
                  
                    s
                    -
                    
                      
                        y
                        ^
                      
                      0
                    
                  
                
                ,
                
                  
                    b
                    →
                  
                  >
                  
                    
                      cos
                      ⁢
                      
                        
                      
                      ⁢
                      
                        φ
                        ^
                      
                    
                    -
                    
                      
                        z
                        ^
                      
                      0
                    
                  
                
              
              )
            
            2
          
        
      
    
  
  ,




with Gδ,σmin the set of pixels s of the grid G whereof the distance in relation to the projection LXY of the axis of rotation exceeds a certain threshold δ and such that the empirical variance of the gray levels assumed by the pixel s in the sequence of cross-section images X0, . . . , Xm exceeds a certain threshold σmin2.

In practice, δ is chosen between 4 and 20 pixels, σmin2 is for example the q-quantile of the empirical variances of the gray levels, calculated for each pixel s of the grid G, assumed by the pixel s in the sequence of cross-section images X0, . . . , Xm (which means that the proportion of the pixels s of the grid G, where the calculated empirical variance of the gray levels is less than or equal to σmin2, is substantially equal to q), q is between 60% and 95% in general.
Limiting the preceding regression to the pixels in Gδ,σmin makes it possible to improve the robustness of the estimate by only using the symmetrization times for the pixels in Gδ,σmin which are generally more reliable.
Step 334 for determining the direction {right arrow over (a)} also comprises a step 358 for determining the direction {right arrow over (a)}, using the projection LXY and the angle φ. The direction {right arrow over (a)} of the axis of rotation L is for example given by the spherical coordinates vector (1,θ,φ), with θ the angle between the axis X and the direction {right arrow over (b)} of the projection LXY, such that

{right arrow over (b)}=πXY({right arrow over (a)})/∥πXY({right arrow over (a)})∥.

Optional Determination of the Passage Point of the Axis of Rotation
As an additional option, step 300 for determining the axis of rotation L and the sequence of disrupting translations T1 . . . Tm also comprises a step 360 for determining a passage point u0, of the axis of rotation L, using the symmetrization time(s) t(s). Preferably, the passage point u0, is chosen as the point of intersection between the axis of rotation L and the straight line perpendicular to plane XY and pasting through point y0. In this case, the passage point u0 is given by the following relation:




  
    
      u
      0
    
    =
    
      
        y
        0
      
      +
      
        
          
            z
            0
          
          
            sin
            ⁢
            
              
            
            ⁢
            φ
          
        
        ⁢
        
          
            e
            →
          
          3
        
      
    
  
  ,




Determining the Sequence of Disrupting Translations
Step 300 for determining the axis of rotation L and the sequence of disrupting translations T1 . . . Tm, also comprises a step 362 for determining the sequence of disrupting translations T1 . . . Tm, done in the same way as in step 330.
More precisely, step 362 for determining the sequence of disrupting translations T1 . . . Tm comprises a step 364 for determining, in the same way as in step 324, reference translation vectors vk,k′ of the affine transformations Fk,k′, the values of which minimize the cost function U, with the direction {right arrow over (a)} of the axis of rotation L and the angular speed τ known.
Owing to the knowledge of the direction {right arrow over (a)} of the axis of rotation, the minimization of U is greatly simplified, by independently minimizing each term χ(Aff(Qk,k′({right arrow over (a)},τ),vk,k′,Xk), Xk′) to determine each reference translation vector vk,k′ individually.
Step 362 for determining the sequence of disrupting translations T1 . . . Tm also comprises a step 366 for determining the sequence of disrupting translations T1 . . . Tm and the passage point u0 having the values that minimize the cost function K of step 332, knowing the reference translation vectors vk,k′.
Alternatively, when the optional step 360 is carried out, the minimization of the cost function K is simplified by using the passage point u0, determined in that step. Moreover, step 366 amounts to determining the sequence of disrupting translations T1 . . . Tm, having the value that minimizes the cost function K, the passage point u0 being known.
Determining a Three-Dimensional Representation of the Object
In reference to FIG. 10, step 400 for determining the three-dimensional representation V of the object O comprises a step 402 for determining a volume D, included in the volume U of the chamber 20, on which the object O will be “reconstructed”. The volume D comprises the object O (at least in part) and its close environment. For example, the volume D is a rhomb.
Step 400 for determining the three-dimensional representation V of the object O also comprises a step 404 for determining a set Ω of points u of the volume D and a value X(u) of each of those points u at a reference moment. The initial moment t0 will subsequently be chosen as reference moment. The set Ω points u in particular comprises points o of the object O in its position O0 at the initial moment t0.
It will be noted that, in a very simple embodiment, the set Ω of points u already forms a three-dimensional representation of the object, the latter being given by the constellation of the points of the set Ω.
Step 404 for determining the set of points Ω is carried out using the positions of the object O in relation to the cross-section plane P at each image capture moment t0 . . . tm, and the sequence of cross-section images X0 . . . Xm.
More precisely, step 404 for determining the set of points Ω comprises a step 406 for calculating, for each cross-section image Xk, the position Ck(s) of each pixel s of the cross-section image Xk at the initial moment t0, assuming that said pixel belongs to the object O. The calculation is done using parameters of the movement of the object O (angular speed τ, axis of rotation L, sequence of disrupting translations T1 . . . Tm) determined beforehand. The position Ck(s) of each point u of the set Ω is given by a respective original pixel s of an image Xk, by:

Ck(s)=R{right arrow over (a)},τ(tk−t0)t(π3s−u0− Tk)+u0 

with Tk=R{right arrow over (a)},τ(tk−t1)T1+R{right arrow over (a)},τ(tk−t2)T2+ . . . +R{right arrow over (a)},τ(tk−tk−1)Tk−1+Tk, with T0=0, R{right arrow over (a)},τ(tk−t0)t is the transposed matrix of R{right arrow over (a)},τ(tk−t0). The set Ω is made up of points u of the volume D of position Ck(s) for each pixel s and each cross-section image Xk. We thus say that s and Xk are the pixel and the original image of the point u, respectively. It is noted that Ck(s) is also the position in the reference of the chosen object O of the pixel s of the cross-section plane at the moment tk the cross-section image Xk is captured.

Each of the points u of the set Ω is associated with the value Xk(s) of the original pixel s of the original cross-section image Xk: X(u)=Xk(s).
It will be noted that the set Ω can thus comprise a same point u several times, each time associated with a respective value, these different values coming from different cross-section images.
Step 400 for determining the three-dimensional representation V of the object O also comprises a step 408 for choosing a three-dimensional representation function Vβ parameterizable with parameters β, and an operation Op giving, from the three-dimensional representation function Vβ, an estimation function {tilde over (X)}=Op(Vβ) of the value of each point u of the set Ω.
Once the parameters β are determined, the three-dimensional representation V is given by the three-dimensional representation function Vβ, preferably at all points of the volume D.
First Alternative of Three-Dimensional Representation Function Vβ
In reference to FIG. 11, in a first alternative, the three-dimensional representation function Vβ, chosen during step 408, is a decomposition into B-spline functions having degree r with nodes w spaced equally in the space:

w=b+(a1k1,a2k2,a3k3),

with bεR3, and a1, a2 and a3 respectively the sampling interval in directions X, Y, and Z, and k1, k2, k3 integers. Each parameter of the three-dimensional representation function Vβ is associated with a respective node.

The three-dimensional representation function Vβ is then written:




  
    
      ∀
      u
    
    =
    
      
        (
        
          
            u
            1
          
          ,
          
            u
            2
          
          ,
          
            u
            3
          
        
        )
      
      ∈
      D
    
  
  ,
  
    

  
  ⁢
  
    
      V
      ⁡
      
        (
        u
        )
      
    
    =
    
      
        ∑
        
          w
          =
          
            
              (
              
                
                  w
                  1
                
                ,
                
                  w
                  2
                
                ,
                
                  w
                  3
                
              
              )
            
            ∈
            W
          
        
      
      ⁢
      
        
          η
          ⁡
          
            (
            
              
                
                  u
                  1
                
                -
                
                  w
                  1
                
              
              
                a
                1
              
            
            )
          
        
        ⁢
        
          η
          ⁡
          
            (
            
              
                
                  u
                  2
                
                -
                
                  w
                  2
                
              
              
                a
                2
              
            
            )
          
        
        ⁢
        
          η
          ⁡
          
            (
            
              
                
                  u
                  3
                
                -
                
                  w
                  3
                
              
              
                a
                3
              
            
            )
          
        
        ⁢
        
          β
          ⁡
          
            (
            w
            )
          
        
      
    
  
  ,




with η the central B-spline function having degrees r defined on a set of real numbers R, W the set of nodes in the volume D. η is for example the indicator function on the interval





  
    [
    
      
        -
        
          1
          2
        
      
      ,
      
        1
        2
      
    
    ]
  
  ,




convoluted r times with itself:





  
    η
    ⁡
    
      (
      x
      )
    
  
  =
  
    
      
        
          1
          
            [
            
              
                -
                
                  1
                  2
                
              
              ,
              
                1
                2
              
            
            ]
          
        
        *
        
          1
          
            [
            
              
                -
                
                  1
                  2
                
              
              ,
              
                1
                2
              
            
            ]
          
        
        *
        ...
        *
        
          1
          
            [
            
              
                -
                
                  1
                  2
                
              
              ,
              
                1
                2
              
            
            ]
          
        
      
      
        }
        
          r
          +
          
            1
            ⁢
            fois
          
        
      
    
    ⁢
    
      
        (
        x
        )
      
      .
    
  




In particular, if r=3 (cubic B-spline function):





  
    η
    ⁡
    
      (
      x
      )
    
  
  =
  
    {
    
      
        
          
            
              
                2
                3
              
              -
              
                x
                2
              
              +
              
                
                  
                    
                    x
                    
                  
                  3
                
                2
              
            
            ,
          
        
        
          
            
              if
              ⁢
              
                
              
              ⁢
              0
            
            ≤
            
              
              x
              
            
            <
            1
          
        
      
      
        
          
            
              
                
                  (
                  
                    2
                    -
                    
                      
                      x
                      
                    
                  
                  )
                
                3
              
              6
            
            ,
          
        
        
          
            
              if
              ⁢
              
                
              
              ⁢
              1
            
            ≤
            
              
              x
              
            
            <
            2
          
        
      
      
        
          
            0
            ,
          
        
        
          
            
              if
              ⁢
              
                
              
              ⁢
              
                
                x
                
              
            
            ≥
            2.
          
        
      
    
  




In this first alternative, the operation Op is chosen as the identity function, such that the estimation function {tilde over (X)} of the value of each point u of the set Ω is equal to the three-dimensional representation function Vβ: {tilde over (X)}(u)=Vβ(u).
With this choice of three-dimensional representation function Vβ, step 400 for determining the three-dimensional representation V of the object O also comprises a step 410 for dividing the volume D into a plurality of sub-volumes D, disconnected from each other. Alternatively, the edges of the sub-volumes D, overlap. In this way, the nodes w are also divided into groups {w}i each associated with a respective sub-volume Di. Likewise, the points of the set Ω are divided into groups Ωi each associated with a respective sub-volume Di.
More precisely, each group {w}i comprises the nodes situated in the respective sub-volume Di and the parameters {β}i of those nodes {w}i are thus associated with the sub-volume Di. Likewise, the points of each group Ωi are the points of the set Ω situated in the sub-volume Di.
Step 400 for determining the three-dimensional representation V of the object O also comprises a step 412 for determining parameters β, such that, for each point u of the set Ω, the estimation {tilde over (X)}(u) of la value of the point u substantially gives the value X(u) of the point u.
More precisely, step 412 for determining the parameters β comprises, for each sub-volume Di, a step 414 for determining parameters {β}i associated with that sub-volume Di, such that, for each point u of the group Ωi and of the groups directly adjacent to the group Ωi, the estimation {tilde over (X)}(u) of the value of the point u substantially gives the value X(u) of the point u, the parameters {β}j≠i associated with the other subset Dj≠i being set at a given value.
Step 414 for determining the parameters {β}i is carried out several times, iteratively, each time sweeping all of the sub-volumes Di: during the first iteration, one successively determines each group of parameters {β}i (to determine a parameter group {β}i, the given value of each of the other parameter groups {β}j≠i is set at a predetermined value, for example zero); during later iterations, one successively determines each group of parameters {β}i (to determine a parameter group {β}i, the given value of each of the other parameter groups {β}j≠i is the last result determined previously).
Preferably, the parameters {β}i are determined by minimizing the following cost function:




  
    
      U
      ⁡
      
        (
        β
        )
      
    
    =
    
      
        
          ∑
          
            u
            ∈
            Ω
          
        
        ⁢
        
          
            (
            
              
                
                  X
                  ~
                
                ⁡
                
                  (
                  u
                  )
                
              
              -
              
                X
                ⁡
                
                  (
                  u
                  )
                
              
            
            )
          
          2
        
      
      +
      
        
          λβ
          t
        
        ⁢
        A
        ⁢
        
          
        
        ⁢
        β
      
    
  
  ,




A is a positive definite symmetrical matrix, or more generally semi positive definite, β′Aβ is a measurement of the quadratic irregularity, and λ>0 is a compromise between the appropriateness of the three-dimensional representation function and the regularity.

Second Alternative of Three-Dimensional Representation Function Vβ
In reference to FIG. 12, the three-dimensional representation function Vβ, chosen during step 408, is written in the form of a decomposition into radial basis functions (RBF) φ with the nodes w:




  
    ∀
    
      u
      ∈
      D
    
  
  ,
  
    
      
        V
        β
      
      ⁡
      
        (
        u
        )
      
    
    =
    
      
        ∑
        
          w
          ∈
          W
        
      
      ⁢
      
        
          φ
          ⁡
          
            (
            
              u
              -
              w
            
            )
          
        
        ⁢
        
          β
          ⁡
          
            (
            w
            )
          
        
      
    
  
  ,




with W the set of nodes in the volume D.

The function φ(u−w) depends on the distance between the point u and the node w, but not the direction between the point u and the node w. For example, φ(x)=exp(−c∥x∥2) or φ(x)=η(∥x∥) with η the cubic central B-spline function.
Moreover, the operation Op gives an estimation function {tilde over (X)}=Op(Vβ, ƒR) of the value of each point u of the set Ω using the three-dimensional representation function Vβ and a point spread function ƒR. In the described example, the operation Op is a convolution of the three-dimensional representation function Vβ with the following point spread function ƒR: ƒR(u)=ƒ(Ru), with ƒ the point spread function without rotation, which is known (for example, given by the builder of the imaging system, or determined experimentally), ƒR the point spread function for the rotation R, and Ru the point resulting from the rotation of the point u by the rotation R.
The point spread function ƒR depends on the rotation R between the position of the object O at the moment tk the respective cross-section image Xk is taken associated with the point u, and the position of the object O at the reference moment t0: R=R{right arrow over (a)},τ(tk−t0). Thus, for each point u of the set Ω,

{tilde over (X)}(u)=Op(Vβ,ƒR)(u)=(Vβ*ƒR)(Ck(s)),

with * the convolution operation, the pixel s and the cross-section image Xk (of the capture moment tk) respectively being the pixel and the original image of the point u.

Because of the choice of radial basis functions φ, the property is obtained that, for each point u of the volume D:

Op(φ,ƒR)(u)=Op(φ,ƒ)(Ru)

or:

(φ*ƒR)(u)=(φ*ƒ)(Ru)

with * the convolution operation, R any rotation whatsoever.

Thus, the step 400 for determining the three-dimensional representation V of the object O also comprises a step 416 for determining the parameters β of the three-dimensional representation function Vβ by minimizing the following cost function:




  
    
      E
      ⁡
      
        (
        β
        )
      
    
    =
    
      
        
          ∑
          
            u
            ∈
            Ω
          
        
        ⁢
        
          
            (
            
              
                X
                ⁡
                
                  (
                  u
                  )
                
              
              -
              
                
                  ∑
                  
                    w
                    ∈
                    W
                  
                
                ⁢
                
                  
                    (
                    
                      φ
                      *
                      f
                    
                    )
                  
                  ⁢
                  
                    (
                    
                      
                        
                          π
                          3
                        
                        ⁢
                        s
                      
                      -
                      
                        ς
                        i
                      
                      -
                      
                        
                          R
                          i
                        
                        ⁢
                        w
                      
                    
                    )
                  
                  ⁢
                  
                    β
                    ⁡
                    
                      (
                      w
                      )
                    
                  
                
              
            
            )
          
          2
        
      
      +
      
        
          λβ
          t
        
        ⁢
        A
        ⁢
        
          
        
        ⁢
        β
      
    
  
  ,




where s and the cross-section image Xi (of the image capture moment ti) are respectively the pixel and the original image of the point u, Ri=R{right arrow over (a)},τ(ti−t0) and ζi= Ti+u0−R{right arrow over (a)},τ(ti−t0)u0.

Advantageously, this minimization is done by calculating a unique convolution (the convolution φ*f) and solving the linear system, which follows from the calculation of the convolution φ*f, on the parameters β. The parameters β are thus easily determined.
Alternative for Determining the Set Ω
FIG. 13 illustrates an alternative embodiment of step 404 for determining a set Ω of points u of the volume D and a value X(u) for each of those points u at a reference moment, to determine the set of points Ω using several sequences of cross-section images, noted Sl with l=1 . . . I. Each sequence Sl is made up of cross-section images X0l, X1l, . . . , Xmll, at respective image capture moments t0l, t1l, . . . tmll.
In this alternative, step 404 for determining the set of points Ω comprises a step 420 for determining a three-dimensional representation function Vl, on a respective sub-volume Dl, for each sequence Sl. This step 420 is for example carried out following the preceding steps 408 to 416 for each sequence Sl.
Each three-dimensional representation function Vl gives a representation of the object O in a respective position, noted Ol.
Step 404 for determining the set of points Ω also comprises a step 422 for discretizing each sub-volume Di according to a three-dimensional grid G3. One thus obtains, for each sequence Sl, a discretized sub-volume {tilde over (D)}l, grouping together the points of the sub-volume Dl situated on the three-dimensional grid G3. One thus obtains:

∀uε{tilde over (D)}l,{tilde over (V)}l(u)=Vl(u).

Preferably, the three-dimensional grid G3 has a pitch less than or equal to that of the grid G of the cross-section plane P.
Step 404 for determining the set of points Ω also comprises a step 424 for determining, for each sequence Sl, a rotation Ql and a translation hl making it possible to substantially place all of the positions Ol of the representations of the object O in a same reference position.
Preferably, the reference position is that of one of the sequences Sl. Subsequently, the position Ol of the object of the first sequence Sl will be the reference position. Thus, a point o of the object O of the position u in Ol is at position Qlu+hl in Ol, with l≠1.
Step 424 for determining, for each sequence Sl, a rotation Ql and a translation hl comprises a step 426 for determining a quantile level q, such that the light points of the object O substantially have values {tilde over (V)}l(u) greater than or equal to the q-quantile ρl(q), for each discretized subset {tilde over (D)}l. For example, the level of quantile q is between 60% and 95%.
Step 424 for determining, for each sequence Sl, a rotation Ql and a translation hl also comprises a step 428 for selecting at least three groups g1 . . . gk, preferably four or more, of points of the discretized subset {tilde over (D)}l, according to a selection criterion. The selection criterion is the same for all of the sequences S1. The selection criterion is applicable in particular on the value of the points of the discretized subset {tilde over (D)}l.
By thus applying the same selection criterion for all of the sequences Sl, it is possible to obtain substantially the same points of the object O for all of the sequences Sl, even if these points do not have the same position for the different sequences Sl.
In a first alternative illustrated in FIG. 13, step 428 for selecting groups g1 . . . gk, comprises a step 430 for selecting a same number n of the brightest points (having the highest values). Preferably, the number n is the integer closest to min(qn1, . . . , qnl) (q is expressed as a percentage), with n1, . . . , nl the numbers of points in the discretized subsets {tilde over (D)}1, . . . {tilde over (D)}l.
Step 428 for selecting groups g1 . . . gk also comprises a step 432 for ranking, for each sequence Sl, the n brightest points (selected beforehand) by their values, in decreasing order, for example.
Step 428 for selecting groups g1 . . . gk also comprises a step 434 for dividing the n ranked brightest points into k groups g1 . . . gk of substantially equal size: the points in group g1 are brighter than those in group g2, which are brighter than those in g3, etc.
In a second alternative illustrated in FIG. 14, step 428 for selecting the groups g1 . . . gk, comprises a step 436 for calculating, for each sequence Sl, the barycenter bl of the points u of the subset {tilde over (D)}l, the values of which are greater than the q-quantile ρl(q). The set of these points will be noted {tilde over (D)}l,q hereinafter. To calculate the barycenter bl, all of the points are weighted the same.
Step 428 for selecting the groups g1 . . . gk, also comprises a step 437 for determining the highest distance between the barycenter and the points of all of the sets {tilde over (D)}l,q. This distance will subsequently be called radius r.
Preferably, the radius r is given by:

r=min(max{∥u−b1∥:uε{tilde over (D)}1,q}, . . . ,max{∥u−bl∥:uεDl,q})

Step 428 for selecting the groups g1 . . . gk, also comprises a step 438 for dividing the radius r (i.e. the segment [0,r]) into k segments of equal sizes (segi=[(k−i)r/k,(k−i+1)r/k], with 1≦i≦k), each group g1 . . . gk being associated with a respective segment.
More precisely, for each sequence Sl, each group gi comprises the points u of {tilde over (D)}l,q for which the distance to the barycenter bl is comprised in the associated segment segi.
Returning to FIG. 13, step 424 for determining, for each sequence Sl, a rotation Ql and a translation hl also comprises a step 440 for determining, for each sequence Sl, the barycenter ωl,i of the points u of each of the groups gi. To calculate the barycenters ωl,i, the points are weighted identically.
Step 424 for determining, for each sequence Sl, a rotation Ql and a translation hl also comprises a step 442 for determining the rotation Ql and the translation hl using the barycenters ωl,i.
Preferably, the rotations Ql and translations hl are determined by minimizing a cost function:




  
    
      (
      
        
          Q
          l
        
        ,
        
          h
          l
        
      
      )
    
    =
    
      arg
      ⁢
      
        
      
      ⁢
      
        
          min
          
            
              Q
              ∈
              
                O
                ⁡
                
                  (
                  3
                  )
                
              
            
            ,
            
              h
              ∈
              
                R
                3
              
            
          
        
        ⁢
        
          
            ∑
            
              i
              =
              1
            
            k
          
          ⁢
          
            
              
              
                
                  Q
                  ⁢
                  
                    
                  
                  ⁢
                  
                    ω
                    
                      1
                      ,
                      i
                    
                  
                
                +
                h
                -
                
                  ω
                  
                    l
                    ,
                    i
                  
                
              
              
            
            2
          
        
      
    
  
  ,




with O(3) the set of 3-by-3 orthogonal matrices.

Preferably, the preceding minimization solution is obtained by calculating:




  {
  
    
      
        
          
            
              
                
                  Q
                  l
                
                =
                
                  
                    P
                    l
                  
                  ⁢
                  
                    P
                    1
                    t
                  
                
              
            
          
          
            
              
                
                  h
                  l
                
                =
                
                  
                    
                      ω
                      _
                    
                    l
                  
                  -
                  
                    
                      Q
                      l
                    
                    ⁢
                    
                      
                        
                          ω
                          _
                        
                        1
                      
                      .
                    
                  
                
              
            
          
        
        ⁢
        
          
        
        ⁢
        with
        ⁢
        
          
        
        ⁢
        
          
            ω
            _
          
          1
        
      
      =
      
        
          
            1
            k
          
          ⁢
          
            
              ∑
              
                i
                =
                1
              
              k
            
            ⁢
            
              
                ω
                
                  1
                  ,
                  i
                
              
              ⁢
              
                
              
              ⁢
              and
              ⁢
              
                
              
              ⁢
              
                
                  ω
                  _
                
                l
              
            
          
        
        =
        
          
            1
            k
          
          ⁢
          
            
              ∑
              
                i
                =
                1
              
              k
            
            ⁢
            
              ω
              
                l
                ,
                i
              
            
          
        
      
    
    ,
  




and P1 and Pl obtained through singular value decomposition (SVD) of the matrix M=Σi=1k(ω1,i− ω1)(ωl,i− ωl)t, i.e. M=P1ΛPlt with P1 and Pl 3-by-3 orthogonal matrices, and Λ a diagonal matrix with non-negative numbers on the diagonal.

Step 404 for determining the set of points Ω, also comprises a step 444 for determining points u of Ω and a value X(u) of each of those points u, using cross-section image sequences Sl with l=1 . . . I, the rotation matrices Ql and the translations hl being known: the set Ω is made up of the points u of the volume D having position

QltR{right arrow over (a)}t,τl(tkl−t0l)t(π3s−u0l− Tkli)+Qlt(u0l−hl),

for each pixel s of each cross-section image Xkl of each sequence Sl, with τl, {right arrow over (a)}l, u0I and Tkl that are, respectively, the parameters of the movement for the sequence Sl, Ql=I3 the 3-by-3 identity matrix, h1=0, Tkl=R{right arrow over (a)}l,τl(tkl−t1l)T1l+R{right arrow over (a)}l,τl(tkl−t2l)T2l+ . . . +R{right arrow over (a)}l,τl(tkl−tk−1l)Tk−11+Tkl and T0l=0. We thus say that s, Xkl and Sl are the original pixel, image, and sequence of the point u, respectively.

Each of the points u of the set Ω is associated with the value Xkl (s) of the pixel s of the cross-section image Xkl of the sequence Sl: X(u)=Xkl(s) with s, Xkl and Sl respectively being the original pixel, image and sequence of the point u.
Taking the Acquisition Time into Account
During the capture of cross-section images, the object O moves during the exposure time to create an image on the focal plane. This creates non-negligible additional kinetic blurring (in addition to the blurring defined by the point spread function).
This kinetic blurring is, in one embodiment of the invention, taken into account in step 400 for the three-dimensional representation V of the object O, preferably when the exposure time is relatively long in relation to the time necessary to go from one cross-section image Xk to the next cross-section image Xk+1.
Thus, in reference to FIG. 15, in a third alternative, step 400 for determining the three-dimensional representation V of the object O comprises a step 450 for determining an acquisition interval Ik=[tk+δ0; tk+δ0+δ] for each cross-section image Xk, defined by an acquisition start moment tk+δ0 and an acquisition end moment tk+δ0+δ, the difference giving the acquisition time δ. One assumes, for 0≦k<m: tk≦tk+δ0<tk+δ0+δ≦tk+1.
In the described example, the acquisition start moment δ0 and the acquisition time δ are the same for all of the cross-section images. Nevertheless, the method easily extends to the case where the acquisition start moment δ0 and the acquisition time δ vary depending on the cross-section images.
During the exposure, the brightness level of each point of the focal plane evolves progressively by adding the luminosity level of the points on the path it crosses.
The position of the object O at each moment tk being known (or estimated), step 400 for the three-dimensional representation V of the object O also comprises a step 452 for determining the continuous position of the object O, as a function of time, between the successive moments tk and tk+1, and more particularly during the acquisition intervals Ik.
Thus, step 452 for determining the continuous position of the object O comprises a step 454 for determining a relation between the position ψk,t(u) at moment tk+t (between successive moments tk and tk+1) of a point of the object, which is situated at position u at moment tk, and the movement parameters of the object: angular speed τ, axis of rotation L and sequence of disrupting translations T1 . . . Tm.
Given that the time between moments tk and tk+1 is short, we consider that the movement is approximately made up of the stable rotation and a translation, linear fraction (proportional to the time) of the disrupting translation between moments tk and tk+1.
Thus, the position ψk,t(u) is given by:




  
    
      
        ψ
        
          k
          ,
          t
        
      
      ⁡
      
        (
        u
        )
      
    
    =
    
      
        
          R
          
            
              a
              →
            
            ,
            
              τ
              ⁢
              
                
              
              ⁢
              t
            
          
        
        ⁡
        
          (
          
            u
            -
            
              u
              0
            
          
          )
        
      
      +
      
        u
        0
      
      +
      
        
          t
          
            
              t
              
                k
                +
                1
              
            
            -
            
              t
              k
            
          
        
        ⁢
        
          T
          
            k
            +
            1
          
        
      
    
  
  ,
  
    
      for
      ⁢
      
        
      
      ⁢
      0
    
    ≤
    t
    ⁢
    
      ≤
      
        k
        +
        1
      
    
    ⁢
    
      -
      
        t
        k
      
    
  
  ,




with Tm+1=0.

One will note that ψk,0(u)=u and ψk,t+1−tk (u) is the position of the point of the object at moment tk+1 starting from position u at moment tk.
Step 452 for determining the continuous position of the object O also comprises a step 456 for determining a relation between the position Ck,t(s) at the initial moment t0 of the point o of the object, whereof the position at moment tk+t is on the pixel s of the cross-section plane P, using the preceding function ψk′. This position Ck,t(s) is given by:




  
    
      C
      
        k
        ,
        t
      
    
    ⁡
    
      (
      s
      )
    
  
  =
  
    
      
        R
        
          
            a
            →
          
          ,
          
            τ
            ⁡
            
              (
              
                
                  t
                  k
                
                -
                
                  t
                  0
                
                +
                t
              
              )
            
          
        
        t
      
      ⁡
      
        (
        
          
            
              π
              3
            
            ⁢
            s
          
          -
          
            u
            0
          
          -
          
            
              t
              
                
                  t
                  
                    k
                    +
                    1
                  
                
                -
                
                  t
                  k
                
              
            
            ⁢
            
              T
              
                k
                +
                1
              
            
          
          -
          
            
              R
              
                
                  a
                  _
                
                ,
                
                  τ
                  ⁢
                  
                    
                  
                  ⁢
                  t
                
              
            
            ⁢
            
              
                T
                _
              
              k
            
          
        
        )
      
    
    +
    
      
        u
        0
      
      .
    
  




Step 400 for determining the three-dimensional representation V of the object O also comprises a step 460 for choosing the operator Op to be integral on the acquisition interval Ik of the convolution of the three-dimensional representation function Vβ with the point spread function ƒR.
Step 400 for determining the three-dimensional representation V of the object O also comprises a step 462 for calculating, for each point u of the set Ω, an estimate {tilde over (X)}(u) of the value of that point u, using the operator Op. The estimate {tilde over (X)}(u) is given by:




  
    
      
        X
        ~
      
      ⁡
      
        (
        u
        )
      
    
    =
    
      
        ∫
        
          =
          
            δ
            0
          
        
        
          
            δ
            0
          
          +
          δ
        
      
      ⁢
      
        
          V
          β
        
        *
        
          
            f
            
              R
              
                
                  a
                  →
                
                ,
                
                  τ
                  ⁡
                  
                    (
                    
                      
                        t
                        k
                      
                      -
                      
                        t
                        0
                      
                      +
                      t
                    
                    )
                  
                
              
            
          
          ⁡
          
            (
            
              
                C
                
                  k
                  ,
                  t
                
              
              ⁡
              
                (
                s
                )
              
            
            )
          
        
        ⁢
        
          
        
        ⁢
        
          d
          t
        
      
    
  
  ,




where s and the cross-section image Xk (of image capture moment tk) are the original pixel and image of the point u, respectively.

As for the second alternative, step 400 for determining the three-dimensional representation V of the object O also comprises a step 464 for choosing the three-dimensional representation Vβ as a decomposition into radial basis functions φ:




  
    ∀
    
      u
      ∈
      D
    
  
  ,
  
    
      
        V
        β
      
      ⁡
      
        (
        u
        )
      
    
    =
    
      
        ∑
        
          w
          ∈
          W
        
      
      ⁢
      
        
          φ
          ⁡
          
            (
            
              u
              -
              w
            
            )
          
        
        ⁢
        
          β
          ⁡
          
            (
            w
            )
          
        
      
    
  
  ,




with W the set of nodes in the volume D.

Step 400 for determining the three-dimensional representation V of the object O also comprises a step 466 for determining the coefficients β(w).
More precisely, by replacing Vβ with its decomposition in radial basis functions, one obtains, for each point u of the set Ω:




  
    
      
        X
        ~
      
      ⁡
      
        (
        u
        )
      
    
    =
    
      
        ∑
        
          w
          ∈
          W
        
      
      ⁢
      
        
          ∫
          
            =
            
              δ
              0
            
          
          
            
              δ
              0
            
            +
            δ
          
        
        ⁢
        
          
            γ
            ⁡
            
              (
              
                
                  
                    π
                    3
                  
                  ⁢
                  s
                
                +
                
                  
                    R
                    
                      
                        a
                        →
                      
                      ,
                      
                        τ
                        ⁡
                        
                          (
                          
                            
                              t
                              k
                            
                            -
                            
                              t
                              0
                            
                            +
                            t
                          
                          )
                        
                      
                    
                  
                  ⁢
                  
                    (
                    
                      
                        u
                        0
                      
                      -
                      w
                    
                    )
                  
                
                -
                
                  
                    t
                    
                      
                        t
                        
                          k
                          +
                          1
                        
                      
                      -
                      
                        t
                        k
                      
                    
                  
                  ⁢
                  
                    T
                    
                      k
                      +
                      1
                    
                  
                
                -
                
                  
                    R
                    
                      
                        a
                        →
                      
                      ,
                      
                        τ
                        ⁢
                        
                          
                        
                        ⁢
                        t
                      
                    
                  
                  ⁢
                  
                    
                      T
                      _
                    
                    k
                  
                
                -
                
                  u
                  0
                
              
              )
            
          
          ⁢
          
            d
            t
          
          ⁢
          
            
          
          ⁢
          
            β
            ⁡
            
              (
              w
              )
            
          
        
      
    
  
  ,




where s and the cross-section image Xk (of image capture moment tk) are the original pixel and image, respectively, of the point u, and γ=φ*f.

Let us note




  
    
      γ
      
        k
        ,
        w
      
    
    ⁡
    
      (
      x
      )
    
  
  =
  
    
      ∫
      
        =
        
          δ
          0
        
      
      
        
          δ
          0
        
        +
        δ
      
    
    ⁢
    
      
        γ
        ⁡
        
          (
          
            x
            +
            
              
                R
                
                  
                    a
                    →
                  
                  ,
                  
                    τ
                    ⁡
                    
                      (
                      
                        
                          t
                          k
                        
                        -
                        
                          t
                          0
                        
                        +
                        t
                      
                      )
                    
                  
                
              
              ⁡
              
                (
                
                  
                    u
                    0
                  
                  -
                  w
                
                )
              
            
            -
            
              
                t
                
                  
                    t
                    
                      k
                      +
                      1
                    
                  
                  -
                  
                    t
                    k
                  
                
              
              ⁢
              
                T
                
                  k
                  +
                  1
                
              
            
            -
            
              
                R
                
                  
                    a
                    _
                  
                  ,
                  
                    τ
                    ⁢
                    
                      
                    
                    ⁢
                    t
                  
                
              
              ⁢
              
                
                  T
                  _
                
                k
              
            
            -
            
              u
              0
            
          
          )
        
      
      ⁢
      
        
          d
          t
        
        .
      
    
  




We therefore have, for each point u of the set Ω:




  
    
      
        X
        ~
      
      ⁡
      
        (
        u
        )
      
    
    =
    
      
        ∑
        
          w
          ∈
          W
        
      
      ⁢
      
        
          
            γ
            
              k
              ,
              w
            
          
          ⁡
          
            (
            
              
                π
                3
              
              ⁢
              s
            
            )
          
        
        ⁢
        
          β
          ⁡
          
            (
            w
            )
          
        
      
    
  
  ,




where s and the cross-section image Xk (of image capture moment tk) are the original pixel and image, respectively, of the point u.
As for the second alternative, γ=φ*ƒ is preferably calculated analytically or numerically by approximate calculations.
Thus, step 466 for determining the coefficients β(w) comprises a step for calculating γk,w.
If the analytical calculation of γk,w is heavy or impossible, γk,w is approximated by a discrete sum, for example the Riemann sum:




  
    
      
        γ
        
          k
          ,
          w
        
      
      ⁡
      
        (
        x
        )
      
    
    ≈
    
      
        δ
        J
      
      ⁢
      
        
          ∑
          
            j
            =
            0
          
          
            J
            -
            1
          
        
        ⁢
        
          γ
          (
          
            
          
          ⁢
          
            x
            +
            
              
                R
                
                  
                    a
                    →
                  
                  ,
                  
                    τ
                    ⁡
                    
                      (
                      
                        
                          t
                          k
                        
                        -
                        
                          t
                          0
                        
                        +
                        
                          δ
                          0
                        
                        +
                        
                          δ
                          ⁢
                          
                            j
                            J
                          
                        
                      
                      )
                    
                  
                
              
              ⁡
              
                (
                
                  
                    u
                    0
                  
                  -
                  w
                
                )
              
            
            -
            
              
                
                  
                    δ
                    0
                  
                  +
                  
                    δ
                    ⁢
                    
                      j
                      J
                    
                  
                
                
                  
                    t
                    
                      k
                      +
                      1
                    
                  
                  -
                  
                    t
                    k
                  
                
              
              ⁢
              
                T
                
                  k
                  +
                  1
                
              
            
            -
            
              
                R
                
                  
                    a
                    _
                  
                  ,
                  
                    τ
                    ⁡
                    
                      (
                      
                        
                          δ
                          0
                        
                        +
                        
                          δ
                          ⁢
                          
                            j
                            J
                          
                        
                      
                      )
                    
                  
                
              
              ⁢
              
                
                  T
                  _
                
                k
              
            
            -
            
              u
              0
            
          
          )
        
      
    
  
  ,




with J a fairly large integer, for example J≈20.

More generally, in the case where several sequences Sl, l=1, . . . , I are used, the rotations Ql and translations hl being known, step 466 for determining the coefficients β(w) comprises a step consisting of substantially placing all of the positions Ol of the representations of the object O in a same reference position. Let us recall that by choosing, for example, the first sequence as the reference sequence, Ql=I3 (with I3 the 3-by-3 identity matrix) and h1=0, one can decompose, in the reference point of the reference sequence, for each point u of the set Ω, the estimate {tilde over (X)}(u) of the value of that point u in a linear combination associated with the coefficients β(w):



  
    
      
    
    ⁢
    
      
        
          
            X
            ~
          
          ⁡
          
            (
            u
            )
          
        
        =
        
          
            ∑
            
              w
              ∈
              W
            
          
          ⁢
          
            
              
                γ
                
                  k
                  ,
                  w
                
                l
              
              ⁡
              
                (
                
                  
                    π
                    3
                  
                  ⁢
                  s
                
                )
              
            
            ⁢
            
              β
              ⁡
              
                (
                w
                )
              
            
          
        
      
      ,
      with
    
  




  
    
      
        γ
        
          k
          ,
          w
        
        l
      
      ⁡
      
        (
        x
        )
      
    
    =
    
       
      
        
          ∫
          
            t
            =
            
              δ
              0
            
          
          
            
              δ
              0
            
            +
            δ
          
        
        ⁢
        
          γ
          (
          
            
          
          ⁢
          
            
              R
              
                
                  
                    a
                    →
                  
                  l
                
                ,
                
                  
                    τ
                    l
                  
                  ⁡
                  
                    (
                    
                      
                        t
                        k
                        l
                      
                      -
                      
                        t
                        0
                        l
                      
                      +
                      t
                    
                    )
                  
                
              
            
            ⁢
            
              Q
              l
              t
            
            ⁢
            
              
                R
                
                  
                    
                      a
                      _
                    
                    l
                  
                  ,
                  
                    
                      τ
                      l
                    
                    ⁡
                    
                      (
                      
                        
                          t
                          k
                          l
                        
                        -
                        
                          t
                          0
                          l
                        
                        +
                        t
                      
                      )
                    
                  
                
                t
              
              ⁡
              
                (
                
                  x
                  -
                  
                    u
                    0
                    l
                  
                  -
                  
                    
                      t
                      
                        
                          t
                          
                            k
                            +
                            1
                          
                          l
                        
                        -
                        
                          t
                          k
                          l
                        
                      
                    
                    ⁢
                    
                      T
                      
                        k
                        +
                        1
                      
                      l
                    
                  
                  -
                  
                     
                    
                      
                        R
                        
                          
                            
                              a
                              →
                            
                            l
                          
                          ,
                          
                            
                              τ
                              l
                            
                            ⁢
                            t
                          
                        
                      
                      ⁢
                      
                        
                          T
                          _
                        
                        k
                        l
                      
                    
                    )
                  
                  +
                  
                    
                      R
                      
                        
                          
                            a
                            →
                          
                          l
                        
                        ,
                        
                          
                            τ
                            l
                          
                          ⁡
                          
                            (
                            
                              
                                t
                                k
                                t
                              
                              -
                              
                                t
                                0
                                t
                              
                              +
                              t
                            
                            )
                          
                        
                      
                    
                    ⁡
                    
                      (
                      
                        
                          
                            Q
                            l
                            t
                          
                          ⁡
                          
                            (
                            
                              
                                u
                                0
                                l
                              
                              -
                              
                                h
                                l
                              
                            
                            )
                          
                        
                        -
                        w
                      
                      )
                    
                  
                
                )
              
            
            ⁢
            
              d
              t
            
          
        
      
    
  



where s, the cross-section image Xkl (of image capture moment tkl) and Sl are the original pixel, image and sequence, respectively, of the point u, τl, {right arrow over (a)}l, u0l and Tkl are the parameters of the movement for the sequence Sl, Tml+1l=0 and γ=φ*ƒ.
Thus, in this case, step 466 comprises a step for calculating γk,wl.
Like γk,w(x), γk,w(x) can also be approximated by a discrete sum (e.g. Riemann sum).
Preferably, D is chosen to be large enough that D contains all of the points u of position

QltR{right arrow over (a)}t,τl(tkl−t0l)t(π3s−u0l− Tkl)+Qlt(u0l−hl),

for each pixel s of each cross-section image Xkl of each sequence Sl, with l=1 . . . I. The parameters β of the three-dimensional representation function Vβ are determined by minimizing the following quadratic cost function:





  
    
      E
      ⁡
      
        (
        β
        )
      
    
    =
    
      
        
          ∑
          
            l
            =
            1
          
          I
        
        ⁢
        
          
            ∑
            
              k
              =
              1
            
            
              m
              t
            
          
          ⁢
          
            
              ∑
              s
            
            ⁢
            
              
                (
                
                  
                    
                      X
                      k
                      l
                    
                    ⁡
                    
                      (
                      s
                      )
                    
                  
                  -
                  
                    
                      ∑
                      
                        w
                        ∈
                        W
                      
                    
                    ⁢
                    
                      
                        
                          γ
                          
                            k
                            ,
                            w
                          
                          l
                        
                        ⁡
                        
                          (
                          
                            
                              π
                              3
                            
                            ⁢
                            s
                          
                          )
                        
                      
                      ⁢
                      
                        β
                        ⁡
                        
                          (
                          w
                          )
                        
                      
                    
                  
                
                )
              
              2
            
          
        
      
      +
      
        
          λβ
          
            
              
            
            ⁢
            t
          
        
        ⁢
        A
        ⁢
        
          
        
        ⁢
        β
      
    
  
  ,




where A is a semi positive definite matrix, βmtAβ measures the irregularity of the coefficients β(w), for example A is chosen in the same way as before, λ>0 is the compromise parameter between the appropriateness of the three-dimensional representation function and the regularity of the coefficients. E(β) is written in the follow matrix form:

E(β)=∥y−Kβ∥2+λβtAβ, 

where β is the vector of the coefficients, the elements of the vector y are the Xkl(s), and the matrix K is made up of the elements yk,wl(π3s).

The result of minimizing E(β) is the solution of the following linear system:

(KtK+λA))β=Kty. 

To calculate the coefficients β(w), numerical optimization methods can be used, for example the conjugate gradient method, or the block optimization method previously presented.
Adjusting the Position of the Cross-Section Plane
It is desirable for the largest part of the object O to pass in the cross-section plane P. In one optimal configuration the axis of rotation L is substantially contained in the cross-section plane P. Thus, in one embodiment of the invention, the method comprises a step 53, inserted between steps 52 and 54, for adjusting the position of the cross-section plane P in relation to the axis of rotation L. Of course, the adjustment step can also be carried out independently of steps 50, 52, 54 and 56.
There are three ways to carry out adjustment step 53, in order to obtain substantially the optimal configuration, or, failing that, a suitable configuration.
First Adjustment Alternative
In a first alternative, adjustment step 53 comprises the movement and tilting of the optical microscope in order to place the axis of rotation L in the focal plane P. This alternative makes it possible to obtain the optical configuration.
This alternative comprises determining the axis of rotation L in the same way as in steps 54 and 56.
Second Adjustment Alternative
In a second alternative, adjustment step 53 comprises modifying the electric or electromagnetic field in order to place the axis of rotation L in the focal plane P. This alternative also makes it possible to obtain the optimal configuration.
This alternative comprises determining the axis of rotation L in the same way as in steps 54 and 56.
Third Adjustment Alternative
In reference to FIG. 16, in a third alternative, adjustment step 53 comprises the movement of the focal plane P along the direction Z (the direction perpendicular to the focal plane P) so that the focal plane P cuts the object O substantially in the middle thereof. The movement of the focal plane P is obtained by translating the optical microscope 12, and/or by acting on the lens 14 to change the focal distance.
In the described example, the middle of the object O is used as the barycenter of the object O.
Thus, in reference to FIG. 16, adjustment step 53 first comprises a step 53A for capturing a sequence of cross-section images X0 . . . Xm, then a step 53B for determining the angular speed τ, the axis of rotation L and the sequence of disrupting translations T1 . . . Tm. Step 53B is carried out in the same way as step 56.
Adjustment step 53 also comprises a step 53C for determining the set of points Ω in the same way as in step 404.
Adjustment step 53 also comprises a step 53D for determining a barycenter b of the light points of the set Ω. The barycenter b is preferably determined by the relation:




  
    b
    =
    
      
        
          ∑
          
            i
            =
            0
          
          n
        
        ⁢
        
          
            ∑
            s
          
          ⁢
          
            
              1
              
                
                  
                    X
                    i
                  
                  ⁡
                  
                    (
                    s
                    )
                  
                
                >
                α
              
            
            ⁢
            
              
                C
                k
              
              ⁡
              
                (
                s
                )
              
            
          
        
      
      
        
          ∑
          
            i
            =
            0
          
          n
        
        ⁢
        
          
            ∑
            s
          
          ⁢
          
            1
            
              
                
                  X
                  i
                
                ⁡
                
                  (
                  s
                  )
                
              
              >
              α
            
          
        
      
    
  
  ,




with X0 . . . Xn, n≦m, the part of the cross-section images X0 . . . Xm, that were captured during the time interval during which the object performs a maximum number of complete revolutions by rotating around the axis of rotation, 1B>A=1 when B is greater than A and 0 otherwise, and α is for example the gray level q-quantile of the cross-section images X0 . . . Xn (which means that the proportion of the pixels that have a gray level less than or equal to α is substantially equal to q). Generally, q is between 60% and 95%.

Adjustment step 53 also comprises a step 53E for calculating the projection b of the barycenter b on the axis of rotation L.
Adjustment step 53 also comprises a step 53F for adjusting the imaging system 10 to bring the cross-section plane P on the projection b of the barycenter b.
Annex
Let φ be a radial basis function. For any rotation matrix R, φ(Rx)=φ(x). We then have:




  
    
      
        
          φ
          *
          
            
              f
              R
            
            ⁡
            
              (
              x
              )
            
          
        
        =
        
          ∫
          
            
              φ
              ⁡
              
                (
                
                  x
                  -
                  u
                
                )
              
            
            ⁢
            
              f
              ⁡
              
                (
                
                  R
                  ⁢
                  
                    
                  
                  ⁢
                  u
                
                )
              
            
            ⁢
            
              d
              u
            
          
        
      
    
  
  
    
      
        =
        
          ∫
          
            
              φ
              ⁡
              
                (
                
                  x
                  -
                  
                    
                      R
                      t
                    
                    ⁢
                    y
                  
                
                )
              
            
            ⁢
            
              f
              ⁡
              
                (
                y
                )
              
            
            ⁢
            
              d
              y
            
          
        
      
    
  
  
    
      
        =
        
          ∫
          
            
              φ
              ⁡
              
                (
                
                  
                    R
                    ⁢
                    
                      
                    
                    ⁢
                    x
                  
                  -
                  y
                
                )
              
            
            ⁢
            
              f
              ⁡
              
                (
                y
                )
              
            
            ⁢
            
              d
              y
            
          
        
      
    
  
  
    
      
        =
        
          φ
          *
          
            
              f
              ⁡
              
                (
                
                  R
                  ⁢
                  
                    
                  
                  ⁢
                  x
                
                )
              
            
            .

## Claims

1. A method for determining a three-dimensional representation (V) of an object (O) using a sequence of cross-section images (X0 . . . Xm) of the object (O) in a cross-section plane(P), each cross-section image (X0 . . . Xm) having been taken at a respective moment (t0 . . . tm) of image capture, of the type comprising:
for each cross-section image (Xk), determining theposition of the object (O) in relation to the cross-section plane(P) at the moment (tk) the cross-section image (Xk) is captured,
determining a three-dimensional representation (V) of the object (O) using cross-section images (X0 . . . Xm) and positions of the object (O) in relation to the cross-section plane (P) at each moment of image capture (t0 . . . tm), the positions of the object(O) in relation to the cross-section plane (P) being determined using the cross-section images (X0 . . . Xm),
the method being characterized in that the object (O) has a movement in relation to the cross-section plane (P) made up, on one hand, of a regular movement defined by regular movement parameters (L, τ) and, on the other hand, a sequence of displacements (T1 . . . Tm) disrupting the object between the moment (tk) the cross-section image (Xk) is captured and the following moment(tk+1) of image (Xk+1) capture, the disrupting displacements (T1 . . . Tm) being disrupting translations, of variable direction and value, the determination of the positions of the object (O) in relation to the cross-section plane (P) comprises the determination of the sequence of disrupting displacements (T1 . . . Tm) using the cross-section images (X0 . . . Xm), the regular movement being a rotational movement defined by the following movement parameters: a fixed axis of rotation (L) and a constant angular speed (τ) determining the positions of the object (O) in relation to the cross-section plane (P) comprises determining at least part of the regular movement parameters (L, τ) using the cross-section images (X0 . . . Xm).

2. The method according to claim 1, also characterized in that determining the positions of the object (O) in relation to the cross-section plane (P) comprises determining the absolute value of the angular speed (|τ51 ) from the determination of a period of revolution (p) such that two cross-section images (Xk, Xk′) taken at respective moments (tk, tk′) separated from each other by a time substantially equal to a non-zero multiple of the period (p), are substantially similar.

3. The method according to claim 2, also characterized in that determining the period (p) comprises:
determining an initial group of candidate periods (p1 . . . pn),
selecting (240), among the candidate periods (p1 . . . pn), for each cross-section image (Xk), a first subset (pj(k,1), . . . , pj(k,e)) grouping together candidate periods p1 . . . pn having the highest levels of periodic similarity, the periodic similarity level (sim(Xk,pj)) characterizing, for each candidate period (pj), the similarity level of the cross-section image (Xk) with substantially periodic cross-section images (Xk′), for the candidate period (pj),
determining a law of probability that a candidate period will be selected, the law of probability having the period (p) as a parameter,
selecting the period (p) from the candidate periods (p1 . . . pn) of the initial group, as being that which maximizes the likelihood, according to the law of probability, knowing the selected candidate periods.

4. The method according to claim 1, also characterized in that determining positions of the object (O) in relation to the cross-section plane (P) comprises, the absolute value of the angular speed (|τ|) being known:
for each cross-section image (Xk), determining cross-section images (Xk′) spatially neighboring on the cross-section image (Xk),
for each cross-section image (Xk), calculating an orthogonal projection (Fk,k′(Xk)) of the cross-section image (Xk) on a support plane of each spatially neighboring cross-section image (Xk), the orthogonal projection (Fk,k′(Xk)) being expressed as a function of an axis of rotation value and a disrupting translation sequence value,
comparing the value of each of a set of pixels of the orthogonal projection (Fk,k′(Xk)) with the value of said pixel of the spatially neighboring cross-section image (Xk′) on which the projection was done,
determining the axis of rotation (L) and sequence of disrupting translations (T1 . . . Tm) whereof the values yield the values of close pixels during the comparison.

5. The method according to claim 4, also characterized in that the orthogonal projection (Fk,k′(Xk)) of the cross-section image (Xk) is done by calculating an affine transformation (Aff(Qk,k′;vk,k′,Xk)) of the cross-section image (Xk), the affine transformation (Aff(Qk,k′,vk,k′,Xk) having a linear transformation component (Qk,k′({right arrow over (a)},τ)) and a translation component (vk,k′(L,τ,T1 . . . Tm)) depending on, respectively, the direction ({right arrow over (a)}) of the axis of rotation value and the rotational speed (τ) on one hand, and, on the other hand, the axis of rotation value (L), the rotational speed (τ) and the disrupting translation sequence value (T1 . . . Tm).

6. The method according to claim 1, also characterized in that determining positions of the object (O) in relation to the cross-section plane (P) comprises, the absolute value of the angular speed (|τ|) being known:
for each cross-section image (Xk), determining cross-section images (Xk′) spatially neighboring on the cross-section image (Xk),
for each cross-section image (Xk), calculating an orthogonal projection (Fk,k′(Xk)) of the cross-section image (Xk) on a support plane (Xk′) of each spatially neighboring cross-section image (Xk′), the orthogonal projection (Fk,k′(Xk)) being expressed as a function of an axis of rotation value and a disrupting translation sequence value, the orthogonal projection (Fk,k′(Xk)) being done by calculating an affine transformation (Aff(Qk,k′,vk,k′,Xk)) of the cross-section image (Xk), the affine transformation (Aff(Qk,k′,vk,k′,Xk)) having a linear transformation component (Qk,k′({right arrow over (a)},τ)) and a translation component (vk,k′), the linear transformation component (Qk,k′({right arrow over (a)},τ)) depending on the direction of the axis ({right arrow over (a)}) of rotation value and the speed of rotation (τ),
comparing the value of each of a set of pixels of the orthogonal projection (Fk,k′(Xk)) with the value of that pixel of the spatially neighboring cross-section image (Xk′) on which the projection was done,
determining the direction of the axis of rotation (L) and the so-called reference translation components (vk,k′), whereof the values yield close pixel values during the comparison.

7. The method according to claim 6, also characterized in that determining an estimate of the position of the object (O) in relation to the cross-section plane (P) comprises, the direction ({right arrow over (a)}) of the axis of rotation (L) and the reference translation components (vk,k′) having been determined:
the expression of affine transformation translation components (vk,k′(u0,{right arrow over (a)},τ,T1 . . . Tm)) as a function of a disrupting translation sequence value (T1 . . . Tm) and a passage point value (u0) of the axis of rotation (L),
determining the sequence of disrupting translations (T1 . . . Tm) and a passage point (u0) of the axis of rotation (L) whereof the values yield translation components (vk,k′(u0,{right arrow over (a)},τ,T1 . . . Tm)) close to reference translation components (vk,k′).

8. The method according to claim 1, also characterized in that determining positions of the object (O) in relation to the cross-section plane (P) comprises determining the axis of rotation(L).

9. The method according to claim 8, also characterized in that determining the axis of rotation (L) comprises determining a projection (LXY) of the axis of rotation (L) on the cross-section plane (P), using cross-section images (X0 . . . Xm).

10. The method according to claim 9, also characterized in that each cross-section image (Xk) comprising a grid of pixels, identical for all of the cross-section images (X0 . . . Xm), the method comprises:
for each pixel (s) of the grid, determining the histogram of the values said pixel (s) assumes in at least part of the sequence of cross-section images (X0 . . . Xm) and the histogram of the values assumed by the point symmetrical to the pixel in relation to a projection value (LXY) of the axis of rotation (L),
determining (344) the projection (LXY) of the axis of rotation (L) whereof the value yields close histograms between each pixel and its point of symmetry.

11. The method according to claim 10, also characterized in that it comprises:
determining, for at least one pixel (s) of the grid, a symmetrization time (t(s)), necessary for a point of the object to move from said pixel (s) to the point of symmetry (LXY(s)) in relation to the projection (LXY) of the axis of rotation (L),
determining an angle between the axis of rotation (L) and its projection (LXY) on the cross-section plane (P), using the symmetrization time(s) (t(s)).

12. The method according to claim 8, also characterized in that determining positions of the object (O) in relation to the cross-section plane (P) comprises, the angular speed (τ) and axis of rotation (L) being known:
for each cross-section image (Xk), determining cross-section images (Xk′) spatially neighboring on the cross-section image (Xk),
for each cross-section image (Xk), calculating an orthogonal projection (Fk,k′(Xk) of the cross-section image (Xk) on a support plane (Xk′) of each spatially neighboring cross-section image (Xk′), the orthogonal projection (Fk,k′(Xk)) being expressed as a function of an axis of rotation value and a disrupting translation sequence value, the orthogonal projection (Fk,k′(Xk)) being done by calculating an affine transformation (Aff(Qk,k′,vk,k′,Xk)) of the cross-section image (Xk), the affine transformation (Aff(Qk,k′,vk,k′,Xk)) having a linear transformation component Qk,k′({right arrow over (a)},τ)) and a translation component (vk,k′), the linear transformation component (Qk,k′({right arrow over (a)},τ) depending on the direction ({right arrow over (a)}) of the axis of rotation value and the angular speed of rotation (τ),
comparing the value of each of a set of pixels of the orthogonal projection (Fk,k′(Xk)) with the value of said pixel of the spatially neighboring cross-section image (Xk′) on which the projection was done,
determining, the axis of rotation (L) being known, so-called reference translation components (vk,k′), the values of which yield close pixel values during the comparison,
expressing affine transformation translation components (vk,k′(u0,{right arrow over (a)},τ,T1 . . . Tm)) as a function of a disrupting translation sequence value (T1 . . . Tm) and a passage point value (u0) of the axis of rotation (L),
determining the sequence of disrupting translations (T1 . . . Tm) and a passage point (u0) of the axis of rotation (L) whereof the values yield translation components (vk,k′(u0,{right arrow over (a)},τ,T1 . . . Tm)) close to the reference translation components (vk,k′).

13. The computer program stored on a non-transitory readable medium, characterized in that it is designed, when implemented on a computer, to carry out the method according to claim 1.

14. The method for analyzing an object (O) to determine a three-dimensional representation (V) thereof, characterized in that it comprises:
introducing the object (O) into an imaging system defining a focal plane forming a cross-section plane (P),
configuring the imaging system to rotate the object (O) around a fixed axis of rotation (L) and with a fixed angular speed (τ),
capturing a sequence of cross-section images (X0 . . . Xm) of the object (O),
determining a three-dimensional representation (V) of the object (O) from the sequence of cross-section images (X0 . . . Xm), using a method according to claim 1.

15. The method according to claim 14, also characterized in that it comprises:
capturing a sequence of cross-section images of the object (O), called adjustment cross-section images,
determining the axis of rotation (L) of the object (O), using the adjustment cross-section images.

16. The method according to claim 15, also characterized in that it comprises adjusting the imaging system to move the cross-section plane (P) in relation to the object (O) such that the cross-section plane (P) passes through a middle ( b) of the object (O).

17. The method according to claim 16, also characterized in that the adjustment of the imaging system comprises determining a barycenter (b) of the object, and the middle of the object (O) is the projection of the barycenter (b) on the axis of rotation (L).

18. An imaging system characterized in that it comprises:
means making it possible to obtain images in a focal plane P,
a receptacle for receiving an object (O),
means for setting the object (O) in motion,
means for receiving cross-section images taken in the focal plane, adapted to implement a method according to claim 1.