1

I am working with climate data in R, I have 3 folders of raster maps: temperature (tmax,tmin) and precip (RF) - 9131 rasters in each folder, 27393 total.

I need to open each raster and extract data for every Spatial Point (284479 pts) and save the data in a new data.frame.

The goal is to have 284479 files, one for each Spatial Point with the data from each raster in columns (nrows = 9131, ncol = 3).

#Set spatial coordinates to create a Spatial object:
coordinates(GridPoint) = ~x + y
Coord <- coordinates(GridPoint)

#Create a loop to generate a input file for each station
for (i in 1:9131) {

#Open each raster per day of record
Tmax <- raster(Tmax_asc_files[i])
Tmin <- raster(Tmin_asc_files[i])
RF   <- raster(RF_asc_files[i])

#extract Values from each raster for that day
Tmax_E <- extract(Tmax, Coord)
tmax.df<- data.frame(Tmax_E)

Tmin_E <- extract(Tmin, Coord)
tmin.df<- data.frame(Tmin_E)

Prcp_E   <- extract(RF, Coord)
prcp.df  <- data.frame(Prcp_E)

for (j in 1:length(tmax.df))  {

Cell.matrix <- matrix(nrow = 9131,ncol = 4)
Cell.matrix <- as.data.frame(Cell.matrix)
colnames(Cell.matrix) <- c("yearday","tmax","tmin","prcp")

tmax.row <- tmax.df[j,]
tmin.row <- tmin.df[j,]
prcp.row <- prcp.df[j,]

#Add date to first column
Cell.matrix[1] <- Date
Cell.matrix[i,2:4] <- c(tmax.row, tmin.row, prcp.row)


filename <- paste("GridPt_", j, sep = "")
write.table(Cell.matrix, file = paste(path, filename, ".mtcin", sep = ""),
            row.names = F, quote = F) 
    }
  }
John Paul
  • 12,196
  • 6
  • 55
  • 75
Heidi Needham
  • 11
  • 2
  • Using extract, especially in combination with loops and large dataframes can be dangerous in terms of performance. 284479 points seems a lot. Are you extracting the enitre raster (every pixel)? – maRtin Dec 09 '16 at 01:00
  • Even though I am not extracting data at every pixel, the overlay grid with spatial points is very close to the entire raster. – Heidi Needham Dec 09 '16 at 01:15
  • If it would be the entire raster you could convert your rasters to matrices and use or use standard subsetting `raster[i]`or `matrix[i]` to access the values. this should be much faster than the `extract` function. – maRtin Dec 09 '16 at 01:18
  • I think you'be better off creating a RasterStack beforehand for each variable by `Tmax_stack <- stack(Tmax_asc_files, quick = TRUE)` instead than cycling on dates. That way, "extract" will at least result already in time series for each "position" you wish to extract, which should be faster. – lbusett Dec 10 '16 at 15:52

0 Answers0