Segfault when using project()

[PMassicotte]

I have one terra raster:

r$> sic
class       : SpatRaster 
dimensions  : 448, 304, 12  (nrow, ncol, nlyr)
resolution  : 24917.24, 24943.67  (x, y)
extent      : -3837420, 3737422, -5337387, 5837377  (xmin, xmax, ymin, ymax)
coord. ref. : WGS 84 / NSIDC Sea Ice Polar Stereographic North (EPSG:3413) 
source      : SIC_month_3413_in.nc:sic 
varname     : sic 
names       : sic_1, sic_2, sic_3, sic_4, sic_5, sic_6, ... 

which I would like to reproject into another grid from another terraraster:

r$> dest_grid
class       : SpatRaster 
dimensions  : 361, 361, 1  (nrow, ncol, nlyr)
resolution  : 25073.32, 25073.32  (x, y)
extent      : -4535762, 4515707, -4515707, 4535762  (xmin, xmax, ymin, ymax)
coord. ref. : WGS 84 / NSIDC EASE-Grid 2.0 North (EPSG:6931) 
source      : AVHRR_albedo_grid.nc:cdr_surface_albedo 
varname     : cdr_surface_albedo (NOAA CDR of surface broadband albedo) 
name        : cdr_surface_albedo_Time=106.1666666666667 
unit        :                                         1 

If I use project(), it crashes my R session:

r$> project(sic, dest_grid)

 *** caught segfault ***
address 0xa0, cause 'memory not mapped'

Traceback:
 1: .External(list(name = "CppMethod__invoke_notvoid", address = <pointer: 0x55962dcb5f20>,     dll = list(name = "Rcpp", path = "/home/xxx/R/x86_64-pc-linux-gnu-library/4.2/Rcpp/libs/Rcpp.so",         dynamicLookup = TRUE, handle = <pointer: 0x55962eceff30>,         info = <pointer: 0x55962b691bc0>), numParameters = -1L),     <pointer: 0x5596366251a0>, <pointer: 0x55963662f960>, .pointer,     ...)
 2: x@ptr$warp(y@ptr, "", method, mask[1], align[1], FALSE, opt)
 3: .local(x, ...)
 4: project(sic, dest_grid)
 5: project(sic, dest_grid)

Possible actions:
1: abort (with core dump, if enabled)
2: normal R exit
3: exit R without saving workspace
4: exit R saving workspace

Selection:

However, if I convert to raster, and use projectRaster(), it works fine:

r$> raster::projectRaster(raster::raster(sic), raster::raster(dest_grid))
class      : RasterLayer 
dimensions : 361, 361, 130321  (nrow, ncol, ncell)
resolution : 25073.32, 25073.32  (x, y)
extent     : -4535762, 4515707, -4515707, 4535762  (xmin, xmax, ymin, ymax)
crs        : +proj=laea +lat_0=90 +lon_0=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs 
source     : memory
names      : sic_1 
values     : -2, 1  (min, max)

Any ideas what could cause the problem?

r$> sessionInfo()
R version 4.2.1 (2022-06-23)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Linux Mint 21

Matrix products: default
BLAS:   /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3
LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.20.so

locale:
 [1] LC_CTYPE=en_CA.UTF-8       LC_NUMERIC=C               LC_TIME=en_CA.UTF-8       
 [4] LC_COLLATE=en_CA.UTF-8     LC_MONETARY=en_CA.UTF-8    LC_MESSAGES=en_CA.UTF-8   
 [7] LC_PAPER=en_CA.UTF-8       LC_NAME=C                  LC_ADDRESS=C              
[10] LC_TELEPHONE=C             LC_MEASUREMENT=en_CA.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] terra_1.6-9

loaded via a namespace (and not attached):
[1] compiler_4.2.1   cli_3.3.0        Rcpp_1.0.9       codetools_0.2-18 jsonlite_1.8.0  
[6] rlang_1.0.4 

[rhijmans]

I cannot reproduce that. I see

library(terra)
sic <- rast(nrow=448, ncol=304, nlyr=12, ext=c(-3837420, 3737422, -5337387, 5837377), crs="EPSG:3413", vals=1)
dst <- rast(nrow=361, ncol=361, ext=c(-4535762, 4515707, -4515707, 4535762), crs="EPSG:6931")
 
project(sic, dst)
#class       : SpatRaster 
#dimensions  : 361, 361, 12  (nrow, ncol, nlyr)
#resolution  : 25073.32, 25073.32  (x, y)
#extent      : -4535762, 4515707, -4515707, 4535762  (xmin, xmax, ymin, ymax)
#coord. ref. : WGS 84 / NSIDC EASE-Grid 2.0 North (EPSG:6931) 
#source      : memory 
#names       : lyr.1, lyr.2, lyr.3, lyr.4, lyr.5, lyr.6, ... 
#min values  :     1,     1,     1,     1,     1,     1, ... 
#max values  :     1,     1,     1,     1,     1,     1, ... 

Does the above also work for you? If it does, it may be specific to the file (perhaps you can share it with me?) or the the versions of GDAL/PROJ

The above works for me with these versions of GDAL/PROJ

gdal(lib="")
#   gdal    proj    geos 
#"3.4.3" "7.2.1" "3.9.1" 

And also with

gdal(lib="")
#    gdal     proj     geos
# "3.4.0"  "8.2.0" "3.10.1" 

[PMassicotte]

Yes, this does work.

> gdal(lib="")
    gdal     proj     geos 
 "3.4.3"  "8.2.1" "3.10.2" 

I do not know what is going on. These two rasters are created within the code. If I try to project directly, it crashes. However, it works if I save these to files on disk, reopen them and then project.