Traitement du signal hyperspectral

De nombreux minéraux, ainsi que les assemblages minéralogiques correspondants, ont été caractérisés en laboratoire. L’intérêt de la spectroscopie de réflectance est de pouvoir effectuer la démarche inverse, à savoir retrouver la composition et les propriétés de la roche à partir d’un spectre. Comme cela a été décrit dans la partie précédente, il existe une grande variété de composition chimique pour un minéral donné et d’autant plus pour les possibilités d’assemblages minéralogiques dans les roches. Il est donc nécessaire d’utiliser des techniques de traitement pour déconvoluer les données. Quelques-unes de ces différentes techniques seront brièvement présentées dans cette dernière partie pour donner un aperçu des possibilités existantes.

1 Différentes techniques de déconvolution spectrale

1.1 Première approche : les traitements ” exploratoires ”

Le développement des techniques de déconvolution (synthétisées notamment par Mustard et Sunshine, 1999) est étroitement lié au développement des résolutions spectrales des appareillages de mesures. Le nombre de canaux disponibles, c’est-à-dire le nombre de longueurs d’onde analysées, va fortement conditionner les traitements possibles pour un spectre.

Classifications statistiques

L’objectif principal des traitements statistiques est de réduire la dimension des données hyperspectrales afin d’accéder efficacement à la partie intéressante de l’information. L’analyse en composante principale (ACP, ou Principal Component Analysis en anglais), permet d’extraire les spectres les plus extrêmes dans un jeu de données mais ne donne aucunement accès à l’information minéralogique contenue. Cette technique a par exemple été utilisée pour Mars sur les données ISM par Murchie et al. (2000), sur les données TIGER par Martin et al. (1996) ou sur Terre pour le massif péridotitique de Ronda Chabrillat et al. (2000).

Rapport de bandes

Une seconde technique simple, qui nécessite peu de canaux pour être mise en oeuvre, consiste à analyser la forme des spectres à partir de rapports de bandes. Dans un premier temps, deux canaux peuvent être comparés l’un à l’autre. Une longueur d’onde non comprise dans l’absorption sert de référence par rapport à une seconde située dans le domaine de l’absorption. Dans l’exemple de la figure 22, le point gris sert de référence et les points bleus et rouges lui sont comparés. Pour le spectre d’olivine (en bleu), si la longueur d’onde correspondant à l’absorption maximum (point bleu) est comparée à la référence, alors la valeur trouvée est importante. Si au contraire la longueur d’onde comparée correspond à l’absorption maximum de l’orthopyroxène (point rouge) alors la valeur du rapport est plus faible. A l’inverse, sur le spectre de l’orthopyroxène (en rouge), la valeur de rapport la plus élevée correspond à la longueur d’onde de l’absorption maximum de l’enstatite. Dans le cas d’un spectre de mélange (en violet), les absorptions de chacun des minéraux sont présentes et les valeurs de rapport pour les deux constituants sont du même ordre.

Fig. 22 – Principe du rapport de bandes : la valeur pour une longueur d’onde hors de l’absorption est comparée à celle pour une longueur d’onde dans l’absorption. Spectres de forstérite, d’enstatite et de mélange (75% Ol – 25% Opx) acquis au RELAB (Brown University).

Grâce à cette technique, il est possible d’obtenir une estimation de la composition d’une roche. De plus, l’utilisation d’un rapport entre deux bandes permet de s’affranchir au premier ordre des variations photométriques (e.g. relief, éclairement), celles-ci étant considérées comme identiques pour les longueurs d’onde concernées (Chevrel et Pinet, 1990, 1992).

Il existe cependant des limitations à cette technique. En effet, il est nécessaire de connaître les positions respectives des absorptions pour chacun des minéraux potentiellement présents dans le spectre et être certain qu’il n’existe pas de superposition de ces absorptions. De plus, les valeurs déterminées à partir des rapports ne correspondent en aucun cas à une représentation physique de la composition et il n’est donc pas possible de quantifier les teneurs respectives de chacun des constituants (e.g. les valeurs de rapport de bandes similaires dans le cas du spectre de mélange alors que l’olivine est présente en proportion trois fois supérieure à l’orthopyroxène).

Pour améliorer l’efficacité de ce type de méthode, il est possible d’utiliser des indices, c’est-à-dire appliquer des coefficients aux valeurs des différentes longueurs d’onde pour faire ressortir l’information liée à une absorption caractérisée sur des données de laboratoire. Poulet et al. (2007) ont ainsi défini des indices pour déterminer la présence de constituants tels que l’olivine, les pyroxènes ou les oxydes ferriques à la surface de Mars.

Forme des absorptions

Pour caractériser la forme des absorptions, il est nécessaire de prendre en compte le continuum. Une fois normalisée, la forme de la bande d’absorption peut être étudiée en utilisant des paramètres simples que sont la profondeur, la largeur à mi-hauteur et la symétrie. La largeur à mi-hauteur (en µm ou nm) correspond à la distance entre les deux bords de l’absorption à hauteur de la demi-profondeur. La symétrie est quant à elle définie par les demi-surfaces de l’aire associée à l’absorption, c’est-à-dire les surfaces situées à gauche et à droite du maximum de l’absorption. L’intérêt de ce type de méthode est de réduire considérablement la dimension des données hyperspectrales.

Comparaison spectrale

Lorsque les roches ont une minéralogie plus complexe ou que les absorptions sont moins bien définies, les techniques détaillées précédemment restent limitées. Dans le cas où les données sont suffisamment résolues spectralement, il devient possible de réaliser une étude détaillée à partir de la forme générale d’un spectre, et de comparer celui-ci aux spectres existants dans les différentes librairies spectrales accessibles.

Deux méthodes peuvent être utilisées pour comparer deux spectres entre eux. La première consiste à mesurer la distance Euclidienne dans un espace à n dimensions, n étant le nombre de longueurs d’onde analysées, c’est-à-dire que pour une longueur d’onde donnée la distance entre le vecteur correspondant au spectre connu et le vecteur correspondant au spectre inconnu est mesurée. La distance spectrale (DS) calculée permet alors de définir un degré de ressemblance, c’est-à-dire que plus la distance spectrale est faible plus les spectres sont similaires. La seconde méthode, correspond à l’utilisation de l’angle spectral (SAM : Spectral Angle Mapper en anglais). Elle est similaire à la première si ce n’est qu’au lieu de mesurer la distance entre les deux vecteurs, c’est l’angle [latex]\alpha[/latex] entre ces derniers qui est mesuré. Launeau et al. (2002) et Combe et al. (2006) ont notamment mis en pratique ce genre de technique respectivement sur le massif péridotitique de Ronda (Espagne) et sur le massif ophiolitique de Sumail.

L’inconvénient de ces méthodes est la nécéssité d’établir une librairie spectrale de référence (de laboratoire ou de terrain) en tenant compte de l’ensemble des mélanges possibles pour chacune des conditions d’observations possibles. Le nombre de combinaison étant très important, il est difficile de mettre en oeuvre une telle procédure directe de façon automatique. Il existe toutefois une méthode améliorée, l’algorithme Tetracorder développé par Clark et al. (2003), qui permet de comparer les spectres uniquement dans la gamme de longueur d’onde affectée par une absorption.

1.2 Vers l’aspect quantitatif

Les techniques de traitement précédentes permettent d’obtenir une estimation des composants d’un spectre. Ils restent toutefois limités en termes de quantification des différents minéraux, surtout dans le cas de mélanges non-linéaires pour lesquels les librairies spectrales manquent de données. Il a donc été nécessaire de développer des techniques plus complexes pour accéder à ces informations.

Théorie du transfert radiatif

La première approche consiste à quantifier les processus de diffusion via l’équation du transfert radiatif. Cette équation permet de calculer l’interaction entre une onde lumineuse et un milieu constitué d’un empilement de couches homogènes en prenant en compte l’absorption et la diffusion. Hapke (1981a), complété par une série successive de publications (Hapke et Wells, 1981b, Hapke, 1984, 1986, 2002), a proposé une solution analytique au problème de transfert radiatif appliqué aux surfaces particulaires telles que les surfaces planétaires. Sa théorie prend en compte la diffusion anisotrope de la lumière par les particules, les effets d’ombrage et les phénomènes de réflexion à la surface des grains. Elle peut, de plus, être utilisée dans le cas de mélanges.

Bien que ce genre de théorie permette de calculer des spectres en réflectance, il faut tenir compte de l’ensemble des paramètres physiques impliqués pour inverser des données (Poulet et Erard, 2004). Il s’agit notamment de la taille des grains, de l’agencement des grains entre eux et de la géométrie de surface. Cette méthode n’est donc pas la plus aisée à mettre en oeuvre dans le cas de traitements d’images hyperspectrales de surfaces planétaires où ces paramètres sont mal contraints. Elle est de plus également soumises aux informations de laboratoire puisque les constantes optiques, fonction de la composition chimique, doivent être connus auparavant.

Mélanges linéaires

Il existe deux types possibles de mélanges pour une surface. Les mélanges linéaires correspondent d’un point de vue physique à la juxtaposition de surfaces homogènes au sein d’un même pixel. Le but de cette approche est donc de trouver, pour chaque pixel, la meilleure combinaison linéaire d’une série de spectres de minéraux purs enregistrés en laboratoire, et permettant de reproduire les données acquises. Cette méthode, initiée par Adams et al. (1986), permet de décrire les variations spatiales de l’information spectrale contenue dans une image. Cependant, en raison des problèmes de non-linéarité selon les types de mélanges, le résultat est utile pour une information au premier ordre mais il apparaît généralement des différences entre les proportions réelles des minéraux et celles calculées. Cette méthode a été très utilisée au cours des vingt dernières années (e.g. Pinet et al., 1993, Head et al., 1993, Tompkins et al., 1994, Martin et al., 1996) et a donné lieu à des améliorations successives concernant la détermination et l’optimisation du choix des pôles de mélange, avec l’introduction de la notion de mélanges itératifs (e.g., Roberts et al., 1998, Chabrillat et al., 2000, Pinet et al., 2000, Adams, 2006, Chevrel et al., 2008).

Un algorithme basé sur les mélanges linéaires de spectres a notamment été développé pour une utilisation en contexte martien (Combe et al., 2008). Cette approche utilise des spectres artificiels plats et à pentes pures dans la librairie de départ permettant de prendre en compte au premier ordre des effets liés à la diffusion atmosphérique, à la taille des grains, à l’ombrage et à la photométrie.

2 Le Modèle Gaussien Modifié (MGM)

Un dernier modèle de déconvolution couramment utilisé est le Modèle Gaussien Modifié, développé par Sunshine et al. (1990). L’intérêt de ce modèle est de prendre en compte directement les processus de transition électronique pour obtenir une quantification des minéraux constitutifs de la surface observée. Il s’applique donc parfaitement dans le cas des olivines et pyroxènes dont les spectres sont marqués par la présence du fer (Sunshine et Pieters, 1993, 1998).

Comme nous l’avons vu précédemment, chacune des techniques de traitement des données hyperspectrales possède ses avantages et ses inconvénients, leur intérêt respectif résidant dans l’objectif fixé par l’utilisateur. Les différences, en termes de détection des unités lithologiques, ont par exemple été mises en évidence par Gendrin (2004) qui a comparé les résultats obtenus à partir de certaines de ces techniques (SAM, Tetracorder, MGM et transformée en ondelettes) sur des données AVIRIS et ISM (Gendrin, 2004, Gendrin et al., 2006).

Le Modèle Gaussien Modifié est a priori tout à fait adapté pour caractériser à la fois les variations de composition modale de surfaces planétaires mafiques et ultra-mafiques et les compositions chimiques des minéraux associés, car il peut couvrir l’ensemble de la gamme des compositions envisagées. Les paragraphes suivants ont donc pour objectif de présenter le principe et les applications du Modèle Gaussien Modifié, tels que décrits dans la littérature.

2.1 Principe

Modèle Gaussien classique et Modifié

Le travail de Sunshine et al. (1990) repose sur le fait que les variations du spectre dans le visible et proche-infrarouge sont composées de bandes d’absorptions assimilables à des gaussiennes (Clark, 1984). Le modèle gaussien classique considère que l’énergie est la variable aléatoire pour toutes les absorptions.

Cependant, dans le cas des absorptions dues aux transitions électroniques, l’énergie est fonction de la distorsion des sites cristallographiques et de la longueur moyenne de la liaison ion-ligand (Burns, 1970). Dans les minéraux naturels, il existe des millions de mailles élémentaires et donc des millions de sites cristallographiques. Bien que ceux-ci présentent des caractères similaires, ils ne vont pas être parfaitement identiques à cause des défauts, des vides ou des substitutions qui peuvent se produire. Ces irrégularités, couplées aux vibrations thermiques, vont mener à une distribution statistique des longueurs moyennes de liaison ion-ligand (Sunshine et Pieters, 1993). La variable aléatoire devient alors cette longueur de liaison (r) et elle peut être reliée à l’énergie de l’absorption (e) par une loi de puissance :

[latex]e{\varpropto}r^{n}[/latex]

A partir de cette loi, Sunshine et al. (1990) suggère qu’une distribution gaussienne des longueurs moyennes de liaison peut être transformée en une distribution gaussienne ” modifiée ” des énergies d’absorption, suivant la relation :

[latex]m(x)=S.\exp\left\{ \frac{-(x^{n}-\nu^{n})^{2}}{2\sigma^{2}}\right\}[/latex]

avec la distribution gaussienne modifiée m(x) qui s’exprime en fonction de son intensité (S), de son centre ([latex]\nu[/latex]) et de sa largeur ([latex]\sigma[/latex]), x étant alors l’énergie. Modifier l’exposant n revient à modifier la symétrie de la distribution, c’est-à-dire la pente relative des ailes gauche et droite de la gaussienne.

Sunshine et al. (1990) ont pu déterminer la valeur du coefficient n empiriquement : dans le cas d’une gaussienne classique n=1, tandis que pour une gaussienne modifiée [latex]n\simeq-1[/latex]. Cette modification permet alors de modéliser une absorption électronique simple avec seulement une distribution gaussienne. La solution est donc plus proche de la réalité physique. Un exemple de différence entre les deux modèles est donné dans la figure 23. Le Modèle Gaussien Modifié va donc permettre d’analyser directement les absorptions dues aux transitions électroniques.

Fig. 23 – Différence entre la modélisation effectuée par le Modèle Gaussien classique (à gauche) et le Modèle Gaussien Modifié (à droite). Exemple sur une forstérite, d’après Sunshine et al. (1990)

Procédure mathématique

Après avoir présenté la particularité du Modèle Gaussien Modifié, ce paragraphe a pour but de décrire le fonctionnement mathématique général du modèle. Celui-ci est présenté plus en détail dans Sunshine et al. (1990) et Sunshine et Pieters (1993).

Les bandes d’absorption suivent la loi de Beer-Lambert, le spectre doit donc être converti en son logarithme népérien. Il est alors modélisé comme étant une superposition de gaussiennes modifiées et d’un continuum de forme polynomiale. La formulation mathématique s’exprime de la manière suivante :

[latex]\ln R(x)=\sum_{i=1}^{m}S_{i}\exp\left\{ \frac{-\left(x^{-1}-\nu_{i}^{-1}\right)^{2}}{2\sigma_{i}^{2}}\right\} +\ln\left(p_{0}+p_{1}.x+p_{2}.x^{2}\right)[/latex]

où R représente la réflectance, [latex]S_{i}[/latex], [latex]\nu_{i}[/latex] et [latex]\sigma_{i}[/latex] respectivement l’intensité, le centre et la largeur de chaque gaussienne, [latex]p_{0}[/latex], [latex]p_{1}[/latex], [latex]p_{2}[/latex] et [latex]p_{3}[/latex] les paramètres du polynôme et x la longueur d’onde. Chaque gaussienne modifiée est effectivement une fonction gaussienne de la longueur d’onde et la profondeur de bande [latex]S_{i}[/latex] doit toujours être négative. D’après cette équation, lors de l’inversion, les paramètres associés aux gaussiennes et au polynôme vont évoluer de façon conjuguée. La largeur à mi-hauteur, qui sera utilisé par la suite dès que nous parlerons de largeur de gaussienne, peut être obtenue à partir de [latex]\sigma[/latex] :

[latex]FWHM=2(2\ln2\sigma)^{\frac{1}{2}}[/latex]

Les paramètres [latex]S_{i}[/latex], [latex]\nu_{i}[/latex] et [latex]\sigma_{i}[/latex] sont optimisés par un algorithme d’inversion non-linéaire basé sur une méthode de moindres carrés appliquée de manière itérative. Sunshine et Pieters (1993) ont choisi d’utiliser la méthode d’inversion stochastique modifiée de Tarantola et Valette (1982). Cette méthode présente l’intérêt d’appliquer des incertitudes aux paramètres initiaux des gaussiennes. En optimisant au préalable la définition de ces incertitudes en fonction des absorptions présentes, le modèle devient à la fois plus réaliste et plus performant du point de vue du temps de calcul. La qualité de la modélisation peut être mesurée par l’erreur résiduelle (root mean square ou rms) :

[latex]rms=\left(\frac{\sum\varepsilon_{\lambda}}{n}\right)^{\frac{1}{2}}[/latex]

Le problème posé est alors de déterminer le nombre de gaussiennes et les paramètres de chacune d’entres elles pour que le modèle reproduise au mieux le spectre observé.