Projet Air

Précédent : Grands domaines d'application Remonter : Projet AIR, Traitement d'image et Suivant : Actions industrielles

• Fonctions implicites définies par systèmes de particules
• Comparaison de méthodes de mouvement adaptées aux images météorologiques
• Estimation du mouvement par hypothèse de conservation de la matière
• Evolution de la végétation
• Validation d'une chaîne de traitement interférométrique radar
• Détection d'effets de phasimétriques SAR associés à des changements d'état de surface
• Chemins géodésiques pour la mise en correspondance de structures météorologiques
• Etude des déformations de surface à partir de courbes d'iso-élévation
• Caractérisation de la température radiative en imagerie NOAA
• Définition d'un système de gestion pour l'environnement littoral

Résultats nouveaux

Fonctions implicites définies par systèmes de particules

Participants : Hussein Yahia , Jean-Paul Berroir , Gilles Mazars


Le modèle de fonction implicite étudié les anneés précédentes a été étendu et généralisé pour obtenir un système de contours actifs spécifique à la segmentation des structures déformables dans les séquences d'images satellite.

Formulation énergétique
Nous nous intéressons à des fonctions implicites définies par des points de contrôle appelés particules. Une énergie interne du type potentiel de Lennard-Jones : $${\cal{U}} = \displaystyle \frac{1}{2} \sum_{i \neq j} U_{ij}$$ avec $$U_{ij}=\displaystyle\frac{(r_{i}+r_{j})^2}{8d_{ij}^{4}}-\displaystyle\displaystyle\frac{1}{d_{ij}^{2}}$$ (où $r_i$ sont les rayons d'influence des particules et $d_{ij}$les distances entre particules $i$ et $j$)permet de modéliser un comportement visco-élastique du contour implicite. Les énergies externes sont au nombre de deux. La première $$E_{contour}=\sum_{\omega \in {P}}(\phi(\omega)-c)^2$$ utilise les données points ${P}$ fournies par un détecteur de contour pour forcer l'iso-contour à s'approcher de ces points.

Cette énergie est insuffisante car l'iso-contour peut être attiré par des points de contour correspondant à du bruit. La solution consiste alors à minimiser l'énergie sur un voisinage ``collier'' de l'iso-contour. Pour cela on introduit un masque:

 $$\chi_{\sigma} (P) = \exp - \displaystyle \frac{(\varphi (P)-c)^2}{\sigma^2}$$

dont les valeurs maximales sont atteintes pour les points $P$ situés sur une «bande» proche de l'iso-contour $\varphi ^{-1}(c)$. La largeur de cette bande dépend du paramètre $\sigma$. Nous considérons alors l'énergie de collier: $${E_{collier}} (...,x_i, y_i, r_i,...) = \displaystyle \frac{1}{K} \int \int_{I} D(P) \chi_{\sigma} (P)d(P)$$ où $D$ est la carte des distances. L'utilisation de $E_{collier}$permet d'éliminer les petites composantes connexes indésirables proches de l'iso-contour, liées au bruit de l'image.

L'énergie totale est minimisée par la méthode de la descente de gradient, qui est robuste et assez souple pour fournir une méthode d'ajout automatique de particules quand c'est nécessAire. Les résultats sont présentés sur les figures (,)
  Figure:  Segmentation d'une théière avec ajoût hiérarchique de particules: (a) initialisation, (b) résultat après minimisation et ajoût automatique de particules.

  Figure:  (a) initialisation des particules, résultat du seuillage; (b) (c) et (d): résultat de la minimisation sur trois images consécutives.

Comparaison de méthodes de mouvement adaptées aux images météorologiques

Participants : Dominique Béréziat , Isabelle Herlin , Laurent Younes


Nous étudions et comparons deux méthodes d'estimation de mouvement en les adaptant aux images météorologiques afin de résoudre certains problèmes inhérents à ce type de données.

La première méthode est une méthode de flot optique incorporant une contrainte de régularisation de type $L^1$ et une formulation variationelle, développée dans le projet par Isaac Cohen, permet de réaliser la minimisation. C'est une méthodologie classique qui nous sert d'étalon. La seconde méthode utilise une hypothèse de variation fixe des niveaux de gris et un modèle de mouvement supposé affine. La méthode de résolution est alors réduite à une estimation paramétrique.

  Figure: Première méthode appliquée à une image de vortex

  Figure: Seconde méthode

Ces deux méthodes ont des approches différentes et l'intérêt réside dans leur complémentarité. La première fournit un champ de vecteur dense estimé localement, la seconde fournit les paramètres d'un mouvement affine estimé globalement sur l'image. On constate que la première approche est bien adaptée pour l'estimation du mouvement de points isolés mais se révèle décevante pour estimer le mouvement global de structures importantes. Au contrAire, la seconde méthode donne une mauvaise estimation des points isolés mais permet une bonne estimation du mouvement global des structures.

Estimation du mouvement par hypothèse de conservation de la matière

Participants : Dominique Béréziat , Laurent Younes , Isabelle Herlin

De nombreuses méthodes d'estimation de mouvement utilisent l'hypothèse du flot optique. Or, la détection du mouvement de structures nuageuses ne vérifie que partiellement cette hypothèse. Nous essayons donc d'établir un critère de détection de mouvement mieux adapté aux structures étudiées. Nous proposons d'utiliser l'hypothèse suivante: ``entre deux instants suffisamment proches, une structure nuageuse garde le même volume dans un espace particulier''. En effet, les images infrarouges mesurent une température, laquelle est directement reliée à l'altitude du nuage observé. Un volume peut donc s'écrire comme un élément de surface de nuage multiplié par le niveau de gris (qui fournit la hauteur).

Soit $I(x,y,t)$ le niveau de gris du pixel de coordonnée $(x,y)$ au temps $t$. Un élément de volume au temps $t=0$ s'écrit:

$${d\cal V}_0 = I(x_0,y_0, 0)dS \mbox{ ,}$$ où $dS=dx_0dy_0$ est l'élément de surface.

Considérons la fonction $\varphi(x_0,y_0,t)$ qui représente au temps $t$ la position du pixel qui était à la position $(x_0,y_0)$ au temps $t=0$. Au temps $t$ l'élément de volume s'écrit donc:

$${d\cal V}_t = I(\varphi(x_0,y_0,t),t) d\varphi(S,t) \mbox{ ,}$$ soit encore: $${d\cal V}_t = I(\varphi(x_0,y_0,t),t) J(\varphi(x_0,y_0,t)dS \mbox{ .}$$ où $J$ est l'opérateur jacobien résultant du changement de variable.

La condition de conservation de la matière peut alors s'écrire:

$$\frac{d}{dt} {\cal V}_t = \frac{d}{dt} (I(\varphi(x,y,t),0)\varphi(x,y,t)) = 0$$

En différenciant, nous obtenons la condition:

 $$w \nabla I + I_t + div(w) I=0 \mbox{ ,}$$

où $w=\varphi_t$.

Une formulation variationelle de () est envisageable avec une contrainte sur le jacobien pour résoudre en $\varphi$ ou une contrainte sur $w$ pour résoudre en $w$.

$$E(w)=\int_\Omega ( w \nabla I + I_t + div(w)I)^2 dxdy + \int_\Omega (\vert\nabla u\vert^2+\vert\nabla v\vert^2)dxdy$$

Evolution de la végétation

Participants : Sonia Bouzidi , Jean-Paul Berroir , Isabelle Herlin

Cette étude fait depuis deux ans l'objet d'une convention de recherche avec le CNES, visant à mettre au point des outils de modélisation de l'évolution de la végétation pour le futur capteur VEGETATION, prochainement embarqué sur SPOT-4. Les objectifs sont l'estimation des rendements agricoles et l'étude du stress hydrique à l'échelle européenne. Une telle étude nécessite une estimation à haute fréquence temporelle de la réflectance des principaux types de culture, dans les bandes rouge et proche infra-rouge, pour en déduire l'indice de végétation normalisé (NDVI). Pour accéder à cette information avec précision, il est nécessAire de disposer de capteurs complémentaires:

1.

à haute résolution spatiale (SPOT), fournissant des mesures à l'échelle des parcelles; ces capteurs ont par contre une résolution temporelle insuffisante (26 jours);

2.

à haute résolution temporelle: NOAA-AVHRR (VEGETATION à l'avenir) fournit des données quotidiennes, mais avec une résolution spatiale de 1,1km, bien supérieure à la taille typique des exploitations en Europe.

La première année d'étude visait à calculer les réflectances quotidiennes des couverts végétaux étant données des relevés terrain à haute précision spatiale et une séquence (144 images) NOAA.

La deuxième année correspond à la génération de ces données terrain à partir d'une séquence SPOT couvrant la période d'activité de la végétation, soient 6 images de Mars à Juillet. En d'autres termes, il s'agit d'effectuer une segmentation spatio-temporelle des données SPOT et une reconnaissance supervisée de l'occupation du sol. Le site d'étude de Chartres est divisé en une zone d'apprentissage et une zone de test. Sur la zone d'apprentissage, des projections optimales sont calculées à chaque date suite à une ACP. Un processus de segmentation markovienne est alors défini dans l'espace des deux premières composantes principales, séparant au mieux les différents couverts. La levée des confusions est effectuée en fusionnant les résultats correspondant aux différentes images SPOT par des opérations logiques. On aboutit à un taux de bonne classification de 79% sur la zone de test sachant que les erreurs de reconnaissance sont très faibles sur les cultures et sont surtout importantes dans les zones urbaines (figure ). Par conséquent, la donnée de 6 images SPOT consécutives sur un site de pratiques agricoles similAires à celles du site d'étude permet de générer des données terrain suffisamment précises.

Les travaux de la dernière année du contrat (en cours) concernent un scénario différent, pour lequel l'utilisateur ne dispose pas de l'ensemble de la séquence SPOT, mais d'une seule image, et bien sûr de la séquence NOAA. L'information extraite de l'image SPOT est purement spatiale: segmentation du site en parcelles (figure ). L'image est alors dégradée à la résolution NOAA. Si l'on affecte à chaque parcelle un type de couvert, en utilisant les réflectances typiques des cultures du site (1ère année d'étude), on peut reconstituer un signal NOAA. Les types de couverts choisis minimisent la différence entre cette prédiction et le signal réel observé.

Ces travaux sont en relation avec les recherches menés par le projet ``Integration of VEGETATION and HRVIR data into yield estimation approach'' soutenu par le programme VEGETATION de l'International User Comitty et à ce titre nous avons bénéficié des données mises à disposition par le Centre Commun de Recherche d'Ispra (projet MARS).

  Figure:  Extrait d'une image SPOT-XS3, données terrain correspondantes et segmentation spatio-temporelle (79% de bonne classification).

  Figure:  Image de NDVI obtenue à partir des canaux 2 et 3 de SPOT, superposée avec la grille de pixels NOAA (gauche), et segmentation spatiale par croissance de région.

  
Validation d'une chaîne de traitement interférométrique radar

Participants : Laurent Da Silva, Etienne Huot, Isabelle Herlin


Mots-clés : interférométrie, radar

L'éruption du Volcan Grimsvötn en Islande a suscité de nombreuses observations et mesures in situ. Les données acquises par les satellites radar ERS1 et ERS2 avant et après les activités sismiques ont permis de générer des interférogrammes de cette région.

Grâce aux mesures réalisées sur le terrain et aux données fournies par l'agence spatiale européenne (ESA), la chaîne de déroulement de phase réalisée les années précédentes au sein du projet AIR a pu être éprouvée. Cette méthode est basée sur un déroulement local acceptant des contraintes globales. Il s'agit d'une technique itérative en trois étapes (déroulement local, amélioration par une segmentation markovienne, détection puis correction des zones résiduelles). Les résultats on été comparés à la fois aux mesures réelles et à ceux obtenus par une autre méthode de déroulement de phase, appliquée à l'ESA.

Il résulte de cette étude que le déroulement obtenu par notre méthode est acceptable tant que le seuil de cohérence est supérieur à 0,35. Les discontinuités dans la phase déroulée apparaissent à 89,9$\%$ dans les régions incohérentes. Dans ces régions, les erreurs peuvent atteindre quelques dizaines de mètres. Rappelons que ces discontinuités peuvent également être liées à des discontinuités réelles du terrain.

On peut donc conclure que la chaîne de déroulement que nous avons proposé offre des résultats dont la validité est fonction de la cohérence interférométrique. Elle doit donc être améliorée en prenant en compte un modèle de discontinuité.

Détection d'effets de phasimétriques SAR associés à des changements d'état de surface

Participants : Etienne Huot, Isabelle Herlin, Jean-Paul Rudant


Mots-clés : radar, multi-temporel, phasimétrie

L'analyse d'une séquence d'images multi-temporelles issues de l'imagerie radar montre que des effets de phase importants peuvent accompagner l'évolution des couverts végétaux à la suite d'épisodes pluvieux. Cette étude, menée en collaboration avec le laboratoire de géologie structurale de Paris VII, est consacrée à la détection de forts effets phasimétriques. Elle s'inscrit dans le cadre d'une action du GdR-PRC ISIS, qui a pour but d'amorcer, sur un nombre réduit de thèmes, une concertation approfondie entre les spécialistes des sciences de la terre utilisateurs de données radar et les traiteurs du signal et de l'image. Financée par le Programme National de Télédétection Spatiale, cette étude appartient aux travaux du groupe de travail lié aux «spécificités du multi-temporel en imagerie SAR».

La région test choisie pour cette étude correspond au site de Naizin, en Bretagne. Cette zone étant stable et les interférogrammes, générés par le CNES, corrigés des effets topographiques grâce à un modèle numérique de terrain, les perturbations phasimétriques ne peuvent être liées qu'à des modifications d'état de surface. Les parcelles cultivées ayant des comportements différents suivant leur nature, certaines structures géométriques correspondant aux champs apparaissent nettement (cf. flèche sur la figure ).

Plutôt que de chercher à localiser chacune de ces parcelles à fort effet phasimétrique, on segmente la région qui les entoure. Cette région de fond se caractérise par trois propriétés complémentAires:

1.

variation homogène de la phase hors de ces régions à fort effet phasimétrique,

2.

cohérence localement homogène,

3.

régularité de la région résultant de la segmentation.

Ces propriétés sont traduites en une formulation markovienne. Les paramètres du modèle sont estimés sur une fenêtre de taille fixe grâce à une pré-segmentation utilisant une détection de contour sur l'interférogramme.

  Figure:   Différentes étapes de la segmentation : (a) partie d'un interférogramme déroulé, (b) cohérence correspondante, (c) résultat de la présegmentation et résultat final projeté sur l'interférogramme non déroulé.

Les images d'amplitude radar n'ont pas été utilisées lors de la segmentation, elles devraient intervenir dans une étape de classification.

Chemins géodésiques pour la mise en correspondance de structures météorologiques

Participants : Isaac Cohen , Isabelle Herlin

  Le suivi de structures nuageuses telles que les vortex dans les images météorologiques satellitAires est souvent limité par le faible échantillonnage temporel (une image toutes les demi-heures). La caractérisation de l'évolution nécessite donc l'utilisation d'un modèle permettant d'apparier des courbes ayant subies de larges déformations et de topologie variable.

Notre approche est basée sur le calcul d'un ensemble de chemins reliant les courbes à mettre en correspondance. Ces courbes, définies par deux ensembles source ${\cal S}$ et destination ${\cal D}$, sont représentées à l'aide de fonctions implicites. Les trajectoires minimisent une fonction de coût mesurant la similarité locale: elles représentent des géodésiques de la surface représentant cette fonction. Ces géodésiques sont calculées à l'aide de l'équation:

 $$\varphi_{t} = \sqrt{a\varphi_{x}^{2} + b\varphi_{y}^{2}-c\varphi_{x}\varphi_{y}}$$


où $\varphi$ est la représentation implicite de la courbe et $a,b$ et $c$ sont des scalAires fonction de la surface sur laquelle sont calculées les géodésiques. Cette surface est définie en coordonnées cartésiennes par: $$Z=(x,y,z(x,y)) = (x,y,min(\vert\varphi_{0}\vert,\vert\psi_{0}\vert))$$ où $\varphi_{0}$ et $\psi_{0}$ sont les représentations implicites des courbes à apparier. L'équation  permet de calculer les cartes de distances géodésiques des courbes à apparier sur la surface $Z$.

La caractérisation des trajectoires consiste à définir pour chaque point $X_S$ appartenant à la source, la courbe de coût minimal le reliant à un point de la destination. La fonction de coût est:

 $$f(x,y) = \varphi(x,y) + \psi(x,y)$$

Les trajectoires sont définies par le point de départ: $p(0)=X_{S}$, donné et par: $$\frac{\partial p}{\partial s} = - \nabla({\cal D}_{S} + {\cal D}_{D})$$ où ${\cal D}_{S}$ et ${\cal D}_{D}$ sont les solutions de l'équation ().

La figure  représente deux images météosat consécutives d'une tempête tropicale. Nous remarquons que la structure correspondant au nuage a subi une déformation de forte amplitude et un changement de topologie. Les chemins d'appariement entre ces deux structures, obtenus par la méthode présentée, sont illustrés par la figure .
  Figure:  Segmentation des régions à apparier à l'aide d'un modèle de contours actifs.

  Figure:   Tracés des chemins ou trajectoires d'appariement entre les deux structures nuageuses.

Etude des déformations de surface à partir de courbes d'iso-élévation

Participants : Christophe Cubier , Isaac Cohen

L'interférométrie d'images radar à synthèse d'ouverture permet de détecter des déformations au sol de faible amplitude. L'étude de ces déformations se fait par rapport à un modèle numérique de terrain et plus particulièrement par rapport à un ensemble de courbes de niveaux. Pour cela, il est nécessAire d'être en mesure de mettre en correspondance des courbes assez complexes et de topologie variable.

Dans le cadre de cette étude nous proposons une méthode de mise en correspondance de courbes basée sur le calcul de chemins géodesiques. Ces chemins de mise en correspondance sont des géodésiques d'une surface calculée à partir des courbes à mettre en correspondance. Cette partie est réalisée en utilisant la méthode décrite dans la section . Mais cette étude spécifique prend en compte la géométrie des courbe lors de la définition de la surface sur laquelle sont calculées les géodésiques: dans le cas de déformations de faible amplitude, l'utilisation de la courbure de la courbe est pertinente. Pour cela, nous avons redéfini la surface sur laquelle sont calculées les géodésiques. Nous utilisons la surface cartésienne suivante:

$$Z=\left(x,y,min\left( \left(1-\frac{(H_S-H_D)^2}{1+(H_S-H_D)^2\... ...-\frac{(H_S-H_D)^2}{1+(H_S-H_D)^2\psi ^2}\right)\psi\right)\right)$$

où $H$ est la courbure moyenne des contours source ou destination definies respectivement par $\varphi ^{-1}(0)$ et $\psi ^{-1}(0)$.

Nous avons choisi, pour illustrer cette méthode, les données interférométriques issues du glacier Vatnajökull présenté section . La figure  représente une courbe de niveau avant et après la déformation de surface. La mise en correspondance de ces courbes est illustrée par la figure .

  Figure:  Deux ensembles de courbes de niveaux à mettre en correspondance.

  Figure:  Représentation des chemins de mise en correspondance entre les deux ensembles.

Caractérisation de la température radiative en imagerie NOAA

Participants : Yannick Marrec , Isabelle Herlin , Isaac Cohen

La différence entre la température de surface d'un couvert végétal et la température de l'Air est signifiante de l'état de ce couvert. L'écart reste faible dans le cas d'un couvert en bonne santé et augmente dans un état de stress hydrique.

Le but de ce travail est donc d'étudier et de caractériser la variabilité de la température de surface à partir de séquences d'images quotidiennes NOAA en fonction de différents paramètres tels que le type de végétation, la topographie, la vapeur d'eau,...

Dans un premier temps, une étude expérimentale est menée sur un site de faible relief afin de caractériser le lien entre les variations de température et le type de végétation présent au sol. Les variations de températures sont caractérisées à l'aide de la composante normale du flot optique permettant ainsi de fAire apparaître les régions pour lesquelles des variations en temps et en espace sont observées.

La figure met en évidence les corrélations existants entre le couvert végétal et la température de surface.
  Figure:  Concordances entre la température de surface calibrée et le NDVI.

Les variations de cette température sont illustrées dans la figure  et sont étroitement liées aux variations du NDVI et donc à l'évolution phénologique du couvert végétal.
  Figure:  Correspondances entre l'évolution de la température et les valeurs de NDVI.

Ce travail s'effectue dans le cadre de la proposition européenne INCO-PED IWRMS en collaboration avec l'Afrique du Sud.

Définition d'un système de gestion pour l'environnement littoral

Participants : Isabelle Herlin , Eric Simon

Le but de cette étude, qui débute fin 1997, est de répondre aux besoins des scientifiques et décideurs qui sont quotidiennement obligés de rechercher, traiter et visualiser des données stockées sous différents formats dans différents lieux.

Nous avons donc, pour l'application de gestion du littoral méditerranéen, défini une architecture système pour répondre aux besoins identifiés: accès par une interface WEB à des données SIG, à des images satellites ou à des cartes, à des mesures in-situ, à des documents et à des outils d'interprétation, de simulation et de visualisation.

Cette architecture utilise les fonctionnalités de deux systèmes différents: un moteur de recherche distribué, et un médiateur permettant des requêtes dans des bases de données distribuées.

La contribution spécifique du projet AIR dans ce système concerne l'indexation et la description des données satellites et des modèles sous-jacents pour permettre l'accès par requêtes. Ce travail s'effectue dans le cadre de la proposition européenne Telematics THETIS.

Précédent : Grands domaines d'application Remonter : Projet AIR, Traitement d'image et Suivant : Actions industrielles