---
title: "World maps"
output:
html_document:
df_print: paged
---
```{r echo = FALSE, message = FALSE}
library(tidyverse)
library(sf)
```
## Longitude--latitude
```{r world-longlat, fig.width = 8.5}
world_sf <- sf::st_as_sf(rworldmap::getMap(resolution = "low"))
crs_longlat <- "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"
ggplot(world_sf) +
geom_sf(fill = "#E69F00B0", color = "black", size = 0.5/.pt) +
coord_sf(expand = FALSE, crs = crs_longlat) +
scale_x_continuous(name = "longitude", breaks = seq(-120, 120, by = 60)) +
scale_y_continuous(name = "latitude", breaks = seq(-60, 60, by = 30)) +
ggtitle("Cartesian longitude and latitude") +
theme_minimal() +
theme(
panel.background = element_rect(fill = "#56B4E950"),
panel.grid.major = element_line(color = "gray30", size = 0.25),
axis.ticks = element_line(color = "gray30", size = 0.5/.pt)
)
```
## Mercator
```{r world-mercator, fig.width = 8.5}
crs_mercator <- "+proj=merc"
# calculate bounding box in transformed coordinates
mercator_bbox <-
rbind(c(-180, -85), c(180, 85)) %>%
st_multipoint() %>%
st_sfc(
crs = "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"
) %>%
st_transform(crs = crs_mercator)
ggplot(world_sf) +
geom_sf(fill = "#E69F00B0", color = "black", size = 0.5/.pt) +
scale_x_continuous(name = NULL, breaks = seq(-120, 120, by = 60)) +
scale_y_continuous(name = NULL, breaks = seq(-60, 60, by = 30)) +
coord_sf(
xlim = mercator_bbox[[1]][, 1],
ylim = mercator_bbox[[1]][, 2],
expand = FALSE,
crs = crs_mercator
) +
ggtitle("Mercator") +
theme_minimal() +
theme(
panel.background = element_rect(fill = "#56B4E950", color = "white", size = 1),
panel.grid.major = element_line(color = "gray30", size = 0.25)
)
```
## Robinson
```{r world-robin, fig.asp = 0.6}
crs_robin <- "+proj=robin +lat_0=0 +lon_0=0 +x0=0 +y0=0"
# projection outline in long-lat coordinates
lats <- c(90:-90, -90:90, 90)
longs <- c(rep(c(180, -180), each = 181), 180)
robin_outline <-
list(cbind(longs, lats)) %>%
st_polygon() %>%
st_sfc(
crs = "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"
) %>%
st_transform(crs = crs_robin)
# bounding box in transformed coordinates
xlim <- c(-18494733, 18613795)
ylim <- c(-9473396, 9188587)
robin_bbox <-
list(
cbind(
c(xlim[1], xlim[2], xlim[2], xlim[1], xlim[1]),
c(ylim[1], ylim[1], ylim[2], ylim[2], ylim[1])
)
) %>%
st_polygon() %>%
st_sfc(crs = crs_robin)
# area outside the earth outline
robin_without <- st_difference(robin_bbox, robin_outline)
ggplot(world_sf) +
geom_sf(fill = "#E69F00B0", color = "black", size = 0.5/.pt) +
geom_sf(data = robin_without, fill = "white", color = NA) +
geom_sf(data = robin_outline, fill = NA, color = "grey30", size = 0.5/.pt) +
scale_x_continuous(name = NULL, breaks = seq(-120, 120, by = 60)) +
scale_y_continuous(name = NULL, breaks = seq(-60, 60, by = 30)) +
coord_sf(xlim = 0.95*xlim, ylim = 0.95*ylim, expand = FALSE, crs = crs_robin, ndiscr = 1000) +
ggtitle("Robinson") +
theme_minimal() +
theme(
panel.background = element_rect(fill = "#56B4E950", color = "white", size = 1),
panel.grid.major = element_line(color = "gray30", size = 0.25)
)
```
## Interrupted Goode homolosine
```{r world-goode, fig.width = 8.5}
crs_goode <- "+proj=igh"
# projection outline in long-lat coordinates
lats <- c(
90:-90, # right side down
-90:0, 0:-90, # third cut bottom
-90:0, 0:-90, # second cut bottom
-90:0, 0:-90, # first cut bottom
-90:90, # left side up
90:0, 0:90, # cut top
90 # close
)
longs <- c(
rep(180, 181), # right side down
rep(c(80.01, 79.99), each = 91), # third cut bottom
rep(c(-19.99, -20.01), each = 91), # second cut bottom
rep(c(-99.99, -100.01), each = 91), # first cut bottom
rep(-180, 181), # left side up
rep(c(-40.01, -39.99), each = 91), # cut top
180 # close
)
goode_outline <-
list(cbind(longs, lats)) %>%
st_polygon() %>%
st_sfc(
crs = "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"
) %>%
st_transform(crs = crs_goode)
# bounding box in transformed coordinates
xlim <- c(-21945470, 21963330)
ylim <- c(-9538022, 9266738)
goode_bbox <-
list(
cbind(
c(xlim[1], xlim[2], xlim[2], xlim[1], xlim[1]),
c(ylim[1], ylim[1], ylim[2], ylim[2], ylim[1])
)
) %>%
st_polygon() %>%
st_sfc(crs = crs_goode)
# area outside the earth outline
goode_without <- st_difference(goode_bbox, goode_outline)
ggplot(world_sf) +
geom_sf(fill = "#E69F00B0", color = "black", size = 0.5/.pt) +
geom_sf(data = goode_without, fill = "white", color = NA) +
geom_sf(data = goode_outline, fill = NA, color = "grey30", size = 0.5/.pt) +
scale_x_continuous(name = NULL, breaks = seq(-120, 120, by = 60)) +
scale_y_continuous(name = NULL, breaks = seq(-60, 60, by = 30)) +
coord_sf(xlim = 0.95*xlim, ylim = 0.95*ylim, expand = FALSE, crs = crs_goode, ndiscr = 1000) +
ggtitle("Interrupted Goode homolosine") +
theme_minimal() +
theme(
panel.background = element_rect(fill = "#56B4E950", color = "white", size = 1),
panel.grid.major = element_line(color = "gray30", size = 0.25)
)
```
## Winkel tripel
```{r world-winkel-tri, fig.width = 8.5}
# Winkel tripel projection needs to be done manually, is not
# supported by sf.
crs_wintri <- "+proj=wintri +datum=WGS84 +no_defs +over"
# world
world_wintri <- lwgeom::st_transform_proj(world_sf, crs = crs_wintri)
# graticule
grat_wintri <- sf::st_graticule(lat = c(-89.9, seq(-80, 80, 20), 89.9))
grat_wintri <- lwgeom::st_transform_proj(grat_wintri, crs = crs_wintri)
# earth outline
lats <- c(90:-90, -90:90, 90)
longs <- c(rep(c(180, -180), each = 181), 180)
wintri_outline <-
list(cbind(longs, lats)) %>%
st_polygon() %>%
st_sfc(
crs = "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"
) %>%
lwgeom::st_transform_proj(crs = crs_wintri)
ggplot() +
geom_sf(data = wintri_outline, fill = "#56B4E950", color = NA) +
geom_sf(data = grat_wintri, color = "gray30", size = 0.25/.pt) +
geom_sf(data = world_wintri, fill = "#E69F00B0", color = "black", size = 0.5/.pt) +
geom_sf(data = wintri_outline, fill = NA, color = "grey30", size = 0.5/.pt) +
ggtitle("Winkel tripel") +
coord_sf(datum = NA) +
theme_minimal()
```