批量提取逐日tif文件的像元值——基于arcgis and R

代码语言居然没有R 
## arcgis准备提取点的shp文件
# 研究区域的shp文件
# 使用arcgis栅格重采样工具确定DEM的分辨率
#  extract by mask 工具确定研究区域的栅格范围
#  raster to point 工具实现栅格转shp
#  栅格、shp统一的地理坐标系，所有的栅格（tif）影像要统一投影坐标系
#  在shp文件的属性表中，添加经度和纬度字段，获取每一个点的经纬度（单位m）。

## R 中实现批量提取数据
# 定义函数
extr_tif <- function(Tifpath,SHPpath,savepath){
     #  Tifpath tif文件所在的文件夹路径
     #  SHPpath  sho文件所在的文件夹路径
     #  savepath  数据保存到本地的路径
library(raster)
library(rgdal)
library(MODIS)
#定义数据框和列表
data1 <- list()
b <- data.frame()
# 变量的初始化
pathshp <- SHPpath
tifpath <- Tifpath
filepath <- savepath
setwd(tifpath)
# # 读取shp文件
ogrInfo(pathshp)
p <- readOGR(dsn = pathshp)# shp文件
# 读取tif文件
file <- list.files(tifpath,pattern=".tif$", full.names = FALSE)
# modis转tif的数据，所以要确定每一期影像的日期
#  modis的日期转化成年—月—日
a <- extractDate(file, asDate = TRUE)
date <- a$inputLayerDates 
a <- p@data
# 提取tif数据
for (j in 1:length(file)) {
    data1[[file[j]]] <- extract(raster(file[j]),p)
    b <- as.data.frame(data1[[j]])
    e <- data.frame(a$lat,a$lon,a$grid_code,b)
    #  保存数据
    filename <- paste(date[j],'.csv',sep = '')
    file1 <- paste(filepath,filename,sep = '/')
    write.csv(assign(paste('data_',date[j],sep='_'),e),file = file1,row.names = F)
   }
}
#  并行计算
library(parallel)
system.time({
    cl <- detectCores(logical = F)-1
    c2 <- makeCluster(c1)
    clusterEvalQ(c2, library(xts))
    clusterMap(c2, extr_tif, Tifpath <- '',# 添加路径，
               MoreArgs=list(SHPpath <- '',#  添加路径
                             savepath <- ''))#  添加路径
                             
    stopCluster(c2)
})
                                                                              2019-2-22