Outils de traitement des images satellites dans le cadre de

l'estimation des émissions de carbone dues à la

déforestation et à la dégradation des forêts​.

Préparé par Bruno COMBAL, Klara Janouskova et Hugh EVA

Centre Commun de Recherche, Commission Européenne, Ispra, Italie
 Activités de l'atelier
​1.​ Les images de satellites d'observation de la Terre 3

​2.​ Le logiciel Impact Toolbox 4

​3.​ Importer des images satellites 4

​4.​ Traitements d'images 6

​5.1.​ Combinaisons de bandes spectrales 6
​5.1.1.​ ND(V,W,S,I) tool 7
​5.1.2.​ Raster calculator 8
​5.2.​ Classifications 13
​5.3.​ Classification automatique 13

​6.​ K-Means (K-moyennes) et PCA (Analyse en Composante Principale) 13

​6.1.​ Principe 13
​6.2.​ Recoder le k-means 15

​7.​ Linear Spectral unmixing (séparation spectrale linéaire) 17

​8.​ Détection de changements 18

​8.1.​ Détection de coupes forestières 19
​8.1.1.​ Détection de changement 19
​8.1.2.​ Masque de forêt 19
​8.2.​ Détection d'exploitations minières 22

​9.​ CarbEF: Estimation des émissions de carbones dues à la déforestation et à la

dégradation 26
​9.1.​ La carte de changement 27
​9.2.​ Création d'une carte de la biomasse aérienne 28
​9.3.​ La carte d'exception 29
​9.3.​ Estimation des émissions de carbone 29

Note sur les exercices

Ce guide utilise des données, disponibles dans 4 répertoires distincts. Ils contiennent les
images satellites nécessaires à la réalisation des exercices : vous trouverez au début de
chaque exercice le nom des données à utiliser. Ces données sont disponibles depuis notre
serveur, sous deux formes : image brutes complètes (à importer avec la fonction
d'importation), ou une zone plus petite (et donc un fichier plus petit à télécharger)
découpée par nos soin sur l'image initiale. Notez que les images extraites n'ont pas besoin
d'être importées dans Impact au moyen de la fonction d'importation, il vous suffit de
décompresser dans le répertoire DATA d'Impact.

1/31
Note: cet atelier a été réalisé avec la version 4.2.2 de Impact toolbox.

​1.​ Les images de satellites d'observation de la Terre

Nous nous intéressons ici aux satellites d'orbite basses polaires (par opposition aux
satellites géostationnaires). Ces satellites passent au dessus du même point après un
certain nombre d'orbites, ce qui définit la fréquence de revisite. Pour Landsat 8 elle est de 16
jours, et pour Sentinel elle est de 5 jours avec deux satellites (S2-A et S2-B).

Les propriétés comme la résolution au sol des images et le nombres de bandes spectrales
est une propriété des capteurs embarqués.
Pour Landsat 8, le capteur OLI offre 11 bandes, avec des résolutions au sol de 30m pour le
visible et le "SWIR" (Short wave infrarouge) et de 100m pour les bandes dans l'infrarouge
thermique, et 15m pour la bande panchromatique.
Pour Sentinel 2, les capteurs MSI fournissent des images hyperspectrales (13 bandes), avec
différentes résolutions en fonction de la bande spectrale, soit 10m, 20m et 60m.
Nous verrons à la section 5 plus de détails sur les propriétés spectrales.

Landsat OLI
Band µm Résolution
Number
1 0.433–0.453, bleu profond, UV 30 m
2 0.450–0.515, bleu 30 m
3 0.525–0.600, vert 30 m
4 0.630–0.680, rouge 30 m
5 0.845–0.885, PIR 30 m
6 1.560–1.660, SWIR 1 30 m
7 2.100–2.300, SWIR 2 30 m
8 0.500–0.680, pan-chromatique 15 m
9 1.360–1.390, nuages 30 m
10 10.6-11.2, Infrarouge therm 1 100 m
11 11.5-12.5, Infrarouge therm 2 100 m

Sentinel 2 MSI
Resolution Band number Central Band width (nm)
(m) wavelength (nm)

10 2 496.6 (bleu) 98
3 560 (vert) 45
4 664.5 (rouge) 38
8 835.1 (PIR) 145
20 5 703.9 (red edge) 19
6 740.2 (red edge) 18
7 782.5 (red edge) 28
8a 864.8 (vegetation red edge) 33
11 1613.7 (SWIR) 143
12 2202.4 (SWIR) 242
60 1 443.9 (aerosols) 27

2/31
 9 945 (vapeur d'eau) 26
10 1373.5 (cirrus) 75

​2.​ Le logiciel Impact Toolbox

Le logiciel IMPACT Toolbox a été développé au JRC, dans le cadre du projet ReCaREDD. Il
fonctionne sous Windows. Le logiciel est conçu de manière à faciliter le traitement des
données Landsat et Sentinel-2, et reconnaît différentes versions de leurs formats.
Vous pouvez le télécharger ici:
[Link]
Vous trouverez sur cette même page une série d'exercices avec leur données, prêts au
téléchargement.

Au lancement, le logiciel contacte notre serveur pour vérifier si des mises à jour sont
disponibles. Une mise à jour ne requiert pas de télécharger l'intégralité du logiciel, mais
seulement les éléments qui ont été modifiés depuis la dernière version, afin de ne pas
consommer une trop grande bande passante internet.
NOTE​: le logiciel ouvre une fenêtre MS-DOS. Le logiciel requiert que cette fenêtre reste
ouverte (sa fermeture met fins aux fonctions d'Impact).

En cas de problème logiciel, ou de suggestion, vous pouvez contacter le développeur

principal: ​[Link]@[Link]​, ou Hugh Eva (​[Link]@[Link]​), leader du
groupe des activités pour l'Afrique.
Nous allons procéder à l'installation du logiciel sur vos PC, puis à une visite guidée de ses
fonctionnalités, avant de débuter la revue des fonctionnalités et les cas d'étude.

​3.​ Importer des images satellites

Dans les exercices suivants, les images satellites vous seront fournies. Si vous souhaitez en
obtenir de nouvelles, vous pouvez utiliser les sites suivants:
- EarthExplorer, pour les images Landsat: [Link]
- JRC CID portal, pour les images sentinel:
[Link]

Ce qu'il faut savoir :

Les fichiers d'imagerie satellite sont généralement des fichiers conteneurs compressés
([Link], zip, etc.) qui contiennent plusieurs fichiers : les fichier images et des
méta-informations relatives à la prise de vue et à l'image elle-même, et autres fichiers
auxiliaires.
Tous ces fichiers ont des formats spécifiques, c'est-à-dire que l'information qu'ils contiennent
est représentée de différentes manières. De plus, ces formats ont évolué dans le temps.
Importer une image satellite brute dans son SIG peut se révéler relativement complexe.
Impact Toolbox fait se travail pour vous. Pour chaque type d'image, il applique les fonctions
de décodages ad-hoc et transforme la donnée brute en image geotiff projetée.
- LANDSAT: plus de 15 formats de métadonnées et leur paramètres de conversion
des comptes numériques en réflectances (en particulier pour Landsat 4 et 5 qui ont
connu une dégradation des capteurs en cours de mission), qui ont changé au cours
des missions;
- Sentinel-2: 3 formats de métadonnées et la structure du fichier conteneur de l'image.

Le logiciel est régulièrement mis à jour, afin de suivre l'évolution des formats. Si le format de
vos données Landsat ou Sentinel-2 n'est pas reconnu, n'hésitez pas à nous contacter.

3/31
Ce qu'Impact fait pour vous :
- Interprète le format des fichiers;
- ainsi que le format des données (Byte, Int, float);
- Offre de convertir les images Landsat en TOA, reflectances, CN;
- Permet de clipper les images (découper une zone d'intérêt), de sélectionner des
bandes, de changer la projection;
- Réalise de nombreuses opérations d'amélioration, classification, segmentation;
- Offre le module CarbEF d'estimation des émissions de carbone.

ATTENTION
Le processus d'importation offre de ne retenir que certaines bandes. Avant de combiner les
bandes spectrales, il est essentiel de vérifier quelles bandes ont été retenues par le
processus d'importation. Pour cela, ouvrez le panneau d'information du fichier. Au bas du
panneau, vous trouverez la liste des bandes importées à la fin, sous le nom : Impact_bands.

Image 1: L'outil d'information d'un fichier raster vous donne la liste des bandes du raster et leur nom d'origine.

Exemple d'importation de données :

Données.1_ ​ACTIVITÉ ​à télécharger depuis
[Link]

- Importons la série d'image satellites Landsat, dans le répertoire activité

(​LC08_L1TP_182061_20160618_20170323_01_T1.[Link]​ et
LC08_L1TP_182061_20161211_20170316_01_T1.[Link]​). Impact vous
permet de sélectionner plusieurs images dans son interface d'importation;
- Définissez vos choix dans l'interface : DN, TOA;
- Je vous recommande de toujours utiliser le format de date YYYYMMDD. Quel est
son intérêt par rapport à DDMMYYYY?;
- Affichez les images: changez l'affichage (visualisation options);

4/31
 - Ces images vous sont-elles familières, savez-vous les interpréter?
- Vous pouvez créer un sous-répertoire dans IMPACT/DATA, par exemple
/IMPACT/DATA/activite, puis déplacer vos images à l'intérieur soit à l'aide de
l'interface d'IMPACT, soit depuis Windows.
- les bandes individuelles sont disponibles sous forme compressée (.gz), ainsi que la
bande "pansharpening". Vous pouvez supprimer ces fichiers, si vous le souhaitez et
ne conserver que le fichier tif de l'import et le fichier de méta-information MTL.

Notez que chaque fois qu'une image est crée dans Impact ou déplacée depuis Windows,
Impact analyse cette image, en particulier pour en extraire ses statistiques par bandes, et
faciliter son affichage.

​4.​ Traitements d'images

​4.1.​ Combinaisons de bandes spectrales
La figure ci-dessous représente une signature spectrale typique de la végétation.
Pouvez-vous commenter cette signature?

Image 2: Signature spectrale typique d'un couvert végétal.

Image 3: Signatures spectrales typiques d'un couvert végétal, d'un sol nu et de l'eau.

5/31
La signature spectrale change avec le type de cible (sol nu, végétation, eau, nuages,
surfaces artificielles), avec l'évolution de la cible, mais aussi avec les conditions
d'observation (angle d'éclairement du soleil, de la prise de vue, turbidité de l'atmosphère).
Ne perdez jamais de vue que la signature spectrale d'un pixel est la combinaison de tous les
éléments composant la scène couverte par un pixel. Par exemple, un pixel de couvert
végétal combine la signature du fond (le sol, nu ou en herbe) et de la végétation qui pousse
dessus. Suivant l'état de développement de la végétation, la contribution du sol dans la
signature du pixel est plus ou moins importante. Un pixel sur un lac a la signature spectrale
de l'eau, mais aussi du fond si l'eau est peu profonde, et des matériaux qui se trouvent dans
l'eau (turbidité).

Image 4: Principe du transfert radiatif

Il est possible de normaliser ces signatures, notamment en corrigeant les angles

d'illumination et d'observation (réflectance TOA). On peut également corriger les effets de
l'atmosphère, mais cela nécessite de connaître les propriétés locale de l'atmosphère, en
particulier le contenu en vapeur d'eau, en aérosol, et les pressions partielles des différents
gaz. Ces informations ne sont pas toujours accessibles. Pour minimiser les effets de
l'atmosphère, on utilise des indices spectraux, qui sont des valeurs relatives entre les
bandes spectrales, plutôt que des valeurs absolues.

​4.1.1.​ ND(V,W,S,I) tool

Cet outil permet de calculer des indices spectraux de la forme:
(B1-B2)/(B1+B2)

Par exemple, des indices célèbres basés sur cette formule sont :
NDVI: B1=proche infrarouge, B2=rouge, indice de végétation
NDWI: B1=proche infrarouge, B2=moyen infrarouge (SWIR), indice de contenu en eau des
feuilles
NDSI: B1=vert, B2=moyen infrarouge, indice de neige
NBR: B1=proche infrarouge, B2=moyen infrarouge (SWIR2) : normalised burnt ratio
etc.
Des indices plus complexes peuvent etres calculés a l'aide du raster calculator (voir le
paragraphe suivant).

6/31
Données.1_ ​ACTIVITÉ
Calculer le NDVI et NDWI pour les images Landsat que vous avez importé dans l'exemple
précédent. Attention: les numéros de bandes que vous spécifiez correspondent à l'ordre des
bandes du fichier: vérifiez à quelle bande cela correspondait dans le fichier Landsat ou
Sentinel d'origine (utiliser la fonction d'information du fichier).

4.1.2.​ Raster calculator

→​
​

Cet outil permet d'aller plus loin dans l'analyse spectrale. Nous verrons également le cas
particulier des opérateurs logiques qui permettent de combiner des bandes, sous certaines
conditions, et de créer des masques.

Imaginons que nous souhaitons créer un masque de couvert végétal, sur la base de
l'hypothèse que nous conservons les pixels pour un NDVI > 0.4, c'est-à-dire :
- Le résultat est une image masque;
- Si le pixel dans l'image d'origine a un NDVI>0.4, le pixel du masque vaut 1, sinon il
vaut 0.

Données.1_ ​ACTIVITÉ
Utilisez le Raster Calculator, avec l'image de NDVI précédente comme donnée d'entrée,
représentée par la variable 'A'. L'opération est A>0.4. Les données en sortie sont des bytes
(savez-vous pourquoi?). Donnez un nom de sortie comme mask_ndvi_gt_0dot4.tif.

7/31
VISUALISEZ

Supposons maintenant que nous souhaitons créer une copie de l'image d'origine et ne
conserver que les valeurs de NDVI >0.4, les autres valant 0. On profite ici du fait que 0 est le
nul de la multiplication, et 1 le neutre, c'est à dire que l'on souhaite multiplier le masque
précédent par l'image de NDVI. Cela peut être fait en une seule opération.

8/31
Remarque: vous multipliez votre image d'entrée (NDVI) par un masque (NDVI>0.4). Le
résultat est sauvegardé dans un fichier dont le type de données doit être défini (Byte, Int16,
Float32). Dans le cas de notre exemple, le fichier de NDVI était codé sur des Float32, nous
sauvegardons le résultat en ​Float32​.
- Pouvez-vous commenter? Que se passe t-il si on sauvegarde le résultat en Byte?
- Vous pouvez créer une table de couleur pour visualiser ce résultat. Par exemple noir
pour les valeurs de 0 à 0.4, bleu à brun pour 0.4 à 0.5 et brun à vert de 0.5 à 0.9;
- J'ai un NDVI en float32, avec des valeurs entre 0 et 0.82, mais je veux sauvegarder
le résultat en Byte. Que faire? Je peux multiplier les valeurs par 255, ainsi les valeurs
codées seront entre 0 et 0.82*255 = 209. La résolution de la donnée en Byte est de
1/255, soit dans ce cas 0.82/255 = 0.0032 (1 représente 0.0032 de l'échelle
originale).

9/31
VISUALISEZ
Question: pouvez-vous mettre les valeurs rejetées à -1, plutôt que 0? Quelle précaution
prendre pour le fichier de sortie?

10/31
Décomposition des calculs :
A A>0.4 -(A<=0.4) Résultat: A(A>0.4)-(A<=0.4)

0 0 -1 -1

-0.5 0 -1 -1

0.2 0 -1 -1

0.5 1 0 0.5

0.9 1 0 0.9

Disons maintenant que nous souhaitons créer une image de moyen infrarouge (SWIR), pour
laquelle on garde la valeur de SWIR si NDVI>0.4, -1 sinon.
CODER
VISUALISER

11/31
​4.2.​ Classifications

Données​.2_CLASSIF_BRAZA
Nous allons commenter différents exemples de classifications d'une image Sentinel-2
couvrant la région de Brazzaville. L'intérêt de cette zone est de présenter une grande variété
de types d'occupation du sol, artificiels et naturels.

Image de référence, Sentinel 2 du 27 juillet 2017 :

S2B_MSIL1C_20170727T085559_N0205_R007_T33MWR_20170727T091039
Elle a été obtenue à partir du site du JRC:
[Link]

Image 5: Image Sentinel 2 de la région de Brazzaville et du Pool, depuis le site du JRC

​4.3.​ Classification automatique

Impact possède un outil de classification automatique. Celui-ci se base sur les signatures
spectrales des pixels, et sur des seuils permettant de faire la différence entre les différentes
classes. L'outil possède une option permettant de réajuster les réflectances de l'image de
fonction du spectre de référence des forêts. Cela permet d'atténuer les problèmes de
modification des réflectances par les conditions d'observation et par l'état de l'atmosphère.
1. Utilisez le classifieur automatique ​sans ​normalisation des forêts permanentes
(evergreen forest normalisation);
2. Puis réalisez la même classification ​avec ​l'option de normalisation des forêts
permanentes.
Que constatez-vous, comment interprétez-vous cette différence ?

5.​ K-Means (K-moyennes) et PCA (Analyse en Composante

Principale)
​5.1.​ Principe
Afin de classer les différentes couvertures du sol, on fait l'hypothèse que celles-ci ont des
signatures typiques. L'idée de l'algorithme du k-means (ou k-moyennes) est d'identifier des
groupes (ou clusters) de points dans l'espace à n-dimensions des n-bandes de l'image.

12/31
 Image 6: Représentation du concept des k-moyennes.

Il est bien évident que plus les nuages de points seront séparés, plus il sera simple de les
identifier, et meilleur sera le résultat de la classification. Il est donc assez courant d'utiliser
au préalable un algorithme d'analyse en composante principale PCA (Principal Component
Analysis) pour réaliser des combinaison de bandes permettant d'accroître la séparation
entre les différents clusters. La PCA crée des combinaisons de bandes dont le but est
d'augmenter la variance le long des nouveaux axes créés. On conservera par la suite
uniquement les axes qui expliquent (au sens statistique du terme) le plus la variance du jeux
de données.

Image 7: Le principe de l'analyse en composante principale est de créer de nouveaux axes maximisant la
variance des grandeurs représentées sur ces axes.

13/31
Données​.2_CLASSIF_BRAZA
Image de référence : Sentinel-2, région de Brazzaville:
S2B_MSIL1C_20170727T085559_N0205_R007_T33MWR_20170727T091039

- Faire une classification par K-means, avec les paramètres par défaut. Que
constatez-vous?
- Faire un PCA, puis un K-means. Que constatez-vous?

k-means PCA suivi de k-means

5.2.​ Recoder le k-means

Le k-means est un outil généraliste : il n'a aucune notion de la couverture du sol à laquelle
correspond chaque classe, il se contente de trier les objets en classes dont les signatures
spectrales sont distinctes.
Nous devons donc recoder le k-means. Pour cela, nous devons décider de classes de
couverture du sol que nous souhaitons garder, et de leur code de classe. En passant en
revue la classification, nous en déduirons à quelle couverture du sol correspond chaque
classe de k-means. Pour cela, nous préparons un table de couverture du sol, inspirée de
celle défini pour la classification automatique. Notez que nous ne pouvons pas attribuer plus
de couverture du sol qu'il n'y a de classes utilisées pour le k-means, soit 10 dans notre cas.
Par convention, 0 correspond à l'absence de données (no-data).

Notez que le k-means utilise un générateur aléatoire pour débuter son cycle itératif de
classification: vos valeurs de classes ne correspondent pas forcément à celle de
votre voisin, et peuvent changer si vous recommencez.

Nom de la classe Valeur cible Valeur(s) d'origine (​exemples​)

No data 0 2

Nuage 1 -

Eau 2 5

Routes/surfaces artificialisées 3 8

Surface artificialisée/sol nus 4 6

Sol nu 5 10, 3

Herbacées 6 1, 4, 9

14/31
 Foret fermee 9 7

Données​.2 CLASSIF_BRAZA
Vous pouvez à présent recoder votre classification, en utilisant la fonction "Start Raster
Recode" du menu contextuel du raster. Utilisez la fonction "use unique value" et reportez les
valeurs de votre tableau.

ATTENTION ​: nous vous conseillons de sauvegarder le recodage dans un nouveau fichier,

plutôt que de recoder le même fichier, la moindre erreur risquant de provoquer une
confusion dans les classes.

Maintenant, vous pouvez créer une nouvelle légende colorée, et l'appliquer à l'image. Faites
un clic droit sur l'image, choisissez "define new style" puis choisissez un nom pour ce style.
Vous allez ensuite ajouter des classes définies par une valeur de départ puis de fin, et une
couleur de départ puis de fin (les couleurs sont interpolées). Vous retrouverez toutes les
légendes dans la partie "settings" du logiciel.

15/31
 Image 8: Exemple de classification colorée et la classification recolorée correspondante.

6.​ Linear Spectral unmixing (séparation spectrale linéaire)

Comme nous l'avons déjà évoqué, un pixel correspondant à un certain type de couverture
(forêt, sol nus, eau, couvert agricole, zone construite) n'est pas uniquement fait de ce type
de couverture. Un pixel de zone forestière peut-être plus ou moins fermé, et la contribution
du sous bois et du sol peut être plus ou moins importante sur la signature spectrale totale du
pixel.
L'idée de la séparation spectrale est d'exprimer le spectre du pixel comme une combinaison
linéaire de signatures typiques de différents objets, à savoir: sol nu, eau, et végétation (Les
valeurs pour Landsat sont les mêmes que pour Sentinel-2). Il est recommandé d'utiliser la
fonction de normalisation des forêts permanentes (sur les régions où elles sont présentes).

Données.1 ​ACTIVITÉ
Nous utiliserons ici l'image Landsat: ​oli_182-061_20160618
Importez l'image, si ce n'est pas déjà fait, puis lancez l'opération de séparation spectrale, en
choisissant la fonction de normalisation des forêts permanentes.
Le résultat de l'opération de séparation spectrale est un nombre de 0 à 255 indiquant la
proportion de chacun des "endmember" dans le spectre du pixel (255 indique 100%, soit un
pixel ne contenant que ce "endmember").
Les proportions sont organisées comme suit:
Bande 1 (rouge): proportion de sol nu
Bande 2 (vert): proportion de végétation
Bande 3 (bleu): proportion d'eau
Ainsi, en affichant l'image de séparation spectrale, on a immédiatement une idée du type de
couverture du sol. Bien entendu, les couleurs correspondant à des mélanges de Rouge, Vert
et Bleu correspondent à des combinaisons de sol, végétation et eau.

16/31
 Image 9: Image Sentinel après une séparation spectrale linéaire.

De manière pratique, la décomposition de la signature de chaque pixel en 3 signatures

typiques ne peut pas toujours donner 3 proportions dont la somme totale est 100%.
La principale raison est qu'il n'est pas toujours possible de décomposer le signal d'un pixel
sur cet ensemble de trois signatures typiques, soit parce que il y a autre chose dans le pixel,
soit parce que la signature de chaque élément du pixel n'est pas une transformation linéaire
de la signature du end-member correspondant. Il peut en résulter une proportion négative (il
faut retirer une fraction du signal d'un des end-members pour que la somme totale se
rapproche du signal du pixel).
Le fichier image résultat code donc les valeurs de proportion attendues, etre 0 et 100%,
ainsi que les valeurs hors bornes, c'est-à-dire excédant 100% ou inférieures à 0%.
La représentation des valeurs est donc la suivante:

Proportion de signature du end-member Codage dans le fichier résultat

]-∞, 0] [0, 100]

]0, 1] ]100, 200]

]1, +∞[ ]200, 255]

En résumé, si vous affichez votre image avec un rendu RGB entre 0 et 255, toutes les
valeurs possibles seront affichées. En limitant le rendu sur chacun des canaux (R, G et B)
entre 100 et 200, vous n'afficher que les valeurs correspondant à une proportion entre 0 et
100%.

​7.​ Détection de changements

Une des activités principales requise dans l'estimation des émissions de carbones dues aux
coupes forestières nécessite de détecter les zones déboisées entre deux dates
d'observations.

17/31
Nous bénéficions de séries temporelles qui nous permettent d'observer des changements de
l'occupation du sol. Nous allons examiner différentes approches permettant de capturer de
tels changements.

​7.1.​ Détection de coupes forestières

Données.1 ​ACTIVITÉ
Nous utiliserons ici le ΔNBR. Cette approche met en valeur des perturbations dans le
couvert arboré, qui peuvent être liées à des coupes sélectives, en calculant la différence de
NBR (voir section 5.1.1) entre deux dates. Nous nous proposons d'observer l'évolution des
coupes forestières, dans la zone correspondant au fichier '[Link]', entre juin et
décembre 2016.

Region: "activite"
Images:​ oli_182-061_20160618, oli_182-061_20161211

​7.1.1.​ Détection de changement

Le module ΔNBR requiert que la paire d'images soit superposable, pixel à pixel. Pour cela,
nous allons découper nos deux images sur une zone commune, en utilisant la fonction "clip"
et le shapefile "activite".
Affichez successivement les deux images : voyez-vous les zones qui ont été déboisées?
Vous pouvez créer votre propre zone d'intérêt avec "draw rectangle".

On calcule le ΔNBR. Le NBR utilise le SWIR2 (band 7 de Landsat) et le PIR (bande 5), qui
sont les bandes 6 et 4 dans notre fichier. Comment le savons nous?
Le résultat du calcul est disponible sous la forme de deux fichiers: un fichier contenant les
valeurs du ΔNBR (nombre code entre -1 et 0, sur un Float32), et un fichier de classes, avec
les classes suivantes:
Classe 10: pas de changement, ΔNBR >0.04
Classe 30: changement, ΔNBR <= -0.04
On change l'affichage (-1 à -0.05). Quels sont vos commentaires?
On peut également afficher l'image de classes.

​7.1.2.​ Masque de forêt

On observe également que des détections ont lieu sur des zones d'activités agricoles, entre
autres, ce qui est normal car les couverts agricoles sont généralement des plantes annuelles
qui ont vocations à être coupées une fois l'an. Nous devons donc nous concentrer sur les
zones forestières.
Il existe de nombreuses manières de discriminer la forêt dans un paysage, tout dépend des
données dont vous disposez, des conditions d'observations, du paysage et de vos
connaissances du terrain.
Utilisons par exemple le moyen infrarouge (SWIR, bande 5 de notre exemple). Les forêts
apparaissent très sombres dans cette longueur d'onde. Appliquons un k-means sur 5
niveaux à la bande SWIR (ne prenons pas trop peu de classes pour éviter les confusions de
couverture du sol) :
- Extraire la bande SWIR avec raster conversion, sur l'image de la première date;
- Appliquer le k-means avec 5 niveaux.

18/31
Prenons note des valeurs de classes qui font partie de la forêt (notons les '1') et celles qui
n'en font pas partie (notons les '0').

Gardez bien à l'esprit que le fond d'image Google ne correspond pas à la date de
l'image que vous analysez. Vous pouvez l'utiliser pour comprendre le contexte, mais
les valeurs de la date d'observation sont celles de l'image d'origine.
Nous constatons que notre masque présente des "trous", car sa création est basée sur la
valeur du SWIR de chaque pixel. L'outil de segmentation peut nous aider à détourer les
formes géographiques de la forêt.

19/31
L'outil de segmentation nous permet d'affecter la classe de chaque polygone en fonction de
la classification que nous avons réalisé en utilisant le k-means sur la bande SWIR :
- Baatz;
- input : la bande SWIR;
- classification: le masque de forêt;
- Bande: 1;
- MMU: 20.

Ainsi, nous obtenons un domaine forestier fermé et continu, sans effet de pixelisation.

20/31
Nous pouvons à présent exporter ce shapefile sous forme de raster afin de l'utiliser comme
masque du domaine forestier pour les résultats du ΔNBR. Utilisez vector to raster
conversion, et imposez le fichier ΔNBR comme grille de référence (savez-vous pourquoi?).
Pour masquer le ΔNBR, il suffit pour cela d'utiliser l'outil band calculator avec le fichier de
ΔNBR et le masque de forêt que nous venons de créer.

Afin de filtrer les résultats de ΔNBR, il suffit d'utiliser band calculator et de multiplier le ΔNBR
par le masque. Si vous masquez le ΔNBR, le résultat est sous forme de float32, si vous
masquez les classes de ΔNBR, le résultat est un Byte.
La fonction show histogram statistics vous donnera la proportion de pixels dans chaque
classe, vous pourrez en déduire la surface de dégradation.

​7.2.​ Détection d'exploitations minières

Données.3 ​MINES
L'exploitation des mines a débuté aux environs de mars 2016, au nord de la ville de Kelle,
pour cela les environs ont été déboisés. Les exploitations se trouvent près de rivières et leur
impact est visible sur les images satellites.

21/31
Les images que nous utiliserons seront celles de février 2016 (avant les activités) et de la fin
de l'année.D'autres images vous montrent que les activités se poursuivent. Les images sont
des imports de Sentinel-2, avec les 6 bandes: B2, B3, B4, B8, B11 et B12 (bleu, vert, rouge,
proche infrarouge, et les deux moyens infrarouge):
Nous traiterons: ​S2A_MSIL1C_20160201mines.tif et
S2A_MSIL1C_20161217mines.tif
Les autres images nous permettent d'apprécier visuellement l'évolution:
S2A_MSIL1C_20170814mines.tif
S2A_MSIL1C_20180121mines.tif

Image du 1er février 2016 et du 17 décembre 2016: les mines suivent le cours d'eau
(bandes 5, 4, 3).

22/31
La bande de moyen infrarouge (bande 6) vous montre également les routes d'exploitation
forestière.
Nous vous proposons de vous consacrer à une zone plus réduite, en "clippant" les images
avec le fichier shapefile "[Link]".
Nous allons maintenant estimer la dégradation en utilisant l'outil de ΔNBR sur les deux dates
de 2016. Par défaut, le kernel est de 7 pixels pour Landsat (soit 7*30m=210m), nous le
fixons donc a 21 (soit 21*10m=210m) pour conserver les mêmes propriétés.
Le résultat est une image pour laquelle le code 10 représente une absence de changement
(ΔNBR<0.04), et le code 30 un changement (ΔNBR>=0.04). Le code 40 est réservé pour les
valeurs hors masque.

Le ΔNBR repère les mines et les aires de déboisement. Malheureusement, les nuages de la
seconde image apparaissent également. Impact vous permet d'éditer le ΔNBR à la main, et
ainsi d'enlever les nuages en les encodant en forêts. Pour cela nous :
23/31
 - définissons une zone d'intérêt avec l'outil "draw polygone";
- faisons un clic droit sur le résultat du ΔNBR, choisissons "start raster recode", puis
recodons l'intervalle 30-31 en 10. Lorsque Impact vous demande s'il faut effacer les
éléments temporaires (delete temporary features), répondez non, de manière à ce
que le polygone que vous venez de tracer rester en mémoire et soit disponible pour
recoder uniquement dans la zone tracée;
- Dans l'interface de recodage, choisissez la bande 1, puis cliquez sur "use unique
values", ce qui fera apparaître toutes les valeurs disponible dans l'image (dans notre
cas 10 et 30);
- Indique que 30 est recodé en 10, et n'oubliez pas de cliquer sur "Recode Only Within
AOI: YES" de manière à ce que le recodages des valeurs de la classe 30 en 10 n'ait
lieu que dans la zone que vous avez dessiné.
- Choisissez soit de recoder directement le fichier, soir de sauvegarder le résultat du
recodage dans un nouveau fichier (pour être certain de ne pas écraser les valeurs
d'origine par une mauvaise manipulation): indique donc soit "Save recoded band to
new file: YES or NO" suivant que vous souhaitez écrire directement dans le fichier
d'origine (YES) ou dans un nouveau fichier (NO).

Il reste quelques éléments liés à des fluctuations du signal, que nous allons éliminer.
Nous lançons une segmentation de l'image, avec le but d'isoler les plus gros éléments
cohérents, et en remplissant les polygones résultants avec la majorité des pixels de l'image
de ΔNBR. Ainsi, le bruit sera "noyé" dans la forêt, à l'échelle des polygones. Nous lançons
une segmentation avec les paramètres suivants :
- en entrée l'image Sentinel correspondant à la seconde date;
- La classification ΔNBR;
- Bandes: R=5, G=4, B=3;
- Algorithm: Baatz;
- MMU: 10 (une valeur plus grande permettra des calculs plus rapide, avec des
résultats plus grossiers);
- similarity: 0.9
NOTE : afin de gagner du temps, vous pouvez découper de nouveau l'image de ΔNBR sur
une région plus petite.

24/31
​8.​ CarbEF: Estimation des émissions de carbones dues à la
déforestation et à la dégradation
Données.3 ​EMISSIONS

→​
Dans le cadre de ReCaREDD, nous avons développé un module, CarbEF, permettant
d'estimer les émissions de C0​2 dues aux processus de déforestation et de dégradation des
forêts.

Ce module exploite des cartes de changement de l'état des forêts, sur deux périodes de
temps, codés de la manière suivante :

25/31
 0: pas de données;
1: T-T: le pixel est resté un couvert arboré sur les deux périodes;
2: NT-NT: sur les deux périodes de temps, le pixel n'était pas un couvert arboré;
3: T-NT1: le couvert arboré a disparu lors de la première période;
4: T-NT2: le couvert arboré a disparu lors de la seconde période.

Ce genre de carte est établi par l'analyse de données de télédétection, par des approches
similaires à celles que nous avons vues précédemment.

Un élément de forêt, ou Unité Minimale de Forêt (UMF) est définie suivant les règles IPCC :
- par une surface minimale (généralement 0.5ha ou 1 ha);
- un taux de couverture du sol par les couronnes d'arbres (souvent 30%);
- une hauteur minimale pour les arbres.
Si ces critères ne sont pas remplis, l'UMF est déclarée "déboisée", sinon elle est simplement
degradée. Notez que l'on ne peut pas apprécier par télédétection une variation de la hauteur
des arbres.
Le logiciel CarbEF analyse l'ensemble des données de changement, et en déduit les
surfaces affectées par une dégradation et celles affectées par une déforestation. Le logiciel
en tire un rapport de synthèse, ainsi qu'une estimation des émissions de carbone associées.
La distinction entre les UMF dégradées et déboisées est fondamentale, car si une UMF
déboisée libère la totalité de carbone qui était séquestrée par la forêt, une entité dégradée
n'en libère qu'une partie.
Il s'ensuit deux hypothèses de calcul (appelées option A et option B dans le rapport de
sortie). L'option A réalise des comptes d'aires et d'émission à l'échelle des UFM, alors que
l'option B fait les même comptages à l'échelle d'un pixel. Soit la dégradation d'une UMF
dégradée est directement proportionnelle à la surface perdue (option B), soit elle correspond
à une proportion du total de l'UMF (et ce, quelque soit la surface perdue, avant d'atteindre le
seuil limite pour être déclarée déboisée). Cette seconde approche (option B) est justifiée par
le fait que la résolution des images satellites ne permet pas d'accéder à une mesure réaliste
de la surface réellement dégradée, c'est l'ensemble de l'UMF qui est concernée.
Le logiciel associe ensuite les surfaces déboisées et dégradées à une donnée de quantité
de carbone disponible. Celle-ci peut soit être une valeur constante pour l'ensemble de la
zone d'étude, soit se présenter sous la forme de carte.
De façon optionnelle, CarbEF permet également à l'utilisateur de définir une carte de zone
d'exception, qui ne doivent pas être prises en comptes dans le calcul des surfaces et des
émissions.
D'autre part, en fournissant un fichier de type shapefile au logiciel, les valeurs de surfaces
déboisées et dégradées peut être réparti en fonctions des différentes classes de fichier
shapefile.

​8.1.​ La carte de changement

La création d'une carte de changement requiert d'analyser plusieurs années de données, et
de procéder à de laborieuses étapes de validations et corrections.
Nous utiliserons ici un extrait du produit ​Roadless​, développé au JRC, qui se base sur
l'analyse de données Landsat.
Le produit Roadless analyse l'archive Landsat et détaille les différents événements qui ont
eu lieu pour chaque pixel de l'image. Ce produit nous a été fourni déjà recode pour deux
périodes, comme indiqués ci-dessous, et correspond à l'image emissions\1_cartes de
changements\Sangha_recentchangesCongo_eq_area

26/31
 Code Avant periode 1 Période 1 Période 2
2000 2013-2014 2015-2016
0 no data no data no data

1: T-T couvert arboré (T) Couvert arboré (T) Couvert arboré (T)

2: NT-NT Pas un couvert arboré (NT) Pas un couvert arboré (NT) Pas un couvert arboré (NT)

3: T-NT1 Couvert arboré (T) Pas un couvert arboré (NT) Pas un couvert arboré (NT)

4: T-NT2 Couvert arboré (T) Couvert arboré (T) Pas un couvert arboré (NT)

Créez un style de couleur pour le fichier image et visualisez la carte.

​8.2.​ Création d'une carte de la biomasse aérienne

On dispose d'une carte de stratification du pays.
Pour chaque strate, nous avons une idée de la biomasse aérienne disponible. Nous allons
croiser ces deux informations. Le taux en carbone de chaque couche peut être attribué à
partir du tableau 1.

Table 1: Estimation du carbone des arbres dans les forêts. (Source: Table 7, Niveau Des Émissions De
Référence Pour Les Forêts De La République Du Congo, 2016)

Utilisez l'outil de ‘raster recode’ pour créer un nouveau fichier avec la biomasse attribuée à
chaque strate.

27/31
​8.3.​ La carte d'exception
On peut également visualiser les exceptions correspondants aux feux, dans le département
de la Sangha. Le travail d'identification des feux a deja ete realise, et se trouve dans le
fichier emissions\3_exceptions\burnt_area_verhegghenetal2016_roadless_class4.tif
Le feux n'a pas consume les arbres, mais seulement le sous bois. La biomasse
cartographiée à cet endroit ne doit donc pas être ajoutée au compte total des émissions.

​8.3.​ Estimation des émissions de carbone

On passe en revue les données dans Impact.
On les ingère dans CarbEF.

28/31
Les résultats sont :
- une image des classes de dégradations/déforestation
- un fichier résumé au format html;
- un fichier excel avec tous les détails de calculs
Les chemins vers les fichiers statistiques sont dans le "log monitor".

Visualisation des résultats : affichez l'image des classes de dégradation et déforestations

issues de CarbEF (voir ci-dessous), le fichier HTML de résumé et le fichier Excel des détails
de calculs (un exemple de résultat est disponible avec le jeux de données).

29/31
Image 10: résultat des classes d'activités des Unités Minimales de Forêt, calculées par CarbEF

--- Fin du document ---

30/31