Fundamental misunderstanding with crs and map projections in R

Question

I think I have a fundamental misunderstanding with how using crs and map projections work in R. First let me show what my final goal is here. I would like to create something like the following: Using rnaturalearth. Specifically I am interested in the curved longitude lines and the diagonal latitude lines.

So far, I have used

library(rnaturalearth)
library(rnaturalearthdata)
library(ggplot2)
library(sf)

gl <- ne_countries(scale = "large", country = "Greenland",returnclass = "sf")
ggplot() + 
  geom_sf(data = gl)

to create the following

The part I am struggling with is indicating these latitude and longitude lines. I know this has something to do with how I project Earth onto the 2D plane (by changing the crs) but I am struggling to tie it all together in R. I believe I have a fundamental misunderstanding of how all this works and could do with some explanation.

I should note that my final goal is to indicate (through the use of dashed lines) areas bounded by longitudes 60,62,64,66 and 68.

Answer 1

The CRS used in your picture looks like EPSG:3184, which is the appropriate longitudinal slice of the Universal Transverse Mercator system. As long as you tell coord_sf this is the co-ordinate reference system that you want to replicate, it will automatically convert your map to that system for plotting. You can do this as follows:

library(rnaturalearth)
library(rnaturalearthdata)
library(ggplot2)
library(sf)

north <- 72.5; south <- 59; east <- -52; west <- -45

extent <- st_sf(a = 1:2, crs = 4326,
            geom = st_sfc(st_point(c(east, south)), 
                          st_point(c(west, north)))) |>
            st_transform(crs = 'EPSG:3184') |>
            st_bbox()

lat <- lapply(c(60, 62, 64, 66, 68), 
            function(x) cbind(seq(-70, -40, 0.1), x)) |>
  st_multilinestring() |> st_sfc(crs = 4326) |> st_sf() |>
  st_transform(crs = 'EPSG:3184') 

lon <- lapply(c(-60, -55, -50, -45, -40), 
              function(x) cbind(x, seq(56, 74, 0.2))) |>
  st_multilinestring() |> st_sfc(crs = 4326) |> st_sf() |>
  st_transform(crs = 'EPSG:3184') 
  
ne_countries(scale = "large", country = "Greenland",returnclass = "sf") |>
  ggplot() + 
  geom_sf(fill = 'white') +
  geom_sf(data = lat, linetype = 2, color = 'gray50') +
  geom_sf(data = lon, color = 'gray', linewidth = 0.2) +
  coord_sf(crs = 'EPSG:3184', xlim = extent[c(1, 3)], ylim = extent[c(2, 4)]) +
  scale_x_continuous(breaks = -5 * 1:10) +
  theme_bw() +
  theme(panel.background = element_rect(fill = '#f4f6f5'),
        panel.grid = element_blank())