以下栅格数据计算有错误,请修改library(terra) path1 <- dir('G:/ERA5hr/P/1/', pattern = '*.tif', full.names = TRUE) path2 <- dir('G:/ERA5hr/ZUST/1/', pattern = '*.tif', full.names = TRUE) path3 <- dir('G:/ERA5hr/ut_750/1/', pattern = '*.tif', full.names = TRUE) for (i in 1:length(path1)) { print(path1[i]) # 读取三个.tif文件 p <- stack(path1[i]) zust <- stack(path2[i]) ut <- stack(path3[i]) # 进行计算,使用if语句替代ifelse函数 if (zust > ut) { file_data <- (p/9.8) * zust^3 * (1 - ut/zust) * (1 + (ut/zust)^2) } else { file_data <- 0 } # 输出为新的.tif文件 output_path <- paste('G:/ERA5hr/Q/Q750/1/', substr(path1[i], 15, 28), '.tif', sep = '') writeRaster(file_data, output_path, overwrite = TRUE, format = "GTiff") }
时间: 2024-04-07 22:31:58 浏览: 63
在这段代码中,if语句不能直接用于两个raster对象之间的比较,因此需要使用相关函数来进行计算。以下是修改后的代码:
library(terra)
path1 <- dir('G:/ERA5hr/P/1/', pattern = '*.tif', full.names = TRUE)
path2 <- dir('G:/ERA5hr/ZUST/1/', pattern = '*.tif', full.names = TRUE)
path3 <- dir('G:/ERA5hr/ut_750/1/', pattern = '*.tif', full.names = TRUE)
for (i in 1:length(path1)) {
print(path1[i])
# 读取三个.tif文件
p <- rast(path1[i])
zust <- rast(path2[i])
ut <- rast(path3[i])
# 进行计算,使用ifelse函数替代if语句
file_data <- ifelse(zust > ut,
(p/9.8) * zust^3 * (1 - ut/zust) * (1 + (ut/zust)^2),
0)
# 输出为新的.tif文件
output_path <- paste('G:/ERA5hr/Q/Q750/1/', substr(path1[i], 15, 28), '.tif', sep = '')
writeRaster(file_data, output_path, overwrite = TRUE, format = "GTiff")
}
在修改后的代码中,使用了rast()函数来读取raster对象,使用了ifelse()函数来进行条件判断。
阅读全文