我尝试使用 getData(来自栅格包)读取 DEM 栅格，然后将 RasterLayer 转换为 SpatRaster(terra 包)。第一步成功了，但第二步失败了。

library(raster)
library(terra)
 
(alt <- getData('alt', country='DEU', mask=T))
#class      : RasterLayer 
#dimensions : 960, 1116, 1071360  (nrow, ncol, ncell)
#resolution : 0.008333333, 0.008333333  (x, y)
#extent     : 5.8, 15.1, 47.2, 55.2  (xmin, xmax, ymin, ymax)
#crs        : +proj=longlat +datum=WGS84 +no_defs 
#source     : D:/dummy/DEU_msk_alt.grd 
#names      : DEU_msk_alt 
#values     : -179, 2701  (min, max)

altT <- rast(alt)            
# rast is supposed to be able to read RasterLayer, but it fails. Why?
# Error : [rast] cannot read from D:/dummy/DEU_msk_alt.grd

一些提示？ :

rast("DEU_msk_alt.grd")
# Error : [rast] cannot read from D:/dummy/DEU_msk_alt.grd

rast("DEU_msk_alt.vrt")
#class       : SpatRaster 
#dimensions  : 960, 1116, 1  (nrow, ncol, nlyr)
#resolution  : 0.008333333, 0.008333333  (x, y)
#extent      : 5.8, 15.1, 47.2, 55.2  (xmin, xmax, ymin, ymax)
#coord. ref. : +proj=longlat +ellps=WGS84 +no_defs 
#data source : DEU_msk_alt.vrt 
#names       : DEU_msk_alt

看起来 rast 函数正在寻找 .vrt 文件，而 getData 将栅格关联到 grd 文件。无论如何，根据文档，当应用于 RasterLayer 时，rast 应该可以工作。

有什么想法吗？如何将这样的 RasterLayer 对象转换为 terra 对象？我想念什么？提前致谢，

JL

最佳答案

发生这种情况是因为 GDAL 驱动程序没有检测到正确的数据类型，因为 .grd 文件中的数据类型描述中有尾随空格:"INT2S " 而不是 "INT2S"

raster 包使用它自己的代码来读取这些文件。 terra 对 GDAL 的依赖性更强，并使用它来读取所有文件类型。由于这些文件相对较小，您也可以这样传输它们:

library(terra)
alt <- raster::getData('alt', country='DEU', mask=TRUE)
x <- rast(alt * 1)

或者你的解决方法

y <- rast( gsub("grd$", "vrt", filename(alt)) )

terra 的当前开发版本(版本 1.1-18；2021 年 4 月)现在可以读取这些文件，即使它仍然发出警告

#Unhandled datatype=INT2S  (GDAL error 1) 

https://stackoverflow.com/questions/65224689/