Login| Sign Up| Help| Contact|

Patent Searching and Data


Title:
METHOD FOR OBTAINING A PHASE IMAGE FROM AN INTENSITY IMAGE
Document Type and Number:
WIPO Patent Application WO/2006/082327
Kind Code:
A1
Abstract:
The invention concerns the field of image processing. The invention concerns more particularly a method for obtaining a phase image from an intensity image. The method is characterized in that it consists in using filters which are measurements of a representation of a Green's function in a space of functions pertaining to Helmholtz equation, so that it is possible to apply integral transforms in a dimension of the space of the radiation wave field, thereby considerably improving the image processing time.

Inventors:
LYUBOSHENKO IGOR (FR)
Application Number:
PCT/FR2006/000254
Publication Date:
August 10, 2006
Filing Date:
February 03, 2006
Export Citation:
Click for automatic bibliography generation   Help
Assignee:
LYUBOSHENKO IGOR (FR)
International Classes:
G01J9/00
Domestic Patent References:
WO2000026622A12000-05-11
Other References:
LYUBOSHENKO I ET AL: "Phase unwrapping for interferometric synthetic aperture radar by use of Helmholtz equation eigenfunctions and the first Green's identity", JOURNAL OF THE OPTICAL SOCIETY OF AMERICA A (OPTICS, IMAGE SCIENCE AND VISION) OPT. SOC. AMERICA USA, vol. 16, no. 2, February 1999 (1999-02-01), pages 378 - 395, XP002347238, ISSN: 0740-3232
Attorney, Agent or Firm:
Breese, Derambure Majerowicz (PARIS, FR)
Download PDF:
Claims:
REVENDICATIONS
1. Procédé de récupération de phase de champ d'ondes de rayonnement comprenant les étapes suivantes : calculer une mesure représentative de la variation d' intensité dudit champ d'ondes de rayonnement sur une surface sélectionnée s ' étendant globalement d'un bout à l ' autre du champ d'ondes ; calculer une mesure représentative de l ' intensité dudit champ d' ondes de rayonnement sur ladite surface sélectionnée ; transformer ladite mesure de la variation d' intensité pour produire une première représentation de transformée intégrale et appliquer à ladite première représentation de transformée intégrale un premier filtre , pour produire une première représentation de transformée intégrale modifiée ; appliquer une première fonction de ladite première transformée intégrale à ladite première représentation de transformée intégrale modifiée pour produire une représentation non transformée ; appliquer une correction basée sur ladite mesure d' intensité sur ladite surface sélectionnée à ladite représentation non transformée ; transformer ladite représentation non transformée corrigée pour produire une seconde représentation de transformée intégrale et appliquer à ladite seconde représentation de transformée intégrale un second filtre ; appliquer une seconde fonction de ladite seconde transformée intégrale à ladite seconde représentation de transformée intégrale modifiée pour produire une mesure de phase dudit champ d' ondes de rayonnement d' un bout à l ' autre dudit plan sélectionné , caractérisé en ce que lesdits premier et second filtres sont des mesures d'une représentation d'une fonction de Green dans un espace de fonctions propres de l ' équation de Helmholtz à coordonnées séparables , de sorte de que lesdites première et seconde transformées intégrales sont des transformées intégrales en une dimension de l ' espace dudit champ d'ondes de rayonnement.
2. Procédé selon la revendication 1 , caractérisé en ce que lesdites première et seconde transformées intégrales sont produites en utilisant une transformée de Fourier.
3. Procédé selon la revendication 2 , caractérisé en ce que ladite transformée de Fourier est une transformée de Fourier rapide .
4. Procédé selon la revendication 1 dans lequel ladite première fonction utilise une fonction trigonométrique de type cosinus , et dans lequel ladite seconde fonction utilise une fonction trigonométrique de type sinus .
5. Procédé selon la revendication 1 dans lequel ladite première fonction utilise une fonction trigonométrique de type sinus , et dans lequel ladite seconde fonction utilise une fonction trigonométrique de type cosinus .
6. Procédé selon la revendication 1 , caractérisé en ce que lesdits premier et second filtres sont sensiblement les mêmes .
7. Procédé selon la revendication 1 , caractérisé en ce que lesdites surfaces sont sensiblement planes .
8. Programme informatique pour exécuter les étapes selon l 'une quelconques des revendications 1 à 7.
9. Programme informatique stocké sur au moins un support de stockage lisible par ordinateur, caractérisé en ce qu' il comprend un moyen pour exécuter les étapes selon l 'une des revendications 1 à 7.
10. Appareil pour la récupération de la phase de champ d' ondes de rayonnement comprenant: un moyen d' acquisition d'une mesure représentative de la variation d ' intensité dudit champ d'ondes de rayonnement sur une surface sélectionnée s ' étendant globalement d'un bout à l ' autre du champ d'ondes ; un moyen d' acquisition d'une mesure représentative de l ' intensité dudit champ d' ondes de rayonnement sur ladite surface sélectionnée ; un moyen de traitement pour séquentiellement réaliser les étapes suivantes : ( i ) transformer ladite mesure de la variation d' intensité pour produire une première représentation de transformée intégrale et appliquer à ladite première représentation de transformée intégrale un premier filtre, pour produire une première représentation de transformée intégrale modifiée ; ( ii ) appliquer une première fonction de ladite première transformée intégrale à ladite première représentation de transformée intégrale modifiée pour produire une représentation non transformée ; ( iii) appliquer une correction basée sur ladite mesure d' intensité sur ladite surface sélectionnée à ladite représentation non transformée ; ( iv) transformer ladite représentation non transformée corrigée pour produire une seconde représentation de transformée intégrale et appliquer à ladite seconde représentation de transformée intégrale un second filtre ; (v) appliquer une seconde fonction de ladite seconde transformée intégrale à ladite seconde représentation de transformée intégrale modifiée pour produire une mesure de phase dudit champ d' ondes de rayonnement d'un bout à l ' autre dudit plan sélectionné, caractérisé en ce que lesdits premier et second filtres sont des mesures d'une représentation d 'une fonction de Green dans un espace de fonctions propres de l ' équation de Helmholtz à coordonnées séparables , de sorte de que lesdites premières et seconde transformées intégrales sont des transformées intégrales en une dimension de l 'espace dudit champ d' ondes de rayonnement .
11. Appareil selon la revendication 8 , caractérisé en ce que ledit moyen d ' acquisition d 'une mesure représentative de la variation d ' intensité dudit champ d'ondes de rayonnement comprend au moins une plateforme mobile .
12. Appareil selon la revendication 8 , caractérisé en ce que ledit moyen d' acquisition d'une mesure représentative de la variation d' intensité dudit champ d 'ondes de rayonnement comprend au moins une lentille à focale variable .
Description:
OBTENTION D 'UNE IMAGE DE PHASE À PARTIR D 'UNE IMAGE D' INTENSITÉ

La présente invention se rapporte au domaine du traitement d' image .

La présente invention se rapporte plus particulièrement à un procédé d ' obtention d'une image de phase à partir d'une image d' intensité .

On s ' intéresse depuis une vingtaine d' années à l ' obtention d'une image de phase à partir d'une image d' intensité . Dans la publication « Relationship between two- dimensional intensity and phase in Fresnel Diffraction zone » , Abramoshkin et al . ( 1989 ) fournissent la valeur du gradient de la phase en fonction de l ' intensité, du gradient de l ' intensité et des fonctions de Green dans l ' espace infini .

Plus précisément, dans la publication de Fornaro et al . « Interferometric SAR phase unwrapping Using Green' s Formulation » , on divulgue l 'utilisation de transformées de Fourier rapides pour calculer la phase en fonction du gradient de phase par les fonctions de Green.

Cette formulation par transformée de Fourier est utilisée dans la demande PCT WO 00/26622 (Nugent et Al . ) qui décrit un procédé d'obtention d'une telle image de phase . Ce procédé concerne la récupération de phase de champ d ' ondes de rayonnement par résolution de l ' équation de transfert d ' intensité . On détermine d ' abord le taux de variation de l ' intensité de l ' image, orthogonalement par rapport à une surface qui s ' étend dans le champ d ' ondes (c ' est-à-dire en mesurant l ' intensité dans deux plans séparés ) . Ce taux est

ensuite soumis au processus de calcul suivant : établissement d ' une transformée intégrale , multiplication par un filtre correspondant à l ' inversion d ' un opérateur différentiel, et établissement d ' une transformée intégrale inverse . On multiplie le résultat par une fonction de l ' intensité par rapport à la surface . Les filtres ont une forme fondée sur les caractéristiques du système optique utilisé pour l ' acquisition des données d ' intensité, du type ouvertures numériques et fréquences spatiales .

Plus spécifiquement, la transformée de Fourier utilisée dans ce document de l ' art antérieur est une transformée de Fourier à deux dimensions .

Or, un des problèmes majeurs dans le domaine du traitement d' image , est le temps de calcul des applications représentant les résolutions des équations caractéristiques de l ' image .

Dans le présent cas , les temps de calculs proposés par le document de l ' art antérieur ne sont pas satisfaisants pour l ' homme du métier désirant effectuer ces calculs sur un grand nombre d' images .

La présente invention résout donc ce problème technique en fournissant une nouvelle méthode d'obtention d'une image de phase à partir d' une image d' intensité .

Par ailleurs , on notera que la demande PCT précitée ne mentionne pas l ' intervention de contraintes aux limites , en supposant seulement que la phase est nulle à l ' infini .

Un objet de la présente invention est donc de prendre en compte ces contraintes , tout en gardant un temps de calcul efficace .

Pour ce faire , la présente invention est du type décrit ci-dessus et elle est remarquable , dans son acception la plus large, en ce qu ' elle concerne un procédé de récupération de phase de champ d' ondes de rayonnement comprenant les étapes suivantes : - calculer une mesure représentative de la variation d' intensité dudit champ d' ondes de rayonnement sur une surface sélectionnée s ' étendant globalement d'un bout à l ' autre du champ d' ondes ;

- calculer une mesure représentative de l ' intensité dudit champ d'ondes de rayonnement sur ladite surface sélectionnée ;

- transformer ladite mesure de la variation d' intensité pour produire une première représentation de transformée intégrale et appliquer à ladite première représentation de transformée intégrale un premier filtre , pour produire une première représentation de transformée intégrale modifiée ;

- appliquer une première fonction de ladite première transformée intégrale à ladite première représentation de transformée intégrale modifiée pour produire une représentation non transformée ;

- appliquer une correction basée sur ladite mesure d' intensité sur ladite surface sélectionnée à ladite représentation non transformée ; - transformer ladite représentation non transformée corrigée pour produire une seconde représentation de transformée intégrale et appliquer à ladite seconde représentation de transformée intégrale un second filtre ;

- appliquer une seconde fonction de ladite seconde transformée intégrale à ladite seconde représentation de transformée intégrale modifiée pour produire une mesure de phase dudit champ d' ondes de rayonnement d'un bout à l ' autre dudit plan sélectionné ,

caractérisé en ce que lesdits premier et second filtres sont des mesures d'une représentation d'une fonction de Green dans un espace de fonctions propres de l ' équation de Helmholtz à coordonnées séparables , de sorte de que lesdites première et seconde transformées intégrales sont des transformées intégrales en une dimension de l ' espace dudit champ d'ondes de rayonnement.

De préférence , ladite première fonction utilise une fonction trigonométrique de type cosinus , et dans lequel ladite seconde fonction utilise une fonction trigonométrique de type sinus .

De la même façon, il est possible que ladite première fonction utilise une fonction trigonométrique de type sinus , et dans lequel ladite seconde fonction utilise une fonction trigonométrique de type cosinus .

Avantageusement, lesdites première et seconde transformées intégrales sont produites en utilisant une transformée de Fourier, par exemple une transformée rapide .

Selon une variante , lesdits premier et second filtre sont sensiblement les mêmes .

Selon un mode de réalisation, lesdites surfaces sont sensiblement planes .

L' invention concerne également un programme informatique, éventuellement stocké sur un support, pour exécuter les différentes étapes du procédé selon l ' invention.

L ' invention concerne également un appareil pour la récupération de la phase de champ d'ondes de rayonnement comprenant:

- un moyen d' acquisition d'une mesure représentative de la variation d' intensité dudit champ d'ondes de rayonnement sur une surface sélectionnée s ' étendant globalement d'un bout à l ' autre du champ d'ondes ;

- un moyen d' acquisition d'une mesure représentative de l ' intensité dudit champ d' ondes de rayonnement sur ladite surface sélectionnée ;

- un moyen de traitement pour séquentiellement réaliser les étapes suivantes :

( i ) transformer ladite mesure de la variation d' intensité pour produire une première représentation de transformée intégrale et appliquer à ladite première représentation de transformée intégrale un premier filtre, pour produire une première représentation de transformée intégrale modifiée ; ( ii ) appliquer une première fonction de ladite première transformée intégrale à ladite première représentation de transformée intégrale modifiée pour produire une représentation non transformée ; ( iii ) appliquer une correction basée sur ladite mesure d' intensité sur ladite surface sélectionnée à ladite représentation non transformée ;

( iv) transformer ladite représentation non transformée corrigée pour produire une seconde représentation de transformée intégrale et appliquer à

ladite seconde représentation de transformée intégrale un second filtre ;

(v) appliquer une seconde fonction de ladite seconde transformée intégrale à ladite seconde représentation de transformée intégrale modifiée pour produire une mesure de phase dudit champ d' ondes de rayonnement d'un bout à l ' autre dudit plan sélectionné , caractérisé en ce que lesdits premier et second filtres sont des mesures d'une représentation d'une fonction de Green dans un espace de fonctions propres de l ' équation de Helmholtz à coordonnées séparables , de sorte de que lesdites première et seconde transformées intégrales sont des transformées intégrales en une dimension de l ' espace dudit champ d' ondes de rayonnement .

Selon un premier mode de réalisation, ledit moyen d' acquisition d'une mesure représentative de la variation d' intensité dudit champ d 'ondes de rayonnement comprend au moins une plateforme mobile .

Selon un second mode de réalisation, ledit moyen d' acquisition d'une mesure représentative de la variation d' intensité dudit champ d' ondes de rayonnement comprend au moins une lentille à focale variable .

On comprendra mieux l ' invention à l ' aide de la description, faite ci-après à titre purement explicatif , d'un mode de réalisation de l ' invention, en référence aux figures annexées où : - la figure 1 représente les données d' image selon 1 ' invention ; la figure 2 illustre la comparaison des multiplications complexes entre la présente invention et les solutions de l ' art antérieur, en particulier celle du

document WO 00/26622 pour une image traitée allant de 64*64 pixels à 128*128 pixels ; la figure 3 illustre la comparaison des multiplications complexes entre la présente invention et les solutions de l ' art antérieur, en particulier celle du document WO 00/26622 pour une image traitée allant de

1024*1024 pixels à 2047*2047 pixels ; la figure 4 illustre les données d' images initiales à l ' origine de l ' algorithme selon la présente invention ; les figures 5A, 5B et 5C, qui doivent se lire de façon successive , représentent un mode d' implémentation exemplaire associé à la présente invention ; les figures 6 à 9 représentent des dispositifs pour la mise en œuvre de l ' invention.

La présente invention est décrite en détail ci-dessous dans un système de coordonnées cartésiennes et une image de type rectangulaire . Il est entendu que l ' homme du métier est apte à adapter l ' enseignement de la présente demande pour d' autres types d' images et/ou d' autres types de coordonnées , à la lecture de la présente description.

NOTATIONS ET PRINCIPE GENERAL DE L' INVENTION

La présente invention concerne de façon générale une méthode de reconstruction de la fonction de phase en un point du plan 2D à partir d'échantillons de la différence du premier ordre aux points (x,y)= (nΔx,mΔy). Nous calculons la phase seulement dans les points pour lesquels x =nΔx et y = mΔy où les lignes d'une grille bi-dimensionnelle ( 2D) se croisent, c ' est-à-dire dans les centres des points où se situent les éléments (pixels ) d'un capteur de

l ' énergie lumineuse , par exemple une caméra CCD. Les différences de phase du premier ordre aux points (x,y)= (nΔx,mAy) sont obtenues par le calcul du gradient de la phase à partir de la dérivée de l ' intensité δl/δz selon l ' axe optique z. Ce calcul sera détaillé ultérieurement . Afin de faciliter la transition de la notation discrète vers la notation continue , le descriptif des termes utilisés dans la description du procédé selon l' invention est présenté dans le Tableau I annexé à la présente demande .

Dans le document mentionné ci-dessus , Fornaro et al . ont utilisé la première identité de Green pour l 'estimation de la phase à partir des différences de premier ordre présentées sous forme continue comme le vecteur du gradient

. La phase est alors calculée comme la somme d 'une intégrale sur la région du support S de l ' image de phase et sur la frontière C de cette région :

Équation ( 1) :

où h s est un vecteur unitaire, perpendiculaire à C, orienté vers la surface extérieure de la région S comme illustré figure 1 , est la fonction de Green. Dans l ' équation ( 1 ) , la phase Φ entre tant dans la partie gauche que dans l ' intégrale sur le contour C. Cette double introduction de la même valeur dans deux parties de l 'équation ( 1 ) ne permet pas d'obtenir la phase si elle n 'est pas égale à zéro ou connue a priori sur le contour C.

Ainsi, l ' intégration de la phase sur le contour C signifie que la phase doit soit être connue , soit être égale

à zéro, soit être pré-calculée par un autre dispositif optique, par exemple des réseaux de lentilles de Shack- Hartmann avant la reconstruction de la phase à l ' intérieur de la surface S. L' intégration sur le contour dans l 'équation ( 1 ) peut également être supprimée si la fonction de Green choisie satisfait les conditions aux limites de Neumann sur la frontière C de la surface S :

Équation 2 :

Or, l 'équation ( 2 ) est valable pour la fonction de Green, donc la phase reconstruite selon l ' équation ( 1 ) satisfait les conditions aux limites de Neumann. L' intégration sur le contour dans l 'équation ( 1 ) peut aussi être supprimée si la fonction de Green choisie satisfait sur le contour C les conditions aux limites de Dirichlet et s ' il est connu que la phase φ(r) est nulle sur C, ce qui s 'écrit :

Équation 3 :

Or, comme précédemment, l ' équation ( 3 ) est valable pour la fonction de Green et pour la phase sur le contour C, donc la phase calculée selon ( 1 ) en utilisant la fonction de

Green de l 'équation ( 3 ) satisfait les conditions aux limites de Dirichlet.

Par ailleurs, il est connu de l 'homme du métier que la fonction de Green satisfait l ' équation :

équation 4 :

où est la fonction de Dirac .

Selon un aspect de l ' invention, pour satisfaire les équations précédentes , on choisit une fonction de Green exprimée sous forme d'une série de fonctions propres de l 'équation de Helmholtz comme suit. L 'homme du métier remarquera qu 'elle satisfera automatiquement les conditions aux limites ( 2 ) ou ( 3 ) selon la forme choisie ainsi que l ' identité ( 4 ) .

Pour une région de support rectangulaire S, ces fonctions spéciales sont :

équation 5 :

pour les conditions aux limites de Neumann, et

équation 6 :

pour les conditions aux limites de Dirichlet.

Pour un support S par exemple de dimension axb la constante A m n dans ( 5 ) et ( 6 ) est égale à si m ou n

sont nuls et dans les autres cas , où a est la longueur, h est la hauteur de la région rectangulaire de support S comme sur la figure 1.

On utilise aussi selon l ' invention, les valeurs propres associées aux fonctions ( 5 ) et ( 6 ) qui sont :

Équation 7 :

La fonction de Green pour un domaine S s ' écrit donc comme une série infinie aux fonctions propres de l ' équation de Helmholtz :

Équation 8 :

La dérivation de l ' expression générale ( 8 ) et la vérification des conditions aux limites ( 2 ) et ( 3 ) sont facilement réalisables par l ' homme du métier pour une telle fonction.

Ainsi, si la fonction de Green ( 8 ) est utilisée à la place de dans l 'équation ( 1 ) l ' intégration sur le contour dans l 'équation ( 1 ) peut être supprimée et la phase peut être exprimé sous forme compacte comme :

Équation 9 :

dans laquelle on notera simplement que le gradient est calculé par rapport à la variable r c 'est-à-dire que

V = Vr .

CALCULS NUMÉRIQUES PAR TRANSFORMÉE DE FOURIER SELON L' INVENTION

Un des objets de la présente invention est de calculer numériquement l ' équation ( 9 ) , à l ' aide des fonctions de

Green telles que définies pour obtenir la phase cherchée ,

La solution tout en obtenant un temps de calcul acceptable . est obtenue par la possibilité d'utilisation d ' une transformée de Fourier rapide uni-dimensionnelle ( IDFFT ou « One Dimensional Fast Fourier Transform » en anglais ) .

L 'utilisation de cette transformée en une dimension ne s ' impose pas du tout à l ' homme du métier au regard des équations générales ci-dessus , et nous fournissons ci- dessous les différentes étapes nécessaires à la compréhension du fait qu ' il est possible d'utiliser une telle transformée en une dimension.

Plusieurs transformations préliminaires sont nécessaires . Afin de faciliter les calculs de l ' homme du métier désirant réaliser l ' invention, nous mentionnons différentes étapes intermédiaires du calcul qui ne doivent pas être considérées comme des composantes essentielles de l ' invention, mais simplement des moyens d' obtenir une forme d' équation sur laquelle il est possible de mettre en œuvre les transformations en une dimension adéquate .

Pour plus de clarté , on considère d' abord les conditions aux limites de Neumann, avec γ = a/b et la

fonction de Green dans l 'expression de l 'équation ( 8 ) avec les fonctions spéciales ( 5 ) :

Équation 10 :

On sépare alors les cas n = 0, m ≠ O et n ≠ O, m = 0 pour obtenir :

équation 11 :

Or , on sait que :

Équation 12 :

Donc, Équation 13 :

avec Équation 14 :

Or,

Équation 15 :

Donc la fonction de Green est donnée par selon cette variante de l ' invention par :

Équation 16 :

Mais on a aussi

Équation 17 :

avec

équation 18

En dérivant ces expressions de la fonction de Green par rapport à x et y, il vient :

Équation 19

et

Équation 20

et en les substituant l 'équation ( 9 ) , il vient donc

Équation 21 :

où M' et N' sont les nombres entiers . Ils peuvent être choisis pour atteindre la précision numérique requise . En particulier, nous choisissons ces deux nombres comme les puissances de deux, plus grands que M +1 et N+ l , respectivement, M et N représentant les dimensions de l ' image à traiter.

Or, la présence des facteurs sin(η m y) et sin(ξ n x) dans les équations ( 19 ) et ( 20 ) , l ' intégration (par rapport à y dans la première somme de la partie droite de l 'équation ( 21 ) pour obtenir D x m ) , et par rapport à x dans la seconde somme de la partie droite pour obtenir D y n )) permet de calculer les séquences D x m ) et D y n ) comme les parties imaginaires (prises avec le signe négatif ) de la transformée de Fourier (FT) des dérivés (par rapport à y et par rapport à x, respectivement) de la phase .

Ces deux séquences peuvent être calculées numériquement en utilisant l ' algorithme de la transformée de Fourier Discrète, par l ' algorithme de la transformée de Fourier rapide ( IDFFT ou FFT ID indifféremment) appliqué aux différences de phase de premier ordre F y et F x , respectivement :

équations 22 et 23

On notera que «≈const signifie que le premier indice (i.e. n) est fixé dans F{n,m) tandis que la FFT ID de la

séquence résultante, dans laquelle le second indice ( i . e . m) varie, est calculé . Inversement, m = const signifie que le second indice ( i . e . m) est fixé dans F x (n,m) tandis que la FFT ID de la séquence résultante, dans laquelle le premier indice ( i . e . n) varie , est calculé . Les indices η m et ξ n indiquent les fréquences discrètes pour lesquelles les FFT ID doivent être évaluées .

L ' homme du métier comprendra par ailleurs facilement au regard des calculs ci-dessus que dans le cas des conditions aux limites de Dirichlet, la phase est donnée de la même façon par :

équation 24 :

avec les parties réelles des transformées de Fourier ID.

Equation 25 :

CALCUL A PARTIR DE L' INTENSITE

En référence encore une fois à la publication de Abramoshkin et Al. précédemment mentionnée, on connaît le gradient de la phase en fonction de la fonction de Green

Equation 26 :

où la fonction de Green utilisée est la fonction de Green pour l 'espace infini :

Equation 27 :

Donc , en réutilisant l 'équation ( 1 ) ,

Équation 28 :

En utilisant les fonctions de Green selon l ' invention (qui annulent la deuxième intégrale sur le contour) , il vient donc :

Équation 29 :

avec , pour les conditions de Neumann

Équation 30 :

avec

Équation 31

On obtient également les fonctions équivalentes pour les conditions aux limites de Dirichlet .

L ' homme du métier comprendra que l ' équation 28 est substantiellement celle utilisée par Νugent et al . dans la demande PCT WO 0026622. Dans cette demande, le deuxième terme est annulé en supposant que la phase est nulle sur le contour, alors que dans la présente demande, la phase peut être non nulle , le calcul est quand même réalisable par les conditions aux limites de Νeumann et de Dirichlet annulant ce terme . Par ailleurs , dans cette demande , l ' équation 28 est exprimée dans le domaine de Fourier et calculée par des

transformées de Fourier en deux dimensions . Au contraire , la présente invention entend effectuer ce calcul de phase en utilisant des transformées de Fourier à une dimension, grâce à l 'usage de fonctions de Green spécifiques , ce qui diminue fortement le nombre d' opérations complexes , et donc le temps de calcul des algorithmes associés au procédé selon l' invention.

Illustré figure 2 et 3 , on a calculé le nombre de multiplications complexes entre la solution de la présente invention par transformée de Fourier ID, et la solution de l ' art antérieur du document WO 00/26622 en fonction de la taille de l ' image traitée. Il est clair pour l ' homme du métier que la solution de la présente invention économise grandement le temps de traitement de l ' image pour la reconstruction de la phase .

ALGORITHME DE TRAITEMENT POUR LA RÉCUPÉRATION DE LA PHASE SELON L ' INVENTION.

Nous détaillons maintenant un exemple d' implémentation de la présente invention afin de compléter la description pour l 'homme du métier.

Pour cela, on utilise les notations suivantes :

Abréviations ;

- CL — Conditions aux limites ; - ID FFT — Transformée de Fourier Directe uni- dimensionnelle ;

ID IFFT — Transformée de Fourier Inverse uni- dimensionnelle ;

- — Partie réelle d'un nombre complexe entre parenthèses ;

- — Partie imaginaire d'un nombre complexe entre parenthèses .

Les paramètres d ' images en entrée suivants ;

- M x N, où M est le nombre de lignes de l ' image, indexées par i = 1 ,..., M ; et N est le nombre de colonnes , indexées par j = 1 ,..., N.

Constantes :

La figure 4 illustre une première étape de l ' algorithme selon l ' invention dans lequel, à partir d 'une image, on calcule dl/dz par différence entre les images de l ' intensité dans les plans z 1 et z 2 . Cette différence est bien sûr calculée pixel par pixel sur les images I 2 (m, n) et I 1 (m,n) . L' homme du métier comprendra que, dans le cas de coordonnées quelconques , cette mesure de dl/dz correspond en fait à toute mesure de la variation d' intensité sur une surface .

De façon générale , l ' homme du métier remarquera, en lien avec la description de l ' invention ci-dessus , que l ' algorithme donne la possibilité à chaque fois que cela est

nécessaire, de choisir entre les conditions aux limites de Neumann et de Dirichlet.

Illustré figure 5A, 5B et 5C, l ' homme du métier reconnaîtra une implémentation des équations données plus haut selon la présente invention. On notera aussi que l ' aspect discret des données réelles sur les images ( les valeurs en x et en y étant obtenues par des données de pixels ) , permettent un calcul discret des intégrales données dans les équations et une interversion des calculs dans les sommes discrètes finies . Les modes de calculs de ces sommes discrètes sont données dans le présent algorithme à titre d' exemple et peuvent bien sûr varier selon les choix de l ' homme du métier sans limiter la portée de la présente invention.

Illustré figure 5A, les hachures obliques de l ' image initiale représente la fenêtre de FFT utilisée pour les calculs . Cette fenêtre dépasse l ' image utile et peut être remplie par des 0 , selon la méthode connue du « zéro padding » , ou « remplissage par des zéros » , ou par tout autre moyen bien connu dans le domaine du traitement des images .

Les lignes verticales de la figure 1 indiquent alors le sens de la FFT ID utilisée selon l ' invention dans cette partie illustrée de l ' algorithme . Dans ce mode de réalisation de l ' invention, on calcule la FFT en fixant une colonne et en faisant varier l ' indice de ligne .

On notera par ailleurs figure 5A l ' introduction du paramètre M n ^ x utilisé de préférence selon ce mode d' implémentation pour éviter la divergence des termes en sinus hyperbolique et un calcul spécifique du coefficient

c y (m, n) au cas où l ' indice représentant les coordonnées dépasse cette valeur seuil .

On reconnaîtra à l ' opération A l ' obtention de la dérivée de la phase selon y tel qu' il est décrit dans la première partie de l ' équation 30. On notera que l ' opération

A nécessite la division par une mesure de l ' intensité qui peut faire diverger la solution si cette intensité est nulle. Dans ce cas , il est possible d' effectuer préalablement une détermination des zéros de l ' intensité pour corriger cette divergence, et/ou d' affecter une valeur minimale à l ' intensité en cas de nullité en un point de l ' espace considéré .

Selon l 'équation 22 , on applique ensuite le procédé de calcul selon l ' invention utilisant une transformée de Fourier ID au gradient de la phase selon y pour obtenir le premier membre de la phase comme dans l ' équation 21. Ce premier membre φ 2 (m, n) est obtenu à l ' opération B de la figure 5C .

De la même façon ( non illustré ) , on obtient la dérivée de la phase selon x tel qu' elle est décrite dans la première partie de l 'équation 30 , et selon l ' équation 23 , on applique ensuite le procédé de calcul selon l ' invention utilisant une transformée de Fourier ID au gradient de la phase selon y pour obtenir le second membre de la phase comme dans l ' équation 21 , φ 1 (m, n) .

La phase reconstruite selon l ' invention est alors calculée par φ(m, n) = φ 1 (m, n)+ φ 2 (m, n) .

De façon générale, illustré figure 5A, 5B et 5C, l ' algorithme selon l ' invention comprend donc , suite au

calcul d'une mesure de la variation d' intensité sur une surface (δl/δz dans le mode de réalisation illustré ) , l ' application d'une première transformation intégrale , de préférence sous la forme d'une transformée de Fourier rapide , pour obtenir la variable nommée (r x (m r n) dépendant des conditions aux limites utilisées . Cette variable est ensuite multipliée ( dans l ' espace de Fourier ou plus généralement dans l 'espace de projection) par un filtre correspondant aux fonctions et modifiées éventuellement par les contraintes de calcul de type M 103x .

Conformément aux équations 19 et 20 , ces filtres et sont des représentations d'une fonction de Green dans un espace de fonctions propres de l ' équation de Helmholtz , par exemple les fonctions associées aux équations 5 et 6.

On applique ensuite une transformée intégrale inverse, par exemple de type transformée de Fourier inverse , en prenant les parties réelles ou imaginaires selon les types de conditions aux limites .

Comme l ' indique par exemple l ' équation 26 , la division par une mesure de l ' intensité, par exemple celle prise sur un plan de calcul de la variation de l ' intensité, permet d'obtenir le gradient de la phase par rapport à une des coordonnées .

En conformité avec les équations 21 , 22 , et 23 , on réapplique ensuite une transformée intégrale à une dimension, par exemple de type Fourier, on réapplique un second filtre à la variable obtenue ( D x (m,n) dans le mode de réalisation proposé,) , et on reprend une transformée

intégrale inverse pour obtenir sensiblement un des termes de l 'expression définissant la phase recherchée .

Il est entendu qu' au regard de l ' exemple de mise en œuvre de l ' algorithme donné dans les figures et décrit ci- dessus , l ' homme du métier est apte à apporter des modifications d' implémentation, par exemple selon les capacités de calcul des machines utilisées , ou selon une caractéristique spécifique de l ' image , sans toutefois sortir du cadre de l ' invention.

AUTRES MODES DE RÉALISATION ASSOCIÉS À L' INVENTION

Au regard de la description faite ci-dessus , d'une part des équations nécessaires à la mise en œuvre de l ' algorithme de l ' invention, et d' autre part de l ' algorithme de traitement lui-même, l ' homme du métier comprendra que le mode de réalisation décrit n'est qu'un exemple non limitatif de la mise en œuvre de l ' invention.

En particulier, on comprendra que dans le cas d'un support non rectangulaire, l ' homme du métier est apte à déterminer les fonctions propres de l 'équation de Helmholtz pour une forme de support spéciale . À titre indicatif seulement, on donne la forme de ces fonctions dans le cas d'un support circulaire

Équation 5 ' :

où J m est la fonction de Bessel d' ordre m avec J m ' (πα mn )= 0, m,n = 0,K ,∞ , i . e . πα nm est le n-ième zéro de la

dérivée de la fonction de Bessel d'ordre m. L'équation 8 permet alors d'obtenir la fonction de Green associée .

L'utilisation de la fonction de Green associée à cette fonction propre est particulièrement utile dans le cas où la reconstruction de la phase dans une région circulaire donne les renseignements suffisants sur la topographie de l ' objet de géométrie circulaire, comme par exemple dans le cas de la visualisation et de la mesure de la face avant d'une fibre optique .

L' homme du métier est enfin apte à déterminer ces fonctions propres dans le cas de supports triangulaires ou semi-infinis .

Par ailleurs , il est entendu que nous avons décrit un mode de réalisation faisant intervenir des transformées de Fourier discrètes en une dimension, mais que l ' homme du métier est apte à mettre en œuvre des transformées de Bessel discrètes en une dimension, en tant qu ' analogue curviligne des transformées de Fourier. L ' homme du métier peut également utiliser tout type de transformation intégrale, en tant que projection d'une fonction sur un espace de fonction, par exemple déterminé par le type de fonction utilisé dans les équations 5 ou 5 ' sans pour autant sortir du cadre de la présente invention.

Dans tous les cas , l 'homme du métier aura compris que si les fonctions propres de l ' équation de Helmholtz s ' expriment en des variables séparées , en particulier si le domaine de l ' image n ' est pas trop distordu, on pourra utiliser la présente invention et n'effectuer des transformations intégrales qu 'en une seule dimension.

Pour des domaines très distordus , on pourra utiliser des transformations conformes pour ramener ces domaines à des domaines plus simples , ou bien définir ce domaine comme une composition de domaine plus simples .

De plus , il est entendu que le mode de réalisation décrit pour une image est également applicable pour tout type d'onde de rayonnement selon l 'ensemble du spectre lumineux. Afin d' appliquer l ' enseignement de la présente invention aux différents spectres , l ' homme du métier est apte à modifier les accessoires de captation tels que décrits dans la partie suivante consacrée aux accessoires associées à l ' invention.

DISPOSITIFS DE MESURE POUR LA MISE EN ŒUVRE DE L' INVENTION

Afin de mettre en œuvre l ' invention, et conformément par exemple à la figure 4 , il convient d'utiliser des dispositifs de mesure de l ' intensité dans des plans distincts .

Selon un mode de réalisation, ces mesures sont réalisées dans des plans perpendiculaires à l ' axe optique .

Pour ce faire, illustré figure 6 , on réalise une prise d' image sur un échantillon transparent. L' homme du métier comprendra que pour effectuer des mesures sur deux plans perpendiculaires , on peut soit modifier la position de l 'échantillon par rapport au dispositif optique fixe , soit déplacer le dispositif optique .

Sur la figure 6 , pour un échantillon transparent, le dispositif selon un mode de réalisation de l ' invention

comprend donc un capteur d' images 61 , par exemple sous la forme d'un capteur CCD, un système optique de projection 62 pour projeter la lumière issue de l ' échantillon 63 vers le capteur. La lumière est émise par un système optique d'éclairage 64. L' échantillon est mobile selon l ' axe des z grâce à un dispositif de déplacement (non représenté ) .

Illustré figure 7 , selon un autre mode de réalisation de l ' invention, l ' échantillon est fixé à une plateforme immobile 75. Le dispositif optique comprend toujours un capteur 71 , un système de projection 72 et un système d' éclairage 73 , et il comprend en outre une plateforme 76 déplaçable à laquelle sont fixés les éléments 71 , 72 et 73. Cette plateforme est actionnée par un système de déplacement 77 sous la forme par exemple d'un moteur ou d'un élément piézo-électrique .

Illustré figure 8 , l ' invention comprend également un dispositif de mesure utilisant une ou plusieurs lentilles à focale variable . Dans ce cas , on dispose par exemple une lentille à focale variable 88 au niveau du système de projection et selon une variante une lentille à focale variable 89 au niveau du système d'éclairage pour des mesures sur un échantillon transparent fixe . L' homme du métier comprendra que dans ce cas le déplacement des mesures d' intensité correspond à la différence de focale, c ' est-à- dire Δz= F 2 -F 1 . Dans le cas de deux lentilles , les modifications de la focale se feront simultanément pour les deux lentilles . On comprendra également que ce système est également réalisable avec une seule lentille à focale variable .

Illustré figure 9 , pour un échantillon diffusant et non plus transparent, on peut, selon un mode de réalisation,

utiliser une séparatrice du faisceau 90 destiné à envoyer un rayon lumineux issu du système d' éclairage 94. On peut alors utiliser l 'ensemble des éléments précédemment décrits tels que le capteur, le système de projection, l ' ensemble étant soit fixe soit mobile , ainsi qu'une ou plusieurs lentilles à focale variable tel que précédemment décrit.

Il est entendu que dans tous ces modes de réalisation, le capteur d' image est relié à un dispositif de traitement de l ' image comprenant un logiciel pour la mise en œuvre de l ' algorithme selon l ' invention.

Les images de phase obtenues par ce logiciel peuvent alors être visualisées sur un écran d' ordinateur.

L ' invention est décrite dans ce qui précède à titre d' exemple . Il est entendu que l ' homme du métier est à même de réaliser différentes variantes de l ' invention sans pour autant sortir du cadre du brevet.

U)

O

Tableau I