Move the Canary Islands and zoom in on the peninsula on a ggplot map with r and sf

Viewed 135

I have tried to move the geometries of the Canary Islands towards the peninsula so that the map of Spain is seen on a larger scale to represent information by provinces.

I have started from the shapefile of this ESRI location

https://www.arcgis.com/home/item.html?id=83d81d9336c745fd839465beab885ab7

The code used is this, you need ggplot2 and sf


geoprovincias <- st_read("ruta/Provincias_ETRS89_30N.shp")

geopenínsula <- geoprovincias[!geoprovincias$CCAA=="Canarias",]
geocanarias <- geoprovincias[geoprovincias$CCAA=="Canarias",]

crs <- st_crs(geopenínsula)

geocanariasco <- geocanarias %>% 
 st_transform(crs) %>%
 st_geometry()
geocanariasco <- geocanariasco + c(5, 7)

geocanarias <- geocanarias %>% 
 st_set_geometry(geocanariasco) %>%
 st_set_crs(crs)

geoprovincias <- rbind(geopenínsula, geocanarias)

mprovincias <- geoprovincias %>% 
 ggplot() +
  geom_sf(
   fill="white",
   color="#1D1D1D",
   size=0.5
  ) +
 theme_void()
mprovincias

2 Answers

I have managed to solve from a different source,

http://centrodedescargas.cnig.es/CentroDescargas/index.jsp

Following the links: Reference geographic information - Municipal, provincial and autonomous limits - Download: lineas_limite.zip.

The path in the uncompressed folder: SIGLIM_Publico_INSPIRE:

  • SHP_ETRS89 - recintos_provinciales_inspire_peninbal_etrs89 - recintos_provinciales_inspire_peninbal_etrs89.shp
  • SHP_WGS84 - recintos_provinciales_inspire_canarias_wgs84 - recintos_provinciales_inspire_canarias_wgs84.shp

The code:

library(ggplot2)
library(sf)
library(rmapshaper)

geopropenínsula <- st_read("ruta/recintos_provinciales_inspire_peninbal_etrs89.shp")
geoprocanarias <- st_read("ruta/recintos_provinciales_inspire_canarias_wgs84.shp")

crs <- st_crs(geopropenínsula) #CRS Península, Baleares, Ceuta y Melilla

geoprocanariasco <- geoprocanarias %>% 
 st_transform(crs) %>% #We transform the geometries of the Canary Islands
 st_geometry() #We extract the geometries.
geoprocanariasco <- geoprocanariasco + c(5, 7) #Move Canarias

geoprocanarias <- geoprocanarias %>% 
 st_set_geometry(geoprocanariasco) %>%  #We change the geometries of the Canary Islands.
 st_set_crs(crs) # We will specify the CRS of the Península for the Islas Canarias again.

geoprovincias <- rbind(geopropenínsula, geoprocanarias) #Join
geoprovincias$CódigoProvincia <- substring(geoprovincias$NATCODE, 5, 6)

# We eliminate the territory not associated with any province.
geoprovincias <- geoprovincias[!geoprovincias$CódigoProvincia=="54",]

# We simplify the provincial boundaries and smooth the lines.
geosprovincias <- ms_simplify(geoprovincias, keep = 0.00050, weighting = 6)

mapapro <- geosprovincias %>%
 ggplot() +
  geom_sf() +
  theme_void()

Here is a solution using some coordinate limits and plot them separately on two plots then to combine them into one. Still need some work to make it more elegant.

library(ggplot2)
library(sf)
library(rnaturalearth)
library(gridExtra)

# my location of the shape file
geoprovincias <- st_read("/S/sinh/Documents/Professional/SO/data/Provincias_ETRS89_30N/Provincias_ETRS89_30N.shp")

geopenínsula <- geoprovincias[!geoprovincias$CCAA=="Canarias",]
geocanarias <- geoprovincias[geoprovincias$CCAA=="Canarias",]

crs <- st_crs(geopenínsula)

geocanariasco <- geocanarias %>% 
  st_transform(crs) %>%
  st_geometry()
geocanariasco <- geocanariasco + c(5, 7)

geocanarias <- geocanarias %>% 
  st_set_geometry(geocanariasco) %>%
  st_set_crs(crs)

geoprovincias <- rbind(geopenínsula, geocanarias)

# This part is for calculate the limit for sf cordinate #####
{
  # all countries at scale 10m
  ctrys <- ne_countries(scale = 10, type = "countries", returnclass = "sf")
  
  crsLONGLAT <- "+proj=longlat +datum=WGS84 +no_defs"
  crsLAEA <- "+proj=laea +lat_0=12 +lon_0=5 +x_0=4321000 +y_0=3210000 +datum=WGS84 +units=m +no_defs"
  
  bb <- st_sfc(
    st_polygon(list(cbind(
      c(-10, 5, 5, -10, -10), # x-coordinates (longitudes) of points A,B,C,D
      c(35, 35, 45, 45, 35)     # y-coordinates (latitudes) of points A,B,C,D
    ))),
    crs = crsLONGLAT)
  
  # now in in LAEA projection
  laeabb <- st_transform(bb, crs = crsLAEA)
  
  # the extent of the bounding box in the new projection
  b <- st_bbox(laeabb)
  big_part <- b
}

# the big part plot
mprovincias_1 <- geoprovincias %>% 
  ggplot() +
  geom_sf(
    fill="white",
    color="#1D1D1D",
    size=0.5
  ) + coord_sf(crs = crsLAEA, xlim = c(big_part["xmin"], big_part["xmax"]),
               ylim = c(big_part["ymin"], big_part["ymax"])) +
  theme_void()

# This part is for calculate the coordlimit for the small island
{
  bb <- st_sfc(
    st_polygon(list(cbind(
      c(-19, -13, -13, -19, -19), # x-coordinates (longitudes) of points A,B,C,D
      c(27, 27, 30, 30, 27)     # y-coordinates (latitudes) of points A,B,C,D
    ))),
    crs = crsLONGLAT)
  
  # now in in LAEA projection
  laeabb <- st_transform(bb, crs = crsLAEA)
  
  # the extent of the bounding box in the new projection
  b <- st_bbox(laeabb)
  small_island <- b
  }

# the island plot
mprovincias_2 <- geoprovincias %>% 
  ggplot() +
  geom_sf(
    fill="white",
    color="#1D1D1D",
    size=0.5
  ) + coord_sf(crs = crsLAEA, xlim = c(small_island["xmin"], small_island["xmax"]),
               ylim = c(small_island["ymin"], small_island["ymax"])) +
  theme_void() + 
  theme(panel.border = element_rect(colour = "black", size = 1, fill = NA))


# combine the two in one a bit clumsy
grid.arrange(arrangeGrob(mprovincias_1, mprovincias_2, ncol = 1, nrow = 2,
                         heights = c(4, 1)))

Map with zoom

Related