Diffusion de données LANDSAT sur la France par THEIA

=>

Nous avons le plaisir de vous annoncer que les données LANDSAT de niveau 2A, produites par le pôle THEIA sont disponibles à l'adresse suivante :

http://theia.cnes.fr/

Cet article vous décrit les données et le fonctionnement du serveur de distribution. Un autre article vous décrit la méthode pour télécharger des séries entières de données.

Les données de Theia sont traitées jusqu'au Niveau 2A : elles sont donc fournies en réflectance de surface après correction atmosphérique, et sont accompagnées de masques de nuages. La façon dont nous les avons produites est décrite ici pour Landsat 8, et ici pour LANDSAT 5 et 7. Le serveur, mis en place par mon collègue du CNES Jérôme Gaspéri, aidé par Rémi Mourembles de CAP Gemini, est doté d'une interface épurée et résolument moderne, avec une seule petite fenêtre pour saisir les requêtes, qui peuvent être rédigées sous forme de phrases. L'outil pratique donc une analyse sémantique (waouh!) de vos requêtes. Il fonctionne aussi sur votre tablette ou téléphone, et bientôt sur votre réfrigérateur.

 

Exemples de requêtes :

1) Recherche de lieu et de date :
images LANDSAT 8 à Paris entre janvier et juin 2013
images LANDSAT 7 à Toulouse acquises en mai 2013
2) Recherche de caractéristiques :
zone herbacée en 2013
images avec forêts en octobre 2009
images sans forêts en octobre 2013
3) Combinaison de caractéristiques :
Images avec zones cultivées et forêts près de Paris entre mars et juin 2013
zone cultivées à Bordeaux en 2013.
4) Style télégraphique
Landsat8 Juillet 2013
LANDSAT5 2011
Arcachon 2009

 

Pour le choix de la zone géographique, vous pouvez aussi zoomer sur la carte pour définir votre zone d'intérêt.

 

Enfin pour télécharger les données, il faut, dans un premier temps créer un compte (ça se passe sur l'icône orange, en haut à gauche), puis vous identifier. Cette inscription nous sera utile pour établir des statistiques d'utilisation. Il ne vous reste plus qu'à télécharger la donnée en cliquant sur l'icône avec le petit nuage. Chaque image est aussi téléchargeable directement à partir de son url définie logiquement avec le nom du produit. Je vais essayer de faire  rapidement un script de téléchargement plus automatique, mais déjà, vous pouvez utiliser le plugin DownThemAll de Firefox, en lui spécifiant, après vous être authentifiés, de télécharger tous les liens du type "*download" présents dans la page (mode d'emploi ici)

 

La publication de ces données est le résultat de plusieurs années de travail, au CESBIO et au CNES, et avec nos sous traitants, mais c'est aussi la première version. Merci de nous faire parvenir vos constats et remarques, surtout si vos avis sont favorables, mais pas seulement, les critiques nous seront très utiles pour améliorer le service, avant l'arrivée de Sentinel 2 dont le lancement est prévu dans un an !

 

Enfin, nous souhaitons vraiment remercier nos collègues de la NASA et de l'USGS pour LANDSAT 7 et 8, et de l'ESA pour LANDSAT 5 qui mettent ces données à disposition sans aucune restriction, ce qui nous permet d'en faire de même au pôle THEIA. Les données sont donc disponibles pour tout utilisateur et pour toute utilisation. Utilisez les donc autant que vous voudrez, mais n'oubliez pas de nous faire part de ce que vous en faites, c'est très important pour nos futures demandes budgétaires.

Completion of the processing of LANDSAT 8 Level 2A products taken above France in 2013.

=>

That was fast ! The processing of all the LANDSAT 8 images taken above France in 2013 took less than 15 days. The first LANDSAT 8 images were taken in April 2013. The MUSCATE team processed the data for the THEIA land data center, using CNES computing center.

 

A few more days will be necessary to upload the data on the THEIA website and to check that the data are correct. Finally, the longest part in the processing is the downloading of the input Level 1T products from USGS earthexplorer website (equivalent to the Level 1C in THEIA's nomenclature).

 

The Level 2A data quality is quite good, as may be seen on the browse products on the right, as shown by images on the right, which come from the times series obtained on the tile of Paris. As usual on this blog, the clouds are circled in green, the shadows in black, the snow in pink and the water in blue. A few clouds are sometimes missed by our multi-temporal method, when the repetitivity of cloud free acquisition is too low, as in the image on the right which was acquired during a cloudy spring. The following images in the time series are not affected by this kind of defect.

 

This paper aims at describing the main steps of the processing.

 

Starting point : Level 1T

We download the input data from the earthexplorer website, using an enhanced version of the script described here. These products are ortho-rectified by USGS, using a global data base of ground control points.

 

The location requirement for LANDSAT 8 is 50 m, which seems to be met by the L1T products. We found location errors around 1.5 pixels near Toulouse, but most regions seem to have better performances. USGS confirmed a 38m bias Southward near Toulouse and will try to correct them.  Our processing does not correct for these small errors, and the next version of the USGS LANDSAT 8 processing only wil lcorrect for this bias.

 

Regarding LANDSAT 8 absolute radiometric calibration, we use the coefficient values recommended by LANDSAT 8 and provided with the L1T products.

 

Resampling to Lambert'93 projection


Level 1T data are provided with the UTM projection. This projections uses three different zones over France, for which the registration of data is not direct. We decided to resample the data on a Lambert'93 projection, which is the official French projection.

Tiling of products

We chose to tile the data in 110*110 km tiles spaced with a 100 km interval, as it will be done for Sentinel-2. The (1,1) tile is in the SouthWest corner of France. The tile of Toulouse is the 5th to the West, and the 2nd to the North. It is named D0005H0002 (D for "droite", H for "Haut")

 

For Corsica, a different tiling made of 2 tiles was defined.

 

For each tile, we provide the whole set of dates for which a LANDSAT 8 image intersects the tile. A few date may be missing, for several reasons, in general related to the cloud cover :

  • The image was not acquired by LANDSAT 8 (when a 100% cloud cover is forecast, the image is not acquired).
  • The image was acquired but not processed to level 1T by LANDSAT8, because the cloud cover prevented from using a sufficient number of ground control points
  • The Level 2A processing rejects images with more than 90% of cloud cover.

 

Level 2A processing (atmospheric correction and cloud screening)

First of all, we would like to outline that our processor does not process the themal bands of LANDSAT 8.

For the visible, near and short wave infrared bands, we use the same method as for SPOT4(Take5). It involves also the MACCS processor, developed and maintained by Mireille Huc at CESBIO. It is based on multi-temporal methods for cloud screening, cloud shadow detection, water detection as well as for the estimation of the aerosol optical thickness.

 

However, thanks to LANDSAT 8 spectral bands, our processing was enriched compared to SPOT4 (Take5) : LANDSAT8's 1.38µm band enables an enhanced detection of high and thin clouds. And thanks to the blue band, we have an additional criterion to detect the aerosols, thanks to the quasi constant relationship between the surface reflectances in the blue and in the red above vegetation. The precision gain due to this criterion compensates for the precision loss due the lower repetitivity of  LANDSAT8 images.

Level 2A images from Paris's tile, from 3 different LANDSAT 8 tracks (From left to right, tracks 200, 199, 198). The viewing angle differs as the image is from the west on the left image, at nadir in the center and from the east for the right image.

 

To enhance the cloud screening accuracy, we decided to use the data from adjacent satellite tracks within the same time series. These data are not acquired under exactly the same angle (+/- 7 degrees), which is the assumed by the multi-temporal method, but the difference is small enough to allow a large accuracy gain due to the enhanced repetitivity. However, because of this approximation, a few artefacts may be observed.

 

For a greater enhancement, we might also use LANDSAT 7 and LANDSAT 8 data in the same time series, but we will implement that later on...

 

Data Format

For LANDSAT 8, we used the same format as for SPOT4 (Take5), excepted a few details, that I will describe soon...

 

 

 

Le traitement des données acquises par LANDSAT8 sur la France en 2013 est terminé.

=>

Ça torche* ! J'annonçais le démarrage du traitement dans un précédent article: il aura fallu moins de 15 jours pour traiter l'ensemble des données acquises par LANDSAT 8 en 2013 sur la France métropolitaine (y compris la corse). La période traitée s'étend d'avril 2013 à fin 2013. Les traitements ont été effectués par l'équipe MUSCATE, pour le compte du pôle THEIA, en utilisant les moyens du centre informatique du CNES.

 

Il suffira de quelques jours de plus pour mettre les données en ligne sur le site du pôle et vérifier que les données peuvent être distribuées. Finalement, le plus long dans ce traitement est le téléchargement des produits de Niveau 1T (l'équivalent du Niveau 1C dans la nomenclature de THEIA) sur le site earthexplorer de l'USGS.

 

Les données obtenues sont de bonne qualité, comme le montrent les exemples de quicklooks de niveau 2A, à droite, tirés de la série temporelle obtenue sur la tuile de la région parisienne. Comme d'habitude, les nuages détectés sont entourés en vert, les ombres en noir, la neige en rose et l'eau en bleu. Quelques nuages sont parfois oubliés par notre méthode multi-temporelle de détection nuageuse, lorsque la répétitivité des données est insuffisante, comme, sur l'image ci contre, lors de ce printemps 2013 extrêmement nuageux. Les images suivantes de la série ne présentent pas ce genre d'erreurs.

 

Cet article décrit brièvement les principales étapes du traitement qui ont permis d'obtenir ce résultat.

 

Point de départ : Niveau 1T

Nous téléchargeons les données depuis le site earthexplorer, en utilisant une version améliorée du script de téléchargement décrit ici. Ces données sont orthorectifiées par l'USGS, à partir d'une base de données mondiale de points d'appuis.

 

Les spécifications de localisation pour les données LANDSAT de l'USGS sont de 50 mètres. Nous avons constaté, sur la région de Toulouse, des erreurs de localisation égales à 1,5 pixels, mais d'autres régions sont mieux loties. L'USGS a confirmé les erreurs que nous observons du côté de Toulouse (38m vers le sud) et promet de les corriger. Notre traitement ne corrige pas ces erreurs qui restent modérées, et pour obtenir une performance améliorée, il faudra attendre le retraitement de ces données par l'USGS, ou rattraper le décalage éventuel par vous mêmes.

 

En ce qui concerne l'étalonnage absolu du satellite LANDSAT 8, nous utilisons celui recommandé par la NASA (et fourni avec les produits LANDSAT 8).

 

Reprojection en Lambert 93


Les données de l'USGS sont fournies en projection UTM. Cette projection utilise trois fuseaux différents au dessus de la France, qui se divise donc en 3 zones différentes, l'Ouest de la France, le Centre et l'Est. Comme les données de deux fuseaux différents ne se superposent pas directement, nous avons donc décidé de reprojeter les données en Lambert 93, qui est la projection officielle pour la France.

Découpage des produits en tuiles

Nous avons pris le parti de suivre la même logique que celle utilisée par Sentinel-2, et de découper les données en tuiles de 110*110 km décalées de 100 km les unes par rapport aux autres. La tuile 1x1 se trouve au sud-ouest de la France, lorsqu'on va vers l'est (vers la droite), on incrémente la première coordonnée de D0001 à D0010 (D pour Droite), lorsqu'on va vers le Nord (vers le haut), on incrémente la seconde coordonnée, de H0001 à H0010 (H pour Haut). La tuile de Toulouse s'appelle donc D0005H0002  Le découpage en tuiles est visible sur l'image ci-jointe.

 

Pour la Corse, nous avons défini une autre grille de tuiles, composée de deux tuiles.

 

Pour chacune des tuiles, nous fournissons l'ensemble des dates pour lesquelles une image LANDSAT8 a une intersection non nulle avec la tuile. Quelques dates peuvent manquer, pour plusieurs raisons, liées en général à la couverture nuageuse :

  • l'image n'a pas été acquise par LANDSAT 8 (quand les prévisions météo indiquent un temps très couvert, les images ne sont pas acquises).
  • l'image a été acquise mais s'est avérée trop nuageuse pour être traitée au niveau 1T par l'USGS
  • l'image est trop nuageuse pour être traitée par la chaîne de Niveau 2A

 

Traitement de Niveau 2A (correction atmosphériques et détection des nuages)

Il est important de noter que notre chaîne ne traite pas les bandes thermiques pour l'instant. Une correction est à l'étude, mais celle-ci ne sera pas opérationnelle avant un ou deux ans.

 

La méthode utilisée pour les bandes visible, proche et moyen infra rouge est quasiment la même que pour SPOT4(Take5). Le traitement a été effectué avec la même chaîne, le prototype de MACCS, développé et maintenu au CESBIO par Mireille Huc. Notre méthode de base est une méthode multi-temporelle à la fois pour la détection des nuages, des ombres de nuages, de l'eau et pour l'estimation de l'épaisseur optique des aérosols.

 

Cependant, grâce à la richesse spectrale de LANDSAT, nous avons pu enrichir nos méthodes par rapport à la version utilisée pour SPOT4 (Take5) :  la bande 1.38µm de LANDSAT8 permet de détecter les nuages hauts. Et grâce à la bande bleue, nous pouvons utiliser un critère complémentaire pour détecter les aérosols, grâce à la relation quasi constante observée entre les réflectances des bandes bleues et rouges au dessus de la végétation. Le gain de précision dû à la présence de cette bande permet de compenser la perte de précision de la méthode multi-temporelle due à la faible répétitivité de LANDSAT8.

Images de la tuile obtenue sur la région parisienne, en provenance de trois traces différentes de LANDSAT (de gauche à droite, les traces 200, 199 et 198). Les angles de visée sont légèrement différents sur chacune des traces (visée depuis l'ouest sur l'image de gauche, visée au nadir au centre, visée depuis l'est sur l'image de droite).

 

Pour augmenter la précision de la détection des nuages, nous avons décidé d'utiliser les données issues de traces adjacentes de LANDSAT8 dans les séries temporelles de niveau 2A. Ces données ne sont pas acquises exactement sous le même angle (+/- 7 degrés), mais la différence d'angle est suffisamment petite pour qu'il y ait un vrai gain de précision sur les zones d'intersection entre traces. En raison de cette approximation, quelques artefacts peuvent être observés.

 

Pour gagner davantage de répétitivité, nous pourrions aussi entrelacer des données LANDSAT 7 et LANDSAT 8 dans des séries temporelles communes, mais ceci ne sera implémenté que dans une prochaine version.

 

Format des données

Le format des produits de Niveau 2A de LANDSAT 8 est le même que celui des données SPOT4 (Take5), à quelques détails près que je fournirai dans une nouvelle page, dès que j'ai un moment...

 

* Ca va vite, oui, c'est comme ça qu'on parle chez les anciens jeunes du CNES...

 

 

The V2.0 of SPOT4 (Take5) data set is available.

Voilà ! The new version (V2.0) of SPOT4-Take5 data set is available, for the 45 sites. I would like to thank the development and processing teams of MUSCATE center in CNES, who work for THEIA, the image quality teams at CNES (SI/QI and SI/MO), and of course Mireille Huc at CESBIO, for the production of this new version, which finally required a lot of work.

The product version number is not included in the filenames, but you can recognise a V2.0 product by looking into the xml metadata file :

<METADATA>
  <HEADER>
    <VERSION>2.0</VERSION>

 

This reprocessing brings the following new features :

    Quicklooks are now provided with the images. The clouds are circled in green, the shadows in black, water in blue and snow in pink.

  • We provide quicklooks on which you can see the cloud and shadows masks
  • We enhanced the quality of the ortho-rectification :
    • By changing the référence ortho-image for the sites in France (GEOSUD, processing done by the french institute for geography IGN)
    • by replacing the LANDSAT 5 otho-images by LANDSAT 8 images for most other sites outside France. LANDSAT 8 geometric performances are enhanced compared to  LANDSAT 5.
    • however, for a few sites (Borneo, Gabon, Congo (1,CNES), CCRS, Cameroun), no clear LANDSAT 8 images was available yet and we had to keep the LANDSAT 5 reference.
      • It's not too bad for Congo, CCRS et Cameroon, as LANDSAT 5 references where quite good, for Gabon, we used a reference made with the cloud free image obtained with SPOT4-Take5, and finally, we just have Borneo site for which the level 1C obtained are quite bad with large registration errors (I am sorry Jukka)
      • A large enhancement of the performances has been observed for Sumatra, Gabon and Congo (2,ESA), for which the first version was quite bad.
  • SPOT4 radiometric calibration updated
    • A the end of SPOT4's life, my CNES colleagues updated its absolute calibration. Spot calibration is obtained using desert sites, using another satellite as reference. Up to now, it was POLDER, but now it is MERIS/ENVISAT. Moreover, the calibration coefficients we used in the first version had been extrapolated from older measurements, while now recent measurements have been used. The differences are not too big, except for the near infra-red band which varied by 4%..
  • The level 2A have been reprocessed with a new version of the aerosol model, with larger aerosols. The previous model had been tuned for sites in France, but we found that  the larger particles fitted better the in situ data on all the sites.
  • For users of mountain sites, we added a few flags about the correction of terrain effects. If the slope is in the shade, or nearly in the shade, the correction we have to do is infinite ! We limited the value of the correction and flagged the pixels for which we had to limit it in the .DIV files.
  • And at last, the Maricopa site was finally processed. This site was acquired under two angles, one from the East, one from the West. It has therefore been observed twice every 5 days under different viewing angles. Such a case was not anticipated in our prototype, and we had to correct it. The site has been divided in 2 sites Maricopa_J1 for observations from the West, and Maricopa_J5 for observations from the East. This site, which benefits from New Mexico blue skies, is a very interesting one for remote sensing geeks, as it combines multi angular and multi-temporal observations at constant angles !

Les données SPOT4 (Take5) V2.0 sont disponibles.

=>

Voilà, les données SPOT4 (Take 5) V2.0 sont en ligne pour les 45 sites. Un grand merci à l'équipe de développement et d'exploitation de MUSCATE, pour le compte de THEIA et du CNES, aux équipes de qualité image au CNES (SI/QI et SI/MO) au CNES, et à Mireille Huc, au CESBIO, pour l'obtention de cette nouvelle version des données, qui a finalement demandé beaucoup de travail.

Nous n'avons pas inclus la version du produit dans le nom du fichier dans ce système prototype, mais vous pouvez reconnaître un fichier de la version 2.0 de la manière suivante :
Il suffit d'aller chercher le tag <VERSION> dans le fichier .xml qui porte le même nom que le produit:

<METADATA>
  <HEADER>
    <VERSION>2.0</VERSION>

 

Ce retraitement apporte les nouveautés suivantes :

    Des "quicklooks" sont maintenant fournis avec les images. Les nuages y sont entourés en vert et les ombres en noir, la neige en rose et l'eau en bleu

  • fourniture de quicklooks avec affichage du masque de nuages
  • Amélioration de l'ortho-rectification des données :
    • changement de référence géométrique pour les sites en France (GEOSUD, traitement IGN)
    • remplacement des ortho-images de référence pour l'ortho-rectification issues de LANDSAT 5 pour la plupart des autres sites dans le monde. LANDSAT 8 a des performances géométriques bien meilleures que celles de LANDSAT 5.
    • Cependant, pour quelques sites (Borneo, Gabon, Congo (1,CNES), CCRS, Cameroun), nous n'avions pas trouvé d'image LANDSAT 8 sans nuage, nous avons gardé l'image LANDSAT 5 comme référence, sauf pour le Gabon.
      • Pour Congo, CCRS et Cameroon, les références LANDSAT 5 étaient assez bonnes, pour le Gabon, nous avons utilisé une bonne référence issue de SPOT, il n'y a que Borneo pour lequel l'ortho-rectification obtenue est restée assez mauvaise, avec beaucoup d'images éliminées et de grosses erreurs de superposition
      • Une très forte amélioration de la qualité est notamment observée sur les sites Sumatra, Gabon et Congo(2 (ESA))
  • Mise à jour de l'étalonnage radiométrique de SPOT4
    • A la fin de la vie de SPOT4, les collègues du CNES ont mis à jour son étalonnage. L'étalonnage de SPOT est calculé sur des sites désertiques, par rapport à un capteur de référence. Jusqu'ici, c'était POLDER qui servait de référence aux étalonnage de SPOT, maintenant c'est MERIS/ENVISAT. De plus, la variation temporelle des coefficients d'étalonnage était extrapolée à partir de mesures plus anciennes, alors que cette fois des mesures récentes ont été intégrées. Les différences ne sont pas très importantes sauf pour la bande proche infra-rouge, où elles atteignent 4%.
  • Nouvelle version du traitement N2A avec un nouveau modèle d'aérosols, composé de particules plus grosses. Le modèle précédent avait été optimisé pour des sites en France. Nous avons trouvé que le modèle utilisé maintenant fournissait de meilleurs résultats en général.
  • Pour les sites montagneux surtout, ajout de flags sur le traitement de la correction des effets du relief. Si la pente est à l'ombre ou presque à l'ombre, la correction à appliquer s'approche de l'infini ! Nous avons donc limité la correction à apporter, et les pixels pour lesquels cette limitation s'exerce vous sont signalés dans les fichiers .DIV
  • Le site Maricopa a enfin été traité. Ce site a été acquis sous deux angles de prise de vue, l'un depuis l'est, l'autre depuis l'Ouest. Il a donc été vu deux fois tous les 5 jours, mais sous des angles différents. Ce cas n'était pas prévu dans notre système, mais ce problème est maintenant corrigé. Le site a été divisé en deux sites, Maricopa_J1 pour les observations depuis l'ouest, et Maricopa_J5 pour les observations depuis l'est. Ce site, situé au nouveau Mexique est très intéressant pour les geeks de la télédétection puisqu'il combine observations multi-temporelles à angles constants et observations sous des angles différents. Sans compter qu'il y fait toujours très beau temps !

The level 3A products

Among the products prepared to be processed by the THEIA land data center, the level 3A product was not yet described in this blog. The level 3A products provide a monthly synthesis of the level 2A. These products should be very useful for the following reasons :

  • The level 3A, produced once a month, uses up to six times less volume than the level 2A products acquires during a month.
  • The level 3A provides a regular time sampling of the reflectances variation, while the level 2A sampling is dependent on the cloud cover
  • Several processing methods and applications are hindered by the data gaps due to cloud cover. The level 3A product aims at minimizing the residual gaps.

 

Thanks to SPOT4 (Take5) data set, we were able to try and test several methods to produce level 3A products on varous types of landscapes and climates. This work, suprvised by Mireille Huc and myself, is performed by Mohamed Kadiri, at CESBIO, and is funded by the CNES budget of THEIA Land Data Center. Our method consists in computing, foe each pixel, a weighted average of the surface reflectances of the cloud free observations, obtained within a N day distance frome the central date TO of the level 3A product. For instance, the example below was obtained with N=21, for the 15 th of each month. As a result, the level 2A used in the average for the level 3A product dated on March the 15th, were acquired from Feruary the 24th to April the 4th.

 

 

The weighted average gives more weight to

  • the cloud free images
  • the pixels which are far from clouds
  • the images with a low aerosol content
  • the images acquired near the level 3A product date.

Les values of the weight and of the duration N, have a large influence on the product data quality. To tune their values, we set up three quality criteria :

  • The percentage of residual data gaps for which all the observations were cloudy
  • The difference of the level 3A reflectances with the values of a selected level 2A product acquired near the central date T0.
  • A measurement of the artefacts standard deviation. The artefacts appear near the borders of data gaps that affect one of the dates used in the level 3A synthesis.

 

For instance, here are the performances obtained on the Versailles site, which was heavily clouded in the spring of 2013. For this site, one can note, that the residual gap percentage is very low despite the bad weather, confirming that Sentinel-2 should be able to provide cloud free Level 3A products each month. For this site, the optimal duration of the synthesis is somewhere between 2* 21 and 2*28 days.

 

Performances obtained for Versailles SPOT4(Take5) site, for several values of the half-period N. In red, the residual percentage of data gaps (scale on the right), in yellow and green, the maximum value of the difference of the level 3A to the central level 2A, for resp the best 70% and 95% of pixels. In blue, the residual error standard deviation.

 

 

 

For Sentinel-2, the level 3A will have to include a correction for directional effects, in order to use in the same level 3A product, the data acquired from different satellite tracks, from different viewing angles. Finally, as an option, we might include a gap-filling method to fill the residual gaps.

In short, we still have work to do. A comparison with the classical NDVI Maximum Value Composite is provided in this post.

Les produits de Niveau 3A

=>

Parmi les produits que nous préparons pour les traitements de données Sentinel-2 du pôle THEIA, les produits de Niveau 3A n'ont pas encore été abordés dans ce blog. Ceux-ci fournissent une synthèse mensuelle des produits de niveau 2A. Les produits de Niveau 3A devraient être bien utiles pour les raisons suivantes :

  • Le produit de Niveau 3A, fourni une fois par mois, représente un volume de données jusqu'à 6 fois inférieur à celui des niveaux 2A acquis pendant un mois.
  • Le produit de niveau 3 fournit un échantillonnage temporel régulier de l'évolution des réflectances, alors que l'échantillonnage du Niveau 2A est dépendant de la couverture nuageuse
  • Beaucoup de méthodes de traitement sont perturbées par la présence de trous dans les images dus à la couverture nuageuse. Le Niveau 3A vise à minimiser les pixels manquants.

 

Grâce au jeu de données SPOT4(Take5), nous avons pu expérimenter différentes méthodes de création des produits de Niveau 3A sur différents types de sites. C'est Mohamed Kadiri, au CESBIO, qui a pris en charge ce travail, financé par le budget CNES du projet THEIA, et avec le soutien de Mirelle Huc. La méthode consiste à calculer, pour chaque pixel, la moyenne des réflectances de surface des observations non nuageuses, obtenues pendant une demi période de N jours autour de la date T0 du produit de niveau 3A. Par exemple, dans l'exemple ci-dessous obtenu avec N=21, pour le produit du 15 mars (T0), les données utilisées couvrent la période du 24 février au 4 avril.

 

 

La moyenne calculée est une moyenne pondérée, qui attribue plus de poids :

  • aux images sans nuages
  • aux pixels situés loin de nuages
  • aux images avec peu d'aérosols
  • aux images proches de la date du produit de Niveau 3

 

Les poids et surtout la demi-période N ont une forte influence sur la qualité des résultats. Pour déterminer leurs valeurs, nous avons mis au point 3 critères de qualité :

  • le pourcentage de pixels dont la réflectance est manquante en raison des nuages
  • la différence du produit de Niveau 3A avec un produit de Niveau 2A faiblement nuageux acquis à quelques jours de la date centrale T0.
  • une mesure des artefacts. Les artefacts apparaissent en bordure des trous (dus aux nuages ou aux ombres) présents sur chacune des images.

 

Voici par exemple les performances obtenues sur le site de Versailles, site fortement nuageux lors du printemps 2013. Cette courbe confirme que malgré le mauvais  temps, Sentinel-2 devrait pouvoir fournir des produits quasiment sans nuages chaque mois sur ce genre de sites. Sur ce site, l'optimum de la durée de synthèse est compris entre 2*21 et 2*28 jours :

Performances obtenues sur le site de Versailles, pour plusieurs valeurs de la longueur de période N. En rouge, les trous résiduels (échelle à gauche), en jaune et vert, l'écart maximal à l'image centrale exprimé en réflectance, pour les 70% et 90% de pixels les meilleurs, et en bleu l'écart-type moyen des artefacts

 

Les produits de niveau 3A de Sentinel-2 devront aussi inclure une correction des effets directionnels, afin de pouvoir inclure dans un même produit de Niveau 3A des données issues de traces orbitales différentes, et donc acquises avec des angles de prise de vue différents. Enfin, en option, nous pourrons proposer une opération de bouchage des trous résiduels par interpolation temporelle ou en utilisant le comportement moyen de pixels similaires. Bref, il nous reste du pain sur la planche. Une comparaison avec la méthode classique de maximum de NDVI est fournie dans cet article.

SMAC prend en charge de nouveaux satellites

=>

De nouveaux coefficients ont été ajoutés au site du CESBIO . Les nouveaux satellites pris en compte sont :

  • Landsat 8, RapidEye
  • Quickbird, Worldview2, Ikonos
  • Pléiades1A (PHR1A)
  • DMC-DEIMOS1

 

Le Simplifié Modèle d'Atmosphérique Correction (SMAC) est parfaitement adapté à l'implémentation rapide et approchée de corrections atmosphériques. Il s'agit de fonctions analytiques dérivées du modèle 5S. Les 49 coefficients de ce modèle sont ajustés à partir de simulations de transfert radiatif obtenues avec le modèle 6S (l'ancienne version, pas la récente version vectorielle). SMAC n'est pas un modèle très précis (beaucoup moins que MACCS), et il faut lui fournir des données auxiliaires pour l'épaisseur optique des aérosols ou pour les contenus atmosphériques en ozone et vapeur d'eau. Quand ces données sont précisément connues, la précision des simulations est en général meilleure que deux à trois pour cent, sauf parfois pour les grands angles (au dessus de 70°) ou dans de fortes bandes d'absorption et si on ne prend pas en compte les effets d'environnement et les effets de pente.

 

SMAC a été conçu pour être facile à utiliser :

#lecture des 49 cofficients
nom_smac ='COEFS/coef_FORMOSAT2_B1_CONT.dat'
coefs=coeff(nom_smac)
 
#Lire la réflectance TOA de la bande à traiter
#(celà va dépendre du format de l'image)
 
#Lire les valeurs des angles dans les métadonnées (ici, on les fixe)
theta_s=30
phi_s=180
theta_v=0
phi_v=0
# calculer la pression atmosphérique à l'altitude du pixel
pressure=PdeZ(1300)
 
#trouver les valeurs des variables atmsophériques (ici, on les fixe arbitrairement)
AOT550=0.1
UO3=0.3
UH2O=3
 
#calculer la correction atmosphérique
r_surf=smac_inv(r_toa,theta_s,phi_s,theta_v,phi_v,pressure,AOT,UO3,UH2O,coefs)

Dans la dernière ligne ci-dessus :

  • theta_s, phi_s sont resp. l'angle solaire zenithal et azimuthal
  • theta_v, phi_v sont resp. l'angle de visée zenithal et azimutha
  • AOT est l'épaisseur optique à 550 nm qui peut provenir d'une station Aeronet, ou fixée au jugé, ou égale à 0.1 (si on se contente d'une correction très approchée).
  • UO3 est le contenu en ozone, en cm.atm (0.3 est souvent suffisamment précis)
  • UH2O est le contenu en vapeur d'eau, en kg/m². J'utilise souvent une valeur égale à 3, quand je ne cherche pas une grande précision, dans des bandes ou l'absorption est faible
References

[1] Rahman, H., & Dedieu, G. (1994). SMAC: a simplified method for the atmospheric correction of satellite measurements in the solar spectrum. REMOTE SENSING, 15(1), 123-143.
"[2]"Tanré, D., Deroo, C., Duhaut, P., Herman, M., Morcrette, J. J., Perbos, J., & Deschamps, P. Y. (1990). Technical note Description of a computer code to simulate the satellite signal in the solar spectrum: the 5S code. International Journal of Remote Sensing, 11(4), 659-668.

"[3]"Vermote, E. F., Tanré, D., Deuze, J. L., Herman, M., & Morcette, J. J. (1997). Second simulation of the satellite signal in the solar spectrum, 6S: An overview. Geoscience and Remote Sensing, IEEE Transactions on, 35(3), 675-686.>
"[4]"Kotchenova, S. Y., Vermote, E. F., Matarrese, R., & Klemm Jr, F. J. (2006). Validation of a vector version of the 6S radiative transfer code for atmospheric correction of satellite data. Part I: Path radiance. Applied Optics, 45(26), 6762-6774.
"[5]"Kotchenova, S. Y., & Vermote, E. F. (2007). Validation of a vector version of the 6S radiative transfer code for atmospheric correction of satellite data. Part II. Homogeneous Lambertian and anisotropic surfaces. Applied Optics, 46(20), 4455-4464.

 

New satellites added to SMAC atmospheric correction

=>

New coefficients have been added to the CESBIO repository for SMAC coefficients. The new sensors taken into account are :

  • Landsat 8, RapidEye
  • Quickbird, Worldview2, Ikonos
  • Pléiades1A (PHR1A)
  • DMC-DEIMOS1

 

The Simplified Model for Atmospheric Correction (SMAC) is the perfect model to perform easy, quick and not too dirty atmospheric corrections. It is based on very simple analytic formulas, based on the 5S model. The 49 coefficients of this model are fitted using a large number of radiative transfer simulations with the 6S model (the old historic version, not the recent vector version). This software is not very accurate (much less than MACCS), and it requires in-situ measurements for the aerosol optical thickness, and weather analyses for ozone and water vapour. If these data are available,  in most cases, its accuracy is within 2 and 3 percent, if we do not account for adjacency effects and slope effects, and it may be worse for large viewing and solar angles (above 70°) or within strong absorption bands.

 

SMAC is very easy to use:

#read the 49 coefficients in smac_soefs table
nom_smac ='COEFS/coef_FORMOSAT2_B1_CONT.dat'
coefs=coeff(nom_smac)
 
#read the TOA reflectance image in r_toa variable
#depends on the file format
 
#read the angle values in the image metadata
theta_s=30
phi_s=180
theta_v=0
phi_v=0
# compute pressure at pixel altitude
pressure=PdeZ(1300)
 
#find the values of AOT, UO3, UH2O
AOT550=0.1
UO3=0.3
UH2O=3
 
#compute the atmospheric correction
r_surf=smac_inv(r_toa,theta_s,phi_s,theta_v,phi_v,pressure,AOT,UO3,UH2O,coefs)

where :

  • theta_s, phi_s are the solar zenith and azimuth angles
  • theta_v, phi_v are the viewing zenith and azimuth angles
  • AOT is the aerosol optical thickness at 550 nm which may be obtained from an Aeronet stations, or guessed, or equal to 0.1 for a really dirty atmospheric correction.
  • UO3 is the ozone content in cm.atm (0.3 is OK)
  • UH2O is the water vapour integrated content in kg/m². When I do quick and dirty atmospheric correction, I often use a value equal to 3, but I do not process spectral bands with strong water vapour absorption bands.
References

[1] Rahman, H., & Dedieu, G. (1994). SMAC: a simplified method for the atmospheric correction of satellite measurements in the solar spectrum. REMOTE SENSING, 15(1), 123-143.
"[2]"Tanré, D., Deroo, C., Duhaut, P., Herman, M., Morcrette, J. J., Perbos, J., & Deschamps, P. Y. (1990). Technical note Description of a computer code to simulate the satellite signal in the solar spectrum: the 5S code. International Journal of Remote Sensing, 11(4), 659-668.

"[3]"Vermote, E. F., Tanré, D., Deuze, J. L., Herman, M., & Morcette, J. J. (1997). Second simulation of the satellite signal in the solar spectrum, 6S: An overview. Geoscience and Remote Sensing, IEEE Transactions on, 35(3), 675-686.>
"[4]"Kotchenova, S. Y., Vermote, E. F., Matarrese, R., & Klemm Jr, F. J. (2006). Validation of a vector version of the 6S radiative transfer code for atmospheric correction of satellite data. Part I: Path radiance. Applied Optics, 45(26), 6762-6774.
"[5]"Kotchenova, S. Y., & Vermote, E. F. (2007). Validation of a vector version of the 6S radiative transfer code for atmospheric correction of satellite data. Part II. Homogeneous Lambertian and anisotropic surfaces. Applied Optics, 46(20), 4455-4464.

 

Schedule for SPOT4 (Take5) reprocessing

We should be able to start reprocessing  SPOT4(Take5) data in January, and to release them in February (it is always risky to announce a date, we have already announced November, then December). This new version will also be produced on MUSCATE prototype processing center, implemented at CNES, for the THEIA Land data center.

 

Why does it take so long ?

 

  • The longest part was the negociation between CNES and CAP GEMINI about the modifications of MUSCATE.
  • The modification itself took a while, as it includes :
    • A major version change of SIGMA, which is CNES software for ortho-rectification. This new version will allow us to improve geometric accuracy using new reference images to take ground control points (GCPs), while the previous version only allowed us to work with LANDSAT 5 or 7, while the other solutions did not work because of a strange bug. Instead of LANDSAT 5 and 7, we will use the GEOSUD mosaic over France, which is a high quality mosaic produced by IGN, that covers the whole France. Elsewhere, we will use LANDSAT 8, which also provides a large gain in accuracy. And for the 4 sites, whose ortho-rectification was most difficult, we will the most cloudfree SPOT4(Take5) image as référence, which should provide us better quality GCPs.
    • We will be able to tell normal Level 1C failures due to cloud cover from possibly abnormal ones.
    • We will at last be able to process NASA's maricopa site, y separating it in two sites, one taken from the East, and the other from the West.
    • The level 2A chain was updated, and we will use a new aerosol model, with slightly larger particle sizes, which provide better results.
    • The Level 2 will provide 2 new flags that indicate the pixels for which the slope correction cannot be performed, for instance because these pixels lie in the shadow. We had forgotten these flags in our first version.
    • I do not know whether you are using the saturation masks or not. But these ones had a few defects that will be corrected in this new version.
    • The level 2A products will also contain nice quicklooks, on which clouds and shadow will be circled in colors. We will use them to perform a quicker validation of the processing, and we hope they will be useful for you too.
  • Until the end of the year, we will install the new version of MUSCATE prototype, then the processors, the parameters, and the new image references. The preparation of 45 reference images is quite some work.
  • Then, banzaï !,  we will validate all this stuff and start the reprocessing of version 2.0 of  SPOT4 (Take5) data.