使用循环提取坐标，匹配它们并写入文件

Question

我正在尝试使用for 循环（或apply 函数作为替代）从data.frame 中提取坐标，在E-OBS gridded dataset 中搜索最近的点，提取温度-时间 x1-x2 的数据并将其写入另一个 excel 文件。

虽然代码可以提取单个数据点，但我似乎无法将此代码包含在循环中并将结果可能添加到输入坐标旁边。

library(sp)
library(raster)
library(ncdf4)

#Coordinates
    df
       ID    site                 E        N
1       1   site_place_date1  7.558758 47.81004
2       2   site_place_date2  7.582749 47.63411
3       3   site_place_date3  7.607968 48.01475
4       4   site_place_date4  7.644660 47.67139
       ...     ...   ...              ...`

手动设置目标点坐标：

lon <- 7.558758  # longitude of location                
lat <- 47.81004 # latitude  of location

#Mean daily temperature
    ncin <- nc_open("tg_0.25deg_reg_v17.0.nc")
      print(ncin)
      t <- ncvar_get(ncin,"time")
      tunits <- ncatt_get(ncin,"time","units")nt <- dim(t)
      nt
      obsoutput <- ncvar_get(ncin, 
                       start= c(which.min(abs(ncin$dim$longitude$vals -   lon)), # look for closest long
                                which.min(abs(ncin$dim$latitude$vals -  lat)),  # look for closest lat
                                1),
                       count=c(1,1,-1))
      DataMeanT <- data.frame(DateN= t, MeanDailyT = obsoutput)
      nc_close(ncin)
      head(DataMeanT)


#check if there are NAs =999
    summary(DataMeanT)

    Data = DataMeanT
    Data$Date = as.Date(Data$DateN,origin="20000-01-01")
    Data$Year = format(Data$Date,"%Y")
    Data$Month = format(Data$Date,"%m")
    head(Data)
    Data$YearMonth = format(Data$Date, format="%Y-%b")

    Data_annual = aggregate(("T_AnnualMean" = MeanDailyT) ~ Year,data = Data, FUN = mean,na.action = na.pass)
    names(Data_annual)[2] <- "AirT"
    head(Data_annual)

#Export table
    write.table(Data_annual, "Site_AirTemp.csv", row.names = FALSE, append = FALSE, col.names = TRUE, sep = ", ", quote = TRUE)

目的是将脚本作为 df 中所有坐标的循环的一部分运行，并将温度数据写入包含站点 ID 信息的新数据表，或者写入 df 的下一列。

Answer 1

只需将整个过程包装在一个定义的方法中，并使用 apply 函数来传递 lon/lat 坐标。一个很好的候选者是mapply 或其包装器Map 在df$E 和df$N 的两个向量之间逐元素迭代。此外，第三个参数 df$site 被传递给唯一 CSV 名称的方法，因为现在同一个文件将被覆盖。

下面的一些非赋值行如head 或summary 被删除，因为它们在方法中什么都不做。此外，上下文管理器within 和with 用于避免Data$ 的重复，以实现更简化的数据操作。 Map 调用写入文件并构建聚合数据帧列表以供以后使用。

功能

my_function <- function(lon, lat, site) {    
    # Mean daily temperature
    ncin <- nc_open("tg_0.25deg_reg_v17.0.nc")
      print(ncin)
      t <- ncvar_get(ncin,"time")
      tunits <- ncatt_get(ncin,"time","units")nt <- dim(t)

      # look for closest lon and lat
      obsoutput <- ncvar_get(ncin, 
                             start = c(which.min(abs(ncin$dim$longitude$vals - lon)),
                                      which.min(abs(ncin$dim$latitude$vals - lat)),
                                      1),
                             count = c(1,1,-1))

      DataMeanT <- data.frame(DateN = t, MeanDailyT = obsoutput)
    nc_close(ncin)    

    Data <- within(DataMeanT, {
               Date <- as.Date(DateN, origin="2000-01-01")
               Year <- format(Date,"%Y")
               Month <- format(Date,"%m")
               YearMonth <- format(Date, format="%Y-%b")
            })

    Data_annual <- with(Data, aggregate(list("AirT" = MeanDailyT), list(Year=Year),
                                        FUN = mean, na.action = na.pass))    
    # Export table
    write.table(Data_annual, paste0("Site_AirTemp_", site, "_.csv"), row.names=FALSE,
                append = FALSE, col.names = TRUE, sep = ", ", quote = TRUE)

    # SAVE AGGREGATED DATA FRAME
    return(Data_annual)
}

致电

# ITERATE THROUGH EACH LON/LAT PAIR ELEMENTWISE
df_list <- Map(my_function, df$E, df$N, df$site)

# df_list <- mapply(my_function, df$E, df$N, df$site, SIMPLIFY=FALSE)    # EQUIVALENT CALL

Answer 2

很难回答您的问题，因为它不可重现。但你可能可以这样做：

library(raster)
b <- brick("tg_0.25deg_reg_v17.0.nc")
e <- extract(b, df[, c('E', 'N')])