本文以TIF格式的气象数据为例,演示如何根据excel表格中的经纬度信息(坐标点),来提取TIF格式的气象数据集中的各气象因素。(前提:TIF格式数据中包含excel中的坐标点)。
library(terra)
library(stringr)
library(dplyr)
library(readxl)
# 设置工作目录
setwd('项目工作路径')
# 列出工作路径中所有TIFF文件
fls <- list.files(pattern = '.tif$', full.names = TRUE)
# 读取Excel文件中的经纬度数据
excel_file <- '项目工作路径/excel文件名.xlsx' # 替换为你的Excel文件路径
sheet_name <- 'Sheet1' # 替换为你的工作表名称
dt1 <- read_excel(excel_file, sheet = sheet_name)
# 使用“地址”列作为ID
#选择1
dt1$ID <- dt1$地址 #可能会有存在多对多的情况,即地址相同的情况,最好是用唯一ID号(可以提前在excel中生成)
#选择2
dt1$ID <- dt1$ID #建议用这个ID
# 创建矢量点
pits <- vect(dt1, geom=c("Longitude", "Latitude"), crs="+proj=longlat")
# 显示生成的点
plot(pits)
# 提取气象因子类型
factors <- unique(str_extract(basename(fls), "^[a-z]+"))
# 创建一个存放结果的表格
rst <- data.frame(ID = unique(dt1$ID))
# 创建进度条
total_iterations <- length(factors) * nrow(dt1)
progress_bar <- txtProgressBar(min = 0, max = total_iterations, style = 3)
iteration_count <- 0
# 循环处理每个气象因子类型
for (factor in factors) {
factor_files <- fls[grep(paste0("^", factor, "_"), basename(fls))]
r <- rast(factor_files)
plot(r, main = paste(factor, " data"))
# 确保矢量点和栅格数据具有相同的CRS
pits_proj <- project(pits, crs(r))
pit_bf <- buffer(pits_proj, 500) #设置缓冲区XXX米(500米)
# 初始化一个空的data.frame来存储提取的值
pit_value <- data.frame(ID = dt1$ID)
for (i in 1:nrow(dt1)) {
if (!is.na(dt1$Longitude[i]) && !is.na(dt1$Latitude[i])) {
# 提取缓冲区内的平均值
value <- extract(r, pit_bf[i], fun = 'mean', na.rm = TRUE)
pit_value[i, factor] <- value[1, 2] # 假设extract返回的data.frame有两列:ID和值
} else {
pit_value[i, factor] <- NA # 如果经纬度有缺失值,则填充NA
}
# 更新进度条
iteration_count <- iteration_count + 1
setTxtProgressBar(progress_bar, iteration_count)
}
rst <- left_join(rst, pit_value, by = 'ID')
}
# 关闭进度条
close(progress_bar)
# 保存结果
write.csv(rst, '保存的文件名.csv', row.names = FALSE)
数据以下图为例,获取以下11一个城市(有缺失值)的相应气象因素值。
最终结果大致形式如下:(结果仅供参考,这个几个湿度和温度的值是我随便编的,不是程序跑出来的)