Binding several sf together only plots first layer

Question

I'm trying to create a script to generate generic maps with bathymetry, and I'm struggling to get it to work. My issue is that depending on the map I want to make, I will call different bathymetric layers, but when I bind them together, only the first one is plotted. I don't get where the issue comes from. Any idea?

I created a function to load my bathymetric shapefiles from the 'rnaturalearth' package:

load_shapefile <- function(file, url, shapefile){ 
  if (!file.exists(file.path(tempdir(), file))) {
    url <- url
    download.file(url, file.path(tempdir(), file))
    unzip(file.path(tempdir(), file), exdir = tempdir())
  }
  st_read(dsn = tempdir(), layer = shapefile,
          quiet = TRUE)
}

I then load a few shapefiles:

bathy_200 <- load_shapefile("ne_10m_bathymetry_K_200.zip",
                            "https://www.naturalearthdata.com/http//www.naturalearthdata.com/download/10m/physical/ne_10m_bathymetry_K_200.zip",
                            "ne_10m_bathymetry_K_200")
bathy_1000 <- load_shapefile("ne_10m_bathymetry_J_1000.zip",
                             "https://www.naturalearthdata.com/http//www.naturalearthdata.com/download/10m/physical/ne_10m_bathymetry_J_1000.zip",
                             "ne_10m_bathymetry_J_1000")
bathy_2000 <- load_shapefile("ne_10m_bathymetry_I_2000.zip",
                             "https://www.naturalearthdata.com/http//www.naturalearthdata.com/download/10m/physical/ne_10m_bathymetry_I_2000.zip",
                             "ne_10m_bathymetry_I_2000")

Let's say I only want to plot the 200 and 2000 isobaths, but still want to create a generic function that would also easily allow me to plot 1000 and 2000, or 200 and 1000 instead:

plot_map <- function(bathy){
ggplot() +
  geom_sf(data = bathy %>%
            left_join(tibble(depth = c(200, 1000, 2000), 
                             fill = c("#E2EFF6", "#B7D7EA", "#8DBEDC")),
                      by = "depth"),
          aes(fill = fill),
          color = NA) +
  coord_sf(xlim = c(-20, 20),
           ylim = c(-20, 20),
           expand = FALSE) +
  scale_fill_identity()
}

When I plot my data, only the 200m shapefile is plotted...

plot_map(bind_rows(bathy_200, bathy_2000))

What dit I miss? I've tried several solutions, e.g. arranging data etc. but the outcome remains the same...

Answer 1

I have found that the git repository is corrupt for this file. When I look at the master repo this file has two identical copies and another file from the same directory is missing. I propose to copy them and then delete the two files involved. After committing the deletes I will re-add them.