在R中循环从.nc文件提取指定经纬度温度时间序列求助

步骤1：准备观测点数据 #

先整理你的生物观测经纬度数据，示例如下（替换为你的实际数据）：

obs_points <- data.frame(
  site_id = c("site_01", "site_02"),  # 站点标识
  lon = c(116.4, 120.2),              # 站点经度
  lat = c(39.9, 30.3)                 # 站点纬度
)

步骤2：获取所有.nc文件并排序 #

先批量获取目标文件夹下的所有.nc文件，并从文件名解析日期、按时间排序：

folder <- "D:/TMP"
# 匹配符合命名规则的.nc文件
nc_files <- list.files(folder, pattern = "globaltemp_\\d{8}\\.nc$", full.names = TRUE, recursive = TRUE)

# 从文件名提取日期（格式：YYYYMMDD）
file_dates <- as.Date(sub("globaltemp_(\\d{8})\\.nc", "\\1", basename(nc_files)), format = "%Y%m%d")

# 按日期排序文件，保证时间序列顺序正确
file_order <- order(file_dates)
nc_files <- nc_files[file_order]
file_dates <- file_dates[file_order]

步骤3：编写单文件数据提取函数 #

封装单个.nc文件的处理逻辑，重点实现最近邻网格匹配和温度提取：

library(ncdf4)

extract_temp <- function(file_path, obs_df, lon_grid, lat_grid) {
  # 打开nc文件
  nc <- nc_open(file_path)
  
  # 读取温度数据并替换填充值为NA
  temp_array <- ncvar_get(nc, "temperature")
  fill_val <- ncatt_get(nc, "temperature", "_FillValue")$value
  temp_array[temp_array == fill_val] <- NA
  
  nc_close(nc)
  
  # 为每个观测点匹配最近的网格索引
  obs_df$lon_idx <- sapply(obs_df$lon, function(x) which.min(abs(lon_grid - x)))
  obs_df$lat_idx <- sapply(obs_df$lat, function(x) which.min(abs(lat_grid - x)))
  
  # 提取对应网格的温度值
  obs_df$temperature <- sapply(1:nrow(obs_df), function(i) {
    temp_array[obs_df$lon_idx[i], obs_df$lat_idx[i]]
  })
  
  # 添加日期列
  obs_df$date <- as.Date(sub("globaltemp_(\\d{8})\\.nc", "\\1", basename(file_path)), format = "%Y%m%d")
  
  # 返回关键列
  return(obs_df[, c("site_id", "date", "lon", "lat", "temperature")])
}

步骤4：批量处理所有文件 #

利用lapply批量处理文件，合并结果得到完整时间序列：

# 提前读取一次网格信息（所有文件网格一致时，避免重复读取）
sample_nc <- nc_open(nc_files[1])
lon_grid <- ncvar_get(sample_nc, "lon")
lat_grid <- ncvar_get(sample_nc, "lat")
nc_close(sample_nc)

# 批量处理所有文件
temp_results <- lapply(nc_files, extract_temp, obs_df = obs_points, lon_grid = lon_grid, lat_grid = lat_grid)

# 合并所有结果为一个数据框
temp_time_series <- do.call(rbind, temp_results)

# 按站点和日期排序
temp_time_series <- temp_time_series[order(temp_time_series$site_id, temp_time_series$date), ]

关于插值的说明 #

• 最近邻插值合理性：你推测的"匹配最近网格中心"完全可行。5km分辨率网格下，观测点与最近网格中心的最大距离约3.5km（网格对角线的一半），这个误差在气候数据的尺度下可以接受，且计算效率最高。
• 高精度插值可选：如果需要更精确的结果，可以使用双线性插值，示例代码如下（需先安装fields包）：

# library(fields)
# 将温度数组转换为矩阵格式
# temp_matrix <- matrix(temp_array, nrow = length(lon_grid), ncol = length(lat_grid))
# 双线性插值
# obs_df$temperature <- interp.surface(list(x = lon_grid, y = lat_grid, z = temp_matrix), obs_df[, c("lon", "lat")])

效率优化建议 #

• 并行处理：针对10年共3650个文件，可使用并行计算加速：

library(parallel)
num_cores <- detectCores() - 1  # 保留一个核心给系统
temp_results <- mclapply(nc_files, extract_temp, obs_df = obs_points, lon_grid = lon_grid, lat_grid = lat_grid, mc.cores = num_cores)

内容的提问来源于stack exchange，提问作者RGR_288