Missing crs info when using `terra::writeCDF` followed by `raster::raster()` [and vice versa]

[elgabbas]

Hello,

I create a NetCDF file using terra::writeCDF(). When reading the NetCDF file using terra::rast(), it reads the CRS information correctly; but when loading the same file using raster::raster(), it does not read the CRS information.

Here is a reprex:

require(raster)
require(terra)

# Creating spatRast
Map1 <- terra::rast(
  nrows = 404, ncols = 390, nlyrs = 1, xmin = 2630000, xmax = 6530000, ymin = 1380000, ymax = 5420000,
  crs = "epsg:3035", vals = 1, names = "Map1")

Map1
# class       : SpatRaster 
# dimensions  : 404, 390, 1  (nrow, ncol, nlyr)
# resolution  : 10000, 10000  (x, y)
# extent      : 2630000, 6530000, 1380000, 5420000  (xmin, xmax, ymin, ymax)
# coord. ref. : ETRS89-extended / LAEA Europe (EPSG:3035) 
# source(s)   : memory
# name        : Map1 
# min value   :    1 
# max value   :    1 

# Save as netCDF
terra::writeCDF(Map1, filename = "Map1.nc", overwrite = TRUE)

Reading using terra

# CRS is read correctly
terra::rast("Map1.nc")
# class       : SpatRaster 
# dimensions  : 404, 390, 1  (nrow, ncol, nlyr)
# resolution  : 10000, 10000  (x, y)
# extent      : 2630000, 6530000, 1380000, 5420000  (xmin, xmax, ymin, ymax)
# coord. ref. : ETRS89-extended / LAEA Europe (EPSG:3035) 
# source      : Map1.nc 
# name        : Map1 

Reading using raster

# CRS is NA
raster::raster("Map1.nc")
# class      : RasterLayer 
# dimensions : 404, 390, 157560  (nrow, ncol, ncell)
# resolution : 10000, 10000  (x, y)
# extent     : 2630000, 6530000, 1380000, 5420000  (xmin, xmax, ymin, ymax)
# crs        : NA 
# source     : Map1.nc 
# names      : Map1 
# zvar       : Map1 

---

Similar results occurred when the NetCDF file was created with raster::writeRaster() then opened with terra::rast()

Map2 <- raster::raster(
    nrows = 404, ncols = 390,  xmn = 2630000, xmx = 6530000, ymn = 1380000, 
    ymx = 5420000, crs = "epsg:3035", vals = 1)
Map2
# class      : RasterLayer 
# dimensions : 404, 390, 157560  (nrow, ncol, ncell)
# resolution : 10000, 10000  (x, y)
# extent     : 2630000, 6530000, 1380000, 5420000  (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)

writeRaster(Map2, "Map2.nc", overwrite=TRUE, format="CDF", varname = "Map2")

# reading using `terra::rast`: missing CRS
terra::rast("Map2.nc")
# class       : SpatRaster 
# dimensions  : 404, 390, 1  (nrow, ncol, nlyr)
# resolution  : 10000, 10000  (x, y)
# extent      : 2630000, 6530000, 1380000, 5420000  (xmin, xmax, ymin, ymax)
# coord. ref. :  
#   source      : Map2.nc 
# name        :       Map2 
# unit        : m +no_defs

# reading using `raster::raster`: correct CRS
raster::raster("Map2.nc")
# class      : RasterLayer 
# dimensions : 404, 390, 157560  (nrow, ncol, ncell)
# resolution : 10000, 10000  (x, y)
# extent     : 2630000, 6530000, 1380000, 5420000  (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     : Map2.nc 
# names      : Map2 
# values     : 1, 1  (min, max)
# zvar       : Map2 

Is it possible to avoid this?