terra
terra copied to clipboard
rspatial
writeRaster: invalid GDAL options cause R to abort when using a large raster
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
──────────────────────────────────────────────────────────────────────
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?
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"
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
