-
-
Save DrRoad/8025b798b5f7ef3e12de23c1c84f89bd to your computer and use it in GitHub Desktop.
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| library(raster) | |
| library(cartography) | |
| library(sf) | |
| library(SpatialPosition) | |
| mtq <- st_read(system.file("shape/martinique.shp", package="cartography")) | |
| # use WGS84 proj | |
| mtq_latlon <- st_transform(mtq, 4326) | |
| # this call throw an error but seems to work... | |
| getData('SRTM', lon=-61, lat=14) | |
| # import raster | |
| ras <- raster("srtm_24_10.tif") | |
| # crop on martinique area | |
| mtq_ras <- crop(ras, st_bbox(mtq_latlon)[c(1,3,2,4)]) | |
| # aggregate the raster | |
| mtq_ras <- aggregate(mtq_ras, fact=4,fun=mean) | |
| mtq_ras <- projectRaster(mtq_ras, crs = "+proj=utm +zone=20 +datum=WGS84 +units=m +no_defs") | |
| # break values | |
| bv <- c(seq(0,1300,100),1339) | |
| # contour extraction | |
| mtq_cont <- rasterToContourPoly(r = mtq_ras, breaks = bv, mask = as(mtq, "Spatial")) | |
| # custom palette | |
| pal <- c("#5D9D52", "#8DBC80", "#B8D9A9", "#FDEBBE", "#F7E0AC", "#F2D69B", | |
| "#EDCC8A", "#E8C279", "#E2B563", "#DBA84C", "#D49B36", "#BA8428", | |
| "#9A6A1E", "#7B5114") | |
| # sp to sf | |
| k <- st_as_sf(mtq_cont) | |
| # order the sf | |
| k <- k[order(k$center),] | |
| png(filename = "mtq.png", width = 700, height = 800, res = 100, bg = NA) | |
| par(mar = c(0,0,1.2,0)) | |
| plot(st_geometry(mtq), col = NA, border = NA, bg = "lightblue") | |
| for(i in 1:nrow(k)){ | |
| p <- st_geometry(k[i,]) | |
| plot(p + c(-50, 50), add=T, border = "#ffffff90",col = "#ffffff90") | |
| plot(p + c(100, -100), col = "#00000090", add=T, border = "#00000090") | |
| plot(p, col = pal[i], border = "NA", add=T) | |
| } | |
| legendChoro(pos = "left", breaks = bv, col = pal, nodata = F, | |
| title.txt = "Elevation\n(metres)", cex = 0.8) | |
| layoutLayer(title = "Martinique Relief", north = T, | |
| sources = 'T. Giraud, 2018', author = "SRTM, 2018", col = "lightblue", | |
| tabtitle = T, coltitle = "black") | |
| dev.off() | |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment