writeRaster: invalid GDAL options cause R to abort when using a large raster

[dfriend21]

I've been playing around with setting GDAL options when using writeRaster. When using invalid GDAL settings with a relatively small raster, terra outputs a warning that prints out the GDAL error(s) that occurred, which is as expected. However, when doing the same thing with a larger raster, R aborts. I've run across this a couple times, with different GDAL options. Here's some reproducible code demonstrating the problem.

library(terra)

#========================================
# USING A SMALL RASTER
#========================================

dims_small <- c(1000,1000)
r_small <- rast(matrix(runif(prod(dims_small)), nrow = dims_small[1], ncol = dims_small[2]))

# try to write it out with invalid GDAL options
path_small <- paste0(tempfile(), ".tif")
writeRaster(r_small, 
            path_small,
            datatype = "FLT4S",
            gdal = c("NBITS=1",
                     "COMPRESS=CCITTFAX3"))
# the above code triggers a GDAL error, but the code runs and a file is created that consists entirely of NaN
r_small_read <- rast(path_small) 
head(values(r_small_read))
unlink(path_small)

#========================================
# USING A LARGE RASTER
#========================================

dims_large <- c(15000, 15000)
r_large <- rast(matrix(runif(prod(dims_large)), nrow = dims_large[1], ncol = dims_large[2]))

# try to write it out with invalid GDAL options
path_large <- paste0(tempfile(), ".tif")
# this next line causes R to abort
writeRaster(r_large, 
            path_large,
            datatype = "FLT4S",
            gdal = c("NBITS=1",
                     "COMPRESS=CCITTFAX3")) 
r_large_read <- rast(path_large)
head(values(r_large_read))
unlink(path_large)

Note that this also occurred with different GDAL options (for example, setting the gdal parameter to COMPRESS=JPEG in the above code produced the same result).

Here is the output of sessioninfo::session_info():

─ Session info ───────────────────────────────────────────────────────
 setting  value
 version  R version 4.2.2 (2022-10-31)
 os       macOS Monterey 12.4
 system   aarch64, darwin20
 ui       RStudio
 language (EN)
 collate  en_US.UTF-8
 ctype    en_US.UTF-8
 tz       America/Chicago
 date     2022-11-29
 rstudio  2022.07.2+576 Spotted Wakerobin (desktop)
 pandoc   NA

─ Packages ───────────────────────────────────────────────────────────
 package     * version date (UTC) lib source
 cli           3.4.1   2022-09-23 [1] CRAN (R 4.2.0)
 codetools     0.2-18  2020-11-04 [1] CRAN (R 4.2.0)
 Rcpp          1.0.9   2022-07-08 [1] CRAN (R 4.2.0)
 rstudioapi    0.14    2022-08-22 [1] CRAN (R 4.2.0)
 sessioninfo   1.2.2   2021-12-06 [1] CRAN (R 4.2.0)
 terra       * 1.6-41  2022-11-18 [1] CRAN (R 4.2.0)

 [1] /Library/Frameworks/R.framework/Versions/4.2-arm64/Resources/library

──────────────────────────────────────────────────────────────────────

[rhijmans]

I do not get a crash with the large raster on windows with GDAL 3.5.2 nor on Ubuntu with GDAL 3.4.0 or GDAL 2.2.3. I do not have access to a Mac right now.

What is the version of GDAL you are using (as returned by terra::gdal())? And how did you install it terra? That is, do you use a binary install from CRAN or did you install from source?

[kadyb]

I can confirm that this causes the crash. Large raster probably means that needs to be processed in blocks. I tested this example with 8 GB RAM, so maybe you should increase the number of pixels to reproduce it. I installed {terra} from CRAN.

Session info

``` R version 4.2.1 (2022-06-23 ucrt) Platform: x86_64-w64-mingw32/x64 (64-bit) Running under: Windows 8.1 x64 (build 9600) Matrix products: default attached base packages: [1] stats graphics grDevices utils datasets methods [7] base other attached packages: [1] terra_1.6-41 loaded via a namespace (and not attached): [1] compiler_4.2.1 tools_4.2.1 Rcpp_1.0.9.3 [4] codetools_0.2-18 gdal proj geos "3.5.2" "8.2.1" "3.9.3" ```

[dfriend21]

I installed terra from CRAN, and running terra::gdal() outputs "3.5.3".

[valentingar]

Maybe I can add to this? I experience the same problem with very large datasets. I can create these objects in memory, but cannot read or write them using terra.

I tested this using the CRAN version, as well as binary and source versions of 1.7.30.

Here is a reprex:

Sys.info()
#>          sysname          release          version         nodename 
#>        "Windows"         "10 x64"    "build 19045" "MN030124-C0014" 
#>          machine            login             user   effective_user 
#>         "x86-64"       "Valentin"       "Valentin"       "Valentin"
packageVersion("terra")
#> [1] '1.7.30'
terra::gdal()
#> [1] "3.5.2"

r <- terra::rast(system.file("external/test.grd", package="raster"))

tfile <- tempfile(fileext = ".grd")

# stacking to make a large raster
# this works
terra::rast(lapply(1:100, function(i) r)) 
#> class       : SpatRaster 
#> dimensions  : 115, 80, 100  (nrow, ncol, nlyr)
#> resolution  : 40, 40  (x, y)
#> extent      : 178400, 181600, 329400, 334000  (xmin, xmax, ymin, ymax)
#> coord. ref. : +proj=sterea +lat_0=52.1561605555556 +lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000 +datum=WGS84 +units=m +no_defs
#> Warning in (grepl("Windows", osVersion)) && (nchar(x) > 256): 'length(x) = 100
#> > 1' in coercion to 'logical(1)'
#> sources     : test.grd  
#>               test.grd  
#>               test.grd  
#>               ... and 97 more source(s)
#> names       :      test,      test,      test,      test,      test,      test, ... 
#> min values  :  138.7071,  138.7071,  138.7071,  138.7071,  138.7071,  138.7071, ... 
#> max values  : 1736.0579, 1736.0579, 1736.0579, 1736.0579, 1736.0579, 1736.0579, ...

# this as well
terra::writeRaster(terra::rast(lapply(1:50, function(i) r)),
                   filename = tfile, 
                   overwrite = TRUE)

# but not this
terra::writeRaster(terra::rast(lapply(1:100, function(i) r)),
                   filename = tfile, 
                   overwrite = TRUE) 
#> Warning: `C:/Users/Valentin/AppData/Local/Temp/RtmpiuIm9J/file10383cc1aeb.grd'
#> not recognized as a supported file format. (GDAL error 4)
#> Error: [rast] cannot open this file as a SpatRaster: C:/Users/Valentin/AppData/Local/Temp/RtmpiuIm9J/file10383cc1aeb.grd

# Using raster package to write:
r <- raster::raster(system.file("external/test.grd", package="raster"))
raster::writeRaster(raster::stack(lapply(1:100, function(i) r)),
                   filename = tfile, 
                   overwrite = TRUE) 

raster::stack(tfile)
#> class      : RasterStack 
#> dimensions : 115, 80, 9200, 100  (nrow, ncol, ncell, nlayers)
#> resolution : 40, 40  (x, y)
#> extent     : 178400, 181600, 329400, 334000  (xmin, xmax, ymin, ymax)
#> crs        : +proj=sterea +lat_0=52.1561605555556 +lon_0=5.38763888888889 +k=0.9999079 +x_0=155000 +y_0=463000 +datum=WGS84 +units=m +no_defs 
#> names      :   test.1,   test.2,   test.3,   test.4,   test.5,   test.6,   test.7,   test.8,   test.9,  test.10,  test.11,  test.12,  test.13,  test.14,  test.15, ... 
#> min values : 138.7071, 138.7071, 138.7071, 138.7071, 138.7071, 138.7071, 138.7071, 138.7071, 138.7071, 138.7071, 138.7071, 138.7071, 138.7071, 138.7071, 138.7071, ... 
#> max values : 1736.058, 1736.058, 1736.058, 1736.058, 1736.058, 1736.058, 1736.058, 1736.058, 1736.058, 1736.058, 1736.058, 1736.058, 1736.058, 1736.058, 1736.058, ...

# but reading still does not work
terra::rast(tfile)
#> Warning: `C:/Users/Valentin/AppData/Local/Temp/RtmpiuIm9J/file10383cc1aeb.grd'
#> not recognized as a supported file format. (GDAL error 4)
#> Error: [rast] cannot open this file as a SpatRaster: C:/Users/Valentin/AppData/Local/Temp/RtmpiuIm9J/file10383cc1aeb.grd

^{Created on 2023-04-26 with reprex v2.0.2}