Auteur : Paul Passy Licence : |cc_by_nc_sa| .. |cc_by_nc_sa| image:: figures/Cc-by-nc-sa_icone.png :width: 80 px .. _outils-raster: Outils raster =================== Il est très souvent nécessaire d'appliquer quelques traitements sur les données rasters avant de les interpréter ou pour en tirer des informations supplémentaires. Il existe de nombreux outils de géomatique dédiés à l'exploitation des données rasters. Nous verrons ici quelques outils très utiles et disponibles dans différents logiciels. Nous détaillerons par exemple comment découper ou masquer un raster ou comment appliquer des opérations algébriques sur une ou plusieurs couches rasters. .. contents:: Table des matières :local: .. _histogramme-raster: Histogramme d'un raster --------------------------- Il est souvent utile de regarder la forme de la distribution des valeurs des pixels d'un raster donné en traçant ce que nous appelons l'histogramme du raster considéré. Sur ce type d'histogramme, il est possible de connaître le nombre de pixels se trouvant dans chaque tranche de valeurs. Ces tranches de valeurs sont en général définies par l'utilisateur au moment de la construction de l'histogramme. La construction de ce genre d'histogramme est possible dans tous les logiciels de géomatique et dans les différents langages de programmation. Dans les différents exemples ci-après, nous tracerons l'histogramme du MNT du `bassin de la Roya`_ (Alpes-Maritime) issu du SRTM et nommé *srtm_roya_L93.tif*. .. _histogramme-raster-qgis: Histogramme d'un raster avec QGIS ************************************ Version de QGIS : 3.26.1 Tout d'abord nous chargeons le raster à étudier dans QGIS. Il est possible d'avoir un premier aperçu de l'histogramme du raster en allant dans ses ``Propriétés`` comme présenté dans la partie dédiée à l'import de raster avec QGIS (:ref:`explorer-raster-qgis`). Cependant, cet histogramme est assez sommaire et il existe une autre fonctionnalité de QGIS pour tracer un histogramme plus ergonomique. Nous choisissons le menu :menuselection:`Boîte à outils de traitements --> Points --> Histogramme d'une couche raster`. La fenêtre suivante s'affiche (:numref:`histo-raster-qgis`). .. figure:: figures/fen_qgis_histo_raster.png :width: 30em :align: center :alt: alternate text :figclass: align-center :name: histo-raster-qgis Menu de création d'un histogramme raster dans QGIS. À la ligne ``Couche source`` nous indiquons le raster à analyser, ici *srtm_roya_L93.tif*. Dans le menu ``Numéro de bande`` nous indiquons la bande à analyser, dans le cas d'un raster multi-bandes. Pour le menu ``Nombre de classes``, nous indiquons le nombre d'intervalles à considérer dans l'histogramme, nous pouvons mettre *100* ici. Enfin, à la ligne ``Histogramme`` nous indiquons un chemin et un nom pour l'histogramme qui sera produit. Remarquons que l'histogramme produit est au format *html*. En effet, l'histogramme sera tracé en utilisant la technologie `Plotly`_ qui permet de créer des graphiques dynamiques. Enfin, nous cliquons sur :guilabel:`Exécuter`. À l'exécution rien n'apparaît mais dans l'onglet ``Journal`` de la fenêtre précédente nous trouvons le chemin du fichier *html* produit (:numref:`histo-raster-chemin-qgis`). .. figure:: figures/fen_qgis_histo_raster_chemin.png :width: 35em :align: center :alt: alternate text :figclass: align-center :name: histo-raster-chemin-qgis Le chemin de l'histogramme html produit, à la ligne Output. Il suffit alors de copier le chemin, d'ouvrir un navigateur web comme *Firefox* ou autre et de coller ce chemin dans la barre d'adresse. L'histogramme s'ouvre alors dans le navigateur. Grâce à l'utilisation de Plotly, nous pouvons survoler les différents intervalles de l'histogramme pour faire apparaître les informations quantitatives sous forme d'infobulles. Il est également possible de zoomer, de se déplacer et d'exporter le graphique en *png*. .. _histogramme-raster-R: Histogramme d'un raster avec R ************************************ Version de R : 4.2.1 Il est tout à fait possible de tracer l'histogramme d'un raster avec R. Ici, nous le tracerons sur un raster stocké au format *SpatRaster* en utilisant le package *terra*. Après avoir chargé le raster dans le code R, nous utilisons simplement la fonction ``hist`` du package *terra*. Nous pouvons y ajouter des options pour augmenter le nombre d'intervalles à prendre en compte et améliorer l'habillage. Le tout est résumé dans les lignes suivantes. .. code-block:: R # chargement de la librairie terra library(terra) # chargement du raster à analyser mnt_roya <- terra::rast("srtm_roya_L93.tif") # construction de l'histogramme avec un habillage pertinent hist(mnt_roya, breaks = 100, main='Histogramme du MNT de la Roya', col='lightblue', xlab='Altitudes (m)', ylab='Fréquence') Ici, nous avons utilisé les options suivantes pour améliorer le rendu de l'histogramme : * *breaks* : le nombre d'intervalles à prendre en compte * *main* : définition du titre de la figure * *col* : définition de la `couleur des bâtons`_ * *xlab* : le titre de l'axe des abscisses * *ylab* : le titre de l'axe des ordonnées Nous obtenons au final le résultat suivant (:numref:`histo-raster-R`). .. figure:: figures/fig_R_histo_srtm_Roya.png :width: 35em :align: center :alt: alternate text :figclass: align-center :name: histo-raster-R Histogramme d'un raster calculé avec R. .. _decouper-raster: Découper un raster -------------------- Il est souvent utile de *découper* un raster afin d'appliquer des traitements seulement sur une zone en particulier. Par *découper*, nous entendons extraire une sous partie (rectangulaire) d'un raster "trop grand" par rapport à notre zone d'études. Par exemple, si nous disposons d'un modèle numérique de terrain (MNT) à l'échelle de la France et que nous souhaitons simplement étudier les Alpes, il sera pertinent d'extraire simplement la région des Alpes de notre MNT général. Faire ce type d'extraction nous permet de gagner du temps de calcul, d'alléger nos données d'entrée et de cibler spécifiquement notre méthode sur notre zone d'intérêt. Il s'agît d'une opération extrêmement courante qu'il est possible de faire avec de nombreux outils différents. Dans les exemples suivants, nous allons extraire d'une (ou plusieurs) bande spectrale de la région parisienne, un zoom sur la ville de Paris. .. _decouper-raster-dans-qgis: Découper un raster dans QGIS ****************************** Version de QGIS : 3.18.1 Nous commençons par charger dans QGIS le raster que nous souhaitons découper et nous repérons visuellement la zone à extraire (:numref:`raster-a-decouper`). .. figure:: figures/fig_raster_a_decouper.png :width: 35em :align: center :alt: alternate text :figclass: align-center :name: raster-a-decouper Le raster à découper selon le carré rouge. L'opération se fait simplement en allant dans le menu :menuselection:`Raster --> Extraction --> Découper un raster selon une emprise...` La fenêtre suivante s'affiche (:numref:`decouper-raster-qgis`). .. figure:: figures/fen_decouper_raster_qgis.png :width: 25em :align: center :alt: alternate text :figclass: align-center :name: decouper-raster-qgis Découper un raster dans QGIS. Dans le champ ``Couche source``, nous indiquons le raster à découper. Dans le champ ``Étendue de découpage``, soit nous écrivons à la main les coordonnées du rectangle à découper (ce qui est plutôt rare...), soit nous récupérons ces coordonnées selon trois façons possibles en cliquant sur l'icône |icone_browse|. Lors du clic sur cette icône, il nous est proposé soit de calculer ces coordonnées selon l'étendue d'une autre couche (``Calculer depuis la couche``), soit de faire coïncider ces coordonnées avec la vue actuelle dans QGIS (``Use Map Canvas Extent``) soit de récupérer ces coordonnées en traçant à la main un rectangle dans la fenêtre principale de QGIS (``Dessiner sur le canevas``). C'est très souvent cette dernière option qu'on utilise. Lors de la sélection de cette option, la fenêtre se cache et le curseur se transforme en croix. Il suffit alors de tracer le rectangle correspondant à la zone à extraire. Une fois ce triangle tracé, les coordonnées s'affichent dans le champ correspondant. Enfin, dans le champ ``Découpé``, nous spécifions le chemin et le nom du raster découpé qui sera créé. Nous pouvons voir dans le panneau ``Console GDAL/OGR``, qu'en fait ce module ne fait qu'interfacer une commande GDAL. Après avoir cliqué sur :guilabel:`Exécuter`, le raster découpé apparaît dans la fenêtre principale de QGIS. .. tip:: Il est tout à fait possible de découper de la même manière un raster multi-bandes, qu'il soit sous forme de raster virtuel (*.vrt*) ou de raster "en dur" (*.tif*). Dans ce cas, le raster découpé est également multi-bandes. Il est même possible de convertir un raster virtuel en tif multi-bandes et vice versa. .. _decouper-raster-lot-qgis: Découper des rasters par lot dans QGIS ****************************************** Version de QGIS : 3.18.1 Il est tout à fait possible de découper plusieurs rasters se superposant en une seule fois dans QGIS. C'est un cas qui arrive fréquemment lorsque nous travaillons avec des images de télédétection multi-spectrales par exemple. Comme dit précédemment, si nous travaillons avec un raster virtuel ou un raster multi-bandes, le découpage par lot est implicite. Ici, nous verrons comment découper plusieurs rasters indépendants mais qui se superposent. Dans cet exemple, nous allons découper plusieurs bandes spectrales Landsat centrées sur la région parisienne afin de n'avoir que des zooms sur l'étendue de la ville de Paris. Une fois les rasters chargés dans QGIS, comme pour découper un seul raster, nous allons dans le menu :menuselection:`Raster --> Extraction --> Découper un raster selon une emprise...` La fenêtre dédiée de paramétrage du découpage s'ouvre. À ce niveau, nous cliquons sur :guilabel:`Exécuter comme processus de lot...` dans le bas de la fenêtre. La fenêtre suivante apparaît (:numref:`decouper-lot-init-qgis`). .. figure:: figures/fen_decouper_par_lot_init_qgis.png :width: 45em :align: center :alt: alternate text :figclass: align-center :name: decouper-lot-init-qgis Découper des rasters par lot dans QGIS, initialisation. Chaque ligne correspondra à un raster associé aux paramètres de découpage. Il faut donc ajouter autant de lignes que de rasters à découper. Soit on ajoute les lignes une par une en cliquant sur l'icône |icone_plus| soit on sélectionne les rasters à découper parmi ceux déjà chargés dans QGIS. Nous cliquons alors sur *Auto-remplissage...* dans la colonne ``Couche source``, puis nous sélectionnons ``Sélectionner à partir des couches chargées``. Il est également possible d'ajouter des rasters depuis un répertoire sans les avoir chargés au préalable dans QGIS. Une fois cette sélection faite, autant de lignes que de rasters sont apparues, et chaque ligne correspond à un raster d'entrée. Ensuite il est nécessaire de définir l'``Étendue de découpage``. Comme dans le cas d'un seul raster, nous pouvons définir l'étendue à la main en traçant un rectangle via l'option ``Dessiner sur le canevas`` accessible via l'icône |icone_browse|. Afin d'avoir rigoureusement la même étendue de découpage pour chacun des rasters nous faisons un copier-coller des étendues de découpage de ligne en ligne. Les deux champs suivants peuvent être laissées par défaut et nous renseignons finalement le chemin et le nom des rasters découpés qui seront créés dans la colonne ``Découpé (étendue)``. Il est possible de faire des copier-coller des chemins en changeant simplement le nom du raster créé (:numref:`decouper-lot-qgis`). .. figure:: figures/fen_decouper_par_lot_qgis.png :width: 45em :align: center :alt: alternate text :figclass: align-center :name: decouper-lot-qgis Découper des rasters par lot dans QGIS. Si nous avons coché la case ``Charger les couches`` dans la fenêtre précédente, les rasters découpés apparaissent automatiquement, sinon il faut les charger manuellement. .. _decouper-raster-dans-R: Découper un raster dans R ****************************** Version de R : 3.18.1 Version de terra : 1.7.29 Il est possible de découper un raster selon une l'emprise d'une couche veceur (ou raster). Pour ce faire nous pouvons passer par la librairie *terra* et charger un raster au format *SpatRaster* et un vecteur au format *SpatVector*. Ici nous découperons une bande spectrale Landsat selon l'emprise du département des Bouches-du-Rhône. .. code-block:: R # chargement du raster et du vectuer b4 <- terra::rast('Landsat_13/LC08_L2SP_196030_20190613_20200828_02_T1_SR_B4.TIF') dep_13 <- terra::vect('dep_13_buffer_32631.gpkg') # découpage du raster selon l'emprise du vecteur b4_emprise_dep13 <- terra::crop(b4, dep13) Si maintenant le but est de vraiment masquer le raster selon les limites précises du vecteur, dans ce cas là il s'agît de découper un raster selon une couche vecteur, comme plus loin (:ref:`masquer-raster-R`). .. _mosaique-raster: Créer une mosaïque de rasters ---------------------------------- Avant de se lancer dans des traitements sur des rasters il peut être utile de créer des mosaïques de rasters. Par exemple, si vous travaillez sur un bassin-versant qui se trouve *à cheval* sur deux dalles de modèles numériques de terrain, il peut être pertinent de mosaïquer les deux dalles avant d'effectuer des traitements. C'est notamment indispensable si vous souhaitez extraire un réseau hydrographique à partir de ce MNT. Un autre cas peut être l'étude par télédétection d'une région qui se trouve à cheval sur deux images satellites. Dans ce cas une mosaïque des images peut être intéressante avant les traitements. Il faudra tout de même faire attention à ce que les images à mosaïquer soit similaires (saison, couverture nuageuse, niveau de traitement...). Pour illustrer nous allons créer une mosaïque de deux rasters multi-bandes contenant 6 bandes Landsat 8. Ces deux images Landsat 8 recouvrent respectivement le département du Val-de-Marne et la ville de Paris. Notons qu'il est plus prudent que les rasters à mosaïquer soient tous dans le même système de coordonnées. .. _mosaique-raster-qgis: Mosaïque de rasters dans QGIS ****************************** Version de QGIS : 3.18.3 Nous commençons par charger dans QGIS les deux rasters correspondant à nos deux images Landsat 8 (:numref:`raster-a-mosaiquer`). .. figure:: figures/fig_qgis_raster_a_mosaiquer.png :width: 25em :align: center :alt: alternate text :figclass: align-center :name: raster-a-mosaiquer Les deux rasters à mosaïquer. Nous voyons tout de suite que les cadres noirs vont gêner le processus de mosaïque. Ces cadres noirs sont bien sûr les pixels en *no data*. La première étape consiste à relever la valeur numérique de ces pixels. Cela nous permettra de les ignorer par la suite. Nous relevons les valeurs de ces pixels via l'outil ``Identifier des entités`` |icone_identifier| de la barre d'outils principale de QGIS. Dans notre cas nous constatons que ces *no data* ont une valeur de 0 sur toutes les bandes (ce qui n'est pas des plus judicieux mais bon...). Pour créer la mosaïque, nous allons dans le menu :menuselection:`Raster --> Divers --> Fusion...` Le menu suivant apparaît (:numref:`mosaique-param`). .. figure:: figures/fen_qgis_mosaique.png :width: 25em :align: center :alt: alternate text :figclass: align-center :name: mosaique-param Paramétrage de la mosaïque. Dans le champ ``Couches en entrée`` nous sélectionnons les rasters à mosaïquer via l'icône |icone_browse|. Il est possible, mais pas forcément utile, de cocher la case ``Placer chaque fichier en entrée dans une bande séparée``. Le point intéressant ici est de spécifier la valeur de *no data* à utiliser. Ce réglage se fait dans le champ ``Valeur de pixel à considérer comme NoData`` dans le bloc ``Paramètres avancés``. Ici nous la paramétrons à *0*. Enfin, dans le champ ``Fusionné`` nous indiquons le chemin et le nom du raster de mosaïque qui sera créé. Nous pouvons retrouver la commande GDAL sous-jacente à cette manipulation dans le bloc ``Console GDAL/OGR``. Nous cliquons ensuite sur :guilabel:`Exécuter`. La mosaïque apparaît automatiquement (:numref:`mosaique-L8`). .. figure:: figures/fig_qgis_mosaique_L8.png :width: 35em :align: center :alt: alternate text :figclass: align-center :name: mosaique-L8 Résultat de la mosaïque de deux images Landsat 8. .. tip:: Sur la figure :numref:`mosaique-L8`, dans la symbologie nous avons réglé les *no data* à 0 qui apparaissent alors en transparent, donc en blanc. Le point intéressant est qu'avec cet outil les rasters multi-bandes sont mosaïqués bande à bande. La mosaïque résultante est donc également multi-bande. .. _mosaique-raster-R: Mosaïque de rasters dans R ****************************** Version de R : 4.3.1 Version de terra : 1.7.9 Lorsque deux rasters adjacents ou se superposant sont chargés en tant qu'objets *SpatRaster*, il est possible de les mosaïquer grâce à la fonction ``merge`` de la librairie *terra*. Le mot clef *first=TRUE* signifie que là où les deux rasters se superposent, les valeurs du premier raster sont prises pour la mosaïque. Il suffit de mettre ce mot clef à *False* pour ce soient les valeurs du second raster. .. code-block:: R # chargement des deux rasters rast1 <- terra::rast('Landsat_Paris1.tif') rast2 <- terra::rast('Landsat_Paris2.tif') # mosaïque des deux rasters mos <- terra::merge(rast1, rast2, first=TRUE) # affichage de la mosaïque terra::plot(mos) .. tip:: La mosaïque avec *terra* fonctionne aussi sur des rasters multi-bandes. .. _rasters-multi-bandes: Raster multi-bandes ----------------------- La création d'un raster multi-bandes est souvent nécessaire pour différents traitements géomatiques, notamment en télédétection. Un raster multi-bandes est un fichier raster qui contient plusieurs sous rasters superposés les uns aux autres. C'est typiquement utilisé pour la manipulation d'images de télédétection multi-spectrales comme les images Landsat ou Sentinel par exemple. Selon le logiciel utilisé, la création d'un tel raster est nommé *concaténation* ou *fusion*. En franglais il est également possible de parler de *stack* (pile). Une telle concaténation est souvent nécessaire pour effectuer une classification d'image satellite, une composition colorée ou extraire des signatures spectrales. Cette manipulation est faisable avec la plupart des outils de géomatique. .. _rasters-multi-bandes-qgis: Création d'un raster multi-bandes avec QGIS ********************************************* Version de QGIS : 3.26.1 Dans QGIS, il existe trois possibilités de concaténation de rasters. Les deux premières consistent à créer un vrai raster multi-bandes *en dur* et la troisième consiste à créer un *raster virtuel* pointant vers les rasters initiaux. La différence entre les deux premières est que dans un cas nous laissons QGIS simplement numéroter les bandes et dans le second cas nous utilisons un module supplémentaire qui permet de nommer les bandes à l'intérieur du raster fusionner. .. _rasters-multi-bandes-dur: Raster multi-bandes en dur (natif) +++++++++++++++++++++++++++++++++++ .. warning:: Si vous souhaitez pouvoir nommer explicitement les bandes de votre raster fusionné, allez directement à la section :ref:`rasters-multi-bandes-dur-noms`. Dans cet exemple, nous allons concaténer six bandes spectrales Landsat en un seul raster multi-bandes. Nous commençons par charger dans QGIS nos rasters de bandes spectrales. Nous avons bien six rasters correspondant à nos six bandes spectrales. Pour concaténer ces rasters, nous allons dans le menu :menuselection:`Raster --> Divers --> Fusionner ...`. Le menu de fusion s'ouvre (:numref:`fusion-qgis-dur`). .. figure:: figures/fen_fusionner_qgis.png :width: 35em :align: center :alt: alternate text :figclass: align-center :name: fusion-qgis-dur Concaténation de rasters dans QGIS. Dans le champ ``Couches en entrée`` nous spécifions les rasters à concaténer. Il est possible de les sélectionner en cliquant sur l'icône |icone_browse| (:numref:`fusion-select`). .. figure:: figures/fen_fusionner_selection_qgis.png :width: 35em :align: center :alt: alternate text :figclass: align-center :name: fusion-select Sélection des rasters à fusionner. Ici, il faut bien faire attention à l'ordre des rasters car ce même ordre sera repris dans le raster concaténé. Nous sélectionnons les rasters à fusionner, ici nous pouvons cliquer sur :guilabel:`Sélectionner tout`. La fenêtre précédente se met à jour et nous voyons bien que dans le champ ``Couches en entrée`` nous avons *6 inputs selected*. Il est ensuite nécessaire de bien sélectionner ``Placer chaque fichier en entrée dans une bande séparée``, sinon le résultat ne sera pas celui attendu. Dans le champ ``Fusionné`` nous spécifions un chemin et un nom pour le raster concaténé résultant. Nous pouvons noter que ce menu n'est finalement qu'un interface à une commande GDAL. Nous avons accès à cette commande dans le panneau ``Commande GDAL/OGR``. En cliquant sur :guilabel:`Exécuter`, le processus se lance et le raster concaténé apparaît automatiquement dans QGIS. Nous pouvons remarquer qu'une composition colorée par défaut est automatiquement générée. .. warning:: Le numéro des bandes est un point important à bien saisir et qui prête à confusion. Le numéro des bandes du raster concaténé est totalement indépendant du "numéro" des bandes en entrée. Par exemple, si nous souhaitons concaténer les bandes 2 (Bleu), 3 (Vert) et 4 (Rouge) d'une image Landsat 8, nous obtiendrons un raster concaténé avec des bandes numérotées de 1 à 3. Mais il faudra bien faire attention au fait que la bande 1 du raster concaténé ne correspondra pas à la bande 1 de Landsat, mais à la bande 2 et ainsi de suite. .. _rasters-multi-bandes-dur-noms: Raster multi-bandes en dur (avec noms) +++++++++++++++++++++++++++++++++++++++ La solution native de QGIS (:ref:`rasters-multi-bandes-dur`) fonctionne parfaitement mais les bandes du raster fusionné ne sont pas nommées, ce qui peut s'avérer fastidueux pour, par exemple, construire des compositions colorées avec QGIS (:ref:`composition-coloree-qgis`). Il faut se souvenir à quoi correspond la *bande 1*, la *bande 2*, etc du raster multi-bandes. Ça peut être source de confusion. Pour palier cette limitation, il est possible d'utiliser le module supplémentaire nommé ``Rename Bands``. Il s'installe comme toute extension de QGIS (:ref:`extension-qgis`). La seule chose est que cette extension est considérée comme étant *expérimentale* (malgré qu'elle fonctionne très bien). Pour y avoir accès, il est donc nécessaire d'aller dans le menu :menuselection:`Extensions --> Installer/Gérer les extensions --> Onglet Paramètres`. Dans cet onglet, il faut cocher ``Afficher les extensions expérimentales``. Une fois cette option cochée, retourner dans ``Toutes`` et l'extension se trouve maintenant bien dans la liste. Concrètement, ce module va créer un raster fusionné en dur classique mais il va l'accompagner d'un fichier texte au format *.xml* qui portera le même nom que le raster fusionné et dans lequel les bandes seront nommées. Il faudra bien faire attention à toujours conserver ces deux fichiers ensemble si nous souhaitons conserver les noms des bandes. Une fois le module installé, nous commençons par charger les rasters à fusionner dans QGIS. Ensuite, le module est accessible dans la :menuselection:`Boîte à outils de traitements --> Rename bands --> Merge and rename bands`. Une fois le module lancé, la fenêtre suivante s'affiche (:numref:`fusion-qgis-noms`). .. figure:: figures/fen_qgis_fusion_noms.png :width: 25em :align: center :alt: alternate text :figclass: align-center :name: fusion-qgis-noms Création d'un raster virtuel dans QGIS en nommant les bandes. Dans le panneau ``Raster layers`` la liste des rasters chargés dans QGIS s'affiche. Nous devons sélectionner les rasters à fusionner dans ce panneau. Une fois sélectionnés ils apparaissent surlignés en bleu. Les rasters sélectionnés apparaissent dans le panneau suivant, à la colonne ``Band name``. À ce niveau il est possible de changer l'ordre des bandes à l'aide des icônes |icone_bas| et |icone_haut| et d'ajouter ou enlever un raster avec les icônes |icone_plus| et |icone_moins|. Dans ce même panneau, à la colonne ``Band future name`` nous indiquons le nom explicite à donner à chaque bande. Il suffit de cliquer sur la ligne et une fenêtre d'entrée de texte ``Band name`` surgit. Nous pouvons y inscrire ce que nous souhaitons. Ici, nous travaillons sur une image Landsat 8 (:ref:`data_Landsat_OLI`), nous choisissons donc les noms suivants *B1 côtier*, *B2 Bleu*, *B3 Vert*, etc Si nous le souhaitons, à la ligne ``No data value`` nous pouvons indiquer une valeur qui correspondra aux pixels sans données. Enfin, à la ligne ``Raster output``, nous indiquons un chemin et le nom du raster multi-bandes qui sera créé. Nous terminons en cliquant sur :guilabel:`Exécuter`. Le raster multi-bandes apparaît dans QGIS et nous constatons que les bandes sont bien nommées comme souhaité (:numref:`qgis-raster-fusion-noms`). .. figure:: figures/fig_qgis_raster_fusion_noms.png :width: 25em :align: center :alt: alternate text :figclass: align-center :name: qgis-raster-fusion-noms Raster fusionné avec noms de bandes explicites. Notons que dans le répertoire d'export, nous trouvons bien notre raster fusionné au format *.tif* ainsi que le fichier texte *.aux.xml* qui contient le nom des bandes. Par curiosité il est possible d'ouvrir ce fichier XML avec un éditeur de texte. Nous y retrouvons bien le nom de chaque bande à chaque balise ``Description`` (:numref:`qgis-raster-fusion-xml`). .. figure:: figures/fig_raster_fusion_xml.png :width: 18em :align: center :alt: alternate text :figclass: align-center :name: qgis-raster-fusion-xml Le fichier XML accompagnant le raster fusionné. .. _rasters-multi-bandes-virtuel: Raster virtuel +++++++++++++++++ Avec QGIS il est également possible de créer un *raster virtuel*. Un tel raster est qualifié de *virtuel* car le fichier généré ne sera pas un fichier raster mais un simple fichier texte pointant vers les rasters initiaux. Il s'agît d'une sorte de raccourci. L'avantage d'un raster virtuel est de ne pas dupliquer les données comme c'est le cas avec une concaténation en dur. .. warning:: Un raster virtuel ne sera souvent pas accepté comme format d'entrée pour une classification d'images. Il faudra alors préférer un vrai raster concaténé. Dans cet exemple, nous allons de nouveau concaténer nos six bandes spectrales Landsat. Une fois ces bandes chargées, nous allons dans le menu :menuselection:`Raster --> Divers --> Construire un raster virtuel`. La fenêtre *Construire un raster virtuel* s’affiche alors (:numref:`fusion-qgis-virtuel`). .. figure:: figures/fen_raster_virtuel_qgis.png :width: 25em :align: center :alt: alternate text :figclass: align-center :name: fusion-qgis-virtuel Création d'un raster virtuel dans QGIS. La procédure est très similaire à celle de la concaténation de rasters en dur. Dans le champ ``Inuput layers`` nous spécifions les rasters à concaténer de la même manière que précédemment. L'ordre des couches en entrée est important car l'ordre sera le même dans raster virtuel final. Le champ ``Resolution`` est important. Le raster virtuel homogénéisera obligatoirement les résolutions spatiales des différents rasters en entrée. Il est possible de prendre la résolution moyenne (``Average``), la plus haute (``Highest``) ou la plus basse (``Lowest``). Lorsque toutes les images ont la même résolution, comme dans cet exemple, l'option choisie n'a pas d'importance. Par contre si un des rasters a une résolution plus fine, le choix sera important. De plus, lors de l'homogénéisation des résolutions, une interpolation spatiale sera effectuée de façon implicite. Il sera possible de choisir la nature de cette interpolation dans les ``Paramètres avancés``. Il ne faut pas oublier de cocher la case ``Place each input file into a separate band``. Enfin, dans le champ ``Virtuel`` nous spécifions un chemin et un nom pour le raster virtuel résultat. Il s'agît encore d'une commande GDAL interfacée, comme nous pouvons le voir dans le panneau ``Console GDAL/OGR``. À la fin de l'exécution du menu, le raster virtuel apparaît automatiquement dans QGIS associée à une composition colorée par défaut. Comme pour le raster concaténé, il faut être prudent sur les numéros de bandes, voir l'avertissement à ce sujet dans la section précédente. Si nous regardons dans notre explorateur de fichiers, nous constatons bien que le fichier créé est en fait un simple fichier texte avec l'extension *.vrt*. Il est possible de l'ouvrir dans un éditeur de texte de type bloc-note (:numref:`vrt-texte`). .. figure:: figures/raster_virtuel_text.png :width: 55em :align: center :alt: alternate text :figclass: align-center :name: vrt-texte Raster virtuel dans un éditeur de texte. Ce fichier texte est un fichier *XML*. Nous y retrouvons les chemins vers les rasters initiaux. Par exemple, ici nous retrouvons le chemin du premier raster au niveau de la première balise *SourceFilename*. Nous retrouvons également d'autres informations relatives aux dimensions du raster, à son encodage, à la méthode de rééchantillonnage... .. warning:: Il faut être prudent lorsque nous déplaçons ou faisons un copier-coller d'un raster virtuel via l'explorateur de fichiers. Comme nous l'avons vu, il s'agît d'un raccourci pointant vers les rasters initiaux. Lorsque nous déplaçons un fichier *.vrt* il faut également déplacer les rasters initiaux. Le plus simple est de créer le raster virtuel dans le même répertoire que les rasters initiaux et déplacer tous les fichiers en même temps lorsque c'est nécessaire. .. _rasters-multi-bandes-R: Création d'un raster multi-bandes avec R ********************************************* Version de R : 4.2.1 Il est très simple de créer un raster multi-bandes avec R. Ici nous présenterons le cas de la création d'un raster contenant quatre bandes spectrales Landsat avec le package *terra*. Après avoir chargé les rasters mono bandes, nous les combinons grâce à la fonction *c()* qui signifie *concatenate*. Le tout se fait avec les lignes suivantes. .. code-block:: R # chargement du package terra library(terra) # chargement dess rasters mono-bandes à fusionner b2 <- terra::rast("Landsat_13/LC08_L2SP_196030_20190613_20200828_02_T1_SR_B2.TIF") b3 <- terra::rast("Landsat_13/LC08_L2SP_196030_20190613_20200828_02_T1_SR_B3.TIF") b4 <- terra::rast("Landsat_13/LC08_L2SP_196030_20190613_20200828_02_T1_SR_B4.TIF") b5 <- terra::rast("Landsat_13/LC08_L2SP_196030_20190613_20200828_02_T1_SR_B5.TIF") # création du raster multi-bandes multi_bandes <- c(b2, b3, b4, b5) # affichage du résultat plot(multi_bandes, col=gray.colors(100, start = 0.1, end = 0.9)) Nous obtenons bien un objet *SpatRaster* à 4 bandes (:numref:`multi-bandes-R`). .. figure:: figures/fig_R_creation_raster_multi_bandes.png :width: 40em :align: center :alt: alternate text :figclass: align-center :name: multi-bandes-R Création d'un raster multi-bandes avec R. Une fois ce raster multi-bandes créé, il est possible de changer le nom des bandes, grâce à la fonction ``names()``, comme indiqué ci-après. .. code-block:: R # renommage des bandes selon leur numéro names(multi_bandes[[1]]) <- 'Bleu' names(multi_bandes[[2]]) <- 'Vert' names(multi_bandes[[3]]) <- 'Rouge' names(multi_bandes[[4]]) <- 'Proche infrarouge' # affichage du résultat plot(multi_bandes, col=gray.colors(100, start = 0.1, end = 0.9)) Les bandes du raster multi-bandes ont bien été renommées comme souhaité (:numref:`renommer-bandes-R`). .. figure:: figures/fig_R_renommer_bandes.png :width: 40em :align: center :alt: alternate text :figclass: align-center :name: renommer-bandes-R Renommer les bandes d'un raster multi-bandes dans R. .. _defusion-raster: Défusionner un raster multi-bandes avec SCP (QGIS) *************************************************** Il est tout à fait possible de "défusionner" un raster multi-bandes afin de récupérer chaque bande dans un fichier *.tif* indépendant. Une solution consiste à passer par le module SCP. Une fois le raster multi-bandes au format .tif chargé dans QGIS, il suffit d'aller dans le menu :menuselection:`SCP --> Pré-traitement --> Séparer les bandes`. Le menu suivant s'affiche (:numref:`separer-bandes`). .. figure:: figures/fen_scp_separer_bandes.png :width: 40em :align: center :alt: alternate text :figclass: align-center :name: separer-bandes Séparer un raster multi-bandes avec SCP. Dans le champ ``Sélectionner un raster multibandes``, nous sélectionnons le raster à séparer. Si nécessaire nous actualisons la liste des rasters en cliquant sur l'icône |icone_refresh|. Nous pouvons spécifier un préfixe pour les rasters indépendants qui seront créés, mais ce champ peut rester vide. Enfin nous cliquons sur :guilabel:`Lancer`. Il nous est alors demander un répertoire pour sauver les rasters séparés. Les rasters indépendants s'affichent automatiquement. .. _raster-rapport-valeurs-uniques: Rapport de valeurs uniques ------------------------------ Dans de nombreux cas il est utile de compter le nombre d’occurrences de chaque valeur de pixels sur un raster. C'est-à-dire de lister toutes les valeurs que prennent les pixels d'un raster donné et, pour chaque valeur rencontrée, compter le nombre de fois où cette valeur apparaît. Cette manipulation est souvent appelée *Rapport de valeurs uniques*. Évidemment, ce rapport de valeurs uniques se calcule sur des rasters discrets, qui ne prennent qu'une certaine gammes de valeurs. C'est typiquement employé sur les rasters d'usage du sol. Pour chaque classe d'usage du sol, nous pouvons obtenir le nombre de pixels. Et connaissant la résolution spatiale du raster (donc la superficie d'un pixel), il est possible de calculer la superficie occupée par chaque classe d'usage du sol. Cette manipulation peut s'effectuer de différentes manières, avec différents outils. Dans les exemples présentés ci-dessous, nous calculerons le rapport de valeurs uniques sur un raster d'occupation du sol de la région du Caire. Ce raster a été obtenu par classification supervisée et présente 7 classes d'occupations du sol comme présentées dans le tableau suivant. .. list-table:: Classification visée :widths: 15 25 :header-rows: 1 * - ID classe - Label classe * - 1 - Eau * - 2 - Végétation inondée * - 3 - Végétation dense * - 4 - Sol nu sableux * - 5 - Sol nu rocheux * - 6 - Bâti dense * - 7 - Bâti diffus .. _raster-rapport-valeurs-uniques-qgis: Rapport de valeurs uniques avec QGIS ************************************** Version de QGIS : 3.18.3 Nous commençons par charger dans QGIS le raster sur lequel nous souhaitons calculer le rapport de valeurs uniques. Ensuite, nous allons dans la :menuselection:`Boîte à outils de traitements --> Analyse raster --> Rapport sur les valeurs uniques de la couche raster`. Le menu suivant apparaît (:numref:`rapport-qgis`). .. figure:: figures/fen_qgis_rapport_valeurs_uniques.png :width: 28em :align: center :alt: alternate text :figclass: align-center :name: rapport-qgis Calcul du rapport de valeurs uniques dans QGIS. Dans le champ ``Couche source`` nous indiquons le raster sur lequel nous travaillons. Dans le champ ``Numéro de bande`` nous sélectionnons la bande du raster à utiliser (seulement utile dans le cas d'un raster multi-bandes). Dans le champ ``Rapport de valeurs uniques`` nous indiquons un chemin et un nom pour le rapport qui sera créé. Ce rapport sera au format HTML et sera donc consultable avec un navigateur Web. Puis dans le champ ``Table de valeurs uniques`` nous indiquons un chemin et un nom pour la table contenant le rapport qui sera créé. Ces deux sorties sont optionnelles. Puis nous cliquons sur :guilabel:`Exécuter`. La table de rapport créée apparaît automatiquement et est nommée par défaut *report*. Il s'agît d'un fichier GeoPackage sans attribut géographique. Nous pouvons ouvrir ses attributs en faisant un clic droit dessus dans le panneau des couches et en sélectionnant ``Ouvrir la Table d'Attributs`` (:numref:`table-rapport-qgis`). .. figure:: figures/fen_qgis_table_rapport.png :width: 20em :align: center :alt: alternate text :figclass: align-center :name: table-rapport-qgis Table de rapport des valeurs uniques. La première colonne *fid* est juste une colonne *technique*. La colonne *value* contient toutes les valeurs uniques trouvées dans le raster, soit, ici, les identifiants des 7 classes d'occupations du sol. La valeur *0* correspond aux *no data*. La colonne *count* contient le nombre de pixels trouvés pour chaque valeur. La colonne *m2* est la conversion du nombre de pixels en superficie exprimée en unités du raster. Cette superficie est calculée selon la résolution spatiale du raster. Si nous l'avons exporté, nous trouvons également un fichier HTML de rapport nommé *report.html* contenant ces mêmes informations plus quelques autres comme le SCR du raster, son emprise... .. note:: Le fait d'exporter le rapport en une table GeoPackage peut être utile si nous souhaitons ensuite faire une jointure attributaire avec ce rapport. .. _raster-rapport-valeurs-uniques-scp: Rapport de valeurs uniques avec SCP ************************************** Version de QGIS : 3.18.3 Version de SCP : 7.8.21 Le module supplémentaire SCP de QGIS offre aussi une fonction pour calculer un rapport de valeurs uniques. Nous commençons par charger le raster à analyser dans QGIS. Ensuite nous allons dans le menu :menuselection:`SCP --> Post-traitement --> Rapport de classification`. Bien que la fonction soit ici nommée "Rapport de classification", elle est applicable à tout raster discret même non issu de classification. Le menu suivant s'affiche (:numref:`scp-rapport`). .. figure:: figures/fen_scp_rapport.png :width: 38em :align: center :alt: alternate text :figclass: align-center :name: scp-rapport Rapport de valeurs uniques avec SCP. Dans le champ ``Sélectionner la classification``, nous sélectionnons la classification à analyser depuis la liste déroulante. Il peut être nécessaire de rafraîchir la liste en cliquant sur l'icône. Il est possible de définir une valeur de *no data* à ne pas compter dans le rapport. Ici ce sera la valeur *0*. Il n'y a plus qu'à cliquer sur :guilabel:`Lancer`. Il est demandé un chemin et un nom pour le fichier de rapport au format csv qui sera créé. Puis le rapport apparaît sous forme de texte dans le panneau ``Sortie`` du menu précédent (:numref:`scp-rapport-sortie`). .. figure:: figures/fen_scp_rapport_sortie.png :width: 38em :align: center :alt: alternate text :figclass: align-center :name: scp-rapport-sortie Sortie du rapport de valeurs uniques avec SCP. La colonne *class* répertorie toutes les valeurs de pixels identifiées sur le raster. La colonne *PixelSum* présente le nombre d'occurrences de chaque valeur de pixel. La colonne *Percentage %* indique le pourcentage que représente chaque valeur de pixel par rapport au raster total. Enfin, la colonne *Area* traduit le nombre de pixels en superficie en mètres carrés. Ce rapport texte se retrouve exactement dans le fichier csv exporté. .. note:: Le petit plus de cette façon de procéder est de sortir un fichier csv facilement chargeable dans un tableur de type Excel ou LibreOffice Calc ou dans un script R ou Python. Le fait d'avoir directement les pourcentages est également intéressant. .. _raster-rapport-valeurs-uniques-R: Rapport de valeurs uniques avec R ************************************** Version de R : 4.8.1 Il est possible de sortir un rapport de valeurs uniques dans R à partir d'un raster. Dans cet exemple, nous allons générer un rapport de valeurs uniques pour un raster d'occupation du sol d'une région de la ville du Caire. Ce raster d'occupation du sol sera chargé sous forme d'un objet *SpatRaster* via *terra*. La fonction à utiliser pour générer ce rapport est ``freq()``. La procédure à suivre est décrite ci-après. .. code-block:: R # chargement de la librairie terra library(terra) # chargement du raster d'occupation du sol usol <- terra::rast("xtr_usol_2011.tif") # rapport de valeurs uniques val_uniques <- freq(usol) # rapport de valeurs uniques pour une valeur donnée val_uniques <- freq(usol, value=5) Le rapport est ici stocké sous forme d'un dataframe dans la variable *val_uniques*. Comme présenté dans le code ci-dessus, il est possible de ne générer ce rapport que pour une valeur donnée grâce à l'ajout de ``value=ma valeur``. .. _raster-calc: Calcul raster mono-bande --------------------------- Fréquemment il est utile d'effectuer des calculs sur une ou plusieurs couches rasters. Un tel calcul peut par exemple être utile pour isoler les altitudes supérieures à un certain seuil sur un MNT ou sur tout autre raster ou pour calculer un indice radiométrique à partir de rasters de bandes spectrales comme le NDVI. En général, il y a trois formes de calculs rasters possibles. La première consiste à combiner un raster et un scalaire par une opération arithmétique. La deuxième consiste à appliquer un opérateur logique sur un raster afin d'identifier les pixels respectant la condition. Enfin, la troisième consiste à combiner plusieurs rasters par un ou plusieurs opérateurs arithmétiques. Dans ce dernier cas, les rasters en entrée doivent, dans la plupart des cas, présenter la même résolution spatiale, la même étendue et être géoréférencés dans le même système de coordonnées. Bien sûr, ces trois formes peuvent faire intervenir plusieurs rasters ou plusieurs scalaires en même temps, ainsi que plusieurs opérateurs arithmétiques et logiques. Dans tous les cas, la logique est que le calcul s'effectue pixel par pixel comme schématisé sur la figure suivante (:numref:`calc_raster_theorie`). Sur cette figure, nous considérons le cas d'un raster composé de 9 pixels, en 3 lignes et 3 colonnes. .. figure:: figures/fig_raster_calc_theorie.png :width: 30em :align: center :alt: alternate text :figclass: align-center :name: calc_raster_theorie Les trois principaux cas de figures du calcul raster. Dans l'exemple précédent, pour le premier cas de figure, chaque pixel est multiplié par 10. Le raster résultat est géographiquement identique au raster d'entrée mais chaque pixel vaut maintenant 10 fois plus que le pixel d'entrée. Dans le deuxième cas, nous avons extrait les pixels dont la valeur est strictement supérieure à 3. Pour chaque pixel, nous avons comparé sa valeur au chiffre 3. Si le pixel était supérieur, alors on le recodait en 1 et dans le cas contraire en 0. En effet, en informatique il n'est pas possible de stocker directement les booléens *Vrai* et *Faux* sur les pixels. Les pixels sont obligatoirement de type numérique. Par défaut, *Vrai* se code en *1* et *Faux* en *0*. Mais notons bien que ces valeurs de *1* et *0* n'ont aucune signification quantitative. Enfin, par défaut (souvent mais pas toujours) à l'affichage, les pixels à *1* apparaissent en blanc et les pixels à *0* en noir. Enfin, dans le troisième cas, nous avons additionné deux rasters géographiquement identiques. Le raster produit contient pour chaque pixel la somme des pixels des deux rasters d'entrée. .. _raster-calc-qgis: Calcul raster avec QGIS *********************** Version de QGIS : 3.20.1 Le calcul raster peut se faire très simplement avec QGIS en utilisant la calculatrice raster intégrée. Elle permet de faire des calculs numériques et logiques. .. _raster-calc-logique-qgis: Calcul logique ++++++++++++++++ Dans cet exemple, nous allons extraire les zones du `bassin de la Roya`_ dont l'altitude est supérieure à 2000 m. Nous commençons par charger le MNT du bassin dans QGIS. Puis l'extraction va se faire grâce à l'opérateur logique " *>* ". Pour ouvrir la calculatrice raster nous allons dans le menu :menuselection:`Raster --> Calculatrice Raster`. La calculatrice s'ouvre (:numref:`qgis-rast-calc-logique`). Dans le panneau ``Bandes raster`` nous trouvons la liste des rasters chargés dans QGIS. Nous remarquons qu'après le nom du raster nous trouvons un *@* suivi d'un numéro (*srtm_roya@1*). Le numéro correspond au numéro de bande dans le raster. Dans le cas d'un raster multi-bandes, nous aurions *mon_raster@1*, *mon_raster@2* et ainsi de suite. Dans le panneau ``Opérateur`` nous trouvons les opérateurs disponibles dans la calculatrice et dans le panneau ``Couche résultat`` les réglages pour le raster qui sera créé. Nous allons maintenant entrer l'expression qui va nous permettre d'extraire les altitudes supérieures à 2000 m dans le panneau ``Expression de la calculatrice raster``. Nous commençons par sélectionner le raster qui rentre dans le calcul en double-cliquant sur son nom dans le panneau ``Bandes raster``. Le raster apparaît alors dans le panneau d'expression. Nous entrons ensuite la suite de l'expression soit au clavier, soit en nous servant des opérateurs du panneau dédié. "srtm_roya@1" > 2000 À la ligne ``Couche en sortie``, nous spécifions un chemin et un nom vers le raster qui sera créé. Nous pouvons également régler le SCR du raster résultat à la ligne ``SCR en sortie``, généralement nous laissons le même SCR que le raster en entrée. Les réglages de l'emprise du résultat peuvent garder leurs valeurs par défaut (:numref:`qgis-rast-calc-logique`). Nous cliquons ensuite sur :guilabel:`OK`. .. figure:: figures/fen_qgis_raster_calc_logique.png :width: 30em :align: center :alt: alternate text :figclass: align-center :name: qgis-rast-calc-logique Extraction des altitudes supérieures à 2000 m du bassin de la Roya via la calculatrice raster de QGIS. Avec ce type de calcul logique, pour chaque pixel du raster nous regardons si son altitude est supérieure ou non à 2000. La réponse est *Oui* ou *Non*, c'est binaire. Par convention, en logique numérique, *Oui* est transcrit en *1* et *Non* est transcrit en *0*. Ainsi, sur le raster résultat, les pixels supérieurs à 2000 m seront codés en *1* et ceux inférieurs à 2000 m seront codés en *0*. Par défaut, dans la plupart des logiciels, les pixels à *1* sont colorés en *blanc* et les pixels à *0* le sont en *noir*. Nous obtenons ainsi le raster binaire suivant (:numref:`raster-alti` B). .. figure:: figures/fig_calc_raster_alti.png :width: 38em :align: center :alt: alternate text :figclass: align-center :name: raster-alti MNT initial (A), raster binaire des altitudes supérieures à 2000 m en blanc et inférieures en noir (B), altitudes supérieures à 2000 m en vert et inférieures en transparent (C). Il est tout à fait possible de complexifier ce genre d'extraction à l'aide d'opérateurs booléens comme le *ET* (*AND*) et le *OU* (*OR*) Ainsi, si nous souhaitons extraire à la fois les altitudes comprises entre 500 m et 700 m et celles comprises entre 1000 m et 1200 m, il nous suffira d'entrer l'expression suivante dans le panneau d'expression de la calculatrice raster. ("srtm_roya@1" > 500 AND "srtm_roya@1"< 700) OR ( "srtm_roya@1" > 1000 AND "srtm_roya@1" < 1200) Les *AND* servent à définir les bandes d'altitude et le *OR* sert à sélectionner les deux bandes. Nous obtenons bien un raster binaire correspondant à ces deux gammes d'altitudes (:numref:`raster-alti` B). .. figure:: figures/fig_calc_raster_alti2.png :width: 15em :align: center :alt: alternate text :figclass: align-center :name: raster-alti2 Raster binaire des deux bandes d'altitudes extraites. .. _raster-calc-multi-qgis: Calcul sur plusieurs rasters ++++++++++++++++++++++++++++++ Il est également possible d'effectuer des calculs combinant plusieurs rasters. Dans l'exemple suivant nous allons calculer un NDVI sur une image Landsat 8 prise au-dessus de Shangaï en août 2020. Après avoir chargé les deux rasters dans QGIS, nous ouvrons la calculatrice raster et nous entrons simplement la formule du NDVI avec les bonnes bandes (:ref:`data_Landsat_OLI`) dans le panneau d'expression en faisant bien attention aux parenthèses (:numref:`qgis-rast-calc-ndvi`). .. warning:: Les différents rasters qui entrent dans le calcul doivent être dans le même SCR. .. figure:: figures/fen_qgis_raster_calc_ndvi.png :width: 30em :align: center :alt: alternate text :figclass: align-center :name: qgis-rast-calc-ndvi Calcul de NDVI avec la calculatrice raster de QGIS. Nous obtenons bien le raster de NDVI de la zone (:numref:`raster-ndvi`). .. figure:: figures/fig_qgis_raster_calc_ndvi.png :width: 20em :align: center :alt: alternate text :figclass: align-center :name: raster-ndvi NDVI calculé pour la zone via la calculatrice raster de QGIS. .. warning:: Le numéro qui suit le *@* dans le nom des rasters dans la calculatrice raster correspond à la position de la bande dans un raster multi-bandes. Il ne correspond en aucun cas au numéro de la bande spectrale d'une image satellite. Par exemple, si vous travaillez avec un raster multi-bandes contenant les bandes *4* et *5* d'une image Landsat 8, ces bandes seront connues comme étant les bandes *1* et *2* du raster multi-bandes. .. _raster-calc-scp: Calcul raster avec SCP ************************ Version de QGIS : 3.20.1 Version de SCP : 7.9.5 Le module SCP (:ref:`logiciels-SCP`) possède une calculatrice raster similaire à celle de QGIS mais avec quelques fonctionnalités en plus. Il est par exemple possible d'utiliser des opérateurs statistiques et des clauses conditionnelles. .. _raster-calc-stats-scp: Opérateurs statistiques ++++++++++++++++++++++++++ Ici nous allons voir un exemple de calcul avec un opérateur statistique appliqué sur le MNT de la Roya. Pour l'exercice nous souhaitons savoir de combien de mètres s'écarte chaque pixel de l'altitude moyenne de la zone. Nous allons ainsi faire la différence entre la valeur de chaque pixel et la valeur moyenne du raster, ce qui est possible en une seule opération avec SCP. Nous allons dans le menu :menuselection:`SCP --> Calcul de bande` et la calculatrice de bande de SCP s'ouvre (:numref:`scp-calc-mean`). Il peut être nécessaire de rafraîchir la liste des bandes dans le panneau ``Liste de bandes`` en cliquant sur l'icône d'actualisation |icone_refresh|. Nous entrerons ensuite l'expression dans le panneau ``Expression`` en nous servant des fonctions proposées dans le panneau ``Fonctions``. Pour répondre à notre question nous entrons ainsi la formule suivante. "srtm_roya" - mean("srtm_roya") .. note:: Nous pouvons récupérer le nom du raster sans l'entrer à la main en double-cliquant sur son nom dans le panneau ``Liste de bandes``. De même, il suffit de double-cliquer sur la fonction *mean* dans le panneau ``Fonctions``. Nous évitons ainsi des fautes de frappe. .. figure:: figures/fen_SCP_calc_raster_mean.png :width: 40em :align: center :alt: alternate text :figclass: align-center :name: scp-calc-mean Calcul avec opérateur statistique dans SCP. Il suffit alors de cliquer sur :guilabel:`Lancer`. Il nous est demandé d'indiquer un chemin et un nom pour le raster résultat, puis le raster calculé apparaît (:numref:`scp-mean`). .. figure:: figures/fig_scp_raster_calc_mean.png :width: 15em :align: center :alt: alternate text :figclass: align-center :name: scp-mean Différence d'altitude entre chaque pixel et l'altitude moyenne de la région. Sur le raster résultat (:numref:`scp-mean`), les pixels allant vers le bleu sont les pixels dont l'altitude est inférieure à l'altitude moyenne de la région, et ceux tirant vers le rouge ceux dont l'altitude est supérieure à la moyenne. .. _raster-calc-scp-condition: Opérateurs logiques ++++++++++++++++++++++++++ Dans ce module SCP, nous trouvons également des opérateurs de conditions de type *Si telle condition est remplie, alors ceci, sinon cela*. Cet opérateur est identifié dans SCP sous le nom de *where*. La syntaxe est *where(condition, opération si Vrai, opération si Faux)*. Par exemple, si nous souhaitons recoder notre MNT en mettant à *10* les pixels dont l'altitude est supérieure à 1000 m et à *20* les pixels dont l'altitude est inférieure à 1000 m, nous entrons la formule suivante dans le panneau d'expression : *where("srtm_roya" > 1000 , 10, 20)* (:numref:`scp-calc-where`). .. figure:: figures/fen_SCP_calc_raster_where.png :width: 40em :align: center :alt: alternate text :figclass: align-center :name: scp-calc-where Calcul avec opérateur conditionnel dans SCP. Le raster résultat est similaire à un raster binaire sauf que, grâce à cet outil de condition, nous pouvons le coder en d'autres valeurs que *0* et *1*. Il est même possible d'emboîter les conditions afin d'obtenir plus de deux valeurs. .. _raster-calc-otb: Calcul raster avec OTB ************************ Version de OTB : 9.0.0 Le logiciel OrfeoToolbox (:ref:`logiciels-OTB`) propose également un module de calcul raster. Ce module est basé sur une librairie écrite en C++ nommée `muparser`_ optimisée pour le calcul matriciel. Ainsi, la calculatrice raster de OTB est performante en terme de puissance. L'autre avantage de OTB est de pouvoir être utilisé via QGIS ou via la console, ce qui est intéressant dans une optique d'automatisation. Dans les exemples suivants, nous calculerons un NDVI à partir d'une image Landsat 8. .. _raster-calc-otb-qgis: Avec OTB via QGIS +++++++++++++++++++ Le module permettant les calculs rasters avec OTB se nomme ``BandMath`` et se trouve dans la :menuselection:`Boîte à outils de traitements --> OTB --> Image Manipulation`. Une description détaillée est disponible sur `la documentation en ligne`_. Nous pouvons également chercher ce module via la barre de recherche de la ``Boîte à outils de traitements``. Une fois le module sélectionné, la fenêtre suivante s'affiche (:numref:`otb-calc-ndvi`). .. figure:: figures/fen_otb_bandmath_ndvi.png :width: 40em :align: center :alt: alternate text :figclass: align-center :name: otb-calc-ndvi Calcul raster avec l'application *BandMath* de OTB. À la ligne ``Input image list`` nous sélectionnons les rasters à utiliser dans le calcul. Le premier raster de la liste sera reconnu sous le nom de *im1*, le deuxième sous le nom de *im2* et ainsi de suite. Puis à la ligne ``Expression``, nous entrons l'expression à calculer en utilisant la syntaxe adéquate. Les rasters sont numérotés selon leur ordre dans la liste du panneau ``Input image list``. Ensuite, nous devons renseigner la bande à considérer dans chacun de ces rasters. Pour faire un calcul sur la première bande du premier raster nous indiquons *im1b1* et ainsi de suite. Dans notre cas nous avons des rasters mono-bande, la formule est donc : (im2b1 - im1b1) / (im2b1 + im1b1) où *im2b1* correspond à la première bande (*b1*) du deuxième raster chargé (*im2*), soit bien le raster de la bande du proche infrarouge *L8_Shangai_B5.tif* et *im1b1* correspond à la première bande (*b1*) du premier raster chargé (*im1*), soit bien le raster de la bande du rouge *L8_Shangai_B4.tif* Nous renseignons ensuite le chemin et le nom du raster résultat à la ligne ``Ouput image``. Nous cliquons enfin sur :guilabel:`Exécuter` pour lancer le calcul. Avec la commande otbcli +++++++++++++++++++++++++ En ligne de commande, l'opération est simple, il suffit de bien respecter la syntaxe et de bien se placer dans le répertoire qui contient nos rasters. Nous commençons par appeler le module *BandMath*, nous définissons ensuite les rasters en entrée via le paramètre *-il*, puis nous définissons le nom du raster résultat via le paramètre *-out* et enfin nous entrons la même expression que précédemment mais placée entre guillemets via le paramètre *-exp*. Dans notre cas, la ligne de commande est la suivante : .. code-block:: sh otbcli_BandMath -il L8_Shangai_B4.tif L8_Shangai_B5.tif -out ndvi_otbcli.tif -exp "(im2b1 - im1b1) / (im2b1 + im1b1)" Le point important est de bien faire attention à l'ordre des rasters. Le premier raster de la commande sera le *im1*, le deuxième *im2* et ainsi de suite. Dans tous les cas, le raster calculé est exporté en *.tif* et peut être ouvert dans QGIS ou tout autre logiciel de géomatique. Avec un opérateur de condition ++++++++++++++++++++++++++++++++ Ce module *BandMath* propose aussi un opérateur de condition, comme vu dans la partie dédiée au calcul raster avec SCP (:ref:`raster-calc-scp-condition`). La syntaxe est *(condition ? valeur si Vraie : valeur si Faux)*. Ainsi, dans notre exemple de coder les altitudes supérieures à 1000 m en 10 et celles inférieures à 1000 m en 20 sur le MNT de la Roya, l'expression à entrer à la ligne ``Expression`` (ou après *-exp*) est : (im1b1 > 1000 ? 10 : 20) où im1b1 correspont au raster de MNT. .. _raster-calc-R: Calcul raster avec R ************************ Version de R : 4.8.1 Lorsqu'un ou plusieur rasters sont chargés dans R, il est possible d'effectuer des opérations mathématiques dessus. Le calcul raster avec R est très simple et très utile. La façon de calculer est légèrement différente selon les packages utilisés. Ici, nous montrerons la marche à suivre pour un raster chargé sous forme d'un objet *SpatRaster* via *terra*. Nous travaillerons sur des bandes spectrales Landsat 8 prises au-dessus des Bouches-du-Rhône. Raster et scalaire ++++++++++++++++++++ Il est possible d'appliquer une opération sur un raster, i.e. sur chacun des pixels de ce raster. Dans notre exemple nous multiplierons un raster, de type *SpatRaster* stocké dans une variable nommé *proche_infrarouge*, par 1000. Concrètement, nous allons multiplier chacun des pixels constitutifs du raster *proche_infrarouge* par 1000. Nous obtiendrons ainsi en sortie un nouveau raster aux même propriétés géométriques que le raster initial mais avec des valeurs 1000 fois plus grandes pour chacun des pixels du raster. Ce raster sera aussi un objet *SpatRaster* que nous stockerons dans une variable nommée *pir_1000*. .. code-block:: R # chargement de la librairie terra library(terra) # chargement du raster pir <- terra::rast("Landsat_13/dep_13_L8_20190613_B5.tif") # multiplication de ce raster par 1000 pir_1000 <- pir * 1000 # affichage du résultat plot(pir_1000, col=gray.colors(100, start = 0.1, end = 0.9)) Nous constatons que l'opération est très simple, et nous obtenons bien notre raster initial où chacun des pixels ont été multipliés par 1000 (:numref:`calc-raster-R-1000`). .. figure:: figures/fig_R_calcul_raster_1000.png :width: 30em :align: center :alt: alternate text :figclass: align-center :name: calc-raster-R-1000 Résultat de la bande spectrale initiale multipliée par 1000 (remarquer le facteur 1000 dans la légende). Entre plusieurs rasters ++++++++++++++++++++++++ Il est également possible d'effectuer des opérations entre deux, ou plus de deux, rasters. Il est conseillé que les rasters en entrée aient les mêmes dimensions et soient dans les mêmes systèmes de coordonnées. Ici, nous calculerons un :ref:`ndvi` à partir de deux rasters de bandes spectrales du proche infrarouge (*pir*) et du rouge (*rouge*). Ces rasters sont au format *SpatRaster*. Nous stockerons le résultat dans une variable nommée *ndvi*. .. code-block:: R # chargement de la librairie terra library(terra) # chargement des rasters des bandes spectrales Rouge et PIR rouge <- terra::rast("Landsat_13/dep_13_L8_20190613_B4.tif") pir <- terra::rast("Landsat_13/dep_13_L8_20190613_B5.tif") # calcul du NDVI ndvi <- (pir - rouge) / (pir + rouge) La formulation est très simple. Le raster résultat sera également un objet *SpatRaster*. Opération logique ++++++++++++++++++ Dans certains cas, nous sommes amenés à appliquer des opérations logiques sur des rasters. À l'issue d'une telle opération logique nous obtenons un nouveau raster dans lequel les pixels répondant *Vrai* à la condition sont codés en *1* et les pixels répondant *Faux* sont codés à *0*. Le raster obtenu est dit *binaire* car exprimé en *0* et *1*. Ce type d'opération est courante lors de l'application de seuils. Dans l'exemple ici, nous appliquerons un seuil sur le MNT du `bassin de la Roya`_ afin d'extraire les pixels dont l'altitude est spérieure à 1500 m. .. code-block:: R # chargement de la librairie terra library(terra) # chargement du MNT mnt_roya <- terra::rast("srtm_roya_L93.tif") # application du seuil sur les altitudes mnt_roya_1500 <- mnt_roya > 1500 # affichage du résultat plot(mnt_roya_1500) Aucune difficulté dans la formulation, et le raster correspondant est bien un raster binaire dans lequel les pixels à *1* correspondent aux pixels ayant une altitude supérieure à 1500 m (:numref:`calc-raster-R-seuil`). .. figure:: figures/fig_R_calcul_raster_seuil.png :width: 30em :align: center :alt: alternate text :figclass: align-center :name: calc-raster-R-seuil Raster binaire après un calcul raster de type seuil dans R. .. _raster-calc-gdal: Calcul raster avec GDAL ************************* Version de GDAL : 3.0.4 GDAL a une fonctionnalité permettant de faire du calcul raster. Comme GDAL s'exécute en lignes de commandes, son usage peut être intéressant pour faire un calcul rapide ou bien un calcul raster dans un script. La commande dédiée est *gdal_calc.py*. C'est une extension Python de GDAl installée par défaut avec GDAL. Les détails de la commande et la liste des options possibles se trouvent sur la `page du manuel dédiée`_. Nous verrons ici les fonctionnalités les plus fréquentes. Raster et scalaire +++++++++++++++++++ Il est très simple d'associer un raster à un scalaire par un calcul. Par exemple, si nous souhaitons multiplier les altitudes d'un MNT, nommé *srtm_Roya.tif*, par 10, la commande à entrer est la suivante. .. code-block:: sh gdal_calc.py -A srtm_roya.tif --outfile=mnt_10.tif --calc="A*10" Avec : * *gdal_calc.py* : appel à la commande GDAL * *-A srtm_roya.tif* : le raster à prendre en entrée qu'on stocke dans la variable *A* * *--outfile=mnt_10.tif* : le nom du raster résultat en sortie * *--calc="A*10"* : l'expression à calculer. Bien noter que cette expression fait appel à la variable *A* et pas directement au raster et que cette expression est entre guillemets et sans espaces En sortie, nous obtenons bien un raster au format GeoTiff que nous pouvons charger dans QGIS. Sur ce raster, les pixels correspondent aux pixels du MNT de départ mais multipliés par 10. Raster et calcul logique ++++++++++++++++++++++++++ Sur le même principe, il est facile d'appliquer une opération logique sur un raster. Par exemple, si nous souhaitons extraire les pixels dont l'altitude est comprise entre 1000 m et 2000 m à partir d'un MNT nommé *srtm_Roya.tif*, la commande à entrer est la suivante. .. code-block:: sh gdal_calc.py -A srtm_roya.tif --outfile=alti_1000-2000.tif --calc="(A>1000)&(A<2000)" Notez que l'opérateur à employer pour le *ET* est l'esperluette *&*. Le résultat est un raster binaire où les pixels à *1* correspondent aux pixels respectant la condition et les pixels à *0* correspondent à tous les autres pixels. Cette expression peut également prendre la forme suivante. .. code-block:: sh gdal_calc.py -A srtm_roya.tif --outfile=alti_1000-2000.tif --calc="logical_and(A>1000, A<2000)" Ici nous employons le mot clef *logical_and()* avec la condition à l'intérieur des parenthèses. Si le calcul fait appel à une condition *OU*, il faut alors employer le mot clef *logical_or()*. Par exemple, si nous disposons d'un raster d'occupation du sol de type Corine Land Cover et que nous souhaitons extraire les occupations du sol de classe *211* et *511*, la commande GDAL sera la suivante. .. code-block:: sh gdal_calc.py -A CLC.tif --outfile=CLC_extract.tif --calc="logical_or(A==211,A==511)" Notez l'emploi de deux signes *=* pour l'égalité. Nous obtenons bien un nouveau raster binaire avec les pixels à *1* qui correspondent aux pixels satisfaisant la condition. Calcul sur plusieurs rasters ++++++++++++++++++++++++++++++ Il est également tout à fait possible de combiner plusieurs rasters lors d'un calcul. Les rasters d'entrée doivent par contre être de mêmes dimensions et dans le même système de projection. Par exemple, pour calculer un :ref:`NDVI` à partir d'une bande de proche infrarouge nommée *L8_B5.TIF* et d'une bande rouge nommée *L8_B4.TIF*, nous entrons la commande suivante. .. code-block:: sh gdal_calc.py -A L8_B4.TIF -B L8_B5.TIF --outfile=ndvi.tif --calc="(B-A)/(B+A)" La différence avec précédemment est que nous spécifions deux rasters en entrée stockés respectivement dans les variables *A* et *B*. Ce sont ces deux variables qui sont utilisées dans l'expression de calcul que nous retrouvons dans *--calc="(B-A)/(B+A)"*. En sortie, nous obtenons bien notre raster de NDVI. .. warning:: Le raster de sortie d'un calcul raster avec GDAL est du même type que les rasters en entrée. Ainsi, si les rasters en entrée sont de type *Integer* (Entiers), le résultat sera aussi de type entier, ce qui peut amener à des surprises inattendues et non souhaitées. Il est censé être possible de forcer le type souhaité en spécifiant par exemple *A.astype(numpy.float64)* dans le calcul mais ça n'a pas l'air de toujours fonctionner... Si les rasters à combiner dans le calcul ne sont pas constitués de fichiers séparés mais proviennent de bandes différentes se trouvant au sein d'un raster multi-bandes, la syntaxe doit être modifiée. Par exemple, ici nous disposons d'un raster multi-bandes contenant des bandes spectrales Landsat 8 de 1 à 7. La bande du *rouge* est donc la bande 4 et celle du *proche infrarouge* la bande 5. Pour calculer le NDVI à partir de ce raster multi-bandes, la syntaxe est la suivante. .. code-block:: sh gdal_calc.py -A L8_stack.tif --A_band=5 -B L8_stack.tif --B_band=4 --outfile=ndvi_stack.tif --calc="(A-B)/(A+B)" Notez que nous stockons le raster multi-bandes une première fois dans la variable *A* et que nous spécifions que nous y prenons simplement la bande *5* via le mot clef *--A_band=5*. Puis nous stockons le raster multi-bandes une seconde fois dans la variable *B* et nous spécifions que nous y prenons cette fois-ci la bande *4* via le mot clef *--B_band=4*. Ainsi les *A* et *B* appelés dans l'expression du calcul correspondent bien aux bandes désirées. .. _raster-calc-stats: Calcul raster avec statistiques -------------------------------- En télédétection il est fréquent de manipuler des rasters multi-bandes. Dans la plupart des cas, ces rasters multi-bandes sont des regroupements de réflectances dans différentes bandes spectrales pour une date donnée ou bien des regroupements d'une même réflectance ou d'un indice pour plusieurs dates. Dans ce second cas nous parlons de *séries temporelles*. Pour effectuer des calculs, ou des statistiques pixel à pixel, sur ces rasters multi-bandes, il n'est pas possible d'utiliser les outils basiques de calculs rasters. Il est nécessaire d'employer des outils dédiés. Dans tous les cas, le raster résultat sera géographiquement identique aux rasters en entrée mais les valeurs des pixels auront pris la valeur de la statistique utilisée en fonction des rasters en entrée, comme présenté sur la figure suivante (:numref:`raster_stats_calc`). .. figure:: figures/fig_raster_calc_stats.png :width: 30em :align: center :alt: alternate text :figclass: align-center :name: raster_stats_calc Calcul de statistiques sur un stack de rasters. Dans cet exemple théorique nous avons en entrée 4 rasters de 9 pixels chacun. Dans un premier temps, nous calculons le raster issu de la statistique *Maximum* appliquée au stack de ces 4 rasters. Nous obtenons un nouveau raster dans lequel chaque pixel prend la valeur maximale rencontrée dans les 4 pixels d'entrée. Dans un second temps, nous calculons le raster issu de la statistique *Moyenne* appliquée au stack des 4 rasters. Nous obtenons un nouveau raster dans lequel chaque pixel prend la valeur moyenne des 4 pixels d'entrée. .. _raster-calc-stats-avec-scp: Calcul raster et stats avec QGIS et SCP ******************************************** Version de QGIS : 3.22.0 Version de SCP : 7.10.5 Le plugin SCP de QGIS propose une fonctionnalité intéressante pour calculer des statistiques pixel à pixel sur un raster multi-bandes. L'idée est de créer un nouveau raster sur lequel les valeurs des pixels correspondront à un résumé statistique pour chacun des pixels. Dans l'exemple ci-dessous, nous allons calculer la réflectance moyenne dans le proche-infrarouge de 4 images Sentinel-2 prises sur la même zone mais à 4 dates différentes. Nous commençons par charger dans QGIS ces 4 rasters. Puis nous allons dans le menu :menuselection:`SCP --> Jeu de bandes` pour créer un jeu de bandes qui contiendra nos 4 rasters de réflectance. Une fois le jeu de bandes (*Band set*) défini, nous allons dans le menu :menuselection:`SCP --> Calcul de bandes`. Le menu suivant s'affiche (:numref:`multi-bandes-scp`). .. figure:: figures/fen_moy_multi_bandes_scp.png :width: 35em :align: center :alt: alternate text :figclass: align-center :name: multi-bandes-scp Moyenne pixel à pixel sur un raster multi-bandes avec SCP dans QGIS. Dans le panneau ``Expression`` nous entrons l'expression suivante : ``mean("bandset#b*")``. Par défaut, le *bandset* est le Jeu de bandes actif. Cette syntaxe signifie que nous travaillons sur le *bandset* actif et sur toutes les bandes le composant *#b**. Et sur toutes ces bandes nous demandons la moyenne *mean()*. Nous cliquons ensuite sur :guilabel:`Exécuter` et nous spécifions un chemin et un nom pour le raster résultat. Chaque pixel de ce raster résultat correspondra à la moyenne des réflectances sur les 4 dates. .. _raster-calc-stats-otb: Calcul raster et stats avec OTB ********************************** Version de OTB : 9.0.0 OrfeoToolbox propose une fonctionnalité qui permet de faire du calcul multibandes combinant des statistiques. Dans l'exemple ci-après, nous allons combiner 4 rasters de NDVI calculés sur une même région mais à 4 dates différentes au sein d'une année. Nous allons construire un nouveau raster, qui sera géographiquement identique aux rasters de NDVI initiaux, mais dont les pixels correspondront au NDVI maximum des 4 dates. Avec OTB via QGIS ++++++++++++++++++ Le module à utiliser se nomme ``BandMathX``, nous le trouvons dans la :menuselection:`Boîte à outils de traitements --> OTB --> Image Manipulation --> BandMathX`. La fenêtre suivante s'ouvre (:numref:`raster-calc-stats-monteverdi`). .. figure:: figures/fen_otb_BandMathX.png :width: 35em :align: center :alt: alternate text :figclass: align-center :name: raster-calc-stats-monteverdi Calcul raster avec statistiques sur plusieurs bandes. Dans le panneau ``Input image list`` nous renseignons les rasters sur lesquels nous souhaitons calculer des statistiques. Ces rasters doivent avoir exactement les mêmes dimensions et les mêmes systèmes de projection. Ici nous ajoutons nos 4 NDVI calculés pour nos 4 dates de l'année 2016. À la ligne ``Output Image``, nous spécifions un chemin et un nom pour le raster de sortie. Enfin, le plus important est l'expression à entrer à la ligne ``Expressions``. Dans le cas du calcul d'un maximum, l'expression est : .. code-block:: sh max(im1b1,im2b1,im3b1,im4b1) où *max* est la statistique à calculer et *im1b1* signifie que nous prenons la *bande 1* (*b1*) de *l'image 1* (*im1*) correspondant au raster chargé en premier, et ainsi de suite. .. note:: Dans le cas où nous ne travaillons pas sur des rasters séparés mais sur des bandes au sein d'un raster multi-bandes, l'expression devient simplement : max(im1b1,im1b2,im1b3,im1b4). Nous ne travaillons qu'avec l'image 1 *im1* mais nous y appelons les différentes bandes de *b1* à *b4*. Une fois cliqué sur :guilabel:`Exécuter`, le raster contenant la valeur maximums des NDVI est calculé. Il est possible de le charger dans QGIS. Avec la commande otbcli ++++++++++++++++++++++++ Dans une optique d'automatisation de tâches, il est possible d'appeler ce module *BandMathX* en ligne de commandes *otbcli*. Cette ligne s'écrit alors simplement comme ci-dessous. .. code-block:: sh otbcli_BandMathX -il ndvi_2016-01-25.tif ndvi_2016-03-15.tif ndvi_2016-05-04.tif ndvi_2016-08-22.tif -out max_ndvi2.tif -exp "max(im1b1,im2b1,im3b1,im4b1)" Où : * *otbcli_BandMathX* est l'appel au module *BandMathX* * *-il* est le mot clef pour spécifier les rasters en entrée * *-out* est le mot clef pour spécifier le nom du raster en sortie * *-exp* est le mot clef pour entrer l'expression au même format que vu dans le paragraphe précédent .. note:: Il peut être laborieux d'entrer le nom à la main de tous les rasters en entrée. Le plus simple est de créer au préalable un :ref:`rasters-multi-bandes` puis d'appeler les différentes bandes de ce raster dans l'expression du calcul. .. _reclassifier-raster: Reclassifier un raster ----------------------- Il est souvent utile de reclassifier un raster, c'est-à-dire d'assigner de nouvelles valeurs aux pixels en fonction de règles. Par exemple, il est possible de changer tous les pixels dont les valeurs sont comprises entre une valeur *a* et une valeur *b* par une valeur *c*. Il est également possible de changer tous les pixels qui ont exactement la valeur *a* par la valeur *c*. Dans les exemples qui suivent, nous reclassifierons un raster d'usage du sol du Caire en 7 classes, obtenu par classification supervisée, pour obtenir un nouveau raster présentant de nouveaux identifiants de classe. Les 7 classes initiales sont rappelées dans le tableau suivant. .. list-table:: Occupation du sol initiale :widths: 15 25 :header-rows: 1 * - ID classe - Label classe * - 1 - Eau * - 2 - Végétation inondée * - 3 - Végétation dense * - 4 - Sol nu sableux * - 5 - Sol nu rocheux * - 6 - Bâti dense * - 7 - Bâti diffus La reclassification nous permettra d'assigner de nouveaux numéros de classes aux différentes classes ou de regrouper plusieurs classes en une seule classe. Pour l'exemple, dans un premier temps nous assignerons l'identifiant 11 à la classe 1, 22 à la classe 2 et ainsi de suite jusqu'à 77 à la classe 7. Dans un second temps, nous regrouperons les classes par grandes catégories comme indiqué dans le tableau suivant. .. list-table:: Reclassification (regroupement) de l'occupation du sol :widths: 15 25 15 25 :header-rows: 1 * - ID classe old - Label old - ID classe new - Label new * - 1 - Eau - 1 - Eau * - 2 - Végétation inondée - 2 - Végétation * - 3 - Végétation dense - 2 - Végétation * - 4 - Sol nu sableux - 3 - Sol nu * - 5 - Sol nu rocheux - 3 - Sol nu * - 6 - Bâti dense - 4 - Bâti * - 7 - Bâti diffus - 4 - Bâti .. _reclassifier-raster-qgis: Reclassifier un raster avec QGIS ********************************* Version de QGIS : 3.20.3 Une fois le raster à reclassifier chargé, nous ouvrons le menu ``Reclassification par table`` dans le panneau :menuselection:`Boîte à outils de traitements --> Analyse raster --> Reclassification par table`. Le menu s'ouvre (:numref:`qgis-reclass`). .. figure:: figures/fen_qgis_reclass.png :width: 25em :align: center :alt: alternate text :figclass: align-center :name: qgis-reclass Reclassification de raster avec QGIS. **Reclassification classe par classe** Nous commençons par reclassifier le raster classe par classe, sans regroupement. Nous changeons simplement l'identifiant de chaque classe par le nouvel identifiant (1 par 11, 2 par 22...). À la ligne ``Couche raster``, nous sélectionnons le raster à reclassifier, ici ``classif_7``. Dans le cas d'un raster multi-bandes, nous précisons la bande à reclassifier à la ligne ``Numéro de bande``. Ensuite, le travail de reclassification s'effectue via une *table de reclassification* à paramétrer à la ligne ``Table de reclassification`` en cliquant sur l'icône |icone_browse| en bout de ligne. Le menu de paramétrage de la table de reclassification s'ouvre (:numref:`qgis-reclass-table1`). Nous commençons par ajouter autant de lignes que de classes à reclassifier en cliquant autant de fois que nécessaire sur ``Ajouter une ligne``. Ensuite, nous indiquons la borne inférieure et supérieure des gammes de valeurs à changer dans les colonnes ``Minimum`` et ``Maximum`` et nous renseignons la nouvelle valeur dans la colonne ``Valeur``. Ici, comme nous changeons les valeurs de toutes les classes, classe par classe, les bornes minimum et maximum sont les mêmes et correspondent aux classes à changer (:numref:`qgis-reclass-table1`). .. figure:: figures/fen_qgis_reclass_table1.png :width: 25em :align: center :alt: alternate text :figclass: align-center :name: qgis-reclass-table1 Table de reclassification de raster classe par classe. Une fois ce paramétrage effectué, nous cliquons sur ``OK`` pour revenir à la fenêtre principale du module. Il faut maintenant régler la règle de reclassification en dépliant le panneau ``Paramètres avancés``. Dans ce panneau, à la ligne ``Limites de plages``, nous spécifions la règle à suivre pour la reclassification. Ici, les valeurs maximum et minimum correspondent directement à nos classes. Nous choisissons donc la règle *min <= valeur <= max*. Nous renseignons ensuite le chemin d'export du raster qui sera créé à la ligne ``Raster reclassifié``. Puis nous cliquons sur ``Exécuter``. Le raster reclassifié apparaît automatiquement. Les classes ont bien été réattribuées. **Reclassification pour regroupement de classes** Dans ce cas de figure, la table de reclassification doit être paramétrée différemment. Les bornes minimum et maximum doivent être renseignées de telle sorte que ce soient des gammes de valeurs qui vont être réassignées. Ici, pour opérer un regroupement de classes tel que présenté dans l'introduction de cette section (:ref:`reclassifier-raster`), nous paramétrons la table tel que présenté ci-dessous (:numref:`qgis-reclass-table2`). .. figure:: figures/fen_qgis_reclass_table2.png :width: 25em :align: center :alt: alternate text :figclass: align-center :name: qgis-reclass-table2 Table de reclassification de raster pour regroupement de classes. Puis nous réglons la règle de reclassification dans le panneau ``Paramètres avancés``. Tel que nous avons paramétré nos bornes de classes, la règle appropriée est alors *min < valeur <= max*. Nous renseignons le chemin du raster de regroupement de classes et nous cliquons sur ``Exécuter``. Le raster aux classes regroupées apparaît (:numref:`usol-caire-reclass`). .. figure:: figures/fig_usol_caire_reclass.png :width: 35em :align: center :alt: alternate text :figclass: align-center :name: usol-caire-reclass Occupation du sol initiale en 7 classes (A) et occupation du sol après regroupement en 4 classes (B). .. _reclassifier-raster-R: Reclassifier un raster avec R ********************************* Version de R : 4.3.1 Version de terra : 1.7.29 Nous allons reproduire les manipulations de reclassification présentées kdans l'introduction de cette partie (:ref:`reclassifier-raster`) avec R en utilisant la librairie *terra*. **Reclassification classe par classe** Nous commencerons par reclassifier un raster d'occupation du sol classe par classe. À chaque classe initiale nous allons attribuer une nouveau numéro de classe. Cela se fait avec la méthode *classify* de *terra*, comme présenté ci-après. .. code-block:: R # import de l'occupation du sol à reclasser usol <- terra::rast('classif_C_2011.tif') # reclass classe par classe # on créé une matrice de reclassification m <- c(1, 11, 2, 22, 3, 33, 4, 44, 5, 55, 6, 66, 7, 77) rclmat <- matrix(m, ncol=2, byrow=TRUE) # on applique cette matrice de reclassification usol_reclass <- terra::classify(usol, rclmat) Comme nous pouvons le voir, nous commençons par créer une matrice de reclassification de type *ancienne valeur, nouvelle valeur* puis nous appliquons cette matrice au raster à reclassifier. **Reclassification pour regroupement de classes** Il est également possible de reclassifier en regroupant des classes. Dans cet exemple, nous allons simplifier notre occupation du sol initialement en 7 classes en une occupation du sol en 4 classes (:numref:`qgis-reclass-table2`). Nous utilisons également la fonction *classify*. .. code-block:: R # import de l'occupation du sol à reclasser usol <- terra::rast('classif_C_2011.tif') # reclass classe par classe # on créé la matrice de reclassification m <- c(1, 1, 1, 1, 3, 2, 3, 5, 3, 5, 7, 4) rclmat <- matrix(m, ncol=3, byrow=TRUE) # on applique la matrice de reclassification usol_reclass <- terra::classify(usol, rclmat) Cette fois-ci la matrice s'écrit sur 3 colonnes de la façon *de telle valeur, à telle valeur, nouvelle valeur*. Nous pouvons toutefois noter que par défaut, la fonction *classify* ne prend pas la valeur minimum de chaque intervalle mais prend la valeur maximale. C'est pourquoi nous lui indiquons que les pixels allant de *1* (non inclus) à *3* (inclus) doivent prendre la nouvelle valeur de *2*. Il est également possible de reclassifier spécifiquement certaines valeurs et de mettre toutes les autres valeurs à telle autre valeur. Par exemple, si nous souhaitons reclassifier notre occupation du sol en deux classes *Végéation* et *Pas végétation* en mettant les deux classes de végétation à *1* et tout le reste à *0*, nous pouvons utiliser le mot clef *others*, comme présenté ci-dessous. .. code-block:: R # import de l'occupation du sol à reclasser usol <- terra::rast('classif_C_2011.tif') # reclass en deux classes # on créé la matrice de reclassification simplement pour la classe d'intérêt m <- c(1, 3, 1) rclmat <- matrix(m, ncol=3, byrow=TRUE) # on applique la matrice de reclassification en spécifiant "tout le reste à 0" usol_reclass <- terra::classify(usol, rclmat, others=0) .. note:: La fonction *classify* de *terra* présente d'autres options intéressantes pour des usages plus précis, à `regarder de plus près`_ en cas de besoin. .. _compresser-raster: Compresser un raster ----------------------- Il arrive d'avoir à faire à des rasters qui deviennent trop volumineux pour un ordinateur donné ou pour un logiciel. Une façon de faire peut être de le découper en tuiles, comme expliqué dans la rubrique dédiée (:ref:`decouper-raster-grille`). Mais il est également possible de réduite le poids d'un raster en utilisant des algorithmes de compression. Une présentation détaillée de ces `algorithmes peut être trouvée en ligne`_. Ici, nous nous contenterons d'en utiliser un sans regarder les processus sous-jacents. Néanmoins, le point clé à bien retenir est que cette compression n'altère pas les données. Ainsi, les valeurs des pixels sont bien conservées ce qui permet de garder l'intégrité de l'information du raster. Dans les exemples ci après nous verrons comment compresser un raster d'une bande spectrale Landsat 8. .. _compresser-raster-qgis: Compresser un raster avec QGIS ********************************* Version de QGIS : 3.26.0 Tout d'abord nous chargeons le raster à compresser dans QGIS. Une fois le raster chargé, il suffit de faire un clic droit sur la couche dans le ``panneau des couches``. Puis nous sélectionnons :menuselection:`Exporter --> Enregistrer sous...` Le menu suivant s'affiche (:numref:`compress-raster-qgis`). .. figure:: figures/fen_qgis_compresser_raster.png :width: 35em :align: center :alt: alternate text :figclass: align-center :name: compress-raster-qgis Compression d'un raster avec QGIS. Dans le champ ``Format`` nous laissons ``GeoTIFF``. Dans le champ ``Nom de fichier``, nous précisons un chemin et un nom pour le raster compressé qui sera créé. Nous pouvons laisser le champ ``SCR`` à sa valeur pour conserver le même SCR. Les options de compression s'obtiennent en dépliant le menu ``Options de création``. Nous choisissons ``Compression élevée`` dans le menu ``Profil`` et nous pouvons laisser les valeurs proposées par défaut pour les variables ``Compress``, ``Predictor`` et ``Zlevel`` qui produisent de bons résultats. Pour plus de détails techniques et des comparaisons de différentes compressions vous pouvez vous référer `à cette page`_. Il ne reste qu'à cliquer sur :guilabel:`OK`. Le raster passe ainsi d'un poids de 228.6 MiO à 85.2 MiO, soit une diminution de poids de près de 60 %, ce qui n'est pas négligeable du tout. De plus, comme dit en préambule, les valeurs des pixels sont conservées à l'identique. .. _compresser-raster-gdal: Compresser un raster avec GDAL ********************************* Version de GDAL : 3.0.4 Il peut être intéressant de passer par GDAL en ligne de commandes pour compresser un raster. Ce type d'usage permet d'aller vite ou de chaîner des traitements facilement. Pour cela, il suffit que GDAL soit bien renseigné dans le *PATH* de votre système d'exploitation. La ligne de commande à utiliser est la suivante. .. code-block:: sh gdal_translate -of GTiff -co "COMPRESS=DEFLATE" -co "PREDICTOR=2" -co "ZLEVEL=9" LC08_L2SP_196030_20190613_20200828_02_T1_SR_B5.TIF gdal-deflate_2-9.tiff Cette ligne de commande se comprend de la façon suivante : * *gdal_translate* : la fonction GDAL dédié à la conversion de format des rasters * *-of GTiff* : spécification du format de raster en sortie, ici *GeoTIFF* * *-co "COMPRESS=DEFLATE"* : l'option qui stipule qu'on va compresser le raster avec l'option *DEFLATE*. D'autres `options sont possibles`_. * *-co "PREDICTOR=2" -co "ZLEVEL=9"* : les deux options qui précisent le niveau de compression *DEFLATE*. Il est possible de changer ces valeurs mais celles-ci donnent de bons résultats. * *LC08_L2SP_196030_20190613_20200828_02_T1_SR_B5.TIF* : le raster à compresser * *gdal-deflate_2-9.tiff* : le nom du raster compressé qui sera créé Le raster passe ainsi d'un poids de 228.6 MiO à 85.2 MiO. .. _stats-zonales-raster: Statistiques zonales raster ----------------------------- Il est possible de calculer des statistiques zonales d'un raster en fonction d'un raster classifié. C'est utile lorsque nous disposons d'un raster avec des valeurs continues (des altitudes, des réflectances, des températures, ...), qui peut être multi-bandes, et que nous souhaitons avoir des statistiques selon des zones définies sur un second raster en valeurs discrètes. Ces valeurs discrètes peuvent être des labels de classes, aussi appelées *zones*. Sauf cas contraire, lans les exemples ci-après nous montrerons comment calculer des signatures spectrales moyennes d'une image Landsat TM en six bandes prise au-dessus du Caire le 13 juillet 2011 pour 10 classes définies par classification non supervisée (:ref:`classification-non-supervisee`). Le raster multi-bandes contient ainsi les réflectances de la scènes dans six bandes spectrales et le raster labellisé contient 10 valeurs différentes, correspondant à 10 classes. .. warning:: Pour la plupart des outils permettant ce calcul de statistiques zonales, il est nécessaire d'avoir les deux rasters dans le même système de coordonnées de référence. Si nécessaire, penser à reprojeter l'un des deux rasters (:ref:`reprojection-raster`). .. _stats-zonales-rasters-QGIS: Statistiques zonales raster avec QGIS *************************************** Version de QGIS : 3.38.1 Pour cet exemple, nous calculerons l'altitude moyenne sous différentes occupations du sol sur le bassin-versant de la Roya (Alpes-Maritimes). Le MNT utilisé est le SRTM et l'occupation du sol est le raster Corine Land Cover de 2018. .. tip:: QGIS permet de prendre en entrée deux rasters dans des systèmes de coordonnées différents et/ou de résolutions différentes. Le calcul de statistiques zonales entre deux rasters se fait facilement dans QGIS. Une fois le raster *continu* et le raster labellisé chargés, il suffit d'aller dans :menuselection:`Boîte à outils de traitements --> Analyse raster --> Statistiques zonales de la couche raster`. Le menu suivant s'affiche (:numref:`stats-zonales-raster-qgis`). .. figure:: figures/fen_qgis_stats_zonales_rasters.png :width: 35em :align: center :alt: alternate text :figclass: align-center :name: stats-zonales-raster-qgis Statistiques zonales entre deux rasters avec QGIS. Dans le champ ``Couche source`` nous renseignons le raster sur lequel calculer les statistiques zonales, ici le modèle numérique de terrain *srtm_roya_L93*. À la ligne, ``Couche zonale`` nous spécifions le raster contenant les labels, ici l'occupation du sol *CLC_2018_Roya*. Enfin, à la ligne ``Statistiques`` nous renseignons un chemin et un nom pour le fichier qui contiendra les statistiques. Il s'agira d'un GeoPackage sans géométrie. .. warning:: Comme il est possible de le constater dans le menu précédent, QGIS ne propose pas de calculer des statistiques zonales sur un raster multi-bandes. Si besoin, privilégier un des outils présentés plus lion dans cette rubrique. Le résultat apparaît dans le ``panneau des couches`` de QGIS et peut s'ouvrir comme une table attributaire (:numref:`tab-stats-zonales`). .. figure:: figures/fig_qgis_tab_stats_zonales_rasters.png :width: 35em :align: center :alt: alternate text :figclass: align-center :name: tab-stats-zonales Table des statistiques zonales entre deux rasters. Dans le champ ``Zone`` nous retrouvons les numéros des zones, qui correspondent ici aux classes d'occupation du sol Corine Land Cover, puis les différentes statistiques d'altitudes dans les colonnes suivantes. Ici, seule la colonne ``mean`` (*moyenne*) a du sens. Par exemple, l'occupation du sol numéro *23* a une altitude moyenne de 903.4 m. .. _stats-zonales-raster-OTB: Statistiques zonales raster avec OTB ********************************************* Version de OTB : 9.0.0 OTB propose un outil pour calculer des statistiques zonales entre deux rasters, en permettant d'utiliser un raster multi-bandes. Par contre, les deux rasters doivent être dans le même système de coordonnées de référence et de la même résolution spatiale. .. _stats-zonales-raster-OTB-QGIS: Avec OTB via QGIS +++++++++++++++++++++++++++++ Dans cet exemple, nous allons calculer les signatures spectrales moyennes d'un raster multi-bandes suite à une classification non supervisée. Une fois les deux rasters chargés dans QGIS, le menu OTB dédié se trouve dans la :menuselection:`Boîte à outils de traitements --> OTB --> Image manipulation --> ZonalStatistics`. Le menu suivant s'affiche (:numref:`otb-stats-zonales_rasters`). .. figure:: figures/fen_otb_stats_zonales_rasters.png :width: 35em :align: center :alt: alternate text :figclass: align-center :name: otb-stats-zonales_rasters Statistiques zonales entre deux rasters avec OTB. À la ligne ``Input image`` nous sélectionnons le raster sur lequel calculer les statistiques. Ce raster peut être multi-bandes. Ici, nous choisissons le stack multi-bandes de l'image Landsat du Caire. À la ligne ``Type of input for the zone definitions`` nous choisissons ``labelimage`` et nous choisissons le raster contenant les zones (i.e. les classes) à la ligne ``Input label image``. Dans le champ ``Format of the output stats`` nous pouvons prendre ``vector`` et nous stipulons où sauvegarder le fichier vecteur qui sera produit à la ligne ``Filename for the output vector data``. Le fichier produit sera une vectorisation du raster de classes avec une table attributaire associée qui contiendra les statistiques zonales pour chaque bande spectrale du raster initial. Il est également possible de sélectionner ``xml`` à la ligne ``Format of the output stats``. Auquel cas il faut écrire un chemin où stocker le fichier xml produit à la ligne ``Filename for the output xml data``. Un fichier XML contenant les statistiques zonales pour chaque classe et chaque bande du raster initial sera produit et lisible avec un éditeur de texte. .. note:: L'ergonomie des sorties fournies par OTB n'est pas très pratique... .. _aligner-rasters: Aligner deux rasters ----------------------- Pour certains traitements, comme le calcul entre plusieurs rasters, il est nécessaire d'aligner des rasters lorsque ceux-ci ne sont pas de la même résolution et donc pas alignés sur une même grille. C'est par exemple le cas lorsque nous calculons certains indices radiométriques faisant appel à des bandes spectrales de différentes résolutions spatiales. Il faudra alors choisir un raster qui servira de référence et aligner le, ou les autres rasters, sur ce raster de référence. Il existe différents outils pour aligner un raster selon un autre raster. .. warning:: Certains logiciels, comme QGIS, permettent d'emblée de faire des calculs rasters entre rasters de différentes résolutions. Dans ce cas, le logiciel réalise par lui même, de façon cachée, un alignement entre les différents rasters. L'avantage est que l'utilisateur n'a pas à se soucier de cette étape, mais le désavantage est que l'utilisateur n'a pas la main sur l'alignement. Quel raster utilisé comme référence ? Quelle méthode de rééchantillonnage pour l'alignement ? Le logiciel décide sans rien dire... .. warning:: Un alignement entre deux rasters entraîne obligatoirement un ré-échantillonnage de l'un des rasters. Or la méthode de ré-échantillonnage a son importance. Si nous souhaitons aligner un raster *discret* (une occupation du sol par exemple) il faut choisir le plus proche voisin, sinon une interpolation linéaire entre les valeurs de classes sera faite, ce qui n'a pas de sens. Par contre, sur un raster *continu* (comme un MNT ou une bande spectrale) il est possible de choisir une méthode qui interpole. .. _aligner-rasters-QGIS: Aligner deux rasters avec QGIS ******************************** Version de QGIS : 3.22.1 Même si QGIS permet de faire des calculs entre rasters de différentes résolutions, il est possible d'aligner des rasters manuellement en choisissant le raster de référence et la méthode de ré-échantillonnage. Pour l'exemple nous allons aligner un raster correspondant à la bande 3 à 10 m d'une image Sentinel-2 sur un raster correspondant à la bande 11 à 20 m. Une fois les deux rasters chargés dans QGIS, nous allons dans le menu :menuselection:`Raster --> Aligner les rasters...`. Dans le menu qui s'ouvre, nous commençons par choisir les rasters à aligner entre eux en cliquant sur le |icone_plus|. La fenêtre suivante s'affiche (:numref:`aligner-rasters-qgis-1`). .. figure:: figures/fen_qgis_aligner_rasters_1.png :width: 35em :align: center :alt: alternate text :figclass: align-center :name: aligner-rasters-qgis-1 Sélection du premier raster à aligner dans QGIS. À ``Couche d'entrée`` nous choisissons au moins les deux rasters à aligner. Nous pouvons commencer par le raster à aligner, à savoir celui à 10 m nommé *S2_B03_10m*. À la ligne ``Nom de du fichier de raster de sortie`` nous indiquons un chemin et un nom pour le raster qui sera aligné, par exemple *S2_B03_10m_resample.tif*. Puis nous devons choisir une ``Méthode de ré-échantillonnage`` parmi celles proposées dans le menu déroulant. Nous pouvons choisir ``Bilinéaire`` dans notre cas. Enfin nous cliquons sur :guilabel:`OK`. Nous recliquons sur le |icone_plus| pour ajouter le second raster, ici *S2_B11_20m*. Nous indiquons ce raster à la ligne ``Nom de du fichier de raster de sortie`` et nous choisissons aussi un nom de sortie et une méthode de ré-échantillonnage pour ce raster. La fenêtre précédente s'est mise à jour comme suit (:numref:`aligner-rasters-qgis-2`). .. figure:: figures/fen_qgis_aligner_rasters_2.png :width: 35em :align: center :alt: alternate text :figclass: align-center :name: aligner-rasters-qgis-2 Paramétrage de l'alignement des rasters. Nos deux rasters apparaissent dans le cadre ``Couches rasters à aligner``. Il faut ensuite choisir une ``Couche de référence``, ici *S2_B11_20m*. Nous pouvons laisser les autres options vides. Si nos deux rasters ne se superposaient pas totalement, nous pourrions régler l'emprise du résultat dans le cadre ``Découper selon l'emprise actuelle``. Enfin nous cliquons sur :guilabel:`OK`. Bien que deux fichiers rasters ont été générés, celui qui a servi de référence n'est en fait pas modifié. .. _aligner-rasters-OTB: Aligner deux rasters avec OTB ******************************** Version de OTB : 9.0.0 Le logiciel OrfeoToolbox (:ref:`logiciels-OTB`) propose un module pour aligner un raster sur un autre. Ce module se nomme *Superimpose*. Comme tout module de OTB, il est possible de l'utiliser via l'interface graphique QGIS ou via la ligne de commande. Avec OTB via QGIS +++++++++++++++++++ Nous trouvons l'application *Superimpose* dans la :menuselection:`Boîte à outils de traitements --> OTB --> Geometry --> Superimpose`. Une fois l'application lancée, la fenêtre suivante apparaît (:numref:`superimpose`). Dans l'exemple, nous alignons une bande spectrale Sentinel-2 à 10 m (*B3*) sur une bande spectrale à 20 m (*B11*). .. figure:: figures/fen_otb_monteverdi_aligner_rasters.png :width: 50em :align: center :alt: alternate text :figclass: align-center :name: superimpose Aligner deux rasters dans OTB. À la ligne ``Reference input`` nous pointons vers le raster à utiliser comme référence, à savoir ici *S2_B11_20m.tiff*. À la ligne, ``The image to reproject`` nous pointons vers le raster à aligner, ici *S2_B03_10m.tiff*. Dans le menu déroulant ``Interpolation``, nous choisissons une méthode de rééchantillonnage. Comme nous modifions un raster à valeurs continues, nous pouvons choisir *linear* ou *bco*. Si nous avions travaillé avec un raster à valeurs discrètes, comme une occupation du sol, nous aurions choisi *nn* (*Nearest Neighboor*). Puis, nous renseignons le chemin et le nom du raster résultat à la ligne ``Output image``. Les autres paramètres peuvent conserver leurs valeurs par défaut. Une fois tout paramétré, il ne reste plus qu'à cliquer sur :guilabel:`Exécuter`. Avec la commande otbcli +++++++++++++++++++++++++++++ Il est possible d'exécuter l'application *Superimpose* via la ligne de commande otbcli. La synataxe est la suivante : .. code-block:: sh otbcli_Superimpose -inr S2_B11_20m.tiff -inm S2_B03_10m.tiff -interpolator bco -out S2_B03_resample_otbcli.tiff Où : * *otbcli_Superimpose* est l'appel au module *Superimpose* * *-inr* est le mot clef pour spécifier le raster à utiliser comme référence * *-inm* est le mot clef pour spécifier le raster à aligner * *-interpolator* est le mot clef pour choisir la méthode d'interpolation (*bco* pour *Bicubic interpolation* et *nn* pour *Nearest neighbor interpolation*) * *-out* est le mot clef pour le raster de sortie aligné .. _aligner-rasters-R: Aligner deux rasters avec R ****************************** Version de R : 4.2.2 Version de Terra : 1.5.34 À l'aide du package *terra* de R, il est facile d'aligner un raster sur un autre. Il suffit d'utiliser la fonction *resample*, comme présenté ci-après, où la bande 3 à 10 m d'une image Sentinel-2 est alignée sur la bande 11 à 20 m. .. code-block:: R # chargement de la librairie terra library(terra) # chargement du raster qui devra être aligné sur un autre s2_b3_10m <- terra::rast('S2_B03_10m.tiff') # chargement du raster sur lequel le premier sera aligné s2_b11_20m <- terra::rast('S2_B11_20m.tiff') # alignement de la bande à 10m sur la bande à 20m b3_20m <- terra::resample(s2_b3_10m, s2_b11_20m) À la fin de ce code, un nouveau raster de type *SpatRaster* est créé dans la variable *b3_20m*. Et ce nouveau raster est identique au premier mais qui a été aligné sur le second. .. |icone_browse| image:: figures/icone_browse.png :width: 20 px .. |icone_plus| image:: figures/icone_export_3D_add_frame.png :width: 18 px .. |icone_refresh| image:: figures/icone_actualiser_SCP.png :width: 20 px .. |icone_identifier| image:: figures/icone_identifier.png :width: 20 px .. |icone_otb_plus| image:: figures/icone_otb_plus.png :width: 20 px .. |icone_bas| image:: figures/icone_qgis_deplacer_champ_bas.png :width: 20 px .. |icone_haut| image:: figures/icone_qgis_deplacer_champ_haut.png :width: 20 px .. |icone_moins| image:: figures/icone_qgis_moins.png :width: 20 px .. _bassin de la Roya: https://fr.wikipedia.org/wiki/Roya .. _muparser: https://beltoforion.de/en/muparser/index.php .. _la documentation en ligne: https://www.orfeo-toolbox.org/CookBook/Applications/app_BandMath.html .. _page du manuel dédiée: https://gdal.org/programs/gdal_calc.html .. _algorithmes peut être trouvée en ligne: https://kokoalberti.com/articles/geotiff-compression-optimization-guide/?utm_source=pocket_mylist .. _à cette page: https://kokoalberti.com/articles/geotiff-compression-optimization-guide/?utm_source=pocket_mylist .. _options sont possibles: https://gdal.org/drivers/raster/gtiff.html#creation-options .. _Plotly: https://plotly.com/python/ .. _couleur des bâtons: https://r-charts.com/colors/ .. _regarder de plus près: https://rdrr.io/cran/terra/man/classify.html