Brevet FR2984503A1 — publié le 20130621

Procédé de reconstruction de propriétés optiques d'un milieu diffusant à l'aide d'une combinaison d'une pluralité de transformées de Mellin-Laplace d'une grandeur comportant une distribution temporelle d'un signal reçu, et système de reconstruction associé La présente invention concerne un procédé de reconstruction de propriétés optiques d'un milieu diffusant à l'aide d'un système de reconstruction comportant au moins une source de rayonnement propre à éclairer un milieu et au moins un détecteur propre à recevoir un signal émis par le milieu, le procédé comprenant: - l'éclairement du milieu par la ou chaque source de rayonnement, - la réception, par le ou chaque détecteur, d'un signal émis par le milieu, et - l'élaboration, pour au moins un couple source-détecteur, d'une première distribution du signal reçu par le détecteur correspondant, - le calcul de la transformée de Mellin-Laplace, pour un ordre et une variable donnés, d'une grandeur comportant la première distribution, l'ordre étant un nombre entier, la variable étant un nombre réel, et - la reconstruction de propriétés optiques du milieu à l'aide de la transformée de Mellin-Laplace de ladite grandeur. L'invention concerne également un tel système de reconstruction. L'invention se situe notamment dans le domaine de l'instrumentation médicale, et en particulier dans le domaine de l'imagerie optique de diffusion appliquée à des tissus biologiques. L'invention s'applique à la détection de différentes pathologies, telles que des cancers du sein, de la prostate ou encore du tube digestif, sein ou cerveau. L'invention s'applique également à la dermatologie, notamment à la caractérisation d'une réaction d'inflammation cutanée. On connaît l'article « Optical tomography in medical imaging » de S.R. Arridge publiée dans la revue « Inverse Problems » en 1999, citant l'utilisation possible d'une transformée de Mellin-Laplace à l'ordre n pour des applications en tomographie optique. Cependant, cet article ne détaille pas précisément comment cette transformée peut être mise en oeuvre dans le cadre de la reconstruction des propriétés optiques d'un milieu diffusant, avec les avantages décrits ci-après. Le but de l'invention est donc de proposer un procédé de reconstruction des propriétés optiques d'un milieu diffusant ainsi qu'un système de reconstruction utilisant avantageusement la transformée de Mellin Laplace, permettant notamment l'accès à des propriétés optiques de diffusion usuellement difficiles à établir, et par exemple la distribution spatiale des coefficients d'absorption ou de diffusion d'un milieu. A cet effet, l'invention a pour objet un procédé de reconstruction du type précité, caractérisé en ce que l'étape de calcul comporte le calcul d'une pluralité de transformées de Mellin-Laplace de ladite grandeur pour des valeurs distinctes de l'ordre, et en ce que l'étape de reconstruction est effectuée à partir d'une combinaison de la pluralité de transformées de Mellin-Laplace. Par reconstruction des propriétés optiques, on entend par exemple : - la reconstruction des propriétés d'absorption, ces dernières étant notamment caractérisées par le coefficient d'absorption, noté pta, - la reconstruction des propriétés de diffusion, ces dernières étant notamment caractérisées par le coefficient de diffusion réduit il', ou par le coefficient de diffusion D, et - la reconstruction des propriétés de fluorescence, ces dernières étant notamment caractérisées par une fonction de réponse (F) d'un fluorophore, ou par une concentration c d'un fluorophore, ou encore par toute autre grandeur exprimant une quantité q d'un fluorophore, ce dernier étant par exemple endogène ou exogène. Suivant d'autres aspects avantageux de l'invention, le procédé de reconstruction comporte une ou plusieurs des caractéristiques suivantes, prises isolément ou suivant toutes les combinaisons techniquement possibles : - l'étape de calcul comporte le calcul d'au moins une transformée de MellinLaplace de ladite grandeur pour une valeur de l'ordre supérieure ou égale à 5, de préférence supérieure ou égale à 8, - l'étape de calcul comporte le calcul de nmax+l transformées de Mellin-Laplace, l'ordre prenant successivement toutes les valeurs comprises entre 0 et nmax, avec nmax supérieur ou égal à 5, de préférence supérieur ou égal à 8, - la valeur de la variable est comprise entre 1 ns-1 et 20 ns-1, - le procédé comprend en outre : + la détermination, pour au moins un couple source-détecteur, d'une première fonction de modélisation d'un signal de diffusion de la lumière entre la source et le détecteur dans le milieu, + l'élaboration, pour ledit au moins un couple source-détecteur, d'une deuxième distribution d'un signal reçu par le détecteur pour un milieu de référence, le signal reçu étant émis par le milieu de référence suite à l'éclairement dudit milieu par la source, et + la détermination, pour ledit au moins un couple source-détecteur, d'une deuxième fonction de modélisation d'un signal de diffusion de la lumière entre la source et le détecteur dans le milieu de référence, et dans lequel ladite grandeur dépend de la première distribution, de la première fonction de modélisation, de la deuxième distribution et de la deuxième fonction de modélisation, - ladite grandeur est obtenue par soustraction du produit de la deuxième distribution et de la première fonction de modélisation du produit de la première distribution et de la deuxième fonction de modélisation, - ladite grandeur est le ratio entre le produit de la première distribution et de la deuxième fonction de modélisation et le produit de la deuxième distribution et de la première fonction de modélisation, - la ou chaque source de rayonnement comporte une source de lumière pulsée, - le ou chaque détecteur est un détecteur résolu en temps, et - les propriétés optiques comprennent au moins un élément parmi le groupe consistant en : + les propriétés d'absorption de la lumière, caractérisées notamment par le coefficient d'absorption, + les propriétés de diffusion, caractérisées notamment par le coefficient de diffusion réduit ou par le coefficient de diffusion, et + les propriétés de fluorescence, caractérisées notamment par une fonction de réponse d'un fluorophore, ou par une concentration du fluorophore, ou encore par une grandeur fonction d'une quantité du fluorophore. L'invention a également pour objet un système de reconstruction comprenant - au moins une source de rayonnement propre à éclairer le milieu, - au moins un détecteur propre à recevoir un signal émis par le milieu, - des moyens d'élaboration, pour au moins un couple source-détecteur, d'une distribution du signal reçu par le détecteur correspondant, - des moyens de calcul de la transformée de Mellin-Laplace, pour un ordre et une variable donnés, d'une grandeur comportant la distribution, l'ordre étant un nombre entier, la variable étant un nombre réel, et - des moyens de reconstruction de propriétés optiques du milieu à l'aide de la transformée de Mellin-Laplace de ladite grandeur, caractérisé en ce que les moyens de calcul sont propres à calculer une pluralité de transformées de Mellin-Laplace de ladite grandeur pour des valeurs distinctes de l'ordre, et en ce que les moyens de reconstruction sont propres à reconstruire les propriétés optiques du milieu à partir d'une combinaison de la pluralité de transformées de Mellin- Laplace. Suivant un autre aspect avantageux de l'invention, le système de reconstruction comprend en outre des moyens de calcul propres à calculer au moins une transformée de Mellin-Laplace de ladite grandeur pour une valeur de l'ordre supérieure ou égale à 5, de préférence supérieure ou égale à 8. Ces caractéristiques et avantages de l'invention apparaîtront à la lecture de la description qui va suivre, donnée uniquement à titre d'exemple, et faite en référence aux dessins annexés sur lesquels : - la figure 1 est une représentation schématique d'un système de reconstruction selon l'invention, propre à reconstruire une image d'un milieu, - la figure 2 est un ensemble de courbes représentant une distribution temporelle d'un signal reçu par le système de la figure 1 et les transformées de Mellin-Laplace de cette distribution temporelle calculées pour des valeurs d'ordre comprises entre 0 et 9, - la figure 3 est un organigramme d'un procédé de reconstruction selon l'invention, - la figure 4 est un ensemble d'images représentant des cartes de sensibilité pour les transformées de Mellin-Laplace calculées, chaque image correspondant à une valeur d'ordre respective et comprise entre 0 et 8, - la figure 5 est une vue d'une image du milieu reconstruite avec un système de reconstruction et un procédé de reconstruction de l'état de la technique, et - la figure 6 est une vue analogue à celle de la figure 5 avec un système de reconstruction et un procédé de reconstruction selon l'invention. Sur la figure 1, un système de reconstruction 10 est destiné à réaliser un examen d'un milieu 12 via l'acquisition d'une image du milieu 12, puis le traitement de cette image. Le système de reconstruction 10 comprend une pluralité de sources de rayonnement 14 aptes à éclairer le milieu 12 et une pluralité de détecteurs 16 propres à recevoir un signal émis par le milieu 12. Les sources de rayonnement 14 et les détecteurs 16, sont, par exemple, agencés dans une sonde, non représentée. Le système 10 comprend une unité de traitement d'informations 18 et un écran 20 d'affichage d'une image du milieu. Dans le mode de réalisation décrit, le système 10 est un système temporel, et chaque source de rayonnement 14 est une source de rayonnement pulsée. Chaque détecteur 16 est un détecteur résolu en temps, c'est-à-dire un détecteur permettant de mesurer la distribution temporelle du temps d'arrivée des photons, également appelés TCSPC (de l'anglais Time-Correlated Single Photon Counting). Selon l'exemple représenté, le milieu 12 présente une surface d'observation 22, en forme d'un plan parallèle à un axe longitudinal X et un axe transversal Y, et au contact de laquelle les sources 14 et les détecteurs 16 sont propres à être appliqués. Le milieu 12 présente une direction d'observation 24 s'étendant selon un axe vertical Z et sensiblement perpendiculairement à la surface d'observation 22. Dans l'exemple de réalisation de la figure 1, le milieu 12 est un milieu à caractériser, noté MC et visible sur la figure 1. Le milieu à caractériser MC comporte des tissus biologiques d'un animal ou de l'homme. Le milieu à caractériser MC est, par exemple, une zone d'un organe tel que le cerveau, un sein, la prostate, le tube digestif ou encore d'autres organes dans lesquels des fluorophores sont susceptibles d'être injectés. Plus généralement, le milieu à caractériser MC est un milieu biologique, par exemple un organe, dont on souhaite déterminer des propriétés optiques de diffusion et notamment la distribution spatiale de coefficients tels que le coefficient d'absorption ou le coefficient de diffusion. Le milieu à caractériser MC est un milieu diffusant contenant des inclusions dont les propriétés optiques d'absorption ou de diffusion sont différentes de celles du milieu. 26, Une seule inclusion 26 est représentée sur la figure 1 pour la clarté du dessin. Chaque source de rayonnement pulsée 14 comporte une source lumineuse pulsée 28 et une fibre optique d'excitation 30 reliée à la source pulsée 28 pour la transmission de l'impulsion lumineuse jusqu'au milieu 12. Lorsqu'une telle source fibrée est utilisée, on assimile l'extrémité libre de la fibre d'excitation 30 à la source s. En variante non représentée, chaque source de rayonnement pulsée 14 comporte une fibre optique d'excitation 30 reliée à une unique source lumineuse pulsée commune à la pluralité de sources de rayonnement 14. Selon cette variante, le système 10 comprend en outre un commutateur optique ou un multiplexeur pour sélectionner la fibre d'excitation 30 dans laquelle le faisceau lumineux est envoyé. En variante encore, l'ensemble des sources de rayonnement pulsée 14 est formé d'une unique source lumineuse pulsée et d'un dispositif à miroir de type MEMS (de l'anglais MicroElectroMechanical Systems), non représenté, pour balayer le milieu 12 avec la lumière issue de la source lumineuse pulsée. En variante encore, la source 14 est une source multiplexée, délivrant des impulsions lumineuses à différentes longueur d'onde. Dans le cas d'un examen de la surface 22 du milieu ou d'un examen à une faible profondeur du milieu 12, c'est-à-dire une profondeur de quelques millimètres à quelques centimètres, la longueur d'ondes de chaque source de rayonnement pulsée 14 est, de préférence, dans le domaine du visible ou du proche infrarouge, c'est-à-dire comprise entre 500 nm et 1300 nm. Le taux de répétition est d'environ 50 MHz, et en général compris entre 10 MHz et 100 MHz. Les impulsions émises par chaque source de rayonnement pulsée 14 présentent une durée comprise entre 10 picosecondes et 1 nanoseconde, de préférence comprise 10 picosecondes et 100 picosecondes, chaque source lumineuse pulsée 28 étant propre à générer une impulsion de largeur temporelle de durée comprise entre 10 picosecondes et 1 nanoseconde, de préférence comprise 10 picosecondes et 100 picosecondes. Chaque détecteur résolu en temps 16, noté également d, comporte une fibre optique de détection 32 reliée à un module de détection résolu en temps 34. Lorsqu'un tel détecteur fibré est utilisé, l'extrémité libre de cette fibre optique de détection 32 est assimilée au détecteur d. Dans l'exemple de réalisation de la figure 1, le module de détection 34 comporte un organe de détection 36 pour chaque détecteur 16. En variante non représentée, le module de détection 34 comporte un organe de détection commun à plusieurs détecteurs 16, en particulier un unique organe de détection pour l'ensemble des détecteurs 16. L'organe de détection 36 est, par exemple, un photomultiplicateur, ou une photodiode à avalanche, également appelée APD (de l'anglais Avalanche Photodiode), ou une photodiode à avalanche à photon unique, également SPAD (de l'anglais Single-Photon Avalanche Diode) ou encore un tube intensificateur d'images ayant une ou plusieurs anodes (de l'anglais Multi-Channel Plate). La sonde est, dans cet exemple, une sonde compacte pour le diagnostic de certains cancers, tels qu'une sonde portative pour le diagnostic du cancer du sein, une sonde endorectale pour le diagnostic du cancer de la prostate, ou encore une sonde de dermatologie. En variante, la sonde est une sonde d'endoscope telle qu'une sonde souple pour le diagnostic du cancer digestif. En variante, la sonde est un dispositif applicable à d'autres organes. L'unité de traitement de l'information 18 comporte un processeur de données 38 et une mémoire 40 associée au processeur 38. L'unité de traitement 18 comporte des moyens 42 d'élaboration, pour au moins un couple source 14 - détecteur 16, d'une première distribution Bsd d'un signal reçu par le détecteur 16 correspondant pour le milieu à caractériser MC, le signal reçu étant émis par le milieu diffusant à caractériser MC suite à l'éclairement dudit milieu MC par la source 14 correspondante. La source 14, et respectivement le détecteur 16, correspondants sont repérés par les indices s, et respectivement d. Dans le mode de réalisation décrit, la première distribution est une distribution temporelle respectivement notée Bsd(t). Les moyens d'élaboration 42 sont, de préférence, réalisés sous forme d'une ou plusieurs cartes électroniques connectées au détecteur 16 et permettant de mesurer la distribution temporelle du temps d'arrivée des photons. Chaque source de lumière pulsée 28 comporte un laser impulsionnel. En variante, chaque source lumineuse pulsée 28 comporte une diode laser. En variante encore, chaque source lumineuse pulsée 28 comporte une source de lumière constante dont l'intensité lumineuse est modulée en des impulsions de durée équivalente par un dispositif à obturation rapide. Ainsi, la source est apte à émettre un rayonnement lumineux prenant la forme d'une impulsion temporelle. Selon ce mode particulier, l'extrémité de chaque fibre d'excitation 30 s'étend perpendiculairement à la surface d'observation 22 c'est-à-dire selon l'axe Z, ou encore de manière oblique par rapport à l'axe Z afin d'émettre la lumière en biais par rapport à la surface d'observation 22. La mémoire 40 est apte à stocker un logiciel 44 de calcul d'une pluralité de transformées de Mellin-Laplace Ysd(P'n), pour une variable p, également appelée largeur, et des valeurs distinctes d'un ordre n donnés, d'une grandeur Ysd(t) comportant la première distribution 13,d(t), l'ordre n étant un nombre entier et la variable p étant un nombre réel. La transformée de Mellin-Laplace d'ordre n et de largeur p d'une fonction f est notée f(P'n), et vérifie l'équation suivante : 1 f (") = - f (pt)n .exp(-pt).f (t).dt (1) n! où 1/n! représente un terme de normalisation. La transformée de Mellin-Laplace f(P'11) de la fonction f s'écrit également : +00 f (") = f W(") f (t).dt (2) -00 où w(P'n) représente une fenêtre temporelle définie par l'équation suivante : w (p,n) H (t).(pt)n .exp(- pt) (3) n! avec H(t) représentant la fonction connue de Heaviside. Sur la figure 2, la fonction f est représentée par la courbe 70 en trait gras, et les fenêtres temporelles W(P'n) sont représentées par les courbes 71 à 80, pour p égal à 2 ns-1 et pour l'ordre n prenant respectivement les valeurs entières successives comprises entre 0 et nmax, avec nmax égal à 9. Les positions ((n+1)/p, pf(P'n)) sont repérés par les points 81 à 90, pour les mêmes valeurs de la largeur p et de l'ordre n. L'aire des fenêtres temporelles w(P'n) est égale à 1/p, et la position moyenne de ces fenêtres temporelles W(P'n) est égale à (n+1)/p. De manière générale, la transformée de Mellin-Laplace f(P'11) d"ordre n et de largeur p de la fonction f sera considérée comme toute fonction comprenant le terme +o. f (pt)n .exp(- pt). f (t).dt à un coefficient multiplicatif près. 0 Dans le mode de réalisation décrit, au moins une transformée de Mellin-Laplace Ysd(P'n) est calculée pour une valeur de l'ordre n supérieure ou égale à 5, de préférence supérieure ou égale à 8. Dans le mode de réalisation décrit, l'ordre n présente une valeur comprise entre 0 et une valeur maximale nmax. L'ordre n prend, par exemple, successivement toutes les valeurs comprises entre 0 et nmax. Autrement dit, nmax+l transformées de Mellin-Laplace Ysd(P'n) sont alors calculées par le logiciel de calcul 44. La valeur maximale nmax est supérieure ou égale à 5, de préférence supérieure ou égale à 8. La combinaison de transformées de Laplace d'ordre différents, et notamment des faibles ordres, c'est-à-dire n proche de 0, avec des ordres importants, c'est-à-dire n 5, permet la prise en compte de la contribution, au signal détecté, de photons émis à des instants très proches après l'impulsion d'excitation (ordres faibles) ainsi qu'à des instants plus éloignés (ordres importants), ces derniers étant appelés photons retardés, car ils sont détectés quelques ns après l'impulsion lumineuse d'excitation. Dans le mode de réalisation décrit, la valeur de la largeur p est comprise entre 1 ns-1 et 20 ns-1, où ns est le symbole de la nanoseconde. La mémoire 40 est également apte à stocker un logiciel 46 de reconstruction de propriétés optiques du milieu 12 à partir d'une combinaison de la pluralité de transformées de Mellin-Laplace Ysd(P'n) de ladite grandeur Ysd (t). Le logiciel de reconstruction 46 est propre à reconstruire de la distribution spatiale de coefficients optiques, tels que le coefficient d'absorption optique pl, ou le coefficient de diffusion D. En variante, le logiciel de reconstruction 46 est propre à reconstruire une carte de fluorescence du milieu 12. Dans le mode de réalisation décrit, la grandeur Ysd(t) est un signal corrigé permettant de prendre en compte la réponse instrument, également appelé distribution corrigée, et est décrite plus en détail ci-après. En variante, la grandeur sur laquelle est calculée la pluralité de transformées de Mellin-Laplace est la première distribution Bsd(t). En variante, les moyens de calcul 44 et les moyens de reconstruction 46 sont réalisés sous forme de composants logiques programmables ou encore sous forme de circuits intégrés dédiés. Le fonctionnement du système d'imagerie 10 selon l'invention est décrit ci-après à l'aide de la figure 3 représentant un organigramme du procédé de traitement d'image selon l'invention. Lors de l'étape 100, le milieu à caractériser MC est éclairé par la ou les sources 14 correspondantes. En cas d'une pluralité de sources 14, l'éclairement du milieu 12 est effectué successivement pour chacune des sources 14. Lorsqu'il est éclairé, le milieu à caractériser MC émet en retour un signal optique de diffusion, susceptible d'être détecté par le ou chaque détecteur 16, et la première distribution du signal reçu 13,d(t) est élaborée, lors de l'étape 110, pour chaque couple source 14 - détecteur 16, également noté (s, d), à l'aide des moyens d'élaboration 42. La première distribution temporelle Bsd(t) vérifie l'équation suivante : B sd (t) = IRF sd (t) * G sd (t) (4) où IRF,d(t) représente la réponse instrument, c'est-à-dire l'influence de la source s et du détecteur d sur la première distribution élaborée, et Gsd(t) représente une première fonction de modélisation d'un signal de diffusion de la lumière entre la source 14 et le détecteur 16 dans le milieu à caractériser MC. On note Ss(t) la réponse temporelle de la source 14, également notée s. Dans cet exemple, la source 14 est une brève impulsion lumineuse, souvent modélisée par une distribution de type Dirac ponctuelle. Ss(t) représente l'intensité temporelle du signal émis. Lorsque la source 14 est fibrée (source de lumière couplée à une fibre optique), c'est l'extrémité de la fibre optique d'excitation 30, qui est alors considérée comme la source, comme précédemment indiqué. On note Dd(t) la réponse temporelle du détecteur 16, également noté d. Cette réponse traduit notamment le délai entre l'arrivée d'un photon sur le détecteur 16 et sa détection effective. Lorsque le détecteur 16 est fibré (détecteur couplé à une fibre optique), c'est l'extrémité de la fibre de détection 32, qui est considérée comme le détecteur, comme précédemment indiqué. On estime généralement la réponse instrument, notée IRF,d(t) l'indice sd corrrespondant à un couple source d'excitation s et détecteur d. Cette réponse instrument est obtenue en combinant Ss(t) et Dd(t) suivant : IRF,d(t) = Ss(t)*Dd(t), où * désigne le produit de convolution. Lorsqu'on a plusieurs sources et plusieurs détecteurs, on détermine IRF,d(t) pour chaque couple source-détecteur, car chaque source et chaque détecteur a sa propre réponse. En pratique, cela se fait en plaçant une source (ou une extrémité de fibre optique constituant la source) en face d'un détecteur (ou d'une extrémité de fibre optique constituant le détecteur). Cela suppose un alignement soigné entre la source et le détecteur, et cela prend du temps, notamment lorsqu'on augmente le nombre de sources et de détecteurs. A titre d'exemple, 100 fonctions de réponse instrument IRF,d sont à déterminer lorsqu'on dispose de 10 sources et 10 détecteurs. Lors de l'étape 120, la réponse instrument IRF,d(t) est prise en compte à l'aide d'une deuxième distribution temporelle Asd(t) et d'une deuxième fonction de modélisation Gsd°(t). La deuxième distribution temporelle Asd(t) est élaborée pour un milieu de référence MR et dans les mêmes conditions expérimentales que pour la première distribution temporelle 13,d(t). Le milieu de référence MR, également appelé fantôme, est connu en soi et est, par exemple, décrit dans le document FR 2 950 241 A1. Le fantôme MR a des caractéristiques optiques connues. Le fantôme MR est, par exemple, un fantôme ayant des coefficients d'absorption pa° et de diffusion D° de valeurs médianes pour les tissus humains, par exemple égales à 0,2 cm-1, et respectivement à 10 cm-1. L'élaboration de la deuxième distribution temporelle Asd(t) est de préférence effectuée de manière préliminaire, avant l'étape 100, et le milieu de référence MR est alors remplacé par le milieu à caractériser MC au début de l'étape 100. En variante, l'élaboration de la deuxième distribution temporelle Asd(t) est effectuée lors de l'étape 120, ce qui nécessite de changer le milieu 12, et de remplacer le milieu à caractériser MC par le milieu de référence MR. La deuxième distribution temporelle Asd(t) vérifie alors l'équation suivante : A $d (t) = IRF sd (t)* G1 (t) (5) Les première et deuxième fonctions de modélisation Gsd(t), Gsd°(t) sont, par exemple, des fonctions de Green, bien connues dans le domaine la diffusion optique. Chaque fonction de Green Gsd(t), Gsd°(t) représente la densité de photons au niveau du détecteur 16 lorsque le milieu 12, c'est-à-dire le milieu de référence MR pour la première fonction de modélisation et respectivement le milieu à caractériser MC pour la deuxième fonction de modélisation, est éclairé par la source 14, où s et d sont les indices désignant respectivement la source 14 et le détecteur 16. Le signal corrigé Ysd(t) est alors calculé en fonction de la première distribution 13,d(t), de la première fonction de modélisation Gsd(t), de la deuxième distribution Asd(t) et de la deuxième fonction de modélisation Gsd°(t). Le signal corrigé Ysd(t) est le résultat d'une opération de comparaison du produit de la deuxième distribution Asd et de la première fonction de modélisation Gsd avec le produit de la première distribution Bsd et de la deuxième fonction de modélisation Gsd°. L'opération de comparaison comporte une opération arithmétique du produit de la deuxième distribution Asd et de la première fonction de modélisation Gsd avec le produit de la première distribution Bsd et de la deuxième fonction de modélisation Gsd°. L'opération arithmétique est, par exemple, une soustraction du produit de la deuxième distribution Asd et de la première fonction de modélisation Gsd du produit de la première distribution Bsd et de la deuxième fonction de modélisation Gsd°. Dans le mode de réalisation décrit, le produit de la deuxième distribution Asd(t) et de la première fonction de modélisation Gsd(t) et le produit de la première distribution Bsd(t) et de la deuxième fonction de modélisation Gsd°(t) sont des produits de convolution, les première Bsd(t) et deuxième Asd(t) distributions étant des distributions temporelles. Le signal corrig é Ysd(t) vérifie l'équation suivante : Y sd (t) - Bsd (t) * G1 (t) Asd (t) * Gsd (t) (6) En variante, l'opération arithmétique est un ratio entre le produit de la première distribution temporelle Bsd(t) et de la deuxième fonction de modélisation Gsd°(t) et le produit de la deuxième distribution temporelle Asd(t) et de la première fonction de modélisation Gsd(t). Selon cette variante, le signal corrigé Ysd(t) vérifie alors l'équation suivante : Bsd (t) * Gsd (t) Y sd (t) = (7) Asd (t) * Gsd (t) En complément, la première fonction de modélisation Gsd(t) est déterminée de manière approchée à l'aide de la fonction de Green Gsdl(t) pour le milieu à caractériser MC, et on obtient alors les équations suivantes : Y sd (t) = Bsd (t)* Gsd (t) - Asd (t)* Gsd (t) = IRF sd (t)* G sd (t)* G(s)d (t) - IRF sd (t)* G(s)d (t)* Gisd (t) = I RF sd (t) * G`s'd (t) * (Gsd (t) - Gsd (t)) (8) = Asd (t) * (Gsd (t) - Gsd (t)) (I. Gis (i- ,t) * cîu a (i- ,t) * Gid (i- , t)cli + I r> Gis (i- ,t)* SD (i- ,t) * r Gid (i- ,t)cli) = - A sd (t)* où pa et D représentent respectivement le coefficient d'absorption et le coefficient de diffusion du milieu à caractériser MC. IRFsd(t) représente la réponse instrument, Gs(r,t) = G(rs,r,t) est la fonction de Green représentant la densité de photons à un endroit r du milieu MC lorsque le milieu MC est éclairé par la source s située en rs, Gd(r,t) = G(rd,r,t) est la fonction de Green représentant la densité de photons à l'endroit r du milieu MC lorsque le milieu à caractériser MC est éclairé par la source s située en rd. On écrit également Gd(r,t) = G(r,rd,t). Ainsi, cette fonction de Green représente également la densité de photons au niveau du détecteur (rd) lorsque le milieu à caractériser MC est éclairé par une source située en r. Autrement dit, Gsd(t) correspond aux vraies valeurs pa(r, t) et D(r, t), c'est-à-dire aux valeurs recherchées, et Gsdl(t) correspond à des valeurs approchées pla(r, t) et D1(r, t), avec pa(r, t) = pla(r, t) + Elpa(r, t) et D(r, t) = D1(r, t) + ÔD(r, t). Le processus de reconstruction des propriétés optiques est généralement itératif, de telle sorte qu'au cours de chaque itération, on obtient des valeurs pla(r) et D1(r), auxquelles correspondent des fonctions de Green G si d (t) pour chaque couple source-détecteur sd. L'objectif du procédé de reconstruction est alors de minimiser Elpa(r, t) et ÔD(r, t), ce qui permet de faire tendre pla(r, t) et D1(r, t) respectivement vers pa(r, t) et D(r, t). Les coefficients pa et D sont constants dans le temps, et les valeurs pa(r, t) = pa(r)a(t) et D(r,t) = D(r)3(t), 3(t) représentent une distribution de Dirac à l'instant t. Le calcul du signal corrigé Ysd(t), également appelé distribution temporelle corrigée, s'affranchit donc de la connaissance précise de la réponse instrument IRF,d(t) d'après les équations (6), (7) ou encore (8). Ceci permet d'éviter des mesures expérimentales relativement longues et fastidieuses afin de déterminer la réponse instrument du système d'imagerie 10. En variante, la réponse instrument IRFsd(t), représentant la répartition temporelle de l'impulsion générée par la source de rayonnement et détectée par le détecteur en l'absence de milieu diffusant, est mesurée, pour chaque couple source-détecteur (s, d), lors d'une opération de calibrage au cours de laquelle la source s est placée face au détecteur d en l'absence du milieu 12. Lors de l'étape 130, les transformées de Mellin-Laplace Ysd(P'n) de la distribution corrigée sont alors calculées pour différentes valeurs de l'ordre n comprises entre 0 et nmax, de préférence pour toutes les valeurs successives comprises entre 0 et n'x, suivant l'équation suivante : Y(sdP') (Bs(l:ii).Gs°:/ P..1) - A( "'il) .G(';i1)) (9) + j=n Or, la distribution corrigée Ysd(t) vérifie l'équation suivante : Y sd (t) = - A sd (t) * Gis , t) * (5,1 , t) * Gd t)cli + J OGi ,t) * e(i-,t)*Gld(i-,t)e)(io) Les équations (9) et (10) permettent alors d'obtenir l'équation suivante, étant donné que la transformée de Mellin-Laplace d'un produit de convolution de deux termes est égale au produit des transformées de Mellin-Laplace de chacun de ces deux termes : Y(") = IA(i:ik).(.1.Gsl(P) (i).4,(7.).Gdi(")(i)cli + f V>"Gsl(P.I) (i).5D(i).V"Gd1(")(i)di) (1 1 ) j+k=n L'équation (11) permet ainsi de calculer la transformée de Mellin-Laplace à l'ordre n de la distribution corrigée Ysd(P'n) des transformées de Mellin-Laplace de la deuxième distribution temporelle Asd(t) et des fonctions de Green Gsl(r), Gd1(r) approchant la première fonction de modélisation Gsd(t). Le milieu 12 est discrétisé en une pluralité M de voxels référencés m, avec m compris entre 1 et M, afin de calculer les transformées de Mellin-Laplace des fonctions de Green Gsl(r), Gd1(r). Les fonctions de Green Gsl(r), Gd1(r) discrétisées pour chaque voxel m sont notées Gsl(rm), Gdi (rm). Selon la variante dans laquelle la réponse instrument IRFsd(t) est mesurée pour chaque couple source-détecteur (s, d), lors d'une opération de calibrage, les transformées de Mellin-Laplace de la première distribution Bsd(P'n) et de la réponse instrument IRFsd(P'n) sont calculées pour différentes valeurs de l'ordre n comprises entre 0 et nmax, de préférence pour toutes les valeurs successives comprises entre 0 et nmax. La reconstruction de propriétés optiques du milieu 12 est alors effectuée lors de _ l'étape 140, afin de déterminer par exemple les vecteurs fia et D , dont les termes pa(m) et D(m), discrétisés pour chaque voxel m, c'est-à-dire ma (i- = m) et D (i- = m) , constituent les cartes des propriétés optiques recherchées. L'étape de reconstruction vise alors à résoudre le système matriciel suivant : Y =W X , (12) où Y représente un vecteur des observations, W représente une matrice de passage, et X représente un vecteur des inconnues. Le vecteur des observations Y contient les Nn transformées de Mellin-Laplace de la distribution corrigée Ysd(P'n), avec Nn égal à nmax+l, calculées pour toutes les valeurs successives de l'ordre n comprises entre 0 et nmax et pour les couples source-détecteur (s, d) considérés. On définit, pour chaque triplet (s, d, n) des indices s, d et de l'ordre n, un indice I de la manière suivante : / = (s -1)x Nd xNn +(d -1)xNn+n (13) où Nd est le nombre de détecteurs 16 considérés. La valeur maximale de l'indice I est égale à Imax, tel que : I. = NsxNdxNn (14) Le vecteur des observations Y comporte alors Imax lignes. La matrice de passage W comprend une première partie WP comportant des premiers termes notés WP (/,m) et une deuxième partie WD comportant des deuxièmes termes notés WD(/,m). Cette matrice est appelée matrice de sensibilité des mesures Y aux propriétés optiques X. Les premiers termes WP (/,m) vérifient l'équation suivante : (I,m) IA(d,k).Gsi(p,i)(i.m).Gdi(pd)(i.m)vm (15) i+ j+k=n et les deuxièmes termes W D (I, M) vérifient l'équation suivante WD(/,m)= -~A~d'k).OGS(P,i)(Ym).OGd(P,j)(~m)~m (16) j+k=n où V, représente un élément de volume entourant la maille m, tel que le volume de la cellule de Voronoi. Autrement dit, V, représente le volume du voxel m. Selon ce mode de réalisation, au cours de chaque itération, les transformées de Mellin Laplace, de différents ordres, de Gsl(t) et de Gd1(t) sont déterminées, notamment en utilisant l'expression suivante : ( 5(r) si n = 0 -VD(i)VG,l(P')(i)+ ma (i)-E Gsl(P')(i)P (17) c s Gi(p,no(i) sin> 0 Ainsi, lors de chaque itération, on détermine Gsl(P.n-°)puis on détermine les transformées d'ordre supérieur de façon itérative, en utilisant l'expression ci-dessus lorsque n > O. Lors de la première itération, on utilise D = D° et pl, = Ita°, et lors des itérations suivantes, on utilise D = D1 et pl, = pal. Le même raisonnement s'applique naturellement pour la détermination des transformées de Mellin Laplace de la fonctionGd1(t). La matrice de passage W comporte alors Imax lignes et 2M colonnes. Chaque terme WP (/,m), WD(/,m) de la matrice de passage W est déterminé par modélisation, à chaque itération, en fonction des propriétés optiques déterminées lors de l'itération précédente, ou lors de la première itération, en fonction de propriétés optiques initialisées par l'opérateur. Cette initialisation est réalisée en considérant par exemple un milieu 12 homogène. Le vecteur des inconnues X comporte 2M lignes et une colonne, et contient les inconnues pa(m) et D(m) pour chacun des M voxels. Les M premières lignes correspondent aux inconnues pa(m) et les M lignes suivantes du vecteur X correspondent aux inconnues D(m). L'inversion de la matrice de passage W est réalisée selon des algorithmes d'inversion bien connus de l'homme du métier, tels qu'un algorithme de descente de gradient, une technique de reconstruction algébrique (ART), une décomposition en valeurs singulières (SVD), ou encore une méthode des gradients conjugués. Le processus s'arrête lorsqu'un critère de convergence est atteint, par exemple lorsque la distance entre deux vecteurs X successifs est inférieure à un seuil prédéterminé. La figure 4 est un ensemble d'images 200 à 208 représentant une sensibilité Jsd(P'n) définie par l'équation : "sd (p,n) ,^.) am ci) pour les transformées de Mellin-Laplace de l'ensemble des distributions temporelles corrigées Ysd(P'n), avec 21 sources 14 et 21 détecteurs 16, symbolisés respectivement par une flèche descendante et une flèche montante. L'espace entre une source 14 et un détecteur 16 successifs est d'environ un millimètre. Chaque image 200 à 208 correspondant à une valeur respective et croissante de l'ordre n comprise entre 0 et 8, l'image 200 étant associée à l'ordre 0 et l'image 208 étant associée à l'ordre 8. Sur chaque image 200 à 208, l'ensemble des distributions temporelles corrigées Ysd(P'n) correspond à une zone claire 210 en forme de banane. L'homme du métier remarquera alors que plus la valeur de l'ordre n est élevée, meilleure est la sensibilité dans une zone profonde du milieu 12, par exemple à une profondeur supérieure à un centimètre. Autrement dit, plus la valeur de l'ordre n est élevée, plus le système de reconstruction 10 suivant le procédé selon l'invention est apte à observer la zone profonde du milieu 12. Le système de reconstruction 10 et le procédé de reconstruction selon l'invention permettent de localiser plus précisément des absorbeurs disposés en profondeur qu'un système d'imagerie et un procédé de reconstruction de l'état de la technique, comme le montrent les figures 5 et 6. La figure 5 est une vue d'une image d'un milieu comportant deux absorbeurs disposés à une profondeur de 1 cm, et respectivement de 2 cm, repérés par les cercles 220, et respectivement 230, reconstruite avec un système d'imagerie et un procédé de reconstruction de l'état de la technique, alors que la figure 6 est une vue d'une image du même milieu reconstruite avec un système d'imagerie et un procédé de reconstruction selon l'invention. En comparant les figures 5 et 6, l'homme du métier notera que les absorbeurs sont localisés avec une meilleure précision dans le cas du système d'imagerie et du procédé de reconstruction selon l'invention, notamment pour l'absorbeur disposé à 2 cm de profondeur. (18) La reconstruction des propriétés optiques à l'aide de la combinaison d'une pluralité de transformées de Mellin-Laplace de la grandeur comportant la distribution élaborée, en particulier à l'aide de la combinaison des nmax+l transformées de Mellin-Laplace Ysd(P'n), permet de converger plus rapidement vers la solution, afin de réduire le temps de traitement nécessaire à la reconstruction des propriétés optiques du milieu. Un autre avantage d'une telle combinaison est qu'elle permet de prendre en compte les photons diffusés détectés peu après l'impulsion lumineuse, c'est-à-dire jusqu'à 1 ou 2 ns après l'impulsion lumineuse, également appelés photons prompts, en utilisant une transformée de Mellin Laplace aux faibles ordres, telles que des valeurs de l'ordre n comprises entre 0 et 5, tout en prenant également en compte les photons dits retardés, c'est-à-dire détectés quelques ns, voire plus de 5ns, après l'impulsion lumineuse. En utilisant des transformées de Mellin Laplace de différents ordres pour réaliser la reconstruction, on prend en compte à la fois les photons dits prompts et les photons dits retardés, ces derniers provenant généralement de profondeurs plus importantes. Il est donc avantageux d'effectuer la reconstruction en utilisant des données d'entrées combinant à la fois des transformées de Mellin Laplace de différents ordres. Par données d'entrée, on entend les données observées, c'est-à-dire les termes constituant le vecteur des observations Y . On décrit à présent une variante selon laquelle on détermine expérimentalement la réponse instrument IRF,d(t) correspondant chaque couple source-détecteur. La reconstruction s'effectue selon un procédé itératif, en réalisant une mesure 13,d(t) pour un ou différents couple(s) source-détecteur. Pour chaque couple source-détecteur, on obtient alors l'équation suivante : Bsd(t)- Bsd(t)= S s(t)*(f Gis ,t)* êtta(r, t) * Gid(r, + Gis(r, t)* , t) * Gld(r, t)cli)* Dd (t) (19) où : - B si d représente une mesure modélisée en considérant que le milieu présente les coefficients d'absorption pal(r) et de diffusion Dal(r) déterminés lors de l'itération précédente, les valeurs pal(r) et Dal(r) étant initialisées lors de la première itération par une première estimation réalisée par l'opérateur, - les coefficients 4u,(i,t) et ÔD(r,t) représentent la différence entre les propriétés optiques du milieu (pa, D) et les propriétés optiques (pal, D1) initiales, ou bien résultant de l'itération précédente. De manière analogue au mode de réalisation précédent, l'objectif est de minimiser Spa(i,t) et ÔD(r,t) . Si 13,d(t) est une mesure réalisée sur le milieu à caractériser MC pour un couple source détecteur sd donné, on détermine les transformées de Mellin Laplace d'une fonction Ysd , ces transformées se présentant sous la forme : B(Pr) -In 1 IRF(P'n-i)B(P y (p,n) s,d s," t=0 s,d s,d I(P'()) s,d Ainsi, en fonction des transformées de Mellin-Laplace de la réponse instrument /RFs(d°k), avec 0 k n, et de la mesure Bs(':/k) avec 0 k n, on détermine 17:,7) : ( f IGs1(")(i).41,(i).G1(dP'k)(i-)di- + IGsl(")(i).5D(i-).G1(dP'k)(i-)di- (21) ~çZj+k=n j+k=n Après discrétisation du milieu en M voxels m, cette expression devient : ( IGs1(")(im).,u,(im).G1(dP'k)(i'm)Vm + /Gs1(")(im).bp(im).G1(dP'k)(im)Vm (22) m=1 \, j+k=n j+k=n On met alors en oeuvre un algorithme d'inversion itératif pour définir les propriétés optiques minimisant les vecteurs Spa (im) et bD(rm), à partir d'une pluralité de transformées de Mellin Laplace Y(d°n) , avec des ordres différents, obtenues à l'aide de mesures réalisées avec une pluralité de couples sources détecteurs. Selon un deuxième mode de réalisation, on utilise la transformée de Mellin Laplace pour établir une carte de fluorescence F(r) , avec F (i) a(i-)t ,17 étant le rendement de fluorescence du fluorophore localisé dans le volume di- . On définit une fonction de réponse du fluorophore F(r,t), telle que F(i,t)- F(i)exp(--t ) (23) avec i désignant le temps de vie du fluorophore. A l'aide d'une source d'excitation s et d'un détecteur d tels que ceux définis dans le premier mode de réalisation, on réalise une mesure temporelle Bsd(t) sur le milieu à caractériser MC. Cette mesure s'écrit : B sd (t) /iffsd (t) * L G s(i , t) * F (i- ,t)* Gd (i- ,t).di (24) où Gs(i,t) et Gd (i" ,t) sont des fonctions de Green représentant la densité de photons respectivement entre la source s et le point r du milieu, et entre le point r du milieu et le détecteur d. (20) v(,n) sdp 175(1;1'n) -A M 18 (1,,n)- pn n exp(-t12) = t .dt =2 (p n! Comme (25) nous obtenons : F(") F (i)(pr)n (26) La transformée de Mellin-Laplace de la distribution temporelle mesurée Bsd est donc : B (n) = MF s( d P .(p .LGs(p,k) ) F(r) . G s( '1) ) . d (27) i+ j+k+1=n Après discrétisation du milieu en M voxels m, on obtient : M B5(7) = IRF .(p ( s(P'k) (7,').G s(,1'1) (28) j+k+1=n m=1 On inverse enfin ce système après discrétisation du milieu en voxels. On met alors en oeuvre un algorithme d'inversion, tel qu'un des algorithmes connus mentionnés précédemment, pour définir le vecteur F (im) correspondant à la carte de fluorescence recherchée, à partir de la pluralité de transformées de Mellin Laplace Y»'), avec des ordres différents, obtenues à l'aide de mesures réalisées avec une pluralité de couples sources détecteurs sd. De manière analogue au premier mode de réalisation, la combinaison de différentes transformées de Mellin Laplace d'une grandeur Ys(7)correspondant à un couple source détecteur, et notamment le recours à différentes valeurs de l'ordre n, permet d'obtenir une meilleure prise en compte les contributions des photons prompts et des photons retardés, ces derniers provenant essentiellement de profondeurs plus élevées. REVENDICATIONS1.- Procédé de reconstruction de propriétés optiques d'un milieu (12), à l'aide d'un système de reconstruction (10) comportant au moins une source de rayonnement (14, s) propre à éclairer le milieu (12) et au moins un détecteur (16, d) propre à recevoir un signal émis par le milieu (12), le procédé comprenant les étapes suivantes : - l'éclairement (100) du milieu (12) par la ou chaque source de rayonnement (14, s), - la réception, par le ou chaque détecteur (16, d), d'un signal émis par le milieu (12), et - l'élaboration (110), pour au moins un couple source-détecteur (14, s, 16, d), d'une première distribution (Bsd(t)) du signal reçu par le détecteur (16, d) correspondant, - le calcul (130) de la transformée de Mellin-Laplace (Yscl(P'n) ; Bsd(P'n)), pour un ordre (n) et une variable (p) donnés, d'une grandeur (Ysd(t) ; Bsd(t)) comportant la première distribution (Bsd(t)), l'ordre (n) étant un nombre entier, la variable (p) étant un nombre réel, et - la reconstruction (140) de propriétés optiques du milieu (12) à l'aide de la transformée de Mellin-Laplace (Ysd(P'n) ; Bsd(P'n)) de ladite grandeur (Ysd(t) ; Bsd(t)), caractérisé en ce que l'étape de calcul (130) comporte le calcul d'une pluralité de transformées de Mellin-Laplace (Ysd(P'n) ; Bsd(P'n) de ladite grandeur pour des valeurs distinctes de l'ordre (n), et en ce que l'étape de reconstruction (140) est effectuée à partir d'une combinaison de la pluralité de transformées de Mellin-Laplace (Ysd(P'n) ; Bsd(P'n)). 2.- Procédé selon la revendication 1, dans lequel l'étape de calcul (130) comporte le calcul d'au moins une transformée de Mellin-Laplace (Ysd(P'n) ; Bsd(P'n)) de ladite grandeur pour une valeur de l'ordre (n) supérieure ou égale à 5, de préférence supérieure ou égale à 8. 3.- Procédé selon la revendication 2, dans lequel l'étape de calcul (130) comporte le calcul de n,'+1 transformées de Mellin-Laplace (Ysd(P'n) ; Bsd(P'n), l'ordre (n) prenant successivement toutes les valeurs comprises entre 0 et nmax, avec rima, supérieur ou égal à 5, de préférence supérieur ou égal à 8. 4.- Procédé selon l'une quelconque des revendications précédentes, dans lequel la valeur de la variable (p) est comprise entre 1 ns-1 et 20 ns-1. 5.- Procédé selon l'une quelconque des revendications précédentes, dans lequel le procédé comprend en outre :- la détermination, pour au moins un couple source-détecteur (14, s, 16, d), d'une première fonction (Gsd(t)) de modélisation d'un signal de diffusion de la lumière entre la source (14, s) et le détecteur (16, d) dans le milieu (12), - l'élaboration, pour ledit au moins un couple source-détecteur (14, s, 16, d), d'une deuxième distribution (Asd(t)) d'un signal reçu par le détecteur (16, d) pour un milieu de référence (MR), le signal reçu étant émis par le milieu de référence (MR) suite à l'éclairement dudit milieu (MR) par la source (14, s), et - la détermination, pour ledit au moins un couple source-détecteur (14, s, 16, d), d'une deuxième fonction (Gsd°(t)) de modélisation d'un signal de diffusion de la lumière entre la source (14, s) et le détecteur (16, d) dans le milieu de référence (MR), et dans lequel ladite grandeur (Ysd(t)) dépend de la première distribution (Bsd(t)), de la première fonction de modélisation (Gsd(t)), de la deuxième distribution (Asd(t)) et de la deuxième fonction de modélisation (Gsd°(t)). 6.- Procédé selon la revendication 5, dans lequel ladite grandeur (Ysd(t)) est obtenue par soustraction du produit de la deuxième distribution (Asd(t)) et de la première fonction de modélisation (Gsd(t)) du produit de la première distribution (Bsd(t)) et de la deuxième fonction de modélisation (Gs2(t)). 7.- Procédé selon la revendication 5, dans lequel ladite grandeur (Ysd(t)) est le ratio entre le produit de la première distribution (Bsd(t)) et de la deuxième fonction de modélisation (Gs2(t)) et le produit de la deuxième distribution (Asd(t)) et de la première fonction de modélisation (Gsd(t)). 8.- Procédé selon l'une quelconque des revendications précédentes, dans lequel la ou chaque source de rayonnement (14, s) comporte une source de lumière pulsée. 9.- Procédé selon l'une quelconque des revendications précédentes, dans lequel le ou chaque détecteur (16, d) est un détecteur résolu en temps. 10.- Procédé selon l'une quelconque des revendications précédentes, dans lequel les propriétés optiques comprennent au moins un élément parmi le groupe consistant en : - les propriétés d'absorption de la lumière, caractérisées notamment par le coefficient d'absorption (14), - les propriétés de diffusion, caractérisées notamment par le coefficient de diffusion réduit (ii's) ou par le coefficient de diffusion (D), et - les propriétés de fluorescence, caractérisées notamment par une fonction de réponse (F) d'un fluorophore, ou par une concentration (c) du fluorophore, ou encore par une grandeur fonction d'une quantité (q) du fluorophore. 11.- Système (10) de reconstruction de propriétés optiques d'un milieu (12), comprenant :- au moins une source (14, s) de rayonnement propre à éclairer le milieu (12), - au moins un détecteur (16, d) propre à recevoir un signal émis par le milieu (12), - des moyens (42) d'élaboration, pour au moins un couple source-détecteur (14, s, 16, d), d'une distribution (Bsd(t)) du signal reçu par le détecteur (16, d) correspondant, - des moyens (44) de calcul de la transformée de Mellin-Laplace (Ysd(P'n) ; Bsd(P'n)), pour un ordre (n) et une variable (p) donnés, d'une grandeur (Ysd(t) ; Bsd(t)) comportant la distribution (Bsd(t)); l'ordre (n) étant un nombre entier, la variable (p) étant un nombre réel, et - des moyens (46) de reconstruction de propriétés optiques du milieu (12) à l'aide de la transformée de Mellin-Laplace (Ysd(P'n) ; Bsd(P'n)) de ladite grandeur (Ysd(t) ; Bsd(t)), caractérisé en ce que les moyens de calcul (44) sont propres à calculer une pluralité de transformées de Mellin-Laplace (Ysd(P'n) ; Bsd(P'n) de ladite grandeur pour des valeurs distinctes de l'ordre (n), et en ce que les moyens de reconstruction (46) sont propres à reconstruire les propriétés optiques du milieu (12) à partir d'une combinaison de la pluralité de transformées de Mellin-Laplace (Ysd(P'n) ; Bsd(P'n). 12.- Système (10) selon la revendication 11, dans lequel les moyens de calcul (44) sont propres à calculer au moins une transformée de Mellin-Laplace (Ysd(P'n) ; Bsd(P'n)) de ladite grandeur pour une valeur de l'ordre (n) supérieure ou égale à 5, de préférence supérieure ou égale à 8.20