.Rhistory

text(14.4*k - d, y = 57.4*k, "Discontinuities",
cex = 0.6, pos = 4, font = 2, col="#404040")
# Intitulé classe 1
text(15.5*k - d, y = 56.6*k, "Between two regions of the same country",
cex = 0.4, pos = 4, font = 1, col="#404040")
# Intitulé classe 2
text(15.5*k - d, y = 55.7*k, "Between two regions of two different countries",
cex = 0.4, pos = 4, font = 1, col="#404040")
# Description
text(15.5*k - d, y = 55.3*k,
"(The height is proportional to the value of the dicontinuity)",
cex = 0.4, pos = 4, font = 1, col="#404040")
# NB
text(10.7*k, y = 54.4*k,
"NB: only the 10% highest discontinuities are represented on the map.",
cex = 0.4, pos = 4, font = 3, col="#404040")
# SYMBOLE pour la légende
# Récupération d'un ligne infra-nationale presque horizontale
myline <- disc[disc$id == "TR21_BG341",]
# On duplique a ligne
myline2 <- myline
# Translation de la ligne pour la positionner dans la légende (Classe 1)
st_geometry(myline) <- st_geometry(myline) + c(4.2*k, 5*k)
# Translation de la ligne pour la positionner dans la légende (Classe 2)
myline2 <- myline
st_geometry(myline2) <- st_geometry(myline2) + c(0, 1.5*k)
# Construction du symbole pour les discontinuités entre deux pays
plot(myline, col= "#66666690",lwd = 0.5 ,add= T)
# 40 itération de la ligne
for (i in 1:40){
myline <- st_geometry(myline) + c(0,delta)
plot(st_geometry(myline), col= "#66666690",lwd = 0.5 ,add= T)
}
# Dernière ligne en rouge
plot(myline, col= "#cf0e00",lwd = 1.2 ,add= T)
# Construction du symbole pour les discontinuités infra-nationale
plot(myline2, col= "#66666690",lwd = 0.5 ,add= T)
# 8 itération de la ligne
for (i in 1:8){
myline2 <- st_geometry(myline2) + c(0,delta)
plot(st_geometry(myline2), col= "#66666690",lwd = 0.5 ,add= T)
}
plot(myline2, col= "#66666690",lwd = 0.5 ,add= T)
# Ne pas oublier
dev.off()
# Chunk 1: setup
## Global options
knitr::opts_chunk$set(echo=TRUE,
cache=FALSE,
prompt=FALSE,
comment=NA,
message=FALSE,
warning=FALSE,
class.source="bg-info",
class.output="bg-warning")
# Chunk 2
library("sf")
library("mapsf")
library("eurostat")
library("rnaturalearth")
library("scales")
# Chunk 3
# NUTS 2016
nuts2016 <- get_eurostat_geospatial(
output_class = "sf",
resolution = "20",
nuts_level = "all",
year = "2016")
# Couche géographique des NUTS3 et NUTS2 pour 2016
nuts2016_3 <- nuts2016[nuts2016$LEVL_CODE == 3, ]
nuts2016_2 <- nuts2016[nuts2016$LEVL_CODE == 2, ]
# Liste des pays (code ISO) sans découpage NUTS3
N2 <- c("AT", "BE", "CH", "DE", "EL", "NL", "UK", "TR", "IE", "IS", "NO")
# Couche géographique NUTS2/3 en fonction de la liste de pays N2
nuts <- rbind(nuts2016_2[nuts2016_2$CNTR_CODE %in% N2, ],
nuts2016_3[!nuts2016_3$CNTR_CODE %in% N2, ])
# Suppression des départements d'Outre-mer français
nuts <- nuts[!nuts$id %in% c("FRY10", "FRY20", "FRY30", "FRY40", "FRY50"), ]
# Suppression des NUTS 2016 du Royaume-Uni
nuts <- nuts[nuts$CNTR_CODE != "UK", ]
# Sélection et renommage de colonnes
nuts <- nuts[,c("id","NUTS_NAME","geometry")]
colnames(nuts) <- c("id","name","geometry")
# Chunk 4
# Récupération du découpage NUTS2 2013
nuts2013 <- get_eurostat_geospatial(
output_class = "sf",
resolution = "20",
nuts_level = "2",
year = "2013")
# Sélection des NUTS2 du Royaume-Uni
uk <- nuts2013[nuts2013$CNTR_CODE == "UK",]
# Sélection et renommage de colonnes
uk <- uk[,c("id","NUTS_NAME","geometry")]
colnames(uk) <- c("id","name","geometry")
# Fusion des NUTS2 2013 du Royaume-Uni avec le fond de carte NUTS2/3 2016
nuts <- rbind(nuts, uk)
# Chunk 5
par(mar = c(0, 0, 0, 0), mfrow = c(1, 3))
mf_map(nuts, col = "#CCCCCC", border = NA)
mf_map(nuts[substr(nuts$id,1,2) %in%  c("AT", "BE", "CH", "DE", "EL", "NL", "TR", "IE", "IS", "NO"),],
col = "#6eb1db",
border = "white",
lwd = 0.2, add = TRUE)
mf_title("NUTS 2 (version 2016)", bg = "#6eb1db")
mf_map(nuts, col = "#CCCCCC", border = NA)
mf_map(nuts[!substr(nuts$id,1,2) %in%  c("AT", "BE", "CH", "DE", "EL", "NL", "UK", "TR", "IE", "IS", "NO", "UK"),],
col = "#6eb1db",
border = "white",
lwd = 0.2, add = TRUE)
mf_title("NUTS 3 (version 2016)", bg = "#6eb1db")
mf_map(nuts, col = "#CCCCCC", border = NA)
mf_map(nuts[substr(nuts$id,1,2) == "UK",],
col = "#6eb1db",
border = "white",
lwd = 0.2, add = TRUE)
mf_title("NUTS 2 (version 2013)", bg = "#6eb1db")
# Chunk 6
# Récupération des données (PIB/hab)
var <- "nama_10r_3gdp"
gdpinh <- get_eurostat(var, time_format = "num")
# Sélection de l'unité de mesure et de la date
gdpinh <- gdpinh[gdpinh$unit == "EUR_HAB",]
gdpinh <- gdpinh[gdpinh$time == 2016, c("geo","values")]
colnames(gdpinh) <- c("id","GDPINH_2016")
# Chunk 7
# Import des donnèes
missing <- read.csv("data/missing.csv")
# Fusion des deux tableaux
gdpinh <- rbind(gdpinh, missing)
# Chunk 8
write.csv(gdpinh, "data/gdpinh.csv")
# Chunk 9
nuts <- merge(
x = nuts,
y = gdpinh,
by = "id",
all.x = TRUE)
# Chunk 10
# Surface terrestre, à petite échelle
land <- ne_download(
scale = 110,
type = "land",
category = "physical",
returnclass = "sf")
# Chunk 11
mf_map(land, border = NA, col = "#6eb1db")
# Chunk 12
# Mers et océans, à petite échelle
ocean <- ne_download(
scale = 110,
type = "ocean",
category = "physical",
returnclass = "sf")
# Chunk 13
mf_map(ocean, border = NA, col = "#6eb1db")
# Chunk 14
graticule = st_graticule(
crs = st_crs(4326),
ndiscr = 100,
lon = seq(-180, 180, by = 2),
lat = seq(-90, 90, by = 1),
margin = 0.01)
# Chunk 15
mf_map(graticule, col = "#6eb1db")
# Chunk 16
# Construction du rectangle (bounding box)
bb <- st_as_sfc(x = st_bbox(c(xmin = -50 , xmax = 70, ymin = 20, ymax = 80),
crs = st_crs(4326)))
# Chunk 17
sf_use_s2(FALSE)
ocean <- st_intersection(ocean, bb)
ocean <-  st_segmentize(ocean, 100)
land <- st_intersection(land, bb)
land <-  st_segmentize(land, 100)
graticule <- st_intersection(graticule, bb)
sf_use_s2(TRUE)
# Chunk 18
# Définition de la projection en format "proj-strings" (PROJ.4)
ortho <- "+proj=ortho +lat_0=-10 +lon_0=15 +x_0=0 +y_0=0
+ellps=WGS84 +units=m +no_defs"
# Projection des couches géographiques
ocean <- st_transform(ocean, ortho)
land <- st_transform(land, ortho)
graticule <-  st_transform(graticule, ortho)
nuts <- st_transform(nuts, ortho)
# Chunk 19
par(mar = c(0, 0, 0, 0), mfrow = c(2, 2))
mf_map(land, col = "#6eb1db", border = NA)
mf_title("land", bg = "#6eb1db")
mf_map(ocean, col = "#6eb1db", border = NA)
mf_title("ocean", bg = "#6eb1db")
mf_map(graticule, col = "#6eb1db", lwd = 1)
mf_title("graticule", bg = "#6eb1db")
mf_map(nuts, col = "#6eb1db", border = "white", lwd = 0.2)
mf_title("NUTS2/3", bg = "#6eb1db")
# Chunk 20
# Union de l'ensemble des entités de la couche géographique 'nuts'
fr <- st_union(nuts[substr(nuts$id, 1, 2) == "FR", ])
# Paramétrage des marges de la fenêtre graphique
par(mar = c(0, 0, 0, 0))
# Affichage multiple de la couche, en appliquant à chaque fois un léger décalage
mf_map(fr + c(5000,-5000), col = "#827e6c40", border = NA)
mf_map(fr + c(10000,-10000), col = "#827e6c40", border = NA, add = TRUE)
mf_map(fr + c(15000,-15000), col = "#827e6c40", border = NA, add = TRUE)
mf_map(fr + c(20000,-20000), col = "#827e6c40", border = NA, add = TRUE)
mf_map(fr + c(25000,-25000), col = "#827e6c40", border = NA, add = TRUE)
mf_map(fr, col = "#6eb1db", border = "white", lwd = 0.1, add = TRUE)
# Chunk 21
# Création de la variable k
k <- 100000
# Coordonnées de l'emprise * k
extent <- c(-20, 42, 24.5, 63) * k
# Construction de l'emprise (objet sfc)
bb <- st_as_sfc(x= st_bbox(c(xmin = extent[1], xmax = extent[3], ymin = extent[2], ymax = extent[4]),
crs = st_crs(nuts)))
# Chunk 22
# Création de la fonction
template = function(file) {
# Création d'un thème
theme <- mf_theme(x = "default",
bg = "#f2efe6",
fg = "#f2efe6",
mar = c(0, 0, 0, 0),
tab = TRUE,
pos = "left",
inner = FALSE,
line = 2,
cex = 1.9,
font = 3)
# Paramétrage de l'export de la carte en format png
mf_export(bb,
export = "png",
width = 2000,
filename = file,
res = 150,
theme = theme,
expandBB = c(-.02, 0, 0.05, 0))
# Affichage de la couche 'ocean'
mf_map(ocean, col = "#9acbe3", border = "#9acbe3", lwd = 5, add = TRUE)
# Affichage de la couche 'graticule'
mf_map(graticule, col = "#FFFFFF80", lwd = 1.5, lty = 3, add = TRUE)
# Affichage d'un effet d'ombrage pour la couche 'nuts'
ue <- st_union(nuts)
mf_map(ue + c(5000,-5000), col = "#827e6c40", border = NA, add = TRUE)
mf_map(ue + c(10000,-10000), col = "#827e6c40", border = NA, add = TRUE)
mf_map(ue + c(15000,-15000), col = "#827e6c40", border = NA, add = TRUE)
mf_map(ue + c(20000,-20000), col = "#827e6c40", border = NA, add = TRUE)
mf_map(ue + c(25000,-25000), col = "#827e6c40", border = NA, add = TRUE)
# Affichage de la couche 'nuts'
mf_map(nuts, col = "#dbccb6", border = "white", lwd = 0.3, add = TRUE)
# Affichage d'un bandeau bleu transparent sous le titre
rect(-22 * k,
41.3 * k,
28 * k,
41.3 * k + 250000,
border = NA, col = "#2369bd80")
# Affichage du titre
text(x = -21.5 * k,
y = 42.4 * k,
labels = "30 YEARS LATER, THE IRON CURTAIN IS STILL VISIBLE",
cex = 2.14,
pos = 4,
font = 2,
col = "#FFFFFF80")
# Affichage des sources
text(x = -21.75 * k,
y = 44.25 * k,
labels = "Map 100% designed in the R language by Nicolas Lambert, 2019.
Code source available here: https://github.com/neocarto/ironcurtain). Data sources: Eurostat & Natural Earth, 2019",
cex = 0.5,
pos = 4,
font = 1,
col = "#3f4654")
# Affichage d'une échelle
mf_scale(size = 700,
lwd = 0.6,
cex = 0.5,
col = "#3f4654",
pos = c(19 * k, y = 44 * k))
}
# Chunk 23
# Indiquez le nom du fichier png crée
# Utilisez dev.off() pour clôturer l'enregistrement du fichier
template("figures/fig1.png")
dev.off()
# Chunk 24
# Discrétisation par quantile (calcul des bornes de classes)
bks <- mf_get_breaks(x = nuts$GDPINH_2016,
nbreaks = 6,
breaks = "quantile")
# Vecteur de couleurs (1 par classe)
cols <- c("#50b160",
"#98c17e",
"#cce3c4",
"#fbf5bd",
"#fcc34f",
"#e97d40")
# Chunk 25
# Utilisation du modèle de mise en page (+ export de la carte)
template("figures/fig2.png")
# Cartographie du PIB par habitant
mf_map(x = nuts,
var = "GDPINH_2016",
type = "choro",
breaks = bks,
pal = cols,
lwd = 0.2,
leg_pos = "n",
add = TRUE)
# Ajout d'une légende
mf_legend(type = "choro",
pos = c(11 * k, 59.05 * k),
title = "",
val = bks,
val_cex = 0.4,
pal = cols,
fg = "#333333",
cex = 0.85,
border = "red",
val_rnd = 0,
no_data = FALSE,
frame = FALSE)
# Ajout d'un titre de légende
text(10.5 * k,
y = 59.1 * k,
labels = "Gross Domestic Product",
cex = 0.75,
pos = 4,
font = 2,
col = "#404040")
# Ajout d'un sous-titre de légende
text(10.5 * k,
y = 58.7 * k,
labels = "(in € per inh. in 2016)",
cex = 0.55,
pos = 4,
font = 1,
col = "#404040")
# NE PAS OUBLIER
dev.off()
# Chunk 26
# Auto-intersection de la couche NUTS (avec zone tampon de 5m)
nuts.borders  <- st_intersection(st_buffer(nuts, 5), st_buffer(nuts, 5))
# Transformation des géométries en 'MULTILINESTRING'
nuts.borders  <- st_cast(nuts.borders ,"MULTILINESTRING")
# Chunk 27
# Suppression des intersections entre un même polygone
nuts.borders  <- nuts.borders [nuts.borders $id != nuts.borders $id.1, ]
# Construction d'un identifiant unique pour chaque frontière
nuts.borders$id1 <- nuts.borders$id
nuts.borders$id2 <- nuts.borders$id.1
nuts.borders$id <-  paste0(nuts.borders$id1, "_", nuts.borders$id2)
rownames(nuts.borders) <- nuts.borders$id
nuts.borders  <- nuts.borders [,c("id","id1","id2","geometry")]
# Chunk 28
# Récupération des données de PIB par habitant, en supprimant la géométrie associée
vals <- st_set_geometry(x = nuts[, c("id","GDPINH_2016")],
value = NULL)
# Double jointure pour récupérer les valeurs des NUTS limitrophes
nuts.borders <- merge (x = nuts.borders, y = vals, by.x = "id1", by.y = "id", all.x = T)
nuts.borders <- merge (x = nuts.borders, y = vals, by.x = "id2", by.y = "id", all.x = T)
# Chunk 29
nuts.borders$disc <- nuts.borders$GDPINH_2016.x / nuts.borders$GDPINH_2016.y
# Chunk 30
threshold <- 0.95
disc <- nuts.borders[nuts.borders$disc >= quantile(nuts.borders$disc,threshold),]
# Chunk 31
template("figures/fig3.png")
mf_map(x = disc,
col = "#d92e94",
lwd = 3,
add = TRUE)
dev.off()
# Chunk 32
# On sélectionne une ligne (frontière) au hasard
line <- st_geometry(disc[5, ])
# Choix d'un nombre d'itérations
nb <- 15
# Choix d'une valeur de translation
delta <- 200
# Chunk 33
# Paramétrage des marges et de la couleur de fond
par(mar = c(0, 0, 0, 0), bg = "#f2efe6")
# Affichage de la frontière
mf_map(line, col = "#66666690", lwd = 0.5)
# 'nb' affichage de la frontière, décalée de 'delta' à chaque fois
for (j in 1:nb) {
line <- line + c(0, delta)
mf_map(line,
col = "#66666690",
lwd = 0.5 ,
add = TRUE)
}
# Affichage de la dernière ligne calculée en rouge, avec une épaisseur de 1.2
mf_map(line,
col = "#cf0e00",
lwd = 1.2,
add = TRUE)
# Chunk 34
# Transformation du type de géométrie
disc <- st_cast(disc,"LINESTRING")
# Calcul des coordonnées des centroïdes de chaque lignes
c <- as.data.frame(st_coordinates(st_centroid(disc)))
# Récupération de la valeur y de chaque centroïde
disc$Y <- c$Y
# Tri des discontinuités en fonction de la valeur de y
disc <- disc[order(disc$Y, decreasing = TRUE), ]
# Chunk 35
# Rééchelonnage des valeurs de discontinuité, de 30 à 70
disc$height <- round(rescale(disc$disc, to=c(30,70)),0)
# Chunk 36
disc$col <-"#cf0e00"
disc$thickness <- 1.2
# Chunk 37
disc$c1 <- substr(disc$id1,1,2)
disc$c2 <- substr(disc$id2,1,2)
# Chunk 38
# Boucle qui parcours toutes les lignes
for (i in 1:length(disc$disc)){
# Pays 1 == Pays 2 ?
if (disc$c1[i]== disc$c2[i]) {
# 8 itérations
disc$height[i] <- 8
# Couleur Grise
disc$col[i] <-"#66666690"
# Petite épaisseur
disc$thickness[i] <- 0.5
}
}
# Chunk 39
# Valeur de translation
delta <- 2500
# Création de la fonction extrude()
extrude <- function(id){
line <- st_geometry(disc[id,])
mf_map(line, col= "#66666690",lwd = 0.5 ,add= TRUE)
nb <- as.numeric(disc[id,"height"])[1]
for (j in 1:nb){
line <- st_geometry(line) + c(0,delta)
mf_map(st_geometry(line), col= "#66666690",lwd = 0.5 ,add= TRUE)
}
mf_map(line, col= disc$col[id],lwd = disc$thickness[id] ,add= TRUE)
}
template("figures/fig4.png")
# BOUCLE for pour appliquer la fonction extrude() à toutes les lignes
for (i in 1:length(disc$height))
{
extrude(i)
}
# TEXTE de la légende pour les discontinuités
# Utilisation du facteur k pour un placement précis des éléments
d = 0.75 * k
# Titre légende
text(14.4*k - d, y = 57.4*k, "Discontinuities",
cex = 0.6, pos = 4, font = 2, col="#404040")
# Intitulé classe 1
text(15.5*k - d, y = 56.6*k, "Between two regions of the same country",
cex = 0.4, pos = 4, font = 1, col="#404040")
# Intitulé classe 2
text(15.5*k - d, y = 55.7*k, "Between two regions of two different countries",
cex = 0.4, pos = 4, font = 1, col="#404040")
# Description
text(15.5*k - d, y = 55.3*k,
"(The height is proportional to the value of the dicontinuity)",
cex = 0.4, pos = 4, font = 1, col="#404040")
# NB
text(10.7*k, y = 54.4*k,
"NB: only the 10% highest discontinuities are represented on the map.",
cex = 0.4, pos = 4, font = 3, col="#404040")
# SYMBOLE pour la légende des discontinuités
# Récupération d'un ligne infra-nationale presque horizontale
myline1 <- disc[disc$id == "TR21_BG341",]
# Translation de la ligne 1 pour la positionner dans la légende (Classe 1)
st_geometry(myline1) <- st_geometry(myline1) + c(4.2*k, 5*k)
# On duplique a ligne
myline2 <- myline1
# Translation de la ligne 2 pour la positionner dans la légende (Classe 2)
st_geometry(myline2) <- st_geometry(myline2) + c(0, 1.5*k)
# Symbole de légende pour les discontinuités entre deux pays
plot(myline1, col= "#66666690",lwd = 0.5 ,add= T)
# 40 itération de la ligne
for (i in 1:40){
myline1 <- st_geometry(myline1) + c(0,delta)
plot(st_geometry(myline1), col= "#66666690",lwd = 0.5 ,add= T)
}
# Dernière ligne en rouge
plot(myline1, col= "#cf0e00",lwd = 1.2 ,add= T)
# Symbole de légende pour les discontinuités infra-nationale
plot(myline2, col= "#66666690",lwd = 0.5 ,add= T)
# 8 itération de la ligne
for (i in 1:8){
myline2 <- st_geometry(myline2) + c(0,delta)
plot(st_geometry(myline2), col= "#66666690",lwd = 0.5 ,add= T)
}
plot(myline2, col= "#66666690",lwd = 0.5 ,add= T)
# Ne pas oublier
dev.off()
citation(sf)
citation("sf")
rref <- bibentry(
bibtype = "misc",
title = "Le nouveau rideau de fer",
subtitle = "Un exemple de carte en 2.5D",
author = c("Nicolas Lambert"),
doi = "10.48645/a4ra-yr11",
url = "https://rzine.fr/publication_rzine/20191125_ironcurtain/",
keywords ="FOS: Other social sciences",
language = "fr",
publisher = "FR2007 CIST",
year = 2021,
copyright = "Creative Commons Attribution Share Alike 4.0 International")