R: цикл для сохранения нескольких шейп-файлов

Когда я пытаюсь использовать writeOGR плюс цикл для сохранения моих шейп-файлов, он ничего не делает, кроме как выдает сообщение об ошибке:

Ошибка в 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")
}

Помощь будет очень признательна. Я, наверное, делаю что-то очень глупое.


r shapefile rgdal ncdf4
person GIS_newb    schedule 14.04.2018    source источник

Ответы (2)


arrow_upward
0
arrow_downward

Я почти уверен, что ваша проблема в строке, в которой вы задаете имя шейп-файла. Использование gsub(".csv","cru_drought_",m) заменит ".csv" на "cru_daught_" в строке m, которая является просто целым числом, которое вы перебираете в цикле. Эта строка не найдена, поэтому она просто присваивает целое число m dsn, так что вы в конечном итоге пытаетесь записать шейп-файлы с именами типа «1», «2» и т. Д. Я думаю, что у вас просто несколько зашифрованных аргументов 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"

Создано 14 апреля 2018 г. пакетом REPEX (v0.2.0).

person camille    schedule 15.04.2018
comment
Спасибо большое, это работает! У меня еще один вопрос: похоже, что этот процесс создает отдельные папки для каждого шейп-файла. Есть ли способ просто сделать это и просто сохранить шейп-файлы в папке рабочего каталога? - person GIS_newb; 15.04.2018
comment
Аргумент dsn для writeOGR - это каталог, в котором вы сохраняете шейп-файл. Попробуйте изменить на что-то вроде writeOGR(shapefile, "folder_of_shapes", layer1, driver="ESRI Shapefile") - person camille; 15.04.2018

arrow_upward
0
arrow_downward

Немного проще использовать raster::shapefile вместо writeOGR

shapefile(plot.locationsSP_DROUGHT, dsn)
person Robert Hijmans    schedule 15.04.2018