How to add european countries borders to a Raster in R?

Question

I have this raster file:

fdiff_north:

class      : RasterLayer 
dimensions : 5500, 12000, 6.6e+07  (nrow, ncol, ncell)
resolution : 100, 100  (x, y)
extent     : 4e+06, 5200000, 2250000, 2800000  (xmin, xmax, ymin, ymax)
crs        : +proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +units=m +no_defs 
source     : memory
names      : layer 
values     : -1, 1  (min, max)

to get a plot:

cl <- colorRampPalette(c("red","white","blue"))(100)
plot(fdiff_north, col=cl)

I need boundaries to make the plot more usable.

I also tried to change the crs of my raster:

rfdiff_north:

class      : RasterLayer 
dimensions : 12571, 6151, 77324221  (nrow, ncol, ncell)
resolution : 9e-04, 0.00128  (x, y)
extent     : 42.7836, 48.3195, 5.672145, 21.76303  (xmin, xmax, ymin, ymax)
crs        : +proj=longlat +datum=WGS84 +no_defs 
source     : r_tmp_2022-02-07_151752_9872_54603.grd 
names      : layer 
values     : -1.466443, 1.366832  (min, max)

I want to add european countries boundaries to my raster, but I can't get it. I've tried with many different functions but I achieved no results. Can you help me to do that ?

Answer 1

The boundaries of European countries you can get from spMaps or any other source (OpenStreetMap data for example)

# devtools::install_github("rte-antares-rpackage/spMaps", ref = "master")
library(spMaps)
eu <- getEuropeCountries(mergeCountry = FALSE)

Trying to create similar (in size and CRS) raster:

library(terra)
#> terra 1.5.17
x <- rast(nrows=5500, ncols = 12000, nlyrs = 1,
     xmin = 4e+06, xmax = 5200000, ymin = 2250000, ymax = 2800000,
     crs = "+proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +units=m +no_defs",
     resolution = 100)
x <- setValues(x,runif(ncell(x)))
x
#> class       : SpatRaster 
#> dimensions  : 5500, 12000, 1  (nrow, ncol, nlyr)
#> resolution  : 100, 100  (x, y)
#> extent      : 4e+06, 5200000, 2250000, 2800000  (xmin, xmax, ymin, ymax)
#> coord. ref. : +proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +units=m +no_defs 
#> source      : memory 
#> name        :        lyr.1 
#> min value   : 4.889444e-09 
#> max value   :            1

Now it's time to convert SpatialPolygons to SimpleFeature and transform the coordinates to yours CRS:

library(sf)
eu <- st_as_sf(eu)

plot(x)
plot(st_transform(eu$geometry, 
    crs = "+proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +units=m +no_defs"),
    add = TRUE)

^{Created on 2022-02-09 by the reprex package (v2.0.1)}

You can crop boundaries to your raster bbox with st_crop()

Regards, Grzegorz