当我尝试使用writeOGR加上循环来保存我的shapefile时,它除了给我一条错误消息外什么也没做:
writeOGR(plot.locationsSP_DROUGHT, dsn, layer1, driver = "ESRI Shapefile"( 中的错误:图层存在,使用新的图层名称
从本质上讲,我正在将每个对象转换为CSV文件,然后转换为形状文件,并希望同时保存CSV文件和形状文件。 这是我的代码片段:
for (m in 1:500){
#First I want to save my CSV files:
drought.slice <- rotate(drought.array[m,,])
drought.vec <- as.vector(drought.slice)
length(drought.vec)
drought.df01 <- data.frame(cbind(lonlat, drought.vec))
names(drought.df01) <- c("lon", "lat", paste(dname, as.character(m), sep = "_"))
head(na.omit(drought.df01))
csvfile<-paste0("cru_drought_",m,".csv")
#Next I want to create shapefiles from the CSV files:
plot.locations_DROUGHT <- read.csv(paste0("cru_drought_",m,".csv"), stringsAsFactors = FALSE)
plot.locationsSP_DROUGHT <- SpatialPointsDataFrame(plot.locations_DROUGHT[,1:2], plot.locations_DROUGHT)
proj4string(plot.locationsSP_DROUGHT) <- CRS("+init=epsg:4326")
dsn <- layer1 <- gsub(".csv","cru_drought_",m)
writeOGR(plot.locationsSP_DROUGHT, dsn, layer1, driver="ESRI Shapefile")
}
这是我正在使用的完整代码:
#Open and read the NCDF file, along with longitude and latitude
rm(list=ls())
library(lattice)
library(ncdf4)
library(chron)
library(rgdal)
library(sp)
library(raster)
library(RColorBrewer)
setwd('/Users/Neil/Dropbox/Drought Maps')
ncname <- "owda-orig"
ncfname <- paste(ncname,".nc",sep="")
dname <- "pdsi"
ncin <- nc_open(ncfname)
print(ncin)
lon <- ncvar_get(ncin, "lon")
nlon <- dim(lon)
head(lon)
lat <- ncvar_get(ncin, "lat", verbose = F)
nlat <- dim(lat)
head(lat)
print(c(nlon, nlat))
t <- ncvar_get(ncin, "time")
nt <- dim(t)
head(t)
drought.array <- ncvar_get(ncin, dname)
dlname <- ncatt_get(ncin, dname, "long_name")
dunits <- ncatt_get(ncin, dname, "units")
#fillvalue <- ncatt_get(ncin, dname, "_FillValue")
dim(drought.array)
creation_date <- ncatt_get(ncin, 0, "creation_date")
Description <- ncatt_get(ncin, 0, "Description")
nc_close(ncin)
rotate <- function(x) t(apply(x, c(1, 2), rev))
m <- 333
drought.slice <- rotate(drought.array[m,,])
image(lon, lat, drought.slice, col = brewer.pal(10, "BrBG"))
lonlat <- expand.grid(lon, lat)
drought.vec <- as.vector(drought.slice)
length(drought.vec)
drought.df01 <- data.frame(cbind(lonlat, drought.vec))
names(drought.df01) <- c("lon", "lat", paste(dname, as.character(m), sep = "_"))
head(na.omit(drought.df01))
for (m in 1:500){
drought.slice <- rotate(drought.array[m,,])
drought.vec <- as.vector(drought.slice)
length(drought.vec)
drought.df01 <- data.frame(cbind(lonlat, drought.vec))
names(drought.df01) <- c("lon", "lat", paste(dname, as.character(m), sep = "_"))
head(na.omit(drought.df01))
csvfile<-paste0("cru_drought_",m,".csv")
plot.locations_DROUGHT <- read.csv(paste0("cru_drought_",m,".csv"), stringsAsFactors = FALSE)
plot.locationsSP_DROUGHT <- SpatialPointsDataFrame(plot.locations_DROUGHT[,1:2], plot.locations_DROUGHT)
proj4string(plot.locationsSP_DROUGHT) <- CRS("+init=epsg:4326")
dsn <- layer1 <- gsub(".csv","cru_drought_",m)
writeOGR(plot.locationsSP_DROUGHT, dsn, layer1, driver="ESRI Shapefile")
}
帮助将不胜感激。 我可能正在做一些非常愚蠢的事情。
确定你的问题出在你制作形状文件名的行中。使用 gsub(".csv","cru_drought_",m)
会将字符串中的".csv"替换为"cru_drought_"m
,这只是您循环的整数。找不到该字符串,因此它只是将整数m
分配给dsn
,因此您最终会尝试编写名称为"1"、"2"等的 shapefile。我认为你只是把gsub
的论点弄乱了。
我注释掉了特定于您的文件的代码部分,并尝试将看起来像您所追求的文件名放在一起。
for (m in 1:5){
# drought.slice <- rotate(drought.array[m,,])
# drought.vec <- as.vector(drought.slice)
# length(drought.vec)
# drought.df01 <- data.frame(cbind(lonlat, drought.vec))
# names(drought.df01) <- c("lon", "lat", paste(dname, as.character(m), sep = "_"))
# head(na.omit(drought.df01))
csvfile<-paste0("cru_drought_",m,".csv")
# plot.locations_DROUGHT <- read.csv(paste0("cru_drought_",m,".csv"), stringsAsFactors = FALSE)
# plot.locationsSP_DROUGHT <- SpatialPointsDataFrame(plot.locations_DROUGHT[,1:2], plot.locations_DROUGHT)
# proj4string(plot.locationsSP_DROUGHT) <- CRS("+init=epsg:4326")
dsn <- layer1 <- gsub(".csv","_shape",csvfile)
# writeOGR(plot.locationsSP_DROUGHT, dsn, layer1, driver="ESRI Shapefile")
print(dsn)
}
#> [1] "cru_drought_1_shape"
#> [1] "cru_drought_2_shape"
#> [1] "cru_drought_3_shape"
#> [1] "cru_drought_4_shape"
#> [1] "cru_drought_5_shape"
创建于 2018-04-14 由 reprex 包 (v0.2.0(.
使用 raster::shapefile
writeOGR
更容易一些 shapefile(plot.locationsSP_DROUGHT, dsn)