3

我尝试使用 getData（来自 raster 包）读取 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

4

1 回答 1

2

发生这种情况是因为 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)
于 2020-12-10T03:00:33.263 回答