These data are not in a netcdf format but, rather a space delimited ASCII format. This data is a bit difficult to deal with because of the lack of headers and any type of unique station identifier.
Returns a vector of file names on disk and pulls associate years from file name.
f <- list.files(getwd())
y <- as.numeric(unlist(lapply(strsplit(f, "[.]"), function(x) {x[2]})))
Reads files, adds associated year and combines into single data.frame.
stations <- data.frame()
for(i in 1:length(f)) {
d <- read.table(f[1], sep = "" , header = FALSE, na.strings = "")
d[,"year"] <- y[i]
names(d) <- c("long","lat", format(ISOdate(2000, 1:12, 1), "%b"), "year")
stations <- rbind(stations, d)
}
stations$ID <- paste(stations[,1], stations[,2], sep="_")
str(stations)
Since there is no unique identifier for a given station, it is difficult to aggregate the duplicate coordinates to each given unique station location. I had to get creative and use a unique concatenated vector of the coordinate pairs. Once you have the data subset to the desired region, you can use merge() or which(x %in% y) to create a subset index, with the ID column, containing the unique coordinate pairs to relate back to the data. However, since the data is by year, in a wide format, you cannot store it as a spatial class object without duplicating the feature geometry. Please give some though on how you are going to analyze this data before just jumping in.
library(sp)
xy <- unique( paste(stations[,1], stations[,2], sep="_") )
stations.xy <- data.frame(long = as.numeric(unlist(lapply(strsplit(xy, "_"), function(x) {x[1]}))),
lat = as.numeric(unlist(lapply(strsplit(xy, "_"), function(x) {x[2]}))))
stations.xy$ID <- xy
coordinates(stations.xy) <- ~long+lat
str(stations.xy@data)
plot(stations.xy)
I am skipping any overlay/subset analysis because it is covered in your previous question. However, this can easily be accomplished using raster::intersect, spatialEco:point.in.poly, sp::over, rgeos::gIntersects, etc...
What happens if you run data.df <- sfSapply(list(veg.r), extract, y=data.sp)
? The question you link to is set up to extract data from a list of rasters, not just one.
Otherwise, how's performance if you rasterize your points layer and use it to extract values? Pseudocode:
presabs <- rasterize(points, vegraster, field = 'ID')
presabs[!is.na(values(presabs))] <- 0
xtrct <- presabs + vegraster
vals <- rasterToPoints(xtrct)
Then you could use this spatial join method to append the extracted data back to your original points. Suspect it might not be foolproof (multiple points in the same cell might cause issues) but maybe worth a shot.
Best Answer
getData()
only downloads tmean, tmin, tmax, prec, bio and alt from WorldClim v1.4. Check lines 252 to 325.You can download files with
utils::download.file()
function using WorlClim data URLs, but it will download a full file, not a specific tiles asgetData()
does.