2.1.3 SRTM

Les variables topographiques telles que l’élévation, la pente et l’aspect ne sont actuellement utilisées ni pour l’évaluation des forestières ni pour l’évaluation de la biomasse, car il a été constaté que ces variables ont peu d’influence sur le pouvoir explicatif des modèles de classification ou de régression. Néanmoins, les données sont conservées et disponibles pour autres analyses.

2.1.3.1 Acquisition des données

Les données topographiques du SRTM ont été obtenues à partir de deux sources:

• Les données originales avec une résolution spatiale de 1 seconde d’arc (environ 30 mètres sur l’équateur) ont été obtenues du jeu de données SRTM 1 Arc-Second Global disponible sur USGS Earthexplorer. Ces données à haute résolution ont des fois des lacunes (manque des données).

• Les données SRTM ajustées avec une résolution de 3 secondes d’arc ont été obtenues à partir de CGIAR SRTM 90m Digital Elevation Database v4.1 pour combler les lacunes des données à haute résolution.

2.1.3.2 Prétraitement des données

Les données de 1 seconde d’arc sont lues et les lacunes éventuelles sont comblées avec les données ajustées de 3 secondes d’arc. Enfin, les données sont reprojetées sur le raster des images Landsat (UTM 31, résolution de 30 mètres) et coupées à la taille.

Script R: 01_SSTS/01_data/_src/prep-SRTM.R

###############################################################################
# prep-SRTM.R: lire, nettoyer et empiler des données topographiques SRTM
# -----------------------------------------------------------------------------
# Bern University of Applied Sciences
# Oliver Gardi, <oliver.gardi@bfh.ch>
# 13 Mai 2020

# Définitions des variables ===================================================
IN.DIR  <- paste0(DIR.RAW.DAT, "/SRTM")
OUT.DIR <- paste0(DIR.SST.DAT, "/SRTM")
if(!dir.exists(OUT.DIR)) dir.create(OUT.DIR)


# Préparation du modèle numérique d'élévation SRTM ============================

# Lire 90m SRTM MNE sans vides (source: http://srtm.csi.cgiar.org/srtmdata/),
dem.90 <- do.call(merge, lapply(as.list(dir(paste0(IN.DIR, "/3arcsecond"), 
                                            pattern=".*[.]tif$", 
                                            full.names=TRUE)), raster))


# Lire 30m SRTM MNE (source: USGS Earthexplorer) et remplir vides avec 90m SRTM 
dem.30 <- foreach(tile=dir(paste0(IN.DIR, "/1arcsecond"), pattern=".*[.]tif$", full.names=TRUE),
                  .combine=merge, .multicombine=TRUE) %dopar% {
                    dem.30.t <- raster(tile)
                    dem.90.t <- round(projectRaster(dem.90, dem.30.t))
                    merge(dem.30.t, dem.90.t)
                  }

# Reprojection d'images MNE vers Landsat (résolution 30m, UTM 31, ...) 
writeRaster(dem.30, paste0(OUT.DIR, "/SRTM-1arcsec_raw.tif"), 
            datatype="INT2S", 
            overwrite=TRUE)
system(paste("gdalwarp",
             paste0(OUT.DIR, "/SRTM-1arcsec_raw.tif"),
             "-t_srs '+proj=utm +zone=31 +datum=WGS84'",
             "-tr 30 30",
             paste("-te", TGO.EXT@xmin, TGO.EXT@ymin, TGO.EXT@xmax, TGO.EXT@ymax),
             paste0(OUT.DIR, "/SRTM-1arcsec.tif"),
             "-ot 'Int16'",
             "-co COMPRESS='LZW'",
             "-co INTERLEAVE='BAND'",
             "-overwrite"))
file.remove(paste0(OUT.DIR, "/SRTM-1arcsec_raw.tif"))

# Calculer les statistiques GDAL (min, max, ...)
system(paste("gdalinfo -stats", paste0(OUT.DIR, "/SRTM-1arcsec.tif")))


# Créer une vignettes du MNE ==================================================
dem <- raster(paste0(OUT.DIR, "/SRTM-1arcsec.tif"))
jpeg(paste0(OUT.DIR, "/SRTM-1arcsec.jpeg"), width=1350, height=3000)
par(plt=c(0,1,0,1))
plot(dem)
plot(mask(dem, TGO, inverse=TRUE), col="#FFFFFF66", legend=FALSE, add=TRUE)
plot(TGO, add=TRUE, lwd=3)
dev.off()