R语言提取站点的nc文件时间序列数据

时间:2024-04-02 20:19:39
# 安装和加载所需的包 #install.packages("ncdf4") library(ncdf4) library(readxl) library(openxlsx) library(raster) library(geosphere) rm(list=ls()) setwd("E:/Program1/input/surfout-beijing-202007") # 读取包含经纬度的NetCDF文件 #nc_open 将数据读入我称为 nc_data 的数据结构中 nc_coords <- nc_open('../wrfout_d03_2020-08-02_00_00_00.nc') crs(nc_coords) names(nc_coords $var) # 读取经纬度变量lat"通常用来表示纬度(Latitude) lat <- ncvar_get(nc_coords, "XLAT") lon <- ncvar_get(nc_coords, "XLONG") dim(lat) # 关闭文件 nc_close(nc_coords) #读取站点信息 Stations=read_excel("../Beijing_staions.xlsx",sheet = 1, col_names = c("beijing","stations","city","lat","lon"), col_types = NULL, na = "", skip = 0) Stations <- as.data.frame(Stations) indexs=c("ALBEDO","T2","Q2","U10","V10","PSFC") # 提取WRDF对应站点模拟结果 for (station in 1:18){ city=Stations[station,2] station_lat=Stations[station,]$lat station_lon=Stations[station,]$lon # 计算站点与每个网格点的距离 #但是经纬度不等同于米的计算所以不可以 #计算栈距离最近的网格 distance=data.frame(matrix( nrow = 150, ncol = 150)) for(i in 1:150){ for(j in 1:150){ #distance[i,j] <- distGeo(c(lon[i,j],lat[i,j]), c(station_lon,station_lat)) distance[i,j] <- distm(c(lon[i,j],lat[i,j]), c(station_lon,station_lat)) } } #distance是米的距离 # 找到最小距离的索引按列平铺 #min_index <- which.min(distances) #返回最小值的行列号 #arr.ind = TRUE 则告诉函数返回的索引是数组形式的,而不是向量形式的 value_hl=which(distance==min(distance),arr.ind=T)#71 30 # 打开每日的nc文件,提取数据 station_city <- data.frame(matrix(nrow = 31 * 24, ncol = length(indexs) + 2)) colnames(station_city) <- c("day", "t","ALBEDO","T2","Q2","U10","V10","PSFC") #将nc读取合并 for (day in 1:31) { # 假设站点数据存储在另一个文件中 julyday <- sprintf("%02d", day) nc=paste("surfout_d03_2020-07-", julyday, "_00_00_00", sep = "") ncfile <- nc_open(nc) for(dex in 1:6){ index=indexs[dex] index_data=ncvar_get(ncfile,index) for(t in 1:24){ index_time=index_data[,,t] index_value=index_time[70,30] station_city[(day - 1) * 24 + t, "t"] <- t station_city[(day - 1) * 24 + t, "day"] <- day station_city[(day - 1) * 24 + t, index] <- index_value } } } nc_close(ncfile) # 将结果写入Excel文件 output_folder <- "../nc/" station_city_file <- paste(output_folder,"station", city, "_data.xlsx", sep = "") write.xlsx(station_city,station_city_file , rowNames= FALSE) }