为什么我可以绘制我的 stars_proxy 对象但不能写出来

Question

我是 Star 新手，所以希望这是一个简单的答案，只是我无法正确理解 Star 的工作流程。
R 版本：4.1.1
星星版本：0.5-5

library(stars)
library(starsdata) #install.packages("starsdata", repos = "http://gis-bigdata.uni-muenster.de", type = "source") 
#Create the rasters to read in as proxy
granule = system.file("sentinel/S2A_MSIL1C_20180220T105051_N0206_R051_T32ULE_20180221T134037.zip", package = "starsdata")
s2 = paste0("SENTINEL2_L1C:/vsizip/", granule, "/S2A_MSIL1C_20180220T105051_N0206_R051_T32ULE_20180221T134037.SAFE/MTD_MSIL1C.xml:10m:EPSG_32632")
r1<-read_stars(s2,,RasterIO=list(bands=1),proxy=T)
r2<-read_stars(s2,,RasterIO=list(bands=2),proxy=T)
r3<-read_stars(s2,,RasterIO=list(bands=3),proxy=T)
write_stars(r1,dsn="r1.tif")
write_stars(r2,dsn="r2.tif")
write_stars(r3,dsn="r3.tif")

然后我从我的环境中清除对象并重新启动 R 会话。

#I clear all the objects and restart my R session here.
library(stars)
foo<-read_stars(c("r1.tif","r2.tif","r3.tif"),proxy=T)
r1<- foo[1]*0
r1[foo[1] > 4000 & foo[2] < 3000] <- 1
r1[foo[1] > 4000 & foo[2] >= 3000 & foo[2] <= 8000] <- 2  
r1[foo[1] > 4000 & foo[2] > 8000 & foo[3] < 2000] <- 4
r1[foo[1] > 4000 & foo[2] > 8000 & foo[3] >= 2000] <- 2
# plot(r1) #this works just fine if you run it
#why doesn't the below work?
write_stars(r1,dsn="out.tif")

尝试写出文件会导致以下错误：

Error in st_as_stars.list(mapply(fun, x, i, value = value, SIMPLIFY = FALSE),  : 
  !is.null(dx) is not TRUE

如果我不是写出文件，而是绘制栅格，它就可以正常工作/按预期工作。 也许问题只是我不明白这个答案也适用于我： How to reassign cell/pixel values in R stars objects

Answer 1

首先感谢您为提供最小可重复示例所做的努力。不幸的是，您使用的图像很重......而且我的电脑很旧！ ;-) 所以我选择将您的示例与另一个图像（stars 库的测试图像）一起使用，这对于我的旧计算机来说更容易处理。

所以请在下面找到一个描述逐步方法的代表。

Reprex

• 第 1 步：从 stars 库的测试图像创建三个虚拟 stars proxy 对象

library(stars)

r <- system.file("tif/L7_ETMs.tif", package = "stars")

r1 <- read_stars(r, RasterIO = list(bands=1), proxy = TRUE)
r2 <- read_stars(r, RasterIO = list(bands=2), proxy = TRUE)
r3 <- read_stars(r, RasterIO = list(bands=3), proxy = TRUE)

• 第 2 步：将每个 stars proxy 对象作为 .tif 文件写入磁盘

write_stars(r1,dsn="r1.tif")
write_stars(r2,dsn="r2.tif")
write_stars(r3,dsn="r3.tif")

• 第 3 步：在 stars proxy 对象 foo 中合并 r1、r2 和 r3

foo <- read_stars(c("r1.tif","r2.tif","r3.tif"), proxy = TRUE)

foo <- merge(foo)

• 第 4 步：foo stars 代理对象的可视化

plot(foo)

如果要显示特定波段，请按以下步骤操作（此处显示波段 3）：

plot(foo[,,,3], main = st_dimensions(foo)["band"]$band$values[3])

• 第 5 步：运行您的代码块

r1 <- foo[,,,3]*0 #0 create a proxy with 0s that we will replace using rules below
r1[foo[,,,1] > 40 & foo[,,,2] < 30] <- 1
r1[foo[,,,1] > 40 & foo[,,,2] >= 30 & foo[,,,2] <= 70] <- 2  
r1[foo[,,,1] > 40 & foo[,,,2] > 70 & foo[,,,3] < 7] <- 4
r1[foo[,,,1] > 40 & foo[,,,2] > 70 & foo[,,,3] >= 7] <- 2

• 第 6 步：输出的可视化

plot(r1)

注意：我故意不在这里包含输出光栅，因为在您的代码块执行结束时，测试图像的所有像素的值都为 2。因此输出图像是单色光栅，没有任何兴趣[此结果与原始测试图像的像素值一致]。

• 第 7 步：保存输出

write_stars(r1, dsn = "out.tif")

• 第 8 步：检查文件是否已成功写入磁盘

file.exists("out.tif")
#> [1] TRUE

^{由reprex package (v2.0.1) 于 2021 年 12 月 10 日创建}