使用 R 中的 terra 包对不同范围的栅格进行马赛克时出错?

2024-01-02

我有一个文件夹,其中包含大约 191 个 GeoTIFF 文件(每个文件都是更大区域的不同 DEM(高程)图块)。我想将所有图块合并到一个光栅文件中。我正在使用terra包并成功加载每个栅格并将其从 2 米分辨率聚合到 30 米分辨率。然而,当运行mosaic函数将它们全部合并,我遇到错误(请参阅下面的错误消息)。我已经能够在仅三个图块的较小子集上运行马赛克功能,但是当我扩展到所有文件时,这就成为一个问题。

通过调用栅格摘要(见下文),聚合确实会稍微改变范围 - 这可能是问题所在吗?resample可能是一个选项,但每个单独的栅格都有不同的范围,我不确定如何实现此修复。

不确定示例数据集是否有帮助,因为我知道这些函数可以工作。我在高性能集群上运行此代码,因此运行小批量代码的效率不是很高。

library(terra)

files <- as.list(list.files("./DEM_tiles", full.names = TRUE))

raster.list <- lapply(files, rast)
 
for(i in 1:length(raster.list)){
 raster.list[[i]] <- aggregate(raster.list[[i]], fact = 15)
 }

raster.mosaic <- do.call(mosaic, raster.list)

> Error: [mosaic] internal error: extents do not match ()
> Execution halted

下面是两个图块的示例:

### Before Aggregation
[[1]]
class       : SpatRaster 
dimensions  : 25000, 25000, 1  (nrow, ncol, nlyr)
resolution  : 2, 2  (x, y)
extent      : -1800000, -1750000, -6e+05, -550000  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=stere +lat_0=90 +lat_ts=70 +lon_0=-45 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs 
source      : 35_23_1_1_2m_v3.0_reg_dem.tif 
name        : 35_23_1_1_2m_v3.0_reg_dem 

[[2]]
class       : SpatRaster 
dimensions  : 25000, 25000, 1  (nrow, ncol, nlyr)
resolution  : 2, 2  (x, y)
extent      : -1800000, -1750000, -550000, -5e+05  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=stere +lat_0=90 +lat_ts=70 +lon_0=-45 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs 
source      : 35_23_1_2_2m_v3.0_reg_dem.tif 
name        : 35_23_1_2_2m_v3.0_reg_dem 


### After Aggregation
[[1]]
class       : SpatRaster 
dimensions  : 1667, 1667, 1  (nrow, ncol, nlyr)
resolution  : 30, 30  (x, y)
extent      : -1800000, -1749990, -600010, -550000  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=stere +lat_0=90 +lat_ts=70 +lon_0=-45 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs 
source      : memory 
name        : 35_23_1_1_2m_v3.0_reg_dem 
min value   :                 -15.62178 
max value   :                  233.6489 

[[2]]
class       : SpatRaster 
dimensions  : 1667, 1667, 1  (nrow, ncol, nlyr)
resolution  : 30, 30  (x, y)
extent      : -1800000, -1749990, -550010, -5e+05  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=stere +lat_0=90 +lat_ts=70 +lon_0=-45 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs 
source      : memory 
name        : 35_23_1_2_2m_v3.0_reg_dem 
min value   :                 -15.27713 
max value   :                  243.0772 

该错误是因为光栅未对齐并且由于错误而导致。我现在得到

library(terra)
#terra version 1.2.1
crs <- "+proj=stere +lat_0=90 +lat_ts=70 +lon_0=-45 +x_0=0 +y_0=0 +datum=WGS84 +units=m"
r1 <- rast(nrow=1667, ncol=1667, ext=c(-1800000, -1749990, -600010, -550000), crs=crs)
r2 <- rast(nrow=1667, ncol=1667, ext=c(-1800000, -1749990, -550010, -5e+05), crs=crs)
values(r1) <- 1:ncell(r1)
values(r2) <- 1:ncell(r2)
m <- mosaic(r1, r2)
#Warning message:
#[mosaic] rasters did not align and were resampled

m
#class       : SpatRaster 
#dimensions  : 3334, 1667, 1  (nrow, ncol, nlyr)
#resolution  : 30, 30  (x, y)
#extent      : -1800000, -1749990, -600010, -499990  (xmin, xmax, ymin, ymax)
#coord. ref. : +proj=stere +lat_0=90 +lat_ts=70 +lon_0=-45 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs 
#source      : memory 
#name        :   lyr.1 
#min value   :       1 
#max value   : 2778889 

Also

@Elia 的建议是一个很好的解决方法:

r1 <- writeRaster(r1, "test1.tif", overwrite=TRUE)
r2 <- writeRaster(r2, "test2.tif", overwrite=TRUE)
v <- vrt(c("test1.tif", "test2.tif"), "test.vrt", overwrite=TRUE)

v
#class       : SpatRaster 
#dimensions  : 3334, 1667, 1  (nrow, ncol, nlyr)
#resolution  : 30, 30  (x, y)
#extent      : -1800000, -1749990, -600020, -5e+05  (xmin, xmax, ymin, ymax)
#coord. ref. : +proj=stere +lat_0=90 +lat_ts=70 +lon_0=-45 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs 
#source      : test.vrt 
#name        :    test 
#min value   :       1 
#max value   : 2778889 
本文内容由网友自发贡献,版权归原作者所有,本站不承担相应法律责任。如您发现有涉嫌抄袭侵权的内容,请联系:hwhale#tublm.com(使用前将#替换为@)

使用 R 中的 terra 包对不同范围的栅格进行马赛克时出错? 的相关文章

  • R 在列中按分隔符分割字符串

    我有一个包含几行的文件 例如 A B C awer ttp net Code 554 abcd ttp net Code 747 asdf ttp net Part 554 xyz ttp net Part 747 我想使用 R 仅拆分表的
  • 如何使用核心 R 操作/访问“dist”类实例的元素?

    R 中的基本 公共类称为 dist 并且是对称距离矩阵的相对有效的表示 不像一个 matrix 对象 但是 似乎不支持操纵 dist 使用索引对实例 操作员 例如 以下代码不返回任何内容 NULL 或出现错误 First create an
  • 在数据框中使用 Ifelse

    我正在使用的数据框是 gt df lt data frame Name c Joy Jane Jack Jad M1 c 10 40 55 90 gt df Name M1 1 Joy 10 2 Jane 40 3 Jack 55 4 Ja
  • 张量流:RStudio 中的 [NOT FOUND] 错误

    我尝试在中运行以下代码RStudio library tensorflow x data lt runif 100 min 0 max 1 y data lt x data 0 1 0 3 W lt tf Variable tf rando
  • 如何从 R 中的列表列表中提取元素?

    我有一堆列表 其中包含列表 广义线性模型输出 我想编写一个函数 该函数将从每个列表中提取多个元素 然后将结果组合到数据框中 我想提取modelset 1 likelihood modelset 1 fixef modelset 2 like
  • R:读取多个Excel文件,提取第一个工作表名称,并创建新列

    我有多个 Excel 文件 并且它们具有唯一的工作表名称 在我的情况下是文件创建日期 我批量阅读它们 需要将工作表名称分配给新列 id 中的每个文件 我知道如何制作数字 id 或 id 文件名 但找不到将工作表名称获取为 id 的方法 li
  • int NA 的内部表示[重复]

    这个问题在这里已经有答案了 这是关于 R 内部结构的问题 R 中如何表示整数 NA 值 与浮点不同 没有神奇的位序列来表示 NaN Create big array newer versions of R won t allocate me
  • 在 R 中创建多维 NetCDF

    我正在尝试使用 R 包创建多维 NetCDF 文件ncdf http cran r project org web packages ncdf index html 我正在对一组 1500 个点进行气候日常观测 每个点的观测数量约为 182
  • 尽管包在本地构建并通过了所有检查,但 CRAN 上的自动包提交错误

    我正在尝试向 CRAN 提交包 但它未通过一些自动检查 https win builder r project org incoming pretest influenceR 0 1 3 20230517 194638 Debian 00c
  • 将多个ggplot2图保存为列表中的R对象并在网格中重新显示

    我想在大型 for 循环期间将多个绘图 使用 ggplot2 保存到列表中 然后随后在网格中显示图像 使用 grid arrange 我已经尝试了两种解决方案 1 将其存储在列表中 如下所示 pltlist qplot lt qplot 然
  • 用 R 求解非平方线性系统

    如何用 R 求解非平方线性系统 A X B 系统无解或有无穷多个解的情况 例子 A matrix c 0 1 2 3 5 3 1 2 5 2 1 1 3 4 T B matrix c 17 28 11 3 1 T A 1 2 3 4 1 0
  • R 中的波形符(~) 运算符

    根据 R 文档 运算符在公式中用于分隔公式的右侧和左侧 右侧是自变量 左侧是因变量 我了解 lm 包中何时使用 然而以下是什么意思呢 x 1 右边是1 什么意思 可以是除 1 之外的任何其他数字吗 From lm 拟合线性模型时 y x 1
  • 使用 sf 与多多边形几何体进行分组(使用 R)

    我有一个放在一起的自定义形状文件 当我一次绘制所有内容时 效果很好 但我想按某些变量进行分组来绘制特定形状的区域 例如 county region sales washoe 1 5 carson city 1 10 clark 2 15 h
  • RQuantLib 的安装

    我尝试从 RStudio 安装 RQuantLib 但它给我带来了问题 我将 R 版本更新到 3 3 1 并尝试使用安装包的常用方法 install packages RQuantLib 按照作者网页上的推荐 http dirk eddel
  • GitHub Actions 中针对 Oauth 令牌的 Web 浏览器身份验证

    我正在尝试使用 GitHub Actions 通过 R 脚本使用 OAuth 2 0 令牌查询 GitHub API 当我运行代码时 它在我的本地计算机上运行良好 其中弹出一个浏览器窗口 指示 等待浏览器中的身份验证 我可以手动关闭该窗口
  • 从字符串列表中,识别哪些是人名,哪些不是

    我有一个如下所示的向量 想确定列表中的哪些元素是人名 哪些不是 我找到了 humaniformat 包 它可以格式化名称 但不幸的是它无法确定字符串是否实际上是名称 我还发现了一些用于实体提取的包 但它们似乎需要实际文本来进行词性标记 而不
  • 将函数参数传递给公式

    我试图理解为什么 foo function d y x fit with d lm y x foo myData Y X 不起作用 例如 myData data frame Y rnorm 50 X runif 50 对我来说似乎棘手的一点
  • R hdf5数据集写错了?

    当我执行以下命令时 我的 预测器 数据集已正确填充 library rhdf5 library forecast library sltl library tseries fid lt H5Fcreate output file TODO
  • R 闪亮动态输入

    我想构建一个 R 闪亮应用程序 它具有动态输入 要求用户输入数字 然后根据该输入生成另外 4 个输入字段 这就是我的想法 library shiny Define UI for random distribution application
  • 查找条目多于 R 中某个值的行总和

    我有以下矩阵 m structure 1 20 Dim 4 5 m 1 2 3 4 5 1 1 5 9 13 17 2 2 6 10 14 18 3 3 7 11 15 19 4 4 8 12 16 20 gt 我想找到每行中条目值大于 5

随机推荐